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...
 
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 {
44  // Check that the correct FluidProperties UserObjects have been provided
45  if (_co2_fp.fluidName() != "co2")
46  paramError("co2_fp", "A valid CO2 FluidProperties UserObject must be provided");
47 
48  if (_brine_fp.fluidName() != "brine")
49  paramError("brine_fp", "A valid Brine FluidProperties UserObject must be provided");
50 
51  // Set the number of phases and components, and their indexes
52  _num_phases = 2;
53  _num_components = 3;
56 
57  // Check that _aqueous_phase_number is <= total number of phases
59  paramError("liquid_phase_number",
60  "This value is larger than the possible number of phases ",
61  _num_phases);
62 
63  // Check that _aqueous_fluid_component is <= total number of fluid components
65  paramError("liquid_fluid_component",
66  "This value is larger than the possible number of fluid components",
68 
69  // Check that the salt component index is not identical to the liquid fluid component
71  paramError(
72  "salt_component",
73  "The value provided must be different from the value entered in liquid_fluid_component");
74 
75  // Check that _salt_component is <= total number of fluid components
77  paramError("salt_component",
78  "The value provided is larger than the possible number of fluid components",
80 }
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.
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 835 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractionsLowTemp().

842 {
843  // Need pressure in bar
844  const Real pbar = pressure * 1.0e-5;
845  // Need NaCl molality (mol/kg)
846  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
847  const Real dmnacl_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
848 
849  const Real lambda = -0.411370585 + 6.07632013e-4 * temperature + 97.5347708 / temperature -
850  0.0237622469 * pbar / temperature +
851  0.0170656236 * pbar / (630.0 - temperature) +
852  1.41335834e-5 * temperature * std::log(pbar);
853 
854  const Real xi = 3.36389723e-4 - 1.9829898e-5 * temperature + 2.12220830e-3 * pbar / temperature -
855  5.24873303e-3 * pbar / (630.0 - temperature);
856 
857  gamma = std::exp(2.0 * lambda * mnacl + xi * mnacl * mnacl);
858 
859  // Derivative wrt pressure
860  const Real dlambda_dp = -0.0237622469 / temperature + 0.0170656236 / (630.0 - temperature) +
861  1.41335834e-5 * temperature / pbar;
862  const Real dxi_dp = 2.12220830e-3 / temperature - 5.24873303e-3 / (630.0 - temperature);
863  dgamma_dp = (2.0 * mnacl * dlambda_dp + mnacl * mnacl * dxi_dp) * gamma * 1.0e-5;
864 
865  // Derivative wrt temperature
866  const Real dlambda_dT = 6.07632013e-4 - 97.5347708 / temperature / temperature +
867  0.0237622469 * pbar / temperature / temperature +
868  0.0170656236 * pbar / (630.0 - temperature) / (630.0 - temperature) +
869  1.41335834e-5 * std::log(pbar);
870  const Real dxi_dT = -1.9829898e-5 - 2.12220830e-3 * pbar / temperature / temperature -
871  5.24873303e-3 * pbar / (630.0 - temperature) / (630.0 - temperature);
872  dgamma_dT = (2.0 * mnacl * dlambda_dT + mnacl * mnacl * dxi_dT) * gamma;
873 
874  // Derivative wrt salt mass fraction
875  dgamma_dX = (2.0 * lambda * dmnacl_dX + 2.0 * xi * mnacl * dmnacl_dX) * gamma;
876 }
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 820 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp().

821 {
822  if (temperature <= 373.15)
823  return 1.0;
824  else
825  {
826  const Real Tref = temperature - 373.15;
827  const Real xh2o = 1.0 - xco2;
828  const Real Am = -3.084e-2 * Tref + 1.927e-5 * Tref * Tref;
829 
830  return std::exp(2.0 * Am * xco2 * xh2o * xh2o);
831  }
832 }
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 805 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp().

806 {
807  if (temperature <= 373.15)
808  return 1.0;
809  else
810  {
811  const Real Tref = temperature - 373.15;
812  const Real xh2o = 1.0 - xco2;
813  const Real Am = -3.084e-2 * Tref + 1.927e-5 * Tref * Tref;
814 
815  return std::exp((Am - 2.0 * Am * xh2o) * xco2 * xco2);
816  }
817 }
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 879 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp().

881 {
882  // Need NaCl molality (mol/kg)
883  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
884  const Real dmnacl_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
885 
886  const Real T2 = temperature * temperature;
887  const Real T3 = temperature * T2;
888 
889  const Real lambda = 2.217e-4 * temperature + 1.074 / temperature + 2648.0 / T2;
890  const Real xi = 1.3e-5 * temperature - 20.12 / temperature + 5259.0 / T2;
891 
892  gamma = (1.0 - mnacl / _invMh2o) * std::exp(2.0 * lambda * mnacl + xi * mnacl * mnacl);
893 
894  // Derivative wrt temperature
895  const Real dlambda_dT = 2.217e-4 - 1.074 / T2 - 5296.0 / T3;
896  const Real dxi_dT = 1.3e-5 + 20.12 / T2 - 10518.0 / T3;
897  dgamma_dT = (2.0 * mnacl * dlambda_dT + mnacl * mnacl * dxi_dT) * gamma;
898 
899  // Derivative wrt salt mass fraction
900  dgamma_dX = (2.0 * lambda * dmnacl_dX + 2.0 * xi * mnacl * dmnacl_dX) * gamma -
901  dmnacl_dX / _invMh2o * std::exp(2.0 * lambda * mnacl + xi * mnacl * mnacl);
902 }
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 1541 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties(), and totalMassFraction().

1542 {
1543  // The calculation of mass fractions is valid from 12C <= T <= 300C, and
1544  // pressure less than 60 MPa
1545  if (temperature < 285.15 || temperature > 573.15)
1546  mooseException(name() + ": temperature " + Moose::stringify(temperature) +
1547  " is outside range 285.15 K <= T <= 573.15 K");
1548 
1549  if (pressure > 6.0e7)
1550  mooseException(name() + ": pressure " + Moose::stringify(pressure) +
1551  " must be less than 60 MPa");
1552 }
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 1511 of file PorousFlowBrineCO2.C.

Referenced by liquidProperties().

1512 {
1513  // Linear fit to model of Duan and Sun (2003) (in kJ/mol)
1514  const Real delta_h = -58.3533 + 0.134519 * temperature;
1515 
1516  // Convert to J/kg
1517  hdis = delta_h * 1000.0 / _Mco2;
1518  dhdis_dT = 134.519 / _Mco2;
1519 }
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 1486 of file PorousFlowBrineCO2.C.

1488 {
1489  // Henry's constant
1490  Real Kh, dKh_dT, dKh_dX;
1491  henryConstant(temperature, Xnacl, Kh, dKh_dT, dKh_dX);
1492 
1493  hdis = -_R * temperature * temperature * dKh_dT / Kh / _Mco2;
1494 
1495  // Derivative of enthalpy of dissolution wrt temperature and xnacl requires the second
1496  // derivatives of Henry's constant. For simplicity, approximate these numerically
1497  const Real dT = temperature * 1.0e-8;
1498  const Real T2 = temperature + dT;
1499  henryConstant(T2, Xnacl, Kh, dKh_dT, dKh_dX);
1500 
1501  dhdis_dT = (-_R * T2 * T2 * dKh_dT / Kh / _Mco2 - hdis) / dT;
1502 
1503  const Real dX = Xnacl * 1.0e-8;
1504  const Real X2 = Xnacl + dX;
1505  henryConstant(temperature, X2, Kh, dKh_dT, dKh_dX);
1506 
1507  dhdis_dX = (-_R * temperature * temperature * dKh_dT / Kh / _Mco2 - hdis) / dX;
1508 }
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 938 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp(), and funcABLowTemp().

939 {
940  // Uses temperature in Celcius
941  const Real Tc = temperature - _T_c2k;
942  const Real Tc2 = Tc * Tc;
943  const Real Tc3 = Tc2 * Tc;
944 
945  Real logK0CO2, dlogK0CO2;
946 
947  if (Tc <= 99.0)
948  {
949  logK0CO2 = 1.189 + 1.304e-2 * Tc - 5.446e-5 * Tc2;
950  dlogK0CO2 = 1.304e-2 - 1.0892e-4 * Tc;
951  }
952  else if (Tc > 99.0 && Tc < 109.0)
953  {
954  const Real Tint = (Tc - 99.0) / 10.0;
955  const Real Tint2 = Tint * Tint;
956  logK0CO2 = 1.9462 + 2.25692e-2 * Tint - 9.49577e-3 * Tint2 - 6.77721e-3 * Tint * Tint2;
957  dlogK0CO2 = 2.25692e-2 - 1.899154e-2 * Tint + 2.033163e-2 * Tint2;
958  }
959  else // 109 <= Tc <= 300
960  {
961  logK0CO2 = 1.668 + 3.992e-3 * Tc - 1.156e-5 * Tc2 + 1.593e-9 * Tc3;
962  dlogK0CO2 = 3.992e-3 - 2.312e-5 * Tc + 4.779e-9 * Tc2;
963  }
964 
965  kco2 = std::pow(10.0, logK0CO2);
966  dkco2_dT = std::log(10.0) * dlogK0CO2 * kco2;
967 }
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 905 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp(), and funcABLowTemp().

906 {
907  // Uses temperature in Celcius
908  const Real Tc = temperature - _T_c2k;
909  const Real Tc2 = Tc * Tc;
910  const Real Tc3 = Tc2 * Tc;
911  const Real Tc4 = Tc3 * Tc;
912 
913  Real logK0H2O, dlogK0H2O;
914 
915  if (Tc <= 99.0)
916  {
917  logK0H2O = -2.209 + 3.097e-2 * Tc - 1.098e-4 * Tc2 + 2.048e-7 * Tc3;
918  dlogK0H2O = 3.097e-2 - 2.196e-4 * Tc + 6.144e-7 * Tc2;
919  }
920  else if (Tc > 99.0 && Tc < 109.0)
921  {
922  const Real Tint = (Tc - 99.0) / 10.0;
923  const Real Tint2 = Tint * Tint;
924  logK0H2O = -0.0204026 + 0.0152513 * Tint + 0.417565 * Tint2 - 0.278636 * Tint * Tint2;
925  dlogK0H2O = 0.0152513 + 0.83513 * Tint - 0.835908 * Tint2;
926  }
927  else // 109 <= Tc <= 300
928  {
929  logK0H2O = -2.1077 + 2.8127e-2 * Tc - 8.4298e-5 * Tc2 + 1.4969e-7 * Tc3 - 1.1812e-10 * Tc4;
930  dlogK0H2O = 2.8127e-2 - 1.68596e-4 * Tc + 4.4907e-7 * Tc2 - 4.7248e-10 * Tc3;
931  }
932 
933  kh2o = std::pow(10.0, logK0H2O);
934  dkh2o_dT = std::log(10.0) * dlogK0H2O * kh2o;
935 }
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 524 of file PorousFlowBrineCO2.C.

Referenced by massFractions(), and totalMassFraction().

535 {
536  Real co2_density, dco2_density_dp, dco2_density_dT;
537  _co2_fp.rho_from_p_T(pressure, temperature, co2_density, dco2_density_dp, dco2_density_dT);
538 
539  // Mole fractions at equilibrium
540  Real xco2, dxco2_dp, dxco2_dT, dxco2_dX, yh2o, dyh2o_dp, dyh2o_dT, dyh2o_dX;
542  temperature,
543  Xnacl,
544  xco2,
545  dxco2_dp,
546  dxco2_dT,
547  dxco2_dX,
548  yh2o,
549  dyh2o_dp,
550  dyh2o_dT,
551  dyh2o_dX);
552 
553  // The mass fraction of H2O in gas (assume no salt in gas phase) and derivatives
554  // wrt p, T, and X
555  Yh2o = yh2o * _Mh2o / (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2);
556  dYh2o_dp = _Mco2 * _Mh2o * dyh2o_dp / (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2) /
557  (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2);
558  dYh2o_dT = _Mco2 * _Mh2o * dyh2o_dT / (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2) /
559  (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2);
560  dYh2o_dX = dyh2o_dX * _Mh2o * _Mco2 / (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2) /
561  (yh2o * _Mh2o + (1.0 - yh2o) * _Mco2);
562 
563  // NaCl molality (mol/kg)
564  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
565  const Real dmnacl_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
566 
567  // The molality of CO2 in 1kg of H2O
568  const Real mco2 = xco2 * (2.0 * mnacl + _invMh2o) / (1.0 - xco2);
569  // The mass fraction of CO2 in brine is then
570  const Real denominator = (1.0 + mnacl * _Mnacl + mco2 * _Mco2);
571  Xco2 = mco2 * _Mco2 / denominator;
572 
573  // Derivatives of Xco2 wrt p, T and X
574  const Real denominator2 = denominator * denominator;
575  const Real dmco2_dp = dxco2_dp * (2.0 * mnacl + _invMh2o) / (1.0 - xco2) / (1.0 - xco2);
576  dXco2_dp = (1.0 + mnacl * _Mnacl) * _Mco2 * dmco2_dp / denominator2;
577 
578  const Real dmco2_dT = dxco2_dT * (2.0 * mnacl + _invMh2o) / (1.0 - xco2) / (1.0 - xco2);
579  dXco2_dT = (1.0 + mnacl * _Mnacl) * _Mco2 * dmco2_dT / denominator2;
580 
581  const Real dmco2_dX =
582  (dxco2_dX * (2.0 * mnacl + _invMh2o) + 2.0 * (1.0 - xco2) * xco2 * dmnacl_dX) / (1.0 - xco2) /
583  (1.0 - xco2);
584  dXco2_dX = _Mco2 * ((1.0 + mnacl * _Mnacl) * dmco2_dX - dmnacl_dX * mco2 * _Mnacl) / denominator2;
585 }
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
virtual Real rho_from_p_T(Real p, Real T) const
Density from pressure and temperature.
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 970 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMassFractions().

981 {
982  // CO2 density and derivatives wrt pressure and temperature
983  Real co2_density, dco2_density_dp, dco2_density_dT;
984  _co2_fp.rho_from_p_T(pressure, temperature, co2_density, dco2_density_dp, dco2_density_dT);
985 
986  if (temperature <= _Tlower)
987  {
989  temperature,
990  Xnacl,
991  xco2,
992  dxco2_dp,
993  dxco2_dT,
994  dxco2_dX,
995  yh2o,
996  dyh2o_dp,
997  dyh2o_dT,
998  dyh2o_dX);
999  }
1000  else if (temperature > _Tlower && temperature < _Tupper)
1001  {
1002  // Cubic polynomial in this regime
1003  const Real Tint = (temperature - _Tlower) / 10.0;
1004 
1005  // Equilibrium mole fractions and derivatives at the lower temperature
1006  Real xco2_lower, dxco2_dT_lower, yh2o_lower, dyh2o_dT_lower;
1008  _Tlower,
1009  Xnacl,
1010  xco2_lower,
1011  dxco2_dp,
1012  dxco2_dT_lower,
1013  dxco2_dX,
1014  yh2o_lower,
1015  dyh2o_dp,
1016  dyh2o_dT_lower,
1017  dyh2o_dX);
1018 
1019  // Equilibrium mole fractions and derivatives at the upper temperature
1020  Real xco2_upper, yh2o_upper;
1021  Real co2_density_upper = _co2_fp.rho_from_p_T(pressure, _Tupper);
1022 
1024  pressure, _Tupper, Xnacl, co2_density_upper, xco2_upper, yh2o_upper);
1025 
1026  Real A, dA_dp, dA_dT, B, dB_dp, dB_dT, dB_dX;
1028  _Tupper,
1029  Xnacl,
1030  co2_density,
1031  xco2_upper,
1032  yh2o_upper,
1033  A,
1034  dA_dp,
1035  dA_dT,
1036  B,
1037  dB_dp,
1038  dB_dT,
1039  dB_dX);
1040 
1041  const Real dyh2o_dT_upper =
1042  ((1.0 - B) * dA_dT + (A - 1.0) * A * dB_dT) / (1.0 - A * B) / (1.0 - A * B);
1043  const Real dxco2_dT_upper = dB_dT * (1.0 - yh2o_upper) - B * dyh2o_dT_upper;
1044 
1045  // The mole fractions in this regime are then found by interpolation
1046  Real dxco2_dT, dyh2o_dT;
1048  Tint, xco2_lower, dxco2_dT_lower, xco2_upper, dxco2_dT_upper, xco2, dxco2_dT);
1050  Tint, yh2o_lower, dyh2o_dT_lower, yh2o_upper, dyh2o_dT_upper, yh2o, dyh2o_dT);
1051  }
1052  else
1053  {
1054  // Equilibrium mole fractions solved using iteration in this regime
1055  solveEquilibriumMoleFractionHighTemp(pressure, temperature, Xnacl, co2_density, xco2, yh2o);
1056 
1057  // Can use these in funcABHighTemp() to get derivatives analyticall rather than by iteration
1058  Real A, dA_dp, dA_dT, B, dB_dp, dB_dT, dB_dX;
1060  temperature,
1061  Xnacl,
1062  co2_density,
1063  xco2,
1064  yh2o,
1065  A,
1066  dA_dp,
1067  dA_dT,
1068  B,
1069  dB_dp,
1070  dB_dT,
1071  dB_dX);
1072 
1073  dyh2o_dp = ((1.0 - B) * dA_dp + (A - 1.0) * A * dB_dp) / (1.0 - A * B) / (1.0 - A * B);
1074  dxco2_dp = dB_dp * (1.0 - yh2o) - B * dyh2o_dp;
1075 
1076  dyh2o_dT = ((1.0 - B) * dA_dT + (A - 1.0) * A * dB_dT) / (1.0 - A * B) / (1.0 - A * B);
1077  dxco2_dT = dB_dT * (1.0 - yh2o) - B * dyh2o_dT;
1078 
1079  dyh2o_dX = ((A - 1.0) * A * dB_dX) / (1.0 - A * B) / (1.0 - A * B);
1080  dxco2_dX = dB_dX * (1.0 - yh2o) - B * dyh2o_dX;
1081  }
1082 }
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) ...
virtual Real rho_from_p_T(Real p, Real T) const
Density from pressure and temperature.
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 1085 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractions().

1096 {
1097  if (temperature > 373.15)
1098  mooseError(name(),
1099  ": equilibriumMoleFractionsLowTemp() is not valid for T > 373.15K. Use "
1100  "equilibriumMoleFractions() instead");
1101 
1102  // CO2 density and derivatives wrt pressure and temperature
1103  Real co2_density, dco2_density_dp, dco2_density_dT;
1104  _co2_fp.rho_from_p_T(pressure, temperature, co2_density, dco2_density_dp, dco2_density_dT);
1105 
1106  // Assume infinite dilution (yh20 = 0 and xco2 = 0) in low temperature regime
1107  Real A, dA_dp, dA_dT, B, dB_dp, dB_dT;
1109  temperature,
1110  co2_density,
1111  dco2_density_dp,
1112  dco2_density_dT,
1113  A,
1114  dA_dp,
1115  dA_dT,
1116  B,
1117  dB_dp,
1118  dB_dT);
1119 
1120  // As the activity coefficient for CO2 in brine used in this regime isn't a 'true'
1121  // activity coefficient, we instead calculate the molality of CO2 in water, then
1122  // correct it for brine, and then calculate the mole fractions.
1123  // The mole fraction in pure water is
1124  const Real yh2ow = (1.0 - B) / (1.0 / A - B);
1125  const Real xco2w = B * (1.0 - yh2ow);
1126 
1127  // Molality of CO2 in pure water
1128  const Real mco2w = xco2w * _invMh2o / (1.0 - xco2w);
1129  // Molality of CO2 in brine is then calculated using gamma
1130  Real gamma, dgamma_dp, dgamma_dT, dgamma_dX;
1131  activityCoefficient(pressure, temperature, Xnacl, gamma, dgamma_dp, dgamma_dT, dgamma_dX);
1132  const Real mco2 = mco2w / gamma;
1133 
1134  // Need NaCl molality (mol/kg)
1135  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
1136  const Real dmnacl_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
1137 
1138  // Mole fractions of CO2 and H2O in liquid and gas phases
1139  const Real total_moles = 2.0 * mnacl + _invMh2o + mco2;
1140  xco2 = mco2 / total_moles;
1141  yh2o = A * (1.0 - xco2 - 2.0 * mnacl / total_moles);
1142 
1143  // The derivatives of the mole fractions wrt pressure
1144  const Real dyh2ow_dp =
1145  ((1.0 - B) * dA_dp + (A - 1.0) * A * dB_dp) / (1.0 - A * B) / (1.0 - A * B);
1146  const Real dxco2w_dp = dB_dp * (1.0 - yh2ow) - B * dyh2ow_dp;
1147 
1148  const Real dmco2w_dp = _invMh2o * dxco2w_dp / (1.0 - xco2w) / (1.0 - xco2w);
1149  const Real dmco2_dp = dmco2w_dp / gamma - mco2w * dgamma_dp / gamma / gamma;
1150  dxco2_dp = (2.0 * mnacl + _invMh2o) * dmco2_dp / total_moles / total_moles;
1151 
1152  dyh2o_dp = (1.0 - xco2 - 2.0 * mnacl / total_moles) * dA_dp - A * dxco2_dp +
1153  2.0 * A * mnacl * dmco2_dp / total_moles / total_moles;
1154 
1155  // The derivatives of the mole fractions wrt temperature
1156  const Real dyh2ow_dT =
1157  ((1.0 - B) * dA_dT + (A - 1.0) * A * dB_dT) / (1.0 - A * B) / (1.0 - A * B);
1158  const Real dxco2w_dT = dB_dT * (1.0 - yh2ow) - B * dyh2ow_dT;
1159 
1160  const Real dmco2w_dT = _invMh2o * dxco2w_dT / (1.0 - xco2w) / (1.0 - xco2w);
1161  const Real dmco2_dT = dmco2w_dT / gamma - mco2w * dgamma_dT / gamma / gamma;
1162  dxco2_dT = (2.0 * mnacl + _invMh2o) * dmco2_dT / total_moles / total_moles;
1163 
1164  dyh2o_dT = (1.0 - xco2 - 2.0 * mnacl / total_moles) * dA_dT - A * dxco2_dT +
1165  2.0 * A * mnacl * dmco2_dT / total_moles / total_moles;
1166 
1167  // The derivatives of the mole fractions wrt salt mass fraction
1168  const Real dmco2_dX = -mco2w * dgamma_dX / gamma / gamma;
1169  dxco2_dX =
1170  ((2.0 * mnacl + _invMh2o) * dmco2_dX - 2.0 * mco2 * dmnacl_dX) / total_moles / total_moles;
1171  dyh2o_dX =
1172  A * (2.0 * (mnacl * dmco2_dX - (mco2 + _invMh2o) * dmnacl_dX) / total_moles / total_moles -
1173  dxco2_dX);
1174 }
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.
virtual Real rho_from_p_T(Real p, Real T) const
Density from pressure and temperature.
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 83 of file PorousFlowBrineCO2.C.

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

◆ fugacityCoefficientCO2HighTemp()

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

Definition at line 762 of file PorousFlowBrineCO2.C.

Referenced by fugacityCoefficientsHighTemp().

764 {
765  // Need pressure in bar
766  const Real pbar = pressure * 1.0e-5;
767  // Molar volume in cm^3/mol
768  const Real V = _Mco2 / co2_density * 1.0e6;
769 
770  // Redlich-Kwong parameters
771  const Real yco2 = 1.0 - yh2o;
772  const Real xh2o = 1.0 - xco2;
773 
774  const Real aCO2 = 8.008e7 - 4.984e4 * temperature;
775  const Real aH2O = 1.337e8 - 1.4e4 * temperature;
776  const Real bCO2 = 28.25;
777  const Real bH2O = 15.7;
778  const Real KH2OCO2 = 1.427e-2 - 4.037e-4 * temperature;
779  const Real KCO2H2O = 0.4228 - 7.422e-4 * temperature;
780  const Real kH2OCO2 = KH2OCO2 * yh2o + KCO2H2O * yco2;
781  const Real kCO2H2O = KCO2H2O * yh2o + KH2OCO2 * yco2;
782 
783  const Real aH2OCO2 = std::sqrt(aCO2 * aH2O) * (1.0 - kH2OCO2);
784  const Real aCO2H2O = std::sqrt(aCO2 * aH2O) * (1.0 - kCO2H2O);
785 
786  const Real amix = yh2o * yh2o * aH2O + yh2o * yco2 * (aH2OCO2 + aCO2H2O) + yco2 * yco2 * aCO2;
787  const Real bmix = yh2o * bH2O + yco2 * bCO2;
788 
789  const Real t15 = std::pow(temperature, 1.5);
790 
791  Real lnPhiCO2 = bCO2 / bmix * (pbar * V / (_Rbar * temperature) - 1.0) -
792  std::log(pbar * (V - bmix) / (_Rbar * temperature));
793 
794  Real term3 = (2.0 * yco2 * aCO2 + yh2o * (aH2OCO2 + aCO2H2O) -
795  yh2o * yco2 * std::sqrt(aH2O * aCO2) * (kH2OCO2 - kCO2H2O) * (yh2o - yco2) +
796  xh2o * xco2 * std::sqrt(aH2O * aCO2) * (kCO2H2O - kH2OCO2)) /
797  amix;
798 
799  lnPhiCO2 += (term3 - bCO2 / bmix) * amix / (bmix * _Rbar * t15) * std::log(V / (V + bmix));
800 
801  return std::exp(lnPhiCO2);
802 }
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 719 of file PorousFlowBrineCO2.C.

Referenced by fugacityCoefficientsHighTemp().

721 {
722  // Need pressure in bar
723  const Real pbar = pressure * 1.0e-5;
724  // Molar volume in cm^3/mol
725  const Real V = _Mco2 / co2_density * 1.0e6;
726 
727  // Redlich-Kwong parameters
728  const Real yco2 = 1.0 - yh2o;
729  const Real xh2o = 1.0 - xco2;
730 
731  const Real aCO2 = 8.008e7 - 4.984e4 * temperature;
732  const Real aH2O = 1.337e8 - 1.4e4 * temperature;
733  const Real bCO2 = 28.25;
734  const Real bH2O = 15.7;
735  const Real KH2OCO2 = 1.427e-2 - 4.037e-4 * temperature;
736  const Real KCO2H2O = 0.4228 - 7.422e-4 * temperature;
737  const Real kH2OCO2 = KH2OCO2 * yh2o + KCO2H2O * yco2;
738  const Real kCO2H2O = KCO2H2O * yh2o + KH2OCO2 * yco2;
739 
740  const Real aH2OCO2 = std::sqrt(aCO2 * aH2O) * (1.0 - kH2OCO2);
741  const Real aCO2H2O = std::sqrt(aCO2 * aH2O) * (1.0 - kCO2H2O);
742 
743  const Real amix = yh2o * yh2o * aH2O + yh2o * yco2 * (aH2OCO2 + aCO2H2O) + yco2 * yco2 * aCO2;
744  const Real bmix = yh2o * bH2O + yco2 * bCO2;
745 
746  const Real t15 = std::pow(temperature, 1.5);
747 
748  Real lnPhiH2O = bH2O / bmix * (pbar * V / (_Rbar * temperature) - 1.0) -
749  std::log(pbar * (V - bmix) / (_Rbar * temperature));
750  Real term3 = (2.0 * yh2o * aH2O + yco2 * (aH2OCO2 + aCO2H2O) -
751  yh2o * yco2 * std::sqrt(aH2O * aCO2) * (kH2OCO2 - kCO2H2O) * (yh2o - yco2) +
752  xh2o * xco2 * std::sqrt(aH2O * aCO2) * (kH2OCO2 - kCO2H2O)) /
753  amix;
754  term3 -= bH2O / bmix;
755  term3 *= amix / (bmix * _Rbar * t15) * std::log(V / (V + bmix));
756  lnPhiH2O += term3;
757 
758  return std::exp(lnPhiH2O);
759 }
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 663 of file PorousFlowBrineCO2.C.

Referenced by funcABHighTemp().

670 {
671  if (temperature <= 373.15)
672  mooseError(name(),
673  ": fugacityCoefficientsHighTemp() is not valid for T <= 373.15K. Use "
674  "fugacityCoefficientsLowTemp() instead");
675 
676  fh2o = fugacityCoefficientH2OHighTemp(pressure, temperature, co2_density, xco2, yh2o);
677  fco2 = fugacityCoefficientCO2HighTemp(pressure, temperature, co2_density, xco2, yh2o);
678 }
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 681 of file PorousFlowBrineCO2.C.

692 {
693  if (temperature <= 373.15)
694  mooseError(name(),
695  ": fugacityCoefficientsHighTemp() is not valid for T <= 373.15K. Use "
696  "fugacityCoefficientsLowTemp() instead");
697 
698  fh2o = fugacityCoefficientH2OHighTemp(pressure, temperature, co2_density, xco2, yh2o);
699  fco2 = fugacityCoefficientCO2HighTemp(pressure, temperature, co2_density, xco2, yh2o);
700 
701  // Derivatives by finite difference
702  const Real dp = 1.0e-2;
703  const Real dT = 1.0e-4;
704 
705  Real fh2o2 = fugacityCoefficientH2OHighTemp(pressure + dp, temperature, co2_density, xco2, yh2o);
706  Real fco22 = fugacityCoefficientCO2HighTemp(pressure + dp, temperature, xco2, co2_density, yh2o);
707 
708  dfh2o_dp = (fh2o2 - fh2o) / dp;
709  dfco2_dp = (fco22 - fco2) / dp;
710 
711  fh2o2 = fugacityCoefficientH2OHighTemp(pressure, temperature + dT, co2_density, xco2, yh2o);
712  fco22 = fugacityCoefficientCO2HighTemp(pressure, temperature + dT, co2_density, xco2, yh2o);
713 
714  dfh2o_dT = (fh2o2 - fh2o) / dT;
715  dfco2_dT = (fco22 - fco2) / dT;
716 }
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 588 of file PorousFlowBrineCO2.C.

Referenced by funcABLowTemp().

599 {
600  if (temperature > 373.15)
601  mooseError(name(),
602  ": fugacityCoefficientsLowTemp() is not valid for T > 373.15K. Use "
603  "fugacityCoefficientsHighTemp() instead");
604 
605  // Need pressure in bar
606  const Real pbar = pressure * 1.0e-5;
607 
608  // Molar volume in cm^3/mol
609  const Real V = _Mco2 / co2_density * 1.0e6;
610  const Real dV_dp = -V / co2_density * dco2_density_dp;
611  const Real dV_dT = -V / co2_density * dco2_density_dT;
612 
613  // Redlich-Kwong parameters
614  const Real aCO2 = 7.54e7 - 4.13e4 * temperature;
615  const Real daCO2_dT = -4.13e4;
616  const Real bCO2 = 27.8;
617  const Real aCO2H2O = 7.89e7;
618  const Real bH2O = 18.18;
619 
620  const Real t15 = std::pow(temperature, 1.5);
621 
622  // The fugacity coefficients for H2O and CO2
623  auto lnPhi = [V, aCO2, bCO2, t15, this](Real a, Real b) {
624  return std::log(V / (V - bCO2)) + b / (V - bCO2) -
625  2.0 * a / (_Rbar * t15 * bCO2) * std::log((V + bCO2) / V) +
626  aCO2 * b / (_Rbar * t15 * bCO2 * bCO2) * (std::log((V + bCO2) / V) - bCO2 / (V + bCO2));
627  };
628 
629  const Real lnPhiH2O = lnPhi(aCO2H2O, bH2O) - std::log(pbar * V / (_Rbar * temperature));
630  const Real lnPhiCO2 = lnPhi(aCO2, bCO2) - std::log(pbar * V / (_Rbar * temperature));
631 
632  fh2o = std::exp(lnPhiH2O);
633  fco2 = std::exp(lnPhiCO2);
634 
635  // The derivative of the fugacity coefficients wrt pressure
636  auto dlnPhi_dV = [V, aCO2, bCO2, t15, this](Real a, Real b) {
637  return (bCO2 * bCO2 - (bCO2 + b) * V) / V / (V - bCO2) / (V - bCO2) +
638  (2.0 * a * (V + bCO2) - aCO2 * b) / (_Rbar * t15 * V * (V + bCO2) * (V + bCO2));
639  };
640 
641  dfh2o_dp = (dlnPhi_dV(aCO2H2O, bH2O) * dV_dp - 1.0e-5 / pbar - dV_dp / V) * fh2o;
642  dfco2_dp = (dlnPhi_dV(aCO2, bCO2) * dV_dp - 1.0e-5 / pbar - dV_dp / V) * fco2;
643 
644  // The derivative of the fugacity coefficient wrt temperature
645  auto dlnPhi_dT = [V, aCO2, daCO2_dT, bCO2, t15, temperature, this](Real a, Real b, Real da_dT) {
646  return (3.0 * bCO2 * a * std::log((V + bCO2) / V) +
647  1.5 * aCO2 * b * (bCO2 / (V + bCO2) - std::log((V + bCO2) / V)) -
648  (2.0 * da_dT * bCO2 * std::log((V + bCO2) / V) -
649  daCO2_dT * b * (std::log((V + bCO2) / V) - bCO2 / (V + bCO2))) *
650  temperature) /
651  (temperature * t15 * _Rbar * bCO2 * bCO2);
652  };
653 
654  dfh2o_dT = (dlnPhi_dT(aCO2H2O, bH2O, 0) + dlnPhi_dV(aCO2H2O, bH2O) * dV_dT - dV_dT / V +
655  1.0 / temperature) *
656  fh2o;
657  dfco2_dT = (dlnPhi_dT(aCO2, bCO2, daCO2_dT) + dlnPhi_dV(aCO2, bCO2) * dV_dT - dV_dT / V +
658  1.0 / temperature) *
659  fco2;
660 }
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 1242 of file PorousFlowBrineCO2.C.

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

1250 {
1251  if (temperature <= 373.15)
1252  mooseError(name(),
1253  ": funcABHighTemp() is not valid for T <= 373.15K. Use funcABLowTemp() instead");
1254 
1255  // Pressure in bar
1256  const Real pbar = pressure * 1.0e-5;
1257  // Temperature in C
1258  const Real Tc = temperature - _T_c2k;
1259 
1260  // Reference pressure and partial molar volumes
1261  const Real pref = -1.9906e-1 + 2.0471e-3 * Tc + 1.0152e-4 * Tc * Tc - 1.4234e-6 * Tc * Tc * Tc +
1262  1.4168e-8 * Tc * Tc * Tc * Tc;
1263  const Real vCO2 = 32.6 + 3.413e-2 * (Tc - 100.0);
1264  const Real vH2O = 18.1 + 3.137e-2 * (Tc - 100.0);
1265 
1266  const Real delta_pbar = pbar - pref;
1267  const Real Rt = _Rbar * temperature;
1268 
1269  // Equilibrium constants
1270  Real K0H2O, dK0H2O_dT, K0CO2, dK0CO2_dT;
1271  equilibriumConstantH2O(temperature, K0H2O, dK0H2O_dT);
1272  equilibriumConstantCO2(temperature, K0CO2, dK0CO2_dT);
1273 
1274  // Fugacity coefficients
1275  Real phiH2O, phiCO2;
1276  fugacityCoefficientsHighTemp(pressure, temperature, co2_density, xco2, yh2o, phiCO2, phiH2O);
1277 
1278  // Activity coefficients
1279  const Real gammaH2O = activityCoefficientH2O(temperature, xco2);
1280  const Real gammaCO2 = activityCoefficientCO2(temperature, xco2);
1281 
1282  // Activity coefficient for CO2 in brine
1283  Real gamma, dgamma_dT, dgamma_dX;
1284  activityCoefficientHighTemp(temperature, Xnacl, gamma, dgamma_dT, dgamma_dX);
1285 
1286  A = K0H2O * gammaH2O / (phiH2O * pbar) * std::exp(delta_pbar * vH2O / Rt);
1287  B = phiCO2 * pbar / (_invMh2o * K0CO2 * gamma * gammaCO2) * std::exp(-delta_pbar * vCO2 / Rt);
1288 }
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 1291 of file PorousFlowBrineCO2.C.

1304 {
1305  funcABHighTemp(pressure, temperature, Xnacl, co2_density, xco2, yh2o, A, B);
1306 
1307  // Use finite differences for derivatives in the high temperature regime
1308  const Real dp = 1.0e-2;
1309  const Real dT = 1.0e-6;
1310  const Real dX = 1.0e-8;
1311 
1312  Real A2, B2;
1313  funcABHighTemp(pressure + dp, temperature, Xnacl, co2_density, xco2, yh2o, A2, B2);
1314  dA_dp = (A2 - A) / dp;
1315  dB_dp = (B2 - B) / dp;
1316 
1317  funcABHighTemp(pressure, temperature + dT, Xnacl, co2_density, xco2, yh2o, A2, B2);
1318  dA_dT = (A2 - A) / dT;
1319  dB_dT = (B2 - B) / dT;
1320 
1321  funcABHighTemp(pressure, temperature, Xnacl + dX, co2_density, xco2, yh2o, A2, B2);
1322  dB_dX = (B2 - B) / dX;
1323 }
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 1177 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractionsLowTemp().

1188 {
1189  if (temperature > 373.15)
1190  mooseError(name(),
1191  ": funcABLowTemp() is not valid for T > 373.15K. Use funcABHighTemp() instead");
1192 
1193  // Pressure in bar
1194  const Real pbar = pressure * 1.0e-5;
1195 
1196  // Reference pressure and partial molar volumes
1197  const Real pref = 1.0;
1198  const Real vCO2 = 32.6;
1199  const Real vH2O = 18.1;
1200 
1201  const Real delta_pbar = pbar - pref;
1202  const Real Rt = _Rbar * temperature;
1203 
1204  // Equilibrium constants
1205  Real K0H2O, dK0H2O_dT, K0CO2, dK0CO2_dT;
1206  equilibriumConstantH2O(temperature, K0H2O, dK0H2O_dT);
1207  equilibriumConstantCO2(temperature, K0CO2, dK0CO2_dT);
1208 
1209  // Fugacity coefficients
1210  Real phiH2O, dphiH2O_dp, dphiH2O_dT;
1211  Real phiCO2, dphiCO2_dp, dphiCO2_dT;
1213  temperature,
1214  co2_density,
1215  dco2_density_dp,
1216  dco2_density_dT,
1217  phiCO2,
1218  dphiCO2_dp,
1219  dphiCO2_dT,
1220  phiH2O,
1221  dphiH2O_dp,
1222  dphiH2O_dT);
1223 
1224  A = K0H2O / (phiH2O * pbar) * std::exp(delta_pbar * vH2O / Rt);
1225  B = phiCO2 * pbar / (_invMh2o * K0CO2) * std::exp(-delta_pbar * vCO2 / Rt);
1226 
1227  dA_dp = (-1.0e-5 / pbar + 1.0e-5 * vH2O / Rt - dphiH2O_dp / phiH2O) * A;
1228 
1229  dB_dp = (1.0e-5 * phiCO2 + pbar * dphiCO2_dp - 1.0e-5 * vCO2 * pbar * phiCO2 / Rt) *
1230  std::exp(-delta_pbar * vCO2 / Rt) / (_invMh2o * K0CO2);
1231 
1232  dA_dT =
1233  (dK0H2O_dT - dphiH2O_dT * K0H2O / phiH2O - delta_pbar * vH2O * K0H2O / (Rt * temperature)) *
1234  std::exp(delta_pbar * vH2O / Rt) / (pbar * phiH2O);
1235 
1236  dB_dT = (-pbar * phiCO2 * dK0CO2_dT / K0CO2 + pbar * dphiCO2_dT +
1237  delta_pbar * vCO2 * pbar * phiCO2 / (Rt * temperature)) *
1238  std::exp(-delta_pbar * vCO2 / Rt) / (_invMh2o * K0CO2);
1239 }
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 267 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties(), and totalMassFraction().

270 {
272 
273  // Gas density, viscosity and enthalpy are approximated with pure CO2 - no correction due
274  // to the small amount of water vapor is made
275  Real co2_density, dco2_density_dp, dco2_density_dT;
276  Real co2_viscosity, dco2_viscosity_dp, dco2_viscosity_dT;
277  Real co2_enthalpy, dco2_enthalpy_dp, dco2_enthalpy_dT;
279  temperature,
280  co2_density,
281  dco2_density_dp,
282  dco2_density_dT,
283  co2_viscosity,
284  dco2_viscosity_dp,
285  dco2_viscosity_dT);
286 
287  _co2_fp.h_from_p_T(pressure, temperature, co2_enthalpy, dco2_enthalpy_dp, dco2_enthalpy_dT);
288 
289  // Save the values to the FluidStateProperties object. Note that derivatives wrt z are 0
290  gas.density = co2_density;
291  gas.ddensity_dp = dco2_density_dp;
292  gas.ddensity_dT = dco2_density_dT;
293  gas.ddensity_dZ = 0.0;
294 
295  gas.viscosity = co2_viscosity;
296  gas.dviscosity_dp = dco2_viscosity_dp;
297  gas.dviscosity_dT = dco2_viscosity_dT;
298  gas.dviscosity_dZ = 0.0;
299 
300  gas.enthalpy = co2_enthalpy;
301  gas.denthalpy_dp = dco2_enthalpy_dp;
302  gas.denthalpy_dT = dco2_enthalpy_dT;
303  gas.denthalpy_dZ = 0.0;
304 }
virtual void rho_mu_dpT(Real pressure, Real temperature, Real &rho, Real &drho_dp, Real &drho_dT, Real &mu, Real &dmu_dp, Real &dmu_dT) const
Density and viscosity and their derivatives wrt pressure and temperature.
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 Real h_from_p_T(Real p, Real T) const
Specific enthalpy from pressure and temperature.

◆ 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 1455 of file PorousFlowBrineCO2.C.

Referenced by enthalpyOfDissolutionGas().

1457 {
1458  // Henry's constant for dissolution in water
1459  Real Kh_h2o, dKh_h2o_dT;
1460  _co2_fp.henryConstant_dT(temperature, Kh_h2o, dKh_h2o_dT);
1461 
1462  // The correction to salt is obtained through the salting out coefficient
1463  const std::vector<Real> b{1.19784e-1, -7.17823e-4, 4.93854e-6, -1.03826e-8, 1.08233e-11};
1464 
1465  // Need temperature in Celcius
1466  const Real Tc = temperature - _T_c2k;
1467 
1468  Real kb = 0.0;
1469  for (unsigned int i = 0; i < b.size(); ++i)
1470  kb += b[i] * MathUtils::pow(Tc, i);
1471 
1472  Real dkb_dT = 0.0;
1473  for (unsigned int i = 1; i < b.size(); ++i)
1474  dkb_dT += i * b[i] * MathUtils::pow(Tc, i - 1);
1475 
1476  // Need salt mass fraction in molality
1477  const Real xmol = Xnacl / (1.0 - Xnacl) / _Mnacl;
1478  const Real dxmol_dX = 1.0 / (1.0 - Xnacl) / (1.0 - Xnacl) / _Mnacl;
1479  // Henry's constant and its derivative wrt temperature and salt mass fraction
1480  Kh = Kh_h2o * std::pow(10.0, xmol * kb);
1481  dKh_dT = (xmol * std::log(10.0) * dkb_dT + dKh_h2o_dT / Kh_h2o) * Kh;
1482  dKh_dX = std::log(10.0) * kb * dxmol_dX * Kh;
1483 }
const Real _T_c2k
Conversion from C to K.
const Real _Mnacl
Molar mass of NaCL.
virtual void henryConstant_dT(Real temperature, Real &Kh, Real &dKh_dT) const
Henry&#39;s law constant for dissolution in water and derivative wrt temperature.
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)

◆ 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 307 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties(), and totalMassFraction().

311 {
313 
314  // The liquid density includes the density increase due to dissolved CO2
315  Real brine_density, dbrine_density_dp, dbrine_density_dT, dbrine_density_dX;
317  temperature,
318  Xnacl,
319  brine_density,
320  dbrine_density_dp,
321  dbrine_density_dT,
322  dbrine_density_dX);
323 
324  // Mass fraction of CO2 in liquid phase
325  const Real Xco2 = liquid.mass_fraction[_gas_fluid_component];
326  const Real dXco2_dp = liquid.dmass_fraction_dp[_gas_fluid_component];
327  const Real dXco2_dT = liquid.dmass_fraction_dT[_gas_fluid_component];
328  const Real dXco2_dZ = liquid.dmass_fraction_dZ[_gas_fluid_component];
329  const Real dXco2_dX = liquid.dmass_fraction_dX[_gas_fluid_component];
330 
331  // The liquid density
332  Real co2_partial_density, dco2_partial_density_dT;
333  partialDensityCO2(temperature, co2_partial_density, dco2_partial_density_dT);
334 
335  const Real liquid_density = 1.0 / (Xco2 / co2_partial_density + (1.0 - Xco2) / brine_density);
336 
337  const Real dliquid_density_dp =
338  (dXco2_dp / brine_density + (1.0 - Xco2) * dbrine_density_dp / brine_density / brine_density -
339  dXco2_dp / co2_partial_density) *
340  liquid_density * liquid_density;
341 
342  const Real dliquid_density_dT =
343  (dXco2_dT / brine_density + (1.0 - Xco2) * dbrine_density_dT / brine_density / brine_density -
344  dXco2_dT / co2_partial_density +
345  Xco2 * dco2_partial_density_dT / co2_partial_density / co2_partial_density) *
346  liquid_density * liquid_density;
347 
348  const Real dliquid_density_dZ =
349  (dXco2_dZ / brine_density - dXco2_dZ / co2_partial_density) * liquid_density * liquid_density;
350 
351  const Real dliquid_density_dX =
352  (dXco2_dX / brine_density + (1.0 - Xco2) * dbrine_density_dX / brine_density / brine_density -
353  dXco2_dX / co2_partial_density) *
354  liquid_density * liquid_density;
355 
356  // Assume that liquid viscosity is just the brine viscosity
357  Real liquid_viscosity, dliquid_viscosity_dp, dliquid_viscosity_dT, dliquid_viscosity_dX;
359  temperature,
360  Xnacl,
361  liquid_viscosity,
362  dliquid_viscosity_dp,
363  dliquid_viscosity_dT,
364  dliquid_viscosity_dX);
365 
366  // Liquid enthalpy (including contribution due to the enthalpy of dissolution)
367  Real brine_enthalpy, dbrine_enthalpy_dp, dbrine_enthalpy_dT, dbrine_enthalpy_dX;
369  temperature,
370  Xnacl,
371  brine_enthalpy,
372  dbrine_enthalpy_dp,
373  dbrine_enthalpy_dT,
374  dbrine_enthalpy_dX);
375 
376  // Enthalpy of CO2
377  Real co2_enthalpy, dco2_enthalpy_dp, dco2_enthalpy_dT;
378  _co2_fp.h_from_p_T(pressure, temperature, co2_enthalpy, dco2_enthalpy_dp, dco2_enthalpy_dT);
379 
380  // Enthalpy of dissolution
381  Real hdis, dhdis_dT;
382  enthalpyOfDissolution(temperature, hdis, dhdis_dT);
383 
384  const Real liquid_enthalpy = (1.0 - Xco2) * brine_enthalpy + Xco2 * (co2_enthalpy + hdis);
385  const Real dliquid_enthalpy_dp = (1.0 - Xco2) * dbrine_enthalpy_dp + Xco2 * dco2_enthalpy_dp +
386  dXco2_dp * (co2_enthalpy + hdis - brine_enthalpy);
387  const Real dliquid_enthalpy_dT = (1.0 - Xco2) * dbrine_enthalpy_dT +
388  Xco2 * (dco2_enthalpy_dT + dhdis_dT) +
389  dXco2_dT * (co2_enthalpy + hdis - brine_enthalpy);
390  const Real dliquid_enthalpy_dZ = dXco2_dZ * (co2_enthalpy + hdis - brine_enthalpy);
391  const Real dliquid_enthalpy_dX =
392  (1.0 - Xco2) * dbrine_enthalpy_dX + dXco2_dX * (co2_enthalpy + hdis - brine_enthalpy);
393 
394  // Save the values to the FluidStateProperties object
395  liquid.density = liquid_density;
396  liquid.ddensity_dp = dliquid_density_dp;
397  liquid.ddensity_dT = dliquid_density_dT;
398  liquid.ddensity_dZ = dliquid_density_dZ;
399  liquid.ddensity_dX = dliquid_density_dX;
400 
401  liquid.viscosity = liquid_viscosity;
402  liquid.dviscosity_dp = dliquid_viscosity_dp;
403  liquid.dviscosity_dT = dliquid_viscosity_dT;
404  liquid.dviscosity_dZ = 0.0;
405  liquid.dviscosity_dX = dliquid_viscosity_dX;
406 
407  liquid.enthalpy = liquid_enthalpy;
408  liquid.denthalpy_dp = dliquid_enthalpy_dp;
409  liquid.denthalpy_dT = dliquid_enthalpy_dT;
410  liquid.denthalpy_dZ = dliquid_enthalpy_dZ;
411  liquid.denthalpy_dX = dliquid_enthalpy_dX;
412 }
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
virtual void rho_dpTx(Real pressure, Real temperature, Real xnacl, Real &rho, Real &drho_dp, Real &drho_dT, Real &drho_dx) const override
Density and its derivatives wrt pressure, temperature and mass fraction.
const std::string temperature
Definition: NS.h:27
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.
virtual void h_dpTx(Real pressure, Real temperature, Real xnacl, Real &h, Real &dh_dp, Real &dh_dT, Real &dh_dx) const override
Enthalpy and derivatives wrt pressure, temperature and mass fraction.
std::vector< Real > dmass_fraction_dp
const std::string pressure
Definition: NS.h:26
virtual void mu_dpTx(Real pressure, Real temperature, Real xnacl, Real &mu, Real &dmu_dp, Real &dmu_dT, Real &dmu_dx) const override
virtual Real h_from_p_T(Real p, Real T) const
Specific enthalpy from pressure and temperature.
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)

◆ 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 159 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties().

165 {
168 
169  // Equilibrium mass fraction of CO2 in liquid and H2O in gas phases
170  Real Xco2, dXco2_dp, dXco2_dT, dXco2_dX, Yh2o, dYh2o_dp, dYh2o_dT, dYh2o_dX;
172  temperature,
173  Xnacl,
174  Xco2,
175  dXco2_dp,
176  dXco2_dT,
177  dXco2_dX,
178  Yh2o,
179  dYh2o_dp,
180  dYh2o_dT,
181  dYh2o_dX);
182 
183  Real Yco2 = 1.0 - Yh2o;
184  Real dYco2_dp = -dYh2o_dp;
185  Real dYco2_dT = -dYh2o_dT;
186  Real dYco2_dX = -dYh2o_dX;
187 
188  // Determine which phases are present based on the value of z
189  phaseState(Z, Xco2, Yco2, phase_state);
190 
191  // The equilibrium mass fractions calculated above are only correct in the two phase
192  // state. If only liquid or gas phases are present, the mass fractions are given by
193  // the total mass fraction z
194  Real Xh2o = 0.0;
195  Real dXco2_dZ = 0.0, dYco2_dZ = 0.0;
196 
197  switch (phase_state)
198  {
200  {
201  Xco2 = Z;
202  Yco2 = 0.0;
203  Xh2o = 1.0 - Z;
204  Yh2o = 0.0;
205  dXco2_dp = 0.0;
206  dXco2_dT = 0.0;
207  dXco2_dX = 0.0;
208  dXco2_dZ = 1.0;
209  dYco2_dp = 0.0;
210  dYco2_dT = 0.0;
211  dYco2_dX = 0.0;
212  break;
213  }
214 
216  {
217  Xco2 = 0.0;
218  Yco2 = Z;
219  Yh2o = 1.0 - Z;
220  dXco2_dp = 0.0;
221  dXco2_dT = 0.0;
222  dXco2_dX = 0.0;
223  dYco2_dZ = 1.0;
224  dYco2_dp = 0.0;
225  dYco2_dT = 0.0;
226  dYco2_dX = 0.0;
227  break;
228  }
229 
231  {
232  // Keep equilibrium mass fractions
233  Xh2o = 1.0 - Xco2;
234  break;
235  }
236  }
237 
238  // Save the mass fractions in the FluidStateProperties object
240  liquid.mass_fraction[_gas_fluid_component] = Xco2;
241  liquid.mass_fraction[_salt_component] = Xnacl;
244 
245  // Save the derivatives wrt PorousFlow variables
246  liquid.dmass_fraction_dp[_aqueous_fluid_component] = -dXco2_dp;
247  liquid.dmass_fraction_dp[_gas_fluid_component] = dXco2_dp;
248  liquid.dmass_fraction_dT[_aqueous_fluid_component] = -dXco2_dT;
249  liquid.dmass_fraction_dT[_gas_fluid_component] = dXco2_dT;
250  liquid.dmass_fraction_dX[_aqueous_fluid_component] = -dXco2_dX;
251  liquid.dmass_fraction_dX[_gas_fluid_component] = dXco2_dX;
252  liquid.dmass_fraction_dX[_salt_component] = 1.0;
253  liquid.dmass_fraction_dZ[_aqueous_fluid_component] = -dXco2_dZ;
254  liquid.dmass_fraction_dZ[_gas_fluid_component] = dXco2_dZ;
255 
257  gas.dmass_fraction_dp[_gas_fluid_component] = dYco2_dp;
259  gas.dmass_fraction_dT[_gas_fluid_component] = dYco2_dT;
261  gas.dmass_fraction_dX[_gas_fluid_component] = dYco2_dX;
263  gas.dmass_fraction_dZ[_gas_fluid_component] = dYco2_dZ;
264 }
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

◆ 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 1392 of file PorousFlowBrineCO2.C.

Referenced by liquidProperties(), and saturationTwoPhase().

1395 {
1396  // This correlation uses temperature in C
1397  const Real Tc = temperature - _T_c2k;
1398  // The parial molar volume
1399  const Real V = 37.51 - 9.585e-2 * Tc + 8.74e-4 * Tc * Tc - 5.044e-7 * Tc * Tc * Tc;
1400  const Real dV_dT = -9.585e-2 + 1.748e-3 * Tc - 1.5132e-6 * Tc * Tc;
1401 
1402  partial_density = 1.0e6 * _Mco2 / V;
1403  dpartial_density_dT = -1.0e6 * _Mco2 * dV_dT / V / V;
1404 }
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 415 of file PorousFlowBrineCO2.C.

Referenced by thermophysicalProperties().

420 {
423 
424  // Approximate liquid density as saturation isn't known yet
425  Real brine_density, dbrine_density_dp, dbrine_density_dT, dbrine_density_dX;
427  temperature,
428  Xnacl,
429  brine_density,
430  dbrine_density_dp,
431  dbrine_density_dT,
432  dbrine_density_dX);
433 
434  // Mass fraction of CO2 in liquid phase
435  const Real Xco2 = liquid.mass_fraction[_gas_fluid_component];
436  const Real dXco2_dp = liquid.dmass_fraction_dp[_gas_fluid_component];
437  const Real dXco2_dT = liquid.dmass_fraction_dT[_gas_fluid_component];
438  const Real dXco2_dX = liquid.dmass_fraction_dX[_gas_fluid_component];
439 
440  // The liquid density
441  Real co2_partial_density, dco2_partial_density_dT;
442  partialDensityCO2(temperature, co2_partial_density, dco2_partial_density_dT);
443 
444  const Real liquid_density = 1.0 / (Xco2 / co2_partial_density + (1.0 - Xco2) / brine_density);
445 
446  const Real dliquid_density_dp =
447  (dXco2_dp / brine_density + (1.0 - Xco2) * dbrine_density_dp / brine_density / brine_density -
448  dXco2_dp / co2_partial_density) *
449  liquid_density * liquid_density;
450 
451  const Real dliquid_density_dT =
452  (dXco2_dT / brine_density + (1.0 - Xco2) * dbrine_density_dT / brine_density / brine_density -
453  dXco2_dT / co2_partial_density +
454  Xco2 * dco2_partial_density_dT / co2_partial_density / co2_partial_density) *
455  liquid_density * liquid_density;
456 
457  const Real dliquid_density_dX =
458  (dXco2_dX / brine_density + (1.0 - Xco2) * dbrine_density_dX / brine_density / brine_density -
459  dXco2_dX / co2_partial_density) *
460  liquid_density * liquid_density;
461 
462  const Real Yco2 = gas.mass_fraction[_gas_fluid_component];
463  const Real dYco2_dp = gas.dmass_fraction_dp[_gas_fluid_component];
464  const Real dYco2_dT = gas.dmass_fraction_dT[_gas_fluid_component];
465  const Real dYco2_dX = gas.dmass_fraction_dX[_gas_fluid_component];
466 
467  // Set mass equilibrium constants used in the calculation of vapor mass fraction
468  const Real K0 = Yco2 / Xco2;
469  const Real K1 = (1.0 - Yco2) / (1.0 - Xco2);
470  const Real vapor_mass_fraction = vaporMassFraction(Z, K0, K1);
471 
472  // The gas saturation in the two phase case
473  gas.saturation = vapor_mass_fraction * liquid_density /
474  (gas.density + vapor_mass_fraction * (liquid_density - gas.density));
475 
476  const Real dv_dZ = (K1 - K0) / ((K0 - 1.0) * (K1 - 1.0));
477  const Real denominator = (gas.density + vapor_mass_fraction * (liquid_density - gas.density)) *
478  (gas.density + vapor_mass_fraction * (liquid_density - gas.density));
479 
480  const Real ds_dZ = gas.density * liquid_density * dv_dZ / denominator;
481 
482  const Real dK0_dp = (Xco2 * dYco2_dp - Yco2 * dXco2_dp) / Xco2 / Xco2;
483  const Real dK0_dT = (Xco2 * dYco2_dT - Yco2 * dXco2_dT) / Xco2 / Xco2;
484  const Real dK0_dX = (Xco2 * dYco2_dX - Yco2 * dXco2_dX) / Xco2 / Xco2;
485 
486  const Real dK1_dp =
487  ((1.0 - Yco2) * dXco2_dp - (1.0 - Xco2) * dYco2_dp) / (1.0 - Xco2) / (1.0 - Xco2);
488  const Real dK1_dT =
489  ((1.0 - Yco2) * dXco2_dT - (1.0 - Xco2) * dYco2_dT) / (1.0 - Xco2) / (1.0 - Xco2);
490  const Real dK1_dX =
491  ((1.0 - Yco2) * dXco2_dX - (1.0 - Xco2) * dYco2_dX) / (1.0 - Xco2) / (1.0 - Xco2);
492 
493  const Real dv_dp =
494  Z * dK1_dp / (K1 - 1.0) / (K1 - 1.0) + (1.0 - Z) * dK0_dp / (K0 - 1.0) / (K0 - 1.0);
495 
496  Real ds_dp = gas.density * liquid_density * dv_dp +
497  vapor_mass_fraction * (1.0 - vapor_mass_fraction) *
498  (gas.density * dliquid_density_dp - gas.ddensity_dp * liquid_density);
499  ds_dp /= denominator;
500 
501  const Real dv_dT =
502  Z * dK1_dT / (K1 - 1.0) / (K1 - 1.0) + (1.0 - Z) * dK0_dT / (K0 - 1.0) / (K0 - 1.0);
503 
504  Real ds_dT = gas.density * liquid_density * dv_dT +
505  vapor_mass_fraction * (1.0 - vapor_mass_fraction) *
506  (gas.density * dliquid_density_dT - gas.ddensity_dT * liquid_density);
507  ds_dT /= denominator;
508 
509  const Real dv_dX =
510  Z * dK1_dX / (K1 - 1.0) / (K1 - 1.0) + (1.0 - Z) * dK0_dX / (K0 - 1.0) / (K0 - 1.0);
511 
512  Real ds_dX = gas.density * liquid_density * dv_dX + vapor_mass_fraction *
513  (1.0 - vapor_mass_fraction) *
514  (gas.density * dliquid_density_dX);
515  ds_dX /= denominator;
516 
517  gas.dsaturation_dp = ds_dp;
518  gas.dsaturation_dT = ds_dT;
519  gas.dsaturation_dZ = ds_dZ;
520  gas.dsaturation_dX = ds_dX;
521 }
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
virtual void rho_dpTx(Real pressure, Real temperature, Real xnacl, Real &rho, Real &drho_dp, Real &drho_dT, Real &drho_dx) const override
Density and its derivatives wrt pressure, temperature and mass fraction.
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 1522 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractions().

1524 {
1525  // Coefficients of cubic polynomial
1526  const Real dT = _Tupper - _Tlower;
1527 
1528  const Real a = f0;
1529  const Real b = df0 * dT;
1530  const Real c = 3.0 * (f1 - f0) - (2.0 * df0 + df1) * dT;
1531  const Real d = 2.0 * (f0 - f1) + (df0 + df1) * dT;
1532 
1533  const Real t2 = temperature * temperature;
1534  const Real t3 = temperature * t2;
1535 
1536  value = a + b * temperature + c * t2 + d * t3;
1537  deriv = b + 2.0 * c * temperature + 3.0 * d * t2;
1538 }
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 1326 of file PorousFlowBrineCO2.C.

Referenced by equilibriumMoleFractions().

1328 {
1329  // Initial guess for yh2o and xco2 (from Spycher and Pruess (2010))
1330  Real y = _brine_fp.vaporPressure(temperature, 0.0) / pressure;
1331  Real x = 0.009;
1332 
1333  // Need salt mass fraction in molality
1334  const Real mnacl = Xnacl / (1.0 - Xnacl) / _Mnacl;
1335 
1336  // If y > 1, then just use y = 1, x = 0 (only a gas phase)
1337  if (y >= 1.0)
1338  {
1339  y = 1.0;
1340  x = 0.0;
1341  }
1342  else
1343  {
1344  // Residual function for Newton-Raphson
1345  auto fy = [mnacl, this](Real y, Real A, Real B) {
1346  return y -
1347  (1.0 - B) * _invMh2o / ((1.0 / A - B) * (2.0 * mnacl + _invMh2o) + 2.0 * mnacl * B);
1348  };
1349 
1350  // Derivative of fy wrt y
1351  auto dfy = [mnacl, this](Real A, Real B, Real dA, Real dB) {
1352  const Real denominator = (1.0 / A - B) * (2.0 * mnacl + _invMh2o) + 2.0 * mnacl * B;
1353  return 1.0 + _invMh2o * dB / denominator +
1354  (1.0 - B) * _invMh2o *
1355  (2.0 * mnacl * dB - (2.0 * mnacl + _invMh2o) * (dB + dA / A / A)) / denominator /
1356  denominator;
1357  };
1358 
1359  Real A, B;
1360  Real dA, dB;
1361  const Real dy = 1.0e-8;
1362 
1363  // Solve for yh2o using Newton-Raphson method
1364  unsigned int iter = 0;
1365  const Real tol = 1.0e-12;
1366  const unsigned int max_its = 10;
1367  funcABHighTemp(pressure, temperature, Xnacl, co2_density, x, y, A, B);
1368 
1369  while (std::abs(fy(y, A, B)) > tol)
1370  {
1371  funcABHighTemp(pressure, temperature, Xnacl, co2_density, x, y, A, B);
1372  // Finite difference derivatives of A and B wrt y
1373  funcABHighTemp(pressure, temperature, Xnacl, co2_density, x, y + dy, dA, dB);
1374  dA = (dA - A) / dy;
1375  dB = (dB - B) / dy;
1376 
1377  y = y - fy(y, A, B) / dfy(A, B, dA, dB);
1378 
1379  x = B * (1.0 - y);
1380 
1381  // Break if not converged and just use the value
1382  if (iter > max_its)
1383  break;
1384  }
1385  }
1386 
1387  yh2o = y;
1388  xco2 = x;
1389 }
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 89 of file PorousFlowBrineCO2.C.

95 {
98 
99  // Check whether the input temperature is within the region of validity
101 
102  // Clear all of the FluidStateProperties data
104 
105  FluidStatePhaseEnum phase_state;
106  massFractions(pressure, temperature, Xnacl, Z, phase_state, fsp);
107 
108  switch (phase_state)
109  {
111  {
112  // Set the gas saturations
113  gas.saturation = 1.0;
114 
115  // Calculate gas properties
117 
118  break;
119  }
120 
122  {
123  // Calculate the liquid properties
124  Real liquid_pressure = pressure - _pc.capillaryPressure(1.0, qp);
125  liquidProperties(liquid_pressure, temperature, Xnacl, fsp);
126 
127  break;
128  }
129 
131  {
132  // Calculate the gas properties
134 
135  // Calculate the saturation
136  saturationTwoPhase(pressure, temperature, Xnacl, Z, fsp);
137 
138  // Calculate the liquid properties
139  Real liquid_pressure = pressure - _pc.capillaryPressure(1.0 - gas.saturation, qp);
140  liquidProperties(liquid_pressure, temperature, Xnacl, fsp);
141 
142  break;
143  }
144  }
145 
146  // Liquid saturations can now be set
147  liquid.saturation = 1.0 - gas.saturation;
148  liquid.dsaturation_dp = -gas.dsaturation_dp;
149  liquid.dsaturation_dT = -gas.dsaturation_dT;
150  liquid.dsaturation_dX = -gas.dsaturation_dX;
151  liquid.dsaturation_dZ = -gas.dsaturation_dZ;
152 
153  // Save pressures to FluidStateProperties object
154  gas.pressure = pressure;
155  liquid.pressure = pressure - _pc.capillaryPressure(liquid.saturation, qp);
156 }
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 1407 of file PorousFlowBrineCO2.C.

1409 {
1410  // Check whether the input pressure and temperature are within the region of validity
1412 
1413  // FluidStateProperties data structure
1414  std::vector<FluidStateProperties> fsp(_num_phases, FluidStateProperties(_num_components));
1417 
1418  // Calculate equilibrium mass fractions in the two-phase state
1419  Real Xco2, dXco2_dp, dXco2_dT, dXco2_dX, Yh2o, dYh2o_dp, dYh2o_dT, dYh2o_dX;
1421  temperature,
1422  Xnacl,
1423  Xco2,
1424  dXco2_dp,
1425  dXco2_dT,
1426  dXco2_dX,
1427  Yh2o,
1428  dYh2o_dp,
1429  dYh2o_dT,
1430  dYh2o_dX);
1431 
1432  // Save the mass fractions in the FluidStateMassFractions object
1433  const Real Yco2 = 1.0 - Yh2o;
1434  liquid.mass_fraction[_aqueous_fluid_component] = 1.0 - Xco2;
1435  liquid.mass_fraction[_gas_fluid_component] = Xco2;
1437  gas.mass_fraction[_gas_fluid_component] = Yco2;
1438 
1439  // Gas properties
1441 
1442  // Liquid properties
1443  const Real liquid_saturation = 1.0 - saturation;
1444  const Real liquid_pressure = pressure - _pc.capillaryPressure(liquid_saturation, qp);
1445  liquidProperties(liquid_pressure, temperature, Xnacl, fsp);
1446 
1447  // The total mass fraction of ncg (z) can now be calculated
1448  const Real Z = (saturation * gas.density * Yco2 + liquid_saturation * liquid.density * Xco2) /
1449  (saturation * gas.density + liquid_saturation * liquid.density);
1450 
1451  return Z;
1452 }
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.


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