27 #ifndef OPM_ECL_MULTIPLEXER_MATERIAL_PARAMS_HPP
28 #define OPM_ECL_MULTIPLEXER_MATERIAL_PARAMS_HPP
37 #include <type_traits>
43 enum class EclMultiplexerApproach {
58 template<
class Traits,
class GasOilMaterialLawT,
class OilWaterMaterialLawT,
class GasWaterMaterialLawT>
61 using Scalar =
typename Traits::Scalar;
62 enum { numPhases = 3 };
69 using Stone1Params =
typename Stone1Material::Params;
70 using Stone2Params =
typename Stone2Material::Params;
71 using DefaultParams =
typename DefaultMaterial::Params;
72 using TwoPhaseParams =
typename TwoPhaseMaterial::Params;
74 template <
class ParamT>
77 inline void operator () (
void* ptr )
79 delete static_cast< ParamT*
> (ptr);
83 using ParamPointerType = std::shared_ptr<void>;
98 setApproach( other.approach() );
104 setApproach( other.approach() );
108 void setApproach(EclMultiplexerApproach newApproach)
110 assert(realParams_ == 0);
111 approach_ = newApproach;
113 switch (approach()) {
114 case EclMultiplexerApproach::EclStone1Approach:
115 realParams_ = ParamPointerType(
new Stone1Params, Deleter< Stone1Params > () );
118 case EclMultiplexerApproach::EclStone2Approach:
119 realParams_ = ParamPointerType(
new Stone2Params, Deleter< Stone2Params > () );
122 case EclMultiplexerApproach::EclDefaultApproach:
123 realParams_ = ParamPointerType(
new DefaultParams, Deleter< DefaultParams > () );
126 case EclMultiplexerApproach::EclTwoPhaseApproach:
127 realParams_ = ParamPointerType(
new TwoPhaseParams, Deleter< TwoPhaseParams > () );
130 case EclMultiplexerApproach::EclOnePhaseApproach:
136 EclMultiplexerApproach approach()
const
137 {
return approach_; }
140 template <EclMultiplexerApproach approachV>
141 typename std::enable_if<approachV == EclMultiplexerApproach::EclStone1Approach, Stone1Params>::type&
144 assert(approach() == approachV);
145 return this->
template castTo<Stone1Params>();
148 template <EclMultiplexerApproach approachV>
149 typename std::enable_if<approachV == EclMultiplexerApproach::EclStone1Approach, const Stone1Params>::type&
150 getRealParams()
const
152 assert(approach() == approachV);
153 return this->
template castTo<Stone1Params>();
157 template <EclMultiplexerApproach approachV>
158 typename std::enable_if<approachV == EclMultiplexerApproach::EclStone2Approach, Stone2Params>::type&
161 assert(approach() == approachV);
162 return this->
template castTo<Stone2Params>();
165 template <EclMultiplexerApproach approachV>
166 typename std::enable_if<approachV == EclMultiplexerApproach::EclStone2Approach, const Stone2Params>::type&
167 getRealParams()
const
169 assert(approach() == approachV);
170 return this->
template castTo<Stone2Params>();
174 template <EclMultiplexerApproach approachV>
175 typename std::enable_if<approachV == EclMultiplexerApproach::EclDefaultApproach, DefaultParams>::type&
178 assert(approach() == approachV);
179 return this->
template castTo<DefaultParams>();
182 template <EclMultiplexerApproach approachV>
183 typename std::enable_if<approachV == EclMultiplexerApproach::EclDefaultApproach, const DefaultParams>::type&
184 getRealParams()
const
186 assert(approach() == approachV);
187 return this->
template castTo<DefaultParams>();
191 template <EclMultiplexerApproach approachV>
192 typename std::enable_if<approachV == EclMultiplexerApproach::EclTwoPhaseApproach, TwoPhaseParams>::type&
195 assert(approach() == approachV);
196 return this->
template castTo<TwoPhaseParams>();
199 template <EclMultiplexerApproach approachV>
200 typename std::enable_if<approachV == EclMultiplexerApproach::EclTwoPhaseApproach, const TwoPhaseParams>::type&
201 getRealParams()
const
203 assert(approach() == approachV);
204 return this->
template castTo<TwoPhaseParams>();
208 template <
class ParamT>
211 return *(
static_cast<ParamT *
> (realParams_.operator->()));
214 template <
class ParamT>
215 const ParamT& castTo()
const
217 return *(
static_cast<const ParamT *
> (realParams_.operator->()));
220 EclMultiplexerApproach approach_;
221 ParamPointerType realParams_;
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Default implementation for asserting finalization of parameter objects.
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclDefaultMaterial.hpp:61
Multiplexer implementation for the parameters required by the multiplexed three-phase material law.
Definition: EclMultiplexerMaterialParams.hpp:60
EclMultiplexerMaterialParams()
The multiplexer constructor.
Definition: EclMultiplexerMaterialParams.hpp:91
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition: EclStone1Material.hpp:60
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition: EclStone2Material.hpp:61
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Definition: EclTwoPhaseMaterial.hpp:57
Default implementation for asserting finalization of parameter objects.
Definition: EnsureFinalized.hpp:47
void finalize()
Mark the object as finalized.
Definition: EnsureFinalized.hpp:75