Exercice: Calcul De Pi Par La Méthode De Monte-Carlo - Normale Sup

Home » Approximation De Pi Par La Méthode De Monte Carlo » Exercice: Calcul De Pi Par La Méthode De Monte-Carlo - Normale Sup

Maybe your like

Lecture notes > Python> Exercice: Calcul de pi par la méthode de Monte-Carlo [Download notebook notebook without solution] Exercice: Calcul de pi par la méthode de Monte-Carlo¶

Les méthodes de Monte-Carlo sont une famille d'algorithmes datant des années 40 permettant la résolution de problèmes d'optimisation, d'intégration ou d'échantillonage. Elles reposent sur le tirage de nombres (pseudo)-aléatoires dont l'espérance de la distribution permet de calculer la quantité recherchée.

In [1]: import numpy as np import matplotlib.pylab as plt %matplotlib inline

Le calcul de $\pi$ par la méthode de Monte-Carlo repose sur l'intégration numérique du quart de cercle. L'aire d'un quart de disque unité est:

\begin{equation} A = \frac{\pi}{4} \end{equation}

Soit $(X_n)_{n\geq1}$, une suite de variables aléatoires indépendantes et uniformes sur le pavé $[0,1]^2$. L'aire de ce pavé est $S=1$

La probabilité que $X_i$ tombe dans le disque unité est de $\frac{A}{S}=\frac{\pi}{4}$.

Soit $(Y_n)_{n\geq1}$, une suite de variables aléatoires de Bernouilli. Pour tout $i$, $Y_i$ vaut 1 si $X_i$ est dans le disque unité. $4\mathbb E(Y_0) = \pi$.

Par la loi des grands nombres, la moyenne empirique de $Y_n$ tend vers l'espérance $\mathbb E(Y_0)$:

\begin{equation} \lim_{n \to \infty} \frac{4}{n} \sum_{i=1}^n Y_i = 4 \mathbb E(Y_0) = \pi \end{equation}

In [ ]: # Calculez pi par la méthode de Monte-Carlo In [2]: N = 10000 # On tire un grand nombre de points uniformément sur [0,1]^2 x = np.random.random(size=N) y = np.random.random(size=N) # On créer un array de booléens indiquant si le point se trouve dans le quart de cercle. z = (x**2+y**2)<1 # On estime pi pi_est = 4*z.sum()/N print('Estimation de pi: {}'.format(pi_est)) Estimation de pi: 3.1252 In [6]: # plt.subplots() permet de créer un objet figure des objets axes. fig, ax = plt.subplots(1,1,figsize=(5,5)) # L'indexage booléen fait partie des "indexages avancés" de numpy. # Il est possible de filter un array en l'indexant par un autre array contenant # des booléens. On obtient un array plus petit ne contenant que les elements aux emplacements "True". plt.scatter(x[z], y[z], marker='.') # Logical not permet d'inverser la valeur logique des éléments d'un np.array. nz = np.logical_not(z) plt.scatter(x[nz], y[nz], marker='.') # La méthode set de l'objet `axe` permet de spécifier le titre. ax.set(title=r"Calcul de $\pi$ par la méthode de Monte Carlo.", xlabel='x', ylabel='y') # La méthode savefig de l'objet `figure` permet de sauvegarder l'image dans un fichier. # Le format est déduit de l'extension. fig.savefig('monte_carlo.png') In [7]: fig, ax = plt.subplots(1,1,figsize=(10,5)) # Notez l'utilisation de np.cumsum(z) qui renvoie un array de la même taille # que z mais avec la somme cumulée le long de l'array. conv = 4*np.cumsum(z)/(np.arange(N)+1) plt.plot(conv, label=r'$4X_n/n$') plt.hlines(np.pi, 0, N, label=r'$\pi$') ax.set(title=r"Convergence du calcul de $\pi$ par la méthode de Monte Carlo.", xlabel='N') plt.legend();

Tag » Approximation De Pi Par La Méthode De Monte Carlo