# Monte Carlo simulation

A question that I read today: “How to estimate \(\pi\)?”

You’d find the thorough overview in the above link (and maybe a lot more blogs available online). I also gave my brain a little exercise today:

#### We can leverage MC simulation!

This is a very straight-forward method that everyone might come up with. As shown below, we know that the relative ratio of the areas of a circle with radius (\(r\)) and a square of width (\(2r\)) can be expressed like:

\[\text{circle}:\text{square} = \pi R^2 : (2R)^2 = \pi : 4\]Therefore, on the Cartesian coordinate, if we can sample N samples of (x, y) pairs uniformly from the square and to keep track of the common samples that also fell into the circle (threshold: \(x^2 + y^2 < R^2\)), we can estimate \(\pi\) as:

\(\hat{\pi} = 4 \times \frac{\text{# of common samples also appears in the circle}}{\text{# of samples in the square}}\). Here is part of my sample codes:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23

R = 1
N = 1000
def estimate_pi(R, N):
"""
Estimate pi from simulation from circle
N: Number of simulation/samples
R: radius of the circle
"""
x = np.random.uniform(-R, R, N)
y = np.random.uniform(-R, R, N)
in_circle = x**2 + y**2
circle_count = 0
for i in range(N):
if in_circle[i] < R**2:
circle_count += 1
pi_hat = 4 * circle_count / N
return pi_hat

Example notebook with above example can be found here.