QZ-Algorithmus

Ziel dieses Schrittes ist es, die Matrix ${\displaystyle A}$  durch orthogonale Transformationen auf obere Hessenbergform und die Matrix ${\displaystyle B}$  auf obere Dreiecksform zu bringen. Durch ${\displaystyle n-1}$  Householder-Spiegelungen von links wird ${\displaystyle B}$  auf obere Dreiecksform transformiert. Wendet man die gleichen Transformationen gleichzeitig auf ${\displaystyle A}$  an, ergibt sich (Veranschaulichung an einem Beispiel der Größe (4,4)): ${\displaystyle A={\begin{pmatrix}*&*&*&*\\{*}&*&*&*\\{*}&*&*&*\\{*}&*&*&*\end{pmatrix}},B={\begin{pmatrix}*&*&*&*\\0&*&*&*\\0&0&*&*\\0&0&0&*\end{pmatrix}}}$ .

Man finde nun eine Givens-Rotation, die von links angewendet auf A folgende Matrix ergibt: ${\displaystyle A={\begin{pmatrix}*&*&*&*\\{*}&*&*&*\\{*}&*&*&*\\0&*&*&*\end{pmatrix}}}$ . Damit erhält man für ${\displaystyle B={\begin{pmatrix}*&*&*&*\\0&*&*&*\\0&0&*&*\\0&0&*&*\end{pmatrix}}}$ .

Durch Anwendung einer Givens-Rotation von rechts kann die obere Dreiecksform von ${\displaystyle B}$  wiederhergestellt werden, ohne die Null an der linken unteren Position von A zu zerstören: ${\displaystyle A={\begin{pmatrix}*&*&*&*\\{*}&*&*&*\\{*}&*&*&*\\0&*&*&*\end{pmatrix}},B={\begin{pmatrix}*&*&*&*\\0&*&*&*\\0&0&*&*\\0&0&0&*\end{pmatrix}}}$ .

Durch analoges spaltenweises Erzeugen von Nullen in ${\displaystyle A}$  erhält man eine obere Hessenbergmatrix:

1. ${\displaystyle A={\begin{pmatrix}*&*&*&*\\{*}&*&*&*\\0&*&*&*\\0&*&*&*\end{pmatrix}},B={\begin{pmatrix}*&*&*&*\\0&*&*&*\\0&*&*&*\\0&0&0&*\end{pmatrix}}}$ 
2. ${\displaystyle A={\begin{pmatrix}*&*&*&*\\{*}&*&*&*\\0&*&*&*\\0&*&*&*\end{pmatrix}},B={\begin{pmatrix}*&*&*&*\\0&*&*&*\\0&0&*&*\\0&0&0&*\end{pmatrix}}}$ 
3. ${\displaystyle A={\begin{pmatrix}*&*&*&*\\{*}&*&*&*\\0&*&*&*\\0&0&*&*\end{pmatrix}},B={\begin{pmatrix}*&*&*&*\\0&*&*&*\\0&0&*&*\\0&0&*&*\end{pmatrix}}}$ 
4. ${\displaystyle A={\begin{pmatrix}*&*&*&*\\{*}&*&*&*\\0&*&*&*\\0&0&*&*\end{pmatrix}},B={\begin{pmatrix}*&*&*&*\\0&*&*&*\\0&0&*&*\\0&0&0&*\end{pmatrix}}}$ .

Falls ${\displaystyle (A,B)}$ -invariante Unterräume berechnet werden sollen, so ist es notwendig, das Produkt der Transformationsmatrizen, die jeweils von links auf ${\displaystyle A}$  und ${\displaystyle B}$  angewendet werden, in einer Matrix ${\displaystyle Q}$  und das Produkt der Transformationsmatrizen, die von rechts angewendet werden, in einer Matrix ${\displaystyle Z}$  zu speichern.

1. ${\displaystyle q:=0}$ 

2. while ${\displaystyle q<n}$  do

3. Bestimme alle ${\displaystyle j\in \{1,\cdots ,n-1\}}$  mit ${\displaystyle |a_{j+1,j}|\leq \varepsilon (|a_{j,j}|+|a_{j+1,j+1}|)}$ . Für diese ${\displaystyle j}$  setze ${\displaystyle a_{j,j+1}=0}$ .

4. Deflation: Finde minimales ${\displaystyle p}$  und maximales ${\displaystyle q}$  mit ${\displaystyle p,q\in \{1,\cdots ,n\}}$  und definiere ${\displaystyle m:=n-p-q}$ , so dass gilt: ${\displaystyle A={\begin{pmatrix}A_{11}&A_{12}&A_{13}\\0&A_{22}&A_{23}\\0&0&A_{33}\end{pmatrix}}}$ , wobei ${\displaystyle A_{11}\in \mathbb {R} ^{p\times p},A_{22}\in \mathbb {R} ^{m\times m},A_{33}\in \mathbb {R} ^{q\times q}}$  und ${\displaystyle A_{11}}$  von oberer Hessenbergform, ${\displaystyle A_{22}}$  von unreduzierter oberer Hessenbergform und ${\displaystyle A_{33}}$  von quasi-oberer Dreiecksform ist.

5. Partitioniere ${\displaystyle B}$  wie ${\displaystyle A}$ , d. h. ${\displaystyle B={\begin{pmatrix}B_{11}&B_{12}&B_{13}\\0&B_{22}&B_{23}\\0&0&B_{33}\end{pmatrix}}}$ , wobei ${\displaystyle B_{11}\in \mathbb {R} ^{p\times p},B_{22}\in \mathbb {R} ^{m\times m},B_{33}\in \mathbb {R} ^{q\times q}}$  obere Dreiecksmatrizen sind.

6. Bringe ${\displaystyle A_{33}}$  in obere Schurform: Finde orthogonale ${\displaystyle Q_{33},Z_{33}}$  so, dass ${\displaystyle A_{33}:=Q_{33}^{T}A_{33}Z_{33}}$  in Schurform und ${\displaystyle B_{33}:=Q_{33}^{T}B_{33}Z_{33}}$  obere Dreiecksmatrix ist.

Falls erforderlich: Aufdatieren von ${\displaystyle Q}$  und ${\displaystyle Z}$ : ${\displaystyle Q:=Q\mathrm {diag} (I_{p},I_{m},Q_{33})}$ , ${\displaystyle Z:=Z\mathrm {diag} (I_{p},I_{m},Z_{33})}$ .

7. if ${\displaystyle q<n}$ :

if ${\displaystyle det(B_{22})=0}$ 

Transformiere mithilfe einer Givens-Rotation von rechts ${\displaystyle a_{n-q,n-q-1}=0}$ , um die Rang-Defizienz von ${\displaystyle B_{22}}$  auf ${\displaystyle B_{33}}$  zu verschieben. Durch die Annullierung von ${\displaystyle a_{n-q,n-q-1}}$  ist ${\displaystyle A_{22}}$  keine unreduzierte Hessenbergmatrix mehr, somit wird ${\displaystyle q}$  erhöht und es besteht die Möglichkeit, dass ${\displaystyle B_{22}}$  in der neuen Partitionierung regulär ist.

else

Führe einen impliziten QZ-Schritt für ${\displaystyle A_{22},B_{22}}$  aus: ${\displaystyle A_{22}:=Q_{22}^{T}A_{22}Z_{22},\quad {B_{22}}:=Q_{22}^{T}{B_{22}}Z_{22}}$ .

end if

8. end if

9. Bestimme Shifts ${\displaystyle a,b}$  als Eigenwerte von ${\displaystyle {\begin{pmatrix}a_{m-1,m-1}&a_{m-1,m}\\a_{m,m-1}&a_{m,m}\end{pmatrix}}{\begin{pmatrix}b_{m-1,m-1}&b_{m-1,m}\\0&b_{m,m}\end{pmatrix}}^{-1}}$ 

10. Bestimme ${\displaystyle (A_{22}B_{22}^{-1}-aI)(A_{22}B_{22}^{-1}-bI)e_{1}={\begin{pmatrix}\alpha \\\beta \\\gamma \\0\\\vdots \\0\end{pmatrix}}}$ 

11. Finde orthogonales ${\displaystyle Q_{1}}$  mit ${\displaystyle Q_{1}^{T}{\begin{pmatrix}\alpha \\\beta \\\gamma \end{pmatrix}}={\begin{pmatrix}*\\0\\0\end{pmatrix}}}$ 

Für ${\displaystyle B_{22}}$  folgt nun: ${\displaystyle {\begin{pmatrix}Q_{1}^{T}&0\\0&I_{m-3}\end{pmatrix}}B_{22}={\begin{pmatrix}*&*&*&\cdots &\cdots &*\\{*}&*&*&\cdots &\cdots &*\\{*}&*&*&\cdots &\cdots &*\\0&0&0&\ddots &&\vdots \\\vdots &\vdots &\vdots &&\ddots &\vdots \\0&0&0&\cdots &\cdots &*\end{pmatrix}}}$ .

Ziel ist es nun, die Dreiecksgestalt von ${\displaystyle B_{22}}$  durch orthogonale Transformationen (Householder-Spiegelungen) von rechts wiederherzustellen:

12. Finde orthogonales ${\displaystyle Z_{1}\in \mathbb {R} ^{3\times 3}}$  mit ${\displaystyle B_{22}\mathrm {diag} (Z_{1},I_{m-3})={\begin{pmatrix}*&*&*&\cdots &\cdots &*\\{*}&*&*&\cdots &\cdots &*\\0&0&*&\cdots &\cdots &*\\0&0&0&\ddots &&\vdots \\\vdots &\vdots &\vdots &&\ddots &\vdots \\0&0&0&\cdots &\cdots &*\end{pmatrix}}}$ . Finde dann orthogonales ${\displaystyle Z_{1}'\in \mathbb {R} ^{2\times 2}}$ , so dass

${\displaystyle B_{22}\mathrm {diag} (Z_{1},I_{m-3})\mathrm {diag} (Z_{1}',I_{m-2})={\begin{pmatrix}*&*&*&\cdots &\cdots &*\\0&*&*&\cdots &\cdots &*\\0&0&*&\cdots &\cdots &*\\0&0&0&\ddots &&\vdots \\\vdots &\vdots &\vdots &&\ddots &\vdots \\0&0&0&\cdots &\cdots &*\end{pmatrix}}}$ .

Für ${\displaystyle A_{22}}$  ergibt sich nun: ${\displaystyle {\tilde {A}}_{22}:={A_{22}}\mathrm {diag} (Z_{1},I_{m-3})\mathrm {diag} (Z'_{1},I_{m-2})={\begin{pmatrix}*&*&*&\cdots &\cdots &\cdots &*\\{*}&*&*&\cdots &\cdots &\cdots &*\\{*}&*&*&\cdots &\cdots &\cdots &*\\{*}&*&*&\cdots &\cdots &\cdots &*\\0&0&0&\ddots &&&\vdots \\\vdots &\vdots &\vdots &&\ddots &&\vdots \\0&0&0&\cdots &0&*&*\end{pmatrix}}}$ . D.h., die Hessenbergstruktur von ${\displaystyle A_{22}}$  ist durch einen unerwünschten 2x2 „Buckel“ zerstört.

13. Dieser Buckel kann durch elementäre, orthogonale Transformationen (z. B. Householder-Spiegelungen) von links eliminiert werden. Finde also orthogonales ${\displaystyle Q''_{1}\in \mathbb {R} ^{3\times 3}}$ , ${\displaystyle Q_{1}'\in \mathbb {R} ^{3\times 3}}$  mit

${\displaystyle \mathrm {diag} (1,Q'_{1},I_{m-4})^{T}\mathrm {diag} (I_{2},Q''_{1},I_{m-5})^{T}{{\tilde {A}}_{22}}={\begin{pmatrix}*&*&*&\cdots &\cdots &\cdots &*\\*&*&*&\cdots &\cdots &\cdots &*\\0&*&*&\ddots &\cdots &\cdots &*\\0&0&*&\ddots &\cdots &\cdots &*\\0&0&0&*&&\vdots \\\vdots &\vdots &\vdots &&\ddots &&\vdots \\0&0&0&\cdots &0&*&*\end{pmatrix}}}$ . Es werden also nacheinander die Vektoren ${\displaystyle {\begin{pmatrix}a_{21}\\a_{31}\\a_{41}\end{pmatrix}}}$  und ${\displaystyle {\begin{pmatrix}a_{32}\\a_{42}\\a_{52}\end{pmatrix}}}$ auf ${\displaystyle {\begin{pmatrix}*\\0\\0\end{pmatrix}}}$ transformiert.

Die Anwendung der Transformation auf ${\displaystyle {\tilde {B}}_{22}}$ von links ergibt jedoch

${\displaystyle \mathrm {diag} (1,Q'_{1},I_{m-4})^{T}\mathrm {diag} (I_{2},Q''_{1},I_{m-5})^{T}{{\tilde {B}}_{22}}={\begin{pmatrix}*&*&*&\cdots &\cdots &\cdots &*\\0&*&*&\cdots &\cdots &\cdots &*\\0&*&*&\ddots &\cdots &\cdots &*\\0&*&*&*&\ddots &\cdots &*\\0&0&0&0&*&\vdots \\\vdots &\vdots &\vdots &&\ddots &&\vdots \\0&0&0&\cdots &0&0&*\end{pmatrix}}}$ , d. h. ein Buckel ist jetzt eine Position tiefer entlang der Diagonalen entstanden.

14. Man wiederhole die Schritte 11–13 so lange, bis ${\displaystyle A_{22}}$  wieder in oberer Hessenberg- und ${\displaystyle B_{22}}$  wieder in oberer Dreieckstruktur vorliegt. Diesen Prozess bezeichnet man, analog zum QR-Algorithmus, auch als „Buckel-Jagen“ oder „Bulge-Chasing“. Die Eliminierung eines Buckels in ${\displaystyle B_{22}}$  an der Diagonalposition j mit Transformationen von links führt zu einem Buckel an der entsprechenden Position in ${\displaystyle A_{22}}$ . Wird dieser Buckel mit Transformationen von rechts eliminiert, führt das zu einem Buckel in ${\displaystyle B_{22}}$  an der Diagonalposition j+1 usw.

15. Nach ${\displaystyle m-2}$  Schritten wird das Ziel erreicht und es ergibt sich ${\displaystyle Q_{22}^{T}=\mathrm {diag} (Q_{1},I_{m-3})^{T}\mathrm {diag} (1,Q'_{1},I_{m-4})^{T}\mathrm {diag} (I_{2},Q_{1}'',I_{m-5})^{T}\cdots \mathrm {diag} (I_{m-3},Q_{m-2})^{T}}$ . Analog erhält man

${\displaystyle Z_{22}=\mathrm {diag} (Z_{1},I_{m-3})\mathrm {diag} (Z_{1}',I_{m-2})\cdots \mathrm {diag} (I_{m-2},Z_{m-2})\mathrm {diag} (I_{m-2},Z_{m-2}')}$ .

Falls ${\displaystyle (A,B)}$ -invarianten Unterräume benötigt werden, ist es notwendig die Matrizen ${\displaystyle {Q}}$  und ${\displaystyle Z}$  aufzudatieren: ${\displaystyle Q:=Q\mathrm {diag} (I_{p},Q_{22},I_{q})}$ , ${\displaystyle Z:=Z\mathrm {diag} (I_{p},Z_{22},I_{q})}$ 

16. end while

In den meisten Fällen konvergiert ${\displaystyle (A,B)}$  im QZ-Algorithmus gegen seine verallgemeinerte, reelle Schur-Form. Für skalare Diagonalblöcke in A gilt ${\displaystyle \lambda _{i}={\frac {a_{ii}}{b_{ii}}}:b_{ii}\neq 0}$  und ${\displaystyle \lambda _{i}=\infty }$  falls ${\displaystyle b_{ii}=0}$ . Falls ein ${\displaystyle i}$  existiert, für das ${\displaystyle a_{ii}=b_{ii}=0}$  ist, so ist ${\displaystyle \Lambda (A,B)=\mathbb {C} }$ . ${\displaystyle 2\times 2}$  Diagonalblöcke von ${\displaystyle A}$  beziehen sich (analog zum QR-Algorithmus) auf Paare komplex konjugierter Eigenwerte ${\displaystyle \lambda ,{\overline {\lambda }}=\Lambda \left({\begin{pmatrix}a_{ii}&a_{i,i+1}\\a_{i+1,i}&a_{i+1,i+1}\end{pmatrix}},{\begin{pmatrix}b_{ii}&b_{i,i+1}\\0&b_{i+1,i+1}\end{pmatrix}}\right)}$ .

• Gene H. Golub, Charles F. Van Loan: Matrix Computations. Johns Hopkins University Press, 1996, ISBN 0-8018-5414-8.
• G. W. Stewart: Matrix Algorithms. Band II: Eigensystems. SIAM 2001, ISBN 0-89871-503-2.