Python中的Gauss-Legendre算法

1. 您忘记了以下括号

   4*t

   ：

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

2. 您可以

   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