25 #ifndef OPM_DEAD_OIL_PVT_HPP
26 #define OPM_DEAD_OIL_PVT_HPP
33 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
37 template <
class Scalar>
38 class GasPvtMultiplexer;
44 template <
class Scalar>
50 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
57 void initFromDeck(DeckConstPtr deck, EclipseStateConstPtr eclState)
59 const auto& pvdoTables = eclState->getTableManager()->getPvdoTables();
60 DeckKeywordConstPtr densityKeyword = deck->getKeyword(
"DENSITY");
62 assert(pvdoTables.size() == densityKeyword->size());
64 size_t numRegions = pvdoTables.size();
67 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
68 Scalar rhoRefO = densityKeyword->getRecord(regionIdx)->getItem(
"OIL")->getSIDouble(0);
69 Scalar rhoRefG = densityKeyword->getRecord(regionIdx)->getItem(
"GAS")->getSIDouble(0);
70 Scalar rhoRefW = densityKeyword->getRecord(regionIdx)->getItem(
"WATER")->getSIDouble(0);
74 const auto& pvdoTable = pvdoTables.getTable<PvdoTable>(regionIdx);
76 const auto& BColumn(pvdoTable.getFormationFactorColumn());
77 std::vector<Scalar> invBColumn(BColumn.size());
78 for (
unsigned i = 0; i < invBColumn.size(); ++i)
79 invBColumn[i] = 1/BColumn[i];
81 inverseOilB_[regionIdx].setXYArrays(pvdoTable.numRows(),
82 pvdoTable.getPressureColumn(),
84 oilMu_[regionIdx].setXYArrays(pvdoTable.numRows(),
85 pvdoTable.getPressureColumn(),
86 pvdoTable.getViscosityColumn());
89 #endif // HAVE_OPM_PARSER
93 oilReferenceDensity_.resize(numRegions);
94 inverseOilB_.resize(numRegions);
95 inverseOilBMu_.resize(numRegions);
96 oilMu_.resize(numRegions);
107 oilReferenceDensity_[regionIdx] = rhoRefOil;
121 { inverseOilB_[regionIdx] = invBo; }
129 { oilMu_[regionIdx] = muo; }
137 size_t numRegions = oilMu_.size();
138 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
141 const auto& oilMu = oilMu_[regionIdx];
142 const auto& invOilB = inverseOilB_[regionIdx];
143 assert(oilMu.numSamples() == invOilB.numSamples());
145 std::vector<Scalar> invBMuColumn;
146 std::vector<Scalar> pressureColumn;
147 invBMuColumn.resize(oilMu.numSamples());
148 pressureColumn.resize(oilMu.numSamples());
150 for (
unsigned pIdx = 0; pIdx < oilMu.numSamples(); ++pIdx) {
151 pressureColumn[pIdx] = invOilB.xAt(pIdx);
152 invBMuColumn[pIdx] = invOilB.valueAt(pIdx)*1/oilMu.valueAt(pIdx);
155 inverseOilBMu_[regionIdx].setXYArrays(pressureColumn.size(),
164 template <
class Evaluation>
167 const Evaluation& pressure,
168 const Evaluation& )
const
170 const Evaluation& invBo = inverseOilB_[regionIdx].eval(pressure,
true);
171 const Evaluation& invMuoBo = inverseOilBMu_[regionIdx].eval(pressure,
true);
173 return invBo/invMuoBo;
179 template <
class Evaluation>
181 const Evaluation& temperature,
182 const Evaluation& pressure,
183 const Evaluation& XoG)
const
185 Scalar rhooRef = oilReferenceDensity_[regionIdx];
194 template <
class Evaluation>
197 const Evaluation& pressure,
198 const Evaluation& )
const
199 {
return 1.0 / inverseOilB_[regionIdx].eval(pressure,
true); }
205 template <
class Evaluation>
208 const Evaluation& pressure)
const
213 return 20e3/pressure;
216 template <
class Evaluation>
218 const Evaluation& temperature,
219 const Evaluation& pressure)
const
227 template <
class Evaluation>
229 const Evaluation& temperature,
230 const Evaluation& pressure)
const
239 template <
class Evaluation>
242 const Evaluation& )
const
251 template <
class Evaluation>
254 const Evaluation& )
const
257 template <
class Evaluation>
260 const Evaluation& )
const
263 template <
class Evaluation>
266 const Evaluation& )
const
270 std::vector<Scalar> oilReferenceDensity_;
271 std::vector<TabulatedOneDFunction> inverseOilB_;
272 std::vector<TabulatedOneDFunction> oilMu_;
273 std::vector<TabulatedOneDFunction> inverseOilBMu_;
Evaluation gasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: DeadOilPvt.hpp:240
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: DeadOilPvt.hpp:165
Evaluation saturatedOilGasMoleFraction(unsigned, const Evaluation &, const Evaluation &) const
Definition: DeadOilPvt.hpp:264
Definition: Air_Mesitylene.hpp:31
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition: ConstantCompressibilityOilPvt.hpp:39
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
Definition: DeadOilPvt.hpp:45
void setNumRegions(size_t numRegions)
Definition: DeadOilPvt.hpp:91
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: DeadOilPvt.hpp:102
Class implementing cubic splines.
void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedOneDFunction &invBo)
Initialize the function for the oil formation volume factor.
Definition: DeadOilPvt.hpp:120
Evaluation fugacityCoefficientOil(unsigned, const Evaluation &, const Evaluation &pressure) const
Returns the fugacity coefficient [Pa] of a component in the fluid phase given a set of parameters...
Definition: DeadOilPvt.hpp:206
Implements a linearly interpolated scalar function that depends on one variable.
Evaluation saturatedOilGasMassFraction(unsigned, const Evaluation &, const Evaluation &) const
Definition: DeadOilPvt.hpp:258
void initEnd(const GasPvtMultiplexer *)
Finish initializing the oil phase PVT properties.
Definition: DeadOilPvt.hpp:134
Evaluation formationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &) const
Returns the formation volume factor [-] of the fluid phase.
Definition: DeadOilPvt.hpp:195
void setOilViscosity(unsigned regionIdx, const TabulatedOneDFunction &muo)
Initialize the viscosity of the oil phase.
Definition: DeadOilPvt.hpp:128
Evaluation oilSaturationPressure(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:252
Evaluation fugacityCoefficientWater(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Definition: DeadOilPvt.hpp:217
Evaluation fugacityCoefficientGas(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Definition: DeadOilPvt.hpp:228
Evaluation density(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &XoG) const
Returns the density [kg/m^3] of the fluid phase given a set of parameters.
Definition: DeadOilPvt.hpp:180
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:44