00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00027 #ifndef OPM_ECL_EPS_TWO_PHASE_LAW_HPP
00028 #define OPM_ECL_EPS_TWO_PHASE_LAW_HPP
00029
00030 #include "EclEpsTwoPhaseLawParams.hpp"
00031
00032 #include <opm/material/fluidstates/SaturationOverlayFluidState.hpp>
00033 #include <opm/common/ErrorMacros.hpp>
00034 #include <opm/common/Exceptions.hpp>
00035
00036 namespace Opm {
00048 template <class EffLawT,
00049 class ParamsT = EclEpsTwoPhaseLawParams<EffLawT> >
00050 class EclEpsTwoPhaseLaw : public EffLawT::Traits
00051 {
00052 typedef EffLawT EffLaw;
00053
00054 public:
00055 typedef typename EffLaw::Traits Traits;
00056 typedef ParamsT Params;
00057 typedef typename EffLaw::Scalar Scalar;
00058
00059 enum { wettingPhaseIdx = Traits::wettingPhaseIdx };
00060 enum { nonWettingPhaseIdx = Traits::nonWettingPhaseIdx };
00061
00063 static const int numPhases = EffLaw::numPhases;
00064 static_assert(numPhases == 2,
00065 "The endpoint scaling applies to the nested twophase laws, not to "
00066 "the threephase one!");
00067
00070 static const bool implementsTwoPhaseApi = true;
00071
00072 static_assert(EffLaw::implementsTwoPhaseApi,
00073 "The material laws put into EclEpsTwoPhaseLaw must implement the "
00074 "two-phase material law API!");
00075
00078 static const bool implementsTwoPhaseSatApi = true;
00079
00080 static_assert(EffLaw::implementsTwoPhaseSatApi,
00081 "The material laws put into EclEpsTwoPhaseLaw must implement the "
00082 "two-phase material law saturation API!");
00083
00086 static const bool isSaturationDependent = true;
00087
00090 static const bool isPressureDependent = false;
00091
00094 static const bool isTemperatureDependent = false;
00095
00098 static const bool isCompositionDependent = false;
00099
00110 template <class Container, class FluidState>
00111 static void capillaryPressures(Container& , const Params& , const FluidState& )
00112 {
00113 OPM_THROW(NotImplemented,
00114 "The capillaryPressures(fs) method is not yet implemented");
00115 }
00116
00127 template <class Container, class FluidState>
00128 static void relativePermeabilities(Container& , const Params& , const FluidState& )
00129 {
00130 OPM_THROW(NotImplemented,
00131 "The pcnw(fs) method is not yet implemented");
00132 }
00133
00145 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00146 static Evaluation pcnw(const Params& , const FluidState& )
00147 {
00148 OPM_THROW(NotImplemented,
00149 "The pcnw(fs) method is not yet implemented");
00150 }
00151
00152 template <class Evaluation>
00153 static Evaluation twoPhaseSatPcnw(const Params& params, const Evaluation& SwScaled)
00154 {
00155 const Evaluation& SwUnscaled = scaledToUnscaledSatPc(params, SwScaled);
00156 const Evaluation& pcUnscaled = EffLaw::twoPhaseSatPcnw(params.effectiveLawParams(), SwUnscaled);
00157 return unscaledToScaledPcnw_(params, pcUnscaled);
00158 }
00159
00160 template <class Evaluation>
00161 static Evaluation twoPhaseSatPcnwInv(const Params& params, const Evaluation& pcnwScaled)
00162 {
00163 Evaluation pcnwUnscaled = scaledToUnscaledPcnw_(params, pcnwScaled);
00164 Evaluation SwUnscaled = EffLaw::twoPhaseSatPcnwInv(params.effectiveLawParams(), pcnwUnscaled);
00165 return unscaledToScaledSatPc(params, SwUnscaled);
00166 }
00167
00171 template <class Container, class FluidState>
00172 static void saturations(Container& , const Params& , const FluidState& )
00173 {
00174 OPM_THROW(NotImplemented,
00175 "The saturations(fs) method is not yet implemented");
00176 }
00177
00182 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00183 static Evaluation Sw(const Params& , const FluidState& )
00184 {
00185 OPM_THROW(NotImplemented,
00186 "The Sw(fs) method is not yet implemented");
00187 }
00188
00189 template <class Evaluation>
00190 static Evaluation twoPhaseSatSw(const Params& , const Evaluation& )
00191 {
00192 OPM_THROW(NotImplemented,
00193 "The twoPhaseSatSw(pc) method is not yet implemented");
00194 }
00195
00200 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00201 static Evaluation Sn(const Params& , const FluidState& )
00202 {
00203 OPM_THROW(NotImplemented,
00204 "The Sn(pc) method is not yet implemented");
00205 }
00206
00207 template <class Evaluation>
00208 static Evaluation twoPhaseSatSn(const Params& , const Evaluation& )
00209 {
00210 OPM_THROW(NotImplemented,
00211 "The twoPhaseSatSn(pc) method is not yet implemented");
00212 }
00213
00223 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00224 static Evaluation krw(const Params& , const FluidState& )
00225 {
00226 OPM_THROW(NotImplemented,
00227 "The krw(fs) method is not yet implemented");
00228 }
00229
00230 template <class Evaluation>
00231 static Evaluation twoPhaseSatKrw(const Params& params, const Evaluation& SwScaled)
00232 {
00233 const Evaluation& SwUnscaled = scaledToUnscaledSatKrw(params, SwScaled);
00234 const Evaluation& krwUnscaled = EffLaw::twoPhaseSatKrw(params.effectiveLawParams(), SwUnscaled);
00235 return unscaledToScaledKrw_(params, krwUnscaled);
00236 }
00237
00238 template <class Evaluation>
00239 static Evaluation twoPhaseSatKrwInv(const Params& params, const Evaluation& krwScaled)
00240 {
00241 Evaluation krwUnscaled = scaledToUnscaledKrw_(params, krwScaled);
00242 Evaluation SwUnscaled = EffLaw::twoPhaseSatKrwInv(params.effectiveLawParams(), krwUnscaled);
00243 return unscaledToScaledSatKrw(params, SwUnscaled);
00244 }
00245
00249 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00250 static Evaluation krn(const Params& , const FluidState& )
00251 {
00252 OPM_THROW(NotImplemented,
00253 "The krn(fs) method is not yet implemented");
00254 }
00255
00256 template <class Evaluation>
00257 static Evaluation twoPhaseSatKrn(const Params& params, const Evaluation& SwScaled)
00258 {
00259 const Evaluation& SwUnscaled = scaledToUnscaledSatKrn(params, SwScaled);
00260 const Evaluation& krnUnscaled = EffLaw::twoPhaseSatKrn(params.effectiveLawParams(), SwUnscaled);
00261 return unscaledToScaledKrn_(params, krnUnscaled);
00262 }
00263
00264 template <class Evaluation>
00265 static Evaluation twoPhaseSatKrnInv(const Params& params, const Evaluation& krnScaled)
00266 {
00267 Evaluation krnUnscaled = scaledToUnscaledKrn_(params, krnScaled);
00268 Evaluation SwUnscaled = EffLaw::twoPhaseSatKrnInv(params.effectiveLawParams(), krnUnscaled);
00269 return unscaledToScaledSatKrn(params, SwUnscaled);
00270 }
00271
00277 template <class Evaluation>
00278 static Evaluation scaledToUnscaledSatPc(const Params& params, const Evaluation& SwScaled)
00279 {
00280 if (!params.config().enableSatScaling())
00281 return SwScaled;
00282
00283
00284
00285 return scaledToUnscaledSatTwoPoint_(SwScaled,
00286 params.unscaledPoints().saturationPcPoints(),
00287 params.scaledPoints().saturationPcPoints());
00288 }
00289
00290 template <class Evaluation>
00291 static Evaluation unscaledToScaledSatPc(const Params& params, const Evaluation& SwUnscaled)
00292 {
00293 if (!params.config().enableSatScaling())
00294 return SwUnscaled;
00295
00296
00297
00298 return unscaledToScaledSatTwoPoint_(SwUnscaled,
00299 params.unscaledPoints().saturationPcPoints(),
00300 params.scaledPoints().saturationPcPoints());
00301 }
00302
00307 template <class Evaluation>
00308 static Evaluation scaledToUnscaledSatKrw(const Params& params, const Evaluation& SwScaled)
00309 {
00310 if (!params.config().enableSatScaling())
00311 return SwScaled;
00312
00313 if (params.config().enableThreePointKrSatScaling()) {
00314 return scaledToUnscaledSatThreePoint_(SwScaled,
00315 params.unscaledPoints().saturationKrwPoints(),
00316 params.scaledPoints().saturationKrwPoints());
00317 }
00318 else {
00319 return scaledToUnscaledSatTwoPoint_(SwScaled,
00320 params.unscaledPoints().saturationKrwPoints(),
00321 params.scaledPoints().saturationKrwPoints());
00322 }
00323 }
00324
00325 template <class Evaluation>
00326 static Evaluation unscaledToScaledSatKrw(const Params& params, const Evaluation& SwUnscaled)
00327 {
00328 if (!params.config().enableSatScaling())
00329 return SwUnscaled;
00330
00331 if (params.config().enableThreePointKrSatScaling()) {
00332 return unscaledToScaledSatThreePoint_(SwUnscaled,
00333 params.unscaledPoints().saturationKrwPoints(),
00334 params.scaledPoints().saturationKrwPoints());
00335 }
00336 else {
00337 return unscaledToScaledSatTwoPoint_(SwUnscaled,
00338 params.unscaledPoints().saturationKrwPoints(),
00339 params.scaledPoints().saturationKrwPoints());
00340 }
00341 }
00342
00347 template <class Evaluation>
00348 static Evaluation scaledToUnscaledSatKrn(const Params& params, const Evaluation& SwScaled)
00349 {
00350 if (!params.config().enableSatScaling())
00351 return SwScaled;
00352
00353 if (params.config().enableThreePointKrSatScaling())
00354 return scaledToUnscaledSatThreePoint_(SwScaled,
00355 params.unscaledPoints().saturationKrnPoints(),
00356 params.scaledPoints().saturationKrnPoints());
00357 else
00358 return scaledToUnscaledSatTwoPoint_(SwScaled,
00359 params.unscaledPoints().saturationKrnPoints(),
00360 params.scaledPoints().saturationKrnPoints());
00361 }
00362
00363
00364 template <class Evaluation>
00365 static Evaluation unscaledToScaledSatKrn(const Params& params, const Evaluation& SwUnscaled)
00366 {
00367 if (!params.config().enableSatScaling())
00368 return SwUnscaled;
00369
00370 if (params.config().enableThreePointKrSatScaling()) {
00371 return unscaledToScaledSatThreePoint_(SwUnscaled,
00372 params.unscaledPoints().saturationKrnPoints(),
00373 params.scaledPoints().saturationKrnPoints());
00374 }
00375 else {
00376 return unscaledToScaledSatTwoPoint_(SwUnscaled,
00377 params.unscaledPoints().saturationKrnPoints(),
00378 params.scaledPoints().saturationKrnPoints());
00379 }
00380 }
00381
00382 private:
00383 template <class Evaluation, class PointsContainer>
00384 static Evaluation scaledToUnscaledSatTwoPoint_(const Evaluation& scaledSat,
00385 const PointsContainer& unscaledSats,
00386 const PointsContainer& scaledSats)
00387 {
00388 return
00389 unscaledSats[0]
00390 +
00391 (scaledSat - scaledSats[0])*((unscaledSats[1] - unscaledSats[0])/(scaledSats[1] - scaledSats[0]));
00392 }
00393
00394 template <class Evaluation, class PointsContainer>
00395 static Evaluation unscaledToScaledSatTwoPoint_(const Evaluation& unscaledSat,
00396 const PointsContainer& unscaledSats,
00397 const PointsContainer& scaledSats)
00398 {
00399 return
00400 scaledSats[0]
00401 +
00402 (unscaledSat - unscaledSats[0])*((scaledSats[1] - scaledSats[0])/(unscaledSats[1] - unscaledSats[0]));
00403 }
00404
00405 template <class Evaluation, class PointsContainer>
00406 static Evaluation scaledToUnscaledSatThreePoint_(const Evaluation& scaledSat,
00407 const PointsContainer& unscaledSats,
00408 const PointsContainer& scaledSats)
00409 {
00410 if (unscaledSats[1] >= unscaledSats[2])
00411 return scaledToUnscaledSatTwoPoint_(scaledSat, unscaledSats, scaledSats);
00412
00413 if (scaledSat < scaledSats[1]) {
00414 Scalar delta = scaledSats[1] - scaledSats[0];
00415 if (delta <= 1e-20)
00416 delta = 1.0;
00417
00418 return
00419 unscaledSats[0]
00420 +
00421 (scaledSat - scaledSats[0])*((unscaledSats[1] - unscaledSats[0])/delta);
00422 }
00423 else {
00424 Scalar delta = scaledSats[2] - scaledSats[1];
00425 if (delta <= 1e-20)
00426 delta = 1.0;
00427
00428 return
00429 unscaledSats[1]
00430 +
00431 (scaledSat - scaledSats[1])*((unscaledSats[2] - unscaledSats[1])/delta);
00432 }
00433 }
00434
00435 template <class Evaluation, class PointsContainer>
00436 static Evaluation unscaledToScaledSatThreePoint_(const Evaluation& unscaledSat,
00437 const PointsContainer& unscaledSats,
00438 const PointsContainer& scaledSats)
00439 {
00440 if (unscaledSats[1] >= unscaledSats[2])
00441 return unscaledToScaledSatTwoPoint_(unscaledSat, unscaledSats, scaledSats);
00442
00443 if (unscaledSat < unscaledSats[1]) {
00444 Scalar delta = unscaledSats[1] - unscaledSats[0];
00445 if (delta <= 1e-20)
00446 delta = 1.0;
00447
00448 return
00449 scaledSats[0]
00450 +
00451 (unscaledSat - unscaledSats[0])*((scaledSats[1] - scaledSats[0])/delta);
00452 }
00453 else {
00454 Scalar delta = unscaledSats[2] - unscaledSats[1];
00455 if (delta <= 1e-20)
00456 delta = 1.0;
00457
00458 return
00459 scaledSats[1]
00460 +
00461 (unscaledSat - unscaledSats[1])*((scaledSats[2] - scaledSats[1])/delta);
00462 }
00463 }
00464
00468 template <class Evaluation>
00469 static Evaluation unscaledToScaledPcnw_(const Params& params, const Evaluation& unscaledPcnw)
00470 {
00471 if (params.config().enableLeverettScaling()) {
00472 Scalar alpha = params.scaledPoints().leverettFactor();
00473 return unscaledPcnw*alpha;
00474 }
00475 else if (params.config().enablePcScaling()) {
00476 Scalar alpha = params.scaledPoints().maxPcnw()/params.unscaledPoints().maxPcnw();
00477 return unscaledPcnw*alpha;
00478 }
00479
00480 return unscaledPcnw;
00481 }
00482
00483 template <class Evaluation>
00484 static Evaluation scaledToUnscaledPcnw_(const Params& params, const Evaluation& scaledPcnw)
00485 {
00486 if (params.config().enableLeverettScaling()) {
00487 Scalar alpha = params.scaledPoints().leverettFactor();
00488 return scaledPcnw/alpha;
00489 }
00490 else if (params.config().enablePcScaling()) {
00491 Scalar alpha = params.unscaledPoints().maxPcnw()/params.scaledPoints().maxPcnw();
00492 return scaledPcnw/alpha;
00493 }
00494
00495 return scaledPcnw;
00496 }
00497
00501 template <class Evaluation>
00502 static Evaluation unscaledToScaledKrw_(const Params& params, const Evaluation& unscaledKrw)
00503 {
00504 if (!params.config().enableKrwScaling())
00505 return unscaledKrw;
00506
00507
00508 Scalar alpha = params.scaledPoints().maxKrw()/params.unscaledPoints().maxKrw();
00509 return unscaledKrw*alpha;
00510 }
00511
00512 template <class Evaluation>
00513 static Evaluation scaledToUnscaledKrw_(const Params& params, const Evaluation& scaledKrw)
00514 {
00515 if (!params.config().enableKrwScaling())
00516 return scaledKrw;
00517
00518 Scalar alpha = params.unscaledPoints().maxKrw()/params.scaledPoints().maxKrw();
00519 return scaledKrw*alpha;
00520 }
00521
00525 template <class Evaluation>
00526 static Evaluation unscaledToScaledKrn_(const Params& params, const Evaluation& unscaledKrn)
00527 {
00528 if (!params.config().enableKrnScaling())
00529 return unscaledKrn;
00530
00531
00532 Scalar alpha = params.scaledPoints().maxKrn()/params.unscaledPoints().maxKrn();
00533 return unscaledKrn*alpha;
00534 }
00535
00536 template <class Evaluation>
00537 static Evaluation scaledToUnscaledKrn_(const Params& params, const Evaluation& scaledKrn)
00538 {
00539 if (!params.config().enableKrnScaling())
00540 return scaledKrn;
00541
00542 Scalar alpha = params.unscaledPoints().maxKrn()/params.scaledPoints().maxKrn();
00543 return scaledKrn*alpha;
00544 }
00545 };
00546 }
00547
00548 #endif