*2.3. Turbulence Models*

Three two-phase flow versions of Reynolds-averaged Navier–Stokes (RANS) turbulence models were used in the simulations: a *k* − *ε* model [9,16], a *k* − *ω*2006 model [15], and a *k* − *ε* turbulence model written in terms of the specific dissipation rate of TKE *ω*, denoted hereafter as the modified *k* − *ε* model. The last model was only used to compare the *k* − *ε* and *k* − *ω*2006 behaviors.

In the framework of two-equation RANS turbulence models, transport equations for the dissipation rate and for the TKE need to be solved to compute the turbulent viscosity *νf t* . The general expression for the TKE transport equation is given by:

$$\begin{split} \frac{\partial k}{\partial t} + \boldsymbol{\nu}\_{j}^{f} \frac{\partial k}{\partial \mathbf{x}\_{j}} &= \frac{R\_{ij}}{\rho^{f}} \frac{\partial \boldsymbol{u}\_{i}^{f}}{\partial \mathbf{x}\_{j}} + \frac{\partial}{\partial \mathbf{x}\_{j}} \left[ \left( \boldsymbol{\nu}^{f} + \sigma\_{k} \boldsymbol{\nu}\_{t}^{f} \right) \frac{\partial k}{\partial \mathbf{x}\_{j}} \right] - \varepsilon \\ &- \frac{2K(1 - a)\phi k}{\rho^{f}} - \frac{\boldsymbol{\nu}\_{t}^{f}}{\sigma\_{\varepsilon}(1 - \phi)} \frac{\partial \phi}{\partial \mathbf{x}\_{j}} (s - 1) \boldsymbol{g}\_{i} \end{split} \tag{25}$$

with *Rfij*the Reynolds stress tensor and *σk* an empirical coefficient.

The transport equation for the dissipation rate and the expression of the turbulent viscosity differ for the different turbulence models.

## 2.3.1. *k* − *ε* Model

For the *k* − *ε* model, the turbulent viscosity *νf t*is calculated as:

$$\nu\_t^f = \mathbb{C}\_{\mu} \frac{k^2}{\varepsilon},\tag{26}$$

and the following transport equation for the dissipation rate *ε* is solved:

$$\begin{split} \frac{\partial \varepsilon}{\partial t} + \boldsymbol{\mu}\_{j}^{f} \frac{\partial \varepsilon}{\partial \boldsymbol{x}\_{j}} &= \mathbb{C}\_{1\varepsilon} \frac{\varepsilon}{k} \frac{R\_{ij}}{\rho^{f}} \frac{\partial \boldsymbol{u}\_{i}^{f}}{\partial \boldsymbol{x}\_{j}} + \frac{\partial}{\partial \boldsymbol{x}\_{j}} \left[ \left( \boldsymbol{\nu}^{f} + \sigma\_{\varepsilon} \boldsymbol{\nu}\_{i}^{f} \right) \frac{\partial \varepsilon}{\partial \boldsymbol{x}\_{j}} \right] - \mathbb{C}\_{2\varepsilon} \frac{\varepsilon^{2}}{k} \\ &- \mathbb{C}\_{3\varepsilon} \frac{\varepsilon}{k} \frac{2K(1-\alpha)\phi k}{\rho^{f}} - \mathbb{C}\_{4\varepsilon} \frac{\varepsilon}{k} \frac{\boldsymbol{\nu}\_{i}^{f}}{\sigma\_{\varepsilon}(1-\phi)} \frac{\partial \phi}{\partial \boldsymbol{x}\_{j}} (\boldsymbol{s} - 1) \boldsymbol{g}\_{i}. \end{split} \tag{27}$$

The values of the empirical coefficients *σk*, *σε*, *C*1*ε*, *C*2*ε*, *C*3*<sup>ε</sup>*, *C*4*ε*, and *Cμ* are listed in Table 1.

**Table 1.** Empirical coefficients for the *k* − *ε* turbulence model from Chauchat et al. (2017) [17].


## 2.3.2. *k* − *ω*2006 Model

The dissipation rate *ε* can be expressed in term of specific dissipation rate *ω* following the expression *ε* = *<sup>C</sup>μk<sup>ω</sup>*. The *k* − *ω*2006 turbulence model uses *ω* and the norm of the deviatoric part of the strain rate tensor || *Sf* || to compute the eddy viscosity:

$$\nu\_t^f = \frac{k}{\max\left[\omega\_\nu \text{ \(C\_{lim} \stackrel{\parallel}{\mid} \frac{\mathbb{S}^f \parallel \text{\(}}{\cdot}\)}{\sqrt{\mathbb{C}\_\mu}}\right]} \tag{28}$$

Compared with the *k* − *ε* or the standard *k* − *ω* model, a stress limiter is incorporated and adjusted by the coefficient *Clim* [14].

The transport equation for the specific dissipation rate *ω* reads:

$$\begin{split} \frac{\partial \boldsymbol{\omega}}{\partial t} + \boldsymbol{u}\_{j}^{f} \frac{\partial \boldsymbol{\omega}}{\partial \mathbf{x}\_{j}} &= \mathbb{C}\_{1\omega} \frac{\omega}{k} \frac{\text{R}\_{ij}}{\rho^{f}} \frac{\partial \boldsymbol{u}\_{i}^{f}}{\partial \mathbf{x}\_{j}} + \frac{\partial}{\partial \mathbf{x}\_{j}} \left[ \left( \boldsymbol{\nu}^{f} + \sigma\_{\omega} \boldsymbol{\nu}\_{t}^{f} \right) \frac{\partial \boldsymbol{\omega}}{\partial \mathbf{x}\_{j}} \right] - \mathbb{C}\_{2\omega} \boldsymbol{\omega}^{2} - \mathbb{C}\_{3\omega} \boldsymbol{\omega} \frac{2\mathcal{K}(1-\boldsymbol{a})\boldsymbol{\phi}}{\rho^{f}} \\ &- \mathbb{C}\_{4\omega} \frac{\omega}{k} \frac{\boldsymbol{\nu}\_{t}^{f}}{\sigma\_{\boldsymbol{\epsilon}}(1-\boldsymbol{\phi})} \frac{\partial \boldsymbol{\phi}}{\partial \mathbf{x}\_{j}} (\boldsymbol{s} - 1) \boldsymbol{g}\_{i} + \sigma\_{d} \frac{1}{\omega} \frac{\partial \boldsymbol{k}}{\partial \mathbf{x}\_{j}} \frac{\partial \boldsymbol{\omega}}{\partial \mathbf{x}\_{j}}. \end{split} \tag{29}$$

The empirical coefficients for this turbulence model are presented in Table 2, and the coefficient before the cross-diffusion term *σd* is given by:

$$\sigma\_d = \begin{cases} 0 & \text{for} \quad \frac{\partial k}{\partial x\_j} \frac{\partial \omega}{\partial x\_j} < 0 \\ 1 & \text{for} \quad \frac{\partial k}{\partial x\_j} \frac{\partial \omega}{\partial x\_j} \ge 0. \end{cases} \tag{30}$$

**Table 2.** Empirical coefficients for the *k* − *ω*2006 turbulence model.


2.3.3. Modified *k* − *ε* Model

The modified *k* − *ε* model is obtained by substituting the dissipation rate *ε* by *<sup>C</sup>μk<sup>ω</sup>* in Equations (26) and (27). The new expression for the eddy viscosity is given by:

$$\nu\_t^f = \frac{k}{\omega} \tag{31}$$

and the new dissipation rate transport equation is the same as Equation (29) from the *k* − *ω*2006 model with different coefficients (see the coefficients in Table 3). The coefficient *σd* in front of the cross-diffusion term allows switching from the *k* − *ε* and the *k* − *ω*2006 models and investigating the sensitivity of the model to this cross-diffusion term.

$$\sigma\_d = \begin{cases} 0 & \text{for} \quad \frac{\partial k}{\partial x\_j} \frac{\partial \omega}{\partial x\_j} < 0 \\\\ 1.712 & \text{for} \quad \frac{\partial k}{\partial x\_j} \frac{\partial \omega}{\partial x\_j} \ge 0. \end{cases} \tag{32}$$

**Table 3.** Empirical coefficients for the modified *k* − *ε* turbulence model.

