계면을 어떻게 다룰 것인가

압축성 다상유동 시뮬레이션에서 가장 근본적인 질문 중 하나는 "두 유체 사이의 **계면(interface)**을 수치적으로 어떻게 표현할 것인가"입니다.

크게 세 가지 계열의 방법이 존재합니다:

1. Volume of Fluid (VOF)
2. Level Set
3. Diffuse Interface (확산 계면)

1. Volume of Fluid (VOF)

원리

각 셀에서 특정 유체가 차지하는 체적 분율(volume fraction) $\alpha$를 추적합니다.

$$\frac{\partial \alpha}{\partial t} + \mathbf{u} \cdot \nabla \alpha = 0$$

$\alpha = 1$이면 해당 셀이 유체 A로 가득 차 있고, $\alpha = 0$이면 유체 B, $0 < \alpha < 1$이면 계면이 셀 내부를 지나갑니다.

계면 재구성 (Interface Reconstruction)

VOF의 핵심은 $\alpha$ 값으로부터 셀 내부의 계면 위치를 재구성하는 것입니다:

• SLIC (Simple Line Interface Calculation): 계면을 좌표축에 평행하게 근사
• PLIC (Piecewise Linear Interface Calculation): 계면을 임의 기울기의 직선으로 근사

PLIC에서 계면의 법선 벡터 $\mathbf{n}$은 $\alpha$ 필드의 기울기로 추정합니다:

$$\mathbf{n} = \frac{\nabla \alpha}{|\nabla \alpha|}$$

장단점

• 장점: 질량이 정확히 보존됨 (보존형)
• 단점: 계면 재구성이 복잡하고 비용이 큼. 3D 확장이 어려움. 곡률 계산 부정확

2. Level Set Method

원리

계면을 부호 거리 함수(signed distance function) $\phi$의 영점 등고선으로 정의합니다.

$$\frac{\partial \phi}{\partial t} + \mathbf{u} \cdot \nabla \phi = 0$$

• $\phi > 0$: 유체 A
• $\phi < 0$: 유체 B
• $\phi = 0$: 계면

Re-initialization

이류(advection) 과정에서 $\phi$는 부호 거리 함수의 성질($|\nabla\phi| = 1$)을 잃어갑니다. 이를 복원하기 위해 re-initialization 방정식을 풉니다:

$$\frac{\partial \phi}{\partial \tau} + \text{sign}(\phi_0)(|\nabla\phi| - 1) = 0$$

여기서 $\tau$는 가상 시간(pseudo-time)입니다.

곡률과 표면장력

Level Set의 큰 장점은 기하학적 양의 계산이 자연스럽다는 것입니다:

$$\kappa = -\nabla \cdot \left(\frac{\nabla\phi}{|\nabla\phi|}\right)$$

표면장력은 CSF(Continuum Surface Force) 모델로 체적력으로 변환됩니다:

$$\mathbf{f}_s = \sigma \kappa \delta(\phi) \nabla\phi$$

장단점

• 장점: 곡률/법선 계산 용이. 위상 변화(topology change) 자연스럽게 처리
• 단점: 질량 보존 실패. Re-initialization 과정에서 계면 위치가 이동 가능

CLSVOF: 두 방법의 결합

VOF의 질량 보존 + Level Set의 기하학적 장점을 합친 Coupled Level Set and VOF (CLSVOF) 방법이 많이 사용됩니다.

3. Diffuse Interface Method (확산 계면)

철학의 전환

앞의 두 방법이 계면을 날카로운(sharp) 것으로 보는 반면, 확산 계면법은 계면이 본질적으로 유한한 두께를 갖는다고 봅니다.

이 접근에서는 별도의 계면 추적 방정식이 필요 없고, 체적 분율이 보존 법칙 시스템의 일부로 포함됩니다.

5-equation Model (Allaire et al., 2002; Kapila et al., 2001)

가장 널리 사용되는 reduced model입니다:

$$\frac{\partial (\alpha_1 \rho_1)}{\partial t} + \nabla \cdot (\alpha_1 \rho_1 \mathbf{u}) = 0$$ $$\frac{\partial (\alpha_2 \rho_2)}{\partial t} + \nabla \cdot (\alpha_2 \rho_2 \mathbf{u}) = 0$$ $$\frac{\partial (\rho \mathbf{u})}{\partial t} + \nabla \cdot (\rho \mathbf{u} \otimes \mathbf{u} + p\mathbf{I}) = 0$$ $$\frac{\partial E}{\partial t} + \nabla \cdot ((E + p)\mathbf{u}) = 0$$ $$\frac{\partial \alpha_1}{\partial t} + \mathbf{u} \cdot \nabla \alpha_1 = 0$$

마지막 방정식이 핵심입니다. $\alpha_1$의 이류 방정식이 **비보존형(non-conservative)**이라는 점. 이것이 바로 계면에서의 압력 진동을 방지하는 열쇠입니다.

왜 비보존형이 필요한가?

$\alpha_1$을 보존형으로 풀면:

$$\frac{\partial \alpha_1}{\partial t} + \nabla \cdot (\alpha_1 \mathbf{u}) = \alpha_1 \nabla \cdot \mathbf{u}$$

우변의 $\alpha_1 \nabla \cdot \mathbf{u}$ 항이 이산화 과정에서 적절히 처리되지 않으면 EOS mixing에 의한 비물리적 압력 진동이 발생합니다.

비보존형 이류는 $\alpha_1$이 계면에서 sharp하게 유지되면서도 압력/속도 평형을 만족시킵니다.

7-equation Model (Baer-Nunziato type)

각 상이 독립적인 압력과 속도를 가질 수 있는 완전 비평형 모델입니다:

$$\frac{\partial \alpha_k \rho_k}{\partial t} + \nabla \cdot (\alpha_k \rho_k \mathbf{u}_k) = 0, \quad k = 1, 2$$

이완(relaxation) 과정을 통해 압력/속도 평형으로 수렴시킵니다. 수학적으로 가장 완전하지만 계산 비용이 높습니다.

비교 요약

| 기법 | 질량 보존 | 계면 선명도 | 위상 변화 | 구현 복잡도 | 압축성 적합도 | |------|:---------:|:----------:|:---------:|:----------:|:------------:| | VOF (PLIC) | O | 높음 | 어려움 | 높음 | 보통 | | Level Set | X | 높음 | 자연스러움 | 보통 | 보통 | | CLSVOF | O | 높음 | 자연스러움 | 매우 높음 | 보통 | | 5-eq Diffuse | O | 보통 | 자연스러움 | 낮음 | 매우 좋음 | | 7-eq Diffuse | O | 보통 | 자연스러움 | 높음 | 매우 좋음 |

결론

압축성 다상유동에서는 **Diffuse Interface (특히 5-equation model)**이 가장 자연스럽고 robust한 선택입니다. 충격파-계면 상호작용, 캐비테이션, 수중 폭발 등의 문제에서 탁월한 성능을 보입니다.

다음 글에서는 이러한 기법들을 실제로 코드로 구현하는 방법을 다룹니다.