My Project
BlackOilFluidSystem.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_BLACK_OIL_FLUID_SYSTEM_HPP
28 #define OPM_BLACK_OIL_FLUID_SYSTEM_HPP
29 
35 
38 
43 
44 #if HAVE_ECL_INPUT
45 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
46 #include <opm/input/eclipse/EclipseState/Tables/FlatTable.hpp>
47 #include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
48 #endif
49 
50 #include <memory>
51 #include <vector>
52 #include <array>
53 
54 namespace Opm {
55 namespace BlackOil {
56 OPM_GENERATE_HAS_MEMBER(Rs, ) // Creates 'HasMember_Rs<T>'.
57 OPM_GENERATE_HAS_MEMBER(Rv, ) // Creates 'HasMember_Rv<T>'.
58 OPM_GENERATE_HAS_MEMBER(Rvw, ) // Creates 'HasMember_Rvw<T>'.
59 OPM_GENERATE_HAS_MEMBER(saltConcentration, )
60 OPM_GENERATE_HAS_MEMBER(saltSaturation, )
61 
62 template <class FluidSystem, class FluidState, class LhsEval>
63 LhsEval getRs_(typename std::enable_if<!HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
64  unsigned regionIdx)
65 {
66  const auto& XoG =
67  decay<LhsEval>(fluidState.massFraction(FluidSystem::oilPhaseIdx, FluidSystem::gasCompIdx));
68  return FluidSystem::convertXoGToRs(XoG, regionIdx);
69 }
70 
71 template <class FluidSystem, class FluidState, class LhsEval>
72 auto getRs_(typename std::enable_if<HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
73  unsigned)
74  -> decltype(decay<LhsEval>(fluidState.Rs()))
75 { return decay<LhsEval>(fluidState.Rs()); }
76 
77 template <class FluidSystem, class FluidState, class LhsEval>
78 LhsEval getRv_(typename std::enable_if<!HasMember_Rv<FluidState>::value, const FluidState&>::type fluidState,
79  unsigned regionIdx)
80 {
81  const auto& XgO =
82  decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::oilCompIdx));
83  return FluidSystem::convertXgOToRv(XgO, regionIdx);
84 }
85 
86 template <class FluidSystem, class FluidState, class LhsEval>
87 auto getRv_(typename std::enable_if<HasMember_Rv<FluidState>::value, const FluidState&>::type fluidState,
88  unsigned)
89  -> decltype(decay<LhsEval>(fluidState.Rv()))
90 { return decay<LhsEval>(fluidState.Rv()); }
91 
92 template <class FluidSystem, class FluidState, class LhsEval>
93 LhsEval getRvw_(typename std::enable_if<!HasMember_Rvw<FluidState>::value, const FluidState&>::type fluidState,
94  unsigned regionIdx)
95 {
96  const auto& XgW =
97  decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::waterCompIdx));
98  return FluidSystem::convertXgWToRvw(XgW, regionIdx);
99 }
100 
101 template <class FluidSystem, class FluidState, class LhsEval>
102 auto getRvw_(typename std::enable_if<HasMember_Rvw<FluidState>::value, const FluidState&>::type fluidState,
103  unsigned)
104  -> decltype(decay<LhsEval>(fluidState.Rvw()))
105 { return decay<LhsEval>(fluidState.Rvw()); }
106 
107 template <class FluidSystem, class FluidState, class LhsEval>
108 LhsEval getSaltConcentration_(typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
109  const FluidState&>::type,
110  unsigned)
111 {return 0.0;}
112 
113 template <class FluidSystem, class FluidState, class LhsEval>
114 auto getSaltConcentration_(typename std::enable_if<HasMember_saltConcentration<FluidState>::value, const FluidState&>::type fluidState,
115  unsigned)
116  -> decltype(decay<LhsEval>(fluidState.saltConcentration()))
117 { return decay<LhsEval>(fluidState.saltConcentration()); }
118 
119 template <class FluidSystem, class FluidState, class LhsEval>
120 LhsEval getSaltSaturation_(typename std::enable_if<!HasMember_saltSaturation<FluidState>::value,
121  const FluidState&>::type,
122  unsigned)
123 {return 0.0;}
124 
125 template <class FluidSystem, class FluidState, class LhsEval>
126 auto getSaltSaturation_(typename std::enable_if<HasMember_saltSaturation<FluidState>::value, const FluidState&>::type fluidState,
127  unsigned)
128  -> decltype(decay<LhsEval>(fluidState.saltSaturation()))
129 { return decay<LhsEval>(fluidState.saltSaturation()); }
130 
131 }
132 
139 template <class Scalar, class IndexTraits = BlackOilDefaultIndexTraits>
140 class BlackOilFluidSystem : public BaseFluidSystem<Scalar, BlackOilFluidSystem<Scalar, IndexTraits> >
141 {
143 
144 public:
148 
150  template <class EvaluationT>
151  struct ParameterCache : public NullParameterCache<EvaluationT>
152  {
153  typedef EvaluationT Evaluation;
154 
155  public:
156  ParameterCache(Scalar maxOilSat = 1.0, unsigned regionIdx=0)
157  {
158  maxOilSat_ = maxOilSat;
159  regionIdx_ = regionIdx;
160  }
161 
169  template <class OtherCache>
170  void assignPersistentData(const OtherCache& other)
171  {
172  regionIdx_ = other.regionIndex();
173  maxOilSat_ = other.maxOilSat();
174  }
175 
183  unsigned regionIndex() const
184  { return regionIdx_; }
185 
193  void setRegionIndex(unsigned val)
194  { regionIdx_ = val; }
195 
196  const Evaluation& maxOilSat() const
197  { return maxOilSat_; }
198 
199  void setMaxOilSat(const Evaluation& val)
200  { maxOilSat_ = val; }
201 
202  private:
203  Evaluation maxOilSat_;
204  unsigned regionIdx_;
205  };
206 
207  /****************************************
208  * Initialization
209  ****************************************/
210 #if HAVE_ECL_INPUT
214  static void initFromState(const EclipseState& eclState, const Schedule& schedule)
215  {
216  size_t numRegions = eclState.runspec().tabdims().getNumPVTTables();
218 
219  numActivePhases_ = 0;
220  std::fill_n(&phaseIsActive_[0], numPhases, false);
221 
222 
223  if (eclState.runspec().phases().active(Phase::OIL)) {
224  phaseIsActive_[oilPhaseIdx] = true;
225  ++ numActivePhases_;
226  }
227 
228  if (eclState.runspec().phases().active(Phase::GAS)) {
229  phaseIsActive_[gasPhaseIdx] = true;
230  ++ numActivePhases_;
231  }
232 
233  if (eclState.runspec().phases().active(Phase::WATER)) {
234  phaseIsActive_[waterPhaseIdx] = true;
235  ++ numActivePhases_;
236  }
237 
238  // set the surface conditions using the STCOND keyword
239  surfaceTemperature = eclState.getTableManager().stCond().temperature;
240  surfacePressure = eclState.getTableManager().stCond().pressure;
241 
242  // The reservoir temperature does not really belong into the table manager. TODO:
243  // change this in opm-parser
244  setReservoirTemperature(eclState.getTableManager().rtemp());
245 
246  // this fluidsystem only supports two or three phases
247  assert(numActivePhases_ >= 1 && numActivePhases_ <= 3);
248 
249  setEnableDissolvedGas(eclState.getSimulationConfig().hasDISGAS());
250  setEnableVaporizedOil(eclState.getSimulationConfig().hasVAPOIL());
251  setEnableVaporizedWater(eclState.getSimulationConfig().hasVAPWAT());
252 
253  if (phaseIsActive(gasPhaseIdx)) {
254  gasPvt_ = std::make_shared<GasPvt>();
255  gasPvt_->initFromState(eclState, schedule);
256  }
257 
258  if (phaseIsActive(oilPhaseIdx)) {
259  oilPvt_ = std::make_shared<OilPvt>();
260  oilPvt_->initFromState(eclState, schedule);
261  }
262 
264  waterPvt_ = std::make_shared<WaterPvt>();
265  waterPvt_->initFromState(eclState, schedule);
266  }
267 
268  // set the reference densities of all PVT regions
269  for (unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
270  setReferenceDensities(phaseIsActive(oilPhaseIdx)? oilPvt_->oilReferenceDensity(regionIdx):700.,
271  phaseIsActive(waterPhaseIdx)? waterPvt_->waterReferenceDensity(regionIdx):1000.,
272  phaseIsActive(gasPhaseIdx)? gasPvt_->gasReferenceDensity(regionIdx):2.,
273  regionIdx);
274  }
275 
276  // set default molarMass and mappings
277  initEnd();
278 
279  // use molarMass of CO2 and Brine as default
280  // when we are using the the CO2STORE option
281  // NB the oil component is used internally for
282  // brine
283  if (eclState.runspec().co2Storage()) {
284  for (unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
285  molarMass_[regionIdx][oilCompIdx] = BrineCo2Pvt<Scalar>::Brine::molarMass();
286  molarMass_[regionIdx][gasCompIdx] = BrineCo2Pvt<Scalar>::CO2::molarMass();
287  }
288  }
289 
290  setEnableDiffusion(eclState.getSimulationConfig().isDiffusive());
291  if(enableDiffusion()) {
292  const auto& diffCoeffTables = eclState.getTableManager().getDiffusionCoefficientTable();
293  if(!diffCoeffTables.empty()) {
294  // if diffusion coefficient table is empty we relay on the PVT model to
295  // to give us the coefficients.
296  diffusionCoefficients_.resize(numRegions,{0,0,0,0,0,0,0,0,0});
297  assert(diffCoeffTables.size() == numRegions);
298  for (unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
299  const auto& diffCoeffTable = diffCoeffTables[regionIdx];
300  molarMass_[regionIdx][oilCompIdx] = diffCoeffTable.oil_mw;
301  molarMass_[regionIdx][gasCompIdx] = diffCoeffTable.gas_mw;
302  setDiffusionCoefficient(diffCoeffTable.gas_in_gas, gasCompIdx, gasPhaseIdx, regionIdx);
303  setDiffusionCoefficient(diffCoeffTable.oil_in_gas, oilCompIdx, gasPhaseIdx, regionIdx);
304  setDiffusionCoefficient(diffCoeffTable.gas_in_oil, gasCompIdx, oilPhaseIdx, regionIdx);
305  setDiffusionCoefficient(diffCoeffTable.oil_in_oil, oilCompIdx, oilPhaseIdx, regionIdx);
306  if(diffCoeffTable.gas_in_oil_cross_phase > 0 || diffCoeffTable.oil_in_oil_cross_phase > 0) {
307  throw std::runtime_error("Cross phase diffusion is set in the deck, but not implemented in Flow. "
308  "Please default DIFFC item 7 and item 8 or set it to zero.");
309  }
310  }
311  }
312  }
313  }
314 #endif // HAVE_ECL_INPUT
315 
324  static void initBegin(size_t numPvtRegions)
325  {
326  isInitialized_ = false;
327 
328  enableDissolvedGas_ = true;
329  enableVaporizedOil_ = false;
330  enableVaporizedWater_ = false;
331  enableDiffusion_ = false;
332 
333  oilPvt_ = nullptr;
334  gasPvt_ = nullptr;
335  waterPvt_ = nullptr;
336 
337  surfaceTemperature = 273.15 + 15.56; // [K]
338  surfacePressure = 1.01325e5; // [Pa]
340 
341  numActivePhases_ = numPhases;
342  std::fill_n(&phaseIsActive_[0], numPhases, true);
343 
344  resizeArrays_(numPvtRegions);
345  }
346 
353  static void setEnableDissolvedGas(bool yesno)
354  { enableDissolvedGas_ = yesno; }
355 
362  static void setEnableVaporizedOil(bool yesno)
363  { enableVaporizedOil_ = yesno; }
364 
371  static void setEnableVaporizedWater(bool yesno)
372  { enableVaporizedWater_ = yesno; }
373 
379  static void setEnableDiffusion(bool yesno)
380  { enableDiffusion_ = yesno; }
381 
382 
386  static void setGasPvt(std::shared_ptr<GasPvt> pvtObj)
387  { gasPvt_ = pvtObj; }
388 
392  static void setOilPvt(std::shared_ptr<OilPvt> pvtObj)
393  { oilPvt_ = pvtObj; }
394 
398  static void setWaterPvt(std::shared_ptr<WaterPvt> pvtObj)
399  { waterPvt_ = pvtObj; }
400 
408  static void setReferenceDensities(Scalar rhoOil,
409  Scalar rhoWater,
410  Scalar rhoGas,
411  unsigned regionIdx)
412  {
413  referenceDensity_[regionIdx][oilPhaseIdx] = rhoOil;
414  referenceDensity_[regionIdx][waterPhaseIdx] = rhoWater;
415  referenceDensity_[regionIdx][gasPhaseIdx] = rhoGas;
416  }
417 
418 
422  static void initEnd()
423  {
424  // calculate the final 2D functions which are used for interpolation.
425  size_t numRegions = molarMass_.size();
426  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
427  // calculate molar masses
428 
429  // water is simple: 18 g/mol
430  molarMass_[regionIdx][waterCompIdx] = 18e-3;
431 
432  if (phaseIsActive(gasPhaseIdx)) {
433  // for gas, we take the density at standard conditions and assume it to be ideal
436  Scalar rho_g = referenceDensity_[/*regionIdx=*/0][gasPhaseIdx];
437  molarMass_[regionIdx][gasCompIdx] = Constants<Scalar>::R*T*rho_g / p;
438  }
439  else
440  // hydrogen gas. we just set this do avoid NaNs later
441  molarMass_[regionIdx][gasCompIdx] = 2e-3;
442 
443  // finally, for oil phase, we take the molar mass from the spe9 paper
444  molarMass_[regionIdx][oilCompIdx] = 175e-3; // kg/mol
445  }
446 
447 
448  int activePhaseIdx = 0;
449  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
450  if(phaseIsActive(phaseIdx)){
451  canonicalToActivePhaseIdx_[phaseIdx] = activePhaseIdx;
452  activeToCanonicalPhaseIdx_[activePhaseIdx] = phaseIdx;
453  activePhaseIdx++;
454  }
455  }
456  isInitialized_ = true;
457  }
458 
459  static bool isInitialized()
460  { return isInitialized_; }
461 
462  /****************************************
463  * Generic phase properties
464  ****************************************/
465 
467  static const unsigned numPhases = 3;
468 
470  static const unsigned waterPhaseIdx = IndexTraits::waterPhaseIdx;
472  static const unsigned oilPhaseIdx = IndexTraits::oilPhaseIdx;
474  static const unsigned gasPhaseIdx = IndexTraits::gasPhaseIdx;
475 
478 
481 
483  static const char* phaseName(unsigned phaseIdx)
484  {
485  switch (phaseIdx) {
486  case waterPhaseIdx:
487  return "water";
488  case oilPhaseIdx:
489  return "oil";
490  case gasPhaseIdx:
491  return "gas";
492 
493  default:
494  throw std::logic_error("Phase index " + std::to_string(phaseIdx) + " is unknown");
495  }
496  }
497 
499  static bool isLiquid(unsigned phaseIdx)
500  {
501  assert(phaseIdx < numPhases);
502  return phaseIdx != gasPhaseIdx;
503  }
504 
505  /****************************************
506  * Generic component related properties
507  ****************************************/
508 
510  static const unsigned numComponents = 3;
511 
513  static const unsigned oilCompIdx = IndexTraits::oilCompIdx;
515  static const unsigned waterCompIdx = IndexTraits::waterCompIdx;
517  static const unsigned gasCompIdx = IndexTraits::gasCompIdx;
518 
519 protected:
520  static unsigned char numActivePhases_;
521  static std::array<bool,numPhases> phaseIsActive_;
522 
523 public:
525  static unsigned numActivePhases()
526  { return numActivePhases_; }
527 
529  static unsigned phaseIsActive(unsigned phaseIdx)
530  {
531  assert(phaseIdx < numPhases);
532  return phaseIsActive_[phaseIdx];
533  }
534 
536  static constexpr unsigned solventComponentIndex(unsigned phaseIdx)
537  {
538  switch (phaseIdx) {
539  case waterPhaseIdx:
540  return waterCompIdx;
541  case oilPhaseIdx:
542  return oilCompIdx;
543  case gasPhaseIdx:
544  return gasCompIdx;
545 
546  default:
547  throw std::logic_error("Phase index " + std::to_string(phaseIdx) + " is unknown");
548  }
549  }
550 
552  static constexpr unsigned soluteComponentIndex(unsigned phaseIdx)
553  {
554  switch (phaseIdx) {
555  case waterPhaseIdx:
556  throw std::logic_error("The water phase does not have any solutes in the black oil model!");
557  case oilPhaseIdx:
558  return gasCompIdx;
559  case gasPhaseIdx:
560  return oilCompIdx;
561 
562  default:
563  throw std::logic_error("Phase index " + std::to_string(phaseIdx) + " is unknown");
564  }
565  }
566 
568  static const char* componentName(unsigned compIdx)
569  {
570  switch (compIdx) {
571  case waterCompIdx:
572  return "Water";
573  case oilCompIdx:
574  return "Oil";
575  case gasCompIdx:
576  return "Gas";
577 
578  default:
579  throw std::logic_error("Component index " + std::to_string(compIdx) + " is unknown");
580  }
581  }
582 
584  static Scalar molarMass(unsigned compIdx, unsigned regionIdx = 0)
585  { return molarMass_[regionIdx][compIdx]; }
586 
588  static bool isIdealMixture(unsigned /*phaseIdx*/)
589  {
590  // fugacity coefficients are only pressure dependent -> we
591  // have an ideal mixture
592  return true;
593  }
594 
596  static bool isCompressible(unsigned /*phaseIdx*/)
597  { return true; /* all phases are compressible */ }
598 
600  static bool isIdealGas(unsigned /*phaseIdx*/)
601  { return false; }
602 
603 
604  /****************************************
605  * Black-oil specific properties
606  ****************************************/
612  static size_t numRegions()
613  { return molarMass_.size(); }
614 
621  static bool enableDissolvedGas()
622  { return enableDissolvedGas_; }
623 
630  static bool enableVaporizedOil()
631  { return enableVaporizedOil_; }
632 
639  static bool enableVaporizedWater()
640  { return enableVaporizedWater_; }
641 
647  static bool enableDiffusion()
648  { return enableDiffusion_; }
649 
655  static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
656  { return referenceDensity_[regionIdx][phaseIdx]; }
657 
658  /****************************************
659  * thermodynamic quantities (generic version)
660  ****************************************/
662  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
663  static LhsEval density(const FluidState& fluidState,
664  const ParameterCache<ParamCacheEval>& paramCache,
665  unsigned phaseIdx)
666  { return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
667 
669  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
670  static LhsEval fugacityCoefficient(const FluidState& fluidState,
671  const ParameterCache<ParamCacheEval>& paramCache,
672  unsigned phaseIdx,
673  unsigned compIdx)
674  {
675  return fugacityCoefficient<FluidState, LhsEval>(fluidState,
676  phaseIdx,
677  compIdx,
678  paramCache.regionIndex());
679  }
680 
682  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
683  static LhsEval viscosity(const FluidState& fluidState,
684  const ParameterCache<ParamCacheEval>& paramCache,
685  unsigned phaseIdx)
686  { return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
687 
689  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
690  static LhsEval enthalpy(const FluidState& fluidState,
691  const ParameterCache<ParamCacheEval>& paramCache,
692  unsigned phaseIdx)
693  { return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
694 
695  /****************************************
696  * thermodynamic quantities (black-oil specific version: Note that the PVT region
697  * index is explicitly passed instead of a parameter cache object)
698  ****************************************/
700  template <class FluidState, class LhsEval = typename FluidState::Scalar>
701  static LhsEval density(const FluidState& fluidState,
702  unsigned phaseIdx,
703  unsigned regionIdx)
704  {
705  assert(phaseIdx <= numPhases);
706  assert(regionIdx <= numRegions());
707 
708  const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
709  const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
710  const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
711 
712  switch (phaseIdx) {
713  case oilPhaseIdx: {
714  if (enableDissolvedGas()) {
715  // miscible oil
716  const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
717  const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
718 
719  return
720  bo*referenceDensity(oilPhaseIdx, regionIdx)
721  + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
722  }
723 
724  // immiscible oil
725  const LhsEval Rs(0.0);
726  const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
727 
728  return referenceDensity(phaseIdx, regionIdx)*bo;
729  }
730 
731  case gasPhaseIdx: {
732  if (enableVaporizedOil()) {
733  // miscible gas
734  const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
735  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
736 
737  return
738  bg*referenceDensity(gasPhaseIdx, regionIdx)
739  + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
740  }
741  //vaporized oil simultaneous with vaporized water is not yet supported
742  if (enableVaporizedWater()) {
743  // gas containing vaporized water
744  const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
745  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rvw);
746 
747  return
748  bg*referenceDensity(gasPhaseIdx, regionIdx)
749  + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
750  }
751 
752  // immiscible gas
753  const LhsEval Rv(0.0);
754  const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
755  return bg*referenceDensity(phaseIdx, regionIdx);
756  }
757 
758  case waterPhaseIdx:
759  return
760  referenceDensity(waterPhaseIdx, regionIdx)
761  * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
762  }
763 
764  throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
765  }
766 
774  template <class FluidState, class LhsEval = typename FluidState::Scalar>
775  static LhsEval saturatedDensity(const FluidState& fluidState,
776  unsigned phaseIdx,
777  unsigned regionIdx)
778  {
779  assert(phaseIdx <= numPhases);
780  assert(regionIdx <= numRegions());
781 
782  const auto& p = fluidState.pressure(phaseIdx);
783  const auto& T = fluidState.temperature(phaseIdx);
784 
785  switch (phaseIdx) {
786  case oilPhaseIdx: {
787  if (enableDissolvedGas()) {
788  // miscible oil
789  const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
790  const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
791 
792  return
793  bo*referenceDensity(oilPhaseIdx, regionIdx)
794  + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
795  }
796 
797  // immiscible oil
798  const LhsEval Rs(0.0);
799  const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
800  return referenceDensity(phaseIdx, regionIdx)*bo;
801  }
802 
803  case gasPhaseIdx: {
804  if (enableVaporizedOil()) {
805  // miscible gas
806  const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
807  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
808 
809  return
810  bg*referenceDensity(gasPhaseIdx, regionIdx)
811  + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
812  }
813  //vaporized oil simultaneous with vaporized water is not yet supported
814  if (enableVaporizedWater()) {
815  // gas containing vaporized water
816  const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
817  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rvw);
818 
819  return
820  bg*referenceDensity(gasPhaseIdx, regionIdx)
821  + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
822  }
823 
824  // immiscible gas
825  const LhsEval Rv(0.0);
826  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
827 
828  return referenceDensity(phaseIdx, regionIdx)*bg;
829 
830  }
831 
832  case waterPhaseIdx:
833  return
834  referenceDensity(waterPhaseIdx, regionIdx)
835  *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
836  }
837 
838  throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
839  }
840 
849  template <class FluidState, class LhsEval = typename FluidState::Scalar>
850  static LhsEval inverseFormationVolumeFactor(const FluidState& fluidState,
851  unsigned phaseIdx,
852  unsigned regionIdx)
853  {
854  assert(phaseIdx <= numPhases);
855  assert(regionIdx <= numRegions());
856 
857  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
858  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
859  const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
860 
861  switch (phaseIdx) {
862  case oilPhaseIdx: {
863  if (enableDissolvedGas()) {
864  const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
865  if (fluidState.saturation(gasPhaseIdx) > 0.0
866  && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
867  {
868  return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
869  } else {
870  return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
871  }
872  }
873 
874  const LhsEval Rs(0.0);
875  return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
876  }
877  case gasPhaseIdx: {
878  if (enableVaporizedOil()) {
879  const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
880  if (fluidState.saturation(oilPhaseIdx) > 0.0
881  && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
882  {
883  return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
884  } else {
885  return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
886  }
887  }
888 
889  const LhsEval Rv(0.0);
890  return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
891  }
892  case waterPhaseIdx:
893  return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
894  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
895  }
896  }
897 
905  template <class FluidState, class LhsEval = typename FluidState::Scalar>
906  static LhsEval saturatedInverseFormationVolumeFactor(const FluidState& fluidState,
907  unsigned phaseIdx,
908  unsigned regionIdx)
909  {
910  assert(phaseIdx <= numPhases);
911  assert(regionIdx <= numRegions());
912 
913  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
914  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
915  const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
916 
917  switch (phaseIdx) {
918  case oilPhaseIdx: return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
919  case gasPhaseIdx: return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
920  case waterPhaseIdx: return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
921  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
922  }
923  }
924 
926  template <class FluidState, class LhsEval = typename FluidState::Scalar>
927  static LhsEval fugacityCoefficient(const FluidState& fluidState,
928  unsigned phaseIdx,
929  unsigned compIdx,
930  unsigned regionIdx)
931  {
932  assert(phaseIdx <= numPhases);
933  assert(compIdx <= numComponents);
934  assert(regionIdx <= numRegions());
935 
936  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
937  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
938 
939  // for the fugacity coefficient of the oil component in the oil phase, we use
940  // some pseudo-realistic value for the vapor pressure to ease physical
941  // interpretation of the results
942  const LhsEval phi_oO = 20e3/p;
943 
944  // for the gas component in the gas phase, assume it to be an ideal gas
945  const Scalar phi_gG = 1.0;
946 
947  // for the fugacity coefficient of the water component in the water phase, we use
948  // the same approach as for the oil component in the oil phase
949  const LhsEval phi_wW = 30e3/p;
950 
951  switch (phaseIdx) {
952  case gasPhaseIdx: // fugacity coefficients for all components in the gas phase
953  switch (compIdx) {
954  case gasCompIdx:
955  return phi_gG;
956 
957  // for the oil component, we calculate the Rv value for saturated gas and Rs
958  // for saturated oil, and compute the fugacity coefficients at the
959  // equilibrium. for this, we assume that in equilibrium the fugacities of the
960  // oil component is the same in both phases.
961  case oilCompIdx: {
962  if (!enableVaporizedOil())
963  // if there's no vaporized oil, the gas phase is assumed to be
964  // immiscible with the oil component
965  return phi_gG*1e6;
966 
967  const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
968  const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
969  const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
970 
971  const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
972  const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
973  const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
974  const auto& x_oOSat = 1.0 - x_oGSat;
975 
976  const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
977  const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
978 
979  return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
980  }
981 
982  case waterCompIdx:
983  // the water component is assumed to be never miscible with the gas phase
984  return phi_gG*1e6;
985 
986  default:
987  throw std::logic_error("Invalid component index "+std::to_string(compIdx));
988  }
989 
990  case oilPhaseIdx: // fugacity coefficients for all components in the oil phase
991  switch (compIdx) {
992  case oilCompIdx:
993  return phi_oO;
994 
995  // for the oil and water components, we proceed analogous to the gas and
996  // water components in the gas phase
997  case gasCompIdx: {
998  if (!enableDissolvedGas())
999  // if there's no dissolved gas, the oil phase is assumed to be
1000  // immiscible with the gas component
1001  return phi_oO*1e6;
1002 
1003  const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1004  const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
1005  const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
1006  const auto& x_gGSat = 1.0 - x_gOSat;
1007 
1008  const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1009  const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
1010  const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
1011 
1012  const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
1013  const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
1014 
1015  return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
1016  }
1017 
1018  case waterCompIdx:
1019  return phi_oO*1e6;
1020 
1021  default:
1022  throw std::logic_error("Invalid component index "+std::to_string(compIdx));
1023  }
1024 
1025  case waterPhaseIdx: // fugacity coefficients for all components in the water phase
1026  // the water phase fugacity coefficients are pretty simple: because the water
1027  // phase is assumed to consist entirely from the water component, we just
1028  // need to make sure that the fugacity coefficients for the other components
1029  // are a few orders of magnitude larger than the one of the water
1030  // component. (i.e., the affinity of the gas and oil components to the water
1031  // phase is lower by a few orders of magnitude)
1032  switch (compIdx) {
1033  case waterCompIdx: return phi_wW;
1034  case oilCompIdx: return 1.1e6*phi_wW;
1035  case gasCompIdx: return 1e6*phi_wW;
1036  default:
1037  throw std::logic_error("Invalid component index "+std::to_string(compIdx));
1038  }
1039 
1040  default:
1041  throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
1042  }
1043 
1044  throw std::logic_error("Unhandled phase or component index");
1045  }
1046 
1048  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1049  static LhsEval viscosity(const FluidState& fluidState,
1050  unsigned phaseIdx,
1051  unsigned regionIdx)
1052  {
1053  assert(phaseIdx <= numPhases);
1054  assert(regionIdx <= numRegions());
1055 
1056  const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1057  const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1058  const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1059 
1060  switch (phaseIdx) {
1061  case oilPhaseIdx: {
1062  if (enableDissolvedGas()) {
1063  const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1064  if (fluidState.saturation(gasPhaseIdx) > 0.0
1065  && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
1066  {
1067  return oilPvt_->saturatedViscosity(regionIdx, T, p);
1068  } else {
1069  return oilPvt_->viscosity(regionIdx, T, p, Rs);
1070  }
1071  }
1072 
1073  const LhsEval Rs(0.0);
1074  return oilPvt_->viscosity(regionIdx, T, p, Rs);
1075  }
1076 
1077  case gasPhaseIdx: {
1078  if (enableVaporizedOil()) {
1079  const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1080  if (fluidState.saturation(oilPhaseIdx) > 0.0
1081  && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1082  {
1083  return gasPvt_->saturatedViscosity(regionIdx, T, p);
1084  } else {
1085  return gasPvt_->viscosity(regionIdx, T, p, Rv);
1086  }
1087  }
1088 
1089  const LhsEval Rv(0.0);
1090  return gasPvt_->viscosity(regionIdx, T, p, Rv);
1091  }
1092 
1093  case waterPhaseIdx:
1094  // since water is always assumed to be immiscible in the black-oil model,
1095  // there is no "saturated water"
1096  return waterPvt_->viscosity(regionIdx, T, p, saltConcentration);
1097  }
1098 
1099  throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1100  }
1101 
1103  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1104  static LhsEval enthalpy(const FluidState& fluidState,
1105  unsigned phaseIdx,
1106  unsigned regionIdx)
1107  {
1108  assert(phaseIdx <= numPhases);
1109  assert(regionIdx <= numRegions());
1110 
1111  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1112  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1113 
1114  switch (phaseIdx) {
1115  case oilPhaseIdx:
1116  return
1117  oilPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1118  + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1119 
1120  case gasPhaseIdx:
1121  return
1122  gasPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1123  + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1124 
1125  case waterPhaseIdx:
1126  return
1127  waterPvt_->internalEnergy(regionIdx, T, p)
1128  + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1129 
1130  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1131  }
1132 
1133  throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1134  }
1135 
1142  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1143  static LhsEval saturatedVaporizationFactor(const FluidState& fluidState,
1144  unsigned phaseIdx,
1145  unsigned regionIdx)
1146  {
1147  assert(phaseIdx <= numPhases);
1148  assert(regionIdx <= numRegions());
1149 
1150  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1151  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1152 
1153  switch (phaseIdx) {
1154  case oilPhaseIdx: return 0.0;
1155  case gasPhaseIdx: return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p);
1156  case waterPhaseIdx: return 0.0;
1157  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1158  }
1159  }
1160 
1167  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1168  static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1169  unsigned phaseIdx,
1170  unsigned regionIdx,
1171  const LhsEval& maxOilSaturation)
1172  {
1173  assert(phaseIdx <= numPhases);
1174  assert(regionIdx <= numRegions());
1175 
1176  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1177  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1178  const auto& So = decay<LhsEval>(fluidState.saturation(oilPhaseIdx));
1179 
1180  switch (phaseIdx) {
1181  case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1182  case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1183  case waterPhaseIdx: return 0.0;
1184  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1185  }
1186  }
1187 
1196  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1197  static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1198  unsigned phaseIdx,
1199  unsigned regionIdx)
1200  {
1201  assert(phaseIdx <= numPhases);
1202  assert(regionIdx <= numRegions());
1203 
1204  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1205  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1206 
1207  switch (phaseIdx) {
1208  case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1209  case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1210  case waterPhaseIdx: return 0.0;
1211  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1212  }
1213  }
1214 
1218  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1219  static LhsEval bubblePointPressure(const FluidState& fluidState,
1220  unsigned regionIdx)
1221  {
1222  return saturationPressure(fluidState, oilPhaseIdx, regionIdx);
1223  }
1224 
1225 
1229  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1230  static LhsEval dewPointPressure(const FluidState& fluidState,
1231  unsigned regionIdx)
1232  {
1233  return saturationPressure(fluidState, gasPhaseIdx, regionIdx);
1234  }
1235 
1246  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1247  static LhsEval saturationPressure(const FluidState& fluidState,
1248  unsigned phaseIdx,
1249  unsigned regionIdx)
1250  {
1251  assert(phaseIdx <= numPhases);
1252  assert(regionIdx <= numRegions());
1253 
1254  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1255 
1256  switch (phaseIdx) {
1257  case oilPhaseIdx: return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1258  case gasPhaseIdx: return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1259  case waterPhaseIdx: return 0.0;
1260  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1261  }
1262  }
1263 
1264  /****************************************
1265  * Auxiliary and convenience methods for the black-oil model
1266  ****************************************/
1271  template <class LhsEval>
1272  static LhsEval convertXoGToRs(const LhsEval& XoG, unsigned regionIdx)
1273  {
1274  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1275  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1276 
1277  return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1278  }
1279 
1284  template <class LhsEval>
1285  static LhsEval convertXgOToRv(const LhsEval& XgO, unsigned regionIdx)
1286  {
1287  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1288  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1289 
1290  return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1291  }
1292 
1297  template <class LhsEval>
1298  static LhsEval convertXgWToRvw(const LhsEval& XgW, unsigned regionIdx)
1299  {
1300  Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1301  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1302 
1303  return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1304  }
1305 
1306 
1311  template <class LhsEval>
1312  static LhsEval convertRsToXoG(const LhsEval& Rs, unsigned regionIdx)
1313  {
1314  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1315  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1316 
1317  const LhsEval& rho_oG = Rs*rho_gRef;
1318  return rho_oG/(rho_oRef + rho_oG);
1319  }
1320 
1325  template <class LhsEval>
1326  static LhsEval convertRvToXgO(const LhsEval& Rv, unsigned regionIdx)
1327  {
1328  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1329  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1330 
1331  const LhsEval& rho_gO = Rv*rho_oRef;
1332  return rho_gO/(rho_gRef + rho_gO);
1333  }
1334 
1339  template <class LhsEval>
1340  static LhsEval convertRvwToXgW(const LhsEval& Rvw, unsigned regionIdx)
1341  {
1342  Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1343  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1344 
1345  const LhsEval& rho_gW = Rvw*rho_wRef;
1346  return rho_gW/(rho_gRef + rho_gW);
1347  }
1348 
1352  template <class LhsEval>
1353  static LhsEval convertXoGToxoG(const LhsEval& XoG, unsigned regionIdx)
1354  {
1355  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1356  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1357 
1358  return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1359  }
1360 
1364  template <class LhsEval>
1365  static LhsEval convertxoGToXoG(const LhsEval& xoG, unsigned regionIdx)
1366  {
1367  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1368  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1369 
1370  return xoG*MG / (xoG*(MG - MO) + MO);
1371  }
1372 
1376  template <class LhsEval>
1377  static LhsEval convertXgOToxgO(const LhsEval& XgO, unsigned regionIdx)
1378  {
1379  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1380  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1381 
1382  return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1383  }
1384 
1388  template <class LhsEval>
1389  static LhsEval convertxgOToXgO(const LhsEval& xgO, unsigned regionIdx)
1390  {
1391  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1392  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1393 
1394  return xgO*MO / (xgO*(MO - MG) + MG);
1395  }
1396 
1404  static const GasPvt& gasPvt()
1405  { return *gasPvt_; }
1406 
1414  static const OilPvt& oilPvt()
1415  { return *oilPvt_; }
1416 
1424  static const WaterPvt& waterPvt()
1425  { return *waterPvt_; }
1426 
1432  static Scalar reservoirTemperature(unsigned = 0)
1433  { return reservoirTemperature_; }
1434 
1440  static void setReservoirTemperature(Scalar value)
1441  { reservoirTemperature_ = value; }
1442 
1443  static short activeToCanonicalPhaseIdx(unsigned activePhaseIdx) {
1444  assert(activePhaseIdx<numActivePhases());
1445  return activeToCanonicalPhaseIdx_[activePhaseIdx];
1446  }
1447 
1448  static short canonicalToActivePhaseIdx(unsigned phaseIdx) {
1449  assert(phaseIdx<numPhases);
1450  assert(phaseIsActive(phaseIdx));
1451  return canonicalToActivePhaseIdx_[phaseIdx];
1452  }
1453 
1455  static Scalar diffusionCoefficient(unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1456  { return diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx]; }
1457 
1459  static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1460  { diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx] = coefficient ; }
1461 
1465  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
1466  static LhsEval diffusionCoefficient(const FluidState& fluidState,
1467  const ParameterCache<ParamCacheEval>& paramCache,
1468  unsigned phaseIdx,
1469  unsigned compIdx)
1470  {
1471  // diffusion is disabled by the user
1472  if(!enableDiffusion())
1473  return 0.0;
1474 
1475  // diffusion coefficients are set, and we use them
1476  if(!diffusionCoefficients_.empty()) {
1477  return diffusionCoefficient(compIdx, phaseIdx, paramCache.regionIndex());
1478  }
1479 
1480  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1481  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1482 
1483  switch (phaseIdx) {
1484  case oilPhaseIdx: return oilPvt().diffusionCoefficient(T, p, compIdx);
1485  case gasPhaseIdx: return gasPvt().diffusionCoefficient(T, p, compIdx);
1486  case waterPhaseIdx: return 0.0;
1487  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1488  }
1489  }
1490 
1491 private:
1492  static void resizeArrays_(size_t numRegions)
1493  {
1494  molarMass_.resize(numRegions);
1495  referenceDensity_.resize(numRegions);
1496  }
1497 
1498  static Scalar reservoirTemperature_;
1499 
1500  static std::shared_ptr<GasPvt> gasPvt_;
1501  static std::shared_ptr<OilPvt> oilPvt_;
1502  static std::shared_ptr<WaterPvt> waterPvt_;
1503 
1504  static bool enableDissolvedGas_;
1505  static bool enableVaporizedOil_;
1506  static bool enableVaporizedWater_;
1507  static bool enableDiffusion_;
1508 
1509  // HACK for GCC 4.4: the array size has to be specified using the literal value '3'
1510  // here, because GCC 4.4 seems to be unable to determine the number of phases from
1511  // the BlackOil fluid system in the attribute declaration below...
1512  static std::vector<std::array<Scalar, /*numPhases=*/3> > referenceDensity_;
1513  static std::vector<std::array<Scalar, /*numComponents=*/3> > molarMass_;
1514  static std::vector<std::array<Scalar, /*numComponents=*/3 * /*numPhases=*/3> > diffusionCoefficients_;
1515 
1516  static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1517  static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1518 
1519  static bool isInitialized_;
1520 };
1521 
1522 template <class Scalar, class IndexTraits>
1523 unsigned char BlackOilFluidSystem<Scalar, IndexTraits>::numActivePhases_;
1524 
1525 template <class Scalar, class IndexTraits>
1526 std::array<bool, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::phaseIsActive_;
1527 
1528 template <class Scalar, class IndexTraits>
1529 std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::activeToCanonicalPhaseIdx_;
1530 
1531 template <class Scalar, class IndexTraits>
1532 std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::canonicalToActivePhaseIdx_;
1533 
1534 template <class Scalar, class IndexTraits>
1535 Scalar
1537 
1538 template <class Scalar, class IndexTraits>
1539 Scalar
1541 
1542 template <class Scalar, class IndexTraits>
1543 Scalar
1544 BlackOilFluidSystem<Scalar, IndexTraits>::reservoirTemperature_;
1545 
1546 template <class Scalar, class IndexTraits>
1547 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGas_;
1548 
1549 template <class Scalar, class IndexTraits>
1550 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedOil_;
1551 
1552 template <class Scalar, class IndexTraits>
1553 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedWater_;
1554 
1555 template <class Scalar, class IndexTraits>
1556 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDiffusion_;
1557 
1558 template <class Scalar, class IndexTraits>
1559 std::shared_ptr<OilPvtMultiplexer<Scalar> >
1560 BlackOilFluidSystem<Scalar, IndexTraits>::oilPvt_;
1561 
1562 template <class Scalar, class IndexTraits>
1563 std::shared_ptr<GasPvtMultiplexer<Scalar> >
1564 BlackOilFluidSystem<Scalar, IndexTraits>::gasPvt_;
1565 
1566 template <class Scalar, class IndexTraits>
1567 std::shared_ptr<WaterPvtMultiplexer<Scalar> >
1568 BlackOilFluidSystem<Scalar, IndexTraits>::waterPvt_;
1569 
1570 template <class Scalar, class IndexTraits>
1571 std::vector<std::array<Scalar, 3> >
1572 BlackOilFluidSystem<Scalar, IndexTraits>::referenceDensity_;
1573 
1574 template <class Scalar, class IndexTraits>
1575 std::vector<std::array<Scalar, 3> >
1576 BlackOilFluidSystem<Scalar, IndexTraits>::molarMass_;
1577 
1578 template <class Scalar, class IndexTraits>
1579 std::vector<std::array<Scalar, 9> >
1580 BlackOilFluidSystem<Scalar, IndexTraits>::diffusionCoefficients_;
1581 
1582 template <class Scalar, class IndexTraits>
1583 bool BlackOilFluidSystem<Scalar, IndexTraits>::isInitialized_ = false;
1584 
1585 } // namespace Opm
1586 
1587 #endif
The base class for all fluid systems.
The class which specifies the default phase and component indices for the black-oil fluid system.
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a CO2-Brine s...
A central place for various physical constants occuring in some equations.
Provides the opm-material specific exception classes.
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
This macro generates a class HasMember_${MEMBER_NAME} which can be used for template specialization.
#define OPM_GENERATE_HAS_MEMBER(MEMBER_NAME,...)
This macro generates a class HasMember_${MEMBER_NAME} which can be used for template specialization.
Definition: HasMemberGeneratorMacros.hpp:49
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
Some templates to wrap the valgrind client request macros.
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
The base class for all fluid systems.
Definition: BaseFluidSystem.hpp:44
Scalar Scalar
The type used for scalar quantities.
Definition: BaseFluidSystem.hpp:49
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
Definition: BlackOilFluidSystem.hpp:141
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition: BlackOilFluidSystem.hpp:525
static Scalar diffusionCoefficient(unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition: BlackOilFluidSystem.hpp:1455
static void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition: BlackOilFluidSystem.hpp:408
static LhsEval convertXoGToRs(const LhsEval &XoG, unsigned regionIdx)
Convert the mass fraction of the gas component in the oil phase to the corresponding gas dissolution ...
Definition: BlackOilFluidSystem.hpp:1272
static LhsEval convertRsToXoG(const LhsEval &Rs, unsigned regionIdx)
Convert a gas dissolution factor to the the corresponding mass fraction of the gas component in the o...
Definition: BlackOilFluidSystem.hpp:1312
static void initEnd()
Finish initializing the black oil fluid system.
Definition: BlackOilFluidSystem.hpp:422
static LhsEval convertRvwToXgW(const LhsEval &Rvw, unsigned regionIdx)
Convert an water vaporization factor to the corresponding mass fraction of the water component in the...
Definition: BlackOilFluidSystem.hpp:1340
static constexpr unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition: BlackOilFluidSystem.hpp:552
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:906
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1353
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: BlackOilFluidSystem.hpp:1049
static const unsigned oilPhaseIdx
Index of the oil phase.
Definition: BlackOilFluidSystem.hpp:472
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1432
static void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:379
static const unsigned waterPhaseIdx
Index of the water phase.
Definition: BlackOilFluidSystem.hpp:470
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition: BlackOilFluidSystem.hpp:1219
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: BlackOilFluidSystem.hpp:701
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:775
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: BlackOilFluidSystem.hpp:683
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: BlackOilFluidSystem.hpp:588
static bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition: BlackOilFluidSystem.hpp:639
static const unsigned numComponents
Number of chemical species in the fluid system.
Definition: BlackOilFluidSystem.hpp:510
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition: BlackOilFluidSystem.hpp:568
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition: BlackOilFluidSystem.hpp:1404
static constexpr unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition: BlackOilFluidSystem.hpp:536
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition: BlackOilFluidSystem.hpp:398
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1440
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx, const LhsEval &maxOilSaturation)
Returns the dissolution factor of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:1168
static size_t numRegions()
Returns the number of PVT regions which are considered.
Definition: BlackOilFluidSystem.hpp:612
static LhsEval enthalpy(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy [J/kg...
Definition: BlackOilFluidSystem.hpp:1104
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition: BlackOilFluidSystem.hpp:499
static LhsEval enthalpy(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy [J/kg...
Definition: BlackOilFluidSystem.hpp:690
static LhsEval fugacityCoefficient(const FluidState &fluidState, unsigned phaseIdx, unsigned compIdx, unsigned regionIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition: BlackOilFluidSystem.hpp:927
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1365
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition: BlackOilFluidSystem.hpp:584
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:1197
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition: BlackOilFluidSystem.hpp:1230
static LhsEval saturationPressure(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the saturation pressure of a given phase [Pa] depending on its composition.
Definition: BlackOilFluidSystem.hpp:1247
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1377
static LhsEval diffusionCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx, unsigned compIdx)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition: BlackOilFluidSystem.hpp:1466
static LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition: BlackOilFluidSystem.hpp:1143
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition: BlackOilFluidSystem.hpp:621
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition: BlackOilFluidSystem.hpp:655
static Scalar surfaceTemperature
The temperature at the surface.
Definition: BlackOilFluidSystem.hpp:480
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition: BlackOilFluidSystem.hpp:392
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1389
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:850
static void setEnableDissolvedGas(bool yesno)
Specify whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition: BlackOilFluidSystem.hpp:353
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition: BlackOilFluidSystem.hpp:483
static const unsigned numPhases
Number of fluid phases in the fluid system.
Definition: BlackOilFluidSystem.hpp:467
static LhsEval convertXgOToRv(const LhsEval &XgO, unsigned regionIdx)
Convert the mass fraction of the oil component in the gas phase to the corresponding oil vaporization...
Definition: BlackOilFluidSystem.hpp:1285
static LhsEval fugacityCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition: BlackOilFluidSystem.hpp:670
static const unsigned gasPhaseIdx
Index of the gas phase.
Definition: BlackOilFluidSystem.hpp:474
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: BlackOilFluidSystem.hpp:663
static LhsEval convertRvToXgO(const LhsEval &Rv, unsigned regionIdx)
Convert an oil vaporization factor to the corresponding mass fraction of the oil component in the gas...
Definition: BlackOilFluidSystem.hpp:1326
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition: BlackOilFluidSystem.hpp:630
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: BlackOilFluidSystem.hpp:600
static unsigned phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition: BlackOilFluidSystem.hpp:529
static Scalar surfacePressure
The pressure at the surface.
Definition: BlackOilFluidSystem.hpp:477
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition: BlackOilFluidSystem.hpp:1424
static void setEnableVaporizedWater(bool yesno)
Specify whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition: BlackOilFluidSystem.hpp:371
static void initBegin(size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition: BlackOilFluidSystem.hpp:324
static bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:647
static LhsEval convertXgWToRvw(const LhsEval &XgW, unsigned regionIdx)
Convert the mass fraction of the water component in the gas phase to the corresponding water vaporiza...
Definition: BlackOilFluidSystem.hpp:1298
static const unsigned gasCompIdx
Index of the gas component.
Definition: BlackOilFluidSystem.hpp:517
static const unsigned oilCompIdx
Index of the oil component.
Definition: BlackOilFluidSystem.hpp:513
static void setEnableVaporizedOil(bool yesno)
Specify whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition: BlackOilFluidSystem.hpp:362
static const unsigned waterCompIdx
Index of the water component.
Definition: BlackOilFluidSystem.hpp:515
static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition: BlackOilFluidSystem.hpp:1459
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition: BlackOilFluidSystem.hpp:1414
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: BlackOilFluidSystem.hpp:596
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition: BlackOilFluidSystem.hpp:386
static Scalar molarMass()
The molar mass in of the component.
Definition: Brine.hpp:80
static Scalar molarMass()
The mass in [kg] of one mole of CO2.
Definition: CO2.hpp:66
A central place for various physical constants occuring in some equations.
Definition: Constants.hpp:41
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition: GasPvtMultiplexer.hpp:93
Evaluation diffusionCoefficient(const Evaluation &temperature, const Evaluation &pressure, unsigned compIdx) const
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition: GasPvtMultiplexer.hpp:302
A parameter cache which does nothing.
Definition: NullParameterCache.hpp:40
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
Definition: OilPvtMultiplexer.hpp:96
Evaluation diffusionCoefficient(const Evaluation &temperature, const Evaluation &pressure, unsigned compIdx) const
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition: OilPvtMultiplexer.hpp:271
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition: WaterPvtMultiplexer.hpp:75
The type of the fluid system's parameter cache.
Definition: BlackOilFluidSystem.hpp:152
void assignPersistentData(const OtherCache &other)
Copy the data which is not dependent on the type of the Scalars from another parameter cache.
Definition: BlackOilFluidSystem.hpp:170
void setRegionIndex(unsigned val)
Set the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:193
unsigned regionIndex() const
Return the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:183