하우스홀더 변환

1. 개요[편집]

하우스홀더 변환 Householder transformation
정의	$H = I - 2\,vv^{\top}/(v^{\top}v)$
성질	대칭 · 직교 · 대합($H^2=I$) · $\det H = -1$
기하	$v$ 에 수직인 초평면에 대한 반사
QR 비용	$2n^2(m - n/3)$ flops ($m\times n$, $m\ge n$)
안정성	후진 안정 (Wilkinson)
LAPACK	LAPACK dgeqrf / dlarf / dlarfb

회전은 한 번에 성분 하나를 지운다. 반사는 한 번에 열 하나를 통째로 지운다.

하우스홀더 변환(Householder transformation, 하우스홀더 반사)은 0이 아닌 벡터 $v$ 에 대해

$$H = I - \frac{2\,v v^{\top}}{v^{\top} v}$$

로 정의되는 행렬로, 기하적으로는 $v$ 를 법선으로 하는 초평면에 대한 거울 반사다. 임의의 $x$ 를 $v$ 방향 성분과 초평면 성분으로 쪼개면 $Hx$ 는 앞쪽만 부호를 뒤집는다. 이름은 1958년 이 반사를 써서 행렬을 삼각화하는 절차를 발표한 앨스턴 하우스홀더에게서 왔다.¹

성질이 지나칠 정도로 좋다. $H^{\top} = H$ 이고 $H^{\top}H = H^2 = I$ 이므로 대칭이면서 직교이고 자기 자신이 역행렬이다. 고유값은 $-1$ 이 한 번, $+1$ 이 $n-1$ 번이라 $\det H = -1$ — 반사이니 방향을 뒤집는 것이 당연하다. 직교행렬이라 2-노름을 정확히 보존하고 조건수가 1이며, 그래서 반올림 오차를 증폭시키지 않는다.

2. 열 하나를 통째로 0으로[편집]

실전 가치는 하나의 사실에서 나온다. 임의의 벡터 $x$ 에 대해 $v$ 를 적절히 고르면 $Hx$ 를 $e_1$ 축 위의 한 점으로 한 방에 보낼 수 있다. 노름이 보존되므로 목적지는 $\mp\|x\|_2 e_1$ 뿐이고,

$$v = x + \operatorname{sign}(x_1)\,\|x\|_2\, e_1, \qquad Hx = -\operatorname{sign}(x_1)\,\|x\|_2\, e_1$$

로 잡으면 된다. 첫 성분 아래가 전부 0이 되므로, QR 분해에서 한 열의 대각 아래 전체를 반사 한 번으로 청소할 수 있다. 기븐스 회전이 성분을 하나씩 지우는 것과 대비되는 지점이다.

2.1. 부호를 왜 $\operatorname{sign}(x_1)$ 로 고르는가[편집]

수학적으로는 $v = x \pm \|x\|e_1$ 둘 다 옳다. 수치적으로는 아니다. $x_1 > 0$ 인데 $v = x - \|x\|e_1$ 을 쓰면 $v_1 = x_1 - \|x\|$ 인데, $x$ 가 이미 $e_1$ 에 거의 나란한 경우 $x_1 \approx \|x\|$ 라서 거의 같은 두 수의 뺄셈이 된다. 유효숫자가 통째로 날아가는 전형적인 상쇄(cancellation)이고, 그렇게 오염된 $v$ 로 만든 $H$ 는 직교성을 잃는다. $x_1$ 과 같은 부호를 더하면 $|v_1| = |x_1| + \|x\|$ 라 상쇄가 원천적으로 일어나지 않는다.

부득이 반대 부호가 필요하면(예컨대 $R$ 의 대각을 양수로 강제하고 싶을 때) 분자를 유리화해

$$v_1 = \frac{x_1^2 - \|x\|_2^2}{x_1 + \|x\|_2} = \frac{-(x_2^2 + \cdots + x_m^2)}{x_1 + \|x\|_2}$$

로 계산하면 상쇄를 피할 수 있다. 부동소수점 연산에서 “대수적으로 같은 두 식이 수치적으로 전혀 다르다”는 원칙의 교과서적 사례다.

또 하나의 실무 규약: $H$ 를 절대 명시적으로 만들지 않는다. $HB = B - \frac{2}{v^{\top}v} v(v^{\top}B)$ 이므로 행렬-벡터 곱 두 번이면 끝이고, $m\times k$ 행렬에 적용하는 데 $4mk$ flops다. $H$ 를 만들어 곱하면 $2m^2k$ — 차원 하나만큼 손해다. LAPACK은 $v_1 = 1$ 로 정규화한 $v$ 의 나머지 성분을 소거된 자리(하삼각)에 그대로 저장하고, 스칼라 $\tau = 2/(v^{\top}v)$ 만 따로 배열에 담아 $H = I - \tau vv^{\top}$ 로 쓴다. 추가 메모리가 사실상 0이다.

3. QR 분해와 축약[편집]

$m\times n$ ($m \ge n$) 행렬에 반사를 $n$ 번 순차 적용하면 $H_n\cdots H_1 A = R$, 즉 $A = QR$ 이다. 비용은

$$2n^2\!\left(m - \frac{n}{3}\right) \ \text{flops}$$

이고 정사각행렬이면 $\tfrac{4}{3}n^3$ 이다. 여기서 $Q$ 는 반사 $n$ 개의 곱으로만 존재하며, 필요한 경우에만 축약 $Q$($m\times n$, 같은 $2n^2(m-n/3)$)나 완전 $Q$($m\times m$, $4(m^2n - mn^2 + n^3/3)$)를 조립한다. 최소자승법에서는 $Q^{\top}b$ 만 있으면 되므로 $Q$ 를 만들 일 자체가 없다.

같은 도구가 축약(reduction)에도 그대로 쓰인다. 양쪽에서 상사변환으로 치되 첫 부대각은 남겨 두는 방식이다.

• 상헤센베르크 축약 — 비대칭 $A$ 를 $Q^{\top}AQ = H$ 로. $\tfrac{10}{3}n^3$($Q$ 누적 시 $\tfrac{14}{3}n^3$). 슈어 분해와 QR 알고리즘의 필수 전처리다.
• 삼중대각화 — 대칭 $A$ 를 삼중대각으로. $\tfrac{4}{3}n^3$. 대칭 고유값 solver의 입구.
• 쌍대각화 — 좌우에서 서로 다른 반사를 걸어 상쌍대각으로(Golub–Kahan). $4mn^2 - \tfrac{4}{3}n^3$. 특이값 분해의 1단계.

복소행렬에서는 반사가 유니터리여야 하므로 $H = I - 2vv^{H}/(v^{H}v)$ 로 쓰고, 부호 대신 $x_1$ 의 위상 $x_1/|x_1|$ 을 곱해 같은 상쇄 회피를 구현한다. 이때 $H$ 는 에르미트·유니터리를 유지하지만 대각 성분이 복소수가 되므로, ”$R$ 의 대각을 실수 양수로”라는 관례를 지키려면 마지막에 대각 위상 행렬을 한 번 더 곱해줘야 한다.

왜 전부 완전 삼각화가 아니라 “한 칸 남기는” 형태에서 멈추는가? 상사변환은 양쪽에서 곱해야 고유값을 보존하는데, 왼쪽에서 지운 성분을 오른쪽 곱이 되살려 놓기 때문이다. 유한 단계로 도달 가능한 한계가 정확히 헤센베르크/삼중대각이고, 나머지는 반복법의 몫이다.

4. 안정성 — 그람-슈미트와의 대비[편집]

하우스홀더 QR은 후진 안정하다. 계산된 $\hat R$ 에 대해 $A + \Delta A = \tilde{Q}\hat{R}$ 를 만족하는 정확히 직교인 $\tilde Q$ 가 존재하고 $\|\Delta A\|_F \le c(m,n)\,u\,\|A\|_F$ 다. 곱해지는 것이 전부 노름 보존 연산이라 오차가 커질 자리가 없다.²

그람-슈미트 직교화는 사정이 다르다. 계산된 $\hat Q$ 의 직교성 손실이

$$\|I - \hat{Q}^{\top}\hat{Q}\| \approx \begin{cases} u\,\kappa(A)^2 & \text{고전 그람-슈미트} \\ u\,\kappa(A) & \text{수정 그람-슈미트} \end{cases}$$

로 조건수에 비례해 커진다. $\kappa(A) = 10^8$ 이면 고전 GS의 $\hat Q$ 는 배정도에서 직교성을 완전히 잃는다. 하우스홀더는 이 항이 $O(u)$ — $\kappa$ 가 들어오지 않는다.³ 대신 GS는 열을 하나씩 만들어내므로 아놀디 알고리즘처럼 부분공간을 점증적으로 키우는 상황에 맞고, 그래서 두 방법은 대체재가 아니라 용도가 갈린다.

기븐스 회전과의 선택 기준은 더 단순하다. 밀집 행렬 전체를 삼각화할 때 기븐스는 하우스홀더보다 약 50% 비싸다. 그러나 이미 대부분 0인 희소·헤센베르크·띠 구조에서는 실제로 죽여야 할 성분만 골라 칠 수 있고, 행 하나를 덧붙인 QR 갱신처럼 국소 수정에도 유리하며, 서로 겹치지 않는 회전은 병렬로 돌아간다. GMRES 내부에서 헤센베르크 최소자승을 푸는 것이 기븐스인 이유다.

5. 블록화 — WY 표현[편집]

flops만 보면 하우스홀더 QR은 이상적이지만, 반사를 하나씩 적용하는 코드는 전부 행렬-벡터 곱(BLAS-2)이라 현대 CPU의 캐시와 벡터 유닛을 놀린다. 해법은 반사 $b$ 개를 묶어 한 덩어리로 표현하는 것이다.

$$H_1 H_2 \cdots H_b = I - V T V^{\top}$$

여기서 $V$ 는 $v_i$ 들을 열로 모은 $m \times b$ 행렬, $T$ 는 $b\times b$ 상삼각행렬이다(compact WY 표현, Schreiber–Van Loan 1989). 이 형태면 적용이 행렬-행렬 곱(BLAS-3)이 되어 캐시 재사용이 살아난다. LAPACK의 dgeqrf 가 정확히 이 구조로, dgeqr2 로 패널을 처리하고 dlarft 로 $T$ 를 만든 뒤 dlarfb 로 나머지 블록에 한꺼번에 반사를 먹인다. 연산량은 오히려 조금 늘지만 실측 성능은 몇 배 차이가 난다.

6. 관련 문서[편집]

• QR 분해 · 기븐스 회전 · 그람-슈미트 · 최소자승법
• 슈어 분해 · 고유값 문제 · 특이값 분해 · 아놀디 알고리즘
• 조건수 · 후진 오차 해석 · 부동소수점 연산
• 회전행렬 · LAPACK

7. Footnotes[편집]

1. Householder, A. S. (1958). “Unitary Triangularization of a Nonsymmetric Matrix.” JACM 5(4). 다섯 쪽짜리 논문 하나가 이후 60년간 밀집 선형대수 라이브러리의 뼈대가 됐다. 참고로 하우스홀더는 원래 수리생물학자로 커리어를 시작했다. ↩

2. “그럼 $\kappa$ 가 $10^{16}$ 이어도 하우스홀더면 안전한가?” 후진 안정은 “이웃 문제의 정확한 답을 준다”는 약속이지 “정확한 답을 준다”는 약속이 아니다. 전진 오차는 여전히 $\kappa(A)u$ 로 커진다. 알고리즘은 무죄, 행렬이 유죄. ↩

3. 그렇다고 수정 그람-슈미트가 쓰레기라는 뜻은 아니다. MGS로 푼 최소자승 해 자체는 하우스홀더와 같은 수준으로 후진 안정하다 — 망가지는 것은 $\hat Q$ 의 직교성이지 $x$ 가 아니다. 문제는 그 $\hat Q$ 를 직교라고 믿고 다음 계산에 넘길 때 터진다. ↩