indc.py

from numpy import log10

def logfact(n):
    if n==0: 
        return 0
    s=0
    for i in xrange(1,n+1):
        s=s+log10(i)
    return s

def logbinom(p,q):
    return (logfact(p)-logfact(q)-logfact(p-q))

def fact(n):
    if n==0:
        return 1
    else:
        return n*fact(n-1)
        
def binom(p,q):
    return fact(p)/fact(q)/fact(p-q)

N=41
T=[binom(2*N,i) for i in xrange(2*N+1)]
s=[T[0]]
for i,x in enumerate(T[1:]):
    s.append(s[i]+x)

for x in reversed(s[:-1]):
    print round(log10(float(x))-N*log10(4),3),