조각포물선법

1. 개요

조각포물선법 Piecewise Parabolic Method (PPM)
고안	필립 콜렐라 · 폴 우드워드 (1984)
계보	고두노프(상수) → MUSCL(선형) → PPM(포물선)
재구성	셀 평균 보존 포물선, 경계값은 4차 보간
진동 억제	극값 평탄화 + 포물선 단조 제약 + 충격 근처 평탄화
부속 장치	접촉 불연속 가파르게 하기, 특성 추적
대표 코드	FLASH · Athena++ · Enzo · VH-1 · CASTRO

상수(고두노프)에서 직선(판 레이어)으로 왔다. 그럼 다음은 당연히 포물선이다. 문제는 포물선이 더 잘 튄다는 것.

조각포물선법(piecewise parabolic method, PPM)은 유한체적법에서 각 셀 내부의 해를 셀 평균을 보존하는 포물선으로 재구성하고, 그 포물선이 새로운 극값을 만들지 않도록 단조 제약을 거는 고해상도 도식이다. 1984년 로런스 버클리 연구소의 필립 콜렐라와 로런스 리버모어 연구소의 폴 우드워드가 압축성 유동·천체물리 기체역학용으로 발표했다. 고두노프 도식의 조각상수, MUSCL의 조각선형에 이은 재구성 차수의 다음 단계이며, 매끄러운 영역에서 같은 격자로 MUSCL보다 확연히 작은 오차를 낸다.

MUSCL 문서가 재구성-발전-평균(REA) 틀과 기울기 제한을 이미 정리했다. 이 문서는 PPM이 그 틀 위에 무엇을 더 얹었는가 — 4차 경계 보간, 포물선 단조 제약, 특성 추적, 그리고 원 논문의 잘 안 알려진 부속 장치들 — 와, 그 대가인 매끄러운 극값에서의 정확도 손실을 판다.

2. 재구성 — 포물선을 어떻게 고르나

등간격 격자에서 셀 $j$ 의 평균 $\bar a_j$ 가 주어졌다고 하자. PPM은 세 단계로 셀 내부 포물선을 정한다.

1단계: 경계값 보간. 셀 경계 $j+\tfrac12$ 의 값을 이웃 네 셀의 평균으로부터 4차 정확도로 보간한다.

$$a_{j+1/2} = \frac{7}{12}\left(\bar a_j + \bar a_{j+1}\right) - \frac{1}{12}\left(\bar a_{j-1} + \bar a_{j+2}\right).$$

원 논문은 이것을 판 레이어식 제한 기울기 $\delta a_j$ 를 써서 $a_{j+1/2} = \bar a_j + \tfrac12(\bar a_{j+1} - \bar a_j) - \tfrac16(\delta a_{j+1} - \delta a_j)$ 로 쓴다. $\delta a$ 가 제한되지 않은 중심 차분이면 위 식과 같고, 제한되면 경계값이 이웃 평균 사이에 머문다.

2단계: 포물선 맞추기. 셀 왼쪽 경계값 $a_L$, 오른쪽 $a_R$, 그리고 셀 평균 $\bar a_j$ 를 동시에 만족하는 포물선은 유일하다. 셀 안의 정규화 좌표 $\xi \in [0, 1]$ 로

$$a(\xi) = a_L + \xi\left[\Delta a + a_6(1 - \xi)\right], \qquad \Delta a = a_R - a_L, \quad a_6 = 6\left(\bar a_j - \tfrac12(a_L + a_R)\right).$$

$a_6$ 이 곡률 항이다. 셀 평균이 두 경계값의 평균보다 크면 위로 볼록해진다.

3단계: 단조 제약. 포물선이 셀 안에서 $a_L$ 이나 $a_R$ 을 넘는 극값을 가지면 진동의 씨앗이 된다. PPM의 제약은 두 경우로 나뉜다.

• $\bar a_j$ 가 국소 극값이면 — 즉 $(a_R - \bar a_j)(\bar a_j - a_L) \le 0$ — 포물선을 버리고 $a_L = a_R = \bar a_j$ 로 평탄화한다(1차 정확도로 떨어진다).
• 그렇지 않은데 포물선 꼭짓점이 셀 안에 들어오면, 꼭짓점이 셀 경계에 오도록 반대쪽 경계값을 고친다. $\Delta a \cdot a_6 > (\Delta a)^2$ 이면 $a_L \leftarrow 3\bar a_j - 2a_R$, $-(\Delta a)^2 > \Delta a \cdot a_6$ 이면 $a_R \leftarrow 3\bar a_j - 2a_L$.

이 제약 뒤의 포물선은 셀 안에서 단조이고, 새 극값을 만들지 않는다.

3. 발전 — 특성 추적

재구성만으로는 시간 적분이 끝나지 않는다. PPM은 셀 경계에서 리만 솔버에 넣을 좌우 상태를, 경계에 반 스텝 동안 도달할 수 있는 영역 위에서 포물선을 평균해 구한다. 속도 $\lambda > 0$ 인 특성에 대해 $\sigma = \lambda\Delta t/\Delta x$ 라 하면 경계 왼쪽 상태는

$$f_{j+1/2,L} = a_R - \frac{\sigma}{2}\left[\Delta a - \left(1 - \frac{2}{3}\sigma\right)a_6\right].$$

선형 이류에서는 이 값에 속도를 곱한 것이 시간 구간 전체의 정확한 플럭스다. 즉 PPM은 공간 재구성이 정해지면 그 포물선을 정확히 이류시키는 도식이다. 오일러 방정식 같은 계에서는 특성곡선법처럼 특성 변수로 분해해 파마다 따로 추적한 뒤 다시 합친다. 원 논문은 라그랑주 좌표에서 풀고 고정 격자로 재매핑하는 PPMLR과, 고정 격자에서 직접 푸는 PPMDE 두 판을 모두 제시했다.

4. 원 논문의 부속 장치들

교과서판 PPM은 위의 재구성만 보여 주지만, 1984년 논문에는 강한 충격을 다루기 위한 장치가 더 있다.

• 접촉 불연속 가파르게 하기(contact steepening). 접촉 불연속은 충격파와 달리 스스로 가파르게 유지되는 메커니즘이 없어서 수치 확산에 계속 퍼진다. PPM은 밀도 프로파일이 접촉 불연속처럼 보이면 경계값을 이웃 셀의 선형 외삽으로 바꿔 계단을 날카롭게 만든다.
• 평탄화(flattening). 강한 충격파 바로 뒤에서는 포물선 재구성 자체가 진동을 만든다. 압력 점프와 속도 수렴(압축) 정도로 충격 강도를 감지해, 강한 충격 주변에서는 포물선을 1차(상수)쪽으로 섞는다.
• 소량의 인공점성. 다차원 강충격 문제에서 남는 잡음을 잡기 위해 인공점성을 약간 더한다. 주 안정화 장치가 아니라 마무리용이다.

이 장치들은 대부분 경험적 임계값을 갖고 있어서, 코드마다 세부 구현이 조금씩 다르다. “PPM을 쓴다”는 말이 코드 사이에서 완전히 같은 도식을 뜻하지는 않는다는 점은 알아 둘 만하다.¹

5. 검산 — MUSCL과 맞붙이기

주기 경계 선형 이류 $a_t + a_x = 0$ 을 CFL 0.5로 한 주기 돌린 뒤 셀 평균의 $L_1$ 오차를 쟀다. MUSCL도 같은 방식(한 스텝 정확 이류 플럭스)으로 구현했다.

초기 조건	$N$	MUSCL minmod	MUSCL MC	PPM(제약 포함)	PPM(무제약 4차 경계)
사인파	50	1.62e-2	3.48e-3	2.46e-3	4.24e-5
사인파	100	4.59e-3	7.55e-4	4.36e-4	5.20e-6
사인파	200	1.25e-3	1.46e-4	8.27e-5	6.47e-7
사각파	100	4.93e-2	2.86e-2	1.94e-2	2.61e-2

세 가지가 보인다.

• 무제약 PPM은 정확히 3차다. 오차 비가 8.15 → 8.04 → 8.0(400칸에서 8.08e-8)으로 $2^3$ 에 붙는다.
• 제약을 걸면 사인파에서 차수가 약 2.4로 떨어진다. 50→100→200→400칸 오차 비가 5.6·5.3·5.5다. 사인파 꼭짓점에서 극값 평탄화가 1차로 깎기 때문이다. 그래도 같은 격자에서 MC 제한 MUSCL보다 오차가 1.4~2배 작다.
• 사각파에서 제약은 필수다. 무제약 PPM은 6.5% 오버슈트(최대 1.0647, 최소 −0.0647)를 남겼고, 제약 PPM은 $[0, 1]$ 을 정확히 지키면서 $L_1$ 오차도 가장 작았다.

6. 극값 문제와 후속 연구

두 번째 관찰이 PPM의 알려진 약점이다. 단조 제약은 “셀 평균이 극값이면 평탄화”라는 규칙 때문에 매끄러운 극값과 진짜 불연속을 구별하지 못한다. 난류 시뮬레이션처럼 매끄러운 극값이 사방에 있는 문제에서는 이게 상당한 수치 소산이 된다.

콜렐라·세코라(2008)는 극값 근처에서 2차 도함수 추정치 여러 개의 부호와 크기를 비교해, 매끄러운 극값이면 포물선을 살리고 불연속이면 제한하는 극값 보존 제한자를 제안했다. 이후 맥코커데일·콜렐라(2011)의 4차 유한체적법으로 이어졌다. 같은 문제의 다른 해법이 ENO 도식·WENO 도식이다 — 고정 스텐실에 제약을 거는 대신 매끄러운 스텐실을 골라 쓰는 접근으로, 극값에서 차수를 덜 잃는다.

7. 쓰이는 곳

PPM은 특히 천체물리 기체역학의 표준 재구성으로 자리 잡았다. FLASH의 분할 유체 솔버, Enzo, VH-1, CASTRO가 PPM을 주력 또는 기본 선택지로 쓰고, Athena++도 재구성 차수 옵션으로 PPM을 제공한다.² 초신성 폭발, 레일리-테일러 불안정, 강충격 상호작용 같은 문제가 주 무대다. 대기 과학에서는 린·루드 계열의 유한체적 역학 코어가 PPM형 재구성을 이류에 쓴다.

2D 시험 문제로 유명한 이중 마하 반사(double Mach reflection)와 전방 계단 유동(forward-facing step)도 같은 해 우드워드·콜렐라의 비교 논문에서 나왔다. 지금도 새 충격 포착 도식이 나오면 이 두 문제로 PPM과 비교하는 게 관행이다.

8. 관련 문서

• MUSCL · 고두노프 도식 · 리만 솔버 · 유한체적법
• 전변분 감소 · ENO 도식 · WENO 도식 · 인공점성
• 수치 소산과 분산 · CFL 조건 · 특성곡선법 · 적응 격자 세분화
• 압축성 유동 · 오일러 방정식 (비점성 유동) · 충격파 · 초신성

9. Footnotes

1. 같은 이름 아래 구현이 갈라지는 건 이 바닥의 고질병이다. 코드 비교 논문이 “PPM 대 PPM”에서 차이가 나오면, 대개 평탄화 임계값이나 특성 추적 여부가 범인이다. ↩

2. 우드워드는 이후로도 PPM 기반 코드로 별 내부 대류와 난류 혼합 대규모 계산을 계속했다. 한 사람이 40년 가까이 같은 도식을 갈아넣은, 보기 드문 장기 연재다. ↩