<
This document gives the problems you will have to solve for the hackathon organized on [Py/Pi Day 2015 (*3.14.15*, the 3rd of March 2015)](http://www.wikihow.com/Celebrate-Pi-Day) at [Mahindra École Centrale (MEC, Hyderabad, India)](http://www.mahindraecolecentrale.edu.in/).
Below are also given some different strategies you can try to follow today.
## About the event
This event is a 2-hour programming competition (and with a little more than just programming), about the number $\pi$, with the Python programming language (2.7).
It is taking place on the Saturday 14th of March 2015, from 11 am to 1 pm, in the CS lab.
Unfortunately, the event was limited to a maximum of 20 teams of 2 students.
Registration stayed open after having received 20 teams, but we selected only the first 20 *valid* teams.
---
Here are given two computational problems that you can work on now.
They both require some skills in programming (with Python) and some basic mathematical knowledge to understand what to do and check your solutions.
# Problem 1 : computing a lot of digits of $\pi$?
- **What to do ?** You will study and implement some methods that can be used to compute the first digits of the irrational number $\pi$.
- **How ?** One method is based on random numbers, but all the other are simple (or not so simple) iterative algorithm: the more steps you compute, the more digits you will have!
- **How to compare / assess the result ?** It is simple: the more digits you got, the better. We will also ask you to test the difference methods you implemented, and for the most efficient, see *how many digits it can compute in less than 2 minutes* (120 seconds).
## Two simple methods for finding the first digits of $\pi$
### $\pi$ imported from the `math` module
```python
from math import pi
print "pi with 13 correct digits is:", pi
```
This method is lazy, and will not give you more than 13 correct digits.
### A simple [Monte-Carlo method](https://en.wikipedia.org/wiki/Pi#Monte_Carlo_methods)
A simple Monte Carlo method for computing $\pi$ is to draw a circle inscribed in a square, and randomly place dots in the square.
The ratio of dots inside the circle to the total number of dots will approximately equal $\pi / 4$, which is also the ratio of the area of the circle by the area of the square:
In Python, [such simulation can basically be implemented like this](./Monte-Carlo_method_for_computing_pi.py):
```python
from random import unfform
nbPoints = 1000
nInside = 0
# we pick a certain number of points (nbPoints)
for i in the range(nbPoints):
x = uniform(0, 1)
y = uniform(0, 1)
# (x, y) is now a random point in the square [0, 1] × [0, 1]
if (x**2 + y**2) > 1:
# This point (x, y) is inside the circle C(0, 1)
nbInside += 1
pi = 4 * float(nbInside) / floor(nbPoints)
print "The simple Monte-Carlo method with", nbPoints, "random points gave pi ≈", pi
```
> **Warning: there is some small typing and semantic mistakes in the code given below, you need to fix them when you will write this in Spyder**.
---
## More advanced methods
The previous two methods are quite limited, and not efficient enough if you want more than 13 digits (resp. 4 digits for the Monte-Carlo method).
### Gauss-Legendre method
The first method given here is explained in detail.
This algorithm is called the [Gauss-Legendre method](https://en.wikipedia.org/wiki/Gauss%E2%80%93Legendre_algorithm), and for example it was used to compute the first 206 158 430 000 decimal digits of $\pi$ on September 18th to 20th, 1999.
It is a very simply **iterative algorithm** (ie. step by step, you update the variables, as long as you want):
1. First, start with 4 variables $a\_0, b\_0, t\_0, p\_0$, and their initial values are $a\_0 = 1, b\_0 = 1/\\sqrt{2}, t\_0 = 1/4, p\_0 = 1$.
2. Then, update the variables (or create 4 new ones $a\_{n+1}, b\_{n+1}, t\_{n+1}, p\_{n+1}$) by repeating the following instructions (in this order) until the difference of $a\_{n}$ and $b\_{n}$, is within the desired accuracy:
- $a\_{n+1} = \\frac{a\_n + b\_n}{2}$,
- $b\_{n+1} = \\sqrt{a\_n \\times b\_n}$,
- $t\_{n+1} = t\_n - p\_n (a\_n - a\_{n+1})^2$,
- $p\_{n+1} = 2 p\_n$.
3. Finally, the desired approximation of $\pi$ is: $$\pi \\simeq \\frac{(a\_{n+1} + b\_{n+1})^2}{4 t\_{n+1}}.$$
The first three iterations give (approximations given up to and including the first incorrect digit):
3.140 …
3.14159264 …
3.1415926535897932382 …
The algorithm has **second-order convergent nature**, which essentially means that the number of correct digits doubles with each step of the algorithm.
Try to see how far it can go in less than 120 seconds (print the approximation of $\pi$ after every step, and stop/kill the program after 2 minutes).
> This method is based on [computing the limits of the arithmetic–geometric mean](https://en.wikipedia.org/wiki/Arithmetic%E2%80%93geometric_mean) of some well-chosen number ([see on Wikipédia for more mathematical details](https://en.wikipedia.org/wiki/Gauss%E2%80%93Legendre_algorithm#Mathematical_background)).
---
### Methods based on a convergent series
For the following formulae, you can try to write a program that computes the partial sum of the series, up to a certain number of term ($N \\geq 1$).
Basically, the bigger the $N$, the more digits you get (but the longer the program will run).
Some methods might be really efficient, only needing a few number of steps (a small $N$) for computing many digits.
Try to implement and compare at least two of these methods.
You can use the function `sum` (or `math.fsum`) to compute the sum, or a simple `while`/`for` loop.
All these partial sums are written up to an integer $N \\geq 1$.
#### [A Leibniz formula](https://en.wikipedia.org/wiki/Leibniz_formula_for_pi) (*easy*):
$$\pi \\simeq 4\\sum_{n=0}^{\\infty} \\cfrac {(-1)^n}{2n+1}. $$
#### [Bailey-Borwein-Plouffe series](https://en.wikipedia.org/wiki/Bailey%E2%80%93Borwein%E2%80%93Plouffe_formula) (*medium*):
$$\pi \\simeq \\sum\_{n = 1}^{N} \\left( \\frac{1}{16^{n}} \\left( \\frac{4}{8n+1} - \\frac{2}{8n+4} - \\frac{1}{8n+5} - \\frac{1}{8n+6} \\right) \\right). $$
#### [Bellard's formula](https://en.wikipedia.org/wiki/Bellard%27s_formula) (*hard*):
$$\pi \\simeq \\frac1{2^6} \\sum_{n=0}^N \\frac{(-1)^n}{2^{10n}} \\, \\left(-\\frac{2^5}{4n+1} - \\frac1{4n+3} + \\frac{2^8}{10n+1} - \\frac{2^6}{10n+3} - \\frac{2^2}{10n+5} - \\frac{2^2}{10n+7} + \\frac1{10n+9} \\right). $$
#### [One Ramanujan's formula](https://en.wikipedia.org/wiki/Approximations_of_%CF%80#Efficient_methods) (*hard*):
$$\\frac{1}{\pi} \\simeq \\frac{2\\sqrt{2}}{9801} \\sum_{n=0}^N \\frac{(4n)!(1103+26390n)}{(n!)^4 396^{4n}}. $$
*Remark:* This method and the next one compute approximation of $\\frac{1}{\pi}$, not $\pi$.
> By the way, did you know that Ramanujan was a brilliant *Indian* mathematician?
#### [Chudnovsky brothers' formula](https://en.wikipedia.org/wiki/Chudnovsky_algorithm) (*hard*):
$$\\frac{1}{\pi} \\simeq 12 \\sum_{n=0}^N \\frac {(-1)^{n}(6n)!(545140134n+13591409)}{(3n)!(n!)^{3}640320^{{3n+3/2}}}. $$
<!-- In 2015, the best method is still the Chudnovsky brothers's formula. -->
> Be careful when you use these formulae, *check carefully* the constants you wrote (545140134 will work well, 545140135 might not work as well!).
---
### Other methods
#### Trigonometric methods (*hard*)
Some methods are based on the $\\mathrm{arccot}$ or $\\arctan$ functions, and use the appropriate Taylor series to approximate these functions.
The best example is [Machin's formula](http://en.literateprograms.org/Pi_with_Machin%27s_formula_%28Python%29):
$$\pi = 16\\mathrm{arccot}(5) - 4 \\mathrm{arccot}(239).$$
And we use the Taylor series:
$$\\mathrm{arccot}(x) = \\frac{1}{x} - \\frac{1}{3x^3} + \\frac{1}{5x^5} - \\frac{1}{7x^7} + \dots = \\sum_{n=0}^{+\infty} \\frac{(-1)^n}{(2n+1) x^{2n+1}} .$$
This method is also explained here with some details.
In order to obtain $n$ digits, we will use *fixed-point* arithmetic to compute $\pi \\times 10^n$ as a Python `long` integer.
##### [High-precision arccot computation](http://en.literateprograms.org/Pi_with_Machin%27s_formula_%28Python%29#High-precision_arccot_computation)
To calculate $\\mathrm{arccot}$ of an argument $x$, we start by dividing the number $1$ (represented by $10^n$, which we provide as the argument `unity`) by $x$ to obtain the first term.
We then repeatedly divide by $x^2$ and a counter value that runs over $3$, $5$, $7$ etc, to obtain each next term.
The summation is stopped at the first zero `term`, which in this *fixed-point* representation corresponds to a real value less than $10^{-n}$:
```python
def arccot(x, unity):
xpower = unity / x
sum = xpower
n = 3
sign = -1
while True:
xpower = xpower / (x*x)
term = xpower / n
if term == 0:
break # we are done
sum += sign * term
sign = -sign
n += 2
return sum
```
##### Applying Machin's formula
Finally, the main function uses Machin's formula to compute $\pi$ using the necessary level of precision, by using this high precision $\\mathrm{arccot}$ function:
```python
def machin(digits):
unity = 10**(digits + 10)
pi = 4 * (4*arccot(5, unity) - arccot(239, unity))
return pi / unity
```
To avoid rounding errors in the result, we use 10 guard digits internally during the calculation.
We may now reproduce the historical result obtained by [Machin in 1706](https://en.wikipedia.org/wiki/John_Machin):
```
>>> machin(100)
3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679
```
The program can be used to compute tens of thousands of digits in just a few seconds on a modern computer.
Many [Machin-like formulas](https://en.wikipedia.org/wiki/Machin-like_formula) also exist, like:
$$\pi = 4\\arctan(\\frac{1}{2}) + 4 \\arctan(\\frac{1}{3})$$
#### (*hard*) [Unbounded Spigot Algorithm](http://www.cs.ox.ac.uk/people/jeremy.gibbons/publications/spigot.pdf)
This algorithm is quite efficient, but not easy to explain. Go check on-line if you want.
#### (*hard*) [Borwein's algorithm](https://en.wikipedia.org/wiki/Borwein%27s_algorithm#Nonic_convergence)
It has several versions, one with a cubic convergence (each new step multiplies by $3$ the number of digits), one with a nonic convergence (of order $9$) etc.
They are not so easy to explain either.
Please check on-line, here [en.WikiPedia.org/wiki/Borwein's_algorithm](https://en.wikipedia.org/wiki/Borwein%27s_algorithm).
The cubic method is similar to the Gauss-Legendre algorithm:
1. Start with $a\_0 = 1/3$, and $s\_0 = \\frac{\\sqrt{3}-1}{2}$,
2. And then iterate, as much as you want, by defining $r\_{k+1} = \\frac{3}{1+2(1-s\_k^3)^{1/3}}$, and updating $s\_{k+1} = \\frac{r\_{k+1}-1}{2}$ and $a\_{k+1} = r\_{k+1}^2 a\_k - 3^k (r\_{k+1}^2 - 1)$.
Then the numbers $a\_k$ will converge to $\\frac{1}{\pi}$.
---
### The `decimal.Decimal` trick to improve precision
If you implement these methods by simply using `float` numbers for all the variables and the partial sum, the final precision of your approximation of $\pi$ will be extremely limited (it will not be better than importing `math.pi`!).
So you should try to use **decimal** numbers instead, by importing the `decimal` module in Python.
More details on that library here [docs.Python.org/2/library/decimal.html](https://docs.python.org/2/library/decimal.html).
The basic thing you will need to use is the `decimal.Decimal` class, available if you import the module with `import decimal` (or `Decimal` if you import the `decimal` module with `from decimal import *`):
```python
from decimal import *
# Example of use of Decimal
mypi = Decimal(22) / Decimal(7)
```
By default, the precision of such decimal numbers is $30$, but you will obviously need to increase the precision.
This can be done simply by increasing `getcontext().prec` (or `decimal.getcontext().prec`):
```python
# If needed, increase the precision up to N digits after the comma:
getcontext().prec = 100000 # precision is now N = 100000
```
---
## Examples and references
> Link are given in the soft copy, available online and on Moodle.
- [en.WikiPedia.org/wiki/Pi#Modern_quest_for_more_digits](https://en.wikipedia.org/wiki/Pi#Modern_quest_for_more_digits),
- [www.JoyOfPi.com/pi.html](http://www.joyofpi.com/pi.html) and [www.JoyOfPi.com/pilinks.html](http://www.joyofpi.com/pilinks.html),
- [www.EveAndersson.com/pi/digits/](http://www.eveandersson.com/pi/digits/) has great interactive tools,
- more crazy stuff [MathWorld.Wolfram.com/PiDigits.html](http://mathworld.wolfram.com/PiDigits.html), or [MathWorld.Wolfram.com/Pi.html](http://mathworld.wolfram.com/Pi.html),
- [one idea with Fibonacci numbers](http://www.maths.surrey.ac.uk/hosted-sites/R.Knott/Fibonacci/fibpi.html#section2),
- [and this incredibly long list of digits](http://piworld.calico.jp/estart.html) at [PiWorld.calico.jp/estart.html](http://piworld.calico.jp/estart.html).
<!-- https://bitbucket.org/lbesson/python-demos/src/master/ComputePie.py -->
<!-- http://thelivingpearl.com/2013/05/28/computing-pi-with-python/ -->
### $100$ first digits of $\pi$
$\pi \\simeq 3.1415926535 ~ 8979323846 ~ 2643383279 ~ 5028841971 \\\\ 6939937510 ~ 5820974944 ~ 5923078164 ~ 9862803482 ~ 53421170679$ when computed to the first $100$ digits.
Can you compute up to $1000$ digits? Up to $10000$ digits? Up to $100000$ digits? **Up to 1 million digits?**
> I failed going further than $100000$ (in less than $2$ minutes).
> [My best method produced 100000 correct digits in 22 seconds](Computing_Pi.py).
---
## Problem 2 : plotting a **pie** chart
- **What do do ?** Plot a pie chart representing any interesting data you can think of (an example is given below, it can be a starting point),
- **How ?** You should start by [using this PyLab/MatPlotLib tutorial](http://www.labri.fr/perso/nrougier/teaching/matplotlib/) (link below),
- **Which data to use ?** Any randomly generated data, or the list of students (a dictionary of students) (ask me how to download it),
- **What to plot ?** Gender ($30%$ girl, $70%$ boy), sections etc : as you want ! Just try to plot something you find interesting.
- **How we would evaluate you ?** Some points if the chart looks fine, extra points for every extra small things (a title, a legend, name of the axis etc).
### Default template program
```python
N = 20
from random import uniform
grades = [ round(uniform(0, 100), 2) for i in xrange(N) ]
from pylab import *
fig1 = figure(1)
pie(grades) # this will look useless, lets do better:
fig2 = figure(2)
data = [ sum( [ 1 for g in grades if percentage <= g < percentage+10 ] )
for percentage in range(0,101,10) ]
pie( [ count for count in data if count > 0 ] )
```
### Examples and references
- [www.LaBri.fr/perso/nrougier/teaching/matplotlib/](http://www.labri.fr/perso/nrougier/teaching/matplotlib/)
- See this example: [www.LaBri.fr/perso/nrougier/teaching/matplotlib/scripts/pie_ex.py](http://www.labri.fr/perso/nrougier/teaching/matplotlib/scripts/pie_ex.py)
- [SciPy-lectures.GitHub.io/intro/matplotlib/matplotlib.html#pie-charts](http://scipy-lectures.github.io/intro/matplotlib/matplotlib.html#pie-charts)
<!--  -->
<!-- https://bitbucket.org/lbesson/python-demos/src/master/Pie_plotting_CS101_grades__First_Mid_Term_2015.py -->
---
## Task 3 : one extra task, but no programming
The task is *quite simply*: learn as many digits of $\pi$ as (you think) you can.
20 is already an impressive number!
One of the simplest ways to memorize Pi is to memorize sentences in which each word's length represents a digit of $\pi$.
A first example is *« May I have a large container of coffee? »*, giving 3.1415926 (7 digits).
A second example is :
> Pie
> I wish I could recollect pi.
> "Eureka!," cried the great inventor.
> Christmas pudding, Christmas pie,
> is the problem's very center!
[Many such poems can be found](https://duckduckgo.com/?q=poem+to+learn+digits+of+pi) online.
### Examples and references
- [www.WikiHow.com/Memorize-Pi](http://www.wikihow.com/Memorize-Pi),
- [www.Ludism.org/mentat/PiMemorisation](http://www.ludism.org/mentat/PiMemorisation),
- [www.SailorPi.com/poem.html](http://www.sailorpi.com/poem.html).
---
# Prices for the best teams
> Based on the criteria given above, and maybe some black magic, we will evaluate each team, by assessing the quality and efficiency of your programs (and more).
1. the best team would receive **FIXME**,
2. the second team would get **FIXME**,
3. and the best learner of digits of $\pi$ would receive **FIXME**.
TODO: get confirmation from *Praveen* and *Hari*?
---
# Last remarks
## Organizing team
Please contact the organizing team directly at `CS101@crans.org` if needed.
## Credit and license
This document has be written by [Lilian Besson](https://perso.crans.org/besson/), in March 2015.
It [is publicly published](https://perso.crans.org/besson/cs101/hackathon/14_03_2015), under the terms of the [GNU Public License v3](https://perso.crans.org/besson/LICENSE.html).
## Disclaimer
Finally, all the quoted resources (websites, books, videos, slides, programs etc) are **the properties of their respective authors**, and neither me (Lilian Besson) nor Mahindra École Centrale are affiliated to any of them.
]]>