How is the value of $\pi$ ( Pi ) actually calculated?

Solution 1:

There is the famous formula $$\dfrac{\pi}{4} = 1 - \dfrac{1}{3} + \dfrac{1}{5} - \dfrac{1}{7} \pm \cdots.$$ Unfortunately, this converges rather slowly. (If you compute to $n$ terms and multiply by $4$, you get within roughly $1/n$ of the value of $\pi$; so to get 3 decimal places of accuracy, you'll need to sum something like the first thousand terms of the series.)

There is the Machin formula $$\dfrac{\pi}{4} = 4\arctan \dfrac{1}{5} - \arctan \dfrac{1}{239},$$ which can be combined with an infinite series formula for $\arctan$ to provide a much more rapidly convergent infinite series which can be used to compute $\pi$ to many digits of accuracy.

The Wikipedia entry provides more details about these and other methods, both historical and contemporary, for computing $\pi$.

Solution 2:

There is an excellent list of ways that $\pi$ has been computed throught history here. In this answer, I will explain a variant of the method that Archimedes used to compute $\pi$.

Consider an isosceles triangle $ABC$ with $AB = AC = 1$. We will start angle $BAC$ at $60^\circ$ and repeatedly halve it, all the while keeping track of the length $BC$. After we have halved the angle $n$ times:

  • The angle $BAC$ is $60^\circ / 2^n$.
  • Therefore, we can fit $6 \cdot 2^n$ copies of the triangle like slices of cake with $A$ being the centre of the cake.
  • Let us call the length $BC$, which we will compute below, $L_n$.
  • The length of the outside of the cake will therefore be $6 \cdot 2^n L_n$.

The outside of the cake is roughly a circle of radius 1, and as $n$ increases it becomes more and more like a circle. The circumference of the circle is $2\pi$. Therefore, this gives us a way to compute $\pi$. Specifically, $\pi \approx 3 \cdot 2^n L_n$.

So how do we compute $L_n$? The starting case is easy. In this case the triangle is equilateral, so $L_0 = 1$. If we know $L_n$ we can compute $L_{n+1}$ using Pythagoras's theorem twice, as explained below. This allows us to compute $L_n$ for any $n$.

First, we bisect the angle $BAC$ to make two right-angled triangles that meet along a line $AP$. We know $AB = 1$ and $BP = L_n / 2$ so by Pythagoras's theorem $AP = \sqrt{1 - L_n^2/4}$.

Next, we extend $AP$ to $Q$ so that the distance $AQ = 1$. Therefore, the distance $PQ = 1 - \sqrt{1 - L_n^2/4}$.

Next, we join $B$ to $Q$ to make a right-angled triangle $BPQ$. We know $BP = L_n / 2$ and $PQ = 1 - \sqrt{1 - L_n^2/4}$ so by Pythagoras's theorem,

$$BQ = \sqrt{(L_n / 2)^2 + \left(1 - \sqrt{1 - L_n^2/4}\right)^2}$$ $$= \sqrt{L_n^2/4 + 1 - 2\sqrt{1 - L_n^2/4} + (1 - L_n^2/4)}$$ $$= \sqrt{2 - \sqrt{4 - L_n^2}}$$

But now $ABQ$ is an isosceles triangle and angle $BAQ$ is half angle $BAC$, so $BQ = L_{n+1}$.

You will notice that in this algorithm we sometimes take a square-root and then immediately square it again as part of the next step, and we sometimes subtract something from 2 and then immediately subtract the result from 4. We can avoid these inefficiencies by defining $M_n = 4 - L_n^2$.

Using that trick, a summary of the algorithm is as follows:

  • $M_0 = 3$
  • $M_{n+1} = 2 + \sqrt{M_n}$
  • $3 \cdot 2^n \sqrt{4 - M_n}$ approaches $\pi$ as $n$ increases.

Here it is in Python:

m = 3
for n in range(1, 11):
     m = 2 + sqrt(m)
     print 3*(2**n)*sqrt(4-m)

Here is the output:

3.10582854123
3.13262861328
3.13935020305
3.14103195089
3.14145247229
3.14155760791
3.14158389215
3.14159046324
3.14159210604
3.14159251659

Unfortunately the numerical accuracy breaks down soon after this point, but it at least shows that the algorithm works in principle.

When Archimedes did it, he was not very good a computing square-roots, and he had to introduce extra approximations for that reason. He got as far as $n=4$, and his answer was $$3\frac{10}{71} = \frac{223}{71} \approx 3.140845070422535$$ which is slightly less than the result we get for $n=4$, which is $3.14103195089$. Archimedes knew that this answer was too small. He also used a similar method to compute an answer that he knew was too big: $$3\frac{1}{7} = \frac{22}{7} \approx 3.142857142857143$$

That's where it comes from.

Solution 3:

When I was a child I was taught this funny experiment:

Bring a thread of known length $l$ and wrap it as a circle. You can easily get the diameter $d$ (measure it) then $\pi = \frac{l}{d}$

Repeat the experiment and you will get approximately equal values.

Solution 4:

Here's $\pi$ accurate to 48 decimal digits: $$ \mathtt{3.141592653589793238462643383279502884197169399375} $$

Archimedes' method can be implemented in Python to compute 48 digits of $\pi$. Python provides long integers which can be used to represent fixed-point numbers with as much precision as required. For example, the fixed-point number $1.000+(375\times10^{-24})$ can be represented by the long integer $\mathtt{1000000000000000000000375L}$.

$\pi$ is defined as the ratio of a circle's circumference $C$ to its diameter $D$. That is, $\pi=C/D$. Equivalently, $\pi$ may also be defined as the ratio of a circle's hemicircumference $H=C/2$ to its radius $r=D/2$. That is, $\pi=C/D=H/r$. And with $r=1$, $\pi=H$.

Archimedes' method is an iterative squeeze algorithm. The unit circle ($r=1$) is squeezed between an inscribed regular polygon and a circumscribed regular polygon (both initially regular hexagons). The circle has hemicircumference $H$. The $\mathrm{n^{th}}$ inscribed polygon has hemiperimeter $I_n$ which is slightly less than $H$. The $\mathrm{n^{th}}$ circumscribed polygon has hemiperimeter $C_n$ which is slightly greater than $H$. And since $\pi=H$, $$I_n < \pi < C_n$$

On each iteration of the algorithm, the number of sides of each polygon is doubled, so that $I_n$ increases a little bit and $C_n$ decreases a little bit. In effect, the hemicircumference $H$ is squeezed tighter and tighter between the two hemiperimeters $I_n$ and $C_n$. In the limit as $n$ approaches infinity, both $I_n$ and $C_n$ converge to $\pi$. $$ \lim_{n \to \infty} \quad I_n = \pi = C_n $$

$\mathtt{apt1002}$ provided a derivation for the $\mathrm{n^{th}}$ inscribed polygon hemiperimeter $I_n$. The iterative operation $M_n=2+\sqrt{M_{n-1}}$ was defined, with $M_0=3$. Also, $M_n$ was specified in terms of the length $L_n$ of the side $BC$ of the $\mathrm{n^{th}}$ inscribed polygon such that $M_n=4-L_n^2$. Solving for $L_n$ we have $L_n = \sqrt{4 - M_n}$. Therefore, the $\mathrm{n^{th}}$ inscribed polygon hemiperimeter is $$I_n = 3 \cdot 2^n \cdot L_n $$

Now we provide a derivation for the $\mathrm{n^{th}}$ circumscribed polygon hemiperimeter $C_n$. Following $\mathtt{apt1002}$'s narrative, we extend $AC$ to $D$ such that angle $QDA$ equals angle $PCA$. This gives us two congruent right triangles $AQD$ and $APC$, with $AQ=1$ and $AP=\sqrt{1-L_n^2/4}$. Let $K_n$ be the length of a side of the $\mathrm{n^{th}}$ circumscribed polygon. Then $QD=K_n/2$ and $PC=L_n/2$. Since the two triangles are congruent, the ratios of their sides are equal. $$ \frac{QD}{AQ}=\frac{PC}{AP} \quad\to\quad \frac{K_n/2}{1}=\frac{L_n/2}{\sqrt{1-L_n^2/4}} %% \quad\to\quad %% K_n = \frac{L_n}{\frac{1}{2}\sqrt{4-L_n^2}} = \frac{2 \cdot L_n}{\sqrt{M_n}} $$

Solving for $K_n$ we have $K_n = 2 \cdot L_n / \sqrt{M_n}$. Therefore, the $\mathrm{n^{th}}$ circumscribed polygon hemiperimeter is $$ C_n = 3 \cdot 2^n \cdot K_n $$

The Python code below specifies three functions. The first function $\mathtt{huge\_int}$ creates a huge integer. It's used to initialize constants and variables, and to scale other variables. Given integer $n$ and exponent $x$, it returns a long integer equal to the value $n \times 10^x$.

The function $\mathtt{huge\_sqrt}$ calculates the square root of a huge integer. Given a long integer $N$ and an exponent $x$, it returns a long integer $\sqrt{N}$. It uses Newton's method, which is derived from the identity $2r^2 = r^2 + r^2$. We replace the second $r^2$ term on the right with $N$, divide both sides by $2r$, factor out the $1/2$ term, and replace each $r$ on the right with $r'$ (which signifies the previous value of $r$). This gives us the iteration equation for root $r$. $$ r = \frac{1}{2} \left( r' + \frac{N}{r'} \right) $$

The main function $\mathtt{archimedes}$ uses Archimedes' method to calculate $\pi$ accurate to $d$ decimal-digits. First it initializes the working exponent $x$, the long integer constants $\mathtt{i4}$ and $\mathtt{i2}$, the long integer variables $M_n$, $I_n$, and $C_n$, and the iteration index $n$. Then it enters the $\mathtt{while}$ loop and iterates until $I_n$ and $C_n$ converge to $\pi$.

$$ $$

def huge_int(n, x):
    return(n * 10**x)

def huge_sqrt(N, x):
    rp = 0
    r = N // 2
    N = N * huge_int(1, x)
    while r != rp:
        rp = r
        r = (rp + (N // rp)) // 2
    return(r)

def archimedes(d):
    x = 2*d+6
    i4 = huge_int(4, x)
    i2 = huge_int(2, x)
    Mn = huge_int(1, x)
    In = 1
    Cn = 2
    n = 0
    while In != Cn:
        Mn = i2 + huge_sqrt(Mn, x)
        Ln = huge_sqrt(i4-Mn, x)
        Kn = 2*Ln*huge_int(1, x) // huge_sqrt(Mn, x)
        In = 3*(2**n)*Ln // huge_int(1, d+6)
        Cn = 3*(2**n)*Kn // huge_int(1, d+6)
        print n, 6*2**n
        print Cn
        print In
        print
        n = n+1
    print "DONE"

$$ $$

Here's part of the output for $\mathtt{archimedes(48)}$. On each iteration, four lines are printed. The first line is the iteration number $n$ along with the number of sides of both polygons. The second and third lines are the $d$-digit values of $C_n$ and $I_n$. Archimedes originally approximated $\pi$ with a 96-side polygon. In order to get an approximation of $\pi$ accurate to 48 decimal digits, we need a polygon with $6 \times 2^{80} \approx 7.25 \times 10^{24}$ sides.

0 6
3464101615137754587054892683011744733885610507620
3000000000000000000000000000000000000000000000000

1 12
3215390309173472477670643901929531596686336954275
3105828541230249148186786051488579940188826815839

2 24
3159659942097500483316634977833210486227538453357
3132628613281238197161749469491736244649776915481

3 48
3146086215131434971098098794237254156265359587343
3139350203046867207135146821208421189150350589362

4 96
3142714599645368298168859093772123871000969091511
3141031950890509638111352926459660107036412216162

.
.
.

78 1813388729421943762059264
3141592653589793238462643383279502884197169399378
3141592653589793238462643383279502884197169399373

79 3626777458843887524118528
3141592653589793238462643383279502884197169399375
3141592653589793238462643383279502884197169399374

80 7253554917687775048237056
3141592653589793238462643383279502884197169399375
3141592653589793238462643383279502884197169399375

DONE