The constant p ¼

4  Pi and the Arithmetic-Geometric Mean

4.1  The AGM sequence

Studied by Lagrange, Gauss and Legendre the Arithmetic-Geometric Mean (or AGM) of two positive numbers a and b is defined to be the common limit of the two sequences { a_n,b_n} _{n ³ 0} determined by:

We assume a ³ b. From the immediate relation

we have a_n ³ b_n and then a_n ³ (a_n+b_n)/2 = a_n+1 and b_n+1 = Ö{a_nb_n} ³ b_n for every n. Therefore, by induction:

From the first two relations we observe that a_n and b_n converge monotonically to a common limit determined by the initial conditions (a,b). We denote the common limit by:

More, from the third relation, the convergence is of order 2 (or quadratic):

The Arithmetic-Geometric Mean has obvious properties:

Example 1 M(Ö2,1) (this value is related to the lemniscate's constant and the following numerical computations were made by Gauss in 1799)

The rate of convergence is impressive, | a₆-b₆| < 10^-86!!

For more details on the AGM, and for it's generalization on complex numbers, read the very nice article from Cox ([3]).

4.2  The AGM and the complete elliptic integral of the first kind

The main result is the following theorem, where I(a,b) denotes the first kind complete elliptic integral:

Theorem 1 We suppose that a ³ b are positive numbers then:

Proof : One proof may be given by two substitutions. The first one, u = btanq is classic and transforms the elliptic integral into

The second substitution v = (u-ab/u)/2 gives (with some care at u = 0):

and the last integral is trivial so:

This theorem is most important and allows to compute the first kind complete elliptic integral with the Arithmetic-Geometric Mean.

4.3  The Legendre relation

The second kind complete elliptic integral is defined as the length of a quarter of an ellipse with major axis a and minor axis b:

If (b,b¢) are linked by the relation b²+b¢² = 1, we have the famous Legendre's relation:

Theorem 2

Rather technical, it can be proved, like Legendre did, by derivation of the function L(b) with respect to the variable b, the result is the null function. The constant value of this function is given with b = 0 and b^¢ = 1.

Theorem 3 Like for the function I(a,b) we may prove :

Using the AGM iteration's notation we may write the previous relation as:

with some basic manipulations we find:

4.4  Expression for p

We apply the Legendre's relation with b = b^¢ = 1/Ö2 so:

one can now eliminate J(1,1/Ö2) in the Legendre's relation:

and from I(a,b) = p/(2M(a,b)) we obtain:

we have just proved the very important result:

Theorem 4 With a₀ = 1, b₀ = 1/Ö2 and { a_k,b_k} _{k ³ 0} defined by the arithmetic-geometric mean:

This formula for p was discovered in 1976 by Eugene Salamin and also at the same time by Richard P. Brent ([1], [2]). The convergence is of order 2 (the number of decimal places is doubled at each iteration) and if one is able to compute M(1,1/Ö2) with efficient algorithm, this leads to a very nice way to compute many digits of p.

4.5  Description of an algorithm using Brent-Salamin formula

With the previous notation the algorithm starts with the initial values:

with (a,b,d,s,t) being multi-precision numbers and k an integer. Then the following calculations are made until the required precision is reached (a-b < e):

The value for p is then:

In this algorithm, multi-precision squaring, square roots and one division are required. During the last years, many record of computation were achieved with such formulae. Jointly other formulae are used to check the result (for instance the Borwein quartically algorithm is used to compare the result with Brent-Salamin algorithm).

4.6  The Borwein quartic iteration

In [4], the following algorithm is given

starting with

then a_k converge quartically to 1/p. This algorithm provides a very fast and efficient way to compute p and is often used with the Brent-Salamin algorithm to establish new record of computation. Of course it's necessary to compute square roots with high precision and this can be done efficiently thanks to the Fast Fourier Transform.

Here are the first iterates :

and 1/a₃ is about 170 digits, 1/a₄ about 694 digits, ... The number of digits increases by a factor of four at each iteration (see [5]).

In the same book we also find quadratic, cubic, quintic, septic, ... algorithms.

4.7  Some records of computation

1988 Kanada gave 201,326,000 decimal places with the Brent-Salamin formula. Only twenty eight iterations of the Arithmetic-Geometric Mean were needed to achieve such an accuracy. The Borwein iteration was used for the verification.

1996 Kanada gave 17,100,000,000 digits using the same algorithm. The computation was performed on a Hitachi parallel computer.

1999 Kanada gave 68 billions of digits with the same algorithm, then 206 billions digits (Sept. 1999).

References

[1]

R.P. Brent, Fast multiple-Precision evaluation of elementary functions, J. Assoc. Comput. Mach., (1976), vol. 23, p. 242-251

[2]

E. Salamin, Computation of p Using Arithmetic-Geometric Mean, Mathematics of Computation, (1976), vol. 30, p. 565-570

[3]

A. Cox, The Arithmetic-Geometric Mean of Gauss, L'enseignement Mathématique, (1984), vol. 30, p. 275-330

[4]

J.M. Borwein and P.B. Borwein, Pi and the AGM - A study in Analytic Number Theory and Computational Complexity, A Wiley-Interscience Publication, New York, (1987)

[5]

J.M. Borwein and P.B. Borwein, Ramanujan and Pi, Scientific American, (1988), p. 112-117