25 #ifndef OPM_LINEAR_MATERIAL_HPP
26 #define OPM_LINEAR_MATERIAL_HPP
32 #include <opm/common/Exceptions.hpp>
33 #include <opm/common/ErrorMacros.hpp>
36 #include <type_traits>
50 template <
class TraitsT,
class ParamsT = LinearMaterialParams<TraitsT> >
56 typedef typename Traits::Scalar
Scalar;
97 template <
class ContainerT,
class Flu
idState>
100 const FluidState &state)
102 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
105 for (
unsigned phaseIdx = 0; phaseIdx < Traits::numPhases; ++phaseIdx) {
106 const Evaluation& S =
107 FsToolbox::template toLhs<Evaluation>(state.saturation(phaseIdx));
111 S*params.pcMaxSat(phaseIdx) +
112 (1.0 - S)*params.pcMinSat(phaseIdx);
119 template <
class ContainerT,
class Flu
idState>
124 OPM_THROW(std::runtime_error,
"Not implemented: LinearMaterial::saturations()");
130 template <
class ContainerT,
class Flu
idState>
133 const FluidState &state)
135 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
139 for (
unsigned phaseIdx = 0; phaseIdx < Traits::numPhases; ++phaseIdx) {
140 const Evaluation& S =
141 FsToolbox::template toLhs<Evaluation>(state.saturation(phaseIdx));
151 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
152 static Evaluation
pcnw(
const Params ¶ms,
const FluidState &fs)
155 const Evaluation&
Sw =
156 FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
159 const Evaluation& wPhasePressure =
160 Sw*params.pcMaxSat(Traits::wettingPhaseIdx) +
161 (1.0 -
Sw)*params.pcMinSat(Traits::wettingPhaseIdx);
163 const Evaluation&
Sn =
164 FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
167 const Evaluation& nPhasePressure =
168 Sn*params.pcMaxSat(Traits::nonWettingPhaseIdx) +
169 (1.0 -
Sn)*params.pcMinSat(Traits::nonWettingPhaseIdx);
171 return nPhasePressure - wPhasePressure;
175 template <
class Evaluation = Scalar>
176 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
179 const Evaluation& wPhasePressure =
180 Sw*params.pcMaxSat(Traits::wettingPhaseIdx) +
181 (1.0 -
Sw)*params.pcMinSat(Traits::wettingPhaseIdx);
183 const Evaluation& nPhasePressure =
184 (1.0 -
Sw)*params.pcMaxSat(Traits::nonWettingPhaseIdx) +
185 Sw*params.pcMinSat(Traits::nonWettingPhaseIdx);
187 return nPhasePressure - wPhasePressure;
194 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
195 static Evaluation
Sw(
const Params &,
const FluidState &)
196 { OPM_THROW(std::runtime_error,
"Not implemented: Sw()"); }
198 template <
class Evaluation = Scalar>
199 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
201 { OPM_THROW(std::runtime_error,
"Not implemented: twoPhaseSatSw()"); }
207 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
208 static Evaluation
Sn(
const Params &,
const FluidState &)
209 { OPM_THROW(std::runtime_error,
"Not implemented: Sn()"); }
211 template <
class Evaluation = Scalar>
212 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
214 { OPM_THROW(std::runtime_error,
"Not implemented: twoPhaseSatSn()"); }
222 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
223 static typename std::enable_if<Traits::numPhases == 3, Evaluation>::type
224 Sg(
const Params &,
const FluidState &)
225 { OPM_THROW(std::runtime_error,
"Not implemented: Sg()"); }
230 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
231 static Evaluation
krw(
const Params &,
const FluidState &fs)
236 const Evaluation&
Sw =
237 FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
241 template <
class Evaluation = Scalar>
242 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
253 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
254 static Evaluation
krn(
const Params &,
const FluidState &fs)
259 const Evaluation&
Sn =
260 FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
264 template <
class Evaluation = Scalar>
265 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
278 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
279 static typename std::enable_if<Traits::numPhases == 3, Evaluation>::type
280 krg(
const Params &,
const FluidState &fs)
285 const Evaluation&
Sg =
286 FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::gasPhaseIdx));
295 template <
class Flu
idState,
class Evaluation = Scalar>
296 static typename std::enable_if<Traits::numPhases == 3, Evaluation>::type
297 pcgn(
const Params ¶ms,
const FluidState &fs)
301 const Evaluation&
Sn =
302 FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
305 const Evaluation& nPhasePressure =
306 Sn*params.pcMaxSat(Traits::nonWettingPhaseIdx) +
307 (1.0 -
Sn)*params.pcMinSat(Traits::nonWettingPhaseIdx);
309 const Evaluation&
Sg =
310 FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::gasPhaseIdx));
313 const Evaluation& gPhasePressure =
314 Sg*params.pcMaxSat(Traits::gasPhaseIdx) +
315 (1.0 -
Sg)*params.pcMinSat(Traits::gasPhaseIdx);
317 return gPhasePressure - nPhasePressure;
Reference implementation of params for the linear M-phase material material.
bool CheckDefined(const T &value OPM_UNUSED)
Make valgrind complain if any of the memory occupied by an object is undefined.
Definition: Valgrind.hpp:74
static std::enable_if< Traits::numPhases==2, Evaluation >::type twoPhaseSatSw(const Params &, const Evaluation &)
Definition: LinearMaterial.hpp:200
static std::enable_if< Traits::numPhases==2, Evaluation >::type twoPhaseSatKrw(const Params &, const Evaluation &Sw)
Definition: LinearMaterial.hpp:243
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: LinearMaterial.hpp:120
Definition: Air_Mesitylene.hpp:31
static void relativePermeabilities(ContainerT &values, const Params &, const FluidState &state)
The relative permeability of all phases.
Definition: LinearMaterial.hpp:131
static Evaluation Sw(const Params &, const FluidState &)
Calculate wetting phase saturation given that the rest of the fluid state has been initialized...
Definition: LinearMaterial.hpp:195
static void capillaryPressures(ContainerT &values, const Params ¶ms, const FluidState &state)
The linear capillary pressure-saturation curve.
Definition: LinearMaterial.hpp:98
static const bool isSaturationDependent
Definition: LinearMaterial.hpp:71
Some templates to wrap the valgrind client request macros.
Evaluation< Scalar, VarSetTag, numVars > max(const Evaluation< Scalar, VarSetTag, numVars > &x1, const Evaluation< Scalar, VarSetTag, numVars > &x2)
Definition: Math.hpp:114
static std::enable_if< Traits::numPhases==3, Evaluation >::type krg(const Params &, const FluidState &fs)
The relative permability of the gas phase.
Definition: LinearMaterial.hpp:280
static std::enable_if< Traits::numPhases==2, Evaluation >::type twoPhaseSatPcnw(const Params ¶ms, const Evaluation &Sw)
Definition: LinearMaterial.hpp:177
static const bool isPressureDependent
Definition: LinearMaterial.hpp:75
static const bool isTemperatureDependent
Definition: LinearMaterial.hpp:79
Traits::Scalar Scalar
Definition: LinearMaterial.hpp:56
static const bool implementsTwoPhaseSatApi
Definition: LinearMaterial.hpp:67
Implements a linear saturation-capillary pressure relation.
Definition: LinearMaterial.hpp:51
Evaluation< Scalar, VarSetTag, numVars > min(const Evaluation< Scalar, VarSetTag, numVars > &x1, const Evaluation< Scalar, VarSetTag, numVars > &x2)
Definition: Math.hpp:61
ParamsT Params
Definition: LinearMaterial.hpp:55
TraitsT Traits
Definition: LinearMaterial.hpp:54
static const int numPhases
The number of fluid phases.
Definition: LinearMaterial.hpp:59
static const bool implementsTwoPhaseApi
Definition: LinearMaterial.hpp:63
static std::enable_if< Traits::numPhases==3, Evaluation >::type Sg(const Params &, const FluidState &)
Calculate gas phase saturation given that the rest of the fluid state has been initialized.
Definition: LinearMaterial.hpp:224
static std::enable_if< Traits::numPhases==2, Evaluation >::type twoPhaseSatSn(const Params &, const Evaluation &)
Definition: LinearMaterial.hpp:213
static Evaluation krw(const Params &, const FluidState &fs)
The relative permability of the wetting phase.
Definition: LinearMaterial.hpp:231
static Evaluation Sn(const Params &, const FluidState &)
Calculate non-wetting liquid phase saturation given that the rest of the fluid state has been initial...
Definition: LinearMaterial.hpp:208
static Evaluation krn(const Params &, const FluidState &fs)
The relative permability of the liquid non-wetting phase.
Definition: LinearMaterial.hpp:254
static std::enable_if< Traits::numPhases==3, Evaluation >::type pcgn(const Params ¶ms, const FluidState &fs)
The difference between the pressures of the gas and the non-wetting phase.
Definition: LinearMaterial.hpp:297
static const bool isCompositionDependent
Definition: LinearMaterial.hpp:83
static std::enable_if< Traits::numPhases==2, Evaluation >::type twoPhaseSatKrn(const Params &, const Evaluation &Sw)
Definition: LinearMaterial.hpp:266
static Evaluation pcnw(const Params ¶ms, const FluidState &fs)
The difference between the pressures of the non-wetting and wetting phase.
Definition: LinearMaterial.hpp:152