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 
38 
39 #include <algorithm>
40 
41 namespace Opm {
42 
49 template <class TraitsT,
50  class GasOilMaterialLawT,
51  class OilWaterMaterialLawT,
52  class GasWaterMaterialLawT,
53  class ParamsT = EclMultiplexerMaterialParams<TraitsT,
54  GasOilMaterialLawT,
55  OilWaterMaterialLawT,
56  GasWaterMaterialLawT> >
57 class EclMultiplexerMaterial : public TraitsT
58 {
59 public:
60  typedef GasOilMaterialLawT GasOilMaterialLaw;
61  typedef OilWaterMaterialLawT OilWaterMaterialLaw;
62  typedef GasWaterMaterialLawT GasWaterMaterialLaw;
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  typedef TraitsT Traits;
88  typedef ParamsT Params;
89  typedef typename Traits::Scalar Scalar;
90 
91  static const int numPhases = 3;
92  static const int waterPhaseIdx = Traits::wettingPhaseIdx;
93  static const int oilPhaseIdx = Traits::nonWettingPhaseIdx;
94  static const int gasPhaseIdx = Traits::gasPhaseIdx;
95 
98  static const bool implementsTwoPhaseApi = false;
99 
102  static const bool implementsTwoPhaseSatApi = false;
103 
106  static const bool isSaturationDependent = true;
107 
110  static const bool isPressureDependent = false;
111 
114  static const bool isTemperatureDependent = false;
115 
118  static const 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  switch (params.approach()) {
140  case EclMultiplexerApproach::EclStone1Approach:
142  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
143  fluidState);
144  break;
145 
146  case EclMultiplexerApproach::EclStone2Approach:
148  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
149  fluidState);
150  break;
151 
152  case EclMultiplexerApproach::EclDefaultApproach:
154  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
155  fluidState);
156  break;
157 
158  case EclMultiplexerApproach::EclTwoPhaseApproach:
160  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
161  fluidState);
162  break;
163 
164  case EclMultiplexerApproach::EclOnePhaseApproach:
165  values[0] = 0.0;
166  break;
167  }
168  }
169 
170  /*
171  * Hysteresis parameters for oil-water
172  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
173  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
174  * \param params Parameters
175  */
176  static void oilWaterHysteresisParams(Scalar& pcSwMdc,
177  Scalar& krnSwMdc,
178  const Params& params)
179  {
180  switch (params.approach()) {
181  case EclMultiplexerApproach::EclStone1Approach:
182  Stone1Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
183  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
184  break;
185 
186  case EclMultiplexerApproach::EclStone2Approach:
187  Stone2Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
188  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
189  break;
190 
191  case EclMultiplexerApproach::EclDefaultApproach:
192  DefaultMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
193  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
194  break;
195 
196  case EclMultiplexerApproach::EclTwoPhaseApproach:
197  TwoPhaseMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
198  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
199  break;
200 
201  case EclMultiplexerApproach::EclOnePhaseApproach:
202  // Do nothing.
203  break;
204  }
205  }
206 
207  /*
208  * Hysteresis parameters for oil-water
209  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
210  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
211  * \param params Parameters
212  */
213  static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
214  const Scalar& krnSwMdc,
215  Params& params)
216  {
217  switch (params.approach()) {
218  case EclMultiplexerApproach::EclStone1Approach:
219  Stone1Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
220  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
221  break;
222 
223  case EclMultiplexerApproach::EclStone2Approach:
224  Stone2Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
225  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
226  break;
227 
228  case EclMultiplexerApproach::EclDefaultApproach:
229  DefaultMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
230  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
231  break;
232 
233  case EclMultiplexerApproach::EclTwoPhaseApproach:
234  TwoPhaseMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
235  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
236  break;
237 
238  case EclMultiplexerApproach::EclOnePhaseApproach:
239  // Do nothing.
240  break;
241  }
242  }
243 
244  /*
245  * Hysteresis parameters for gas-oil
246  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
247  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
248  * \param params Parameters
249  */
250  static void gasOilHysteresisParams(Scalar& pcSwMdc,
251  Scalar& krnSwMdc,
252  const Params& params)
253  {
254  switch (params.approach()) {
255  case EclMultiplexerApproach::EclStone1Approach:
256  Stone1Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
257  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
258  break;
259 
260  case EclMultiplexerApproach::EclStone2Approach:
261  Stone2Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
262  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
263  break;
264 
265  case EclMultiplexerApproach::EclDefaultApproach:
266  DefaultMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
267  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
268  break;
269 
270  case EclMultiplexerApproach::EclTwoPhaseApproach:
271  TwoPhaseMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
272  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
273  break;
274 
275  case EclMultiplexerApproach::EclOnePhaseApproach:
276  // Do nothing.
277  break;
278  }
279  }
280 
281  /*
282  * Hysteresis parameters for gas-oil
283  * @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
284  * @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
285  * \param params Parameters
286  */
287  static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
288  const Scalar& krnSwMdc,
289  Params& params)
290  {
291  switch (params.approach()) {
292  case EclMultiplexerApproach::EclStone1Approach:
293  Stone1Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
294  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>());
295  break;
296 
297  case EclMultiplexerApproach::EclStone2Approach:
298  Stone2Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
299  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>());
300  break;
301 
302  case EclMultiplexerApproach::EclDefaultApproach:
303  DefaultMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
304  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>());
305  break;
306 
307  case EclMultiplexerApproach::EclTwoPhaseApproach:
308  TwoPhaseMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
309  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>());
310  break;
311 
312  case EclMultiplexerApproach::EclOnePhaseApproach:
313  // Do nothing.
314  break;
315  }
316  }
317 
327  template <class FluidState, class Evaluation = typename FluidState::Scalar>
328  static Evaluation pcgn(const Params& /* params */,
329  const FluidState& /* fs */)
330  {
331  throw std::logic_error("Not implemented: pcgn()");
332  }
333 
343  template <class FluidState, class Evaluation = typename FluidState::Scalar>
344  static Evaluation pcnw(const Params& /* params */,
345  const FluidState& /* fs */)
346  {
347  throw std::logic_error("Not implemented: pcnw()");
348  }
349 
353  template <class ContainerT, class FluidState>
354  static void saturations(ContainerT& /* values */,
355  const Params& /* params */,
356  const FluidState& /* fs */)
357  {
358  throw std::logic_error("Not implemented: saturations()");
359  }
360 
364  template <class FluidState, class Evaluation = typename FluidState::Scalar>
365  static Evaluation Sg(const Params& /* params */,
366  const FluidState& /* fluidState */)
367  {
368  throw std::logic_error("Not implemented: Sg()");
369  }
370 
374  template <class FluidState, class Evaluation = typename FluidState::Scalar>
375  static Evaluation Sn(const Params& /* params */,
376  const FluidState& /* fluidState */)
377  {
378  throw std::logic_error("Not implemented: Sn()");
379  }
380 
384  template <class FluidState, class Evaluation = typename FluidState::Scalar>
385  static Evaluation Sw(const Params& /* params */,
386  const FluidState& /* fluidState */)
387  {
388  throw std::logic_error("Not implemented: Sw()");
389  }
390 
406  template <class ContainerT, class FluidState>
407  static void relativePermeabilities(ContainerT& values,
408  const Params& params,
409  const FluidState& fluidState)
410  {
411  switch (params.approach()) {
412  case EclMultiplexerApproach::EclStone1Approach:
414  params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
415  fluidState);
416  break;
417 
418  case EclMultiplexerApproach::EclStone2Approach:
420  params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
421  fluidState);
422  break;
423 
424  case EclMultiplexerApproach::EclDefaultApproach:
426  params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
427  fluidState);
428  break;
429 
430  case EclMultiplexerApproach::EclTwoPhaseApproach:
432  params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
433  fluidState);
434  break;
435 
436  case EclMultiplexerApproach::EclOnePhaseApproach:
437  values[0] = 1.0;
438  break;
439 
440  default:
441  throw std::logic_error("Not implemented: relativePermeabilities() option for unknown EclMultiplexerApproach (="
442  + std::to_string(static_cast<int>(params.approach())) + ")");
443  }
444  }
445 
449  template <class Evaluation, class FluidState>
450  static Evaluation relpermOilInOilGasSystem(const Params& params,
451  const FluidState& fluidState)
452  {
453  switch (params.approach()) {
454  case EclMultiplexerApproach::EclStone1Approach:
455  return Stone1Material::template relpermOilInOilGasSystem<Evaluation>
456  (params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
457  fluidState);
458 
459  case EclMultiplexerApproach::EclStone2Approach:
460  return Stone2Material::template relpermOilInOilGasSystem<Evaluation>
461  (params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
462  fluidState);
463 
464  case EclMultiplexerApproach::EclDefaultApproach:
465  return DefaultMaterial::template relpermOilInOilGasSystem<Evaluation>
466  (params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
467  fluidState);
468 
469  default:
470  throw std::logic_error {
471  "relpermOilInOilGasSystem() is specific to three phases"
472  };
473  }
474  }
475 
479  template <class Evaluation, class FluidState>
480  static Evaluation relpermOilInOilWaterSystem(const Params& params,
481  const FluidState& fluidState)
482  {
483  switch (params.approach()) {
484  case EclMultiplexerApproach::EclStone1Approach:
485  return Stone1Material::template relpermOilInOilWaterSystem<Evaluation>
486  (params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
487  fluidState);
488 
489  case EclMultiplexerApproach::EclStone2Approach:
490  return Stone2Material::template relpermOilInOilWaterSystem<Evaluation>
491  (params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
492  fluidState);
493 
494  case EclMultiplexerApproach::EclDefaultApproach:
495  return DefaultMaterial::template relpermOilInOilWaterSystem<Evaluation>
496  (params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
497  fluidState);
498 
499  default:
500  throw std::logic_error {
501  "relpermOilInOilWaterSystem() is specific to three phases"
502  };
503  }
504  }
505 
509  template <class FluidState, class Evaluation = typename FluidState::Scalar>
510  static Evaluation krg(const Params& /* params */,
511  const FluidState& /* fluidState */)
512  {
513  throw std::logic_error("Not implemented: krg()");
514  }
515 
519  template <class FluidState, class Evaluation = typename FluidState::Scalar>
520  static Evaluation krw(const Params& /* params */,
521  const FluidState& /* fluidState */)
522  {
523  throw std::logic_error("Not implemented: krw()");
524  }
525 
529  template <class FluidState, class Evaluation = typename FluidState::Scalar>
530  static Evaluation krn(const Params& /* params */,
531  const FluidState& /* fluidState */)
532  {
533  throw std::logic_error("Not implemented: krn()");
534  }
535 
536 
544  template <class FluidState>
545  static void updateHysteresis(Params& params, const FluidState& fluidState)
546  {
547  switch (params.approach()) {
548  case EclMultiplexerApproach::EclStone1Approach:
549  Stone1Material::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::EclStone1Approach>(),
550  fluidState);
551  break;
552 
553  case EclMultiplexerApproach::EclStone2Approach:
554  Stone2Material::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::EclStone2Approach>(),
555  fluidState);
556  break;
557 
558  case EclMultiplexerApproach::EclDefaultApproach:
559  DefaultMaterial::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::EclDefaultApproach>(),
560  fluidState);
561  break;
562 
563  case EclMultiplexerApproach::EclTwoPhaseApproach:
564  TwoPhaseMaterial::updateHysteresis(params.template getRealParams<EclMultiplexerApproach::EclTwoPhaseApproach>(),
565  fluidState);
566  break;
567  case EclMultiplexerApproach::EclOnePhaseApproach:
568  break;
569  }
570  }
571 };
572 } // namespace Opm
573 
574 #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.
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
Some templates to wrap the valgrind client request macros.
Implements the default three phase capillary pressure law used by the ECLipse simulator.
Definition: EclDefaultMaterial.hpp:59
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:134
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclDefaultMaterial.hpp:308
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclDefaultMaterial.hpp:422
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:480
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: EclMultiplexerMaterial.hpp:98
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: EclMultiplexerMaterial.hpp:106
static Evaluation krw(const Params &, const FluidState &)
The relative permeability of the wetting phase.
Definition: EclMultiplexerMaterial.hpp:520
static Evaluation Sg(const Params &, const FluidState &)
The saturation of the gas phase.
Definition: EclMultiplexerMaterial.hpp:365
static Evaluation relpermOilInOilGasSystem(const Params &params, const FluidState &fluidState)
The relative permeability of oil in oil/gas system.
Definition: EclMultiplexerMaterial.hpp:450
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: EclMultiplexerMaterial.hpp:110
static Evaluation pcgn(const Params &, const FluidState &)
Capillary pressure between the gas and the non-wetting liquid (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:328
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:344
static void updateHysteresis(Params &params, const FluidState &fluidState)
Update the hysteresis parameters after a time step.
Definition: EclMultiplexerMaterial.hpp:545
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: EclMultiplexerMaterial.hpp:118
static void relativePermeabilities(ContainerT &values, const Params &params, const FluidState &fluidState)
The relative permeability of all phases.
Definition: EclMultiplexerMaterial.hpp:407
static Evaluation Sn(const Params &, const FluidState &)
The saturation of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:375
static Evaluation krg(const Params &, const FluidState &)
The relative permeability of the gas phase.
Definition: EclMultiplexerMaterial.hpp:510
static void saturations(ContainerT &, const Params &, const FluidState &)
The inverse of the capillary pressure.
Definition: EclMultiplexerMaterial.hpp:354
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: EclMultiplexerMaterial.hpp:102
static Evaluation krn(const Params &, const FluidState &)
The relative permeability of the non-wetting (i.e., oil) phase.
Definition: EclMultiplexerMaterial.hpp:530
static Evaluation Sw(const Params &, const FluidState &)
The saturation of the wetting (i.e., water) phase.
Definition: EclMultiplexerMaterial.hpp:385
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: EclMultiplexerMaterial.hpp:114
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