Random

# Estimating Pi with Monte Carlo method.

## Using the Monte Carlo method to estimate the value of pi.

Monte Carlo is statistical method based on a series of random numbers. Monte Carlos are used in a wide-range of physical systems, finance and other research field.

In this post, we are going to use the Monte Carlo method to estimate the value of pi.

First, consider a circle inscribed in a square (as in the figure). f assume that the radius of the circle is $R$, then the Area of the circle $= Pi * R^2$ and the Area of the squar $= (2 * R)^2 = 4 * R^2$.

If we throw a dart blindly inside of the square, what is the probability (P) the dart will actually land inside the circle?

P = Area of the circle / Area of the square = Pi / 4

So the chances of hitting the circle are Pi / 4.

In other words, $pi = 4 * P$

Let’s try to do this in python:

1 from __future__ import division
2 import numpy
3
4 NBR_POINTS = 1000000
6
7 print (len(filter(lambda x:numpy.sqrt(numpy.random.randint(0, RADIUS)**2+numpy.random.randint(0, RADIUS)**2)<RADIUS,xrange(NBR_POINTS)))/NBR_POINTS)*4

If you have problem with the filter, read this!

To go a little bit further with this we can imagine a map-reduce system to make the estimation faster.

 1 from __future__ import division
2 import collections
3 import itertools
4 import multiprocessing
5 import numpy
6
7 class MapReduce(object):
8     """
9     The map reduce object, should be initialized with:
10         map_fn
11         reduce_fn
12         nbr_workers
13     """
14     def __init__(self, map_fn, reduce_fn, num_workers=None):
15         """
16         initiaize the mapreduce object
17             map_fn : Function to map inputs to intermediate data, takes as
18             input one arg and returns a tuple (key, value)
19             reduce_fn : Function to reduce intermediate data to final result
20             takes as arg keys as produced from the map, and the values associated with it
21         """
22         self.map_fn = map_fn
23         self.reduce_fn = reduce_fn
24         self.pool = multiprocessing.Pool(num_workers)
25
26     def partition(self, mapped_values):
27         """
28         returns the mapped_values organised by their keys. (keys, associated values)
29         """
30         organised_data = collections.defaultdict(list)
31         for key, value in mapped_values:
32             organised_data[key].append(value)
33         return organised_data.items()
34
35     def __call__(self, inputs=None, chunk_size=1):
36         """
37         process the data through the map reduce functions.
38         inputs : iterable
39         chank_size : amount of data to hand to each worker
40         """
41         mapped_data = self.pool.map(self.map_fn, inputs, chunksize = chunk_size)
42         partioned_data = self.partition(itertools.chain(*mapped_data))
43         reduced_data = self.pool.map(self.reduce_fn, partioned_data)
44         return reduced_data
45
46
47 NBR_POINTS = 1000000
49 NBR_WORKERS = 4
50 NBR_PER_WORKER = int(NBR_POINTS/NBR_WORKERS)
51
52
53 def probability_calculation(item):
54     print multiprocessing.current_process().name, 'calculating', item
55     output = []
66 print pi