intégration numérique le match!

Les méthodes d'intégration numériques ça fait parti des sujets  que j'aime bien enseigner.  En effet c'est un sujet qui permet de donner une signification concrète à un concept abstrait (une intégrale) tout en  mêlant de la géométrie élémentaire, de l'analyse  et des bases de programmation.  Comme  pour l'interpolation polynomiale voici un petit comparatif des méthodes les plus communes :

Dans toute la suite on va prendre en exemple  le calcul de l'intégrale de la fonction  $f(x)=\sin(2x/3)$  sur l'intervalle $[0,2\pi]$:
$$ \int_0^{2\pi}f(x) dx= 9/4 $$
dans l'animation le fond passe en cyan quand l'erreur entre la valeur exacte 9/4 et les valeurs approchée d'une méthode passe en deçà de $10^{-5}$.

1. La méthode des rectangles

C'est la méthode est la plus simple, et qui va servir de modèle, elle consiste  a :

• découper l'intervalle d'intégration  en n sous-intervalles à partir des points équidistants $$x_k= a+{b-a\over n}(k-1)$$
•  remplacer l'intégrale sur chaque intervalle $[x_k,x_{k+1}]$  par l'aire du rectangle de hauteur $f(x_k)$ la valeur de $f$ du coté gauche de l'intervalle $${\mathcal A}_k=f(x_k)\times {b-a\over n}$$
• approcher l'intégrale par la somme des aires des rectangles $$\int_a^b f(x) dx \approx \sum_{k=1}^n {\mathcal A}_k$$

Voici un exemple en images :

 En faisant une hypothèse sur la régularité de f on peut majorer la vitesse de convergence vers la valeur exacte de l'intégrale.

Théorème [erreur dans la méthode des rectangles]
Si $f$  est   dérivable sur $[a,b]$ et sa dérivée est bornée sur $[a,b]$  alors
$$\left\vert \int_a^bf(x) dx -\sum_{k=1}^n {\mathcal A}_k\right\vert \leq \sup_{t\in[a,b]}\vert f'(t)\vert{(b-a)^2\over 2n}$$

D'après cette formule l'erreur doit être majorée par
 $$erreur\leq \sup_{t\in[0,2\pi]}\vert (2/3)\cos(2t/3)\vert{(2\pi)^2\over 2n}\approx {13.16\over n}$$

C'est bien le cas dans l'exemple ici puisque l'erreur se comporte asymptotiquement comme $2.7/n$ quand $n\to\infty$ . Le gros défaut de cette méthode est la lenteur avec laquelle on s'approche de la valeur exacte de l'intégrale, avec 10000 rectangles l'erreur est encore de l'ordre de $10^{-3}$.

2. La méthode des Trapèzes

Pour améliorer la méthode précédente on  peut remplacer les rectangles par des trapèzes construits sur les sous-intervalles  $[x_k,x_{k+1}]$ en prenant pour hauteurs $f(x_k)$  et $f(x_{k+1})$. Ceci permet de prendre en compte le sens de variation de la fonction sur chaque  sous-intervalle et devrait réduire l'erreur. Dans la formule d'approximation, l'aire du rectangle est donc remplacée par celle du trapèze :

$${\mathcal A}_k={f(x_k)+f(x_{k+1})\over 2}\times {b-a\over n}$$
en image sur le même exemple que pour la méthode des rectangles :

Avec une hypothèse de régularité sur f on peut majorer la vitesse de convergence vers la valeur exacte de l'intégrale.

Théorème [erreur dans la méthode des trapèzes]
Si $f$  est  deux fois dérivable sur $[a,b]$ et sa dérivée seconde est bornée sur $[a,b]$  alors
$$\left\vert \int_a^bf(x) dx -\sum_{k=1}^n {\mathcal A}_k\right\vert \leq \sup_{t\in[a,b]}\vert f''(t)\vert{(b-a)^3\over 12n^2}$$

D'après cette formule l'erreur doit être majorée par
 $$erreur\leq \sup_{t\in[0,2\pi]}\vert (2/3)^2\sin(2t/3)\vert{(2\pi)^3\over 12n^2}\approx {9.19\over n^2}$$

C'est bien le cas dans l'exemple ici puisque l'erreur se comporte asymptotiquement comme $3.28/n^2$ quand $n\to\infty$ . Concrètement on gagne deux fois plus de décimales exactes  que dans la méthode des rectangles pour la même valeur de n .

3. La méthode du point milieu

Étonnamment il est possible d'obtenir une méthode aussi efficace que celle des Trapèzes en  utilisant des rectangles mais en prenant pour hauteur du rectangle la valeur de f au centre de l'intervalle $[x_k,x_{k+1}]$ ce qui donne pour l'aire du rectangle :

$${\mathcal A}_k=f\left({x_k+x_{k+1}\over 2}\right)\times {b-a\over n}$$
Avec la même hypothèse de régularité sur f on obtient une vitesse de convergence comparable à celle donnée par la méthode des trapèze.

Théorème [erreur dans la méthode du point milieu]
Si $f$  est  deux fois dérivable sur $[a,b]$ et sa dérivée seconde est bornée sur $[a,b]$  alors
$$\left\vert \int_a^bf(x) dx -\sum_{k=1}^n {\mathcal A}_k\right\vert \leq \sup_{t\in[a,b]}\vert f''(t)\vert{(b-a)^3\over 24n^2}$$

ce n'est pas si étonnant car si on fait un développement de Taylor-Lagrange de f autour du point $x_k$, en posant h=(b-a)/n  on a , on obtient :
$$\begin{eqnarray*}
f\left(x_{k+1}\right)&=& f(x_k)+h f'(x_k)+ {h^2\over 2} f''(e_k)\\
{f(x_k)+f(x_{k+1})\over 2}&=& f(x_k)+{h\over 2} f'(x_k)+ {h^2\over 4} f''(e_k) \\
f\left({x_k+x_{k+1}\over 2}\right)&=& f(x_k)+{h\over 2} f'(x_k)+ {h^2\over 8} f''(d_k)\\
\end{eqnarray*}$$

la formule du trapèze et du point milieu sont donc égales à l'ordre 1, et à l'ordre suivant le terme  ${h^2\over 8}$  pour le point milieu est deux fois plus petit que pour les  trapèzes. Dans l'exemple  l'erreur se comporte asymptotiquement comme $1.64/n^2$ quand $n\to\infty$ c'est exactement deux fois moins que ce qu'on a trouvé pour la méthode des trapèzes.

4. La méthode de Simpson

Pour améliorer la méthode on peut chercher une forme géométrique simple qui permette sur chaque intervalle d'approcher le profil de la courbe en tenant compte de la concavité. La solution la plus simple est d'utiliser un morceau de parabole passant par les 3 points $x_k$, $x_{k+1}$ et ${x_k+x_{k+1}\over 2}$, puis de calculer la valeur de l'intégrale  du trinôme correspondant en fonction des valeurs de f  en ces trois points ...  c'est de cette manière que Simpson a obtenu au XVIIIième siècle la formule :

$${\mathcal A}_k=\left(f(x_k)+4f\left({x_k+x_{k+1}\over2}\right)+f(x_{k+1})\right)\times {b-a\over 6n}$$

Théorème [erreur dans la méthode de Simpson]
Si $f$  est  quatre fois dérivable sur $[a,b]$ et sa dérivée quatrième est bornée sur $[a,b]$  alors
$$\left\vert \int_a^bf(x) dx -\sum_{k=1}^n {\mathcal A}_k\right\vert \leq \sup_{t\in[a,b]}\vert f''''(t)\vert{(b-a)^5\over 2880n^4}$$

La méthode est  très rapide comme on peut le voir sur l'animation. La puissance $n^4$  signifie qu'on gagne 4 décimales à chaque fois qu'on multiplie par 10 le nombre d'intervalles. Dans l'exemple  l'erreur se comporte asymptotiquement comme $0.24/n^4$ quand $n\to\infty$ et pour $n=30$ la valeur approchée  est  déjà égale à la valeur exacte pour la précision affichée!

5. Autres méthodes

On pourrait raffiner encore  la méthode mais le gain en précision se ferait au prix d'une complexité  beaucoup plus grandes dans les formules et les calculs à mettre en œuvre. Il existent pourtant d'autres méthodes d'intégration basées sur des approches différentes :

• les méthodes de quadrature , basées sur des approximations polynomiales
• les méthodes probabilistes, basées sur des tirages aléatoires, comme la méthode de Monte-Carlo  (voir aussi mon billet sur la loi des grands nombres )