Pi has always been attractive for humans, especially mathematicians, even from the ancient Greece. 300 years before Christ, Euclides proposed to approximate a circle by regular polygones with increasing number of edges, to prove that the ratio of the surface of a circle to its radius was a constant. Around 50 years later, Archimedes provided the first approximate value of this number, and he gave it as a fraction!

3 + 10/71 < pi < 3 + 1/7

During the Middle Ages, not much improvement was done in Europe, but many mathematicians worked on Pi in Asia such as Aryabhata (476-550) in India or Tsu Chung-Chi (430-501) in China, and Al-Kashi (1390-1450) in Persia.

1500 years ago, Aryabhata also approximated Pi as a fraction : 62832/20000 = 3,1416
Roughly at the same time, Tsu Chung-Chi provided Pi with 6 correct decimal places,
600 years ago, the great scientist Al-Kashi provided an approximation with 16 correct digits : 3,1415926535897932 ! I guess most of us don't know that much...

In France, René Descartes (1596-1650) proposed, with no demonstration, a formula in the form of an infinite series of products and fractions of square roots (that you can see above). For a mathematician, this formula is really nice, as it involves only one number, one operation (square root) and it can be written easily although it is quite complex.

Alas! In 1882, Carl-Louis-Ferdinand von Lindemann, a German mathematician, demonstrated that Pi is a transcendental number. This means it cannot be the root of a polynomial with integer coefficients. Think about it for a while. You've learned about polynomials in high school: a polynomial is a function of a variable 'x' which can be expressed as

P(x) = a x**n + b x**(n-1) + ... + k x + m

'a', 'b' ... 'm' are the coefficients. 'n' is called the degree of the polynomial. For example : P(x) = x**2 - x - 1.

The equation P(x) = 0 sometimes has 2 roots (solutions), sometimes only one and sometimes no solution at all (for real numbers). The root of the example polynomial is the golden ratio (roughly 1,618039), another number that attracts many people.

Lindemann demonstrated that you cannot find any polynomial of any degree (as high as you want) with integer coefficients of which Pi is a root. Such a number is called a transcendental number.

As such, Pi cannot be expressed as a fraction. So the values provided by Archimedes and Aryabhata are, and will always be, approximations.

So what? Should we stay there and change subjects? No, let's talk about another mathematical tool : continued fractions. And see how it relates to Pi...

Step 1: What Are Continued Fractions

Continued fractions are just fractions made of fractions. Every number, rational or irrational, can be written as a continued fraction. I won't go into the theory related to continued fractions (CF), as you can find much more information about them on Wikipedia for instance.

An example of continued fraction is shown in the image above. The CF is merely defined by its coefficients a0, a1, a2 ... an. So it's a compact way of writing a number.

The continued fraction of Pi is also given above : its coefficients are 3, 7, 15, 1, 292, 1, 1, 1, 2, 1, 3, 1... This sequence of numbers is registered in the OEIS (Online Encyclopedia of Integer Sequences) under the number A001203.

I wrote a short Python script to compute the coefficient of the CF for a given number (Pi in our case) :

import mpmathfrom mpmath import *mp.dps = 100 # decimal placesPI = mpmathify(pi)print ('Continued fraction of ')## Compute CFsteps = 30x = PIprint(x)print("=")a = mpmath.floor(x)print (int(a))for i in range (steps):    b = fsub(x, a)    x = fdiv(1, b)    a = mpmath.floor(x)    print(int(a))

It uses the mpmath library, which enables increased accuracy. This ensures that all 30 values are correct.

When running the script, the output is :

Continued fraction of
3.141592653589793238462643383279502884197169399375105820974944592307816406 = 3 7 15 1 292 1 1 1 2 1 3 1 14 2 1 1 2 2 2 2 1 84 2 1 1 15 3 13 1 4 2

All the numbers 3, 7, 15, 1 etc are present. But when examining this sequence, it seems rather random.

Compare this sequence with the sequence of e, Euler's number (change x = PI by x = mpmathify(e) in the code above):

Continued fraction of
2.718281828459045235360287471352662497757247093699959574966967627724076630353547594571382178525166427 = 2 1 2 1 1 4 1 1 6 1 1 8 1 1 10 1 1 12 1 1 14 1 1 16 1 1 18 1 1 20

The pattern can be found easily : 1 2 1, then 1 4 1, then 1 6 1, etc, although it's a transcendental number as well. And it's even simpler with the golden ratio:

1 1 1 1 1 1 1 1 1 1 1 1 1 1...

So, you can see that some numbers can be represented as continued fractions with very simple patterns. The above Python script can help you find some others. A nice one is the square root of 2 : change x = mpmath.sqrt(2) and run the script: you'll see by yourself :-)

But Pi is different... Nobody knows why, but that's what makes its charm...

One concluding remark here: seeing the general shape of a CF on the picture above, we see that the inverse of the number is also present in this representation. Indeed, starting from this CF, we find that the CF of

1 / ( x - a0)

is simply the CF of x starting from the second stage, with a1. We will remember this later.

Step 2: Let's Generalize...

This is where William Brouncker comes in. William Brouncker was an English mathematician, living in the 17th century. Among other things, he invented a generalized continued fraction involving Pi.

Wait a minute : generalized ?

We have seen in the previous step that the standard CF always has a 1 above the fraction. This can be generalized if we use other positive integer numbers as well. We can even generalize one step further by using negative numbers, fractions or even complex numbers, but that's another story...

As far as I know, William Brouncker was the first person to come with a generalized CF, and it involved Pi. He calculated a GCF for Pi / 4, which has coefficients that exhibit a simple pattern:

The number in front of the 'plus' sign is 2 (except the first one which is 1)
The numbers over the fractions are the sequence of the squares of odd numbers: 1, 9, 25, 49, 81...

How can we calculate that? Here is another Python script that will compute this GCF:

#python 3.7.1def CF (n):    # Brouncker continued fraction    # a : number before the plus sign    if n == 0: a = 1    else: a = 2    # b: numerator of the fraction    b = (2 * n - 1) ** 2    return a, bNmax = 500000N = Nmax + 1S = 1for i in range(1,N):  j = Nmax - i  a, b = CF(j)  S = b / (a + S)  # print (str(j)+' a = '+str(a)+' b = '+str(b)) # uncomment if you want to see the coefficientsprint ("Continued Fraction = {0:18.16f}".format(S))

If you run it, the output is:

Continued Fraction = 0.7853976633984483

which is quite close to the real value of Pi/4 ~ 0.785398163397 (5 correct digits, in green)

Increasing Nmax will improve the precision at the cost of computation time.

There are other GCF that can be used to compute Pi. We just need to change the CF(n) function in the script above to compute them. Here is a short collection (from Wikipedia):

Pi / 4 =

The numbers in front of the 'plus' sign are the sequence of odd numbers : 1, 3, 5, 7...
The numbers over the fractions are the sequence of the squared numbers: 1, 1, 4, 9, 16...

The function to use is:

def CF (n):    # number before the plus sign    a = 2 * n + 1    # numerator of the fraction    if n == 0: b = 1    else: b = n ** 2    return a, b

and the result is:

Continued Fraction = 0.7853981633974483

which is more precise than Brouncker's formula.

Pi - 3 =

The numbers in front of the 'plus' sign are always 6
The numbers over the fractions are the sequence of the squared odd numbers: 1, 9, 25, 49, 81...

Here, the function is much simpler:

def CF (n):    # number before the plus sign    a = 6    # numerator of the fraction    b = (2 * n + 1) ** 2    return a, b

and the result is well know now:

Continued Fraction = 0.1415926535897932

with very good precision too.

Pi² / 12 =

The numbers in front of the 'plus' sign are the odd numbers
The numbers over the fractions are the sequence of the 4th power of the numbers: 1, 16, 81, 256, 625...

Here, the function is:

def CF (n):    # number before the plus sign    a = 2*n+1    # numerator of the fraction    if n == 0: b = 1    else: b = n**4    return a, b

and the result is:

Continued Fraction = 0.8224670334221132

which is precise up to 11 digits.

Step 3: Searching for Pi

How can we use this in our endless search of Pi in fractions?

Well, just try and change the CF function and see what comes out:

def CF (n):    # number before the plus sign    if n % 2 == 0: a = 1    else: a = 2    # numerator of the fraction    b = 2 * n + 1    return a, b

Here, I have tried alternating 1 and 2 for 'a' and odd numbers for 'b'. The result is:

Continued Fraction = 0.5792047726384800

But what is this number? How is it related to Pi?

To try to answer this question, I added a very simple search engine to the Python script. It assumes that the CF is the output of a polynomial P(x) for x = Pi, and that the polynomial has integer coefficients in a given range (here between -12 and +12 -- you can change that of course). For simplicity, I supposed that the polynomial is of degree 2 or lower. It's easy to change that for higher degrees, at the cost of longer computing time.

The code is now:

import mathdef CF (n):    # number before the plus sign    if n % 2 == 0: a = 1    else: a = 2    # numerator of the fraction    b = 2 * n + 1    return a, bNmax = 500000N = Nmax + 1S = 1for i in range(1,N):  j = Nmax - i  a, b = CF(j)  S = b / (a + S)  # print (str(j)+' a = '+str(a)+' b = '+str(b))print ("Continued Fraction = {0:18.16f}".format(S))# optimizationdef opti(x, S): # the function representing the polynomial    if x[3] == 0: return 100 # prevents division by 0    return abs((x[0]*math.pi**2 + x[1]*math.pi + x[2]) / x[3] - S)min = -12max = 12thresh = 1e-9# search for a polynomialfound = Falsebest = 100x = (0,1,2,3)for x0 in range(min, max+1):    for x1 in range(min, max+1):        for x2 in range(min, max+1):            for x3 in range(min, max+1):                x = (x0,x1,x2,x3)                y = opti(x,S)                if (y < best):                    best = y                    bestx = x                if (y < thresh):                    found = True                    print("{0:10.7f} = ({1:2d} * pi^2 + {2:2d} * pi + {3:2d}) / {4:2d}, precision  {5:8.5e}".format(S,x[0],x[1],x[2],x[3],y))

Running it provides no answer, because I set the precision at

thresh = 1e-9

Then, my number is not generated by such a polynomial. We can although refine the search by adding the following lines (using scipy's optimization engine):

from scipy import optimizeif (not found):    print ("Best polynomial found: ({0:2d} * pi^2 + {1:2d} * pi + {2:2d}) / {3:2d}, precision  {4:8.5e}".format(bestx[0],bestx[1],bestx[2],bestx[3],best))    print("Expecting {0:4.2e} precision".format(thresh))    # Use this for starting point of optimization    print("Running thorough search...")    optx = optimize.minimize(opti, bestx, args=(S,)) # run optimized search    print("Found: ({0:12.9f} * pi^2 + {1:12.9f} * pi + {2:12.9f}) / {3:12.9f}, precision  {4:8.5e}".format(optx.x[0],optx.x[1],optx.x[2],optx.x[3],optx.fun))

Then, the answer is:

Continued Fraction = 0.5792047726384800
Best polynomial found: (-3 * pi^2 + 8 * pi + 1) / -6, precision 1.40556e-04
Expecting 1.00e-09 precision
Running thorough search...
Found: (-2.999923375 * pi^2 + 8.000024407 * pi + 1.000007788) / -6.000004512, precision 8.55427e-13

A possible approximation would be (3*Pi² - 8*Pi - 1)/6, but it's only precise at 0.00014, so only 3 correct digits. That's not good! Refining the search provides another approximation, but the coefficients are not integer any more.

Just a safety check: if I try the last function of the previous step, it provides the correct result:

Continued Fraction = 0.8224670334221132
0.8224670 = (-1 * pi^2 + 0 * pi + 0) / -12, precision 1.99996e-12
0.8224670 = ( 1 * pi^2 + 0 * pi + 0) / 12, precision 1.99996e-12

My number is correctly recognized as Pi² / 12, with 11 correct digits.

One last remark: as shown earlier, the CF of a number and its inverse are roughly the same. So we don't need to search for the inverse of a polynomial, which eases the search.

Want to search by yourself?

If you'd like to use this script to search for your own GCF, I suggest that you try several patterns and code them in the CF(n) function. But do not use the optimize section, which may be ambiguous and let you think that a GCF is correct when it is not the case.

When you run the script, if you get an answer such as the one just above, with precision less than 1e-10, then Bingo! Bob is your uncle! (*) : you found a new continued fraction involving Pi !

(*) I just learned this expression, and really love it :-)))

You can also decide that a polynomial is not the good fitting function: why not try exponential, logarithms, gamma, Bessel, and so on. I think this field of the Mathematical landscape is rather empty, so fell free to adapt the script accordingly.

In December 2008, a new continued fraction for Pi was proposed, involving radicals:

The numbers in front of the 'plus' sign are the inverse of the integer numbers: 1, 1/2, 1/3, 1/4, etc
The number over the fractions is always 1

If you want to try it, just use:

def CF (n):<br>    # number before the plus sign    a = 1/(n+1)    # numerator of the fraction    b = 1    return a,b

and run the script. The answer is

Continued Fraction = 0.5707963267956820

which is Pi / 2 - 1 with 12 exact decimal places. The search provides:

0.5707963 = ( 0 * pi^2 + -6 * pi + 12) / -12, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + -5 * pi + 10) / -10, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + -4 * pi + 8) / -8, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + -3 * pi + 6) / -6, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + -2 * pi + 4) / -4, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + -1 * pi + 2) / -2, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + 1 * pi + -2) / 2, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + 2 * pi + -4) / 4, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + 3 * pi + -6) / 6, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + 4 * pi + -8) / 8, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + 5 * pi + -10) / 10, precision 7.85483e-13
0.5707963 = ( 0 * pi^2 + 6 * pi + -12) / 12, precision 7.85483e-13

which are 12 ways of confirming that...

Step 4: Conclusion

I hope this will enable you to find new formulas involving Pi, and why not become famous in the Mathematics world!

Oh, I forgot one last formula : use Pi to compute the volume of a thin cylinder. Imagine a very thin colorful cylinder of height 'a' and radius 'z'. Then the volume is simply:

V = Pi z z a

I call it the Yummi formula.

If you would like to know more (and I really mean much more) about Pi, I recommend this website pi314, full of very interesting and intriguing facts and information about your now preferred number...

You can also watch the movie The Life of Pi, which has great CGI effects. Indeed, the life of Pi is very very rich, and we've met in our 'oddissey' famous people ranging from Euclides to Louis Armstrong. Why not join them yourself?