27 #ifndef OPM_DRY_GAS_PVT_HPP
28 #define OPM_DRY_GAS_PVT_HPP
35 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
36 #include <opm/input/eclipse/Schedule/Schedule.hpp>
37 #include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
38 #include <opm/input/eclipse/EclipseState/Tables/PvdgTable.hpp>
48 template <
class Scalar>
51 using SamplingPoints = std::vector<std::pair<Scalar, Scalar>>;
57 DryGasPvt(
const std::vector<Scalar>& gasReferenceDensity,
58 const std::vector<TabulatedOneDFunction>& inverseGasB,
59 const std::vector<TabulatedOneDFunction>& gasMu,
60 const std::vector<TabulatedOneDFunction>& inverseGasBMu)
61 : gasReferenceDensity_(gasReferenceDensity)
62 , inverseGasB_(inverseGasB)
64 , inverseGasBMu_(inverseGasBMu)
73 void initFromState(
const EclipseState& eclState,
const Schedule&)
75 const auto& pvdgTables = eclState.getTableManager().getPvdgTables();
76 const auto& densityTable = eclState.getTableManager().getDensityTable();
78 assert(pvdgTables.size() == densityTable.size());
83 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
84 Scalar rhoRefO = densityTable[regionIdx].oil;
85 Scalar rhoRefG = densityTable[regionIdx].gas;
86 Scalar rhoRefW = densityTable[regionIdx].water;
91 constexpr Scalar p = 1.01325e5;
92 constexpr Scalar T = 273.15 + 15.56;
93 constexpr Scalar MO = 175e-3;
95 constexpr Scalar MW = 18.0e-3;
100 const auto& pvdgTable = pvdgTables.getTable<PvdgTable>(regionIdx);
104 std::vector<Scalar> invB(pvdgTable.numRows());
105 const auto& Bg = pvdgTable.getFormationFactorColumn();
106 for (
unsigned i = 0; i < Bg.size(); ++ i) {
110 size_t numSamples = invB.size();
111 inverseGasB_[regionIdx].setXYArrays(numSamples, pvdgTable.getPressureColumn(), invB);
112 gasMu_[regionIdx].setXYArrays(numSamples, pvdgTable.getPressureColumn(), pvdgTable.getViscosityColumn());
136 gasReferenceDensity_[regionIdx] = rhoRefGas;
154 { gasMu_[regionIdx] = mug; }
163 SamplingPoints tmp(samplePoints);
164 auto it = tmp.begin();
165 const auto& endIt = tmp.end();
166 for (; it != endIt; ++ it)
167 std::get<1>(*it) = 1.0/std::get<1>(*it);
169 inverseGasB_[regionIdx].setContainerOfTuples(tmp);
170 assert(inverseGasB_[regionIdx].monotonic());
180 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
183 const auto& gasMu = gasMu_[regionIdx];
184 const auto& invGasB = inverseGasB_[regionIdx];
185 assert(gasMu.numSamples() == invGasB.numSamples());
187 std::vector<Scalar> pressureValues(gasMu.numSamples());
188 std::vector<Scalar> invGasBMuValues(gasMu.numSamples());
189 for (
unsigned pIdx = 0; pIdx < gasMu.numSamples(); ++pIdx) {
190 pressureValues[pIdx] = invGasB.xAt(pIdx);
191 invGasBMuValues[pIdx] = invGasB.valueAt(pIdx) * (1.0/gasMu.valueAt(pIdx));
194 inverseGasBMu_[regionIdx].setXYContainers(pressureValues, invGasBMuValues);
202 {
return gasReferenceDensity_.size(); }
207 template <
class Evaluation>
211 const Evaluation&)
const
213 throw std::runtime_error(
"Requested the enthalpy of gas but the thermal option is not enabled");
219 template <
class Evaluation>
221 const Evaluation& temperature,
222 const Evaluation& pressure,
224 const Evaluation& )
const
230 template <
class Evaluation>
233 const Evaluation& pressure)
const
235 const Evaluation& invBg = inverseGasB_[regionIdx].eval(pressure,
true);
236 const Evaluation& invMugBg = inverseGasBMu_[regionIdx].eval(pressure,
true);
238 return invBg/invMugBg;
244 template <
class Evaluation>
246 const Evaluation& temperature,
247 const Evaluation& pressure,
249 const Evaluation& )
const
255 template <
class Evaluation>
258 const Evaluation& pressure)
const
259 {
return inverseGasB_[regionIdx].eval(pressure,
true); }
267 template <
class Evaluation>
270 const Evaluation& )
const
276 template <
class Evaluation>
279 const Evaluation& )
const
285 template <
class Evaluation>
290 const Evaluation& )
const
296 template <
class Evaluation>
299 const Evaluation& )
const
302 template <
class Evaluation>
303 Evaluation diffusionCoefficient(
const Evaluation& ,
307 throw std::runtime_error(
"Not implemented: The PVT model does not provide a diffusionCoefficient()");
310 const Scalar gasReferenceDensity(
unsigned regionIdx)
const
311 {
return gasReferenceDensity_[regionIdx]; }
313 const std::vector<TabulatedOneDFunction>& inverseGasB()
const
314 {
return inverseGasB_; }
316 const std::vector<TabulatedOneDFunction>& gasMu()
const
319 const std::vector<TabulatedOneDFunction> inverseGasBMu()
const
320 {
return inverseGasBMu_; }
322 bool operator==(
const DryGasPvt<Scalar>& data)
const
324 return gasReferenceDensity_ == data.gasReferenceDensity_ &&
325 inverseGasB_ == data.inverseGasB_ &&
326 gasMu_ == data.gasMu_ &&
327 inverseGasBMu_ == data.inverseGasBMu_;
331 std::vector<Scalar> gasReferenceDensity_;
332 std::vector<TabulatedOneDFunction> inverseGasB_;
333 std::vector<TabulatedOneDFunction> gasMu_;
334 std::vector<TabulatedOneDFunction> inverseGasBMu_;
A central place for various physical constants occuring in some equations.
Implements a linearly interpolated scalar function that depends on one variable.
A central place for various physical constants occuring in some equations.
Definition: Constants.hpp:41
This class represents the Pressure-Volume-Temperature relations of the gas phase without vaporized oi...
Definition: DryGasPvt.hpp:50
Evaluation saturatedViscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of oil saturated gas at given pressure.
Definition: DryGasPvt.hpp:231
Evaluation saturationPressure(unsigned, const Evaluation &, const Evaluation &) const
Returns the saturation pressure of the gas phase [Pa] depending on its mass fraction of the oil compo...
Definition: DryGasPvt.hpp:268
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: DryGasPvt.hpp:201
void setGasFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the formation volume factor of dry gas.
Definition: DryGasPvt.hpp:161
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &, const Evaluation &) const
Returns the formation volume factor [-] of the fluid phase.
Definition: DryGasPvt.hpp:245
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of gas given a set of parameters.
Definition: DryGasPvt.hpp:208
Evaluation saturatedOilVaporizationFactor(unsigned, const Evaluation &, const Evaluation &) const
Returns the oil vaporization factor [m^3/m^3] of the oil phase.
Definition: DryGasPvt.hpp:297
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &, const Evaluation &) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: DryGasPvt.hpp:220
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of oil saturated gas at given pressure.
Definition: DryGasPvt.hpp:256
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: DryGasPvt.hpp:176
Evaluation saturatedWaterVaporizationFactor(unsigned, const Evaluation &, const Evaluation &) const
Returns the water vaporization factor [m^3/m^3] of the water phase.
Definition: DryGasPvt.hpp:277
void setMolarMasses(unsigned, Scalar, Scalar, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: DryGasPvt.hpp:142
void setGasViscosity(unsigned regionIdx, const TabulatedOneDFunction &mug)
Initialize the viscosity of the gas phase.
Definition: DryGasPvt.hpp:153
void setReferenceDensities(unsigned regionIdx, Scalar, Scalar rhoRefGas, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: DryGasPvt.hpp:131
Evaluation saturatedOilVaporizationFactor(unsigned, const Evaluation &, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the oil vaporization factor [m^3/m^3] of the oil phase.
Definition: DryGasPvt.hpp:286
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:48