#include "cs_defs.h"
Go to the source code of this file.
Macros | |
#define | CS_ISOTROPIC_DIFFUSION (1 << 0) |
#define | CS_ORTHOTROPIC_DIFFUSION (1 << 1) |
#define | CS_ANISOTROPIC_LEFT_DIFFUSION (1 << 2) |
#define | CS_ANISOTROPIC_RIGHT_DIFFUSION (1 << 3) |
#define | CS_ANISOTROPIC_DIFFUSION ((1 << 2) + (1 << 3)) |
Typedefs | |
typedef void() | cs_analytic_func_t(cs_real_t time, cs_lnum_t n_elts, const cs_lnum_t *elt_ids, const cs_real_t *coords, bool dense_output, void *input, cs_real_t *retval) |
Generic function pointer for an evaluation relying on an analytic function elt_ids is optional. If not NULL, it enables to access to the coords array with an indirection. The same indirection can be applied to fill retval if dense_output is set to false. More... | |
typedef void() | cs_dof_func_t(cs_lnum_t n_elts, const cs_lnum_t *elt_ids, bool dense_output, void *input, cs_real_t *retval) |
Generic function pointer for computing a quantity at predefined locations such as degrees of freedom (DoF): cells, faces, edges or or vertices. elt_ids is optional. If not NULL, it enables to fill retval partially with this indirection if dense_output is set to false. More... | |
typedef void() | cs_time_func_t(double time, void *input, cs_real_t *retval) |
Function which defines the evolution of a quantity according to the current time and any structure given as a parameter. More... | |
Enumerations | |
enum | cs_param_space_scheme_t { CS_SPACE_SCHEME_LEGACY, CS_SPACE_SCHEME_CDOVB, CS_SPACE_SCHEME_CDOVCB, CS_SPACE_SCHEME_CDOEB, CS_SPACE_SCHEME_CDOFB, CS_SPACE_SCHEME_HHO_P0, CS_SPACE_SCHEME_HHO_P1, CS_SPACE_SCHEME_HHO_P2, CS_SPACE_N_SCHEMES } |
Type of numerical scheme for the discretization in space. More... | |
enum | cs_param_dof_reduction_t { CS_PARAM_REDUCTION_DERHAM, CS_PARAM_REDUCTION_AVERAGE, CS_PARAM_N_REDUCTIONS } |
Numerical settings for time discretization | |
enum | cs_param_time_scheme_t { CS_TIME_SCHEME_STEADY, CS_TIME_SCHEME_EULER_IMPLICIT, CS_TIME_SCHEME_EULER_EXPLICIT, CS_TIME_SCHEME_CRANKNICO, CS_TIME_SCHEME_THETA, CS_TIME_SCHEME_BDF2, CS_TIME_N_SCHEMES } |
typedef void() cs_analytic_func_t(cs_real_t time, cs_lnum_t n_elts, const cs_lnum_t *elt_ids, const cs_real_t *coords, bool dense_output, void *input, cs_real_t *retval) |
Generic function pointer for an evaluation relying on an analytic function elt_ids is optional. If not NULL, it enables to access to the coords array with an indirection. The same indirection can be applied to fill retval if dense_output is set to false.
[in] | time | when ? |
[in] | n_elts | number of elements to consider |
[in] | elt_ids | list of elements ids (in coords and retval) |
[in] | coords | where ? |
[in] | dense_output | perform an indirection in retval or not |
[in] | input | NULL or pointer to a structure cast on-the-fly |
[in,out] | retval | resulting value(s). Must be allocated. |
typedef void() cs_dof_func_t(cs_lnum_t n_elts, const cs_lnum_t *elt_ids, bool dense_output, void *input, cs_real_t *retval) |
Generic function pointer for computing a quantity at predefined locations such as degrees of freedom (DoF): cells, faces, edges or or vertices. elt_ids is optional. If not NULL, it enables to fill retval partially with this indirection if dense_output is set to false.
[in] | n_elts | number of elements to consider |
[in] | elt_ids | list of elements ids |
[in] | dense_output | perform an indirection in retval or not |
[in] | input | NULL or pointer to a structure cast on-the-fly |
[in,out] | retval | resulting value(s). Must be allocated. |
typedef void() cs_time_func_t(double time, void *input, cs_real_t *retval) |
Function which defines the evolution of a quantity according to the current time and any structure given as a parameter.
[in] | time | value of the time at the end of the last iteration |
[in] | input | NULL or pointer to a structure cast on-the-fly |
[in] | retval | result of the evaluation |
Choice of how to extrapolate the advection field in the advection term.
Type of formulation for the advection term
Type of numerical scheme for the advection term
Enumerator | |
---|---|
CS_PARAM_ADVECTION_SCHEME_CENTERED | centered discretization |
CS_PARAM_ADVECTION_SCHEME_CIP | Continuous Interior Penalty discretization. Only available for CS_SPACE_SCHEME_CDOVCB |
CS_PARAM_ADVECTION_SCHEME_CIP_CW | Continuous Interior Penalty discretization. Only available for CS_SPACE_SCHEME_CDOVCB Variant with a cellwise constant approximation of the advection field |
CS_PARAM_ADVECTION_SCHEME_HYBRID_CENTERED_UPWIND | Centered discretization with a portion between [0,1] of upwinding. The portion is specified thanks to CS_EQKEY_ADV_UPWIND_PORTION If the portion is equal to 0, then one recovers CS_PARAM_ADVECTION_SCHEME_CENTERED. If the portion is equal to 1, then one recovers CS_PARAM_ADVECTION_SCHEME_UPWIND |
CS_PARAM_ADVECTION_SCHEME_SAMARSKII | Weighting between an upwind and a centered discretization relying on the Peclet number. Weighting function = Samarskii |
CS_PARAM_ADVECTION_SCHEME_SG | Weighting between an upwind and a centered discretization relying on the Peclet number. Weighting function = Scharfetter-Gummel |
CS_PARAM_ADVECTION_SCHEME_UPWIND | Low order upwind discretization |
CS_PARAM_N_ADVECTION_SCHEMES |
Choice of how to handle the advection term in an equation.
Enumerator | |
---|---|
CS_PARAM_ADVECTION_IMPLICIT_FULL | The advection term is implicitly treated. The non-linearity stemming from this term has to be solved using a specific algorithm such as the Picard (fixed-point) technique or more elaborated techniques. |
CS_PARAM_ADVECTION_IMPLICIT_LINEARIZED | The advection term is implicitly treated. The non-linearity stemming from this term is simplified. Namely, one assumes a linearized advection. This is equivalent to a one-step Picard technique. |
CS_PARAM_ADVECTION_EXPLICIT | The advection term is treated explicitly. One keeps the non-linearity stemming from this term at the right hand-side. An extrapolation can be used for the advection field (cf. cs_param_advection_extrapol_t) |
CS_PARAM_N_ADVECTION_STRATEGIES |
enum cs_param_amg_type_t |
Type of AMG (Algebraic MultiGrid) algorithm to use (either as a preconditioner with or a solver).
Type of method for enforcing the boundary conditions. According to the type of numerical scheme, some enforcements are not available.
Enumerator | |
---|---|
CS_PARAM_BC_ENFORCE_ALGEBRAIC | Weak enforcement of the boundary conditions (i.e. one keeps the degrees of freedom in the algebraic system) with an algebraic manipulation of the linear system |
CS_PARAM_BC_ENFORCE_PENALIZED | Weak enforcement of the boundary conditions (i.e. one keeps the degrees of freedom in the algebraic system) with a penalization technique using a huge value. |
CS_PARAM_BC_ENFORCE_WEAK_NITSCHE | Weak enforcement of the boundary conditions (i.e. one keeps the degrees of freedom in the algebraic system) with a Nitsche-like penalization technique. This technique does not increase the conditioning number as much as CS_PARAM_BC_ENFORCE_PENALIZED but is more computationally intensive. The computation of the diffusion term should be activated with this choice. |
CS_PARAM_BC_ENFORCE_WEAK_SYM | Weak enforcement of the boundary conditions (i.e. one keeps the degrees of freedom in the algebraic system) with a Nitsche-like penalization technique. This variant enables to keep the symmetry of the algebraic contrary to CS_PARAM_BC_ENFORCE_WEAK_NITSCHE. The computation of the diffusion term should be activated with this choice. |
CS_PARAM_N_BC_ENFORCEMENTS |
enum cs_param_bc_type_t |
Type of boundary condition to enforce.
Type of solver to use to solve a linear system. Some of the mentionned solver are available only if the PETSc library is linked with code_saturne.
Enumerator | |
---|---|
CS_PARAM_ITSOL_NONE | No iterative solver (equivalent to a "preonly" choice in PETSc) |
CS_PARAM_ITSOL_AMG | Algebraic multigrid solver (additional options may be set using cs_param_amg_type_t) |
CS_PARAM_ITSOL_BICG | Bi-Conjuguate gradient (useful for non-symmetric systems) |
CS_PARAM_ITSOL_BICGSTAB2 | Stabilized Bi-Conjuguate gradient (useful for non-symmetric systems) |
CS_PARAM_ITSOL_CG | Conjuguate Gradient (solver of choice for symmetric positive definite systems) |
CS_PARAM_ITSOL_CR3 | 3-layer conjugate residual (can handle non-symmetric systems) |
CS_PARAM_ITSOL_FCG | Flexible Conjuguate Gradient (variant of the CG when the preconditioner may change from one iteration to another. For instance when using an AMG as preconditioner) |
CS_PARAM_ITSOL_FGMRES | Only with PETsc Flexible Generalized Minimal RESidual |
CS_PARAM_ITSOL_GAUSS_SEIDEL | Gauss-Seidel |
CS_PARAM_ITSOL_GKB_CG | Golub-Kahan Bidiagonalization algorithm. Useful for solving saddle-point systems. The inner solver is a (flexible) CG solver. |
CS_PARAM_ITSOL_GKB_GMRES | Golub-Kahan Bidiagonalization algorithm. Useful for solving saddle-point systems. The inner solver is a (flexible) GMRES solver. |
CS_PARAM_ITSOL_GMRES | Only with PETsc Generalized Minimal RESidual |
CS_PARAM_ITSOL_JACOBI | Jacobi (diagonal relaxation) |
CS_PARAM_ITSOL_MINRES | Only with PETsc Mininal residual algorithm |
CS_PARAM_ITSOL_MUMPS | Only with PETsc/MUMPS MUMPS direct solver (LU factorization) |
CS_PARAM_ITSOL_MUMPS_LDLT | Only with PETsc/MUMPS MUMPS direct solver (LDLT factorization also known as Cholesky factorization) |
CS_PARAM_ITSOL_SYM_GAUSS_SEIDEL | Symmetric Gauss-Seidel |
CS_PARAM_N_ITSOL_TYPES |
Type of preconditioner to use with the iterative solver. Some of the mentionned preconditioners are available only if the PETSc library is linked with code_saturne
Enumerator | |
---|---|
CS_PARAM_PRECOND_NONE | No preconditioner |
CS_PARAM_PRECOND_BJACOB_ILU0 | Block Jacobi with an ILU zero fill-in in each block |
CS_PARAM_PRECOND_BJACOB_SGS | Block Jacobi with a symmetric Gauss-Seidel in each block (rely on Eisenstat's trick with PETsc) |
CS_PARAM_PRECOND_AMG | Algebraic multigrid preconditioner (additional options may be set using cs_param_amg_type_t) |
CS_PARAM_PRECOND_AMG_BLOCK | Algebraic multigrid preconditioner by block (useful in case of vector valued variables) |
CS_PARAM_PRECOND_AS | Only with PETSc Additive Schwarz preconditioner |
CS_PARAM_PRECOND_DIAG | Diagonal (also Jacobi) preconditioner. The cheapest one but not the most efficient one. |
CS_PARAM_PRECOND_GKB_CG | Only with PETSc Golub-Kahan Bidiagonalization solver used as a preconditioner. Only useful if one has to solve a saddle-point system (such systems arise when solving Stokes or Navier-Stokes in a fully couple manner). Variant with CG as inner solver. |
CS_PARAM_PRECOND_GKB_GMRES | Only with PETSc Golub-Kahan Bidiagonalization solver used as a preconditioner. Only useful if one has to solve a saddle-point system (such systems arise when solving Stokes or Navier-Stokes in a fully couple manner). Variant with GMRES as inner solver. |
CS_PARAM_PRECOND_ILU0 | Only with PETSc Incomplute LU factorization (fill-in coefficient set to 0) |
CS_PARAM_PRECOND_ICC0 | Only with PETSc Incomplute Cholesky factorization (fill-in coefficient set to 0). This is variant of the ILU0 preconditioner dedicated to symmetric positive definite system |
CS_PARAM_PRECOND_POLY1 | Neumann polynomial preconditioning. Polynoms of order 1. |
CS_PARAM_PRECOND_POLY2 | Neumann polynomial preconditioning. Polynoms of order 2. |
CS_PARAM_PRECOND_SSOR | Symmetric Successive OverRelaxations (can be seen as a symmetric Gauss-Seidel preconditioner) |
CS_PARAM_N_PRECOND_TYPES |
Way to renormalize (or not) the residual arising during the resolution of linear systems
Class of iterative solvers to consider for solver the linear system.
Type of numerical scheme for the discretization in space.
How is defined the degree of freedom.
Type of numerical scheme for the discretization in time
const char* cs_param_get_advection_extrapol_name | ( | cs_param_advection_extrapol_t | extrapol | ) |
Get the label associated to the extrapolation used for the advection field.
[in] | adv_stra | type of extrapolation for the advection field |
const char* cs_param_get_advection_form_name | ( | cs_param_advection_form_t | adv_form | ) |
Get the label associated to the advection formulation.
[in] | adv_form | type of advection formulation |
const char* cs_param_get_advection_scheme_name | ( | cs_param_advection_scheme_t | scheme | ) |
Get the label of the advection scheme.
[in] | scheme | type of advection scheme |
const char* cs_param_get_advection_strategy_name | ( | cs_param_advection_strategy_t | adv_stra | ) |
Get the label associated to the advection strategy.
[in] | adv_stra | type of advection strategy |
const char* cs_param_get_amg_type_name | ( | cs_param_amg_type_t | type | ) |
Get the name of the type of algebraic multigrid (AMG)
[in] | type | type of AMG |
const char* cs_param_get_bc_enforcement_name | ( | cs_param_bc_enforce_t | type | ) |
Get the name of the type of enforcement of the boundary condition.
[in] | type | type of enforcement of boundary conditions |
const char* cs_param_get_bc_name | ( | cs_param_bc_type_t | type | ) |
Get the name of the type of boundary condition.
[in] | type | type of boundary condition |
const char* cs_param_get_precond_name | ( | cs_param_precond_type_t | precond | ) |
Get the name of the preconditioner.
[in] | precond | type of preconditioner |
const char* cs_param_get_solver_name | ( | cs_param_itsol_type_t | solver | ) |
Get the name of the solver.
[in] | solver | type of iterative solver |
const char* cs_param_get_space_scheme_name | ( | cs_param_space_scheme_t | scheme | ) |
Get the name of the space discretization scheme.
[in] | scheme | type of space scheme |
const char* cs_param_get_time_scheme_name | ( | cs_param_time_scheme_t | scheme | ) |
Get the name of the time discretization scheme.
[in] | scheme | type of time scheme |
bool cs_param_space_scheme_is_face_based | ( | cs_param_space_scheme_t | scheme | ) |
Return true if the space scheme has degrees of freedom on faces, otherwise false.
[in] | scheme | type of space scheme |
const char cs_med_sepline[50] |
const char cs_sep_h1[80] |
const char cs_sep_h2[80] |
const char cs_sepline[80] |