27 #ifndef OPM_LINEAR_MATERIAL_HPP
28 #define OPM_LINEAR_MATERIAL_HPP
36 #include <type_traits>
50 template <
class TraitsT,
class ParamsT = LinearMaterialParams<TraitsT> >
54 typedef TraitsT Traits;
55 typedef ParamsT Params;
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;
104 for (
unsigned phaseIdx = 0; phaseIdx < Traits::numPhases; ++phaseIdx) {
105 const Evaluation& S =
106 decay<Evaluation>(state.saturation(phaseIdx));
107 Valgrind::CheckDefined(S);
110 S*params.pcMaxSat(phaseIdx) +
111 (1.0 - S)*params.pcMinSat(phaseIdx);
118 template <
class ContainerT,
class Flu
idState>
123 throw std::runtime_error(
"Not implemented: LinearMaterial::saturations()");
129 template <
class ContainerT,
class Flu
idState>
132 const FluidState& state)
134 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
136 for (
unsigned phaseIdx = 0; phaseIdx < Traits::numPhases; ++phaseIdx) {
137 const Evaluation& S =
138 decay<Evaluation>(state.saturation(phaseIdx));
139 Valgrind::CheckDefined(S);
141 values[phaseIdx] = max(min(S,1.0),0.0);
148 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
149 static Evaluation
pcnw(
const Params& params,
const FluidState& fs)
151 const Evaluation&
Sw =
152 decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
153 Valgrind::CheckDefined(
Sw);
155 const Evaluation& wPhasePressure =
156 Sw*params.pcMaxSat(Traits::wettingPhaseIdx) +
157 (1.0 -
Sw)*params.pcMinSat(Traits::wettingPhaseIdx);
159 const Evaluation&
Sn =
160 decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
161 Valgrind::CheckDefined(
Sn);
163 const Evaluation& nPhasePressure =
164 Sn*params.pcMaxSat(Traits::nonWettingPhaseIdx) +
165 (1.0 -
Sn)*params.pcMinSat(Traits::nonWettingPhaseIdx);
167 return nPhasePressure - wPhasePressure;
171 template <
class Evaluation = Scalar>
172 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
173 twoPhaseSatPcnw(
const Params& params,
const Evaluation&
Sw)
175 const Evaluation& wPhasePressure =
176 Sw*params.pcMaxSat(Traits::wettingPhaseIdx) +
177 (1.0 -
Sw)*params.pcMinSat(Traits::wettingPhaseIdx);
179 const Evaluation& nPhasePressure =
180 (1.0 -
Sw)*params.pcMaxSat(Traits::nonWettingPhaseIdx) +
181 Sw*params.pcMinSat(Traits::nonWettingPhaseIdx);
183 return nPhasePressure - wPhasePressure;
190 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
191 static Evaluation
Sw(
const Params& ,
const FluidState& )
192 {
throw std::runtime_error(
"Not implemented: Sw()"); }
194 template <
class Evaluation = Scalar>
195 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
196 twoPhaseSatSw(
const Params& ,
const Evaluation& )
197 {
throw std::runtime_error(
"Not implemented: twoPhaseSatSw()"); }
203 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
204 static Evaluation
Sn(
const Params& ,
const FluidState& )
205 {
throw std::runtime_error(
"Not implemented: Sn()"); }
207 template <
class Evaluation = Scalar>
208 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
209 twoPhaseSatSn(
const Params& ,
const Evaluation& )
210 {
throw std::runtime_error(
"Not implemented: twoPhaseSatSn()"); }
218 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
219 static typename std::enable_if<Traits::numPhases == 3, Evaluation>::type
220 Sg(
const Params& ,
const FluidState& )
221 {
throw std::runtime_error(
"Not implemented: Sg()"); }
226 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
227 static Evaluation
krw(
const Params& ,
const FluidState& fs)
229 const Evaluation&
Sw =
230 decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
231 return max(0.0, min(1.0,
Sw));
234 template <
class Evaluation = Scalar>
235 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
236 twoPhaseSatKrw(
const Params& ,
const Evaluation&
Sw)
237 {
return max(0.0, min(1.0,
Sw)); }
242 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
243 static Evaluation
krn(
const Params& ,
const FluidState& fs)
245 const Evaluation&
Sn =
246 decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
247 return max(0.0, min(1.0,
Sn));
250 template <
class Evaluation = Scalar>
251 static typename std::enable_if<Traits::numPhases == 2, Evaluation>::type
252 twoPhaseSatKrn(
const Params& ,
const Evaluation&
Sw)
254 return max(0.0, min(1.0,
Sw));
262 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
263 static typename std::enable_if<Traits::numPhases == 3, Evaluation>::type
264 krg(
const Params& ,
const FluidState& fs)
266 const Evaluation&
Sg =
267 decay<Evaluation>(fs.saturation(Traits::gasPhaseIdx));
268 return max(0.0, min(1.0,
Sg));
276 template <
class Flu
idState,
class Evaluation = Scalar>
277 static typename std::enable_if<Traits::numPhases == 3, Evaluation>::type
278 pcgn(
const Params& params,
const FluidState& fs)
280 const Evaluation&
Sn =
281 decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
282 Valgrind::CheckDefined(
Sn);
284 const Evaluation& nPhasePressure =
285 Sn*params.pcMaxSat(Traits::nonWettingPhaseIdx) +
286 (1.0 -
Sn)*params.pcMinSat(Traits::nonWettingPhaseIdx);
288 const Evaluation&
Sg =
289 decay<Evaluation>(fs.saturation(Traits::gasPhaseIdx));
290 Valgrind::CheckDefined(
Sg);
292 const Evaluation& gPhasePressure =
293 Sg*params.pcMaxSat(Traits::gasPhaseIdx) +
294 (1.0 -
Sg)*params.pcMinSat(Traits::gasPhaseIdx);
296 return gPhasePressure - nPhasePressure;
Reference implementation of params for the linear M-phase material material.
Some templates to wrap the valgrind client request macros.
Implements a linear saturation-capillary pressure relation.
Definition: LinearMaterial.hpp:52
static Evaluation pcnw(const Params ¶ms, const FluidState &fs)
The difference between the pressures of the non-wetting and wetting phase.
Definition: LinearMaterial.hpp:149
static Evaluation krn(const Params &, const FluidState &fs)
The relative permability of the liquid non-wetting phase.
Definition: LinearMaterial.hpp:243
static void relativePermeabilities(ContainerT &values, const Params &, const FluidState &state)
The relative permeability of all phases.
Definition: LinearMaterial.hpp:130
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: LinearMaterial.hpp:67
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: LinearMaterial.hpp:75
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:278
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:220
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: LinearMaterial.hpp:71
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:264
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:191
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:204
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: LinearMaterial.hpp:79
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: LinearMaterial.hpp:119
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: LinearMaterial.hpp:83
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: LinearMaterial.hpp:63
static void capillaryPressures(ContainerT &values, const Params ¶ms, const FluidState &state)
The linear capillary pressure-saturation curve.
Definition: LinearMaterial.hpp:98
static const int numPhases
The number of fluid phases.
Definition: LinearMaterial.hpp:59
static Evaluation krw(const Params &, const FluidState &fs)
The relative permability of the wetting phase.
Definition: LinearMaterial.hpp:227