This section details the turbulence models implemented in Caelus 4.10. Caelus 4.10 supports both steady-state and transient turbulence models. From an averaging point of view, Caelus 4.10 supports Reynolds averaged (RANS) and large eddy simulations (LES). The current focus of turbulence modelling in Caelus 4.10 is on incompressible flows. Future versions of Caelus will support compressible flows and appropriate compressible turbulence models.

**Model Description**

The Spalart-Allmaras model uses one equation to compute the turbulent viscosity. The Spalart-Allmaras model used in Caelus 4.10 corresponds to the version of Shur et al. [9] with curvature correction. The transport equation in this version is the same as the original Spalart-Allmaras model [10] but with rotation and curvature effects accounted for. This modification is accomplished by introducing the \(f_{r1}\) function into the transport equation for the viscosity-like variable

\[\partial_t \tilde{\nu} + u_j \partial_j \tilde{\nu} = c_{b1}f_{r1}\tilde{S}\tilde{\nu} - c_{w1}f_{w1}\left(\frac{\tilde{\nu}}{d}\right)^2 + \frac{1}{\sigma}\left[\partial_j(\nu+\tilde{\nu})\partial_j \tilde{\nu} + c_{b2}\partial_i \tilde{\nu} \partial_i \tilde{\nu} \right]\]

Terms appearing in the Spalart-Allmaras turbulence model have the following definitions:

**Model definitions**

Symbol | Definition | Formula |
---|---|---|

\(\tilde{\nu}\) | Viscosity-like variable | |

\(\nu_t\) | Turbulent eddy viscosity | \(\tilde{\nu} f_{v1}\) |

\(f_{v1}\) | Viscous Damping Function | \(\frac{\chi^3}{\chi^3 + c_{v1}}\) |

\(\chi\) | Viscosity ratio | \(\frac{\tilde{\nu}}{\nu}\) |

\(\nu\) | Kinematic Viscosity | |

\(\tilde{S}\) | Modified vorticity | \(\Omega + \frac{\tilde{\nu}}{\kappa^2 d^2}f_{v2}\) |

\(\Omega\) | Vorticity magnitude | \(\Omega = \sqrt{2W_{ij}W_{ij}}\) |

\(W_{ij}\) | Vorticity | \(\frac{1}{2}(\partial_j u_i - \partial_i u_j)\) |

\(d\) | Wall distance | |

\(f_{v2}\) | Function | \(1-\frac{\chi}{1+\chi f_{v1}}\) |

\(f_{w}\) | Function | \(g\left[\frac{1+c_{w3}^6}{g^5+c_{w3}^6}\right]^{1/6}\) |

\(g\) | Function | \(r+c_{w2}(r^6-r)\) |

\(r\) | Variable | \(r=min\left[\frac{\tilde{\nu}}{\tilde{S} \kappa^2 d^2},10\right]\) |

It should be noted that the Caelus 4.10 implementation of the Spalart-Allmaras turbulence model omits the rotational and turbulence tripping terms as the flow is assumed fully turbulent.

**Model constants**

\(c_{b1} = 0.1355\)

\(\sigma = \frac{2}{3}\)

\(c_{b2} = 0.622\)

\(c_{w2} = 0.3\)

\(c_{w3} = 2.0\)

\(c_{v1} = 7.1\)

\(c_{w1} = \frac{c_{b1}}{\kappa^2} + \frac{1+c_{b2}}{\sigma}\)

\(\kappa = 0.41\)

The Spalart-Allmaras model neglects the contribution of the turbulent kinetic energy to the eddy viscosity. This is deemed acceptable for many incompressible and low Mach number compressible flows.

**Boundary conditions**

The Spalart-Allmaras turbulence model was originally developed for external aerodynamics of streamlined bodies. Typically, mesh resolutions used for external aerodynamics have \(y^+\) values close to 1. Therefore, the Spalart-Allmaras model is often referred to as a low Reynolds number model when it comes to near wall treatment. In the case of meshes resolved near the wall, i.e. \(y^+ \sim 1\), the following set of boundary conditions are recommended:

\(\tilde{\nu}|_{wall} = 0\)

\(\nu_t|_{wall} = 0\)

However, in Caelus 4.10 it is possible to use the Spalart-Allmaras on meshes that are not wall-resolved, i.e., meshes that have the first cell height next to the wall \(y^+ \gg 1\). In such cases, a wall function for \(\nu_t\) based on the velocity is recommended.

\(\nu_t|_{wall} =\) nutUWallFunction (\(y^+ > 30\))

Inlet values of \(\tilde{\nu}\) at the far field boundary for external flows is suggested as follows:

\(\tilde{\nu}|_{inlet} = 3 \nu\) to \(5 \nu\).

\(\nu_t\) should be set to “calculated”.

It is also possible to specify \(\tilde{\nu}\) from the turbulence intensity and length scale as follows:

\(\tilde{\nu}|_{inlet} = \sqrt{\frac{3}{2}} |u_i| I l\)

where \(I = \frac{|u_i'|}{|u_i|}\) is turbulence intensity and

\(l\) is a turbulence length scale (measure of largest eddies in the flow field).

In the case of internal flows, depending on the level of turbulence intensity, \(\tilde{\nu}\) can be determined from the turbulence viscosity ratio

\(\tilde{\nu} = 10\nu\) or higher and

\(\nu_t\) should be set to “calculated”.

At the outlet, zero gradient an acceptable choice.

**Model description**

The realizable \(k\)–\(\varepsilon\) turbulence model is an improvement on the standard \(k\)–\(\varepsilon\) model as it is modified to guarantee realizability of normal Reynolds stresses [13]. In the standard \(k\)–\(\varepsilon\) turbulence model it is possible, under certain conditions, to have negative normal Reynolds stresses, i.e. non-realizable. The realizable \(k\)–\(\varepsilon\) turbulence model resolves this problem by introducing a non-constant definition of \(C_{\mu}\) and by modifying the \(\varepsilon\) equation. The transport equations for the realizable \(k\)–\(\varepsilon\) turbulence model are as follows:

\[\partial_t k + \partial_j (u_j k) = \partial_j\left[\mu + \frac{\mu_t}{\sigma_k}\partial_j k\right] + G_k -\varepsilon\]

\[\partial_t \varepsilon + \partial_j (u_j \varepsilon) = \partial_j\left[\mu + \frac{\mu_t}{\sigma_{\varepsilon}}\partial_j \varepsilon \right] + C_1 S \varepsilon - C_2 \frac{\varepsilon^2}{k+ \sqrt{\nu\varepsilon}}\]

The terms appearing in the realizable \(k\)–\(\varepsilon\) turbulence model have the following definitions:

**Model definitions**

Symbol | Definition | Formula |
---|---|---|

\(\nu_t\) | Turbulent viscosity | \(C_{\mu}\frac{k^2}{\varepsilon}\) |

\(C_1\) | Model parameter | \(max \left(0.43,\frac{\eta}{\eta+5}\right)\) |

\(\eta\) | Strain magnitude over time scale | \(S\frac{k}{\varepsilon}\) |

\(S\) | Strain magnitude | \(S = \sqrt{2S_{ij}S_{ij}}\) |

\(S_{ij}\) | Strain tensor | \(\frac{1}{2}(\partial_j u_i + \partial_i u_j)\) |

\(C_{\mu}\) | Model parameter | \(\frac{1}{A_0 + A_S \frac{k U^*}{\varepsilon}}\) |

\(U^*\) | Velocity | \(\sqrt{S_{ij}S_{ij}+\tilde{\Omega_{ij}}\tilde{\Omega_{ij}}}\) |

\(\overline{\Omega}_{ij}\) | Rate of rotation tensor | \(\frac{1}{2}(\partial_j u_i - \partial_i u_j)\) |

\(\omega_k\) | Angular velocity | |

\(\Omega_{ij}\) | Rotation tensor | \(\overline{\Omega}_{ij} - \varepsilon_{ijk} \omega_k\) |

\(\tilde{\Omega}_{ij}\) | Modified rotation tensor | \(\Omega_{ij} - 2 \varepsilon_{ijk} \omega_k\) |

\(A_0\) | Model constant | \(4.04\) |

\(A_S\) | Function | \(\sqrt{6} cos(\phi)\) |

\(\phi\) | Function | \(\frac{1}{3} cos^{-1} (\sqrt{6} W)\) |

\(W\) | Variable | \(\frac{S_ijS_{jk}S{ki}}{S^3}\) |

**Model constants**

\(C_{1\varepsilon} = 1.44\)

\(C_2 = 1.9\)

\(\sigma_k = 1.0\)

\(\sigma_{\varepsilon} = 1.2\)

**Boundary conditions**

The realizable \(k\)–\(\varepsilon\) was developed as a high-Reynolds number model. In other words, the wall treatment is such that wall functions for \(\nu_t\), \(k\), and \(\varepsilon\) fields must be always used. It is recommended that velocity based wall functions for the \(\nu_t\) field is used. On the inlet, the following relations can be used to determine the values of \(k\) and \(\varepsilon\) fields:

\(k|_{inlet} = \frac{3}{2} (U_{inlet} I)^2\),

\(\varepsilon = C_{\mu}^{3/4} \frac{k^{3/2}}{l}\),

\(\nu_t\) should be set to “calculated” or computed from the viscosity ratio.

\(I = \frac{|u_i'|}{|u_i|}\) is the turbulence intensity

and \(l\) is a turbulence length scale.

On outlets, zero gradient is an acceptable choice.

**Model description**

The \(k\)–\(\omega\) shear stress transport (SST) turbulence model was introduce by Menter [16] as an improvement on the original two equation \(k\)–\(\omega\) model of Wilcox [14]. The SST formulation uses the \(k\)–\(\omega\) formulation in the inner parts of the boundary layer and it switches to a \(k\)–\(\varepsilon\) behaviour in the free-stream portions of the flow domain. Thus, the \(k\)–\(\omega\) SST turbulence model does not require a damping function close to the wall enabling its use with near-wall resolved meshes without modification. At the same time, the \(k\)–\(\omega\) SST models reverts to a \(k\)–\(\varepsilon\) model behaviour away from the walls thus removing the sensitivity of the original \(k\)–\(\omega\) model to the \(\omega\) free-stream boundary value.

The implementation of the \(k\)–\(\omega\) SST model in Caelus 4.10 is that of Kuntz et al [3]. The transport equations for turbulence kinetic energy and specific dissipation are:

\[\partial_t k + \partial_j (u_j k) = \frac{P}{\rho} - \beta^* \omega k + \partial_j[(\nu + \sigma_k \nu_t) \partial_j k]\]

\[\partial_t \omega + \partial_j(u_j \omega) = \frac{\gamma P}{\rho \nu_t} - \beta \omega^2 + \partial_j[(\nu + \nu_t) \partial_j \omega] + 2(1-F_1) \frac{\sigma_{\omega2}}{\omega} \partial_j k \partial_j \omega\]

It should be noted that this set of transport equations only applicable for incompressible flows as density was assumed constant. Terms appearing in the \(k\)–\(\omega\) SST model have the following definitions:

**Model definitions**

Symbol | Definition | Formula |
---|---|---|

\(\nu_t\) | Turbulent viscosity | \(\frac{a_1 k}{max (a_1 \omega, S F_2)}\) |

\(S_{ij}\) | Symmetric part of stress tensor | \(\frac{1}{2}(\partial_j u_i + \partial_i u_j)\) |

\(S\) | Strain invariant | \(S=\sqrt{S_{ij}S_{ij}}\) |

\(F_1\) | Blending function | \(tanh(arg_1^4)\) |

\(arg_1\) | Argument for \(F_1\) function | \(min\left[max \left(\frac{\sqrt{k}}{\beta^* \omega d},\frac{500 \nu}{d^2 \omega}\right), \frac{4 \sigma_{\omega 2} k}{CD_{k \omega} d^2}\right]\) |

\(CD_{k \omega}\) | Function | \(max\left(2 \sigma_{\omega 2}\frac{1}{\omega} \frac{\partial k}{\partial x_j} \frac{\partial \omega}{\partial x_j}, 10^{-10}\right)\) |

\(F_2\) | Blending function | \(tanh(arg_2^2)\) |

\(arg_2\) | Argument for \(F_2\) function | \(max \left(2 \frac{\sqrt{k}}{\beta^* \omega d},\frac{500 \nu}{d^2 \omega}\right)\) |

This variant of the \(k\)–\(\omega\) SST turbulence model uses a production limiter in both \(k\) and \(\omega\) equations,

\[P = min(P, 10 \beta^* \omega k).\]

**Model constants**

\(\gamma_1 = \frac{5}{9}\)

\(\gamma2 = 0.44\)

\(\sigma_{k 1} = 0.85\)

\(\sigma_{\omega 1} = 0.5\)

\(\beta_1 = 0.075\)

\(\sigma_{k 2} = 1.0\)

\(\sigma_{\omega 2} = 0.856\)

\(\beta_2 = 0.0828\)

\(\beta^* = 0.09\)

\(\kappa = 0.41\)

\(a_1 = 0.31\)

**Boundary conditions**

As described above, the \(k\)–\(\omega\) SST turbulence model does not require damping functions near the wall. Thus it can be used for both resolved and unresolved near-wall meshes. In the case of meshes resolved near the wall, i.e. \(y^+ \sim 1\), the following set of boundary conditions are recommended:

\(k_{wall} = 0\)

\(\omega_{wall} =\) omegaWallFunction

\(\nu_t|_{wall} = 0\)

For meshes that are not wall-resolved, i.e., meshes that have the first cell height next to the wall \(y^+ \gg 1\), wall functions will be needed. In such cases, a wall function for \(\nu_t\) based on the velocity is recommended.

\(k_{wall} =\) kqRWallFunction

\(\omega_{wall} =\) omegaWallFunction

\(\nu_t|_{wall} =\) nutUSpaldingWallFunction or alternatively nutUWallFunction if \(y^+ \gt 30\)

Inlet values of can be determined as follows:

\(k_{inlet} = \frac{3}{2} (U_{inlet} I)^2\),

\(\omega = C_{\mu}^{-1/4} \frac{\sqrt{k}}{l}\),

\(\nu_t\) should be set to “calculated” or computed from the turbulent viscosity relation.

\(I = \frac{|u_i'|}{|u_i|}\) is the turbulence intensity

and \(l\) is a turbulence length scale.

On outlets, zero gradient is an acceptable choice.