Gauß-Seidel-Verfahren: Konvergenz und Programm

Die Aufgabe

Betrachtet werde das numerische Gauß-Seidel-Verfahren. Mit der Matrixzerlegung $A = L + D + R$ ( $L$ : strikt untere Dreiecksmatrix, $D$ : Diagonale, $R$ : strikt obere Dreiecksmatrix) folgt die Iteration

x_{i+1} = -(D+L)^{-1}\cdot R\cdot x_i + (D+L)^{-1}\cdot b\,.

Hinweis: Beachten Sie, dass die $x_i\in\mathbb{R}^n$ , also Vektoren sind.

a) Recherchieren Sie Konvergenzaussagen des Verfahrens.

b) Programmieren Sie das Verfahren allgemein und testen Sie dieses anhand des Beispiels

\begin{pmatrix} 16 & 3 \\ 7 & -11 \end{pmatrix}\cdot x = \begin{pmatrix} 22 \\ -15 \end{pmatrix}.

c) Lösen Sie das LGS

\begin{bmatrix} 1 & 4 & 3 \\ 2 & 5 & 4 \\ 1 & -3 & -2 \end{bmatrix}\cdot x = \begin{pmatrix} 1 \\ 4 \\ 5 \end{pmatrix}

von Hand und mit dem Programm. Was beobachten Sie bei der Iteration?

Gegeben ist ein Verfahren, das ein Gleichungssystem $Ax = b$ näherungsweise löst. In Teil a) sollen wir nachschlagen, wann es gegen die Lösung läuft. In Teil b) schreiben wir ein Programm und testen es an einem Beispiel mit zwei Unbekannten. In Teil c) rechnen wir ein Beispiel mit drei Unbekannten zweimal: von Hand und mit dem Programm. Dabei sollen wir beobachten, was die Näherungen tun.

Eine Schreibweise ändern wir: In der Aufgabe zählt der Index $i$ die Näherungen. Wir schreiben $x^{(k)}$ für die Näherung nach $k$ Durchläufen. Das Wort „Durchlauf“ meint: einmal alle Komponenten neu berechnen. Die Komponenten eines Vektors $x$ heißen $x_1, x_2, \dots$ . So wird der Hinweis der Aufgabe sichtbar: $x^{(k)}$ ist ein ganzer Vektor.

Merke: Die Musterlösung antwortet in Teil a) kurz. Ob diese Kurzantwort genau genug ist, klären die Schritte 5 und 6.

So gehen wir vor:

Die Zerlegung und die Iterationsformel verstehen (Schritte 1 bis 3).
Die Konvergenzaussagen angeben: Spektralradius, Diagonaldominanz, positive Definitheit (Schritte 4 bis 6).
Das Verfahren in Python programmieren und am Beispiel b) testen (Schritte 7 bis 12).
Das System c) von Hand exakt lösen, dann mit dem Verfahren iterieren und beobachten (Schritte 13 bis 19).

Die Zerlegung A = L + D + R

Die Matrix aus b) wird so zerlegt (der Eintrag 7 steht unter der Diagonale, die 3 darüber):

Zerlegung der Matrix aus Teil b)

\underbrace{\begin{pmatrix} 16 & 3 \\ 7 & -11 \end{pmatrix}}_{A} = \underbrace{\begin{pmatrix} 0 & 0 \\ 7 & 0 \end{pmatrix}}_{L} + \underbrace{\begin{pmatrix} 16 & 0 \\ 0 & -11 \end{pmatrix}}_{D} + \underbrace{\begin{pmatrix} 0 & 3 \\ 0 & 0 \end{pmatrix}}_{R}

Das Verfahren zerlegt die Matrix $A$ in drei Teile, so wie das Bild rechts es zeigt. Für jede Matrix gilt:

$L$ enthält alle Einträge unterhalb der Diagonale, sonst nur Nullen. Das „L“ steht für „lower“, englisch „unten“.
$D$ enthält nur die Einträge auf der Diagonale, sonst Nullen. Das „D“ steht für „Diagonale“.
$R$ enthält alle Einträge oberhalb der Diagonale, sonst Nullen. Das „R“ steht für „rechts oben“.

Addiert man die drei Matrizen Feld für Feld, kommt wieder $A$ heraus, weil jedes Feld von $A$ in genau einer der drei Matrizen steht.

Die Iterationsformel entsteht aus Ax = b

Von Ax = b zur Iteration

\begin{aligned} (L + D + R)\,x &= b \\ (D + L)\,x &= b - R\,x \\ (D + L)\,x^{(k+1)} &= b - R\,x^{(k)} \\ x^{(k+1)} &= -(D+L)^{-1} R\,x^{(k)} + (D+L)^{-1} b \end{aligned}

Kurzschreibweise

x^{(k+1)} = T\,x^{(k)} + c \quad\text{mit}\quad T = -(D+L)^{-1}R,\quad c = (D+L)^{-1}b

Warum hat die Formel aus der Aufgabe diese Gestalt? Wir starten bei dem System $Ax = b$ und setzen $A = L + D + R$ ein: $(L + D + R)\,x = b$ . Jetzt bringen wir $R\,x$ auf die rechte Seite: $(D + L)\,x = b - R\,x$ .

Nun kommt die Idee des Verfahrens. Rechts setzen wir die alte Näherung $x^{(k)}$ ein. Links steht dann die neue Näherung $x^{(k+1)}$ . Um sie allein zu erhalten, multipliziert man von links mit der Inversen $(D+L)^{-1}$ (der Matrix, die $D + L$ rückgängig macht), ergibt sich genau die Formel der Aufgabe.

Das geht nur, wenn $D + L$ eine Inverse hat. $D + L$ hat unterhalb und auf der Diagonale Einträge und darüber nur Nullen. Bei so einer Dreiecksmatrix ist die Determinante das Produkt der Diagonaleinträge. Sie ist also nicht null, wenn kein Diagonaleintrag von $A$ null ist. Das setzen wir ab jetzt voraus.

Aus der Lösung linearer Gleichungssysteme: Eine quadratische Matrix hat genau dann eine Inverse, wenn ihre Determinante nicht null ist.

Die Abkürzungen $T$ und $c$ machen die Formel kürzer: $T$ heißt Iterationsmatrix, $c$ ist der feste Vektor, der nur von $A$ und $b$ abhängt.

Schritt 0Die Aufgabe

Schritt 1Die Zerlegung A = L + D + R

Schritt 2Die Iterationsformel entsteht aus Ax = b

Die Zerlegung A = L + D + R

Abbildung 1Die Zerlegung A = L + D + R

Die Aufgabe

Die Zerlegung A = L + D + R

Die Iterationsformel entsteht aus Ax = b