Prediction of Turbulent Flows in a Square Duct Subjected to Lorentz Force

In-Gyu Jeon; Ho-Yong Park; Dong-Kurl Kwak; Jong-Keun Shin

doi:10.21729/ksds.2025.18.2.1

Preview

Original Article

Journal of Korean Society of Disaster and Security. 30 June 2025. 1-12
https://doi.org/10.21729/ksds.2025.18.2.1

Prediction of Turbulent Flows in a Square Duct Subjected to Lorentz Force

로렌츠 힘이 작용하는 정사각 덕트 내 난류유동 해석

In-Gyu Jeon¹

Ho-Yong Park¹

Dong-Kurl Kwak²^*

Jong-Keun Shin²^*

전 인규¹

박 호용¹

곽 동걸²^*

신 종근²^*

¹Ph.D. Course, Graduate School of Disaster Prevention, Kangwon National University

²Professor, Graduate School of Disaster Prevention, Kangwon National University

¹강원대학교 방재전문대학원 박사과정

²강원대학교 방재전문대학원 교수

^{*Corresponding Author}

License (open-access, https://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 (https://creativecommons.org/licenses/by-nc/4.0/) which permits unrestricted non-commercial use, distribution, and reproduction in any medium, provided the original work is properly cited.

ABSTRACT

This study was conducted to understand the effects of a magnetic field on mean velocity and turbulent characteristics in fully developed turbulent magnetohydrodynamic (MHD) flow through a square duct. Based on the coupled Navier-Stokes MHD equations, an extended Reynolds stress transport model including the Lorentz force effect was applied, and the analysis was performed using an in-house computational fluid dynamics code. In order to directly apply the Lorenz force effect to the turbulent MHD flow, models related to the generation and dissipation rate of the Lorenz force were additionally applied to the Reynolds stress and dissipation rate equations, and the indirect effect was also considered through the modification of the velocity-pressure gradient correlation model of the Reynolds stress equation. In addition, the model coefficients of the models related to the Lorenz force were tuned using the fully developed channel MHD data so that they could be linked to the elliptic-blending Reynolds stress model and then the square duct flow was analyzed. The prediction results were validated through comparison with the direct numerical simulation (DNS) data of Chaudhary et al. (2010), confirming the flattening of the mean velocity distribution and the suppression trends of turbulence intensity and shear stress. The results of the analysis show that the proposed extended model contributes to improving the predictive accuracy of MHD-affected flow fields and suggests its applicability to various engineering applications involving complex boundary conditions.

Keywords

magnetohydrodynamic

Reynolds stress transport equation

turbulence model

elliptic-blending

Lorentz force

본 연구는 정사각형 단면을 갖는 덕트 내 완전 발달된 난류 자기유체역학(magnetohydrodynamic, MHD) 유동에서 자기장이 평균속도와 난류 특성에 미치는 영향을 이해하기 위하여 수행되었다. 결합된 형태의 나비에-스토크스(Navier-Stokes) MHD 방정식들을 기반으로 로렌츠 힘 효과를 포함한 확장된 형태의 레이놀즈응력 수송 모형을 적용하였으며, 이 해석은 전산 유체역학을 위한 자체 개발 코드로 수행되었다. 난류 MHD 유동에 로렌츠 힘 효과를 직접적으로 부여하기 위하여 레이놀즈응력과 소멸률 방정식들에 로렌츠 힘 관련 생성과 소멸에 관한 모형들을 추가하였으며, 레이놀즈응력 방정식의 속도-압력구배 상관관계 모형의 수정을 통해서도 간접적인 효과를 발휘하도록 고려하였다. 또한, 본 연구에서는 타원-혼합 레이놀즈응력 모형과 연계되도록 로렌츠 힘과 관련된 모형계수들을 완전 발달된 평판 MHD 자료들을 활용하여 조율하는 과정을 거친 후에 사각덕트 유동을 해석하였다. 예측된 결과는 Chaudhary et al.(2010)의 직접수치모사(direct numerical simulation, DNS) 자료와 비교하여 타당성을 검증하였으며, 평균속도 분포의 평탄화와 난류강도 및 전단응력의 억제 경향을 확인하였다. 해석 결과는 제안한 확장형 모형이 자기장 영향을 받는 유동장의 예측 정확도 향상에 기여할 수 있음을 보여주고 복잡한 경계조건을 갖는 다양한 공학적 응용 분야에서 활용 가능성을 제시한다.

키워드

자기유체역학

레이놀즈응력 수송방정식

난류 모형

타원-혼합

로렌츠 힘

MAIN

1. 서 론
2. 이론 해석
2.1 MHD 유동 관련 운동량 방정식
2.2 난류 유동 해석을 위한 타원-혼합 모형
2.3 난류 모형에 로렌츠 힘 효과 항 부여
2.4 로렌츠 힘 효과가 고려된 난류 모형의 모형 계수
3. 수치해석
4. 결과 및 고찰
5. 결 론

1. 서 론

다양한 산업분야에서 자기장이 부가된 전도성 유체의 유동특성을 이해하려는 노력이 진행되고 있다. 즉, 자기유체역학(magnetohydrodynamic, MHD)의 적용 분야는 용융된 금속의 생산과정에서 결함을 최소화하는 주조 공정, 고품질의 철강생산을 위한 금속산업, MHD 펌프, MHD 발전기, 자동차산업 및 핵융합로의 플라즈마 등이 있으며, 이러한 분야에서 MHD 유동 특성을 예측하고 설계하는 과정은 이제 필수조건이 되고 있다(Kenjeres and Hanjalic, 2000; Kenjeres et al., 2004; Wilson et al., 2015; Shin and Kwak, 2016; Kenjeres, 2018).

자기장이 유동장에 부가되면 전류가 유도되고 이 전류와 자기장의 상호작용으로 인해 로렌츠 힘이 생성되며, 이 로렌츠 힘은 흐름을 구속하고 속도장을 변화시키게 된다. 또한, 자기장은 난류 유동을 재층류화하고 난류의 구조를 크게 변화시킬 수 있기 때문에 자기장에 노출된 난류의 마찰 특성 및 혼합 현상은 자기장이 없는 흐름에 비해 크게 달라질 수 있게 된다 Kenjeres et al., 2004; Wilson et al., 2015). 이러한 로렌츠 힘에 의한 유동 특성을 이해하고자, 본 연구에서는 완전 발달된 정사각 덕트 내 난류 MHD 유동을 레이놀즈 평균된 Navier-Stokes 방정식 모형으로 해석을 시도하였다.

공학적으로 많이 사용되고 있지만 자기장의 작용 방향에 상관없이 난류 MHD 유동 특성을 예측해 낼 수 있는 $k - ϵ$ 모형은 구현되어 있지 않은 상태이고, 이 모형 보다 발전된 형태의 레이놀즈응력 모형들은 자기장의 방향과는 무관하게 해를 예측해 주고 있지만 벽면 수직벡터를 직접 구해야 하는 것과 벽면을 인지하기 위해 난류 레이놀즈수와 같은 경험함수 등을 사용하고 있기 때문에, 유동 벽면이 복잡한 3차원 유동장에는 그 적용에 한계가 있다(Wilson et al., 2015; Shin and Kwak, 2016).

본 해석에서는 다양한 난류 유동에 적용되어 비교적 타당성 있는 유동 특성을 예측해 주고 있는 Thielen et al.(2004)의 타원-혼합 방정식에 의한 레이놀즈 응력 모형에 Kenjeres et al.(2004)이 제안한 로렌츠 힘 효과 항을 결합시켜 이 모형의 적용 범위를 확장하는 시도를 먼저 하였다. 이후 확장 모형을 정사각 덕트 내 난류 MHD 유동장에 적용하여 그 예측 해를 통해 기존 모형들이 안고 있는 문제점들을 해소할 수 있는지의 여부를 판단하게 된다.

2. 이론 해석

2.1 MHD 유동 관련 운동량 방정식

MHD 유동을 기술하는 방정식들은 Navier-Stokes 운동량 방정식, Maxwell 방정식 그리고 이동 매체(moving media)에 대한 옴(ohm)의 법칙 등으로 구성된다. 이 방정식들의 결합은 운동량 방정식에서 F^L=J×B로 정의되는 로렌츠 힘을 통해서 이루어진다. 여기서, J는 전류밀도, B는 부가된 자기장을 의미한다. 이동 매체에 대한 옴(Ohm)의 법칙 J=𝜎(-▽𝜱+U×B)을 적용하는 것에 의해 로렌츠 힘은 식 (1)과 같이 유도될 수 있다(Kenjeres and Hanjalic, 2000; Kenjeres et al., 2004; Kenjeres, 2018).

(1)

F_{i}^{L} = \frac{σ}{ρ} (- ϵ_{i j k} B_{k} \frac{\partial Φ}{\partial x_{j}} + U_{k} B_{i} B_{k} - U_{i} B_{k}^{2})

식 (1)에서 𝜎는 전기전도율(electrical conductivity), $ϵ_{i j k}$ 는 퍼뮤테이션 텐서이고, 전위(electric potential) 𝜱는 옴의 법칙과 전류보존 조건 ▽•J=0을 사용하는 것에 의해 식 (2)의 푸아송 방정식으로부터 구해진다.

(2)

\nabla^{2} Φ = \nabla \cdot (U \times B)

최종적으로 로렌츠 힘을 고려한 레이놀즈 평균된 운동을 나타내기 위한 운동량 방정식은 식 (3)과 같이 표현된다.

(3)

\frac{D U_{i}}{D t} = - \frac{1}{ρ} \frac{\partial P}{\partial x_{i}} + \frac{\partial}{\partial x_{j}} [ν (\frac{\partial U_{i}}{\partial x_{j}} + \frac{\partial U_{j}}{\partial x_{i}}) - \bar{u_{i} u_{j}}] + F_{i}^{L}

여기서, $\bar{u_{i} u_{j}}$ 는 레이놀즈응력 텐서이다.

2.2 난류 유동 해석을 위한 타원-혼합 모형

난류 MHD 유동장을 해석하기 위한 기본 난류 모형으로 본 연구에서는 타원-혼합 2차 모멘트 모형(elliptic blending second moment model, EBM)을 채택하였는데(Manceau and Hanjalic, 2002; Thielen et al., 2004), 이는 모형의 단순성과 편리성 외에 벽을 인지하기 위한 경험함수 등을 사용하지 않고도 벽 근방 난류 거동에 대한 물리적인 특성을 예측해 낼 수 있는 능력을 보유하고 있기 때문이다. 이 기본 모형에 로렌츠 힘 효과 $S_{i j}^{M}$ 을 포함한 레이놀즈응력 수송방정식은 식 (4)와 같이 정리된다(Kenjeres et al., 2004; Thielen et al., 2004).

(4)

\frac{D \bar{u_{i} u_{j}}}{D t} = \frac{\partial}{\partial x_{k}} ((ν δ_{k l} + C_{s} \bar{u_{k} u_{l}} T) \frac{\partial \bar{u_{i} u_{j}}}{\partial x_{l}}) + P_{i j} + Π_{i j}^{*} + S_{i j}^{M} - ϵ_{i j}

식 (4)에서 우변의 첫 번째 항은 점성 및 난류확산 항이고, 두 번째 항 $P_{i j}$ 는 식 (5)와 같이 난류 생성을 나타내는 항으로써 따로 모형화가 필요하지 않다.

(5)

P_{i j} = - \bar{u_{i} u_{k}} \frac{\partial U_{j}}{\partial x_{k}} - \bar{u_{j} u_{k}} \frac{\partial U_{i}}{\partial x_{k}}

또한 식 (4)에서 $Π_{i j}^{*}$ 는 속도-압력구배 상관관계, $ϵ_{i j}$ 는 소멸항을 나타내며 EBM에서 이들은 각각 식 (6)과 식 (7)로 모형화 된다.

(6)

Π_{i j}^{*} = (1 - ψ^{2}) Π_{i j}^{w} + ψ^{2} Π_{i j}^{h}

(7)

ϵ_{i j} = (1 - ψ^{2}) \frac{\bar{u_{i} u_{j}}}{k} ϵ + ψ^{2} \frac{2}{3} ϵ δ_{i j}

이 모형들에서 스칼라 𝜓를 구하기 위하여 식 (8)과 같은 1개의 타원 미분방정식을 풀기 때문에 모형의 타원성을 유지할 수 있게 된다(Manceau and Hanjalic, 2002; Thielen et al., 2004).

(8)

ψ - L^{2} \nabla^{2} ψ = 1

식 (6)에서 속도-압력구배항에 대한 벽 근방 모형 $Π_{i j}^{w}$ 은 식 (9)와 같이 모형화하여 레이놀즈응력의 벽면 경계조건을 만족하도록 하였다.

(9)

Π_{i j}^{w} = - 5 \frac{ϵ}{k} (\bar{u_{i} u_{k}} n_{j} n_{k} + \bar{u_{j} u_{k}} n_{i} n_{k} - \frac{1}{2} \bar{u_{k} u_{l}} n_{k} n_{l} (n_{i} n_{j} + δ_{i j}))

그러나 이 모형은 벽면에 수직한 벡터가 포함되는 단점을 갖고 있기 때문에 식 (10)과 같이 식 (8)의 타원 방정식으로부터 유동장 내에서 벽면에 수직한 방향을 나타내는 단위벡터를 구하였다(Manceau and Hanjalic, 2002; Thielen et al., 2004).

(10)

n = \frac{\nabla ψ}{‖ \nabla ψ ‖}

또한 식 (8)에서 길이 척도 $L$ 은 식 (11)의 Kolmogorov 길이 척도를 사용하였다.

(11)

L = C_{L} \max (\frac{k^{3 / 2}}{ϵ}, C_{η} \frac{ν^{3 / 4}}{ϵ^{1 / 4}})

식 (6)에서 균질 난류에 대한 속도-압력구배항 $Π_{i j}^{h}$ 은 일반적으로 사용되고 있는 Speziale et al.(1991)이 제안한 SSG모형을 확장하여 적용된다.

2.3 난류 모형에 로렌츠 힘 효과 항 부여

난류 MHD유동장에서 로렌츠 힘은 전위 변화와 자기장 그리고 자기장과 평균속도 간 상호작용을 거쳐 식 (1)을 통해 직접 유동장에 작용하게 된다. 또한, 레이놀즈응력 $\bar{u_{i} u_{j}}$ 수송방정식과 난류운동에너지 소멸률 𝜖 방정식을 통해서도 직접적으로 작용하며, 레이놀즈응력 방정식의 압력-속도구배 상관관계를 통해서는 간접적으로 유동장에 영향을 미치게 된다.

2.3.1 레이놀즈응력 방정식을 통한 로렌츠 힘의 직접적인 기여

식 (1)을 통해 로렌츠 힘이 유동장에 직접 작용하지만, 레이놀즈응력 $\bar{u_{i} u_{j}}$ 에 대한 수송방정식 안으로 요동 로렌츠 힘( $f_{i}^{L}$ )을 고려하면 식 (12)와 같은 추가 항이 생성되며 이 항을 통해서도 직접적으로 로렌츠 힘에 의한 효과가 유동장에 영향을 미치게 된다(Kenjeres et al., 2004). 이 추가 항은 식 (4)의 레이놀즈응력 방정식에서 로렌츠 힘에 의한 생성항 $S_{i j}^{M}$ 이다.

(12)

S_{i j}^{M} = \frac{σ}{ρ} (- ϵ_{i k l} B_{l} \bar{u_{j} \frac{\partial ϕ}{\partial x_{k}}} - ϵ_{j k l} B_{l} \bar{u_{i} \frac{\partial ϕ}{\partial x_{k}}}) + \frac{σ}{ρ} (B_{i} B_{k} \bar{u_{j} u_{k}} + B_{j} B_{k} \bar{u_{i} u_{k}} - 2 B_{k}^{2} \bar{u_{i} u_{j}})

여기서, $B_{i}$ 는 자기장의 세기이고 𝜙는 요동 전위(fluctuating electric potential)이다. 식 (12)에서 요동 속도-전기장 상관관계( $\bar{u_{i} \partial ϕ / \partial x_{j}}$ )는 완전 발달된 채널 MHD 유동에 대한 근사화를 통해 Kenjeres et al.(2004)이 제안한 식 (13)의 모형을 활용하였다.

(13)

\bar{u_{i} \frac{\partial ϕ}{\partial x_{j}}} \approx C_{λ} ϵ_{j k l} B_{l} \bar{u_{i} u_{k}}

따라서, $S_{i j}^{M}$ 에 대한 최종모형은 식 (14)와 같이 정리된다.

(14)

S_{i j}^{M} = (1 - C_{λ}) \frac{σ}{ρ} (B_{i} B_{k} \bar{u_{j} u_{k}} + B_{j} B_{k} \bar{u_{i} u_{k}} - 2 B_{k}^{2} \bar{u_{i} u_{j}})

여기서, 모형 계수 $C_{λ}$ 는 Kenjeres et al.(2004)을 따라 본 해석에서도 0.6으로 채택하였다.

2.3.2 난류운동에너지 소멸률 방정식을 통한 로렌츠 힘의 직접적인 기여

레이놀즈응력 방정식의 추가 항 $S_{i j}^{M}$ 과 유사하게 난류운동에너지 소멸률 𝜖 수송방정식에도 요동 로렌츠 힘( $f_{i}^{L}$ )을 고려하면 식 (15)와 같은 추가 항 $S_{ϵ}^{M}$ 이 생성되며, 이 항도 로렌츠 힘 효과가 유동장에 직접 영향을 미치는 요인이 된다(Kenjeres et al., 2004).

(15)

S_{ϵ}^{M} = \frac{1}{2} S_{i i}^{M} \frac{ϵ}{k} = \frac{ϵ}{k} \frac{σ}{ρ} (- ϵ_{i j k} B_{k} \bar{u_{i} \frac{\partial ϕ}{\partial x_{j}}} - 2 k B_{k}^{2} + B_{i} B_{k} \bar{u_{i} u_{k}})

Kenjeres et al.(2004)은 식 (15)의 $S_{ϵ}^{M}$ 에 대한 모형을 난류운동에너지 소멸률 방정식에 포함시키고 이 항에 대한 예측 분포를 DNS(direct numerical simulation) 자료에 근접시키기 위해 식 (16)과 같은 모형 계수 $C_{ϵ 4}$ 를 도입하였다.

(16)

S_{ϵ}^{M} = C_{ϵ 4} \frac{S_{\ddot{u}}^{M}}{2} \frac{ϵ}{k} = C_{ϵ 4} \frac{S_{i i}^{M}}{2 T}

또한 식 (16)에서 시간척도 $T$ 는 소멸률에 대한 벽면 인접효과를 고려하기 위하여 최초 제안된 $k / ϵ$ 를 대체한 것이며, 이 시간척도 $T$ 는 식 (17)과 같이 Kolmogorov 시간척도로 제한된다(Manceau and Hanjalic, 2002; Thielen et al., 2004).

(17)

T = \max (\frac{k}{ϵ}, C_{T} {(\frac{ν}{ϵ})}^{1 / 2})

식 (16)의 추가 항을 고려하여 최종 얻어지는 난류운동에너지 소멸률 방정식은 식 (18)과 같이 정리된다.

(18)

\frac{D ϵ}{D t} = \frac{\partial}{\partial x_{k}} ((ν δ_{k l} + C_{ϵ} \bar{u_{k} u_{l}}) \frac{\partial ϵ}{\partial x_{l}}) + C_{ϵ 1} \frac{P_{k k}}{2 T} - C_{ϵ 2} \frac{ϵ}{T} + C_{ϵ 3} \frac{S_{k k}^{M}}{2 T}

식 (18)에서 모형 계수 $C_{ϵ 1}$ 은 식 (19)와 같이 표현된다(Thielen et al., 2004).

(19)

C_{ϵ 1} = C_{ϵ 1}^{o} (1 + 0.03 (1 - ψ^{2}) \sqrt{\frac{k}{\bar{u_{i} u_{j}} n_{i} n_{j}}})

한편, 모형 계수 $C_{ϵ 3}$ 은 자기장이 흐름방향 및 채널에 수직인 방향으로 작용하는 채널 MHD 유동에 관한 DNS 자료(Noguchi et al., 1998)를 근간으로 하여 본 연구에서 1.75로 적용하였다.

2.3.3 속도-압력구배 상관관계를 통한 로렌츠 힘의 간접적인 기여

부가된 자기장의 간접적 기여는 요동 압력에 대한 푸아송 방정식에서 요동 로렌츠 힘으로 인한 속도-압력구배(velocity-pressure gradient) 상관관계에서 비롯된다. 균질 난류에서 요동 자기력으로 인한 속도-압력구배항 $Π_{i j}^{h}$ 에 대한 모형화를 위해 자기장이 없는 경우에 대한 SSG모형(Speziale et al., 1991)에 로렌츠 힘 효과를 고려한 식 (20)과 같은 모형을 본 연구에 도입하였다.

(20)

Π_{i j}^{h} = - (C_{1} + C_{1}^{*} \frac{P}{ϵ}) ϵ b_{i j} + C_{2} (b_{i k} b_{k j} - \frac{1}{3} b_{k l} b_{k l} δ_{i j}) + (C_{3} - C_{3}^{*} \sqrt{b_{k l} b_{k l}}) k S_{i j} + C_{4} k (b_{i k} S_{j k} + b_{j k} S_{i k} - \frac{2}{3} b_{l m} S_{l m} δ_{i j}) + C_{5} k (b_{i k} W_{j k} + b_{j k} W_{i k}) - C_{6 s} (S_{i j}^{M} - \frac{1}{3} S_{k k}^{M} δ_{i j})

여기서, $b_{i j}$ 는 레이놀즈응력의 비등방 텐서로서 식 (21)과 같이 표현되고, 평균 변형률 $S_{i j}$ 와 평균 와도(mean vorticity) 텐서 $W_{i j}$ 는 식 (22)로 주어진다.

(21)

b_{i j} = \frac{\bar{u_{i} u_{j}}}{2 k} - \frac{1}{3} δ_{i j}

(22)

S_{i j} = \frac{1}{2} (\frac{\partial U_{i}}{\partial x_{j}} + \frac{\partial U_{j}}{\partial x_{i}}), W_{i j} = \frac{1}{2} (\frac{\partial U_{i}}{\partial x_{j}} - \frac{\partial U_{j}}{\partial x_{i}})

또한 추가 항에 대한 모형 계수 $C_{6 s}$ 는 평판 MHD 유동에 대한 DNS 자료(Noguchi et al., 1998)와 비교하는 과정을 거쳐 본 연구에서 0.55로 선정하여 사용하였다.

2.4 로렌츠 힘 효과가 고려된 난류 모형의 모형 계수

자기장에 의한 로렌츠 힘 효과를 반영할 수 있는, 이른바 확장된 형태의 타원-혼합 레이놀즈응력 난류 모형의 모형 계수들은 Table 1과 같이 요약된다.

Table 1.

Model coefficients considering MHD effect

$C_{s}$	$C_{λ}$	$C_{ϵ}$	$C_{ϵ 1}^{o}$	$C_{ϵ 2}$	$C_{ϵ 3}$	$C_{1}$	$C_{1}^{*}$	$C_{2}$	$C_{3}$	$C_{3}^{*}$	$C_{4}$	$C_{5}$	$C_{6 s}$	$C_{T}$	$C_{L}$	$C_{η}$
0.21	0.6	0.18	1.44	1.83	1.75	3.4	1.8	4.2	0.8	1.3	1.25	0.4	0.55	6.0	0.161	80.0

3. 수치해석

Fig. 1은 본 해석에서 수행된 완전 발달된 정사각 덕트 난류 MHD 유동장이고, 자기장( $B_{o}$ )은 덕트 하단에서 상단으로 작용하는 경우를 고려하였으며, 자기장의 세기는 하트만 수 $H a (\equiv D_{h} B_{o} \sqrt{σ / (ρ ν)})$ 로 나타내었다.

https://cdn.apub.kr/journalsite/sites/ksds/2025-018-02/N0240180201/images/ksds_2025_182_1_F1.jpg

Fig. 1.

MHD flow field in a square duct

본 해석에 사용한 전산 해석 프로그램은 주유동 방향의 대류와 확산을 배제한 이산화 방정식을 바탕으로 자체 작성(in-house code) 하여 사용하였으며, 이러한 방법의 적용은 전진 해법(forward marching procedure)을 한 단면에서 반복 계산법으로 연속방정식이 만족될 때까지 종속 변수 값들을 다시 계산하는 방법이 된다. 유동 단면에서 압력, 레이놀즈응력 그리고 전위와 같은 스칼라 변수는 격자점에 위치시키고, 단면 속도 성분은 스칼라 격자점 사이에 위치시키는 이른바 엇갈림(staggered) 격자계를 활용하였다. 그 잇점으로는 속도가 스칼라 검사체적의 경계에 있으므로 대류항의 계산시 필요한 위치에 존재하게 되며, 또한 압력 격자점이 속도 격자점의 양쪽에 위치하게 되므로 속도에 대한 압력구배를 정확히 구할 수 있게 된다(Shin et al., 2008). 계산 격자는 다양한 격자수에 대한 사전 해석을 통해, 격자수가 예측 결과에 영향을 미치지 않는 100 × 100 격자를 선정하였으며 벽면에서 첫 번째 위치한 격자는 $y^{+}$ ≤0.5가 되게 하였다.

수송방정식의 차분화는 유동 단면의 $U$ , $V$ 방정식의 대류항에 대하여 QUICK 해법을 기타 방정식에 대해서는 HYBRID 해법을 적용하였으며, 압력 수정은 SIMPLE 알고리즘을 사용하였다(Shin et al., 2008). 경계조건으로는 벽면에서 평균속도 및 레이놀즈응력들에 대해서는 No-slip 조건을 𝜖에 대해서는 등방성 소멸률이 0이 되는 조건을 적용하였다. 전위 𝜱에 대해서는 벽면에서 수직 구배가 0이 되도록 하였다(Dey and Zikanov, 2012). 출구 조건에 대해서는 완전 발달된 유동이라고 보아 Newmann 조건을 적용하였다.

4. 결과 및 고찰

Fig. 2는 각각 수평 이등분선과 수직 이등분선을 따라 예측된 흐름방향 무차원 평균속도를 보여주고 있으며, 전반적으로 DNS 결과(Chaudhary et al., 2010)를 잘 따르고 있다. 자기장이 없는 경우에 비해 수직 이등분선을 따른 평균속도는 수평 이등분선을 따른 경우보다 중심부에서 더 평편한 분포로 예측된 것을 볼 수 있는데, 이는 $H a$ =21.2에서 자기장이 하단에서 상단 방향으로 작용할 때 측벽 방향보다는 상 ‧ 하단 방향에서 난류를 더 억제하기 때문인 것으로 판단된다. 즉, 덕트 하단에서 상단 방향으로 자기장 $B_{o}$ 가 작용하는 경우 완전 발달된 덕트 유동에서 식 (1)에 의한 로렌츠 힘은 식 (23)과 같이 전개되므로 자기장 방향에 수직인 속도 $U$ , $W$ 에 로렌츠 힘이 직접적으로 작용하며, 이 로렌츠 힘에 의해 유동속도가 변형된다는 것을 알 수 있다.

(23)

F_{x}^{L} = - \frac{σ}{ρ} B_{o}^{2} U, F_{y}^{L} = 0, F_{z}^{L} = - \frac{σ}{ρ} B_{o}^{2} (W + \frac{1}{B_{o}} \frac{\partial Φ}{\partial x})

특히, 주유동 방향의 로렌츠 힘에는 전위에 대한 변화 효과가 포함되어 있는데, 이 효과가 유동장의 수평 방향보다 수직 방향으로 더 강하게 작용하여 수직 이등분선을 따른 주유동 속도를 보다 평편하게 예측하는 것으로 여겨진다. 로렌츠 힘 효과가 커질수록 유동장 중심부에서의 속도분포가 더 평편해지는 경향을 자기장의 세기가 보다 강해진 $H a$ =30의 예측 결과에서 볼 수 있다.

https://cdn.apub.kr/journalsite/sites/ksds/2025-018-02/N0240180201/images/ksds_2025_182_1_F2.jpg

Fig. 2.

Streamwise mean velocity profiles along (a) horizontal bisector, (b) vertical bisector

$H a$ =30인 경우는 DNS 자료가 제시되어 있지 않아 그 해석 결과의 정확성을 판단할 수 없으나, 자기장의 세기에 따라 예측되어야 할 속도 분포를 비교적 타당성 있게 보여주고 있다.

Fig. 3은 수평 및 수직 이등분선을 따라 예측된 흐름방향 난류강도( $w$ ) 분포를 나타낸다. 이 예측 결과로부터 자기장이 난류를 억제하는 효과를 덕트 중심부와 수직 이등분선을 따라 하단 벽 근처에서 명확하게 볼 수 있다. 즉, 주유동 방향으로 로렌츠 힘의 직접적인 작용에 의한 속도 변형과 로렌츠 힘 생성항 $S_{33}^{M}$ 이 레이놀즈응력 방정식에 직접 영향을 주어 흐름방향 난류강도를 억제하기 때문이다. 특히, 로렌츠 힘에 기인한 평균속도 변형이 유동장의 수평 방향보다는 수직 방향의 하단 벽면 인접부에서 강하게 발생하여 이 부근의 난류강도를 상당히 변화시키는 것으로 판단된다.

https://cdn.apub.kr/journalsite/sites/ksds/2025-018-02/N0240180201/images/ksds_2025_182_1_F3.jpg

Fig. 3.

Streamwise turbulence intensity profiles along (a) horizontal bisector, (b) vertical bisector

$H a$ =21.2에서 하단 벽 근처의 흐름방향 난류강도( $w$ ) 분포는 피크 위치가 덕트 중심부로 약간 이동하면서 자기장이 없는 경우에 비해 약 30% 정도 더 억제되고 있다. 예측 결과들은 두 이등분선을 따라 난류가 동시에 억제되는 모습을 보여주고 있으며, 이는 DNS 결과(Chaudhary et al., 2010)와 비교적 잘 일치하고 있음을 나타내주고 있다. 또한, 자기장이 강해질수록 수직 방향 전체 및 유동 중심부에서 흐름방향 난류강도가 상당히 감소하는 것을 $H a$ =30의 결과로부터 알 수 있다. Fig. 4는 이등분선을 따라 예측된 수평 방향 난류강도( $u$ ) 분포와 수직 방향 난류강도( $v$ ) 분포를 보여준다. 이 경우에는 로렌츠 힘에 의한 속도 변형이 $\bar{u u}$ 및 $\bar{v v}$ 레이놀즈응력 방정식 모두에 영향을 주고, 로렌츠 생성항은 $\bar{u u}$ 레이놀즈응력 방정식에 직접 작용하게 되며, 이 영향들의 상호작용을 통해 수평 및 수직 방향 난류강도 모두를 억제하게 된다. 자기장이 없는 경우에 비해 $H a$ =21.2에서 예측된 본 해석 결과들은 DNS(Chaudhary et al., 2010)와 다소 차이가 있으나, 그 변화 경향은 DNS 분포를 그대로 추종하고 있다.

https://cdn.apub.kr/journalsite/sites/ksds/2025-018-02/N0240180201/images/ksds_2025_182_1_F4.jpg

Fig. 4.

x-direction (a) and y-direction (b) turbulence intensity profiles along horizontal and vertical bisectors

자기장이 부가된 유동장에서 난류강도들( $w$ , $u$ , $v$ )이 DNS 자료보다 다소 크게 예측되고 있는데, 이는 $\bar{u_{i} u_{j}}$ 방정식 모형에 도입한 로렌츠 힘에 의한 생성항 $S_{i j}^{M}$ 모형의 부정확성으로부터 기인하는 것과 본 해석에서 기본 모형으로 채택한 타원-혼합 난류 모형이 현재 계산된 정사각 덕트의 평균 레이놀즈수보다 다소 높은 채널유동장에서 난류 모형 계수들이 조율된 것도 또 다른 원인이 되는 것으로 고려된다(Kenjeres et al., 2004; Thielen et al., 2004). 두 번째 경우는 $H a$ =0인 유동 중심부에서 DNS 자료 보다 크게 예측되는 것으로부터 유추할 수 있다.

Fig. 5(a)는 수평 및 수직 이등분선을 따라 예측된 레이놀즈 전단응력( $\bar{u w}$ 및 $\bar{v w}$ ) 분포를 나타낸다. 로렌츠 힘이 주 방향 난류강도들( $w$ , $u$ , $v$ )에 미치는 영향과 유사하게 레이놀즈 전단응력( $\bar{v w}$ )도 수직 이등분선을 따라 하단 벽 근처에서 상당히 강하게 억제되고 있으나, 측벽 근처에서 $\bar{u w}$ 의 억제는 다소 약하게 작용하는 것을 보여준다.

Fig. 5(b)는 수직 및 수평 이등분선을 따르는 난류운동에너지 분포인데, 유동 방향 난류강도들( $w$ , $u$ , $v$ )의 예측 결과들로부터 이미 인지할 수 있듯이 수평 방향에 비해 수직 방향 하단의 난류운동에너지가 상당히 억제된 것을 보여준다. 난류운동에너지는 각 주유동 방향의 레이놀즈 응력들을 모두 더해서 2로 나눈 값 $k = (\bar{u u} + \bar{v v} + \bar{w w}) / 2$ 으로 정의된다. Fig. 6은 유동 단면을 흐르는 속도벡터 선도를 나타내는데, 자기장이 없는 경우의 Fig. 6(a)에는 유동장 모서리를 향하는 전형적인 8개의 와류가 유동장 전체에 발생하는 것을 보여주고 있다. 이에 반해 유동장의 $y$ 방향(하단에서 상단으로)으로 자기장이 작용하는 경우인 Fig. 6(b)에는 주유동 및 2차 유동의 속도변형에 의해 보다 강한 4개의 와류가 수평 및 수직 이등분선으로 이루어진 각 사분면 내에 형성되고, 나머지 4개의 와류는 양 측벽 모서리 근처에서 상당히 위축된 형태로 유지되고 있는 것을 볼 수 있다.

https://cdn.apub.kr/journalsite/sites/ksds/2025-018-02/N0240180201/images/ksds_2025_182_1_F5.jpg

Fig. 5.

(a) Reynolds shear stress -uw profile along horizontal bisector and Reynolds shear stress vw profile along vertical bisector, (b) Turbulent kinetic energy $k$ profile along horizontal bisector and vertical bisector

https://cdn.apub.kr/journalsite/sites/ksds/2025-018-02/N0240180201/images/ksds_2025_182_1_F6.jpg

Fig. 6.

Cross-stream vector fields for (a) Ha = 0, (b) Ha = 30

Fig. 7은 유동 단면에서 $y$ 방향으로 $H a$ =21.2의 자기장이 작용할 때의 전위(electric potential) 𝜱의 분포도이며 양 측벽 인접부에서 전위의 변화가 있는 것을 보여주는데, 이들 관계는 식 (24)의 전위에 대한 푸아송 방정식으로부터 얻어진다.

(24)

\frac{\partial^{2} Φ}{\partial x^{2}} + \frac{\partial^{2} Φ}{\partial y^{2}} = - B_{y} \frac{\partial W}{\partial x}

전위가 자기장의 작용 방향에 따라 변화를 보이는 것은 정사각 덕트 단면에서 $y$ 방향으로만 자기장( $B_{y} = B_{o}$ )이 작용하는 경우 전위는 $x$ 방향에 따른 주유동속도( $W$ )의 변화에 영향을 받기 때문에 Fig. 7과 같이 양 측벽으로의 변화가 생긴다. 또한, 이 전위는 주유동 방향의 로렌츠 힘에 상당한 영향을 주어 유동장의 평균속도 변형을 통해 난류를 억제시키는 요인이 된다.

https://cdn.apub.kr/journalsite/sites/ksds/2025-018-02/N0240180201/images/ksds_2025_182_1_F7.jpg

Fig. 7.

Contours of electric potential 𝜱 for Ha = 21.2

5. 결 론

자기장이 덕트 하단에서 상단 방향으로 작용하는 완전 발달된 정사각 덕트 내 난류 MHD 유동장에 확장된 형태의 타원-혼합 레이놀즈응력 방정식 난류 모형을 적용하여 다음과 같은 결론을 얻었다.

(1) 평행평판 난류 MHD 유동에서 조율된 로렌츠 힘 효과 관련 난류 모형 계수들을 사용하여 예측된 정사각 덕트 내 유동 특성들은 전반적으로 DNS 분포(Chaudhary et al., 2010)에 상당히 근접한 경향을 보여주고 있어, 본 유동장 해석에 적용한 MHD 유동 난류 모형의 타당성을 입증해 주고 있다.

(2) 자기장이 없는 경우에 비해 수직 이등분선을 따른 평균속도는 수평 이등분선을 따른 경우보다 중심부에서 더 평편한 분포로 예측된 것을 볼 수 있는데, 이는 흐름방향의 로렌츠 힘에 포함된 전위에 대한 변화 효과가 유동장의 수평 방향보다 수직 방향으로 더 강하게 작용하기 때문인 것으로 판단된다.

(3) 유동 방향으로 작용하는 로렌츠 힘과 변형된 속도들의 상호작용 그리고 레이놀즈응력 방정식에 직접 도입된 로렌츠 힘 생성항의 영향으로 난류 전단응력 및 난류강도 분포들이 억제되고 있지만, 전위 변화에 대한 영향으로 인해 수직 방향 분포들이 보다 강하게 억제되고 있다.

(4) 본 해석에 도입된 확장 형태의 타원-혼합 레이놀즈응력 난류 모형은 벽 경계가 복잡한 유동장에 대해서 벽을 인지하기 위한 경험함수 등을 사용하지 않고도 우수한 예측이 가능하기 때문에, 공학적 응용이 필요한 각종 난류 MHD유동장 해석에 대한 활용의 폭을 보다 넓혔다.

References

Chaudhary, R., S. P. Vanka, and B. G. Thomas. (2010). Direct Numerical Simulation of Magnetic Field Effects on Turbulent Flow in a Square Duct. Physics of Fluids. 22: 075102.

10.1063/1.3456724

Dey, P. K. and O. Zikanov. (2012). Turbulence and Transport of Passive Scalar in Magnetohydrodynamic Channel Flows with Different Orientation of Magnetic Field. International Journal of Heat and Fluid Flow. 36: 101-117.

10.1016/j.ijheatfluidflow.2012.01.005

Kenjeres, S. (2018). On Modeling and Eddy-Resolving Simulations of Flow, Turbulence, Mixing and Heat Transfer of Electrically Conducting and Magnetizing Fluids: A Review. International Journal of Heat and Fluid Flow. 73: 270-279.

10.1016/j.ijheatfluidflow.2018.09.003

Kenjeres, S. and K. Hanjalic. (2000). On the Implementation of the Effect of Lorentz Force in Turbulence Closure Models. International Journal of Heat and Fluid Flow. 21(3): 329-337.

10.1016/S0142-727X(00)00017-5

Kenjeres, S., K. Hanjalic, and D. Bal. (2004). A Direct-Numerical-Simulation-Based Second-Moment Closure for Turbulent Magnetohydrodynamic Flows. Physics of Fluids. 16(5): 1229-1241.

10.1063/1.1649335

Manceau, R. and K. Hanjalic. (2002). Elliptic Blending Model: A New Near Wall Reynolds Stress Turbulence Closure. Physics of Fluids. 14(2): 744-754.

10.1063/1.1432693

Noguchi, H., Y. Ohtsubo, and N. Kasagi. (1998). DNS Database of Turbulence and Heat Transfer. https://www.thtlab.jp/

Shin, J. K., J. S. An, Y. D. Choi, Y. C. Kim, and M. S. Kim. (2008). Elliptic Relaxation Second Moment Closure for the Turbulent Heat Fluxes. Journal of Turbulence. 9(3): 1-29.

10.1080/14685240701823101

Shin, Jong Keun and Dong Kurl Kwak. (2016). Evaluation and Improvement of Turbulence Model for Prediction of Turbulent Magnetohydrodynamic Flows. Journal of the Korean Society of Mechanical Technology. 18(5): 707-714.

10.17958/ksmt.18.5.201610.707

Speziale, C. G., S. Sarkar, and T. B. Gatski. (1991). Modeling of the Pressure-Strain Correlation Tensor: An Invariant Dynamical Systems Approach. Journal of Fluid Mechanics. 227: 245-272.

10.1017/S0022112091000101

Thielen, L., K. Hanjalic, H. Jonker, and R. Manceau. (2004). Prediction of Flow and Heat Transfer in Multiple Impinging Jets with an Elliptic-Blending Second-Moment Closure. International Journal of Heat and Mass Transfer. 48(8): 1683-1598.

10.1016/j.ijheatmasstransfer.2004.10.025

Wilson, D. R., T. J. Craft, and H. Iacovides. (2015). Application of Reynolds Stress Transport Turbulence Closure Models to Flows Affected by Lorentz and Buoyancy Forces. International Journal of Heat and Fluid Flow. 55: 180-197.

10.1016/j.ijheatfluidflow.2015.06.007

Korean References Translated from the English

신종근, 곽동걸 (2016). 난류 자기유체역학 유동 예측을 위한 난류모형의 평가와 개선. 한국기계기술학회지. 18(5): 707-714.

10.17958/ksmt.18.5.201610.707

Journal of Korean Society of Disaster and Security ISSN:2466-1147(Print) 2508-285X(Online) 한국방재안전학회 논문집

Preview

Prediction of Turbulent Flows in a Square Duct Subjected to Lorentz Force

ABSTRACT

MAIN

(1)

(2)

(3)

(4)

(5)

(6)

(7)

(8)

(9)

(10)

(11)

(12)

(13)

(14)

(15)

(16)

(17)

(18)

(19)

(20)

(21)

(22)

Table 1.

Model coefficients considering MHD effect

Fig. 1.

MHD flow field in a square duct

(23)

Fig. 2.

Streamwise mean velocity profiles along (a) horizontal bisector, (b) vertical bisector

Fig. 3.

Streamwise turbulence intensity profiles along (a) horizontal bisector, (b) vertical bisector

Fig. 4.

x-direction (a) and y-direction (b) turbulence intensity profiles along horizontal and vertical bisectors

Fig. 5.

(a) Reynolds shear stress -uw profile along horizontal bisector and Reynolds shear stress vw profile along vertical bisector, (b) Turbulent kinetic energy k profile along horizontal bisector and vertical bisector

Fig. 6.

Cross-stream vector fields for (a) Ha = 0, (b) Ha = 30

(24)

Fig. 7.

Contours of electric potential 𝜱 for Ha = 21.2

References

Korean References Translated from the English

(a) Reynolds shear stress -uw profile along horizontal bisector and Reynolds shear stress vw profile along vertical bisector, (b) Turbulent kinetic energy $k$ profile along horizontal bisector and vertical bisector