1 #ifndef __CS_PARAM_TYPES_H__ 2 #define __CS_PARAM_TYPES_H__ 3 4 /*============================================================================ 5 * Manage the definition/setting of a computation 6 *============================================================================*/ 7 8 /* 9 This file is part of Code_Saturne, a general-purpose CFD tool. 10 11 Copyright (C) 1998-2021 EDF S.A. 12 13 This program is free software; you can redistribute it and/or modify it under 14 the terms of the GNU General Public License as published by the Free Software 15 Foundation; either version 2 of the License, or (at your option) any later 16 version. 17 18 This program is distributed in the hope that it will be useful, but WITHOUT 19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS 20 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more 21 details. 22 23 You should have received a copy of the GNU General Public License along with 24 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin 25 Street, Fifth Floor, Boston, MA 02110-1301, USA. 26 */ 27 28 /*----------------------------------------------------------------------------*/ 29 30 /*---------------------------------------------------------------------------- 31 * Local headers 32 *----------------------------------------------------------------------------*/ 33 34 #include "cs_defs.h" 35 36 /*----------------------------------------------------------------------------*/ 37 38 BEGIN_C_DECLS 39 40 /*============================================================================ 41 * Macro definitions 42 *============================================================================*/ 43 44 /*! 45 * @addtogroup scalar_params 46 * 47 * @{ 48 */ 49 50 /* 51 * Field property type 52 */ 53 54 /*! isotropic diffusion */ 55 56 #define CS_ISOTROPIC_DIFFUSION (1 << 0) 57 58 /*! orthotropic diffusion */ 59 60 #define CS_ORTHOTROPIC_DIFFUSION (1 << 1) 61 62 /*! diffusion by a left-multiplied symmetric 3x3 tensor */ 63 64 #define CS_ANISOTROPIC_LEFT_DIFFUSION (1 << 2) 65 66 /*! diffusion by a right-multiplied symmetric 3x3 tensor */ 67 68 #define CS_ANISOTROPIC_RIGHT_DIFFUSION (1 << 3) 69 70 /*! diffusion by a symmetric 3x3 tensor */ 71 72 #define CS_ANISOTROPIC_DIFFUSION ((1 << 2) + (1 << 3)) 73 74 /*! @} */ 75 76 /*============================================================================ 77 * Type definitions 78 *============================================================================*/ 79 80 /*----------------------------------------------------------------------------*/ 81 /*! 82 * \brief Generic function pointer for an evaluation relying on an analytic 83 * function. 84 * 85 * For the calling function, elt_ids is optional. If not NULL, the coords 86 * array should be accessed with an indirection. The same indirection can 87 * be applied to fill retval if dense_output is set to false. 88 * 89 * \param[in] time when ? 90 * \param[in] n_elts number of elements to consider 91 * \param[in] elt_ids list of elements ids (in coords and retval) 92 * \param[in] coords where ? 93 * \param[in] dense_output perform an indirection in retval or not 94 * \param[in] input NULL or pointer to a structure cast on-the-fly 95 * \param[in, out] retval resulting value(s). Must be allocated. 96 */ 97 /*----------------------------------------------------------------------------*/ 98 99 typedef void 100 (cs_analytic_func_t) (cs_real_t time, 101 cs_lnum_t n_elts, 102 const cs_lnum_t *elt_ids, 103 const cs_real_t *coords, 104 bool dense_output, 105 void *input, 106 cs_real_t *retval); 107 108 /*----------------------------------------------------------------------------*/ 109 /*! 110 * \brief Generic function pointer for computing a quantity at predefined 111 * locations such as degrees of freedom (DoF): cells, faces, edges or 112 * or vertices. 113 * 114 * For the calling function, elt_ids is optional. If not NULL, array(s) should 115 * be accessed with an indirection. The same indirection can be applied to fill 116 * retval if dense_output is set to false. 117 * 118 * \param[in] n_elts number of elements to consider 119 * \param[in] elt_ids list of elements ids 120 * \param[in] dense_output perform an indirection in retval or not 121 * \param[in] input NULL or pointer to a structure cast on-the-fly 122 * \param[in, out] retval resulting value(s). Must be allocated. 123 */ 124 /*----------------------------------------------------------------------------*/ 125 126 typedef void 127 (cs_dof_func_t) (cs_lnum_t n_elts, 128 const cs_lnum_t *elt_ids, 129 bool dense_output, 130 void *input, 131 cs_real_t *retval); 132 133 /*----------------------------------------------------------------------------*/ 134 /*! 135 * \brief Function which defines the evolution of a quantity according to the 136 * current time and any structure given as a parameter 137 * 138 * \param[in] time value of the time at the end of the last iteration 139 * \param[in] input NULL or pointer to a structure cast on-the-fly 140 * \param[in] retval result of the evaluation 141 */ 142 /*----------------------------------------------------------------------------*/ 143 144 typedef void 145 (cs_time_func_t) (double time, 146 void *input, 147 cs_real_t *retval); 148 149 /*! \enum cs_param_space_scheme_t 150 * \brief Type of numerical scheme for the discretization in space 151 * 152 * \var CS_SPACE_SCHEME_LEGACY 153 * Cell-centered Finite Volume Two-Point Flux 154 * 155 * \var CS_SPACE_SCHEME_CDOVB 156 * CDO scheme with vertex-based positioning 157 * 158 * \var CS_SPACE_SCHEME_CDOVCB 159 * CDO scheme with vertex+cell-based positioning 160 * 161 * \var CS_SPACE_SCHEME_CDOEB 162 * CDO scheme with edge-based positioning 163 * 164 * \var CS_SPACE_SCHEME_CDOFB 165 * CDO scheme with face-based positioning 166 * 167 * \var CS_SPACE_SCHEME_HHO_P0 168 * Hybrid High Order (HHO) schemes 169 * HHO scheme with face-based positioning (lowest order) 170 * 171 * \var CS_SPACE_SCHEME_HHO_P1 172 * Hybrid High Order (HHO) schemes 173 * HHO scheme with face-based positioning (k=1 up to order 3) 174 * 175 * \var CS_SPACE_SCHEME_HHO_P2 176 * Hybrid High Order (HHO) schemes 177 * HHO scheme with face-based positioning (k=2 up to order 4) 178 */ 179 180 typedef enum { 181 182 CS_SPACE_SCHEME_LEGACY, 183 CS_SPACE_SCHEME_CDOVB, 184 CS_SPACE_SCHEME_CDOVCB, 185 CS_SPACE_SCHEME_CDOEB, 186 CS_SPACE_SCHEME_CDOFB, 187 CS_SPACE_SCHEME_HHO_P0, 188 CS_SPACE_SCHEME_HHO_P1, 189 CS_SPACE_SCHEME_HHO_P2, 190 191 CS_SPACE_N_SCHEMES 192 193 } cs_param_space_scheme_t; 194 195 /*! \enum cs_param_space_scheme_t 196 * \brief How is defined the degree of freedom 197 * 198 * \var CS_PARAM_REDUCTION_DERHAM 199 * Evaluation at vertices for potentials, integral along a line for 200 * circulations, integral across the normal component of a face for fluxes and 201 * integral inside a cell for densities 202 * 203 * \var CS_PARAM_REDUCTION_AVERAGE 204 * Degrees of freedom are defined as the average on the element 205 */ 206 207 typedef enum { 208 209 CS_PARAM_REDUCTION_DERHAM, 210 CS_PARAM_REDUCTION_AVERAGE, 211 212 CS_PARAM_N_REDUCTIONS 213 214 } cs_param_dof_reduction_t; 215 216 /*! 217 * @name Numerical settings for time discretization 218 * @{ 219 * 220 * \enum cs_param_time_scheme_t 221 * Type of numerical scheme for the discretization in time 222 * 223 * \var CS_TIME_SCHEME_STEADY 224 * No time scheme. Steady-state computation. 225 * 226 * \var CS_TIME_SCHEME_EULER_IMPLICIT 227 * fully implicit (forward Euler/theta-scheme = 1) 228 * 229 * \var CS_TIME_SCHEME_EULER_EXPLICIT 230 * fully explicit (backward Euler/theta-scheme = 0) 231 * 232 * \var CS_TIME_SCHEME_CRANKNICO 233 * Crank-Nicolson (theta-scheme = 0.5) 234 * 235 * \var CS_TIME_SCHEME_THETA 236 * theta-scheme 237 * 238 * \var CS_TIME_SCHEME_BDF2 239 * theta-scheme 240 */ 241 242 typedef enum { 243 244 CS_TIME_SCHEME_STEADY, 245 CS_TIME_SCHEME_EULER_IMPLICIT, 246 CS_TIME_SCHEME_EULER_EXPLICIT, 247 CS_TIME_SCHEME_CRANKNICO, 248 CS_TIME_SCHEME_THETA, 249 CS_TIME_SCHEME_BDF2, 250 251 CS_TIME_N_SCHEMES 252 253 } cs_param_time_scheme_t; 254 255 /*! 256 * @} 257 * @name Settings for the advection term 258 * @{ 259 * 260 * \enum cs_param_advection_form_t 261 * Type of formulation for the advection term 262 * 263 * \var CS_PARAM_ADVECTION_FORM_CONSERV 264 * conservative formulation (i.e. divergence of a flux) 265 * 266 * \var CS_PARAM_ADVECTION_FORM_NONCONS 267 * non-conservative formation (i.e advection field times 268 * the gradient) 269 * 270 * \var CS_PARAM_ADVECTION_FORM_SKEWSYM 271 * skew-symmetric form 272 * 273 */ 274 275 typedef enum { 276 277 CS_PARAM_ADVECTION_FORM_CONSERV, 278 CS_PARAM_ADVECTION_FORM_NONCONS, 279 CS_PARAM_ADVECTION_FORM_SKEWSYM, 280 281 CS_PARAM_N_ADVECTION_FORMULATIONS 282 283 } cs_param_advection_form_t; 284 285 /*! \enum cs_param_advection_scheme_t 286 * Type of numerical scheme for the advection term 287 * 288 * \var CS_PARAM_ADVECTION_SCHEME_CENTERED 289 * centered discretization 290 * 291 * \var CS_PARAM_ADVECTION_SCHEME_CIP 292 * Continuous Interior Penalty discretization. Only available for 293 * \ref CS_SPACE_SCHEME_CDOVCB 294 * 295 * \var CS_PARAM_ADVECTION_SCHEME_CIP_CW 296 * Continuous Interior Penalty discretization. Only available for 297 * \ref CS_SPACE_SCHEME_CDOVCB 298 * Variant with a cellwise constant approximation of the advection field 299 * 300 * \var CS_PARAM_ADVECTION_SCHEME_HYBRID_CENTERED_UPWIND 301 * Centered discretization with a portion between [0,1] of upwinding. 302 * The portion is specified thanks to \ref CS_EQKEY_ADV_UPWIND_PORTION 303 * If the portion is equal to 0, then one recovers 304 * \ref CS_PARAM_ADVECTION_SCHEME_CENTERED. If the portion is equal to 1, then 305 * one recovers CS_PARAM_ADVECTION_SCHEME_UPWIND 306 * 307 * \var CS_PARAM_ADVECTION_SCHEME_SAMARSKII 308 * Weighting between an upwind and a centered discretization relying on the 309 * Peclet number. Weighting function = Samarskii 310 * 311 * \var CS_PARAM_ADVECTION_SCHEME_SG 312 * Weighting between an upwind and a centered discretization relying on the 313 * Peclet number. Weighting function = Scharfetter-Gummel 314 * 315 * \var CS_PARAM_ADVECTION_SCHEME_UPWIND 316 * Low order upwind discretization 317 */ 318 319 typedef enum { 320 321 CS_PARAM_ADVECTION_SCHEME_CENTERED, 322 CS_PARAM_ADVECTION_SCHEME_CIP, 323 CS_PARAM_ADVECTION_SCHEME_CIP_CW, 324 CS_PARAM_ADVECTION_SCHEME_HYBRID_CENTERED_UPWIND, 325 CS_PARAM_ADVECTION_SCHEME_SAMARSKII, 326 CS_PARAM_ADVECTION_SCHEME_SG, 327 CS_PARAM_ADVECTION_SCHEME_UPWIND, 328 329 CS_PARAM_N_ADVECTION_SCHEMES 330 331 } cs_param_advection_scheme_t; 332 333 /*! \enum cs_param_advection_strategy_t 334 * \brief Choice of how to handle the advection term in an equation 335 * 336 * \var CS_PARAM_ADVECTION_IMPLICIT_FULL 337 * The advection term is implicitly treated. The non-linearity stemming from 338 * this term has to be solved using a specific algorithm such as the Picard 339 * (fixed-point) technique or more elaborated techniques. 340 * 341 * \var CS_PARAM_ADVECTION_IMPLICIT_LINEARIZED 342 * The advection term is implicitly treated. The non-linearity stemming from 343 * this term is simplified. Namely, one assumes a linearized advection. This is 344 * equivalent to a one-step Picard technique. 345 * 346 * \var CS_PARAM_ADVECTION_EXPLICIT 347 * The advection term is treated explicitly. One keeps the non-linearity 348 * stemming from this term at the right hand-side. An extrapolation can be 349 * used for the advection field (cf. \ref cs_param_advection_extrapol_t) 350 */ 351 352 typedef enum { 353 354 CS_PARAM_ADVECTION_IMPLICIT_FULL, 355 CS_PARAM_ADVECTION_IMPLICIT_LINEARIZED, 356 CS_PARAM_ADVECTION_EXPLICIT, 357 358 CS_PARAM_N_ADVECTION_STRATEGIES 359 360 } cs_param_advection_strategy_t; 361 362 /*! \enum cs_param_advection_extrapol_t 363 * \brief Choice of how to extrapolate the advection field in the advection 364 * term 365 * 366 * \var CS_PARAM_ADVECTION_EXTRAPOL_NONE 367 * The advection field is not extrapolated. The last known advection field 368 * is considered \phi^n if one computes u^(n+1) knowing u^n. In case of a 369 * a non-linearity, this is \phi^(n+1,k) if one computes u^(n+1,k+1) knowing 370 * u^(n+1,k) and u^n. The initial step is such that u^(n+1,0) = u^n 371 * This is a good choice with a first-order forward or backward time scheme. 372 * 373 * \var CS_PARAM_ADVECTION_EXTRAPOL_TAYLOR_2 374 * The advection field is extrapolated with a 2nd order Taylor expansion 375 * yielding \phi^extrap = 2 \phi^n - \phi^(n-1) 376 * This corresponds to an estimation of \phi at n+1. Thus, this is a good 377 * choice when associated with a BDF2 time scheme. 378 * 379 * \var CS_PARAM_ADVECTION_EXTRAPOL_ADAMS_BASHFORTH_2 380 * The advection field is extrapolated with a 2nd order Adams-Bashforth 381 * technique yielding \phi^extrap = 3/2 \phi^n - 1/2 \phi^(n-1) 382 * This corresponds to an estimation of \phi at n+1/2. Thus, this is a good 383 * choice when associated with a Crank-Nilcolson time scheme. 384 */ 385 386 typedef enum { 387 388 CS_PARAM_ADVECTION_EXTRAPOL_NONE, 389 CS_PARAM_ADVECTION_EXTRAPOL_TAYLOR_2, 390 CS_PARAM_ADVECTION_EXTRAPOL_ADAMS_BASHFORTH_2, 391 392 CS_PARAM_N_ADVECTION_EXTRAPOLATIONS 393 394 } cs_param_advection_extrapol_t; 395 396 /*! 397 * @} 398 * @name Settings for the boundary conditions 399 * @{ 400 * 401 * \enum cs_param_bc_type_t 402 * Type of boundary condition to enforce. 403 * 404 * \var CS_PARAM_BC_HMG_DIRICHLET 405 * Homogeneous Dirichlet boundary conditions. The value of a variable is set 406 * to zero. 407 * 408 * \var CS_PARAM_BC_DIRICHLET 409 * Dirichlet boundary conditions. The value of the variable is set to the user 410 * requirements. 411 * 412 * \var CS_PARAM_BC_HMG_NEUMANN 413 * Homogeneous Neumann conditions. The value of the flux of variable is set 414 * to zero. 415 * 416 * \var CS_PARAM_BC_NEUMANN 417 * Neumann conditions (along the face normal). The value of the flux of 418 * variable is set to the user requirements. 419 * 420 * \var CS_PARAM_BC_NEUMANN_FULL 421 * Neumann conditions with full definition relative to all directions, i.e. 422 * not only the normal direction. The value of the flux of variable is set to 423 * the user requirements. 424 * 425 * \var CS_PARAM_BC_ROBIN 426 * Robin conditions. 427 * 428 * \var CS_PARAM_BC_SLIDING 429 * Sliding conditions. Homogeneous Dirichlet for the normal component and 430 * homogeneous Neumann for the tangential components. Only available for 431 * vector-valued equations. 432 * 433 * \var CS_PARAM_BC_CIRCULATION 434 * Set the tangential part of a vector-valued field. This is the part lying on 435 * the boundary of a part of the computatinal domain. Nothing is prescribed for 436 * the normal part of the vector-valued field. 437 */ 438 439 typedef enum { 440 441 CS_PARAM_BC_HMG_DIRICHLET, 442 CS_PARAM_BC_DIRICHLET, 443 CS_PARAM_BC_HMG_NEUMANN, 444 CS_PARAM_BC_NEUMANN, 445 CS_PARAM_BC_NEUMANN_FULL, 446 CS_PARAM_BC_ROBIN, 447 CS_PARAM_BC_SLIDING, 448 CS_PARAM_BC_CIRCULATION, 449 450 CS_PARAM_N_BC_TYPES 451 452 } cs_param_bc_type_t; 453 454 /*! \enum cs_param_bc_enforce_t 455 * Type of method for enforcing the boundary conditions. According to the type 456 * of numerical scheme, some enforcements are not available. 457 * 458 * \var CS_PARAM_BC_ENFORCE_ALGEBRAIC 459 * Weak enforcement of the boundary conditions (i.e. one keeps the degrees of 460 * freedom in the algebraic system) with an algebraic manipulation of the 461 * linear system 462 * 463 * \var CS_PARAM_BC_ENFORCE_PENALIZED 464 * Weak enforcement of the boundary conditions (i.e. one keeps the degrees of 465 * freedom in the algebraic system) with a penalization technique using a huge 466 * value. 467 * 468 * \var CS_PARAM_BC_ENFORCE_WEAK_NITSCHE 469 * Weak enforcement of the boundary conditions (i.e. one keeps the degrees of 470 * freedom in the algebraic system) with a Nitsche-like penalization technique. 471 * This technique does not increase the conditioning number as much as 472 * \ref CS_PARAM_BC_ENFORCE_PENALIZED but is more computationally intensive. 473 * The computation of the diffusion term should be activated with this choice. 474 * 475 * \var CS_PARAM_BC_ENFORCE_WEAK_SYM 476 * Weak enforcement of the boundary conditions (i.e. one keeps the degrees of 477 * freedom in the algebraic system) with a Nitsche-like penalization technique. 478 * This variant enables to keep the symmetry of the algebraic contrary to 479 * \ref CS_PARAM_BC_ENFORCE_WEAK_NITSCHE. 480 * The computation of the diffusion term should be activated with this choice. 481 */ 482 483 typedef enum { 484 485 CS_PARAM_BC_ENFORCE_ALGEBRAIC, 486 CS_PARAM_BC_ENFORCE_PENALIZED, 487 CS_PARAM_BC_ENFORCE_WEAK_NITSCHE, 488 CS_PARAM_BC_ENFORCE_WEAK_SYM, 489 490 CS_PARAM_N_BC_ENFORCEMENTS 491 492 } cs_param_bc_enforce_t; 493 494 /*! 495 * @} 496 * @name Settings for non-linear algorithms 497 * @{ 498 * 499 * \enum cs_param_nl_algo_t 500 * \brief Class of non-linear iterative algorithm 501 * 502 * \var CS_PARAM_NL_ALGO_PICARD 503 * Picard (also called fixed point) algorithm 504 * 505 * \var CS_PARAM_NL_ALGO_ANDERSON 506 * Anderson acceleration 507 * 508 * \var CS_PARAM_N_NL_ALGOS 509 */ 510 511 typedef enum { 512 513 CS_PARAM_NL_ALGO_PICARD, 514 CS_PARAM_NL_ALGO_ANDERSON, 515 CS_PARAM_N_NL_ALGOS 516 517 } cs_param_nl_algo_t; 518 519 /*! 520 * @} 521 * @name Settings for the linear solvers or SLES (Sparse Linear Equation Solver) 522 * @{ 523 * 524 * \enum cs_param_sles_class_t 525 * \brief Class of iterative solvers to consider for solver the linear system 526 * 527 * \var CS_PARAM_SLES_CLASS_CS 528 * Iterative solvers available in Code_Saturne 529 * 530 * \var CS_PARAM_SLES_CLASS_HYPRE 531 * Solvers available in HYPRE through the PETSc library 532 * 533 * \var CS_PARAM_SLES_CLASS_MUMPS 534 * Solvers available with MUMPS (without the PETSc interface) 535 * 536 * \var CS_PARAM_SLES_CLASS_PETSC 537 * Solvers available in PETSc. Please notice that 538 * the MUMPS solver can be handled within PETSc if the installation of PETSc 539 * includes the MUMPS library 540 * 541 * \var CS_PARAM_SLES_N_CLASSES 542 */ 543 544 typedef enum { 545 546 CS_PARAM_SLES_CLASS_CS, 547 CS_PARAM_SLES_CLASS_HYPRE, 548 CS_PARAM_SLES_CLASS_MUMPS, 549 CS_PARAM_SLES_CLASS_PETSC, 550 551 CS_PARAM_SLES_N_CLASSES 552 553 } cs_param_sles_class_t; 554 555 /*! 556 * \enum cs_param_amg_type_t 557 * Type of AMG (Algebraic MultiGrid) algorithm to use (either as a 558 * preconditioner with or a solver). There are different choices of 559 * implementation and of type of cycle 560 */ 561 562 typedef enum { 563 564 CS_PARAM_AMG_NONE, /*!< No specified algorithm */ 565 CS_PARAM_AMG_HYPRE_BOOMER_V, /*!< V-cycle Boomer algorithm (Hypre lib.) */ 566 CS_PARAM_AMG_HYPRE_BOOMER_W, /*!< W-cycle Boomer algorithm (Hypre lib.) */ 567 CS_PARAM_AMG_PETSC_GAMG_V , /*!< V-cycle GAMG algorithm (PETSc lib.) */ 568 CS_PARAM_AMG_PETSC_GAMG_W , /*!< W-cycle GAMG algorithm (PETSc lib.) */ 569 CS_PARAM_AMG_PETSC_PCMG , /*!< preconditioned MG algorithm from PETSc */ 570 CS_PARAM_AMG_HOUSE_V, /*!< In-house algorithm with V-cycle */ 571 CS_PARAM_AMG_HOUSE_K, /*!< In-house algorithm with K-cycle */ 572 573 CS_PARAM_N_AMG_TYPES 574 575 } cs_param_amg_type_t; 576 577 /*! \enum cs_param_schur_approx_t 578 * 579 * \brief Strategy to build the Schur complement approximation. This appears 580 * in block preconditioning or uzawa algorithms when a monolithic (fully 581 * coupled approach) is used. 582 * | A B^t| 583 * | B 0 | 584 * 585 * \var CS_PARAM_SCHUR_NONE 586 * There is no schur complement approximation. 587 * 588 * \var CS_PARAM_SCHUR_DIAG_INVERSE 589 * The schur complement approximation is defined as B.diag(A)^-1.B^t 590 * 591 * \var CS_PARAM_SCHUR_ELMAN 592 * The inverse of the schur complement matrix is approximated by 593 * (BBt)^-1 B.A.B^t (B.Bt)^-1 594 * This formulation is detailed in Elman'99, SIAM J. SCI. COMPUT. 595 * 596 * \var CS_PARAM_SCHUR_IDENTITY 597 * The schur complement approximation is simply the identity matrix 598 * 599 * \var CS_PARAM_SCHUR_LUMPED_INVERSE 600 * The schur complement approximation is defined as B.lumped(A^-1).B^t where 601 * x=lumped(A^-1) results from A.x = 1 (1 is the array fills with 1 in each 602 * entry) 603 * 604 * \var CS_PARAM_SCHUR_MASS_SCALED 605 * The schur complement approximation is simply a scaled diagonal mass matrix 606 * related to the 22 block 607 * 608 * \var CS_PARAM_SCHUR_MASS_SCALED_DIAG_INVERSE 609 * The schur complement approximation is defined as 610 * S \approx alpha.M22 + 1/dt*B.diag(A)^-1.B^t 611 * where M22 is the mass matrix related to the (2,2) block 612 * 613 * \var CS_PARAM_SCHUR_MASS_SCALED_LUMPED_INVERSE 614 * The schur complement approximation is defined as 615 * S \approx alpha.M22 + 1/dt*B.lumped(A^-1).B^t 616 * where M22 is the mass matrix related to the (2,2) block and where 617 * x=lumped(A^-1) results from A.x = 1 (1 is the array fills with 1 in each 618 * entry) 619 */ 620 621 typedef enum { 622 623 CS_PARAM_SCHUR_NONE, 624 625 CS_PARAM_SCHUR_DIAG_INVERSE, 626 CS_PARAM_SCHUR_ELMAN, 627 CS_PARAM_SCHUR_IDENTITY, 628 CS_PARAM_SCHUR_LUMPED_INVERSE, 629 CS_PARAM_SCHUR_MASS_SCALED, 630 CS_PARAM_SCHUR_MASS_SCALED_DIAG_INVERSE, 631 CS_PARAM_SCHUR_MASS_SCALED_LUMPED_INVERSE, 632 633 CS_PARAM_N_SCHUR_APPROX 634 635 } cs_param_schur_approx_t; 636 637 /*! 638 * \enum cs_param_precond_block_t 639 * Type of preconditioning by block 640 * 641 * \var CS_PARAM_PRECOND_BLOCK_NONE 642 * No block preconditioner is requested (default) 643 * 644 * \var CS_PARAM_PRECOND_BLOCK_DIAG 645 * Only the diagonal blocks are considered in the preconditioner 646 * 647 * \var CS_PARAM_PRECOND_BLOCK_FULL_DIAG 648 * Only the diagonal blocks are considered in the preconditioner. 649 * All possible blocks are considered. With simple cases, this is equivalent to 650 * \ref CS_PARAM_PRECOND_BLOCK_DIAG 651 * 652 * \var CS_PARAM_PRECOND_BLOCK_FULL_LOWER_TRIANGULAR 653 * The diagonal blocks and the lower blocks are considered in the 654 * preconditioner. 655 * All possible blocks are considered. With simple cases, this is equivalent to 656 * \ref CS_PARAM_PRECOND_BLOCK_FULL_LOWER_TRIANGULAR 657 * 658 * \var CS_PARAM_PRECOND_BLOCK_FULL_SYM_GAUSS_SEIDEL 659 * A symmetric Gauss-Seidel block preconditioning is considered 660 * (cf. Y. Notay, "A new analysis of block preconditioners for saddle-point 661 * problems" (2014), SIAM J. Matrix. Anal. Appl.) 662 * All possible blocks are considered. With simple cases, this is equivalent to 663 * \ref CS_PARAM_PRECOND_BLOCK_FULL_SYM_GAUSS_SEIDEL 664 * 665 * \var CS_PARAM_PRECOND_BLOCK_FULL_UPPER_TRIANGULAR 666 * The diagonal blocks and the upper blocks are considered in the 667 * preconditioner 668 * All possible blocks are considered. With simple cases, this is equivalent to 669 * \ref CS_PARAM_PRECOND_BLOCK_FULL_UPPER_TRIANGULAR 670 * 671 * \var CS_PARAM_PRECOND_BLOCK_LOWER_TRIANGULAR 672 * The diagonal blocks and the lower blocks are considered in the 673 * preconditioner 674 * 675 * \var CS_PARAM_PRECOND_BLOCK_SYM_GAUSS_SEIDEL 676 * A symmetric Gauss-Seidel block preconditioning is considered 677 * (cf. Y. Notay, "A new analysis of block preconditioners for saddle-point 678 * problems" (2014), SIAM J. Matrix. Anal. Appl.) 679 * 680 * \var CS_PARAM_PRECOND_BLOCK_UPPER_TRIANGULAR 681 * The diagonal blocks and the upper blocks are considered in the 682 * preconditioner 683 * 684 * \var CS_PARAM_PRECOND_BLOCK_UZAWA 685 * One iteration of block Uzawa is performed as preconditioner 686 */ 687 688 typedef enum { 689 690 CS_PARAM_PRECOND_BLOCK_NONE, 691 CS_PARAM_PRECOND_BLOCK_DIAG, 692 CS_PARAM_PRECOND_BLOCK_FULL_DIAG, 693 CS_PARAM_PRECOND_BLOCK_FULL_LOWER_TRIANGULAR, 694 CS_PARAM_PRECOND_BLOCK_FULL_SYM_GAUSS_SEIDEL, 695 CS_PARAM_PRECOND_BLOCK_FULL_UPPER_TRIANGULAR, 696 CS_PARAM_PRECOND_BLOCK_LOWER_TRIANGULAR, 697 CS_PARAM_PRECOND_BLOCK_SYM_GAUSS_SEIDEL, 698 CS_PARAM_PRECOND_BLOCK_UPPER_TRIANGULAR, 699 CS_PARAM_PRECOND_BLOCK_UZAWA, 700 701 CS_PARAM_N_PCD_BLOCK_TYPES 702 703 } cs_param_precond_block_t; 704 705 /*! 706 * \enum cs_param_precond_type_t 707 * Type of preconditioner to use with the iterative solver. 708 * Some of the mentionned preconditioners are available only if the PETSc 709 * library is linked with code_saturne 710 * 711 * \var CS_PARAM_PRECOND_NONE 712 * No preconditioner 713 * 714 * \var CS_PARAM_PRECOND_BJACOB_ILU0 715 * Block Jacobi with an ILU zero fill-in in each block 716 * 717 * \var CS_PARAM_PRECOND_BJACOB_SGS 718 * Block Jacobi with a symmetric Gauss-Seidel in each block (rely on 719 * Eisenstat's trick with PETSc) 720 * 721 * \var CS_PARAM_PRECOND_AMG 722 * Algebraic multigrid preconditioner (additional options may be set using 723 * \ref cs_param_amg_type_t) 724 * 725 * \var CS_PARAM_PRECOND_DIAG 726 * Diagonal (also Jacobi) preconditioner. The cheapest one but not the most 727 * efficient one. 728 * 729 * \var CS_PARAM_PRECOND_GKB_CG 730 * Golub-Kahan Bidiagonalization solver used as a preconditioner. Only 731 * useful if one has to solve a saddle-point system (such systems arise 732 * when solving Stokes or Navier-Stokes in a fully couple manner). 733 * Variant with CG as inner solver. 734 * 735 * \var CS_PARAM_PRECOND_GKB_GMRES 736 * Golub-Kahan Bidiagonalization solver used as a preconditioner. Only 737 * useful if one has to solve a saddle-point system (such systems arise 738 * when solving Stokes or Navier-Stokes in a fully couple manner). 739 * Variant with GMRES as inner solver. 740 * 741 * \var CS_PARAM_PRECOND_LU 742 * LU factorization (direct solver) 743 * 744 * \var CS_PARAM_PRECOND_ILU0 745 * Incomplete LU factorization (fill-in coefficient set to 0) 746 * 747 * \var CS_PARAM_PRECOND_ICC0 748 * Incomplete Cholesky factorization (fill-in coefficient set to 0). This is 749 * variant of the ILU0 preconditioner dedicated to symmetric positive definite 750 * system 751 * 752 * \var CS_PARAM_PRECOND_POLY1 753 * Neumann polynomial preconditioning. Polynoms of order 1. 754 * 755 * \var CS_PARAM_PRECOND_POLY2 756 * Neumann polynomial preconditioning. Polynoms of order 2. 757 * 758 * \var CS_PARAM_PRECOND_SSOR 759 * Symmetric Successive OverRelaxations (can be seen as a symmetric 760 * Gauss-Seidel preconditioner) 761 */ 762 763 typedef enum { 764 765 CS_PARAM_PRECOND_NONE, 766 767 CS_PARAM_PRECOND_BJACOB_ILU0, /*!< Only with PETSc */ 768 CS_PARAM_PRECOND_BJACOB_SGS, /*!< Only with PETSc */ 769 CS_PARAM_PRECOND_AMG, 770 CS_PARAM_PRECOND_DIAG, 771 CS_PARAM_PRECOND_GKB_CG, /*!< Only with PETSc */ 772 CS_PARAM_PRECOND_GKB_GMRES, /*!< Only with PETSc */ 773 CS_PARAM_PRECOND_LU, /*!< Only with PETSc */ 774 CS_PARAM_PRECOND_ILU0, /*!< Only with PETSc */ 775 CS_PARAM_PRECOND_ICC0, /*!< Only with PETSc*/ 776 CS_PARAM_PRECOND_POLY1, 777 CS_PARAM_PRECOND_POLY2, 778 CS_PARAM_PRECOND_SSOR, 779 780 CS_PARAM_N_PRECOND_TYPES 781 782 } cs_param_precond_type_t; 783 784 /*! 785 * \enum cs_param_itsol_type_t 786 * Type of solver to use to solve a linear system. 787 * Some of the mentionned solver are available only if the PETSc library is 788 * linked with code_saturne. 789 * 790 * \var CS_PARAM_ITSOL_NONE 791 * No iterative solver (equivalent to a "preonly" choice in PETSc) 792 * 793 * \var CS_PARAM_ITSOL_AMG 794 * Algebraic multigrid solver (additional options may be set using 795 * \ref cs_param_amg_type_t) 796 * 797 * \var CS_PARAM_ITSOL_BICG 798 * Bi-Conjuguate gradient (useful for non-symmetric systems) 799 * 800 * \var CS_PARAM_ITSOL_BICGSTAB2 801 * Stabilized Bi-Conjuguate gradient (useful for non-symmetric systems) 802 * 803 * \var CS_PARAM_ITSOL_CG 804 * Conjuguate Gradient (solver of choice for symmetric positive definite 805 * systems) 806 * 807 * \var CS_PARAM_ITSOL_CR3 808 * 3-layer conjugate residual (can handle non-symmetric systems) 809 * 810 * \var CS_PARAM_ITSOL_FCG 811 * Flexible Conjuguate Gradient (variant of the CG when the preconditioner 812 * may change from one iteration to another. For instance when using an AMG 813 * as preconditioner) 814 * 815 * \var CS_PARAM_ITSOL_FGMRES 816 * Flexible Generalized Minimal RESidual 817 * 818 * \var CS_PARAM_ITSOL_GAUSS_SEIDEL 819 * Gauss-Seidel 820 * 821 * \var CS_PARAM_ITSOL_GCR 822 * Generalized conjugate residual (flexible iterative solver for symmetric or 823 * non-symmetric system) 824 * 825 * \var CS_PARAM_ITSOL_GKB_CG 826 * Golub-Kahan Bidiagonalization algorithm. Useful for solving saddle-point 827 * systems. The inner solver is a (flexible) CG solver. 828 * 829 * \var CS_PARAM_ITSOL_GKB_GMRES 830 * Golub-Kahan Bidiagonalization algorithm. Useful for solving saddle-point 831 * systems. The inner solver is a (flexible) GMRES solver. 832 * 833 * \var CS_PARAM_ITSOL_GMRES 834 * Generalized Minimal RESidual 835 * 836 * \var CS_PARAM_ITSOL_JACOBI 837 * Jacobi (diagonal relaxation) 838 * 839 * \var CS_PARAM_ITSOL_MINRES 840 * Mininal residual algorithm 841 * 842 * \var CS_PARAM_ITSOL_MUMPS 843 * MUMPS direct solver (LU factorization) 844 * 845 * \var CS_PARAM_ITSOL_MUMPS_FLOAT 846 * MUMPS direct solver (LU factorization with float instead of double). This 847 * enables a saving of the memory consumption which is acceptable if one 848 * considers a preconditioner 849 * 850 * \var CS_PARAM_ITSOL_MUMPS_FLOAT_LDLT 851 * MUMPS direct solver (LDLT factorization also known as Cholesky 852 * factorization). This enables a saving of the memory consumption which is 853 * acceptable if one considers a preconditioner 854 * 855 * \var CS_PARAM_ITSOL_MUMPS_LDLT 856 * MUMPS direct solver (LDLT factorization also known as Cholesky factorization) 857 * 858 * \var CS_PARAM_ITSOL_SYM_GAUSS_SEIDEL 859 * Symmetric Gauss-Seidel 860 * 861 * \var CS_PARAM_ITSOL_USER_DEFINED 862 * User-defined iterative solver. It relies on the implementation of the 863 * the function cs_user_sles_it_solver() 864 */ 865 866 typedef enum { 867 868 CS_PARAM_ITSOL_NONE, 869 870 CS_PARAM_ITSOL_AMG, 871 CS_PARAM_ITSOL_BICG, 872 CS_PARAM_ITSOL_BICGSTAB2, 873 CS_PARAM_ITSOL_CG, 874 CS_PARAM_ITSOL_CR3, 875 CS_PARAM_ITSOL_FCG, 876 CS_PARAM_ITSOL_FGMRES, /*!< Only with PETSc */ 877 CS_PARAM_ITSOL_GAUSS_SEIDEL, 878 CS_PARAM_ITSOL_GCR, 879 CS_PARAM_ITSOL_GKB_CG, 880 CS_PARAM_ITSOL_GKB_GMRES, 881 CS_PARAM_ITSOL_GMRES, /*!< Only with PETSc */ 882 CS_PARAM_ITSOL_JACOBI, 883 CS_PARAM_ITSOL_MINRES, /*!< Only with PETSc */ 884 CS_PARAM_ITSOL_MUMPS, /*!< Only with PETSc or MUMPS */ 885 CS_PARAM_ITSOL_MUMPS_FLOAT, /*!< Only with MUMPS */ 886 CS_PARAM_ITSOL_MUMPS_FLOAT_LDLT, /*!< Only with MUMPS */ 887 CS_PARAM_ITSOL_MUMPS_LDLT, /*!< Only with PETSc or MUMPS */ 888 CS_PARAM_ITSOL_SYM_GAUSS_SEIDEL, 889 CS_PARAM_ITSOL_USER_DEFINED, 890 891 CS_PARAM_N_ITSOL_TYPES 892 893 } cs_param_itsol_type_t; 894 895 /*! 896 * \enum cs_param_resnorm_type_t 897 * Way to renormalize (or not) the residual arising during the resolution of 898 * linear systems 899 */ 900 901 typedef enum { 902 903 CS_PARAM_RESNORM_NONE, /*!< No renormalization */ 904 CS_PARAM_RESNORM_NORM2_RHS, /*!< Renormalization based on the Euclidean 905 norm of the right-hand side */ 906 CS_PARAM_RESNORM_WEIGHTED_RHS, /*!< Renormalization based on a weighted 907 Euclidean norm of the right-hand side */ 908 CS_PARAM_RESNORM_FILTERED_RHS, /*!< Renormalization based on an Euclidean 909 norm of a selection of the right-hand 910 side (penalized terms are filtered) */ 911 CS_PARAM_N_RESNORM_TYPES 912 913 } cs_param_resnorm_type_t; 914 915 /*! 916 * \enum cs_param_dotprod_type_t 917 * Way to compute a dot product or its related norm 918 */ 919 920 typedef enum { 921 922 CS_PARAM_DOTPROD_EUCLIDEAN, /*!< Usual Eucliean 2-norm (no weight) */ 923 CS_PARAM_DOTPROD_CDO, /*!< Weighted 2-norm from CDO quantities */ 924 925 CS_PARAM_N_DOTPROD_TYPES 926 927 } cs_param_dotprod_type_t; 928 929 /*============================================================================ 930 * Global variables 931 *============================================================================*/ 932 933 /* Separation lines: header1, header2 (compatible with markdown), other */ 934 935 extern const char cs_sep_h1[80]; 936 extern const char cs_sep_h2[80]; 937 extern const char cs_sepline[80]; 938 extern const char cs_med_sepline[50]; 939 940 /*============================================================================ 941 * Public function prototypes 942 *============================================================================*/ 943 944 /*----------------------------------------------------------------------------*/ 945 /*! 946 * \brief Return true if the space scheme has degrees of freedom on faces, 947 * otherwise false 948 * 949 * \param[in] scheme type of space scheme 950 * 951 * \return true or false 952 */ 953 /*----------------------------------------------------------------------------*/ 954 955 bool 956 cs_param_space_scheme_is_face_based(cs_param_space_scheme_t scheme); 957 958 /*----------------------------------------------------------------------------*/ 959 /*! 960 * \brief Get the name of the space discretization scheme 961 * 962 * \param[in] scheme type of space scheme 963 * 964 * \return the associated space scheme name 965 */ 966 /*----------------------------------------------------------------------------*/ 967 968 const char * 969 cs_param_get_space_scheme_name(cs_param_space_scheme_t scheme); 970 971 /*----------------------------------------------------------------------------*/ 972 /*! 973 * \brief Get the name of the time discretization scheme 974 * 975 * \param[in] scheme type of time scheme 976 * 977 * \return the associated time scheme name 978 */ 979 /*----------------------------------------------------------------------------*/ 980 981 const char * 982 cs_param_get_time_scheme_name(cs_param_time_scheme_t scheme); 983 984 /*----------------------------------------------------------------------------*/ 985 /*! 986 * \brief Get the label associated to the advection formulation 987 * 988 * \param[in] adv_form type of advection formulation 989 * 990 * \return the associated label 991 */ 992 /*----------------------------------------------------------------------------*/ 993 994 const char * 995 cs_param_get_advection_form_name(cs_param_advection_form_t adv_form); 996 997 /*----------------------------------------------------------------------------*/ 998 /*! 999 * \brief Get the label of the advection scheme 1000 * 1001 * \param[in] scheme type of advection scheme 1002 * 1003 * \return the associated advection scheme label 1004 */ 1005 /*----------------------------------------------------------------------------*/ 1006 1007 const char * 1008 cs_param_get_advection_scheme_name(cs_param_advection_scheme_t scheme); 1009 1010 /*----------------------------------------------------------------------------*/ 1011 /*! 1012 * \brief Get the label associated to the advection strategy 1013 * 1014 * \param[in] adv_stra type of advection strategy 1015 * 1016 * \return the associated label 1017 */ 1018 /*----------------------------------------------------------------------------*/ 1019 1020 const char * 1021 cs_param_get_advection_strategy_name(cs_param_advection_strategy_t adv_stra); 1022 1023 /*----------------------------------------------------------------------------*/ 1024 /*! 1025 * \brief Get the label associated to the extrapolation used for the advection 1026 * field 1027 * 1028 * \param[in] adv_stra type of extrapolation for the advection field 1029 * 1030 * \return the associated label 1031 */ 1032 /*----------------------------------------------------------------------------*/ 1033 1034 const char * 1035 cs_param_get_advection_extrapol_name(cs_param_advection_extrapol_t extrapol); 1036 1037 /*----------------------------------------------------------------------------*/ 1038 /*! 1039 * \brief Get the name of the type of boundary condition 1040 * 1041 * \param[in] type type of boundary condition 1042 * 1043 * \return the associated bc name 1044 */ 1045 /*----------------------------------------------------------------------------*/ 1046 1047 const char * 1048 cs_param_get_bc_name(cs_param_bc_type_t bc); 1049 1050 /*----------------------------------------------------------------------------*/ 1051 /*! 1052 * \brief Get the name of the type of enforcement of the boundary condition 1053 * 1054 * \param[in] type type of enforcement of boundary conditions 1055 * 1056 * \return the associated name 1057 */ 1058 /*----------------------------------------------------------------------------*/ 1059 1060 const char * 1061 cs_param_get_bc_enforcement_name(cs_param_bc_enforce_t type); 1062 1063 /*----------------------------------------------------------------------------*/ 1064 /*! 1065 * \brief Get the name of the non-linear algorithm 1066 * 1067 * \param[in] algo type of algorithm 1068 * 1069 * \return the associated algorithm name 1070 */ 1071 /*----------------------------------------------------------------------------*/ 1072 1073 const char * 1074 cs_param_get_nl_algo_name(cs_param_nl_algo_t algo); 1075 1076 /*----------------------------------------------------------------------------*/ 1077 /*! 1078 * \brief Get the name of the type of dot product to apply 1079 * 1080 * \param[in] dp_type type of dot product 1081 * 1082 * \return the associated type name 1083 */ 1084 /*----------------------------------------------------------------------------*/ 1085 1086 const char * 1087 cs_param_get_dotprod_type_name(cs_param_dotprod_type_t dp_type); 1088 1089 /*----------------------------------------------------------------------------*/ 1090 /*! 1091 * \brief Get the name of the solver 1092 * 1093 * \param[in] solver type of iterative solver 1094 * 1095 * \return the associated solver name 1096 */ 1097 /*----------------------------------------------------------------------------*/ 1098 1099 const char * 1100 cs_param_get_solver_name(cs_param_itsol_type_t solver); 1101 1102 /*----------------------------------------------------------------------------*/ 1103 /*! 1104 * \brief Get the name of the preconditioner 1105 * 1106 * \param[in] precond type of preconditioner 1107 * 1108 * \return the associated preconditioner name 1109 */ 1110 /*----------------------------------------------------------------------------*/ 1111 1112 const char * 1113 cs_param_get_precond_name(cs_param_precond_type_t precond); 1114 1115 /*----------------------------------------------------------------------------*/ 1116 /*! 1117 * \brief Get the name of the type of block preconditioning 1118 * 1119 * \param[in] type type of block preconditioning 1120 * 1121 * \return the associated type name 1122 */ 1123 /*----------------------------------------------------------------------------*/ 1124 1125 const char * 1126 cs_param_get_precond_block_name(cs_param_precond_block_t type); 1127 1128 /*----------------------------------------------------------------------------*/ 1129 /*! 1130 * \brief Get the name of the type of Schur complement approximation 1131 * 1132 * \param[in] type type of Schur complement approximation 1133 * 1134 * \return the associated type name 1135 */ 1136 /*----------------------------------------------------------------------------*/ 1137 1138 const char * 1139 cs_param_get_schur_approx_name(cs_param_schur_approx_t type); 1140 1141 /*----------------------------------------------------------------------------*/ 1142 /*! 1143 * \brief Get the name of the type of algebraic multigrid (AMG) 1144 * 1145 * \param[in] type type of AMG 1146 * 1147 * \return the associated type name 1148 */ 1149 /*----------------------------------------------------------------------------*/ 1150 1151 const char * 1152 cs_param_get_amg_type_name(cs_param_amg_type_t type); 1153 1154 /*----------------------------------------------------------------------------*/ 1155 1156 END_C_DECLS 1157 1158 #endif /* __CS_PARAM_TYPES_H__ */ 1159