뉴턴-크리로프법

1. 개요

뉴턴-크리로프법 Newton–Krylov / JFNK
바깥 반복	(비정확) 뉴턴법
안쪽 반복	크리로프 부분공간법 (주로 GMRES)
야코비안	저장 안 함 — 함수 평가 차분으로 Jv 만 근사
핵심 조율	강제항 ηk (아이젠슈타트–워커), 섭동 크기 ε, 전처리기
대표 구현	PETSc SNES, SUNDIALS KINSOL, Trilinos NOX, MOOSE(PJFNK)

야코비안은 필요 없다. 야코비안이 벡터에 하는 짓만 알면 된다.

뉴턴-크리로프법(Newton–Krylov method)은 비선형 연립방정식 $\mathbf F(\mathbf u) = \mathbf 0$ 을 뉴턴-랩슨법으로 풀되, 매 뉴턴 스텝의 선형계 $J\,\delta\mathbf u = -\mathbf F$ 를 GMRES 같은 크리로프 부분공간법으로 근사적으로 푸는 이중 반복 기법이다. 크리로프 솔버가 행렬 자체가 아니라 행렬-벡터 곱만 요구한다는 점을 이용해 $J\mathbf v$ 를 함수 평가의 차분으로 대신하면 야코비안-프리 뉴턴-크리로프법(Jacobian-free Newton–Krylov, JFNK)이 된다. 크놀과 키스의 2004년 서베이가 이 분야의 표준 참고문헌이다.

이 아이디어 자체는 GMRES·자코비안 행렬·PETSc 문서에서 이미 한 문단씩 소개됐다. 이 문서는 그 한 문단 뒤에 숨은 세 개의 손잡이 — 안쪽 반복을 언제 멈추나(강제항), 차분 간격을 얼마로 잡나(섭동 크기), 크리로프 반복을 어떻게 줄이나(전처리) — 를 숫자와 함께 판다. 실무에서 JFNK가 “되느냐 안 되느냐”는 거의 이 셋에서 갈린다.

2. 비정확 뉴턴 — 안쪽을 대충 푸는 게 정답이다

뉴턴 스텝을 정확히 풀 필요는 없다. 뎀보·아이젠슈타트·스테이하우그(1982)의 비정확 뉴턴법(inexact Newton)은 스텝 $\mathbf s_k$ 가 다음만 만족하면 된다고 본다.¹

$$\|\mathbf F(\mathbf u_k) + J(\mathbf u_k)\,\mathbf s_k\| \le \eta_k\,\|\mathbf F(\mathbf u_k)\|, \qquad 0 \le \eta_k < 1.$$

$\eta_k$ 를 강제항(forcing term)이라 한다. 이론은 깔끔하다.

• $\eta_k \le \eta_{\max} < 1$ 로 고정하면 국소 선형 수렴.
• $\eta_k \to 0$ 이면 초선형 수렴.
• $\eta_k = O(\|\mathbf F(\mathbf u_k)\|)$ 이면 2차 수렴 — 정확한 뉴턴과 같다.

여기서 실무의 함정이 하나 보인다. 해에서 멀 때는 선형화 자체가 부정확하므로, 그 선형계를 기계정밀도까지 푸는 건 순수한 낭비다. 이걸 과잉 풀이(oversolving)라 한다. 아이젠슈타트와 워커(1996)는 바깥 잔차의 줄어드는 속도를 보고 $\eta_k$ 를 자동으로 정하는 규칙을 제안했다.² 그중 널리 쓰이는 두 번째 선택은

$$\eta_k = \gamma \left( \frac{\|\mathbf F(\mathbf u_k)\|}{\|\mathbf F(\mathbf u_{k-1})\|} \right)^{\alpha}, \qquad \gamma = 0.9,\ \alpha = 2,$$

에 $\eta_{k-1}$ 이 너무 급하게 작아지지 않게 하는 안전장치($\gamma\eta_{k-1}^{\alpha} > 0.1$ 이면 그 값과의 최댓값을 취함)를 붙인 것이다. 바깥이 빠르게 수렴할 때만 안쪽을 빡빡하게 푼다.

3. 섭동 크기 ε — 둥근 오차와 절단 오차의 줄다리기

JFNK의 핵심 근사는 방향 도함수의 1차 차분이다.

$$J(\mathbf u)\,\mathbf v \approx \frac{\mathbf F(\mathbf u + \varepsilon\mathbf v) - \mathbf F(\mathbf u)}{\varepsilon}.$$

$\varepsilon$ 이 크면 절단 오차 $O(\varepsilon\|\mathbf v\|^2)$ 가, 작으면 부동소수점 연산의 소거 오차 $O(\epsilon_{\text{mach}}/\varepsilon)$ 가 커진다. 둘을 맞추면 $\varepsilon \sim \sqrt{\epsilon_{\text{mach}}} \approx 1.5\times10^{-8}$ 근처가 나오고, 실무에서는 $\mathbf u$ 와 $\mathbf v$ 의 크기로 정규화해

$$\varepsilon = \frac{\sqrt{\epsilon_{\text{mach}}}\,(1 + \|\mathbf u\|)}{\|\mathbf v\|}$$

같은 형태를 쓴다(크놀–키스 서베이에 여러 변형이 정리돼 있다). 이건 수치미분의 최적 간격 문제와 판박이다. 차분 대신 자동 미분의 전진 모드로 $J\mathbf v$ 를 정확히 계산하면 이 손잡이가 아예 사라진다.

미분 불가능성 주의. $\mathbf F$ 안에 max, if, 플럭스 제한자, 상변화 스위치처럼 매끄럽지 않은 부분이 있으면 $\varepsilon$ 섭동이 분기를 넘나들며 $J\mathbf v$ 가 엉망이 된다. 충격 포착 CFD에서 JFNK가 이유 없이 정체되면 제일 먼저 의심할 곳이다.

4. 전처리 — 행렬 없이 행렬을 쓰는 법

전처리기 없는 GMRES는 조건수가 나쁜 야코비안에서 반복 수가 폭발한다. 그런데 JFNK에는 행렬이 없다. 해결책은 오른쪽 전처리다.

$$(J P^{-1})(P\,\delta\mathbf u) = -\mathbf F \quad\Rightarrow\quad J P^{-1}\mathbf v \approx \frac{\mathbf F(\mathbf u + \varepsilon P^{-1}\mathbf v) - \mathbf F(\mathbf u)}{\varepsilon}.$$

$P$ 는 $J$ 와 비슷하기만 하면 되고 정확할 필요가 없다. 그래서 진짜 야코비안은 행렬 없이 쓰고, 근사 야코비안만 행렬로 조립해 전처리에 쓰는 구성이 가능하다. PETSc의 -snes_mf_operator, MOOSE의 기본 solve_type = PJFNK가 이 방식이다.

크놀–키스가 강조한 물리 기반 전처리(physics-based preconditioning)가 이 자유를 활용한다. 예를 들어 강하게 결합된 다물리 문제에서 기존의 반음해법이나 연산자 분리 솔버 한 스텝을 그대로 $P^{-1}$ 로 쓴다. 그러면 원래 분리 솔버의 결합 오차는 바깥 뉴턴이 지워 주고, 분리 솔버는 크리로프 반복 수를 줄이는 역할만 한다. 레거시 코드를 버리지 않고 완전 음해 결합으로 올라가는 길이라 원자로·핵융합·복사 유체 코드에서 많이 쓰였다.

5. 검산 — 브라투 문제에서 세 손잡이 돌려 보기

2차원 브라투 문제 $-\Delta u - \lambda e^u = 0$ ($\lambda = 6$, 단위 정사각형, 내부 31×31 격자, 5점 차분, 디리클레 0)을 JFNK로 풀었다. GMRES는 재시작 60, 바깥 수렴 기준은 $\|\mathbf F\|$ 를 초기값의 $10^{-8}$ 배로. 최대 해 $u_{\max} = 0.79695$ 는 정확한 뉴턴(희소 직접 솔버, 5회 반복)과 다섯 자리까지 일치했다.

설정	뉴턴 반복	크리로프 반복 합	함수 평가 합
$\eta = 10^{-1}$ 고정	8	216	233
$\eta = 10^{-4}$ 고정	4	152	161
$\eta = 10^{-8}$ 고정	4	11,920	12,125
아이젠슈타트–워커	5	140	151
아이젠슈타트–워커 + 라플라시안 전처리	4	10	26

읽을거리가 셋이다. 과잉 풀이는 실재한다 — $\eta = 10^{-8}$ 은 뉴턴 반복 수는 최소지만 전처리 없는 재시작 GMRES가 그 정밀도에서 정체해 함수 평가가 80배로 뛰었다. 아이젠슈타트–워커가 손 튜닝 없이 최적 근처에 갔다. 그리고 전처리가 압도적이다 — 비선형 항을 뺀 라플라시안 $L$ 을 한 번 분해해 $P$ 로 쓰기만 해도 크리로프 반복이 140에서 10으로 줄었다. 물리 기반 전처리의 축소판이다.

섭동 크기도 확인했다. $\varepsilon = 10^{-2}/\|\mathbf v\|$ 로 크게 잡으면 크리로프 반복이 169로 조금 늘었을 뿐 수렴은 했다(바깥 잔차는 정확히 평가하므로 $J\mathbf v$ 오차가 수렴 자체를 막지는 않는다). 반대로 $\varepsilon = 10^{-14}/\|\mathbf v\|$ 로 기계정밀도 근처까지 줄이자 차분이 소거 오차에 먹혀 GMRES가 붕괴하고 NaN이 나왔다.³

6. 전역화와 실무

전역화. 뉴턴-크리로프도 국소 기법이라 멀리서 출발하면 장치가 필요하다. 잔차 노름에 대한 되돌이 선탐색이 기본이고, 정상 상태 CFD에서는 유사 과도 연속법(pseudo-transient continuation, ΨTC)이 국룰이다.

$$\left(\frac{1}{\Delta\tau_k}I + J\right)\delta\mathbf u = -\mathbf F(\mathbf u_k).$$

가상 시간 간격 $\Delta\tau$ 를 처음엔 작게(안정한 시간 전진), 잔차가 줄면 크게(뉴턴에 접근) 키운다. 흔한 규칙은 $\Delta\tau_k \propto 1/\|\mathbf F(\mathbf u_k)\|$ 다(SER, switched evolution relaxation). 켈리·키스(1998)가 수렴 이론을 정리했다.

구현체. PETSc SNES(-snes_mf, -snes_mf_operator), SUNDIALS의 KINSOL, Trilinos의 NOX가 대표적이다. 시간 적분기 안의 음해 스텝도 같은 기계를 쓴다 — 강성 방정식용 BDF 적분기의 매 스텝 비선형 풀이가 바로 뉴턴-크리로프다.

언제 쓰지 말까. 야코비안을 싸게 정확히 조립할 수 있고 희소 직접 솔버가 메모리에 들어가면 그냥 뉴턴 + 직접 풀이가 더 견고하다. 함수 평가가 매우 비싸고 크기가 작으면 브로이던 방법이 나을 수 있다. JFNK의 진짜 강점은 수백만 자유도 이상, 야코비안 조립이 고통스럽고, 쓸 만한 전처리기가 있는 문제다.

7. 관련 문서

• 뉴턴-랩슨법 · 크리로프 부분공간법 · GMRES · 자코비안 행렬
• 전처리기 · 다중격자법 · 연산자 분리 · 다물리 연성해석
• 브로이던 방법 · 준-뉴턴법 · 수치미분 · 자동 미분
• PETSc · MOOSE · 강성 방정식 · 비정확 뉴턴법

8. Footnotes

1. 1982년 논문 제목이 그냥 Inexact Newton methods 다. 한 단어짜리 제목으로 분야 이름을 선점한, 짧고 강한 논문의 표본이다. ↩

2. 실무자에게 아이젠슈타트–워커는 “튜닝하기 귀찮을 때 켜 두면 대충 맞는 기본값”이다. PETSc에서는 -snes_ksp_ew 한 줄이면 켜진다. ↩

3. 그래서 JFNK 코드가 NaN을 뱉으면 물리를 의심하기 전에 $\varepsilon$ 부터 찍어 보는 게 순서다. 대부분의 “수치 불안정”은 차분 간격 버그다. ↩