채프먼-엔스코그 전개

편집 역사 토론
수치해석 전산유체역학 유체역학 마지막 수정: 2026-07-27 04:07:12

1. 개요[편집]

채프먼-엔스코그 전개(Chapman-Enskog expansion)는 분자 하나하나의 통계를 다루는 볼츠만 수송방정식으로부터, 크누센수를 작은 매개변수로 삼는 다중척도 점근 전개를 통해 오일러 방정식·나비에-스토크스 방정식 같은 연속체 방정식과 그 수송계수를 유도해 내는 방법이다. 1916~17년 시드니 채프먼과 다비트 엔스코그가 각각 독립적으로 완성했으며1, “왜 기체 분자의 아수라장이 매끄러운 유체 방정식으로 보이는가”에 대한 인류의 표준 답변이다.

CFD 하는 사람에게 이 전개가 중요한 이유는 따로 있다. 격자 볼츠만 방법이 “가상 입자 놀이”를 하는데도 진짜 유체 해를 내놓는 근거가 바로 이것이고, LBM 사용자가 매일 손으로 만지는 완화시간 τ\tau 와 점성 ν\nu 의 관계식 역시 이 전개의 1차 항에서 떨어져 나온 물건이기 때문이다. 즉 LBM 코드에서 tau = 3*nu + 0.5 한 줄을 아무 생각 없이 치는 순간, 우리는 이미 채프먼-엔스코그의 결과를 쓰고 있는 셈이다.

화면은 D2Q9 BGK 격자 볼츠만법으로 실제로 푼 원기둥 주위 유동이다. 이 커널이 점성을 지정하는 방법은 딱 하나 — 완화시간 τ 를 정하는 것뿐이고, 그 τ 가 거시 점성으로 번역되는 환율 ν = c_s²(τ − 1/2)Δt 가 바로 채프먼-엔스코그 전개의 크누센수 1차 항에서 나오는 결과다. Re=90 은 후류가 아직 크게 요동치지 않는 매끄러운 영역이라, 격자 간격이 평균자유행로보다 압도적으로 큰 연속체 극한(Kn≪1)에 해당한다 — 전개가 가장 잘 먹히는 바로 그 구간이다. Re 슬라이더를 올리면 τ 가 1/2 쪽으로 붙으면서 전개의 전제와 수치 안정성이 동시에 나빠진다.

2. 크누센수와 척도 분리[편집]

전개의 소전제는 하나다. 분자 하나가 충돌 사이에 날아가는 거리(평균자유행로 λ\lambda)가, 우리가 보고 싶은 유동 구조의 크기 LL 보다 압도적으로 작을 것.

Kn=λL≪1\mathrm{Kn} = \frac{\lambda}{L} \ll 1

이 무차원수가 크누센수다. Kn≲10−2\mathrm{Kn} \lesssim 10^{-2} 면 연속체 영역이라 나비에-스토크스가 잘 맞고, 10−2∼10−110^{-2} \sim 10^{-1} 은 미끄럼 유동, 그 위는 천이·자유분자 영역이라 연속체 방정식이 무너진다. MEMS 채널, 진공 장비, 재진입 기체가 후자에 해당하며, 이때는 아예 볼츠만 쪽에서 직접 푸는 DSMC 같은 도구로 갈아탄다.

전개의 두 번째 아이디어가 시간 척도의 분리다. 충돌은 빠르고 유동은 느리다. 그래서 시간 미분을 척도별로 쪼갠다.

f=f(0)+Kn f(1)+Kn2f(2)+⋯ ,∂t=∂t0+Kn ∂t1+⋯f = f^{(0)} + \mathrm{Kn}\, f^{(1)} + \mathrm{Kn}^2 f^{(2)} + \cdots, \qquad \partial_t = \partial_{t_0} + \mathrm{Kn}\, \partial_{t_1} + \cdots

3. 차수별로 무엇이 나오는가[편집]

무차원화한 볼츠만 방정식 ∂tf+ξ⋅∇f=Kn−1Q(f,f)\partial_t f + \boldsymbol{\xi}\cdot\nabla f = \mathrm{Kn}^{-1} Q(f,f) 에 위 전개를 대입하고 크누센수 차수별로 정리하면, 유체역학 교과서의 계보가 순서대로 튀어나온다.

  • Kn−1\mathrm{Kn}^{-1}: Q(f(0),f(0))=0Q(f^{(0)}, f^{(0)}) = 0. 충돌항의 영공간이라 해는 국소 맥스웰 분포 하나뿐이다. 여기서 밀도·속도·온도라는 다섯 개의 거시 변수가 등장한다.
  • Kn0\mathrm{Kn}^{0}: 압력이 등방적인 오일러 방정식. 점성도 열전도도 없다. 즉 비점성 유동은 볼츠만의 0차 근사다.
  • Kn1\mathrm{Kn}^{1}: 선형화된 충돌 연산자를 가해성 조건(질량·운동량·에너지 모멘트가 0) 아래 뒤집어 f(1)f^{(1)} 을 구하면, 변형률에 비례하는 뉴턴 점성응력과 온도구배에 비례하는 푸리에 열유속이 유도된다. 나비에-스토크스-푸리에가 여기서 완성된다.
  • Kn2\mathrm{Kn}^{2} 이상: 버넷·초버넷 방정식. 미끄럼 유동에서 정확도가 조금 올라가지만, 짧은 파장에서 선형 불안정(보빌레프 불안정)이 알려져 있어 그대로 쓰면 폭발한다. 그래서 실무는 정규화 버넷이나 그라드 13-모멘트, R13 같은 모멘트법 계열로 우회한다.2

1차 결과의 위력은 수송계수가 공짜로 따라 나온다는 데 있다. 지름 dd 인 강체구 기체의 1차 근사 점성은

μ=516d2mkBTπ\mu = \frac{5}{16 d^2}\sqrt{\frac{m k_B T}{\pi}}

로, 밀도에 무관하고 T\sqrt{T} 로 증가한다. 압력을 두 배 올려도 점성이 그대로라는 맥스웰의 반직관적 예측이 여기서 나오고, 실험이 그걸 지지했다. 단원자 기체의 프란틀수가 Pr=2/3\mathrm{Pr} = 2/3 로 떨어지는 것도 같은 계산의 부산물이다.

4. LBM에서 ν = c_s²(τ − 1/2)Δt 가 나오는 자리[편집]

격자 볼츠만 방법의 BGK 갱신식은 완전히 이산적이다.

fi(x+eiΔt, t+Δt)−fi(x,t)=−1τ^(fi−fieq)f_i(\mathbf{x} + \mathbf{e}_i \Delta t,\, t + \Delta t) - f_i(\mathbf{x}, t) = -\frac{1}{\hat\tau}\left(f_i - f_i^{\mathrm{eq}}\right)

좌변을 Δt\Delta t 로 테일러 전개하면 1차 항 Δt(∂t+ei⋅∇)fi\Delta t(\partial_t + \mathbf{e}_i\cdot\nabla)f_i 뒤에 2차 항 Δt22(∂t+ei⋅∇)2fi\tfrac{\Delta t^2}{2}(\partial_t + \mathbf{e}_i\cdot\nabla)^2 f_i 가 남는다. 여기에 위와 똑같이 fi=fieq+εfi(1)+ε2fi(2)f_i = f_i^{\mathrm{eq}} + \varepsilon f_i^{(1)} + \varepsilon^2 f_i^{(2)}, ∂t=ε∂t1+ε2∂t2\partial_t = \varepsilon\partial_{t_1} + \varepsilon^2\partial_{t_2} 를 넣고 차수를 맞추면, ε2\varepsilon^2 방정식의 유속 항 앞에 인자 (1−12τ^)\left(1 - \frac{1}{2\hat\tau}\right) 가 붙는다. 이 인자를 2차 모멘트로 옮기면 그대로

ν=cs2(τ^−12)Δt,cs2=13Δx2Δt2 (D2Q9)\nu = c_s^2\left(\hat\tau - \frac{1}{2}\right)\Delta t, \qquad c_s^2 = \frac{1}{3}\frac{\Delta x^2}{\Delta t^2}\ (\text{D2Q9})

가 된다. 그 악명 높은 −1/2-1/2 는 물리가 아니라 이산 격자의 2차 테일러 잔여항이다. 이걸 흡수해 준 덕에 LBM은 1차가 아니라 2차 정확도를 갖고, 동시에 τ^>1/2\hat\tau > 1/2 라는 하한이 강제된다. 레이놀즈수를 올릴수록 τ^→1/2+\hat\tau \to 1/2^{+} 로 붙으면서 유효 점성이 0에 수렴하고, 그 지점에서 SRT-LBM이 불안정해지는 것도 같은 식이 말해 주는 사실이다.3 SRT에서는 체적점성까지 같은 τ^\hat\tau 하나에 묶여 따로 못 만지는데, MRT가 완화율을 모멘트별로 쪼개 이 족쇄를 푸는 것도 채프먼-엔스코그로 각 모멘트의 역할을 확인했기에 가능한 설계다.

5. 격자 기체 오토마타와 갈릴레이 불변성[편집]

이 전개가 진짜 무기로 쓰인 첫 무대는 LBM이 아니라 격자 기체 오토마타였다. FHP류 LGCA는 점유수 ni∈{0,1}n_i \in \{0,1\} 만 다루는 셀룰러 오토마타인데, 앙상블 평균 후 채프먼-엔스코그를 돌려 보면 이류항 앞에 밉살스러운 인자가 하나 붙는다.

∂tu+g(ρ) (u⋅∇)u=−∇P+ν∇2u,g(ρ)=DD+2⋅1−2d1−d\partial_t \mathbf{u} + g(\rho)\,(\mathbf{u}\cdot\nabla)\mathbf{u} = -\nabla P + \nu\nabla^2\mathbf{u}, \qquad g(\rho) = \frac{D}{D+2}\cdot\frac{1-2d}{1-d}

여기서 dd 는 방향당 점유확률이다. g≠1g \neq 1 이라는 것은 갈릴레이 불변성이 깨졌다는 뜻이다. 등속으로 움직이는 좌표계에서 물리가 달라진다는, 유체 코드가 절대 하면 안 되는 짓을 격자가 하고 있는 것이다. LBM이 페르미-디랙형 점유수를 실수 분포함수로 바꾸고 2차 다항식 평형분포를 도입한 결정적 이유가 이 g(ρ)g(\rho) 를 1로 되돌리기 위해서였다.

그런데 완전히 고쳐진 것은 아니다. D2Q9 격자는 4차 모멘트 텐서까지는 등방적이지만 3차 모멘트 ∑iwieiαeiβeiγ\sum_i w_i e_{i\alpha}e_{i\beta}e_{i\gamma} 를 연속체 값과 일치시키지 못한다. 그 결과 운동량 방정식에 O(u3)O(u^3) 짜리 오차항(입방 속도 결함)이 남고, 이것이 표준 LBM에 저마하수 제약이 붙는 근본 원인이다. 마하수를 0.1~0.2 아래로 유지하라는 국룰은 경험칙이 아니라 이 전개가 계산해 준 오차 항의 크기다.4

6. 한계와 오해[편집]

  • 점근 급수이지 수렴 급수가 아니다. 차수를 올린다고 정확해지는 게 아니라, 어느 지점부터 나빠진다. 버넷이 불안정한 것이 그 대표적 증상이다.
  • 가해성 조건이 본질이다. 각 차수에서 f(n)f^{(n)} 의 보존 모멘트를 0으로 두는 선택이 없으면 해가 유일하지 않다. 이 선택이 “밀도·속도·온도는 전부 f(0)f^{(0)} 이 정의한다”는 물리적 규약이다.
  • 강한 비평형에는 안 통한다. 충격파 내부, 재진입 유동, 진공 장비처럼 국소 Kn\mathrm{Kn} 이 커지는 곳에서는 1차 근사가 무너진다. 압축성 유동 코드가 충격파 두께를 물리적으로 재현하지 못하는 이유이기도 하다.
  • LBM은 볼츠만 방정식을 푸는 게 아니다. 정확히는 “채프먼-엔스코그를 거쳐 나비에-스토크스와 점근적으로 일치하도록 설계된 이산 동역학”을 푸는 것이다. 미시 물리가 진짜인 척하면 곤란하다.

7. 관련 문서[편집]

8. Footnotes[편집]

  1. 두 사람은 서로의 작업을 모른 채 거의 동시에 도달했고, 방법론도 달랐다. 채프먼은 맥스웰식 수송방정식을, 엔스코그는 볼츠만 방정식의 체계적 섭동 해법을 썼다. 지금 교과서에 실린 정돈된 형태는 엔스코그 쪽에 가깝다. ↩

  2. 버넷 방정식은 “차수를 올렸으니 더 정확하겠지”라는 순진한 기대를 배신하는 대표 사례다. 짧은 파장 교란에 대해 선형 불안정이라, 격자를 촘촘히 할수록 더 빨리 터진다. 격자 수렴 연구를 돌리다 발산하는 걸 보고 코드를 의심하면 안 되고 방정식을 의심해야 하는, 흔치 않은 상황. ↩

  3. 그래서 LBM 튜닝의 국룰은 “레이놀즈수를 올리고 싶으면 τ\tau 를 1/2로 밀지 말고 격자를 늘려라”이다. 물론 격자를 늘리면 시간 스텝도 같이 줄어들어 계산량이 순식간에 갈려 나간다. 세상에 공짜는 없다. ↩

  4. 3차 모멘트 결함을 없애려면 속도 방향을 더 준 격자(D2Q17, D2Q37 등)를 써야 한다. 열유동·고마하수용 고차 격자가 존재하는 이유가 이것인데, 방향 수가 늘면 메모리와 대역폭이 그대로 늘어난다. 결국 “어느 모멘트까지 살릴 것인가”가 격자 설계의 전부다. ↩