27#ifndef OPM_WATER_PVT_THERMAL_HPP
28#define OPM_WATER_PVT_THERMAL_HPP
39template <
class Scalar,
bool enableThermal,
bool enableBrine>
40class WaterPvtMultiplexer;
48template <
class Scalar,
bool enableBrine>
57 enableThermalDensity_ =
false;
58 enableThermalViscosity_ =
false;
59 enableInternalEnergy_ =
false;
60 isothermalPvt_ =
nullptr;
64 const std::vector<Scalar>& viscrefPress,
65 const std::vector<Scalar>& watdentRefTemp,
66 const std::vector<Scalar>& watdentCT1,
67 const std::vector<Scalar>& watdentCT2,
68 const std::vector<Scalar>& watJTRefPres,
69 const std::vector<Scalar>& watJTC,
70 const std::vector<Scalar>& pvtwRefPress,
71 const std::vector<Scalar>& pvtwRefB,
72 const std::vector<Scalar>& pvtwCompressibility,
73 const std::vector<Scalar>& pvtwViscosity,
74 const std::vector<Scalar>& pvtwViscosibility,
75 const std::vector<TabulatedOneDFunction>& watvisctCurves,
76 const std::vector<TabulatedOneDFunction>& internalEnergyCurves,
80 bool enableInternalEnergy)
81 : isothermalPvt_(isothermalPvt)
82 , viscrefPress_(viscrefPress)
83 , watdentRefTemp_(watdentRefTemp)
84 , watdentCT1_(watdentCT1)
85 , watdentCT2_(watdentCT2)
86 , watJTRefPres_(watJTRefPres)
88 , pvtwRefPress_(pvtwRefPress)
90 , pvtwCompressibility_(pvtwCompressibility)
91 , pvtwViscosity_(pvtwViscosity)
92 , pvtwViscosibility_(pvtwViscosibility)
93 , watvisctCurves_(watvisctCurves)
94 , internalEnergyCurves_(internalEnergyCurves)
98 , enableInternalEnergy_(enableInternalEnergy)
105 {
delete isothermalPvt_; }
119 pvtwRefPress_.resize(numRegions);
120 pvtwRefB_.resize(numRegions);
121 pvtwCompressibility_.resize(numRegions);
122 pvtwViscosity_.resize(numRegions);
123 pvtwViscosibility_.resize(numRegions);
124 viscrefPress_.resize(numRegions);
125 watvisctCurves_.resize(numRegions);
126 watdentRefTemp_.resize(numRegions);
127 watdentCT1_.resize(numRegions);
128 watdentCT2_.resize(numRegions);
129 watJTRefPres_.resize(numRegions);
130 watJTC_.resize(numRegions);
131 internalEnergyCurves_.resize(numRegions);
144 {
return enableThermalDensity_; }
150 {
return enableJouleThomson_; }
156 {
return enableThermalViscosity_; }
158 size_t numRegions()
const
159 {
return pvtwRefPress_.size(); }
164 template <
class Evaluation>
166 const Evaluation& temperature,
167 const Evaluation& pressure,
168 const Evaluation& Rsw,
169 const Evaluation& saltconcentration)
const
171 if (!enableInternalEnergy_)
172 throw std::runtime_error(
"Requested the internal energy of water but it is disabled");
174 if (!enableJouleThomson_) {
178 return internalEnergyCurves_[regionIdx].eval(temperature,
true);
181 Evaluation Tref = watdentRefTemp_[regionIdx];
182 Evaluation Pref = watJTRefPres_[regionIdx];
183 Scalar JTC =watJTC_[regionIdx];
186 Evaluation Cp = internalEnergyCurves_[regionIdx].eval(temperature,
true)/temperature;
187 Evaluation density = invB * waterReferenceDensity(regionIdx);
189 Evaluation enthalpyPres;
191 enthalpyPres = -Cp * JTC * (pressure -Pref);
193 else if(enableThermalDensity_) {
194 Scalar c1T = watdentCT1_[regionIdx];
195 Scalar c2T = watdentCT2_[regionIdx];
197 Evaluation alpha = (c1T + 2 * c2T * (temperature - Tref)) /
198 (1 + c1T *(temperature - Tref) + c2T * (temperature - Tref) * (temperature - Tref));
200 constexpr const int N = 100;
201 Evaluation deltaP = (pressure - Pref)/N;
202 Evaluation enthalpyPresPrev = 0;
203 for (
size_t i = 0; i < N; ++i) {
204 Evaluation Pnew = Pref + i * deltaP;
206 Evaluation jouleThomsonCoefficient = -(1.0/Cp) * (1.0 - alpha * temperature)/rho;
207 Evaluation deltaEnthalpyPres = -Cp * jouleThomsonCoefficient * deltaP;
208 enthalpyPres = enthalpyPresPrev + deltaEnthalpyPres;
209 enthalpyPresPrev = enthalpyPres;
213 throw std::runtime_error(
"Requested Joule-thomson calculation but thermal water density (WATDENT) is not provided");
216 Evaluation enthalpy = Cp * (temperature - Tref) + enthalpyPres;
218 return enthalpy - pressure/density;
225 template <
class Evaluation>
227 const Evaluation& temperature,
228 const Evaluation& pressure,
229 const Evaluation& Rsw,
230 const Evaluation& saltconcentration)
const
232 const auto& isothermalMu = isothermalPvt_->
viscosity(regionIdx, temperature, pressure, Rsw, saltconcentration);
236 Scalar x = -pvtwViscosibility_[regionIdx]*(viscrefPress_[regionIdx] - pvtwRefPress_[regionIdx]);
237 Scalar muRef = pvtwViscosity_[regionIdx]/(1.0 + x + 0.5*x*x);
240 const auto& muWatvisct = watvisctCurves_[regionIdx].eval(temperature,
true);
241 return isothermalMu * muWatvisct/muRef;
247 template <
class Evaluation>
249 const Evaluation& temperature,
250 const Evaluation& pressure,
251 const Evaluation& saltconcentration)
const
253 const auto& isothermalMu = isothermalPvt_->
saturatedViscosity(regionIdx, temperature, pressure, saltconcentration);
257 Scalar x = -pvtwViscosibility_[regionIdx]*(viscrefPress_[regionIdx] - pvtwRefPress_[regionIdx]);
258 Scalar muRef = pvtwViscosity_[regionIdx]/(1.0 + x + 0.5*x*x);
261 const auto& muWatvisct = watvisctCurves_[regionIdx].eval(temperature,
true);
262 return isothermalMu * muWatvisct/muRef;
268 template <
class Evaluation>
270 const Evaluation& temperature,
271 const Evaluation& pressure,
272 const Evaluation& saltconcentration)
const
274 Evaluation Rsw = 0.0;
280 template <
class Evaluation>
282 const Evaluation& temperature,
283 const Evaluation& pressure,
284 const Evaluation& Rsw,
285 const Evaluation& saltconcentration)
const
290 Scalar BwRef = pvtwRefB_[regionIdx];
291 Scalar TRef = watdentRefTemp_[regionIdx];
292 const Evaluation& X = pvtwCompressibility_[regionIdx]*(pressure - pvtwRefPress_[regionIdx]);
293 Scalar cT1 = watdentCT1_[regionIdx];
294 Scalar cT2 = watdentCT2_[regionIdx];
295 const Evaluation& Y = temperature - TRef;
300 return 1.0/(((1 - X)*(1 + cT1*Y + cT2*Y*Y))*BwRef);
309 template <
class Evaluation>
313 const Evaluation& )
const
319 template <
class Evaluation>
323 const Evaluation& )
const
326 template <
class Evaluation>
327 Evaluation diffusionCoefficient(
const Evaluation& ,
331 throw std::runtime_error(
"Not implemented: The PVT model does not provide a diffusionCoefficient()");
334 const IsothermalPvt* isoThermalPvt()
const
335 {
return isothermalPvt_; }
337 const Scalar waterReferenceDensity(
unsigned regionIdx)
const
340 const std::vector<Scalar>& viscrefPress()
const
341 {
return viscrefPress_; }
343 const std::vector<Scalar>& watdentRefTemp()
const
344 {
return watdentRefTemp_; }
346 const std::vector<Scalar>& watdentCT1()
const
347 {
return watdentCT1_; }
349 const std::vector<Scalar>& watdentCT2()
const
350 {
return watdentCT2_; }
352 const std::vector<Scalar>& pvtwRefPress()
const
353 {
return pvtwRefPress_; }
355 const std::vector<Scalar>& pvtwRefB()
const
356 {
return pvtwRefB_; }
358 const std::vector<Scalar>& pvtwCompressibility()
const
359 {
return pvtwCompressibility_; }
361 const std::vector<Scalar>& pvtwViscosity()
const
362 {
return pvtwViscosity_; }
364 const std::vector<Scalar>& pvtwViscosibility()
const
365 {
return pvtwViscosibility_; }
367 const std::vector<TabulatedOneDFunction>& watvisctCurves()
const
368 {
return watvisctCurves_; }
370 const std::vector<TabulatedOneDFunction> internalEnergyCurves()
const
371 {
return internalEnergyCurves_; }
373 bool enableInternalEnergy()
const
374 {
return enableInternalEnergy_; }
376 const std::vector<Scalar>& watJTRefPres()
const
377 {
return watJTRefPres_; }
379 const std::vector<Scalar>& watJTC()
const
383 bool operator==(
const WaterPvtThermal<Scalar, enableBrine>& data)
const
385 if (isothermalPvt_ && !data.isothermalPvt_)
387 if (!isothermalPvt_ && data.isothermalPvt_)
390 return this->viscrefPress() == data.viscrefPress() &&
391 this->watdentRefTemp() == data.watdentRefTemp() &&
392 this->watdentCT1() == data.watdentCT1() &&
393 this->watdentCT2() == data.watdentCT2() &&
394 this->watdentCT2() == data.watdentCT2() &&
395 this->watJTRefPres() == data.watJTRefPres() &&
396 this->watJTC() == data.watJTC() &&
397 this->pvtwRefPress() == data.pvtwRefPress() &&
398 this->pvtwRefB() == data.pvtwRefB() &&
399 this->pvtwCompressibility() == data.pvtwCompressibility() &&
400 this->pvtwViscosity() == data.pvtwViscosity() &&
401 this->pvtwViscosibility() == data.pvtwViscosibility() &&
402 this->watvisctCurves() == data.watvisctCurves() &&
403 this->internalEnergyCurves() == data.internalEnergyCurves() &&
407 this->enableInternalEnergy() == data.enableInternalEnergy();
410 WaterPvtThermal<Scalar, enableBrine>& operator=(
const WaterPvtThermal<Scalar, enableBrine>& data)
412 if (data.isothermalPvt_)
413 isothermalPvt_ =
new IsothermalPvt(*data.isothermalPvt_);
415 isothermalPvt_ =
nullptr;
416 viscrefPress_ = data.viscrefPress_;
417 watdentRefTemp_ = data.watdentRefTemp_;
418 watdentCT1_ = data.watdentCT1_;
419 watdentCT2_ = data.watdentCT2_;
420 watJTRefPres_ = data.watJTRefPres_;
421 watJTC_ = data.watJTC_;
422 pvtwRefPress_ = data.pvtwRefPress_;
423 pvtwRefB_ = data.pvtwRefB_;
424 pvtwCompressibility_ = data.pvtwCompressibility_;
425 pvtwViscosity_ = data.pvtwViscosity_;
426 pvtwViscosibility_ = data.pvtwViscosibility_;
427 watvisctCurves_ = data.watvisctCurves_;
428 internalEnergyCurves_ = data.internalEnergyCurves_;
429 enableThermalDensity_ = data.enableThermalDensity_;
430 enableJouleThomson_ = data.enableJouleThomson_;
431 enableThermalViscosity_ = data.enableThermalViscosity_;
432 enableInternalEnergy_ = data.enableInternalEnergy_;
438 IsothermalPvt* isothermalPvt_;
442 std::vector<Scalar> viscrefPress_;
444 std::vector<Scalar> watdentRefTemp_;
445 std::vector<Scalar> watdentCT1_;
446 std::vector<Scalar> watdentCT2_;
448 std::vector<Scalar> watJTRefPres_;
449 std::vector<Scalar> watJTC_;
451 std::vector<Scalar> pvtwRefPress_;
452 std::vector<Scalar> pvtwRefB_;
453 std::vector<Scalar> pvtwCompressibility_;
454 std::vector<Scalar> pvtwViscosity_;
455 std::vector<Scalar> pvtwViscosibility_;
457 std::vector<TabulatedOneDFunction> watvisctCurves_;
460 std::vector<TabulatedOneDFunction> internalEnergyCurves_;
462 bool enableThermalDensity_;
463 bool enableJouleThomson_;
464 bool enableThermalViscosity_;
465 bool enableInternalEnergy_;
Definition: EclipseState.hpp:55
Definition: Schedule.hpp:130
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:51
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition: WaterPvtMultiplexer.hpp:82
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rsw, const Evaluation &saltconcentration) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WaterPvtMultiplexer.hpp:163
Evaluation saturatedViscosity(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:177
const Scalar waterReferenceDensity(unsigned regionIdx) const
Return the reference density which are considered by this PVT-object.
Definition: WaterPvtMultiplexer.hpp:145
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rsw, const Evaluation &saltconcentration) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WaterPvtMultiplexer.hpp:190
This class implements temperature dependence of the PVT properties of water.
Definition: WaterPvtThermal.hpp:50
bool enableThermalViscosity() const
Returns true iff the viscosity of the water phase is temperature dependent.
Definition: WaterPvtThermal.hpp:155
Evaluation saturationPressure(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the saturation pressure of the water phase [Pa] depending on its mass fraction of the gas com...
Definition: WaterPvtThermal.hpp:310
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rsw, const Evaluation &saltconcentration) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WaterPvtThermal.hpp:226
void initEnd()
Finish initializing the thermal part of the water phase PVT properties.
Definition: WaterPvtThermal.hpp:137
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rsw, const Evaluation &saltconcentration) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WaterPvtThermal.hpp:281
bool enableThermalDensity() const
Returns true iff the density of the water phase is temperature dependent.
Definition: WaterPvtThermal.hpp:143
Evaluation saturatedGasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the gas dissolution factor [m^3/m^3] of the water phase.
Definition: WaterPvtThermal.hpp:320
Evaluation internalEnergy(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rsw, const Evaluation &saltconcentration) const
Returns the specific internal energy [J/kg] of water given a set of parameters.
Definition: WaterPvtThermal.hpp:165
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &saltconcentration) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WaterPvtThermal.hpp:269
bool enableJouleThomson() const
Returns true iff Joule-Thomson effect for the water phase is active.
Definition: WaterPvtThermal.hpp:149
Evaluation saturatedViscosity(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:248
void setNumRegions(size_t numRegions)
Set the number of PVT-regions considered by this object.
Definition: WaterPvtThermal.hpp:117
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:30