#include "cs_defs.h"
Go to the source code of this file.
Data Structures | |
struct | cs_turb_model_t |
Turbulence model general options descriptor. More... | |
struct | cs_turb_ref_values_t |
struct | cs_turb_rans_model_t |
RANS turbulence model descriptor. More... | |
struct | cs_turb_les_model_t |
LES turbulence model descriptor. More... | |
struct | cs_turb_hybrid_model_t |
Hybrid turbulence model descriptor. More... | |
Enumerations | |
enum | cs_turb_model_type_t { CS_TURB_NONE = 0 , CS_TURB_MIXING_LENGTH = 10 , CS_TURB_K_EPSILON = 20 , CS_TURB_K_EPSILON_LIN_PROD = 21 , CS_TURB_K_EPSILON_LS = 22 , CS_TURB_K_EPSILON_QUAD = 23 , CS_TURB_RIJ_EPSILON_LRR = 30 , CS_TURB_RIJ_EPSILON_SSG = 31 , CS_TURB_RIJ_EPSILON_EBRSM = 32 , CS_TURB_LES_SMAGO_CONST = 40 , CS_TURB_LES_SMAGO_DYN = 41 , CS_TURB_LES_WALE = 42 , CS_TURB_V2F_PHI = 50 , CS_TURB_V2F_BL_V2K = 51 , CS_TURB_K_OMEGA = 60 , CS_TURB_SPALART_ALLMARAS = 70 } |
enum | { CS_TURB_TYPE_NONE = 0 , CS_TURB_RANS = 1 , CS_TURB_LES = 2 , CS_TURB_HYBRID = 3 } |
enum | { CS_TURB_ALGEBRAIC = 0 , CS_TURB_FIRST_ORDER = 1 , CS_TURB_SECOND_ORDER = 2 } |
enum | { CS_HYBRID_NONE = 0 , CS_HYBRID_DES = 1 , CS_HYBRID_DDES = 2 , CS_HYBRID_SAS = 3 , CS_HYBRID_HTLES = 4 } |
Functions | |
void | cs_turb_model_init (void) |
Initialize turbulence model structures. More... | |
void | cs_set_type_order_turbulence_model (void) |
Initialize type and order members of turbulence model structure. More... | |
void | cs_set_glob_turb_model (void) |
Set global pointer to turbulence model structure. More... | |
cs_turb_model_t * | cs_get_glob_turb_model (void) |
Provide write access to turbulence model structure. More... | |
void | cs_turb_compute_constants (void) |
Compute turbulence model constants, some of which may depend on the model choice. More... | |
cs_turb_ref_values_t * | cs_get_glob_turb_ref_values (void) |
Provide access to cs_glob_turb_ref_values. More... | |
cs_turb_rans_model_t * | cs_get_glob_turb_rans_model (void) |
Provide access to cs_glob_turb_rans_model. More... | |
cs_turb_les_model_t * | cs_get_glob_turb_les_model (void) |
Provide access to cs_glob_turb_les_model. More... | |
cs_turb_hybrid_model_t * | cs_get_glob_turb_hybrid_model (void) |
Provide access to cs_glob_turb_hybrid_model. More... | |
void | cs_turb_model_log_setup (void) |
Print the turbulence model parameters to setup.log. More... | |
void | cs_turb_constants_log_setup (void) |
Print the turbulent constants to setup.log. More... | |
void | cs_turb_init_ref_quantities (void) |
Compute characteristic length for turbulence if not already done. More... | |
void | cs_clip_turbulent_fluxes (int flux_id, int ivartt) |
Clipping for the turbulence flux vector. More... | |
void | cs_turbulence_function_k (int location_id, cs_lnum_t n_elts, const cs_lnum_t *elt_ids, void *input, void *vals) |
Return or estimate the value of the turbulent kinetic energy over specified elements. More... | |
void | cs_turbulence_function_eps (int location_id, cs_lnum_t n_elts, const cs_lnum_t *elt_ids, void *input, void *vals) |
Return or estimate the value of the turbulent dissipation over specified elements. More... | |
void | cs_turbulence_function_rij (int location_id, cs_lnum_t n_elts, const cs_lnum_t *elt_ids, void *input, void *vals) |
Return or estimate the value of the Reynolds stresses over specified elements. More... | |
enum cs_turb_model_type_t |
void cs_clip_turbulent_fluxes | ( | int | flux_id, |
int | variance_id | ||
) |
Clipping for the turbulence flux vector.
[in] | flux_id | turbulent flux index |
[in] | variance_id | scalar variance index |
cs_turb_hybrid_model_t* cs_get_glob_turb_hybrid_model | ( | void | ) |
Provide access to cs_glob_turb_hybrid_model.
needed to initialize structure with GUI
cs_turb_les_model_t* cs_get_glob_turb_les_model | ( | void | ) |
Provide access to cs_glob_turb_les_model.
needed to initialize structure with GUI
cs_turb_model_t* cs_get_glob_turb_model | ( | void | ) |
Provide write access to turbulence model structure.
cs_turb_rans_model_t* cs_get_glob_turb_rans_model | ( | void | ) |
Provide access to cs_glob_turb_rans_model.
needed to initialize structure with GUI
cs_turb_ref_values_t* cs_get_glob_turb_ref_values | ( | void | ) |
Provide access to cs_glob_turb_ref_values.
needed to initialize structure with GUI
void cs_set_glob_turb_model | ( | void | ) |
Set global pointer to turbulence model structure.
This global pointer provides a read-only access to the structure.
void cs_set_type_order_turbulence_model | ( | void | ) |
Initialize type and order members of turbulence model structure.
void cs_turb_compute_constants | ( | void | ) |
Compute turbulence model constants, some of which may depend on the model choice.
void cs_turb_constants_log_setup | ( | void | ) |
Print the turbulent constants to setup.log.
void cs_turb_init_ref_quantities | ( | void | ) |
Compute characteristic length for turbulence if not already done.
void cs_turb_model_init | ( | void | ) |
Initialize turbulence model structures.
void cs_turb_model_log_setup | ( | void | ) |
Print the turbulence model parameters to setup.log.
void cs_turbulence_function_eps | ( | int | location_id, |
cs_lnum_t | n_elts, | ||
const cs_lnum_t * | elt_ids, | ||
void * | input, | ||
void * | vals | ||
) |
Return or estimate the value of the turbulent dissipation over specified elements.
Returned values are zero for turbulence models other than RANS.
This function matches the cs_eval_at_location_t function profile.
[in] | location_id | base associated mesh location id |
[in] | n_elts | number of associated elements |
[in] | elt_ids | ids of associated elements, or NULL if no filtering is required |
[in,out] | input | ignored |
[in,out] | vals | pointer to output values (size: n_elts*dimension) |
void cs_turbulence_function_k | ( | int | location_id, |
cs_lnum_t | n_elts, | ||
const cs_lnum_t * | elt_ids, | ||
void * | input, | ||
void * | vals | ||
) |
Return or estimate the value of the turbulent kinetic energy over specified elements.
Returned values are zero for turbulence models other than RANS.
This function matches the cs_eval_at_location_t function profile.
[in] | location_id | base associated mesh location id |
[in] | n_elts | number of associated elements |
[in] | elt_ids | ids of associated elements, or NULL if no filtering is required |
[in,out] | input | ignored |
[in,out] | vals | pointer to output values (size: n_elts*dimension) |
void cs_turbulence_function_rij | ( | int | location_id, |
cs_lnum_t | n_elts, | ||
const cs_lnum_t * | elt_ids, | ||
void * | input, | ||
void * | vals | ||
) |
Return or estimate the value of the Reynolds stresses over specified elements.
Returned values are zero for turbulence models other than RANS.
This function matches the cs_eval_at_location_t function profile.
[in] | location_id | base associated mesh location id |
[in] | n_elts | number of associated elements |
[in] | elt_ids | ids of associated elements, or NULL if no filtering is required |
[in,out] | input | ignored |
[in,out] | vals | pointer to output values (size: n_elts*dimension) |
|
extern |
|
extern |
|
extern |
|
extern |
|
extern |
|
extern |
Constant used to define, for each cell \(\Omega_i\), the width of the (implicit) filter:
Useful if and only if iturb = 40 or 41.
|
extern |
Werner and Wengle coefficient
|
extern |
Constant used to define, for each cell \(Omega_i\), the width of the (implicit) filter:
Useful if and only if iturb = 40 or 41.
|
extern |
Werner and Wengle coefficient
|
extern |
Coefficient of turbulent DFM flow model.
|
extern |
Coefficient of turbulent DFM flow model.
|
extern |
Coefficient of turbulent DFM flow model.
|
extern |
Coefficient of turbulent DFM flow model.
|
extern |
Constants of the Cazalbou rotation/curvature correction.
|
extern |
Constants of the Cazalbou rotation/curvature correction.
|
extern |
Constants of the Cazalbou rotation/curvature correction.
|
extern |
Constants of the Cazalbou rotation/curvature correction.
|
extern |
Constants of the Cazalbou rotation/curvature correction.
|
extern |
Constants of the Cazalbou rotation/curvature correction.
|
extern |
Constant \( C_{DDES} \) for the \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST) and hybrid_turb=1.
|
extern |
Van Driest constant appearing in the van Driest damping function applied to the Smagorinsky constant:
Useful if and only if iturb = 40 or 41.
|
extern |
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\)).
|
extern |
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).
|
extern |
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.
|
extern |
Constant of the Rij-epsilon EBRSM.
|
extern |
Constant of the Rij-epsilon EBRSM.
|
extern |
|
extern |
|
extern |
|
extern |
|
extern |
|
extern |
Constant of the Rij-epsilon EBRSM.
|
extern |
Constant of the Rij-epsilon EBRSM.
|
extern |
|
extern |
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).
|
extern |
Constant \(\beta_1\) for the \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST).
|
extern |
Constant \(\beta_2\) for the \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST).
|
extern |
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.
|
extern |
\(\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).
|
extern |
\(\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).
|
extern |
Constant \(\sigma_{k1}\) for the \(k-\omega\) SST model. Useful if and only if iturb=60.
|
extern |
Constant \(\sigma_{k2}\) for the \(k-\omega\) SST model. Useful if and only if iturb=60.
|
extern |
Constant \(\sigma_{\omega 1}\) for the \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST).
|
extern |
Constant \(\sigma_{\omega 2}\) for the \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST).
|
extern |
Constant \(C_\mu\) for all the RANS turbulence models. Warning: different value for v2f models. Useful only for RANS models ( \(k-\varepsilon\), \(R_{ij}-\varepsilon\) or \(k-\omega\)).
|
extern |
\( C_\mu^\frac{1}{4} \)
|
extern |
Constants for the Baglietto et al. quadratic k-epsilon model. Useful if and only if iturb = CS_TURB_K_EPSILON_QUAD
|
extern |
|
extern |
|
extern |
|
extern |
|
extern |
Specific constant of v2f "BL-v2k" (or phi-alpha).
|
extern |
Specific constant of v2f "BL-v2k" (or phi-alpha).
|
extern |
Specific constant of v2f "BL-v2k" (or phi-alpha).
|
extern |
Specific constant of v2f "BL-v2k" (or phi-alpha).
|
extern |
Specific constant of v2f "BL-v2k" (or phi-alpha).
|
extern |
Specific constant of v2f "BL-v2k" (or phi-alpha).
|
extern |
Specific constant of v2f "BL-v2k" (or phi-alpha).
|
extern |
Specific constant of v2f "BL-v2k" (or phi-alpha).
|
extern |
Specific constant of v2f "BL-v2k" (or phi-alpha).
|
extern |
Constant \(C_1\) for the \(R_{ij}-\varepsilon\) LRR model. Useful if and only if iturb=30 ( \(R_{ij}-\varepsilon\) LRR).
|
extern |
|
extern |
Constant \(C_3\) for the \(R_{ij}-\varepsilon\) models. Value is 0.55 for SSG and LRR, 0.6 for EBRSM.
|
extern |
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 cs_turb_rans_model_t::irijec=1 ( \(R_{ij}-\varepsilon\) LRR).
|
extern |
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 cs_turb_rans_model_t::irijec=1 ( \(R_{ij}-\varepsilon\) LRR).
|
extern |
Specific constant of Spalart-Allmaras.
|
extern |
Specific constant of Spalart-Allmaras.
|
extern |
Constant \( C_{SAS}\) for the hybrid \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST) and hybrid_turb=3.
|
extern |
constant \( C_{DDES}\) for the hybrid \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST) and hybrid_turb=3.
|
extern |
Specific constant of Spalart-Allmaras.
|
extern |
Specific constant of Spalart-Allmaras.
|
extern |
Specific constant of Spalart-Allmaras.
|
extern |
Specific constant of Spalart-Allmaras.
|
extern |
Specific constant of Spalart-Allmaras.
|
extern |
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.
|
extern |
Maximum allowed value for the variable \(C\) appearing in the LES dynamic model. Any larger value yielded by the calculation procedure of the dynamic model will be clipped to \( smagmx\).
Useful if and only if iturb = 41.
|
extern |
Minimum allowed value for the variable \(C\) appearing in the LES dynamic model. Any smaller value yielded by the calculation procedure of the dynamic model will be clipped to \( smagmn\).
Useful if and only if iturb = 41.
|
extern |
Constant \(C_s\) for the \(R_{ij}-\varepsilon\) LRR model. Useful if and only if iturb=30 ( \(R_{ij}-\varepsilon\) LRR).
|
extern |
Constant \(C_{\varepsilon 2}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
|
extern |
Constant \(C_{r1}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
|
extern |
Constant \(C_{r2}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
|
extern |
Constant \(C_{r3}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
|
extern |
constant \(C_{r4}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
|
extern |
Constant \(C_{r1}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
|
extern |
Constant \(C_{s1}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
|
extern |
Constant \(C_{s2}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
|
extern |
Constant of the Spalart-Shur rotation/curvature correction.
|
extern |
Constant of the Spalart-Shur rotation/curvature correction.
|
extern |
Constant of the Spalart-Shur rotation/curvature correction.
|
extern |
Constant of logarithmic smooth 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\)).
|
extern |
Constant \( \alpha \) for logarithmic law function switching from rough to smooth: \( \dfrac{1}{\kappa} \ln(y u_k/(\nu + \alpha \xi u_k)) + cstlog \) ( \( \alpha = \exp \left( -\kappa (8.5 - 5.2) \right) \)).
Useful if and only if iturb >= 10 (mixing length, \(k-\varepsilon\), \(R_{ij}-\varepsilon\), LES, v2f or \(k-\omega\)).
|
extern |
Constant of logarithmic rough law function: \( \dfrac{1}{\kappa} \ln(y/\xi) + cstlog_{rough} \) ( \( cstlog_{rough} = 8.5 \)).
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\)).
|
extern |
Constant of GGDH and AFM on the thermal scalar.
|
extern |
Constant of GGDH and AFM on the thermal scalar.
|
extern |
|
extern |
Constant \(a_1\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
|
extern |
Constant \(C_1\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
|
extern |
Constant \(C_2\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
|
extern |
Constant \(C_L\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
|
extern |
Constant \(C_T\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
|
extern |
Constant \(C_{\varepsilon 2}\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
|
extern |
Constant \(C_\eta\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
|
extern |
Constant of the WALE LES method.
|
extern |
Werner and Wengle coefficient
|
extern |
Coefficient of turbulent AFM flow model.
|
extern |
Van Driest constant. (= 25.6)
Useful if and only if cs_glob_wall_functions::iwallf = CS_WALL_F_2SCALES_VDRIEST. (Two scales log law at the wall using Van Driest mixing length expression).
|
extern |
Constant in the expression of Ce1' for the Rij-epsilon EBRSM.
|
extern |
Constant of the Rij-epsilon EBRSM.
|
extern |
Constant of the Rij-epsilon EBRSM.
|
extern |
Constant of the Rij-epsilon EBRSM.
|
extern |
Coefficient of turbulent AFM flow model.
|
extern |
Karman constant. (= 0.42)
Useful if and only if iturb >= 10. (mixing length, \(k-\varepsilon\), \(R_{ij}-\varepsilon\), LES, v2f or \(k-\omega\)).
|
extern |
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.
|
extern |
Constant used in the definition of LES filtering diameter: \( \delta = \text{xlesfl} . (\text{ales} . volume)^{\text{bles}}\) cs_turb_xlesfl is a 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