Monte Carlo pi

Charpter 2, page 38


In [10]:
import matplotlib.pyplot as plt
import numpy as np
from scipy import stats
import seaborn as sns
import tqdm

%matplotlib inline
sns.set(font_scale=1.5)

In [2]:
N = 10000
x, y = np.random.uniform(-1, 1, size=(2, N))
inside = (x**2 + y**2) <= 1
pi = inside.sum()*4/N
error = abs((pi - np.pi)/pi)* 100
outside = np.invert(inside)

In [3]:
plt.plot(x[inside], y[inside], 'b.')
plt.plot(x[outside], y[outside], 'r.')
plt.plot(0, 0, label='$\hat \pi$ = {:4.3f}\nerror = {:4.3f}%'.
format(pi, error), alpha=0)
plt.axis('square')
plt.legend(frameon=True, framealpha=0.9, fontsize=16);



In [11]:
N = np.logspace(1, 7, 50)
 
ans = []
for nn in tqdm.tqdm_notebook(N):
    x, y = np.random.uniform(-1, 1, size=(2, int(nn)))
    inside = (x**2 + y**2) <= 1
    pi = inside.sum()*4/int(nn)
    error = abs((pi - np.pi)/pi)* 100
    outside = np.invert(inside)
    ans.append((nn, pi, error))
ans = np.asarray(ans)




In [21]:
plt.errorbar(ans[:,0], ans[:,1], yerr=ans[:,2])
plt.xscale('log')
plt.ylim((2.5, 3.5))
plt.axhline(np.pi, lw=1, c='r', ls='--')


Out[21]:
<matplotlib.lines.Line2D at 0x7fe4c05fde90>

In [ ]: