My Project
OilPvtThermal.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_OIL_PVT_THERMAL_HPP
28#define OPM_OIL_PVT_THERMAL_HPP
29
31
32namespace Opm {
33
34#if HAVE_ECL_INPUT
35class EclipseState;
36class Schedule;
37#endif
38
39template <class Scalar, bool enableThermal>
40class OilPvtMultiplexer;
41
48template <class Scalar>
50{
51public:
53 using IsothermalPvt = OilPvtMultiplexer<Scalar, /*enableThermal=*/false>;
54
56 {
57 enableThermalDensity_ = false;
58 enableJouleThomson_ = false;
59 enableThermalViscosity_ = false;
60 enableInternalEnergy_ = false;
61 isothermalPvt_ = nullptr;
62 }
63
64 OilPvtThermal(IsothermalPvt* isothermalPvt,
65 const std::vector<TabulatedOneDFunction>& oilvisctCurves,
66 const std::vector<Scalar>& viscrefPress,
67 const std::vector<Scalar>& viscrefRs,
68 const std::vector<Scalar>& viscRef,
69 const std::vector<Scalar>& oildentRefTemp,
70 const std::vector<Scalar>& oildentCT1,
71 const std::vector<Scalar>& oildentCT2,
72 const std::vector<Scalar>& oilJTRefPres,
73 const std::vector<Scalar>& oilJTC,
74 const std::vector<TabulatedOneDFunction>& internalEnergyCurves,
78 bool enableInternalEnergy)
79 : isothermalPvt_(isothermalPvt)
80 , oilvisctCurves_(oilvisctCurves)
81 , viscrefPress_(viscrefPress)
82 , viscrefRs_(viscrefRs)
83 , viscRef_(viscRef)
84 , oildentRefTemp_(oildentRefTemp)
85 , oildentCT1_(oildentCT1)
86 , oildentCT2_(oildentCT2)
87 , oilJTRefPres_(oilJTRefPres)
88 , oilJTC_(oilJTC)
89 , internalEnergyCurves_(internalEnergyCurves)
90 , enableThermalDensity_(enableThermalDensity)
91 , enableJouleThomson_(enableJouleThomson)
92 , enableThermalViscosity_(enableThermalViscosity)
93 , enableInternalEnergy_(enableInternalEnergy)
94 { }
95
96 OilPvtThermal(const OilPvtThermal& data)
97 { *this = data; }
98
100 { delete isothermalPvt_; }
101
102#if HAVE_ECL_INPUT
106 void initFromState(const EclipseState& eclState, const Schedule& schedule);
107#endif // HAVE_ECL_INPUT
108
112 void setNumRegions(size_t numRegions)
113 {
114 oilvisctCurves_.resize(numRegions);
115 viscrefPress_.resize(numRegions);
116 viscrefRs_.resize(numRegions);
117 viscRef_.resize(numRegions);
118 internalEnergyCurves_.resize(numRegions);
119 oildentRefTemp_.resize(numRegions);
120 oildentCT1_.resize(numRegions);
121 oildentCT2_.resize(numRegions);
122 oilJTRefPres_.resize(numRegions);
123 oilJTC_.resize(numRegions);
124 rhoRefG_.resize(numRegions);
125 }
126
130 void initEnd()
131 { }
132
137 { return enableThermalDensity_; }
138
143 { return enableJouleThomson_; }
144
149 { return enableThermalViscosity_; }
150
151 size_t numRegions() const
152 { return viscrefRs_.size(); }
153
157 template <class Evaluation>
158 Evaluation internalEnergy(unsigned regionIdx,
159 const Evaluation& temperature,
160 const Evaluation& pressure,
161 const Evaluation& Rs) const
162 {
163 if (!enableInternalEnergy_)
164 throw std::runtime_error("Requested the internal energy of oil but it is disabled");
165
166 if (!enableJouleThomson_) {
167 // compute the specific internal energy for the specified tempature. We use linear
168 // interpolation here despite the fact that the underlying heat capacities are
169 // piecewise linear (which leads to a quadratic function)
170 return internalEnergyCurves_[regionIdx].eval(temperature, /*extrapolate=*/true);
171 }
172 else {
173 Evaluation Tref = oildentRefTemp_[regionIdx];
174 Evaluation Pref = oilJTRefPres_[regionIdx];
175 Scalar JTC = oilJTC_[regionIdx]; // if JTC is default then JTC is calculated
176
177 Evaluation invB = inverseFormationVolumeFactor(regionIdx, temperature, pressure, Rs);
178 Evaluation Cp = internalEnergyCurves_[regionIdx].eval(temperature, /*extrapolate=*/true)/temperature;
179 Evaluation density = invB * (oilReferenceDensity(regionIdx) + Rs * rhoRefG_[regionIdx]);
180
181 Evaluation enthalpyPres;
182 if (JTC != 0) {
183 enthalpyPres = -Cp * JTC * (pressure -Pref);
184 }
185 else if(enableThermalDensity_) {
186 Scalar c1T = oildentCT1_[regionIdx];
187 Scalar c2T = oildentCT2_[regionIdx];
188
189 Evaluation alpha = (c1T + 2 * c2T * (temperature - Tref)) /
190 (1 + c1T *(temperature - Tref) + c2T * (temperature - Tref) * (temperature - Tref));
191
192 const int N = 100; // value is experimental
193 Evaluation deltaP = (pressure - Pref)/N;
194 Evaluation enthalpyPresPrev = 0;
195 for (size_t i = 0; i < N; ++i) {
196 Evaluation Pnew = Pref + i * deltaP;
197 Evaluation rho = inverseFormationVolumeFactor(regionIdx, temperature, Pnew, Rs) *
198 (oilReferenceDensity(regionIdx) + Rs * rhoRefG_[regionIdx]) ;
199 // see e.g.https://en.wikipedia.org/wiki/Joule-Thomson_effect for a derivation of the Joule-Thomson coeff.
200 Evaluation jouleThomsonCoefficient = -(1.0/Cp) * (1.0 - alpha * temperature)/rho;
201 Evaluation deltaEnthalpyPres = -Cp * jouleThomsonCoefficient * deltaP;
202 enthalpyPres = enthalpyPresPrev + deltaEnthalpyPres;
203 enthalpyPresPrev = enthalpyPres;
204 }
205 }
206 else {
207 throw std::runtime_error("Requested Joule-thomson calculation but thermal oil density (OILDENT) is not provided");
208 }
209
210 Evaluation enthalpy = Cp * (temperature - Tref) + enthalpyPres;
211
212 return enthalpy - pressure/density;
213 }
214 }
215
219 template <class Evaluation>
220 Evaluation viscosity(unsigned regionIdx,
221 const Evaluation& temperature,
222 const Evaluation& pressure,
223 const Evaluation& Rs) const
224 {
225 const auto& isothermalMu = isothermalPvt_->viscosity(regionIdx, temperature, pressure, Rs);
227 return isothermalMu;
228
229 // compute the viscosity deviation due to temperature
230 const auto& muOilvisct = oilvisctCurves_[regionIdx].eval(temperature, /*extrapolate=*/true);
231 return muOilvisct/viscRef_[regionIdx]*isothermalMu;
232 }
233
237 template <class Evaluation>
238 Evaluation saturatedViscosity(unsigned regionIdx,
239 const Evaluation& temperature,
240 const Evaluation& pressure) const
241 {
242 const auto& isothermalMu = isothermalPvt_->saturatedViscosity(regionIdx, temperature, pressure);
244 return isothermalMu;
245
246 // compute the viscosity deviation due to temperature
247 const auto& muOilvisct = oilvisctCurves_[regionIdx].eval(temperature, /*extrapolate=*/true);
248 return muOilvisct/viscRef_[regionIdx]*isothermalMu;
249 }
250
251
255 template <class Evaluation>
256 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
257 const Evaluation& temperature,
258 const Evaluation& pressure,
259 const Evaluation& Rs) const
260 {
261 const auto& b =
262 isothermalPvt_->inverseFormationVolumeFactor(regionIdx, temperature, pressure, Rs);
263
265 return b;
266
267 // we use the same approach as for the for water here, but with the OPM-specific
268 // OILDENT keyword.
269 Scalar TRef = oildentRefTemp_[regionIdx];
270 Scalar cT1 = oildentCT1_[regionIdx];
271 Scalar cT2 = oildentCT2_[regionIdx];
272 const Evaluation& Y = temperature - TRef;
273
274 return b/(1 + (cT1 + cT2*Y)*Y);
275 }
276
280 template <class Evaluation>
281 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
282 const Evaluation& temperature,
283 const Evaluation& pressure) const
284 {
285 const auto& b =
286 isothermalPvt_->saturatedInverseFormationVolumeFactor(regionIdx, temperature, pressure);
287
289 return b;
290
291 // we use the same approach as for the for water here, but with the OPM-specific
292 // OILDENT keyword.
293 Scalar TRef = oildentRefTemp_[regionIdx];
294 Scalar cT1 = oildentCT1_[regionIdx];
295 Scalar cT2 = oildentCT2_[regionIdx];
296 const Evaluation& Y = temperature - TRef;
297
298 return b/(1 + (cT1 + cT2*Y)*Y);
299 }
300
308 template <class Evaluation>
309 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
310 const Evaluation& temperature,
311 const Evaluation& pressure) const
312 { return isothermalPvt_->saturatedGasDissolutionFactor(regionIdx, temperature, pressure); }
313
321 template <class Evaluation>
322 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
323 const Evaluation& temperature,
324 const Evaluation& pressure,
325 const Evaluation& oilSaturation,
326 const Evaluation& maxOilSaturation) const
327 { return isothermalPvt_->saturatedGasDissolutionFactor(regionIdx, temperature, pressure, oilSaturation, maxOilSaturation); }
328
336 template <class Evaluation>
337 Evaluation saturationPressure(unsigned regionIdx,
338 const Evaluation& temperature,
339 const Evaluation& pressure) const
340 { return isothermalPvt_->saturationPressure(regionIdx, temperature, pressure); }
341
342 template <class Evaluation>
343 Evaluation diffusionCoefficient(const Evaluation& temperature,
344 const Evaluation& pressure,
345 unsigned compIdx) const
346 {
347 return isothermalPvt_->diffusionCoefficient(temperature, pressure, compIdx);
348 }
349
350 const IsothermalPvt* isoThermalPvt() const
351 { return isothermalPvt_; }
352
353 const Scalar oilReferenceDensity(unsigned regionIdx) const
354 { return isothermalPvt_->oilReferenceDensity(regionIdx); }
355
356 const std::vector<TabulatedOneDFunction>& oilvisctCurves() const
357 { return oilvisctCurves_; }
358
359 const std::vector<Scalar>& viscrefPress() const
360 { return viscrefPress_; }
361
362 const std::vector<Scalar>& viscrefRs() const
363 { return viscrefRs_; }
364
365 const std::vector<Scalar>& viscRef() const
366 { return viscRef_; }
367
368 const std::vector<Scalar>& oildentRefTemp() const
369 { return oildentRefTemp_; }
370
371 const std::vector<Scalar>& oildentCT1() const
372 { return oildentCT1_; }
373
374 const std::vector<Scalar>& oildentCT2() const
375 { return oildentCT2_; }
376
377 const std::vector<TabulatedOneDFunction> internalEnergyCurves() const
378 { return internalEnergyCurves_; }
379
380 bool enableInternalEnergy() const
381 { return enableInternalEnergy_; }
382
383 const std::vector<Scalar>& oilJTRefPres() const
384 { return oilJTRefPres_; }
385
386 const std::vector<Scalar>& oilJTC() const
387 { return oilJTC_; }
388
389 bool operator==(const OilPvtThermal<Scalar>& data) const
390 {
391 if (isothermalPvt_ && !data.isothermalPvt_)
392 return false;
393 if (!isothermalPvt_ && data.isothermalPvt_)
394 return false;
395
396 return this->oilvisctCurves() == data.oilvisctCurves() &&
397 this->viscrefPress() == data.viscrefPress() &&
398 this->viscrefRs() == data.viscrefRs() &&
399 this->viscRef() == data.viscRef() &&
400 this->oildentRefTemp() == data.oildentRefTemp() &&
401 this->oildentCT1() == data.oildentCT1() &&
402 this->oildentCT2() == data.oildentCT2() &&
403 this->oilJTRefPres() == data.oilJTRefPres() &&
404 this->oilJTC() == data.oilJTC() &&
405 this->internalEnergyCurves() == data.internalEnergyCurves() &&
406 this->enableThermalDensity() == data.enableThermalDensity() &&
407 this->enableJouleThomson() == data.enableJouleThomson() &&
408 this->enableThermalViscosity() == data.enableThermalViscosity() &&
409 this->enableInternalEnergy() == data.enableInternalEnergy();
410 }
411
412 OilPvtThermal<Scalar>& operator=(const OilPvtThermal<Scalar>& data)
413 {
414 if (data.isothermalPvt_)
415 isothermalPvt_ = new IsothermalPvt(*data.isothermalPvt_);
416 else
417 isothermalPvt_ = nullptr;
418 oilvisctCurves_ = data.oilvisctCurves_;
419 viscrefPress_ = data.viscrefPress_;
420 viscrefRs_ = data.viscrefRs_;
421 viscRef_ = data.viscRef_;
422 oildentRefTemp_ = data.oildentRefTemp_;
423 oildentCT1_ = data.oildentCT1_;
424 oildentCT2_ = data.oildentCT2_;
425 oilJTRefPres_ = data.oilJTRefPres_;
426 oilJTC_ = data.oilJTC_;
427 internalEnergyCurves_ = data.internalEnergyCurves_;
428 enableThermalDensity_ = data.enableThermalDensity_;
429 enableJouleThomson_ = data.enableJouleThomson_;
430 enableThermalViscosity_ = data.enableThermalViscosity_;
431 enableInternalEnergy_ = data.enableInternalEnergy_;
432
433 return *this;
434 }
435
436private:
437 IsothermalPvt* isothermalPvt_;
438
439 // The PVT properties needed for temperature dependence of the viscosity. We need
440 // to store one value per PVT region.
441 std::vector<TabulatedOneDFunction> oilvisctCurves_;
442 std::vector<Scalar> viscrefPress_;
443 std::vector<Scalar> viscrefRs_;
444 std::vector<Scalar> viscRef_;
445
446 // The PVT properties needed for temperature dependence of the density.
447 std::vector<Scalar> oildentRefTemp_;
448 std::vector<Scalar> oildentCT1_;
449 std::vector<Scalar> oildentCT2_;
450
451 std::vector<Scalar> oilJTRefPres_;
452 std::vector<Scalar> oilJTC_;
453
454 std::vector<Scalar> rhoRefG_;
455
456 // piecewise linear curve representing the internal energy of oil
457 std::vector<TabulatedOneDFunction> internalEnergyCurves_;
458
459 bool enableThermalDensity_;
460 bool enableJouleThomson_;
461 bool enableThermalViscosity_;
462 bool enableInternalEnergy_;
463};
464
465} // namespace Opm
466
467#endif
Definition: EclipseState.hpp:55
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
Definition: OilPvtMultiplexer.hpp:97
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:202
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:193
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: OilPvtMultiplexer.hpp:212
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:221
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:183
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:245
const Scalar oilReferenceDensity(unsigned regionIdx) const
Return the reference density which are considered by this PVT-object.
Definition: OilPvtMultiplexer.hpp:166
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:254
This class implements temperature dependence of the PVT properties of oil.
Definition: OilPvtThermal.hpp:50
void initEnd()
Finish initializing the thermal part of the oil phase PVT properties.
Definition: OilPvtThermal.hpp:130
bool enableThermalDensity() const
Returns true iff the density of the oil phase is temperature dependent.
Definition: OilPvtThermal.hpp:136
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:309
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:256
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:238
Evaluation internalEnergy(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rs) const
Returns the specific internal energy [J/kg] of oil given a set of parameters.
Definition: OilPvtThermal.hpp:158
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:220
void setNumRegions(size_t numRegions)
Set the number of PVT-regions considered by this object.
Definition: OilPvtThermal.hpp:112
bool enableThermalViscosity() const
Returns true iff the viscosity of the oil phase is temperature dependent.
Definition: OilPvtThermal.hpp:148
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the saturation pressure of the oil phase [Pa].
Definition: OilPvtThermal.hpp:337
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:322
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of gas-saturated oil phase.
Definition: OilPvtThermal.hpp:281
bool enableJouleThomson() const
Returns true iff Joule-Thomson effect for the oil phase is active.
Definition: OilPvtThermal.hpp:142
Definition: Schedule.hpp:130
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:51
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:30