En analyse numérique, lε-algorithme est un algorithme non linéaire d'accélération de la convergence de suites numériques. Cet algorithme a été proposé par Peter Wynn P. Wynn, « On a device for computing the em(Sn) transformation », MTAC, 10 (1956), . en 1956 pour calculer la transformation de Shanks. C'est l'un des algorithmes d'accélération de la convergence les mieux connus et les plus utilisés. C'est une généralisation de l'algorithme Delta-2 d'Aitken.
Présentation
Soit une suite numérique dont on cherche à connaitre la limite . l'ε-algorithme consiste à calculer un tableau en initialisant la première colonne par la suite , et en calculant les autres cases à l'aide des relations suivantes :
Les différentes colonnes d'indice pair du tableau fournissent d'autres suites convergeant en général plus vite que la suite d'origine. Les colonnes d'indice impair peuvent être considérées comme des intermédiaires de calcul.
On présente traditionnellement les résultats de l'ε-algorithme sous forme d'un tableau aux lignes décalées : la formule de calcul de l'ε-algorithme relie ainsi les termes formant un losange dans le tableau.
Propriétés
Le tableau obtenu par l'ε-algorithme est une des méthodes utilisées pour calculer la transformée de Shanks d'une suite. En effet :
où est la transformée de Shanks d'ordre de la suite . L'ε-algorithme partage donc toutes les propriétés de la transformée de Shanks. Notamment, l'ε-algorithme est un moyen pour calculer la table de Padé d'une série entière. De même, la suite correspond au Delta-2 d'Aitken. Les approximants de Padé à plusieurs variables peuvent aussi très simplement être calculés à l'aide de l'ε-algorithme. Par exemple, pour une fonction à deux variables, dont le développement de Taylor est :
Sa table de Padé peut être obtenue à l'aide de l'ε-algorithme, en initialisant la première colonne par :
Variantes
ε-algorithme vectoriel
Il existe plusieurs variantes de l'ε-algorithme pouvant être utilisées avec des suites vectorielles (ε-algorithme vectoriel, topologique ou scalaire appliqué à chacune des composantes du vecteur). L'ε-algorithme vectoriel est analogue dans sa forme à l'ε-algorithme classique :
où les suites et les 0 sont des vecteurs.
Reste à définir ce qu'est l'« inverse » d'un vecteur. On pourra prendre par exemple comme inverse d'un vecteur :
où la norme du vecteur est la norme euclidienne. L'ε-algorithme vectoriel a de nombreuses applications pratiques mais est très sensible aux erreurs d'arrondi. Il permet notamment de généraliser la méthode de Aitken-Steffensen à des systèmes d'équations non linéaires. Il fournissant ainsi un algorithme à convergence quadratique ne nécessitant pas de calcul du jacobien, à l'inverse de la méthode de Newton.
L'ε-algorithme matriciel existe pour les matrices carrées, dans ce cas l'inverse apparaissant dans l'algorithme est tout simplement celui d'une matrice carrée classique.
ε-algorithme confluent
L'ε-algorithme classique accélère une suite numérique, c'est-à-dire une fonction d'une variable discrète . La version confluente de cet algorithme accélère la convergence d'une fonction d'une variable continue . On l'obtient en faisant un changement de variable à l'algorithme d'origine, et en faisant tendre vers 0. on obtient :
Cet algorithme est utilisé par exemple pour l'évaluation des intégrales impropres lorsque les dérivées de la fonction sont accessibles, ou à l'élaboration de séries asymptotiques pour certaines fonctions.
La dérivée des termes apparaissant dans la formule peuvent se ramener aux dérivées de la fonction de base par calcul formel, où en utilisant les relations des déterminants de Henkel :
avec
et
On trouve, en explicitant les premières expressions :
ε-algorithme dépendant d'une suite auxiliaire
Il arrive fréquemment que la suite à accélérer dépende d'une suite auxiliaire , celle-ci tendant vers l'infini (par exemple des approximations par discrétisation avec un maillage de plus en plus fin). Il est possible d'utiliser l'ε-algorithme de base dans ces cas mais la façon dont évolue la suite associée (par exemple l'évolution du nombre de mailles entre chaque ) est un renseignement précieux qui pourrait aider l'accélération de la convergence. Deux variantes de l'ε-algorithme utilisent cette information supplémentaire et donnent souvent de meilleurs résultats.
Et la deuxième variante :
Ces algorithmes s'apparentent aux méthodes par extrapolation (extrapolation de Richardson ou ρ-algorithme). Ils peuvent accélérer des suites à convergence logarithmique, que l'ε-algorithme classique n'accélère pas.
ε-algorithme d'interpolation
L'ε-algorithme classique permet de calculer la table de Padé d'une fonction lorsque les termes de la suite sont les sommes partielles du développement limité de cette fonction. Une variante. permet de construire la table des fractions rationnelles d'interpolation d'une fonction en partant de la suite des polynômes d'interpolation.
Le polynôme étant la valeur en x du polynôme d'interpolation passant par les points , que l'on pourra calculer par la méthode de Lagrange, de Newton où de Neville, par exemple. est la fraction rationnelle de degré au numérateur et au dénominateur passant par les points .
Cet algorithme peut se présenter comme un complément aux méthodes de Thiele ou Bulirsch-Stoer, avec la particularité que l'ε-algorithme permet de calculer tout le tableau des fractions rationnelles d'interpolation au lieu d'une seule diagonale, ou un escalier. De plus, si le polynôme utilisé est un polynôme d'interpolation d'Hermite (imposant d'interpoler la fonction, mais aussi sa ou ses dérivées), l'ε-algorithme d'interpolation fournira les fraction rationnelle d'interpolation vérifiant les mêmes critères que le polynôme d'Hermite. La méthode d'interpolation de Newton avec confluence des abscisses (pour les points où les dérivées sont à renseigner) est toute indiquée pour calculer le polynôme d'Hermite, d'autant que la suite des abscisses utilisée (avec les confluences) sera nécessaire pour l'ε-algorithme d'interpolation.
Règles particulières
Il peut arriver lors du calcul du tableau de l'ε-algorithme, qu'à certains endroits se produisent des divisions par 0 ou presque. Ceci peut générer des dégradations de précision pour le calcul des zones du tableau dépendant de la case suspecte, voir des plantages de programme. Plusieurs auteurs ont développé des algorithmes spéciaux. pour contourner les cases problématiques en cas de détection de risque de division par 0. Ces algorithmes ne fonctionnent que si les cases singulières ne sont pas trop proches les unes des autres.
Exemples
Voir aussi les exemples de la transformée de Shanks.
Accélération de la convergence d'une suite
La méthode de Bernoulli permet, sous certaines conditions, d'évaluer la plus grande (ou la plus petite) racine d'un polynôme donné. Mais la suite générée par la méthode peut converger lentement vers la valeur de la racine : l'ε-algorithme permet d'accélérer cette convergence. Par exemple, avec le polynôme dont les racines sont le nombre d'or φ et 1/φ, la suite générée par la méthode de Bernoulli donne la suite de Fibonacci dont le ratio des termes successifs converge vers φ=1,6180339887499...
Dans cet exemple, seuls les termes d'indice haut positif ont été calculés.
Nous constatons dans cet exemple qu'effectivement seules les colonnes paires convergent vers la suite d'origine, et ceci plus rapidement (12 chiffres exacts au lieu de 3 dans la suite initiale).
Il est à noter que l'ε-algorithme possède une propriété particulière vis-à-vis de la méthode de Bernoulli qui permet de l'utiliser aussi dans un autre but que l'accélération de la convergence : le calcul simultané de toutes les racines du polynôme au lieu d'une seule à la fois.
Soit le polynôme :
dont les racines sont réelles, telles que |λ1|> |λ2|>...>|λp|
la méthode de Bernoulli génère la suite définie par :
dont le rapport des termes consécutifs converge vers λ1.
La suite est initialisée avec , ,... arbitraires, non tous nuls.
En appliquant l'ε-algorithme, non pas sur le ratio , mais directement sur la suite des , on a :
On constate dans ce cas que les colonnes impaires s'avèrent aussi intéressantes.
La vitesse de convergence, pour fixé, est :
On peut donc utiliser à nouveau l'ε-algorithme, mais cette fois pour accélérer la convergence de chacun de ces ratios, comme on l'a fait pour la suite .
Résolution d'un système d'équation non linéaire
L'une des utilisations de l'ε-algorithme est de fournir une généralisation de la méthode d'Aitken-Steffensen pour des systèmes d'équations non linéaires.
La méthode consiste à réécrire le système d'équation en un système de la forme (, et étant des vecteurs représentant les inconnues ou le système d'équations). En partant d'un vecteur arbitraire (de préférence proche de la solution du système), on génère la suite de vecteurs de la manière suivante :
calculer l'ε-algorithme de cette suite de vecteurs
avec étant le nombre d'équations du système.
La suite des générés par cette méthode est à convergence quadratique.
Plusieurs choix pour étant possibles, il est préférable, bien que non nécessaire, de choisir de telle sorte que la suite générée converge vers la solution du système.
L'ε-algorithme vectoriel parait le plus adapté à cet algorithme, mais l'ε-algorithme scalaire appliqué à chacune des composantes du vecteur fonctionne aussi.
Fraction rationnelle d'interpolation
Exemple: interpolation polynomiale de la fonction logarithme népérien et calcul des fraction rationnelles correspondantes à l'aide de l'ε-algorithme d'interpolation
La fonction logarithme sera interpolée en trois points: 1/2, 1 et e. Pour tester la mise en œuvre de l'algorithme dans le cas d'une interpolation d'Hermite, l'interpolation en 1/2 sera double (interpolation du logarithme et de sa dérivée), et triple en 1 (interpolation du logarithme, et de ses dérivées premières et secondes). Il s'agit d'un polynôme de degré 5 qui doit donc vérifier :
Calcul des différences divisées du polynôme d'interpolation :
+Tableau des différences divisées pour l'interpolation polynomiale de ln(x) en 3 points i xi yi / d0 d1 d2 d3 d4 d5 0 2,71828182845904 1 0,581976706869326 -0,243279819530861 0,149405165216329 -0,178413885100506 0,24817470934455 1 1 0 1 -0,5 0,545177444479562 -0,728935333122626 2 1 0 1 -0,772588722239781 0,909645111040875 3 1 0 1,38629436111989 -1,22741127776022 4 0.5 -0,693147180559945 2 5 0.5 -0,693147180559945
Les cases en italiques correspondent à des indéterminées des différences divisées (du fait des points doubles ou triples), et sont remplacés par les valeurs de la dérivée en ce point (ou de f"(x)/2! pour la différence seconde). Les cases en gras dans le tableau sont celles utiles au calcul de la forme de Newton du polynôme d'interpolation. Le polynôme s'écrit, à l'aide de la formule de Newton : .
Calcul de P5(2), et application de l'ε-algorithme d'interpolation.
La suite des Pi(2) servant à initialiser le tableau de l'ε-algorithme d'interpolation, est la suite de polynômes de degré 0, 1, 2... obtenue en retenant 1,2,3... termes du polynôme P5(2) écrit sous la forme de Newton ci dessus. Cette suite de polynôme passent respectivement par les points d'abscisse (x0,y0); (x0,y0),(x1,y1); (x0,y0),(x1,y1),(x2,y2); ...
+ Interpolation polynomiale de ln(2)=0,693147180559945..., et calcul de l'ε-algorithme d'interpolation i xi ε(i)-1 Pi(x)= ε(i)0 ε(i)1 ε(i)2 ε(i)3 ε(i)4 0 2,71828182845904 0 1 -2,39221119117733 0,705207033512282 -61,0702755830021 0,693879569827545 1 1 0 0,581976706869326 5,72267438319407 0,69023539353372 121,870014050176 0,692833802968095 2 1 0 0,756720180469139 -9,31836050756016 0,695317144633344 -146,585461084686 3 1 0 0,649405165216329 5,20217803449192 0,690925043467957 4 0.5 0 0,777556616828803 -2,49324571003365 5 0.5 0 0,51016754082086