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  using Evaluation = EvaluationT;
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 constexpr unsigned numPhases = 3;
468 
470  static constexpr unsigned waterPhaseIdx = IndexTraits::waterPhaseIdx;
472  static constexpr unsigned oilPhaseIdx = IndexTraits::oilPhaseIdx;
474  static constexpr 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 constexpr unsigned numComponents = 3;
511 
513  static constexpr unsigned oilCompIdx = IndexTraits::oilCompIdx;
515  static constexpr unsigned waterCompIdx = IndexTraits::waterCompIdx;
517  static constexpr 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: {
733  // gas containing vaporized oil and vaporized water
734  const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
735  const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
736  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
737 
738  return
739  bg*referenceDensity(gasPhaseIdx, regionIdx)
740  + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
741  + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
742  }
743  if (enableVaporizedOil()) {
744  // miscible gas
745  const LhsEval Rvw(0.0);
746  const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
747  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
748 
749  return
750  bg*referenceDensity(gasPhaseIdx, regionIdx)
751  + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
752  }
753  if (enableVaporizedWater()) {
754  // gas containing vaporized water
755  const LhsEval Rv(0.0);
756  const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
757  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
758 
759  return
760  bg*referenceDensity(gasPhaseIdx, regionIdx)
761  + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
762  }
763 
764  // immiscible gas
765  const LhsEval Rv(0.0);
766  const LhsEval Rvw(0.0);
767  const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
768  return bg*referenceDensity(phaseIdx, regionIdx);
769  }
770 
771  case waterPhaseIdx:
772  return
773  referenceDensity(waterPhaseIdx, regionIdx)
774  * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
775  }
776 
777  throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
778  }
779 
787  template <class FluidState, class LhsEval = typename FluidState::Scalar>
788  static LhsEval saturatedDensity(const FluidState& fluidState,
789  unsigned phaseIdx,
790  unsigned regionIdx)
791  {
792  assert(phaseIdx <= numPhases);
793  assert(regionIdx <= numRegions());
794 
795  const auto& p = fluidState.pressure(phaseIdx);
796  const auto& T = fluidState.temperature(phaseIdx);
797 
798  switch (phaseIdx) {
799  case oilPhaseIdx: {
800  if (enableDissolvedGas()) {
801  // miscible oil
802  const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
803  const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
804 
805  return
806  bo*referenceDensity(oilPhaseIdx, regionIdx)
807  + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
808  }
809 
810  // immiscible oil
811  const LhsEval Rs(0.0);
812  const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
813  return referenceDensity(phaseIdx, regionIdx)*bo;
814  }
815 
816  case gasPhaseIdx: {
818  // gas containing vaporized oil and vaporized water
819  const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
820  const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
821  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
822 
823  return
824  bg*referenceDensity(gasPhaseIdx, regionIdx)
825  + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
826  + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx) ;
827  }
828 
829  if (enableVaporizedOil()) {
830  // miscible gas
831  const LhsEval Rvw(0.0);
832  const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
833  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
834 
835  return
836  bg*referenceDensity(gasPhaseIdx, regionIdx)
837  + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
838  }
839 
840  if (enableVaporizedWater()) {
841  // gas containing vaporized water
842  const LhsEval Rv(0.0);
843  const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
844  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
845 
846  return
847  bg*referenceDensity(gasPhaseIdx, regionIdx)
848  + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
849  }
850 
851  // immiscible gas
852  const LhsEval Rv(0.0);
853  const LhsEval Rvw(0.0);
854  const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
855 
856  return referenceDensity(phaseIdx, regionIdx)*bg;
857 
858  }
859 
860  case waterPhaseIdx:
861  return
862  referenceDensity(waterPhaseIdx, regionIdx)
863  *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
864  }
865 
866  throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
867  }
868 
877  template <class FluidState, class LhsEval = typename FluidState::Scalar>
878  static LhsEval inverseFormationVolumeFactor(const FluidState& fluidState,
879  unsigned phaseIdx,
880  unsigned regionIdx)
881  {
882  assert(phaseIdx <= numPhases);
883  assert(regionIdx <= numRegions());
884 
885  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
886  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
887  const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
888 
889  switch (phaseIdx) {
890  case oilPhaseIdx: {
891  if (enableDissolvedGas()) {
892  const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
893  if (fluidState.saturation(gasPhaseIdx) > 0.0
894  && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
895  {
896  return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
897  } else {
898  return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
899  }
900  }
901 
902  const LhsEval Rs(0.0);
903  return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
904  }
905  case gasPhaseIdx: {
907  const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
908  const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
909  if (fluidState.saturation(waterPhaseIdx) > 0.0
910  && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
911  && fluidState.saturation(oilPhaseIdx) > 0.0
912  && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
913  {
914  return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
915  } else {
916  return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
917  }
918  }
919 
920  if (enableVaporizedOil()) {
921  const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
922  if (fluidState.saturation(oilPhaseIdx) > 0.0
923  && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
924  {
925  return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
926  } else {
927  const LhsEval Rvw(0.0);
928  return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
929  }
930  }
931 
932  if (enableVaporizedWater()) {
933  const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
934  if (fluidState.saturation(waterPhaseIdx) > 0.0
935  && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
936  {
937  return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
938  } else {
939  const LhsEval Rv(0.0);
940  return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
941  }
942  }
943 
944  const LhsEval Rv(0.0);
945  const LhsEval Rvw(0.0);
946  return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
947  }
948  case waterPhaseIdx:
949  return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
950  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
951  }
952  }
953 
961  template <class FluidState, class LhsEval = typename FluidState::Scalar>
962  static LhsEval saturatedInverseFormationVolumeFactor(const FluidState& fluidState,
963  unsigned phaseIdx,
964  unsigned regionIdx)
965  {
966  assert(phaseIdx <= numPhases);
967  assert(regionIdx <= numRegions());
968 
969  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
970  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
971  const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
972 
973  switch (phaseIdx) {
974  case oilPhaseIdx: return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
975  case gasPhaseIdx: return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
976  case waterPhaseIdx: return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
977  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
978  }
979  }
980 
982  template <class FluidState, class LhsEval = typename FluidState::Scalar>
983  static LhsEval fugacityCoefficient(const FluidState& fluidState,
984  unsigned phaseIdx,
985  unsigned compIdx,
986  unsigned regionIdx)
987  {
988  assert(phaseIdx <= numPhases);
989  assert(compIdx <= numComponents);
990  assert(regionIdx <= numRegions());
991 
992  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
993  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
994 
995  // for the fugacity coefficient of the oil component in the oil phase, we use
996  // some pseudo-realistic value for the vapor pressure to ease physical
997  // interpretation of the results
998  const LhsEval phi_oO = 20e3/p;
999 
1000  // for the gas component in the gas phase, assume it to be an ideal gas
1001  constexpr const Scalar phi_gG = 1.0;
1002 
1003  // for the fugacity coefficient of the water component in the water phase, we use
1004  // the same approach as for the oil component in the oil phase
1005  const LhsEval phi_wW = 30e3/p;
1006 
1007  switch (phaseIdx) {
1008  case gasPhaseIdx: // fugacity coefficients for all components in the gas phase
1009  switch (compIdx) {
1010  case gasCompIdx:
1011  return phi_gG;
1012 
1013  // for the oil component, we calculate the Rv value for saturated gas and Rs
1014  // for saturated oil, and compute the fugacity coefficients at the
1015  // equilibrium. for this, we assume that in equilibrium the fugacities of the
1016  // oil component is the same in both phases.
1017  case oilCompIdx: {
1018  if (!enableVaporizedOil())
1019  // if there's no vaporized oil, the gas phase is assumed to be
1020  // immiscible with the oil component
1021  return phi_gG*1e6;
1022 
1023  const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1024  const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
1025  const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
1026 
1027  const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1028  const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
1029  const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
1030  const auto& x_oOSat = 1.0 - x_oGSat;
1031 
1032  const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
1033  const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
1034 
1035  return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
1036  }
1037 
1038  case waterCompIdx:
1039  // the water component is assumed to be never miscible with the gas phase
1040  return phi_gG*1e6;
1041 
1042  default:
1043  throw std::logic_error("Invalid component index "+std::to_string(compIdx));
1044  }
1045 
1046  case oilPhaseIdx: // fugacity coefficients for all components in the oil phase
1047  switch (compIdx) {
1048  case oilCompIdx:
1049  return phi_oO;
1050 
1051  // for the oil and water components, we proceed analogous to the gas and
1052  // water components in the gas phase
1053  case gasCompIdx: {
1054  if (!enableDissolvedGas())
1055  // if there's no dissolved gas, the oil phase is assumed to be
1056  // immiscible with the gas component
1057  return phi_oO*1e6;
1058 
1059  const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1060  const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
1061  const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
1062  const auto& x_gGSat = 1.0 - x_gOSat;
1063 
1064  const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1065  const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
1066  const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
1067 
1068  const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
1069  const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
1070 
1071  return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
1072  }
1073 
1074  case waterCompIdx:
1075  return phi_oO*1e6;
1076 
1077  default:
1078  throw std::logic_error("Invalid component index "+std::to_string(compIdx));
1079  }
1080 
1081  case waterPhaseIdx: // fugacity coefficients for all components in the water phase
1082  // the water phase fugacity coefficients are pretty simple: because the water
1083  // phase is assumed to consist entirely from the water component, we just
1084  // need to make sure that the fugacity coefficients for the other components
1085  // are a few orders of magnitude larger than the one of the water
1086  // component. (i.e., the affinity of the gas and oil components to the water
1087  // phase is lower by a few orders of magnitude)
1088  switch (compIdx) {
1089  case waterCompIdx: return phi_wW;
1090  case oilCompIdx: return 1.1e6*phi_wW;
1091  case gasCompIdx: return 1e6*phi_wW;
1092  default:
1093  throw std::logic_error("Invalid component index "+std::to_string(compIdx));
1094  }
1095 
1096  default:
1097  throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
1098  }
1099 
1100  throw std::logic_error("Unhandled phase or component index");
1101  }
1102 
1104  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1105  static LhsEval viscosity(const FluidState& fluidState,
1106  unsigned phaseIdx,
1107  unsigned regionIdx)
1108  {
1109  assert(phaseIdx <= numPhases);
1110  assert(regionIdx <= numRegions());
1111 
1112  const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1113  const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1114  const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1115 
1116  switch (phaseIdx) {
1117  case oilPhaseIdx: {
1118  if (enableDissolvedGas()) {
1119  const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1120  if (fluidState.saturation(gasPhaseIdx) > 0.0
1121  && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
1122  {
1123  return oilPvt_->saturatedViscosity(regionIdx, T, p);
1124  } else {
1125  return oilPvt_->viscosity(regionIdx, T, p, Rs);
1126  }
1127  }
1128 
1129  const LhsEval Rs(0.0);
1130  return oilPvt_->viscosity(regionIdx, T, p, Rs);
1131  }
1132 
1133  case gasPhaseIdx: {
1135  const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1136  const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1137  if (fluidState.saturation(waterPhaseIdx) > 0.0
1138  && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1139  && fluidState.saturation(oilPhaseIdx) > 0.0
1140  && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1141  {
1142  return gasPvt_->saturatedViscosity(regionIdx, T, p);
1143  } else {
1144  return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1145  }
1146  }
1147  if (enableVaporizedOil()) {
1148  const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1149  if (fluidState.saturation(oilPhaseIdx) > 0.0
1150  && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1151  {
1152  return gasPvt_->saturatedViscosity(regionIdx, T, p);
1153  } else {
1154  const LhsEval Rvw(0.0);
1155  return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1156  }
1157  }
1158  if (enableVaporizedWater()) {
1159  const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1160  if (fluidState.saturation(waterPhaseIdx) > 0.0
1161  && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1162  {
1163  return gasPvt_->saturatedViscosity(regionIdx, T, p);
1164  } else {
1165  const LhsEval Rv(0.0);
1166  return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1167  }
1168  }
1169 
1170  const LhsEval Rv(0.0);
1171  const LhsEval Rvw(0.0);
1172  return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1173  }
1174 
1175  case waterPhaseIdx:
1176  // since water is always assumed to be immiscible in the black-oil model,
1177  // there is no "saturated water"
1178  return waterPvt_->viscosity(regionIdx, T, p, saltConcentration);
1179  }
1180 
1181  throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1182  }
1183 
1185  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1186  static LhsEval enthalpy(const FluidState& fluidState,
1187  unsigned phaseIdx,
1188  unsigned regionIdx)
1189  {
1190  assert(phaseIdx <= numPhases);
1191  assert(regionIdx <= numRegions());
1192 
1193  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1194  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1195 
1196  switch (phaseIdx) {
1197  case oilPhaseIdx:
1198  return
1199  oilPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1200  + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1201 
1202  case gasPhaseIdx:
1203  return
1204  gasPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1205  + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1206 
1207  case waterPhaseIdx:
1208  return
1209  waterPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1210  + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1211 
1212  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1213  }
1214 
1215  throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1216  }
1217 
1224  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1225  static LhsEval saturatedVaporizationFactor(const FluidState& fluidState,
1226  unsigned phaseIdx,
1227  unsigned regionIdx)
1228  {
1229  assert(phaseIdx <= numPhases);
1230  assert(regionIdx <= numRegions());
1231 
1232  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1233  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1234 
1235  switch (phaseIdx) {
1236  case oilPhaseIdx: return 0.0;
1237  case gasPhaseIdx: return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p);
1238  case waterPhaseIdx: return 0.0;
1239  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1240  }
1241  }
1242 
1249  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1250  static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1251  unsigned phaseIdx,
1252  unsigned regionIdx,
1253  const LhsEval& maxOilSaturation)
1254  {
1255  assert(phaseIdx <= numPhases);
1256  assert(regionIdx <= numRegions());
1257 
1258  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1259  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1260  const auto& So = decay<LhsEval>(fluidState.saturation(oilPhaseIdx));
1261 
1262  switch (phaseIdx) {
1263  case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1264  case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1265  case waterPhaseIdx: return 0.0;
1266  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1267  }
1268  }
1269 
1278  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1279  static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1280  unsigned phaseIdx,
1281  unsigned regionIdx)
1282  {
1283  assert(phaseIdx <= numPhases);
1284  assert(regionIdx <= numRegions());
1285 
1286  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1287  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1288 
1289  switch (phaseIdx) {
1290  case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1291  case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1292  case waterPhaseIdx: return 0.0;
1293  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1294  }
1295  }
1296 
1300  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1301  static LhsEval bubblePointPressure(const FluidState& fluidState,
1302  unsigned regionIdx)
1303  {
1304  return saturationPressure(fluidState, oilPhaseIdx, regionIdx);
1305  }
1306 
1307 
1311  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1312  static LhsEval dewPointPressure(const FluidState& fluidState,
1313  unsigned regionIdx)
1314  {
1315  return saturationPressure(fluidState, gasPhaseIdx, regionIdx);
1316  }
1317 
1328  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1329  static LhsEval saturationPressure(const FluidState& fluidState,
1330  unsigned phaseIdx,
1331  unsigned regionIdx)
1332  {
1333  assert(phaseIdx <= numPhases);
1334  assert(regionIdx <= numRegions());
1335 
1336  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1337 
1338  switch (phaseIdx) {
1339  case oilPhaseIdx: return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1340  case gasPhaseIdx: return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1341  case waterPhaseIdx: return 0.0;
1342  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1343  }
1344  }
1345 
1346  /****************************************
1347  * Auxiliary and convenience methods for the black-oil model
1348  ****************************************/
1353  template <class LhsEval>
1354  static LhsEval convertXoGToRs(const LhsEval& XoG, unsigned regionIdx)
1355  {
1356  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1357  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1358 
1359  return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1360  }
1361 
1366  template <class LhsEval>
1367  static LhsEval convertXgOToRv(const LhsEval& XgO, unsigned regionIdx)
1368  {
1369  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1370  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1371 
1372  return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1373  }
1374 
1379  template <class LhsEval>
1380  static LhsEval convertXgWToRvw(const LhsEval& XgW, unsigned regionIdx)
1381  {
1382  Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1383  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1384 
1385  return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1386  }
1387 
1388 
1393  template <class LhsEval>
1394  static LhsEval convertRsToXoG(const LhsEval& Rs, unsigned regionIdx)
1395  {
1396  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1397  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1398 
1399  const LhsEval& rho_oG = Rs*rho_gRef;
1400  return rho_oG/(rho_oRef + rho_oG);
1401  }
1402 
1407  template <class LhsEval>
1408  static LhsEval convertRvToXgO(const LhsEval& Rv, unsigned regionIdx)
1409  {
1410  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1411  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1412 
1413  const LhsEval& rho_gO = Rv*rho_oRef;
1414  return rho_gO/(rho_gRef + rho_gO);
1415  }
1416 
1421  template <class LhsEval>
1422  static LhsEval convertRvwToXgW(const LhsEval& Rvw, unsigned regionIdx)
1423  {
1424  Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1425  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1426 
1427  const LhsEval& rho_gW = Rvw*rho_wRef;
1428  return rho_gW/(rho_gRef + rho_gW);
1429  }
1430 
1434  template <class LhsEval>
1435  static LhsEval convertXoGToxoG(const LhsEval& XoG, unsigned regionIdx)
1436  {
1437  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1438  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1439 
1440  return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1441  }
1442 
1446  template <class LhsEval>
1447  static LhsEval convertxoGToXoG(const LhsEval& xoG, unsigned regionIdx)
1448  {
1449  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1450  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1451 
1452  return xoG*MG / (xoG*(MG - MO) + MO);
1453  }
1454 
1458  template <class LhsEval>
1459  static LhsEval convertXgOToxgO(const LhsEval& XgO, unsigned regionIdx)
1460  {
1461  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1462  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1463 
1464  return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1465  }
1466 
1470  template <class LhsEval>
1471  static LhsEval convertxgOToXgO(const LhsEval& xgO, unsigned regionIdx)
1472  {
1473  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1474  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1475 
1476  return xgO*MO / (xgO*(MO - MG) + MG);
1477  }
1478 
1486  static const GasPvt& gasPvt()
1487  { return *gasPvt_; }
1488 
1496  static const OilPvt& oilPvt()
1497  { return *oilPvt_; }
1498 
1506  static const WaterPvt& waterPvt()
1507  { return *waterPvt_; }
1508 
1514  static Scalar reservoirTemperature(unsigned = 0)
1515  { return reservoirTemperature_; }
1516 
1522  static void setReservoirTemperature(Scalar value)
1523  { reservoirTemperature_ = value; }
1524 
1525  static short activeToCanonicalPhaseIdx(unsigned activePhaseIdx) {
1526  assert(activePhaseIdx<numActivePhases());
1527  return activeToCanonicalPhaseIdx_[activePhaseIdx];
1528  }
1529 
1530  static short canonicalToActivePhaseIdx(unsigned phaseIdx) {
1531  assert(phaseIdx<numPhases);
1532  assert(phaseIsActive(phaseIdx));
1533  return canonicalToActivePhaseIdx_[phaseIdx];
1534  }
1535 
1537  static Scalar diffusionCoefficient(unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1538  { return diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx]; }
1539 
1541  static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1542  { diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx] = coefficient ; }
1543 
1547  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
1548  static LhsEval diffusionCoefficient(const FluidState& fluidState,
1549  const ParameterCache<ParamCacheEval>& paramCache,
1550  unsigned phaseIdx,
1551  unsigned compIdx)
1552  {
1553  // diffusion is disabled by the user
1554  if(!enableDiffusion())
1555  return 0.0;
1556 
1557  // diffusion coefficients are set, and we use them
1558  if(!diffusionCoefficients_.empty()) {
1559  return diffusionCoefficient(compIdx, phaseIdx, paramCache.regionIndex());
1560  }
1561 
1562  const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1563  const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1564 
1565  switch (phaseIdx) {
1566  case oilPhaseIdx: return oilPvt().diffusionCoefficient(T, p, compIdx);
1567  case gasPhaseIdx: return gasPvt().diffusionCoefficient(T, p, compIdx);
1568  case waterPhaseIdx: return 0.0;
1569  default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1570  }
1571  }
1572 
1573 private:
1574  static void resizeArrays_(size_t numRegions)
1575  {
1576  molarMass_.resize(numRegions);
1577  referenceDensity_.resize(numRegions);
1578  }
1579 
1580  static Scalar reservoirTemperature_;
1581 
1582  static std::shared_ptr<GasPvt> gasPvt_;
1583  static std::shared_ptr<OilPvt> oilPvt_;
1584  static std::shared_ptr<WaterPvt> waterPvt_;
1585 
1586  static bool enableDissolvedGas_;
1587  static bool enableVaporizedOil_;
1588  static bool enableVaporizedWater_;
1589  static bool enableDiffusion_;
1590 
1591  // HACK for GCC 4.4: the array size has to be specified using the literal value '3'
1592  // here, because GCC 4.4 seems to be unable to determine the number of phases from
1593  // the BlackOil fluid system in the attribute declaration below...
1594  static std::vector<std::array<Scalar, /*numPhases=*/3> > referenceDensity_;
1595  static std::vector<std::array<Scalar, /*numComponents=*/3> > molarMass_;
1596  static std::vector<std::array<Scalar, /*numComponents=*/3 * /*numPhases=*/3> > diffusionCoefficients_;
1597 
1598  static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1599  static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1600 
1601  static bool isInitialized_;
1602 };
1603 
1604 template <class Scalar, class IndexTraits>
1605 unsigned char BlackOilFluidSystem<Scalar, IndexTraits>::numActivePhases_;
1606 
1607 template <class Scalar, class IndexTraits>
1608 std::array<bool, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::phaseIsActive_;
1609 
1610 template <class Scalar, class IndexTraits>
1611 std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::activeToCanonicalPhaseIdx_;
1612 
1613 template <class Scalar, class IndexTraits>
1614 std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::canonicalToActivePhaseIdx_;
1615 
1616 template <class Scalar, class IndexTraits>
1617 Scalar
1619 
1620 template <class Scalar, class IndexTraits>
1621 Scalar
1623 
1624 template <class Scalar, class IndexTraits>
1625 Scalar
1626 BlackOilFluidSystem<Scalar, IndexTraits>::reservoirTemperature_;
1627 
1628 template <class Scalar, class IndexTraits>
1629 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGas_;
1630 
1631 template <class Scalar, class IndexTraits>
1632 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedOil_;
1633 
1634 template <class Scalar, class IndexTraits>
1635 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedWater_;
1636 
1637 template <class Scalar, class IndexTraits>
1638 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDiffusion_;
1639 
1640 template <class Scalar, class IndexTraits>
1641 std::shared_ptr<OilPvtMultiplexer<Scalar> >
1642 BlackOilFluidSystem<Scalar, IndexTraits>::oilPvt_;
1643 
1644 template <class Scalar, class IndexTraits>
1645 std::shared_ptr<GasPvtMultiplexer<Scalar> >
1646 BlackOilFluidSystem<Scalar, IndexTraits>::gasPvt_;
1647 
1648 template <class Scalar, class IndexTraits>
1649 std::shared_ptr<WaterPvtMultiplexer<Scalar> >
1650 BlackOilFluidSystem<Scalar, IndexTraits>::waterPvt_;
1651 
1652 template <class Scalar, class IndexTraits>
1653 std::vector<std::array<Scalar, 3> >
1654 BlackOilFluidSystem<Scalar, IndexTraits>::referenceDensity_;
1655 
1656 template <class Scalar, class IndexTraits>
1657 std::vector<std::array<Scalar, 3> >
1658 BlackOilFluidSystem<Scalar, IndexTraits>::molarMass_;
1659 
1660 template <class Scalar, class IndexTraits>
1661 std::vector<std::array<Scalar, 9> >
1662 BlackOilFluidSystem<Scalar, IndexTraits>::diffusionCoefficients_;
1663 
1664 template <class Scalar, class IndexTraits>
1665 bool BlackOilFluidSystem<Scalar, IndexTraits>::isInitialized_ = false;
1666 
1667 } // namespace Opm
1668 
1669 #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 constexpr unsigned waterCompIdx
Index of the water component.
Definition: BlackOilFluidSystem.hpp:515
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:1537
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:1354
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:1394
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition: BlackOilFluidSystem.hpp:510
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:1422
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:962
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1435
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: BlackOilFluidSystem.hpp:1105
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1514
static void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:379
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition: BlackOilFluidSystem.hpp:1301
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:788
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 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:1486
static constexpr unsigned numPhases
Number of fluid phases in the fluid system.
Definition: BlackOilFluidSystem.hpp:467
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:1522
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:1250
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:1186
static constexpr unsigned oilCompIdx
Index of the oil component.
Definition: BlackOilFluidSystem.hpp:513
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:983
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1447
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:1279
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition: BlackOilFluidSystem.hpp:1312
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:1329
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1459
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:1548
static LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition: BlackOilFluidSystem.hpp:1225
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 constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition: BlackOilFluidSystem.hpp:470
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1471
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:878
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition: BlackOilFluidSystem.hpp:474
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 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:1367
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 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 constexpr unsigned gasCompIdx
Index of the gas component.
Definition: BlackOilFluidSystem.hpp:517
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:1408
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:1506
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 constexpr unsigned oilPhaseIdx
Index of the oil phase.
Definition: BlackOilFluidSystem.hpp:472
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:1380
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 void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition: BlackOilFluidSystem.hpp:1541
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition: BlackOilFluidSystem.hpp:1496
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:100
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:322
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