ADI 법

편집 역사 토론
수치해석 전산유체역학 마지막 수정: 2026-08-02 04:55:41

1. 개요[편집]

ADI 법
Alternating Direction Implicit method
한국어교대방향 암시법
원 논문Peaceman & Rachford (1955)
핵심한 스텝을 방향별로 쪼개 매 반스텝마다 삼중대각계만 푼다
스텝당 비용$O(N)$ per line — 토머스 알고리즘
정확도 / 안정성2차 · 2D에서 무조건 안정 (3D는 별도 도식 필요)

2차원을 한 번에 암시적으로 푸는 건 비싸다. 그러면 1차원을 두 번 풀면 되잖아?

ADI 법(Alternating Direction Implicit method, 교대방향 암시법)은 다차원 확산형 방정식을 암시적으로 적분할 때 한 시간 스텝을 방향 수만큼의 부분 스텝으로 쪼개, 각 부분 스텝에서 한 방향만 암시적으로 처리해 삼중대각 연립방정식만 풀도록 만드는 시간 전진 기법이다. 1955년 피스먼과 랙포드가 석유 저류층 유동 계산을 위해 제안했고1, 이후 30년 동안 구조격자 기반 전산유체역학·열전달 해석 코드의 사실상 표준 선형해법기였다.

발상의 출발점은 지극히 실용적이다. 2D 확산을 크랭크-니콜슨법으로 풀면 계수행렬이 **오대각(pentadiagonal)**이 된다. N×NN \times N 격자면 미지수 N2N^2개에 대역폭 NN이라, 대역 LU 분해로 직접 풀면 O(N4)O(N^4)이 든다. 반면 방향별로 쪼개면 각 격자선이 독립적인 삼중대각계가 되고, 이건 토머스 알고리즘으로 선당 O(N)O(N), 반스텝당 O(N2)O(N^2) — 즉 미지수 개수에 선형이다. 무려 N2N^2배의 차이다.

2. 피스먼-랙포드 도식[편집]

2D 열방정식 ut=α(uxx+uyy)u_t = \alpha(u_{xx} + u_{yy})를 균일 격자에서 이산화한다. δx2ui,j=ui+1,j−2ui,j+ui−1,j\delta_x^2 u_{i,j} = u_{i+1,j} - 2u_{i,j} + u_{i-1,j}, μ=αΔt/h2\mu = \alpha\Delta t/h^2로 두면 피스먼-랙포드 도식은 두 반스텝으로 구성된다.

un+1/2−unΔt/2=αh2(δx2un+1/2+δy2un)\frac{u^{n+1/2} - u^{n}}{\Delta t/2} = \frac{\alpha}{h^2}\left(\delta_x^2 u^{n+1/2} + \delta_y^2 u^{n}\right) un+1−un+1/2Δt/2=αh2(δx2un+1/2+δy2un+1)\frac{u^{n+1} - u^{n+1/2}}{\Delta t/2} = \frac{\alpha}{h^2}\left(\delta_x^2 u^{n+1/2} + \delta_y^2 u^{n+1}\right)

첫 반스텝은 xx 암시·yy 명시, 둘째 반스텝은 정확히 반대다. 연산자 형태로 정리하면

(I−μ2δx2)un+1/2=(I+μ2δy2)un,(I−μ2δy2)un+1=(I+μ2δx2)un+1/2\left(I - \tfrac{\mu}{2}\delta_x^2\right)u^{n+1/2} = \left(I + \tfrac{\mu}{2}\delta_y^2\right)u^{n}, \qquad \left(I - \tfrac{\mu}{2}\delta_y^2\right)u^{n+1} = \left(I + \tfrac{\mu}{2}\delta_x^2\right)u^{n+1/2}

좌변 행렬이 각각 xx-선, yy-선 방향의 삼중대각 행렬이다. 첫 반스텝에서는 y=y = const인 행마다, 둘째 반스텝에서는 x=x = const인 열마다 독립적으로 풀면 된다. 각 선이 서로 독립이므로 병렬화도 자연스럽다 — 다만 방향을 바꿀 때 **데이터 전치(transpose)**가 필요하고, 이게 분산 메모리 환경에서 ADI의 발목을 잡는다.2

3. 왜 무조건 안정인가[편집]

두 반스텝을 소거하면(상수계수라 두 연산자가 교환한다) 다음 인수분해 형태가 나온다.

(I−μ2δx2)(I−μ2δy2)un+1=(I+μ2δx2)(I+μ2δy2)un\left(I - \tfrac{\mu}{2}\delta_x^2\right)\left(I - \tfrac{\mu}{2}\delta_y^2\right)u^{n+1} = \left(I + \tfrac{\mu}{2}\delta_x^2\right)\left(I + \tfrac{\mu}{2}\delta_y^2\right)u^{n}

이것은 2D 크랭크-니콜슨에 μ24δx2δy2(un+1−un)\frac{\mu^2}{4}\delta_x^2\delta_y^2(u^{n+1} - u^{n})라는 여분 항이 붙은 것과 같다. 이 항은 O(Δt3)O(\Delta t^3)이므로 국소 절단오차의 차수를 떨어뜨리지 않는다. 2차 정확도가 유지된다.

폰 노이만 안정성 해석을 적용해 보자. u∼Gnei(kxx+kyy)u \sim G^n e^{\mathrm{i}(k_x x + k_y y)}를 대입하면 δx2→−4sin⁡2(θx/2)\delta_x^2 \to -4\sin^2(\theta_x/2) (θx=kxh\theta_x = k_x h)이므로, X=2μsin⁡2(θx/2)X = 2\mu\sin^2(\theta_x/2), Y=2μsin⁡2(θy/2)Y = 2\mu\sin^2(\theta_y/2)로 두면 증폭인자가 두 방향의 곱으로 깔끔하게 갈라진다.

G=1−X1+X⋅1−Y1+YG = \frac{1 - X}{1 + X}\cdot\frac{1 - Y}{1 + Y}

μ>0\mu > 0이면 X,Y≥0X, Y \ge 0이고, 각 인자는 항상 절댓값이 1 이하다. 따라서 ∣G∣≤1\lvert G \rvert \le 1이 Δt\Delta t와 무관하게 성립한다 — 무조건 안정. 이 “곱으로 분리된다”는 구조가 ADI의 전부라고 해도 과언이 아니다.3

문제는 3D다. 같은 방식으로 세 방향으로 확장한 피스먼-랙포드는 증폭인자가 위처럼 예쁘게 인수분해되지 않고, ∣G∣>1\lvert G \rvert > 1인 파수 조합이 존재한다. 즉 3D 피스먼-랙포드는 조건부 안정이다. 그래서 3D에서는 더글라스-랙포드(1956) 계열을 쓴다. 이쪽은 무조건 안정이지만 시간 1차이고, 2차 정확도까지 원하면 브라이언(1961)·더글라스(1962)의 변형이 필요하다.

더글라스-건(Douglas-Gunn, 1964) 형식은 이 난립을 정리한 일반 틀이다. 도식을 증분(delta) 형태로 쓴다는 것이 핵심으로, 미지수를 Δu=un+1−un\Delta u = u^{n+1} - u^n로 놓으면 우변이 곧 잔차가 되어 (1) 정상해에 수렴할 때 분할오차가 자동으로 사라지고, (2) 경계조건 처리가 훨씬 자연스러워지며, (3) 변수계수·혼합미분항으로도 확장이 쉽다.

4. 연산자 분리와의 관계[편집]

ADI는 연산자 분리의 한 종류지만, 스트랑 분할처럼 ”xx 문제를 풀고 그 결과를 yy 문제의 초기값으로 넘기는” 순차 분할이 아니라 **근사 인수분해(approximate factorization)**형 분할이다. 즉 (I−μ2(δx2+δy2))(I - \frac{\mu}{2}(\delta_x^2 + \delta_y^2))를 (I−μ2δx2)(I−μ2δy2)(I - \frac{\mu}{2}\delta_x^2)(I - \frac{\mu}{2}\delta_y^2)로 근사한 것이고, 두 표현의 차이가 곧 분할오차다.

상수계수 확산에서는 δx2\delta_x^2와 δy2\delta_y^2가 교환하므로 이 오차가 위에서 본 O(Δt3)O(\Delta t^3) 항 하나로 정리되고 차수 손실이 없다. 그러나 두 연산자가 교환하지 않으면 분할오차가 교환자 [A,B]=AB−BA[A, B] = AB - BA에 비례하고, 이게 사라지지 않으면 시간 차수가 떨어진다. 실전에서 교환자가 살아나는 경우는 다음과 같다.

  • 변수계수 확산(비균질 물성, 온도의존 열전도도)
  • 혼합미분항 uxyu_{xy}가 있는 경우(비직교 곡선격자에서 항상 생긴다)
  • 방향별로 다른 경계조건이 강하게 작용하는 경우

그래서 곡선격자 CFD 코드는 ADI를 그대로 쓰지 않고 혼합항을 우변으로 넘겨 명시적으로 처리하거나, 부분 스텝을 더 얹는 식으로 보정한다.

5. 경계조건 — 제일 많이 틀리는 곳[편집]

ADI에서 가장 악명 높은 함정. 중간값 un+1/2u^{n+1/2}는 시각 tn+Δt/2t^n + \Delta t/2의 물리적 해가 아니다. 그건 인수분해 과정에서 튀어나온 순수한 중간 변수다. 그런데 대부분의 사람이 첫 반스텝의 경계에 ubn+1/2=g(tn+1/2)u^{n+1/2}_b = g(t^{n+1/2})를 그냥 넣는다. 그 순간 전체 도식의 시간 정확도가 2차에서 1차로 떨어진다. 격자를 아무리 조여도 오차가 Δt\Delta t에 선형으로만 줄어드는, 원인을 찾기 극도로 어려운 종류의 버그다.4

올바른 중간 경계조건은 두 반스텝 관계식에서 직접 유도된다. 둘째 식을 (I+μ2δx2)un+1/2=(I−μ2δy2)un+1\left(I + \frac{\mu}{2}\delta_x^2\right)u^{n+1/2} = \left(I - \frac{\mu}{2}\delta_y^2\right)u^{n+1}로 다시 쓰고 첫 식과 변끼리 더하면 좌변에서 δx2\delta_x^2 항이 상쇄되어

un+1/2∣경계=12[(I+μ2δy2)un+(I−μ2δy2)un+1]경계u^{n+1/2}\Big|_{\text{경계}} = \frac{1}{2}\left[\left(I + \tfrac{\mu}{2}\delta_y^2\right)u^{n} + \left(I - \tfrac{\mu}{2}\delta_y^2\right)u^{n+1}\right]_{\text{경계}}

디리클레 경계라면 우변의 unu^n, un+1u^{n+1}이 모두 알려진 값이므로 그대로 계산하면 된다. 경계를 따라 접선 방향 2차 미분이 들어간다는 것 — 즉 중간 경계값은 시간 보간이 아니라 접선 방향 확산이 섞인 값이라는 것이 핵심이다. 노이만 경계나 시간의존 경계에서는 유도가 더 지저분해지는데, 더글라스-건 증분 형식을 쓰면 경계에서 Δu=0\Delta u = 0(디리클레 고정) 같은 자연스러운 조건으로 바뀌어 이 문제가 상당 부분 사라진다.

6. 실전 구현 메모[편집]

교과서 도식과 실제 코드 사이에는 늘 몇 개의 함정이 더 있다.

  • 대류항이 섞이면 삼중대각이 비대칭이 된다. ut+aux=αuxxu_t + a u_x = \alpha u_{xx}를 중심차분하면 부대각 원소가 −μ/2∓σ/2-\mu/2 \mp \sigma/2 꼴이 되고, 셀 페클레 수 Peh=ah/α>2Pe_h = ah/\alpha > 2에서 대각우세가 깨진다. 토머스 알고리즘은 피벗팅을 하지 않으므로 이때 수치가 진동하거나 아예 망가진다. 풍상 차분으로 대각우세를 회복시키든지, 대류항을 우변으로 넘겨 명시적으로 처리(IMEX)하는 것이 표준 대응이다.
  • 주기 경계는 순환 삼중대각계가 된다. 모서리에 원소 두 개가 더 붙어 순수 삼중대각이 아니게 되므로, 셔먼-모리슨 보정이나 순환 토머스(Sherman–Morrison / Thomas-with-correction)를 써야 한다. 이걸 모르고 일반 토머스를 그대로 쓰면 경계 한 줄만 조용히 틀린다.
  • LOD와 혼동하지 말 것. 국소 1차원(locally one-dimensional, LOD) 도식은 각 부분 스텝을 완결된 1D 문제로 취급해 순차적으로 푸는 진짜 연산자 분리이고, ADI는 하나의 다차원 도식을 인수분해한 것이다. LOD는 구현이 더 단순하지만 중간해가 물리적이라는 착각을 주기 쉽고, 일반적으로 시간 1차에 그친다.
  • 정상해를 원한다면 Δt\Delta t는 그냥 완화 파라미터다. 시간 정확도를 포기하고 정상상태만 노린다면 Δt\Delta t를 반복마다 바꾸는 것이 이득이다. 큰 값은 저주파 오차를, 작은 값은 고주파 오차를 잘 죽이므로, 값을 주기적으로 순환시키는 파라미터열이 수렴을 극적으로 가속한다. 다중격자법이 격자 레벨로 하는 일을 ADI는 시간 스텝 크기로 흉내 내는 셈이다.

7. 지금은 어디에 쓰나[편집]

1980년대까지 ADI는 SIMPLE 알고리즘 계열 CFD 코드에서 운동량 방정식과 포아송 방정식 계열 보정식을 푸는 기본 도구였다. 압축성 코드 쪽에서는 빔-워밍(Beam-Warming) 근사 인수분해가 같은 아이디어를 오일러/나비에-스토크스 방정식에 적용한 형태로, NASA의 ARC3D·OVERFLOW 계보를 만들었다. 엇갈림 격자 위에서의 압력 보정도 대개 라인 단위로 풀렸다.

오늘날 범용 해법기로서의 ADI는 다중격자법과 크리로프 부분공간법에 밀렸다. 이유는 명확하다. ADI의 수렴은 격자 종횡비와 조건수에 민감하고, 비정렬 격자에는 “방향”이라는 개념 자체가 없으며, 방향 전환마다 전치 통신이 필요하기 때문이다. 그럼에도 다음 자리에서는 여전히 현역이다.

  • 구조격자 + 강한 방향성. 경계층 격자처럼 종횡비 1000:1짜리 얇은 셀이 깔린 영역에서는 점 완화(point relaxation)가 벽 수직 방향 결합을 못 이겨 기어간다. 그 방향만 암시적으로 처리하는 **라인 완화(line relaxation)**가 정확히 ADI 반스텝 하나이며, 다중격자의 평활자(smoother)로 얹으면 궁합이 매우 좋다.
  • 전처리기로서의 ADI. 자기 자신이 해법기일 필요는 없다. ADI 스텝 하나를 GMRES/BiCGSTAB의 전처리로 쓰면 싸고 효과적이다.
  • 분리 가능한 포아송 방정식. 가속 파라미터열(Wachspress 파라미터)을 쓴 ADI 반복은 분리 가능한 문제에서 반복 횟수가 O(log⁡(1/ϵ))O(\log(1/\epsilon)) 수준으로 떨어져, 지금 기준으로도 빠르다.
  • 강성 방정식 맥락의 IMEX. 확산항만 방향별 암시로 처리하고 대류는 명시적으로 두는 조합.

요약하면 ADI는 “죽은 방법”이 아니라 범용 자리에서 물러나 특수 자리에 정착한 방법이다. 그리고 유한한 시간 안에 O(N)O(N)짜리 삼중대각 풀이로 무조건 안정성을 사는 그 아이디어는, 이후 등장한 거의 모든 근사 인수분해 도식의 원형으로 남아 있다.

8. 관련 문서[편집]

9. Footnotes[편집]

  1. 피스먼과 랙포드는 정유회사(Humble Oil) 연구소 소속이었다. 즉 ADI는 대학 수학과가 아니라 석유 시추 현장의 예산 압박에서 태어난 알고리즘이다. 저류층 시뮬레이션은 지금도 ADI 계열 전처리의 주요 고객이다. ↩

  2. 라인 solve 자체는 부끄러울 정도로 병렬적인데, xx 방향을 다 풀고 yy 방향으로 넘어가려면 도메인 전체를 뒤집어야 한다. MPI 관점에서는 전면적인 all-to-all 통신이고, 노드 수를 늘릴수록 통신이 계산을 잡아먹는다. ADI가 대규모 병렬 시대에 밀려난 결정적 이유 중 하나다. ↩

  3. 그래서 “ADI를 쓰면 무조건 안정하니까 Δt\Delta t를 무한히 키워도 된다”는 말은 반만 맞다. 안정하다고 정확한 것은 아니고, 특히 크랭크-니콜슨 계열이라 큰 Δt\Delta t에서 고주파 성분이 부호를 바꾸며 진동한다. L-안정이 아닌 것의 대가는 여기서도 똑같이 청구된다. ↩

  4. 페어웨더와 미첼이 1960년대에 이미 지적한 문제인데, 60년이 지난 지금도 학위논문 코드에서 꾸준히 재발견된다. 증상이 “왜인지 모르게 시간 수렴 차수가 1.0 근처”라서, 대부분은 경계조건을 의심하지 않고 격자를 더 조이다가 학기를 보낸다. ↩