programmer's documentation
cs_convection_diffusion.h
Go to the documentation of this file.
1 #ifndef __CS_CONVECTION_DIFFUSION_H__
2 #define __CS_CONVECTION_DIFFUSION_H__
3 
4 /*============================================================================
5  * Convection-diffusion operators.
6  *============================================================================*/
7 
8 /*
9  This file is part of Code_Saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2015 EDF S.A.
12 
13  This program is free software; you can redistribute it and/or modify it under
14  the terms of the GNU General Public License as published by the Free Software
15  Foundation; either version 2 of the License, or (at your option) any later
16  version.
17 
18  This program is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20  FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
21  details.
22 
23  You should have received a copy of the GNU General Public License along with
24  this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
25  Street, Fifth Floor, Boston, MA 02110-1301, USA.
26 */
27 
28 /*----------------------------------------------------------------------------*/
29 
30 #include "cs_defs.h"
31 
32 /*----------------------------------------------------------------------------
33  * Local headers
34  *----------------------------------------------------------------------------*/
35 
36 #include "cs_base.h"
37 #include "cs_halo.h"
38 
39 /*----------------------------------------------------------------------------*/
40 
42 
43 /*=============================================================================
44  * Local Macro definitions
45  *============================================================================*/
46 
47 /*============================================================================
48  * Type definition
49  *============================================================================*/
50 
51 /*============================================================================
52  * Global variables
53  *============================================================================*/
54 
55 /*============================================================================
56  * Private function definitions
57  *============================================================================*/
58 
59 /*----------------------------------------------------------------------------*/
76 /*----------------------------------------------------------------------------*/
77 
78 inline static void
80  const cs_real_t pj,
81  const cs_real_t distf,
82  const cs_real_t srfan,
83  const cs_real_3_t i_face_normal,
84  const cs_real_3_t gradi,
85  const cs_real_3_t gradj,
86  const cs_real_3_t grdpai,
87  const cs_real_3_t grdpaj,
88  const cs_real_t i_massflux,
89  double *testij,
90  double *tesqck)
91 {
92  double testi, testj;
93  double dcc, ddi, ddj;
94 
95  /* Slope test
96  ----------*/
97 
98  testi = grdpai[0]*i_face_normal[0]
99  + grdpai[1]*i_face_normal[1]
100  + grdpai[2]*i_face_normal[2];
101  testj = grdpaj[0]*i_face_normal[0]
102  + grdpaj[1]*i_face_normal[1]
103  + grdpaj[2]*i_face_normal[2];
104  *testij = grdpai[0]*grdpaj[0]
105  + grdpai[1]*grdpaj[1]
106  + grdpai[2]*grdpaj[2];
107 
108  if (i_massflux>0.) {
109  dcc = gradi[0]*i_face_normal[0]
110  + gradi[1]*i_face_normal[1]
111  + gradi[2]*i_face_normal[2];
112  ddi = testi;
113  ddj = (pj-pi)/distf *srfan;
114  } else {
115  dcc = gradj[0]*i_face_normal[0]
116  + gradj[1]*i_face_normal[1]
117  + gradj[2]*i_face_normal[2];
118  ddi = (pj-pi)/distf *srfan;
119  ddj = testj;
120  }
121  *tesqck = pow(dcc, 2.) - pow(ddi-ddj, 2.);
122 }
123 /*----------------------------------------------------------------------------*/
140 /*----------------------------------------------------------------------------*/
141 
142 inline static void
144  const cs_real_6_t pj,
145  const cs_real_t distf,
146  const cs_real_t srfan,
147  const cs_real_3_t i_face_normal,
148  const cs_real_63_t gradi,
149  const cs_real_63_t gradj,
150  const cs_real_63_t grdpai,
151  const cs_real_63_t grdpaj,
152  const cs_real_t i_massflux,
153  cs_real_t testij[6],
154  cs_real_t tesqck[6])
155 {
156  double testi[6], testj[6];
157  double dcc[6], ddi[6], ddj[6];
158 
159  /* Slope test
160  ----------*/
161  for (int isou = 0; isou < 6; isou++) {
162  testi[isou] = grdpai[isou][0]*i_face_normal[0]
163  + grdpai[isou][1]*i_face_normal[1]
164  + grdpai[isou][2]*i_face_normal[2];
165  testj[isou] = grdpaj[isou][0]*i_face_normal[0]
166  + grdpaj[isou][1]*i_face_normal[1]
167  + grdpaj[isou][2]*i_face_normal[2];
168  testij[isou] = grdpai[isou][0]*grdpaj[isou][0]
169  + grdpai[isou][1]*grdpaj[isou][1]
170  + grdpai[isou][2]*grdpaj[isou][2];
171 
172  if (i_massflux>0.) {
173  dcc[isou] = gradi[isou][0]*i_face_normal[0]
174  + gradi[isou][1]*i_face_normal[1]
175  + gradi[isou][2]*i_face_normal[2];
176  ddi[isou] = testi[isou];
177  ddj[isou] = (pj[isou]-pi[isou])/distf *srfan;
178  } else {
179  dcc[isou] = gradj[isou][0]*i_face_normal[0]
180  + gradj[isou][1]*i_face_normal[1]
181  + gradj[isou][2]*i_face_normal[2];
182  ddi[isou] = (pj[isou]-pi[isou])/distf *srfan;
183  ddj[isou] = testj[isou];
184  }
185  tesqck[isou] = pow(dcc[isou], 2.) - pow(ddi[isou]-ddj[isou], 2.);
186  }
187 }
188 
189 /*----------------------------------------------------------------------------*/
208 /*----------------------------------------------------------------------------*/
209 
210 inline static void
211 cs_i_compute_quantities(const int ircflp,
212  const double pnd,
213  const cs_real_3_t cell_ceni,
214  const cs_real_3_t cell_cenj,
215  const cs_real_3_t i_face_cog,
216  const cs_real_3_t dijpf,
217  const cs_real_3_t gradi,
218  const cs_real_3_t gradj,
219  const cs_real_t pi,
220  const cs_real_t pj,
221  cs_real_t *recoi,
222  cs_real_t *recoj,
223  cs_real_t *pip,
224  cs_real_t *pjp)
225 {
226  cs_real_t diipfx, diipfy, diipfz, djjpfx, djjpfy, djjpfz;
227  cs_real_t dpxf, dpyf, dpzf;
228 
229  /* Recompute II' and JJ' at this level */
230 
231  diipfx = i_face_cog[0] - (cell_ceni[0] + (1.-pnd) * dijpf[0]);
232  diipfy = i_face_cog[1] - (cell_ceni[1] + (1.-pnd) * dijpf[1]);
233  diipfz = i_face_cog[2] - (cell_ceni[2] + (1.-pnd) * dijpf[2]);
234 
235  djjpfx = i_face_cog[0] - cell_cenj[0] + pnd * dijpf[0];
236  djjpfy = i_face_cog[1] - cell_cenj[1] + pnd * dijpf[1];
237  djjpfz = i_face_cog[2] - cell_cenj[2] + pnd * dijpf[2];
238 
239  dpxf = 0.5*(gradi[0] + gradj[0]);
240  dpyf = 0.5*(gradi[1] + gradj[1]);
241  dpzf = 0.5*(gradi[2] + gradj[2]);
242 
243  /* reconstruction only if IRCFLP = 1 */
244  *recoi = ircflp*(dpxf*diipfx+dpyf*diipfy+dpzf*diipfz);
245  *recoj = ircflp*(dpxf*djjpfx+dpyf*djjpfy+dpzf*djjpfz);
246  *pip = pi + *recoi;
247  *pjp = pj + *recoj;
248 }
249 
250 /*----------------------------------------------------------------------------*/
269 /*----------------------------------------------------------------------------*/
270 
271 inline static void
273  const double pnd,
274  const cs_real_3_t cell_ceni,
275  const cs_real_3_t cell_cenj,
276  const cs_real_3_t i_face_cog,
277  const cs_real_3_t dijpf,
278  const cs_real_63_t gradi,
279  const cs_real_63_t gradj,
280  const cs_real_6_t pi,
281  const cs_real_6_t pj,
282  cs_real_t recoi[6],
283  cs_real_t recoj[6],
284  cs_real_t pip[6],
285  cs_real_t pjp[6])
286 {
287  cs_real_3_t dijpfv, diipfv, djjpfv;
288  cs_real_3_t dpvf;
289 
290 
291  for (int jsou = 0; jsou < 3; jsou++) {
292  dijpfv[jsou] = dijpf[jsou];
293  }
294 
295  /* Recompute II' and JJ' at this level */
296  for (int jsou = 0; jsou < 3; jsou++) {
297  diipfv[jsou] = i_face_cog[jsou]
298  - (cell_ceni[jsou] + (1.-pnd) * dijpfv[jsou]);
299  djjpfv[jsou] = i_face_cog[jsou]
300  - cell_cenj[jsou] + pnd * dijpfv[jsou];
301  }
302 
303 
304  /*-----------------
305  X-Y-Z components, p=u, v, w */
306  for (int isou = 0; isou < 6; isou++) {
307 
308  for (int jsou = 0; jsou < 3; jsou++)
309  dpvf[jsou] = 0.5*( gradi[isou][jsou]
310  + gradj[isou][jsou]);
311 
312  /* reconstruction only if IRCFLP = 1 */
313 
314  recoi[isou] = ircflp*( dpvf[0]*diipfv[0]
315  + dpvf[1]*diipfv[1]
316  + dpvf[2]*diipfv[2]);
317 
318 
319  recoj[isou] = ircflp*( dpvf[0]*djjpfv[0]
320  + dpvf[1]*djjpfv[1]
321  + dpvf[2]*djjpfv[2]);
322 
323 
324  pip[isou] = pi[isou] + recoi[isou];
325 
326  pjp[isou] = pj[isou] + recoj[isou];
327 
328  }
329 }
330 
331 
332 /*----------------------------------------------------------------------------*/
348 /*----------------------------------------------------------------------------*/
349 
350 inline static void
351 cs_i_relax_c_val(const double relaxp,
352  const cs_real_t pia,
353  const cs_real_t pja,
354  const cs_real_t recoi,
355  const cs_real_t recoj,
356  const cs_real_t pi,
357  const cs_real_t pj,
358  cs_real_t *pir,
359  cs_real_t *pjr,
360  cs_real_t *pipr,
361  cs_real_t *pjpr)
362 {
363  *pir = pi/relaxp - (1.-relaxp)/relaxp * pia;
364  *pjr = pj/relaxp - (1.-relaxp)/relaxp * pja;
365 
366  *pipr = *pir + recoi;
367  *pjpr = *pjr + recoj;
368 }
369 
370 /*----------------------------------------------------------------------------*/
386 /*----------------------------------------------------------------------------*/
387 
388 inline static void
389 cs_i_relax_c_val_tensor(const double relaxp,
390  const cs_real_6_t pia,
391  const cs_real_6_t pja,
392  const cs_real_6_t recoi,
393  const cs_real_6_t recoj,
394  const cs_real_6_t pi,
395  const cs_real_6_t pj,
396  cs_real_t pir[6],
397  cs_real_t pjr[6],
398  cs_real_t pipr[6],
399  cs_real_t pjpr[6])
400 {
401  for (int isou = 0; isou < 6; isou++){
402 
403  pir[isou] = pi[isou] /relaxp - (1.-relaxp)/relaxp * pia[isou];
404  pjr[isou] = pj[isou] /relaxp - (1.-relaxp)/relaxp * pja[isou];
405 
406  pipr[isou] = pir[isou] + recoi[isou];
407  pjpr[isou] = pjr[isou] + recoj[isou];
408  }
409 }
410 
411 /*----------------------------------------------------------------------------*/
425 /*----------------------------------------------------------------------------*/
426 
427 inline static void
429  const cs_real_t pj,
430  const cs_real_t pip,
431  const cs_real_t pjp,
432  cs_real_t *pir,
433  cs_real_t *pjr,
434  cs_real_t *pipr,
435  cs_real_t *pjpr)
436 {
437  *pir = pi;
438  *pjr = pj;
439 
440  *pipr = pip;
441  *pjpr = pjp;
442 }
443 
444 /*----------------------------------------------------------------------------*/
458 /*----------------------------------------------------------------------------*/
459 
460 inline static void
462  const cs_real_6_t pj,
463  const cs_real_6_t pip,
464  const cs_real_6_t pjp,
465  cs_real_t pir[6],
466  cs_real_t pjr[6],
467  cs_real_t pipr[6],
468  cs_real_t pjpr[6])
469 {
470 
471  for (int isou = 0 ; isou < 6 ; isou++){
472  pir[isou] = pi[isou];
473  pjr[isou] = pj[isou];
474 
475  pipr[isou] = pip[isou];
476  pjpr[isou] = pjp[isou];
477  }
478 }
479 
480 /*----------------------------------------------------------------------------*/
493 /*----------------------------------------------------------------------------*/
494 
495 inline static void
497  const cs_real_t pj,
498  const cs_real_t pir,
499  const cs_real_t pjr,
500  cs_real_t *pifri,
501  cs_real_t *pifrj,
502  cs_real_t *pjfri,
503  cs_real_t *pjfrj)
504 {
505  *pifri = pir;
506  *pifrj = pi;
507  *pjfri = pj;
508  *pjfrj = pjr;
509 }
510 
511 /*----------------------------------------------------------------------------*/
524 /*----------------------------------------------------------------------------*/
525 
526 inline static void
528  const cs_real_6_t pj,
529  const cs_real_6_t pir,
530  const cs_real_6_t pjr,
531  cs_real_t pifri[6],
532  cs_real_t pifrj[6],
533  cs_real_t pjfri[6],
534  cs_real_t pjfrj[6])
535 {
536  int isou;
537 
538  for (isou = 0; isou < 6; isou++){
539  pifri[isou] = pir[isou];
540  pifrj[isou] = pi[isou];
541  pjfri[isou] = pj[isou];
542  pjfrj[isou] = pjr[isou];
543  }
544 }
545 
546 /*----------------------------------------------------------------------------*/
560 /*----------------------------------------------------------------------------*/
561 
562 inline static void
563 cs_centered_f_val(const double pnd,
564  const cs_real_t pip,
565  const cs_real_t pjp,
566  const cs_real_t pipr,
567  const cs_real_t pjpr,
568  cs_real_t *pifri,
569  cs_real_t *pifrj,
570  cs_real_t *pjfri,
571  cs_real_t *pjfrj)
572 {
573  *pifri = pnd*pipr + (1.-pnd)*pjp;
574  *pjfri = *pifri;
575  *pifrj = pnd*pip + (1.-pnd)*pjpr;
576  *pjfrj = *pifrj;
577 }
578 
579 /*----------------------------------------------------------------------------*/
593 /*----------------------------------------------------------------------------*/
594 
595 inline static void
596 cs_centered_f_val_tensor(const double pnd,
597  const cs_real_6_t pip,
598  const cs_real_6_t pjp,
599  const cs_real_6_t pipr,
600  const cs_real_6_t pjpr,
601  cs_real_t pifri[6],
602  cs_real_t pifrj[6],
603  cs_real_t pjfri[6],
604  cs_real_t pjfrj[6])
605 {
606  for (int isou = 0 ; isou < 6 ; isou++){
607  pifri[isou] = pnd*pipr[isou] + (1.-pnd)*pjp[isou];
608  pjfri[isou] = pifri[isou];
609  pifrj[isou] = pnd*pip[isou] + (1.-pnd)*pjpr[isou];
610  pjfrj[isou] = pifrj[isou];
611  }
612 }
613 
614 /*----------------------------------------------------------------------------*/
632 /*----------------------------------------------------------------------------*/
633 
634 inline static void
635 cs_solu_f_val(const cs_real_3_t cell_ceni,
636  const cs_real_3_t cell_cenj,
637  const cs_real_3_t i_face_cog,
638  const cs_real_3_t gradi,
639  const cs_real_3_t gradj,
640  const cs_real_t pi,
641  const cs_real_t pj,
642  const cs_real_t pir,
643  const cs_real_t pjr,
644  cs_real_t *pifri,
645  cs_real_t *pifrj,
646  cs_real_t *pjfri,
647  cs_real_t *pjfrj)
648 {
649  cs_real_t difx, dify, difz, djfx, djfy, djfz;
650 
651  difx = i_face_cog[0] - cell_ceni[0];
652  dify = i_face_cog[1] - cell_ceni[1];
653  difz = i_face_cog[2] - cell_ceni[2];
654  djfx = i_face_cog[0] - cell_cenj[0];
655  djfy = i_face_cog[1] - cell_cenj[1];
656  djfz = i_face_cog[2] - cell_cenj[2];
657 
658  /* leave reconstruction of PIF and PJF even if IRCFLP=0
659  otherwise, it is the same as using upwind */
660  *pifri = pir + difx*gradi[0]+dify*gradi[1]+difz*gradi[2];
661  *pifrj = pi + difx*gradi[0]+dify*gradi[1]+difz*gradi[2];
662  *pjfrj = pjr + djfx*gradj[0]+djfy*gradj[1]+djfz*gradj[2];
663  *pjfri = pj + djfx*gradj[0]+djfy*gradj[1]+djfz*gradj[2];
664 }
665 
666 /*----------------------------------------------------------------------------*/
684 /*----------------------------------------------------------------------------*/
685 
686 inline static void
688  const cs_real_3_t cell_cenj,
689  const cs_real_3_t i_face_cog,
690  const cs_real_63_t gradi,
691  const cs_real_63_t gradj,
692  const cs_real_6_t pi,
693  const cs_real_6_t pj,
694  const cs_real_6_t pir,
695  const cs_real_6_t pjr,
696  cs_real_t pifri[6],
697  cs_real_t pifrj[6],
698  cs_real_t pjfri[6],
699  cs_real_t pjfrj[6])
700 {
701  cs_real_3_t difv, djfv;
702 
703  for (int jsou = 0; jsou < 3; jsou++) {
704  difv[jsou] = i_face_cog[jsou] - cell_ceni[jsou];
705  djfv[jsou] = i_face_cog[jsou] - cell_cenj[jsou];
706  }
707 
708 
709  /* leave reconstruction of PIF and PJF even if IRCFLP=0
710  otherwise, it is the same as using upwind */
711 
712  for (int isou = 0; isou < 6; isou++) {
713  pifri[isou] = pir[isou]+ difv[0]*gradi[isou][0]
714  + difv[1]*gradi[isou][1]
715  + difv[2]*gradi[isou][2];
716  pifrj[isou] = pi[isou] + difv[0]*gradi[isou][0]
717  + difv[1]*gradi[isou][1]
718  + difv[2]*gradi[isou][2];
719 
720  pjfrj[isou] = pjr[isou]+ djfv[0]*gradj[isou][0]
721  + djfv[1]*gradj[isou][1]
722  + djfv[2]*gradj[isou][2];
723  pjfri[isou] = pj[isou] + djfv[0]*gradj[isou][0]
724  + djfv[1]*gradj[isou][1]
725  + djfv[2]*gradj[isou][2];
726  }
727 }
728 
729 
730 /*----------------------------------------------------------------------------*/
746 /*----------------------------------------------------------------------------*/
747 
748 inline static void
749 cs_blend_f_val(const double blencp,
750  const cs_real_t pi,
751  const cs_real_t pj,
752  const cs_real_t pir,
753  const cs_real_t pjr,
754  cs_real_t *pifri,
755  cs_real_t *pifrj,
756  cs_real_t *pjfri,
757  cs_real_t *pjfrj)
758 {
759  *pifri = blencp*(*pifri)+(1.-blencp)*pir;
760  *pifrj = blencp*(*pifrj)+(1.-blencp)*pi;
761  *pjfri = blencp*(*pjfri)+(1.-blencp)*pj;
762  *pjfrj = blencp*(*pjfrj)+(1.-blencp)*pjr;
763 }
764 
765 /*----------------------------------------------------------------------------*/
781 /*----------------------------------------------------------------------------*/
782 
783 inline static void
784 cs_blend_f_val_tensor(const double blencp,
785  const cs_real_6_t pi,
786  const cs_real_6_t pj,
787  const cs_real_6_t pir,
788  const cs_real_6_t pjr,
789  cs_real_t pifri[6],
790  cs_real_t pifrj[6],
791  cs_real_t pjfri[6],
792  cs_real_t pjfrj[6])
793 {
794  for (int isou = 0; isou < 6; isou++) {
795  pifri[isou] = blencp*(pifri[isou])+(1.-blencp)*pir[isou];
796  pifrj[isou] = blencp*(pifrj[isou])+(1.-blencp)*pi[isou];
797  pjfri[isou] = blencp*(pjfri[isou])+(1.-blencp)*pj[isou];
798  pjfrj[isou] = blencp*(pjfrj[isou])+(1.-blencp)*pjr[isou];
799  }
800 }
801 
802 
803 /*----------------------------------------------------------------------------*/
818 /*----------------------------------------------------------------------------*/
819 
820 inline static void
821 cs_i_conv_flux(const int iconvp,
822  const cs_real_t pi,
823  const cs_real_t pj,
824  const cs_real_t pifri,
825  const cs_real_t pifrj,
826  const cs_real_t pjfri,
827  const cs_real_t pjfrj,
828  const cs_real_t i_massflux,
829  cs_real_2_t fluxij)
830 {
831  cs_real_t flui, fluj;
832 
833  flui = 0.5*(i_massflux + fabs(i_massflux));
834  fluj = 0.5*(i_massflux - fabs(i_massflux));
835 
836  fluxij[0] += iconvp*(flui*pifri + fluj*pjfri - i_massflux*pi);
837  fluxij[1] += iconvp*(flui*pifrj + fluj*pjfrj - i_massflux*pj);
838 }
839 
840 /*----------------------------------------------------------------------------*/
856 /*----------------------------------------------------------------------------*/
857 
858 inline static void
859 cs_i_conv_flux_tensor(const int iconvp,
860  const cs_real_6_t pi,
861  const cs_real_6_t pj,
862  const cs_real_6_t pifri,
863  const cs_real_6_t pifrj,
864  const cs_real_6_t pjfri,
865  const cs_real_6_t pjfrj,
866  const cs_real_t i_massflux,
867  cs_real_t fluxi[6],
868  cs_real_t fluxj[6])
869 {
870  cs_real_t flui, fluj;
871 
872  flui = 0.5*(i_massflux + fabs(i_massflux));
873  fluj = 0.5*(i_massflux - fabs(i_massflux));
874 
875  for (int isou = 0; isou < 6; isou++){
876 
877  fluxi[isou] += iconvp*(flui*pifri[isou] + fluj*pjfri[isou] - i_massflux*pi[isou]);
878  fluxj[isou] += iconvp*(flui*pifrj[isou] + fluj*pjfrj[isou] - i_massflux*pj[isou]);
879  }
880 }
881 
882 /*----------------------------------------------------------------------------*/
898 /*----------------------------------------------------------------------------*/
899 
900 inline static void
901 cs_i_conv_flux_cons(const int iconvp,
902  const cs_real_t pifri,
903  const cs_real_t pifrj,
904  const cs_real_t pjfri,
905  const cs_real_t pjfrj,
906  const cs_real_t i_massflux,
907  const cs_real_t xcppi,
908  const cs_real_t xcppj,
909  cs_real_2_t fluxij)
910 {
911  cs_real_t flui, fluj;
912 
913  flui = 0.5*(i_massflux +fabs(i_massflux));
914  fluj = 0.5*(i_massflux -fabs(i_massflux));
915 
916  fluxij[0] += iconvp*xcppi*(flui*pifri + fluj*pjfri);
917  fluxij[1] += iconvp*xcppj*(flui*pifrj + fluj*pjfrj);
918 }
919 
920 /*----------------------------------------------------------------------------*/
932 /*----------------------------------------------------------------------------*/
933 
934 inline static void
935 cs_i_diff_flux(const int idiffp,
936  const cs_real_t pip,
937  const cs_real_t pjp,
938  const cs_real_t pipr,
939  const cs_real_t pjpr,
940  const cs_real_t i_visc,
941  cs_real_2_t fluxij)
942 {
943  fluxij[0] += idiffp*i_visc*(pipr -pjp);
944  fluxij[1] += idiffp*i_visc*(pip -pjpr);
945 }
946 
947 /*----------------------------------------------------------------------------*/
960 /*----------------------------------------------------------------------------*/
961 
962 inline static void
963 cs_i_diff_flux_tensor(const int idiffp,
964  const cs_real_6_t pip,
965  const cs_real_6_t pjp,
966  const cs_real_6_t pipr,
967  const cs_real_6_t pjpr,
968  const cs_real_t i_visc,
969  cs_real_t fluxi[6],
970  cs_real_t fluxj[6])
971 {
972  for (int isou = 0; isou < 6; isou++){
973  fluxi[isou] += idiffp*i_visc*(pipr[isou] -pjp[isou]);
974  fluxj[isou] += idiffp*i_visc*(pip[isou] -pjpr[isou]);
975  }
976 }
977 
978 /*----------------------------------------------------------------------------*/
1005 /*----------------------------------------------------------------------------*/
1006 
1007 inline static void
1008 cs_i_cd_steady_upwind(const int ircflp,
1009  const cs_real_t relaxp,
1010  const cs_real_t weight,
1011  const cs_real_3_t cell_ceni,
1012  const cs_real_3_t cell_cenj,
1013  const cs_real_3_t i_face_cog,
1014  const cs_real_3_t dijpf,
1015  const cs_real_3_t gradi,
1016  const cs_real_3_t gradj,
1017  const cs_real_t pi,
1018  const cs_real_t pj,
1019  const cs_real_t pia,
1020  const cs_real_t pja,
1021  cs_real_t *pifri,
1022  cs_real_t *pifrj,
1023  cs_real_t *pjfri,
1024  cs_real_t *pjfrj,
1025  cs_real_t *pip,
1026  cs_real_t *pjp,
1027  cs_real_t *pipr,
1028  cs_real_t *pjpr)
1029 {
1030  cs_real_t pir, pjr;
1031  cs_real_t recoi, recoj;
1032 
1033  cs_i_compute_quantities(ircflp,
1034  weight,
1035  cell_ceni,
1036  cell_cenj,
1037  i_face_cog,
1038  dijpf,
1039  gradi,
1040  gradj,
1041  pi,
1042  pj,
1043  &recoi,
1044  &recoj,
1045  pip,
1046  pjp);
1047 
1048  cs_i_relax_c_val(relaxp,
1049  pia,
1050  pja,
1051  recoi,
1052  recoj,
1053  pi,
1054  pj,
1055  &pir,
1056  &pjr,
1057  pipr,
1058  pjpr);
1059 
1060  cs_upwind_f_val(pi,
1061  pj,
1062  pir,
1063  pjr,
1064  pifri,
1065  pifrj,
1066  pjfri,
1067  pjfrj);
1068 }
1069 
1070 /*----------------------------------------------------------------------------*/
1097 /*----------------------------------------------------------------------------*/
1098 
1099 inline static void
1101  const cs_real_t relaxp,
1102  const cs_real_t weight,
1103  const cs_real_3_t cell_ceni,
1104  const cs_real_3_t cell_cenj,
1105  const cs_real_3_t i_face_cog,
1106  const cs_real_3_t dijpf,
1107  const cs_real_63_t gradi,
1108  const cs_real_63_t gradj,
1109  const cs_real_6_t pi,
1110  const cs_real_6_t pj,
1111  const cs_real_6_t pia,
1112  const cs_real_6_t pja,
1113  cs_real_t pifri[6],
1114  cs_real_t pifrj[6],
1115  cs_real_t pjfri[6],
1116  cs_real_t pjfrj[6],
1117  cs_real_t pip[6],
1118  cs_real_t pjp[6],
1119  cs_real_t pipr[6],
1120  cs_real_t pjpr[6])
1121 {
1122  cs_real_6_t pir, pjr;
1123  cs_real_6_t recoi, recoj;
1124 
1126  weight,
1127  cell_ceni,
1128  cell_cenj,
1129  i_face_cog,
1130  dijpf,
1131  gradi,
1132  gradj,
1133  pi,
1134  pj,
1135  recoi,
1136  recoj,
1137  pip,
1138  pjp);
1139 
1140  cs_i_relax_c_val_tensor(relaxp,
1141  pia,
1142  pja,
1143  recoi,
1144  recoj,
1145  pi,
1146  pj,
1147  pir,
1148  pjr,
1149  pipr,
1150  pjpr);
1151 
1153  pj,
1154  pir,
1155  pjr,
1156  pifri,
1157  pifrj,
1158  pjfri,
1159  pjfrj);
1160 }
1161 
1162 /*----------------------------------------------------------------------------*/
1186 /*----------------------------------------------------------------------------*/
1187 
1188 inline static void
1189 cs_i_cd_unsteady_upwind(const int ircflp,
1190  const cs_real_t weight,
1191  const cs_real_3_t cell_ceni,
1192  const cs_real_3_t cell_cenj,
1193  const cs_real_3_t i_face_cog,
1194  const cs_real_3_t dijpf,
1195  const cs_real_3_t gradi,
1196  const cs_real_3_t gradj,
1197  const cs_real_t pi,
1198  const cs_real_t pj,
1199  cs_real_t *pifri,
1200  cs_real_t *pifrj,
1201  cs_real_t *pjfri,
1202  cs_real_t *pjfrj,
1203  cs_real_t *pip,
1204  cs_real_t *pjp,
1205  cs_real_t *pipr,
1206  cs_real_t *pjpr)
1207 {
1208  cs_real_t pir, pjr;
1209  cs_real_t recoi, recoj;
1210 
1211  cs_i_compute_quantities(ircflp,
1212  weight,
1213  cell_ceni,
1214  cell_cenj,
1215  i_face_cog,
1216  dijpf,
1217  gradi,
1218  gradj,
1219  pi,
1220  pj,
1221  &recoi,
1222  &recoj,
1223  pip,
1224  pjp);
1225 
1227  pj,
1228  *pip,
1229  *pjp,
1230  &pir,
1231  &pjr,
1232  pipr,
1233  pjpr);
1234 
1235  cs_upwind_f_val(pi,
1236  pj,
1237  pir,
1238  pjr,
1239  pifri,
1240  pifrj,
1241  pjfri,
1242  pjfrj);
1243 }
1244 
1245 /*----------------------------------------------------------------------------*/
1269 /*----------------------------------------------------------------------------*/
1270 inline static void
1272  const cs_real_t weight,
1273  const cs_real_3_t cell_ceni,
1274  const cs_real_3_t cell_cenj,
1275  const cs_real_3_t i_face_cog,
1276  const cs_real_3_t dijpf,
1277  const cs_real_63_t gradi,
1278  const cs_real_63_t gradj,
1279  const cs_real_6_t pi,
1280  const cs_real_6_t pj,
1281  cs_real_t pifri[6],
1282  cs_real_t pifrj[6],
1283  cs_real_t pjfri[6],
1284  cs_real_t pjfrj[6],
1285  cs_real_t pip[6],
1286  cs_real_t pjp[6],
1287  cs_real_t pipr[6],
1288  cs_real_t pjpr[6])
1289 {
1290  cs_real_6_t pir, pjr;
1291  cs_real_6_t recoi, recoj;
1292 
1294  weight,
1295  cell_ceni,
1296  cell_cenj,
1297  i_face_cog,
1298  dijpf,
1299  gradi,
1300  gradj,
1301  pi,
1302  pj,
1303  recoi,
1304  recoj,
1305  pip,
1306  pjp);
1307 
1309  pj,
1310  pip,
1311  pjp,
1312  pir,
1313  pjr,
1314  pipr,
1315  pjpr);
1316 
1318  pj,
1319  pir,
1320  pjr,
1321  pifri,
1322  pifrj,
1323  pjfri,
1324  pjfrj);
1325 }
1326 
1327 /*----------------------------------------------------------------------------*/
1357 /*----------------------------------------------------------------------------*/
1358 
1359 inline static void
1360 cs_i_cd_steady(const int ircflp,
1361  const int ischcp,
1362  const double relaxp,
1363  const double blencp,
1364  const cs_real_t weight,
1365  const cs_real_3_t cell_ceni,
1366  const cs_real_3_t cell_cenj,
1367  const cs_real_3_t i_face_cog,
1368  const cs_real_3_t dijpf,
1369  const cs_real_3_t gradi,
1370  const cs_real_3_t gradj,
1371  const cs_real_t pi,
1372  const cs_real_t pj,
1373  const cs_real_t pia,
1374  const cs_real_t pja,
1375  cs_real_t *pifri,
1376  cs_real_t *pifrj,
1377  cs_real_t *pjfri,
1378  cs_real_t *pjfrj,
1379  cs_real_t *pip,
1380  cs_real_t *pjp,
1381  cs_real_t *pipr,
1382  cs_real_t *pjpr)
1383 {
1384  cs_real_t pir, pjr;
1385  cs_real_t recoi, recoj;
1386 
1387  cs_i_compute_quantities(ircflp,
1388  weight,
1389  cell_ceni,
1390  cell_cenj,
1391  i_face_cog,
1392  dijpf,
1393  gradi,
1394  gradj,
1395  pi,
1396  pj,
1397  &recoi,
1398  &recoj,
1399  pip,
1400  pjp);
1401 
1402  cs_i_relax_c_val(relaxp,
1403  pia,
1404  pja,
1405  recoi,
1406  recoj,
1407  pi,
1408  pj,
1409  &pir,
1410  &pjr,
1411  pipr,
1412  pjpr);
1413 
1414  if (ischcp == 1) {
1415 
1416  /* Centered
1417  --------*/
1418 
1419  cs_centered_f_val(weight,
1420  *pip,
1421  *pjp,
1422  *pipr,
1423  *pjpr,
1424  pifri,
1425  pifrj,
1426  pjfri,
1427  pjfrj);
1428 
1429  } else {
1430 
1431  /* Second order
1432  ------------*/
1433 
1434  cs_solu_f_val(cell_ceni,
1435  cell_cenj,
1436  i_face_cog,
1437  gradi,
1438  gradj,
1439  pi,
1440  pj,
1441  pir,
1442  pjr,
1443  pifri,
1444  pifrj,
1445  pjfri,
1446  pjfrj);
1447 
1448  }
1449 
1450  /* Blending
1451  --------*/
1452 
1453  cs_blend_f_val(blencp,
1454  pi,
1455  pj,
1456  pir,
1457  pjr,
1458  pifri,
1459  pifrj,
1460  pjfri,
1461  pjfrj);
1462 }
1463 
1464 /*----------------------------------------------------------------------------*/
1494 /*----------------------------------------------------------------------------*/
1495 
1496 inline static void
1497 cs_i_cd_steady_tensor(const int ircflp,
1498  const int ischcp,
1499  const double relaxp,
1500  const double blencp,
1501  const cs_real_t weight,
1502  const cs_real_3_t cell_ceni,
1503  const cs_real_3_t cell_cenj,
1504  const cs_real_3_t i_face_cog,
1505  const cs_real_3_t dijpf,
1506  const cs_real_63_t gradi,
1507  const cs_real_63_t gradj,
1508  const cs_real_6_t pi,
1509  const cs_real_6_t pj,
1510  const cs_real_6_t pia,
1511  const cs_real_6_t pja,
1512  cs_real_t pifri[6],
1513  cs_real_t pifrj[6],
1514  cs_real_t pjfri[6],
1515  cs_real_t pjfrj[6],
1516  cs_real_t pip[6],
1517  cs_real_t pjp[6],
1518  cs_real_t pipr[6],
1519  cs_real_t pjpr[6])
1520 
1521 {
1522  cs_real_6_t pir, pjr;
1523  cs_real_6_t recoi, recoj;
1524 
1526  weight,
1527  cell_ceni,
1528  cell_cenj,
1529  i_face_cog,
1530  dijpf,
1531  gradi,
1532  gradj,
1533  pi,
1534  pj,
1535  recoi,
1536  recoj,
1537  pip,
1538  pjp);
1539 
1540  cs_i_relax_c_val_tensor(relaxp,
1541  pia,
1542  pja,
1543  recoi,
1544  recoj,
1545  pi,
1546  pj,
1547  pir,
1548  pjr,
1549  pipr,
1550  pjpr);
1551 
1552  if (ischcp == 1) {
1553 
1554  /* Centered
1555  --------*/
1556 
1557  cs_centered_f_val_tensor(weight,
1558  pip,
1559  pjp,
1560  pipr,
1561  pjpr,
1562  pifri,
1563  pifrj,
1564  pjfri,
1565  pjfrj);
1566 
1567  } else {
1568 
1569  /* Second order
1570  ------------*/
1571 
1572  cs_solu_f_val_tensor(cell_ceni,
1573  cell_cenj,
1574  i_face_cog,
1575  gradi,
1576  gradj,
1577  pi,
1578  pj,
1579  pir,
1580  pjr,
1581  pifri,
1582  pifrj,
1583  pjfri,
1584  pjfrj);
1585 
1586  }
1587 
1588  /* Blending
1589  --------*/
1590 
1591  cs_blend_f_val_tensor(blencp,
1592  pi,
1593  pj,
1594  pir,
1595  pjr,
1596  pifri,
1597  pifrj,
1598  pjfri,
1599  pjfrj);
1600 }
1601 
1602 
1603 /*----------------------------------------------------------------------------*/
1630 /*----------------------------------------------------------------------------*/
1631 
1632 inline static void
1633 cs_i_cd_unsteady(const int ircflp,
1634  const int ischcp,
1635  const double blencp,
1636  const cs_real_t weight,
1637  const cs_real_3_t cell_ceni,
1638  const cs_real_3_t cell_cenj,
1639  const cs_real_3_t i_face_cog,
1640  const cs_real_3_t dijpf,
1641  const cs_real_3_t gradi,
1642  const cs_real_3_t gradj,
1643  const cs_real_t pi,
1644  const cs_real_t pj,
1645  cs_real_t *pifri,
1646  cs_real_t *pifrj,
1647  cs_real_t *pjfri,
1648  cs_real_t *pjfrj,
1649  cs_real_t *pip,
1650  cs_real_t *pjp,
1651  cs_real_t *pipr,
1652  cs_real_t *pjpr)
1653 {
1654  cs_real_t pir, pjr;
1655  cs_real_t recoi, recoj;
1656 
1657  cs_i_compute_quantities(ircflp,
1658  weight,
1659  cell_ceni,
1660  cell_cenj,
1661  i_face_cog,
1662  dijpf,
1663  gradi,
1664  gradj,
1665  pi,
1666  pj,
1667  &recoi,
1668  &recoj,
1669  pip,
1670  pjp);
1671 
1673  pj,
1674  *pip,
1675  *pjp,
1676  &pir,
1677  &pjr,
1678  pipr,
1679  pjpr);
1680 
1681  if (ischcp == 1) {
1682 
1683  /* Centered
1684  --------*/
1685 
1686  cs_centered_f_val(weight,
1687  *pip,
1688  *pjp,
1689  *pipr,
1690  *pjpr,
1691  pifri,
1692  pifrj,
1693  pjfri,
1694  pjfrj);
1695 
1696  } else {
1697 
1698  /* Second order
1699  ------------*/
1700 
1701  cs_solu_f_val(cell_ceni,
1702  cell_cenj,
1703  i_face_cog,
1704  gradi,
1705  gradj,
1706  pi,
1707  pj,
1708  pir,
1709  pjr,
1710  pifri,
1711  pifrj,
1712  pjfri,
1713  pjfrj);
1714 
1715  }
1716 
1717  /* Blending
1718  --------*/
1719 
1720  cs_blend_f_val(blencp,
1721  pi,
1722  pj,
1723  pir,
1724  pjr,
1725  pifri,
1726  pifrj,
1727  pjfri,
1728  pjfrj);
1729 }
1730 
1731 /*----------------------------------------------------------------------------*/
1758 /*----------------------------------------------------------------------------*/
1759 
1760 inline static void
1761 cs_i_cd_unsteady_tensor(const int ircflp,
1762  const int ischcp,
1763  const double blencp,
1764  const cs_real_t weight,
1765  const cs_real_3_t cell_ceni,
1766  const cs_real_3_t cell_cenj,
1767  const cs_real_3_t i_face_cog,
1768  const cs_real_3_t dijpf,
1769  const cs_real_63_t gradi,
1770  const cs_real_63_t gradj,
1771  const cs_real_6_t pi,
1772  const cs_real_6_t pj,
1773  cs_real_t pifri[6],
1774  cs_real_t pifrj[6],
1775  cs_real_t pjfri[6],
1776  cs_real_t pjfrj[6],
1777  cs_real_t pip[6],
1778  cs_real_t pjp[6],
1779  cs_real_t pipr[6],
1780  cs_real_t pjpr[6])
1781 
1782 {
1783  cs_real_6_t pir, pjr;
1784  cs_real_6_t recoi, recoj;
1785 
1787  weight,
1788  cell_ceni,
1789  cell_cenj,
1790  i_face_cog,
1791  dijpf,
1792  gradi,
1793  gradj,
1794  pi,
1795  pj,
1796  recoi,
1797  recoj,
1798  pip,
1799  pjp);
1800 
1802  pj,
1803  pip,
1804  pjp,
1805  pir,
1806  pjr,
1807  pipr,
1808  pjpr);
1809 
1810  if (ischcp == 1) {
1811 
1812  /* Centered
1813  --------*/
1814 
1815  cs_centered_f_val_tensor(weight,
1816  pip,
1817  pjp,
1818  pipr,
1819  pjpr,
1820  pifri,
1821  pifrj,
1822  pjfri,
1823  pjfrj);
1824 
1825  } else {
1826 
1827  /* Second order
1828  ------------*/
1829 
1830  cs_solu_f_val_tensor(cell_ceni,
1831  cell_cenj,
1832  i_face_cog,
1833  gradi,
1834  gradj,
1835  pi,
1836  pj,
1837  pir,
1838  pjr,
1839  pifri,
1840  pifrj,
1841  pjfri,
1842  pjfrj);
1843 
1844  }
1845 
1846  /* Blending
1847  --------*/
1848 
1849  cs_blend_f_val_tensor(blencp,
1850  pi,
1851  pj,
1852  pir,
1853  pjr,
1854  pifri,
1855  pifrj,
1856  pjfri,
1857  pjfrj);
1858 }
1859 /*----------------------------------------------------------------------------*/
1896 /*----------------------------------------------------------------------------*/
1897 
1898 inline static void
1899 cs_i_cd_steady_slope_test(bool *upwind_switch,
1900  const int ircflp,
1901  const int ischcp,
1902  const double relaxp,
1903  const double blencp,
1904  const cs_real_t weight,
1905  const cs_real_t i_dist,
1906  const cs_real_t i_face_surf,
1907  const cs_real_3_t cell_ceni,
1908  const cs_real_3_t cell_cenj,
1909  const cs_real_3_t i_face_normal,
1910  const cs_real_3_t i_face_cog,
1911  const cs_real_3_t dijpf,
1912  const cs_real_t i_massflux,
1913  const cs_real_3_t gradi,
1914  const cs_real_3_t gradj,
1915  const cs_real_3_t grdpai,
1916  const cs_real_3_t grdpaj,
1917  const cs_real_t pi,
1918  const cs_real_t pj,
1919  const cs_real_t pia,
1920  const cs_real_t pja,
1921  cs_real_t *pifri,
1922  cs_real_t *pifrj,
1923  cs_real_t *pjfri,
1924  cs_real_t *pjfrj,
1925  cs_real_t *pip,
1926  cs_real_t *pjp,
1927  cs_real_t *pipr,
1928  cs_real_t *pjpr)
1929 {
1930  cs_real_t pir, pjr;
1931  cs_real_t recoi, recoj;
1932  cs_real_t distf, srfan, testij, tesqck;
1933 
1934  distf = i_dist;
1935  srfan = i_face_surf;
1936 
1937  *upwind_switch = false;
1938 
1939  cs_i_compute_quantities(ircflp,
1940  weight,
1941  cell_ceni,
1942  cell_cenj,
1943  i_face_cog,
1944  dijpf,
1945  gradi,
1946  gradj,
1947  pi,
1948  pj,
1949  &recoi,
1950  &recoj,
1951  pip,
1952  pjp);
1953 
1954  cs_i_relax_c_val(relaxp,
1955  pia,
1956  pja,
1957  recoi,
1958  recoj,
1959  pi,
1960  pj,
1961  &pir,
1962  &pjr,
1963  pipr,
1964  pjpr);
1965 
1966  cs_slope_test(pi,
1967  pj,
1968  distf,
1969  srfan,
1970  i_face_normal,
1971  gradi,
1972  gradj,
1973  grdpai,
1974  grdpaj,
1975  i_massflux,
1976  &testij,
1977  &tesqck);
1978 
1979  if (tesqck<=0. || testij<=0.) {
1980 
1981  /* Upwind
1982  --------*/
1983 
1984  cs_upwind_f_val(pi,
1985  pj,
1986  pir,
1987  pjr,
1988  pifri,
1989  pifrj,
1990  pjfri,
1991  pjfrj);
1992 
1993  *upwind_switch = true;
1994 
1995  } else {
1996 
1997  if (ischcp==1) {
1998 
1999  /* Centered
2000  --------*/
2001 
2002  cs_centered_f_val(weight,
2003  *pip,
2004  *pjp,
2005  *pipr,
2006  *pjpr,
2007  pifri,
2008  pifrj,
2009  pjfri,
2010  pjfrj);
2011 
2012  } else {
2013 
2014  /* Second order
2015  ------------*/
2016 
2017  cs_solu_f_val(cell_ceni,
2018  cell_cenj,
2019  i_face_cog,
2020  gradi,
2021  gradj,
2022  pi,
2023  pj,
2024  pir,
2025  pjr,
2026  pifri,
2027  pifrj,
2028  pjfri,
2029  pjfrj);
2030 
2031  }
2032  }
2033 
2034  /* Blending
2035  --------*/
2036 
2037  cs_blend_f_val(blencp,
2038  pi,
2039  pj,
2040  pir,
2041  pjr,
2042  pifri,
2043  pifrj,
2044  pjfri,
2045  pjfrj);
2046 }
2047 /*----------------------------------------------------------------------------*/
2084 /*----------------------------------------------------------------------------*/
2085 
2086 inline static void
2087 cs_i_cd_steady_slope_test_tensor(bool upwind_switch[6],
2088  const int ircflp,
2089  const int ischcp,
2090  const double relaxp,
2091  const double blencp,
2092  const cs_real_t weight,
2093  const cs_real_t i_dist,
2094  const cs_real_t i_face_surf,
2095  const cs_real_3_t cell_ceni,
2096  const cs_real_3_t cell_cenj,
2097  const cs_real_3_t i_face_normal,
2098  const cs_real_3_t i_face_cog,
2099  const cs_real_3_t dijpf,
2100  const cs_real_t i_massflux,
2101  const cs_real_63_t gradi,
2102  const cs_real_63_t gradj,
2103  const cs_real_63_t grdpai,
2104  const cs_real_63_t grdpaj,
2105  const cs_real_6_t pi,
2106  const cs_real_6_t pj,
2107  const cs_real_6_t pia,
2108  const cs_real_6_t pja,
2109  cs_real_t pifri[6],
2110  cs_real_t pifrj[6],
2111  cs_real_t pjfri[6],
2112  cs_real_t pjfrj[6],
2113  cs_real_t pip[6],
2114  cs_real_t pjp[6],
2115  cs_real_t pipr[6],
2116  cs_real_t pjpr[6])
2117 {
2118  cs_real_6_t pir, pjr;
2119  cs_real_6_t recoi, recoj;
2120  cs_real_t distf, srfan;
2121  cs_real_6_t testij, tesqck;
2122  int isou;
2123 
2124  distf = i_dist;
2125  srfan = i_face_surf;
2126 
2127  for (isou = 0; isou < 6; isou++) {
2128  upwind_switch[isou] = false;
2129  }
2130 
2132  weight,
2133  cell_ceni,
2134  cell_cenj,
2135  i_face_cog,
2136  dijpf,
2137  gradi,
2138  gradj,
2139  pi,
2140  pj,
2141  recoi,
2142  recoj,
2143  pip,
2144  pjp);
2145 
2146  cs_i_relax_c_val_tensor(relaxp,
2147  pia,
2148  pja,
2149  recoi,
2150  recoj,
2151  pi,
2152  pj,
2153  pir,
2154  pjr,
2155  pipr,
2156  pjpr);
2157 
2159  pj,
2160  distf,
2161  srfan,
2162  i_face_normal,
2163  gradi,
2164  gradj,
2165  grdpai,
2166  grdpaj,
2167  i_massflux,
2168  testij,
2169  tesqck);
2170 
2171  for (isou = 0 ; isou < 6; isou++) {
2172  if (tesqck[isou]<=0. || testij[isou]<=0.) {
2173 
2174  /* Upwind
2175  --------*/
2176 
2177  cs_upwind_f_val(pi[isou],
2178  pj[isou],
2179  pir[isou],
2180  pjr[isou],
2181  &pifri[isou],
2182  &pifrj[isou],
2183  &pjfri[isou],
2184  &pjfrj[isou]);
2185 
2186  upwind_switch[isou] = true;
2187 
2188  } else {
2189 
2190  if (ischcp==1) {
2191 
2192  /* Centered
2193  --------*/
2194 
2195  cs_centered_f_val(weight,
2196  pip[isou],
2197  pjp[isou],
2198  pipr[isou],
2199  pjpr[isou],
2200  &pifri[isou],
2201  &pifrj[isou],
2202  &pjfri[isou],
2203  &pjfrj[isou]);
2204 
2205  } else {
2206 
2207  /* Second order
2208  ------------*/
2209 
2210  cs_solu_f_val(cell_ceni,
2211  cell_cenj,
2212  i_face_cog,
2213  gradi[isou],
2214  gradj[isou],
2215  pi[isou],
2216  pj[isou],
2217  pir[isou],
2218  pjr[isou],
2219  &pifri[isou],
2220  &pifrj[isou],
2221  &pjfri[isou],
2222  &pjfrj[isou]);
2223 
2224  }
2225  }
2226  }
2227 
2228  /* Blending
2229  --------*/
2230 
2231  cs_blend_f_val_tensor(blencp,
2232  pi,
2233  pj,
2234  pir,
2235  pjr,
2236  pifri,
2237  pifrj,
2238  pjfri,
2239  pjfrj);
2240 }
2241 
2242 /*----------------------------------------------------------------------------*/
2276 /*----------------------------------------------------------------------------*/
2277 
2278 inline static void
2279 cs_i_cd_unsteady_slope_test(bool *upwind_switch,
2280  const int ircflp,
2281  const int ischcp,
2282  const double blencp,
2283  const cs_real_t weight,
2284  const cs_real_t i_dist,
2285  const cs_real_t i_face_surf,
2286  const cs_real_3_t cell_ceni,
2287  const cs_real_3_t cell_cenj,
2288  const cs_real_3_t i_face_normal,
2289  const cs_real_3_t i_face_cog,
2290  const cs_real_3_t dijpf,
2291  const cs_real_t i_massflux,
2292  const cs_real_3_t gradi,
2293  const cs_real_3_t gradj,
2294  const cs_real_3_t grdpai,
2295  const cs_real_3_t grdpaj,
2296  const cs_real_t pi,
2297  const cs_real_t pj,
2298  cs_real_t *pifri,
2299  cs_real_t *pifrj,
2300  cs_real_t *pjfri,
2301  cs_real_t *pjfrj,
2302  cs_real_t *pip,
2303  cs_real_t *pjp,
2304  cs_real_t *pipr,
2305  cs_real_t *pjpr)
2306 {
2307  cs_real_t pir, pjr;
2308  cs_real_t recoi, recoj;
2309  cs_real_t distf, srfan, testij, tesqck;
2310 
2311  distf = i_dist;
2312  srfan = i_face_surf;
2313 
2314  *upwind_switch = false;
2315 
2316  cs_i_compute_quantities(ircflp,
2317  weight,
2318  cell_ceni,
2319  cell_cenj,
2320  i_face_cog,
2321  dijpf,
2322  gradi,
2323  gradj,
2324  pi,
2325  pj,
2326  &recoi,
2327  &recoj,
2328  pip,
2329  pjp);
2330 
2332  pj,
2333  *pip,
2334  *pjp,
2335  &pir,
2336  &pjr,
2337  pipr,
2338  pjpr);
2339 
2340  cs_slope_test(pi,
2341  pj,
2342  distf,
2343  srfan,
2344  i_face_normal,
2345  gradi,
2346  gradj,
2347  grdpai,
2348  grdpaj,
2349  i_massflux,
2350  &testij,
2351  &tesqck);
2352 
2353  if (tesqck<=0. || testij<=0.) {
2354 
2355  /* Upwind
2356  --------*/
2357 
2358  cs_upwind_f_val(pi,
2359  pj,
2360  pir,
2361  pjr,
2362  pifri,
2363  pifrj,
2364  pjfri,
2365  pjfrj);
2366 
2367  *upwind_switch = true;
2368 
2369  } else {
2370 
2371  if (ischcp==1) {
2372 
2373  /* Centered
2374  --------*/
2375 
2376  cs_centered_f_val(weight,
2377  *pip,
2378  *pjp,
2379  *pipr,
2380  *pjpr,
2381  pifri,
2382  pifrj,
2383  pjfri,
2384  pjfrj);
2385 
2386  } else {
2387 
2388  /* Second order
2389  ------------*/
2390 
2391  cs_solu_f_val(cell_ceni,
2392  cell_cenj,
2393  i_face_cog,
2394  gradi,
2395  gradj,
2396  pi,
2397  pj,
2398  pir,
2399  pjr,
2400  pifri,
2401  pifrj,
2402  pjfri,
2403  pjfrj);
2404 
2405  }
2406  }
2407 
2408  /* Blending
2409  --------*/
2410 
2411  cs_blend_f_val(blencp,
2412  pi,
2413  pj,
2414  pir,
2415  pjr,
2416  pifri,
2417  pifrj,
2418  pjfri,
2419  pjfrj);
2420 }
2421 
2422 /*----------------------------------------------------------------------------*/
2456 /*----------------------------------------------------------------------------*/
2457 
2458 inline static void
2460  const int ircflp,
2461  const int ischcp,
2462  const double blencp,
2463  const cs_real_t weight,
2464  const cs_real_t i_dist,
2465  const cs_real_t i_face_surf,
2466  const cs_real_3_t cell_ceni,
2467  const cs_real_3_t cell_cenj,
2468  const cs_real_3_t i_face_normal,
2469  const cs_real_3_t i_face_cog,
2470  const cs_real_3_t dijpf,
2471  const cs_real_t i_massflux,
2472  const cs_real_63_t gradi,
2473  const cs_real_63_t gradj,
2474  const cs_real_63_t grdpai,
2475  const cs_real_63_t grdpaj,
2476  const cs_real_6_t pi,
2477  const cs_real_6_t pj,
2478  cs_real_t pifri[6],
2479  cs_real_t pifrj[6],
2480  cs_real_t pjfri[6],
2481  cs_real_t pjfrj[6],
2482  cs_real_t pip[6],
2483  cs_real_t pjp[6],
2484  cs_real_t pipr[6],
2485  cs_real_t pjpr[6])
2486 {
2487  cs_real_6_t pir, pjr;
2488  cs_real_6_t recoi, recoj;
2489  cs_real_t distf, srfan;
2490  cs_real_6_t testij, tesqck;
2491  int isou;
2492 
2493  distf = i_dist;
2494  srfan = i_face_surf;
2495 
2496  for (isou = 0; isou < 6; isou++) {
2497  upwind_switch[isou] = false;
2498  }
2499 
2501  weight,
2502  cell_ceni,
2503  cell_cenj,
2504  i_face_cog,
2505  dijpf,
2506  gradi,
2507  gradj,
2508  pi,
2509  pj,
2510  recoi,
2511  recoj,
2512  pip,
2513  pjp);
2514 
2516  pj,
2517  pip,
2518  pjp,
2519  pir,
2520  pjr,
2521  pipr,
2522  pjpr);
2523 
2525  pj,
2526  distf,
2527  srfan,
2528  i_face_normal,
2529  gradi,
2530  gradj,
2531  grdpai,
2532  grdpaj,
2533  i_massflux,
2534  testij,
2535  tesqck);
2536  for (isou = 0 ; isou < 6; isou++) {
2537  if (tesqck[isou]<=0. || testij[isou]<=0.) {
2538 
2539  /* Upwind
2540  --------*/
2541 
2542  cs_upwind_f_val(pi[isou],
2543  pj[isou],
2544  pir[isou],
2545  pjr[isou],
2546  &pifri[isou],
2547  &pifrj[isou],
2548  &pjfri[isou],
2549  &pjfrj[isou]);
2550 
2551  upwind_switch[isou] = true;
2552 
2553  } else {
2554 
2555  if (ischcp==1) {
2556 
2557  /* Centered
2558  --------*/
2559 
2560  cs_centered_f_val(weight,
2561  pip[isou],
2562  pjp[isou],
2563  pipr[isou],
2564  pjpr[isou],
2565  &pifri[isou],
2566  &pifrj[isou],
2567  &pjfri[isou],
2568  &pjfrj[isou]);
2569 
2570  } else {
2571 
2572  /* Second order
2573  ------------*/
2574 
2575  cs_solu_f_val(cell_ceni,
2576  cell_cenj,
2577  i_face_cog,
2578  gradi[isou],
2579  gradj[isou],
2580  pi[isou],
2581  pj[isou],
2582  pir[isou],
2583  pjr[isou],
2584  &pifri[isou],
2585  &pifrj[isou],
2586  &pjfri[isou],
2587  &pjfrj[isou]);
2588 
2589  }
2590  }
2591  }
2592 
2593  /* Blending
2594  --------*/
2595 
2596  cs_blend_f_val_tensor(blencp,
2597  pi,
2598  pj,
2599  pir,
2600  pjr,
2601  pifri,
2602  pifrj,
2603  pjfri,
2604  pjfrj);
2605 }
2606 
2607 
2608 /*----------------------------------------------------------------------------*/
2617 /*----------------------------------------------------------------------------*/
2618 
2619 inline static void
2621  const cs_real_3_t gradi,
2622  const int ircflp,
2623  cs_real_t *recoi)
2624 {
2625  *recoi = ircflp * ( gradi[0]*diipb[0]
2626  + gradi[1]*diipb[1]
2627  + gradi[2]*diipb[2]);
2628 }
2629 
2630 /*----------------------------------------------------------------------------*/
2639 /*----------------------------------------------------------------------------*/
2640 
2641 inline static void
2643  const cs_real_63_t gradi,
2644  const int ircflp,
2645  cs_real_t recoi[6])
2646 {
2647  for (int isou = 0; isou < 6; isou++) {
2648  recoi[isou] = ircflp * (gradi[isou][0]*diipb[0]
2649  + gradi[isou][1]*diipb[1]
2650  + gradi[isou][2]*diipb[2]);
2651  }
2652 }
2653 /*----------------------------------------------------------------------------*/
2664 /*----------------------------------------------------------------------------*/
2665 
2666 inline static void
2667 cs_b_relax_c_val(const double relaxp,
2668  const cs_real_t pi,
2669  const cs_real_t pia,
2670  const cs_real_t recoi,
2671  cs_real_t *pir,
2672  cs_real_t *pipr)
2673 {
2674  *pir = pi/relaxp - (1.-relaxp)/relaxp*pia;
2675  *pipr = *pir + recoi;
2676 }
2677 /*----------------------------------------------------------------------------*/
2688 /*----------------------------------------------------------------------------*/
2689 
2690 inline static void
2691 cs_b_relax_c_val_tensor(const double relaxp,
2692  const cs_real_6_t pi,
2693  const cs_real_6_t pia,
2694  const cs_real_6_t recoi,
2695  cs_real_t pir[6],
2696  cs_real_t pipr[6])
2697 {
2698  for (int isou = 0; isou < 6; isou++) {
2699  pir[isou] = pi[isou]/relaxp - (1.-relaxp)/relaxp*pia[isou];
2700  pipr[isou] = pir[isou] + recoi[isou];
2701  }
2702 }
2703 /*----------------------------------------------------------------------------*/
2712 /*----------------------------------------------------------------------------*/
2713 
2714 inline static void
2716  const cs_real_t recoi,
2717  cs_real_t *pir,
2718  cs_real_t *pipr)
2719 {
2720  *pir = pi;
2721  *pipr = pi + recoi;
2722 }
2723 /*----------------------------------------------------------------------------*/
2732 /*----------------------------------------------------------------------------*/
2733 
2734 inline static void
2736  const cs_real_6_t recoi,
2737  cs_real_t pir[6],
2738  cs_real_t pipr[6])
2739 {
2740  for(int isou = 0; isou< 6; isou++) {
2741  pir[isou] = pi[isou];
2742  pipr[isou] = pi[isou] + recoi[isou];
2743  }
2744 }
2745 
2746 
2747 /*----------------------------------------------------------------------------*/
2770 /*----------------------------------------------------------------------------*/
2771 
2772 inline static void
2773 cs_b_imposed_conv_flux(const int iconvp,
2774  const int inc,
2775  const int ifaccp,
2776  const cs_int_t bc_type,
2777  const int icvfli,
2778  const cs_real_t pi,
2779  const cs_real_t pir,
2780  const cs_real_t pipr,
2781  const cs_real_t coefap,
2782  const cs_real_t coefbp,
2783  const cs_real_t coface,
2784  const cs_real_t cofbce,
2785  const cs_real_t b_massflux,
2786  cs_real_t *flux)
2787 {
2788  cs_real_t flui, fluj, pfac;
2789 
2790  /* Computed convective flux */
2791 
2792  if (icvfli == 0) {
2793 
2794  /* Remove decentering for coupled faces */
2795  if (ifaccp==1 && bc_type==CS_COUPLED) {
2796  flui = 0.0;
2797  fluj = b_massflux;
2798  } else {
2799  flui = 0.5*(b_massflux +fabs(b_massflux));
2800  fluj = 0.5*(b_massflux -fabs(b_massflux));
2801  }
2802 
2803  pfac = inc*coefap + coefbp*pipr;
2804  *flux += iconvp*(flui*pir + fluj*pfac - b_massflux*pi);
2805 
2806  /* Imposed convective flux */
2807 
2808  } else {
2809 
2810  pfac = inc*coface + cofbce*pipr;
2811  *flux += iconvp*(-b_massflux*pi + pfac);
2812 
2813  }
2814 }
2815 
2816 /*----------------------------------------------------------------------------*/
2839 /*----------------------------------------------------------------------------*/
2840 
2841 inline static void
2843  const int inc,
2844  const int ifaccp,
2845  const cs_int_t bc_type,
2846  const int icvfli,
2847  const cs_real_t pir,
2848  const cs_real_t pipr,
2849  const cs_real_t coefap,
2850  const cs_real_t coefbp,
2851  const cs_real_t coface,
2852  const cs_real_t cofbce,
2853  const cs_real_t xcpp,
2854  const cs_real_t b_massflux,
2855  cs_real_t *flux)
2856 {
2857  cs_real_t flui, fluj, pfac;
2858 
2859  /* Computed convective flux */
2860 
2861  if (icvfli == 0) {
2862 
2863  /* Remove decentering for coupled faces */
2864  if (ifaccp==1 && bc_type==CS_COUPLED) {
2865  flui = 0.0;
2866  fluj = b_massflux;
2867  } else {
2868  flui = 0.5*(b_massflux +fabs(b_massflux));
2869  fluj = 0.5*(b_massflux -fabs(b_massflux));
2870  }
2871 
2872  pfac = inc*coefap + coefbp*pipr;
2873  *flux += iconvp*xcpp*(flui*pir + fluj*pfac);
2874 
2875  /* Imposed convective flux */
2876 
2877  } else {
2878 
2879  pfac = inc*coface + cofbce*pipr;
2880  *flux += iconvp*pfac;
2881 
2882  }
2883 }
2884 
2885 /*----------------------------------------------------------------------------*/
2904 /*----------------------------------------------------------------------------*/
2905 
2906 inline static void
2907 cs_b_upwind_flux(const int iconvp,
2908  const int inc,
2909  const int ifaccp,
2910  const int bc_type,
2911  const cs_real_t pi,
2912  const cs_real_t pir,
2913  const cs_real_t pipr,
2914  const cs_real_t coefap,
2915  const cs_real_t coefbp,
2916  const cs_real_t b_massflux,
2917  cs_real_t *flux)
2918 {
2919  cs_real_t flui, fluj, pfac;
2920 
2921  /* Remove decentering for coupled faces */
2922  if (ifaccp==1 && bc_type==CS_COUPLED) {
2923  flui = 0.0;
2924  fluj = b_massflux;
2925  } else {
2926  flui = 0.5*(b_massflux +fabs(b_massflux));
2927  fluj = 0.5*(b_massflux -fabs(b_massflux));
2928  }
2929 
2930  pfac = inc*coefap + coefbp*pipr;
2931  *flux += iconvp*(flui*pir + fluj*pfac - b_massflux*pi);
2932 }
2933 
2934 /*----------------------------------------------------------------------------*/
2953 /*----------------------------------------------------------------------------*/
2954 
2955 inline static void
2956 cs_b_upwind_flux_tensor(const int iconvp,
2957  const int inc,
2958  const int ifaccp,
2959  const int bc_type,
2960  const cs_real_6_t pi,
2961  const cs_real_6_t pir,
2962  const cs_real_6_t pipr,
2963  const cs_real_6_t coefa,
2964  const cs_real_66_t coefb,
2965  const cs_real_t b_massflux,
2966  cs_real_t flux[6])
2967 {
2968  cs_real_t flui, fluj, pfac;
2969 
2970  /* Remove decentering for coupled faces */
2971  if (ifaccp == 1 && bc_type == CS_COUPLED) {
2972  flui = 0.0;
2973  fluj = b_massflux;
2974  } else {
2975  flui = 0.5*(b_massflux +fabs(b_massflux));
2976  fluj = 0.5*(b_massflux -fabs(b_massflux));
2977  }
2978  for (int isou = 0; isou < 6; isou++) {
2979  pfac = inc*coefa[isou];
2980  for (int jsou = 0; jsou < 6; jsou++) {
2981  pfac += coefb[isou][jsou]*pipr[jsou];
2982  }
2983  flux[isou] += iconvp*(flui*pir[isou]
2984  + fluj*pfac
2985  - b_massflux*pi[isou]);
2986  }
2987 }
2988 /*----------------------------------------------------------------------------*/
3010 /*----------------------------------------------------------------------------*/
3011 
3012 inline static void
3013 cs_b_upwind_flux_cons(const int iconvp,
3014  const int inc,
3015  const int ifaccp,
3016  const cs_int_t bc_type,
3017  const cs_real_t pir,
3018  const cs_real_t pipr,
3019  const cs_real_t coefap,
3020  const cs_real_t coefbp,
3021  const cs_real_t b_massflux,
3022  const cs_real_t xcpp,
3023  cs_real_t *flux)
3024 {
3025  cs_real_t flui, fluj, pfac;
3026 
3027  /* Remove decentering for coupled faces */
3028  if (ifaccp==1 && bc_type==CS_COUPLED) {
3029  flui = 0.0;
3030  fluj = b_massflux;
3031  } else {
3032  flui = 0.5*(b_massflux +fabs(b_massflux));
3033  fluj = 0.5*(b_massflux -fabs(b_massflux));
3034  }
3035 
3036  pfac = inc*coefap + coefbp*pipr;
3037  *flux += iconvp*xcpp*(flui*pir + fluj*pfac);
3038 }
3039 
3040 /*----------------------------------------------------------------------------*/
3052 /*----------------------------------------------------------------------------*/
3053 
3054 inline static void
3055 cs_b_diff_flux(const int idiffp,
3056  const int inc,
3057  const cs_real_t pipr,
3058  const cs_real_t cofafp,
3059  const cs_real_t cofbfp,
3060  const cs_real_t b_visc,
3061  cs_real_t *flux)
3062 {
3063  cs_real_t pfacd = inc*cofafp + cofbfp*pipr;
3064  *flux += idiffp*b_visc*pfacd;
3065 }
3066 /*----------------------------------------------------------------------------*/
3078 /*----------------------------------------------------------------------------*/
3079 
3080 inline static void
3081 cs_b_diff_flux_tensor(const int idiffp,
3082  const int inc,
3083  const cs_real_6_t pipr,
3084  const cs_real_6_t cofaf,
3085  const cs_real_66_t cofbf,
3086  const cs_real_t b_visc,
3087  cs_real_t flux[6])
3088 {
3089  cs_real_t pfacd ;
3090  for (int isou = 0; isou < 6; isou++) {
3091  pfacd = inc*cofaf[isou];
3092  for (int jsou = 0; jsou < 6; jsou++) {
3093  pfacd += cofbf[isou][jsou]*pipr[jsou];
3094  }
3095  flux[isou] += idiffp*b_visc*pfacd;
3096  }
3097 }
3098 /*----------------------------------------------------------------------------*/
3112 /*----------------------------------------------------------------------------*/
3113 
3114 inline static void
3115 cs_b_cd_steady(const int ircflp,
3116  const double relaxp,
3117  const cs_real_3_t diipb,
3118  const cs_real_3_t gradi,
3119  const cs_real_t pi,
3120  const cs_real_t pia,
3121  cs_real_t *pir,
3122  cs_real_t *pipr)
3123 {
3124  cs_real_t recoi;
3125 
3127  gradi,
3128  ircflp,
3129  &recoi);
3130 
3131  cs_b_relax_c_val(relaxp,
3132  pi,
3133  pia,
3134  recoi,
3135  pir,
3136  pipr);
3137 }
3138 
3139 /*----------------------------------------------------------------------------*/
3153 /*----------------------------------------------------------------------------*/
3154 
3155 inline static void
3156 cs_b_cd_steady_tensor(const int ircflp,
3157  const double relaxp,
3158  const cs_real_3_t diipb,
3159  const cs_real_63_t gradi,
3160  const cs_real_6_t pi,
3161  const cs_real_6_t pia,
3162  cs_real_t pir[6],
3163  cs_real_t pipr[6])
3164 {
3165  cs_real_6_t recoi;
3166 
3168  gradi,
3169  ircflp,
3170  recoi);
3171 
3172  cs_b_relax_c_val_tensor(relaxp,
3173  pi,
3174  pia,
3175  recoi,
3176  pir,
3177  pipr);
3178 }
3179 /*----------------------------------------------------------------------------*/
3191 /*----------------------------------------------------------------------------*/
3192 
3193 inline static void
3194 cs_b_cd_unsteady(const int ircflp,
3195  const cs_real_3_t diipb,
3196  const cs_real_3_t gradi,
3197  const cs_real_t pi,
3198  cs_real_t *pir,
3199  cs_real_t *pipr)
3200 {
3201  cs_real_t recoi;
3202 
3204  gradi,
3205  ircflp,
3206  &recoi);
3207 
3209  recoi,
3210  pir,
3211  pipr);
3212 }
3213 /*----------------------------------------------------------------------------*/
3225 /*----------------------------------------------------------------------------*/
3226 
3227 inline static void
3228 cs_b_cd_unsteady_tensor(const int ircflp,
3229  const cs_real_3_t diipb,
3230  const cs_real_63_t gradi,
3231  const cs_real_6_t pi,
3232  cs_real_t pir[6],
3233  cs_real_t pipr[6])
3234 {
3235  cs_real_6_t recoi;
3236 
3238  gradi,
3239  ircflp,
3240  recoi);
3241 
3243  recoi,
3244  pir,
3245  pipr);
3246 }
3247 /*----------------------------------------------------------------------------*/
3263 /*----------------------------------------------------------------------------*/
3264 
3265 inline static void
3266 cs_slope_test_gradient(const int f_id,
3267  const int inc,
3268  const cs_halo_type_t halo_type,
3269  cs_real_3_t *grad,
3270  cs_real_3_t *grdpa,
3271  cs_real_t *pvar,
3272  const cs_real_t *coefap,
3273  const cs_real_t *coefbp,
3274  const cs_real_t *i_massflux)
3275 {
3276  const cs_mesh_t *m = cs_glob_mesh;
3277  const cs_halo_t *halo = m->halo;
3279 
3280  const cs_lnum_t n_cells = m->n_cells;
3281 
3282  const cs_lnum_2_t *restrict i_face_cells
3283  = (const cs_lnum_2_t *restrict)m->i_face_cells;
3284  const cs_lnum_t *restrict b_face_cells
3285  = (const cs_lnum_t *restrict)m->b_face_cells;
3286  const cs_real_t *restrict cell_vol = fvq->cell_vol;
3287  const cs_real_3_t *restrict cell_cen
3288  = (const cs_real_3_t *restrict)fvq->cell_cen;
3289  const cs_real_3_t *restrict i_face_normal
3290  = (const cs_real_3_t *restrict)fvq->i_face_normal;
3291  const cs_real_3_t *restrict b_face_normal
3292  = (const cs_real_3_t *restrict)fvq->b_face_normal;
3293  const cs_real_3_t *restrict i_face_cog
3294  = (const cs_real_3_t *restrict)fvq->i_face_cog;
3295  const cs_real_3_t *restrict diipb
3296  = (const cs_real_3_t *restrict)fvq->diipb;
3297 
3298  const int n_i_groups = m->i_face_numbering->n_groups;
3299  const int n_i_threads = m->i_face_numbering->n_threads;
3300  const int n_b_groups = m->b_face_numbering->n_groups;
3301  const int n_b_threads = m->b_face_numbering->n_threads;
3302  const cs_lnum_t *restrict i_group_index = m->i_face_numbering->group_index;
3303  const cs_lnum_t *restrict b_group_index = m->b_face_numbering->group_index;
3304 
3305  cs_lnum_t cell_id, face_id, ii, jj;
3306  int g_id, t_id;
3307 
3308  cs_real_t difx,dify,difz,djfx,djfy,djfz;
3309  cs_real_t diipbx, diipby, diipbz;
3310  cs_real_t pfac, pfac1, pfac2, pfac3, unsvol;
3311  cs_real_t pif,pjf;
3312 
3313  for (g_id = 0; g_id < n_i_groups; g_id++) {
3314 # pragma omp parallel for private(face_id, ii, jj, difx, dify, difz, \
3315  djfx, djfy, djfz, pif, pjf, pfac, \
3316  pfac1, pfac2, pfac3)
3317  for (t_id = 0; t_id < n_i_threads; t_id++) {
3318  for (face_id = i_group_index[(t_id*n_i_groups + g_id)*2];
3319  face_id < i_group_index[(t_id*n_i_groups + g_id)*2 + 1];
3320  face_id++) {
3321 
3322  ii = i_face_cells[face_id][0];
3323  jj = i_face_cells[face_id][1];
3324 
3325  difx = i_face_cog[face_id][0] - cell_cen[ii][0];
3326  dify = i_face_cog[face_id][1] - cell_cen[ii][1];
3327  difz = i_face_cog[face_id][2] - cell_cen[ii][2];
3328  djfx = i_face_cog[face_id][0] - cell_cen[jj][0];
3329  djfy = i_face_cog[face_id][1] - cell_cen[jj][1];
3330  djfz = i_face_cog[face_id][2] - cell_cen[jj][2];
3331 
3332  pif = pvar[ii]
3333  + difx*grad[ii][0]+dify*grad[ii][1]+difz*grad[ii][2];
3334  pjf = pvar[jj]
3335  + djfx*grad[jj][0]+djfy*grad[jj][1]+djfz*grad[jj][2];
3336 
3337  pfac = pjf;
3338  if (i_massflux[face_id] > 0.) pfac = pif;
3339 
3340  pfac1 = pfac*i_face_normal[face_id][0];
3341  pfac2 = pfac*i_face_normal[face_id][1];
3342  pfac3 = pfac*i_face_normal[face_id][2];
3343 
3344  grdpa[ii][0] = grdpa[ii][0] + pfac1;
3345  grdpa[ii][1] = grdpa[ii][1] + pfac2;
3346  grdpa[ii][2] = grdpa[ii][2] + pfac3;
3347 
3348  grdpa[jj][0] = grdpa[jj][0] - pfac1;
3349  grdpa[jj][1] = grdpa[jj][1] - pfac2;
3350  grdpa[jj][2] = grdpa[jj][2] - pfac3;
3351 
3352  }
3353  }
3354  }
3355 
3356  for (g_id = 0; g_id < n_b_groups; g_id++) {
3357 # pragma omp parallel for private(face_id, ii, diipbx, diipby, diipbz, pfac) \
3358  if(m->n_b_faces > CS_THR_MIN)
3359  for (t_id = 0; t_id < n_b_threads; t_id++) {
3360  for (face_id = b_group_index[(t_id*n_b_groups + g_id)*2];
3361  face_id < b_group_index[(t_id*n_b_groups + g_id)*2 + 1];
3362  face_id++) {
3363 
3364  ii = b_face_cells[face_id];
3365 
3366  diipbx = diipb[face_id][0];
3367  diipby = diipb[face_id][1];
3368  diipbz = diipb[face_id][2];
3369  pfac = inc*coefap[face_id]
3370  + coefbp[face_id] * (pvar[ii] + diipbx*grad[ii][0]
3371  + diipby*grad[ii][1]
3372  + diipbz*grad[ii][2]);
3373  grdpa[ii][0] = grdpa[ii][0] + pfac*b_face_normal[face_id][0];
3374  grdpa[ii][1] = grdpa[ii][1] + pfac*b_face_normal[face_id][1];
3375  grdpa[ii][2] = grdpa[ii][2] + pfac*b_face_normal[face_id][2];
3376 
3377  }
3378  }
3379  }
3380 
3381 # pragma omp parallel for private(unsvol)
3382  for (cell_id = 0; cell_id < n_cells; cell_id++) {
3383 
3384  unsvol = 1./cell_vol[cell_id];
3385 
3386  grdpa[cell_id][0] = grdpa[cell_id][0]*unsvol;
3387  grdpa[cell_id][1] = grdpa[cell_id][1]*unsvol;
3388  grdpa[cell_id][2] = grdpa[cell_id][2]*unsvol;
3389 
3390  }
3391 
3392  /* Synchronization for parallelism or periodicity */
3393 
3394  if (halo != NULL) {
3395  cs_halo_sync_var_strided(halo, halo_type, (cs_real_t *)grdpa, 3);
3396  if (cs_glob_mesh->n_init_perio > 0)
3397  cs_halo_perio_sync_var_vect(halo, halo_type, (cs_real_t *)grdpa, 3);
3398 
3399  /* Gradient periodicity of rotation for Reynolds stress components */
3400  if (cs_glob_mesh->have_rotation_perio > 0 && f_id != -1)
3401  cs_gradient_perio_process_rij(&f_id, grdpa);
3402  }
3403 
3404 }
3405 
3406 /*----------------------------------------------------------------------------*/
3422 /*----------------------------------------------------------------------------*/
3423 
3424 inline static void
3426  const int inc,
3427  const cs_halo_type_t halo_type,
3428  cs_real_63_t *grad,
3429  cs_real_63_t *grdpa,
3430  cs_real_6_t *pvar,
3431  const cs_real_6_t *coefa,
3432  const cs_real_66_t *coefb,
3433  const cs_real_t *i_massflux)
3434 {
3435  const cs_mesh_t *m = cs_glob_mesh;
3436  const cs_halo_t *halo = m->halo;
3438 
3439  const cs_lnum_t n_cells = m->n_cells;
3440 
3441  const cs_lnum_2_t *restrict i_face_cells
3442  = (const cs_lnum_2_t *restrict)m->i_face_cells;
3443  const cs_lnum_t *restrict b_face_cells
3444  = (const cs_lnum_t *restrict)m->b_face_cells;
3445  const cs_real_t *restrict cell_vol = fvq->cell_vol;
3446  const cs_real_3_t *restrict cell_cen
3447  = (const cs_real_3_t *restrict)fvq->cell_cen;
3448  const cs_real_3_t *restrict i_face_normal
3449  = (const cs_real_3_t *restrict)fvq->i_face_normal;
3450  const cs_real_3_t *restrict b_face_normal
3451  = (const cs_real_3_t *restrict)fvq->b_face_normal;
3452  const cs_real_3_t *restrict i_face_cog
3453  = (const cs_real_3_t *restrict)fvq->i_face_cog;
3454  const cs_real_3_t *restrict diipb
3455  = (const cs_real_3_t *restrict)fvq->diipb;
3456 
3457  const int n_i_groups = m->i_face_numbering->n_groups;
3458  const int n_i_threads = m->i_face_numbering->n_threads;
3459  const int n_b_groups = m->b_face_numbering->n_groups;
3460  const int n_b_threads = m->b_face_numbering->n_threads;
3461  const cs_lnum_t *restrict i_group_index = m->i_face_numbering->group_index;
3462  const cs_lnum_t *restrict b_group_index = m->b_face_numbering->group_index;
3463 
3464  cs_lnum_t cell_id, face_id, ii, jj;
3465  int g_id, t_id;
3466  int isou, jsou;
3467  double vfac[3];
3468 
3469  cs_real_t pfac, unsvol;
3470  cs_real_t pif,pjf;
3471  cs_real_6_t difv, djfv, diipbv;
3472 
3473 
3474  for (g_id = 0; g_id < n_i_groups; g_id++) {
3475 # pragma omp parallel for private(face_id, ii, jj, isou, jsou, \
3476  difv, djfv, pif, pjf, pfac, vfac)
3477  for (t_id = 0; t_id < n_i_threads; t_id++) {
3478  for (face_id = i_group_index[(t_id*n_i_groups + g_id)*2];
3479  face_id < i_group_index[(t_id*n_i_groups + g_id)*2 + 1];
3480  face_id++) {
3481 
3482  ii = i_face_cells[face_id][0];
3483  jj = i_face_cells[face_id][1];
3484 
3485  for (jsou = 0; jsou < 3; jsou++) {
3486  difv[jsou] = i_face_cog[face_id][jsou] - cell_cen[ii][jsou];
3487  djfv[jsou] = i_face_cog[face_id][jsou] - cell_cen[jj][jsou];
3488  }
3489 
3490  /*-----------------
3491  X-Y-Z component, p=u, v, w */
3492 
3493  for (isou = 0; isou < 6; isou++) {
3494  pif = pvar[ii][isou];
3495  pjf = pvar[jj][isou];
3496  for (jsou = 0; jsou < 3; jsou++) {
3497  pif = pif + grad[ii][isou][jsou]*difv[jsou];
3498  pjf = pjf + grad[jj][isou][jsou]*djfv[jsou];
3499  }
3500 
3501  pfac = pjf;
3502  if (i_massflux[face_id] > 0.) pfac = pif;
3503 
3504  /* U gradient */
3505 
3506  for (jsou = 0; jsou < 3; jsou++) {
3507  vfac[jsou] = pfac*i_face_normal[face_id][jsou];
3508 
3509  grdpa[ii][isou][jsou] = grdpa[ii][isou][jsou] + vfac[jsou];
3510  grdpa[jj][isou][jsou] = grdpa[jj][isou][jsou] - vfac[jsou];
3511  }
3512  }
3513 
3514  }
3515  }
3516  }
3517 
3518  for (g_id = 0; g_id < n_b_groups; g_id++) {
3519 # pragma omp parallel for private(face_id, ii, isou, jsou, diipbv, pfac) \
3520  if(m->n_b_faces > CS_THR_MIN)
3521  for (t_id = 0; t_id < n_b_threads; t_id++) {
3522  for (face_id = b_group_index[(t_id*n_b_groups + g_id)*2];
3523  face_id < b_group_index[(t_id*n_b_groups + g_id)*2 + 1];
3524  face_id++) {
3525 
3526  ii = b_face_cells[face_id];
3527 
3528  for (jsou = 0; jsou < 3; jsou++)
3529  diipbv[jsou] = diipb[face_id][jsou];
3530 
3531  /*-----------------
3532  X-Y-Z components, p=u, v, w */
3533  for (isou = 0; isou <6; isou++) {
3534  pfac = inc*coefa[face_id][isou];
3535  /*coefu is a matrix */
3536  for (jsou = 0; jsou < 6; jsou++)
3537  pfac += coefb[face_id][jsou][isou]*( pvar[ii][jsou]
3538  + grad[ii][jsou][0]*diipbv[0]
3539  + grad[ii][jsou][1]*diipbv[1]
3540  + grad[ii][jsou][2]*diipbv[2]);
3541 
3542  for (jsou = 0; jsou < 3; jsou++)
3543  grdpa[ii][isou][jsou] += pfac*b_face_normal[face_id][jsou];
3544  }
3545 
3546  }
3547  }
3548  }
3549 
3550 # pragma omp parallel for private(isou, jsou, unsvol)
3551  for (cell_id = 0; cell_id < n_cells; cell_id++) {
3552  unsvol = 1./cell_vol[cell_id];
3553  for (isou = 0; isou < 6; isou++) {
3554  for (jsou = 0; jsou < 3; jsou++)
3555  grdpa[cell_id][isou][jsou] = grdpa[cell_id][isou][jsou]*unsvol;
3556  }
3557  }
3558 
3559  /* Handle parallelism and periodicity */
3560 
3561  if (halo != NULL) {
3562  cs_halo_sync_var_strided(halo, halo_type, (cs_real_t *)grdpa, 18);
3563  if (m->n_init_perio > 0)
3564  cs_halo_perio_sync_var_sym_tens(halo, halo_type, (cs_real_t *)grdpa);
3565  }
3566 
3567 }
3568 
3569 /*============================================================================
3570  * Public function prototypes for Fortran API
3571  *============================================================================*/
3572 
3573 /*----------------------------------------------------------------------------
3574  * Wrapper to cs_convection_diffusion_scalar
3575  *----------------------------------------------------------------------------*/
3576 
3577 void CS_PROCF (bilsc2, BILSC2)
3578 (
3579  const cs_int_t *const idtvar,
3580  const cs_int_t *const f_id,
3581  const cs_var_cal_opt_t *const var_cal_opt,
3582  const cs_int_t *const icvflb,
3583  const cs_int_t *const inc,
3584  const cs_int_t *const iccocg,
3585  const cs_int_t *const ifaccp,
3586  cs_real_t pvar[],
3587  const cs_real_t pvara[],
3588  const cs_int_t bc_type[],
3589  const cs_int_t icvfli[],
3590  const cs_real_t coefap[],
3591  const cs_real_t coefbp[],
3592  const cs_real_t cofafp[],
3593  const cs_real_t cofbfp[],
3594  const cs_real_t i_massflux[],
3595  const cs_real_t b_massflux[],
3596  const cs_real_t i_visc[],
3597  const cs_real_t b_visc[],
3598  cs_real_t rhs[]
3599 );
3600 
3601 /*----------------------------------------------------------------------------
3602  * Wrapper to cs_convection_diffusion_vector
3603  *----------------------------------------------------------------------------*/
3604 
3605 void CS_PROCF (bilsc4, BILSC4)
3606 (
3607  const cs_int_t *const idtvar,
3608  const cs_int_t *const f_id,
3609  const cs_var_cal_opt_t *const var_cal_opt,
3610  const cs_int_t *const icvflb,
3611  const cs_int_t *const inc,
3612  const cs_int_t *const ifaccp,
3613  const cs_int_t *const ivisep,
3614  cs_real_3_t pvar[],
3615  const cs_real_3_t pvara[],
3616  const cs_int_t bc_type[],
3617  const cs_int_t icvfli[],
3618  const cs_real_3_t coefav[],
3619  const cs_real_33_t coefbv[],
3620  const cs_real_3_t cofafv[],
3621  const cs_real_33_t cofbfv[],
3622  const cs_real_t i_massflux[],
3623  const cs_real_t b_massflux[],
3624  const cs_real_t i_visc[],
3625  const cs_real_t b_visc[],
3626  const cs_real_t secvif[],
3627  cs_real_3_t rhs[]
3628 );
3629 /*----------------------------------------------------------------------------
3630  * Wrapper to cs_convection_diffusion_tensor
3631  *----------------------------------------------------------------------------*/
3632 
3633 void CS_PROCF (bilsc6, BILSC6)
3634 (
3635  const cs_int_t *const idtvar,
3636  const cs_int_t *const f_id,
3637  const cs_var_cal_opt_t *const var_cal_opt,
3638  const cs_int_t *const icvflb,
3639  const cs_int_t *const inc,
3640  const cs_int_t *const ifaccp,
3641  cs_real_6_t pvar[],
3642  const cs_real_6_t pvara[],
3643  const cs_int_t bc_type[],
3644  const cs_int_t icvfli[],
3645  const cs_real_6_t coefa[],
3646  const cs_real_66_t coefb[],
3647  const cs_real_6_t cofaf[],
3648  const cs_real_66_t cofbf[],
3649  const cs_real_t i_massflux[],
3650  const cs_real_t b_massflux[],
3651  const cs_real_t i_visc[],
3652  const cs_real_t b_visc[],
3653  cs_real_6_t rhs[]
3654 );
3655 /*----------------------------------------------------------------------------
3656  * Wrapper to cs_convection_diffusion_thermal
3657  *----------------------------------------------------------------------------*/
3658 
3659 void CS_PROCF (bilsct, BILSCT)
3660 (
3661  const cs_int_t *const idtvar,
3662  const cs_int_t *const f_id,
3663  const cs_var_cal_opt_t *const var_cal_opt,
3664  const cs_int_t *const inc,
3665  const cs_int_t *const iccocg,
3666  const cs_int_t *const ifaccp,
3667  cs_real_t pvar[],
3668  const cs_real_t pvara[],
3669  const cs_int_t bc_type[],
3670  const cs_real_t coefap[],
3671  const cs_real_t coefbp[],
3672  const cs_real_t cofafp[],
3673  const cs_real_t cofbfp[],
3674  const cs_real_t i_massflux[],
3675  const cs_real_t b_massflux[],
3676  const cs_real_t i_visc[],
3677  const cs_real_t b_visc[],
3678  const cs_real_t xcpp[],
3679  cs_real_t rhs[]
3680 );
3681 
3682 /*----------------------------------------------------------------------------
3683  * Wrapper to cs_anisotropic_diffusion_scalar
3684  *----------------------------------------------------------------------------*/
3685 
3686 void CS_PROCF (diften, DIFTEN)
3687 (
3688  const cs_int_t *const idtvar,
3689  const cs_int_t *const f_id,
3690  const cs_var_cal_opt_t *const var_cal_opt,
3691  const cs_int_t *const inc,
3692  const cs_int_t *const iccocg,
3693  cs_real_t pvar[],
3694  const cs_real_t pvara[],
3695  const cs_real_t coefap[],
3696  const cs_real_t coefbp[],
3697  const cs_real_t cofafp[],
3698  const cs_real_t cofbfp[],
3699  const cs_real_t i_visc[],
3700  const cs_real_t b_visc[],
3701  cs_real_6_t viscel[],
3702  const cs_real_2_t weighf[],
3703  const cs_real_t weighb[],
3704  cs_real_t rhs[]
3705 );
3706 
3707 /*----------------------------------------------------------------------------
3708  * Wrapper to cs_anisotropic_diffusion_vector
3709  *----------------------------------------------------------------------------*/
3710 
3711 void CS_PROCF (diftnv, DIFTNV)
3712 (
3713  const cs_int_t *const idtvar,
3714  const cs_int_t *const f_id,
3715  const cs_var_cal_opt_t *const var_cal_opt,
3716  const cs_int_t *const inc,
3717  const cs_int_t *const ifaccp,
3718  const cs_int_t *const ivisep,
3719  cs_real_3_t pvar[],
3720  const cs_real_3_t pvara[],
3721  const cs_int_t bc_type[],
3722  const cs_real_3_t coefav[],
3723  const cs_real_33_t coefbv[],
3724  const cs_real_3_t cofafv[],
3725  const cs_real_33_t cofbfv[],
3726  const cs_real_33_t i_visc[],
3727  const cs_real_t b_visc[],
3728  const cs_real_t secvif[],
3729  cs_real_3_t rhs[]
3730 );
3731 /*----------------------------------------------------------------------------
3732  * Wrapper to cs_anisotropic_diffusion_tensor
3733  *----------------------------------------------------------------------------*/
3734 
3735 void CS_PROCF (diftnts, DIFTNTS)
3736 (
3737  const cs_int_t *const idtvar,
3738  const cs_int_t *const f_id,
3739  const cs_var_cal_opt_t *const var_cal_opt,
3740  const cs_int_t *const inc,
3741  cs_real_6_t pvar[],
3742  const cs_real_6_t pvara[],
3743  const cs_int_t bc_type[],
3744  const cs_real_6_t coefa[],
3745  const cs_real_66_t coefb[],
3746  const cs_real_6_t cofaf[],
3747  const cs_real_66_t cofbf[],
3748  const cs_real_t i_visc[],
3749  const cs_real_t b_visc[],
3750  cs_real_6_t viscel[],
3751  const cs_real_2_t weighf[],
3752  const cs_real_t weighb[],
3753  cs_real_6_t rhs[]
3754 );
3755 /*----------------------------------------------------------------------------
3756  * Wrapper to cs_face_diffusion_potential
3757  *----------------------------------------------------------------------------*/
3758 
3759 void CS_PROCF (itrmas, ITRMAS)
3760 (
3761  const cs_int_t *const f_id,
3762  const cs_int_t *const init,
3763  const cs_int_t *const inc,
3764  const cs_int_t *const imrgra,
3765  const cs_int_t *const iccocg,
3766  const cs_int_t *const nswrgp,
3767  const cs_int_t *const imligp,
3768  const cs_int_t *const iphydp,
3769  const cs_int_t *const iwarnp,
3770  const cs_real_t *const epsrgp,
3771  const cs_real_t *const climgp,
3772  const cs_real_t *const extrap,
3773  cs_real_3_t frcxt[],
3774  cs_real_t pvar[],
3775  const cs_real_t coefap[],
3776  const cs_real_t coefbp[],
3777  const cs_real_t cofafp[],
3778  const cs_real_t cofbfp[],
3779  const cs_real_t i_visc[],
3780  const cs_real_t b_visc[],
3781  const cs_real_t viselx[],
3782  const cs_real_t visely[],
3783  const cs_real_t viselz[],
3784  cs_real_t i_massflux[],
3785  cs_real_t b_massflux[]
3786 );
3787 
3788 /*----------------------------------------------------------------------------
3789  * Wrapper to cs_face_anisotropic_diffusion_potential
3790  *----------------------------------------------------------------------------*/
3791 
3792 void CS_PROCF (itrmav, ITRMAV)
3793 (
3794  const cs_int_t *const f_id,
3795  const cs_int_t *const init,
3796  const cs_int_t *const inc,
3797  const cs_int_t *const imrgra,
3798  const cs_int_t *const iccocg,
3799  const cs_int_t *const nswrgp,
3800  const cs_int_t *const imligp,
3801  const cs_int_t *const ircflp,
3802  const cs_int_t *const iphydp,
3803  const cs_int_t *const iwarnp,
3804  const cs_real_t *const epsrgp,
3805  const cs_real_t *const climgp,
3806  const cs_real_t *const extrap,
3807  cs_real_3_t frcxt[],
3808  cs_real_t pvar[],
3809  const cs_real_t coefap[],
3810  const cs_real_t coefbp[],
3811  const cs_real_t cofafp[],
3812  const cs_real_t cofbfp[],
3813  const cs_real_t i_visc[],
3814  const cs_real_t b_visc[],
3815  cs_real_6_t viscel[],
3816  const cs_real_2_t weighf[],
3817  const cs_real_t weighb[],
3818  cs_real_t i_massflux[],
3819  cs_real_t b_massflux[]
3820 );
3821 
3822 /*----------------------------------------------------------------------------
3823  * Wrapper to cs_diffusion_potential
3824  *----------------------------------------------------------------------------*/
3825 
3826 void CS_PROCF (itrgrp, ITRGRP)
3827 (
3828  const cs_int_t *const f_id,
3829  const cs_int_t *const init,
3830  const cs_int_t *const inc,
3831  const cs_int_t *const imrgra,
3832  const cs_int_t *const iccocg,
3833  const cs_int_t *const nswrgp,
3834  const cs_int_t *const imligp,
3835  const cs_int_t *const iphydp,
3836  const cs_int_t *const iwarnp,
3837  const cs_real_t *const epsrgp,
3838  const cs_real_t *const climgp,
3839  const cs_real_t *const extrap,
3840  cs_real_3_t frcxt[],
3841  cs_real_t pvar[],
3842  const cs_real_t coefap[],
3843  const cs_real_t coefbp[],
3844  const cs_real_t cofafp[],
3845  const cs_real_t cofbfp[],
3846  const cs_real_t i_visc[],
3847  const cs_real_t b_visc[],
3848  const cs_real_t viselx[],
3849  const cs_real_t visely[],
3850  const cs_real_t viselz[],
3851  cs_real_t diverg[]
3852 );
3853 
3854 /*----------------------------------------------------------------------------
3855  * Wrapper to cs_anisotropic_diffusion_potential
3856  *----------------------------------------------------------------------------*/
3857 
3858 void CS_PROCF (itrgrv, ITRGRV)
3859 (
3860  const cs_int_t *const f_id,
3861  const cs_int_t *const init,
3862  const cs_int_t *const inc,
3863  const cs_int_t *const imrgra,
3864  const cs_int_t *const iccocg,
3865  const cs_int_t *const nswrgp,
3866  const cs_int_t *const imligp,
3867  const cs_int_t *const ircflp,
3868  const cs_int_t *const iphydp,
3869  const cs_int_t *const iwarnp,
3870  const cs_real_t *const epsrgp,
3871  const cs_real_t *const climgp,
3872  const cs_real_t *const extrap,
3873  cs_real_3_t frcxt[],
3874  cs_real_t pvar[],
3875  const cs_real_t coefap[],
3876  const cs_real_t coefbp[],
3877  const cs_real_t cofafp[],
3878  const cs_real_t cofbfp[],
3879  const cs_real_t i_visc[],
3880  const cs_real_t b_visc[],
3881  cs_real_6_t viscel[],
3882  const cs_real_2_t weighf[],
3883  const cs_real_t weighb[],
3884  cs_real_t diverg[]
3885 );
3886 
3887 /*=============================================================================
3888  * Public function prototypes
3889  *============================================================================*/
3890 
3891 /*----------------------------------------------------------------------------*/
3946 /*----------------------------------------------------------------------------*/
3947 
3948 void
3950  int f_id,
3951  const cs_var_cal_opt_t var_cal_opt,
3952  int icvflb,
3953  int inc,
3954  int iccocg,
3955  int ifaccp,
3956  cs_real_t *restrict pvar,
3957  const cs_real_t *restrict pvara,
3958  const cs_int_t bc_type[],
3959  const cs_int_t icvfli[],
3960  const cs_real_t coefap[],
3961  const cs_real_t coefbp[],
3962  const cs_real_t cofafp[],
3963  const cs_real_t cofbfp[],
3964  const cs_real_t i_massflux[],
3965  const cs_real_t b_massflux[],
3966  const cs_real_t i_visc[],
3967  const cs_real_t b_visc[],
3968  cs_real_t *restrict rhs);
3969 
3970 /*----------------------------------------------------------------------------*/
4032 /*----------------------------------------------------------------------------*/
4033 
4034 void
4036  int f_id,
4037  const cs_var_cal_opt_t var_cal_opt,
4038  int icvflb,
4039  int inc,
4040  int ifaccp,
4041  int ivisep,
4042  cs_real_3_t *restrict pvar,
4043  const cs_real_3_t *restrict pvara,
4044  const cs_int_t bc_type[],
4045  const cs_int_t icvfli[],
4046  const cs_real_3_t coefav[],
4047  const cs_real_33_t coefbv[],
4048  const cs_real_3_t cofafv[],
4049  const cs_real_33_t cofbfv[],
4050  const cs_real_t i_massflux[],
4051  const cs_real_t b_massflux[],
4052  const cs_real_t i_visc[],
4053  const cs_real_t b_visc[],
4054  const cs_real_t secvif[],
4055  cs_real_3_t *restrict rhs);
4056 /*----------------------------------------------------------------------------*/
4107 /*----------------------------------------------------------------------------*/
4108 
4109 void
4111  int f_id,
4112  const cs_var_cal_opt_t var_cal_opt,
4113  int icvflb,
4114  int inc,
4115  int ifaccp,
4116  cs_real_6_t *restrict pvar,
4117  const cs_real_6_t *restrict pvara,
4118  const cs_int_t bc_type[],
4119  const cs_int_t icvfli[],
4120  const cs_real_6_t coefa[],
4121  const cs_real_66_t coefb[],
4122  const cs_real_6_t cofaf[],
4123  const cs_real_66_t cofbf[],
4124  const cs_real_t i_massflux[],
4125  const cs_real_t b_massflux[],
4126  const cs_real_t i_visc[],
4127  const cs_real_t b_visc[],
4128  cs_real_6_t *restrict rhs);
4129 
4130 /*----------------------------------------------------------------------------*/
4179 /*----------------------------------------------------------------------------*/
4180 
4181 void
4183  int f_id,
4184  const cs_var_cal_opt_t var_cal_opt,
4185  int inc,
4186  int iccocg,
4187  int ifaccp,
4188  cs_real_t *restrict pvar,
4189  const cs_real_t *restrict pvara,
4190  const cs_int_t bc_type[],
4191  const cs_real_t coefap[],
4192  const cs_real_t coefbp[],
4193  const cs_real_t cofafp[],
4194  const cs_real_t cofbfp[],
4195  const cs_real_t i_massflux[],
4196  const cs_real_t b_massflux[],
4197  const cs_real_t i_visc[],
4198  const cs_real_t b_visc[],
4199  const cs_real_t xcpp[],
4200  cs_real_t *restrict rhs);
4201 
4202 /*----------------------------------------------------------------------------*/
4250 /*----------------------------------------------------------------------------*/
4251 
4252 void
4254  int f_id,
4255  const cs_var_cal_opt_t var_cal_opt,
4256  int inc,
4257  int iccocg,
4258  cs_real_t *restrict pvar,
4259  const cs_real_t *restrict pvara,
4260  const cs_real_t coefap[],
4261  const cs_real_t coefbp[],
4262  const cs_real_t cofafp[],
4263  const cs_real_t cofbfp[],
4264  const cs_real_t i_visc[],
4265  const cs_real_t b_visc[],
4266  cs_real_6_t *restrict viscel,
4267  const cs_real_2_t weighf[],
4268  const cs_real_t weighb[],
4269  cs_real_t *restrict rhs);
4270 
4271 /*-----------------------------------------------------------------------------*/
4318 /*----------------------------------------------------------------------------*/
4319 
4320 void
4322  int f_id,
4323  const cs_var_cal_opt_t var_cal_opt,
4324  int inc,
4325  int ifaccp,
4326  int ivisep,
4327  cs_real_3_t *restrict pvar,
4328  const cs_real_3_t *restrict pvara,
4329  const cs_int_t bc_type[],
4330  const cs_real_3_t coefav[],
4331  const cs_real_33_t coefbv[],
4332  const cs_real_3_t cofafv[],
4333  const cs_real_33_t cofbfv[],
4334  const cs_real_33_t i_visc[],
4335  const cs_real_t b_visc[],
4336  const cs_real_t secvif[],
4337  cs_real_3_t *restrict rhs);
4338 /*----------------------------------------------------------------------------*/
4383 /*----------------------------------------------------------------------------*/
4384 
4385 void
4387  int f_id,
4388  const cs_var_cal_opt_t var_cal_opt,
4389  int inc,
4390  cs_real_6_t *restrict pvar,
4391  const cs_real_6_t *restrict pvara,
4392  const cs_int_t bc_type[],
4393  const cs_real_6_t coefa[],
4394  const cs_real_66_t coefb[],
4395  const cs_real_6_t cofaf[],
4396  const cs_real_66_t cofbf[],
4397  const cs_real_t i_visc[],
4398  const cs_real_t b_visc[],
4399  cs_real_6_t *restrict viscel,
4400  const cs_real_2_t weighf[],
4401  const cs_real_t weighb[],
4402  cs_real_6_t *restrict rhs);
4403 
4404 /*----------------------------------------------------------------------------*/
4463 /*----------------------------------------------------------------------------*/
4464 
4465 void
4466 cs_face_diffusion_potential(const int f_id,
4467  const cs_mesh_t *m,
4468  cs_mesh_quantities_t *fvq,
4469  int init,
4470  int inc,
4471  int imrgra,
4472  int iccocg,
4473  int nswrgp,
4474  int imligp,
4475  int iphydp,
4476  int iwarnp,
4477  double epsrgp,
4478  double climgp,
4479  double extrap,
4480  cs_real_3_t *restrict frcxt,
4481  cs_real_t *restrict pvar,
4482  const cs_real_t coefap[],
4483  const cs_real_t coefbp[],
4484  const cs_real_t cofafp[],
4485  const cs_real_t cofbfp[],
4486  const cs_real_t i_visc[],
4487  const cs_real_t b_visc[],
4488  const cs_real_t viselx[],
4489  const cs_real_t visely[],
4490  const cs_real_t viselz[],
4491  cs_real_t *restrict i_massflux,
4492  cs_real_t *restrict b_massflux);
4493 
4494 /*----------------------------------------------------------------------------*/
4562 /*----------------------------------------------------------------------------*/
4563 
4564 void
4566  const cs_mesh_t *m,
4567  cs_mesh_quantities_t *fvq,
4568  int init,
4569  int inc,
4570  int imrgra,
4571  int iccocg,
4572  int nswrgp,
4573  int imligp,
4574  int ircflp,
4575  int iphydp,
4576  int iwarnp,
4577  double epsrgp,
4578  double climgp,
4579  double extrap,
4580  cs_real_3_t *restrict frcxt,
4581  cs_real_t *restrict pvar,
4582  const cs_real_t coefap[],
4583  const cs_real_t coefbp[],
4584  const cs_real_t cofafp[],
4585  const cs_real_t cofbfp[],
4586  const cs_real_t i_visc[],
4587  const cs_real_t b_visc[],
4588  cs_real_6_t *restrict viscel,
4589  const cs_real_2_t weighf[],
4590  const cs_real_t weighb[],
4591  cs_real_t *restrict i_massflux,
4592  cs_real_t *restrict b_massflux);
4593 
4594 /*----------------------------------------------------------------------------*/
4652 /*----------------------------------------------------------------------------*/
4653 
4654 void
4655 cs_diffusion_potential(const int f_id,
4656  const cs_mesh_t *m,
4657  cs_mesh_quantities_t *fvq,
4658  int init,
4659  int inc,
4660  int imrgra,
4661  int iccocg,
4662  int nswrgp,
4663  int imligp,
4664  int iphydp,
4665  int iwarnp,
4666  double epsrgp,
4667  double climgp,
4668  double extrap,
4669  cs_real_3_t *restrict frcxt,
4670  cs_real_t *restrict pvar,
4671  const cs_real_t coefap[],
4672  const cs_real_t coefbp[],
4673  const cs_real_t cofafp[],
4674  const cs_real_t cofbfp[],
4675  const cs_real_t i_visc[],
4676  const cs_real_t b_visc[],
4677  const cs_real_t viselx[],
4678  const cs_real_t visely[],
4679  const cs_real_t viselz[],
4680  cs_real_t *restrict diverg);
4681 
4682 /*----------------------------------------------------------------------------*/
4751 /*----------------------------------------------------------------------------*/
4752 
4753 void
4754 cs_anisotropic_diffusion_potential(const int f_id,
4755  const cs_mesh_t *m,
4756  cs_mesh_quantities_t *fvq,
4757  int init,
4758  int inc,
4759  int imrgra,
4760  int iccocg,
4761  int nswrgp,
4762  int imligp,
4763  int ircflp,
4764  int iphydp,
4765  int iwarnp,
4766  double epsrgp,
4767  double climgp,
4768  double extrap,
4769  cs_real_3_t *restrict frcxt,
4770  cs_real_t *restrict pvar,
4771  const cs_real_t coefap[],
4772  const cs_real_t coefbp[],
4773  const cs_real_t cofafp[],
4774  const cs_real_t cofbfp[],
4775  const cs_real_t i_visc[],
4776  const cs_real_t b_visc[],
4777  cs_real_6_t *restrict viscel,
4778  const cs_real_2_t weighf[],
4779  const cs_real_t weighb[],
4780  cs_real_t *restrict diverg);
4781 
4782 /*----------------------------------------------------------------------------*/
4783 
4785 
4786 #endif /* __CS_CONVECTION_DIFFUSION_H__ */
static void cs_b_cd_steady(const int ircflp, const double relaxp, const cs_real_3_t diipb, const cs_real_3_t gradi, const cs_real_t pi, const cs_real_t pia, cs_real_t *pir, cs_real_t *pipr)
Handle preparation of boundary face values for the flux computation in case of a steady algorithm...
Definition: cs_convection_diffusion.h:3115
void cs_gradient_perio_process_rij(const cs_int_t *f_id, cs_real_3_t grad[])
Process grad buffers in case of rotation on Reynolds stress tensor.
Definition: cs_gradient_perio.c:688
static void cs_b_no_relax_c_val_tensor(const cs_real_6_t pi, const cs_real_6_t recoi, cs_real_t pir[6], cs_real_t pipr[6])
Copy values at bounadry cell i for consistency with relaxation case.
Definition: cs_convection_diffusion.h:2735
double precision, dimension(:,:), pointer diipb
Definition: mesh.f90:216
#define restrict
Definition: cs_defs.h:122
static void cs_i_diff_flux(const int idiffp, 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_2_t fluxij)
Add diffusive fluxes to fluxes at face ij.
Definition: cs_convection_diffusion.h:935
static void cs_slope_test_gradient(const int f_id, const int inc, const cs_halo_type_t halo_type, cs_real_3_t *grad, cs_real_3_t *grdpa, cs_real_t *pvar, const cs_real_t *coefap, const cs_real_t *coefbp, const cs_real_t *i_massflux)
Compute the upwind gradient used in the slope tests.
Definition: cs_convection_diffusion.h:3266
void cs_halo_sync_var_strided(const cs_halo_t *halo, cs_halo_type_t sync_mode, cs_real_t var[], int stride)
Definition: cs_halo.c:1288
cs_real_t cs_real_2_t[2]
vector of 2 floating-point values
Definition: cs_defs.h:306
static void cs_centered_f_val_tensor(const double pnd, const cs_real_6_t pip, const cs_real_6_t pjp, const cs_real_6_t pipr, const cs_real_6_t pjpr, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6])
Prepare value at face ij by using a centered scheme.
Definition: cs_convection_diffusion.h:596
cs_real_t cs_real_6_t[6]
vector of 6 floating-point values
Definition: cs_defs.h:309
static void cs_blend_f_val_tensor(const double blencp, const cs_real_6_t pi, const cs_real_6_t pj, const cs_real_6_t pir, const cs_real_6_t pjr, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6])
Blend face values for a centered or SOLU scheme with face values for an upwind scheme.
Definition: cs_convection_diffusion.h:784
void cs_halo_perio_sync_var_sym_tens(const cs_halo_t *halo, cs_halo_type_t sync_mode, cs_real_t var[])
Definition: cs_halo_perio.c:1002
void cs_anisotropic_diffusion_vector(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int inc, int ifaccp, int ivisep, cs_real_3_t *restrict pvar, const cs_real_3_t *restrict pvara, const cs_int_t bc_type[], const cs_real_3_t coefav[], const cs_real_33_t coefbv[], const cs_real_3_t cofafv[], const cs_real_33_t cofbfv[], const cs_real_33_t i_visc[], const cs_real_t b_visc[], const cs_real_t secvif[], cs_real_3_t *restrict rhs)
Add the explicit part of the diffusion terms with a symmetric tensorial diffusivity for a transport e...
Definition: cs_convection_diffusion.c:6088
static void cs_i_conv_flux(const int iconvp, 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, cs_real_2_t fluxij)
Add convective fluxes (substracting the mass accumulation from them) to fluxes at face ij...
Definition: cs_convection_diffusion.h:821
void cs_anisotropic_diffusion_scalar(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int inc, int iccocg, cs_real_t *restrict pvar, const cs_real_t *restrict pvara, const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t *restrict viscel, const cs_real_2_t weighf[], const cs_real_t weighb[], cs_real_t *restrict rhs)
Add the explicit part of the diffusion terms with a symmetric tensor diffusivity for a transport equa...
Definition: cs_convection_diffusion.c:5513
integer, dimension(:), allocatable icvfli
Definition: cfpoin.f90:48
static void cs_i_cd_unsteady_slope_test_tensor(bool upwind_switch[6], const int ircflp, const int ischcp, const double blencp, const cs_real_t weight, const cs_real_t i_dist, const cs_real_t i_face_surf, const cs_real_3_t cell_ceni, const cs_real_3_t cell_cenj, const cs_real_3_t i_face_normal, const cs_real_3_t i_face_cog, const cs_real_3_t dijpf, const cs_real_t i_massflux, const cs_real_63_t gradi, const cs_real_63_t gradj, const cs_real_63_t grdpai, const cs_real_63_t grdpaj, const cs_real_6_t pi, const cs_real_6_t pj, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6], cs_real_t pip[6], cs_real_t pjp[6], cs_real_t pipr[6], cs_real_t pjpr[6])
Handle preparation of internal face values for the fluxes computation in case of a unsteady algorithm...
Definition: cs_convection_diffusion.h:2459
void cs_anisotropic_diffusion_tensor(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int inc, cs_real_6_t *restrict pvar, const cs_real_6_t *restrict pvara, const cs_int_t bc_type[], const cs_real_6_t coefa[], const cs_real_66_t coefb[], const cs_real_6_t cofaf[], const cs_real_66_t cofbf[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t *restrict viscel, const cs_real_2_t weighf[], const cs_real_t weighb[], cs_real_6_t *restrict rhs)
Add the explicit part of the diffusion terms with a symmetric tensor diffusivity for a transport equa...
Definition: cs_convection_diffusion.c:6619
static void cs_b_upwind_flux(const int iconvp, const int inc, const int ifaccp, 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, cs_real_t *flux)
Add convective flux (substracting the mass accumulation from it) to flux at boundary face...
Definition: cs_convection_diffusion.h:2907
cs_real_t cs_real_66_t[6][6]
6x6 matrix of floating-point values
Definition: cs_defs.h:312
static void cs_b_cd_unsteady(const int ircflp, const cs_real_3_t diipb, const cs_real_3_t gradi, const cs_real_t pi, cs_real_t *pir, cs_real_t *pipr)
Handle preparation of boundary face values for the flux computation in case of an unsteady algorithm...
Definition: cs_convection_diffusion.h:3194
double precision pi
value with 16 digits
Definition: cstnum.f90:48
#define BEGIN_C_DECLS
Definition: cs_defs.h:429
int cs_int_t
Fortran-compatible integer.
Definition: cs_defs.h:295
static void cs_b_no_relax_c_val(const cs_real_t pi, const cs_real_t recoi, cs_real_t *pir, cs_real_t *pipr)
Copy values at bounadry cell i for consistency with relaxation case.
Definition: cs_convection_diffusion.h:2715
cs_mesh_t * cs_glob_mesh
static void cs_i_cd_steady_tensor(const int ircflp, const int ischcp, const double relaxp, 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_3_t dijpf, const cs_real_63_t gradi, const cs_real_63_t gradj, const cs_real_6_t pi, const cs_real_6_t pj, const cs_real_6_t pia, const cs_real_6_t pja, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6], cs_real_t pip[6], cs_real_t pjp[6], cs_real_t pipr[6], cs_real_t pjpr[6])
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion.h:1497
static void cs_i_relax_c_val_tensor(const double relaxp, const cs_real_6_t pia, const cs_real_6_t pja, const cs_real_6_t recoi, const cs_real_6_t recoj, const cs_real_6_t pi, const cs_real_6_t pj, cs_real_t pir[6], cs_real_t pjr[6], cs_real_t pipr[6], cs_real_t pjpr[6])
Compute relaxed values at cell i and j.
Definition: cs_convection_diffusion.h:389
static 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.h:351
static void cs_b_upwind_flux_cons(const int iconvp, const int inc, const int ifaccp, const cs_int_t bc_type, 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 (conservative formulation) to flux at boundary face. The convective flux is a pur...
Definition: cs_convection_diffusion.h:3013
void diftnv(const cs_int_t *const idtvar, const cs_int_t *const f_id, const cs_var_cal_opt_t *const var_cal_opt, const cs_int_t *const inc, const cs_int_t *const ifaccp, const cs_int_t *const ivisep, cs_real_3_t pvar[], const cs_real_3_t pvara[], const cs_int_t bc_type[], const cs_real_3_t coefav[], const cs_real_33_t coefbv[], const cs_real_3_t cofafv[], const cs_real_33_t cofbfv[], const cs_real_33_t i_visc[], const cs_real_t b_visc[], const cs_real_t secvif[], cs_real_3_t rhs[])
Definition: cs_convection_diffusion.c:405
Definition: cs_halo.h:70
integer(c_int), pointer, save idtvar
option for a variable time step
Definition: optcal.f90:427
static void cs_b_upwind_flux_tensor(const int iconvp, const int inc, const int ifaccp, const int bc_type, const cs_real_6_t pi, const cs_real_6_t pir, const cs_real_6_t pipr, const cs_real_6_t coefa, const cs_real_66_t coefb, const cs_real_t b_massflux, cs_real_t flux[6])
Add convective flux (substracting the mass accumulation from it) to flux at boundary face...
Definition: cs_convection_diffusion.h:2956
static void cs_i_cd_unsteady_upwind_tensor(const int ircflp, 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_3_t dijpf, const cs_real_63_t gradi, const cs_real_63_t gradj, const cs_real_6_t pi, const cs_real_6_t pj, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6], cs_real_t pip[6], cs_real_t pjp[6], cs_real_t pipr[6], cs_real_t pjpr[6])
Handle preparation of internal face values for the fluxes computation in case of an unsteady algorith...
Definition: cs_convection_diffusion.h:1271
static void cs_i_cd_steady_slope_test_tensor(bool upwind_switch[6], const int ircflp, const int ischcp, const double relaxp, const double blencp, const cs_real_t weight, const cs_real_t i_dist, const cs_real_t i_face_surf, const cs_real_3_t cell_ceni, const cs_real_3_t cell_cenj, const cs_real_3_t i_face_normal, const cs_real_3_t i_face_cog, const cs_real_3_t dijpf, const cs_real_t i_massflux, const cs_real_63_t gradi, const cs_real_63_t gradj, const cs_real_63_t grdpai, const cs_real_63_t grdpaj, const cs_real_6_t pi, const cs_real_6_t pj, const cs_real_6_t pia, const cs_real_6_t pja, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6], cs_real_t pip[6], cs_real_t pjp[6], cs_real_t pipr[6], cs_real_t pjpr[6])
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion.h:2087
void cs_convection_diffusion_vector(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int icvflb, int inc, int ifaccp, int ivisep, cs_real_3_t *restrict pvar, const cs_real_3_t *restrict pvara, const cs_int_t bc_type[], const cs_int_t icvfli[], const cs_real_3_t coefav[], const cs_real_33_t coefbv[], const cs_real_3_t cofafv[], const cs_real_33_t cofbfv[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t secvif[], cs_real_3_t *restrict rhs)
Add the explicit part of the convection/diffusion terms of a transport equation of a vector field ...
Definition: cs_convection_diffusion.c:1808
void bilsc2(const cs_int_t *const idtvar, const cs_int_t *const f_id, const cs_var_cal_opt_t *const var_cal_opt, const cs_int_t *const icvflb, const cs_int_t *const inc, const cs_int_t *const iccocg, const cs_int_t *const ifaccp, cs_real_t pvar[], const cs_real_t pvara[], const cs_int_t bc_type[], const cs_int_t icvfli[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_t rhs[])
Definition: cs_convection_diffusion.c:161
static void cs_solu_f_val_tensor(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_63_t gradi, const cs_real_63_t gradj, const cs_real_6_t pi, const cs_real_6_t pj, const cs_real_6_t pir, const cs_real_6_t pjr, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6])
Prepare value at face ij by using a Second Order Linear Upwind scheme.
Definition: cs_convection_diffusion.h:687
void cs_convection_diffusion_tensor(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int icvflb, int inc, int ifaccp, cs_real_6_t *restrict pvar, const cs_real_6_t *restrict pvara, const cs_int_t bc_type[], const cs_int_t icvfli[], const cs_real_6_t coefa[], const cs_real_66_t coefb[], const cs_real_6_t cofaf[], const cs_real_66_t cofbf[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t *restrict rhs)
Add the explicit part of the convection/diffusion terms of a transport equation of a vector field ...
Definition: cs_convection_diffusion.c:3388
Definition: cs_mesh.h:62
static void cs_blend_f_val(const double blencp, const cs_real_t pi, const cs_real_t pj, const cs_real_t pir, const cs_real_t pjr, cs_real_t *pifri, cs_real_t *pifrj, cs_real_t *pjfri, cs_real_t *pjfrj)
Blend face values for a centered or SOLU scheme with face values for an upwind scheme.
Definition: cs_convection_diffusion.h:749
static void cs_i_cd_unsteady_slope_test(bool *upwind_switch, const int ircflp, const int ischcp, const double blencp, const cs_real_t weight, const cs_real_t i_dist, const cs_real_t i_face_surf, const cs_real_3_t cell_ceni, const cs_real_3_t cell_cenj, const cs_real_3_t i_face_normal, const cs_real_3_t i_face_cog, const cs_real_3_t dijpf, const cs_real_t i_massflux, const cs_real_3_t gradi, const cs_real_3_t gradj, const cs_real_3_t grdpai, const cs_real_3_t grdpaj, const cs_real_t pi, const cs_real_t pj, 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 unsteady algorithm...
Definition: cs_convection_diffusion.h:2279
cs_lnum_t n_cells
Definition: cs_mesh.h:72
static void cs_i_compute_quantities_tensor(const int ircflp, const double pnd, 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_3_t dijpf, const cs_real_63_t gradi, const cs_real_63_t gradj, const cs_real_6_t pi, const cs_real_6_t pj, cs_real_t recoi[6], cs_real_t recoj[6], cs_real_t pip[6], cs_real_t pjp[6])
Reconstruct values in I' and J'.
Definition: cs_convection_diffusion.h:272
int n_init_perio
Definition: cs_mesh.h:111
cs_lnum_t * group_index
Definition: cs_numbering.h:93
static void cs_b_relax_c_val_tensor(const double relaxp, const cs_real_6_t pi, const cs_real_6_t pia, const cs_real_6_t recoi, cs_real_t pir[6], cs_real_t pipr[6])
Compute relaxed values at boundary cell i.
Definition: cs_convection_diffusion.h:2691
static 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.h:2667
static void cs_i_cd_steady_upwind_tensor(const int ircflp, const cs_real_t relaxp, 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_3_t dijpf, const cs_real_63_t gradi, const cs_real_63_t gradj, const cs_real_6_t pi, const cs_real_6_t pj, const cs_real_6_t pia, const cs_real_6_t pja, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6], cs_real_t pip[6], cs_real_t pjp[6], cs_real_t pipr[6], cs_real_t pjpr[6])
Handle preparation of internal face values for the fluxes computation in case of a steady algorithm a...
Definition: cs_convection_diffusion.h:1100
cs_real_t * cell_cen
Definition: cs_mesh_quantities.h:53
static void cs_b_imposed_conv_flux_cons(const int iconvp, const int inc, const int ifaccp, const cs_int_t bc_type, const int icvfli, const cs_real_t pir, const cs_real_t pipr, const cs_real_t coefap, const cs_real_t coefbp, const cs_real_t coface, const cs_real_t cofbce, const cs_real_t xcpp, const cs_real_t b_massflux, cs_real_t *flux)
Add convective flux (conservative formulation) to flux at boundary face. The convective flux can be e...
Definition: cs_convection_diffusion.h:2842
void diften(const cs_int_t *const idtvar, const cs_int_t *const f_id, const cs_var_cal_opt_t *const var_cal_opt, const cs_int_t *const inc, const cs_int_t *const iccocg, cs_real_t pvar[], const cs_real_t pvara[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t viscel[], const cs_real_2_t weighf[], const cs_real_t weighb[], cs_real_t rhs[])
Definition: cs_convection_diffusion.c:361
static void cs_i_cd_steady_slope_test(bool *upwind_switch, const int ircflp, const int ischcp, const double relaxp, const double blencp, const cs_real_t weight, const cs_real_t i_dist, const cs_real_t i_face_surf, const cs_real_3_t cell_ceni, const cs_real_3_t cell_cenj, const cs_real_3_t i_face_normal, const cs_real_3_t i_face_cog, const cs_real_3_t dijpf, const cs_real_t i_massflux, const cs_real_3_t gradi, const cs_real_3_t gradj, const cs_real_3_t grdpai, const cs_real_3_t grdpaj, 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.h:1899
cs_mesh_quantities_t * cs_glob_mesh_quantities
Definition: cs_parameters.h:130
static void cs_slope_test_gradient_tensor(const int f_id, const int inc, const cs_halo_type_t halo_type, cs_real_63_t *grad, cs_real_63_t *grdpa, cs_real_6_t *pvar, const cs_real_6_t *coefa, const cs_real_66_t *coefb, const cs_real_t *i_massflux)
Compute the upwind gradient used in the slope tests.
Definition: cs_convection_diffusion.h:3425
cs_real_t * i_face_normal
Definition: cs_mesh_quantities.h:57
void bilsc6(const cs_int_t *const idtvar, const cs_int_t *const f_id, const cs_var_cal_opt_t *const var_cal_opt, const cs_int_t *const icvflb, const cs_int_t *const inc, const cs_int_t *const ifaccp, cs_real_6_t pvar[], const cs_real_6_t pvara[], const cs_int_t bc_type[], const cs_int_t icvfli[], const cs_real_6_t coefa[], const cs_real_66_t coefb[], const cs_real_6_t cofaf[], const cs_real_66_t cofbf[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t rhs[])
Definition: cs_convection_diffusion.c:263
Definition: cs_mesh_quantities.h:51
cs_halo_type_t
Definition: cs_halo.h:49
void diftnts(const cs_int_t *const idtvar, const cs_int_t *const f_id, const cs_var_cal_opt_t *const var_cal_opt, const cs_int_t *const inc, cs_real_6_t pvar[], const cs_real_6_t pvara[], const cs_int_t bc_type[], const cs_real_6_t coefa[], const cs_real_66_t coefb[], const cs_real_6_t cofaf[], const cs_real_66_t cofbf[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t viscel[], const cs_real_2_t weighf[], const cs_real_t weighb[], cs_real_6_t rhs[])
Definition: cs_convection_diffusion.c:448
static void cs_upwind_f_val(const cs_real_t pi, const cs_real_t pj, const cs_real_t pir, const cs_real_t pjr, cs_real_t *pifri, cs_real_t *pifrj, cs_real_t *pjfri, cs_real_t *pjfrj)
Prepare value at face ij by using an upwind scheme.
Definition: cs_convection_diffusion.h:496
void cs_anisotropic_diffusion_potential(const int f_id, const cs_mesh_t *m, cs_mesh_quantities_t *fvq, int init, int inc, int imrgra, int iccocg, int nswrgp, int imligp, int ircflp, int iphydp, int iwarnp, double epsrgp, double climgp, double extrap, cs_real_3_t *restrict frcxt, cs_real_t *restrict pvar, const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t *restrict viscel, const cs_real_2_t weighf[], const cs_real_t weighb[], cs_real_t *restrict diverg)
Add the explicit part of the divergence of the mass flux due to the pressure gradient (routine analog...
Definition: cs_convection_diffusion.c:8391
static void cs_slope_test(const cs_real_t pi, const cs_real_t pj, const cs_real_t distf, const cs_real_t srfan, const cs_real_3_t i_face_normal, const cs_real_3_t gradi, const cs_real_3_t gradj, const cs_real_3_t grdpai, const cs_real_3_t grdpaj, const cs_real_t i_massflux, double *testij, double *tesqck)
Compute slope test criteria at internal face between cell i and j.
Definition: cs_convection_diffusion.h:79
int cs_lnum_2_t[2]
vector of 2 local mesh-entity ids
Definition: cs_defs.h:301
static void cs_i_cd_unsteady_upwind(const int ircflp, 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_3_t dijpf, 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 *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 an unsteady algorith...
Definition: cs_convection_diffusion.h:1189
static void cs_b_compute_quantities_tensor(const cs_real_3_t diipb, const cs_real_63_t gradi, const int ircflp, cs_real_t recoi[6])
Reconstruct values in I' at boundary cell i.
Definition: cs_convection_diffusion.h:2642
static void cs_slope_test_tensor(const cs_real_6_t pi, const cs_real_6_t pj, const cs_real_t distf, const cs_real_t srfan, const cs_real_3_t i_face_normal, const cs_real_63_t gradi, const cs_real_63_t gradj, const cs_real_63_t grdpai, const cs_real_63_t grdpaj, const cs_real_t i_massflux, cs_real_t testij[6], cs_real_t tesqck[6])
Compute slope test criteria at internal face between cell i and j.
Definition: cs_convection_diffusion.h:143
void cs_face_diffusion_potential(const int f_id, const cs_mesh_t *m, cs_mesh_quantities_t *fvq, int init, int inc, int imrgra, int iccocg, int nswrgp, int imligp, int iphydp, int iwarnp, double epsrgp, double climgp, double extrap, cs_real_3_t *restrict frcxt, cs_real_t *restrict pvar, const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t viselx[], const cs_real_t visely[], const cs_real_t viselz[], cs_real_t *restrict i_massflux, cs_real_t *restrict b_massflux)
Update the face mass flux with the face pressure (or pressure increment, or pressure double increment...
Definition: cs_convection_diffusion.c:7202
static void cs_b_cd_unsteady_tensor(const int ircflp, const cs_real_3_t diipb, const cs_real_63_t gradi, const cs_real_6_t pi, cs_real_t pir[6], cs_real_t pipr[6])
Handle preparation of boundary face values for the flux computation in case of a steady algorithm...
Definition: cs_convection_diffusion.h:3228
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:307
void cs_diffusion_potential(const int f_id, const cs_mesh_t *m, cs_mesh_quantities_t *fvq, int init, int inc, int imrgra, int iccocg, int nswrgp, int imligp, int iphydp, int iwarnp, double epsrgp, double climgp, double extrap, cs_real_3_t *restrict frcxt, cs_real_t *restrict pvar, const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t viselx[], const cs_real_t visely[], const cs_real_t viselz[], cs_real_t *restrict diverg)
Update the cell mass flux divergence with the face pressure (or pressure increment, or pressure double increment) gradient.
Definition: cs_convection_diffusion.c:7989
int n_groups
Definition: cs_numbering.h:85
void itrmav(const cs_int_t *const f_id, const cs_int_t *const init, const cs_int_t *const inc, const cs_int_t *const imrgra, const cs_int_t *const iccocg, const cs_int_t *const nswrgp, const cs_int_t *const imligp, const cs_int_t *const ircflp, const cs_int_t *const iphydp, const cs_int_t *const iwarnp, const cs_real_t *const epsrgp, const cs_real_t *const climgp, const cs_real_t *const extrap, cs_real_3_t frcxt[], cs_real_t pvar[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t viscel[], const cs_real_2_t weighf[], const cs_real_t weighb[], cs_real_t i_massflux[], cs_real_t b_massflux[])
Definition: cs_convection_diffusion.c:556
Definition: cs_parameters.h:53
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:292
static void cs_b_compute_quantities(const cs_real_3_t diipb, const cs_real_3_t gradi, const int ircflp, cs_real_t *recoi)
Reconstruct values in I' at boundary cell i.
Definition: cs_convection_diffusion.h:2620
void itrgrv(const cs_int_t *const f_id, const cs_int_t *const init, const cs_int_t *const inc, const cs_int_t *const imrgra, const cs_int_t *const iccocg, const cs_int_t *const nswrgp, const cs_int_t *const imligp, const cs_int_t *const ircflp, const cs_int_t *const iphydp, const cs_int_t *const iwarnp, const cs_real_t *const epsrgp, const cs_real_t *const climgp, const cs_real_t *const extrap, cs_real_3_t frcxt[], cs_real_t pvar[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t viscel[], const cs_real_2_t weighf[], const cs_real_t weighb[], cs_real_t diverg[])
Definition: cs_convection_diffusion.c:686
void itrmas(const cs_int_t *const f_id, const cs_int_t *const init, const cs_int_t *const inc, const cs_int_t *const imrgra, const cs_int_t *const iccocg, const cs_int_t *const nswrgp, const cs_int_t *const imligp, const cs_int_t *const iphydp, const cs_int_t *const iwarnp, const cs_real_t *const epsrgp, const cs_real_t *const climgp, const cs_real_t *const extrap, cs_real_3_t frcxt[], cs_real_t pvar[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t viselx[], const cs_real_t visely[], const cs_real_t viselz[], cs_real_t i_massflux[], cs_real_t b_massflux[])
Definition: cs_convection_diffusion.c:491
static void cs_i_conv_flux_tensor(const int iconvp, const cs_real_6_t pi, const cs_real_6_t pj, const cs_real_6_t pifri, const cs_real_6_t pifrj, const cs_real_6_t pjfri, const cs_real_6_t pjfrj, const cs_real_t i_massflux, cs_real_t fluxi[6], cs_real_t fluxj[6])
Add convective fluxes (substracting the mass accumulation from them) to fluxes at face ij...
Definition: cs_convection_diffusion.h:859
static void cs_b_imposed_conv_flux(const int iconvp, const int inc, const int ifaccp, const cs_int_t bc_type, const int icvfli, 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 coface, const cs_real_t cofbce, const cs_real_t b_massflux, cs_real_t *flux)
Add convective flux (substracting the mass accumulation from it) to flux at boundary face...
Definition: cs_convection_diffusion.h:2773
#define END_C_DECLS
Definition: cs_defs.h:430
static void cs_i_cd_steady_upwind(const int ircflp, const cs_real_t relaxp, 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_3_t dijpf, const cs_real_3_t gradi, const cs_real_3_t gradj, 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.h:1008
static void cs_centered_f_val(const double pnd, const cs_real_t pip, const cs_real_t pjp, const cs_real_t pipr, const cs_real_t pjpr, cs_real_t *pifri, cs_real_t *pifrj, cs_real_t *pjfri, cs_real_t *pjfrj)
Prepare value at face ij by using a centered scheme.
Definition: cs_convection_diffusion.h:563
double cs_real_t
Definition: cs_defs.h:296
static void cs_i_compute_quantities(const int ircflp, const double pnd, 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_3_t dijpf, 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.h:211
cs_numbering_t * b_face_numbering
Definition: cs_mesh.h:131
cs_numbering_t * i_face_numbering
Definition: cs_mesh.h:130
#define CS_PROCF(x, y)
Definition: cs_defs.h:453
void bilsct(const cs_int_t *const idtvar, const cs_int_t *const f_id, const cs_var_cal_opt_t *const var_cal_opt, const cs_int_t *const inc, const cs_int_t *const iccocg, const cs_int_t *const ifaccp, cs_real_t pvar[], const cs_real_t pvara[], const cs_int_t bc_type[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t xcpp[], cs_real_t rhs[])
Definition: cs_convection_diffusion.c:313
int n_threads
Definition: cs_numbering.h:84
void cs_face_anisotropic_diffusion_potential(const int f_id, const cs_mesh_t *m, cs_mesh_quantities_t *fvq, int init, int inc, int imrgra, int iccocg, int nswrgp, int imligp, int ircflp, int iphydp, int iwarnp, double epsrgp, double climgp, double extrap, cs_real_3_t *restrict frcxt, cs_real_t *restrict pvar, const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t *restrict viscel, const cs_real_2_t weighf[], const cs_real_t weighb[], cs_real_t *restrict i_massflux, cs_real_t *restrict b_massflux)
Add the explicit part of the pressure gradient term to the mass flux in case of anisotropic diffusion...
Definition: cs_convection_diffusion.c:7565
static void cs_i_cd_unsteady(const int ircflp, 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_3_t dijpf, 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 *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 unsteady algorithm...
Definition: cs_convection_diffusion.h:1633
cs_real_t cs_real_33_t[3][3]
3x3 matrix of floating-point values
Definition: cs_defs.h:311
cs_real_t * b_face_normal
Definition: cs_mesh_quantities.h:59
int have_rotation_perio
Definition: cs_mesh.h:114
static void cs_i_cd_unsteady_tensor(const int ircflp, 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_3_t dijpf, const cs_real_63_t gradi, const cs_real_63_t gradj, const cs_real_6_t pi, const cs_real_6_t pj, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6], cs_real_t pip[6], cs_real_t pjp[6], cs_real_t pipr[6], cs_real_t pjpr[6])
Handle preparation of internal face values for the fluxes computation in case of an unsteady algorith...
Definition: cs_convection_diffusion.h:1761
cs_real_t * diipb
Definition: cs_mesh_quantities.h:75
void cs_halo_perio_sync_var_vect(const cs_halo_t *halo, cs_halo_type_t sync_mode, cs_real_t var[], int incvar)
Definition: cs_halo_perio.c:640
cs_real_t cs_real_63_t[6][3]
Definition: cs_defs.h:316
static void cs_b_diff_flux_tensor(const int idiffp, const int inc, const cs_real_6_t pipr, const cs_real_6_t cofaf, const cs_real_66_t cofbf, const cs_real_t b_visc, cs_real_t flux[6])
Add diffusive flux to flux at boundary face.
Definition: cs_convection_diffusion.h:3081
cs_real_t * cell_vol
Definition: cs_mesh_quantities.h:54
static void cs_i_cd_steady(const int ircflp, const int ischcp, const double relaxp, 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_3_t dijpf, const cs_real_3_t gradi, const cs_real_3_t gradj, 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.h:1360
void bilsc4(const cs_int_t *const idtvar, const cs_int_t *const f_id, const cs_var_cal_opt_t *const var_cal_opt, const cs_int_t *const icvflb, const cs_int_t *const inc, const cs_int_t *const ifaccp, const cs_int_t *const ivisep, cs_real_3_t pvar[], const cs_real_3_t pvara[], const cs_int_t bc_type[], const cs_int_t icvfli[], const cs_real_3_t coefav[], const cs_real_33_t coefbv[], const cs_real_3_t cofafv[], const cs_real_33_t cofbfv[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t secvif[], cs_real_3_t rhs[])
Definition: cs_convection_diffusion.c:211
cs_lnum_2_t * i_face_cells
Definition: cs_mesh.h:86
static void cs_b_diff_flux(const int idiffp, 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.h:3055
double precision, dimension(:,:), pointer dijpf
Definition: mesh.f90:208
static void cs_b_cd_steady_tensor(const int ircflp, const double relaxp, const cs_real_3_t diipb, const cs_real_63_t gradi, const cs_real_6_t pi, const cs_real_6_t pia, cs_real_t pir[6], cs_real_t pipr[6])
Handle preparation of boundary face values for the flux computation in case of a steady algorithm...
Definition: cs_convection_diffusion.h:3156
static void cs_upwind_f_val_tensor(const cs_real_6_t pi, const cs_real_6_t pj, const cs_real_6_t pir, const cs_real_6_t pjr, cs_real_t pifri[6], cs_real_t pifrj[6], cs_real_t pjfri[6], cs_real_t pjfrj[6])
Prepare value at face ij by using an upwind scheme.
Definition: cs_convection_diffusion.h:527
static void cs_i_no_relax_c_val(const cs_real_t pi, const cs_real_t pj, const cs_real_t pip, const cs_real_t pjp, cs_real_t *pir, cs_real_t *pjr, cs_real_t *pipr, cs_real_t *pjpr)
Copy reconstructed or not cell values for consistency with relaxation case.
Definition: cs_convection_diffusion.h:428
static void cs_solu_f_val(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_3_t gradi, const cs_real_3_t gradj, const cs_real_t pi, const cs_real_t pj, const cs_real_t pir, const cs_real_t pjr, cs_real_t *pifri, cs_real_t *pifrj, cs_real_t *pjfri, cs_real_t *pjfrj)
Prepare value at face ij by using a Second Order Linear Upwind scheme.
Definition: cs_convection_diffusion.h:635
void cs_convection_diffusion_scalar(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int icvflb, int inc, int iccocg, int ifaccp, cs_real_t *restrict pvar, const cs_real_t *restrict pvara, const cs_int_t bc_type[], const cs_int_t icvfli[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_t *restrict rhs)
Add the explicit part of the convection/diffusion terms of a standard transport equation of a scalar ...
Definition: cs_convection_diffusion.c:808
void cs_convection_diffusion_thermal(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int inc, int iccocg, int ifaccp, cs_real_t *restrict pvar, const cs_real_t *restrict pvara, const cs_int_t bc_type[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t xcpp[], cs_real_t *restrict rhs)
Add the explicit part of the convection/diffusion terms of a transport equation of a scalar field su...
Definition: cs_convection_diffusion.c:4294
static void cs_i_diff_flux_tensor(const int idiffp, const cs_real_6_t pip, const cs_real_6_t pjp, const cs_real_6_t pipr, const cs_real_6_t pjpr, const cs_real_t i_visc, cs_real_t fluxi[6], cs_real_t fluxj[6])
Add diffusive fluxes to fluxes at face ij.
Definition: cs_convection_diffusion.h:963
cs_real_t * i_face_cog
Definition: cs_mesh_quantities.h:65
static void cs_i_no_relax_c_val_tensor(const cs_real_6_t pi, const cs_real_6_t pj, const cs_real_6_t pip, const cs_real_6_t pjp, cs_real_t pir[6], cs_real_t pjr[6], cs_real_t pipr[6], cs_real_t pjpr[6])
Copy reconstructed or not cell values for consistency with relaxation case.
Definition: cs_convection_diffusion.h:461
static void cs_i_conv_flux_cons(const int iconvp, 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 (cons. formulation) to fluxes at face ij.
Definition: cs_convection_diffusion.h:901
void itrgrp(const cs_int_t *const f_id, const cs_int_t *const init, const cs_int_t *const inc, const cs_int_t *const imrgra, const cs_int_t *const iccocg, const cs_int_t *const nswrgp, const cs_int_t *const imligp, const cs_int_t *const iphydp, const cs_int_t *const iwarnp, const cs_real_t *const epsrgp, const cs_real_t *const climgp, const cs_real_t *const extrap, cs_real_3_t frcxt[], cs_real_t pvar[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t viselx[], const cs_real_t visely[], const cs_real_t viselz[], cs_real_t diverg[])
Definition: cs_convection_diffusion.c:623
integer(c_int), pointer, save imrgra
type of gradient reconstruction
Definition: optcal.f90:266
cs_halo_t * halo
Definition: cs_mesh.h:127
cs_lnum_t * b_face_cells
Definition: cs_mesh.h:87
integer, save ifaccp
indicator coupling face / face only
Definition: cplsat.f90:45