Computing π

From Wikipedia, the free encyclopedia

(Redirected from Software for calculating π)

The following contains information regarding the computation of pi.

[edit] Standard methods

[edit] Summing a circle's area

Pi can be obtained from a circle if its radius and area are known. Since the area of a circle is given by this formula:

$A = \pi r^2\,\!$

If a circle with radius r is drawn with its center at the point (0, 0), any point whose distance from the origin is less than r will fall inside the circle. The Pythagorean theorem gives the distance from any point (x, y) to the center:

$d=\sqrt{x^2+y^2}.\!$

Mathematical "graph paper" is formed by imagining a 1×1 square centered around each cell (x, y), where x and y are integers between −r and r. Squares whose center resides inside or exactly on the border of the circle can then be counted by testing whether, for each cell (x, y),

$\sqrt{x^2+y^2} \le r.\!$

The total number of cells satisfying that condition thus approximates the area of the circle, which then can be used to calculate an approximation of π. Closer approximations can be produced by using larger values of r.

Mathematically, this formula can be written:

$\pi = \lim_{r \to \infty} \frac{1}{r^2} \sum_{x=-r}^{r} \; \sum_{y=-r}^{r} \begin{cases} 1 & \text{if } \sqrt{x^2+y^2} \le r \\ 0 & \text{if } \sqrt{x^2+y^2} > r \end{cases}$

In other words, begin by choosing a value for r. Consider all cells (x, y) in which both x and y are integers between −r and r. Starting at 0, add 1 for each cell whose distance to the origin (0,0) is less than or equal to r. When finished, divide the sum, representing the area of a circle of radius r, by r² to find the approximation of π. For example, if r is 5, then the cells considered are:

(−5,5)	(−4,5)	(−3,5)	(−2,5)	(−1,5)	(0,5)	(1,5)	(2,5)	(3,5)	(4,5)	(5,5)
(−5,4)	(−4,4)	(−3,4)	(−2,4)	(−1,4)	(0,4)	(1,4)	(2,4)	(3,4)	(4,4)	(5,4)
(−5,3)	(−4,3)	(−3,3)	(−2,3)	(−1,3)	(0,3)	(1,3)	(2,3)	(3,3)	(4,3)	(5,3)
(−5,2)	(−4,2)	(−3,2)	(−2,2)	(−1,2)	(0,2)	(1,2)	(2,2)	(3,2)	(4,2)	(5,2)
(−5,1)	(−4,1)	(−3,1)	(−2,1)	(−1,1)	(0,1)	(1,1)	(2,1)	(3,1)	(4,1)	(5,1)
(−5,0)	(−4,0)	(−3,0)	(−2,0)	(−1,0)	(0,0)	(1,0)	(2,0)	(3,0)	(4,0)	(5,0)
(−5,−1)	(−4,−1)	(−3,−1)	(−2,−1)	(−1,−1)	(0,−1)	(1,−1)	(2,−1)	(3,−1)	(4,−1)	(5,−1)
(−5,−2)	(−4,−2)	(−3,−2)	(−2,−2)	(−1,−2)	(0,−2)	(1,−2)	(2,−2)	(3,−2)	(4,−2)	(5,−2)
(−5,−3)	(−4,−3)	(−3,−3)	(−2,−3)	(−1,−3)	(0,−3)	(1,−3)	(2,−3)	(3,−3)	(4,−3)	(5,−3)
(−5,−4)	(−4,−4)	(−3,−4)	(−2,−4)	(−1,−4)	(0,−4)	(1,−4)	(2,−4)	(3,−4)	(4,−4)	(5,−4)
(−5,−5)	(−4,−5)	(−3,−5)	(−2,−5)	(−1,−5)	(0,−5)	(1,−5)	(2,−5)	(3,−5)	(4,−5)	(5,−5)

This circle as it would be drawn on a Cartesian coordinate graph. The cells (±3, ±4) and (±4, ±3) are labeled.

The 12 cells (0, ±5), (±5, 0), (±3, ±4), (±4, ±3) are exactly on the circle, and 69 cells are completely inside, so the approximate area is 81, and π is calculated to be approximately 3.24 because 81 / 5² = 3.24. Results for some values of r are shown in the table below:

r	area	approximation of π
2	13	3.25
3	29	3.22222
4	49	3.0625
5	81	3.24
10	317	3.17
20	1257	3.1425
100	31417	3.1417
1000	3141549	3.141549

Similarly, the more complex approximations of π given below involve repeated calculations of some sort, yielding closer and closer approximations with increasing numbers of calculations.

[edit] Approximation with a regular polygon

π is defined as the ratio of the circumference of a circle to its diameter. Circles can be approximated as regular polygons with an increasing number of sides, approaching infinity. Archimedes used this method with a 96-sided polygon to show that π is between 223/71 and 22/7.

[edit] Continued fractions

Besides its simple continued fraction representation [3; 7, 15, 1, 292, 1, 1, …], which displays no discernible pattern, π has many generalized continued fraction representations generated by a simple rule, including these two.

$\pi= {3 + \cfrac{1^2}{6 + \cfrac{3^2}{6 + \cfrac{5^2}{6 + \cfrac{7^2}{\ddots\,}}}}}\!$

$\pi = \cfrac{4}{1 + \cfrac{1^2}{3 + \cfrac{2^2}{5 + \cfrac{3^2}{7 + \cfrac{4^2}{\ddots}}}}}\!$

(Other representations are available at The Wolfram Functions Site.)

[edit] Trigonometry

[edit] Gregory-Leibniz series

The Gregory-Leibniz series

$\frac{\pi}{4} = \frac{1}{1} - \frac{1}{3} + \frac{1}{5} - \frac{1}{7} + \cdots + \frac{(-1)^n}{2n+1} + \cdots\!$

is the power series for arctan(x) specialized to x = 1. It converges too slowly to be of practical interest. However, the power series converges much faster for smaller values of $x$ , which leads to formulas where $π$ arises as the sum of small angles with rational tangents, such as these two by John Machin:

$\frac{\pi}{4} = \arctan\frac{1}{2} + \arctan\frac{1}{3}\!$

$\frac{\pi}{4} = 4 \arctan\frac{1}{5} - \arctan\frac{1}{239}\!$

Formulas for π of this type are known as Machin-like formulae.

[edit] Arctangent

Knowing that $\arctan(1) \cdot 4 = \pi\!$ the formula can be simplified to get:

$\cfrac {\pi}{2} = \sum_{n=0}^{\infty} \cfrac {n!}{(2n + 1)!!} = 1 + \cfrac{1}{3} + \cfrac{1\cdot2}{3\cdot5} + \cfrac{1\cdot2\cdot3}{3\cdot5\cdot7}+\cdots\!$

See: Double Factorial

[edit] Arcsine

Observing an equilateral triangle and noting that

$\sin\left (\frac{\pi}{6}\right )=\frac{1}{2}\!$

yields

$\pi = 3 \sum_{n=0}^\infty {2n \choose n} \frac{1}{(2n+1) 16^n} = 3 + \frac{1}{8} + \frac{9}{640} + \frac{15}{7168} + \cdots\!$

with a convergence such that each additional 5 terms yields at least 3 more digits.

[edit] The Salamin–Brent algorithm

The Gauss–Legendre algorithm or Salamin–Brent algorithm was discovered independently by Richard Brent and Eugene Salamin in 1975. This can compute pi to N digits in time proportional to N log(N) log(log(N)), much faster than the trigonometric formulae.

[edit] Digit extraction methods

[edit] BBP formula (base 16)

The Bailey–Borwein–Plouffe formula (BBP) for calculating pi was discovered in 1995 by Simon Plouffe. The formula computes pi in base 16 without needing to compute the previous digits (digit extraction). ^[1]

$\pi=\sum_{n=0}^\infty \left(\frac{4}{8n+1}-\frac{2}{8n+4}-\frac{1}{8n+5}-\frac{1}{8n+6}\right)\left(\frac{1}{16}\right)^n\!$

[edit] Bellard's improvement (base 2)

An alternative formula for computing pi in base 2 was derived by Fabrice Bellard. This O(n²) algorithm is an improvement of the O(n³log(n)³) algorithm, and has been measured to make computing binary digits of pi 43% faster. ^[2]

$\pi=\frac{1}{2^6}\sum_{n=0}^\infty \frac{(-1)^n}{2^{10n}} \left (-\frac{2^5}{4n+1}-\frac{1}{4n+3}+\frac{2^8}{10n+1}-\frac{2^6}{10n+3}-\frac{2^2}{10n+5}-\frac{2^2}{10n+7}+\frac{1}{10n+9}\right )\!$

[edit] Extending to arbitrary bases

In 1996, Simon Plouffe derived an algorithm to extract the nth digit of pi in an arbitrary base in O(n³log(n)³) time. ^[3]

[edit] Improvement using the Gosper formula

In 1997, Fabrice Bellard improved Plouffe's formula for digit-extraction in an arbitrary base to reduce the runtime to O(n²). ^[4]

[edit] Efficient methods

In the early years of the computer, the first expansion of π to 100,000 decimal places was computed by Maryland mathematician Dr. Daniel Shanks and his team at the United States Naval Research Laboratory (N.R.L.) in 1961.

Daniel Shanks and his team used two different power series for calculating the digits of π. For one it was known that any error would produce a value slightly high, and for the other, it was known that any error would produce a value slightly low. And hence, as long as the two series produced the same digits, there was a very high confidence that they were correct. The first 100,000 digits of π were published by the Naval Research Laboratory.

None of the formulæ given above can serve as an efficient way of approximating π. For fast calculations, one may use a formula such as Machin's:

$\frac{\pi}{4} = 4 \arctan\frac{1}{5} - \arctan\frac{1}{239}\!$

together with the Taylor series expansion of the function arctan(x). This formula is most easily verified using polar coordinates of complex numbers, starting with

$(5+i)^4\cdot(-239+i)=-114244-114244i.\!$

Formulae of this kind are known as Machin-like formulae.

Many other expressions for π were developed and published by Indian mathematician Srinivasa Ramanujan. He worked with mathematician Godfrey Harold Hardy in England for a number of years.

Extremely long decimal expansions of π are typically computed with the Gauss-Legendre algorithm and Borwein's algorithm; the Salamin-Brent algorithm which was invented in 1976 has also been used.

The first one million digits of π and ¹/_π are available from Project Gutenberg (see external links below). The record as of December 2002 by Yasumasa Kanada of Tokyo University stands at 1,241,100,000,000 digits, which were computed in September 2002 on a 64-node Hitachi supercomputer with 1 terabyte of main memory, which carries out 2 trillion operations per second, nearly twice as many as the computer used for the previous record (206 billion digits). The following Machin-like formulæ were used for this:

$\frac{\pi}{4} = 12 \arctan\frac{1}{49} + 32 \arctan\frac{1}{57} - 5 \arctan\frac{1}{239} + 12 \arctan\frac{1}{110443}\!$

K. Takano (1982).

$\frac{\pi}{4} = 44 \arctan\frac{1}{57} + 7 \arctan\frac{1}{239} - 12 \arctan\frac{1}{682} + 24 \arctan\frac{1}{12943}\!$

F. C. W. Störmer (1896).

These approximations have so many digits that they are no longer of any practical use, except for testing new supercomputers. (Normality of π will always depend on the infinite string of digits on the end, not on any finite computation.)

In 1997, David H. Bailey, Peter Borwein and Simon Plouffe published a paper (Bailey, 1997) on a new formula for π as an infinite series:

$\pi = \sum_{k = 0}^{\infty} \frac{1}{16^k} \left( \frac{4}{8k + 1} - \frac{2}{8k + 4} - \frac{1}{8k + 5} - \frac{1}{8k + 6}\right).\!$

This formula permits one to fairly readily compute the k^th binary or hexadecimal digit of π, without having to compute the preceding k − 1 digits. Bailey's website contains the derivation as well as implementations in various programming languages. The PiHex project computed 64-bits around the quadrillionth bit of π (which turns out to be 0).

Fabrice Bellard claims to have beaten the efficiency record set by Bailey, Borwein, and Plouffe with his formula to calculate binary digits of π [1]:

$\pi = \frac{1}{2^6} \sum_{n=0}^{\infty} \frac{{(-1)}^n}{2^{10n}} \left( - \frac{2^5}{4n+1} - \frac{1}{4n+3} + \frac{2^8}{10n+1} - \frac{2^6}{10n+3} - \frac{2^2}{10n+5} - \frac{2^2}{10n+7} + \frac{1}{10n+9} \right)\!$

Other formulæ that have been used to compute estimates of π include:

$\frac{\pi}{2}=\sum_{k=0}^\infty\frac{k!}{(2k+1)!!}=1+\frac{1}{3}\left(1+\frac{2}{5}\left(1+\frac{3}{7}\left(1+\frac{4}{9}(1+\cdots)\right)\right)\right)\!$

Newton.

$\frac{1}{\pi} = \frac{2\sqrt{2}}{9801} \sum^\infty_{k=0} \frac{(4k)!(1103+26390k)}{(k!)^4 396^{4k}}\!$

Srinivasa Ramanujan.

This converges extraordinarily rapidly. Ramanujan's work is the basis for the fastest algorithms used, as of the turn of the millennium, to calculate π.

$\frac{1}{\pi} = 12 \sum^\infty_{k=0} \frac{(-1)^k (6k)! (13591409 + 545140134k)}{(3k)!(k!)^3 640320^{3k + 3/2}}\!$

David Chudnovsky and Gregory Chudnovsky.

[edit] Projects

[edit] Pi Hex

Pi Hex was a project to compute three specific bits of π using a distributed network of several hundred computers. In 2000, after two years, the project finished computing the five trillionth (10¹²), the forty trillionth, and the quadrillionth (10¹⁵) bits. All three of them turned out to be 0.

[edit] Background pi

Inspired by Pi Hex and Project Pi, Background Pi seeks to compute decimal digits of pi sequentially. The project has computed over a hundred thousand digits using spare CPU cycles. Background Pi is oriented to be more for an average end user than for a power user by offering an unobtrusive user interface. Research is underway on the efficiency of converting computed hex digits to decimal as computing hex digits is faster than computing decimal. A new version is in development that would manage multiple computation projects in a friendlier interface than BOINC.

[edit] Software for calculating π

Over the years, several programs have been written for calculating pi (π) to many digits on personal computers.

[edit] General-purpose

Most computer algebra systems can calculate π and other common mathematical constants to any desired precision.

Functions for calculating π are also included in many general libraries for arbitrary-precision arithmetic, for instance CLN and MPFR.

An example of a computer script (written in PHP) has been made to calculate pi and is below. It is designed to be able to follow Viète's formula for calculating pi and if this script is given enough time to run, it can display 1000 digits of pi.

<?
set_time_limit(0);
echo "                                                                                           
 
                                   <br>";
$digits=1000;
$pival='1';
$str1='';
while ($pirow<2000)
    {
    $pirow+=1;
    $pisubrow=0;
    $tempval='0';
    while ($pisubrow<$pirow)
        {
        $tempval=bcsqrt(bcadd('2',$tempval,$digits),$digits);
        $pisubrow+=1;
        sleep(1);
        }
    $tempval=bcdiv('2',$tempval,$digits);
    $pival=bcmul($tempval,$pival,$digits);
    unset($tempval);      
    $tempval=bcmul('2',$pival,$digits);
    $common = array();
    $length = (strlen($str1) >= strlen($tempval)) ? strlen($str1) : strlen($tempval);
    for($x=0;$x<$length;$x++){
        if($str1[$x] == $tempval[$x]){
            $common[] = $str1[$x]; 
        }else{
            break;
        }
    }
    if ($pirow>5)
        {
        echo "pi=".substr(join('',$common),0,-2);
        }
    unset($str1);
    $str1=$tempval;    
    unset($tempval);
    echo "<p>";
    flush();
    sleep(1);
    }
$pival=bcmul($pival,'2',$digits);
unset($pirow);
?>

The above PHP script will need to be given at least 2 minutes for a decent result to be displayed and of course, this script can continue to run for days with no errors. Also if you are not worried about cpu usage in this script then you can remove the lines with the sleep(); function for a faster result.

[edit] Special-purpose

Programs designed for the specific purpose of calculating π may have better performance than general-purpose mathematical software. They typically implement checkpointing and efficient disk swapping to facilitate extremely long-running and memory-expensive computations.

PiFast, by Xavier Gourdon was the fastest program for Microsoft Windows in 2003. According to its author, it can compute one million digits in 3.5 seconds on a 2.4 GHz Pentium 4.^[5] PiFast can also compute other irrational numbers like e and √2. It can also work at lesser efficiency with very little memory (down to a few tens of megabytes to compute well over a billion (10⁹) digits). This tool is a popular benchmark in the overclocking community. PiFast 4.4 is available from Stu's Pi page. PiFast 4.3 is available from Gourdon's page.

QuickPi by Steve Pagliarulo for Windows is faster than PiFast for runs of under 400,000,000 digits. Version 4.5 is available on Stu's Pi Page below. Like PiFast, QuickPi can also compute other irrational numbers like e, √2, and √3. The software may be obtained from the Pi-Hacks Yahoo! forum, or from Stu's Pi page.

[edit] See also

[edit] References

^ MathWorld: BBP Formula http://mathworld.wolfram.com/BBPFormula.html
^ Bellard's Website: http://fabrice.bellard.free.fr/pi/pi_bin/pi_bin.html
^ Simon Plouffe, On the computation of the n'th decimal digit of various transcendental numbers, November 1996
^ Bellard's Website: http://fabrice.bellard.free.fr/pi/pi_n2/pi_n2.html
^ PiFast timings

[edit] External links