대칭행렬의 LDU 분해
📂행렬대수 대칭행렬의 LDU 분해 개요 L T L^{T} L T 는 상삼각행렬이므로 A = L U A = LU A = LU 에서의 U U U 를 U : = D L T U:= DL^{T} U := D L T 으로 바꾼다고 보면 된다. 일반적인 LU 분해 보다 조건이 까다로워진만큼 계산량은 많이 줄어든다.
정리 가역대칭행렬 A A A 에 대해 A = L D L T A = LDL^{T} A = L D L T 를 만족하는 하삼각행렬 L L L 과 대각행렬 D D D 가 존재한다.
증명 A A A 는 가역행렬 이므로 A = L U A = L U A = LU 를 만족하는 하삼각행렬 L L L 과 상삼각행렬 U U U 가 존재한다. 한편 A = A T A = A^{T} A = A T 이므로
L U = A = A T = U T L T
L U = A = A^{T} = U^{T} L^{T}
LU = A = A T = U T L T
LU 분해에서 L L L 은 대각선 성분이 모두 1 1 1 인 하삼각행렬로 구해지므로 det L = 1 \det L = 1 det L = 1 , 즉 역행렬 L − 1 L^{-1} L − 1 이 존재한다. L − 1 L^{-1} L − 1 는 하삼각행렬이고
L U = U T L T
L U = U^{T} L^{T}
LU = U T L T
을 U U U 에 대해 정리하면
U = L − 1 U T L T
U = L^{-1} U^{T} L^{T}
U = L − 1 U T L T
양변의 오른쪽에 ( L T ) − 1 (L^{T})^{-1} ( L T ) − 1 을 곱하면
U ( L T ) − 1 = L − 1 U T
U (L^{T})^{-1} = L^{-1} U^{T} U ( L T ) − 1 = L − 1 U T
을 얻는다. U ( L T ) − 1 = L − 1 U T U (L^{T})^{-1} = L^{-1} U^{T} U ( L T ) − 1 = L − 1 U T 의 좌변은 상삼각행렬이고 우변은 하삼각행렬이므로
D : = U ( L T ) − 1 = L − 1 U T
D: = U (L^{T})^{-1} = L^{-1} U^{T}
D := U ( L T ) − 1 = L − 1 U T
는 대각행렬 임을 알 수 있다. 이를 U U U 에 대해서 정리하면 U = D L T U = DL^{T} U = D L T 이고 다음을 얻는다.
A = L U = L D L T
A = LU = LDL^{T}
A = LU = L D L T
■
계산 이제 구체적으로 L L L 과 D D D 를 계산하기 위해 아래와 같이 정의해보자.
L : = [ 1 0 0 ⋯ 0 l 21 1 0 ⋯ 0 l 31 l 32 1 ⋯ 0 ⋮ ⋮ ⋮ ⋱ ⋮ l n 1 l n 2 l n 3 ⋯ 1 ] , D : = [ d 11 0 0 ⋯ 0 0 d 22 0 ⋯ 0 0 0 d 33 ⋯ 0 ⋮ ⋮ ⋮ ⋱ ⋮ 0 0 0 ⋯ d n n ]
L : = \begin{bmatrix}
1 & 0 & 0 & \cdots & 0
\\ l_{21} & 1 & 0 & \cdots & 0
\\ l_{31} & l_{32} & 1 & \cdots & 0
\\ \vdots & \vdots & \vdots & \ddots & \vdots
\\ l_{n1} & l_{n2} & l_{n3} & \cdots & 1
\end{bmatrix},\quad
D : = \begin{bmatrix}
d_{11} & 0 & 0 & \cdots & 0
\\ 0 & d_{22} & 0 & \cdots & 0
\\ 0 & 0 & d_{33} & \cdots & 0
\\ \vdots & \vdots & \vdots & \ddots & \vdots
\\ 0 & 0 & 0 & \cdots & d_{nn}
\end{bmatrix}
L := 1 l 21 l 31 ⋮ l n 1 0 1 l 32 ⋮ l n 2 0 0 1 ⋮ l n 3 ⋯ ⋯ ⋯ ⋱ ⋯ 0 0 0 ⋮ 1 , D := d 11 0 0 ⋮ 0 0 d 22 0 ⋮ 0 0 0 d 33 ⋮ 0 ⋯ ⋯ ⋯ ⋱ ⋯ 0 0 0 ⋮ d nn
A = ( a i j ) = L D L T A = (a_{ij}) = LDL^{T} A = ( a ij ) = L D L T 를 직접 계산해보자.
l i j > 0 l_{ij} > 0 l ij > 0 이고 l i i = 1 l_{ii} = 1 l ii = 1 이므로
a i j = ∑ ν = 1 n l i ν ∑ μ = 1 n d ν μ l μ j T = ∑ ν = 1 n l i ν ∑ μ = 1 n δ ν μ d ν μ l j μ = ∑ ν = 1 n l i ν d ν ν l j ν = ∑ ν = 1 min ( i , j ) l i ν d ν ν l j ν
\begin{align*}
a_{ij} =& \sum_{\nu = 1}^{n} l_{i \nu} \sum_{\mu = 1}^{n} d_{\nu \mu} l_{\mu j}^{T}
\\ =& \sum_{\nu = 1}^{n} l_{i \nu} \sum_{\mu = 1}^{n} \delta_{\nu \mu} d_{\nu \mu} l_{ j \mu}
\\ =& \sum_{\nu = 1}^{n} l_{i \nu} d_{\nu \nu} l_{ j \nu}
\\ =& \sum_{\nu = 1}^{ \min (i,j) } l_{i \nu} d_{\nu \nu} l_{ j \nu}
\end{align*}
a ij = = = = ν = 1 ∑ n l i ν μ = 1 ∑ n d νμ l μ j T ν = 1 ∑ n l i ν μ = 1 ∑ n δ νμ d νμ l j μ ν = 1 ∑ n l i ν d νν l j ν ν = 1 ∑ m i n ( i , j ) l i ν d νν l j ν
우선 j ≤ i j \le i j ≤ i 인 경우로 가정해보자.
a i j = ∑ ν = 1 j l i ν d ν ν l j ν = ∑ ν = 1 j − 1 l i ν d ν ν l j ν + l i j d j j l j j = ∑ ν = 1 j − 1 l i ν d ν ν l j ν + l i j d j j
a_{ij} = \sum_{\nu = 1}^{j} l_{i \nu} d_{\nu \nu} l_{j \nu} = \sum_{\nu = 1}^{j-1} l_{i \nu} d_{\nu \nu} l_{j \nu} + l_{ij} d_{jj} l_{jj} = \sum_{\nu = 1}^{j-1} l_{i \nu} d_{\nu \nu} l_{j \nu} + l_{ij} d_{jj}
a ij = ν = 1 ∑ j l i ν d νν l j ν = ν = 1 ∑ j − 1 l i ν d νν l j ν + l ij d jj l jj = ν = 1 ∑ j − 1 l i ν d νν l j ν + l ij d jj
여기서 i = j = k i = j = k i = j = k 인 경우
a k k = ∑ ν = 1 k − 1 l k ν d ν ν l k ν + d k k d k k
a_{kk} = \sum_{\nu = 1}^{k-1} l_{k \nu} d_{\nu \nu} l_{k \nu} + d_{kk}
d_{kk}
a kk = ν = 1 ∑ k − 1 l k ν d νν l k ν + d kk d kk
에 대해 정리하면
d k k = a k k − ∑ ν = 1 k − 1 d ν ν l k ν 2
d_{kk} = a_{kk} - \sum_{\nu = 1}^{k-1} d_{\nu \nu} l_{k \nu}^2
d kk = a kk − ν = 1 ∑ k − 1 d νν l k ν 2
반대로 j = k j = k j = k 이고 i = k + 1 , ⋯ , n i = k+1, \cdots , n i = k + 1 , ⋯ , n 이라고 하면
a i k = ∑ ν = 1 k − 1 l i ν d ν ν l k ν + l i k d k k
a_{ik} = \sum_{\nu = 1}^{k-1} l_{i \nu} d_{\nu \nu} l_{k \nu} + l_{ik} d_{kk}
a ik = ν = 1 ∑ k − 1 l i ν d νν l k ν + l ik d kk
l i k l_{ik} l ik 에 대해 정리하면
l i k = 1 d k k { a i k − ∑ ν = 1 k − 1 l i ν d ν ν l k ν }
l_{ik} = {{1} \over {d_{kk}}} \left\{ a_{ik} - \sum_{\nu = 1}^{k-1} l_{i \nu} d_{\nu \nu} l_{k \nu} \right\}
l ik = d kk 1 { a ik − ν = 1 ∑ k − 1 l i ν d νν l k ν }
알고리즘 ( a i j ) ∈ R n × n (a_{ij}) \in \mathbb{R}^{n \times n} ( a ij ) ∈ R n × n 가 가역대칭행렬이라고 하자.
**Step 1. k = 1 k = 1 k = 1
d 11 = a 11 d_{11} = a_{11} d 11 = a 11 을 대입하고 l i 1 = 1 d 11 a i 1 \displaystyle l_{i1} = {{1} \over {d_{11}}} a_{i1} l i 1 = d 11 1 a i 1 을 계산한다.
**Step 2. k = 2 , 3 , ⋯ , n − 1 k = 2, 3, \cdots , n-1 k = 2 , 3 , ⋯ , n − 1
Step 2-1. 다음을 계산한다.
d k k = a k k − ∑ ν = 1 k − 1 d ν ν l k ν 2
d_{kk} = a_{kk} - \sum_{\nu = 1}^{k-1} d_{\nu \nu} l_{k \nu}^2
d kk = a kk − ν = 1 ∑ k − 1 d νν l k ν 2
**Step 2-2. i = k + 1 , k + 2 , ⋯ , n − 1 i = k+1 , k+2 , \cdots , n-1 i = k + 1 , k + 2 , ⋯ , n − 1 에 대해 다음을 계산한다.
l i k = 1 d k k { a i k − ∑ ν = 1 k − 1 l i ν d ν ν l k ν }
l_{ik} = {{1} \over {d_{kk}}} \left\{ a_{ik} - \sum_{\nu = 1}^{k-1} l_{i \nu} d_{\nu \nu} l_{k \nu} \right\}
l ik = d kk 1 { a ik − ν = 1 ∑ k − 1 l i ν d νν l k ν }
**Step 3. k = n k = n k = n 에 대해 다음을 계산한다.
d n n = a n n − ∑ ν = 1 n − 1 d ν ν l n ν 2
d_{nn} = a_{nn} - \sum_{\nu = 1}^{n-1} d_{\nu \nu} l_{n \nu}^2
d nn = a nn − ν = 1 ∑ n − 1 d νν l n ν 2
손으로 계산할 땐 복잡해보였지만 알고리즘으로 정리해서 보면 계산량이 확 줄어들었음을 알 수 있다.언제나 D D D 의 대각성분을 먼저 구한 후 L L L 의 값을 한 열 단위로 구해나간다는 것을 기억하자. 코드 이제 LU 분해에서 했던 것과 같은 예제를 가져와서 LDU 분해를 실제로 해보자.
아래의 코드는 위의 알고리즘을 실제로 R 에서 구현한 것이다.
B <- matrix( c ( 2 , 1 , 0 , 1 , 2 , 1 , 0 , 1 , 2 ) , ncol= 3 ) ; B
LDU<- function ( A)
{
n= length ( A[ , 1 ] )
L<- diag( 1 , n)
D<- diag( 0 , n)
k= 1
D[ 1 , 1 ] = A[ 1 , 1 ]
L[ , 1 ] = A[ , 1 ] / D[ 1 , 1 ]
for ( k in 2 : ( n- 1 ) )
{
D[ k, k] = A[ k, k] - sum ( diag( D) [ 1 : ( k- 1 ) ] * L[ k, 1 : ( k- 1 ) ] ^ 2 )
for ( i in ( k+ 1 ) : n)
{
L[ i, k] = ( A[ i, k] - sum ( L[ i, 1 : ( k- 1 ) ] * diag( D) [ 1 : ( k- 1 ) ] * L[ k, 1 : ( k- 1 ) ] ) ) / D[ k, k]
}
}
D[ n, n] = A[ n, n] - sum ( diag( D) [ 1 : ( n- 1 ) ] * L[ n, 1 : ( n- 1 ) ] ^ 2 )
return ( list ( L= L, D= D) )
}
LDU( B)
LDU( B) $ L%*% LDU( B) $ D%*% t( LDU( B) $ L)
실행 결과는 다음과 같다.
R 에는 따로 LDU 분해가 구현되어 있지 않아 검산만 해보았는데, 정확하게 분해가 되었음을 알 수 있다.
B B B 의 사이즈가 그렇게 크지 않고 계산도 그럭저럭 할만한 편이므로 정 의심된다면 직접 손계산을 해봐도 좋다.