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
34
35#include <algorithm>
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::pcSwMdc(...)
186 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
187 * \param params Parameters
188 */
189 static void oilWaterHysteresisParams(Scalar& pcSwMdc,
190 Scalar& krnSwMdc,
191 const Params& params)
192 {
193 OPM_TIMEFUNCTION_LOCAL();
194 pcSwMdc = params.oilWaterParams().pcSwMdc();
195 krnSwMdc = params.oilWaterParams().krnSwMdc();
196
197 Valgrind::CheckDefined(pcSwMdc);
198 Valgrind::CheckDefined(krnSwMdc);
199 }
200
201 /*
202 * Hysteresis parameters for oil-water
203 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
204 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
205 * \param params Parameters
206 */
207 static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
208 const Scalar& krnSwMdc,
209 Params& params)
210 {
211 OPM_TIMEFUNCTION_LOCAL();
212 constexpr const Scalar krwSw = 2.0; //Should not be used
213 params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
214 }
215
216 /*
217 * Hysteresis parameters for gas-oil
218 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
219 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
220 * \param params Parameters
221 */
222 static void gasOilHysteresisParams(Scalar& pcSwMdc,
223 Scalar& krnSwMdc,
224 const Params& params)
225 {
226 pcSwMdc = params.gasOilParams().pcSwMdc();
227 krnSwMdc = params.gasOilParams().krnSwMdc();
228
229 Valgrind::CheckDefined(pcSwMdc);
230 Valgrind::CheckDefined(krnSwMdc);
231 }
232
233 /*
234 * Hysteresis parameters for gas-oil
235 * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
236 * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
237 * \param params Parameters
238 */
239 static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
240 const Scalar& krnSwMdc,
241 Params& params)
242 {
243 constexpr const Scalar krwSw = 2.0; //Should not be used
244 params.gasOilParams().update(pcSwMdc, krwSw, krnSwMdc);
245 }
246
247 static Scalar trappedGasSaturation(const Params& params){
248 if(params.approach() == EclTwoPhaseApproach::GasOil)
249 return params.gasOilParams().SnTrapped();
250 if(params.approach() == EclTwoPhaseApproach::GasWater)
251 return params.gasWaterParams().SnTrapped();
252 return 0.0; // oil-water case
253 }
254
255 static Scalar trappedOilSaturation(const Params& params){
256 if(params.approach() == EclTwoPhaseApproach::GasOil)
257 return params.gasOilParams().SwTrapped();
258 if(params.approach() == EclTwoPhaseApproach::OilWater)
259 return params.oilWaterParams().SnTrapped();
260 return 0.0; // gas-water case
261 }
271 template <class FluidState, class Evaluation = typename FluidState::Scalar>
272 static Evaluation pcgn(const Params& /* params */,
273 const FluidState& /* fs */)
274 {
275 throw std::logic_error("Not implemented: pcgn()");
276 }
277
287 template <class FluidState, class Evaluation = typename FluidState::Scalar>
288 static Evaluation pcnw(const Params& /* params */,
289 const FluidState& /* fs */)
290 {
291 throw std::logic_error("Not implemented: pcnw()");
292 }
293
297 template <class ContainerT, class FluidState>
298 static void saturations(ContainerT& /* values */,
299 const Params& /* params */,
300 const FluidState& /* fs */)
301 {
302 throw std::logic_error("Not implemented: saturations()");
303 }
304
308 template <class FluidState, class Evaluation = typename FluidState::Scalar>
309 static Evaluation Sg(const Params& /* params */,
310 const FluidState& /* fluidState */)
311 {
312 throw std::logic_error("Not implemented: Sg()");
313 }
314
318 template <class FluidState, class Evaluation = typename FluidState::Scalar>
319 static Evaluation Sn(const Params& /* params */,
320 const FluidState& /* fluidState */)
321 {
322 throw std::logic_error("Not implemented: Sn()");
323 }
324
328 template <class FluidState, class Evaluation = typename FluidState::Scalar>
329 static Evaluation Sw(const Params& /* params */,
330 const FluidState& /* fluidState */)
331 {
332 throw std::logic_error("Not implemented: Sw()");
333 }
334
350 template <class ContainerT, class FluidState>
351 static void relativePermeabilities(ContainerT& values,
352 const Params& params,
353 const FluidState& fluidState)
354 {
355 OPM_TIMEFUNCTION_LOCAL();
356 using Evaluation = typename std::remove_reference<decltype(values[0])>::type;
357
358 switch (params.approach()) {
359 case EclTwoPhaseApproach::GasOil: {
360 const Evaluation& So =
361 decay<Evaluation>(fluidState.saturation(oilPhaseIdx));
362
363 values[oilPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrw(params.gasOilParams(), So);
364 values[gasPhaseIdx] = GasOilMaterialLaw::twoPhaseSatKrn(params.gasOilParams(), So);
365 break;
366 }
367
368 case EclTwoPhaseApproach::OilWater: {
369 const Evaluation& Sw =
370 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
371
372 values[waterPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrw(params.oilWaterParams(), Sw);
373 values[oilPhaseIdx] = OilWaterMaterialLaw::twoPhaseSatKrn(params.oilWaterParams(), Sw);
374 break;
375 }
376
377 case EclTwoPhaseApproach::GasWater: {
378 const Evaluation& Sw =
379 decay<Evaluation>(fluidState.saturation(waterPhaseIdx));
380
381 values[waterPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatKrw(params.gasWaterParams(), Sw);
382 values[gasPhaseIdx] = GasWaterMaterialLaw::twoPhaseSatKrn(params.gasWaterParams(), Sw);
383
384 break;
385 }
386 }
387 }
388
392 template <class FluidState, class Evaluation = typename FluidState::Scalar>
393 static Evaluation krg(const Params& /* params */,
394 const FluidState& /* fluidState */)
395 {
396 throw std::logic_error("Not implemented: krg()");
397 }
398
402 template <class FluidState, class Evaluation = typename FluidState::Scalar>
403 static Evaluation krw(const Params& /* params */,
404 const FluidState& /* fluidState */)
405 {
406 throw std::logic_error("Not implemented: krw()");
407 }
408
412 template <class FluidState, class Evaluation = typename FluidState::Scalar>
413 static Evaluation krn(const Params& /* params */,
414 const FluidState& /* fluidState */)
415 {
416 throw std::logic_error("Not implemented: krn()");
417 }
418
419
427 template <class FluidState>
428 static bool updateHysteresis(Params& params, const FluidState& fluidState)
429 {
430 OPM_TIMEFUNCTION_LOCAL();
431 switch (params.approach()) {
432 case EclTwoPhaseApproach::GasOil: {
433 Scalar So = scalarValue(fluidState.saturation(oilPhaseIdx));
434
435 return params.gasOilParams().update(/*pcSw=*/So, /*krwSw=*/So, /*krnSw=*/So);
436 break;
437 }
438
439 case EclTwoPhaseApproach::OilWater: {
440 Scalar Sw = scalarValue(fluidState.saturation(waterPhaseIdx));
441
442 return params.oilWaterParams().update(/*pcSw=*/Sw, /*krwSw=*/Sw, /*krnSw=*/Sw);
443 break;
444 }
445
446 case EclTwoPhaseApproach::GasWater: {
447 Scalar Sw = scalarValue(fluidState.saturation(waterPhaseIdx));
448
449 return params.gasWaterParams().update(/*pcSw=*/Sw, /*krwSw=*/Sw, /*krnSw=*/Sw);
450 break;
451 }
452 }
453
454 // Should not get here...
455 return false;
456 }
457};
458
459} // namespace Opm
460
461#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:351
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:329
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition EclTwoPhaseMaterial.hpp:272
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition EclTwoPhaseMaterial.hpp:298
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition EclTwoPhaseMaterial.hpp:413
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition EclTwoPhaseMaterial.hpp:319
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:393
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition EclTwoPhaseMaterial.hpp:309
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
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:288
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition EclTwoPhaseMaterial.hpp:403
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