1 /*============================================================================
2 * Base turbulence model data.
3 *============================================================================*/
4
5 /*
6 This file is part of Code_Saturne, a general-purpose CFD tool.
7
8 Copyright (C) 1998-2021 EDF S.A.
9
10 This program is free software; you can redistribute it and/or modify it under
11 the terms of the GNU General Public License as published by the Free Software
12 Foundation; either version 2 of the License, or (at your option) any later
13 version.
14
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
18 details.
19
20 You should have received a copy of the GNU General Public License along with
21 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
22 Street, Fifth Floor, Boston, MA 02110-1301, USA.
23 */
24
25 /*----------------------------------------------------------------------------*/
26
27 #include "cs_defs.h"
28
29 /*----------------------------------------------------------------------------*/
30
31 /*----------------------------------------------------------------------------
32 * Standard C library headers
33 *----------------------------------------------------------------------------*/
34
35 #include <assert.h>
36 #include <stdio.h>
37 #include <stdlib.h>
38 #include <string.h>
39 #include <math.h>
40
41 /*----------------------------------------------------------------------------
42 * Local headers
43 *----------------------------------------------------------------------------*/
44
45 #include "bft_mem.h"
46 #include "bft_error.h"
47 #include "bft_printf.h"
48
49 #include "cs_field.h"
50 #include "cs_field_pointer.h"
51 #include "cs_log.h"
52 #include "cs_map.h"
53 #include "cs_parall.h"
54 #include "cs_parameters.h"
55 #include "cs_mesh_location.h"
56 #include "cs_time_step.h"
57 #include "cs_wall_functions.h"
58
59 /*----------------------------------------------------------------------------
60 * Header for the current file
61 *----------------------------------------------------------------------------*/
62
63 #include "cs_turbulence_model.h"
64 #include "cs_math.h"
65 #include "cs_log_iteration.h"
66
67 /*----------------------------------------------------------------------------*/
68
69 BEGIN_C_DECLS
70
71 /*=============================================================================
72 * Additional doxygen documentation
73 *============================================================================*/
74
75 /*!
76 \file cs_turbulence_model.c
77 Base turbulence model data.
78 */
79
80 /*----------------------------------------------------------------------------*/
81
82 /*! \struct cs_turb_model_t
83
84 \brief Turbulence model general options descriptor.
85
86 Members of this turbulence model are publicly accessible, to allow for concise
87 syntax, as it is expected to be used in many places.
88
89 \var cs_turb_model_t::iturb
90 turbulence model
91 - CS_TURB_NONE: no turbulence model (laminar flow)
92 - CS_TURB_MIXING_LENGTH: mixing length model
93 - CS_TURB_K_EPSILON: standard \f$ k-\varepsilon \f$ model
94 - CS_TURB_K_EPSILON_LIN_PROD: \f$ k-\varepsilon \f$ model with
95 Linear Production (LP) correction
96 - CS_TURB_K_EPSILON_LS: Launder-Sharma \f$ k-\varepsilon \f$ model
97 - CS_TURB_K_EPSILON_QUAD: Baglietto et al. quadratic
98 \f$ k-\varepsilon \f$ model
99 - CS_TURB_RIJ_EPSILON_LRR: \f$ R_{ij}-\epsilon \f$ (LRR)
100 - CS_TURB_RIJ_EPSILON_SSG: \f$ R_{ij}-\epsilon \f$ (SSG)
101 - CS_TURB_RIJ_EPSILON_EBRSM: \f$ R_{ij}-\epsilon \f$ (EBRSM)
102 - CS_TURB_LES_SMAGO_CONST: LES (constant Smagorinsky model)
103 - CS_TURB_LES_SMAGO_DYN: LES ("classical" dynamic Smagorisky model)
104 - CS_TURB_LES_WALE: LES (WALE)
105 - CS_TURB_V2F_PHI: v2f phi-model
106 - CS_TURB_V2F_BL_V2K: v2f \f$ BL-v^2-k \f$
107 - CS_TURB_K_OMEGA: \f$ k-\omega \f$ SST
108 - CS_TURB_SPALART_ALLMARAS: Spalart-Allmaras model
109 \var cs_turb_model_t::itytur
110 class of turbulence model (integer value iturb/10, deprecated)
111 \var cs_turb_model_t::hybrid_turb
112 Type of hybrid turbulence model
113 - 0: No model
114 - 1: Detached Eddy Simulation
115 - 2: Delayed Detached Eddy Simulation
116 - 3: Scale Adaptive Model (Menter et al.)
117 \var cs_turb_model_t::type
118 Type of modelling
119 - CS_TURB_NONE: No model
120 - CS_TURB_RANS: RANS
121 - CS_TURB_LES: LES
122 - CS_TURB_HYBRID: Hybrid RANS LES
123 */
124
125 /*----------------------------------------------------------------------------*/
126
127 /*! \struct cs_turb_rans_model_t
128
129 \brief RANS turbulence model descriptor.
130
131 Members of this turbulence model are publicly accessible, to allow for concise
132 syntax, as it is expected to be used in many places.
133
134 \var cs_turb_rans_model_t::irccor
135 activation of rotation/curvature correction for an eddy viscosity
136 turbulence models
137 - 0: false
138 - 1: true
139 \var cs_turb_rans_model_t::itycor
140 type of rotation/curvature correction for an eddy viscosity turbulence
141 models
142 - 1: Cazalbou correction (default when \ref irccor = 1 and
143 \ref cs_turb_model_t::itytur "itytur" = 2 or 5)
144 - 2: Spalart-Shur correction (default when \ref irccor = 1 and
145 \ref iturb = 60 or 70)
146 \var cs_turb_rans_model_t::idirsm
147 turbulent diffusion model for second moment closure
148 - 0: scalar diffusivity (Shir model)
149 - 1: tensorial diffusivity (Daly and Harlow model, default model)
150 \var cs_turb_rans_model_t::iclkep
151 Indicates the clipping method used for \f$k\f$ and
152 \f$\varepsilon\f$, for the \f$k-\epsilon\f$ and v2f models\n
153 - 0: clipping in absolute value
154 - 1: coupled clipping based on physical relationships\n
155 Useful if and only if \ref iturb = 20, 21 or 50 (\f$k-\epsilon\f$ and
156 v2f models). The results obtained with the method corresponding to
157 \ref iclkep =1 showed in some cases a substantial sensitivity to the
158 values of the length scale \ref almax.\n
159 The option \ref iclkep = 1 is therefore not recommended, and,
160 if chosen, must be used cautiously.
161 \var cs_turb_rans_model_t::igrhok
162 Indicates if the term \f$\frac{2}{3}\grad \rho k\f$
163 is taken into account in the velocity equation.
164 - 1: true
165 - 0: false in the velocity
166 Useful if and only if \ref iturb = 20, 21, 50 or 60.\n
167 This term may generate non-physical velocities at the wall.
168 When it is not explicitly taken into account, it is
169 implicitly included into the pressure.
170 \var cs_turb_rans_model_t::igrake
171 Indicates if the terms related to gravity are taken
172 into account in the equations of \f$k-\epsilon\f$.\n
173 - 1: true (default if \f$ \rho \f$ is variable)
174 - 0: false
175 Useful if and only if \ref iturb = 20, 21, 50 or 60 and
176 (\ref cs_physical_constants_t::gravity "gravity")
177 \f$\ne\f$ (0,0,0) and the density is not uniform.
178 \var cs_turb_rans_model_t::igrari
179 Indicates if the terms related to gravity are taken
180 into account in the equations of \f$R_{ij}-\epsilon\f$.\n
181 - 1: true (default if \f$ \rho \f$ is variable)
182 - 0: false
183 Useful if and only if \ref iturb = 30 or 31 and
184 (\ref cs_physical_constants_t::gravity "gravity") \f$\ne\f$
185 (0,0,0) (\f$R_{ij}-\epsilon\f$ model with gravity) and the
186 density is not uniform.
187 \var cs_turb_rans_model_t::ikecou
188 Indicates if the coupling of the source terms of
189 \f$k\f$ and \f$\epsilon\f$ or \f$k\f$ and \f$\omega\f$
190 is taken into account or not.
191 - 1: true,
192 - 0: false\n
193 If \ref ikecou = 0 in \f$k-\epsilon\f$ model, the term
194 in \f$\epsilon\f$ in the equation of \f$k\f$ is made implicit.\n
195 \ref ikecou is initialised to 0 if \ref iturb = 21 or 60, and
196 to 1 if \ref iturb = 20.\n
197 \ref ikecou = 1 is forbidden when using the v2f model (\ref iturb = 50).\n
198 Useful if and only if \ref iturb = 20, 21 or 60 (\f$k-\epsilon\f$ and
199 \f$k-\omega\f$ models)
200 \var cs_turb_rans_model_t::reinit_turb
201 Advanced re-init for EBRSM and k-omega models
202 - 1: true
203 - 0: false (default)
204 \var cs_turb_rans_model_t::irijco
205 coupled solving of Rij
206 - 1: true (default)
207 - 0: false
208 \var cs_turb_rans_model_t::irijnu
209 pseudo eddy viscosity in the matrix of momentum equation to partially
210 implicit \f$ \divv \left( \rho \tens{R} \right) \f$
211 - 1: true
212 - 0: false (default)
213 The goal is to improve the stability of the calculation.
214 The usefulness of \ref irijnu = 1 has however not been
215 clearly demonstrated.\n Since the system is solved in
216 incremental form, this extra turbulent viscosity does
217 not change the final solution for steady flows. However,
218 for unsteady flows, the parameter \ref cs_var_cal_opt_t::nswrsm "nswrsm"
219 should be increased.\n Useful if and only if \ref iturb = 30
220 or 31 (\f$R_{ij}-\epsilon\f$ model).
221 \var cs_turb_rans_model_t::irijrb
222 accurate treatment of \f$ \tens{R} \f$ at the boundary (see \ref condli)
223 - 1: true
224 - 0: false (default)
225 \var cs_turb_rans_model_t::irijec
226 Indicates if the wall echo terms in
227 \f$R_{ij}-\epsilon\f$ LRR model are taken into account:
228 - 1: true,
229 - 0: false (default)\n
230 Useful if and only if \ref iturb = 30 (\f$R_{ij}-\epsilon\f$
231 LRR).\n It is not recommended to take these terms into account:
232 they have an influence only near the walls, their expression is hardly
233 justifiable according to some authors and, in the configurations
234 studied with Code_Saturne, they did not bring any improvement in the
235 results.\n
236 In addition, their use induces an increase in the calculation time.\n
237 The wall echo terms imply the calculation of the distance to the wall
238 for every cell in the domain. See \ref optcal::icdpar "icdpar" for
239 potential
240 restrictions due to this.
241 \var cs_turb_rans_model_t::idifre
242 whole treatment of the diagonal part of the diffusion tensor of \f$
243 \tens{R} \f$ and \f$ \varepsilon \f$
244 - 1: true (default)
245 - 0: simplified treatment
246 \var cs_turb_rans_model_t::iclsyr
247 partial implicitation of symmetry BCs of \f$ \tens{R} \f$
248 - 1: true (default)
249 - 0: false
250 \var cs_turb_rans_model_t::iclptr
251 partial implicitation of wall BCs of \f$ \tens{R} \f$
252 - 1: true
253 - 0: false (default)
254 \var cs_turb_ref_values_t::almax
255 characteristic macroscopic length of the domain, used for the
256 initialization of the turbulence and the potential clipping (with
257 \ref cs_turb_rans_model_t::iclkep "iclkep"= 1)
258 - Negative value: not initialized (the code then uses the cubic root of
259 the domain volume).
260
261 Useful mainly for RANS models.
262 \var cs_turb_ref_values_t::uref
263 characteristic flow velocity, used for the initialization of the
264 turbulence
265 - Negative value: not initialized.
266
267 Useful mainly for RANS models and if
268 the turbulence is not initialized somewhere else (restart file or
269 subroutine \ref cs\_user\_initialization).
270 \var cs_turb_rans_model_t::xlomlg
271 mixing length for the mixing length model
272
273 Useful if and only if \ref iturb= 10 (mixing length).
274 */
275
276 /*----------------------------------------------------------------------------*/
277
278 /*! \struct cs_turb_les_model_t
279
280 \brief LES turbulence model descriptor.
281
282 Members of this turbulence model are publicly accessible, to allow for concise
283 syntax, as it is expected to be used in many places.
284
285 \var cs_turb_les_model_t::idries
286 Activates or the van Driest wall-damping for the
287 Smagorinsky constant (the Smagorinsky constant
288 is multiplied by the damping function
289 \f$1-e^{-y^+/cdries}\f$, where \f$y^+\f$
290 designates the non-dimensional distance to the
291 nearest wall).
292 - 1: true
293 - 0: false
294 The default value is 1 for the Smagorinsky model
295 and 0 for the dynamic model.\n The van Driest
296 wall-damping requires the knowledge of the
297 distance to the nearest wall for each cell
298 in the domain. Refer to keyword \ref optcal::icdpar "icdpar"
299 for potential limitations.\n
300 Useful if and only if \ref iturb = 40 or 41
301 */
302
303 /*----------------------------------------------------------------------------*/
304
305 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
306
307 /*=============================================================================
308 * Macro definitions
309 *============================================================================*/
310
311 /*============================================================================
312 * Type definitions
313 *============================================================================*/
314
315 /* main turbulence model structure and associated pointer */
316
317 static cs_turb_model_t _turb_model =
318 {
319 .iturb = -999,
320 .itytur = -999,
321 .hybrid_turb = 0,
322 .type = -1,
323 .order = -1
324 };
325
326 const cs_turb_model_t *cs_glob_turb_model = NULL;
327
328 /* Reference values for turbulence structure and associated pointer */
329
330 static cs_turb_ref_values_t
331 _turb_ref_values =
332 {
333 .almax = -999,
334 .uref =-1e13
335 };
336
337 const cs_turb_ref_values_t *cs_glob_turb_ref_values = &_turb_ref_values;
338
339 /* RANS turbulence model structure and associated pointer */
340
341 static cs_turb_rans_model_t
342 _turb_rans_model =
343 {
344 .irccor = 0,
345 .itycor = -999,
346 .idirsm = 0,
347 .iclkep = 0,
348 .igrhok = 0,
349 .igrake = 1,
350 .igrari = 1,
351 .ikecou = 0,
352 .reinit_turb= 1,
353 .irijco = 1, /* Coupled version of DRSM models */
354 .irijnu = 0,
355 .irijrb = 0,
356 .irijec = 0,
357 .idifre = 1,
358 .iclsyr = 1,
359 .iclptr = 0,
360 .xlomlg = -1e13
361 };
362
363 const cs_turb_rans_model_t *cs_glob_turb_rans_model = &_turb_rans_model;
364
365 /* LES turbulence model structure and associated pointer */
366
367 static cs_turb_les_model_t _turb_les_model =
368 {
369 .idries = -1,
370 };
371
372 const cs_turb_les_model_t *cs_glob_turb_les_model = &_turb_les_model;
373
374 /*============================================================================
375 * Global variables
376 *============================================================================*/
377
378 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
379
380 /*!
381 * Karman constant. (= 0.42)
382 *
383 * Useful if and only if \ref iturb >= 10.
384 * (mixing length, \f$k-\varepsilon\f$, \f$R_{ij}-\varepsilon\f$,
385 * LES, v2f or \f$k-\omega\f$).
386 */
387 const double cs_turb_xkappa = 0.42;
388
389 /*!
390 * Van Driest constant. (= 25.6)
391 *
392 * Useful if and only if \ref cs_wall_functions_t::iwallf
393 * "cs_glob_wall_functions::iwallf" = CS_WALL_F_2SCALES_VDRIEST.
394 * (Two scales log law at the wall using Van Driest mixing length expression).
395 */
396 const double cs_turb_vdriest = 25.6;
397
398 /*!
399 * Constant of logarithmic smooth law function:
400 * \f$ \dfrac{1}{\kappa} \ln(y^+) + cstlog \f$
401 * (\f$ cstlog = 5.2 \f$).
402 *
403 * Constant of the logarithmic wall function.
404 * Useful if and only if \ref iturb >= 10 (mixing length, \f$k-\varepsilon\f$,
405 * \f$R_{ij}-\varepsilon\f$, LES, v2f or \f$k-\omega\f$).
406 */
407 const double cs_turb_cstlog = 5.2;
408
409 /*!
410 * Constant of logarithmic rough law function:
411 * \f$ \dfrac{1}{\kappa} \ln(y/\xi) + cstlog_{rough} \f$
412 * (\f$ cstlog_{rough} = 8.5 \f$).
413 *
414 * Constant of the logarithmic wall function.
415 * Useful if and only if \ref iturb >= 10 (mixing length, \f$k-\varepsilon\f$,
416 * \f$R_{ij}-\varepsilon\f$, LES, v2f or \f$k-\omega\f$).
417 */
418 const double cs_turb_cstlog_rough = 8.5;
419
420 /*!
421 * Constant \f$ \alpha \f$ for logarithmic law function switching from rough
422 * to smooth:
423 * \f$ \dfrac{1}{\kappa} \ln(y u_k/(\nu + \alpha \xi u_k)) + cstlog \f$
424 * (\f$ \alpha = \exp \left( -\kappa (8.5 - 5.2) \right) \f$).
425 *
426 * Useful if and only if \ref iturb >= 10 (mixing length, \f$k-\varepsilon\f$,
427 * \f$R_{ij}-\varepsilon\f$, LES, v2f or \f$k-\omega\f$).
428 */
429 double cs_turb_cstlog_alpha;
430
431 /*! Werner and Wengle coefficient */
432 const double cs_turb_apow = 8.3;
433
434 /*! Werner and Wengle coefficient */
435 const double cs_turb_bpow = 1.0/7.0;
436
437 /*! Werner and Wengle coefficient */
438 double cs_turb_dpow = -1.;
439
440 /*!
441 * Constant \f$C_\mu\f$ for all the RANS turbulence models except for the
442 * Warning: different value for v2f models. Useful if and only if \ref iturb = 20, 21, 30, 31, 50, 51 or 60
443 * (\f$k-\varepsilon\f$, \f$R_{ij}-\varepsilon\f$ or \f$k-\omega\f$).
444 */
445 double cs_turb_cmu = 0.09;
446
447 /*! \f$ C_\mu^\frac{1}{4} \f$ */
448 double cs_turb_cmu025 = 0.547722557; /* computed more precisely later */
449
450 /*!
451 * Constant \f$C_{\varepsilon 1}\f$ for all the RANS turbulence models except
452 * for the v2f and the \f$k-\omega\f$ models.
453 * Useful if and only if \ref iturb= 20, 21, 30 or 31 (\f$k-\varepsilon\f$
454 * or \f$R_{ij}-\varepsilon\f$).
455 */
456 double cs_turb_ce1 = 1.44;
457
458 /*!
459 * Constant \f$C_{\varepsilon 2}\f$ for the \f$k-\varepsilon\f$ and
460 * \f$R_{ij}-\varepsilon\f$ LRR models.
461 * Useful if and only if \ref iturb = 20, 21 or 30
462 * (\f$k-\varepsilon\f$ or \f$R_{ij}-\varepsilon\f$ LRR).
463 */
464 double cs_turb_ce2 = 1.92;
465
466 /*!
467 * Coefficient of interfacial coefficient in k-eps, used in Lagrange treatment.
468 *
469 * Constant \f$C_{\varepsilon 4}\f$ for the interfacial term (Lagrangian module)
470 * in case of two-way coupling. Useful in case of Lagrangian modelling,
471 * in \f$k-\varepsilon\f$ and \f$R_{ij}-\varepsilon\f$ with two-way coupling.
472 */
473 double cs_turb_ce4 = 1.20;
474
475 /*!
476 * Constant \f$C_1\f$ for the \f$R_{ij}-\varepsilon\f$ LRR model.
477 * Useful if and only if \ref iturb=30 (\f$R_{ij}-\varepsilon\f$ LRR).
478 */
479 double cs_turb_crij1 = 1.80;
480
481 /*
482 * Constant \f$C_2\f$ for the \f$R_{ij}-\varepsilon\f$ LRR model.
483 * Useful if and only if \ref iturb=30 (\f$R_{ij}-\varepsilon\f$ LRR).
484 */
485 double cs_turb_crij2 = 0.60;
486
487 /*!
488 * Constant \f$C_3\f$ for the \f$R_{ij}-\varepsilon\f$ models.
489 * Value is 0.55 for SSG and LRR, 0.6 for EBRSM.
490 */
491 double cs_turb_crij3 = 0.55;
492
493 /*!
494 * Constant \f$C_1^\prime\f$ for the \f$R_{ij}-\varepsilon\f$ LRR model,
495 * corresponding to the wall echo terms.
496 * Useful if and only if \ref iturb=30 and \ref cs_turb_rans_model_t::irijec
497 * "cs_turb_rans_model_t::irijec"=1
498 * (\f$R_{ij}-\varepsilon\f$ LRR).
499 */
500 const double cs_turb_crijp1 = 0.50;
501
502 /*!
503 * Constant \f$C_2^\prime\f$ for the \f$R_{ij}-\varepsilon\f$ LRR model,
504 * corresponding to the wall echo terms.
505 * Useful if and only if \ref iturb=30 and \ref cs_turb_rans_model_t::irijec=1
506 * (\f$R_{ij}-\varepsilon\f$ LRR).
507 */
508 const double cs_turb_crijp2 = 0.30;
509
510 /*!
511 * Constant \f$C_{\varepsilon 2}\f$ for the \f$R_{ij}-\varepsilon\f$ SSG model.
512 * Useful if and only if \ref iturb=31 (\f$R_{ij}-\varepsilon\f$ SSG).
513 */
514 const double cs_turb_cssge2 = 1.83;
515
516 /*!
517 * Constant \f$C_{s1}\f$ for the \f$R_{ij}-\varepsilon\f$ SSG model.
518 * Useful if and only if \ref iturb=31 (\f$R_{ij}-\varepsilon\f$ SSG).
519 */
520 const double cs_turb_cssgs1 = 1.70;
521
522 /*!
523 * Constant \f$C_{s2}\f$ for the \f$R_{ij}-\varepsilon\f$ SSG model.
524 * Useful if and only if \ref iturb=31 (\f$R_{ij}-\varepsilon\f$ SSG).
525 */
526 const double cs_turb_cssgs2 = -1.05;
527
528 /*!
529 * Constant \f$C_{r1}\f$ for the \f$R_{ij}-\varepsilon\f$ SSG model.
530 * Useful if and only if \ref iturb=31 (\f$R_{ij}-\varepsilon\f$ SSG).
531 */
532 const double cs_turb_cssgr1 = 0.90;
533
534 /*!
535 * Constant \f$C_{r2}\f$ for the \f$R_{ij}-\varepsilon\f$ SSG model.
536 * Useful if and only if \ref iturb=31 (\f$R_{ij}-\varepsilon\f$ SSG).
537 */
538 const double cs_turb_cssgr2 = 0.80;
539
540 /*!
541 * Constant \f$C_{r3}\f$ for the \f$R_{ij}-\varepsilon\f$ SSG model.
542 * Useful if and only if \ref iturb=31 (\f$R_{ij}-\varepsilon\f$ SSG).
543 */
544 const double cs_turb_cssgr3 = 0.65;
545
546 /*!
547 * constant \f$C_{r4}\f$ for the \f$R_{ij}-\varepsilon\f$ SSG model.
548 * Useful if and only if \ref iturb=31 (\f$R_{ij}-\varepsilon\f$ SSG).
549 */
550 const double cs_turb_cssgr4 = 0.625;
551
552 /*!
553 * Constant \f$C_{r1}\f$ for the \f$R_{ij}-\varepsilon\f$ SSG model.
554 * Useful if and only if \ref iturb=31 (\f$R_{ij}-\varepsilon\f$ SSG).
555 */
556 const double cs_turb_cssgr5 = 0.20;
557
558 /*! Constant of the Rij-epsilon EBRSM. */
559 const double cs_turb_cebms1 = 1.70;
560
561 /*! Constant of the Rij-epsilon EBRSM. */
562 const double cs_turb_cebms2 = 0.;
563
564 const double cs_turb_cebmr1 = 0.90;
565 const double cs_turb_cebmr2 = 0.80;
566 const double cs_turb_cebmr3 = 0.65;
567 const double cs_turb_cebmr4 = 0.625;
568 const double cs_turb_cebmr5 = 0.20;
569
570 /*!
571 * Constant \f$C_s\f$ for the \f$R_{ij}-\varepsilon\f$ LRR model.
572 * Useful if and only if \ref iturb=30 (\f$R_{ij}-\varepsilon\f$ LRR).
573 */
574 double cs_turb_csrij;
575
576 /*! Constant of the Rij-epsilon EBRSM. */
577 const double cs_turb_cebme2 = 1.83;
578
579 /*! Constant of the Rij-epsilon EBRSM. */
580 const double cs_turb_cebmmu = 0.22;
581
582 /*! Constant of the Rij-epsilon EBRSM. */
583 const double cs_turb_xcl = 0.122;
584
585 /*! Constant in the expression of Ce1' for the Rij-epsilon EBRSM. */
586 const double cs_turb_xa1 = 0.1;
587
588 /*! Constant of the Rij-epsilon EBRSM. */
589 const double cs_turb_xct = 6.0;
590
591 /*! Constant of the Rij-epsilon EBRSM. */
592 const double cs_turb_xceta = 80.0;
593
594 /*! Specific constant of v2f "BL-v2k" (or phi-alpha). */
595 const double cs_turb_cpale1 = 1.44;
596
597 /*! Specific constant of v2f "BL-v2k" (or phi-alpha). */
598 const double cs_turb_cpale2 = 1.83;
599
600 /*! Specific constant of v2f "BL-v2k" (or phi-alpha). */
601 const double cs_turb_cpale3 = 2.3;
602
603 /*! Specific constant of v2f "BL-v2k" (or phi-alpha). */
604 const double cs_turb_cpale4 = 0.4;
605
606 /*! Specific constant of v2f "BL-v2k" (or phi-alpha). */
607 const double cs_turb_cpalc1 = 1.7;
608
609 /*! Specific constant of v2f "BL-v2k" (or phi-alpha). */
610 const double cs_turb_cpalc2 = 0.9;
611
612 /*! Specific constant of v2f "BL-v2k" (or phi-alpha). */
613 const double cs_turb_cpalct = 4.0;
614
615 /*! Specific constant of v2f "BL-v2k" (or phi-alpha). */
616 const double cs_turb_cpalcl = 0.164;
617
618 /*! Specific constant of v2f "BL-v2k" (or phi-alpha). */
619 const double cs_turb_cpalet = 75.0;
620
621 /*!
622 * Constant \f$\sigma_{k1}\f$ for the \f$k-\omega\f$ SST model.
623 * Useful if and only if \ref iturb=60.
624 */
625 const double cs_turb_ckwsk1 = 1.0/0.85;
626
627 /*!
628 * Constant \f$\sigma_{k2}\f$ for the \f$k-\omega\f$ SST model.
629 * Useful if and only if \ref iturb=60.
630 */
631 const double cs_turb_ckwsk2 = 1.0;
632
633 /*!
634 * Constant \f$\sigma_{\omega 1}\f$ for the \f$k-\omega\f$ SST model.
635 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST).
636 */
637 const double cs_turb_ckwsw1 = 2.0;
638
639 /*!
640 * Constant \f$\sigma_{\omega 2}\f$ for the \f$k-\omega\f$ SST model.
641 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST).
642 */
643 const double cs_turb_ckwsw2 = 1.0/0.856;
644
645 /*!
646 * Constant \f$\beta_1\f$ for the \f$k-\omega\f$ SST model.
647 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST).
648 */
649 const double cs_turb_ckwbt1 = 0.075;
650
651 /*!
652 * Constant \f$\beta_2\f$ for the \f$k-\omega\f$ SST model.
653 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST).
654 */
655 const double cs_turb_ckwbt2 = 0.0828;
656
657
658 /*!
659 * \f$\frac{\beta_1}{C_\mu}-\frac{\kappa^2}{\sqrt{C_\mu}\sigma_{\omega 1}}\f$.
660 * Constant \f$\gamma_1\f$ for the \f$k-\omega\f$ SST model.
661 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST).
662 * \warning: \f$\gamma_1\f$ is calculated before the call to
663 * \ref usipsu. Hence, if \f$\beta_1\f$, \f$C_\mu\f$, \f$\kappa\f$ or
664 * \f$\sigma_{\omega 1}\f$ is modified in \ref usipsu,
665 * \ref cs_turb_ckwgm1 must also be modified in accordance.
666 */
667 double cs_turb_ckwgm1 = -1.;
668
669 /*!
670 * \f$\frac{\beta_2}{C_\mu}-\frac{\kappa^2}{\sqrt{C_\mu}\sigma_{\omega 2}}\f$.
671 * Constant \f$\gamma_2\f$ for the \f$k-\omega\f$ SST model.
672 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST).
673 * \warning: \f$\gamma_2\f$ is calculated before the call to \ref usipsu. Hence,
674 * if \f$\beta_2\f$, \f$C_\mu\f$, \f$\kappa\f$ or \f$\sigma_{\omega 2}\f$ is
675 * modified in \ref usipsu, \ref cs_turb_ckwgm2 must also be modified
676 * in accordance.
677 */
678 double cs_turb_ckwgm2 = -1.;
679
680 /*!
681 * Specific constant of k-omega SST.
682 * Constant \f$a_1\f$ for the \f$k-\omega\f$ SST model.
683 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST).
684 */
685 const double cs_turb_ckwa1 = 0.31;
686
687 /*!
688 * Constant \f$ c_1 \f$ for the \f$k-\omega\f$ SST model.
689 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST).
690 * Specific constant of k-omega SST.
691 */
692 const double cs_turb_ckwc1 = 10.0;
693
694 /*!
695 * Constant \f$ C_{DDES} \f$ for the \f$k-\omega\f$ SST model.
696 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST) and hybrid_turb=1.
697 */
698 double cs_turb_cddes = -1.;
699
700 /*!
701 * Constant \f$ C_{SAS}\f$ for the hybrid \f$k-\omega\f$ SST model.
702 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST) and hybrid_turb=3.
703 */
704 const double cs_turb_csas = 0.11;
705
706 /*! constant \f$ C_{DDES}\f$ for the hybrid \f$k-\omega\f$ SST model.
707 * Useful if and only if \ref iturb=60 (\f$k-\omega\f$ SST) and hybrid_turb=3.
708 */
709 const double cs_turb_csas_eta2 = 3.51;
710
711 /*!
712 * Specific constant of Spalart-Allmaras.
713 */
714 const double cs_turb_csab1 = 0.1355;
715
716 /*!
717 * Specific constant of Spalart-Allmaras.
718 */
719 const double cs_turb_csab2 = 0.622;
720
721 /*!
722 * Specific constant of Spalart-Allmaras.
723 */
724 const double cs_turb_csasig = 2.0/3.0;
725
726 /*!
727 * Specific constant of Spalart-Allmaras.
728 */
729 const double cs_turb_csav1 = 7.1;
730
731 /*!
732 * Specific constant of Spalart-Allmaras.
733 */
734 double cs_turb_csaw1 = -1.;
735
736 /*!
737 * Specific constant of Spalart-Allmaras.
738 */
739 const double cs_turb_csaw2 = 0.3;
740
741 /*!
742 * Specific constant of Spalart-Allmaras.
743 */
744 const double cs_turb_csaw3 = 2.0;
745
746 /*!
747 * Constant of the Spalart-Shur rotation/curvature correction.
748 */
749 const double cs_turb_cssr1 = 1.0;
750
751 /*!
752 * Constant of the Spalart-Shur rotation/curvature correction.
753 */
754 const double cs_turb_cssr2 = 2.0;
755
756 /*!
757 * Constant of the Spalart-Shur rotation/curvature correction.
758 */
759 const double cs_turb_cssr3 = 1.0;
760
761 /*!
762 * Constants of the Cazalbou rotation/curvature correction.
763 */
764 const double cs_turb_ccaze2 = 1.83;
765
766 /*!
767 * Constants of the Cazalbou rotation/curvature correction.
768 */
769 const double cs_turb_ccazsc = 0.119;
770
771 /*!
772 * Constants of the Cazalbou rotation/curvature correction.
773 */
774 const double cs_turb_ccaza = 4.3;
775
776 /*!
777 * Constants of the Cazalbou rotation/curvature correction.
778 */
779 const double cs_turb_ccazb = 5.130;
780
781 /*!
782 * Constants of the Cazalbou rotation/curvature correction.
783 */
784 const double cs_turb_ccazc = 0.453;
785
786 /*!
787 * Constants of the Cazalbou rotation/curvature correction.
788 */
789 const double cs_turb_ccazd = 0.682;
790
791 /*!
792 * Constant used in the definition of LES filtering diameter:
793 * \f$ \delta = \text{xlesfl} . (\text{ales} . volume)^{\text{bles}}\f$
794 * \ref cs_turb_xlesfl is a constant used to define, for
795 * each cell \f$\Omega_i\f$, the width of the (implicit) filter:
796 * \f$\overline{\Delta}=xlesfl(ales*|\Omega_i|)^{bles}\f$\n
797 * Useful if and only if \ref iturb = 40 or 41
798 */
799 const double cs_turb_xlesfl = 2.0;
800
801 /*!
802 * Constant used to define, for each cell \f$\Omega_i\f$, the width of
803 * the (implicit) filter:
804 * - \f$\overline{\Delta}=xlesfl(ales*|\Omega_i|)^{bles}\f$
805 *
806 * Useful if and only if \ref iturb = 40 or 41.
807 */
808 const double cs_turb_ales = 1.0;
809
810 /*!
811 * Constant used to define, for each cell \f$Omega_i\f$, the width of
812 * the (implicit) filter:
813 * - \f$\overline{\Delta}=xlesfl(ales*|\Omega_i|)^{bles}\f$
814 *
815 * Useful if and only if \ref iturb = 40 or 41.
816 */
817 const double cs_turb_bles = 1.0/3.0;
818
819 /*!
820 * Smagorinsky constant used in the Smagorinsky model for LES.
821 * The sub-grid scale viscosity is calculated by
822 * \f$\displaystyle\mu_{sg}=
823 * \rho C_{smago}^2\bar{\Delta}^2\sqrt{2\bar{S}_{ij}\bar{S}_{ij}}\f$
824 * where \f$\bar{\Delta}\f$ is the width of the filter
825 * and \f$\bar{S}_{ij}\f$ the filtered strain rate.
826 *
827 * Useful if and only if \ref iturb = 40.
828 * \note In theory Smagorinsky constant is 0.18.
829 * For a channel, 0.065 value is rather taken.
830 */
831 double cs_turb_csmago = 0.065;
832
833 /*!
834 * Ratio between explicit and explicit filter width for a dynamic model.
835 * Constant used to define, for each cell \f$\Omega_i\f$, the width of the
836 * explicit filter used in the framework of the LES dynamic model:
837 * \f$\widetilde{\overline{\Delta}}=xlesfd\overline{\Delta}\f$.
838 *
839 * Useful if and only if \ref iturb = 41.
840 */
841 double cs_turb_xlesfd = 1.5;
842
843 /*!
844 * Maximum allowed value for the variable \f$C\f$ appearing in the LES dynamic
845 * model.
846 * Any larger value yielded by the calculation procedure of the dynamic model
847 * will be clipped to \f$ smagmx\f$.
848 *
849 * Useful if and only if \ref iturb = 41.
850 */
851 double cs_turb_csmago_max = -1.;
852
853 /*!
854 * Minimum allowed value for the variable \f$C\f$ appearing in the LES dynamic
855 * model.
856 * Any smaller value yielded by the calculation procedure of the dynamic model
857 * will be clipped to \f$ smagmn\f$.
858 *
859 * Useful if and only if \ref iturb = 41.
860 */
861 double cs_turb_csmago_min = 0.;
862
863 /*!
864 * Van Driest constant appearing in the van Driest damping function applied to
865 * the Smagorinsky constant:
866 * - \f$ (1-\exp^{(-y^+/cdries}) \f$.
867 *
868 * Useful if and only if \ref iturb = 40 or 41.
869 */
870 double cs_turb_cdries = 26.0;
871
872 /*!
873 * Constant \f$a_1\f$ for the v2f \f$\varphi\f$-model.
874 * Useful if and only if \ref iturb=50 (v2f \f$\varphi\f$-model).
875 */
876 const double cs_turb_cv2fa1 = 0.05;
877
878 /*!
879 * Constant \f$C_{\varepsilon 2}\f$ for the v2f \f$\varphi\f$-model.
880 * Useful if and only if \ref iturb=50 (v2f \f$\varphi\f$-model).
881 */
882 const double cs_turb_cv2fe2 = 1.85;
883
884 /*!
885 * Constant \f$C_1\f$ for the v2f \f$\varphi\f$-model.
886 * Useful if and only if \ref iturb=50 (v2f \f$\varphi\f$-model).
887 */
888 const double cs_turb_cv2fc1 = 1.4;
889
890 /*!
891 * Constant \f$C_2\f$ for the v2f \f$\varphi\f$-model.
892 * Useful if and only if \ref iturb=50 (v2f \f$\varphi\f$-model).
893 */
894 const double cs_turb_cv2fc2 = 0.3;
895
896 /*!
897 * Constant \f$C_T\f$ for the v2f \f$\varphi\f$-model.
898 * Useful if and only if \ref iturb=50 (v2f \f$\varphi\f$-model).
899 */
900 const double cs_turb_cv2fct = 6.0;
901
902 /*!
903 * Constant \f$C_L\f$ for the v2f \f$\varphi\f$-model.
904 * Useful if and only if \ref iturb=50 (v2f \f$\varphi\f$-model).
905 */
906 const double cs_turb_cv2fcl = 0.25;
907
908 /*!
909 * Constant \f$C_\eta\f$ for the v2f \f$\varphi\f$-model.
910 * Useful if and only if \ref iturb=50 (v2f \f$\varphi\f$-model).
911 */
912 const double cs_turb_cv2fet = 110.0;
913
914 /*!
915 * Constants for the Baglietto et al. quadratic k-epsilon model.
916 * Useful if and only if \ref iturb = CS_TURB_K_EPSILON_QUAD
917 */
918 const double cs_turb_cnl1 = 0.8;
919 const double cs_turb_cnl2 = 11.;
920 const double cs_turb_cnl3 = 4.5;
921 const double cs_turb_cnl4 = 1e3;
922 const double cs_turb_cnl5 = 1.;
923
924 /*!
925 * Constant of the WALE LES method.
926 */
927 double cs_turb_cwale = 0.25;
928
929 /*!
930 * Coefficient of turbulent AFM flow model.
931 */
932 const double cs_turb_xiafm = 0.7;
933
934 /*!
935 * Coefficient of turbulent AFM flow model.
936 */
937 const double cs_turb_etaafm = 0.4;
938
939 /*!
940 * Coefficient of turbulent DFM flow model.
941 */
942 const double cs_turb_c1trit = 4.15;
943
944 /*!
945 * Coefficient of turbulent DFM flow model.
946 */
947 const double cs_turb_c2trit = 0.55;
948
949 /*!
950 * Coefficient of turbulent DFM flow model.
951 */
952 const double cs_turb_c3trit= 0.5;
953
954 /*!
955 * Coefficient of turbulent DFM flow model.
956 */
957 const double cs_turb_c4trit = 0.;
958
959 /*!
960 * Constant of GGDH and AFM on the thermal scalar.
961 */
962 const double cs_turb_cthafm = 0.236;
963
964 /*!
965 * Constant of GGDH and AFM on the thermal scalar.
966 */
967 const double cs_turb_cthdfm = 0.31;
968
969 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
970
971 /*============================================================================
972 * Prototypes for functions intended for use only by Fortran wrappers.
973 * (descriptions follow, with function bodies).
974 *============================================================================*/
975
976 void
977 cs_f_turb_model_get_pointers(int **iturb,
978 int **itytur,
979 int **hybrid_turb);
980
981 void
982 cs_f_turb_rans_model_get_pointers(int **irccor,
983 int **itycor,
984 int **idirsm,
985 int **iclkep,
986 int **igrhok,
987 int **igrake,
988 int **igrari,
989 int **ikecou,
990 int **reinit_turb,
991 int **irijco,
992 int **irijnu,
993 int **irijrb,
994 int **irijec,
995 int **idifre,
996 int **iclsyr,
997 int **iclptr);
998
999 void
1000 cs_f_turb_les_model_get_pointers(int **idries);
1001
1002 void
1003 cs_f_turb_reference_values(double **almax,
1004 double **uref,
1005 double **xlomlg);
1006
1007 void
1008 cs_f_turb_model_constants_get_pointers(double **cmu,
1009 double **cmu025,
1010 double **crij1,
1011 double **crij2,
1012 double **csmago,
1013 double **xlesfd,
1014 double **smagmx,
1015 double **smagmn,
1016 double **cwale,
1017 const double **xlesfl,
1018 const double **ales,
1019 const double **bles,
1020 double **cdries
1021 );
1022
1023 /*============================================================================
1024 * Private function definitions
1025 *============================================================================*/
1026
1027 /*============================================================================
1028 * Fortran wrapper function definitions
1029 *============================================================================*/
1030
1031 /*----------------------------------------------------------------------------
1032 * Get pointers to members of the global turbulence model structure.
1033 *
1034 * This function is intended for use by Fortran wrappers, and
1035 * enables mapping to Fortran global pointers.
1036 *
1037 * parameters:
1038 * iturb --> pointer to cs_glob_turb_model->iturb
1039 * itytur --> pointer to cs_glob_turb_model->itytur
1040 * hybrid_turb --> pointer to cs_glob_turb_model->hybrid_turb
1041 *----------------------------------------------------------------------------*/
1042
1043 void
cs_f_turb_model_get_pointers(int ** iturb,int ** itytur,int ** hybrid_turb)1044 cs_f_turb_model_get_pointers(int **iturb,
1045 int **itytur,
1046 int **hybrid_turb)
1047 {
1048 *iturb = &(_turb_model.iturb);
1049 *itytur = &(_turb_model.itytur);
1050 *hybrid_turb = &(_turb_model.hybrid_turb);
1051 }
1052
1053 /*----------------------------------------------------------------------------
1054 * Get pointers to members of the RANS turbulence functions structure.
1055 *
1056 * This function is intended for use by Fortran wrappers, and
1057 * enables mapping to Fortran global pointers.
1058 *
1059 * parameters:
1060 * irccor --> pointer to cs_glob_turb_rans_model->irccor
1061 * itycor --> pointer to cs_glob_turb_rans_model->itycor
1062 * idirsm --> pointer to cs_glob_turb_rans_model->idirsm
1063 * iclkep --> pointer to cs_glob_turb_rans_model->iclkep
1064 * igrhok --> pointer to cs_glob_turb_rans_model->igrhok
1065 * igrake --> pointer to cs_glob_turb_rans_model->igrake
1066 * igrari --> pointer to cs_glob_turb_rans_model->igrari
1067 * ikecou --> pointer to cs_glob_turb_rans_model->ikecou
1068 * reinit_turb --> pointer to cs_glob_turb_rans_model->reinit_turb
1069 * irijco --> pointer to cs_glob_turb_rans_model->irijco
1070 * irijnu --> pointer to cs_glob_turb_rans_model->irijnu
1071 * irijrb --> pointer to cs_glob_turb_rans_model->irijrb
1072 * irijec --> pointer to cs_glob_turb_rans_model->irijec
1073 * idifre --> pointer to cs_glob_turb_rans_model->idifre
1074 * iclsyr --> pointer to cs_glob_turb_rans_model->iclsyr
1075 * iclptr --> pointer to cs_glob_turb_rans_model->iclptr
1076 *----------------------------------------------------------------------------*/
1077
1078 void
cs_f_turb_rans_model_get_pointers(int ** irccor,int ** itycor,int ** idirsm,int ** iclkep,int ** igrhok,int ** igrake,int ** igrari,int ** ikecou,int ** reinit_turb,int ** irijco,int ** irijnu,int ** irijrb,int ** irijec,int ** idifre,int ** iclsyr,int ** iclptr)1079 cs_f_turb_rans_model_get_pointers(int **irccor,
1080 int **itycor,
1081 int **idirsm,
1082 int **iclkep,
1083 int **igrhok,
1084 int **igrake,
1085 int **igrari,
1086 int **ikecou,
1087 int **reinit_turb,
1088 int **irijco,
1089 int **irijnu,
1090 int **irijrb,
1091 int **irijec,
1092 int **idifre,
1093 int **iclsyr,
1094 int **iclptr)
1095 {
1096 *irccor = &(_turb_rans_model.irccor);
1097 *itycor = &(_turb_rans_model.itycor);
1098 *idirsm = &(_turb_rans_model.idirsm);
1099 *iclkep = &(_turb_rans_model.iclkep);
1100 *igrhok = &(_turb_rans_model.igrhok);
1101 *igrake = &(_turb_rans_model.igrake);
1102 *igrari = &(_turb_rans_model.igrari);
1103 *ikecou = &(_turb_rans_model.ikecou);
1104 *reinit_turb= &(_turb_rans_model.reinit_turb);
1105 *irijco = &(_turb_rans_model.irijco);
1106 *irijnu = &(_turb_rans_model.irijnu);
1107 *irijrb = &(_turb_rans_model.irijrb);
1108 *irijec = &(_turb_rans_model.irijec);
1109 *idifre = &(_turb_rans_model.idifre);
1110 *iclsyr = &(_turb_rans_model.iclsyr);
1111 *iclptr = &(_turb_rans_model.iclptr);
1112 }
1113
1114 /*----------------------------------------------------------------------------
1115 * Get pointers to members of the LES turbulence model structure.
1116 *
1117 * This function is intended for use by Fortran wrappers, and
1118 * enables mapping to Fortran global pointers.
1119 *
1120 * parameters:
1121 * idries --> pointer to cs_glob_turb_les_model->idries
1122 *----------------------------------------------------------------------------*/
1123
1124 void
cs_f_turb_les_model_get_pointers(int ** idries)1125 cs_f_turb_les_model_get_pointers(int **idries)
1126 {
1127 *idries = &(_turb_les_model.idries);
1128 }
1129
1130 /*----------------------------------------------------------------------------
1131 * Get pointers to members of the RANS turbulence functions structure.
1132 *
1133 * This function is intended for use by Fortran wrappers, and
1134 * enables mapping to Fortran global pointers.
1135 *
1136 * parameters:
1137 * almax --> pointer to cs_glob_turb_ref_values->almax
1138 * uref --> pointer to cs_glob_turb_ref_values->uref
1139 * xlomlg --> pointer to cs_glob_turb_rans_model->xlomlg
1140 *----------------------------------------------------------------------------*/
1141
1142 void
cs_f_turb_reference_values(double ** almax,double ** uref,double ** xlomlg)1143 cs_f_turb_reference_values(double **almax,
1144 double **uref,
1145 double **xlomlg)
1146 {
1147 *almax = &(_turb_ref_values.almax);
1148 *uref = &(_turb_ref_values.uref);
1149 *xlomlg = &(_turb_rans_model.xlomlg);
1150 }
1151
1152 /*----------------------------------------------------------------------------
1153 * Get pointers to constants for turbulence models.
1154 *
1155 * This function is intended for use by Fortran wrappers, and
1156 * enables mapping to Fortran global pointers.
1157 *
1158 * parameters:
1159 * cmu --> pointer to cs_turb_cmu
1160 * cmu025 --> pointer to cs_turb_cmu025
1161 * crij1 --> pointer to cs_turb_crij1
1162 * crij2 --> pointer to cs_turb_crij2
1163 * csmago --> pointer to cs_turb_csmago
1164 * xlesfd --> pointer to cs_turb_xlesfd
1165 * smagmx --> pointer to cs_turb_smago_max
1166 * smagmn --> pointer to cs_turb_smago_min
1167 * cwale --> pointer to cs_turb_cwale
1168 * xlesfl --> pointer to cs_turb_xlesfl
1169 * ales --> pointer to cs_turb_ales
1170 * bles --> pointer to cs_turb_bles
1171 * cdries --> pointer to cs_turb_cdries
1172 *----------------------------------------------------------------------------*/
1173
1174 void
cs_f_turb_model_constants_get_pointers(double ** cmu,double ** cmu025,double ** crij1,double ** crij2,double ** csmago,double ** xlesfd,double ** smagmx,double ** smagmn,double ** cwale,const double ** xlesfl,const double ** ales,const double ** bles,double ** cdries)1175 cs_f_turb_model_constants_get_pointers(double **cmu,
1176 double **cmu025,
1177 double **crij1,
1178 double **crij2,
1179 double **csmago,
1180 double **xlesfd,
1181 double **smagmx,
1182 double **smagmn,
1183 double **cwale,
1184 const double **xlesfl,
1185 const double **ales,
1186 const double **bles,
1187 double **cdries)
1188 {
1189 *cmu = &cs_turb_cmu;
1190 *cmu025 = &cs_turb_cmu025;
1191 *crij1 = &cs_turb_crij1;
1192 *crij2 = &cs_turb_crij2;
1193 *csmago= &cs_turb_csmago;
1194 *xlesfd= &cs_turb_xlesfd;
1195 *csmago= &cs_turb_csmago;
1196 *smagmx= &cs_turb_csmago_max;
1197 *smagmn= &cs_turb_csmago_min;
1198 *cwale = &cs_turb_cwale;
1199 *xlesfl= &cs_turb_xlesfl;
1200 *ales = &cs_turb_ales;
1201 *bles = &cs_turb_bles;
1202 *cdries= &cs_turb_cdries;
1203 }
1204
1205 /*============================================================================
1206 * Private function definitions
1207 *============================================================================*/
1208
1209 /*----------------------------------------------------------------------------
1210 * Return enumeration name associated with a turbulence model value
1211 *
1212 * parameters:
1213 * id <-- model type
1214 *
1215 * returns:
1216 * pointer to enum name.
1217 *----------------------------------------------------------------------------*/
1218
1219 static const char *
_turbulence_model_enum_name(cs_turb_model_type_t id)1220 _turbulence_model_enum_name(cs_turb_model_type_t id)
1221 {
1222 const char *s = NULL;
1223 switch(id) {
1224 case CS_TURB_NONE:
1225 s = "CS_TURB_NONE";
1226 break;
1227 case CS_TURB_MIXING_LENGTH:
1228 s = "CS_TURB_MIXING_LENGTH";
1229 break;
1230 case CS_TURB_K_EPSILON:
1231 s = "CS_TURB_K_EPSILON";
1232 break;
1233 case CS_TURB_K_EPSILON_LIN_PROD:
1234 s = "CS_TURB_K_EPSILON_LIN_PROD";
1235 break;
1236 case CS_TURB_K_EPSILON_LS:
1237 s = "CS_TURB_K_EPSILON_LS";
1238 break;
1239 case CS_TURB_K_EPSILON_QUAD:
1240 s = "CS_TURB_K_EPSILON_QUAD";
1241 break;
1242 case CS_TURB_RIJ_EPSILON_LRR:
1243 s = "CS_TURB_RIJ_EPSILON_LRR";
1244 break;
1245 case CS_TURB_RIJ_EPSILON_SSG:
1246 s = "CS_TURB_RIJ_EPSILON_SSG";
1247 break;
1248 case CS_TURB_RIJ_EPSILON_EBRSM:
1249 s = "CS_TURB_RIJ_EPSILON_EBRSM";
1250 break;
1251 case CS_TURB_LES_SMAGO_CONST:
1252 s = "CS_TURB_LES_SMAGO_CONST";
1253 break;
1254 case CS_TURB_LES_SMAGO_DYN:
1255 s = "CS_TURB_LES_SMAGO_DYN";
1256 break;
1257 case CS_TURB_LES_WALE:
1258 s = "CS_TURB_LES_WALE";
1259 break;
1260 case CS_TURB_V2F_PHI:
1261 s = "CS_TURB_V2F_PHI";
1262 break;
1263 case CS_TURB_V2F_BL_V2K:
1264 s = "CS_TURB_V2F_BL_V2K";
1265 break;
1266 case CS_TURB_K_OMEGA:
1267 s = "CS_TURB_K_OMEGA";
1268 break;
1269 case CS_TURB_SPALART_ALLMARAS:
1270 s = "CS_TURB_SPALART_ALLMARAS";
1271 break;
1272 default:
1273 bft_error(__FILE__, __LINE__, 0,
1274 _("Unknown cs_turb_model_type_t value: %d"), id);
1275 }
1276
1277 return s;
1278 }
1279
1280 /*----------------------------------------------------------------------------
1281 * Return associated with a turbulence model value
1282 *
1283 * parameters:
1284 * id <-- model type
1285 *
1286 * returns:
1287 * pointer to enum name.
1288 *----------------------------------------------------------------------------*/
1289
1290 static const char *
_turbulence_model_name(cs_turb_model_type_t id)1291 _turbulence_model_name(cs_turb_model_type_t id)
1292 {
1293 const char *s = NULL;
1294 switch(id) {
1295 case CS_TURB_NONE:
1296 s = _("no turbulence model");
1297 break;
1298 case CS_TURB_MIXING_LENGTH:
1299 s = _("mixing length model");
1300 break;
1301 case CS_TURB_K_EPSILON:
1302 s = _("standard k-epsilon model");
1303 break;
1304 case CS_TURB_K_EPSILON_LIN_PROD:
1305 s = _("k-epsilon model with Linear Production (LP) correction");
1306 break;
1307 case CS_TURB_K_EPSILON_LS:
1308 s = _("Launder-Sharma k-epsilon model");
1309 break;
1310 case CS_TURB_K_EPSILON_QUAD:
1311 s = _("Baglietto et al. quadratic k-epsilon model");
1312 break;
1313 case CS_TURB_RIJ_EPSILON_LRR:
1314 s = _("Rij-epsilon (LRR) model");
1315 break;
1316 case CS_TURB_RIJ_EPSILON_SSG:
1317 s = _("Rij-epsilon (SSG)");
1318 break;
1319 case CS_TURB_RIJ_EPSILON_EBRSM:
1320 s = _("Rij-epsilon (EBRSM))");
1321 break;
1322 case CS_TURB_LES_SMAGO_CONST:
1323 s = _("LES (constant Smagorinsky model)");
1324 break;
1325 case CS_TURB_LES_SMAGO_DYN:
1326 s = _("LES (classical dynamic Smagorisky model)");
1327 break;
1328 case CS_TURB_LES_WALE:
1329 s = _("LES (WALE)");
1330 break;
1331 case CS_TURB_V2F_PHI:
1332 s = _("v2f phi-model");
1333 break;
1334 case CS_TURB_V2F_BL_V2K:
1335 s = _("v2f BL-v2-k)");
1336 break;
1337 case CS_TURB_K_OMEGA:
1338 s = _("k-omega SST");
1339 break;
1340 case CS_TURB_SPALART_ALLMARAS:
1341 s = _("Spalart-Allmaras model");
1342 break;
1343 default:
1344 bft_error(__FILE__, __LINE__, 0,
1345 _("Unknown cs_turb_model_type_t value: %d"), id);
1346 }
1347
1348 return s;
1349 }
1350
1351 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
1352
1353 /*=============================================================================
1354 * Public function definitions
1355 *============================================================================*/
1356
1357 /*----------------------------------------------------------------------------*/
1358 /*!
1359 * \brief Initialize turbulence model structures
1360 */
1361 /*----------------------------------------------------------------------------*/
1362
1363 void
cs_turb_model_init(void)1364 cs_turb_model_init(void) {
1365
1366 /* set global pointer to turbulence model */
1367 cs_set_glob_turb_model();
1368
1369 }
1370
1371 /*----------------------------------------------------------------------------*/
1372 /*!
1373 * \brief Initialize type and order members of turbulence model structure
1374 */
1375 /*----------------------------------------------------------------------------*/
1376
1377 void
cs_set_type_order_turbulence_model(void)1378 cs_set_type_order_turbulence_model(void)
1379 {
1380 _turb_model.type = CS_TURB_NONE;
1381 if (_turb_model.iturb == CS_TURB_MIXING_LENGTH) {
1382 _turb_model.type = CS_TURB_RANS;
1383 _turb_model.order = CS_TURB_ALGEBRAIC;
1384 }
1385 else if ( _turb_model.iturb == CS_TURB_K_EPSILON
1386 || _turb_model.iturb == CS_TURB_K_EPSILON_LIN_PROD
1387 || _turb_model.iturb == CS_TURB_K_EPSILON_LS
1388 || _turb_model.iturb == CS_TURB_K_EPSILON_QUAD
1389 || _turb_model.iturb == CS_TURB_V2F_PHI
1390 || _turb_model.iturb == CS_TURB_V2F_BL_V2K
1391 || _turb_model.iturb == CS_TURB_K_OMEGA
1392 || _turb_model.iturb == CS_TURB_SPALART_ALLMARAS) {
1393 _turb_model.type = CS_TURB_RANS;
1394 _turb_model.order = CS_TURB_FIRST_ORDER;
1395 }
1396 else if ( _turb_model.iturb == CS_TURB_RIJ_EPSILON_LRR
1397 || _turb_model.iturb == CS_TURB_RIJ_EPSILON_SSG
1398 || _turb_model.iturb == CS_TURB_RIJ_EPSILON_EBRSM) {
1399 _turb_model.type = CS_TURB_RANS;
1400 _turb_model.order = CS_TURB_SECOND_ORDER;
1401 }
1402 else if ( _turb_model.iturb == CS_TURB_LES_SMAGO_CONST
1403 || _turb_model.iturb == CS_TURB_LES_SMAGO_DYN
1404 || _turb_model.iturb == CS_TURB_LES_WALE) {
1405 _turb_model.type = CS_TURB_LES;
1406 _turb_model.order = CS_TURB_ALGEBRAIC;
1407 }
1408
1409 else {
1410 _turb_model.iturb = 0;
1411 _turb_model.itytur = CS_TURB_TYPE_NONE;
1412 }
1413 }
1414
1415 /*----------------------------------------------------------------------------*/
1416 /*!
1417 * \brief Provide write access to turbulence model structure
1418 */
1419 /*----------------------------------------------------------------------------*/
1420
1421 cs_turb_model_t *
cs_get_glob_turb_model(void)1422 cs_get_glob_turb_model(void)
1423 {
1424 return &_turb_model;
1425 }
1426
1427 /*----------------------------------------------------------------------------*/
1428 /*!
1429 * \brief Set global pointer to turbulence model structure
1430 *
1431 * This global pointer provides a read-only access to the structure.
1432 */
1433 /*----------------------------------------------------------------------------*/
1434
1435 void
cs_set_glob_turb_model(void)1436 cs_set_glob_turb_model(void)
1437 {
1438 /* If not set yet, points to the locally defined structure */
1439 if (cs_glob_turb_model == NULL)
1440 cs_glob_turb_model = &_turb_model;
1441 }
1442
1443 /*----------------------------------------------------------------------------*/
1444 /*!
1445 * \brief Compute turbulence model constants,
1446 * some of which may depend on the model choice.
1447 */
1448 /*----------------------------------------------------------------------------*/
1449
1450 void
cs_turb_compute_constants(void)1451 cs_turb_compute_constants(void)
1452 {
1453 cs_turb_dpow = 1./(1.+cs_turb_bpow);
1454 cs_turb_cmu025 = pow(cs_turb_cmu, 0.25);
1455 cs_turb_cstlog_alpha = exp(-cs_turb_xkappa
1456 * (cs_turb_cstlog_rough - cs_turb_cstlog));
1457
1458 int k_turb_schmidt = cs_field_key_id("turbulent_schmidt");
1459
1460 cs_field_pointer_ensure_init();
1461
1462 if (CS_F_(k) != NULL)
1463 cs_field_set_key_double(CS_F_(k), k_turb_schmidt, 1.);
1464
1465 if (CS_F_(phi) != NULL)
1466 cs_field_set_key_double(CS_F_(phi), k_turb_schmidt, 1.);
1467
1468 if ( cs_glob_turb_model->iturb == CS_TURB_RIJ_EPSILON_LRR
1469 || cs_glob_turb_model->iturb == CS_TURB_RIJ_EPSILON_SSG)
1470 cs_field_set_key_double(CS_F_(eps), k_turb_schmidt, 1.22);
1471 else if (cs_glob_turb_model->iturb == CS_TURB_RIJ_EPSILON_EBRSM)
1472 cs_field_set_key_double(CS_F_(eps), k_turb_schmidt, 1.15);
1473 else if (cs_glob_turb_model->iturb == CS_TURB_V2F_BL_V2K)
1474 cs_field_set_key_double(CS_F_(eps), k_turb_schmidt, 1.5);
1475 else
1476 cs_field_set_key_double(CS_F_(eps), k_turb_schmidt, 1.30);
1477
1478 if (cs_glob_turb_model->iturb == CS_TURB_RIJ_EPSILON_EBRSM)
1479 cs_turb_csrij = 0.21;
1480 else
1481 cs_turb_csrij = 0.22;
1482
1483 if (cs_glob_turb_model->iturb == CS_TURB_K_OMEGA){
1484 /* SST DDES */
1485 if (cs_glob_turb_model->hybrid_turb == 2)
1486 cs_turb_cddes = 0.65;
1487 else if (cs_glob_turb_model->hybrid_turb == 1)
1488 cs_turb_cddes = 0.61;
1489 }
1490 else if (cs_glob_turb_model->iturb == CS_TURB_V2F_BL_V2K) {
1491 cs_turb_cddes = 0.60;
1492 }
1493
1494 double xkappa2 = cs_turb_xkappa*cs_turb_xkappa;
1495 cs_turb_ckwgm1 = cs_turb_ckwbt1/cs_turb_cmu
1496 - xkappa2/(cs_turb_ckwsw1*sqrt(cs_turb_cmu));
1497 cs_turb_ckwgm2 = cs_turb_ckwbt2/cs_turb_cmu
1498 - xkappa2/(cs_turb_ckwsw2*sqrt(cs_turb_cmu));
1499 cs_turb_csaw1 = cs_turb_csab1/xkappa2
1500 + 1./cs_turb_csasig*(1. + cs_turb_csab2);
1501 cs_turb_csmago_max = cs_turb_csmago*cs_turb_csmago;
1502 cs_turb_csmago_min = 0.;
1503
1504 /* LRR constants */
1505 cs_turb_crij1 = 1.80;
1506 cs_turb_crij2 = 0.60;
1507 }
1508
1509 /*----------------------------------------------------------------------------*/
1510 /*!
1511 * \brief Provide access to cs_glob_turb_ref_values
1512 *
1513 * needed to initialize structure with GUI
1514 */
1515 /*----------------------------------------------------------------------------*/
1516
1517 cs_turb_ref_values_t *
cs_get_glob_turb_ref_values(void)1518 cs_get_glob_turb_ref_values(void)
1519 {
1520 return &_turb_ref_values;
1521 }
1522
1523 /*----------------------------------------------------------------------------*/
1524 /*!
1525 * \brief Provide access to cs_glob_turb_rans_model
1526 *
1527 * needed to initialize structure with GUI
1528 */
1529 /*----------------------------------------------------------------------------*/
1530
1531 cs_turb_rans_model_t *
cs_get_glob_turb_rans_model(void)1532 cs_get_glob_turb_rans_model(void)
1533 {
1534 return &_turb_rans_model;
1535 }
1536
1537 /*----------------------------------------------------------------------------*/
1538 /*!
1539 * \brief Provide access to cs_glob_turb_les_model
1540 *
1541 * needed to initialize structure with GUI
1542 */
1543 /*----------------------------------------------------------------------------*/
1544
1545 cs_turb_les_model_t *
cs_get_glob_turb_les_model(void)1546 cs_get_glob_turb_les_model(void)
1547 {
1548 return &_turb_les_model;
1549 }
1550
1551 /*----------------------------------------------------------------------------*/
1552 /*!
1553 * \brief Print the turbulence model parameters to setup.log.
1554 */
1555 /*----------------------------------------------------------------------------*/
1556
1557 void
cs_turb_model_log_setup(void)1558 cs_turb_model_log_setup(void)
1559 {
1560 if (cs_glob_turb_model == NULL)
1561 return;
1562
1563 const cs_turb_model_t *turb_model = cs_glob_turb_model;
1564 const cs_wall_functions_t *wall_fns = cs_get_glob_wall_functions();
1565
1566
1567 cs_var_cal_opt_t var_cal_opt;
1568 int key_cal_opt_id = cs_field_key_id("var_cal_opt");
1569
1570 cs_log_printf(CS_LOG_SETUP, _("\n"
1571 "Turbulence model options\n"
1572 "------------------------\n\n"
1573 " Continuous phase:\n\n"));
1574
1575 if (turb_model->type == CS_TURB_RANS)
1576 cs_log_printf
1577 (CS_LOG_SETUP,
1578 _(" RANS model (type = CS_TURB_RANS)\n"));
1579 else if (turb_model->type == CS_TURB_LES)
1580 cs_log_printf
1581 (CS_LOG_SETUP,
1582 _(" LES model (type = CS_TURB_LES)\n"));
1583 else if (turb_model->order == CS_TURB_ALGEBRAIC)
1584 cs_log_printf
1585 (CS_LOG_SETUP,
1586 _(" Algebraic model (order = CS_TURB_ALGEBRAIC)\n"));
1587 else if (turb_model->order == CS_TURB_FIRST_ORDER)
1588 cs_log_printf
1589 (CS_LOG_SETUP,
1590 _(" First order model (order = CS_TURB_FIRST_ORDER)\n"));
1591 else if (turb_model->order == CS_TURB_SECOND_ORDER)
1592 cs_log_printf
1593 (CS_LOG_SETUP,
1594 _(" Second order model (order = CS_TURB_SECOND_ORDER)\n"));
1595
1596 cs_log_printf
1597 (CS_LOG_SETUP,
1598 _("\n %s\n"
1599 " (iturb = %s)\n\n"),
1600 _turbulence_model_name(turb_model->iturb),
1601 _turbulence_model_enum_name(turb_model->iturb));
1602
1603 const char *iwallf_value_str[]
1604 = {N_("Disabled"),
1605 N_("One scale power law, forbidden for k-epsilon"),
1606 N_("One scale log law"),
1607 N_("Two scales log law"),
1608 N_("Scalable wall function"),
1609 N_("Two scales Van Driest"),
1610 N_("Two scales smooth/rough"),
1611 N_("All y+")};
1612
1613 if (wall_fns->iwallf >= 0)
1614 cs_log_printf
1615 (CS_LOG_SETUP,
1616 _(" iwallf (wall function:\n"
1617 " %s)\n"),
1618 iwallf_value_str[wall_fns->iwallf]);
1619
1620 const char *iwalfs_value_str[]
1621 = {N_("Arparci and Larsen"),
1622 N_("Van Driest"),
1623 N_("Louis (atmo flows)"),
1624 N_("Monin Obukhov (atmo flows)"),
1625 N_("smooth/rough")};
1626
1627 if (wall_fns->iwalfs >= 0)
1628 cs_log_printf
1629 (CS_LOG_SETUP,
1630 _(" iwalfs (Scalar wall function:\n"
1631 " %s)\n"),
1632 iwalfs_value_str[wall_fns->iwalfs]);
1633
1634 cs_log_printf(CS_LOG_SETUP,
1635 _(" ypluli: %14.5e (Limit Y+)\n\n"),
1636 wall_fns->ypluli);
1637
1638 const char *igrhok_value_str[]
1639 = {N_("0 (ignore Grad(rho k) in velocity equation)"),
1640 N_("1 (Grad(rho k) in velocity equation)")};
1641
1642 cs_log_printf(CS_LOG_SETUP,
1643 _(" igrhok: %s\n"),
1644 _(igrhok_value_str[cs_glob_turb_rans_model->igrhok]));
1645
1646 if (turb_model->iturb == CS_TURB_MIXING_LENGTH) {
1647
1648 cs_log_printf(CS_LOG_SETUP,
1649 _(" xlomlg: %14.5e (Characteristic length)\n"),
1650 cs_glob_turb_rans_model->xlomlg);
1651
1652 }
1653 else if ( turb_model->iturb == CS_TURB_K_EPSILON
1654 || turb_model->iturb == CS_TURB_K_EPSILON_LIN_PROD
1655 || turb_model->iturb == CS_TURB_K_EPSILON_LS
1656 || turb_model->iturb == CS_TURB_K_EPSILON_QUAD) {
1657
1658 cs_log_printf
1659 (CS_LOG_SETUP,
1660 _(" almax: %14.5e (Characteristic length)\n"
1661 " uref: %14.5e (Characteristic velocity)\n"
1662 " iclkep: %14d (k-epsilon clipping model)\n"
1663 " ikecou: %14d (k-epsilon coupling mode)\n"
1664 " igrake: %14d (Account for gravity)\n"),
1665 cs_glob_turb_ref_values->almax,
1666 cs_glob_turb_ref_values->uref,
1667 cs_glob_turb_rans_model->iclkep,
1668 cs_glob_turb_rans_model->ikecou,
1669 cs_glob_turb_rans_model->igrake);
1670
1671 if ( cs_glob_turb_rans_model->ikecou == 0
1672 && cs_glob_time_step_options->idtvar >= 0) {
1673 cs_field_get_key_struct(CS_F_(k), key_cal_opt_id, &var_cal_opt);
1674 cs_real_t relaxvk = var_cal_opt.relaxv;
1675 cs_field_get_key_struct(CS_F_(eps), key_cal_opt_id, &var_cal_opt);
1676 cs_real_t relaxve = var_cal_opt.relaxv;
1677 cs_log_printf
1678 (CS_LOG_SETUP,
1679 _(" relaxv: %14.5e for k (Relaxation)\n"
1680 " relaxv: %14.5e for epsilon (Relaxation)\n"),
1681 relaxvk, relaxve);
1682 }
1683 else
1684 cs_log_printf(CS_LOG_SETUP, _("\n"));
1685
1686 }
1687 else if ( turb_model->iturb == CS_TURB_RIJ_EPSILON_LRR
1688 || turb_model->iturb == CS_TURB_RIJ_EPSILON_SSG
1689 || turb_model->iturb == CS_TURB_RIJ_EPSILON_EBRSM) {
1690
1691 cs_log_printf(CS_LOG_SETUP,
1692 _(" almax: %14.5e (Characteristic length)\n"
1693 " uref: %14.5e (Characteristic velocity)\n"
1694 " irijco: %14d (Coupled resolution)\n"
1695 " irijnu: %14d (Matrix stabilization)\n"
1696 " irijrb: %14d (Reconstruct at boundaries)\n"
1697 " igrari: %14d (Account for gravity)\n"
1698 " iclsyr: %14d (Symmetry implicitation)\n"
1699 " iclptr: %14d (Wall implicitation)\n"),
1700 cs_glob_turb_ref_values->almax,
1701 cs_glob_turb_ref_values->uref,
1702 cs_glob_turb_rans_model->irijco,
1703 cs_glob_turb_rans_model->irijnu,
1704 cs_glob_turb_rans_model->irijrb,
1705 cs_glob_turb_rans_model->igrari,
1706 cs_glob_turb_rans_model->iclsyr,
1707 cs_glob_turb_rans_model->iclptr);
1708
1709 int idirsm = cs_glob_turb_rans_model->idirsm;
1710 if (idirsm < 0 || idirsm > 1)
1711 idirsm = 2;
1712
1713 const char *s_turb_diff_model[3]
1714 = {"scalar diffusivity: Shir model",
1715 "tensorial diffusivity: Daly-Harlow model",
1716 "<unknown>"};
1717
1718 cs_log_printf(CS_LOG_SETUP,
1719 _(" idirsm: %14d (%s)\n"),
1720 idirsm, s_turb_diff_model[idirsm]);
1721
1722 if (turb_model->iturb == CS_TURB_RIJ_EPSILON_LRR) {
1723 cs_log_printf(CS_LOG_SETUP,
1724 _(" irijec: %14d (Wall echo terms)\n"
1725 " idifre: %14d (Handle diffusion tensor)\n"),
1726 cs_glob_turb_rans_model->irijec,
1727 cs_glob_turb_rans_model->idifre);
1728 }
1729 else if (turb_model->iturb == CS_TURB_RIJ_EPSILON_EBRSM) {
1730 cs_log_printf(CS_LOG_SETUP,
1731 _(" reinit_turb: %14d (turbulence reinitialization)\n"),
1732 cs_glob_turb_rans_model->reinit_turb);
1733 }
1734
1735 }
1736 else if (turb_model->type == CS_TURB_LES) {
1737 cs_log_printf(CS_LOG_SETUP,
1738 _(" csmago: %14.5e (Smagorinsky constant)\n"
1739 " cwale: %14.5e (WALE model constant)\n"
1740 " xlesfl: %14.5e (Filter with in a cell is)\n"
1741 " ales: %14.5e (written as)\n"
1742 " bles: %14.5e (xlesfl*(ales*volume)**(bles))\n"
1743 " idries: %14d (=1 Van Driest damping)\n"
1744 " cdries: %14.5e (Van Driest constant)\n"
1745 " xlesfd: %14.5e (Ratio between the explicit)\n"
1746 " (filter and LES filter)\n"
1747 " (recommended value: 1.5)\n"
1748 " smagmx: %14.5e (Max Smagorinsky in the)\n"
1749 " (dynamic model case)\n"),
1750 cs_turb_csmago, cs_turb_cwale, cs_turb_xlesfl,
1751 cs_turb_ales, cs_turb_bles, cs_glob_turb_les_model->idries,
1752 cs_turb_cdries, cs_turb_xlesfd, cs_turb_csmago_max);
1753
1754 }
1755 else if ( turb_model->iturb == CS_TURB_V2F_PHI
1756 || turb_model->iturb == CS_TURB_V2F_BL_V2K) {
1757
1758 cs_log_printf(CS_LOG_SETUP,
1759 _(" almax: %14.5e (Characteristic length)\n"
1760 " uref: %14.5e (Characteristic velocity)\n"
1761 " iclkep: %14d (k-epsilon clipping model)\n"
1762 " ikecou: %14d (k-epsilon coupling mode)\n"
1763 " igrake: %14d (Account for gravity)\n"),
1764 cs_glob_turb_ref_values->almax,
1765 cs_glob_turb_ref_values->uref,
1766 cs_glob_turb_rans_model->iclkep,
1767 cs_glob_turb_rans_model->ikecou,
1768 cs_glob_turb_rans_model->igrake);
1769
1770 if ( cs_glob_turb_rans_model->ikecou == 0
1771 && cs_glob_time_step_options->idtvar >= 0) {
1772
1773 cs_real_t relaxvk, relaxve;
1774 cs_field_get_key_struct(CS_F_(k), key_cal_opt_id, &var_cal_opt);
1775 relaxvk = var_cal_opt.relaxv;
1776 cs_field_get_key_struct(CS_F_(eps), key_cal_opt_id, &var_cal_opt);
1777 relaxve = var_cal_opt.relaxv;
1778 cs_log_printf
1779 (CS_LOG_SETUP,
1780 _(" relaxv: %14.5e for k (Relaxation)\n"
1781 " relaxv: %14.5e for epsilon (Relaxation)\n"),
1782 relaxvk, relaxve);
1783
1784 }
1785 else
1786 cs_log_printf(CS_LOG_SETUP,_("\n"));
1787
1788 }
1789 else if (turb_model->iturb == CS_TURB_K_OMEGA) {
1790
1791 const char *hybrid_turb_value_str[]
1792 = {N_("CS_HYBRID_NONE (no RANS-LES hybrid model)"),
1793 N_("CS_HYBRID_DES (RANS-LES hybrid model)"),
1794 N_("CS_HYBRID_DDES (RANS-LES hybrid model)"),
1795 N_("CS_HYBRID_SAS (Scale Adpative Model)")};
1796
1797 cs_log_printf(CS_LOG_SETUP,
1798 _(" almax: %14.5e (Characteristic length)\n"
1799 " uref: %14.5e (Characteristic velocity)\n"
1800 " ikecou: %14d (k-epsilon coupling mode)\n"
1801 " hybrid_turb: %s\n"
1802 " igrake: %14d (Account for gravity)\n"),
1803 cs_glob_turb_ref_values->almax,
1804 cs_glob_turb_ref_values->uref,
1805 cs_glob_turb_rans_model->ikecou,
1806 hybrid_turb_value_str[turb_model->hybrid_turb],
1807 cs_glob_turb_rans_model->igrake);
1808
1809 if ( cs_glob_turb_rans_model->ikecou == 0
1810 && cs_glob_time_step_options->idtvar >= 0) {
1811
1812 cs_field_get_key_struct(CS_F_(k), key_cal_opt_id, &var_cal_opt);
1813 cs_real_t relaxvk = var_cal_opt.relaxv;
1814 cs_field_get_key_struct(CS_F_(omg), key_cal_opt_id, &var_cal_opt);
1815 cs_real_t relaxvo = var_cal_opt.relaxv;
1816 cs_log_printf(CS_LOG_SETUP,
1817 _(" relaxv: %14.5e for k (Relaxation)\n"
1818 " relaxv: %14.5e for omega (Relaxation)\n"),
1819 relaxvk, relaxvo);
1820 }
1821 else
1822 cs_log_printf(CS_LOG_SETUP,_("\n"));
1823
1824 }
1825 else if (turb_model->iturb == CS_TURB_SPALART_ALLMARAS) {
1826
1827 cs_field_get_key_struct(CS_F_(nusa), key_cal_opt_id, &var_cal_opt);
1828 cs_log_printf(CS_LOG_SETUP,
1829 _(" almax: %14.5e (Characteristic length)\n"
1830 " uref: %14.5e (Characteristic velocity)\n"
1831 " relaxv: %14.5e for nu (Relaxation)\n"),
1832 cs_glob_turb_ref_values->almax,
1833 cs_glob_turb_ref_values->uref,
1834 var_cal_opt.relaxv);
1835
1836 }
1837
1838 if ( turb_model->type == CS_TURB_RANS
1839 && turb_model->order == CS_TURB_FIRST_ORDER){
1840
1841 const char *irccor_value_str[]
1842 = {N_("0 (no rotation/curvature correction)"),
1843 N_("1 (rotation/curvature correction)")};
1844 cs_log_printf(CS_LOG_SETUP,
1845 _(" irccor: %s\n"),
1846 irccor_value_str[cs_glob_turb_rans_model->irccor]);
1847 }
1848 }
1849
1850 /*----------------------------------------------------------------------------*/
1851 /*!
1852 * \brief Print the turbulent constants to setup.log.
1853 *
1854 *----------------------------------------------------------------------------*/
1855
1856 void
cs_turb_constants_log_setup(void)1857 cs_turb_constants_log_setup(void)
1858 {
1859 if (cs_glob_turb_model == NULL)
1860 return;
1861
1862 const cs_turb_model_t *turb_model = cs_glob_turb_model;
1863
1864 cs_log_printf(CS_LOG_SETUP,
1865 _("\nConstants:\n\n"
1866 " xkappa: %14.5e (Von Karman constant)\n"
1867 " cstlog: %14.5e (U+=Log(y+)/kappa +cstlog)\n"
1868 " apow: %14.5e (U+=apow (y+)**bpow (W&W law))\n"
1869 " bpow: %14.5e (U+=apow (y+)**bpow (W&W law))\n\n"),
1870 cs_turb_xkappa, cs_turb_cstlog, cs_turb_apow, cs_turb_bpow);
1871
1872 if (turb_model->iturb != CS_TURB_NONE)
1873 cs_log_printf(CS_LOG_SETUP,
1874 _(" %s constants:\n"),
1875 _turbulence_model_name(turb_model->iturb));
1876
1877 if ( turb_model->iturb == CS_TURB_K_EPSILON
1878 || turb_model->iturb == CS_TURB_K_EPSILON_LIN_PROD
1879 || turb_model->iturb == CS_TURB_K_EPSILON_LS
1880 || turb_model->iturb == CS_TURB_K_EPSILON_QUAD)
1881 cs_log_printf
1882 (CS_LOG_SETUP,
1883 _(" ce1: %14.5e (Cepsilon 1: production coef.)\n"
1884 " ce2: %14.5e (Cepsilon 2: dissipat. coef.)\n"
1885 " cmu: %14.5e (Cmu constant)\n"),
1886 cs_turb_ce1, cs_turb_ce2, cs_turb_cmu);
1887
1888 else if (turb_model->iturb == CS_TURB_RIJ_EPSILON_LRR)
1889 cs_log_printf
1890 (CS_LOG_SETUP,
1891 _(" ce1: %14.5e (Cepsilon 1: production coef.)\n"
1892 " ce2: %14.5e (Cepsilon 2: dissipat. coef.)\n"
1893 " crij1: %14.5e (Slow term coefficient)\n"
1894 " crij2: %14.5e (Fast term coefficient)\n"
1895 " crij3: %14.5e (Gravity term coefficient)\n"
1896 " csrij: %14.5e (Rij diffusion coeff.)\n"
1897 " crijp1: %14.5e (Slow coeff. for wall echo)\n"
1898 " crijp2: %14.5e (Fast coeff. for wall echo)\n"
1899 " cmu: %14.5e (Cmu constant)\n"),
1900 cs_turb_ce1, cs_turb_ce2, cs_turb_crij1, cs_turb_crij2,
1901 cs_turb_crij3, cs_turb_csrij, cs_turb_crijp1,
1902 cs_turb_crijp2, cs_turb_cmu);
1903
1904 else if (turb_model->iturb == CS_TURB_RIJ_EPSILON_SSG)
1905 cs_log_printf
1906 (CS_LOG_SETUP,
1907 _(" cssgs1: %14.5e (Cs1 coeff.)\n"
1908 " cssgs2: %14.5e (Cs2 coeff.)\n"
1909 " cssgr1: %14.5e (Cr1 coeff.)\n"
1910 " cssgr2: %14.5e (Cr2 coeff.)\n"
1911 " cssgr3: %14.5e (Cr3 coeff.)\n"
1912 " cssgr4: %14.5e (Cr4 coeff.)\n"
1913 " cssgr5: %14.5e (Cr5 coeff.)\n"
1914 " csrij: %14.5e (Rij Cs diffusion coeff.)\n"
1915 " crij3: %14.5e (Gravity term coeff.)\n"
1916 " ce1: %14.5e (Ceps1 coeff.)\n"
1917 " cssge2: %14.5e (Ceps2 coeff.)\n"
1918 " cmu: %14.5e (Cmu constant)\n"),
1919 cs_turb_cssgs1, cs_turb_cssgs2, cs_turb_cssgr1,
1920 cs_turb_cssgr2, cs_turb_cssgr3, cs_turb_cssgr4,
1921 cs_turb_cssgr5, cs_turb_csrij, cs_turb_crij3,
1922 cs_turb_ce1, cs_turb_cssge2,
1923 cs_turb_cmu);
1924
1925 else if (turb_model->iturb == CS_TURB_RIJ_EPSILON_EBRSM) {
1926 cs_turb_crij3 = 0.6;
1927 cs_log_printf
1928 (CS_LOG_SETUP,
1929 _(" cebms1: %14.5e (Cs1 coeff.)\n"
1930 " cebmr1: %14.5e (Cr1 coeff.)\n"
1931 " cebmr2: %14.5e (Cr2 coeff.)\n"
1932 " cebmr3: %14.5e (Cr3 coeff.)\n"
1933 " cebmr4: %14.5e (Cr4 coeff.)\n"
1934 " cebmr5: %14.5e (Cr5 coeff.)\n"
1935 " csrij: %14.5e (Rij Cs diffusion coeff.)\n"
1936 " crij3: %14.5e (Gravity term coeff.)\n"
1937 " cebme2: %14.5e (Coef Ceps2)\n"
1938 " ce1: %14.5e (Coef Ceps1)\n"
1939 " xa1: %14.5e (Coef A1)\n"
1940 " xceta: %14.5e (Coef Ceta)\n"
1941 " xct: %14.5e (Coef CT)\n"),
1942 cs_turb_cebms1, cs_turb_cebmr1, cs_turb_cebmr2,
1943 cs_turb_cebmr3, cs_turb_cebmr4, cs_turb_cebmr5,
1944 cs_turb_csrij, cs_turb_crij3, cs_turb_cebme2,
1945 cs_turb_ce1, cs_turb_xa1,
1946 cs_turb_xceta, cs_turb_xct);
1947
1948 }
1949
1950 else if (turb_model->iturb == CS_TURB_V2F_PHI)
1951 cs_log_printf
1952 (CS_LOG_SETUP,
1953 _(" cv2fa1: %14.5e (a1 to calculate Cepsilon1)\n"
1954 " cv2fe2: %14.5e (Cepsilon 2: dissip. coeff.)\n"
1955 " cmu : %14.5e (Cmu constant)\n"
1956 " cv2fct: %14.5e (CT constant)\n"
1957 " cv2fcl: %14.5e (CL constant)\n"
1958 " cv2fet: %14.5e (C_eta constant)\n"
1959 " cv2fc1: %14.5e (C1 constant)\n"
1960 " cv2fc2: %14.5e (C2 constant)\n"),
1961 cs_turb_cv2fa1, cs_turb_cv2fe2,
1962 cs_turb_cmu, cs_turb_cv2fct,
1963 cs_turb_cv2fcl, cs_turb_cv2fet, cs_turb_cv2fc1,
1964 cs_turb_cv2fc2);
1965
1966 else if (turb_model->iturb == CS_TURB_V2F_BL_V2K)
1967 cs_log_printf
1968 (CS_LOG_SETUP,
1969 _(" cpale1: %14.5e (Cepsilon 1 : Prod. coeff.)\n"
1970 " cpale2: %14.5e (Cepsilon 2 : Diss. coeff.)\n"
1971 " cpale3: %14.5e (Cepsilon 3 : E term coeff.)\n"
1972 " cpale4: %14.5e (Cepsilon 4 : Mod Diss. coef.)\n"
1973 " cmu : %14.5e (Cmu constant)\n"
1974 " cpalct: %14.5e (CT constant)\n"
1975 " cpalcl: %14.5e (CL constant)\n"
1976 " cpalet: %14.5e (C_eta constant)\n"
1977 " cpalc1: %14.5e (C1 constant)\n"
1978 " cpalc2: %14.5e (C2 constant)\n"),
1979 cs_turb_cpale1, cs_turb_cpale2, cs_turb_cpale3,
1980 cs_turb_cpale4,
1981 cs_turb_cmu, cs_turb_cpalct, cs_turb_cpalcl,
1982 cs_turb_cpalet, cs_turb_cpalc1, cs_turb_cpalc2);
1983
1984 else if (turb_model->iturb == CS_TURB_K_OMEGA)
1985 cs_log_printf
1986 (CS_LOG_SETUP,
1987 _(" ckwsk1: %14.5e (sigma_k1 constant)\n"
1988 " ckwsk2: %14.5e (sigma_k2 constant)\n"
1989 " ckwsw1: %14.5e (sigma_omega1 constant)\n"
1990 " ckwsw2: %14.5e (sigma_omega2 constant)\n"
1991 " ckwbt1: %14.5e (beta1 constant)\n"
1992 " ckwbt2: %14.5e (beta2 constant)\n"
1993 " ckwgm1: %14.5e (gamma1 constant)\n"
1994 " ckwgm2: %14.5e (gamma2 constant)\n"
1995 " ckwa1: %14.5e (a1 constant to compute mu_t)\n"
1996 " ckwc1: %14.5e (c1 const. for prod. limiter)\n"
1997 " cmu: %14.5e (Cmu (or Beta*) constant for)\n"
1998 " omega/epsilon conversion)\n"),
1999 cs_turb_ckwsk1, cs_turb_ckwsk2, cs_turb_ckwsw1,
2000 cs_turb_ckwsw2, cs_turb_ckwbt1, cs_turb_ckwbt2,
2001 cs_turb_ckwgm1, cs_turb_ckwgm2, cs_turb_ckwa1,
2002 cs_turb_ckwc1, cs_turb_cmu);
2003
2004 else if (turb_model->iturb == CS_TURB_SPALART_ALLMARAS)
2005 cs_log_printf(CS_LOG_SETUP,
2006 _(" csab1: %14.5e (b1 constant)\n"
2007 " csab2: %14.5e (b2 constant)\n"
2008 " csasig: %14.5e (sigma constant)\n"
2009 " csav1: %14.5e (v1 constant)\n"
2010 " csaw1: %14.5e (w1 constant)\n"
2011 " csaw2: %14.5e (w2 constant)\n"
2012 " csaw3: %14.5e (w3 constant)\n"),
2013 cs_turb_csab1, cs_turb_csab2, cs_turb_csasig,
2014 cs_turb_csav1, cs_turb_csaw1, cs_turb_csaw2,
2015 cs_turb_csaw3);
2016
2017 int iokss = 0, iokcaz = 0;
2018
2019 if (cs_glob_turb_rans_model != NULL) {
2020 if (cs_glob_turb_rans_model->irccor == 1) {
2021 if (cs_glob_turb_rans_model->itycor == 1)
2022 iokcaz = 1;
2023 else if (cs_glob_turb_rans_model->itycor == 2)
2024 iokss = 1;
2025 }
2026 }
2027
2028 if (iokcaz > 0)
2029 cs_log_printf(CS_LOG_SETUP,
2030 _(" Rotation/curvature correction (Cazalbou)\n"
2031 " ccaze2: %14.5e (Coef Ce2^0)\n"
2032 " ccazsc: %14.5e (Coef Csc)\n"
2033 " ccaza: %14.5e (Coef a)\n"
2034 " ccazb: %14.5e (Coef b)\n"
2035 " ccazc: %14.5e (Coef c)\n"
2036 " ccazd: %14.5e (Coef d)\n"),
2037 cs_turb_ccaze2, cs_turb_ccazsc, cs_turb_ccaza,
2038 cs_turb_ccazb, cs_turb_ccazc, cs_turb_ccazd);
2039
2040 if (iokss > 0)
2041 cs_log_printf(CS_LOG_SETUP,
2042 _(" Rotation/curvature correction (Spalart-Shur)\n"
2043 " cssr1: %14.5e (Coef c_r1)\n"
2044 " cssr2: %14.5e (Coef c_r2)\n"
2045 " cssr3: %14.5e (Coef c_r3)\n"),
2046 cs_turb_cssr1, cs_turb_cssr2, cs_turb_cssr3);
2047 }
2048
2049 /*----------------------------------------------------------------------------*/
2050 /*!
2051 * \brief Clipping for the turbulence flux vector
2052 *
2053 * \param[in] flux_id turbulent flux index
2054 * \param[in] variance_id scalar variance index
2055 *
2056 *----------------------------------------------------------------------------*/
2057
2058 void
cs_clip_turbulent_fluxes(int flux_id,int variance_id)2059 cs_clip_turbulent_fluxes(int flux_id,
2060 int variance_id)
2061 {
2062 cs_lnum_t n_cells = cs_glob_mesh->n_cells;
2063
2064 const cs_real_6_t *cvar_rij = (const cs_real_6_t *)CS_F_(rij)->val;
2065 const cs_real_t *cvar_tt
2066 = (const cs_real_t *)cs_field_by_id(variance_id)->val;
2067 cs_real_3_t *cvar_rit = (cs_real_3_t *) cs_field_by_id(flux_id)->val;
2068 cs_real_3_t *cvar_clip_rit = NULL;
2069
2070 cs_field_t *field_rit = cs_field_by_id(flux_id);
2071
2072 /* Local variables */
2073 const cs_real_t tol_jacobi = 1.0e-12;
2074 const cs_real_t l_threshold = 1.0e12;
2075 cs_lnum_t iclip = 0;
2076 cs_real_t flsq, maj;
2077 cs_real_33_t rij;
2078 cs_real_33_t eigvect = {{1.0,0.0,0.0},{0.0,1.0,0.0},{0.0,0.0,1.0}};
2079 cs_real_3_t eigval;
2080 cs_real_3_t rot_rit;
2081 cs_real_3_t rit;
2082
2083 /* Get clippings field for DFM */
2084 cs_lnum_t kclipp, clip_rit_id;
2085 kclipp = cs_field_key_id("clipping_id");
2086 clip_rit_id = cs_field_get_key_int(field_rit, kclipp);
2087 if (clip_rit_id >= 0) {
2088 cvar_clip_rit = (cs_real_3_t *) cs_field_by_id(clip_rit_id)->val;
2089 for (cs_lnum_t cell_id = 0; cell_id < n_cells; cell_id++) {
2090 cvar_clip_rit[cell_id][0] = 0;
2091 cvar_clip_rit[cell_id][1] = 0;
2092 cvar_clip_rit[cell_id][2] = 0;
2093 }
2094 }
2095
2096 cs_real_t rit_min_prcoord[3] = {l_threshold,l_threshold,l_threshold};
2097 cs_real_t rit_max_prcoord[3] = {-l_threshold,-l_threshold,-l_threshold};
2098
2099 cs_lnum_t iclip_tab_max[3] = {0,0,0};
2100 cs_lnum_t iclip_tab_min[3] = {0,0,0};
2101
2102 for (int cell_id = 0; cell_id < n_cells; cell_id++) {
2103 rij[0][0] = cvar_rij[cell_id][0];
2104 rij[1][1] = cvar_rij[cell_id][1];
2105 rij[2][2] = cvar_rij[cell_id][2];
2106 rij[0][1] = cvar_rij[cell_id][3];
2107 rij[1][0] = cvar_rij[cell_id][3];
2108 rij[0][2] = cvar_rij[cell_id][5];
2109 rij[2][0] = cvar_rij[cell_id][5];
2110 rij[1][2] = cvar_rij[cell_id][4];
2111 rij[2][1] = cvar_rij[cell_id][4];
2112
2113 rit[0] = cvar_rit[cell_id][0];
2114 rit[1] = cvar_rit[cell_id][1];
2115 rit[2] = cvar_rit[cell_id][2];
2116
2117 cs_math_33_eig_val_vec(rij,tol_jacobi,eigval,eigvect);
2118 cs_math_33_3_product(eigvect,cvar_rit[cell_id],rot_rit);
2119
2120 for (int i = 0; i < 3; i++) {
2121 rit_min_prcoord[i] = CS_MIN(rit_min_prcoord[i],rot_rit[i]);
2122 rit_max_prcoord[i] = CS_MAX(rit_max_prcoord[i],rot_rit[i]);
2123 }
2124
2125 for (int i = 0; i < 3; i++) {
2126 flsq = pow(rot_rit[i],2);
2127 maj = eigval[i]*cvar_tt[cell_id];
2128 if ((flsq > 1.0e-12) && (flsq > maj)) {
2129 rot_rit[i] = rot_rit[i]*sqrt(maj/flsq);
2130 iclip = 1;
2131 if (rot_rit[i] > 0)
2132 iclip_tab_max[i] += 1;
2133 else
2134 iclip_tab_min[i] += 1;
2135 }
2136 }
2137
2138 if (iclip > 0) {
2139 cs_math_33t_3_product(eigvect,rot_rit,cvar_rit[cell_id]);
2140 cvar_clip_rit[cell_id][0] = cvar_rit[cell_id][0] - rit[0];
2141 cvar_clip_rit[cell_id][1] = cvar_rit[cell_id][1] - rit[1];
2142 cvar_clip_rit[cell_id][2] = cvar_rit[cell_id][2] - rit[2];
2143 }
2144 }
2145
2146 cs_lnum_t iclip_max = iclip_tab_max[0] + iclip_tab_max[1]
2147 + iclip_tab_max[2];
2148 cs_lnum_t iclip_min = iclip_tab_min[0] + iclip_tab_min[1]
2149 + iclip_tab_min[2];
2150
2151 /* Save clippings for log */
2152 cs_log_iteration_clipping_field(flux_id,
2153 iclip_min,
2154 iclip_max,
2155 rit_min_prcoord,
2156 rit_max_prcoord,
2157 iclip_tab_min,
2158 iclip_tab_max);
2159 return;
2160 }
2161
2162 /*----------------------------------------------------------------------------*/
2163
2164 END_C_DECLS
2165