programmer's documentation
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
Variables
cs_turbulence_model.c File Reference
#include "cs_defs.h"
#include <assert.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "bft_mem.h"
#include "bft_error.h"
#include "bft_printf.h"
#include "cs_log.h"
#include "cs_map.h"
#include "cs_parall.h"
#include "cs_mesh_location.h"
#include "cs_turbulence_model.h"
Include dependency graph for cs_turbulence_model.c:

Variables

const double xkappa = 0.42
 
const double cstlog = 5.2
 
const double apow = 8.3
 
const double bpow = 1.0/7.0
 
double dpow
 
const double cmu = 0.09
 
double cmu025
 
const double ce1 = 1.44
 
const double ce2 = 1.92
 
const double ce4 = 1.20
 
const double sigmak = 1.0
 
double sigmae
 
const double crij1 = 1.80
 
const double crij2 = 0.60
 
const double crij3 = 0.55
 
const double crijp1 = 0.50
 
const double crijp2 = 0.30
 
const double cssge2 = 1.83
 
const double cssgs1 = 1.70
 
const double cssgs2 = -1.05
 
const double cssgr1 = 0.90
 
const double cssgr2 = 0.80
 
const double cssgr3 = 0.65
 
const double cssgr4 = 0.625
 
const double cssgr5 = 0.20
 
const double cebms1 = 1.70
 
const double cebms2 = 0.
 
const double cebmr1 = 0.90
 
const double cebmr2 = 0.80
 
const double cebmr3 = 0.65
 
const double cebmr4 = 0.625
 
const double cebmr5 = 0.20
 
const double cebmr6 = 0.6
 
double csrij
 
const double cebme2 = 1.83
 
const double cebmmu = 0.22
 
const double xcl = 0.122
 
const double xa1 = 0.1
 
const double xct = 6.0
 
const double xceta = 80.0
 
const double cpale1 = 1.44
 
const double cpale2 = 1.83
 
const double cpale3 = 2.3
 
const double cpale4 = 0.4
 
const double cpalse = 1.5
 
const double cpalmu = 0.22
 
const double cpalc1 = 1.7
 
const double cpalc2 = 0.9
 
const double cpalct = 4.0
 
const double cpalcl = 0.164
 
const double cpalet = 75.0
 
const double ckwsk1 = 1.0/0.85
 
const double ckwsk2 = 1.0
 
const double ckwsw1 = 2.0
 
const double ckwsw2 = 1.0/0.856
 
const double ckwbt1 = 0.075
 
const double ckwbt2 = 0.0828
 
double ckwgm1
 
double ckwgm2
 
const double ckwa1 = 0.31
 
const double ckwc1 = 10.0
 
const double csab1 = 0.1355
 
const double csab2 = 0.622
 
const double csasig = 2.0/3.0
 
const double csav1 = 7.1
 
double csaw1
 
const double csaw2 = 0.3
 
const double csaw3 = 2.0
 
const double cssr1 = 1.0
 
const double cssr2 = 12.0
 
const double cssr3 = 1.0
 
const double ccaze2 = 1.83
 
const double ccazsc = 0.119
 
const double ccaza = 4.3
 
const double ccazb = 5.130
 
const double ccazc = 0.453
 
const double ccazd = 0.682
 
const double xlesfl = 2.0
 
const double ales = 1.0
 
const double bles = 1.0/3.0
 
const double csmago = 0.065
 
const double xlesfd = 1.5
 
double smagmx
 
const double cdries = 26.0
 
const double cv2fa1 = 0.05
 
const double cv2fe2 = 1.85
 
const double cv2fmu = 0.22
 
const double cv2fc1 = 1.4
 
const double cv2fc2 = 0.3
 
const double cv2fct = 6.0
 
const double cv2fcl = 0.25
 
const double cv2fet = 110.0
 
const double cwale = 0.25
 
const double xiafm = 0.7
 
const double etaafm = 0.4
 
const double c1trit = 4.15
 
const double c2trit = 0.55
 
const double c3trit = 0.5
 
const double c4trit = 0.
 
const double cthafm = 0.236
 
const double cthdfm = 0.31
 

Detailed Description

Base turbulence model data.

Variable Documentation

const double ales = 1.0

Constant used to define, for each cell $\Omega_i$, the width of the (implicit) filter:

  • $\overline{\Delta}=xlesfl(ales*|\Omega_i|)^{bles}$

Useful if and only if iturb = 40 or 41.

const double apow = 8.3

Werner and Wengle coefficient

const double bles = 1.0/3.0

Constant used to define, for each cell $$, the width of the (implicit) filter:

  • $\overline{\Delta}=xlesfl(ales*|\Omega_i|)^{bles}$

Useful if and only if iturb = 40 or 41.

const double bpow = 1.0/7.0

Werner and Wengle coefficient

const double c1trit = 4.15

Coefficient of turbulent DFM flow model.

const double c2trit = 0.55

Coefficient of turbulent DFM flow model.

const double c3trit = 0.5

Coefficient of turbulent DFM flow model.

const double c4trit = 0.

Coefficient of turbulent DFM flow model.

const double ccaza = 4.3

Constants of the Cazalbou rotation/curvature correction.

const double ccazb = 5.130

Constants of the Cazalbou rotation/curvature correction.

const double ccazc = 0.453

Constants of the Cazalbou rotation/curvature correction.

const double ccazd = 0.682

Constants of the Cazalbou rotation/curvature correction.

const double ccaze2 = 1.83

Constants of the Cazalbou rotation/curvature correction.

const double ccazsc = 0.119

Constants of the Cazalbou rotation/curvature correction.

const double cdries = 26.0

Van Driest constant appearing in the van Driest damping function applied to the Smagorinsky constant:

  • $ (1-\exp^{(-y^+/cdries}) $.

Useful if and only if iturb = 40 or 41.

const double ce1 = 1.44

Constant $C_{\varepsilon 1}$ for all the RANS turbulence models except for the v2f and the $k-\omega$ models. Useful if and only if iturb= 20, 21, 30 or 31 ( $k-\varepsilon$ or $R_{ij}-\varepsilon$).

const double ce2 = 1.92

Constant $C_{\varepsilon 2}$ for the $k-\varepsilon$ and $R_{ij}-\varepsilon$ LRR models. Useful if and only if { iturb}= 20, 21 or 30 ( $k-\varepsilon$ or $R_{ij}-\varepsilon$ LRR).

const double ce4 = 1.20

Coefficient of interfacial coefficient in k-eps, used in Lagrange treatment.

Constant $C_{\varepsilon 4}$ for the interfacial term (Lagrangian module) in case of two-way coupling. Useful in case of Lagrangian modelling, in $k-\varepsilon$ and $R_{ij}-\varepsilon$ with two-way coupling.

const double cebme2 = 1.83

Constant of the Rij-epsilon EBRSM.

const double cebmmu = 0.22

Constant of the Rij-epsilon EBRSM.

const double cebmr1 = 0.90
const double cebmr2 = 0.80
const double cebmr3 = 0.65
const double cebmr4 = 0.625
const double cebmr5 = 0.20
const double cebmr6 = 0.6
const double cebms1 = 1.70

Constant of the Rij-epsilon EBRSM.

const double cebms2 = 0.

Constant of the Rij-epsilon EBRSM.

const double ckwa1 = 0.31

Specific constant of k-omega SST. Constant $a_1$ for the $k-\omega$ SST model. Useful if and only if iturb=60 ( $k-\omega$ SST).

const double ckwbt1 = 0.075

Constant $\beta_1$ for the $k-\omega$ SST model. Useful if and only if iturb=60 ( $k-\omega$ SST).

const double ckwbt2 = 0.0828

Constant $\beta_2$ for the $k-\omega$ SST model. Useful if and only if iturb=60 ( $k-\omega$ SST).

const double ckwc1 = 10.0

Constant $ c_1 $ for the $k-\omega$ SST model. Useful if and only if iturb=60 ( $k-\omega$ SST). Specific constant of k-omega SST.

double ckwgm1

$\frac{\beta_1}{C_\mu}-\frac{\kappa^2}{\sqrt{C_\mu}\sigma_{\omega 1}}$. Constant $\gamma_1$ for the $k-\omega$ SST model. Useful if and only if iturb=60 ( $k-\omega$ SST).

Warning
: $\gamma_1$ is calculated before the call to usipsu. Hence, if $\beta_1$, $C_\mu$, $\kappa$ or $\sigma_{\omega 1}$ is modified in usipsu, ckwgm1 must also be modified in accordance.
double ckwgm2

$\frac{\beta_2}{C_\mu}-\frac{\kappa^2}{\sqrt{C_\mu}\sigma_{\omega 2}}$. Constant $\gamma_2$ for the $k-\omega$ SST model. Useful if and only if iturb=60 ( $k-\omega$ SST).

Warning
: $\gamma_2$ is calculated before the call to usipsu. Hence, if $\beta_2$, $C_\mu$, $\kappa$ or $\sigma_{\omega 2}$ is modified in usipsu, ckwgm2 must also be modified in accordance.
const double ckwsk1 = 1.0/0.85

Constant $\sigma_{k1}$ for the $k-\omega$ SST model. Useful if and only if iturb=60.

const double ckwsk2 = 1.0

Constant $\sigma_{k2}$ for the $k-\omega$ SST model. Useful if and only if iturb=60.

const double ckwsw1 = 2.0

Constant $\sigma_{\omega 1}$ for the $k-\omega$ SST model. Useful if and only if iturb=60 ( $k-\omega$ SST).

const double ckwsw2 = 1.0/0.856

Constant $\sigma_{\omega 2}$ for the $k-\omega$ SST model. Useful if and only if iturb=60 ( $k-\omega$ SST).

const double cmu = 0.09

Constant $C_\mu$ for all the RANS turbulence models except for the v2f model (see cv2fmu for the value of $C_\mu$ in case of v2f modelling). Useful if and only if iturb = 20, 21, 30, 31 or 60 ( $k-\varepsilon$, $R_{ij}-\varepsilon$ or $k-\omega$).

double cmu025

$ C_\mu^\frac{1}{4} $

const double cpalc1 = 1.7

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpalc2 = 0.9

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpalcl = 0.164

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpalct = 4.0

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpale1 = 1.44

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpale2 = 1.83

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpale3 = 2.3

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpale4 = 0.4

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpalet = 75.0

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpalmu = 0.22

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double cpalse = 1.5

Specific constant of v2f "BL-v2k" (or phi-alpha).

const double crij1 = 1.80

Constant $C_1$ for the $R_{ij}-\varepsilon$ LRR model. Useful if and only if iturb=30 ( $R_{ij}-\varepsilon$ LRR).

const double crij2 = 0.60
const double crij3 = 0.55

Constant $C_3$ for the $R_{ij}-\varepsilon$ LRR model. Useful if and only if iturb=30 ( $R_{ij}-\varepsilon$ LRR).

const double crijp1 = 0.50

Constant $C_1^\prime$ for the $R_{ij}-\varepsilon$ LRR model, corresponding to the wall echo terms. Useful if and only if iturb=30 and irijec=1 ( $R_{ij}-\varepsilon$ LRR).

const double crijp2 = 0.30

Constant $C_2^\prime$ for the $R_{ij}-\varepsilon$ LRR model, corresponding to the wall echo terms. Useful if and only if iturb=30 and irijec=1 ( $R_{ij}-\varepsilon$ LRR).

const double csab1 = 0.1355

Specific constant of Spalart-Allmaras.

const double csab2 = 0.622

Specific constant of Spalart-Allmaras.

const double csasig = 2.0/3.0

Specific constant of Spalart-Allmaras.

const double csav1 = 7.1

Specific constant of Spalart-Allmaras.

double csaw1

Specific constant of Spalart-Allmaras.

const double csaw2 = 0.3

Specific constant of Spalart-Allmaras.

const double csaw3 = 2.0

Specific constant of Spalart-Allmaras.

const double csmago = 0.065

Smagorinsky constant used in the Smagorinsky model for LES. The sub-grid scale viscosity is calculated by $\displaystyle\mu_{sg}= \rho C_{smago}^2\bar{\Delta}^2\sqrt{2\bar{S}_{ij}\bar{S}_{ij}}$ where $\bar{\Delta}$ is the width of the filter and $\bar{S}_{ij}$ the filtered strain rate.

Useful if and only if iturb = 40.

Note
In theory Smagorinsky constant is 0.18. For a channel, 0.065 value is rather taken.
double csrij

Constant $C_s$ for the $R_{ij}-\varepsilon$ LRR model. Useful if and only if iturb=30 ( $R_{ij}-\varepsilon$ LRR).

const double cssge2 = 1.83

Constant $C_{\varepsilon 2}$ for the $R_{ij}-\varepsilon$ SSG model. Useful if and only if iturb=31 ( $R_{ij}-\varepsilon$ SSG).

const double cssgr1 = 0.90

Constant $C_{r1}$ for the $R_{ij}-\varepsilon$ SSG model. Useful if and only if iturb=31 ( $R_{ij}-\varepsilon$ SSG).

const double cssgr2 = 0.80

Constant $C_{r2}$ for the $R_{ij}-\varepsilon$ SSG model. Useful if and only if iturb=31 ( $R_{ij}-\varepsilon$ SSG).

const double cssgr3 = 0.65

Constant $C_{r3}$ for the $R_{ij}-\varepsilon$ SSG model. Useful if and only if iturb=31 ( $R_{ij}-\varepsilon$ SSG).

const double cssgr4 = 0.625

constant $C_{r4}$ for the $R_{ij}-\varepsilon$ SSG model. Useful if and only if iturb=31 ( $R_{ij}-\varepsilon$ SSG).

const double cssgr5 = 0.20

Constant $C_{r1}$ for the $R_{ij}-\varepsilon$ SSG model. Useful if and only if iturb=31 ( $R_{ij}-\varepsilon$ SSG).

const double cssgs1 = 1.70

Constant $C_{s1}$ for the $R_{ij}-\varepsilon$ SSG model. Useful if and only if iturb=31 ( $R_{ij}-\varepsilon$ SSG).

const double cssgs2 = -1.05

Constant $C_{s2}$ for the $R_{ij}-\varepsilon$ SSG model. Useful if and only if iturb=31 ( $R_{ij}-\varepsilon$ SSG).

const double cssr1 = 1.0

Constant of the Spalart-Shur rotation/curvature correction.

const double cssr2 = 12.0

Constant of the Spalart-Shur rotation/curvature correction.

const double cssr3 = 1.0

Constant of the Spalart-Shur rotation/curvature correction.

const double cstlog = 5.2

Constant of logarithmic law function: $ \dfrac{1}{\kappa} \ln(y^+) + cstlog $ ( $ cstlog = 5.2 $).

Constant of the logarithmic wall function. Useful if and only if iturb >= 10 (mixing length, $k-\varepsilon$, $R_{ij}-\varepsilon$, LES, v2f or $k-\omega$).

const double cthafm = 0.236

Constant of GGDH and AFM on the thermal scalar.

const double cthdfm = 0.31

Constant of GGDH and AFM on the thermal scalar.

const double cv2fa1 = 0.05

Constant $a_1$ for the v2f $\varphi$-model. Useful if and only if iturb=50 (v2f $\varphi$-model).

const double cv2fc1 = 1.4

Constant $C_1$ for the v2f $\varphi$-model. Useful if and only if iturb=50 (v2f $\varphi$-model).

const double cv2fc2 = 0.3

Constant $C_2$ for the v2f $\varphi$-model. Useful if and only if iturb=50 (v2f $\varphi$-model).

const double cv2fcl = 0.25

Constant $C_L$ for the v2f $\varphi$-model. Useful if and only if iturb=50 (v2f $\varphi$-model).

const double cv2fct = 6.0

Constant $C_T$ for the v2f $\varphi$-model. Useful if and only if iturb=50 (v2f $\varphi$-model).

const double cv2fe2 = 1.85

Constant $C_{\varepsilon 2}$ for the v2f $\varphi$-model. Useful if and only if iturb=50 (v2f $\varphi$-model).

const double cv2fet = 110.0

Constant $C_\eta$ for the v2f $\varphi$-model. Useful if and only if iturb=50 (v2f $\varphi$-model).

const double cv2fmu = 0.22

Constant $C_\mu$ for the v2f $\varphi$-model. Useful if and only if iturb=50 (v2f $\varphi$-model).

const double cwale = 0.25

Constant of the WALE LES method.

double dpow

Werner and Wengle coefficient

const double etaafm = 0.4

Coefficient of turbulent AFM flow model.

double sigmae

Prandtl number for $\varepsilon$. Useful if and only if iturb= 20, 21, 30, 31 or 50 ( $k-\varepsilon$, $R_{ij}-\varepsilon$ or v2f).

const double sigmak = 1.0

Prandtl number for $k$ with $k-\varepsilon$ and v2f models. Useful if and only if iturb=20, 21 or 50 ( $k-\varepsilon$ or v2f).

double smagmx

Maximum allowed value for the variable $C$ appearing in the LES dynamic model (the "square" comes from the fact that the variable of the dynamic model corresponds to the square of the constant of the Smagorinsky model). Any larger value yielded by the calculation procedure of the dynamic model will be clipped to $ smagmx^2$.

Useful if and only if iturb = 41.

const double xa1 = 0.1

Constant in the expression of Ce1' for the Rij-epsilon EBRSM.

const double xceta = 80.0

Constant of the Rij-epsilon EBRSM.

const double xcl = 0.122

Constant of the Rij-epsilon EBRSM.

const double xct = 6.0

Constant of the Rij-epsilon EBRSM.

const double xiafm = 0.7

Coefficient of turbulent AFM flow model.

const double xkappa = 0.42

(end ignore by Doxygen)

Karman constant. (= 0.42)

Useful if and only if iturb >= 10. (mixing length, $k-\varepsilon$, $R_{ij}-\varepsilon$, LES, v2f or $k-\omega$).

const double xlesfd = 1.5

Ratio between explicit and explicit filter width for a dynamic model. Constant used to define, for each cell $\Omega_i$, the width of the explicit filter used in the framework of the LES dynamic model: $\widetilde{\overline{\Delta}}=xlesfd\overline{\Delta}$.

Useful if and only if iturb = 41.

const double xlesfl = 2.0

Constant used in the definition of LES filtering diameter: $ \delta = \text{xlesfl} . (\text{ales} . volume)^{\text{bles}} $.