如何优化Python实现的Chudnovsky算法以更快计算π？

误差分析 #

你的代码实现逻辑是正确的，不存在原理性误差，具体验证如下：

• 算法匹配：代码完全符合Chudnovsky公式的定义，通过级数累加后计算π = C / sum（其中C = 426880*sqrt(10005)）的逻辑正确。
• 收敛次数合理：循环次数ceil(digits / 14.182)符合算法特性——每一项约贡献14.18位有效数字，确保级数足够收敛到目标精度。
• 精度缓冲充足：设置mp.dps = digits + 2提供了额外精度，抵消浮点运算的累积误差，保证最终结果的有效位数。
• 结果截取正确：通过str(...)[:digits + 1]截取字符串，正好保留整数位3加上digits位小数，输出格式符合预期。

只要输入的位数合理，代码输出的前digits位π值是准确的。

性能优化方向 #

要在1分钟内计算200万位π，需要针对性优化计算效率，以下是可行方案：

1. 替换高精度库：用gmpy2替代mpmath #

mpmath是纯Python实现的高精度库，运算速度有限；gmpy2基于C语言的GMP库开发，高精度运算效率远超mpmath，是提升速度的核心优化点。

示例优化代码：

import gmpy2
from math import ceil
from time import time

print("Chudnovsky algorithm (gmpy2 optimized)")

def chudnovsky(digits):
    # 设置足够的精度缓冲
    gmpy2.get_context().precision = digits + 10
    C = gmpy2.mpz(426880) * gmpy2.sqrt(gmpy2.mpz(10005))
    
    # 用整数类型存储递推变量，避免浮点开销
    L = gmpy2.mpz(13591409)
    M = gmpy2.mpz(1)
    X = gmpy2.mpz(1)
    K = gmpy2.mpz(-6)
    # 用分数类型累加，自动约分，减少误差和计算量
    sum_terms = gmpy2.mpq(0)
    
    iterations = ceil(digits / 14.181647462)  # 更精确的每迭代位数
    for q in range(iterations):
        sum_terms += gmpy2.mpq(M * L, X)
        L += 545140134
        X *= -262537412640768000
        K += 12
        # 整数递推M，避免浮点运算
        numerator = K**3 - 16 * K
        denominator = (q + 1)**3
        M = M * numerator // denominator
    
    pi = C / sum_terms
    # 转换为指定位数的字符串
    pi_str = gmpy2.to_str(pi, digits)
    return pi_str

# 用户输入逻辑不变
while True:
    try:
        user_input = int(input("\nHow many digits?: "))
        break
    except ValueError:
        print("Enter an integer")

start = time()
print()
print(chudnovsky(user_input))
print()
print(f"Elapsed time: [{time() - start}s]")

2. 优化递推逻辑，减少重复计算 #

• 避免在循环内重复计算(q+1)**3，可以用递推式推导：(q+1)**3 = q**3 + 3*q*q + 3*q + 1，进一步降低运算量。
• 预计算常数因子（如-262537412640768000）为整数，避免每次循环的类型转换开销。

3. 利用多核CPU并行计算 #

将级数拆分为多个子段，用multiprocessing模块分配到不同CPU核心并行计算，最后合并结果。例如将1000次迭代拆分为4个进程，每个进程计算250次，能显著缩短计算时间。

4. 调整精度设置 #

不要过度设置额外精度（如digits + 20），过多的精度会大幅增加计算时间。gmpy2中设置digits + 10已足够抵消误差。

内容的提问来源于stack exchange，提问作者PrancerFrisco