행렬방정식 - 일반행렬 with Scipy

Posted on | In | 
1
from scipy import linalg

Determinant 구하기

linalg.det(A)

기본 알고리즘 : LU Decomposition

→ A = LU \(det(A) = det(L * U) = det(L) * det(U) = det(U)\) : U의 대각항들만 모두 곱해주는 형태로 구해짐

Lapack 함수 (subroutine)

  1. zgetrf : complex128
  2. dgetrf : float64

​

역행렬 구하기

목적이 행렬방정식 Ax=b를 구하는 거라면 사용하지 말자!

linalg.inv(A) > 역행렬이 없으면 singular matrix 에러 출현! : 수학적으로는 에러가 아니지만, 수치적으로 에러에 기인해 판단될 수도 있음 (구분이 불가함)

기본 알고리즘

LU Decomposition : Solve \(LUA^{-1} = I\) : Backward phase (backsubstitution)

Lapack

  • getrf : LU decomposition
  • getri : inverse from triangular matrix : 트라이앵귤러 매트릭스를 푸는 solver 역할을 겸함

​

Ax = b 풀기

linalg.solve(A, b, assume_a="gen")

assume_a : 입력으로 들어오는 행렬의 특성

설정을 잘못 넣어도 별다른 오류 출력이 없음

  • gen : 행렬의 특성을 모를 때 LU decomposition - gesv
  • sym : (Real or Complex) Symmetric Metrix \(A = A^T\) Diagonal pivoting \(A = LDL^T\) : block diagonal - sysv
  • her : Hermitian matrix \(A = A^*\) Diagonal Pivoting \(A = LDL^*\) - hesv
  • pos : Positive Definite \(x^TAx > 0\) Symmetric or Hermitian Cholesky decomposition : posv

​

Triangular Matrix Solver

Ax = b, A = Upper(lower) triangular matrix

Backward phase (backsubstitution) 과정만이 필요함

linalg.solve_triangular(A, b, {lower = False})

Lapack

trtrs

​

해의 정확성

Ax = b, x = 수치적 계산으로 근사된 값

  • Q1. Ax와 b가 충분히 비슷한가 ?
  • Q2. Ax - b 충분히 작은가? 0에 충분히 가까운가?

엔트리의 오차를 직접 계산할 수 있음 (절대 / 상대 오차 모두)

np.allclose(A@x, b) 변수의 모든 엔트리를 비교해 충분히 비슷하다면 True → np.allclose(A@x - b, np.zeros((b.shape[0], )))