Pi to a Million (or a Billion or a Trillion) Decimal Places

This combines two subjects, other than astronomy, in which I'm very interested: mathematics and computing.

A Pi Pie!

Have you ever read about pi (the ratio of the circumference of a circle to it's diameter) being calculated to so many billion decimal places and wondered how they do it?

Well maybe you haven't but I have!

There's two parts to it: obviously you need some sort of formula from which to work out the value of pi, and you also need to be able to perform arithmetic to a million, or a billion, or however many, decimal places.

A Formula!

The first bit is fairly easy: we use the following infinite series expression for the inverse tangent (or acrtangent) which is thought to have been first discovered by the Indian mathematician Madhava of Sangamagrama:-

For any value z this series will give us the angle (in radians) which has that values as it's tangent (i.e. it will give us the angle x such that tan(x)=z).

Now the tan of 45 degrees is one, but in radians 45 degrees is , so if we let z=1 in the above formula we get:-

As with any infinite series, we don't actually carry on adding terms indefinitely - we only use enough terms to give us the accuracy we require. Unfortunately, this formula isn't very useful because the terms get smaller and smaller very slowly. It  would take 500 thousand terms just to get an accuracy of only 5 decimal places.

A formula which converges more quickly was discovered by the English mathematician and astronomer John Machin in 1706:-

While we need to work out the arctan of two numbers, because they are much smaller than one, the terms converge much more quickly.  For example, for arctan(1/5) each term is less than 1/25th of the previous one.

Modern methods to calculate pi to a large number of decimal places use formula which contain even more terms, such as the following:-

This formula was discovered by the Japanese mathematician Yasumasa Kanada who used it in 2002 to calculate pi to over a trillion decimal places.

Notice that the largest value for which we need to calculate the arctangent is 1/49.  Because this number is almost ten times smaller than 1/5 the series will converge much more quickly. 

Trillion Decimal Place Arithmetic

Any computer can do integer arithmetic on integers up to a certain size.  That size depends on the particular architecture of the computer.  Most computer manufacturers these days produce processors that support 64-bit numbers and arithmetic - that's large enough to hold an 18-digit decimal (whole) number.

(Of course they can also do floating point arithmetic, but that's not what we're going to use.)

We use a large number of integers each holding a fixed number (but less than 18) of decimal digits.  For example, if we want to work out pi to a million decimal places then we can use 100,000 integers, each holding 10 digits.

Of course, we need to store more just one of these numbers.  If you look at Machin's formula, or the last formula I gave above, you can see that we need to be able to add and subtract these numbers and to multiply them by "normal" sized numbers such as 4, 12, 32 etc.

Addition is fairly easy: you just add the decimal integers going from right to left.  If the result of adding two integers is greater than 9,999,999,999 we subtract 10,000,000,000 and carry one:-

Subtraction would be done in a similar way. Multiplication by a "normal" sized number is also done from right to left.

The other thing which we need to do is calculate the arctangent of a number to (say) a million decimal places.  However, all formula involve expressions of the form:-

where n is an integer.  This makes things easier; the series for arctangent I showed earlier can be re-written with z replaced by1/n:-

The only new operation we have here is division by a "normal" sized number.  (Update: if we used the origional formula, we would have to multiply two billion decimal place numbers to gether and that would be much less efficient.)  This division is done a bit like multiplication, only proceeding from right to left.  Any remainder from one division operation being carried to the right.

To actually perform this calculation, we would need to keep both the running total and the current term as 1 million digit numbers.

The next term is more efficiently calculated from the previous term than from scratch.  For example, the fourth term above in the equation above: , can be calculated from the third term: ,  by multiplying by 5, dividing by 7 and then dividing again by n squared.

Trivia

 Here is a visual representation of the first four million digits of pi.

Find your date of birth in pi here.

Buy pi related t-shirts here or here! (No, I don't get any commission.)

Update (20/11/12)

• I wrote the program and got it to work (correctly!) up to 100,000 digits.  It took about 2 minutes on my PC - in 1961 this took 8.7 hours on an IBM 7090.   I got bored with it with after that.

• An alternative formula for evaluating pi is the Chudnovsky algorithm:-

• According to the novel Contact by Carl Sagan, after a huge number of digits of pi, they eventually stop being random.  Instead you get a sequence of zeros and ones which, when plotted on a 121 X 121 array form a circle, but only if you calculate pi in base 11!  This is science fiction.  Although they haven't proved it, mathematiciens believe that pi is a normal number.