Rho algorithme

La formule d'interpolation de Thiele (en)[2] est une méthode de calcul de fraction rationnelle d'interpolation s'apparentant à celle de Newton pour les polynômes d'interpolation. La fraction rationnelle de Thiele passant par les points (x_i, S_i) pour i=n, n+1, ..., n+k se présente sous forme d'une fraction continue généralisée qui s'écrit :

${\displaystyle T_{k}^{(n)}(x)=\rho _{0}^{(n)}+{\cfrac {x-x_{n}}{\rho _{1}^{(n)}-\rho _{-1}^{(n)}+{\cfrac {x-x_{n+1}}{\rho _{2}^{(n)}-\rho _{0}^{(n)}+{\cfrac {x-x_{n+2}}{\rho _{3}^{(n)}-\rho _{1}^{(n)}+{\cfrac {x-x_{n+3}}{\rho _{4}^{(n)}-\rho _{2}^{(n)}+\cdots +{\cfrac {...}{\rho _{k-1}^{(n)}-\rho _{k-3}^{(n)}+{\cfrac {x-x_{n+k-1}}{\rho _{k}^{(n)}-\rho _{k-2}^{(n)}}}}}}}}}}}}}}$

soit, pour tous les points à partir de (x₀, S₀) :

${\displaystyle S(x)=\rho _{0}^{(0)}+{\cfrac {x-x_{0}}{\rho _{1}^{(0)}-\rho _{-1}^{(0)}+{\cfrac {x-x_{1}}{\rho _{2}^{(0)}-\rho _{0}^{(0)}+{\cfrac {x-x_{2}}{\rho _{3}^{(0)}-\rho _{1}^{(0)}+{\cfrac {x-x_{3}}{\rho _{4}^{(0)}-\rho _{2}^{(0)}+\cdots }}}}}}}}}$

où les différences réciproques ρ⁽ⁿ⁾_k sont calculées avec la suite S_i et la suite auxiliaire x_i. La fonction de Thiele est une fraction rationnelle de degrés identiques au numérateur et dénominateur pour k pair et de degré supérieur de 1 au numérateur pour k impair. Plus précisément, on a :

${\displaystyle T_{2k}^{(0)}(x)={\frac {\rho _{2k}^{(0)}x^{k}+...}{x^{k}+...}}}$

${\displaystyle T_{2k+1}^{(0)}(x)={\frac {x^{k+1}+...}{\rho _{2k+1}^{(0)}x^{k}+...}}}$

On en déduit facilement que, lorsque k est pair, le ρ-algorithme correspond à la limite en l'infini de la fraction d'interpolation de Thiele. De plus, on peut en déduire aussi que les différences réciproques ρ⁽ⁿ⁾_k ne dépendent pas de l'ordre des points ayant servi à les calculer (puisque la fraction d'interpolation ne dépend pas de cet ordre et que la différence réciproque est un des coefficients de cette fraction).

La série suivante, série dite du problème de Bâle,

${\displaystyle \sum _{n=1}^{\infty }{\frac {1}{n^{2}}}={\frac {1}{1^{2}}}+{\frac {1}{2^{2}}}+{\frac {1}{3^{2}}}+{\frac {1}{4^{2}}}+\cdots ={\frac {\pi ^{2}}{6}}=1,64493406684822643...}$

converge très lentement. Le ρ-algorithme peut être utilisé pour accélérer la convergence des sommes partielles de cette série. Voici le résultat :

${\displaystyle \scriptstyle {\begin{matrix}&\rho _{0}^{(n)}=S_{n}&&\rho _{2}^{(n)}&&\rho _{4}^{(n)}&&\rho _{6}^{(n)}&&\rho _{8}^{(n)}\\0&&&&&&&\\&1&&&&&&\\0&&4&&&&&\\&1,25&&1,65&&&&\\0&&9&&-936&&\\&1,3611111&&1,64682540&&1,64489489&&\\0&&16&&-3008&&177552,00\\&1,42361111&&1,64583333&&1,64492259&&1,64493438\\0&&25&&-7400&&686475,00&&-29649576,9\\&1,46361111&&1,64542929&&1,64492979&&1,64493438&&1,644934064\\0&&36&&-15408&&2064852,04&&-128422187,3\\&1,49138889&&1,64523504&&1,64493220&&1,64493415\\0&&49&&-28616&&5229818,72\\&1,51179705&&1,64523504&&1,64493315\\0&&64&&-48896\\&1,52742205&&1,64513039\\0&&81\\&1,53976773\\\end{matrix}}}$

Le ρ-algorithme classique a été utilisé pour les calculs précédents. La suite à accélérer est la somme partielle de la série, en ajoutant un terme à la fois. On constate effectivement que seules les colonnes paires convergent vers la limite, et ceci beaucoup plus vite que la suite d'origine (pour obtenir la même précision avec la suite d'origine, il aurait fallu sommer plusieurs milliards de termes au lieu de 9).

Comparaison de l'interpolation polynomiale et rationnelle de la fonction Arctan(x)

Le calcul des différences réciproques (tableau de valeur du ρ-algorithme) permet d'utiliser la formule de Thiele pour le calcul de fraction rationnelle d'interpolation. L'interpolation rationnelle fournit de meilleurs résultats par rapport à l'interpolation polynomiale, lorsque la fonction à interpoler :

• présente des asymptotes ;
• présente des pôles à proximité ou dans la zone d'interpolation ;
• présente des pôles complexes au voisinage de l'intervalle d'interpolation (phénomène de Runge).

Le graphique ci-contre montre le gain que peut apporter l'interpolation rationnelle sur l'interpolation polynomiale dans certains cas. L'exemple porte sur l'interpolation de la fonction Arc-tangente à l'aide de 21 points d'interpolation régulièrement espacés sur l'intervalle [-10;10]. On constate que l'interpolation polynomiale présente un phénomène de Runge. L'interpolation rationnelle est quasiment confondue avec la fonction arc tangente dans l'intervalle d'interpolation, et même légèrement au-delà. La fonction arctangente présente des asymptotes horizontales et des pôles complexes en -i et +i, conditions non favorables pour une interpolation polynomiale. Le phénomène de Runge peut être estompé en resserrant les abscisses d'interpolation aux extrémités de l'intervalle, par exemple avec les abscisses de Chebychev : cependant l'amélioration obtenue ne suffit pas à tenir la comparaison avec l'interpolation rationnelle (erreur max de l'ordre de 0,07 contre 0,00005).

Plusieurs méthodes numériques utilisent l'extrapolation de Richardson pour accélérer la convergence de méthodes d'ordre peu élevé. On trouve par exemple son utilisation dans :

• l'évaluation de la dérivée d'une fonction, par accélération de la convergence de la méthode des différences centrales ;
• l'intégration numérique, en accélérant la méthode des trapèzes (méthode de Romberg) ;
• l'intégration des équations différentielles, en accélérant la méthode du point milieu (méthode de Gragg-Bulirsch-Stoer).

Dans chacune de ces applications, il est possible de remplacer l'extrapolation de Richardson par le ρ-algorithme. La suite associée x_n devant tendre vers l'infini dans le cas du ρ-algorithme, on prendra l'inverse de la suite associée à l'extrapolation de Richardson (ou leur carré, le cas échéant). Les colonnes impaires de l'algorithme ne convergeant pas vers la limite cherchée, la valeur à retenir pour l'extrapolation est alternativement ρ⁽⁰⁾_k pour k pair et ρ⁽¹⁾_k–1 pour k impair. Bulirsch et Stoer[3] ont proposé dans ce but une autre méthode d'interpolation rationnelle qui, pour un nombre pair de points, donne les mêmes résultats que le ρ-algorithme, mais nécessite plus de calculs.

Le ρ-algorithme s'est révélé intéressant pour divers applications pratiques, entre autres :

• l'inversion numérique de la transformation de Laplace[4] ;
• l'accélération de la convergence d'algorithmes de restauration d'images (algorithme de Richardson-Lucy).

• Claude Brezinski, Algorithmes d'accélération de la convergence : étude numérique, éditions Technip, 1978, (ISBN 2-7108-0341-0)

• ε-algorithme
• θ-algorithme