27 #ifndef OPM_DEAD_OIL_PVT_HPP
28 #define OPM_DEAD_OIL_PVT_HPP
35 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
36 #include <opm/parser/eclipse/EclipseState/Tables/PvdoTable.hpp>
37 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
45 template <
class Scalar>
48 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
54 DeadOilPvt(
const std::vector<Scalar>& oilReferenceDensity,
55 const std::vector<TabulatedOneDFunction>& inverseOilB,
56 const std::vector<TabulatedOneDFunction>& oilMu,
57 const std::vector<TabulatedOneDFunction>& inverseOilBMu)
58 : oilReferenceDensity_(oilReferenceDensity)
59 , inverseOilB_(inverseOilB)
61 , inverseOilBMu_(inverseOilBMu)
68 void initFromState(
const EclipseState& eclState,
const Schedule&)
70 const auto& pvdoTables = eclState.getTableManager().getPvdoTables();
71 const auto& densityTable = eclState.getTableManager().getDensityTable();
73 assert(pvdoTables.size() == densityTable.size());
78 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
79 Scalar rhoRefO = densityTable[regionIdx].oil;
80 Scalar rhoRefG = densityTable[regionIdx].gas;
81 Scalar rhoRefW = densityTable[regionIdx].water;
85 const auto& pvdoTable = pvdoTables.getTable<PvdoTable>(regionIdx);
87 const auto& BColumn(pvdoTable.getFormationFactorColumn());
88 std::vector<Scalar> invBColumn(BColumn.size());
89 for (
unsigned i = 0; i < invBColumn.size(); ++i)
90 invBColumn[i] = 1/BColumn[i];
92 inverseOilB_[regionIdx].setXYArrays(pvdoTable.numRows(),
93 pvdoTable.getPressureColumn(),
95 oilMu_[regionIdx].setXYArrays(pvdoTable.numRows(),
96 pvdoTable.getPressureColumn(),
97 pvdoTable.getViscosityColumn());
120 oilReferenceDensity_[regionIdx] = rhoRefOil;
134 { inverseOilB_[regionIdx] = invBo; }
142 { oilMu_[regionIdx] = muo; }
151 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
154 const auto& oilMu = oilMu_[regionIdx];
155 const auto& invOilB = inverseOilB_[regionIdx];
156 assert(oilMu.numSamples() == invOilB.numSamples());
158 std::vector<Scalar> invBMuColumn;
159 std::vector<Scalar> pressureColumn;
160 invBMuColumn.resize(oilMu.numSamples());
161 pressureColumn.resize(oilMu.numSamples());
163 for (
unsigned pIdx = 0; pIdx < oilMu.numSamples(); ++pIdx) {
164 pressureColumn[pIdx] = invOilB.xAt(pIdx);
165 invBMuColumn[pIdx] = invOilB.valueAt(pIdx)*1/oilMu.valueAt(pIdx);
168 inverseOilBMu_[regionIdx].setXYArrays(pressureColumn.size(),
178 {
return inverseOilBMu_.size(); }
183 template <
class Evaluation>
187 const Evaluation&)
const
189 throw std::runtime_error(
"Requested the enthalpy of oil but the thermal option is not enabled");
195 template <
class Evaluation>
197 const Evaluation& temperature,
198 const Evaluation& pressure,
199 const Evaluation& )
const
205 template <
class Evaluation>
208 const Evaluation& pressure)
const
210 const Evaluation& invBo = inverseOilB_[regionIdx].eval(pressure,
true);
211 const Evaluation& invMuoBo = inverseOilBMu_[regionIdx].eval(pressure,
true);
213 return invBo/invMuoBo;
219 template <
class Evaluation>
222 const Evaluation& pressure,
223 const Evaluation& )
const
224 {
return inverseOilB_[regionIdx].eval(pressure,
true); }
231 template <
class Evaluation>
234 const Evaluation& pressure)
const
235 {
return inverseOilB_[regionIdx].eval(pressure,
true); }
240 template <
class Evaluation>
243 const Evaluation& )
const
249 template <
class Evaluation>
254 const Evaluation& )
const
263 template <
class Evaluation>
266 const Evaluation& )
const
269 template <
class Evaluation>
270 Evaluation saturatedGasMassFraction(
unsigned ,
272 const Evaluation& )
const
275 template <
class Evaluation>
276 Evaluation saturatedGasMoleFraction(
unsigned ,
278 const Evaluation& )
const
281 template <
class Evaluation>
282 Evaluation diffusionCoefficient(
const Evaluation& ,
286 throw std::runtime_error(
"Not implemented: The PVT model does not provide a diffusionCoefficient()");
289 const Scalar oilReferenceDensity(
unsigned regionIdx)
const
290 {
return oilReferenceDensity_[regionIdx]; }
292 const std::vector<TabulatedOneDFunction>& inverseOilB()
const
293 {
return inverseOilB_; }
295 const std::vector<TabulatedOneDFunction>& oilMu()
const
298 const std::vector<TabulatedOneDFunction>& inverseOilBMu()
const
299 {
return inverseOilBMu_; }
301 bool operator==(
const DeadOilPvt<Scalar>& data)
const
303 return this->oilReferenceDensity_ == data.oilReferenceDensity_ &&
304 this->inverseOilB() == data.inverseOilB() &&
305 this->oilMu() == data.oilMu() &&
306 this->inverseOilBMu() == data.inverseOilBMu();
310 std::vector<Scalar> oilReferenceDensity_;
311 std::vector<TabulatedOneDFunction> inverseOilB_;
312 std::vector<TabulatedOneDFunction> oilMu_;
313 std::vector<TabulatedOneDFunction> inverseOilBMu_;
Class implementing cubic splines.
Implements a linearly interpolated scalar function that depends on one variable.
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
Definition: DeadOilPvt.hpp:47
void setOilViscosity(unsigned regionIdx, const TabulatedOneDFunction &muo)
Initialize the viscosity of the oil phase.
Definition: DeadOilPvt.hpp:141
Evaluation saturatedGasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: DeadOilPvt.hpp:241
void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedOneDFunction &invBo)
Initialize the function for the oil formation volume factor.
Definition: DeadOilPvt.hpp:133
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: DeadOilPvt.hpp:196
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: DeadOilPvt.hpp:115
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: DeadOilPvt.hpp:147
Evaluation saturatedViscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of gas saturated oil given a pressure.
Definition: DeadOilPvt.hpp:206
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of oil given a set of parameters.
Definition: DeadOilPvt.hpp:184
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of saturated oil.
Definition: DeadOilPvt.hpp:232
Evaluation saturationPressure(unsigned, const Evaluation &, const Evaluation &) const
Returns the saturation pressure of the oil phase [Pa] depending on its mass fraction of the gas compo...
Definition: DeadOilPvt.hpp:264
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &) const
Returns the formation volume factor [-] of the fluid phase.
Definition: DeadOilPvt.hpp:220
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: DeadOilPvt.hpp:177
Evaluation saturatedGasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: DeadOilPvt.hpp:250
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47