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 
33 
36 
38 #include <opm/common/Valgrind.hpp>
40 #include <opm/common/Exceptions.hpp>
41 #include <opm/common/ErrorMacros.hpp>
42 
43 #include <memory>
44 #include <vector>
45 #include <array>
46 
47 namespace Opm {
48 namespace BlackOil {
49 OPM_GENERATE_HAS_MEMBER(Rs, ) // Creates 'HasMember_Rs<T>'.
50 OPM_GENERATE_HAS_MEMBER(Rv, ) // Creates 'HasMember_Rv<T>'.
51 
52 template <class FluidSystem, class LhsEval, class FluidState>
53 LhsEval getRs_(typename std::enable_if<!HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
54  unsigned regionIdx)
55 {
56  const auto& XoG =
57  Opm::decay<LhsEval>(fluidState.massFraction(FluidSystem::oilPhaseIdx,
58  FluidSystem::gasCompIdx));
59  return FluidSystem::convertXoGToRs(XoG, regionIdx);
60 }
61 
62 template <class FluidSystem, class LhsEval, class FluidState>
63 auto getRs_(typename std::enable_if<HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
64  unsigned regionIdx OPM_UNUSED)
66  ::template decay<LhsEval>(fluidState.Rs()))
67 {
68  return Opm::decay<LhsEval>(fluidState.Rs());
69 }
70 
71 template <class FluidSystem, class LhsEval, class FluidState>
72 LhsEval getRv_(typename std::enable_if<!HasMember_Rv<FluidState>::value, const FluidState&>::type fluidState,
73  unsigned regionIdx)
74 {
75  const auto& XgO =
76  Opm::decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx,
77  FluidSystem::oilCompIdx));
78  return FluidSystem::convertXgOToRv(XgO, regionIdx);
79 }
80 
81 template <class FluidSystem, class LhsEval, class FluidState>
82 auto getRv_(typename std::enable_if<HasMember_Rv<FluidState>::value, const FluidState&>::type fluidState,
83  unsigned regionIdx OPM_UNUSED)
85  ::template decay<LhsEval>(fluidState.Rv()))
86 {
87  return Opm::decay<LhsEval>(fluidState.Rv());
88 }
89 }
90 
91 namespace FluidSystems {
92 
99 template <class Scalar>
100 class BlackOil : public BaseFluidSystem<Scalar, BlackOil<Scalar> >
101 {
102  typedef BlackOil ThisType;
103 
104 public:
108 
110  template <class EvaluationT>
111  struct ParameterCache : public Opm::NullParameterCache<EvaluationT>
112  {
113  typedef EvaluationT Evaluation;
114 
115  public:
116  ParameterCache(Scalar maxOilSat = 1.0, unsigned regionIdx=0)
117  {
118  maxOilSat_ = maxOilSat;
119  regionIdx_ = regionIdx;
120  }
121 
129  template <class OtherCache>
130  void assignPersistentData(const OtherCache& other)
131  {
132  regionIdx_ = other.regionIndex();
133  maxOilSat_ = other.maxOilSat();
134  }
135 
143  unsigned regionIndex() const
144  { return regionIdx_; }
145 
153  void setRegionIndex(unsigned val)
154  { regionIdx_ = val; }
155 
156  const Scalar maxOilSat() const
157  { return maxOilSat_; }
158 
159  void setMaxOilSat(Scalar val)
160  { maxOilSat_ = val; }
161 
162  private:
163  Scalar maxOilSat_;
164  unsigned regionIdx_;
165  };
166 
167  /****************************************
168  * Initialization
169  ****************************************/
170 #if HAVE_OPM_PARSER
171 
174  static void initFromDeck(const Deck& deck, const EclipseState& eclState)
175  {
176  auto densityKeyword = deck.getKeyword("DENSITY");
177  size_t numRegions = densityKeyword.size();
179 
180  numActivePhases_ = 0;
181  std::fill(&phaseIsActive_[0], &phaseIsActive_[numPhases], false);
182 
183  if (deck.hasKeyword("OIL")) {
184  phaseIsActive_[oilPhaseIdx] = true;
185  ++ numActivePhases_;
186  }
187 
188  if (deck.hasKeyword("GAS")) {
189  phaseIsActive_[gasPhaseIdx] = true;
190  ++ numActivePhases_;
191  }
192 
193  if (deck.hasKeyword("WATER")) {
194  phaseIsActive_[waterPhaseIdx] = true;
195  ++ numActivePhases_;
196  }
197 
198  // The reservoir temperature does not really belong into the table manager. TODO:
199  // change this in opm-parser
200  setReservoirTemperature(eclState.getTableManager().rtemp());
201 
202  // this fluidsystem only supports two or three phases
203  assert(numActivePhases_ >= 2 && numActivePhases_ <= 3);
204 
205  setEnableDissolvedGas(deck.hasKeyword("DISGAS"));
206  setEnableVaporizedOil(deck.hasKeyword("VAPOIL"));
207 
208  // set the reference densities of all PVT regions
209  for (unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
210  const auto& densityRecord = densityKeyword.getRecord(regionIdx);
211  setReferenceDensities(densityRecord.getItem("OIL").getSIDouble(0),
212  densityRecord.getItem("WATER").getSIDouble(0),
213  densityRecord.getItem("GAS").getSIDouble(0),
214  regionIdx);
215  }
216 
217  if (phaseIsActive(gasPhaseIdx)) {
218  gasPvt_ = std::make_shared<GasPvt>();
219  gasPvt_->initFromDeck(deck, eclState);
220  }
221 
222  if (phaseIsActive(oilPhaseIdx)) {
223  oilPvt_ = std::make_shared<OilPvt>();
224  oilPvt_->initFromDeck(deck, eclState);
225  }
226 
228  waterPvt_ = std::make_shared<WaterPvt>();
229  waterPvt_->initFromDeck(deck, eclState);
230  }
231 
232  initEnd();
233  }
234 #endif // HAVE_OPM_PARSER
235 
244  static void initBegin(size_t numPvtRegions)
245  {
246  enableDissolvedGas_ = true;
247  enableVaporizedOil_ = false;
248 
249  numActivePhases_ = numPhases;
250  std::fill(&phaseIsActive_[0], &phaseIsActive_[numPhases], true);
251 
252  resizeArrays_(numPvtRegions);
253 
255  }
256 
263  static void setEnableDissolvedGas(bool yesno)
264  { enableDissolvedGas_ = yesno; }
265 
272  static void setEnableVaporizedOil(bool yesno)
273  { enableVaporizedOil_ = yesno; }
274 
278  static void setGasPvt(std::shared_ptr<GasPvt> pvtObj)
279  { gasPvt_ = pvtObj; }
280 
284  static void setOilPvt(std::shared_ptr<OilPvt> pvtObj)
285  { oilPvt_ = pvtObj; }
286 
290  static void setWaterPvt(std::shared_ptr<WaterPvt> pvtObj)
291  { waterPvt_ = pvtObj; }
292 
300  static void setReferenceDensities(Scalar rhoOil,
301  Scalar rhoWater,
302  Scalar rhoGas,
303  unsigned regionIdx)
304  {
305  referenceDensity_[regionIdx][oilPhaseIdx] = rhoOil;
306  referenceDensity_[regionIdx][waterPhaseIdx] = rhoWater;
307  referenceDensity_[regionIdx][gasPhaseIdx] = rhoGas;
308  }
309 
313  static void initEnd()
314  {
315  // calculate the final 2D functions which are used for interpolation.
316  size_t numRegions = molarMass_.size();
317  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
318  // calculate molar masses
319 
320  // water is simple: 18 g/mol
321  molarMass_[regionIdx][waterCompIdx] = 18e-3;
322 
323  if (phaseIsActive(gasPhaseIdx)) {
324  // for gas, we take the density at standard conditions and assume it to be ideal
327  Scalar rho_g = referenceDensity_[/*regionIdx=*/0][gasPhaseIdx];
328  molarMass_[regionIdx][gasCompIdx] = Opm::Constants<Scalar>::R*T*rho_g / p;
329  }
330  else
331  // hydrogen gas. we just set this do avoid NaNs later
332  molarMass_[regionIdx][gasCompIdx] = 2e-3;
333 
334  // finally, for oil phase, we take the molar mass from the spe9 paper
335  molarMass_[regionIdx][oilCompIdx] = 175e-3; // kg/mol
336  }
337 
338  isInitialized_ = true;
339  }
340 
341  static bool isInitialized()
342  { return isInitialized_; }
343 
344  /****************************************
345  * Generic phase properties
346  ****************************************/
347 
349  static const unsigned numPhases = 3;
350 
352  static const unsigned waterPhaseIdx = 0;
354  static const unsigned oilPhaseIdx = 1;
356  static const unsigned gasPhaseIdx = 2;
357 
359  static const Scalar surfacePressure;
360 
363 
365  static const char* phaseName(unsigned phaseIdx)
366  {
367  static const char* name[] = { "water", "oil", "gas" };
368 
369  assert(0 <= phaseIdx && phaseIdx < numPhases + 1);
370  return name[phaseIdx];
371  }
372 
374  static bool isLiquid(unsigned phaseIdx)
375  {
376  assert(0 <= phaseIdx && phaseIdx < numPhases);
377  return phaseIdx != gasPhaseIdx;
378  }
379 
380  /****************************************
381  * Generic component related properties
382  ****************************************/
383 
385  static const unsigned numComponents = 3;
386 
388  static const unsigned oilCompIdx = 0;
390  static const unsigned waterCompIdx = 1;
392  static const unsigned gasCompIdx = 2;
393 
394 protected:
395  static const int phaseToSolventCompIdx_[3];
396  static const int phaseToSoluteCompIdx_[3];
397 
398  static unsigned char numActivePhases_;
399  static bool phaseIsActive_[numPhases];
400 
401 public:
403  static unsigned numActivePhases()
404  { return numActivePhases_; }
405 
407  static unsigned phaseIsActive(unsigned phaseIdx)
408  {
409  assert(phaseIdx < numPhases);
410  return phaseIsActive_[phaseIdx];
411  }
412 
414  static constexpr unsigned solventComponentIndex(unsigned phaseIdx)
415  { return static_cast<unsigned>(phaseToSolventCompIdx_[phaseIdx]); }
416 
418  static constexpr unsigned soluteComponentIndex(unsigned phaseIdx)
419  { return static_cast<unsigned>(phaseToSoluteCompIdx_[phaseIdx]); }
420 
422  static const char* componentName(unsigned compIdx)
423  {
424  static const char* name[] = { "Oil", "Water", "Gas" };
425 
426  assert(0 <= compIdx && compIdx < numComponents);
427  return name[compIdx];
428  }
429 
431  static Scalar molarMass(unsigned compIdx, unsigned regionIdx = 0)
432  { return molarMass_[regionIdx][compIdx]; }
433 
435  static bool isIdealMixture(unsigned /*phaseIdx*/)
436  {
437  // fugacity coefficients are only pressure dependent -> we
438  // have an ideal mixture
439  return true;
440  }
441 
443  static bool isCompressible(unsigned /*phaseIdx*/)
444  { return true; /* all phases are compressible */ }
445 
447  static bool isIdealGas(unsigned /*phaseIdx*/)
448  { return false; }
449 
450 
451  /****************************************
452  * Black-oil specific properties
453  ****************************************/
459  static size_t numRegions()
460  { return molarMass_.size(); }
461 
468  static bool enableDissolvedGas()
469  { return enableDissolvedGas_; }
470 
477  static bool enableVaporizedOil()
478  { return enableVaporizedOil_; }
479 
485  static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
486  { return referenceDensity_[regionIdx][phaseIdx]; }
487 
488  /****************************************
489  * thermodynamic quantities (generic version, only isothermal)
490  ****************************************/
492  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
493  static LhsEval density(const FluidState& fluidState,
494  const ParameterCache<ParamCacheEval>& paramCache,
495  unsigned phaseIdx)
496  { return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
497 
499  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
500  static LhsEval fugacityCoefficient(const FluidState& fluidState,
501  const ParameterCache<ParamCacheEval>& paramCache,
502  unsigned phaseIdx,
503  unsigned compIdx)
504  {
505  return fugacityCoefficient<FluidState, LhsEval>(fluidState,
506  phaseIdx,
507  compIdx,
508  paramCache.regionIndex());
509  }
510 
512  template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
513  static LhsEval viscosity(const FluidState& fluidState,
514  const ParameterCache<ParamCacheEval>& paramCache,
515  unsigned phaseIdx)
516  { return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
517 
518 
519  /****************************************
520  * thermodynamic quantities (black-oil specific version: Note that the PVT region
521  * index is explicitly passed instead of a parameter cache object)
522  ****************************************/
524  template <class FluidState, class LhsEval = typename FluidState::Scalar>
525  static LhsEval density(const FluidState& fluidState,
526  unsigned phaseIdx,
527  unsigned regionIdx)
528  {
529  assert(0 <= phaseIdx && phaseIdx <= numPhases);
530  assert(0 <= regionIdx && regionIdx <= numRegions());
531 
532  switch (phaseIdx) {
533  case oilPhaseIdx: {
534  if (!enableDissolvedGas()) {
535  // immiscible oil
536  const auto& bo = inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
537  return referenceDensity(phaseIdx, regionIdx)*bo;
538  }
539 
540  // miscible oil
541  const auto& bo = inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
542  const auto& Rs = Opm::BlackOil::template getRs_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
543 
544  return
545  bo*referenceDensity(oilPhaseIdx, regionIdx)
546  + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
547  }
548 
549  case gasPhaseIdx: {
550  if (!enableVaporizedOil()) {
551  // immiscible gas
552  const auto& bg = inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
553  return bg*referenceDensity(phaseIdx, regionIdx);
554  }
555 
556  // miscible gas
557  const auto& bg = inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
558  const auto& Rv = Opm::BlackOil::template getRv_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
559 
560  return
561  bg*referenceDensity(gasPhaseIdx, regionIdx)
562  + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
563  }
564 
565  case waterPhaseIdx:
566  return
567  referenceDensity(waterPhaseIdx, regionIdx)
568  *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
569  }
570 
571  OPM_THROW(std::logic_error, "Unhandled phase index " << phaseIdx);
572  }
573 
581  template <class FluidState, class LhsEval = typename FluidState::Scalar>
582  static LhsEval saturatedDensity(const FluidState& fluidState,
583  unsigned phaseIdx,
584  unsigned regionIdx)
585  {
586  assert(0 <= phaseIdx && phaseIdx <= numPhases);
587  assert(0 <= regionIdx && regionIdx <= numRegions());
588 
589  switch (phaseIdx) {
590  case oilPhaseIdx: {
591  if (!enableDissolvedGas()) {
592  // immiscible oil
593  const auto& bo = inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
594  return referenceDensity(phaseIdx, regionIdx)*bo;
595  }
596 
597  // miscible oil
598  const auto& bo = saturatedInverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
599  const auto& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
600 
601  return
602  bo*referenceDensity(oilPhaseIdx, regionIdx)
603  + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
604  }
605 
606  case gasPhaseIdx: {
607  if (!enableVaporizedOil()) {
608  // immiscible gas
609  const auto& bg = inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
610  return referenceDensity(phaseIdx, regionIdx)*bg;
611  }
612 
613  // miscible gas
614  const auto& bg = saturatedInverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
615  const auto& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
616 
617  return
618  bg*referenceDensity(gasPhaseIdx, regionIdx)
619  + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
620  }
621 
622  case waterPhaseIdx:
623  return
624  referenceDensity(waterPhaseIdx, regionIdx)
625  *saturatedInverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
626  }
627 
628  OPM_THROW(std::logic_error, "Unhandled phase index " << phaseIdx);
629  }
630 
639  template <class FluidState, class LhsEval = typename FluidState::Scalar>
640  static LhsEval inverseFormationVolumeFactor(const FluidState& fluidState,
641  unsigned phaseIdx,
642  unsigned regionIdx)
643  {
644  assert(0 <= phaseIdx && phaseIdx <= numPhases);
645  assert(0 <= regionIdx && regionIdx <= numRegions());
646 
647  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
648  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
649 
650  switch (phaseIdx) {
651  case oilPhaseIdx: {
652  if (enableDissolvedGas()) {
653  if (fluidState.phaseIsPresent(gasPhaseIdx)) {
654  if (fluidState.saturation(gasPhaseIdx) < 1e-4) {
655  // here comes the relatively expensive case: first calculate and then
656  // interpolate between the saturated and undersaturated quantities to
657  // avoid a discontinuity
658  const auto& Rs = Opm::BlackOil::template getRs_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
659  const auto& alpha = Opm::decay<LhsEval>(fluidState.saturation(gasPhaseIdx))/1e-4;
660  const auto& bSat = oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
661  const auto& bUndersat = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
662  return alpha*bSat + (1.0 - alpha)*bUndersat;
663  }
664 
665  return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
666  }
667 
668  const auto& Rs = Opm::BlackOil::template getRs_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
669  return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
670  }
671 
672  const LhsEval Rs(0.0);
673  return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
674  }
675  case gasPhaseIdx: {
676  if (enableVaporizedOil()) {
677  if (fluidState.phaseIsPresent(oilPhaseIdx)) {
678  if (fluidState.saturation(oilPhaseIdx) < 1e-4) {
679  // here comes the relatively expensive case: first calculate and then
680  // interpolate between the saturated and undersaturated quantities to
681  // avoid a discontinuity
682  const auto& Rv = Opm::BlackOil::template getRv_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
683  const auto& alpha = Opm::decay<LhsEval>(fluidState.saturation(oilPhaseIdx))/1e-4;
684  const auto& bSat = gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
685  const auto& bUndersat = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
686  return alpha*bSat + (1.0 - alpha)*bUndersat;
687  }
688 
689  return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
690  }
691 
692  const auto& Rv = Opm::BlackOil::template getRv_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
693  return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
694  }
695 
696  const LhsEval Rv(0.0);
697  return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
698  }
699  case waterPhaseIdx:
700  return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p);
701  default: OPM_THROW(std::logic_error, "Unhandled phase index " << phaseIdx);
702  }
703  }
704 
712  template <class FluidState, class LhsEval = typename FluidState::Scalar>
713  static LhsEval saturatedInverseFormationVolumeFactor(const FluidState& fluidState,
714  unsigned phaseIdx,
715  unsigned regionIdx)
716  {
717  assert(0 <= phaseIdx && phaseIdx <= numPhases);
718  assert(0 <= regionIdx && regionIdx <= numRegions());
719 
720  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
721  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
722 
723  switch (phaseIdx) {
724  case oilPhaseIdx: return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
725  case gasPhaseIdx: return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
726  case waterPhaseIdx: return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p);
727  default: OPM_THROW(std::logic_error, "Unhandled phase index " << phaseIdx);
728  }
729  }
730 
732  template <class FluidState, class LhsEval = typename FluidState::Scalar>
733  static LhsEval fugacityCoefficient(const FluidState& fluidState,
734  unsigned phaseIdx,
735  unsigned compIdx,
736  unsigned regionIdx)
737  {
738  assert(0 <= phaseIdx && phaseIdx <= numPhases);
739  assert(0 <= compIdx && compIdx <= numComponents);
740  assert(0 <= regionIdx && regionIdx <= numRegions());
741 
742  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
743  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
744 
745  // for the fugacity coefficient of the oil component in the oil phase, we use
746  // some pseudo-realistic value for the vapor pressure to ease physical
747  // interpretation of the results
748  const LhsEval phi_oO = 20e3/p;
749 
750  // for the gas component in the gas phase, assume it to be an ideal gas
751  const Scalar phi_gG = 1.0;
752 
753  // for the fugacity coefficient of the water component in the water phase, we use
754  // the same approach as for the oil component in the oil phase
755  const LhsEval phi_wW = 30e3/p;
756 
757  switch (phaseIdx) {
758  case gasPhaseIdx: // fugacity coefficients for all components in the gas phase
759  switch (compIdx) {
760  case gasCompIdx:
761  return phi_gG;
762 
763  // for the oil component, we calculate the Rv value for saturated gas and Rs
764  // for saturated oil, and compute the fugacity coefficients at the
765  // equilibrium. for this, we assume that in equilibrium the fugacities of the
766  // oil component is the same in both phases.
767  case oilCompIdx: {
768  if (!enableVaporizedOil())
769  // if there's no vaporized oil, the gas phase is assumed to be
770  // immiscible with the oil component
771  return phi_gG*1e6;
772 
773  const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
774  const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
775  const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
776 
777  const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
778  const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
779  const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
780  const auto& x_oOSat = 1.0 - x_oGSat;
781 
782  const auto& p_o = Opm::decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
783  const auto& p_g = Opm::decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
784 
785  return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
786  }
787 
788  case waterCompIdx:
789  // the water component is assumed to be never miscible with the gas phase
790  return phi_gG*1e6;
791 
792  default:
793  OPM_THROW(std::logic_error,
794  "Invalid component index " << compIdx);
795  }
796 
797  case oilPhaseIdx: // fugacity coefficients for all components in the oil phase
798  switch (compIdx) {
799  case oilCompIdx:
800  return phi_oO;
801 
802  // for the oil and water components, we proceed analogous to the gas and
803  // water components in the gas phase
804  case gasCompIdx: {
805  if (!enableDissolvedGas())
806  // if there's no dissolved gas, the oil phase is assumed to be
807  // immiscible with the gas component
808  return phi_oO*1e6;
809 
810  const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
811  const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
812  const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
813  const auto& x_gGSat = 1.0 - x_gOSat;
814 
815  const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
816  const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
817  const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
818 
819  const auto& p_o = Opm::decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
820  const auto& p_g = Opm::decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
821 
822  return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
823  }
824 
825  case waterCompIdx:
826  return phi_oO*1e6;
827 
828  default:
829  OPM_THROW(std::logic_error,
830  "Invalid component index " << compIdx);
831  }
832 
833  case waterPhaseIdx: // fugacity coefficients for all components in the water phase
834  // the water phase fugacity coefficients are pretty simple: because the water
835  // phase is assumed to consist entirely from the water component, we just
836  // need to make sure that the fugacity coefficients for the other components
837  // are a few orders of magnitude larger than the one of the water
838  // component. (i.e., the affinity of the gas and oil components to the water
839  // phase is lower by a few orders of magnitude)
840  switch (compIdx) {
841  case waterCompIdx: return phi_wW;
842  case oilCompIdx: return 1.1e6*phi_wW;
843  case gasCompIdx: return 1e6*phi_wW;
844  default:
845  OPM_THROW(std::logic_error,
846  "Invalid component index " << compIdx);
847  }
848 
849  default:
850  OPM_THROW(std::logic_error,
851  "Invalid phase index " << phaseIdx);
852  }
853 
854  OPM_THROW(std::logic_error, "Unhandled phase or component index");
855  }
856 
858  template <class FluidState, class LhsEval = typename FluidState::Scalar>
859  static LhsEval viscosity(const FluidState& fluidState,
860  unsigned phaseIdx,
861  unsigned regionIdx)
862  {
863  assert(0 <= phaseIdx && phaseIdx <= numPhases);
864  assert(0 <= regionIdx && regionIdx <= numRegions());
865 
866  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
867  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
868 
869  switch (phaseIdx) {
870  case oilPhaseIdx: {
871  if (enableDissolvedGas()) {
872  if (fluidState.phaseIsPresent(gasPhaseIdx)) {
873  if (fluidState.saturation(gasPhaseIdx) < 1e-4) {
874  // here comes the relatively expensive case: first calculate and then
875  // interpolate between the saturated and undersaturated quantities to
876  // avoid a discontinuity
877  const auto& Rs = Opm::BlackOil::template getRs_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
878  const auto& alpha = Opm::decay<LhsEval>(fluidState.saturation(gasPhaseIdx))/1e-4;
879  const auto& muSat = oilPvt_->saturatedViscosity(regionIdx, T, p);
880  const auto& muUndersat = oilPvt_->viscosity(regionIdx, T, p, Rs);
881  return alpha*muSat + (1.0 - alpha)*muUndersat;
882  }
883 
884  return oilPvt_->saturatedViscosity(regionIdx, T, p);
885  }
886 
887  const auto& Rs = Opm::BlackOil::template getRs_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
888  return oilPvt_->viscosity(regionIdx, T, p, Rs);
889  }
890 
891  const LhsEval Rs(0.0);
892  return oilPvt_->viscosity(regionIdx, T, p, Rs);
893  }
894 
895  case gasPhaseIdx: {
896  if (enableVaporizedOil()) {
897  if (fluidState.phaseIsPresent(oilPhaseIdx)) {
898  if (fluidState.saturation(oilPhaseIdx) < 1e-4) {
899  // here comes the relatively expensive case: first calculate and then
900  // interpolate between the saturated and undersaturated quantities to
901  // avoid a discontinuity
902  const auto& Rv = Opm::BlackOil::template getRv_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
903  const auto& alpha = Opm::decay<LhsEval>(fluidState.saturation(oilPhaseIdx))/1e-4;
904  const auto& muSat = gasPvt_->saturatedViscosity(regionIdx, T, p);
905  const auto& muUndersat = gasPvt_->viscosity(regionIdx, T, p, Rv);
906  return alpha*muSat + (1.0 - alpha)*muUndersat;
907  }
908 
909  return gasPvt_->saturatedViscosity(regionIdx, T, p);
910  }
911 
912  const auto& Rv = Opm::BlackOil::template getRv_<ThisType, LhsEval, FluidState>(fluidState, regionIdx);
913  return gasPvt_->viscosity(regionIdx, T, p, Rv);
914  }
915 
916  const LhsEval Rv(0.0);
917  return gasPvt_->viscosity(regionIdx, T, p, Rv);
918  }
919 
920  case waterPhaseIdx:
921  // since water is always assumed to be immiscible in the black-oil model,
922  // there is no "saturated water"
923  return waterPvt_->viscosity(regionIdx, T, p);
924  }
925 
926  OPM_THROW(std::logic_error, "Unhandled phase index " << phaseIdx);
927  }
928 
935  template <class FluidState, class LhsEval = typename FluidState::Scalar>
936  static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
937  unsigned phaseIdx,
938  unsigned regionIdx,
939  Scalar maxOilSaturation)
940  {
941  assert(0 <= phaseIdx && phaseIdx <= numPhases);
942  assert(0 <= regionIdx && regionIdx <= numRegions());
943 
944  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
945  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
946  const auto& So = Opm::decay<LhsEval>(fluidState.saturation(oilPhaseIdx));
947 
948  switch (phaseIdx) {
949  case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
950  case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
951  case waterPhaseIdx: return 0.0;
952  default: OPM_THROW(std::logic_error, "Unhandled phase index " << phaseIdx);
953  }
954  }
955 
964  template <class FluidState, class LhsEval = typename FluidState::Scalar>
965  static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
966  unsigned phaseIdx,
967  unsigned regionIdx)
968  {
969  assert(0 <= phaseIdx && phaseIdx <= numPhases);
970  assert(0 <= regionIdx && regionIdx <= numRegions());
971 
972  const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
973  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
974 
975  switch (phaseIdx) {
976  case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
977  case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
978  case waterPhaseIdx: return 0.0;
979  default: OPM_THROW(std::logic_error, "Unhandled phase index " << phaseIdx);
980  }
981  }
982 
986  template <class FluidState, class LhsEval = typename FluidState::Scalar>
987  static LhsEval bubblePointPressure(const FluidState& fluidState,
988  unsigned regionIdx)
989  {
990  return saturationPressure(fluidState, oilPhaseIdx, regionIdx);
991  }
992 
993 
997  template <class FluidState, class LhsEval = typename FluidState::Scalar>
998  static LhsEval dewPointPressure(const FluidState& fluidState,
999  unsigned regionIdx)
1000  {
1001  return saturationPressure(fluidState, gasPhaseIdx, regionIdx);
1002  }
1003 
1014  template <class FluidState, class LhsEval = typename FluidState::Scalar>
1015  static LhsEval saturationPressure(const FluidState& fluidState,
1016  unsigned phaseIdx,
1017  unsigned regionIdx)
1018  {
1019  assert(0 <= phaseIdx && phaseIdx <= numPhases);
1020  assert(0 <= regionIdx && regionIdx <= numRegions());
1021 
1022  const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
1023 
1024  switch (phaseIdx) {
1025  case oilPhaseIdx: return oilPvt_->saturationPressure(regionIdx, T, Opm::BlackOil::template getRs_<ThisType, LhsEval, FluidState>(fluidState, regionIdx));
1026  case gasPhaseIdx: return gasPvt_->saturationPressure(regionIdx, T, Opm::BlackOil::template getRv_<ThisType, LhsEval, FluidState>(fluidState, regionIdx));
1027  case waterPhaseIdx: return 0.0;
1028  default: OPM_THROW(std::logic_error, "Unhandled phase index " << phaseIdx);
1029  }
1030  }
1031 
1032  /****************************************
1033  * Auxiliary and convenience methods for the black-oil model
1034  ****************************************/
1039  template <class LhsEval>
1040  static LhsEval convertXoGToRs(const LhsEval& XoG, unsigned regionIdx)
1041  {
1042  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1043  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1044 
1045  return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1046  }
1047 
1052  template <class LhsEval>
1053  static LhsEval convertXgOToRv(const LhsEval& XgO, unsigned regionIdx)
1054  {
1055  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1056  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1057 
1058  return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1059  }
1060 
1065  template <class LhsEval>
1066  static LhsEval convertRsToXoG(const LhsEval& Rs, unsigned regionIdx)
1067  {
1068  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1069  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1070 
1071  const LhsEval& rho_oG = Rs*rho_gRef;
1072  return rho_oG/(rho_oRef + rho_oG);
1073  }
1074 
1079  template <class LhsEval>
1080  static LhsEval convertRvToXgO(const LhsEval& Rv, unsigned regionIdx)
1081  {
1082  Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1083  Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1084 
1085  const LhsEval& rho_gO = Rv*rho_oRef;
1086  return rho_gO/(rho_gRef + rho_gO);
1087  }
1088 
1092  template <class LhsEval>
1093  static LhsEval convertXoGToxoG(const LhsEval& XoG, unsigned regionIdx)
1094  {
1095  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1096  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1097 
1098  return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1099  }
1100 
1104  template <class LhsEval>
1105  static LhsEval convertxoGToXoG(const LhsEval& xoG, unsigned regionIdx)
1106  {
1107  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1108  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1109 
1110  return xoG*MG / (xoG*(MG - MO) + MO);
1111  }
1112 
1116  template <class LhsEval>
1117  static LhsEval convertXgOToxgO(const LhsEval& XgO, unsigned regionIdx)
1118  {
1119  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1120  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1121 
1122  return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1123  }
1124 
1128  template <class LhsEval>
1129  static LhsEval convertxgOToXgO(const LhsEval& xgO, unsigned regionIdx)
1130  {
1131  Scalar MO = molarMass_[regionIdx][oilCompIdx];
1132  Scalar MG = molarMass_[regionIdx][gasCompIdx];
1133 
1134  return xgO*MO / (xgO*(MO - MG) + MG);
1135  }
1136 
1144  static const GasPvt& gasPvt()
1145  { return *gasPvt_; }
1146 
1154  static const OilPvt& oilPvt()
1155  { return *oilPvt_; }
1156 
1164  static const WaterPvt& waterPvt()
1165  { return *waterPvt_; }
1166 
1172  static Scalar reservoirTemperature(unsigned pvtRegionIdx OPM_UNUSED = 0)
1173  { return reservoirTemperature_; }
1174 
1180  static void setReservoirTemperature(Scalar value)
1181  { reservoirTemperature_ = value; }
1182 
1183 private:
1184  static void resizeArrays_(size_t numRegions)
1185  {
1186  molarMass_.resize(numRegions);
1187  referenceDensity_.resize(numRegions);
1188  }
1189 
1190  static Scalar reservoirTemperature_;
1191 
1192  static std::shared_ptr<GasPvt> gasPvt_;
1193  static std::shared_ptr<OilPvt> oilPvt_;
1194  static std::shared_ptr<WaterPvt> waterPvt_;
1195 
1196  static bool enableDissolvedGas_;
1197  static bool enableVaporizedOil_;
1198 
1199  // HACK for GCC 4.4: the array size has to be specified using the literal value '3'
1200  // here, because GCC 4.4 seems to be unable to determine the number of phases from
1201  // the BlackOil fluid system in the attribute declaration below...
1202  static std::vector<std::array<Scalar, /*numPhases=*/3> > referenceDensity_;
1203  static std::vector<std::array<Scalar, /*numComponents=*/3> > molarMass_;
1204 
1205  static bool isInitialized_;
1206 };
1207 
1208 template <class Scalar>
1209 const int BlackOil<Scalar>::phaseToSolventCompIdx_[3] =
1210 {
1211  waterCompIdx, // water phase
1212  oilCompIdx, // oil phase
1213  gasCompIdx // gas phase
1214 };
1215 
1216 
1217 template <class Scalar>
1218 const int BlackOil<Scalar>::phaseToSoluteCompIdx_[3] =
1219 {
1220  -1, // water phase
1221  gasCompIdx, // oil phase
1222  oilCompIdx // gas phase
1223 };
1224 
1225 template <class Scalar>
1226 unsigned char BlackOil<Scalar>::numActivePhases_;
1227 
1228 template <class Scalar>
1229 bool BlackOil<Scalar>::phaseIsActive_[numPhases];
1230 
1231 template <class Scalar>
1232 const Scalar
1233 BlackOil<Scalar>::surfaceTemperature = 273.15 + 15.56; // [K]
1234 
1235 template <class Scalar>
1236 const Scalar
1237 BlackOil<Scalar>::surfacePressure = 101325.0; // [Pa]
1238 
1239 template <class Scalar>
1240 Scalar
1241 BlackOil<Scalar>::reservoirTemperature_;
1242 
1243 template <class Scalar>
1244 bool BlackOil<Scalar>::enableDissolvedGas_;
1245 
1246 template <class Scalar>
1247 bool BlackOil<Scalar>::enableVaporizedOil_;
1248 
1249 template <class Scalar>
1250 std::shared_ptr<OilPvtMultiplexer<Scalar> >
1251 BlackOil<Scalar>::oilPvt_;
1252 
1253 template <class Scalar>
1254 std::shared_ptr<Opm::GasPvtMultiplexer<Scalar> >
1255 BlackOil<Scalar>::gasPvt_;
1256 
1257 template <class Scalar>
1258 std::shared_ptr<WaterPvtMultiplexer<Scalar> >
1259 BlackOil<Scalar>::waterPvt_;
1260 
1261 template <class Scalar>
1262 std::vector<std::array<Scalar, 3> >
1263 BlackOil<Scalar>::referenceDensity_;
1264 
1265 template <class Scalar>
1266 std::vector<std::array<Scalar, 3> >
1267 BlackOil<Scalar>::molarMass_;
1268 
1269 template <class Scalar>
1270 bool BlackOil<Scalar>::isInitialized_ = false;
1271 
1272 }} // namespace Opm, FluidSystems
1273 
1274 #endif
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx, Scalar maxOilSaturation)
Returns the dissolution factor of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:936
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition: BlackOilFluidSystem.hpp:1154
unsigned regionIndex() const
Return the index of the region which should be used to determine the thermodynamic properties...
Definition: BlackOilFluidSystem.hpp:143
static void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition: BlackOilFluidSystem.hpp:300
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1180
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
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:263
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1117
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:640
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition: BlackOilFluidSystem.hpp:987
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition: BlackOilFluidSystem.hpp:365
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition: BlackOilFluidSystem.hpp:468
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:130
void setRegionIndex(unsigned val)
Set the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:153
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: BlackOilFluidSystem.hpp:859
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition: BlackOilFluidSystem.hpp:374
static const Scalar surfaceTemperature
The temperature at the surface.
Definition: BlackOilFluidSystem.hpp:362
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition: BlackOilFluidSystem.hpp:403
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:493
Definition: Evaluation.hpp:531
static unsigned phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition: BlackOilFluidSystem.hpp:407
#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
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
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:272
Definition: Air_Mesitylene.hpp:33
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1105
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition: GasPvtMultiplexer.hpp:75
static const unsigned gasCompIdx
Index of the gas component.
Definition: BlackOilFluidSystem.hpp:392
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:513
static constexpr unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition: BlackOilFluidSystem.hpp:414
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition: BlackOilFluidSystem.hpp:477
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:1040
static const unsigned oilCompIdx
Index of the oil component.
Definition: BlackOilFluidSystem.hpp:388
A central place for various physical constants occuring in some equations.
Definition: Constants.hpp:40
This macro generates a class HasMember_${MEMBER_NAME} which can be used for template specialization...
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1093
static const unsigned gasPhaseIdx
Index of the gas phase.
Definition: BlackOilFluidSystem.hpp:356
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition: BlackOilFluidSystem.hpp:431
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition: BlackOilFluidSystem.hpp:485
static const unsigned waterPhaseIdx
Index of the water phase.
Definition: BlackOilFluidSystem.hpp:352
static const unsigned numComponents
Number of chemical species in the fluid system.
Definition: BlackOilFluidSystem.hpp:385
Definition: MathToolbox.hpp:48
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: BlackOilFluidSystem.hpp:525
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
Definition: BlackOilFluidSystem.hpp:100
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:582
static Scalar reservoirTemperature(unsigned pvtRegionIdx OPM_UNUSED=0)
Set the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1172
static size_t numRegions()
Returns the number of PVT regions which are considered.
Definition: BlackOilFluidSystem.hpp:459
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1129
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:1053
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:1066
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
Definition: OilPvtMultiplexer.hpp:75
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition: BlackOilFluidSystem.hpp:284
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition: WaterPvtMultiplexer.hpp:55
The type of the fluid system&#39;s parameter cache.
Definition: BlackOilFluidSystem.hpp:111
The base class for all fluid systems.
Definition: BaseFluidSystem.hpp:43
A central place for various physical constants occuring in some equations.
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition: BlackOilFluidSystem.hpp:1144
static void initBegin(size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition: BlackOilFluidSystem.hpp:244
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition: BlackOilFluidSystem.hpp:290
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:713
static void initEnd()
Finish initializing the black oil fluid system.
Definition: BlackOilFluidSystem.hpp:313
static constexpr unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition: BlackOilFluidSystem.hpp:418
static const unsigned oilPhaseIdx
Index of the oil phase.
Definition: BlackOilFluidSystem.hpp:354
static const Scalar surfacePressure
The pressure at the surface.
Definition: BlackOilFluidSystem.hpp:359
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition: BlackOilFluidSystem.hpp:1164
A parameter cache which does nothing.
Definition: NullParameterCache.hpp:39
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: BlackOilFluidSystem.hpp:435
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:1015
Scalar Scalar
The type used for scalar quantities.
Definition: BaseFluidSystem.hpp:49
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:1080
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:733
static const unsigned numPhases
Number of fluid phases in the fluid system.
Definition: BlackOilFluidSystem.hpp:349
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition: BlackOilFluidSystem.hpp:278
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:965
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: BlackOilFluidSystem.hpp:447
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:500
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: BlackOilFluidSystem.hpp:443
The base class for all fluid systems.
static const unsigned waterCompIdx
Index of the water component.
Definition: BlackOilFluidSystem.hpp:390
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition: BlackOilFluidSystem.hpp:422
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition: BlackOilFluidSystem.hpp:998