1 #ifndef __CS_EQUATION_PARAM_H__
2 #define __CS_EQUATION_PARAM_H__
3
4 /*============================================================================
5 * Functions related to the structure cs_equation_param_t storing the settings
6 * related to an equation.
7 *============================================================================*/
8
9 /*
10 This file is part of Code_Saturne, a general-purpose CFD tool.
11
12 Copyright (C) 1998-2021 EDF S.A.
13
14 This program is free software; you can redistribute it and/or modify it under
15 the terms of the GNU General Public License as published by the Free Software
16 Foundation; either version 2 of the License, or (at your option) any later
17 version.
18
19 This program is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
21 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
22 details.
23
24 You should have received a copy of the GNU General Public License along with
25 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
26 Street, Fifth Floor, Boston, MA 02110-1301, USA.
27 */
28
29 /*----------------------------------------------------------------------------
30 * Local headers
31 *----------------------------------------------------------------------------*/
32
33 #include "cs_advection_field.h"
34 #include "cs_enforcement.h"
35 #include "cs_hodge.h"
36 #include "cs_param_cdo.h"
37 #include "cs_param_sles.h"
38 #include "cs_property.h"
39 #include "cs_xdef.h"
40
41 /*----------------------------------------------------------------------------*/
42
43 BEGIN_C_DECLS
44
45 /*============================================================================
46 * Macro definitions
47 *============================================================================*/
48
49 /*!
50 * @name Flags specifying which term is needed for an equation.
51 * @{
52 *
53 * \def CS_EQUATION_LOCKED
54 * \brief Parameters for setting an equation are not modifiable anymore
55 *
56 * \def CS_EQUATION_UNSTEADY
57 * \brief Unsteady term is needed
58 *
59 * \def CS_EQUATION_CONVECTION
60 * \brief Convection term is needed
61 *
62 * \def CS_EQUATION_DIFFUSION
63 * \brief Diffusion term is needed. A scalar-/vector-valued Laplacian
64 * with div .grad
65 *
66 * \def CS_EQUATION_CURLCURL
67 * \brief The term related to the curl-curl operator is needed
68 *
69 * \def CS_EQUATION_GRADDIV
70 * \brief The term related to the grad-div operator is needed
71 *
72 * \def CS_EQUATION_REACTION
73 * \brief Reaction term is needed
74 *
75 * \def CS_EQUATION_FORCE_VALUES
76 * \brief Add an algebraic manipulation to set the value of a given set of
77 * interior degrees of freedom
78 *
79 * \def CS_EQUATION_USER_HOOK
80 * \brief Activate a user hook to get a fine control of the discretization
81 * process during the cellwise building of the linear system
82 * Need to match the cs_equation_user_hook_t prototype
83 *
84 */
85
86 #define CS_EQUATION_LOCKED (1 << 0) /* 1 */
87 #define CS_EQUATION_UNSTEADY (1 << 1) /* 2 */
88 #define CS_EQUATION_CONVECTION (1 << 2) /* 4 */
89 #define CS_EQUATION_DIFFUSION (1 << 3) /* 8 */
90 #define CS_EQUATION_CURLCURL (1 << 4) /* 16 */
91 #define CS_EQUATION_GRADDIV (1 << 5) /* 32 */
92 #define CS_EQUATION_REACTION (1 << 6) /* 64 */
93 #define CS_EQUATION_FORCE_VALUES (1 << 7) /* 128 */
94 #define CS_EQUATION_USER_HOOK (1 << 8) /* 256 */
95 #define CS_EQUATION_USER_TRIGGERED (1 << 9) /* 512 */
96
97
98 /*!
99 * @}
100 * @name Flags specifying which extra operation is needed for an equation.
101 * @{
102 *
103 * \def CS_EQUATION_POST_BALANCE
104 * \brief Compute and postprocess the equation balance
105 *
106 * \def CS_EQUATION_POST_PECLET
107 * \brief Compute and postprocess the Peclet number
108 *
109 * \def CS_EQUATION_POST_UPWIND_COEF
110 * \brief Postprocess the value of the upwinding coefficient
111 *
112 * \def CS_EQUATION_POST_NORMAL_FLUX
113 * \brief Postprocess the value of the normal flux across the boundary faces
114 *
115 */
116
117 #define CS_EQUATION_POST_BALANCE (1 << 0) /* 1 */
118 #define CS_EQUATION_POST_PECLET (1 << 1) /* 2 */
119 #define CS_EQUATION_POST_UPWIND_COEF (1 << 2) /* 4 */
120 #define CS_EQUATION_POST_NORMAL_FLUX (1 << 3) /* 8 */
121
122 /*! @} */
123
124 /*============================================================================
125 * Type definitions
126 *============================================================================*/
127
128 /*! \enum cs_equation_type_t
129 * \brief Type of equations managed by the solver
130 *
131 * \var CS_EQUATION_TYPE_GROUNDWATER
132 * Equation related to the groundwater flow module
133 *
134 * \var CS_EQUATION_TYPE_MAXWELL
135 * Equation related to the Maxwell module
136 *
137 * \var CS_EQUATION_TYPE_NAVSTO
138 * Equation related to the resolution of the Navier-Stokes system
139 * - Example: momentum, prediction, correction, energy...
140 *
141 * \var CS_EQUATION_TYPE_PREDEFINED
142 * Predefined equation (most part of the setting is already done)
143 * - Example: equation for the wall distance or ALE
144 *
145 * \var CS_EQUATION_TYPE_THERMAL
146 * Equation related to the heat transfer
147 *
148 * \var CS_EQUATION_TYPE_SOLIDIFICATION
149 * Equation related to the solidification module
150 *
151 * \var CS_EQUATION_TYPE_USER
152 * User-defined equation (no interaction with energy and momentum equations)
153 * Resolution is performed after all predefined equations have been solved.
154 */
155
156 typedef enum {
157
158 CS_EQUATION_TYPE_GROUNDWATER,
159 CS_EQUATION_TYPE_MAXWELL,
160 CS_EQUATION_TYPE_NAVSTO,
161 CS_EQUATION_TYPE_PREDEFINED,
162 CS_EQUATION_TYPE_THERMAL,
163 CS_EQUATION_TYPE_SOLIDIFICATION,
164 CS_EQUATION_TYPE_USER,
165
166 CS_EQUATION_N_TYPES
167
168 } cs_equation_type_t;
169
170 /*----------------------------------------------------------------------------*/
171
172 /*! \struct cs_equation_param_t
173 * \brief Set of parameters to handle an unsteady convection-diffusion-reaction
174 * equation with term sources
175 */
176
177 typedef struct {
178
179 /*!
180 * @name General settings
181 * @{
182 */
183
184 char *restrict name; /*!< name of the equation */
185 cs_equation_type_t type; /*!< type of equation: predefined... */
186 int dim; /*!< Dimension of the unknown */
187
188 /*! \var verbosity
189 * Verbosity for the resolution (0 or 1 for a reasonable log size, 2 or more
190 * for troubleshooting).
191 *
192 * \var iwarni
193 * \deprecated use verbosity instead (iwarni is an alias to verbosity)
194 */
195
196 union {
197 int verbosity;
198 int iwarni;
199 };
200
201 /*! \var flag
202 * Flag to know if unsteady or diffusion or convection or reaction
203 * or source terms are activated or not
204 */
205
206 cs_flag_t flag;
207
208 /*! \var process_flag
209 * Flag to determine if predefined post-treatments such as Peclet,
210 * are requested
211 */
212
213 cs_flag_t process_flag;
214
215 /* Numerical settings */
216
217 cs_param_space_scheme_t space_scheme; /*!< Space discretization scheme */
218 cs_param_dof_reduction_t dof_reduction; /*!< How is defined DoF */
219
220 /*! \var space_poly_degree
221 * Maximum degree of the polynomial basis
222 */
223
224 int space_poly_degree;
225
226 /*!
227 * @}
228 * @name Legacy Settings
229 * @{
230 *
231 * \var iconv
232 * Indicate if the convection is taken into account (1) or not (0). By default,
233 * 0 for the pressure or f in v2f model, 1 for the other unknowns.
234 *
235 * \var istat
236 * Indicate whether unsteady terms are present (1) or not (0) in the matrices.
237 * By default, 0 for the pressure or f in v2f model, 1 for the other unknowns.
238 *
239 * \var idircl
240 * Indicate whether the diagonal of the matrix should be slightly shifted if
241 * there is no Dirichlet boundary condition and if \ref istat = 0
242 * (0: false / 1: true). Indeed, in such a case, the matrix for the general
243 * advection/diffusion equation is singular. A slight shift in the diagonal
244 * will make it invertible again.
245 * By default, \ref idircl is set to 1 for all the unknowns, except
246 * \f$\overline{f}\f$ in the v2f model (whose equation already contain another
247 * diagonal term).
248 * \remark this code is defined automatically based on the
249 * presence of Dirichlet BCs.
250 *
251 * \var ndircl
252 * Number of Dirichlet BCs
253 *
254 * \var idiff
255 * Indifcate if diffusion is taken into account (1) or not (0).
256 *
257 * \var idifft
258 * When diffusion is taken into account (\ref idiff = 1), indicate if the
259 * turbulent diffusion is taken into account (\ref idifft = 1) or not (0).
260 *
261 * \var idften
262 * Type of diffusivity flag (sum of mask constants defining if diffusivity is
263 * isotropic, anisotropic, ... Masks are defined in \ref scalar_params).
264 *
265 * \var iswdyn
266 * Dynamic relaxation type:
267 * - 0 no dynamic relaxation
268 * - 1 dynamic relaxation depending on \f$ \delta \varia^k \f$
269 * - 2 dynamic relaxation depending on \f$ \delta \varia^k \f$
270 * and \f$ \delta \varia^{k-1} \f$.
271 *
272 * \var ischcv
273 * Indicate the type of second-order convective scheme
274 * - 0: Second Order Linear Upwind
275 * - 1: Centered
276 * - 2: Second Order with upwind-gradient reconstruction (SOLU)
277 * - 3: Blending between SOLU and Centered scheme
278 * - 4: NVD/TVD Scheme
279 * Then "limiter_choice" keyword must be set:
280 * - 0: Gamma
281 * - 1: SMART
282 * - 2: CUBISTA
283 * - 3: SUPERBEE
284 * - 4: MUSCL
285 * - 5: MINMOD
286 * - 6: CLAM
287 * - 7: STOIC
288 * - 8: OSHER
289 * - 9: WASEB
290 * --- VOF scheme ---
291 * - 10: M-HRIC
292 * - 11: M-CICSAM
293 *
294 * \var ibdtso
295 * Backward differential scheme in time order.
296 *
297 * \var isstpc
298 * Indicate whether a slope test should be used to switch from a second-order
299 * to an upwind convective scheme under certain conditions, to ensure
300 * stability.
301 * - 0: slope test activated
302 * - 1: slope test deactivated
303 * - 2: continuous limiter ensuring boundedness (beta limiter)
304 * The use of the slope test stabilises the calculation but may reduce the
305 * spatial convergence order.
306 *
307 * \var nswrgr
308 * Iteration limit for the iterative gradient reconstruction (\ref imrgra = 0).
309 * If \ref imrgra = 0 and \ref nswrgr <= 1, gradients are not reconstructed.
310 *
311 * \var nswrsm
312 * Iteration limit for the reconstruction of the right-hand sides of the
313 * equations with a first-order scheme in time (standard case), the default
314 * values are 2 for pressure and 1 for the other variables. With a
315 * second-order scheme in time (\ref optcal::ischtp "ischtp" = 2) or LES, the
316 * default values are 5 for pressure and 10 for the other variables.
317 *
318 * \var imvisf
319 * Face viscosity field interpolation
320 * - 0: aritmetic (default)
321 * - 1: harmonic
322 *
323 * \var imrgra
324 * Type of gradient reconstruction
325 * - 0: iterative reconstruction of the non-orthogonalities
326 * - 1: least squares method based on the first neighbor cells
327 * (those which share a face with the treated cell)
328 * - 2, 3: least squares method using the extended neighborhood
329 * - 4: Green-Gauss based using the least squares method
330 * (first neighbors) to compute face values
331 * - 5, 6: Green-Gauss based using the least squares method with
332 * and extended neighborhood to compute face values\n
333 * If the computation fails due to mesh quality aspects,
334 * it is usually effective to use \ref imrgra = 3, 5, or 6.
335 *
336 * \var imligr
337 * Type of gradient limiter
338 * - -1 (CS_GRADIENT_LIMIT_NONE): no limitation
339 * - 0 (CS_GRADIENT_LIMIT_CELL): based on the neighbors
340 * - 1 (CS_GRADIENT_LIMIT_FACE): superior order\n
341 * \ref imligr is applied only to least-squares gradients.
342 * In the case of the Green-Gauss gradient with least-squares
343 * based face gradients, applied to the least-squares step.
344 *
345 * \var ircflu
346 * Indicate whether the convective and diffusive fluxes at the faces should be
347 * reconstructed:
348 * - 0: no reconstruction
349 * - 1: reconstruction
350 * Deactivating the reconstruction of the fluxes can have a stabilizing
351 * effect on the calculation. It is sometimes useful with the \f$ k-\epsilon
352 * \f$ model, if the mesh is strongly non-orthogonal in the near-wall region,
353 * where the gradients of k and \f$ \epsilon \f$ are strong. In such a case,
354 * setting \ref ircflu = 0 will probably help (switching to a first order
355 * convective scheme, \ref blencv = 0, for k and \f$ \epsilon \f$ might also
356 * help in that case).
357 *
358 * \var iwgrec
359 * Gradient calculation weighting
360 * - 0: standard
361 * - 1: weighted
362 *
363 * \var icoupl
364 * Internal coupling indicator
365 * - -1: not coupled (default)
366 * - 1: coupled
367 *
368 * \var thetav
369 * Value of \f$ \theta \f$ used to express at the second order the terms of
370 * convection, diffusion and the source terms which are linear functions of
371 * the solved variable, according to the formula
372 * \f$ \phi^{n+\theta} = (1-\theta) \phi^n + \theta \phi^{n+1} \f$.
373 * Generally, only the values 1 and 0.5 are used. The user is not allowed to
374 * modify this variable.
375 * - 1: first-order
376 * - 0.5: second-order \n
377 * For the pressure, \ref thetav is always 1. For the other variables,
378 * \ref thetav = 0.5 is used when the second-order time scheme is activated
379 * (\ref optcal::ischtp "ischtp = 2", standard for LES calculations),
380 * otherwise \ref thetav = 1.
381 *
382 * \var blencv
383 * Proportion of second-order convective scheme (0 corresponds to an upwind
384 * first-order scheme); in case of LES calculation, a second-order scheme is
385 * recommended and activated by default (\ref blencv = 1).\n
386 * Relevant where \ref iconv = 1.
387 *
388 * \var blend_st
389 * Proportion of second-order convective scheme (0 corresponds to an upwind
390 * first-order scheme) after the slope test is activated; in case of LES
391 * calculation, a second-order scheme is recommended and activated by
392 * default (\ref blend_st = 1).\n
393 * Relevant where\ref iconv = 1.
394 *
395 * \var epsilo
396 * Relative precision for the solution of the linear system. The default is
397 * \ref epsilo = \f$ 10^-8 \f$ . When there are enough iterations on the
398 * reconstruction of the right-hand side of the equation, the value may be
399 * increased (by default, in case of second-order in time, with \ref nswrsm =
400 * 5 or 10, \ref epsilo is increased to \f$ 10^-5 \f$.
401 *
402 * \var epsrsm
403 * Relative precision on the reconstruction of the right hand-side. The
404 * default is \ref epsrsm = \f$ 10^-8 \f$. When there are not enough
405 * iterations on the reconstruction of the right-hand side of the equation,
406 * the value may be increased (by default, in case of second-order in time,
407 * with \ref nswrsm = 5 or 10, \ref epsrsm is increased to \f$ 10^-5 \f$ ).
408 *
409 * \var epsrgr
410 * Relative precision for the iterative gradient reconstruction.
411 * (when \ref imrgra = 0).
412 *
413 * \var climgr
414 * For least squares gradients, factor of gradient limitation
415 * (high value means little limitation).\n
416 * Relevant for all the variables using least-squares gradientsfor which
417 * \ref imligr > CS_GRADIENT_LIMIT_NONE.
418 *
419 * \var relaxv
420 * Relaxation coefficient for the associated variable. This relaxation
421 * parameter is only useful for the pressure with the unsteady algorithm (so
422 * as to improve the convergence in case of meshes of insufficient quality or
423 * of some turbulent models (k-epsilon, v2f, k-omega) and
424 * \ref cs_turb_rans_model_t::ikecou "ikecou" = 0; if
425 * \ref cs_turb_rans_model_t::ikecou "ikecou" = 1, \ref relaxv is ignored.\n
426 * Default values are 0.7 for turbulent variables and 1. for pressure.
427 * \ref relaxv also stores the value of the relaxation coefficient when using
428 * the steady algorithm, deduced from the value of
429 * \ref cs_time_step_options_t::relxst "relxst" (defaulting to
430 * \ref relaxv = 1. - relxst).\n
431 * Used only for the pressure and for turbulent variables
432 * (\f$ k-\epsilon \f$, v2f or \f$ k-\omega \f$ models without coupling) with
433 * the unsteady algorithm. Always used with the steady algorithm.
434 *
435 */
436
437 int iconv;
438 int istat;
439 int idircl;
440 int ndircl;
441 int idiff;
442 int idifft;
443 int idften;
444 int iswdyn;
445 int ischcv;
446 int ibdtso;
447 int isstpc;
448 int nswrgr;
449 int nswrsm;
450 int imvisf;
451 int imrgra;
452 int imligr;
453 int ircflu;
454 int iwgrec;
455 int icoupl;
456 double thetav; /* TODO: merge with theta */
457 double blencv;
458 double blend_st;
459 double epsilo;
460 double epsrsm;
461 double epsrgr;
462 double climgr;
463 double relaxv;
464
465 /*!
466 * @}
467 * @name Settings for the boundary conditions
468 * @{
469 *
470 * \var default_bc
471 * Default boundary condition related to this equation. Valid choices:
472 * - \ref CS_PARAM_BC_HMG_NEUMANN
473 * - \ref CS_PARAM_BC_HMG_DIRICHLET
474 *
475 * \var n_bc_defs
476 * Number of boundary conditions which are defined for this equation
477 *
478 * \var bc_defs
479 * Pointers to the definitions of the boundary conditions
480 *
481 * \var default_enforcement
482 * Type of strategy to enforce an essential boundary conditions (Dirichlet for
483 * instance) when no predefined strategy is prescribed. See \ref
484 * cs_param_bc_enforce_t for more details.
485 *
486 * \var strong_pena_bc_coeff
487 * Value of the penalization coefficient used to enforce the Dirichlet
488 * boundary conditions when \ref CS_PARAM_BC_ENFORCE_PENALIZED is set. This
489 * value should be sufficiently large in order to neglect off-diagonal terms.
490 *
491 * \var weak_pena_bc_coeff
492 * Value of the penalization coefficient used to enforce the Dirichlet
493 * boundary condition when \ref CS_PARAM_BC_ENFORCE_WEAK_NITSCHE or \ref
494 * CS_PARAM_BC_ENFORCE_WEAK_SYM is set. This two latter strategies have a
495 * lesser influence on the conditioning number of the linear system than the
496 * choice \ref CS_PARAM_BC_ENFORCE_PENALIZED
497 *
498 */
499
500 cs_param_bc_type_t default_bc;
501 int n_bc_defs;
502 cs_xdef_t **bc_defs;
503
504 cs_param_bc_enforce_t default_enforcement;
505 cs_real_t strong_pena_bc_coeff;
506 cs_real_t weak_pena_bc_coeff;
507
508 /*!
509 * @}
510 * @name Numerical settings for the time-dependent parameters
511 * @{
512 *
513 * \var n_ic_defs
514 * Number of definitions for setting the intial condition
515 *
516 * \var ic_defs
517 * List of pointers to the definition of the inititial condition
518 */
519
520 int n_ic_defs;
521 cs_xdef_t **ic_defs;
522
523 /*!
524 * @}
525 * @name Numerical settings for the time-dependent parameters
526 * @{
527 *
528 * \var do_lumping
529 * Set to true or false. Activate several actions:
530 * Perform a mass lumping of the matrices related to the time and reaction
531 * discretization. All source terms are evaluated using a barycentric
532 * quadrature.
533 */
534
535 bool do_lumping;
536
537 /*!
538 * \var time_hodgep
539 * Set of parameters for the discrete Hodge operator related to the unsteady
540 * term
541 *
542 * \var time_property
543 * Pointer to the \ref cs_property_t structure related to the unsteady term
544 *
545 * \var time_scheme
546 * Type of numerical scheme used for the time discretization
547 *
548 * \var theta
549 * Value of the coefficient for a theta scheme (between 0 and 1)
550 *
551 */
552
553 cs_hodge_param_t time_hodgep;
554 cs_property_t *time_property;
555 cs_param_time_scheme_t time_scheme;
556 cs_real_t theta;
557
558 /*!
559 * @}
560 * @name Numerical settings for the diffusion (Laplacian div-grad) term
561 * @{
562 *
563 * \var diffusion_hodgep
564 * Set of parameters for the discrete Hodge operator related to the diffusion
565 *
566 * \var diffusion_property
567 * Pointer to the property related to the diffusion term
568 */
569
570 cs_hodge_param_t diffusion_hodgep;
571 cs_property_t *diffusion_property;
572
573 /*!
574 * @}
575 * @name Numerical settings for the curl-curl term
576 * @{
577 *
578 * \var curlcurl_hodgep
579 * Set of parameters for the discrete Hodge operator related to the
580 * curl-curl operator
581 *
582 * \var curlcurl_property
583 * Pointer to the property related to the curl-curl term
584 */
585
586 cs_hodge_param_t curlcurl_hodgep;
587 cs_property_t *curlcurl_property;
588
589 /*!
590 * @}
591 * @name Numerical settings for the grad-div term
592 * @{
593 *
594 * \var graddiv_hodgep
595 * Set of parameters for the discrete Hodge operator related to the grad-div
596 * operator
597 *
598 * \var graddiv_property
599 * Pointer to the property related to the grad-div term
600 */
601
602 cs_hodge_param_t graddiv_hodgep;
603 cs_property_t *graddiv_property;
604
605 /*!
606 * @}
607 * @name Numerical settings for the advection term
608 * @{
609 *
610 * \var adv_formulation
611 * Type of formulation (conservative, non-conservative...) for the advective
612 * term
613 *
614 * \var adv_scheme
615 * Numerical scheme used for the discretization of the advection term
616 *
617 * \var adv_strategy
618 * Strategy used to handle the advection term (please refer to \ref
619 * cs_param_advection_strategy_t)
620 *
621 * \var adv_extrapol
622 * Extrapolation used to estimate the advection field (please refer to \ref
623 * cs_param_advection_extrapol_t)
624 *
625 * \var upwind_portion
626 * Value between 0. and 1. (0: centered scheme, 1: pure upwind scheme)
627 * Introduce a constant portion of upwinding in a centered scheme
628 * Only useful if the advection scheme is set to
629 * CS_PARAM_ADVECTION_SCHEME_HYBRID_CENTERED_UPWIND
630 *
631 * \var adv_field
632 * Pointer to the \ref cs_adv_field_t structure associated to the advection
633 * term
634 *
635 * \var adv_scaling_property
636 * May be set to NULL even if the advection term is activated. The value of
637 * this property in each cell is multiplicative coefficient in front of the
638 * advection term (boundary terms are also considered) This is useful to treat
639 * the thermal module using the variable temperature instead of the enthalpy
640 * for instance or in the solidification module.
641 */
642
643 cs_param_advection_form_t adv_formulation;
644 cs_param_advection_scheme_t adv_scheme;
645 cs_param_advection_strategy_t adv_strategy;
646 cs_param_advection_extrapol_t adv_extrapol;
647 cs_real_t upwind_portion;
648 cs_adv_field_t *adv_field;
649 cs_property_t *adv_scaling_property;
650
651 /*!
652 * @}
653 * @name Numerical settings for the reaction term
654 * The contribution of a reaction term to the algebraic system is either at
655 * the left-hand and/or right-hand side according to the choice of time
656 * scheme
657 * @{
658 *
659 * \var reaction_hodgep
660 * Set of parameters for the discrete Hodge operator related to the reaction
661 * terms
662 *
663 * \var n_reaction_terms
664 * Number of reaction terms to consider.
665 *
666 * \var reaction_properties
667 * List of properties associated to each reaction term
668 */
669
670 cs_hodge_param_t reaction_hodgep;
671 int n_reaction_terms;
672 cs_property_t **reaction_properties;
673
674 /*!
675 * @}
676 * @name Definition of the related source terms
677 * The contribution of a source term to the algebraic system is always at
678 * right-hand side whatever is the choice of time scheme
679 * @{
680 *
681 * \var n_source_terms
682 * Number of source terms to consider.
683 *
684 * \var source_terms
685 * List of definition of each source term
686 */
687
688 int n_source_terms;
689 cs_xdef_t **source_terms;
690
691 /*!
692 * @}
693 * @name Definition of the related volume mass injection
694 * The contribution of a volume mass injection to the algebraic system
695 * is always at right-hand side whatever is the choice of time scheme,
696 * and is defined in a manner similar to boundary conditions. For
697 * variables whose injection value matches the "ambient" value, no term
698 * needs to be added.
699 * @{
700 *
701 * \var n_volume_mass_injections
702 * Number of volume injections to consider.
703 *
704 * \var volume_mass_injections
705 * List of definitions of injection values
706 */
707
708 int n_volume_mass_injections;
709 cs_xdef_t **volume_mass_injections;
710
711 /*!
712 * @}
713 * @name Enforcement of values inside the computational domain
714 *
715 * This is different from the enforcement of boundary conditions but rely on
716 * the same algebraic manipulation. Only an algebraic manipulation is used to
717 * enforce interior/border DoFs.
718 *
719 * @{
720 *
721 * \var n_enforcements
722 * Number of enforcements which have been specified
723 *
724 * \var enforcement_params
725 * Array of \ref cs_enforcement_param_t structures storing the settings of
726 * each enforcement
727 */
728
729 int n_enforcements;
730 cs_enforcement_param_t **enforcement_params;
731
732 /*!
733 * @}
734 * @name Settings related to the resolution of the algebraic system
735 * @{
736 *
737 * \var sles_param
738 * Set of parameters to specify how to to solve the algebraic system
739 * - iterative solver
740 * - preconditioner
741 * - tolerance...
742 */
743
744 cs_param_sles_t *sles_param;
745
746 /*!
747 * @}
748 * @name Settings related to the optimization of the performance
749 * @{
750 *
751 * \var omp_assembly_choice
752 * When OpenMP is active, choice of parallel reduction for the assembly
753 */
754
755 cs_param_assemble_omp_strategy_t omp_assembly_choice;
756
757 /*! @} */
758
759 } cs_equation_param_t;
760
761 /*----------------------------------------------------------------------------*/
762
763 /*! \enum cs_equation_key_t
764 * \brief List of available keys for setting the parameters of an equation
765 *
766 * \var CS_EQKEY_ADV_EXTRAPOL
767 * Choice in the way to extrapolate the advection field when building the
768 * advection operator
769 * - "none" (default)
770 * - "taylor"
771 * - "adams_bashforth"
772 * Please refer to \ref cs_param_advection_extrapol_t for more details.
773 *
774 * \var CS_EQKEY_ADV_FORMULATION
775 * Kind of formulation of the advective term. Available choices are:
776 * - "conservative"
777 * - "non_conservative"
778 *
779 * \var CS_EQKEY_ADV_SCHEME
780 * Type of numerical scheme for the advective term. The available choices
781 * depend on the space discretization scheme.
782 * - "upwind" (cf. \ref CS_PARAM_ADVECTION_SCHEME_UPWIND)
783 * - "centered" (cf. \ref CS_PARAM_ADVECTION_SCHEME_CENTERED)
784 * - "mix_centered_upwind" or "hybrid_centered_upwind"
785 * (cf. \ref CS_PARAM_ADVECTION_SCHEME_HYBRID_CENTERED_UPWIND)
786 * - "samarskii" --> switch smoothly between an upwind and a centered scheme
787 * thanks to a weight depending on the Peclet number. (cf.
788 * \ref CS_PARAM_ADVECTION_SCHEME_SAMARSKII). Only for CDO-Vb schemes.
789 * - "sg" --> closely related to "samarskii" but with a different definition of
790 * the weight (cf. \ref CS_PARAM_ADVECTION_SCHEME_SG). Only for CDO-Vb schemes
791 * - "cip" --> means "continuous interior penalty" (only for CDOVCB schemes).
792 * Enable a better accuracy. (cf. \ref CS_PARAM_ADVECTION_SCHEME_CIP)
793 * - "cip_cw" --> means "continuous interior penalty" (only for CDOVCB schemes).
794 * Enable a better accuracy.
795 * Consider a cellwise approximation of the advection field.
796 * (cf. \ref CS_PARAM_ADVECTION_SCHEME_CIP_CW)
797 *
798 * \var CS_EQKEY_ADV_STRATEGY
799 * Strategy used to handle the advection term
800 * - "fully_implicit" or "implicit" (default choice)
801 * - "linearized" or "implicit_linear"
802 * - "explicit"
803 * There is a difference between the two first choices when the advection term
804 * induces a non-linearity. In this situation, the first choice implies that a
805 * non-linear algorithm has to be used.
806 *
807 * \var CS_EQKEY_ADV_UPWIND_PORTION
808 * Value between 0 and 1 specifying the portion of upwind added to a centered
809 * discretization.
810 *
811 * \var CS_EQKEY_AMG_TYPE
812 * Specify which type of algebraic multigrid (AMG) to choose.
813 * Available choices are:
814 * - "none" --> (default) No predefined AMG solver
815 * - "boomer"/"boomer_v" --> Boomer AMG V-cycle multigrid from the Hypre library
816 * - "boomer_w" --> Boomer AMG W-cycle multigrid from the Hypre library
817 * - "gamg"/"gamg_v" --> GAMG V-cycle multigrid from the PETSc library
818 * - "gamg_w" --> GAMG W-cycle multigrid from the PETSc library
819 * - "pcmg" --> MG multigrid as preconditioner from the PETSc library
820 * - "v_cycle" --> code_saturne's built-in multigrid with a V-cycle strategy
821 * - "k_cycle" --> code_saturne's built-in multigrid with a K-cycle strategy
822 * WARNING: When selecting "boomer", "boomer_w" and "gamg"/"gamg_w", one needs
823 * to install code_saturne with the PETSc library and PETSc with Hypre if one
824 * wants to use boomeramg
825 *
826 * \var CS_EQKEY_BC_ENFORCEMENT
827 * Set the type of enforcement of the boundary conditions.
828 * Available choices are:
829 * - "algebraic": Modify the linear system so as to add the contribution of the
830 * Dirichlet in the right-hand side and replace the line associated to a
831 * Dirichlet BC by identity. This is a good choice for pure diffusinon or pure
832 * convection problem.
833 * - "penalization": Add a huge penalization coefficient on the diagonal term
834 * of the line related to DoFs associated a Dirichlet BC. The right-hand side is
835 * also modified to take into account this penalization. Be aware that it may
836 * worsen the matrix conditioning.
837 * - "weak": weak enforcement using the Nitsche method (there is also
838 * penalization term but the scaling is such that a moderate value (1-100) of
839 * the penalization coefficient is sufficient). This a good choice for
840 * convection/diffusion problem.
841 * - "weak_sym": Same as the "weak" option but with a modification so as to
842 * add a symmetric contribution to the system. If the problem yields a symmetric
843 * matrix. This choice is more relevant than "weak". This a good choice for a
844 * diffusion problem.
845 *
846 * For HHO and CDO-Face based schemes, only the "penalization" and "algebraic"
847 * technique is available up to now.
848 *
849 * \var CS_EQKEY_BC_QUADRATURE
850 * Set the quadrature algorithm used for evaluating integral quantities on
851 * faces or volumes. Available choices are:
852 * - "bary" used the barycenter approximation
853 * - "higher" used 4 Gauss points for approximating the integral
854 * - "highest" used 5 Gauss points for approximating the integral
855 *
856 * Remark: "higher" and "highest" implies automatically a subdivision into
857 * tetrahedra of each cell.
858 *
859 * \var CS_EQKEY_BC_STRONG_PENA_COEFF
860 * Set the value of the penalization coefficient when "penalization" is
861 * activated The default value is 1e12.
862 * cf. \ref CS_PARAM_BC_ENFORCE_PENALIZED
863 *
864 * \var CS_EQKEY_BC_WEAK_PENA_COEFF
865 * Set the value of the penalization coefficient when "weak" or "weak_sym" is
866 * activated. The default value is 100.
867 * cf. \ref CS_PARAM_BC_ENFORCE_WEAK_NITSCHE
868 * or \ref CS_PARAM_BC_ENFORCE_WEAK_SYM
869 *
870 * \var CS_EQKEY_DOF_REDUCTION
871 * Set how is defined each degree of freedom (DoF).
872 * - "de_rham" (default): Evaluation at vertices for potentials, integral
873 * along a line for circulations, integral across the normal component of a
874 * face for fluxes and integral inside a cell for densities
875 * - "average": DoF are defined as the average on the element
876 *
877 * \var CS_EQKEY_EXTRA_OP
878 * Set the additional post-processing to perform. Available choices are:
879 * - "balance" post-process the balance result in each control volume for
880 * each main term of an equation (diffusion, convection, time...)
881 * - "peclet" post-process an estimation of the Peclet number in each cell
882 * - "upwind_coef" post-process an estimation of the upwinding coefficient
883 * - "normal_flux" post-process the normal flux across boundary faces
884 *
885 * \var CS_EQKEY_HODGE_DIFF_ALGO
886 * Set the algorithm used for building the discrete Hodge operator used
887 * in the diffusion term. Available choices are:
888 * - "voronoi" --> leads to a diagonal discrete Hodge operator but is not
889 * consistent for all meshes. Require an "orthogonal" (or admissible) mesh;
890 * - "cost" --> (default for diffusion) is more robust (i.e. it handles more
891 * general meshes but is is less efficient)
892 * - "wbs" --> is robust and accurate but is limited to the reconstruction of
893 * potential-like degrees of freedom and needs a correct computation of the
894 * cell barycenter
895 *
896 * \var CS_EQKEY_HODGE_DIFF_COEF
897 * This key is only useful if CS_EQKEY_HODGE_{TIME, DIFF, REAC}_ALGO is set to
898 * "cost".
899 * keyval is either a name or a value:
900 * - "dga" corresponds to the value \f$ 1./3. \f$
901 * - "sushi" corresponds to the value \f$1./\sqrt(3.)\f$
902 * - "gcr" corresponds to the value \f$1\f$.
903 * - or "1.5", "9" for instance
904 *
905 * \var CS_EQKEY_HODGE_TIME_ALGO
906 * Set the algorithm used for building the discrete Hodge operator used
907 * in the unsteady term. Available choices are:
908 * - "voronoi" --> (default) leads to diagonal discrete Hodge operator. It acts
909 * as a mass lumping.
910 * - "wbs" is robust and accurate but is less efficient. It needs a correct
911 * computation of the cell barycenter
912 *
913 * \var CS_EQKEY_HODGE_REAC_ALGO
914 * Set the algorithm used for building the discrete Hodge operator used in the
915 * reaction term. Available choices are:
916 * - "voronoi" --> leads to diagonal discrete Hodge operator (similar to a
917 * lumping).
918 * - "wbs" --> (default) is robust and accurate but is limited to the
919 * reconstruction of potential-like degrees of freedom and needs a correct
920 * computation of the cell barycenter
921 *
922 * \var CS_EQKEY_ITSOL
923 * Specify the iterative solver for solving the linear system related to an
924 * equation. Avalaible choices are:
925
926 * - "amg" --> Algebraic MultiGrid iterative solver.
927 * Good choice for a scalable solver
928 * related to symmetric positive definite
929 * system.
930 * - "jacobi","diag" or "diagonal" --> simpliest iterative solver
931 * - "gauss_seidel" or "gs" --> Gauss-Seidel algorithm
932 * - "cg" --> conjuguate gradient algorithm
933 * - "cr3" --> a 3-layer conjugate residual solver (when "cs" is
934 * chosen as the solver family)
935 * - "fcg" --> flexible version of the conjuguate gradient
936 * algorithm used when the preconditioner can change
937 * iteration by iteration
938 * - "bicg" --> Bi-CG algorithm (for non-symmetric linear systems)
939 * - "bicgstab2" --> BiCG-Stab2 algorithm (for non-symmetric linear
940 * systems)
941 * - "gcr" --> robust and flexible iterative solver. Not the best
942 * choice if the system is easy to solve
943 * - "gmres" --> robust iterative solver. Not the best choice if the
944 * system is easy to solve
945 * - "fgmres" --> Flexible gmres (only with PETSc installation up to
946 * now). An evolutive preconditioner can be used with
947 * this solver. This is a very robust iterative solver.
948 * Not the best choice if the system is easy to solve
949 * - "minres" --> Solver of choice for symmetric indefinite systems
950 * - "mumps" --> Direct solver (very robust but memory consumming)
951 * LU factorization (through PETSc or MUMPS)
952 * - "mumps_float" --> Direct solver (very robust but memory consumming)
953 * LU factorization (through PETSc or MUMPS).
954 * Reduction of the memory consumption thanks to the usage
955 * of float instead of double. This is a recommended
956 * choice when used inside a preconditioner
957 * - "mumps_float_ldlt" --> Direct solver (very robust but memory consumming)
958 * LDLT factorization (through PETSc or MUMPS).
959 * Reduction of the memory consumption thanks to the
960 * usage of float instead of double. This is a
961 * recommended choice when used inside a preconditioner
962 * - "mumps_ldlt" --> Direct solver (very robust but memory consumming)
963 * LDLT factorization (through PETSc or MUMPS).
964 * - "sym_gauss_seidel" or "sgs" --> Symmetric Gauss-Seidel algorithm
965 * - "user" --> User-defined iterative solver (rely on the function
966 * cs_user_sles_it_solver())
967 * - "none" --> No solver.
968 *
969 * \var CS_EQKEY_ITSOL_EPS
970 * Tolerance factor for stopping the iterative processus for solving the
971 * linear system related to an equation
972 * - Example: "1e-10"
973 *
974 * \var CS_EQKEY_ITSOL_MAX_ITER
975 * Maximum number of iterations for solving the linear system
976 * - Example: "2000"
977 *
978 * \var CS_EQKEY_ITSOL_RESNORM_TYPE
979 * Normalized or not the residual before testing if one continues iterating
980 * for solving the linear system. This normalization is performed before
981 * applying the boundary conditions to avoid handling the penalization of
982 * Dirichlet boundary conditions. If the RHS norm is equal to zero, then
983 * the "vol_tot" option is used for rescaling the residual.
984 *
985 * Available choices are:
986 * "false" or "none" (default)
987 * "rhs"
988 * "weighted_rhs" or "weighted"
989 * "filtered_rhs" or "filtered_rhs"
990 *
991 * \var CS_EQKEY_ITSOL_RESTART
992 * Maximum number of iterations before restarting a Krylov solver
993 * Only useful with GMRES, flexible GMRES or GCR solvers.
994 * This value has an impact on the memory footprint since one has to store a
995 * buffer of double with a size equal to restart*sizeof(solution array)
996 * - Example: "20"
997 *
998 * \var CS_EQKEY_OMP_ASSEMBLY_STRATEGY
999 * Choice of the way to perform the assembly when OpenMP is active
1000 * Available choices are:
1001 * - "atomic" or "critical"
1002 *
1003 * \var CS_EQKEY_PRECOND
1004 * Specify the preconditioner associated to an iterative solver. Be careful
1005 * some options are only available with a given solver class. Be sure that your
1006 * installation has been installed with the appropriate library.
1007 * Available choices are:
1008 * - "none": no preconditioner is used
1009 * - "jacobi" or "diag": diagonal preconditoner
1010 * - "block_jacobi"/"bjacobi": Block Jacobi with ILU(0) in each block. A block
1011 * is the matrix associated to a rank. (Only with
1012 * PETSc)
1013 * - "bjacobi_sgs"/"bjacobi_ssor": Block Jacobi with Symm. Gauss-Seidel in each
1014 * block. A block is the matrix associated to a
1015 * rank. (Only with PETSc)
1016 * - "poly1": Neumann polynomial of order 1 (only with Code_Saturne)
1017 * - "poly2": Neumann polynomial of order 2 (only with Code_Saturne)
1018 * - "ssor": symmetric successive over-relaxation (only with PETSC)
1019 * - "ilu0": incomplete LU factorization (only with PETSc)
1020 * - "icc0": incomplete Cholesky factorization (for symmetric matrices and
1021 * only with PETSc)
1022 * - "lu": LU factorization (only with PETSc). It may use MUMPS ig PETSc is
1023 * built with MUMPS
1024 * - "amg": algebraic multigrid technique (see \ref CS_EQKEY_AMG_TYPE for
1025 additional settings)
1026 * - "amg_block"/"block_amg: algebraic multigrid by block (useful for
1027 * vector-valued equations). By default, a diagonal block
1028 * preconditioning is used if nothing else is set.
1029 *
1030 * \var CS_EQKEY_PRECOND_BLOCK_TYPE
1031 * Specify the type of block preconditioner associated to a preconditioner.
1032 * When "full" is specified. That means that the smallest possible blocks are
1033 * considered. Be careful: Most of these options are available only with
1034 * PETSc/HYPRE.
1035 * Available choices are:
1036 * - "none": no block preconditioner (default choice)
1037 * - "diag": diagonal (or additive) block preconditoner
1038 * - "lower": lower triangular (multiplicative) block preconditioner
1039 * - "upper": upper triangular (multiplicative) block preconditioner
1040 * - "symm": symmetric Gauss-Seidel block preconditioner
1041 * - "full_diag": full diagonal (or additive) block preconditoner.
1042 * - "full_lower": full lower triangular (multiplicative) block preconditioner
1043 * - "full_upper": upper triangular (multiplicative) block preconditioner
1044 * - "full_symm": symmetric Gauss-Seidel block preconditioner
1045 *
1046 * \var CS_EQKEY_SLES_VERBOSITY
1047 * Level of details written by the code for the resolution of the linear system
1048 * - Examples: "0", "1", "2" or higher
1049 *
1050 * \var CS_EQKEY_SPACE_SCHEME
1051 * Set the space discretization scheme. Available choices are:
1052 * - "cdo_vb" or "cdovb" for CDO vertex-based scheme
1053 * - "cdo_vcb" or "cdovcb" for CDO vertex+cell-based scheme
1054 * - "cdo_fb" or "cdofb" for CDO face-based scheme
1055 * - "cdo_eb" or "cdoeb" for CDO edge-based scheme
1056 * - "hho_p1" for HHO schemes with \f$\mathbb{P}_1\f$ polynomial approximation
1057 * - "hho_p2" for HHO schemes with \f$\mathbb{P}_2\f$ polynomial approximation
1058 *
1059 * \var CS_EQKEY_TIME_SCHEME
1060 * Set the scheme for the temporal discretization. Available choices are:
1061 * - "euler_implicit" or "forward_euler" (1st order)
1062 * - "euler_explicit" or "backward_euler" (1st order, conditional stability)
1063 * - "theta_scheme": Time scheme encompassing several other schemes according
1064 * to the value of theta. One recovers for instance "euler_implicit" with
1065 * theta equal to "1", "euler_explicit" with "0" or "crank_nicolson" with
1066 * theta equal to "0.5"
1067 * - "crank_nicolson": Shortcut to set a theta scheme with theta equal to
1068 * 0.5. This time discretization is second_order in time accurate.
1069 * - "bdf2": 2nd order, implicit time scheme
1070 *
1071 * \var CS_EQKEY_TIME_THETA
1072 * Set a value betwwen 0 and 1 for the theta parameter when a "theta_scheme" is
1073 * set for the time discretization. Only useful if CS_EQKEY_TIME_SCHEME is set
1074 * to "theta_scheme".
1075 * - Example: "0.75"
1076 *
1077 * \var CS_EQKEY_VERBOSITY
1078 * Set the level of details written by the code for an equation.
1079 * The higher the more detailed information is displayed.
1080 * - "0" (default)
1081 * - "1" detailed setup resume and coarse grain timer stats
1082 * - "2" fine grain for timer stats
1083 *
1084 */
1085
1086 typedef enum {
1087
1088 CS_EQKEY_ADV_EXTRAPOL,
1089 CS_EQKEY_ADV_FORMULATION,
1090 CS_EQKEY_ADV_SCHEME,
1091 CS_EQKEY_ADV_STRATEGY,
1092 CS_EQKEY_ADV_UPWIND_PORTION,
1093 CS_EQKEY_AMG_TYPE,
1094 CS_EQKEY_BC_ENFORCEMENT,
1095 CS_EQKEY_BC_QUADRATURE,
1096 CS_EQKEY_BC_STRONG_PENA_COEFF,
1097 CS_EQKEY_BC_WEAK_PENA_COEFF,
1098 CS_EQKEY_DO_LUMPING,
1099 CS_EQKEY_DOF_REDUCTION,
1100 CS_EQKEY_EXTRA_OP,
1101 CS_EQKEY_HODGE_DIFF_ALGO,
1102 CS_EQKEY_HODGE_DIFF_COEF,
1103 CS_EQKEY_HODGE_TIME_ALGO,
1104 CS_EQKEY_HODGE_REAC_ALGO,
1105 CS_EQKEY_ITSOL,
1106 CS_EQKEY_ITSOL_EPS,
1107 CS_EQKEY_ITSOL_MAX_ITER,
1108 CS_EQKEY_ITSOL_RESNORM_TYPE,
1109 CS_EQKEY_ITSOL_RESTART,
1110 CS_EQKEY_OMP_ASSEMBLY_STRATEGY,
1111 CS_EQKEY_PRECOND,
1112 CS_EQKEY_PRECOND_BLOCK_TYPE,
1113 CS_EQKEY_SLES_VERBOSITY,
1114 CS_EQKEY_SOLVER_FAMILY,
1115 CS_EQKEY_SPACE_SCHEME,
1116 CS_EQKEY_TIME_SCHEME,
1117 CS_EQKEY_TIME_THETA,
1118 CS_EQKEY_VERBOSITY,
1119
1120 CS_EQKEY_N_KEYS
1121
1122 } cs_equation_key_t;
1123
1124 /*============================================================================
1125 * Static inline public function prototypes
1126 *============================================================================*/
1127
1128 /*----------------------------------------------------------------------------*/
1129 /*!
1130 * \brief Update the flag related to a cs_equation_param_t structure
1131 *
1132 * \param[in, out] eqp pointer to a \ref cs_equation_param_t
1133 * \param[in] flag flag to add
1134 */
1135 /*----------------------------------------------------------------------------*/
1136
1137 static inline void
cs_equation_param_set_flag(cs_equation_param_t * eqp,cs_flag_t flag)1138 cs_equation_param_set_flag(cs_equation_param_t *eqp,
1139 cs_flag_t flag)
1140 {
1141 assert(eqp != NULL);
1142 eqp->flag |= flag;
1143 }
1144
1145 /*----------------------------------------------------------------------------*/
1146 /*!
1147 * \brief Ask if the parameters of the equation needs a diffusion term
1148 *
1149 * \param[in] eqp pointer to a \ref cs_equation_param_t
1150 *
1151 * \return true or false
1152 */
1153 /*----------------------------------------------------------------------------*/
1154
1155 static inline bool
cs_equation_param_has_diffusion(const cs_equation_param_t * eqp)1156 cs_equation_param_has_diffusion(const cs_equation_param_t *eqp)
1157 {
1158 assert(eqp != NULL);
1159 if (eqp->flag & CS_EQUATION_DIFFUSION)
1160 return true;
1161 else
1162 return false;
1163 }
1164
1165 /*----------------------------------------------------------------------------*/
1166 /*!
1167 * \brief Ask if the parameters of the equation needs a curl-curl term
1168 *
1169 * \param[in] eqp pointer to a \ref cs_equation_param_t
1170 *
1171 * \return true or false
1172 */
1173 /*----------------------------------------------------------------------------*/
1174
1175 static inline bool
cs_equation_param_has_curlcurl(const cs_equation_param_t * eqp)1176 cs_equation_param_has_curlcurl(const cs_equation_param_t *eqp)
1177 {
1178 assert(eqp != NULL);
1179 if (eqp->flag & CS_EQUATION_CURLCURL)
1180 return true;
1181 else
1182 return false;
1183 }
1184
1185 /*----------------------------------------------------------------------------*/
1186 /*!
1187 * \brief Ask if the parameters of the equation needs a grad-div term
1188 *
1189 * \param[in] eqp pointer to a \ref cs_equation_param_t
1190 *
1191 * \return true or false
1192 */
1193 /*----------------------------------------------------------------------------*/
1194
1195 static inline bool
cs_equation_param_has_graddiv(const cs_equation_param_t * eqp)1196 cs_equation_param_has_graddiv(const cs_equation_param_t *eqp)
1197 {
1198 assert(eqp != NULL);
1199 if (eqp->flag & CS_EQUATION_GRADDIV)
1200 return true;
1201 else
1202 return false;
1203 }
1204
1205 /*----------------------------------------------------------------------------*/
1206 /*!
1207 * \brief Ask if the parameters of the equation needs a convection term
1208 *
1209 * \param[in] eqp pointer to a \ref cs_equation_param_t
1210 *
1211 * \return true or false
1212 */
1213 /*----------------------------------------------------------------------------*/
1214
1215 static inline bool
cs_equation_param_has_convection(const cs_equation_param_t * eqp)1216 cs_equation_param_has_convection(const cs_equation_param_t *eqp)
1217 {
1218 assert(eqp != NULL);
1219 if (eqp->flag & CS_EQUATION_CONVECTION)
1220 return true;
1221 else
1222 return false;
1223 }
1224
1225 /*----------------------------------------------------------------------------*/
1226 /*!
1227 * \brief Ask if the parameters of the equation needs a reaction term
1228 *
1229 * \param[in] eqp pointer to a \ref cs_equation_param_t
1230 *
1231 * \return true or false
1232 */
1233 /*----------------------------------------------------------------------------*/
1234
1235 static inline bool
cs_equation_param_has_reaction(const cs_equation_param_t * eqp)1236 cs_equation_param_has_reaction(const cs_equation_param_t *eqp)
1237 {
1238 assert(eqp != NULL);
1239 if (eqp->flag & CS_EQUATION_REACTION)
1240 return true;
1241 else
1242 return false;
1243 }
1244
1245 /*----------------------------------------------------------------------------*/
1246 /*!
1247 * \brief Ask if the parameters of the equation needs an unsteady term
1248 *
1249 * \param[in] eqp pointer to a \ref cs_equation_param_t
1250 *
1251 * \return true or false
1252 */
1253 /*----------------------------------------------------------------------------*/
1254
1255 static inline bool
cs_equation_param_has_time(const cs_equation_param_t * eqp)1256 cs_equation_param_has_time(const cs_equation_param_t *eqp)
1257 {
1258 assert(eqp != NULL);
1259 if (eqp->flag & CS_EQUATION_UNSTEADY)
1260 return true;
1261 else
1262 return false;
1263 }
1264
1265 /*----------------------------------------------------------------------------*/
1266 /*!
1267 * \brief Ask if the parameters of the equation needs a source term
1268 *
1269 * \param[in] eqp pointer to a \ref cs_equation_param_t
1270 *
1271 * \return true or false
1272 */
1273 /*----------------------------------------------------------------------------*/
1274
1275 static inline bool
cs_equation_param_has_sourceterm(const cs_equation_param_t * eqp)1276 cs_equation_param_has_sourceterm(const cs_equation_param_t *eqp)
1277 {
1278 assert(eqp != NULL);
1279 if (eqp->n_source_terms > 0)
1280 return true;
1281 else
1282 return false;
1283 }
1284
1285 /*----------------------------------------------------------------------------*/
1286 /*!
1287 * \brief Ask if the parameters of the equation has an internal enforcement
1288 * of the degrees of freedom
1289 *
1290 * \param[in] eqp pointer to a \ref cs_equation_param_t
1291 *
1292 * \return true or false
1293 */
1294 /*----------------------------------------------------------------------------*/
1295
1296 static inline bool
cs_equation_param_has_internal_enforcement(const cs_equation_param_t * eqp)1297 cs_equation_param_has_internal_enforcement(const cs_equation_param_t *eqp)
1298 {
1299 assert(eqp != NULL);
1300 if (eqp->flag & CS_EQUATION_FORCE_VALUES)
1301 return true;
1302 else
1303 return false;
1304 }
1305
1306 /*----------------------------------------------------------------------------*/
1307 /*!
1308 * \brief Ask if the parameters of the equation induces an implicit treatment
1309 * of the advection term
1310 *
1311 * \param[in] eqp pointer to a \ref cs_equation_param_t
1312 *
1313 * \return true or false
1314 */
1315 /*----------------------------------------------------------------------------*/
1316
1317 static inline bool
cs_equation_param_has_implicit_advection(const cs_equation_param_t * eqp)1318 cs_equation_param_has_implicit_advection(const cs_equation_param_t *eqp)
1319 {
1320 assert(eqp != NULL);
1321 if (eqp->flag & CS_EQUATION_CONVECTION) {
1322 if (eqp->adv_strategy == CS_PARAM_ADVECTION_IMPLICIT_FULL ||
1323 eqp->adv_strategy == CS_PARAM_ADVECTION_IMPLICIT_LINEARIZED)
1324 return true;
1325 else
1326 return false;
1327 }
1328 else
1329 return false;
1330 }
1331
1332 /*----------------------------------------------------------------------------*/
1333 /*!
1334 * \brief Ask if the parameters of the equation has activated a user hook
1335 * to get a fine tuning of the cellwise system building
1336 *
1337 * \param[in] eqp pointer to a \ref cs_equation_param_t
1338 *
1339 * \return true or false
1340 */
1341 /*----------------------------------------------------------------------------*/
1342
1343 static inline bool
cs_equation_param_has_user_hook(const cs_equation_param_t * eqp)1344 cs_equation_param_has_user_hook(const cs_equation_param_t *eqp)
1345 {
1346 assert(eqp != NULL);
1347 if (eqp->flag & CS_EQUATION_USER_HOOK)
1348 return true;
1349 else
1350 return false;
1351 }
1352
1353 /*----------------------------------------------------------------------------*/
1354 /*!
1355 * \brief Check if a \ref cs_equation_param_t structure has its name member
1356 * equal to the given name
1357 *
1358 * \param[in] eqp pointer to a \ref cs_equation_param_t structure
1359 * \param[in] name name of the equation
1360 *
1361 * \return true if the given eqp has the same name the one given as parameter
1362 * otherwise false
1363 */
1364 /*----------------------------------------------------------------------------*/
1365
1366 static inline bool
cs_equation_param_has_name(cs_equation_param_t * eqp,const char * name)1367 cs_equation_param_has_name(cs_equation_param_t *eqp,
1368 const char *name)
1369 {
1370 if (eqp == NULL)
1371 return false;
1372 if (eqp->name == NULL)
1373 return false;
1374 if (strcmp(eqp->name, name) == 0)
1375 return true;
1376 else
1377 return false;
1378 }
1379
1380 /*============================================================================
1381 * Public function prototypes
1382 *============================================================================*/
1383
1384 /*----------------------------------------------------------------------------*/
1385 /*!
1386 * \brief Create a \ref cs_equation_param_t structure
1387 *
1388 * \param[in] name name of the equation
1389 * \param[in] type type of equation
1390 * \param[in] dim dim of the variable associated to this equation
1391 * \param[in] default_bc type of boundary condition set by default
1392 *
1393 * \return a pointer to a new allocated \ref cs_equation_param_t structure
1394 */
1395 /*----------------------------------------------------------------------------*/
1396
1397 cs_equation_param_t *
1398 cs_equation_param_create(const char *name,
1399 cs_equation_type_t type,
1400 int dim,
1401 cs_param_bc_type_t default_bc);
1402
1403 /*----------------------------------------------------------------------------*/
1404 /*!
1405 * \brief Create a \ref cs_equation_param_t structure
1406 *
1407 * \deprecated Renamed to\ref cs_equation_param_create.
1408 *
1409 * \param[in] name name of the equation
1410 * \param[in] type type of equation
1411 * \param[in] dim dim of the variable associated to this equation
1412 * \param[in] default_bc type of boundary condition set by default
1413 *
1414 * \return a pointer to a new allocated \ref cs_equation_param_t structure
1415 */
1416 /*----------------------------------------------------------------------------*/
1417
1418 inline static cs_equation_param_t *
cs_equation_create_param(const char * name,cs_equation_type_t type,int dim,cs_param_bc_type_t default_bc)1419 cs_equation_create_param(const char *name,
1420 cs_equation_type_t type,
1421 int dim,
1422 cs_param_bc_type_t default_bc)
1423 {
1424 return cs_equation_param_create(name, type, dim, default_bc);
1425 }
1426
1427 /*----------------------------------------------------------------------------*/
1428 /*!
1429 * \brief Copy the settings from one \ref cs_equation_param_t structure to
1430 * another one. The name is not copied.
1431 *
1432 * \param[in] ref pointer to the reference \ref cs_equation_param_t
1433 * \param[in, out] dst pointer to the \ref cs_equation_param_t to update
1434 * \param[in] copy_fid copy also the field id or not
1435 */
1436 /*----------------------------------------------------------------------------*/
1437
1438 void
1439 cs_equation_param_copy_from(const cs_equation_param_t *ref,
1440 cs_equation_param_t *dst,
1441 bool copy_fid);
1442
1443 /*----------------------------------------------------------------------------*/
1444 /*!
1445 * \brief Copy the settings from one \ref cs_equation_param_t structure to
1446 * another one. The name is not copied.
1447 *
1448 * \deprecated Renamed to\ref cs_equation_param_copy_from.
1449 *
1450 * \param[in] ref pointer to the reference \ref cs_equation_param_t
1451 * \param[in, out] dst pointer to the \ref cs_equation_param_t to update
1452 * \param[in] copy_fid copy also the field id or not
1453 */
1454 /*----------------------------------------------------------------------------*/
1455
1456 inline static void
cs_equation_copy_param_from(const cs_equation_param_t * ref,cs_equation_param_t * dst,bool copy_fid)1457 cs_equation_copy_param_from(const cs_equation_param_t *ref,
1458 cs_equation_param_t *dst,
1459 bool copy_fid)
1460 {
1461 cs_equation_param_copy_from(ref, dst, copy_fid);
1462 }
1463
1464 /*----------------------------------------------------------------------------*/
1465 /*!
1466 * \brief Free the contents of a \ref cs_equation_param_t
1467 *
1468 * The cs_equation_param_t structure itself is not freed, but all its
1469 * sub-structures are freed.
1470 *
1471 * This is useful for equation parameters which are accessed through
1472 * field keywords.
1473 *
1474 * \param[in, out] eqp pointer to a \ref cs_equation_param_t
1475 */
1476 /*----------------------------------------------------------------------------*/
1477
1478 void
1479 cs_equation_param_clear(cs_equation_param_t *eqp);
1480
1481 /*----------------------------------------------------------------------------*/
1482 /*!
1483 * \brief Free the contents of a \ref cs_equation_param_t
1484 *
1485 * \deprecated Renamed to\ref cs_equation_param_clear.
1486 *
1487 * The cs_equation_param_t structure itself is not freed, but all its
1488 * sub-structures are freed.
1489 *
1490 * This is useful for equation parameters which are accessed through
1491 * field keywords.
1492 *
1493 * \param[in, out] eqp pointer to a \ref cs_equation_param_t
1494 */
1495 /*----------------------------------------------------------------------------*/
1496
1497 inline static void
cs_equation_clear_param(cs_equation_param_t * eqp)1498 cs_equation_clear_param(cs_equation_param_t *eqp)
1499 {
1500 cs_equation_param_clear(eqp);
1501 }
1502
1503 /*----------------------------------------------------------------------------*/
1504 /*!
1505 * \brief Free a \ref cs_equation_param_t
1506 *
1507 * \param[in] eqp pointer to a \ref cs_equation_param_t
1508 *
1509 * \return a NULL pointer
1510 */
1511 /*----------------------------------------------------------------------------*/
1512
1513 cs_equation_param_t *
1514 cs_equation_param_free(cs_equation_param_t *eqp);
1515
1516 /*----------------------------------------------------------------------------*/
1517 /*!
1518 * \brief Free a \ref cs_equation_param_t
1519 *
1520 * \deprecated Renamed to\ref cs_equation_param_free.
1521 *
1522 * \param[in] eqp pointer to a \ref cs_equation_param_t
1523 *
1524 * \return a NULL pointer
1525 */
1526 /*----------------------------------------------------------------------------*/
1527
1528 inline static cs_equation_param_t *
cs_equation_free_param(cs_equation_param_t * eqp)1529 cs_equation_free_param(cs_equation_param_t *eqp)
1530 {
1531 return cs_equation_param_free(eqp);
1532 }
1533
1534 /*----------------------------------------------------------------------------*/
1535 /*!
1536 * \brief Set a parameter attached to a keyname in a \ref cs_equation_param_t
1537 * structure
1538 *
1539 * \param[in, out] eqp pointer to a \ref cs_equation_param_t structure
1540 * \param[in] key key related to the member of eq to set
1541 * \param[in] keyval accessor to the value to set
1542 */
1543 /*----------------------------------------------------------------------------*/
1544
1545 void
1546 cs_equation_param_set(cs_equation_param_t *eqp,
1547 cs_equation_key_t key,
1548 const char *keyval);
1549
1550 /*----------------------------------------------------------------------------*/
1551 /*!
1552 * \brief Set a parameter attached to a keyname in a \ref cs_equation_param_t
1553 * structure
1554 *
1555 * \deprecated Renamed to\ref cs_equation_param_set.
1556 *
1557 * \param[in, out] eqp pointer to a \ref cs_equation_param_t structure
1558 * \param[in] key key related to the member of eq to set
1559 * \param[in] keyval accessor to the value to set
1560 */
1561 /*----------------------------------------------------------------------------*/
1562
1563 inline static void
cs_equation_set_param(cs_equation_param_t * eqp,cs_equation_key_t key,const char * keyval)1564 cs_equation_set_param(cs_equation_param_t *eqp,
1565 cs_equation_key_t key,
1566 const char *keyval)
1567 {
1568 cs_equation_param_set(eqp, key, keyval);
1569 }
1570
1571 /*----------------------------------------------------------------------------*/
1572 /*!
1573 * \brief Set parameters for initializing SLES structures used for the
1574 * resolution of the linear system.
1575 * Settings are related to this equation.
1576 *
1577 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1578 */
1579 /*----------------------------------------------------------------------------*/
1580
1581 void
1582 cs_equation_param_set_sles(cs_equation_param_t *eqp);
1583
1584 /*----------------------------------------------------------------------------*/
1585 /*!
1586 * \brief Last modification of the cs_equation_param_t structure before
1587 * launching the computation
1588 *
1589 * \param[in, out] eqp pointer to a \ref cs_equation_param_t structure
1590 */
1591 /*----------------------------------------------------------------------------*/
1592
1593 void
1594 cs_equation_param_last_stage(cs_equation_param_t *eqp);
1595
1596 /*----------------------------------------------------------------------------*/
1597 /*!
1598 * \brief Print the detail of a \ref cs_equation_param_t structure
1599 *
1600 * \param[in] eqp pointer to a \ref cs_equation_param_t structure
1601 */
1602 /*----------------------------------------------------------------------------*/
1603
1604 void
1605 cs_equation_param_log(const cs_equation_param_t *eqp);
1606
1607 /*----------------------------------------------------------------------------*/
1608 /*!
1609 * \brief Ask if the parameter settings of the equation has requested the
1610 * treatment of Robin boundary conditions
1611 *
1612 * \param[in] eqp pointer to a \ref cs_equation_param_t
1613 *
1614 * \return true or false
1615 */
1616 /*----------------------------------------------------------------------------*/
1617
1618 bool
1619 cs_equation_param_has_robin_bc(const cs_equation_param_t *eqp);
1620
1621 /*----------------------------------------------------------------------------*/
1622 /*!
1623 * \brief Define the initial condition for the unknown related to this equation
1624 * This definition can be done on a specified mesh location.
1625 * By default, the unknown is set to zero everywhere.
1626 * Here a constant value is set to all the entities belonging to the
1627 * given mesh location
1628 *
1629 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1630 * \param[in] z_name name of the associated zone (if NULL or
1631 * "" all cells are considered)
1632 * \param[in] val pointer to the value
1633 *
1634 * \return a pointer to the new \ref cs_xdef_t structure
1635 */
1636 /*----------------------------------------------------------------------------*/
1637
1638 cs_xdef_t *
1639 cs_equation_add_ic_by_value(cs_equation_param_t *eqp,
1640 const char *z_name,
1641 cs_real_t *val);
1642
1643 /*----------------------------------------------------------------------------*/
1644 /*!
1645 * \brief Define the initial condition for the unknown related to this equation
1646 * This definition can be done on a specified mesh location.
1647 * By default, the unknown is set to zero everywhere.
1648 * Here the value related to all the entities belonging to the
1649 * given mesh location is such that the integral over these cells
1650 * returns the requested quantity
1651 *
1652 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1653 * \param[in] z_name name of the associated zone (if NULL or
1654 * "" all cells are considered)
1655 * \param[in] quantity quantity to distribute over the mesh location
1656 *
1657 * \return a pointer to the new \ref cs_xdef_t structure
1658 */
1659 /*----------------------------------------------------------------------------*/
1660
1661 cs_xdef_t *
1662 cs_equation_add_ic_by_qov(cs_equation_param_t *eqp,
1663 const char *z_name,
1664 double quantity);
1665
1666 /*----------------------------------------------------------------------------*/
1667 /*!
1668 * \brief Define the initial condition for the unknown related to this
1669 * equation. This definition can be done on a specified mesh location.
1670 * By default, the unknown is set to zero everywhere.
1671 * Here the initial value is set according to an analytical function
1672 *
1673 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1674 * \param[in] z_name name of the associated zone (if NULL or "" if
1675 * all cells are considered)
1676 * \param[in] analytic pointer to an analytic function
1677 * \param[in] input NULL or pointer to a structure cast on-the-fly
1678 *
1679 * \return a pointer to the new \ref cs_xdef_t structure
1680 */
1681 /*----------------------------------------------------------------------------*/
1682
1683 cs_xdef_t *
1684 cs_equation_add_ic_by_analytic(cs_equation_param_t *eqp,
1685 const char *z_name,
1686 cs_analytic_func_t *analytic,
1687 void *input);
1688
1689 /*----------------------------------------------------------------------------*/
1690 /*!
1691 * \brief Set a boundary condition from an existing \ref cs_xdef_t structure
1692 * The lifecycle of the cs_xdef_t structure is now managed by the
1693 * current \ref cs_equation_param_t structure.
1694 *
1695 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1696 * \param[in] xdef pointer to the \ref cs_xdef_t structure to transfer
1697 */
1698 /*----------------------------------------------------------------------------*/
1699
1700 void
1701 cs_equation_add_xdef_bc(cs_equation_param_t *eqp,
1702 cs_xdef_t *xdef);
1703
1704 /*----------------------------------------------------------------------------*/
1705 /*!
1706 * \brief Define and initialize a new structure to set a boundary condition
1707 * related to the given equation structure
1708 * z_name corresponds to the name of a pre-existing cs_zone_t
1709 *
1710 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1711 * \param[in] bc_type type of boundary condition to add
1712 * \param[in] z_name name of the related boundary zone
1713 * \param[in] values pointer to a array storing the values
1714 *
1715 * \return a pointer to the new \ref cs_xdef_t structure
1716 */
1717 /*----------------------------------------------------------------------------*/
1718
1719 cs_xdef_t *
1720 cs_equation_add_bc_by_value(cs_equation_param_t *eqp,
1721 const cs_param_bc_type_t bc_type,
1722 const char *z_name,
1723 cs_real_t *values);
1724
1725 /*----------------------------------------------------------------------------*/
1726 /*!
1727 * \brief Define and initialize a new structure to set a boundary condition
1728 * related to the given equation structure
1729 * z_name corresponds to the name of a pre-existing cs_zone_t
1730 *
1731 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1732 * \param[in] bc_type type of boundary condition to add
1733 * \param[in] z_name name of the related boundary zone
1734 * \param[in] loc information to know where are located values
1735 * \param[in] array pointer to an array
1736 * \param[in] is_owner transfer the lifecycle to the cs_xdef_t structure
1737 * (true or false)
1738 * \param[in] index optional pointer to the array index
1739 *
1740 * \return a pointer to the new allocated \ref cs_xdef_t structure
1741 */
1742 /*----------------------------------------------------------------------------*/
1743
1744 cs_xdef_t *
1745 cs_equation_add_bc_by_array(cs_equation_param_t *eqp,
1746 const cs_param_bc_type_t bc_type,
1747 const char *z_name,
1748 cs_flag_t loc,
1749 cs_real_t *array,
1750 bool is_owner,
1751 cs_lnum_t *index);
1752
1753 /*----------------------------------------------------------------------------*/
1754 /*!
1755 * \brief Define and initialize a new structure to set a boundary condition
1756 * related to the given equation structure
1757 * z_name corresponds to the name of a pre-existing cs_zone_t
1758 *
1759 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1760 * \param[in] bc_type type of boundary condition to add
1761 * \param[in] z_name name of the related boundary zone
1762 * \param[in] field pointer to a cs_field_t structure
1763 *
1764 * \return a pointer to the new allocated \ref cs_xdef_t structure
1765 */
1766 /*----------------------------------------------------------------------------*/
1767
1768 cs_xdef_t *
1769 cs_equation_add_bc_by_field(cs_equation_param_t *eqp,
1770 const cs_param_bc_type_t bc_type,
1771 const char *z_name,
1772 cs_field_t *field);
1773
1774 /*----------------------------------------------------------------------------*/
1775 /*!
1776 * \brief Define and initialize a new structure to set a boundary condition
1777 * related to the given equation param structure
1778 * ml_name corresponds to the name of a pre-existing cs_mesh_location_t
1779 *
1780 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1781 * \param[in] bc_type type of boundary condition to add
1782 * \param[in] z_name name of the associated zone (if NULL or "" if
1783 * all cells are considered)
1784 * \param[in] analytic pointer to an analytic function defining the value
1785 * \param[in] input NULL or pointer to a structure cast on-the-fly
1786 *
1787 * \return a pointer to the new \ref cs_xdef_t structure
1788 */
1789 /*----------------------------------------------------------------------------*/
1790
1791 cs_xdef_t *
1792 cs_equation_add_bc_by_analytic(cs_equation_param_t *eqp,
1793 const cs_param_bc_type_t bc_type,
1794 const char *z_name,
1795 cs_analytic_func_t *analytic,
1796 void *input);
1797
1798 /*----------------------------------------------------------------------------*/
1799 /*!
1800 * \brief Define and initialize a new structure to set a boundary condition
1801 * related to the given cs_equation_param_t structure
1802 * ml_name corresponds to the name of a pre-existing cs_mesh_location_t
1803 *
1804 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1805 * \param[in] bc_type type of boundary condition to add
1806 * \param[in] z_name name of the associated zone (if NULL or "" if
1807 * all cells are considered)
1808 * \param[in] loc_flag location where values are computed
1809 * \param[in] func pointer to cs_dof_func_t function
1810 * \param[in] input NULL or pointer to a structure cast on-the-fly
1811 *
1812 * \return a pointer to the new \ref cs_xdef_t structure
1813 */
1814 /*----------------------------------------------------------------------------*/
1815
1816 cs_xdef_t *
1817 cs_equation_add_bc_by_dof_func(cs_equation_param_t *eqp,
1818 const cs_param_bc_type_t bc_type,
1819 const char *z_name,
1820 cs_flag_t loc_flag,
1821 cs_dof_func_t *func,
1822 void *input);
1823
1824 /*----------------------------------------------------------------------------*/
1825 /*!
1826 * \brief Return pointer to existing boundary condition definition structure
1827 * for the given equation param structure and zone.
1828 *
1829 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1830 * \param[in] z_name name of the associated zone (if NULL or "" if
1831 * all cells are considered)
1832 *
1833 * \return a pointer to the \ref cs_xdef_t structure if present, or NULL
1834 */
1835 /*----------------------------------------------------------------------------*/
1836
1837 cs_xdef_t *
1838 cs_equation_find_bc(cs_equation_param_t *eqp,
1839 const char *z_name);
1840
1841 /*----------------------------------------------------------------------------*/
1842 /*!
1843 * \brief Remove boundary condition from the given equation param structure
1844 * for a given zone.
1845 *
1846 * If no matching boundary condition is found, the function returns
1847 * silently.
1848 *
1849 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1850 * \param[in] z_name name of the associated zone (if NULL or "" if
1851 * all cells are considered)
1852 */
1853 /*----------------------------------------------------------------------------*/
1854
1855 void
1856 cs_equation_remove_bc(cs_equation_param_t *eqp,
1857 const char *z_name);
1858
1859 /*----------------------------------------------------------------------------*/
1860 /*!
1861 * \brief Define and initialize a new structure to set a sliding boundary
1862 * condition related to the given equation structure
1863 * z_name corresponds to the name of a pre-existing cs_zone_t
1864 *
1865 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1866 * \param[in] z_name name of the related boundary zone
1867 */
1868 /*----------------------------------------------------------------------------*/
1869
1870 void
1871 cs_equation_add_sliding_condition(cs_equation_param_t *eqp,
1872 const char *z_name);
1873
1874 /*----------------------------------------------------------------------------*/
1875 /*!
1876 * \brief Associate a new term related to the Laplacian operator for the
1877 * equation associated to the given \ref cs_equation_param_t structure
1878 * Laplacian means div-grad (either for vector-valued or scalar-valued
1879 * equations)
1880 *
1881 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1882 * \param[in] property pointer to a cs_property_t structure
1883 */
1884 /*----------------------------------------------------------------------------*/
1885
1886 void
1887 cs_equation_add_diffusion(cs_equation_param_t *eqp,
1888 cs_property_t *property);
1889
1890 /*----------------------------------------------------------------------------*/
1891 /*!
1892 * \brief Associate a new term related to the curl-curl operator for the
1893 * equation associated to the given \ref cs_equation_param_t structure
1894 *
1895 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1896 * \param[in] property pointer to a cs_property_t structure
1897 * \param[in] inversion 1: true, false otherwise
1898 */
1899 /*----------------------------------------------------------------------------*/
1900
1901 void
1902 cs_equation_add_curlcurl(cs_equation_param_t *eqp,
1903 cs_property_t *property,
1904 int inversion);
1905
1906 /*----------------------------------------------------------------------------*/
1907 /*!
1908 * \brief Associate a new term related to the grad-div operator for the
1909 * equation associated to the given \ref cs_equation_param_t structure
1910 *
1911 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1912 * \param[in] property pointer to a cs_property_t structure
1913 */
1914 /*----------------------------------------------------------------------------*/
1915
1916 void
1917 cs_equation_add_graddiv(cs_equation_param_t *eqp,
1918 cs_property_t *property);
1919
1920 /*----------------------------------------------------------------------------*/
1921 /*!
1922 * \brief Associate a new term related to the time derivative operator for
1923 * the equation associated to the given \ref cs_equation_param_t
1924 * structure
1925 *
1926 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1927 * \param[in] property pointer to a cs_property_t structure
1928 */
1929 /*----------------------------------------------------------------------------*/
1930
1931 void
1932 cs_equation_add_time(cs_equation_param_t *eqp,
1933 cs_property_t *property);
1934
1935 /*----------------------------------------------------------------------------*/
1936 /*!
1937 * \brief Associate a new term related to the advection operator for the
1938 * equation associated to the given \ref cs_equation_param_t structure
1939 *
1940 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1941 * \param[in] adv_field pointer to a cs_adv_field_t structure
1942 */
1943 /*----------------------------------------------------------------------------*/
1944
1945 void
1946 cs_equation_add_advection(cs_equation_param_t *eqp,
1947 cs_adv_field_t *adv_field);
1948
1949 /*----------------------------------------------------------------------------*/
1950 /*!
1951 * \brief Associate a scaling property to the advection
1952 *
1953 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1954 * \param[in] property pointer to a cs_property_t structure
1955 */
1956 /*----------------------------------------------------------------------------*/
1957
1958 void
1959 cs_equation_add_advection_scaling_property(cs_equation_param_t *eqp,
1960 cs_property_t *property);
1961
1962 /*----------------------------------------------------------------------------*/
1963 /*!
1964 * \brief Associate a new term related to the reaction operator for the
1965 * equation associated to the given \ref cs_equation_param_t structure
1966 *
1967 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1968 * \param[in] property pointer to a cs_property_t structure
1969 *
1970 * \return the id related to the reaction term
1971 */
1972 /*----------------------------------------------------------------------------*/
1973
1974 int
1975 cs_equation_add_reaction(cs_equation_param_t *eqp,
1976 cs_property_t *property);
1977
1978 /*----------------------------------------------------------------------------*/
1979 /*!
1980 * \brief Add a new source term by initializing a cs_xdef_t structure.
1981 * Case of a definition by a constant value
1982 *
1983 * \param[in, out] eqp pointer to a cs_equation_param_t structure
1984 * \param[in] z_name name of the associated zone (if NULL or
1985 * "" all cells are considered)
1986 * \param[in] val pointer to the value
1987 *
1988 * \return a pointer to the new \ref cs_xdef_t structure
1989 */
1990 /*----------------------------------------------------------------------------*/
1991
1992 cs_xdef_t *
1993 cs_equation_add_source_term_by_val(cs_equation_param_t *eqp,
1994 const char *z_name,
1995 cs_real_t *val);
1996
1997 /*----------------------------------------------------------------------------*/
1998 /*!
1999 * \brief Add a new source term by initializing a cs_xdef_t structure.
2000 * Case of a definition by an analytical function
2001 *
2002 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2003 * \param[in] z_name name of the associated zone (if NULL or "" if
2004 * all cells are considered)
2005 * \param[in] func pointer to an analytical function
2006 * \param[in] input NULL or pointer to a structure cast on-the-fly
2007 *
2008 * \return a pointer to the new \ref cs_xdef_t structure
2009 */
2010 /*----------------------------------------------------------------------------*/
2011
2012 cs_xdef_t *
2013 cs_equation_add_source_term_by_analytic(cs_equation_param_t *eqp,
2014 const char *z_name,
2015 cs_analytic_func_t *func,
2016 void *input);
2017
2018 /*----------------------------------------------------------------------------*/
2019 /*!
2020 * \brief Add a new source term by initializing a cs_xdef_t structure.
2021 * Case of a definition by a DoF function
2022 *
2023 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2024 * \param[in] z_name name of the associated zone (if NULL or "" if
2025 * all cells are considered)
2026 * \param[in] loc_flag location of element ids given as parameter
2027 * \param[in] func pointer to a DoF function
2028 * \param[in] input NULL or pointer to a structure cast on-the-fly
2029 *
2030 * \return a pointer to the new \ref cs_xdef_t structure
2031 */
2032 /*----------------------------------------------------------------------------*/
2033
2034 cs_xdef_t *
2035 cs_equation_add_source_term_by_dof_func(cs_equation_param_t *eqp,
2036 const char *z_name,
2037 cs_flag_t loc_flag,
2038 cs_dof_func_t *func,
2039 void *input);
2040
2041 /*----------------------------------------------------------------------------*/
2042 /*!
2043 * \brief Add a new source term by initializing a cs_xdef_t structure.
2044 * Case of a definition by an array.
2045 *
2046 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2047 * \param[in] z_name name of the associated zone (if NULL or "" if
2048 * all cells are considered)
2049 * \param[in] loc information to know where are located values
2050 * \param[in] array pointer to an array
2051 * \param[in] is_owner transfer the lifecycle to the cs_xdef_t structure
2052 * (true or false)
2053 * \param[in] index optional pointer to the array index
2054 *
2055 * \return a pointer to the new \ref cs_xdef_t structure
2056 */
2057 /*----------------------------------------------------------------------------*/
2058
2059 cs_xdef_t *
2060 cs_equation_add_source_term_by_array(cs_equation_param_t *eqp,
2061 const char *z_name,
2062 cs_flag_t loc,
2063 cs_real_t *array,
2064 bool is_owner,
2065 cs_lnum_t *index);
2066
2067 /*----------------------------------------------------------------------------*/
2068 /*!
2069 * \brief Add a new volume mass injection definition source term by
2070 * initializing a cs_xdef_t structure, using a constant value.
2071 *
2072 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2073 * \param[in] z_name name of the associated zone (if NULL or "" if
2074 * all cells are considered)
2075 * \param[in] val pointer to the value
2076 *
2077 * \return a pointer to the new \ref cs_xdef_t structure
2078 */
2079 /*----------------------------------------------------------------------------*/
2080
2081 cs_xdef_t *
2082 cs_equation_add_volume_mass_injection_by_value(cs_equation_param_t *eqp,
2083 const char *z_name,
2084 double *val);
2085
2086 /*----------------------------------------------------------------------------*/
2087 /*!
2088 * \brief Add a new volume mass injection definition source term by
2089 * initializing a cs_xdef_t structure, using a constant quantity
2090 * distributed over the associated zone's volume.
2091 *
2092 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2093 * \param[in] z_name name of the associated zone (if NULL or "" if
2094 * all cells are considered)
2095 * \param[in] quantity pointer to quantity to distribute over the zone
2096 *
2097 * \return a pointer to the new \ref cs_xdef_t structure
2098 */
2099 /*----------------------------------------------------------------------------*/
2100
2101 cs_xdef_t *
2102 cs_equation_add_volume_mass_injection_by_qov(cs_equation_param_t *eqp,
2103 const char *z_name,
2104 double *quantity);
2105
2106 /*----------------------------------------------------------------------------*/
2107 /*!
2108 * \brief Add a new volume mass injection definition source term by
2109 * initializing a cs_xdef_t structure, using an analytical function.
2110 *
2111 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2112 * \param[in] z_name name of the associated zone (if NULL or "" if
2113 * all cells are considered)
2114 * \param[in] func pointer to an analytical function
2115 * \param[in] input NULL or pointer to a structure cast on-the-fly
2116 *
2117 * \return a pointer to the new \ref cs_xdef_t structure
2118 */
2119 /*----------------------------------------------------------------------------*/
2120
2121 cs_xdef_t *
2122 cs_equation_add_volume_mass_injection_by_analytic(cs_equation_param_t *eqp,
2123 const char *z_name,
2124 cs_analytic_func_t *func,
2125 void *input);
2126
2127 /*----------------------------------------------------------------------------*/
2128 /*!
2129 * \brief Add an enforcement of the value of degrees of freedom located at
2130 * the mesh vertices.
2131 * The spatial discretization scheme for the given equation has to be
2132 * CDO vertex-based or CDO vertex+cell-based schemes.
2133 *
2134 * One assumes that values are interlaced if eqp->dim > 1
2135 * ref_value or elt_values has to be defined. If both parameters are
2136 * defined, one keeps the definition in elt_values
2137 *
2138 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2139 * \param[in] n_vertices number of vertices to enforce
2140 * \param[in] vertex_ids list of vertices
2141 * \param[in] ref_value default values or ignored (may be NULL)
2142 * \param[in] vtx_values list of associated values, ignored if NULL
2143 *
2144 * \return a pointer to a cs_enforcement_param_t structure
2145 */
2146 /*----------------------------------------------------------------------------*/
2147
2148 cs_enforcement_param_t *
2149 cs_equation_add_vertex_dof_enforcement(cs_equation_param_t *eqp,
2150 cs_lnum_t n_vertices,
2151 const cs_lnum_t vertex_ids[],
2152 const cs_real_t ref_value[],
2153 const cs_real_t vtx_values[]);
2154
2155 /*----------------------------------------------------------------------------*/
2156 /*!
2157 * \brief Add an enforcement of the value of degrees of freedom located at
2158 * the mesh edges.
2159 * The spatial discretization scheme for the given equation has to be
2160 * CDO edge-based schemes.
2161 *
2162 * One assumes that values are interlaced if eqp->dim > 1
2163 * ref_value or elt_values has to be defined. If both parameters are
2164 * defined, one keeps the definition in elt_values
2165 *
2166 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2167 * \param[in] n_edges number of edges to enforce
2168 * \param[in] edge_ids list of edges
2169 * \param[in] ref_value default values or ignored (may be NULL)
2170 * \param[in] edge_values list of associated values, ignored if NULL
2171 *
2172 * \return a pointer to a cs_enforcement_param_t structure
2173 */
2174 /*----------------------------------------------------------------------------*/
2175
2176 cs_enforcement_param_t *
2177 cs_equation_add_edge_dof_enforcement(cs_equation_param_t *eqp,
2178 cs_lnum_t n_edges,
2179 const cs_lnum_t edge_ids[],
2180 const cs_real_t ref_value[],
2181 const cs_real_t edge_values[]);
2182
2183 /*----------------------------------------------------------------------------*/
2184 /*!
2185 * \brief Add an enforcement of the value of degrees of freedom located at
2186 * the mesh faces.
2187 * The spatial discretization scheme for the given equation has to be
2188 * CDO face-based schemes.
2189 *
2190 * One assumes that values are interlaced if eqp->dim > 1
2191 * ref_value or elt_values has to be defined. If both parameters are
2192 * defined, one keeps the definition in elt_values
2193 *
2194 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2195 * \param[in] n_faces number of faces to enforce
2196 * \param[in] face_ids list of faces
2197 * \param[in] ref_value default values or ignored (may be NULL)
2198 * \param[in] face_values list of associated values, ignored if NULL
2199 *
2200 * \return a pointer to a cs_enforcement_param_t structure
2201 */
2202 /*----------------------------------------------------------------------------*/
2203
2204 cs_enforcement_param_t *
2205 cs_equation_add_face_dof_enforcement(cs_equation_param_t *eqp,
2206 cs_lnum_t n_faces,
2207 const cs_lnum_t face_ids[],
2208 const cs_real_t ref_value[],
2209 const cs_real_t face_values[]);
2210
2211 /*----------------------------------------------------------------------------*/
2212 /*!
2213 * \brief Add an enforcement of the value related to the degrees of freedom
2214 * associated to the list of selected cells.
2215 *
2216 * One assumes that values are interlaced if eqp->dim > 1
2217 * ref_value or elt_values has to be defined. If both parameters are
2218 * defined, one keeps the definition in elt_values
2219 *
2220 * \param[in, out] eqp pointer to a cs_equation_param_t structure
2221 * \param[in] n_cells number of selected cells
2222 * \param[in] cell_ids list of cell ids
2223 * \param[in] ref_value ignored if NULL
2224 * \param[in] cell_values list of associated values, ignored if NULL
2225 *
2226 * \return a pointer to a cs_enforcement_param_t structure
2227 */
2228 /*----------------------------------------------------------------------------*/
2229
2230 cs_enforcement_param_t *
2231 cs_equation_add_cell_enforcement(cs_equation_param_t *eqp,
2232 cs_lnum_t n_cells,
2233 const cs_lnum_t elt_ids[],
2234 const cs_real_t ref_value[],
2235 const cs_real_t cell_values[]);
2236
2237 /*----------------------------------------------------------------------------*/
2238
2239 END_C_DECLS
2240
2241 #endif /* __CS_EQUATION_PARAM_H__ */
2242