A = LU 분해 (2)

2023.01.08 - [선형대수학/implement_linearAlgebra] - 소거법을 이용한 Ax = b 의 풀이

소거법을 이용한 Ax = b 의 풀이

원래 행렬 A 는 마지막으로 얻은 행렬 U 와 어떤 관계가 있을까? 

1 단계에서 4*4 크기의 문제는 1행의 배수를 다른 행에서 제거함으로써 3*3 크기의 문제로 축소했다.

위 식에서 우변의 첫 번째 행렬이 행렬 A 에서 제거되었다. 이것이 랭크 1인 행렬이다. 

A = np.array(([1,2,3,4],[6,7,8,9],[7,13,12,11],[13,17,4,3]))
A
>>>
array([[ 1,  2,  3,  4],
       [ 6,  7,  8,  9],
       [ 7, 13, 12, 11],
       [13, 17,  4,  3]])
       
for i in range(A.shape[0]):
  exec(f"l{i}u{i} = []")
  a = A[i][i]
  exec(f"l{i}u{i}.append(A[i])")
  for j in range(i+1, A.shape[0]):
    exec(f"l{i}u{i}.append((A[j][i] / a) * A[i])")
    A[j] = A[j] - (A[j][i] / a) * A[i]
    
l0u0
>>>
[array([1, 2, 3, 4]),
 array([ 6., 12., 18., 24.]),
 array([ 7., 14., 21., 28.]),
 array([13., 26., 39., 52.])]
 
l1u1
>>>
[array([  0,  -5, -10, -15]),
 array([ 0., -1., -2., -3.]),
 array([  0.,  -9., -18., -27.])]
 
l2u2
>>>
[array([  0,   0,  -7, -14]),
 array([  0.,   0., -17., -34.])]
 
l3u3
>>>
[array([ 0,  0,  0, 12])]

A = l0u0 + l1u1 + l2u2 + l3u3

numpy array 로의 변환,

list_a = []

for i in range(A.shape[0]):
  a = []
  for j in range(A.shape[0]-i):
    for k in range(A.shape[1]):
      exec(f"a.append(l{i}u{i}[j][k])")
  list_a.append(a)
  
for i in range(len(list_a)):
  exec(f"l{i}u{i} = list_a[i]")
  exec(f"l{i}u{i} = np.array(l{i}u{i}).reshape(len(list_a)-i,-1)")
  
l0u0
>>>
array([[ 1.,  2.,  3.,  4.],
       [ 6., 12., 18., 24.],
       [ 7., 14., 21., 28.],
       [13., 26., 39., 52.]])
       
l1u1
>>>
array([[  0.,  -5., -10., -15.],
       [  0.,  -1.,  -2.,  -3.],
       [  0.,  -9., -18., -27.]])
       
l2u2
>>>
array([[  0.,   0.,  -7., -14.],
       [  0.,   0., -17., -34.]])
       
l3u3
>>>
array([[ 0,  0,  0, 12]])

네 np array 를 합치기 위한 과정 수정,

for i in range(len(list_a)):
  exec(f"l{i}u{i} = list_a[i]")
  exec(f"l{i}u{i} = np.array(l{i}u{i}).reshape(len(list_a)-i,-1)")
  if i != 0:
    exec(f"l{i}u{i} = np.concatenate((np.zeros((i,4)), l{i}u{i}))")
    
l0u0
>>>
array([[ 1.,  2.,  3.,  4.],
       [ 6., 12., 18., 24.],
       [ 7., 14., 21., 28.],
       [13., 26., 39., 52.]])
       
l1u1
>>>
array([[  0.,   0.,   0.,   0.],
       [  0.,  -5., -10., -15.],
       [  0.,  -1.,  -2.,  -3.],
       [  0.,  -9., -18., -27.]])
       
l2u2
>>>
array([[  0.,   0.,   0.,   0.],
       [  0.,   0.,   0.,   0.],
       [  0.,   0.,  -7., -14.],
       [  0.,   0., -17., -34.]])
       
l3u3
>>>
array([[ 0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.],
       [ 0.,  0.,  0., 12.]])
       
new_A = np.zeros((4,4))
for i in range(len(list_a)):
  exec(f"A = A + l{i}u{i}")
  
new_A
>>>
array([[ 1.,  2.,  3.,  4.],
       [ 6.,  7.,  8.,  9.],
       [ 7., 13., 12., 11.],
       [13., 17.,  4.,  3.]])
       
A
>>>
array([[ 1.,  2.,  3.,  4.],
       [ 6.,  7.,  8.,  9.],
       [ 7., 13., 12., 11.],
       [13., 17.,  4.,  3.]])

합이 동일한 것을 확인할 수 있다.

A = LU 를 얻을 수 있다.

위 코드의 수정

L = np.zeros((A.shape[0], A.shape[1]))
L
>>>
array([[0., 0., 0., 0.],
       [0., 0., 0., 0.],
       [0., 0., 0., 0.],
       [0., 0., 0., 0.]])
       
U = np.zeros((A.shape[0], A.shape[1]))
U[0] = A[0]
U
>>>
array([[1., 2., 3., 4.],
       [0., 0., 0., 0.],
       [0., 0., 0., 0.],
       [0., 0., 0., 0.]])
       
for i in range(A.shape[0]-1):
  a = A[i][i]
  if(a == 0):
    continue
  for j in range(i+1, A.shape[0]):
    L[j][i] = (A[j][i] / a)
    A[j] = A[j] - (A[j][i] / a) * A[i]
    U[j] = A[j]
    
for i in range(A.shape[0]):
  L[i][i] = 1
 
L
>>>
array([[ 1.        ,  0.        ,  0.        ,  0.        ],
       [ 6.        ,  1.        ,  0.        ,  0.        ],
       [ 7.        ,  0.2       ,  1.        ,  0.        ],
       [13.        ,  1.8       ,  2.42857143,  1.        ]])
       
U
>>>
array([[  1.,   2.,   3.,   4.],
       [  0.,  -5., -10., -15.],
       [  0.,   0.,  -7., -14.],
       [  0.,   0.,   0.,  12.]])
       
L.dot(U)
>>>
array([[ 1.,  2.,  3.,  4.],
       [ 6.,  7.,  8.,  9.],
       [ 7., 13., 12., 11.],
       [13., 17.,  4.,  3.]])
       
A
>>>
array([[ 1,  2,  3,  4],
       [ 6,  7,  8,  9],
       [ 7, 13, 12, 11],
       [13, 17,  4,  3]])

LU 의 값이 A 와 동일한 것을 확인할 수 있다.

 

LU 분해 시 소거법의 아이디어를 이용해서 A = LU 를 이끌어 냈다.

즉, 나머지 n-1 개의 방정식에서 x_1 를 제거하여 문제의 크기를 n 에서 n-1 로 줄인다. 이 과정에서는 1 행의 배수를 뺏고, 제거한 행렬의 랭크는 1이다. 

아이디언느 A 의 열을 다루는 대신 U 의 행을 이용하는 것이다.

핵심은 뺄셈하는 행이 추축행이고, 이는 행렬 U 안에 이미 존재하는 것,

이로써 행 교환 없이 A = LU 를 얻었다.