1. History
The methods used for approximating $ \pi $ span millennia and severely vary in complexity; we will describe a few of the popular methods. The reference here is [1]. Around 200 BC, Archimedes approximated the circumference $ C $ and radius $ r $ of a circle by inscribing it in a polygon with $ n $ sides. It is not hard to see that as $ n $ approaches infinity, we get $ \pi $ using the circumference formula $ C = 2 \pi r $. Then during the invention of calculus in the 1600s, Newton and others used integrals and power series expansions to calculate $ \pi $. For instance, the identity
\[\arctan(x) = \sum_{n=0}^{\infty} (-1)^{n} \frac{x^{2n+1}}{2n+1}\]evaluated at $ x = 1 $ gives us
\[\frac{\pi}{4} = 1-\frac{1}{3}+\frac{1}{5}-\frac{1}{7}+\cdots.\]Furthermore, in the 1700s, Euler calculated the values of the Riemann zeta function. Famously, for $ x=2 $, we get that
\[\zeta(2) = \sum_{n=1}^{\infty} \frac{1}{n^2} = \frac{\pi^2}{6}.\]The Gregory series obtained from $ x=1 $ converges especially slowly, making it unsuitable for high-precision computation without further acceleration or a better arctangent identity.
In contrast, the Gauss-Legendre algorithm has quadratic convergence, so the number of correct digits roughly doubles at each step. For example, about 9 iterations suffice to calculate 512 decimal places of $ \pi $. As an aside, Ramanujan's equation for $ \pi $
\[\frac{1}{\pi} = \frac{2 \sqrt{2}}{9801} \sum_{n=0}^{\infty} \frac{(4n)! (1103+26390n)}{(n!)^4 396^{4n}}.\]needs around $ 65 $ terms to reach the same precision [3]. Ramanujan's work also inspired later rapidly convergent series used in high-precision computations of $ \pi $.
2. Elliptic integrals
Let $ (a \cos \theta, b \sin \theta) $ be an ellipse parameterized by $ \theta \in [0, 2 \pi] $. Then its arc-length is given by
\[\int_{0}^{2 \pi} \sqrt{a^2 \sin^2 \theta + b^2 \cos^2 \theta} d \theta. \tag{1}\]This integral is not easy to compute in itself, and it was generalized to the study of so-called elliptic integrals in the early 1700s. Now notice that in (1) all of the data concerning our specific ellipse's arc-length is contained within the interval $ 0 \leq \theta \leq \pi/2 $. Hence, it would make sense to reduce our study to complete elliptic integrals, meaning those with amplitude $ \pi/2 $. Furthermore, we will restrict ourselves to complete elliptic integrals of the first and second kind. For the proof that every elliptic integral is of the first, second, or third kind, see [2].
Definition. Let
\[\begin{aligned} F(k) & = \int_0^{\pi/2} \frac{1}{\sqrt{1 - k^2 \sin^2 \theta}} d \theta, \\ E(k) & = \int_0^{\pi/2} \sqrt{1 - k^2 \sin^2 \theta} d \theta, \\ F_S(a, b) & = \int_0^{\pi/2} \frac{1}{\sqrt{a^2 \cos^2 \theta + b^2 \sin^2 \theta}} d \theta, \\ E_S(a, b) & = \int_0^{\pi/2} \sqrt{a^2 \cos^2 \theta + b^2 \sin^2 \theta} d \theta. \end{aligned}\]Then $ F(k) $ and $ E(k) $ are called complete elliptic integrals of the first and second kinds, respectively, with symmetric forms $ F_S(a, b) $ and $ E_S(a, b) $. We refer to $ k $ as the modulus of our integral.
From here on we will assume that by elliptic integral we mean a complete elliptic integral of the first or second kind. We see that the ordinary and symmetric forms of elliptic integrals are related by the following proposition.
Proposition 2.1. Let $ a>0 $, $ 0<b\leq a $, and $ k^2 = 1 - b^2/a^2 $. Then the identities
\[\begin{aligned} F(k) & = a F_S(a, b) \\ E(k) & = \frac{1}{a} E_S(a, b). \end{aligned}\]are true.
Proof. Let us prove the second equation with the first following similarly. Substituting, we get
\[\begin{aligned} E(k) & = \int_0^{\pi/2} \sqrt{1 - k^2 \sin^2 \theta} d \theta \\ & = \int_0^{\pi/2} \sqrt{1 - (1-b^2/a^2) \sin^2 \theta} d \theta \\ & = \int_0^{\pi/2} \sqrt{\cos^2 \theta + b^2/a^2 \sin^2 \theta} d \theta \\ & = \frac{1}{a} \int_0^{\pi/2} \sqrt{a^2 \cos^2 \theta + b^2 \sin^2 \theta} d \theta \\ & = \frac{1}{a} E_S(a, b), \end{aligned}\]and hence our result. $ \blacksquare $
We also recall a result regarding the symmetry of these integrals originally discovered by Legendre.
Proposition 2.2. (Legendre's Identity). Suppose $ 0<k_1,k_2<1 $ and $ k_1^2 + k_2^2 = 1 $. Then
\[F(k_1) E(k_2) + F(k_2) E(k_1) - F(k_1) F(k_2) = \frac{\pi}{2}\]holds.
Proof. We leave the details of this proof to the reader. Taking the derivative with respect to $ k_1 $ shows that the L.H.S. is constant. To see that this value is $ \pi/2 $, we take the limit as $ k_1 $ goes to $ 0 $. $ \blacksquare $
3. Main results
We now prove the Gauss-Legendre algorithm. We will not discuss error analysis, which is done in [4].
Definition. Let $ a_0, b_0>0 $. Let $ a_{n+1} = (a_n + b_n)/2 $ and $ b_{n+1} = (a_n b_n)^{1/2} $ be the arithmetic and geometric means, respectively, of the $ n $th terms. After possibly interchanging $ a_0 $ and $ b_0 $, the sequence $ (a_n) $ is decreasing, $ (b_n) $ is increasing, and $ b_n\leq a_n $. Their limits agree because the recurrence gives $ 2\lim a_n=\lim a_n+\lim b_n $. We call their common limit
\[\lim_{n \to \infty} a_n = \lim_{n \to \infty} b_n = M(a_0, b_0)\]the arithmetic-geometric (AM-GM) mean of $ a_0 $ and $ b_0 $.
Theorem (Gauss-Legendre). Set $ a_0 = 1 $ and $ b_0 = \sqrt{2}/2 $. Then the series in the denominator below converges, and
\[\pi = \frac{2 M (1, \sqrt{2}/2)^2}{1 - \displaystyle{\sum_{n=0}^{\infty} 2^n (a_n^2 - b_n^2)}}\]More explicitly, if
\[\pi_N= \frac{2a_{N+1}^2} {1-\displaystyle\sum_{n=0}^{N}2^n(a_n^2-b_n^2)},\]then $\pi_N\to\pi$ quadratically. This is the computational version analyzed in [4].
We need two lemmas in order to prove our result.
Lemma 3.1. Let $ a_0 = a\geq b=b_0>0 $ as in the definition above. Set
\[S = a^2 - \sum_{n=0}^{\infty} 2^{n-1} \left( a_n^2 - b_n^2 \right).\]Then $ E_S(a, b) = S F_S(a, b) $.
Proof. This proof is taken from [5]. Consider the integral
\[L(a, b) = a^2 F_S(a, b) - E_S(a, b) \tag{2}\]Explicitly, this expands to
\[L(a, b) = (a^2 - b^2) \int_0^{\pi/2} \frac{\sin^2 \theta}{\sqrt{a^2 \cos^2 \theta + b^2 \sin^2 \theta}} d \theta.\]Substituting $ x^2 = a^2 \cos^2 \theta + b^2 \sin^2 \theta $,
\[L(a, b) = \int_b^a \sqrt{\frac{a^2-x^2}{x^2-b^2}} dx.\]The substitution used above also gives
\[F_S(a,b) =\int_b^a\frac{dx}{\sqrt{(a^2-x^2)(x^2-b^2)}}.\]Splitting this integral at $ \sqrt{ab} $ and applying $ y=(x+ab/x)/2 $ on both pieces gives Gauss's invariance
\[F_S(a,b)=F_S(a_1,b_1).\]Using the associated substitution $ y = (x+ab/x)/2 $ in $ L(a,b) $,
\[\begin{aligned} L(a, b) & = \frac{1}{2} \int_{b_1}^{a_1} \frac{(a^2 - b^2) + 4(a_1^2 - y^2)}{\sqrt{(a_1^2 - y^2)(y^2 - b_1^2)}} dy \\ & = \frac{1}{2} (a^2 - b^2) F_S(a_1, b_1) + 2 L(a_1, b_1) \\ & = \frac{1}{2} (a^2 - b^2) F_S(a, b) + 2 L(a_1, b_1), \end{aligned}\]and thus
\[\frac{L(a, b)}{F_S(a, b)} = \frac{1}{2}(a_0^2 - b_0^2) + 2 \frac{L(a_1, b_1)}{F_S(a_1, b_1)}.\]Repeatedly applying this identity gives, for every $ N $,
\[\frac{L(a,b)}{F_S(a,b)} =\sum_{n=0}^{N}2^{n-1}(a_n^2-b_n^2) +2^{N+1}\frac{L(a_{N+1},b_{N+1})}{F_S(a_{N+1},b_{N+1})}.\]The identity
\[a_{n+1}^2-b_{n+1}^2=\frac{(a_n-b_n)^2}{4}\]shows the quadratic convergence to zero. Since $ b_n $ remains bounded away from zero, the integral formula for $ L $ in (2) shows that the final remainder tends to zero. Therefore,
\[L(a, b) = \sum_{n=0}^{\infty} 2^{n-1} (a_n^2 - b_n^2) F_S(a, b),\]or equivalently
\[L(a, b) = - S F_S(a, b) + a^2 F_S(a, b).\]Combining this equation with (2) gives us the result. $ \blacksquare $
Lemma 3.2. For $ a,b>0 $, we have
\[F_S(a, b) = \frac{\pi}{2 M(a, b)}\]In particular,
\[F_S(1, \sqrt{2}/2) = \frac{\pi}{2 M(1, \sqrt{2}/2)}.\]Proof. The substitution used in the proof of lemma 3.1 gives
\[F_S(a_n,b_n)=F_S(a_{n+1},b_{n+1})\]at every step of the AM-GM iteration. Since $ a_n $ and $ b_n $ converge to $ M(a,b)>0 $, dominated convergence gives
\[\begin{aligned} F_S(a,b) &=\lim_{n\to\infty}F_S(a_n,b_n)\\ &=\int_0^{\pi/2}\frac{1}{M(a,b)}\,d\theta\\ &=\frac{\pi}{2M(a,b)}. \end{aligned}\]This proves the lemma. $ \blacksquare $
We finally have the tools we need to prove the Gauss-Legendre Algorithm does indeed converge to $ \pi $.
Proof (Gauss-Legendre). Set $ k = \sqrt{2}/2 $ to be our modulus. Then we notice $ 2k^2 = 1 $, and hence we can apply proposition 2.2 to get
\[2 F (k) E (k) - F(k)^2 = \frac{\pi}{2}.\]Now let us evaluate these integrals by first converting them into symmetric form then applying our lemmas. Since $ k^2 = 1-k^2 $, we can set $ a = 1 $ and $ b = k $ in proposition 2.1 to get
\[2 F_S(1, k) E_S(1, k) - F_S(1, k)^2 = \frac{\pi}{2}.\]Then applying lemma 3.1 gives us an equation only dependent on $ F_S(1, k) $
\[(2 S - 1) F_S(1, k)^2 = \frac{\pi}{2}.\]Finally, we apply lemma 3.2 to write our equation in terms of the AM-GM mean
\[\begin{aligned} \frac{\pi}{2} & = (2S-1) \left( \frac{\pi}{2 M(1, k)} \right)^2 \\ \pi & = \frac{2M(1, k)^2}{2S-1}. \end{aligned}\]Plugging in $ S $ gives us our result. $ \blacksquare $
4. References
- David H. Bailey, Simon M. Plouffe, Peter B. Borwein, and Jonathan M. Borwein, The quest for $ \pi $, The Mathematical Intelligencer 19 (1997), no. 1, 50–56.
- Gosta Mittag-Leffler, An introduction to the theory of elliptic functions, Annals of Mathematics 24 (1923), no. 4, 271–351.
- Srinivasa Ramanujan, Modular equations and approximations to $ \pi $, The Quarterly Journal of Pure and Applied Mathematics 45 (1914), 350–372.
- Eugene Salamin, Computation of $ \pi $ using arithmetic-geometric mean, Mathematics of Computation 30 (1976), no. 135, 565–570.
- Paramanand Singh, $ \pi $ and the AGM: Evaluating elliptic integrals, 2009. Link.