Skip to content
cfd-lab:~/ko/posts/2026-06-30-van-der-waals…● online
NOTE #088DAY TUE 유체역학DATE 2026.06.30UPDATED 2026.10.19READ 4 min read#Fluid-Mechanics#Equation-of-State#Van-der-Waals#Compressibility#Supercritical

이상기체가 거짓말을 시작하는 곳 — 압축인자 Z와 van der Waals 상태방정식

같은 속도에서 운동량 플럭스는 밀도 비의 제곱이다. 26%의 밀도 오차는 ρu²에서 약 60%가 된다.

100기압, 200 K 질소를 이상기체로 넣으면 밀도가 실제보다 4분의 1 가까이 빗나간다. 같은 압력, 같은 온도다. 분자는 서로를 당기고, 자기 부피도 차지한다.

그 어긋남은 압축인자 ZZ 하나다. van der Waals는 상수 두 개로 기체를 흉내 내고, 임계점 아래에서는 S자 루프가 액체와 기체를 가른다.

환산 온도 슬라이더를 1 아래로 내리고 Maxwell 체크를 켜면, 루프가 자라면서 면적이 같은 수평 포화선이 생긴다.

saturation pressure P_r,sat = 0.647
Z(sat. vapor) = 0.633  ·   Z(sat. liquid) = 0.163

0.85 근처에서 포화 액체의 ZZ는 0.1 아래로 떨어진다. 슬라이더를 1에 붙이면 액체와 기체의 ZZ가 한 점으로 모인다. 구분이 사라지는 자리다.

100기압의 질소는 이상기체보다 무겁다#

Pv=RTPv = RT는 분자를 점으로 보고 사이 힘을 버린다. 희박하면 맞다. 분자가 멀리 떨어져 있다.

압력을 올리면 분자 사이가 좁아진다. 약한 인력(van der Waals 힘)이 부피를 줄인다. 동시에 분자 자신의 부피가 남는다.

두 효과는 반대다. 인력은 부피를 줄이고, 유한 부피는 부피를 키운다. 어느 쪽이 이기는지는 온도와 압력이 정한다. 그래서 실제 기체가 이상기체보다 무거운지 가벼운지는 조건마다 뒤집힌다. 그 밀도 순위가 뒤집히면, 바닥의 가벼운 유체가 대류를 켜는 문턱도 같이 움직인다.

압축인자 Z — 어긋남을 하나의 숫자로#

어긋남을 한 수로 묶으면 압축인자다.

Z=PvRT=vactualvidealZ = \frac{Pv}{RT} = \frac{v_\text{actual}}{v_\text{ideal}}

PP는 압력, vv는 몰부피, RR은 기체상수, TT는 온도다. Z=1Z=1이면 이상기체다. 1에서 멀수록 더 많이 어긋난다.

Z<1Z<1은 인력 우세다. 실제 부피가 이상값보다 작다. Z>1Z>1은 유한 부피 우세다. 분자가 밀어 부피가 커진다. 질소 예의 Z≈0.79Z \approx 0.79는 실제 밀도가 이상값의 1/0.79≈1.271/0.79 \approx 1.27배라는 뜻이다.

대응 상태 원리: 모든 기체가 닮는 좌표#

기체마다 ZZ 곡선은 다르다. 압력과 온도를 각자 임계값으로 나누면 곡선이 거의 포개진다.

PR=PPc,TR=TTcP_R = \frac{P}{P_c}, \qquad T_R = \frac{T}{T_c}

PcP_c, TcT_c는 임계 압력·온도다. 대응 상태 원리(principle of corresponding states)는 환산 압력 PRP_R과 환산 온도 TRT_R이 같으면 ZZ도 거의 같다는 말이다.

질소든 메탄이든 이산화탄소든, TR=1.5T_R=1.5, PR=2P_R=2에서는 비슷한 ZZ다. 분자의 이름이 아니라 임계점에서의 상대 위치가 행동을 정한다. 그래서 차트 한 장으로 많은 기체를 덮는다.

코드의 질소는 그 예시와 다른 칸에 있다. T=200T = 200, Tc=126.2T_c = 126.2라 TR=200/126.2≈1.585T_R = 200/126.2 \approx 1.585다. P=107P = 10^7, Pc=3.39×106P_c = 3.39\times 10^6이라 PR=10/3.39≈2.95P_R = 10/3.39 \approx 2.95다. 출력 Z=0.7905Z = 0.7905를 환산 부피로 되돌리면 vR=Z (T/Tc) (Pc/P) (8/3)≈1.13v_R = Z\,(T/T_c)\,(P_c/P)\,(8/3) \approx 1.13이다. 이 점을 환산 등온선에 다시 넣으면 PRP_R이 2.95로 돌아온다. 인쇄된 ZZ는 그 초임계 등온선 위의 점이다.

van der Waals: 두 상수로 분자를 흉내내다#

1873년 van der Waals는 이상기체에 보정 두 개를 넣었다.

(P+av2)(v−b)=RT\left(P + \frac{a}{v^2}\right)(v - b) = RT

a/v2a/v^2는 인력이 깎는 압력이다. 부피가 작을수록 크다. bb는 분자가 차지하는 부피라 자유 공간은 v−bv-b다.

두 상수는 임계점에서 정한다. 거기에서 PP–vv 등온선은 변곡점이고, 기울기와 곡률이 동시에 0이다.

(∂P∂v)Tc=0,(∂2P∂v2)Tc=0\left(\frac{\partial P}{\partial v}\right)_{T_c} = 0, \qquad \left(\frac{\partial^2 P}{\partial v^2}\right)_{T_c} = 0

풀면 상수와 임계 압축인자가 나온다.

a=27R2Tc264Pc,b=RTc8Pc,Zc=PcvcRTc=38a = \frac{27 R^2 T_c^2}{64 P_c}, \qquad b = \frac{R T_c}{8 P_c}, \qquad Z_c = \frac{P_c v_c}{R T_c} = \frac{3}{8}

Zc=3/8=0.375Z_c = 3/8 = 0.375는 기체 종류와 무관하다. 실제 기체의 ZcZ_c는 0.27~0.29다. 0.375는 0.28 근처에 비해 약 3분의 1 크다. 정량적으로는 과대하고, 루프와 임계점의 그림은 맞다.

환산 변수로 쓰면 aa, bb가 사라진다.

PR=8TR3vR−1−3vR2P_R = \frac{8 T_R}{3 v_R - 1} - \frac{3}{v_R^2}

이 한 줄이 대응 상태 원리의 식이다. 임계점은 (vR,PR,TR)=(1,1,1)(v_R, P_R, T_R) = (1, 1, 1)에 고정된다. 분모 3vR−13v_R - 1이 0인 곳은 vR=1/3v_R = 1/3이다. 그 왼쪽은 분자 부피 bb보다 작아 등온선을 그리지 않는다.

임계점 아래의 S자 곡선과 Maxwell 작도#

TR>1T_R > 1이면 등온선은 단조롭게 감소한다. 압력을 올리면 부피가 줄어든다.

TR<1T_R < 1이면 S자 루프가 생긴다. 가운데는 (∂P/∂v)T>0(\partial P/\partial v)_T > 0이다. 압력을 올렸는데 부피가 커진다는 뜻이라 역학적으로 불안정하다. 실제 물질은 그 구간을 건너뛴다.

대신 수평선 하나가 액체와 기체를 잇는다. 높이가 포화 압력이다. Maxwell 등면적 작도(equal-area rule)는 루프에서 잘린 위·아래 면적이 같아지는 높이를 고른다. 액체와 기체의 Gibbs 자유에너지가 같다는 조건과 동치다. 두 상이 같이 있을 수 있는 압력은 그 하나다.

Python으로 Z를 푼다#

van der Waals를 ZZ의 다항식으로 쓰면 3차다. 임계 상수만 있으면 PP, TT에서 ZZ가 나온다.

import numpy as np
 
R = 8.314  # J/(mol·K)
 
def vdw_constants(Tc, Pc):
    """임계 상수로부터 van der Waals a, b."""
    a = 27.0 * R**2 * Tc**2 / (64.0 * Pc)
    b = R * Tc / (8.0 * Pc)
    return a, b
 
def vdw_compressibility(P, T, Tc, Pc):
    """Z^3 - (1+B)Z^2 + A Z - A B = 0 의 물리적 근."""
    a, b = vdw_constants(Tc, Pc)
    A = a * P / (R * T)**2          # 인력 항의 무차원 크기
    B = b * P / (R * T)            # 유한 부피 항의 무차원 크기
    coeffs = [1.0, -(1.0 + B), A, -A * B]
    roots = np.roots(coeffs)
    real = roots[np.abs(roots.imag) < 1e-8].real
    real = real[real > B]          # v > b  ⇒  Z > B 만 물리적
    return real.max(), real.min(), real.size
 
# 질소: Tc = 126.2 K, Pc = 3.39 MPa, M = 28 g/mol
Tc, Pc, M = 126.2, 3.39e6, 0.028
P, T = 1.0e7, 200.0                # 100 bar, 200 K
Z_gas, Z_liq, n = vdw_compressibility(P, T, Tc, Pc)
 
rho_ideal = P * M / (R * T)        # 이상기체 밀도
rho_real = rho_ideal / Z_gas       # 실제 밀도
print(f"근 개수 = {n},  Z = {Z_gas:.4f}")
print(f"이상기체 밀도 = {rho_ideal:6.1f} kg/m^3")
print(f"실제   밀도 = {rho_real:6.1f} kg/m^3  ({100*(rho_real/rho_ideal-1):+.0f}%)")

실행 결과는 세 줄이다.

근 개수 = 1,  Z = 0.7905
이상기체 밀도 =  168.4 kg/m^3
실제   밀도 =  213.0 kg/m^3  (+27%)

213.0/168.4=1.265213.0/168.4 = 1.265다. 증가분은 26.5%이고, 출력의 +27%+27\%는 그 값을 정수로 올린 것이다. TR≈1.59T_R \approx 1.59라 초임계다. 실근이 하나라 액체·기체 구분이 없다. 그래도 이상기체 밀도는 실제보다 약 4분의 1 가볍다.

초임계 인젝션과 CFD — 왜 상태방정식이 격자를 흔드는가#

액체 로켓 연소실은 보통 추진제 임계 압력보다 위에서 돈다. LOX(액체 산소)나 극저온 질소가 그 챔버로 들어가면 익숙한 물방울이 안 생긴다.

표면장력이 사라진다. 임계점을 넘으면 액체와 기체의 경계가 무너진다. 제트는 방울로 쪼개지지 않고, 밀도 구배가 빗살처럼 풀리며 주변과 섞인다. 초임계(또는 천이임계) 인젝션이다.

상태방정식이 격자를 흔드는 이유는 밀도가 압력·온도와 비선형으로 묶여서다. 이상기체를 쓰면 제트 밀도가 수십 % 틀리고, 운동량 플럭스 ρu2\rho u^2와 혼합 길이가 같이 밀린다. 실제 코드는 cubic EOS(Peng–Robinson, SRK)나 NASG(Noble–Abel Stiffened Gas)를 쓴다. ZZ는 그 모델이 칸마다 푸는 양이다.

밀도 오차를 운동량에도 27%라고 적으면 실패한다#

같은 속도에서 ρu2\rho u^2는 밀도 비의 제곱이다. 1.2652≈1.601.265^2 \approx 1.60이라, 26.5%의 밀도 오차는 운동량 플럭스에서 약 60%가 된다. 밀도가 27% 틀렸으니 플럭스도 27%라고 적으면 제곱을 빠뜨린 것이다.

그 압력을 힘에 그대로 넣으면, Maxwell이 고른 수평선 위의 계면은 가만히 있을까. 압력을 직접 넣었더니 정지한 계면이 떨린 비이상 격자볼츠만은 그 힘이 자유에너지로 바뀌어야 잠잠한지, 아직 열려 있다.

도움이 됐다면 공유해주세요.