====== HYPOFE2 ======
===== Description =====
Multiscale law for water-air seepage, pollutant diffusion and advection. Inspired from WAVAT and ADVEC.

Can be parallelized with ELEMB (macroscale) or at the perturbation loop (microscale).

Takes into account the hysteresis in the water retention law when used with FKRSAT.
==== The model ====
This law is only  used for water seepage - air seepage- pollutant diffusion and advection (coupled) for non linear analysis in 2D porous media.

=== Mass conservation of liquid water ===
\[
\underbrace{\frac{\partial}{\partial t} (\rho_s . n . S_{r,w}) + div(\rho_w \vec{q_l})}_{\text{Liquide}} = 0
\]

=== Liquid flow ===
Starting from Darcy's law, the liquid water velocity is:
\[
\vec{q_l} = - \frac{k_w}{\mu_w}\left[ \vec{grad}(p_w) \right]\ \text{where}\ k_w = K_w \frac{\mu_w}{\rho_w g}\left[ m^2\right]
\]

=== Liquid State Equations ===

  - Density: $\rho_w$: \[\rho_w (p_w) = \rho_{wo}\left[ 1+\frac{p_w-p_{w0}}{\chi_w}\right]\]
  - Intrinsic Permeability $k_w$: \\ Depending on the water saturation degree $S_w$ : $k_{r,w} = f(S_w)$ with $k_{w,eff} = k_f k_{r,w}$
  - Saturation degree $S_w$: \\ Depending on succion $s = p_a - p_w : S_w = f(s)$

=== Saturation degree equation (with FKRSAT) ===
ISR = 53 Van Genuchten model (ISR=5) with hysteresis implemented.

The main water retention curves (d=drying, w=wetting) are, according to the Van Genuchten model:
\[S_{ed} = S_{res} + (S_{max}-S_{res}) \left[1 + \left(\frac{s}{a_d}\right)^{n_d}\right]^{-m_d}\] 
\[S_{ew} = S_{res} + (S_{max}-S_{res}) \left[1 + \left(\frac{s}{a_w}\right)^{n_w}\right]^{-m_w}\]

The hysteresis is then defined by:
\[\frac{\partial S_{es}}{\partial s} (\text{wetting}) = \left(\frac{s_w}{s}\right)^b\left(\frac{\partial S_{ew}}{\partial s}\right) \text{ with } s_w = a_w \left(S_e^{-1/m_w}\right)^{1/n_w}\] 
\[\frac{\partial S_{es}}{\partial s} (\text{drying}) = \left(\frac{s_d}{s}\right)^{-b}\left(\frac{\partial S_{ed}}{\partial s}\right) \text{ with } s_d = a_d \left(S_e^{-1/m_d}\right)^{1/n_d}\]

And therefore:
\[S_e^{t+1} = S_e^t + \left(\frac{\partial S_{es}}{\partial s}\right)\times ds\]

The ISR=53 parameters are: CSRW1=$a_d$, CSRW2=$n_d$, CSRW3=$a_w$, CSRW4=$n_w$ and CSRW5=$b$

=== Mass conservation of dry air ===
\[\frac{\partial}{\partial t} (\rho_a . n . S_{r,g}) + div(\rho_a \vec{q_g}) = 0\]

=== Gas flows ===
Starting from Darcy's law, the gas velocity is:
\[
\vec{q_g} = - \frac{k_g}{\mu_g}\left[ \vec{grad}(p_g) + \right]\ \text{où}\ k_g = K_g \frac{\mu_g}{\rho_g g}\left[ m^2\right]
\]

=== Gas State Equation ===

  - Density $\rho_a$ :\\ //Hypothesis// : The air is supposed to be a perfect gas. \[\rho_a (p_a) = \rho_{a,0}\frac{p_a}{p_{a,0}} \]
  - Intrinsic Permeability $k_g$: \\ Depending on the saturation degree $S_g$ : $k_{r,g} = f(S_g)$ with $k_{g,effectif} = k_{g, intrinsic}k_{a,w}$

=== Balance Equation of Pollutant ===
\[\frac{\partial}{\partial x_i} (v_i^p) = 0\]

=== Pollutant flows ===
\[ v_i^p = v_i^{advection} + v_i^{diffusion+dispersion} = C_M v_i^w - D \frac{\partial C_m}{\partial x_i} \]\\
With C_M and C_m [-] the concentration in pollutant at the macroscale and subscale, respectively. $v_i^w$ is the water velocity obtained from Darcy's law and $D$ [m$^2$/s] is the diffusion and dispersion coefficient.
==== Files ====
Prepro: LHYPOFE2.F \\
Lagamine: HYPOFE2.F \\
===== Availability =====
|Plane stress state| NO |
|Plane strain state| YES |
|Axisymmetric state| NO|
|3D state| YES |
|Generalized plane state| NO |
===== Input file =====
==== Parameters defining the type of constitutive law ====
^ Line 1 (2I5, 60A1)^^
|IL|Law number|
|ITYPE| 629 |
|COMMENT| Any comment (up to 60 characters) that will be reproduced on the output listing|
==== Integer parameters ====
^ Line 1 (3I10,2G10.0) ^^
|NLAWFEM2|Number of constitutive laws at the subscale|
|KFLU|Number of DOF: 1=Pw, 2=Pw+C, 3=Pw+Pg, 4=Pw+C+Pg with C the concentration in pollutant|
|MITER|Maximum number of iterations at the subscale|
|CNORM|Norm for the solver of the subscale|
|FACONV|Units of conversion of the RVE (it has a size of 1[-])|



==== Real parameters ====
^ Line 1 (3E10.2,2G10.0) ^^
|VISCW0|Liquid dynamic viscosity $(=\mu_{w,0})\ \left[ Pa.s \right]$|
|RHOW0|Liquid density $(=\rho_{w,0})\ \left[ kg.m^{-3}\right]$|
|UXHIW|Liquid compressibility coefficient $(=1/ \chi_{w})\ \left[ Pa^{-1}\right]$|
|PW0|Initial water pressure $\left[ Pa\right]$|
|T0|Initial temperature $\left[ K\right]$|
^ Line 2 (1G10.0) ^^
|CPINI|Initial pollutant concentration $\left[ -\right]$|
^ Line 3 (3E10.2,2G10.0) ^^
|VISCA0|Gas dynamic viscosity $(=\mu_{a,0})\ \left[Pa.s \right]$|
|RHOA0|Gaz density $(=\rho_{a,0})\ \left[kg.m^{-3}\right]$|
|PMGAS|Gas molar mass $[g/mol]$|
|PA0|Initial gas pressure $\left[ Pa\right]$|
|PHENRY|Henry coefficient|
^ Line 4 (1I10) ^^
|IVAP|= 1 for vapour, = 0 if liquid water only (VAPOUR NOT IMPLEMENTED YET)|
^ Line 5 (3I10) ^^
|ISR|Retention curve (=53 for Van Genuchten with hysteresis)|
|IKW|Water relative permeability curve (=7 for Van Genuchten)|
|IKA|Gas relative permeability curve (=6 for Van Genuchten)|
^ Line 6 (3G10.0)^^
|CKW1|First parameter of IKW|
|CKW2|Second paremeter of IKW|
|CKW3|Third parameter of IKW|
^ Line 7 (2G10.0)^^
|CKA1|First parameter of IKA|
|CKA2|Second parameter of IKA|
^ Line 8 (5G10.0)^^
|CSR1|First parameter of ISR|
|CSR2|Second parameter of ISR|
|CSR3|Third parameter of ISR|
|CSR4|Fourth parameter of ISR|
|CSR5|Fifth parameter of ISR|
^ Line 9 (5G10.0)^^
|SRES|Residual saturation degree $(=S_{res})$|
|SRFIELD|Field saturation degree $(=S_{r, field})$|
|AIREV|Air entry pressure $\left[Pa\right]$|
|AKRMIN|Minimum value of relative permeabikity|
|SRINI|Initial saturation degree|

==== Subscale parameters ====
To be repeated as many time as NLAWFEM2.
^ Line 1 (2I5) ^^
|ILAW2|Number of the subscale constitutive law (=1:NLAWFEM2)|
|ITYPE2|Type of subscale law (=1 for Hydraulic pollutant microscale law)|
^ Line 2 (4G10.0) ^^
|POROS|Material porosity ($=n$)|
|PERMINT|Material intrinsic permeability ($=k_{int}$) $[m^2]$|
|DIFFC|Material diffusion coefficient of the pollutant ($D_{app}$) $[m^2/s]$|
|TORTU|Material tortuosity ($=\tau$)|

===== Stresses =====
==== Number of stresses ====
28
==== Meaning ====
__In 2D state :__
|SIG(1)|$\sigma_x$ (unused)|
|SIG(2)|$\sigma_y$ (unused)|
|SIG(3)|$\sigma_{xy}$ (unused)|
|SIG(4)|$\sigma_z$ (unused)|
|SIG(5)|Homogenised liquid flow along $x$ $(=f_{wx})$|
|SIG(6)|Homogenised liquid flow along $y$ $(=f_{wy})$|
|SIG(7)|Homogenised liquid flow stored $(=f_{we})$|
|SIG(8)|Homogenised mean flow of the pollutant along $x$ $(=(f_{px,a}+f_{px,b})/2)$|
|SIG(9)|Homogenised mean flow of the pollutant along $y$ $(=(f_{py,a}+f_{py,b})/2)$|
|SIG(10)|Homogenised pollutant flow stored (takes advection into account) $(=f_{pe})$|
|SIG(11)|Homogenised diffusive flow of the pollutant along $x$ for the current step $(=f_{px,b})$|
|SIG(12)|Homogenised diffusive flow of the pollutant along $y$ for the current step $(=f_{py,b})$|
|SIG(13)|Homogenised gas flow along $x$ $(=f_{gx})$|
|SIG(14)|Homogenised gas flow along $y$ $(=f_{gy})$|
|SIG(15)|Homogenised gas flow stored $(=f_{ge})$|
|SIG(16)|Advective flow of dissolved gas along $x$ (unused)|
|SIG(17)|Advective flow of dissolved gas along $y$ (unused)|
|SIG(18)|Unused|
|SIG(19)|Unused|
|SIG(20)|Unused|
|SIG(21)|Unused|
|SIG(22)|Unused|
|SIG(23)|Unused|
|SIG(24)|Unused|
|SIG(25)|Unused|
|SIG(26)|Unused|
|SIG(27)|Unused|
|SIG(28)|Unused|

===== State variables =====
==== Number of state variables ====
10 + 5*(Number of Subscale Nodes)\\
/!\ The state variables vector also contains the following information for each subscale node: X,Y,Pw,C,Pg
==== List of state variables ====
|Q(1)|Liquid water mass at the RVE|
|Q(2)|Pollutant mass at the RVE|
|Q(3)|Gaseous air mass at the RVE|
|Q(4)|Homogenised macroscale porosity|
|Q(5)|Water saturation degree|
|Q(6)|Homogenised water relative permeability|
|Q(7)|Homogenised gas relative permeability|
|Q(8)|Homogenised macroscale tortuosity|
|Q(9)|Vapour mass at the RVE (unused)|
|Q(10)|Homogenised succion|
|Q(11 + (i-1)*5)|$X_i$|
|Q(11 + (i-1)*5 +1)|$Y_i$|
|Q(11 + (i-1)*5 +2)|$P_{w,i}$|
|Q(11 + (i-1)*5 +3)|$C_i$|
|Q(11 + (i-1)*5 +4)|$P_{g,i}$|