Given an array 'D' with n elements: d[0], d[1], ..., d[n-1], you can perform the following two steps on the array.  

Randomly choose two indexes (l, r) with l < r, swap (d[l], d[r])   
Randomly choose two indexes (l, r) with l < r, reverse (d[l...r]) (both inclusive)   

After you perform the first operation a times and the second operation b times, you randomly choose two indices l & r with l < r and calculate the S = sum(d[l...r]) (both inclusive).  

Now, you are to find the expected value of S.

Input Format 

The first line of the input contains 3 space separated integers - n, a and b. 

The next line contains n space separated integers which are the elements of the array d.   

n a b
d[0] d[1] ... d[n-1]

Output Format 

Print the expected value of S.  

E(S)

Constraints  

2 <= n <= 1000 

1 <= a <= 10^{9} 

1 <= b <= 10  

The answer will be considered correct, if the absolute or relative error doesn't exceed 10^{-4}.  

Sample Input #00:

3 1 1 
1 2 3

Sample Output #00:

4.666667

Explanation #00:

At step 1): 

You have three choices: 

1. swap(0, 1), 2 1 3 

2. swap(0, 2), 3 2 1 

3. swap(1, 2), 1 3 2  

At step 2): 

For every result you have three choices for reversing: 

1. [2 1 3] -> [1 2 3] [3 1 2] [2 3 1] 

2. [3 2 1] -> [2 3 1] [1 2 3] [3 1 2] 

3. [1 3 2] -> [3 1 2] [2 3 1] [1 2 3]  

So you have 9 possible arrays with each having a 1/9 probability.

For the last step: 

Each of the 9 arrays will have 3 possible sums with equal probability. For [1 2 3], you can get 1+2, 2+3 and 1+2+3. 

Since there will be 27 outcome for this input, one can calculate the expected value by finding sum of all 27 S and dividing it by 27.

Source: hackerrank | Difficulty: MEDIUM_HARD
Tags: Probability


Solutions:
def apply_Pa(x, a):
	n = len(x)
	lambda_1 = 1
	lambda_2 = ((n - 3) / (n - 1)) ** a
	diag_on = (lambda_1 * 1 + lambda_2 * (n - 1)) / n
	diag_off = (lambda_1 * 1 + lambda_2 * -1) / n
	x_sum = sum(x)
	x = [(diag_on - diag_off) * v + diag_off * x_sum for v in x]
	return x

def apply_Pb(x, b):
	n = len(x)
	mat = [[None] * n for _ in range(n)]
	for i in range(n):
		for j in range(i):
			mat[i][j] = min(n - i, j + 1)
		mat[i][i] = i * (i - 1) // 2 + (n - i - 1) * (n - i - 2) // 2 + min(i, n - i - 1)
		for j in range(i + 1, n):
			mat[i][j] = min(i + 1, n - j)
	res = x
	for _ in range(b):
		prev_res = res
		res = [0] * n
		for i in range(n):
			for j in range(n):
				res[i] += mat[i][j] * prev_res[j]
	norm = n * (n - 1) // 2
	norm = pow(norm, b)
	for i in range(n):
		res[i] = res[i] / norm
	return res

def expect_sum(x):

	def sum_k(l):
		return l * (l - 1) // 2
	n = len(x)
	number_of_sums = n * (n - 1) // 2
	e_sum = 0
	for i in range(len(x)):
		e_sum += x[i] * (number_of_sums - sum_k(i) - sum_k(n - i - 1))
	e_sum = e_sum / number_of_sums
	return e_sum
(n, a, b) = [int(x) for x in input().split()]
d = [int(x) for x in input().split()]
d1 = apply_Pa(d, a)
d2 = apply_Pb(d1, b)
d3 = expect_sum(d2)
print(d3)


def apply_Pa(x, a):
	n = len(x)
	lambda_1 = 1
	lambda_2 = ((n - 3) / (n - 1)) ** a
	diag_on = (lambda_1 * 1 + lambda_2 * (n - 1)) / n
	diag_off = (lambda_1 * 1 + lambda_2 * -1) / n
	x_sum = sum(x)
	x = [(diag_on - diag_off) * v + diag_off * x_sum for v in x]
	return x

def apply_Pb(x, b):
	n = len(x)
	mat = [[None] * n for _ in range(n)]
	for i in range(n):
		for j in range(i):
			mat[i][j] = min(n - i, j + 1)
		mat[i][i] = i * (i - 1) // 2 + (n - i - 1) * (n - i - 2) // 2 + min(i, n - i - 1)
		for j in range(i + 1, n):
			mat[i][j] = min(i + 1, n - j)
	res = x
	for _ in range(b):
		prev_res = res
		res = [0] * n
		for i in range(n):
			for j in range(n):
				res[i] += mat[i][j] * prev_res[j]
	norm = n * (n - 1) // 2
	norm = pow(norm, b)
	for i in range(n):
		res[i] = res[i] / norm
	return res

def expect_sum(x):

	def sum_k(l):
		return l * (l - 1) // 2
	n = len(x)
	number_of_sums = n * (n - 1) // 2
	e_sum = 0
	for i in range(len(x)):
		e_sum += x[i] * (number_of_sums - sum_k(i) - sum_k(n - i - 1))
	e_sum = e_sum / number_of_sums
	return e_sum
(n, a, b) = [int(x) for x in input().split()]
d = [int(x) for x in input().split()]
d1 = apply_Pa(d, a)
d2 = apply_Pb(d1, b)
d3 = expect_sum(d2)
print(d3)


import math

def center(x):
	mean = sum(x) / len(x)
	return (mean, [xi - mean for xi in x])

def smul(a, x):
	return [a * xi for xi in x]

def dot(x, y):
	return sum((xi * yi for (xi, yi) in zip(x, y)))

def mul(x, y):
	return [xi * yi for (xi, yi) in zip(x, y)]

def add(x, y):
	return [xi + yi for (xi, yi) in zip(x, y)]

def matmul(A, x):
	return [dot(ai, x) for ai in A]

def step1(n, a, d):
	(mean, d) = center(d)
	alpha = 2 / (n - 1)
	if alpha < 1:
		log1alpha = math.log1p(-alpha)
		beta = math.exp(a * log1alpha)
	else:
		beta = 0
	return (mean, smul(beta, d))

def step2(n, b):
	P = calc_P(n)
	nchoose2 = P[-1]
	(_, rho) = center([((i + 1) * (n - i) - 1) / nchoose2 for i in range(n)])
	return Pbx(P, b, rho)

def calc_P(n):
	A = [[min(min(i, j), n + 1 - max(i, j)) for j in range(1, n + 1)] for i in range(1, n + 1)]
	e = [1 for _ in range(n)]
	nchoose2 = n * (n - 1) // 2
	diag = [nchoose2 - di for di in matmul(A, e)]
	return (A, diag, nchoose2)

def Pbx(P, b, x):
	(A, diag, nchoose2) = P
	scale = 1 / nchoose2
	for _ in range(b):
		x = smul(scale, add(matmul(A, x), mul(diag, x)))
	return x

def get_ES(n, a, b, d):
	if n == 2:
		return sum(d)
	(mean, diff) = step1(n, a, d)
	rho = step2(n, b)
	return mean * ((n + 4) / 3) + dot(diff, rho)
(n, a, b) = map(int, input().split())
d = list(map(int, input().split()))
print(get_ES(n, a, b, d))


def apply_Pa(x, a):
	n = len(x)
	lambda_1 = 1
	lambda_2 = ((n - 3) / (n - 1)) ** a
	diag_on = (lambda_1 * 1 + lambda_2 * (n - 1)) / n
	diag_off = (lambda_1 * 1 + lambda_2 * -1) / n
	x_sum = sum(x)
	x = [(diag_on - diag_off) * v + diag_off * x_sum for v in x]
	return x

def apply_Pb(x, b):
	n = len(x)
	mat = [[None] * n for _ in range(n)]
	for i in range(n):
		for j in range(i):
			mat[i][j] = min(n - i, j + 1)
		mat[i][i] = i * (i - 1) // 2 + (n - i - 1) * (n - i - 2) // 2 + min(i, n - i - 1)
		for j in range(i + 1, n):
			mat[i][j] = min(i + 1, n - j)
	res = x
	for _ in range(b):
		prev_res = res
		res = [0] * n
		for i in range(n):
			for j in range(n):
				res[i] += mat[i][j] * prev_res[j]
	norm = n * (n - 1) // 2
	norm = pow(norm, b)
	for i in range(n):
		res[i] = res[i] / norm
	return res

def expect_sum(x):

	def sum_k(l):
		return l * (l - 1) // 2
	n = len(x)
	number_of_sums = n * (n - 1) // 2
	e_sum = 0
	for i in range(len(x)):
		e_sum += x[i] * (number_of_sums - sum_k(i) - sum_k(n - i - 1))
	e_sum = e_sum / number_of_sums
	return e_sum
(n, a, b) = [int(x) for x in input().split()]
d = [int(x) for x in input().split()]
d1 = apply_Pa(d, a)
d2 = apply_Pb(d1, b)
d3 = expect_sum(d2)
print(d3)