25 #ifndef OPM_ECL_EPS_TWO_PHASE_LAW_HPP
26 #define OPM_ECL_EPS_TWO_PHASE_LAW_HPP
31 #include <opm/common/ErrorMacros.hpp>
32 #include <opm/common/Exceptions.hpp>
46 template <
class EffLawT,
47 class ParamsT = EclEpsTwoPhaseLawParams<EffLawT> >
50 typedef EffLawT EffLaw;
53 typedef typename EffLaw::Traits
Traits;
55 typedef typename EffLaw::Scalar
Scalar;
62 static_assert(numPhases == 2,
63 "The endpoint scaling applies to the nested twophase laws, not to "
64 "the threephase one!");
68 static const bool implementsTwoPhaseApi =
true;
70 static_assert(EffLaw::implementsTwoPhaseApi,
71 "The material laws put into EclEpsTwoPhaseLaw must implement the "
72 "two-phase material law API!");
76 static const bool implementsTwoPhaseSatApi =
true;
78 static_assert(EffLaw::implementsTwoPhaseSatApi,
79 "The material laws put into EclEpsTwoPhaseLaw must implement the "
80 "two-phase material law saturation API!");
84 static const bool isSaturationDependent =
true;
88 static const bool isPressureDependent =
false;
92 static const bool isTemperatureDependent =
false;
96 static const bool isCompositionDependent =
false;
108 template <
class Container,
class Flu
idState>
109 static void capillaryPressures(Container& ,
const Params& ,
const FluidState& )
111 OPM_THROW(NotImplemented,
112 "The capillaryPressures(fs) method is not yet implemented");
125 template <
class Container,
class Flu
idState>
126 static void relativePermeabilities(Container& ,
const Params& ,
const FluidState& )
128 OPM_THROW(NotImplemented,
129 "The pcnw(fs) method is not yet implemented");
143 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
144 static Evaluation pcnw(
const Params& ,
const FluidState& )
146 OPM_THROW(NotImplemented,
147 "The pcnw(fs) method is not yet implemented");
150 template <
class Evaluation>
151 static Evaluation twoPhaseSatPcnw(
const Params ¶ms,
const Evaluation& SwScaled)
153 const Evaluation& SwUnscaled = scaledToUnscaledSatPc(params, SwScaled);
154 const Evaluation& pcUnscaled = EffLaw::twoPhaseSatPcnw(params.effectiveLawParams(), SwUnscaled);
155 return unscaledToScaledPcnw_(params, pcUnscaled);
158 template <
class Evaluation>
159 static Evaluation twoPhaseSatPcnwInv(
const Params ¶ms,
const Evaluation& pcnwScaled)
161 Evaluation pcnwUnscaled = scaledToUnscaledPcnw_(params, pcnwScaled);
162 Evaluation SwUnscaled = EffLaw::twoPhaseSatPcnwInv(params.effectiveLawParams(), pcnwUnscaled);
163 return unscaledToScaledSatPc(params, SwUnscaled);
169 template <
class Container,
class Flu
idState>
170 static void saturations(Container& ,
const Params& ,
const FluidState& )
172 OPM_THROW(NotImplemented,
173 "The saturations(fs) method is not yet implemented");
180 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
181 static Evaluation Sw(
const Params& ,
const FluidState& )
183 OPM_THROW(NotImplemented,
184 "The Sw(fs) method is not yet implemented");
187 template <
class Evaluation>
188 static Evaluation twoPhaseSatSw(
const Params& ,
const Evaluation& )
190 OPM_THROW(NotImplemented,
191 "The twoPhaseSatSw(pc) method is not yet implemented");
198 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
199 static Evaluation Sn(
const Params& ,
const FluidState& )
201 OPM_THROW(NotImplemented,
202 "The Sn(pc) method is not yet implemented");
205 template <
class Evaluation>
206 static Evaluation twoPhaseSatSn(
const Params& ,
const Evaluation& )
208 OPM_THROW(NotImplemented,
209 "The twoPhaseSatSn(pc) method is not yet implemented");
221 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
222 static Evaluation krw(
const Params& ,
const FluidState& )
224 OPM_THROW(NotImplemented,
225 "The krw(fs) method is not yet implemented");
228 template <
class Evaluation>
229 static Evaluation twoPhaseSatKrw(
const Params ¶ms,
const Evaluation& SwScaled)
231 const Evaluation& SwUnscaled = scaledToUnscaledSatKrw(params, SwScaled);
232 const Evaluation& krwUnscaled = EffLaw::twoPhaseSatKrw(params.effectiveLawParams(), SwUnscaled);
233 return unscaledToScaledKrw_(params, krwUnscaled);
236 template <
class Evaluation>
237 static Evaluation twoPhaseSatKrwInv(
const Params ¶ms,
const Evaluation& krwScaled)
239 Evaluation krwUnscaled = scaledToUnscaledKrw_(params, krwScaled);
240 Evaluation SwUnscaled = EffLaw::twoPhaseSatKrwInv(params.effectiveLawParams(), krwUnscaled);
241 return unscaledToScaledSatKrw(params, SwUnscaled);
247 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
248 static Evaluation krn(
const Params& ,
const FluidState& )
250 OPM_THROW(NotImplemented,
251 "The krn(fs) method is not yet implemented");
254 template <
class Evaluation>
255 static Evaluation twoPhaseSatKrn(
const Params ¶ms,
const Evaluation& SwScaled)
257 const Evaluation& SwUnscaled = scaledToUnscaledSatKrn(params, SwScaled);
258 const Evaluation& krnUnscaled = EffLaw::twoPhaseSatKrn(params.effectiveLawParams(), SwUnscaled);
259 return unscaledToScaledKrn_(params, krnUnscaled);
262 template <
class Evaluation>
263 static Evaluation twoPhaseSatKrnInv(
const Params ¶ms,
const Evaluation& krnScaled)
265 Evaluation krnUnscaled = scaledToUnscaledKrn_(params, krnScaled);
266 Evaluation SwUnscaled = EffLaw::twoPhaseSatKrnInv(params.effectiveLawParams(), krnUnscaled);
267 return unscaledToScaledSatKrn(params, SwUnscaled);
275 template <
class Evaluation>
276 static Evaluation scaledToUnscaledSatPc(
const Params ¶ms,
const Evaluation& SwScaled)
278 if (!params.config().enableSatScaling())
283 return scaledToUnscaledSatTwoPoint_(SwScaled,
284 params.unscaledPoints().saturationPcPoints(),
285 params.scaledPoints().saturationPcPoints());
288 template <
class Evaluation>
289 static Evaluation unscaledToScaledSatPc(
const Params ¶ms,
const Evaluation& SwUnscaled)
291 if (!params.config().enableSatScaling())
296 return unscaledToScaledSatTwoPoint_(SwUnscaled,
297 params.unscaledPoints().saturationPcPoints(),
298 params.scaledPoints().saturationPcPoints());
305 template <
class Evaluation>
306 static Evaluation scaledToUnscaledSatKrw(
const Params ¶ms,
const Evaluation& SwScaled)
308 if (!params.config().enableSatScaling())
311 if (params.config().enableThreePointKrSatScaling()) {
312 return scaledToUnscaledSatThreePoint_(SwScaled,
313 params.unscaledPoints().saturationKrwPoints(),
314 params.scaledPoints().saturationKrwPoints());
317 return scaledToUnscaledSatTwoPoint_(SwScaled,
318 params.unscaledPoints().saturationKrwPoints(),
319 params.scaledPoints().saturationKrwPoints());
323 template <
class Evaluation>
324 static Evaluation unscaledToScaledSatKrw(
const Params ¶ms,
const Evaluation& SwUnscaled)
326 if (!params.config().enableSatScaling())
329 if (params.config().enableThreePointKrSatScaling()) {
330 return unscaledToScaledSatThreePoint_(SwUnscaled,
331 params.unscaledPoints().saturationKrwPoints(),
332 params.scaledPoints().saturationKrwPoints());
335 return unscaledToScaledSatTwoPoint_(SwUnscaled,
336 params.unscaledPoints().saturationKrwPoints(),
337 params.scaledPoints().saturationKrwPoints());
345 template <
class Evaluation>
346 static Evaluation scaledToUnscaledSatKrn(
const Params ¶ms,
const Evaluation& SwScaled)
348 if (!params.config().enableSatScaling())
351 if (params.config().enableThreePointKrSatScaling())
352 return scaledToUnscaledSatThreePoint_(SwScaled,
353 params.unscaledPoints().saturationKrnPoints(),
354 params.scaledPoints().saturationKrnPoints());
356 return scaledToUnscaledSatTwoPoint_(SwScaled,
357 params.unscaledPoints().saturationKrnPoints(),
358 params.scaledPoints().saturationKrnPoints());
362 template <
class Evaluation>
363 static Evaluation unscaledToScaledSatKrn(
const Params ¶ms,
const Evaluation& SwUnscaled)
365 if (!params.config().enableSatScaling())
368 if (params.config().enableThreePointKrSatScaling()) {
369 return unscaledToScaledSatThreePoint_(SwUnscaled,
370 params.unscaledPoints().saturationKrnPoints(),
371 params.scaledPoints().saturationKrnPoints());
374 return unscaledToScaledSatTwoPoint_(SwUnscaled,
375 params.unscaledPoints().saturationKrnPoints(),
376 params.scaledPoints().saturationKrnPoints());
381 template <
class Evaluation,
class Po
intsContainer>
382 static Evaluation scaledToUnscaledSatTwoPoint_(
const Evaluation& scaledSat,
383 const PointsContainer& unscaledSats,
384 const PointsContainer& scaledSats)
389 (scaledSat - scaledSats[0])*((unscaledSats[1] - unscaledSats[0])/(scaledSats[1] - scaledSats[0]));
392 template <
class Evaluation,
class Po
intsContainer>
393 static Evaluation unscaledToScaledSatTwoPoint_(
const Evaluation& unscaledSat,
394 const PointsContainer& unscaledSats,
395 const PointsContainer& scaledSats)
400 (unscaledSat - unscaledSats[0])*((scaledSats[1] - scaledSats[0])/(unscaledSats[1] - unscaledSats[0]));
403 template <
class Evaluation,
class Po
intsContainer>
404 static Evaluation scaledToUnscaledSatThreePoint_(
const Evaluation& scaledSat,
405 const PointsContainer& unscaledSats,
406 const PointsContainer& scaledSats)
408 if (unscaledSats[1] >= unscaledSats[2])
409 return scaledToUnscaledSatTwoPoint_(scaledSat, unscaledSats, scaledSats);
411 if (scaledSat < scaledSats[1])
415 (scaledSat - scaledSats[0])*((unscaledSats[1] - unscaledSats[0])/(scaledSats[1] - scaledSats[0]));
420 (scaledSat - scaledSats[1])*((unscaledSats[2] - unscaledSats[1])/(scaledSats[2] - scaledSats[1]));
423 template <
class Evaluation,
class Po
intsContainer>
424 static Evaluation unscaledToScaledSatThreePoint_(
const Evaluation& unscaledSat,
425 const PointsContainer& unscaledSats,
426 const PointsContainer& scaledSats)
428 if (unscaledSats[1] >= unscaledSats[2])
429 return unscaledToScaledSatTwoPoint_(unscaledSat, unscaledSats, scaledSats);
431 if (unscaledSat < unscaledSats[1])
435 (unscaledSat - unscaledSats[0])*((scaledSats[1] - scaledSats[0])/(unscaledSats[1] - unscaledSats[0]));
440 (unscaledSat - unscaledSats[1])*((scaledSats[2] - scaledSats[1])/(unscaledSats[2] - unscaledSats[1]));
446 template <
class Evaluation>
447 static Evaluation unscaledToScaledPcnw_(
const Params ¶ms,
const Evaluation& unscaledPcnw)
449 if (!params.config().enablePcScaling())
452 Scalar alpha = params.scaledPoints().maxPcnw()/params.unscaledPoints().maxPcnw();
453 return unscaledPcnw*alpha;
456 template <
class Evaluation>
457 static Evaluation scaledToUnscaledPcnw_(
const Params ¶ms,
const Evaluation& scaledPcnw)
459 if (!params.config().enablePcScaling())
462 Scalar alpha = params.unscaledPoints().maxPcnw()/params.scaledPoints().maxPcnw();
463 return scaledPcnw/alpha;
469 template <
class Evaluation>
470 static Evaluation unscaledToScaledKrw_(
const Params ¶ms,
const Evaluation& unscaledKrw)
472 if (!params.config().enableKrwScaling())
476 Scalar alpha = params.scaledPoints().maxKrw()/params.unscaledPoints().maxKrw();
477 return unscaledKrw*alpha;
480 template <
class Evaluation>
481 static Evaluation scaledToUnscaledKrw_(
const Params ¶ms,
const Evaluation& scaledKrw)
483 if (!params.config().enableKrwScaling())
486 Scalar alpha = params.unscaledPoints().maxKrw()/params.scaledPoints().maxKrw();
487 return scaledKrw*alpha;
493 template <
class Evaluation>
494 static Evaluation unscaledToScaledKrn_(
const Params ¶ms,
const Evaluation& unscaledKrn)
496 if (!params.config().enableKrnScaling())
500 Scalar alpha = params.scaledPoints().maxKrn()/params.unscaledPoints().maxKrn();
501 return unscaledKrn*alpha;
504 template <
class Evaluation>
505 static Evaluation scaledToUnscaledKrn_(
const Params ¶ms,
const Evaluation& scaledKrn)
507 if (!params.config().enableKrnScaling())
510 Scalar alpha = params.unscaledPoints().maxKrn()/params.scaledPoints().maxKrn();
511 return scaledKrn*alpha;
Definition: Air_Mesitylene.hpp:31
Definition: EclEpsTwoPhaseLaw.hpp:58
Definition: EclEpsTwoPhaseLaw.hpp:57
A default implementation of the parameters for the material law adapter class which implements ECL en...
ParamsT Params
Definition: EclEpsTwoPhaseLaw.hpp:54
EffLaw::Scalar Scalar
Definition: EclEpsTwoPhaseLaw.hpp:55
EffLaw::Traits Traits
Definition: EclEpsTwoPhaseLaw.hpp:53
This material law takes a material law defined for unscaled saturation and converts it to a material ...
Definition: EclEpsTwoPhaseLaw.hpp:48
This is a fluid state which allows to set the fluid saturations and takes all other quantities from a...
static const int numPhases
The number of fluid phases.
Definition: EclEpsTwoPhaseLaw.hpp:61