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 
31 #include "EclDefaultMaterial.hpp"
32 #include "EclStone1Material.hpp"
33 #include "EclStone2Material.hpp"
34 #include "EclTwoPhaseMaterial.hpp"
35 
36 #include <algorithm>
37 #include <stdexcept>
38 
39 namespace Opm {
40 
47 template <class TraitsT,
48  class GasOilMaterialLawT,
49  class OilWaterMaterialLawT,
50  class GasWaterMaterialLawT,
51  class ParamsT = EclMultiplexerMaterialParams<TraitsT,
52  GasOilMaterialLawT,
53  OilWaterMaterialLawT,
54  GasWaterMaterialLawT> >
55 class EclMultiplexerMaterial : public TraitsT
56 {
57 public:
58  using GasOilMaterialLaw = GasOilMaterialLawT;
59  using OilWaterMaterialLaw = OilWaterMaterialLawT;
60  using GasWaterMaterialLaw = GasWaterMaterialLawT;
61 
66 
67  // some safety checks
68  static_assert(TraitsT::numPhases == 3,
69  "The number of phases considered by this capillary pressure "
70  "law is always three!");
71  static_assert(GasOilMaterialLaw::numPhases == 2,
72  "The number of phases considered by the gas-oil capillary "
73  "pressure law must be two!");
74  static_assert(OilWaterMaterialLaw::numPhases == 2,
75  "The number of phases considered by the oil-water capillary "
76  "pressure law must be two!");
77  static_assert(GasWaterMaterialLaw::numPhases == 2,
78  "The number of phases considered by the gas-water capillary "
79  "pressure law must be two!");
80  static_assert(std::is_same<typename GasOilMaterialLaw::Scalar,
81  typename OilWaterMaterialLaw::Scalar>::value,
82  "The two two-phase capillary pressure laws must use the same "
83  "type of floating point values.");
84 
85  using Traits = TraitsT;
86  using Params = ParamsT;
87  using Scalar = typename Traits::Scalar;
88 
89  static constexpr int numPhases = 3;
90  static constexpr int waterPhaseIdx = Traits::wettingPhaseIdx;
91  static constexpr int oilPhaseIdx = Traits::nonWettingPhaseIdx;
92  static constexpr int gasPhaseIdx = Traits::gasPhaseIdx;
93 
96  static constexpr bool implementsTwoPhaseApi = false;
97 
100  static constexpr bool implementsTwoPhaseSatApi = false;
101 
104  static constexpr bool isSaturationDependent = true;
105 
108  static constexpr bool isPressureDependent = false;
109 
112  static constexpr bool isTemperatureDependent = false;
113 
116  static constexpr bool isCompositionDependent = false;
117 
132  template <class ContainerT, class FluidState>
133  static void capillaryPressures(ContainerT& values,
134  const Params& params,
135  const FluidState& fluidState)
136  {
137  switch (params.approach()) {
138  case EclMultiplexerApproach::EclStone1Approach:
140  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
141  fluidState);
142  break;
143 
144  case EclMultiplexerApproach::EclStone2Approach:
146  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
147  fluidState);
148  break;
149 
150  case EclMultiplexerApproach::EclDefaultApproach:
152  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
153  fluidState);
154  break;
155 
156  case EclMultiplexerApproach::EclTwoPhaseApproach:
158  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
159  fluidState);
160  break;
161 
162  case EclMultiplexerApproach::EclOnePhaseApproach:
163  values[0] = 0.0;
164  break;
165  }
166  }
167 
168  /*
169  * Hysteresis parameters for oil-water
170  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
171  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
172  * \param params Parameters
173  */
174  static void oilWaterHysteresisParams(Scalar& pcSwMdc,
175  Scalar& krnSwMdc,
176  const Params& params)
177  {
178  switch (params.approach()) {
179  case EclMultiplexerApproach::EclStone1Approach:
180  Stone1Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
181  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
182  break;
183 
184  case EclMultiplexerApproach::EclStone2Approach:
185  Stone2Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
186  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
187  break;
188 
189  case EclMultiplexerApproach::EclDefaultApproach:
190  DefaultMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
191  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
192  break;
193 
194  case EclMultiplexerApproach::EclTwoPhaseApproach:
195  TwoPhaseMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
196  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
197  break;
198 
199  case EclMultiplexerApproach::EclOnePhaseApproach:
200  // Do nothing.
201  break;
202  }
203  }
204 
205  /*
206  * Hysteresis parameters for oil-water
207  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
208  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
209  * \param params Parameters
210  */
211  static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
212  const Scalar& krnSwMdc,
213  Params& params)
214  {
215  switch (params.approach()) {
216  case EclMultiplexerApproach::EclStone1Approach:
217  Stone1Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
218  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
219  break;
220 
221  case EclMultiplexerApproach::EclStone2Approach:
222  Stone2Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
223  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
224  break;
225 
226  case EclMultiplexerApproach::EclDefaultApproach:
227  DefaultMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
228  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
229  break;
230 
231  case EclMultiplexerApproach::EclTwoPhaseApproach:
232  TwoPhaseMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
233  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
234  break;
235 
236  case EclMultiplexerApproach::EclOnePhaseApproach:
237  // Do nothing.
238  break;
239  }
240  }
241 
242  /*
243  * Hysteresis parameters for gas-oil
244  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
245  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
246  * \param params Parameters
247  */
248  static void gasOilHysteresisParams(Scalar& pcSwMdc,
249  Scalar& krnSwMdc,
250  const Params& params)
251  {
252  switch (params.approach()) {
253  case EclMultiplexerApproach::EclStone1Approach:
254  Stone1Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
255  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
256  break;
257 
258  case EclMultiplexerApproach::EclStone2Approach:
259  Stone2Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
260  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
261  break;
262 
263  case EclMultiplexerApproach::EclDefaultApproach:
264  DefaultMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
265  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
266  break;
267 
268  case EclMultiplexerApproach::EclTwoPhaseApproach:
269  TwoPhaseMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
270  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
271  break;
272 
273  case EclMultiplexerApproach::EclOnePhaseApproach:
274  // Do nothing.
275  break;
276  }
277  }
278 
279  /*
280  * Hysteresis parameters for gas-oil
281  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
282  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
283  * \param params Parameters
284  */
285  static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
286  const Scalar& krnSwMdc,
287  Params& params)
288  {
289  switch (params.approach()) {
290  case EclMultiplexerApproach::EclStone1Approach:
291  Stone1Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
292  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
293  break;
294 
295  case EclMultiplexerApproach::EclStone2Approach:
296  Stone2Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
297  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
298  break;
299 
300  case EclMultiplexerApproach::EclDefaultApproach:
301  DefaultMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
302  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
303  break;
304 
305  case EclMultiplexerApproach::EclTwoPhaseApproach:
306  TwoPhaseMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
307  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
308  break;
309 
310  case EclMultiplexerApproach::EclOnePhaseApproach:
311  // Do nothing.
312  break;
313  }
314  }
315 
325  template <class FluidState, class Evaluation = typename FluidState::Scalar>
326  static Evaluation pcgn(const Params& /* params */,
327  const FluidState& /* fs */)
328  {
329  throw std::logic_error("Not implemented: pcgn()");
330  }
331 
341  template <class FluidState, class Evaluation = typename FluidState::Scalar>
342  static Evaluation pcnw(const Params& /* params */,
343  const FluidState& /* fs */)
344  {
345  throw std::logic_error("Not implemented: pcnw()");
346  }
347 
351  template <class ContainerT, class FluidState>
352  static void saturations(ContainerT& /* values */,
353  const Params& /* params */,
354  const FluidState& /* fs */)
355  {
356  throw std::logic_error("Not implemented: saturations()");
357  }
358 
362  template <class FluidState, class Evaluation = typename FluidState::Scalar>
363  static Evaluation Sg(const Params& /* params */,
364  const FluidState& /* fluidState */)
365  {
366  throw std::logic_error("Not implemented: Sg()");
367  }
368 
372  template <class FluidState, class Evaluation = typename FluidState::Scalar>
373  static Evaluation Sn(const Params& /* params */,
374  const FluidState& /* fluidState */)
375  {
376  throw std::logic_error("Not implemented: Sn()");
377  }
378 
382  template <class FluidState, class Evaluation = typename FluidState::Scalar>
383  static Evaluation Sw(const Params& /* params */,
384  const FluidState& /* fluidState */)
385  {
386  throw std::logic_error("Not implemented: Sw()");
387  }
388 
404  template <class ContainerT, class FluidState>
405  static void relativePermeabilities(ContainerT& values,
406  const Params& params,
407  const FluidState& fluidState)
408  {
409  switch (params.approach()) {
410  case EclMultiplexerApproach::EclStone1Approach:
412  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
413  fluidState);
414  break;
415 
416  case EclMultiplexerApproach::EclStone2Approach:
418  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
419  fluidState);
420  break;
421 
422  case EclMultiplexerApproach::EclDefaultApproach:
424  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
425  fluidState);
426  break;
427 
428  case EclMultiplexerApproach::EclTwoPhaseApproach:
430  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
431  fluidState);
432  break;
433 
434  case EclMultiplexerApproach::EclOnePhaseApproach:
435  values[0] = 1.0;
436  break;
437 
438  default:
439  throw std::logic_error("Not implemented: relativePermeabilities() option for unknown EclMultiplexerApproach (="
440  + std::to_string(static_cast<int>(params.approach())) + ")");
441  }
442  }
443 
447  template <class Evaluation, class FluidState>
448  static Evaluation relpermOilInOilGasSystem(const Params& params,
449  const FluidState& fluidState)
450  {
451  switch (params.approach()) {
452  case EclMultiplexerApproach::EclStone1Approach:
453  return Stone1Material::template relpermOilInOilGasSystem<Evaluation>
454  (params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
455  fluidState);
456 
457  case EclMultiplexerApproach::EclStone2Approach:
458  return Stone2Material::template relpermOilInOilGasSystem<Evaluation>
459  (params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
460  fluidState);
461 
462  case EclMultiplexerApproach::EclDefaultApproach:
463  return DefaultMaterial::template relpermOilInOilGasSystem<Evaluation>
464  (params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
465  fluidState);
466 
467  default:
468  throw std::logic_error {
469  "relpermOilInOilGasSystem() is specific to three phases"
470  };
471  }
472  }
473 
477  template <class Evaluation, class FluidState>
478  static Evaluation relpermOilInOilWaterSystem(const Params& params,
479  const FluidState& fluidState)
480  {
481  switch (params.approach()) {
482  case EclMultiplexerApproach::EclStone1Approach:
483  return Stone1Material::template relpermOilInOilWaterSystem<Evaluation>
484  (params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
485  fluidState);
486 
487  case EclMultiplexerApproach::EclStone2Approach:
488  return Stone2Material::template relpermOilInOilWaterSystem<Evaluation>
489  (params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
490  fluidState);
491 
492  case EclMultiplexerApproach::EclDefaultApproach:
493  return DefaultMaterial::template relpermOilInOilWaterSystem<Evaluation>
494  (params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
495  fluidState);
496 
497  default:
498  throw std::logic_error {
499  "relpermOilInOilWaterSystem() is specific to three phases"
500  };
501  }
502  }
503 
507  template <class FluidState, class Evaluation = typename FluidState::Scalar>
508  static Evaluation krg(const Params& /* params */,
509  const FluidState& /* fluidState */)
510  {
511  throw std::logic_error("Not implemented: krg()");
512  }
513 
517  template <class FluidState, class Evaluation = typename FluidState::Scalar>
518  static Evaluation krw(const Params& /* params */,
519  const FluidState& /* fluidState */)
520  {
521  throw std::logic_error("Not implemented: krw()");
522  }
523 
527  template <class FluidState, class Evaluation = typename FluidState::Scalar>
528  static Evaluation krn(const Params& /* params */,
529  const FluidState& /* fluidState */)
530  {
531  throw std::logic_error("Not implemented: krn()");
532  }
533 
534 
542  template <class FluidState>
543  static void updateHysteresis(Params& params, const FluidState& fluidState)
544  {
545  switch (params.approach()) {
546  case EclMultiplexerApproach::EclStone1Approach:
547  Stone1Material::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
548  fluidState);
549  break;
550 
551  case EclMultiplexerApproach::EclStone2Approach:
552  Stone2Material::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
553  fluidState);
554  break;
555 
556  case EclMultiplexerApproach::EclDefaultApproach:
557  DefaultMaterial::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
558  fluidState);
559  break;
560 
561  case EclMultiplexerApproach::EclTwoPhaseApproach:
562  TwoPhaseMaterial::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
563  fluidState);
564  break;
565  case EclMultiplexerApproach::EclOnePhaseApproach:
566  break;
567  }
568  }
569 };
570 
571 } // namespace Opm
572 
573 #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 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:56
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:133
static Evaluation relpermOilInOilWaterSystem(const Params &params, const FluidState &fluidState)
The relative permeability of oil in oil/water system.
Definition: EclMultiplexerMaterial.hpp:478
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition: EclMultiplexerMaterial.hpp:518
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition: EclMultiplexerMaterial.hpp:363
static Evaluation relpermOilInOilGasSystem(const Params &params, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition: EclMultiplexerMaterial.hpp:448
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:326
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:342
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclMultiplexerMaterial.hpp:543
static constexpr bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: EclMultiplexerMaterial.hpp:108
static constexpr bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclMultiplexerMaterial.hpp:100
static constexpr bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclMultiplexerMaterial.hpp:104
static constexpr bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclMultiplexerMaterial.hpp:112
static constexpr bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclMultiplexerMaterial.hpp:96
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclMultiplexerMaterial.hpp:405
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:373
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition: EclMultiplexerMaterial.hpp:508
static constexpr bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: EclMultiplexerMaterial.hpp:116
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: EclMultiplexerMaterial.hpp:352
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:528
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition: EclMultiplexerMaterial.hpp:383
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