27 #ifndef OPM_OIL_PVT_THERMAL_HPP
28 #define OPM_OIL_PVT_THERMAL_HPP
38 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
39 #include <opm/parser/eclipse/EclipseState/Tables/SimpleTable.hpp>
40 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
44 template <
class Scalar,
bool enableThermal>
45 class OilPvtMultiplexer;
53 template <
class Scalar>
62 enableThermalDensity_ =
false;
63 enableThermalViscosity_ =
false;
64 enableInternalEnergy_ =
false;
65 isothermalPvt_ =
nullptr;
69 const std::vector<TabulatedOneDFunction>& oilvisctCurves,
70 const std::vector<Scalar>& viscrefPress,
71 const std::vector<Scalar>& viscrefRs,
72 const std::vector<Scalar>& viscRef,
73 const std::vector<Scalar>& oildentRefTemp,
74 const std::vector<Scalar>& oildentCT1,
75 const std::vector<Scalar>& oildentCT2,
76 const std::vector<TabulatedOneDFunction>& internalEnergyCurves,
79 bool enableInternalEnergy)
80 : isothermalPvt_(isothermalPvt)
81 , oilvisctCurves_(oilvisctCurves)
82 , viscrefPress_(viscrefPress)
83 , viscrefRs_(viscrefRs)
85 , oildentRefTemp_(oildentRefTemp)
86 , oildentCT1_(oildentCT1)
87 , oildentCT2_(oildentCT2)
88 , internalEnergyCurves_(internalEnergyCurves)
91 , enableInternalEnergy_(enableInternalEnergy)
98 {
delete isothermalPvt_; }
104 void initFromState(
const EclipseState& eclState,
const Schedule& schedule)
110 isothermalPvt_->initFromState(eclState, schedule);
115 const auto& tables = eclState.getTableManager();
117 enableThermalDensity_ = tables.OilDenT().size() > 0;
118 enableThermalViscosity_ = tables.hasTables(
"OILVISCT");
119 enableInternalEnergy_ = tables.hasTables(
"SPECHEAT");
121 unsigned numRegions = isothermalPvt_->
numRegions();
125 if (enableThermalViscosity_) {
126 if (tables.getViscrefTable().empty())
127 throw std::runtime_error(
"VISCREF is required when OILVISCT is present");
129 const auto& oilvisctTables = tables.getOilvisctTables();
130 const auto& viscrefTable = tables.getViscrefTable();
132 assert(oilvisctTables.size() == numRegions);
133 assert(viscrefTable.size() == numRegions);
135 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
136 const auto& TCol = oilvisctTables[regionIdx].getColumn(
"Temperature").vectorCopy();
137 const auto& muCol = oilvisctTables[regionIdx].getColumn(
"Viscosity").vectorCopy();
138 oilvisctCurves_[regionIdx].setXYContainers(TCol, muCol);
140 viscrefPress_[regionIdx] = viscrefTable[regionIdx].reference_pressure;
141 viscrefRs_[regionIdx] = viscrefTable[regionIdx].reference_rs;
146 Scalar Tref = 273.15 + 20;
149 viscRef_[regionIdx] =
152 viscrefPress_[regionIdx],
153 viscrefRs_[regionIdx]);
158 const auto& oilDenT = tables.OilDenT();
159 if (oilDenT.size() > 0) {
160 assert(oilDenT.size() == numRegions);
161 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
162 const auto& record = oilDenT[regionIdx];
164 oildentRefTemp_[regionIdx] = record.T0;
165 oildentCT1_[regionIdx] = record.C1;
166 oildentCT2_[regionIdx] = record.C2;
170 if (enableInternalEnergy_) {
174 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
175 const auto& specheatTable = tables.getSpecheatTables()[regionIdx];
176 const auto& temperatureColumn = specheatTable.getColumn(
"TEMPERATURE");
177 const auto& cvOilColumn = specheatTable.getColumn(
"CV_OIL");
179 std::vector<double> uSamples(temperatureColumn.size());
181 Scalar u = temperatureColumn[0]*cvOilColumn[0];
182 for (
size_t i = 0;; ++i) {
185 if (i >= temperatureColumn.size() - 1)
190 Scalar c_v0 = cvOilColumn[i];
191 Scalar c_v1 = cvOilColumn[i + 1];
192 Scalar T0 = temperatureColumn[i];
193 Scalar T1 = temperatureColumn[i + 1];
194 u += 0.5*(c_v0 + c_v1)*(T1 - T0);
197 internalEnergyCurves_[regionIdx].setXYContainers(temperatureColumn.vectorCopy(), uSamples);
208 oilvisctCurves_.resize(numRegions);
209 viscrefPress_.resize(numRegions);
210 viscrefRs_.resize(numRegions);
211 viscRef_.resize(numRegions);
212 internalEnergyCurves_.resize(numRegions);
225 {
return enableThermalDensity_; }
231 {
return enableThermalViscosity_; }
233 size_t numRegions()
const
234 {
return viscrefRs_.size(); }
239 template <
class Evaluation>
241 const Evaluation& temperature,
243 const Evaluation&)
const
245 if (!enableInternalEnergy_)
246 throw std::runtime_error(
"Requested the internal energy of oil but it is disabled");
251 return internalEnergyCurves_[regionIdx].eval(temperature,
true);
257 template <
class Evaluation>
259 const Evaluation& temperature,
260 const Evaluation& pressure,
261 const Evaluation& Rs)
const
263 const auto& isothermalMu = isothermalPvt_->
viscosity(regionIdx, temperature, pressure, Rs);
268 const auto& muOilvisct = oilvisctCurves_[regionIdx].eval(temperature,
true);
269 return muOilvisct/viscRef_[regionIdx]*isothermalMu;
275 template <
class Evaluation>
277 const Evaluation& temperature,
278 const Evaluation& pressure)
const
280 const auto& isothermalMu = isothermalPvt_->
saturatedViscosity(regionIdx, temperature, pressure);
285 const auto& muOilvisct = oilvisctCurves_[regionIdx].eval(temperature,
true);
286 return muOilvisct/viscRef_[regionIdx]*isothermalMu;
293 template <
class Evaluation>
295 const Evaluation& temperature,
296 const Evaluation& pressure,
297 const Evaluation& Rs)
const
307 Scalar TRef = oildentRefTemp_[regionIdx];
308 Scalar cT1 = oildentCT1_[regionIdx];
309 Scalar cT2 = oildentCT2_[regionIdx];
310 const Evaluation& Y = temperature - TRef;
312 return b/(1 + (cT1 + cT2*Y)*Y);
318 template <
class Evaluation>
320 const Evaluation& temperature,
321 const Evaluation& pressure)
const
331 Scalar TRef = oildentRefTemp_[regionIdx];
332 Scalar cT1 = oildentCT1_[regionIdx];
333 Scalar cT2 = oildentCT2_[regionIdx];
334 const Evaluation& Y = temperature - TRef;
336 return b/(1 + (cT1 + cT2*Y)*Y);
346 template <
class Evaluation>
348 const Evaluation& temperature,
349 const Evaluation& pressure)
const
359 template <
class Evaluation>
361 const Evaluation& temperature,
362 const Evaluation& pressure,
363 const Evaluation& oilSaturation,
364 const Evaluation& maxOilSaturation)
const
374 template <
class Evaluation>
376 const Evaluation& temperature,
377 const Evaluation& pressure)
const
380 template <
class Evaluation>
381 Evaluation diffusionCoefficient(
const Evaluation& temperature,
382 const Evaluation& pressure,
383 unsigned compIdx)
const
388 const IsothermalPvt* isoThermalPvt()
const
389 {
return isothermalPvt_; }
391 const Scalar oilReferenceDensity(
unsigned regionIdx)
const
394 const std::vector<TabulatedOneDFunction>& oilvisctCurves()
const
395 {
return oilvisctCurves_; }
397 const std::vector<Scalar>& viscrefPress()
const
398 {
return viscrefPress_; }
400 const std::vector<Scalar>& viscrefRs()
const
401 {
return viscrefRs_; }
403 const std::vector<Scalar>& viscRef()
const
406 const std::vector<Scalar>& oildentRefTemp()
const
407 {
return oildentRefTemp_; }
409 const std::vector<Scalar>& oildentCT1()
const
410 {
return oildentCT1_; }
412 const std::vector<Scalar>& oildentCT2()
const
413 {
return oildentCT2_; }
415 const std::vector<TabulatedOneDFunction> internalEnergyCurves()
const
416 {
return internalEnergyCurves_; }
418 bool enableInternalEnergy()
const
419 {
return enableInternalEnergy_; }
421 bool operator==(
const OilPvtThermal<Scalar>& data)
const
423 if (isothermalPvt_ && !data.isothermalPvt_)
425 if (!isothermalPvt_ && data.isothermalPvt_)
428 return (!this->isoThermalPvt() ||
429 (*this->isoThermalPvt() == *data.isoThermalPvt())) &&
430 this->oilvisctCurves() == data.oilvisctCurves() &&
431 this->viscrefPress() == data.viscrefPress() &&
432 this->viscrefRs() == data.viscrefRs() &&
433 this->viscRef() == data.viscRef() &&
434 this->oildentRefTemp() == data.oildentRefTemp() &&
435 this->oildentCT1() == data.oildentCT1() &&
436 this->oildentCT2() == data.oildentCT2() &&
437 this->internalEnergyCurves() == data.internalEnergyCurves() &&
440 this->enableInternalEnergy() == data.enableInternalEnergy();
443 OilPvtThermal<Scalar>& operator=(
const OilPvtThermal<Scalar>& data)
445 if (data.isothermalPvt_)
446 isothermalPvt_ =
new IsothermalPvt(*data.isothermalPvt_);
448 isothermalPvt_ =
nullptr;
449 oilvisctCurves_ = data.oilvisctCurves_;
450 viscrefPress_ = data.viscrefPress_;
451 viscrefRs_ = data.viscrefRs_;
452 viscRef_ = data.viscRef_;
453 oildentRefTemp_ = data.oildentRefTemp_;
454 oildentCT1_ = data.oildentCT1_;
455 oildentCT2_ = data.oildentCT2_;
456 internalEnergyCurves_ = data.internalEnergyCurves_;
457 enableThermalDensity_ = data.enableThermalDensity_;
458 enableThermalViscosity_ = data.enableThermalViscosity_;
459 enableInternalEnergy_ = data.enableInternalEnergy_;
465 IsothermalPvt* isothermalPvt_;
469 std::vector<TabulatedOneDFunction> oilvisctCurves_;
470 std::vector<Scalar> viscrefPress_;
471 std::vector<Scalar> viscrefRs_;
472 std::vector<Scalar> viscRef_;
475 std::vector<Scalar> oildentRefTemp_;
476 std::vector<Scalar> oildentCT1_;
477 std::vector<Scalar> oildentCT2_;
480 std::vector<TabulatedOneDFunction> internalEnergyCurves_;
482 bool enableThermalDensity_;
483 bool enableThermalViscosity_;
484 bool enableInternalEnergy_;
A central place for various physical constants occuring in some equations.
This file provides a wrapper around the "final" C++-2011 statement.
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 in the black-oil mod...
Definition: OilPvtMultiplexer.hpp:96
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: OilPvtMultiplexer.hpp:177
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rs) const
Returns the formation volume factor [-] of the fluid phase.
Definition: OilPvtMultiplexer.hpp:219
Evaluation saturatedViscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: OilPvtMultiplexer.hpp:210
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: OilPvtMultiplexer.hpp:229
Evaluation saturatedGasDissolutionFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the gas dissolution factor [m^3/m^3] of saturated oil.
Definition: OilPvtMultiplexer.hpp:238
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rs) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: OilPvtMultiplexer.hpp:200
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &temperature, const Evaluation &Rs) const
Returns the saturation pressure [Pa] of oil given the mass fraction of the gas component in the oil p...
Definition: OilPvtMultiplexer.hpp:262
const Scalar oilReferenceDensity(unsigned regionIdx) const
Return the reference density which are considered by this PVT-object.
Definition: OilPvtMultiplexer.hpp:183
Evaluation diffusionCoefficient(const Evaluation &temperature, const Evaluation &pressure, unsigned compIdx) const
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition: OilPvtMultiplexer.hpp:271
This class implements temperature dependence of the PVT properties of oil.
Definition: OilPvtThermal.hpp:55
void initEnd()
Finish initializing the thermal part of the oil phase PVT properties.
Definition: OilPvtThermal.hpp:218
bool enableThermalDensity() const
Returns true iff the density of the oil phase is temperature dependent.
Definition: OilPvtThermal.hpp:224
Evaluation internalEnergy(unsigned regionIdx, const Evaluation &temperature, const Evaluation &, const Evaluation &) const
Returns the specific internal energy [J/kg] of oil given a set of parameters.
Definition: OilPvtThermal.hpp:240
Evaluation saturatedGasDissolutionFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: OilPvtThermal.hpp:347
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rs) const
Returns the formation volume factor [-] of the fluid phase.
Definition: OilPvtThermal.hpp:294
Evaluation saturatedViscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: OilPvtThermal.hpp:276
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rs) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: OilPvtThermal.hpp:258
void setNumRegions(size_t numRegions)
Set the number of PVT-regions considered by this object.
Definition: OilPvtThermal.hpp:206
bool enableThermalViscosity() const
Returns true iff the viscosity of the oil phase is temperature dependent.
Definition: OilPvtThermal.hpp:230
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the saturation pressure of the oil phase [Pa].
Definition: OilPvtThermal.hpp:375
Evaluation saturatedGasDissolutionFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &oilSaturation, const Evaluation &maxOilSaturation) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: OilPvtThermal.hpp:360
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of gas-saturated oil phase.
Definition: OilPvtThermal.hpp:319
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47