PYTHON 尝试用 Chudnovsky 算法准确得到圆周率的 100 位

Question

我正在尝试使用 Python 和此维基百科页面自学 Chudnovsky 算法：

https://en.wikipedia.org/wiki/Chudnovsky_algorithm

在 wiki 上，我专注于“高性能迭代实现，[that] 可以简化为”：

我尝试编写最右边使用 Sigma 符号的方程式。我熟悉 Python 但数学不是很好。我给自己定下的目标是看能不能准确打印出至少100位圆周率。

公式中有 5 组括号，因此我尝试对 5 个不同的组件分别进行编码。我还编写了一个执行阶乘的函数，因为阶乘用于 5 个 components/parentheses.

中的 3 个

这是我的 23 行工作代码，有人可以帮我理解为什么它不能准确地达到 100 位吗？它精确到第 28 位：3.1415926535897932384626433832。然后对于第 29 个数字，它说 8 但它应该是 7...

import math
from decimal import *

def factorial(n):
    if n == 0:
        return 1
    memory = n
    
    for i in range(1, n):
        memory *= i
    
    return memory

iterations = 500
_sum = 0

#here's the Sigma part
for q in range(0, iterations):
    a = factorial(6*q)
    b = (545140134*q) + 13591409
       
    c = factorial(3*q)
    d = (factorial(q))**3
    e = (-262537412640768000)**q
    
    numerator = (a*b)
    denominator = (c*d*e)

    _sum += Decimal(numerator / denominator)

#ensures that you get 100 digits for pi
getcontext().prec = 100

sq = Decimal(10005).sqrt()
overPI = Decimal(426880 * sq)

pi = (overPI) * (Decimal(1 / _sum))
print("Pi is", pi)

感谢您提供的任何帮助。

Answer 1

在文件的开始设置小数精度，以便应用于所有后续操作：
```
from decimal import *
getcontext().prec = 100 
```

在除法或乘法之前将一个或两个操作数转换为 Decimal：

# bad, float to decimal -> loss of precision
_sum += Decimal(numerator / denominator)

# better, precision preserved
_sum += Decimal(numerator) / Decimal(denominator)

结果 - 精确到 98 d.p。（100 s.f。减去舍入误差）：

#     3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679
Pi is 3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117069

PYTHON 尝试用 Chudnovsky 算法准确得到圆周率的 100 位

Trying to accurately get 100 digits of pi with Chudnovsky's algorithm in PYTHON

python

python-3.x

math

symbolic-math

pi