3230.sagews

from operator import mul
def W(n):
    return reduce(mul, range(1,2*n,2))
def L(n):
    if (n==1):
        return 1
    else:
        return reduce(mul, range(2,2*n,2))

def bp(n, k, v):
  result = 0
  for p in Partitions(n, length=k).list():
    factorialProduct = 1
    powerFactorialProduct = 1
    for index in set(p):
      count = list(p).count(index)
      factorialProduct = factorialProduct * factorial(count)
      powerFactorialProduct = powerFactorialProduct * factorial(index)^count
    coefficient = factorial(n) / (factorialProduct * powerFactorialProduct)
    result = result + coefficient * reduce(mul, [v[i-1] for i in p])
  return result

def shoelaceExpectationNumerator(n):
    # build vector of L1...Ln
    v = []
    for i in range(1,n+1):
        v.append(L(i))
    sum = 0
    for i in range(1,n+1):
        sum = sum + i*bp(n,i,v)
    return sum

def shoelaceExpectation(n):
    return float(shoelaceExpectationNumerator(n))/W(n)

print shoelaceExpectation(2)
print shoelaceExpectation(3)
print shoelaceExpectation(4)
print shoelaceExpectation(30)

1.33333333333 1.53333333333 1.67619047619 2.68237684749