00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00027 #ifndef OPM_GAS_PVT_MULTIPLEXER_HPP
00028 #define OPM_GAS_PVT_MULTIPLEXER_HPP
00029
00030 #include "DryGasPvt.hpp"
00031 #include "WetGasPvt.hpp"
00032 #include "GasPvtThermal.hpp"
00033
00034 #if HAVE_OPM_PARSER
00035 #include <opm/parser/eclipse/Deck/Deck.hpp>
00036 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
00037 #include <opm/parser/eclipse/Deck/DeckKeyword.hpp>
00038 #include <opm/parser/eclipse/Deck/DeckRecord.hpp>
00039 #endif
00040
00041 namespace Opm {
00042 #define OPM_GAS_PVT_MULTIPLEXER_CALL(codeToCall) \
00043 switch (gasPvtApproach_) { \
00044 case DryGasPvt: { \
00045 auto& pvtImpl = getRealPvt<DryGasPvt>(); \
00046 codeToCall; \
00047 break; \
00048 } \
00049 case WetGasPvt: { \
00050 auto& pvtImpl = getRealPvt<WetGasPvt>(); \
00051 codeToCall; \
00052 break; \
00053 } \
00054 case ThermalGasPvt: { \
00055 auto& pvtImpl = getRealPvt<ThermalGasPvt>(); \
00056 codeToCall; \
00057 break; \
00058 } \
00059 case NoGasPvt: \
00060 OPM_THROW(std::logic_error, "Not implemented: Gas PVT of this deck!"); \
00061 }
00062
00063
00074 template <class Scalar, bool enableThermal = true>
00075 class GasPvtMultiplexer
00076 {
00077 public:
00078 typedef Opm::GasPvtThermal<Scalar> GasPvtThermal;
00079
00080 enum GasPvtApproach {
00081 NoGasPvt,
00082 DryGasPvt,
00083 WetGasPvt,
00084 ThermalGasPvt
00085 };
00086
00087 GasPvtMultiplexer()
00088 {
00089 gasPvtApproach_ = NoGasPvt;
00090 }
00091
00092 ~GasPvtMultiplexer()
00093 {
00094 switch (gasPvtApproach_) {
00095 case DryGasPvt: {
00096 delete &getRealPvt<DryGasPvt>();
00097 break;
00098 }
00099 case WetGasPvt: {
00100 delete &getRealPvt<WetGasPvt>();
00101 break;
00102 }
00103 case ThermalGasPvt: {
00104 delete &getRealPvt<ThermalGasPvt>();
00105 break;
00106 }
00107 case NoGasPvt:
00108 break;
00109 }
00110 }
00111
00112 #if HAVE_OPM_PARSER
00113
00118 void initFromDeck(const Deck& deck, const EclipseState& eclState)
00119 {
00120 bool enableGas = deck.hasKeyword("GAS");
00121 if (!enableGas)
00122 return;
00123
00124 if (enableThermal
00125 && (deck.hasKeyword("TREF")
00126 || deck.hasKeyword("GASVISCT")))
00127 setApproach(ThermalGasPvt);
00128 else if (deck.hasKeyword("PVTG"))
00129 setApproach(WetGasPvt);
00130 else if (deck.hasKeyword("PVDG"))
00131 setApproach(DryGasPvt);
00132
00133 OPM_GAS_PVT_MULTIPLEXER_CALL(pvtImpl.initFromDeck(deck, eclState));
00134 }
00135 #endif // HAVE_OPM_PARSER
00136
00137 void setApproach(GasPvtApproach gasPvtAppr)
00138 {
00139 switch (gasPvtAppr) {
00140 case DryGasPvt:
00141 realGasPvt_ = new Opm::DryGasPvt<Scalar>;
00142 break;
00143
00144 case WetGasPvt:
00145 realGasPvt_ = new Opm::WetGasPvt<Scalar>;
00146 break;
00147
00148 case ThermalGasPvt:
00149 realGasPvt_ = new Opm::GasPvtThermal<Scalar>;
00150 break;
00151
00152 case NoGasPvt:
00153 OPM_THROW(std::logic_error, "Not implemented: Gas PVT of this deck!");
00154 }
00155
00156 gasPvtApproach_ = gasPvtAppr;
00157 }
00158
00159 void initEnd()
00160 { OPM_GAS_PVT_MULTIPLEXER_CALL(pvtImpl.initEnd()); }
00161
00165 unsigned numRegions() const
00166 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.numRegions()); return 1; }
00167
00171 template <class Evaluation = Scalar>
00172 Evaluation viscosity(unsigned regionIdx,
00173 const Evaluation& temperature,
00174 const Evaluation& pressure,
00175 const Evaluation& Rv) const
00176 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.viscosity(regionIdx, temperature, pressure, Rv)); return 0; }
00177
00181 template <class Evaluation = Scalar>
00182 Evaluation saturatedViscosity(unsigned regionIdx,
00183 const Evaluation& temperature,
00184 const Evaluation& pressure) const
00185 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedViscosity(regionIdx, temperature, pressure)); return 0; }
00186
00190 template <class Evaluation = Scalar>
00191 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
00192 const Evaluation& temperature,
00193 const Evaluation& pressure,
00194 const Evaluation& Rv) const
00195 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.inverseFormationVolumeFactor(regionIdx, temperature, pressure, Rv)); return 0; }
00196
00200 template <class Evaluation = Scalar>
00201 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
00202 const Evaluation& temperature,
00203 const Evaluation& pressure) const
00204 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedInverseFormationVolumeFactor(regionIdx, temperature, pressure)); return 0; }
00205
00209 template <class Evaluation = Scalar>
00210 Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
00211 const Evaluation& temperature,
00212 const Evaluation& pressure) const
00213 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedOilVaporizationFactor(regionIdx, temperature, pressure)); return 0; }
00214
00218 template <class Evaluation = Scalar>
00219 Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
00220 const Evaluation& temperature,
00221 const Evaluation& pressure,
00222 const Evaluation& oilSaturation,
00223 Scalar maxOilSaturation) const
00224 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturatedOilVaporizationFactor(regionIdx, temperature, pressure, oilSaturation, maxOilSaturation)); return 0; }
00225
00232 template <class Evaluation = Scalar>
00233 Evaluation saturationPressure(unsigned regionIdx,
00234 const Evaluation& temperature,
00235 const Evaluation& Rv) const
00236 { OPM_GAS_PVT_MULTIPLEXER_CALL(return pvtImpl.saturationPressure(regionIdx, temperature, Rv)); return 0; }
00237
00243 GasPvtApproach gasPvtApproach() const
00244 { return gasPvtApproach_; }
00245
00246
00247 template <GasPvtApproach approachV>
00248 typename std::enable_if<approachV == DryGasPvt, Opm::DryGasPvt<Scalar> >::type& getRealPvt()
00249 {
00250 assert(gasPvtApproach() == approachV);
00251 return *static_cast<Opm::DryGasPvt<Scalar>* >(realGasPvt_);
00252 }
00253
00254 template <GasPvtApproach approachV>
00255 typename std::enable_if<approachV == DryGasPvt, const Opm::DryGasPvt<Scalar> >::type& getRealPvt() const
00256 {
00257 assert(gasPvtApproach() == approachV);
00258 return *static_cast<const Opm::DryGasPvt<Scalar>* >(realGasPvt_);
00259 }
00260
00261
00262 template <GasPvtApproach approachV>
00263 typename std::enable_if<approachV == WetGasPvt, Opm::WetGasPvt<Scalar> >::type& getRealPvt()
00264 {
00265 assert(gasPvtApproach() == approachV);
00266 return *static_cast<Opm::WetGasPvt<Scalar>* >(realGasPvt_);
00267 }
00268
00269 template <GasPvtApproach approachV>
00270 typename std::enable_if<approachV == WetGasPvt, const Opm::WetGasPvt<Scalar> >::type& getRealPvt() const
00271 {
00272 assert(gasPvtApproach() == approachV);
00273 return *static_cast<const Opm::WetGasPvt<Scalar>* >(realGasPvt_);
00274 }
00275
00276
00277 template <GasPvtApproach approachV>
00278 typename std::enable_if<approachV == ThermalGasPvt, Opm::GasPvtThermal<Scalar> >::type& getRealPvt()
00279 {
00280 assert(gasPvtApproach() == approachV);
00281 return *static_cast<Opm::GasPvtThermal<Scalar>* >(realGasPvt_);
00282 }
00283
00284 template <GasPvtApproach approachV>
00285 typename std::enable_if<approachV == ThermalGasPvt, const Opm::GasPvtThermal<Scalar> >::type& getRealPvt() const
00286 {
00287 assert(gasPvtApproach() == approachV);
00288 return *static_cast<const Opm::GasPvtThermal<Scalar>* >(realGasPvt_);
00289 }
00290
00291 private:
00292 GasPvtApproach gasPvtApproach_;
00293 void* realGasPvt_;
00294 };
00295
00296 #undef OPM_GAS_PVT_MULTIPLEXER_CALL
00297
00298 }
00299
00300 #endif