9.0
general documentation
Loading...
Searching...
No Matches
cs_balance.h
Go to the documentation of this file.
1#ifndef __CS_BALANCE_H__
2#define __CS_BALANCE_H__
3
4/*============================================================================
5 * Building of the right hand side for a transport of a field.
6 *============================================================================*/
7
8/*
9 This file is part of code_saturne, a general-purpose CFD tool.
10
11 Copyright (C) 1998-2025 EDF S.A.
12
13 This program is free software; you can redistribute it and/or modify it under
14 the terms of the GNU General Public License as published by the Free Software
15 Foundation; either version 2 of the License, or (at your option) any later
16 version.
17
18 This program is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
21 details.
22
23 You should have received a copy of the GNU General Public License along with
24 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
25 Street, Fifth Floor, Boston, MA 02110-1301, USA.
26*/
27
28/*----------------------------------------------------------------------------*/
29
30#include "base/cs_defs.h"
31
32/*----------------------------------------------------------------------------
33 * Standard C library headers
34 *----------------------------------------------------------------------------*/
35
36/*----------------------------------------------------------------------------
37 * Local headers
38 *----------------------------------------------------------------------------*/
39
40#include "base/cs_parameters.h"
41
42/*----------------------------------------------------------------------------*/
44
45/*=============================================================================
46 * Public function prototypes
47 *============================================================================*/
48
49/*----------------------------------------------------------------------------*/
53/*----------------------------------------------------------------------------*/
54
55void
57
58/*----------------------------------------------------------------------------*/
59/*
60 * \brief Wrapper to the function which adds the explicit part of the
61 * convection/diffusion terms of a transport equation of
62 * a scalar field \f$ \varia \f$.
63 *
64 * More precisely, the right hand side \f$ Rhs \f$ is updated as
65 * follows:
66 * \f[
67 * Rhs = Rhs - \sum_{\fij \in \Facei{\celli}} \left(
68 * \dot{m}_\ij \left( \varia_\fij - \varia_\celli \right)
69 * - \mu_\fij \gradv_\fij \varia \cdot \vect{S}_\ij \right)
70 * \f]
71 *
72 * Warning:
73 * - \f$ Rhs \f$ has already been initialized
74 * before calling cs_balance_scalar!
75 * - mind the minus sign
76 *
77 * Options for the convective scheme:
78 * - blencp = 0: upwind scheme for the advection
79 * - blencp = 1: no upwind scheme except in the slope test
80 * - ischcp = 0: SOLU
81 * - ischcp = 1: centered
82 * - ischcp = 2: SOLU with upwind gradient reconstruction
83 * - ischcp = 3: blending SOLU centered
84 * - ischcp = 4: NVD-TVD
85 * - imucpp = 0: do not multiply the convective part by \f$ C_p \f$
86 * - imucpp = 1: multiply the convective part by \f$ C_p \f$
87 *
88 * \param[in] idtvar indicator of the temporal scheme
89 * \param[in] f_id field id (or -1)
90 * \param[in] imucpp indicator
91 * - 0 do not multiply the convectiv term by Cp
92 * - 1 do multiply the convectiv term by Cp
93 * \param[in] imasac take mass accumulation into account?
94 * \param[in] inc indicator
95 * - 0 when solving an increment
96 * - 1 otherwise
97 * \param[in] eqp pointer to a cs_equation_param_t structure which
98 * contains variable calculation options
99 * \param[in] pvar solved variable (current time step)
100 * may be NULL if pvara != NULL
101 * \param[in] pvara solved variable (previous time step)
102 * may be NULL if pvar != NULL
103 * \param[in] bc_coeffs boundary condition structure for the variable
104 * \param[in] i_massflux mass flux at interior faces
105 * \param[in] b_massflux mass flux at boundary faces
106 * \param[in] i_visc \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
107 * at interior faces for the r.h.s.
108 * \param[in] b_visc \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
109 * at boundary faces for the r.h.s.
110 * \param[in] viscel symmetric cell tensor \f$ \tens{\mu}_\celli \f$
111 * \param[in] xcpp array of specific heat (Cp)
112 * \param[in] weighf internal face weight between cells i j in case
113 * of tensor diffusion
114 * \param[in] weighb boundary face weight for cells i in case
115 * of tensor diffusion
116 * \param[in] icvflb global indicator of boundary convection flux
117 * - 0 upwind scheme at all boundary faces
118 * - 1 imposed flux at some boundary faces
119 * \param[in] icvfli boundary face indicator array of convection flux
120 * - 0 upwind scheme
121 * - 1 imposed flux
122 * \param[in,out] smbrp right hand side \f$ \vect{Rhs} \f$
123 */
124/*----------------------------------------------------------------------------*/
125
126void
127cs_balance_scalar(int idtvar,
128 int f_id,
129 int imucpp,
130 int imasac,
131 int inc,
133 cs_real_t pvar[],
134 const cs_real_t pvara[],
135 const cs_field_bc_coeffs_t *bc_coeffs,
136 const cs_real_t i_massflux[],
137 const cs_real_t b_massflux[],
138 const cs_real_t i_visc[],
139 const cs_real_t b_visc[],
140 cs_real_6_t viscel[],
141 const cs_real_t xcpp[],
142 const cs_real_2_t weighf[],
143 const cs_real_t weighb[],
144 int icvflb,
145 const int icvfli[],
146 cs_real_t smbrp[]);
147
149
150#if defined(__cplusplus)
151
152/*----------------------------------------------------------------------------*/
153/*
154 * \brief Wrapper to the function which adds the explicit part of the
155 * convection/diffusion terms of a transport equation of
156 * a scalar field \f$ \varia \f$.
157 *
158 * More precisely, the right hand side \f$ Rhs \f$ is updated as
159 * follows:
160 * \f[
161 * Rhs = Rhs - \sum_{\fij \in \Facei{\celli}} \left(
162 * \dot{m}_\ij \left( \varia_\fij - \varia_\celli \right)
163 * - \mu_\fij \gradv_\fij \varia \cdot \vect{S}_\ij \right)
164 * \f]
165 *
166 * Warning:
167 * - \f$ Rhs \f$ has already been initialized
168 * before calling cs_balance_scalar!
169 * - mind the minus sign
170 *
171 * Options for the convective scheme:
172 * - blencp = 0: upwind scheme for the advection
173 * - blencp = 1: no upwind scheme except in the slope test
174 * - ischcp = 0: SOLU
175 * - ischcp = 1: centered
176 * - ischcp = 2: SOLU with upwind gradient reconstruction
177 * - ischcp = 3: blending SOLU centered
178 * - ischcp = 4: NVD-TVD
179 * - imucpp = 0: do not multiply the convective part by \f$ C_p \f$
180 * - imucpp = 1: multiply the convective part by \f$ C_p \f$
181 *
182 * \param[in] idtvar indicator of the temporal scheme
183 * \param[in] f_id field id (or -1)
184 * \param[in] imucpp indicator
185 * - 0 do not multiply the convectiv term by Cp
186 * - 1 do multiply the convectiv term by Cp
187 * \param[in] imasac take mass accumulation into account?
188 * \param[in] inc indicator
189 * - 0 when solving an increment
190 * - 1 otherwise
191 * \param[in] eqp pointer to a cs_equation_param_t structure which
192 * contains variable calculation options
193 * \param[in] pvar solved variable (current time step)
194 * may be NULL if pvara != NULL
195 * \param[in] pvara solved variable (previous time step)
196 * may be NULL if pvar != NULL
197 * \param[in] bc_coeffs boundary condition structure for the variable
198 * \param[in] i_massflux mass flux at interior faces
199 * \param[in] b_massflux mass flux at boundary faces
200 * \param[in] i_visc \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
201 * at interior faces for the r.h.s.
202 * \param[in] b_visc \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
203 * at boundary faces for the r.h.s.
204 * \param[in] viscel symmetric cell tensor \f$ \tens{\mu}_\celli \f$
205 * \param[in] xcpp array of specific heat (Cp)
206 * \param[in] weighf internal face weight between cells i j in case
207 * of tensor diffusion
208 * \param[in] weighb boundary face weight for cells i in case
209 * of tensor diffusion
210 * \param[in] icvflb global indicator of boundary convection flux
211 * - 0 upwind scheme at all boundary faces
212 * - 1 imposed flux at some boundary faces
213 * \param[in] icvfli boundary face indicator array of convection flux
214 * - 0 upwind scheme
215 * - 1 imposed flux
216 * \param[in,out] smbrp right hand side \f$ \vect{Rhs} \f$
217 * \param[in,out] i_flux interior flux (or nullptr)
218 * \param[in,out] b_flux boundary flux (or nullptr)
219 */
220/*----------------------------------------------------------------------------*/
221
222void
223cs_balance_scalar(int idtvar,
224 int f_id,
225 int imucpp,
226 int imasac,
227 int inc,
229 cs_real_t pvar[],
230 const cs_real_t pvara[],
231 const cs_field_bc_coeffs_t *bc_coeffs,
232 const cs_real_t i_massflux[],
233 const cs_real_t b_massflux[],
234 const cs_real_t i_visc[],
235 const cs_real_t b_visc[],
236 cs_real_6_t viscel[],
237 const cs_real_t xcpp[],
238 const cs_real_2_t weighf[],
239 const cs_real_t weighb[],
240 int icvflb,
241 const int icvfli[],
242 cs_real_t smbrp[],
243 cs_real_2_t i_flux[],
244 cs_real_t b_flux[]);
245
246#endif // defined(__cplusplus)
247
249
250/*----------------------------------------------------------------------------*/
251/*
252 * \brief Wrapper to the function which adds the explicit part of the
253 * convection/diffusion
254 * terms of a transport equation of a vector field \f$ \vect{\varia} \f$.
255 *
256 * More precisely, the right hand side \f$ \vect{Rhs} \f$ is updated as
257 * follows:
258 * \f[
259 * \vect{Rhs} = \vect{Rhs} - \sum_{\fij \in \Facei{\celli}} \left(
260 * \dot{m}_\ij \left( \vect{\varia}_\fij - \vect{\varia}_\celli \right)
261 * - \mu_\fij \gradt_\fij \vect{\varia} \cdot \vect{S}_\ij \right)
262 * \f]
263 *
264 * Remark:
265 * if ivisep = 1, then we also take \f$ \mu \transpose{\gradt\vect{\varia}}
266 * + \lambda \trace{\gradt\vect{\varia}} \f$, where \f$ \lambda \f$ is
267 * the secondary viscosity, i.e. usually \f$ -\frac{2}{3} \mu \f$.
268 *
269 * Warning:
270 * - \f$ \vect{Rhs} \f$ has already been initialized
271 * before calling cs_balance_vector!
272 * - mind the sign minus
273 *
274 * Options for the convective scheme:
275 * - blencp = 0: upwind scheme for the advection
276 * - blencp = 1: no upwind scheme except in the slope test
277 * - ischcp = 0: SOLU
278 * - ischcp = 1: centered
279 * - ischcp = 2: SOLU with upwind gradient reconstruction
280 * - ischcp = 3: blending SOLU centered
281 * - ischcp = 4: NVD-TVD
282 *
283 * \param[in] idtvar indicator of the temporal scheme
284 * \param[in] f_id field id (or -1)
285 * \param[in] imasac take mass accumulation into account?
286 * \param[in] inc indicator
287 * - 0 when solving an increment
288 * - 1 otherwise
289 * \param[in] ivisep indicator to take \f$ \divv
290 * \left(\mu \gradt \transpose{\vect{a}} \right)
291 * -2/3 \grad\left( \mu \dive \vect{a} \right)\f$
292 * - 1 take into account,
293 * - 0 otherwise
294 * \param[in] eqp pointer to a cs_equation_param_t structure which
295 * contains variable calculation options
296 * \param[in] pvar solved velocity (current time step)
297 * \param[in] pvara solved velocity (previous time step)
298 * \param[in] bc_coeffs_v boundary condition structure for the variable
299 * \param[in] bc_coeffs_solve_v sweep loop boundary conditions structure
300 * \param[in] i_massflux mass flux at interior faces
301 * \param[in] b_massflux mass flux at boundary faces
302 * \param[in] i_visc \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
303 * at interior faces for the r.h.s.
304 * \param[in] b_visc \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
305 * at boundary faces for the r.h.s.
306 * \param[in] secvif secondary viscosity at interior faces
307 * \param[in] secvib secondary viscosity at boundary faces
308 * \param[in] viscel symmetric cell tensor \f$ \tens{\mu}_\celli \f$
309 * \param[in] weighf internal face weight between cells i j in case
310 * of tensor diffusion
311 * \param[in] weighb boundary face weight for cells i in case
312 * of tensor diffusion
313 * \param[in] icvflb global indicator of boundary convection flux
314 * - 0 upwind scheme at all boundary faces
315 * - 1 imposed flux at some boundary faces
316 * \param[in] icvfli boundary face indicator array of convection flux
317 * - 0 upwind scheme
318 * - 1 imposed flux
319 * \param[in,out] smbr right hand side \f$ \vect{Rhs} \f$
320 */
321/*----------------------------------------------------------------------------*/
322
323void
324cs_balance_vector(int idtvar,
325 int f_id,
326 int imasac,
327 int inc,
328 int ivisep,
330 cs_real_3_t pvar[],
331 const cs_real_3_t pvara[],
332 const cs_field_bc_coeffs_t *bc_coeffs_v,
333 const cs_bc_coeffs_solve_t *bc_coeffs_solve_v,
334 const cs_real_t i_massflux[],
335 const cs_real_t b_massflux[],
336 const cs_real_t i_visc[],
337 const cs_real_t b_visc[],
338 const cs_real_t secvif[],
339 const cs_real_t secvib[],
340 cs_real_6_t viscel[],
341 const cs_real_2_t weighf[],
342 const cs_real_t weighb[],
343 int icvflb,
344 const int icvfli[],
345 cs_real_3_t i_pvar[],
346 cs_real_3_t b_pvar[],
347 cs_real_3_t smbr[]);
348
349/*----------------------------------------------------------------------------*/
350/*
351 * \brief Wrapper to the function which adds the explicit part of the
352 * convection/diffusion
353 * terms of a transport equation of a tensor field \f$ \tens{\varia} \f$.
354 *
355 * More precisely, the right hand side \f$ \vect{Rhs} \f$ is updated as
356 * follows:
357 * \f[
358 * \tens{Rhs} = \tens{Rhs} - \sum_{\fij \in \Facei{\celli}} \left(
359 * \dot{m}_\ij \left( \tens{\varia}_\fij - \tens{\varia}_\celli \right)
360 * - \mu_\fij \gradt_\fij \tens{\varia} \cdot \tens{S}_\ij \right)
361 * \f]
362 *
363 * Warning:
364 * - \f$ \tens{Rhs} \f$ has already been initialized before calling bilscts!
365 * - mind the sign minus
366 *
367 * Options for the convective scheme:
368 * - blencp = 0: upwind scheme for the advection
369 * - blencp = 1: no upwind scheme except in the slope test
370 * - ischcp = 0: SOLU
371 * - ischcp = 1: centered
372 * - ischcp = 2: SOLU with upwind gradient reconstruction
373 * - ischcp = 3: blending SOLU centered
374 * - ischcp = 4: NVD-TVD
375 *
376 * \param[in] idtvar indicator of the temporal scheme
377 * \param[in] f_id field id (or -1)
378 * \param[in] imasac take mass accumulation into account?
379 * \param[in] inc indicator
380 * \param[in] eqp pointer to a cs_equation_param_t structure which
381 * contains variable calculation options
382 * \param[in] pvar solved velocity (current time step)
383 * \param[in] pvara solved velocity (previous time step)
384 * \param[in] bc_coeffs_ts boundary condition structure for the variable
385 * \param[in] bc_coeffs_solve_ts sweep loop boundary conditions structure
386 * \param[in] i_massflux mass flux at interior faces
387 * \param[in] b_massflux mass flux at boundary faces
388 * \param[in] i_visc \f$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} \f$
389 * at interior faces for the r.h.s.
390 * \param[in] b_visc \f$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} \f$
391 * at boundary faces for the r.h.s.
392 * \param[in] viscel symmetric cell tensor \f$ \tens{\mu}_\celli \f$
393 * \param[in] weighf internal face weight between cells i j in case
394 * of tensor diffusion
395 * \param[in] weighb boundary face weight for cells i in case
396 * of tensor diffusion
397 * \param[in] icvflb global indicator of boundary convection flux
398 * - 0 upwind scheme at all boundary faces
399 * - 1 imposed flux at some boundary faces
400 * \param[in] icvfli boundary face indicator array of convection flux
401 * - 0 upwind scheme
402 * - 1 imposed flux
403 * \param[in,out] smbrp right hand side \f$ \vect{Rhs} \f$
404 */
405/*----------------------------------------------------------------------------*/
406
407void
408cs_balance_tensor(int idtvar,
409 int f_id,
410 int imasac,
411 int inc,
413 cs_real_6_t pvar[],
414 const cs_real_6_t pvara[],
415 const cs_field_bc_coeffs_t *bc_coeffs_ts,
416 const cs_bc_coeffs_solve_t *bc_coeffs_solve_ts,
417 const cs_real_t i_massflux[],
418 const cs_real_t b_massflux[],
419 const cs_real_t i_visc[],
420 const cs_real_t b_visc[],
421 cs_real_6_t viscel[],
422 const cs_real_2_t weighf[],
423 const cs_real_t weighb[],
424 int icvflb,
425 const int icvfli[],
426 cs_real_6_t smbrp[]);
427
428/*----------------------------------------------------------------------------*/
429
431
432#endif /* __CS_BALANCE_H__ */
void cs_balance_initialize(void)
Initialize balance timers.
Definition cs_balance.cpp:116
void cs_balance_vector(int idtvar, int f_id, int imasac, int inc, int ivisep, cs_equation_param_t *eqp, cs_real_3_t pvar[], const cs_real_3_t pvara[], const cs_field_bc_coeffs_t *bc_coeffs_v, const cs_bc_coeffs_solve_t *bc_coeffs_solve_v, const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t secvif[], const cs_real_t secvib[], cs_real_6_t viscel[], const cs_real_2_t weighf[], const cs_real_t weighb[], int icvflb, const int icvfli[], cs_real_3_t i_pvar[], cs_real_3_t b_pvar[], cs_real_3_t smbr[])
Wrapper to the function which adds the explicit part of the convection/diffusion terms of a transport...
Definition cs_balance.cpp:551
void cs_balance_tensor(int idtvar, int f_id, int imasac, int inc, cs_equation_param_t *eqp, cs_real_6_t pvar[], const cs_real_6_t pvara[], const cs_field_bc_coeffs_t *bc_coeffs_ts, const cs_bc_coeffs_solve_t *bc_coeffs_solve_ts, 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 viscel[], const cs_real_2_t weighf[], const cs_real_t weighb[], int icvflb, const int icvfli[], cs_real_6_t smbrp[])
Wrapper to the function which adds the explicit part of the convection/diffusion terms of a transport...
Definition cs_balance.cpp:776
void cs_balance_scalar(int idtvar, int f_id, int imucpp, int imasac, int inc, cs_equation_param_t *eqp, cs_real_t pvar[], const cs_real_t pvara[], const cs_field_bc_coeffs_t *bc_coeffs, 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 viscel[], const cs_real_t xcpp[], const cs_real_2_t weighf[], const cs_real_t weighb[], int icvflb, const int icvfli[], cs_real_t smbrp[])
Wrapper to the function which adds the explicit part of the convection/diffusion terms of a transport...
Definition cs_balance.cpp:196
#define BEGIN_C_DECLS
Definition cs_defs.h:542
double cs_real_t
Floating-point value.
Definition cs_defs.h:342
cs_real_t cs_real_2_t[2]
vector of 2 floating-point values
Definition cs_defs.h:358
cs_real_t cs_real_6_t[6]
vector of 6 floating-point values
Definition cs_defs.h:361
#define END_C_DECLS
Definition cs_defs.h:543
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition cs_defs.h:359
Definition cs_field.h:104
Set of parameters to handle an unsteady convection-diffusion-reaction equation with term sources.
Definition cs_equation_param.h:193
Field boundary condition descriptor (for variables)
Definition cs_field.h:121