www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
PorousFlowBrineCO2 Class Reference

Specialized class for brine and CO2 including calculation of mutual solubility of the two fluids using a high-accuracy fugacity-based formulation. More...

#include <PorousFlowBrineCO2.h>

Inheritance diagram for PorousFlowBrineCO2:
[legend]

Public Member Functions

 PorousFlowBrineCO2 (const InputParameters &parameters)
 
virtual std::string fluidStateName () const override
 Name of FluidState. More...
 
virtual void thermophysicalProperties (Real pressure, Real temperature, Real Xnacl, Real Z, unsigned int qp, std::vector< FluidStateProperties > &fsp) const override
 Determines the complete thermophysical state of the system for a given set of primary variables. More...
 
void equilibriumMoleFractions (Real pressure, Real temperature, Real Xnacl, Real &xco2, Real &dxco2_dp, Real &dxco2_dT, Real &dxco2_dX, Real &yh2o, Real &dyh2o_dp, Real &dyh2o_dT, Real &dyh2o_dX) const
 Mole fractions of CO2 in brine and water vapor in CO2 at equilibrium. More...
 
void equilibriumMassFractions (Real pressure, Real temperature, Real Xnacl, Real &Xco2, Real &dXco2_dp, Real &dXco2_dT, Real &dXco2_dX, Real &Yh2o, Real &dYh2o_dp, Real &dYh2o_dT, Real &dYh2o_dX) const
 Mass fractions of CO2 in brine and water vapor in CO2 at equilibrium. More...
 
void massFractions (Real pressure, Real temperature, Real Xnacl, Real Z, FluidStatePhaseEnum &phase_state, std::vector< FluidStateProperties > &fsp) const
 Mass fractions of CO2 and H2O in both phases, as well as derivatives wrt PorousFlow variables. More...
 
void gasProperties (Real pressure, Real temperature, std::vector< FluidStateProperties > &fsp) const
 Thermophysical properties of the gaseous state. More...
 
void liquidProperties (Real pressure, Real temperature, Real Xnacl, std::vector< FluidStateProperties > &fsp) const
 Thermophysical properties of the liquid state. More...
 
void saturationTwoPhase (Real pressure, Real temperature, Real Xnacl, Real Z, std::vector< FluidStateProperties > &fsp) const
 Gas and liquid saturations for the two-phase region. More...
 
void fugacityCoefficientsLowTemp (Real pressure, Real temperature, Real co2_density, Real dco2_density_dp, Real dco2_density_dT, Real &fco2, Real &dfco2_dp, Real &dfco2_dT, Real &fh2o, Real &dfh2o_dp, Real &dfh2o_dT) const
 Fugacity coefficients for H2O and CO2 for T <= 100C Eq. More...
 
Real activityCoefficientH2O (Real temperature, Real xco2) const
 Activity coefficient of H2O Eq. More...
 
Real activityCoefficientCO2 (Real temperature, Real xco2) const
 Activity coefficient of CO2 Eq. More...
 
void activityCoefficient (Real pressure, Real temperature, Real Xnacl, Real &gamma, Real &dgamma_dp, Real &dgamma_dT, Real &dgamma_dX) const
 Activity coefficient for CO2 in brine. More...
 
void activityCoefficientHighTemp (Real temperature, Real Xnacl, Real &gamma, Real &dgamma_dT, Real &dgamma_dX) const
 Activity coefficient for CO2 in brine used in the elevated temperature formulation. More...
 
void equilibriumConstantH2O (Real temperature, Real &kh2o, Real &dkh2o_dT) const
 Equilibrium constant for H2O For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2003) For temperatures 109C <= T <= 300C, uses Spycher and Pruess (2010) For temperatures 99C < T < 109C, the value is calculated by smoothly interpolating the two formulations. More...
 
void equilibriumConstantCO2 (Real temperature, Real &kco2, Real &dkco2_dT) const
 Equilibrium constant for CO2 For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2003) For temperatures 109C <= T <= 300C, uses Spycher and Pruess (2010) For temperatures 99C < T < 109C, the value is calculated by smoothly interpolating the two formulations. More...
 
void partialDensityCO2 (Real temperature, Real &partial_density, Real &dpartial_density_dT) const
 Partial density of dissolved CO2 From Garcia, Density of aqueous solutions of CO2, LBNL-49023 (2001) More...
 
virtual Real totalMassFraction (Real pressure, Real temperature, Real Xnacl, Real saturation, unsigned int qp) const override
 Total mass fraction of fluid component summed over all phases in the two-phase state for a specified gas saturation. More...
 
unsigned int saltComponentIndex () const
 The index of the salt component. More...
 
void henryConstant (Real temperature, Real Xnacl, Real &Kh, Real &dKh_dT, Real &dKh_dx) const
 Henry's constant of dissolution of gas phase CO2 in brine. More...
 
void enthalpyOfDissolutionGas (Real temperature, Real Xnacl, Real &hdis, Real &dhdis_dT, Real &dhdis_dx) const
 Enthalpy of dissolution of gas phase CO2 in brine calculated using Henry's constant From Himmelblau, Partial molal heats and entropies of solution for gases dissolved in water from the freezing to the near critical point, J. More...
 
void enthalpyOfDissolution (Real temperature, Real &hdis, Real &dhdis_dT) const
 Enthalpy of dissolution of CO2 in brine calculated using linear fit to model of Duan and Sun, An improved model calculating CO2 solubility in pure water and aqueous NaCl solutions from 273 to 533 K and from 0 to 2000 bar, Chemical geology, 193, 257–271 (2003). More...
 
void equilibriumMoleFractionsLowTemp (Real pressure, Real temperature, Real Xnacl, Real &xco2, Real &dxco2_dp, Real &dxco2_dT, Real &dxco2_dX, Real &yh2o, Real &dyh2o_dp, Real &dyh2o_dT, Real &dyh2o_dX) const
 Mole fractions of CO2 in brine and water vapor in CO2 at equilibrium in the low temperature regime (T <= 99C). More...
 
void solveEquilibriumMoleFractionHighTemp (Real pressure, Real temperature, Real Xnacl, Real co2_density, Real &xco2, Real &yh2o) const
 Function to solve for yh2o and xco2 iteratively in the elevated temperature regime (T > 100C) More...
 
unsigned int numPhases () const
 The maximum number of phases in this model. More...
 
unsigned int numComponents () const
 The maximum number of components in this model. More...
 
unsigned int aqueousPhaseIndex () const
 The index of the aqueous phase. More...
 
unsigned int gasPhaseIndex () const
 The index of the gas phase. More...
 
unsigned int aqueousComponentIndex () const
 The index of the aqueous fluid component. More...
 
unsigned int gasComponentIndex () const
 The index of the gas fluid component. More...
 
void clearFluidStateProperties (std::vector< FluidStateProperties > &fsp) const
 Clears the contents of the FluidStateProperties data structure. More...
 
void initialize () final
 
void execute () final
 
void finalize () final
 
Real rachfordRice (Real vf, std::vector< Real > &Zi, std::vector< Real > &Ki) const
 Rachford-Rice equation for vapor fraction. More...
 
Real rachfordRiceDeriv (Real vf, std::vector< Real > &Zi, std::vector< Real > &Ki) const
 Derivative of Rachford-Rice equation wrt vapor fraction. More...
 
Real vaporMassFraction (Real Z0, Real K0, Real K1) const
 Solves Rachford-Rice equation to provide vapor mass fraction. More...
 
Real vaporMassFraction (std::vector< Real > &Zi, std::vector< Real > &Ki) const
 
void fugacityCoefficientsHighTemp (Real pressure, Real temperature, Real co2_density, Real xco2, Real yh2o, Real &fco2, Real &fh2o) const
 Fugacity coefficients for H2O and CO2 at elevated temperatures (100C < T <= 300C). More...
 
void fugacityCoefficientsHighTemp (Real pressure, Real temperature, Real co2_density, Real xco2, Real yh2o, Real &fco2, Real &dfco2_dp, Real &dfco2_dT, Real &fh2o, Real &dfh2o_dp, Real &dfh2o_dT) const
 
Real fugacityCoefficientH2OHighTemp (Real pressure, Real temperature, Real co2_density, Real xco2, Real yh2o) const
 
Real fugacityCoefficientCO2HighTemp (Real pressure, Real temperature, Real co2_density, Real xco2, Real yh2o) const
 

Protected Member Functions

void checkVariables (Real pressure, Real temperature) const
 Check the input variables. More...
 
void smoothCubicInterpolation (Real temperature, Real f0, Real df0, Real f1, Real df1, Real &value, Real &deriv) const
 Cubic function to smoothly interpolate between the low temperature and elevated temperature models for 99C < T < 109C. More...
 
void phaseState (Real Zi, Real Xi, Real Yi, FluidStatePhaseEnum &phase_state) const
 Determines the phase state gven the total mass fraction and equilibrium mass fractions. More...
 
void funcABHighTemp (Real pressure, Real temperature, Real Xnacl, Real co2_density, Real xco2, Real yh2o, Real &A, Real &B) const
 The function A (Eq. More...
 
void funcABHighTemp (Real pressure, Real temperature, Real Xnacl, Real co2_density, Real xco2, Real yh2o, Real &A, Real &dA_dp, Real &dA_dT, Real &B, Real &dB_dp, Real &dB_dT, Real &dB_dX) const
 
void funcABLowTemp (Real pressure, Real temperature, Real co2_density, Real dco2_density_dp, Real dco2_density_dT, Real &A, Real &dA_dp, Real &dA_dT, Real &B, Real &dB_dp, Real &dB_dT) const
 

Protected Attributes

const unsigned int _salt_component
 Salt component index. More...
 
const BrineFluidProperties_brine_fp
 Fluid properties UserObject for water. More...
 
const SinglePhaseFluidProperties_co2_fp
 Fluid properties UserObject for the CO2. More...
 
const SinglePhaseFluidProperties_water_fp
 Fluid properties UserObject for H20. More...
 
const Real _Mh2o
 Molar mass of water (kg/mol) More...
 
const Real _invMh2o
 Inverse of molar mass of H2O (mol/kg) More...
 
const Real _Mco2
 Molar mass of CO2 (kg/mol) More...
 
const Real _Mnacl
 Molar mass of NaCL. More...
 
const Real _Rbar
 Molar gas constant in bar cm^3 /(K mol) More...
 
const Real _Tlower
 Temperature below which the Spycher, Pruess & Ennis-King (2003) model is used (K) More...
 
const Real _Tupper
 Temperature above which the Spycher & Pruess (2010) model is used (K) More...
 
const Real _Zmin
 Minimum Z - below this value all CO2 will be dissolved. More...
 
unsigned int _num_phases
 Number of phases. More...
 
unsigned int _num_components
 Number of components. More...
 
const unsigned int _aqueous_phase_number
 Phase number of the aqueous phase. More...
 
unsigned int _gas_phase_number
 Phase number of the gas phase. More...
 
const unsigned int _aqueous_fluid_component
 Fluid component number of the aqueous component. More...
 
unsigned int _gas_fluid_component
 Fluid component number of the gas phase. More...
 
const Real _R
 Universal gas constant (J/mol/K) More...
 
const Real _T_c2k
 Conversion from C to K. More...
 
const Real _nr_max_its
 Maximum number of iterations for the Newton-Raphson iterations. More...
 
const Real _nr_tol
 Tolerance for Newton-Raphson iterations. More...
 
const PorousFlowCapillaryPressure_pc
 Capillary pressure UserObject. More...
 

Detailed Description

Specialized class for brine and CO2 including calculation of mutual solubility of the two fluids using a high-accuracy fugacity-based formulation.

For temperatures 12C <= T <= 99C, the formulation is based on Spycher, Pruess and Ennis-King, CO2-H2O mixtures in the geological sequestration of CO2. I. Assessment and calculation of mutual solubilities from 12 to 100C and up to 600 bar, Geochimica et Cosmochimica Acta, 67, 3015-3031 (2003), and Spycher and Pruess, CO2-H2O mixtures in the geological sequestration of CO2. II. Partitioning in chloride brine at 12-100C and up to 600 bar, Geochimica et Cosmochimica Acta, 69, 3309-3320 (2005).

For temperatures 109C <= T <= 300C, the formulation is based on Spycher and Pruess, A Phase-Partitioning Model for CO2-Brine Mixtures at Elevated Temperatures and Pressures: Application to CO2-Enhanced Geothermal Systems, Transport in Porous Media, 82, 173-196 (2010)

As the two formulations do not coincide at temperatures near 100C, a cubic polynomial is used in the intermediate temperature range 99C < T < 109C to provide a smooth transition from the two formulations in this region.

Notation convention Throughout this class, both mole fractions and mass fractions will be used. The following notation will be used: yk: mole fraction of component k in the gas phase xk: mole fraction of component k in the liquid phase Yk: mass fraction of component k in the gas phase Xk: mass fraction of component k in the liquid phase

Definition at line 51 of file PorousFlowBrineCO2.h.

Constructor & Destructor Documentation

◆ PorousFlowBrineCO2()

PorousFlowBrineCO2::PorousFlowBrineCO2 ( const InputParameters &  parameters)

Definition at line 30 of file PorousFlowBrineCO2.C.

31  : PorousFlowFluidStateBase(parameters),
32  _salt_component(getParam<unsigned int>("salt_component")),
33  _brine_fp(getUserObject<BrineFluidProperties>("brine_fp")),
34  _co2_fp(getUserObject<SinglePhaseFluidProperties>("co2_fp")),
37  _invMh2o(1.0 / _Mh2o),
40  _Rbar(_R * 10.0),
41  _Tlower(372.15),
42  _Tupper(382.15),
43  _Zmin(1.0e-4)
44 {
45  // Check that the correct FluidProperties UserObjects have been provided
46  if (_co2_fp.fluidName() != "co2")
47  paramError("co2_fp", "A valid CO2 FluidProperties UserObject must be provided");
48 
49  if (_brine_fp.fluidName() != "brine")
50  paramError("brine_fp", "A valid Brine FluidProperties UserObject must be provided");
51 
52  // Set the number of phases and components, and their indexes
53  _num_phases = 2;
54  _num_components = 3;
57 
58  // Check that _aqueous_phase_number is <= total number of phases
60  paramError("liquid_phase_number",
61  "This value is larger than the possible number of phases ",
62  _num_phases);
63 
64  // Check that _aqueous_fluid_component is <= total number of fluid components
66  paramError("liquid_fluid_component",
67  "This value is larger than the possible number of fluid components",
69 
70  // Check that the salt component index is not identical to the liquid fluid component
72  paramError(
73  "salt_component",
74  "The value provided must be different from the value entered in liquid_fluid_component");
75 
76  // Check that _salt_component is <= total number of fluid components
78  paramError("salt_component",
79  "The value provided is larger than the possible number of fluid components",
81 }
const Real _Tlower
Temperature below which the Spycher, Pruess & Ennis-King (2003) model is used (K) ...
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
unsigned int _num_phases
Number of phases.
const unsigned int _salt_component
Salt component index.
const Real _Mh2o
Molar mass of water (kg/mol)
const SinglePhaseFluidProperties & _water_fp
Fluid properties UserObject for H20.
const Real _Mnacl
Molar mass of NaCL.
const Real _invMh2o
Inverse of molar mass of H2O (mol/kg)
virtual std::string fluidName() const override
Fluid name.
virtual Real molarMass() const
Molar mass [kg/mol].
unsigned int _gas_phase_number
Phase number of the gas phase.
unsigned int _num_components
Number of components.
const unsigned int _aqueous_fluid_component
Fluid component number of the aqueous component.
const Real _Mco2
Molar mass of CO2 (kg/mol)
Real molarMassH2O() const
H2O molar mass.
const BrineFluidProperties & _brine_fp
Fluid properties UserObject for water.
const SinglePhaseFluidProperties & _co2_fp
Fluid properties UserObject for the CO2.
Real molarMassNaCl() const
NaCl molar mass.
virtual const SinglePhaseFluidProperties & getComponent(unsigned int component) const override
Get UserObject for specified component.
const Real _R
Universal gas constant (J/mol/K)
unsigned int _gas_fluid_component
Fluid component number of the gas phase.
PorousFlowFluidStateBase(const InputParameters &parameters)
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
virtual std::string fluidName() const
Fluid name.
const Real _Zmin
Minimum Z - below this value all CO2 will be dissolved.
static const unsigned int WATER
Fluid component numbers for water and NaCl.
const Real _Tupper
Temperature above which the Spycher & Pruess (2010) model is used (K)

Member Function Documentation

◆ activityCoefficient()

void PorousFlowBrineCO2::activityCoefficient ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real &  gamma,
Real &  dgamma_dp,
Real &  dgamma_dT,
Real &  dgamma_dX 
) const

Activity coefficient for CO2 in brine.

From Duan and Sun, An improved model calculating CO2 solubility in pure water and aqueous NaCl solutions from 257 to 533 K and from 0 to 2000 bar, Chem. Geol., 193, 257-271 (2003)

Note: this is not a 'true' activity coefficient, and is instead related to the molality of CO2 in water and brine. Nevertheless, Spycher and Pruess (2005) refer to it as an activity coefficient, so this notation is followed here.

Parameters
pressurephase pressure (Pa)
temperaturephase temperature (K)
Xnaclsalt mass fraction (kg/kg)
[out]gammaactivity coefficient for CO2 in brine
[out]dgamma_dpderivative of activity coefficient wrt pressure
[out]dgamma_dTderivative of activity coefficient wrt temperature
[out]dgamma_dXderivative of activity coefficient wrt salt mass fraction

Definition at line 847 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractionsLowTemp().

854 {
855  // Need pressure in bar
856  const Real pbar = pressure * 1.0e-5;
857  // Need NaCl molality (mol/kg)
858  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
859  const Real dmnacl_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
860 
861  const Real lambda = -0.411370585 + 6.07632013e-4 * temperature + 97.5347708 / temperature -
862  0.0237622469 * pbar / temperature +
863  0.0170656236 * pbar / (630.0 - temperature) +
864  1.41335834e-5 * temperature * std::log(pbar);
865 
866  const Real xi = 3.36389723e-4 - 1.9829898e-5 * temperature + 2.12220830e-3 * pbar / temperature -
867  5.24873303e-3 * pbar / (630.0 - temperature);
868 
869  gamma = std::exp(2.0 * lambda * mnacl + xi * mnacl * mnacl);
870 
871  // Derivative wrt pressure
872  const Real dlambda_dp = -0.0237622469 / temperature + 0.0170656236 / (630.0 - temperature) +
873  1.41335834e-5 * temperature / pbar;
874  const Real dxi_dp = 2.12220830e-3 / temperature - 5.24873303e-3 / (630.0 - temperature);
875  dgamma_dp = (2.0 * mnacl * dlambda_dp + mnacl * mnacl * dxi_dp) * gamma * 1.0e-5;
876 
877  // Derivative wrt temperature
878  const Real dlambda_dT = 6.07632013e-4 - 97.5347708 / temperature / temperature +
879  0.0237622469 * pbar / temperature / temperature +
880  0.0170656236 * pbar / (630.0 - temperature) / (630.0 - temperature) +
881  1.41335834e-5 * std::log(pbar);
882  const Real dxi_dT = -1.9829898e-5 - 2.12220830e-3 * pbar / temperature / temperature -
883  5.24873303e-3 * pbar / (630.0 - temperature) / (630.0 - temperature);
884  dgamma_dT = (2.0 * mnacl * dlambda_dT + mnacl * mnacl * dxi_dT) * gamma;
885 
886  // Derivative wrt salt mass fraction
887  dgamma_dX = (2.0 * lambda * dmnacl_dX + 2.0 * xi * mnacl * dmnacl_dX) * gamma;
888 }
const Real _Mnacl
Molar mass of NaCL.
const std::string temperature
Definition: NS.h:27
const std::string pressure
Definition: NS.h:26

◆ activityCoefficientCO2()

Real PorousFlowBrineCO2::activityCoefficientCO2 ( Real  temperature,
Real  xco2 
) const

Activity coefficient of CO2 Eq.

(13) from Spycher & Pruess (2010)

Parameters
temperaturetemperature (K)
xco2mole fraction of CO2 in liquid phase (-)
Returns
activity coefficient

Definition at line 832 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp().

833 {
834  if (temperature <= 373.15)
835  return 1.0;
836  else
837  {
838  const Real Tref = temperature - 373.15;
839  const Real xh2o = 1.0 - xco2;
840  const Real Am = -3.084e-2 * Tref + 1.927e-5 * Tref * Tref;
841 
842  return std::exp(2.0 * Am * xco2 * xh2o * xh2o);
843  }
844 }
const std::string temperature
Definition: NS.h:27

◆ activityCoefficientH2O()

Real PorousFlowBrineCO2::activityCoefficientH2O ( Real  temperature,
Real  xco2 
) const

Activity coefficient of H2O Eq.

(12) from Spycher & Pruess (2010)

Parameters
temperaturetemperature (K)
xco2mole fraction of CO2 in liquid phase (-)
Returns
activity coefficient

Definition at line 817 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp().

818 {
819  if (temperature <= 373.15)
820  return 1.0;
821  else
822  {
823  const Real Tref = temperature - 373.15;
824  const Real xh2o = 1.0 - xco2;
825  const Real Am = -3.084e-2 * Tref + 1.927e-5 * Tref * Tref;
826 
827  return std::exp((Am - 2.0 * Am * xh2o) * xco2 * xco2);
828  }
829 }
const std::string temperature
Definition: NS.h:27

◆ activityCoefficientHighTemp()

void PorousFlowBrineCO2::activityCoefficientHighTemp ( Real  temperature,
Real  Xnacl,
Real &  gamma,
Real &  dgamma_dT,
Real &  dgamma_dX 
) const

Activity coefficient for CO2 in brine used in the elevated temperature formulation.

Eq. (18) from Spycher and Pruess (2010).

Note: unlike activityCoefficient(), no pressure dependence is included in this formulation

Parameters
temperaturephase temperature (K)
Xnaclsalt mass fraction (kg/kg)
[out]gammaactivity coefficient for CO2 in brine
[out]dgamma_dTderivative of activity coefficient wrt temperature
[out]dgamma_dXderivative of activity coefficient wrt salt mass fraction

Definition at line 891 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp().

893 {
894  // Need NaCl molality (mol/kg)
895  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
896  const Real dmnacl_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
897 
898  const Real T2 = temperature * temperature;
899  const Real T3 = temperature * T2;
900 
901  const Real lambda = 2.217e-4 * temperature + 1.074 / temperature + 2648.0 / T2;
902  const Real xi = 1.3e-5 * temperature - 20.12 / temperature + 5259.0 / T2;
903 
904  gamma = (1.0 - mnacl / _invMh2o) * std::exp(2.0 * lambda * mnacl + xi * mnacl * mnacl);
905 
906  // Derivative wrt temperature
907  const Real dlambda_dT = 2.217e-4 - 1.074 / T2 - 5296.0 / T3;
908  const Real dxi_dT = 1.3e-5 + 20.12 / T2 - 10518.0 / T3;
909  dgamma_dT = (2.0 * mnacl * dlambda_dT + mnacl * mnacl * dxi_dT) * gamma;
910 
911  // Derivative wrt salt mass fraction
912  dgamma_dX = (2.0 * lambda * dmnacl_dX + 2.0 * xi * mnacl * dmnacl_dX) * gamma -
913  dmnacl_dX / _invMh2o * std::exp(2.0 * lambda * mnacl + xi * mnacl * mnacl);
914 }
const Real _Mnacl
Molar mass of NaCL.
const Real _invMh2o
Inverse of molar mass of H2O (mol/kg)
const std::string temperature
Definition: NS.h:27

◆ aqueousComponentIndex()

unsigned int PorousFlowFluidStateBase::aqueousComponentIndex ( ) const
inlineinherited

The index of the aqueous fluid component.

Returns
aqueous fluid component number

Definition at line 117 of file PorousFlowFluidStateBase.h.

117 { return _aqueous_fluid_component; };
const unsigned int _aqueous_fluid_component
Fluid component number of the aqueous component.

◆ aqueousPhaseIndex()

unsigned int PorousFlowFluidStateBase::aqueousPhaseIndex ( ) const
inlineinherited

The index of the aqueous phase.

Returns
aqueous phase number

Definition at line 105 of file PorousFlowFluidStateBase.h.

105 { return _aqueous_phase_number; };
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.

◆ checkVariables()

void PorousFlowBrineCO2::checkVariables ( Real  pressure,
Real  temperature 
) const
protected

Check the input variables.

Definition at line 1553 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties(), and totalMassFraction().

1554 {
1555  // The calculation of mass fractions is valid from 12C <= T <= 300C, and
1556  // pressure less than 60 MPa
1557  if (temperature < 285.15 || temperature > 573.15)
1558  mooseException(name() + ": temperature " + Moose::stringify(temperature) +
1559  " is outside range 285.15 K <= T <= 573.15 K");
1560 
1561  if (pressure > 6.0e7)
1562  mooseException(name() + ": pressure " + Moose::stringify(pressure) +
1563  " must be less than 60 MPa");
1564 }
const std::string temperature
Definition: NS.h:27
const std::string name
Definition: Setup.h:22
const std::string pressure
Definition: NS.h:26

◆ clearFluidStateProperties()

void PorousFlowFluidStateBase::clearFluidStateProperties ( std::vector< FluidStateProperties > &  fsp) const
inherited

Clears the contents of the FluidStateProperties data structure.

Parameters
[out]fspFluidStateProperties data structure with all data initialized to 0

Definition at line 41 of file PorousFlowFluidStateBase.C.

Referenced by PorousFlowWaterNCG::thermophysicalProperties(), and thermophysicalProperties().

42 {
43  std::fill(fsp.begin(), fsp.end(), FluidStateProperties(_num_components));
44 }
unsigned int _num_components
Number of components.
Data structure to pass calculated thermophysical properties.

◆ enthalpyOfDissolution()

void PorousFlowBrineCO2::enthalpyOfDissolution ( Real  temperature,
Real &  hdis,
Real &  dhdis_dT 
) const

Enthalpy of dissolution of CO2 in brine calculated using linear fit to model of Duan and Sun, An improved model calculating CO2 solubility in pure water and aqueous NaCl solutions from 273 to 533 K and from 0 to 2000 bar, Chemical geology, 193, 257–271 (2003).

In the region of interest, the more complicated model given in Eq. (8) of Duan and Sun is well approximated by a simple linear fit (R^2 = 0.9997).

Note: as the effect of salt mass fraction is small, it is not included in this function.

Parameters
temperaturefluid temperature (K)
[out]hdisenthalpy of dissolution (J/kg)
[out]dhdis_dTderivative of enthalpy of dissolution wrt temperature

Definition at line 1523 of file PorousFlowBrineCO2.C.

Referenced by liquidProperties().

1524 {
1525  // Linear fit to model of Duan and Sun (2003) (in kJ/mol)
1526  const Real delta_h = -58.3533 + 0.134519 * temperature;
1527 
1528  // Convert to J/kg
1529  hdis = delta_h * 1000.0 / _Mco2;
1530  dhdis_dT = 134.519 / _Mco2;
1531 }
const Real _Mco2
Molar mass of CO2 (kg/mol)
const std::string temperature
Definition: NS.h:27

◆ enthalpyOfDissolutionGas()

void PorousFlowBrineCO2::enthalpyOfDissolutionGas ( Real  temperature,
Real  Xnacl,
Real &  hdis,
Real &  dhdis_dT,
Real &  dhdis_dx 
) const

Enthalpy of dissolution of gas phase CO2 in brine calculated using Henry's constant From Himmelblau, Partial molal heats and entropies of solution for gases dissolved in water from the freezing to the near critical point, J.

Phys. Chem. 63 (1959). Correction due to salinity from Battistelli et al, A fluid property module for the TOUGH2 simulator for saline brines with non-condensible gas, Proc. Eighteenth Workshop on Geothermal Reservoir Engineering (1993).

Parameters
temperaturefluid temperature (K)
XnaclNaCl mass fraction (kg/kg)
[out]hdisenthalpy of dissolution (J/kg)
[out]dhdis_dTderivative of enthalpy of dissolution wrt temperature
[out]dhdis_dxderivative of enthalpy of dissolution wrt salt mass fraction

Definition at line 1498 of file PorousFlowBrineCO2.C.

1500 {
1501  // Henry's constant
1502  Real Kh, dKh_dT, dKh_dX;
1503  henryConstant(temperature, Xnacl, Kh, dKh_dT, dKh_dX);
1504 
1505  hdis = -_R * temperature * temperature * dKh_dT / Kh / _Mco2;
1506 
1507  // Derivative of enthalpy of dissolution wrt temperature and xnacl requires the second
1508  // derivatives of Henry's constant. For simplicity, approximate these numerically
1509  const Real dT = temperature * 1.0e-8;
1510  const Real T2 = temperature + dT;
1511  henryConstant(T2, Xnacl, Kh, dKh_dT, dKh_dX);
1512 
1513  dhdis_dT = (-_R * T2 * T2 * dKh_dT / Kh / _Mco2 - hdis) / dT;
1514 
1515  const Real dX = Xnacl * 1.0e-8;
1516  const Real X2 = Xnacl + dX;
1517  henryConstant(temperature, X2, Kh, dKh_dT, dKh_dX);
1518 
1519  dhdis_dX = (-_R * temperature * temperature * dKh_dT / Kh / _Mco2 - hdis) / dX;
1520 }
const Real _Mco2
Molar mass of CO2 (kg/mol)
const std::string temperature
Definition: NS.h:27
const Real _R
Universal gas constant (J/mol/K)
void henryConstant(Real temperature, Real Xnacl, Real &Kh, Real &dKh_dT, Real &dKh_dx) const
Henry&#39;s constant of dissolution of gas phase CO2 in brine.

◆ equilibriumConstantCO2()

void PorousFlowBrineCO2::equilibriumConstantCO2 ( Real  temperature,
Real &  kco2,
Real &  dkco2_dT 
) const

Equilibrium constant for CO2 For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2003) For temperatures 109C <= T <= 300C, uses Spycher and Pruess (2010) For temperatures 99C < T < 109C, the value is calculated by smoothly interpolating the two formulations.

Parameters
temperaturetemperature (K)
[out]kco2equilibrium constant for CO2
[out]dkco2_dTderivative of equilibrium constant wrt temperature

Definition at line 950 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp(), and funcABLowTemp().

951 {
952  // Uses temperature in Celcius
953  const Real Tc = temperature - _T_c2k;
954  const Real Tc2 = Tc * Tc;
955  const Real Tc3 = Tc2 * Tc;
956 
957  Real logK0CO2, dlogK0CO2;
958 
959  if (Tc <= 99.0)
960  {
961  logK0CO2 = 1.189 + 1.304e-2 * Tc - 5.446e-5 * Tc2;
962  dlogK0CO2 = 1.304e-2 - 1.0892e-4 * Tc;
963  }
964  else if (Tc > 99.0 && Tc < 109.0)
965  {
966  const Real Tint = (Tc - 99.0) / 10.0;
967  const Real Tint2 = Tint * Tint;
968  logK0CO2 = 1.9462 + 2.25692e-2 * Tint - 9.49577e-3 * Tint2 - 6.77721e-3 * Tint * Tint2;
969  dlogK0CO2 = 2.25692e-2 - 1.899154e-2 * Tint + 2.033163e-2 * Tint2;
970  }
971  else // 109 <= Tc <= 300
972  {
973  logK0CO2 = 1.668 + 3.992e-3 * Tc - 1.156e-5 * Tc2 + 1.593e-9 * Tc3;
974  dlogK0CO2 = 3.992e-3 - 2.312e-5 * Tc + 4.779e-9 * Tc2;
975  }
976 
977  kco2 = std::pow(10.0, logK0CO2);
978  dkco2_dT = std::log(10.0) * dlogK0CO2 * kco2;
979 }
const Real _T_c2k
Conversion from C to K.
const std::string temperature
Definition: NS.h:27
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)

◆ equilibriumConstantH2O()

void PorousFlowBrineCO2::equilibriumConstantH2O ( Real  temperature,
Real &  kh2o,
Real &  dkh2o_dT 
) const

Equilibrium constant for H2O For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2003) For temperatures 109C <= T <= 300C, uses Spycher and Pruess (2010) For temperatures 99C < T < 109C, the value is calculated by smoothly interpolating the two formulations.

Parameters
temperaturetemperature (K)
[out]kh2oequilibrium constant for H2O
[out]dkh2o_dTderivative of equilibrium constant wrt temperature

Definition at line 917 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp(), and funcABLowTemp().

918 {
919  // Uses temperature in Celcius
920  const Real Tc = temperature - _T_c2k;
921  const Real Tc2 = Tc * Tc;
922  const Real Tc3 = Tc2 * Tc;
923  const Real Tc4 = Tc3 * Tc;
924 
925  Real logK0H2O, dlogK0H2O;
926 
927  if (Tc <= 99.0)
928  {
929  logK0H2O = -2.209 + 3.097e-2 * Tc - 1.098e-4 * Tc2 + 2.048e-7 * Tc3;
930  dlogK0H2O = 3.097e-2 - 2.196e-4 * Tc + 6.144e-7 * Tc2;
931  }
932  else if (Tc > 99.0 && Tc < 109.0)
933  {
934  const Real Tint = (Tc - 99.0) / 10.0;
935  const Real Tint2 = Tint * Tint;
936  logK0H2O = -0.0204026 + 0.0152513 * Tint + 0.417565 * Tint2 - 0.278636 * Tint * Tint2;
937  dlogK0H2O = 0.0152513 + 0.83513 * Tint - 0.835908 * Tint2;
938  }
939  else // 109 <= Tc <= 300
940  {
941  logK0H2O = -2.1077 + 2.8127e-2 * Tc - 8.4298e-5 * Tc2 + 1.4969e-7 * Tc3 - 1.1812e-10 * Tc4;
942  dlogK0H2O = 2.8127e-2 - 1.68596e-4 * Tc + 4.4907e-7 * Tc2 - 4.7248e-10 * Tc3;
943  }
944 
945  kh2o = std::pow(10.0, logK0H2O);
946  dkh2o_dT = std::log(10.0) * dlogK0H2O * kh2o;
947 }
const Real _T_c2k
Conversion from C to K.
const std::string temperature
Definition: NS.h:27
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)

◆ equilibriumMassFractions()

void PorousFlowBrineCO2::equilibriumMassFractions ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real &  Xco2,
Real &  dXco2_dp,
Real &  dXco2_dT,
Real &  dXco2_dX,
Real &  Yh2o,
Real &  dYh2o_dp,
Real &  dYh2o_dT,
Real &  dYh2o_dX 
) const

Mass fractions of CO2 in brine and water vapor in CO2 at equilibrium.

Parameters
pressurephase pressure (Pa)
temperaturephase temperature (K)
XnaclNaCl mass fraction (kg/kg)
[out]Xco2mass fraction of CO2 in liquid (kg/kg)
[out]dXco2_dpderivative of mass fraction of CO2 in liquid wrt pressure
[out]dXco2_dTderivative of mass fraction of CO2 in liqiud wrt temperature
[out]dXco2_dXderivative of mass fraction of CO2 in liqiud wrt salt mass fraction
[out]Yh2omass fraction of H2O in gas (kg/kg)
[out]dYh2og_dpderivative of mass fraction of H2O in gas wrt pressure
[out]dYh2o_dTderivative of mass fraction of H2O in gas wrt temperature
[out]dYh2o_dXderivative of mass fraction of H2O in gas wrt salt mass fraction

Definition at line 536 of file PorousFlowBrineCO2.C.

Referenced by massFractions(), and totalMassFraction().

547 {
548  Real co2_density, dco2_density_dp, dco2_density_dT;
549  _co2_fp.rho_from_p_T(pressure, temperature, co2_density, dco2_density_dp, dco2_density_dT);
550 
551  // Mole fractions at equilibrium
552  Real xco2, dxco2_dp, dxco2_dT, dxco2_dX, yh2o, dyh2o_dp, dyh2o_dT, dyh2o_dX;
554  temperature,
555  Xnacl,
556  xco2,
557  dxco2_dp,
558  dxco2_dT,
559  dxco2_dX,
560  yh2o,
561  dyh2o_dp,
562  dyh2o_dT,
563  dyh2o_dX);
564 
565  // The mass fraction of H2O in gas (assume no salt in gas phase) and derivatives
566  // wrt p, T, and X
567  Yh2o = yh2o * _Mh2o / (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2);
568  dYh2o_dp = _Mco2 * _Mh2o * dyh2o_dp / (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2) /
569  (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2);
570  dYh2o_dT = _Mco2 * _Mh2o * dyh2o_dT / (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2) /
571  (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2);
572  dYh2o_dX = dyh2o_dX * _Mh2o * _Mco2 / (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2) /
573  (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2);
574 
575  // NaCl molality (mol/kg)
576  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
577  const Real dmnacl_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
578 
579  // The molality of CO2 in 1kg of H2O
580  const Real mco2 = xco2 * (2.0 * mnacl + _invMh2o) / (1.0 - xco2);
581  // The mass fraction of CO2 in brine is then
582  const Real denominator = (1.0 + mnacl * _Mnacl + mco2 * _Mco2);
583  Xco2 = mco2 * _Mco2 / denominator;
584 
585  // Derivatives of Xco2 wrt p, T and X
586  const Real denominator2 = denominator * denominator;
587  const Real dmco2_dp = dxco2_dp * (2.0 * mnacl + _invMh2o) / (1.0 - xco2) / (1.0 - xco2);
588  dXco2_dp = (1.0 + mnacl * _Mnacl) * _Mco2 * dmco2_dp / denominator2;
589 
590  const Real dmco2_dT = dxco2_dT * (2.0 * mnacl + _invMh2o) / (1.0 - xco2) / (1.0 - xco2);
591  dXco2_dT = (1.0 + mnacl * _Mnacl) * _Mco2 * dmco2_dT / denominator2;
592 
593  const Real dmco2_dX =
594  (dxco2_dX * (2.0 * mnacl + _invMh2o) + 2.0 * (1.0 - xco2) * xco2 * dmnacl_dX) / (1.0 - xco2) /
595  (1.0 - xco2);
596  dXco2_dX = _Mco2 * ((1.0 + mnacl * _Mnacl) * dmco2_dX - dmnacl_dX * mco2 * _Mnacl) / denominator2;
597 }
const Real _Mh2o
Molar mass of water (kg/mol)
const Real _Mnacl
Molar mass of NaCL.
const Real _invMh2o
Inverse of molar mass of H2O (mol/kg)
const Real _Mco2
Molar mass of CO2 (kg/mol)
const std::string temperature
Definition: NS.h:27
const SinglePhaseFluidProperties & _co2_fp
Fluid properties UserObject for the CO2.
const std::string pressure
Definition: NS.h:26
void equilibriumMoleFractions(Real pressure, Real temperature, Real Xnacl, Real &xco2, Real &dxco2_dp, Real &dxco2_dT, Real &dxco2_dX, Real &yh2o, Real &dyh2o_dp, Real &dyh2o_dT, Real &dyh2o_dX) const
Mole fractions of CO2 in brine and water vapor in CO2 at equilibrium.

◆ equilibriumMoleFractions()

void PorousFlowBrineCO2::equilibriumMoleFractions ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real &  xco2,
Real &  dxco2_dp,
Real &  dxco2_dT,
Real &  dxco2_dX,
Real &  yh2o,
Real &  dyh2o_dp,
Real &  dyh2o_dT,
Real &  dyh2o_dX 
) const

Mole fractions of CO2 in brine and water vapor in CO2 at equilibrium.

In the low temperature regime (T <= 99C), the mole fractions are calculated directly, while in the elevated temperature regime (T >= 109C), they are calculated iteratively. In the intermediate regime (99C < T < 109C), a cubic polynomial is used to smoothly connect the low and elevated temperature regimes.

Parameters
pressurephase pressure (Pa)
temperaturephase temperature (K)
XnaclNaCl mass fraction (kg/kg)
[out]xco2mole fraction of CO2 in liquid
[out]dxco2_dpderivative of mole fraction of CO2 in liquid wrt pressure
[out]dxco2_dTderivative of mole fraction of CO2 in liqiud wrt temperature
[out]dxco2_dXderivative of mole fraction of CO2 in liqiud wrt salt mass fraction
[out]yh2omass fraction of mole in gas
[out]dyh2og_dpderivative of mole fraction of H2O in gas wrt pressure
[out]dyh2o_dTderivative of mole fraction of H2O in gas wrt temperature
[out]dyh2o_dXderivative of mole fraction of H2O in gas wrt salt mass fraction

Definition at line 982 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMassFractions().

993 {
994  // CO2 density and derivatives wrt pressure and temperature
995  Real co2_density, dco2_density_dp, dco2_density_dT;
996  _co2_fp.rho_from_p_T(pressure, temperature, co2_density, dco2_density_dp, dco2_density_dT);
997 
998  if (temperature <= _Tlower)
999  {
1001  temperature,
1002  Xnacl,
1003  xco2,
1004  dxco2_dp,
1005  dxco2_dT,
1006  dxco2_dX,
1007  yh2o,
1008  dyh2o_dp,
1009  dyh2o_dT,
1010  dyh2o_dX);
1011  }
1012  else if (temperature > _Tlower && temperature < _Tupper)
1013  {
1014  // Cubic polynomial in this regime
1015  const Real Tint = (temperature - _Tlower) / 10.0;
1016 
1017  // Equilibrium mole fractions and derivatives at the lower temperature
1018  Real xco2_lower, dxco2_dT_lower, yh2o_lower, dyh2o_dT_lower;
1020  _Tlower,
1021  Xnacl,
1022  xco2_lower,
1023  dxco2_dp,
1024  dxco2_dT_lower,
1025  dxco2_dX,
1026  yh2o_lower,
1027  dyh2o_dp,
1028  dyh2o_dT_lower,
1029  dyh2o_dX);
1030 
1031  // Equilibrium mole fractions and derivatives at the upper temperature
1032  Real xco2_upper, yh2o_upper;
1033  Real co2_density_upper = _co2_fp.rho_from_p_T(pressure, _Tupper);
1034 
1036  pressure, _Tupper, Xnacl, co2_density_upper, xco2_upper, yh2o_upper);
1037 
1038  Real A, dA_dp, dA_dT, B, dB_dp, dB_dT, dB_dX;
1040  _Tupper,
1041  Xnacl,
1042  co2_density,
1043  xco2_upper,
1044  yh2o_upper,
1045  A,
1046  dA_dp,
1047  dA_dT,
1048  B,
1049  dB_dp,
1050  dB_dT,
1051  dB_dX);
1052 
1053  const Real dyh2o_dT_upper =
1054  ((1.0 - B) * dA_dT + (A - 1.0) * A * dB_dT) / (1.0 - A * B) / (1.0 - A * B);
1055  const Real dxco2_dT_upper = dB_dT * (1.0 - yh2o_upper) - B * dyh2o_dT_upper;
1056 
1057  // The mole fractions in this regime are then found by interpolation
1058  Real dxco2_dT, dyh2o_dT;
1060  Tint, xco2_lower, dxco2_dT_lower, xco2_upper, dxco2_dT_upper, xco2, dxco2_dT);
1062  Tint, yh2o_lower, dyh2o_dT_lower, yh2o_upper, dyh2o_dT_upper, yh2o, dyh2o_dT);
1063  }
1064  else
1065  {
1066  // Equilibrium mole fractions solved using iteration in this regime
1067  solveEquilibriumMoleFractionHighTemp(pressure, temperature, Xnacl, co2_density, xco2, yh2o);
1068 
1069  // Can use these in funcABHighTemp() to get derivatives analyticall rather than by iteration
1070  Real A, dA_dp, dA_dT, B, dB_dp, dB_dT, dB_dX;
1072  temperature,
1073  Xnacl,
1074  co2_density,
1075  xco2,
1076  yh2o,
1077  A,
1078  dA_dp,
1079  dA_dT,
1080  B,
1081  dB_dp,
1082  dB_dT,
1083  dB_dX);
1084 
1085  dyh2o_dp = ((1.0 - B) * dA_dp + (A - 1.0) * A * dB_dp) / (1.0 - A * B) / (1.0 - A * B);
1086  dxco2_dp = dB_dp * (1.0 - yh2o) - B * dyh2o_dp;
1087 
1088  dyh2o_dT = ((1.0 - B) * dA_dT + (A - 1.0) * A * dB_dT) / (1.0 - A * B) / (1.0 - A * B);
1089  dxco2_dT = dB_dT * (1.0 - yh2o) - B * dyh2o_dT;
1090 
1091  dyh2o_dX = ((A - 1.0) * A * dB_dX) / (1.0 - A * B) / (1.0 - A * B);
1092  dxco2_dX = dB_dX * (1.0 - yh2o) - B * dyh2o_dX;
1093  }
1094 }
const Real _Tlower
Temperature below which the Spycher, Pruess & Ennis-King (2003) model is used (K) ...
void equilibriumMoleFractionsLowTemp(Real pressure, Real temperature, Real Xnacl, Real &xco2, Real &dxco2_dp, Real &dxco2_dT, Real &dxco2_dX, Real &yh2o, Real &dyh2o_dp, Real &dyh2o_dT, Real &dyh2o_dX) const
Mole fractions of CO2 in brine and water vapor in CO2 at equilibrium in the low temperature regime (T...
const std::string temperature
Definition: NS.h:27
void solveEquilibriumMoleFractionHighTemp(Real pressure, Real temperature, Real Xnacl, Real co2_density, Real &xco2, Real &yh2o) const
Function to solve for yh2o and xco2 iteratively in the elevated temperature regime (T > 100C) ...
const SinglePhaseFluidProperties & _co2_fp
Fluid properties UserObject for the CO2.
void funcABHighTemp(Real pressure, Real temperature, Real Xnacl, Real co2_density, Real xco2, Real yh2o, Real &A, Real &B) const
The function A (Eq.
const std::string pressure
Definition: NS.h:26
void smoothCubicInterpolation(Real temperature, Real f0, Real df0, Real f1, Real df1, Real &value, Real &deriv) const
Cubic function to smoothly interpolate between the low temperature and elevated temperature models fo...
const Real _Tupper
Temperature above which the Spycher & Pruess (2010) model is used (K)

◆ equilibriumMoleFractionsLowTemp()

void PorousFlowBrineCO2::equilibriumMoleFractionsLowTemp ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real &  xco2,
Real &  dxco2_dp,
Real &  dxco2_dT,
Real &  dxco2_dX,
Real &  yh2o,
Real &  dyh2o_dp,
Real &  dyh2o_dT,
Real &  dyh2o_dX 
) const

Mole fractions of CO2 in brine and water vapor in CO2 at equilibrium in the low temperature regime (T <= 99C).

Parameters
pressurephase pressure (Pa)
temperaturephase temperature (K)
XnaclNaCl mass fraction (kg/kg)
[out]xco2mole fraction of CO2 in liquid
[out]dxco2_dpderivative of mole fraction of CO2 in liquid wrt pressure
[out]dxco2_dTderivative of mole fraction of CO2 in liqiud wrt temperature
[out]dxco2_dXderivative of mole fraction of CO2 in liqiud wrt salt mass fraction
[out]yh2omass fraction of mole in gas
[out]dyh2og_dpderivative of mole fraction of H2O in gas wrt pressure
[out]dyh2o_dTderivative of mole fraction of H2O in gas wrt temperature
[out]dyh2o_dXderivative of mole fraction of H2O in gas wrt salt mass fraction

Definition at line 1097 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractions().

1108 {
1109  if (temperature > 373.15)
1110  mooseError(name(),
1111  ": equilibriumMoleFractionsLowTemp() is not valid for T > 373.15K. Use "
1112  "equilibriumMoleFractions() instead");
1113 
1114  // CO2 density and derivatives wrt pressure and temperature
1115  Real co2_density, dco2_density_dp, dco2_density_dT;
1116  _co2_fp.rho_from_p_T(pressure, temperature, co2_density, dco2_density_dp, dco2_density_dT);
1117 
1118  // Assume infinite dilution (yh20 = 0 and xco2 = 0) in low temperature regime
1119  Real A, dA_dp, dA_dT, B, dB_dp, dB_dT;
1121  temperature,
1122  co2_density,
1123  dco2_density_dp,
1124  dco2_density_dT,
1125  A,
1126  dA_dp,
1127  dA_dT,
1128  B,
1129  dB_dp,
1130  dB_dT);
1131 
1132  // As the activity coefficient for CO2 in brine used in this regime isn't a 'true'
1133  // activity coefficient, we instead calculate the molality of CO2 in water, then
1134  // correct it for brine, and then calculate the mole fractions.
1135  // The mole fraction in pure water is
1136  const Real yh2ow = (1.0 - B) / (1.0 / A - B);
1137  const Real xco2w = B * (1.0 - yh2ow);
1138 
1139  // Molality of CO2 in pure water
1140  const Real mco2w = xco2w * _invMh2o / (1.0 - xco2w);
1141  // Molality of CO2 in brine is then calculated using gamma
1142  Real gamma, dgamma_dp, dgamma_dT, dgamma_dX;
1143  activityCoefficient(pressure, temperature, Xnacl, gamma, dgamma_dp, dgamma_dT, dgamma_dX);
1144  const Real mco2 = mco2w / gamma;
1145 
1146  // Need NaCl molality (mol/kg)
1147  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
1148  const Real dmnacl_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
1149 
1150  // Mole fractions of CO2 and H2O in liquid and gas phases
1151  const Real total_moles = 2.0 * mnacl + _invMh2o + mco2;
1152  xco2 = mco2 / total_moles;
1153  yh2o = A * (1.0 - xco2 - 2.0 * mnacl / total_moles);
1154 
1155  // The derivatives of the mole fractions wrt pressure
1156  const Real dyh2ow_dp =
1157  ((1.0 - B) * dA_dp + (A - 1.0) * A * dB_dp) / (1.0 - A * B) / (1.0 - A * B);
1158  const Real dxco2w_dp = dB_dp * (1.0 - yh2ow) - B * dyh2ow_dp;
1159 
1160  const Real dmco2w_dp = _invMh2o * dxco2w_dp / (1.0 - xco2w) / (1.0 - xco2w);
1161  const Real dmco2_dp = dmco2w_dp / gamma - mco2w * dgamma_dp / gamma / gamma;
1162  dxco2_dp = (2.0 * mnacl + _invMh2o) * dmco2_dp / total_moles / total_moles;
1163 
1164  dyh2o_dp = (1.0 - xco2 - 2.0 * mnacl / total_moles) * dA_dp - A * dxco2_dp +
1165  2.0 * A * mnacl * dmco2_dp / total_moles / total_moles;
1166 
1167  // The derivatives of the mole fractions wrt temperature
1168  const Real dyh2ow_dT =
1169  ((1.0 - B) * dA_dT + (A - 1.0) * A * dB_dT) / (1.0 - A * B) / (1.0 - A * B);
1170  const Real dxco2w_dT = dB_dT * (1.0 - yh2ow) - B * dyh2ow_dT;
1171 
1172  const Real dmco2w_dT = _invMh2o * dxco2w_dT / (1.0 - xco2w) / (1.0 - xco2w);
1173  const Real dmco2_dT = dmco2w_dT / gamma - mco2w * dgamma_dT / gamma / gamma;
1174  dxco2_dT = (2.0 * mnacl + _invMh2o) * dmco2_dT / total_moles / total_moles;
1175 
1176  dyh2o_dT = (1.0 - xco2 - 2.0 * mnacl / total_moles) * dA_dT - A * dxco2_dT +
1177  2.0 * A * mnacl * dmco2_dT / total_moles / total_moles;
1178 
1179  // The derivatives of the mole fractions wrt salt mass fraction
1180  const Real dmco2_dX = -mco2w * dgamma_dX / gamma / gamma;
1181  dxco2_dX =
1182  ((2.0 * mnacl + _invMh2o) * dmco2_dX - 2.0 * mco2 * dmnacl_dX) / total_moles / total_moles;
1183  dyh2o_dX =
1184  A * (2.0 * (mnacl * dmco2_dX - (mco2 + _invMh2o) * dmnacl_dX) / total_moles / total_moles -
1185  dxco2_dX);
1186 }
const Real _Mnacl
Molar mass of NaCL.
const Real _invMh2o
Inverse of molar mass of H2O (mol/kg)
const std::string temperature
Definition: NS.h:27
void activityCoefficient(Real pressure, Real temperature, Real Xnacl, Real &gamma, Real &dgamma_dp, Real &dgamma_dT, Real &dgamma_dX) const
Activity coefficient for CO2 in brine.
const SinglePhaseFluidProperties & _co2_fp
Fluid properties UserObject for the CO2.
const std::string name
Definition: Setup.h:22
const std::string pressure
Definition: NS.h:26
void funcABLowTemp(Real pressure, Real temperature, Real co2_density, Real dco2_density_dp, Real dco2_density_dT, Real &A, Real &dA_dp, Real &dA_dT, Real &B, Real &dB_dp, Real &dB_dT) const

◆ execute()

void PorousFlowFluidStateFlash::execute ( )
inlinefinalinherited

Definition at line 37 of file PorousFlowFluidStateFlash.h.

37 {};

◆ finalize()

void PorousFlowFluidStateFlash::finalize ( )
inlinefinalinherited

Definition at line 38 of file PorousFlowFluidStateFlash.h.

38 {};

◆ fluidStateName()

std::string PorousFlowBrineCO2::fluidStateName ( ) const
overridevirtual

Name of FluidState.

Implements PorousFlowFluidStateBase.

Definition at line 84 of file PorousFlowBrineCO2.C.

85 {
86  return "brine-co2";
87 }

◆ fugacityCoefficientCO2HighTemp()

Real PorousFlowBrineCO2::fugacityCoefficientCO2HighTemp ( Real  pressure,
Real  temperature,
Real  co2_density,
Real  xco2,
Real  yh2o 
) const

Definition at line 774 of file PorousFlowBrineCO2.C.

Referenced by fugacityCoefficientsHighTemp().

776 {
777  // Need pressure in bar
778  const Real pbar = pressure * 1.0e-5;
779  // Molar volume in cm^3/mol
780  const Real V = _Mco2 / co2_density * 1.0e6;
781 
782  // Redlich-Kwong parameters
783  const Real yco2 = 1.0 - yh2o;
784  const Real xh2o = 1.0 - xco2;
785 
786  const Real aCO2 = 8.008e7 - 4.984e4 * temperature;
787  const Real aH2O = 1.337e8 - 1.4e4 * temperature;
788  const Real bCO2 = 28.25;
789  const Real bH2O = 15.7;
790  const Real KH2OCO2 = 1.427e-2 - 4.037e-4 * temperature;
791  const Real KCO2H2O = 0.4228 - 7.422e-4 * temperature;
792  const Real kH2OCO2 = KH2OCO2 * yh2o + KCO2H2O * yco2;
793  const Real kCO2H2O = KCO2H2O * yh2o + KH2OCO2 * yco2;
794 
795  const Real aH2OCO2 = std::sqrt(aCO2 * aH2O) * (1.0 - kH2OCO2);
796  const Real aCO2H2O = std::sqrt(aCO2 * aH2O) * (1.0 - kCO2H2O);
797 
798  const Real amix = yh2o * yh2o * aH2O + yh2o * yco2 * (aH2OCO2 + aCO2H2O) + yco2 * yco2 * aCO2;
799  const Real bmix = yh2o * bH2O + yco2 * bCO2;
800 
801  const Real t15 = std::pow(temperature, 1.5);
802 
803  Real lnPhiCO2 = bCO2 / bmix * (pbar * V / (_Rbar * temperature) - 1.0) -
804  std::log(pbar * (V - bmix) / (_Rbar * temperature));
805 
806  Real term3 = (2.0 * yco2 * aCO2 + yh2o * (aH2OCO2 + aCO2H2O) -
807  yh2o * yco2 * std::sqrt(aH2O * aCO2) * (kH2OCO2 - kCO2H2O) * (yh2o - yco2) +
808  xh2o * xco2 * std::sqrt(aH2O * aCO2) * (kCO2H2O - kH2OCO2)) /
809  amix;
810 
811  lnPhiCO2 += (term3 - bCO2 / bmix) * amix / (bmix * _Rbar * t15) * std::log(V / (V + bmix));
812 
813  return std::exp(lnPhiCO2);
814 }
const Real _Mco2
Molar mass of CO2 (kg/mol)
const std::string temperature
Definition: NS.h:27
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
const std::string pressure
Definition: NS.h:26

◆ fugacityCoefficientH2OHighTemp()

Real PorousFlowBrineCO2::fugacityCoefficientH2OHighTemp ( Real  pressure,
Real  temperature,
Real  co2_density,
Real  xco2,
Real  yh2o 
) const

Definition at line 731 of file PorousFlowBrineCO2.C.

Referenced by fugacityCoefficientsHighTemp().

733 {
734  // Need pressure in bar
735  const Real pbar = pressure * 1.0e-5;
736  // Molar volume in cm^3/mol
737  const Real V = _Mco2 / co2_density * 1.0e6;
738 
739  // Redlich-Kwong parameters
740  const Real yco2 = 1.0 - yh2o;
741  const Real xh2o = 1.0 - xco2;
742 
743  const Real aCO2 = 8.008e7 - 4.984e4 * temperature;
744  const Real aH2O = 1.337e8 - 1.4e4 * temperature;
745  const Real bCO2 = 28.25;
746  const Real bH2O = 15.7;
747  const Real KH2OCO2 = 1.427e-2 - 4.037e-4 * temperature;
748  const Real KCO2H2O = 0.4228 - 7.422e-4 * temperature;
749  const Real kH2OCO2 = KH2OCO2 * yh2o + KCO2H2O * yco2;
750  const Real kCO2H2O = KCO2H2O * yh2o + KH2OCO2 * yco2;
751 
752  const Real aH2OCO2 = std::sqrt(aCO2 * aH2O) * (1.0 - kH2OCO2);
753  const Real aCO2H2O = std::sqrt(aCO2 * aH2O) * (1.0 - kCO2H2O);
754 
755  const Real amix = yh2o * yh2o * aH2O + yh2o * yco2 * (aH2OCO2 + aCO2H2O) + yco2 * yco2 * aCO2;
756  const Real bmix = yh2o * bH2O + yco2 * bCO2;
757 
758  const Real t15 = std::pow(temperature, 1.5);
759 
760  Real lnPhiH2O = bH2O / bmix * (pbar * V / (_Rbar * temperature) - 1.0) -
761  std::log(pbar * (V - bmix) / (_Rbar * temperature));
762  Real term3 = (2.0 * yh2o * aH2O + yco2 * (aH2OCO2 + aCO2H2O) -
763  yh2o * yco2 * std::sqrt(aH2O * aCO2) * (kH2OCO2 - kCO2H2O) * (yh2o - yco2) +
764  xh2o * xco2 * std::sqrt(aH2O * aCO2) * (kH2OCO2 - kCO2H2O)) /
765  amix;
766  term3 -= bH2O / bmix;
767  term3 *= amix / (bmix * _Rbar * t15) * std::log(V / (V + bmix));
768  lnPhiH2O += term3;
769 
770  return std::exp(lnPhiH2O);
771 }
const Real _Mco2
Molar mass of CO2 (kg/mol)
const std::string temperature
Definition: NS.h:27
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
const std::string pressure
Definition: NS.h:26

◆ fugacityCoefficientsHighTemp() [1/2]

void PorousFlowBrineCO2::fugacityCoefficientsHighTemp ( Real  pressure,
Real  temperature,
Real  co2_density,
Real  xco2,
Real  yh2o,
Real &  fco2,
Real &  fh2o 
) const

Fugacity coefficients for H2O and CO2 at elevated temperatures (100C < T <= 300C).

Eq. (A8) from Spycher & Pruess (2010)

Parameters
pressuregas pressure (Pa)
temperaturetemperature (K)
co2_densityCO2 density (kg/m^3)
xco2mole fraction of CO2 in liquid phase (-)
yh2omole fraction of H2O in gas phase (-)
[out]fco2fugacity coefficient for CO2
[out]dfco2_dpderivative of fugacity coefficient wrt pressure
[out]dfco2_dTderivative of fugacity coefficient wrt temperature
[out]fh2ofugacity coefficient for H2O
[out]dfh2o_dpderivative of fugacity coefficient wrt pressure
[out]dfh2o_dTderivative of fugacity coefficient wrt temperature

Definition at line 675 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp().

682 {
683  if (temperature <= 373.15)
684  mooseError(name(),
685  ": fugacityCoefficientsHighTemp() is not valid for T <= 373.15K. Use "
686  "fugacityCoefficientsLowTemp() instead");
687 
688  fh2o = fugacityCoefficientH2OHighTemp(pressure, temperature, co2_density, xco2, yh2o);
689  fco2 = fugacityCoefficientCO2HighTemp(pressure, temperature, co2_density, xco2, yh2o);
690 }
Real fugacityCoefficientH2OHighTemp(Real pressure, Real temperature, Real co2_density, Real xco2, Real yh2o) const
Real fugacityCoefficientCO2HighTemp(Real pressure, Real temperature, Real co2_density, Real xco2, Real yh2o) const
const std::string temperature
Definition: NS.h:27
const std::string name
Definition: Setup.h:22
const std::string pressure
Definition: NS.h:26

◆ fugacityCoefficientsHighTemp() [2/2]

void PorousFlowBrineCO2::fugacityCoefficientsHighTemp ( Real  pressure,
Real  temperature,
Real  co2_density,
Real  xco2,
Real  yh2o,
Real &  fco2,
Real &  dfco2_dp,
Real &  dfco2_dT,
Real &  fh2o,
Real &  dfh2o_dp,
Real &  dfh2o_dT 
) const

Definition at line 693 of file PorousFlowBrineCO2.C.

704 {
705  if (temperature <= 373.15)
706  mooseError(name(),
707  ": fugacityCoefficientsHighTemp() is not valid for T <= 373.15K. Use "
708  "fugacityCoefficientsLowTemp() instead");
709 
710  fh2o = fugacityCoefficientH2OHighTemp(pressure, temperature, co2_density, xco2, yh2o);
711  fco2 = fugacityCoefficientCO2HighTemp(pressure, temperature, co2_density, xco2, yh2o);
712 
713  // Derivatives by finite difference
714  const Real dp = 1.0e-2;
715  const Real dT = 1.0e-4;
716 
717  Real fh2o2 = fugacityCoefficientH2OHighTemp(pressure + dp, temperature, co2_density, xco2, yh2o);
718  Real fco22 = fugacityCoefficientCO2HighTemp(pressure + dp, temperature, xco2, co2_density, yh2o);
719 
720  dfh2o_dp = (fh2o2 - fh2o) / dp;
721  dfco2_dp = (fco22 - fco2) / dp;
722 
723  fh2o2 = fugacityCoefficientH2OHighTemp(pressure, temperature + dT, co2_density, xco2, yh2o);
724  fco22 = fugacityCoefficientCO2HighTemp(pressure, temperature + dT, co2_density, xco2, yh2o);
725 
726  dfh2o_dT = (fh2o2 - fh2o) / dT;
727  dfco2_dT = (fco22 - fco2) / dT;
728 }
Real fugacityCoefficientH2OHighTemp(Real pressure, Real temperature, Real co2_density, Real xco2, Real yh2o) const
Real fugacityCoefficientCO2HighTemp(Real pressure, Real temperature, Real co2_density, Real xco2, Real yh2o) const
const std::string temperature
Definition: NS.h:27
const std::string name
Definition: Setup.h:22
const std::string pressure
Definition: NS.h:26

◆ fugacityCoefficientsLowTemp()

void PorousFlowBrineCO2::fugacityCoefficientsLowTemp ( Real  pressure,
Real  temperature,
Real  co2_density,
Real  dco2_density_dp,
Real  dco2_density_dT,
Real &  fco2,
Real &  dfco2_dp,
Real &  dfco2_dT,
Real &  fh2o,
Real &  dfh2o_dp,
Real &  dfh2o_dT 
) const

Fugacity coefficients for H2O and CO2 for T <= 100C Eq.

(B7) from Spycher, Pruess and Ennis-King (2003)

Parameters
pressuregas pressure (Pa)
temperaturetemperature (K)
co2_densityCO2 density (kg/m^3)
dco2_density_dpderivative of CO2 density wrt pressure
dco2_density_dTderivative of CO2 density wrt temperature
[out]fco2fugacity coefficient for CO2
[out]dfco2_dpderivative of fugacity coefficient wrt pressure
[out]dfco2_dTderivative of fugacity coefficient wrt temperature
[out]fh2ofugacity coefficient for H2O
[out]dfh2o_dpderivative of fugacity coefficient wrt pressure
[out]dfh2o_dTderivative of fugacity coefficient wrt temperature

Definition at line 600 of file PorousFlowBrineCO2.C.

Referenced by funcABLowTemp().

611 {
612  if (temperature > 373.15)
613  mooseError(name(),
614  ": fugacityCoefficientsLowTemp() is not valid for T > 373.15K. Use "
615  "fugacityCoefficientsHighTemp() instead");
616 
617  // Need pressure in bar
618  const Real pbar = pressure * 1.0e-5;
619 
620  // Molar volume in cm^3/mol
621  const Real V = _Mco2 / co2_density * 1.0e6;
622  const Real dV_dp = -V / co2_density * dco2_density_dp;
623  const Real dV_dT = -V / co2_density * dco2_density_dT;
624 
625  // Redlich-Kwong parameters
626  const Real aCO2 = 7.54e7 - 4.13e4 * temperature;
627  const Real daCO2_dT = -4.13e4;
628  const Real bCO2 = 27.8;
629  const Real aCO2H2O = 7.89e7;
630  const Real bH2O = 18.18;
631 
632  const Real t15 = std::pow(temperature, 1.5);
633 
634  // The fugacity coefficients for H2O and CO2
635  auto lnPhi = [V, aCO2, bCO2, t15, this](Real a, Real b) {
636  return std::log(V / (V - bCO2)) + b / (V - bCO2) -
637  2.0 * a / (_Rbar * t15 * bCO2) * std::log((V + bCO2) / V) +
638  aCO2 * b / (_Rbar * t15 * bCO2 * bCO2) * (std::log((V + bCO2) / V) - bCO2 / (V + bCO2));
639  };
640 
641  const Real lnPhiH2O = lnPhi(aCO2H2O, bH2O) - std::log(pbar * V / (_Rbar * temperature));
642  const Real lnPhiCO2 = lnPhi(aCO2, bCO2) - std::log(pbar * V / (_Rbar * temperature));
643 
644  fh2o = std::exp(lnPhiH2O);
645  fco2 = std::exp(lnPhiCO2);
646 
647  // The derivative of the fugacity coefficients wrt pressure
648  auto dlnPhi_dV = [V, aCO2, bCO2, t15, this](Real a, Real b) {
649  return (bCO2 * bCO2 - (bCO2 + b) * V) / V / (V - bCO2) / (V - bCO2) +
650  (2.0 * a * (V + bCO2) - aCO2 * b) / (_Rbar * t15 * V * (V + bCO2) * (V + bCO2));
651  };
652 
653  dfh2o_dp = (dlnPhi_dV(aCO2H2O, bH2O) * dV_dp - 1.0e-5 / pbar - dV_dp / V) * fh2o;
654  dfco2_dp = (dlnPhi_dV(aCO2, bCO2) * dV_dp - 1.0e-5 / pbar - dV_dp / V) * fco2;
655 
656  // The derivative of the fugacity coefficient wrt temperature
657  auto dlnPhi_dT = [V, aCO2, daCO2_dT, bCO2, t15, temperature, this](Real a, Real b, Real da_dT) {
658  return (3.0 * bCO2 * a * std::log((V + bCO2) / V) +
659  1.5 * aCO2 * b * (bCO2 / (V + bCO2) - std::log((V + bCO2) / V)) -
660  (2.0 * da_dT * bCO2 * std::log((V + bCO2) / V) -
661  daCO2_dT * b * (std::log((V + bCO2) / V) - bCO2 / (V + bCO2))) *
662  temperature) /
663  (temperature * t15 * _Rbar * bCO2 * bCO2);
664  };
665 
666  dfh2o_dT = (dlnPhi_dT(aCO2H2O, bH2O, 0) + dlnPhi_dV(aCO2H2O, bH2O) * dV_dT - dV_dT / V +
667  1.0 / temperature) *
668  fh2o;
669  dfco2_dT = (dlnPhi_dT(aCO2, bCO2, daCO2_dT) + dlnPhi_dV(aCO2, bCO2) * dV_dT - dV_dT / V +
670  1.0 / temperature) *
671  fco2;
672 }
const Real _Mco2
Molar mass of CO2 (kg/mol)
const std::string temperature
Definition: NS.h:27
const std::string name
Definition: Setup.h:22
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
const std::string pressure
Definition: NS.h:26

◆ funcABHighTemp() [1/2]

void PorousFlowBrineCO2::funcABHighTemp ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real  co2_density,
Real  xco2,
Real  yh2o,
Real &  A,
Real &  B 
) const
protected

The function A (Eq.

(11) of Spycher, Pruess and Ennis-King (2003) for T <= 100C, and Eqs. (10) and (17) of Spycher and Pruess (2010) for T > 100C)

Parameters
pressuregas pressure (Pa)
temperaturefluid temperature (K)
XnaclNaCl mass fraction (kg/kg)
co2_densityCO2 density (kg/m^3)
xco2mole fraction of CO2 in liquid phase (-)
yh2omole fraction of H2O in gas phase (-)
[out]Athe function A
[out]Bthe function B

Definition at line 1254 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractions(), funcABHighTemp(), and solveEquilibriumMoleFractionHighTemp().

1262 {
1263  if (temperature <= 373.15)
1264  mooseError(name(),
1265  ": funcABHighTemp() is not valid for T <= 373.15K. Use funcABLowTemp() instead");
1266 
1267  // Pressure in bar
1268  const Real pbar = pressure * 1.0e-5;
1269  // Temperature in C
1270  const Real Tc = temperature - _T_c2k;
1271 
1272  // Reference pressure and partial molar volumes
1273  const Real pref = -1.9906e-1 + 2.0471e-3 * Tc + 1.0152e-4 * Tc * Tc - 1.4234e-6 * Tc * Tc * Tc +
1274  1.4168e-8 * Tc * Tc * Tc * Tc;
1275  const Real vCO2 = 32.6 + 3.413e-2 * (Tc - 100.0);
1276  const Real vH2O = 18.1 + 3.137e-2 * (Tc - 100.0);
1277 
1278  const Real delta_pbar = pbar - pref;
1279  const Real Rt = _Rbar * temperature;
1280 
1281  // Equilibrium constants
1282  Real K0H2O, dK0H2O_dT, K0CO2, dK0CO2_dT;
1283  equilibriumConstantH2O(temperature, K0H2O, dK0H2O_dT);
1284  equilibriumConstantCO2(temperature, K0CO2, dK0CO2_dT);
1285 
1286  // Fugacity coefficients
1287  Real phiH2O, phiCO2;
1288  fugacityCoefficientsHighTemp(pressure, temperature, co2_density, xco2, yh2o, phiCO2, phiH2O);
1289 
1290  // Activity coefficients
1291  const Real gammaH2O = activityCoefficientH2O(temperature, xco2);
1292  const Real gammaCO2 = activityCoefficientCO2(temperature, xco2);
1293 
1294  // Activity coefficient for CO2 in brine
1295  Real gamma, dgamma_dT, dgamma_dX;
1296  activityCoefficientHighTemp(temperature, Xnacl, gamma, dgamma_dT, dgamma_dX);
1297 
1298  A = K0H2O * gammaH2O / (phiH2O * pbar) * std::exp(delta_pbar * vH2O / Rt);
1299  B = phiCO2 * pbar / (_invMh2o * K0CO2 * gamma * gammaCO2) * std::exp(-delta_pbar * vCO2 / Rt);
1300 }
const Real _T_c2k
Conversion from C to K.
const Real _invMh2o
Inverse of molar mass of H2O (mol/kg)
const std::string temperature
Definition: NS.h:27
void equilibriumConstantH2O(Real temperature, Real &kh2o, Real &dkh2o_dT) const
Equilibrium constant for H2O For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2...
void activityCoefficientHighTemp(Real temperature, Real Xnacl, Real &gamma, Real &dgamma_dT, Real &dgamma_dX) const
Activity coefficient for CO2 in brine used in the elevated temperature formulation.
const std::string name
Definition: Setup.h:22
void equilibriumConstantCO2(Real temperature, Real &kco2, Real &dkco2_dT) const
Equilibrium constant for CO2 For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2...
Real activityCoefficientCO2(Real temperature, Real xco2) const
Activity coefficient of CO2 Eq.
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
Real activityCoefficientH2O(Real temperature, Real xco2) const
Activity coefficient of H2O Eq.
const std::string pressure
Definition: NS.h:26
void fugacityCoefficientsHighTemp(Real pressure, Real temperature, Real co2_density, Real xco2, Real yh2o, Real &fco2, Real &fh2o) const
Fugacity coefficients for H2O and CO2 at elevated temperatures (100C < T <= 300C).

◆ funcABHighTemp() [2/2]

void PorousFlowBrineCO2::funcABHighTemp ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real  co2_density,
Real  xco2,
Real  yh2o,
Real &  A,
Real &  dA_dp,
Real &  dA_dT,
Real &  B,
Real &  dB_dp,
Real &  dB_dT,
Real &  dB_dX 
) const
protected

Definition at line 1303 of file PorousFlowBrineCO2.C.

1316 {
1317  funcABHighTemp(pressure, temperature, Xnacl, co2_density, xco2, yh2o, A, B);
1318 
1319  // Use finite differences for derivatives in the high temperature regime
1320  const Real dp = 1.0e-2;
1321  const Real dT = 1.0e-6;
1322  const Real dX = 1.0e-8;
1323 
1324  Real A2, B2;
1325  funcABHighTemp(pressure + dp, temperature, Xnacl, co2_density, xco2, yh2o, A2, B2);
1326  dA_dp = (A2 - A) / dp;
1327  dB_dp = (B2 - B) / dp;
1328 
1329  funcABHighTemp(pressure, temperature + dT, Xnacl, co2_density, xco2, yh2o, A2, B2);
1330  dA_dT = (A2 - A) / dT;
1331  dB_dT = (B2 - B) / dT;
1332 
1333  funcABHighTemp(pressure, temperature, Xnacl + dX, co2_density, xco2, yh2o, A2, B2);
1334  dB_dX = (B2 - B) / dX;
1335 }
const std::string temperature
Definition: NS.h:27
void funcABHighTemp(Real pressure, Real temperature, Real Xnacl, Real co2_density, Real xco2, Real yh2o, Real &A, Real &B) const
The function A (Eq.
const std::string pressure
Definition: NS.h:26

◆ funcABLowTemp()

void PorousFlowBrineCO2::funcABLowTemp ( Real  pressure,
Real  temperature,
Real  co2_density,
Real  dco2_density_dp,
Real  dco2_density_dT,
Real &  A,
Real &  dA_dp,
Real &  dA_dT,
Real &  B,
Real &  dB_dp,
Real &  dB_dT 
) const
protected

Definition at line 1189 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractionsLowTemp().

1200 {
1201  if (temperature > 373.15)
1202  mooseError(name(),
1203  ": funcABLowTemp() is not valid for T > 373.15K. Use funcABHighTemp() instead");
1204 
1205  // Pressure in bar
1206  const Real pbar = pressure * 1.0e-5;
1207 
1208  // Reference pressure and partial molar volumes
1209  const Real pref = 1.0;
1210  const Real vCO2 = 32.6;
1211  const Real vH2O = 18.1;
1212 
1213  const Real delta_pbar = pbar - pref;
1214  const Real Rt = _Rbar * temperature;
1215 
1216  // Equilibrium constants
1217  Real K0H2O, dK0H2O_dT, K0CO2, dK0CO2_dT;
1218  equilibriumConstantH2O(temperature, K0H2O, dK0H2O_dT);
1219  equilibriumConstantCO2(temperature, K0CO2, dK0CO2_dT);
1220 
1221  // Fugacity coefficients
1222  Real phiH2O, dphiH2O_dp, dphiH2O_dT;
1223  Real phiCO2, dphiCO2_dp, dphiCO2_dT;
1225  temperature,
1226  co2_density,
1227  dco2_density_dp,
1228  dco2_density_dT,
1229  phiCO2,
1230  dphiCO2_dp,
1231  dphiCO2_dT,
1232  phiH2O,
1233  dphiH2O_dp,
1234  dphiH2O_dT);
1235 
1236  A = K0H2O / (phiH2O * pbar) * std::exp(delta_pbar * vH2O / Rt);
1237  B = phiCO2 * pbar / (_invMh2o * K0CO2) * std::exp(-delta_pbar * vCO2 / Rt);
1238 
1239  dA_dp = (-1.0e-5 / pbar + 1.0e-5 * vH2O / Rt - dphiH2O_dp / phiH2O) * A;
1240 
1241  dB_dp = (1.0e-5 * phiCO2 + pbar * dphiCO2_dp - 1.0e-5 * vCO2 * pbar * phiCO2 / Rt) *
1242  std::exp(-delta_pbar * vCO2 / Rt) / (_invMh2o * K0CO2);
1243 
1244  dA_dT =
1245  (dK0H2O_dT - dphiH2O_dT * K0H2O / phiH2O - delta_pbar * vH2O * K0H2O / (Rt * temperature)) *
1246  std::exp(delta_pbar * vH2O / Rt) / (pbar * phiH2O);
1247 
1248  dB_dT = (-pbar * phiCO2 * dK0CO2_dT / K0CO2 + pbar * dphiCO2_dT +
1249  delta_pbar * vCO2 * pbar * phiCO2 / (Rt * temperature)) *
1250  std::exp(-delta_pbar * vCO2 / Rt) / (_invMh2o * K0CO2);
1251 }
const Real _invMh2o
Inverse of molar mass of H2O (mol/kg)
const std::string temperature
Definition: NS.h:27
void fugacityCoefficientsLowTemp(Real pressure, Real temperature, Real co2_density, Real dco2_density_dp, Real dco2_density_dT, Real &fco2, Real &dfco2_dp, Real &dfco2_dT, Real &fh2o, Real &dfh2o_dp, Real &dfh2o_dT) const
Fugacity coefficients for H2O and CO2 for T <= 100C Eq.
void equilibriumConstantH2O(Real temperature, Real &kh2o, Real &dkh2o_dT) const
Equilibrium constant for H2O For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2...
const std::string name
Definition: Setup.h:22
void equilibriumConstantCO2(Real temperature, Real &kco2, Real &dkco2_dT) const
Equilibrium constant for CO2 For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2...
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
const std::string pressure
Definition: NS.h:26

◆ gasComponentIndex()

unsigned int PorousFlowFluidStateBase::gasComponentIndex ( ) const
inlineinherited

The index of the gas fluid component.

Returns
gas fluid component number

Definition at line 123 of file PorousFlowFluidStateBase.h.

123 { return _gas_fluid_component; };
unsigned int _gas_fluid_component
Fluid component number of the gas phase.

◆ gasPhaseIndex()

unsigned int PorousFlowFluidStateBase::gasPhaseIndex ( ) const
inlineinherited

The index of the gas phase.

Returns
gas phase number

Definition at line 111 of file PorousFlowFluidStateBase.h.

111 { return _gas_phase_number; };
unsigned int _gas_phase_number
Phase number of the gas phase.

◆ gasProperties()

void PorousFlowBrineCO2::gasProperties ( Real  pressure,
Real  temperature,
std::vector< FluidStateProperties > &  fsp 
) const

Thermophysical properties of the gaseous state.

Parameters
pressuregas pressure (Pa)
temperaturetemperature (K)
[out]FluidStateDensitydata structure

Definition at line 279 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties(), and totalMassFraction().

282 {
284 
285  // Gas density, viscosity and enthalpy are approximated with pure CO2 - no correction due
286  // to the small amount of water vapor is made
287  Real co2_density, dco2_density_dp, dco2_density_dT;
288  Real co2_viscosity, dco2_viscosity_dp, dco2_viscosity_dT;
289  Real co2_enthalpy, dco2_enthalpy_dp, dco2_enthalpy_dT;
291  temperature,
292  co2_density,
293  dco2_density_dp,
294  dco2_density_dT,
295  co2_viscosity,
296  dco2_viscosity_dp,
297  dco2_viscosity_dT);
298 
299  _co2_fp.h_from_p_T(pressure, temperature, co2_enthalpy, dco2_enthalpy_dp, dco2_enthalpy_dT);
300 
301  // Save the values to the FluidStateProperties object. Note that derivatives wrt z are 0
302  gas.density = co2_density;
303  gas.ddensity_dp = dco2_density_dp;
304  gas.ddensity_dT = dco2_density_dT;
305  gas.ddensity_dZ = 0.0;
306 
307  gas.viscosity = co2_viscosity;
308  gas.dviscosity_dp = dco2_viscosity_dp;
309  gas.dviscosity_dT = dco2_viscosity_dT;
310  gas.dviscosity_dZ = 0.0;
311 
312  gas.enthalpy = co2_enthalpy;
313  gas.denthalpy_dp = dco2_enthalpy_dp;
314  gas.denthalpy_dT = dco2_enthalpy_dT;
315  gas.denthalpy_dZ = 0.0;
316 }
unsigned int _gas_phase_number
Phase number of the gas phase.
const std::string temperature
Definition: NS.h:27
const SinglePhaseFluidProperties & _co2_fp
Fluid properties UserObject for the CO2.
Data structure to pass calculated thermophysical properties.
const std::string pressure
Definition: NS.h:26
virtual void rho_mu_from_p_T(Real pressure, Real temperature, Real &rho, Real &mu) const

◆ henryConstant()

void PorousFlowBrineCO2::henryConstant ( Real  temperature,
Real  Xnacl,
Real &  Kh,
Real &  dKh_dT,
Real &  dKh_dx 
) const

Henry's constant of dissolution of gas phase CO2 in brine.

From Battistelli et al, A fluid property module for the TOUGH2 simulator for saline brines with non-condensible gas, Proc. Eighteenth Workshop on Geothermal Reservoir Engineering (1993)

Parameters
temperaturefluid temperature (K)
XnaclNaCl mass fraction (kg/kg)
[out]KhHenry's constant (Pa)
[out]dKh_dTderivative of Henry's constant wrt temperature
[out]dKh_dxderivative of Henry's constant wrt salt mass fraction

Definition at line 1467 of file PorousFlowBrineCO2.C.

Referenced by enthalpyOfDissolutionGas().

1469 {
1470  // Henry's constant for dissolution in water
1471  Real Kh_h2o, dKh_h2o_dT;
1472  _co2_fp.henryConstant(temperature, Kh_h2o, dKh_h2o_dT);
1473 
1474  // The correction to salt is obtained through the salting out coefficient
1475  const std::vector<Real> b{1.19784e-1, -7.17823e-4, 4.93854e-6, -1.03826e-8, 1.08233e-11};
1476 
1477  // Need temperature in Celcius
1478  const Real Tc = temperature - _T_c2k;
1479 
1480  Real kb = 0.0;
1481  for (unsigned int i = 0; i < b.size(); ++i)
1482  kb += b[i] * MathUtils::pow(Tc, i);
1483 
1484  Real dkb_dT = 0.0;
1485  for (unsigned int i = 1; i < b.size(); ++i)
1486  dkb_dT += i * b[i] * MathUtils::pow(Tc, i - 1);
1487 
1488  // Need salt mass fraction in molality
1489  const Real xmol = Xnacl / (1.0 - Xnacl) / _Mnacl;
1490  const Real dxmol_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
1491  // Henry's constant and its derivative wrt temperature and salt mass fraction
1492  Kh = Kh_h2o * std::pow(10.0, xmol * kb);
1493  dKh_dT = (xmol * std::log(10.0) * dkb_dT + dKh_h2o_dT / Kh_h2o) * Kh;
1494  dKh_dX = std::log(10.0) * kb * dxmol_dX * Kh;
1495 }
const Real _T_c2k
Conversion from C to K.
const Real _Mnacl
Molar mass of NaCL.
const std::string temperature
Definition: NS.h:27
const SinglePhaseFluidProperties & _co2_fp
Fluid properties UserObject for the CO2.
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
virtual Real henryConstant(Real temperature) const
Henry&#39;s law constant for dissolution in water.

◆ initialize()

void PorousFlowFluidStateFlash::initialize ( )
inlinefinalinherited

Definition at line 36 of file PorousFlowFluidStateFlash.h.

36 {};

◆ liquidProperties()

void PorousFlowBrineCO2::liquidProperties ( Real  pressure,
Real  temperature,
Real  Xnacl,
std::vector< FluidStateProperties > &  fsp 
) const

Thermophysical properties of the liquid state.

Parameters
pressureliquid pressure (Pa)
temperaturetemperature (K)
XnaclNaCl mass fraction (kg/kg)
[out]FluidStateDensitydata structure

Definition at line 319 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties(), and totalMassFraction().

323 {
325 
326  // The liquid density includes the density increase due to dissolved CO2
327  Real brine_density, dbrine_density_dp, dbrine_density_dT, dbrine_density_dX;
329  temperature,
330  Xnacl,
331  brine_density,
332  dbrine_density_dp,
333  dbrine_density_dT,
334  dbrine_density_dX);
335 
336  // Mass fraction of CO2 in liquid phase
337  const Real Xco2 = liquid.mass_fraction[_gas_fluid_component];
338  const Real dXco2_dp = liquid.dmass_fraction_dp[_gas_fluid_component];
339  const Real dXco2_dT = liquid.dmass_fraction_dT[_gas_fluid_component];
340  const Real dXco2_dZ = liquid.dmass_fraction_dZ[_gas_fluid_component];
341  const Real dXco2_dX = liquid.dmass_fraction_dX[_gas_fluid_component];
342 
343  // The liquid density
344  Real co2_partial_density, dco2_partial_density_dT;
345  partialDensityCO2(temperature, co2_partial_density, dco2_partial_density_dT);
346 
347  const Real liquid_density = 1.0 / (Xco2 / co2_partial_density + (1.0 - Xco2) / brine_density);
348 
349  const Real dliquid_density_dp =
350  (dXco2_dp / brine_density + (1.0 - Xco2) * dbrine_density_dp / brine_density / brine_density -
351  dXco2_dp / co2_partial_density) *
352  liquid_density * liquid_density;
353 
354  const Real dliquid_density_dT =
355  (dXco2_dT / brine_density + (1.0 - Xco2) * dbrine_density_dT / brine_density / brine_density -
356  dXco2_dT / co2_partial_density +
357  Xco2 * dco2_partial_density_dT / co2_partial_density / co2_partial_density) *
358  liquid_density * liquid_density;
359 
360  const Real dliquid_density_dZ =
361  (dXco2_dZ / brine_density - dXco2_dZ / co2_partial_density) * liquid_density * liquid_density;
362 
363  const Real dliquid_density_dX =
364  (dXco2_dX / brine_density + (1.0 - Xco2) * dbrine_density_dX / brine_density / brine_density -
365  dXco2_dX / co2_partial_density) *
366  liquid_density * liquid_density;
367 
368  // Assume that liquid viscosity is just the brine viscosity
369  Real liquid_viscosity, dliquid_viscosity_dp, dliquid_viscosity_dT, dliquid_viscosity_dX;
371  temperature,
372  Xnacl,
373  liquid_viscosity,
374  dliquid_viscosity_dp,
375  dliquid_viscosity_dT,
376  dliquid_viscosity_dX);
377 
378  // Liquid enthalpy (including contribution due to the enthalpy of dissolution)
379  Real brine_enthalpy, dbrine_enthalpy_dp, dbrine_enthalpy_dT, dbrine_enthalpy_dX;
381  temperature,
382  Xnacl,
383  brine_enthalpy,
384  dbrine_enthalpy_dp,
385  dbrine_enthalpy_dT,
386  dbrine_enthalpy_dX);
387 
388  // Enthalpy of CO2
389  Real co2_enthalpy, dco2_enthalpy_dp, dco2_enthalpy_dT;
390  _co2_fp.h_from_p_T(pressure, temperature, co2_enthalpy, dco2_enthalpy_dp, dco2_enthalpy_dT);
391 
392  // Enthalpy of dissolution
393  Real hdis, dhdis_dT;
394  enthalpyOfDissolution(temperature, hdis, dhdis_dT);
395 
396  const Real liquid_enthalpy = (1.0 - Xco2) * brine_enthalpy + Xco2 * (co2_enthalpy + hdis);
397  const Real dliquid_enthalpy_dp = (1.0 - Xco2) * dbrine_enthalpy_dp + Xco2 * dco2_enthalpy_dp +
398  dXco2_dp * (co2_enthalpy + hdis - brine_enthalpy);
399  const Real dliquid_enthalpy_dT = (1.0 - Xco2) * dbrine_enthalpy_dT +
400  Xco2 * (dco2_enthalpy_dT + dhdis_dT) +
401  dXco2_dT * (co2_enthalpy + hdis - brine_enthalpy);
402  const Real dliquid_enthalpy_dZ = dXco2_dZ * (co2_enthalpy + hdis - brine_enthalpy);
403  const Real dliquid_enthalpy_dX =
404  (1.0 - Xco2) * dbrine_enthalpy_dX + dXco2_dX * (co2_enthalpy + hdis - brine_enthalpy);
405 
406  // Save the values to the FluidStateProperties object
407  liquid.density = liquid_density;
408  liquid.ddensity_dp = dliquid_density_dp;
409  liquid.ddensity_dT = dliquid_density_dT;
410  liquid.ddensity_dZ = dliquid_density_dZ;
411  liquid.ddensity_dX = dliquid_density_dX;
412 
413  liquid.viscosity = liquid_viscosity;
414  liquid.dviscosity_dp = dliquid_viscosity_dp;
415  liquid.dviscosity_dT = dliquid_viscosity_dT;
416  liquid.dviscosity_dZ = 0.0;
417  liquid.dviscosity_dX = dliquid_viscosity_dX;
418 
419  liquid.enthalpy = liquid_enthalpy;
420  liquid.denthalpy_dp = dliquid_enthalpy_dp;
421  liquid.denthalpy_dT = dliquid_enthalpy_dT;
422  liquid.denthalpy_dZ = dliquid_enthalpy_dZ;
423  liquid.denthalpy_dX = dliquid_enthalpy_dX;
424 }
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
virtual Real rho_from_p_T_X(Real pressure, Real temperature, Real xnacl) const override
const std::string temperature
Definition: NS.h:27
virtual Real mu_from_p_T_X(Real pressure, Real temperature, Real xnacl) const override
std::vector< Real > dmass_fraction_dT
void enthalpyOfDissolution(Real temperature, Real &hdis, Real &dhdis_dT) const
Enthalpy of dissolution of CO2 in brine calculated using linear fit to model of Duan and Sun...
const BrineFluidProperties & _brine_fp
Fluid properties UserObject for water.
const SinglePhaseFluidProperties & _co2_fp
Fluid properties UserObject for the CO2.
std::vector< Real > mass_fraction
std::vector< Real > dmass_fraction_dX
Data structure to pass calculated thermophysical properties.
std::vector< Real > dmass_fraction_dZ
unsigned int _gas_fluid_component
Fluid component number of the gas phase.
std::vector< Real > dmass_fraction_dp
const std::string pressure
Definition: NS.h:26
void partialDensityCO2(Real temperature, Real &partial_density, Real &dpartial_density_dT) const
Partial density of dissolved CO2 From Garcia, Density of aqueous solutions of CO2, LBNL-49023 (2001)
FPDualReal h_from_p_T_X(const FPDualReal &pressure, const FPDualReal &temperature, const FPDualReal &xnacl) const override

◆ massFractions()

void PorousFlowBrineCO2::massFractions ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real  Z,
FluidStatePhaseEnum phase_state,
std::vector< FluidStateProperties > &  fsp 
) const

Mass fractions of CO2 and H2O in both phases, as well as derivatives wrt PorousFlow variables.

Values depend on the phase state (liquid, gas or two phase)

Parameters
pressurephase pressure (Pa)
temperaturephase temperature (K)
XnaclNaCl mass fraction (kg/kg)
Ztotal mass fraction of CO2 component
[out]PhaseStateEnumcurrent phase state
[out]FluidStateMassFractionsdata structure

Definition at line 160 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties().

166 {
169 
170  Real Xco2 = 0.0, dXco2_dp = 0.0, dXco2_dT = 0.0, dXco2_dX = 0.0;
171  Real Yh2o = 0.0, dYh2o_dp = 0.0, dYh2o_dT = 0.0, dYh2o_dX = 0.0;
172  Real Yco2 = 0.0, dYco2_dp = 0.0, dYco2_dT = 0.0, dYco2_dX = 0.0;
173 
174  // If the amount of CO2 is less than the smallest solubility, then all CO2 will
175  // be dissolved, and the equilibrium mass fractions do not need to be computed
176  if (Z < _Zmin)
177  phase_state = FluidStatePhaseEnum::LIQUID;
178 
179  else
180  {
181  // Equilibrium mass fraction of CO2 in liquid and H2O in gas phases
183  temperature,
184  Xnacl,
185  Xco2,
186  dXco2_dp,
187  dXco2_dT,
188  dXco2_dX,
189  Yh2o,
190  dYh2o_dp,
191  dYh2o_dT,
192  dYh2o_dX);
193 
194  Yco2 = 1.0 - Yh2o;
195  dYco2_dp = -dYh2o_dp;
196  dYco2_dT = -dYh2o_dT;
197  dYco2_dX = -dYh2o_dX;
198 
199  // Determine which phases are present based on the value of z
200  phaseState(Z, Xco2, Yco2, phase_state);
201  }
202 
203  // The equilibrium mass fractions calculated above are only correct in the two phase
204  // state. If only liquid or gas phases are present, the mass fractions are given by
205  // the total mass fraction z
206  Real Xh2o = 0.0;
207  Real dXco2_dZ = 0.0, dYco2_dZ = 0.0;
208 
209  switch (phase_state)
210  {
212  {
213  Xco2 = Z;
214  Yco2 = 0.0;
215  Xh2o = 1.0 - Z;
216  Yh2o = 0.0;
217  dXco2_dp = 0.0;
218  dXco2_dT = 0.0;
219  dXco2_dX = 0.0;
220  dXco2_dZ = 1.0;
221  dYco2_dp = 0.0;
222  dYco2_dT = 0.0;
223  dYco2_dX = 0.0;
224  break;
225  }
226 
228  {
229  Xco2 = 0.0;
230  Yco2 = Z;
231  Yh2o = 1.0 - Z;
232  dXco2_dp = 0.0;
233  dXco2_dT = 0.0;
234  dXco2_dX = 0.0;
235  dYco2_dZ = 1.0;
236  dYco2_dp = 0.0;
237  dYco2_dT = 0.0;
238  dYco2_dX = 0.0;
239  break;
240  }
241 
243  {
244  // Keep equilibrium mass fractions
245  Xh2o = 1.0 - Xco2;
246  break;
247  }
248  }
249 
250  // Save the mass fractions in the FluidStateProperties object
252  liquid.mass_fraction[_gas_fluid_component] = Xco2;
253  liquid.mass_fraction[_salt_component] = Xnacl;
256 
257  // Save the derivatives wrt PorousFlow variables
258  liquid.dmass_fraction_dp[_aqueous_fluid_component] = -dXco2_dp;
259  liquid.dmass_fraction_dp[_gas_fluid_component] = dXco2_dp;
260  liquid.dmass_fraction_dT[_aqueous_fluid_component] = -dXco2_dT;
261  liquid.dmass_fraction_dT[_gas_fluid_component] = dXco2_dT;
262  liquid.dmass_fraction_dX[_aqueous_fluid_component] = -dXco2_dX;
263  liquid.dmass_fraction_dX[_gas_fluid_component] = dXco2_dX;
264  liquid.dmass_fraction_dX[_salt_component] = 1.0;
265  liquid.dmass_fraction_dZ[_aqueous_fluid_component] = -dXco2_dZ;
266  liquid.dmass_fraction_dZ[_gas_fluid_component] = dXco2_dZ;
267 
269  gas.dmass_fraction_dp[_gas_fluid_component] = dYco2_dp;
271  gas.dmass_fraction_dT[_gas_fluid_component] = dYco2_dT;
273  gas.dmass_fraction_dX[_gas_fluid_component] = dYco2_dX;
275  gas.dmass_fraction_dZ[_gas_fluid_component] = dYco2_dZ;
276 }
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
const unsigned int _salt_component
Salt component index.
unsigned int _gas_phase_number
Phase number of the gas phase.
const unsigned int _aqueous_fluid_component
Fluid component number of the aqueous component.
void phaseState(Real Zi, Real Xi, Real Yi, FluidStatePhaseEnum &phase_state) const
Determines the phase state gven the total mass fraction and equilibrium mass fractions.
const std::string temperature
Definition: NS.h:27
std::vector< Real > dmass_fraction_dT
std::vector< Real > mass_fraction
std::vector< Real > dmass_fraction_dX
Data structure to pass calculated thermophysical properties.
std::vector< Real > dmass_fraction_dZ
unsigned int _gas_fluid_component
Fluid component number of the gas phase.
void equilibriumMassFractions(Real pressure, Real temperature, Real Xnacl, Real &Xco2, Real &dXco2_dp, Real &dXco2_dT, Real &dXco2_dX, Real &Yh2o, Real &dYh2o_dp, Real &dYh2o_dT, Real &dYh2o_dX) const
Mass fractions of CO2 in brine and water vapor in CO2 at equilibrium.
std::vector< Real > dmass_fraction_dp
const std::string pressure
Definition: NS.h:26
const Real _Zmin
Minimum Z - below this value all CO2 will be dissolved.

◆ numComponents()

unsigned int PorousFlowFluidStateBase::numComponents ( ) const
inlineinherited

The maximum number of components in this model.

Returns
number of components

Definition at line 99 of file PorousFlowFluidStateBase.h.

99 { return _num_components; };
unsigned int _num_components
Number of components.

◆ numPhases()

unsigned int PorousFlowFluidStateBase::numPhases ( ) const
inlineinherited

The maximum number of phases in this model.

Returns
number of phases

Definition at line 93 of file PorousFlowFluidStateBase.h.

Referenced by PorousFlowFluidState::PorousFlowFluidState().

93 { return _num_phases; };
unsigned int _num_phases
Number of phases.

◆ partialDensityCO2()

void PorousFlowBrineCO2::partialDensityCO2 ( Real  temperature,
Real &  partial_density,
Real &  dpartial_density_dT 
) const

Partial density of dissolved CO2 From Garcia, Density of aqueous solutions of CO2, LBNL-49023 (2001)

Parameters
temperaturefluid temperature (K)
[out]partialmolar density (kg/m^3)
[out]derivativeof partial molar density wrt temperature

Definition at line 1404 of file PorousFlowBrineCO2.C.

Referenced by liquidProperties(), and saturationTwoPhase().

1407 {
1408  // This correlation uses temperature in C
1409  const Real Tc = temperature - _T_c2k;
1410  // The parial molar volume
1411  const Real V = 37.51 - 9.585e-2 * Tc + 8.74e-4 * Tc * Tc - 5.044e-7 * Tc * Tc * Tc;
1412  const Real dV_dT = -9.585e-2 + 1.748e-3 * Tc - 1.5132e-6 * Tc * Tc;
1413 
1414  partial_density = 1.0e6 * _Mco2 / V;
1415  dpartial_density_dT = -1.0e6 * _Mco2 * dV_dT / V / V;
1416 }
const Real _T_c2k
Conversion from C to K.
const Real _Mco2
Molar mass of CO2 (kg/mol)
const std::string temperature
Definition: NS.h:27

◆ phaseState()

void PorousFlowFluidStateFlash::phaseState ( Real  Zi,
Real  Xi,
Real  Yi,
FluidStatePhaseEnum phase_state 
) const
protectedinherited

Determines the phase state gven the total mass fraction and equilibrium mass fractions.

Parameters
Zitotal mass fraction
Xiequilibrium mass fraction in liquid
Yiequilibrium mass fraction in gas
[out]phase_statethe phase state (gas, liquid, two phase)

Definition at line 29 of file PorousFlowFluidStateFlash.C.

Referenced by PorousFlowWaterNCG::massFractions(), and massFractions().

33 {
34  if (Zi <= Xi)
35  {
36  // In this case, there is not enough component i to form a gas phase,
37  // so only a liquid phase is present
38  phase_state = FluidStatePhaseEnum::LIQUID;
39  }
40  else if (Zi > Xi && Zi < Yi)
41  {
42  // Two phases are present
43  phase_state = FluidStatePhaseEnum::TWOPHASE;
44  }
45  else // (Zi >= Yi)
46  {
47  // In this case, there is not enough water to form a liquid
48  // phase, so only a gas phase is present
49  phase_state = FluidStatePhaseEnum::GAS;
50  }
51 }

◆ rachfordRice()

Real PorousFlowFluidStateFlash::rachfordRice ( Real  vf,
std::vector< Real > &  Zi,
std::vector< Real > &  Ki 
) const
inherited

Rachford-Rice equation for vapor fraction.

Can be solved analytically for two components in two phases, but must be solved iteratively using a root finding algorithm for more components. This equation has the nice property that it is monotonic in the interval [0,1], so that only a small number of iterations are typically required to find the root.

The Rachford-Rice equation can also be used to check whether the phase state is two phase, single phase gas, or single phase liquid. Evaluate f(v), the Rachford-Rice equation evaluated at the vapor mass fraction.

If f(0) < 0, then the mixture is below the bubble point, and only a single phase liquid can exist

If f(1) > 0, then the mixture is above the dew point, and only a single phase gas exists.

If f(0) >= 0 and f(1) <= 0, the mixture is between the bubble and dew points, and both gas and liquid phases exist.

Parameters
vfvapor fraction
Zimass fractions
Kiequilibrium constants
Returns
f(x)

Definition at line 54 of file PorousFlowFluidStateFlash.C.

Referenced by PorousFlowFluidStateFlash::vaporMassFraction().

57 {
58  const std::size_t num_z = Zi.size();
59  // Check that the sizes of the mass fractions and equilibrium constant vectors are correct
60  if (Ki.size() != num_z + 1)
61  mooseError("The number of mass fractions or equilibrium components passed to rachfordRice is "
62  "not correct");
63 
64  Real f = 0.0;
65  Real Z_total = 0.0;
66 
67  for (std::size_t i = 0; i < num_z; ++i)
68  {
69  f += Zi[i] * (Ki[i] - 1.0) / (1.0 + x * (Ki[i] - 1.0));
70  Z_total += Zi[i];
71  }
72 
73  // Add the last component (with total mass fraction = 1 - z_total)
74  f += (1.0 - Z_total) * (Ki[num_z] - 1.0) / (1.0 + x * (Ki[num_z] - 1.0));
75 
76  return f;
77 }

◆ rachfordRiceDeriv()

Real PorousFlowFluidStateFlash::rachfordRiceDeriv ( Real  vf,
std::vector< Real > &  Zi,
std::vector< Real > &  Ki 
) const
inherited

Derivative of Rachford-Rice equation wrt vapor fraction.

Has the nice property that it is strictly negative in the interval [0,1]

Parameters
vfvapor fraction
Zimass fractions
Kiequilibrium constants
Returns
f'(x)

Definition at line 80 of file PorousFlowFluidStateFlash.C.

Referenced by PorousFlowFluidStateFlash::vaporMassFraction().

83 {
84  const std::size_t num_Z = Zi.size();
85  // Check that the sizes of the mass fractions and equilibrium constant vectors are correct
86  if (Ki.size() != num_Z + 1)
87  mooseError("The number of mass fractions or equilibrium components passed to rachfordRice is "
88  "not correct");
89 
90  Real df = 0.0;
91  Real Z_total = 0.0;
92 
93  for (std::size_t i = 0; i < num_Z; ++i)
94  {
95  df -= Zi[i] * (Ki[i] - 1.0) * (Ki[i] - 1.0) / (1.0 + x * (Ki[i] - 1.0)) /
96  (1.0 + x * (Ki[i] - 1.0));
97  Z_total += Zi[i];
98  }
99 
100  // Add the last component (with total mass fraction = 1 - z_total)
101  df -= (1.0 - Z_total) * (Ki[num_Z] - 1.0) * (Ki[num_Z] - 1.0) / (1.0 + x * (Ki[num_Z] - 1.0)) /
102  (1.0 + x * (Ki[num_Z] - 1.0));
103 
104  return df;
105 }

◆ saltComponentIndex()

unsigned int PorousFlowBrineCO2::saltComponentIndex ( ) const
inline

The index of the salt component.

Returns
salt component number

Definition at line 357 of file PorousFlowBrineCO2.h.

357 { return _salt_component; };
const unsigned int _salt_component
Salt component index.

◆ saturationTwoPhase()

void PorousFlowBrineCO2::saturationTwoPhase ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real  Z,
std::vector< FluidStateProperties > &  fsp 
) const

Gas and liquid saturations for the two-phase region.

Parameters
pressuregas pressure (Pa)
temperaturephase temperature (K)
XnaclNaCl mass fraction (kg/kg)
Ztotal mass fraction of CO2 component
[out]FluidStateSaturationdata structure

Definition at line 427 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties().

432 {
435 
436  // Approximate liquid density as saturation isn't known yet
437  Real brine_density, dbrine_density_dp, dbrine_density_dT, dbrine_density_dX;
439  temperature,
440  Xnacl,
441  brine_density,
442  dbrine_density_dp,
443  dbrine_density_dT,
444  dbrine_density_dX);
445 
446  // Mass fraction of CO2 in liquid phase
447  const Real Xco2 = liquid.mass_fraction[_gas_fluid_component];
448  const Real dXco2_dp = liquid.dmass_fraction_dp[_gas_fluid_component];
449  const Real dXco2_dT = liquid.dmass_fraction_dT[_gas_fluid_component];
450  const Real dXco2_dX = liquid.dmass_fraction_dX[_gas_fluid_component];
451 
452  // The liquid density
453  Real co2_partial_density, dco2_partial_density_dT;
454  partialDensityCO2(temperature, co2_partial_density, dco2_partial_density_dT);
455 
456  const Real liquid_density = 1.0 / (Xco2 / co2_partial_density + (1.0 - Xco2) / brine_density);
457 
458  const Real dliquid_density_dp =
459  (dXco2_dp / brine_density + (1.0 - Xco2) * dbrine_density_dp / brine_density / brine_density -
460  dXco2_dp / co2_partial_density) *
461  liquid_density * liquid_density;
462 
463  const Real dliquid_density_dT =
464  (dXco2_dT / brine_density + (1.0 - Xco2) * dbrine_density_dT / brine_density / brine_density -
465  dXco2_dT / co2_partial_density +
466  Xco2 * dco2_partial_density_dT / co2_partial_density / co2_partial_density) *
467  liquid_density * liquid_density;
468 
469  const Real dliquid_density_dX =
470  (dXco2_dX / brine_density + (1.0 - Xco2) * dbrine_density_dX / brine_density / brine_density -
471  dXco2_dX / co2_partial_density) *
472  liquid_density * liquid_density;
473 
474  const Real Yco2 = gas.mass_fraction[_gas_fluid_component];
475  const Real dYco2_dp = gas.dmass_fraction_dp[_gas_fluid_component];
476  const Real dYco2_dT = gas.dmass_fraction_dT[_gas_fluid_component];
477  const Real dYco2_dX = gas.dmass_fraction_dX[_gas_fluid_component];
478 
479  // Set mass equilibrium constants used in the calculation of vapor mass fraction
480  const Real K0 = Yco2 / Xco2;
481  const Real K1 = (1.0 - Yco2) / (1.0 - Xco2);
482  const Real vapor_mass_fraction = vaporMassFraction(Z, K0, K1);
483 
484  // The gas saturation in the two phase case
485  gas.saturation = vapor_mass_fraction * liquid_density /
486  (gas.density + vapor_mass_fraction * (liquid_density - gas.density));
487 
488  const Real dv_dZ = (K1 - K0) / ((K0 - 1.0) * (K1 - 1.0));
489  const Real denominator = (gas.density + vapor_mass_fraction * (liquid_density - gas.density)) *
490  (gas.density + vapor_mass_fraction * (liquid_density - gas.density));
491 
492  const Real ds_dZ = gas.density * liquid_density * dv_dZ / denominator;
493 
494  const Real dK0_dp = (Xco2 * dYco2_dp - Yco2 * dXco2_dp) / Xco2 / Xco2;
495  const Real dK0_dT = (Xco2 * dYco2_dT - Yco2 * dXco2_dT) / Xco2 / Xco2;
496  const Real dK0_dX = (Xco2 * dYco2_dX - Yco2 * dXco2_dX) / Xco2 / Xco2;
497 
498  const Real dK1_dp =
499  ((1.0 - Yco2) * dXco2_dp - (1.0 - Xco2) * dYco2_dp) / (1.0 - Xco2) / (1.0 - Xco2);
500  const Real dK1_dT =
501  ((1.0 - Yco2) * dXco2_dT - (1.0 - Xco2) * dYco2_dT) / (1.0 - Xco2) / (1.0 - Xco2);
502  const Real dK1_dX =
503  ((1.0 - Yco2) * dXco2_dX - (1.0 - Xco2) * dYco2_dX) / (1.0 - Xco2) / (1.0 - Xco2);
504 
505  const Real dv_dp =
506  Z * dK1_dp / (K1 - 1.0) / (K1 - 1.0) + (1.0 - Z) * dK0_dp / (K0 - 1.0) / (K0 - 1.0);
507 
508  Real ds_dp = gas.density * liquid_density * dv_dp +
509  vapor_mass_fraction * (1.0 - vapor_mass_fraction) *
510  (gas.density * dliquid_density_dp - gas.ddensity_dp * liquid_density);
511  ds_dp /= denominator;
512 
513  const Real dv_dT =
514  Z * dK1_dT / (K1 - 1.0) / (K1 - 1.0) + (1.0 - Z) * dK0_dT / (K0 - 1.0) / (K0 - 1.0);
515 
516  Real ds_dT = gas.density * liquid_density * dv_dT +
517  vapor_mass_fraction * (1.0 - vapor_mass_fraction) *
518  (gas.density * dliquid_density_dT - gas.ddensity_dT * liquid_density);
519  ds_dT /= denominator;
520 
521  const Real dv_dX =
522  Z * dK1_dX / (K1 - 1.0) / (K1 - 1.0) + (1.0 - Z) * dK0_dX / (K0 - 1.0) / (K0 - 1.0);
523 
524  Real ds_dX = gas.density * liquid_density * dv_dX + vapor_mass_fraction *
525  (1.0 - vapor_mass_fraction) *
526  (gas.density * dliquid_density_dX);
527  ds_dX /= denominator;
528 
529  gas.dsaturation_dp = ds_dp;
530  gas.dsaturation_dT = ds_dT;
531  gas.dsaturation_dZ = ds_dZ;
532  gas.dsaturation_dX = ds_dX;
533 }
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
virtual Real rho_from_p_T_X(Real pressure, Real temperature, Real xnacl) const override
unsigned int _gas_phase_number
Phase number of the gas phase.
const std::string temperature
Definition: NS.h:27
std::vector< Real > dmass_fraction_dT
const BrineFluidProperties & _brine_fp
Fluid properties UserObject for water.
std::vector< Real > mass_fraction
std::vector< Real > dmass_fraction_dX
Data structure to pass calculated thermophysical properties.
unsigned int _gas_fluid_component
Fluid component number of the gas phase.
std::vector< Real > dmass_fraction_dp
Real vaporMassFraction(Real Z0, Real K0, Real K1) const
Solves Rachford-Rice equation to provide vapor mass fraction.
const std::string pressure
Definition: NS.h:26
void partialDensityCO2(Real temperature, Real &partial_density, Real &dpartial_density_dT) const
Partial density of dissolved CO2 From Garcia, Density of aqueous solutions of CO2, LBNL-49023 (2001)

◆ smoothCubicInterpolation()

void PorousFlowBrineCO2::smoothCubicInterpolation ( Real  temperature,
Real  f0,
Real  df0,
Real  f1,
Real  df1,
Real &  value,
Real &  deriv 
) const
protected

Cubic function to smoothly interpolate between the low temperature and elevated temperature models for 99C < T < 109C.

Parameters
temperaturetemperature (K)
f0function value at T = 372K (99C)
df0derivative of function at T = 372K (99C)
f1function value at T = 382K (109C)
df1derivative of function at T = 382K (109C)
[out]valuevalue at the given temperature
[out]derivderivative at the given temperature

Definition at line 1534 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractions().

1536 {
1537  // Coefficients of cubic polynomial
1538  const Real dT = _Tupper - _Tlower;
1539 
1540  const Real a = f0;
1541  const Real b = df0 * dT;
1542  const Real c = 3.0 * (f1 - f0) - (2.0 * df0 + df1) * dT;
1543  const Real d = 2.0 * (f0 - f1) + (df0 + df1) * dT;
1544 
1545  const Real t2 = temperature * temperature;
1546  const Real t3 = temperature * t2;
1547 
1548  value = a + b * temperature + c * t2 + d * t3;
1549  deriv = b + 2.0 * c * temperature + 3.0 * d * t2;
1550 }
const Real _Tlower
Temperature below which the Spycher, Pruess & Ennis-King (2003) model is used (K) ...
const std::string temperature
Definition: NS.h:27
const Real _Tupper
Temperature above which the Spycher & Pruess (2010) model is used (K)

◆ solveEquilibriumMoleFractionHighTemp()

void PorousFlowBrineCO2::solveEquilibriumMoleFractionHighTemp ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real  co2_density,
Real &  xco2,
Real &  yh2o 
) const

Function to solve for yh2o and xco2 iteratively in the elevated temperature regime (T > 100C)

Parameters
pressuregas pressure (Pa)
temperaturefluid temperature (K)
XnaclNaCl mass fraction (kg/kg)
co2_densityCO2 density (kg/m^3)
[out]xco2mole fraction of CO2 in liquid phase (-)
[out]yh2omole fraction of H2O in gas phase (-)

Definition at line 1338 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractions().

1340 {
1341  // Initial guess for yh2o and xco2 (from Spycher and Pruess (2010))
1342  Real y = _brine_fp.vaporPressure(temperature, 0.0) / pressure;
1343  Real x = 0.009;
1344 
1345  // Need salt mass fraction in molality
1346  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
1347 
1348  // If y > 1, then just use y = 1, x = 0 (only a gas phase)
1349  if (y >= 1.0)
1350  {
1351  y = 1.0;
1352  x = 0.0;
1353  }
1354  else
1355  {
1356  // Residual function for Newton-Raphson
1357  auto fy = [mnacl, this](Real y, Real A, Real B) {
1358  return y -
1359  (1.0 - B) * _invMh2o / ((1.0 / A - B) * (2.0 * mnacl + _invMh2o) + 2.0 * mnacl * B);
1360  };
1361 
1362  // Derivative of fy wrt y
1363  auto dfy = [mnacl, this](Real A, Real B, Real dA, Real dB) {
1364  const Real denominator = (1.0 / A - B) * (2.0 * mnacl + _invMh2o) + 2.0 * mnacl * B;
1365  return 1.0 + _invMh2o * dB / denominator +
1366  (1.0 - B) * _invMh2o *
1367  (2.0 * mnacl * dB - (2.0 * mnacl + _invMh2o) * (dB + dA / A / A)) / denominator /
1368  denominator;
1369  };
1370 
1371  Real A, B;
1372  Real dA, dB;
1373  const Real dy = 1.0e-8;
1374 
1375  // Solve for yh2o using Newton-Raphson method
1376  unsigned int iter = 0;
1377  const Real tol = 1.0e-12;
1378  const unsigned int max_its = 10;
1379  funcABHighTemp(pressure, temperature, Xnacl, co2_density, x, y, A, B);
1380 
1381  while (std::abs(fy(y, A, B)) > tol)
1382  {
1383  funcABHighTemp(pressure, temperature, Xnacl, co2_density, x, y, A, B);
1384  // Finite difference derivatives of A and B wrt y
1385  funcABHighTemp(pressure, temperature, Xnacl, co2_density, x, y + dy, dA, dB);
1386  dA = (dA - A) / dy;
1387  dB = (dB - B) / dy;
1388 
1389  y = y - fy(y, A, B) / dfy(A, B, dA, dB);
1390 
1391  x = B * (1.0 - y);
1392 
1393  // Break if not converged and just use the value
1394  if (iter > max_its)
1395  break;
1396  }
1397  }
1398 
1399  yh2o = y;
1400  xco2 = x;
1401 }
const double tol
Definition: Setup.h:19
const Real _Mnacl
Molar mass of NaCL.
const Real _invMh2o
Inverse of molar mass of H2O (mol/kg)
const std::string temperature
Definition: NS.h:27
const BrineFluidProperties & _brine_fp
Fluid properties UserObject for water.
Real vaporPressure(Real temperature, Real xnacl) const
Brine vapour pressure From Haas, Physical properties of the coexisting phases and thermochemical prop...
void funcABHighTemp(Real pressure, Real temperature, Real Xnacl, Real co2_density, Real xco2, Real yh2o, Real &A, Real &B) const
The function A (Eq.
const std::string pressure
Definition: NS.h:26

◆ thermophysicalProperties()

void PorousFlowBrineCO2::thermophysicalProperties ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real  Z,
unsigned int  qp,
std::vector< FluidStateProperties > &  fsp 
) const
overridevirtual

Determines the complete thermophysical state of the system for a given set of primary variables.

Parameters
pressuregas phase pressure (Pa)
temperaturefluid temperature (K)
Xnaclmass fraction of NaCl
Ztotal mass fraction of fluid component
qpquadpoint index
[out]fspthe FluidStateProperties struct containing all properties

Implements PorousFlowFluidStateBase.

Definition at line 90 of file PorousFlowBrineCO2.C.

96 {
99 
100  // Check whether the input temperature is within the region of validity
102 
103  // Clear all of the FluidStateProperties data
105 
106  FluidStatePhaseEnum phase_state;
107  massFractions(pressure, temperature, Xnacl, Z, phase_state, fsp);
108 
109  switch (phase_state)
110  {
112  {
113  // Set the gas saturations
114  gas.saturation = 1.0;
115 
116  // Calculate gas properties
118 
119  break;
120  }
121 
123  {
124  // Calculate the liquid properties
125  Real liquid_pressure = pressure - _pc.capillaryPressure(1.0, qp);
126  liquidProperties(liquid_pressure, temperature, Xnacl, fsp);
127 
128  break;
129  }
130 
132  {
133  // Calculate the gas properties
135 
136  // Calculate the saturation
137  saturationTwoPhase(pressure, temperature, Xnacl, Z, fsp);
138 
139  // Calculate the liquid properties
140  Real liquid_pressure = pressure - _pc.capillaryPressure(1.0 - gas.saturation, qp);
141  liquidProperties(liquid_pressure, temperature, Xnacl, fsp);
142 
143  break;
144  }
145  }
146 
147  // Liquid saturations can now be set
148  liquid.saturation = 1.0 - gas.saturation;
149  liquid.dsaturation_dp = -gas.dsaturation_dp;
150  liquid.dsaturation_dT = -gas.dsaturation_dT;
151  liquid.dsaturation_dX = -gas.dsaturation_dX;
152  liquid.dsaturation_dZ = -gas.dsaturation_dZ;
153 
154  // Save pressures to FluidStateProperties object
155  gas.pressure = pressure;
156  liquid.pressure = pressure - _pc.capillaryPressure(liquid.saturation, qp);
157 }
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
FluidStatePhaseEnum
Phase state enum.
unsigned int _gas_phase_number
Phase number of the gas phase.
void clearFluidStateProperties(std::vector< FluidStateProperties > &fsp) const
Clears the contents of the FluidStateProperties data structure.
const std::string temperature
Definition: NS.h:27
virtual Real capillaryPressure(Real saturation, unsigned qp=0) const
Capillary pressure is calculated as a function of true saturation.
Data structure to pass calculated thermophysical properties.
void gasProperties(Real pressure, Real temperature, std::vector< FluidStateProperties > &fsp) const
Thermophysical properties of the gaseous state.
const PorousFlowCapillaryPressure & _pc
Capillary pressure UserObject.
const std::string pressure
Definition: NS.h:26
void liquidProperties(Real pressure, Real temperature, Real Xnacl, std::vector< FluidStateProperties > &fsp) const
Thermophysical properties of the liquid state.
void checkVariables(Real pressure, Real temperature) const
Check the input variables.
void saturationTwoPhase(Real pressure, Real temperature, Real Xnacl, Real Z, std::vector< FluidStateProperties > &fsp) const
Gas and liquid saturations for the two-phase region.
void massFractions(Real pressure, Real temperature, Real Xnacl, Real Z, FluidStatePhaseEnum &phase_state, std::vector< FluidStateProperties > &fsp) const
Mass fractions of CO2 and H2O in both phases, as well as derivatives wrt PorousFlow variables...

◆ totalMassFraction()

Real PorousFlowBrineCO2::totalMassFraction ( Real  pressure,
Real  temperature,
Real  Xnacl,
Real  saturation,
unsigned int  qp 
) const
overridevirtual

Total mass fraction of fluid component summed over all phases in the two-phase state for a specified gas saturation.

Parameters
pressuregas pressure (Pa)
temperaturetemperature (K)
XnaclNaCl mass fraction (kg/kg)
saturationgas saturation (-)
qpquadpoint index
Returns
total mass fraction Z (-)

Implements PorousFlowFluidStateBase.

Definition at line 1419 of file PorousFlowBrineCO2.C.

1421 {
1422  // Check whether the input pressure and temperature are within the region of validity
1424 
1425  // FluidStateProperties data structure
1426  std::vector<FluidStateProperties> fsp(_num_phases, FluidStateProperties(_num_components));
1429 
1430  // Calculate equilibrium mass fractions in the two-phase state
1431  Real Xco2, dXco2_dp, dXco2_dT, dXco2_dX, Yh2o, dYh2o_dp, dYh2o_dT, dYh2o_dX;
1433  temperature,
1434  Xnacl,
1435  Xco2,
1436  dXco2_dp,
1437  dXco2_dT,
1438  dXco2_dX,
1439  Yh2o,
1440  dYh2o_dp,
1441  dYh2o_dT,
1442  dYh2o_dX);
1443 
1444  // Save the mass fractions in the FluidStateMassFractions object
1445  const Real Yco2 = 1.0 - Yh2o;
1446  liquid.mass_fraction[_aqueous_fluid_component] = 1.0 - Xco2;
1447  liquid.mass_fraction[_gas_fluid_component] = Xco2;
1449  gas.mass_fraction[_gas_fluid_component] = Yco2;
1450 
1451  // Gas properties
1453 
1454  // Liquid properties
1455  const Real liquid_saturation = 1.0 - saturation;
1456  const Real liquid_pressure = pressure - _pc.capillaryPressure(liquid_saturation, qp);
1457  liquidProperties(liquid_pressure, temperature, Xnacl, fsp);
1458 
1459  // The total mass fraction of ncg (z) can now be calculated
1460  const Real Z = (saturation * gas.density * Yco2 + liquid_saturation * liquid.density * Xco2) /
1461  (saturation * gas.density + liquid_saturation * liquid.density);
1462 
1463  return Z;
1464 }
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
unsigned int _num_phases
Number of phases.
unsigned int _gas_phase_number
Phase number of the gas phase.
unsigned int _num_components
Number of components.
const unsigned int _aqueous_fluid_component
Fluid component number of the aqueous component.
const std::string temperature
Definition: NS.h:27
virtual Real capillaryPressure(Real saturation, unsigned qp=0) const
Capillary pressure is calculated as a function of true saturation.
std::vector< Real > mass_fraction
Data structure to pass calculated thermophysical properties.
unsigned int _gas_fluid_component
Fluid component number of the gas phase.
void equilibriumMassFractions(Real pressure, Real temperature, Real Xnacl, Real &Xco2, Real &dXco2_dp, Real &dXco2_dT, Real &dXco2_dX, Real &Yh2o, Real &dYh2o_dp, Real &dYh2o_dT, Real &dYh2o_dX) const
Mass fractions of CO2 in brine and water vapor in CO2 at equilibrium.
void gasProperties(Real pressure, Real temperature, std::vector< FluidStateProperties > &fsp) const
Thermophysical properties of the gaseous state.
const PorousFlowCapillaryPressure & _pc
Capillary pressure UserObject.
const std::string pressure
Definition: NS.h:26
void liquidProperties(Real pressure, Real temperature, Real Xnacl, std::vector< FluidStateProperties > &fsp) const
Thermophysical properties of the liquid state.
void checkVariables(Real pressure, Real temperature) const
Check the input variables.

◆ vaporMassFraction() [1/2]

Real PorousFlowFluidStateFlash::vaporMassFraction ( Real  Z0,
Real  K0,
Real  K1 
) const
inherited

Solves Rachford-Rice equation to provide vapor mass fraction.

For two components, the analytical solution is used, while for cases with more than two components, a Newton-Raphson iterative solution is calculated.

Parameters
Zitotal mass fraction(s)
Kiequilibrium constant(s)
Returns
vapor mass fraction

Definition at line 108 of file PorousFlowFluidStateFlash.C.

Referenced by PorousFlowWaterNCG::saturationTwoPhase(), saturationTwoPhase(), and PorousFlowFluidStateFlash::vaporMassFraction().

109 {
110  return (Z0 * (K1 - K0) - (K1 - 1.0)) / ((K0 - 1.0) * (K1 - 1.0));
111 }

◆ vaporMassFraction() [2/2]

Real PorousFlowFluidStateFlash::vaporMassFraction ( std::vector< Real > &  Zi,
std::vector< Real > &  Ki 
) const
inherited

Definition at line 114 of file PorousFlowFluidStateFlash.C.

115 {
116  // Check that the sizes of the mass fractions and equilibrium constant vectors are correct
117  if (Ki.size() != Zi.size() + 1)
118  mooseError("The number of mass fractions or equilibrium components passed to rachfordRice is "
119  "not correct");
120  Real v;
121 
122  // If there are only two components, an analytical solution is possible
123  if (Ki.size() == 2)
124  v = vaporMassFraction(Zi[0], Ki[0], Ki[1]);
125  else
126  {
127  // More than two components - solve the Rachford-Rice equation using
128  // Newton-Raphson method
129  // Initial guess for vapor mass fraction
130  Real v0 = 0.5;
131  unsigned int iter = 0;
132 
133  while (std::abs(rachfordRice(v0, Zi, Ki)) > _nr_tol)
134  {
135  v0 = v0 - rachfordRice(v0, Zi, Ki) / rachfordRiceDeriv(v0, Zi, Ki);
136  iter++;
137 
138  if (iter > _nr_max_its)
139  break;
140  }
141  v = v0;
142  }
143  return v;
144 }
Real rachfordRice(Real vf, std::vector< Real > &Zi, std::vector< Real > &Ki) const
Rachford-Rice equation for vapor fraction.
Real vaporMassFraction(Real Z0, Real K0, Real K1) const
Solves Rachford-Rice equation to provide vapor mass fraction.
const Real _nr_tol
Tolerance for Newton-Raphson iterations.
Real rachfordRiceDeriv(Real vf, std::vector< Real > &Zi, std::vector< Real > &Ki) const
Derivative of Rachford-Rice equation wrt vapor fraction.
const Real _nr_max_its
Maximum number of iterations for the Newton-Raphson routine.

Member Data Documentation

◆ _aqueous_fluid_component

const unsigned int PorousFlowFluidStateBase::_aqueous_fluid_component
protectedinherited

◆ _aqueous_phase_number

const unsigned int PorousFlowFluidStateBase::_aqueous_phase_number
protectedinherited

◆ _brine_fp

const BrineFluidProperties& PorousFlowBrineCO2::_brine_fp
protected

Fluid properties UserObject for water.

Definition at line 522 of file PorousFlowBrineCO2.h.

Referenced by liquidProperties(), PorousFlowBrineCO2(), saturationTwoPhase(), and solveEquilibriumMoleFractionHighTemp().

◆ _co2_fp

const SinglePhaseFluidProperties& PorousFlowBrineCO2::_co2_fp
protected

◆ _gas_fluid_component

unsigned int PorousFlowFluidStateBase::_gas_fluid_component
protectedinherited

◆ _gas_phase_number

unsigned int PorousFlowFluidStateBase::_gas_phase_number
protectedinherited

◆ _invMh2o

const Real PorousFlowBrineCO2::_invMh2o
protected

◆ _Mco2

const Real PorousFlowBrineCO2::_Mco2
protected

◆ _Mh2o

const Real PorousFlowBrineCO2::_Mh2o
protected

Molar mass of water (kg/mol)

Definition at line 528 of file PorousFlowBrineCO2.h.

Referenced by equilibriumMassFractions().

◆ _Mnacl

const Real PorousFlowBrineCO2::_Mnacl
protected

◆ _nr_max_its

const Real PorousFlowFluidStateBase::_nr_max_its
protectedinherited

Maximum number of iterations for the Newton-Raphson iterations.

Definition at line 194 of file PorousFlowFluidStateBase.h.

◆ _nr_tol

const Real PorousFlowFluidStateBase::_nr_tol
protectedinherited

Tolerance for Newton-Raphson iterations.

Definition at line 196 of file PorousFlowFluidStateBase.h.

◆ _num_components

unsigned int PorousFlowFluidStateBase::_num_components
protectedinherited

◆ _num_phases

unsigned int PorousFlowFluidStateBase::_num_phases
protectedinherited

◆ _pc

const PorousFlowCapillaryPressure& PorousFlowFluidStateBase::_pc
protectedinherited

◆ _R

const Real PorousFlowFluidStateBase::_R
protectedinherited

Universal gas constant (J/mol/K)

Definition at line 190 of file PorousFlowFluidStateBase.h.

Referenced by PorousFlowWaterNCG::enthalpyOfDissolution(), and enthalpyOfDissolutionGas().

◆ _Rbar

const Real PorousFlowBrineCO2::_Rbar
protected

Molar gas constant in bar cm^3 /(K mol)

Definition at line 536 of file PorousFlowBrineCO2.h.

Referenced by fugacityCoefficientCO2HighTemp(), fugacityCoefficientH2OHighTemp(), fugacityCoefficientsLowTemp(), funcABHighTemp(), and funcABLowTemp().

◆ _salt_component

const unsigned int PorousFlowBrineCO2::_salt_component
protected

Salt component index.

Definition at line 520 of file PorousFlowBrineCO2.h.

Referenced by massFractions(), PorousFlowBrineCO2(), and saltComponentIndex().

◆ _T_c2k

const Real PorousFlowFluidStateBase::_T_c2k
protectedinherited

◆ _Tlower

const Real PorousFlowBrineCO2::_Tlower
protected

Temperature below which the Spycher, Pruess & Ennis-King (2003) model is used (K)

Definition at line 538 of file PorousFlowBrineCO2.h.

Referenced by equilibriumMoleFractions(), and smoothCubicInterpolation().

◆ _Tupper

const Real PorousFlowBrineCO2::_Tupper
protected

Temperature above which the Spycher & Pruess (2010) model is used (K)

Definition at line 540 of file PorousFlowBrineCO2.h.

Referenced by equilibriumMoleFractions(), and smoothCubicInterpolation().

◆ _water_fp

const SinglePhaseFluidProperties& PorousFlowBrineCO2::_water_fp
protected

Fluid properties UserObject for H20.

Definition at line 526 of file PorousFlowBrineCO2.h.

◆ _Zmin

const Real PorousFlowBrineCO2::_Zmin
protected

Minimum Z - below this value all CO2 will be dissolved.

This reduces the computational burden when small values of Z are present

Definition at line 543 of file PorousFlowBrineCO2.h.

Referenced by massFractions().


The documentation for this class was generated from the following files: