In der numerischen Mathematik ist das Gauß-Seidel-Verfahren oder Einzelschrittverfahren (nach Carl Friedrich Gauß und Ludwig Seidel) ein Algorithmus zur näherungsweisen Lösung von linearen Gleichungssystemen. Es ist, wie das Jacobi-Verfahren und das SOR-Verfahren, ein spezielles Splitting-Verfahren. Das Verfahren wurde zuerst von Gauß entwickelt, aber nicht veröffentlicht, sondern nur in einem Brief im Jahr 1823 an Gerling erwähnt.[1] Erst 1874 wurde es, bevor seine Anwendung durch Gauß bekannt war, von Seidel veröffentlicht.
Entwickelt wurde das Verfahren, da das Gaußsche Eliminationsverfahren, ein exakter Löser, bei Handrechnung für Rechenfehler sehr anfällig ist. Eine iterative Vorgehensweise hat diesen Nachteil nicht.
Beschreibung des Verfahrens
Gegeben ist ein lineares Gleichungssystem in
Variablen mit den
Gleichungen
![{\displaystyle {\begin{matrix}a_{11}\cdot x_{1}+\dotsb +a_{1n}\cdot x_{n}&=&b_{1}\\a_{21}\cdot x_{1}+\dotsb +a_{2n}\cdot x_{n}&=&b_{2}\\&\vdots &\\a_{n1}\cdot x_{1}+\dotsb +a_{nn}\cdot x_{n}&=&b_{n}\\\end{matrix}}.}](https://wikimedia.org/api/rest_v1/media/math/render/svg/48290aae24e9b4d1ea6187aff49594b7ca562dc7)
Um dieses zu lösen, wird ein Iterationsverfahren durchgeführt. Gegeben sei ein Näherungsvektor
an die exakte Lösung. Nun wird die
-te Gleichung nach der
-ten Variablen
aufgelöst, wobei die vorher berechneten Werte des aktuellen Iterationsschritts mit verwendet werden, im Gegensatz zum Jacobi-Verfahren, bei dem nur die Werte des letzten Iterationsschrittes verwendet werden. Das heißt für den
-ten Iterationsschritt:
![{\displaystyle x_{k}^{(m+1)}:={\frac {1}{a_{kk}}}\left(b_{k}-\sum _{i=1}^{k-1}a_{ki}\cdot x_{i}^{(m+1)}-\sum _{i=k+1}^{n}a_{ki}\cdot x_{i}^{(m)}\right),\,k=1,\dotsc ,n}](https://wikimedia.org/api/rest_v1/media/math/render/svg/5ad19b90158178cbb957c435f49601e640fbd0f5)
Das Ergebnis der Rechnung ist ein neuer Näherungsvektor
für den gesuchten Lösungsvektor
. Wiederholt man diesen Vorgang, gewinnt man eine Folge von Werten, die sich dem Lösungsvektor im Falle der Konvergenz immer mehr annähern:
![{\displaystyle x^{(0)},x^{(1)},x^{(2)},\dotsc \rightarrow x.}](https://wikimedia.org/api/rest_v1/media/math/render/svg/c3d4007551bb4c34c36d9bc84806adfdd5d621dd)
Das Gauß-Seidel-Verfahren ist inhärent sequentiell, das heißt bevor eine Gleichung aufgelöst werden kann, müssen die Ergebnisse der vorherigen Gleichungen vorliegen. Es ist damit wie das Vorwärts- und Rückwärtseinsetzen mit einer LR-Zerlegung nur eingeschränkt zur Nutzung auf Parallelrechnern geeignet, insbesondere nur für dünnbesetzte Matrizen.
Als Algorithmusskizze mit Abbruchbedingung ergibt sich:
- wiederhole
![{\displaystyle {\text{fehler}}:=0}](https://wikimedia.org/api/rest_v1/media/math/render/svg/d90a478741296c04763d863f0e7e03690cb96386)
- für
bis
![{\displaystyle x_{k}^{(m+1)}:={\frac {1}{a_{kk}}}\left(b_{k}-\sum _{i=1}^{k-1}a_{ki}\cdot x_{i}^{(m+1)}-\sum _{i=k+1}^{n}a_{ki}\cdot x_{i}^{(m)}\right)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/1c241536b74e77de5979a8bb05ce058b5ebded62)
![{\displaystyle {\text{fehler}}:=\max({\text{fehler}},|x_{k}^{(m+1)}-x_{k}^{(m)}|)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/482a783c590f221ca77886b958119fbcd1f6a74f)
- nächstes
![{\displaystyle k}](https://wikimedia.org/api/rest_v1/media/math/render/svg/c3c9a2c7b599b37105512c5d570edc034056dd40)
![{\displaystyle m:=m+1}](https://wikimedia.org/api/rest_v1/media/math/render/svg/70912e69169cb8f7b8c97960c6a80486474de8c9)
- bis
![{\displaystyle {\text{fehler}}<{\text{fehlerschranke}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/36e8dbb7fe143e3cc95d13fb25efacaab5905903)
Dabei wurden die Erstbelegung des Variablenvektors, die willkürlich gewählt werden kann, und eine Fehlerschranke als Eingangsgrößen des Algorithmus angenommen, die Näherungslösung ist die vektorielle Rückgabegröße. Die Fehlerschranke misst hier, welche Größe die letzte Änderung des Variablenvektors hatte. Als Bedingung für die Durchführbarkeit des Algorithmus lässt sich festhalten, dass die Hauptdiagonalelemente
von Null verschieden sein müssen.
Bei dünnbesetzten Matrizen reduziert sich der Aufwand des Verfahrens pro Iteration deutlich.
Beschreibung des Verfahrens in Matrixschreibweise
Die Matrix
des linearen Gleichungssystems
wird hierzu in eine Diagonalmatrix
, eine strikte obere Dreiecksmatrix
und eine strikte untere Dreiecksmatrix
zerlegt, so dass gilt:
.
In jedem Iterationsschritt gilt dann
. Nach Umstellen ergibt sich formal
und daraus
.
Man legt dann einen Startvektor
fest und setzt ihn in die Iterationsvorschrift ein:
.
Beispiel
Es ist das lineare Gleichungssystem
mit
und ![{\displaystyle b={\begin{pmatrix}11\\13\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/cd9558f9ba40c36cd8bdff92439738223364458c)
gegeben. Dazu wird die Gleichung
in der Form
verwendet, wobei
und
ist. Dafür muss die Matrix
als Summe einer unteren Dreiecksmatrix
und einer oberen Dreiecksmatrix
dargestellt werden:
![{\displaystyle A=L_{*}+U={\begin{pmatrix}16&0\\7&-11\\\end{pmatrix}}+{\begin{pmatrix}0&3\\0&0\\\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/84db1384db3227f99af44de3f76a3c330979df74)
Die inverse Matrix von
ist
![{\displaystyle L_{*}^{-1}={\begin{pmatrix}16&0\\7&-11\end{pmatrix}}^{-1}\approx {\begin{pmatrix}0{,}0625&0{,}0000\\0{,}0398&-0{,}0909\\\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/0b990d31c991071abe724bca87afdeef722372b6)
Daraus ergibt sich
![{\displaystyle T=-L_{*}^{-1}U=-{\begin{pmatrix}0{,}0625&0{,}0000\\0{,}0398&-0{,}0909\end{pmatrix}}\times {\begin{pmatrix}0&3\\0&0\end{pmatrix}}={\begin{pmatrix}0{,}000&-0{,}1875\\0{,}000&-0{,}1194\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/cee2de5d5e8f17a479b8087c382528d2fba90924)
![{\displaystyle C=L_{*}^{-1}b={\begin{pmatrix}0{,}0625&0{,}0000\\0{,}0398&-0{,}0909\end{pmatrix}}\times {\begin{pmatrix}11\\13\end{pmatrix}}={\begin{pmatrix}0{,}6875\\-0{,}7439\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/0b56b7843c354674a85b804307e965f659d66288)
Daraus können die Vektoren
mithilfe der Iterationsvorschrift
iterativ berechnet werden.
Als erstes wird der Startvektor gesetzt, z. B.
![{\displaystyle x^{(0)}={\begin{pmatrix}1{,}0\\1{,}0\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f24d70a9dcf629730d75fc1149698d522102092c)
Daraus ergibt sich schrittweise
![{\displaystyle x^{(1)}={\begin{pmatrix}0{,}000&-0{,}1875\\0{,}000&-0{,}1193\end{pmatrix}}\times {\begin{pmatrix}1{,}0\\1{,}0\end{pmatrix}}+{\begin{pmatrix}0{,}6875\\-0{,}7443\end{pmatrix}}={\begin{pmatrix}0{,}5000\\-0{,}8636\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/89ad9557ebe935600d791596e79107f89bf562a7)
![{\displaystyle x^{(2)}={\begin{pmatrix}0{,}000&-0{,}1875\\0{,}000&-0{,}1193\end{pmatrix}}\times {\begin{pmatrix}0{,}5000\\-0{,}8636\end{pmatrix}}+{\begin{pmatrix}0{,}6875\\-0{,}7443\end{pmatrix}}={\begin{pmatrix}0{,}8494\\-0{,}6413\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a4d99c085726fca309ca452bde8d76ca13bd20af)
![{\displaystyle x^{(3)}={\begin{pmatrix}0{,}000&-0{,}1875\\0{,}000&-0{,}1193\end{pmatrix}}\times {\begin{pmatrix}0{,}8494\\-0{,}6413\\\end{pmatrix}}+{\begin{pmatrix}0{,}6875\\-0{,}7443\end{pmatrix}}={\begin{pmatrix}0{,}8077\\-0{,}6678\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/6444eced7f663336cbf3e6654bbe7e4f055e7d7e)
![{\displaystyle x^{(4)}={\begin{pmatrix}0{,}000&-0{,}1875\\0{,}000&-0{,}1193\end{pmatrix}}\times {\begin{pmatrix}0{,}8077\\-0{,}6678\end{pmatrix}}+{\begin{pmatrix}0{,}6875\\-0{,}7443\end{pmatrix}}={\begin{pmatrix}0{,}8127\\-0{,}6646\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/991253bd081eea61d91e8b05ec0103158f1f1727)
![{\displaystyle x^{(5)}={\begin{pmatrix}0{,}000&-0{,}1875\\0{,}000&-0{,}1193\end{pmatrix}}\times {\begin{pmatrix}0{,}8127\\-0{,}6646\end{pmatrix}}+{\begin{pmatrix}0{,}6875\\-0{,}7443\end{pmatrix}}={\begin{pmatrix}0{,}8121\\-0{,}6650\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/ffbd8567855f50c02a4bef51ff2418922c7f1bf9)
![{\displaystyle x^{(6)}={\begin{pmatrix}0{,}000&-0{,}1875\\0{,}000&-0{,}1193\end{pmatrix}}\times {\begin{pmatrix}0{,}8121\\-0{,}6650\end{pmatrix}}+{\begin{pmatrix}0{,}6875\\-0{,}7443\end{pmatrix}}={\begin{pmatrix}0{,}8122\\-0{,}6650\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/344b60812f69f396049ac556cc7175271eab4279)
![{\displaystyle x^{(7)}={\begin{pmatrix}0{,}000&-0{,}1875\\0{,}000&-0{,}1193\end{pmatrix}}\times {\begin{pmatrix}0{,}8122\\-0{,}6650\end{pmatrix}}+{\begin{pmatrix}0{,}6875\\-0{,}7443\end{pmatrix}}={\begin{pmatrix}0{,}8122\\-0{,}6650\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/d50afbb74cf8342ef92616eeca10f2aefbbac2e9)
Der Algorithmus konvergiert also gegen die Lösung
![{\displaystyle x=A^{-1}b\approx {\begin{pmatrix}0{,}8122\\-0{,}6650\end{pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/b7094928cb86c0c48523e898e339305da74b458a)
Diagonaldominanz und Konvergenz
Das Verfahren konvergiert linear, wenn der Spektralradius der Iterationsmatrix
kleiner 1 ist. In diesem Falle ist der Fixpunktsatz von Banach bzw. der Konvergenzsatz der Neumann-Reihe (auf eine hinreichend große Potenz von
) anwendbar und das Verfahren konvergiert. Im gegenteiligen Fall divergiert das Verfahren, wenn die rechte Seite des Gleichungssystems einen Anteil eines Eigenvektors zu einem Eigenwert mit Betrag größer als 1 beinhaltet. Je geringer der Spektralradius, desto schneller konvergiert das Verfahren.
Die Bestimmung des Spektralradius ist für den praktischen Einsatz meist ungeeignet, weswegen über die hinreichende Bedingung, dass eine Matrixnorm der Verfahrensmatrix
kleiner als 1 sein muss, bequemere Kriterien gefunden werden können. Diese Matrixnorm ist gleichzeitig die Kontraktionskonstante im Sinne des banachschen Fixpunktsatzes.
Im Falle, dass sowohl
als auch
„kleine“ Matrizen bzgl. der gewählten Matrixnorm sind, gibt es die folgende Abschätzung der Matrixnorm für
(siehe Neumann-Reihe für den ersten Faktor):
![{\displaystyle {\begin{aligned}\|T\|&\leq \|(I+D^{-1}L)^{-1}\|\cdot \|D^{-1}U\|\leq {\frac {\|D^{-1}U\|}{1-\|D^{-1}L\|}}\\&=1-{\frac {1-\left(\|D^{-1}L\|+\|D^{-1}U\|\right)}{1-\|D^{-1}L\|}}\end{aligned}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/9f56ec3465c134edb459bf04651f645f5e44cef1)
Der letzte Ausdruck ist für
ebenfalls kleiner als 1. Obwohl die Konvergenzbedingung diejenige des Jacobi-Verfahrens ist, ist die so erhaltene Abschätzung der Kontraktionskonstante
des Gauß-Seidel-Verfahrens immer kleinergleich der entsprechenden Abschätzung der Kontraktionskonstante des Jacobi-Verfahrens.
Das einfachste und gebräuchlichste hinreichende Konvergenzkriterium der Diagonaldominanz ergibt sich für die Supremumsnorm der Vektoren und die Zeilensummennorm als zugehörige induzierte Matrixnorm. Es verlangt die Erfüllung des Zeilensummenkriteriums, also der Ungleichung
für
.
Je größer die kleinste Differenz zwischen rechten und linken Seiten der Ungleichung ist, desto schneller konvergiert das Verfahren. Man kann versuchen, diese Differenz mittels Zeilen- und Spaltenvertauschungen zu vergrößern, d. h. durch Umnummerieren der Zeilen und Spalten. Dabei kann man beispielsweise anstreben, die betragsgrößten Elemente der Matrix auf die Diagonale zu bringen. Die Zeilenvertauschungen müssen natürlich auch auf der rechten Seite des Gleichungssystems vollzogen werden.
Anwendungen
Für moderne Anwendungen wie die Lösung großer dünnbesetzter Gleichungssysteme, die aus der Diskretisierung partieller Differentialgleichungen stammen, ist das Verfahren ungeeignet. Es wird jedoch mit Erfolg als Vorkonditionierer in Krylow-Unterraum-Verfahren oder als Glätter in Mehrgitterverfahren eingesetzt.
Erweiterung auf nichtlineare Gleichungssysteme
Die Idee des Gauß-Seidel-Verfahrens lässt sich auf nichtlineare Gleichungssysteme
mit einer mehrdimensionalen nichtlinearen Funktion
erweitern. Wie im linearen Fall wird im
-ten Schritt die
-te Gleichung bezüglich der
-ten Variablen gelöst, wobei für die weiteren Variablen der bisherige Näherungswert und für die vorherigen Variablen der vorher berechnete neue Wert genommen wird:
- Für k=1, … bis zur Konvergenz
- Für i=1, …, n:
- Löse
nach
auf.
Hierbei ist das Lösen in der Regel als die Anwendung eines weiteren iterativen Verfahrens zur Lösung nichtlinearer Gleichungen zu verstehen. Um dieses Verfahren vom Gauß-Seidel-Verfahren für lineare Gleichungssysteme zu unterscheiden, wird häufig vom Gauß-Seidel-Prozess gesprochen. Die Konvergenz des Prozesses folgt aus dem Banachschen Fixpunktsatz wieder als linear.
Literatur
- A. Meister: Numerik linearer Gleichungssysteme, 2. Auflage, Vieweg 2005, ISBN 3-528-13135-7
- R. Barrett et al.: Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods, 2nd Edition, SIAM Philadelphia, 1994
- W. C. Rheinboldt: Methods for Solving Systems of Nonlinear Equations, 2. Auflage, SIAM, 1998, ISBN 0-89871-415-X
Weblinks
- Eric W. Weisstein et al.: Gauss-Seidel Method. MathWorld
- C code (Beispiel)
Einzelnachweise
- ↑ http://gdz.sub.uni-goettingen.de/en/dms/loader/img/?PPN=PPN23601515X&DMDID=DMDLOG_0112&LOGID=LOG_0112&PHYSID=PHYS_0286