Unter der iterativen Lösung linearer Gleichungssysteme Ax = b mit \(A\in {{\mathbb{R}}}^{n\times n},b\in {{\mathbb{R}}}^{n}\) versteht man Verfahren, welche die Lösung als Grenzwert einer unendlichen Folge von Näherungslösungen berechnen; im Gegensatz dazu versteht man unter einer direkten Lösung linearer Gleichungssysteme Verfahren, welche (bei rundungsfehlerfreier Rechnung) die exakte Lösung nach endlich vielen Schritten berechnen. Dabei verändern direkte Verfahren typischerweise die gegebene Koeffizientenmatrix und transformieren das ursprüngliche Problem in ein einfacher zu lösendes. Iterative Verfahren berechnen sukzessive Näherungslösungen an die gesuchte Lösung. Ein Iterationsschritt besteht häufig in der Ausführung einer Matrix-VektorMultiplikation. Während direkte Verfahren zumindest theoretisch die exakte Lösung des Problems berechnen, sind bei iterativen Verfahren Fragen der Konvergenz und Konvergenzgeschwindigkeit von Bedeutung.

Die klassischen Iterationsverfahren beruhen auf dem Umschreiben des linearen Gleichungssystems Ax = b in eine Fixpunktgleichung \begin{eqnarray}x=Tx+f\end{eqnarray} mit einer geeigneten Matrix \(T\in {{\mathbb{R}}}^{n\times n}\) und einem Vektor \(f\in {{\mathbb{R}}}^{n}\). Ist die Fixpunktgleichung eindeutig lösbar, so wählt man einen Startvektor \({x}^{(0)}\in {{\mathbb{R}}}^{n}\) und bildet die Folge \begin{eqnarray}{x}^{(k+1)}=T{x}^{(k)}+f\end{eqnarray} für k = 1, 2,…. Sind alle Eigenwerte von T betragsmäßig kleiner als 1, dann konvergiert diese Folge gegen die Lösung x der Gleichung Ax = b. Eine Fixpunktgleichung der Form x = Tx + f erhältman dabei aus Ax = b typischerweise aus einerZerlegung von A in \begin{eqnarray}A=B-C,\end{eqnarray}\(B,C\in {{\mathbb{R}}}^{n\times n}\), B nichtsingulär, denn dann ist (B – C)x = b, bzw. \begin{eqnarray}x={B}^{-1}Cx+{B}^{-1}b.\end{eqnarray} Man erhält also ein System x = Tx + f mit \(T={B}^{-1}C\) und \(f={B}^{-1}b\). Die prominentesten Vertreter dieser Iterationsverfahren sind das Jacobi-Verfahren, das Gauß-Seidel-Verfahren und dasSOR-Verfahren.

Die Konvergenzgeschwindigkeit einer Fixpunktiteration hängt von den Eigenwerten von T ab. Nurwenn alle Eigenwerte von T betragsmäßig kleinerals 1 sind, d. h. wenn \begin{eqnarray}\varrho (T)=\text {max}\{|\lambda |,\lambda\, \text{Eigenwert von}\,T\}\lt 1,\end{eqnarray} tritt Konvergenz ein. Je kleiner \(\varrho (T)\) ist, destoschneller ist die Konvergenz. Falls \(\varrho (T)\) nahe bei1 ist, so ist die Konvergenz recht langsam und mankann versuchen, mittels Relaxation oder poly-nomieller Konvergenzbeschleunigung aus der gege-benen Fixpunktiteration eine schneller konvergie-rende herzuleiten.

Moderne Iterationsverfahren sind häufig Krylow-Raum-basierte Verfahren. Dabei wird, ausgehendvon einem (beliebigen) Startvektor \({x}^{(0)}\), eine Folgevon Näherungsvektoren \({x}^{(k)}\) an die gesuchte Lö-sung x gebildet. \({x}^{(k)}\) wird dazu aus einem verscho-benen Krylow-Raum \begin{eqnarray}\{{x}^{(0)}\}+{{\mathcal{K}}}_{k}(B,{r}^{(0)})\end{eqnarray} mit \begin{eqnarray}{{\mathcal{K}}}_{k}(B,{r}^{(0)})=\{{r}^{(0)},B{r}^{(0)},{B}^{2}{r}^{(0)},\mathrm{\ldots },{B}^{k-1}{r}^{(0)}\}\end{eqnarray} und r⁽⁰⁾ = b − Ax⁽⁰⁾ so gewählt, daß eine Bedingung der Art \begin{eqnarray}b-A{x}^{(k)}\perp { {\mathcal L} }_{k}\end{eqnarray} erfüllt ist. Die verschiedenen Versionen von Krylow-Raum-Verfahren unterscheiden sich in der Wahl der Matrix B und des k-dimensionalen Raums \({ {\mathcal L} }_{k}\).

Das älteste Verfahren dieser Klasse ist das konjugierte Gradientenverfahren, welches lineare Gleichungssysteme Ax = b mit symmetrisch positiv definiten Matrizen \(A\in {{\mathbb{R}}}^{n\times n}\) löst. Die nächste Näherung \({x}^{(k)}\) wird so gewählt, daß \({x}^{(k)}\) den Ausdruck \begin{eqnarray}(x-{x}^{(k)})A(x-{x}^{(k)})\end{eqnarray} über dem verschobenen Krylow-Raum \begin{eqnarray}\{{x}^{(0)}\}+{{\mathcal{K}}}_{k}(A,{r}^{(0)})\end{eqnarray} minimiert. Pro Iterationsschritt wird lediglich eine Matrix-Vektor-Multiplikation benötigt; die Matrix A selbst bleibt unverändert. Das Verfahren ist daher insbesondere für große, sparse Koeffizientenmatrizen A geeignet. Zur Berechnung der k-ten Iterier- ten wird lediglich Information aus dem (k − 1)-ten Schritt benötigt. Informationen aus den Schritten 1, 2,…,k−2 muß nicht gespeichert werden. Theoretisch ist bei exakter Rechnung spätestens nach n Schritten die gesuchte Lösung \({x}^{(n)}\) berechnet. Doch aufgrund von Rundungsfehlern ist dies in der Praxis i. a. nicht der Fall. Man setzt das Verfahren einfach solange fort, bis \({x}^{(n)}\) eine genügend gute Näherung an die gesuchte Lösung x ist, d. h. bis \(b-A{x}^{(k)}\) klein genug ist. Das Konvergenzverhalten wird wesentlich durch die Konditionszahl (Kondition) der Matrix A bestimmt. Hat A eine kleine Kondition, so wird das konjugierte Gradientenverfahren in der Regel schnell konvergieren. Hat A aber eine große Kondition, tritt Konvergenz häufig nur sehr langsam ein. Um dieses Problem zu beheben, kann man versuchen, das Gleichungssystem Ax = b in ein äquivalentes Gleichungssystem \(\mathop{\tilde A}\mathop{\tilde x}=\mathop{\tilde b}\) zu überführen, für welches das Iterationsverfahren ein besseres Konvergenzverhalten hat. Diese Technik wird Vorkonditionierung genannt und kann die Konvergenz des konjugierten Gradientenverfahrens erheblich verbessern.

Ein wesentlicher Punkt bei der Herleitung des konjugierten Gradientenverfahrens ist die Restriktion auf symmetrisch positiv definite Matrizen. Um das konjugierte Gradientenverfahren für beliebige (unsymmetrische oder symmetrische, nicht positiv definite) Matrizen zu verallgemeinern, wurden zahlreiche Ansätze vorgeschlagen. So ist z. B. das lineare Gleichungssystem Ax = b mit beliebiger (n × n)-Matrix A äquivalent zu dem linearen Gleichungssystem \begin{eqnarray}{A}^{T}Ax={A}^{T}b,\end{eqnarray} dessen Koeffizientenmatrix A^TA symmetrisch positiv definit ist. Wendet man nun das konjugierte Gradientenverfahren auf A^TAx = A^Tb an, so erhält man das CGNR-Verfahren (Conjugate Gradient applied to Normal equations minimizing the Residual), welches die Iterierte \({x}^{(k)}\) so berechnet, daß \begin{eqnarray}{(b-A{x}^{(k)})}^{T}(b-A{x}^{(k)})\end{eqnarray} über \begin{eqnarray}\{{x}^{(0)}\}+{{\mathcal{K}}}_{k}({A},{r}^{(0)}).\end{eqnarray} minimiert wird. Bei diesem Verfahren bestimmt die Kondition von A^TA im wesentlichen die Konvergenzrate. Da die Kondition von A^TA das Quadrat der Kondition von A sein kann, wird das CGNR- Verfahren schon für Matrizen A mit nicht allzu großer Kondition nicht gut konvergieren.

Weitere Varianten unterscheiden sich hauptsächlich in der Wahl des Krylow-Raums und der Minimierungsaufgabe. So minimiert das GCR- Verfahren (Generalized Conjugate Residual) \begin{eqnarray}{(b-A{x}^{(k)})}^{T}(b-A{x}^{(k)})\end{eqnarray} über \begin{eqnarray}\{{x}^{(0)}\}+{{\mathcal{K}}}_{k}(A,{r}^{(0)}).\end{eqnarray} Hier hat man nicht das Problem der eventuell großen Kondition von A^TA, dafür benötigt man aber zur Berechnung der k-ten Iterierten Informationen aus allen vorangegangenen Schritten, d. h. mit wachsendem k wird mehr und mehr Speicherplatz benötigt. Darüberhinaus kann das GCR-Verfahren zusammenbrechen, ohne eine Lösung des linearen Gleichungssystems zu berechnen. Dieses Problem kann man durch Wahl einer geeigneten Basis für den Krylow-Raum \({{\mathcal{K}}}_{k}(A,{r}^{(0)})\) beheben.

Dies führt dann auf das GMRES-Verfahren (Generalized Minimal RESidual), welches nicht zusammenbrechen kann und die minimale Anzahl von Matrix-Vektor-Multiplikationen zur Minimierung von \begin{eqnarray}{(b-A{x}^{(k)})}^{T}(b-A{x}^{(k)})\end{eqnarray} über einem gegebenen Krylow-Raum benötigt. Die benötigte orthogonale Basis des \({{\mathcal{K}}}_{k}(A,{r}^{(0)})\) wird mittels des Arnoldi-Verfahrens berechnet.

Zwei weitere, häufig verwendete Varianten beruhen statt auf dem Arnoldi-Verfahren zur Berechnung einer orthogonalen Basis eines Krylow- Raums auf dem unsymmetrischen Lanczos- Verfahren. Das BiCG-Verfahren (BiConjugate Gradient-Verfahren) und das QMR-Verfahren (Quasi Minimal Residual-Verfahren) können allerdings zusammenbrechen, sodaß hier sogenannte look-ahead-Methoden angewendet werden sollten, um diese Verfahren stets durchführen zu können.

In der Literatur existieren zahlreiche weitere Varianten von Krylow-Raum-basierten Verfahren. Jedes dieser Verfahren hat gewisse Vor- und Nachteile, für jedes Verfahren lassen sich Beispiele finden, für welche das jeweilige Verfahren besonders gut oder besonders schlecht geeignet ist. Eine allgemeine Regel, welches Verfahren wann angewendet werden sollte, gibt es (noch) nicht.

[3] Hackbusch, W.: Iterative Lösung großer schwachbesetzer Gleichungssysteme. B.G. Teubner Stuttgart, 1993.

[4] Saad, Y.: Iterative Methods for Sparse Linear Systems. The Pws Series in Computer Science, 1996.

[5] Schwarz, H.R.: Numerische Mathematik. B.G. Teubner Stuttgart, 1993.

[6] Stoer, J. und Bulirsch, R.: Numerische Mathematik I und II. Springer Heidelberg/Berlin, 1994/1991.