GENERAL Mode

Governing Equations

The GENERAL mode involves two phase liquid water-gas flow coupled to the reactive transport mode. Mass conservation equations have the form

(1)\[\frac{{{\partial}}}{{{\partial}}t} \porosity \Big(s_l^{} \density_l^{} x_i^l + \saturation_g^{} \density_g^{} x_i^g \Big) + {\boldsymbol{\nabla}}\cdot\Big({\boldsymbol{q}}_l^{} \density_l^{} x_i^l + {\boldsymbol{q}}_g \density_g^{} x_i^g -\porosity \saturation_l^{} D_l^{} \density_l^{} {\boldsymbol{\nabla}}x_i^l -\porosity \saturation_g^{} D_g^{} \density_g^{} {\boldsymbol{\nabla}}x_i^g \Big) = Q_i^{},\]

for liquid and gas saturation \(s_{l,\,g}^{}\), density \(\density_{l,\,g}^{}\), diffusivity \(D_{l,\,g}^{}\), Darcy velocity \({\boldsymbol{q}}_{l,\,g}^{}\) and mole fraction \(x_i^{l,\,g}\). The energy conservation equation can be written in the form

(2)\[\sum_{{{\alpha}}=l,\,g}\left\{\frac{{{\partial}}}{{{\partial}}t} \big(\porosity \saturation_{{\alpha}}\density_{{\alpha}}U_{{\alpha}}\big) + {\boldsymbol{\nabla}}\cdot\big({\boldsymbol{q}}_{{\alpha}}\density_{{\alpha}}H_{{\alpha}}\big) \right\} + \frac{{{\partial}}}{{{\partial}}t}\big( (1-\porosity)\density_r C_p T \big) - {\boldsymbol{\nabla}}\cdot (\kappa{\boldsymbol{\nabla}}T) = Q,\]

as the sum of contributions from liquid and gas fluid phases and rock, with internal energy \(U_{{\alpha}}\) and enthalpy \(H_{{\alpha}}\) of fluid phase \({{\alpha}}\), rock heat capacity \(C_p\) and thermal conductivity \(\kappa\). Note that

(3)\[U_{{\alpha}}= H_{{\alpha}}-\frac{P_{{\alpha}}}{\density_{{\alpha}}}.\]

Thermal conductivity \(\kappa\) is determined from THERMAL_CHARACTERISTIC_CURVES, where the DEFAULT option uses the equation from Somerton et al., 1974:

(4)\[\kappa = \kappa_{\rm dry} + \sqrt{s_l^{}} (\kappa_{\rm sat} - \kappa_{\rm dry}),\]

where \(\kappa_{\rm dry}\) and \(\kappa_{\rm sat}\) are dry and fully saturated rock thermal conductivities, respectively.

The Darcy velocity of the \(\alpha^{th}\) phase is equal to

(5)\[\boldsymbol{q}_\alpha = -\frac{k k^{r}_{\alpha}}{\mu_\alpha} \boldsymbol{\nabla} (p_\alpha - \gamma_\alpha \boldsymbol{g} z), \ \ \ (\alpha=l,g),\]

where \(\boldsymbol{g}\) denotes the acceleration of gravity, \(k\) denotes the saturated permeability, \(k^{r}_{\alpha}\) the relative permeability, \(\mu_\alpha\) the viscosity, \(p_\alpha\) the pressure of the \(\alpha^{th}\) fluid phase, and

(6)\[\gamma_\alpha^{} = W_\alpha^{} \density_\alpha^{},\]

with \(W_\alpha\) the gram formula weight of the \(\alpha^{th}\) phase

(7)\[W_\alpha = \sum_{i=w,\,a} W_i^{} x_i^\alpha,\]

where \(W_i\) refers to the formula weight of the \(i^{th}\) component.

Capillary Pressure - Saturation Functions

Capillary pressure is related to effective liquid saturation by the van Genuchten and Brooks-Corey relations, as described under the sections van Genuchten Saturation Function and Brooks-Corey Saturation Function under RICHARDS Mode. Because both a liquid (wetting) and gas (non-wetting) phase are considered, the effective saturation \(s_e\) in the van Genuchten and Brooks-Corey relations under RICHARDS Mode becomes the effective liquid saturation \(s_{el}\) in the multiphase formulation. Liquid saturation \(s_l\) is obtained from the effective liquid saturation by

(8)\[\saturation_{l} = \saturation_{el}s_0 - \saturation_{el}s_{rl} + \saturation_{rl},\]

where \(s_{rl}\) denotes the liquid residual saturation, and \(s_0\) denotes the maximum liquid saturation. The gas saturation can be obtained from the relation

(9)\[\saturation_l + \saturation_g = 1\]

The effective gas saturation \(s_{eg}\) is defined by the relation

(10)\[\saturation_{eg} = 1 - \frac{s_l-s_{rl}}{1-s_{rl}-s_{rg}}\]

Additionally, a linear relationship between capillary pressure \(p_c\) and effective liquid saturation can be described as

(11)\[\saturation_{el} = {{p_c-p_c^{max}}\over{\frac{1}{\alpha}-p_c^{max}}}\]

where \(\alpha\) is a fitting parameter representing the air entry pressure [Pa]. The inverse relationship for capillary pressure is

(12)\[p_c = \left({\frac{1}{\alpha}-p_c^{max}}\right)s_{el} + p_c^{max}\]

Relative Permeability Functions

Two forms of each relative permeability function are implemented based on the Mualem and Burdine formulations as in RICHARDS Mode, but the effective liquid saturation \(s_{el}\) and the effective gas saturation \(s_{eg}\) are used. A summary of the relationships used can be found in Chen et al. (1999), where the tortuosity \(\eta\) is set to \(1/2\). The implemented relative permeability functions include: Mualem-van Genuchten, Mualem-Brooks-Corey, Mualem-linear, Burdine-van Genuchten, Burdine-Brooks-Corey, and Burdine-linear. For each relationship, the following definitions apply:

\[ \begin{align}\begin{aligned}S_{el} = \frac{S_{l}-S_{rl}}{1-S_{rl}}\\S_{eg} = \frac{S_{l}-S_{rl}}{1-S_{rl}-S_{rg}}\end{aligned}\end{align} \]

For the Mualem relative permeability function based on the van Genuchten saturation function, the liquid and gas relative permeability functions are given by the expressions

(13)\[ \begin{align}\begin{aligned}k^{r}_{l} =& \sqrt{s_{el}} \left\{1 - \left[1- \left( \saturation_{el} \right)^{1/m} \right]^m \right\}^2\\k^{r}_{g} =& \sqrt{1-s_{eg}} \left\{1 - \left( \saturation_{eg} \right)^{1/m} \right\}^{2m}.\end{aligned}\end{align} \]

For the Mualem relative permeability function based on the Brooks-Corey saturation function, the liquid and gas relative permeability functions are given by the expressions

(14)\[ \begin{align}\begin{aligned}k^{r}_{l} =& \big(s_{el}\big)^{5/2+2/\lambda}\\k^{r}_{g} =& \sqrt{1-s_{eg}}\left({1-s_{eg}^{1+1/\lambda}}\right)^{2}.\end{aligned}\end{align} \]

For the Mualem relative permeability function based on the linear saturation functions, the liquid and gas relative permeability functions are given by the expressions

(15)\[ \begin{align}\begin{aligned}k^{r}_{l} =& \sqrt{s_{el}}\frac{\ln\left({p_c/p_c^{max}}\right)}{\ln\left({\frac{1}{\alpha}/p_c^{max}}\right)}\\k^{r}_{g} =& \sqrt{1-s_{eg}}\left({1-\frac{k^{r}_{l}}{\sqrt{s_{eg}}}}\right)\end{aligned}\end{align} \]

For the Burdine relative permeability function based on the van Genuchten saturation function, the liquid and gas relative permeability functions are given by the expressions

(16)\[ \begin{align}\begin{aligned}k^{r}_{l} =& \saturation_{el}^2 \left\{1 - \left[1- \left( \saturation_{el} \right)^{1/m} \right]^m \right\}\\k^{r}_{g} =& (1-s_{eg})^2 \left\{1 - \left( \saturation_{eg} \right)^{1/m} \right\}^{m}.\end{aligned}\end{align} \]

For the Burdine relative permeability function based on the Brooks-Corey saturation function, the liquid and gas relative permeability functions have the form

(17)\[ \begin{align}\begin{aligned}k^{r}_{l} =& \big(s_{el}\big)^{3+2/\lambda}\\k^{r}_{g} =& (1-s_{eg})^2\left[{1-(s_{eg})^{1+2/\lambda}}\right].\end{aligned}\end{align} \]

For the Burdine relative permeability function based on the linear saturation functions, the liquid and gas relative permeability functions are given by the expressions

(18)\[ \begin{align}\begin{aligned}k^{r}_{l} =& \saturation_{el}\\k^{r}_{g} =& 1 - \saturation_{eg}.\end{aligned}\end{align} \]

Kelvin’s Equation for Vapor Pressure Lowering

Vapor pressure lowering resulting from capillary suction is described by Kelvin’s equation given by

(19)\[p_v = p_{\rm sat} (T) e^{-p_c/\density_l RT},\]

where \(p_v\) represents the vapor pressure, \(p_{\rm sat}\) the saturation pressure of pure water, \(p_c\) capillary pressure, \(\density_l\) liquid mole density, \(T\) denotes the temperature, and \(R\) the gas constant.

Fully Implicit Solute Coupling

A third mass conservation equation for a solute can additionally be solved:

(20)\[\frac{{{\partial}}}{{{\partial}}t} \porosity \Big(s_l \density_l x_i^l + \saturation_p \density_p x_i^p\Big) + {\boldsymbol{\nabla}}\cdot\Big({\boldsymbol{q}}_l\density_l x_i^l -\porosity \saturation_l D_l \density_l {\boldsymbol{\nabla}}x_i^l\Big) = Q_i\]

for liquid and solid precipitate saturation \(s_{l,\,p}\), density \(\density_{l,\,p}\), diffusivity \(D_{l}\), Darcy velocity \({\boldsymbol{q}}_{l}\) and mole fraction \(x_i^{l,\,p}\).

At solute concentrations above solubility, a solid precipitate phase forms. Transport of the solute is only possible through advection or diffusion within the liquid phase.

Alternatively, the rock matrix can be soluble, in which case the mass conservation equation becomes

(21)\[\frac{{{\partial}}}{{{\partial}}t} \Big(\porosity (s_l \density_l x_i^l)+(1-\porosity)\Big) + {\boldsymbol{\nabla}}\cdot\Big({\boldsymbol{q}}_l \density_l x_i^l -\porosity \saturation_l D_l \density_l {\boldsymbol{\nabla}}x_i^l\Big) = Q_i\]