Quickly converging computation of π
The Gauss–Legendre algorithm is an algorithm to compute the digits of π . It is notable for being rapidly convergent, with only 25 iterations producing 45 million correct digits of π . However, it has some drawbacks (for example, it is computer memory -intensive) and therefore all record-breaking calculations for many years have used other methods, almost always the Chudnovsky algorithm . For details, see Chronology of computation of π .
The method is based on the individual work of Carl Friedrich Gauss (1777–1855) and Adrien-Marie Legendre (1752–1833) combined with modern algorithms for multiplication and square roots . It repeatedly replaces two numbers by their arithmetic and geometric mean , in order to approximate their arithmetic-geometric mean .
The version presented below is also known as the Gauss–Euler , Brent–Salamin (or Salamin–Brent ) algorithm ;[ 1] it was independently discovered in 1975 by Richard Brent and Eugene Salamin . It was used to compute the first 206,158,430,000 decimal digits of π on September 18 to 20, 1999, and the results were checked with Borwein's algorithm .
Initial value setting:
a
0
=
1
b
0
=
1
2
p
0
=
1
t
0
=
1
4
.
{\displaystyle a_{0}=1\qquad b_{0}={\frac {1}{\sqrt {2}}}\qquad p_{0}=1\qquad t_{0}={\frac {1}{4}}.}
Repeat the following instructions until the difference between
a
n
+
1
{\displaystyle a_{n+1}}
and
b
n
+
1
{\displaystyle b_{n+1}}
is within the desired accuracy:
a
n
+
1
=
a
n
+
b
n
2
,
b
n
+
1
=
a
n
b
n
,
p
n
+
1
=
2
p
n
,
t
n
+
1
=
t
n
−
p
n
(
a
n
+
1
−
a
n
)
2
.
{\displaystyle {\begin{aligned}a_{n+1}&={\frac {a_{n}+b_{n}}{2}},\\\\b_{n+1}&={\sqrt {a_{n}b_{n}}},\\\\p_{n+1}&=2p_{n},\\\\t_{n+1}&=t_{n}-p_{n}(a_{n+1}-a_{n})^{2}.\\\end{aligned}}}
π is then approximated as:
π
≈
(
a
n
+
1
+
b
n
+
1
)
2
4
t
n
+
1
.
{\displaystyle \pi \approx {\frac {(a_{n+1}+b_{n+1})^{2}}{4t_{n+1}}}.}
The first three iterations give (approximations given up to and including the first incorrect digit):
3.140
…
{\displaystyle 3.140\dots }
3.14159264
…
{\displaystyle 3.14159264\dots }
3.1415926535897932382
…
{\displaystyle 3.1415926535897932382\dots }
3.14159265358979323846264338327950288419711
…
{\displaystyle 3.14159265358979323846264338327950288419711\dots }
3.141592653589793238462643383279502884197169399375105820974944592307816406286208998625
…
{\displaystyle 3.141592653589793238462643383279502884197169399375105820974944592307816406286208998625\dots }
The algorithm has quadratic convergence , which essentially means that the number of correct digits doubles with each iteration of the algorithm.
Mathematical background [ edit ]
Limits of the arithmetic–geometric mean[ edit ]
The arithmetic–geometric mean of two numbers, a0 and b0 , is found by calculating the limit of the sequences
a
n
+
1
=
a
n
+
b
n
2
,
b
n
+
1
=
a
n
b
n
,
{\displaystyle {\begin{aligned}a_{n+1}&={\frac {a_{n}+b_{n}}{2}},\\[6pt]b_{n+1}&={\sqrt {a_{n}b_{n}}},\end{aligned}}}
which both converge to the same limit.
If
a
0
=
1
{\displaystyle a_{0}=1}
and
b
0
=
cos
φ
{\displaystyle b_{0}=\cos \varphi }
then the limit is
π
2
K
(
sin
φ
)
{\textstyle {\pi \over 2K(\sin \varphi )}}
where
K
(
k
)
{\displaystyle K(k)}
is the complete elliptic integral of the first kind
K
(
k
)
=
∫
0
π
/
2
d
θ
1
−
k
2
sin
2
θ
.
{\displaystyle K(k)=\int _{0}^{\pi /2}{\frac {d\theta }{\sqrt {1-k^{2}\sin ^{2}\theta }}}.}
If
c
0
=
sin
φ
{\displaystyle c_{0}=\sin \varphi }
,
c
i
+
1
=
a
i
−
a
i
+
1
{\displaystyle c_{i+1}=a_{i}-a_{i+1}}
, then
∑
i
=
0
∞
2
i
−
1
c
i
2
=
1
−
E
(
sin
φ
)
K
(
sin
φ
)
{\displaystyle \sum _{i=0}^{\infty }2^{i-1}c_{i}^{2}=1-{E(\sin \varphi ) \over K(\sin \varphi )}}
where
E
(
k
)
{\displaystyle E(k)}
is the complete elliptic integral of the second kind :
E
(
k
)
=
∫
0
π
/
2
1
−
k
2
sin
2
θ
d
θ
{\displaystyle E(k)=\int _{0}^{\pi /2}{\sqrt {1-k^{2}\sin ^{2}\theta }}\;d\theta }
Gauss knew of these two results.[ 2]
[ 3]
[ 4]
Legendre’s identity[ edit ]
Legendre proved the following identity:
K
(
cos
θ
)
E
(
sin
θ
)
+
K
(
sin
θ
)
E
(
cos
θ
)
−
K
(
cos
θ
)
K
(
sin
θ
)
=
π
2
,
{\displaystyle K(\cos \theta )E(\sin \theta )+K(\sin \theta )E(\cos \theta )-K(\cos \theta )K(\sin \theta )={\pi \over 2},}
for all
θ
{\displaystyle \theta }
.[ 2]
Elementary proof with integral calculus [ edit ]
The Gauss-Legendre algorithm can be proven to give results converging to
π
{\displaystyle \pi }
using only integral calculus. This is done here[ 5] and here.[ 6]
^ Brent, Richard , Old and New Algorithms for pi , Letters to the Editor, Notices of the AMS 60(1), p. 7
^ a b Brent, Richard (1975), Traub, J F (ed.), "Multiple-precision zero-finding methods and the complexity of elementary function evaluation" , Analytic Computational Complexity , New York: Academic Press, pp. 151– 176, archived from the original on 23 July 2008, retrieved 8 September 2007
^ Salamin, Eugene , Computation of pi , Charles Stark Draper Laboratory ISS memo 74–19, 30 January 1974, Cambridge, Massachusetts
^ Salamin, Eugene (1976), "Computation of pi Using Arithmetic–Geometric Mean", Mathematics of Computation , vol. 30, no. 135, pp. 565– 570, doi :10.2307/2005327 , ISSN 0025-5718 , JSTOR 2005327
^ Lord, Nick (1992), "Recent Calculations of π: The Gauss-Salamin Algorithm", The Mathematical Gazette , 76 (476): 231– 242, doi :10.2307/3619132 , JSTOR 3619132 , S2CID 125865215
^ Milla, Lorenz (2019), Easy Proof of Three Recursive π-Algorithms , arXiv :1907.04110