Application numérique : la méthode des trapèzes

Leçon optionnelle — Illustration concrète de l'intégration numérique

Contexte

Dans la leçon « Pourquoi certaines intégrales n'ont pas de solution analytique », nous avons vu que l'intégrale

n'admet aucune primitive en forme close. C'est une intégrale elliptique, prouvablement non élémentaire par le théorème de Liouville.

Pourtant, cette intégrale existe bel et bien — c'est l'aire sous une courbe parfaitement définie. Comment la calculer? C'est ici qu'intervient l'analyse numérique.

La méthode des trapèzes

Idée géométrique

La méthode des trapèzes approche l'aire sous la courbe par une somme d'aires de trapèzes. Sur chaque sous-intervalle, on relie les deux points de la courbe par un segment de droite.

Formule

On découpe l'intervalle $[a, b]$ en $n$ sous-intervalles de largeur $h = \frac{b-a}{n}$. Les points de subdivision sont :

L'approximation de l'intégrale est alors :

Ou, sous forme compacte :

Implémentation

import numpy as np

def f(x):
  """Notre intégrale elliptique."""
  return np.sqrt(1 + np.cos(x)**2)

def trapezes(a, b, n):
  """
  Méthode des trapèzes avec n subdivisions.

  Paramètres :
      a, b : bornes de l'intervalle
      n    : nombre de subdivisions

  Retourne :
      Approximation de l'intégrale
  """
  h = (b - a) / n
  x = np.linspace(a, b, n + 1)
  y = f(x)

  # Formule des trapèzes
  return h * (0.5 * y[0] + np.sum(y[1:-1]) + 0.5 * y[-1])

# Calcul pour différentes valeurs de n
print("Méthode des trapèzes pour J = ∫₀^π √(1 + cos²(x)) dx")
print("-" * 50)

for n in [2, 4, 6, 8, 10, 12, 16]:
  I = trapezes(0, np.pi, n)
  print(f"n = {n:3d} : I ≈ {I:.10f}")

Résultats

Ajouter deux subdivisions divise l'erreur par un facteur qui reste compris entre 40 et 65 tout au long du tableau. Au-delà de $n = 16$, l'erreur atteint la précision machine et cesse de décroître : ce n'est plus la méthode qu'on observe, mais la limite du calcul en double précision.

Analyse de la convergence

Ordre de convergence

On observe une convergence bien plus rapide que $O(h^2)$ : dès $n = 10$ l'erreur vaut environ $1{,}5 \times 10^{-9}$, et à partir de $n = 20$ on atteint la précision machine.

C'est un cas particulièrement favorable : $f(x) = \sqrt{1 + \cos^2(x)}$ est analytique et $\pi$-périodique, et $[0, \pi]$ couvre exactement une période. Toutes les dérivées coïncident aux deux bornes, les termes de bord de la formule d'Euler-Maclaurin s'annulent, et la méthode des trapèzes converge alors plus vite que n'importe quelle puissance de $h$.

Borne théorique de l'erreur

Pour une fonction $f$ deux fois dérivable, l'erreur de la méthode des trapèzes est bornée par :

Dans le cas général, cette borne garantit une convergence en $O(h^2)$ : doubler $n$ divise l'erreur au moins par 4. Ce n'est qu'une majoration, et l'exemple ci-dessus fait bien mieux.

Le cas général, sur la même fonction

Pour voir le $O(h^2)$ annoncé, il suffit de garder la même fonction et de changer l'intervalle. Sur $[0, 1]$, qui ne correspond à aucune période, les termes de bord ne s'annulent plus :

La valeur exacte vaut ici $1{,}3114424982$, et le rapport se tient sur 4.00 à chaque doublement : c'est exactement le comportement d'ordre 2. Deux intervalles, une seule fonction, deux régimes de convergence — ce n'est donc pas l'intégrande qui décide, mais la façon dont ses dérivées se comportent aux bornes.

La formule d'Euler-Maclaurin fait apparaître des termes de bord en $f'(b) - f'(a)$, puis en $f'''(b) - f'''(a)$, et ainsi de suite. Sur $[0, \pi]$, les dérivées impaires prennent la même valeur aux deux bornes : toutes ces différences s'annulent, et il ne reste aucun terme, à aucun ordre fini. Sur $[0, 1]$, $f'(0) = 0$ mais $f'(1) \approx -0{,}40$ : le premier terme survit et impose son $O(h^2)$.

Le cas général reste donc celui-ci : la méthode des trapèzes est d'ordre 2. La convergence spectrale observée plus haut tient à la périodicité de l'exemple, ce n'est pas la règle.

Ce que montre cet exemple

Lien avec les fonctions spéciales

Cette valeur correspond exactement à :

où $E(k)$ est l'intégrale elliptique complète de seconde espèce :

Les tables de fonctions spéciales (comme le Handbook of Mathematical Functions d'Abramowitz et Stegun) donnent $E(1/\sqrt{2}) \approx 1.3506...$.

Pour aller plus loin

La méthode des trapèzes est la plus simple des méthodes d'intégration numérique, mais pas la plus efficace. Dans les chapitres suivants, nous verrons :

• La méthode de Simpson (convergence en $O(h^4)$)
• Les méthodes de Gauss (précision optimale pour un nombre donné d'évaluations)
• Les méthodes adaptatives (qui ajustent automatiquement le pas)