My Project
BlackOilFluidSystem.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_BLACK_OIL_FLUID_SYSTEM_HPP
28#define OPM_BLACK_OIL_FLUID_SYSTEM_HPP
29
35
36#include <opm/common/TimingMacros.hpp>
37
40
44
45#include <array>
46#include <cstddef>
47#include <memory>
48#include <stdexcept>
49#include <vector>
50
51namespace Opm {
52
53#if HAVE_ECL_INPUT
54class EclipseState;
55class Schedule;
56#endif
57
58namespace BlackOil {
59OPM_GENERATE_HAS_MEMBER(Rs, ) // Creates 'HasMember_Rs<T>'.
60OPM_GENERATE_HAS_MEMBER(Rv, ) // Creates 'HasMember_Rv<T>'.
61OPM_GENERATE_HAS_MEMBER(Rvw, ) // Creates 'HasMember_Rvw<T>'.
62OPM_GENERATE_HAS_MEMBER(Rsw, ) // Creates 'HasMember_Rsw<T>'.
63OPM_GENERATE_HAS_MEMBER(saltConcentration, )
64OPM_GENERATE_HAS_MEMBER(saltSaturation, )
65
66template <class FluidSystem, class FluidState, class LhsEval>
67LhsEval getRs_(typename std::enable_if<!HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
68 unsigned regionIdx)
69{
70 const auto& XoG =
71 decay<LhsEval>(fluidState.massFraction(FluidSystem::oilPhaseIdx, FluidSystem::gasCompIdx));
72 return FluidSystem::convertXoGToRs(XoG, regionIdx);
73}
74
75template <class FluidSystem, class FluidState, class LhsEval>
76auto getRs_(typename std::enable_if<HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
77 unsigned)
78 -> decltype(decay<LhsEval>(fluidState.Rs()))
79{ return decay<LhsEval>(fluidState.Rs()); }
80
81template <class FluidSystem, class FluidState, class LhsEval>
82LhsEval getRv_(typename std::enable_if<!HasMember_Rv<FluidState>::value, const FluidState&>::type fluidState,
83 unsigned regionIdx)
84{
85 const auto& XgO =
86 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::oilCompIdx));
87 return FluidSystem::convertXgOToRv(XgO, regionIdx);
88}
89
90template <class FluidSystem, class FluidState, class LhsEval>
91auto getRv_(typename std::enable_if<HasMember_Rv<FluidState>::value, const FluidState&>::type fluidState,
92 unsigned)
93 -> decltype(decay<LhsEval>(fluidState.Rv()))
94{ return decay<LhsEval>(fluidState.Rv()); }
95
96template <class FluidSystem, class FluidState, class LhsEval>
97LhsEval getRvw_(typename std::enable_if<!HasMember_Rvw<FluidState>::value, const FluidState&>::type fluidState,
98 unsigned regionIdx)
99{
100 const auto& XgW =
101 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::waterCompIdx));
102 return FluidSystem::convertXgWToRvw(XgW, regionIdx);
103}
104
105template <class FluidSystem, class FluidState, class LhsEval>
106auto getRvw_(typename std::enable_if<HasMember_Rvw<FluidState>::value, const FluidState&>::type fluidState,
107 unsigned)
108 -> decltype(decay<LhsEval>(fluidState.Rvw()))
109{ return decay<LhsEval>(fluidState.Rvw()); }
110
111template <class FluidSystem, class FluidState, class LhsEval>
112LhsEval getRsw_(typename std::enable_if<!HasMember_Rsw<FluidState>::value, const FluidState&>::type fluidState,
113 unsigned regionIdx)
114{
115 const auto& XwG =
116 decay<LhsEval>(fluidState.massFraction(FluidSystem::waterPhaseIdx, FluidSystem::gasCompIdx));
117 return FluidSystem::convertXwGToRsw(XwG, regionIdx);
118}
119
120template <class FluidSystem, class FluidState, class LhsEval>
121auto getRsw_(typename std::enable_if<HasMember_Rsw<FluidState>::value, const FluidState&>::type fluidState,
122 unsigned)
123 -> decltype(decay<LhsEval>(fluidState.Rsw()))
124{ return decay<LhsEval>(fluidState.Rsw()); }
125
126template <class FluidSystem, class FluidState, class LhsEval>
127LhsEval getSaltConcentration_(typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
128 const FluidState&>::type,
129 unsigned)
130{return 0.0;}
131
132template <class FluidSystem, class FluidState, class LhsEval>
133auto getSaltConcentration_(typename std::enable_if<HasMember_saltConcentration<FluidState>::value, const FluidState&>::type fluidState,
134 unsigned)
135 -> decltype(decay<LhsEval>(fluidState.saltConcentration()))
136{ return decay<LhsEval>(fluidState.saltConcentration()); }
137
138template <class FluidSystem, class FluidState, class LhsEval>
139LhsEval getSaltSaturation_(typename std::enable_if<!HasMember_saltSaturation<FluidState>::value,
140 const FluidState&>::type,
141 unsigned)
142{return 0.0;}
143
144template <class FluidSystem, class FluidState, class LhsEval>
145auto getSaltSaturation_(typename std::enable_if<HasMember_saltSaturation<FluidState>::value, const FluidState&>::type fluidState,
146 unsigned)
147 -> decltype(decay<LhsEval>(fluidState.saltSaturation()))
148{ return decay<LhsEval>(fluidState.saltSaturation()); }
149
150}
151
158template <class Scalar, class IndexTraits = BlackOilDefaultIndexTraits>
159class BlackOilFluidSystem : public BaseFluidSystem<Scalar, BlackOilFluidSystem<Scalar, IndexTraits> >
160{
162
163public:
167
169 template <class EvaluationT>
170 struct ParameterCache : public NullParameterCache<EvaluationT>
171 {
172 using Evaluation = EvaluationT;
173
174 public:
175 ParameterCache(Scalar maxOilSat = 1.0, unsigned regionIdx=0)
176 {
177 maxOilSat_ = maxOilSat;
178 regionIdx_ = regionIdx;
179 }
180
188 template <class OtherCache>
189 void assignPersistentData(const OtherCache& other)
190 {
191 regionIdx_ = other.regionIndex();
192 maxOilSat_ = other.maxOilSat();
193 }
194
202 unsigned regionIndex() const
203 { return regionIdx_; }
204
212 void setRegionIndex(unsigned val)
213 { regionIdx_ = val; }
214
215 const Evaluation& maxOilSat() const
216 { return maxOilSat_; }
217
218 void setMaxOilSat(const Evaluation& val)
219 { maxOilSat_ = val; }
220
221 private:
222 Evaluation maxOilSat_;
223 unsigned regionIdx_;
224 };
225
226 /****************************************
227 * Initialization
228 ****************************************/
229#if HAVE_ECL_INPUT
233 static void initFromState(const EclipseState& eclState, const Schedule& schedule);
234#endif // HAVE_ECL_INPUT
235
244 static void initBegin(std::size_t numPvtRegions);
245
252 static void setEnableDissolvedGas(bool yesno)
253 { enableDissolvedGas_ = yesno; }
254
261 static void setEnableVaporizedOil(bool yesno)
262 { enableVaporizedOil_ = yesno; }
263
270 static void setEnableVaporizedWater(bool yesno)
271 { enableVaporizedWater_ = yesno; }
272
279 static void setEnableDissolvedGasInWater(bool yesno)
280 { enableDissolvedGasInWater_ = yesno; }
286 static void setEnableDiffusion(bool yesno)
287 { enableDiffusion_ = yesno; }
288
289
293 static void setGasPvt(std::shared_ptr<GasPvt> pvtObj)
294 { gasPvt_ = pvtObj; }
295
299 static void setOilPvt(std::shared_ptr<OilPvt> pvtObj)
300 { oilPvt_ = pvtObj; }
301
305 static void setWaterPvt(std::shared_ptr<WaterPvt> pvtObj)
306 { waterPvt_ = pvtObj; }
307
315 static void setReferenceDensities(Scalar rhoOil,
316 Scalar rhoWater,
317 Scalar rhoGas,
318 unsigned regionIdx);
319
323 static void initEnd();
324
325 static bool isInitialized()
326 { return isInitialized_; }
327
328 /****************************************
329 * Generic phase properties
330 ****************************************/
331
333 static constexpr unsigned numPhases = 3;
334
336 static constexpr unsigned waterPhaseIdx = IndexTraits::waterPhaseIdx;
338 static constexpr unsigned oilPhaseIdx = IndexTraits::oilPhaseIdx;
340 static constexpr unsigned gasPhaseIdx = IndexTraits::gasPhaseIdx;
341
344
347
349 static const char* phaseName(unsigned phaseIdx);
350
352 static bool isLiquid(unsigned phaseIdx)
353 {
354 assert(phaseIdx < numPhases);
355 return phaseIdx != gasPhaseIdx;
356 }
357
358 /****************************************
359 * Generic component related properties
360 ****************************************/
361
363 static constexpr unsigned numComponents = 3;
364
366 static constexpr unsigned oilCompIdx = IndexTraits::oilCompIdx;
368 static constexpr unsigned waterCompIdx = IndexTraits::waterCompIdx;
370 static constexpr unsigned gasCompIdx = IndexTraits::gasCompIdx;
371
372protected:
373 static unsigned char numActivePhases_;
374 static std::array<bool,numPhases> phaseIsActive_;
375
376public:
378 static unsigned numActivePhases()
379 { return numActivePhases_; }
380
382 static unsigned phaseIsActive(unsigned phaseIdx)
383 {
384 assert(phaseIdx < numPhases);
385 return phaseIsActive_[phaseIdx];
386 }
387
389 static unsigned solventComponentIndex(unsigned phaseIdx);
390
392 static unsigned soluteComponentIndex(unsigned phaseIdx);
393
395 static const char* componentName(unsigned compIdx);
396
398 static Scalar molarMass(unsigned compIdx, unsigned regionIdx = 0)
399 { return molarMass_[regionIdx][compIdx]; }
400
402 static bool isIdealMixture(unsigned /*phaseIdx*/)
403 {
404 // fugacity coefficients are only pressure dependent -> we
405 // have an ideal mixture
406 return true;
407 }
408
410 static bool isCompressible(unsigned /*phaseIdx*/)
411 { return true; /* all phases are compressible */ }
412
414 static bool isIdealGas(unsigned /*phaseIdx*/)
415 { return false; }
416
417
418 /****************************************
419 * Black-oil specific properties
420 ****************************************/
426 static std::size_t numRegions()
427 { return molarMass_.size(); }
428
435 static bool enableDissolvedGas()
436 { return enableDissolvedGas_; }
437
438
446 { return enableDissolvedGasInWater_; }
447
454 static bool enableVaporizedOil()
455 { return enableVaporizedOil_; }
456
464 { return enableVaporizedWater_; }
465
471 static bool enableDiffusion()
472 { return enableDiffusion_; }
473
479 static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
480 { return referenceDensity_[regionIdx][phaseIdx]; }
481
482 /****************************************
483 * thermodynamic quantities (generic version)
484 ****************************************/
486 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
487 static LhsEval density(const FluidState& fluidState,
488 const ParameterCache<ParamCacheEval>& paramCache,
489 unsigned phaseIdx)
490 { return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
491
493 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
494 static LhsEval fugacityCoefficient(const FluidState& fluidState,
495 const ParameterCache<ParamCacheEval>& paramCache,
496 unsigned phaseIdx,
497 unsigned compIdx)
498 {
499 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
500 phaseIdx,
501 compIdx,
502 paramCache.regionIndex());
503 }
504
506 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
507 static LhsEval viscosity(const FluidState& fluidState,
508 const ParameterCache<ParamCacheEval>& paramCache,
509 unsigned phaseIdx)
510 { return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
511
513 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
514 static LhsEval enthalpy(const FluidState& fluidState,
515 const ParameterCache<ParamCacheEval>& paramCache,
516 unsigned phaseIdx)
517 { return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
518
519 /****************************************
520 * thermodynamic quantities (black-oil specific version: Note that the PVT region
521 * index is explicitly passed instead of a parameter cache object)
522 ****************************************/
524 template <class FluidState, class LhsEval = typename FluidState::Scalar>
525 static LhsEval density(const FluidState& fluidState,
526 unsigned phaseIdx,
527 unsigned regionIdx)
528 {
529 assert(phaseIdx <= numPhases);
530 assert(regionIdx <= numRegions());
531
532 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
533 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
534 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
535
536 switch (phaseIdx) {
537 case oilPhaseIdx: {
538 if (enableDissolvedGas()) {
539 // miscible oil
540 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
541 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
542
543 return
544 bo*referenceDensity(oilPhaseIdx, regionIdx)
545 + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
546 }
547
548 // immiscible oil
549 const LhsEval Rs(0.0);
550 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
551
552 return referenceDensity(phaseIdx, regionIdx)*bo;
553 }
554
555 case gasPhaseIdx: {
557 // gas containing vaporized oil and vaporized water
558 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
559 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
560 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
561
562 return
563 bg*referenceDensity(gasPhaseIdx, regionIdx)
564 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
565 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
566 }
567 if (enableVaporizedOil()) {
568 // miscible gas
569 const LhsEval Rvw(0.0);
570 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
571 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
572
573 return
574 bg*referenceDensity(gasPhaseIdx, regionIdx)
575 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
576 }
577 if (enableVaporizedWater()) {
578 // gas containing vaporized water
579 const LhsEval Rv(0.0);
580 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
581 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
582
583 return
584 bg*referenceDensity(gasPhaseIdx, regionIdx)
585 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
586 }
587
588 // immiscible gas
589 const LhsEval Rv(0.0);
590 const LhsEval Rvw(0.0);
591 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
592 return bg*referenceDensity(phaseIdx, regionIdx);
593 }
594
595 case waterPhaseIdx:
597 // gas miscible in water
598 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
599 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
600 return
601 bw*referenceDensity(waterPhaseIdx, regionIdx)
602 + Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
603 }
604 const LhsEval Rsw(0.0);
605 return
607 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
608 }
609
610 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
611 }
612
622 template <class FluidState, class LhsEval = typename FluidState::Scalar>
623 static LhsEval saturatedDensity(const FluidState& fluidState,
624 unsigned phaseIdx,
625 unsigned regionIdx)
626 {
627 assert(phaseIdx <= numPhases);
628 assert(regionIdx <= numRegions());
629
630 const auto& p = fluidState.pressure(phaseIdx);
631 const auto& T = fluidState.temperature(phaseIdx);
632
633 switch (phaseIdx) {
634 case oilPhaseIdx: {
635 if (enableDissolvedGas()) {
636 // miscible oil
637 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
638 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
639
640 return
641 bo*referenceDensity(oilPhaseIdx, regionIdx)
642 + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
643 }
644
645 // immiscible oil
646 const LhsEval Rs(0.0);
647 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
648 return referenceDensity(phaseIdx, regionIdx)*bo;
649 }
650
651 case gasPhaseIdx: {
653 // gas containing vaporized oil and vaporized water
654 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
655 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
656 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
657
658 return
659 bg*referenceDensity(gasPhaseIdx, regionIdx)
660 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
661 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx) ;
662 }
663
664 if (enableVaporizedOil()) {
665 // miscible gas
666 const LhsEval Rvw(0.0);
667 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
668 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
669
670 return
671 bg*referenceDensity(gasPhaseIdx, regionIdx)
672 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
673 }
674
675 if (enableVaporizedWater()) {
676 // gas containing vaporized water
677 const LhsEval Rv(0.0);
678 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
679 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
680
681 return
682 bg*referenceDensity(gasPhaseIdx, regionIdx)
683 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
684 }
685
686 // immiscible gas
687 const LhsEval Rv(0.0);
688 const LhsEval Rvw(0.0);
689 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
690
691 return referenceDensity(phaseIdx, regionIdx)*bg;
692
693 }
694
695 case waterPhaseIdx:
696 {
698 // miscible in water
699 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
700 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
701 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
702 return
703 bw*referenceDensity(waterPhaseIdx, regionIdx)
704 + Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
705 }
706 return
708 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
709 }
710 }
711
712 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
713 }
714
723 template <class FluidState, class LhsEval = typename FluidState::Scalar>
724 static LhsEval inverseFormationVolumeFactor(const FluidState& fluidState,
725 unsigned phaseIdx,
726 unsigned regionIdx)
727 {
728 OPM_TIMEBLOCK_LOCAL(inverseFormationVolumeFactor);
729 assert(phaseIdx <= numPhases);
730 assert(regionIdx <= numRegions());
731
732 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
733 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
734
735 switch (phaseIdx) {
736 case oilPhaseIdx: {
737 if (enableDissolvedGas()) {
738 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
739 if (fluidState.saturation(gasPhaseIdx) > 0.0
740 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
741 {
742 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
743 } else {
744 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
745 }
746 }
747
748 const LhsEval Rs(0.0);
749 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
750 }
751 case gasPhaseIdx: {
753 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
754 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
755 if (fluidState.saturation(waterPhaseIdx) > 0.0
756 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
757 && fluidState.saturation(oilPhaseIdx) > 0.0
758 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
759 {
760 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
761 } else {
762 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
763 }
764 }
765
766 if (enableVaporizedOil()) {
767 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
768 if (fluidState.saturation(oilPhaseIdx) > 0.0
769 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
770 {
771 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
772 } else {
773 const LhsEval Rvw(0.0);
774 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
775 }
776 }
777
778 if (enableVaporizedWater()) {
779 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
780 if (fluidState.saturation(waterPhaseIdx) > 0.0
781 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
782 {
783 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
784 } else {
785 const LhsEval Rv(0.0);
786 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
787 }
788 }
789
790 const LhsEval Rv(0.0);
791 const LhsEval Rvw(0.0);
792 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
793 }
794 case waterPhaseIdx:
795 {
796 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
798 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
799 if (fluidState.saturation(gasPhaseIdx) > 0.0
800 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
801 {
802 return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
803 } else {
804 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
805 }
806 }
807 const LhsEval Rsw(0.0);
808 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
809 }
810 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
811 }
812 }
813
823 template <class FluidState, class LhsEval = typename FluidState::Scalar>
824 static LhsEval saturatedInverseFormationVolumeFactor(const FluidState& fluidState,
825 unsigned phaseIdx,
826 unsigned regionIdx)
827 {
828 OPM_TIMEBLOCK_LOCAL(saturatedInverseFormationVolumeFactor);
829 assert(phaseIdx <= numPhases);
830 assert(regionIdx <= numRegions());
831
832 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
833 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
834 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
835
836 switch (phaseIdx) {
837 case oilPhaseIdx: return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
838 case gasPhaseIdx: return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
839 case waterPhaseIdx: return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
840 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
841 }
842 }
843
845 template <class FluidState, class LhsEval = typename FluidState::Scalar>
846 static LhsEval fugacityCoefficient(const FluidState& fluidState,
847 unsigned phaseIdx,
848 unsigned compIdx,
849 unsigned regionIdx)
850 {
851 assert(phaseIdx <= numPhases);
852 assert(compIdx <= numComponents);
853 assert(regionIdx <= numRegions());
854
855 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
856 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
857
858 // for the fugacity coefficient of the oil component in the oil phase, we use
859 // some pseudo-realistic value for the vapor pressure to ease physical
860 // interpretation of the results
861 const LhsEval phi_oO = 20e3/p;
862
863 // for the gas component in the gas phase, assume it to be an ideal gas
864 constexpr const Scalar phi_gG = 1.0;
865
866 // for the fugacity coefficient of the water component in the water phase, we use
867 // the same approach as for the oil component in the oil phase
868 const LhsEval phi_wW = 30e3/p;
869
870 switch (phaseIdx) {
871 case gasPhaseIdx: // fugacity coefficients for all components in the gas phase
872 switch (compIdx) {
873 case gasCompIdx:
874 return phi_gG;
875
876 // for the oil component, we calculate the Rv value for saturated gas and Rs
877 // for saturated oil, and compute the fugacity coefficients at the
878 // equilibrium. for this, we assume that in equilibrium the fugacities of the
879 // oil component is the same in both phases.
880 case oilCompIdx: {
881 if (!enableVaporizedOil())
882 // if there's no vaporized oil, the gas phase is assumed to be
883 // immiscible with the oil component
884 return phi_gG*1e6;
885
886 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
887 const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
888 const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
889
890 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
891 const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
892 const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
893 const auto& x_oOSat = 1.0 - x_oGSat;
894
895 const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
896 const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
897
898 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
899 }
900
901 case waterCompIdx:
902 // the water component is assumed to be never miscible with the gas phase
903 return phi_gG*1e6;
904
905 default:
906 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
907 }
908
909 case oilPhaseIdx: // fugacity coefficients for all components in the oil phase
910 switch (compIdx) {
911 case oilCompIdx:
912 return phi_oO;
913
914 // for the oil and water components, we proceed analogous to the gas and
915 // water components in the gas phase
916 case gasCompIdx: {
917 if (!enableDissolvedGas())
918 // if there's no dissolved gas, the oil phase is assumed to be
919 // immiscible with the gas component
920 return phi_oO*1e6;
921
922 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
923 const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
924 const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
925 const auto& x_gGSat = 1.0 - x_gOSat;
926
927 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
928 const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
929 const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
930
931 const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
932 const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
933
934 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
935 }
936
937 case waterCompIdx:
938 return phi_oO*1e6;
939
940 default:
941 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
942 }
943
944 case waterPhaseIdx: // fugacity coefficients for all components in the water phase
945 // the water phase fugacity coefficients are pretty simple: because the water
946 // phase is assumed to consist entirely from the water component, we just
947 // need to make sure that the fugacity coefficients for the other components
948 // are a few orders of magnitude larger than the one of the water
949 // component. (i.e., the affinity of the gas and oil components to the water
950 // phase is lower by a few orders of magnitude)
951 switch (compIdx) {
952 case waterCompIdx: return phi_wW;
953 case oilCompIdx: return 1.1e6*phi_wW;
954 case gasCompIdx: return 1e6*phi_wW;
955 default:
956 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
957 }
958
959 default:
960 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
961 }
962
963 throw std::logic_error("Unhandled phase or component index");
964 }
965
967 template <class FluidState, class LhsEval = typename FluidState::Scalar>
968 static LhsEval viscosity(const FluidState& fluidState,
969 unsigned phaseIdx,
970 unsigned regionIdx)
971 {
972 OPM_TIMEBLOCK_LOCAL(viscosity);
973 assert(phaseIdx <= numPhases);
974 assert(regionIdx <= numRegions());
975
976 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
977 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
978
979 switch (phaseIdx) {
980 case oilPhaseIdx: {
981 if (enableDissolvedGas()) {
982 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
983 if (fluidState.saturation(gasPhaseIdx) > 0.0
984 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
985 {
986 return oilPvt_->saturatedViscosity(regionIdx, T, p);
987 } else {
988 return oilPvt_->viscosity(regionIdx, T, p, Rs);
989 }
990 }
991
992 const LhsEval Rs(0.0);
993 return oilPvt_->viscosity(regionIdx, T, p, Rs);
994 }
995
996 case gasPhaseIdx: {
998 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
999 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1000 if (fluidState.saturation(waterPhaseIdx) > 0.0
1001 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1002 && fluidState.saturation(oilPhaseIdx) > 0.0
1003 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1004 {
1005 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1006 } else {
1007 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1008 }
1009 }
1010 if (enableVaporizedOil()) {
1011 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1012 if (fluidState.saturation(oilPhaseIdx) > 0.0
1013 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1014 {
1015 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1016 } else {
1017 const LhsEval Rvw(0.0);
1018 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1019 }
1020 }
1021 if (enableVaporizedWater()) {
1022 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1023 if (fluidState.saturation(waterPhaseIdx) > 0.0
1024 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1025 {
1026 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1027 } else {
1028 const LhsEval Rv(0.0);
1029 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1030 }
1031 }
1032
1033 const LhsEval Rv(0.0);
1034 const LhsEval Rvw(0.0);
1035 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1036 }
1037
1038 case waterPhaseIdx:
1039 {
1040 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1042 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1043 if (fluidState.saturation(gasPhaseIdx) > 0.0
1044 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
1045 {
1046 return waterPvt_->saturatedViscosity(regionIdx, T, p, saltConcentration);
1047 } else {
1048 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1049 }
1050 }
1051 const LhsEval Rsw(0.0);
1052 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1053 }
1054 }
1055
1056 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1057 }
1058
1060 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1061 static LhsEval enthalpy(const FluidState& fluidState,
1062 unsigned phaseIdx,
1063 unsigned regionIdx)
1064 {
1065 assert(phaseIdx <= numPhases);
1066 assert(regionIdx <= numRegions());
1067
1068 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1069 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1070
1071 switch (phaseIdx) {
1072 case oilPhaseIdx:
1073 return
1074 oilPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1075 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1076
1077 case gasPhaseIdx:
1078 return
1079 gasPvt_->internalEnergy(regionIdx, T, p,
1080 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1081 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1082 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1083
1084 case waterPhaseIdx:
1085 return
1086 waterPvt_->internalEnergy(regionIdx, T, p,
1087 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1088 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1089 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1090
1091 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1092 }
1093
1094 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1095 }
1096
1103 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1104 static LhsEval saturatedVaporizationFactor(const FluidState& fluidState,
1105 unsigned phaseIdx,
1106 unsigned regionIdx)
1107 {
1108 assert(phaseIdx <= numPhases);
1109 assert(regionIdx <= numRegions());
1110
1111 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1112 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1113 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
1114
1115 switch (phaseIdx) {
1116 case oilPhaseIdx: return 0.0;
1117 case gasPhaseIdx: return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p, saltConcentration);
1118 case waterPhaseIdx: return 0.0;
1119 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1120 }
1121 }
1122
1129 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1130 static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1131 unsigned phaseIdx,
1132 unsigned regionIdx,
1133 const LhsEval& maxOilSaturation)
1134 {
1135 OPM_TIMEBLOCK_LOCAL(saturatedDissolutionFactor);
1136 assert(phaseIdx <= numPhases);
1137 assert(regionIdx <= numRegions());
1138
1139 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1140 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1141 const auto& So = decay<LhsEval>(fluidState.saturation(oilPhaseIdx));
1142
1143 switch (phaseIdx) {
1144 case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1145 case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1146 case waterPhaseIdx: return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1147 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1148 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1149 }
1150 }
1151
1160 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1161 static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1162 unsigned phaseIdx,
1163 unsigned regionIdx)
1164 {
1165 OPM_TIMEBLOCK_LOCAL(saturatedDissolutionFactor);
1166 assert(phaseIdx <= numPhases);
1167 assert(regionIdx <= numRegions());
1168
1169 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1170 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1171
1172 switch (phaseIdx) {
1173 case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1174 case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1175 case waterPhaseIdx: return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1176 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1177 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1178 }
1179 }
1180
1184 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1185 static LhsEval bubblePointPressure(const FluidState& fluidState,
1186 unsigned regionIdx)
1187 {
1188 return saturationPressure(fluidState, oilPhaseIdx, regionIdx);
1189 }
1190
1191
1195 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1196 static LhsEval dewPointPressure(const FluidState& fluidState,
1197 unsigned regionIdx)
1198 {
1199 return saturationPressure(fluidState, gasPhaseIdx, regionIdx);
1200 }
1201
1212 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1213 static LhsEval saturationPressure(const FluidState& fluidState,
1214 unsigned phaseIdx,
1215 unsigned regionIdx)
1216 {
1217 assert(phaseIdx <= numPhases);
1218 assert(regionIdx <= numRegions());
1219
1220 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1221
1222 switch (phaseIdx) {
1223 case oilPhaseIdx: return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1224 case gasPhaseIdx: return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1225 case waterPhaseIdx: return waterPvt_->saturationPressure(regionIdx, T,
1226 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1227 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1228 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1229 }
1230 }
1231
1232 /****************************************
1233 * Auxiliary and convenience methods for the black-oil model
1234 ****************************************/
1239 template <class LhsEval>
1240 static LhsEval convertXoGToRs(const LhsEval& XoG, unsigned regionIdx)
1241 {
1242 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1243 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1244
1245 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1246 }
1247
1252 template <class LhsEval>
1253 static LhsEval convertXwGToRsw(const LhsEval& XwG, unsigned regionIdx)
1254 {
1255 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1256 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1257
1258 return XwG/(1.0 - XwG)*(rho_wRef/rho_gRef);
1259 }
1260
1265 template <class LhsEval>
1266 static LhsEval convertXgOToRv(const LhsEval& XgO, unsigned regionIdx)
1267 {
1268 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1269 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1270
1271 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1272 }
1273
1278 template <class LhsEval>
1279 static LhsEval convertXgWToRvw(const LhsEval& XgW, unsigned regionIdx)
1280 {
1281 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1282 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1283
1284 return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1285 }
1286
1287
1292 template <class LhsEval>
1293 static LhsEval convertRsToXoG(const LhsEval& Rs, unsigned regionIdx)
1294 {
1295 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1296 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1297
1298 const LhsEval& rho_oG = Rs*rho_gRef;
1299 return rho_oG/(rho_oRef + rho_oG);
1300 }
1301
1306 template <class LhsEval>
1307 static LhsEval convertRswToXwG(const LhsEval& Rsw, unsigned regionIdx)
1308 {
1309 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1310 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1311
1312 const LhsEval& rho_wG = Rsw*rho_gRef;
1313 return rho_wG/(rho_wRef + rho_wG);
1314 }
1315
1320 template <class LhsEval>
1321 static LhsEval convertRvToXgO(const LhsEval& Rv, unsigned regionIdx)
1322 {
1323 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1324 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1325
1326 const LhsEval& rho_gO = Rv*rho_oRef;
1327 return rho_gO/(rho_gRef + rho_gO);
1328 }
1329
1334 template <class LhsEval>
1335 static LhsEval convertRvwToXgW(const LhsEval& Rvw, unsigned regionIdx)
1336 {
1337 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1338 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1339
1340 const LhsEval& rho_gW = Rvw*rho_wRef;
1341 return rho_gW/(rho_gRef + rho_gW);
1342 }
1343
1347 template <class LhsEval>
1348 static LhsEval convertXgWToxgW(const LhsEval& XgW, unsigned regionIdx)
1349 {
1350 Scalar MW = molarMass_[regionIdx][waterCompIdx];
1351 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1352
1353 return XgW*MG / (MW*(1 - XgW) + XgW*MG);
1354 }
1355
1359 template <class LhsEval>
1360 static LhsEval convertXwGToxwG(const LhsEval& XwG, unsigned regionIdx)
1361 {
1362 Scalar MW = molarMass_[regionIdx][waterCompIdx];
1363 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1364
1365 return XwG*MW / (MG*(1 - XwG) + XwG*MW);
1366 }
1367
1371 template <class LhsEval>
1372 static LhsEval convertXoGToxoG(const LhsEval& XoG, unsigned regionIdx)
1373 {
1374 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1375 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1376
1377 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1378 }
1379
1383 template <class LhsEval>
1384 static LhsEval convertxoGToXoG(const LhsEval& xoG, unsigned regionIdx)
1385 {
1386 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1387 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1388
1389 return xoG*MG / (xoG*(MG - MO) + MO);
1390 }
1391
1395 template <class LhsEval>
1396 static LhsEval convertXgOToxgO(const LhsEval& XgO, unsigned regionIdx)
1397 {
1398 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1399 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1400
1401 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1402 }
1403
1407 template <class LhsEval>
1408 static LhsEval convertxgOToXgO(const LhsEval& xgO, unsigned regionIdx)
1409 {
1410 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1411 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1412
1413 return xgO*MO / (xgO*(MO - MG) + MG);
1414 }
1415
1423 static const GasPvt& gasPvt()
1424 { return *gasPvt_; }
1425
1433 static const OilPvt& oilPvt()
1434 { return *oilPvt_; }
1435
1443 static const WaterPvt& waterPvt()
1444 { return *waterPvt_; }
1445
1451 static Scalar reservoirTemperature(unsigned = 0)
1452 { return reservoirTemperature_; }
1453
1460 { reservoirTemperature_ = value; }
1461
1462 static short activeToCanonicalPhaseIdx(unsigned activePhaseIdx);
1463
1464 static short canonicalToActivePhaseIdx(unsigned phaseIdx);
1465
1467 static Scalar diffusionCoefficient(unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1468 { return diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx]; }
1469
1471 static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1472 { diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx] = coefficient ; }
1473
1477 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
1478 static LhsEval diffusionCoefficient(const FluidState& fluidState,
1479 const ParameterCache<ParamCacheEval>& paramCache,
1480 unsigned phaseIdx,
1481 unsigned compIdx)
1482 {
1483 // diffusion is disabled by the user
1484 if(!enableDiffusion())
1485 return 0.0;
1486
1487 // diffusion coefficients are set, and we use them
1488 if(!diffusionCoefficients_.empty()) {
1489 return diffusionCoefficient(compIdx, phaseIdx, paramCache.regionIndex());
1490 }
1491
1492 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1493 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1494
1495 switch (phaseIdx) {
1496 case oilPhaseIdx: return oilPvt().diffusionCoefficient(T, p, compIdx);
1497 case gasPhaseIdx: return gasPvt().diffusionCoefficient(T, p, compIdx);
1498 case waterPhaseIdx: return waterPvt().diffusionCoefficient(T, p, compIdx);
1499 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1500 }
1501 }
1502
1503private:
1504 static void resizeArrays_(std::size_t numRegions);
1505
1506 static Scalar reservoirTemperature_;
1507
1508 static std::shared_ptr<GasPvt> gasPvt_;
1509 static std::shared_ptr<OilPvt> oilPvt_;
1510 static std::shared_ptr<WaterPvt> waterPvt_;
1511
1512 static bool enableDissolvedGas_;
1513 static bool enableDissolvedGasInWater_;
1514 static bool enableVaporizedOil_;
1515 static bool enableVaporizedWater_;
1516 static bool enableDiffusion_;
1517
1518 // HACK for GCC 4.4: the array size has to be specified using the literal value '3'
1519 // here, because GCC 4.4 seems to be unable to determine the number of phases from
1520 // the BlackOil fluid system in the attribute declaration below...
1521 static std::vector<std::array<Scalar, /*numPhases=*/3> > referenceDensity_;
1522 static std::vector<std::array<Scalar, /*numComponents=*/3> > molarMass_;
1523 static std::vector<std::array<Scalar, /*numComponents=*/3 * /*numPhases=*/3> > diffusionCoefficients_;
1524
1525 static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1526 static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1527
1528 static bool isInitialized_;
1529};
1530
1531template <class Scalar, class IndexTraits>
1532unsigned char BlackOilFluidSystem<Scalar, IndexTraits>::numActivePhases_;
1533
1534template <class Scalar, class IndexTraits>
1535std::array<bool, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::phaseIsActive_;
1536
1537template <class Scalar, class IndexTraits>
1538std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::activeToCanonicalPhaseIdx_;
1539
1540template <class Scalar, class IndexTraits>
1541std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::canonicalToActivePhaseIdx_;
1542
1543template <class Scalar, class IndexTraits>
1544Scalar
1546
1547template <class Scalar, class IndexTraits>
1548Scalar
1550
1551template <class Scalar, class IndexTraits>
1552Scalar
1553BlackOilFluidSystem<Scalar, IndexTraits>::reservoirTemperature_;
1554
1555template <class Scalar, class IndexTraits>
1556bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGas_;
1557
1558template <class Scalar, class IndexTraits>
1559bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGasInWater_;
1560
1561
1562template <class Scalar, class IndexTraits>
1563bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedOil_;
1564
1565template <class Scalar, class IndexTraits>
1566bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedWater_;
1567
1568template <class Scalar, class IndexTraits>
1569bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDiffusion_;
1570
1571template <class Scalar, class IndexTraits>
1572std::shared_ptr<OilPvtMultiplexer<Scalar> >
1573BlackOilFluidSystem<Scalar, IndexTraits>::oilPvt_;
1574
1575template <class Scalar, class IndexTraits>
1576std::shared_ptr<GasPvtMultiplexer<Scalar> >
1577BlackOilFluidSystem<Scalar, IndexTraits>::gasPvt_;
1578
1579template <class Scalar, class IndexTraits>
1580std::shared_ptr<WaterPvtMultiplexer<Scalar> >
1581BlackOilFluidSystem<Scalar, IndexTraits>::waterPvt_;
1582
1583template <class Scalar, class IndexTraits>
1584std::vector<std::array<Scalar, 3> >
1585BlackOilFluidSystem<Scalar, IndexTraits>::referenceDensity_;
1586
1587template <class Scalar, class IndexTraits>
1588std::vector<std::array<Scalar, 3> >
1589BlackOilFluidSystem<Scalar, IndexTraits>::molarMass_;
1590
1591template <class Scalar, class IndexTraits>
1592std::vector<std::array<Scalar, 9> >
1593BlackOilFluidSystem<Scalar, IndexTraits>::diffusionCoefficients_;
1594
1595template <class Scalar, class IndexTraits>
1596bool BlackOilFluidSystem<Scalar, IndexTraits>::isInitialized_ = false;
1597
1598} // namespace Opm
1599
1600#endif
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...
Some templates to wrap the valgrind client request macros.
The base class for all fluid systems.
Definition: BaseFluidSystem.hpp:46
Scalar 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:160
static constexpr unsigned waterCompIdx
Index of the water component.
Definition: BlackOilFluidSystem.hpp:368
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition: BlackOilFluidSystem.hpp:378
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:1467
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:1240
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:1293
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition: BlackOilFluidSystem.hpp:363
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:1335
static LhsEval convertXwGToxwG(const LhsEval &XwG, unsigned regionIdx)
Convert a gas mass fraction in the water phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1360
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:824
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:1372
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: BlackOilFluidSystem.hpp:968
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1451
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:1307
static void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:286
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition: BlackOilFluidSystem.hpp:1185
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: BlackOilFluidSystem.hpp:525
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:623
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:507
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: BlackOilFluidSystem.hpp:402
static bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition: BlackOilFluidSystem.hpp:463
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:333
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition: BlackOilFluidSystem.hpp:305
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1459
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:279
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:1130
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:1061
static constexpr unsigned oilCompIdx
Index of the oil component.
Definition: BlackOilFluidSystem.hpp:366
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:352
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:514
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:846
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1384
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition: BlackOilFluidSystem.hpp:398
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:1161
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition: BlackOilFluidSystem.hpp:1196
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:1213
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:1253
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1396
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:1478
static LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition: BlackOilFluidSystem.hpp:1104
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition: BlackOilFluidSystem.hpp:435
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition: BlackOilFluidSystem.hpp:479
static Scalar surfaceTemperature
The temperature at the surface.
Definition: BlackOilFluidSystem.hpp:346
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition: BlackOilFluidSystem.hpp:299
static constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition: BlackOilFluidSystem.hpp:336
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1408
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:724
static LhsEval convertXgWToxgW(const LhsEval &XgW, unsigned regionIdx)
Convert a water mass fraction in the gas phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1348
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition: BlackOilFluidSystem.hpp:340
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:252
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:1266
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:494
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:487
static constexpr unsigned gasCompIdx
Index of the gas component.
Definition: BlackOilFluidSystem.hpp:370
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:1321
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition: BlackOilFluidSystem.hpp:1443
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition: BlackOilFluidSystem.hpp:454
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: BlackOilFluidSystem.hpp:414
static unsigned phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition: BlackOilFluidSystem.hpp:382
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:343
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:270
static constexpr unsigned oilPhaseIdx
Index of the oil phase.
Definition: BlackOilFluidSystem.hpp:338
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition: BlackOilFluidSystem.hpp:1433
static bool enableDissolvedGasInWater()
Returns whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition: BlackOilFluidSystem.hpp:445
static bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:471
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:1279
static std::size_t numRegions()
Returns the number of PVT regions which are considered.
Definition: BlackOilFluidSystem.hpp:426
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:261
static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition: BlackOilFluidSystem.hpp:1471
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition: BlackOilFluidSystem.hpp:1423
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: BlackOilFluidSystem.hpp:410
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition: BlackOilFluidSystem.hpp:293
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition: GasPvtMultiplexer.hpp:102
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:316
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:97
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:254
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition: WaterPvtMultiplexer.hpp:82
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:245
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:171
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:189
void setRegionIndex(unsigned val)
Set the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:212
unsigned regionIndex() const
Return the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:202