Lineare Gleichungssysteme haben jeweils die Form $A \cdot x = b$ wobei $A$ und $b$ gegeben und $x$ gesucht ist:
$$A = \left(
\begin{matrix}
a_{11} & a_{12} & \cdots & a_{1n} \\
a_{21} & a_{22} & \cdots & a_{2n} \\
\vdots & \vdots && \vdots \\
a_{n1} & a_{n2} & \cdots & a_{nn}
\end{matrix}
\right),
x = \left(
\begin{matrix}
x_1 \\
x_2 \\
\vdots \\
x_n
\end{matrix}
\right),
b = \left(
\begin{matrix}
b_1 \\
b_2 \\
\vdots \\
b_n
\end{matrix}
\right)$$
</div>
### Eigenschaften
- Gleich viele gesuchte Variablen $x_n$ wie Gleichungen $n$. Folglich:
- Die Matrix $A$ ist eine quadratische Matrix mit Dimensionen $n \times n$
- $A$ ist invertierbar
- $A$ hat eine Determinante $\det(A)$
### Dreiecks-Matrizen
<divclass="letters">
***$L$: Untere Dreiecksmatrix***
Eine Matrix, die in der oberen rechten Ecke nur den Wert $0$ und auf der Diagonale nur den Wert $1$ hat. Eine Untere Dreiecksmatrix hat also folgende Form:
$$L = \left(
\begin{matrix}
1 & 0 & 0 & \cdots & 0 \\
l_{21} & 1 & 0 & \cdots & 0 \\
l_{31} & l_{32} & 1 & \ddots & 0 \\
\vdots & \vdots & \ddots & \ddots & 0 \\
l_{n1} & l_{n2} & \cdots & l_{nn - 1} & 1
\end{matrix}
\right)$$
***$R$: Obere Dreiecksmatrix***
Eine Matrix, die unten links von der Diagonale nur den Wert $0$ beinhaltet. Eine Obere Dreiecksmatrix hat dementsprechend folgende Form:
$$R = \left(
\begin{matrix}
r_{11} & r_{12} & r_{13} & \cdots & r_{1n} \\
0 & r_{22} & r_{23} & \cdots & r_{2n} \\
0 & 0 & r_{33} & \cdots & r_{3n} \\
\vdots & \vdots & \ddots & \ddots & \vdots \\
0 & 0 & \cdots & 0 & r_{nn}
\end{matrix}
\right)$$
</div>
***Code-Beispiele:***
_Umwandlung in $R$-Matrix:_
```py
for i in range(n):
if A[i, i] == 0:
index = -1
for j in range(i + 1, n):
if A[j, i] > 0:
index = j
if index == -1:
raise Exception("Invalid Matrix")
else:
# Swap lines
A[[i, index]] = A[[index, i]]
for j in range(i + 1, n):
factor = A[j, i] / A[i, i]
A[j] = A[j] - (factor * A[i])
```
### Der Gauss-Algorithmus
Der Gauss-Algorithmus basiert darauf, dass ein lineares Gleichungssystem leicht lösbar ist, falls $A$ eine obere Dreiecksmatrix ist. $A$ muss also hierfür die Form einer oberen Dreiecksmatrix $R$ haben.
<divclass="formula">
***Gauss-Algorithmus:***
$$x_i = \frac{b_i - \sum_{j = i + 1}^n{a_{ij} \cdot x_j}}{a_{ii}}, i = n, n - 1, \dots, 1$$
</div>
Um den Gauss-Algorithmus anzuwenden, muss die Matrix $A$ erst in eine $R$-Matrix umgewandelt werden. Dies funktioniert wie folgt:
1. Mit $i$ von $1$ bis $n$
2. Falls $a_{ii}$ den Wert $0$ hat:
1. Mit $j$ von $i + 1$ bis $n$
2. Prüfen, ob $a_{ji}$ einen höheren Wert als $0$ hat
- Falls Zeile gefunden wurde:
- $a_{i}$ mit $a_{j}$ tauschen
- $b_{i}$ mit $b_{j}$ tauschen
- Sonst beenden: ungültige Matrix
3. Mit $j$ von $i + 1$ bis $n$
1. $a_k = a_k - \frac{a_{ki}}{a_{ii}} \cdot a_i$
2. $b_k = b_k - \frac{a_{ki}}{a_{ii}} \cdot b_i$
***Code-Beispiel:***
```py
from numpy import array, zeros
def gaussMethod(A, b):
A = array(A)
n = A.shape[0]
A = A.reshape((n, n))
b = array(b).reshape((n))
result = zeros(n)
# Convert to R-Matrix
for i in range(n):
maxIndex = i
for j in range(i + 1, n):
if A[j, i] > A[maxIndex, i]:
maxIndex = j
# Swap lines
A[[i, maxIndex]] = A[[maxIndex, i]]
b[[i, maxIndex]] = b[[maxIndex, i]]
for j in range(i + 1, n):
factor = A[j, i] / A[i, i]
A[j] = A[j] - (factor * A[i])
b[j] = b[j] - (factor * b[i])
# Calculate result
for index in range(n, 0, -1):
i = index - 1
value = b[i]
for j in range(i, n):
value = value - A[i, j] * result[j]
result[i] = value / A[i, i]
return result.reshape((n, 1))
```
### Fehlerfortpflanzung und Pivotisierung
- Da beim Umwandeln einer Matrix $A$ in die $R$-Form Zeilen in jedem Schritt mit dem Faktor $\lambda = \frac{a_{ji}}{a_{ii}}$ multipliziert werden, vergrössert sich der Schritt immer um $|\lambda|$
- $\lambda$ kann klein gehalten werden, indem Zeilen der Grösse nach sortiert werden
- In den Code-Beispielen ist dies bereits berücksichtigt
### Determinanten-Bestimmung
Die Determinante einer Matrix $A$ lässt sich einfach berechnen, sobald sie in die $R$-Form gebracht wurde mit folgender Formel:
<divclass="formula">
Determinanten-Bestimmung mit Matrix $\tilde{A}$ (die Matrix $A$ in der $R$-Form):
$$\det(A) =
(-1)^l \cdot \det(\tilde{A}) =
(-1)^l \cdot \prod_{i = 1}^n{\tilde{a_{ii}}}$$
</div>
***Code-Beispiel:***
```py
from numpy import array
def det(A):
l = 0
n = A.shape[0]
A = A.reshape((n, n))
# Convert to R-Matrix
for i in range(n):
maxIndex = i
for j in range(i + 1, n):
if A[j, i] > A[maxIndex, i]:
maxIndex = j
# Swap lines
A[[i, maxIndex]] = A[[maxIndex, i]]
l = l + 1
for j in range(i + 1, n):
factor = A[j, i] / A[i, i]
A[j] = A[j] - (factor * A[i])
result = 1
for i in range(n):
result = result * A[i, i]
return (-1 ** l) * result
```
### Die $LR$-Zerlegung
In der $LR$-Zerlegung wird die Matrix $A$ in die Matrizen $L$ und $R$ aufgeteilt, sodass $A = L \cdot R$ gilt.
Alternative Namen dieses Vorgangs sind ***$LR$-Faktorisierung*** und $LU$-decomposition.
<divclass="formula">
Für in $L$ und $R$ zerlegte Matrizen gilt:
$$A \cdot x = b$$
und
$$A \cdot x = L \cdot R \cdot x = L \cdot y = b$$
Aufwand: Berechnung der $LR$-Zerlegung mit Gauss-Algorithmus benötigt ca. $\frac{2}{3}n^3$ Punktoperationen.
</div>
Falls Zeilenvertauschungen stattfinden, entsteht bei der $LR$-Zerlegung eine zusätzliche Permutations-Matrix $P$.
<divclass="formula">
Für $L$ und $R$ zerlegte Matrizen mit Permutation $P$ gilt:
$$P \cdot A = L \cdot R$$
$$L \cdot y = P \cdot b$$
$$R \cdot x = y$$
</div>
Das Verfahren für die $LR$-Zerlegung ist identisch zu den Schritten bei der Umwandlung in eine $R$-Matrix. Jedoch wird jeweils der Wert $l_{ji}$ in der (zu Beginn) leeren Matrix $L$ mit dem im aktuellen Eliminationsschritt gesetzt. Zudem muss bei Vertauschungen die Permutations-Matrix $P$ entsprechend angepasst werden:
***Code-Beispiel:***
```py
from numpy import array, identity, zeros
def decomposite(A):
l = 0
n = A.shape[0]
R = A.reshape((n, n))
L = zeros((n, n))
P = identity((n, n))
# Convert to LR-Matrix
for i in range(n):
maxIndex = i
for j in range(i + 1, n):
if A[j, i] > A[maxIndex, i]:
maxIndex = j
# Swap lines
Pn = identity((n, n))
A[[i, maxIndex]] = A[[maxIndex, i]]
Pn[[i, maxIndex]] = Pn[[maxIndex, i]]
P = P * Pn
for j in range(i + 1, n):
factor = R[j, i] / R[i, i]
L[j, i] = factor
R[j] = R[j] - (factor * R[i])
result = 1
for i in range(n):
result = result * R[i, i]
return [L, R, P]
```
Wenn die $LR$-Zerlegung, wie in diesem Code, Zeilenaustausch und das Berechnen von $P$ involviert, spricht man von einer $LR$-Zerlegung mit **Spaltenmaximum-Strategie**.
***Vorgang:***
1. Gemäss vorhergehender Beschreibung und Code-Beispiel die Matrizen $L$ und $R$ berechnen
2. Mit Hilfe des Gauss-Algorithmus $L \cdot y = P \cdot b$ nach $y$ auflösen
3. Mit Hilfe des Gauss-Algorithmus $R \cdot x = y$ nach $x$ auflösen
