1 #ifndef __CS_GWF_PRIV_H__
2 #define __CS_GWF_PRIV_H__
3 
4 /*============================================================================
5  * Structures/types related to the groundwater flow module
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  *  Local headers
30  *----------------------------------------------------------------------------*/
31 
32 #include "cs_advection_field.h"
33 #include "cs_gwf.h"
34 
35 /*----------------------------------------------------------------------------*/
36 
37 BEGIN_C_DECLS
38 
39 /*============================================================================
40  * Macro definitions
41  *============================================================================*/
42 
43 /*============================================================================
44  * Type definitions
45  *============================================================================*/
46 
47 /*! \struct cs_gwf_darcy_flux_t
48  *
49  * \brief Structure to handle the Darcy flux
50  */
51 
52 typedef struct {
53 
54   /*!
55    * \var adv_field
56    * Pointer to a \ref cs_adv_field_t structure. Darcy advective flux in the
57    * liquid phase. This structure is used to define the advective term in
58    * tracer equations.
59    *
60    * \var flux_location
61    * Indicate where the arrays defining the Darcy fluxes are located
62    *
63    * \var flux_val
64    * Array storing the liquid Darcian flux in each location (for instance the
65    * dual faces associated to each cell)
66    *
67    * \var boundary_flux_val
68    * Array storing the normal Darcian flux across the boundary of the
69    * computational domain for the liquid phase. This is an optional array.
70    */
71 
72   cs_adv_field_t               *adv_field;
73   cs_flag_t                     flux_location;
74   cs_real_t                    *flux_val;
75   cs_real_t                    *boundary_flux_val;
76 
77 } cs_gwf_darcy_flux_t;
78 
79 
80 /*! \struct cs_gwf_saturated_single_phase_t
81  *
82  * \brief Structure to handle the modelling of a single-phase flows in a porous
83  *        media considered as saturated.
84  *
85  * Several simplifications are operated in this context. Only the liquid phase
86  * is taken into account.
87  */
88 
89 typedef struct {
90 
91   /*!
92    * @name Equation
93    * @{
94    *
95    * \var richards
96 
97    * The Richards equation is the governing equation which corresponds to the
98    * mass conservation of water. "hydraulic_head" is the associated variable
99    * "permeability" is the diffusion property related to this equation.
100    */
101 
102   cs_equation_t                *richards;
103 
104   /*!
105    * @}
106    * @name Darcy (advection) field
107    * @{
108    *
109    * \var darcy
110    * Pointer to a \ref cs_gwf_darcy_flux_t structure. Darcy advective flux in
111    * the liquid phase. This structure is used to define the advective term in
112    * tracer equations.
113    */
114 
115   cs_gwf_darcy_flux_t          *darcy;
116 
117   /*!
118    * @}
119    * @name Properties related to the model
120    * @{
121    *
122    * \var moisture_content
123    * This quantity describes the level of saturation in a soil. This is a
124    * constant value in case of a satured soil and variable one in case of a
125    * unsaturated soil.
126    */
127 
128   cs_property_t                *moisture_content;
129 
130   /*!
131    * @}
132    * @name Additional fields/arrays
133    * @{
134    *
135    * \var pressure_head
136    * Pointer to a \ref cs_field_t structure. Allocated only if the gravitation
137    * effect is active. Location of this field depends on the discretization
138    * scheme used to solve the Richards equation.
139    * The pressure head is denoted by h, hydraulic head (unknowns solved in the
140    * Richards eq.) is denoted by H and there are linked by:
141    * h = H - gravity_potential
142    */
143 
144   cs_field_t                   *pressure_head;
145 
146   /*!
147    * @}
148    */
149 
150 } cs_gwf_saturated_single_phase_t;
151 
152 /*! \struct cs_gwf_unsaturated_single_phase_t
153  *
154  * \brief Structure to handle the modelling of a single-phase flows in a porous
155  *        media considered as saturated or not. Several simplifications can be
156  *        be operated in this context. Only the liquid phase is taken into
157  *        account.
158  */
159 
160 typedef struct {
161 
162     /*!
163    * @name Equation
164    * @{
165    *
166    * \var richards
167    * The Richards equation is the governing equation which corresponds to the
168    * mass conservation of water. "hydraulic_head" is the associated variable
169    * "permeability" is the diffusion property related to this equation.
170    */
171 
172   cs_equation_t                *richards;
173 
174   /*!
175    * @}
176    * @name Darcy (advection) field
177    * @{
178    *
179    * \var darcy
180    * Pointer to a \ref cs_gwf_darcy_flux_t structure. Darcy advective flux in
181    * the liquid phase. This structure is used to define the advective term in
182    * tracer equations.
183    */
184 
185   cs_gwf_darcy_flux_t          *darcy;
186 
187   /*!
188    * @}
189    * @name Properties related to the model
190    * @{
191    *
192    * \var permeability
193    * This quantity the product of the absolute permeability and the relative
194    * permeability
195    *
196    * \var moisture_content
197    * This quantity describes the level of saturation in a soil. This is a
198    * constant value in case of a satured soil and variable one in case of a
199    * unsaturated soil.
200    *
201    * \var soil_capacity
202    * property attached to the unsteady term in the Richards equation
203    */
204 
205   cs_property_t                *permeability;
206   cs_property_t                *moisture_content;
207   cs_property_t                *soil_capacity;
208 
209   /*!
210    * @}
211    * @name Additional fields/arrays
212    * @{
213    *
214    * \var permeability_field
215    * Pointer to a \ref cs_field_t structure. May be not allocated according to
216    * the set of options. Store the value of the full permeability field in each
217    * cell (k = k_abs * k_rel)
218    *
219    * \var moisture_field
220    * Pointer to a \ref cs_field_t structure. Structure storing the value of the
221    * moisture content in each cell. This is an optional structure i.e. it may
222    * be set to NULL (for instance in case of a satured soil on the full
223    * domain).
224    *
225    * \var capacity_field
226    * Pointer to a \ref cs_field_t structure. Structure storing the value of the
227    * soil capacity in each cell. This is an optional structure i.e. it is
228    * set to NULL in case of a satured soil on the full domain.
229    *
230    * \var pressure_head
231    * Pointer to a \ref cs_field_t structure. Allocated only if the gravitation
232    * effect is active. Location of this field depends on the discretization
233    * scheme used to solve the Richards equation.
234    * The pressure head is denoted by h, hydraulic head (unknowns solved in the
235    * Richards eq.) is denoted by H and there are linked by:
236    * h = H - gravity_potential
237    *
238    * \var head_in_law
239    * Array of values located at the same location as the hydraulic head solved
240    * in the Richards equation. The values stored in this array are used to
241    * update related quantities with law such as Van Genuchten/Mualen.
242    */
243 
244   cs_field_t                   *permeability_field;
245   cs_field_t                   *moisture_field;
246   cs_field_t                   *capacity_field;
247   cs_field_t                   *pressure_head;
248 
249   cs_real_t                    *head_in_law;
250 
251   /*!
252    * @}
253    */
254 
255 } cs_gwf_unsaturated_single_phase_t;
256 
257 /*! \struct cs_gwf_miscible_two_phase_t
258  *
259  * \brief Structure to handle the modelling of miscible two-phase flows in a
260  *        porous media.
261 
262  * The model follows what is depicted in "Finite volume approximation of a
263  * diffusion-dissolution model and application to nuclear waste storage"
264  * O. Angelini, C. Chavant, E. Chénier, R. Eymard and S. Granet in Mathematics
265  * and Computers in Simulation (2011), 81 (10), pp. 2001--2017
266  *
267  * Main assumptions are:
268  *   - No water in the gas phase
269  *   - Incompressibility of the liquid phase
270  *   - Hydrogen pressure is given by the "perfect gas" law in the gas phase and
271  *     the Henry's law in the liquid phase
272  *
273  * The two primitive variables are the liquid and gas pressures with a specific
274  * treatment in the saturated case to handle the gaz pressure (cf. the cited
275  * article or Angelini's PhD thesis
276  *
277  * Notations are the following :
278  * - Two phases: Liquid phase denoted by "l" and gaseous phase denoted by g
279  * - Two components: water denoted by w and a gaseous component (let's say
280  *   hydrogen) denoted by h. The gaseous component is present in the two phases
281  *   whereas water is only considered in the liquid phase.
282  */
283 
284 typedef struct {
285 
286   /* Set of equations associated to this modelling */
287 
288   /*!
289    * @name Equations
290    * @{
291    *
292    * \var w_eq
293    * Equation of conservation for the water component. Only the liquid phase is
294    * considered. One assumes to wapter vapour in the gas phase.
295    *
296    */
297 
298   cs_equation_t                *w_eq;
299 
300   /*! \var h_eq
301    * Equation of conservation for the (di)hydrogen. Hydrogen can be present in
302    * the liquid or in the gas phase.
303    */
304 
305   cs_equation_t                *h_eq;
306 
307   /*!
308    * @}
309    * @name Darcy (advection) fields
310    * @{
311    *
312    * \var l_darcy
313    * Pointer to a \ref cs_gwf_darcy_flux_t structure. Darcy advective flux in
314    * the liquid phase
315    *
316    * \var g_darcy
317    * Pointer to a \ref cs_gwf_darcy_flux_t structure. Darcy advective flux in
318    * the gas phase
319    */
320 
321   cs_gwf_darcy_flux_t          *l_darcy;
322   cs_gwf_darcy_flux_t          *g_darcy;
323 
324   /*!
325    * @}
326    * @name Properties related to the model
327    * @{
328    *
329    * \var time_w_eq_pty
330    * Property related to the unsteady term of the water conservation equation
331    *
332    * \var diff_w_eq_pty
333    * Property related to the diffusion term of the water conservation equation
334    *
335    * \var time_h_eq_pty
336    * Property related to the unsteady term of the hydrogen conservation equation
337    *
338    * \var diff_hl_eq_pty
339    * Property related to the diffusion term of the hydrogen conservation
340    * equation (part related to the liquid phase)
341    *
342    * \var diff_hg_eq_pty
343    * Property related to the diffusion term of the hydrogen conservation
344    * equation (part related to the gas phase)
345    */
346 
347   cs_property_t                *time_w_eq_pty;
348   cs_property_t                *diff_w_eq_pty;
349   cs_property_t                *time_h_eq_pty;
350   cs_property_t                *diff_hl_eq_pty;
351   cs_property_t                *diff_hg_eq_pty;
352 
353   /*!
354    * @}
355    * @name Additional fields
356    * @{
357    *
358    * \var l_saturation
359    * Pointer to a \ref cs_field_t structure. Liquid saturation at cells.
360    *
361    * \var c_pressure
362    * Pointer to a \ref cs_field_t structure named "capillarity_pressure".
363    * Capillarity pressure P_c = P_g - P_l
364    *
365    * \var l_pressure
366    * Pointer to a \ref cs_field_t structure named "liquid_pressure".
367    * Pressure in the liquid phase is denoted by P_l.
368    *
369    * \var g_pressure
370    * Pointer to a \ref cs_field_t structure named "gas_pressure".
371    * Pressure in the gas phase is denoted by P_g.
372    */
373 
374   cs_field_t                   *l_saturation;
375   cs_field_t                   *c_pressure;
376   cs_field_t                   *l_pressure;
377   cs_field_t                   *g_pressure;
378 
379   /*!
380    * @}
381    * @name Additional arrays
382    * @{
383    *
384    * \var time_w_eq_array
385    * Values in each cell of the coefficient appearing in front of the unsteady
386    * term in the water conservation equation. This array is linked to the
387    * \ref time_w_eq_pty (size = n_cells)
388    *
389    * \var diff_w_eq_array
390    * Values in each cell of the coefficient appearing in the diffusion
391    * term in the water conservation equation. This array is linked to the
392    * \ref diff_w_eq_pty (size = n_cells)
393    *
394    * \var time_h_eq_array
395    * Values in each cell of the coefficient appearing in front of the unsteady
396    * term in the hydrogen conservation equation. This array is linked to the
397    * \ref time_h_eq_pty (size = n_cells)
398    *
399    * \var diff_hl_eq_array
400    * Values in each cell of the coefficient appearing in the diffusion
401    * term for the liquid phase in the hydrogen conservation equation.
402    * This array is linked to the \ref diff_hl_eq_pty (size = n_cells)
403    *
404    * \var diff_hg_eq_array
405    * Values in each cell of the coefficient appearing in the diffusion
406    * term for the gas phase in the hydrogen conservation equation.
407    * This array is linked to the \ref diff_hl_eq_pty (size = n_cells)
408    *
409    * \var l_rel_permeability
410    * Values in each cell of the relative permeability in the liquid phase.
411    * This quantity is used either in the water conservation or in the hydrogen
412    * conservation. This enables also to recover the (full) permeability since
413    * permeability = abs_permeability * rel_l_permeability
414    */
415 
416   cs_real_t                    *time_w_eq_array;
417   cs_real_t                    *diff_w_eq_array;
418   cs_real_t                    *time_h_eq_array;
419   cs_real_t                    *diff_hl_eq_array;
420   cs_real_t                    *diff_hg_eq_array;
421   cs_real_t                    *l_rel_permeability;
422 
423   /*!
424    * @}
425    * @name Model parameters
426    * @{
427    *
428    * \var l_mass_density
429    * Mass density in the liquid phase. With the model assumptions, this
430    * corresponds to the mass density of the main component in the liquid phase
431    * (e.g. water) in kg.m^-3
432    *
433    * \var l_viscosity
434    * Viscosity in the liquid phase (assumed to be constant) in Pa.s
435    *
436    * \var g_viscosity
437    * Viscosity in the gas phase (assumed to be constant) in Pa.s
438    *
439    * \var w_molar_mass
440    * Molar mass of the main component in the liquid phase (e.g. water) in
441    * kg.mol^-1
442    *
443    * \var h_molar_mass
444    * Molar mass of the main component in the gas phase (e.g. hydrogen) in
445    * kg.mol^-1
446    *
447    * \var l_diffusivity_h
448    * Molecular diffusivity of the hydrogen in the liquid phase in m^2.s^-1
449    *
450    * \var ref_temperature
451    * Reference temperature used in the "perfect gas" law (this is used when no
452    * thermal equation is solved). One expects a temperature in Kelvin.
453    *
454    * \henry_constant
455    * Value of the Henry constant used in the Henry's law. Setting a very low
456    * value for this constant enables the model to degenerate into an immiscible
457    * model.
458    */
459 
460   cs_real_t                     l_mass_density;
461   cs_real_t                     l_viscosity;
462   cs_real_t                     g_viscosity;
463   cs_real_t                     l_diffusivity_h;
464   cs_real_t                     w_molar_mass;
465   cs_real_t                     h_molar_mass;
466   cs_real_t                     ref_temperature;
467   cs_real_t                     henry_constant;
468 
469   /*!
470    * @}
471    */
472 
473 } cs_gwf_miscible_two_phase_t;
474 
475 
476 /*! \struct _gwf_t
477  *
478  * \brief Main set of parameters/structures to manage the groundwater flow
479  *        (GWF) module. This is an explicit definition of the structure
480  *        \ref cs_gwf_t
481  */
482 
483 struct _gwf_t {
484 
485   /*!
486    * @name Metadata
487    * @{
488    *
489    * \var model
490    * Model used to describe the behavior of the flow in the GWF module (system
491    * of equations related to the chosen physical modelling). See \ref
492    * cs_gwf_model_type_t for more details on each model
493    *
494    * \var flag
495    * Flag dedicated to general options to handle the GWF module
496    *
497    * \var post_flag
498    * Flag dedicated to the (automatic) post-processing of the GWF module
499    */
500 
501   cs_gwf_model_type_t           model;
502   cs_gwf_option_flag_t          flag;
503   cs_flag_t                     post_flag;
504 
505   /*!
506    * @}
507    * @name Properties
508    * @{
509    *
510    * \var abs_permeability
511    * Absolute (or intrinsic) permeability which characterizes the behavior of a
512    * soil. According to the model of soil, this absolute permeability can be
513    * weigthed by a relative (scalar-valued) permeability. In the simplest case
514    * (saturated soil) the relative permeability is useless since this is equal
515    * to 1 (no weigth).
516    */
517 
518   cs_property_t                *abs_permeability;
519 
520   /*!
521    * @}
522    * @name Modelling context
523    * @{
524    *
525    * \var model_context
526    * Pointer to a structure cast on-the-fly which depends on the choice of
527    * model (for instance single-phase or two-phase flows in porous media)
528    */
529 
530   void                         *model_context;
531 
532   /*!
533    * @}
534    * @name Associated tracers
535    * @{
536    *
537    * \var n_tracers
538    * Number of tracers to consider. Each tracer is related to a scalar-valued
539    * transport equation. There is at least an unsteady term and an advection
540    * term. The advection term is linked to the Darcy flux in the liquid phase.
541    *
542    * \var tracers
543    * Array of pointers to the \ref cs_gwf_tracer_t structure which manages each
544    * tracer equation and its related quantities/metadata.
545    *
546    * \var finalize_tracer_setup
547    * This is a function pointer to finalize the setup of a tracer
548    * equation. There is a default pointer but this can be overloaded by a
549    * user-defined function in the case of a user-defined tracer.
550    *
551    * \var add_tracer_terms
552    * This is a function pointer to add non-standard terms in a tracer
553    * equation. There is a default pointer but this can be overloaded by a
554    * user-defined function in the case of a user-defined tracer.
555    */
556 
557   int                            n_tracers;
558   cs_gwf_tracer_t              **tracers;
559   cs_gwf_tracer_setup_t        **finalize_tracer_setup;
560   cs_gwf_tracer_add_terms_t    **add_tracer_terms;
561 
562   /*!
563    * @}
564    */
565 
566 };
567 
568 /*============================================================================
569  * Public function prototypes
570  *============================================================================*/
571 
572 /*----------------------------------------------------------------------------*/
573 /*!
574  * \brief  Allocate and initialize a \ref cs_gwf_darcy_flux_t structure
575  *
576  * \param[in]      loc_flag   flag to define where the flux is defined
577  *
578  * \return a pointer to the newly allocated structure
579  */
580 /*----------------------------------------------------------------------------*/
581 
582 cs_gwf_darcy_flux_t *
583 cs_gwf_darcy_flux_create(cs_flag_t         loc_flag);
584 
585 /*----------------------------------------------------------------------------*/
586 /*!
587  * \brief  Free a \ref cs_gwf_darcy_flux_t structure
588  *
589  * \param[in, out] p_darcy   pointer of pointer to the darcy structure
590  */
591 /*----------------------------------------------------------------------------*/
592 
593 void
594 cs_gwf_darcy_flux_free(cs_gwf_darcy_flux_t  **p_darcy);
595 
596 /*----------------------------------------------------------------------------*/
597 /*!
598  * \brief  Log a \ref cs_gwf_darcy_flux_t structure
599  *
600  * \param[in, out] darcy   pointer to the darcy structure
601  */
602 /*----------------------------------------------------------------------------*/
603 
604 void
605 cs_gwf_darcy_flux_log(cs_gwf_darcy_flux_t    *darcy);
606 
607 /*----------------------------------------------------------------------------*/
608 /*!
609  * \brief  Set the definition of the advection field attached to a
610  *         \ref cs_gwf_darcy_flux_t structure
611  *
612  * \param[in]       connect       pointer to a cs_cdo_connect_t structure
613  * \param[in]       quant         pointer to a cs_cdo_quantities_t structure
614  * \param[in]       space_scheme  space discretization using this structure
615  * \param[in, out]  darcy         pointer to the darcy structure
616  */
617 /*----------------------------------------------------------------------------*/
618 
619 void
620 cs_gwf_darcy_flux_define(const cs_cdo_connect_t       *connect,
621                          const cs_cdo_quantities_t    *quant,
622                          cs_param_space_scheme_t       space_scheme,
623                          cs_gwf_darcy_flux_t          *darcy);
624 
625 /*----------------------------------------------------------------------------*/
626 /*!
627  * \brief  Update the advection field/arrays related to the Darcy flux.
628  *
629  * \param[in]      t_eval    time at which one performs the evaluation
630  * \param[in]      eq        pointer to the equation related to this Darcy flux
631  * \param[in]      cur2prev  true or false
632  * \param[in, out] darcy     pointer to the darcy structure
633  */
634 /*----------------------------------------------------------------------------*/
635 
636 void
637 cs_gwf_darcy_flux_update(const cs_real_t              t_eval,
638                          const cs_equation_t         *eq,
639                          bool                         cur2prev,
640                          cs_gwf_darcy_flux_t         *darcy);
641 
642 /*----------------------------------------------------------------------------*/
643 /*!
644  * \brief  Operate the balance by zone (relying on the splitting arising from
645  *         the boundary settings) for the advection field attached to a \ref
646  *         cs_gwf_darcy_flux_t structure
647  *
648  * \param[in]       connect       pointer to a cs_cdo_connect_t structure
649  * \param[in]       quant         pointer to a cs_cdo_quantities_t structure
650  * \param[in]       eqp           pointer to the set of equation parameters
651  * \param[in, out]  darcy         pointer to the darcy structure
652  */
653 /*----------------------------------------------------------------------------*/
654 
655 void
656 cs_gwf_darcy_flux_balance(const cs_cdo_connect_t       *connect,
657                           const cs_cdo_quantities_t    *quant,
658                           const cs_equation_param_t    *eqp,
659                           cs_gwf_darcy_flux_t          *darcy);
660 
661 /*----------------------------------------------------------------------------*/
662 
663 END_C_DECLS
664 
665 #endif /* __CS_GWF_PRIV_H__ */
666