Algoritmo de Gauss-Legendre en python

Question

Necesito ayuda para calcular Pi. Estoy tratando de escribir un programa de Python que calculará Pi a X dígitos. He intentado varios de la lista de correo de python, y es lento para mi uso. He leído sobre Gauss-Legendre Algorithm, y he tratado de portarlo a Python sin éxito.

Estoy leyendo desde Here, y agradecería cualquier opinión sobre dónde me estoy equivocando!

Genera: 0,163991276262

from __future__ import division 
import math 
def square(x):return x*x 
a = 1 
b = 1/math.sqrt(2) 
t = 1/4 
x = 1 
for i in range(1000): 
    y = a 
    a = (a+b)/2 
    b = math.sqrt(b*y) 
    t = t - x * square((y-a)) 
    x = 2* x 

pi = (square((a+b)))/4*t 
print pi 
raw_input() 

Answer 1

1. se le olvidó paréntesis alrededor 4*t:

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

2. Puede utilizar decimal para realizar el cálculo con mayor precisión.

   #!/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() 
   

Salida:

3.141592653589793238462643383279502884197169399375105820974944592307816406286208\ 
    998628034825342117068 

Answer 2

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

debería ser

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

Answer 3

1. Si desea calcular Pi a 1000 dígitos es necesario utilizar un tipo de datos que se admite 1000 dígitos de precisión (por ejemplo, mxNumber)
2. Necesita calcular a, b, t y x hasta | a-b | < 10 ** - dígitos, no iterar dígitos veces.
3. Calcula cuadrado y pi como @JF. sugiere.