My Project
Loading...
Searching...
No Matches
GasPvtMultiplexer.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_GAS_PVT_MULTIPLEXER_HPP
28#define OPM_GAS_PVT_MULTIPLEXER_HPP
29
30#include "DryGasPvt.hpp"
31#include "DryHumidGasPvt.hpp"
32#include "WetHumidGasPvt.hpp"
33#include "WetGasPvt.hpp"
34#include "GasPvtThermal.hpp"
35#include "Co2GasPvt.hpp"
36#include "H2GasPvt.hpp"
37
38namespace Opm {
39
40#if HAVE_ECL_INPUT
41class EclipseState;
42class Schedule;
43#endif
44
45#define OPM_GAS_PVT_MULTIPLEXER_CALL(codeToCall) \
46 switch (gasPvtApproach_) { \
47 case GasPvtApproach::DryGas: { \
48 auto& pvtImpl = getRealPvt<GasPvtApproach::DryGas>(); \
49 codeToCall; \
50 break; \
51 } \
52 case GasPvtApproach::DryHumidGas: { \
53 auto& pvtImpl = getRealPvt<GasPvtApproach::DryHumidGas>(); \
54 codeToCall; \
55 break; \
56 } \
57 case GasPvtApproach::WetHumidGas: { \
58 auto& pvtImpl = getRealPvt<GasPvtApproach::WetHumidGas>(); \
59 codeToCall; \
60 break; \
61 } \
62 case GasPvtApproach::WetGas: { \
63 auto& pvtImpl = getRealPvt<GasPvtApproach::WetGas>(); \
64 codeToCall; \
65 break; \
66 } \
67 case GasPvtApproach::ThermalGas: { \
68 auto& pvtImpl = getRealPvt<GasPvtApproach::ThermalGas>(); \
69 codeToCall; \
70 break; \
71 } \
72 case GasPvtApproach::Co2Gas: { \
73 auto& pvtImpl = getRealPvt<GasPvtApproach::Co2Gas>(); \
74 codeToCall; \
75 break; \
76 } \
77 case GasPvtApproach::H2Gas: { \
78 auto& pvtImpl = getRealPvt<GasPvtApproach::H2Gas>(); \
79 codeToCall; \
80 break; \
81 } \
82 case GasPvtApproach::NoGas: \
83 throw std::logic_error("Not implemented: Gas PVT of this deck!"); \
84 }
85
86enum class GasPvtApproach {
87 NoGas,
88 DryGas,
89 DryHumidGas,
90 WetHumidGas,
91 WetGas,
92 ThermalGas,
93 Co2Gas,
94 H2Gas
95};
96
107template <class Scalar, bool enableThermal = true>
109{
110public:
112 {
113 gasPvtApproach_ = GasPvtApproach::NoGas;
114 realGasPvt_ = nullptr;
115 }
116
117 GasPvtMultiplexer(GasPvtApproach approach, void* realGasPvt)
118 : gasPvtApproach_(approach)
119 , realGasPvt_(realGasPvt)
120 { }
121
123 {
124 *this = data;
125 }
126
128 {
129 switch (gasPvtApproach_) {
130 case GasPvtApproach::DryGas: {
131 delete &getRealPvt<GasPvtApproach::DryGas>();
132 break;
133 }
134 case GasPvtApproach::DryHumidGas: {
135 delete &getRealPvt<GasPvtApproach::DryHumidGas>();
136 break;
137 }
138 case GasPvtApproach::WetHumidGas: {
139 delete &getRealPvt<GasPvtApproach::WetHumidGas>();
140 break;
141 }
142 case GasPvtApproach::WetGas: {
143 delete &getRealPvt<GasPvtApproach::WetGas>();
144 break;
145 }
146 case GasPvtApproach::ThermalGas: {
147 delete &getRealPvt<GasPvtApproach::ThermalGas>();
148 break;
149 }
150 case GasPvtApproach::Co2Gas: {
151 delete &getRealPvt<GasPvtApproach::Co2Gas>();
152 break;
153 }
154 case GasPvtApproach::H2Gas: {
155 delete &getRealPvt<GasPvtApproach::H2Gas>();
156 break;
157 }
158 case GasPvtApproach::NoGas:
159 break;
160 }
161 }
162
163#if HAVE_ECL_INPUT
169 void initFromState(const EclipseState& eclState, const Schedule& schedule);
170#endif // HAVE_ECL_INPUT
171
172 void setApproach(GasPvtApproach gasPvtAppr)
173 {
174 switch (gasPvtAppr) {
175 case GasPvtApproach::DryGas:
176 realGasPvt_ = new DryGasPvt<Scalar>;
177 break;
178
179 case GasPvtApproach::DryHumidGas:
180 realGasPvt_ = new DryHumidGasPvt<Scalar>;
181 break;
182
183 case GasPvtApproach::WetHumidGas:
184 realGasPvt_ = new WetHumidGasPvt<Scalar>;
185 break;
186
187 case GasPvtApproach::WetGas:
188 realGasPvt_ = new WetGasPvt<Scalar>;
189 break;
190
191 case GasPvtApproach::ThermalGas:
192 realGasPvt_ = new GasPvtThermal<Scalar>;
193 break;
194
195 case GasPvtApproach::Co2Gas:
196 realGasPvt_ = new Co2GasPvt<Scalar>;
197 break;
198
199 case GasPvtApproach::H2Gas:
200 realGasPvt_ = new H2GasPvt<Scalar>;
201 break;
202
203 case GasPvtApproach::NoGas:
204 throw std::logic_error("Not implemented: Gas PVT of this deck!");
205 }
206
207 gasPvtApproach_ = gasPvtAppr;
208 }
209
210 void initEnd()
211 { OPM_GAS_PVT_MULTIPLEXER_CALL(pvtImpl.initEnd()); }
212
216 unsigned numRegions() const
217 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.numRegions()); return 1; }
218
219 void setVapPars(const Scalar par1, const Scalar par2)
220 {
221 OPM_GAS_PVT_MULTIPLEXER_CALL(pvtImpl.setVapPars(par1, par2));
222 }
223
227 const Scalar gasReferenceDensity(unsigned regionIdx)
228 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.gasReferenceDensity(regionIdx)); return 2.; }
229
233 template <class Evaluation>
234 Evaluation internalEnergy(unsigned regionIdx,
235 const Evaluation& temperature,
236 const Evaluation& pressure,
237 const Evaluation& Rv,
238 const Evaluation& Rvw) const
239 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.internalEnergy(regionIdx, temperature, pressure, Rv, Rvw)); return 0; }
240
244 template <class Evaluation = Scalar>
245 Evaluation viscosity(unsigned regionIdx,
246 const Evaluation& temperature,
247 const Evaluation& pressure,
248 const Evaluation& Rv,
249 const Evaluation& Rvw ) const
250 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.viscosity(regionIdx, temperature, pressure, Rv, Rvw)); return 0; }
251
255 template <class Evaluation = Scalar>
256 Evaluation saturatedViscosity(unsigned regionIdx,
257 const Evaluation& temperature,
258 const Evaluation& pressure) const
259 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedViscosity(regionIdx, temperature, pressure)); return 0; }
260
264 template <class Evaluation = Scalar>
265 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
266 const Evaluation& temperature,
267 const Evaluation& pressure,
268 const Evaluation& Rv,
269 const Evaluation& Rvw) const
270 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.inverseFormationVolumeFactor(regionIdx, temperature, pressure, Rv, Rvw)); return 0; }
271
275 template <class Evaluation = Scalar>
276 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
277 const Evaluation& temperature,
278 const Evaluation& pressure) const
279 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedInverseFormationVolumeFactor(regionIdx, temperature, pressure)); return 0; }
280
284 template <class Evaluation = Scalar>
285 Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
286 const Evaluation& temperature,
287 const Evaluation& pressure) const
288 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedOilVaporizationFactor(regionIdx, temperature, pressure)); return 0; }
289
293 template <class Evaluation = Scalar>
294 Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
295 const Evaluation& temperature,
296 const Evaluation& pressure,
297 const Evaluation& oilSaturation,
298 const Evaluation& maxOilSaturation) const
299 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedOilVaporizationFactor(regionIdx, temperature, pressure, oilSaturation, maxOilSaturation)); return 0; }
300
304 template <class Evaluation = Scalar>
305 Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx,
306 const Evaluation& temperature,
307 const Evaluation& pressure) const
308 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedWaterVaporizationFactor(regionIdx, temperature, pressure)); return 0; }
309
313 template <class Evaluation = Scalar>
314 Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx,
315 const Evaluation& temperature,
316 const Evaluation& pressure,
317 const Evaluation& saltConcentration) const
318 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedWaterVaporizationFactor(regionIdx, temperature, pressure, saltConcentration)); return 0; }
319
326 template <class Evaluation = Scalar>
327 Evaluation saturationPressure(unsigned regionIdx,
328 const Evaluation& temperature,
329 const Evaluation& Rv) const
330 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturationPressure(regionIdx, temperature, Rv)); return 0; }
331
335 template <class Evaluation>
336 Evaluation diffusionCoefficient(const Evaluation& temperature,
337 const Evaluation& pressure,
338 unsigned compIdx) const
339 {
340 OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.diffusionCoefficient(temperature, pressure, compIdx)); return 0;
341 }
342
348 GasPvtApproach gasPvtApproach() const
349 { return gasPvtApproach_; }
350
351 // get the parameter object for the dry gas case
352 template <GasPvtApproach approachV>
353 typename std::enable_if<approachV == GasPvtApproach::DryGas, DryGasPvt<Scalar> >::type& getRealPvt()
354 {
355 assert(gasPvtApproach() == approachV);
356 return *static_cast<DryGasPvt<Scalar>* >(realGasPvt_);
357 }
358
359 template <GasPvtApproach approachV>
360 typename std::enable_if<approachV == GasPvtApproach::DryGas, const DryGasPvt<Scalar> >::type& getRealPvt() const
361 {
362 assert(gasPvtApproach() == approachV);
363 return *static_cast<const DryGasPvt<Scalar>* >(realGasPvt_);
364 }
365
366 // get the parameter object for the dry humid gas case
367 template <GasPvtApproach approachV>
368 typename std::enable_if<approachV == GasPvtApproach::DryHumidGas, DryHumidGasPvt<Scalar> >::type& getRealPvt()
369 {
370 assert(gasPvtApproach() == approachV);
371 return *static_cast<DryHumidGasPvt<Scalar>* >(realGasPvt_);
372 }
373
374 template <GasPvtApproach approachV>
375 typename std::enable_if<approachV == GasPvtApproach::DryHumidGas, const DryHumidGasPvt<Scalar> >::type& getRealPvt() const
376 {
377 assert(gasPvtApproach() == approachV);
378 return *static_cast<const DryHumidGasPvt<Scalar>* >(realGasPvt_);
379 }
380
381 // get the parameter object for the wet humid gas case
382 template <GasPvtApproach approachV>
383 typename std::enable_if<approachV == GasPvtApproach::WetHumidGas, WetHumidGasPvt<Scalar> >::type& getRealPvt()
384 {
385 assert(gasPvtApproach() == approachV);
386 return *static_cast<WetHumidGasPvt<Scalar>* >(realGasPvt_);
387 }
388
389 template <GasPvtApproach approachV>
390 typename std::enable_if<approachV == GasPvtApproach::WetHumidGas, const WetHumidGasPvt<Scalar> >::type& getRealPvt() const
391 {
392 assert(gasPvtApproach() == approachV);
393 return *static_cast<const WetHumidGasPvt<Scalar>* >(realGasPvt_);
394 }
395
396 // get the parameter object for the wet gas case
397 template <GasPvtApproach approachV>
398 typename std::enable_if<approachV == GasPvtApproach::WetGas, WetGasPvt<Scalar> >::type& getRealPvt()
399 {
400 assert(gasPvtApproach() == approachV);
401 return *static_cast<WetGasPvt<Scalar>* >(realGasPvt_);
402 }
403
404 template <GasPvtApproach approachV>
405 typename std::enable_if<approachV == GasPvtApproach::WetGas, const WetGasPvt<Scalar> >::type& getRealPvt() const
406 {
407 assert(gasPvtApproach() == approachV);
408 return *static_cast<const WetGasPvt<Scalar>* >(realGasPvt_);
409 }
410
411 // get the parameter object for the thermal gas case
412 template <GasPvtApproach approachV>
413 typename std::enable_if<approachV == GasPvtApproach::ThermalGas, GasPvtThermal<Scalar> >::type& getRealPvt()
414 {
415 assert(gasPvtApproach() == approachV);
416 return *static_cast<GasPvtThermal<Scalar>* >(realGasPvt_);
417 }
418 template <GasPvtApproach approachV>
419 typename std::enable_if<approachV == GasPvtApproach::ThermalGas, const GasPvtThermal<Scalar> >::type& getRealPvt() const
420 {
421 assert(gasPvtApproach() == approachV);
422 return *static_cast<const GasPvtThermal<Scalar>* >(realGasPvt_);
423 }
424
425 template <GasPvtApproach approachV>
426 typename std::enable_if<approachV == GasPvtApproach::Co2Gas, Co2GasPvt<Scalar> >::type& getRealPvt()
427 {
428 assert(gasPvtApproach() == approachV);
429 return *static_cast<Co2GasPvt<Scalar>* >(realGasPvt_);
430 }
431
432 template <GasPvtApproach approachV>
433 typename std::enable_if<approachV == GasPvtApproach::Co2Gas, const Co2GasPvt<Scalar> >::type& getRealPvt() const
434 {
435 assert(gasPvtApproach() == approachV);
436 return *static_cast<const Co2GasPvt<Scalar>* >(realGasPvt_);
437 }
438
439 template <GasPvtApproach approachV>
440 typename std::enable_if<approachV == GasPvtApproach::H2Gas, H2GasPvt<Scalar> >::type& getRealPvt()
441 {
442 assert(gasPvtApproach() == approachV);
443 return *static_cast<H2GasPvt<Scalar>* >(realGasPvt_);
444 }
445
446 template <GasPvtApproach approachV>
447 typename std::enable_if<approachV == GasPvtApproach::H2Gas, const H2GasPvt<Scalar> >::type& getRealPvt() const
448 {
449 assert(gasPvtApproach() == approachV);
450 return *static_cast<const H2GasPvt<Scalar>* >(realGasPvt_);
451 }
452
453 const void* realGasPvt() const { return realGasPvt_; }
454
455 GasPvtMultiplexer<Scalar,enableThermal>& operator=(const GasPvtMultiplexer<Scalar,enableThermal>& data)
456 {
457 gasPvtApproach_ = data.gasPvtApproach_;
458 switch (gasPvtApproach_) {
459 case GasPvtApproach::DryGas:
460 realGasPvt_ = new DryGasPvt<Scalar>(*static_cast<const DryGasPvt<Scalar>*>(data.realGasPvt_));
461 break;
462 case GasPvtApproach::DryHumidGas:
463 realGasPvt_ = new DryHumidGasPvt<Scalar>(*static_cast<const DryHumidGasPvt<Scalar>*>(data.realGasPvt_));
464 break;
465 case GasPvtApproach::WetHumidGas:
466 realGasPvt_ = new WetHumidGasPvt<Scalar>(*static_cast<const WetHumidGasPvt<Scalar>*>(data.realGasPvt_));
467 break;
468 case GasPvtApproach::WetGas:
469 realGasPvt_ = new WetGasPvt<Scalar>(*static_cast<const WetGasPvt<Scalar>*>(data.realGasPvt_));
470 break;
471 case GasPvtApproach::ThermalGas:
472 realGasPvt_ = new GasPvtThermal<Scalar>(*static_cast<const GasPvtThermal<Scalar>*>(data.realGasPvt_));
473 break;
474 case GasPvtApproach::Co2Gas:
475 realGasPvt_ = new Co2GasPvt<Scalar>(*static_cast<const Co2GasPvt<Scalar>*>(data.realGasPvt_));
476 break;
477 case GasPvtApproach::H2Gas:
478 realGasPvt_ = new H2GasPvt<Scalar>(*static_cast<const H2GasPvt<Scalar>*>(data.realGasPvt_));
479 break;
480 default:
481 break;
482 }
483
484 return *this;
485 }
486
487private:
488 GasPvtApproach gasPvtApproach_;
489 void* realGasPvt_;
490};
491
492} // namespace Opm
493
494#endif
This class represents the Pressure-Volume-Temperature relations of the gas phase for CO2.
This class represents the Pressure-Volume-Temperature relations of the gas phase without vaporized oi...
This class represents the Pressure-Volume-Temperature relations of the gas phase with vaporized water...
This class implements temperature dependence of the PVT properties of gas.
This class represents the Pressure-Volume-Temperature relations of the gas phase for H2.
This class represents the Pressure-Volume-Temperature relations of the gas phas with vaporized oil.
This class represents the Pressure-Volume-Temperature relations of the gas phase with vaporized oil a...
This class represents the Pressure-Volume-Temperature relations of the gas phase for CO2.
Definition Co2GasPvt.hpp:54
This class represents the Pressure-Volume-Temperature relations of the gas phase without vaporized oi...
Definition DryGasPvt.hpp:49
This class represents the Pressure-Volume-Temperature relations of the gas phase with vaporized water...
Definition DryHumidGasPvt.hpp:52
Definition EclipseState.hpp:63
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition GasPvtMultiplexer.hpp:109
GasPvtApproach gasPvtApproach() const
Returns the concrete approach for calculating the PVT relations.
Definition GasPvtMultiplexer.hpp:348
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of oil saturated gas given a set of parameters.
Definition GasPvtMultiplexer.hpp:276
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition GasPvtMultiplexer.hpp:216
Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the water vaporization factor [m^3/m^3] of water saturated gas.
Definition GasPvtMultiplexer.hpp:305
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rv, const Evaluation &Rvw) const
Returns the formation volume factor [-] of the fluid phase.
Definition GasPvtMultiplexer.hpp:265
Evaluation internalEnergy(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rv, const Evaluation &Rvw) const
Returns the specific enthalpy [J/kg] of gas given a set of parameters.
Definition GasPvtMultiplexer.hpp:234
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rv, const Evaluation &Rvw) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition GasPvtMultiplexer.hpp:245
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 GasPvtMultiplexer.hpp:336
Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &saltConcentration) const
Returns the water vaporization factor [m^3/m^3] of water saturated gas.
Definition GasPvtMultiplexer.hpp:314
Evaluation saturatedOilVaporizationFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the oil vaporization factor [m^3/m^3] of oil saturated gas.
Definition GasPvtMultiplexer.hpp:285
const Scalar gasReferenceDensity(unsigned regionIdx)
Return the reference density which are considered by this PVT-object.
Definition GasPvtMultiplexer.hpp:227
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &temperature, const Evaluation &Rv) const
Returns the saturation pressure of the gas phase [Pa] depending on its mass fraction of the oil compo...
Definition GasPvtMultiplexer.hpp:327
Evaluation saturatedOilVaporizationFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &oilSaturation, const Evaluation &maxOilSaturation) const
Returns the oil vaporization factor [m^3/m^3] of oil saturated gas.
Definition GasPvtMultiplexer.hpp:294
Evaluation saturatedViscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of oil saturated gas given a set of parameters.
Definition GasPvtMultiplexer.hpp:256
This class implements temperature dependence of the PVT properties of gas.
Definition GasPvtThermal.hpp:50
This class represents the Pressure-Volume-Temperature relations of the gas phase for H2.
Definition H2GasPvt.hpp:50
Definition Schedule.hpp:88
This class represents the Pressure-Volume-Temperature relations of the gas phas with vaporized oil.
Definition WetGasPvt.hpp:51
This class represents the Pressure-Volume-Temperature relations of the gas phase with vaporized oil a...
Definition WetHumidGasPvt.hpp:51
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30