7.1
general 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-2021 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 #include "cs_math.h"
39 #include "cs_mesh_quantities.h"
40 #include "cs_parameters.h"
41 #include "cs_field_pointer.h"
42 
43 /*----------------------------------------------------------------------------*/
44 
45 BEGIN_C_DECLS
46 
47 /*=============================================================================
48  * Macro definitions
49  *============================================================================*/
50 
51 /*============================================================================
52  * Type definition
53  *============================================================================*/
54 
55 /*----------------------------------------------------------------------------
56  * NVD/TVD Advection Scheme
57  *----------------------------------------------------------------------------*/
58 
59 typedef enum {
60 
61  CS_NVD_GAMMA = 0, /* GAMMA */
62  CS_NVD_SMART = 1, /* SMART */
63  CS_NVD_CUBISTA = 2, /* CUBISTA */
64  CS_NVD_SUPERBEE = 3, /* SUPERBEE */
65  CS_NVD_MUSCL = 4, /* MUSCL */
66  CS_NVD_MINMOD = 5, /* MINMOD */
67  CS_NVD_CLAM = 6, /* CLAM */
68  CS_NVD_STOIC = 7, /* STOIC */
69  CS_NVD_OSHER = 8, /* OSHER */
70  CS_NVD_WASEB = 9, /* WASEB */
71  CS_NVD_VOF_HRIC = 10, /* M-HRIC for VOF */
72  CS_NVD_VOF_CICSAM = 11, /* M-CICSAM for VOF */
73  CS_NVD_VOF_STACS = 12, /* STACS for VOF */
74  CS_NVD_N_TYPES = 13 /* number of NVD schemes */
75 
76 } cs_nvd_type_t;
77 
78 /*============================================================================
79  * Global variables
80  *============================================================================*/
81 
82 /*============================================================================
83  * Public inlined function
84  *============================================================================*/
85 
86 /*----------------------------------------------------------------------------
87  * Synchronize halos for scalar variables.
88  *
89  * parameters:
90  * m <-- pointer to associated mesh structure
91  * pvar <-> variable
92  *----------------------------------------------------------------------------*/
93 
94 inline static void
95 cs_sync_scalar_halo(const cs_mesh_t *m,
96  cs_real_t pvar[])
97 {
98  if (m->halo != NULL)
99  cs_halo_sync_var(m->halo, CS_HALO_STANDARD, pvar);
100 }
101 
102 /*----------------------------------------------------------------------------
103  * Return pointer to slope test indicator field values if active.
104  *
105  * parameters:
106  * f_id <-- field id (or -1)
107  * var_cal_opt <-- variable calculation options
108  *
109  * return:
110  * pointer to local values array, or NULL;
111  *----------------------------------------------------------------------------*/
112 
113 inline static cs_real_t *
114 cs_get_v_slope_test(int f_id,
115  const cs_var_cal_opt_t var_cal_opt)
116 {
117  const int iconvp = var_cal_opt.iconv;
118  const int isstpp = var_cal_opt.isstpc;
119  const double blencp = var_cal_opt.blencv;
120 
121  cs_real_t *v_slope_test = NULL;
122 
123  if (f_id > -1 && iconvp > 0 && blencp > 0. && isstpp == 0) {
124 
125  static int _k_slope_test_f_id = -1;
126 
127  cs_field_t *f = cs_field_by_id(f_id);
128 
129  int f_track_slope_test_id = -1;
130 
131  if (_k_slope_test_f_id < 0)
132  _k_slope_test_f_id = cs_field_key_id_try("slope_test_upwind_id");
133  if (_k_slope_test_f_id > -1 && isstpp == 0)
134  f_track_slope_test_id = cs_field_get_key_int(f, _k_slope_test_f_id);
135 
136  if (f_track_slope_test_id > -1)
137  v_slope_test = (cs_field_by_id(f_track_slope_test_id))->val;
138 
139  if (v_slope_test != NULL) {
140  const cs_lnum_t n_cells_ext = cs_glob_mesh->n_cells_with_ghosts;
141 # pragma omp parallel for
142  for (cs_lnum_t cell_id = 0; cell_id < n_cells_ext; cell_id++)
143  v_slope_test[cell_id] = 0.;
144  }
145 
146  }
147 
148  return v_slope_test;
149 }
150 
151 /*----------------------------------------------------------------------------
152  * Compute the local cell Courant number as the maximum of all cell face based
153  * Courant number at each cell.
154  *
155  * parameters:
156  * f_id <-- field id (or -1)
157  * courant --> cell Courant number
158  */
159 /*----------------------------------------------------------------------------*/
160 
161 inline static void
162 cs_cell_courant_number(const int f_id,
163  cs_real_t *courant)
164 {
165  const cs_mesh_t *m = cs_glob_mesh;
166  const cs_mesh_quantities_t *fvq = cs_glob_mesh_quantities;
167 
168  const cs_lnum_t n_cells_ext = m->n_cells_with_ghosts;
169  const int n_i_groups = m->i_face_numbering->n_groups;
170  const int n_i_threads = m->i_face_numbering->n_threads;
171  const int n_b_groups = m->b_face_numbering->n_groups;
172  const int n_b_threads = m->b_face_numbering->n_threads;
173  const cs_lnum_t *restrict i_group_index = m->i_face_numbering->group_index;
174  const cs_lnum_t *restrict b_group_index = m->b_face_numbering->group_index;
175 
176  const cs_lnum_2_t *restrict i_face_cells
177  = (const cs_lnum_2_t *restrict)m->i_face_cells;
178  const cs_lnum_t *restrict b_face_cells
179  = (const cs_lnum_t *restrict)m->b_face_cells;
180 
181  const cs_real_t *restrict vol
182  = (cs_real_t *restrict)fvq->cell_vol;
183 
184  cs_field_t *f = cs_field_by_id(f_id);
185  const int kimasf = cs_field_key_id("inner_mass_flux_id");
186  const int kbmasf = cs_field_key_id("boundary_mass_flux_id");
187  const cs_real_t *restrict i_massflux
188  = cs_field_by_id( cs_field_get_key_int(f, kimasf) )->val;
189  const cs_real_t *restrict b_massflux
190  = cs_field_by_id( cs_field_get_key_int(f, kbmasf) )->val;
191 
192  const cs_real_t *restrict dt
193  = (const cs_real_t *restrict)CS_F_(dt)->val;
194 
195  /* Initialisation */
196 
197 # pragma omp parallel for
198  for (cs_lnum_t ii = 0; ii < n_cells_ext; ii++) {
199  courant[ii] = 0.;
200  }
201 
202  /* ---> Contribution from interior faces */
203 
204  cs_real_t cnt;
205 
206  for (int g_id = 0; g_id < n_i_groups; g_id++) {
207 # pragma omp parallel for
208  for (int t_id = 0; t_id < n_i_threads; t_id++) {
209  for (cs_lnum_t face_id = i_group_index[(t_id*n_i_groups + g_id)*2];
210  face_id < i_group_index[(t_id*n_i_groups + g_id)*2 + 1];
211  face_id++) {
212  cs_lnum_t ii = i_face_cells[face_id][0];
213  cs_lnum_t jj = i_face_cells[face_id][1];
214 
215  cnt = CS_ABS(i_massflux[face_id])*dt[ii]/vol[ii];
216  courant[ii] = CS_MAX(courant[ii], cnt);
217 
218  cnt = CS_ABS(i_massflux[face_id])*dt[jj]/vol[jj];
219  courant[jj] = CS_MAX(courant[jj], cnt);
220  }
221  }
222  }
223 
224  /* ---> Contribution from boundary faces */
225 
226  for (int g_id = 0; g_id < n_b_groups; g_id++) {
227 # pragma omp parallel for
228  for (int t_id = 0; t_id < n_b_threads; t_id++) {
229  for (cs_lnum_t face_id = b_group_index[(t_id*n_b_groups + g_id)*2];
230  face_id < b_group_index[(t_id*n_b_groups + g_id)*2 + 1];
231  face_id++) {
232  cs_lnum_t ii = b_face_cells[face_id];
233 
234  cnt = CS_ABS(b_massflux[face_id])*dt[ii]/vol[ii];
235  courant[ii] = CS_MAX(courant[ii], cnt);
236  }
237  }
238  }
239 }
240 
241 /*----------------------------------------------------------------------------*/
252 /*----------------------------------------------------------------------------*/
253 
254 inline static cs_real_t
255 cs_nvd_scheme_scalar(const cs_nvd_type_t limiter,
256  const cs_real_t nvf_p_c,
257  const cs_real_t nvf_r_f,
258  const cs_real_t nvf_r_c)
259 {
260  cs_real_t nvf_p_f;
261 
262  cs_real_t beta_m, rfc, r1f, r1, r2, r3, b1, b2;
263 
264  switch (limiter) {
265  case CS_NVD_GAMMA: /* Gamma scheme */
266  beta_m = 0.1; /* in [0.1, 0.5] */
267  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
268 
269  if (nvf_p_c < beta_m) {
270  nvf_p_f = nvf_p_c*(1.+rfc*(1.-nvf_p_c)/beta_m);
271  } else {
272  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
273 
274  nvf_p_f = r1f*nvf_p_c+rfc;
275  }
276 
277  break;
278 
279  case CS_NVD_SMART: /* SMART scheme */
280  if (nvf_p_c < (nvf_r_c/3.)) {
281  r1 = nvf_r_f*(1.-3.*nvf_r_c+2.*nvf_r_f);
282  r2 = nvf_r_c*(1.-nvf_r_c);
283 
284  nvf_p_f = nvf_p_c*r1/r2;
285  } else if (nvf_p_c <= (nvf_r_c*(1.+nvf_r_f-nvf_r_c)/nvf_r_f)) {
286  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
287  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
288 
289  nvf_p_f = nvf_r_f*(r1f*nvf_p_c/nvf_r_c + rfc);
290  } else {
291  nvf_p_f = 1.;
292  }
293 
294  break;
295 
296  case CS_NVD_CUBISTA: /* CUBISTA scheme */
297  if (nvf_p_c < (3.*nvf_r_c/4.)) {
298  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
299 
300  nvf_p_f = nvf_r_f*(1.+rfc/3.)*nvf_p_c/nvf_r_c;
301  } else if (nvf_p_c <= (nvf_r_c*(1.+2.*(nvf_r_f-nvf_r_c))/(2.*nvf_r_f-nvf_r_c))) {
302  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
303  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
304 
305  nvf_p_f = nvf_r_f*(r1f*nvf_p_c/nvf_r_c+rfc);
306  } else {
307  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
308 
309  nvf_p_f = 1.-.5*r1f*(1.-nvf_p_c);
310  }
311 
312  break;
313 
314  case CS_NVD_SUPERBEE: /* SuperBee scheme */
315  if (nvf_p_c < (nvf_r_c/(2.-nvf_r_c))) {
316  nvf_p_f = (2.*nvf_r_f-nvf_r_c)*nvf_p_c/nvf_r_c;
317  } else if (nvf_p_c < nvf_r_c) {
318  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
319  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
320 
321  nvf_p_f = r1f*nvf_p_c+rfc;
322  } else if (nvf_p_c < (nvf_r_c/nvf_r_f)) {
323  nvf_p_f = nvf_r_f*nvf_p_c/nvf_r_c;
324  } else {
325  nvf_p_f = 1.;
326  }
327 
328  break;
329 
330  case CS_NVD_MUSCL: /* MUSCL scheme */
331  if (nvf_p_c < (.5*nvf_r_c)) {
332  nvf_p_f = (2.*nvf_r_f-nvf_r_c)*nvf_p_c/nvf_r_c;
333  } else if (nvf_p_c < (1.+nvf_r_c-nvf_r_f)) {
334  nvf_p_f = nvf_p_c+nvf_r_f-nvf_r_c;
335  } else {
336  nvf_p_f = 1.;
337  }
338 
339  break;
340 
341  case CS_NVD_MINMOD: /* MINMOD scheme */
342  if (nvf_p_c < nvf_r_c) {
343  nvf_p_f = nvf_r_f*nvf_p_c/nvf_r_c;
344  } else {
345  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
346  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
347 
348  nvf_p_f = r1f*nvf_p_c+rfc;
349  }
350 
351  break;
352 
353  case CS_NVD_CLAM: /* CLAM scheme */
354  r1 = nvf_r_c*nvf_r_c-nvf_r_f;
355  r2 = nvf_r_c*(nvf_r_c-1.);
356  r3 = nvf_r_f-nvf_r_c;
357 
358  nvf_p_f = nvf_p_c*(r1+r3*nvf_p_c)/r2;
359 
360  break;
361 
362  case CS_NVD_STOIC: /* STOIC scheme */
363  b1 = (nvf_r_c-nvf_r_f)*nvf_r_c;
364  b2 = nvf_r_c+nvf_r_f+2.*nvf_r_f*nvf_r_f-4.*nvf_r_f*nvf_r_c;
365 
366  if (nvf_p_c < (b1/b2)) {
367  r1 = -nvf_r_f*(1.-3.*nvf_r_c+2.*nvf_r_f);
368  r2 = nvf_r_c*(nvf_r_c-1.);
369 
370  nvf_p_f = nvf_p_c*r1/r2;
371  } else if (nvf_p_c < nvf_r_c) {
372  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
373  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
374 
375  nvf_p_f = r1f*nvf_p_c+rfc;
376  } else if (nvf_p_c < (nvf_r_c*(1.+nvf_r_f-nvf_r_c)/nvf_r_f)) {
377  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
378  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
379 
380  nvf_p_f = nvf_r_f*(nvf_p_c*r1f/nvf_r_c+rfc);
381  } else {
382  nvf_p_f = 1.;
383  }
384 
385  break;
386 
387  case CS_NVD_OSHER: /* OSHER scheme */
388  if (nvf_p_c < (nvf_r_c/nvf_r_f)) {
389  nvf_p_f = nvf_r_f*nvf_p_c/nvf_r_c;
390  } else {
391  nvf_p_f = 1.;
392  }
393 
394  break;
395 
396  case CS_NVD_WASEB: /* WASEB scheme */
397  r1 = nvf_r_c*nvf_r_f*(nvf_r_f-nvf_r_c);
398  r2 = 2.*nvf_r_c*(1.-nvf_r_c)-nvf_r_f*(1.-nvf_r_f);
399 
400  if (nvf_p_c < (r1/r2)) {
401  nvf_p_f = 2.*nvf_p_c;
402  } else if (nvf_p_c <= (nvf_r_c*(1.+nvf_r_f-nvf_r_c)/nvf_r_f)) {
403  rfc = (nvf_r_f-nvf_r_c)/(1.-nvf_r_c);
404  r1f = (1.-nvf_r_f)/(1.-nvf_r_c);
405 
406  nvf_p_f = nvf_r_f*(nvf_p_c*r1f/nvf_r_c+rfc);
407  } else {
408  nvf_p_f = 1.;
409  }
410 
411  break;
412 
413  default: /* Upwinding */
414  nvf_p_f = nvf_p_c;
415  break;
416  }
417 
418  return nvf_p_f;
419 }
420 
421 /*----------------------------------------------------------------------------*/
437 /*----------------------------------------------------------------------------*/
438 
439 inline static cs_real_t
440 cs_nvd_vof_scheme_scalar(const cs_nvd_type_t limiter,
441  const cs_real_3_t i_face_normal,
442  const cs_real_t nvf_p_c,
443  const cs_real_t nvf_r_f,
444  const cs_real_t nvf_r_c,
445  const cs_real_3_t gradv_c,
446  const cs_real_t c_courant)
447 {
448  assert(limiter >= CS_NVD_VOF_HRIC);
449 
450  cs_real_t nvf_p_f;
451  cs_real_t blend, high_order, low_order, ratio;
452 
453  /* Compute gradient angle indicator */
454  cs_real_t dotp = CS_ABS(cs_math_3_dot_product(gradv_c, i_face_normal));
455  cs_real_t sgrad = cs_math_3_norm(gradv_c);
456  cs_real_t snorm = cs_math_3_norm(i_face_normal);
457  cs_real_t denom = snorm*sgrad;
458 
459  if (limiter == CS_NVD_VOF_HRIC) { /* M-HRIC scheme */
460  /* High order scheme : Bounded Downwind */
461  if (nvf_p_c <= .5) {
462  high_order = 2.*nvf_p_c;
463  } else {
464  high_order = 1.;
465  }
466 
467  /* Low order scheme : MUSCL */
468  low_order = cs_nvd_scheme_scalar(CS_NVD_MUSCL,
469  nvf_p_c,
470  nvf_r_f,
471  nvf_r_c);
472 
473  /* Compute the blending factor */
474  if (denom <= (cs_math_epzero*dotp)) {
475  blend = 1.;
476  } else {
477  ratio = dotp/denom;
478  blend = CS_MIN(1., pow(ratio, .5));
479  }
480 
481  /* Blending */
482  nvf_p_f = blend*high_order + (1.-blend)*low_order;
483 
484  /* Extra blending due to the cell Courant number */
485  if (c_courant < .7 && c_courant > .3) {
486  nvf_p_f = nvf_p_f + (nvf_p_f - low_order)*(.7 - c_courant )/.4;
487  } else if (c_courant >= .7) {
488  nvf_p_f = low_order;
489  }
490  } else if (limiter == CS_NVD_VOF_CICSAM) { /* M-CICSAM scheme */
491  /* High order scheme : HYPER-C + SUPERBEE */
492  if (c_courant <= .3) {
493  high_order = CS_MIN(1., nvf_p_c/(c_courant+cs_math_epzero));
494  } else if (c_courant <= .6) {
495  high_order = CS_MIN(1., nvf_p_c/.3);
496  } else if (c_courant <= .7) {
497  cs_real_t superbee = cs_nvd_scheme_scalar(CS_NVD_SUPERBEE,
498  nvf_p_c,
499  nvf_r_f,
500  nvf_r_c);
501  high_order = 10.*( (.7-c_courant)*CS_MIN(1., nvf_p_c/.3)
502  + (c_courant-.6)*superbee);
503  }
504  else {
505  high_order = cs_nvd_scheme_scalar(CS_NVD_SUPERBEE,
506  nvf_p_c,
507  nvf_r_f,
508  nvf_r_c);
509  }
510 
511  /* Low order scheme : MUSCL */
512  low_order = cs_nvd_scheme_scalar(CS_NVD_MUSCL,
513  nvf_p_c,
514  nvf_r_f,
515  nvf_r_c);
516 
517  /* Compute the blending factor */
518  if (denom <= (cs_math_epzero*dotp)) {
519  blend = 1.;
520  } else {
521  ratio = dotp/denom;
522  blend = CS_MIN(1., pow(ratio, 2.));
523  }
524 
525  /* Blending */
526  nvf_p_f = blend*high_order + (1.-blend)*low_order;
527  } else { /* STACS scheme */
528  /* High order scheme : SUPERBEE */
529  high_order = cs_nvd_scheme_scalar(CS_NVD_SUPERBEE,
530  nvf_p_c,
531  nvf_r_f,
532  nvf_r_c);
533 
534  /* Low order scheme : STOIC */
535  low_order = cs_nvd_scheme_scalar(CS_NVD_STOIC,
536  nvf_p_c,
537  nvf_r_f,
538  nvf_r_c);
539 
540  /* Compute the blending factor */
541  if (denom <= (cs_math_epzero*dotp)) {
542  blend = 1.;
543  } else {
544  ratio = dotp/denom;
545  blend = CS_MIN(1., pow(ratio, 4.));
546  }
547 
548  /* Blending */
549  nvf_p_f = blend*high_order + (1.-blend)*low_order;
550  }
551 
552  return nvf_p_f;
553 }
554 
555 /*----------------------------------------------------------------------------*/
572 /*----------------------------------------------------------------------------*/
573 
574 inline static void
575 cs_slope_test(const cs_real_t pi,
576  const cs_real_t pj,
577  const cs_real_t distf,
578  const cs_real_t srfan,
579  const cs_real_t i_face_normal[3],
580  const cs_real_t gradi[3],
581  const cs_real_t gradj[3],
582  const cs_real_t grdpai[3],
583  const cs_real_t grdpaj[3],
584  const cs_real_t i_massflux,
585  cs_real_t *testij,
586  cs_real_t *tesqck)
587 {
588  cs_real_t testi, testj;
589  cs_real_t dcc, ddi, ddj;
590 
591  /* Slope test */
592 
593  testi = grdpai[0]*i_face_normal[0]
594  + grdpai[1]*i_face_normal[1]
595  + grdpai[2]*i_face_normal[2];
596  testj = grdpaj[0]*i_face_normal[0]
597  + grdpaj[1]*i_face_normal[1]
598  + grdpaj[2]*i_face_normal[2];
599 
600  *testij = grdpai[0]*grdpaj[0]
601  + grdpai[1]*grdpaj[1]
602  + grdpai[2]*grdpaj[2];
603 
604  if (i_massflux>0.) {
605  dcc = gradi[0]*i_face_normal[0]
606  + gradi[1]*i_face_normal[1]
607  + gradi[2]*i_face_normal[2];
608  ddi = testi;
609  ddj = (pj-pi)/distf *srfan;
610  } else {
611  dcc = gradj[0]*i_face_normal[0]
612  + gradj[1]*i_face_normal[1]
613  + gradj[2]*i_face_normal[2];
614  ddi = (pj-pi)/distf *srfan;
615  ddj = testj;
616  }
617 
618  *tesqck = cs_math_sq(dcc) - cs_math_sq(ddi-ddj);
619 }
620 
621 /*----------------------------------------------------------------------------*/
638 /*----------------------------------------------------------------------------*/
639 
640 inline static void
641 cs_slope_test_vector(const cs_real_t pi[3],
642  const cs_real_t pj[3],
643  const cs_real_t distf,
644  const cs_real_t srfan,
645  const cs_real_t i_face_normal[3],
646  const cs_real_t gradi[3][3],
647  const cs_real_t gradj[3][3],
648  const cs_real_t gradsti[3][3],
649  const cs_real_t gradstj[3][3],
650  const cs_real_t i_massflux,
651  cs_real_t *testij,
652  cs_real_t *tesqck)
653 {
654  cs_real_t testi[3], testj[3];
655  cs_real_t dcc[3], ddi[3], ddj[3];
656  *testij = 0.;
657  *tesqck = 0.;
658 
659  /* Slope test */
660 
661  for (int i = 0; i < 3; i++) {
662  *testij += gradsti[i][0]*gradstj[i][0]
663  + gradsti[i][1]*gradstj[i][1]
664  + gradsti[i][2]*gradstj[i][2];
665 
666  testi[i] = gradsti[i][0]*i_face_normal[0]
667  + gradsti[i][1]*i_face_normal[1]
668  + gradsti[i][2]*i_face_normal[2];
669  testj[i] = gradstj[i][0]*i_face_normal[0]
670  + gradstj[i][1]*i_face_normal[1]
671  + gradstj[i][2]*i_face_normal[2];
672 
673  if (i_massflux > 0.) {
674  dcc[i] = gradi[i][0]*i_face_normal[0]
675  + gradi[i][1]*i_face_normal[1]
676  + gradi[i][2]*i_face_normal[2];
677  ddi[i] = testi[i];
678  ddj[i] = (pj[i]-pi[i])/distf *srfan;
679  } else {
680  dcc[i] = gradj[i][0]*i_face_normal[0]
681  + gradj[i][1]*i_face_normal[1]
682  + gradj[i][2]*i_face_normal[2];
683  ddi[i] = (pj[i]-pi[i])/distf *srfan;
684  ddj[i] = testj[i];
685  }
686  }
687 
688  *tesqck = cs_math_3_square_norm(dcc) - cs_math_3_square_distance(ddi, ddj);
689 }
690 
691 /*----------------------------------------------------------------------------*/
708 /*----------------------------------------------------------------------------*/
709 
710 inline static void
711 cs_slope_test_tensor(const cs_real_t pi[6],
712  const cs_real_t pj[6],
713  const cs_real_t distf,
714  const cs_real_t srfan,
715  const cs_real_t i_face_normal[3],
716  const cs_real_t gradi[6][3],
717  const cs_real_t gradj[6][3],
718  const cs_real_t gradsti[6][3],
719  const cs_real_t gradstj[6][3],
720  const cs_real_t i_massflux,
721  cs_real_t *testij,
722  cs_real_t *tesqck)
723 {
724  cs_real_t testi[6], testj[6];
725  cs_real_t dcc[6], ddi[6], ddj[6];
726 
727  *testij = 0.;
728  *tesqck = 0.;
729 
730  /* Slope test */
731 
732  for (int ij = 0; ij < 6; ij++) {
733  *testij += gradsti[ij][0]*gradstj[ij][0]
734  + gradsti[ij][1]*gradstj[ij][1]
735  + gradsti[ij][2]*gradstj[ij][2];
736  testi[ij] = gradsti[ij][0]*i_face_normal[0]
737  + gradsti[ij][1]*i_face_normal[1]
738  + gradsti[ij][2]*i_face_normal[2];
739  testj[ij] = gradstj[ij][0]*i_face_normal[0]
740  + gradstj[ij][1]*i_face_normal[1]
741  + gradstj[ij][2]*i_face_normal[2];
742 
743  if (i_massflux > 0.) {
744  dcc[ij] = gradi[ij][0]*i_face_normal[0]
745  + gradi[ij][1]*i_face_normal[1]
746  + gradi[ij][2]*i_face_normal[2];
747  ddi[ij] = testi[ij];
748  ddj[ij] = (pj[ij]-pi[ij])/distf *srfan;
749  }
750  else {
751  dcc[ij] = gradj[ij][0]*i_face_normal[0]
752  + gradj[ij][1]*i_face_normal[1]
753  + gradj[ij][2]*i_face_normal[2];
754  ddi[ij] = (pj[ij]-pi[ij])/distf *srfan;
755  ddj[ij] = testj[ij];
756  }
757 
758  *tesqck += cs_math_sq(dcc[ij]) - cs_math_sq(ddi[ij]-ddj[ij]);
759  }
760 }
761 
762 /*----------------------------------------------------------------------------*/
778 /*----------------------------------------------------------------------------*/
779 
780 inline static void
781 cs_i_compute_quantities(const cs_real_t bldfrp,
782  const cs_real_3_t diipf,
783  const cs_real_3_t djjpf,
784  const cs_real_3_t gradi,
785  const cs_real_3_t gradj,
786  const cs_real_t pi,
787  const cs_real_t pj,
788  cs_real_t *recoi,
789  cs_real_t *recoj,
790  cs_real_t *pip,
791  cs_real_t *pjp)
792 {
793  cs_real_t gradpf[3] = {0.5*(gradi[0] + gradj[0]),
794  0.5*(gradi[1] + gradj[1]),
795  0.5*(gradi[2] + gradj[2])};
796 
797  /* reconstruction only if IRCFLP = 1 */
798  *recoi = bldfrp*(cs_math_3_dot_product(gradpf, diipf));
799  *recoj = bldfrp*(cs_math_3_dot_product(gradpf, djjpf));
800  *pip = pi + *recoi;
801  *pjp = pj + *recoj;
802 }
803 
804 /*----------------------------------------------------------------------------*/
820 /*----------------------------------------------------------------------------*/
821 
822 inline static void
823 cs_i_compute_quantities_vector(const cs_real_t bldfrp,
824  const cs_real_3_t diipf,
825  const cs_real_3_t djjpf,
826  const cs_real_33_t gradi,
827  const cs_real_33_t gradj,
828  const cs_real_3_t pi,
829  const cs_real_3_t pj,
830  cs_real_t recoi[3],
831  cs_real_t recoj[3],
832  cs_real_t pip[3],
833  cs_real_t pjp[3])
834 {
835  cs_real_3_t dpvf;
836 
837  /* x-y-z components, p = u, v, w */
838 
839  for (int isou = 0; isou < 3; isou++) {
840 
841  for (int jsou = 0; jsou < 3; jsou++)
842  dpvf[jsou] = 0.5*( gradi[isou][jsou]
843  + gradj[isou][jsou]);
844 
845  /* reconstruction only if IRCFLP = 1 */
846 
847  recoi[isou] = bldfrp*(cs_math_3_dot_product(dpvf, diipf));
848  recoj[isou] = bldfrp*(cs_math_3_dot_product(dpvf, djjpf));
849 
850  pip[isou] = pi[isou] + recoi[isou];
851  pjp[isou] = pj[isou] + recoj[isou];
852 
853  }
854 }
855 
856 /*----------------------------------------------------------------------------*/
872 /*----------------------------------------------------------------------------*/
873 
874 inline static void
875 cs_i_compute_quantities_tensor(const cs_real_t bldfrp,
876  const cs_real_3_t diipf,
877  const cs_real_3_t djjpf,
878  const cs_real_63_t gradi,
879  const cs_real_63_t gradj,
880  const cs_real_6_t pi,
881  const cs_real_6_t pj,
882  cs_real_t recoi[6],
883  cs_real_t recoj[6],
884  cs_real_t pip[6],
885  cs_real_t pjp[6])
886 {
887  cs_real_3_t dpvf;
888 
889  /* x-y-z components, p = u, v, w */
890 
891  for (int isou = 0; isou < 6; isou++) {
892 
893  for (int jsou = 0; jsou < 3; jsou++)
894  dpvf[jsou] = 0.5*( gradi[isou][jsou]
895  + gradj[isou][jsou]);
896 
897  /* reconstruction only if IRCFLP = 1 */
898 
899  recoi[isou] = bldfrp*(cs_math_3_dot_product(dpvf, diipf));
900  recoj[isou] = bldfrp*(cs_math_3_dot_product(dpvf, djjpf));
901 
902  pip[isou] = pi[isou] + recoi[isou];
903  pjp[isou] = pj[isou] + recoj[isou];
904 
905  }
906 }
907 
908 /*----------------------------------------------------------------------------*/
924 /*----------------------------------------------------------------------------*/
925 
926 inline static void
927 cs_i_relax_c_val(const double relaxp,
928  const cs_real_t pia,
929  const cs_real_t pja,
930  const cs_real_t recoi,
931  const cs_real_t recoj,
932  const cs_real_t pi,
933  const cs_real_t pj,
934  cs_real_t *pir,
935  cs_real_t *pjr,
936  cs_real_t *pipr,
937  cs_real_t *pjpr)
938 {
939  *pir = pi/relaxp - (1.-relaxp)/relaxp * pia;
940  *pjr = pj/relaxp - (1.-relaxp)/relaxp * pja;
941 
942  *pipr = *pir + recoi;
943  *pjpr = *pjr + recoj;
944 }
945 
946 /*----------------------------------------------------------------------------*/
962 /*----------------------------------------------------------------------------*/
963 
964 inline static void
965 cs_i_relax_c_val_vector(const double relaxp,
966  const cs_real_3_t pia,
967  const cs_real_3_t pja,
968  const cs_real_3_t recoi,
969  const cs_real_3_t recoj,
970  const cs_real_3_t pi,
971  const cs_real_3_t pj,
972  cs_real_t pir[3],
973  cs_real_t pjr[3],
974  cs_real_t pipr[3],
975  cs_real_t pjpr[3])
976 {
977  for (int isou = 0; isou < 3; isou++) {
978  pir[isou] = pi[isou] /relaxp - (1.-relaxp)/relaxp * pia[isou];
979  pjr[isou] = pj[isou] /relaxp - (1.-relaxp)/relaxp * pja[isou];
980 
981  pipr[isou] = pir[isou] + recoi[isou];
982  pjpr[isou] = pjr[isou] + recoj[isou];
983  }
984 }
985 
986 /*----------------------------------------------------------------------------*/
1002 /*----------------------------------------------------------------------------*/
1003 
1004 inline static void
1005 cs_i_relax_c_val_tensor(const cs_real_t relaxp,
1006  const cs_real_t pia[6],
1007  const cs_real_t pja[6],
1008  const cs_real_t recoi[6],
1009  const cs_real_t recoj[6],
1010  const cs_real_t pi[6],
1011  const cs_real_t pj[6],
1012  cs_real_t pir[6],
1013  cs_real_t pjr[6],
1014  cs_real_t pipr[6],
1015  cs_real_t pjpr[6])
1016 {
1017  for (int isou = 0; isou < 6; isou++) {
1018  pir[isou] = pi[isou] /relaxp - (1.-relaxp)/relaxp * pia[isou];
1019  pjr[isou] = pj[isou] /relaxp - (1.-relaxp)/relaxp * pja[isou];
1020 
1021  pipr[isou] = pir[isou] + recoi[isou];
1022  pjpr[isou] = pjr[isou] + recoj[isou];
1023  }
1024 }
1025 
1026 /*----------------------------------------------------------------------------*/
1033 /*----------------------------------------------------------------------------*/
1034 
1035 inline static void
1036 cs_upwind_f_val(const cs_real_t p,
1037  cs_real_t *pf)
1038 {
1039  *pf = p;
1040 }
1041 
1042 /*----------------------------------------------------------------------------*/
1049 /*----------------------------------------------------------------------------*/
1050 
1051 inline static void
1052 cs_upwind_f_val_vector(const cs_real_3_t p,
1053  cs_real_t pf[3])
1054 {
1055  for (int isou = 0; isou < 3; isou++)
1056  pf[isou] = p[isou];
1057 }
1058 
1059 /*----------------------------------------------------------------------------*/
1066 /*----------------------------------------------------------------------------*/
1067 
1068 inline static void
1069 cs_upwind_f_val_tensor(const cs_real_6_t p,
1070  cs_real_t pf[6])
1071 {
1072  for (int isou = 0; isou < 6; isou++)
1073  pf[isou] = p[isou];
1074 }
1075 
1076 /*----------------------------------------------------------------------------*/
1085 /*----------------------------------------------------------------------------*/
1086 
1087 inline static void
1088 cs_centered_f_val(const double pnd,
1089  const cs_real_t pip,
1090  const cs_real_t pjp,
1091  cs_real_t *pf)
1092 {
1093  *pf = pnd*pip + (1.-pnd)*pjp;
1094 }
1095 
1096 /*----------------------------------------------------------------------------*/
1105 /*----------------------------------------------------------------------------*/
1106 
1107 inline static void
1108 cs_centered_f_val_vector(const double pnd,
1109  const cs_real_3_t pip,
1110  const cs_real_3_t pjp,
1111  cs_real_t pf[3])
1112 {
1113  for (int isou = 0; isou < 3; isou++)
1114  pf[isou] = pnd*pip[isou] + (1.-pnd)*pjp[isou];
1115 }
1116 
1117 /*----------------------------------------------------------------------------*/
1126 /*----------------------------------------------------------------------------*/
1127 
1128 inline static void
1129 cs_centered_f_val_tensor(const double pnd,
1130  const cs_real_6_t pip,
1131  const cs_real_6_t pjp,
1132  cs_real_t pf[6])
1133 {
1134  for (int isou = 0; isou < 6; isou++)
1135  pf[isou] = pnd*pip[isou] + (1.-pnd)*pjp[isou];
1136 }
1137 
1138 /*----------------------------------------------------------------------------*/
1148 /*----------------------------------------------------------------------------*/
1149 
1150 inline static void
1151 cs_solu_f_val(const cs_real_3_t cell_cen,
1152  const cs_real_3_t i_face_cog,
1153  const cs_real_3_t grad,
1154  const cs_real_t p,
1155  cs_real_t *pf)
1156 {
1157  cs_real_t df[3];
1158 
1159  df[0] = i_face_cog[0] - cell_cen[0];
1160  df[1] = i_face_cog[1] - cell_cen[1];
1161  df[2] = i_face_cog[2] - cell_cen[2];
1162 
1163  *pf = p + cs_math_3_dot_product(df, grad);
1164 }
1165 
1166 /*----------------------------------------------------------------------------*/
1176 /*----------------------------------------------------------------------------*/
1177 
1178 inline static void
1179 cs_solu_f_val_vector(const cs_real_3_t cell_cen,
1180  const cs_real_3_t i_face_cog,
1181  const cs_real_33_t grad,
1182  const cs_real_3_t p,
1183  cs_real_t pf[3])
1184 {
1185  cs_real_t df[3];
1186 
1187  for (int jsou = 0; jsou < 3; jsou++)
1188  df[jsou] = i_face_cog[jsou] - cell_cen[jsou];
1189 
1190  for (int isou = 0; isou < 3; isou++) {
1191  pf[isou] = p[isou] + df[0]*grad[isou][0]
1192  + df[1]*grad[isou][1]
1193  + df[2]*grad[isou][2];
1194 
1195  }
1196 }
1197 
1198 /*----------------------------------------------------------------------------*/
1208 /*----------------------------------------------------------------------------*/
1209 
1210 inline static void
1211 cs_solu_f_val_tensor(const cs_real_3_t cell_cen,
1212  const cs_real_3_t i_face_cog,
1213  const cs_real_63_t grad,
1214  const cs_real_6_t p,
1215  cs_real_t pf[6])
1216 {
1217  cs_real_t df[3];
1218 
1219  for (int jsou = 0; jsou < 3; jsou++)
1220  df[jsou] = i_face_cog[jsou] - cell_cen[jsou];
1221 
1222  for (int isou = 0; isou < 6; isou++) {
1223  pf[isou] = p[isou] + df[0]*grad[isou][0]
1224  + df[1]*grad[isou][1]
1225  + df[2]*grad[isou][2];
1226  }
1227 }
1228 
1229 /*----------------------------------------------------------------------------*/
1239 /*----------------------------------------------------------------------------*/
1240 
1241 inline static void
1242 cs_blend_f_val(const double blencp,
1243  const cs_real_t p,
1244  cs_real_t *pf)
1245 {
1246  *pf = blencp * (*pf) + (1. - blencp) * p;
1247 }
1248 
1249 /*----------------------------------------------------------------------------*/
1259 /*----------------------------------------------------------------------------*/
1260 
1261 inline static void
1262 cs_blend_f_val_vector(const double blencp,
1263  const cs_real_3_t p,
1264  cs_real_t pf[3])
1265 {
1266  for (int isou = 0; isou < 3; isou++)
1267  pf[isou] = blencp*(pf[isou])+(1.-blencp)*p[isou];
1268 }
1269 
1270 /*----------------------------------------------------------------------------*/
1280 /*----------------------------------------------------------------------------*/
1281 
1282 inline static void
1283 cs_blend_f_val_tensor(const double blencp,
1284  const cs_real_t p[6],
1285  cs_real_t pf[6])
1286 {
1287  for (int isou = 0; isou < 6; isou++)
1288  pf[isou] = blencp*(pf[isou])+(1.-blencp)*p[isou];
1289 }
1290 
1291 /*----------------------------------------------------------------------------*/
1312 /*----------------------------------------------------------------------------*/
1313 
1314 inline static void
1315 cs_i_conv_flux(const int iconvp,
1316  const cs_real_t thetap,
1317  const int imasac,
1318  const cs_real_t pi,
1319  const cs_real_t pj,
1320  const cs_real_t pifri,
1321  const cs_real_t pifrj,
1322  const cs_real_t pjfri,
1323  const cs_real_t pjfrj,
1324  const cs_real_t i_massflux,
1325  const cs_real_t xcppi,
1326  const cs_real_t xcppj,
1327  cs_real_2_t fluxij)
1328 {
1329  cs_real_t flui, fluj;
1330 
1331  flui = 0.5*(i_massflux + fabs(i_massflux));
1332  fluj = 0.5*(i_massflux - fabs(i_massflux));
1333 
1334  fluxij[0] += iconvp*xcppi*(thetap*(flui*pifri + fluj*pjfri) - imasac*i_massflux*pi);
1335  fluxij[1] += iconvp*xcppj*(thetap*(flui*pifrj + fluj*pjfrj) - imasac*i_massflux*pj);
1336 }
1337 
1338 /*----------------------------------------------------------------------------*/
1356 /*----------------------------------------------------------------------------*/
1357 
1358 inline static void
1359 cs_i_conv_flux_vector(const int iconvp,
1360  const cs_real_t thetap,
1361  const int imasac,
1362  const cs_real_t pi[3],
1363  const cs_real_t pj[3],
1364  const cs_real_t pifri[3],
1365  const cs_real_t pifrj[3],
1366  const cs_real_t pjfri[3],
1367  const cs_real_t pjfrj[3],
1368  const cs_real_t i_massflux,
1369  cs_real_t fluxi[3],
1370  cs_real_t fluxj[3])
1371 {
1372  cs_real_t flui, fluj;
1373 
1374  flui = 0.5*(i_massflux + fabs(i_massflux));
1375  fluj = 0.5*(i_massflux - fabs(i_massflux));
1376 
1377  for (int isou = 0; isou < 3; isou++) {
1378 
1379  fluxi[isou] += iconvp*( thetap*(flui*pifri[isou] + fluj*pjfri[isou])
1380  - imasac*i_massflux*pi[isou]);
1381  fluxj[isou] += iconvp*( thetap*(flui*pifrj[isou] + fluj*pjfrj[isou])
1382  - imasac*i_massflux*pj[isou]);
1383  }
1384 }
1385 
1386 /*----------------------------------------------------------------------------*/
1404 /*----------------------------------------------------------------------------*/
1405 
1406 inline static void
1407 cs_i_conv_flux_tensor(const int iconvp,
1408  const cs_real_t thetap,
1409  const int imasac,
1410  const cs_real_t pi[6],
1411  const cs_real_t pj[6],
1412  const cs_real_t pifri[6],
1413  const cs_real_t pifrj[6],
1414  const cs_real_t pjfri[6],
1415  const cs_real_t pjfrj[6],
1416  const cs_real_t i_massflux,
1417  cs_real_t fluxi[6],
1418  cs_real_t fluxj[6])
1419 {
1420  cs_real_t flui, fluj;
1421 
1422  flui = 0.5*(i_massflux + fabs(i_massflux));
1423  fluj = 0.5*(i_massflux - fabs(i_massflux));
1424 
1425  for (int isou = 0; isou < 6; isou++) {
1426  fluxi[isou] += iconvp*( thetap*(flui*pifri[isou] + fluj*pjfri[isou])
1427  - imasac*i_massflux*pi[isou]);
1428  fluxj[isou] += iconvp*( thetap*(flui*pifrj[isou] + fluj*pjfrj[isou])
1429  - imasac*i_massflux*pj[isou]);
1430  }
1431 }
1432 
1433 /*----------------------------------------------------------------------------*/
1446 /*----------------------------------------------------------------------------*/
1447 
1448 inline static void
1449 cs_i_diff_flux(const int idiffp,
1450  const cs_real_t thetap,
1451  const cs_real_t pip,
1452  const cs_real_t pjp,
1453  const cs_real_t pipr,
1454  const cs_real_t pjpr,
1455  const cs_real_t i_visc,
1456  cs_real_t fluxij[2])
1457 {
1458  fluxij[0] += idiffp*thetap*i_visc*(pipr -pjp);
1459  fluxij[1] += idiffp*thetap*i_visc*(pip -pjpr);
1460 }
1461 
1462 /*----------------------------------------------------------------------------*/
1476 /*----------------------------------------------------------------------------*/
1477 
1478 inline static void
1479 cs_i_diff_flux_vector(const int idiffp,
1480  const cs_real_t thetap,
1481  const cs_real_t pip[3],
1482  const cs_real_t pjp[3],
1483  const cs_real_t pipr[3],
1484  const cs_real_t pjpr[3],
1485  const cs_real_t i_visc,
1486  cs_real_t fluxi[3],
1487  cs_real_t fluxj[3])
1488 {
1489  for (int isou = 0; isou < 3; isou++) {
1490  fluxi[isou] += idiffp*thetap*i_visc*(pipr[isou] -pjp[isou]);
1491  fluxj[isou] += idiffp*thetap*i_visc*(pip[isou] -pjpr[isou]);
1492  }
1493 }
1494 
1495 /*----------------------------------------------------------------------------*/
1509 /*----------------------------------------------------------------------------*/
1510 
1511 inline static void
1512 cs_i_diff_flux_tensor(const int idiffp,
1513  const cs_real_t thetap,
1514  const cs_real_t pip[6],
1515  const cs_real_t pjp[6],
1516  const cs_real_t pipr[6],
1517  const cs_real_t pjpr[6],
1518  const cs_real_t i_visc,
1519  cs_real_t fluxi[6],
1520  cs_real_t fluxj[6])
1521 {
1522  for (int isou = 0; isou < 6; isou++) {
1523  fluxi[isou] += idiffp*thetap*i_visc*(pipr[isou] -pjp[isou]);
1524  fluxj[isou] += idiffp*thetap*i_visc*(pip[isou] -pjpr[isou]);
1525  }
1526 }
1527 
1528 /*----------------------------------------------------------------------------*/
1552 /*----------------------------------------------------------------------------*/
1553 
1554 inline static void
1555 cs_i_cd_steady_upwind(const cs_real_t bldfrp,
1556  const cs_real_t relaxp,
1557  const cs_real_t diipf[3],
1558  const cs_real_t djjpf[3],
1559  const cs_real_t gradi[3],
1560  const cs_real_t gradj[3],
1561  const cs_real_t pi,
1562  const cs_real_t pj,
1563  const cs_real_t pia,
1564  const cs_real_t pja,
1565  cs_real_t *pifri,
1566  cs_real_t *pifrj,
1567  cs_real_t *pjfri,
1568  cs_real_t *pjfrj,
1569  cs_real_t *pip,
1570  cs_real_t *pjp,
1571  cs_real_t *pipr,
1572  cs_real_t *pjpr)
1573 {
1574  cs_real_t pir, pjr;
1575  cs_real_t recoi, recoj;
1576 
1577  cs_i_compute_quantities(bldfrp,
1578  diipf,
1579  djjpf,
1580  gradi,
1581  gradj,
1582  pi,
1583  pj,
1584  &recoi,
1585  &recoj,
1586  pip,
1587  pjp);
1588 
1589  cs_i_relax_c_val(relaxp,
1590  pia,
1591  pja,
1592  recoi,
1593  recoj,
1594  pi,
1595  pj,
1596  &pir,
1597  &pjr,
1598  pipr,
1599  pjpr);
1600 
1601  cs_upwind_f_val(pi,
1602  pifrj);
1603  cs_upwind_f_val(pir,
1604  pifri);
1605  cs_upwind_f_val(pj,
1606  pjfri);
1607  cs_upwind_f_val(pjr,
1608  pjfrj);
1609 }
1610 
1611 /*----------------------------------------------------------------------------*/
1635 /*----------------------------------------------------------------------------*/
1636 
1637 inline static void
1638 cs_i_cd_steady_upwind_vector(const cs_real_t bldfrp,
1639  const cs_real_t relaxp,
1640  const cs_real_t diipf[3],
1641  const cs_real_t djjpf[3],
1642  const cs_real_t gradi[3][3],
1643  const cs_real_t gradj[3][3],
1644  const cs_real_t pi[3],
1645  const cs_real_t pj[3],
1646  const cs_real_t pia[3],
1647  const cs_real_t pja[3],
1648  cs_real_t pifri[3],
1649  cs_real_t pifrj[3],
1650  cs_real_t pjfri[3],
1651  cs_real_t pjfrj[3],
1652  cs_real_t pip[3],
1653  cs_real_t pjp[3],
1654  cs_real_t pipr[3],
1655  cs_real_t pjpr[3])
1656 {
1657  cs_real_t pir[3], pjr[3];
1658  cs_real_t recoi[3], recoj[3];
1659 
1660  cs_i_compute_quantities_vector(bldfrp,
1661  diipf,
1662  djjpf,
1663  gradi,
1664  gradj,
1665  pi,
1666  pj,
1667  recoi,
1668  recoj,
1669  pip,
1670  pjp);
1671 
1672  cs_i_relax_c_val_vector(relaxp,
1673  pia,
1674  pja,
1675  recoi,
1676  recoj,
1677  pi,
1678  pj,
1679  pir,
1680  pjr,
1681  pipr,
1682  pjpr);
1683 
1684  cs_upwind_f_val_vector(pi,
1685  pifrj);
1686  cs_upwind_f_val_vector(pir,
1687  pifri);
1688  cs_upwind_f_val_vector(pj,
1689  pjfri);
1690  cs_upwind_f_val_vector(pjr,
1691  pjfrj);
1692 }
1693 
1694 /*----------------------------------------------------------------------------*/
1718 /*----------------------------------------------------------------------------*/
1719 
1720 inline static void
1721 cs_i_cd_steady_upwind_tensor(const cs_real_t bldfrp,
1722  const cs_real_t relaxp,
1723  const cs_real_t diipf[3],
1724  const cs_real_t djjpf[3],
1725  const cs_real_t gradi[6][3],
1726  const cs_real_t gradj[6][3],
1727  const cs_real_t pi[6],
1728  const cs_real_t pj[6],
1729  const cs_real_t pia[6],
1730  const cs_real_t pja[6],
1731  cs_real_t pifri[6],
1732  cs_real_t pifrj[6],
1733  cs_real_t pjfri[6],
1734  cs_real_t pjfrj[6],
1735  cs_real_t pip[6],
1736  cs_real_t pjp[6],
1737  cs_real_t pipr[6],
1738  cs_real_t pjpr[6])
1739 {
1740  cs_real_t pir[6], pjr[6];
1741  cs_real_t recoi[6], recoj[6];
1742 
1743  cs_i_compute_quantities_tensor(bldfrp,
1744  diipf,
1745  djjpf,
1746  gradi,
1747  gradj,
1748  pi,
1749  pj,
1750  recoi,
1751  recoj,
1752  pip,
1753  pjp);
1754 
1755  cs_i_relax_c_val_tensor(relaxp,
1756  pia,
1757  pja,
1758  recoi,
1759  recoj,
1760  pi,
1761  pj,
1762  pir,
1763  pjr,
1764  pipr,
1765  pjpr);
1766 
1767  cs_upwind_f_val_tensor(pi,
1768  pifrj);
1769  cs_upwind_f_val_tensor(pir,
1770  pifri);
1771  cs_upwind_f_val_tensor(pj,
1772  pjfri);
1773  cs_upwind_f_val_tensor(pjr,
1774  pjfrj);
1775 }
1776 
1777 /*----------------------------------------------------------------------------*/
1794 /*----------------------------------------------------------------------------*/
1795 
1796 inline static void
1797 cs_i_cd_unsteady_upwind(const cs_real_t bldfrp,
1798  const cs_real_t diipf[3],
1799  const cs_real_t djjpf[3],
1800  const cs_real_t gradi[3],
1801  const cs_real_t gradj[3],
1802  const cs_real_t pi,
1803  const cs_real_t pj,
1804  cs_real_t *pif,
1805  cs_real_t *pjf,
1806  cs_real_t *pip,
1807  cs_real_t *pjp)
1808 {
1809  cs_real_t recoi, recoj;
1810 
1811  cs_i_compute_quantities(bldfrp,
1812  diipf,
1813  djjpf,
1814  gradi,
1815  gradj,
1816  pi,
1817  pj,
1818  &recoi,
1819  &recoj,
1820  pip,
1821  pjp);
1822 
1823  cs_upwind_f_val(pi, pif);
1824  cs_upwind_f_val(pj, pjf);
1825 }
1826 
1827 /*----------------------------------------------------------------------------*/
1844 /*----------------------------------------------------------------------------*/
1845 
1846 inline static void
1847 cs_i_cd_unsteady_upwind_vector(const cs_real_t bldfrp,
1848  const cs_real_t diipf[3],
1849  const cs_real_t djjpf[3],
1850  const cs_real_t gradi[3][3],
1851  const cs_real_t gradj[3][3],
1852  const cs_real_t pi[3],
1853  const cs_real_t pj[3],
1854  cs_real_t pif[3],
1855  cs_real_t pjf[3],
1856  cs_real_t pip[3],
1857  cs_real_t pjp[3])
1858 {
1859  cs_real_t recoi[3], recoj[3];
1860 
1861  cs_i_compute_quantities_vector(bldfrp,
1862  diipf,
1863  djjpf,
1864  gradi,
1865  gradj,
1866  pi,
1867  pj,
1868  recoi,
1869  recoj,
1870  pip,
1871  pjp);
1872 
1873  cs_upwind_f_val_vector(pi, pif);
1874  cs_upwind_f_val_vector(pj, pjf);
1875 }
1876 
1877 /*----------------------------------------------------------------------------*/
1894 /*----------------------------------------------------------------------------*/
1895 
1896 inline static void
1897 cs_i_cd_unsteady_upwind_tensor(const cs_real_t bldfrp,
1898  const cs_real_t diipf[3],
1899  const cs_real_t djjpf[3],
1900  const cs_real_t gradi[6][3],
1901  const cs_real_t gradj[6][3],
1902  const cs_real_t pi[6],
1903  const cs_real_t pj[6],
1904  cs_real_t pif[6],
1905  cs_real_t pjf[6],
1906  cs_real_t pip[6],
1907  cs_real_t pjp[6])
1908 {
1909  cs_real_t recoi[6], recoj[6];
1910 
1911  cs_i_compute_quantities_tensor(bldfrp,
1912  diipf,
1913  djjpf,
1914  gradi,
1915  gradj,
1916  pi,
1917  pj,
1918  recoi,
1919  recoj,
1920  pip,
1921  pjp);
1922 
1923  cs_upwind_f_val_tensor(pi, pif);
1924  cs_upwind_f_val_tensor(pj, pjf);
1925 
1926 }
1927 
1928 /*----------------------------------------------------------------------------*/
1961 /*----------------------------------------------------------------------------*/
1962 
1963 inline static void
1964 cs_i_cd_steady(const cs_real_t bldfrp,
1965  const int ischcp,
1966  const double relaxp,
1967  const double blencp,
1968  const cs_real_t weight,
1969  const cs_real_t cell_ceni[3],
1970  const cs_real_t cell_cenj[3],
1971  const cs_real_t i_face_cog[3],
1972  const cs_real_t diipf[3],
1973  const cs_real_t djjpf[3],
1974  const cs_real_t gradi[3],
1975  const cs_real_t gradj[3],
1976  const cs_real_t gradupi[3],
1977  const cs_real_t gradupj[3],
1978  const cs_real_t pi,
1979  const cs_real_t pj,
1980  const cs_real_t pia,
1981  const cs_real_t pja,
1982  cs_real_t *pifri,
1983  cs_real_t *pifrj,
1984  cs_real_t *pjfri,
1985  cs_real_t *pjfrj,
1986  cs_real_t *pip,
1987  cs_real_t *pjp,
1988  cs_real_t *pipr,
1989  cs_real_t *pjpr)
1990 {
1991  cs_real_t pir, pjr;
1992  cs_real_t recoi, recoj;
1993 
1994  cs_i_compute_quantities(bldfrp,
1995  diipf,
1996  djjpf,
1997  gradi,
1998  gradj,
1999  pi,
2000  pj,
2001  &recoi,
2002  &recoj,
2003  pip,
2004  pjp);
2005 
2006  cs_i_relax_c_val(relaxp,
2007  pia,
2008  pja,
2009  recoi,
2010  recoj,
2011  pi,
2012  pj,
2013  &pir,
2014  &pjr,
2015  pipr,
2016  pjpr);
2017 
2018  if (ischcp == 1) {
2019 
2020  /* Centered
2021  --------*/
2022 
2023  cs_centered_f_val(weight,
2024  *pip,
2025  *pjpr,
2026  pifrj);
2027  cs_centered_f_val(weight,
2028  *pipr,
2029  *pjp,
2030  pifri);
2031  cs_centered_f_val(weight,
2032  *pipr,
2033  *pjp,
2034  pjfri);
2035  cs_centered_f_val(weight,
2036  *pip,
2037  *pjpr,
2038  pjfrj);
2039 
2040  } else if (ischcp == 0) {
2041 
2042  /* Original SOLU
2043  --------------*/
2044 
2045  cs_solu_f_val(cell_ceni,
2046  i_face_cog,
2047  gradi,
2048  pi,
2049  pifrj);
2050  cs_solu_f_val(cell_ceni,
2051  i_face_cog,
2052  gradi,
2053  pir,
2054  pifri);
2055  cs_solu_f_val(cell_cenj,
2056  i_face_cog,
2057  gradj,
2058  pj,
2059  pjfri);
2060  cs_solu_f_val(cell_cenj,
2061  i_face_cog,
2062  gradj,
2063  pjr,
2064  pjfrj);
2065 
2066  } else {
2067 
2068  /* SOLU
2069  ----*/
2070 
2071  cs_solu_f_val(cell_ceni,
2072  i_face_cog,
2073  gradupi,
2074  pi,
2075  pifrj);
2076  cs_solu_f_val(cell_ceni,
2077  i_face_cog,
2078  gradupi,
2079  pir,
2080  pifri);
2081  cs_solu_f_val(cell_cenj,
2082  i_face_cog,
2083  gradupj,
2084  pj,
2085  pjfri);
2086  cs_solu_f_val(cell_cenj,
2087  i_face_cog,
2088  gradupj,
2089  pjr,
2090  pjfrj);
2091 
2092  }
2093 
2094  /* Blending
2095  --------*/
2096 
2097  cs_blend_f_val(blencp,
2098  pi,
2099  pifrj);
2100  cs_blend_f_val(blencp,
2101  pir,
2102  pifri);
2103  cs_blend_f_val(blencp,
2104  pj,
2105  pjfri);
2106  cs_blend_f_val(blencp,
2107  pjr,
2108  pjfrj);
2109 }
2110 
2111 /*----------------------------------------------------------------------------*/
2142 /*----------------------------------------------------------------------------*/
2143 
2144 inline static void
2145 cs_i_cd_steady_vector(const cs_real_t bldfrp,
2146  const int ischcp,
2147  const double relaxp,
2148  const double blencp,
2149  const cs_real_t weight,
2150  const cs_real_3_t cell_ceni,
2151  const cs_real_3_t cell_cenj,
2152  const cs_real_3_t i_face_cog,
2153  const cs_real_3_t diipf,
2154  const cs_real_3_t djjpf,
2155  const cs_real_33_t gradi,
2156  const cs_real_33_t gradj,
2157  const cs_real_3_t pi,
2158  const cs_real_3_t pj,
2159  const cs_real_3_t pia,
2160  const cs_real_3_t pja,
2161  cs_real_t pifri[3],
2162  cs_real_t pifrj[3],
2163  cs_real_t pjfri[3],
2164  cs_real_t pjfrj[3],
2165  cs_real_t pip[3],
2166  cs_real_t pjp[3],
2167  cs_real_t pipr[3],
2168  cs_real_t pjpr[3])
2169 {
2170  cs_real_t pir[3], pjr[3];
2171  cs_real_t recoi[3], recoj[3];
2172 
2173  cs_i_compute_quantities_vector(bldfrp,
2174  diipf,
2175  djjpf,
2176  gradi,
2177  gradj,
2178  pi,
2179  pj,
2180  recoi,
2181  recoj,
2182  pip,
2183  pjp);
2184 
2185  cs_i_relax_c_val_vector(relaxp,
2186  pia,
2187  pja,
2188  recoi,
2189  recoj,
2190  pi,
2191  pj,
2192  pir,
2193  pjr,
2194  pipr,
2195  pjpr);
2196 
2197  if (ischcp == 1) {
2198 
2199  /* Centered
2200  --------*/
2201 
2202  cs_centered_f_val_vector(weight,
2203  pip,
2204  pjpr,
2205  pifrj);
2206  cs_centered_f_val_vector(weight,
2207  pipr,
2208  pjp,
2209  pifri);
2210  cs_centered_f_val_vector(weight,
2211  pipr,
2212  pjp,
2213  pjfri);
2214  cs_centered_f_val_vector(weight,
2215  pip,
2216  pjpr,
2217  pjfrj);
2218 
2219  } else {
2220 
2221  /* Second order
2222  ------------*/
2223 
2224  cs_solu_f_val_vector(cell_ceni,
2225  i_face_cog,
2226  gradi,
2227  pi,
2228  pifrj);
2229  cs_solu_f_val_vector(cell_ceni,
2230  i_face_cog,
2231  gradi,
2232  pir,
2233  pifri);
2234  cs_solu_f_val_vector(cell_cenj,
2235  i_face_cog,
2236  gradj,
2237  pj,
2238  pjfri);
2239  cs_solu_f_val_vector(cell_cenj,
2240  i_face_cog,
2241  gradj,
2242  pjr,
2243  pjfrj);
2244 
2245  }
2246 
2247  /* Blending
2248  --------*/
2249  cs_blend_f_val_vector(blencp,
2250  pi,
2251  pifrj);
2252  cs_blend_f_val_vector(blencp,
2253  pir,
2254  pifri);
2255  cs_blend_f_val_vector(blencp,
2256  pj,
2257  pjfri);
2258  cs_blend_f_val_vector(blencp,
2259  pjr,
2260  pjfrj);
2261 
2262 }
2263 
2264 /*----------------------------------------------------------------------------*/
2295 /*----------------------------------------------------------------------------*/
2296 
2297 inline static void
2298 cs_i_cd_steady_tensor(const cs_real_t bldfrp,
2299  const int ischcp,
2300  const double relaxp,
2301  const double blencp,
2302  const cs_real_t weight,
2303  const cs_real_3_t cell_ceni,
2304  const cs_real_3_t cell_cenj,
2305  const cs_real_3_t i_face_cog,
2306  const cs_real_3_t diipf,
2307  const cs_real_3_t djjpf,
2308  const cs_real_63_t gradi,
2309  const cs_real_63_t gradj,
2310  const cs_real_6_t pi,
2311  const cs_real_6_t pj,
2312  const cs_real_6_t pia,
2313  const cs_real_6_t pja,
2314  cs_real_t pifri[6],
2315  cs_real_t pifrj[6],
2316  cs_real_t pjfri[6],
2317  cs_real_t pjfrj[6],
2318  cs_real_t pip[6],
2319  cs_real_t pjp[6],
2320  cs_real_t pipr[6],
2321  cs_real_t pjpr[6])
2322 
2323 {
2324  cs_real_t pir[6], pjr[6];
2325  cs_real_t recoi[6], recoj[6];
2326 
2327  cs_i_compute_quantities_tensor(bldfrp,
2328  diipf,
2329  djjpf,
2330  gradi,
2331  gradj,
2332  pi,
2333  pj,
2334  recoi,
2335  recoj,
2336  pip,
2337  pjp);
2338 
2339  cs_i_relax_c_val_tensor(relaxp,
2340  pia,
2341  pja,
2342  recoi,
2343  recoj,
2344  pi,
2345  pj,
2346  pir,
2347  pjr,
2348  pipr,
2349  pjpr);
2350 
2351  if (ischcp == 1) {
2352 
2353  /* Centered
2354  --------*/
2355 
2356  cs_centered_f_val_tensor(weight,
2357  pip,
2358  pjpr,
2359  pifrj);
2360  cs_centered_f_val_tensor(weight,
2361  pipr,
2362  pjp,
2363  pifri);
2364  cs_centered_f_val_tensor(weight,
2365  pipr,
2366  pjp,
2367  pjfri);
2368  cs_centered_f_val_tensor(weight,
2369  pip,
2370  pjpr,
2371  pjfrj);
2372 
2373  } else {
2374 
2375  /* Second order
2376  ------------*/
2377 
2378  cs_solu_f_val_tensor(cell_ceni,
2379  i_face_cog,
2380  gradi,
2381  pi,
2382  pifrj);
2383  cs_solu_f_val_tensor(cell_ceni,
2384  i_face_cog,
2385  gradi,
2386  pir,
2387  pifri);
2388  cs_solu_f_val_tensor(cell_cenj,
2389  i_face_cog,
2390  gradj,
2391  pj,
2392  pjfri);
2393  cs_solu_f_val_tensor(cell_cenj,
2394  i_face_cog,
2395  gradj,
2396  pjr,
2397  pjfrj);
2398 
2399  }
2400 
2401  /* Blending
2402  --------*/
2403 
2404  cs_blend_f_val_tensor(blencp,
2405  pi,
2406  pifrj);
2407  cs_blend_f_val_tensor(blencp,
2408  pir,
2409  pifri);
2410  cs_blend_f_val_tensor(blencp,
2411  pj,
2412  pjfri);
2413  cs_blend_f_val_tensor(blencp,
2414  pjr,
2415  pjfrj);
2416 }
2417 
2418 /*----------------------------------------------------------------------------*/
2448 /*----------------------------------------------------------------------------*/
2449 
2450 inline static void
2451 cs_i_cd_unsteady(const cs_real_t bldfrp,
2452  const int ischcp,
2453  const double blencp,
2454  const cs_real_t weight,
2455  const cs_real_3_t cell_ceni,
2456  const cs_real_3_t cell_cenj,
2457  const cs_real_3_t i_face_cog,
2458  const cs_real_t hybrid_blend_i,
2459  const cs_real_t hybrid_blend_j,
2460  const cs_real_3_t diipf,
2461  const cs_real_3_t djjpf,
2462  const cs_real_3_t gradi,
2463  const cs_real_3_t gradj,
2464  const cs_real_3_t gradupi,
2465  const cs_real_3_t gradupj,
2466  const cs_real_t pi,
2467  const cs_real_t pj,
2468  cs_real_t *pif,
2469  cs_real_t *pjf,
2470  cs_real_t *pip,
2471  cs_real_t *pjp)
2472 {
2473  cs_real_t recoi, recoj;
2474 
2475  cs_i_compute_quantities(bldfrp,
2476  diipf,
2477  djjpf,
2478  gradi,
2479  gradj,
2480  pi,
2481  pj,
2482  &recoi,
2483  &recoj,
2484  pip,
2485  pjp);
2486 
2487 
2488  if (ischcp == 1) {
2489 
2490  /* Centered
2491  --------*/
2492 
2493  cs_centered_f_val(weight,
2494  *pip,
2495  *pjp,
2496  pif);
2497  cs_centered_f_val(weight,
2498  *pip,
2499  *pjp,
2500  pjf);
2501 
2502  } else if (ischcp == 0) {
2503 
2504  /* Legacy SOLU
2505  -----------*/
2506 
2507  cs_solu_f_val(cell_ceni,
2508  i_face_cog,
2509  gradi,
2510  pi,
2511  pif);
2512  cs_solu_f_val(cell_cenj,
2513  i_face_cog,
2514  gradj,
2515  pj,
2516  pjf);
2517 
2518  } else if (ischcp == 3) {
2519 
2520  /* Centered
2521  --------*/
2522 
2523  cs_centered_f_val(weight,
2524  *pip,
2525  *pjp,
2526  pif);
2527  cs_centered_f_val(weight,
2528  *pip,
2529  *pjp,
2530  pjf);
2531 
2532  /* Legacy SOLU
2533  -----------*/
2534  cs_real_t pif_up, pjf_up;
2535  cs_real_t hybrid_blend_interp;
2536 
2537  cs_solu_f_val(cell_ceni,
2538  i_face_cog,
2539  gradi,
2540  pi,
2541  &pif_up);
2542  cs_solu_f_val(cell_cenj,
2543  i_face_cog,
2544  gradj,
2545  pj,
2546  &pjf_up);
2547 
2548  hybrid_blend_interp = fmin(hybrid_blend_i,hybrid_blend_j);
2549  *pif = hybrid_blend_interp*(*pif) + (1. - hybrid_blend_interp)*pif_up;
2550  *pjf = hybrid_blend_interp*(*pjf) + (1. - hybrid_blend_interp)*pjf_up;
2551 
2552  } else {
2553 
2554  /* SOLU
2555  ----*/
2556 
2557  cs_solu_f_val(cell_ceni,
2558  i_face_cog,
2559  gradupi,
2560  pi,
2561  pif);
2562  cs_solu_f_val(cell_cenj,
2563  i_face_cog,
2564  gradupj,
2565  pj,
2566  pjf);
2567 
2568  }
2569 
2570 
2571  /* Blending
2572  --------*/
2573 
2574  cs_blend_f_val(blencp,
2575  pi,
2576  pif);
2577  cs_blend_f_val(blencp,
2578  pj,
2579  pjf);
2580 }
2581 
2582 /*----------------------------------------------------------------------------*/
2608 /*----------------------------------------------------------------------------*/
2609 
2610 inline static void
2611 cs_i_cd_unsteady_vector(const cs_real_t bldfrp,
2612  const int ischcp,
2613  const double blencp,
2614  const cs_real_t weight,
2615  const cs_real_3_t cell_ceni,
2616  const cs_real_3_t cell_cenj,
2617  const cs_real_3_t i_face_cog,
2618  const cs_real_t hybrid_blend_i,
2619  const cs_real_t hybrid_blend_j,
2620  const cs_real_3_t diipf,
2621  const cs_real_3_t djjpf,
2622  const cs_real_33_t gradi,
2623  const cs_real_33_t gradj,
2624  const cs_real_3_t pi,
2625  const cs_real_3_t pj,
2626  cs_real_t pif[3],
2627  cs_real_t pjf[3],
2628  cs_real_t pip[3],
2629  cs_real_t pjp[3])
2630 
2631 {
2632  cs_real_t recoi[3], recoj[3];
2633 
2634  cs_i_compute_quantities_vector(bldfrp,
2635  diipf,
2636  djjpf,
2637  gradi,
2638  gradj,
2639  pi,
2640  pj,
2641  recoi,
2642  recoj,
2643  pip,
2644  pjp);
2645 
2646  if (ischcp == 1) {
2647 
2648  /* Centered
2649  --------*/
2650 
2651  cs_centered_f_val_vector(weight,
2652  pip,
2653  pjp,
2654  pif);
2655  cs_centered_f_val_vector(weight,
2656  pip,
2657  pjp,
2658  pjf);
2659  } else if (ischcp == 3) {
2660 
2661  /* Centered
2662  --------*/
2663 
2664  cs_centered_f_val_vector(weight,
2665  pip,
2666  pjp,
2667  pif);
2668  cs_centered_f_val_vector(weight,
2669  pip,
2670  pjp,
2671  pjf);
2672 
2673  /* SOLU
2674  -----*/
2675  cs_real_t pif_up[3], pjf_up[3];
2676  cs_real_t hybrid_blend_interp;
2677 
2678  cs_solu_f_val_vector(cell_ceni,
2679  i_face_cog,
2680  gradi,
2681  pi,
2682  pif_up);
2683  cs_solu_f_val_vector(cell_cenj,
2684  i_face_cog,
2685  gradj,
2686  pj,
2687  pjf_up);
2688 
2689  hybrid_blend_interp = fmin(hybrid_blend_i,hybrid_blend_j);
2690  for (int isou = 0; isou < 3; isou++) {
2691  pif[isou] = hybrid_blend_interp *pif[isou]
2692  + (1. - hybrid_blend_interp)*pif_up[isou];
2693  pjf[isou] = hybrid_blend_interp *pjf[isou]
2694  + (1. - hybrid_blend_interp)*pjf_up[isou];
2695  }
2696  } else {
2697 
2698  /* Second order
2699  ------------*/
2700 
2701  cs_solu_f_val_vector(cell_ceni,
2702  i_face_cog,
2703  gradi,
2704  pi,
2705  pif);
2706  cs_solu_f_val_vector(cell_cenj,
2707  i_face_cog,
2708  gradj,
2709  pj,
2710  pjf);
2711 
2712  }
2713 
2714  /* Blending
2715  --------*/
2716 
2717  cs_blend_f_val_vector(blencp,
2718  pi,
2719  pif);
2720  cs_blend_f_val_vector(blencp,
2721  pj,
2722  pjf);
2723 
2724 }
2725 
2726 /*----------------------------------------------------------------------------*/
2750 /*----------------------------------------------------------------------------*/
2751 
2752 inline static void
2753 cs_i_cd_unsteady_tensor(const cs_real_t bldfrp,
2754  const int ischcp,
2755  const double blencp,
2756  const cs_real_t weight,
2757  const cs_real_3_t cell_ceni,
2758  const cs_real_3_t cell_cenj,
2759  const cs_real_3_t i_face_cog,
2760  const cs_real_3_t diipf,
2761  const cs_real_3_t djjpf,
2762  const cs_real_63_t gradi,
2763  const cs_real_63_t gradj,
2764  const cs_real_6_t pi,
2765  const cs_real_6_t pj,
2766  cs_real_t pif[6],
2767  cs_real_t pjf[6],
2768  cs_real_t pip[6],
2769  cs_real_t pjp[6])
2770 
2771 {
2772  cs_real_t recoi[6], recoj[6];
2773 
2774  cs_i_compute_quantities_tensor(bldfrp,
2775  diipf,
2776  djjpf,
2777  gradi,
2778  gradj,
2779  pi,
2780  pj,
2781  recoi,
2782  recoj,
2783  pip,
2784  pjp);
2785 
2786  if (ischcp == 1) {
2787 
2788  /* Centered
2789  --------*/
2790 
2791  cs_centered_f_val_tensor(weight,
2792  pip,
2793  pjp,
2794  pif);
2795  cs_centered_f_val_tensor(weight,
2796  pip,
2797  pjp,
2798  pjf);
2799 
2800  } else {
2801 
2802  /* Second order
2803  ------------*/
2804 
2805  cs_solu_f_val_tensor(cell_ceni,
2806  i_face_cog,
2807  gradi,
2808  pi,
2809  pif);
2810  cs_solu_f_val_tensor(cell_cenj,
2811  i_face_cog,
2812  gradj,
2813  pj,
2814  pjf);
2815 
2816  }
2817 
2818  /* Blending
2819  --------*/
2820 
2821  cs_blend_f_val_tensor(blencp,
2822  pi,
2823  pif);
2824  cs_blend_f_val_tensor(blencp,
2825  pj,
2826  pjf);
2827 
2828 }
2829 
2830 /*----------------------------------------------------------------------------*/
2874 /*----------------------------------------------------------------------------*/
2875 
2876 inline static void
2877 cs_i_cd_steady_slope_test(bool *upwind_switch,
2878  const int iconvp,
2879  const cs_real_t bldfrp,
2880  const int ischcp,
2881  const double relaxp,
2882  const double blencp,
2883  const double blend_st,
2884  const cs_real_t weight,
2885  const cs_real_t i_dist,
2886  const cs_real_t i_face_surf,
2887  const cs_real_3_t cell_ceni,
2888  const cs_real_3_t cell_cenj,
2889  const cs_real_3_t i_face_normal,
2890  const cs_real_3_t i_face_cog,
2891  const cs_real_3_t diipf,
2892  const cs_real_3_t djjpf,
2893  const cs_real_t i_massflux,
2894  const cs_real_3_t gradi,
2895  const cs_real_3_t gradj,
2896  const cs_real_3_t gradupi,
2897  const cs_real_3_t gradupj,
2898  const cs_real_3_t gradsti,
2899  const cs_real_3_t gradstj,
2900  const cs_real_t pi,
2901  const cs_real_t pj,
2902  const cs_real_t pia,
2903  const cs_real_t pja,
2904  cs_real_t *pifri,
2905  cs_real_t *pifrj,
2906  cs_real_t *pjfri,
2907  cs_real_t *pjfrj,
2908  cs_real_t *pip,
2909  cs_real_t *pjp,
2910  cs_real_t *pipr,
2911  cs_real_t *pjpr)
2912 {
2913  cs_real_t pir, pjr;
2914  cs_real_t recoi, recoj;
2915  cs_real_t testij, tesqck;
2916 
2917  *upwind_switch = false;
2918 
2919  cs_i_compute_quantities(bldfrp,
2920  diipf,
2921  djjpf,
2922  gradi,
2923  gradj,
2924  pi,
2925  pj,
2926  &recoi,
2927  &recoj,
2928  pip,
2929  pjp);
2930 
2931  cs_i_relax_c_val(relaxp,
2932  pia,
2933  pja,
2934  recoi,
2935  recoj,
2936  pi,
2937  pj,
2938  &pir,
2939  &pjr,
2940  pipr,
2941  pjpr);
2942 
2943  /* Convection slope test is needed only when iconv >0 */
2944  if (iconvp > 0) {
2945  cs_slope_test(pi,
2946  pj,
2947  i_dist,
2948  i_face_surf,
2949  i_face_normal,
2950  gradi,
2951  gradj,
2952  gradsti,
2953  gradstj,
2954  i_massflux,
2955  &testij,
2956  &tesqck);
2957 
2958  if (ischcp==1) {
2959 
2960  /* Centered
2961  --------*/
2962 
2963  cs_centered_f_val(weight,
2964  *pip,
2965  *pjpr,
2966  pifrj);
2967  cs_centered_f_val(weight,
2968  *pipr,
2969  *pjp,
2970  pifri);
2971  cs_centered_f_val(weight,
2972  *pipr,
2973  *pjp,
2974  pjfri);
2975  cs_centered_f_val(weight,
2976  *pip,
2977  *pjpr,
2978  pjfrj);
2979 
2980  } else if (ischcp == 0) {
2981 
2982  /* Second order
2983  ------------*/
2984 
2985  cs_solu_f_val(cell_ceni,
2986  i_face_cog,
2987  gradi,
2988  pi,
2989  pifrj);
2990  cs_solu_f_val(cell_ceni,
2991  i_face_cog,
2992  gradi,
2993  pir,
2994  pifri);
2995  cs_solu_f_val(cell_cenj,
2996  i_face_cog,
2997  gradj,
2998  pj,
2999  pjfri);
3000  cs_solu_f_val(cell_cenj,
3001  i_face_cog,
3002  gradj,
3003  pjr,
3004  pjfrj);
3005 
3006  } else {
3007 
3008  /* SOLU
3009  -----*/
3010 
3011  cs_solu_f_val(cell_ceni,
3012  i_face_cog,
3013  gradupi,
3014  pi,
3015  pifrj);
3016  cs_solu_f_val(cell_ceni,
3017  i_face_cog,
3018  gradupi,
3019  pir,
3020  pifri);
3021  cs_solu_f_val(cell_cenj,
3022  i_face_cog,
3023  gradupj,
3024  pj,
3025  pjfri);
3026  cs_solu_f_val(cell_cenj,
3027  i_face_cog,
3028  gradupj,
3029  pjr,
3030  pjfrj);
3031  }
3032 
3033 
3034  /* Slope test: Pourcentage of upwind
3035  ----------------------------------*/
3036 
3037  if (tesqck <= 0. || testij <= 0.) {
3038 
3039  cs_blend_f_val(blend_st,
3040  pi,
3041  pifrj);
3042  cs_blend_f_val(blend_st,
3043  pir,
3044  pifri);
3045  cs_blend_f_val(blend_st,
3046  pj,
3047  pjfri);
3048  cs_blend_f_val(blend_st,
3049  pjr,
3050  pjfrj);
3051 
3052  *upwind_switch = true;
3053 
3054  }
3055 
3056 
3057  /* Blending
3058  --------*/
3059 
3060  cs_blend_f_val(blencp,
3061  pi,
3062  pifrj);
3063  cs_blend_f_val(blencp,
3064  pir,
3065  pifri);
3066  cs_blend_f_val(blencp,
3067  pj,
3068  pjfri);
3069  cs_blend_f_val(blencp,
3070  pjr,
3071  pjfrj);
3072 
3073  /* If iconv=0 p*fr* are useless */
3074  } else {
3075  cs_upwind_f_val(pi,
3076  pifrj);
3077  cs_upwind_f_val(pir,
3078  pifri);
3079  cs_upwind_f_val(pj,
3080  pjfri);
3081  cs_upwind_f_val(pjr,
3082  pjfrj);
3083  }
3084 
3085 }
3086 
3087 /*----------------------------------------------------------------------------*/
3129 /*----------------------------------------------------------------------------*/
3130 
3131 inline static void
3132 cs_i_cd_steady_slope_test_vector(bool *upwind_switch,
3133  const int iconvp,
3134  const cs_real_t bldfrp,
3135  const int ischcp,
3136  const double relaxp,
3137  const double blencp,
3138  const double blend_st,
3139  const cs_real_t weight,
3140  const cs_real_t i_dist,
3141  const cs_real_t i_face_surf,
3142  const cs_real_3_t cell_ceni,
3143  const cs_real_3_t cell_cenj,
3144  const cs_real_3_t i_face_normal,
3145  const cs_real_3_t i_face_cog,
3146  const cs_real_3_t diipf,
3147  const cs_real_3_t djjpf,
3148  const cs_real_t i_massflux,
3149  const cs_real_33_t gradi,
3150  const cs_real_33_t gradj,
3151  const cs_real_33_t grdpai,
3152  const cs_real_33_t grdpaj,
3153  const cs_real_3_t pi,
3154  const cs_real_3_t pj,
3155  const cs_real_3_t pia,
3156  const cs_real_3_t pja,
3157  cs_real_t pifri[3],
3158  cs_real_t pifrj[3],
3159  cs_real_t pjfri[3],
3160  cs_real_t pjfrj[3],
3161  cs_real_t pip[3],
3162  cs_real_t pjp[3],
3163  cs_real_t pipr[3],
3164  cs_real_t pjpr[3])
3165 {
3166  cs_real_t pir[3], pjr[3];
3167  cs_real_t recoi[3], recoj[3];
3168  cs_real_t testij, tesqck;
3169  int isou;
3170 
3171  cs_i_compute_quantities_vector(bldfrp,
3172  diipf,
3173  djjpf,
3174  gradi,
3175  gradj,
3176  pi,
3177  pj,
3178  recoi,
3179  recoj,
3180  pip,
3181  pjp);
3182 
3183  cs_i_relax_c_val_vector(relaxp,
3184  pia,
3185  pja,
3186  recoi,
3187  recoj,
3188  pi,
3189  pj,
3190  pir,
3191  pjr,
3192  pipr,
3193  pjpr);
3194 
3195  /* Convection slope test is needed only when iconv >0 */
3196  if (iconvp > 0) {
3197  cs_slope_test_vector(pi,
3198  pj,
3199  i_dist,
3200  i_face_surf,
3201  i_face_normal,
3202  gradi,
3203  gradj,
3204  grdpai,
3205  grdpaj,
3206  i_massflux,
3207  &testij,
3208  &tesqck);
3209 
3210  for (isou = 0; isou < 3; isou++) {
3211  if (ischcp==1) {
3212 
3213  /* Centered
3214  --------*/
3215 
3216  cs_centered_f_val(weight,
3217  pip[isou],
3218  pjpr[isou],
3219  &pifrj[isou]);
3220  cs_centered_f_val(weight,
3221  pipr[isou],
3222  pjp[isou],
3223  &pifri[isou]);
3224  cs_centered_f_val(weight,
3225  pipr[isou],
3226  pjp[isou],
3227  &pjfri[isou]);
3228  cs_centered_f_val(weight,
3229  pip[isou],
3230  pjpr[isou],
3231  &pjfrj[isou]);
3232 
3233  } else {
3234 
3235  /* Second order
3236  ------------*/
3237 
3238  cs_solu_f_val(cell_ceni,
3239  i_face_cog,
3240  gradi[isou],
3241  pi[isou],
3242  &pifrj[isou]);
3243  cs_solu_f_val(cell_ceni,
3244  i_face_cog,
3245  gradi[isou],
3246  pir[isou],
3247  &pifri[isou]);
3248  cs_solu_f_val(cell_cenj,
3249  i_face_cog,
3250  gradj[isou],
3251  pj[isou],
3252  &pjfri[isou]);
3253  cs_solu_f_val(cell_cenj,
3254  i_face_cog,
3255  gradj[isou],
3256  pjr[isou],
3257  &pjfrj[isou]);
3258 
3259  }
3260 
3261  }
3262 
3263  /* Slope test: Pourcentage of upwind
3264  ----------------------------------*/
3265 
3266  if (tesqck <= 0. || testij <= 0.) {
3267  cs_blend_f_val_vector(blend_st,
3268  pi,
3269  pifrj);
3270  cs_blend_f_val_vector(blend_st,
3271  pir,
3272  pifri);
3273  cs_blend_f_val_vector(blend_st,
3274  pj,
3275  pjfri);
3276  cs_blend_f_val_vector(blend_st,
3277  pjr,
3278  pjfrj);
3279 
3280  *upwind_switch = true;
3281  }
3282 
3283 
3284  /* Blending
3285  --------*/
3286  cs_blend_f_val_vector(blencp,
3287  pi,
3288  pifrj);
3289  cs_blend_f_val_vector(blencp,
3290  pir,
3291  pifri);
3292  cs_blend_f_val_vector(blencp,
3293  pj,
3294  pjfri);
3295  cs_blend_f_val_vector(blencp,
3296  pjr,
3297  pjfrj);
3298 
3299  /* If iconv=0 p*fr* are useless */
3300  } else {
3301  for (isou = 0; isou < 3; isou++) {
3302  cs_upwind_f_val(pi[isou],
3303  &pifrj[isou]);
3304  cs_upwind_f_val(pir[isou],
3305  &pifri[isou]);
3306  cs_upwind_f_val(pj[isou],
3307  &pjfri[isou]);
3308  cs_upwind_f_val(pjr[isou],
3309  &pjfrj[isou]);
3310  }
3311  }
3312 
3313 }
3314 
3315 /*----------------------------------------------------------------------------*/
3357 /*----------------------------------------------------------------------------*/
3358 
3359 inline static void
3360 cs_i_cd_steady_slope_test_tensor(bool *upwind_switch,
3361  const int iconvp,
3362  const cs_real_t bldfrp,
3363  const int ischcp,
3364  const double relaxp,
3365  const double blencp,
3366  const double blend_st,
3367  const cs_real_t weight,
3368  const cs_real_t i_dist,
3369  const cs_real_t i_face_surf,
3370  const cs_real_3_t cell_ceni,
3371  const cs_real_3_t cell_cenj,
3372  const cs_real_3_t i_face_normal,
3373  const cs_real_3_t i_face_cog,
3374  const cs_real_3_t diipf,
3375  const cs_real_3_t djjpf,
3376  const cs_real_t i_massflux,
3377  const cs_real_63_t gradi,
3378  const cs_real_63_t gradj,
3379  const cs_real_63_t grdpai,
3380  const cs_real_63_t grdpaj,
3381  const cs_real_6_t pi,
3382  const cs_real_6_t pj,
3383  const cs_real_6_t pia,
3384  const cs_real_6_t pja,
3385  cs_real_t pifri[6],
3386  cs_real_t pifrj[6],
3387  cs_real_t pjfri[6],
3388  cs_real_t pjfrj[6],
3389  cs_real_t pip[6],
3390  cs_real_t pjp[6],
3391  cs_real_t pipr[6],
3392  cs_real_t pjpr[6])
3393 {
3394  cs_real_t pir[6], pjr[6];
3395  cs_real_t recoi[6], recoj[6];
3396  cs_real_t testij, tesqck;
3397  int isou;
3398 
3399  cs_i_compute_quantities_tensor(bldfrp,
3400  diipf,
3401  djjpf,
3402  gradi,
3403  gradj,
3404  pi,
3405  pj,
3406  recoi,
3407  recoj,
3408  pip,
3409  pjp);
3410 
3411  cs_i_relax_c_val_tensor(relaxp,
3412  pia,
3413  pja,
3414  recoi,
3415  recoj,
3416  pi,
3417  pj,
3418  pir,
3419  pjr,
3420  pipr,
3421  pjpr);
3422 
3423  /* Convection slope test is needed only when iconv >0 */
3424  if (iconvp > 0) {
3425  cs_slope_test_tensor(pi,
3426  pj,
3427  i_dist,
3428  i_face_surf,
3429  i_face_normal,
3430  gradi,
3431  gradj,
3432  grdpai,
3433  grdpaj,
3434  i_massflux,
3435  &testij,
3436  &tesqck);
3437 
3438  for (isou = 0; isou < 6; isou++) {
3439  if (ischcp==1) {
3440 
3441  /* Centered
3442  --------*/
3443 
3444  cs_centered_f_val(weight,
3445  pip[isou],
3446  pjpr[isou],
3447  &pifrj[isou]);
3448  cs_centered_f_val(weight,
3449  pipr[isou],
3450  pjp[isou],
3451  &pifri[isou]);
3452  cs_centered_f_val(weight,
3453  pipr[isou],
3454  pjp[isou],
3455  &pjfri[isou]);
3456  cs_centered_f_val(weight,
3457  pip[isou],
3458  pjpr[isou],
3459  &pjfrj[isou]);
3460 
3461  } else {
3462 
3463  /* Second order
3464  ------------*/
3465 
3466  cs_solu_f_val(cell_ceni,
3467  i_face_cog,
3468  gradi[isou],
3469  pi[isou],
3470  &pifrj[isou]);
3471  cs_solu_f_val(cell_ceni,
3472  i_face_cog,
3473  gradi[isou],
3474  pir[isou],
3475  &pifri[isou]);
3476  cs_solu_f_val(cell_cenj,
3477  i_face_cog,
3478  gradj[isou],
3479  pj[isou],
3480  &pjfri[isou]);
3481  cs_solu_f_val(cell_cenj,
3482  i_face_cog,
3483  gradj[isou],
3484  pjr[isou],
3485  &pjfrj[isou]);
3486 
3487  }
3488 
3489  }
3490 
3491  /* Slope test: Pourcentage of upwind
3492  ----------------------------------*/
3493 
3494  if (tesqck <= 0. || testij <= 0.) {
3495 
3496  cs_blend_f_val_tensor(blend_st,
3497  pi,
3498  pifrj);
3499  cs_blend_f_val_tensor(blend_st,
3500  pir,
3501  pifri);
3502  cs_blend_f_val_tensor(blend_st,
3503  pj,
3504  pjfri);
3505  cs_blend_f_val_tensor(blend_st,
3506  pjr,
3507  pjfrj);
3508 
3509  *upwind_switch = true;
3510 
3511  }
3512 
3513 
3514  /* Blending
3515  --------*/
3516 
3517  cs_blend_f_val_tensor(blencp,
3518  pi,
3519  pifrj);
3520  cs_blend_f_val_tensor(blencp,
3521  pir,
3522  pifri);
3523  cs_blend_f_val_tensor(blencp,
3524  pj,
3525  pjfri);
3526  cs_blend_f_val_tensor(blencp,
3527  pjr,
3528  pjfrj);
3529 
3530  /* If iconv=0 p*fr* are useless */
3531  } else {
3532  for (isou = 0; isou < 6; isou++) {
3533  cs_upwind_f_val(pi[isou],
3534  &pifrj[isou]);
3535  cs_upwind_f_val(pir[isou],
3536  &pifri[isou]);
3537  cs_upwind_f_val(pj[isou],
3538  &pjfri[isou]);
3539  cs_upwind_f_val(pjr[isou],
3540  &pjfrj[isou]);
3541  }
3542  }
3543 
3544 }
3545 
3546 /*----------------------------------------------------------------------------*/
3583 /*----------------------------------------------------------------------------*/
3584 
3585 inline static void
3586 cs_i_cd_unsteady_slope_test(bool *upwind_switch,
3587  const int iconvp,
3588  const cs_real_t bldfrp,
3589  const int ischcp,
3590  const double blencp,
3591  const double blend_st,
3592  const cs_real_t weight,
3593  const cs_real_t i_dist,
3594  const cs_real_t i_face_surf,
3595  const cs_real_3_t cell_ceni,
3596  const cs_real_3_t cell_cenj,
3597  const cs_real_3_t i_face_normal,
3598  const cs_real_3_t i_face_cog,
3599  const cs_real_3_t diipf,
3600  const cs_real_3_t djjpf,
3601  const cs_real_t i_massflux,
3602  const cs_real_3_t gradi,
3603  const cs_real_3_t gradj,
3604  const cs_real_3_t gradupi,
3605  const cs_real_3_t gradupj,
3606  const cs_real_3_t gradsti,
3607  const cs_real_3_t gradstj,
3608  const cs_real_t pi,
3609  const cs_real_t pj,
3610  cs_real_t *pif,
3611  cs_real_t *pjf,
3612  cs_real_t *pip,
3613  cs_real_t *pjp)
3614 {
3615  CS_UNUSED(blend_st);
3616 
3617  cs_real_t recoi, recoj;
3618  cs_real_t testij, tesqck;
3619 
3620  *upwind_switch = false;
3621 
3622  cs_i_compute_quantities(bldfrp,
3623  diipf,
3624  djjpf,
3625  gradi,
3626  gradj,
3627  pi,
3628  pj,
3629  &recoi,
3630  &recoj,
3631  pip,
3632  pjp);
3633 
3634  /* Convection slope test is needed only when iconv >0 */
3635  if (iconvp > 0) {
3636  cs_slope_test(pi,
3637  pj,
3638  i_dist,
3639  i_face_surf,
3640  i_face_normal,
3641  gradi,
3642  gradj,
3643  gradsti,
3644  gradstj,
3645  i_massflux,
3646  &testij,
3647  &tesqck);
3648 
3649  if (ischcp==1) {
3650 
3651  /* Centered
3652  --------*/
3653 
3654  cs_centered_f_val(weight,
3655  *pip,
3656  *pjp,
3657  pif);
3658  cs_centered_f_val(weight,
3659  *pip,
3660  *pjp,
3661  pjf);
3662 
3663  } else if (ischcp == 0) {
3664 
3665  /* Original SOLU
3666  --------------*/
3667 
3668  cs_solu_f_val(cell_ceni,
3669  i_face_cog,
3670  gradi,
3671  pi,
3672  pif);
3673  cs_solu_f_val(cell_cenj,
3674  i_face_cog,
3675  gradj,
3676  pj,
3677  pjf);
3678 
3679  } else {
3680 
3681  /* SOLU
3682  -----*/
3683 
3684  cs_solu_f_val(cell_ceni,
3685  i_face_cog,
3686  gradupi,
3687  pi,
3688  pif);
3689  cs_solu_f_val(cell_cenj,
3690  i_face_cog,
3691  gradupj,
3692  pj,
3693  pjf);
3694 
3695  }
3696 
3697  /* Slope test: Pourcentage of upwind
3698  ----------------------------------*/
3699 
3700  if (tesqck<=0. || testij<=0.) {
3701 
3702  cs_blend_f_val(blend_st,
3703  pi,
3704  pif);
3705  cs_blend_f_val(blend_st,
3706  pj,
3707  pjf);
3708 
3709  *upwind_switch = true;
3710 
3711  }
3712 
3713  /* Blending
3714  --------*/
3715 
3716  cs_blend_f_val(blencp,
3717  pi,
3718  pif);
3719  cs_blend_f_val(blencp,
3720  pj,
3721  pjf);
3722 
3723  /* If iconv=0 p*f are useless */
3724  } else {
3725  cs_upwind_f_val(pi,
3726  pif);
3727  cs_upwind_f_val(pj,
3728  pjf);
3729  }
3730 
3731 }
3732 
3733 /*----------------------------------------------------------------------------*/
3744 /*----------------------------------------------------------------------------*/
3745 
3746 inline static void
3747 cs_central_downwind_cells(const cs_lnum_t ii,
3748  const cs_lnum_t jj,
3749  const cs_real_t i_massflux,
3750  cs_lnum_t *ic,
3751  cs_lnum_t *id)
3752 {
3753  if (i_massflux >= 0.) {
3754  *ic = ii;
3755  *id = jj;
3756  } else {
3757  *ic = jj;
3758  *id = ii;
3759  }
3760 }
3761 
3762 /*----------------------------------------------------------------------------*/
3784 /*----------------------------------------------------------------------------*/
3785 
3786 inline static void
3787 cs_i_cd_unsteady_nvd(const cs_nvd_type_t limiter,
3788  const double beta,
3789  const cs_real_3_t cell_cen_c,
3790  const cs_real_3_t cell_cen_d,
3791  const cs_real_3_t i_face_normal,
3792  const cs_real_3_t i_face_cog,
3793  const cs_real_3_t gradv_c,
3794  const cs_real_t p_c,
3795  const cs_real_t p_d,
3796  const cs_real_t local_max_c,
3797  const cs_real_t local_min_c,
3798  const cs_real_t courant_c,
3799  cs_real_t *pif,
3800  cs_real_t *pjf)
3801 {
3802  /* distance between face center and central cell center */
3803  cs_real_t dist_fc;
3804  cs_real_3_t nfc;
3805  cs_math_3_length_unitv(cell_cen_c, i_face_cog, &dist_fc, nfc);
3806 
3807  /* unit vector and distance between central and downwind cells centers */
3808  cs_real_t dist_dc;
3809  cs_real_3_t ndc;
3810  cs_math_3_length_unitv(cell_cen_c, cell_cen_d, &dist_dc, ndc);
3811 
3812  /* Place the upwind point on the line that joins
3813  the two cells on the upwind side and the same
3814  distance as that between the two cells */
3815  const cs_real_t dist_cu = dist_dc;
3816  const cs_real_t dist_du = dist_dc + dist_cu;
3817 
3818  /* Compute the property on the upwind assuming a parabolic
3819  variation of the property between the two cells */
3820  const cs_real_t gradc = cs_math_3_dot_product(gradv_c, ndc);
3821 
3822  const cs_real_t grad2c = ((p_d - p_c)/dist_dc - gradc)/dist_dc;
3823 
3824  cs_real_t p_u = p_c + (grad2c*dist_cu - gradc)*dist_cu;
3825  p_u = CS_MAX(CS_MIN(p_u, local_max_c), local_min_c);
3826 
3827  /* Compute the normalised distances */
3828  const cs_real_t nvf_r_f = (dist_fc+dist_cu)/dist_du;
3829  const cs_real_t nvf_r_c = dist_cu/dist_du;
3830 
3831  /* Check for the bounds of the NVD diagram and compute the face
3832  property according to the selected NVD scheme */
3833  const cs_real_t _small = cs_math_epzero
3834  * (CS_ABS(p_u)+CS_ABS(p_c)+CS_ABS(p_d));
3835 
3836  if (CS_ABS(p_d-p_u) <= _small) {
3837  *pif = p_c;
3838  *pjf = p_c;
3839  } else {
3840  const cs_real_t nvf_p_c = (p_c - p_u)/(p_d - p_u);
3841 
3842  if (nvf_p_c <= 0. || nvf_p_c >= 1.) {
3843  *pif = p_c;
3844  *pjf = p_c;
3845  } else {
3846  cs_real_t nvf_p_f;
3847 
3848  /* Highly compressive NVD scheme for VOF */
3849  if (limiter >= CS_NVD_VOF_HRIC) {
3850  nvf_p_f = cs_nvd_vof_scheme_scalar(limiter,
3851  i_face_normal,
3852  nvf_p_c,
3853  nvf_r_f,
3854  nvf_r_c,
3855  gradv_c,
3856  courant_c);
3857  } else { /* Regular NVD scheme */
3858  nvf_p_f = cs_nvd_scheme_scalar(limiter,
3859  nvf_p_c,
3860  nvf_r_f,
3861  nvf_r_c);
3862  }
3863 
3864  *pif = p_u + nvf_p_f*(p_d - p_u);
3865  *pif = CS_MAX(CS_MIN(*pif, local_max_c), local_min_c);
3866 
3867  cs_blend_f_val(beta,
3868  p_c,
3869  pif);
3870 
3871  *pjf = *pif;
3872  }
3873  }
3874 }
3875 
3876 /*----------------------------------------------------------------------------*/
3911 /*----------------------------------------------------------------------------*/
3912 
3913 inline static void
3914 cs_i_cd_unsteady_slope_test_vector(bool *upwind_switch,
3915  const int iconvp,
3916  const cs_real_t bldfrp,
3917  const int ischcp,
3918  const double blencp,
3919  const double blend_st,
3920  const cs_real_t weight,
3921  const cs_real_t i_dist,
3922  const cs_real_t i_face_surf,
3923  const cs_real_3_t cell_ceni,
3924  const cs_real_3_t cell_cenj,
3925  const cs_real_3_t i_face_normal,
3926  const cs_real_3_t i_face_cog,
3927  const cs_real_3_t diipf,
3928  const cs_real_3_t djjpf,
3929  const cs_real_t i_massflux,
3930  const cs_real_33_t gradi,
3931  const cs_real_33_t gradj,
3932  const cs_real_33_t grdpai,
3933  const cs_real_33_t grdpaj,
3934  const cs_real_3_t pi,
3935  const cs_real_3_t pj,
3936  cs_real_t pif[3],
3937  cs_real_t pjf[3],
3938  cs_real_t pip[3],
3939  cs_real_t pjp[3])
3940 {
3941  cs_real_t recoi[3], recoj[3];
3942  cs_real_t testij, tesqck;
3943 
3944  cs_i_compute_quantities_vector(bldfrp,
3945  diipf,
3946  djjpf,
3947  gradi,
3948  gradj,
3949  pi,
3950  pj,
3951  recoi,
3952  recoj,
3953  pip,
3954  pjp);
3955 
3956  /* Convection slope test is needed only when iconv >0 */
3957  if (iconvp > 0) {
3958  cs_slope_test_vector(pi,
3959  pj,
3960  i_dist,
3961  i_face_surf,
3962  i_face_normal,
3963  gradi,
3964  gradj,
3965  grdpai,
3966  grdpaj,
3967  i_massflux,
3968  &testij,
3969  &tesqck);
3970 
3971  for (int isou = 0; isou < 3; isou++) {
3972  if (ischcp == 1) {
3973 
3974  /* Centered
3975  --------*/
3976 
3977  cs_centered_f_val(weight,
3978  pip[isou],
3979  pjp[isou],
3980  &pif[isou]);
3981  cs_centered_f_val(weight,
3982  pip[isou],
3983  pjp[isou],
3984  &pjf[isou]);
3985 
3986  } else {
3987 
3988  /* Second order
3989  ------------*/
3990 
3991  cs_solu_f_val(cell_ceni,
3992  i_face_cog,
3993  gradi[isou],
3994  pi[isou],
3995  &pif[isou]);
3996  cs_solu_f_val(cell_cenj,
3997  i_face_cog,
3998  gradj[isou],
3999  pj[isou],
4000  &pjf[isou]);
4001 
4002  }
4003 
4004  }
4005 
4006  /* Slope test: Pourcentage of upwind
4007  ----------------------------------*/
4008 
4009  if (tesqck <= 0. || testij <= 0.) {
4010 
4011  cs_blend_f_val_vector(blend_st,
4012  pi,
4013  pif);
4014  cs_blend_f_val_vector(blend_st,
4015  pj,
4016  pjf);
4017 
4018  *upwind_switch = true;
4019 
4020  }
4021 
4022 
4023  /* Blending
4024  --------*/
4025  cs_blend_f_val_vector(blencp,
4026  pi,
4027  pif);
4028  cs_blend_f_val_vector(blencp,
4029  pj,
4030  pjf);
4031 
4032  /* If iconv=0 p*f are useless */
4033  } else {
4034 
4035  for (int isou = 0; isou < 3; isou++) {
4036  cs_upwind_f_val(pi[isou],
4037  &pif[isou]);
4038  cs_upwind_f_val(pj[isou],
4039  &pjf[isou]);
4040 
4041  }
4042  }
4043 
4044 }
4045 
4046 /*----------------------------------------------------------------------------*/
4081 /*----------------------------------------------------------------------------*/
4082 
4083 inline static void
4084 cs_i_cd_unsteady_slope_test_tensor(bool *upwind_switch,
4085  const int iconvp,
4086  const cs_real_t bldfrp,
4087  const int ischcp,
4088  const double blencp,
4089  const double blend_st,
4090  const cs_real_t weight,
4091  const cs_real_t i_dist,
4092  const cs_real_t i_face_surf,
4093  const cs_real_3_t cell_ceni,
4094  const cs_real_3_t cell_cenj,
4095  const cs_real_3_t i_face_normal,
4096  const cs_real_3_t i_face_cog,
4097  const cs_real_3_t diipf,
4098  const cs_real_3_t djjpf,
4099  const cs_real_t i_massflux,
4100  const cs_real_63_t gradi,
4101  const cs_real_63_t gradj,
4102  const cs_real_63_t grdpai,
4103  const cs_real_63_t grdpaj,
4104  const cs_real_6_t pi,
4105  const cs_real_6_t pj,
4106  cs_real_t pif[6],
4107  cs_real_t pjf[6],
4108  cs_real_t pip[6],
4109  cs_real_t pjp[6])
4110 {
4111  cs_real_t recoi[6], recoj[6];
4112  cs_real_t testij, tesqck;
4113  int isou;
4114 
4115  cs_i_compute_quantities_tensor(bldfrp,
4116  diipf,
4117  djjpf,
4118  gradi,
4119  gradj,
4120  pi,
4121  pj,
4122  recoi,
4123  recoj,
4124  pip,
4125  pjp);
4126 
4127  /* Convection slope test is needed only when iconv >0 */
4128  if (iconvp > 0) {
4129  cs_slope_test_tensor(pi,
4130  pj,
4131  i_dist,
4132  i_face_surf,
4133  i_face_normal,
4134  gradi,
4135  gradj,
4136  grdpai,
4137  grdpaj,
4138  i_massflux,
4139  &testij,
4140  &tesqck);
4141 
4142  for (isou = 0; isou < 6; isou++) {
4143 
4144  if (ischcp==1) {
4145 
4146  /* Centered
4147  --------*/
4148 
4149  cs_centered_f_val(weight,
4150  pip[isou],
4151  pjp[isou],
4152  &pif[isou]);
4153  cs_centered_f_val(weight,
4154  pip[isou],
4155  pjp[isou],
4156  &pjf[isou]);
4157 
4158  } else {
4159 
4160  /* Second order
4161  ------------*/
4162 
4163  cs_solu_f_val(cell_ceni,
4164  i_face_cog,
4165  gradi[isou],
4166  pi[isou],
4167  &pif[isou]);
4168  cs_solu_f_val(cell_cenj,
4169  i_face_cog,
4170  gradj[isou],
4171  pj[isou],
4172  &pjf[isou]);
4173  }
4174 
4175  }
4176 
4177  /* Slope test activated: poucentage of upwind */
4178  if (tesqck <= 0. || testij <= 0.) {
4179 
4180  /* Upwind
4181  --------*/
4182 
4183  cs_blend_f_val_tensor(blend_st,
4184  pi,
4185  pif);
4186  cs_blend_f_val_tensor(blend_st,
4187  pj,
4188  pjf);
4189 
4190  *upwind_switch = true;
4191  }
4192 
4193 
4194  /* Blending
4195  --------*/
4196 
4197  cs_blend_f_val_tensor(blencp,
4198  pi,
4199  pif);
4200  cs_blend_f_val_tensor(blencp,
4201  pj,
4202  pjf);
4203 
4204  /* If iconv=0 p*fr* are useless */
4205  } else {
4206 
4207  for (isou = 0; isou < 6; isou++) {
4208  cs_upwind_f_val(pi[isou],
4209  &pif[isou]);
4210  cs_upwind_f_val(pj[isou],
4211  &pjf[isou]);
4212  }
4213  }
4214 }
4215 
4216 /*----------------------------------------------------------------------------*/
4225 /*----------------------------------------------------------------------------*/
4226 
4227 inline static void
4228 cs_b_compute_quantities(const cs_real_3_t diipb,
4229  const cs_real_3_t gradi,
4230  const cs_real_t bldfrp,
4231  cs_real_t *recoi)
4232 {
4233  *recoi = bldfrp * ( gradi[0]*diipb[0]
4234  + gradi[1]*diipb[1]
4235  + gradi[2]*diipb[2]);
4236 }
4237 
4238 /*----------------------------------------------------------------------------*/
4247 /*----------------------------------------------------------------------------*/
4248 
4249 inline static void
4250 cs_b_compute_quantities_vector(const cs_real_3_t diipb,
4251  const cs_real_33_t gradi,
4252  const cs_real_t bldfrp,
4253  cs_real_t recoi[3])
4254 {
4255  for (int isou = 0; isou < 3; isou++) {
4256  recoi[isou] = bldfrp * ( gradi[isou][0]*diipb[0]
4257  + gradi[isou][1]*diipb[1]
4258  + gradi[isou][2]*diipb[2]);
4259  }
4260 }
4261 
4262 /*----------------------------------------------------------------------------*/
4271 /*----------------------------------------------------------------------------*/
4272 
4273 inline static void
4274 cs_b_compute_quantities_tensor(const cs_real_3_t diipb,
4275  const cs_real_63_t gradi,
4276  const cs_real_t bldfrp,
4277  cs_real_t recoi[6])
4278 {
4279  for (int isou = 0; isou < 6; isou++) {
4280  recoi[isou] = bldfrp * ( gradi[isou][0]*diipb[0]
4281  + gradi[isou][1]*diipb[1]
4282  + gradi[isou][2]*diipb[2]);
4283  }
4284 }
4285 
4286 /*----------------------------------------------------------------------------*/
4297 /*----------------------------------------------------------------------------*/
4298 
4299 inline static void
4300 cs_b_relax_c_val(const double relaxp,
4301  const cs_real_t pi,
4302  const cs_real_t pia,
4303  const cs_real_t recoi,
4304  cs_real_t *pir,
4305  cs_real_t *pipr)
4306 {
4307  *pir = pi/relaxp - (1.-relaxp)/relaxp*pia;
4308  *pipr = *pir + recoi;
4309 }
4310 
4311 /*----------------------------------------------------------------------------*/
4322 /*----------------------------------------------------------------------------*/
4323 
4324 inline static void
4325 cs_b_relax_c_val_vector(const double relaxp,
4326  const cs_real_3_t pi,
4327  const cs_real_3_t pia,
4328  const cs_real_3_t recoi,
4329  cs_real_t pir[3],
4330  cs_real_t pipr[3])
4331 {
4332  for (int isou = 0; isou < 3; isou++) {
4333  pir[isou] = pi[isou]/relaxp - (1.-relaxp)/relaxp*pia[isou];
4334  pipr[isou] = pir[isou] + recoi[isou];
4335  }
4336 }
4337 
4338 /*----------------------------------------------------------------------------*/
4349 /*----------------------------------------------------------------------------*/
4350 
4351 inline static void
4352 cs_b_relax_c_val_tensor(const double relaxp,
4353  const cs_real_6_t pi,
4354  const cs_real_6_t pia,
4355  const cs_real_6_t recoi,
4356  cs_real_t pir[6],
4357  cs_real_t pipr[6])
4358 {
4359  for (int isou = 0; isou < 6; isou++) {
4360  pir[isou] = pi[isou]/relaxp - (1.-relaxp)/relaxp*pia[isou];
4361  pipr[isou] = pir[isou] + recoi[isou];
4362  }
4363 }
4364 
4365 /*----------------------------------------------------------------------------*/
4389 /*----------------------------------------------------------------------------*/
4390 
4391 inline static void
4392 cs_b_imposed_conv_flux(int iconvp,
4393  cs_real_t thetap,
4394  int imasac,
4395  int inc,
4396  int bc_type,
4397  int icvfli,
4398  cs_real_t pi,
4399  cs_real_t pir,
4400  cs_real_t pipr,
4401  cs_real_t coefap,
4402  cs_real_t coefbp,
4403  cs_real_t coface,
4404  cs_real_t cofbce,
4405  cs_real_t b_massflux,
4406  cs_real_t xcpp,
4407  cs_real_t *flux)
4408 {
4409  cs_real_t flui, fluj, pfac;
4410 
4411  /* Computed convective flux */
4412 
4413  if (icvfli == 0) {
4414 
4415  /* Remove decentering for coupled faces */
4416  if (bc_type == CS_COUPLED_FD) {
4417  flui = 0.0;
4418  fluj = b_massflux;
4419  } else {
4420  flui = 0.5*(b_massflux +fabs(b_massflux));
4421  fluj = 0.5*(b_massflux -fabs(b_massflux));
4422  }
4423 
4424  pfac = inc*coefap + coefbp*pipr;
4425  *flux += iconvp*xcpp*(thetap*(flui*pir + fluj*pfac) -imasac*( b_massflux*pi));
4426 
4427  /* Imposed convective flux */
4428 
4429  } else {
4430 
4431  pfac = inc*coface + cofbce*pipr;
4432  *flux += iconvp*xcpp*(-imasac*(b_massflux*pi) + thetap*(pfac));
4433 
4434  }
4435 }
4436 
4437 /*----------------------------------------------------------------------------*/
4459 /*----------------------------------------------------------------------------*/
4460 
4461 inline static void
4462 cs_b_imposed_conv_flux_vector(int iconvp,
4463  cs_real_t thetap,
4464  int imasac,
4465  int inc,
4466  int bc_type,
4467  int icvfli,
4468  const cs_real_t pi[restrict 3],
4469  const cs_real_t pir[restrict 3],
4470  const cs_real_t pipr[restrict 3],
4471  const cs_real_t coefap[restrict 3],
4472  const cs_real_t coefbp[restrict 3][3],
4473  const cs_real_t coface[restrict 3],
4474  const cs_real_t cofbce[restrict 3][3],
4475  cs_real_t b_massflux,
4476  cs_real_t flux[restrict 3])
4477 {
4478  cs_real_t flui, fluj, pfac;
4479 
4480  /* Computed convective flux */
4481 
4482  if (icvfli == 0) {
4483 
4484  /* Remove decentering for coupled faces */
4485  if (bc_type == CS_COUPLED_FD) {
4486  flui = 0.0;
4487  fluj = b_massflux;
4488  } else {
4489  flui = 0.5*(b_massflux +fabs(b_massflux));
4490  fluj = 0.5*(b_massflux -fabs(b_massflux));
4491  }
4492  for (int isou = 0; isou < 3; isou++) {
4493  pfac = inc*coefap[isou];
4494  for (int jsou = 0; jsou < 3; jsou++) {
4495  pfac += coefbp[isou][jsou]*pipr[jsou];
4496  }
4497  flux[isou] += iconvp*( thetap*(flui*pir[isou] + fluj*pfac)
4498  - imasac*b_massflux*pi[isou]);
4499  }
4500 
4501  /* Imposed convective flux */
4502 
4503  } else {
4504 
4505  for (int isou = 0; isou < 3; isou++) {
4506  pfac = inc*coface[isou];
4507  for (int jsou = 0; jsou < 3; jsou++) {
4508  pfac += cofbce[isou][jsou]*pipr[jsou];
4509  }
4510  flux[isou] += iconvp*( thetap*pfac
4511  - imasac*b_massflux*pi[isou]);
4512  }
4513 
4514  }
4515 }
4516 
4517 /*----------------------------------------------------------------------------*/
4537 /*----------------------------------------------------------------------------*/
4538 
4539 inline static void
4540 cs_b_upwind_flux(const int iconvp,
4541  const cs_real_t thetap,
4542  const int imasac,
4543  const int inc,
4544  const int bc_type,
4545  const cs_real_t pi,
4546  const cs_real_t pir,
4547  const cs_real_t pipr,
4548  const cs_real_t coefap,
4549  const cs_real_t coefbp,
4550  const cs_real_t b_massflux,
4551  const cs_real_t xcpp,
4552  cs_real_t *flux)
4553 {
4554  cs_real_t flui, fluj, pfac;
4555 
4556  /* Remove decentering for coupled faces */
4557  if (bc_type == CS_COUPLED_FD) {
4558  flui = 0.0;
4559  fluj = b_massflux;
4560  } else {
4561  flui = 0.5*(b_massflux +fabs(b_massflux));
4562  fluj = 0.5*(b_massflux -fabs(b_massflux));
4563  }
4564 
4565  pfac = inc*coefap + coefbp*pipr;
4566  *flux += iconvp*xcpp*(thetap*(flui*pir + fluj*pfac) -imasac*( b_massflux*pi));
4567 }
4568 
4569 /*----------------------------------------------------------------------------*/
4589 /*----------------------------------------------------------------------------*/
4590 
4591 inline static void
4592 cs_b_upwind_flux_vector(const int iconvp,
4593  const cs_real_t thetap,
4594  const int imasac,
4595  const int inc,
4596  const int bc_type,
4597  const cs_real_3_t pi,
4598  const cs_real_3_t pir,
4599  const cs_real_3_t pipr,
4600  const cs_real_3_t coefa,
4601  const cs_real_33_t coefb,
4602  const cs_real_t b_massflux,
4603  cs_real_t flux[3])
4604 {
4605  cs_real_t flui, fluj, pfac;
4606 
4607  /* Remove decentering for coupled faces */
4608  if (bc_type == CS_COUPLED_FD) {
4609  flui = 0.0;
4610  fluj = b_massflux;
4611  } else {
4612  flui = 0.5*(b_massflux +fabs(b_massflux));
4613  fluj = 0.5*(b_massflux -fabs(b_massflux));
4614  }
4615  for (int isou = 0; isou < 3; isou++) {
4616  pfac = inc*coefa[isou];
4617  for (int jsou = 0; jsou < 3; jsou++) {
4618  pfac += coefb[isou][jsou]*pipr[jsou];
4619  }
4620  flux[isou] += iconvp*( thetap*(flui*pir[isou] + fluj*pfac)
4621  - imasac*b_massflux*pi[isou]);
4622  }
4623 }
4624 
4625 /*----------------------------------------------------------------------------*/
4645 /*----------------------------------------------------------------------------*/
4646 
4647 inline static void
4648 cs_b_upwind_flux_tensor(const int iconvp,
4649  const cs_real_t thetap,
4650  const int imasac,
4651  const int inc,
4652  const int bc_type,
4653  const cs_real_6_t pi,
4654  const cs_real_6_t pir,
4655  const cs_real_6_t pipr,
4656  const cs_real_6_t coefa,
4657  const cs_real_66_t coefb,
4658  const cs_real_t b_massflux,
4659  cs_real_t flux[6])
4660 {
4661  cs_real_t flui, fluj, pfac;
4662 
4663  /* Remove decentering for coupled faces */
4664  if (bc_type == CS_COUPLED_FD) {
4665  flui = 0.0;
4666  fluj = b_massflux;
4667  } else {
4668  flui = 0.5*(b_massflux +fabs(b_massflux));
4669  fluj = 0.5*(b_massflux -fabs(b_massflux));
4670  }
4671  for (int isou = 0; isou < 6; isou++) {
4672  pfac = inc*coefa[isou];
4673  for (int jsou = 0; jsou < 6; jsou++) {
4674  pfac += coefb[isou][jsou]*pipr[jsou];
4675  }
4676  flux[isou] += iconvp*( thetap*(flui*pir[isou] + fluj*pfac)
4677  - imasac*b_massflux*pi[isou]);
4678  }
4679 }
4680 
4681 /*----------------------------------------------------------------------------*/
4694 /*----------------------------------------------------------------------------*/
4695 
4696 inline static void
4697 cs_b_diff_flux(const int idiffp,
4698  const cs_real_t thetap,
4699  const int inc,
4700  const cs_real_t pipr,
4701  const cs_real_t cofafp,
4702  const cs_real_t cofbfp,
4703  const cs_real_t b_visc,
4704  cs_real_t *flux)
4705 {
4706  cs_real_t pfacd = inc*cofafp + cofbfp*pipr;
4707  *flux += idiffp*thetap*b_visc*pfacd;
4708 }
4709 
4710 /*----------------------------------------------------------------------------*/
4723 /*----------------------------------------------------------------------------*/
4724 
4725 inline static void
4726 cs_b_diff_flux_vector(const int idiffp,
4727  const cs_real_t thetap,
4728  const int inc,
4729  const cs_real_3_t pipr,
4730  const cs_real_3_t cofaf,
4731  const cs_real_33_t cofbf,
4732  const cs_real_t b_visc,
4733  cs_real_t flux[3])
4734 {
4735  cs_real_t pfacd ;
4736  for (int isou = 0; isou < 3; isou++) {
4737  pfacd = inc*cofaf[isou];
4738  for (int jsou = 0; jsou < 3; jsou++) {
4739  pfacd += cofbf[isou][jsou]*pipr[jsou];
4740  }
4741  flux[isou] += idiffp*thetap*b_visc*pfacd;
4742  }
4743 }
4744 
4745 /*----------------------------------------------------------------------------*/
4758 /*----------------------------------------------------------------------------*/
4759 
4760 inline static void
4761 cs_b_diff_flux_tensor(const int idiffp,
4762  const cs_real_t thetap,
4763  const int inc,
4764  const cs_real_6_t pipr,
4765  const cs_real_6_t cofaf,
4766  const cs_real_66_t cofbf,
4767  const cs_real_t b_visc,
4768  cs_real_t flux[6])
4769 {
4770  cs_real_t pfacd ;
4771  for (int isou = 0; isou < 6; isou++) {
4772  pfacd = inc*cofaf[isou];
4773  for (int jsou = 0; jsou < 6; jsou++) {
4774  pfacd += cofbf[isou][jsou]*pipr[jsou];
4775  }
4776  flux[isou] += idiffp*thetap*b_visc*pfacd;
4777  }
4778 }
4779 
4780 /*----------------------------------------------------------------------------*/
4794 /*----------------------------------------------------------------------------*/
4795 
4796 inline static void
4797 cs_b_cd_steady(const cs_real_t bldfrp,
4798  const double relaxp,
4799  const cs_real_3_t diipb,
4800  const cs_real_3_t gradi,
4801  const cs_real_t pi,
4802  const cs_real_t pia,
4803  cs_real_t *pir,
4804  cs_real_t *pipr)
4805 {
4806  cs_real_t recoi;
4807 
4808  cs_b_compute_quantities(diipb,
4809  gradi,
4810  bldfrp,
4811  &recoi);
4812 
4813  cs_b_relax_c_val(relaxp,
4814  pi,
4815  pia,
4816  recoi,
4817  pir,
4818  pipr);
4819 }
4820 
4821 /*----------------------------------------------------------------------------*/
4835 /*----------------------------------------------------------------------------*/
4836 
4837 inline static void
4838 cs_b_cd_steady_vector(const cs_real_t bldfrp,
4839  const double relaxp,
4840  const cs_real_3_t diipb,
4841  const cs_real_33_t gradi,
4842  const cs_real_3_t pi,
4843  const cs_real_3_t pia,
4844  cs_real_t pir[3],
4845  cs_real_t pipr[3])
4846 {
4847  cs_real_3_t recoi;
4848 
4849  cs_b_compute_quantities_vector(diipb,
4850  gradi,
4851  bldfrp,
4852  recoi);
4853 
4854  cs_b_relax_c_val_vector(relaxp,
4855  pi,
4856  pia,
4857  recoi,
4858  pir,
4859  pipr);
4860 }
4861 
4862 /*----------------------------------------------------------------------------*/
4876 /*----------------------------------------------------------------------------*/
4877 
4878 inline static void
4879 cs_b_cd_steady_tensor(const cs_real_t bldfrp,
4880  const double relaxp,
4881  const cs_real_3_t diipb,
4882  const cs_real_63_t gradi,
4883  const cs_real_6_t pi,
4884  const cs_real_6_t pia,
4885  cs_real_t pir[6],
4886  cs_real_t pipr[6])
4887 {
4888  cs_real_6_t recoi;
4889 
4890  cs_b_compute_quantities_tensor(diipb,
4891  gradi,
4892  bldfrp,
4893  recoi);
4894 
4895  cs_b_relax_c_val_tensor(relaxp,
4896  pi,
4897  pia,
4898  recoi,
4899  pir,
4900  pipr);
4901 }
4902 
4903 /*----------------------------------------------------------------------------*/
4914 /*----------------------------------------------------------------------------*/
4915 
4916 inline static void
4917 cs_b_cd_unsteady(const cs_real_t bldfrp,
4918  const cs_real_3_t diipb,
4919  const cs_real_3_t gradi,
4920  const cs_real_t pi,
4921  cs_real_t *pip)
4922 {
4923  cs_real_t recoi;
4924 
4925  cs_b_compute_quantities(diipb,
4926  gradi,
4927  bldfrp,
4928  &recoi);
4929 
4930  *pip = pi + recoi;
4931 }
4932 
4933 /*----------------------------------------------------------------------------*/
4944 /*----------------------------------------------------------------------------*/
4945 
4946 inline static void
4947 cs_b_cd_unsteady_vector(const cs_real_t bldfrp,
4948  const cs_real_3_t diipb,
4949  const cs_real_33_t gradi,
4950  const cs_real_3_t pi,
4951  cs_real_t pip[3])
4952 {
4953  cs_real_3_t recoi;
4954 
4955  cs_b_compute_quantities_vector(diipb,
4956  gradi,
4957  bldfrp,
4958  recoi);
4959 
4960  for (int isou = 0; isou < 3; isou++)
4961  pip[isou] = pi[isou] + recoi[isou];
4962 }
4963 
4964 /*----------------------------------------------------------------------------*/
4975 /*----------------------------------------------------------------------------*/
4976 
4977 inline static void
4978 cs_b_cd_unsteady_tensor(const cs_real_t bldfrp,
4979  const cs_real_3_t diipb,
4980  const cs_real_63_t gradi,
4981  const cs_real_6_t pi,
4982  cs_real_t pip[6])
4983 {
4984  cs_real_6_t recoi;
4985 
4986  cs_b_compute_quantities_tensor(diipb,
4987  gradi,
4988  bldfrp,
4989  recoi);
4990 
4991  for(int isou = 0; isou< 6; isou++)
4992  pip[isou] = pi[isou] + recoi[isou];
4993 }
4994 
4995 /*----------------------------------------------------------------------------*/
5006 /*----------------------------------------------------------------------------*/
5007 
5008 inline static void
5009 cs_b_diff_flux_coupling(int idiffp,
5010  cs_real_t pi,
5011  cs_real_t pj,
5012  cs_real_t b_visc,
5013  cs_real_t *fluxi)
5014 {
5015  *fluxi += idiffp*b_visc*(pi - pj);
5016 }
5017 
5018 /*----------------------------------------------------------------------------*/
5029 /*----------------------------------------------------------------------------*/
5030 
5031 inline static void
5032 cs_b_diff_flux_coupling_vector(int idiffp,
5033  const cs_real_t pi[3],
5034  const cs_real_t pj[3],
5035  cs_real_t b_visc,
5036  cs_real_t fluxi[3])
5037 {
5038  for (int k = 0; k < 3; k++)
5039  fluxi[k] += idiffp*b_visc*(pi[k] - pj[k]);
5040 }
5041 
5042 
5043 /*============================================================================
5044  * Public function prototypes for Fortran API
5045  *============================================================================*/
5046 
5047 /*----------------------------------------------------------------------------
5048  * Wrapper to cs_face_diffusion_potential
5049  *----------------------------------------------------------------------------*/
5050 
5051 void CS_PROCF (itrmas, ITRMAS)
5052 (
5053  const int *const f_id,
5054  const int *const init,
5055  const int *const inc,
5056  const int *const imrgra,
5057  const int *const iccocg,
5058  const int *const nswrgp,
5059  const int *const imligp,
5060  const int *const iphydp,
5061  const int *const iwgrp,
5062  const int *const iwarnp,
5063  const cs_real_t *const epsrgp,
5064  const cs_real_t *const climgp,
5065  const cs_real_t *const extrap,
5066  cs_real_3_t frcxt[],
5067  cs_real_t pvar[],
5068  const cs_real_t coefap[],
5069  const cs_real_t coefbp[],
5070  const cs_real_t cofafp[],
5071  const cs_real_t cofbfp[],
5072  const cs_real_t i_visc[],
5073  const cs_real_t b_visc[],
5074  cs_real_t visel[],
5075  cs_real_t i_massflux[],
5076  cs_real_t b_massflux[]
5077 );
5078 
5079 /*----------------------------------------------------------------------------
5080  * Wrapper to cs_face_anisotropic_diffusion_potential
5081  *----------------------------------------------------------------------------*/
5082 
5083 void CS_PROCF (itrmav, ITRMAV)
5084 (
5085  const int *const f_id,
5086  const int *const init,
5087  const int *const inc,
5088  const int *const imrgra,
5089  const int *const iccocg,
5090  const int *const nswrgp,
5091  const int *const imligp,
5092  const int *const ircflp,
5093  const int *const iphydp,
5094  const int *const iwgrp,
5095  const int *const iwarnp,
5096  const cs_real_t *const epsrgp,
5097  const cs_real_t *const climgp,
5098  const cs_real_t *const extrap,
5099  cs_real_3_t frcxt[],
5100  cs_real_t pvar[],
5101  const cs_real_t coefap[],
5102  const cs_real_t coefbp[],
5103  const cs_real_t cofafp[],
5104  const cs_real_t cofbfp[],
5105  const cs_real_t i_visc[],
5106  const cs_real_t b_visc[],
5107  cs_real_6_t viscel[],
5108  const cs_real_2_t weighf[],
5109  const cs_real_t weighb[],
5110  cs_real_t i_massflux[],
5111  cs_real_t b_massflux[]
5112 );
5113 
5114 /*----------------------------------------------------------------------------
5115  * Wrapper to cs_diffusion_potential
5116  *----------------------------------------------------------------------------*/
5117 
5118 void CS_PROCF (itrgrp, ITRGRP)
5119 (
5120  const int *const f_id,
5121  const int *const init,
5122  const int *const inc,
5123  const int *const imrgra,
5124  const int *const iccocg,
5125  const int *const nswrgp,
5126  const int *const imligp,
5127  const int *const iphydp,
5128  const int *const iwgrp,
5129  const int *const iwarnp,
5130  const cs_real_t *const epsrgp,
5131  const cs_real_t *const climgp,
5132  const cs_real_t *const extrap,
5133  cs_real_3_t frcxt[],
5134  cs_real_t pvar[],
5135  const cs_real_t coefap[],
5136  const cs_real_t coefbp[],
5137  const cs_real_t cofafp[],
5138  const cs_real_t cofbfp[],
5139  const cs_real_t i_visc[],
5140  const cs_real_t b_visc[],
5141  cs_real_t visel[],
5142  cs_real_t diverg[]
5143 );
5144 
5145 /*----------------------------------------------------------------------------
5146  * Wrapper to cs_anisotropic_diffusion_potential
5147  *----------------------------------------------------------------------------*/
5148 
5149 void CS_PROCF (itrgrv, ITRGRV)
5150 (
5151  const int *const f_id,
5152  const int *const init,
5153  const int *const inc,
5154  const int *const imrgra,
5155  const int *const iccocg,
5156  const int *const nswrgp,
5157  const int *const imligp,
5158  const int *const ircflp,
5159  const int *const iphydp,
5160  const int *const iwgrp,
5161  const int *const iwarnp,
5162  const cs_real_t *const epsrgp,
5163  const cs_real_t *const climgp,
5164  const cs_real_t *const extrap,
5165  cs_real_3_t frcxt[],
5166  cs_real_t pvar[],
5167  const cs_real_t coefap[],
5168  const cs_real_t coefbp[],
5169  const cs_real_t cofafp[],
5170  const cs_real_t cofbfp[],
5171  const cs_real_t i_visc[],
5172  const cs_real_t b_visc[],
5173  cs_real_6_t viscel[],
5174  const cs_real_2_t weighf[],
5175  const cs_real_t weighb[],
5176  cs_real_t diverg[]
5177 );
5178 
5179 /*=============================================================================
5180  * Public function prototypes
5181  *============================================================================*/
5182 
5183 /*----------------------------------------------------------------------------*/
5202 /*----------------------------------------------------------------------------*/
5203 
5204 void
5205 cs_slope_test_gradient(int f_id,
5206  int inc,
5207  cs_halo_type_t halo_type,
5208  const cs_real_3_t *grad,
5209  cs_real_3_t *grdpa,
5210  const cs_real_t *pvar,
5211  const cs_real_t *coefap,
5212  const cs_real_t *coefbp,
5213  const cs_real_t *i_massflux);
5214 
5215 /*----------------------------------------------------------------------------*/
5231 /*----------------------------------------------------------------------------*/
5232 
5233 void
5234 cs_upwind_gradient(const int f_id,
5235  const int inc,
5236  const cs_halo_type_t halo_type,
5237  const cs_real_t coefap[],
5238  const cs_real_t coefbp[],
5239  const cs_real_t i_massflux[],
5240  const cs_real_t b_massflux[],
5241  const cs_real_t *restrict pvar,
5242  cs_real_3_t *restrict grdpa);
5243 
5244 /*----------------------------------------------------------------------------*/
5262 /*----------------------------------------------------------------------------*/
5263 
5264 void
5265 cs_slope_test_gradient_vector(const int inc,
5266  const cs_halo_type_t halo_type,
5267  const cs_real_33_t *grad,
5268  cs_real_33_t *grdpa,
5269  const cs_real_3_t *pvar,
5270  const cs_real_3_t *coefa,
5271  const cs_real_33_t *coefb,
5272  const cs_real_t *i_massflux);
5273 
5274 /*----------------------------------------------------------------------------*/
5292 /*----------------------------------------------------------------------------*/
5293 
5294 void
5295 cs_slope_test_gradient_tensor(const int inc,
5296  const cs_halo_type_t halo_type,
5297  const cs_real_63_t *grad,
5298  cs_real_63_t *grdpa,
5299  const cs_real_6_t *pvar,
5300  const cs_real_6_t *coefa,
5301  const cs_real_66_t *coefb,
5302  const cs_real_t *i_massflux);
5303 
5304 /*----------------------------------------------------------------------------*/
5313 /*----------------------------------------------------------------------------*/
5314 
5315 void
5316 cs_beta_limiter_building(int f_id,
5317  int inc,
5318  const cs_real_t rovsdt[]);
5319 
5320 /*----------------------------------------------------------------------------*/
5372 /*----------------------------------------------------------------------------*/
5373 
5374 void
5375 cs_convection_diffusion_scalar(int idtvar,
5376  int f_id,
5377  const cs_var_cal_opt_t var_cal_opt,
5378  int icvflb,
5379  int inc,
5380  int iccocg,
5381  int imasac,
5382  cs_real_t *restrict pvar,
5383  const cs_real_t *restrict pvara,
5384  const int icvfli[],
5385  const cs_real_t coefap[],
5386  const cs_real_t coefbp[],
5387  const cs_real_t cofafp[],
5388  const cs_real_t cofbfp[],
5389  const cs_real_t i_massflux[],
5390  const cs_real_t b_massflux[],
5391  const cs_real_t i_visc[],
5392  const cs_real_t b_visc[],
5393  cs_real_t *restrict rhs);
5394 
5395 /*----------------------------------------------------------------------------*/
5434 /*----------------------------------------------------------------------------*/
5435 
5436 void
5437 cs_face_convection_scalar(int idtvar,
5438  int f_id,
5439  const cs_var_cal_opt_t var_cal_opt,
5440  int icvflb,
5441  int inc,
5442  int iccocg,
5443  int imasac,
5444  cs_real_t *restrict pvar,
5445  const cs_real_t *restrict pvara,
5446  const int icvfli[],
5447  const cs_real_t coefap[],
5448  const cs_real_t coefbp[],
5449  const cs_real_t i_massflux[],
5450  const cs_real_t b_massflux[],
5451  cs_real_2_t i_conv_flux[],
5452  cs_real_t b_conv_flux[]);
5453 
5454 /*----------------------------------------------------------------------------*/
5514 /*----------------------------------------------------------------------------*/
5515 
5516 void
5517 cs_convection_diffusion_vector(int idtvar,
5518  int f_id,
5519  const cs_var_cal_opt_t var_cal_opt,
5520  int icvflb,
5521  int inc,
5522  int ivisep,
5523  int imasac,
5524  cs_real_3_t *restrict pvar,
5525  const cs_real_3_t *restrict pvara,
5526  const int icvfli[],
5527  const cs_real_3_t coefav[],
5528  const cs_real_33_t coefbv[],
5529  const cs_real_3_t cofafv[],
5530  const cs_real_33_t cofbfv[],
5531  const cs_real_t i_massflux[],
5532  const cs_real_t b_massflux[],
5533  const cs_real_t i_visc[],
5534  const cs_real_t b_visc[],
5535  const cs_real_t i_secvis[],
5536  const cs_real_t b_secvis[],
5537  cs_real_3_t *restrict rhs);
5538 
5539 /*----------------------------------------------------------------------------*/
5584 /*----------------------------------------------------------------------------*/
5585 
5586 void
5587 cs_convection_diffusion_tensor(int idtvar,
5588  int f_id,
5589  const cs_var_cal_opt_t var_cal_opt,
5590  int icvflb,
5591  int inc,
5592  int imasac,
5593  cs_real_6_t *restrict pvar,
5594  const cs_real_6_t *restrict pvara,
5595  const cs_real_6_t coefa[],
5596  const cs_real_66_t coefb[],
5597  const cs_real_6_t cofaf[],
5598  const cs_real_66_t cofbf[],
5599  const cs_real_t i_massflux[],
5600  const cs_real_t b_massflux[],
5601  const cs_real_t i_visc[],
5602  const cs_real_t b_visc[],
5603  cs_real_6_t *restrict rhs);
5604 
5605 /*----------------------------------------------------------------------------*/
5651 /*----------------------------------------------------------------------------*/
5652 
5653 void
5654 cs_convection_diffusion_thermal(int idtvar,
5655  int f_id,
5656  const cs_var_cal_opt_t var_cal_opt,
5657  int inc,
5658  int iccocg,
5659  int imasac,
5660  cs_real_t *restrict pvar,
5661  const cs_real_t *restrict pvara,
5662  const cs_real_t coefap[],
5663  const cs_real_t coefbp[],
5664  const cs_real_t cofafp[],
5665  const cs_real_t cofbfp[],
5666  const cs_real_t i_massflux[],
5667  const cs_real_t b_massflux[],
5668  const cs_real_t i_visc[],
5669  const cs_real_t b_visc[],
5670  const cs_real_t xcpp[],
5671  cs_real_t *restrict rhs);
5672 
5673 /*----------------------------------------------------------------------------*/
5721 /*----------------------------------------------------------------------------*/
5722 
5723 void
5724 cs_anisotropic_diffusion_scalar(int idtvar,
5725  int f_id,
5726  const cs_var_cal_opt_t var_cal_opt,
5727  int inc,
5728  int iccocg,
5729  cs_real_t *restrict pvar,
5730  const cs_real_t *restrict pvara,
5731  const cs_real_t coefap[],
5732  const cs_real_t coefbp[],
5733  const cs_real_t cofafp[],
5734  const cs_real_t cofbfp[],
5735  const cs_real_t i_visc[],
5736  const cs_real_t b_visc[],
5737  cs_real_6_t *restrict viscel,
5738  const cs_real_2_t weighf[],
5739  const cs_real_t weighb[],
5740  cs_real_t *restrict rhs);
5741 
5742 /*-----------------------------------------------------------------------------*/
5792 /*----------------------------------------------------------------------------*/
5793 
5794 void
5795 cs_anisotropic_left_diffusion_vector(int idtvar,
5796  int f_id,
5797  const cs_var_cal_opt_t var_cal_opt,
5798  int inc,
5799  int ivisep,
5800  cs_real_3_t *restrict pvar,
5801  const cs_real_3_t *restrict pvara,
5802  const cs_real_3_t coefav[],
5803  const cs_real_33_t coefbv[],
5804  const cs_real_3_t cofafv[],
5805  const cs_real_33_t cofbfv[],
5806  const cs_real_33_t i_visc[],
5807  const cs_real_t b_visc[],
5808  const cs_real_t i_secvis[],
5809  cs_real_3_t *restrict rhs);
5810 
5811 /*-----------------------------------------------------------------------------*/
5856 /*----------------------------------------------------------------------------*/
5857 
5858 void
5859 cs_anisotropic_right_diffusion_vector(int idtvar,
5860  int f_id,
5861  const cs_var_cal_opt_t var_cal_opt,
5862  int inc,
5863  cs_real_3_t *restrict pvar,
5864  const cs_real_3_t *restrict pvara,
5865  const cs_real_3_t coefav[],
5866  const cs_real_33_t coefbv[],
5867  const cs_real_3_t cofafv[],
5868  const cs_real_33_t cofbfv[],
5869  const cs_real_t i_visc[],
5870  const cs_real_t b_visc[],
5871  cs_real_6_t *restrict viscel,
5872  const cs_real_2_t weighf[],
5873  const cs_real_t weighb[],
5874  cs_real_3_t *restrict rhs);
5875 
5876 /*----------------------------------------------------------------------------*/
5920 /*----------------------------------------------------------------------------*/
5921 
5922 void
5923 cs_anisotropic_diffusion_tensor(int idtvar,
5924  int f_id,
5925  const cs_var_cal_opt_t var_cal_opt,
5926  int inc,
5927  cs_real_6_t *restrict pvar,
5928  const cs_real_6_t *restrict pvara,
5929  const cs_real_6_t coefa[],
5930  const cs_real_66_t coefb[],
5931  const cs_real_6_t cofaf[],
5932  const cs_real_66_t cofbf[],
5933  const cs_real_t i_visc[],
5934  const cs_real_t b_visc[],
5935  cs_real_6_t *restrict viscel,
5936  const cs_real_2_t weighf[],
5937  const cs_real_t weighb[],
5938  cs_real_6_t *restrict rhs);
5939 
5940 /*----------------------------------------------------------------------------*/
5999 /*----------------------------------------------------------------------------*/
6000 
6001 void
6002 cs_face_diffusion_potential(const int f_id,
6003  const cs_mesh_t *m,
6004  cs_mesh_quantities_t *fvq,
6005  int init,
6006  int inc,
6007  int imrgra,
6008  int iccocg,
6009  int nswrgp,
6010  int imligp,
6011  int iphydp,
6012  int iwgrp,
6013  int iwarnp,
6014  double epsrgp,
6015  double climgp,
6016  cs_real_3_t *restrict frcxt,
6017  cs_real_t *restrict pvar,
6018  const cs_real_t coefap[],
6019  const cs_real_t coefbp[],
6020  const cs_real_t cofafp[],
6021  const cs_real_t cofbfp[],
6022  const cs_real_t i_visc[],
6023  const cs_real_t b_visc[],
6024  cs_real_t *restrict visel,
6025  cs_real_t *restrict i_massflux,
6026  cs_real_t *restrict b_massflux);
6027 
6028 /*----------------------------------------------------------------------------*/
6098 /*----------------------------------------------------------------------------*/
6099 
6100 void
6101 cs_face_anisotropic_diffusion_potential(const int f_id,
6102  const cs_mesh_t *m,
6103  cs_mesh_quantities_t *fvq,
6104  int init,
6105  int inc,
6106  int imrgra,
6107  int iccocg,
6108  int nswrgp,
6109  int imligp,
6110  int ircflp,
6111  int iphydp,
6112  int iwgrp,
6113  int iwarnp,
6114  double epsrgp,
6115  double climgp,
6116  cs_real_3_t *restrict frcxt,
6117  cs_real_t *restrict pvar,
6118  const cs_real_t coefap[],
6119  const cs_real_t coefbp[],
6120  const cs_real_t cofafp[],
6121  const cs_real_t cofbfp[],
6122  const cs_real_t i_visc[],
6123  const cs_real_t b_visc[],
6124  cs_real_6_t *restrict viscel,
6125  const cs_real_2_t weighf[],
6126  const cs_real_t weighb[],
6127  cs_real_t *restrict i_massflux,
6128  cs_real_t *restrict b_massflux);
6129 
6130 /*----------------------------------------------------------------------------*/
6188 /*----------------------------------------------------------------------------*/
6189 
6190 void
6191 cs_diffusion_potential(const int f_id,
6192  const cs_mesh_t *m,
6193  cs_mesh_quantities_t *fvq,
6194  int init,
6195  int inc,
6196  int imrgra,
6197  int iccocg,
6198  int nswrgp,
6199  int imligp,
6200  int iphydp,
6201  int iwgrp,
6202  int iwarnp,
6203  double epsrgp,
6204  double climgp,
6205  cs_real_3_t *restrict frcxt,
6206  cs_real_t *restrict pvar,
6207  const cs_real_t coefap[],
6208  const cs_real_t coefbp[],
6209  const cs_real_t cofafp[],
6210  const cs_real_t cofbfp[],
6211  const cs_real_t i_visc[],
6212  const cs_real_t b_visc[],
6213  cs_real_t visel[],
6214  cs_real_t *restrict diverg);
6215 
6216 /*----------------------------------------------------------------------------*/
6287 /*----------------------------------------------------------------------------*/
6288 
6289 void
6290 cs_anisotropic_diffusion_potential(const int f_id,
6291  const cs_mesh_t *m,
6292  cs_mesh_quantities_t *fvq,
6293  int init,
6294  int inc,
6295  int imrgra,
6296  int iccocg,
6297  int nswrgp,
6298  int imligp,
6299  int ircflp,
6300  int iphydp,
6301  int iwgrp,
6302  int iwarnp,
6303  double epsrgp,
6304  double climgp,
6305  cs_real_3_t *restrict frcxt,
6306  cs_real_t *restrict pvar,
6307  const cs_real_t coefap[],
6308  const cs_real_t coefbp[],
6309  const cs_real_t cofafp[],
6310  const cs_real_t cofbfp[],
6311  const cs_real_t i_visc[],
6312  const cs_real_t b_visc[],
6313  cs_real_6_t *restrict viscel,
6314  const cs_real_2_t weighf[],
6315  const cs_real_t weighb[],
6316  cs_real_t *restrict diverg);
6317 
6318 /*----------------------------------------------------------------------------*/
6319 
6320 END_C_DECLS
6321 
6322 #endif /* __CS_CONVECTION_DIFFUSION_H__ */
k
Definition: cs_field_pointer.h:70
itrgrv
void itrgrv(const int *const f_id, const int *const init, const int *const inc, const int *const imrgra, const int *const iccocg, const int *const nswrgp, const int *const imligp, const int *const ircflp, const int *const iphydp, const int *const iwgrp, const int *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:812
cs_convection_diffusion_vector
void cs_convection_diffusion_vector(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int icvflb, int inc, int ivisep, int imasac, cs_real_3_t *restrict pvar, const cs_real_3_t *restrict pvara, const int 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 i_secvis[], const cs_real_t b_secvis[], 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:4169
mesh::diipb
double precision, dimension(:,:), pointer diipb
Definition: mesh.f90:212
restrict
#define restrict
Definition: cs_defs.h:142
cs_real_2_t
cs_real_t cs_real_2_t[2]
vector of 2 floating-point values
Definition: cs_defs.h:334
cs_base.h
cs_real_6_t
cs_real_t cs_real_6_t[6]
vector of 6 floating-point values
Definition: cs_defs.h:337
cs_field_key_id_try
int cs_field_key_id_try(const char *name)
Return an id associated with a given key name if present.
Definition: cs_field.c:2524
cs_anisotropic_diffusion_tensor
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_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:10164
cs_nvd_type_t
cs_nvd_type_t
Definition: cs_convection_diffusion.h:59
cs_face_diffusion_potential
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 iwgrp, int iwarnp, double epsrgp, double climgp, 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_t *restrict visel, 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:10781
cs_anisotropic_diffusion_scalar
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:8154
cs_field_t
Field descriptor.
Definition: cs_field.h:125
cs_equation_param_t
Set of parameters to handle an unsteady convection-diffusion-reaction equation with term sources...
Definition: cs_equation_param.h:177
CS_ABS
#define CS_ABS(a)
Definition: cs_defs.h:472
cs_field_by_id
cs_field_t * cs_field_by_id(int id)
Return a pointer to a field based on its id.
Definition: cs_field.c:2314
cs_real_66_t
cs_real_t cs_real_66_t[6][6]
6x6 matrix of floating-point values
Definition: cs_defs.h:343
itrgrp
void itrgrp(const int *const f_id, const int *const init, const int *const inc, const int *const imrgra, const int *const iccocg, const int *const nswrgp, const int *const imligp, const int *const iphydp, const int *const iwgrp, const int *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_t visel[], cs_real_t diverg[])
Definition: cs_convection_diffusion.c:750
CS_NVD_STOIC
Definition: cs_convection_diffusion.h:68
cs_slope_test_gradient_tensor
void cs_slope_test_gradient_tensor(const int inc, const cs_halo_type_t halo_type, const cs_real_63_t *grad, cs_real_63_t *grdpa, const 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.c:1328
cstnum::pi
double precision pi
value with 16 digits
Definition: cstnum.f90:48
cs_defs.h
BEGIN_C_DECLS
#define BEGIN_C_DECLS
Definition: cs_defs.h:510
cs_face_convection_scalar
void cs_face_convection_scalar(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int icvflb, int inc, int iccocg, int imasac, cs_real_t *restrict pvar, const cs_real_t *restrict pvara, const int icvfli[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], cs_real_2_t i_conv_flux[], cs_real_t b_conv_flux[])
Update face flux with convection contribution of a standard transport equation of a scalar field ...
Definition: cs_convection_diffusion.c:3058
cs_math_epzero
const cs_real_t cs_math_epzero
cs_glob_mesh
cs_mesh_t * cs_glob_mesh
CS_UNUSED
#define CS_UNUSED(x)
Definition: cs_defs.h:496
cs_field_pointer.h
cs_slope_test_gradient_vector
void cs_slope_test_gradient_vector(const int inc, const cs_halo_type_t halo_type, const cs_real_33_t *grad, cs_real_33_t *grdpa, const cs_real_3_t *pvar, const cs_real_3_t *coefa, const cs_real_33_t *coefb, const cs_real_t *i_massflux)
Compute the upwind gradient used in the slope tests.
Definition: cs_convection_diffusion.c:1173
cpincl::beta
double precision, dimension(ncharm), save beta
Definition: cpincl.f90:99
optcal::idtvar
integer(c_int), pointer, save idtvar
option for a variable time step
Definition: optcal.f90:383
cs_anisotropic_left_diffusion_vector
void cs_anisotropic_left_diffusion_vector(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int inc, int ivisep, cs_real_3_t *restrict pvar, const cs_real_3_t *restrict pvara, 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 i_secvis[], cs_real_3_t *restrict rhs)
Add explicit part of the terms of diffusion by a left-multiplying symmetric tensorial diffusivity for...
Definition: cs_convection_diffusion.c:8907
p
Definition: cs_field_pointer.h:67
cs_math_3_length_unitv
void cs_math_3_length_unitv(const cs_real_t xa[3], const cs_real_t xb[3], cs_real_t *len, cs_real_3_t unitv)
Compute the length (Euclidean norm) between two points xa and xb in a Cartesian coordinate system of ...
Definition: cs_math.c:438
CS_NVD_CUBISTA
Definition: cs_convection_diffusion.h:63
cs_field_t::val
cs_real_t * val
Definition: cs_field.h:146
CS_F_
#define CS_F_(e)
Macro used to return a field pointer by its enumerated value.
Definition: cs_field_pointer.h:51
cs_parameters.h
CS_NVD_MINMOD
Definition: cs_convection_diffusion.h:66
cfpoin::icvfli
integer, dimension(:), allocatable icvfli
boundary convection flux indicator of a Rusanov or an analytical flux (some boundary contributions of...
Definition: cfpoin.f90:52
cs_real_t
double cs_real_t
Floating-point value.
Definition: cs_defs.h:322
itrmas
void itrmas(const int *const f_id, const int *const init, const int *const inc, const int *const imrgra, const int *const iccocg, const int *const nswrgp, const int *const imligp, const int *const iphydp, const int *const iwgrp, const int *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_t visel[], cs_real_t i_massflux[], cs_real_t b_massflux[])
Definition: cs_convection_diffusion.c:616
cs_mesh_t
Definition: cs_mesh.h:84
cs_field_get_key_int
int cs_field_get_key_int(const cs_field_t *f, int key_id)
Return a integer value for a given key associated with a field.
Definition: cs_field.c:2991
cs_field_key_id
int cs_field_key_id(const char *name)
Return an id associated with a given key name.
Definition: cs_field.c:2497
var_cal_opt
cs_equation_param_t * var_cal_opt
Definition: keywords.h:135
cs_halo_sync_var
void cs_halo_sync_var(const cs_halo_t *halo, cs_halo_type_t sync_mode, cs_real_t var[])
Definition: cs_halo.c:1903
cs_mesh_quantities.h
cpincl::b2
double precision, dimension(ncharm), save b2
Definition: cpincl.f90:233
cs_equation_param_t::iconv
int iconv
Definition: cs_equation_param.h:437
CS_NVD_MUSCL
Definition: cs_convection_diffusion.h:65
cs_numbering_t::group_index
cs_lnum_t * group_index
Definition: cs_numbering.h:98
cs_anisotropic_right_diffusion_vector
void cs_anisotropic_right_diffusion_vector(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int inc, cs_real_3_t *restrict pvar, const cs_real_3_t *restrict pvara, 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_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_3_t *restrict rhs)
Add explicit part of the terms of diffusion by a right-multiplying symmetric tensorial diffusivity fo...
Definition: cs_convection_diffusion.c:9453
dt
Definition: cs_field_pointer.h:65
cs_glob_mesh_quantities
cs_mesh_quantities_t * cs_glob_mesh_quantities
cs_mesh_quantities_t
Definition: cs_mesh_quantities.h:89
cs_equation_param_t::blencv
double blencv
Definition: cs_equation_param.h:457
CS_MIN
#define CS_MIN(a, b)
Definition: cs_defs.h:473
cs_halo_type_t
cs_halo_type_t
Definition: cs_halo.h:56
cs_halo.h
cs_face_anisotropic_diffusion_potential
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 iwgrp, int iwarnp, double epsrgp, double climgp, 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:11130
CS_NVD_SMART
Definition: cs_convection_diffusion.h:62
cs_lnum_2_t
cs_lnum_t cs_lnum_2_t[2]
vector of 2 local mesh-entity ids
Definition: cs_defs.h:328
cs_real_3_t
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:335
cs_numbering_t::n_groups
int n_groups
Definition: cs_numbering.h:90
cs_equation_param_t::isstpc
int isstpc
Definition: cs_equation_param.h:447
cs_anisotropic_diffusion_potential
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 iwgrp, int iwarnp, double epsrgp, double climgp, 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:11993
CS_NVD_N_TYPES
Definition: cs_convection_diffusion.h:74
cs_slope_test_gradient
void cs_slope_test_gradient(int f_id, int inc, cs_halo_type_t halo_type, const cs_real_3_t *grad, cs_real_3_t *grdpa, const 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.c:901
cs_convection_diffusion_tensor
void cs_convection_diffusion_tensor(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int icvflb, int inc, int imasac, cs_real_6_t *restrict pvar, const cs_real_6_t *restrict pvara, 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:5857
coincl::fmin
double precision, save fmin
Definition: coincl.f90:133
cs_lnum_t
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:316
cs_mesh_t::n_cells_with_ghosts
cs_lnum_t n_cells_with_ghosts
Definition: cs_mesh.h:150
CS_NVD_WASEB
Definition: cs_convection_diffusion.h:70
CS_HALO_STANDARD
Definition: cs_halo.h:58
cs_math.h
END_C_DECLS
#define END_C_DECLS
Definition: cs_defs.h:511
CS_NVD_CLAM
Definition: cs_convection_diffusion.h:67
CS_NVD_VOF_CICSAM
Definition: cs_convection_diffusion.h:72
cs_mesh_t::b_face_numbering
cs_numbering_t * b_face_numbering
Definition: cs_mesh.h:162
cs_mesh_t::i_face_numbering
cs_numbering_t * i_face_numbering
Definition: cs_mesh.h:161
CS_PROCF
#define CS_PROCF(x, y)
Definition: cs_defs.h:524
CS_MAX
#define CS_MAX(a, b)
Definition: cs_defs.h:474
CS_COUPLED_FD
Definition: cs_parameters.h:93
cs_numbering_t::n_threads
int n_threads
Definition: cs_numbering.h:89
cs_diffusion_potential
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 iwgrp, int iwarnp, double epsrgp, double climgp, 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_t visel[], 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:11615
numvar::kimasf
integer, save kimasf
interior and boundary convective mass flux key ids of the variables
Definition: numvar.f90:187
cs_real_33_t
cs_real_t cs_real_33_t[3][3]
3x3 matrix of floating-point values
Definition: cs_defs.h:342
CS_NVD_VOF_HRIC
Definition: cs_convection_diffusion.h:71
CS_NVD_GAMMA
Definition: cs_convection_diffusion.h:61
CS_NVD_SUPERBEE
Definition: cs_convection_diffusion.h:64
cs_convection_diffusion_scalar
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 imasac, cs_real_t *restrict pvar, const cs_real_t *restrict pvara, const int 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:1625
cs_real_63_t
cs_real_t cs_real_63_t[6][3]
Definition: cs_defs.h:350
cs_mesh_quantities_t::cell_vol
cs_real_t * cell_vol
Definition: cs_mesh_quantities.h:92
numvar::kbmasf
integer, save kbmasf
Definition: numvar.f90:187
cs_mesh_t::i_face_cells
cs_lnum_2_t * i_face_cells
Definition: cs_mesh.h:110
itrmav
void itrmav(const int *const f_id, const int *const init, const int *const inc, const int *const imrgra, const int *const iccocg, const int *const nswrgp, const int *const imligp, const int *const ircflp, const int *const iphydp, const int *const iwgrp, const int *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:680
cs_beta_limiter_building
void cs_beta_limiter_building(int f_id, int inc, const cs_real_t rovsdt[])
Compute the beta blending coefficient of the beta limiter (ensuring preservation of a given min/max p...
Definition: cs_convection_diffusion.c:1474
cs_convection_diffusion_thermal
void cs_convection_diffusion_thermal(int idtvar, int f_id, const cs_var_cal_opt_t var_cal_opt, int inc, int iccocg, int imasac, 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_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:6820
cpincl::b1
double precision, dimension(ncharm), save b1
Definition: cpincl.f90:233
cs_upwind_gradient
void cs_upwind_gradient(const int f_id, const int inc, const cs_halo_type_t halo_type, const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t *restrict pvar, cs_real_3_t *restrict grdpa)
Compute the upwind gradient in order to cope with SOLU schemes observed in the litterature.
Definition: cs_convection_diffusion.c:1043
CS_NVD_OSHER
Definition: cs_convection_diffusion.h:69
optcal::imrgra
integer(c_int), pointer, save imrgra
type of gradient reconstruction
Definition: optcal.f90:260
cs_mesh_t::halo
cs_halo_t * halo
Definition: cs_mesh.h:155
cs_mesh_t::b_face_cells
cs_lnum_t * b_face_cells
Definition: cs_mesh.h:111
CS_NVD_VOF_STACS
Definition: cs_convection_diffusion.h:73