고유값 조건수

1. 개요

고유값 조건수 Eigenvalue Condition Number
정의	$\kappa(\lambda) = \lVert\mathbf x\rVert\lVert\mathbf y\rVert/\lvert\mathbf y^*\mathbf x\rvert$
뜻	$\lvert\delta\lambda\rvert \lesssim \kappa(\lambda)\lVert E\rVert$ (1차)
기하학	$1/\cos\angle(\mathbf x, \mathbf y)$
최솟값	1 (정규 행렬은 모든 고유값이 1)
무한대	결함 고유값(조르당 블록)
계산	LAPACK xTRSNA·xGEEVX, MATLAB condeig
제안	윌킨슨 (1965)

행렬 전체의 조건수가 멀쩡해도, 고유값 하나하나는 제각각 예민할 수 있다.

고유값 조건수(eigenvalue condition number) $\kappa(\lambda)$ 는 행렬 $A$ 를 $E$ 만큼 흔들었을 때 단순 고유값 $\lambda$ 가 얼마나 움직이는지를 재는 수로,

$$\kappa(\lambda) = \frac{\lVert\mathbf x\rVert_2\,\lVert\mathbf y\rVert_2}{\lvert\mathbf y^*\mathbf x\rvert}$$

로 정의된다. 여기서 $\mathbf x$ 는 오른쪽 고유벡터($A\mathbf x = \lambda\mathbf x$), $\mathbf y$ 는 왼쪽 고유벡터($\mathbf y^*A = \lambda\mathbf y^*$)다. 윌킨슨이 1965년 『대수적 고유값 문제』에서 체계화했다.

선형계의 조건수 $\kappa(A) = \lVert A\rVert\lVert A^{-1}\rVert$ 는 행렬 하나에 숫자 하나지만, 고유값 문제는 고유값마다 따로 조건수가 있다. 바우어-파이크 정리가 모든 고유값을 고유벡터 행렬의 조건수 $\kappa(V)$ 하나로 묶어 버리는 데 비해, 고유값 조건수는 「어느 고유값을 믿어도 되는가」 를 개별로 알려 준다. 수치 고유값 해법은 후진 안정이므로, 계산된 고유값의 실제 오차는 대략 (기계 엡실론) × $\lVert A\rVert$ × $\kappa(\lambda)$ 다. 이 문서는 그 공식을 유도하고, 윌킨슨의 악명 높은 예제로 실측한다.

2. 유도 — 1차 섭동

단순 고유값 $\lambda$ 와 고유벡터 $\mathbf x$ 가 $A + \varepsilon E$ 에서 $\lambda(\varepsilon)$, $\mathbf x(\varepsilon)$ 로 매끄럽게 변한다고 하자. $(A + \varepsilon E)\mathbf x(\varepsilon) = \lambda(\varepsilon)\mathbf x(\varepsilon)$ 을 $\varepsilon$ 으로 미분하고 $\varepsilon = 0$ 에서

$$A\mathbf x' + E\mathbf x = \lambda'\mathbf x + \lambda\mathbf x'$$

이다. 왼쪽에서 $\mathbf y^*$ 를 곱하면 $\mathbf y^*A = \lambda\mathbf y^*$ 이므로 $\mathbf x'$ 항이 지워지고

$$\lambda'(0) = \frac{\mathbf y^*E\,\mathbf x}{\mathbf y^*\mathbf x}$$

가 남는다. 코시-슈바르츠로 $\lvert\mathbf y^*E\mathbf x\rvert \le \lVert\mathbf y\rVert\lVert E\rVert_2\lVert\mathbf x\rVert$ 이므로

$$\lvert\delta\lambda\rvert \le \kappa(\lambda)\,\lVert\delta A\rVert_2 + O(\lVert\delta A\rVert^2)$$

이고, 등호는 $E = \mathbf y\mathbf x^*/(\lVert\mathbf y\rVert\lVert\mathbf x\rVert)$ 에서 달성된다. 즉 $\kappa(\lambda)$ 는 1차 섭동에서 정확한(최악의) 증폭률이다. 같은 계산에서 $\lambda$ 가 단순하다는 가정이 $\mathbf y^*\mathbf x \ne 0$ 을 보장한다.

기하학적 뜻. $\kappa(\lambda) = 1/\cos\theta$, $\theta$ 는 왼쪽·오른쪽 고유벡터 사이 각이다. 정규 행렬(대칭·에르미트·유니터리 포함, 반대 개념은 비정규 행렬)은 $\mathbf y = \mathbf x$ 라서 모든 고유값이 $\kappa = 1$ 이고, 고유값이 $\lVert E\rVert_2$ 이상 움직이지 못한다(에르미트의 경우 바일 부등식). 비정규성이 커져 두 벡터가 거의 직교할수록 조건수가 커진다.

결함 고유값. $\lambda$ 가 크기 $m$ 의 조르당 표준형 블록에 속하면 $\mathbf y^*\mathbf x = 0$ 이라 $\kappa = \infty$ 다. 이때 섭동은 1차가 아니라 $\lvert\delta\lambda\rvert \sim \lVert E\rVert^{1/m}$ 로 움직인다(횔더 연속). 기계 엡실론 $10^{-16}$ 이 $m = 2$ 에서는 $10^{-8}$, $m = 16$ 에서는 $0.1$ 이 된다.

간단한 예. $A = \begin{pmatrix}1 & t\\ 0 & 2\end{pmatrix}$ 이면 두 고유값 모두 $\kappa = \sqrt{1+t^2}$ 이다. $t = 10^4$ 일 때 $(2,1)$ 성분에 $10^{-8}$ 을 더하면 고유값이 각각 $\mp 0.99990\times10^{-4}$ 만큼, 즉 섭동의 약 1만 배 움직인다. 1차 예측 $\kappa\cdot\varepsilon = 10^{-4}$ 과 맞는다.

3. 윌킨슨의 행렬 — 무엇이 실제로 깨지나

윌킨슨이 즐겨 쓴 20×20 상삼각 이중대각 행렬

$$W = \begin{pmatrix} 20 & 20 & & & \\ & 19 & 20 & & \\ & & \ddots & \ddots & \\ & & & 2 & 20\\ & & & & 1 \end{pmatrix}$$

은 고유값이 대각성분 $1, 2, \dots, 20$ 이고, 서로 1 씩 떨어진 멀쩡한 정수다. $\lVert W\rVert_2 = 36.8$ 이다. 이 행렬의 고유값 조건수를 계산하면

$\lambda$	1, 20	5	10, 11	15
$\kappa(\lambda)$	$8.4\times10^{7}$	$2.4\times10^{11}$	$5.1\times10^{12}$	$6.9\times10^{11}$

가운데 고유값이 가장 예민하다. 실제로 확인해 보자.

한 성분만 흔들기. 왼쪽 아래 구석 $(20, 1)$ 성분에 $\varepsilon$ 을 넣으면 특성다항식이 $\varepsilon\cdot 20^{19} \approx 5.2\times10^{24}\,\varepsilon$ 만큼 바뀐다. $\prod(\lambda - i)$ 의 10~11 사이 극값이 $4.09\times10^{11}$ 이므로 $\varepsilon \approx 7.8\times10^{-14}$ 에서 고유값 10 과 11 이 충돌해 복소쌍이 된다. 행렬 원소의 상대 크기로 따지면 $10^{-15}$ 급의 섭동, 즉 배정밀도 반올림 몇 번에 해당한다. $\varepsilon = 10^{-10}$ 이면 20개 중 14개가 복소수가 되고 허수부가 최대 2.73 까지 커진다. 행렬 원소는 한 자리도 안 바뀌었는데 고유값은 이미 정수 근처에도 없다.

후진 안정 해법의 정직한 오차. 직교 닮음 변환 $QWQ^T$ 로 섞은 뒤(삼각 구조를 숨겨야 해법이 공짜로 정답을 내지 못한다) 표준 QR 고유값 해법을 돌리면, 무작위 $Q$ 세 개에서

• 고유값 1, 20 의 오차: $4\times10^{-8} \sim 2\times10^{-7}$ (예측 $u\lVert W\rVert\kappa \approx 3\times10^{-7}$)
• 가장 나쁜 고유값의 오차: $0.007 \sim 0.013$ (예측 $\approx 0.02$)

로, 조건수가 예측한 자릿수가 정확히 빠진다. 해법이 잘못한 것이 아니다. 해법은 반올림 수준의 후진 오차만 냈고(후진 오차 해석), 문제 자체가 그 후진 오차를 $10^{12}$ 배로 키운 것이다¹.

가까운 다중 고유값. 윌킨슨(1972)은 $\kappa(\lambda) > 1$ 이면 $\lVert E\rVert_2 \le \lVert A\rVert_2/\sqrt{\kappa(\lambda)^2 - 1}$ 인 섭동으로 $\lambda$ 를 다중 고유값으로 만들 수 있음을 보였다. 위 행렬에 대입하면 $7.3\times10^{-12}$ 이하의 섭동이 존재한다는 뜻이고, 실제로 성분 하나에 $7.8\times10^{-14}$ 만 넣어도 충돌했다(정리는 상한만 준다). 직관적으로, 조건수가 크다 = 결함 행렬이 바로 옆에 있다 는 뜻이다. 조건수가 큰 고유값 근처의 의사스펙트럼이 반지름 대략 $\kappa(\lambda)\varepsilon$ 의 원판으로 부풀어 이웃과 합쳐지는 것이 같은 현상의 그림이다.

4. 고유벡터는 또 다르다

고유값이 멀쩡해도 고유벡터는 예민할 수 있다. 1차 섭동에서 고유벡터의 변화는 다른 고유값과의 간격에 반비례한다. 에르미트 행렬에서는 데이비스-카한 정리가 $\sin\angle(\mathbf x, \tilde{\mathbf x}) \lesssim \lVert E\rVert/\mathrm{gap}$ 을 준다. 그래서 대칭 행렬이라도 고유값이 거의 겹치면 고유값은 $\lVert E\rVert$ 만큼만 움직이는데 고유벡터는 제멋대로 돈다. 이때 개별 고유벡터 대신 군집이 펼치는 불변 부분공간을 계산하는 것이 표준 처방이다. 비대칭 행렬에서는 간격 대신 블록 사이의 분리량 $\mathrm{sep}(\lambda, T_{22})$ 가 그 역할을 하고, 이것은 간격보다 훨씬 작을 수 있다.

5. 실무에서 계산하기

• LAPACK. xGEEVX 가 고유값과 함께 고유값 조건수의 역수 RCONDE 와 고유벡터 쪽 RCONDV 를 돌려준다. 내부적으로 슈어 분해에서 xTRSNA 를 부른다. 역수를 주는 이유는 결함 고유값에서 0 이 되어 오버플로가 없기 때문이다².
• MATLAB·SciPy. MATLAB condeig(A). SciPy 에는 전용 함수가 없지만 scipy.linalg.eig(A, left=True, right=True) 로 양쪽 고유벡터를 받아 정의대로 계산하면 된다(위 표가 그렇게 나왔다).
• 오차 막대. 계산된 고유값 $\hat\lambda$ 에 $\epsilon_{\text{mach}}\lVert A\rVert/\mathrm{rcond}$ 를 오차 막대로 붙이는 것이 LAPACK 사용자 지침의 권장이다. 고유값을 설계 판정(안정성·공진 회피)에 쓴다면 이 막대가 판정 경계를 넘는지부터 본다.

시뮬레이션에서 조건수가 큰 고유값이 나오는 대표적인 곳은 비정규 연산자다. 대류가 지배하는 이산화 행렬, 유동 안정성 해석의 오어-좀머펠트 연산자, 감쇠가 큰 구조 모드가 그렇다. 이런 문제에서 「고유값의 실수부가 음수니 안정」 이라는 판정은 조건수를 같이 보지 않으면 믿을 수 없고, 그래서 의사스펙트럼과 과도 성장 해석이 따로 발달했다³.

6. 관련 문서

• 바우어-파이크 정리 · 조건수 · 후진 오차 해석 · 거슈고린 원판 정리
• 고유값 문제 · QR 알고리즘 · 슈어 분해 · LAPACK
• 비정규 행렬 · 의사스펙트럼 · 조르당 표준형
• 바일 부등식 · 섭동 이론 · 윌킨슨 다항식

7. Footnotes

1. 윌킨슨은 다항식 쪽에서도 $\prod_{i=1}^{20}(x-i)$ 의 $x^{19}$ 계수를 $2^{-23}$ 만큼 바꾸면 근 절반이 복소평면으로 날아가는 「배신의 다항식」 을 보였다. 이쪽은 고유값이 아니라 근이지만, 동반 행렬로 바꾸면 같은 이야기다. 정수 1~20 은 그에게 거의 저주받은 숫자였던 셈. ↩

2. LAPACK 이 조건수 대신 그 역수를 주는 설계는 「무한대를 숫자로 표현하지 말고 0 으로 표현하라」 는 수치해석 국룰의 교과서 사례다. ↩

3. 고유값만 보고 「안정」 도장을 찍었는데 시간적분을 돌려 보니 에너지가 $10^{5}$ 배 뛰었다가 줄어든다. 고유값은 거짓말을 하지 않았다, 다만 다 말하지 않았을 뿐. ↩