目錄
- Gauss 消去法的矩陣意義
- 矩陣的 LU 分解
- 其它的三角分解
-
- LDU 分解
- 對稱正定矩陣
- 帶狀矩陣的分解
-
- 三對角陣的分解
- 矩陣分解的應用
-
- 1.求解方程組 $Ax=b$ 的解
- 2.求矩陣行列式
- 3.求矩陣的逆
- 4.求解成組方程組
Gauss 消去法的矩陣意義
注意：這裡面讨論的不是列主元高斯消去法，是以算法可能不穩定。
根據《線性代數》我們知道，對于矩陣 A 行初等變換相當于左乘相應初等矩陣，對于矩陣 A 列初等變換相當于右乘相應初等矩陣。
由此可知：
A = L U A=LU A=LU
即 A 可以表示成為一個機關下三角陣 L 與一個上三角陣 U 的乘積
這就是矩陣的 LU 分解,也稱為三角分解或 Doolittle 分解
矩陣的 LU 分解
那是不是所有的方陣都能夠進行 LU 分解呢，如果不是，那又是哪些能夠進行 LU 分解呢？
聯系前面我們知道 Gauss 消元法順利進行的條件是各個消元的主元不為 0 （即 a k k ( k − 1 ) ≠ 0 a_{kk}^{(k-1)} \neq 0 akk(k−1)̸=0 ）。我們知道，能夠進行 Gauss 消元與能夠進行 LU 分解是等價的， Gauss 消元法也可以用來求行列式。
若将矩陣 A 的左上方各階順序主子矩陣記為 A 1 , A 2 , . . . , A n A_1,A_2,...,A_n A1,A2,...,An 。使用高斯消元法使得第 k 個消元的主元如果為 0 ，那麼 d e t ( A k ) = 0 det(A_k)=0 det(Ak)=0 ，是以有：
若 n 階方陣的各階順序主子行列式不為 0 ，則存在唯一的機關下三角矩陣 L 和上三角矩陣 U ，滿足 A=LU
由此定理可以得到 LU 分解的計算公式：
對于 A 中的每一項，有
是以得到 LU 分解的每項的計算公式：
{ μ 1 j = α 1 j , j = 1 , 2 , . . . , n l i 1 = α i 1 μ 11 , i = 2 , 3 , . . . , n μ i j = α i j − ∑ k = 1 i − 1 l i k u k j , j = i , i + 1 , . . . , n , i = 2 , 3 , . . , n l k i = 1 μ i i ( α k i − ∑ t = 1 i − 1 l k t μ t i ) , k = i + 1 , . . . , n , i = 2 , 3 , . . , n \left\{ \begin{array}{l} {\mu _{1j}} = {\alpha _{1j}}{\rm{ }},j = 1,2,...,n\\ {l_{i1}} = \frac{{{\alpha _{i1}}}}{{{\mu _{11}}}}{\rm{ }},i = 2,3,...,n\\ {\mu _{ij}} = {\alpha _{ij}} - \sum\limits_{k = 1}^{i - 1} {{l_{ik}}{u_{kj}}} {\rm{ }},j = i,i + 1,...,n,i = 2,3,..,n\\ {l_{ki}} = \frac{1}{{{\mu _{ii}}}}({\alpha _{ki}} - \sum\limits_{t = 1}^{i - 1} {{l_{kt}}{\mu _{ti}}} ){\rm{ }},k = i + 1,...,n,i = 2,3,..,n \end{array} \right. ⎩⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎧μ1j=α1j,j=1,2,...,nli1=μ11αi1,i=2,3,...,nμij=αij−k=1∑i−1likukj,j=i,i+1,...,n,i=2,3,..,nlki=μii1(αki−t=1∑i−1lktμti),k=i+1,...,n,i=2,3,..,n
是以要先計算 U 中的第 i 行元素，再計算 L 中的第 i 列元素（因為計算 L 中的第 i 列元素需要知道 u i i u_{ii} uii ）。總的時間複雜度為 O ( n 3 ) O(n^3) O(n3) 。
我們在存儲時，可以将 LU 存儲在一個矩陣中，減小空間消耗：
( a 11 ( 0 ) a 12 ( 0 ) a 13 ( 0 ) ⋯ a 1 n ( 0 ) β 1 ( 0 ) l 21 a 22 ( 1 ) a 23 ( 1 ) ⋯ a 2 n ( 1 ) β 2 ( 1 ) l 31 l 32 a 33 ( 2 ) ⋯ a 3 n ( 2 ) β 3 ( 2 ) ⋱ ⋮ ⋮ l n 1 l n 2 l n 3 ⋯ a n n ( n − 1 ) β n ( n − 1 ) ) {\begin{pmatrix} {a_{11}^{(0)}}&{a_{12}^{(0)}}&{a_{13}^{(0)}}& \cdots &{a_{1n}^{(0)}}&{\beta _1^{(0)}}\\ {{l_{21}}}&{a_{22}^{(1)}}&{a_{23}^{(1)}}& \cdots &{a_{2n}^{(1)}}&{\beta _2^{(1)}}\\ {{l_{31}}}&{{l_{32}}}&{a_{33}^{(2)}}& \cdots &{a_{3n}^{(2)}}&{\beta _3^{(2)}}\\ {}&{}&{}& \ddots & \vdots & \vdots \\ {{l_{n1}}}&{{l_{n2}}}&{{l_{n3}}}& \cdots &{a_{nn}^{(n - 1)}}&{\beta _n^{(n - 1)}} \end{pmatrix}} ⎝⎜⎜⎜⎜⎜⎜⎛a11(0)l21l31ln1a12(0)a22(1)l32ln2a13(0)a23(1)a33(2)ln3⋯⋯⋯⋱⋯a1n(0)a2n(1)a3n(2)⋮ann(n−1)β1(0)β2(1)β3(2)⋮βn(n−1)⎠⎟⎟⎟⎟⎟⎟⎞
故而空間複雜度為 O ( n 2 ) O(n^2) O(n2) 。
其它的三角分解
LDU 分解
對于上三角矩陣 U ，如果 U 的對角元存在 0 ，則 d e t ( U ) = 0 det(U)=0 det(U)=0 。是以當 A 的三角分解存在時，上三角矩陣 U 的全體對角元 u k k ≠ 0 u_{kk}\neq 0 ukk̸=0,若令
D = d i a g ( u 11 , u 22 , . . . , u n n ) D=diag(u_{11},u_{22},...,u_{nn}) D=diag(u11,u22,...,unn)
則可以将 LU 分解轉化為：
A = L U = L D D − 1 U = L D M T A=LU=LDD^{-1}U=LDM^T A=LU=LDD−1U=LDMT
我們知道 D − 1 = d i a g ( 1 u 11 , 1 u 22 , . . . , 1 u n n ) D^{-1}=diag(\frac 1 {u_{11}},\frac 1 {u_{22}},...,\frac 1 {u_{nn}}) D−1=diag(u111,u221,...,unn1), D − 1 U D^{-1}U D−1U 為機關上三角矩陣。
M T M^T MT 為機關上三角矩陣
M = ( 1 m 21 1 m 31 m 32 1 ⋮ ⋮ ⋮ ⋱ m n 1 m n 2 m n 3 ⋯ 1 ) , m j i = μ i j μ i i , ( i = 1 , 2 , . . . , n − 1 ; j = n + 1 , . . . , n ) M = {\begin{pmatrix} 1&{}&{}&{}&{}\\ {{m_{21}}}&1&{}&{}&{}\\ {{m_{31}}}&{{m_{32}}}&1&{}&{}\\ \vdots & \vdots & \vdots & \ddots &{}\\ {{m_{n1}}}&{{m_{n2}}}&{{m_{n3}}}& \cdots &1 \end{pmatrix}} ,{m_{ji}} = \frac{{{\mu _{ij}}}}{{{\mu _{ii}}}},\begin{array}{l} (i = 1,2,...,n - 1;\\ j = n + 1,...,n) \end{array} M=⎝⎜⎜⎜⎜⎜⎛1m21m31⋮mn11m32⋮mn21⋮mn3⋱⋯1⎠⎟⎟⎟⎟⎟⎞,mji=μiiμij,(i=1,2,...,n−1;j=n+1,...,n)
其中 LDU 分解和 LU 分解的條件是一樣的。
我們知道，矩陣的 LU 分解能夠通過計算式計算出來，那麼 L 、 D 、 M 的元的公式也可以通過 LU 的元的公式計算出來。
{ d k = α k k − ∑ p = 1 k − 1 l k p ( d p m k p ) l i k = 1 d k [ α i k − ∑ p = 1 k − 1 l i p ( d p m k p ) ] m i k = 1 d k [ α k i − ∑ p = 1 k − 1 ( l k p d p ) m i p ) ] k = 1 , 2 , . . . , n , i = k + 1 , k + 2 , . . . , n \left\{ \begin{matrix} {d_k} = {\alpha _{kk}} - \sum\limits_{p = 1}^{k - 1} {{l_{kp}}({d_p}{m_{kp}})} {\rm{ }}\\ {l_{ik}} = \frac{1}{{{d_k}}}[{\alpha _{ik}} - \sum\limits_{p = 1}^{k - 1} {{l_{ip}}({d_p}{m_{kp}})]} \\ {m_{ik}} = \frac{1}{{{d_k}}}[{\alpha _{ki}} - \sum\limits_{p = 1}^{k - 1} {({l_{kp}}{d_p}){m_{ip}})]} {\rm{ }}\\ {\rm{ }}k = 1,2,...,n,i = k + 1,k + 2,...,n \end{matrix} \right. ⎩⎪⎪⎪⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎧dk=αkk−p=1∑k−1lkp(dpmkp)lik=dk1[αik−p=1∑k−1lip(dpmkp)]mik=dk1[αki−p=1∑k−1(lkpdp)mip)]k=1,2,...,n,i=k+1,k+2,...,n
時間複雜度為 O ( n 3 ) O(n^3) O(n3) 。
當 A 為對稱矩陣時（即 A T = A A^T=A AT=A)，有 A = L D M T A=LDM^T A=LDMT ，則：
A = A T = M D L T A=A^T=MDL^T A=AT=MDLT,又因為 A 的 LDU 分解具有唯一性，是以：
M = L M=L M=L
是以，當 A 對稱時，有分解式：
A = L D L T A=LDL^T A=LDLT
此時，隻需要計算 L 和 D ，計算量減少約一半。
對稱正定矩陣
任何對稱正定矩陣 A 都存在唯一的下三角矩陣 G,滿足 A = G G T A=GG^T A=GGT
理由如下：
在實際應用中，對稱正定矩陣非常重要和常見。
線上性代數中，正定矩陣等價于順序主子式全大于 0 ，即可以進行 LU 分解或者 LDU 分解。
隻要 A 對稱且可以分解，那麼 A = L D L T A=LDL^T A=LDLT 成立，令 D 1 2 = d i a g ( d 1 , d 2 , . . . , d n , ) D^{\frac 1 2}=diag(\sqrt {d_1},\sqrt {d_2},...,\sqrt {d_n},) D21=diag(d1
,d2
,...,dn
,)
令 G = L D 1 2 G=LD^{\frac 1 2} G=LD21
則 A = L D L T = L D 1 2 D 1 2 L T = G G T A=LDL^T=LD^{\frac 1 2}D^{\frac 1 2}L^T=GG^T A=LDLT=LD21D21LT=GGT
稱為 G 的Cholesky 分解
帶狀矩陣的分解
帶狀矩陣在有限元分析等實際應用場景非常普遍，是以非常重要。
定義：設 A = ( α i j ) A=(\alpha_{ij}) A=(αij) 是 n 階方陣，對于小于 n 的正整數 p 和 q ，當 j>i+p 或 i>j+p 時，有 α i j = 0 \alpha_{ij}=0 αij=0 ，則稱 A 是具有上帶寬 q 和下帶寬 p 的帶狀矩陣。
分解結果：
設 n 階矩陣 A 有 LU 分解， A=LU,如果 A 是上下帶寬分别為 q 和 p 的帶狀矩陣，則 L 是帶寬為 p 的下三角矩陣， U 是帶寬為 q 的上三角矩陣。
計算公式：
{ μ i j = α i j − ∑ t = max { 1 , i − p , j − q } i − 1 l i t μ t j , j = i , i + 1 , . . . , min { n , i + q } l k i = 1 μ i i ( α k i − ∑ t = max { 1 , k − p , i − q } i − 1 l i t μ t j ) , k = i + 1 , . . . , min { n , i + p } \left\{ \begin{array}{l} {\mu _{{\rm{ij}}}} = {\alpha _{ij}} - \sum\limits_{t = \max \{ 1,i - p,j - q\} }^{i - 1} {{l_{it}}{\mu _{tj}},j = i,i + 1,...,\min \{ n,i + q\} } \\ {l_{ki}} = \frac{1}{{{\mu _{ii}}}}({\alpha _{ki}} - \sum\limits_{t = \max \{ 1,k - p,i - q\} }^{i - 1} {{l_{it}}{\mu _{tj}}),k = i + 1,...,\min \{ n,i + p\} } \end{array} \right. ⎩⎪⎪⎪⎨⎪⎪⎪⎧μij=αij−t=max{1,i−p,j−q}∑i−1litμtj,j=i,i+1,...,min{n,i+q}lki=μii1(αki−t=max{1,k−p,i−q}∑i−1litμtj),k=i+1,...,min{n,i+p}
隻需要計算帶狀矩陣中不為 0 的元素，是以運算量比 O ( n 3 ) O(n^3) O(n3) 小得多， p 和 q 固定時，運算量是 O ( n ) O(n) O(n)
三對角陣的分解
當 p=q=1 時，此時的帶狀矩陣稱為三對角陣
T = ( b 1 c 1 a 2 b 2 c 2 a 3 b 3 c 3 ⋱ ⋱ ⋱ a n − 1 b n − 1 c n − 1 a n b n ) {\rm{ }}T = {\begin{pmatrix} {{b_1}}&{{c_1}}&{}&{}&{}&{}\\ {{a_2}}&{{b_2}}&{{c_2}}&{}&{}&{}\\ {}&{{a_3}}&{{b_3}}&{{c_3}}&{}&{}\\ {}&{}& \ddots & \ddots & \ddots &{}\\ {}&{}&{}&{{a_{n - 1}}}&{{b_{n - 1}}}&{{c_{n - 1}}}\\ {}&{}&{}&{}&{{a_n}}&{{b_n}} \end{pmatrix}} {\rm{ }} T=⎝⎜⎜⎜⎜⎜⎜⎛b1a2c1b2a3c2b3⋱c3⋱an−1⋱bn−1ancn−1bn⎠⎟⎟⎟⎟⎟⎟⎞
根據帶狀矩陣分解的定理,上面帶狀矩陣可以被分解為：
L = ( 1 l 2 1 l 3 1 ⋱ ⋱ l n 1 ) U = ( μ 1 r 1 μ 2 r 2 ⋱ ⋱ μ n − 1 r n − 1 μ n ) L = \left( {\begin{matrix} 1&{}&{}&{}&{}\\ {{l_2}}&1&{}&{}&{}\\ {}&{{l_3}}&1&{}&{}\\ {}&{}& \ddots & \ddots &{}\\ {}&{}&{}&{{l_n}}&1 \end{matrix}} \right) U = \left( {\begin{matrix} {{\mu _1}}&{{r_1}}&{}&{}&{}\\ {}&{{\mu _2}}&{{r_2}}&{}&{}\\ {}&{}& \ddots & \ddots &{}\\ {}&{}&{}&{{\mu _{n - 1}}}&{{r_{n - 1}}}\\ {}&{}&{}&{}&{{\mu _n}} \end{matrix}} \right) L=⎝⎜⎜⎜⎜⎛1l21l31⋱⋱ln1⎠⎟⎟⎟⎟⎞U=⎝⎜⎜⎜⎜⎛μ1r1μ2r2⋱⋱μn−1rn−1μn⎠⎟⎟⎟⎟⎞
是以得到計算公式：
{ μ 1 = b 1 , r 1 = c 1 l k = α k μ k − 1 , μ k = b k − l k r k − 1 , r k = c k . k = 2 , 3 , . . . , n \left\{ \begin{matrix} {\mu _1} = {b_1},{r_1} = {c_1}\\ {l_k} = \frac{{{\alpha _k}}}{{{\mu _{k - 1}}}},{\mu _k} = {b_k} - {l_k}{r_{k - 1}},{r_k} = {c_k}.k = 2,3,...,n \end{matrix} \right. {μ1=b1,r1=c1lk=μk−1αk,μk=bk−lkrk−1,rk=ck.k=2,3,...,n
求解三對角方程組： T x = d Tx=d Tx=d,其中 d = ( d 1 , d 2 , . . . , d n ) T d={(d_1,d_2,...,d_n)}^T d=(d1,d2,...,dn)T ，将 T 進行三角分解 T=LU ，則 LUx=d,令 y=Ux ，則 Ly=d,得到計算 x 和 y 的公式：
可以推導出：
追的過程（從前往後）：
{ y 1 = d 1 y k = d k − l k y k − 1 , k = 2 , 3 , . . . , n \left\{ \begin{array}{l} {y_1} = {d_1}\\ {y_k} = {d_k} - {l_k}{y_{k - 1}},k = 2,3,...,n \end{array} \right. {y1=d1yk=dk−lkyk−1,k=2,3,...,n
趕的過程（從後往前）：
{ x n = y n μ n x k = 1 μ k ( y k − r k x k + 1 ) , k = n − 1 , n − 2 , . . . , 2 , 1 \left\{ \begin{array}{l} {x_n} = \frac {y_n} {\mu _n}\\ {x_k} = \frac{1}{{{\mu _k}}}({y_k} - {r_k}{x_{k + 1}}),k = n - 1,n - 2,...,2,1 \end{array} \right. {xn=μnynxk=μk1(yk−rkxk+1),k=n−1,n−2,...,2,1
是以該方法被稱為追趕法。
矩陣分解的應用
若對 A 進行 LU 分解
1.求解方程組 A x = b Ax=b Ax=b 的解
原方程等價于：
{ L y = b U x = y \left\{ \begin{array}{l} Ly = b\\ Ux = y \end{array} \right. {Ly=bUx=y
則可以通過類似三對角方程組追趕法求解。 LU 分解的時間複雜度為 O ( n 3 ) O(n^3) O(n3) ，之後的求解過程的時間複雜度為 O ( n 2 ) O(n^2) O(n2)
2.求矩陣行列式
det ( A ) = det ( L ) det ( U ) \det (A) = \det (L)\det (U) det(A)=det(L)det(U)
因為 L 是機關下三角矩陣，是以
det ( A ) = det ( L ) det ( u ) = det ( u ) = μ 11 μ 22 ⋯ μ n n \det (A) = \det (L)\det (u) = \det (u) = {\mu _{11}}{\mu _{22}} \cdots {\mu _{nn}} det(A)=det(L)det(u)=det(u)=μ11μ22⋯μnn
3.求矩陣的逆
對稱矩陣 A 進行 LDU 分解有 A = L D L T A=LDL^T A=LDLT ， A − 1 = L − T D − 1 L − 1 {A^{ - 1}} = {L^{ - T}}{D^{ - 1}}{L^{ - 1}} A−1=L−TD−1L−1
D − 1 = d i a g ( d 1 − 1 , d 2 − 1 , . . . , d n − 1 ) {D^{ - 1}} = diag(d_1^{ - 1},d_2^{ - 1},...,d_n^{ - 1}) D−1=diag(d1−1,d2−1,...,dn−1)
L − 1 L^{-1} L−1 是機關下三角矩陣，由追趕法可以解出 L − 1 L^{-1} L−1 ，進而求出 A − 1 A^{-1} A−1
4.求解成組方程組
若有一系列方程組
A x = B k ( k = 1 , . . . , m , m ≥ 1 ) Ax = {B_k}{\rm{ (}}k = 1,...,m,m \ge 1) Ax=Bk(k=1,...,m,m≥1)
若 A 存在 LU 分解， A=LU ，
{ L y = B k U x = y \left\{ \begin{array}{l} Ly = {B_k}\\ Ux = y \end{array} \right. {Ly=BkUx=y
對每個 B k {B_k} Bk 分别按照求解方程組的方法就可以求解成組方程組。