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