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