My Project
Loading...
Searching...
No Matches
EclTwoPhaseMaterial.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_ECL_TWO_PHASE_MATERIAL_HPP
28#define OPM_ECL_TWO_PHASE_MATERIAL_HPP
29
31
32#include <opm/common/TimingMacros.hpp>
33
36
37namespace Opm {
38
48template <class TraitsT,
49 class GasOilMaterialLawT,
50 class OilWaterMaterialLawT,
51 class GasWaterMaterialLawT,
52 class ParamsT = EclTwoPhaseMaterialParams<TraitsT,
53 typename GasOilMaterialLawT::Params,
54 typename OilWaterMaterialLawT::Params,
55 typename GasWaterMaterialLawT::Params> >
56class EclTwoPhaseMaterial : public TraitsT
57{
58public:
59 using GasOilMaterialLaw = GasOilMaterialLawT;
60 using OilWaterMaterialLaw = OilWaterMaterialLawT;
61 using GasWaterMaterialLaw = GasWaterMaterialLawT;
62
63 // some safety checks
64 static_assert(TraitsT::numPhases == 3,
65 "The number of phases considered by this capillary pressure "
66 "law is always three!");
67 static_assert(GasOilMaterialLaw::numPhases == 2,
68 "The number of phases considered by the gas-oil capillary "
69 "pressure law must be two!");
70 static_assert(OilWaterMaterialLaw::numPhases == 2,
71 "The number of phases considered by the oil-water capillary "
72 "pressure law must be two!");
73 static_assert(GasWaterMaterialLaw::numPhases == 2,
74 "The number of phases considered by the gas-water capillary "
75 "pressure law must be two!");
76 static_assert(std::is_same<typename GasOilMaterialLaw::Scalar,
77 typename OilWaterMaterialLaw::Scalar>::value,
78 "The two two-phase capillary pressure laws must use the same "
79 "type of floating point values.");
80
81 using Traits = TraitsT;
82 using Params = ParamsT;
83 using Scalar = typename Traits::Scalar;
84
85 static constexpr int numPhases = 3;
86 static constexpr int waterPhaseIdx = Traits::wettingPhaseIdx;
87 static constexpr int oilPhaseIdx = Traits::nonWettingPhaseIdx;
88 static constexpr int gasPhaseIdx = Traits::gasPhaseIdx;
89
92 static constexpr bool implementsTwoPhaseApi = false;
93
96 static constexpr bool implementsTwoPhaseSatApi = false;
97
100 static constexpr bool isSaturationDependent = true;
101
104 static constexpr bool isPressureDependent = false;
105
108 static constexpr bool isTemperatureDependent = false;
109
112 static constexpr bool isCompositionDependent = false;
113
114 template <class ContainerT, class FluidState>
115 static Scalar relpermOilInOilGasSystem(const Params& /*params*/,
116 const FluidState& /*fluidState*/) {
117 throw std::logic_error {
118 "relpermOilInOilGasSystem() is specific to three phases"
119 };
120 }
121 template <class ContainerT, class FluidState>
122 static Scalar relpermOilInOilWaterSystem(const Params& /*params*/,
123 const FluidState& /*fluidState*/) {
124 throw std::logic_error {
125 "relpermOilInOilWaterSystem() is specific to three phases"
126 };
127 }
128
144 template <class ContainerT, class FluidState>
145 static void capillaryPressures(ContainerT& values,
146 const Params& params,
147 const FluidState& fluidState)
148 {
149 OPM_TIMEFUNCTION_LOCAL();
150 using Evaluation = typename std::remove_reference<decltype(values[0])>::type;
151
152 switch (params.approach()) {
153 case EclTwoPhaseApproach::GasOil: {
154 const Evaluation& So =
155 decay<Evaluation>(fluidState.saturation(oilPhaseIdx));
156
157 values[oilPhaseIdx] = 0.0;
158 values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatPcnw(params.gasOilParams(), So);
159 break;
160 }
161
162 case EclTwoPhaseApproach::OilWater: {
163 const Evaluation& Sw =
164 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
165
166 values[waterPhaseIdx] = 0.0;
167 values[oilPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatPcnw(params.oilWaterParams(), Sw);
168 break;
169 }
170
171 case EclTwoPhaseApproach::GasWater: {
172 const Evaluation& Sw =
173 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
174
175 values[waterPhaseIdx] = 0.0;
176 values[gasPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatPcnw(params.gasWaterParams(), Sw);
177 break;
178 }
179
180 }
181 }
182
183 /*
184 * Hysteresis parameters for oil-water
185 * @see EclHysteresisTwoPhaseLawParams::soMax(...)
186 * @see EclHysteresisTwoPhaseLawParams::swMax(...)
187 * @see EclHysteresisTwoPhaseLawParams::swMin(...)
188 * \param params Parameters
189 */
190 static void oilWaterHysteresisParams(Scalar& soMax,
191 Scalar& swMax,
192 Scalar& swMin,
193 const Params& params)
194 {
195 soMax = 1.0 - params.oilWaterParams().krnSwMdc();
196 swMax = params.oilWaterParams().krwSwMdc();
197 swMin = params.oilWaterParams().pcSwMdc();
198 Valgrind::CheckDefined(soMax);
199 Valgrind::CheckDefined(swMax);
200 Valgrind::CheckDefined(swMin);
201 }
202
203 /*
204 * Hysteresis parameters for oil-water
205 * @see EclHysteresisTwoPhaseLawParams::soMax(...)
206 * @see EclHysteresisTwoPhaseLawParams::swMax(...)
207 * @see EclHysteresisTwoPhaseLawParams::swMin(...)
208 * \param params Parameters
209 */
210 static void setOilWaterHysteresisParams(const Scalar& soMax,
211 const Scalar& swMax,
212 const Scalar& swMin,
213 Params& params)
214 {
215 params.oilWaterParams().update(swMin, swMax, 1.0 - soMax);
216 }
217
218
219 /*
220 * Hysteresis parameters for gas-oil
221 * @see EclHysteresisTwoPhaseLawParams::sgmax(...)
222 * @see EclHysteresisTwoPhaseLawParams::shmax(...)
223 * @see EclHysteresisTwoPhaseLawParams::somin(...)
224 * \param params Parameters
225 */
226 static void gasOilHysteresisParams(Scalar& sgmax,
227 Scalar& shmax,
228 Scalar& somin,
229 const Params& params)
230 {
231 sgmax = 1.0 - params.gasOilParams().krnSwMdc();
232 shmax = params.gasOilParams().krwSwMdc();
233 somin = params.gasOilParams().pcSwMdc();
234 Valgrind::CheckDefined(sgmax);
235 Valgrind::CheckDefined(shmax);
236 Valgrind::CheckDefined(somin);
237 }
238
239 /*
240 * Hysteresis parameters for gas-oil
241 * @see EclHysteresisTwoPhaseLawParams::sgmax(...)
242 * @see EclHysteresisTwoPhaseLawParams::shmax(...)
243 * \param params Parameters
244 */
245 static void setGasOilHysteresisParams(const Scalar& sgmax,
246 const Scalar& shmax,
247 const Scalar& somin,
248 Params& params)
249 {
250 params.gasOilParams().update(somin , shmax, 1.0 - sgmax);
251 }
252
253 static Scalar trappedGasSaturation(const Params& params, bool maximumTrapping){
254 if(params.approach() == EclTwoPhaseApproach::GasOil)
255 return params.gasOilParams().SnTrapped(maximumTrapping);
256 if(params.approach() == EclTwoPhaseApproach::GasWater)
257 return params.gasWaterParams().SnTrapped(maximumTrapping);
258 return 0.0; // oil-water case
259 }
260
261 static Scalar strandedGasSaturation(const Params& params, Scalar Sg, Scalar Kg){
262 if(params.approach() == EclTwoPhaseApproach::GasOil)
263 return params.gasOilParams().SnStranded(Sg, Kg);
264 if(params.approach() == EclTwoPhaseApproach::GasWater)
265 return params.gasWaterParams().SnStranded(Sg, Kg);
266 return 0.0; // oil-water case
267 }
268
269 static Scalar trappedOilSaturation(const Params& params, bool maximumTrapping){
270 if(params.approach() == EclTwoPhaseApproach::GasOil)
271 return params.gasOilParams().SwTrapped();
272 if(params.approach() == EclTwoPhaseApproach::OilWater)
273 return params.oilWaterParams().SnTrapped(maximumTrapping);
274 return 0.0; // gas-water case
275 }
276
277 static Scalar trappedWaterSaturation(const Params& params){
278 if(params.approach() == EclTwoPhaseApproach::GasWater)
279 return params.gasWaterParams().SwTrapped();
280 if(params.approach() == EclTwoPhaseApproach::OilWater)
281 return params.oilWaterParams().SwTrapped();
282 return 0.0; // gas-oil case
283 }
284
285
295 template <class FluidState, class Evaluation = typename FluidState::Scalar>
296 static Evaluation pcgn(const Params& /* params */,
297 const FluidState& /* fs */)
298 {
299 throw std::logic_error("Not implemented: pcgn()");
300 }
301
311 template <class FluidState, class Evaluation = typename FluidState::Scalar>
312 static Evaluation pcnw(const Params& /* params */,
313 const FluidState& /* fs */)
314 {
315 throw std::logic_error("Not implemented: pcnw()");
316 }
317
321 template <class ContainerT, class FluidState>
322 static void saturations(ContainerT& /* values */,
323 const Params& /* params */,
324 const FluidState& /* fs */)
325 {
326 throw std::logic_error("Not implemented: saturations()");
327 }
328
332 template <class FluidState, class Evaluation = typename FluidState::Scalar>
333 static Evaluation Sg(const Params& /* params */,
334 const FluidState& /* fluidState */)
335 {
336 throw std::logic_error("Not implemented: Sg()");
337 }
338
342 template <class FluidState, class Evaluation = typename FluidState::Scalar>
343 static Evaluation Sn(const Params& /* params */,
344 const FluidState& /* fluidState */)
345 {
346 throw std::logic_error("Not implemented: Sn()");
347 }
348
352 template <class FluidState, class Evaluation = typename FluidState::Scalar>
353 static Evaluation Sw(const Params& /* params */,
354 const FluidState& /* fluidState */)
355 {
356 throw std::logic_error("Not implemented: Sw()");
357 }
358
374 template <class ContainerT, class FluidState>
375 static void relativePermeabilities(ContainerT& values,
376 const Params& params,
377 const FluidState& fluidState)
378 {
379 OPM_TIMEFUNCTION_LOCAL();
380 using Evaluation = typename std::remove_reference<decltype(values[0])>::type;
381
382 switch (params.approach()) {
383 case EclTwoPhaseApproach::GasOil: {
384 const Evaluation& So =
385 decay<Evaluation>(fluidState.saturation(oilPhaseIdx));
386
387 values[oilPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), So);
388 values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), So);
389 break;
390 }
391
392 case EclTwoPhaseApproach::OilWater: {
393 const Evaluation& sw =
394 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
395
396 values[waterPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), sw);
397 values[oilPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), sw);
398 break;
399 }
400
401 case EclTwoPhaseApproach::GasWater: {
402 const Evaluation& sw =
403 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
404
405 values[waterPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatKrw(params.gasWaterParams(), sw);
406 values[gasPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatKrn(params.gasWaterParams(), sw);
407
408 break;
409 }
410 }
411 }
412
416 template <class FluidState, class Evaluation = typename FluidState::Scalar>
417 static Evaluation krg(const Params& /* params */,
418 const FluidState& /* fluidState */)
419 {
420 throw std::logic_error("Not implemented: krg()");
421 }
422
426 template <class FluidState, class Evaluation = typename FluidState::Scalar>
427 static Evaluation krw(const Params& /* params */,
428 const FluidState& /* fluidState */)
429 {
430 throw std::logic_error("Not implemented: krw()");
431 }
432
436 template <class FluidState, class Evaluation = typename FluidState::Scalar>
437 static Evaluation krn(const Params& /* params */,
438 const FluidState& /* fluidState */)
439 {
440 throw std::logic_error("Not implemented: krn()");
441 }
442
443
451 template <class FluidState>
452 static bool updateHysteresis(Params& params, const FluidState& fluidState)
453 {
454 OPM_TIMEFUNCTION_LOCAL();
455 switch (params.approach()) {
456 case EclTwoPhaseApproach::GasOil: {
457 Scalar So = scalarValue(fluidState.saturation(oilPhaseIdx));
458 return params.gasOilParams().update(/*pcSw=*/So, /*krwSw=*/So, /*krnSw=*/So);
459 }
460
461 case EclTwoPhaseApproach::OilWater: {
462 Scalar sw = scalarValue(fluidState.saturation(waterPhaseIdx));
463 return params.oilWaterParams().update(/*pcSw=*/sw, /*krwSw=*/sw, /*krnSw=*/sw);
464 }
465
466 case EclTwoPhaseApproach::GasWater: {
467 Scalar sw = scalarValue(fluidState.saturation(waterPhaseIdx));
468 return params.gasWaterParams().update(/*pcSw=*/sw, /*krwSw=*/sw, /*krnSw=*/sw);
469 }
470 }
471
472 // Should not get here...
473 return false;
474 }
475};
476
477} // namespace Opm
478
479#endif
Implementation for the parameters required by the material law for two-phase simulations.
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
Some templates to wrap the valgrind client request macros.
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Definition EclTwoPhaseMaterial.hpp:57
static constexpr bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition EclTwoPhaseMaterial.hpp:108
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition EclTwoPhaseMaterial.hpp:375
static constexpr bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition EclTwoPhaseMaterial.hpp:92
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition EclTwoPhaseMaterial.hpp:353
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition EclTwoPhaseMaterial.hpp:296
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition EclTwoPhaseMaterial.hpp:322
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition EclTwoPhaseMaterial.hpp:437
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition EclTwoPhaseMaterial.hpp:343
static constexpr bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition EclTwoPhaseMaterial.hpp:100
static constexpr bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition EclTwoPhaseMaterial.hpp:104
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition EclTwoPhaseMaterial.hpp:417
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition EclTwoPhaseMaterial.hpp:333
static bool updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclTwoPhaseMaterial.hpp:452
static void capillaryPressures(ContainerT &values, const Params &params, const FluidState &fluidState)
Implements the multiplexer three phase capillary pressure law used by the ECLipse simulator.
Definition EclTwoPhaseMaterial.hpp:145
static constexpr bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition EclTwoPhaseMaterial.hpp:112
static Evaluation pcnw(const Params &, const FluidState &)
Capillary pressure between the non-wetting liquid (i.e., oil) and the wetting liquid (i....
Definition EclTwoPhaseMaterial.hpp:312
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition EclTwoPhaseMaterial.hpp:427
static constexpr bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition EclTwoPhaseMaterial.hpp:96
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30