27 #ifndef OPM_LINEAR_MATERIAL_HPP
28 #define OPM_LINEAR_MATERIAL_HPP
33 #include <opm/common/Valgrind.hpp>
34 #include <opm/common/Exceptions.hpp>
35 #include <opm/common/ErrorMacros.hpp>
38 #include <type_traits>
52 template <
class TraitsT,
class ParamsT = LinearMaterialParams<TraitsT> >
56 typedef TraitsT Traits;
57 typedef ParamsT Params;
58 typedef typename Traits::Scalar Scalar;
99 template <
class ContainerT,
class Flu
idState>
101 const Params& params,
102 const FluidState& state)
104 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
106 for (
unsigned phaseIdx = 0; phaseIdx < Traits::numPhases; ++phaseIdx) {
107 const Evaluation& S =
108 Opm::decay<Evaluation>(state.saturation(phaseIdx));
109 Valgrind::CheckDefined(S);
112 S*params.pcMaxSat(phaseIdx) +
113 (1.0 - S)*params.pcMinSat(phaseIdx);
120 template <
class ContainerT,
class Flu
idState>
125 OPM_THROW(std::runtime_error,
"Not implemented: LinearMaterial::saturations()");
131 template <
class ContainerT,
class Flu
idState>
134 const FluidState& state)
136 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
138 for (
unsigned phaseIdx = 0; phaseIdx < Traits::numPhases; ++phaseIdx) {
139 const Evaluation& S =
140 Opm::decay<Evaluation>(state.saturation(phaseIdx));
141 Valgrind::CheckDefined(S);
143 values[phaseIdx] = Opm::max(Opm::min(S,1.0),0.0);
150 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
151 static Evaluation
pcnw(
const Params& params,
const FluidState& fs)
153 const Evaluation&
Sw =
154 Opm::decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
155 Valgrind::CheckDefined(Sw);
157 const Evaluation& wPhasePressure =
158 Sw*params.pcMaxSat(Traits::wettingPhaseIdx) +
159 (1.0 -
Sw)*params.pcMinSat(Traits::wettingPhaseIdx);
161 const Evaluation&
Sn =
162 Opm::decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
163 Valgrind::CheckDefined(Sn);
165 const Evaluation& nPhasePressure =
166 Sn*params.pcMaxSat(Traits::nonWettingPhaseIdx) +
167 (1.0 -
Sn)*params.pcMinSat(Traits::nonWettingPhaseIdx);
169 return nPhasePressure - wPhasePressure;
173 template <
class Evaluation = Scalar>
174 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
175 twoPhaseSatPcnw(
const Params& params,
const Evaluation&
Sw)
177 const Evaluation& wPhasePressure =
178 Sw*params.pcMaxSat(Traits::wettingPhaseIdx) +
179 (1.0 -
Sw)*params.pcMinSat(Traits::wettingPhaseIdx);
181 const Evaluation& nPhasePressure =
182 (1.0 -
Sw)*params.pcMaxSat(Traits::nonWettingPhaseIdx) +
183 Sw*params.pcMinSat(Traits::nonWettingPhaseIdx);
185 return nPhasePressure - wPhasePressure;
192 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
193 static Evaluation
Sw(
const Params& ,
const FluidState& )
194 { OPM_THROW(std::runtime_error,
"Not implemented: Sw()"); }
196 template <
class Evaluation = Scalar>
197 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
198 twoPhaseSatSw(
const Params& ,
const Evaluation& )
199 { OPM_THROW(std::runtime_error,
"Not implemented: twoPhaseSatSw()"); }
205 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
206 static Evaluation
Sn(
const Params& ,
const FluidState& )
207 { OPM_THROW(std::runtime_error,
"Not implemented: Sn()"); }
209 template <
class Evaluation = Scalar>
210 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
211 twoPhaseSatSn(
const Params& ,
const Evaluation& )
212 { OPM_THROW(std::runtime_error,
"Not implemented: twoPhaseSatSn()"); }
220 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
221 static typename std::enable_if<Traits::numPhases == 3, Evaluation>::type
222 Sg(
const Params& ,
const FluidState& )
223 { OPM_THROW(std::runtime_error,
"Not implemented: Sg()"); }
228 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
229 static Evaluation
krw(
const Params& ,
const FluidState& fs)
231 const Evaluation& Sw =
232 Opm::decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
233 return Opm::max(0.0, Opm::min(1.0, Sw));
236 template <
class Evaluation = Scalar>
237 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
238 twoPhaseSatKrw(
const Params& ,
const Evaluation& Sw)
239 {
return Opm::max(0.0, Opm::min(1.0, Sw)); }
244 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
245 static Evaluation
krn(
const Params& ,
const FluidState& fs)
247 const Evaluation&
Sn =
248 Opm::decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
249 return Opm::max(0.0, Opm::min(1.0, Sn));
252 template <
class Evaluation = Scalar>
253 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
254 twoPhaseSatKrn(
const Params& ,
const Evaluation& Sw)
256 return Opm::max(0.0, Opm::min(1.0, Sw));
264 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
265 static typename std::enable_if<Traits::numPhases == 3, Evaluation>::type
266 krg(
const Params& ,
const FluidState& fs)
268 const Evaluation&
Sg =
269 Opm::decay<Evaluation>(fs.saturation(Traits::gasPhaseIdx));
270 return Opm::max(0.0, Opm::min(1.0, Sg));
278 template <
class Flu
idState,
class Evaluation = Scalar>
279 static typename std::enable_if<Traits::numPhases == 3, Evaluation>::type
280 pcgn(
const Params& params,
const FluidState& fs)
282 const Evaluation&
Sn =
283 Opm::decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
284 Valgrind::CheckDefined(Sn);
286 const Evaluation& nPhasePressure =
287 Sn*params.pcMaxSat(Traits::nonWettingPhaseIdx) +
288 (1.0 -
Sn)*params.pcMinSat(Traits::nonWettingPhaseIdx);
290 const Evaluation&
Sg =
291 Opm::decay<Evaluation>(fs.saturation(Traits::gasPhaseIdx));
292 Valgrind::CheckDefined(Sg);
294 const Evaluation& gPhasePressure =
295 Sg*params.pcMaxSat(Traits::gasPhaseIdx) +
296 (1.0 -
Sg)*params.pcMinSat(Traits::gasPhaseIdx);
298 return gPhasePressure - nPhasePressure;
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition...
Definition: LinearMaterial.hpp:85
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: LinearMaterial.hpp:73
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: LinearMaterial.hpp:69
static const int numPhases
The number of fluid phases.
Definition: LinearMaterial.hpp:61
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: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:151
static Evaluation krw(const Params &, const FluidState &fs)
The relative permability of the wetting phase.
Definition: LinearMaterial.hpp:229
static void relativePermeabilities(ContainerT &values, const Params &, const FluidState &state)
The relative permeability of all phases.
Definition: LinearMaterial.hpp:132
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: LinearMaterial.hpp:121
Reference implementation of params for the linear M-phase material material.
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:193
static Evaluation krn(const Params &, const FluidState &fs)
The relative permability of the liquid non-wetting phase.
Definition: LinearMaterial.hpp:245
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:280
static void capillaryPressures(ContainerT &values, const Params ¶ms, const FluidState &state)
The linear capillary pressure-saturation curve.
Definition: LinearMaterial.hpp:100
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: LinearMaterial.hpp:81
Implements a linear saturation-capillary pressure relation.
Definition: LinearMaterial.hpp:53
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure...
Definition: LinearMaterial.hpp:77
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:222
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:206
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: LinearMaterial.hpp:65