1. Macierze Symetryczne
Problem pozostaje taki sam — mamy do policzenia wynik równania A x = b ,
Ax = b,
A x = b , A T = A .
A^T = A.
A T = A .
1.1. Zmodyfikowany rozkład L U LU L U
Stosujemy rozkład L U LU L U :
A = L × U
A = L \times U
A = L × U
przy czym teraz zapisujemy to samo nieco inaczej:
L U = L D U ‾
L U = L D \overline{U}
L U = L D U
gdzie
L = [ 1 l 21 1 l 31 l 32 1 ⋮ ⋱ l n 1 l n 2 ⋯ l n , n − 1 1 ]
L =
\begin{bmatrix}
1\\
l_{21} & 1\\
l_{31} & l_{32} & 1\\
\vdots &&&& \ddots\\
l_{n1} & l_{n2} & \dotsb && l_{n,n-1} & 1
\end{bmatrix}
L = ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ 1 l 2 1 l 3 1 ⋮ l n 1 1 l 3 2 l n 2 1 ⋯ ⋱ l n , n − 1 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤
U ‾ = [ 1 u 12 … … u 1 n 1 u 23 … u 2 n ⋱ ⋮ 1 u n − 1 , n 1 ]
\overline{U} =
\begin{bmatrix}
1 & u_{12} & \dots & \dots & u_{1n}\\
& 1 & u_{23} & \dots & u_{2n}\\
&& \ddots && \vdots\\
&&& 1 & u_{n-1,n}\\
&&&& 1
\end{bmatrix}
U = ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ 1 u 1 2 1 … u 2 3 ⋱ … … 1 u 1 n u 2 n ⋮ u n − 1 , n 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤
D = diag ( d i )
D = \operatorname{diag}(d_i)
D = d i a g ( d i )
Warto zauważyć, że w takiej konfiguracji mamy: A T = U ‾ T D L T = L D U ‾ = A
A^T = \overline{U}^T D L^T = L D \overline{U} = A
A T = U T D L T = L D U = A A A A
1.2. Sposób liczenia wektora x x x
Jesteśmy uzbrojeni w rozkład macierzy: A = L D L T
A = L D L^T
A = L D L T
taki, że
$$ LD =
\[10pt]
L^T =
\[10pt]
A =
$$
Czyli liczymy po kolei:
mamy d 1 = a 11 d_{1} = a_{11} d 1 = a 1 1 d 1 d_{1} d 1
mamy:
d 1 l 21 = a 21 d_1 l_{21} = a_{21} d 1 l 2 1 = a 2 1 l 21 l_{21} l 2 1 d 1 l 21 ⋅ l 21 + d 2 = a 22 d_1 l_{21} \cdot l_{21} + d_2 = a_{22} d 1 l 2 1 ⋅ l 2 1 + d 2 = a 2 2 d 2 d_2 d 2 mamy:
d 1 l 31 = a 31 d_1 l_{31} = a_{31} d 1 l 3 1 = a 3 1 l 31 l_{31} l 3 1 d 1 l 31 ⋅ l 21 + d 2 l 32 = a 32 d_1 l_{31} \cdot l_{21} + d_2 l_{32} = a_{32} d 1 l 3 1 ⋅ l 2 1 + d 2 l 3 2 = a 3 2 l 32 l_{32} l 3 2 d 1 l 31 ⋅ l 31 + d 2 l 32 ⋅ l 32 + d 3 = a 33 d_1 l_{31} \cdot l_{31} + d_2 l_{32} \cdot l_{32} + d_3 = a_{33} d 1 l 3 1 ⋅ l 3 1 + d 2 l 3 2 ⋅ l 3 2 + d 3 = a 3 3 d 3 d_3 d 3 itd.
aż w końcu:
d 1 l n 1 = a n 1 d_1 l_{n1} = a_{n1} d 1 l n 1 = a n 1 d 1 l n 1 ⋅ l 21 + d 2 l n 2 = a n 2 d_1 l_{n1} \cdot l_{21} + d_2 l_{n2} = a_{n2} d 1 l n 1 ⋅ l 2 1 + d 2 l n 2 = a n 2 …
d 1 l n 1 ⋅ l n 1 + d 2 l n 2 ⋅ l n 2 + ⋯ + d n − 1 l n − 1 , n ⋅ l n − 1 , n + d n = a n n d_1 l_{n1} \cdot l_{n1} + d_2 l_{n2} \cdot l_{n2} + \dots + d_{n-1} l_{n-1, n} \cdot l_{n-1, n} + d_n = a_{nn} d 1 l n 1 ⋅ l n 1 + d 2 l n 2 ⋅ l n 2 + ⋯ + d n − 1 l n − 1 , n ⋅ l n − 1 , n + d n = a n n
Możemy dodatkowo określić ciąg r j = d i l n , j r_j = d_i l_{n,j} r j = d i l n , j
Powyższy algorytm w formie programu:
for k : = 1 k := 1 k : = 1 to n n n
for j : = 1 j := 1 j : = 1 to k − 1 k-1 k − 1
r j : = a k , j − ∑ q = 1 j − 1 l j , q ⋅ r q r_{j} := a_{k,j} - \sum_{q=1}^{j-1} l_{j,q} \cdot r_{q} r j : = a k , j − ∑ q = 1 j − 1 l j , q ⋅ r q l k , j : = r j d j l_{k,j} := \frac{r_j}{d_j} l k , j : = d j r j d k : = a k , k − ∑ q = 1 k − 1 l k , q ⋅ r q d_k := a_{k,k} - \sum_{q=1}^{k-1} l_{k,q} \cdot r_q d k : = a k , k − ∑ q = 1 k − 1 l k , q ⋅ r q
W faktycznej implementacji przechowywanie danych wyglądałoby następująco:
macierze trójkątne L , A L, A L , A n × n n\times n n × n
macierz diagonalną D D D
ciąg r j r_j r j
1.3. Ostatecznie
Czyli mając A x = b
Ax = b
A x = b L D L T x = b
L D L^T x = b
L D L T x = b { z = L T x O ( 1 2 n 2 ) D z = y O ( n ) L y = b O ( 1 2 n 2 )
\begin{cases}
z = L^T x & O\left( \frac{1}{2} n^2 \right)\\
D z = y & O\left( n \right)\\
L y = b & O\left( \frac{1}{2} n^2 \right)\\
\end{cases}
⎩ ⎪ ⎪ ⎨ ⎪ ⎪ ⎧ z = L T x D z = y L y = b O ( 2 1 n 2 ) O ( n ) O ( 2 1 n 2 )
Złożoność obliczeniowa: O ( 1 6 n 3 ) + O ( n 2 ) + O ( n ) O\left( \frac{1}{6} n^3 \right) + O\left( n^2 \right) + O(n) O ( 6 1 n 3 ) + O ( n 2 ) + O ( n )
2. Rozkład Choleskyego
Mamy
A x = b
Ax = b
A x = b
i tak jak wcześniej macierz A A A A A A dodatnio określona .
2.1. DEF: Macierz symetryczna dodatnio określona
Macierzą symetryczną A T = A ∈ R n × n A^T = A \in \mathbb{R}^{n\times n} A T = A ∈ R n × n x ∈ R n x \in \mathbb{R}^n x ∈ R n x T A x > 0 x^T A x > 0 x T A x > 0
2.2. Kolejny zmodyfikowany rozkład L U LU L U
Wprowadzamy małą modyfikację do wcześniej określonego rozkładu : A = L U = L D L T = L D ⏞ L ^ D L T ⏞ L ^ T = L ^ L ^ T
A = L U = L D L^T = \overbrace{L \sqrt{D}}^{\hat{L}} \overbrace{\sqrt{D} L^T}^{\hat{L}^T} = \hat{L} \hat{L}^T
A = L U = L D L T = L D L ^ D L T L ^ T = L ^ L ^ T
gdzie
D = diag ( d i ) D = \operatorname{diag}(d_i) D = d i a g ( d i ) D = diag ( d i ) \sqrt{D} = \operatorname{diag}(\sqrt{d_i}) D = d i a g ( d i )
(robimy pierwiastki, więc potrzebujemy macierzy dodatnio określonych)
Implementacja
Mamy: [ l ^ 11 l ^ 21 l ^ 22 ⋮ ⋱ l ^ n 1 l ^ n 2 … l ^ n n ] × [ l ^ 11 l ^ 21 … l ^ n 1 l ^ 22 … l ^ n 2 ⋱ ⋮ l ^ n n ] = [ a 11 a 21 … a n 1 a 21 a 22 … a n 2 ⋮ ⋮ a n 1 a n 2 … a n n ]
\begin{bmatrix}
\hat{l}_{11}\\
\hat{l}_{21} & \hat{l}_{22}\\
\vdots && \ddots\\
\hat{l}_{n1} & \hat{l}_{n2} & \dots &\hat{l}_{nn}
\end{bmatrix}
\times
\begin{bmatrix}
\hat{l}_{11} & \hat{l}_{21} & \dots & \hat{l}_{n1}\\
& \hat{l}_{22} & \dots & \hat{l}_{n2}\\
&& \ddots & \vdots\\
&&& \hat{l}_{nn}\\
\end{bmatrix}
=
\begin{bmatrix}
a_{11} & a_{21} & \dots & a_{n1}\\
a_{21} & a_{22} & \dots & a_{n2}\\
\vdots &&& \vdots\\
a_{n1} & a_{n2} & \dots & a_{nn}\\
\end{bmatrix}
⎣ ⎢ ⎢ ⎢ ⎢ ⎡ l ^ 1 1 l ^ 2 1 ⋮ l ^ n 1 l ^ 2 2 l ^ n 2 ⋱ … l ^ n n ⎦ ⎥ ⎥ ⎥ ⎥ ⎤ × ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ l ^ 1 1 l ^ 2 1 l ^ 2 2 … … ⋱ l ^ n 1 l ^ n 2 ⋮ l ^ n n ⎦ ⎥ ⎥ ⎥ ⎥ ⎤ = ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ a 1 1 a 2 1 ⋮ a n 1 a 2 1 a 2 2 a n 2 … … … a n 1 a n 2 ⋮ a n n ⎦ ⎥ ⎥ ⎥ ⎥ ⎤
I podobnie jak wcześniej mamy równania postaci:
l ^ 11 ⋅ l ^ 11 = a 11 \hat{l}_{11} \cdot \hat{l}_{11} = a_{11} l ^ 1 1 ⋅ l ^ 1 1 = a 1 1 l ^ i 1 l ^ 11 = a i 1 i = 2 , … , n \hat{l}_{i1} \hat{l}_{11} = a_{i1} \quad i = 2,\dots,n l ^ i 1 l ^ 1 1 = a i 1 i = 2 , … , n
l ^ 21 2 + l ^ 22 2 = a 22 \hat{l}_{21}^2 + \hat{l}_{22}^2 = a_{22} l ^ 2 1 2 + l ^ 2 2 2 = a 2 2 l ^ i 1 ⋅ l ^ 21 + l ^ i 2 ⋅ l ^ 22 = a i 2 i = 3 , … , n \hat{l}_{i1} \cdot \hat{l}_{21} + \hat{l}_{i2} \cdot \hat{l}_{22} = a_{i2} \quad i = 3,\dots,n l ^ i 1 ⋅ l ^ 2 1 + l ^ i 2 ⋅ l ^ 2 2 = a i 2 i = 3 , … , n …
(j j j
l ^ j 1 2 + l ^ j 2 2 + ⋯ + l ^ j j 2 = a j j \hat{l}_{j1}^2 + \hat{l}_{j2}^2 + \dots + \hat{l}_{jj}^2 = a_{jj} l ^ j 1 2 + l ^ j 2 2 + ⋯ + l ^ j j 2 = a j j l ^ i 1 ⋅ l ^ j 1 + l ^ i 2 ⋅ l ^ j 2 + ⋯ + l ^ i j ⋅ l ^ j j = a i j i = j + 1 , … , n \hat{l}_{i1} \cdot \hat{l}_{j1} + \hat{l}_{i2} \cdot \hat{l}_{j2} + \dots + \hat{l}_{ij} \cdot \hat{l}_{jj} = a_{ij} \quad i = j+1,\dots,n l ^ i 1 ⋅ l ^ j 1 + l ^ i 2 ⋅ l ^ j 2 + ⋯ + l ^ i j ⋅ l ^ j j = a i j i = j + 1 , … , n
Program:
for j : = 1 j := 1 j : = 1 to n n n
l ^ j j : = a j j − ∑ k = 1 j − 1 l ^ j k 2 \hat{l}_{jj} := \sqrt{a_{jj} - \sum_{k=1}^{j-1} \hat{l}_{jk}^2} l ^ j j : = a j j − ∑ k = 1 j − 1 l ^ j k 2 for i : = j + 1 i := j+1 i : = j + 1 to n n n
l ^ i j : = a i j − ∑ k = 1 j − 1 l ^ i k − l ^ j k l ^ j j \hat{l}_{ij} := \frac{a_{ij} - \sum_{k=1}^{j-1} \hat{l}_{ik} - \hat{l}_{jk}}{\hat{l}_{jj}} l ^ i j : = l ^ j j a i j − ∑ k = 1 j − 1 l ^ i k − l ^ j k
W pamięci przechowujemy:
macierze trójkątne A A A L ^ \hat{L} L ^ n × n n \times n n × n
Dostajemy układ równań:
{ L ^ y = b O ( 1 2 n 2 ) L ^ x = y O ( 1 2 n 2 )
\begin{cases}
\hat{L} y = b & O\left( \frac{1}{2}n^2 \right)\\
\hat{L} x = y & O\left( \frac{1}{2}n^2 \right)\\
\end{cases}
{ L ^ y = b L ^ x = y O ( 2 1 n 2 ) O ( 2 1 n 2 )
Złożoność obliczeniowa: O ( 1 6 n 3 ) + O ( n 2 ) O\left(\frac{1}{6} n^3\right) + O(n^2) O ( 6 1 n 3 ) + O ( n 2 )