LR Zerlegung - Max Planck Societypeople.mpi-inf.mpg.de/~msagralo/LR_Zerlegung.pdfAx+b = x ()(Id3 A)x...
48
LR Zerlegung Michael Sagraloff
Transcript of LR Zerlegung - Max Planck Societypeople.mpi-inf.mpg.de/~msagralo/LR_Zerlegung.pdfAx+b = x ()(Id3 A)x...
LR ZerlegungMichael Sagraloff
Beispiel eines linearen Gleichungssystems in derÖkonomie (Input-Output Analyse)
Wir nehmen an, dass es 3 Güter G1, G2, und G3 gibt. Dannentspricht der Eintrag ai,j der sogenannten Input-Output Matrix
A =
(a1,1 a1,2 a1,3
a2,1 a2,2 a2,3
a3,1 a3,2 a3,4
):=
(0 0.2 0.4
0.2 0 00.4 0.4 0
)der Anzahl der Einheiten des Gutes Gi , die zur Produktioneiner Einheit des Gutes Gj benötigt werden.Der Endbedarf des Marktes an Gut Gi sei bi , mit
b := (b1,b2,b3)t := (24,64,36)t .
Wieviele Einheiten xi von Gi müssen insgesamt hergestelltwerden?Die Lösung ergibt sich aus dem Gleichungssystem
A · x + b = x⇐⇒ (Id3−A) · x = b
LR Zerlegung Michael Sagraloff 15.06.2016 2
Beispiel eines linearen Gleichungssystems in derÖkonomie (Input-Output Analyse)
Wir nehmen an, dass es 3 Güter G1, G2, und G3 gibt. Dannentspricht der Eintrag ai,j der sogenannten Input-Output Matrix
A =
(a1,1 a1,2 a1,3
a2,1 a2,2 a2,3
a3,1 a3,2 a3,4
):=
(0 0.2 0.4
0.2 0 00.4 0.4 0
)der Anzahl der Einheiten des Gutes Gi , die zur Produktioneiner Einheit des Gutes Gj benötigt werden.Der Endbedarf des Marktes an Gut Gi sei bi , mit
b := (b1,b2,b3)t := (24,64,36)t .
Wieviele Einheiten xi von Gi müssen insgesamt hergestelltwerden?Die Lösung ergibt sich aus dem Gleichungssystem
A · x + b = x⇐⇒ (Id3−A) · x = b
LR Zerlegung Michael Sagraloff 15.06.2016 2
Input-Output Analyse
(1 −0.2 −0.4−0.2 1 0−0.4 −0.4 1
)·
(x1
x2
x3
)=
(246436
)⇔
x1 − 0.2x2 − 0.4x3 = 24−0.2x1 + x2 = 64
−0.4x1 − 0.4x2 + x3 = 36
Elimination von Variablen durch Abziehen entsprechender Viel-fache einer Gleichung von einer anderen (Gauß Elimination): 1 −0.2 −0.4
−0.2 1 0−0.4 −0.4 1
· x1
x2
x3
=
246436
⇔
1 −0.2 −0.40 (−0.2) 0.96 −0.080 (−0.4) −0.48 0.84
· x1
x2
x3
=
2468.845.6
⇔
1 −0.2 −0.40 (−0.2) 0.96 −0.080 (−0.4) 0 (−0.5) 0.8
· x1
x2
x3
=
2468.880
Rückwärtseinsetzen: x3 = 100, x2 = − (68.8+0.08·x3)
0.96 = 80, undx1 = 24 + 0.2 · x2 + 0.4 · x3 = 80.
LR Zerlegung Michael Sagraloff 15.06.2016 3
Input-Output Analyse
(1 −0.2 −0.4−0.2 1 0−0.4 −0.4 1
)·
(x1
x2
x3
)=
(246436
)⇔
x1 − 0.2x2 − 0.4x3 = 24−0.2x1 + x2 = 64
−0.4x1 − 0.4x2 + x3 = 36
Elimination von Variablen durch Abziehen entsprechender Viel-fache einer Gleichung von einer anderen (Gauß Elimination): 1 −0.2 −0.4
−0.2 1 0−0.4 −0.4 1
· x1
x2
x3
=
246436
⇔
1 −0.2 −0.40 (−0.2) 0.96 −0.080 (−0.4) −0.48 0.84
· x1
x2
x3
=
2468.845.6
⇔
1 −0.2 −0.40 (−0.2) 0.96 −0.080 (−0.4) 0 (−0.5) 0.8
· x1
x2
x3
=
2468.880
Rückwärtseinsetzen: x3 = 100, x2 = − (68.8+0.08·x3)
0.96 = 80, undx1 = 24 + 0.2 · x2 + 0.4 · x3 = 80.
LR Zerlegung Michael Sagraloff 15.06.2016 3
Input-Output Analyse
(1 −0.2 −0.4−0.2 1 0−0.4 −0.4 1
)·
(x1
x2
x3
)=
(246436
)⇔
x1 − 0.2x2 − 0.4x3 = 24−0.2x1 + x2 = 64
−0.4x1 − 0.4x2 + x3 = 36
Elimination von Variablen durch Abziehen entsprechender Viel-fache einer Gleichung von einer anderen (Gauß Elimination): 1 −0.2 −0.4
−0.2 1 0−0.4 −0.4 1
· x1
x2
x3
=
246436
⇔
1 −0.2 −0.40 (−0.2) 0.96 −0.080 (−0.4) −0.48 0.84
· x1
x2
x3
=
2468.845.6
⇔
1 −0.2 −0.40 (−0.2) 0.96 −0.080 (−0.4) 0 (−0.5) 0.8
· x1
x2
x3
=
2468.880
Rückwärtseinsetzen: x3 = 100, x2 = − (68.8+0.08·x3)
0.96 = 80, undx1 = 24 + 0.2 · x2 + 0.4 · x3 = 80.
LR Zerlegung Michael Sagraloff 15.06.2016 3
Input-Output Analyse
(1 −0.2 −0.4−0.2 1 0−0.4 −0.4 1
)·
(x1
x2
x3
)=
(246436
)⇔
x1 − 0.2x2 − 0.4x3 = 24−0.2x1 + x2 = 64
−0.4x1 − 0.4x2 + x3 = 36
Elimination von Variablen durch Abziehen entsprechender Viel-fache einer Gleichung von einer anderen (Gauß Elimination): 1 −0.2 −0.4
−0.2 1 0−0.4 −0.4 1
· x1
x2
x3
=
246436
⇔
1 −0.2 −0.40 (−0.2) 0.96 −0.080 (−0.4) −0.48 0.84
· x1
x2
x3
=
2468.845.6
⇔
1 −0.2 −0.40 (−0.2) 0.96 −0.080 (−0.4) 0 (−0.5) 0.8
· x1
x2
x3
=
2468.880
Rückwärtseinsetzen: x3 = 100, x2 = − (68.8+0.08·x3)
0.96 = 80, undx1 = 24 + 0.2 · x2 + 0.4 · x3 = 80.
LR Zerlegung Michael Sagraloff 15.06.2016 3
Input-Output Analyse
Man rechnet leicht nach, dass
A =
1 0 0−0.2 1 0−0.4 −0.5 1
︸ ︷︷ ︸
L
· 1 −0.2 −0.4
0 0.96 −0.080 0 0.8
︸ ︷︷ ︸
R
,
d.h. A lässt sich als Produkt zweier Dreiecksmatrizen L und Rschreiben. Das ist kein Zufall!Beachte: Die Matrix A ist üblicherweise viel größer!Beispiel für das Auftreten eines sehr großen linearenGleichungssystems in der Medizin: Computertomographie
LR Zerlegung Michael Sagraloff 15.06.2016 4
Input-Output Analyse
Man rechnet leicht nach, dass
A =
1 0 0−0.2 1 0−0.4 −0.5 1
︸ ︷︷ ︸
L
· 1 −0.2 −0.4
0 0.96 −0.080 0 0.8
︸ ︷︷ ︸
R
,
d.h. A lässt sich als Produkt zweier Dreiecksmatrizen L und Rschreiben. Das ist kein Zufall!Beachte: Die Matrix A ist üblicherweise viel größer!Beispiel für das Auftreten eines sehr großen linearenGleichungssystems in der Medizin: Computertomographie
LR Zerlegung Michael Sagraloff 15.06.2016 4
LR Zerlegung (bzw. LU Zerlegung)
Problemstellung 1Gegeben eine n × n Matrix
A =
a1,1 a1,2 · · · a1,n
a2,1 a2,2 · · · a2,n...
......
...an,1 an,2 · · · an,n
∈ Kn×n
über einem Körper K (üblicherweise K = R oder K = C),bestimme Links- bzw. Rechtsdreiecksmatrizen
L =
l1,1 0 · · · 0l2,1 l2,2 · · · 0...
......
...ln,1 ln,2 · · · ln,n
und R =
r1,1 r1,2 · · · r1,n
0 r2,2 · · · r2,n...
......
...0 0 · · · rn,n
mit Einträgen li,j ∈ K und ri,j ∈ K, so dass A = L · R.
LR Zerlegung Michael Sagraloff 15.06.2016 5
Sinnvolle Fragestellung?
(a) Eindeutigkeit? 3 2 16 6 39 10 3
︸ ︷︷ ︸
A
=
1 0 02 1 03 2 1
︸ ︷︷ ︸
L
·
3 2 10 2 10 0 1
︸ ︷︷ ︸
R
=
1 0 02 2 03 4 1
︸ ︷︷ ︸
L′
·
3 2 10 1 1
20 0 −2
︸ ︷︷ ︸
R′
(b) Existenz?
A =
(3 2 16 4 39 10 3
)=
(l1,1 0 0l2,1 l2,2 0l3,1 l3,2 l3,3
)·
(r1,1 r1,2 r1,3
0 r2,2 r2,3
0 0 r3,3
)
Es gilt(
3 26 4
)=
(l1,1 0l2,1 l2,2
)·(
r1,1 r1,2
0 r2,2
), und daher ist
einer der beiden 2× 2-Matrizen auf der rechten Seite derGleichung singulär. Somit gilt auch det A = 0
LR Zerlegung Michael Sagraloff 15.06.2016 6
Sinnvolle Fragestellung?
(a) Eindeutigkeit? 3 2 16 6 39 10 3
︸ ︷︷ ︸
A
=
1 0 02 1 03 2 1
︸ ︷︷ ︸
L
·
3 2 10 2 10 0 1
︸ ︷︷ ︸
R
=
1 0 02 2 03 4 1
︸ ︷︷ ︸
L′
·
3 2 10 1 1
20 0 −2
︸ ︷︷ ︸
R′
(b) Existenz?
A =
(3 2 16 4 39 10 3
)=
(l1,1 0 0l2,1 l2,2 0l3,1 l3,2 l3,3
)·
(r1,1 r1,2 r1,3
0 r2,2 r2,3
0 0 r3,3
)
Es gilt(
3 26 4
)=
(l1,1 0l2,1 l2,2
)·(
r1,1 r1,2
0 r2,2
), und daher ist
einer der beiden 2× 2-Matrizen auf der rechten Seite derGleichung singulär. Somit gilt auch det A = 0
LR Zerlegung Michael Sagraloff 15.06.2016 6
Sinnvolle Fragestellung?
(a) Eindeutigkeit? 3 2 16 6 39 10 3
︸ ︷︷ ︸
A
=
1 0 02 1 03 2 1
︸ ︷︷ ︸
L
·
3 2 10 2 10 0 1
︸ ︷︷ ︸
R
=
1 0 02 2 03 4 1
︸ ︷︷ ︸
L′
·
3 2 10 1 1
20 0 −2
︸ ︷︷ ︸
R′
(b) Existenz?
A =
(3 2 16 4 39 10 3
)=
(l1,1 0 0l2,1 l2,2 0l3,1 l3,2 l3,3
)·
(r1,1 r1,2 r1,3
0 r2,2 r2,3
0 0 r3,3
)
Es gilt(
3 26 4
)=
(l1,1 0l2,1 l2,2
)·(
r1,1 r1,2
0 r2,2
), und daher ist
einer der beiden 2× 2-Matrizen auf der rechten Seite derGleichung singulär. Somit gilt auch det A = 0
LR Zerlegung Michael Sagraloff 15.06.2016 6
Sinnvolle Fragestellung?
(a) Eindeutigkeit? 3 2 16 6 39 10 3
︸ ︷︷ ︸
A
=
1 0 02 1 03 2 1
︸ ︷︷ ︸
L
·
3 2 10 2 10 0 1
︸ ︷︷ ︸
R
=
1 0 02 2 03 4 1
︸ ︷︷ ︸
L′
·
3 2 10 1 1
20 0 −2
︸ ︷︷ ︸
R′
(b) Existenz?
A =
(3 2 16 4 39 10 3
)=
(l1,1 0 0l2,1 l2,2 0l3,1 l3,2 l3,3
)·
(r1,1 r1,2 r1,3
0 r2,2 r2,3
0 0 r3,3
)Es gilt
(3 26 4
)=
(l1,1 0l2,1 l2,2
)·(
r1,1 r1,2
0 r2,2
), und daher ist
einer der beiden 2× 2-Matrizen auf der rechten Seite derGleichung singulär. Somit gilt auch det A = 0
LR Zerlegung Michael Sagraloff 15.06.2016 6
Was ist hinsichtlich der Existenz schiefgegangen?
Auch wenn A regulär ist, kann eine der i × i-Untermatrizen
A[i] :=
a1,1 a2,2 · · · a1,i
a2,1 a2,2 · · · a2,i...
......
...ai,1 ai,2 · · · ai,i
singulär sein, d.h. der i -te Hauptminor det A[i] verschwindet.Dann existiert keine LR-Zerlegung wie in Problemstellung 1.
geeignete Umsortierung der Zeilen: det A[i] 6= 0 für alle i .
HauptsatzFalls A regulär ist, gibt es eine (Zeilen-) PermutationsmatrixP ∈ Kn×n und eindeutig bestimmte Matrizen L und R wie inProblemstellung 1 angegeben, so dass
li,i = 1 für alle i = 1, . . . ,n, undP · A = L · R.
LR Zerlegung Michael Sagraloff 15.06.2016 7
Was ist hinsichtlich der Existenz schiefgegangen?
Auch wenn A regulär ist, kann eine der i × i-Untermatrizen
A[i] :=
a1,1 a2,2 · · · a1,i
a2,1 a2,2 · · · a2,i...
......
...ai,1 ai,2 · · · ai,i
singulär sein, d.h. der i -te Hauptminor det A[i] verschwindet.Dann existiert keine LR-Zerlegung wie in Problemstellung 1.geeignete Umsortierung der Zeilen: det A[i] 6= 0 für alle i .
HauptsatzFalls A regulär ist, gibt es eine (Zeilen-) PermutationsmatrixP ∈ Kn×n und eindeutig bestimmte Matrizen L und R wie inProblemstellung 1 angegeben, so dass
li,i = 1 für alle i = 1, . . . ,n, undP · A = L · R.
LR Zerlegung Michael Sagraloff 15.06.2016 7
Was ist hinsichtlich der Existenz schiefgegangen?
Auch wenn A regulär ist, kann eine der i × i-Untermatrizen
A[i] :=
a1,1 a2,2 · · · a1,i
a2,1 a2,2 · · · a2,i...
......
...ai,1 ai,2 · · · ai,i
singulär sein, d.h. der i -te Hauptminor det A[i] verschwindet.Dann existiert keine LR-Zerlegung wie in Problemstellung 1.geeignete Umsortierung der Zeilen: det A[i] 6= 0 für alle i .
HauptsatzFalls A regulär ist, gibt es eine (Zeilen-) PermutationsmatrixP ∈ Kn×n und eindeutig bestimmte Matrizen L und R wie inProblemstellung 1 angegeben, so dass
li,i = 1 für alle i = 1, . . . ,n, undP · A = L · R.
LR Zerlegung Michael Sagraloff 15.06.2016 7
Hauptsatz
Noch einmal zurück zum Beispiel: 0 0 10 1 01 0 0
︸ ︷︷ ︸
P
·
(3 2 16 4 39 10 3
)=
(1 0 03 1 02 0 1
)·
(3 2 10 4 00 0 1
)
Der Beweis des Hauptsatzes ist konstruktiv und liefert einenAlgorithmus zur Bestimmung von P, L, und R.
(A) Im ersten Schritt behandeln wir nur den Spezialfall, in dem alleHauptminoren von A nicht verschwinden. In diesem Fallkönnen wir P = Idn wählen.
(B) Im zweiten Schritt zeigen wir, wie der allgemeine Fall durcheine entsprechende Umsortierung der Zeilen auf (1) zurück-geführt werden kann; Verfahren wird nur anhand eines Bsp.verdeutlicht; Beweis findet sich im Anhang.
LR Zerlegung Michael Sagraloff 15.06.2016 8
Warum will man eine LR-Zerlegung berechnen?
Berechnung der Determinante einer Matrix:
det(A) = det(P)−1 · det(L) · det(R) =
wobei v die Anzahl der Zeilenwechsel bezeichnet, die beiAnwendung mit P herbeigeführt werden.
Lösen eines linearen Gleichungssystems: Bestimme x ∈ Kn,so dass
A · x = b ,
wobei b ∈ Kn ein gegebener Vektor ist;
LR Zerlegung Michael Sagraloff 15.06.2016 9
Warum will man eine LR-Zerlegung berechnen?
Berechnung der Determinante einer Matrix:
det(A) = det(P)−1 · det(L) · det(R) = (−1)v ·n∏
i=1
ri,i ,
wobei v die Anzahl der Zeilenwechsel bezeichnet, die beiAnwendung mit P herbeigeführt werden.
Lösen eines linearen Gleichungssystems: Bestimme x ∈ Kn,so dass
A · x = b ,
wobei b ∈ Kn ein gegebener Vektor ist;
LR Zerlegung Michael Sagraloff 15.06.2016 9
Warum will man eine LR-Zerlegung berechnen?
Berechnung der Determinante einer Matrix:
det(A) = det(P)−1 · det(L) · det(R) = (−1)v ·n∏
i=1
ri,i ,
wobei v die Anzahl der Zeilenwechsel bezeichnet, die beiAnwendung mit P herbeigeführt werden.
Lösen eines linearen Gleichungssystems: Bestimme x ∈ Kn,so dass
A · x = b ,
wobei b ∈ Kn ein gegebener Vektor ist;
LR Zerlegung Michael Sagraloff 15.06.2016 9
Warum will man eine LR-Zerlegung berechnen?
Berechnung der Determinante einer Matrix:
det(A) = det(P)−1 · det(L) · det(R) = (−1)v ·n∏
i=1
ri,i ,
wobei v die Anzahl der Zeilenwechsel bezeichnet, die beiAnwendung mit P herbeigeführt werden.
Lösen eines linearen Gleichungssystems: Bestimme x ∈ Kn,so dass
A · x = b ⇔ P · A · x = P · b⇔ L · (R · x)︸ ︷︷ ︸y
= P · b︸︷︷︸b′
,
wobei b ∈ Kn ein gegebener Vektor ist; Das Lösen diesesSystems reduziert sich dann auf das iterative Lösen derGleichungssyteme L · y = b′ und R · x = y in Dreiecksform.
LR Zerlegung Michael Sagraloff 15.06.2016 9
Beweis des Hauptsatzes, Teil A (Induktion über n)
Wir nehmen an, dass det A[i] 6= 0 für alle i = 1, . . . ,n.Insbesondere gilt also auch a1,1 6= 0.Für alle i = 2, . . . ,n, subtrahiere das li,1-fachen der erstenZeile von der i-ten Zeile, wobei li,1 :=
ai,1a1,1
:
Z1Z2...
Zn
a1,1 a1,2 · · · a1,na2,1 a2,2 · · · a2,n
.
.
....
.
.
....
an,1 an,2 · · · an,n
Z ′i =Zi−li,1·Z1−→
Z ′1Z ′2...
Z ′n
a1,1 a1,2 · · · a1,n0 a′2,2 · · · a′2,n...
.
.
....
.
.
.0 a′n,2 · · · a′n,n
,
Die Zeilenoperationen lassen sich als MatrixmultiplikationL1 · A darstellen, wobei
L1 :=
1 0 · · · 0−l2,1 1 · · · 0
.
.
....
.
.
....
−ln,1 0 · · · 1
LR Zerlegung Michael Sagraloff 15.06.2016 10
Beweis des Hauptsatzes, Teil A
Alle Hauptminoren der (n − 1)× (n − 1)-Matrix
A′ :=
(a′2,2 a′2,3 · · · a′2,n
.
.
....
.
.
....
a′n,2 a′n,3 · · · a′n,n
)sind von 0 verschieden da die Zeilenoperationen dieDeterminante einer beliebigen k × k -Untermatrix von A nichtverändert und daher det A[i+1] = a1,1 · det(A′)[i] gilt.Nach Induktionsannahme gibt es DreiecksmatrizenL′ = (l ′i,j)i,j=1,...,n−1 und R′ = (r ′i,j)i,j=1,...,n−1 mit A′ = L′ ·R′ undl ′i,i = 1 für alle i = 1, . . . ,n − 1, und somit gilt:
L1 · A :=
1 0 · · · 00 l′1,1 · · · 0
.
.
....
.
.
....
0 l′n−1,1 · · · l′n−1,n−1
︸ ︷︷ ︸
=:L̂
·
a1,1 a1,2 · · · an,n
0 r ′1,1 · · · r ′1,n−1...
.
.
....
.
.
.0 0 · · · r ′n−1,n−1
︸ ︷︷ ︸
=:R
LR Zerlegung Michael Sagraloff 15.06.2016 11
Beweis des Hauptsatzes, Teil A
Alle Hauptminoren der (n − 1)× (n − 1)-Matrix
A′ :=
(a′2,2 a′2,3 · · · a′2,n
.
.
....
.
.
....
a′n,2 a′n,3 · · · a′n,n
)sind von 0 verschieden da die Zeilenoperationen dieDeterminante einer beliebigen k × k -Untermatrix von A nichtverändert und daher det A[i+1] = a1,1 · det(A′)[i] gilt.Nach Induktionsannahme gibt es DreiecksmatrizenL′ = (l ′i,j)i,j=1,...,n−1 und R′ = (r ′i,j)i,j=1,...,n−1 mit A′ = L′ ·R′ undl ′i,i = 1 für alle i = 1, . . . ,n − 1, und somit gilt:
L1 · A :=
1 0 · · · 00 l′1,1 · · · 0
.
.
....
.
.
....
0 l′n−1,1 · · · l′n−1,n−1
︸ ︷︷ ︸
=:L̂
·
a1,1 a1,2 · · · an,n
0 r ′1,1 · · · r ′1,n−1...
.
.
....
.
.
.0 0 · · · r ′n−1,n−1
︸ ︷︷ ︸
=:R
LR Zerlegung Michael Sagraloff 15.06.2016 11
Beweis des Hauptsatzes, Teil A
Wir erhalten somit die Zerlegung A = L · R, wobei
L = L−11 · L̂ =
1 0 · · · 0
l2,1 1 · · · 0...
......
...ln,1 0 · · · 1
·
1 0 · · · 00 l ′1,1 · · · 0...
......
...0 l ′n−1,1 · · · l ′n−1,n−1
=
1 0 · · · 0
l2,1 l ′1,1 · · · 0...
......
...ln,1 l ′n−1,1 · · · l ′n−1,n−1
und l ′i,i = 1 für alle i nach Induktionsvoraussetzung.
Übungsaufgabe: Zeige die Eindeutigkeit der Zerlegung!Beachte, dass bei iterativer Anwendung des Verfahrens dieMatrizen L und R direkt bestimmt werden −→ Beispiel
LR Zerlegung Michael Sagraloff 15.06.2016 12
Beispiel einer LR Zerlegung
5 2 −1 02 5 0 1−1 0 5 −20 1 −2 5
︸ ︷︷ ︸
A0 :=A
L1·→
5 2 −1 00 21
525 1
0 25
245 −2
0 1 −2 5
︸ ︷︷ ︸
A1
L2·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 − 44
2110021
︸ ︷︷ ︸
A2
L3·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 0 96
25
︸ ︷︷ ︸
A3=:R
1 0 0 025 1 0 0− 1
5 0 1 00 0 0 1
︸ ︷︷ ︸
L̂1=L−11
→
1 0 0 025 1 0 0− 1
52
21 1 00 5
21 0 1
︸ ︷︷ ︸
L̂2=L̂1·L−12
→
1 0 0 025 1 0 0− 1
5221 1 0
0 521 − 44
100 1
︸ ︷︷ ︸
L̂3=:L̂2·L−13
Beachte, dass wir die Einträge der Matrix L (bzw. L̂i ) direkt indie freiwerdenden Stellen der Matrizen Ai speichern können.
LR Zerlegung Michael Sagraloff 15.06.2016 13
Beispiel einer LR Zerlegung
5 2 −1 02 5 0 1−1 0 5 −20 1 −2 5
︸ ︷︷ ︸
A0 :=A
L1·→
5 2 −1 00 21
525 1
0 25
245 −2
0 1 −2 5
︸ ︷︷ ︸
A1
L2·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 − 44
2110021
︸ ︷︷ ︸
A2
L3·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 0 96
25
︸ ︷︷ ︸
A3=:R
1 0 0 025 1 0 0− 1
5 0 1 00 0 0 1
︸ ︷︷ ︸
L̂1=L−11
→
1 0 0 025 1 0 0− 1
52
21 1 00 5
21 0 1
︸ ︷︷ ︸
L̂2=L̂1·L−12
→
1 0 0 025 1 0 0− 1
5221 1 0
0 521 − 44
100 1
︸ ︷︷ ︸
L̂3=:L̂2·L−13
Beachte, dass wir die Einträge der Matrix L (bzw. L̂i ) direkt indie freiwerdenden Stellen der Matrizen Ai speichern können.
LR Zerlegung Michael Sagraloff 15.06.2016 13
Beispiel einer LR Zerlegung
5 2 −1 02 5 0 1−1 0 5 −20 1 −2 5
︸ ︷︷ ︸
A0 :=A
L1·→
5 2 −1 00 21
525 1
0 25
245 −2
0 1 −2 5
︸ ︷︷ ︸
A1
L2·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 − 44
2110021
︸ ︷︷ ︸
A2
L3·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 0 96
25
︸ ︷︷ ︸
A3=:R
1 0 0 025 1 0 0− 1
5 0 1 00 0 0 1
︸ ︷︷ ︸
L̂1=L−11
→
1 0 0 025 1 0 0− 1
52
21 1 00 5
21 0 1
︸ ︷︷ ︸
L̂2=L̂1·L−12
→
1 0 0 025 1 0 0− 1
5221 1 0
0 521 − 44
100 1
︸ ︷︷ ︸
L̂3=:L̂2·L−13
Beachte, dass wir die Einträge der Matrix L (bzw. L̂i ) direkt indie freiwerdenden Stellen der Matrizen Ai speichern können.
LR Zerlegung Michael Sagraloff 15.06.2016 13
Beispiel einer LR Zerlegung
5 2 −1 02 5 0 1−1 0 5 −20 1 −2 5
︸ ︷︷ ︸
A0 :=A
L1·→
5 2 −1 00 21
525 1
0 25
245 −2
0 1 −2 5
︸ ︷︷ ︸
A1
L2·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 − 44
2110021
︸ ︷︷ ︸
A2
L3·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 0 96
25
︸ ︷︷ ︸
A3=:R
1 0 0 025 1 0 0− 1
5 0 1 00 0 0 1
︸ ︷︷ ︸
L̂1=L−11
→
1 0 0 025 1 0 0− 1
52
21 1 00 5
21 0 1
︸ ︷︷ ︸
L̂2=L̂1·L−12
→
1 0 0 025 1 0 0− 1
5221 1 0
0 521 − 44
100 1
︸ ︷︷ ︸
L̂3=:L̂2·L−13
Beachte, dass wir die Einträge der Matrix L (bzw. L̂i ) direkt indie freiwerdenden Stellen der Matrizen Ai speichern können.
LR Zerlegung Michael Sagraloff 15.06.2016 13
Beispiel einer LR Zerlegung
5 2 −1 02 5 0 1−1 0 5 −20 1 −2 5
︸ ︷︷ ︸
A0 :=A
L1·→
5 2 −1 00 21
525 1
0 25
245 −2
0 1 −2 5
︸ ︷︷ ︸
A1
L2·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 − 44
2110021
︸ ︷︷ ︸
A2
L3·→
5 2 −1 00 21
525 1
0 0 10021 − 44
210 0 0 96
25
︸ ︷︷ ︸
A3=:R
1 0 0 025 1 0 0− 1
5 0 1 00 0 0 1
︸ ︷︷ ︸
L̂1=L−11
→
1 0 0 025 1 0 0− 1
52
21 1 00 5
21 0 1
︸ ︷︷ ︸
L̂2=L̂1·L−12
→
1 0 0 025 1 0 0− 1
5221 1 0
0 521 − 44
100 1
︸ ︷︷ ︸
L̂3=:L̂2·L−13
Beachte, dass wir die Einträge der Matrix L (bzw. L̂i ) direkt indie freiwerdenden Stellen der Matrizen Ai speichern können.
LR Zerlegung Michael Sagraloff 15.06.2016 13
Algorithmus in Pseudocode
Algorithm 1: LREingabe :Matrix A = (ai,j)i,j=1,...,n ∈ Kn×n mit det A[i] 6= 0 für alle i = 1, . . . , n.Ausgabe :Untere Dreiecksmatrix L ∈ Kn×n mit Einseinträgen auf der
Diagonalen und obere Dreiecksmatrix R ∈ Kn×n, so dass A = L · R.L := (li,j)i,j=1,...,n := Idn
R := (ri,j)i,j=1,...,n := Afor j ← 1 to n − 1 do
for i ← j + 1 to n doli,j :=
ai,jaj,j
for k ← j to n dori,k := ri,k − li,j · ri,j
endend
endreturn (L,R) := ((li,j)i,j=1,...,n, (ri,j)i,j=1,...,n)
LR Zerlegung Michael Sagraloff 15.06.2016 14
Algorithmische Komplexität
LR-Zerlegung:In der j-ten Iteration benötigen wir n − j Divisionen zurBerechnung der Quotienten li,j für i = j + 1, . . . ,n, sowiejeweils (n − j)2 Additionen und Multiplikationen für dieZeilenoperationen.Summation über j = 1, . . . ,n liefert
2·n∑
j=1
(n− j)·(n− j+1) < 2·n∑
j=1
j2 = 2· n(n + 1)(2n + 1)6
≈ 2n3
3.
Lösen eine Gleichungssystems durch Vorwärts- undRückwärtseinsetzen unter Verwendung der LR-Zerlegung:
Eine Division und jeweils n − j Additionen und Multiplikationenim j-ten SchrittArithmetische Komplexität is also beschränkt durch
2 ·∑n
j=1(2(n − j) + 1) = 2n(n − 1) + 2n = 2n2
LR Zerlegung Michael Sagraloff 15.06.2016 15
Algorithmische Komplexität
LR-Zerlegung:In der j-ten Iteration benötigen wir n − j Divisionen zurBerechnung der Quotienten li,j für i = j + 1, . . . ,n, sowiejeweils (n − j)2 Additionen und Multiplikationen für dieZeilenoperationen.Summation über j = 1, . . . ,n liefert
2·n∑
j=1
(n− j)·(n− j+1) < 2·n∑
j=1
j2 = 2· n(n + 1)(2n + 1)6
≈ 2n3
3.
Lösen eine Gleichungssystems durch Vorwärts- undRückwärtseinsetzen unter Verwendung der LR-Zerlegung:
Eine Division und jeweils n − j Additionen und Multiplikationenim j-ten SchrittArithmetische Komplexität is also beschränkt durch
2 ·∑n
j=1(2(n − j) + 1) = 2n(n − 1) + 2n = 2n2
LR Zerlegung Michael Sagraloff 15.06.2016 15
Zusammenhang mit dem Gauß-Algorithmus
Beachte, dass das Gleichungssytem A · x = b auch direkt gelöstwerden kann; siehe dazu auch das Beispiel vom Anfang:
Alle Operationen, die auf A ausgeführt werden, werdenparallel auf dem Vektor b ausgeführtLösung des so erhaltenen Dreieckssystems bleibt invariantund kann durch Rückwärtseinsetzen gewonnen werden.Vorteil der direkten Berechnung einer LR Zerlegung:Bei gegebener LR Zerlegung kann man die Lösung für eingegebenes b in quadratischer Laufzeit bestimmen,wohingegen der Gaußalgorithmus kubische Laufzeit hat.
LR Zerlegung Michael Sagraloff 15.06.2016 16
Zur Implementierung
Exakte Berechnung über Q (eher theoretisch interessant):polynomielle Laufzeit des Algorithmus ist nicht trivial!Zwischenergebnisse lassen sich als Quotienten vonDeterminanten geeigneter Untermatrizen schreiben.
Bei der Implementierung des Algorithmus wird üblicherweiseFließkommaarithmetik verwendet (d.h. single oder doubleprecision in IEEE754 Format)⇒ Rechenfehler!?Wahl des Pivotelements ist entscheidend für die numerischeStabilität des VerfahrensAufgrund der kubischen Laufzeit wird der Algorithmus nur fürkleine und mittelgroße Matrizen (d.h. n ≤ 10000) verwendet.Speicherbedarf bei n = 104: 108 Einträgen a 64 bit ≈ 1 GBintegriert in den LA Bibliotheken NAG, IMSL, und LAPACK.
LR Zerlegung Michael Sagraloff 15.06.2016 17
Zur Implementierung
Exakte Berechnung über Q (eher theoretisch interessant):polynomielle Laufzeit des Algorithmus ist nicht trivial!Zwischenergebnisse lassen sich als Quotienten vonDeterminanten geeigneter Untermatrizen schreiben.
Bei der Implementierung des Algorithmus wird üblicherweiseFließkommaarithmetik verwendet (d.h. single oder doubleprecision in IEEE754 Format)⇒ Rechenfehler!?Wahl des Pivotelements ist entscheidend für die numerischeStabilität des VerfahrensAufgrund der kubischen Laufzeit wird der Algorithmus nur fürkleine und mittelgroße Matrizen (d.h. n ≤ 10000) verwendet.Speicherbedarf bei n = 104: 108 Einträgen a 64 bit ≈ 1 GBintegriert in den LA Bibliotheken NAG, IMSL, und LAPACK.
LR Zerlegung Michael Sagraloff 15.06.2016 17
Zur Implementierung
Exakte Berechnung über Q (eher theoretisch interessant):polynomielle Laufzeit des Algorithmus ist nicht trivial!Zwischenergebnisse lassen sich als Quotienten vonDeterminanten geeigneter Untermatrizen schreiben.
Bei der Implementierung des Algorithmus wird üblicherweiseFließkommaarithmetik verwendet (d.h. single oder doubleprecision in IEEE754 Format)⇒ Rechenfehler!?Wahl des Pivotelements ist entscheidend für die numerischeStabilität des VerfahrensAufgrund der kubischen Laufzeit wird der Algorithmus nur fürkleine und mittelgroße Matrizen (d.h. n ≤ 10000) verwendet.Speicherbedarf bei n = 104: 108 Einträgen a 64 bit ≈ 1 GBintegriert in den LA Bibliotheken NAG, IMSL, und LAPACK.
LR Zerlegung Michael Sagraloff 15.06.2016 17
Zur Implementierung
Exakte Berechnung über Q (eher theoretisch interessant):polynomielle Laufzeit des Algorithmus ist nicht trivial!Zwischenergebnisse lassen sich als Quotienten vonDeterminanten geeigneter Untermatrizen schreiben.
Bei der Implementierung des Algorithmus wird üblicherweiseFließkommaarithmetik verwendet (d.h. single oder doubleprecision in IEEE754 Format)⇒ Rechenfehler!?Wahl des Pivotelements ist entscheidend für die numerischeStabilität des VerfahrensAufgrund der kubischen Laufzeit wird der Algorithmus nur fürkleine und mittelgroße Matrizen (d.h. n ≤ 10000) verwendet.Speicherbedarf bei n = 104: 108 Einträgen a 64 bit ≈ 1 GBintegriert in den LA Bibliotheken NAG, IMSL, und LAPACK.
LR Zerlegung Michael Sagraloff 15.06.2016 17
Zur Implementierung
Exakte Berechnung über Q (eher theoretisch interessant):polynomielle Laufzeit des Algorithmus ist nicht trivial!Zwischenergebnisse lassen sich als Quotienten vonDeterminanten geeigneter Untermatrizen schreiben.
Bei der Implementierung des Algorithmus wird üblicherweiseFließkommaarithmetik verwendet (d.h. single oder doubleprecision in IEEE754 Format)⇒ Rechenfehler!?Wahl des Pivotelements ist entscheidend für die numerischeStabilität des VerfahrensAufgrund der kubischen Laufzeit wird der Algorithmus nur fürkleine und mittelgroße Matrizen (d.h. n ≤ 10000) verwendet.Speicherbedarf bei n = 104: 108 Einträgen a 64 bit ≈ 1 GBintegriert in den LA Bibliotheken NAG, IMSL, und LAPACK.
LR Zerlegung Michael Sagraloff 15.06.2016 17
Zur Implementierung
Exakte Berechnung über Q (eher theoretisch interessant):polynomielle Laufzeit des Algorithmus ist nicht trivial!Zwischenergebnisse lassen sich als Quotienten vonDeterminanten geeigneter Untermatrizen schreiben.
Bei der Implementierung des Algorithmus wird üblicherweiseFließkommaarithmetik verwendet (d.h. single oder doubleprecision in IEEE754 Format)⇒ Rechenfehler!?Wahl des Pivotelements ist entscheidend für die numerischeStabilität des VerfahrensAufgrund der kubischen Laufzeit wird der Algorithmus nur fürkleine und mittelgroße Matrizen (d.h. n ≤ 10000) verwendet.Speicherbedarf bei n = 104: 108 Einträgen a 64 bit ≈ 1 GBintegriert in den LA Bibliotheken NAG, IMSL, und LAPACK.
LR Zerlegung Michael Sagraloff 15.06.2016 17
Beispiel zu unterschiedlichen Pivotstrategien
Gegeben sei die Matrix:
A :=
(2.3 1.8 11.4 1.1 −0.70.8 4.3 2.1
)≈ 1 0 0
0.609 1 00.348 845 1
︸ ︷︷ ︸
L
· 2.3 1.8 1
0 0.04 1.3090 0 1108
︸ ︷︷ ︸
R
und das Gleichungssytem A · x = b := (1.2,−2.1,0.6)t mit der"exakten" Lösung x ≈ (0.3,−0.980,2.160).
Berechnung der LR Zerlegung mit dreistelliger Genauigkeit(ohne Pivotstrategie): 2.3 1.8 1.0
1.4 1.1 −0.70.8 4.3 2.1
→ 2.3 1.8 1.0
0.609 0.0038 −1.310.348 3.67 1.75
→
2.3 1.8 1.00.609 0.0038 −1.310.348 966 1270
Rückwärtseinsetzen liefert die Lösung x̃ = (2.37,−3.55,2.15).
LR Zerlegung Michael Sagraloff 15.06.2016 18
Beispiel zu unterschiedlichen Pivotstrategien
Wahl des Pivotelements 3.67 an Stelle von 0.0038 im zweitenIterationschritt liefert die LR Zerlegung: 1 0 0
0 0 10 1 0
︸ ︷︷ ︸
P
·A ≈ 1 0 0
0.348 1 00.609 0.001 1
︸ ︷︷ ︸
L̃
· 2.3 1.8 1
0 3.67 1.750 0 −1.31
︸ ︷︷ ︸
R̃
,
also eine gute Approximation der "exakten" LR Zerlegung:
P · A =
1 0 00.349 1 00.609 0.001 1
︸ ︷︷ ︸
≈L
· 2.3 1.8 1
0 3.674 1.7520 0 −1.311
︸ ︷︷ ︸
≈R
,
Rückwärtseinsetzen liefert eine gute Approximation derLösung:
x̃ = (0.35,−0.98,2.16).
LR Zerlegung Michael Sagraloff 15.06.2016 19
Anhang: Beweis des Hauptsatzes, Teil B
Vertausche die erste Zeile von A mit einer beliebigen Zeile i0für die ai0,1 6= 0; Den Eintrag ai0,1 nennt man Pivotelement.
Z1...
Zi0...
Zn
a1,1 a1,2 · · · a1,n
.
.
....
.
.
....
ai0,1· · · · · · ai0,n
.
.
....
.
.
....
an,1 an,2 · · · an,n
Z1↔Zi0−→
Zi0...
Z1...
Zn
ai0,1ai0,2
· · · ai0,n
.
.
....
.
.
....
a1,1 a1,2 · · · a1,n...
.
.
....
.
.
.0 an,2 · · · an,n
︸ ︷︷ ︸
A′=(a′i,j )i,j=1,...,n
,
Die Zeilenvertauschung der i-ten und j-ten Zeile ist als Matrix-multiplikation P[i , j] · A darstellbar (d.h. A′ = P[1, i0] · A), wobei
P[i , j] :=
1 · · · · · · · · · · · · · · · 0...
.
.
....
.
.
....
.
.
....
0 · · · 0 · · · 1 · · · 0...
.
.
....
.
.
....
.
.
....
0 · · · 1 · · · 0 · · · 0...
.
.
....
.
.
....
.
.
....
0 · · · · · · · · · · · · · · · 1
1...i...j...n
LR Zerlegung Michael Sagraloff 15.06.2016 20
Beweis des Hauptsatzes, Teil B
Für alle i 6= 1, subtrahiere das li,1-fache der ersten Zeile von A′
von der i-ten Zeile von A′, wobei li,1 :=a′1,ia′1,1
:
Z ′1Z ′2...
Z ′n
a′1,1 a′1,2 · · · a′1,na′2,1 a′2,2 · · · a′2,n
.
.
....
.
.
....
a′n,1 a′n,2 · · · a′n,n
Z ′′i =Z ′i −li,1·Z ′1−→
Z ′′1Z ′′2
.
.
.Z ′′n
ai0,1
ai0,2· · · ai0,n
0 a′′2,2 · · · a′′2,n...
.
.
....
.
.
.0 a′′n,2 · · · a′′n,n
︸ ︷︷ ︸
A′′=(a′′i,j )i,j=1,...,n
,
Wie im Beweis zu Teil A stellen wir die Zeilenoperationen alsMatrixmultiplikation dar, d.h.
1 0 · · · 0−l2,1 1 · · · 0
.
.
....
.
.
....
−ln,1 0 · · · 1
︸ ︷︷ ︸
L1
·A′ = A′′,
und somit L1 · P[1, i0] · A = A′′.
LR Zerlegung Michael Sagraloff 15.06.2016 21
Beweis des Hauptsatzes, Teil B
Nach Induktionsannahme gibt es eine (n − 1)× (n − 1)-Permutationsmatrix P ′′ = (p′′i,j)
′i,j=1,...,n−1 und Dreiecksmatrizen
L′′ = (l ′′i,j)′i,j=1,...,n−1 und R′′ = (r ′′i,j)
′i,j=1,...,n−1 mit
1 0 · · · 00 p′′1,1 · · · p′′1,n−1...
.
.
....
.
.
.0 p′′n−1,1 · · · p′′n−1,n−1
︸ ︷︷ ︸
P∗
·
ai0,1
ai0,2· · · ai0,n
0 a′′2,2 · · · a′′2,n...
.
.
....
.
.
.0 a′′n,2 · · · a′′n,n
︸ ︷︷ ︸
A′′
=
1 0 · · · 00 l′′1,1 · · · 0
.
.
....
.
.
....
0 l′′n−1,1 · · · l′′n−1,n−1
︸ ︷︷ ︸
L∗
·
ai0,1
ai0,2· · · ai0,n
0 r ′′1,1 · · · r ′′1,n−1...
.
.
....
.
.
.0 0 · · · r ′′n−1,n−1
︸ ︷︷ ︸
R∗
,
und somit P∗ · L1 · P[1, i0] · A = L∗ · R∗.
Leider kommutieren P∗ und L1 nicht, aber ...
LR Zerlegung Michael Sagraloff 15.06.2016 22
Beweis des Hauptsatzes, Teil B
Nach Induktionsannahme gibt es eine (n − 1)× (n − 1)-Permutationsmatrix P ′′ = (p′′i,j)
′i,j=1,...,n−1 und Dreiecksmatrizen
L′′ = (l ′′i,j)′i,j=1,...,n−1 und R′′ = (r ′′i,j)
′i,j=1,...,n−1 mit
1 0 · · · 00 p′′1,1 · · · p′′1,n−1...
.
.
....
.
.
.0 p′′n−1,1 · · · p′′n−1,n−1
︸ ︷︷ ︸
P∗
·
ai0,1
ai0,2· · · ai0,n
0 a′′2,2 · · · a′′2,n...
.
.
....
.
.
.0 a′′n,2 · · · a′′n,n
︸ ︷︷ ︸
A′′
=
1 0 · · · 00 l′′1,1 · · · 0
.
.
....
.
.
....
0 l′′n−1,1 · · · l′′n−1,n−1
︸ ︷︷ ︸
L∗
·
ai0,1
ai0,2· · · ai0,n
0 r ′′1,1 · · · r ′′1,n−1...
.
.
....
.
.
.0 0 · · · r ′′n−1,n−1
︸ ︷︷ ︸
R∗
,
und somit P∗ · L1 · P[1, i0] · A = L∗ · R∗.
Leider kommutieren P∗ und L1 nicht, aber ...
LR Zerlegung Michael Sagraloff 15.06.2016 22
Beweis des Hauptsatzes, Teil B
mit (1, l∗2 , . . . , l∗n)
t := P∗ · (1, l2,1, . . . , ln,1)t gilt:
P∗ · L1 =
1 0 · · · 00 p′′1,1 · · · p′′1,n−1...
......
...0 p′′n−1,1 · · · p′′n−1,n−1
·
1 0 · · · 0−l2,1 1 · · · 0
......
......
−ln,1 0 · · · 1
=
1 0 · · · 0−l∗2 p′′1,1 · · · p′′1,n−1
......
......
−l∗n p′′n−1,1 · · · p′′n−1,n−1
=
1 0 · · · 0−l∗2 1 · · · 0
......
......
−l∗n 0 · · · 1
︸ ︷︷ ︸
L∗1
·
1 0 · · · 00 p′′1,1 · · · p′′1,n−1...
......
...0 p′′n−1,1 · · · p′′n−1,n−1
,
und somit L∗1 · P∗ · P[1, i0]︸ ︷︷ ︸
=:P
·A = L∗ ·R∗ ⇔ P · A = (L∗1)−1 · L∗︸ ︷︷ ︸=:L
·R.
LR Zerlegung Michael Sagraloff 15.06.2016 23
