역반복법

수치해석 구조해석 마지막 수정: 2026-07-30 04:17:41

1. 개요[편집]

역반복법 Inverse Iteration
제안	Helmut Wielandt, 1944
반복식	$(A - \sigma I)\mathbf{y}_{k+1} = \mathbf{x}_k$
구하는 것	$\sigma$에 가장 가까운 고유쌍
수렴 속도	선형 (RQI로 확장 시 세제곱)
비용	LU/LDL^T 1회 + 반복당 전후치환
주 무대	모드 해석, 좌굴, QR 후 고유벡터 복원

“행렬이 거의 특이한데 괜찮나요?” — 괜찮다. 그게 목적이다.

역반복법(inverse iteration, 또는 빌란트 반복 Wielandt iteration)은 시프트 $\sigma$ 를 하나 잡고 $(A - \sigma I)^{-1}$ 에 거듭제곱법을 적용해, $\sigma$ 에 가장 가까운 고유값의 고유벡터를 뽑아내는 반복 알고리즘이다. 1944년 헬무트 빌란트가 제안했고, 오늘날에도 모드 해석과 좌굴 해석의 내부 엔진으로, 그리고 QR 알고리즘이 고유값만 뱉었을 때 고유벡터를 복원하는 표준 수단으로 살아 있다.

역행렬을 실제로 만들지는 않는다. 매 반복은 선형계

(A - \sigma I)\,\mathbf{y}_{k+1} = \mathbf{x}_k, \qquad \mathbf{x}_{k+1} = \frac{\mathbf{y}_{k+1}}{\|\mathbf{y}_{k+1}\|}

를 푸는 것이며, $\sigma$ 가 고정이면 $A - \sigma I$ 의 LU 분해(대칭이면 $LDL^{\mathsf{T}}$ )를 한 번만 해두고 반복마다 전진·후진 대입만 반복한다. 즉 초기 비용은 분해 한 방, 반복 비용은 삼각계 풀이 두 방이라 매우 싸다.

2. 왜 시프트가 가까울수록 빠른가[편집]

$A$ 의 고유값이 $\lambda_1, \dots, \lambda_n$ 이면 $(A - \sigma I)^{-1}$ 의 고유값은 $1/(\lambda_i - \sigma)$ 이고, 고유벡터는 그대로다. 거듭제곱법은 절댓값이 가장 큰 고유값 방향으로 수렴하므로, 결국 $|\lambda_i - \sigma|$ 가 가장 작은 고유벡터가 이긴다.

수렴은 선형이며 한 스텝당 축소율이

\rho = \frac{|\lambda_{j} - \sigma|}{|\lambda_{j'} - \sigma|}

이다( $\lambda_j$ 가 $\sigma$ 에 가장 가까운 고유값, $\lambda_{j'}$ 가 그다음). 시프트를 참 고유값 바로 옆에 붙이면 $\rho$ 가 $10^{-6}$ 같은 값이 되어 반복 한두 번에 기계정밀도에 도달한다. 실무에서 역반복이 “반복법”이라기보다 “한 방짜리 후처리”처럼 쓰이는 이유다.¹

시프트가 없는( $\sigma = 0$ ) 경우는 그냥 최소 절댓값 고유값을 찾는 것이 되어, 거듭제곱법의 정확한 반대편 도구가 된다. 스펙트럼 한복판을 노릴 수 있다는 점이 거듭제곱법과 결정적으로 다르며, 이 발상을 크릴로프 부분공간으로 확장한 것이 이동-역변환(shift-invert) 란초스·아놀디이고 ARPACK의 모드 3이 그것이다.

3. 거의 특이한 계를 일부러 푼다는 역설[편집]

$\sigma \to \lambda_j$ 이면 $A - \sigma I$ 는 특이행렬에 한없이 가까워진다. 조건수가 $10^{14}$ 쯤 되는 계를 매 반복 푸는 셈이고, 수치해석 수업에서 배운 대로라면 결과는 쓰레기여야 한다. 그런데 역반복법은 바로 그 상황에서 가장 잘 작동한다. 이것이 역반복법의 고전적 역설이며, 윌킨슨(J. H. Wilkinson)이 정리한 답은 다음과 같다.

부분 피벗팅 가우스 소거는 후진 안정하므로, 계산된 해 $\hat{\mathbf{y}}$ 는 약간 섭동된 계의 정확한 해다.

(A - \sigma I + E)\,\hat{\mathbf{y}} = \mathbf{x}_k, \qquad \|E\| \approx u\,\|A - \sigma I\|

여기서 $u$ 는 단위 반올림. 그런데 $\sigma$ 가 고유값에 가까우므로 $\|\hat{\mathbf{y}}\|$ 는 어마어마하게 크고, 양변을 $\|\hat{\mathbf{y}}\|$ 로 나누면

\frac{\|(A - \sigma I)\hat{\mathbf{y}}\|}{\|\hat{\mathbf{y}}\|} = \frac{\|\mathbf{x}_k - E\hat{\mathbf{y}}\|}{\|\hat{\mathbf{y}}\|} \;\lesssim\; \frac{1}{\|\hat{\mathbf{y}}\|} + u\|A\|

즉 정규화된 해의 고유잔차는 반올림 수준으로 작다. 오차가 폭발하긴 하는데, 그 폭발이 정확히 “원하는 고유벡터 방향”으로 일어난다. 계가 나쁜 조건이라는 사실은 여기서 병이 아니라 증폭기다.²

이 논점에는 정확히 짚어야 할 단서가 셋 붙는다.

보증되는 것은 작은 잔차이지, 작은 고유벡터 오차가 아니다. 잔차를 고유벡터 오차로 바꾸려면 고유값 간격 $\delta$ 로 나눠야 하고( $\varepsilon \lesssim \|\mathbf{r}\|/\delta$ ), 고유값이 뭉쳐 있으면 이 변환이 무의미해진다. 비정규·결함 행렬에서는 애초에 좋은 시프트로도 좋은 고유벡터를 못 얻는다.
정확도를 “개선”하려 들면 안 된다. 반복 개선(iterative refinement)으로 선형계 해를 더 정확히 만들면, 그건 잔차 $\mathbf{x}_k - E\hat{\mathbf{y}}$ 를 없애는 방향이라 오히려 도움이 안 된다. 선형계의 잔차와 고유쌍의 잔차는 다른 물건이므로, 수렴 판정은 반드시 $\|A\mathbf{x} - \theta\mathbf{x}\|$ 로 한다.
$\sigma$ 가 정확히 고유값이면 $A - \sigma I$ 가 부동소수점에서도 특이가 되어 분해가 0으로 나눌 수 있다. 관행적 처방은 피벗이 0에 가까우면 $u\|A\|$ 수준의 값으로 갈아끼우는 것이다. 정확히 특이한 계를 만나는 것은 사고가 아니라 목표에 도달했다는 신호에 가깝다.

4. 뭉친 고유값 — 진짜 어려운 부분[편집]

역반복법이 실제로 고생하는 지점은 조건수가 아니라 축퇴와 클러스터다. 서로 다른 시작 벡터로 독립적으로 역반복을 돌려 얻은 두 고유벡터는, 대응 고유값이 가까우면 서로 직교하지 않는다. 각자 “그 클러스터가 펼치는 불변 부분공간 안의 아무 방향”으로 수렴해버리기 때문이다.

전통적 처방은 그람-슈미트 재직교화다. LAPACK의 대칭 삼중대각 고유벡터 루틴 xSTEIN이 정확히 이 방식이고, 클러스터가 크면 재직교화 비용이 $\mathcal{O}(n^3)$ 까지 치솟는 데다 결과의 직교성이 완전히 보증되지도 않는다. 이 결함을 정면으로 해결하려고 나온 것이 MRRR(Multiple Relatively Robust Representations, xSTEMR)로, 클러스터마다 시프트를 다르게 잡은 $LDL^{\mathsf{T}}$ 표현을 따로 두어 재직교화 없이 $\mathcal{O}(n^2)$ 에 직교 고유벡터를 얻는다. “역반복법이 잘 안 되는 유일한 경우를 위해 만들어진 알고리즘”이라 봐도 된다.³

5. 레일리 몫 반복으로의 확장[편집]

시프트를 고정하지 않고 매 반복 현재 근사에서 다시 뽑으면 레일리 몫 반복(Rayleigh Quotient Iteration, RQI)이 된다.

\sigma_k = R(\mathbf{x}_k) = \frac{\mathbf{x}_k^{\mathsf{T}}A\mathbf{x}_k}{\mathbf{x}_k^{\mathsf{T}}\mathbf{x}_k}, \qquad (A - \sigma_k I)\,\mathbf{y}_{k+1} = \mathbf{x}_k

대칭 행렬에서 RQI는 세제곱 수렴한다. 근거는 두 조각의 곱이다 — 레일리 몫이 벡터 오차 $\varepsilon$ 에서 고유값 오차 $\mathcal{O}(\varepsilon^2)$ 을 주고, 역반복 한 스텝이 벡터 오차를 시프트 오차에 비례해 줄인다. 합치면 $\varepsilon_{k+1} = \mathcal{O}(\varepsilon_k^3)$ . 대칭 행렬에서는 거의 모든 시작 벡터에 대해 전역 수렴한다는 파를렛의 결과까지 있어서, 이론적으로 매우 강력하다. 비대칭에서는 제곱 수렴으로 떨어진다.

대가는 명확하다. 시프트가 매번 바뀌므로 분해를 재사용할 수 없다. 반복마다 $\mathcal{O}(n^3)$ (희소면 재분해 비용)를 새로 물어야 해서, 실무에서는 고정 시프트 역반복으로 대충 수렴시킨 뒤 마지막 한두 번만 RQI로 마무리하거나, 안쪽 선형계를 반복 솔버로 근사해 푸는 야코비-데이비드슨 계열로 갈아탄다.

6. 구조해석 현장에서[편집]

모드 해석: $K\boldsymbol{\phi} = \omega^2 M\boldsymbol{\phi}$ 에서 관심 주파수 근처에 $\sigma$ 를 놓고 $K - \sigma M$ 을 한 번 분해한 뒤, 벡터 여러 개를 동시에 밀어 넣는 부분공간 반복(subspace iteration, Bathe)을 돌린다. 블록 역반복 + 레일리 몫 사영의 조합이며, 상용 코드에서 블록 란초스와 함께 오랫동안 양대 기본 해석기였다.
모드 누락 검사: $K - \sigma M$ 의 $LDL^{\mathsf{T}}$ 분해에서 음수 피벗의 개수 = $\sigma$ 보다 작은 고유값의 개수다(스투름 수열 성질). 어차피 역반복을 위해 분해를 했으니 이 카운트는 공짜로 나오며, “요구한 20개 모드가 정말 최저 20개인가”를 검증하는 표준 수단이다. 실린더 껍질처럼 모드가 촘촘한 구조에서 이 체크를 빼먹으면 모드를 통째로 건너뛴 것도 모른 채 해석이 끝난다.
좌굴 고유치: $(K + \lambda K_G)\boldsymbol{\phi} = \mathbf{0}$ 에서 예상 임계하중 배수 근처로 시프트를 잡는다. 좌굴 모드는 종종 근접 중근이라, 위에서 말한 클러스터 재직교화가 그대로 문제가 된다.
QR 후처리: 조밀 행렬에서 QR 분해 기반 알고리즘이 고유값을 먼저 다 구한 뒤, 필요한 고유벡터만 그 값을 시프트로 삼아 역반복 한두 번으로 복원한다. LAPACK xHSEIN이 이 일을 한다.

7. 관련 문서[편집]

그래서 역반복 코드에 정교한 수렴 감시 로직을 넣으면 대부분 놀고 있다. 시프트가 좋으면 2회, 나빠도 5회 안에 끝난다. 반복이 20회를 넘어가면 알고리즘을 의심할 게 아니라 시프트나 행렬을 의심해야 한다. ↩
처음 이 논리를 들으면 사기당한 기분이 든다. 실제로 1950~60년대에는 “거의 특이한 계를 푸는 것은 미친 짓”이라는 이유로 역반복법을 기피한 문헌이 있었고, 윌킨슨의 해명 이후에야 표준 도구가 됐다. 나쁜 조건수는 죄가 아니라 어느 방향으로 나쁜지가 문제라는 교훈. ↩
참고로 MRRR은 이론적으로 깔끔하지만 구현이 악명 높게 까다로워서, LAPACK의 dstemr가 특정 입력에서 직교성을 잃는 버그 리포트가 몇 년 주기로 올라온다. 그래서 안전을 중시하는 코드는 여전히 dstein(역반복 + 재직교화)이나 분할정복 dstedc를 쓰기도 한다. 새 알고리즘이 항상 이기는 건 아니다. ↩