1 #ifndef __CS_GWF_H__ 2 #define __CS_GWF_H__ 3 4 /*============================================================================ 5 * Set of main functions to handle the groundwater flow module with CDO 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_base.h" 35 #include "cs_equation.h" 36 #include "cs_gwf_soil.h" 37 #include "cs_gwf_tracer.h" 38 39 /*----------------------------------------------------------------------------*/ 40 41 BEGIN_C_DECLS 42 43 /*============================================================================ 44 * Macro definitions 45 *============================================================================*/ 46 47 /*! 48 * \enum cs_gwf_model_type_t 49 * \brief Type of system of equation(s) to consider for the physical modelling 50 */ 51 52 typedef enum { 53 54 /*! 55 * \brief Single phase (liquid phase) modelling in a porous media. 56 * 57 * All soils are assumed to be saturated. This yields several simplifications 58 * in the Richards equation governing the water conservation. The Richards 59 * equation is steady. The saturation is constant and there is no relative 60 * permeability. 61 */ 62 63 CS_GWF_MODEL_SATURATED_SINGLE_PHASE, 64 65 /*! 66 * \brief Single phase (liquid phase) modelling in a porous media. 67 * 68 * Some soils are not saturated and are described by a more complex model 69 * such as the Van Genuchten-Mualen model. Simplifications made in the case 70 * of \ref CS_GWF_MODEL_SATURATED_SINGLE_PHASE do not hold anymore. Richards 71 * equation is unsteady and there may be a non-linearity to handle according 72 * to the type of soil model. Soil properties such as permeability, soil 73 * capacity and liquid saturation (also called moisture content) are neither 74 * uniform nor steady. 75 */ 76 77 CS_GWF_MODEL_UNSATURATED_SINGLE_PHASE, 78 79 /*! 80 * \brief Two phase flow modelling (gaz and liquid phases) in porous media. 81 * 82 * A Richards-like equation is considered in each phase to take into account 83 * the mass conservation of water and one other component. The component can 84 * be disolved in the liquid phase. No water vapour is taken into 85 * account. Please refer to \ref cs_gwf_miscible_two_phase_t for more details. 86 */ 87 88 CS_GWF_MODEL_TWO_PHASE, 89 90 CS_GWF_N_MODEL_TYPES /*!< Number of predefined models (not a model) */ 91 92 } cs_gwf_model_type_t; 93 94 typedef cs_flag_t cs_gwf_option_flag_t; 95 96 /*! 97 * \enum cs_gwf_model_bit_t 98 * \brief Elemental modelling choice either from the physical viewpoint or the 99 * numerical viewpoint 100 */ 101 102 typedef enum { 103 104 /* -------- Main physical modelling */ 105 106 /*! 107 * \brief Gravitation effects are taken into account in the Richards equation 108 */ 109 110 CS_GWF_GRAVITATION = 1<< 0, /* = 1 */ 111 112 /* --------- Main numerical options */ 113 114 /*! 115 * \brief Even if the Richards equation is steady-state, this equation is 116 * solved at each iteration. 117 */ 118 119 CS_GWF_FORCE_RICHARDS_ITERATIONS = 1<< 6, /* = 64 */ 120 121 /*! 122 * \brief Compute the mean-value of the hydraulic head field and subtract 123 * this mean-value to get a field with zero mean-value. It's important 124 * to set this flag if no boundary condition is given. 125 */ 126 127 CS_GWF_RESCALE_HEAD_TO_ZERO_MEAN_VALUE = 1<< 7, /* = 128 */ 128 129 /*! 130 * \brief Activate a treatment to enforce a Darcy flux to be divergence-free 131 */ 132 133 CS_GWF_ENFORCE_DIVERGENCE_FREE = 1<< 8 /* = 256 */ 134 135 /*! 136 * @} 137 */ 138 139 } cs_gwf_model_bit_t; 140 141 /*! 142 * @name Flags specifying the kind of post-processing to perform in 143 * the groundwater flow module 144 * @{ 145 * 146 * \def CS_GWF_POST_SOIL_CAPACITY 147 * \brief Activate the post-processing of the soil capacity (property in front 148 * of the unsteady term in Richards equation) 149 * 150 * \def CS_GWF_POST_LIQUID_SATURATION 151 * \brief Activate the post-processing of the liquid saturation (also nammed 152 * "moisture content" in case of single phase flow) 153 * 154 * \def CS_GWF_POST_PERMEABILITY 155 * \brief Activate the post-processing of the permeability field 156 * 157 * \def CS_GWF_POST_DARCY_FLUX_BALANCE 158 * \brief Compute the overall balance at the different boundaries of 159 * the Darcy flux 160 * 161 * \def CS_GWF_POST_DARCY_FLUX_DIVERGENCE 162 * \brief Compute in each control volume (vertices or cells w.r.t the space 163 * scheme) the divergence of the Darcy flux 164 * 165 * \def CS_GWF_POST_DARCY_FLUX_AT_BOUNDARY 166 * \brief Define a field at boundary faces for the Darcy flux and activate the 167 * post-processing 168 */ 169 170 #define CS_GWF_POST_SOIL_CAPACITY (1 << 0) 171 #define CS_GWF_POST_LIQUID_SATURATION (1 << 1) 172 #define CS_GWF_POST_PERMEABILITY (1 << 2) 173 #define CS_GWF_POST_DARCY_FLUX_BALANCE (1 << 3) 174 #define CS_GWF_POST_DARCY_FLUX_DIVERGENCE (1 << 4) 175 #define CS_GWF_POST_DARCY_FLUX_AT_BOUNDARY (1 << 5) 176 177 /*! 178 * @} 179 */ 180 181 /*============================================================================ 182 * Type definitions 183 *============================================================================*/ 184 185 typedef struct _gwf_t cs_gwf_t; 186 187 /*============================================================================ 188 * Public function prototypes 189 *============================================================================*/ 190 191 /*----------------------------------------------------------------------------*/ 192 /*! 193 * \brief Check if the groundwater flow module has been activated 194 * 195 * \return true or false 196 */ 197 /*----------------------------------------------------------------------------*/ 198 199 bool 200 cs_gwf_is_activated(void); 201 202 /*----------------------------------------------------------------------------*/ 203 /*! 204 * \brief Initialize the module dedicated to groundwater flows 205 * 206 * \param[in] pty_type type of permeability (iso, ortho...) 207 * \param[in] model type of physical modelling 208 * \param[in] option_flag optional flag to specify this module 209 * 210 * \return a pointer to a new allocated groundwater flow structure 211 */ 212 /*----------------------------------------------------------------------------*/ 213 214 cs_gwf_t * 215 cs_gwf_activate(cs_property_type_t pty_type, 216 cs_gwf_model_type_t model, 217 cs_gwf_option_flag_t option_flag); 218 219 /*----------------------------------------------------------------------------*/ 220 /*! 221 * \brief Free all structures related to groundwater flows 222 * 223 * \return a NULL pointer 224 */ 225 /*----------------------------------------------------------------------------*/ 226 227 cs_gwf_t * 228 cs_gwf_destroy_all(void); 229 230 /*----------------------------------------------------------------------------*/ 231 /*! 232 * \brief Set the parameters defining the two-phase flow model. 233 * Use SI unit if not prescribed otherwise. 234 * 235 * \param[in] l_mass_density mass density of the main liquid component 236 * \param[in] l_viscosity viscosity in the liquid phase (Pa.s) 237 * \param[in] g_viscosity viscosity in the gas phase (Pa.s) 238 * \param[in] l_diffusivity_h diffusivity of the main gas component in the 239 * liquid phase 240 * \param[in] w_molar_mass molar mass of the main liquid component 241 * \param[in] h_molar_mass molar mass of the main gas component 242 * \param[in] ref_temperature reference temperature 243 * \param[in] henry_constant constant in the Henry law 244 */ 245 /*----------------------------------------------------------------------------*/ 246 247 void 248 cs_gwf_set_two_phase_model(cs_real_t l_mass_density, 249 cs_real_t l_viscosity, 250 cs_real_t g_viscosity, 251 cs_real_t l_diffusivity_h, 252 cs_real_t w_molar_mass, 253 cs_real_t h_molar_mass, 254 cs_real_t ref_temperature, 255 cs_real_t henry_constant); 256 257 /*----------------------------------------------------------------------------*/ 258 /*! 259 * \brief Summary of the main cs_gwf_t structure 260 */ 261 /*----------------------------------------------------------------------------*/ 262 263 void 264 cs_gwf_log_setup(void); 265 266 /*----------------------------------------------------------------------------*/ 267 /*! 268 * \brief Set the flag dedicated to the post-processing of the GWF module 269 * 270 * \param[in] post_flag flag to set 271 * \param[in] reset reset post flag before 272 */ 273 /*----------------------------------------------------------------------------*/ 274 275 void 276 cs_gwf_set_post_options(cs_flag_t post_flag, 277 bool reset); 278 279 /*----------------------------------------------------------------------------*/ 280 /*! 281 * \brief Retrieve the advection field related to the Darcy flux in the liquid 282 * phase 283 * 284 * \return a pointer to a cs_adv_field_t structure or NULL 285 */ 286 /*----------------------------------------------------------------------------*/ 287 288 cs_adv_field_t * 289 cs_gwf_get_adv_field(void); 290 291 /*----------------------------------------------------------------------------*/ 292 /*! 293 * \brief Retrieve the head used in soil updates when an unsaturated 294 * single-phase flow is considered. These values are located at cells. 295 * 296 * \return a pointer to the requested array of values or NULL if not defined 297 */ 298 /*----------------------------------------------------------------------------*/ 299 300 cs_real_t * 301 cs_gwf_get_uspf_head_in_law(void); 302 303 /*----------------------------------------------------------------------------*/ 304 /*! 305 * \brief Create and add a new cs_gwf_soil_t structure. An initialization by 306 * default of all members is performed. 307 * 308 * \param[in] z_name name of the volume zone corresponding to the soil 309 * \param[in] bulk_density value of the mass density 310 * \param[in] sat_moisture value of the saturated moisture content 311 * \param[in] soil_model type of modelling for the hydraulic behavior 312 * 313 * \return a pointer to the new allocated soil structure 314 */ 315 /*----------------------------------------------------------------------------*/ 316 317 cs_gwf_soil_t * 318 cs_gwf_add_soil(const char *z_name, 319 double bulk_density, 320 double sat_moisture, 321 cs_gwf_soil_hydraulic_model_t soil_model); 322 323 /*----------------------------------------------------------------------------*/ 324 /*! 325 * \brief Add a new equation related to the groundwater flow module 326 * This equation is a particular type of unsteady advection-diffusion 327 * reaction eq. 328 * Tracer is advected thanks to the darcian velocity and 329 * diffusion/reaction parameters result from a physical modelling. 330 * Terms solved in the equation are activated according to the settings. 331 * The advection field corresponds to that of the liquid phase. 332 * 333 * \param[in] tr_model physical modelling to consider (0 = default settings) 334 * \param[in] eq_name name of the tracer equation 335 * \param[in] var_name name of the related variable 336 */ 337 /*----------------------------------------------------------------------------*/ 338 339 cs_gwf_tracer_t * 340 cs_gwf_add_tracer(cs_gwf_tracer_model_t tr_model, 341 const char *eq_name, 342 const char *var_name); 343 344 /*----------------------------------------------------------------------------*/ 345 /*! 346 * \brief Add a new equation related to the groundwater flow module 347 * This equation is a particular type of unsteady advection-diffusion 348 * reaction eq. 349 * Tracer is advected thanks to the darcian velocity and 350 * diffusion/reaction parameters result from a physical modelling. 351 * Terms are activated according to the settings. 352 * Modelling of the tracer parameters are left to the user 353 * 354 * \param[in] eq_name name of the tracer equation 355 * \param[in] var_name name of the related variable 356 * \param[in] setup function pointer (predefined prototype) 357 * \param[in] add_terms function pointer (predefined prototype) 358 */ 359 /*----------------------------------------------------------------------------*/ 360 361 cs_gwf_tracer_t * 362 cs_gwf_add_user_tracer(const char *eq_name, 363 const char *var_name, 364 cs_gwf_tracer_setup_t *setup, 365 cs_gwf_tracer_add_terms_t *add_terms); 366 367 /*----------------------------------------------------------------------------*/ 368 /*! 369 * \brief Retrieve the pointer to the cs_gwf_tracer_t structure associated to 370 * the name given as parameter 371 * 372 * \param[in] eq_name name of the tracer equation 373 * 374 * \return the pointer to a cs_gwf_tracer_t structure 375 */ 376 /*----------------------------------------------------------------------------*/ 377 378 cs_gwf_tracer_t * 379 cs_gwf_tracer_by_name(const char *eq_name); 380 381 /*----------------------------------------------------------------------------*/ 382 /*! 383 * \brief Add new terms if needed (such as diffusion or reaction) to tracer 384 * equations according to the settings 385 */ 386 /*----------------------------------------------------------------------------*/ 387 388 void 389 cs_gwf_add_tracer_terms(void); 390 391 /*----------------------------------------------------------------------------*/ 392 /*! 393 * \brief Predefined settings for the Richards equation and the related 394 * equations defining the groundwater flow module. 395 * At this stage, all soils have been defined and equatino parameters 396 * are set. Create new cs_field_t structures according to the setting. 397 */ 398 /*----------------------------------------------------------------------------*/ 399 400 void 401 cs_gwf_init_setup(void); 402 403 /*----------------------------------------------------------------------------*/ 404 /*! 405 * \brief Last initialization step of the groundwater flow module. At this 406 * stage, the mesh quantities are defined. 407 * 408 * \param[in] connect pointer to a cs_cdo_connect_t structure 409 * \param[in] quant pointer to a cs_cdo_quantities_t structure 410 */ 411 /*----------------------------------------------------------------------------*/ 412 413 void 414 cs_gwf_finalize_setup(const cs_cdo_connect_t *connect, 415 const cs_cdo_quantities_t *quant); 416 417 /*----------------------------------------------------------------------------*/ 418 /*! 419 * \brief Update the groundwater system (pressure head, head in law, moisture 420 * content, darcian velocity, soil capacity or permeability if needed) 421 * 422 * \param[in] mesh pointer to a cs_mesh_t structure 423 * \param[in] connect pointer to a cs_cdo_connect_t structure 424 * \param[in] quant pointer to a cs_cdo_quantities_t structure 425 * \param[in] ts pointer to a cs_time_step_t structure 426 * \param[in] cur2prev true or false 427 */ 428 /*----------------------------------------------------------------------------*/ 429 430 void 431 cs_gwf_update(const cs_mesh_t *mesh, 432 const cs_cdo_connect_t *connect, 433 const cs_cdo_quantities_t *quant, 434 const cs_time_step_t *ts, 435 bool cur2prev); 436 437 /*----------------------------------------------------------------------------*/ 438 /*! 439 * \brief Compute the steady-state of the groundwater flows module. 440 * Nothing is done if all equations are unsteady. 441 * 442 * \param[in] mesh pointer to a cs_mesh_t structure 443 * \param[in] time_step pointer to a cs_time_step_t structure 444 * \param[in] connect pointer to a cs_cdo_connect_t structure 445 * \param[in] cdoq pointer to a cs_cdo_quantities_t structure 446 */ 447 /*----------------------------------------------------------------------------*/ 448 449 void 450 cs_gwf_compute_steady_state(const cs_mesh_t *mesh, 451 const cs_time_step_t *time_step, 452 const cs_cdo_connect_t *connect, 453 const cs_cdo_quantities_t *cdoq); 454 455 /*----------------------------------------------------------------------------*/ 456 /*! 457 * \brief Compute the system related to groundwater flows module 458 * 459 * \param[in] mesh pointer to a cs_mesh_t structure 460 * \param[in] time_step pointer to a cs_time_step_t structure 461 * \param[in] connect pointer to a cs_cdo_connect_t structure 462 * \param[in] cdoq pointer to a cs_cdo_quantities_t structure 463 */ 464 /*----------------------------------------------------------------------------*/ 465 466 void 467 cs_gwf_compute(const cs_mesh_t *mesh, 468 const cs_time_step_t *time_step, 469 const cs_cdo_connect_t *connect, 470 const cs_cdo_quantities_t *cdoq); 471 472 /*----------------------------------------------------------------------------*/ 473 /*! 474 * \brief Compute the integral over a given set of cells of the field related 475 * to a tracer equation. This integral turns out to be exact for linear 476 * functions. 477 * 478 * \param[in] connect pointer to a \ref cs_cdo_connect_t structure 479 * \param[in] cdoq pointer to a \ref cs_cdo_quantities_t structure 480 * \param[in] tracer pointer to a \ref cs_gwf_tracer_t structure 481 * \param[in] z_name name of the volumic zone where the integral is done 482 * (if NULL or "" all cells are considered) 483 * 484 * \return the value of the integral 485 */ 486 /*----------------------------------------------------------------------------*/ 487 488 cs_real_t 489 cs_gwf_integrate_tracer(const cs_cdo_connect_t *connect, 490 const cs_cdo_quantities_t *cdoq, 491 const cs_gwf_tracer_t *tracer, 492 const char *z_name); 493 494 /*----------------------------------------------------------------------------*/ 495 /*! 496 * \brief Predefined extra-operations for the groundwater flow module 497 * 498 * \param[in] connect pointer to a cs_cdo_connect_t structure 499 * \param[in] cdoq pointer to a cs_cdo_quantities_t structure 500 */ 501 /*----------------------------------------------------------------------------*/ 502 503 void 504 cs_gwf_extra_op(const cs_cdo_connect_t *connect, 505 const cs_cdo_quantities_t *cdoq); 506 507 /*----------------------------------------------------------------------------*/ 508 /*! 509 * \brief Predefined post-processing output for the groundwater flow module 510 * in case of saturated single-phase flows (sspf) in porous media. 511 * Prototype of this function is given since it is a function pointer 512 * defined in cs_post.h (\ref cs_post_time_mesh_dep_output_t) 513 * 514 * \param[in, out] input pointer to a optional structure (here a 515 * cs_gwf_t structure) 516 * \param[in] mesh_id id of the output mesh for the current call 517 * \param[in] cat_id category id of the output mesh for this call 518 * \param[in] ent_flag indicate global presence of cells (ent_flag[0]), 519 * interior faces (ent_flag[1]), boundary faces 520 * (ent_flag[2]), particles (ent_flag[3]) or probes 521 * (ent_flag[4]) 522 * \param[in] n_cells local number of cells of post_mesh 523 * \param[in] n_i_faces local number of interior faces of post_mesh 524 * \param[in] n_b_faces local number of boundary faces of post_mesh 525 * \param[in] cell_ids list of cells (0 to n-1) 526 * \param[in] i_face_ids list of interior faces (0 to n-1) 527 * \param[in] b_face_ids list of boundary faces (0 to n-1) 528 * \param[in] time_step pointer to a cs_time_step_t struct. 529 */ 530 /*----------------------------------------------------------------------------*/ 531 532 void 533 cs_gwf_extra_post_sspf(void *input, 534 int mesh_id, 535 int cat_id, 536 int ent_flag[5], 537 cs_lnum_t n_cells, 538 cs_lnum_t n_i_faces, 539 cs_lnum_t n_b_faces, 540 const cs_lnum_t cell_ids[], 541 const cs_lnum_t i_face_ids[], 542 const cs_lnum_t b_face_ids[], 543 const cs_time_step_t *time_step); 544 545 /*----------------------------------------------------------------------------*/ 546 /*! 547 * \brief Predefined post-processing output for the groundwater flow module 548 * in case of unsaturated single-phase flows (uspf) in porous media. 549 * Prototype of this function is given since it is a function pointer 550 * defined in cs_post.h (\ref cs_post_time_mesh_dep_output_t) 551 * 552 * \param[in, out] input pointer to a optional structure (here a 553 * cs_gwf_t structure) 554 * \param[in] mesh_id id of the output mesh for the current call 555 * \param[in] cat_id category id of the output mesh for this call 556 * \param[in] ent_flag indicate global presence of cells (ent_flag[0]), 557 * interior faces (ent_flag[1]), boundary faces 558 * (ent_flag[2]), particles (ent_flag[3]) or probes 559 * (ent_flag[4]) 560 * \param[in] n_cells local number of cells of post_mesh 561 * \param[in] n_i_faces local number of interior faces of post_mesh 562 * \param[in] n_b_faces local number of boundary faces of post_mesh 563 * \param[in] cell_ids list of cells (0 to n-1) 564 * \param[in] i_face_ids list of interior faces (0 to n-1) 565 * \param[in] b_face_ids list of boundary faces (0 to n-1) 566 * \param[in] time_step pointer to a cs_time_step_t struct. 567 */ 568 /*----------------------------------------------------------------------------*/ 569 570 void 571 cs_gwf_extra_post_uspf(void *input, 572 int mesh_id, 573 int cat_id, 574 int ent_flag[5], 575 cs_lnum_t n_cells, 576 cs_lnum_t n_i_faces, 577 cs_lnum_t n_b_faces, 578 const cs_lnum_t cell_ids[], 579 const cs_lnum_t i_face_ids[], 580 const cs_lnum_t b_face_ids[], 581 const cs_time_step_t *time_step); 582 583 /*----------------------------------------------------------------------------*/ 584 /*! 585 * \brief Predefined post-processing output for the groundwater flow module 586 * in case of miscible two-phase flows (mtpf) in porous media. 587 * Prototype of this function is given since it is a function pointer 588 * defined in cs_post.h (\ref cs_post_time_mesh_dep_output_t) 589 * 590 * \param[in, out] input pointer to a optional structure (here a 591 * cs_gwf_t structure) 592 * \param[in] mesh_id id of the output mesh for the current call 593 * \param[in] cat_id category id of the output mesh for this call 594 * \param[in] ent_flag indicate global presence of cells (ent_flag[0]), 595 * interior faces (ent_flag[1]), boundary faces 596 * (ent_flag[2]), particles (ent_flag[3]) or probes 597 * (ent_flag[4]) 598 * \param[in] n_cells local number of cells of post_mesh 599 * \param[in] n_i_faces local number of interior faces of post_mesh 600 * \param[in] n_b_faces local number of boundary faces of post_mesh 601 * \param[in] cell_ids list of cells (0 to n-1) 602 * \param[in] i_face_ids list of interior faces (0 to n-1) 603 * \param[in] b_face_ids list of boundary faces (0 to n-1) 604 * \param[in] time_step pointer to a cs_time_step_t struct. 605 */ 606 /*----------------------------------------------------------------------------*/ 607 608 void 609 cs_gwf_extra_post_mtpf(void *input, 610 int mesh_id, 611 int cat_id, 612 int ent_flag[5], 613 cs_lnum_t n_cells, 614 cs_lnum_t n_i_faces, 615 cs_lnum_t n_b_faces, 616 const cs_lnum_t cell_ids[], 617 const cs_lnum_t i_face_ids[], 618 const cs_lnum_t b_face_ids[], 619 const cs_time_step_t *time_step); 620 621 /*----------------------------------------------------------------------------*/ 622 623 END_C_DECLS 624 625 #endif /* __CS_GWF_H__ */ 626