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_EFF_TO_ABS_LAW_HPP
00028 #define OPM_EFF_TO_ABS_LAW_HPP
00029
00030 #include "EffToAbsLawParams.hpp"
00031
00032 #include <opm/material/fluidstates/SaturationOverlayFluidState.hpp>
00033
00034 namespace Opm {
00069 template <class EffLawT, class ParamsT = EffToAbsLawParams<typename EffLawT::Params, EffLawT::numPhases> >
00070 class EffToAbsLaw : public EffLawT::Traits
00071 {
00072 typedef EffLawT EffLaw;
00073
00074 public:
00075 typedef typename EffLaw::Traits Traits;
00076 typedef ParamsT Params;
00077 typedef typename EffLaw::Scalar Scalar;
00078
00080 static const int numPhases = EffLaw::numPhases;
00081
00084 static const bool implementsTwoPhaseApi = EffLaw::implementsTwoPhaseApi;
00085
00088 static const bool implementsTwoPhaseSatApi = EffLaw::implementsTwoPhaseSatApi;
00089
00092 static const bool isSaturationDependent = EffLaw::isSaturationDependent;
00093
00096 static const bool isPressureDependent = EffLaw::isPressureDependent;
00097
00100 static const bool isTemperatureDependent = EffLaw::isTemperatureDependent;
00101
00104 static const bool isCompositionDependent = EffLaw::isCompositionDependent;
00105
00116 template <class Container, class FluidState>
00117 static void capillaryPressures(Container& values, const Params& params, const FluidState& fs)
00118 {
00119 typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
00120
00121 OverlayFluidState overlayFs(fs);
00122 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
00123 overlayFs.setSaturation(phaseIdx,
00124 effectiveSaturation(params,
00125 fs.saturation(phaseIdx),
00126 phaseIdx));
00127 }
00128
00129 EffLaw::template capillaryPressures<Container, OverlayFluidState>(values, params, overlayFs);
00130 }
00131
00142 template <class Container, class FluidState>
00143 static void relativePermeabilities(Container& values, const Params& params, const FluidState& fs)
00144 {
00145 typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
00146
00147 OverlayFluidState overlayFs(fs);
00148 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
00149 overlayFs.setSaturation(phaseIdx,
00150 effectiveSaturation(params,
00151 fs.saturation(phaseIdx),
00152 phaseIdx));
00153 }
00154
00155 EffLaw::template relativePermeabilities<Container, OverlayFluidState>(values, params, overlayFs);
00156 }
00157
00169 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00170 static Evaluation pcnw(const Params& params, const FluidState& fs)
00171 {
00172 typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
00173
00174 static_assert(FluidState::numPhases == numPhases,
00175 "The fluid state and the material law must exhibit the same "
00176 "number of phases!");
00177
00178 OverlayFluidState overlayFs(fs);
00179 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
00180 overlayFs.setSaturation(phaseIdx,
00181 effectiveSaturation(params,
00182 fs.saturation(phaseIdx),
00183 phaseIdx));
00184 }
00185
00186 return EffLaw::template pcnw<OverlayFluidState, Evaluation>(params, overlayFs);
00187 }
00188
00189 template <class Evaluation>
00190 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
00191 twoPhaseSatPcnw(const Params& params, const Evaluation& SwAbs)
00192 {
00193 const Evaluation& SwEff = effectiveSaturation(params, SwAbs, Traits::wettingPhaseIdx);
00194
00195 return EffLaw::twoPhaseSatPcnw(params, SwEff);
00196 }
00197
00201 template <class Container, class FluidState>
00202 static void saturations(Container& values, const Params& params, const FluidState& fs)
00203 {
00204 EffLaw::template saturations<Container, FluidState>(values, params, fs);
00205 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
00206 values[phaseIdx] = absoluteSaturation(params, values[phaseIdx], phaseIdx);
00207 }
00208 }
00209
00214 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00215 static Evaluation Sw(const Params& params, const FluidState& fs)
00216 {
00217 return absoluteSaturation(params,
00218 EffLaw::template Sw<FluidState, Evaluation>(params, fs),
00219 Traits::wettingPhaseIdx);
00220 }
00221
00222 template <class Evaluation>
00223 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
00224 twoPhaseSatSw(const Params& params, const Evaluation& Sw)
00225 { return absoluteSaturation(params,
00226 EffLaw::twoPhaseSatSw(params, Sw),
00227 Traits::wettingPhaseIdx); }
00228
00233 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00234 static Evaluation Sn(const Params& params, const FluidState& fs)
00235 {
00236 return absoluteSaturation(params,
00237 EffLaw::template Sn<FluidState, Evaluation>(params, fs),
00238 Traits::nonWettingPhaseIdx);
00239 }
00240
00241 template <class Evaluation>
00242 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
00243 twoPhaseSatSn(const Params& params, const Evaluation& Sw)
00244 {
00245 return absoluteSaturation(params,
00246 EffLaw::twoPhaseSatSn(params, Sw),
00247 Traits::nonWettingPhaseIdx);
00248 }
00249
00256 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00257 static typename std::enable_if< (Traits::numPhases > 2), Evaluation>::type
00258 Sg(const Params& params, const FluidState& fs)
00259 {
00260 return absoluteSaturation(params,
00261 EffLaw::template Sg<FluidState, Evaluation>(params, fs),
00262 Traits::gasPhaseIdx);
00263 }
00264
00274 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00275 static Evaluation krw(const Params& params, const FluidState& fs)
00276 {
00277 typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
00278
00279 static_assert(FluidState::numPhases == numPhases,
00280 "The fluid state and the material law must exhibit the same "
00281 "number of phases!");
00282
00283 OverlayFluidState overlayFs(fs);
00284 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
00285 overlayFs.setSaturation(phaseIdx,
00286 effectiveSaturation(params,
00287 fs.saturation(phaseIdx),
00288 phaseIdx));
00289 }
00290
00291 return EffLaw::template krw<OverlayFluidState, Evaluation>(params, overlayFs);
00292 }
00293
00294 template <class Evaluation>
00295 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
00296 twoPhaseSatKrw(const Params& params, const Evaluation& Sw)
00297 { return EffLaw::twoPhaseSatKrw(params, effectiveSaturation(params, Sw, Traits::nonWettingPhaseIdx)); }
00298
00302 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00303 static Evaluation krn(const Params& params, const FluidState& fs)
00304 {
00305 typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
00306
00307 static_assert(FluidState::numPhases == numPhases,
00308 "The fluid state and the material law must exhibit the same "
00309 "number of phases!");
00310
00311 OverlayFluidState overlayFs(fs);
00312 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
00313 overlayFs.setSaturation(phaseIdx,
00314 effectiveSaturation(params,
00315 fs.saturation(phaseIdx),
00316 phaseIdx));
00317 }
00318
00319 return EffLaw::template krn<OverlayFluidState, Evaluation>(params, overlayFs);
00320 }
00321
00322 template <class Evaluation>
00323 static typename std::enable_if<implementsTwoPhaseSatApi, Evaluation>::type
00324 twoPhaseSatKrn(const Params& params, const Evaluation& Sw)
00325 { return EffLaw::twoPhaseSatKrn(params, effectiveSaturation(params, Sw, Traits::nonWettingPhaseIdx)); }
00326
00332 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00333 static typename std::enable_if< (Traits::numPhases > 2), Evaluation>::type
00334 krg(const Params& params, const FluidState& fs)
00335 {
00336 typedef Opm::SaturationOverlayFluidState<FluidState> OverlayFluidState;
00337
00338 static_assert(FluidState::numPhases == numPhases,
00339 "The fluid state and the material law must exhibit the same "
00340 "number of phases!");
00341
00342 OverlayFluidState overlayFs(fs);
00343 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
00344 overlayFs.setSaturation(phaseIdx,
00345 effectiveSaturation(params,
00346 fs.saturation(phaseIdx),
00347 phaseIdx));
00348 }
00349
00350 return EffLaw::template krg<OverlayFluidState, Evaluation>(params, overlayFs);
00351 }
00352
00356 template <class Evaluation>
00357 static Evaluation effectiveSaturation(const Params& params, const Evaluation& S, unsigned phaseIdx)
00358 { return (S - params.residualSaturation(phaseIdx))/(1.0 - params.sumResidualSaturations()); }
00359
00363 template <class Evaluation>
00364 static Evaluation absoluteSaturation(const Params& params, const Evaluation& S, unsigned phaseIdx)
00365 { return S*(1.0 - params.sumResidualSaturations()) + params.residualSaturation(phaseIdx); }
00366
00367 private:
00376 static Scalar dSeff_dSabs_(const Params& params, int )
00377 { return 1.0/(1 - params.sumResidualSaturations()); }
00378
00387 static Scalar dSabs_dSeff_(const Params& params, int )
00388 { return 1 - params.sumResidualSaturations(); }
00389 };
00390 }
00391
00392 #endif