거듭제곱법

수치해석 마지막 수정: 2026-07-31 05:14:26

1. 개요[편집]

거듭제곱법 Power Method
별칭	멱승법, power iteration
제안	von Mises · Geiringer, 1929
반복식	$\mathbf{x}_{k+1} = A\mathbf{x}_k / \lVert A\mathbf{x}_k \rVert$
구하는 것	절댓값 최대 고유쌍 딱 1개
수렴	선형, 축소율 $\|\lambda_2/\lambda_1\|$
반복당 비용	행렬-벡터 곱 1회 + 정규화

곱하고, 길이로 나누고, 또 곱한다. 이게 전부다.

거듭제곱법(power method, 멱승법)은 행렬 $A$ 에 벡터를 반복해서 곱하고 매번 정규화하는 것만으로 절댓값이 가장 큰 고유값(지배 고유값)과 그 고유벡터를 얻는 반복 알고리즘이다.¹ 고유값 문제를 푸는 방법 중 가장 단순하며, 란초스 알고리즘·아놀디 알고리즘·역반복법이 전부 이것의 변형이거나 이것을 부품으로 쓴다.

\mathbf{x}_{k+1} = \frac{A\mathbf{x}_k}{\lVert A\mathbf{x}_k \rVert}, \qquad \mu_k = \frac{\mathbf{x}_k^{\mathsf T} A \mathbf{x}_k}{\mathbf{x}_k^{\mathsf T}\mathbf{x}_k}

$A$ 를 명시적으로 저장할 필요조차 없다. 벡터에 곱하는 방법만 있으면 되므로 희소행렬이나 행렬 없는(matrix-free) 연산자에 그대로 붙는다. 고유값 추정치는 레일리 몫 $\mu_k$ 로 뽑는다.

주는 정보가 딱 하나(지배 고유쌍)뿐이고 수렴이 느릴 수 있다는 약점이 분명한데도 이 방법이 90년 넘게 살아남은 이유는 세 가지다. 구현이 자명하고, 필요한 것이 행렬-벡터 곱뿐이며, 지배 고유쌍 하나만 있으면 되는 문제가 세상에 아주 많기 때문이다. 웹 랭킹, 마르코프 연쇄의 정상상태, 스펙트럼 반지름 추정이 전부 그런 문제다.

2. 왜 수렴하는가[편집]

$A$ 가 대각화 가능하고 고유벡터 $\mathbf{v}_i$ 가 기저를 이룬다고 하자. 초기 벡터를 $\mathbf{x}_0 = \sum_i c_i \mathbf{v}_i$ 로 전개하면

A^k \mathbf{x}_0 = \lambda_1^k \left( c_1 \mathbf{v}_1 + \sum_{i \ge 2} c_i \Big(\frac{\lambda_i}{\lambda_1}\Big)^{k} \mathbf{v}_i \right)

이다. $|\lambda_1| > |\lambda_2| \ge \cdots$ 이면 괄호 안의 잔여항이 전부 0으로 죽고 방향은 $\mathbf{v}_1$ 만 남는다. 즉 각오차는 $|\lambda_2/\lambda_1|^k$ 로 줄어드는 선형 수렴이고, 축소율은 오직 스펙트럼 갭이 결정한다. 갭이 넉넉하면 열 번이면 끝나고, 비율이 0.99면 자릿수 하나 얻는 데 230회가 필요하다. 대칭 행렬이면 고유벡터 오차가 $\varepsilon$ 일 때 레일리 몫 오차는 $O(\varepsilon^2)$ 이라, 고유값만 필요하면 체감 속도가 제곱으로 좋아진다.

3. 언제 죽는가[편집]

갭이 없을 때. $|\lambda_1| = |\lambda_2|$ 면 수렴 자체가 없다. 실행렬의 복소 켤레쌍이 지배 고유값이면 반복 벡터가 2차원 부분공간 안에서 영원히 회전한다.
$\lambda_2 = -\lambda_1$ . 짝수·홀수 반복이 서로 다른 벡터로 진동한다. 두 스텝을 묶으면( $A^2$ ) 갭이 생기지만 부호는 따로 복원해야 한다.
$c_1 = 0$ . 이론상 사망이지만 실제로는 반올림 오차가 $\mathbf{v}_1$ 성분을 다시 심어 줘서 결국 수렴한다.²
오버플로/언더플로. 정규화를 빼먹으면 $\lambda_1^k$ 가 그대로 자라거나 0으로 가라앉는다. 정규화는 장식이 아니다.
결함 행렬. 대각화가 안 되는(조던 블록이 있는) 경우에도 방향은 결국 수렴하지만, 속도가 기하급수가 아니라 $1/k$ 꼴로 떨어져 체감상 정지한 것처럼 보인다.

여기서 중요한 것은 이 실패들이 버그가 아니라 알고리즘의 정의역 밖이라는 점이다. 거듭제곱법은 “지배 고유값이 유일하게 존재한다”는 가정을 대놓고 깔고 있고, 그 가정이 깨지면 조용히 틀린 답을 주는 게 아니라 아예 수렴을 멈춘다. 이 점만큼은 정직한 알고리즘이다.

4. 시프트, 디플레이션, 블록화[편집]

축소율이 마음에 안 들면 스펙트럼을 옮기면 된다. $A - \sigma I$ 의 고유값은 $\lambda_i - \sigma$ 이므로 $\sigma$ 를 잘 고르면 비율 $|\lambda_2 - \sigma| / |\lambda_1 - \sigma|$ 를 줄일 수 있다. 이 발상을 역행렬 쪽으로 돌린 것이 $(A-\sigma I)^{-1}$ 에 거듭제곱법을 거는 역반복법이고, $\sigma$ 를 매 반복 레일리 몫으로 갱신하면 대칭계에서 세제곱 수렴하는 레일리 몫 반복이 된다. 그쪽 이야기는 해당 문서에 위임한다.

두 번째 고유쌍이 필요하면 디플레이션을 쓴다. 대칭 행렬이면 호텔링 수축 $A' = A - \lambda_1 \mathbf{v}_1\mathbf{v}_1^{\mathsf T}$ 로 $\lambda_1$ 을 0으로 눌러 버리면 된다.³ 비대칭이면 좌고유벡터가 필요하고, 수축을 거듭할수록 오차가 누적돼 서너 개가 실용 한계다.

애초에 여러 개를 한꺼번에 원한다면 벡터 하나 대신 $n \times p$ 블록을 곱하고 매 스텝 QR 분해로 직교화하는 블록 거듭제곱법, 곧 부분공간(직교) 반복이 정석이다. $j$ 번째 벡터의 수렴률은 $|\lambda_{p+1}/\lambda_j|$ 가 되어 블록을 키울수록 빨라진다. 여기에 전처리기를 물린 현대적 후예가 LOBPCG다. 여기서 한 걸음만 더 나가 지금까지 만든 벡터를 전부 기저로 재활용하면 크릴로프 계열(크리로프 부분공간법)이 되고, 실무 구현이 ARPACK이다.

5. 정지 조건과 실전 감각[편집]

반복을 언제 끊을지는 벡터 변화량이 아니라 잔차 노름 $r_k = \lVert A\mathbf{x}_k - \mu_k \mathbf{x}_k \rVert$ 로 판정하는 것이 정석이다. 대칭 행렬이면 $r_k$ 가 그대로 후진 오차 해석에서 말하는 섭동 크기이고, 실제 고유값과의 거리는 $|\mu_k - \lambda_1| \le r_k$ 로 즉시 묶인다. 비대칭이면 이 보장이 깨져서 고유값 조건수만큼 부풀 수 있으니 잔차만 보고 안심하면 안 된다.

세 가지를 기억하면 실무에서 판단이 빨라진다.

상황	증상	대처
갭이 좁다	잔차가 완만한 직선으로 기어감	시프트, 또는 크릴로프 계열로 갈아탐
복소 켤레쌍 지배	잔차가 주기적으로 요동	2차원 블록 반복
고유값만 필요	잔차는 크나 $\mu_k$ 는 이미 안정	조기 종료

반복당 비용이 행렬-벡터 곱 하나라는 점은 이 알고리즘의 최대 장점이자 존재 이유다. 코드 열 줄이면 되고, 메모리는 벡터 두 개면 되고, 병렬 컴퓨팅 환경에서는 곱셈 커널 하나만 병렬화하면 그대로 스케일한다. 정확도 요구가 낮은 스펙트럼 반지름 추정(예: 명시적 시간적분의 안정 시간간격 산정)에서는 아직도 이게 최선의 선택이다.

6. 페이지랭크와 정상분포[편집]

구글의 초기 페이지랭크는 사실상 거듭제곱법 그 자체다. 링크 구조를 확률행렬 $S$ 로 만들고 감쇠계수 $\alpha$ 로 균등 점프를 섞은 구글 행렬

G = \alpha S + (1-\alpha)\tfrac{1}{n}\mathbf{e}\mathbf{e}^{\mathsf T}

에 반복을 돌리면 $\lambda_1 = 1$ 의 고유벡터, 즉 정상분포가 나온다. 핵심은 $|\lambda_2(G)| \le \alpha$ 가 증명되어 있다는 것 — $\alpha = 0.85$ 면 축소율이 웹 그래프의 생김새와 무관하게 0.85 이하로 보장되고, 그래서 수십억 노드에서도 100회 남짓이면 8자리가 맞는다.⁴ 감쇠계수를 1에 붙일수록 원래 링크 구조에는 충실해지지만 수렴은 끝없이 느려진다. 즉 여기서는 감쇠계수라는 모델링 파라미터가 곧 수치적 축소율이라, 모델과 알고리즘이 한 몸으로 묶여 있는 보기 드문 구조다.

같은 논리로 마르코프 연쇄의 정상분포(마르코프 연쇄 몬테카를로에서 목표분포 수렴을 논할 때의 그 스펙트럼 갭)도 전이행렬에 거듭제곱법을 건 것이며, $|\lambda_2|$ 가 곧 혼합 시간을 지배한다. 구조해석에서는 방향이 반대다 — 필요한 건 최저차 모드 해석 고유쌍인데 거듭제곱법은 최대 쪽으로 가므로, 이동-역변환을 씌워 원하는 대역을 지배 고유값으로 만들어 준 뒤에야 쓸 수 있다.

7. 관련 문서[편집]

리하르트 폰 미제스와 힐다 가이링거가 1929년 ZAMM에 정리해 실은 것이 표준 출처로 꼽힌다. 통계에서 폰 미제스 분포로 유명한 그 폰 미제스가 아니라 그의 형제 쪽이냐고 묻는 사람이 많은데, 폰 미제스 응력의 리하르트 본인 맞다. 응력 항복조건과 고유값 반복을 둘 다 남긴 사람이다. ↩
반올림 오차가 알고리즘을 살려 주는 흔치 않은 사례. 물론 “그러니까 초기 벡터는 아무거나 넣어도 된다”가 아니라 “난수 벡터를 넣어라”가 실무 결론이다. 특별한 구조를 가진 초기 벡터(예: 전부 1)가 하필 대칭성 때문에 $c_1 = 0$ 인 부분공간에 정확히 놓이는 일은 생각보다 자주 일어난다. ↩
이름은 심리학자 해럴드 호텔링에서 왔다. 주성분 분석을 정리한 그 사람이 맞고, 실제로 초기 PCA 계산 절차가 “거듭제곱법으로 1주성분 뽑고 수축하고 반복”이었다. ↩
Haveliwala & Kamvar(2003)의 결과다. 브린과 페이지가 0.85를 쓴 이유를 두고 온갖 해석이 붙지만, 결과적으로 이 값은 “수렴 속도”와 “원 그래프 충실도”의 타협점이고 그 트레이드오프가 $|\lambda_2| \le \alpha$ 한 줄에 전부 들어 있다. ↩