opm-common
Loading...
Searching...
No Matches
EclMultiplexerMaterial.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_MULTIPLEXER_MATERIAL_HPP
28#define OPM_ECL_MULTIPLEXER_MATERIAL_HPP
29
32#include "EclStone1Material.hpp"
33#include "EclStone2Material.hpp"
35
36#include <opm/common/TimingMacros.hpp>
37
38#include <algorithm>
39#include <stdexcept>
40
41namespace Opm {
42
43#define OPM_ECL_MULTIPLEXER_MATERIAL_CALL(codeToCall, onePhaseCode) \
44 switch (params.approach()) { \
45 case EclMultiplexerApproach::Stone1: { \
46 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::Stone1; \
47 auto& realParams = params.template getRealParams<approach>(); \
48 using ActualLaw = Stone1Material; \
49 codeToCall; \
50 break; \
51 } \
52 case EclMultiplexerApproach::Stone2: { \
53 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::Stone2; \
54 auto& realParams = params.template getRealParams<approach>(); \
55 using ActualLaw = Stone2Material; \
56 codeToCall; \
57 break; \
58 } \
59 case EclMultiplexerApproach::Default: { \
60 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::Default; \
61 auto& realParams = params.template getRealParams<approach>(); \
62 using ActualLaw = DefaultMaterial; \
63 codeToCall; \
64 break; \
65 } \
66 case EclMultiplexerApproach::TwoPhase: { \
67 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::TwoPhase; \
68 auto& realParams = params.template getRealParams<approach>(); \
69 using ActualLaw = TwoPhaseMaterial; \
70 codeToCall; \
71 break; \
72 } \
73 case EclMultiplexerApproach::OnePhase: { \
74 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::OnePhase; \
75 onePhaseCode; \
76 break; \
77 } \
78 }
79
80// The static_assert does not compile with gcc 12 and earlier when placed in the multiplexer below.
81#if defined(__GNUC__) && (__GNUC__ < 13)
82 #define STATIC_ASSERT_ECL_MULTIPLEXER_UNLESS_GCC_LT_13 throw std::logic_error("Unhandled EclMultiplexerApproach")
83#else
84 #define STATIC_ASSERT_ECL_MULTIPLEXER_UNLESS_GCC_LT_13 static_assert(false, "Unhandled EclMultiplexerApproach")
85#endif
86
87#define OPM_ECL_MULTIPLEXER_MATERIAL_CALL_COMPILETIME(codeToCall, onePhaseCode) \
88 if constexpr (Head::approach == EclMultiplexerApproach::Stone1) { \
89 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::Stone1; \
90 auto& realParams = params.template getRealParams<approach>(); \
91 using ActualLaw = Stone1Material; \
92 codeToCall; \
93 } else if constexpr (Head::approach == EclMultiplexerApproach::Stone2) { \
94 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::Stone2; \
95 auto& realParams = params.template getRealParams<approach>(); \
96 using ActualLaw = Stone2Material; \
97 codeToCall; \
98 } else if constexpr (Head::approach == EclMultiplexerApproach::Default) { \
99 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::Default; \
100 auto& realParams = params.template getRealParams<approach>(); \
101 using ActualLaw = DefaultMaterial; \
102 codeToCall; \
103 } else if constexpr (Head::approach == EclMultiplexerApproach::TwoPhase) { \
104 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::TwoPhase; \
105 auto& realParams = params.template getRealParams<approach>(); \
106 using ActualLaw = TwoPhaseMaterial; \
107 codeToCall; \
108 } else if constexpr (Head::approach == EclMultiplexerApproach::OnePhase) { \
109 [[maybe_unused]] constexpr EclMultiplexerApproach approach = EclMultiplexerApproach::OnePhase; \
110 onePhaseCode; \
111 } else { \
112 STATIC_ASSERT_ECL_MULTIPLEXER_UNLESS_GCC_LT_13; \
113 }
114
115// Pass this for the onePhaseCode argument if nothing is to be done.
116inline void doNothing() { };
117
124template <class TraitsT,
125 class GasOilMaterialLawT,
126 class OilWaterMaterialLawT,
127 class GasWaterMaterialLawT,
128 class ParamsT = EclMultiplexerMaterialParams<TraitsT,
129 GasOilMaterialLawT,
130 OilWaterMaterialLawT,
131 GasWaterMaterialLawT> >
132class EclMultiplexerMaterial : public TraitsT
133{
134public:
135 using GasOilMaterialLaw = GasOilMaterialLawT;
136 using OilWaterMaterialLaw = OilWaterMaterialLawT;
137 using GasWaterMaterialLaw = GasWaterMaterialLawT;
138
143
144 // some safety checks
145 static_assert(TraitsT::numPhases == 3,
146 "The number of phases considered by this capillary pressure "
147 "law is always three!");
148 static_assert(GasOilMaterialLaw::numPhases == 2,
149 "The number of phases considered by the gas-oil capillary "
150 "pressure law must be two!");
151 static_assert(OilWaterMaterialLaw::numPhases == 2,
152 "The number of phases considered by the oil-water capillary "
153 "pressure law must be two!");
154 static_assert(GasWaterMaterialLaw::numPhases == 2,
155 "The number of phases considered by the gas-water capillary "
156 "pressure law must be two!");
157 static_assert(std::is_same<typename GasOilMaterialLaw::Scalar,
158 typename OilWaterMaterialLaw::Scalar>::value,
159 "The two two-phase capillary pressure laws must use the same "
160 "type of floating point values.");
161
162 using Traits = TraitsT;
163 using Params = ParamsT;
164 using Scalar = typename Traits::Scalar;
165
166 static constexpr int numPhases = 3;
167 static constexpr int waterPhaseIdx = Traits::wettingPhaseIdx;
168 static constexpr int oilPhaseIdx = Traits::nonWettingPhaseIdx;
169 static constexpr int gasPhaseIdx = Traits::gasPhaseIdx;
170
173 static constexpr bool implementsTwoPhaseApi = false;
174
177 static constexpr bool implementsTwoPhaseSatApi = false;
178
181 static constexpr bool isSaturationDependent = true;
182
185 static constexpr bool isPressureDependent = false;
186
189 static constexpr bool isTemperatureDependent = false;
190
193 static constexpr bool isCompositionDependent = false;
194
209 template <class ContainerT, class FluidState, class ...Args>
210 static void capillaryPressures(ContainerT& values,
211 const Params& params,
212 const FluidState& fluidState)
213 {
214 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
215 if constexpr (FrontIsEclMultiplexerDispatchV<Args...>) {
216 capillaryPressuresT<ContainerT, FluidState, Args...>(values, params, fluidState);
217 return;
218 }
219 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(ActualLaw::capillaryPressures(values, realParams, fluidState),
220 values[0] = 0.0);
221 }
222
223
224 template <class ContainerT, class FluidState, class Head, class ...Args>
225 static void capillaryPressuresT(ContainerT& values,
226 const Params& params,
227 const FluidState& fluidState)
228 {
229#define OPM_LOCAL_TEMPLATE_ARGS ContainerT, FluidState, Args...
230 OPM_ECL_MULTIPLEXER_MATERIAL_CALL_COMPILETIME(
231 ActualLaw::template capillaryPressures<OPM_LOCAL_TEMPLATE_ARGS>(values, realParams, fluidState),
232 values[0] = 0.0
233 );
234#undef OPM_LOCAL_TEMPLATE_ARGS
235 }
236
237 /*
238 * Hysteresis parameters for oil-water
239 * @see EclHysteresisTwoPhaseLawParams::soMax(...)
240 * @see EclHysteresisTwoPhaseLawParams::swMax(...)
241 * @see EclHysteresisTwoPhaseLawParams::swMin(...)
242 * \param params Parameters
243 */
244 static void oilWaterHysteresisParams(Scalar& soMax,
245 Scalar& swMax,
246 Scalar& swMin,
247 const Params& params)
248 {
249 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
250 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(ActualLaw::oilWaterHysteresisParams(soMax, swMax, swMin, realParams),
251 doNothing());
252 }
253
254 /*
255 * Hysteresis parameters for oil-water
256 * @see EclHysteresisTwoPhaseLawParams::soMax(...)
257 * @see EclHysteresisTwoPhaseLawParams::swMax(...)
258 * @see EclHysteresisTwoPhaseLawParams::swMin(...)
259 * \param params Parameters
260 */
261 static void setOilWaterHysteresisParams(const Scalar& soMax,
262 const Scalar& swMax,
263 const Scalar& swMin,
264 Params& params)
265 {
266 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
267 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(ActualLaw::setOilWaterHysteresisParams(soMax, swMax, swMin, realParams),
268 doNothing());
269 }
270
271 /*
272 * Hysteresis parameters for gas-oil
273 * @see EclHysteresisTwoPhaseLawParams::sgmax(...)
274 * @see EclHysteresisTwoPhaseLawParams::shmax(...)
275 * @see EclHysteresisTwoPhaseLawParams::somin(...)
276 * \param params Parameters
277 */
278 static void gasOilHysteresisParams(Scalar& sgmax,
279 Scalar& shmax,
280 Scalar& somin,
281 const Params& params)
282 {
283 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
284 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(ActualLaw::gasOilHysteresisParams(sgmax, shmax, somin, realParams),
285 doNothing());
286 }
287
288 static Scalar trappedGasSaturation(const Params& params, bool maximumTrapping)
289 {
290 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
291 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(return ActualLaw::trappedGasSaturation(realParams, maximumTrapping),
292 return 0.0);
293 return 0.0;
294 }
295
296 static Scalar strandedGasSaturation(const Params& params, Scalar Sg, Scalar Kg)
297 {
298 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
299 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(return ActualLaw::strandedGasSaturation(realParams, Sg, Kg),
300 return 0.0);
301 return 0.0;
302 }
303
304 static Scalar trappedOilSaturation(const Params& params, bool maximumTrapping)
305 {
306 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
307 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(return ActualLaw::trappedOilSaturation(realParams, maximumTrapping),
308 return 0.0);
309 return 0.0;
310 }
311
312 static Scalar trappedWaterSaturation(const Params& params)
313 {
314 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
315 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(return ActualLaw::trappedWaterSaturation(realParams),
316 return 0.0);
317 return 0.0;
318 }
319 /*
320 * Hysteresis parameters for gas-oil
321 * @see EclHysteresisTwoPhaseLawParams::sgmax(...)
322 * @see EclHysteresisTwoPhaseLawParams::shmax(...)
323 * @see EclHysteresisTwoPhaseLawParams::somin(...)
324 * \param params Parameters
325 */
326 static void setGasOilHysteresisParams(const Scalar& sgmax,
327 const Scalar& shmax,
328 const Scalar& somin,
329 Params& params)
330 {
331 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
332 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(ActualLaw::setGasOilHysteresisParams(sgmax, shmax, somin, realParams),
333 doNothing());
334 }
335
345 template <class FluidState, class Evaluation = typename FluidState::Scalar>
346 static Evaluation pcgn(const Params& /* params */,
347 const FluidState& /* fs */)
348 {
349 throw std::logic_error("Not implemented: pcgn()");
350 }
351
361 template <class FluidState, class Evaluation = typename FluidState::Scalar>
362 static Evaluation pcnw(const Params& /* params */,
363 const FluidState& /* fs */)
364 {
365 throw std::logic_error("Not implemented: pcnw()");
366 }
367
371 template <class ContainerT, class FluidState>
372 static void saturations(ContainerT& /* values */,
373 const Params& /* params */,
374 const FluidState& /* fs */)
375 {
376 throw std::logic_error("Not implemented: saturations()");
377 }
378
382 template <class FluidState, class Evaluation = typename FluidState::Scalar>
383 static Evaluation Sg(const Params& /* params */,
384 const FluidState& /* fluidState */)
385 {
386 throw std::logic_error("Not implemented: Sg()");
387 }
388
392 template <class FluidState, class Evaluation = typename FluidState::Scalar>
393 static Evaluation Sn(const Params& /* params */,
394 const FluidState& /* fluidState */)
395 {
396 throw std::logic_error("Not implemented: Sn()");
397 }
398
402 template <class FluidState, class Evaluation = typename FluidState::Scalar>
403 static Evaluation Sw(const Params& /* params */,
404 const FluidState& /* fluidState */)
405 {
406 throw std::logic_error("Not implemented: Sw()");
407 }
408
424 template <class ContainerT, class FluidState, class ...Args>
425 static void relativePermeabilities(ContainerT& values,
426 const Params& params,
427 const FluidState& fluidState)
428 {
429 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
430 if constexpr (FrontIsEclMultiplexerDispatchV<Args...>) {
431 relativePermeabilitiesT<ContainerT, FluidState, Args...>(values, params, fluidState);
432 return;
433 }
434 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(ActualLaw::relativePermeabilities(values, realParams, fluidState),
435 values[0] = 1.0);
436 }
437
438 template <class ContainerT, class FluidState, class Head, class ...Args>
439 static void relativePermeabilitiesT(ContainerT& values,
440 const Params& params,
441 const FluidState& fluidState)
442 {
443#define OPM_LOCAL_TEMPLATE_ARGS ContainerT, FluidState, Args...
444 OPM_ECL_MULTIPLEXER_MATERIAL_CALL_COMPILETIME(
445 ActualLaw::template relativePermeabilities<OPM_LOCAL_TEMPLATE_ARGS>(values, realParams, fluidState),
446 values[0] = 1.0
447 );
448#undef OPM_LOCAL_TEMPLATE_ARGS
449 }
450
451
452
453
457 template <class Evaluation, class FluidState>
458 static Evaluation relpermOilInOilGasSystem(const Params& params,
459 const FluidState& fluidState)
460 {
461 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
462 switch (params.approach()) {
463 case EclMultiplexerApproach::Stone1:
464 return Stone1Material::template relpermOilInOilGasSystem<Evaluation>
465 (params.template getRealParams<EclMultiplexerApproach::Stone1>(),
466 fluidState);
467
468 case EclMultiplexerApproach::Stone2:
469 return Stone2Material::template relpermOilInOilGasSystem<Evaluation>
470 (params.template getRealParams<EclMultiplexerApproach::Stone2>(),
471 fluidState);
472
473 case EclMultiplexerApproach::Default:
474 return DefaultMaterial::template relpermOilInOilGasSystem<Evaluation>
475 (params.template getRealParams<EclMultiplexerApproach::Default>(),
476 fluidState);
477
478 default:
479 throw std::logic_error {
480 "relpermOilInOilGasSystem() is specific to three phases"
481 };
482 }
483 }
484
488 template <class Evaluation, class FluidState>
489 static Evaluation relpermOilInOilWaterSystem(const Params& params,
490 const FluidState& fluidState)
491 {
492 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
493 switch (params.approach()) {
494 case EclMultiplexerApproach::Stone1:
495 return Stone1Material::template relpermOilInOilWaterSystem<Evaluation>
496 (params.template getRealParams<EclMultiplexerApproach::Stone1>(),
497 fluidState);
498
499 case EclMultiplexerApproach::Stone2:
500 return Stone2Material::template relpermOilInOilWaterSystem<Evaluation>
501 (params.template getRealParams<EclMultiplexerApproach::Stone2>(),
502 fluidState);
503
504 case EclMultiplexerApproach::Default:
505 return DefaultMaterial::template relpermOilInOilWaterSystem<Evaluation>
506 (params.template getRealParams<EclMultiplexerApproach::Default>(),
507 fluidState);
508
509 default:
510 throw std::logic_error {
511 "relpermOilInOilWaterSystem() is specific to three phases"
512 };
513 }
514 }
515
519 template <class FluidState, class Evaluation = typename FluidState::Scalar>
520 static Evaluation krg(const Params& /* params */,
521 const FluidState& /* fluidState */)
522 {
523 throw std::logic_error("Not implemented: krg()");
524 }
525
529 template <class FluidState, class Evaluation = typename FluidState::Scalar>
530 static Evaluation krw(const Params& /* params */,
531 const FluidState& /* fluidState */)
532 {
533 throw std::logic_error("Not implemented: krw()");
534 }
535
539 template <class FluidState, class Evaluation = typename FluidState::Scalar>
540 static Evaluation krn(const Params& /* params */,
541 const FluidState& /* fluidState */)
542 {
543 throw std::logic_error("Not implemented: krn()");
544 }
545
546
554 template <class FluidState>
555 static bool updateHysteresis(Params& params, const FluidState& fluidState)
556 {
557 OPM_TIMEFUNCTION_LOCAL(Subsystem::SatProps);
558 OPM_ECL_MULTIPLEXER_MATERIAL_CALL(return ActualLaw::updateHysteresis(realParams, fluidState),
559 return false);
560 return false;
561 }
562};
563
564} // namespace Opm
565
566#endif
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Multiplexer implementation for the parameters required by the multiplexed three-phase material law.
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition EclDefaultMaterial.hpp:62
Multiplexer implementation for the parameters required by the multiplexed three-phase material law.
Definition EclMultiplexerMaterialParams.hpp:75
Implements a multiplexer class that provides all three phase capillary pressure laws used by the ECLi...
Definition EclMultiplexerMaterial.hpp:133
static Evaluation relpermOilInOilWaterSystem(const Params &params, const FluidState &fluidState)
The relative permeability of oil in oil/water system.
Definition EclMultiplexerMaterial.hpp:489
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition EclMultiplexerMaterial.hpp:530
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition EclMultiplexerMaterial.hpp:383
static Evaluation relpermOilInOilGasSystem(const Params &params, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition EclMultiplexerMaterial.hpp:458
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition EclMultiplexerMaterial.hpp:346
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition EclMultiplexerMaterial.hpp:425
static Evaluation pcnw(const Params &, const FluidState &)
Capillary pressure between the non-wetting liquid (i.e., oil) and the wetting liquid (i....
Definition EclMultiplexerMaterial.hpp:362
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition EclMultiplexerMaterial.hpp:393
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition EclMultiplexerMaterial.hpp:520
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 EclMultiplexerMaterial.hpp:210
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition EclMultiplexerMaterial.hpp:372
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition EclMultiplexerMaterial.hpp:540
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition EclMultiplexerMaterial.hpp:403
static bool updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclMultiplexerMaterial.hpp:555
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition EclStone1Material.hpp:61
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition EclStone2Material.hpp:61
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Definition EclTwoPhaseMaterial.hpp:57
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30