My Project
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
49template <class TraitsT,
50 class GasOilMaterialLawT,
51 class OilWaterMaterialLawT,
52 class GasWaterMaterialLawT,
53 class ParamsT = EclMultiplexerMaterialParams<TraitsT,
54 GasOilMaterialLawT,
55 OilWaterMaterialLawT,
56 GasWaterMaterialLawT> >
57class EclMultiplexerMaterial : public TraitsT
58{
59public:
60 using GasOilMaterialLaw = GasOilMaterialLawT;
61 using OilWaterMaterialLaw = OilWaterMaterialLawT;
62 using GasWaterMaterialLaw = GasWaterMaterialLawT;
63
68
69 // some safety checks
70 static_assert(TraitsT::numPhases == 3,
71 "The number of phases considered by this capillary pressure "
72 "law is always three!");
73 static_assert(GasOilMaterialLaw::numPhases == 2,
74 "The number of phases considered by the gas-oil capillary "
75 "pressure law must be two!");
76 static_assert(OilWaterMaterialLaw::numPhases == 2,
77 "The number of phases considered by the oil-water capillary "
78 "pressure law must be two!");
79 static_assert(GasWaterMaterialLaw::numPhases == 2,
80 "The number of phases considered by the gas-water capillary "
81 "pressure law must be two!");
82 static_assert(std::is_same<typename GasOilMaterialLaw::Scalar,
83 typename OilWaterMaterialLaw::Scalar>::value,
84 "The two two-phase capillary pressure laws must use the same "
85 "type of floating point values.");
86
87 using Traits = TraitsT;
88 using Params = ParamsT;
89 using Scalar = typename Traits::Scalar;
90
91 static constexpr int numPhases = 3;
92 static constexpr int waterPhaseIdx = Traits::wettingPhaseIdx;
93 static constexpr int oilPhaseIdx = Traits::nonWettingPhaseIdx;
94 static constexpr int gasPhaseIdx = Traits::gasPhaseIdx;
95
98 static constexpr bool implementsTwoPhaseApi = false;
99
102 static constexpr bool implementsTwoPhaseSatApi = false;
103
106 static constexpr bool isSaturationDependent = true;
107
110 static constexpr bool isPressureDependent = false;
111
114 static constexpr bool isTemperatureDependent = false;
115
118 static constexpr bool isCompositionDependent = false;
119
134 template <class ContainerT, class FluidState>
135 static void capillaryPressures(ContainerT& values,
136 const Params& params,
137 const FluidState& fluidState)
138 {
139 OPM_TIMEBLOCK_LOCAL(capillaryPressures);
140 switch (params.approach()) {
141 case EclMultiplexerApproach::Stone1:
143 params.template getRealParams<EclMultiplexerApproach::Stone1>(),
144 fluidState);
145 break;
146
147 case EclMultiplexerApproach::Stone2:
149 params.template getRealParams<EclMultiplexerApproach::Stone2>(),
150 fluidState);
151 break;
152
153 case EclMultiplexerApproach::Default:
155 params.template getRealParams<EclMultiplexerApproach::Default>(),
156 fluidState);
157 break;
158
159 case EclMultiplexerApproach::TwoPhase:
161 params.template getRealParams<EclMultiplexerApproach::TwoPhase>(),
162 fluidState);
163 break;
164
165 case EclMultiplexerApproach::OnePhase:
166 values[0] = 0.0;
167 break;
168 }
169 }
170
171 /*
172 * Hysteresis parameters for oil-water
173 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
174 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
175 * \param params Parameters
176 */
177 static void oilWaterHysteresisParams(Scalar& pcSwMdc,
178 Scalar& krnSwMdc,
179 const Params& params)
180 {
181 switch (params.approach()) {
182 case EclMultiplexerApproach::Stone1:
183 Stone1Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
184 params.template getRealParams<EclMultiplexerApproach::Stone1>());
185 break;
186
187 case EclMultiplexerApproach::Stone2:
188 Stone2Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
189 params.template getRealParams<EclMultiplexerApproach::Stone2>());
190 break;
191
192 case EclMultiplexerApproach::Default:
193 DefaultMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
194 params.template getRealParams<EclMultiplexerApproach::Default>());
195 break;
196
197 case EclMultiplexerApproach::TwoPhase:
198 TwoPhaseMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
199 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
200 break;
201
202 case EclMultiplexerApproach::OnePhase:
203 // Do nothing.
204 break;
205 }
206 }
207
208 /*
209 * Hysteresis parameters for oil-water
210 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
211 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
212 * \param params Parameters
213 */
214 static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
215 const Scalar& krnSwMdc,
216 Params& params)
217 {
218 switch (params.approach()) {
219 case EclMultiplexerApproach::Stone1:
220 Stone1Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
221 params.template getRealParams<EclMultiplexerApproach::Stone1>());
222 break;
223
224 case EclMultiplexerApproach::Stone2:
225 Stone2Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
226 params.template getRealParams<EclMultiplexerApproach::Stone2>());
227 break;
228
229 case EclMultiplexerApproach::Default:
230 DefaultMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
231 params.template getRealParams<EclMultiplexerApproach::Default>());
232 break;
233
234 case EclMultiplexerApproach::TwoPhase:
235 TwoPhaseMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
236 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
237 break;
238
239 case EclMultiplexerApproach::OnePhase:
240 // Do nothing.
241 break;
242 }
243 }
244
245 /*
246 * Hysteresis parameters for gas-oil
247 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
248 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
249 * \param params Parameters
250 */
251 static void gasOilHysteresisParams(Scalar& pcSwMdc,
252 Scalar& krnSwMdc,
253 const Params& params)
254 {
255 switch (params.approach()) {
256 case EclMultiplexerApproach::Stone1:
257 Stone1Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
258 params.template getRealParams<EclMultiplexerApproach::Stone1>());
259 break;
260
261 case EclMultiplexerApproach::Stone2:
262 Stone2Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
263 params.template getRealParams<EclMultiplexerApproach::Stone2>());
264 break;
265
266 case EclMultiplexerApproach::Default:
267 DefaultMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
268 params.template getRealParams<EclMultiplexerApproach::Default>());
269 break;
270
271 case EclMultiplexerApproach::TwoPhase:
272 TwoPhaseMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
273 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
274 break;
275
276 case EclMultiplexerApproach::OnePhase:
277 // Do nothing.
278 break;
279 }
280 }
281
282 /*
283 * Hysteresis parameters for gas-oil
284 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
285 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
286 * \param params Parameters
287 */
288 static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
289 const Scalar& krnSwMdc,
290 Params& params)
291 {
292 switch (params.approach()) {
293 case EclMultiplexerApproach::Stone1:
294 Stone1Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
295 params.template getRealParams<EclMultiplexerApproach::Stone1>());
296 break;
297
298 case EclMultiplexerApproach::Stone2:
299 Stone2Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
300 params.template getRealParams<EclMultiplexerApproach::Stone2>());
301 break;
302
303 case EclMultiplexerApproach::Default:
304 DefaultMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
305 params.template getRealParams<EclMultiplexerApproach::Default>());
306 break;
307
308 case EclMultiplexerApproach::TwoPhase:
309 TwoPhaseMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
310 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
311 break;
312
313 case EclMultiplexerApproach::OnePhase:
314 // Do nothing.
315 break;
316 }
317 }
318
328 template <class FluidState, class Evaluation = typename FluidState::Scalar>
329 static Evaluation pcgn(const Params& /* params */,
330 const FluidState& /* fs */)
331 {
332 throw std::logic_error("Not implemented: pcgn()");
333 }
334
344 template <class FluidState, class Evaluation = typename FluidState::Scalar>
345 static Evaluation pcnw(const Params& /* params */,
346 const FluidState& /* fs */)
347 {
348 throw std::logic_error("Not implemented: pcnw()");
349 }
350
354 template <class ContainerT, class FluidState>
355 static void saturations(ContainerT& /* values */,
356 const Params& /* params */,
357 const FluidState& /* fs */)
358 {
359 throw std::logic_error("Not implemented: saturations()");
360 }
361
365 template <class FluidState, class Evaluation = typename FluidState::Scalar>
366 static Evaluation Sg(const Params& /* params */,
367 const FluidState& /* fluidState */)
368 {
369 throw std::logic_error("Not implemented: Sg()");
370 }
371
375 template <class FluidState, class Evaluation = typename FluidState::Scalar>
376 static Evaluation Sn(const Params& /* params */,
377 const FluidState& /* fluidState */)
378 {
379 throw std::logic_error("Not implemented: Sn()");
380 }
381
385 template <class FluidState, class Evaluation = typename FluidState::Scalar>
386 static Evaluation Sw(const Params& /* params */,
387 const FluidState& /* fluidState */)
388 {
389 throw std::logic_error("Not implemented: Sw()");
390 }
391
407 template <class ContainerT, class FluidState>
408 static void relativePermeabilities(ContainerT& values,
409 const Params& params,
410 const FluidState& fluidState)
411 {
412 OPM_TIMEBLOCK_LOCAL(relativePermeabilities);
413 switch (params.approach()) {
414 case EclMultiplexerApproach::Stone1:
416 params.template getRealParams<EclMultiplexerApproach::Stone1>(),
417 fluidState);
418 break;
419
420 case EclMultiplexerApproach::Stone2:
422 params.template getRealParams<EclMultiplexerApproach::Stone2>(),
423 fluidState);
424 break;
425
426 case EclMultiplexerApproach::Default:
428 params.template getRealParams<EclMultiplexerApproach::Default>(),
429 fluidState);
430 break;
431
432 case EclMultiplexerApproach::TwoPhase:
434 params.template getRealParams<EclMultiplexerApproach::TwoPhase>(),
435 fluidState);
436 break;
437
438 case EclMultiplexerApproach::OnePhase:
439 values[0] = 1.0;
440 break;
441
442 default:
443 throw std::logic_error("Not implemented: relativePermeabilities() option for unknown EclMultiplexerApproach (="
444 + std::to_string(static_cast<int>(params.approach())) + ")");
445 }
446 }
447
451 template <class Evaluation, class FluidState>
452 static Evaluation relpermOilInOilGasSystem(const Params& params,
453 const FluidState& fluidState)
454 {
455 switch (params.approach()) {
456 case EclMultiplexerApproach::Stone1:
457 return Stone1Material::template relpermOilInOilGasSystem<Evaluation>
458 (params.template getRealParams<EclMultiplexerApproach::Stone1>(),
459 fluidState);
460
461 case EclMultiplexerApproach::Stone2:
462 return Stone2Material::template relpermOilInOilGasSystem<Evaluation>
463 (params.template getRealParams<EclMultiplexerApproach::Stone2>(),
464 fluidState);
465
466 case EclMultiplexerApproach::Default:
467 return DefaultMaterial::template relpermOilInOilGasSystem<Evaluation>
468 (params.template getRealParams<EclMultiplexerApproach::Default>(),
469 fluidState);
470
471 default:
472 throw std::logic_error {
473 "relpermOilInOilGasSystem() is specific to three phases"
474 };
475 }
476 }
477
481 template <class Evaluation, class FluidState>
482 static Evaluation relpermOilInOilWaterSystem(const Params& params,
483 const FluidState& fluidState)
484 {
485 switch (params.approach()) {
486 case EclMultiplexerApproach::Stone1:
487 return Stone1Material::template relpermOilInOilWaterSystem<Evaluation>
488 (params.template getRealParams<EclMultiplexerApproach::Stone1>(),
489 fluidState);
490
491 case EclMultiplexerApproach::Stone2:
492 return Stone2Material::template relpermOilInOilWaterSystem<Evaluation>
493 (params.template getRealParams<EclMultiplexerApproach::Stone2>(),
494 fluidState);
495
496 case EclMultiplexerApproach::Default:
497 return DefaultMaterial::template relpermOilInOilWaterSystem<Evaluation>
498 (params.template getRealParams<EclMultiplexerApproach::Default>(),
499 fluidState);
500
501 default:
502 throw std::logic_error {
503 "relpermOilInOilWaterSystem() is specific to three phases"
504 };
505 }
506 }
507
511 template <class FluidState, class Evaluation = typename FluidState::Scalar>
512 static Evaluation krg(const Params& /* params */,
513 const FluidState& /* fluidState */)
514 {
515 throw std::logic_error("Not implemented: krg()");
516 }
517
521 template <class FluidState, class Evaluation = typename FluidState::Scalar>
522 static Evaluation krw(const Params& /* params */,
523 const FluidState& /* fluidState */)
524 {
525 throw std::logic_error("Not implemented: krw()");
526 }
527
531 template <class FluidState, class Evaluation = typename FluidState::Scalar>
532 static Evaluation krn(const Params& /* params */,
533 const FluidState& /* fluidState */)
534 {
535 throw std::logic_error("Not implemented: krn()");
536 }
537
538
546 template <class FluidState>
547 static void updateHysteresis(Params& params, const FluidState& fluidState)
548 {
549 switch (params.approach()) {
550 case EclMultiplexerApproach::Stone1:
551 Stone1Material::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::Stone1>(),
552 fluidState);
553 break;
554
555 case EclMultiplexerApproach::Stone2:
556 Stone2Material::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::Stone2>(),
557 fluidState);
558 break;
559
560 case EclMultiplexerApproach::Default:
561 DefaultMaterial::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::Default>(),
562 fluidState);
563 break;
564
565 case EclMultiplexerApproach::TwoPhase:
566 TwoPhaseMaterial::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::TwoPhase>(),
567 fluidState);
568 break;
569 case EclMultiplexerApproach::OnePhase:
570 break;
571 }
572 }
573};
574
575} // namespace Opm
576
577#endif
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclDefaultMaterial.hpp:61
static void capillaryPressures(ContainerT &values, const Params &params, const FluidState &state)
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclDefaultMaterial.hpp:136
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclDefaultMaterial.hpp:310
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclDefaultMaterial.hpp:424
Implements a multiplexer class that provides all three phase capillary pressure laws used by the ECLi...
Definition: EclMultiplexerMaterial.hpp:58
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:135
static Evaluation relpermOilInOilWaterSystem(const Params &params, const FluidState &fluidState)
The relative permeability of oil in oil/water system.
Definition: EclMultiplexerMaterial.hpp:482
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition: EclMultiplexerMaterial.hpp:522
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition: EclMultiplexerMaterial.hpp:366
static Evaluation relpermOilInOilGasSystem(const Params &params, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition: EclMultiplexerMaterial.hpp:452
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:329
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:345
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclMultiplexerMaterial.hpp:547
static constexpr bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: EclMultiplexerMaterial.hpp:110
static constexpr bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclMultiplexerMaterial.hpp:102
static constexpr bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclMultiplexerMaterial.hpp:106
static constexpr bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclMultiplexerMaterial.hpp:114
static constexpr bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclMultiplexerMaterial.hpp:98
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclMultiplexerMaterial.hpp:408
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:376
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition: EclMultiplexerMaterial.hpp:512
static constexpr bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: EclMultiplexerMaterial.hpp:118
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: EclMultiplexerMaterial.hpp:355
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:532
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition: EclMultiplexerMaterial.hpp:386
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition: EclStone1Material.hpp:60
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclStone1Material.hpp:313
static void capillaryPressures(ContainerT &values, const Params &params, const FluidState &state)
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclStone1Material.hpp:135
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclStone1Material.hpp:420
Implements the second phase capillary pressure/relperm law suggested by Stone as used by the ECLipse ...
Definition: EclStone2Material.hpp:61
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclStone2Material.hpp:404
static void capillaryPressures(ContainerT &values, const Params &params, const FluidState &state)
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclStone2Material.hpp:136
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclStone2Material.hpp:314
Implements a multiplexer class that provides ECL saturation functions for twophase simulations.
Definition: EclTwoPhaseMaterial.hpp:57
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclTwoPhaseMaterial.hpp:317
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclTwoPhaseMaterial.hpp:393
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:129
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:30