8.2
general documentation
cs_convection_diffusion_priv.h
Go to the documentation of this file.
1 #ifndef __CS_CONVECTION_DIFFUSION_PRIV_H__
2 #define __CS_CONVECTION_DIFFUSION_PRIV_H__
3 
4 /*============================================================================
5  * Private functions for convection-diffusion operators.
6  *============================================================================*/
7 
8 /*
9  This file is part of code_saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2024 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 #include "cs_defs.h"
31 
32 /*----------------------------------------------------------------------------
33  * Local headers
34  *----------------------------------------------------------------------------*/
35 
36 #include "cs_math.h"
37 #include "cs_parameters.h" // for BC types
38 
40 
41 /*=============================================================================
42  * Macro definitions
43  *============================================================================*/
44 
45 /*============================================================================
46  * Type definition
47  *============================================================================*/
48 
49 /*============================================================================
50  * Global variables
51  *============================================================================*/
52 
53 /*============================================================================
54  * Public inlined function
55  *============================================================================*/
56 
57 /*----------------------------------------------------------------------------
58  * Synchronize halos for scalar variables.
59  *
60  * parameters:
61  * m <-- pointer to associated mesh structure
62  * halo_type <-> halo type
63  * pvar <-> variable
64  *----------------------------------------------------------------------------*/
65 
66 inline static void
68  cs_halo_type_t halo_type,
69  cs_real_t pvar[])
70 {
71  if (m->halo != NULL)
72  cs_halo_sync_var(m->halo, halo_type, pvar);
73 }
74 
75 /*----------------------------------------------------------------------------*/
86 /*----------------------------------------------------------------------------*/
87 
88 CS_F_HOST_DEVICE inline static cs_real_t
90  const cs_real_t nvf_p_c,
91  const cs_real_t nvf_r_f,
92  const cs_real_t nvf_r_c)
93 {
94  cs_real_t nvf_p_f;
95 
96  cs_real_t beta_m, rfc, r1f, r1, r2, r3, b1, b2;
97 
98  switch (limiter) {
99  case CS_NVD_GAMMA: /* Gamma scheme */
100  beta_m = 0.1; /* in [0.1, 0.5] */
101  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
102 
103  if (nvf_p_c < beta_m) {
104  nvf_p_f = nvf_p_c*(1.+rfc*(1.-nvf_p_c)/beta_m);
105  } else {
106  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
107 
108  nvf_p_f = r1f*nvf_p_c+rfc;
109  }
110 
111  break;
112 
113  case CS_NVD_SMART: /* SMART scheme */
114  if (nvf_p_c < (nvf_r_c/3.)) {
115  r1 = nvf_r_f*(1.-3.*nvf_r_c+2.*nvf_r_f);
116  r2 = nvf_r_c*(1.-nvf_r_c);
117 
118  nvf_p_f = nvf_p_c*r1/r2;
119  } else if (nvf_p_c <= (nvf_r_c*(1.+nvf_r_f-nvf_r_c)/nvf_r_f)) {
120  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
121  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
122 
123  nvf_p_f = nvf_r_f*(r1f*nvf_p_c/nvf_r_c + rfc);
124  } else {
125  nvf_p_f = 1.;
126  }
127 
128  break;
129 
130  case CS_NVD_CUBISTA: /* CUBISTA scheme */
131  if (nvf_p_c < (3.*nvf_r_c/4.)) {
132  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
133 
134  nvf_p_f = nvf_r_f*(1.+rfc/3.)*nvf_p_c/nvf_r_c;
135  } else if (nvf_p_c <= (nvf_r_c*(1.+2.*(nvf_r_f-nvf_r_c))/(2.*nvf_r_f-nvf_r_c))) {
136  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
137  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
138 
139  nvf_p_f = nvf_r_f*(r1f*nvf_p_c/nvf_r_c+rfc);
140  } else {
141  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
142 
143  nvf_p_f = 1.-.5*r1f*(1.-nvf_p_c);
144  }
145 
146  break;
147 
148  case CS_NVD_SUPERBEE: /* SuperBee scheme */
149  if (nvf_p_c < (nvf_r_c/(2.-nvf_r_c))) {
150  nvf_p_f = (2.*nvf_r_f-nvf_r_c)*nvf_p_c/nvf_r_c;
151  } else if (nvf_p_c < nvf_r_c) {
152  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
153  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
154 
155  nvf_p_f = r1f*nvf_p_c+rfc;
156  } else if (nvf_p_c < (nvf_r_c/nvf_r_f)) {
157  nvf_p_f = nvf_r_f*nvf_p_c/nvf_r_c;
158  } else {
159  nvf_p_f = 1.;
160  }
161 
162  break;
163 
164  case CS_NVD_MUSCL: /* MUSCL scheme */
165  if (nvf_p_c < (.5*nvf_r_c)) {
166  nvf_p_f = (2.*nvf_r_f-nvf_r_c)*nvf_p_c/nvf_r_c;
167  } else if (nvf_p_c < (1.+nvf_r_c-nvf_r_f)) {
168  nvf_p_f = nvf_p_c+nvf_r_f-nvf_r_c;
169  } else {
170  nvf_p_f = 1.;
171  }
172 
173  break;
174 
175  case CS_NVD_MINMOD: /* MINMOD scheme */
176  if (nvf_p_c < nvf_r_c) {
177  nvf_p_f = nvf_r_f*nvf_p_c/nvf_r_c;
178  } else {
179  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
180  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
181 
182  nvf_p_f = r1f*nvf_p_c+rfc;
183  }
184 
185  break;
186 
187  case CS_NVD_CLAM: /* CLAM scheme */
188  r1 = nvf_r_c*nvf_r_c-nvf_r_f;
189  r2 = nvf_r_c*(nvf_r_c-1.);
190  r3 = nvf_r_f-nvf_r_c;
191 
192  nvf_p_f = nvf_p_c*(r1+r3*nvf_p_c)/r2;
193 
194  break;
195 
196  case CS_NVD_STOIC: /* STOIC scheme */
197  b1 = (nvf_r_c-nvf_r_f)*nvf_r_c;
198  b2 = nvf_r_c+nvf_r_f+2.*nvf_r_f*nvf_r_f-4.*nvf_r_f*nvf_r_c;
199 
200  if (nvf_p_c < (b1/b2)) {
201  r1 = -nvf_r_f*(1.-3.*nvf_r_c+2.*nvf_r_f);
202  r2 = nvf_r_c*(nvf_r_c-1.);
203 
204  nvf_p_f = nvf_p_c*r1/r2;
205  } else if (nvf_p_c < nvf_r_c) {
206  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
207  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
208 
209  nvf_p_f = r1f*nvf_p_c+rfc;
210  } else if (nvf_p_c < (nvf_r_c*(1.+nvf_r_f-nvf_r_c)/nvf_r_f)) {
211  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
212  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
213 
214  nvf_p_f = nvf_r_f*(nvf_p_c*r1f/nvf_r_c+rfc);
215  } else {
216  nvf_p_f = 1.;
217  }
218 
219  break;
220 
221  case CS_NVD_OSHER: /* OSHER scheme */
222  if (nvf_p_c < (nvf_r_c/nvf_r_f)) {
223  nvf_p_f = nvf_r_f*nvf_p_c/nvf_r_c;
224  } else {
225  nvf_p_f = 1.;
226  }
227 
228  break;
229 
230  case CS_NVD_WASEB: /* WASEB scheme */
231  r1 = nvf_r_c*nvf_r_f*(nvf_r_f-nvf_r_c);
232  r2 = 2.*nvf_r_c*(1.-nvf_r_c)-nvf_r_f*(1.-nvf_r_f);
233 
234  if (nvf_p_c < (r1/r2)) {
235  nvf_p_f = 2.*nvf_p_c;
236  } else if (nvf_p_c <= (nvf_r_c*(1.+nvf_r_f-nvf_r_c)/nvf_r_f)) {
237  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
238  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
239 
240  nvf_p_f = nvf_r_f*(nvf_p_c*r1f/nvf_r_c+rfc);
241  } else {
242  nvf_p_f = 1.;
243  }
244 
245  break;
246 
247  default: /* Upwinding */
248  nvf_p_f = nvf_p_c;
249  break;
250  }
251 
252  return nvf_p_f;
253 }
254 
255 /*----------------------------------------------------------------------------*/
271 /*----------------------------------------------------------------------------*/
272 
273 CS_F_HOST_DEVICE inline static cs_real_t
275  const cs_real_t i_face_u_normal[3],
276  cs_real_t nvf_p_c,
277  cs_real_t nvf_r_f,
278  cs_real_t nvf_r_c,
279  const cs_real_t gradv_c[3],
280  const cs_real_t c_courant)
281 {
282  assert(limiter >= CS_NVD_VOF_HRIC);
283 
284  cs_real_t nvf_p_f;
285  cs_real_t blend, high_order, low_order, ratio;
286 
287  /* Compute gradient angle indicator */
288  cs_real_t dotp = cs_math_fabs(cs_math_3_dot_product(gradv_c, i_face_u_normal));
289  cs_real_t sgrad = cs_math_3_norm(gradv_c);
290  cs_real_t denom = sgrad;
291 
292  if (limiter == CS_NVD_VOF_HRIC) { /* M-HRIC scheme */
293  /* High order scheme : Bounded Downwind */
294  if (nvf_p_c <= .5) {
295  high_order = 2.*nvf_p_c;
296  } else {
297  high_order = 1.;
298  }
299 
300  /* Low order scheme : MUSCL */
302  nvf_p_c,
303  nvf_r_f,
304  nvf_r_c);
305 
306  /* Compute the blending factor */
307  if (denom <= (cs_math_epzero*dotp)) {
308  blend = 1.;
309  } else {
310  ratio = dotp/denom;
311  blend = cs_math_fmin(1., pow(ratio, .5));
312  }
313 
314  /* Blending */
315  nvf_p_f = blend*high_order + (1.-blend)*low_order;
316 
317  /* Extra blending due to the cell Courant number */
318  if (c_courant < .7 && c_courant > .3) {
319  nvf_p_f = nvf_p_f + (nvf_p_f - low_order)*(.7 - c_courant )/.4;
320  } else if (c_courant >= .7) {
321  nvf_p_f = low_order;
322  }
323  } else if (limiter == CS_NVD_VOF_CICSAM) { /* M-CICSAM scheme */
324  /* High order scheme : HYPER-C + SUPERBEE */
325  if (c_courant <= .3) {
326  high_order = cs_math_fmin(1., nvf_p_c/(c_courant+cs_math_epzero));
327  } else if (c_courant <= .6) {
328  high_order = cs_math_fmin(1., nvf_p_c/.3);
329  } else if (c_courant <= .7) {
331  nvf_p_c,
332  nvf_r_f,
333  nvf_r_c);
334  high_order = 10.*( (.7-c_courant)*cs_math_fmin(1., nvf_p_c/.3)
335  + (c_courant-.6)*superbee);
336  }
337  else {
339  nvf_p_c,
340  nvf_r_f,
341  nvf_r_c);
342  }
343 
344  /* Low order scheme : MUSCL */
346  nvf_p_c,
347  nvf_r_f,
348  nvf_r_c);
349 
350  /* Compute the blending factor */
351  if (denom <= (cs_math_epzero*dotp)) {
352  blend = 1.;
353  } else {
354  ratio = dotp/denom;
355  blend = cs_math_fmin(1., pow(ratio, 2.));
356  }
357 
358  /* Blending */
359  nvf_p_f = blend*high_order + (1.-blend)*low_order;
360  } else { /* STACS scheme */
361  /* High order scheme : SUPERBEE */
363  nvf_p_c,
364  nvf_r_f,
365  nvf_r_c);
366 
367  /* Low order scheme : STOIC */
369  nvf_p_c,
370  nvf_r_f,
371  nvf_r_c);
372 
373  /* Compute the blending factor */
374  if (denom <= (cs_math_epzero*dotp)) {
375  blend = 1.;
376  } else {
377  ratio = dotp/denom;
378  blend = cs_math_fmin(1., pow(ratio, 4.));
379  }
380 
381  /* Blending */
382  nvf_p_f = blend*high_order + (1.-blend)*low_order;
383  }
384 
385  return nvf_p_f;
386 }
387 
388 /*----------------------------------------------------------------------------*/
404 /*----------------------------------------------------------------------------*/
405 
406 CS_F_HOST_DEVICE inline static void
408  const cs_real_t pj,
409  const cs_real_t distf,
410  const cs_real_t i_face_u_normal[3],
411  const cs_real_t gradi[3],
412  const cs_real_t gradj[3],
413  const cs_real_t grdpai[3],
414  const cs_real_t grdpaj[3],
415  const cs_real_t i_massflux,
416  cs_real_t *testij,
417  cs_real_t *tesqck)
418 {
419  cs_real_t dcc, ddi, ddj;
420 
421  /* Slope test */
422 
423  *testij = grdpai[0]*grdpaj[0]
424  + grdpai[1]*grdpaj[1]
425  + grdpai[2]*grdpaj[2];
426 
427  if (i_massflux > 0.) {
428  dcc = gradi[0]*i_face_u_normal[0]
429  + gradi[1]*i_face_u_normal[1]
430  + gradi[2]*i_face_u_normal[2];
431  ddi = grdpai[0]*i_face_u_normal[0]
432  + grdpai[1]*i_face_u_normal[1]
433  + grdpai[2]*i_face_u_normal[2];
434  ddj = (pj-pi)/distf;
435  }
436  else {
437  dcc = gradj[0]*i_face_u_normal[0]
438  + gradj[1]*i_face_u_normal[1]
439  + gradj[2]*i_face_u_normal[2];
440  ddi = (pj-pi)/distf;
441  ddj = grdpaj[0]*i_face_u_normal[0]
442  + grdpaj[1]*i_face_u_normal[1]
443  + grdpaj[2]*i_face_u_normal[2];
444  }
445 
446  *tesqck = cs_math_sq(dcc) - cs_math_sq(ddi-ddj);
447 }
448 
449 /*----------------------------------------------------------------------------*/
468 /*----------------------------------------------------------------------------*/
469 
470 template <cs_lnum_t stride>
471 CS_F_HOST_DEVICE inline static void
473  const cs_real_t pj[stride],
474  const cs_real_t distf,
475  const cs_real_t i_face_u_normal[3],
476  const cs_real_t gradi[stride][3],
477  const cs_real_t gradj[stride][3],
478  const cs_real_t gradsti[stride][3],
479  const cs_real_t gradstj[stride][3],
480  const cs_real_t i_massflux,
481  cs_real_t *testij,
482  cs_real_t *tesqck)
483 {
484  cs_real_t dcc[stride], ddi[stride], ddj[stride];
485 
486  *testij = 0.;
487  *tesqck = 0.;
488 
489  /* Slope test */
490 
491  for (int i = 0; i < stride; i++) {
492  *testij += gradsti[i][0]*gradstj[i][0]
493  + gradsti[i][1]*gradstj[i][1]
494  + gradsti[i][2]*gradstj[i][2];
495 
496  if (i_massflux > 0.) {
497  dcc[i] = gradi[i][0]*i_face_u_normal[0]
498  + gradi[i][1]*i_face_u_normal[1]
499  + gradi[i][2]*i_face_u_normal[2];
500  ddi[i] = gradsti[i][0]*i_face_u_normal[0]
501  + gradsti[i][1]*i_face_u_normal[1]
502  + gradsti[i][2]*i_face_u_normal[2];
503  ddj[i] = (pj[i]-pi[i])/distf;
504  }
505  else {
506  dcc[i] = gradj[i][0]*i_face_u_normal[0]
507  + gradj[i][1]*i_face_u_normal[1]
508  + gradj[i][2]*i_face_u_normal[2];
509  ddi[i] = (pj[i]-pi[i])/distf;
510  ddj[i] = gradstj[i][0]*i_face_u_normal[0]
511  + gradstj[i][1]*i_face_u_normal[1]
512  + gradstj[i][2]*i_face_u_normal[2];
513  }
514 
515  *tesqck += cs_math_sq(dcc[i]) - cs_math_sq(ddi[i]-ddj[i]);
516  }
517 }
518 
519 /*----------------------------------------------------------------------------*/
535 /*----------------------------------------------------------------------------*/
536 
537 CS_F_HOST_DEVICE inline static void
539  const cs_real_3_t diipf,
540  const cs_real_3_t djjpf,
541  const cs_real_3_t gradi,
542  const cs_real_3_t gradj,
543  const cs_real_t pi,
544  const cs_real_t pj,
545  cs_real_t *recoi,
546  cs_real_t *recoj,
547  cs_real_t *pip,
548  cs_real_t *pjp)
549 {
550  cs_real_t gradpf[3] = {0.5*(gradi[0] + gradj[0]),
551  0.5*(gradi[1] + gradj[1]),
552  0.5*(gradi[2] + gradj[2])};
553 
554  /* reconstruction only if IRCFLP = 1 */
555  *recoi = bldfrp*(cs_math_3_dot_product(gradpf, diipf));
556  *recoj = bldfrp*(cs_math_3_dot_product(gradpf, djjpf));
557  *pip = pi + *recoi;
558  *pjp = pj + *recoj;
559 }
560 
561 /*----------------------------------------------------------------------------*/
580 /*----------------------------------------------------------------------------*/
581 
582 template <cs_lnum_t stride>
583 CS_F_HOST_DEVICE inline static void
585  const cs_real_t diipf[3],
586  const cs_real_t djjpf[3],
587  const cs_real_t gradi[stride][3],
588  const cs_real_t gradj[stride][3],
589  const cs_real_t pi[stride],
590  const cs_real_t pj[stride],
591  cs_real_t recoi[stride],
592  cs_real_t recoj[stride],
593  cs_real_t pip[stride],
594  cs_real_t pjp[stride])
595 {
596  cs_real_t dpvf[3];
597 
598  /* x-y-z components, p = u, v, w */
599 
600  for (int isou = 0; isou < stride; isou++) {
601 
602  for (int jsou = 0; jsou < 3; jsou++)
603  dpvf[jsou] = 0.5*( gradi[isou][jsou]
604  + gradj[isou][jsou]);
605 
606  /* reconstruction only if IRCFLP = 1 */
607 
608  recoi[isou] = bldfrp*(cs_math_3_dot_product(dpvf, diipf));
609  recoj[isou] = bldfrp*(cs_math_3_dot_product(dpvf, djjpf));
610 
611  pip[isou] = pi[isou] + recoi[isou];
612  pjp[isou] = pj[isou] + recoj[isou];
613 
614  }
615 }
616 
617 /*----------------------------------------------------------------------------*/
633 /*----------------------------------------------------------------------------*/
634 
635 CS_F_HOST_DEVICE inline static void
636 cs_i_relax_c_val(const double relaxp,
637  const cs_real_t pia,
638  const cs_real_t pja,
639  const cs_real_t recoi,
640  const cs_real_t recoj,
641  const cs_real_t pi,
642  const cs_real_t pj,
643  cs_real_t *pir,
644  cs_real_t *pjr,
645  cs_real_t *pipr,
646  cs_real_t *pjpr)
647 {
648  *pir = pi/relaxp - (1.-relaxp)/relaxp * pia;
649  *pjr = pj/relaxp - (1.-relaxp)/relaxp * pja;
650 
651  *pipr = *pir + recoi;
652  *pjpr = *pjr + recoj;
653 }
654 
655 /*----------------------------------------------------------------------------*/
674 /*----------------------------------------------------------------------------*/
675 
676 template <cs_lnum_t stride>
677 CS_F_HOST_DEVICE inline static void
679  const cs_real_t pia[stride],
680  const cs_real_t pja[stride],
681  const cs_real_t recoi[stride],
682  const cs_real_t recoj[stride],
683  const cs_real_t pi[stride],
684  const cs_real_t pj[stride],
685  cs_real_t pir[stride],
686  cs_real_t pjr[stride],
687  cs_real_t pipr[stride],
688  cs_real_t pjpr[stride])
689 {
690  for (int isou = 0; isou < stride; isou++) {
691  pir[isou] = pi[isou] /relaxp - (1.-relaxp)/relaxp * pia[isou];
692  pjr[isou] = pj[isou] /relaxp - (1.-relaxp)/relaxp * pja[isou];
693 
694  pipr[isou] = pir[isou] + recoi[isou];
695  pjpr[isou] = pjr[isou] + recoj[isou];
696  }
697 }
698 
699 /*----------------------------------------------------------------------------*/
706 /*----------------------------------------------------------------------------*/
707 
708 CS_F_HOST_DEVICE inline static void
710  cs_real_t *pf)
711 {
712  *pf = p;
713 }
714 
715 /*----------------------------------------------------------------------------*/
725 /*----------------------------------------------------------------------------*/
726 
727 template <cs_lnum_t stride>
728 CS_F_HOST_DEVICE inline static void
730  cs_real_t pf[stride])
731 {
732  for (int isou = 0; isou < stride; isou++)
733  pf[isou] = p[isou];
734 }
735 
736 /*----------------------------------------------------------------------------*/
745 /*----------------------------------------------------------------------------*/
746 
747 CS_F_HOST_DEVICE inline static void
748 cs_centered_f_val(double pnd,
749  cs_real_t pip,
750  cs_real_t pjp,
751  cs_real_t *pf)
752 {
753  *pf = pnd*pip + (1.-pnd)*pjp;
754 }
755 
756 /*----------------------------------------------------------------------------*/
768 /*----------------------------------------------------------------------------*/
769 
770 template <cs_lnum_t stride>
771 CS_F_HOST_DEVICE inline static void
773  const cs_real_t pip[stride],
774  const cs_real_t pjp[stride],
775  cs_real_t pf[stride])
776 {
777  for (int isou = 0; isou < stride; isou++)
778  pf[isou] = pnd*pip[isou] + (1.-pnd)*pjp[isou];
779 }
780 
781 /*----------------------------------------------------------------------------*/
791 /*----------------------------------------------------------------------------*/
792 
793 CS_F_HOST_DEVICE inline static void
794 cs_solu_f_val(const cs_real_t cell_cen[3],
795  const cs_real_t i_face_cog[3],
796  const cs_real_t grad[3],
797  const cs_real_t p,
798  cs_real_t *pf)
799 {
800  cs_real_t df[3];
801 
802  df[0] = i_face_cog[0] - cell_cen[0];
803  df[1] = i_face_cog[1] - cell_cen[1];
804  df[2] = i_face_cog[2] - cell_cen[2];
805 
806  *pf = p + cs_math_3_dot_product(df, grad);
807 }
808 
809 /*----------------------------------------------------------------------------*/
822 /*----------------------------------------------------------------------------*/
823 
824 template <cs_lnum_t stride>
825 CS_F_HOST_DEVICE inline static void
827  const cs_real_t i_face_cog[3],
828  const cs_real_t grad[stride][3],
829  const cs_real_t p[stride],
830  cs_real_t pf[stride])
831 {
832  cs_real_t df[3];
833 
834  for (cs_lnum_t jsou = 0; jsou < 3; jsou++)
835  df[jsou] = i_face_cog[jsou] - cell_cen[jsou];
836 
837  for (cs_lnum_t isou = 0; isou < stride; isou++) {
838  pf[isou] = p[isou] + df[0]*grad[isou][0]
839  + df[1]*grad[isou][1]
840  + df[2]*grad[isou][2];
841  }
842 }
843 
844 /*----------------------------------------------------------------------------*/
854 /*----------------------------------------------------------------------------*/
855 
856 CS_F_HOST_DEVICE inline static void
857 cs_blend_f_val(const double blencp,
858  const cs_real_t p,
859  cs_real_t *pf)
860 {
861  *pf = blencp * (*pf) + (1. - blencp) * p;
862 }
863 
864 /*----------------------------------------------------------------------------*/
877 /*----------------------------------------------------------------------------*/
878 
879 template <cs_lnum_t stride>
880 CS_F_HOST_DEVICE inline static void
881 cs_blend_f_val_strided(const double blencp,
882  const cs_real_t p[stride],
883  cs_real_t pf[stride])
884 {
885  for (int isou = 0; isou < stride; isou++)
886  pf[isou] = blencp*(pf[isou])+(1.-blencp)*p[isou];
887 }
888 
889 /*----------------------------------------------------------------------------*/
910 /*----------------------------------------------------------------------------*/
911 
912 CS_F_HOST_DEVICE inline static void
913 cs_i_conv_flux(const int iconvp,
914  const cs_real_t thetap,
915  const int imasac,
916  const cs_real_t pi,
917  const cs_real_t pj,
918  const cs_real_t pifri,
919  const cs_real_t pifrj,
920  const cs_real_t pjfri,
921  const cs_real_t pjfrj,
922  const cs_real_t i_massflux,
923  const cs_real_t xcppi,
924  const cs_real_t xcppj,
925  cs_real_2_t fluxij)
926 {
927  cs_real_t flui, fluj;
928 
929  flui = 0.5*(i_massflux + fabs(i_massflux));
930  fluj = 0.5*(i_massflux - fabs(i_massflux));
931 
932  fluxij[0] += iconvp*xcppi*(thetap*(flui*pifri + fluj*pjfri) - imasac*i_massflux*pi);
933  fluxij[1] += iconvp*xcppj*(thetap*(flui*pifrj + fluj*pjfrj) - imasac*i_massflux*pj);
934 }
935 
936 /*----------------------------------------------------------------------------*/
957 /*----------------------------------------------------------------------------*/
958 
959 template <cs_lnum_t stride>
960 CS_F_HOST_DEVICE inline static void
962  cs_real_t thetap,
963  int imasac,
964  const cs_real_t pi[stride],
965  const cs_real_t pj[stride],
966  const cs_real_t pifri[stride],
967  const cs_real_t pifrj[stride],
968  const cs_real_t pjfri[stride],
969  const cs_real_t pjfrj[stride],
970  cs_real_t i_massflux,
971  cs_real_t fluxi[stride],
972  cs_real_t fluxj[stride])
973 {
974  cs_real_t flui, fluj;
975 
976  flui = 0.5*(i_massflux + fabs(i_massflux));
977  fluj = 0.5*(i_massflux - fabs(i_massflux));
978 
979  for (int isou = 0; isou < stride; isou++) {
980  fluxi[isou] += iconvp*( thetap*(flui*pifri[isou] + fluj*pjfri[isou])
981  - imasac*i_massflux*pi[isou]);
982  fluxj[isou] += iconvp*( thetap*(flui*pifrj[isou] + fluj*pjfrj[isou])
983  - imasac*i_massflux*pj[isou]);
984  }
985 }
986 
987 /*----------------------------------------------------------------------------*/
1000 /*----------------------------------------------------------------------------*/
1001 
1002 CS_F_HOST_DEVICE inline static void
1003 cs_i_diff_flux(const int idiffp,
1004  const cs_real_t thetap,
1005  const cs_real_t pip,
1006  const cs_real_t pjp,
1007  const cs_real_t pipr,
1008  const cs_real_t pjpr,
1009  const cs_real_t i_visc,
1010  cs_real_t fluxij[2])
1011 {
1012  fluxij[0] += idiffp*thetap*i_visc*(pipr -pjp);
1013  fluxij[1] += idiffp*thetap*i_visc*(pip -pjpr);
1014 }
1015 
1016 /*----------------------------------------------------------------------------*/
1033 /*----------------------------------------------------------------------------*/
1034 
1035 template <cs_lnum_t stride>
1036 CS_F_HOST_DEVICE inline static void
1038  cs_real_t thetap,
1039  const cs_real_t pip[stride],
1040  const cs_real_t pjp[stride],
1041  const cs_real_t pipr[stride],
1042  const cs_real_t pjpr[stride],
1043  const cs_real_t i_visc,
1044  cs_real_t fluxi[stride],
1045  cs_real_t fluxj[stride])
1046 {
1047  for (int isou = 0; isou < stride; isou++) {
1048  fluxi[isou] += idiffp*thetap*i_visc*(pipr[isou] -pjp[isou]);
1049  fluxj[isou] += idiffp*thetap*i_visc*(pip[isou] -pjpr[isou]);
1050  }
1051 }
1052 
1053 /*----------------------------------------------------------------------------*/
1077 /*----------------------------------------------------------------------------*/
1078 
1079 CS_F_HOST_DEVICE inline static void
1081  const cs_real_t relaxp,
1082  const cs_real_t diipf[3],
1083  const cs_real_t djjpf[3],
1084  const cs_real_t gradi[3],
1085  const cs_real_t gradj[3],
1086  const cs_real_t pi,
1087  const cs_real_t pj,
1088  const cs_real_t pia,
1089  const cs_real_t pja,
1090  cs_real_t *pifri,
1091  cs_real_t *pifrj,
1092  cs_real_t *pjfri,
1093  cs_real_t *pjfrj,
1094  cs_real_t *pip,
1095  cs_real_t *pjp,
1096  cs_real_t *pipr,
1097  cs_real_t *pjpr)
1098 {
1099  cs_real_t pir, pjr;
1100  cs_real_t recoi, recoj;
1101 
1102  cs_i_compute_quantities(bldfrp,
1103  diipf,
1104  djjpf,
1105  gradi,
1106  gradj,
1107  pi,
1108  pj,
1109  &recoi,
1110  &recoj,
1111  pip,
1112  pjp);
1113 
1114  cs_i_relax_c_val(relaxp,
1115  pia,
1116  pja,
1117  recoi,
1118  recoj,
1119  pi,
1120  pj,
1121  &pir,
1122  &pjr,
1123  pipr,
1124  pjpr);
1125 
1127  pifrj);
1128  cs_upwind_f_val(pir,
1129  pifri);
1130  cs_upwind_f_val(pj,
1131  pjfri);
1132  cs_upwind_f_val(pjr,
1133  pjfrj);
1134 }
1135 
1136 /*----------------------------------------------------------------------------*/
1163 /*----------------------------------------------------------------------------*/
1164 
1165 template <cs_lnum_t stride>
1166 CS_F_HOST_DEVICE inline static void
1168  const cs_real_t relaxp,
1169  const cs_real_t diipf[3],
1170  const cs_real_t djjpf[3],
1171  const cs_real_t gradi[stride][3],
1172  const cs_real_t gradj[stride][3],
1173  const cs_real_t pi[stride],
1174  const cs_real_t pj[stride],
1175  const cs_real_t pia[stride],
1176  const cs_real_t pja[stride],
1177  cs_real_t pifri[stride],
1178  cs_real_t pifrj[stride],
1179  cs_real_t pjfri[stride],
1180  cs_real_t pjfrj[stride],
1181  cs_real_t pip[stride],
1182  cs_real_t pjp[stride],
1183  cs_real_t pipr[stride],
1184  cs_real_t pjpr[stride])
1185 {
1186  cs_real_t pir[stride], pjr[stride];
1187  cs_real_t recoi[stride], recoj[stride];
1188 
1189  cs_i_compute_quantities_strided<stride>(bldfrp,
1190  diipf,
1191  djjpf,
1192  gradi,
1193  gradj,
1194  pi,
1195  pj,
1196  recoi,
1197  recoj,
1198  pip,
1199  pjp);
1200 
1201  cs_i_relax_c_val_strided<stride>(relaxp,
1202  pia,
1203  pja,
1204  recoi,
1205  recoj,
1206  pi,
1207  pj,
1208  pir,
1209  pjr,
1210  pipr,
1211  pjpr);
1212 
1213  cs_upwind_f_val_strided<stride>(pi, pifrj);
1214  cs_upwind_f_val_strided<stride>(pir, pifri);
1215  cs_upwind_f_val_strided<stride>(pj, pjfri);
1216  cs_upwind_f_val_strided<stride>(pjr, pjfrj);
1217 }
1218 
1219 /*----------------------------------------------------------------------------*/
1236 /*----------------------------------------------------------------------------*/
1237 
1238 CS_F_HOST_DEVICE inline static void
1240  const cs_real_t diipf[3],
1241  const cs_real_t djjpf[3],
1242  const cs_real_t gradi[3],
1243  const cs_real_t gradj[3],
1244  const cs_real_t pi,
1245  const cs_real_t pj,
1246  cs_real_t *pif,
1247  cs_real_t *pjf,
1248  cs_real_t *pip,
1249  cs_real_t *pjp)
1250 {
1251  cs_real_t recoi, recoj;
1252 
1253  cs_i_compute_quantities(bldfrp,
1254  diipf,
1255  djjpf,
1256  gradi,
1257  gradj,
1258  pi,
1259  pj,
1260  &recoi,
1261  &recoj,
1262  pip,
1263  pjp);
1264 
1265  cs_upwind_f_val(pi, pif);
1266  cs_upwind_f_val(pj, pjf);
1267 }
1268 
1269 /*----------------------------------------------------------------------------*/
1289 /*----------------------------------------------------------------------------*/
1290 
1291 template <cs_lnum_t stride>
1292 CS_F_HOST_DEVICE inline static void
1294  const cs_real_t diipf[3],
1295  const cs_real_t djjpf[3],
1296  const cs_real_t gradi[stride][3],
1297  const cs_real_t gradj[stride][3],
1298  const cs_real_t pi[stride],
1299  const cs_real_t pj[stride],
1300  cs_real_t pif[stride],
1301  cs_real_t pjf[stride],
1302  cs_real_t pip[stride],
1303  cs_real_t pjp[stride])
1304 {
1305  cs_real_t recoi[stride], recoj[stride];
1306 
1307  cs_i_compute_quantities_strided<stride>(bldfrp,
1308  diipf,
1309  djjpf,
1310  gradi,
1311  gradj,
1312  pi,
1313  pj,
1314  recoi,
1315  recoj,
1316  pip,
1317  pjp);
1318 
1319  cs_upwind_f_val_strided<stride>(pi, pif);
1320  cs_upwind_f_val_strided<stride>(pj, pjf);
1321 }
1322 
1323 /*----------------------------------------------------------------------------*/
1356 /*----------------------------------------------------------------------------*/
1357 
1358 CS_F_HOST_DEVICE inline static void
1360  const int ischcp,
1361  const double relaxp,
1362  const double blencp,
1363  const cs_real_t weight,
1364  const cs_real_t cell_ceni[3],
1365  const cs_real_t cell_cenj[3],
1366  const cs_real_t i_face_cog[3],
1367  const cs_real_t diipf[3],
1368  const cs_real_t djjpf[3],
1369  const cs_real_t gradi[3],
1370  const cs_real_t gradj[3],
1371  const cs_real_t gradupi[3],
1372  const cs_real_t gradupj[3],
1373  const cs_real_t pi,
1374  const cs_real_t pj,
1375  const cs_real_t pia,
1376  const cs_real_t pja,
1377  cs_real_t *pifri,
1378  cs_real_t *pifrj,
1379  cs_real_t *pjfri,
1380  cs_real_t *pjfrj,
1381  cs_real_t *pip,
1382  cs_real_t *pjp,
1383  cs_real_t *pipr,
1384  cs_real_t *pjpr)
1385 {
1386  cs_real_t pir, pjr;
1387  cs_real_t recoi, recoj;
1388 
1389  cs_i_compute_quantities(bldfrp,
1390  diipf,
1391  djjpf,
1392  gradi,
1393  gradj,
1394  pi,
1395  pj,
1396  &recoi,
1397  &recoj,
1398  pip,
1399  pjp);
1400 
1401  cs_i_relax_c_val(relaxp,
1402  pia,
1403  pja,
1404  recoi,
1405  recoj,
1406  pi,
1407  pj,
1408  &pir,
1409  &pjr,
1410  pipr,
1411  pjpr);
1412 
1413  if (ischcp == 1) {
1414 
1415  /* Centered
1416  --------*/
1417 
1418  cs_centered_f_val(weight,
1419  *pip,
1420  *pjpr,
1421  pifrj);
1422  cs_centered_f_val(weight,
1423  *pipr,
1424  *pjp,
1425  pifri);
1426  cs_centered_f_val(weight,
1427  *pipr,
1428  *pjp,
1429  pjfri);
1430  cs_centered_f_val(weight,
1431  *pip,
1432  *pjpr,
1433  pjfrj);
1434 
1435  } else if (ischcp == 0) {
1436 
1437  /* Original SOLU
1438  --------------*/
1439 
1440  cs_solu_f_val(cell_ceni,
1441  i_face_cog,
1442  gradi,
1443  pi,
1444  pifrj);
1445  cs_solu_f_val(cell_ceni,
1446  i_face_cog,
1447  gradi,
1448  pir,
1449  pifri);
1450  cs_solu_f_val(cell_cenj,
1451  i_face_cog,
1452  gradj,
1453  pj,
1454  pjfri);
1455  cs_solu_f_val(cell_cenj,
1456  i_face_cog,
1457  gradj,
1458  pjr,
1459  pjfrj);
1460 
1461  } else {
1462 
1463  /* SOLU
1464  ----*/
1465 
1466  cs_solu_f_val(cell_ceni,
1467  i_face_cog,
1468  gradupi,
1469  pi,
1470  pifrj);
1471  cs_solu_f_val(cell_ceni,
1472  i_face_cog,
1473  gradupi,
1474  pir,
1475  pifri);
1476  cs_solu_f_val(cell_cenj,
1477  i_face_cog,
1478  gradupj,
1479  pj,
1480  pjfri);
1481  cs_solu_f_val(cell_cenj,
1482  i_face_cog,
1483  gradupj,
1484  pjr,
1485  pjfrj);
1486 
1487  }
1488 
1489  /* Blending
1490  --------*/
1491 
1492  cs_blend_f_val(blencp,
1493  pi,
1494  pifrj);
1495  cs_blend_f_val(blencp,
1496  pir,
1497  pifri);
1498  cs_blend_f_val(blencp,
1499  pj,
1500  pjfri);
1501  cs_blend_f_val(blencp,
1502  pjr,
1503  pjfrj);
1504 }
1505 
1506 /*----------------------------------------------------------------------------*/
1540 /*----------------------------------------------------------------------------*/
1541 
1542 template <cs_lnum_t stride>
1543 CS_F_HOST_DEVICE inline static void
1545  int ischcp,
1546  double relaxp,
1547  double blencp,
1548  cs_real_t weight,
1549  const cs_real_t cell_ceni[3],
1550  const cs_real_t cell_cenj[3],
1551  const cs_real_t i_face_cog[3],
1552  const cs_real_t diipf[3],
1553  const cs_real_t djjpf[3],
1554  const cs_real_t gradi[stride][3],
1555  const cs_real_t gradj[stride][3],
1556  const cs_real_t pi[stride],
1557  const cs_real_t pj[stride],
1558  const cs_real_t pia[stride],
1559  const cs_real_t pja[stride],
1560  cs_real_t pifri[stride],
1561  cs_real_t pifrj[stride],
1562  cs_real_t pjfri[stride],
1563  cs_real_t pjfrj[stride],
1564  cs_real_t pip[stride],
1565  cs_real_t pjp[stride],
1566  cs_real_t pipr[stride],
1567  cs_real_t pjpr[stride])
1568 
1569 {
1570  cs_real_t pir[stride], pjr[stride];
1571  cs_real_t recoi[stride], recoj[stride];
1572 
1573  cs_i_compute_quantities_strided<stride>(bldfrp,
1574  diipf,
1575  djjpf,
1576  gradi,
1577  gradj,
1578  pi,
1579  pj,
1580  recoi,
1581  recoj,
1582  pip,
1583  pjp);
1584 
1585  cs_i_relax_c_val_strided<stride>(relaxp,
1586  pia,
1587  pja,
1588  recoi,
1589  recoj,
1590  pi,
1591  pj,
1592  pir,
1593  pjr,
1594  pipr,
1595  pjpr);
1596 
1597  if (ischcp == 1) {
1598 
1599  /* Centered
1600  --------*/
1601 
1602  cs_centered_f_val_strided<stride>(weight, pip, pjpr, pifrj);
1603  cs_centered_f_val_strided<stride>(weight, pipr, pjp, pifri);
1604  cs_centered_f_val_strided<stride>(weight, pipr, pjp, pjfri);
1605  cs_centered_f_val_strided<stride>(weight, pip, pjpr, pjfrj);
1606 
1607  }
1608  else {
1609 
1610  /* Second order
1611  ------------*/
1612 
1613  cs_solu_f_val_strided<stride>(cell_ceni,
1614  i_face_cog,
1615  gradi,
1616  pi,
1617  pifrj);
1618  cs_solu_f_val_strided<stride>(cell_ceni,
1619  i_face_cog,
1620  gradi,
1621  pir,
1622  pifri);
1623  cs_solu_f_val_strided<stride>(cell_cenj,
1624  i_face_cog,
1625  gradj,
1626  pj,
1627  pjfri);
1628  cs_solu_f_val_strided<stride>(cell_cenj,
1629  i_face_cog,
1630  gradj,
1631  pjr,
1632  pjfrj);
1633 
1634  }
1635 
1636  /* Blending
1637  --------*/
1638 
1639  cs_blend_f_val_strided<stride>(blencp, pi, pifrj);
1640  cs_blend_f_val_strided<stride>(blencp, pir, pifri);
1641  cs_blend_f_val_strided<stride>(blencp, pj, pjfri);
1642  cs_blend_f_val_strided<stride>(blencp, pjr, pjfrj);
1643 }
1644 
1645 /*----------------------------------------------------------------------------*/
1675 /*----------------------------------------------------------------------------*/
1676 
1677 CS_F_HOST_DEVICE inline static void
1679  const int ischcp,
1680  const double blencp,
1681  const cs_real_t weight,
1682  const cs_real_3_t cell_ceni,
1683  const cs_real_3_t cell_cenj,
1684  const cs_real_3_t i_face_cog,
1685  const cs_real_t hybrid_blend_i,
1686  const cs_real_t hybrid_blend_j,
1687  const cs_real_3_t diipf,
1688  const cs_real_3_t djjpf,
1689  const cs_real_3_t gradi,
1690  const cs_real_3_t gradj,
1691  const cs_real_3_t gradupi,
1692  const cs_real_3_t gradupj,
1693  const cs_real_t pi,
1694  const cs_real_t pj,
1695  cs_real_t *pif,
1696  cs_real_t *pjf,
1697  cs_real_t *pip,
1698  cs_real_t *pjp)
1699 {
1700  cs_real_t recoi, recoj;
1701 
1702  cs_i_compute_quantities(bldfrp,
1703  diipf,
1704  djjpf,
1705  gradi,
1706  gradj,
1707  pi,
1708  pj,
1709  &recoi,
1710  &recoj,
1711  pip,
1712  pjp);
1713 
1714 
1715  if (ischcp == 1) {
1716 
1717  /* Centered
1718  --------*/
1719 
1720  cs_centered_f_val(weight,
1721  *pip,
1722  *pjp,
1723  pif);
1724  cs_centered_f_val(weight,
1725  *pip,
1726  *pjp,
1727  pjf);
1728 
1729  } else if (ischcp == 0) {
1730 
1731  /* Legacy SOLU
1732  -----------*/
1733 
1734  cs_solu_f_val(cell_ceni,
1735  i_face_cog,
1736  gradi,
1737  pi,
1738  pif);
1739  cs_solu_f_val(cell_cenj,
1740  i_face_cog,
1741  gradj,
1742  pj,
1743  pjf);
1744 
1745  } else if (ischcp == 3) {
1746 
1747  /* Centered
1748  --------*/
1749 
1750  cs_centered_f_val(weight,
1751  *pip,
1752  *pjp,
1753  pif);
1754  cs_centered_f_val(weight,
1755  *pip,
1756  *pjp,
1757  pjf);
1758 
1759  /* Legacy SOLU
1760  -----------*/
1761  cs_real_t pif_up, pjf_up;
1762  cs_real_t hybrid_blend_interp;
1763 
1764  cs_solu_f_val(cell_ceni,
1765  i_face_cog,
1766  gradi,
1767  pi,
1768  &pif_up);
1769  cs_solu_f_val(cell_cenj,
1770  i_face_cog,
1771  gradj,
1772  pj,
1773  &pjf_up);
1774 
1775  hybrid_blend_interp = fmin(hybrid_blend_i,hybrid_blend_j);
1776  *pif = hybrid_blend_interp*(*pif) + (1. - hybrid_blend_interp)*pif_up;
1777  *pjf = hybrid_blend_interp*(*pjf) + (1. - hybrid_blend_interp)*pjf_up;
1778 
1779  } else {
1780 
1781  /* SOLU
1782  ----*/
1783 
1784  cs_solu_f_val(cell_ceni,
1785  i_face_cog,
1786  gradupi,
1787  pi,
1788  pif);
1789  cs_solu_f_val(cell_cenj,
1790  i_face_cog,
1791  gradupj,
1792  pj,
1793  pjf);
1794 
1795  }
1796 
1797  /* Blending
1798  --------*/
1799 
1800  cs_blend_f_val(blencp,
1801  pi,
1802  pif);
1803  cs_blend_f_val(blencp,
1804  pj,
1805  pjf);
1806 }
1807 
1808 /*----------------------------------------------------------------------------*/
1837 /*----------------------------------------------------------------------------*/
1838 
1839 template <cs_lnum_t stride>
1840 CS_F_HOST_DEVICE inline static void
1842  int ischcp,
1843  double blencp,
1844  cs_real_t weight,
1845  const cs_real_t cell_ceni[3],
1846  const cs_real_t cell_cenj[3],
1847  const cs_real_t i_face_cog[3],
1848  const cs_real_t hybrid_blend_i,
1849  const cs_real_t hybrid_blend_j,
1850  const cs_real_t diipf[3],
1851  const cs_real_t djjpf[3],
1852  const cs_real_t gradi[stride][3],
1853  const cs_real_t gradj[stride][3],
1854  const cs_real_t pi[stride],
1855  const cs_real_t pj[stride],
1856  cs_real_t pif[stride],
1857  cs_real_t pjf[stride],
1858  cs_real_t pip[stride],
1859  cs_real_t pjp[stride])
1860 
1861 {
1862  cs_real_t recoi[stride], recoj[stride];
1863 
1864  cs_i_compute_quantities_strided<stride>(bldfrp,
1865  diipf,
1866  djjpf,
1867  gradi,
1868  gradj,
1869  pi,
1870  pj,
1871  recoi,
1872  recoj,
1873  pip,
1874  pjp);
1875 
1876  if (ischcp == 1) {
1877 
1878  /* Centered
1879  --------*/
1880 
1881  cs_centered_f_val_strided<stride>(weight, pip, pjp, pif);
1882  cs_centered_f_val_strided<stride>(weight, pip, pjp, pjf);
1883 
1884  }
1885  else if (ischcp == 3) {
1886 
1887  /* Centered
1888  --------*/
1889 
1890  cs_centered_f_val_strided<stride>(weight, pip, pjp, pif);
1891  cs_centered_f_val_strided<stride>(weight, pip, pjp, pjf);
1892 
1893  /* SOLU
1894  -----*/
1895  cs_real_t pif_up[stride], pjf_up[stride];
1896  cs_real_t hybrid_blend_interp;
1897 
1898  cs_solu_f_val_strided<stride>(cell_ceni,
1899  i_face_cog,
1900  gradi,
1901  pi,
1902  pif_up);
1903  cs_solu_f_val_strided<stride>(cell_cenj,
1904  i_face_cog,
1905  gradj,
1906  pj,
1907  pjf_up);
1908 
1909  hybrid_blend_interp = fmin(hybrid_blend_i, hybrid_blend_j);
1910  for (int isou = 0; isou < stride; isou++) {
1911  pif[isou] = hybrid_blend_interp *pif[isou]
1912  + (1. - hybrid_blend_interp)*pif_up[isou];
1913  pjf[isou] = hybrid_blend_interp *pjf[isou]
1914  + (1. - hybrid_blend_interp)*pjf_up[isou];
1915  }
1916  }
1917  else {
1918 
1919  /* Second order
1920  ------------*/
1921 
1922  cs_solu_f_val_strided<stride>(cell_ceni,
1923  i_face_cog,
1924  gradi,
1925  pi,
1926  pif);
1927  cs_solu_f_val_strided<stride>(cell_cenj,
1928  i_face_cog,
1929  gradj,
1930  pj,
1931  pjf);
1932 
1933  }
1934 
1935  /* Blending
1936  --------*/
1937 
1938  cs_blend_f_val_strided<stride>(blencp, pi, pif);
1939  cs_blend_f_val_strided<stride>(blencp, pj, pjf);
1940 }
1941 
1942 /*----------------------------------------------------------------------------*/
1985 /*----------------------------------------------------------------------------*/
1986 
1987 CS_F_HOST_DEVICE inline static void
1988 cs_i_cd_steady_slope_test(bool *upwind_switch,
1989  const int iconvp,
1990  const cs_real_t bldfrp,
1991  const int ischcp,
1992  const double relaxp,
1993  const double blencp,
1994  const double blend_st,
1995  const cs_real_t weight,
1996  const cs_real_t i_dist,
1997  const cs_real_t cell_ceni[3],
1998  const cs_real_t cell_cenj[3],
1999  const cs_real_t i_face_u_normal[3],
2000  const cs_real_t i_face_cog[3],
2001  const cs_real_t diipf[3],
2002  const cs_real_t djjpf[3],
2003  const cs_real_t i_massflux,
2004  const cs_real_t gradi[3],
2005  const cs_real_t gradj[3],
2006  const cs_real_t gradupi[3],
2007  const cs_real_t gradupj[3],
2008  const cs_real_t gradsti[3],
2009  const cs_real_t gradstj[3],
2010  const cs_real_t pi,
2011  const cs_real_t pj,
2012  const cs_real_t pia,
2013  const cs_real_t pja,
2014  cs_real_t *pifri,
2015  cs_real_t *pifrj,
2016  cs_real_t *pjfri,
2017  cs_real_t *pjfrj,
2018  cs_real_t *pip,
2019  cs_real_t *pjp,
2020  cs_real_t *pipr,
2021  cs_real_t *pjpr)
2022 {
2023  cs_real_t pir, pjr;
2024  cs_real_t recoi, recoj;
2025  cs_real_t testij, tesqck;
2026 
2027  *upwind_switch = false;
2028 
2029  cs_i_compute_quantities(bldfrp,
2030  diipf,
2031  djjpf,
2032  gradi,
2033  gradj,
2034  pi,
2035  pj,
2036  &recoi,
2037  &recoj,
2038  pip,
2039  pjp);
2040 
2041  cs_i_relax_c_val(relaxp,
2042  pia,
2043  pja,
2044  recoi,
2045  recoj,
2046  pi,
2047  pj,
2048  &pir,
2049  &pjr,
2050  pipr,
2051  pjpr);
2052 
2053  /* Convection slope test is needed only when iconv >0 */
2054  if (iconvp > 0) {
2055  cs_slope_test(pi,
2056  pj,
2057  i_dist,
2058  i_face_u_normal,
2059  gradi,
2060  gradj,
2061  gradsti,
2062  gradstj,
2063  i_massflux,
2064  &testij,
2065  &tesqck);
2066 
2067  if (ischcp==1) {
2068 
2069  /* Centered
2070  --------*/
2071 
2072  cs_centered_f_val(weight,
2073  *pip,
2074  *pjpr,
2075  pifrj);
2076  cs_centered_f_val(weight,
2077  *pipr,
2078  *pjp,
2079  pifri);
2080  cs_centered_f_val(weight,
2081  *pipr,
2082  *pjp,
2083  pjfri);
2084  cs_centered_f_val(weight,
2085  *pip,
2086  *pjpr,
2087  pjfrj);
2088 
2089  } else if (ischcp == 0) {
2090 
2091  /* Second order
2092  ------------*/
2093 
2094  cs_solu_f_val(cell_ceni,
2095  i_face_cog,
2096  gradi,
2097  pi,
2098  pifrj);
2099  cs_solu_f_val(cell_ceni,
2100  i_face_cog,
2101  gradi,
2102  pir,
2103  pifri);
2104  cs_solu_f_val(cell_cenj,
2105  i_face_cog,
2106  gradj,
2107  pj,
2108  pjfri);
2109  cs_solu_f_val(cell_cenj,
2110  i_face_cog,
2111  gradj,
2112  pjr,
2113  pjfrj);
2114 
2115  } else {
2116 
2117  /* SOLU
2118  -----*/
2119 
2120  cs_solu_f_val(cell_ceni,
2121  i_face_cog,
2122  gradupi,
2123  pi,
2124  pifrj);
2125  cs_solu_f_val(cell_ceni,
2126  i_face_cog,
2127  gradupi,
2128  pir,
2129  pifri);
2130  cs_solu_f_val(cell_cenj,
2131  i_face_cog,
2132  gradupj,
2133  pj,
2134  pjfri);
2135  cs_solu_f_val(cell_cenj,
2136  i_face_cog,
2137  gradupj,
2138  pjr,
2139  pjfrj);
2140  }
2141 
2142 
2143  /* Slope test: percentage of upwind
2144  ----------------------------------*/
2145 
2146  if (tesqck <= 0. || testij <= 0.) {
2147 
2148  cs_blend_f_val(blend_st,
2149  pi,
2150  pifrj);
2151  cs_blend_f_val(blend_st,
2152  pir,
2153  pifri);
2154  cs_blend_f_val(blend_st,
2155  pj,
2156  pjfri);
2157  cs_blend_f_val(blend_st,
2158  pjr,
2159  pjfrj);
2160 
2161  *upwind_switch = true;
2162 
2163  }
2164 
2165 
2166  /* Blending
2167  --------*/
2168 
2169  cs_blend_f_val(blencp,
2170  pi,
2171  pifrj);
2172  cs_blend_f_val(blencp,
2173  pir,
2174  pifri);
2175  cs_blend_f_val(blencp,
2176  pj,
2177  pjfri);
2178  cs_blend_f_val(blencp,
2179  pjr,
2180  pjfrj);
2181 
2182  /* If iconv=0 p*fr* are useless */
2183  }
2184  else {
2185  cs_upwind_f_val(pi, pifrj);
2186  cs_upwind_f_val(pir, pifri);
2187  cs_upwind_f_val(pj, pjfri);
2188  cs_upwind_f_val(pjr, pjfrj);
2189  }
2190 
2191 }
2192 
2193 /*----------------------------------------------------------------------------*/
2237 /*----------------------------------------------------------------------------*/
2238 
2239 template <cs_lnum_t stride>
2240 CS_F_HOST_DEVICE inline static void
2242  int iconvp,
2243  cs_real_t bldfrp,
2244  int ischcp,
2245  double relaxp,
2246  double blencp,
2247  double blend_st,
2248  cs_real_t weight,
2249  cs_real_t i_dist,
2250  const cs_real_t cell_ceni[3],
2251  const cs_real_t cell_cenj[3],
2252  const cs_real_t i_face_u_normal[3],
2253  const cs_real_t i_face_cog[3],
2254  const cs_real_t diipf[3],
2255  const cs_real_t djjpf[3],
2256  cs_real_t i_massflux,
2257  const cs_real_t gradi[stride][3],
2258  const cs_real_t gradj[stride][3],
2259  const cs_real_t grdpai[stride][3],
2260  const cs_real_t grdpaj[stride][3],
2261  const cs_real_t pi[stride],
2262  const cs_real_t pj[stride],
2263  const cs_real_t pia[stride],
2264  const cs_real_t pja[stride],
2265  cs_real_t pifri[stride],
2266  cs_real_t pifrj[stride],
2267  cs_real_t pjfri[stride],
2268  cs_real_t pjfrj[stride],
2269  cs_real_t pip[stride],
2270  cs_real_t pjp[stride],
2271  cs_real_t pipr[stride],
2272  cs_real_t pjpr[stride])
2273 {
2274  cs_real_t pir[stride], pjr[stride];
2275  cs_real_t recoi[stride], recoj[stride];
2276  cs_real_t testij, tesqck;
2277  int isou;
2278 
2279  cs_i_compute_quantities_strided<stride>(bldfrp,
2280  diipf,
2281  djjpf,
2282  gradi,
2283  gradj,
2284  pi,
2285  pj,
2286  recoi,
2287  recoj,
2288  pip,
2289  pjp);
2290 
2291  cs_i_relax_c_val_strided<stride>(relaxp,
2292  pia,
2293  pja,
2294  recoi,
2295  recoj,
2296  pi,
2297  pj,
2298  pir,
2299  pjr,
2300  pipr,
2301  pjpr);
2302 
2303  /* Convection slope test is needed only when iconv >0 */
2304  if (iconvp > 0) {
2305  cs_slope_test_strided<stride>(pi,
2306  pj,
2307  i_dist,
2308  i_face_u_normal,
2309  gradi,
2310  gradj,
2311  grdpai,
2312  grdpaj,
2313  i_massflux,
2314  &testij,
2315  &tesqck);
2316 
2317  for (isou = 0; isou < stride; isou++) {
2318  if (ischcp==1) {
2319 
2320  /* Centered
2321  --------*/
2322 
2323  cs_centered_f_val(weight, pip[isou], pjpr[isou], &pifrj[isou]);
2324  cs_centered_f_val(weight, pipr[isou], pjp[isou], &pifri[isou]);
2325  cs_centered_f_val(weight, pipr[isou], pjp[isou], &pjfri[isou]);
2326  cs_centered_f_val(weight, pip[isou], pjpr[isou], &pjfrj[isou]);
2327 
2328  }
2329  else {
2330 
2331  /* Second order
2332  ------------*/
2333 
2334  cs_solu_f_val(cell_ceni, i_face_cog, gradi[isou], pi[isou],
2335  &pifrj[isou]);
2336  cs_solu_f_val(cell_ceni, i_face_cog, gradi[isou], pir[isou],
2337  &pifri[isou]);
2338  cs_solu_f_val(cell_cenj, i_face_cog, gradj[isou], pj[isou],
2339  &pjfri[isou]);
2340  cs_solu_f_val(cell_cenj, i_face_cog, gradj[isou], pjr[isou],
2341  &pjfrj[isou]);
2342 
2343  }
2344 
2345  }
2346 
2347  /* Slope test: percentage of upwind
2348  ----------------------------------*/
2349 
2350  if (tesqck <= 0. || testij <= 0.) {
2351 
2352  cs_blend_f_val_strided<stride>(blend_st, pi, pifrj);
2353  cs_blend_f_val_strided<stride>(blend_st, pir, pifri);
2354  cs_blend_f_val_strided<stride>(blend_st, pj, pjfri);
2355  cs_blend_f_val_strided<stride>(blend_st, pjr, pjfrj);
2356 
2357  *upwind_switch = true;
2358 
2359  }
2360 
2361  /* Blending
2362  --------*/
2363 
2364  cs_blend_f_val_strided<stride>(blencp, pi, pifrj);
2365  cs_blend_f_val_strided<stride>(blencp, pir, pifri);
2366  cs_blend_f_val_strided<stride>(blencp, pj, pjfri);
2367  cs_blend_f_val_strided<stride>(blencp, pjr, pjfrj);
2368 
2369  /* If iconv=0 p*fr* are useless */
2370  }
2371  else {
2372  for (isou = 0; isou < stride; isou++) {
2373  cs_upwind_f_val(pi[isou], &pifrj[isou]);
2374  cs_upwind_f_val(pir[isou], &pifri[isou]);
2375  cs_upwind_f_val(pj[isou], &pjfri[isou]);
2376  cs_upwind_f_val(pjr[isou], &pjfrj[isou]);
2377  }
2378  }
2379 }
2380 
2381 /*----------------------------------------------------------------------------*/
2417 /*----------------------------------------------------------------------------*/
2418 
2419 CS_F_HOST_DEVICE inline static void
2420 cs_i_cd_unsteady_slope_test(bool *upwind_switch,
2421  const int iconvp,
2422  const cs_real_t bldfrp,
2423  const int ischcp,
2424  const double blencp,
2425  const double blend_st,
2426  const cs_real_t weight,
2427  const cs_real_t i_dist,
2428  const cs_real_3_t cell_ceni,
2429  const cs_real_3_t cell_cenj,
2430  const cs_real_3_t i_face_u_normal,
2431  const cs_real_3_t i_face_cog,
2432  const cs_real_3_t diipf,
2433  const cs_real_3_t djjpf,
2434  const cs_real_t i_massflux,
2435  const cs_real_3_t gradi,
2436  const cs_real_3_t gradj,
2437  const cs_real_3_t gradupi,
2438  const cs_real_3_t gradupj,
2439  const cs_real_3_t gradsti,
2440  const cs_real_3_t gradstj,
2441  const cs_real_t pi,
2442  const cs_real_t pj,
2443  cs_real_t *pif,
2444  cs_real_t *pjf,
2445  cs_real_t *pip,
2446  cs_real_t *pjp)
2447 {
2448  CS_UNUSED(blend_st);
2449 
2450  cs_real_t recoi, recoj;
2451  cs_real_t testij, tesqck;
2452 
2453  *upwind_switch = false;
2454 
2455  cs_i_compute_quantities(bldfrp,
2456  diipf,
2457  djjpf,
2458  gradi,
2459  gradj,
2460  pi,
2461  pj,
2462  &recoi,
2463  &recoj,
2464  pip,
2465  pjp);
2466 
2467  /* Convection slope test is needed only when iconv >0 */
2468  if (iconvp > 0) {
2469  cs_slope_test(pi,
2470  pj,
2471  i_dist,
2472  i_face_u_normal,
2473  gradi,
2474  gradj,
2475  gradsti,
2476  gradstj,
2477  i_massflux,
2478  &testij,
2479  &tesqck);
2480 
2481  if (ischcp==1) {
2482 
2483  /* Centered
2484  --------*/
2485 
2486  cs_centered_f_val(weight,
2487  *pip,
2488  *pjp,
2489  pif);
2490  cs_centered_f_val(weight,
2491  *pip,
2492  *pjp,
2493  pjf);
2494 
2495  } else if (ischcp == 0) {
2496 
2497  /* Original SOLU
2498  --------------*/
2499 
2500  cs_solu_f_val(cell_ceni,
2501  i_face_cog,
2502  gradi,
2503  pi,
2504  pif);
2505  cs_solu_f_val(cell_cenj,
2506  i_face_cog,
2507  gradj,
2508  pj,
2509  pjf);
2510 
2511  } else {
2512 
2513  /* SOLU
2514  -----*/
2515 
2516  cs_solu_f_val(cell_ceni,
2517  i_face_cog,
2518  gradupi,
2519  pi,
2520  pif);
2521  cs_solu_f_val(cell_cenj,
2522  i_face_cog,
2523  gradupj,
2524  pj,
2525  pjf);
2526 
2527  }
2528 
2529  /* Slope test: percentage of upwind
2530  -------------------------------- */
2531 
2532  if (tesqck<=0. || testij<=0.) {
2533 
2534  cs_blend_f_val(blend_st,
2535  pi,
2536  pif);
2537  cs_blend_f_val(blend_st,
2538  pj,
2539  pjf);
2540 
2541  *upwind_switch = true;
2542 
2543  }
2544 
2545  /* Blending
2546  --------*/
2547 
2548  cs_blend_f_val(blencp,
2549  pi,
2550  pif);
2551  cs_blend_f_val(blencp,
2552  pj,
2553  pjf);
2554 
2555  /* If iconv=0 p*f are useless */
2556  } else {
2558  pif);
2559  cs_upwind_f_val(pj,
2560  pjf);
2561  }
2562 
2563 }
2564 
2565 /*----------------------------------------------------------------------------*/
2576 /*----------------------------------------------------------------------------*/
2577 
2578 CS_F_HOST_DEVICE inline static void
2580  const cs_lnum_t jj,
2581  const cs_real_t i_massflux,
2582  cs_lnum_t *ic,
2583  cs_lnum_t *id)
2584 {
2585  if (i_massflux >= 0.) {
2586  *ic = ii;
2587  *id = jj;
2588  } else {
2589  *ic = jj;
2590  *id = ii;
2591  }
2592 }
2593 
2594 /*----------------------------------------------------------------------------*/
2615 /*----------------------------------------------------------------------------*/
2616 
2617 CS_F_HOST_DEVICE inline static void
2619  const double beta,
2620  const cs_real_3_t cell_cen_c,
2621  const cs_real_3_t cell_cen_d,
2622  const cs_real_3_t i_face_u_normal,
2623  const cs_real_3_t i_face_cog,
2624  const cs_real_3_t gradv_c,
2625  const cs_real_t p_c,
2626  const cs_real_t p_d,
2627  const cs_real_t local_max_c,
2628  const cs_real_t local_min_c,
2629  const cs_real_t courant_c,
2630  cs_real_t *pif,
2631  cs_real_t *pjf)
2632 {
2633  /* distance between face center and central cell center */
2634  /* Distance between face center and central cell center */
2635  cs_real_t dist_fc = cs_math_3_distance(cell_cen_c, i_face_cog);
2636 
2637  /* Unit vector and distance between central and downwind cells centers */
2638  cs_real_t diff[3] = {cell_cen_d[0] - cell_cen_c[0],
2639  cell_cen_d[1] - cell_cen_c[1],
2640  cell_cen_d[2] - cell_cen_c[2]};
2641 
2642  cs_real_t dist_dc = cs_math_3_norm(diff);
2643  cs_real_t invl = 1./dist_dc;
2644 
2645  cs_real_t ndc[3] = {invl*diff[0],
2646  invl*diff[1],
2647  invl*diff[2]};
2648 
2649  /* Place the upwind point on the line that joins
2650  the two cells on the upwind side and the same
2651  distance as that between the two cells */
2652  const cs_real_t dist_cu = dist_dc;
2653  const cs_real_t dist_du = dist_dc + dist_cu;
2654 
2655  /* Compute the property on the upwind assuming a parabolic
2656  variation of the property between the two cells */
2657  const cs_real_t gradc = cs_math_3_dot_product(gradv_c, ndc);
2658 
2659  const cs_real_t grad2c = ((p_d - p_c)/dist_dc - gradc)/dist_dc;
2660 
2661  cs_real_t p_u = p_c + (grad2c*dist_cu - gradc)*dist_cu;
2662  p_u = cs_math_fmax(cs_math_fmin(p_u, local_max_c), local_min_c);
2663 
2664  /* Compute the normalised distances */
2665  const cs_real_t nvf_r_f = (dist_fc+dist_cu)/dist_du;
2666  const cs_real_t nvf_r_c = dist_cu/dist_du;
2667 
2668  /* Check for the bounds of the NVD diagram and compute the face
2669  property according to the selected NVD scheme */
2670  const cs_real_t _small = cs_math_epzero
2671  * (CS_ABS(p_u)+CS_ABS(p_c)+CS_ABS(p_d));
2672 
2673  if (CS_ABS(p_d-p_u) <= _small) {
2674  *pif = p_c;
2675  *pjf = p_c;
2676  } else {
2677  const cs_real_t nvf_p_c = (p_c - p_u)/(p_d - p_u);
2678 
2679  if (nvf_p_c <= 0. || nvf_p_c >= 1.) {
2680  *pif = p_c;
2681  *pjf = p_c;
2682  } else {
2683  cs_real_t nvf_p_f;
2684 
2685  /* Highly compressive NVD scheme for VOF */
2686  if (limiter >= CS_NVD_VOF_HRIC) {
2687  nvf_p_f = cs_nvd_vof_scheme_scalar(limiter,
2688  i_face_u_normal,
2689  nvf_p_c,
2690  nvf_r_f,
2691  nvf_r_c,
2692  gradv_c,
2693  courant_c);
2694  } else { /* Regular NVD scheme */
2695  nvf_p_f = cs_nvd_scheme_scalar(limiter,
2696  nvf_p_c,
2697  nvf_r_f,
2698  nvf_r_c);
2699  }
2700 
2701  *pif = p_u + nvf_p_f*(p_d - p_u);
2702  *pif = cs_math_fmax(cs_math_fmin(*pif, local_max_c), local_min_c);
2703 
2704  cs_blend_f_val(beta,
2705  p_c,
2706  pif);
2707 
2708  *pjf = *pif;
2709  }
2710  }
2711 }
2712 
2713 /*----------------------------------------------------------------------------*/
2750 /*----------------------------------------------------------------------------*/
2751 
2752 template <cs_lnum_t stride>
2753 CS_F_HOST_DEVICE inline static void
2755  int iconvp,
2756  cs_real_t bldfrp,
2757  int ischcp,
2758  double blencp,
2759  double blend_st,
2760  cs_real_t weight,
2761  cs_real_t i_dist,
2762  const cs_real_t cell_ceni[3],
2763  const cs_real_t cell_cenj[3],
2764  const cs_real_t i_face_u_normal[3],
2765  const cs_real_t i_face_cog[3],
2766  const cs_real_t diipf[3],
2767  const cs_real_t djjpf[3],
2768  cs_real_t i_massflux,
2769  const cs_real_t gradi[stride][3],
2770  const cs_real_t gradj[stride][3],
2771  const cs_real_t grdpai[stride][3],
2772  const cs_real_t grdpaj[stride][3],
2773  const cs_real_t pi[stride],
2774  const cs_real_t pj[stride],
2775  cs_real_t pif[stride],
2776  cs_real_t pjf[stride],
2777  cs_real_t pip[stride],
2778  cs_real_t pjp[stride])
2779 {
2780  cs_real_t recoi[stride], recoj[stride];
2781  cs_real_t testij, tesqck;
2782  int isou;
2783 
2784  cs_i_compute_quantities_strided<stride>(bldfrp,
2785  diipf,
2786  djjpf,
2787  gradi,
2788  gradj,
2789  pi,
2790  pj,
2791  recoi,
2792  recoj,
2793  pip,
2794  pjp);
2795 
2796  /* Convection slope test is needed only when iconv >0 */
2797  if (iconvp > 0) {
2798  cs_slope_test_strided<stride>(pi,
2799  pj,
2800  i_dist,
2801  i_face_u_normal,
2802  gradi,
2803  gradj,
2804  grdpai,
2805  grdpaj,
2806  i_massflux,
2807  &testij,
2808  &tesqck);
2809 
2810  for (isou = 0; isou < stride; isou++) {
2811 
2812  if (ischcp==1) {
2813 
2814  /* Centered
2815  --------*/
2816 
2817  cs_centered_f_val(weight,
2818  pip[isou],
2819  pjp[isou],
2820  &pif[isou]);
2821  cs_centered_f_val(weight,
2822  pip[isou],
2823  pjp[isou],
2824  &pjf[isou]);
2825 
2826  }
2827  else {
2828 
2829  /* Second order
2830  ------------*/
2831 
2832  cs_solu_f_val(cell_ceni,
2833  i_face_cog,
2834  gradi[isou],
2835  pi[isou],
2836  &pif[isou]);
2837  cs_solu_f_val(cell_cenj,
2838  i_face_cog,
2839  gradj[isou],
2840  pj[isou],
2841  &pjf[isou]);
2842  }
2843 
2844  }
2845 
2846  /* Slope test activated: percentage of upwind */
2847  if (tesqck <= 0. || testij <= 0.) {
2848 
2849  /* Upwind
2850  --------*/
2851 
2852  cs_blend_f_val_strided<stride>(blend_st, pi, pif);
2853  cs_blend_f_val_strided<stride>(blend_st, pj, pjf);
2854 
2855  *upwind_switch = true;
2856  }
2857 
2858  /* Blending
2859  --------*/
2860 
2861  cs_blend_f_val_strided<stride>(blencp, pi, pif);
2862  cs_blend_f_val_strided<stride>(blencp, pj, pjf);
2863 
2864  /* If iconv=0 p*fr* are useless */
2865  }
2866  else {
2867 
2868  for (isou = 0; isou < stride; isou++) {
2869  cs_upwind_f_val(pi[isou], &pif[isou]);
2870  cs_upwind_f_val(pj[isou], &pjf[isou]);
2871  }
2872  }
2873 }
2874 
2875 /*----------------------------------------------------------------------------*/
2884 /*----------------------------------------------------------------------------*/
2885 
2886 CS_F_HOST_DEVICE inline static void
2888  const cs_real_3_t gradi,
2889  const cs_real_t bldfrp,
2890  cs_real_t *recoi)
2891 {
2892  *recoi = bldfrp * ( gradi[0]*diipb[0]
2893  + gradi[1]*diipb[1]
2894  + gradi[2]*diipb[2]);
2895 }
2896 
2897 /*----------------------------------------------------------------------------*/
2909 /*----------------------------------------------------------------------------*/
2910 
2911 template <cs_lnum_t stride>
2912 CS_F_HOST_DEVICE inline static void
2914  const cs_real_t gradi[stride][3],
2915  const cs_real_t bldfrp,
2916  cs_real_t recoi[stride])
2917 {
2918  for (int isou = 0; isou < stride; isou++) {
2919  recoi[isou] = bldfrp * ( gradi[isou][0]*diipb[0]
2920  + gradi[isou][1]*diipb[1]
2921  + gradi[isou][2]*diipb[2]);
2922  }
2923 }
2924 
2925 /*----------------------------------------------------------------------------*/
2936 /*----------------------------------------------------------------------------*/
2937 
2938 CS_F_HOST_DEVICE inline static void
2939 cs_b_relax_c_val(const double relaxp,
2940  const cs_real_t pi,
2941  const cs_real_t pia,
2942  const cs_real_t recoi,
2943  cs_real_t *pir,
2944  cs_real_t *pipr)
2945 {
2946  *pir = pi/relaxp - (1.-relaxp)/relaxp*pia;
2947  *pipr = *pir + recoi;
2948 }
2949 
2950 /*----------------------------------------------------------------------------*/
2964 /*----------------------------------------------------------------------------*/
2965 
2966 template <cs_lnum_t stride>
2967 CS_F_HOST_DEVICE inline static void
2968 cs_b_relax_c_val_strided(const double relaxp,
2969  const cs_real_t pi[stride],
2970  const cs_real_t pia[stride],
2971  const cs_real_t recoi[stride],
2972  cs_real_t pir[stride],
2973  cs_real_t pipr[stride])
2974 {
2975  for (int isou = 0; isou < stride; isou++) {
2976  pir[isou] = pi[isou]/relaxp - (1.-relaxp)/relaxp*pia[isou];
2977  pipr[isou] = pir[isou] + recoi[isou];
2978  }
2979 }
2980 
2981 /*----------------------------------------------------------------------------*/
3005 /*----------------------------------------------------------------------------*/
3006 
3007 CS_F_HOST_DEVICE inline static void
3009  cs_real_t thetap,
3010  int imasac,
3011  int inc,
3012  int bc_type,
3013  int icvfli,
3014  cs_real_t pi,
3015  cs_real_t pir,
3016  cs_real_t pipr,
3017  cs_real_t coefap,
3018  cs_real_t coefbp,
3019  cs_real_t coface,
3020  cs_real_t cofbce,
3021  cs_real_t b_massflux,
3022  cs_real_t xcpp,
3023  cs_real_t *flux)
3024 {
3025  cs_real_t flui, fluj, pfac;
3026 
3027  /* Computed convective flux */
3028 
3029  if (icvfli == 0) {
3030 
3031  /* Remove decentering for coupled faces */
3032  if (bc_type == CS_COUPLED_FD) {
3033  flui = 0.0;
3034  fluj = b_massflux;
3035  } else {
3036  flui = 0.5*(b_massflux +fabs(b_massflux));
3037  fluj = 0.5*(b_massflux -fabs(b_massflux));
3038  }
3039 
3040  pfac = inc*coefap + coefbp*pipr;
3041  *flux += iconvp*xcpp*(thetap*(flui*pir + fluj*pfac) -imasac*( b_massflux*pi));
3042 
3043  /* Imposed convective flux */
3044 
3045  } else {
3046 
3047  pfac = inc*coface + cofbce*pipr;
3048  *flux += iconvp*xcpp*(-imasac*(b_massflux*pi) + thetap*(pfac));
3049 
3050  }
3051 }
3052 
3053 /*----------------------------------------------------------------------------*/
3078 /*----------------------------------------------------------------------------*/
3079 
3080 template <cs_lnum_t stride>
3081 CS_F_HOST_DEVICE inline static void
3083  cs_real_t thetap,
3084  int imasac,
3085  int inc,
3086  int bc_type,
3087  int icvfli,
3088  const cs_real_t pi[stride],
3089  const cs_real_t pir[stride],
3090  const cs_real_t pipr[stride],
3091  const cs_real_t coefap[stride],
3092  const cs_real_t coefbp[stride][stride],
3093  const cs_real_t coface[stride],
3094  const cs_real_t cofbce[stride][stride],
3095  cs_real_t b_massflux,
3096  cs_real_t pfac[stride],
3097  cs_real_t flux[stride])
3098 {
3099  cs_real_t flui, fluj;
3100 
3101  /* Computed convective flux */
3102 
3103  if (icvfli == 0) {
3104 
3105  /* Remove decentering for coupled faces */
3106  if (bc_type == CS_COUPLED_FD) {
3107  flui = 0.0;
3108  fluj = b_massflux;
3109  }
3110  else {
3111  flui = 0.5*(b_massflux +fabs(b_massflux));
3112  fluj = 0.5*(b_massflux -fabs(b_massflux));
3113  }
3114  for (int isou = 0; isou < stride; isou++) {
3115  pfac[isou] = inc*coefap[isou];
3116  for (int jsou = 0; jsou < stride; jsou++) {
3117  pfac[isou] += coefbp[isou][jsou]*pipr[jsou];
3118  }
3119  flux[isou] += iconvp*( thetap*(flui*pir[isou] + fluj*pfac[isou])
3120  - imasac*b_massflux*pi[isou]);
3121  }
3122 
3123  /* Imposed convective flux */
3124 
3125  }
3126  else {
3127 
3128  for (int isou = 0; isou < stride; isou++) {
3129  pfac[isou] = inc*coface[isou];
3130  for (int jsou = 0; jsou < stride; jsou++) {
3131  pfac[isou] += cofbce[isou][jsou]*pipr[jsou];
3132  }
3133  flux[isou] += iconvp*( thetap*pfac[isou]
3134  - imasac*b_massflux*pi[isou]);
3135  }
3136 
3137  }
3138 }
3139 
3140 /*----------------------------------------------------------------------------*/
3160 /*----------------------------------------------------------------------------*/
3161 
3162 CS_F_HOST_DEVICE inline static void
3163 cs_b_upwind_flux(const int iconvp,
3164  const cs_real_t thetap,
3165  const int imasac,
3166  const int inc,
3167  const int bc_type,
3168  const cs_real_t pi,
3169  const cs_real_t pir,
3170  const cs_real_t pipr,
3171  const cs_real_t coefap,
3172  const cs_real_t coefbp,
3173  const cs_real_t b_massflux,
3174  const cs_real_t xcpp,
3175  cs_real_t *flux)
3176 {
3177  cs_real_t flui, fluj, pfac;
3178 
3179  /* Remove decentering for coupled faces */
3180  if (bc_type == CS_COUPLED_FD) {
3181  flui = 0.0;
3182  fluj = b_massflux;
3183  } else {
3184  flui = 0.5*(b_massflux +fabs(b_massflux));
3185  fluj = 0.5*(b_massflux -fabs(b_massflux));
3186  }
3187 
3188  pfac = inc*coefap + coefbp*pipr;
3189  *flux += iconvp*xcpp*(thetap*(flui*pir + fluj*pfac) -imasac*( b_massflux*pi));
3190 }
3191 
3192 /*----------------------------------------------------------------------------*/
3216 /*----------------------------------------------------------------------------*/
3217 
3218 template <cs_lnum_t stride>
3219 CS_F_HOST_DEVICE inline static void
3221  cs_real_t thetap,
3222  int imasac,
3223  int inc,
3224  int bc_type,
3225  const cs_real_t pi[stride],
3226  const cs_real_t pir[stride],
3227  const cs_real_t pipr[stride],
3228  const cs_real_t coefa[stride],
3229  const cs_real_t coefb[stride][stride],
3230  cs_real_t b_massflux,
3231  cs_real_t pfac[stride],
3232  cs_real_t flux[stride])
3233 {
3234  cs_real_t flui, fluj;
3235 
3236  /* Remove decentering for coupled faces */
3237  if (bc_type == CS_COUPLED_FD) {
3238  flui = 0.0;
3239  fluj = b_massflux;
3240  } else {
3241  flui = 0.5*(b_massflux +fabs(b_massflux));
3242  fluj = 0.5*(b_massflux -fabs(b_massflux));
3243  }
3244  for (int isou = 0; isou < stride; isou++) {
3245  pfac[isou] = inc*coefa[isou];
3246  for (int jsou = 0; jsou < stride; jsou++) {
3247  pfac[isou] += coefb[isou][jsou]*pipr[jsou];
3248  }
3249  flux[isou] += iconvp*( thetap*(flui*pir[isou] + fluj*pfac[isou])
3250  - imasac*b_massflux*pi[isou]);
3251  }
3252 }
3253 
3254 /*----------------------------------------------------------------------------*/
3267 /*----------------------------------------------------------------------------*/
3268 
3269 CS_F_HOST_DEVICE inline static void
3270 cs_b_diff_flux(const int idiffp,
3271  const cs_real_t thetap,
3272  const int inc,
3273  const cs_real_t pipr,
3274  const cs_real_t cofafp,
3275  const cs_real_t cofbfp,
3276  const cs_real_t b_visc,
3277  cs_real_t *flux)
3278 {
3279  cs_real_t pfacd = inc*cofafp + cofbfp*pipr;
3280  *flux += idiffp*thetap*b_visc*pfacd;
3281 }
3282 
3283 /*----------------------------------------------------------------------------*/
3299 /*----------------------------------------------------------------------------*/
3300 
3301 template <cs_lnum_t stride>
3302 CS_F_HOST_DEVICE inline static void
3304  cs_real_t thetap,
3305  int inc,
3306  const cs_real_t pipr[stride],
3307  const cs_real_t cofaf[stride],
3308  const cs_real_t cofbf[stride][stride],
3309  cs_real_t b_visc,
3310  cs_real_t flux[stride])
3311 {
3312  cs_real_t pfacd;
3313  for (int isou = 0; isou < stride; isou++) {
3314  pfacd = inc*cofaf[isou];
3315  for (int jsou = 0; jsou < stride; jsou++) {
3316  pfacd += cofbf[isou][jsou]*pipr[jsou];
3317  }
3318  flux[isou] += idiffp*thetap*b_visc*pfacd;
3319  }
3320 }
3321 
3322 /*----------------------------------------------------------------------------*/
3336 /*----------------------------------------------------------------------------*/
3337 
3338 CS_F_HOST_DEVICE inline static void
3340  const double relaxp,
3341  const cs_real_3_t diipb,
3342  const cs_real_3_t gradi,
3343  const cs_real_t pi,
3344  const cs_real_t pia,
3345  cs_real_t *pir,
3346  cs_real_t *pipr)
3347 {
3348  cs_real_t recoi;
3349 
3351  gradi,
3352  bldfrp,
3353  &recoi);
3354 
3355  cs_b_relax_c_val(relaxp,
3356  pi,
3357  pia,
3358  recoi,
3359  pir,
3360  pipr);
3361 }
3362 
3363 /*----------------------------------------------------------------------------*/
3380 /*----------------------------------------------------------------------------*/
3381 
3382 template <cs_lnum_t stride>
3383 CS_F_HOST_DEVICE inline static void
3385  double relaxp,
3386  const cs_real_t diipb[3],
3387  const cs_real_t gradi[stride][3],
3388  const cs_real_t pi[stride],
3389  const cs_real_t pia[stride],
3390  cs_real_t pir[stride],
3391  cs_real_t pipr[stride])
3392 {
3393  cs_real_t recoi[stride];
3394 
3395  cs_b_compute_quantities_strided<stride>(diipb,
3396  gradi,
3397  bldfrp,
3398  recoi);
3399 
3400  cs_b_relax_c_val_strided<stride>(relaxp,
3401  pi,
3402  pia,
3403  recoi,
3404  pir,
3405  pipr);
3406 }
3407 
3408 /*----------------------------------------------------------------------------*/
3419 /*----------------------------------------------------------------------------*/
3420 
3421 CS_F_HOST_DEVICE inline static void
3423  const cs_real_3_t diipb,
3424  const cs_real_3_t gradi,
3425  const cs_real_t pi,
3426  cs_real_t *pip)
3427 {
3428  cs_real_t recoi;
3429 
3431  gradi,
3432  bldfrp,
3433  &recoi);
3434 
3435  *pip = pi + recoi;
3436 }
3437 
3438 /*----------------------------------------------------------------------------*/
3452 /*----------------------------------------------------------------------------*/
3453 
3454 template <cs_lnum_t stride>
3455 CS_F_HOST_DEVICE inline static void
3457  const cs_real_t diipb[3],
3458  const cs_real_t gradi[stride][3],
3459  const cs_real_t pi[stride],
3460  cs_real_t pip[stride])
3461 {
3462  cs_real_t recoi[stride];
3463 
3464  cs_b_compute_quantities_strided<stride>(diipb,
3465  gradi,
3466  bldfrp,
3467  recoi);
3468 
3469  for (int isou = 0; isou< stride; isou++)
3470  pip[isou] = pi[isou] + recoi[isou];
3471 }
3472 
3473 /*----------------------------------------------------------------------------*/
3484 /*----------------------------------------------------------------------------*/
3485 
3486 CS_F_HOST_DEVICE inline static void
3488  cs_real_t pi,
3489  cs_real_t pj,
3490  cs_real_t b_visc,
3491  cs_real_t *fluxi)
3492 {
3493  *fluxi += idiffp*b_visc*(pi - pj);
3494 }
3495 
3496 /*----------------------------------------------------------------------------*/
3510 /*----------------------------------------------------------------------------*/
3511 
3512 template <cs_lnum_t stride>
3513 CS_F_HOST_DEVICE inline static void
3515  const cs_real_t pi[stride],
3516  const cs_real_t pj[stride],
3517  cs_real_t b_visc,
3518  cs_real_t fluxi[stride])
3519 {
3520  for (int k = 0; k < stride; k++)
3521  fluxi[k] += idiffp*b_visc*(pi[k] - pj[k]);
3522 }
3523 
3524 /*----------------------------------------------------------------------------*/
3525 
3526 #endif /* __CS_CONVECTION_DIFFUSION_PRIV_H__ */
cs_nvd_type_t
Definition: cs_convection_diffusion.h:123
@ CS_NVD_SUPERBEE
Definition: cs_convection_diffusion.h:128
@ CS_NVD_SMART
Definition: cs_convection_diffusion.h:126
@ CS_NVD_CUBISTA
Definition: cs_convection_diffusion.h:127
@ CS_NVD_STOIC
Definition: cs_convection_diffusion.h:132
@ CS_NVD_WASEB
Definition: cs_convection_diffusion.h:134
@ CS_NVD_CLAM
Definition: cs_convection_diffusion.h:131
@ CS_NVD_OSHER
Definition: cs_convection_diffusion.h:133
@ CS_NVD_VOF_CICSAM
Definition: cs_convection_diffusion.h:136
@ CS_NVD_GAMMA
Definition: cs_convection_diffusion.h:125
@ CS_NVD_MINMOD
Definition: cs_convection_diffusion.h:130
@ CS_NVD_VOF_HRIC
Definition: cs_convection_diffusion.h:135
@ CS_NVD_MUSCL
Definition: cs_convection_diffusion.h:129
static CS_F_HOST_DEVICE void cs_b_diff_flux_coupling(int idiffp, cs_real_t pi, cs_real_t pj, cs_real_t b_visc, cs_real_t *fluxi)
Add diffusive flux to flux at an internal coupling face.
Definition: cs_convection_diffusion_priv.h:3487
static CS_F_HOST_DEVICE void cs_i_diff_flux_strided(int idiffp, cs_real_t thetap, const cs_real_t pip[stride], const cs_real_t pjp[stride], const cs_real_t pipr[stride], const cs_real_t pjpr[stride], const cs_real_t i_visc, cs_real_t fluxi[stride], cs_real_t fluxj[stride])
Add diffusive fluxes to fluxes at face ij.
Definition: cs_convection_diffusion_priv.h:1037
static CS_F_HOST_DEVICE void cs_centered_f_val_strided(double pnd, const cs_real_t pip[stride], const cs_real_t pjp[stride], cs_real_t pf[stride])
Prepare value at face ij by using a centered scheme.
Definition: cs_convection_diffusion_priv.h:772
static CS_F_HOST_DEVICE void cs_i_cd_steady_slope_test(bool *upwind_switch, const int iconvp, const cs_real_t bldfrp, const int ischcp, const double relaxp, const double blencp, const double blend_st, const cs_real_t weight, const cs_real_t i_dist, const cs_real_t cell_ceni[3], const cs_real_t cell_cenj[3], const cs_real_t i_face_u_normal[3], const cs_real_t i_face_cog[3], const cs_real_t diipf[3], const cs_real_t djjpf[3], const cs_real_t i_massflux, const cs_real_t gradi[3], const cs_real_t gradj[3], const cs_real_t gradupi[3], const cs_real_t gradupj[3], const cs_real_t gradsti[3], const cs_real_t gradstj[3], const cs_real_t pi, const cs_real_t pj, const cs_real_t pia, const cs_real_t pja, cs_real_t *pifri, cs_real_t *pifrj, cs_real_t *pjfri, cs_real_t *pjfrj, cs_real_t *pip, cs_real_t *pjp, cs_real_t *pipr, cs_real_t *pjpr)
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion_priv.h:1988
static CS_F_HOST_DEVICE void cs_i_cd_unsteady_strided(cs_real_t bldfrp, int ischcp, double blencp, cs_real_t weight, const cs_real_t cell_ceni[3], const cs_real_t cell_cenj[3], const cs_real_t i_face_cog[3], const cs_real_t hybrid_blend_i, const cs_real_t hybrid_blend_j, const cs_real_t diipf[3], const cs_real_t djjpf[3], const cs_real_t gradi[stride][3], const cs_real_t gradj[stride][3], const cs_real_t pi[stride], const cs_real_t pj[stride], cs_real_t pif[stride], cs_real_t pjf[stride], cs_real_t pip[stride], cs_real_t pjp[stride])
Handle preparation of internal face values for the fluxes computation in case of an unsteady algorith...
Definition: cs_convection_diffusion_priv.h:1841
static CS_F_HOST_DEVICE void cs_b_diff_flux_strided(int idiffp, cs_real_t thetap, int inc, const cs_real_t pipr[stride], const cs_real_t cofaf[stride], const cs_real_t cofbf[stride][stride], cs_real_t b_visc, cs_real_t flux[stride])
Add diffusive flux to flux at boundary face.
Definition: cs_convection_diffusion_priv.h:3303
static CS_F_HOST_DEVICE void cs_b_cd_steady(const cs_real_t bldfrp, const double relaxp, const cs_real_3_t diipb, const cs_real_3_t gradi, const cs_real_t pi, const cs_real_t pia, cs_real_t *pir, cs_real_t *pipr)
Handle preparation of boundary face values for the flux computation in case of a steady algorithm.
Definition: cs_convection_diffusion_priv.h:3339
static CS_F_HOST_DEVICE void cs_i_cd_unsteady(const cs_real_t bldfrp, const int ischcp, const double blencp, const cs_real_t weight, const cs_real_3_t cell_ceni, const cs_real_3_t cell_cenj, const cs_real_3_t i_face_cog, const cs_real_t hybrid_blend_i, const cs_real_t hybrid_blend_j, const cs_real_3_t diipf, const cs_real_3_t djjpf, const cs_real_3_t gradi, const cs_real_3_t gradj, const cs_real_3_t gradupi, const cs_real_3_t gradupj, const cs_real_t pi, const cs_real_t pj, cs_real_t *pif, cs_real_t *pjf, cs_real_t *pip, cs_real_t *pjp)
Handle preparation of internal face values for the fluxes computation in case of a unsteady algorithm...
Definition: cs_convection_diffusion_priv.h:1678
static CS_F_HOST_DEVICE void cs_i_relax_c_val_strided(cs_real_t relaxp, const cs_real_t pia[stride], const cs_real_t pja[stride], const cs_real_t recoi[stride], const cs_real_t recoj[stride], const cs_real_t pi[stride], const cs_real_t pj[stride], cs_real_t pir[stride], cs_real_t pjr[stride], cs_real_t pipr[stride], cs_real_t pjpr[stride])
Compute relaxed values at cell i and j.
Definition: cs_convection_diffusion_priv.h:678
static CS_F_HOST_DEVICE void cs_i_cd_steady_upwind_strided(const cs_real_t bldfrp, const cs_real_t relaxp, const cs_real_t diipf[3], const cs_real_t djjpf[3], const cs_real_t gradi[stride][3], const cs_real_t gradj[stride][3], const cs_real_t pi[stride], const cs_real_t pj[stride], const cs_real_t pia[stride], const cs_real_t pja[stride], cs_real_t pifri[stride], cs_real_t pifrj[stride], cs_real_t pjfri[stride], cs_real_t pjfrj[stride], cs_real_t pip[stride], cs_real_t pjp[stride], cs_real_t pipr[stride], cs_real_t pjpr[stride])
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion_priv.h:1167
static CS_F_HOST_DEVICE void cs_i_cd_unsteady_slope_test(bool *upwind_switch, const int iconvp, const cs_real_t bldfrp, const int ischcp, const double blencp, const double blend_st, const cs_real_t weight, const cs_real_t i_dist, const cs_real_3_t cell_ceni, const cs_real_3_t cell_cenj, const cs_real_3_t i_face_u_normal, const cs_real_3_t i_face_cog, const cs_real_3_t diipf, const cs_real_3_t djjpf, const cs_real_t i_massflux, const cs_real_3_t gradi, const cs_real_3_t gradj, const cs_real_3_t gradupi, const cs_real_3_t gradupj, const cs_real_3_t gradsti, const cs_real_3_t gradstj, const cs_real_t pi, const cs_real_t pj, cs_real_t *pif, cs_real_t *pjf, cs_real_t *pip, cs_real_t *pjp)
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion_priv.h:2420
static CS_F_HOST_DEVICE void cs_b_relax_c_val_strided(const double relaxp, const cs_real_t pi[stride], const cs_real_t pia[stride], const cs_real_t recoi[stride], cs_real_t pir[stride], cs_real_t pipr[stride])
Compute relaxed values at boundary cell i.
Definition: cs_convection_diffusion_priv.h:2968
static CS_F_HOST_DEVICE void cs_i_compute_quantities_strided(const cs_real_t bldfrp, const cs_real_t diipf[3], const cs_real_t djjpf[3], const cs_real_t gradi[stride][3], const cs_real_t gradj[stride][3], const cs_real_t pi[stride], const cs_real_t pj[stride], cs_real_t recoi[stride], cs_real_t recoj[stride], cs_real_t pip[stride], cs_real_t pjp[stride])
Reconstruct values in I' and J'.
Definition: cs_convection_diffusion_priv.h:584
static CS_F_HOST_DEVICE void cs_b_cd_unsteady(const cs_real_t bldfrp, const cs_real_3_t diipb, const cs_real_3_t gradi, const cs_real_t pi, cs_real_t *pip)
Handle preparation of boundary face values for the flux computation in case of an unsteady algorithm.
Definition: cs_convection_diffusion_priv.h:3422
static CS_F_HOST_DEVICE void cs_b_upwind_flux(const int iconvp, const cs_real_t thetap, const int imasac, const int inc, const int bc_type, const cs_real_t pi, const cs_real_t pir, const cs_real_t pipr, const cs_real_t coefap, const cs_real_t coefbp, const cs_real_t b_massflux, const cs_real_t xcpp, cs_real_t *flux)
Add convective flux (substracting the mass accumulation from it) to flux at boundary face....
Definition: cs_convection_diffusion_priv.h:3163
static CS_F_HOST_DEVICE void cs_solu_f_val(const cs_real_t cell_cen[3], const cs_real_t i_face_cog[3], const cs_real_t grad[3], const cs_real_t p, cs_real_t *pf)
Prepare value at face ij by using a Second Order Linear Upwind scheme.
Definition: cs_convection_diffusion_priv.h:794
static void cs_sync_scalar_halo(const cs_mesh_t *m, cs_halo_type_t halo_type, cs_real_t pvar[])
Definition: cs_convection_diffusion_priv.h:67
static CS_F_HOST_DEVICE void cs_i_relax_c_val(const double relaxp, const cs_real_t pia, const cs_real_t pja, const cs_real_t recoi, const cs_real_t recoj, const cs_real_t pi, const cs_real_t pj, cs_real_t *pir, cs_real_t *pjr, cs_real_t *pipr, cs_real_t *pjpr)
Compute relaxed values at cell i and j.
Definition: cs_convection_diffusion_priv.h:636
static CS_F_HOST_DEVICE cs_real_t cs_nvd_vof_scheme_scalar(cs_nvd_type_t limiter, const cs_real_t i_face_u_normal[3], cs_real_t nvf_p_c, cs_real_t nvf_r_f, cs_real_t nvf_r_c, const cs_real_t gradv_c[3], const cs_real_t c_courant)
Compute the normalised face scalar using the specified NVD scheme for the case of a Volume-of-Fluid (...
Definition: cs_convection_diffusion_priv.h:274
static CS_F_HOST_DEVICE void cs_b_cd_unsteady_strided(cs_real_t bldfrp, const cs_real_t diipb[3], const cs_real_t gradi[stride][3], const cs_real_t pi[stride], cs_real_t pip[stride])
Handle preparation of boundary face values for the flux computation in case of a steady algorithm.
Definition: cs_convection_diffusion_priv.h:3456
static CS_F_HOST_DEVICE void cs_b_upwind_flux_strided(int iconvp, cs_real_t thetap, int imasac, int inc, int bc_type, const cs_real_t pi[stride], const cs_real_t pir[stride], const cs_real_t pipr[stride], const cs_real_t coefa[stride], const cs_real_t coefb[stride][stride], cs_real_t b_massflux, cs_real_t pfac[stride], cs_real_t flux[stride])
Add convective flux (substracting the mass accumulation from it) to flux at boundary face....
Definition: cs_convection_diffusion_priv.h:3220
static CS_F_HOST_DEVICE void cs_b_diff_flux_coupling_strided(int idiffp, const cs_real_t pi[stride], const cs_real_t pj[stride], cs_real_t b_visc, cs_real_t fluxi[stride])
Add diffusive flux to flux at an internal coupling face for a vector.
Definition: cs_convection_diffusion_priv.h:3514
static CS_F_HOST_DEVICE void cs_slope_test(const cs_real_t pi, const cs_real_t pj, const cs_real_t distf, const cs_real_t i_face_u_normal[3], const cs_real_t gradi[3], const cs_real_t gradj[3], const cs_real_t grdpai[3], const cs_real_t grdpaj[3], const cs_real_t i_massflux, cs_real_t *testij, cs_real_t *tesqck)
Compute slope test criteria at internal face between cell i and j.
Definition: cs_convection_diffusion_priv.h:407
static CS_F_HOST_DEVICE void cs_upwind_f_val_strided(const cs_real_t p[stride], cs_real_t pf[stride])
Prepare value at face ij by using an upwind scheme.
Definition: cs_convection_diffusion_priv.h:729
static CS_F_HOST_DEVICE void cs_i_cd_unsteady_slope_test_strided(bool *upwind_switch, int iconvp, cs_real_t bldfrp, int ischcp, double blencp, double blend_st, cs_real_t weight, cs_real_t i_dist, const cs_real_t cell_ceni[3], const cs_real_t cell_cenj[3], const cs_real_t i_face_u_normal[3], const cs_real_t i_face_cog[3], const cs_real_t diipf[3], const cs_real_t djjpf[3], cs_real_t i_massflux, const cs_real_t gradi[stride][3], const cs_real_t gradj[stride][3], const cs_real_t grdpai[stride][3], const cs_real_t grdpaj[stride][3], const cs_real_t pi[stride], const cs_real_t pj[stride], cs_real_t pif[stride], cs_real_t pjf[stride], cs_real_t pip[stride], cs_real_t pjp[stride])
Handle preparation of internal face values for the fluxes computation in case of a unsteady algorithm...
Definition: cs_convection_diffusion_priv.h:2754
static CS_F_HOST_DEVICE void cs_b_relax_c_val(const double relaxp, const cs_real_t pi, const cs_real_t pia, const cs_real_t recoi, cs_real_t *pir, cs_real_t *pipr)
Compute relaxed values at boundary cell i.
Definition: cs_convection_diffusion_priv.h:2939
static CS_F_HOST_DEVICE void cs_b_imposed_conv_flux_strided(int iconvp, cs_real_t thetap, int imasac, int inc, int bc_type, int icvfli, const cs_real_t pi[stride], const cs_real_t pir[stride], const cs_real_t pipr[stride], const cs_real_t coefap[stride], const cs_real_t coefbp[stride][stride], const cs_real_t coface[stride], const cs_real_t cofbce[stride][stride], cs_real_t b_massflux, cs_real_t pfac[stride], cs_real_t flux[stride])
Add convective flux (substracting the mass accumulation from it) to flux at boundary face....
Definition: cs_convection_diffusion_priv.h:3082
static CS_F_HOST_DEVICE void cs_i_cd_unsteady_upwind_strided(const cs_real_t bldfrp, const cs_real_t diipf[3], const cs_real_t djjpf[3], const cs_real_t gradi[stride][3], const cs_real_t gradj[stride][3], const cs_real_t pi[stride], const cs_real_t pj[stride], cs_real_t pif[stride], cs_real_t pjf[stride], cs_real_t pip[stride], cs_real_t pjp[stride])
Handle preparation of internal face values for the fluxes computation in case of an unsteady algorith...
Definition: cs_convection_diffusion_priv.h:1293
static CS_F_HOST_DEVICE void cs_centered_f_val(double pnd, cs_real_t pip, cs_real_t pjp, cs_real_t *pf)
Prepare value at face ij by using a centered scheme.
Definition: cs_convection_diffusion_priv.h:748
static CS_F_HOST_DEVICE void cs_b_compute_quantities_strided(const cs_real_t diipb[3], const cs_real_t gradi[stride][3], const cs_real_t bldfrp, cs_real_t recoi[stride])
Reconstruct values in I' at boundary cell i.
Definition: cs_convection_diffusion_priv.h:2913
static CS_F_HOST_DEVICE void cs_blend_f_val_strided(const double blencp, const cs_real_t p[stride], cs_real_t pf[stride])
Blend face values for a centered or SOLU scheme with face values for an upwind scheme.
Definition: cs_convection_diffusion_priv.h:881
static CS_F_HOST_DEVICE void cs_i_compute_quantities(const cs_real_t bldfrp, const cs_real_3_t diipf, const cs_real_3_t djjpf, const cs_real_3_t gradi, const cs_real_3_t gradj, const cs_real_t pi, const cs_real_t pj, cs_real_t *recoi, cs_real_t *recoj, cs_real_t *pip, cs_real_t *pjp)
Reconstruct values in I' and J'.
Definition: cs_convection_diffusion_priv.h:538
static CS_F_HOST_DEVICE void cs_b_compute_quantities(const cs_real_3_t diipb, const cs_real_3_t gradi, const cs_real_t bldfrp, cs_real_t *recoi)
Reconstruct values in I' at boundary cell i.
Definition: cs_convection_diffusion_priv.h:2887
static CS_F_HOST_DEVICE void cs_i_conv_flux_strided(int iconvp, cs_real_t thetap, int imasac, const cs_real_t pi[stride], const cs_real_t pj[stride], const cs_real_t pifri[stride], const cs_real_t pifrj[stride], const cs_real_t pjfri[stride], const cs_real_t pjfrj[stride], cs_real_t i_massflux, cs_real_t fluxi[stride], cs_real_t fluxj[stride])
Add convective fluxes (substracting the mass accumulation from them) to fluxes at face ij.
Definition: cs_convection_diffusion_priv.h:961
static CS_F_HOST_DEVICE void cs_upwind_f_val(const cs_real_t p, cs_real_t *pf)
Prepare value at face ij by using an upwind scheme.
Definition: cs_convection_diffusion_priv.h:709
static CS_F_HOST_DEVICE void cs_b_diff_flux(const int idiffp, const cs_real_t thetap, const int inc, const cs_real_t pipr, const cs_real_t cofafp, const cs_real_t cofbfp, const cs_real_t b_visc, cs_real_t *flux)
Add diffusive flux to flux at boundary face.
Definition: cs_convection_diffusion_priv.h:3270
static CS_F_HOST_DEVICE void cs_i_cd_unsteady_nvd(const cs_nvd_type_t limiter, const double beta, const cs_real_3_t cell_cen_c, const cs_real_3_t cell_cen_d, const cs_real_3_t i_face_u_normal, const cs_real_3_t i_face_cog, const cs_real_3_t gradv_c, const cs_real_t p_c, const cs_real_t p_d, const cs_real_t local_max_c, const cs_real_t local_min_c, const cs_real_t courant_c, cs_real_t *pif, cs_real_t *pjf)
Handle preparation of internal face values for the convection flux computation in case of an unsteady...
Definition: cs_convection_diffusion_priv.h:2618
static CS_F_HOST_DEVICE void cs_i_cd_steady_strided(cs_real_t bldfrp, int ischcp, double relaxp, double blencp, cs_real_t weight, const cs_real_t cell_ceni[3], const cs_real_t cell_cenj[3], const cs_real_t i_face_cog[3], const cs_real_t diipf[3], const cs_real_t djjpf[3], const cs_real_t gradi[stride][3], const cs_real_t gradj[stride][3], const cs_real_t pi[stride], const cs_real_t pj[stride], const cs_real_t pia[stride], const cs_real_t pja[stride], cs_real_t pifri[stride], cs_real_t pifrj[stride], cs_real_t pjfri[stride], cs_real_t pjfrj[stride], cs_real_t pip[stride], cs_real_t pjp[stride], cs_real_t pipr[stride], cs_real_t pjpr[stride])
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion_priv.h:1544
static CS_F_HOST_DEVICE void cs_i_diff_flux(const int idiffp, const cs_real_t thetap, const cs_real_t pip, const cs_real_t pjp, const cs_real_t pipr, const cs_real_t pjpr, const cs_real_t i_visc, cs_real_t fluxij[2])
Add diffusive fluxes to fluxes at face ij.
Definition: cs_convection_diffusion_priv.h:1003
static CS_F_HOST_DEVICE void cs_i_cd_unsteady_upwind(const cs_real_t bldfrp, const cs_real_t diipf[3], const cs_real_t djjpf[3], const cs_real_t gradi[3], const cs_real_t gradj[3], const cs_real_t pi, const cs_real_t pj, cs_real_t *pif, cs_real_t *pjf, cs_real_t *pip, cs_real_t *pjp)
Handle preparation of internal face values for the fluxes computation in case of an unsteady algorith...
Definition: cs_convection_diffusion_priv.h:1239
static CS_F_HOST_DEVICE void cs_solu_f_val_strided(const cs_real_t cell_cen[3], const cs_real_t i_face_cog[3], const cs_real_t grad[stride][3], const cs_real_t p[stride], cs_real_t pf[stride])
Prepare value at face ij by using a Second Order Linear Upwind scheme.
Definition: cs_convection_diffusion_priv.h:826
static CS_F_HOST_DEVICE void cs_b_cd_steady_strided(cs_real_t bldfrp, double relaxp, const cs_real_t diipb[3], const cs_real_t gradi[stride][3], const cs_real_t pi[stride], const cs_real_t pia[stride], cs_real_t pir[stride], cs_real_t pipr[stride])
Handle preparation of boundary face values for the flux computation in case of a steady algorithm.
Definition: cs_convection_diffusion_priv.h:3384
static CS_F_HOST_DEVICE void cs_central_downwind_cells(const cs_lnum_t ii, const cs_lnum_t jj, const cs_real_t i_massflux, cs_lnum_t *ic, cs_lnum_t *id)
Determine the upwind and downwind sides of an internal face and matching cell indices.
Definition: cs_convection_diffusion_priv.h:2579
static CS_F_HOST_DEVICE void cs_i_conv_flux(const int iconvp, const cs_real_t thetap, const int imasac, const cs_real_t pi, const cs_real_t pj, const cs_real_t pifri, const cs_real_t pifrj, const cs_real_t pjfri, const cs_real_t pjfrj, const cs_real_t i_massflux, const cs_real_t xcppi, const cs_real_t xcppj, cs_real_2_t fluxij)
Add convective fluxes (substracting the mass accumulation from them) to fluxes at face ij.
Definition: cs_convection_diffusion_priv.h:913
static CS_F_HOST_DEVICE void cs_i_cd_steady(const cs_real_t bldfrp, const int ischcp, const double relaxp, const double blencp, const cs_real_t weight, const cs_real_t cell_ceni[3], const cs_real_t cell_cenj[3], const cs_real_t i_face_cog[3], const cs_real_t diipf[3], const cs_real_t djjpf[3], const cs_real_t gradi[3], const cs_real_t gradj[3], const cs_real_t gradupi[3], const cs_real_t gradupj[3], const cs_real_t pi, const cs_real_t pj, const cs_real_t pia, const cs_real_t pja, cs_real_t *pifri, cs_real_t *pifrj, cs_real_t *pjfri, cs_real_t *pjfrj, cs_real_t *pip, cs_real_t *pjp, cs_real_t *pipr, cs_real_t *pjpr)
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion_priv.h:1359
static CS_F_HOST_DEVICE cs_real_t cs_nvd_scheme_scalar(const cs_nvd_type_t limiter, const cs_real_t nvf_p_c, const cs_real_t nvf_r_f, const cs_real_t nvf_r_c)
Compute the normalized face scalar using the specified NVD scheme.
Definition: cs_convection_diffusion_priv.h:89
static CS_F_HOST_DEVICE void cs_blend_f_val(const double blencp, const cs_real_t p, cs_real_t *pf)
Blend face values for a centered or SOLU scheme with face values for an upwind scheme.
Definition: cs_convection_diffusion_priv.h:857
static CS_F_HOST_DEVICE void cs_i_cd_steady_slope_test_strided(bool *upwind_switch, int iconvp, cs_real_t bldfrp, int ischcp, double relaxp, double blencp, double blend_st, cs_real_t weight, cs_real_t i_dist, const cs_real_t cell_ceni[3], const cs_real_t cell_cenj[3], const cs_real_t i_face_u_normal[3], const cs_real_t i_face_cog[3], const cs_real_t diipf[3], const cs_real_t djjpf[3], cs_real_t i_massflux, const cs_real_t gradi[stride][3], const cs_real_t gradj[stride][3], const cs_real_t grdpai[stride][3], const cs_real_t grdpaj[stride][3], const cs_real_t pi[stride], const cs_real_t pj[stride], const cs_real_t pia[stride], const cs_real_t pja[stride], cs_real_t pifri[stride], cs_real_t pifrj[stride], cs_real_t pjfri[stride], cs_real_t pjfrj[stride], cs_real_t pip[stride], cs_real_t pjp[stride], cs_real_t pipr[stride], cs_real_t pjpr[stride])
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion_priv.h:2241
static CS_F_HOST_DEVICE void cs_slope_test_strided(const cs_real_t pi[stride], const cs_real_t pj[stride], const cs_real_t distf, const cs_real_t i_face_u_normal[3], const cs_real_t gradi[stride][3], const cs_real_t gradj[stride][3], const cs_real_t gradsti[stride][3], const cs_real_t gradstj[stride][3], const cs_real_t i_massflux, cs_real_t *testij, cs_real_t *tesqck)
Compute slope test criteria at internal face between cell i and j.
Definition: cs_convection_diffusion_priv.h:472
static CS_F_HOST_DEVICE void cs_b_imposed_conv_flux(int iconvp, cs_real_t thetap, int imasac, int inc, int bc_type, int icvfli, cs_real_t pi, cs_real_t pir, cs_real_t pipr, cs_real_t coefap, cs_real_t coefbp, cs_real_t coface, cs_real_t cofbce, cs_real_t b_massflux, cs_real_t xcpp, cs_real_t *flux)
Add convective flux (substracting the mass accumulation from it) to flux at boundary face....
Definition: cs_convection_diffusion_priv.h:3008
static CS_F_HOST_DEVICE void cs_i_cd_steady_upwind(const cs_real_t bldfrp, const cs_real_t relaxp, const cs_real_t diipf[3], const cs_real_t djjpf[3], const cs_real_t gradi[3], const cs_real_t gradj[3], const cs_real_t pi, const cs_real_t pj, const cs_real_t pia, const cs_real_t pja, cs_real_t *pifri, cs_real_t *pifrj, cs_real_t *pjfri, cs_real_t *pjfrj, cs_real_t *pip, cs_real_t *pjp, cs_real_t *pipr, cs_real_t *pjpr)
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion_priv.h:1080
#define CS_F_HOST_DEVICE
Definition: cs_defs.h:546
double cs_real_t
Floating-point value.
Definition: cs_defs.h:332
#define CS_ABS(a)
Definition: cs_defs.h:490
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:347
cs_real_t cs_real_2_t[2]
vector of 2 floating-point values
Definition: cs_defs.h:346
#define CS_UNUSED(x)
Definition: cs_defs.h:514
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:325
@ p
Definition: cs_field_pointer.h:67
@ k
Definition: cs_field_pointer.h:70
void cs_halo_sync_var(const cs_halo_t *halo, cs_halo_type_t sync_mode, cs_real_t var[])
Definition: cs_halo.cpp:2071
cs_halo_type_t
Definition: cs_halo.h:56
static CS_F_HOST_DEVICE cs_real_t cs_math_3_distance(const cs_real_t xa[3], const cs_real_t xb[3])
Compute the (euclidean) distance between two points xa and xb in a Cartesian coordinate system of dim...
Definition: cs_math.h:324
static CS_F_HOST_DEVICE cs_real_t cs_math_3_dot_product(const cs_real_t u[3], const cs_real_t v[3])
Compute the dot product of two vectors of 3 real values.
Definition: cs_math.h:391
const cs_real_t cs_math_epzero
static CS_F_HOST_DEVICE cs_real_t cs_math_3_norm(const cs_real_t v[3])
Compute the euclidean norm of a vector of dimension 3.
Definition: cs_math.h:460
static CS_F_HOST_DEVICE cs_real_t cs_math_fabs(cs_real_t x)
Compute the absolute value of a real value.
Definition: cs_math.h:163
static CS_F_HOST_DEVICE cs_real_t cs_math_fmin(cs_real_t x, cs_real_t y)
Compute the min value of two real values.
Definition: cs_math.h:181
static CS_F_HOST_DEVICE cs_real_t cs_math_sq(cs_real_t x)
Compute the square of a real value.
Definition: cs_math.h:242
static CS_F_HOST_DEVICE cs_real_t cs_math_fmax(cs_real_t x, cs_real_t y)
Compute the max value of two real values.
Definition: cs_math.h:200
@ CS_COUPLED_FD
Definition: cs_parameters.h:101
integer, dimension(:), allocatable, target icvfli
boundary convection flux indicator of a Rusanov or an analytical flux (some boundary contributions of...
Definition: cfpoin.f90:55
double precision pi
value with 16 digits
Definition: cstnum.f90:48
double precision, dimension(:,:), pointer diipb
vector II' for interior faces for every boundary face, the three components of the vector ....
Definition: mesh.f90:168
real(c_double), pointer, save fmin
Definition: coincl.f90:185
real(c_double), dimension(:), pointer, save b2
Definition: cpincl.f90:109
real(c_double), dimension(:), pointer, save b1
Definition: cpincl.f90:109
Definition: cs_mesh.h:85
cs_halo_t * halo
Definition: cs_mesh.h:156