General Turbulence Guidelines

Turbulence is probably the most challenging area in fluid dynamics and the most limiting factor in accurate computer simulation of engineering flows. Typical CFD software solve the Navier Stokes and conservation equations, but as long as we are not resolving flows directly, some other mechanism have to be involved to capture Reynolds number. Most typical approach is called a Reynolds Averaged Navier Stokes (RANS). More advanced option – Eddy Simulation technique resolves the larger eddies in the flow – required for separation or large recirculating regions.

The RANS models apply a Reynolds decomposition technique to the Navier Stokes equations which breaks the velocity down into its mean and fluctuating components. This decomposition leaves us with one unknown value, which is termed the Reynolds Stress. We use turbulence models to resolve the Reynolds Stress and close the equation set. There are two ways we can go about resolving this, the first (and most common) is to use an isotropic value for the turbulent viscosity which is called the an Eddy Viscosity Model, the other way is to solve the Reynolds Stress Model (RSM) for the 6 separate Reynolds Stresses equations, which results in an anisotropic solution. This problem also also applies to Eddy Viscosity Models using an isotropic value may not be appropriate and could increase the numerical diffusion. Solving additional 6 Reynolds Stresses equations and dissipation will be more accurate, but run time will increase significantly (and stability decrease proportionally).

The common options are:

  • Spalart-Allmaras – One equation model for attached aerodynamic analysis.
  • k-epsilon – Two equation model for free shear and non-wall bounded flow behavior. Was the previous industrial standard.
  • k-omega – Two equation model for wall bounded flows, not commonly used.
  • SST (Shear Stress Transport) – Two equation model blending the free-stream advantages of the k-epsilon model with the wall bounded advantages of the k-omega model. This is the new industrial standard and should be the default choice for most applications.

Some additional features like Curvature Correction are very useful in tackling problems where previously RSM model were required. The Transition Models are a good solution depending on the problem  (they may give much better drag prediction as they will maintain laminar flow along a body and develop natural transition points, as well as calculate regions of bypass transition).

For example to simulate the flow in a cyclone separator with high swirl number Reynolds Stress model should be used. However SST model with the Curvature Correction model should be a good option as well (more stable).

The more advanced options for resolving turbulence:

  • Scale Adaptive Simulation
  • Detached Eddy Simulation
  • Large Eddy Simulation

The Scale Adaptive Simulation is the best point to start when stepping up from a RANS modelling analysis. The Scale Adaptive Simulation is based on an unsteady SST RANS model but calculates the local length scales and resolves accordingly. The advantage of this is that it can be run on a good quality RANS mesh without any additional meshing requirements. For some flow behaviours this has been shown to give very similar results to Large Eddy Simulation (LES) and Detached Eddy Simulation (DES). This is the best recommendation as a starting point if unsteady simulations is required (highly separated or recirculating flows).

The Detached Eddy Simulation resolves the boundary layer and the smaller eddies by using a RANS solution in those regions. This approach can produce quite good results, but this method is very sensitive for the mesh resolution. High quality mesh in the boundary layer and a transition region from the boundary layer to the free-stream is required to help the solver switch from the RANS solution in the boundary layer to an LES solution for the larger eddies. This adds complexity and has been known to cause erroneous results.

The Large Eddy Simulation is the most expensive of the three simulation options and requires a high resolution mesh to accurately resolve the eddies, especially in the boundary layer. Once the eddies are too small to be captured by the mesh, they fall back into a sub-grid model. As this sub-grid model is used only to resolve the smallest eddies which tend to have more universal properties, this yields very accurate results, but because of the large mesh requirements this is not a commonly used approach in industrial settings. LES an be used only in specified region of model this means that an unsteady RANS in the bulk of the flow is used and a full LES calculation in only activated in separated region – it is more practical for industrial use.

 

 

For calculations of turbulent flows, a variety of turbulence models are available, which could be categorized into three general groups, namely, RANS (Reynolds–Averaged Navier–Stokes), DES (Detached Eddy Simulation) and hybrid RANSLES, as well as LES (Large Eddy Simulation).

Turbulence Models

1. Spalart & Allmaras One-equation Model (1994)

This model solves one transport equation for a quantity \tilde{\nu} which is equivalent to the eddy viscosity. Since turbulence is characterised by two scales, e.g. velocity and length scales, and the model only solves for one propoerty additional information is needed. The Spalart & Allmaras model uses the wall distance, that would be active through the complete boundary layer, not only in the viscous sub layer. The detailed definition of the model is given in the reference, and is not repeted here. The model is integrated all the way to the wall which requires a good resolution normal to the wall (y^+\sim{1}).

2. Eddy Viscosity Two-Equation Models

In these turbulence models two additional transport equations for the turbulent kinetic energy {k} and some auxiliary quantity (\varepsilon, \omega or \tau) is solved. The models below are integrated all the way to the wall which requires a good resolution normal to the wall (y^+\sim{1}).

2.1. The Wilcox Standard k-\omega Turbulence Model (1988)
The additional equations can be put in the form.

\displaystyle U = \begin{pmatrix} \rho k\\ \rho \omega \end{pmatrix},\quad f_{I_i} = \begin{pmatrix} \rho k w_i\\ \rho \omega w_i \end{pmatrix},\quad -f_{V_i} = \begin{pmatrix} (\mu + \mu_t \sigma^*) k_{,i}\\ (\mu + \mu_t \sigma) \omega_{,i}\\ \end{pmatrix}, \ \ \ \ \ (1)

and the source term is defined as:

\displaystyle Q = \begin{pmatrix} P_k - \beta^*\rho k \omega \\ \displaystyle\gamma\frac{\omega}{k}P_k - \beta\rho \omega^2 \\ \end{pmatrix}, \ \ \ \ \ (2)

where {P_k} is the production of turbulent kinetic energy given by

\displaystyle P_k = \frac{\partial w_i}{\partial x_j}\left( \mu\left[ \frac{\partial w_i}{\partial x_j} + \frac{\partial w_j}{\partial x_i} - \frac{2}{3}(\nabla w)\delta_{ij}\right] -\frac{2}{3}\rho k \delta_{ij} \right). \ \ \ \ \ (3)

The eddy viscosity is given by

\displaystyle \mu_T = \rho\frac{k}{\omega}. \ \ \ \ \ (4)

Finally, the closure coefficients are as follow

\displaystyle \gamma=0.556, \quad \beta=0.075, \quad \beta^*=0.09, \quad \sigma=0.5, \quad \text{and} \quad \sigma^*=0.5. \ \ \ \ \ (5)

One of the major problems with the Wilcox {k-\omega} model is that the model has an unphysical free stream dependency.

2.3. Low Reynolds Number (LRN) {k-\omega} Models
Low Reynolds number (LRN) models are developed to account for the viscous and wall-damping effects. For linear two-equation models, this has usually been achieved by means of some empirical damping functions, which also help in many modelling variants to attain correct asymptotic properties when integrated to the wall surface. The LRN {k-\omega} model by Wilcox (1994) has been derived on the basis of the Wilcox standard {k-\omega} model. This suggests that the LRN Wilcox model returns to the standard version, when all the damping functions are set to a constant value of unity.

The LRN model by Peng (1997) is a modified version of the {k-\omega} model. It is cast in a transformed form from {k-\varepsilon} type model by some approximations and, consequently, has a cross diffusion term in the {\omega} equation, which helps to suppress the freestream sensitivity. The cross diffusion term takes the following form

\displaystyle C_\omega \frac{\mu_t}{k} \frac{\partial{k}}{\partial{x_k}}\frac{\partial{\omega}}{\partial{x_k}}. \ \ \ \ \ (6)

Unlike the following Kok TNT {k-\omega} model, this LRN model employs the cross diffusion term over the whole boundary layer. The model constants and damping functions are calibrated for internal flows characterized by recirculation with separation and reattachment.

2.4. The Kok TNT {k-\omega} Turbulence Model
The {k-\omega} model by Kok (2000} is derived for overcoming the unphysical free stream dependency. The equations are on the same form as the {k-\omega} model in (1) but with an added cross diffusion term. Also the model coefficients are modified. The cross diffusion term is written as

\displaystyle \sigma_d\frac{\rho}{\omega}\max\left(0, \frac{\partial{k}}{\partial{x_k}}\frac{\partial{\omega}}{\partial{x_k}}\right). \ \ \ \ \ (7)

2.5. The Menter BSL and SST {k-\omega} Turbulence Models
The {k-\omega} model by
Menter 1994 is also derived for overcoming the unphysical free stream dependency. Also here the equations are on the same form as the {k-\omega} model in (1) but with an added cross diffusion term as in (7). The model is a blending of the standard {k-\varepsilon} model in the outer part of the boundary layer and the Wilcox {k-\omega} model in the near wall part of the boundary layer. The {k-\varepsilon} model is transformed into {k-\omega} resulting in a cross diffusion term and modified coefficients which are blended according to a blending function.

In the SST model an additional modification is made, compared to the BSL model. A limiter on the eddy viscosity is added so that {\overline{uv}\le{a_1k}}, with {a_1=0.31}. This is in accordance with the Bradshaw assumption, and will reduce the turbulence level around stagnation regions and in adverse pressure gradient boundary layers compared to a standard eddy-viscosity model and the prediction of separated flows are improved.

3. Explicit Algebraic Reynolds Stress Models (EARSM)

The EARSM by Wallin & Johansson (2000) is a rational approximation of a full Reynolds stress transport model in the weak equilibrium limit where the Reynolds stress anisotropy may be considered constant in time and space. The Reynolds stress tensor is explicitly expressed in terms of the velocity gradient and the turbulence scales.

3.1 The Wallin & Johansson EARSM
The EARSM may be written in a common way where the anisotropy tensor {\mathbf{a}} is written in terms of the strain- and rotation rate tensors {\mathbf{S}} and {\mathbf{\Omega}} as

\displaystyle \mathbf{a} = \sum_{\lambda = 1}^{10} \beta_{\lambda} \mathbf{T}^{(\lambda)}, \ \ \ \ \ (8)

where the {\beta} coefficients are functions of the five invariants of {\mathbf{S}} and {\mathbf{\Omega}}. The {\mathbf{T}}‘s are

\displaystyle \begin{aligned} \mathbf{T}^{(1)} & = \mathbf{S},\\ \mathbf{T}^{(2)} & = \mathbf{S}^2 - \frac{1}{3}\textit{II}_S\mathbf{I},\\ \mathbf{T}^{(3)} & = \mathbf{\Omega}^2 - \frac{1}{3}\textit{II}_{\Omega}\mathbf{I},\\ \mathbf{T}^{(4)} & = \mathbf{S}\mathbf{\Omega} - \mathbf{\Omega}\mathbf{S},\\ \mathbf{T}^{(5)} & = \mathbf{S}^2\mathbf{\Omega} - \mathbf{\Omega}\mathbf{S}^2,\\ \mathbf{T}^{(6)} & = \mathbf{S}\mathbf{\Omega}^2 + \mathbf{\Omega}^2\mathbf{S} - \frac{2}{3}\textit{IV}\mathbf{I},\\ \mathbf{T}^{(7)} & = \mathbf{S}^2\mathbf{\Omega}^2 + \mathbf{\Omega}^2\mathbf{S}^2 - \frac{2}{3}\textit{V}\mathbf{I},\\ \mathbf{T}^{(8)} & = \mathbf{S}\mathbf{\Omega}\mathbf{S}^2 - \mathbf{S}^2\mathbf{\Omega}\mathbf{S}^2,\\ \mathbf{T}^{(9)} & = \mathbf{\Omega}\mathbf{S}\mathbf{\Omega}^2 - \mathbf{\Omega}^2\mathbf{S}\mathbf{\Omega},\\ \mathbf{T}^{(10)} & = \mathbf{\Omega}\mathbf{S}^2\mathbf{\Omega}^2 - \mathbf{\Omega}^2\mathbf{S}^2\mathbf{\Omega}. \end{aligned} \ \ \ \ \ (9)

The Reynolds stress tensor is then related to the anisotropy as

\displaystyle \overline{\rho{u_iu_j}} = \rho{k}\left(a_{ij}+\frac{2}{3}\delta_{ij}\right), \ \ \ \ \ (10)

which can be rewritten as an effective eddy viscosity, {\mu_t}, and an extra contribution to the expression for the stress tensor (!xx!), which reads

\displaystyle \rho k a_{ij}^{\mathrm extra}. \ \ \ \ \ (11)

3.2. The length-scale determining ({\omega}) models
This constitutive relation for the Reynolds stress tensor can, in principle, be coupled to any two-equation model platform. EARSM is usually available with the standard {k-\omega} model by Wilcox (1988), the {k-\omega} model by Kok (2000) and the {k-\omega} model by Menter (1994).

These {k-\omega} models are basically derived for a standard eddy-viscosity relation but the requirements from an EARSM are different. A new {k-\omega} model derived with the Wallin & Johansson (2000) EARSM and, thus, completely consistent with it is the model by Hellsten (2004). This model was extensively tested in different flows and is the recommended as default model.

3.3. Curvature Correction
In strong rotation and curvature the EARSM approximation is not perfectly valid and a full Reynolds stress model is to prefere. However, an approximation of the missing terms is given by Wallin & Johansson (2002). The approximation is given in terms of gradients of the strain- and rotation rate tensors. The model has some stability problems, mostly because of the dependency on second spatial derivatives, and should be used with caution. However, in cases with strong swirl the model is expected to give improvements. The curvature correction is available together with the standard {k-\omega} model by Wilcox (1988), the {k-\omega} model by Menter (1994) and the {k-\omega} model by Hellsten (2004).

4. Differential Reynolds Stress Models

The transport equation for the Reynolds stress tensor may be derived from the Navier-Stokes equations

\displaystyle \frac{{\mathrm D}\overline{u_iu_j}}{{\mathrm D} t} = {\mathcal P}_{ij}-\varepsilon_{ij}+\Pi_{ij}+{\mathcal D}_{ij}. \ \ \ \ \ (12)

The terms represent production, dissipation, pressure-strain rate, and diffusion (molecular and turbulent), respectively. One should particularly note that the production term is explicit in the Reynolds stresses,

\displaystyle {\mathcal P}_{ij} = -\overline{u_iu_k}\frac{\partial U_j}{\partial x_k} -\overline{u_ju_k}\frac{\partial U_i}{\partial x_k} \ \ \ \ \ (13)

whereas the other terms need to be modelled.

In differential Reynolds stress models (DRSM), or Reynolds stress transport (RST) models, all different terms in (12) are kept or modelled which results in a transport equation for every individual Reynolds stress component. In general three-dimensional mean-flows this implies six equations due to symmetry in the Reynolds stress tensor.

The dissipation rate tensor {\varepsilon_{ij}} is usually decomposed into an isotropic part and a deviation from that, {\varepsilon_{ij}=\varepsilon(e_{ij}+2\delta_{ij}/3)}. First, the total dissipation rate {\varepsilon} is modelled through a transport equation, similar to the {\varepsilon} equation in the {K}—[/latex]{\varepsilon}[/latex] models. Also other alternatives to {\varepsilon}, such as {\omega} or {\tau} exist. The dissipation rate anisotropy {e_{ij}} is typically explicitly modelled in terms of the Reynolds stress anisotropy or included into the modelling of the pressure strain rate.

The standard Wilcox {\omega} model, used together with the DRSM model reads

\displaystyle \frac{{\mathrm D}\omega}{{\mathrm D} t} = \alpha\frac{\omega}{K}{\mathcal P}-\beta\omega^2 + \frac{\partial}{\partial x_l}\left[ \left(\nu+\frac{\nu_T^{\mathrm(eff)}}{\sigma_\omega}\right) \frac{\partial \omega}{\partial x_l}\right], \ \ \ \ \ (14)

from which the dissipation {\varepsilon=\beta^*\omega{K}}.

The diffusion is modelled using a simple gradient diffusion model

\displaystyle {\mathcal D}_{ij} = \frac{\partial}{\partial x_l}\left[ \left(\nu+\frac{\nu_T^{\mathrm(eff)}}{\sigma_K}\right) \frac{\partial \overline{u_iu_j}}{\partial x_l}\right], \ \ \ \ \ (15)

with the effective viscosity defined as {\nu_T^{\mathrm(eff)}=K/\omega}.

The modelling of the pressure strain rate is given in the following form that includes many of the linear or quasi-linear models in literature, such as e.g. the SSG model \cite{speziale-sarkar-gatski:1991}. A general form for the pressure-strain rate and dissipation rate anisotropy {e_{ij}} lumped together reads \frac{\Pi_{ij}}{\varepsilon}-e_{ij} = -\frac{1}{2}\left(C_1^0+C_1^1\frac{\mathcal P}{\varepsilon}\right)a_{ij} +\left(C_2-\frac{C_2^*}{2}\sqrt{a_{kl}a_{lk}}\right){\tau}S_{ij} \nonumber

+{C_3\over{2}}{\tau}\left(a_{ik}S_{kj}+S_{ik}a_{kj} -{2\over{3}}a_{kl}S_{lk}\delta_{ij}\right) -{C_4\over{2}}{\tau}\left(a_{ik}\Omega_{kj}-\Omega_{ik}a_{kj}\right)

+\frac{C_5}{4}\left(a_{ik}a_{kj}-\frac{1}{3}a_{kl}a_{lk}\delta_{ij}\right) \nonumber where {\tau=K/\varepsilon} is the turbulent timescale. The {C} coefficients for the models are given in table below.

{C_1^0} {C_1^1} {C_2} {C_2^*} {C_3} {C_4} {C_5}
W&J {4.6} {1.24} {0.47} {0} {2} {0.56} {0}
SSG {3.4} {1.8} {0.8} {1.30} {1.25} {0.40} {4.2}

The two differential Reynolds stress models are available together with the standard {k-\omega} model by Wilcox (1988) and the {k-\omega} model by Hellsten (2004).

5. DES and Hybrid RANS-LES Models

DES (Detached Eddy Simulation) and other hybrid methods, combining RANS (Reynolds–Averaged Navier–Stokes) and LES (Large Eddy Simulation) in the turbulence model, are emerging turbulence modelling approaches. Note that DES is one type of “hybrid RANS-LES models”. The development of such methods has been originally motivated due to aeronautical applications, where the turbulent flow is characterized by unsteadiness, massive separation and vortical motions. In order to avoid using huge near-wall grid resolution as required in well-resolved LES for high-Reynolds number flows, these models often use RANS mode in the wall boundary layer being coupled with LES mode in the off-wall region and in regions where the flow is “detached” from wall surfaces (e.g. afterbody flows). The DES and other hybrid RANS-LES models are often viewed as a promising compromise between RANS modelling (in terms of computational efficiency) and LES (in terms of computational accuracy).

The use of LES mode in hybrid RANS-LES methods implies that instantaneous, large-scale flow structures are resolved in the LES region. When DES or other RANS-LES models are selected, the simulation should be carried out being time-dependent and usually in a three-dimensional computational domain.

5.1. The Spalart-Allmaras DES Model
The DES modelling approach was first proposed in 1997 by the pioneering work of Spalart-Allmaras (1997). This approach employs the Spalart-Allmaras (S-A) one-equation model in the wall boundary layer in combination with its SGS (SubGrid Scale) modelling variant in the LES region away from the wall. The same turbulence transport equation is thus invoked in both RANS and LES regions through a switch of turbulence length scales.

The S-A model solves a transport equation for a working eddy viscosity, {\tilde{\nu}_t}. For details, the readers should refer to Spalart-Allmaras (19947) for the S-A RANS model and Spalart-Allmaras (1997) for its DES version. Nevertheless, without repeating the details about the model constants and empirical functions, the equation for {\tilde{\nu}_t} is defined as \frac{\partial \tilde{\nu}_t}{\partial t} = C_{b1} \tilde{S} \tilde{\nu}_t + \frac{1}{\sigma} \left\{\frac{\partial}{\partial x_j} \left[ \left( \nu + \tilde{\nu}_t \right) \frac{\partial \tilde{\nu}_t}{\partial x_j} \right] + C_{b2} \frac{\partial \tilde{\nu}_t}{\partial x_j} \frac{\partial \tilde{\nu}_t}{\partial x_j}\right\} -C_{w1} f_w \left[ \frac{\tilde{\nu}_t}{l}\right]^{2} with \tilde{S} = \Omega + \frac{\tilde{\nu}_t}{\kappa^2 l^2} f_{\nu 2} where {\Omega} is the magnitude of vorticity, {l} is the local turbulent length scale, {\kappa} is the von Karman constant and {\sigma} is a model constant. The turbulence eddy viscosity is computed by {\mu_t = f_{\nu 1} \rho \tilde{\nu}_t}. The {C}‘s and {f}‘s with different subscripts are model constants and empirical model functions.

When the turbulent length scale, {l}, is taken as being the local wall distance, {d}, namely, {l=d}, the model takes the form of the original S-A RANS model. The SGS model in the off-wall LES region is adjusted from the RANS model by switching the length scale {l} from the local wall distance to a SGS turbulence length scale in association to the local cell size. This is done by the formulation of {l= \min(d, C_{des} \Delta)}, where {\Delta} is the maximum size of each cell in the three directions for structured grid and {C_{des}} is a model constant. With unstructured grid, we have estimated {\Delta} by taking the maximum edge of each cell. The model constant, {C_{des}}, has been calibrated in LES for decaying, homogeneous, isotropic turbulence, which gives {C_{des} =0.65} for best simulation of turbulence energy decay. Note that, when the {\tilde{\nu}_t}-equation functions as a transport equation for SGS turbulence, under the assumption of local equilibrium the model is similar to the well-known SGS model by Smagorinsky (1963) (but with a linear alignment between {\tilde{\nu}_t} and {\Omega}).

5.2 The Peng Hybrid RANS-LES Model
The hybrid RANS-LES model by Peng (2005) is based on an algebraic formulation for the turbulent eddy viscosity for both the RANS and LES modes. Simple algebraic (zero-equation) RANS models have been proved very robust in modelling attached boundary layers. The use of an algebraic RANS mode in hybrid RANS-LES modelling is thus an effective choice to pave the near-wall attached boundary layer. On the LES side, a presumable argumentation is that the unresolved SGS turbulence resulting from a filtering process (with a sufficiently small filter width) is more isotropic and thus its modelling may also be accomplished with a relatively simple and amenable approximation.

The Peng hybrid RANS-LES model combines an algebraic mixing-length type RANS mode in the wall layer with the SGS model by Smagorinsky (1963). The RANS eddy viscosity is formulated by \tilde{\mu}_t = \rho \tilde{l}^2_{\mu}\vert S \vert where the length scale {\tilde{l}_{\mu} = f_{\mu}\kappa d} and {f_{\mu}} is damping function. In the off-wall LES region, the SGS eddy viscosity with the Smagorinsky model reads \mu_ {sgs} = \rho (C_s \Delta)^2 \vert S \vert with {C_s = 0.12} and {\Delta} is the filter width. The matching between the RANS and LES modes is accomplished by modifying the RANS turbulent length scale over the RANS-LES interface into {l_\mu = \tilde{l}_{\mu} f_s} so that {\mu_t = \rho l_{\mu}^2 \vert S \vert} in the RANS region, where {f_s} is a matching function. The eddy viscosity, {\mu_h}, in the hybrid RANS-LES model is computed by \mu_h = \left\{ {ll} \mu_t &amp \textrm{if {\tilde{l}_{\mu} < \Delta}}

\mu_{sgs} & \textrm{if {\tilde{l}_{\mu} \geq \Delta}} \right.

Apart from its use for external aerodynamic flows, the use of the model for internal flows with moderate separation has also been demonstrated by Peng (2006). This hybrid model does not solve for any additional trubulence transport equation. It is thus more computationally efficient than the S-A DES model with a reduction of CPU time by about {20\%}.