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_TWO_PHASE_MATERIAL_HPP
00028 #define OPM_ECL_TWO_PHASE_MATERIAL_HPP
00029
00030 #include "EclTwoPhaseMaterialParams.hpp"
00031
00032 #include <opm/common/Valgrind.hpp>
00033 #include <opm/material/common/MathToolbox.hpp>
00034
00035 #include <opm/common/Exceptions.hpp>
00036 #include <opm/common/ErrorMacros.hpp>
00037
00038 #include <algorithm>
00039
00040 namespace Opm {
00041
00051 template <class TraitsT,
00052 class GasOilMaterialLawT,
00053 class OilWaterMaterialLawT,
00054 class ParamsT = EclTwoPhaseMaterialParams<TraitsT,
00055 typename GasOilMaterialLawT::Params,
00056 typename OilWaterMaterialLawT::Params> >
00057 class EclTwoPhaseMaterial : public TraitsT
00058 {
00059 public:
00060 typedef GasOilMaterialLawT GasOilMaterialLaw;
00061 typedef OilWaterMaterialLawT OilWaterMaterialLaw;
00062
00063
00064 static_assert(TraitsT::numPhases == 3,
00065 "The number of phases considered by this capillary pressure "
00066 "law is always three!");
00067 static_assert(GasOilMaterialLaw::numPhases == 2,
00068 "The number of phases considered by the gas-oil capillary "
00069 "pressure law must be two!");
00070 static_assert(OilWaterMaterialLaw::numPhases == 2,
00071 "The number of phases considered by the oil-water capillary "
00072 "pressure law must be two!");
00073 static_assert(std::is_same<typename GasOilMaterialLaw::Scalar,
00074 typename OilWaterMaterialLaw::Scalar>::value,
00075 "The two two-phase capillary pressure laws must use the same "
00076 "type of floating point values.");
00077
00078 typedef TraitsT Traits;
00079 typedef ParamsT Params;
00080 typedef typename Traits::Scalar Scalar;
00081
00082 static const int numPhases = 3;
00083 static const int waterPhaseIdx = Traits::wettingPhaseIdx;
00084 static const int oilPhaseIdx = Traits::nonWettingPhaseIdx;
00085 static const int gasPhaseIdx = Traits::gasPhaseIdx;
00086
00089 static const bool implementsTwoPhaseApi = false;
00090
00093 static const bool implementsTwoPhaseSatApi = false;
00094
00097 static const bool isSaturationDependent = true;
00098
00101 static const bool isPressureDependent = false;
00102
00105 static const bool isTemperatureDependent = false;
00106
00109 static const bool isCompositionDependent = false;
00110
00125 template <class ContainerT, class FluidState>
00126 static void capillaryPressures(ContainerT& values,
00127 const Params& params,
00128 const FluidState& fluidState)
00129 {
00130 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
00131
00132 switch (params.approach()) {
00133 case EclTwoPhaseGasOil: {
00134 const Evaluation& So =
00135 Opm::decay<Evaluation>(fluidState.saturation(oilPhaseIdx));
00136
00137 values[oilPhaseIdx] = 0.0;
00138 values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(), So);
00139 break;
00140 }
00141
00142 case EclTwoPhaseOilWater: {
00143 const Evaluation& Sw =
00144 Opm::decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
00145
00146 values[waterPhaseIdx] = 0.0;
00147 values[oilPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(), Sw);
00148 break;
00149 }
00150
00151 case EclTwoPhaseGasWater: {
00152 const Evaluation& Sw =
00153 Opm::decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
00154
00155 values[waterPhaseIdx] = 0.0;
00156 values[gasPhaseIdx] =
00157 OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(), Sw)
00158 + GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(), 0.0);
00159 break;
00160 }
00161
00162 }
00163 }
00164
00165
00166
00167
00168
00169
00170
00171 static void oilWaterHysteresisParams(Scalar& pcSwMdc,
00172 Scalar& krnSwMdc,
00173 const Params& params)
00174 {
00175 pcSwMdc = params.oilWaterParams().pcSwMdc();
00176 krnSwMdc = params.oilWaterParams().krnSwMdc();
00177
00178 Valgrind::CheckDefined(pcSwMdc);
00179 Valgrind::CheckDefined(krnSwMdc);
00180 }
00181
00182
00183
00184
00185
00186
00187
00188 static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
00189 const Scalar& krnSwMdc,
00190 Params& params)
00191 {
00192 const Scalar krwSw = 2.0;
00193 params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
00194 }
00195
00196
00197
00198
00199
00200
00201
00202 static void gasOilHysteresisParams(Scalar& pcSwMdc,
00203 Scalar& krnSwMdc,
00204 const Params& params)
00205 {
00206 pcSwMdc = params.gasOilParams().pcSwMdc();
00207 krnSwMdc = params.gasOilParams().krnSwMdc();
00208
00209 Valgrind::CheckDefined(pcSwMdc);
00210 Valgrind::CheckDefined(krnSwMdc);
00211 }
00212
00213
00214
00215
00216
00217
00218
00219 static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
00220 const Scalar& krnSwMdc,
00221 Params& params)
00222 {
00223 const Scalar krwSw = 2.0;
00224 params.gasOilParams().update(pcSwMdc, krwSw, krnSwMdc);
00225 }
00226
00236 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00237 static Evaluation pcgn(const Params& ,
00238 const FluidState& )
00239 {
00240 OPM_THROW(std::logic_error, "Not implemented: pcgn()");
00241 }
00242
00252 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00253 static Evaluation pcnw(const Params& ,
00254 const FluidState& )
00255 {
00256 OPM_THROW(std::logic_error, "Not implemented: pcnw()");
00257 }
00258
00262 template <class ContainerT, class FluidState>
00263 static void saturations(ContainerT& ,
00264 const Params& ,
00265 const FluidState& )
00266 {
00267 OPM_THROW(std::logic_error, "Not implemented: saturations()");
00268 }
00269
00273 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00274 static Evaluation Sg(const Params& ,
00275 const FluidState& )
00276 {
00277 OPM_THROW(std::logic_error, "Not implemented: Sg()");
00278 }
00279
00283 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00284 static Evaluation Sn(const Params& ,
00285 const FluidState& )
00286 {
00287 OPM_THROW(std::logic_error, "Not implemented: Sn()");
00288 }
00289
00293 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00294 static Evaluation Sw(const Params& ,
00295 const FluidState& )
00296 {
00297 OPM_THROW(std::logic_error, "Not implemented: Sw()");
00298 }
00299
00315 template <class ContainerT, class FluidState>
00316 static void relativePermeabilities(ContainerT& values,
00317 const Params& params,
00318 const FluidState& fluidState)
00319 {
00320 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
00321
00322 switch (params.approach()) {
00323 case EclTwoPhaseGasOil: {
00324 const Evaluation& So =
00325 Opm::decay<Evaluation>(fluidState.saturation(oilPhaseIdx));
00326
00327 values[oilPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), So);
00328 values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), So);
00329 break;
00330 }
00331
00332 case EclTwoPhaseOilWater: {
00333 const Evaluation& Sw =
00334 Opm::decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
00335
00336 values[waterPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), Sw);
00337 values[oilPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), Sw);
00338 break;
00339 }
00340
00341 case EclTwoPhaseGasWater: {
00342 const Evaluation& Sw =
00343 Opm::decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
00344
00345 values[waterPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), Sw);
00346 values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), Sw);
00347 break;
00348 }
00349 }
00350 }
00351
00355 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00356 static Evaluation krg(const Params& ,
00357 const FluidState& )
00358 {
00359 OPM_THROW(std::logic_error, "Not implemented: krg()");
00360 }
00361
00365 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00366 static Evaluation krw(const Params& ,
00367 const FluidState& )
00368 {
00369 OPM_THROW(std::logic_error, "Not implemented: krw()");
00370 }
00371
00375 template <class FluidState, class Evaluation = typename FluidState::Scalar>
00376 static Evaluation krn(const Params& ,
00377 const FluidState& )
00378 {
00379 OPM_THROW(std::logic_error, "Not implemented: krn()");
00380 }
00381
00382
00390 template <class FluidState>
00391 static void updateHysteresis(Params& params, const FluidState& fluidState)
00392 {
00393 switch (params.approach()) {
00394 case EclTwoPhaseGasOil: {
00395 Scalar So = Opm::scalarValue(fluidState.saturation(oilPhaseIdx));
00396
00397 params.gasOilParams().update(So, So, So);
00398 break;
00399 }
00400
00401 case EclTwoPhaseOilWater: {
00402 Scalar Sw = Opm::scalarValue(fluidState.saturation(waterPhaseIdx));
00403
00404 params.oilWaterParams().update(Sw, Sw, Sw);
00405 break;
00406 }
00407
00408 case EclTwoPhaseGasWater: {
00409 Scalar Sw = Opm::scalarValue(fluidState.saturation(waterPhaseIdx));
00410
00411 params.oilWaterParams().update(Sw, Sw, 0);
00412 params.gasOilParams().update(1.0, 0.0, Sw);
00413 break;
00414 }
00415 }
00416 }
00417 };
00418 }
00419
00420 #endif