TD H2 : Matplotlib, calculs d'aire sous une courbe par la méthode de Monte-Carlo (correction)¶
!!! info Info Pour calculer l'aire sous une courbe, une technique statistique consiste à encadrer cette courbe dans un rectangle et tirer aléatoirement des points dans ce rectangle. Le rapport du nombre de points tombés sous la courbe par rapport aux nombres de points tirés au total et le même que le rapport de l'aire sous la courbe par rapport à l'aire du rectangle...
Par exemple, sur l'image ci-dessous, il y a environ 35% des points tirés dans le rectangle qui sont situés en-dessous de la courbe, or le rectangle a pour aire 1, donc l'aire sous la courbe est approximativment de $35\%\times 1=0,35$. !!!
On pourrrait tout faire d'un coup, mais comme souvent en informatique, on préfaire découper le travail en tâches élémentaires. On séparera ici deux parties. Dans la première partie on s'occupera du côté calculs, puis dans un deuxième temps on gèrera les représentatiosn graphiques.
Parie 1 : les calculs¶
Q1. On commence par définir la zone de travail. Les coordonnées de la zone de travail seront stockées dans un dictionnaire. Compléter le script suivant, sachant que l'on veut estimer l'aire sous la courbe de la fonction carrée entre 0 et 1.
# On définit le cadre à l'aide d'un dictionnaire
zone = {
'xmin' : 0,
'ymin' : 0,
'xmax' : 1,
'ymax' : 1
}
# fonction carrée
def f1(x) :
return x*x
Q2. Pour faire des tirages aléatoires de points, on va utiliser la fonction unniform du module random. Utilsier l'aide intégrée de python pour savoir comment elle fonctionne.
import random
help(random.uniform)
Help on method uniform in module random:
uniform(a, b) method of random.Random instance
Get a random number in the range [a, b) or [a, b] depending on rounding.
Q3. Pour la partie calculs. Compléter le script ci-dessous en suivant les étapes suivantes:
- Écrire une fonction python
points_aleatoires(n,zone)qui reçoit les coordonnées d'une zone (via le dictionnaire du début du TD) et un entier n et qui génrère aléatoirement n points uniforméments répartis dans la zone. La fonction va renvoyer deux listes. La première est la liste des abscisses, la deuxième la liste des ordonnées. - Écrire une fonction python
monte_carlo(f,abscisses,ordonnees,zone)qui renvoie une estimation de l'aire sous la courbe de la fonction f, calculée par la méthode de monte-carlo, décrite en introduction.
# partie calculs
from random import uniform
def points_aleatoires(n,zone):
""" Cette fonction reçoit un entier n et un dictionnaire délimitant
une zone rectangulaire. Elle tire n points aléatoires dans cette zone.
Elle renvoie deux liste : celle des abscisses et celle des ordonnées
"""
x = []
y = []
for i in range(n):
x.append(uniform(zone['xmin'],zone['xmax']))
y.append(uniform(zone['ymin'],zone['ymax']))
return x,y
def monte_carlo(f,abscisses,ordonnees,zone):
""" Cette fonction reçoit une fonction f, deux listes donnant les coordonnées
de points et une zone. Elle renvoie l'aire sous la courbe de f, estimée par
la méthode de Monte-Carlo"""
n_inf = 0
n_tot = len(abscisses)
for i in range(n_tot):
if ordonnees[i] < f(abscisses[i]) :
n_inf += 1
aire_zone = (zone['xmax']-zone['xmin'])*(zone['ymax']-zone['ymin'])
proportion_inf = n_inf/n_tot
return proportion_inf*aire_zone
X,Y = points_aleatoires(100,zone)
aire = monte_carlo(f1,X,Y,zone)
print(f"Estiamtion de l'aire : {aire}")
Estiamtion de l'aire : 0.35
Partie 2 : représentations graphiques¶
Q4. On s'occupe maintenant de la partie dessin. Compléter le script ci-dessous en suivant les étapes ci-dessous :
- Écrire une fonction python
dessine_fonction(f,zone)qui trace la représentation graphique de la fonction f pour x entre xmin et xmax à l'aide du module matplotlib. - Écrire une fonction python
dessine_zone(zone)qui trace un rectangle dont les coins on pour coordonnées (xmin,ymin) et (xmax,ymax) représentant la zone de tirage des points (on utilisera également la fonction plot du module matplotlib). - Écrire une fonction
dessine_monte_carlo(X,Y,f)qui dessine les points donnés par les listes X et Y (en rouge ceux en-dessous de la courbe de la fonction f et en bleu au dessus).
# partie dessin
import matplotlib.pyplot as plt
def dessine_fonction(f,zone):
"""Cette fonction reçoit une fonction f et une zone (sous forme d'un dictionnaire)
et elle trace la courbe représentative de cette fonction"""
pas = 0.01
abscisses = []
ordonnees = []
x = zone['xmin']
while x <= zone['xmax']:
abscisses.append(x)
ordonnees.append(f(x))
x += pas
plt.plot(abscisses,ordonnees)
def dessine_zone(zone):
""" Cette fonction reçoit une zone (sous forme d'un dictionnaire) et
trace le rectangle qui la représente"""
x = [zone['xmin'], zone['xmin'], zone['xmax'], zone['xmax'], zone['xmin']]
y = [zone['ymin'], zone['ymax'], zone['ymax'], zone['ymin'], zone['ymin']]
plt.plot(x,y,'-')
def dessine_monte_carlo(f,X,Y):
"""Cette fonction reçoit deux listes (abscisses,ordonnees) représentant des points et
une fonction f. Elle affiche les points situés sous la courbe en rouge et ceux
au-dessus de la courbe en bleu"""
x_sub = []
y_sub = []
x_sup = []
y_sup = []
for i in range(len(X)) :
if Y[i] < f(X[i]):
x_sub.append(X[i])
y_sub.append(Y[i])
else :
x_sup.append(X[i])
y_sup.append(Y[i])
plt.plot(x_sub,y_sub,'r*')
plt.plot(x_sup,y_sup,'b+')
#On teste :
plt.figure(1)
plt.clf()
dessine_zone(zone)
dessine_fonction(f1,zone)
X,Y = points_aleatoires(100,zone)
dessine_monte_carlo(f1,X,Y)
plt.show()
Q5. Savez-vous déterminer l'aire exacte sous cette courbe? Combien vaut-elle ?
Réponse : c'est $\displaystyle\int_0^1 x^2 dx=\left[\dfrac{x^3}{3}\right]_0^1=\dfrac{1}{3}$
Partie 3 : avec une autre fonction¶
Pour voir si on a compris, on considère maintenant la fonction $f(x)=\sqrt{4-x^2}$ pour $x\in[0,2]$.
Q6. Compléter le script ci-dessous pour qu'il trace comme dans la partie précédente, cette courbe et calcule l'aire située sous la courbe.
import math
def f2(x):
return math.sqrt(4-x*x)
zone = {
'xmin' : 0,
'ymin' : 0,
'xmax' : 2,
'ymax' : 2
}
plt.figure(2)
plt.clf()
X,Y = points_aleatoires(100,zone)
aire = monte_carlo(f2,X,Y,zone)
print(f"Estimation de l'aire : {aire}")
dessine_zone(zone)
dessine_fonction(f2,zone)
dessine_monte_carlo(f2,X,Y)
plt.axis('equal') # pour rendre le repère orthonormé
plt.show()
Estimation de l'aire : 2.92
Q7. Pour contrôler la réponse estimée, on va faire le calcul exact de l'aire sous la courbe. Montrer que cette courbe est un quart de cercle de rayon 2 et donner la valeur exacte de l'aire estimée à la question précédente.
Réponse : La distance du point (0,0) au point (x,f(x)) vaut toujours 2 (car $x^2+\sqrt(4-x^2)^2=4$), donc cette courbe est bien une partie du cercle de centre O et de rayon 2. Son aire vaut un quart de l'aire du disque, d'où aire = $\pi\times 4/4=\pi$.
Partie 4 : pour les plus rapides¶
Q8. Pour ceux qui sont plus à l'aise, écrire une fonction python determine_zone(f,xmin,xmax) qui renvoie le dictionnaire zone en calculant les ymin et ymax.
def determine_zone(f,xmin,xmax):
""" calcule une estimation du minimum et du maximum
de la fonction f, sur l'intervalle [xmin,xmax].
Elle teste les valeurs de f(x) pour x variant entre
xmin et xmax avec un pas de 0.01.
Elle renvoie la zone sous forme d'un dictionnaire"""
dx = 0.01 # le pas d'incrémentation de x
x = xmin
ymin = f(x)
ymax = f(x)
while x<=xmax :
y = f(x)
if y<ymin :
ymin = y
if y>ymax :
ymax = y
x += dx
zone = {
'xmin' : xmin,
'ymin' : ymin,
'xmax' : xmax,
'ymax' : ymax
}
return zone
# On teste avec f1 par exemple
print(determine_zone(f1,0,2))
{'xmin': 0, 'ymin': 0, 'xmax': 2, 'ymax': 3.960100000000006}
# Encore un exemple complet
def f3(x):
return math.log(1+x)
zone = determine_zone(f3,0,5)
X,Y = points_aleatoires(100,zone)
aire = monte_carlo(f3,X,Y,zone)
print(f"Estimation de l'aire : {aire}")
plt.figure(2)
plt.clf()
dessine_zone(zone)
dessine_fonction(f3,zone)
dessine_monte_carlo(f3,X,Y)
plt.axis('equal') # pour rendre le repère orthonormé
plt.show()
Estimation de l'aire : 6.091982195375352