목록으로
Chapter 7: 선형시스템
진행률: 7 / 8
예상 시간: 25분난이도: 중급
선형시스템
학습 목표
- 연립방정식의 다양한 해법 마스터
- LU 분해와 효율적인 계산
- 반복법과 수치적 안정성
선형시스템 Ax = b
유일해
det(A) ≠ 0
rank(A) = n
무수히 많은 해
rank(A) = rank([A|b]) < n
해 없음
rank(A) < rank([A|b])
가우스 소거법
1
전진 소거
상삼각 행렬로 변환
2
후진 대입
아래에서 위로 해 계산
계산 복잡도: O(n³)
LU 분해
PA = LU
- • P: 순열 행렬 (피벗팅)
- • L: 하삼각 행렬
- • U: 상삼각 행렬
장점
- • 여러 b에 대해 효율적
- • 행렬식 계산 용이
- • 역행렬 계산 효율적
해법 과정
- 1. Ly = Pb (전진 대입)
- 2. Ux = y (후진 대입)
반복법
야코비 방법
x_i^(k+1) = (b_i - Σ a_ij x_j^(k)) / a_ii
가우스-자이델
새 값을 즉시 사용
더 빠른 수렴
SOR 방법
과이완 계수 ω 사용
최적 수렴
선형시스템 솔버 구현
import numpy as np
from scipy.linalg import lu, solve_triangular
# 가우스 소거법
def gaussian_elimination(A, b):
"""가우스 소거법으로 Ax = b 해결"""
n = len(A)
Ab = np.column_stack([A.copy(), b.copy()])
# 전진 소거
for i in range(n):
# 피벗팅
max_row = i + np.argmax(np.abs(Ab[i:, i]))
Ab[[i, max_row]] = Ab[[max_row, i]]
# 소거
for j in range(i+1, n):
factor = Ab[j, i] / Ab[i, i]
Ab[j, i:] -= factor * Ab[i, i:]
# 후진 대입
x = np.zeros(n)
for i in range(n-1, -1, -1):
x[i] = (Ab[i, -1] - np.dot(Ab[i, i+1:n], x[i+1:])) / Ab[i, i]
return x
# LU 분해 솔버
def lu_solver(A, b):
"""LU 분해를 이용한 선형시스템 해결"""
P, L, U = lu(A)
# Ly = Pb
y = solve_triangular(L, P @ b, lower=True)
# Ux = y
x = solve_triangular(U, y)
return x
# 야코비 반복법
def jacobi_method(A, b, x0=None, tol=1e-10, max_iter=1000):
"""야코비 반복법"""
n = len(A)
x = x0 if x0 is not None else np.zeros(n)
for iteration in range(max_iter):
x_new = np.zeros(n)
for i in range(n):
sum_ax = sum(A[i, j] * x[j] for j in range(n) if j != i)
x_new[i] = (b[i] - sum_ax) / A[i, i]
# 수렴 확인
if np.linalg.norm(x_new - x) < tol:
print(f"수렴 (반복: {iteration+1})")
return x_new
x = x_new
print("최대 반복 횟수 도달")
return x
# 가우스-자이델 방법
def gauss_seidel(A, b, x0=None, tol=1e-10, max_iter=1000):
"""가우스-자이델 반복법"""
n = len(A)
x = x0 if x0 is not None else np.zeros(n)
for iteration in range(max_iter):
x_old = x.copy()
for i in range(n):
sum1 = sum(A[i, j] * x[j] for j in range(i))
sum2 = sum(A[i, j] * x_old[j] for j in range(i+1, n))
x[i] = (b[i] - sum1 - sum2) / A[i, i]
# 수렴 확인
if np.linalg.norm(x - x_old) < tol:
print(f"수렴 (반복: {iteration+1})")
return x
print("최대 반복 횟수 도달")
return x
# 조건수 분석
def condition_analysis(A):
"""행렬의 조건수 분석"""
cond = np.linalg.cond(A)
print(f"조건수: {cond:.2e}")
if cond < 10:
print("매우 안정적")
elif cond < 100:
print("안정적")
elif cond < 1000:
print("주의 필요")
else:
print("불안정 (ill-conditioned)")
return cond
# 예제
A = np.array([[10, -1, 2],
[-1, 11, -1],
[2, -1, 10]])
b = np.array([6, 25, -11])
print("가우스 소거법:", gaussian_elimination(A, b))
print("LU 분해:", lu_solver(A, b))
print("야코비 방법:", jacobi_method(A, b))
print("가우스-자이델:", gauss_seidel(A, b))
print("\n조건수 분석:")
condition_analysis(A)