处理Python中大二项式的总和

作者: ℘ 宋唐 ℒº ѵ ℯ楚厦
发布时间: 2024-07-19 02:48:20 (17天前)
转自：

3 条回复

0#
回复此人
烏鴉喝酒 | 2019-08-31 10-32

<div class =“post-text”itemprop =“text”> 我自己发现了这个问题： 执行 <code> scipy.special.binom(99,50) </code> 给 <pre> <code> 5.044567227278209e+28 </code> </pre> 而在WolframAlpha上计算二项式（99,50）给出了 <pre> <code> 5.0445672272782096667406248628e+28 </code> </pre> 存在绝对差异，其数量级为10 ^ 12。 这就是为什么，对于高k值，python函数的结果肯定是不可靠的。所以我需要改变二项式的计算方式。 </DIV>

编辑
1#
回复此人
Hey ou | 2019-08-31 10-32

<div class =“post-text”itemprop =“text”> 我不明白你为什么要参与 <code> numpy </code> 功能在这里，以及为什么要转换为 <code> float </code> 对象。真的，对于这个公式，如果你的输入总是整数，那么简单地坚持下去 <code> int </code> 和 <code> fractions.Fraction </code> 你的答案永远都是 精确 。很容易实现自己的 <code> binom </code> 功能： <pre> <code> In [8]: def binom(n, k): ...: return ( ...: factorial(n) ...: // (factorial(k)*factorial(n-k)) ...: ) ...: </code> </pre> 注意，我使用整数除法： <code> // </code> 。最后，你的总和： <pre> <code> In [9]: from fractions import Fraction ...: def F(k, a, s): ...: result = Fraction(0, 1) ...: for i in range(1, k+1): ...: b = binom(k, i)*pow(-1, i+1) ...: x = Fraction(s, (a-1)*i - 1) ...: result += b*x ...: return result ...: </code> </pre> 结果如下： <pre> <code> In [10]: F(99, 3, 2) Out[10]: Fraction(47372953498165579239913571130715220654368322523193013011418, 1421930192463933435386372127473055337225260516259631545875) </code> </pre> 哪个基于wolfram-alpha似乎正确... 请注意，如果说， <code> alpha </code> 可以是非整数，你可以使用 <code> decimal.Decimal </code> 对于任意精度浮点运算： <pre> <code> In [17]: from decimal import Decimal ...: def F(k, a, s): ...: result = Decimal('0') ...: for i in range(1, k+1): ...: b = binom(k, i)*pow(-1, i+1) ...: x = Decimal(s) / Decimal((a-1)*i - 1) ...: result += b*x ...: return result ...: In [18]: F(99, 3, 2) Out[18]: Decimal('33.72169506311642881389682714') </code> </pre> 让我们提高精度： <pre> <code> In [20]: import decimal In [21]: decimal.getcontext().prec Out[21]: 28 In [22]: decimal.getcontext().prec = 100 In [23]: F(99, 3, 2) Out[23]: Decimal('33.31594880623309576443774363783112352607484321721989160481537847749994248174570647797323789728798446') </code> </pre> </DIV>

编辑

登录后才能参与评论