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 $AX-XB=C$.

Historique

Développé par Richard Harold Bartels et Gilbert Wright Stewart en 1972, 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 $A$ et $B$ pour transformer l'équation $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 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 $X$ est de taille petite ou moyenne.

L'algorithme

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

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

   $R=U^{T}AU,$

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

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

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

   $RY-Y{S}^T=F$,

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

   $(R-s_{k,k}I)y_k=f_k+\sum _{j=k+1}^n{s}_{k,j}{y}_j,$

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

4. Calculer $X=UY{V}^T.$

Coût du calcul

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

Simplifications et cas particuliers

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

L'algorithme de Hessenberg-Schur

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

$H=Q^{T}AQ$,

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

$HY-Y{S}^T=F$

qui peut être résolu en utilisant la substitution directe. L'avantage de cette approche est que $H=Q^{T}AQ$ peut être calculée en utilisant les réflexions de Householder au prix de $(5/3)m^3$ flops, par rapport aux $10m^3$ flops nécessaires pour calculer la décomposition réelle de Schur de $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 $𝒪(m^3+n^3)$ le coût de l'algorithme de Bartels-Stewart peut être prohibitif. Quand $A$ et $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. Des méthodes itératives peuvent également être utilisées pour construire directement des approximations de rang inférieur à $X$ lors de la résolution $AX-XB=C$.