*2.2. Governing Equations and Turbulence Model*

To analyze the three-dimensional compressible fluid flow in a gas turbine under unsteady-state conditions, we used the continuity equation, momentum equation, and energy equation, which can be expressed as shown below:

Continuity equation:

$$\frac{\partial \rho}{\partial t} + \frac{\partial}{\partial x\_i} (\rho u\_i) = 0.\tag{1}$$

Momentum equation:

$$\frac{\partial}{\partial t}(\rho u\_i) + \frac{\partial \{\rho u\_i u\_j\}}{\partial \mathbf{x}\_i} = -\frac{\partial \mathbf{P}}{\partial \mathbf{x}\_i} + \frac{\partial}{\partial \mathbf{x}\_j} \left[ \mu \left( \frac{\partial u\_i}{\partial \mathbf{x}\_j} + \frac{\partial u\_j}{\partial \mathbf{x}\_i} - \frac{2}{3} \delta\_{ij} \frac{\partial u\_k}{\partial \mathbf{x}\_k} \right) \right] + \frac{\partial \left( -\rho \left| u\_i' u\_j' \right| \right)}{\partial \mathbf{x}\_j}.\tag{2}$$

Energy equation:

$$\frac{\partial}{\partial t}(\rho E) + \frac{\partial}{\partial \mathbf{x}\_j} \{u\_j(\rho E + P)\} = \frac{\partial}{\partial \mathbf{x}\_j} \left[ (\mathbf{k}\_{eff}) \frac{\partial T}{\partial \mathbf{x}\_j} \right] + \frac{\partial}{\partial \mathbf{x}\_j} \left[ u\_i \mu\_{eff} \left( \frac{\partial u\_i}{\partial \mathbf{x}\_j} + \frac{\partial u\_j}{\partial \mathbf{x}\_i} - \frac{2}{3} \delta\_{ij} \frac{\partial u\_k}{\partial \mathbf{x}\_k} \right) \right] \tag{3}$$
 
$$i, j, k = 1, 2, \text{ and } 3,$$

where ρ is the fluid density, *u* is the fluid velocity, *P* is the fluid pressure, and μ is the fluid viscosity. In the energy equation, *E* is the specific internal energy, *keff* is the effective thermal conductivity, and μ*eff* is the effective dynamic viscosity. A finite volume method (FVM)-based commercial CFD software, ANSYS CFX [18], was used to solve the governing equations.

 =

Accurate prediction of the complex heat flow in a gas turbine requires an appropriate turbulence model for the simulation. Various methods, such as direct numerical simulation (DNS), large eddy simulation (LES), and the Reynolds-average Navier–Stokes (RANS) method, have been introduced for this purpose. DNS and LES can provide details of turbulence statistics but are high-cost methods [19–21]. The RANS method is usually used for CFD simulations due to its lower computational cost [22,23], especially for turbomachinery simulation. Moreover, previous studies have confirmed that the SST γ model and SST γ−θ model are the most suitable for analysis of transitional flows [24,25]. Furthermore, Choi and Ryu [9] concluded that results obtained using the *k* − ω SST γ turbulence model were in agreemen<sup>t</sup> with the corresponding experimental results [26]. Therefore, to accurately predict the complex fluid flow and heat transfer characteristics in a gas turbine, we used the *k* − ω SST γ turbulence model in this study.

The *k* − ω SST model was combined with free stream formulations and the *k* − ω formulation in the near wall using a blending function proposed by Menter [27,28]. The corresponding continuity, momentum, turbulence kinetic energy (*k*) equation, and eddy dissipation (ω) equations were formulated to express the *k* − ω baseline (BSL) model:

$$\frac{\partial}{\partial \mathbf{x}\_i}(\rho u\_i) = 0,\tag{4}$$

$$\frac{\partial \left(\rho u\_i u\_j\right)}{\partial \mathbf{x}\_i} = -\frac{\partial P^\*}{\partial \mathbf{x}\_i} + \frac{\partial}{\partial \mathbf{x}\_j} \left[\mu\_{eff} \left(\frac{\partial u\_i}{\partial \mathbf{x}\_j} + \frac{\partial u\_j}{\partial \mathbf{x}\_i}\right)\right] \tag{5}$$

$$\frac{\partial(\rho k)}{\partial t} + \frac{\partial(\rho u\_{\rangle} k)}{\partial \mathbf{x}\_{\rangle}} = P\_K - 0.09 \rho a k + \frac{\partial}{\partial \mathbf{x}\_{\rangle}} \left[ \left( \mu + \frac{\mu\_l}{\sigma\_k} \right) \frac{\partial k}{\partial \mathbf{x}\_{\rangle}} \right] \tag{6}$$

$$\frac{\partial(\rho\omega)}{\partial t} + \frac{\partial(\rho u\_f\omega)}{\partial \mathbf{x}\_j} = \frac{\mathcal{V}}{\nu\_t}\mathcal{P}\mathbf{x} - \beta\rho\omega^2 + \frac{\partial}{\partial \mathbf{x}\_j} \left[ \left(\mu + \frac{\mu\_t}{\sigma\_\omega} \right) \frac{\partial \mathbf{k}}{\partial \mathbf{x}\_j} \right] + 2\rho(1 - F\_1) \frac{1}{\sigma\_{\omega/2}} \frac{1}{\omega} \frac{\partial \mathbf{k}}{\partial \mathbf{x}\_j} \frac{\partial \omega}{\partial \mathbf{x}\_j},\tag{7}$$

where:

$$P^\* = P + \frac{2}{3} (\rho k + (\mu + \mu\_l) \frac{\partial u\_k}{\partial \mathbf{x}\_k}) \tag{8}$$

$$P\_K = \mu l \left(\frac{\partial u\_i}{\partial \mathbf{x}\_j} + \frac{\partial u\_j}{\partial \mathbf{x}\_i}\right) \frac{\partial u\_i}{\partial \mathbf{x}\_j} - \frac{2}{3} \frac{\partial u\_k}{\partial \mathbf{x}\_k} \left(3\mu\_l \frac{\partial u\_k}{\partial \mathbf{x}\_k} + \rho k\right), \tag{9}$$

μ*t*(*kg*/*ms*) is the turbulence viscosity calculated using the following equation:

$$\frac{1}{\rho}\mu\_t = \frac{0.31k}{\max(0.31\omega, SF\_2)},\tag{10}$$

The blending functions *F*1 and *F*2 are defined by the following variables:

$$\arg z\_1 = \min \left[ \max \left( \frac{\sqrt{k}}{0.09\omega y} \frac{500\nu}{y^2 \alpha} \right) \frac{4\rho \sigma\_{\omega,2} k}{CD\_{k\omega} y^2} \right] \tag{11}$$

$$\arg g\_2 = \max \left( 2 \frac{\sqrt{k}}{0.09\omega y'}, \frac{500\nu}{y^2 \omega} \right) \tag{12}$$

as follows:

$$F\_1 = \tanh(\arg\_{1}^{2}),\tag{13}$$

$$F\_2 = \tanh(\arg\_2^2),\tag{14}$$

where σ*k* and σω are the turbulent Prandtl numbers for *k* and ω, respectively. The formulations for these equations are expressed below:

$$
\sigma\_k = \frac{1}{F\_1/\sigma\_{k,1} + (1 - F\_1)/\sigma\_{k,2}} \, ^\prime \tag{15}
$$

$$
\sigma\_{\omega} = \frac{1}{F\_1/\sigma\_{\omega,1} + (1 - F\_1)/\sigma\_{\omega,2}},
\tag{16}
$$

$$\text{C}D\_{k\omega} = \max\left(2\rho \frac{1}{\sigma\_{\omega,2}} \frac{1}{\omega} \frac{\partial k}{\partial \mathbf{x}\_i} \frac{\partial \omega}{\partial \mathbf{x}\_j}, 10^{-10}\right). \tag{17}$$

Moreover, it is necessary to define the transport equation for the intermittency (γ) to obtain the complete expression for the *k* − ω SST γ turbulence model. The transport equation of γ can be defined as:

$$\frac{\partial(\rho \mathbf{y}\boldsymbol{\gamma})}{\partial t} + \frac{\partial(\rho u\_{\boldsymbol{\gamma}} \boldsymbol{\gamma})}{\partial \mathbf{x}\_{\boldsymbol{\gamma}}} = P\_{\boldsymbol{\gamma}1} + P\_{\boldsymbol{\gamma}2} - (E\_{\boldsymbol{\gamma}1} + E\_{\boldsymbol{\gamma}2}) + \frac{\partial}{\partial \mathbf{x}\_{\boldsymbol{\gamma}}} \Big[ \left( \mu + \frac{\mu\_{t}}{\sigma\_{\boldsymbol{\gamma}}} \right) \frac{\partial \mathbf{k}}{\partial \mathbf{x}\_{\boldsymbol{\gamma}}} \Big]. \tag{18}$$

The transition sources are defined as follows:

$$P\_{\uparrow\uparrow} = 2F\_{\text{length}} \rho S(\chi F\_{\text{orset}})^{c\_{\uparrow\uparrow}},\tag{19}$$

$$E\_{\uparrow\uparrow} = \gamma P\_{\uparrow\uparrow}.\tag{20}$$

The destruction/ relaminarization sources are defined as follows:

$$P\_{\gamma2} = 2c\_{\gamma1} \rho \Omega \gamma F\_{\text{turb}\prime} \tag{21}$$

$$E\_{\gamma 2} = \gamma P\_{\gamma 2} c\_{\gamma 2}.\tag{22}$$

The transition onset is controlled by the following functions:

$$R\_T = \frac{\rho k}{\mu \nu'}\tag{23}$$

$$F\_{onset1} = \frac{\rho y^2 S}{2.193 \mu Re\_{0c}} \,\text{\,\,\,}\tag{24}$$

$$F\_{onset2} = \min\{\max\{F\_{onset1}, F\_{onset1}4}, 2\},\tag{25}$$

$$F\_{onset3} = \max\left(1 - \left(\frac{R\_T}{2.5}\right)^3, 0\right),\tag{26}$$

*Fonset* = max(*Fonset*2 − *Fonset*3 , <sup>0</sup>), (27)

*Fturb* = *e*<sup>−</sup>(0.25*RT*)<sup>4</sup> , (28)

where *S* is the strain rate magnitude, *<sup>F</sup>*length is an empirical correlation, Ω is the vorticity magnitude, and *Re*θ*c* is the critical Reynolds number, at which the intermittency first starts to increase in the boundary layer. The other constant coefficients for the equations above are as follows [5,29]:

$$
\sigma\_{k,1} = 1.176, \ \sigma\_{k,2} = 1.0, \ c\_{\gamma^1} = 0.03, \ c\_{\gamma^2} = 50,
$$

$$
\sigma\_{\alpha;1} = 2.0, \ \sigma\_{\alpha;2} = 1.168, \ c\_{\gamma^3} = 0.5, \ \sigma\_{\gamma} = 1.
$$
