3D CFD Analysis of Slug Flow Frequency in Liquid-Gas Two-Phase Flow

YongJun Kwon; Hyung-Suk Kim

doi:10.17820/eri.2025.12.4.266

Preview

Original Article

Ecology and Resilient Infrastructure. 31 December 2025. 266-278
https://doi.org/10.17820/eri.2025.12.4.266

3D CFD Analysis of Slug Flow Frequency in Liquid-Gas Two-Phase Flow

슬러그 흐름 및 빈도 특성 분석을 위한 액체-기체 2상 흐름 3차원 수치모의

YongJun Kwon¹

Hyung-Suk Kim²^*

권 용준¹

김 형석²^*

¹Ph.D. Student, Department of Civil and Environmental Engineering, Kunsan National University

²Associate Professor, Department of Civil Engineering, Kunsan National University

¹국립군산대학교 토목환경공학부 박사과정

²국립군산대학교 토목공학과 부교수

^{*Corresponding Author}

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

This is an open-access article distributed under the terms of the Creative Commons Attribution-Non Commercial-No Derives (http://creativecommons.org/licenses/by-nc-nd/4.0/) which permits unrestricted non-commercial use, distribution, and reproduction in any medium, provided the original work is properly cited.

ABSTRACT

There is a need for numerical model to evaluate and reflect the unsteady characteristics of two-phase flow in the design of underground stormwater tunnels. In this study, the interFoam solver in OpenFOAM, capable of three-dimensional two-phase flow analysis, was employed to simulate slug flow and perform a quantitative evaluation of its characteristics. Numerical simulations were conducted under various flow conditions defined by different superficial Reynolds numbers. The results demonstrated that the process of slug formation and development could be reproduced through detailed analysis of slug generation frequency, spatial occurrence, and propagation behavior. For both Re_SG=1,377 and Re_SG=4,320, the total number of slugs generally increased with higher Re_SL values. At Re_SG=4,320, slug generation was concentrated during the initial inflow stage. The slug frequency tended to increase with Re_SL at Re_SG=1,377, whereas at Re_SG=4,320, it showed a temporary increase near the inlet but decreased along the downstream direction. The dimensionless slug frequency increased as Re_SG rose from 1,377 to 4,320 under Re_SL=7,070 and 12,730, whereas no distinct trend was observed with respect to Re_SG when Re_SL=16,500. The Strouhal number showed some deviations from the experimental results when Re_SG≥4,000, but exhibited good agreement with previous experiments when Re_SG≤1,377. The three-dimensional numerical model used in this study is expected to be capable of predicting the behavior of complex slug flow.

Keywords

Slug frequency

3D Numerical analysis

Slug flow

Turbulence model

OpenFOAM

지하방수로에서 발생하는 2상(two-phase) 유동의 비정상적 거동을 예측하고 이를 설계에 반영하기 위한 기술이 필요하다. 본 연구는 3차원 2상 유동해석이 가능한 OpenFOAM의 interFoam 솔버를 활용하여 슬러그 흐름을 모의하고 그 특성에 대한 정량적 분석을 수행하였다. 다양한 겉보기 레이놀즈수(Superficial Reynolds number)에 따라 총 12가지 유동 조건에서 수치모의를 수행하였으며, 그 결과 슬러그 형성과 이동 과정 및 발생 횟수와 위치 분석을 통해 슬러그 발달 과정을 재현할 수 있음을 확인하였다. 슬러그 발생 횟수는 Re_SG=1,377 및 4,320의 경우 Re_SL이 증가함에 따라 전반적으로 많아지는 경향이 나타났으며, Re_SG=4,320일 때 유입 초기에 슬러그 발생 횟수가 집중되는 것으로 모의되었다. 슬러그 빈도는 Re_SG=1,377인 경우 Re_SL이 커지면서 증가하는 경향을 보였으나, Re_SG=4,320인 경우 Re_SL이 커지면서 유입부에서 일시적으로 증가하였지만 하류로 갈수록 감소하였다. 무차원 슬러그 빈도는 Re_SL=7,070과 12,730인 경우 Re_SG=1,377에서 4,320으로 커질 때 증가하는 것으로 나타났으나, Re_SL=16,500인 경우에는 Re_SG 조건에 따라 특정 경향성이 확인되지 않았다. Strouhal number는 Re_SG≥4,000인 조건에서 수리실험과 다소 차이가 발생하였으나, Re_SG≤1,377인 조건에서 기존 수리실험 결과와 잘 일치하는 것으로 나타났다. 본 연구의 3차원 수치모의를 활용하여 복잡한 슬러그 흐름에 대한 거동을 예측할 수 있는 것으로 판단된다.

키워드

슬러그 빈도

3차원 수치해석

슬러그 흐름

난류모델

OpenFOAM

MAIN

1. 서 론
2. 수치모델
2.1 지배방정식
2.2 수치기법
3. 수치해석
3.1 수치모델 검증
3.2 수치모델 구성 및 조건
4. 수치해석 결과
4.1 슬러그 흐름 형성 과정
4.2 슬러그 빈도 분석
4.3 고찰
5. 결 론

1. 서 론

전 지구적 기후 변화로 인해 단기간에 집중되는 강우가 빈번하게 발생하고 있으며, 그 강도와 패턴 또한 과거보다 더욱 불규칙하고 변동성이 커지고 있다. 특히 도심 지역에서는 도로, 보도, 지붕 등 불투수 면적의 증가로 인해 우수가 빠르게 유출되어 배수시설의 용량이 초과되는 사례가 발생하고 있으며, 이는 지하 인프라의 훼손, 수질 악화, 도시 홍수 등 복합적인 문제를 초래하고 있다. 이러한 문제를 완화하기 위해 국내에서는 집중호우 시 대량의 초기 유출수를 신속히 지하로 유도하여 배수할 수 있는 신월 빗물배수터널이 설치되었으며, 추가적인 건설이 계획되어 있다.

지하 방수로는 도심지 침수 피해를 직접적으로 줄일 수 있는 핵심 기반시설로, 지상의 우수가 기존 관거를 통해 방수로 유입부로 집수된 뒤 수직구를 따라 지하 수십 미터 깊이의 본선으로 낙하되어 일시적으로 저류된다. 유입 초기에는 자유수면을 가진 개수로 유동이 나타나지만, 시간이 지남에 따라 관수로 흐름으로 전이될 수 있다. 이 과정에서 액체와 기체가 함께 흐르는 2상 유동(two-phase flow)이 발생하며, 공기주머니(air pocket) 형성, 급격한 압력 변동 및 공기 분출 등의 현상이 나타날 수 있다. 이러한 현상은 지하방수로 시스템 안정성에 부정적인 영향을 미치며, 심한 경우 방수로의 진동이나 누수, 파손으로 이어질 위험이 있다.

현재 국내의 도시 침수 대응 연구는 주로 단상(one-phase) 1차원 모형(SWMM 등)을 이용하여 수행되고 있으나, 이러한 접근은 2상 유동의 복잡한 물리적 거동을 정밀하게 재현하기에는 한계가 있다(Yang et al. 2022). 반면, 2상 유동에 대한 연구는 석유·가스 산업 분야에서 활발히 이루어졌으며, 주로 액체–기체 혼합 유동의 특성 해석에 초점이 맞추어져 왔다(Thaker and Banerjee 2015). 관내 액체와 기체가 동시에 흐를 때 특정 유속 조건에서 슬러그 유동(slug flow)이 형성될 수 있다. 슬러그 유동은 관을 채운 고속의 액체 슬러그와 그 뒤를 따르는 긴 기포가 반복적으로 나타나는 간헐적 유동 형태로, 일반적으로 층상흐름(stratified flow)이나 파형흐름(wavy flow) 상태에서 기체 유속 증가에 따라 액체 표면의 파동이 불안정하게 발달할 때 형성된다. 이러한 불안정성은 Kelvin-Helmholtz 이론으로 설명되며(Taitel and Dukler 1977), 파고가 관 상부에 도달해 처오름(run-up)이 발생하면 액체 슬러그가 형성되고, 뒤따르는 공기층이 슬러그를 하류로 가속 시킨다. 슬러그는 이동 과정에서 전방 유체를 포획하며 성장하거나 일정한 크기를 유지한 채 이동하며, 이러한 유동은 비정상적이고 무작위적인 특성을 보인다(Kushiyama et al. 2003). 슬러그 유동은 층상·파동·플러그 등 다른 유동 형태와 연속적인 관계를 가지며, 빈도, 길이, 속도, 압력 변동 등 다양한 특성을 나타낸다(Mandhane et al. 1974, Taitel and Dukler 1977, Al-Safran 2009). 특히 슬러그 빈도는 관내 불안정성, 압력 변동 그리고 관 부식 등과도 밀접한 관련이 있는 주요 현상이다(Hill et al. 1996). 이를 규명하기 위해 유속과 레이놀즈수 조건을 바탕으로 주로 수리실험이 수행되어 왔다(Gokcal et al. 2010, Thaker and Banerjee 2015, Archibong-Eso et al. 2018, Abdulkadir et al. 2020, Sassi et al. 2022, Cao et al. 2024).

최근에는 이를 기반으로 슬러그 유동에 대한 수치모의 연구도 활발히 이루어지고 있다. Valluri et al. (2008)은 자유수면 추적 기법인 Level Set Method (LSM)를 이용하여 2차원 수치해석을 수행하고 층상 흐름에서 슬러그 발생 시점까지의 발달 과정을 모의하였으며, Vallée et al. (2008)은 PIV 실험과 ANSYS CFX 기반 2차원 해석 결과를 비교하여 슬러그 발달을 재현하였다. De Schepper et al. (2008)은 PLIC (piecewise linear interface calculation) 기법을 적용하여 Baker 차트의 주요 유동 영역을 성공적으로 재현하였다. Lakehal et al. (2012)은 URANS와 LES를 결합한 V-LES (Very Large Eddy Simulation) 기법을 활용하여 3차원 해석을 수행하였으며, 슬러그 발생 주기를 합리적으로 재현할 수 있음을 보였다. Deendarlianto et al. (2016)은 ANSYS Fluent의 VOF (Volume of Fluid) 기법을 이용해 3차원 2상 유동을 모의하고 긴 기포의 형상이 실험과 잘 부합함을 확인하였다. López et al. (2016)은 STAR‑CCM+와 VOF를 적용하여 수평관 내 액체–기체 2상 유동의 부피분율 및 유동 패턴을 실험과 비교하였고, Akhlaghi et al. (2019)은 OpenFOAM 기반 Multi-Fluid VOF 기법을 통해 복잡한 슬러그 거동을 재현하였다. 또한 Elfaki et al. (2021)은 LSM 기법으로 가스–원유 조건에서 겉보기 액체 유속 증가에 따른 슬러그 특성 변화를 분석하였고, Wu et al. (2021)은 OpenFOAM을 이용하여 슬러그의 형성과 이동 과정을 수치적으로 분석하였다. He et al. (2022)은 ANSYS Fluent를 이용한 3차원 CFD 해석을 통해 오리피스가 있는 수평관 내 압력 변동 특성을 검토하였다.

선행 연구들은 1차원 혹은 2차원 해석에 한정되어 실제 발생하는 3차원적 슬러그 유동의 발생과 전이 거동을 충분히 반영하지 못하였다. 3차원 수치해석(V-LES, Multi-Fluid VOF 등)은 보다 상세한 해석이 가능하지만, 격자 세분화와 계산 비용이 매우 크며 난류 특성을 명확히 재현하는데 한계를 보였다. 특히 3차원 수치해석을 활용하여 슬러그 유동을 정량적으로 분석한 연구는 국내외에서도 제한적으로 이루어지고 있으며, 기존 연구의 대부분은 정성적 유동 패턴 재현에 머물러 있다. 국내 수공학 분야에서는 3차원 슬러그 유동에 대한 정량적 해석 연구는 찾아보기 어렵다. 이로 인해 실제 대규모 지하 배수시설 혹은 방수로에서 발생할 수 있는 2상 유동의 비정상적 거동을 예측하고 설계에 반영하기 위한 기반이 부족한 것이 현실이다.

본 연구에서는 OpenFOAM을 활용하여 수평관 내 액체–기체 2상 유동의 슬러그 형성과 발달 과정에 대하여 3차원 수치모의를 수행하였다. 겉보기 레이놀즈수 조건에 따른 슬러그 이동 및 발생과정을 검토하였으며, 슬러그 횟수 및 발생 빈도를 정량적으로 분석하고, 기존 수리실험에서 제시된 무차원 슬러그 빈도와 비교하여 수치모형의 슬러그 흐름에 대한 적용성을 검토하였다.

2. 수치모델

2.1 지배방정식

본 연구에서는 3차원 이상류를 적절하게 모의할 수 있는 OpenFOAM v2412 toolbox (OpenFOAM 2024)를 활용하였다. OpenFOAM에서 액체와 기체 이상류를 모의하기 위해 interFoam 솔버를 활용하였으며, 자유수면 표현에 VOF (volume of fluid) 기법이 적용된다. 유동 해석에 Navier-Stokes 방정식을 레이놀즈 평균한 비압축성 RANS (Reynolds-Averaged Naiver- Stokes) 방정식이 적용된다. 유동 해석을 위해 사용한 연속방정식과 운동량방정식은 Eqs. 1 and 2와 같다.

(Eq. 1)

\frac{\partial {\bar{u}}_{i}}{\partial x_{i}} = 0

(Eq. 2)

ρ \frac{\partial {\bar{u}}_{i}}{\partial t} + ρ {\bar{u}}_{j} \frac{\partial {\bar{u}}_{i}}{\partial x_{j}} = - \frac{\partial \bar{P}}{\partial x_{i}} + \frac{\partial}{\partial x_{j}} [(μ + μ_{t}) (\frac{\partial {\bar{u}}_{i}}{\partial x_{j}} + \frac{\partial {\bar{u}}_{j}}{\partial x_{i}})] + ρ g_{j} + σ κ \frac{\partial α}{\partial x_{i}}

여기서, t는 시간, $u_{i}$ 는 레이놀즈 평균 유속, $P$ 는 압력, 𝜌는 밀도, 𝜇는 점성계수, $μ_{t}$ 는 난류 와점성계수 $g_{j}$ 는 중력가속도 벡터, 𝜎는 표면장력, 𝜅는 두 유체사이의 평균 곡률, 𝛼는 체적분율(volume fraction)이다. 두 유체의 경계면 포착을 위해 적용된 VOF 기법은 스칼라 이송방정식을 해석하여 구한 체적분율 𝛼로 유체의 성질을 다음의 Eqs. 3 and 4로 결정하게 된다.

(Eq. 3)

ρ = α ρ_{w} + (1 - α) ρ_{a}

(Eq. 4)

μ = α μ_{w} + (1 - α) μ_{a}

여기서, $ρ_{w}$ , $ρ_{a}$ 는 각각 액체와 기체의 밀도, $μ_{w}$ , $μ_{a}$ 는 각각 액체와 기체의 점성계수이다. 𝛼가 0이면 기체, 1이면 액체, 0과 1사이에 자유수면이 표현된다. 체적분율의 이송방정식에는 가상 압축항이 포함되며, 수치확산(numerical diffusion)을 방지하고 두 유체의 경계면을 선명하게 유지할 수 있도록 한다. 체적분율 𝛼에 대한 지배방정식은 아래 Eq. 5와 같이 표현된다(Weller 2008).

(Eq. 5)

\frac{\partial α}{\partial t} + \frac{\partial (α u)}{\partial x_{i}} + \frac{\partial (α (1 - α) u_{r})}{\partial x_{j}} = 0

여기서, 마지막 항은 가상 압축항이며, 가상 압축항의 $u_{r}$ 은 인공 압축 유속(artificial compression velocity)으로 유속장 최대값을 반영하고, 자유수면에만 영향을 준다.

2.2 수치기법

지배방정식 계산과정에 유한체적법(finite volume method)을 사용하였으며, 유체 유동 해석을 위한 운동량방정식 이송항의 이산화를 위해 2차 정확도의 중앙차분(central differencing)에 제한자(limiter)를 적용하여 2차 정확도의 특성은 유지하면서 수치적 불안정을 완화하였다. 난류수송방정식에서 이송항은 1차 정확도의 풍상차분(upwind) 기법을 적용하였으며, VOF 방정식의 이송항은 van Leer의 TVD (total bariation diminishing) 기법을 활용하여 이산화 하였다. 유속과 압력 조합으로 PIMPLE 기법을 적용하였으며, 이는 PISO와 SIMPLE 알고리즘을 결합한 방법으로 Courant Number가 급증할 가능성이 있는 불안정한 계산에서 유용한 알고리즘이다. OpenFOAM에서 난류모형은 k-ε, k-ω, LES 모델 등을 적용할 수 있다.

3. 수치해석

3.1 수치모델 검증

액체–기체 2상 유동에 대한 3차원 수치모형 및 난류모델 선정을 위하여 Adouni et al. (2021)의 수리실험 결과와 비교하였다. 해당 실험은 길이 X=5.0 m, 폭 Y=0.075 m, 높이 Z=0.15 m의 사각형 관로에서 수행되었으며, 겉보기 액체 유속은 0.10 m/s, 겉보기 기체 유속은 2.50 m/s 조건이었다(Fig. 1). 수치모의는 동일한 형상과 유입 조건을 구현하였고, 난류 모형은 k-ε 유형의 standard k-ε, RNG k-ε, realizable k-ε 난류모델과 k-ω SST 난류모델을 적용하여 검토하였다. 격자 구성은 X×Y×Z= 1,400×30×75로 총 3.15×10⁶개의 셀을 사용하였다. 초기 유입 조건으로 하부 0.215Z에서 액체가 유입하는 것으로 설정하였으며 나머지 상부 영역은 기체가 유입되는 조건을 적용하였다. 계산은 총 T=10 s 동안 수행하였으며, X= 2.83 m, Y=0.0375 m 지점에서 결과를 추출하여 비교하였다. Fig. 2의 유속 프로파일 비교 결과, 관로 바닥과 상부 부근의 액체와 기체 유속에 대한 수치해석 결과가 실측 실험값과 잘 일치하였다. 반면 액체와 기체 2상 흐름의 경계면 인근에서는 유속을 다소 과소평가하는 경향이 관찰되었다. 그럼에도 전반적인 액체-기체 2상 흐름에 대한 유속을 정교하게 예측하는 것으로 판단되며 특히 난류 모형 중 RNG k–ε 모형이 전반적으로 가장 양호한 재현성을 보였다. 따라서 이후 해석에서는 RNG k–ε 모형을 채택하여 수치모의에 적용하였다.

https://cdn.apub.kr/journalsite/sites/kseie/2025-012-04/N0190120407/images/kseie_2025_124_266_F1.jpg

Fig. 1.

Computational domain for model validation.

https://cdn.apub.kr/journalsite/sites/kseie/2025-012-04/N0190120407/images/kseie_2025_124_266_F2.jpg

Fig. 2.

Comparison of the velocity between the experiment (Adouni et al. 2021) and numerical results.

3.2 수치모델 구성 및 조건

3.2.1 수치모델 구성

Fig. 3에 3차원 수치해석에 사용한 제원과 경계 구성을 제시하였으며, 원형 관수로의 제원은 직경(D) 0.025 m, 길이(x) 8 m이다. 원형 관로의 격자는 1.7 - 3.0 mm 간격으로 구성하였으며 총 격자 개수는 2.15×10⁶개를 사용하였다. 경계조건으로 유입부(inlet)는 유체의 평균 유속 값을 주는 fixed value 활용하였으며, 하부 0.5D는 액체의 유속, 상부 0.5D는 기체의 유속을 적용하였다. 벽면 유속 경계조건으로 유체와 벽 경계부 유속을 0으로 가정하는 no-slip 조건을 활용하였으며, 유출부(outlet)는 zero gradient를 적용하였다.

https://cdn.apub.kr/journalsite/sites/kseie/2025-012-04/N0190120407/images/kseie_2025_124_266_F3.jpg

Fig. 3.

Computational domain of circular pipe.

3.2.2 수치모의 조건

슬러그 흐름은 액체와 기체 경계면 사이의 속도 차이로 발생하며, 다양한 흐름 패턴이 존재한다. Fig. 4에 Thaker and Banerjee (2015)가 액체와 기체의 겉보기 레이놀즈수(Re_SL, Re_SG) 조건에 따라 슬러그 흐름 패턴을 분류한 결과를 제시하였다. 슬러그 흐름은 크게 액체와 기체 경계면 사이의 변화가 없는 층상흐름(stratified flow), 관 상부에 액체가 닿아 처오름 현상이 발생하는 slug formation, 액체와 기체 경계면에 약간의 물결이 형성되는 파형흐름(wavy flow), 기체 함유 여부에 따른 high 또는 less aerated 슬러그 흐름, 그리고 플러그(Plug) 흐름으로 나뉜다. 여기서 플러그 흐름은 액체 덩어리 안에 기체 기포가 거의 없거나 매우 적은 상태로 이동하는 형태의 유동이다(Nydal et al. 1992, Hurlburt and Hanratty 2002, Kadri et al. 2009). 겉보기 레이놀즈수는 관로에 여러 유체가 흐를 때 한 유체가 관로를 완전히 채우고 있다고 가정하여 계산된 레이놀즈수이다. 두 유체에 대한 겉보기 레이놀즈수와 겉보기 유속 식은 아래 Eqs. 6, 7, 8, 9과 같다.

(Eq. 6)

R e_{S L} = \frac{ρ_{L} V_{S L} D}{μ_{L}}

(Eq. 7)

R e_{S G} = \frac{ρ_{g} V_{S G} D}{μ_{g}}

(Eq. 8)

V_{S L} = \frac{Q_{L}}{A}

(Eq. 9)

V_{S G} = \frac{Q_{g}}{A}

여기서, $V_{S L}$ 은 겉보기 액체 유속, $V_{S G}$ 는 겉보기 기체 유속, $Q_{L}$ 은 액체 유량, $Q_{g}$ 는 기체 유량, $A$ 는 원형 관수로 단면적이다.

https://cdn.apub.kr/journalsite/sites/kseie/2025-012-04/N0190120407/images/kseie_2025_124_266_F4.jpg

Fig. 4.

Flow regime map (Thaker and Banerjee 2015) and simulation conditions.

수치해석을 통하여 겉보기 레이놀즈수에 따른 다양한 흐름 패턴과 슬러그 빈도 분석을 위해 총 12가지 조건에 대하여 모의를 수행하였으며, 수행한 액체 레이놀즈수(Re_SL)와 기체 레이놀즈수(Re_SG) 조건은 Table 1에 제시하였다. Fig. 5은 액체와 기체 겉보기 유속 조건에 대한 슬러그 흐름 발생 시작 조건을 나타낸 그래프이며, viscous 켈빈-헬름홀츠(VKH) 이론과 Lin (1985) 및 Andritsos et al. (1989)의 수리실험을 통해 검증되었다. 여기서 VKH 이론에 대한 공식은 Shimada et al. (2025)에서 확인할 수 있다. Fig. 5의 결과는 슬러그 흐름 발생 최소 조건을 의미하며, 액체와 기체의 겉보기 유속 조건이 해당 기준보다 위에 존재한다는 것은 슬러그 흐름이 발생한다는 것을 의미한다. 본 연구에서 수치모의 조건은 슬러그 발생 최소조건을 상회하는 영역에 해당하는 조건으로 하였다.

Table 1.

Numerical cases and flow pattern

Cases	Re_SL	Re_SG	Flow pattern
Case 1	1,420	3,540	Slug Formation (SF)
Case 2	16,500	765	Plug Flow (PF)
Case 3	16,500	4,000	High Aerated Slug (HAS)
Case 4	6,500	1,377	Less Aerated Slug (LAS)
Case 5	7,070	1,377	Less Aerated Slug (LAS)
Case 6	9,200	1,377	Slug and Plug (SAP)
Case 7	12,730	1,377	Plug Flow (PF)
Case 8	16,324	1,377	Plug Flow (PF)
Case 9	2,830	4,320	Slug Formation (SF)
Case 10	7,070	4,320	High Aerated Slug (HAS)
Case 11	12,730	4,320	High Aerated Slug (HAS)
Case 12	16,500	4,320	High Aerated Slug (HAS)

https://cdn.apub.kr/journalsite/sites/kseie/2025-012-04/N0190120407/images/kseie_2025_124_266_F5.jpg

Fig. 5.

Conditions for onset of slug from stratified flow.

4. 수치해석 결과

4.1 슬러그 흐름 형성 과정

슬러그 흐름의 발달 과정은 다음과 같다. 먼저 액체와 기체의 상대적 유속 차이 및 유체와 관 내벽 사이의 전단 효과로 인해 액체 표면에 작은 파동이 발생하게 되고, 액체에 의한 단면적 감소로 기체의 유속이 높아지면서 압력이 낮아지게 된다. 이때 파동의 발달은 중력에 의해 상쇄되려고 하며, 기체의 유속이 중력 효과를 초과할 만큼 큰 경우 관 상부에 액체가 닿으면서 슬러그가 형성된다(Thaker and Banerjee 2015). Fig. 6은 Case 2 조건에서 수치모의 결과로 나타낸 시간에 따른 체적분율(𝛼, volume rate)을 Thaker and Banerjee (2015)의 수리실험으로 제시된 관측 결과와 비교한 그림이며, 𝛼가 1이면 액체, 0이면 기체이다. 유입 초기 층상 흐름 상태에서 수표면은 액체와 기체 사이의 유속 차이로 인해 약간 상승하는 부분이 발생하며, 기체가 통과하는 단면적을 감소시키면서 기체 가속 및 압력 감소를 유발한다(Figs. 6(a) - (d)). 그 다음 Fig. 6(e)와 같이 발생한 파동은 점차 진폭이 증가하여 기체가 지나가는 통로를 차단하게 된다. 이렇게 형성된 파동은 가속된 기체에 의해 액체의 진행 방향 쪽으로 쓸어내려지며, Figs. 6(f) - (i)와 같이 슬러그를 형성되게 된다. Fig. 6의 수리실험 사진에서 동그라미 점선은 슬러그 상류 쪽에서 수위가 평형 수위(equilibrium level)보다 낮게 형성된 곳이다. Fig. 6(j)에서 정량적으로 차이가 발생하지만 슬러그 상류 수위가 평형 수위보다 낮게 형성되는 영역을 확인할 수 있으며, 슬러그 중 일부는 중첩되기도 하고, 일부는 상쇄되어 층상 흐름 상태로 회복되는 과정에서 발달 및 생성과 소멸을 반복하였다. 슬러그 흐름의 발달 과정에 대한 내용은 Taitel and Dukler (1977)의 관측 결과로부터 확인할 수 있으며, 본 연구의 수치모의 결과도 유사한 현상이 재현되는 것을 확인할 수 있었다. 전반적으로 슬러그 흐름 발생 이후 형상은 다소 차이가 나타나는 것을 볼 수 있지만, 슬러그 이동 및 흐름 발달 과정에 대한 현상은 정교하게 모의되고 있는 것으로 판단되다.

https://cdn.apub.kr/journalsite/sites/kseie/2025-012-04/N0190120407/images/kseie_2025_124_266_F6.jpg

Fig. 6.

Comparison of the onset of slug flow.

4.2 슬러그 빈도 분석

슬러그 횟수는 특정 시간 동안 관내 특정 지점을 통과하는 슬러그 수로 정의된다(Al-Safran 2009). 수치해석을 통하여 슬러그를 측정 하는 방법에는 압력 또는 체적분율 데이터의 PSD (power spectral density) 혹은 pwelch를 활용하거나 시간에 따른 체적분율 데이터와 임계값을 활용하여 측정하는 방법 등이 있으며 각 측정 방법마다 한계가 존재한다. PSD 방법은 파동과 슬러그를 구분하는 데 한계가 있으며, 체적분율과 임계값 활용은 높은 기체량이 함유되어 있는 경우 슬러그와 기포를 구분하기 어려운 단점이 있다(Webner et al. 2023). 정확한 슬러그 측정을 위해 Webner et al. (2023)에서 제시된 시간에 따른 체적분율의 line-averaged 데이터( $h_{l}^{l i n e} (x, t)$ )와 임계값을 기반한 측정 방법을 활용하였으며, 슬러그에 기체가 다량 함유된 경우 이를 구분하기 위하여 가시화 도구인 paraview를 활용하여 체적분율 데이터와 비교하면서 임계값을 선정하고 슬러그를 측정하였다. 체적분율의 line-averaged 데이터는 단면에서 등간격의 z방향 수직선을 따라 존재하는 격자 개수를 충분히 포함할 수 있도록 측점(N) 100개를 활용하여 계산되었으며, 식은 Eq. 10과 같다.

(Eq. 10)

h_{l}^{l i n e} (x, t) = \frac{1}{N + 1} (1.5 (α (x, 0, z_{1}, t) + α (x, 0, z_{N}, t) + \sum_{i = 2}^{N - 1} α (x, 0, z_{i}, t))

여기서, x는 길이방향 측점, 0은 폭(y)방향의 중앙 지점, $z_{1}$ 과 $z_{N}$ 은 각각 관의 최하부와 최상부에 있는 측점 위치이며, 상세한 내용은 Webner et al. (2023)을 참고할 수 있다.

Fig. 7은 총 8 m 관의 1 m 간격 측점 지점(=40D, 80D, 120D, 160D, 200D 240D 280D) 마다 시간에 따른 슬러그 발생 횟수를 나타낸 결과이다. Figs. 7(a) - (e)는 Re_SG=1,377 조건인 경우를 나타낸 결과이며, Figs. 7(f) - (i)는 Re_SG=4,320인 조건의 결과가 제시되어 있다. Fig. 7(a)에서 Re_SG=1,377로 일정할 때 Re_SL=6,500인 경우 처음 슬러그는 3 m (120D) 지점에서 발생하기 시작하는 것을 볼 수 있으며, 하류 방향으로 갈수록 10초 이후에 더 많은 슬러그가 발생하는 것으로 나타났다. Fig. 7(b)와 Fig. 7(c)에서 Re_SL=7,070 및 Re_SL=9,200으로 겉보기 액체 레이놀즈수가 증가하면서 슬러그 횟수가 증가함과 동시에 10초 이전에 슬러그가 나타나는 것을 확인할 수 있다. Fig. 7(d)에서 Re_SL=12,730으로 증가하면서 2 m (80D) 지점에서 슬러그가 발생하기 시작하였고, 전반적으로 슬러그 발생 횟수가 증가하는 것으로 나타났다. Fig. 7(e)에서 Re_SL=16,324로 증가한 경우 슬러그 발생 횟수가 점차 많아지면서 10초 이전의 슬러그 발생 횟수가 증가하는 특성을 보였다. 종합적으로 Re_SG=1,377의 경우 Re_SL=6,500에서 16,324로 점차 증가하면서 슬러그 발생 횟수가 증가하는 경향이 나타났다. Fig. 7(f)를 보면, Re_SG=4,320으로 일정할 때 Re_SL=2,830 조건은 1 m (40D) 지점에서 슬러그가 발생하기 시작하였으며, 대부분의 슬러그가 10초 이전, 특히 2초 내에 발생하는 것으로 나타났다. Fig. 7(g)에서 겉보기 액체 레이놀즈수 Re_SL=7,070으로 증가하면서 2초 내에 발생하는 슬러그 횟수가 증가하였지만, Fig. 7(h)와 같이 Re_SL=12,730으로 증가하면서 슬러그는 하류로 갈수록 시간에 따라 고른 분포를 보여주었다. Fig. 7(i)에서 Re_SL=16,500으로 증가하는 경우 더 많은 슬러그가 발생하였고, 겉보기 액체 레이놀즈수가 높아지면서 슬러그 횟수가 증가하는 경향을 볼 수 있었다.

https://cdn.apub.kr/journalsite/sites/kseie/2025-012-04/N0190120407/images/kseie_2025_124_266_F7.jpg

Fig. 7.

Number of slug along the pipeline.

Re_SG=1,377 및 4,320의 경우 모두 Re_SL이 증가함에 따라 슬러그 발생 횟수가 증가하는 경향이 나타났다. Re_SG가 증가한 경우 1-2 m (40D-80D)에서 슬러그 발생 횟수가 증가하였고, 유입 초반 1-2초 내에 발생하는 슬러그가 증가하는 경향을 확인할 수 있었다. 이러한 현상은 겉보기 액체 및 기체 유속이 증가할수록 슬러그가 상류에서 주로 발생한다는 수리실험 내용과 일치하며, 관 길이에 따른 슬러그 형성 위치(30D-60D) 및 횟수는 정성적으로 모의가 가능함을 보였다.

Figs. 8(a) - (b)는 거리에 따른 슬러그 빈도(f_s)를 나타낸 결과이다. Fig. 8(a)에서 슬러그 빈도는 Re_SG=1,377로 일정할 때 Re_SL이 증가함에 따라 점차 높은 값이 나타나는 것을 볼 수 있으며, 전반적으로 하류로 갈수록 점차 감소하는 경향이 발생하였다. Re_SL=9,200까지 슬러그 빈도는 비교적 유사하게 나타났지만, Re_SL=12,730으로 증가하면서 슬러그 빈도가 급격하게 상승하였다. 특히, Re_SL이 증가한 경우 슬러그는 보다 더 상류 쪽에서 발생하기 시작하였으며, 이 현상은 Cao et al. (2024)의 수리실험 결과와 일치한다. Re_SG=1,377인 경우 Re_SL=12,730으로 높아지면서 슬러그 빈도는 x=3 m 지점에서 최댓값이 1.73으로 급격한 상승이 나타났다. Fig. 8(b)에서 Re_SG=4,320로 일정할 때 Re_SL 4가지 조건의 슬러그 빈도는 유입부 주변 x=2-3 m 까지 급격한 상승을 보여주었지만, x=3-4 m 지점에서 급격하게 감소하는 것으로 나타났으며, 하류에서 약간 상승하다가 다시 감소하는 경향이 발생하였다. Re_SG=4,320인 경우 슬러그 빈도 최대값은 Re_SL=7,070에서 3.95로 나타났다. Re_SG=4,320으로 높은 경우 슬러그 빈도는 Re_SL이 증가함에 따라 일관적인 특성을 보이지 않고 무작위적으로 나타나면서 다른 양상을 보여주었으며, 이전의 수리실험 경향과 다소 차이를 보였다(Andritsos et al. 1989, Archibong-Eso et al. 2018). Fig. 8(c)는 거리에 따른 무차원 슬러그 빈도를 제시한 결과이다. 무차원 슬러그 빈도는 Strouhal number (f_sD/V_M)와 Froude number $(V_{M}^{2} / g D)$ 의 곱으로 표현되며, 여기서 V_M=V_SL+V_SG을 나타낸다. Fig. 8(c)에서 Re_SL=7,070일 때 Re_SG가 1,377에서 4,320으로 증가하면 무차원 슬러그 빈도는 크게 상승하는 것으로 나타났으며, Re_SL=12,730인 경우도 유사한 경향이 발생하는 것으로 나타났다. Re_SL=16,500으로 증가한 경우 Re_SG가 765에서 4,000으로 증가하였음에도 무차원 슬러그 빈도는 감소하는 것으로 나타났으며, Re_SG=4,320으로 증가하면서 무차원 슬러그 빈도는 x=2 m (80D)에서 크게 증가하고 이후 감소하는 것으로 나타났다.

https://cdn.apub.kr/journalsite/sites/kseie/2025-012-04/N0190120407/images/kseie_2025_124_266_F8.jpg

Fig. 8.

Slug frequency and non-dimensional slug frequency for distance of x.

일반적으로 슬러그 빈도는 Re_SL 값이 증가함에 따라 커지고, 하류로 갈수록 점차 감소하는 경향이 있음을 여러 수리실험 연구를 통해 확인되어 왔다(Woods and Hanratty 1999, Woods et al 2006, Wang et al. 2007, Dinaryanto et al. 2017). 수치해석 결과도 Re_SG=1,377인 경우 Re_SL 값이 커질수록 슬러그 빈도가 증가하는 것으로 나타나면서 기존 수리실험 결과와 동일한 경향을 보였다. 하지만 Re_SG=4,320인 경우 Re_SL이 커질수록 슬러그 빈도가 증가하는 경향이 발생하지 않았으며, 같은 Re_SL 조건에서 Re_SG가 높아지는 경우도 슬러그 빈도는 상이한 결과를 보이면서 수리실험 결과와 다소 차이가 나타났다.

Fig. 9는 수치모의를 통한 무차원 슬러그 빈도에 대한 결과를 Fossa et al. (2003)과 Wang et al. (2007)의 수리실험 결과에서 제안한 상관관계와 비교하였다. 여기서 가로축에 혼합 흐름에서 액체 비율(V_L/V_M)을 세로축에는 Strouhal Number (f_s/D/V_SG)를 제시하였다. Fig. 8의 결과에서도 볼 수 있지만 슬러그 빈도는 거리와 조건에 따라 급격한 증가 또는 감소를 보였다. 12가지 조건들에서 수행한 수치해석 결과로 40D-280D 구간에 대하여 Strouhal number의 중앙값(median)을 제시하였으며, 거리에 따라 발생한 Strouhal number의 표준편차를 함께 도시하였다. Fig. 9에서 수치해석 결과의 Strouhal number는 수리실험 결과와 비교하여 전반적으로 잘 일치하는 것을 볼 수 있다. 겉보기 기체 레이놀즈수에 따라 비교해 보면, Strouhal number는 특히 Re_SG=1,377 이하인 조건에서 대부분 잘 일치하는 것으로 나타났다. Re_SG=4,000 이상으로 높아지는 경우 거리에 따른 Strouhal number 편차가 큰 경향을 보이고, 수리실험과 다소 차이를 보이는 것으로 나타났다.

https://cdn.apub.kr/journalsite/sites/kseie/2025-012-04/N0190120407/images/kseie_2025_124_266_F9.jpg

Fig. 9.

Comparison of non-dimensional slug frequency between experiments data and simulation results.

슬러그 형성은 본질적으로 복잡하고 불안정한 특성을 내포하고 있기 때문에 이를 정확하게 수치모의하기는 매우 어렵다. 본 연구의 OpenFOAM을 활용한 3차원 수치해석 결과로 슬러그 형성 과정이 잘 재현되는 것을 확인하였으며, 특히, Re_SG=1,377 이하의 낮은 겉보기 기체 유속 조건에서 발생한 슬러그 특성과 현상이 수리실험에 나타난 경향과 잘 일치하는 것을 볼 수 있었다.

4.3 고찰

슬러그 흐름의 본질적인 불안정하고 무작위적인 특성으로 인해 수치모델을 활용하여 슬러그 발생 빈도와 개별 슬러그의 거동을 정밀하게 재현하는 데 한계가 존재한다. 자유수면 처리 방법, 격자 해상도 등이 슬러그 특성을 해석하는데 영향을 미치며, 특히, V-LES나 Multi-VOF 등 보다 정밀한 수치기법을 적용할 수 있으나 계산 시간과 물리적 하드웨어를 고려한 현실적인 제약을 간과하지 않을 수 없다. 또한, 슬러그 빈도 정량화에 대한 연구는 많은 수리실험을 통해 검토되었지만 여전히 명확하게 정립되지 않았다. 실제로 수평관 내 액체-기체 2상 흐름에서 슬러그 빈도는 유속, 관경, 관 길이, 유체의 점성과 밀도, 유입조건 등의 여러 변수에 민감하게 반응하며, 이 변수들 간의 상관관계도 실험 마다 다양하게 나타난다(Cao et al. 2024).

본 연구에서 Re_SG에 따른 슬러그 빈도 분석 결과 역시 이러한 불규칙성을 반영하고 있는 것으로 판단된다. Thaker and Banerjee (2015)에서 같은 Re_SL조건에서 Re_SG가 높아질수록 슬러그 빈도는 서서히 증가하는 경향을 보여주었으며, Re_SG가 높은 슬러그 빈도 형성에 미치는 영향은 미미하다고 보고되었다. 하지만 Al-Kayiem et al. (2017)은 동일한 Re_SL에서 Re_SG가 증가하면서 슬러그 빈도가 감소한다고 제안하였으며, Archibong-Eso et al. (2018)는 기체의 겉보기 레이놀즈수에 따라 Re_SG≤2,500일 때, 슬러그 빈도는 V_SG가 높아지면서 증가하였고, Re_SG≥2,500일 때, 슬러그 빈도는 V_SG증가에 따라 감소하는 것을 확인하였다. 또한 Cao et al. (2024)의 수리실험 결과로 같은 겉보기 액체 유속 조건에서 V_SG≤1.0 m/s일 때 슬러그 빈도는 V_SG가 높아질수록 증가하다가, V_SG≥1.0 m/s일 때 V_SG가 높아질수록 슬러그 빈도는 감소하는 경향을 보인다고 제시하였다. 본 수치모의 결과는 동일한 Re_SL 조건에서 Re_SG가 1,377에서 4,320으로 높아지는 경우 슬러그 빈도가 증가하는 경향을 보였다. 그러나 Re_SG가 4,000에서 4,320으로 높아지면서 슬러그 빈도가 감소하는 실험결과와 달리 오히려 증가하는 경우도 발생하면서 기존 수리실험 결과들과 상이함을 나타냈다. Andritsos et al. (1989)는 V_SG> 1.56 m/s로 상대적으로 높은 겉보기 기체 유속의 경우 관 입구로부터 먼 거리인 260D 이후 슬러그가 발생하기 시작하는 것을 관측하였으며, Cao et al. (2024)는 V_SG≥1.0 m/s인 경우 슬러그의 준발달(quasi-developed) 거리를 약 1,200D로 제시하였다. 일반적으로 슬러그 빈도는 하류로 갈수록 기체로 인한 점성감쇠에 의해 감소하는 경향을 보인다(Thaker and Banerjee 2015). 수치해석 결과의 차이는 Re_SG=4,000 이상의 조건인 경우 V_SG≥1.0 m/s으로 높아지면서 슬러그가 상당히 발달할 정도로 관 길이가 충분하지 않기 때문에 나타난 것으로 추측할 수 있다. 이러한 경우 충분한 길이의 관에서 수치모의를 수행하거나 Schmelter et al. (2021)이 제안한 유입부 유속을 사인파 형태로 적용하는 인위적인 교란을 활용한다면 더 나은 결과를 도출할 수 있을 것으로 판단된다.

5. 결 론

본 연구는 OpenFOAM을 활용하여 액체와 기체 2상 유체에 대한 3차원 수치모의를 수행하여 슬러그 흐름에 대한 분석을 수행하였다. 2상 유동 해석을 위해 자유수면 해석 기법인 VOF 기반의 interFoam 솔버를 사용하였다. 여러 난류모델 비교를 통해 RNG k–ε 난류모델을 선정하고 수평관 직경 0.025 m, 길이 8 m에서 겉보기 액체 및 기체 레이놀즈수에 따라 수치모의를 수행하였다. 겉보기 액체 및 기체 레이놀즈수 조건에 따른 슬러그 이동 및 발달 과정을 검토하였고, 슬러그 빈도를 측정하여 거리에 따른 특성을 분석하여 수리실험과 비교하였다. 3차원 수치해석을 통해 도출된 결론은 다음과 같다.

(1)3D 수치모의를 활용하여 파동의 발달, 액체로 인한 단면적 감소, 처오름 등 일련의 슬러그 이동과 발달 과정 및 현상이 합리적으로 재현 되었음을 확인하였다.

(2)슬러그 발생 횟수는 Re_SG=1,377 및 4,320의 경우 모두 Re_SL가 증가함에 따라 전반적으로 많아지는 경향을 보였다. 특히, Re_SG=4,320인 경우 유입 초반 1-2초 내에 유입부 근처에 슬러그가 집중되는 것으로 나타났으며, 이는 기존 수리실험 경향과 일치하였다.

(3)슬러그 빈도는 Re_SG=1,377의 경우 Re_SL가 6,500-16,324 조건에서 점차 커질수록 전반적으로 증가하였지만, Re_SG=4,320의 경우 Re_SL가 2,830-16,500로 커지면서 유입부에서 일시적으로 증가하였으나 하류로 갈수록 감소하였다.

(4)무차원 슬러그 빈도는 Re_SL=7,070 및 Re_SL=12,730 조건인 경우 Re_SG가 1,377에서 4,320으로 높아지지면서 증가하는 것으로 나타났지만, Re_SL=16,500인 경우 Re_SG 조건에 따라 특정 경향성이 확인되지 않았다.

(5)액체 비율에 대한 Strouhal number는 Re_SG≤1,377의 조건에서 수리실험 결과와 잘 일치하는 것으로 나타났지만, Re_SG≥4,000의 조건인 경우 수리실험과 다소 차이를 보여주었다. 이러한 차이는 슬러그 조건에서의 흐름이 불안정성과 무작위성 및 기체에 의한 점성감쇠 현상으로 인해 발생한 것으로 판단된다.

향후 본 연구에서 구축한 3차원 수치모델을 지하방수로 해석에 확장 적용하여 다양한 관경과 비정상 유입 조건을 반영한 공기 연행 및 공기 갇힘 등의 현상 규명 연구에 활용할 수 있을 것으로 기대된다.

Acknowledgements

This work is financially supported by Korea Ministry of Climate, Energy, Environment(MCEE) as 「Technology development project to optimize planning, operation, and maintenance of urban flood control facilities)(RS-2024-00397821)」.

References

Abdulkadir, M., Zhao, D., Abdulkareem, L.A., Asikolaye, N.O., and Hernandez-Perez, V. 2020. Insights into the transition from plug to slug flow in a horizontal pipe: An experimental study. Chemical Engineering Research and Design 163: 85-95.

10.1016/j.cherd.2020.08.025

Adouni, H., Chouari, Y., Kriaa, W., and Bournot, H. 2021. The effect of air velocity on slugs in a confined channel. Flow Measurement and Instrumentation 81: 1-11.

10.1016/j.flowmeasinst.2021.102013

Akhlaghi, M., Mohammadi, V., Nouri, N.M., Taherkhani, M., and Karimi, M. 2019. Multi-Fluid VoF model assessment to simulate the horizontal air–water intermittent flow. Chemical Engineering Research and Design 152: 48-59.

10.1016/j.cherd.2019.09.031

Al-Kayiem, H.H., Mohmmed, A.O., Al-Hashimy, Z.I., and Time, R.W. 2017. Statistical assessment of experimental observation on the slug body length and slug translational velocity in a horizontal pipe. International Journal of Heat and Mass Transfer 105: 252-260.

10.1016/j.ijheatmasstransfer.2016.09.105

Al-Safran, E. 2009. Investigation and prediction of slug frequency in gas/liquid horizontal pipe flow. Journal of Petroleum Science and Engineering 69(1-2): 143-155.

10.1016/j.petrol.2009.08.009

Andritsos, N., Williams, L., and Hanratty, T.J. 1989. Effect of liquid viscosity on the stratified-slug transition in horizontal pipe flow. International Journal of Multiphase Flow 15(6): 877-892.

10.1016/0301-9322(89)90017-7

Archibong-Eso, A., Baba, Y., Aliyu, A., Zhao, Y., Yan, W., and Yeung, H. 2018. On slug frequency in concurrent high viscosity liquid and gas flow. Journal of Petroleum Science and Engineering 163: 600-610.

10.1016/j.petrol.2017.12.071

Cao, Y., Xu, Q., Yu, H., Huang, B., Liu, T., and Guo L. 2024. Influence of pipeline diameters and fluid properties on slug frequency in horizontal pipelines. Chemical Engineering Science 283: 1-20.

10.1016/j.ces.2023.119397

De Schepper, S.C.K., Heynderickx, G.J., and Marin, G.B. 2008. CFD modeling of all gas–liquid and vapor–liquid flow regimes predicted by the Baker chart. Chemical Engineering Journal 138(1-3): 349-357.

10.1016/j.cej.2007.06.007

Deendarlianto, Andrianto, M., Widyaparaga, A., Dinaryanto, O., Khasani, and Indarto 2016. CFD Studies on the gas-liquid plug two-phase flow in a horizontal pipe. Journal of Petroleum Science and Engineering 147: 779-787.

10.1016/j.petrol.2016.09.019

Dinaryanto, O., Prayitno, Y., Majid, A.I., Hudaya, A.Z., Nusirwan, Y.A., Widyaparaga, A., Indarto, and Deenarlianto. 2017. Experimental investigation on the initiation and flow development of gas-liquid slug two-phase flow in a horizontal pipe. Experimental Thermal and Fluid Science 81: 93-108.

10.1016/j.expthermflusci.2016.10.013

Elfaki, M., Nasif, M.S., and Muhammad, M. 2021. Effect of Changing Crude Oil Grade on Slug Characteristics and Flow Induced Mechanical Stresses in Pipes. Applied Sciences 11(11): 1-19.

10.3390/app11115215

Fossa, M., Guglielmini, G., Marchitto, A. 2003. Intermittent flow parameters from void fraction analysis. Flow Measurement and Instrumentation 14(4-5): 161-168.

10.1016/S0955-5986(03)00021-9

Gokcal, B., Al-Sarkhi, A.S., Sarica, C., and Al-Safran, E.M. 2010. Prediction of slug frequency for high-viscosity oils in horizontal pipes. SPE Projects Facilities and Construction 5(3): 136-144.

10.2118/124057-PA

He, J., Hou, Q., Lian, J., Tijsseling, A.S., Bozkus, Z., and Laanearu J. 2022. Three-dimensional CFD analysis of liquid slug acceleration and impact in a voided pipe line with end orifice. Engineering Applications of Computational Fluid Mechanics 16(1): 1444-1463.

10.1080/19942060.2022.2095440

Hill, T.J., Fairhurst, C.P., Nelson, C.J., Becerra, H., and Bailey, R.S. 1996. Multiphase Production Through Hilly Terrain Pipelines in Cusiana Oilfield, Colombia. Proceedings of the 1996 SPE Annual Technical Conference and Exhibition, SPE, Colorado, USA, pp. 319-334.

10.2118/36606-MS

Hurlburt, E.T. and Hanratty, T.J. 2002. Prediction of the transition from stratified to slug and plug flow for long pipes. International Journal of Multiphase Flow 28(5): 707-729.

10.1016/S0301-9322(02)00009-5

Kadri, U., Mudde, R.F., Oliemans, R.V.A., Bonizzi, M., and Andreussi P. 2009. Prediction of the transition from stratified to slug flow or roll-waves in gas-liquid horizontal pipes. International Journal of Multiphase Flow 35(11): 1001-1010.

10.1016/j.ijmultiphaseflow.2009.07.002

Kushiyama, K., Kamioka, S., and Yamada, T. 2003. Development of numerical model concerning the behavior of water and air in a long inverted siphon. -A proposal of unsteady numerical model based on hydraulic experiment model. Journal of Japan Society Hydrology and Resources 16(5): 527-540. (in Japan)

10.3178/jjshwr.16.527

Lakehal, D., Labois, M., and Narayanan, C. 2012. Advances in the Large-Eddy and Interface Simulation (LEIS) of interfacial multiphase flows in pipes. Progress in Computational Fluid Dynamics 12(2-3): 153-163.

10.1504/PCFD.2012.047458

Lin, P.Y. 1985, Flow Regime Transitions in Horizontal Gas-liquid Flow. Ph.D. dissertation, University of Illinois, Urbana, USA, pp. 77-155.

López, J., Pineda, H., Bello, D., and Ratkovich, N. 2016. Study of liquid–gas two-phase flow in horizontal pipes using high speed filming and computational fluid dynamics. Experimental Thermal and Fluid Science 76: 126-134.

10.1016/j.expthermflusci.2016.02.013

Mandhane, J.M., Gregory, G.A., and AzizA, K. 1974. Flow pattern map for gas—liquid flow in horizontal pipes. International Journal of Multiphase Flow 1(4): 537-553.

10.1016/0301-9322(74)90006-8

Nydal, O.J., Pintus, S., and Andreussi, P. 1992. Statistical characterization of slug flow in horizontal pipes. International Journal of Multiphase Flow 18(3): 439-453.

10.1016/0301-9322(92)90027-E

OpenFOAM 2024. OpenFOAM - The open source CFD toolbox 2024 User’s Guide.

Sassi, P., Fernández, G., Stiriba, Y., and Pallarès, J. 2022. Effect of solid particles on the slug frequency, bubble velocity and bubble length of intermittent gas–liquid two-phase flows in horizontal pipelines. International Journal of Multiphase Flow 149: 1-22.

10.1016/j.ijmultiphaseflow.2022.103985

Schmelter, S., Knotek, S., Olbrich, M., Fiebach, A., and Bär, M. 2021. On the influence of inlet perturbations on slug dynamics in horizontal multiphase flow—a computational study. Metrologia 58(1): 1-14.

10.1088/1681-7575/abd1c9

Shimada, N., Sato, M., Yaegashi, Y., Koch, T., and Le, T.T. 2025. An Application of Modified S-CLSVOF Method to Kelvin-Helmholtz Instability and Comparison with Theoretical Result. Journal of Chemical Engineering of Japan 58(1): 1-8.

10.1080/00219592.2025.2451953

Taitel, Y. and Dukler, A.E. 1977. A model for slug frequency during gas-liquid flow in horizontal and near horizontal pipes. International Journal of Multiphase Flow 3(6): 585-596.

10.1016/0301-9322(77)90031-3

Thaker, J. and Banerjee, J. 2015. Characterization of two-phase slug flow sub-regimes using flow visualization. Journal of Petroleum Science and Engineering 135: 561-576.

10.1016/j.petrol.2015.10.018

Vallée, C., Höhne, T., Prasser, H., and Sühnel T. 2008. Experimental investigation and CFD simulation of horizontal stratified two-phase flow phenomena. Nuclear Engineering and Design 238(3) 637-646.

10.1016/j.nucengdes.2007.02.051

Valluri, P., Spelt, P.D.M., Lawrence, C.J., and Hewitt, G.F. 2008. Numerical simulation of the onset of slug initiation in laminar horizontal channel flow. International Journal of Multiphase Flow 34(2): 206-225.

10.1016/j.ijmultiphaseflow.2007.09.001

Wang, X., Guo, L., and Zhang X. 2007. An experimental study of the statistical parameters of gas–liquid two-phase slug flow in horizontal pipeline. International Journal of Heat and Mass Transfer 50(11-12): 2439-2443.

10.1016/j.ijheatmasstransfer.2006.12.011

Webner, F., Polansky, J., Knotek, S., and Schmelter, S. 2023. Sensitivity analysis of threshold parameters in slug detection algorithms. International Journal of Multiphase Flow 158: 1-14.

10.1016/j.ijmultiphaseflow.2022.104278

Weller, H.G. 2008. A new approach to VOF-based interface capturing methods for incompressible and compressible flow. OpenCFD Ltd., Report TR/HGW 4: 35.

Woods, B.D. and Hanratty T.J. 1999. Influence of Froude number on physical processes determining frequency of slugging in horizontal gas liquid flows. International Journal of Multiphase Flow 25(6-7): 1195-1223.

10.1016/S0301-9322(99)00058-0

Woods, B.D., Fan, Z., and Hanratty, T.J. 2006. Frequency and development of slugs in a horizontal pipe at large liquid flows. International Journal of Multiphase Flow 32(8): 902-925.

10.1016/j.ijmultiphaseflow.2006.02.020

Wu, X., Wang, Z., Dong, M., Ge, Q., and Dong, L. 2021. Simulation Study on the Development Process and Phase Interface Structure of Gas-Liquid Slug Flow in a Horizontal Pipe. Frontiers in Energy Research 9: 1-13.

10.3389/fenrg.2021.762471

Yang, Q., Ma, Z., and Zhang, S. 2022. Urban Pluvial Flood Modeling by Coupling Raster-Based Two-Dimensional Hydrodynamic Model and SWMM. Water 14(11): 1-16.

10.3390/w14111760

Ecology and Resilient Infrastructure ISSN:2288-8527(Online) 응용생태공학회 논문집

Preview

3D CFD Analysis of Slug Flow Frequency in Liquid-Gas Two-Phase Flow

ABSTRACT

MAIN

(Eq. 1)

(Eq. 2)

(Eq. 3)

(Eq. 4)

(Eq. 5)

Fig. 1.

Computational domain for model validation.

Fig. 2.

Comparison of the velocity between the experiment (Adouni et al. 2021) and numerical results.

Fig. 3.

Computational domain of circular pipe.

(Eq. 6)

(Eq. 7)

(Eq. 8)

(Eq. 9)

Fig. 4.

Flow regime map (Thaker and Banerjee 2015) and simulation conditions.

Table 1.

Numerical cases and flow pattern

Fig. 5.

Conditions for onset of slug from stratified flow.

Fig. 6.

Comparison of the onset of slug flow.

(Eq. 10)

Fig. 7.

Number of slug along the pipeline.

Fig. 8.

Slug frequency and non-dimensional slug frequency for distance of x.

Fig. 9.

Comparison of non-dimensional slug frequency between experiments data and simulation results.

Acknowledgements

References