LUP-Zerlegung
3. Umformung in eine Iteration
Dieser Teil folgt in Kürze.
Wenn wir die Matrix A so in das Produkt L·U zerlegen,
dass L eine untere Dreiecksmatrix mit Einsen auf der Hauptdiagonalen ist
und U eine obere Dreiecksmatrix, können wir die relevaten Informationen
in eine einzige Matrix packen:
\def\pink{\bbox[#ffe0e0,2pt]}
\def\green{\bbox[#c0ffc0,2pt]}
\pmatrix{\green{u_{11}} & \green{u_{12}} & \green{u_{13}} & \ldots & \green{u_{1n}} \cr
\pink{l_{21}} & \green{u_{22}} & \green{u_{23}} & \ldots & \green{u_{2n}} \cr
\pink{l_{31}} & \pink{l_{32}} & \green{u_{33}} & \ldots & \green{u_{2n}} \cr
\vdots & \vdots & \vdots & \ddots & \vdots \cr
\pink{l_{n1}} & \pink{l_{n2}} & \pink{l_{n3}} & \ldots & \green{u_{nn}}}
Daraus lassen sich L und U gewinnen:
- L, indem die Diagonalelemente (u11
bis unn) durch 1 und die
übrigen grün unterlegten Elemente durch Null ersetzt werden,
- U, indem alle rötlich unterlegten
Elemente durch Null ersetzt werden.
Wir entwerfen nun einen Algorithmus, der die Ausgangsmatrix A
direkt in diese Kompaktdarstellung umformt.
Im letzten Abschnitt haben wir die Matrix so zerlegt:
\def\yellow{\bbox[#ffffa0,2pt]}
\def\pink{\bbox[#ffe0e0,2pt]}
\def\green{\bbox[#c0ffc0,2pt]}
\def\blue{\bbox[#e8f0ff,2pt;border:1px #8090ff solid]}
\def\white{\bbox[border:1px #8090ff solid]}
A = \pmatrix{a_{11} & \textbf w \cr
\textbf v & A'}
= \pmatrix{ 1 & \textbf 0 \cr
\pink{\frac{1}{a_{11}} \textbf{v}} & I_{n-1}}
\pmatrix{\green {a_{11}} & \green{\textbf w} \cr
\textbf 0 & \yellow{A' - \frac{1}{a_{11}} \textbf{v w}}}
Wenn wir die relevanten Informationen zusammenpacken, erhalten wir:
\def\yellow{\bbox[#ffffa0,2pt]}
\def\pink{\bbox[#ffe0e0,2pt]}
\def\green{\bbox[#c0ffc0,2pt]}
\pmatrix{\green {a_{11}} & \green{\textbf w} \cr
\pink{\frac{1}{a_{11}} \textbf{v}} & \yellow{A' - \frac{1}{a_{11}} \textbf{v w}}}
Im den rötlichen und grünen Bereichen stehen schon Teile der Dreicksmatrizen,
der gelb unterlegte Bereich muss weiter zerlegt werden.
Dies bekommen wir durch folgende Operationen aus der Ausgangsmatrix
(machen Sie sich dazu bitte klar, wie das Produkt
vw gebildet wird):
- Die erste Zeile bleibt unverändert
- Der Rest der ersten Spalte wird durch
a11 dividiert
- Im Rest der Matrix wird von jedem Element
das Produkt aus dem links davon stehenden Element der ersten Spalte
(das schon durch a11 dividiert ist!)
und dem darüber stehenden Element der ersten Zeile subtrahiert.
In mathGUIde (Matrizen werden hier ab Null indiziert!) sieht das so aus:
for i in fromTo(1, n-1):
A[i,0] /= A[0,0]
for i in fromTo(1, n-1):
for j in fromTo(1, n-1):
A[i,j] -= A[i,0] * A[0,j]
Jetzt ist es nur noch ein kleiner Schritt, dieses Verfahren
auf den Rest der Matrix zu erweitern.
Dazu führen wir eine äußere Schleifenvariable k
für
den linken (und oberen) Index der zu bearbeitenden
Restmatrix ein:
def luIter(A):
n = A.height()
for k in fromTo(0, n-2):
assert A[k,k] != 0
for i in fromTo(k+1, n-1):
A[i,k] /= A[k,k]
for j in fromTo(k+1,n-1):
A[i,j] -= A[i,k] * A[k,j]
Wir testen die Funktion am Beispiel aus dem ersten Abschnitt:
\def\yellow{\bbox[#ffffa0,2pt]}
\def\pink{\bbox[#ffe0e0,2pt]}
\def\green{\bbox[#c0ffc0,2pt]}
\yellow{\pmatrix{1 & 2 & 3 \cr
2 & 9 & 10 \cr
3 & 26 & 24 }}
= \pink{\pmatrix{1 & 0 & 0 \cr
2 & 1 & 0 \cr
3 & 4 & 1 }}
\green{\pmatrix{1 & 2 & 3 \cr
0 & 5 & 4 \cr
0 & 0 & -1 }}
A = Matrix.fromString("1, 2, 3; 2, 9, 10; 3, 26, 24")
luIter(A)
print(A)
Bitte überzeugen Sie sich, dass nun in A die gesuchte linke und die rechte
Dreiecksmatrix stecken!
Fortsetzung: 4. Beseitigung des Schönheitsfehlers: LUP-Zerlegung