[WS] Lineare Gleichungssysteme 2 - direkte Verfahren

Neue Frage »

19.07.2007, 03:38	tigerbine	Auf diesen Beitrag antworten »
[WS] Lineare Gleichungssysteme 2 - direkte Verfahren Gliederung Lösbarkeit von Linearen Gleichungssystemen 1a. Satz von Hurwitz (Sylvester Kriterium) 1b. Digonaldominante Matrizen 1c. Strikt Diagonaldominante Matrizen 1d. Dreiecksmatrizen 1e. Frobeniusmatrizen Lösen einfacher Gleichungssysteme 2a. Diagonalmatrix 2b. untere Dreiecksmatrix 2c. obere Dreiecksmatrix LR-Zerlegung ohne Pivotisierung 3a. Berechnung 3b. Existenzsatz 1 3c. Existenzsatz 2 3d. Eindeutigkeit 3e. Algorithmus zur Lösung des LGS Ax = b LR-Zerlegung mit Pivotisierung 4a. Berechnung 4b. Existenz 4c. Eindeutigkeit - Pivotwahl 4d. Numerische Hinweise weitere Möglichkeiten 5a. Determinantenberechnung 5b. Rangebstimmung 5c. Inversenbestimmung 1 - Gauss-Algorithmus 5d. Inversenbestimmung 2 - Gauss-Jordan Algorithmus 5e. Inversenbestimmung 3 - Austauschauschalgorithmus Direkte LR-Zerlegung 6a. Verfahren von Doolittle und Crout Bandmatrizen 7a. Tridiagonalmatrizen 7b. Allgemein zum Speichern Symmetrische positiv definite Matrizen (SPD) 8a. LDL^T-Zerlegungen 8b.Cholesky-Zerlegung QR-Zerlegungen 9a. Eindeutigkeit der QR-Zerlegung 9b. Lösen regulärer LGS Ax=b 9c. Existenz von QR-Zerlegungen 9d. Gram-Schmidt-Orthogonalisierung 9e. Householder-Spiegelungen 9f. Givens-Rotationen
19.07.2007, 03:40	tigerbine	Auf diesen Beitrag antworten »
1. Lösbarkeit von Linearen Gleichungssystemen Betrachtet werden lineare Gleichungssysteme über reellen endlich-dimensionalen Vektorräumen. $\begin{eqnarray} \text{A ist eine reguläre (n}\times \text{n) Matrix } \Leftrightarrow \text{Das LGS Ax = b besitzt }\forall~ b \in \mathbb R^n \text{ genau eine Lösung } x \in \mathbb R^n \end{eqnarray}$ $\begin{eqnarray} A \in \mathbb R^{m \times n}, \operatorname{Rang}(A) = m \Leftrightarrow \text{ Das LGS } Ax=b \text{ ist universell lösbar.} \end{eqnarray}$ $\begin{eqnarray} b \notin \operatorname{Im}(A) \Leftrightarrow \text{ Das LGS besitzt keine Lösung.} \end{eqnarray}$ In den weiteren Betrachtungen wird meist A als regulär angesehen. Ziel ist es die eindeutig existierende Lösung von Ax=b zu bestimmen.
19.07.2007, 04:33	tigerbine	Auf diesen Beitrag antworten »
1a. Satz von Hurwitz (Sylvester Kriterium) $\begin{eqnarray} \text{Ist }A \in \mathbb R^{n \times n} \text{ symmetrisch, so gilt: } \end{eqnarray}$ ******************************** $\begin{eqnarray} \text{A ist positiv definit } \Leftrightarrow \text{ det(A[k]) }>\text{ 0 für alle k=1,...,n } \end{eqnarray}$ Beweis: Vgl. z.B. Fischer - Lineare Algebra
19.07.2007, 04:37	tigerbine	Auf diesen Beitrag antworten »
1b. Diagonaldominate Matrizen Eine (n x n)-Matrix $\begin{eqnarray} A = (a_{ij}) \end{eqnarray}$ heißt schwach diagonaldominant, falls die Beträge ihrer Diagonalelemente $\begin{eqnarray} a_{ii} \end{eqnarray}$ jeweils größer oder gleich der Summe der Beträge der restlichen jeweiligen Zeileneinträge $\begin{eqnarray} a_{ij} \end{eqnarray}$ sind, d.h. wenn für alle $\begin{eqnarray} i\in\{1,\ldots, n\} \end{eqnarray}$ gilt $\begin{eqnarray} \sum_{j=1 \atop j\ne i}^n\|a_{ij}\|\leq\|a_{ii}\| \end{eqnarray}$
19.07.2007, 04:39	tigerbine	Auf diesen Beitrag antworten »
1c. Strikt Diagonaldominate Matrizen und SPD-Matrizen Eine (n x n)-Matrix $\begin{eqnarray} A = (a_{ij}) \end{eqnarray}$ heißt strikt diagonaldominant, falls die Beträge ihrer Diagonalelemente $\begin{eqnarray} a_{ii} \end{eqnarray}$ jeweils größer sind als die Summe der Beträge der restlichen jeweiligen Zeileneinträge $\begin{eqnarray} a_{ij} \end{eqnarray}$ , d.h. wenn für alle $\begin{eqnarray} i\in\{1,\ldots, n\} \end{eqnarray}$ gilt $\begin{eqnarray} \sum_{j=1 \atop j\ne i}^n\|a_{ij}\|<\|a_{ii}\| \end{eqnarray}$ $\begin{eqnarray} \text{A ist strikt diagonaldominant }\Rightarrow \text{ A ist regulär} \end{eqnarray}$ ***************************************** Beweis:** Es wird hier gezeigt, dass A dann eine LR-Zerlegung (ohne Pivotisierung) besitzt. Daraus folgt direkt die Regularität. $\begin{eqnarray} \text{Sei A ist strikt diagonaldominant und symmetrisch. Dann gilt:} \end{eqnarray}$ ******************************************************** $\begin{eqnarray} \text{Die Diagonalelemente von A sind positiv }\Rightarrow \text{ A ist positiv definit} \end{eqnarray}$ Beweis: Eine symmetrische Matrix ist z.B. dann pos. definit, wenn alle EW nicht negativ sind. Sie ist nach dem Spektralsatz diagonalisierbar, d.h. es gibt eine ONB aus Eigenvektoren. Es gilt dann: $\begin{eqnarray} x^TAx = x^T(Q^TDQ)x =(Qx)^TD(Qx) = y^TDy =\sum_{j=1}^n~d_{jj}\cdot y_j^2 ~> 0 \quad \forall x \neq 0 \end{eqnarray}$ Wie sehen nun die EW der vorliegenden Matrix A aus? Sei $\begin{eqnarray} \lambda \end{eqnarray}$ ein EW von A mit EV x, mit $\begin{eqnarray} \|\|x\|\|_{\infty} =1 \end{eqnarray}$ . So gibt es ein maximales Element $\begin{eqnarray} \|x_i\|=1 \end{eqnarray}$ in x. Dann gilt: $\begin{eqnarray} \lambda\cdot x_i = (Ax)_i = \sum_{j=1}^{n}~a_{ij}\cdot x_j = \sum_{j=1 \atop j \neq i}^n~a_{ij}\cdot x_j + a_{ii} \cdot x_i \end{eqnarray}$ $\begin{eqnarray} \|(\lambda - a_{ii})\|\underbrace{\|x_i\| }_{=1}&&= \|\sum_{j=1\atop j\neq i}^n~a_{ij}\cdot x_j\| \\&&\leq \sum_{j=1\atop j\neq i}^n~\|a_{ij}\| \cdot \underbrace{\|x_j\|}_{\leq 1} \end{eqnarray}$ Also gilt wegen der strikten Diagonaldominanz und der pos. Diagonalelemente: $\begin{eqnarray} \|\lambda - a_{ii}\| < \|a_{ii}\| \end{eqnarray}$ $\begin{eqnarray} \Rightarrow \lambda \in (0, 2a_{ii}) \Rightarrow \lambda >0 \end{eqnarray}$ $\begin{eqnarray} \text{Sei A symmetrisch und positiv definit. Dann gilt:} \end{eqnarray}$ ************************************************ $\begin{eqnarray} \text{die Diagonalelemenete sind positive, d.h. }a_{ii} > 0 \end{eqnarray}$ $\begin{eqnarray} \text{Jedes betragsmaximale Element von A ist ein Diagonalement} \end{eqnarray}$ Beweis Da A positiv definit ist, gilt speziell für die Einheitsvektoren $\begin{eqnarray} e_i^TAe_i > 0 \end{eqnarray}$ . Damit folgt die Behauptung. Angenommen, $\begin{eqnarray} a_{ij},~j \neq i \end{eqnarray}$ sei das betragsmaximale Element von A. Es ergibt sich folgende Fallunterscheidung: $\begin{eqnarray} a_{ij}>0 \end{eqnarray}$ Sei $\begin{eqnarray} x \in \mathbb R^n \end{eqnarray}$ so gewählt, dass gilt $\begin{eqnarray} x_k = 0~\forall~k=\{1,...,n\}\setminus \{i,j\},~x_i = -1,~x_j =1 \end{eqnarray}$ . Dann folgt mit $\begin{eqnarray} y_k = -a_{ki} + a_{kj},~\forall k =1,...,n \end{eqnarray}$ $\begin{eqnarray} x^TAx = x^Ty = a_{ii}-a_{ij}-a_{ji}+a_{j,j} = (a_{ii}+a_{jj}) - 2a_{ij} < 0 \end{eqnarray}$ $\begin{eqnarray} a_{ij}<0 \end{eqnarray}$ Sei $\begin{eqnarray} x \in \mathbb R^n \end{eqnarray}$ so gewählt, dass gilt $\begin{eqnarray} x_k = 0~\forall~k=\{1,...,n\}\setminus \{i,j\},~x_i = 1,~x_j =1 \end{eqnarray}$ . Dann folgt mit $\begin{eqnarray} y_k = a_{ki} + a_{kj},~\forall k =1,...,n \end{eqnarray}$ $\begin{eqnarray} x^TAx = x^Ty = a_{ii}+a_{ij}+-a_{ji}+a_{j,j} = (a_{ii}+a_{jj}) + 2a_{ij} < 0 \end{eqnarray}$ Dann ist aber A nicht positiv definit. somit ist die Behauptung falsch und das betragsmaximale Element ist ein Diagonalelement.
19.07.2007, 13:59	tigerbine	Auf diesen Beitrag antworten »
1d. Dreiecksmatrizen Produkt von Dreiecksmatrizen Sind $\begin{eqnarray} L^{(1)},~L^{(2)} \end{eqnarray}$ zwei untere (normierte) Dreiecksmatrizen, so ist auch $\begin{eqnarray} L^{(1)}L^{(2)} \end{eqnarray}$ eine untere (normierte) Dreiecksmatrix. Sind $\begin{eqnarray} R^{(1)},~R^{(2)} \end{eqnarray}$ zwei obere (normierte) Dreiecksmatrizen, so ist auch $\begin{eqnarray} R^{(1)}R^{(2)} \end{eqnarray}$ eine obere (normierte) Dreiecksmatrix. Inverse von Dreiecksmatrizen Ist $\begin{eqnarray} L \end{eqnarray}$ eine reguläre (normierte) untere Dreiecksmatrix, so ist auch $\begin{eqnarray} L^{-1} \end{eqnarray}$ eine untere (normierte) Dreiecksmatrix. Ist $\begin{eqnarray} R \end{eqnarray}$ eine reguläre (normierte) obere Dreiecksmatrix, so ist auch $\begin{eqnarray} L^{-1} \end{eqnarray}$ eine obere (normierte) Dreiecksmatrix. Beweis 1: Den Beweis führt man durch Aufschreiben der Matrizenmultiplikation. Es ist für $\begin{eqnarray} L:=L^{(1)}L^{(2)},~\quad i > j \end{eqnarray}$ : $\begin{eqnarray} l_{ij} &&= \sum_{m=1}^n~l_{im}^{(1)}l_{mj}^{(2)} = \\ && = \sum_{m=1}^{j-1}~l_{im}^{(1)}l_{mj}^{(2)} + \sum_{m=j}^n~l_{im}^{(1)}l_{mj}^{(2)} = 0 \end{eqnarray}$ Handelt es sich um normierte untere Dreiecksmatrizen, so folgt analog: $\begin{eqnarray} l_{ii} &&= \sum_{m=1}^n~l_{im}^{(1)}l_{mi}^{(2)} = \\ && = \sum_{m=1}^{i-1}~l_{im}^{(1)}l_{mi}^{(2)} + l_{ii}^{(1)}l_{ii}^{(2)} +\sum_{m=i+1}^n~l_{im}^{(1)}l_{mi}^{(2)} = 1 \end{eqnarray}$ Der Beweis für B folgt bereits aus A, da er sich durch Transponieren auf diesen Fall zurückführen lässt. Beweis 2: Setze $\begin{eqnarray} X:=L^{-1} \end{eqnarray}$ , dann gilt $\begin{eqnarray} LX=I \end{eqnarray}$ und Spaltenweise genommen $\begin{eqnarray} Lx_i = e_i, ~\forall i=1,...n \end{eqnarray}$ . Nun schreibt man die Matrizen als Blockmatrizen, derart: $\begin{eqnarray} L =\begin{pmatrix} L_{11}^{i} &0 \\ L_{21}^{i}& L_{22}^{i} \end{pmatrix},~x_i = \begin{pmatrix} x_{1}^{i} \\ x_{2}^{i} \end{pmatrix},~e_i= \begin{pmatrix} 0 \\ e_1 \end{pmatrix} \end{eqnarray}$ mit $\begin{eqnarray} L_{11}^{(i)} \in \mathbb R^{(i-1) \times (i-1)} \end{eqnarray}$ $\begin{eqnarray} L_{21}^{(i)} \in \mathbb R^{n-(i-1) \times (i-1)} \end{eqnarray}$ $\begin{eqnarray} L_{22}^{(i)} \in \mathbb R^{n-(i-1) \times n-(i-1)} \end{eqnarray}$ $\begin{eqnarray} x_1^{(i)} \in\mathbb R^{(i-1)} \end{eqnarray}$ $\begin{eqnarray} x_2^{(i)} \in\mathbb R^{n-(i-1)} \end{eqnarray}$ $\begin{eqnarray} e_1 \in\mathbb R^{n-(i-1)} \end{eqnarray}$ so kann man schreiben: $\begin{eqnarray} \begin{pmatrix} L_{11}^{i} &0 \\ L_{21}^{i}& L_{22}^{i} \end{pmatrix} \begin{pmatrix} x_{1}^{i} \\ x_{2}^{i} \end{pmatrix} = \begin{pmatrix} 0 \\ e_1 \end{pmatrix} \end{eqnarray}$ Hieraus folgt: $\begin{eqnarray} L_{11}^{i}x_1^{i} = 0 \Rightarrow x_1^{i} = 0 \end{eqnarray}$ Somit ist die Inverse Matrix eine untere Dreiecksmatrix. Ist L normiert, so folgt $\begin{eqnarray} L_{22}^{(i)}x_2^{(i)} = e_1 \end{eqnarray}$ Da $\begin{eqnarray} L_{22}^{(i)} \end{eqnarray}$ eine normierte untere Dreiecksmatrix ist, folgt, dass das i-te Element des Spaltenvektors $\begin{eqnarray} x_i \end{eqnarray}$ eine Eins sein muss. Da dies für alle i=1,...,n gilt, ist Die Inverse ebenfalls normiert.
Anzeige

19.07.2007, 14:29	tigerbine	Auf diesen Beitrag antworten »
1e. Frobeniusmatrizen Betrachtet wird hier nu das Rechnung mit Frobeniusmatrizen. Inversion Da $\begin{eqnarray} L_i \end{eqnarray}$ eine normierte untere Dreiecksmatrix ist, existiert die Inverse und ist ebenfalls normiert. Für ihre Gestalt sei die Behauptung: $\begin{eqnarray} L_i=I - l_i\cdot e_i^T = \begin{pmatrix}1&&&&& \\ &\ddots &&&& \\ &&1&&& \\ &&-l_{i+1,i} &\ddots && \\ &&\vdots &&\ddots & \\&&-l_{n,i}&&&1 \end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} L_i^{(-1)}=I + l_i\cdot e_i^T = \begin{pmatrix}1&&&&& \\ &\ddots &&&& \\ &&1&&& \\ &&+l_{i+1,i} &\ddots && \\ &&\vdots &&\ddots & \\&&+l_{n,i}&&&1 \end{pmatrix} \end{eqnarray}$ mit $\begin{eqnarray} l_i=(0,...,0,l_{i+1,i},...,l_{n,i})^T \end{eqnarray}$ und $\begin{eqnarray} e_i \end{eqnarray}$ als i-tem Einheitsvektor. Es gilt dann offensichtlich $\begin{eqnarray} e_i^Tl_i = 0 \end{eqnarray}$ und somit ergibt sich: $\begin{eqnarray} (I-l_ie_i^T)(I+l_ie_i^T) = I -l_ie_i^T + l_ie_i^T - l_ie_i^Tl_ie_i^T = I \end{eqnarray}$ Produkt spezieller Frobenius-Matrizen Auch hier soll zunächst die Behauptung allgemein formuliert werden: $\begin{eqnarray} L_1^{-1}L_2^{-1}...L_{k}^{-1} = I + \sum_{i=1}^{k}~l_ie_i^T = \begin{pmatrix} 1&&&&&\\ l_{2,1} & & \ddots &&&& \\ \vdots & \ddots &1 &&& \\ l_{k+1,1} & \hdots & l_{k+1,k} & 1 && \\ \vdots &\ddots &\vdots &&\ddots & \\ l_{n,1} & \vdots &\l_{n,k} &&&1 \end{pmatrix} \end{eqnarray}$ Der Beweis wird dann über Induktion geführt. Der Induktionsanfang folgt aus dem vorherigen Beweis. sei die Aussage nun richtig für ein beligeiges $\begin{eqnarray} 1 \leq k < n-1 \end{eqnarray}$ . Dann folgt wegen $\begin{eqnarray} e_i^Tl_{k+1} = 0\quad \forall i =1,...,k \end{eqnarray}$ sowie mit den Darstellungen der Inversen Matrizen aus dem vorherigen Beweis und mit der IV: $\begin{eqnarray} L_1^{-1}L_2^{-1}...L_{k}^{-1} =(I + \sum_{i=1}^{k}~l_ie_i^T) (I+l_{k+1}e_{k+1}^T) = I+\sum_{i=1}^{k}~l_ie_i^T \end{eqnarray}$
19.07.2007, 14:53	tigerbine	Auf diesen Beitrag antworten »
2a. Diagonalmatrix Im Falle einer Diagonalmatrix gestaltet sich das Lösen des LGS ziemlich einfach. $\begin{eqnarray} \begin{pmatrix} d_{11}&&0 \\ &\ddots&\\0&&d_{nn}\end{pmatrix} \begin{pmatrix} x_1\\ \vdots \\ x_n\end{pmatrix} = \begin{pmatrix} b_1 \\ \vdots \\ b_n \end{pmatrix} \end{eqnarray}$ Somit gilt für die Lösung: $\begin{eqnarray} x_i = \frac{b_i}{d_{ii}} \end{eqnarray}$ Der daraus resultierende Algorithmus ist wohl definiert, da D regulär ist. Die Laufzeit beträgt $\begin{eqnarray} \mathcal O(n) \end{eqnarray}$ .
19.07.2007, 14:55	tigerbine	Auf diesen Beitrag antworten »
2b. untere Dreiecksmatrix Im Fall einer unteren Dreiecksmatrix hat das LGS die Gestalt Lx=b: $\begin{eqnarray} \begin{pmatrix} l_{11}&&0 \\ \vdots &\ddots&\\l_{n1} & \hdots &l_{nn}\end{pmatrix} \begin{pmatrix} x_1\\ \vdots \\ x_n\end{pmatrix} = \begin{pmatrix} b_1 \\ \vdots \\ b_n \end{pmatrix} \end{eqnarray}$ Der daraus resultierende Lösungsalgorithmus heißt Vorwärtssubstitution: $\begin{eqnarray} x_1 = \frac{b_1}{l_{11}} \end{eqnarray}$ FOR i = 2:1:n $\begin{eqnarray} x_i=\frac{(b_i-\sum_{j=1}^{i-1}l_{ij}\cdot x_j)}{l_{ii}} \end{eqnarray}$ END Auch hier ist die Durchfürhrbarkeit an die Regularität der Matrix L geknüpft. Die Anzahl der Rechenoperarationen ist n². Die Laufzeit beträgt $\begin{eqnarray} \mathcal O(n^2) \end{eqnarray}$ .
19.07.2007, 17:12	tigerbine	Auf diesen Beitrag antworten »
2c. obere Dreiecksmatrix Im Fall einer oberen Dreiecksmatrix hat das LGS die Gestalt Rx=b: $\begin{eqnarray} \begin{pmatrix} r_{11}& \hdots & r_{1n} \\ &\ddots& \vdots\\0 & &r_{nn}\end{pmatrix} \begin{pmatrix} x_1\\ \vdots \\ x_n\end{pmatrix} = \begin{pmatrix} b_1 \\ \vdots \\ b_n \end{pmatrix} \end{eqnarray}$ Der daraus resultierende Lösungsalgorithmus heißt Rückwärtssubstitution: $\begin{eqnarray} x_n = \frac{b_n}{r_{nn}} \end{eqnarray}$ FOR i = (n-1):-1:1 $\begin{eqnarray} x_i=(b_i-\sum_{j=i+1}^{n}r_{ij}\cdot x_j)/r_{ii} \end{eqnarray}$ END Auch hier ist die Durchfürhrbarkeit an die Regularität der Matrix R geknüpft. Die Anzahl der Rechenoperarationen ist n². Die Laufzeit beträgt $\begin{eqnarray} \mathcal O(n^2) \end{eqnarray}$ .
19.07.2007, 17:45	tigerbine	Auf diesen Beitrag antworten »
3. LR-Zerlegung ohne Pivotisierung Unter einer LR-Zerlegung (ohne Pivotisierung) einer Matrix $\begin{eqnarray} A \in \mathbb R^{n \times n} \end{eqnarray}$ versteht man eine Zerlegung in der Form $\begin{eqnarray} A = LR \end{eqnarray}$ , wobei $\begin{eqnarray} L \in \mathbb R^{n \times n} \end{eqnarray}$ eine normierte untere Dreiecksmatrix und $\begin{eqnarray} R \in \mathbb R^{n \times n} \end{eqnarray}$ eine obere Dreiecksmatrix ist.
19.07.2007, 17:46	tigerbine	Auf diesen Beitrag antworten »
3a. Berechnung Wie berechnet man eine LR-Zerlegung? Aus der Schulzeit erinnert man sich vielleicht noch an das Lösen von Linearen Gleichungssystemen. Dazu wurde die Matrix A mittels Elementarer Zeileinumformungen auf Dreiecksgestalt R gebracht (Änderungen wirken sich auch auf b aus) und das neue LGS wurde durch Rückwärtssubstitution gelöst. Würde man die an A,b durchgeführten Änderungen in einer Matrix $\begin{eqnarray} L^{-1} \end{eqnarray}$ L speichern, so ergäbe sich: $\begin{eqnarray} L^{-1}Ax = L^{-1}b \Leftrightarrow Rx = y \end{eqnarray}$ Damit sei die Motivation erklärt, wie man mittels Frobeniusmatrizen letztendlich die LR-Zerlegung der Matrix konstruiert. Die Konstruktion des Algorithmus soll hier dabei an einem Beispiel erläutert werden: $\begin{eqnarray} A= \begin{pmatrix} 2& 1&1 &0 \\ 4&3&3&1 \\ 8&7&9&5\\6&7&9&8 \end{pmatrix} \end{eqnarray}$ Nun sollen die Eintrage in der ersten Spalte unterhalb der Diagonale annulliert werden. Dies geschieht (Nachrechnen) durch folgende Matrizenmultiplikation: $\begin{eqnarray} L_1A= \begin{pmatrix}1&0&0&0 \\ -2&1 &0 &0 \\ -4 &0&1&0 \\ -3&0&0&1 \end{pmatrix} \begin{pmatrix} 2& 1&1 &0 \\ 4&3&3&1 \\ 8&7&9&5\\6&7&9&8 \end{pmatrix} = \begin{pmatrix} 2& 1&1 &0 \\ 0&1&1&1 \\ 0&3&5&5\\0&4&6&8 \end{pmatrix} = A^{(2)} \end{eqnarray}$ Nun ist die zweite Spalte an der Reihe: $\begin{eqnarray} L_2A^{(2)} = \begin{pmatrix}1&0&0&0 \\ 0&1 &0 &0 \\ 0 &-3&1&0 \\ 0&-4&0&1 \end{pmatrix} \begin{pmatrix} 2& 1&1 &0 \\ 0&1&1&1 \\ 0&3&5&5\\0&4&6&8 \end{pmatrix} = \begin{pmatrix} 2& 1&1 &0 \\ 0&1&1&1 \\ 0&0&2&2\\0&0&2&4 \end{pmatrix} = A^{(3)} \end{eqnarray}$ Bleibt noch die dritte Spalte: $\begin{eqnarray} L_3A^{(3)} = \begin{pmatrix}1&0&0&0 \\ 0&1 &0 &0 \\ 0 &0&1&0 \\ 0&0&-1&1 \end{pmatrix} \begin{pmatrix} 2& 1&1 &0 \\ 0&1&1&1 \\ 0&0&2&2\\0&0&2&4 \end{pmatrix} =\begin{pmatrix} 2& 1&1 &0 \\ 0&1&1&1 \\ 0&0&2&2\\0&0&0&2 \end{pmatrix} =R \end{eqnarray}$ Nun kommen 2 Interessante Details. Man überprüfe durch nachrechnen, dass gilt (2),(3): $\begin{eqnarray} L_1^{-1} = \begin{pmatrix}1&0&0&0 \\ 2&1 &0 &0 \\ 4 &0&1&0 \\ 3&0&0&1 \end{pmatrix},~L_2^{-1} = \begin{pmatrix}1&0&0&0 \\ 0&1 &0 &0 \\ 0 &3&1&0 \\ 0&4&0&1 \end{pmatrix},~L_3^{-1}=\begin{pmatrix}1&0&0&0 \\ 0&1 &0 &0 \\ 0 &0&1&0 \\ 0&0&1&1 \end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} L:=L_1^{-1}L_2^{-1}L_3^{-1} = \begin{pmatrix}1&0&0&0 \\ 2&1 &0 &0 \\ 4 &3&1&0 \\ 3&4&1&1 \end{pmatrix} \end{eqnarray}$ Dieses einfache Übertragen der Einträge funktioniert nicht bei $\begin{eqnarray} L_3L_2L_1 \end{eqnarray}$ , denn: $\begin{eqnarray} L_3L_2L_1 = \begin{pmatrix} 1&0&0&0 \\-2&1&0&0 \\2&-3&1&0\\3&-1&-1&1 \end{pmatrix} \end{eqnarray}$
19.07.2007, 17:58	tigerbine	Auf diesen Beitrag antworten »
3b. Existenzsatz 1 Die reguläre Matrix A besitzt eine LR-Zerlegung ohne Pivotisierung $\begin{eqnarray} \Leftrightarrow \end{eqnarray}$ Alle Hauptabschnittsmatrizen $\begin{eqnarray} A[k] \end{eqnarray}$ regulär sind. Eine Hauptabschnittsmatrix ist genau dann regulär, wenn gilt $\begin{eqnarray} \det(A[k]) \neq 0 \end{eqnarray}$ . Beweis: $\begin{eqnarray} \Rightarrow \end{eqnarray}$ Besitze A eine LR-Zerlegung. Dann existiert eine normierte untere Dreiecksmatrix $\begin{eqnarray} L \in \mathbb R^{n \times n} \end{eqnarray}$ und eine obere Dreiecksmatrix $\begin{eqnarray} R \in \mathbb R^{n\times n} \end{eqnarray}$ mit $\begin{eqnarray} A=LR \end{eqnarray}$ . Aus $\begin{eqnarray} \det(A)=\det(L)\cdot \det(R) \end{eqnarray}$ und der Regularität von A folgt $\begin{eqnarray} \det(L)\neq 0,~\det(R) \neq 0 \end{eqnarray}$ . Da es sich um Dreiecksmatirzen handelt folgt sofort, dass auch alle ihre Hauptabschnittsmatrizen regulär sind. Durch "Nachrechnen" zeigt man die Gültigkeit von $\begin{eqnarray} (LA)[k] = L[k]A[k] \quad \forall~k=1,...,n~() \end{eqnarray}$ für jede (nicht notwendigerweise normierte) Dreiecksmatrix $\begin{eqnarray} L \in \mathbb R^{n \times n} \end{eqnarray}$ . Damit folgt: $\begin{eqnarray} \det(A[k]) &&= \det(LR[k]) = \det(L[k]R[k])= \\&&=\det(L[k])\cdot\det(R[k]) \neq 0\qquad \forall~k=1,...,n \end{eqnarray}$ $\begin{eqnarray} \Leftarrow \end{eqnarray}$ Seien alle Hauptabschnittsmatrizen positiv. Die Existenz der LR-Zerlegung weißt man durch die in diesem Falle Wohldefiniertheit des Algorithmus ihrer Berechnung nach. Es muss also gezeigt werden, dass alle Pivotelemente $\begin{eqnarray} a_{ii}^{(i)},~i = 1,...,n-1 \end{eqnarray}$ von Null verschieden sind. Der Beweis wird durch Induktion geführt: $\begin{eqnarray} \text{IA: }&& i=1,~a_{11}^{(1)} = \det(A[1]) \neq 0 \text{ nach Voraussetzung}\\\\\text{IV: } && \text{ Sei nun }a_{ii}^{i} \neq 0 \text{ für ein } 1 \leq i < n-1\\\\\text{IS: }&& \text{Aus der IV folgt die Existenz einer Matrix} \\ && A^{(i+1)} = L_i \cdot L_1 \cdot A \\\\&&\text{Durch sukzessive Anwendung von () folgt} \\&&A^{(i+1)}[i+1] = L_i[i+1]...L_1[i+1]A[i+1]\\\\&&\text{Mit der Determinantenmultiplikation folgt dann} \\&& \det(A^{(i+1)})[i+1]) = \det(L_i[i+1])\cdot ...\cdot \det(L_1[i+1]) \cdot \det(A[i+1]) \\\\ &&\text{Da } A^{(i+1)}[i+1] \text{offenbar eine obere Dreiecksmatrix darstellt mit} \\&& a_{i+1,i+1}^{(i+1)} \neq 0 \text{ als letzten Diagonalelement, folgt hieraus wegen } \\&& \det(A^{(i+1)}[i+1]) \neq 0 \text{ sofort } a_{i+1,i+1}^{(i+1)} \neq 0 \end{eqnarray}$
19.07.2007, 18:11	tigerbine	Auf diesen Beitrag antworten »
3c. Existenzsatz 2 Die reguläre Matrix A ist diagonaldominant $\begin{eqnarray} \Rightarrow \end{eqnarray}$ Es existiert eine LR Zerlegung A=LR ohne Pivotisierung. Beweis $\begin{eqnarray} &&\text{Da A regulär und diagonaldominant ist, muss } a_{11} \neq 0\text{ gelten} \\&& \text{Schritt 1 ist daher ohne Pivotisierung möglich} \\ \\&&\text{Es gilt dann:} \\ && a_{jk}^{(1)} = \begin{cases} a_{1k} & j=1,~\qquad k=1,...,n \\ a_{jk}-\frac{a_{j1}}{a_{11}}\cdot a_{1k} & j=2,...,n,~k=1,...,n \end{cases} \\\\ && \end{eqnarray}$ Die Matrix hat dann die Gestalt: $\begin{eqnarray} A^{(1)} = \begin{pmatrix} a_{11} & \hdots & \hdots &a_{1n} \\ 0& & \hdots &* \\ 0 &* & \hdots &* \\ 0&* & \hdots &* \end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} \text{Setzt man nun } q:= \frac{a_{j1}}{a_{11}} \text{, so gilt:} \end{eqnarray}$ $\begin{eqnarray} \sum_{k=2 \atop k \neq j}^n~\|a_{jk}^{(1)}\| && \leq \sum_{k=2 \atop k \neq j}^n~\|a_{jk}\| + \|q\|\cdot \sum_{k=2 \atop k \neq j}^n~\|a_{1k}\| \\ && = \sum_{k=1 \atop k \neq j}^n~\|a_{jk}\| -\|a_{j1}\| + \|q\|\cdot \sum_{k=1 \atop k \neq j}^n~\|a_{1k}\|- \|q\|\cdot \|a_{1j}\| \\ && \leq \|a_{jj}\| - \|a_{j1}\| + \| \frac{a_{j1}}{a_{11}}\| \cdot \|a_{11}\| - \|q\| \cdot \|a_{1j}\| \\&& = \|a_{jj}\| - \|q \cdot a_{1j}\| \\&& \leq \|a_{jj}^{(1)}\| \end{eqnarray}$ $\begin{eqnarray} \text{Somit ist auch } A^{(1)} \text{diagonaldominant und da auch regulär (Blockdiagonalmatrix )} \end{eqnarray}$ $\begin{eqnarray} \Rightarrow a_{22} \neq 0 \end{eqnarray*}$ Man kann also alle Schritte ohne Pivotisierung durchführen. Im Falle der strikten Diagonaldominanz folgt hiermit sogar direkt, dass A regulär ist.
19.07.2007, 20:27	tigerbine	Auf diesen Beitrag antworten »
3d. Eindeutigkeit Seien $\begin{eqnarray} A=L_1R_1 \end{eqnarray}$ und $\begin{eqnarray} A=L_2R_2 \end{eqnarray}$ zwei LR-Zerlegungen von A. Dann gilt: $\begin{eqnarray} L_2R_2 = L_1R_1 \end{eqnarray}$ Aufgrund der Regularität von A folgt dann: $\begin{eqnarray} R_2R_1^{-1} = L_2^{-1}L_1 \end{eqnarray}$ Damit folgt wegen der Gestalt (L ist normiert), dass gelten muss: $\begin{eqnarray} L_2^{-1}L_1 = I = R_2R_1^{-1} \end{eqnarray}$ Somit folgt: $\begin{eqnarray} L_1=L_2,~R_1 = R_2 \end{eqnarray}$ Und die Zerlegung ist Eindeutig.
19.07.2007, 21:14	tigerbine	Auf diesen Beitrag antworten »
3e. Algorithmus zur Lösung des LGS Ax = b Gemäßt 2.6 bedeutet gerade die Gestalt einer Frobeniusmatrix. Bei der Durchführung der LR-Zerlegung wird man die Elemente der Matrix A sukzessive überschreiben. Dabei geht die Diagonale von L verloren, aber die ist ja definiert durch lauter Einsen. Der Algorithmus sieht dann wie folgt aus: Nun soll ja das LGS Ax = b, bzw LRx = b gelöst werden. Dazu löst man zunächst Ly = b durch Vorwärtssubstitution und speichert y durch überschreiben von b. Danach löst man Rx=y durch Rückwärtssubstitution und speichert wieder die Lösung in b. Die Anzahl der Rechenoperationen der Gaußeleiminationen ist $\begin{eqnarray} \approx \frac{4}{6}n^3 \end{eqnarray}$ , so dass man ein LGS in der Laufzeit $\begin{eqnarray} \mathcal O(n^3) \end{eqnarray}$ lösen kann.
19.07.2007, 21:27	tigerbine	Auf diesen Beitrag antworten »
4. LR-Zerlegung mit Pivotisierung In der Praxis wird jedoch häufiger der Fall einer LR-Zerlegung mit Pivotisierung auftreten, d.h. die Matrix A muss erst noch "aufbereitet" werden. Man schreibt die Zerlegung dann in der Form $\begin{eqnarray} PA = LR \end{eqnarray}$ , wobei P eine Permutationsmatrix ist. Selbst wenn eine Pivotisierung nicht nötig ist, ist es am Ende günstiger die Matrix gleich durch diesen Algorithmus zu schicken. Denn bei der LR-Zerlegung ohne Pivotisierung muss man ja noch den Aufwand hinzunehmen, nachzuprüfen, ob diese Überhaupt möglich ist. Ein weiterer Vorteil ist, dass man hier noch nicht einmal die Regularität von A prüfen muss. Bricht der Algorithmus ab, so ist A numerisch singulär.
19.07.2007, 21:29	tigerbine	Auf diesen Beitrag antworten »
4a. Berechnung An einem Beispiel soll wieder die Idee verdeutlicht werden. $\begin{eqnarray} A = \begin{pmatrix} 0&2&-1&-2 \\ 2& -2 & 4 &-1 \\ 1&1&1&1 \\-2 & 1 &-2 & 1 \end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} P_{12}A = \begin{pmatrix} 2& -2 & 4 &-1 \\ 0&2&-1&-2 \\ 1&1&1&1 \\-2 & 1 &-2 & 1 \end{pmatrix} \end{eqnarray}$ Darauf läßt sich nun ein Schritt der LR-Zerlegung (3a) durchführen. $\begin{eqnarray} L_1P_{12}A = \begin{pmatrix} 2&-2&4&-1 \\ 0&2&-1&-2\\ 0 &2& -1&1.5 \\ 0&-1&2&0 \end{pmatrix} \end{eqnarray}$ mit $\begin{eqnarray} L_1:=\begin{pmatrix} 1&0&0&0 \\ 0&1&0&0\\ -0.5&0&1&0 \\ 1&0&0&1 \end{pmatrix} \end{eqnarray}$ Im nächsten Schritt ist keine Vertauschung nötig. Dennoch schreibt man P=I hin. Warum werden wir später sehen. $\begin{eqnarray} L_2P_{22}L_1P_{12}A = \begin{pmatrix} 2&-2&4&-1 \\ 0&2&-1&-2\\ 0 &0& 0&3.5 \\ 0&0&1.5& -1 \end{pmatrix} \end{eqnarray}$ mit $\begin{eqnarray} L_2:=\begin{pmatrix} 1&0&0&0 \\ 0&1&0&0\\ 0 &-1&1&0 \\ 0&0.5&0&1 \end{pmatrix} \end{eqnarray}$ Im letzten Schritt ist nun wieder ein Vertauschen nötig. $\begin{eqnarray} P_{34}L_2P_{22}L_1P_{12}A = \begin{pmatrix} 2&-2&4&-1 \\ 0&2&-1&-2\\ 0&0&1.5& -1 \\ 0 &0& 0&3.5 \end{pmatrix} \end{eqnarray}$ Damit ist schon die gewünschte Dreiecksgestalt erreicht. Dennoch notiert man auch hier: $\begin{eqnarray} L_3P_{34}L_2P_{22}L_1P_{12}A = \begin{pmatrix} 2&-2&4&-1 \\ 0&2&-1&-2\\ 0&0&1.5& -1 \\ 0 &0& 0&3.5 \end{pmatrix} \end{eqnarray}$ mit $\begin{eqnarray} L_3=I \end{eqnarray}$ Wie hängen nun die theoretisch existierende Darstellung PA=LR und die eben berechnete zusammen? Dazu bedient man sich folgendes "Tricks", der am Beispiel illustriert werden soll: $\begin{eqnarray} &&L_3':=L_3,\\ &&L_2':=P_3L_2P_3^{-1},\\ &&L_1':=P_3P_2L_1P_2^{-1}P_3^{-1} \end{eqnarray}$ Diese Matrizen haben sind wieder untere Dreiecksmatrizen mit der gleichen Struktur (Beweis unten). Die Permutationsmatrizen geben in ihrem Index hier nur an, welche Zeile getauscht werden soll. Den Tauschpartner läßt man hier nun einmal außen vor(siehe auch hier). Somit ergibt sich: $\begin{eqnarray} L_3'L_2'L_1'P_3P_2P_1 && = L_3(P_3L_2P_3^{-1}) (P_3P_2L_1P_2^{-1}P_3^{-1})P_3P_2P_1\\ && =L_3P_3L_2(P_3^{-1}P_3)P_2L_1[P_2^{-1}(P_3^{-1}P_3)P_2]P_1 \\&&=L_3P_3L_2P_2L_1P_1 \end{eqnarray}$ Und damit die Zerlegung: $\begin{eqnarray} PA=LR \text{ mit } P:=P_3P_2P_1,~L:=(L_3'L_2'L_1')^{-1} \end{eqnarray}$ Dabei ist P wieder eine Permutationsmatrix und L wieder eine untere normierte Dreiecksmatrix. Im Beispiel also: $\begin{eqnarray} P_1=\begin{pmatrix} 0&1&0&0 \\ 1&0&0&0 \\0&0&1&0 \\0 &0&0&1\end{pmatrix},~P_2=\begin{pmatrix} 1&0&0&0 \\ 0&1&0&0 \\0&0&1&0 \\0 &0&0&1\end{pmatrix},~P_3=\begin{pmatrix} 1&0&0&0 \\ 0&1&0&0 \\0&0&0&1 \\0 &0&1&0\end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} \Rightarrow P =\begin{pmatrix} 0&1&0&0 \\ 1&0&0&0 \\0&0&0&1 \\0 &0&1&0\end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} L_3' = \begin{pmatrix} 1&0&0&0 \\ 0&1&0&0 \\0&0&1&0 \\0 &0&0&1\end{pmatrix},~L_2'=\begin{pmatrix} 1&0&0&0 \\ 0&1&0&0\\ 0 &0.5&1&0 \\ 0&-1&0&1 \end{pmatrix},~L_1'=\begin{pmatrix} 1&0&0&0 \\ 0&1&0&0\\ 1&0&1&0 \\ -0.5&0&0&1 \end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} \Rightarrow L =\begin{pmatrix} 1&0&0&0 \\ 0&1&0&0\\ -1 &-0.5&1&0 \\ 0.5&1&0&1 \end{pmatrix} \end{eqnarray}$ ********************************************************* Beweis:** Hier soll einmal die Struktur der konstruierten Matrizen betrachtet werden. Eine Permutationsmatrix ist zu sich selbst invers. Allgemein definiert ist nun: $\begin{eqnarray} L_i':=P_{n-1}...P_{i+1}~L_i~P_{i+1}^{-1}...P_{n-1}^{-1} \end{eqnarray}$ Betrachten wir erstmal nur die Struktur von: $\begin{eqnarray} P_{i+1}~L_i~P_{i+1}^{-1} \end{eqnarray}$ mit der Rechtsmultiplikation mit $\begin{eqnarray} P_{i+1}^{-1}=P_{i+1} \end{eqnarray}$ werden die Spalte (i+1) und eine Spalte $\begin{eqnarray} i+1\leq j \leq n \end{eqnarray}$ vertauscht. Für diese gilt als Spaltenvektoren: $\begin{eqnarray} l_{i+1}=(0,...,0,\underbrace{1}_{i+1},0,..0)^T,~\text{d.h. }l_{i+1,i+1}=1,~\text{ sonst 0} \end{eqnarray}$ $\begin{eqnarray} l_{j}=(0,...,0,\underbrace{1}_{j},0,..0)^T,~\text{d.h. }l_{j,j}=1,~\text{ sonst 0} \end{eqnarray}$ Nach dem ersten Tausch gilt also: $\begin{eqnarray} l_{i+1,k} = \begin{cases} 1& k=j \\0&\text{sonst} \end{cases},~l_{j,k} = \begin{cases} 1& k=i+1 \\0&\text{sonst} \end{cases} \end{eqnarray}$ Nun werden die Zeilen (i+1) und j getauscht, somit ergibt sich: $\begin{eqnarray} l_{i+1,k} =\begin{cases} 1& k=i+1 \\0&\text{sonst} \end{cases},~l_{j,k} =\begin{cases} 1& k=j \\0&\text{sonst} \end{cases} \end{eqnarray}$ Und somit ist die Struktur der Matrix erhalten geblieben, die Einträge $\begin{eqnarray} l_{i,i+1},~l_{i,j} \end{eqnarray}$ sind jedoch vertauscht worden.
19.07.2007, 21:31	tigerbine	Auf diesen Beitrag antworten »
4b. Existenz Jede reguläre Matrix $\begin{eqnarray} A \in \mathbb R^{n \times n} \end{eqnarray}$ besitzt eine Zerlegung in der Gestalt PA = LR mit mit einer(n x n) Permutationsmatrix P, einer (n x n) normierten unteren Dreiecksmatrix L sowie einer (n x n) oberen Dreiecksmatrix R. Beweis: Der Beweis wird durch Induktion geführt. $\begin{eqnarray} \text{IA:~}&&\text{Da A regulär ist, ergibt sich für n=1 offensichtlich P=L=I, R=a} \end{eqnarray}$ $\begin{eqnarray} \text{IV: } &&\text{Gelte die Beh. für alle reg. A }\in \mathbb R^{(n-1) \times (n-1)} \end{eqnarray}$ $\begin{eqnarray} \text{IS:~}&&\text{Sei A }\in \mathbb R^{n \times n} \text{ regulär} \\&&\text{Dann gibt es in der ersten Spalte mind. ein von 0 verschiedenes Element.} \\&&\text{Somit gibt es eine Permutationsmatrix } P_n \text{ derart:} \\&& P_nA = \begin{pmatrix} \alpha & r^T\\ l &C \end{pmatrix},~\alpha \in \mathbb R \setminus\{0\},~l,r \in \mathbb R^{(n-1)},~C\in \mathbb R^{n \times n} \\&&\text{Setzt man nun } B:=C-\frac{1}{\alpha} \text{ dann folgt:}\\&&P_nA = \begin{pmatrix}1&0\\\frac{1}{\alpha}l&I \end{pmatrix} \begin{pmatrix} \alpha & r^T\\ 0 & B \end{pmatrix} \\&&\text{Mit Determinantenrechnung folgt aus: }\\&&\det(P_nA)=\det(P_n)\det(A) \neq 0 \\&&\det(P_nA)=\alpha \cdot \det(B)\\&&\text{dass B regulär ist. Somt kann man die IV auf B anwenden: }\\&&P_{n-1}B = L_{n-1}R_{n-1} \end{eqnarray}$ $\begin{eqnarray} &&\text{Somit folgt durch Blockmatrixrechnung:}\\&&P_nA = \begin{pmatrix} 1&0\\0 & P_{n-1}^T\end{pmatrix} \begin{pmatrix} 1&0 \\\frac{1}{\alpha} P_{n-1 }l&L_{n-1}\end{pmatrix} \begin{pmatrix}\alpha & r^T\\0 & R_{n-1} \end{pmatrix}\\&& \text{Und damit die Zerlegung PA=LR.} \end{eqnarray}$
19.07.2007, 21:31	tigerbine	Auf diesen Beitrag antworten »
4c. Eindeutigkeit - Pivotwahl In den Überlegungen von (4b) wurde bislang nur festgestellt, dass man bei einer regulären Matrix A immer einen "Tauschpartner" für den Fall einer Null auf der Diagonale findet. Steht einem nicht nur ein "Tauschpartner" zur Verfügung, so stellt sich die Frage nach der "besten" Wahl. Den geringsten Aufwand wird wohl eine Suche nach dem ersten von 0 verschiedenen Eintrag in der betreffenden Spalte unterhalb der Diagonale machen. Ist diese "blinde" Wahl jedoch numerisch sinnvoll? Durch die Implementierung erhält man eine Zerlegung $\begin{eqnarray} \hat L \hat R \end{eqnarray}$ und löst das LGS der Form: $\begin{eqnarray} (A+E)\hat x = b \end{eqnarray}$ mit einer unbekannten Matrix E, für die folgende Abschätzung gilt (Golub & VanLoan[1996, S. 106]): $\begin{eqnarray} \fbox{$ \displaystyle \|E\| \leq n\cdot m (3\|A\| + 5\|\hat L\| \|\hat R\|) + \mathcal O(m^2) $} \end{eqnarray}$ m: Maschinengenauigkeit n: Ordnung von A, also ist A eine (n x n) Matrix \|.\|: Matrix der betragsmäßigen Einträgen Somit ist es erstrebenswert die Einträge von $\begin{eqnarray} \hat L \hat R \end{eqnarray}$ möglichst klein zu halten. 2 LR-Zerlegungen: $\begin{eqnarray} A:=\begin{pmatrix} \varepsilon & 1 \\1&1\end{pmatrix}=\begin{pmatrix}1 & 0\\\frac{1}{\varepsilon}&1 \end{pmatrix} \begin{pmatrix} \varepsilon&1\\ 0 & 1-\frac{1}{\varepsilon} \end{pmatrix}=LR \end{eqnarray}$ $\begin{eqnarray} A:=\begin{pmatrix} 1 & 1 \\\varepsilon&1\end{pmatrix}=\begin{pmatrix}1 & 0\\\varepsilon &1 \end{pmatrix} \begin{pmatrix} 1&1\\ 0 & 1-\varepsilon\end{pmatrix}=LR \end{eqnarray}$ $\begin{eqnarray} \text{mit }0<\varepsilon <<1 \end{eqnarray}$ Im ersten Fall treten also, im Gegensatz zum Zweiten, sehr große Werte auf. Man kann nun zeigen, dass allgemein für sehr kleine Diagonalelemente in der Matrix die Einträge in den Matrizen L und R sehr groß werden. Somit sollte man Kehrschluss die Diagonalelemente durch vertauschen möglichst groß wählen. Somit ergeben sich: weitere Pivotsuchstrategien: Spaltenpivotsuche Unter den Elementen $\begin{eqnarray} a_{j,i+1}^{(i)},~j=i+1,...,n \end{eqnarray}$ der j-ten Spalte wähle man ein dem Betrage nach größtes als neues Diagonalelement. Suchaufwand: $\begin{eqnarray} \mathcal O(n^2) \end{eqnarray}$ Totale Pivotsuche Unter den Elementen $\begin{eqnarray} a_{j,k}^{(i)},~j,k \leq i+1 \end{eqnarray}$ wähle man ein dem Betrag nach größtes als neues Diagonalelement aus. Hier müssen dann ggf. Zeilen und Spalten getauscht werden. Suchaufwand: $\begin{eqnarray} \mathcal O(n^3) \end{eqnarray}$ Bei größeren Matrizen kann es sinnvoll sein, nur so lange zu suchen, bis ein Schwellenwert s überschritten wird. Mit diesem Wissen kann man nun festhalten, dass auch im Falle der Existenz einer LR-Zerlegung ohne Pivotisierung es sinnvoll seien kann, selbige durchzuführen. Hierzu ein Beispiel: $\begin{eqnarray} A=\begin{pmatrix} 11 & 44 & 1 \\ 0.1 & 0.4 & 3 \\ 0 & 1 & -1 \end{pmatrix},~b = \begin{pmatrix} 1\\ 1 \\ 1\end{pmatrix} \end{eqnarray}$ Exakte Lösung $\begin{eqnarray} x^T=\frac{1}{329}(-1732,~438,~109) \approx (-5.26444,~1.33131,~0.331307) \end{eqnarray}$ Ohne Pivotisierung $\begin{eqnarray} A^{(1)} = \begin{pmatrix} 11 & 44 & 1 \\ 0 & -4.00000E-08 & 2.99091 \\ 0 & 1 & -1 \end{pmatrix},~b^{(1)} = \begin{pmatrix} 1\\ 0.990909 \\ 1\end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} A^{(2)} = \begin{pmatrix} 11 & 44 & 1 \\ 0 & -4.00000E-08 & 2.99091 \\ 0 & 0 & 7.47727E+07 \end{pmatrix},~b^{(1)} = \begin{pmatrix} 1\\ 0.990909 \\ 2.47727E+07\end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} \Rightarrow x^T=(-41.8765,~10.4843,~0.331307) \end{eqnarray}$ Mit Spalten-Pivotisierung $\begin{eqnarray} A^{(1)} = \begin{pmatrix} 11 & 44 & 1 \\ 0 & -4.00000E-08 & 2.99091 \\ 0 & 1 & -1 \end{pmatrix},~b^{(1)} = \begin{pmatrix} 1\\ 0.990909 \\ 1\end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} A^{(2)} = \begin{pmatrix} 11 & 44 & 1 \\ 0 & 1 & -1 \\ 0 & 0 & 2.99091 \end{pmatrix},~b^{(1)} = \begin{pmatrix} 1\\ 1\\ 0.990909 \end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} \Rightarrow x^T=(-5.26445,~1.331331,~0.331307) \end{eqnarray}$
19.07.2007, 22:03	tigerbine	Auf diesen Beitrag antworten »
4d. Numerische Hinweise Ein singuläres LGS wird vom Gaußverfahren dadurch erkannt, dass an einer bestimmten Stelle alle Pivotsuchen erfolglos verlaufen. Durch Rundungsfehler kann es aber passieren, dass das System trotzdem als unlösbar erscheint. Eine hieraus berechnete Lösung ist jedoch wenig sinnvoll und man wird in der Praxis Systeme mit Pivotelementen kleiner $\begin{eqnarray} \varepsilon =10^{-8} \end{eqnarray}$ als nicht behandelbar erklären. Aufgrund von Rechenungenaugkeiten wird in der Regel die berechnete Lösung nicht mit der exakten übereinstimmen. Vergleiche hierzu die Überlegungen aus dem [Workshop - lineare Gleichungssysteme 1] Ferner kann es trotzdem noch zu Problemen kommen, wenn die verschiedenen Gleichungen stark unterschiedlich gewichtet sind. Es empfiehlt sich eine Equilibrierung. (siehe Anhang)
20.07.2007, 00:28	tigerbine	Auf diesen Beitrag antworten »
5. Weitere Möglichkeiten Vielleicht erinnert man sich noch an die Lineare Algebra I, und wie man dort die Inverse einer Matrix berechnet hat. Dies und anderes was man neben der Lösung des LGS Ax=b noch aus der LR Zerlegung erhält, soll hier betrachtet werden.
20.07.2007, 01:02	tigerbine	Auf diesen Beitrag antworten »
5a. Determinantenberechnung Für quadratische Matrizen A und B gilt: $\begin{eqnarray} \det(AB) = \det(A) \cdot \det(B) \end{eqnarray}$ Für die Determinante einer Dreiecksmatrix gilt: $\begin{eqnarray} \det(L) = l_{11} \cdot ... \cdot l_{nn} \end{eqnarray}$ Determinante einer Permutationsmatrix: $\begin{eqnarray} \det(P_{ij}) =\det(P_{ij}^{-1})= -1,~\quad i\neq j \end{eqnarray}$ Determinante einer Frobeniusmatrix: $\begin{eqnarray} \det{L_i} = \det(L_i^{-1}) = 1 \end{eqnarray}$ Gestalt der LR-Zerlegung mit Pivotisierung: $\begin{eqnarray} L = (L_{n-1}\cdot P_{n-1}\cdot ...\cdot L_1\cdot P_1)^{-1} \end{eqnarray}$ Somit folgt: $\begin{eqnarray} \det(A) =\begin{cases} -\det(R) & \text{ungerade Anzahl an Zeilenvertauschungen} \\\det(R) & \text{gerade Anzahl an Zeilenvertauschungen} \\0 & \text{A ist singulär} \end{cases} \end{eqnarray}$ Die Dreieckszerlegung bietet so eine einfache Möglichkeit die Determinante zu berechnen und ist im Vergleich zur Leibnizformel auch kostengünstig. $\begin{eqnarray} \mathcal O(n³) << \mathcal O(n!) \end{eqnarray}$
20.07.2007, 01:32	tigerbine	Auf diesen Beitrag antworten »
5b. Rangebstimmung Falls es immer ein von Null verschiedenes Pivotelement gibt, so ist $\begin{eqnarray} \det(A) \neq 0 \end{eqnarray}$ und damit hat A den vollen Rang. Ist im i-ten Eliminationsschritt für alle Elemente der i-ten Spalte $\begin{eqnarray} a_{ji}^{(i-1)}=0, \quad j=i,...,n \end{eqnarray}$ , dann ist A singulär. $\begin{eqnarray} A^{i}=\begin{pmatrix}&&&\\ 0&0&0&0 \\0&&& \\0 &* &* &* \end{pmatrix} \end{eqnarray}$ Um nun den Vollen Rang bestimmen zu können, nimmt man eine Totalpivotisierung vor. Gilt dann nach dem i-ten Eliminationsschritt $\begin{eqnarray} a_{jk}^{(i)} = 0,~\quad j,k = i+1,...,n \end{eqnarray*}$ So ist Rang(A) = i. Gilt auch für nicht quadratische Matizen.
20.07.2007, 18:12	tigerbine	Auf diesen Beitrag antworten »
5c. Inversesbestimmung 1 $\begin{eqnarray} PA=LR \end{eqnarray}$ Gesucht ist $\begin{eqnarray} AA^{-1} = I \end{eqnarray}$ $\begin{eqnarray} PAA^{-1} =PI \end{eqnarray}$ Dieses kann man nun mit Gaußalgorithmus für jeden Spaltenvektor von $\begin{eqnarray} x^{(i)}:=(A^{-1})_i,~b^{i}=Pe_i \end{eqnarray}$ lösen Berechne die LR Zerlegung PA=LR Löse $\begin{eqnarray} Ly^{(i)}=Pe{(i)} \end{eqnarray}$ Löse $\begin{eqnarray} Rx^{i} =y^{(i)} \end{eqnarray}$ i=1,...,n $\begin{eqnarray} A^{-1} = \begin{pmatrix} \vdots & \vdots & \vdots \\x^{1} & \hdots & x^{(n)} \\ \vdots &\vdots&\vdots \end{pmatrix} \end{eqnarray}$
20.07.2007, 20:36	tigerbine	Auf diesen Beitrag antworten »
5d. Inversesbestimmung 2 - Gauß Jordan Algorithmus Man kann ein LGS auch mit dem Gauß-Jordan Algorithmus lösen. Der erste Teil ist Gauß-Algorithmus gleich, danach wird aus der oberen Dreiecksmatrix R jedoch noch eine Diagonalmatrix gemacht. Dazu multipliziert man z.B. von Links mit Matrizen (ähnlich Frobenius). Dabei ändern sich nur die Einträge in der entsprechenden Spalte oberhalb der Diagonalen. Man kann also die Matrizen $\begin{eqnarray} \hat{L} \end{eqnarray}$ direkt aus der Dreiecksmatrix angeben. Ebenso bleiben die Rechenregeln (Inversion, Multiplikation der Inversen) wie bei den Frobenismatrizen erhalten. $\begin{eqnarray} \hat L = \begin{pmatrix} 1 &&&& \\ &1&&&\\ &&1&&\\&&&1\end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} AA^{-1}&&=I \\\\ \hat L_{n-1}... \hat L_1\cdot L_{n-1}P_{n-1} \cdot...\cdot L_1P_1 AA^{-1} &&= \hat L_{n-1}... \hat L_1\cdot L_{n-1}P_{n-1} \cdot...\cdot L_1P_1I \\\\ DA^{-1} &&= \hat L_{n-1}... \hat L_1\cdot L_{n-1}P_{n-1} \cdot...\cdot L_1P_1 I \\\\ A^{-1} &&= D^{-1}\hat L_{n-1}... \hat L_1\cdot L_{n-1}P_{n-1} \cdot...\cdot L_1P_1 I \end{eqnarray}$ Beispiel: $\begin{eqnarray} A = \begin{pmatrix} 3&1&6 \\ 2&1&3 \\ 1&1&1 \end{pmatrix} \end{eqnarray}$ Schreibe $\begin{eqnarray} \begin{pmatrix} 3&1&6 \\ 2&1&3 \\ 1&1&1 \end{pmatrix} \quad \begin{pmatrix} 1&0&0 \\ 0&1&0 \\ 0&0&1 \end{pmatrix} \end{eqnarray}$ Vorwärtselimination 1 (Obere Dreiecksmatrix links erzeugen) $\begin{eqnarray} \begin{pmatrix} 3&1&6 \\ 0&\frac{1}{3}&-1 \\ 0& \frac{2}{3}&-1 \end{pmatrix} \quad \begin{pmatrix} 1&0&0 \\ -\frac{2}{3}&1&0 \\ -\frac{1}{3}&0&1 \end{pmatrix} \end{eqnarray}$ Pivotisierung 1 (größtes Spaltenelement) $\begin{eqnarray} \begin{pmatrix} 3&1&6 \\ 0& \frac{2}{3}&-1 \\ 0&\frac{1}{3}&-1\end{pmatrix} \quad \begin{pmatrix} 1&0&0 \\ -\frac{1}{3}&0&1 \\ -\frac{2}{3}&1&0 \end{pmatrix} \end{eqnarray}$ Vorwärtselimination 2 (Obere Dreiecksmatrix links erzeugen) $\begin{eqnarray} \begin{pmatrix} 3&1&6 \\ 0& \frac{2}{3}&-1 \\ 0&0&-\frac{1}{2}\end{pmatrix} \quad \begin{pmatrix} 1&0&0 \\ -\frac{1}{3}&0&1 \\ -\frac{1}{2}&1&-\frac{1}{2} \end{pmatrix} \end{eqnarray}$ Rückwärtselimination 1 (Diagonalmatrix erzeugen) $\begin{eqnarray} \begin{pmatrix} 3&1&0\\0&\frac{2}{3}&0 \\0&0&-\frac{1}{2} \end{pmatrix} \quad \begin{pmatrix} -5&12&-6 \\ \frac{2}{3} &-2&2 \\ -\frac{1}{2} &1 & -\frac{1}{2} \end{pmatrix} \end{eqnarray}$ Rückwärtselimination 2 (Diagonalmatrix erzeugen) $\begin{eqnarray} \begin{pmatrix} 3&0&0\\0&\frac{2}{3}&0 \\0&0&-\frac{1}{2} \end{pmatrix} \quad \begin{pmatrix} -6&15&-9 \\ -\frac{2}{3} &2&-2 \\ -\frac{1}{2} &1 & -\frac{1}{2} \end{pmatrix} \end{eqnarray}$ Skalierung (auf Einheitsmatrix) $\begin{eqnarray} \begin{pmatrix} 1&0&0\\0&1&0 \\0&0&1 \end{pmatrix} \quad \begin{pmatrix} -2&5&-3 \\ 1 &-3&3 \\ 1 &-2 & 1 \end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} A^{-1}=\begin{pmatrix} -2&5&-3 \\ 1 &-3&3 \\ 1 &-2 & 1 \end{pmatrix} \end{eqnarray*}$
21.07.2007, 00:00	tigerbine	Auf diesen Beitrag antworten »
5e. Austauschalgorithmus Auch hier ist die Grundlage ein Algorithmus zur Lösung des LGS Ax=y. Wieder gibt es Varianten mit und ohne Pivotisierung. Die zweite Variante soll anhand eines Beispiels erläutert werden. Mehr kann man dem angehängten PDF-File entnehmen. Die Laufzeit liegt hier bei $\begin{eqnarray} \mathcal O (n^3) \end{eqnarray}$ Wir betrachten also das LGS: $\begin{eqnarray} 1\cdot x_1+ 1 \cdot x_2 + 1 \cdot x_3 = y_1 \end{eqnarray}$ $\begin{eqnarray} 1\cdot x_1+ 2 \cdot x_2 + 3 \cdot x_3 = y_2 \end{eqnarray}$ $\begin{eqnarray} 1\cdot x_1+ 3 \cdot x_2 + 6 \cdot x_3 = y_3 \end{eqnarray}$ Die erste Spalte enthält nur die 1 als Koeffizient, also ist keine Zeilenvertauschung nötig. Man stellt nun die erste Zeile nach $\begin{eqnarray} x_1 \end{eqnarray}$ um: $\begin{eqnarray} x_1 = 1\cdot y_1 -1 \cdot x_2 - 1 \cdot x_3 \end{eqnarray}$ Die erste Zeile wird nun damit getauscht, in den anderen wird $\begin{eqnarray} x_1 \end{eqnarray}$ ersetzt. $\begin{eqnarray} 1\cdot y_1 - 1 \cdot x_2 - 1 \cdot x_3 = x_1 \end{eqnarray}$ $\begin{eqnarray} 1\cdot y_1+ 1 \cdot x_2 + 2 \cdot x_3 = y_2 \end{eqnarray}$ $\begin{eqnarray} 1\cdot y_1+ 2 \cdot x_2 + 5 \cdot x_3 = y_3 \end{eqnarray}$ Nun schaut man in die zweite Spalte, dort ist ein Tausch nötig. $\begin{eqnarray} 1\cdot y_1 - 1 \cdot x_2 - 1 \cdot x_3 = x_1 \end{eqnarray}$ $\begin{eqnarray} 1\cdot y_1+ 2 \cdot x_2 + 5 \cdot x_3 = y_3 \end{eqnarray}$ $\begin{eqnarray} 1\cdot y_1+ 1 \cdot x_2 + 2 \cdot x_3 = y_2 \end{eqnarray}$ Nun stellt man die Zweite Zeile nach $\begin{eqnarray} x_2 \end{eqnarray}$ $\begin{eqnarray} x_2 = -\frac{1}{2} \cdot y_1 + \frac{1}{2}y_3 -\frac{5}{2} \cdot x_3 \end{eqnarray}$ Nun wieder das LGS aufschreiben $\begin{eqnarray} \frac{3}{2} \cdot y_1 -\frac{1}{2} \cdot y_3 + \frac{3}{2} \cdot x_3 = x_1 \end{eqnarray}$ $\begin{eqnarray} -\frac{1}{2} \cdot y_1 + \frac{1}{2}y_3 -\frac{5}{2} \cdot x_3 = x_2 \end{eqnarray}$ $\begin{eqnarray} \frac{1}{2} \cdot y_1 + \frac{1}{2} \cdot y_3 - \frac{1}{2}\cdot x_3 = y_2 \end{eqnarray}$ Nun noch die Letze Zeile umstellen: $\begin{eqnarray} x_3 = y_1 + y_3 - 2y_2 \end{eqnarray}$ Und im LGS Austauschen und einsetzen: $\begin{eqnarray} 3\cdot y_1 + 1 \cdot y_3 - 3 \cdot y_2 = x_1 \end{eqnarray}$ $\begin{eqnarray} -3 \cdot y_1 - 2\cdot y_3 + 5\cdot y_2 = x_2 \end{eqnarray}$ $\begin{eqnarray} 1\cdot y_1 + 1\cdot y_3 - 2\cdot y_2 = x_3 \end{eqnarray}$ Wollte man nun nur die Lösung des LGS, wäre man an dieser Stellt bereits fertig. Um die Inverse von A zu bestimmen, müssen noch die Spalten 2 und 3 vertauscht werden. Dazu multipliziert man von rechts mit einer Permutationsmatrix. Alles in Matrixschreibweise ergibt dann: $\begin{eqnarray} A^{-1} = \begin{pmatrix} 3&1&-3 \\ -3 &-2&5 \\ 1&1&-2 \end{pmatrix} \begin{pmatrix} 1&0&0 \\ 0&0&1 \\ 0&1&0 \end{pmatrix} = \begin{pmatrix} 3&-3&1 \\ -3 &5&-2 \\ 1&-2 &1 \end{pmatrix} \end{eqnarray}$
21.07.2007, 00:25	tigerbine	Auf diesen Beitrag antworten »
6. Direkte LR-Zerlegung In Abschnitt 3,4 haben wir festgestellt, dass man eine reguläre Matrix A so permutieren kann, dass alle Hauptabschnittsmatrizen A[k] regulär sind. Im folgenden sei A schon in dieser Weise sortiert, so dass es eine eindeutige Zerlegung A = LR gibt. Wie kann man eine Dreieckszerlegung allgemein erzeugt werden? $\begin{eqnarray} L=(l_{jk}),~l_{jk}=0 \text{ für } j < k \end{eqnarray}$ $\begin{eqnarray} R=(r_{kl}),~l_{kl}=0 \text{ für } k > l \end{eqnarray}$ Damit folgt für die Matrix A: $\begin{eqnarray} (a_{jl})=\sum_{k=1}^{\min(j,l)}~l_{jk}r_{rl},\quad j,l=^,1,...,n \end{eqnarray}$ Das ist ein Gleichungssystem mit n² Gleichungen (Einträge von a) und n²+n Unbekannten, nämlich den von Null verschiedenen Einträgen von L und R (hier wird L noch nicht als normiert vorausgesetzt!).
21.07.2007, 00:34	tigerbine	Auf diesen Beitrag antworten »
6a. Verfahren von Doolittle und Crout Um nun eine eindeutige Lösung zu erhalten, könnte man n Einträge, z.B. die Diagonaleinträge von L oder R, vorgeben. Insbesondere sind diese aufgrund der Regularität von A von Null verschieden. Sind nun die Diagonalelemente von L vorgegeben, d.h. es gilt: $\begin{eqnarray} l_{jj} \neq 0 \end{eqnarray}$ ergibt sich das Verfahren von Doolitle $\begin{eqnarray} K:=\min(j,l),~K=1,2,...,n \end{eqnarray}$ $\begin{eqnarray} r_{Kl}=\frac{(a_{Kl}-\sum_{k=1}^{K-1}~l_{Kk}r_{kl})}{l_{KK}},~l \geq K \end{eqnarray}$ $\begin{eqnarray} l_{jK}=\frac{(a_{jK}-\sum_{k=1}^{K-1}~l_{jk}r_{kK})}{r_{KK}},~j \geq K+1 \end{eqnarray}$ Bildlich werden die Einträge dann wie folgt berechnet: $\begin{eqnarray} R=\begin{pmatrix}&+&+&+\\&&&\\&&&\\&&&* \end{pmatrix},~R=\begin{pmatrix}&+&+&+\\&&+&+\\&&&\\&&& \end{pmatrix},R=\begin{pmatrix}&+&+&+\\&&+&+\\&&&+\\&&& \end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} L=\begin{pmatrix}&&&\\+&&&\\+&&&\\+&&& \end{pmatrix},~L=\begin{pmatrix}&&&\\+&&&\\+&+&&\\+&+&& \end{pmatrix},L=\begin{pmatrix}&&&\\+&&&\\+&+&&\\+&+&+& \end{pmatrix} \end{eqnarray}$ Gibt man hingegen die Diagonalelemente von R vor, so berechnet man andersherum mit dem Verfahren von Crout* $\begin{eqnarray} K:=\min(j,l),~K=1,2,...,n \end{eqnarray}$ $\begin{eqnarray} l_{jK}=\frac{(a_{jK}-\sum_{k=1}^{K-1}~l_{jk}r_{kK})}{r_{KK}},~j \geq K \end{eqnarray}$ $\begin{eqnarray} r_{Kl}=\frac{(a_{Kl}-\sum_{k=1}^{K-1}~l_{Kk}r_{kl})}{l_{KK}},~l \geq K+1 \end{eqnarray}$ Bildlich werden die Einträge dann wie folgt berechnet: $\begin{eqnarray} L=\begin{pmatrix}&&&\\+&&&\\+&&&\\+&&&* \end{pmatrix},~L=\begin{pmatrix}&&&\\+&&&\\+&+&&\\+&+&& \end{pmatrix},L=\begin{pmatrix}&&&\\+&&&\\+&+&&\\+&+&+& \end{pmatrix} \end{eqnarray}$ $\begin{eqnarray} R=\begin{pmatrix}&+&+&+\\&&&\\&&&\\&&& \end{pmatrix},~R=\begin{pmatrix}&+&+&+\\&&+&+\\&&&\\&&& \end{pmatrix},R=\begin{pmatrix}&+&+&+\\&&+&+\\&&&+\\&&& \end{pmatrix} \end{eqnarray}$ Anmerkung:* Im Falle K=1 ist die Summe durch 0 zu ersetzen. Wie sollte man nun die Diagonalelemente vorgeben? Am einfachsten ist die Wahl des Wertes 1, da so die Division entfallen kann. Aus vorherigen Überlegungen ist aber auch hier i.A. eine Pivotisierung vorzunehmen.
21.07.2007, 00:42	tigerbine	Auf diesen Beitrag antworten »
7a. Tridiagonalmatrizen Besonders einfach gestalten sich die ebenen genannten Algorithmen im Fall einer Tridiagonalmatrix. Im Falle von Crout ergibt sich dann mit $\begin{eqnarray} r_{kk}=1 \end{eqnarray}$ : $\begin{eqnarray} l_{11} &&= a_{11} \\ r_{j,j+1}&&= \frac{a_{j,j+1}}{l_{jj}} \\l_{j+1,j+1}&&= a_{j+1,j+1}-a_{j+1,j}\cdot r_{j,j+1}\\l_{j+1,j}&&=a_{j+1,j} \end{eqnarray}$ $\begin{eqnarray} j=1,...,n-1 \end{eqnarray}$ Aufwand: $\begin{eqnarray} 2n,~\mathcal O(n) \end{eqnarray}$ Wann ist der Algorithmus wohl definiert? Falls die Matrix A strikt diagonaldominant ist, sind alle Hauptabschnittsmatrizen regulär, insbesondere gilt $\begin{eqnarray} l_{jj} \neq 0 \end{eqnarray}$ . Beweis siehe hier, [Workshop - Spline-Interpolation]
21.07.2007, 01:48	tigerbine	Auf diesen Beitrag antworten »
7b. Allgemein zum Speichern Mit steigendem n kann es zu Speicherproblemen von Matrizen kommen, haben diese doch n² Einträge. So ist es zu überlegen, ob es bei speziellen Matrizen, möglich ist mit weniger Speicherplatz aus zukommen. Dazu betrachtet man Bandmatrizen Die Speicherplatz Ersparnis beruht darauf, dass man nur das Band speichert. Das hat aber auch zur Folge, dass man keine Zeilenvertauschungen mehr durchführen kann. Ist nun bei Bandmatrix A das Gauß'sche Eliminationsverfahren ohne Pivotiserung durchführbar ist, dann sind auch alle reduzierten Matrizen Bandmatrizen desselben Typs und die Faktoren der LR-Zerlegung sind Bandmatrizen vom Typ (p,0), (0,q). Eine dünne Besetzung der Bandmatrix bringt keine Speicherplatzersparnis.
21.07.2007, 01:53	tigerbine	Auf diesen Beitrag antworten »
8a. LDL^T Zerlegungen Sei nun die reguläre Matrix A symmetrisch, d.h. es gilt $\begin{eqnarray} A^T=A \end{eqnarray}$ . Zunächst macht man den verallgemeinerten Ansatz. $\begin{eqnarray} A=LDL^T \end{eqnarray}$ Wobei D eine Diagonalmatrix und L eine normierte untere Dreiecksmatrix ist. Die rechte Seite ergibt offensichtlich eine symmetrische Matrix. Bleibt noch offen, ob es eine solche Zerlegung für jede reguläre symmetrische Matrix A gibt. Zunächst einmal erhält man: $\begin{eqnarray} a_{jl}=\sum_{k=1}^l~d_{kk}\cdot l_{jk} \cdot l_{lk} = \sum_{k=1}^{l-1}~d_{kk} \cdot l_{jk}\cdot l_{lk}+d_{ll}l_{jl} \end{eqnarray}$ $\begin{eqnarray} 1 \leq l \leq j \leq n \end{eqnarray}$ Daraus ergibt sich dann für K=1,2,...,n: $\begin{eqnarray} d_{KK} = a_{KK} - \sum_{k=1}^{K-1}~d_{kk}\cdot l_{Kk}^2 \end{eqnarray}$ $\begin{eqnarray} l_{jK} =\frac{(a_{jK}-\sum_{k=1}^{K-1}~d_{kk}\cdot l_{jk}\cdot l_{Kk})}{d_{KK}},~j \geq K+1 \end{eqnarray}$ Existiert diese Zerlegung für jede reguläre symmetrische Matrix? $\begin{eqnarray} A=\begin{pmatrix} 0&1\\1&0\end{pmatrix} \end{eqnarray}$ Hier ist keine $\begin{eqnarray} A=LDL^T \end{eqnarray}$ Zerlegung möglich, denn $\begin{eqnarray} d_{11} = a_{11} - 0 = 0 \end{eqnarray}$ und damit $\begin{eqnarray} l_{21} \text{ nicht definiert } \end{eqnarray}$
21.07.2007, 02:34	tigerbine	Auf diesen Beitrag antworten »
8b. Cholesky-Zerlegung Ist nun die Matrix A auch noch positiv definit, so sind alle Diagonalelemente von D positiv (betrachte die pos. Definitheit für die Einheitsvektoren). In diesem Fall existiert: $\begin{eqnarray} D^{\frac{1}{2}} = diag(\sqrt{d_{11}},...,\sqrt{d_{nn}}) \end{eqnarray}$ Und es ist: $\begin{eqnarray} A = (LD^{\frac{1}{2}})(LD^{\frac{1}{2}})^T=\hat L\hat L^T \end{eqnarray}$ eine Dreieckszerlegung von A, die Cholesky-Zerlegung. Der Algorithmus folgt aus voherigem: $\begin{eqnarray} l_{KK} = (a_{KK} - \sum_{k=1}^{K-1}~ l_{Kk}^2)^{\frac{1}{2}} \end{eqnarray}$ $\begin{eqnarray} l_{jK} =\frac{(a_{jK}-\sum_{k=1}^{K-1}~ l_{jk}\cdot l_{Kk})}{l_{KK}},~j \geq K+1 \end{eqnarray}$ Weiter ist aufgrund der SPD-Gestalt eine Zerlegung immer ohne Pivotisierung möglich (alle Hauptabschnittsmatrizen sind regulär, ferner sind die Diagonalelemente von A positiv) Beweis Sei $\begin{eqnarray} e_i \end{eqnarray}$ der i-te Einheitsvektor des $\begin{eqnarray} \mathbb R^{n} \end{eqnarray}$ . Aus der positiven Definitheit von A folgt: $\begin{eqnarray} a_{ii} = e_i^TAe_i >0 \end{eqnarray}$ Die Symmetrie der der A[k] folgt direkt aus der von A. Für $\begin{eqnarray} k \in \{1,...,n\} \end{eqnarray}$ sei $\begin{eqnarray} y \in \mathbb R^k,~y \neq 0 \end{eqnarray}$ beliebig gegeben. Mit $\begin{eqnarray} x:=(y^T,0,...,0)^T \in \mathbb R^n \end{eqnarray}$ . Es folgt dann $\begin{eqnarray} y^TA[k]y = x^Tx > 0 \end{eqnarray}$ somit dass A[k] positiv definit ist. Aufwand $\begin{eqnarray} \approx \frac{n^3}{6} + n,~\mathcal O (n^3) \end{eqnarray}$ Leider weichen hier die Angaben über die Laufzeiten von Buch zu Buch ab. Festzuhalten bleibt, das für SPD Matrizen die Cholesky Zerlegung halb so aufwendig ist wie die Gaußelimination und man ohne Pivotisierung auskommt. Der Cholesky-Algorithmus ist die effizienteste allgemeine Testmethode auf positive Definitheit. Man kann auf das Rechnen mit den Quadratwurzeln verzichten, indem man unter Berücksichtigung der Symmetrie den Gaußalgorithmus wie gewohnt durchführt und nur die Pivots in einer Diagonalmatrix (als Vektor) speichert. Dies gibt dann die Eingang erwähnte $\begin{eqnarray} LDL^T \end{eqnarray}$ Zerlegung.
01.01.2009, 17:24	tigerbine	Auf diesen Beitrag antworten »
9. QR-Zerlegung In diesem Workshop soll sich mir der QR-Zerlegung auseinander gesetzt werden. D.h. die Matrix A soll Produkt einer orthogonalen Matrix Q und einer oberen Dreiecksmatrix R dargestellt werden. Als Anwendung soll dann nur das Lösen Linearer Gleichungssysteme im Vergleich zur LR-Zerlegung betrachtet werden. Daher sei A im folgenden auch als regulär angenommen. Die Bedeutung bei der Lösung linearer Ausgleichsprobleme und Eigenwertprobleme wird in eigenen Workshops behandelt.
01.01.2009, 18:08	tigerbine	Auf diesen Beitrag antworten »
9a. Eindeutigkeit der QR-Zerlegung (bis auf Diagonalmatrix) Gegeben sind n² Einträge der Matrix A und zu bestimmen sind [n² + n(n+1)/2] Einträge der Matrizen Q und R. So ohne weiteres ist also (im Falle der Existenz) nicht mit einer Eindeutigen Zerlegung A=QR zu rechnen. Die Spalten der Matrix Q müssen orthonomiert sein. Sollte nun eine QR-Zerlegung existieren, so ist sie im Wesentlichen eindeutig. Betrachtet man also 2 QR-Zerlegungen, so existiert eine Diagonalmatrix D (mit Werten +/- 1), so dass gilt: $\begin{eqnarray} A = Q_1 R_1,~A = Q_2R_2,\quad Q_1 = Q_2D,~R_2 = DR_1 \end{eqnarray}$ Beweis: Es nun nach Voraussetzung $\begin{eqnarray} Q_1 R_1 =Q_2 R_2 \end{eqnarray}$ . Da es sich bei den Qs um orthogonale Matrizen handelt, folgt sofort: $\begin{eqnarray} Q_1^T Q_1 = I_n, \quad \hat Q_2^TQ_2 = I_n \end{eqnarray}$ Damit erhält man dann aus der ersten Gleichheit: $\begin{eqnarray} Q_2^T (Q_1 R_1) &&= Q_2^T (Q_2 R_2) \\ (Q_2^TQ_1) R_1 &&= R_2 \end{eqnarray}$ und andersherum folgt: $\begin{eqnarray} (Q_1^TQ_2)R_2 &&= R_1 \end{eqnarray}$ Es sind nun die Rs regulär, sonst wäre Rang(A) <n, und somit existieren die Inversen: $\begin{eqnarray} Q_2^T Q_1 &&= R_2 R_1^{-1} \end{eqnarray}$ $\begin{eqnarray} Q_1^T Q_2 &&= R_1 R_2^{-1} \end{eqnarray}$ Nun definiert man $\begin{eqnarray} D:=R_2 R_1^{-1} \end{eqnarray}$ und weist nach, dass es sich um eine orthogonale Diagonalmatrix handelt. Mit dem Wissen über Dreiecksmatrizen folgt, dass auch die Inversen obere Dreiecksmatrizen sind. Es ist also auch D eine obere Dreiecksmatrix und ihre Transponierte somit eine untere Dreiecksmatrix. Es gilt weiter: $\begin{eqnarray} D^T = (R_2R_1^{-1})^T = (Q_2^T Q_1)^T = Q_1^T Q_2 = R_1 R_2^{-1} \end{eqnarray}$ Somit muss D eine Diagonalmatrix sein. Nun folgt weiter $\begin{eqnarray} D^T = R_1 R_2^{-1} = (R_2 R_1^{-1})^{-1} = D^{-1} \end{eqnarray}$ Womit nun auch die Orthogonalität von D gezeigt wurde, d.h. die Diagonalwerte sind +/- 1. Mit $\begin{eqnarray} Q_1 = Q_2 (R_2 R_1^{-1}) = Q_2D \end{eqnarray}$ folgt dann die Behauptung.
01.01.2009, 20:40	tigerbine	Auf diesen Beitrag antworten »
9b. Lösen eines reg. LGS Ax=b Sollte einer Zerlegung existieren, so ist Q einfach durch Transponieren invertierbar und man kann das LGS am Ende durch Rückwärtssubstitution lösen. $\begin{eqnarray} Ax&&=b \\ QRx&&=b \\Rx&&=Q^Tb \end{eqnarray}$
01.01.2009, 20:59	tigerbine	Auf diesen Beitrag antworten »
9c. Existenz von QR-Zerlegungen Mit der Angabe eines durchführbaren Algorithmus, der zu einer Matrix A eine Zerlegung QR=A liefert dürfte die Frage nach der Existenz einer solchen Zerlegung beantwortet sein.
01.01.2009, 21:18	tigerbine	Auf diesen Beitrag antworten »
9d. Gram-Schmidt-Orthogonalisierung Gesucht ist nun die QR-Zerlegung. Zunächst notiert man: $\begin{eqnarray} A = Q R \Leftrightarrow \begin{pmatrix} \vdots & & \vdots \\ a^1& \hdots & a^n \\ \vdots&&\vdots \end{pmatrix} =\begin{pmatrix} \vdots & & \vdots \\ q^1& \hdots & q^n \\ \vdots&&\vdots \end{pmatrix}\begin{pmatrix} r_{11} & &r_{1n} \\ & \ddots & \vdots \\ &&r_{nn} \end{pmatrix} \end{eqnarray}$ Die Potenzen Kennzeichen hier die Spaltenvektoren ausgeschrieben bedeutet das: $\begin{eqnarray} a^1=r_{11}\cdot q^1 \end{eqnarray}$ $\begin{eqnarray} a^2=r_{12}\cdot q^1 + r_{22} \cdot q^2 \end{eqnarray}$ $\begin{eqnarray} a^3=r_{13}\cdot q^1 + r_{23} \cdot q^2 + r_{33} \cdot q^3 \end{eqnarray}$ $\begin{eqnarray} \vdots \end{eqnarray}$ $\begin{eqnarray} a^n=r_{1n}\cdot q^1 + r_{2n}\cdot q^2 + \hdots + r_{nn} \cdot q^n \end{eqnarray}$ allgemein: $\begin{eqnarray} a^j = \sum_{i=1}^j~r_{ij} \cdot q^i \quad j=1,...,n \qquad () \end{eqnarray}$ Damit erhält man nun eine Berechnungsvorschrift für die orthonormalen Vektoren "q" und die zugehörigen Koeffizienten "r". Für j=1 ergibt sich: $\begin{eqnarray} \|\|a^1\|\|_2 = \|r_{11}\| \cdot \|\|q^1\|\|_2 \end{eqnarray}$ Da Q orthogonal ist, gilt $\begin{eqnarray} \|\|q^{j}\|\|_2 = 1 ,~j=1,...,n \end{eqnarray}$ $\begin{eqnarray} \Rightarrow r_{11} = \pm \|\|a^1\|\|_2,\quad q^1 = \frac{1}{r_{11}}\cdot a^1 \end{eqnarray}$ Um nun eine eindeutige Lösung zu erhalten, sollen die Diagonalelemente von R positiv sein. Weiter erhält man durch Umstellen zunächst: $\begin{eqnarray} r_{22} = \|\|a^2-r_{12}\cdot q^1\|\|_2 \quad \text{Normierung von } q^2 \end{eqnarray}$ $\begin{eqnarray} q^2 = \frac{1}{r_{22}}\cdot (a^2-r_{12}\cdot q^1)\quad \text{Umstellen von ()} \end{eqnarray}$ Nun die Orthogonalitätsbedingung: $\begin{eqnarray} 0 \stackrel{!}=(q^1)^Tq^2 \end{eqnarray}$ $\begin{eqnarray} 0 \stackrel{!}=(q^1)^T\frac{1}{r_{22}} \cdot (a^2-r_{12} \cdot q^1) \end{eqnarray}$ $\begin{eqnarray} 0 \stackrel{!}=\frac{1}{r_{22}} \cdot [(q^1)^Ta^2-r_{12} \cdot (q^1)^Tq^1] \end{eqnarray}$ $\begin{eqnarray} 0 \stackrel{!}=\frac{1}{r_{22}} \cdot [(q^1)^Ta^2-r_{12} ] \end{eqnarray}$ $\begin{eqnarray} r_{12} = (q^1)^Ta^2 \end{eqnarray}$ Allgemein ergibt sich dann: $\begin{eqnarray} r_{jj} = \|\|a^j-\sum_{i=1}^{j-1}~r_{ij} \cdot q^{i}\|\|_2 \end{eqnarray}$ $\begin{eqnarray} q^j:=\frac{1}{r_{jj}} \cdot \Bigl( a^j - \sum_{k=1}^{j-1}~r_{kj} \cdot q^k \Bigr) \end{eqnarray}$ $\begin{eqnarray} r_{ij} = (q^{i})^Ta^j~\quad \forall ~ i=1,...,j-1 \end{eqnarray}$ Mit dieser Vorschrift (Gram-Schmidt-Verfahren) ist die Existenz einer red. QR-Zerlegung konstruktiv (fast) bewiesen worden. Es ist noch zu zeigen, dass das Verfahren wohldefiniert ist. Dies ist dann der Fall, wenn für die Nenner gilt $\begin{eqnarray} r_{jj} \neq 0 \end{eqnarray}$ . Wäre nun $\begin{eqnarray} r_{jj} = 0 \end{eqnarray}$ für ein kleinstes $\begin{eqnarray} j \in \{1,...,n \} \end{eqnarray}$ , so folgt aus () $\begin{eqnarray} a^j \in \text{span} \{q^1,...,q^{j-1}\} = \text{span} \{a^1,...,a^{j-1}\} \end{eqnarray*}$ Dies steht aber im Widerspruch zu Rang(A) = n. Numerisch ist das Verfahren nicht stabil. Ein modifizierter Algorithmus befindet sich im Anhang.
02.01.2009, 00:41	tigerbine	Auf diesen Beitrag antworten »
9e. Householder-Spiegelungen Bei diesem Ansatz geht man ähnlich der Berechnung der LR-Zerlegung vor. Durch Linksmultiplikation mit Orthogonalen Matrizen soll A auf die Form einer Dreiecksmatrix gebracht werden. $\begin{eqnarray} Q^{-1}A = R \Leftrightarrow Q^{T}A = R \Leftrightarrow A = QR \end{eqnarray}$ Dabei annulliert $\begin{eqnarray} Q_j \end{eqnarray}$ die Eintrage der j-ten Spalte unterhalb der Diagonalen. $\begin{eqnarray} Q^T = Q_l \cdot .... \cdot Q_1,\quad l:=\min\{m-1,n\} \end{eqnarray}$ Householder-Spiegelung (Definition) $\begin{eqnarray} H:=I_n - 2\cdot \underbrace{uu^T}_{Dyade}\quad \text{mit } u \in \mathbb R^n,~\|\|u\|\|_2=1 \end{eqnarray}$ Konstruktion von u x sei nun der Teil des Spaltenvektors von "A", den es zu annullieren gilt, zzgl. des Diagonalelements. $\begin{eqnarray} u:=x + \text{sign}(x_1)\cdot \|\|x\|\|_2 \cdot e_1 \quad \text{(Bedingungen für die zu erzeugenden Nullen)} \end{eqnarray}$ $\begin{eqnarray} \beta:=\frac{2}{u^Tu} \quad \text{(Normierung)} \end{eqnarray}$ $\begin{eqnarray} H=I-\beta\cdot uu^T \end{eqnarray}$
02.01.2009, 00:43	tigerbine	Auf diesen Beitrag antworten »
9f. Givens-Rotatioinen Bei diesem Ansatz geht man ähnlich der Berechnung der LR-Zerlegung vor. Durch Linksmultiplikation mit Orthogonalen Matrizen soll A auf die Form einer Dreiecksmatrix gebracht werden. $\begin{eqnarray} Q^{-1}A = R \Leftrightarrow Q^{T}A = R \Leftrightarrow A = QR \end{eqnarray}$ Dabei annulliert $\begin{eqnarray} G_{ij} \end{eqnarray}$ den Eintrag an der Stelle (i,j). $\begin{eqnarray} Q^T = G_{ml} \cdot .... \cdot G_{21},\quad l:=\min\{m-1,n\} \end{eqnarray}$ Givens-Rotation (Definition) $\begin{eqnarray} G_{ij}:=\begin{pmatrix}I_{j-1} & & & \\ &c&0&\hdots &0&s \\ &0&1&&0&0\\& \vdots & \vdots & \ddots &\vdots &\vdots \\&0&0&&1&0\\ &-s&0&\hdots &0&c \\ &&&&&&I_{n-i} \end{pmatrix} \end{eqnarray}$ Konstruktion von c und s $\begin{eqnarray} p:=\pm \sqrt{a_{jj}^2+a_{ij}^2} \end{eqnarray}$ $\begin{eqnarray} c:=\frac{a_{jj}}{p},~s:=\frac{a_{ij}}{p} \end{eqnarray}$

12 nächste Seite »

Neue Frage »

Antworten »

[WS] Lineare Gleichungssysteme 2 - direkte Verfahren

Verwandte Themen