Gauss Legendre Algorithm in Python

Question

I need help computing Pi. I am trying to write a python program that will compute Pi numbers to X. I tried several from the python mailing list and it should slow down for my use. I read about the Gauss-Legendre Algorithm , and I tried porting it to Python without success.

I read from here and I would appreciate any data on where I am going wrong!

Output: 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((ya)) x = 2* x pi = (square((a+b)))/4*t print pi raw_input() 

Answer 1

• You forgot the brackets around 4*t :

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

• You can use decimal to do the calculation with more precision.

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

Output:

 3.141592653589793238462643383279502884197169399375105820974944592307816406286208\ 998628034825342117068 

Answer 2

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

it should be

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

Answer 3

• If you want to calculate PI up to 1000 digits, you need to use a data type that supports 1000 precision digits (e.g. mxNumber )
• You need to calculate a, b, t and x, while | ab | <10 ** - numbers, not iteration of numbers.
• Compute the square and pi as @JF suggests.