Algorithme de Bartels-Stewart

En algèbre linéaire numérique, l'algorithme de Bartels-Stewart est un algorithme utilisé pour résoudre numériquement l'équation de Sylvester ${\displaystyle AX-XB=C}$.

Historique

Développé par Richard Harold Bartels et Gilbert Wright Stewart en 1972^[1], l'algorithme était la première méthode numériquement stable qui pouvait être systématiquement appliquée pour résoudre de telles équations. L'algorithme opère en utilisant les décompositions de Schur réelles de ${\displaystyle A}$ et ${\displaystyle B}$ pour transformer l'équation ${\displaystyle AX-XB=C}$ en un système triangulaire qui peut ensuite être résolu en utilisant la substitution avant ou arrière. En 1979, Gene H. Golub, Charles F. Van Loan et Stephen Nash ont introduit une version améliorée de l'algorithme^[2] connue sous le nom d'algorithme de Hessenberg-Schur. La méthode de Bartels et Stewartt reste une approche standard pour résoudre les équations de Sylvester lorsque ${\displaystyle X}$ est de taille petite ou moyenne.

L'algorithme

Soient ${\displaystyle X,C\in \mathbb {R} ^{m\times n}}$, et supposons que les valeurs propres de ${\displaystyle A}$ sont distinctes des valeurs propres de ${\displaystyle B}$. L'équation matricielle ${\displaystyle AX-XB=C}$ a alors une solution unique. L'algorithme de Bartels-Stewart calcule ${\displaystyle X}$ en appliquant les étapes suivantes^[1] :

1. Calculer les décompositions de Schur réelles :

   ${\displaystyle R=U^{T}AU,}$

   ${\displaystyle S=V^{T}B^{T}V.}$

   Les matrices ${\displaystyle R}$ et ${\displaystyle S}$ sont des matrices triangulaires par blocs, avec des blocs diagonaux de taille 1 ou 2.

2. Calculer ${\displaystyle F=U^{T}CV.}$
3. Résoudre le système simplifié

   ${\displaystyle RY-YS^{T}=F}$,

   où ${\displaystyle Y=U^{T}XV}$. Pour cela, on peut utiliser la substitution avant par blocs. Concrètement, si ${\displaystyle s_{k-1,k}=0}$, alors on a :

   ${\displaystyle (R-s_{k,k}I)y_{k}=f_{k}+\sum _{j=k+1}^{n}s_{k,j}y_{j},}$

   où ${\displaystyle y_{k}}$ est la ${\displaystyle k}$-ième colonne de ${\displaystyle Y}$. Quand ${\displaystyle s_{k-1,k}\neq 0}$, les colonnes ${\displaystyle y_{k-1}}$ et ${\displaystyle y_{k}}$ doivent être concaténées et résolues simultanément.

4. Calculer ${\displaystyle X=UYV^{T}.}$

Coût du calcul

En utilisant l'algorithme QR, les décompositions réelles de Schur à l'étape 1 nécessitent environ ${\displaystyle 10(m^{3}+n^{3})}$ opérations flottantes, de sorte que le coût global du calcul global est : ${\displaystyle 10(m^{3}+n^{3})+2,5(mn^{2}+nm^{2})}$^[2].

Simplifications et cas particuliers

Dans le cas particulier où ${\displaystyle B=-A^{T}}$ et ${\displaystyle C}$ est une matrice symétrique, la solution ${\displaystyle X}$ est également symétrique. Cette symétrie peut être exploitée pour calculer ${\displaystyle Y}$ plus efficacement à l'étape 3 de l'algorithme^[1].

L'algorithme de Hessenberg-Schur

L'algorithme de Hessenberg–Schur^[2] remplace la décomposition ${\displaystyle R=U^{T}AU}$ de l'étape 1 par la décomposition

${\displaystyle H=Q^{T}AQ}$,

où ${\displaystyle H}$ est une matrice de Hessenberg supérieure. Cela conduit à un système de la forme

${\displaystyle HY-YS^{T}=F}$

qui peut être résolu en utilisant la substitution directe. L'avantage de cette approche est que ${\displaystyle H=Q^{T}AQ}$ peut être calculée en utilisant les réflexions de Householder au prix de ${\displaystyle (5/3)m^{3}}$ flops, par rapport aux ${\displaystyle 10m^{3}}$ flops nécessaires pour calculer la décomposition réelle de Schur de ${\displaystyle A}$.

Logiciel et implémentation

Les sous-programmes requis pour la variante Hessenberg-Schur de l'algorithme Bartels-Stewart sont implémentés dans la bibliothèque SLICOT. Ceux-ci sont utilisés dans la boîte à outils du système de contrôle MATLAB.

Approches alternatives

Pour les systèmes de grande taille, le coût en ${\displaystyle {\mathcal {O}}(m^{3}+n^{3})}$ le coût de l'algorithme de Bartels-Stewart peut être prohibitif. Quand ${\displaystyle A}$ et ${\displaystyle B}$ sont creuses ou structurées, de sorte que les résolutions linéaires et les multiplications vectorielles matricielles sont efficaces, les algorithmes itératifs peuvent être plus rapides. Ceux-ci incluent des méthodes basées sur la projection, qui utilisent des itérations Méthode itérative de sous-espace de Krylov, des méthodes basées sur l'itération implicite de direction alternée (ADI) et des hybridations qui impliquent à la fois la projection et l'ADI^[3]. Des méthodes itératives peuvent également être utilisées pour construire directement des approximations de rang inférieur à ${\displaystyle X}$ lors de la résolution ${\displaystyle AX-XB=C}$.

Notes et références

• (en) Cet article est partiellement ou en totalité issu de l’article de Wikipédia en anglais intitulé « Bartels–Stewart algorithm » (voir la liste des auteurs).

• Portail des mathématiques