My Project
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
36namespace Opm {
37
38#if HAVE_ECL_INPUT
39class EclipseState;
40class Schedule;
41#endif
42
43#define OPM_OIL_PVT_MULTIPLEXER_CALL(codeToCall) \
44 switch (approach_) { \
45 case OilPvtApproach::ConstantCompressibilityOil: { \
46 auto& pvtImpl = getRealPvt<OilPvtApproach::ConstantCompressibilityOil>(); \
47 codeToCall; \
48 break; \
49 } \
50 case OilPvtApproach::DeadOil: { \
51 auto& pvtImpl = getRealPvt<OilPvtApproach::DeadOil>(); \
52 codeToCall; \
53 break; \
54 } \
55 case OilPvtApproach::LiveOil: { \
56 auto& pvtImpl = getRealPvt<OilPvtApproach::LiveOil>(); \
57 codeToCall; \
58 break; \
59 } \
60 case OilPvtApproach::ThermalOil: { \
61 auto& pvtImpl = getRealPvt<OilPvtApproach::ThermalOil>(); \
62 codeToCall; \
63 break; \
64 } \
65 case OilPvtApproach::BrineCo2: { \
66 auto& pvtImpl = getRealPvt<OilPvtApproach::BrineCo2>(); \
67 codeToCall; \
68 break; \
69 } \
70 case OilPvtApproach::NoOil: \
71 throw std::logic_error("Not implemented: Oil PVT of this deck!"); \
72 }
73
74enum class OilPvtApproach {
75 NoOil,
76 LiveOil,
77 DeadOil,
78 ConstantCompressibilityOil,
79 ThermalOil,
80 BrineCo2
81};
82
95template <class Scalar, bool enableThermal = true>
97{
98public:
100 {
101 approach_ = OilPvtApproach::NoOil;
102 realOilPvt_ = nullptr;
103 }
104
105 OilPvtMultiplexer(OilPvtApproach approach, void* realOilPvt)
106 : approach_(approach)
107 , realOilPvt_(realOilPvt)
108 { }
109
111 {
112 *this = data;
113 }
114
116 {
117 switch (approach_) {
118 case OilPvtApproach::LiveOil: {
119 delete &getRealPvt<OilPvtApproach::LiveOil>();
120 break;
121 }
122 case OilPvtApproach::DeadOil: {
123 delete &getRealPvt<OilPvtApproach::DeadOil>();
124 break;
125 }
126 case OilPvtApproach::ConstantCompressibilityOil: {
127 delete &getRealPvt<OilPvtApproach::ConstantCompressibilityOil>();
128 break;
129 }
130 case OilPvtApproach::ThermalOil: {
131 delete &getRealPvt<OilPvtApproach::ThermalOil>();
132 break;
133 }
134 case OilPvtApproach::BrineCo2: {
135 delete &getRealPvt<OilPvtApproach::BrineCo2>();
136 break;
137 }
138
139 case OilPvtApproach::NoOil:
140 break;
141 }
142 }
143
144#if HAVE_ECL_INPUT
150 void initFromState(const EclipseState& eclState, const Schedule& schedule);
151#endif // HAVE_ECL_INPUT
152
153
154 void initEnd()
155 { OPM_OIL_PVT_MULTIPLEXER_CALL(pvtImpl.initEnd()); }
156
160 unsigned numRegions() const
161 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.numRegions()); return 1; }
162
166 const Scalar oilReferenceDensity(unsigned regionIdx) const
167 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.oilReferenceDensity(regionIdx)); return 700.; }
168
172 template <class Evaluation>
173 Evaluation internalEnergy(unsigned regionIdx,
174 const Evaluation& temperature,
175 const Evaluation& pressure,
176 const Evaluation& Rs) const
177 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.internalEnergy(regionIdx, temperature, pressure, Rs)); return 0; }
178
182 template <class Evaluation>
183 Evaluation viscosity(unsigned regionIdx,
184 const Evaluation& temperature,
185 const Evaluation& pressure,
186 const Evaluation& Rs) const
187 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.viscosity(regionIdx, temperature, pressure, Rs)); return 0; }
188
192 template <class Evaluation>
193 Evaluation saturatedViscosity(unsigned regionIdx,
194 const Evaluation& temperature,
195 const Evaluation& pressure) const
196 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedViscosity(regionIdx, temperature, pressure)); return 0; }
197
201 template <class Evaluation>
202 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
203 const Evaluation& temperature,
204 const Evaluation& pressure,
205 const Evaluation& Rs) const
206 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.inverseFormationVolumeFactor(regionIdx, temperature, pressure, Rs)); return 0; }
207
211 template <class Evaluation>
212 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
213 const Evaluation& temperature,
214 const Evaluation& pressure) const
215 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedInverseFormationVolumeFactor(regionIdx, temperature, pressure)); return 0; }
216
220 template <class Evaluation>
221 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
222 const Evaluation& temperature,
223 const Evaluation& pressure) const
224 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedGasDissolutionFactor(regionIdx, temperature, pressure)); return 0; }
225
229 template <class Evaluation>
230 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
231 const Evaluation& temperature,
232 const Evaluation& pressure,
233 const Evaluation& oilSaturation,
234 const Evaluation& maxOilSaturation) const
235 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedGasDissolutionFactor(regionIdx, temperature, pressure, oilSaturation, maxOilSaturation)); return 0; }
236
244 template <class Evaluation>
245 Evaluation saturationPressure(unsigned regionIdx,
246 const Evaluation& temperature,
247 const Evaluation& Rs) const
248 { OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.saturationPressure(regionIdx, temperature, Rs)); return 0; }
249
253 template <class Evaluation>
254 Evaluation diffusionCoefficient(const Evaluation& temperature,
255 const Evaluation& pressure,
256 unsigned compIdx) const
257 {
258 OPM_OIL_PVT_MULTIPLEXER_CALL(return pvtImpl.diffusionCoefficient(temperature, pressure, compIdx)); return 0;
259 }
260
261 void setApproach(OilPvtApproach appr)
262 {
263 switch (appr) {
264 case OilPvtApproach::LiveOil:
265 realOilPvt_ = new LiveOilPvt<Scalar>;
266 break;
267
268 case OilPvtApproach::DeadOil:
269 realOilPvt_ = new DeadOilPvt<Scalar>;
270 break;
271
272 case OilPvtApproach::ConstantCompressibilityOil:
274 break;
275
276 case OilPvtApproach::ThermalOil:
277 realOilPvt_ = new OilPvtThermal<Scalar>;
278 break;
279
280 case OilPvtApproach::BrineCo2:
281 realOilPvt_ = new BrineCo2Pvt<Scalar>;
282 break;
283
284 case OilPvtApproach::NoOil:
285 throw std::logic_error("Not implemented: Oil PVT of this deck!");
286 }
287
288 approach_ = appr;
289 }
290
296 OilPvtApproach approach() const
297 { return approach_; }
298
299 // get the concrete parameter object for the oil phase
300 template <OilPvtApproach approachV>
301 typename std::enable_if<approachV == OilPvtApproach::LiveOil, LiveOilPvt<Scalar> >::type& getRealPvt()
302 {
303 assert(approach() == approachV);
304 return *static_cast<LiveOilPvt<Scalar>* >(realOilPvt_);
305 }
306
307 template <OilPvtApproach approachV>
308 typename std::enable_if<approachV == OilPvtApproach::LiveOil, const LiveOilPvt<Scalar> >::type& getRealPvt() const
309 {
310 assert(approach() == approachV);
311 return *static_cast<LiveOilPvt<Scalar>* >(realOilPvt_);
312 }
313
314 template <OilPvtApproach approachV>
315 typename std::enable_if<approachV == OilPvtApproach::DeadOil, DeadOilPvt<Scalar> >::type& getRealPvt()
316 {
317 assert(approach() == approachV);
318 return *static_cast<DeadOilPvt<Scalar>* >(realOilPvt_);
319 }
320
321 template <OilPvtApproach approachV>
322 typename std::enable_if<approachV == OilPvtApproach::DeadOil, const DeadOilPvt<Scalar> >::type& getRealPvt() const
323 {
324 assert(approach() == approachV);
325 return *static_cast<DeadOilPvt<Scalar>* >(realOilPvt_);
326 }
327
328 template <OilPvtApproach approachV>
329 typename std::enable_if<approachV == OilPvtApproach::ConstantCompressibilityOil, ConstantCompressibilityOilPvt<Scalar> >::type& getRealPvt()
330 {
331 assert(approach() == approachV);
332 return *static_cast<ConstantCompressibilityOilPvt<Scalar>* >(realOilPvt_);
333 }
334
335 template <OilPvtApproach approachV>
336 typename std::enable_if<approachV == OilPvtApproach::ConstantCompressibilityOil, const ConstantCompressibilityOilPvt<Scalar> >::type& getRealPvt() const
337 {
338 assert(approach() == approachV);
339 return *static_cast<ConstantCompressibilityOilPvt<Scalar>* >(realOilPvt_);
340 }
341
342 template <OilPvtApproach approachV>
343 typename std::enable_if<approachV == OilPvtApproach::ThermalOil, OilPvtThermal<Scalar> >::type& getRealPvt()
344 {
345 assert(approach() == approachV);
346 return *static_cast<OilPvtThermal<Scalar>* >(realOilPvt_);
347 }
348
349 template <OilPvtApproach approachV>
350 typename std::enable_if<approachV == OilPvtApproach::ThermalOil, const OilPvtThermal<Scalar> >::type& getRealPvt() const
351 {
352 assert(approach() == approachV);
353 return *static_cast<const OilPvtThermal<Scalar>* >(realOilPvt_);
354 }
355
356 template <OilPvtApproach approachV>
357 typename std::enable_if<approachV == OilPvtApproach::BrineCo2, BrineCo2Pvt<Scalar> >::type& getRealPvt()
358 {
359 assert(approach() == approachV);
360 return *static_cast<BrineCo2Pvt<Scalar>* >(realOilPvt_);
361 }
362
363 template <OilPvtApproach approachV>
364 typename std::enable_if<approachV == OilPvtApproach::BrineCo2, const BrineCo2Pvt<Scalar> >::type& getRealPvt() const
365 {
366 assert(approach() == approachV);
367 return *static_cast<const BrineCo2Pvt<Scalar>* >(realOilPvt_);
368 }
369
370 const void* realOilPvt() const { return realOilPvt_; }
371
372 OilPvtMultiplexer<Scalar,enableThermal>& operator=(const OilPvtMultiplexer<Scalar,enableThermal>& data)
373 {
374 approach_ = data.approach_;
375 switch (approach_) {
376 case OilPvtApproach::ConstantCompressibilityOil:
377 realOilPvt_ = new ConstantCompressibilityOilPvt<Scalar>(*static_cast<const ConstantCompressibilityOilPvt<Scalar>*>(data.realOilPvt_));
378 break;
379 case OilPvtApproach::DeadOil:
380 realOilPvt_ = new DeadOilPvt<Scalar>(*static_cast<const DeadOilPvt<Scalar>*>(data.realOilPvt_));
381 break;
382 case OilPvtApproach::LiveOil:
383 realOilPvt_ = new LiveOilPvt<Scalar>(*static_cast<const LiveOilPvt<Scalar>*>(data.realOilPvt_));
384 break;
385 case OilPvtApproach::ThermalOil:
386 realOilPvt_ = new OilPvtThermal<Scalar>(*static_cast<const OilPvtThermal<Scalar>*>(data.realOilPvt_));
387 break;
388 case OilPvtApproach::BrineCo2:
389 realOilPvt_ = new BrineCo2Pvt<Scalar>(*static_cast<const BrineCo2Pvt<Scalar>*>(data.realOilPvt_));
390 break;
391 default:
392 break;
393 }
394
395 return *this;
396 }
397
398private:
399 OilPvtApproach approach_;
400 void* realOilPvt_;
401};
402
403} // namespace Opm
404
405#endif
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a CO2-Brine s...
Definition: BrineCo2Pvt.hpp:57
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:55
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:97
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: OilPvtMultiplexer.hpp:160
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:173
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
OilPvtApproach approach() const
Returns the concrete approach for calculating the PVT relations.
Definition: OilPvtMultiplexer.hpp:296
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
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:230
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
Definition: Schedule.hpp:130
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:30