27#ifndef OPM_EFF_TO_ABS_LAW_HPP
28#define OPM_EFF_TO_ABS_LAW_HPP
69template <
class EffLawT,
class ParamsT = EffToAbsLawParams<
typename EffLawT::Params, EffLawT::numPhases> >
72 typedef EffLawT EffLaw;
75 typedef typename EffLaw::Traits Traits;
76 typedef ParamsT Params;
77 typedef typename EffLaw::Scalar Scalar;
116 template <
class Container,
class Flu
idState>
121 OverlayFluidState overlayFs(fs);
122 for (
unsigned phaseIdx = 0; phaseIdx <
numPhases; ++ phaseIdx) {
123 overlayFs.setSaturation(phaseIdx,
125 fs.saturation(phaseIdx),
129 EffLaw::template capillaryPressures<Container, OverlayFluidState>(values, params, overlayFs);
142 template <
class Container,
class Flu
idState>
147 OverlayFluidState overlayFs(fs);
148 for (
unsigned phaseIdx = 0; phaseIdx <
numPhases; ++ phaseIdx) {
149 overlayFs.setSaturation(phaseIdx,
151 fs.saturation(phaseIdx),
155 EffLaw::template relativePermeabilities<Container, OverlayFluidState>(values, params, overlayFs);
169 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
170 static Evaluation
pcnw(
const Params& params,
const FluidState& fs)
174 static_assert(FluidState::numPhases ==
numPhases,
175 "The fluid state and the material law must exhibit the same "
176 "number of phases!");
178 OverlayFluidState overlayFs(fs);
179 for (
unsigned phaseIdx = 0; phaseIdx <
numPhases; ++ phaseIdx) {
180 overlayFs.setSaturation(phaseIdx,
182 fs.saturation(phaseIdx),
186 return EffLaw::template pcnw<OverlayFluidState, Evaluation>(params, overlayFs);
189 template <
class Evaluation>
190 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
191 twoPhaseSatPcnw(
const Params& params,
const Evaluation& SwAbs)
195 return EffLaw::twoPhaseSatPcnw(params, SwEff);
201 template <
class Container,
class Flu
idState>
202 static void saturations(Container& values,
const Params& params,
const FluidState& fs)
204 EffLaw::template saturations<Container, FluidState>(values, params, fs);
205 for (
unsigned phaseIdx = 0; phaseIdx <
numPhases; ++phaseIdx) {
214 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
215 static Evaluation
Sw(
const Params& params,
const FluidState& fs)
218 EffLaw::template Sw<FluidState, Evaluation>(params, fs),
219 Traits::wettingPhaseIdx);
222 template <
class Evaluation>
223 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
224 twoPhaseSatSw(
const Params& params,
const Evaluation&
Sw)
226 EffLaw::twoPhaseSatSw(params,
Sw),
227 Traits::wettingPhaseIdx); }
233 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
234 static Evaluation
Sn(
const Params& params,
const FluidState& fs)
237 EffLaw::template Sn<FluidState, Evaluation>(params, fs),
238 Traits::nonWettingPhaseIdx);
241 template <
class Evaluation>
242 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
243 twoPhaseSatSn(
const Params& params,
const Evaluation&
Sw)
246 EffLaw::twoPhaseSatSn(params,
Sw),
247 Traits::nonWettingPhaseIdx);
256 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
257 static typename std::enable_if< (Traits::numPhases > 2), Evaluation>::type
258 Sg(
const Params& params,
const FluidState& fs)
261 EffLaw::template Sg<FluidState, Evaluation>(params, fs),
262 Traits::gasPhaseIdx);
274 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
275 static Evaluation
krw(
const Params& params,
const FluidState& fs)
279 static_assert(FluidState::numPhases ==
numPhases,
280 "The fluid state and the material law must exhibit the same "
281 "number of phases!");
283 OverlayFluidState overlayFs(fs);
284 for (
unsigned phaseIdx = 0; phaseIdx <
numPhases; ++ phaseIdx) {
285 overlayFs.setSaturation(phaseIdx,
287 fs.saturation(phaseIdx),
291 return EffLaw::template krw<OverlayFluidState, Evaluation>(params, overlayFs);
294 template <
class Evaluation>
295 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
296 twoPhaseSatKrw(
const Params& params,
const Evaluation&
Sw)
302 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
303 static Evaluation
krn(
const Params& params,
const FluidState& fs)
307 static_assert(FluidState::numPhases ==
numPhases,
308 "The fluid state and the material law must exhibit the same "
309 "number of phases!");
311 OverlayFluidState overlayFs(fs);
312 for (
unsigned phaseIdx = 0; phaseIdx <
numPhases; ++ phaseIdx) {
313 overlayFs.setSaturation(phaseIdx,
315 fs.saturation(phaseIdx),
319 return EffLaw::template krn<OverlayFluidState, Evaluation>(params, overlayFs);
322 template <
class Evaluation>
323 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
324 twoPhaseSatKrn(
const Params& params,
const Evaluation&
Sw)
332 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
333 static typename std::enable_if< (Traits::numPhases > 2), Evaluation>::type
334 krg(
const Params& params,
const FluidState& fs)
338 static_assert(FluidState::numPhases ==
numPhases,
339 "The fluid state and the material law must exhibit the same "
340 "number of phases!");
342 OverlayFluidState overlayFs(fs);
343 for (
unsigned phaseIdx = 0; phaseIdx <
numPhases; ++ phaseIdx) {
344 overlayFs.setSaturation(phaseIdx,
346 fs.saturation(phaseIdx),
350 return EffLaw::template krg<OverlayFluidState, Evaluation>(params, overlayFs);
356 template <
class Evaluation>
358 {
return (S - params.residualSaturation(phaseIdx))/(1.0 - params.sumResidualSaturations()); }
363 template <
class Evaluation>
365 {
return S*(1.0 - params.sumResidualSaturations()) + params.residualSaturation(phaseIdx); }
376 static Scalar dSeff_dSabs_(
const Params& params,
int )
377 {
return 1.0/(1 - params.sumResidualSaturations()); }
387 static Scalar dSabs_dSeff_(
const Params& params,
int )
388 {
return 1 - params.sumResidualSaturations(); }
A default implementation of the parameters for the adapter class to convert material laws from effect...
This is a fluid state which allows to set the fluid saturations and takes all other quantities from a...
This material law takes a material law defined for effective saturations and converts it to a materia...
Definition EffToAbsLaw.hpp:71
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition EffToAbsLaw.hpp:84
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition EffToAbsLaw.hpp:96
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition EffToAbsLaw.hpp:100
static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fs)
The capillary pressure-saturation curves depending on absolute saturations.
Definition EffToAbsLaw.hpp:117
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition EffToAbsLaw.hpp:104
static Evaluation pcnw(const Params ¶ms, const FluidState &fs)
The capillary pressure-saturation curve.
Definition EffToAbsLaw.hpp:170
static Evaluation krn(const Params ¶ms, const FluidState &fs)
The relative permeability of the non-wetting phase.
Definition EffToAbsLaw.hpp:303
static Evaluation Sn(const Params ¶ms, const FluidState &fs)
Calculate non-wetting liquid phase saturation given that the rest of the fluid state has been initial...
Definition EffToAbsLaw.hpp:234
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition EffToAbsLaw.hpp:92
static Evaluation krw(const Params ¶ms, const FluidState &fs)
The relative permeability for the wetting phase.
Definition EffToAbsLaw.hpp:275
static std::enable_if<(Traits::numPhases >2), Evaluation >::type Sg(const Params ¶ms, const FluidState &fs)
Calculate gas phase saturation given that the rest of the fluid state has been initialized.
Definition EffToAbsLaw.hpp:258
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition EffToAbsLaw.hpp:88
static Evaluation absoluteSaturation(const Params ¶ms, const Evaluation &S, unsigned phaseIdx)
Convert an effective saturation to an absolute one.
Definition EffToAbsLaw.hpp:364
static const int numPhases
The number of fluid phases.
Definition EffToAbsLaw.hpp:80
static std::enable_if<(Traits::numPhases >2), Evaluation >::type krg(const Params ¶ms, const FluidState &fs)
The relative permability of the gas phase.
Definition EffToAbsLaw.hpp:334
static Evaluation Sw(const Params ¶ms, const FluidState &fs)
Calculate wetting liquid phase saturation given that the rest of the fluid state has been initialized...
Definition EffToAbsLaw.hpp:215
static Evaluation effectiveSaturation(const Params ¶ms, const Evaluation &S, unsigned phaseIdx)
Convert an absolute saturation to an effective one.
Definition EffToAbsLaw.hpp:357
static void relativePermeabilities(Container &values, const Params ¶ms, const FluidState &fs)
The relative permeability-saturation curves depending on absolute saturations.
Definition EffToAbsLaw.hpp:143
static void saturations(Container &values, const Params ¶ms, const FluidState &fs)
The saturation-capillary pressure curves.
Definition EffToAbsLaw.hpp:202
This is a fluid state which allows to set the fluid saturations and takes all other quantities from a...
Definition SaturationOverlayFluidState.hpp:44
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30