QZ-Algorithmus

Der QZ-Algorithmus oder die QZ-Faktorisierung ist ein numerisches Verfahren zur Lösung des verallgemeinerten Eigenwertproblems.

$Ax=\lambda B{x}^{}$ , mit $A,B\in {ℝ}^{n\times n}$ bzw. $A,B\in {ℂ}^{n\times n}$

Das verallgemeinerte Eigenwertproblem ist äquivalent zum Eigenwertproblem $A{B}^{-1}y=\lambda y$, wobei $y=Bx$ und $B$ invertierbar sein muss. Es wird jedoch nicht explizit die Matrix $B^{-1}$ berechnet, um die Kondition des Problems nicht zu verschlechtern, sondern $A$ und $B$ werden simultan durch Ähnlichkeitstransformationen (Givens-Rotationen und Householder-Spiegelungen) in verallgemeinerte Schurform gebracht.

Gegeben ist ein Matrixbüschel $A-\lambda B$. Gesucht sind orthogonale Matrizen $Q$ und $Z$, so dass $Q^T(A-\lambda B)Z=T-\lambda S$ von verallgemeinerter Schurform ist, d. h. $T$ ist von quasi-oberer Dreiecksform und $S$ ist von oberer Dreiecksform. Im Fall $A,B\in {ℂ}^{n\times n}$ ist $T$ stets von oberer Dreiecksform. Aus der verallgemeinerten Schurform lassen sich dann die Eigenwerte und aus $Q$ und $Z$ $(A,B)$-invariante Unterräume des Matrixbüschels $A-\lambda B$ bestimmen.

Vortransformation

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

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

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

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

1. $A=\left(\begin{array}{cccc}*& *& *& *\\ *& *& *& *\\ 0& *& *& *\\ 0& *& *& *\end{array}\right),B=\left(\begin{array}{cccc}*& *& *& *\\ 0& *& *& *\\ 0& *& *& *\\ 0& 0& 0& *\end{array}\right)$
2. $A=\left(\begin{array}{cccc}*& *& *& *\\ *& *& *& *\\ 0& *& *& *\\ 0& *& *& *\end{array}\right),B=\left(\begin{array}{cccc}*& *& *& *\\ 0& *& *& *\\ 0& 0& *& *\\ 0& 0& 0& *\end{array}\right)$
3. $A=\left(\begin{array}{cccc}*& *& *& *\\ *& *& *& *\\ 0& *& *& *\\ 0& 0& *& *\end{array}\right),B=\left(\begin{array}{cccc}*& *& *& *\\ 0& *& *& *\\ 0& 0& *& *\\ 0& 0& *& *\end{array}\right)$
4. $A=\left(\begin{array}{cccc}*& *& *& *\\ *& *& *& *\\ 0& *& *& *\\ 0& 0& *& *\end{array}\right),B=\left(\begin{array}{cccc}*& *& *& *\\ 0& *& *& *\\ 0& 0& *& *\\ 0& 0& 0& *\end{array}\right)$.

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

QZ-Algorithmus mit impliziten Shifts

1. $q:=0$

2. while $q<n$ do

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

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

5. Partitioniere $B$ wie $A$, d. h. $B=\left(\begin{array}{ccc}{B}_{11}& B_{12}& B_{13}\\ 0& B_{22}& B_{23}\\ 0& 0& B_{33}\end{array}\right)$, wobei $B_{11}\in {ℝ}^{p\times p},B_{22}\in {ℝ}^{m\times m},B_{33}\in {ℝ}^{q\times q}$ obere Dreiecksmatrizen sind.

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

Falls erforderlich: Aufdatieren von $Q$ und $Z$: $Q:=Q\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_p,I_m,Q_{33})$, $Z:=Z\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_p,I_m,Z_{33})$.

7. if $q<n$:

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

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

else

Führe einen impliziten QZ-Schritt für $A_{22},B_{22}$ aus: $A_{22}:=Q_{22}^T{A}_{22}{Z}_{22},B_{22}:=Q_{22}^T{B}_{22}{Z}_{22}$.

end if

8. end if

Wahl der Shifts

9. Bestimme Shifts $a,b$ als Eigenwerte von $\left(\begin{array}{cc}{a}_{m-1,m-1}& a_{m-1,m}\\ a_{m,m-1}& a_{m,m}\end{array}\right){\left(\begin{array}{cc}{b}_{m-1,m-1}& b_{m-1,m}\\ 0& b_{m,m}\end{array}\right)}^{-1}$

10. Bestimme $(A_{22}{B}_{22}^{-1}-aI)(A_{22}{B}_{22}^{-1}-bI)e_1=\left(\begin{array}{c}\alpha \\ \beta \\ \gamma \\ 0\\ ⋮\\ 0\end{array}\right)$

Der implizite QZ-Schritt

11. Finde orthogonales $Q_1$ mit $Q_1^T\left(\begin{array}{c}\alpha \\ \beta \\ \gamma \end{array}\right)=\left(\begin{array}{c}*\\ 0\\ 0\end{array}\right)$

Für $B_{22}$ folgt nun: $\left(\begin{array}{cc}{Q}_1^T& 0\\ 0& I_{m-3}\end{array}\right)B_{22}=\left(\begin{array}{cccccc}*& *& *& \cdots & \cdots & *\\ *& *& *& \cdots & \cdots & *\\ *& *& *& \cdots & \cdots & *\\ 0& 0& 0& \ddots & & ⋮\\ ⋮& ⋮& ⋮& & \ddots & ⋮\\ 0& 0& 0& \cdots & \cdots & *\end{array}\right)$.

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

12. Finde orthogonales $Z_1\in {ℝ}^{3\times 3}$ mit $B_{22}\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(Z_1,I_{m-3})=\left(\begin{array}{cccccc}*& *& *& \cdots & \cdots & *\\ *& *& *& \cdots & \cdots & *\\ 0& 0& *& \cdots & \cdots & *\\ 0& 0& 0& \ddots & & ⋮\\ ⋮& ⋮& ⋮& & \ddots & ⋮\\ 0& 0& 0& \cdots & \cdots & *\end{array}\right)$. Finde dann orthogonales ${Z_1}^{\prime }\in {ℝ}^{2\times 2}$, so dass

$B_{22}\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(Z_1,I_{m-3})\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}({Z_1}^{\prime },I_{m-2})=\left(\begin{array}{cccccc}*& *& *& \cdots & \cdots & *\\ 0& *& *& \cdots & \cdots & *\\ 0& 0& *& \cdots & \cdots & *\\ 0& 0& 0& \ddots & & ⋮\\ ⋮& ⋮& ⋮& & \ddots & ⋮\\ 0& 0& 0& \cdots & \cdots & *\end{array}\right)$.

Für $A_{22}$ ergibt sich nun: ${\tilde{A}}_{22}:=A_{22}\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(Z_1,I_{m-3})\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(Z{\text{'}}_1,I_{m-2})=\left(\begin{array}{ccccccc}*& *& *& \cdots & \cdots & \cdots & *\\ *& *& *& \cdots & \cdots & \cdots & *\\ *& *& *& \cdots & \cdots & \cdots & *\\ *& *& *& \cdots & \cdots & \cdots & *\\ 0& 0& 0& \ddots & & & ⋮\\ ⋮& ⋮& ⋮& & \ddots & & ⋮\\ 0& 0& 0& \cdots & 0& *& *\end{array}\right)$. D.h., die Hessenbergstruktur von $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 $Q^{\prime }{\text{'}}_1\in {ℝ}^{3\times 3}$, ${Q_1}^{\prime }\in {ℝ}^{3\times 3}$ mit

$\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(1,Q{\text{'}}_1,I_{m-4}{)}^T\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_2,Q^{\prime }{\text{'}}_1,I_{m-5}{)}^T{\tilde{A}}_{22}=\left(\begin{array}{ccccccc}*& *& *& \cdots & \cdots & \cdots & *\\ *& *& *& \cdots & \cdots & \cdots & *\\ 0& *& *& \ddots & \cdots & \cdots & *\\ 0& 0& *& \ddots & \cdots & \cdots & *\\ 0& 0& 0& *& & ⋮\\ ⋮& ⋮& ⋮& & \ddots & & ⋮\\ 0& 0& 0& \cdots & 0& *& *\end{array}\right)$. Es werden also nacheinander die Vektoren $\left(\begin{array}{c}{a}_{21}\\ a_{31}\\ a_{41}\end{array}\right)$ und $\left(\begin{array}{c}{a}_{32}\\ a_{42}\\ a_{52}\end{array}\right)$auf $\left(\begin{array}{c}*\\ 0\\ 0\end{array}\right)$transformiert.

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

$\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(1,Q{\text{'}}_1,I_{m-4}{)}^T\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_2,Q^{\prime }{\text{'}}_1,I_{m-5}{)}^T{\tilde{B}}_{22}=\left(\begin{array}{ccccccc}*& *& *& \cdots & \cdots & \cdots & *\\ 0& *& *& \cdots & \cdots & \cdots & *\\ 0& *& *& \ddots & \cdots & \cdots & *\\ 0& *& *& *& \ddots & \cdots & *\\ 0& 0& 0& 0& *& ⋮\\ ⋮& ⋮& ⋮& & \ddots & & ⋮\\ 0& 0& 0& \cdots & 0& 0& *\end{array}\right)$, d. h. ein Buckel ist jetzt eine Position tiefer entlang der Diagonalen entstanden.

14. Man wiederhole die Schritte 11–13 so lange, bis $A_{22}$ wieder in oberer Hessenberg- und $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 $B_{22}$ an der Diagonalposition j mit Transformationen von links führt zu einem Buckel an der entsprechenden Position in $A_{22}$. Wird dieser Buckel mit Transformationen von rechts eliminiert, führt das zu einem Buckel in $B_{22}$ an der Diagonalposition j+1 usw.

15. Nach $m-2$ Schritten wird das Ziel erreicht und es ergibt sich $Q_{22}^T=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(Q_1,I_{m-3}{)}^T\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(1,Q{\text{'}}_1,I_{m-4}{)}^T\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_2,{Q_1}^{″},I_{m-5}{)}^T\cdots \mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_{m-3},Q_{m-2}{)}^T$. Analog erhält man

$Z_{22}=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(Z_1,I_{m-3})\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}({Z_1}^{\prime },I_{m-2})\cdots \mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_{m-2},Z_{m-2})\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_{m-2},Z_{m^{\prime }-2})$.

Falls $(A,B)$-invarianten Unterräume benötigt werden, ist es notwendig die Matrizen $Q$ und $Z$ aufzudatieren: $Q:=Q\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_p,Q_{22},I_q)$, $Z:=Z\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(I_p,Z_{22},I_q)$

16. end while

Bestimmung der Eigenwerte

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

Literatur

• 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.