APPLICATION OF ENTROPY FIX TO COMPUTATIONS OF HYPERSONIC LAMINAR AND TURBULENT BOUNDARY LAYERS ON HIFIRE-1

J. Park; M. Jeong; S. Jee

doi:10.6112/kscfe.2025.30.3.136

Preview

Original Article

Journal of Computational Fluids Engineering. 30 September 2025. 136-150
https://doi.org/10.6112/kscfe.2025.30.3.136

APPLICATION OF ENTROPY FIX TO COMPUTATIONS OF HYPERSONIC LAMINAR AND TURBULENT BOUNDARY LAYERS ON HIFIRE-1

엔트로피 수정 기법을 이용한 HIFiRE-1 모델의 극초음속 층류 및 난류 경계층 수치해석

J. Park¹

M. Jeong¹

S. Jee¹^†

박 재현¹

정 민재¹

지 솔근¹^†

¹Department of Mechanical and Robotics Engineering, Gwangju Institute of Science and Technology

¹광주과학기술원 기계로봇공학과

^{*Corresponding Author}

License (open-access, http://creativecommons.org/licenses/by-nc/4.0):

This is an Open-Access article distributed under the terms of the Creative Commons Attribution Non-Commercial License (http://creativecommons.org/licenses/by-nc/4.0) which permits unrestricted noncommercial use, distribution, and reproduction in anymedium, provided the original work is properly cited.

ABSTRACT

This study presents a case study of hypersonic laminar and turbulent boundary layers on the HIFiRE-1 model, proposed as a benchmark in the 2024 High-Fidelity CFD Verification Workshop. Computations were initially performed with the Roe flux-differencing scheme, but non-physical carbuncle phenomena appeared near the cone nose. To improve numerical stability, an entropy fix was applied by modifying the eigenvalues comprising the dissipation term. The improved simulation captured key shock structures and showed good agreement with reference data including surface pressure and heat flux. The turbulent heat flux is up to four times higher than in the laminar one, highlighting the impact of the boundary layer state on the aerothermodynamics of a high-speed vehicle.

Keywords

Hypersonic boundary layer

laminar and turbulent flow

Entropy fix

HIFiRE-1

Computational fluid dynamics

키워드

극초음속 경계층 유동

층류와 난류 유동

엔트로피 수정

하이파이어-1

전산유체역학

MAIN

1. 서 론
2. 해석 방법(Methodology)
2.1 지배방정식 및 난류 모델
2.2 수치 해석법 및 해석자
2.3 격자 및 경계 조건
3. 해석 결과(Results)
3.1 유동장(Flow Domain)
3.2 Residual history 및 격자 수렴성 검사
3.3 벽면 압력 및 열유속 분포
4. 결 론(Conclusion)

1. 서 론

최근 극초음속 항공기, 상업용 위성, 우주 재진입체 등 차세대 우주 수송 기술의 개발이 활발히 이루어지면서, 극초음속 유동에 대한 수치해석 기술의 중요성 또한 함께 부각되고 있다. 극초음속 유동(hypersonic flow)은 마하수 5 이상에서 발생하는 고압축성 유동으로, 비행체 주위에 강한 충격파가 형성되며, 난류 경계층에서는 공력 가열과 마찰 저항의 증가로 인해 비행체의 생존성(survivability)에 중대한 영향을 미친다. 이러한 유동 현상을 정밀하게 예측하는 것은 극초음속 비행체의 설계에 있어 필수적인 요소이며, 이를 위해 다양한 대표 비행체 모델들이 실험 및 수치해석 연구에 활용되었다. 대표적인 예로 NASA의 X-43A[1], DLR의 SHEFEX 시리즈[2], 그리고 미국-호주 공동 연구 프로그램인 HIFiRE(High Fidelity Flight Research Experiment) 시리즈[3] 등이 있다. 이들 비행체 모델은 다양한 기하 구조와 유동 조건에서 극초음속 유동의 특성을 포괄적으로 관찰할 수 있도록 설계되었다.

그 중 HIFiRE-1 모델은 원뿔체(cone), 실린더(cylinder), 그리고 플레어(flare)로 구성된 실험체로, 무딘 선두부(blunt nose)에서 형성되는 bow 충격파(shock wave), 원뿔체 위에서의 경사(oblique) 충격파, 실린더 구간의 팽창파(expansion wave), 그리고 플레어에서의 재압축(recompression) 충격파 등 다양한 유동 구조를 포착할 수 있도록 설계되었다(Fig. 1)[4]. 해당 모델은 고도 30km에서 Mach 수 6.58 및 7.16 조건으로 AFRL 연구소에 의해 비행 실험이 수행되었고[5], 이후 유사 조건에서 미국 CUBRC 사에 의해 지상 풍동 실험[4]도 이루어졌다.이러한 특성으로 인해 HIFiRE-1은 극초음속 유동 해석에서의 난류 모델 성능 평가 및 수치 해석자의 해석 정확도 비교를 위한 기준 모델로 널리 활용되고 있다[6,7].

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F1.jpg

Fig. 1.

Geometry of HIFiRE-1 consistings of cone, cylinder, and flare[4]

특히, 2024년에 개최된 High-Fidelity CFD Verification Workshop[8]에서는 수치해석 결과의 상호 검증을 위한 표준 모델 중 하나로 HIFiRE-1 모델이 선정되었으며, 층류 및 난류 경계층 해석이 주요 과제로 제시되었다 . 본 워크숍에는 미국 Sandia National Lab, 광주과학기술원 등 여러 기관이 참여하였고, 각 기관은 동일한 형상 조건과 격자 조건에서 HIFiRE-1 모델에 대한 수치해석을 수행하였다. 자유류 조건은 마하수 $M_{\infty} = 7.18$ , 받음각 $α = 2^{°}$ , 그리고 단위 레이놀즈 수는 $R e_{\infty} = 10.123 \times 10^{6} / m$ 로 설정되었다. 해석 결과 비교는 벽면 압력과 열유속 분포를 기준으로 수행되었으며, 난류 경계층 해석에는 RANS의 SA-negative 모델이 공통적으로 적용되었다. 수치 해석자 및 수치 기법의 선택은 각 기관의 자율에 맡겨졌으며, 기관별 수치 해석자와 수치 해석법의 정확성을 비교 및 검증하는 데 목적이 있다.

HIFiRE-1 기체의 형상과 워크숍에서 제시된 유동 조건으로 인해 기체 주위에는 강한 충격파가 형성되며, 이로 인한 수치적 불안정성을 억제하는 것이 해석의 핵심 과제이다. 충격파 근처에서는 유동 특성이 급변하면서 수치 진동이 발생하기 쉬운데, 이는 수치해의 정확성과 수렴성 모두에 부정적인 영향을 미친다. 이러한 문제를 해결하기 위해 흔히 활용되는 Roe flux 차분 기법[9]은, 셀 경계에서 정의된 평균 상태를 기반으로 고유값(eigenvalue) 분해를 수행하고, 이를 통해 특성파를 분리함으로써 불연속성을 안정적으로 처리할 수 있다. 이 방식은 급격한 유동 변화 구간에서도 수치적 발산을 억제하고, 강건한 해를 도출하는 데 유리하다. 특히, TVD(total variation diminishing) 조건을 만족하는 보정 기법을 병행하면 충격파나 경계층처럼 유동이 급변하는 영역에서도 높은 정확도와 안정적인 수렴을 동시에 확보할 수 있어, 충격파 지배 유동(shock-dominated flow) 해석에서 효과적인 수단으로 사용된다[10].

그러나 Roe flux 차분 기법이 효과적인 충격파 포착 능력을 갖추고 있음에도, 극초음속 조건에서는 비물리적인 해를 유발할 수 있는 한계를 지닌다. 특히 무딘 물체 주위에 형성되는 정상 충격파 영역에서는 일부 고유값이 0에 가까워지면서 수치적 안정성이 저하되고, 충격파를 경계로 엔트로피 조건을 만족하지 못하는 현상이 발생할 수 있다[11]. 이로 인해 충격파가 파단되는 carbuncle 현상이 나타나며, 그 결과 충격파 주변의 수치 진동이 하류로 갈수록 증폭되어 해의 수렴성을 저해하게 된다[12]. 이러한 문제를 완화하기 위한 방안이 엔트로피 수정 기법(entropy fix)으로, 고유값이 0에 가까운 영역에서 인위적으로 수치적 소산을 추가함으로써 Roe 기법의 수렴성과 계산 안정성을 향상시킬 수 있다[13].

본 연구는 2024년 High-Fidelity CFD Verification Workshop[8]을 통해 HIFiRE-1 형상을 대상으로 한 수치해석 사례 분석(case study)이며, 이를 통해 극초음속 층류 및 난류 경계층 유동 해석에 사용되는 수치 기법의 타당성을 검토하고자 하였다. 해석에는 Roe flux 차분 기법을 적용하였으며, 이 과정에서 발생하는 수치적 불안정성을 확인하고, 이를 완화하는 방안으로 엔트로피 수정 기법을 도입하였다. 해당 기법이 해석의 정확도와 계산 안정성에 미치는 영향을 정량적으로 평가하였으며, 결과는 워크숍에서 제공된 문헌 데이터[14]와 비교하였다. 또한, 층류와 난류 해석에서 벽면 압력 및 열 유속 분포의 차이를 비교하였다. 이를 통해 HIFiRE-1과 같은 대표적 극초음속 비행체에 대해, 층류 및 난류 경계층 해석을 위한 수치해석 방법의 타당성과 해석 결과에 나타나는 유동 특성을 고찰하였다.

2. 해석 방법(Methodology)

2.1 지배방정식 및 난류 모델

본 연구에서는 극초음속 층류 및 난류 경계층 해석을 위하여 정상(steady) 상태의 압축성 유동 Navier-Stokes 방정식을 이용하였으며, 그 형태는 식 (1)과 같다.

(1)

\frac{\partial (ρ u_{i})}{\partial x_{i}} = 0 \frac{\partial}{\partial x_{i}} (ρ u_{i} u_{j} + p δ_{i j}) = \frac{\partial}{\partial x_{i}} τ_{i j} \frac{\partial}{\partial x_{i}} (ρ u_{i} H) = \frac{\partial}{\partial x_{i}} (τ_{i j} u_{j} - q_{i})

식 (1)에 있어 $ρ$ 는 밀도, $u_{i}$ 는 속도, $p$ 는 압력을 나타내며 $δ_{i j}$ 는 Kronecker delta 함수를 의미한다. 또한, 에너지 보존 방정식에서 $H$ 는 단위 질량당 총 엔탈피(total enthalpy)로, $H = E + p / ρ$ 로 정의되며, $E$ 는 단위 질량당 총에너지(total energy)로서 $E = e + u_{k} u_{k} / 2$ 로 나타난다. 여기서 $e$ 는 단위 질량당 내부에너지(internal energy)를 의미한다. 전단 응력 텐서(shear stress tensor)는 $τ_{i j}$ 이며 식 (2)와 같이 유효 점성(effective viscosity, $μ_{e f f}$ )과 변형률 텐서 $S_{i j}$ 로 정의된다.

(2)

τ_{i j} = 2 μ_{e f f} S_{i j} = μ_{e f f} (\frac{\partial u_{i}}{\partial x_{j}} + \frac{\partial u_{j}}{\partial x_{i}}) μ_{e f f} = μ + μ_{t}

층류 해석에서 난류 점성 $μ_{t} = 0$ 이며, 기체 점성 $μ$ 는 Sutherland law에 따라 정의된다.

(3)

μ = μ_{0} {(\frac{T}{T_{0}})}^{1.5} \frac{T_{0} + T_{s}}{T + T_{s}}

Sutherland law에서 $μ_{0} = 1.717 \times 10^{- 5} k g / (m s), T_{0} = 273.15 K$ 그리고 $T_{S} = 111 K$ 로 주어진다.

난류 해석에서는 정상 RANS(Reynolds averaged Navier-Stokes) 해석을 수행하였으며, $μ_{t}$ 는 Spalart-Allmaras(SA) negative 모델을 사용하여 계산된다. SA-negative 모델은 SA 모델 변수가 $\hat{ν} ≧ 0$ 인 경우 standard SA 모델[15]을 적용하고, $\hat{ν} < 0$ 인 경우 수정된 수송 방정식[16]을 적용하여 난류 점성을 계산한다. SA-negative 모델의 수송 방정식은 난류 점성이 0보다 작아지는 비물리적인 해석 결과를 방지하며, 국부적으로 $\hat{ν} < 0$ 인 영역에서 적절히 $μ_{t} = 0$ 으로 보정하는 효과를 가진다.

식 (1)의 에너지 보존 방정식에 있어 $q_{i}$ 는 열 유속(heat flux)이며, $K_{e f f}$ 는 유효 열전도율(thermal conductivity)로 층류 및 난류 해석에 따라 아래와 같이 정의된다.

(4)

q_{i} = - K_{e f f} \frac{\partial T}{\partial x_{i}} K_{e f f} = K + K_{t} K = \frac{μ C_{p}}{\Pr} (\Pr = 0.72), K_{t} = \{\begin{cases} 0 & (l a m i n a r) \\ \frac{μ_{t}}{\Pr_{t}} & (t u r b u l e n t, \Pr_{t} = 0.9) \end{cases}

2.2 수치 해석법 및 해석자

본 연구에서는 NASA에서 개발한 정렬 격자 기반 압축성 유동 해석 코드인 CFL3D를 이용하여 수치 해석을 수행하였다[17]. 점성 항은 2차 중심 차분(central differencing) 기법을 적용하였으며, 대류 항은 2차 정확도의 Roe flux 차분 기법(Roe flux differencing scheme, FDS)을 통해 공간 차분화되었다. Roe FDS는 셀 경계에서 Euler flux를 근사하는 방식으로, 유동이 오른쪽 방향으로 하류 진행하는 경우 셀 경계의 좌측(L)과 우측(R) 셀의 정보를 이용하여 식 (5)와 같이 수치 플럭스를 계산한다[9]. 이 과정은 중앙 항(central term)과 소산 항(dissipation term)으로 구분된다. 소산 항의 경우 3차원 압축성 유동 Euler 방정식의 오른쪽 고유벡터 행렬(right eigenvector matrix, $R$ )과 고유값 대각 행렬(digonal matrix of eigenvalues, $Λ$ )로 구성된다.

(5)

F = \frac{1}{2} (F_{L} + F_{R}) + F_{d} F_{d} = - \frac{1}{2} R |Λ| R^{- 1} (Q_{R} - Q_{L}) F = [\begin{matrix} ρ u_{i} \\ ρ u_{i} u_{j} + p δ_{i j} \\ ρ u_{i} H \end{matrix}], Q = [\begin{matrix} ρ \\ ρ u_{i} \\ ρ E \end{matrix}]

$x$ -방향 유동에 대한 오른쪽 고유 행렬 $R$ 과 고유값 대각 행렬 $Λ$ 는 아래와 같이 정의된다.

(6)

R = [\begin{matrix} 1 & 0 & 0 & 0 & 1 \\ \tilde{u} - \tilde{a} & 0 & 0 & 0 & \tilde{u} + \tilde{a} \\ \tilde{v} & 1 & 0 & 0 & \tilde{v} \\ \tilde{w} & 0 & 1 & 0 & \tilde{w} \\ \tilde{H} - \tilde{u} \tilde{a} & \tilde{v} & \tilde{w} & 1 & \tilde{H} + \tilde{u} \tilde{a} \end{matrix}], Λ = d i a g (λ_{1}, λ_{2}, λ_{3}, λ_{4}, λ_{5}) λ_{1} = \tilde{u} - \tilde{a}, λ_{2} = λ_{3} = λ_{4} = \tilde{u}, λ_{5} = \tilde{u} + \tilde{a} \tilde{u} = \frac{\sqrt{ρ_{L}} u_{L} + \sqrt{ρ_{R}} u_{R}}{\sqrt{ρ_{L}} + \sqrt{ρ_{R}}} \tilde{H} = \frac{\sqrt{ρ_{L}} H_{L} + \sqrt{ρ_{R}} H_{R}}{\sqrt{ρ_{L}} + \sqrt{ρ_{R}}} \tilde{a} = \sqrt{(γ - 1) (\tilde{H} - \frac{1}{2} ({\tilde{u_{i}}}^{2}))}

여기서, $\tilde{u}, \tilde{v}, \tilde{w}$ 는 각각 Roe 평균 속도의 $x, y, z$ 성분을 의미하며, $\tilde{H}$ 는 Roe 평균 총 엔탈피, $\tilde{a}$ 는 Roe 평균 음속이다.

극초음속 유동에서 충격파를 경계로 유동장 내 물리량이 급격하고 불연속적으로 변화한다. Roe FDS는 인접한 셀 간의 상태 차이를 기반으로 수치 플럭스를 계산하므로, 충격파 전후의 유동 조건을 동시에 반영할 수 있다는 점에서 극초음속 유동 해석에 유리한 특성을 가진다[9]. 그러나 구성된 고유값 중 일부가 0에 근접할 경우, 소산항( $F_{d}$ )의 효과가 감소하여 수치적 안정성이 저하되는 문제가 발생한다. 또한, 충격파를 경계로 엔트로피가 감소하는 비물리적인 해(non-physical solution)가 도출되는 한계도 존재한다[18].

이러한 고 마하수 유동에서 Roe FDS의 충격파 포착 한계를 보완하기 위해 엔트로피 수정 기법(entropy fix)이 도입되었다. 이 기법은 고유값 $λ_{n}$ 이 일정 임계값 $(h)$ 보다 작아지는 경우, 이를 보정하여 수치적 소산을 인위적으로 증가시키는 방식으로 수치적 불안정을 완화한다 . 본 연구에서는 CFL3D에 내장된 엔트로피 수정 기법[19]으로 Yee[20]가 제안한 보정 함수를 사용하였으며, 그 형태는 식 (7)에 제시되어 있다. 구체적으로는, 식 (5)에 포함된 Roe FDS의 소산항에서 고유치 대각 행렬( $Λ$ )를 보정된 형태인 $\hat{Λ}$ 로 대체하여 적용한다. 엔트로피 수정 기법에 사용되는 임계값은 계산 안정성과 정확성 간의 균형을 고려하여 $ε_{λ}$ 를 $0.05 \leq ε_{λ} \leq 0.2$ 로 설정한다[21]. 본 연구에서는 선행 연구[22]를 참조하여 $ε_{λ} = 0.1$ 을 적용하여 계산을 수행하였다. Gnoffo[22]은 $0.05 \leq ε_{λ} \leq 0.3$ 에 대해 둔체(blunt body) 모델을 대상으로 검토하였으며, $ε_{λ} = 0.1$ 에서 수치적 소산을 최소화하면서 안정적인 계산을 수행할 수 있음을 확인하였다.

(7)

\hat{Λ} = d i a g ({\hat{λ}}_{1}, {\hat{λ}}_{2}, {\hat{λ}}_{3}, {\hat{λ}}_{4}, {\hat{λ}}_{5}) {\hat{λ}}_{n} = \{\begin{cases} λ_{n} & λ_{n} \geq h \\ \frac{1}{2} (\frac{λ_{n}^{2}}{h} + h) & λ_{n} < h \end{cases} h = ε_{λ} \max (λ_{n})

2.3 격자 및 경계 조건

Fig. 2는 HIFiRE-1 비행체 모델의 표면의 층류 및 난류 경계층 유동 해석을 위한 격자 구성과 경계 조건을 나타낸 것이다. 본 연구에서는 2024 High-Fidelity CFD Verification Workshop에서 제공된 기준 격자[8,23]를 활용하였으며, 이는 기관 간 해석자 및 수치해석 기법에 따른 결과의 정확도와 일관성을 비교·평가하기 위한 목적을 갖는다. 본 연구에서는 층류 및 난류 유동 해석을 개별적으로 수행하였다. 층류 유동은 선두부의 원뿔(cone) 형상만을 대상으로 수행하였고, 난류 유동 해석은 원뿔-실린더-플레어(cone-cylinder-flare)를 포함한 전체 모델 형상에 대해 수행되었다. 해석에는 서로 다른 해상도의 Medium, Fine, Finer 세 종류의 격자가 사용되었으며, 각 격자는 모든 방향에서 격자 수가 2 배씩 증가하도록 구성되었다. 이를 바탕으로 격자 해상도에 따른 경계층 유동 해석 결과를 비교하였다. 각 격자에 대한 상세 정보는 Table 1에 제시되어 있으며, $∆ y^{+}$ 는 무차원 벽 거리(non-dimensional wall distance)를 나타내고 식 (8)과 같이 정의된다.

(8)

∆ y^{+} = \frac{∆ y_{w a l l} u_{τ}}{ν}

여기서 $∆ y_{w a l l}$ 는 벽면 첫 번째 격자의 수직 거리, $u_{τ} = \sqrt{τ_{w} / ρ}$ 로 마찰 속도(friction velocity)를 나타내며, $τ_{w}$ 는 벽 전단응력(wall shear stress)을 나타낸다. $y = 1.4 m$ 에서 실린더 표면에서의 무차원 벽 수직 거리를 나타낸다. Table 1의 무차원 벽 거리 정보는 경계층 두께가 얇은 풍상(windward)측을 기준으로 평가된 값으로, 다른 영역에서도 경계층 내 충분한 격자 해상도가 확보됨을 의미한다.

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F2.jpg

Fig. 2.

Computational grid configuration corresponding to Table 1, including boundary conditions and target regions for laminar(cone only) and RANS(whole body: cone, cylinder and flare) simulations

Table 1.

Grid resolution of three grids Medium, Fine, and Finer, in laminar and RANS(turbulence) simulations[8], with non-dimensional wall distance $∆ y_{c o n e}^{+}$ at $x = 0.5 m$ (conical surface, windward) and $∆ y_{c y l i n d e r}^{+}$ at $x = 1.4 m$ (cylinder surface, windward)

Case	$N_{n x} \times N_{n y} \times N_{n θ}$	$N_{x} \times N_{y} \times N_{θ}$		Total	$∆ y_{c o n e}^{+}$	$∆ y_{c y l i n d e r}^{+}$
Medium	$32 \times 128 \times 16$	Laminar	$112 \times 128 \times 64$	$9.83 \times 10^{5}$	0.66	-
Medium	$32 \times 128 \times 16$	RANS	$240 \times 128 \times 64$	$2.03 \times 10^{6}$	0.76	0.31
Fine	$64 \times 256 \times 32$	Laminar	$224 \times 256 \times 128$	$7.86 \times 10^{6}$	0.32	-
Fine	$64 \times 256 \times 32$	RANS	$480 \times 256 \times 128$	$16.3 \times 10^{6}$	0.37	0.16
Finer	$128 \times 512 \times 64$	Laminar	$448 \times 512 \times 256$	$62.9 \times 10^{6}$	0.17	-
Finer	$128 \times 512 \times 64$	RANS	$960 \times 512 \times 256$	$130 \times 10^{6}$	0.19	0.079

해석에서 자유류 조건은 마하수 $M_{\infty} = 7.18$ , 받음각 $α = 2^{°}$ , 정온도(static temperature) $T_{\infty} = 231.9 K$ 로 설정되었으며, 단위 Reynolds 수는 $R e_{\infty} = 10.123 \times 10^{6} / m$ 이다. 이러한 조건은 입구(inlet) 경계에 비반사 경계 조건(non- reflecting boundary condition, NRBC)과 함께 적용되었다[24]. NRBC는 벽면이나 충격파에서 반사된 수치적 파동(numerical wave)이 유입 경계를 통해 계산 영역 내로 다시 들어오는 것을 억제하여, 해석의 안정성과 정확도를 향상시키는 데 기여한다. 벽면에서는 no-slip condition $(u_{i} = 0)$ 을 적용하였으며, 등온 조건으로 $T_{w} = 301.5 K$ 를 적용하여 자유류 온도와 $T_{w} / T_{\infty} = 1.3$ 의 비를 갖는다. $x - y$ 평면에 대하여 대칭(symmetric) 조건을 부여하였으며, 출구(outflow) 경계에서는 내부 셀의 값을 그대로 연장하는 방식으로 외삽(extrapolation)하였다. 난류 모델에 있어 난류 변수는 자유류에서 ${\hat{ν}}_{\infty} / ν = 4.0$ 그리고 벽면에서 ${\hat{ν}}_{w a l l} = 0$ 의 값을 가지며, $ν = μ / ρ$ 는 기체 동점성(kinematic viscosity)이다.

3. 해석 결과(Results)

3.1 유동장(Flow Domain)

Fig. 3은 엔트로피 수정 기법이 적용된 RANS 해석 결과로부터 얻어진 밀도 구배(density gradient) 기반의 수치적 Schlieren 이미지로, 자유류 마하수 $M_{\infty} = 7.18$ 조건에서 HIFiRE-1 기체 주위에 형성되는 주요 유동 구조를 가시화한 것이다. 먼저, 기체 선두부에서는 강한 bow 충격파가 형성된다. 이어서 원뿔체의 경사진 표면을 따라 경사 충격파가 발생하며, 이 두 충격파는 연속적으로 연결되어 있다. 이어지는 실린더 구간에서는 유동이 확장되면서 팽창파가 발생하고, 이에 따라 유속은 증가하며 압력과 온도는 감소한다. 마지막으로, 후방의 flare에서는 유동이 다시 압축되며, 재압축 충격파가 발생하고, 이로 인해 온도와 압력이 다시 상승하게 된다. 이렇듯 극초음속 유동 조건에서 발생하는 다양한 압축 및 팽창 현상의 위치와 상대적인 세기가 뚜렷하게 나타나며, 이는 해석에서의 수치적 안정성을 확보하기 위해 도입된 엔트로피 수정 기법의 효과를 보여주는 결과이다.

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F3.jpg

Fig. 3.

Numerical Schlieren image of flow structures around HIFiRE-1 at $M_{\infty} = 7.18$ and $α = 2^{°}$

Fig. 4는 현 RANS 해석에서, 엔트로피 수정법의 적용 여부에 따른 마하 수(Mach number) 등고선을 비교한 것이다. 동일한 수치 기법을 사용했음에도 불구하고, 엔트로피 수정 기법을 적용하지 않은 경우(좌측)에는 충격파 전면이 파단되는 carbuncle 현상[25]이 발생하여, 충격파가 매끄럽지 않고 비물리적으로 왜곡된 해가 나타났다. 반면, 엔트로피 수정 기법을 적용한 경우(우측)에는 충격파가 안정되고 일관된 형태로 유지되며, carbuncle 현상이 억제되어 물리적으로 타당한 해가 도출되었다. 이러한 경향은 층류 해석에서도 유사하게 관찰되었으며, 엔트로피 수정 기법의 적용이 해석 조건에 관계없이 수치적 안정성 확보에 기여함을 시사한다. 이는 Roe FDS 기법에서도 엔트로피 수정 기법이 충격파 처리의 안정성과 정확성 향상에 중요한 역할을 한다는 점을 뒷받침한다.

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F4.jpg

Fig. 4.

Comparison of Mach number isolines around the nose in RANS with or without the entropy fix method

3.2 Residual history 및 격자 수렴성 검사

질량 보존 방정식의 밀도 항을 기준으로 계산된 잔차(residual) 이력을 통해, 층류 및 RANS 해석에서 엔트로피 수정 기법의 적용 여부가 수렴성에 어떠한 영향을 미치는지 비교하였다(Fig. 5). 모든 해석은 mesh sequencing 기법을 통해 수행되었으며, 이는 상대적으로 조밀하지 않은 격자(coarse grid) 해석 결과를 초기 조건으로 사용하여 점진적으로 더 미세한 격자(fine grid) 해석으로 전이하는 방식이다. 엔트로피 수정을 적용하지 않은 경우, 잔차가 $10^{- 9}$ 수준으로 충분히 감소하지 못하거나 진동 또는 발산하는 경향을 나타냈다. 반면, 엔트로피 수정 기법을 적용한 경우에는 모든 격자 수준에서 잔차가 기계 오차 수준까지 안정적으로 감소하며 수렴하는 양상을 보였다. 이러한 차이는 엔트로피 수정 기법을 적용하지 않은 해석에서 충격파 전면에서 발생한 carbuncle 현상이 수치적 교란을 유발하여 시뮬레이션의 수렴성을 저해한데 기인한다. 이후 제시되는 수치해석 결과는 모두 엔트로피 수정 기법이 적용된 경우에 해당한다.

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F5.jpg

Fig. 5.

Effect of entropy fix on residual drop in laminar and RANS cases with various grid resolutions

Fig. 6와 7은 각각 층류와 난류 해석에서 다양한 격자 해상도에 따른 표면 압력과 열유속 분포를 풍상(windward), 중앙(middle), 풍하(leeward) 영역을 따라 나타낸 것이다. 무차원화된 벽면 압력 $({\tilde{p}}_{w})$ 과 열유속 $({\tilde{q^{''}}}_{w})$ 은 식 (9)과 같다. 여기서 $p_{\infty}$ 는 자유류의 압력이며, $K_{\infty}$ 는 자유류의 열전도율, $r = 0.1375 m$ 은 실린더의 반경을 나타낸다.

(9)

\tilde{p_{w}} = \frac{p_{w}}{p_{\infty}}, \tilde{q_{w}^{''}} = \frac{q_{w}^{''} r}{K_{\infty} T_{\infty}}

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F6.jpg

Fig. 6.

Wall pressure and heat flux distributions on the conical surface in laminar cases with various grid resolutions

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F7.jpg

Fig. 7.

Wall pressure and heat flux distributions in RANS cases with various grid resolutions(legend identical to Fig. 6)

격자 수렴성 분석 결과, 전 계산 영역에서 격자 수준(Medium, Fine, Finer)에 관계없이 압력 및 열유속 분포가 전반적으로 잘 일치하는 것으로 나타났다. 특히, Medium 격자에서의 해석 결과를 Finer 격자와 비교했을 때 정량적인 차이가 거의 없어, 해당 격자 수준에서도 해의 수렴성이 충분히 확보되었음을 확인할 수 있다.

받음각 $α = 2^{°}$ 조건 하에서, 전체 형상 중 바람이 풍상 측에서 풍하 측 대비 상대적으로 높은 압력 및 열유속이 나타났다. 이는 유동이 풍하 측에 직접 충돌하면서 강한 압축이 발생하고, 이로 인해 유동의 압력과 온도가 상대적으로 높은 것과 관련이 있다. 이러한 특성은 모든 $x$ 위치에서 일관되게 나타난다. 선두부를 지난 이후 원뿔체 구간에서는 압력 분포는 비교적 일정하게 유지되며, 열유속은 점진적으로 감소하는 경향을 보인다. 실린더 시작점에서는 유동의 팽창으로 인해 압력과 열유속이 모두 급격히 감소하고, 플레어 구간에서는 재압축으로 인해 다시 증가하는 양상이 나타난다.

3.3 벽면 압력 및 열유속 분포

Finer 격자 수준에서 수행된 층류 및 RANS 해석 결과의 벽면 압력 분포(Fig. 8)와 열유속 분포(Fig. 9)를, 2024 High-Fidelity CFD Verification Workshop[8]에 참가한 미국 Sandia National Lab의 수치해석 결과와 비교하였다. Sandia National Lab에서는 SPARC(Sandia parallel aerodynamics and reentry code)[26]를 사용하여 해석을 수행하였으며, 이 코드는 정렬 격자 기반의 셀 중심 유한 체적법(cell-centered finite volume method)을 기반으로 한다. SPARC에서 엔트로피 수정 기법에 대한 명시적 언급은 없지만, Larson shock sensor와 TVD(total variation diminishing) 기법을 적용하여 충격파 영역에서 발생할 수 있는 수치적 불안정성을 효과적으로 억제하였다[14]. 또한, RANS 해석에서는 본 연구와 동일하게 SA-negative 모델이 사용되었다. Fig. 8과 9에 제시된 벽면 압력과 열 유속 분포 비교 결과, 본 연구에서 사용한 CFL3D와 SPARC 해석 간에 매우 유사한 경향을 보였다. 이를 통해 본 연구에서 적용한 수치해석 기법의 정확성과 타당성을 검증하였다.

Fig. 8의 벽면 압력 분포에 있어 층류 해석과 RANS 해석 결과는 선두부부터 원뿔체에 이르기까지 위치와 관계없이 거의 동일한 분포를 보여, 두 해석 간의 차이는 사실상 존재하지 않았다. Fig. 9의 열유속의 경우, 선두부에서 최대 무차원 열유속은 ${\tilde{q^{''}}}_{w} = 26.7 \times 10^{4}$ , 차원값으로는 $q_{w}^{''} = 15.1 M W / m^{2}$ 로 계산되었으며, 층류 해석과 RANS 해석 간에 큰 차이는 나타나지 않았다. 원뿔체 표면에서는 층류 해석 결과, 유동이 하류로 진행함에 따라 열유속이 점진적으로 감소하는 경향을 보인다. 반면, RANS 해석에서는 선두부에서부터 난류 점성 $μ_{t}$ 의 발달로 열유속이 증가하여 전반적으로 층류 해석보다 높은 열유속을 나타냈다.

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F8.jpg

Fig. 8.

Comparison of wall pressure distributions in laminar and RANS simulation between the current and the reference[14]

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F9.jpg

Fig. 9.

Comparison of wall heat flux distribution in laminar and RANS between the current and the reference[14] (legend identical to Fig. 8)

특히 RANS 해석에서의 최대 열유속은 풍상 측에서 ${\tilde{q^{''}}}_{w} = 2.05 \times 10^{4}$ , $q_{w}^{''} = 1.16 M W / m^{2}$ 로 나타났다. 후방의 플레어 구간에서는 재압축 충격파가 발생하는 위치에서 열유속이 급격히 증가하였으며, 풍상 측 기준으로 평균 무차원 열유속은 ${\tilde{q^{''}}}_{w} = 7.01 \times 10^{4}$ , 차원값으로는 $q_{w}^{''} = 3.97 M W / m^{2}$ 로 계산되었다.

Fig. 10은 풍상 측에서 층류 해석 대비 RANS 해석의 열유속 비율 $R_{h e a t f l u x}$ 과, 원뿔체 형상 전체에서의 평균 열유속 $q_{w_{a v g}}^{''}$ 를 비교한 것이다. 열유속 분포의 비 $R_{h e a t f l u x}$ 와 평균 열유속 $q_{w_{a v g}}^{''}$ 에 대한 정의는 식 (10)과 같으며, $A_{c o n e}$ 은 원뿔체 표면의 넓이를 나타낸다.

https://cdn.apub.kr/journalsite/sites/kscfe/2025-030-03/N0500300309/images/jkscfe_2025_303_136_F10.jpg

Fig. 10.

Heat flux ratio between laminar and RANS results on windward side and surface-average heat flux comparison on the conical surface between laminar and RANS results

(10)

R_{h e a t f l u x} = \frac{q_{w_{R A N S}}^{''}}{q_{W_{l a m i n a r}}^{''}}, q_{w_{a v g}}^{''} = \frac{∯ q'' d A}{A_{c o n e}}

풍상 측에서 층류 해석과 비교한 RANS 해석의 열유속 비율 $R_{h e a t f l u x}$ 은 유동이 하류로 진행할수록 점차 증가하며, 최대 약 4.8배까지 높아진다. 평균 열유속 또한 RANS 해석에서 더 높게 나타났는데, 층류의 경우 평균 열유속은 $q_{w_{a v g}}^{''} = 1.44 M W / m^{2}$ 인 반면, RANS 해석에서는 $q_{w_{a v g}}^{''} = 5.85 M W / m^{2}$ 로 약 4배 높은 수준이다. 이는 난류 유동이 층류보다 더 높은 열유속을 유도함을 뒷받침하며, 특히 하류로 갈수록 난류 효과가 커짐에 따라 열유속 차이가 더욱 두드러짐을 보여준다.

4. 결 론(Conclusion)

본 연구는 2024년 High-Fidelity CFD Verification Wroskhop의 일환으로, 검증 모델로 제시된 HIFiRE-1 모델을 대상으로 극초음속 층류 및 난류 경계층 유동 해석을 수행한 사례 분석이다. 해당 워크숍은 충격파 지배 유동 조건에서의 수치 해석자와 수치해석법의 정확도를 상호 비교하고 검증하기 위해 개최되었으며, 모든 참가자는 동일한 형상, 격자, 유동 조건(자유류 마하수 $M_{\infty} = 7.18$ , 받음각 $α = 2^{°}$ ) 하에서 해석을 수행하였으며, 난류 경계층 해석의 경우 RANS 기반의 SA-negative 모델을 공통적으로 사용하였다. 본 연구에서는 Roe flux 차분 기법을 활용하여 수치해석을 진행하였으나, 모델의 선두부에서 발생하는 고강도의 충격파 주변 비물리적인 수치해로 인하여 carbuncle 현상이 발생하며 수치적 불안정성이 나타나는 한계가 확인되었다. 이를 해결하기 위해 엔트로피 수정 기법을 도입하여 소산항을 보정함으로써 Roe flux 차분 기법의 안정성과 수렴성을 개선하였고, 그 효과는 충격파 전면의 형태 안정성과 잔차 감소 이력을 통해 입증되었다.

해석 결과, 격자 해상도에 관계없이 벽면 압력과 열유속 분포는 전반적으로 잘 수렴되었으며, 중간 수준(Medium) 격자에서도 Finer 격자와 유사한 정량적 정확도를 확보할 수 있었다. 또한, 선두부에서 bow 충격파, 원뿔체에서의 경사 충격파, 실린더 구간의 팽창파, 플레어 구간의 재압축 충격파 등 HIFiRE-1 주위에 형성되는 주요 고속 유동 현상을 성공적으로 모사함으로써 해석의 물리적 타당성을 확보하였다. 본 연구에서는 워크숍에서 제공된 Sandia National Lab의 SPARC 코드와 비교하여 매우 유사한 수치해석 결과를 얻었다. 이를 통해 사용된 수치 기법의 정확성과 신뢰성을 검증하였다. 층류와 RANS 간의 비교에서는 RANS 해석에서 원뿔체 표면에서의 열 유속이 평균적으로 약 4배 이상 높게 나타났으며, 특히 유동이 하류로 갈수록 그 차이는 더욱 두드러졌다. 이러한 경향은 난류 경계층의 공력가열 특성이 극초음속 비행체의 열적 설계에 결정적인 영향을 미칠 수 있음을 시사한다. 향후에는 난류와 천이 해석 기법과의 연계를 통해 극초음속 경계층 유동 해석의 예측 정밀도를 한층 향상할 수 있을 것으로 기대된다.

Acknowledgements

이 연구는 2025년 정부(방위사업청)의 재원으로 국방과학연구소의 지원을 받아 수행된 미래도전국방기술 연구개발 사업임(No. 915067201).

References

2005, Bahm, C., Baumann E., Martin, J., Bose, D., Beck, R. and Strovers, B., “The X-43A Hyper-X Mach 7 flight 2 guidance, navigation, and control overview and flight test results,” in AIAA/CIRA 13^thInternational Space Planes and Hypersonics Systems and Technologies Conference, AIAA 2005-3275.

10.2514/6.2005-3275

2006, Turner, J., Horschgen, M., Jung, W., Stamminger, A. and Turner, P., “SHEFEX-hypersonic re-entry flight experiment vehicle and subsystem design, flight performance and prospects,” in 14^thAIAA/AHI Space Planes and Hypersonic Systems and Technologies Conference, AIAA 2006-8115..

10.2514/6.2006-8115

2012, Bowcutt, K., Paull, A., Dolvin, D. and Smart, M., “HIFiRE: An international collaboration to advance the science and technology of hypersonic flight,” in Proceedings of the 28^thInternational Congress of the Aeronautical Sciences, Vol.1, pp.65-76.

2008, Wadhams, T.P., Mundy, E., MacLean, M.G. and Holden, M.S.“Ground test studies of the HIFiRE-1 transition experiment part 1: experimental results,” J. Spacecr. Rockets, Vol.45, No.6, pp.1134-1148.

10.2514/1.38338

2007, Kimmel, R., Adamczak, D., Gaitonde, D., Rougeux, A. and Hayes, J., “HIFiRE-1 boundary layer transition experiment design,” in 45^thAIAA Aerospace Sciences Meeting and Exhibit, AIAA 2007-0534.

10.2514/6.2007-534

2008, MacLean, M., Wadhams, T., Holden, M. and Johnson, H., “Ground test studies of the HIFiRE-1 transition experiment Part 2: computational analysis,” J. Spacecr. Rockets, Vol.45, No.6, pp.1149-1164.

10.2514/1.37693

2014, Yentsch, R.J., Gaitonde, D.V. and Kimmel, R., “Performance of turbulence modeling in simulation of the HIFiRE-1 flight test,” J. Spacecr. Rockets, Vol.51, No.1, pp.117-127.

10.2514/1.A32535

2023, Holst, K.R., Kim, C. and Zahr, M.J., “High-Fidelity CFD Verification Workshop 2024: Shock-Dominated Flows,” in AIAA Scitech 2023 Forum, AIAA 2023-1242.

10.2514/6.2023-1242

1981, Roe, P.L., “Approximate Riemann solvers, parameter vectors, and difference schemes,” J. Comput. Phys., Vol.43, No.2, pp.357-372.

10.1016/0021-9991(81)90128-5

1998, Gottlieb, S. and Shu, C.W., “Total variation diminishing Runge-Kutta schemes,” Math. Comput., Vol.67, No.221, pp.73-85.

10.1090/S0025-5718-98-00913-2

2000, Yee H.C., Vinokur, M. and Djomehri, M.J., “Entropy splitting and numerical dissipation,” J. Comput. Phys., Vol.162, No.1, pp.33-81.

10.1006/jcph.2000.6517

2013, MacCormack, R.W., “Carbuncle computational fluid dynamics problem for blunt-body flows,” J. Aerosp. Inf. Syst., Vol.10, No.5, pp.229-239.

10.2514/1.53684

2003, Kim, S.S., Kim, C., Rho, O.H. and Hong, S.K., “Cures for the shock instability: development of a shock-stable Roe scheme,” J. Comput. Phys., Vol.185, No.2, pp.342-374.

10.1016/S0021-9991(02)00037-2

2024, Hwang, T., Zahr, M.J., Holst, K., Freels, J., Baker, A.J., Batten, P., ... and Johnsen, E., “High-Fidelity CFD Verification Workshop 2024 summary: Shock-Dominated flows,” in AIAA Aviation Forum and Ascend, AIAA 2024-3695.

10.2514/6.2024-3695

1992, Spalart, P. and Allmaras S.R., “A one-equation turbulence model for aerodynamic flows,” in 30^thAerospace Sciences Meeting and Exhibit, AIAA 1992-0439.

10.2514/6.1992-439

2012, Allmaras, S.R. and Johnson, F.T., “Modifications and clarifications for the implementation of the Spalart-Allmaras turbulence model,” in 7^thInternational Conference on Computational Fluid Dynamics(ICCFD), Vol.1902, pp.1-11.

1998, Krist, S.L., Biedron, R.T. and Rumsey, C.L., “CFL3D user’s manual (version 5.0),” NASA Technical Memorandum 1998-208444.

1976, Harten, A., Hyman, J.M., Lax, P.D. and Keyfitz, B.,“On finite-difference approximations and entropy conditions for shocks,” Commun. Pure Appl. Math., Vol 29, No.3, pp.297-322.

10.1002/cpa.3160290305

2006, Bartels, R.E., Rumsey, C.L. and Biedron, R.T., “CFL3D Version 6.4-General Usage and Aeroelastic Analysis,” NASA Technical Memornadum 2006-214301

1987, Yee, H.C., “Upwind and symmetric shock-capturing shcemes,” NASA Technical Memornadum 1987-89464

1997, Harten, A., “High resolution schemes for hyperbolic conservation laws,” J. Comput. Phys., Vol.135, No.2, pp.260-278.

10.1006/jcph.1997.5713

1990, Gnoffo, P.A., “An upwind-biased, point-implicit relaxation algorithm for viscous, compressible perfect-gas flows,” NASA Technical Paper 1990-2953.

10.2514/6.1989-1972

https://acdl.mit.edu/HFW2022/HiFire1/Pointwise/.

1991, Givoli, D., “Non-reflecting boundary conditions,” J. Comput. Phys., Vol.94, No.1, pp.1-29.

10.1016/0021-9991(91)90135-8

2016, Garicano-Mena, J., Lami, A. and Deconinck, H., “An energy-dissipative remedy agains carbuncle: Application to hypersonic flows around blunt bodies,” Comput. Fluids, Vol.133, pp.43-54.

10.1016/j.compfluid.2016.03.028

2020, Smith, J., “Sandia Parallel Aerodynamics Reentry Code (SPARC) - the Future of Production and Research Aerodynamics,” Sandia National Lab, Technical Report SAND-2020-8483R.

10.2172/1648402

Journal of Computational Fluids Engineering ISSN:1598-6071(Print) 3022-6252(Online) 한국전산유체공학회지

Preview

APPLICATION OF ENTROPY FIX TO COMPUTATIONS OF HYPERSONIC LAMINAR AND TURBULENT BOUNDARY LAYERS ON HIFIRE-1

ABSTRACT

MAIN

Fig. 1.

Geometry of HIFiRE-1 consistings of cone, cylinder, and flare[4]

(1)

(2)

(3)

(4)

(5)

(6)

(7)

(8)

Fig. 2.

Computational grid configuration corresponding to Table 1, including boundary conditions and target regions for laminar(cone only) and RANS(whole body: cone, cylinder and flare) simulations

Table 1.

Grid resolution of three grids Medium, Fine, and Finer, in laminar and RANS(turbulence) simulations[8], with non-dimensional wall distance ∆ycone+ at x=0.5m(conical surface, windward) and ∆ycylinder+ at x=1.4m (cylinder surface, windward)

Fig. 3.

Numerical Schlieren image of flow structures around HIFiRE-1 at M∞=7.18 and α=2°

Fig. 4.

Comparison of Mach number isolines around the nose in RANS with or without the entropy fix method

Fig. 5.

Effect of entropy fix on residual drop in laminar and RANS cases with various grid resolutions

(9)

Fig. 6.

Wall pressure and heat flux distributions on the conical surface in laminar cases with various grid resolutions

Fig. 7.

Wall pressure and heat flux distributions in RANS cases with various grid resolutions(legend identical to Fig. 6)

Fig. 8.

Comparison of wall pressure distributions in laminar and RANS simulation between the current and the reference[14]

Fig. 9.

Comparison of wall heat flux distribution in laminar and RANS between the current and the reference[14] (legend identical to Fig. 8)

Fig. 10.

Heat flux ratio between laminar and RANS results on windward side and surface-average heat flux comparison on the conical surface between laminar and RANS results

(10)

Acknowledgements

References

Grid resolution of three grids Medium, Fine, and Finer, in laminar and RANS(turbulence) simulations[8], with non-dimensional wall distance $∆ y_{c o n e}^{+}$ at $x = 0.5 m$ (conical surface, windward) and $∆ y_{c y l i n d e r}^{+}$ at $x = 1.4 m$ (cylinder surface, windward)

Numerical Schlieren image of flow structures around HIFiRE-1 at $M_{\infty} = 7.18$ and $α = 2^{°}$