My Project
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
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_TIMEFUNCTION_LOCAL();
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 OPM_TIMEFUNCTION_LOCAL();
182 switch (params.approach()) {
183 case EclMultiplexerApproach::Stone1:
184 Stone1Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
185 params.template getRealParams<EclMultiplexerApproach::Stone1>());
186 break;
187
188 case EclMultiplexerApproach::Stone2:
189 Stone2Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
190 params.template getRealParams<EclMultiplexerApproach::Stone2>());
191 break;
192
193 case EclMultiplexerApproach::Default:
194 DefaultMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
195 params.template getRealParams<EclMultiplexerApproach::Default>());
196 break;
197
198 case EclMultiplexerApproach::TwoPhase:
199 TwoPhaseMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
200 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
201 break;
202
203 case EclMultiplexerApproach::OnePhase:
204 // Do nothing.
205 break;
206 }
207 }
208
209 /*
210 * Hysteresis parameters for oil-water
211 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
212 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
213 * \param params Parameters
214 */
215 static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
216 const Scalar& krnSwMdc,
217 Params& params)
218 {
219 OPM_TIMEFUNCTION_LOCAL();
220 switch (params.approach()) {
221 case EclMultiplexerApproach::Stone1:
222 Stone1Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
223 params.template getRealParams<EclMultiplexerApproach::Stone1>());
224 break;
225
226 case EclMultiplexerApproach::Stone2:
227 Stone2Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
228 params.template getRealParams<EclMultiplexerApproach::Stone2>());
229 break;
230
231 case EclMultiplexerApproach::Default:
232 DefaultMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
233 params.template getRealParams<EclMultiplexerApproach::Default>());
234 break;
235
236 case EclMultiplexerApproach::TwoPhase:
237 TwoPhaseMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
238 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
239 break;
240
241 case EclMultiplexerApproach::OnePhase:
242 // Do nothing.
243 break;
244 }
245 }
246
247 /*
248 * Hysteresis parameters for gas-oil
249 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
250 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
251 * \param params Parameters
252 */
253 static void gasOilHysteresisParams(Scalar& pcSwMdc,
254 Scalar& krnSwMdc,
255 const Params& params)
256 {
257 OPM_TIMEFUNCTION_LOCAL();
258 switch (params.approach()) {
259 case EclMultiplexerApproach::Stone1:
260 Stone1Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
261 params.template getRealParams<EclMultiplexerApproach::Stone1>());
262 break;
263
264 case EclMultiplexerApproach::Stone2:
265 Stone2Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
266 params.template getRealParams<EclMultiplexerApproach::Stone2>());
267 break;
268
269 case EclMultiplexerApproach::Default:
270 DefaultMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
271 params.template getRealParams<EclMultiplexerApproach::Default>());
272 break;
273
274 case EclMultiplexerApproach::TwoPhase:
275 TwoPhaseMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
276 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
277 break;
278
279 case EclMultiplexerApproach::OnePhase:
280 // Do nothing.
281 break;
282 }
283 }
284
285 static Scalar trappedGasSaturation(const Params& params)
286 {
287 OPM_TIMEFUNCTION_LOCAL();
288 switch (params.approach()) {
289 case EclMultiplexerApproach::Stone1:
290 return params.template getRealParams<EclMultiplexerApproach::Stone1>().gasOilParams().SnTrapped();
291 case EclMultiplexerApproach::Stone2:
292 return params.template getRealParams<EclMultiplexerApproach::Stone2>().gasOilParams().SnTrapped();
293 case EclMultiplexerApproach::Default:
294 return params.template getRealParams<EclMultiplexerApproach::Default>().gasOilParams().SnTrapped();
295 case EclMultiplexerApproach::TwoPhase:
296 if(params.template getRealParams<EclMultiplexerApproach::TwoPhase>().approach() == EclTwoPhaseApproach::GasOil)
297 return params.template getRealParams<EclMultiplexerApproach::TwoPhase>().gasOilParams().SnTrapped();
298 if(params.template getRealParams<EclMultiplexerApproach::TwoPhase>().approach() == EclTwoPhaseApproach::GasWater)
299 return params.template getRealParams<EclMultiplexerApproach::TwoPhase>().gasWaterParams().SnTrapped();
300 return 0.0; // oil-water case
301 case EclMultiplexerApproach::OnePhase:
302 return 0.0;
303 }
304 return 0.0;
305 }
306
307 static Scalar trappedOilSaturation(const Params& params)
308 {
309 OPM_TIMEFUNCTION_LOCAL();
310 switch (params.approach()) {
311 case EclMultiplexerApproach::Stone1:
312 return params.template getRealParams<EclMultiplexerApproach::Stone1>().oilWaterParams().SnTrapped();
313 case EclMultiplexerApproach::Stone2:
314 return params.template getRealParams<EclMultiplexerApproach::Stone2>().oilWaterParams().SnTrapped();
315 case EclMultiplexerApproach::Default:
316 return params.template getRealParams<EclMultiplexerApproach::Default>().oilWaterParams().SnTrapped();
317 case EclMultiplexerApproach::TwoPhase:
318 if(params.template getRealParams<EclMultiplexerApproach::TwoPhase>().approach() == EclTwoPhaseApproach::GasOil)
319 return params.template getRealParams<EclMultiplexerApproach::TwoPhase>().gasOilParams().SwTrapped();
320 if(params.template getRealParams<EclMultiplexerApproach::TwoPhase>().approach() == EclTwoPhaseApproach::OilWater)
321 return params.template getRealParams<EclMultiplexerApproach::TwoPhase>().oilWaterParams().SnTrapped();
322 return 0.0; // gas-water case
323 case EclMultiplexerApproach::OnePhase:
324 return 0.0;
325 }
326 return 0.0;
327 }
328
329 static Scalar trappedWaterSaturation(const Params& params)
330 {
331 OPM_TIMEFUNCTION_LOCAL();
332 switch (params.approach()) {
333 case EclMultiplexerApproach::Stone1:
334 return params.template getRealParams<EclMultiplexerApproach::Stone1>().oilWaterParams().SwTrapped();
335 case EclMultiplexerApproach::Stone2:
336 return params.template getRealParams<EclMultiplexerApproach::Stone2>().oilWaterParams().SwTrapped();
337 case EclMultiplexerApproach::Default:
338 return params.template getRealParams<EclMultiplexerApproach::Default>().oilWaterParams().SwTrapped();
339 case EclMultiplexerApproach::TwoPhase:
340 if(params.template getRealParams<EclMultiplexerApproach::TwoPhase>().approach() == EclTwoPhaseApproach::GasWater)
341 return params.template getRealParams<EclMultiplexerApproach::TwoPhase>().gasWaterParams().SwTrapped();
342 if(params.template getRealParams<EclMultiplexerApproach::TwoPhase>().approach() == EclTwoPhaseApproach::OilWater)
343 return params.template getRealParams<EclMultiplexerApproach::TwoPhase>().oilWaterParams().SwTrapped();
344 return 0.0; // gas-oil case
345 case EclMultiplexerApproach::OnePhase:
346 return 0.0;
347 }
348 return 0.0;
349 }
350 /*
351 * Hysteresis parameters for gas-oil
352 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
353 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
354 * \param params Parameters
355 */
356 static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
357 const Scalar& krnSwMdc,
358 Params& params)
359 {
360 OPM_TIMEFUNCTION_LOCAL();
361 switch (params.approach()) {
362 case EclMultiplexerApproach::Stone1:
363 Stone1Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
364 params.template getRealParams<EclMultiplexerApproach::Stone1>());
365 break;
366
367 case EclMultiplexerApproach::Stone2:
368 Stone2Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
369 params.template getRealParams<EclMultiplexerApproach::Stone2>());
370 break;
371
372 case EclMultiplexerApproach::Default:
373 DefaultMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
374 params.template getRealParams<EclMultiplexerApproach::Default>());
375 break;
376
377 case EclMultiplexerApproach::TwoPhase:
378 TwoPhaseMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
379 params.template getRealParams<EclMultiplexerApproach::TwoPhase>());
380 break;
381
382 case EclMultiplexerApproach::OnePhase:
383 // Do nothing.
384 break;
385 }
386 }
387
397 template <class FluidState, class Evaluation = typename FluidState::Scalar>
398 static Evaluation pcgn(const Params& /* params */,
399 const FluidState& /* fs */)
400 {
401 throw std::logic_error("Not implemented: pcgn()");
402 }
403
413 template <class FluidState, class Evaluation = typename FluidState::Scalar>
414 static Evaluation pcnw(const Params& /* params */,
415 const FluidState& /* fs */)
416 {
417 throw std::logic_error("Not implemented: pcnw()");
418 }
419
423 template <class ContainerT, class FluidState>
424 static void saturations(ContainerT& /* values */,
425 const Params& /* params */,
426 const FluidState& /* fs */)
427 {
428 throw std::logic_error("Not implemented: saturations()");
429 }
430
434 template <class FluidState, class Evaluation = typename FluidState::Scalar>
435 static Evaluation Sg(const Params& /* params */,
436 const FluidState& /* fluidState */)
437 {
438 throw std::logic_error("Not implemented: Sg()");
439 }
440
444 template <class FluidState, class Evaluation = typename FluidState::Scalar>
445 static Evaluation Sn(const Params& /* params */,
446 const FluidState& /* fluidState */)
447 {
448 throw std::logic_error("Not implemented: Sn()");
449 }
450
454 template <class FluidState, class Evaluation = typename FluidState::Scalar>
455 static Evaluation Sw(const Params& /* params */,
456 const FluidState& /* fluidState */)
457 {
458 throw std::logic_error("Not implemented: Sw()");
459 }
460
476 template <class ContainerT, class FluidState>
477 static void relativePermeabilities(ContainerT& values,
478 const Params& params,
479 const FluidState& fluidState)
480 {
481 OPM_TIMEFUNCTION_LOCAL();
482 switch (params.approach()) {
483 case EclMultiplexerApproach::Stone1:
485 params.template getRealParams<EclMultiplexerApproach::Stone1>(),
486 fluidState);
487 break;
488
489 case EclMultiplexerApproach::Stone2:
491 params.template getRealParams<EclMultiplexerApproach::Stone2>(),
492 fluidState);
493 break;
494
495 case EclMultiplexerApproach::Default:
497 params.template getRealParams<EclMultiplexerApproach::Default>(),
498 fluidState);
499 break;
500
501 case EclMultiplexerApproach::TwoPhase:
503 params.template getRealParams<EclMultiplexerApproach::TwoPhase>(),
504 fluidState);
505 break;
506
507 case EclMultiplexerApproach::OnePhase:
508 values[0] = 1.0;
509 break;
510
511 default:
512 throw std::logic_error("Not implemented: relativePermeabilities() option for unknown EclMultiplexerApproach (="
513 + std::to_string(static_cast<int>(params.approach())) + ")");
514 }
515 }
516
520 template <class Evaluation, class FluidState>
521 static Evaluation relpermOilInOilGasSystem(const Params& params,
522 const FluidState& fluidState)
523 {
524 OPM_TIMEFUNCTION_LOCAL();
525 switch (params.approach()) {
526 case EclMultiplexerApproach::Stone1:
527 return Stone1Material::template relpermOilInOilGasSystem<Evaluation>
528 (params.template getRealParams<EclMultiplexerApproach::Stone1>(),
529 fluidState);
530
531 case EclMultiplexerApproach::Stone2:
532 return Stone2Material::template relpermOilInOilGasSystem<Evaluation>
533 (params.template getRealParams<EclMultiplexerApproach::Stone2>(),
534 fluidState);
535
536 case EclMultiplexerApproach::Default:
537 return DefaultMaterial::template relpermOilInOilGasSystem<Evaluation>
538 (params.template getRealParams<EclMultiplexerApproach::Default>(),
539 fluidState);
540
541 default:
542 throw std::logic_error {
543 "relpermOilInOilGasSystem() is specific to three phases"
544 };
545 }
546 }
547
551 template <class Evaluation, class FluidState>
552 static Evaluation relpermOilInOilWaterSystem(const Params& params,
553 const FluidState& fluidState)
554 {
555 OPM_TIMEFUNCTION_LOCAL();
556 switch (params.approach()) {
557 case EclMultiplexerApproach::Stone1:
558 return Stone1Material::template relpermOilInOilWaterSystem<Evaluation>
559 (params.template getRealParams<EclMultiplexerApproach::Stone1>(),
560 fluidState);
561
562 case EclMultiplexerApproach::Stone2:
563 return Stone2Material::template relpermOilInOilWaterSystem<Evaluation>
564 (params.template getRealParams<EclMultiplexerApproach::Stone2>(),
565 fluidState);
566
567 case EclMultiplexerApproach::Default:
568 return DefaultMaterial::template relpermOilInOilWaterSystem<Evaluation>
569 (params.template getRealParams<EclMultiplexerApproach::Default>(),
570 fluidState);
571
572 default:
573 throw std::logic_error {
574 "relpermOilInOilWaterSystem() is specific to three phases"
575 };
576 }
577 }
578
582 template <class FluidState, class Evaluation = typename FluidState::Scalar>
583 static Evaluation krg(const Params& /* params */,
584 const FluidState& /* fluidState */)
585 {
586 throw std::logic_error("Not implemented: krg()");
587 }
588
592 template <class FluidState, class Evaluation = typename FluidState::Scalar>
593 static Evaluation krw(const Params& /* params */,
594 const FluidState& /* fluidState */)
595 {
596 throw std::logic_error("Not implemented: krw()");
597 }
598
602 template <class FluidState, class Evaluation = typename FluidState::Scalar>
603 static Evaluation krn(const Params& /* params */,
604 const FluidState& /* fluidState */)
605 {
606 throw std::logic_error("Not implemented: krn()");
607 }
608
609
617 template <class FluidState>
618 static bool updateHysteresis(Params& params, const FluidState& fluidState)
619 {
620 OPM_TIMEFUNCTION_LOCAL();
621 switch (params.approach()) {
622 case EclMultiplexerApproach::Stone1:
623 return Stone1Material::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::Stone1>(),
624 fluidState);
625 break;
626
627 case EclMultiplexerApproach::Stone2:
628 return Stone2Material::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::Stone2>(),
629 fluidState);
630 break;
631
632 case EclMultiplexerApproach::Default:
633 return DefaultMaterial::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::Default>(),
634 fluidState);
635 break;
636
637 case EclMultiplexerApproach::TwoPhase:
638 return TwoPhaseMaterial::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::TwoPhase>(),
639 fluidState);
640 break;
641 case EclMultiplexerApproach::OnePhase:
642 return false;
643 break;
644 }
645 return false;
646 }
647};
648
649} // namespace Opm
650
651#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: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 bool updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclDefaultMaterial.hpp:439
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition EclDefaultMaterial.hpp:319
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:552
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition EclMultiplexerMaterial.hpp:593
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition EclMultiplexerMaterial.hpp:435
static Evaluation relpermOilInOilGasSystem(const Params &params, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition EclMultiplexerMaterial.hpp:521
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition EclMultiplexerMaterial.hpp:398
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:414
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:477
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition EclMultiplexerMaterial.hpp:445
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition EclMultiplexerMaterial.hpp:583
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:424
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition EclMultiplexerMaterial.hpp:603
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition EclMultiplexerMaterial.hpp:455
static bool updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclMultiplexerMaterial.hpp:618
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 bool 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 bool 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:351
static bool updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition EclTwoPhaseMaterial.hpp:428
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
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30