My Project
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
37namespace Opm {
38
39#if HAVE_ECL_INPUT
40class EclipseState;
41class Schedule;
42#endif
43
44#define OPM_GAS_PVT_MULTIPLEXER_CALL(codeToCall) \
45 switch (gasPvtApproach_) { \
46 case GasPvtApproach::DryGas: { \
47 auto& pvtImpl = getRealPvt<GasPvtApproach::DryGas>(); \
48 codeToCall; \
49 break; \
50 } \
51 case GasPvtApproach::DryHumidGas: { \
52 auto& pvtImpl = getRealPvt<GasPvtApproach::DryHumidGas>(); \
53 codeToCall; \
54 break; \
55 } \
56 case GasPvtApproach::WetHumidGas: { \
57 auto& pvtImpl = getRealPvt<GasPvtApproach::WetHumidGas>(); \
58 codeToCall; \
59 break; \
60 } \
61 case GasPvtApproach::WetGas: { \
62 auto& pvtImpl = getRealPvt<GasPvtApproach::WetGas>(); \
63 codeToCall; \
64 break; \
65 } \
66 case GasPvtApproach::ThermalGas: { \
67 auto& pvtImpl = getRealPvt<GasPvtApproach::ThermalGas>(); \
68 codeToCall; \
69 break; \
70 } \
71 case GasPvtApproach::Co2Gas: { \
72 auto& pvtImpl = getRealPvt<GasPvtApproach::Co2Gas>(); \
73 codeToCall; \
74 break; \
75 } \
76 case GasPvtApproach::NoGas: \
77 throw std::logic_error("Not implemented: Gas PVT of this deck!"); \
78 }
79
80enum class GasPvtApproach {
81 NoGas,
82 DryGas,
83 DryHumidGas,
84 WetHumidGas,
85 WetGas,
86 ThermalGas,
87 Co2Gas
88};
89
100template <class Scalar, bool enableThermal = true>
102{
103public:
105 {
106 gasPvtApproach_ = GasPvtApproach::NoGas;
107 realGasPvt_ = nullptr;
108 }
109
110 GasPvtMultiplexer(GasPvtApproach approach, void* realGasPvt)
111 : gasPvtApproach_(approach)
112 , realGasPvt_(realGasPvt)
113 { }
114
116 {
117 *this = data;
118 }
119
121 {
122 switch (gasPvtApproach_) {
123 case GasPvtApproach::DryGas: {
124 delete &getRealPvt<GasPvtApproach::DryGas>();
125 break;
126 }
127 case GasPvtApproach::DryHumidGas: {
128 delete &getRealPvt<GasPvtApproach::DryHumidGas>();
129 break;
130 }
131 case GasPvtApproach::WetHumidGas: {
132 delete &getRealPvt<GasPvtApproach::WetHumidGas>();
133 break;
134 }
135 case GasPvtApproach::WetGas: {
136 delete &getRealPvt<GasPvtApproach::WetGas>();
137 break;
138 }
139 case GasPvtApproach::ThermalGas: {
140 delete &getRealPvt<GasPvtApproach::ThermalGas>();
141 break;
142 }
143 case GasPvtApproach::Co2Gas: {
144 delete &getRealPvt<GasPvtApproach::Co2Gas>();
145 break;
146 }
147 case GasPvtApproach::NoGas:
148 break;
149 }
150 }
151
152#if HAVE_ECL_INPUT
158 void initFromState(const EclipseState& eclState, const Schedule& schedule);
159#endif // HAVE_ECL_INPUT
160
161 void setApproach(GasPvtApproach gasPvtAppr)
162 {
163 switch (gasPvtAppr) {
164 case GasPvtApproach::DryGas:
165 realGasPvt_ = new DryGasPvt<Scalar>;
166 break;
167
168 case GasPvtApproach::DryHumidGas:
169 realGasPvt_ = new DryHumidGasPvt<Scalar>;
170 break;
171
172 case GasPvtApproach::WetHumidGas:
173 realGasPvt_ = new WetHumidGasPvt<Scalar>;
174 break;
175
176 case GasPvtApproach::WetGas:
177 realGasPvt_ = new WetGasPvt<Scalar>;
178 break;
179
180 case GasPvtApproach::ThermalGas:
181 realGasPvt_ = new GasPvtThermal<Scalar>;
182 break;
183
184 case GasPvtApproach::Co2Gas:
185 realGasPvt_ = new Co2GasPvt<Scalar>;
186 break;
187
188 case GasPvtApproach::NoGas:
189 throw std::logic_error("Not implemented: Gas PVT of this deck!");
190 }
191
192 gasPvtApproach_ = gasPvtAppr;
193 }
194
195 void initEnd()
196 { OPM_GAS_PVT_MULTIPLEXER_CALL(pvtImpl.initEnd()); }
197
201 unsigned numRegions() const
202 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.numRegions()); return 1; }
203
207 const Scalar gasReferenceDensity(unsigned regionIdx)
208 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.gasReferenceDensity(regionIdx)); return 2.; }
209
213 template <class Evaluation>
214 Evaluation internalEnergy(unsigned regionIdx,
215 const Evaluation& temperature,
216 const Evaluation& pressure,
217 const Evaluation& Rv,
218 const Evaluation& Rvw) const
219 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.internalEnergy(regionIdx, temperature, pressure, Rv, Rvw)); return 0; }
220
224 template <class Evaluation = Scalar>
225 Evaluation viscosity(unsigned regionIdx,
226 const Evaluation& temperature,
227 const Evaluation& pressure,
228 const Evaluation& Rv,
229 const Evaluation& Rvw ) const
230 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.viscosity(regionIdx, temperature, pressure, Rv, Rvw)); return 0; }
231
235 template <class Evaluation = Scalar>
236 Evaluation saturatedViscosity(unsigned regionIdx,
237 const Evaluation& temperature,
238 const Evaluation& pressure) const
239 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedViscosity(regionIdx, temperature, pressure)); return 0; }
240
244 template <class Evaluation = Scalar>
245 Evaluation inverseFormationVolumeFactor(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.inverseFormationVolumeFactor(regionIdx, temperature, pressure, Rv, Rvw)); return 0; }
251
255 template <class Evaluation = Scalar>
256 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
257 const Evaluation& temperature,
258 const Evaluation& pressure) const
259 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedInverseFormationVolumeFactor(regionIdx, temperature, pressure)); return 0; }
260
264 template <class Evaluation = Scalar>
265 Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
266 const Evaluation& temperature,
267 const Evaluation& pressure) const
268 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedOilVaporizationFactor(regionIdx, temperature, pressure)); return 0; }
269
273 template <class Evaluation = Scalar>
274 Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
275 const Evaluation& temperature,
276 const Evaluation& pressure,
277 const Evaluation& oilSaturation,
278 const Evaluation& maxOilSaturation) const
279 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedOilVaporizationFactor(regionIdx, temperature, pressure, oilSaturation, maxOilSaturation)); return 0; }
280
284 template <class Evaluation = Scalar>
285 Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx,
286 const Evaluation& temperature,
287 const Evaluation& pressure) const
288 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedWaterVaporizationFactor(regionIdx, temperature, pressure)); return 0; }
289
293 template <class Evaluation = Scalar>
294 Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx,
295 const Evaluation& temperature,
296 const Evaluation& pressure,
297 const Evaluation& saltConcentration) const
298 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedWaterVaporizationFactor(regionIdx, temperature, pressure, saltConcentration)); return 0; }
299
306 template <class Evaluation = Scalar>
307 Evaluation saturationPressure(unsigned regionIdx,
308 const Evaluation& temperature,
309 const Evaluation& Rv) const
310 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturationPressure(regionIdx, temperature, Rv)); return 0; }
311
315 template <class Evaluation>
316 Evaluation diffusionCoefficient(const Evaluation& temperature,
317 const Evaluation& pressure,
318 unsigned compIdx) const
319 {
320 OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.diffusionCoefficient(temperature, pressure, compIdx)); return 0;
321 }
322
328 GasPvtApproach gasPvtApproach() const
329 { return gasPvtApproach_; }
330
331 // get the parameter object for the dry gas case
332 template <GasPvtApproach approachV>
333 typename std::enable_if<approachV == GasPvtApproach::DryGas, DryGasPvt<Scalar> >::type& getRealPvt()
334 {
335 assert(gasPvtApproach() == approachV);
336 return *static_cast<DryGasPvt<Scalar>* >(realGasPvt_);
337 }
338
339 template <GasPvtApproach approachV>
340 typename std::enable_if<approachV == GasPvtApproach::DryGas, const DryGasPvt<Scalar> >::type& getRealPvt() const
341 {
342 assert(gasPvtApproach() == approachV);
343 return *static_cast<const DryGasPvt<Scalar>* >(realGasPvt_);
344 }
345
346 // get the parameter object for the dry humid gas case
347 template <GasPvtApproach approachV>
348 typename std::enable_if<approachV == GasPvtApproach::DryHumidGas, DryHumidGasPvt<Scalar> >::type& getRealPvt()
349 {
350 assert(gasPvtApproach() == approachV);
351 return *static_cast<DryHumidGasPvt<Scalar>* >(realGasPvt_);
352 }
353
354 template <GasPvtApproach approachV>
355 typename std::enable_if<approachV == GasPvtApproach::DryHumidGas, const DryHumidGasPvt<Scalar> >::type& getRealPvt() const
356 {
357 assert(gasPvtApproach() == approachV);
358 return *static_cast<const DryHumidGasPvt<Scalar>* >(realGasPvt_);
359 }
360
361 // get the parameter object for the wet humid gas case
362 template <GasPvtApproach approachV>
363 typename std::enable_if<approachV == GasPvtApproach::WetHumidGas, WetHumidGasPvt<Scalar> >::type& getRealPvt()
364 {
365 assert(gasPvtApproach() == approachV);
366 return *static_cast<WetHumidGasPvt<Scalar>* >(realGasPvt_);
367 }
368
369 template <GasPvtApproach approachV>
370 typename std::enable_if<approachV == GasPvtApproach::WetHumidGas, const WetHumidGasPvt<Scalar> >::type& getRealPvt() const
371 {
372 assert(gasPvtApproach() == approachV);
373 return *static_cast<const WetHumidGasPvt<Scalar>* >(realGasPvt_);
374 }
375
376 // get the parameter object for the wet gas case
377 template <GasPvtApproach approachV>
378 typename std::enable_if<approachV == GasPvtApproach::WetGas, WetGasPvt<Scalar> >::type& getRealPvt()
379 {
380 assert(gasPvtApproach() == approachV);
381 return *static_cast<WetGasPvt<Scalar>* >(realGasPvt_);
382 }
383
384 template <GasPvtApproach approachV>
385 typename std::enable_if<approachV == GasPvtApproach::WetGas, const WetGasPvt<Scalar> >::type& getRealPvt() const
386 {
387 assert(gasPvtApproach() == approachV);
388 return *static_cast<const WetGasPvt<Scalar>* >(realGasPvt_);
389 }
390
391 // get the parameter object for the thermal gas case
392 template <GasPvtApproach approachV>
393 typename std::enable_if<approachV == GasPvtApproach::ThermalGas, GasPvtThermal<Scalar> >::type& getRealPvt()
394 {
395 assert(gasPvtApproach() == approachV);
396 return *static_cast<GasPvtThermal<Scalar>* >(realGasPvt_);
397 }
398 template <GasPvtApproach approachV>
399 typename std::enable_if<approachV == GasPvtApproach::ThermalGas, const GasPvtThermal<Scalar> >::type& getRealPvt() const
400 {
401 assert(gasPvtApproach() == approachV);
402 return *static_cast<const GasPvtThermal<Scalar>* >(realGasPvt_);
403 }
404
405 template <GasPvtApproach approachV>
406 typename std::enable_if<approachV == GasPvtApproach::Co2Gas, Co2GasPvt<Scalar> >::type& getRealPvt()
407 {
408 assert(gasPvtApproach() == approachV);
409 return *static_cast<Co2GasPvt<Scalar>* >(realGasPvt_);
410 }
411
412 template <GasPvtApproach approachV>
413 typename std::enable_if<approachV == GasPvtApproach::Co2Gas, const Co2GasPvt<Scalar> >::type& getRealPvt() const
414 {
415 assert(gasPvtApproach() == approachV);
416 return *static_cast<const Co2GasPvt<Scalar>* >(realGasPvt_);
417 }
418
419 const void* realGasPvt() const { return realGasPvt_; }
420
421 GasPvtMultiplexer<Scalar,enableThermal>& operator=(const GasPvtMultiplexer<Scalar,enableThermal>& data)
422 {
423 gasPvtApproach_ = data.gasPvtApproach_;
424 switch (gasPvtApproach_) {
425 case GasPvtApproach::DryGas:
426 realGasPvt_ = new DryGasPvt<Scalar>(*static_cast<const DryGasPvt<Scalar>*>(data.realGasPvt_));
427 break;
428 case GasPvtApproach::DryHumidGas:
429 realGasPvt_ = new DryHumidGasPvt<Scalar>(*static_cast<const DryHumidGasPvt<Scalar>*>(data.realGasPvt_));
430 break;
431 case GasPvtApproach::WetHumidGas:
432 realGasPvt_ = new WetHumidGasPvt<Scalar>(*static_cast<const WetHumidGasPvt<Scalar>*>(data.realGasPvt_));
433 break;
434 case GasPvtApproach::WetGas:
435 realGasPvt_ = new WetGasPvt<Scalar>(*static_cast<const WetGasPvt<Scalar>*>(data.realGasPvt_));
436 break;
437 case GasPvtApproach::ThermalGas:
438 realGasPvt_ = new GasPvtThermal<Scalar>(*static_cast<const GasPvtThermal<Scalar>*>(data.realGasPvt_));
439 break;
440 case GasPvtApproach::Co2Gas:
441 realGasPvt_ = new Co2GasPvt<Scalar>(*static_cast<const Co2GasPvt<Scalar>*>(data.realGasPvt_));
442 break;
443 default:
444 break;
445 }
446
447 return *this;
448 }
449
450private:
451 GasPvtApproach gasPvtApproach_;
452 void* realGasPvt_;
453};
454
455} // namespace Opm
456
457#endif
This class represents the Pressure-Volume-Temperature relations of the gas phase for CO2.
Definition: Co2GasPvt.hpp:53
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:55
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition: GasPvtMultiplexer.hpp:102
GasPvtApproach gasPvtApproach() const
Returns the concrete approach for calculating the PVT relations.
Definition: GasPvtMultiplexer.hpp:328
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:256
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: GasPvtMultiplexer.hpp:201
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:285
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:245
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:214
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:225
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:316
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:294
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:265
const Scalar gasReferenceDensity(unsigned regionIdx)
Return the reference density which are considered by this PVT-object.
Definition: GasPvtMultiplexer.hpp:207
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:307
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:274
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:236
This class implements temperature dependence of the PVT properties of gas.
Definition: GasPvtThermal.hpp:50
Definition: Schedule.hpp:130
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