안녕하세요! CFD(Computational Fluid Dynamics, 전산 유체 역학)를 공부하다 보면 반드시 마주치게 되는 알고리즘 중 하나가 바로 SIMPLE (Semi-Implicit Method for Pressure Linked Equations) 알고리즘입니다. 특히 불압축성 유동(Incompressible Flow) 문제를 풀 때 핵심적인 역할을 하죠.
이번 포스팅에서는 제공된 동영상 강의 내용을 바탕으로 SIMPLE 알고리즘이 왜 필요하며, 어떤 원리로 작동하는지 차근차근 알아보겠습니다.
왜 불압축성 Navier-Stokes 방정식은 풀기 어려울까?
CFD의 목표는 유체의 움직임을 지배하는 Navier-Stokes 방정식을 수치적으로 푸는 것입니다. 불압축성 유동의 경우, 주요 방정식은 다음과 같습니다.
- 연속 방정식 (Continuity Equation): 질량 보존
∇ ⋅ U = 0여기서∇ ⋅는 발산(divergence) 연산자이고,U는 속도 벡터입니다. 이 식은 유체가 흘러 들어오고 나가는 양이 같다는, 즉 밀도가 변하지 않는다는 것을 의미합니다. - 운동량 방정식 (Momentum Equation): 뉴턴 제2법칙 (F=ma)
(동영상에서는 정상 상태(steady-state) 및 외력이 없는 경우를 가정)U ⋅ ∇U - ∇ ⋅ (ν∇U) = -∇pU ⋅ ∇U: 대류항 (Convection term) – 유체 자체가 운동량을 가지고 이동하는 효과∇ ⋅ (ν∇U): 확산항 (Diffusion term) – 점성(viscosity)에 의한 운동량 확산 효과 (ν: 동점성 계수)-∇p: 압력 구배항 (Pressure gradient term) – 압력 차이에 의해 발생하는 힘 (p: Kinematic Pressure, 즉 압력/밀도)
문제점:
총 4개의 방정식(연속 방정식 1개, 운동량 방정식 3개)과 4개의 미지수(속도 성분 Ux, Uy, Uz, 압력 p)가 있어 풀 수 있을 것 같지만, 몇 가지 어려움이 있습니다.
- 압력 방정식의 부재: 운동량 방정식은 속도(U)에 대한 식이지, 압력(p)을 직접 구하는 식이 아닙니다.
- 연속 방정식의 역할: 연속 방정식(
∇ ⋅ U = 0)은 미지수를 푸는 방정식이라기보다는, 계산된 속도장(U)이 반드시 만족해야 하는 제약 조건(constraint)으로 작용합니다. 즉, 운동량 방정식으로 구한 속도가 질량 보존 법칙을 만족하도록 만들어야 합니다. 이것이 압력-속도 연계(Pressure-Velocity Coupling) 문제의 핵심입니다. - 비선형성: 운동량 방정식의 대류항(
U ⋅ ∇U)은 속도 자신이 포함된 비선형 항이라 해석적으로 풀기 어렵습니다. - 상태 방정식의 한계: 이상 기체 상태 방정식(
p = ρRT) 같은 것을 압력 계산에 바로 쓰기 어렵습니다. 불압축성 유동에서는 밀도(ρ)가 일정하고, 등온 유동이라면 온도(T)도 일정하여 압력을 직접 구할 수 없습니다.
SIMPLE 알고리즘: 압력-속도 연계 문제 해결하기
SIMPLE 알고리즘은 이러한 압력-속도 연계 문제를 해결하기 위한 반복적(iterative) 방법입니다. 핵심 아이디어는 다음과 같습니다.
- 운동량 방정식과 연속 방정식을 조합하여 압력(또는 압력 보정치)에 대한 방정식을 유도합니다.
- 계산된 속도장이 연속 방정식을 만족하도록 속도 보정 단계를 포함합니다.
이 과정은 예측(Predictor) 단계와 보정(Corrector) 단계로 나눌 수 있습니다.
방정식 유도 과정 (간략화)
1. 운동량 방정식의 행렬 형태:
먼저, 공간에 대해 이산화된(discretized) 운동량 방정식을 행렬 형태로 표현합니다. (예: 유한 체적법 사용)
MU = -∇p (식 5)U: 각 격자 셀 중심에서의 미지 속도 벡터∇p: 각 셀에서의 압력 구배 벡터 (이산화된 형태)M: 계수 행렬 (Coefficient Matrix). 이산화 과정에서 생성되며, 격자 정보, 유체 물성치, 이전 스텝의 속도값 등을 포함하는 알려진(known) 값으로 취급됩니다.
2. 행렬 분해 (Matrix Decomposition):
SIMPLE 알고리즘의 핵심 단계 중 하나는 계수 행렬 M을 대각 성분(diagonal component) A와 비대각 성분(off-diagonal component)의 기여도 H로 분리하는 것입니다.
AU - H = -∇p (식 9)A: 행렬M의 대각 성분만을 포함하는 대각 행렬 (Diagonal Matrix). 역행렬A⁻¹을 매우 쉽게 계산할 수 있습니다. (A⁻¹의 대각 성분은A의 대각 성분의 역수).H: 행렬M의 비대각 성분과 해당 이웃 셀 속도의 곱으로 이루어진 항. 현재 계산 단계에서는 이전 반복 스텝(previous iteration)의 속도값을 사용하여 계산하므로,H역시 알려진(known) 값으로 취급됩니다.
3. 압력 푸아송 방정식 (Pressure Poisson Equation) 유도:
이제 압력을 구하기 위한 방정식을 유도합니다.
- 분해된 운동량 방정식을 속도
U에 대해 정리합니다:U = A⁻¹H - A⁻¹∇p (식 16) - 이 속도
U에 대한 표현식을 연속 방정식∇ ⋅ U = 0에 대입합니다:∇ ⋅ (A⁻¹H - A⁻¹∇p) = 0 (식 17) - 압력 항을 좌변으로, 나머지 항을 우변으로 정리하면 압력에 대한 푸아송 방정식 형태를 얻습니다:
∇ ⋅ (A⁻¹∇p) = ∇ ⋅ (A⁻¹H) (Boxed Equation / 식 20)
이것이 바로 SIMPLE 알고리즘에서 압력(또는 압력 보정치)을 계산하는 데 사용되는 핵심 방정식입니다.
SIMPLE 알고리즘의 반복 계산 절차
이제 유도된 방정식을 사용하여 SIMPLE 알고리즘의 반복 계산 과정을 살펴봅니다.
- 초기 추측: 압력장(p)과 속도장(U)의 초기값을 추측합니다. (예: 이전 시간 스텝의 값 또는 0)
- 운동량 예측 (Momentum Predictor): 이전 스텝의 압력장(
p_old)을 사용하여 이산화된 운동량 방정식(MU* = -∇p_old)을 풀어 예측 속도장(predicted velocity field, U*)을 계산합니다.- 주의: 이 예측 속도장
U*는 일반적으로 연속 방정식(∇ ⋅ U* ≠ 0)을 만족하지 않습니다.
- 주의: 이 예측 속도장
- 압력 보정 (Pressure Correction): 예측 속도장
U*(정확히는U*로 계산된H*항)을 사용하여 압력 푸아송 방정식 (∇ ⋅ (A⁻¹∇p) = ∇ ⋅ (A⁻¹H*))을 풀어 새로운 압력장(p)을 계산합니다. (실제로는 압력 보정치p'를 계산하는 경우가 많지만, 개념적으로는 새로운 압력을 구하는 것입니다.) - 속도 보정 (Velocity Corrector): 새로 계산된 압력장
p를 사용하여 속도 보정식 (U = A⁻¹H* - A⁻¹∇p)으로 보정된 속도장(corrected velocity field, U)을 계산합니다. 이 보정된 속도장U는 연속 방정식(∇ ⋅ U = 0)을 만족하도록 만들어집니다. - (선택) 다른 스칼라 방정식 풀이: 만약 에너지, 난류 모델(k-ε 등), 화학종 수송 방정식 등이 있다면, 이 단계에서 해당 스칼라 변수들을 계산합니다. 계산된 스칼라 값들은 다음 반복 스텝에서 유체 물성치(점성 등) 업데이트에 사용될 수 있습니다.
- 수렴 확인 및 반복: 계산된 해가 수렴 조건을 만족하는지 확인합니다. (예: 잔차(residual)가 충분히 작은지). 만족하지 않으면, 업데이트된 압력(p)과 속도(U)를 가지고 단계 2부터 다시 반복합니다.
이 예측-보정 과정을 해가 수렴할 때까지 반복하는 것이 SIMPLE 알고리즘의 핵심입니다.
추가 방정식 통합 (에너지, 난류 등)
에너지 방정식, 난류 모델 관련 방정식(k, ε, ω 등), 화학종 수송 방정식 등 다른 스칼라 수송 방정식들은 일반적으로 위 SIMPLE 루프 내의 속도 보정 단계(4단계) 이후에 풀게 됩니다. 이 방정식들도 운동량 방정식과 유사하게 M_φ φ = S_φ 형태의 행렬 방정식으로 이산화하여 풀 수 있습니다 (φ는 해당 스칼라 변수, S_φ는 소스항).
이렇게 계산된 스칼라 값(온도, 난류 점성 등)은 다음 반복 계산에서 운동량 방정식의 계수 행렬 M과 소스항 H를 업데이트하는 데 사용되어 상호작용을 반영합니다.
마무리 생각: 압력 기반 알고리즘
SIMPLE 알고리즘은 압력에 대한 푸아송 방정식을 푸는 방식 때문에 압력 기반(Pressure-Based) 알고리즘이라고 불립니다. 불압축성 유동 해석에 널리 사용되지만, 비등온(non-isothermal) 유동의 경우, 상태 방정식(예: 이상 기체 법칙)을 사용하여 압력과 온도 계산 후 밀도를 업데이트합니다.
참고로, 주로 압축성 유동(Compressible Flow)에 사용되는 밀도 기반(Density-Based) 알고리즘도 있습니다. 이 방식은 연속 방정식에서 직접 밀도를 계산하고, 상태 방정식을 사용하여 밀도와 온도(에너지)로부터 압력을 계산하는, 반대 방향의 접근법을 사용합니다.
유용한 참고 자료
더 깊이 있는 내용을 원하시면 다음 자료들을 참고해보세요.
- H. Jasak, “Error Estimation and Analysis for the Finite Volume Method with Application to Fluid Flows”, PhD Thesis, Imperial College London, 1996: 유한 체적법 및 CFD 알고리즘에 대한 고전적인 자료입니다.
- OpenFOAM Wiki (OpenFOAM guide/The SIMPLE algorithm in OpenFOAM): 오픈 소스 CFD 코드인 OpenFOAM에서 SIMPLE 알고리즘이 어떻게 구현되어 있는지 볼 수 있습니다. (동영상에서 언급된 링크를 찾아보세요)
이 포스팅이 SIMPLE 알고리즘을 이해하는 데 도움이 되었기를 바랍니다. CFD는 깊이 있는 분야이므로 꾸준히 학습하는 것이 중요합니다. 궁금한 점이나 더 알고 싶은 내용이 있다면 언제든지 댓글로 남겨주세요!
