python中的高斯 - 勒让德算法

你忘记了括号4*t:

pi = (a+b)**2 / (4*t)

您可以使用decimal更高精度执行计算.

#!/usr/bin/env python
from __future__ import with_statement
import decimal

def pi_gauss_legendre():
    D = decimal.Decimal
    with decimal.localcontext() as ctx:
        ctx.prec += 2                
        a, b, t, p = 1, 1/D(2).sqrt(), 1/D(4), 1                
        pi = None
        while 1:
            an    = (a + b) / 2
            b     = (a * b).sqrt()
            t    -= p * (a - an) * (a - an)
            a, p  = an, 2*p
            piold = pi
            pi    = (a + b) * (a + b) / (4 * t)
            if pi == piold:  # equal within given precision
                break
    return +pi

decimal.getcontext().prec = 100
print pi_gauss_legendre()

输出:

3.141592653589793238462643383279502884197169399375105820974944592307816406286208\
    998628034825342117068