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>
49 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
56 void initFromDeck(
const Deck& deck,
const EclipseState& eclState)
58 const auto& pvdoTables = eclState.getTableManager().getPvdoTables();
59 const auto& densityKeyword = deck.getKeyword(
"DENSITY");
61 assert(pvdoTables.size() == densityKeyword.size());
64 setNumRegions(numRegions);
66 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
67 Scalar rhoRefO = densityKeyword.getRecord(regionIdx).getItem(
"OIL").getSIDouble(0);
68 Scalar rhoRefG = densityKeyword.getRecord(regionIdx).getItem(
"GAS").getSIDouble(0);
69 Scalar rhoRefW = densityKeyword.getRecord(regionIdx).getItem(
"WATER").getSIDouble(0);
73 const auto& pvdoTable = pvdoTables.getTable<PvdoTable>(regionIdx);
75 const auto& BColumn(pvdoTable.getFormationFactorColumn());
76 std::vector<Scalar> invBColumn(BColumn.size());
77 for (
unsigned i = 0; i < invBColumn.size(); ++i)
78 invBColumn[i] = 1/BColumn[i];
80 inverseOilB_[regionIdx].setXYArrays(pvdoTable.numRows(),
81 pvdoTable.getPressureColumn(),
83 oilMu_[regionIdx].setXYArrays(pvdoTable.numRows(),
84 pvdoTable.getPressureColumn(),
85 pvdoTable.getViscosityColumn());
90 #endif // HAVE_OPM_PARSER
94 oilReferenceDensity_.resize(numRegions);
95 inverseOilB_.resize(numRegions);
96 inverseOilBMu_.resize(numRegions);
97 oilMu_.resize(numRegions);
108 oilReferenceDensity_[regionIdx] = rhoRefOil;
122 { inverseOilB_[regionIdx] = invBo; }
130 { oilMu_[regionIdx] = muo; }
138 size_t numRegions = oilMu_.size();
139 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
142 const auto& oilMu = oilMu_[regionIdx];
143 const auto& invOilB = inverseOilB_[regionIdx];
144 assert(oilMu.numSamples() == invOilB.numSamples());
146 std::vector<Scalar> invBMuColumn;
147 std::vector<Scalar> pressureColumn;
148 invBMuColumn.resize(oilMu.numSamples());
149 pressureColumn.resize(oilMu.numSamples());
151 for (
unsigned pIdx = 0; pIdx < oilMu.numSamples(); ++pIdx) {
152 pressureColumn[pIdx] = invOilB.xAt(pIdx);
153 invBMuColumn[pIdx] = invOilB.valueAt(pIdx)*1/oilMu.valueAt(pIdx);
156 inverseOilBMu_[regionIdx].setXYArrays(pressureColumn.size(),
166 {
return inverseOilBMu_.size(); }
171 template <
class Evaluation>
173 const Evaluation& temperature,
174 const Evaluation& pressure,
175 const Evaluation& )
const
181 template <
class Evaluation>
184 const Evaluation& pressure)
const
186 const Evaluation& invBo = inverseOilB_[regionIdx].eval(pressure,
true);
187 const Evaluation& invMuoBo = inverseOilBMu_[regionIdx].eval(pressure,
true);
189 return invBo/invMuoBo;
195 template <
class Evaluation>
198 const Evaluation& pressure,
199 const Evaluation& )
const
200 {
return inverseOilB_[regionIdx].eval(pressure,
true); }
207 template <
class Evaluation>
210 const Evaluation& pressure)
const
211 {
return inverseOilB_[regionIdx].eval(pressure,
true); }
216 template <
class Evaluation>
219 const Evaluation& )
const
225 template <
class Evaluation>
239 template <
class Evaluation>
242 const Evaluation& )
const
245 template <
class Evaluation>
246 Evaluation saturatedGasMassFraction(
unsigned ,
248 const Evaluation& )
const
251 template <
class Evaluation>
252 Evaluation saturatedGasMoleFraction(
unsigned ,
254 const Evaluation& )
const
258 std::vector<Scalar> oilReferenceDensity_;
259 std::vector<TabulatedOneDFunction> inverseOilB_;
260 std::vector<TabulatedOneDFunction> oilMu_;
261 std::vector<TabulatedOneDFunction> inverseOilBMu_;
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
Definition: DeadOilPvt.hpp:46
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: DeadOilPvt.hpp:165
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &) const
Returns the formation volume factor [-] of the fluid phase.
Definition: DeadOilPvt.hpp:196
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:172
void setOilViscosity(unsigned regionIdx, const TabulatedOneDFunction &muo)
Initialize the viscosity of the oil phase.
Definition: DeadOilPvt.hpp:129
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:240
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: DeadOilPvt.hpp:135
Class implementing cubic splines.
Evaluation saturatedGasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: DeadOilPvt.hpp:217
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of saturated oil.
Definition: DeadOilPvt.hpp:208
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: DeadOilPvt.hpp:103
Evaluation saturatedGasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &, const Evaluation &, Scalar) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: DeadOilPvt.hpp:226
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:182
Implements a linearly interpolated scalar function that depends on one variable.
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47
void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedOneDFunction &invBo)
Initialize the function for the oil formation volume factor.
Definition: DeadOilPvt.hpp:121