Sparse Stencil·Jacobi·Conjugate Gradient

0을 저장하는 대신 연산을 표현한다

경계값이 0인 일차원 Poisson 연산자는 대각 2, 이웃 성분 -1의 삼중대각 행렬입니다. Reaction c를 더하면 대각이 2+c가 됩니다. 내부 미지수가 N개이면 nonzero는 3N-2개뿐입니다. 실험은 반복 곱셈에서 dense 행렬을 만들지 않고 세 점 stencil을 직접 적용합니다.

우변은 모두 1입니다. 별도의 pivoted 직접 해법으로 기준 해를 계산합니다. c=0에서는 1부터 시작하는 index i에 대해 i(N+1-i)/2라는 독립적인 해도 있습니다. 복잡한 시뮬레이션 전에 검산하기 좋은 이산 포물선입니다.

Jacobi의 오차를 고유값으로 읽는다

0에서 시작하는 weighted Jacobi는 다음과 같습니다.

오차 전파 행렬의 spectral radius는 max_j |1-omega lambda_j/(2+c)|입니다. 1보다 작으면 이 대칭 모형의 모든 초기 오차가 수렴합니다. 1보다 크면 증가하는 mode가 가능하지만 특정 초기 오차가 그 mode를 곧바로 크게 드러낸다는 뜻은 아닙니다.

c=0에서 N이 커지면 최소 고유값이 작아지고 매끄러운 오차 mode가 천천히 줄어듭니다. 각 sparse matvec의 비용은 N에 비례해도 충분한 정확도를 얻기 위한 반복 수가 커질 수 있습니다.

CG는 실제 conjugate 방향을 만든다

CG는 SPD 행렬에서 잔차와 이전 방향으로 새로운 탐색 방향을 구성합니다. 실험은 실제 step length·잔차·방향 계수를 계산합니다. 정해진 감쇠 곡선을 대신 그리지 않습니다. 새 stencil 곱으로 다시 계산한 상대 잔차가 10^-12 이하이면 갱신을 멈추고 실제 update 수를 표시합니다.

N=8, c=0, 대칭 우변에서는 이 예의 CG가 네 번에 수렴합니다. 대칭성이 활성 부분공간을 줄이기 때문입니다. 정확 산술에서 N회 이내 종료한다는 결과를 모든 부동소수점 시스템의 무조건적 보장으로 읽으면 안 됩니다.

잔차가 오차 인증이 되려면

오차 e=x-xstar, 잔차 r=b-Ax에 대해 SPD 성질로 다음 상한을 얻습니다.

작은 최소 고유값의 역수가 중요합니다. 잔차가 작다는 사실만으로 해가 충분히 정확하다고 단정할 수 없습니다. 실험은 실제 해 오차·다시 계산한 true residual·그 상한을 함께 보여줍니다. 재귀적으로 갱신한 CG 잔차는 반올림으로 달라질 수 있어 중단 판단에는 새 곱셈을 사용합니다.

직접 풀기: Stencil과 비용 확인

N=4, c=0의 해는 (2,3,3,2)입니다. 삼중대각 행렬을 곱하면 모두 1이 되는지 확인하세요. Nonzero는 10개로 dense의 16개보다 적습니다. N을 늘려 최소 고유값과 Jacobi 오차를 비교하고 omega=1.2에서 초기 오차가 작아 보여도 수렴 보장을 주장할 수 없는 이유를 설명하세요.

더 읽기: SciPy CG 문서.