#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_assert.h"
#include "cs_field.h"
#include "cs_field_default.h"
#include "cs_field_pointer.h"
#include "cs_function.h"
#include "cs_gradient.h"
#include "cs_log.h"
#include "cs_log_iteration.h"
#include "cs_math.h"
#include "cs_mesh.h"
#include "cs_parall.h"
#include "cs_mesh_location.h"
#include "cs_time_step.h"
#include "cs_wall_functions.h"
#include "cs_turbulence_model.h"
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... | |
cs_turb_model_t * | cs_get_glob_turb_model (void) |
Provide write access to turbulence model structure. More... | |
void | cs_set_glob_turb_model (void) |
Set global pointer to turbulence model structure. More... | |
void | cs_turb_compute_constants (int phase_id) |
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 variance_id) |
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... | |
Variables | |
double | cs_turb_xkappa = 0.42 |
double | cs_turb_vdriest = 25.6 |
double | cs_turb_cstlog = 5.2 |
double | cs_turb_cstlog_rough = 8.5 |
double | cs_turb_cstlog_alpha |
double | cs_turb_apow = 8.3 |
double | cs_turb_bpow = 1.0/7.0 |
double | cs_turb_dpow = -1. |
double | cs_turb_cmu = 0.09 |
double | cs_turb_cmu025 = 0.547722557 |
double | cs_turb_ce1 = 1.44 |
double | cs_turb_ce2 = 1.92 |
double | cs_turb_ce4 = 1.20 |
double | cs_turb_crij1 = 1.80 |
double | cs_turb_crij2 = 0.60 |
double | cs_turb_crij_c0 = 3.5 |
double | cs_turb_crij3 = 0.55 |
double | cs_turb_crijp1 = 0.50 |
double | cs_turb_crijp2 = 0.30 |
double | cs_turb_cssge2 = 1.83 |
double | cs_turb_cssgs1 = 1.70 |
double | cs_turb_cssgs2 = -1.05 |
double | cs_turb_cssgr1 = 0.90 |
double | cs_turb_cssgr2 = 0.80 |
double | cs_turb_cssgr3 = 0.65 |
double | cs_turb_cssgr4 = 0.625 |
double | cs_turb_cssgr5 = 0.20 |
double | cs_turb_cebms1 = 1.70 |
double | cs_turb_cebms2 = 0. |
double | cs_turb_cebmr1 = 0.90 |
double | cs_turb_cebmr2 = 0.80 |
double | cs_turb_cebmr3 = 0.65 |
double | cs_turb_cebmr4 = 0.625 |
double | cs_turb_cebmr5 = 0.20 |
double | cs_turb_csrij |
double | cs_turb_cebme2 = 1.83 |
double | cs_turb_cebmmu = 0.22 |
double | cs_turb_xcl = 0.122 |
double | cs_turb_xa1 = 0.1 |
double | cs_turb_xct = 6.0 |
double | cs_turb_xceta = 80.0 |
double | cs_turb_cpale1 = 1.44 |
double | cs_turb_cpale2 = 1.83 |
double | cs_turb_cpale3 = 2.3 |
double | cs_turb_cpale4 = 0.4 |
double | cs_turb_cpalc1 = 1.7 |
double | cs_turb_cpalc2 = 0.9 |
double | cs_turb_cpalct = 4.0 |
double | cs_turb_cpalcl = 0.164 |
double | cs_turb_cpalet = 75.0 |
double | cs_turb_ckwsk1 = 1.0/0.85 |
double | cs_turb_ckwsk2 = 1.0 |
double | cs_turb_ckwsw1 = 2.0 |
double | cs_turb_ckwsw2 = 1.0/0.856 |
double | cs_turb_ckwbt1 = 0.075 |
double | cs_turb_ckwbt2 = 0.0828 |
double | cs_turb_ckwgm1 = -1. |
double | cs_turb_ckwgm2 = -1. |
double | cs_turb_ckwa1 = 0.31 |
double | cs_turb_ckwc1 = 10.0 |
double | cs_turb_cddes = 0.65 |
double | cs_turb_csas = 0.11 |
double | cs_turb_csas_eta2 = 3.51 |
double | cs_turb_chtles_bt0 = 0.48 |
double | cs_turb_csab1 = 0.1355 |
double | cs_turb_csab2 = 0.622 |
double | cs_turb_csasig = 2.0/3.0 |
double | cs_turb_csav1 = 7.1 |
double | cs_turb_csaw1 = -1. |
double | cs_turb_csaw2 = 0.3 |
double | cs_turb_csaw3 = 2.0 |
double | cs_turb_cssr1 = 1.0 |
double | cs_turb_cssr2 = 2.0 |
double | cs_turb_cssr3 = 1.0 |
double | cs_turb_ccaze2 = 1.83 |
double | cs_turb_ccazsc = 0.119 |
double | cs_turb_ccaza = 4.3 |
double | cs_turb_ccazb = 5.130 |
double | cs_turb_ccazc = 0.453 |
double | cs_turb_ccazd = 0.682 |
double | cs_turb_xlesfl = 2.0 |
double | cs_turb_ales = 1.0 |
double | cs_turb_bles = 1.0/3.0 |
double | cs_turb_csmago = 0.065 |
double | cs_turb_xlesfd = 1.5 |
double | cs_turb_csmago_max = -1. |
double | cs_turb_csmago_min = 0. |
double | cs_turb_cdries = 26.0 |
double | cs_turb_cv2fa1 = 0.05 |
double | cs_turb_cv2fe2 = 1.85 |
double | cs_turb_cv2fc1 = 1.4 |
double | cs_turb_cv2fc2 = 0.3 |
double | cs_turb_cv2fct = 6.0 |
double | cs_turb_cv2fcl = 0.25 |
double | cs_turb_cv2fet = 110.0 |
double | cs_turb_cnl1 = 0.8 |
double | cs_turb_cnl2 = 11. |
double | cs_turb_cnl3 = 4.5 |
double | cs_turb_cnl4 = 1e3 |
double | cs_turb_cnl5 = 1. |
double | cs_turb_cwale = 0.25 |
double | cs_turb_xiafm = 0.7 |
double | cs_turb_etaafm = 0.4 |
double | cs_turb_c1trit = 4.15 |
double | cs_turb_c2trit = 0.55 |
double | cs_turb_c3trit = 0.5 |
double | cs_turb_c4trit = 0. |
double | cs_turb_cthafm = 0.236 |
double | cs_turb_cthdfm = 0.31 |
double | cs_turb_cthebdfm = 0.22 |
double | cs_turb_xclt = 0.305 |
Base turbulence model data.
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 | ( | int | phase_id | ) |
Compute turbulence model constants, some of which may depend on the model choice.
[in] | phase_id | turbulent phase id (-1 for single phase flow) |
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) |
double cs_turb_ales = 1.0 |
Constant used to define, for each cell \(\Omega_i\), the width of the (implicit) filter:
Useful if and only if iturb = 40 or 41.
double cs_turb_apow = 8.3 |
Werner and Wengle coefficient
double cs_turb_bles = 1.0/3.0 |
Constant used to define, for each cell \(Omega_i\), the width of the (implicit) filter:
Useful if and only if iturb = 40 or 41.
double cs_turb_bpow = 1.0/7.0 |
Werner and Wengle coefficient
double cs_turb_c1trit = 4.15 |
Coefficient of turbulent DFM flow model.
double cs_turb_c2trit = 0.55 |
Coefficient of turbulent DFM flow model.
double cs_turb_c3trit = 0.5 |
Coefficient of turbulent DFM flow model.
double cs_turb_c4trit = 0. |
Coefficient of turbulent DFM flow model.
double cs_turb_ccaza = 4.3 |
Constants of the Cazalbou rotation/curvature correction.
double cs_turb_ccazb = 5.130 |
Constants of the Cazalbou rotation/curvature correction.
double cs_turb_ccazc = 0.453 |
Constants of the Cazalbou rotation/curvature correction.
double cs_turb_ccazd = 0.682 |
Constants of the Cazalbou rotation/curvature correction.
double cs_turb_ccaze2 = 1.83 |
Constants of the Cazalbou rotation/curvature correction.
double cs_turb_ccazsc = 0.119 |
Constants of the Cazalbou rotation/curvature correction.
double cs_turb_cddes = 0.65 |
Constant \( C_{DDES} \) for the \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST) and hybrid_turb=1.
double cs_turb_cdries = 26.0 |
Van Driest constant appearing in the van Driest damping function applied to the Smagorinsky constant:
Useful if and only if iturb = 40 or 41.
double cs_turb_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\)).
double cs_turb_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).
double cs_turb_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.
double cs_turb_cebme2 = 1.83 |
Constant of the Rij-epsilon EBRSM.
double cs_turb_cebmmu = 0.22 |
Constant of the Rij-epsilon EBRSM.
double cs_turb_cebmr1 = 0.90 |
double cs_turb_cebmr2 = 0.80 |
double cs_turb_cebmr3 = 0.65 |
double cs_turb_cebmr4 = 0.625 |
double cs_turb_cebmr5 = 0.20 |
double cs_turb_cebms1 = 1.70 |
Constant of the Rij-epsilon EBRSM.
double cs_turb_cebms2 = 0. |
Constant of the Rij-epsilon EBRSM.
double cs_turb_chtles_bt0 = 0.48 |
double cs_turb_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).
double cs_turb_ckwbt1 = 0.075 |
Constant \(\beta_1\) for the \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST).
double cs_turb_ckwbt2 = 0.0828 |
Constant \(\beta_2\) for the \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST).
double cs_turb_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 cs_turb_ckwgm1 = -1. |
\(\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).
double cs_turb_ckwgm2 = -1. |
\(\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).
double cs_turb_ckwsk1 = 1.0/0.85 |
Constant \(\sigma_{k1}\) for the \(k-\omega\) SST model. Useful if and only if iturb=60.
double cs_turb_ckwsk2 = 1.0 |
Constant \(\sigma_{k2}\) for the \(k-\omega\) SST model. Useful if and only if iturb=60.
double cs_turb_ckwsw1 = 2.0 |
Constant \(\sigma_{\omega 1}\) for the \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST).
double cs_turb_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).
double cs_turb_cmu = 0.09 |
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\)).
double cs_turb_cmu025 = 0.547722557 |
\( C_\mu^\frac{1}{4} \)
double cs_turb_cnl1 = 0.8 |
Constants for the Baglietto et al. quadratic k-epsilon model. Useful if and only if iturb = CS_TURB_K_EPSILON_QUAD
double cs_turb_cnl2 = 11. |
double cs_turb_cnl3 = 4.5 |
double cs_turb_cnl4 = 1e3 |
double cs_turb_cnl5 = 1. |
double cs_turb_cpalc1 = 1.7 |
Specific constant of v2f "BL-v2k" (or phi-alpha).
double cs_turb_cpalc2 = 0.9 |
Specific constant of v2f "BL-v2k" (or phi-alpha).
double cs_turb_cpalcl = 0.164 |
Specific constant of v2f "BL-v2k" (or phi-alpha).
double cs_turb_cpalct = 4.0 |
Specific constant of v2f "BL-v2k" (or phi-alpha).
double cs_turb_cpale1 = 1.44 |
Specific constant of v2f "BL-v2k" (or phi-alpha).
double cs_turb_cpale2 = 1.83 |
Specific constant of v2f "BL-v2k" (or phi-alpha).
double cs_turb_cpale3 = 2.3 |
Specific constant of v2f "BL-v2k" (or phi-alpha).
double cs_turb_cpale4 = 0.4 |
Specific constant of v2f "BL-v2k" (or phi-alpha).
double cs_turb_cpalet = 75.0 |
Specific constant of v2f "BL-v2k" (or phi-alpha).
double cs_turb_crij1 = 1.80 |
Constant \(C_1\) for the \(R_{ij}-\varepsilon\) LRR model. Useful if and only if iturb=30 ( \(R_{ij}-\varepsilon\) LRR).
double cs_turb_crij2 = 0.60 |
Constant \(C_2\) for the \(R_{ij}-\varepsilon\) LRR model. Useful if and only if iturb=30 ( \(R_{ij}-\varepsilon\) LRR).
double cs_turb_crij3 = 0.55 |
Constant \(C_3\) for the \(R_{ij}-\varepsilon\) models. Value is 0.55 for SSG and LRR, 0.6 for EBRSM.
double cs_turb_crij_c0 = 3.5 |
Rotta constant \(C_0\) for the \(R_{ij}-\varepsilon\) model. Useful for the Lagrangian model. The value is set from \(C_1\) if and only if iturb=CS_TURB_RIJ_EPSILON_LRR ( \(R_{ij}-\varepsilon\) LRR) and \(C_2=0\).
double cs_turb_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 cs_turb_rans_model_t::irijec=1 ( \(R_{ij}-\varepsilon\) LRR).
double cs_turb_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 cs_turb_rans_model_t::irijec=1 ( \(R_{ij}-\varepsilon\) LRR).
double cs_turb_csab1 = 0.1355 |
Specific constant of Spalart-Allmaras.
double cs_turb_csab2 = 0.622 |
Specific constant of Spalart-Allmaras.
double cs_turb_csas = 0.11 |
Constant \( C_{SAS}\) for the hybrid \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST) and hybrid_turb=3.
double cs_turb_csas_eta2 = 3.51 |
constant \( C_{DDES}\) for the hybrid \(k-\omega\) SST model. Useful if and only if iturb=60 ( \(k-\omega\) SST) and hybrid_turb=3.
double cs_turb_csasig = 2.0/3.0 |
Specific constant of Spalart-Allmaras.
double cs_turb_csav1 = 7.1 |
Specific constant of Spalart-Allmaras.
double cs_turb_csaw1 = -1. |
Specific constant of Spalart-Allmaras.
double cs_turb_csaw2 = 0.3 |
Specific constant of Spalart-Allmaras.
double cs_turb_csaw3 = 2.0 |
Specific constant of Spalart-Allmaras.
double cs_turb_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.
double cs_turb_csmago_max = -1. |
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.
double cs_turb_csmago_min = 0. |
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.
double cs_turb_csrij |
Constant \(C_s\) for the \(R_{ij}-\varepsilon\) LRR model. Useful if and only if iturb=30 ( \(R_{ij}-\varepsilon\) LRR).
double cs_turb_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).
double cs_turb_cssgr1 = 0.90 |
Constant \(C_{r1}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
double cs_turb_cssgr2 = 0.80 |
Constant \(C_{r2}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
double cs_turb_cssgr3 = 0.65 |
Constant \(C_{r3}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
double cs_turb_cssgr4 = 0.625 |
constant \(C_{r4}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
double cs_turb_cssgr5 = 0.20 |
Constant \(C_{r1}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
double cs_turb_cssgs1 = 1.70 |
Constant \(C_{s1}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
double cs_turb_cssgs2 = -1.05 |
Constant \(C_{s2}\) for the \(R_{ij}-\varepsilon\) SSG model. Useful if and only if iturb=31 ( \(R_{ij}-\varepsilon\) SSG).
double cs_turb_cssr1 = 1.0 |
Constant of the Spalart-Shur rotation/curvature correction.
double cs_turb_cssr2 = 2.0 |
Constant of the Spalart-Shur rotation/curvature correction.
double cs_turb_cssr3 = 1.0 |
Constant of the Spalart-Shur rotation/curvature correction.
double cs_turb_cstlog = 5.2 |
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\)).
double cs_turb_cstlog_alpha |
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\)).
double cs_turb_cstlog_rough = 8.5 |
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\)).
double cs_turb_cthafm = 0.236 |
Constant of GGDH and AFM on the thermal scalar.
double cs_turb_cthdfm = 0.31 |
Constant of GGDH and AFM on the thermal scalar.
double cs_turb_cthebdfm = 0.22 |
double cs_turb_cv2fa1 = 0.05 |
Constant \(a_1\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
double cs_turb_cv2fc1 = 1.4 |
Constant \(C_1\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
double cs_turb_cv2fc2 = 0.3 |
Constant \(C_2\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
double cs_turb_cv2fcl = 0.25 |
Constant \(C_L\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
double cs_turb_cv2fct = 6.0 |
Constant \(C_T\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
double cs_turb_cv2fe2 = 1.85 |
Constant \(C_{\varepsilon 2}\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
double cs_turb_cv2fet = 110.0 |
Constant \(C_\eta\) for the v2f \(\varphi\)-model. Useful if and only if iturb=50 (v2f \(\varphi\)-model).
double cs_turb_cwale = 0.25 |
Constant of the WALE LES method.
double cs_turb_dpow = -1. |
Werner and Wengle coefficient
double cs_turb_etaafm = 0.4 |
Coefficient of turbulent AFM flow model.
double cs_turb_vdriest = 25.6 |
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).
double cs_turb_xa1 = 0.1 |
Constant in the expression of Ce1' for the Rij-epsilon EBRSM.
double cs_turb_xceta = 80.0 |
Constant of the Rij-epsilon EBRSM.
double cs_turb_xcl = 0.122 |
Constant of the Rij-epsilon EBRSM.
double cs_turb_xclt = 0.305 |
constant of EB-AFM and EB-DFM (0.122*2.5, See F. Dehoux thesis)
double cs_turb_xct = 6.0 |
Constant of the Rij-epsilon EBRSM.
double cs_turb_xiafm = 0.7 |
Coefficient of turbulent AFM flow model.
double cs_turb_xkappa = 0.42 |
Karman constant. (= 0.42)
Useful if and only if iturb >= 10. (mixing length, \(k-\varepsilon\), \(R_{ij}-\varepsilon\), LES, v2f or \(k-\omega\)).
double cs_turb_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.
double cs_turb_xlesfl = 2.0 |
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