LA-3) LU분해 feat. 조건수,특이값분해 [R]

 

아래는 위의 예시의 R코드입니다. 패키지 자체에 대한 함수로도 존재합니다.

# LU분해
# sol1)
A <- matrix(c(1,-2,0,
              -1,3,-1,
              0,-1,2),byrow = T,ncol=3)
E1 <- matrix(c(1,0,0,
               1,1,0,
               0,0,1),byrow = T,ncol=3)
E2 <- matrix(c(1,0,0,
               0,1,0,
               0,1,1),byrow = T,ncol=3)
U <- E2%*%E1%*%A
L <- solve(E1)%*%solve(E2)
L%*%U

# sol2)
#library(matlib)
LU(A)
> LU(A)
$P
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1

$L
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]   -1    1    0
[3,]    0   -1    1

$U
     [,1] [,2] [,3]
[1,]    1   -2    0
[2,]    0    1   -1
[3,]    0    0    1

LU분해를 이용한 연립방정식 풀이1

LU분해를 이용한 연립방정식 풀이 2

LU분해를 이용한 연립방정식 풀이 3

아래는 ex2)에 대한 R코드입니다.

B <- matrix(c(1,-2,0,
              -1,3,-1,
              0,-1,2),byrow = T,ncol=3)
E1 <- matrix(c(1,0,0,
               1,1,0,
               0,0,1),byrow = T,ncol=3)
E2 <- matrix(c(1,0,0,
               0,1,0,
               0,1,1),byrow = T,ncol=3)
U <- E2%*%E1%*%B
L <- solve(E1)%*%solve(E2)
#LY=B의 Y값
gaussianElimination(A=L,B = matrix(c(3,-6,5),byrow=T,ncol=1))

#UX =Y의 X값
gaussianElimination(A=U,B = matrix(c(2,-1,1),byrow=T,ncol=1))

추가적으로 LU분해 속도가 확실히 더 빠릅니다.

특이값 분해

조건수 구하기 위한 특이값 R코드입니다.

# 행렬 생성####
A <- matrix(c(1, 2, 3, 4), nrow = 2, byrow = TRUE)

# SVD
svd_result <- svd(A)

> svd(A)
$d
[1] 5.4649857 0.3659662

$u
           [,1]       [,2]
[1,] -0.4045536 -0.9145143
[2,] -0.9145143  0.4045536

$v
           [,1]       [,2]
[1,] -0.5760484  0.8174156
[2,] -0.8174156 -0.5760484

# 특이값
singular_values <- svd_result$d

#실제로 특이값 식이 맞는지 확인
svd_result$u%*%(matrix(c(singular_values[1],0,0,singular_values[2]),byrow = T,ncol = 2))%*%t(svd_result$v)
> svd_result$u%*%(matrix(c(singular_values[1],0,0,singular_values[2]),byrow = T,ncol = 2))%*%t(svd_result$v)
     [,1] [,2]
[1,]    1    2
[2,]    3    4

여기서 우리가 중요하게 볼 것은 특이값이므로,

> singular_values
[1] 5.4649857 0.3659662

 

> singular_values[1]*svd(A_inverse)$d[1]
[1] 14.93303

참고로 역행렬 조건수도 일반행렬 조건수와 동일합니다.

 

조건수가 1에 근접하면 오차에 민감한 변화가 작고 조건수가 크면 작은 오차에 민감한 변화가 큰데

이에 대한 대표적인 행렬들이 조건수가 1인 단위행렬과 조건수가 매우 큰 힐버트행렬입니다.

힐버트 행렬은 1,1 요소자리에서 부터 우측으로 가면서 각 요소자리 합-1분의 1만큼 요소를 가집니다.

#library(pracma)
B <- c(1,1,1 ,1) %>% as.matrix()
B_prime <- c(1,1,1 ,1+0.0001) %>% as.matrix()

H <- matrix(0.0001,ncol=4,nrow=4)

gaussianElimination(A=hilb(4),B = B)
gaussianElimination(A=hilb(4)+H,B=B)
gaussianElimination(A=hilb(4)+diag(0.0001,4),B=B)

#조건수 체크 : 대각행렬에 0.0001을 더한 조건수가 기존 조건수와 매우 달라짐.
> kappa(hilb(4))
[1] 20093.93
> kappa(hilb(4)+diag(0.0001,4))
[1] 9677.346
> kappa(hilb(4)+H)
[1] 20019.34
> gaussianElimination(A=hilb(4),B = B)
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    0    0    0   -4
[2,]    0    1    0    0   60
[3,]    0    0    1    0 -180
[4,]    0    0    0    1  140
> gaussianElimination(A=hilb(4)+H,B=B)
     [,1] [,2] [,3] [,4]       [,5]
[1,]    1    0    0    0   -3.99361
[2,]    0    1    0    0   59.90415
[3,]    0    0    1    0 -179.71246
[4,]    0    0    0    1  139.77636
> gaussianElimination(A=hilb(4)+diag(0.0001,4),B=B)
     [,1] [,2] [,3] [,4]        [,5]
[1,]    1    0    0    0  -0.5889767
[2,]    0    1    0    0  21.1225671
[3,]    0    0    1    0 -85.7591250
[4,]    0    0    0    1  78.4565082

 

조건수가 커야 오차가 가우스소거법과 LU분해의 오차 차이가 확 크기 때문에 

큰 행렬을 이용했다.

library(matlib)
set.seed(8457)

n <- 100
A <- matrix(rnorm(n * n), n, n)
b <- rnorm(n)

# 행렬의 조건수 계산
cond_A <- kappa(A)
paste0("행렬 A의 조건수:", cond_A)

# 가우스 소거법을 사용하여 해를 구함
x_gaussian <- gaussianElimination(A = A,B = b)[,101]

# LU 분해를 사용하여 Ax = b를 푸는 함수
solve_lu <- function(A, b) {
  lu <- lu(A)
  L <- lu$L
  U <- lu$U
  y <- forwardsolve(L, b) #LY = b
  x <- backsolve(U, y) # UX = Y
  x
}

# LU 분해를 사용하여 해를 구함
x_lu <- solve_lu(A, b)


# 실제 해 계산
x_true <- solve(A, b)

# 오차 계산
error_gaussian <- norm(x_true - x_gaussian, type = "2")
error_lu <- norm(x_true - x_lu, type = "2")

paste0("가우스 소거법에 의한 해의 오차:", error_gaussian)
paste0("LU 분해에 의한 해의 오차:", error_lu)
> paste0("가우스 소거법에 의한 해의 오차:", error_gaussian)
[1] "가우스 소거법에 의한 해의 오차:2.904426349253e-08"
> paste0("LU 분해에 의한 해의 오차:", error_lu)
[1] "LU 분해에 의한 해의 오차:6.55434647136738e-12"

 

 

 

행렬 수 200x200으로 가우스소거법을 진행하니 R이 터져버렸다.. (가우스 소거법의 문제점)

LU는 연산이 매우 간단하기에 빠르게 가능했고 이에 대한 오차도 작습니다.

 

##LU분해와 가우스소거법들의 오차검증이 제대로 확인되지 않았습니다.

그렇기에 아래의 오차내용은 무시해주세요.

 

##확실한 것은 LU분해가 행렬 사이즈가 클 수록 명백하게 가우스소거법보다 해를 찾아내는 속도가 빠른다는 점입니다.