크랭크-니콜슨법

1. 개요

크랭크-니콜슨법(Crank-Nicolson method)은 시간에 대한 편미분방정식을 이산화할 때 현재 시간 단계와 다음 시간 단계의 공간 항을 평균(사다리꼴 적분) 하여 시간 전진시키는 반음해(semi-implicit) 시간 적분 도식이다. 1947년 존 크랭크(John Crank)와 필리스 니콜슨(Phyllis Nicolson)이 열전달 해석의 수치 계산을 위해 제안했으며¹, 시간에 대해 2차 정확도를 가지면서 무조건 안정(A-안정)이라는 두 마리 토끼를 동시에 잡은 도식으로 확산 문제의 사실상 표준이 되었다.

명시적 오일러의 가벼운 스텝과 음해 오일러의 안정성을 반씩 섞은 도식이다. 공짜는 없어서, 시간 간격을 키우면 발산 대신 진동이 온다.

그 진동은 이렇게 생겼다. 주기 경계 32칸, 한 칸만 1 인 초기조건을 $r = \alpha \Delta t / \Delta x^2 = 10$ 으로 한 스텝 밀면 뜨거웠던 칸이 −0.564 가 되고 이웃이 +0.280 이 된다. 가장 잘게 출렁이는 모드의 증폭인자는 $(1-2r)/(1+2r) = -0.905$ 라, 부호는 매 스텝 뒤집히면서 크기는 거의 안 죽는다. 같은 모드를 전진 오일러로 $r=0.6$ 에 밀면 증폭인자가 −1.4 가 되어 진동하다 터진다. 크랭크-니콜슨은 터지지 않을 뿐, 큰 $r$ 에서 체스판 무늬를 얌전히 지워 주지는 않는다. 격자를 바꿔 가며 보는 전진 오일러 열확산은 메시 독립성 연구 쪽에 있다.

2. 세타 방법과의 관계

크랭크-니콜슨법은 더 일반적인 세타(θ) 방법의 특수한 경우다. 확산 방정식 $\partial u/\partial t = \alpha\, \partial^2 u/\partial x^2$ 을 공간 이산화한 뒤 시간 항을 다음처럼 쓰는다.

$$\frac{u^{n+1} - u^n}{\Delta t} = \alpha\left[(1-\theta)\,L u^n + \theta\, L u^{n+1}\right]$$

여기서 $L$은 공간 2계 미분 연산자(예: 중심차분)이고 $\theta \in [0,1]$가 음해성의 정도를 조절한다.

• $\theta = 0$ — 명시적 오일러(전진 오일러). 1차 정확도, 조건부 안정.
• $\theta = 1$ — 음해 오일러(후진 오일러). 1차 정확도, 무조건 안정.
• $\theta = 1/2$ — 크랭크-니콜슨법. 시간 2차 정확도, 무조건 안정.

즉 크랭크-니콜슨법은 현재와 미래의 기울기를 정확히 반씩 신뢰한다. 이 “반반”이 사다리꼴 공식(trapezoidal rule)과 수학적으로 같은 구조이기 때문에 시간 오차가 $O(\Delta t^2)$로 떨어지는 것.

3. 정확도와 안정성

크랭크-니콜슨법의 국소 절단 오차는 공간·시간 모두 2차, 즉 $O(\Delta t^2 + \Delta x^2)$이다. 명시적/음해 오일러가 시간에 대해 1차인 것과 대비되는 결정적 장점이다.

안정성은 폰 노이만 안정성 해석으로 확인한다. 파수 $k$에 대한 증폭 인자(amplification factor)는

$$G = \frac{1 - 2r\sin^2(k\Delta x/2)}{1 + 2r\sin^2(k\Delta x/2)}, \qquad r = \frac{\alpha \Delta t}{\Delta x^2}$$

로 주어진다. 분자·분모를 보면 임의의 $r>0$에 대해 $|G| \le 1$이 항상 성립한다. 즉 CFL 조건 같은 시간 간격 제약 없이 무조건 안정하다. 선형 안정성 이론에서 실수 음의 고유값 전체에 대해 안정한 성질을 A-안정(A-stability)이라 부르는데, 크랭크-니콜슨은 이를 만족한다.

4. 스퓨리어스 오실레이션

문제는 여기서 시작된다. “무조건 안정 = 아무 $\Delta t$나 써도 됨”이 절대 아니다. $r$이 커지면(대략 $r > 1/2$) 위 증폭 인자 $G$가 음수가 되어, 매 시간 단계마다 부호가 뒤집히는 톱니 모양 진동이 생긴다. 초기 조건에 급격한 불연속(예: 스텝 함수)이 있을 때 특히 심하다.²

발산하지는 않으니 계산이 폭발하진 않지만, 물리적으로 존재할 수 없는 언더슈트/오버슈트(예: 음의 농도, 절대영도 이하 온도)가 나타난다. 이는 A-안정보다 강한 성질인 **L-안정(L-stability)**을 크랭크-니콜슨이 만족하지 못하기 때문이다. 고주파 성분에 대해 $|G| \to 1$이라 감쇠가 거의 없어서 진동이 오래 살아남는다.

실무적 해결책은 다음과 같다.

• 진동이 문제되는 초반 몇 스텝만 음해 오일러($\theta=1$)로 돌려 고주파를 죽인 뒤 전환하는 란초스(Rannacher) 시동.
• $\theta$를 $1/2$보다 살짝 크게(예: $0.5 + \epsilon$) 잡아 약간의 수치 감쇠를 도입.
• 적응 격자 세분화나 시간 간격 축소로 셀 단위 $r$을 낮추기.

5. 구현과 활용

시간 항이 음해이므로 각 시간 단계마다 선형 연립방정식을 풀어야 한다. 1차원 확산 문제에서는 삼중대각(tridiagonal) 희소행렬이 나오고, 이는 토마스 알고리즘(LU 분해의 삼중대각 특수형)으로 $O(N)$에 풀린다. 2차원 이상에서는 크리로프 부분공간법이나 다중격자법이 동원된다.

적용 분야는 확산이 지배하는 거의 모든 곳이다. 열전달 해석의 비정상 전도, 대류-확산 방정식의 확산 항, 옵션 가격 결정의 블랙-숌즈 방정식(금융공학에서 크랭크-니콜슨은 사실상 디폴트다), 슈뢰딩거 방정식의 시간 전개(이 경우 유니터리성을 보존하는 성질이 특히 유용하다) 등. 순수 대류 문제에서는 오히려 진동 때문에 별로 인기가 없다는 점은 알아두자.

6. 실무 팁

• 초기 불연속이 있으면 Rannacher 시동으로 몇 스텝 음해 오일러를 쓰라.
• 큰 디를타 t에서 진동이 보이면 θ를 약간 크게 조정하라.
• 1D에서는 토마스 알고리즘으로 빠르게 풀고, 2D+에서는 적절한 선형 솔버를 선택하라.
• 순수 대류 문제에는 다른 도식을 고려하라.

7. 관련 문서

• 차분 도식
• 유한차분법
• 폰 노이만 안정성 해석
• CFL 조건
• 대류-확산 방정식
• 열전달 해석
• 룽게-쿠타법
• 뉴마크법
• 수치적분

8. Footnotes

1. Crank, J., & Nicolson, P. (1947). A practical method for numerical evaluation of solutions of partial differential equations of the heat-conduction type. 참고로 필리스 니콜슨은 이 논문 하나로 이름이 영원히 남았다. 수치해석 도식에 자기 이름이 붙는 것만큼 확실한 불멸은 없다. ↩

2. 급격한 초기 불연속을 “함부로 밟으면 빠걱대는 낡은 마루” 정도로 생각하면 된다. A-안정이라 무너진 않지만, 밟을 때마다 빠걱빠걱(진동) 소리가 난다. ↩