Add vapporized water to gas in co2brine model

This commit is contained in:
Tor Harald Sandve 2023-02-17 15:39:19 +01:00
parent 7526a7bb77
commit c84c1d54dd
3 changed files with 135 additions and 33 deletions

View File

@ -93,9 +93,6 @@ int main(int argc, char **argv)
if (argc > 6)
rs = atof(argv[6]);
size_t num_regions = 1;
Opm::Co2GasPvt<double> co2Pvt(num_regions);
const double MmNaCl = 58e-3; // molar mass of NaCl [kg/mol]
// convert to mass fraction
std::vector<double> salinity = {0.0};
@ -103,6 +100,8 @@ int main(int argc, char **argv)
salinity[0] = 1 / ( 1 + 1 / (molality*MmNaCl));
Opm::BrineCo2Pvt<double> brineCo2Pvt(salinity);
Opm::Co2GasPvt<double> co2Pvt(salinity);
double value;
if (prop == "density") {
if (phase == "CO2") {
@ -146,9 +145,19 @@ int main(int argc, char **argv)
throw std::runtime_error("phase " + phase + " not recognized. Use either CO2 or brine");
}
} else if (prop == "rsSat") {
if (phase == "CO2") {
value = co2Pvt.saturatedWaterVaporizationFactor(/*regionIdx=*/0,
T,
p);
} else if (phase == "brine") {
value = brineCo2Pvt.saturatedGasDissolutionFactor(/*regionIdx=*/0,
T,
p);
} else {
throw std::runtime_error("phase " + phase + " not recognized. Use either CO2 or brine");
}
} else if (prop == "diffusionCoefficient") {
size_t comp_idx = 0; // not used
if (phase == "CO2") {

View File

@ -30,6 +30,7 @@
#include <opm/material/Constants.hpp>
#include <opm/material/components/CO2.hpp>
#include <opm/material/components/Brine.hpp>
#include <opm/material/components/SimpleHuDuanH2O.hpp>
#include <opm/material/common/UniformTabulated2DFunction.hpp>
#include <opm/material/binarycoefficients/Brine_CO2.hpp>
@ -52,6 +53,7 @@ class Co2GasPvt
{
using CO2 = ::Opm::CO2<Scalar>;
using H2O = SimpleHuDuanH2O<Scalar>;
using Brine = ::Opm::Brine<Scalar, H2O>;
static constexpr bool extrapolate = true;
public:
@ -60,18 +62,22 @@ public:
explicit Co2GasPvt() = default;
Co2GasPvt(size_t numRegions,
Co2GasPvt(const std::vector<Scalar>& salinity,
Scalar T_ref = 288.71, //(273.15 + 15.56)
Scalar P_ref = 101325)
: salinity_(salinity)
{
setNumRegions(numRegions);
for (size_t i = 0; i < numRegions; ++i) {
int num_regions = salinity_.size();
setNumRegions(num_regions);
Brine::salinity = salinity[0];
for (int i = 0; i < num_regions; ++i) {
gasReferenceDensity_[i] = CO2::gasDensity(T_ref, P_ref, extrapolate);
brineReferenceDensity_[i] = Brine::liquidDensity(T_ref, P_ref, extrapolate);
}
}
#if HAVE_ECL_INPUT
/*!
* \brief Initialize the parameters for co2 gas using an ECL deck.
* \brief Initialize the parameters for CO2 gas using an ECL deck.
*/
void initFromState(const EclipseState& eclState, const Schedule&);
#endif
@ -79,6 +85,8 @@ public:
void setNumRegions(size_t numRegions)
{
gasReferenceDensity_.resize(numRegions);
brineReferenceDensity_.resize(numRegions);
salinity_.resize(numRegions);
}
@ -86,15 +94,16 @@ public:
* \brief Initialize the reference densities of all fluids for a given PVT region
*/
void setReferenceDensities(unsigned regionIdx,
Scalar /*rhoRefOil*/,
Scalar rhoRefBrine,
Scalar rhoRefGas,
Scalar /*rhoRefWater*/)
{
gasReferenceDensity_[regionIdx] = rhoRefGas;
brineReferenceDensity_[regionIdx] = rhoRefBrine;;
}
/*!
* \brief Finish initializing the oil phase PVT properties.
* \brief Finish initializing the co2 phase PVT properties.
*/
void initEnd()
{
@ -111,12 +120,16 @@ public:
* \brief Returns the specific enthalpy [J/kg] of gas given a set of parameters.
*/
template <class Evaluation>
Evaluation internalEnergy(unsigned,
Evaluation internalEnergy(unsigned regionIdx,
const Evaluation& temperature,
const Evaluation& pressure,
const Evaluation&) const
const Evaluation& rv) const
{
return CO2::gasInternalEnergy(temperature, pressure, extrapolate);
Evaluation result = 0;
const Evaluation xBrine = convertRvwToXgW_(rv,regionIdx);
result += xBrine * Brine::gasEnthalpy(temperature, pressure);
result += (1 - xBrine) * CO2::gasEnthalpy(temperature, pressure, extrapolate);
return result;
}
/*!
@ -131,13 +144,14 @@ public:
{ return saturatedViscosity(regionIdx, temperature, pressure); }
/*!
* \brief Returns the dynamic viscosity [Pa s] of oil saturated gas at given pressure.
* \brief Returns the dynamic viscosity [Pa s] of fluid phase at saturated conditions.
*/
template <class Evaluation>
Evaluation saturatedViscosity(unsigned /*regionIdx*/,
const Evaluation& temperature,
const Evaluation& pressure) const
{
// Neglects impact of vapporized water on the visosity
return CO2::gasViscosity(temperature, pressure, extrapolate);
}
@ -153,13 +167,14 @@ public:
{ return saturatedInverseFormationVolumeFactor(regionIdx, temperature, pressure); }
/*!
* \brief Returns the formation volume factor [-] of oil saturated gas at given pressure.
* \brief Returns the formation volume factor [-] of water saturated gas at given pressure.
*/
template <class Evaluation>
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
const Evaluation& temperature,
const Evaluation& pressure) const
{
// Neglects impact of vapporized water on the density
return CO2::gasDensity(temperature, pressure, extrapolate)/gasReferenceDensity_[regionIdx];
}
@ -173,46 +188,46 @@ public:
Evaluation saturationPressure(unsigned /*regionIdx*/,
const Evaluation& /*temperature*/,
const Evaluation& /*Rv*/) const
{ return 0.0; /* this is dry gas! */ }
{ return 0.0; /* not implemented */ }
/*!
* \brief Returns the water vaporization factor \f$R_vw\f$ [m^3/m^3] of the water phase.
*/
template <class Evaluation>
Evaluation saturatedWaterVaporizationFactor(unsigned /*regionIdx*/,
const Evaluation& /*temperature*/,
const Evaluation& /*pressure*/) const
{ return 0.0; /* this is non-humid gas! */ }
Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx,
const Evaluation& temperature,
const Evaluation& pressure) const
{ return rvwSat_(regionIdx, temperature, pressure); }
/*!
* \brief Returns the water vaporization factor \f$R_vw\f$ [m^3/m^3] of water saturated gas.
* \brief Returns the water vaporization factor \f$R_vw\f$ [m^3/m^3] of water phase.
*/
template <class Evaluation = Scalar>
Evaluation saturatedWaterVaporizationFactor(unsigned /*regionIdx*/,
const Evaluation& /*temperature*/,
const Evaluation& /*pressure*/,
Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx,
const Evaluation& temperature,
const Evaluation& pressure,
const Evaluation& /*saltConcentration*/) const
{ return 0.0; }
{ return rvwSat_(regionIdx, temperature, pressure); }
/*!
* \brief Returns the oil vaporization factor \f$R_v\f$ [m^3/m^3] of the oil phase.
*/
template <class Evaluation>
Evaluation saturatedOilVaporizationFactor(unsigned /*regionIdx*/,
const Evaluation& /*temperature*/,
const Evaluation& /*pressure*/,
Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
const Evaluation& temperature,
const Evaluation& pressure,
const Evaluation& /*oilSaturation*/,
const Evaluation& /*maxOilSaturation*/) const
{ return 0.0; /* this is dry gas! */ }
{ return rvwSat_(regionIdx, temperature, pressure); }
/*!
* \brief Returns the oil vaporization factor \f$R_v\f$ [m^3/m^3] of the oil phase.
*/
template <class Evaluation>
Evaluation saturatedOilVaporizationFactor(unsigned /*regionIdx*/,
const Evaluation& /*temperature*/,
const Evaluation& /*pressure*/) const
{ return 0.0; /* this is dry gas! */ }
Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
const Evaluation& temperature,
const Evaluation& pressure) const
{ return rvwSat_(regionIdx, temperature, pressure); }
template <class Evaluation>
Evaluation diffusionCoefficient(const Evaluation& temperature,
@ -222,11 +237,88 @@ public:
return BinaryCoeffBrineCO2::gasDiffCoeff(temperature, pressure, extrapolate);
}
Scalar gasReferenceDensity(unsigned regionIdx) const
const Scalar gasReferenceDensity(unsigned regionIdx) const
{ return gasReferenceDensity_[regionIdx]; }
const Scalar oilReferenceDensity(unsigned regionIdx) const
{ return brineReferenceDensity_[regionIdx]; }
const Scalar waterReferenceDensity(unsigned regionIdx) const
{ return brineReferenceDensity_[regionIdx]; }
const Scalar salinity(unsigned regionIdx) const
{ return salinity_[regionIdx]; }
private:
template <class LhsEval>
LhsEval rvwSat_(unsigned regionIdx,
const LhsEval& temperature,
const LhsEval& pressure) const
{
if (!enableDissolution_)
return 0.0;
// calulate the equilibrium composition for the given
// temperature and pressure.
LhsEval xgH2O;
LhsEval xlCO2;
BinaryCoeffBrineCO2::calculateMoleFractions(temperature,
pressure,
salinity_[regionIdx],
/*knownPhaseIdx=*/-1,
xlCO2,
xgH2O,
extrapolate);
// normalize the phase compositions
xgH2O = max(0.0, min(1.0, xgH2O));
return convertXgOToRvw(convertxgOToXgO(xgH2O), regionIdx);
}
/*!
* \brief Convert the mass fraction of the water component in the gas phase to the
* corresponding water vaporization factor.
*/
template <class LhsEval>
LhsEval convertXgOToRvw(const LhsEval& XgO, unsigned regionIdx) const
{
Scalar rho_oRef = brineReferenceDensity_[regionIdx];
Scalar rho_gRef = gasReferenceDensity_[regionIdx];
return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
}
/*!
* \brief Convert a water vapporization factor to the the corresponding mass fraction
* of the water component in the gas phase.
*/
template <class LhsEval>
LhsEval convertRvwToXgW_(const LhsEval& Rvw, unsigned regionIdx) const
{
Scalar rho_wRef = brineReferenceDensity_[regionIdx];
Scalar rho_gRef = gasReferenceDensity_[regionIdx];
const LhsEval& rho_wG = Rvw*rho_gRef;
return rho_wG/(rho_wRef + rho_wG);
}
/*!
* \brief Convert a water mole fraction in the gas phase the corresponding mass fraction.
*/
template <class LhsEval>
LhsEval convertxgOToXgO(const LhsEval& xgW) const
{
Scalar M_CO2 = CO2::molarMass();
Scalar M_Brine = Brine::molarMass();
return xgW*M_Brine / (xgW*(M_Brine - M_CO2) + M_CO2);
}
std::vector<Scalar> brineReferenceDensity_;
std::vector<Scalar> gasReferenceDensity_;
std::vector<Scalar> salinity_;
bool enableDissolution_ = true;
};
} // namespace Opm

View File

@ -55,6 +55,7 @@ initFromState(const EclipseState& eclState, const Schedule&)
Scalar T_ref = eclState.getTableManager().stCond().temperature;
Scalar P_ref = eclState.getTableManager().stCond().pressure;
gasReferenceDensity_[regionIdx] = CO2::gasDensity(T_ref, P_ref, extrapolate);
brineReferenceDensity_[regionIdx] = Brine::liquidDensity(T_ref, P_ref, extrapolate);
initEnd();
}