# Chudnovsky algorithm

The Chudnovsky algorithm is a fast method for calculating the digits of π. It was published by the Chudnovsky brothers in 1989,[1] and was used in the world record calculations of 2.7 trillion digits of π in December 2009,[2] 5 trillion digits of π in August 2010,[3] 10 trillion digits of π in October 2011,[4][5] and 12.1 trillion digits in December 2013.[6]

The algorithm is based on the negated Heegner number d = −163, the j-function j(1+−163/2) = −6403203, and on the following rapidly convergent generalized hypergeometric series:[2]

${\displaystyle {\frac {1}{\pi }}=12\sum _{k=0}^{\infty }{\frac {(-1)^{k}(6k)!(545\,140\,134k+13\,591\,409)}{(3k)!(k!)^{3}\left(640\,320^{3}\right)^{k+{\frac {1}{2}}}}}.\!}$

For high performance iterative implementation, this can be simplified to -

${\displaystyle {\frac {\sqrt {640320^{3}}}{12\pi }}={\frac {426880{\sqrt {10005}}}{\pi }}=\sum _{k=0}^{\infty }{\frac {(6k)!(545\,140\,134k+13\,591\,409)}{(3k)!(k!)^{3}\left(-262\,537\,412\,640\,768\,000\right)^{k}}}.\!}$

There are 3 big integer terms (multinomial term Mk, linear term Lk and exponential term Xk) that make up the series and π is just the constant (C) divided by the sum of the series, as below.

${\displaystyle \pi =C/\sum _{k=0}^{\infty }{\frac {Mk*Lk}{Xk}}}$, where -
C = 426880√10005Mk+1 = Mk * (12k+2) * (12k+6) * (12k+10) / (k+1)^3 and M0 = 1Lk+1 = Lk + 545,140,134 and L0 = 13,591,409Xk+1 = Xk * -262,537,412,640,768,000 and X0 = 1Mk can be optimized further -Kk+1 = Kk + 12 and K0 = 6Mk+1 = Mk * (K^3 - (K<<4)) / (k+1)^3 and M0 = 1


Note that 545140134 = 163 × 3344418 and,

${\displaystyle e^{\pi {\sqrt {163}}}\approx 640\,320^{3}+743.999\,999\,999\,999\,25\dots }$ and
${\displaystyle 640\,320^{3}=262\,537\,412\,640\,768\,000}$

This identity is similar to some of Ramanujan's formulas involving π,[2] and is an example of a Ramanujan–Sato series.

## Example: Python Implementation

π can be computed to any precision using the above algorithm in any environment that supports Big Integer and configurable precision arithmetic (such as Python, Ruby, Perl, PHP, C#, Java, C++ and so many other languages, libraries, frameworks and runtimes). As an example, here is a Python implementation.

#!/usr/bin/env python -i

from decimal import getcontext, Decimal as Dec

def PI(maxK=70, prec=1008, disp=1007): # parameter defaults chosen to gain 1000+ digits within a few seconds
getcontext().prec = prec
k12_6,M,L,X,S = 6, 1, 13591409, 1, 13591409 # At k=0: 12k+6=6, (6k)!/(3k)!(k!)^3=1, 545140134k+13591409, (-640320^3)^k=1, 0+Dec(M*L)/X=L
for k in xrange(1,maxK+1):
M = (k12_6**3 - (k12_6<<4))*M / k**3; k12_6 += 12
L += 545140134
X *= -262537412640768000 # -640320**3;
S += Dec(M*L) / X
pi = 426880 * 10005**Dec(".5") / S
pi = Dec(("%s"%pi)[:disp]) # drop few digits of precision for accuracy
print "PI(maxK=%d iterations, getcontext().prec=%d, disp=%d digits) =\n%s" % (maxK, prec, disp, pi)
return pi

Pi = PI()
print "\nFor greater precision and more digits (takes a few extra seconds) - Try"
print "Pi = PI(317,4501,4500)"
print "Pi = PI(353,5022,5020)"


## References

1. ^ Chudnovsky, David V.; Chudnovsky, Gregory V. (1989), "The Computation of Classical Constants", Proceedings of the National Academy of Sciences of the United States of America, 86 (21): 8178–8182, doi:10.1073/pnas.86.21.8178, ISSN 0027-8424, JSTOR 34831, PMC , PMID 16594075
2. ^ a b c Baruah, Nayandeep Deka; Berndt, Bruce C.; Chan, Heng Huat (2009), "Ramanujan's series for 1/π: a survey", American Mathematical Monthly, 116 (7): 567–587, doi:10.4169/193009709X458555, JSTOR 40391165, MR 2549375
3. ^
4. ^ Yee, Alexander; Kondo, Shigeru (2011), 10 Trillion Digits of Pi: A Case Study of summing Hypergeometric Series to high precision on Multicore Systems, Technical Report, Computer Science Department, University of Illinois
5. ^ Aron, Jacob (March 14, 2012), "Constants clash on pi day", NewScientist
6. ^ Alexander J. Yee; Shigeru Kondo (28 December 2013). "12.1 Trillion Digits of Pi".