27 #ifndef OPM_WATER_PVT_THERMAL_HPP
28 #define OPM_WATER_PVT_THERMAL_HPP
38 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
39 #include <opm/input/eclipse/EclipseState/Tables/SimpleTable.hpp>
40 #include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
44 template <
class Scalar,
bool enableThermal,
bool enableBrine>
45 class WaterPvtMultiplexer;
53 template <
class Scalar,
bool enableBrine>
62 enableThermalDensity_ =
false;
63 enableThermalViscosity_ =
false;
64 enableInternalEnergy_ =
false;
65 isothermalPvt_ =
nullptr;
69 const std::vector<Scalar>& viscrefPress,
70 const std::vector<Scalar>& watdentRefTemp,
71 const std::vector<Scalar>& watdentCT1,
72 const std::vector<Scalar>& watdentCT2,
73 const std::vector<Scalar>& pvtwRefPress,
74 const std::vector<Scalar>& pvtwRefB,
75 const std::vector<Scalar>& pvtwCompressibility,
76 const std::vector<Scalar>& pvtwViscosity,
77 const std::vector<Scalar>& pvtwViscosibility,
78 const std::vector<TabulatedOneDFunction>& watvisctCurves,
79 const std::vector<TabulatedOneDFunction>& internalEnergyCurves,
82 bool enableInternalEnergy)
83 : isothermalPvt_(isothermalPvt)
84 , viscrefPress_(viscrefPress)
85 , watdentRefTemp_(watdentRefTemp)
86 , watdentCT1_(watdentCT1)
87 , watdentCT2_(watdentCT2)
88 , pvtwRefPress_(pvtwRefPress)
90 , pvtwCompressibility_(pvtwCompressibility)
91 , pvtwViscosity_(pvtwViscosity)
92 , pvtwViscosibility_(pvtwViscosibility)
93 , watvisctCurves_(watvisctCurves)
94 , internalEnergyCurves_(internalEnergyCurves)
97 , enableInternalEnergy_(enableInternalEnergy)
104 {
delete isothermalPvt_; }
110 void initFromState(
const EclipseState& eclState,
const Schedule& schedule)
116 isothermalPvt_->initFromState(eclState, schedule);
121 const auto& tables = eclState.getTableManager();
123 enableThermalDensity_ = tables.WatDenT().size() > 0;
124 enableThermalViscosity_ = tables.hasTables(
"WATVISCT");
125 enableInternalEnergy_ = tables.hasTables(
"SPECHEAT");
127 unsigned numRegions = isothermalPvt_->
numRegions();
130 if (enableThermalDensity_) {
131 const auto& watDenT = tables.WatDenT();
133 assert(watDenT.size() == numRegions);
134 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
135 const auto& record = watDenT[regionIdx];
137 watdentRefTemp_[regionIdx] = record.T0;
138 watdentCT1_[regionIdx] = record.C1;
139 watdentCT2_[regionIdx] = record.C2;
142 const auto& pvtwTables = tables.getPvtwTable();
144 assert(pvtwTables.size() == numRegions);
145 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
146 pvtwRefPress_[regionIdx] = pvtwTables[regionIdx].reference_pressure;
147 pvtwRefB_[regionIdx] = pvtwTables[regionIdx].volume_factor;
151 if (enableThermalViscosity_) {
152 if (tables.getViscrefTable().empty())
153 throw std::runtime_error(
"VISCREF is required when WATVISCT is present");
155 const auto& watvisctTables = tables.getWatvisctTables();
156 const auto& viscrefTables = tables.getViscrefTable();
158 const auto& pvtwTables = tables.getPvtwTable();
160 assert(pvtwTables.size() == numRegions);
161 assert(watvisctTables.size() == numRegions);
162 assert(viscrefTables.size() == numRegions);
164 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
165 const auto& T = watvisctTables[regionIdx].getColumn(
"Temperature").vectorCopy();
166 const auto& mu = watvisctTables[regionIdx].getColumn(
"Viscosity").vectorCopy();
167 watvisctCurves_[regionIdx].setXYContainers(T, mu);
169 viscrefPress_[regionIdx] = viscrefTables[regionIdx].reference_pressure;
172 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
173 pvtwViscosity_[regionIdx] = pvtwTables[regionIdx].viscosity;
174 pvtwViscosibility_[regionIdx] = pvtwTables[regionIdx].viscosibility;
178 if (enableInternalEnergy_) {
182 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
183 const auto& specHeatTable = tables.getSpecheatTables()[regionIdx];
184 const auto& temperatureColumn = specHeatTable.getColumn(
"TEMPERATURE");
185 const auto& cvWaterColumn = specHeatTable.getColumn(
"CV_WATER");
187 std::vector<double> uSamples(temperatureColumn.size());
189 Scalar u = temperatureColumn[0]*cvWaterColumn[0];
190 for (
size_t i = 0;; ++i) {
193 if (i >= temperatureColumn.size() - 1)
198 Scalar c_v0 = cvWaterColumn[i];
199 Scalar c_v1 = cvWaterColumn[i + 1];
200 Scalar T0 = temperatureColumn[i];
201 Scalar T1 = temperatureColumn[i + 1];
202 u += 0.5*(c_v0 + c_v1)*(T1 - T0);
205 internalEnergyCurves_[regionIdx].setXYContainers(temperatureColumn.vectorCopy(), uSamples);
216 pvtwRefPress_.resize(numRegions);
217 pvtwRefB_.resize(numRegions);
218 pvtwCompressibility_.resize(numRegions);
219 pvtwViscosity_.resize(numRegions);
220 pvtwViscosibility_.resize(numRegions);
221 viscrefPress_.resize(numRegions);
222 watvisctCurves_.resize(numRegions);
223 watdentRefTemp_.resize(numRegions);
224 watdentCT1_.resize(numRegions);
225 watdentCT2_.resize(numRegions);
226 internalEnergyCurves_.resize(numRegions);
239 {
return enableThermalDensity_; }
245 {
return enableThermalViscosity_; }
247 size_t numRegions()
const
248 {
return pvtwRefPress_.size(); }
253 template <
class Evaluation>
255 const Evaluation& temperature,
256 const Evaluation&)
const
258 if (!enableInternalEnergy_)
259 throw std::runtime_error(
"Requested the internal energy of oil but it is disabled");
264 return internalEnergyCurves_[regionIdx].eval(temperature,
true);
270 template <
class Evaluation>
272 const Evaluation& temperature,
273 const Evaluation& pressure,
274 const Evaluation& saltconcentration)
const
276 const auto& isothermalMu = isothermalPvt_->
viscosity(regionIdx, temperature, pressure, saltconcentration);
280 Scalar x = -pvtwViscosibility_[regionIdx]*(viscrefPress_[regionIdx] - pvtwRefPress_[regionIdx]);
281 Scalar muRef = pvtwViscosity_[regionIdx]/(1.0 + x + 0.5*x*x);
284 const auto& muWatvisct = watvisctCurves_[regionIdx].eval(temperature,
true);
285 return isothermalMu * muWatvisct/muRef;
291 template <
class Evaluation>
293 const Evaluation& temperature,
294 const Evaluation& pressure,
295 const Evaluation& saltconcentration)
const
300 Scalar BwRef = pvtwRefB_[regionIdx];
301 Scalar TRef = watdentRefTemp_[regionIdx];
302 const Evaluation& X = pvtwCompressibility_[regionIdx]*(pressure - pvtwRefPress_[regionIdx]);
303 Scalar cT1 = watdentCT1_[regionIdx];
304 Scalar cT2 = watdentCT2_[regionIdx];
305 const Evaluation& Y = temperature - TRef;
310 return 1.0/(((1 - X)*(1 + cT1*Y + cT2*Y*Y))*BwRef);
313 const IsothermalPvt* isoThermalPvt()
const
314 {
return isothermalPvt_; }
316 const Scalar waterReferenceDensity(
unsigned regionIdx)
const
319 const std::vector<Scalar>& viscrefPress()
const
320 {
return viscrefPress_; }
322 const std::vector<Scalar>& watdentRefTemp()
const
323 {
return watdentRefTemp_; }
325 const std::vector<Scalar>& watdentCT1()
const
326 {
return watdentCT1_; }
328 const std::vector<Scalar>& watdentCT2()
const
329 {
return watdentCT2_; }
331 const std::vector<Scalar>& pvtwRefPress()
const
332 {
return pvtwRefPress_; }
334 const std::vector<Scalar>& pvtwRefB()
const
335 {
return pvtwRefB_; }
337 const std::vector<Scalar>& pvtwCompressibility()
const
338 {
return pvtwCompressibility_; }
340 const std::vector<Scalar>& pvtwViscosity()
const
341 {
return pvtwViscosity_; }
343 const std::vector<Scalar>& pvtwViscosibility()
const
344 {
return pvtwViscosibility_; }
346 const std::vector<TabulatedOneDFunction>& watvisctCurves()
const
347 {
return watvisctCurves_; }
349 const std::vector<TabulatedOneDFunction> internalEnergyCurves()
const
350 {
return internalEnergyCurves_; }
352 bool enableInternalEnergy()
const
353 {
return enableInternalEnergy_; }
355 bool operator==(
const WaterPvtThermal<Scalar, enableBrine>& data)
const
357 if (isothermalPvt_ && !data.isothermalPvt_)
359 if (!isothermalPvt_ && data.isothermalPvt_)
362 return (!this->isoThermalPvt() ||
363 (*this->isoThermalPvt() == *data.isoThermalPvt())) &&
364 this->viscrefPress() == data.viscrefPress() &&
365 this->watdentRefTemp() == data.watdentRefTemp() &&
366 this->watdentCT1() == data.watdentCT1() &&
367 this->watdentCT2() == data.watdentCT2() &&
368 this->pvtwRefPress() == data.pvtwRefPress() &&
369 this->pvtwRefB() == data.pvtwRefB() &&
370 this->pvtwCompressibility() == data.pvtwCompressibility() &&
371 this->pvtwViscosity() == data.pvtwViscosity() &&
372 this->pvtwViscosibility() == data.pvtwViscosibility() &&
373 this->watvisctCurves() == data.watvisctCurves() &&
374 this->internalEnergyCurves() == data.internalEnergyCurves() &&
377 this->enableInternalEnergy() == data.enableInternalEnergy();
380 WaterPvtThermal<Scalar, enableBrine>& operator=(
const WaterPvtThermal<Scalar, enableBrine>& data)
382 if (data.isothermalPvt_)
383 isothermalPvt_ =
new IsothermalPvt(*data.isothermalPvt_);
385 isothermalPvt_ =
nullptr;
386 viscrefPress_ = data.viscrefPress_;
387 watdentRefTemp_ = data.watdentRefTemp_;
388 watdentCT1_ = data.watdentCT1_;
389 watdentCT2_ = data.watdentCT2_;
390 pvtwRefPress_ = data.pvtwRefPress_;
391 pvtwRefB_ = data.pvtwRefB_;
392 pvtwCompressibility_ = data.pvtwCompressibility_;
393 pvtwViscosity_ = data.pvtwViscosity_;
394 pvtwViscosibility_ = data.pvtwViscosibility_;
395 watvisctCurves_ = data.watvisctCurves_;
396 internalEnergyCurves_ = data.internalEnergyCurves_;
397 enableThermalDensity_ = data.enableThermalDensity_;
398 enableThermalViscosity_ = data.enableThermalViscosity_;
399 enableInternalEnergy_ = data.enableInternalEnergy_;
405 IsothermalPvt* isothermalPvt_;
409 std::vector<Scalar> viscrefPress_;
411 std::vector<Scalar> watdentRefTemp_;
412 std::vector<Scalar> watdentCT1_;
413 std::vector<Scalar> watdentCT2_;
415 std::vector<Scalar> pvtwRefPress_;
416 std::vector<Scalar> pvtwRefB_;
417 std::vector<Scalar> pvtwCompressibility_;
418 std::vector<Scalar> pvtwViscosity_;
419 std::vector<Scalar> pvtwViscosibility_;
421 std::vector<TabulatedOneDFunction> watvisctCurves_;
424 std::vector<TabulatedOneDFunction> internalEnergyCurves_;
426 bool enableThermalDensity_;
427 bool enableThermalViscosity_;
428 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.
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition: WaterPvtMultiplexer.hpp:75
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: WaterPvtMultiplexer.hpp:141
const Scalar waterReferenceDensity(unsigned regionIdx)
Return the reference density which are considered by this PVT-object.
Definition: WaterPvtMultiplexer.hpp:147
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &saltconcentration) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WaterPvtMultiplexer.hpp:176
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &saltconcentration) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WaterPvtMultiplexer.hpp:163
This class implements temperature dependence of the PVT properties of water.
Definition: WaterPvtThermal.hpp:55
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &saltconcentration) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WaterPvtThermal.hpp:292
bool enableThermalViscosity() const
Returns true iff the viscosity of the water phase is temperature dependent.
Definition: WaterPvtThermal.hpp:244
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &saltconcentration) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WaterPvtThermal.hpp:271
Evaluation internalEnergy(unsigned regionIdx, const Evaluation &temperature, const Evaluation &) const
Returns the specific internal energy [J/kg] of water given a set of parameters.
Definition: WaterPvtThermal.hpp:254
void initEnd()
Finish initializing the thermal part of the water phase PVT properties.
Definition: WaterPvtThermal.hpp:232
bool enableThermalDensity() const
Returns true iff the density of the water phase is temperature dependent.
Definition: WaterPvtThermal.hpp:238
void setNumRegions(size_t numRegions)
Set the number of PVT-regions considered by this object.
Definition: WaterPvtThermal.hpp:214