My Project
Loading...
Searching...
No Matches
OilPvtMultiplexer.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_MULTIPLEXER_HPP
28#define OPM_OIL_PVT_MULTIPLEXER_HPP
29
31#include "DeadOilPvt.hpp"
32#include "LiveOilPvt.hpp"
33#include "OilPvtThermal.hpp"
34#include "BrineCo2Pvt.hpp"
35#include "BrineH2Pvt.hpp"
36
37namespace Opm {
38
39#if HAVE_ECL_INPUT
40class EclipseState;
41class Schedule;
42#endif
43
44#define OPM_OIL_PVT_MULTIPLEXER_CALL(codeToCall) \
45 switch (approach_) { \
46 case OilPvtApproach::ConstantCompressibilityOil: { \
47 auto& pvtImpl = getRealPvt<OilPvtApproach::ConstantCompressibilityOil>(); \
48 codeToCall; \
49 break; \
50 } \
51 case OilPvtApproach::DeadOil: { \
52 auto& pvtImpl = getRealPvt<OilPvtApproach::DeadOil>(); \
53 codeToCall; \
54 break; \
55 } \
56 case OilPvtApproach::LiveOil: { \
57 auto& pvtImpl = getRealPvt<OilPvtApproach::LiveOil>(); \
58 codeToCall; \
59 break; \
60 } \
61 case OilPvtApproach::ThermalOil: { \
62 auto& pvtImpl = getRealPvt<OilPvtApproach::ThermalOil>(); \
63 codeToCall; \
64 break; \
65 } \
66 case OilPvtApproach::BrineCo2: { \
67 auto& pvtImpl = getRealPvt<OilPvtApproach::BrineCo2>(); \
68 codeToCall; \
69 break; \
70 } \
71 case OilPvtApproach::BrineH2: { \
72 auto& pvtImpl = getRealPvt<OilPvtApproach::BrineH2>(); \
73 codeToCall; \
74 break; \
75 } \
76 case OilPvtApproach::NoOil: \
77 throw std::logic_error("Not implemented: Oil PVT of this deck!"); \
78 }
79
80enum class OilPvtApproach {
81 NoOil,
82 LiveOil,
83 DeadOil,
84 ConstantCompressibilityOil,
85 ThermalOil,
86 BrineCo2,
87 BrineH2
88};
89
102template <class Scalar, bool enableThermal = true>
104{
105public:
107 {
108 approach_ = OilPvtApproach::NoOil;
109 realOilPvt_ = nullptr;
110 }
111
112 OilPvtMultiplexer(OilPvtApproach approach, void* realOilPvt)
113 : approach_(approach)
114 , realOilPvt_(realOilPvt)
115 { }
116
118 {
119 *this = data;
120 }
121
123 {
124 switch (approach_) {
125 case OilPvtApproach::LiveOil: {
126 delete &getRealPvt<OilPvtApproach::LiveOil>();
127 break;
128 }
129 case OilPvtApproach::DeadOil: {
130 delete &getRealPvt<OilPvtApproach::DeadOil>();
131 break;
132 }
133 case OilPvtApproach::ConstantCompressibilityOil: {
134 delete &getRealPvt<OilPvtApproach::ConstantCompressibilityOil>();
135 break;
136 }
137 case OilPvtApproach::ThermalOil: {
138 delete &getRealPvt<OilPvtApproach::ThermalOil>();
139 break;
140 }
141 case OilPvtApproach::BrineCo2: {
142 delete &getRealPvt<OilPvtApproach::BrineCo2>();
143 break;
144 }
145 case OilPvtApproach::BrineH2: {
146 delete &getRealPvt<OilPvtApproach::BrineH2>();
147 break;
148 }
149 case OilPvtApproach::NoOil:
150 break;
151 }
152 }
153
154#if HAVE_ECL_INPUT
160 void initFromState(const EclipseState& eclState, const Schedule& schedule);
161#endif // HAVE_ECL_INPUT
162
163
164 void initEnd()
165 { OPM_OIL_PVT_MULTIPLEXER_CALL(pvtImpl.initEnd()); }
166
170 unsigned numRegions() const
171 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.numRegions()); return 1; }
172
173 void setVapPars(const Scalar par1, const Scalar par2)
174 {
175 OPM_OIL_PVT_MULTIPLEXER_CALL(pvtImpl.setVapPars(par1, par2));
176 }
177
181 const Scalar oilReferenceDensity(unsigned regionIdx) const
182 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.oilReferenceDensity(regionIdx)); return 700.; }
183
187 template <class Evaluation>
188 Evaluation internalEnergy(unsigned regionIdx,
189 const Evaluation& temperature,
190 const Evaluation& pressure,
191 const Evaluation& Rs) const
192 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.internalEnergy(regionIdx, temperature, pressure, Rs)); return 0; }
193
197 template <class Evaluation>
198 Evaluation viscosity(unsigned regionIdx,
199 const Evaluation& temperature,
200 const Evaluation& pressure,
201 const Evaluation& Rs) const
202 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.viscosity(regionIdx, temperature, pressure, Rs)); return 0; }
203
207 template <class Evaluation>
208 Evaluation saturatedViscosity(unsigned regionIdx,
209 const Evaluation& temperature,
210 const Evaluation& pressure) const
211 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedViscosity(regionIdx, temperature, pressure)); return 0; }
212
216 template <class Evaluation>
217 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
218 const Evaluation& temperature,
219 const Evaluation& pressure,
220 const Evaluation& Rs) const
221 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.inverseFormationVolumeFactor(regionIdx, temperature, pressure, Rs)); return 0; }
222
226 template <class Evaluation>
227 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
228 const Evaluation& temperature,
229 const Evaluation& pressure) const
230 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedInverseFormationVolumeFactor(regionIdx, temperature, pressure)); return 0; }
231
235 template <class Evaluation>
236 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
237 const Evaluation& temperature,
238 const Evaluation& pressure) const
239 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedGasDissolutionFactor(regionIdx, temperature, pressure)); return 0; }
240
244 template <class Evaluation>
245 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
246 const Evaluation& temperature,
247 const Evaluation& pressure,
248 const Evaluation& oilSaturation,
249 const Evaluation& maxOilSaturation) const
250 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedGasDissolutionFactor(regionIdx, temperature, pressure, oilSaturation, maxOilSaturation)); return 0; }
251
259 template <class Evaluation>
260 Evaluation saturationPressure(unsigned regionIdx,
261 const Evaluation& temperature,
262 const Evaluation& Rs) const
263 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturationPressure(regionIdx, temperature, Rs)); return 0; }
264
268 template <class Evaluation>
269 Evaluation diffusionCoefficient(const Evaluation& temperature,
270 const Evaluation& pressure,
271 unsigned compIdx) const
272 {
273 OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.diffusionCoefficient(temperature, pressure, compIdx)); return 0;
274 }
275
276 void setApproach(OilPvtApproach appr)
277 {
278 switch (appr) {
279 case OilPvtApproach::LiveOil:
280 realOilPvt_ = new LiveOilPvt<Scalar>;
281 break;
282
283 case OilPvtApproach::DeadOil:
284 realOilPvt_ = new DeadOilPvt<Scalar>;
285 break;
286
287 case OilPvtApproach::ConstantCompressibilityOil:
289 break;
290
291 case OilPvtApproach::ThermalOil:
292 realOilPvt_ = new OilPvtThermal<Scalar>;
293 break;
294
295 case OilPvtApproach::BrineCo2:
296 realOilPvt_ = new BrineCo2Pvt<Scalar>;
297 break;
298
299 case OilPvtApproach::BrineH2:
300 realOilPvt_ = new BrineH2Pvt<Scalar>;
301 break;
302
303 case OilPvtApproach::NoOil:
304 throw std::logic_error("Not implemented: Oil PVT of this deck!");
305 }
306
307 approach_ = appr;
308 }
309
315 OilPvtApproach approach() const
316 { return approach_; }
317
318 // get the concrete parameter object for the oil phase
319 template <OilPvtApproach approachV>
320 typename std::enable_if<approachV == OilPvtApproach::LiveOil, LiveOilPvt<Scalar> >::type& getRealPvt()
321 {
322 assert(approach() == approachV);
323 return *static_cast<LiveOilPvt<Scalar>* >(realOilPvt_);
324 }
325
326 template <OilPvtApproach approachV>
327 typename std::enable_if<approachV == OilPvtApproach::LiveOil, const LiveOilPvt<Scalar> >::type& getRealPvt() const
328 {
329 assert(approach() == approachV);
330 return *static_cast<LiveOilPvt<Scalar>* >(realOilPvt_);
331 }
332
333 template <OilPvtApproach approachV>
334 typename std::enable_if<approachV == OilPvtApproach::DeadOil, DeadOilPvt<Scalar> >::type& getRealPvt()
335 {
336 assert(approach() == approachV);
337 return *static_cast<DeadOilPvt<Scalar>* >(realOilPvt_);
338 }
339
340 template <OilPvtApproach approachV>
341 typename std::enable_if<approachV == OilPvtApproach::DeadOil, const DeadOilPvt<Scalar> >::type& getRealPvt() const
342 {
343 assert(approach() == approachV);
344 return *static_cast<DeadOilPvt<Scalar>* >(realOilPvt_);
345 }
346
347 template <OilPvtApproach approachV>
348 typename std::enable_if<approachV == OilPvtApproach::ConstantCompressibilityOil, ConstantCompressibilityOilPvt<Scalar> >::type& getRealPvt()
349 {
350 assert(approach() == approachV);
351 return *static_cast<ConstantCompressibilityOilPvt<Scalar>* >(realOilPvt_);
352 }
353
354 template <OilPvtApproach approachV>
355 typename std::enable_if<approachV == OilPvtApproach::ConstantCompressibilityOil, const ConstantCompressibilityOilPvt<Scalar> >::type& getRealPvt() const
356 {
357 assert(approach() == approachV);
358 return *static_cast<ConstantCompressibilityOilPvt<Scalar>* >(realOilPvt_);
359 }
360
361 template <OilPvtApproach approachV>
362 typename std::enable_if<approachV == OilPvtApproach::ThermalOil, OilPvtThermal<Scalar> >::type& getRealPvt()
363 {
364 assert(approach() == approachV);
365 return *static_cast<OilPvtThermal<Scalar>* >(realOilPvt_);
366 }
367
368 template <OilPvtApproach approachV>
369 typename std::enable_if<approachV == OilPvtApproach::ThermalOil, const OilPvtThermal<Scalar> >::type& getRealPvt() const
370 {
371 assert(approach() == approachV);
372 return *static_cast<const OilPvtThermal<Scalar>* >(realOilPvt_);
373 }
374
375 template <OilPvtApproach approachV>
376 typename std::enable_if<approachV == OilPvtApproach::BrineCo2, BrineCo2Pvt<Scalar> >::type& getRealPvt()
377 {
378 assert(approach() == approachV);
379 return *static_cast<BrineCo2Pvt<Scalar>* >(realOilPvt_);
380 }
381
382 template <OilPvtApproach approachV>
383 typename std::enable_if<approachV == OilPvtApproach::BrineCo2, const BrineCo2Pvt<Scalar> >::type& getRealPvt() const
384 {
385 assert(approach() == approachV);
386 return *static_cast<const BrineCo2Pvt<Scalar>* >(realOilPvt_);
387 }
388
389 const void* realOilPvt() const { return realOilPvt_; }
390
391 template <OilPvtApproach approachV>
392 typename std::enable_if<approachV == OilPvtApproach::BrineH2, BrineH2Pvt<Scalar> >::type& getRealPvt()
393 {
394 assert(approach() == approachV);
395 return *static_cast<BrineH2Pvt<Scalar>* >(realOilPvt_);
396 }
397
398 template <OilPvtApproach approachV>
399 typename std::enable_if<approachV == OilPvtApproach::BrineH2, const BrineH2Pvt<Scalar> >::type& getRealPvt() const
400 {
401 assert(approach() == approachV);
402 return *static_cast<const BrineH2Pvt<Scalar>* >(realOilPvt_);
403 }
404
405 OilPvtMultiplexer<Scalar,enableThermal>& operator=(const OilPvtMultiplexer<Scalar,enableThermal>& data)
406 {
407 approach_ = data.approach_;
408 switch (approach_) {
409 case OilPvtApproach::ConstantCompressibilityOil:
410 realOilPvt_ = new ConstantCompressibilityOilPvt<Scalar>(*static_cast<const ConstantCompressibilityOilPvt<Scalar>*>(data.realOilPvt_));
411 break;
412 case OilPvtApproach::DeadOil:
413 realOilPvt_ = new DeadOilPvt<Scalar>(*static_cast<const DeadOilPvt<Scalar>*>(data.realOilPvt_));
414 break;
415 case OilPvtApproach::LiveOil:
416 realOilPvt_ = new LiveOilPvt<Scalar>(*static_cast<const LiveOilPvt<Scalar>*>(data.realOilPvt_));
417 break;
418 case OilPvtApproach::ThermalOil:
419 realOilPvt_ = new OilPvtThermal<Scalar>(*static_cast<const OilPvtThermal<Scalar>*>(data.realOilPvt_));
420 break;
421 case OilPvtApproach::BrineCo2:
422 realOilPvt_ = new BrineCo2Pvt<Scalar>(*static_cast<const BrineCo2Pvt<Scalar>*>(data.realOilPvt_));
423 break;
424 case OilPvtApproach::BrineH2:
425 realOilPvt_ = new BrineH2Pvt<Scalar>(*static_cast<const BrineH2Pvt<Scalar>*>(data.realOilPvt_));
426 break;
427 default:
428 break;
429 }
430
431 return *this;
432 }
433
434private:
435 OilPvtApproach approach_;
436 void* realOilPvt_;
437};
438
439} // namespace Opm
440
441#endif
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a CO2-Brine s...
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a H2-Brine sy...
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
This class represents the Pressure-Volume-Temperature relations of the oil phas with dissolved gas.
This class implements temperature dependence of the PVT properties of oil.
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a CO2-Brine s...
Definition BrineCo2Pvt.hpp:61
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a H2-Brine sy...
Definition BrineH2Pvt.hpp:53
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
Definition ConstantCompressibilityOilPvt.hpp:46
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
Definition DeadOilPvt.hpp:45
Definition EclipseState.hpp:63
This class represents the Pressure-Volume-Temperature relations of the oil phas with dissolved gas.
Definition LiveOilPvt.hpp:51
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
Definition OilPvtMultiplexer.hpp:104
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition OilPvtMultiplexer.hpp:170
Evaluation internalEnergy(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rs) const
Returns the specific enthalpy [J/kg] oil given a set of parameters.
Definition OilPvtMultiplexer.hpp:188
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:217
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:208
OilPvtApproach approach() const
Returns the concrete approach for calculating the PVT relations.
Definition OilPvtMultiplexer.hpp:315
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition OilPvtMultiplexer.hpp:227
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:236
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:198
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:260
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 saturated oil.
Definition OilPvtMultiplexer.hpp:245
const Scalar oilReferenceDensity(unsigned regionIdx) const
Return the reference density which are considered by this PVT-object.
Definition OilPvtMultiplexer.hpp:181
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:269
This class implements temperature dependence of the PVT properties of oil.
Definition OilPvtThermal.hpp:50
Definition Schedule.hpp:88
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30