Calculating Pi in Python

Question

I have implemented a calculation of Pi in Python, just for the lulz. But now I wonder if I could improve this more. Could you help me out? I do not program a lot in Python.

from time import time
from decimal import *
def faculty(m):
 k=m
 n=1
 while(n<m):
 k=k*n
 n+=1
 if(k>0):
 return k
 else:
 return 1
def ramanujan(b, n):
 return b+Decimal(faculty(4*n)*(1103+26390*n))/Decimal((faculty(n)**4)*396**(4*n))
def chudnovsky(b, n):
 return b+(((-1)**n)*(faculty(6*n))*(13591409+(545140134*n)))/(faculty(Decimal(3)*n)*(faculty(n)**3)*640320**(3*n+(Decimal(3)/Decimal(2))))
def calc(x, a, func):
 b=Decimal(0)
 n=Decimal(0)
 while(n<x):
 b = func(b, n)
 n=n+1
 return ((a*b)**(-1))
def calcrama():
 return calc(20, Decimal((2*Decimal(2).sqrt())/9801), ramanujan)
def calcchud():
 return calc(20, 12, chudnovsky)
def save(name, func):
 fobj = open(name, "w")
 t = time()
 pi = func()
 t = time() - t
 fobj.write("Time: "+str(t)+"\nPi: "+str(pi))
 fobj.close()
getcontext().prec = 1000
save("rama.txt", calcrama)
save("chud.txt", calcchud)

Answer 1

I find that your program is obfuscated and therefore incomprehensible to anyone except for perhaps an expert in those calculations methods. You should pick more informative variable names than a, b, k, m, n, and x.

faculty(m) should be renamed factorial(m). Ideally, you wouldn't need it — see below.

If you could verify that the denominator is always a divisor of the numerator, you should use // integer division.

Use more whitespace for readability. Here, I would break with Python guidelines and be more generous with whitespace than usual.

Verbose comments, especially docstrings, would also be appreciated for complicated mathematics:

def chudnovsky(b, n):
 """
 Based on the identity
 n
 1 inf (-1) (6n)! (13591409 + 545140134n)
 ----- = SUM -----------------------------------
 12 PI n=0 3 3n + 1.5
 (3n)! (n!) 640320
 Returns the nth term of the sum.
 """
 return b + (
 ((-1)**n) * (factorial(6 * n)) * 
 (13591409 + (545140134 * n))
 ) // (
 factorial(Decimal(3) * n) *
 (factorial(n)**3) *
 640320**(3 * n + (Decimal(3) / Decimal(2)))
 )

Many parts of that expression would be more efficiently computed based on their counterparts in the previous term, instead of starting from scratch. You might want to build an generator that caches those intermediate results.

Furthermore, you have very large numbers in both your numerators and denominators, especially factorials. You can help keep those small by cancelling out the factorials in the numerator and denominator as soon as possible.

def ChudnovskyTerms():
 """
 Based on the identity
 n
 1 inf 12 (-1) (6n)! (13591409 + 545140134n)
 ---- = SUM --------------------------------------
 PI n=0 3 3n + 1.5
 (3n)! (n!) 640320
 yields successive terms of the sum.
 """
 term = 0
 sign = 1
 factorial_6_3 = 1 # = factorial(6 * n) / factorial(3 * n)
 numerator_sum = 13591409
 n_factorial_cubed = 1
 denominator = 640320**1.5
 while True:
 yield 12 * sign * factorial_6_3 * ...
 term += 1
 sign *= -1
 for i in range(6):
 factorial_6_3 *= (6 * term - i)
 for i in range(3):
 factorial_6_3 //= (3 * term - i)
 numerator_sum += 545140134
 # etc.

• NaCl: So wrapping it in a class increases the speed significantly?
• 200_success: Keeping intermediate results rather than discarding them increases the speed. Packaging the code in a generator gives it a reasonable interface, and the overhead wouldn't be worse than what you had before.
• Gareth Rees: self.term looks wrong to me: what's self?
• 200_success: @GarethRees Thanks for spotting that. It was left over from an old version before I rewrote it as a generator.