≡ Menu

The Gauss-Legendre Algorithm

The Gauss–Legendre algorithm is an algorithm to compute the digits of \(\pi\). It is notable for being rapidly convergent, with only 25 iterations producing 45 million correct digits of \(\pi\).

It repeatedly replaces two numbers (a e b, say) by their arithmetic and geometric means,

\(\frac{{a + b}}{2}{\text{ e }}\sqrt {ab}\)

in order to find their arithmetic-geometric mean.

The algorithm is quite simple, and is described here. Below you can find a Python code for implementing the algorithm, making use of the decimal module:

#!/usr/bin/env python
from __future__ import with_statement
import decimal

def pi_gauss_legendre():
    D = decimal.Decimal
    with decimal.localcontext() as ctx:
        ctx.prec += 2                
        a, b, t, p = 1, 1/D(2).sqrt(), 1/D(4), 1                
        pi = None
        while 1:
            an    = (a + b) / 2
            b     = (a * b).sqrt()
            t    -= p * (a - an) * (a - an)
            a, p  = an, 2*p
            piold = pi
            pi    = (a + b) * (a + b) / (4 * t)
            if pi == piold:  # equal within given precision
                break
    return +pi

decimal.getcontext().prec = 4500
print (pi_gauss_legendre())

This is a screen of output:

{ 0 comments… add one }

Rispondi

Next post:

Previous post: