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

Fluid state class for brine and CO2. More...

#include <PorousFlowFluidStateBrineCO2.h>

Inheritance diagram for PorousFlowFluidStateBrineCO2:
[legend]

Public Member Functions

 PorousFlowFluidStateBrineCO2 (const InputParameters &parameters)
 

Protected Member Functions

virtual void thermophysicalProperties () override
 Calculates all required thermophysical properties and derivatives for each phase and fluid component. More...
 
void massFractions (Real pressure, Real temperature, Real xnacl, Real &xco2l, Real &dxco2l_dp, Real &dxco2l_dT, Real &xh2og, Real &dxh2og_dp, Real &dxh2og_dT) const
 Mass fractions of CO2 and brine calculated using mutual solubilities given by Spycher, Pruess and Ennis-King, CO2-H2O mixtures in the geological sequestration of CO2. More...
 
void fugacityCoefficientCO2 (Real pressure, Real temperature, Real &fco2, Real &dfco2_dp, Real &dfco2_dT) const
 Fugacity coefficient for CO2. More...
 
void fugacityCoefficientH2O (Real pressure, Real temperature, Real &fh2o, Real &dfh2o_dp, Real &dfh2o_dT) const
 Fugacity coefficient for H2O. More...
 
void activityCoefficient (Real pressure, Real temperature, Real xnacl, Real &gamma, Real &dgamma_dp, Real &dgamma_dT) const
 Activity coefficient for CO2 in brine. More...
 
void equilibriumConstantH2O (Real temperature, Real &kh2o, Real &dkh2o_dT) const
 Equilibrium constant for H2O from Spycher, Pruess and Ennis-King, CO2-H2O mixtures in the geological sequestration of CO2. More...
 
void equilibriumConstantCO2 (Real temperature, Real &kco2, Real &dkco2_dT) const
 Equilibrium constant for CO2 from Spycher, Pruess and Ennis-King, CO2-H2O mixtures in the geological sequestration of CO2. 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 void initQpStatefulProperties () override
 
virtual void computeQpProperties () override
 
void setMaterialVectorSize () const
 Size material property vectors and initialise with zeros. More...
 
Real rachfordRice (Real x, std::vector< Real > &Ki) const
 Rachford-Rice equation for vapor fraction. More...
 
Real rachfordRiceDeriv (Real x, std::vector< Real > &Ki) const
 Derivative of Rachford-Rice equation wrt vapor fraction. More...
 
Real vaporMassFraction (std::vector< Real > &Ki) const
 Solves Rachford-Rice equation to provide vapor mass fraction. More...
 

Protected Attributes

const VariableValue & _xnacl
 Salt mass fraction (kg/kg) More...
 
const BrineFluidProperties_brine_fp
 Fluid properties UserObject for brine. More...
 
const SinglePhaseFluidPropertiesPT_co2_fp
 Fluid properties UserObject for CO2. More...
 
const SinglePhaseFluidPropertiesPT_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 VariableValue & _gas_porepressure
 Porepressure. More...
 
const VariableGradient & _gas_gradp_qp
 Gradient of porepressure (only defined at the qps) More...
 
const unsigned int _gas_porepressure_varnum
 Moose variable number of the gas porepressure. More...
 
const unsigned int _pvar
 PorousFlow variable number of the gas porepressure. More...
 
std::vector< const VariableValue * > _z
 Total mass fraction(s) of the gas component(s) summed over all phases. More...
 
std::vector< const VariableGradient * > _gradz_qp
 Gradient(s) of total mass fraction(s) of the gas component(s) (only defined at the qps) More...
 
std::vector< unsigned int > _z_varnum
 Moose variable number of z. More...
 
std::vector< unsigned int > _zvar
 PorousFlow variable number of z. More...
 
const unsigned int _num_z_vars
 Number of coupled total mass fractions. Should be _num_phases - 1. More...
 
const unsigned int _aqueous_phase_number
 Phase number of the aqueous phase. More...
 
const 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...
 
const unsigned int _gas_fluid_component
 Fluid component number of the gas phase. More...
 
const MaterialProperty< Real > & _temperature
 Temperature. More...
 
const MaterialProperty< RealGradient > & _gradT_qp
 Gradient of temperature (only defined at the qps) More...
 
const MaterialProperty< std::vector< Real > > & _dtemperature_dvar
 Derivative of temperature wrt PorousFlow variables. More...
 
MaterialProperty< std::vector< std::vector< Real > > > & _mass_frac
 Mass fraction matrix. More...
 
MaterialProperty< std::vector< std::vector< RealGradient > > > * _grad_mass_frac_qp
 Gradient of the mass fraction matrix (only defined at the qps) More...
 
MaterialProperty< std::vector< std::vector< std::vector< Real > > > > & _dmass_frac_dvar
 Derivative of the mass fraction matrix with respect to the Porous Flow variables. More...
 
const MaterialProperty< std::vector< Real > > & _saturation_old
 Old value of saturation. More...
 
MaterialProperty< std::vector< Real > > & _fluid_density
 Fluid density of each phase. More...
 
MaterialProperty< std::vector< std::vector< Real > > > & _dfluid_density_dvar
 Derivative of the fluid density for each phase wrt PorousFlow variables. More...
 
MaterialProperty< std::vector< Real > > & _fluid_viscosity
 Viscosity of each phase. More...
 
MaterialProperty< std::vector< std::vector< Real > > > & _dfluid_viscosity_dvar
 Derivative of the fluid viscosity for each phase wrt PorousFlow variables. More...
 
const Real _T_c2k
 Conversion from degrees Celsius to degrees Kelvin. More...
 
const Real _R
 Universal gas constant (J/mol/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...
 
bool _is_initqp
 Flag to indicate whether to calculate stateful properties. More...
 
std::vector< FluidStateProperties_fsp
 FluidStateProperties data structure. More...
 
const PorousFlowCapillaryPressure_pc_uo
 Capillary pressure UserObject. More...
 
const unsigned int _num_phases
 Number of phases. More...
 
const unsigned int _num_components
 Number of components. More...
 
const unsigned int _num_pf_vars
 Number of PorousFlow variables. More...
 
MaterialProperty< std::vector< Real > > & _porepressure
 Computed nodal or quadpoint values of porepressure of the phases. More...
 
MaterialProperty< std::vector< std::vector< Real > > > & _dporepressure_dvar
 d(porepressure)/d(PorousFlow variable) More...
 
MaterialProperty< std::vector< RealGradient > > *const _gradp_qp
 Grad(p) at the quadpoints. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dgradp_qp_dgradv
 d(grad porepressure)/d(grad PorousFlow variable) at the quadpoints More...
 
MaterialProperty< std::vector< std::vector< RealGradient > > > *const _dgradp_qp_dv
 d(grad porepressure)/d(PorousFlow variable) at the quadpoints More...
 
MaterialProperty< std::vector< Real > > & _saturation
 Computed nodal or qp saturation of the phases. More...
 
MaterialProperty< std::vector< std::vector< Real > > > & _dsaturation_dvar
 d(saturation)/d(PorousFlow variable) More...
 
MaterialProperty< std::vector< RealGradient > > *const _grads_qp
 Grad(s) at the quadpoints. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dgrads_qp_dgradv
 d(grad saturation)/d(grad PorousFlow variable) at the quadpoints More...
 
MaterialProperty< std::vector< std::vector< RealGradient > > > *const _dgrads_qp_dv
 d(grad saturation)/d(PorousFlow variable) at the quadpoints More...
 

Detailed Description

Fluid state class for brine and CO2.

Includes mutual solubility of CO2 and brine using model of 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)

Definition at line 30 of file PorousFlowFluidStateBrineCO2.h.

Constructor & Destructor Documentation

PorousFlowFluidStateBrineCO2::PorousFlowFluidStateBrineCO2 ( const InputParameters &  parameters)

Definition at line 25 of file PorousFlowFluidStateBrineCO2.C.

26  : PorousFlowFluidStateFlashBase(parameters),
27 
28  _xnacl(_nodal_material ? coupledNodalValue("xnacl") : coupledValue("xnacl")),
29  _brine_fp(getUserObject<BrineFluidProperties>("brine_fp")),
30  _co2_fp(getUserObject<SinglePhaseFluidPropertiesPT>("co2_fp")),
33  _invMh2o(1.0 / _Mh2o),
36  _Rbar(_R * 10.0)
37 {
38  // Check that the correct FluidProperties UserObjects have been provided
39  if (_brine_fp.fluidName() != "brine")
40  mooseError("Only a valid Brine FluidProperties UserObject can be provided in brine_fp");
41 
42  if (_co2_fp.fluidName() != "co2")
43  mooseError("Only a valid CO2 FluidProperties UserObject can be provided in co2_fp");
44 }
virtual std::string fluidName() const =0
Fluid name.
PorousFlowFluidStateFlashBase(const InputParameters &parameters)
const Real _Mnacl
Molar mass of NaCL.
const SinglePhaseFluidPropertiesPT & _water_fp
Fluid properties UserObject for H20.
const BrineFluidProperties & _brine_fp
Fluid properties UserObject for brine.
virtual std::string fluidName() const override
Fluid name.
const VariableValue & _xnacl
Salt mass fraction (kg/kg)
virtual Real molarMass() const =0
Molar mass.
const Real _Mh2o
Molar mass of water (kg/mol)
const SinglePhaseFluidPropertiesPT & _co2_fp
Fluid properties UserObject for CO2.
const Real _R
Universal gas constant (J/mol/K)
const Real _invMh2o
Inverse of molar mass of H2O (mol/kg)
const Real _Mco2
Molar mass of CO2 (kg/mol)
Real molarMassNaCl() const
NaCl molar mass.
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
virtual const SinglePhaseFluidPropertiesPT & getComponent(unsigned int component) const override
Get UserObject for specified component.
Real molarMassH2O() const
H2O molar mass.
static const unsigned int WATER
Fluid component numbers for water and NaCl.

Member Function Documentation

void PorousFlowFluidStateBrineCO2::activityCoefficient ( Real  pressure,
Real  temperature,
Real  xnacl,
Real &  gamma,
Real &  dgamma_dp,
Real &  dgamma_dT 
) const
protected

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)

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

Definition at line 562 of file PorousFlowFluidStateBrineCO2.C.

Referenced by massFractions().

568 {
569  // Need pressure in bar
570  Real pbar = pressure * 1.0e-5;
571  // Need NaCl molality (mol/kg)
572  Real mnacl = xnacl / (1.0 - xnacl) / _Mnacl;
573 
574  Real lambda = -0.411370585 + 6.07632013e-4 * temperature + 97.5347708 / temperature -
575  0.0237622469 * pbar / temperature + 0.0170656236 * pbar / (630.0 - temperature) +
576  1.41335834e-5 * temperature * std::log(pbar);
577 
578  Real xi = 3.36389723e-4 - 1.9829898e-5 * temperature + 2.12220830e-3 * pbar / temperature -
579  5.24873303e-3 * pbar / (630.0 - temperature);
580 
581  gamma = std::exp(2.0 * lambda * mnacl + xi * mnacl * mnacl);
582 
583  // Derivative wrt pressure
584  Real dlambda_dp, dxi_dp;
585  dlambda_dp = -0.0237622469 / temperature + 0.0170656236 / (630.0 - temperature) +
586  1.41335834e-5 * temperature / pbar;
587  dxi_dp = 2.12220830e-3 / temperature - 5.24873303e-3 / (630.0 - temperature);
588  dgamma_dp = (2.0 * mnacl * dlambda_dp + mnacl * mnacl * dxi_dp) *
589  std::exp(2.0 * lambda * mnacl + xi * mnacl * mnacl) * 1.0e-5;
590 
591  // Derivative wrt temperature
592  Real dlambda_dT, dxi_dT;
593  dlambda_dT = 6.07632013e-4 - 97.5347708 / temperature / temperature +
594  0.0237622469 * pbar / temperature / temperature +
595  0.0170656236 * pbar / (630.0 - temperature) / (630.0 - temperature) +
596  1.41335834e-5 * std::log(pbar);
597  dxi_dT = -1.9829898e-5 - 2.12220830e-3 * pbar / temperature / temperature -
598  5.24873303e-3 * pbar / (630.0 - temperature) / (630.0 - temperature);
599  dgamma_dT = (2.0 * mnacl * dlambda_dT + mnacl * mnacl * dxi_dT) *
600  std::exp(2.0 * lambda * mnacl + xi * mnacl * mnacl);
601 }
const Real _Mnacl
Molar mass of NaCL.
const std::string temperature
Definition: NS.h:25
const std::string pressure
Definition: NS.h:24
void PorousFlowFluidStateFlashBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Reimplemented from PorousFlowVariableBase.

Definition at line 161 of file PorousFlowFluidStateFlashBase.C.

162 {
163  _is_initqp = false;
164  // Prepare the derivative vectors
166 
167  // Set the size of all other vectors
169 
170  // Calculate all required thermophysical properties
172 
173  for (unsigned int ph = 0; ph < _num_phases; ++ph)
174  {
175  _saturation[_qp][ph] = _fsp[ph].saturation;
176  _porepressure[_qp][ph] = _fsp[ph].pressure;
177  _fluid_density[_qp][ph] = _fsp[ph].fluid_density;
178  _fluid_viscosity[_qp][ph] = _fsp[ph].fluid_viscosity;
179  _mass_frac[_qp][ph] = _fsp[ph].mass_fraction;
180  }
181 
182  // Derivative of saturation wrt variables
183  for (unsigned int ph = 0; ph < _num_phases; ++ph)
184  {
185  _dsaturation_dvar[_qp][ph][_zvar[0]] = _fsp[ph].dsaturation_dz;
186  _dsaturation_dvar[_qp][ph][_pvar] = _fsp[ph].dsaturation_dp;
187  }
188  // Derivative of capillary pressure
190 
191  // Derivative of porepressure wrt variables
192  if (_dictator.isPorousFlowVariable(_gas_porepressure_varnum))
193  {
194  for (unsigned int ph = 0; ph < _num_phases; ++ph)
195  {
196  _dporepressure_dvar[_qp][ph][_pvar] = 1.0;
197  if (!_nodal_material)
198  (*_dgradp_qp_dgradv)[_qp][ph][_pvar] = 1.0;
199  }
200 
201  if (!_nodal_material)
202  (*_dgradp_qp_dgradv)[_qp][_aqueous_phase_number][_zvar[0]] =
203  -dpc * _fsp[_aqueous_phase_number].dsaturation_dp -
204  dpc * _fsp[_aqueous_phase_number].dsaturation_dz;
205 
206  // The aqueous phase porepressure is also a function of liquid saturation,
207  // which depends on both gas porepressure and z
211  -dpc * _dsaturation_dvar[_qp][_aqueous_phase_number][_zvar[0]];
212  }
213 
214  // Calculate derivatives of material properties wrt primary variables
215  // Derivative of z wrt variables
216  std::vector<Real> dz_dvar;
217  dz_dvar.assign(_num_pf_vars, 0.0);
218  if (_dictator.isPorousFlowVariable(_z_varnum[0]))
219  dz_dvar[_zvar[0]] = 1.0;
220 
221  // Derivatives of properties wrt primary variables
222  for (unsigned int v = 0; v < _num_pf_vars; ++v)
223  {
224  for (unsigned int ph = 0; ph < _num_phases; ++ph)
225  {
226  // Derivative of density in each phase
227  _dfluid_density_dvar[_qp][ph][v] =
228  _fsp[ph].dfluid_density_dp * _dporepressure_dvar[_qp][ph][v];
229  _dfluid_density_dvar[_qp][ph][v] += _fsp[ph].dfluid_density_dT * _dtemperature_dvar[_qp][v];
230  _dfluid_density_dvar[_qp][ph][v] += _fsp[ph].dfluid_density_dz * dz_dvar[v];
231 
232  // Derivative of viscosity in each phase
233  _dfluid_viscosity_dvar[_qp][ph][v] =
234  _fsp[ph].dfluid_viscosity_dp * _dporepressure_dvar[_qp][ph][v];
235  _dfluid_viscosity_dvar[_qp][ph][v] +=
236  _fsp[ph].dfluid_viscosity_dT * _dtemperature_dvar[_qp][v];
237  _dfluid_viscosity_dvar[_qp][ph][v] += _fsp[ph].dfluid_viscosity_dz * dz_dvar[v];
238 
239  // The derivative of the mass fractions for each fluid component in each phase
240  for (unsigned int comp = 0; comp < _num_components; ++comp)
241  {
242  _dmass_frac_dvar[_qp][ph][comp][v] =
243  _fsp[ph].dmass_fraction_dp[comp] * _dporepressure_dvar[_qp][ph][v];
244  _dmass_frac_dvar[_qp][ph][comp][v] +=
245  _fsp[ph].dmass_fraction_dT[comp] * _dtemperature_dvar[_qp][v];
246  _dmass_frac_dvar[_qp][ph][comp][v] += _fsp[ph].dmass_fraction_dz[comp] * dz_dvar[v];
247  }
248  }
249  }
250 
251  // If the material properties are being evaluated at the qps, calculate the
252  // gradients as well. Note: only nodal properties are evaluated in
253  // initQpStatefulProperties(), so no need to check _is_initqp flag for qp
254  // properties
255  if (!_nodal_material)
256  {
257  // Second derivative of capillary pressure
258  Real d2pc = _pc_uo.d2CapillaryPressure(_fsp[_aqueous_phase_number].saturation);
259 
260  (*_grads_qp)[_qp][_gas_phase_number] =
262  _dsaturation_dvar[_qp][_gas_phase_number][_zvar[0]] * (*_gradz_qp[0])[_qp];
263  (*_grads_qp)[_qp][_aqueous_phase_number] = -(*_grads_qp)[_qp][_gas_phase_number];
264 
265  (*_gradp_qp)[_qp][_gas_phase_number] = _gas_gradp_qp[_qp];
266  (*_gradp_qp)[_qp][_aqueous_phase_number] =
267  _gas_gradp_qp[_qp] - dpc * (*_grads_qp)[_qp][_aqueous_phase_number];
268 
269  (*_dgradp_qp_dv)[_qp][_aqueous_phase_number][_zvar[0]] =
270  -d2pc * (*_grads_qp)[_qp][_aqueous_phase_number];
271 
272  (*_grad_mass_frac_qp)[_qp][_aqueous_phase_number][_aqueous_fluid_component] =
274  _gas_gradp_qp[_qp] +
276  (*_gradz_qp[0])[_qp];
277  (*_grad_mass_frac_qp)[_qp][_aqueous_phase_number][_gas_fluid_component] =
278  -(*_grad_mass_frac_qp)[_qp][_aqueous_phase_number][_aqueous_fluid_component];
279  (*_grad_mass_frac_qp)[_qp][_gas_phase_number][_aqueous_fluid_component] =
280  _fsp[_gas_phase_number].dmass_fraction_dp[_aqueous_fluid_component] * _gas_gradp_qp[_qp] +
281  _fsp[_gas_phase_number].dmass_fraction_dz[_aqueous_fluid_component] * (*_gradz_qp[0])[_qp];
282  (*_grad_mass_frac_qp)[_qp][_gas_phase_number][_gas_fluid_component] =
283  -(*_grad_mass_frac_qp)[_qp][_gas_phase_number][_aqueous_fluid_component];
284  }
285 }
std::vector< const VariableGradient * > _gradz_qp
Gradient(s) of total mass fraction(s) of the gas component(s) (only defined at the qps) ...
void setMaterialVectorSize() const
Size material property vectors and initialise with zeros.
MaterialProperty< std::vector< std::vector< Real > > > & _dporepressure_dvar
d(porepressure)/d(PorousFlow variable)
const MaterialProperty< std::vector< Real > > & _dtemperature_dvar
Derivative of temperature wrt PorousFlow variables.
const VariableGradient & _gas_gradp_qp
Gradient of porepressure (only defined at the qps)
MaterialProperty< std::vector< std::vector< Real > > > & _dfluid_viscosity_dvar
Derivative of the fluid viscosity for each phase wrt PorousFlow variables.
const unsigned int _num_pf_vars
Number of PorousFlow variables.
const unsigned int _gas_porepressure_varnum
Moose variable number of the gas porepressure.
std::vector< unsigned int > _z_varnum
Moose variable number of z.
const unsigned int _pvar
PorousFlow variable number of the gas porepressure.
MaterialProperty< std::vector< Real > > & _saturation
Computed nodal or qp saturation of the phases.
MaterialProperty< std::vector< std::vector< Real > > > & _dfluid_density_dvar
Derivative of the fluid density for each phase wrt PorousFlow variables.
MaterialProperty< std::vector< std::vector< Real > > > & _dsaturation_dvar
d(saturation)/d(PorousFlow variable)
virtual Real dCapillaryPressure(Real saturation) const
Derivative of capillary pressure wrt true saturation.
MaterialProperty< std::vector< Real > > & _fluid_density
Fluid density of each phase.
virtual Real d2CapillaryPressure(Real saturation) const
Second derivative of capillary pressure wrt true saturation.
MaterialProperty< std::vector< std::vector< Real > > > & _mass_frac
Mass fraction matrix.
const unsigned int _gas_fluid_component
Fluid component number of the gas phase.
const PorousFlowCapillaryPressure & _pc_uo
Capillary pressure UserObject.
std::vector< FluidStateProperties > _fsp
FluidStateProperties data structure.
std::vector< unsigned int > _zvar
PorousFlow variable number of z.
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
virtual void computeQpProperties() override
const unsigned int _gas_phase_number
Phase number of the gas phase.
bool _is_initqp
Flag to indicate whether to calculate stateful properties.
const unsigned int _aqueous_fluid_component
Fluid component number of the aqueous component.
const unsigned int _num_components
Number of components.
const unsigned int _num_phases
Number of phases.
void FORTRAN_CALL() saturation(double &P, double &T, int &N, int &nerr)
MaterialProperty< std::vector< Real > > & _porepressure
Computed nodal or quadpoint values of porepressure of the phases.
MaterialProperty< std::vector< Real > > & _fluid_viscosity
Viscosity of each phase.
MaterialProperty< std::vector< std::vector< std::vector< Real > > > > & _dmass_frac_dvar
Derivative of the mass fraction matrix with respect to the Porous Flow variables. ...
virtual void thermophysicalProperties()=0
Calculates all required thermophysical properties and derivatives for each phase and fluid component...
void PorousFlowFluidStateBrineCO2::equilibriumConstantCO2 ( Real  temperature,
Real &  kco2,
Real &  dkco2_dT 
) const
protected

Equilibrium constant for CO2 from 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). Note: correlation uses temperature in Celcius

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

Definition at line 619 of file PorousFlowFluidStateBrineCO2.C.

Referenced by massFractions().

622 {
623  // Uses temperature in Celcius
624  Real Tc = temperature - _T_c2k;
625 
626  Real logK0CO2 = 1.189 + 1.304e-2 * Tc - 5.446e-5 * Tc * Tc;
627  Real dlogK0CO2 = 1.304e-2 - 1.0892e-4 * Tc;
628 
629  kco2 = std::pow(10.0, logK0CO2);
630  dkco2_dT = std::log(10.0) * dlogK0CO2 * kco2;
631 }
const std::string temperature
Definition: NS.h:25
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const Real _T_c2k
Conversion from degrees Celsius to degrees Kelvin.
void PorousFlowFluidStateBrineCO2::equilibriumConstantH2O ( Real  temperature,
Real &  kh2o,
Real &  dkh2o_dT 
) const
protected

Equilibrium constant for H2O from 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). Note: correlation uses temperature in Celcius

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

Definition at line 604 of file PorousFlowFluidStateBrineCO2.C.

Referenced by massFractions().

607 {
608  // Uses temperature in Celcius
609  Real Tc = temperature - _T_c2k;
610 
611  Real logK0H2O = -2.209 + 3.097e-2 * Tc - 1.098e-4 * Tc * Tc + 2.048e-7 * Tc * Tc * Tc;
612  Real dlogK0H2O = 3.097e-2 - 2.196e-4 * Tc + 6.144e-7 * Tc * Tc;
613 
614  kh2o = std::pow(10.0, logK0H2O);
615  dkh2o_dT = std::log(10.0) * dlogK0H2O * kh2o;
616 }
const std::string temperature
Definition: NS.h:25
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const Real _T_c2k
Conversion from degrees Celsius to degrees Kelvin.
void PorousFlowFluidStateBrineCO2::fugacityCoefficientCO2 ( Real  pressure,
Real  temperature,
Real &  fco2,
Real &  dfco2_dp,
Real &  dfco2_dT 
) const
protected

Fugacity coefficient for CO2.

Eq. (B7) from 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)

Parameters
pressuregas pressure (Pa)
temperaturetemperature (K)
[out]fco2fugacity coefficient for CO2
[out]dfco2_dpderivative of fugacity coefficient wrt pressure
[out]dfco2_dTderivative of fugacity coefficient wrt temperature

Definition at line 466 of file PorousFlowFluidStateBrineCO2.C.

Referenced by massFractions().

468 {
469  // Need pressure in bar
470  Real pbar = pressure * 1.0e-5;
471  // CO2 density and derivatives wrt pressure and temperature
472  Real gas_density, dgas_density_dp, dgas_density_dT;
473  _co2_fp.rho_dpT(pressure, temperature, gas_density, dgas_density_dp, dgas_density_dT);
474  // Molar volume in cm^3/mol
475  Real V = _Mco2 / gas_density * 1.0e6;
476 
477  // Redlich-Kwong parameters
478  Real aCO2 = 7.54e7 - 4.13e4 * temperature;
479  Real bCO2 = 27.8;
480 
481  Real term1 = std::log(V / (V - bCO2)) + bCO2 / (V - bCO2) -
482  2.0 * aCO2 / (_Rbar * std::pow(temperature, 1.5) * bCO2) * std::log((V + bCO2) / V) +
483  aCO2 / (_Rbar * std::pow(temperature, 1.5) * bCO2) *
484  (std::log((V + bCO2) / V) - bCO2 / (V + bCO2));
485 
486  Real lnPhiCO2 = term1 - std::log(pbar * V / (_Rbar * temperature));
487  fco2 = std::exp(lnPhiCO2);
488 
489  // The derivative of the fugacity coefficient wrt pressure
490  Real dV_dp = -_Mco2 / gas_density / gas_density * dgas_density_dp * 1.0e6;
491  Real dterm1_dV =
492  (bCO2 * (bCO2 - 2.0 * V) / (bCO2 - V) / (bCO2 - V) +
493  aCO2 * (bCO2 + 2.0 * V) / (_Rbar * std::pow(temperature, 1.5) * (bCO2 + V) * (bCO2 + V))) /
494  V;
495  dfco2_dp = (dterm1_dV * dV_dp - 1.0e-5 / pbar - 1.0 / V * dV_dp) * fco2;
496 
497  // The derivative of the fugacity coefficient wrt temperature
498  Real dV_dT = -_Mco2 / gas_density / gas_density * dgas_density_dT * 1.0e6;
499  Real dterm1_dT = 3.0 * aCO2 * _Rbar * std::sqrt(temperature) * bCO2 * std::log((V + bCO2) / V) /
500  (_Rbar * std::pow(temperature, 1.5) * bCO2) /
501  (_Rbar * std::pow(temperature, 1.5) * bCO2) +
502  8.26e4 / (_Rbar * std::pow(temperature, 1.5) * bCO2) * std::log((V + bCO2) / V) -
503  1.5 * aCO2 * _Rbar * std::sqrt(temperature) * bCO2 *
504  (std::log((V + bCO2) / V) - bCO2 / (V + bCO2)) /
505  (_Rbar * std::pow(temperature, 1.5) * bCO2) /
506  (_Rbar * std::pow(temperature, 1.5) * bCO2) -
507  4.13e4 / (_Rbar * std::pow(temperature, 1.5) * bCO2) *
508  (std::log((V + bCO2) / V) - bCO2 / (V + bCO2));
509  dfco2_dT = (dterm1_dT + dterm1_dV * dV_dT - dV_dT / V + 1.0 / temperature) * fco2;
510 }
virtual void rho_dpT(Real pressure, Real temperature, Real &rho, Real &drho_dp, Real &drho_dT) const =0
Density and its derivatives wrt pressure and temperature.
const std::string temperature
Definition: NS.h:25
const SinglePhaseFluidPropertiesPT & _co2_fp
Fluid properties UserObject for CO2.
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const std::string pressure
Definition: NS.h:24
const Real _Mco2
Molar mass of CO2 (kg/mol)
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
void PorousFlowFluidStateBrineCO2::fugacityCoefficientH2O ( Real  pressure,
Real  temperature,
Real &  fh2o,
Real &  dfh2o_dp,
Real &  dfh2o_dT 
) const
protected

Fugacity coefficient for H2O.

Eq. (B7) from 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)

Parameters
pressuregas pressure (Pa)
temperaturetemperature (K)
[out]fh2ofugacity coefficient for H2O
[out]dfh2o_dpderivative of fugacity coefficient wrt pressure
[out]dfh2o_dTderivative of fugacity coefficient wrt temperature

Definition at line 513 of file PorousFlowFluidStateBrineCO2.C.

Referenced by massFractions().

515 {
516  // Need pressure in bar
517  Real pbar = pressure * 1.0e-5;
518  // CO2 density and derivatives wrt pressure and temperature
519  Real gas_density, dgas_density_dp, dgas_density_dT;
520  _co2_fp.rho_dpT(pressure, temperature, gas_density, dgas_density_dp, dgas_density_dT);
521  // Molar volume in cm^3/mol
522  Real V = _Mco2 / gas_density * 1.0e6;
523 
524  // Redlich-Kwong parameters
525  Real aCO2 = 7.54e7 - 4.13e4 * temperature;
526  Real bCO2 = 27.8;
527  Real aCO2H2O = 7.89e7;
528  Real bH2O = 18.18;
529 
530  Real term1 =
531  std::log(V / (V - bCO2)) + bH2O / (V - bCO2) -
532  2.0 * aCO2H2O / (_Rbar * std::pow(temperature, 1.5) * bCO2) * std::log((V + bCO2) / V) +
533  aCO2 * bH2O / (_Rbar * std::pow(temperature, 1.5) * bCO2 * bCO2) *
534  (std::log((V + bCO2) / V) - bCO2 / (V + bCO2));
535 
536  Real lnPhiH2O = term1 - std::log(pbar * V / (_Rbar * temperature));
537  fh2o = std::exp(lnPhiH2O);
538 
539  // The derivative of the fugacity coefficient wrt pressure
540  Real dV_dp = -_Mco2 / gas_density / gas_density * dgas_density_dp * 1.0e6;
541  Real dterm1_dV = ((bCO2 * bCO2 - (bCO2 + bH2O) * V) / (bCO2 - V) / (bCO2 - V) +
542  (2.0 * aCO2H2O * (bCO2 + V) - aCO2 * bH2O) /
543  (_Rbar * std::pow(temperature, 1.5) * (bCO2 + V) * (bCO2 + V))) /
544  V;
545  dfh2o_dp = (dterm1_dV * dV_dp - 1.0e-5 / pbar - dV_dp / V) * fh2o;
546 
547  // The derivative of the fugacity coefficient wrt temperature
548  Real dV_dT = -_Mco2 / gas_density / gas_density * dgas_density_dT * 1.0e6;
549  Real dterm1_dT = 3.0 * _Rbar * std::sqrt(temperature) * bCO2 * aCO2H2O *
550  std::log((V + bCO2) / V) / (_Rbar * std::pow(temperature, 1.5) * bCO2) /
551  (_Rbar * std::pow(temperature, 1.5) * bCO2) -
552  1.5 * aCO2 * bH2O * _Rbar * std::sqrt(temperature) * bCO2 * bCO2 *
553  (std::log((V + bCO2) / V) - bCO2 / (V + bCO2)) /
554  (_Rbar * std::pow(temperature, 1.5) * bCO2 * bCO2) /
555  (_Rbar * std::pow(temperature, 1.5) * bCO2 * bCO2) -
556  4.13e4 * bH2O * (std::log((V + bCO2) / V) - bCO2 / (V + bCO2)) /
557  (_Rbar * std::pow(temperature, 1.5) * bCO2 * bCO2);
558  dfh2o_dT = (dterm1_dT + dterm1_dV * dV_dT - dV_dT / V + 1.0 / temperature) * fh2o;
559 }
virtual void rho_dpT(Real pressure, Real temperature, Real &rho, Real &drho_dp, Real &drho_dT) const =0
Density and its derivatives wrt pressure and temperature.
const std::string temperature
Definition: NS.h:25
const SinglePhaseFluidPropertiesPT & _co2_fp
Fluid properties UserObject for CO2.
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const std::string pressure
Definition: NS.h:24
const Real _Mco2
Molar mass of CO2 (kg/mol)
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
void PorousFlowFluidStateFlashBase::initQpStatefulProperties ( )
overrideprotectedvirtualinherited

Reimplemented from PorousFlowVariableBase.

Definition at line 134 of file PorousFlowFluidStateFlashBase.C.

135 {
136  _is_initqp = true;
137  // Set the size of pressure and saturation vectors
139 
140  // Set the size of all other vectors
142 
143  // Set the initial values of the properties at the nodes.
144  // Note: not required for qp materials as no old values at the qps are requested
145  if (_nodal_material)
146  {
148 
149  for (unsigned int ph = 0; ph < _num_phases; ++ph)
150  {
151  _saturation[_qp][ph] = _fsp[ph].saturation;
152  _porepressure[_qp][ph] = _fsp[ph].pressure;
153  _fluid_density[_qp][ph] = _fsp[ph].fluid_density;
154  _fluid_viscosity[_qp][ph] = _fsp[ph].fluid_viscosity;
155  _mass_frac[_qp][ph] = _fsp[ph].mass_fraction;
156  }
157  }
158 }
void setMaterialVectorSize() const
Size material property vectors and initialise with zeros.
MaterialProperty< std::vector< Real > > & _saturation
Computed nodal or qp saturation of the phases.
MaterialProperty< std::vector< Real > > & _fluid_density
Fluid density of each phase.
MaterialProperty< std::vector< std::vector< Real > > > & _mass_frac
Mass fraction matrix.
virtual void initQpStatefulProperties() override
std::vector< FluidStateProperties > _fsp
FluidStateProperties data structure.
bool _is_initqp
Flag to indicate whether to calculate stateful properties.
const unsigned int _num_phases
Number of phases.
MaterialProperty< std::vector< Real > > & _porepressure
Computed nodal or quadpoint values of porepressure of the phases.
MaterialProperty< std::vector< Real > > & _fluid_viscosity
Viscosity of each phase.
virtual void thermophysicalProperties()=0
Calculates all required thermophysical properties and derivatives for each phase and fluid component...
void PorousFlowFluidStateBrineCO2::massFractions ( Real  pressure,
Real  temperature,
Real  xnacl,
Real &  xco2l,
Real &  dxco2l_dp,
Real &  dxco2l_dT,
Real &  xh2og,
Real &  dxh2og_dp,
Real &  dxh2og_dT 
) const
protected

Mass fractions of CO2 and brine calculated using mutual solubilities given by 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)

Parameters
pressurephase pressure (Pa)
temperaturephase temperature (C)
xnaclNaCl mass fraction (kg/kg)
[out]xco2lmass fraction of CO2 in liquid (kg/kg)
[out]dxco2l_dpderivative of mass fraction of CO2 in liquid wrt pressure
[out]dxco2l_dTderivative of mass fraction of CO2 in liqiud wrt temperature
[out]xh2ogmass fraction of H2O in gas (kg/kg)
[out]dh2ogl_dpderivative of mass fraction of H2O in gas wrt pressure
[out]dh2og_dTderivative of mass fraction of H2O in gas wrt temperature

Definition at line 341 of file PorousFlowFluidStateBrineCO2.C.

Referenced by thermophysicalProperties().

350 {
351  // Pressure in bar
352  Real pbar = pressure * 1.0e-5;
353  // Pressure minus 1 bar
354  Real delta_pbar = pbar - 1.0;
355 
356  // Average partial molar volumes (cm^3/mol) as given by Sypcher, Pruess and Ennis-King (2003)
357  const Real vCO2 = 32.6;
358  const Real vH2O = 18.1;
359 
360  // NaCl molality (mol/kg)
361  Real mnacl = xnacl / (1.0 - xnacl) / _Mnacl;
362 
363  // Equilibrium constants
364  Real K0H2O, dK0H2O_dT, K0CO2, dK0CO2_dT;
365  equilibriumConstantH2O(temperature, K0H2O, dK0H2O_dT);
366  equilibriumConstantCO2(temperature, K0CO2, dK0CO2_dT);
367 
368  // Fugacity coefficients
369  Real phiH2O, dphiH2O_dp, dphiH2O_dT;
370  Real phiCO2, dphiCO2_dp, dphiCO2_dT;
371  fugacityCoefficientCO2(pressure, temperature, phiCO2, dphiCO2_dp, dphiCO2_dT);
372  fugacityCoefficientH2O(pressure, temperature, phiH2O, dphiH2O_dp, dphiH2O_dT);
373 
374  // Activity coefficient
375  Real gamma, dgamma_dp, dgamma_dT;
376  activityCoefficient(pressure, temperature, xnacl, gamma, dgamma_dp, dgamma_dT);
377 
378  Real A = K0H2O / (phiH2O * pbar) * std::exp(delta_pbar * vH2O / (_Rbar * temperature));
379  Real B =
380  phiCO2 * pbar / (_invMh2o * K0CO2) * std::exp(-delta_pbar * vCO2 / (_Rbar * temperature));
381 
382  Real dA_dp = (-1.0e-5 * K0H2O / pbar + 1.0e-5 * vH2O * K0H2O / (_Rbar * temperature) -
383  K0H2O * dphiH2O_dp / phiH2O) *
384  std::exp(delta_pbar * vH2O / (_Rbar * temperature)) / (pbar * phiH2O);
385  Real dB_dp = (1.0e-5 * phiCO2 + pbar * dphiCO2_dp -
386  1.0e-5 * vCO2 * pbar * phiCO2 / (_Rbar * temperature)) *
387  std::exp(-delta_pbar * vCO2 / (_Rbar * temperature)) / (_invMh2o * K0CO2);
388 
389  Real dA_dT = (dK0H2O_dT - dphiH2O_dT * K0H2O / phiH2O -
390  delta_pbar * vH2O * K0H2O / (_Rbar * temperature * temperature)) *
391  std::exp(delta_pbar * vH2O / (_Rbar * temperature)) / (pbar * phiH2O);
392  Real dB_dT = (-pbar * phiCO2 * dK0CO2_dT / K0CO2 + pbar * dphiCO2_dT +
393  delta_pbar * vCO2 * pbar * phiCO2 / (_Rbar * temperature * temperature)) *
394  std::exp(-delta_pbar * vCO2 / (_Rbar * temperature)) / (_invMh2o * K0CO2);
395 
396  // The mole fraction of H2O in the CO2-rich gas phase is then
397  Real yH2O = (1.0 - B) / (1.0 / A - B);
398  // The mole fraction of CO2 in the H2O-rich liquid phase is then (note: no salinty effect)
399  Real xCO2 = B * (1.0 - yH2O);
400  // The molality of CO2 in the H2O-rich liquid phase is (note: no salinty effect)
401  Real mCO2 = xCO2 * _invMh2o / (1.0 - xCO2);
402 
403  // The molality of CO2 in brine is then given by
404  Real mCO2b = mCO2 / gamma;
405  // The mole fraction of CO2 in brine is then
406  Real xCO2b = mCO2b / (2.0 * mnacl + _invMh2o + mCO2b);
407  // The mole fraction of H2O in the CO2-rich gas phase corrected for NaCl mole fraction is
408  Real yH2Ob = A * (1.0 - xCO2b - 2.0 * mnacl / (2.0 * mnacl + _invMh2o + mCO2b));
409 
410  // Convert the mole fractions to mass fractions and then update referenced values
411  // The mass fraction of H2O in gas (assume no salt in gas phase)
412  xh2og = yH2Ob * _Mh2o / (yH2Ob * _Mh2o + (1.0 - yH2Ob) * _Mco2);
413 
414  // The number of moles of CO2 in 1kg of H2O
415  Real nco2 = xCO2b * (2.0 * mnacl + _invMh2o) / (1.0 - xCO2b);
416 
417  // The mass fraction of CO2 in liquid
418  xco2l = nco2 * _Mco2 / (1.0 + mnacl * _Mnacl + nco2 * _Mco2);
419 
420  // The derivatives of the mass fractions wrt pressure
421  Real dyH2O_dp = ((1.0 - B) * dA_dp + (A - 1.0) * A * dB_dp) / (1.0 - A * B) / (1.0 - A * B);
422  Real dxCO2_dp = dB_dp * (1.0 - yH2O) - B * dyH2O_dp;
423 
424  Real dmCO2_dp = _invMh2o * dxCO2_dp / (1.0 - xCO2) / (1.0 - xCO2);
425  Real dmCO2b_dp = dmCO2_dp / gamma - mCO2 * dgamma_dp / gamma / gamma;
426  Real dxCO2b_dp = (2.0 * mnacl + _invMh2o) * dmCO2b_dp / (2.0 * mnacl + _invMh2o + mCO2b) /
427  (2.0 * mnacl + _invMh2o + mCO2b);
428 
429  Real dyH2Ob_dp = (1.0 - xCO2b - 2.0 * mnacl / (2.0 * mnacl + _invMh2o + mCO2b)) * dA_dp -
430  A * dxCO2b_dp +
431  2.0 * A * mnacl * dmCO2b_dp / (2.0 * mnacl + _invMh2o + mCO2b) /
432  (2.0 * mnacl + _invMh2o + mCO2b);
433 
434  dxh2og_dp = _Mco2 * _Mh2o * dyH2Ob_dp / (yH2Ob * _Mh2o + (1.0 - yH2Ob) * _Mco2) /
435  (yH2Ob * _Mh2o + (1.0 - yH2Ob) * _Mco2);
436 
437  Real dnco2_dp = dxCO2b_dp * (2.0 * mnacl + _invMh2o) / (1.0 - xCO2b) / (1.0 - xCO2b);
438 
439  dxco2l_dp = (1.0 + mnacl * _Mnacl) * _Mco2 * dnco2_dp / (1.0 + mnacl * _Mnacl + nco2 * _Mco2) /
440  (1.0 + mnacl * _Mnacl + nco2 * _Mco2);
441 
442  // The derivatives of the mass fractions wrt temperature
443  Real dyH2O_dT = ((1.0 - B) * dA_dT + (A - 1.0) * A * dB_dT) / (1.0 - A * B) / (1.0 - A * B);
444  Real dxCO2_dT = dB_dT * (1.0 - yH2O) - B * dyH2O_dT;
445 
446  Real dmCO2_dT = _invMh2o * dxCO2_dT / (1.0 - xCO2) / (1.0 - xCO2);
447  Real dmCO2b_dT = dmCO2_dT / gamma - mCO2 * dgamma_dT / gamma / gamma;
448  Real dxCO2b_dT = (2.0 * mnacl + _invMh2o) * dmCO2b_dT / (2.0 * mnacl + _invMh2o + mCO2b) /
449  (2.0 * mnacl + _invMh2o + mCO2b);
450 
451  Real dyH2Ob_dT = (1.0 - xCO2b - 2.0 * mnacl / (2.0 * mnacl + _invMh2o + mCO2b)) * dA_dT -
452  A * dxCO2b_dT +
453  2.0 * A * mnacl * dmCO2b_dT / (2.0 * mnacl + _invMh2o + mCO2b) /
454  (2.0 * mnacl + _invMh2o + mCO2b);
455 
456  dxh2og_dT = _Mco2 * _Mh2o * dyH2Ob_dT / (yH2Ob * _Mh2o + (1.0 - yH2Ob) * _Mco2) /
457  (yH2Ob * _Mh2o + (1.0 - yH2Ob) * _Mco2);
458 
459  Real dnco2_dT = dxCO2b_dT * (2.0 * mnacl + _invMh2o) / (1.0 - xCO2b) / (1.0 - xCO2b);
460 
461  dxco2l_dT = (1.0 + mnacl * _Mnacl) * _Mco2 * dnco2_dT / (1.0 + mnacl * _Mnacl + nco2 * _Mco2) /
462  (1.0 + mnacl * _Mnacl + nco2 * _Mco2);
463 }
void fugacityCoefficientCO2(Real pressure, Real temperature, Real &fco2, Real &dfco2_dp, Real &dfco2_dT) const
Fugacity coefficient for CO2.
const Real _Mnacl
Molar mass of NaCL.
void fugacityCoefficientH2O(Real pressure, Real temperature, Real &fh2o, Real &dfh2o_dp, Real &dfh2o_dT) const
Fugacity coefficient for H2O.
const std::string temperature
Definition: NS.h:25
void activityCoefficient(Real pressure, Real temperature, Real xnacl, Real &gamma, Real &dgamma_dp, Real &dgamma_dT) const
Activity coefficient for CO2 in brine.
const Real _Mh2o
Molar mass of water (kg/mol)
const Real _invMh2o
Inverse of molar mass of H2O (mol/kg)
const std::string pressure
Definition: NS.h:24
const Real _Mco2
Molar mass of CO2 (kg/mol)
void equilibriumConstantH2O(Real temperature, Real &kh2o, Real &dkh2o_dT) const
Equilibrium constant for H2O from Spycher, Pruess and Ennis-King, CO2-H2O mixtures in the geological ...
const Real _Rbar
Molar gas constant in bar cm^3 /(K mol)
void equilibriumConstantCO2(Real temperature, Real &kco2, Real &dkco2_dT) const
Equilibrium constant for CO2 from Spycher, Pruess and Ennis-King, CO2-H2O mixtures in the geological ...
void PorousFlowFluidStateBrineCO2::partialDensityCO2 ( Real  temperature,
Real &  partial_density,
Real &  dpartial_density_dT 
) const
protected

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 634 of file PorousFlowFluidStateBrineCO2.C.

Referenced by thermophysicalProperties().

637 {
638  // This correlation uses temperature in C
639  Real Tc = temperature - _T_c2k;
640  // The parial molar volume
641  Real V = 37.51 - 9.585e-2 * Tc + 8.74e-4 * Tc * Tc - 5.044e-7 * Tc * Tc * Tc;
642  Real dV_dT = -9.585e-2 + 1.748e-3 * Tc - 1.5132e-6 * Tc * Tc;
643 
644  partial_density = 1.0e6 * _Mco2 / V;
645  dpartial_density_dT = -1.0e6 * _Mco2 * dV_dT / V / V;
646 }
const std::string temperature
Definition: NS.h:25
const Real _T_c2k
Conversion from degrees Celsius to degrees Kelvin.
const Real _Mco2
Molar mass of CO2 (kg/mol)
Real PorousFlowFluidStateFlashBase::rachfordRice ( Real  x,
std::vector< Real > &  Ki 
) const
protectedinherited

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
xvapor fraction
Kiequilibrium constants
Returns
f(x)

Definition at line 320 of file PorousFlowFluidStateFlashBase.C.

Referenced by PorousFlowFluidStateFlashBase::vaporMassFraction().

321 {
322  // Check that the size of the equilibrium constant vector is correct
323  if (Ki.size() != _num_components)
324  mooseError("The number of equilibrium components passed to rachfordRice is not correct");
325 
326  Real f = 0.0;
327  Real z_total = 0.0;
328 
329  for (unsigned int i = 0; i < _num_z_vars; ++i)
330  {
331  f += (*_z[i])[_qp] * (Ki[i] - 1.0) / (1.0 + x * (Ki[i] - 1.0));
332  z_total += (*_z[i])[_qp];
333  }
334 
335  // Add the last component (with total mass fraction = 1 - z_total)
336  f += (1.0 - z_total) * (Ki[_num_z_vars] - 1.0) / (1.0 + x * (Ki[_num_z_vars] - 1.0));
337 
338  return f;
339 }
const unsigned int _num_z_vars
Number of coupled total mass fractions. Should be _num_phases - 1.
const unsigned int _num_components
Number of components.
std::vector< const VariableValue * > _z
Total mass fraction(s) of the gas component(s) summed over all phases.
Real PorousFlowFluidStateFlashBase::rachfordRiceDeriv ( Real  x,
std::vector< Real > &  Ki 
) const
protectedinherited

Derivative of Rachford-Rice equation wrt vapor fraction.

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

Parameters
xvapor fraction
Kiequilibrium constants
Returns
f'(x)

Definition at line 342 of file PorousFlowFluidStateFlashBase.C.

Referenced by PorousFlowFluidStateFlashBase::vaporMassFraction().

343 {
344  // Check that the size of the equilibrium constant vector is correct
345  if (Ki.size() != _num_components)
346  mooseError("The number of equilibrium components passed to rachfordRiceDeriv is not correct");
347 
348  Real df = 0.0;
349  Real z_total = 0.0;
350 
351  for (unsigned int i = 0; i < _num_z_vars; ++i)
352  {
353  df -= (*_z[i])[_qp] * (Ki[i] - 1.0) * (Ki[i] - 1.0) / (1.0 + x * (Ki[i] - 1.0)) /
354  (1.0 + x * (Ki[i] - 1.0));
355  z_total += (*_z[i])[_qp];
356  }
357 
358  // Add the last component (with total mass fraction = 1 - z_total)
359  df -= (1.0 - z_total) * (Ki[_num_z_vars] - 1.0) * (Ki[_num_z_vars] - 1.0) /
360  (1.0 + x * (Ki[_num_z_vars] - 1.0)) / (1.0 + x * (Ki[_num_z_vars] - 1.0));
361 
362  return df;
363 }
const unsigned int _num_z_vars
Number of coupled total mass fractions. Should be _num_phases - 1.
const unsigned int _num_components
Number of components.
std::vector< const VariableValue * > _z
Total mass fraction(s) of the gas component(s) summed over all phases.
void PorousFlowFluidStateFlashBase::setMaterialVectorSize ( ) const
protectedinherited

Size material property vectors and initialise with zeros.

Definition at line 288 of file PorousFlowFluidStateFlashBase.C.

Referenced by PorousFlowFluidStateFlashBase::computeQpProperties(), and PorousFlowFluidStateFlashBase::initQpStatefulProperties().

289 {
290  _fluid_density[_qp].assign(_num_phases, 0.0);
291  _fluid_viscosity[_qp].assign(_num_phases, 0.0);
292  _mass_frac[_qp].resize(_num_phases);
293 
294  // Derivatives and gradients are not required in initQpStatefulProperties
295  if (!_is_initqp)
296  {
297  _dfluid_density_dvar[_qp].resize(_num_phases);
298  _dfluid_viscosity_dvar[_qp].resize(_num_phases);
299  _dmass_frac_dvar[_qp].resize(_num_phases);
300 
301  if (!_nodal_material)
302  (*_grad_mass_frac_qp)[_qp].resize(_num_phases);
303 
304  for (unsigned int ph = 0; ph < _num_phases; ++ph)
305  {
306  _dfluid_density_dvar[_qp][ph].assign(_num_pf_vars, 0.0);
307  _dfluid_viscosity_dvar[_qp][ph].assign(_num_pf_vars, 0.0);
308  _dmass_frac_dvar[_qp][ph].resize(_num_components);
309 
310  for (unsigned int comp = 0; comp < _num_components; ++comp)
311  _dmass_frac_dvar[_qp][ph][comp].assign(_num_pf_vars, 0.0);
312 
313  if (!_nodal_material)
314  (*_grad_mass_frac_qp)[_qp][ph].assign(_num_components, RealGradient());
315  }
316  }
317 }
MaterialProperty< std::vector< std::vector< Real > > > & _dfluid_viscosity_dvar
Derivative of the fluid viscosity for each phase wrt PorousFlow variables.
const unsigned int _num_pf_vars
Number of PorousFlow variables.
MaterialProperty< std::vector< std::vector< Real > > > & _dfluid_density_dvar
Derivative of the fluid density for each phase wrt PorousFlow variables.
MaterialProperty< std::vector< Real > > & _fluid_density
Fluid density of each phase.
MaterialProperty< std::vector< std::vector< Real > > > & _mass_frac
Mass fraction matrix.
bool _is_initqp
Flag to indicate whether to calculate stateful properties.
const unsigned int _num_components
Number of components.
const unsigned int _num_phases
Number of phases.
MaterialProperty< std::vector< Real > > & _fluid_viscosity
Viscosity of each phase.
MaterialProperty< std::vector< std::vector< std::vector< Real > > > > & _dmass_frac_dvar
Derivative of the mass fraction matrix with respect to the Porous Flow variables. ...
void PorousFlowFluidStateBrineCO2::thermophysicalProperties ( )
overrideprotectedvirtual

Calculates all required thermophysical properties and derivatives for each phase and fluid component.

Must override in all derived classes.

Implements PorousFlowFluidStateFlashBase.

Definition at line 47 of file PorousFlowFluidStateBrineCO2.C.

48 {
49  // The FluidProperty objects use temperature in K
50  Real Tk = _temperature[_qp] + _T_c2k;
51  Real pressure = _gas_porepressure[_qp];
52  Real xnacl = _xnacl[_qp];
53 
54  // Mass fraction of CO2 in liquid and H2O in gas phases
55  Real X0, dX0_dp, dX0_dT, Y0, Y1, dY1_dp, dY1_dT;
56  massFractions(pressure, Tk, xnacl, X0, dX0_dp, dX0_dT, Y1, dY1_dp, dY1_dT);
57  Y0 = 1.0 - Y1;
58 
59  // Determine which phases are present based on the value of z
60  bool is_liquid = false;
61  bool is_gas = false;
62  bool is_twophase = false;
63 
64  if ((*_z[0])[_qp] <= X0)
65  {
66  // In this case, there is not enough CO2 to form a gas phase,
67  // so only a liquid phase is present
68  is_liquid = true;
69  }
70  else if ((*_z[0])[_qp] > X0 && (*_z[0])[_qp] < Y0)
71  {
72  // Two phases are present
73  is_twophase = true;
74  }
75  else // ((*_z[0])[_qp] >= Y0)
76  {
77  // In this case, there is not enough brine to form a liquid
78  // phase, so only a gas phase is present
79  is_gas = true;
80  }
81 
82  // Material must provide the following properties
83  Real vapor_mass_fraction = 0.0;
84  Real gas_saturation, liquid_saturation;
85  Real gas_density, liquid_density;
86  Real gas_viscosity, liquid_viscosity;
87 
88  // And the following derivatives
89  Real dX0_dz = 0.0, dY0_dp = 0.0, dY0_dT = 0.0, dY0_dz = 0.0;
90  Real dgas_density_dp = 0.0, dgas_density_dT = 0.0, dgas_density_dz = 0.0;
91  Real dliquid_density_dp = 0.0, dliquid_density_dT = 0.0, dliquid_density_dz = 0.0;
92  Real dgas_viscosity_dp = 0.0, dgas_viscosity_dT = 0.0, dgas_viscosity_dz = 0.0;
93  Real dliquid_viscosity_dp = 0.0, dliquid_viscosity_dT = 0.0, dliquid_viscosity_dz = 0.0;
94 
95  if (is_liquid)
96  {
97  X0 = (*_z[0])[_qp];
98  Y0 = 0.0;
99  gas_saturation = 0.0;
100  dX0_dp = 0.0;
101  dX0_dT = 0.0;
102  dX0_dz = 1.0;
103  gas_density = 0.0;
104  gas_viscosity = 1.0; // To guard against division by 0
105  }
106  else if (is_twophase)
107  {
108  // Two phase are present. Set mass equilibrium constants used in the
109  // calculation of vapor mass fraction
110  Y0 = 1.0 - Y1;
111  std::vector<Real> Ki(_num_components);
112  Ki[0] = Y0 / X0;
113  Ki[1] = (1.0 - Y0) / (1.0 - X0);
114  vapor_mass_fraction = vaporMassFraction(Ki);
115 
116  // Derivatives of mass fractions wrt PorousFlow variables
117  dY0_dp = -dY1_dp;
118  dY0_dT = -dY1_dT;
119  }
120  else // if (is_gas)
121  {
122  X0 = 0.0;
123  dX0_dp = 0.0;
124  dX0_dT = 0.0;
125  Y0 = (*_z[0])[_qp];
126  dY0_dz = 1.0;
127  vapor_mass_fraction = 1.0;
128  gas_saturation = 1.0;
129  liquid_density = 0.0;
130  liquid_viscosity = 1.0; // To guard against division by 0
131  }
132 
133  // Update all remaining mass fractions and derivatives
134  // Mass fraction of CO2 in gas and H2O in liquid phases
135  Real X1 = 1.0 - X0;
136  Y1 = 1.0 - Y0;
137 
138  // Calculate the gas density and viscosity in the single phase gas or two
139  // phase region
140  if (is_gas || is_twophase)
141  {
142  Real co2_density, dco2_density_dp, dco2_density_dT;
143  _co2_fp.rho_dpT(pressure, Tk, co2_density, dco2_density_dp, dco2_density_dT);
144 
145  // Gas density is given by the CO2 density - no correction due to the small amount of
146  // brine vapor is made
147  gas_density = co2_density;
148  dgas_density_dp = dco2_density_dp;
149  dgas_density_dT = dco2_density_dT;
150  dgas_density_dz = 0.0;
151 
152  Real co2_viscosity, dco2_viscosity_drho, dco2_viscosity_dT;
154  co2_density, Tk, dco2_density_dT, co2_viscosity, dco2_viscosity_drho, dco2_viscosity_dT);
155 
156  // Assume that the viscosity of the gas phase is a weighted sum of the
157  // individual viscosities
158  gas_viscosity = co2_viscosity;
159  dgas_viscosity_dp = dco2_viscosity_drho * dco2_density_dp;
160  dgas_viscosity_dT = dco2_viscosity_dT;
161  dgas_viscosity_dz = 0.0;
162  }
163 
164  // Calculate the saturation in the two phase case using the vapor mass
165  // fraction
166  if (is_twophase)
167  {
168  // Liquid density
169  Real co2_partial_density, dco2_partial_density_dT;
170  partialDensityCO2(Tk, co2_partial_density, dco2_partial_density_dT);
171  // Use old value of gas saturation to estimate liquid saturation
172  Real liqsat = 1.0;
173  if (!_is_initqp)
174  liqsat -= _saturation_old[_qp][_gas_phase_number];
175 
176  liquid_density =
177  1.0 / (X0 / co2_partial_density +
178  X1 / _brine_fp.rho(pressure - _pc_uo.capillaryPressure(liqsat), Tk, xnacl));
179 
180  // The gas saturation in the two phase case
181  gas_saturation = vapor_mass_fraction * liquid_density /
182  (gas_density + vapor_mass_fraction * (liquid_density - gas_density));
183  }
184 
185  // Calculate the saturations and pressures for each phase
186  liquid_saturation = 1.0 - gas_saturation;
187  Real liquid_porepressure = pressure - _pc_uo.capillaryPressure(liquid_saturation);
188 
189  // Calculate liquid density and viscosity if in the two phase or single phase
190  // liquid region, including a density correction due to the presence of dissolved
191  // CO2
192  if (is_liquid || is_twophase)
193  {
194  Real brine_density, dbrine_density_dp, dbrine_density_dT, dbrine_density_dx;
195  Real dliquid_viscosity_drho, dliquid_viscosity_dx;
196 
197  _brine_fp.rho_dpTx(liquid_porepressure,
198  Tk,
199  xnacl,
200  brine_density,
201  dbrine_density_dp,
202  dbrine_density_dT,
203  dbrine_density_dx);
204 
205  // The liquid density
206  Real co2_partial_density, dco2_partial_density_dT;
207  partialDensityCO2(Tk, co2_partial_density, dco2_partial_density_dT);
208  liquid_density = 1.0 / (X0 / co2_partial_density + X1 / brine_density);
209 
210  dliquid_density_dp =
211  (dX0_dp / brine_density + X1 * dbrine_density_dp / brine_density / brine_density -
212  dX0_dp / co2_partial_density) *
213  liquid_density * liquid_density;
214 
215  dliquid_density_dT =
216  (dX0_dT / brine_density + X1 * dbrine_density_dT / brine_density / brine_density -
217  dX0_dT / co2_partial_density +
218  X0 * dco2_partial_density_dT / co2_partial_density / co2_partial_density) *
219  liquid_density * liquid_density;
220 
221  dliquid_density_dz =
222  (dX0_dz / brine_density - dX0_dz / co2_partial_density) * liquid_density * liquid_density;
223 
224  // Liquid viscosity is just the brine viscosity
225  // Note: brine viscosity (and derivatives) requires water density (and derivatives)
226  Real water_density, dwater_density_dp, dwater_density_dT;
227  _water_fp.rho_dpT(liquid_porepressure, Tk, water_density, dwater_density_dp, dwater_density_dT);
228 
229  _brine_fp.mu_drhoTx(water_density,
230  Tk,
231  xnacl,
232  dwater_density_dT,
233  liquid_viscosity,
234  dliquid_viscosity_drho,
235  dliquid_viscosity_dT,
236  dliquid_viscosity_dx);
237 
238  // The derivative of viscosity wrt pressure is given by the chain rule
239  dliquid_viscosity_dp = dliquid_viscosity_drho * dwater_density_dp;
240  }
241 
242  // Save properties in FluidStateProperties vector
243  auto & liquid = _fsp[_aqueous_phase_number];
244  auto & gas = _fsp[_gas_phase_number];
245 
246  liquid.saturation = liquid_saturation;
247  gas.saturation = gas_saturation;
248  liquid.pressure = liquid_porepressure;
249  gas.pressure = _gas_porepressure[_qp];
250 
251  liquid.mass_fraction[_aqueous_fluid_component] = 1.0 - X0;
252  liquid.mass_fraction[_gas_fluid_component] = X0;
253  gas.mass_fraction[_aqueous_fluid_component] = 1.0 - Y0;
254  gas.mass_fraction[_gas_fluid_component] = Y0;
255 
256  liquid.fluid_density = liquid_density;
257  gas.fluid_density = gas_density;
258  liquid.fluid_viscosity = liquid_viscosity;
259  gas.fluid_viscosity = gas_viscosity;
260 
261  // Derivatives wrt PorousFlow variables are required by the kernels
262  // Note: these don't need to be stateful so don't calculate them in
263  // initQpStatefulProperties
264  if (!_is_initqp)
265  {
266  // Derivative of gas saturation wrt variables
267  Real ds_dp = 0.0, ds_dT = 0.0, ds_dz = 0.0;
268  if (is_twophase)
269  {
270  Real K0 = Y0 / X0;
271  Real K1 = (1.0 - Y0) / (1.0 - X0);
272  Real dv_dz = (K1 - K0) / ((K0 - 1.0) * (K1 - 1.0));
273 
274  ds_dz = gas_density * liquid_density * dv_dz;
275  ds_dz /= (gas_density + vapor_mass_fraction * (liquid_density - gas_density)) *
276  (gas_density + vapor_mass_fraction * (liquid_density - gas_density));
277 
278  Real dK0_dp = (-Y0 * dX0_dp - X0 * dY1_dp) / X0 / X0;
279  Real dK0_dT = (-Y0 * dX0_dT - X0 * dY1_dT) / X0 / X0;
280 
281  Real dK1_dp = (Y1 * dX0_dp + (1.0 - X0) * dY1_dp) / (1.0 - X0) / (1.0 - X0);
282  Real dK1_dT = (Y1 * dX0_dT + (1.0 - X0) * dY1_dT) / (1.0 - X0) / (1.0 - X0);
283 
284  Real dv_dp = (*_z[0])[_qp] * dK1_dp / (K1 - 1.0) / (K1 - 1.0) +
285  (1.0 - (*_z[0])[_qp]) * dK0_dp / (K0 - 1.0) / (K0 - 1.0);
286 
287  ds_dp = gas_density * liquid_density * dv_dp +
288  vapor_mass_fraction * (1.0 - vapor_mass_fraction) *
289  (gas_density * dliquid_density_dp - dgas_density_dp * liquid_density);
290  ds_dp /= (gas_density + vapor_mass_fraction * (liquid_density - gas_density)) *
291  (gas_density + vapor_mass_fraction * (liquid_density - gas_density));
292 
293  Real dv_dT = (*_z[0])[_qp] * dK1_dT / (K1 - 1.0) / (K1 - 1.0) +
294  (1.0 - (*_z[0])[_qp]) * dK0_dT / (K0 - 1.0) / (K0 - 1.0);
295 
296  ds_dT = gas_density * liquid_density * dv_dT +
297  vapor_mass_fraction * (1.0 - vapor_mass_fraction) *
298  (gas_density * dliquid_density_dT - dgas_density_dT * liquid_density);
299  ds_dT /= (gas_density + vapor_mass_fraction * (liquid_density - gas_density)) *
300  (gas_density + vapor_mass_fraction * (liquid_density - gas_density));
301  }
302 
303  // Save derivatives in FluidStateProperties vector
304  liquid.dsaturation_dp = -ds_dp;
305  liquid.dsaturation_dT = -ds_dT;
306  liquid.dsaturation_dz = -ds_dz;
307  gas.dsaturation_dp = ds_dp;
308  gas.dsaturation_dT = ds_dT;
309  gas.dsaturation_dz = ds_dz;
310 
311  liquid.dmass_fraction_dp[_aqueous_fluid_component] = -dX0_dp;
312  liquid.dmass_fraction_dp[_gas_fluid_component] = dX0_dp;
313  liquid.dmass_fraction_dT[_aqueous_fluid_component] = -dX0_dT;
314  liquid.dmass_fraction_dT[_gas_fluid_component] = dX0_dT;
315  liquid.dmass_fraction_dz[_aqueous_fluid_component] = -dX0_dz;
316  liquid.dmass_fraction_dz[_gas_fluid_component] = dX0_dz;
317  gas.dmass_fraction_dp[_aqueous_fluid_component] = -dY0_dp;
318  gas.dmass_fraction_dp[_gas_fluid_component] = dY0_dp;
319  gas.dmass_fraction_dT[_aqueous_fluid_component] = -dY0_dT;
320  gas.dmass_fraction_dT[_gas_fluid_component] = dY0_dT;
321  gas.dmass_fraction_dz[_aqueous_fluid_component] = -dY0_dz;
322  gas.dmass_fraction_dz[_gas_fluid_component] = dY0_dz;
323 
324  liquid.dfluid_density_dp = dliquid_density_dp;
325  liquid.dfluid_density_dT = dliquid_density_dT;
326  liquid.dfluid_density_dz = dliquid_density_dz;
327  gas.dfluid_density_dp = dgas_density_dp;
328  gas.dfluid_density_dT = dgas_density_dT;
329  gas.dfluid_density_dz = dgas_density_dz;
330 
331  liquid.dfluid_viscosity_dp = dliquid_viscosity_dp;
332  liquid.dfluid_viscosity_dT = dliquid_viscosity_dT;
333  liquid.dfluid_viscosity_dz = dliquid_viscosity_dz;
334  gas.dfluid_viscosity_dp = dgas_viscosity_dp;
335  gas.dfluid_viscosity_dT = dgas_viscosity_dT;
336  gas.dfluid_viscosity_dz = dgas_viscosity_dz;
337  }
338 }
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 MaterialProperty< Real > & _temperature
Temperature.
virtual void rho_dpT(Real pressure, Real temperature, Real &rho, Real &drho_dp, Real &drho_dT) const =0
Density and its derivatives wrt pressure and temperature.
const SinglePhaseFluidPropertiesPT & _water_fp
Fluid properties UserObject for H20.
const BrineFluidProperties & _brine_fp
Fluid properties UserObject for brine.
const VariableValue & _gas_porepressure
Porepressure.
const MaterialProperty< std::vector< Real > > & _saturation_old
Old value of saturation.
void massFractions(Real pressure, Real temperature, Real xnacl, Real &xco2l, Real &dxco2l_dp, Real &dxco2l_dT, Real &xh2og, Real &dxh2og_dp, Real &dxh2og_dT) const
Mass fractions of CO2 and brine calculated using mutual solubilities given by Spycher, Pruess and Ennis-King, CO2-H2O mixtures in the geological sequestration of CO2.
virtual Real rho(Real pressure, Real temperature, Real xnacl) const override
Density.
const VariableValue & _xnacl
Salt mass fraction (kg/kg)
virtual Real capillaryPressure(Real saturation) const
Capillary pressure is calculated as a function of true saturation.
const unsigned int _gas_fluid_component
Fluid component number of the gas phase.
const PorousFlowCapillaryPressure & _pc_uo
Capillary pressure UserObject.
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)
std::vector< FluidStateProperties > _fsp
FluidStateProperties data structure.
const unsigned int _aqueous_phase_number
Phase number of the aqueous phase.
virtual void mu_drhoT_from_rho_T(Real density, Real temperature, Real ddensity_dT, Real &mu, Real &dmu_drho, Real &dmu_dT) const =0
Dynamic viscosity and its derivatives wrt density and temperature.
const SinglePhaseFluidPropertiesPT & _co2_fp
Fluid properties UserObject for CO2.
const unsigned int _gas_phase_number
Phase number of the gas phase.
bool _is_initqp
Flag to indicate whether to calculate stateful properties.
const unsigned int _aqueous_fluid_component
Fluid component number of the aqueous component.
const Real _T_c2k
Conversion from degrees Celsius to degrees Kelvin.
const unsigned int _num_components
Number of components.
const std::string pressure
Definition: NS.h:24
virtual void mu_drhoTx(Real water_density, Real temperature, Real xnacl, Real dwater_density_dT, Real &mu, Real &dmu_drho, Real &dmu_dT, Real &dmu_dx) const override
Dynamic viscosity and its derivatives wrt density, temperature and mass fraction. ...
std::vector< const VariableValue * > _z
Total mass fraction(s) of the gas component(s) summed over all phases.
Real vaporMassFraction(std::vector< Real > &Ki) const
Solves Rachford-Rice equation to provide vapor mass fraction.
Real PorousFlowFluidStateFlashBase::vaporMassFraction ( std::vector< Real > &  Ki) const
protectedinherited

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
Kiequilibrium constants
Returns
vapor mass fraction

Definition at line 366 of file PorousFlowFluidStateFlashBase.C.

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

367 {
368  // Check that the size of the equilibrium constant vector is correct
369  if (Ki.size() != _num_components)
370  mooseError("The number of equilibrium components passed to vaporMassFraction is not correct");
371 
372  Real v;
373 
374  // If there are only two components, an analytical solution is possible
375  if (_num_components == 2)
376  v = ((*_z[0])[_qp] * (Ki[1] - Ki[0]) - (Ki[1] - 1.0)) / ((Ki[0] - 1.0) * (Ki[1] - 1.0));
377  else
378  {
379  // More than two components - solve the Rachford-Rice equation using
380  // Newton-Raphson method.
381  // Initial guess for vapor mass fraction
382  Real v0 = 0.5;
383  unsigned int iter = 0;
384 
385  while (std::abs(rachfordRice(v0, Ki)) > _nr_tol)
386  {
387  v0 = v0 - rachfordRice(v0, Ki) / rachfordRiceDeriv(v0, Ki);
388  iter++;
389 
390  if (iter > _nr_max_its)
391  break;
392  }
393  v = v0;
394  }
395  return v;
396 }
Real rachfordRiceDeriv(Real x, std::vector< Real > &Ki) const
Derivative of Rachford-Rice equation wrt vapor fraction.
const Real _nr_tol
Tolerance for Newton-Raphson iterations.
Real rachfordRice(Real x, std::vector< Real > &Ki) const
Rachford-Rice equation for vapor fraction.
const unsigned int _num_components
Number of components.
const Real _nr_max_its
Maximum number of iterations for the Newton-Raphson iterations.
std::vector< const VariableValue * > _z
Total mass fraction(s) of the gas component(s) summed over all phases.

Member Data Documentation

const unsigned int PorousFlowFluidStateFlashBase::_aqueous_fluid_component
protectedinherited
const unsigned int PorousFlowFluidStateFlashBase::_aqueous_phase_number
protectedinherited
const BrineFluidProperties& PorousFlowFluidStateBrineCO2::_brine_fp
protected

Fluid properties UserObject for brine.

Definition at line 155 of file PorousFlowFluidStateBrineCO2.h.

Referenced by PorousFlowFluidStateBrineCO2(), and thermophysicalProperties().

const SinglePhaseFluidPropertiesPT& PorousFlowFluidStateBrineCO2::_co2_fp
protected

Fluid properties UserObject for CO2.

Definition at line 157 of file PorousFlowFluidStateBrineCO2.h.

Referenced by fugacityCoefficientCO2(), fugacityCoefficientH2O(), PorousFlowFluidStateBrineCO2(), and thermophysicalProperties().

MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowFluidStateFlashBase::_dfluid_density_dvar
protectedinherited

Derivative of the fluid density for each phase wrt PorousFlow variables.

Definition at line 177 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::computeQpProperties(), and PorousFlowFluidStateFlashBase::setMaterialVectorSize().

MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowFluidStateFlashBase::_dfluid_viscosity_dvar
protectedinherited

Derivative of the fluid viscosity for each phase wrt PorousFlow variables.

Definition at line 181 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::computeQpProperties(), and PorousFlowFluidStateFlashBase::setMaterialVectorSize().

MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowVariableBase::_dgradp_qp_dgradv
protectedinherited

d(grad porepressure)/d(grad PorousFlow variable) at the quadpoints

Definition at line 52 of file PorousFlowVariableBase.h.

MaterialProperty<std::vector<std::vector<RealGradient> > >* const PorousFlowVariableBase::_dgradp_qp_dv
protectedinherited

d(grad porepressure)/d(PorousFlow variable) at the quadpoints

Definition at line 55 of file PorousFlowVariableBase.h.

MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowVariableBase::_dgrads_qp_dgradv
protectedinherited

d(grad saturation)/d(grad PorousFlow variable) at the quadpoints

Definition at line 67 of file PorousFlowVariableBase.h.

MaterialProperty<std::vector<std::vector<RealGradient> > >* const PorousFlowVariableBase::_dgrads_qp_dv
protectedinherited

d(grad saturation)/d(PorousFlow variable) at the quadpoints

Definition at line 70 of file PorousFlowVariableBase.h.

MaterialProperty<std::vector<std::vector<std::vector<Real> > > >& PorousFlowFluidStateFlashBase::_dmass_frac_dvar
protectedinherited

Derivative of the mass fraction matrix with respect to the Porous Flow variables.

Definition at line 170 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::computeQpProperties(), and PorousFlowFluidStateFlashBase::setMaterialVectorSize().

MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowVariableBase::_dporepressure_dvar
protectedinherited
MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowVariableBase::_dsaturation_dvar
protectedinherited
const MaterialProperty<std::vector<Real> >& PorousFlowFluidStateFlashBase::_dtemperature_dvar
protectedinherited

Derivative of temperature wrt PorousFlow variables.

Definition at line 164 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::computeQpProperties().

MaterialProperty<std::vector<Real> >& PorousFlowFluidStateFlashBase::_fluid_density
protectedinherited
MaterialProperty<std::vector<Real> >& PorousFlowFluidStateFlashBase::_fluid_viscosity
protectedinherited
std::vector<FluidStateProperties> PorousFlowFluidStateFlashBase::_fsp
protectedinherited
const unsigned int PorousFlowFluidStateFlashBase::_gas_fluid_component
protectedinherited
const VariableGradient& PorousFlowFluidStateFlashBase::_gas_gradp_qp
protectedinherited

Gradient of porepressure (only defined at the qps)

Definition at line 136 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::computeQpProperties().

const unsigned int PorousFlowFluidStateFlashBase::_gas_phase_number
protectedinherited
const VariableValue& PorousFlowFluidStateFlashBase::_gas_porepressure
protectedinherited
const unsigned int PorousFlowFluidStateFlashBase::_gas_porepressure_varnum
protectedinherited

Moose variable number of the gas porepressure.

Definition at line 138 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::computeQpProperties().

MaterialProperty<std::vector<std::vector<RealGradient> > >* PorousFlowFluidStateFlashBase::_grad_mass_frac_qp
protectedinherited

Gradient of the mass fraction matrix (only defined at the qps)

Definition at line 168 of file PorousFlowFluidStateFlashBase.h.

MaterialProperty<std::vector<RealGradient> >* const PorousFlowVariableBase::_gradp_qp
protectedinherited

Grad(p) at the quadpoints.

Definition at line 49 of file PorousFlowVariableBase.h.

MaterialProperty<std::vector<RealGradient> >* const PorousFlowVariableBase::_grads_qp
protectedinherited

Grad(s) at the quadpoints.

Definition at line 64 of file PorousFlowVariableBase.h.

const MaterialProperty<RealGradient>& PorousFlowFluidStateFlashBase::_gradT_qp
protectedinherited

Gradient of temperature (only defined at the qps)

Definition at line 162 of file PorousFlowFluidStateFlashBase.h.

std::vector<const VariableGradient *> PorousFlowFluidStateFlashBase::_gradz_qp
protectedinherited

Gradient(s) of total mass fraction(s) of the gas component(s) (only defined at the qps)

Definition at line 144 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::computeQpProperties(), and PorousFlowFluidStateFlashBase::PorousFlowFluidStateFlashBase().

const Real PorousFlowFluidStateBrineCO2::_invMh2o
protected

Inverse of molar mass of H2O (mol/kg)

Definition at line 163 of file PorousFlowFluidStateBrineCO2.h.

Referenced by massFractions().

bool PorousFlowFluidStateFlashBase::_is_initqp
protectedinherited
MaterialProperty<std::vector<std::vector<Real> > >& PorousFlowFluidStateFlashBase::_mass_frac
protectedinherited
const Real PorousFlowFluidStateBrineCO2::_Mco2
protected

Molar mass of CO2 (kg/mol)

Definition at line 165 of file PorousFlowFluidStateBrineCO2.h.

Referenced by fugacityCoefficientCO2(), fugacityCoefficientH2O(), massFractions(), and partialDensityCO2().

const Real PorousFlowFluidStateBrineCO2::_Mh2o
protected

Molar mass of water (kg/mol)

Definition at line 161 of file PorousFlowFluidStateBrineCO2.h.

Referenced by massFractions().

const Real PorousFlowFluidStateBrineCO2::_Mnacl
protected

Molar mass of NaCL.

Definition at line 167 of file PorousFlowFluidStateBrineCO2.h.

Referenced by activityCoefficient(), and massFractions().

const Real PorousFlowFluidStateFlashBase::_nr_max_its
protectedinherited

Maximum number of iterations for the Newton-Raphson iterations.

Definition at line 188 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::vaporMassFraction().

const Real PorousFlowFluidStateFlashBase::_nr_tol
protectedinherited

Tolerance for Newton-Raphson iterations.

Definition at line 190 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::vaporMassFraction().

const unsigned int PorousFlowVariableBase::_num_components
protectedinherited
const unsigned int PorousFlowVariableBase::_num_pf_vars
protectedinherited
const unsigned int PorousFlowVariableBase::_num_phases
protectedinherited
const unsigned int PorousFlowFluidStateFlashBase::_num_z_vars
protectedinherited
const PorousFlowCapillaryPressure& PorousFlowFluidStateFlashBase::_pc_uo
protectedinherited
MaterialProperty<std::vector<Real> >& PorousFlowVariableBase::_porepressure
protectedinherited
const unsigned int PorousFlowFluidStateFlashBase::_pvar
protectedinherited

PorousFlow variable number of the gas porepressure.

Definition at line 140 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateFlashBase::computeQpProperties().

const Real PorousFlowFluidStateFlashBase::_R
protectedinherited

Universal gas constant (J/mol/K)

Definition at line 186 of file PorousFlowFluidStateFlashBase.h.

Referenced by PorousFlowFluidStateWaterNCG::enthalpyOfDissolution().

const Real PorousFlowFluidStateBrineCO2::_Rbar
protected

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

Definition at line 169 of file PorousFlowFluidStateBrineCO2.h.

Referenced by fugacityCoefficientCO2(), fugacityCoefficientH2O(), and massFractions().

MaterialProperty<std::vector<Real> >& PorousFlowVariableBase::_saturation
protectedinherited
const MaterialProperty<std::vector<Real> >& PorousFlowFluidStateFlashBase::_saturation_old
protectedinherited

Old value of saturation.

Definition at line 172 of file PorousFlowFluidStateFlashBase.h.

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

const Real PorousFlowFluidStateFlashBase::_T_c2k
protectedinherited
const MaterialProperty<Real>& PorousFlowFluidStateFlashBase::_temperature
protectedinherited
const SinglePhaseFluidPropertiesPT& PorousFlowFluidStateBrineCO2::_water_fp
protected

Fluid properties UserObject for H20.

Definition at line 159 of file PorousFlowFluidStateBrineCO2.h.

Referenced by thermophysicalProperties().

const VariableValue& PorousFlowFluidStateBrineCO2::_xnacl
protected

Salt mass fraction (kg/kg)

Definition at line 153 of file PorousFlowFluidStateBrineCO2.h.

Referenced by thermophysicalProperties().

std::vector<const VariableValue *> PorousFlowFluidStateFlashBase::_z
protectedinherited
std::vector<unsigned int> PorousFlowFluidStateFlashBase::_z_varnum
protectedinherited
std::vector<unsigned int> PorousFlowFluidStateFlashBase::_zvar
protectedinherited

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