计算大二项系数的和

Question

我想计算以下这些和。但是我觉得二项式系数太大了，所以计算的时候出错了。

from __future__ import division  
import numpy as np 
from scipy.special import binom 
print [sum(binom(n,2*k)*np.sqrt(np.pi*k)**(-n/10) for k in xrange(1,int(n/2)+1)) for n in xrange(100000)]

有没有什么方法可以大致估算一下答案呢？

Answer 1

你遇到的问题是 scipy.special.binom 这个函数在计算结果时使用了近似值。你可以尝试用 scipy.misc.comb 这个函数来精确计算，记得加上可选参数 exact=True。

如果你设置 exact=False，它会用一个叫做伽马函数的东西来快速计算，但如果你想要精确的结果，可以把它设置为 exact=True。这样的话，它会返回一个 Python 的 long 类型的结果。

比如说：

In [1]: from scipy.misc import comb
In [2]: comb(1100, 600, exact=1)
Out[2]: 3460566959226705345639127495806232085745599377428662585566293887742644983083368677353972462238094509711079840182716572056521046152741092473183810039372681921994584724384022883591903620756613168264181145704714086085028150718406438428295606240034677372942820551517227766024953527980780035209056864110017856973033878393954438656320L

另外，你也可以尝试用其他的库，比如 gmpy，相关链接在 这里 和 这里。

Answer 2

嗯，二项式的计算很快就会达到上限：

from scipy.special import binom

binom(1019, 509)    # => 1.40313388415e+305
binom(1020, 510)    # => inf

你到底想要进行什么样的计算呢？

---

这里有一个重新整理过的版本，它稍微调整了一下数值；我们可以为每个n找到二项式系列的连续值，而不是每次都从头开始计算，并且我把平方根的幂运算合并成了一个操作。

from math import pi

for n in xrange(100000):
    total = 0.

    binom = 1
    binom_mul = n
    binom_div = 1

    power = -0.05 * n
    for k in xrange(1, n // 2 + 1):
        binom = binom * binom_mul / binom_div
        binom_mul -= 1
        binom_div += 1

        total += binom * (pi * k) ** power

    print(total)