programmer's documentation
cs_wall_functions.h
Go to the documentation of this file.
1 #ifndef __CS_WALL_FUNCTIONS_H__
2 #define __CS_WALL_FUNCTIONS_H__
3 
4 /*============================================================================
5  * Wall functions
6  *============================================================================*/
7 
8 /*
9  This file is part of Code_Saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2015 EDF S.A.
12 
13  This program is free software; you can redistribute it and/or modify it under
14  the terms of the GNU General Public License as published by the Free Software
15  Foundation; either version 2 of the License, or (at your option) any later
16  version.
17 
18  This program is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20  FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
21  details.
22 
23  You should have received a copy of the GNU General Public License along with
24  this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
25  Street, Fifth Floor, Boston, MA 02110-1301, USA.
26 */
27 
28 /*----------------------------------------------------------------------------*/
29 
30 /*----------------------------------------------------------------------------
31  * BFT library headers
32  *----------------------------------------------------------------------------*/
33 
34 #include <bft_printf.h>
35 
36 /*----------------------------------------------------------------------------
37  * Local headers
38  *----------------------------------------------------------------------------*/
39 
40 #include "cs_base.h"
41 #include "cs_turbulence_model.h"
42 
43 /*----------------------------------------------------------------------------*/
44 
46 
47 /*=============================================================================
48  * Local Macro definitions
49  *============================================================================*/
50 
51 /*============================================================================
52  * Type definition
53  *============================================================================*/
54 
55 /* Wall function type */
56 /*--------------------*/
57 
58 typedef enum {
59 
66 
68 
69 typedef enum {
70 
73 
75 
76 /* Wall functions descriptor */
77 /*---------------------------*/
78 
79 typedef struct {
80 
81  cs_wall_f_type_t iwallf; /* wall function type */
82 
83  cs_wall_f_s_type_t iwalfs; /* wall function type for scalars */
84 
85  int iwallt; /* exchange coefficient correlation
86  - 0: not used by default
87  - 1: exchange coefficient computed with a
88  correlation */
89 
90  double ypluli; /* limit value of y+ for the viscous
91  sublayer */
92 
94 
95 /*============================================================================
96  * Global variables
97  *============================================================================*/
98 
99 /* Pointer to wall functions descriptor structure */
100 
102 
103 /*============================================================================
104  * Private function definitions
105  *============================================================================*/
106 
107 /*----------------------------------------------------------------------------*/
124 /*----------------------------------------------------------------------------*/
125 
126 inline static void
128  cs_real_t vel,
129  cs_real_t y,
130  int *iuntur,
131  cs_lnum_t *nsubla,
132  cs_lnum_t *nlogla,
133  cs_real_t *ustar,
134  cs_real_t *uk,
135  cs_real_t *yplus,
136  cs_real_t *ypup,
137  cs_real_t *cofimp)
138 {
139  const double ypluli = cs_glob_wall_functions->ypluli;
140 
141  const double ydvisc = y / l_visc;
142 
143  /* Compute the friction velocity ustar */
144 
145  *ustar = pow((vel/(cs_turb_apow * pow(ydvisc, cs_turb_bpow))), cs_turb_dpow);
146  *uk = *ustar;
147  *yplus = *ustar * ydvisc;
148 
149  /* In the viscous sub-layer: U+ = y+ */
150  if (*yplus <= ypluli) {
151 
152  *ustar = sqrt(vel / ydvisc);
153  *yplus = *ustar * ydvisc;
154  *uk = *ustar;
155  *ypup = 1.;
156  *cofimp = 0.;
157 
158  /* Disable the wall funcion count the cell in the viscous sub-layer */
159  *iuntur = 0;
160  *nsubla += 1;
161 
162  /* In the log layer */
163  } else {
164  *ypup = pow(vel, 2. * cs_turb_dpow-1.)
165  / pow(cs_turb_apow, 2. * cs_turb_dpow);
166  *cofimp = 1. + cs_turb_bpow
167  * pow(*ustar, cs_turb_bpow + 1. - 1./cs_turb_dpow)
168  * (pow(2., cs_turb_bpow - 1.) - 2.);
169 
170  /* Count the cell in the log layer */
171  *nlogla += 1;
172 
173  }
174 }
175 
176 /*----------------------------------------------------------------------------*/
195 /*----------------------------------------------------------------------------*/
196 
197 inline static void
199  cs_real_t l_visc,
200  cs_real_t vel,
201  cs_real_t y,
202  int *iuntur,
203  cs_lnum_t *nsubla,
204  cs_lnum_t *nlogla,
205  cs_real_t *ustar,
206  cs_real_t *uk,
207  cs_real_t *yplus,
208  cs_real_t *ypup,
209  cs_real_t *cofimp)
210 {
211  const double ypluli = cs_glob_wall_functions->ypluli;
212 
213  double ustarwer, ustarmin, ustaro, ydvisc;
214  double eps = 0.001;
215  int niter_max = 100;
216  int iter = 0;
217  double reynolds;
218 
219  /* Compute the local Reynolds number */
220 
221  ydvisc = y / l_visc;
222  reynolds = vel * ydvisc;
223 
224  /*
225  * Compute the friction velocity ustar
226  */
227 
228  /* In the viscous sub-layer: U+ = y+ */
229  if (reynolds <= ypluli * ypluli) {
230 
231  *ustar = sqrt(vel / ydvisc);
232  *yplus = *ustar * ydvisc;
233  *uk = *ustar;
234  *ypup = 1.;
235  *cofimp = 0.;
236 
237  /* Disable the wall funcion count the cell in the viscous sub-layer */
238  *iuntur = 0;
239  *nsubla += 1;
240 
241  /* In the log layer */
242  } else {
243 
244  /* The initial value is Wener or the minimun ustar to ensure convergence */
245  ustarwer = pow(fabs(vel) / cs_turb_apow / pow(ydvisc, cs_turb_bpow),
246  cs_turb_dpow);
247  ustarmin = exp(-cs_turb_cstlog * cs_turb_xkappa)/ydvisc;
248  ustaro = CS_MAX(ustarwer, ustarmin);
249  *ustar = (cs_turb_xkappa * vel + ustaro)
250  / (log(ydvisc * ustaro) + cs_turb_xkappa * cs_turb_cstlog + 1.);
251 
252  /* Iterative solving */
253  for (iter = 0; iter < niter_max
254  && fabs(*ustar - ustaro) >= eps * ustaro; iter++) {
255  ustaro = *ustar;
256  *ustar = (cs_turb_xkappa * vel + ustaro)
257  / (log(ydvisc * ustaro) + cs_turb_xkappa * cs_turb_cstlog + 1.);
258  }
259 
260  if (iter >= niter_max) {
261  bft_printf(_("WARNING: non-convergence in the computation\n"
262  "******** of the friction velocity\n\n"
263  "face number: %d \n"
264  "friction vel: %f \n" ), ifac, *ustar);
265  }
266 
267  *uk = *ustar;
268  *yplus = *ustar * ydvisc;
269  *ypup = *yplus / (log(*yplus) / cs_turb_xkappa + cs_turb_cstlog);
270  *cofimp = 1. - *ypup / cs_turb_xkappa * 1.5 / *yplus;
271 
272  /* Count the cell in the log layer */
273  *nlogla += 1;
274 
275  }
276 
277 }
278 
279 /*----------------------------------------------------------------------------*/
299 /*----------------------------------------------------------------------------*/
300 
301 inline static void
303  cs_real_t t_visc,
304  cs_real_t vel,
305  cs_real_t y,
306  cs_real_t kinetic_en,
307  int *iuntur,
308  cs_lnum_t *nsubla,
309  cs_lnum_t *nlogla,
310  cs_real_t *ustar,
311  cs_real_t *uk,
312  cs_real_t *yplus,
313  cs_real_t *ypup,
314  cs_real_t *cofimp)
315 {
316  const double ypluli = cs_glob_wall_functions->ypluli;
317 
318  double rcprod, ml_visc, Re, g;
319 
320  /* Compute the friction velocity ustar */
321 
322  /* Blending for very low values of k */
323  Re = sqrt(kinetic_en) * y / l_visc;
324  g = exp(-Re/11.);
325 
326  *uk = sqrt( (1.-g) * cs_turb_cmu025 * cs_turb_cmu025 * kinetic_en
327  + g * l_visc * vel / y);
328 
329  *yplus = *uk * y / l_visc;
330 
331  /* log layer */
332  if (*yplus > ypluli) {
333 
334  *ustar = vel / (log(*yplus) / cs_turb_xkappa + cs_turb_cstlog);
335  *ypup = *yplus / (log(*yplus) / cs_turb_xkappa + cs_turb_cstlog);
336  /* Mixing length viscosity */
337  ml_visc = cs_turb_xkappa * l_visc * *yplus;
338  rcprod = CS_MIN(cs_turb_xkappa, CS_MAX(1., sqrt(ml_visc / t_visc)) / *yplus);
339  *cofimp = 1. - *ypup / cs_turb_xkappa * ( 2. * rcprod - 1. / (2. * *yplus));
340 
341  *nlogla += 1;
342 
343  /* viscous sub-layer */
344  } else {
345 
346  if (*yplus > 1.e-12) {
347  *ustar = fabs(vel / *yplus); /* FIXME remove that: its is here only to
348  be fully equivalent to the former code. */
349  } else {
350  *ustar = 0.;
351  }
352  *ypup = 1.;
353  *cofimp = 0.;
354 
355  *iuntur = 0;
356  *nsubla += 1;
357 
358  }
359 }
360 
361 /*----------------------------------------------------------------------------*/
382 /*----------------------------------------------------------------------------*/
383 
384 inline static void
386  cs_real_t t_visc,
387  cs_real_t vel,
388  cs_real_t y,
389  cs_real_t kinetic_en,
390  int *iuntur,
391  cs_lnum_t *nsubla,
392  cs_lnum_t *nlogla,
393  cs_real_t *ustar,
394  cs_real_t *uk,
395  cs_real_t *yplus,
396  cs_real_t *dplus,
397  cs_real_t *ypup,
398  cs_real_t *cofimp)
399 {
400  const double ypluli = cs_glob_wall_functions->ypluli;
401 
402  double rcprod, ml_visc, Re, g;
403  /* Compute the friction velocity ustar */
404 
405  /* Blending for very low values of k */
406  Re = sqrt(kinetic_en) * y / l_visc;
407  g = exp(-Re/11.);
408 
409  *uk = sqrt( (1.-g) * cs_turb_cmu025 * cs_turb_cmu025 * kinetic_en
410  + g * l_visc * vel / y);
411 
412  *yplus = *uk * y / l_visc;
413 
414  /* Compute the friction velocity ustar */
415  *uk = cs_turb_cmu025 * sqrt(kinetic_en);
416  *yplus = *uk * y / l_visc;
417 
418  /* Log layer */
419  if (*yplus > ypluli) {
420 
421  *dplus = 0.;
422 
423  *nlogla += 1;
424 
425  /* Viscous sub-layer and therefore shift */
426  } else {
427 
428  *dplus = ypluli - *yplus;
429  *yplus = ypluli;
430 
431  /* Count the cell as if it was in the viscous sub-layer */
432  *nsubla += 1;
433 
434  }
435 
436  /* Mixing length viscosity */
437  ml_visc = cs_turb_xkappa * l_visc * *yplus;
438  rcprod = CS_MIN(cs_turb_xkappa, CS_MAX(1., sqrt(ml_visc / t_visc)) / *yplus);
439 
440  *ustar = vel / (log(*yplus) / cs_turb_xkappa + cs_turb_cstlog);
441  *ypup = (*yplus - *dplus) / (log(*yplus) / cs_turb_xkappa + cs_turb_cstlog);
442  *cofimp = 1. - *ypup
443  / cs_turb_xkappa * (2. * rcprod - 1. / (2. * *yplus - *dplus));
444 }
445 
446 /*----------------------------------------------------------------------------
447  * Compute u+ for a given yk+ between 0.1 and 200 according to the two
448  * scales wall functions using Van Driest mixing length.
449  * This function holds the coefficients of the polynome fitting log(u+).
450  *
451  * parameters:
452  * yplus <-- dimensionless distance
453  *
454  * returns:
455  * the resulting dimensionless velocity.
456  *----------------------------------------------------------------------------*/
457 
458 inline static cs_real_t
460 {
461  /* Coefficients of the polynome fitting log(u+) for yk < 200 */
462  static double aa[11] = {-0.0091921, 3.9577, 0.031578,
463  -0.51013, -2.3254, -0.72665,
464  2.969, 0.48506, -1.5944,
465  0.087309, 0.1987 };
466 
467  cs_real_t y1,y2,y3,y4,y5,y6,y7,y8,y9,y10, uplus;
468 
469  y1 = 0.25 * log(yplus);
470  y2 = y1 * y1;
471  y3 = y2 * y1;
472  y4 = y3 * y1;
473  y5 = y4 * y1;
474  y6 = y5 * y1;
475  y7 = y6 * y1;
476  y8 = y7 * y1;
477  y9 = y8 * y1;
478  y10 = y9 * y1;
479 
480  uplus = aa[0]
481  + aa[1] * y1
482  + aa[2] * y2
483  + aa[3] * y3
484  + aa[4] * y4
485  + aa[5] * y5
486  + aa[6] * y6
487  + aa[7] * y7
488  + aa[8] * y8
489  + aa[9] * y9
490  + aa[10] * y10;
491 
492  return exp(uplus);
493 }
494 
495 /*----------------------------------------------------------------------------*/
530 /*----------------------------------------------------------------------------*/
531 
532 inline static void
534  cs_real_t l_visc,
535  cs_real_t vel,
536  cs_real_t y,
537  cs_real_t kinetic_en,
538  int *iuntur,
539  cs_lnum_t *nsubla,
540  cs_lnum_t *nlogla,
541  cs_real_t *ustar,
542  cs_real_t *uk,
543  cs_real_t *yplus,
544  cs_real_t *ypup,
545  cs_real_t *cofimp,
546  cs_real_t *lmk,
547  cs_real_t kr,
548  bool wf)
549 {
550  double urplus, dup, lmk15;
551 
552  if (wf)
553  *uk = sqrt(sqrt((1.-cs_turb_crij2)/cs_turb_crij1 * rnnb * kinetic_en));
554 
555  /* Set a low threshold value in case tangential velocity is zero */
556  *yplus = CS_MAX(*uk * y / l_visc, 1.e-4);
557 
558  /* Dimensionless roughness */
559  cs_real_t krp = *uk * kr / l_visc;
560 
561  /* Extension of Van Driest mixing length according to Rotta (1962) with
562  Cebeci & Chang (1978) correlation */
563  cs_real_t dyrp = 0.9 * (sqrt(krp) - krp * exp(-krp / 6.));
564  cs_real_t yrplus = *yplus + dyrp;
565 
566  if (dyrp <= 1.e-1)
567  dup = dyrp;
568  else if (dyrp <= 200.)
569  dup = _vdriest_dupdyp_integral(dyrp);
570  else
571  dup = 16.088739022054590 + log(dyrp/200.) / cs_turb_xkappa;
572 
573  if (yrplus <= 1.e-1) {
574 
575  urplus = yrplus;
576 
577  if (wf) {
578  *iuntur = 0;
579  *nsubla += 1;
580 
581  *lmk = 0.;
582 
583  *ypup = 1.;
584 
585  *cofimp = 0.;
586  }
587 
588  } else if (yrplus <= 200.) {
589 
590  urplus = _vdriest_dupdyp_integral(yrplus);
591 
592  if (wf) {
593  *nlogla += 1;
594 
595  *ypup = *yplus / (urplus-dup);
596 
597  /* Mixing length in y+ */
598  *lmk = cs_turb_xkappa * (*yplus) *(1-exp(- (*yplus) / cs_turb_vdriest));
599 
600  /* Mixing length in 3/2*y+ */
601  lmk15 = cs_turb_xkappa * 1.5 * (*yplus) *(1-exp(- 1.5 * (*yplus)
602  / cs_turb_vdriest));
603 
604  *cofimp = 1. - (2. / (1. + *lmk) - 1. / (1. + lmk15)) * *ypup;
605  }
606 
607  } else {
608 
609  urplus = 16.088739022054590 + log(yrplus/200) / cs_turb_xkappa;
610 
611  if (wf) {
612  *nlogla += 1;
613 
614  *ypup = *yplus / (urplus-dup);
615 
616  /* Mixing length in y+ */
617  *lmk = cs_turb_xkappa * (*yplus) *(1-exp(- (*yplus) / cs_turb_vdriest));
618 
619  /* Mixing length in 3/2*y+ */
620  lmk15 = cs_turb_xkappa * 1.5 * (*yplus) *(1-exp(- 1.5 * (*yplus)
621  / cs_turb_vdriest));
622 
623  *cofimp = 1. - (2. / *lmk - 1. / lmk15) * *ypup;
624  }
625 
626  }
627 
628  *ustar = vel / (urplus-dup);
629 }
630 
631 /*----------------------------------------------------------------------------*/
651 /*----------------------------------------------------------------------------*/
652 
653 inline static void
655  cs_real_t t_visc,
656  cs_real_t vel,
657  cs_real_t y,
658  int *iuntur,
659  cs_lnum_t *nsubla,
660  cs_lnum_t *nlogla,
661  cs_real_t *ustar,
662  cs_real_t *uk,
663  cs_real_t *yplus,
664  cs_real_t *dplus,
665  cs_real_t *ypup,
666  cs_real_t *cofimp)
667 {
668  const double ypluli = cs_glob_wall_functions->ypluli;
669 
670  /* Compute the friction velocity ustar */
671 
672  *ustar = sqrt(vel * l_visc / y);
673  *yplus = *ustar * y / l_visc;
674  *uk = *ustar;
675  *ypup = l_visc / (l_visc + t_visc);
676  *cofimp = 0.;
677  *iuntur = 0;
678 
679  if (*yplus <= ypluli) {
680 
681  /* Disable the wall funcion count the cell in the viscous sub-layer */
682  *nsubla += 1;
683 
684  } else {
685 
686  /* Count the cell as if it was in the viscous sub-layer */
687  *nsubla += 1;
688 
689  }
690 }
691 
692 /*----------------------------------------------------------------------------*/
718 /*----------------------------------------------------------------------------*/
719 
720 inline static void
722  cs_real_t prt,
723  cs_real_t yplus,
724  cs_real_t dplus,
725  cs_real_t *htur,
726  cs_real_t *yplim)
727 {
728  /* Local variables */
729  double tplus;
730  double beta2,a2;
731  double yp2;
732  double prlm1;
733 
734  const double epzero = 1.e-12;
735 
736  /*==========================================================================*/
737 
738  /*==========================================================================
739  1. Initializations
740  ==========================================================================*/
741 
742  /*==========================================================================*/
743 
744  (*htur) = CS_MAX(yplus-dplus,epzero)/CS_MAX(yplus,epzero);
745 
746  prlm1 = 0.1;
747 
748  /*==========================================================================
749  2. Compute htur for small Prandtl numbers
750  ==========================================================================*/
751 
752  if (prl <= prlm1) {
753  (*yplim) = prt/(prl*cs_turb_xkappa);
754  if (yplus > (*yplim)) {
755  tplus = prl*(*yplim) + prt/cs_turb_xkappa * log(yplus/(*yplim));
756  (*htur) = prl*(yplus-dplus)/tplus;
757  }
758 
759  /*========================================================================
760  3. Compute htur for the model with three sub-layers
761  ========================================================================*/
762 
763  } else {
764  yp2 = cs_turb_xkappa*1000./prt;
765  yp2 = sqrt(yp2);
766  (*yplim) = pow(1000./prl,1./3.);
767 
768  a2 = 15.*pow(prl,2./3.);
769  beta2 = a2 - 500./ pow(yp2,2);
770 
771  if (yplus >= (*yplim) && yplus < yp2) {
772  tplus = a2 - 500./(yplus*yplus);
773  (*htur) = prl*(yplus-dplus)/tplus;
774  }
775 
776  if (yplus >= yp2) {
777  tplus = beta2 + prt/cs_turb_xkappa*log(yplus/yp2);
778  (*htur) = prl*(yplus-dplus)/tplus;
779  }
780 
781  }
782 }
783 
784 /*----------------------------------------------------------------------------*/
807 /*----------------------------------------------------------------------------*/
808 
809 inline static void
811  cs_real_t prt,
812  cs_real_t yplus,
813  cs_real_t *htur)
814 {
815  cs_real_t prlrat = prl / prt;
816 
817  /* Parameters of the numerical quadrature */
818  const int ninter_max = 100;
819  const cs_real_t ypmax = 1.e2;
820 
821  /* No correction for very small yplus */
822  if (yplus <= 0.1)
823  *htur = 1.;
824  else {
825  cs_real_t ypint = CS_MIN(yplus, ypmax);
826 
827  /* The number of sub-intervals is taken proportional to yplus and equal to
828  * ninter_max if yplus=ypmax */
829 
830  int npeff = CS_MAX((int)(ypint / ypmax * (double)(ninter_max)), 1);
831 
832  double dy = ypint / (double)(npeff);
833  cs_real_t stplus = 0.;
834  cs_real_t nut1 = 0.;
835  cs_real_t nut2 = 0.;
836 
837  for (int ip = 1; ip <= npeff; ip++) {
838  double yp = ypint * (double)(ip) / (double)(npeff);
839  nut2 = cs_turb_xkappa * yp * (1. - exp(-yp / cs_turb_vdriest));
840  stplus += dy / (1. + prlrat * 0.5 * (nut1 + nut2));
841  nut1 = nut2;
842  }
843 
844  if (yplus > ypint) {
845  cs_real_t r = prlrat * cs_turb_xkappa;
846  stplus += log( (1. + r*yplus) / (1. + r*ypint)) / r;
847  }
848 
849  if (stplus >= 1.e-6)
850  *htur = yplus / stplus;
851  else
852  *htur = 1.;
853  }
854 }
855 
856 /*============================================================================
857  * Public function definitions for Fortran API
858  *============================================================================*/
859 
860 /*----------------------------------------------------------------------------
861  * Wrapper to cs_wall_functions_velocity.
862  *----------------------------------------------------------------------------*/
863 
864 void CS_PROCF (wallfunctions, WALLFUNCTIONS)
865 (
866  const cs_int_t *const iwallf,
867  const cs_lnum_t *const ifac,
868  const cs_real_t *const viscosity,
869  const cs_real_t *const t_visc,
870  const cs_real_t *const vel,
871  const cs_real_t *const y,
872  const cs_real_t *const rnnb,
873  const cs_real_t *const kinetic_en,
874  cs_int_t *iuntur,
875  cs_lnum_t *nsubla,
876  cs_lnum_t *nlogla,
877  cs_real_t *ustar,
878  cs_real_t *uk,
879  cs_real_t *yplus,
880  cs_real_t *ypup,
881  cs_real_t *cofimp,
882  cs_real_t *dplus
883 );
884 
885 /*----------------------------------------------------------------------------
886  * Wrapper to cs_wall_functions_scalar.
887  *----------------------------------------------------------------------------*/
888 
889 void CS_PROCF (hturbp, HTURBP)
890 (
891  const cs_int_t *const iwalfs,
892  const cs_real_t *const prl,
893  const cs_real_t *const prt,
894  const cs_real_t *const yplus,
895  const cs_real_t *const dplus,
896  cs_real_t *htur,
897  cs_real_t *yplim
898 );
899 
900 /*=============================================================================
901  * Public function prototypes
902  *============================================================================*/
903 
907 /*-------------------------------------------------------------------------------
908  Arguments
909  ______________________________________________________________________________.
910  mode name role !
911  ______________________________________________________________________________*/
936 /*-------------------------------------------------------------------------------*/
937 
938 void
940  cs_lnum_t ifac,
941  cs_real_t l_visc,
942  cs_real_t t_visc,
943  cs_real_t vel,
944  cs_real_t y,
945  cs_real_t rnnb,
946  cs_real_t kinetic_en,
947  int *iuntur,
948  cs_lnum_t *nsubla,
949  cs_lnum_t *nlogla,
950  cs_real_t *ustar,
951  cs_real_t *uk,
952  cs_real_t *yplus,
953  cs_real_t *ypup,
954  cs_real_t *cofimp,
955  cs_real_t *dplus);
956 
957 /*-------------------------------------------------------------------------------*/
958 
972 /*-------------------------------------------------------------------------------
973  Arguments
974  ______________________________________________________________________________.
975  mode name role !
976  ______________________________________________________________________________*/
987 /*-------------------------------------------------------------------------------*/
988 
989 void
991  cs_real_t prl,
992  cs_real_t prt,
993  cs_real_t yplus,
994  cs_real_t dplus,
995  cs_real_t *htur,
996  cs_real_t *yplim);
997 
998 /*----------------------------------------------------------------------------*/
999 
1001 
1002 #endif /* __CS_WALL_FUNCTIONS_H__ */
const double cs_turb_xkappa
Definition: cs_turbulence_model.c:272
double precision epzero
epsilon
Definition: cstnum.f90:40
cs_wall_f_type_t iwallf
Definition: cs_wall_functions.h:81
static void cs_wall_functions_2scales_scalable(cs_real_t l_visc, cs_real_t t_visc, cs_real_t vel, cs_real_t y, cs_real_t kinetic_en, int *iuntur, cs_lnum_t *nsubla, cs_lnum_t *nlogla, cs_real_t *ustar, cs_real_t *uk, cs_real_t *yplus, cs_real_t *dplus, cs_real_t *ypup, cs_real_t *cofimp)
Scalable wall function: shift the wall if .
Definition: cs_wall_functions.h:385
void cs_wall_functions_velocity(cs_wall_f_type_t iwallf, cs_lnum_t ifac, cs_real_t l_visc, cs_real_t t_visc, cs_real_t vel, cs_real_t y, cs_real_t rnnb, cs_real_t kinetic_en, int *iuntur, cs_lnum_t *nsubla, cs_lnum_t *nlogla, cs_real_t *ustar, cs_real_t *uk, cs_real_t *yplus, cs_real_t *ypup, cs_real_t *cofimp, cs_real_t *dplus)
Compute the friction velocity and / .
Definition: cs_wall_functions.c:288
#define BEGIN_C_DECLS
Definition: cs_defs.h:419
integer(c_int), pointer, save iwalfs
Wall functions for scalar.
Definition: optcal.f90:562
int cs_int_t
Fortran-compatible integer.
Definition: cs_defs.h:295
double cs_turb_cmu025
Definition: cs_turbulence_model.c:312
static cs_real_t _vdriest_dupdyp_integral(cs_real_t yplus)
Definition: cs_wall_functions.h:459
Definition: cs_wall_functions.h:71
Definition: cs_wall_functions.h:61
double precision, dimension(ncharm), save a2
Definition: cpincl.f90:226
void hturbp(const cs_int_t *const iwalfs, const cs_real_t *const prl, const cs_real_t *const prt, const cs_real_t *const yplus, const cs_real_t *const dplus, cs_real_t *htur, cs_real_t *yplim)
Definition: cs_wall_functions.c:231
Definition: cs_wall_functions.h:60
double cs_turb_dpow
Definition: cs_turbulence_model.c:301
const double cs_turb_bpow
Definition: cs_turbulence_model.c:298
Definition: cs_wall_functions.h:62
static void cs_wall_functions_s_arpaci_larsen(cs_real_t prl, cs_real_t prt, cs_real_t yplus, cs_real_t dplus, cs_real_t *htur, cs_real_t *yplim)
The correction of the exchange coefficient is computed thanks to a similarity model between dynamic v...
Definition: cs_wall_functions.h:721
cs_wall_f_s_type_t
Definition: cs_wall_functions.h:69
Definition: cs_wall_functions.h:64
const double cs_turb_vdriest
Definition: cs_turbulence_model.c:281
static void cs_wall_functions_s_vdriest(cs_real_t prl, cs_real_t prt, cs_real_t yplus, cs_real_t *htur)
The correction of the exchange coefficient is computed thanks to a numerical integration of: with ...
Definition: cs_wall_functions.h:810
#define CS_MIN(a, b)
Definition: cs_defs.h:392
real(c_double), pointer, save ypluli
limit value of for the viscous sublayer. ypluli depends on the chosen wall function: it is initializ...
Definition: cstphy.f90:249
integer(c_int), pointer, save iwallf
Wall functions.
Definition: optcal.f90:557
void wallfunctions(const cs_int_t *const iwallf, const cs_lnum_t *const ifac, const cs_real_t *const viscosity, const cs_real_t *const t_visc, const cs_real_t *const vel, const cs_real_t *const y, const cs_real_t *const rnnb, const cs_real_t *const kinetic_en, cs_int_t *iuntur, cs_lnum_t *nsubla, cs_lnum_t *nlogla, cs_real_t *ustar, cs_real_t *uk, cs_real_t *yplus, cs_real_t *ypup, cs_real_t *cofimp, cs_real_t *dplus)
Definition: cs_wall_functions.c:185
Definition: cs_wall_functions.h:63
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:292
static void cs_wall_functions_1scale_power(cs_real_t l_visc, cs_real_t vel, cs_real_t y, int *iuntur, cs_lnum_t *nsubla, cs_lnum_t *nlogla, cs_real_t *ustar, cs_real_t *uk, cs_real_t *yplus, cs_real_t *ypup, cs_real_t *cofimp)
Power law: Werner & Wengle.
Definition: cs_wall_functions.h:127
static void cs_wall_functions_2scales_vdriest(cs_real_t rnnb, cs_real_t l_visc, cs_real_t vel, cs_real_t y, cs_real_t kinetic_en, int *iuntur, cs_lnum_t *nsubla, cs_lnum_t *nlogla, cs_real_t *ustar, cs_real_t *uk, cs_real_t *yplus, cs_real_t *ypup, cs_real_t *cofimp, cs_real_t *lmk, cs_real_t kr, bool wf)
Two velocity scales wall function using Van Driest mixing length.
Definition: cs_wall_functions.h:533
#define END_C_DECLS
Definition: cs_defs.h:420
wall functions descriptor.
Definition: cs_wall_functions.h:79
#define _(String)
Definition: cs_defs.h:52
double cs_real_t
Definition: cs_defs.h:296
const double cs_turb_cstlog
Definition: cs_turbulence_model.c:292
#define CS_PROCF(x, y)
Definition: cs_defs.h:443
#define CS_MAX(a, b)
Definition: cs_defs.h:393
int iwallt
Definition: cs_wall_functions.h:85
static void cs_wall_functions_2scales_log(cs_real_t l_visc, cs_real_t t_visc, cs_real_t vel, cs_real_t y, cs_real_t kinetic_en, int *iuntur, cs_lnum_t *nsubla, cs_lnum_t *nlogla, cs_real_t *ustar, cs_real_t *uk, cs_real_t *yplus, cs_real_t *ypup, cs_real_t *cofimp)
Log law: piecewise linear and log, with two velocity scales based on the friction and the turbulent k...
Definition: cs_wall_functions.h:302
static void cs_wall_functions_disabled(cs_real_t l_visc, cs_real_t t_visc, cs_real_t vel, cs_real_t y, int *iuntur, cs_lnum_t *nsubla, cs_lnum_t *nlogla, cs_real_t *ustar, cs_real_t *uk, cs_real_t *yplus, cs_real_t *dplus, cs_real_t *ypup, cs_real_t *cofimp)
No wall function.
Definition: cs_wall_functions.h:654
int bft_printf(const char *const format,...)
Replacement for printf() with modifiable behavior.
Definition: bft_printf.c:140
cs_wall_f_s_type_t iwalfs
Definition: cs_wall_functions.h:83
Definition: cs_field_pointer.h:70
double ypluli
Definition: cs_wall_functions.h:90
cs_wall_f_type_t
Definition: cs_wall_functions.h:58
Definition: cs_wall_functions.h:72
const double cs_turb_crij1
Definition: cs_turbulence_model.c:356
const double cs_turb_crij2
Definition: cs_turbulence_model.c:362
static void cs_wall_functions_1scale_log(cs_lnum_t ifac, cs_real_t l_visc, cs_real_t vel, cs_real_t y, int *iuntur, cs_lnum_t *nsubla, cs_lnum_t *nlogla, cs_real_t *ustar, cs_real_t *uk, cs_real_t *yplus, cs_real_t *ypup, cs_real_t *cofimp)
Log law: piecewise linear and log, with one velocity scale based on the friction. ...
Definition: cs_wall_functions.h:198
const cs_wall_functions_t * cs_glob_wall_functions
void cs_wall_functions_scalar(cs_wall_f_s_type_t iwalfs, cs_real_t prl, cs_real_t prt, cs_real_t yplus, cs_real_t dplus, cs_real_t *htur, cs_real_t *yplim)
Compute the correction of the exchange coefficient between the fluid and the wall for a turbulent flo...
Definition: cs_wall_functions.c:446
Definition: cs_wall_functions.h:65
const double cs_turb_apow
Definition: cs_turbulence_model.c:295