Numerische lineare Algebra

Die numerische lineare Algebra ist ein zentrales Teilgebiet der numerischen Mathematik. Sie beschäftigt sich mit der Entwicklung und der Analyse von Rechenmethoden für Problemstellungen der linearen Algebra, insbesondere der Lösung von linearen Gleichungssystemen und Eigenwertproblemen. Solche Probleme treten häufig in den Ingenieurwissenschaften, der Physik oder der Ökonometrie auf.

Einführung in die Problemstellungen

Ein – auch historisch gesehen – zentraler Anfangspunkt der elementaren linearen Algebra sind lineare Gleichungssysteme. Wir betrachten $n$ Gleichungen der Gestalt

$a_{i1}x_{1}+\ldots +a_{in}x_{n}=b_{i}$

für $n$ Unbekannte $x_{1},\dotsc ,x_{n}$. Die Koeffizienten $a_{ij}$ und $b_{i}$ sind gegebene Zahlen; die gesuchten Werte für $x_{1},\dotsc ,x_{n}$ sollen so bestimmt werden, dass alle Gleichungen erfüllt werden. Die Koeffizienten lassen sich zu einer Matrix $A=(a_{ij})$ zusammenfassen; die Zahlen $b_{i}$ und die Unbekannten $x_{j}$ bilden Spaltenvektoren $\mathbf {b} =(b_{i})$ und $\mathbf {x} =(x_{j})$. Auf diese Weise ergibt sich die Matrix-Vektor-Darstellung

$A\cdot \mathbf {x} =\mathbf {b}$

eines linearen Gleichungssystems: Gesucht ist ein Vektor $\mathbf {x}$, der bei der Matrix-Vektor-Multiplikation mit der gegebenen Matrix $A$ den gegebenen Vektor $\mathbf {b}$ ergibt. Als Teilgebiet der Numerik betrachtet auch die numerische lineare Algebra nur sogenannte korrekt gestellte Probleme, also insbesondere nur solche Probleme, die eine Lösung besitzen und bei denen die Lösung eindeutig bestimmt ist. Insbesondere wird im Folgenden stets angenommen, dass die Matrix $A$ regulär ist, also ein Inverse $A^{-1}$ besitzt. Dann gibt es für jede rechte Seite $\mathbf {b}$ eine eindeutig bestimmte Lösung $\mathbf {x}$ des linearen Gleichungssystems, die formal als $\mathbf {x} =A^{-1}\mathbf {b}$ angegeben werden kann.

Viele wichtige Anwendungen führen allerdings auf lineare Gleichungssysteme mit mehr Gleichungen als Unbekannten. In der Matrix-Vektor-Darstellung $A\cdot \mathbf {x} =\mathbf {b}$ hat in diesem Fall die Matrix $A$ mehr Zeilen als Spalten. Solche überbestimmten Systeme haben im Allgemeinen keine Lösung. Man behilft sich deshalb damit, den Vektor $\mathbf {x}$ so zu wählen, dass die Differenz $\mathbf {r} =A\cdot \mathbf {x} -\mathbf {b}$ in einem noch festzulegenden Sinn „möglichst klein“ wird. Beim mit Abstand wichtigsten Fall, dem sogenannten linearen Ausgleichsproblem, wird dazu die Methode der kleinsten Quadrate verwendet: Hierbei wird $\mathbf {x}$ so gewählt, dass die Quadratsumme $r_{1}^{2}+\ldots +r_{m}^{2}$ minimal wird, wobei $r_{1},\dotsc ,r_{m}$ die Komponenten des Differenzvektors $\mathbf {r}$ bezeichnen. Mithilfe der euklidischen Norm lässt sich das auch so schreiben: Man wähle $\mathbf {x}$ so, dass $\|A\cdot \mathbf {x} -\mathbf {b} \|_{2}^{2}$ minimal wird.

Neben den linearen Gleichungen sind die Eigenwertprobleme ein weiteres zentrales Thema der linearen Algebra. Gegeben ist hierbei eine Matrix $A$ mit $n$ Zeilen und $n$ Spalten; gesucht sind Zahlen $\lambda$ und Vektoren $\mathbf {x} \neq \mathbf {0}$, sodass die Gleichung

$A\cdot \mathbf {x} =\lambda \mathbf {x}$

erfüllt ist. Man nennt dann $\mathbf {x}$ einen Eigenvektor von $A$ zum Eigenwert $\lambda$. Das Problem alle Eigenwerte und Eigenvektoren einer Matrix zu bestimmen ist gleichbedeutend damit sie zu diagonalisieren. Das bedeutet: Man finde eine reguläre Matrix $S$ und eine Diagonalmatrix $D$ mit $S^{-1}\cdot A\cdot S=D$. Die Diagonaleinträge von $D$ sind dann die Eigenwerte von $A$ und die Spalten von $S$ die zugehörigen Eigenvektoren.

Grundprinzipien

Fehleranalyse: Vektor- und Matrixnormen

Als Maße für die „Größe“ eines Vektors $\mathbf {x} =(x_{i})\in \mathbb {R} ^{n}$ werden in der Mathematik unterschiedliche Vektornormen verwendet. Am bekanntesten und verbreitetsten ist die euklidische Norm

$\|\mathbf {x} \|_{2}={\sqrt {\sum _{i=1}^{n}x_{i}^{2}}}$,

also die Wurzel aus der Summe der Quadrate aller Vektorkomponenten. Bei der bekannten geometrischen Veranschaulichung von Vektoren als Pfeile im zwei- oder dreidimensionalen Raum entspricht dies gerade der Pfeillänge. Je nach untersuchter Fragestellung können jedoch auch andere Vektornormen wie etwa die Maximumsnorm $\|\mathbf {x} \|_{\infty }$ oder die 1-Norm $\|\mathbf {x} \|_{1}$ geeigneter sein.

Sind $\mathbf {x} ,{\tilde {\mathbf {x} }}\in \mathbb {R} ^{n}$ Vektoren, wobei ${\tilde {\mathbf {x} }}$ als eine Näherung für $\mathbf {x}$ aufgefasst werden soll, so lässt sich mithilfe einer Vektornorm $\|\cdot \|$ die Genauigkeit dieser Näherung quantifizieren. Die Norm des Differenzvektors

$\|{\tilde {\mathbf {x} }}-\mathbf {x} \|$

wird als (normweiser) absoluter Fehler bezeichnet. Betrachtet man den absoluten Fehler im Verhältnis zur Norm des „exakten“ Vektors $\mathbf {x} \neq \mathbf {0}$ erhält man den (normweisen) relativen Fehler

${\frac {\|{\tilde {\mathbf {x} }}-\mathbf {x} \|}{\|\mathbf {x} \|}}$.

Da der relative Fehler nicht durch die Skalierung von $\mathbf {x}$ und ${\tilde {\mathbf {x} }}$ beeinflusst wird, ist dieser das Standardmaß für den Unterschied der beiden Vektoren und wird oft auch vereinfacht nur als „Fehler“ bezeichnet.

Auch die „Größe“ von Matrizen wird mit Normen gemessen, den Matrixnormen. Für die Wahl einer Matrixnorm $\|A\|$ ist es wesentlich, dass sie zur verwendeten Vektornorm „passt“, insbesondere soll die Ungleichung $\|A\mathbf {x} \|\leq \|A\|\|\mathbf {x} \|$ für alle $\mathbf {x}$ erfüllt sein. Definiert man $\|A\|$ für eine gegebene Vektornorm als die kleinste Zahl $L$, sodass $\|A\mathbf {x} \|\leq L\|\mathbf {x} \|$ für alle $\mathbf {x}$ gilt, dann erhält man die sogenannte natürliche Matrixnorm. Für jede Vektornorm gibt es also eine davon induzierte natürliche Matrixnorm: Für die euklidische Norm ist das die Spektralnorm $\|A\|_{2}$, für die Maximumsnorm ist es die Zeilensummennorm $\|A\|_{\infty }$ und für die 1-Norm die Spaltensummennorm $\|A\|_{1}$. Analog zu Vektoren kann mithilfe einer Matrixnorm der relative Fehler

${\frac {\|{\tilde {A}}-A\|}{\|A\|}}$

bei einer Näherung einer Matrix $A$ durch eine Matrix ${\tilde {A}}$ quantifiziert werden.

Kondition und Stabilität

Bei Problemen aus der Praxis sind gegebene Größen meist mit Fehlern behaftet, den Datenfehlern. Zum Beispiel kann bei einem linearen Gleichungssystem $A\mathbf {x} =\mathbf {b}$ die gegebene rechte Seite $\mathbf {b}$ aus einer Messung stammen und daher eine Messabweichung aufweisen. Aber auch bei theoretisch beliebig genau bekannten Größen lassen sich Rundungsfehler bei ihrer Darstellung im Computer als Gleitkommazahlen nicht vermeiden. Es muss also davon ausgegangen werden, dass anstelle des exakten Systems $A\mathbf {x} =\mathbf {b}$ in Wirklichkeit ein System $A{\tilde {\mathbf {x} }}={\tilde {\mathbf {b} }}$ mit einer gestörten rechten Seite ${\tilde {\mathbf {b} }}$ und dementsprechend einer „falschen“ Lösung ${\tilde {\mathbf {x} }}$ vorliegt. Die grundlegende Frage ist nun, wie stark sich Störungen der gegebenen Größen auf Störungen der gesuchten Größen auswirken. Wenn der relative Fehler der Lösung nicht wesentlich größer ist als die relativen Fehler der Eingangsgrößen spricht man von einem gut konditionierten, anderenfalls von einem schlecht konditionierten Problem. Für das Beispiel linearer Gleichungssysteme lässt sich hierzu die Abschätzung

${\frac {\|{\tilde {\mathbf {x} }}-\mathbf {x} \|}{\|\mathbf {x} \|}}\leq \|A\|\|A^{-1}\|\cdot {\frac {\|{\tilde {\mathbf {b} }}-\mathbf {b} \|}{\|\mathbf {b} \|}}$

beweisen. Das Problem ist also gut konditioniert, wenn $\|A\|\|A^{-1}\|$, das Produkt der Norm der Koeffizientenmatrix und der Norm ihrer Inversen, klein ist. Diese wichtige Kenngröße heißt Konditionszahl der Matrix $A$ und wird mit $\kappa (A)$ bezeichnet. In realen Problemen wird meist nicht nur, wie hier dargestellt, die rechte Seite $\mathbf {b}$ fehlerbehaftet sein, sondern auch die Matrix $A$. Dann gilt eine ähnliche, kompliziertere Abschätzung, in der aber ebenfalls $\kappa (A)$ die wesentliche Kennzahl zur Bestimmung der Kondition des Problems bei kleinen Datenfehlern ist.

Während die Kondition eine Eigenschaft des zu lösenden Problems ist, ist Stabilität eine Eigenschaft des dafür verwendeten Verfahrens. Ein numerischer Algorithmus liefert – auch bei exakt gedachten Eingangsdaten – im Allgemeinen nicht die exakte Lösung des Problems. Zum Beispiel muss ein iteratives Verfahren, das eine wahre Lösung schrittweise immer genauer annähert, nach endlich vielen Schritten mit der bis dahin erreichten Näherungslösung abbrechen. Aber auch bei direkten Verfahren, die theoretisch in endlich vielen Rechenschritten die exakte Lösung ergeben, kommt es bei der Umsetzung auf dem Computer bei jeder Rechenoperation zu Rundungsfehlern. In der numerischen Mathematik werden zwei unterschiedliche Stabilitätsbegriffe verwendet, die Vorwärtsstabilität und Rückwärtsstabilität. Sei dazu allgemein $u$ eine Eingabegröße eines Problems und $v=f(u)$ seine exakte Lösung, aufgefasst als Wert einer Funktion $f$ angewendet auf $u$. Auch wenn man die Eingabegröße als exakt vorgegeben betrachtet, wird die Berechnung mit einem Algorithmus ein anderes, „falsches“ Ergebnis ${\tilde {v}}=\operatorname {alg} (u)$ liefern, aufgefasst als Wert einer anderen, „falschen“ Funktion $\operatorname {alg}$ ebenfalls angewendet auf $u$. Ein Algorithmus heißt vorwärtsstabil, wenn sich ${\tilde {v}}$ nicht wesentlich stärker von $v$ unterscheidet, als es aufgrund der Fehler in der Eingangsgröße $u$ und der Kondition des Problems sowieso zu erwarten wäre. Mit einer formalen Definition dieses Begriffs erhält man zwar ein naheliegendes und relativ anschauliches Maß für die Stabilität, aber bei komplizierten Algorithmen ist es oft schwierig, ihre Vorwärtsstabilität zu untersuchen. Daher wird im Allgemeinen zunächst eine sogenannte Rückwärtsanalyse betrachtet: Dazu wird ein ${\tilde {u}}$ bestimmt mit $\operatorname {alg} (u)=f({\tilde {u}})$, das heißt: Der durch das Verfahren berechnete „falsche“ Wert wird aufgefasst als „richtiger“ Wert, der aber mit einem anderen Wert der Eingabegröße berechnet wurde. Ein Algorithmus heißt rückwärtsstabil, wenn sich ${\tilde {u}}$ nicht wesentlich stärker von $u$ unterscheidet, als es aufgrund der Fehler in dieser Eingangsgröße sowieso zu erwarten wäre. Es lässt sich beweisen, dass ein rückwärtsstabiler Algorithmus auch vorwärtsstabil ist.

Orthogonalität und orthogonale Matrizen

Wie die lineare Algebra zeigt, besteht ein enger Zusammenhang zwischen Matrizen und Basen des Vektorraums $\mathbb {R} ^{n}$. Sind $n$ linear unabhängige Vektoren $\mathbf {b} _{1},\dotsc ,\mathbf {b} _{n}$ im $\mathbb {R} ^{n}$ geben, so bilden diese eine Basis. Eine wichtigen Spezialfall bilden die Orthonormalbasen. Hierbei sind die Basisvektoren paarweise orthogonal zueinander („stehen senkrecht aufeinander“) und sind zudem alle auf euklidische Länge 1 normiert, so wie die Standardbasis $(\mathbf {e} _{1},\mathbf {e} _{2},\mathbf {e} _{3})$ im dreidimensionalen Raum. Fasst man die Basisvektoren spaltenweise zu einer Matrix

$(\mathbf {b} _{1}|\cdots |\mathbf {b} _{n})$

zusammen, so erhält man im Fall einer Orthonormalbasis eine sogenannte orthogonale Matrix.

Orthonormalbasen und orthogonale Matrizen besitzen zahlreiche bemerkenswerte Eigenschaften, auf denen die wichtigsten Verfahren der modernen numerischen linearen Algebra basieren. Die Tatsache, dass bei einer orthogonalen Matrix $Q$ die Spalten eine Orthonormalbasis bilden, lässt sich in Matrixschreibweise durch die Gleichung $Q^{T}Q=I$ ausdrücken, wobei $Q^{T}$ die transponierte Matrix und $I$ die Einheitsmatrix bezeichnen. Das zeigt wiederum, dass ein orthogonale Matrix regulär ist und ihre Inverse gleich ihrer Transponierten ist: $Q^{-1}=Q^{T}$. Die Lösung eines linearen Gleichungssystems $Q\mathbf {x} =\mathbf {b}$ lässt sich daher sehr einfach bastimmen, es gilt $\mathbf {x} =Q^{T}\mathbf {b}$. Ein andere grundlegende Eigenschaft ist es, dass eine Multiplikation eines Vektors seine euklidische Norm unverändert lässt

$\|Q\mathbf {x} \|_{2}=\|\mathbf {x} \|_{2}$.

Damit folgt für die Spektralnorm $\|Q\|_{2}=1$ und für die Konditionszahl ebenfalls

$\kappa (Q)=\|Q\|_{2}\|Q^{-1}\|_{2}=1$,

denn $Q^{-1}=Q^{T}$ ist ebenfalls eine orthogonale Matrix. Multiplikationen mit orthogonalen Matrizen bewirken also keine Vergrößerung des relativen Fehlers.

Orthogonale Matrizen spielen auch eine wichtige Rolle in der Theorie und der numerischen Behandlung von Eigenwertproblemen. Nach der einfachsten Version des Spektralsatzes lassen sich symmetrische Matrizen orthogonal diagonalisieren. Damit ist gemeint: Zu einer Matrix $A$, für die $A^{T}=A$ gilt, existiert eine orthogonale Matrix $Q$ und eine Diagonalmatrix $D$ mit

$Q^{T}AQ=D$.

Auf der Diagonale von $D$ stehen die Eigenwerte von $A$ und die Spalten von $Q$ bilden eine Orthonormalbasis aus Eigenvektoren. Mit der sogenannten schurschen Normalform existiert eine Verallgemeinerung dieser orthogonalen Transformation für nicht symmetrische Matrizen.

Literatur

• Wolfgang Bunse, Angelika Bunse-Gerstner: Numerische lineare Algebra. B. G. Teubner, Stuttgart 1985, ISBN 978-3-519-02067-7. 
• Biswa Nath Datta: Numerical Linear Algebra and Applications. 2. Auflage. SIAM, Philadelphia 2010, ISBN 978-0-89871-685-6. 
• James W. Demmel: Applied Numerical Linear Algebra. SIAM, Philadelphia 1997, ISBN 978-0-89871-389-3. 
• Gene H. Golub, Charles F. Van Loan: Matrix Computations. 3. Auflage. The Johns Hopkins University Press, Baltimore 1996, ISBN 0-8018-5413-X. 
• Nicholas J. Higham: Accuracy and Stability of Numerical Algorithms. 2. Auflage. SIAM, Philadelphia 2002, ISBN 0-89871-521-0. 
• Andreas Meister: Numerik linearer Gleichungssysteme. 5. Auflage. Springer Spektrum, Wiesbaden 2015, ISBN 978-3-658-07199-8. 
• Granville Sewell: Computational Methods of Linear Algebra. 3. Auflage. World Scientific, Singapur 2014, ISBN 978-981-4603-85-0. 
• Lloyd N. Trefethen, David Bau, III: Numerical Linear Algebra. SIAM, Philadelphia 1997, ISBN 978-0-89871-361-9. 

Weblinks

• Numerische Lineare Algebra. Institut Computational Mathematics der TU Braunschweig, abgerufen am 2. März 2016 (Überblick über das Fachgebiet). 
• Peter Spellucci: Numerische Lineare Algebra. (PDF) 16. Juli 2009, abgerufen am 2. März 2016 (Vorlesungsskript der TU Darmstadt). 
• Frank Wübbeling: Numerische Lineare Algebra im WS 2012/13. (PDF) 29. März 2013, abgerufen am 2. März 2016 (Vorlesungsskript der Universität Münster).