Méthode de Romberg

On souhaite calculer l'intégrale ${\displaystyle I=\int _{a}^{b}f(x)\,dx}$ d'une fonction f supposée régulière sur [a,b] en subdivisant [a,b] en n sous-intervalles identiques, du type ${\displaystyle [a+kh,a+(k+1)h]}$ pour ${\displaystyle k=0,...,n-1}$ et ${\displaystyle h=(b-a)/n}$. La méthode des trapèzes, notée T(n), est définie à l'aide de cette grille régulière :

${\displaystyle T(n)={\frac {h}{2}}\left[\sum _{k=0}^{n}\omega _{k}f(a+kh)\right]}$

où les pondérations ${\displaystyle \omega _{k}}$ sont égales à 1 pour les points extrêmes, et à 2 pour les autres. La méthode est dite d’ordre 1, car donne un résultat exact pour un polynôme de degré 1. D'après la formule d'Euler-Maclaurin, l'erreur commise, notée ${\displaystyle E=I-T(n)\,}$, vérifie

${\displaystyle E(h)=c_{1}h^{2}+c_{2}h^{4}+c_{3}h^{6}\cdots .}$

où les coefficients c_k ne dépendent pas du pas h choisi, et seules les puissances paires de h apparaissent. Cette propriété (qui est aussi partagée avec la méthode du point médian) sera avantageusement exploitée, pour supprimer, un à un les termes de l'erreur, à l'aide de l'extrapolation de Richardson.

Des manipulations algébriques fournissent une approximation de I en ${\displaystyle O(h^{4})}$. Pour ce faire, il s'agit d'éliminer le terme en h² dans le développement de l'erreur. On y parvient en considérant la quantité [4 T(n) – T(n/2)]/3, qui correspond précisément à la méthode de Simpson. Le même procédé peut être reconduit, permettant ainsi d'annuler successivement les termes en ${\displaystyle h^{4},h^{8},\ldots }$, conduisant ainsi à des approximations de I qui sont de plus en plus précises.

Formalisations de la technique précédente de réduction de l'erreur :

• Initialisations :

${\displaystyle R(0,0)={\frac {1}{2}}(b-a)(f(a)+f(b)).}$

${\displaystyle R(n,0)=T(2^{n}),\,}$

soit le résultat de la méthode des trapèzes basée sur

${\displaystyle h_{n}={\frac {b-a}{2^{n}}}.}$

• Récurrence :

${\displaystyle R(n,m)={\frac {4^{m}R(n,m-1)-R(n-1,m-1)}{4^{m}-1}}}$

On montre par récurrence que l'approximation à l’étape n, soit R(n,n), fournit une approximation de I qui est en ${\displaystyle O(h_{n}^{2n+2})}$, ceci à condition que la fonction soit de classe ${\displaystyle {\mathcal {C}}^{2n+2}}$ (ou ${\displaystyle {\mathcal {C}}_{I}^{2n+2}}$). Le calcul de R(n,n) est exact pour les polynômes de degré 2 n + 1 : elle est d’ordre 2n + 1.

L'algorithme peut être représenté sous la forme d'un tableau. La première diagonale R(n,n) fournit les approximations successives ; la première colonne R(n, 0) correspond (par définition) à la méthode du trapèze, la seconde R(1,n) reflète la méthode de Simpson, la troisième celle de Boole qui est la formule de Newton-Cotes à 5 points appliquée à une décomposition en ${\displaystyle 2^{n-2}+1}$ sous-intervalles. Par contre, les colonnes suivantes correspondent à des formules de quadrature différentes de celles de Newton-Cotes d’ordre supérieur, évitant ainsi les problèmes d’instabilité.

On considère le calcul de l'intégrale

${\displaystyle I=\int _{-1}^{2}{\frac {2}{1+4t^{2}}}\,\mathrm {d} t.}$

On trouve ${\displaystyle I=\arctan(4)+\arctan(2)\simeq 2,432966381462123}$. La matrice de l'algorithme est

d'où l'on tire les approximations consécutives

• Méthode des trapèzes
• Méthode de Simpson
• Extrapolation de Richardson

• (en) Recursive integration formulas from Romberg integration

.