Estimer une aire par tirage au sort

Estimer une aire par tirage au sort.

La méthode de Monte-Carlo consiste à tirer des points au hasard dans un rectangle contenant le domaine, puis à compter la proportion de ceux qui tombent dedans. C'est ainsi que les logiciels d'analyse d'images mesurent la surface d'une feuille ou d'une tache sur une photographie.

Expliquer pourquoi la proportion de points « sous la courbe » approche le rapport de l'aire cherchée à l'aire du rectangle.
Tester la méthode sur $\int_{0}^{1} 1 - x^{2} d x$ , dont on sait qu'elle vaut $\frac{π}{4}$ (aire d'un quart de disque de rayon $1$ ).
L'appliquer à $\int_{0}^{1} e^{- x^{2}} d x$ , dont aucune primitive ne s'exprime avec les fonctions usuelles.

Le principe. Un point tiré au hasard uniformément dans le rectangle $[a; b] \times [0; h_{m a x}]$ a autant de chances de tomber dans n'importe quelle partie de ce rectangle : la proportion de tirages qui tombent dans une région donnée approche donc le rapport de l'aire de cette région à l'aire totale $(b - a) h_{m a x}$ . En multipliant la proportion observée par $(b - a) h_{m a x}$ , on obtient une estimation de l'aire cherchée.
et 3. Programme.


from math import exp, sqrt
from random import random, seed

def monte_carlo(f, a, b, hmax, n):
    dedans = 0
    for _ in range(n):
        x = a + (b - a) * random()
        y = hmax * random()
        if y <= f(x):
            dedans = dedans + 1
    return (b - a) * hmax * dedans / n

demi_cercle = lambda x: sqrt(1 - x * x)
seed(7)
print("quart de disque, n = 10^3 :", round(monte_carlo(demi_cercle, 0, 1, 1, 1000), 5))
seed(7)
print("quart de disque, n = 10^5 :", round(monte_carlo(demi_cercle, 0, 1, 1, 100000), 5))
print("pi / 4                    :", round(3.141592653589793 / 4, 5))

cloche = lambda x: exp(-x * x)
seed(7)
print("cloche, n = 10^4          :", round(monte_carlo(cloche, 0, 1, 1, 10000), 5))
seed(7)
print("cloche, n = 10^6          :", round(monte_carlo(cloche, 0, 1, 1, 1000000), 5))

```text

quart de disque, n = 10^3 : 0.788 quart de disque, n = 10^5 : 0.78556 pi / 4 : 0.7854 cloche, n = 10^4 : 0.7499 cloche, n = 10^6 : 0.7468 ``` Le test réussit : $0, 78556$ contre $0, 785398$ . On accorde donc du crédit à l'estimation $\int_{0}^{1} e^{- x^{2}} d x \approx 0, 7468$ — que les trapèzes avec $n = 1000$ confirment : $0, 746824$ .

Deux enseignements. D'abord, on valide toujours une méthode numérique sur un cas dont on connaît la réponse avant de la lâcher sur un cas inconnu. Ensuite, Monte-Carlo est lent : multiplier $n$ par $100$ ne gagne qu'un facteur $10$ sur la précision. Sa vertu n'est pas la vitesse, c'est de fonctionner là où le découpage régulier est impossible.

Énoncé

Corrigé