My Project
Loading...
Searching...
No Matches
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
36
37#include <opm/common/TimingMacros.hpp>
38
41
45
46#include <array>
47#include <cstddef>
48#include <memory>
49#include <stdexcept>
50#include <string_view>
51#include <vector>
52
53namespace Opm {
54
55#if HAVE_ECL_INPUT
56class EclipseState;
57class Schedule;
58#endif
59
60namespace BlackOil {
61OPM_GENERATE_HAS_MEMBER(Rs, ) // Creates 'HasMember_Rs<T>'.
62OPM_GENERATE_HAS_MEMBER(Rv, ) // Creates 'HasMember_Rv<T>'.
63OPM_GENERATE_HAS_MEMBER(Rvw, ) // Creates 'HasMember_Rvw<T>'.
64OPM_GENERATE_HAS_MEMBER(Rsw, ) // Creates 'HasMember_Rsw<T>'.
65OPM_GENERATE_HAS_MEMBER(saltConcentration, )
66OPM_GENERATE_HAS_MEMBER(saltSaturation, )
67
68template <class FluidSystem, class FluidState, class LhsEval>
69LhsEval getRs_(typename std::enable_if<!HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
70 unsigned regionIdx)
71{
72 const auto& XoG =
73 decay<LhsEval>(fluidState.massFraction(FluidSystem::oilPhaseIdx, FluidSystem::gasCompIdx));
74 return FluidSystem::convertXoGToRs(XoG, regionIdx);
75}
76
77template <class FluidSystem, class FluidState, class LhsEval>
78auto getRs_(typename std::enable_if<HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
79 unsigned)
80 -> decltype(decay<LhsEval>(fluidState.Rs()))
81{ return decay<LhsEval>(fluidState.Rs()); }
82
83template <class FluidSystem, class FluidState, class LhsEval>
84LhsEval getRv_(typename std::enable_if<!HasMember_Rv<FluidState>::value, const FluidState&>::type fluidState,
85 unsigned regionIdx)
86{
87 const auto& XgO =
88 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::oilCompIdx));
89 return FluidSystem::convertXgOToRv(XgO, regionIdx);
90}
91
92template <class FluidSystem, class FluidState, class LhsEval>
93auto getRv_(typename std::enable_if<HasMember_Rv<FluidState>::value, const FluidState&>::type fluidState,
94 unsigned)
95 -> decltype(decay<LhsEval>(fluidState.Rv()))
96{ return decay<LhsEval>(fluidState.Rv()); }
97
98template <class FluidSystem, class FluidState, class LhsEval>
99LhsEval getRvw_(typename std::enable_if<!HasMember_Rvw<FluidState>::value, const FluidState&>::type fluidState,
100 unsigned regionIdx)
101{
102 const auto& XgW =
103 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::waterCompIdx));
104 return FluidSystem::convertXgWToRvw(XgW, regionIdx);
105}
106
107template <class FluidSystem, class FluidState, class LhsEval>
108auto getRvw_(typename std::enable_if<HasMember_Rvw<FluidState>::value, const FluidState&>::type fluidState,
109 unsigned)
110 -> decltype(decay<LhsEval>(fluidState.Rvw()))
111{ return decay<LhsEval>(fluidState.Rvw()); }
112
113template <class FluidSystem, class FluidState, class LhsEval>
114LhsEval getRsw_(typename std::enable_if<!HasMember_Rsw<FluidState>::value, const FluidState&>::type fluidState,
115 unsigned regionIdx)
116{
117 const auto& XwG =
118 decay<LhsEval>(fluidState.massFraction(FluidSystem::waterPhaseIdx, FluidSystem::gasCompIdx));
119 return FluidSystem::convertXwGToRsw(XwG, regionIdx);
120}
121
122template <class FluidSystem, class FluidState, class LhsEval>
123auto getRsw_(typename std::enable_if<HasMember_Rsw<FluidState>::value, const FluidState&>::type fluidState,
124 unsigned)
125 -> decltype(decay<LhsEval>(fluidState.Rsw()))
126{ return decay<LhsEval>(fluidState.Rsw()); }
127
128template <class FluidSystem, class FluidState, class LhsEval>
129LhsEval getSaltConcentration_(typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
130 const FluidState&>::type,
131 unsigned)
132{return 0.0;}
133
134template <class FluidSystem, class FluidState, class LhsEval>
135auto getSaltConcentration_(typename std::enable_if<HasMember_saltConcentration<FluidState>::value, const FluidState&>::type fluidState,
136 unsigned)
137 -> decltype(decay<LhsEval>(fluidState.saltConcentration()))
138{ return decay<LhsEval>(fluidState.saltConcentration()); }
139
140template <class FluidSystem, class FluidState, class LhsEval>
141LhsEval getSaltSaturation_(typename std::enable_if<!HasMember_saltSaturation<FluidState>::value,
142 const FluidState&>::type,
143 unsigned)
144{return 0.0;}
145
146template <class FluidSystem, class FluidState, class LhsEval>
147auto getSaltSaturation_(typename std::enable_if<HasMember_saltSaturation<FluidState>::value, const FluidState&>::type fluidState,
148 unsigned)
149 -> decltype(decay<LhsEval>(fluidState.saltSaturation()))
150{ return decay<LhsEval>(fluidState.saltSaturation()); }
151
152}
153
160template <class Scalar, class IndexTraits = BlackOilDefaultIndexTraits>
161class BlackOilFluidSystem : public BaseFluidSystem<Scalar, BlackOilFluidSystem<Scalar, IndexTraits> >
162{
164
165public:
169
171 template <class EvaluationT>
172 struct ParameterCache : public NullParameterCache<EvaluationT>
173 {
174 using Evaluation = EvaluationT;
175
176 public:
177 ParameterCache(Scalar maxOilSat = 1.0, unsigned regionIdx=0)
178 {
179 maxOilSat_ = maxOilSat;
180 regionIdx_ = regionIdx;
181 }
182
190 template <class OtherCache>
191 void assignPersistentData(const OtherCache& other)
192 {
193 regionIdx_ = other.regionIndex();
194 maxOilSat_ = other.maxOilSat();
195 }
196
204 unsigned regionIndex() const
205 { return regionIdx_; }
206
214 void setRegionIndex(unsigned val)
215 { regionIdx_ = val; }
216
217 const Evaluation& maxOilSat() const
218 { return maxOilSat_; }
219
220 void setMaxOilSat(const Evaluation& val)
221 { maxOilSat_ = val; }
222
223 private:
224 Evaluation maxOilSat_;
225 unsigned regionIdx_;
226 };
227
228 /****************************************
229 * Initialization
230 ****************************************/
231#if HAVE_ECL_INPUT
235 static void initFromState(const EclipseState& eclState, const Schedule& schedule);
236#endif // HAVE_ECL_INPUT
237
246 static void initBegin(std::size_t numPvtRegions);
247
254 static void setEnableDissolvedGas(bool yesno)
255 { enableDissolvedGas_ = yesno; }
256
263 static void setEnableVaporizedOil(bool yesno)
264 { enableVaporizedOil_ = yesno; }
265
272 static void setEnableVaporizedWater(bool yesno)
273 { enableVaporizedWater_ = yesno; }
274
281 static void setEnableDissolvedGasInWater(bool yesno)
282 { enableDissolvedGasInWater_ = yesno; }
288 static void setEnableDiffusion(bool yesno)
289 { enableDiffusion_ = yesno; }
290
291
295 static void setGasPvt(std::shared_ptr<GasPvt> pvtObj)
296 { gasPvt_ = pvtObj; }
297
301 static void setOilPvt(std::shared_ptr<OilPvt> pvtObj)
302 { oilPvt_ = pvtObj; }
303
307 static void setWaterPvt(std::shared_ptr<WaterPvt> pvtObj)
308 { waterPvt_ = pvtObj; }
309
310 static void setVapPars(const Scalar par1, const Scalar par2)
311 {
312 if (gasPvt_) {
313 gasPvt_->setVapPars(par1, par2);
314 }
315 if (oilPvt_) {
316 oilPvt_->setVapPars(par1, par2);
317 }
318 if (waterPvt_) {
319 waterPvt_->setVapPars(par1, par2);
320 }
321 }
322
330 static void setReferenceDensities(Scalar rhoOil,
331 Scalar rhoWater,
332 Scalar rhoGas,
333 unsigned regionIdx);
334
338 static void initEnd();
339
340 static bool isInitialized()
341 { return isInitialized_; }
342
343 /****************************************
344 * Generic phase properties
345 ****************************************/
346
348 static constexpr unsigned numPhases = 3;
349
351 static constexpr unsigned waterPhaseIdx = IndexTraits::waterPhaseIdx;
353 static constexpr unsigned oilPhaseIdx = IndexTraits::oilPhaseIdx;
355 static constexpr unsigned gasPhaseIdx = IndexTraits::gasPhaseIdx;
356
359
362
364 static std::string_view phaseName(unsigned phaseIdx);
365
367 static bool isLiquid(unsigned phaseIdx)
368 {
369 assert(phaseIdx < numPhases);
370 return phaseIdx != gasPhaseIdx;
371 }
372
373 /****************************************
374 * Generic component related properties
375 ****************************************/
376
378 static constexpr unsigned numComponents = 3;
379
381 static constexpr unsigned oilCompIdx = IndexTraits::oilCompIdx;
383 static constexpr unsigned waterCompIdx = IndexTraits::waterCompIdx;
385 static constexpr unsigned gasCompIdx = IndexTraits::gasCompIdx;
386
387protected:
388 static unsigned char numActivePhases_;
389 static std::array<bool,numPhases> phaseIsActive_;
390
391public:
393 static unsigned numActivePhases()
394 { return numActivePhases_; }
395
397 static bool phaseIsActive(unsigned phaseIdx)
398 {
399 assert(phaseIdx < numPhases);
400 return phaseIsActive_[phaseIdx];
401 }
402
404 static unsigned solventComponentIndex(unsigned phaseIdx);
405
407 static unsigned soluteComponentIndex(unsigned phaseIdx);
408
410 static std::string_view componentName(unsigned compIdx);
411
413 static Scalar molarMass(unsigned compIdx, unsigned regionIdx = 0)
414 { return molarMass_[regionIdx][compIdx]; }
415
417 static bool isIdealMixture(unsigned /*phaseIdx*/)
418 {
419 // fugacity coefficients are only pressure dependent -> we
420 // have an ideal mixture
421 return true;
422 }
423
425 static bool isCompressible(unsigned /*phaseIdx*/)
426 { return true; /* all phases are compressible */ }
427
429 static bool isIdealGas(unsigned /*phaseIdx*/)
430 { return false; }
431
432
433 /****************************************
434 * Black-oil specific properties
435 ****************************************/
441 static std::size_t numRegions()
442 { return molarMass_.size(); }
443
450 static bool enableDissolvedGas()
451 { return enableDissolvedGas_; }
452
453
461 { return enableDissolvedGasInWater_; }
462
469 static bool enableVaporizedOil()
470 { return enableVaporizedOil_; }
471
479 { return enableVaporizedWater_; }
480
486 static bool enableDiffusion()
487 { return enableDiffusion_; }
488
494 static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
495 { return referenceDensity_[regionIdx][phaseIdx]; }
496
497 /****************************************
498 * thermodynamic quantities (generic version)
499 ****************************************/
501 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
502 static LhsEval density(const FluidState& fluidState,
503 const ParameterCache<ParamCacheEval>& paramCache,
504 unsigned phaseIdx)
505 { return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
506
508 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
509 static LhsEval fugacityCoefficient(const FluidState& fluidState,
510 const ParameterCache<ParamCacheEval>& paramCache,
511 unsigned phaseIdx,
512 unsigned compIdx)
513 {
514 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
515 phaseIdx,
516 compIdx,
517 paramCache.regionIndex());
518 }
519
521 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
522 static LhsEval viscosity(const FluidState& fluidState,
523 const ParameterCache<ParamCacheEval>& paramCache,
524 unsigned phaseIdx)
525 { return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
526
528 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
529 static LhsEval enthalpy(const FluidState& fluidState,
530 const ParameterCache<ParamCacheEval>& paramCache,
531 unsigned phaseIdx)
532 { return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
533
534 /****************************************
535 * thermodynamic quantities (black-oil specific version: Note that the PVT region
536 * index is explicitly passed instead of a parameter cache object)
537 ****************************************/
539 template <class FluidState, class LhsEval = typename FluidState::Scalar>
540 static LhsEval density(const FluidState& fluidState,
541 unsigned phaseIdx,
542 unsigned regionIdx)
543 {
544 assert(phaseIdx <= numPhases);
545 assert(regionIdx <= numRegions());
546
547 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
548 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
549 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
550
551 switch (phaseIdx) {
552 case oilPhaseIdx: {
553 if (enableDissolvedGas()) {
554 // miscible oil
555 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
556 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
557
558 return
559 bo*referenceDensity(oilPhaseIdx, regionIdx)
560 + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
561 }
562
563 // immiscible oil
564 const LhsEval Rs(0.0);
565 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
566
567 return referenceDensity(phaseIdx, regionIdx)*bo;
568 }
569
570 case gasPhaseIdx: {
572 // gas containing vaporized oil and vaporized water
573 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
574 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
575 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
576
577 return
578 bg*referenceDensity(gasPhaseIdx, regionIdx)
579 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
580 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
581 }
582 if (enableVaporizedOil()) {
583 // miscible gas
584 const LhsEval Rvw(0.0);
585 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
586 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
587
588 return
589 bg*referenceDensity(gasPhaseIdx, regionIdx)
590 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
591 }
592 if (enableVaporizedWater()) {
593 // gas containing vaporized water
594 const LhsEval Rv(0.0);
595 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
596 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
597
598 return
599 bg*referenceDensity(gasPhaseIdx, regionIdx)
600 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
601 }
602
603 // immiscible gas
604 const LhsEval Rv(0.0);
605 const LhsEval Rvw(0.0);
606 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
607 return bg*referenceDensity(phaseIdx, regionIdx);
608 }
609
610 case waterPhaseIdx:
612 // gas miscible in water
613 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
614 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
615 return
616 bw*referenceDensity(waterPhaseIdx, regionIdx)
617 + Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
618 }
619 const LhsEval Rsw(0.0);
620 return
622 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
623 }
624
625 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
626 }
627
637 template <class FluidState, class LhsEval = typename FluidState::Scalar>
638 static LhsEval saturatedDensity(const FluidState& fluidState,
639 unsigned phaseIdx,
640 unsigned regionIdx)
641 {
642 assert(phaseIdx <= numPhases);
643 assert(regionIdx <= numRegions());
644
645 const auto& p = fluidState.pressure(phaseIdx);
646 const auto& T = fluidState.temperature(phaseIdx);
647
648 switch (phaseIdx) {
649 case oilPhaseIdx: {
650 if (enableDissolvedGas()) {
651 // miscible oil
652 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
653 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
654
655 return
656 bo*referenceDensity(oilPhaseIdx, regionIdx)
657 + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
658 }
659
660 // immiscible oil
661 const LhsEval Rs(0.0);
662 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
663 return referenceDensity(phaseIdx, regionIdx)*bo;
664 }
665
666 case gasPhaseIdx: {
668 // gas containing vaporized oil and vaporized water
669 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
670 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
671 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
672
673 return
674 bg*referenceDensity(gasPhaseIdx, regionIdx)
675 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
676 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx) ;
677 }
678
679 if (enableVaporizedOil()) {
680 // miscible gas
681 const LhsEval Rvw(0.0);
682 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
683 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
684
685 return
686 bg*referenceDensity(gasPhaseIdx, regionIdx)
687 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
688 }
689
690 if (enableVaporizedWater()) {
691 // gas containing vaporized water
692 const LhsEval Rv(0.0);
693 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
694 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
695
696 return
697 bg*referenceDensity(gasPhaseIdx, regionIdx)
698 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
699 }
700
701 // immiscible gas
702 const LhsEval Rv(0.0);
703 const LhsEval Rvw(0.0);
704 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
705
706 return referenceDensity(phaseIdx, regionIdx)*bg;
707
708 }
709
710 case waterPhaseIdx:
711 {
713 // miscible in water
714 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
715 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
716 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
717 return
718 bw*referenceDensity(waterPhaseIdx, regionIdx)
719 + Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
720 }
721 return
723 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
724 }
725 }
726
727 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
728 }
729
738 template <class FluidState, class LhsEval = typename FluidState::Scalar>
739 static LhsEval inverseFormationVolumeFactor(const FluidState& fluidState,
740 unsigned phaseIdx,
741 unsigned regionIdx)
742 {
743 OPM_TIMEBLOCK_LOCAL(inverseFormationVolumeFactor);
744 assert(phaseIdx <= numPhases);
745 assert(regionIdx <= numRegions());
746
747 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
748 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
749
750 switch (phaseIdx) {
751 case oilPhaseIdx: {
752 if (enableDissolvedGas()) {
753 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
754 if (fluidState.saturation(gasPhaseIdx) > 0.0
755 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
756 {
757 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
758 } else {
759 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
760 }
761 }
762
763 const LhsEval Rs(0.0);
764 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
765 }
766 case gasPhaseIdx: {
768 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
769 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
770 if (fluidState.saturation(waterPhaseIdx) > 0.0
771 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
772 && fluidState.saturation(oilPhaseIdx) > 0.0
773 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
774 {
775 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
776 } else {
777 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
778 }
779 }
780
781 if (enableVaporizedOil()) {
782 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
783 if (fluidState.saturation(oilPhaseIdx) > 0.0
784 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
785 {
786 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
787 } else {
788 const LhsEval Rvw(0.0);
789 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
790 }
791 }
792
793 if (enableVaporizedWater()) {
794 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
795 if (fluidState.saturation(waterPhaseIdx) > 0.0
796 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
797 {
798 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
799 } else {
800 const LhsEval Rv(0.0);
801 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
802 }
803 }
804
805 const LhsEval Rv(0.0);
806 const LhsEval Rvw(0.0);
807 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
808 }
809 case waterPhaseIdx:
810 {
811 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
813 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
814 if (fluidState.saturation(gasPhaseIdx) > 0.0
815 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
816 {
817 return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
818 } else {
819 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
820 }
821 }
822 const LhsEval Rsw(0.0);
823 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
824 }
825 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
826 }
827 }
828
838 template <class FluidState, class LhsEval = typename FluidState::Scalar>
839 static LhsEval saturatedInverseFormationVolumeFactor(const FluidState& fluidState,
840 unsigned phaseIdx,
841 unsigned regionIdx)
842 {
843 OPM_TIMEBLOCK_LOCAL(saturatedInverseFormationVolumeFactor);
844 assert(phaseIdx <= numPhases);
845 assert(regionIdx <= numRegions());
846
847 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
848 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
849 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
850
851 switch (phaseIdx) {
852 case oilPhaseIdx: return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
853 case gasPhaseIdx: return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
854 case waterPhaseIdx: return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
855 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
856 }
857 }
858
860 template <class FluidState, class LhsEval = typename FluidState::Scalar>
861 static LhsEval fugacityCoefficient(const FluidState& fluidState,
862 unsigned phaseIdx,
863 unsigned compIdx,
864 unsigned regionIdx)
865 {
866 assert(phaseIdx <= numPhases);
867 assert(compIdx <= numComponents);
868 assert(regionIdx <= numRegions());
869
870 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
871 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
872
873 // for the fugacity coefficient of the oil component in the oil phase, we use
874 // some pseudo-realistic value for the vapor pressure to ease physical
875 // interpretation of the results
876 const LhsEval phi_oO = 20e3/p;
877
878 // for the gas component in the gas phase, assume it to be an ideal gas
879 constexpr const Scalar phi_gG = 1.0;
880
881 // for the fugacity coefficient of the water component in the water phase, we use
882 // the same approach as for the oil component in the oil phase
883 const LhsEval phi_wW = 30e3/p;
884
885 switch (phaseIdx) {
886 case gasPhaseIdx: // fugacity coefficients for all components in the gas phase
887 switch (compIdx) {
888 case gasCompIdx:
889 return phi_gG;
890
891 // for the oil component, we calculate the Rv value for saturated gas and Rs
892 // for saturated oil, and compute the fugacity coefficients at the
893 // equilibrium. for this, we assume that in equilibrium the fugacities of the
894 // oil component is the same in both phases.
895 case oilCompIdx: {
896 if (!enableVaporizedOil())
897 // if there's no vaporized oil, the gas phase is assumed to be
898 // immiscible with the oil component
899 return phi_gG*1e6;
900
901 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
902 const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
903 const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
904
905 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
906 const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
907 const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
908 const auto& x_oOSat = 1.0 - x_oGSat;
909
910 const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
911 const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
912
913 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
914 }
915
916 case waterCompIdx:
917 // the water component is assumed to be never miscible with the gas phase
918 return phi_gG*1e6;
919
920 default:
921 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
922 }
923
924 case oilPhaseIdx: // fugacity coefficients for all components in the oil phase
925 switch (compIdx) {
926 case oilCompIdx:
927 return phi_oO;
928
929 // for the oil and water components, we proceed analogous to the gas and
930 // water components in the gas phase
931 case gasCompIdx: {
932 if (!enableDissolvedGas())
933 // if there's no dissolved gas, the oil phase is assumed to be
934 // immiscible with the gas component
935 return phi_oO*1e6;
936
937 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
938 const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
939 const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
940 const auto& x_gGSat = 1.0 - x_gOSat;
941
942 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
943 const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
944 const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
945
946 const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
947 const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
948
949 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
950 }
951
952 case waterCompIdx:
953 return phi_oO*1e6;
954
955 default:
956 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
957 }
958
959 case waterPhaseIdx: // fugacity coefficients for all components in the water phase
960 // the water phase fugacity coefficients are pretty simple: because the water
961 // phase is assumed to consist entirely from the water component, we just
962 // need to make sure that the fugacity coefficients for the other components
963 // are a few orders of magnitude larger than the one of the water
964 // component. (i.e., the affinity of the gas and oil components to the water
965 // phase is lower by a few orders of magnitude)
966 switch (compIdx) {
967 case waterCompIdx: return phi_wW;
968 case oilCompIdx: return 1.1e6*phi_wW;
969 case gasCompIdx: return 1e6*phi_wW;
970 default:
971 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
972 }
973
974 default:
975 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
976 }
977
978 throw std::logic_error("Unhandled phase or component index");
979 }
980
982 template <class FluidState, class LhsEval = typename FluidState::Scalar>
983 static LhsEval viscosity(const FluidState& fluidState,
984 unsigned phaseIdx,
985 unsigned regionIdx)
986 {
987 OPM_TIMEBLOCK_LOCAL(viscosity);
988 assert(phaseIdx <= numPhases);
989 assert(regionIdx <= numRegions());
990
991 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
992 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
993
994 switch (phaseIdx) {
995 case oilPhaseIdx: {
996 if (enableDissolvedGas()) {
997 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
998 if (fluidState.saturation(gasPhaseIdx) > 0.0
999 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
1000 {
1001 return oilPvt_->saturatedViscosity(regionIdx, T, p);
1002 } else {
1003 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1004 }
1005 }
1006
1007 const LhsEval Rs(0.0);
1008 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1009 }
1010
1011 case gasPhaseIdx: {
1013 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1014 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1015 if (fluidState.saturation(waterPhaseIdx) > 0.0
1016 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1017 && fluidState.saturation(oilPhaseIdx) > 0.0
1018 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1019 {
1020 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1021 } else {
1022 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1023 }
1024 }
1025 if (enableVaporizedOil()) {
1026 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1027 if (fluidState.saturation(oilPhaseIdx) > 0.0
1028 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1029 {
1030 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1031 } else {
1032 const LhsEval Rvw(0.0);
1033 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1034 }
1035 }
1036 if (enableVaporizedWater()) {
1037 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1038 if (fluidState.saturation(waterPhaseIdx) > 0.0
1039 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1040 {
1041 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1042 } else {
1043 const LhsEval Rv(0.0);
1044 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1045 }
1046 }
1047
1048 const LhsEval Rv(0.0);
1049 const LhsEval Rvw(0.0);
1050 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1051 }
1052
1053 case waterPhaseIdx:
1054 {
1055 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1057 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1058 if (fluidState.saturation(gasPhaseIdx) > 0.0
1059 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
1060 {
1061 return waterPvt_->saturatedViscosity(regionIdx, T, p, saltConcentration);
1062 } else {
1063 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1064 }
1065 }
1066 const LhsEval Rsw(0.0);
1067 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1068 }
1069 }
1070
1071 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1072 }
1073
1075 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1076 static LhsEval enthalpy(const FluidState& fluidState,
1077 unsigned phaseIdx,
1078 unsigned regionIdx)
1079 {
1080 assert(phaseIdx <= numPhases);
1081 assert(regionIdx <= numRegions());
1082
1083 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1084 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1085
1086 switch (phaseIdx) {
1087 case oilPhaseIdx:
1088 return
1089 oilPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1090 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1091
1092 case gasPhaseIdx:
1093 return
1094 gasPvt_->internalEnergy(regionIdx, T, p,
1095 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1096 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1097 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1098
1099 case waterPhaseIdx:
1100 return
1101 waterPvt_->internalEnergy(regionIdx, T, p,
1102 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1103 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1104 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1105
1106 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1107 }
1108
1109 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1110 }
1111
1118 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1119 static LhsEval saturatedVaporizationFactor(const FluidState& fluidState,
1120 unsigned phaseIdx,
1121 unsigned regionIdx)
1122 {
1123 assert(phaseIdx <= numPhases);
1124 assert(regionIdx <= numRegions());
1125
1126 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1127 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1128 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
1129
1130 switch (phaseIdx) {
1131 case oilPhaseIdx: return 0.0;
1132 case gasPhaseIdx: return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p, saltConcentration);
1133 case waterPhaseIdx: return 0.0;
1134 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1135 }
1136 }
1137
1144 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1145 static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1146 unsigned phaseIdx,
1147 unsigned regionIdx,
1148 const LhsEval& maxOilSaturation)
1149 {
1150 OPM_TIMEBLOCK_LOCAL(saturatedDissolutionFactor);
1151 assert(phaseIdx <= numPhases);
1152 assert(regionIdx <= numRegions());
1153
1154 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1155 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1156 const auto& So = decay<LhsEval>(fluidState.saturation(oilPhaseIdx));
1157
1158 switch (phaseIdx) {
1159 case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1160 case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1161 case waterPhaseIdx: return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1162 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1163 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1164 }
1165 }
1166
1175 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1176 static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1177 unsigned phaseIdx,
1178 unsigned regionIdx)
1179 {
1180 OPM_TIMEBLOCK_LOCAL(saturatedDissolutionFactor);
1181 assert(phaseIdx <= numPhases);
1182 assert(regionIdx <= numRegions());
1183
1184 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1185 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1186
1187 switch (phaseIdx) {
1188 case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1189 case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1190 case waterPhaseIdx: return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1191 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1192 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1193 }
1194 }
1195
1199 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1200 static LhsEval bubblePointPressure(const FluidState& fluidState,
1201 unsigned regionIdx)
1202 {
1203 return saturationPressure(fluidState, oilPhaseIdx, regionIdx);
1204 }
1205
1206
1210 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1211 static LhsEval dewPointPressure(const FluidState& fluidState,
1212 unsigned regionIdx)
1213 {
1214 return saturationPressure(fluidState, gasPhaseIdx, regionIdx);
1215 }
1216
1227 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1228 static LhsEval saturationPressure(const FluidState& fluidState,
1229 unsigned phaseIdx,
1230 unsigned regionIdx)
1231 {
1232 assert(phaseIdx <= numPhases);
1233 assert(regionIdx <= numRegions());
1234
1235 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1236
1237 switch (phaseIdx) {
1238 case oilPhaseIdx: return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1239 case gasPhaseIdx: return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1240 case waterPhaseIdx: return waterPvt_->saturationPressure(regionIdx, T,
1241 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1242 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1243 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1244 }
1245 }
1246
1247 /****************************************
1248 * Auxiliary and convenience methods for the black-oil model
1249 ****************************************/
1254 template <class LhsEval>
1255 static LhsEval convertXoGToRs(const LhsEval& XoG, unsigned regionIdx)
1256 {
1257 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1258 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1259
1260 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1261 }
1262
1267 template <class LhsEval>
1268 static LhsEval convertXwGToRsw(const LhsEval& XwG, unsigned regionIdx)
1269 {
1270 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1271 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1272
1273 return XwG/(1.0 - XwG)*(rho_wRef/rho_gRef);
1274 }
1275
1280 template <class LhsEval>
1281 static LhsEval convertXgOToRv(const LhsEval& XgO, unsigned regionIdx)
1282 {
1283 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1284 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1285
1286 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1287 }
1288
1293 template <class LhsEval>
1294 static LhsEval convertXgWToRvw(const LhsEval& XgW, unsigned regionIdx)
1295 {
1296 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1297 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1298
1299 return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1300 }
1301
1302
1307 template <class LhsEval>
1308 static LhsEval convertRsToXoG(const LhsEval& Rs, unsigned regionIdx)
1309 {
1310 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1311 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1312
1313 const LhsEval& rho_oG = Rs*rho_gRef;
1314 return rho_oG/(rho_oRef + rho_oG);
1315 }
1316
1321 template <class LhsEval>
1322 static LhsEval convertRswToXwG(const LhsEval& Rsw, unsigned regionIdx)
1323 {
1324 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1325 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1326
1327 const LhsEval& rho_wG = Rsw*rho_gRef;
1328 return rho_wG/(rho_wRef + rho_wG);
1329 }
1330
1335 template <class LhsEval>
1336 static LhsEval convertRvToXgO(const LhsEval& Rv, unsigned regionIdx)
1337 {
1338 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1339 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1340
1341 const LhsEval& rho_gO = Rv*rho_oRef;
1342 return rho_gO/(rho_gRef + rho_gO);
1343 }
1344
1349 template <class LhsEval>
1350 static LhsEval convertRvwToXgW(const LhsEval& Rvw, unsigned regionIdx)
1351 {
1352 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1353 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1354
1355 const LhsEval& rho_gW = Rvw*rho_wRef;
1356 return rho_gW/(rho_gRef + rho_gW);
1357 }
1358
1362 template <class LhsEval>
1363 static LhsEval convertXgWToxgW(const LhsEval& XgW, unsigned regionIdx)
1364 {
1365 Scalar MW = molarMass_[regionIdx][waterCompIdx];
1366 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1367
1368 return XgW*MG / (MW*(1 - XgW) + XgW*MG);
1369 }
1370
1374 template <class LhsEval>
1375 static LhsEval convertXwGToxwG(const LhsEval& XwG, unsigned regionIdx)
1376 {
1377 Scalar MW = molarMass_[regionIdx][waterCompIdx];
1378 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1379
1380 return XwG*MW / (MG*(1 - XwG) + XwG*MW);
1381 }
1382
1386 template <class LhsEval>
1387 static LhsEval convertXoGToxoG(const LhsEval& XoG, unsigned regionIdx)
1388 {
1389 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1390 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1391
1392 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1393 }
1394
1398 template <class LhsEval>
1399 static LhsEval convertxoGToXoG(const LhsEval& xoG, unsigned regionIdx)
1400 {
1401 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1402 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1403
1404 return xoG*MG / (xoG*(MG - MO) + MO);
1405 }
1406
1410 template <class LhsEval>
1411 static LhsEval convertXgOToxgO(const LhsEval& XgO, unsigned regionIdx)
1412 {
1413 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1414 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1415
1416 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1417 }
1418
1422 template <class LhsEval>
1423 static LhsEval convertxgOToXgO(const LhsEval& xgO, unsigned regionIdx)
1424 {
1425 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1426 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1427
1428 return xgO*MO / (xgO*(MO - MG) + MG);
1429 }
1430
1438 static const GasPvt& gasPvt()
1439 { return *gasPvt_; }
1440
1448 static const OilPvt& oilPvt()
1449 { return *oilPvt_; }
1450
1458 static const WaterPvt& waterPvt()
1459 { return *waterPvt_; }
1460
1466 static Scalar reservoirTemperature(unsigned = 0)
1467 { return reservoirTemperature_; }
1468
1475 { reservoirTemperature_ = value; }
1476
1477 static short activeToCanonicalPhaseIdx(unsigned activePhaseIdx);
1478
1479 static short canonicalToActivePhaseIdx(unsigned phaseIdx);
1480
1482 static Scalar diffusionCoefficient(unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1483 { return diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx]; }
1484
1486 static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1487 { diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx] = coefficient ; }
1488
1492 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
1493 static LhsEval diffusionCoefficient(const FluidState& fluidState,
1494 const ParameterCache<ParamCacheEval>& paramCache,
1495 unsigned phaseIdx,
1496 unsigned compIdx)
1497 {
1498 // diffusion is disabled by the user
1499 if(!enableDiffusion())
1500 return 0.0;
1501
1502 // diffusion coefficients are set, and we use them
1503 if(!diffusionCoefficients_.empty()) {
1504 return diffusionCoefficient(compIdx, phaseIdx, paramCache.regionIndex());
1505 }
1506
1507 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1508 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1509
1510 switch (phaseIdx) {
1511 case oilPhaseIdx: return oilPvt().diffusionCoefficient(T, p, compIdx);
1512 case gasPhaseIdx: return gasPvt().diffusionCoefficient(T, p, compIdx);
1513 case waterPhaseIdx: return waterPvt().diffusionCoefficient(T, p, compIdx);
1514 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1515 }
1516 }
1517
1518private:
1519 static void resizeArrays_(std::size_t numRegions);
1520
1521 static Scalar reservoirTemperature_;
1522
1523 static std::shared_ptr<GasPvt> gasPvt_;
1524 static std::shared_ptr<OilPvt> oilPvt_;
1525 static std::shared_ptr<WaterPvt> waterPvt_;
1526
1527 static bool enableDissolvedGas_;
1528 static bool enableDissolvedGasInWater_;
1529 static bool enableVaporizedOil_;
1530 static bool enableVaporizedWater_;
1531 static bool enableDiffusion_;
1532
1533 // HACK for GCC 4.4: the array size has to be specified using the literal value '3'
1534 // here, because GCC 4.4 seems to be unable to determine the number of phases from
1535 // the BlackOil fluid system in the attribute declaration below...
1536 static std::vector<std::array<Scalar, /*numPhases=*/3> > referenceDensity_;
1537 static std::vector<std::array<Scalar, /*numComponents=*/3> > molarMass_;
1538 static std::vector<std::array<Scalar, /*numComponents=*/3 * /*numPhases=*/3> > diffusionCoefficients_;
1539
1540 static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1541 static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1542
1543 static bool isInitialized_;
1544};
1545
1546template <class Scalar, class IndexTraits>
1547unsigned char BlackOilFluidSystem<Scalar, IndexTraits>::numActivePhases_;
1548
1549template <class Scalar, class IndexTraits>
1550std::array<bool, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::phaseIsActive_;
1551
1552template <class Scalar, class IndexTraits>
1553std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::activeToCanonicalPhaseIdx_;
1554
1555template <class Scalar, class IndexTraits>
1556std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::canonicalToActivePhaseIdx_;
1557
1558template <class Scalar, class IndexTraits>
1559Scalar
1561
1562template <class Scalar, class IndexTraits>
1563Scalar
1565
1566template <class Scalar, class IndexTraits>
1567Scalar
1568BlackOilFluidSystem<Scalar, IndexTraits>::reservoirTemperature_;
1569
1570template <class Scalar, class IndexTraits>
1571bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGas_;
1572
1573template <class Scalar, class IndexTraits>
1574bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGasInWater_;
1575
1576
1577template <class Scalar, class IndexTraits>
1578bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedOil_;
1579
1580template <class Scalar, class IndexTraits>
1581bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedWater_;
1582
1583template <class Scalar, class IndexTraits>
1584bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDiffusion_;
1585
1586template <class Scalar, class IndexTraits>
1587std::shared_ptr<OilPvtMultiplexer<Scalar> >
1588BlackOilFluidSystem<Scalar, IndexTraits>::oilPvt_;
1589
1590template <class Scalar, class IndexTraits>
1591std::shared_ptr<GasPvtMultiplexer<Scalar> >
1592BlackOilFluidSystem<Scalar, IndexTraits>::gasPvt_;
1593
1594template <class Scalar, class IndexTraits>
1595std::shared_ptr<WaterPvtMultiplexer<Scalar> >
1596BlackOilFluidSystem<Scalar, IndexTraits>::waterPvt_;
1597
1598template <class Scalar, class IndexTraits>
1599std::vector<std::array<Scalar, 3> >
1600BlackOilFluidSystem<Scalar, IndexTraits>::referenceDensity_;
1601
1602template <class Scalar, class IndexTraits>
1603std::vector<std::array<Scalar, 3> >
1604BlackOilFluidSystem<Scalar, IndexTraits>::molarMass_;
1605
1606template <class Scalar, class IndexTraits>
1607std::vector<std::array<Scalar, 9> >
1608BlackOilFluidSystem<Scalar, IndexTraits>::diffusionCoefficients_;
1609
1610template <class Scalar, class IndexTraits>
1611bool BlackOilFluidSystem<Scalar, IndexTraits>::isInitialized_ = false;
1612
1613} // namespace Opm
1614
1615#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...
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a H2-Brine sy...
A central place for various physical constants occuring in some equations.
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:47
ScalarT Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:52
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
Definition BlackOilFluidSystem.hpp:162
static constexpr unsigned waterCompIdx
Index of the water component.
Definition BlackOilFluidSystem.hpp:383
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition BlackOilFluidSystem.hpp:393
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:1482
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:1255
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:1308
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition BlackOilFluidSystem.hpp:378
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:1350
static LhsEval convertXwGToxwG(const LhsEval &XwG, unsigned regionIdx)
Convert a gas mass fraction in the water phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1375
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition BlackOilFluidSystem.hpp:839
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1387
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition BlackOilFluidSystem.hpp:983
static void initEnd()
Finish initializing the black oil fluid system.
Definition BlackOilFluidSystem.cpp:260
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition BlackOilFluidSystem.hpp:1466
static LhsEval convertRswToXwG(const LhsEval &Rsw, unsigned regionIdx)
Convert a gas dissolution factor to the the corresponding mass fraction of the gas component in the w...
Definition BlackOilFluidSystem.hpp:1322
static bool phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition BlackOilFluidSystem.hpp:397
static void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem.hpp:288
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition BlackOilFluidSystem.hpp:1200
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition BlackOilFluidSystem.hpp:540
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition BlackOilFluidSystem.hpp:638
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:522
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition BlackOilFluidSystem.hpp:417
static bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition BlackOilFluidSystem.hpp:478
static constexpr unsigned numPhases
Number of fluid phases in the fluid system.
Definition BlackOilFluidSystem.hpp:348
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition BlackOilFluidSystem.hpp:307
static void initBegin(std::size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition BlackOilFluidSystem.cpp:223
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition BlackOilFluidSystem.hpp:1474
static void setEnableDissolvedGasInWater(bool yesno)
Specify whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem.hpp:281
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:1145
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:1076
static constexpr unsigned oilCompIdx
Index of the oil component.
Definition BlackOilFluidSystem.hpp:381
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition BlackOilFluidSystem.hpp:367
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:529
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:861
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition BlackOilFluidSystem.hpp:1399
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition BlackOilFluidSystem.hpp:413
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition BlackOilFluidSystem.hpp:1176
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition BlackOilFluidSystem.hpp:1211
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:1228
static LhsEval convertXwGToRsw(const LhsEval &XwG, unsigned regionIdx)
Convert the mass fraction of the gas component in the water phase to the corresponding gas dissolutio...
Definition BlackOilFluidSystem.hpp:1268
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition BlackOilFluidSystem.cpp:355
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1411
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:1493
static LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition BlackOilFluidSystem.hpp:1119
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition BlackOilFluidSystem.hpp:450
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition BlackOilFluidSystem.hpp:494
static Scalar surfaceTemperature
The temperature at the surface.
Definition BlackOilFluidSystem.hpp:361
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition BlackOilFluidSystem.hpp:301
static constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition BlackOilFluidSystem.hpp:351
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition BlackOilFluidSystem.hpp:1423
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition BlackOilFluidSystem.hpp:739
static LhsEval convertXgWToxgW(const LhsEval &XgW, unsigned regionIdx)
Convert a water mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1363
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition BlackOilFluidSystem.hpp:355
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:254
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:1281
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:509
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:502
static constexpr unsigned gasCompIdx
Index of the gas component.
Definition BlackOilFluidSystem.hpp:385
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:1336
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition BlackOilFluidSystem.hpp:1458
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition BlackOilFluidSystem.hpp:469
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition BlackOilFluidSystem.hpp:429
static unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition BlackOilFluidSystem.cpp:316
static Scalar surfacePressure
The pressure at the surface.
Definition BlackOilFluidSystem.hpp:358
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:272
static constexpr unsigned oilPhaseIdx
Index of the oil phase.
Definition BlackOilFluidSystem.hpp:353
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition BlackOilFluidSystem.hpp:1448
static bool enableDissolvedGasInWater()
Returns whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem.hpp:460
static unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition BlackOilFluidSystem.cpp:333
static bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem.hpp:486
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:1294
static void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition BlackOilFluidSystem.cpp:249
static std::size_t numRegions()
Returns the number of PVT regions which are considered.
Definition BlackOilFluidSystem.hpp:441
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:263
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition BlackOilFluidSystem.cpp:299
static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition BlackOilFluidSystem.hpp:1486
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition BlackOilFluidSystem.hpp:1438
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition BlackOilFluidSystem.hpp:425
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition BlackOilFluidSystem.hpp:295
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition GasPvtMultiplexer.hpp:109
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:336
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:104
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:269
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition WaterPvtMultiplexer.hpp:89
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 WaterPvtMultiplexer.hpp:261
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
The type of the fluid system's parameter cache.
Definition BlackOilFluidSystem.hpp:173
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:191
void setRegionIndex(unsigned val)
Set the index of the region which should be used to determine the thermodynamic properties.
Definition BlackOilFluidSystem.hpp:214
unsigned regionIndex() const
Return the index of the region which should be used to determine the thermodynamic properties.
Definition BlackOilFluidSystem.hpp:204