Transformation de Shanks
# Transformation de Shanks

En analyse numérique, la transformation de Shanks est une méthode non linéaire d'accélération de la convergence de suites numériques. Cette méthode est nommée d'après Daniel Shanks, qui l'exposa en 1955, bien qu'elle ait été étudiée et publiée par R. J. Schmidt dès 1941. C'est une généralisation de l'algorithme Delta-2 d'Aitken.

## Présentation

Soit une suite numérique $(A_n)$ dont on cherche à connaitre la limite A. La transformée de Shanks d'ordre k de cette suite est donnée par le rapport de deux déterminants :

    $e_{k}(A_{n})={\frac {\begin{vmatrix}A_{n}&\cdots &A_{n+k-1}&A_{n+k}\\\Delta A_{n}&\cdots &\Delta A_{n+k-1}&\Delta A_{n+k}\\\Delta A_{n+1}&\cdots &\Delta A_{n+k}&\Delta A_{n+k+1}\\\vdots &&\vdots &\vdots \\\Delta A_{n+k-1}&\cdots &\Delta A_{n+2k-2}&\Delta A_{n+2k-1}\\\end{vmatrix}}{\begin{vmatrix}1&\cdots &1&1\\\Delta A_{n}&\cdots &\Delta A_{n+k-1}&\Delta A_{n+k}\\\Delta A_{n+1}&\cdots &\Delta A_{n+k}&\Delta A_{n+k+1}\\\vdots &&\vdots &\vdots \\\Delta A_{n+k-1}&\cdots &\Delta A_{n+2k-2}&\Delta A_{n+2k-1}\\\end{vmatrix}}}={\frac {\begin{vmatrix}A_{n}&\cdots &A_{n+k-1}&A_{n+k}\\A_{n+1}&\cdots &A_{n+k}&A_{n+k+1}\\A_{n+2}&\cdots &A_{n+k+1}&A_{n+k+2}\\\vdots &&\vdots &\vdots \\A_{n+k}&\cdots &A_{n+2k-1}&A_{n+2k}\\\end{vmatrix}}{\begin{vmatrix}\Delta ^{2}A_{n}&\cdots &\Delta ^{2}A_{n+k-2}&\Delta ^{2}A_{n+k-1}\\\Delta ^{2}A_{n+1}&\cdots &\Delta ^{2}A_{n+k-1}&\Delta ^{2}A_{n+k}\\\Delta ^{2}A_{n+2}&\cdots &\Delta ^{2}A_{n+k}&\Delta ^{2}A_{n+k+1}\\\vdots &&\vdots &\vdots \\\Delta ^{2}A_{n+k-1}&\cdots &\Delta ^{2}A_{n+2k-3}&\Delta ^{2}A_{n+2k-2}\\\end{vmatrix}}},$

avec $\Delta A_{p}=A_{p+1}-A_{p}$ et $\Delta ^{2}A_{p}=\Delta A_{p+1}-\Delta A_{p}=A_{p+2}-2A_{p+1}+A_{p}.$

## Propriétés

### Modèle de convergence

La transformée de Shanks d'ordre k donne exactement la limite de la suite d'origine $A_n$ si celle-ci vérifie l'équation aux différences linéaire à coefficients constants d'ordre k du type :

    $\sum _{p=0}^{k}a_{p}(A_{n+p}-A)=0\quad \forall n$ avec les constantes arbitraires $a_p$ telles que $\sum _{p=0}^{k}a_{p}\neq 0.$

Ce modèle de comportement de convergence possède 2k + 1 inconnues. En évaluant l'équation ci-dessus pour $A_n,$ $A_{n+1},$ ..., $A_{n+2k},$ on obtient l'expression de la transformée de Shanks $e_k(A_n)$ en résolvant l'inconnue A dans le système obtenu.

En résolvant l'équation aux différences linéaire, on peut expliciter la forme de $A_n$ pour laquelle la transformée de Shanks d'ordre k fournit la limite exacte :

    $A_{n}=A+\sum _{p=1}^{m}Q_{p}^{(k_{p})}(n)\lambda _{p}^{n}.$

Les $λ_p$ étant les m racines du polynôme $P(\lambda )=\sum _{p=0}^{k}a_{p}\lambda ^{p}$, $Q_p^{kp}(n)$ étant un polynôme arbitraire en n, dont le degré $k_p$ est la multiplicité de la racine $λ_p.$ Les racines pouvant être complexes, le comportement de $A_n$ englobe de nombreuses formes : combinaisons linéaires d'exponentielles décroissantes ou croissantes avec des sinusoïdes amorties ou amplifiées.

### Lien avec la table de Padé

Lorsque la suite à accélérer est la somme partielle d'un développement en série entière d'une fonction :

    $A_{n}=\sum _{k=1}^{n}c_{k}x^{k}\quad f(x)=\sum _{k=1}^{\infty }c_{k}x^{k}$

alors, les différentes transformées de Shanks de la suite constituent le tableau de Padé de la série :

$e_{k}(A_{n})=[n+k/k]_{f(x)}\quad k=0,1...\quad n=-k,-k+1...$

avec $[n+k/k]_{f(x)}$ l'approximant de Padé d'indice (k + n, k) de f(x) (la fraction rationnelle de degré n + k au numérateur et k au dénominateur dont le développement limité coïncide avec celui de f(x) jusqu'au degré 2k + n). La transformation de Shanks fournit donc un moyen de calculer les différents approximants de Padé d'une fonction dont on connait le développement limité.

### Cas de convergence démontrée

L'accélération de la convergence de la transformation de Shanks a été démontrée pour certains types de suite, notamment les suites totalement monotones $(A_n$ est totalement monotone si $(-1)^{k}\Delta ^{k}A_{n}\geqslant 0\quad \forall n,k=0,1...$) et les suites totalement oscillantes $(A_n$ est totalement oscillante si $(-1)^{n}A_{n}$ est totalement monotone).

S'il existe des constantes a et b telles que la suite a $A_n$ + b soit totalement monotone ou oscillante, et si

    $\lim _{n\to \infty }{\frac {A_{n+1}-A}{A_{n}-A}}\neq 1$

(la convergence n'est pas de type logarithmique) alors les lignes, colonnes et diagonales du tableau des transformées de Shanks convergent vers la limite de la suite d'origine, et cela plus vite que, respectivement, les lignes, colonnes et diagonales précédentes.

## Algorithmes de calcul

En pratique, le calcul des déterminants étant assez gourmand, la transformée de Shanks n'est effectivement calculée de cette manière que pour de faibles valeurs de k, notamment pour k = 1 où l'on obtient le $Δ^2$ d'Aitken. En effet : $e_{1}(A_{n})={\frac {\begin{vmatrix}A_{n}&A_{n+1}\\\Delta A_{n}&\Delta A_{n+1}\\\end{vmatrix}}{\begin{vmatrix}1&1\\\Delta A_{n}&\Delta A_{n+1}\\\end{vmatrix}}}={\frac {A_{n}A_{n+2}-A_{n+1}^{2}}{A_{n+2}-2A_{n+1}+A_{n}}}$

La deuxième expression donnée plus haut de la transformée est un rapport de deux déterminants possédant une structure particulière, et qui entrent dans la catégorie des déterminants de Hankel. Il existe une relation de récurrence entre les déterminants de Hankel qui permet ainsi de calculer la transformée de Shanks plus économiquement.

Mais surtout, la méthode la plus utilisée pour calculer la transformée de Shanks est l'ε-algorithme de Peter Wynn :

    $\forall n,\ \varepsilon _{-1}^{(n)}=0\quad \varepsilon _{2n}^{(-n-1)}=0\quad \varepsilon _{0}^{(n)}=A_{n}$
    $\forall n,k,\ \varepsilon _{k+1}^{(n)}=\epsilon _{k-1}^{(n+1)}+{\frac {1}{\epsilon _{k}^{(n+1)}-\varepsilon _{k}^{(n)}}}$
    $\forall n,k,\ e_{k}(A_{n})=\epsilon _{2k}^{(n)}.$

L'ε-algorithme se prête particulièrement bien au calcul dans un tableur. Dans l'ε-algorithme, seuls les indices k pairs fournissent des résultats intéressants, les autres servant d'intermédiaires de calcul.

Plusieurs autres algorithmes peuvent servir à calculer les transformées de Shanks (q-d algorithme de Rutishauser, η-algorithme de Bauer...). La règle de la croix peut être aussi utilisée :

    $\forall n,\ e_{-1}(A_{n})=\infty ,\quad e_{n}(A_{-n-1})=0,\quad e_{0}(A_{n})=A_{n},$
    $\forall n,k,\ {\frac {1}{e_{k+1}(A_{n-1})-e_{k}(A_{n})}}+{\frac {1}{e_{k-1}(A_{n+1})-e_{k}(A_{n})}}={\frac {1}{e_{k}(A_{n+1})-e_{k}(A_{n})}}+{\frac {1}{e_{k}(A_{n-1})-e_{k}(A_{n})}}$.

## Exemples

### Extension du domaine de convergence d'un développement en série entière

Le développement en série entière de la fonction arc tangente

    $\arctan x=\sum _{n=0}^{\infty }{\frac {(-1)^{n}x^{2n+1}}{2n+1}}$

a un rayon de convergence égal à 1. Le calcul des transformées de Shanks de ce développement va fournir le tableau de Padé de la fonction arc tangente, dont certains éléments peuvent converger au-delà de 1. Par exemple, pour x = 2, on a :

Application de la transformée de Shanks à arctan(2) = 1,10714871... :

| n | $A_n$ | $e_1(A_n)$ | $e_2(A_n)$ | $e_3(A_n)$ | $e_4(A_n)$ | $e_5(A_n)$ |
| - | - | - | - | - | - | - |
| 0 | 2,000000 | 1,215686 | 1,122699 | 1,109420 | 1,107481 | 1,107197 |
| 1 | –0,666667 | 0,992593 | 1,095083 | 1,105663 | 1,106953 | |
| 2 | 5,733333 | 1,285457 | 1,120862 | 1,108523 | 1,107305 | |
| 3 | –12,552381 | 0,762040 | 1,087293 | 1,105534 | | |
| 4 | 44,336508 | 1,873988 | 1,141097 | 1,109410 | | |
| 5 | –141,845310 | –0,766091 | 1,041663 | | | |
| 6 | 488,308536 | 6,008969 | 1,245533 | | | |
| 7 | –1696,224797 | –12,406001 | | | | |
| 8 | 6013,892850 | 39,911298 | | | | |
| 9 | –21580,212414 | | | | | |
| 10 | 78284,168539 | | | | | |

Malgré la divergence manifeste de la série à accélérer, on constate que sa transformée de Shanks converge vers la valeur de arctan(2). On obtient ainsi une partie de la table de Padé de la fonction arc tangente (l'autre partie peut être calculée en continuant l'ε-algorithme pour les n négatifs). Dans cet exemple, la meilleure estimation étant obtenue pour $e_k(A_0),$ correspondant à la diagonale du tableau de Padé. Il est important de calculer la suite $A_n$ ainsi que les transformées successives en double précision (ou plus) afin de minimiser les propagations d'erreur d'arrondi dont l'algorithme a tendance à amplifier.

### Application répétée et itérée

La transformation de Shanks est une transformation de suite à suite : le résultat est donc une suite (qui contient moins de terme que la suite originale). Il peut venir à l'idée de refaire subir à la suite résultante une nouvelle transformation de Shanks : on appelle ceci une application répétée de la transformation de Shanks.

En reprenant la première ligne $e_k(A_0)$ du tableau précédent comme nouvelle suite initiale, on obtient :

Application répétée de la transformée de Shanks à arctan(2) = 1,10714871... :

| $e_k(A_0)$ | $e_1(e_k(A_0))$ | $e_2(e_k(A_0))$ |
| - | - | - |
| 2,000000 | 1,11019129 | 1,10714844 |
| 1,215686 | 1,10720759 | 1,10714867 |
| 1,122699 | 1,10714937 | |
| 1,109420 | 1,10714868 | |
| 1,107481 | | |
| 1,107197 | | |

On constate une amélioration du résultat, ceci sans calculer de nouveaux termes de la série d'arc tangente (les calculs ont été effectués avec plus de chiffres significatifs que ceux affichés).

Une autre possibilité consiste à limiter l'ordre k des transformées. La dernière colonne calculée $e_k(A_n)$ servira alors de suite initiale pour une nouvelle transformation, poussée jusqu'à l'ordre k, et ainsi de suite. On appelle ceci une application itérée de la transformée de Shanks d'ordre k. On prend le plus souvent k = 1, c'est-à-dire que l'on itère le $Δ^2$ d'Aitken.

En reprenant l'exemple de arctan(2), on trouve :

Application itérée de $e_1(A_n)$ $(Δ^2$ d'Aitken) à arctan(2) = 1,10714871... :

| n | $A_n$ | $e_1(A_n)$ | $e_1^2(A_n)$ | $e_1^3(A_n)$ | $e_1^4(A_n)$ | $e_1^5(A_n)$ |
| - | - | - | - | - | - | - |
| 0 | 2,000000 | 1,215686 | 1,119223 | 1,108112 | 1,107203 | 1,107151 |
| 1 | –0,666667 | 0,992593 | 1,097666 | 1,106475 | 1,107111 | |
| 2 | 5,733333 | 1,285457 | 1,117931 | 1,107786 | 1,107181 | |
| 3 | –12,552381 | 0,762040 | 1,091576 | 1,106394 | | |
| 4 | 44,336508 | 1,873988 | 1,133689 | 1,108205 | | |
| 5 | –141,845310 | –0,766091 | 1,056116 | | | |
| 6 | 488,308536 | 6,008969 | 1,214677 | | | |
| 7 | –1696,224797 | –12,406001 | | | | |
| 8 | 6013,892850 | 39,911298 | | | | |
| 9 | –21580,212414 | | | | | |
| 10 | 78284,168539 | | | | | |

Dans ce cas précis, l'application itérée du $Δ^2$ d'Aitken donne de meilleurs résultats que la transformée de Shanks équivalente. Les applications répétées et itérées sont particulièrement sensibles à la propagation d'erreurs d'arrondi.

### Accélération de la convergence d'une série de Fourier

Les sommes partielles d'une série de Fourier peuvent fournir une suite qui peut être accélérée par la transformée de Shanks. Cependant, une manière plus logique de procéder est de transformer, par un changement de variable, la série de Fourier en une série entière. Les transformées de Shanks de cette série entière en fourniront sa table de Padé.

La série de Fourier :

    $f(\theta )={\frac {a_{0}}{2}}+\sum _{n=1}^{\infty }(a_{n}\cos {n\theta }+b_{n}\sin {n\theta })$

se transforme en une série entière

    $f(z)={\frac {a_{0}}{2}}+\sum _{n=1}^{\infty }(a_{n}-\mathrm {i} b_{n})z^{n}{\text{ avec }}z=\mathrm {e} ^{\mathrm {i} \theta }$

dont la partie réelle est la valeur cherchée.

Par exemple, la fonction en dents de scie sur [–π, +π], dont la série de Fourier vaut :

    $f(\theta )=\sum _{n=1}^{\infty }(-1)^{n+1}{\frac {2}{n}}\sin {n\theta }$

devient après changement de variable la série entière :

    $f(z)=\sum _{n=1}^{\infty }(-1)^{n}{\frac {2\mathrm {i} }{n}}z^{n}{\text{ avec }}z=\mathrm {e} ^{\mathrm {i} \theta }$

dont on pourra accélérer la convergence en utilisant la transformation de Shanks (avec des nombres complexes).

Avec les six premiers termes de la série de Fourier :

    $S_{k}=\sum _{n=1}^{k}(-1)^{n+1}{\frac {2}{n}}\sin {n\theta }$

on trouve, en explicitant la transformée de Shanks :

    $e_{3}(S_{k})={\frac {22z^{3}+120z^{2}+120z}{3\mathrm {i} z^{3}+36\mathrm {i} z^{2}+90\mathrm {i} z+60\mathrm {i} }}$

dont la partie réelle, en revenant à la variable θ, forme une fraction rationnelle trigonométrique :

    ${\frac {440\sin {3\theta }+2940\sin {2\theta }+4704\sin {\theta }}{120\cos {3\theta }+1620\cos {2\theta }+5832\cos {\theta }+4335}}$.

On constate sur le graphique ci-contre une nette accélération de la convergence de la transformée de Shanks par rapport à la série initiale.

