In [2]:
%matplotlib inline

In [4]:
import matplotlib.pyplot as plt
from matplotlib import gridspec
import numpy as np

In [5]:
plt.clf()
circle = plt.Circle((0, 0), radius=0.5, fc='y')
plt.plot(np.random.uniform(-0.5,0.5,1000),np.random.uniform(-0.5,0.5,1000),'b+')
plt.gca().add_patch(circle)
plt.axis('scaled')


Out[5]:
(-0.60000000000000009,
 0.60000000000000009,
 -0.60000000000000009,
 0.60000000000000009)

In [7]:
"""
Using Monte Carlo to find the value of pi.
Intiution: Ratio of area of circle to area of uniq length squre it is inscribed it can be used to approximate pi.
Area of circle = pi/4
Area of square = 1
Ratio = pi/4
Let us generate random points x,y inside the square. The probabilty of the point being inside circle is equal to the above ratio.
Circle centered on origin.
Given: -0.5 <= x,y <= 0.5
A point is inside circle if x**2 + y**2 <= 1/4
"""
import matplotlib.pyplot as plt
from matplotlib import gridspec
import numpy as np

plt.clf()
MAXN = 10000 # Max number of iterations
n = 1 # Current value of iteration
n_arr = [] # Saving values of iterations
pi_arr = [] # Approx values pi for given iteration
n_circle = 0 # Number of points inside circle

fig = plt.figure(figsize=(10,15))
gs = gridspec.GridSpec(2, 1, height_ratios=[2, 1]) 
ax = [plt.subplot(gs[0]),plt.subplot(gs[1])]
# Create a circle just for representation
circle = plt.Circle((0, 0), radius=0.5, fc='y')
ax[0].add_patch(circle)

# Run the simulation
while n <= MAXN:
    p_color = 'b' # If point is outside circle mark it with blue color
    x = np.random.uniform(-0.5,0.5)
    y = np.random.uniform(-0.5,0.5)
    if (x**2 + y**2) <= 0.25:
        n_circle += 1
        p_color = 'r' # If point is outside circle mark it with blue color
    ax[0].plot(x,y,p_color+'+')
    n_arr.append(n)
    pi_arr.append(n_circle*4.0/n) # Value of pi = 4*points in circle/points overall
    n+= 1
    
print n, n_circle, pi_arr[-1], n_arr[-1]    
ax[1].plot(n_arr,pi_arr,"b-",label="monte carlo")
ax[1].plot(n_arr, np.pi*np.ones(len(n_arr)),"r-", label="actual")
ax[1].legend()
ax[1].set_xlabel("Number of iterations")
ax[1].set_ylabel("Value of pi")
plt.title("Monte Carlo approximation of pi")


10001 7848 3.1392 10000
Out[7]:
<matplotlib.text.Text at 0x80fba50>
<matplotlib.figure.Figure at 0x543b2d0>

In [7]: