Cold CO injection into an elastic reservoir - a multi-phase THM problem

Fluid flow: mass

In PorousFlow, the mass of species per volume of rock is written as a sum over all phases present in the system: (1) LaForce et al. make the following assumptions

• The porosity, , is constant. It is independent of time, fluid pressure, temperature and volumetric expansion. However, it is dependent on rock type.

  [porosity_reservoir]
    type = PorousFlowPorosityConst
    porosity = 0.2
  []

• There are two phases, liquid and gas. Hence .

• The liquid is water and the gas is CO. Hence . The water component exists only in the liquid phase and the CO component exists only in the gas phase. Hence , , and .

  [massfrac_ph0_sp0]
    initial_condition = 1 # all H20 in phase=0
  []
  [massfrac_ph1_sp0]
    initial_condition = 0 # no H2O in phase=1
  []

• The density of each phase is constant. This is in order to simplify the fluid equations, which eventually reduce to a well-studied form that models advection of a saturation shock front. However, in order to derive solutions for the pressure distribution, LaForce et al. assume a small fluid compressibility, which is the same for each phase. To satisfy both the incompressibility and small compressibility requirements, the MOOSE model uses a fluid bulk modulus that is larger than any porepressures.

• There is no desorbed fluid.

With these assumptions, the mass formula becomes

The PorousFlow input file uses the porepressure of water and the gas saturation as its primary nonlinear variables describing fluid flow:

  [pwater]
    initial_condition = 18.3e6
  []
  [sgas]
    initial_condition = 0.0
  []

Fluid flow: flux

In PorousFlow, the flux is a sum of advective flux and diffusive-and-dispersive flux: (2) Advective flux is governed by Darcy's law: (3) Diffusion and dispersion are proportional to the gradient of : (4)

LaForce et al. make the following assumptions

• The permeability tensor is diagonal and constant (independent of fluid pressure, temperature and rock stress and strain). Its and components are equal, and its component is zero. It is dependent on rock type. Since the PorousFlow input file is setup in 2D RZ coordinates, the radial permeability enters the first "slot", the component appears in the central slot, and all other components are zero in this situation

  [permeability_reservoir]
    type = PorousFlowPermeabilityConst
    permeability = '2e-12 0 0  0 0 0  0 0 0'
  []

• There is a constant and nonzero residual water saturation and constant and nonzero residual CO saturation .

• The relative permeability functions are functions of the effective saturation and the functions are (17)

  [relperm_liquid]
    type = PorousFlowRelativePermeabilityCorey
    n = 4
    phase = 0
    s_res = 0.200
    sum_s_res = 0.405
  []
  [relperm_gas]
    type = PorousFlowRelativePermeabilityBC
    phase = 1
    s_res = 0.205
    sum_s_res = 0.405
    nw_phase = true
    lambda = 2
  []

• The viscosity is a function of temperature only. LaForce et al. use a particular function of temperature that is not included in PorousFlow. Therefore, here we shall use constant viscosity for each phase that is equal to the average viscosity quoted in LaForce et al., and use LaForce et al.'s solution based on that constant viscosity. This results in only very tiny modifications to LaForce's solutions which would be difficult to see in the finite-resolution of the MOOSE model anyway.

• There is no capillarity. Therefore . Define .

  [pc]
    type = PorousFlowCapillaryPressureConst
    pc = 0

• There is no gravity.

• There is no diffusion and dispersion (which is implicitly implied by the assumptions on ).

With these assumptions, the fluid fluxes are

Fluid flow: mass conservation

In PorousFlow, mass conservation for fluid species is described by the continuity equation (5)

LaForce et al. make the following assumptions:

• The term is ignored.

• There is no radioactive decay, .

• Sources and sinks only occur on the boundary, so except for on the boundary.

Therefore, the fluid-flow Kernels are

  [mass_water_dot]
    type = PorousFlowMassTimeDerivative
    fluid_component = 0
    variable = pwater
  []
  [flux_water]
    type = PorousFlowAdvectiveFlux
    fluid_component = 0
    use_displaced_mesh = false
    variable = pwater
  []
  [mass_co2_dot]
    type = PorousFlowMassTimeDerivative
    fluid_component = 1
    variable = sgas
  []
  [flux_co2]
    type = PorousFlowAdvectiveFlux
    fluid_component = 1
    use_displaced_mesh = false
    variable = sgas
  []

Fluid flow: summary

With all the assumptions so far, the mass conservation equations read The sum of these yields the total flow rate: This total flow rate is dependent on time (controlled by the boundary conditions) but satisfies .

For later use, define the fractional flows (6)

Heat flow: heat energy

In PorousFlow, the heat energy per unit volume is (7)

LaForce et al. make the following assumptions

• The rock grain density, , is constant.

• The specific heat capacity of the rock grains, , is constant.

  [internal_energy_reservoir]
    type = PorousFlowMatrixInternalEnergy
    specific_heat_capacity = 1100
    density = 2350.0
  []

• The internal energy of the liquid phase is , where is constant and independent of all other parameters in the system.

• The internal energy of the gas phase is , where is constant and independent of all other parameters of the system.

With these assumptions the heat energy density is

Heat flow: flux

In PorousFlow, the heat flux is a sum of heat conduction and convection with the fluid: (8)

LaForce et al. make the following assumptions

• The only non-zero component of the thermal conductivity is the component. It varies with water saturation The coefficients and are independent of time, fluid pressure, temperature and mechanical deformation, but they are dependent on rock type.

  [thermal_conductivity_reservoir]
    type = PorousFlowThermalConductivityIdeal
    dry_thermal_conductivity = '0 0 0  0 1.320 0  0 0 0'
    wet_thermal_conductivity = '0 0 0  0 3.083 0  0 0 0'
  []

• The enthalpy of the liquid phase is .

• The enthalpy of the gas phase is .

With the assumptions made so far where the fractional flows of Eq. (6) have been employed.

Heat flow: energy conservation

In PorousFlow, energy conservation for heat is described by the continuity equation (9)

LaForce et al. make the following assumptions

• There is no plastic heating.

• There are no heat sources, save for on the boundary.

Heat flow: summary

With the assumptions made so far, the heat equation reads (10) where the following two variables, and have been defined: (11) Eq. (10) and Eq. (11) are exactly Eq. (A1) in LaForce et al. (2014a).

During the course of their analysis, LaForce et al. replace the diffusive term, with a Lauwerier heat-loss expression: (12) In terms of LaForce et al.'s variables, the new heat-loss coefficient is (13) This is a particularly nice simplification because it means that the overburden and underburden don't need to be considered any further in the analysis. It does mean, however, that the mechanical effects of these surrounding rocks on the aquifer are not modelled by LaForce et al.

Therefore, the heat-flow Kernels are

  [energy_dot]
    type = PorousFlowEnergyTimeDerivative
    variable = temp
  []
  [advection]
    type = PorousFlowHeatAdvection
    use_displaced_mesh = false
    variable = temp
  []
  [conduction]
    type = PorousFlowExponentialDecay
    use_displaced_mesh = false
    variable = temp
    reference = 358
    rate = rate
  []

along with the rate of heat loss:

[Functions<<<{"href": "../../syntax/Functions/index.html"}>>>]
  [decay_rate]
    # Eqn(26) of the first paper of LaForce et al.
    # Ka * (rho C)_a = 10056886.914
    # h = 11
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = 'sqrt(10056886.914/t)/11.0'
  []
[]

Solid mechanics

Denote the effective stress tensor by . In PorousFlow it is defined by (14)

LaForce et al. assume

• The Biot coefficient, , is constant and independent of all other parameters in the theory.

• The effective fluid pressure is

  [eff_fluid_pressure]
    type = PorousFlowEffectiveFluidPressure
  []

In PorousFlow, the solid-mechanical constitutive law reads (15)

LaForce et al. assume

• The drained elasticity tensor, , is isotropic, constant and independent of all other parameters in the theory.

  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    shear_modulus = 6.0E9
    poissons_ratio = 0.2
  []

• There is no plasticity.

• The drained thermal expansion coefficient, , is constant and independent of all other parameters in the theory.

  [strain]
    type = ComputeAxisymmetricRZSmallStrain
    eigenstrain_names = 'thermal_contribution ini_stress'
  []
  [ini_strain]
    type = ComputeEigenstrainFromInitialStress
    initial_stress = '-12.8E6 0 0  0 -51.3E6 0  0 0 -12.8E6'
    eigenstrain_name = ini_stress
  []
  [thermal_contribution]
    type = ComputeThermalExpansionEigenstrain
    temperature = temp
    stress_free_temperature = 358
    thermal_expansion_coeff = 5E-6
    eigenstrain_name = thermal_contribution
  []
  [stress]
    type = ComputeLinearElasticStress
  []

In PorousFlow, conservation of momentum reads (16)

LaForce et al. assume

• The acceleration of the solid skeleton, , is zero.

• There are no body forces, , except for on the boundary.

Further assumptions and implications

LaForce et al. make the following assumptions:

• There is no displacement in the direction, . The solid-mechanics Kernels are therefore

  [grad_stress_r]
    type = StressDivergenceRZTensors
    temperature = temp
    eigenstrain_names = thermal_contribution
    variable = disp_r
    use_displaced_mesh = false
    component = 0

• The reservoir is free to contract and expand radially without being hindered by the impermeable seals.

There is no dependence on the axial coordinates, , because of the use of the Lauwerier heat-loss model (Eq. (12)) that simulates axial heat conduction. Therefore, all independent variables — the water saturation , fluid porepressure , temperature , and radial displacement — are dependent on radius only.

The Variables used in the PorousFlow input file are:

[Variables<<<{"href": "../../syntax/Variables/index.html"}>>>]
  [pwater]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 18.3e6
  []
  [sgas]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 0.0
  []
  [temp]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 358
  []
  [disp_r]
  []
[]