Es gibt viele Artikel über die LU-Zerlegung, aber da es keinen Artikel gab, in dem das Merkmal der LU-Zerlegung erwähnt wurde, dass Berechnungen mit Speichereinsparung durchgeführt werden können, habe ich versucht, die LU-Zerlegung etwas genauer zu erklären und erneut auf die Speichereinsparung zu achten. Überlegen.
Nach dem Studium von React habe ich eine Site erstellt, die leicht eine LU-Zerlegung durchführen kann. Zur Vereinfachung des Renderns ist die nachstehend beschriebene speichersparende Implementierung jedoch nicht implementiert. LU-Zerlegungserfahrung
Bevor wir uns mit der LU-Zerlegung befassen, wollen wir uns mit dreieckigen Matrizen befassen. Das Ziel der LU-Zerlegung ist **, die Koeffizienten simultaner Gleichungen nur in einer Dreiecksmatrix ** darzustellen. Ich werde den Grund erklären. Betrachten Sie zunächst eine einfache simultane Gleichung.
\left\{
\begin{array}{rrrrrrr}
2x_1 & & & & &=& 2 \\
x_1 &+& 4x_2 & & &=& 9 \\
4x_1 &-& 3x_2 &+& 3x_3 &=& -5
\end{array}
\right.
Die Lösung dieser simultanen Gleichung kann leicht erhalten werden, indem sie nacheinander gelöst wird.
\begin{align}
x_1 &= \frac{2}{2} = 1 \\
x_2 &= \frac{9 - x_1}{4} = \frac{9-1}{4} = 2 \\
x_3 &= \frac{-5 - 4x_1 + 3x_2}{3} = \frac{-5 - 4 + 6}{3} = -1
\end{align}
Wenn diese simultane Gleichung durch eine Matrix ausgedrückt wird,
\left(
\begin{matrix}
2 & 0 & 0 \\
1 & 4 & 0 \\
4 & -3 & 3
\end{matrix}
\right)
\left(
\begin{matrix}
x_1 \\
x_2 \\
x_3
\end{matrix}
\right)
=
\left(
\begin{matrix}
2 \\
9 \\
-5
\end{matrix}
\right) \tag{1}
Es sieht aus wie. Sie können sehen, dass die Koeffizienten der simultanen Gleichungen durch eine ** untere Dreiecksmatrix ** dargestellt werden, wobei die Elemente über der diagonalen Komponente $ 0 $ sind. Wenn Sie in diesem Fall $ x_i (i \ in 2, 3) $ finden, kennen Sie bereits $ x_ {i-1} (i \ in 1,2) $, sodass Sie sie der Reihe nach lösen können. Diese Methode wird als ** Vorwärtssubstitution ** bezeichnet. Betrachten Sie als nächstes die folgenden simultanen Gleichungen.
\left\{
\begin{array}{rrrrrrr}
3x_1 &-& 3x_2 &+& 4x_3 &=& -5 \\
& & 4x_2 &+& x_3 &=& 9 \\
& & & & 2x_3 &=& 2
\end{array}
\right. \\
\left(
\begin{matrix}
3 & -3 & 4\\
0 & 4 & 1 \\
0 & 0 & 2
\end{matrix}
\right)
\left(
\begin{matrix}
x_1 \\
x_2 \\
x_3
\end{matrix}
\right)
=
\left(
\begin{matrix}
-5 \\
9 \\
2
\end{matrix}
\right)
In diesem Fall werden die Koeffizienten der simultanen Gleichungen durch eine ** obere Dreiecksmatrix ** dargestellt, wobei die Komponenten unter der diagonalen Komponente $ 0 $ sind. In diesem Fall ist beim Finden von $ x_i (i \ in 2, 1) $ $ x_ {i + 1} (i \ in 3, 2) $ bekannt, also die umgekehrte Reihenfolge im Fall der unteren Dreiecksmatrix Kann gelöst werden. Diese Methode wird als ** Rückwärtssubstitution ** bezeichnet. Auf diese Weise kann, wenn die Koeffizienten simultaner Gleichungen durch eine Dreiecksmatrix dargestellt werden, dies mit einem einfachen Algorithmus gelöst werden. Hier werde ich die Wörter definieren. Der Ausdruck $ (1) $ wird wie folgt ausgedrückt.
A \boldsymbol{x} = \boldsymbol{b} \\
A =
\left(
\begin{matrix}
3 & -3 & 4\\
0 & 4 & 1 \\
0 & 0 & 2
\end{matrix}
\right), \
\boldsymbol{x} =
\left(
\begin{matrix}
x_1 \\
x_2 \\
x_3
\end{matrix}
\right), \
\boldsymbol{b} =
\left(
\begin{matrix}
2 \\
9 \\
-5
\end{matrix}
\right)
Hier wird $ A $ als Systemmatrix bezeichnet, $ \ boldsymbol {x} $ als Lösungsvektor und $ \ boldsymbol {b} $ als Vektor auf der rechten Seite. Im Folgenden wird diese Notation verallgemeinert, wobei $ A $ die Matrix $ n \ times n $ und $ \ boldsymbol {x}, \ \ boldsymbol {b} $ der Dimensionsvektor $ n $ ist.
Wie oben erwähnt, ist es fast so, als würde man die Gleichungen lösen, wenn die Koeffizienten simultaner Gleichungen durch eine dreieckige Matrix dargestellt werden können. Betrachten Sie nun die Koeffizientenmatrix $ A $ als Produkt der unteren Dreiecksmatrix $ L $ und der oberen Dreiecksmatrix $ U $.
LU \boldsymbol{x} = \boldsymbol{b} \\
Wenn $ U \ boldsymbol {x} = \ boldsymbol {y} $, kann $ \ boldsymbol {y} $ durch Vorwärtssubstitution erhalten werden, und unter Verwendung des Ergebnisses kann $ \ boldsymbol {x} $ durch Rückwärtssubstitution erhalten werden. Sie können sehen, dass es erforderlich ist.
\begin{align}
L \boldsymbol{y} &= \boldsymbol{b} \\
U \boldsymbol{x} &= \boldsymbol{y}
\end{align}
Wie finden Sie $ L $ und $ U $? Bevor wir uns dem Hauptthema zuwenden, konzentrieren wir uns auf den Freiheitsgrad. Es gibt $ n ^ 2 $ Elemente für $ A $ und $ n (n + 1) / 2 $ Elemente für $ L bzw. \ U $ für insgesamt $ n ^ 2 + n $ Elemente. .. Daher ist nach der Zerlegung der Freiheitsgrad $ n $ größer, aber dies bedeutet, dass die diagonalen Komponenten von $ L $ alle $ 1 $ sind, so dass der Freiheitsgrad $ n ^ 2 $ wird, was zu $ A $ wird. Andererseits können $ L und \ U $ eindeutig bestimmt werden. (Übrigens spielt es keine Rolle, ob die diagonale Komponente von $ L $ 1 $ oder die diagonale Komponente von $ U $ 1 $ beträgt. Dies hängt von der Person ab.)
Jetzt werde ich erklären, wie man es konkret findet, aber vorher werde ich die kleine Matrix erklären.
Betrachten Sie die folgende Matrix $ M $.
M = \left(
\begin{matrix}
m_{11} & m_{12} &|& m_{13} & m_{14} & m_{15} & m_{16} \\
m_{21} & m_{22} &|& m_{23} & m_{24} & m_{25} & m_{26} \\
- & - & - & - & - & - & - \\
m_{31} & m_{32} &|& m_{33} & m_{34} & m_{35} & m_{36} \\
m_{41} & m_{42} &|& m_{43} & m_{44} & m_{45} & m_{46} \\
m_{51} & m_{52} &|& m_{53} & m_{54} & m_{55} & m_{56} \\
m_{61} & m_{62} &|& m_{63} & m_{64} & m_{65} & m_{66}
\end{matrix}
\right)
Trennen Sie diese Matrix durch eine gepunktete Linie und geben Sie jedem Block (** kleine Matrix **) einen Namen.
M_{11} =
\left(
\begin{matrix}
m_{11} & m_{12} \\
m_{21} & m_{22}
\end{matrix}
\right)
, \
M_{12} =
\left(
\begin{matrix}
m_{13} & m_{14} & m_{15} & m_{16} \\
m_{23} & m_{24} & m_{25} & m_{26}
\end{matrix}
\right) \\
M_{21} =
\left(
\begin{matrix}
m_{31} & m_{32} \\
m_{41} & m_{42} \\
m_{51} & m_{52} \\
m_{61} & m_{62}
\end{matrix}
\right)
, \
M_{22} =
\left(
\begin{matrix}
m_{33} & m_{34} & m_{35} & m_{36} \\
m_{43} & m_{44} & m_{45} & m_{46} \\
m_{53} & m_{54} & m_{55} & m_{56} \\
m_{63} & m_{64} & m_{65} & m_{66}
\end{matrix}
\right)
Dann wird die Matrix $ M $ ausgedrückt als:
M =
\left(
\begin{matrix}
M_{11} & M_{12} \\
M_{21} & M_{22}
\end{matrix}
\right)
Betrachten Sie nun die $ 6 \ times 6 $ -Matrix $ N $ und versuchen Sie, sie auf die gleiche Weise wie $ M $ zu trennen.
N =
\left(
\begin{matrix}
N_{11} & N_{12} \\
N_{21} & N_{22}
\end{matrix}
\right)
Und das Produkt von $ M $ und $ N $ kann tatsächlich ausgedrückt werden als:
\begin{align}
MN &=
\left(
\begin{matrix}
M_{11} & M_{12} \\
M_{21} & M_{22}
\end{matrix}
\right)
\left(
\begin{matrix}
N_{11} & N_{12} \\
N_{21} & N_{22}
\end{matrix}
\right) \\
&=
\left(
\begin{matrix}
N_{11}N_{11} + N_{12}N_{21} & N_{11}N_{12} + N_{12}N_{22} \\
N_{21}N_{11} + N_{22}N_{21} & N_{21}N_{12} + N_{22}N_{22}
\end{matrix}
\right) \\
\end{align}
Mit anderen Worten, Sie können eine kleine Matrix wie ein Element einer Matrix behandeln. (Sie können es sehen, indem Sie es tatsächlich berechnen.) Jede Submatrix muss jedoch in $ M $ und $ N $ unterteilt werden, damit sie konsistent berechnet werden kann.
Die Einführung wurde verlängert, aber hier sind die Schritte, um endlich die Matrizen $ L $ und $ U $ zu finden. Teilen Sie zuerst $ L $ und $ U $ in kleinere Matrizen, wie unten gezeigt. Wobei $ O $ eine Matrix mit allen Elementen $ 0 $ ist.
L =
\left(
\begin{matrix}
1 & O \\
L_{21} & L_{22}
\end{matrix}
\right), \
U =
\left(
\begin{matrix}
u_{11} & U_{12} \\
O & U_{22}
\end{matrix}
\right)
$ L_ {22} und U_ {22} $ sind das untere und obere Dreieck von $ (n-1) \ times (n-1) $. Teilen Sie auch die Koeffizientenmatrix $ A $ auf die gleiche Weise.
A =
\left(
\begin{matrix}
a_{11} & A_{12} \\
A_{21} & A_{22}
\end{matrix}
\right)
Da $ A = LU $,
\begin{align}
\left(
\begin{matrix}
a_{11} & A_{12} \\
A_{21} & A_{22}
\end{matrix}
\right)
&=
\left(
\begin{matrix}
1 & O \\
L_{21} & L_{22}
\end{matrix}
\right)
\left(
\begin{matrix}
u_{11} & U_{12} \\
O & U_{22}
\end{matrix}
\right) \\
&=
\left(
\begin{matrix}
u_{11} & U_{12} \\
L_{21}u_{11} & L_{21}U_{12} + L_{22}U_{22}
\end{matrix}
\right)
\end{align}
Wird ausgedrückt als. Aus dieser Beziehung
u_{11} = a_{11}, \ U_{12} = A_{12}
Bekommen Unter der Annahme von $ a_ {11} \ neq 0 $,
L_{21} = \frac{A_{21}}{a_{11}}
Und $ L_ {21} $ werden ebenfalls bestimmt. Dies gibt Ihnen die erste Spalte von $ L $ und die erste Zeile von $ U $. Nun, ich denke, es ist bisher wahrscheinlich leicht zu verstehen. Die Frage ist, wie man mit $ A_ {22} = L_ {21} U_ {12} + L_ {22} U_ {22} $ umgeht. Wie bereits erwähnt, wurden hier bereits $ L_ {21} $ und $ U_ {12} $ erhalten, sodass der erste Term auf der rechten Seite bereits erhalten wurde. In Bezug auf den zweiten Term auf der rechten Seite, wie man ihn findet Wie am Anfang des Hauptthemas erwähnt, ist $ L_ {22} $ die untere Dreiecksmatrix und $ U_ {22} $ die obere Dreiecksmatrix. Durch Transponieren des ersten Terms
\hat{A}_{22} = A_{22} - L_{21}U_{12}
vorausgesetzt,
A^\prime_{22} = L_{22}U_{22}
Kann ausgedrückt werden als. Dies bedeutet, dass die $ (n-1) \ times (n-1) $ -Matrix $ A ^ \ prime_ {22} $ in eine untere Dreiecksmatrix und eine obere Dreiecksmatrix zerlegt wird. Daher wird durch Anwenden des gleichen Verfahrens wie oben für dieses $ A ^ \ prime_ {22} $ das Problem der LU-Zerlegung der Matrix von $ (n-2) \ mal (n-2) $ verringert. Darüber hinaus kann es sequentiell gelöst werden, was zu dem Problem der LU-Zerlegung der Matrix der $ (n-3) \ times (n-3) $ -Matrix führt ... Sobald Sie dies verstanden haben, müssen Sie es nur noch in Ihren Code einfügen.
Jetzt werde ich den Code schreiben, um die Matrizen $ L $ und $ U $ tatsächlich zu finden, aber wenn ich an nichts denke, das zweidimensionale Array für die Koeffizientenmatrix $ A $ und das zweidimensionale Array für die untere Dreiecksmatrix $ L $ , Möglicherweise möchten Sie drei Variablen im zweidimensionalen Array für die obere Dreiecksmatrix $ U $ haben. Tatsächlich werden in dem Artikel, in dem die LU-Zerlegung eingeführt wird (z. B. Ich habe die LU-Zerlegungsmethode programmgesteuert geschrieben), die Variablen für jede Matrix definiert. tun. Obwohl es sich um eine destruktive Methode handelt, ist nur eine Variable erforderlich, wenn Sie die Methode zum Überschreiben von $ A $ verwenden. Betrachten Sie insbesondere die Tatsache, dass die diagonale Komponente von $ L $ $ 1 $ ist, und verwalten Sie $ L $ und $ U $ in einer Matrix $ V $ wie folgt.
L =
\left(
\begin{matrix}
1 & 0 & 0 \\
l_{11} & 1 & 0 \\
l_{21} & l_{22} & 1 \\
\end{matrix}
\right), \
U =
\left(
\begin{matrix}
u_{11} & u_{12} & u_{13} \\
0 & u_{22} & u_{22} \\
0 & 0 & u_{33} \\
\end{matrix}
\right) \\
V =
\left(
\begin{matrix}
u_{11} & u_{12} & u_{13} \\
l_{11} & u_{22} & u_{22} \\
l_{21} & l_{22} & u_{33} \\
\end{matrix}
\right)
Wenn Sie $ V $ wie folgt betrachten und $ V $ über $ A $ überschreiben, beträgt der zum Speichern des Arrays erforderliche Speicher nur $ A $, wodurch Speicherplatz gespart werden kann.
Die Einführung ist ziemlich umfangreich geworden, aber ich werde ein Implementierungsbeispiel in Python vorstellen. Wie ich ausführlich erklärt habe, ist die Implementierung selbst überraschend einfach.
import numpy as np
def main():
n = int(input())
A = np.array([np.array(list(map(float, input().split()))) for i in range(n)])
for i in range(n):
for j in range(i + 1, n):
A[j][i] /= A[i][i]
for k in range(i + 1, n):
A[j][k] -= A[j][i] * A[i][k]
return A
if __name__ == "__main__":
print(main())
Kommentar wird hinzugefügt. Ich werde nur den Code setzen.
import numpy as np
def main():
n = int(input())
A = np.array([np.array(list(map(float, input().split()))) for i in range(n)])
P = np.identity(n)
for i in range(n):
maxidx = np.argmax(abs(A[i:n, i])) + i
A[i, :], A[maxidx, :] = A[maxidx, :], A[i, :].copy()
P[i, :], P[maxidx, :] = P[maxidx, :], P[i, :].copy()
for j in range(i + 1, n):
A[j][i] /= A[i][i]
for k in range(i + 1, n):
A[j][k] -= A[j][i] * A[i][k]
return A, P
if __name__ == "__main__":
A, P = main()
Tatsächlich ist die Berechnung möglicherweise nicht so stabil wie sie ist. Das liegt daran, dass wir $ a_ {ii} \ neq 0 $ annehmen. Um dieses Problem zu vermeiden, müssen Sie eine ** teilweise Pivot-Auswahl ** durchführen, die die Zeilen vertauscht. Ich werde diese Erklärung hinzufügen, wenn ich Zeit habe.
Numerische Berechnung für das Engineering
Wenn Sie Fragen haben, fragen Sie bitte! Ich bin mit Mathematik nicht sehr vertraut, daher kann ich möglicherweise nicht damit umgehen.
Recommended Posts