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
34
35#include <opm/common/TimingMacros.hpp>
36
39
44
45#include <array>
46#include <cstddef>
47#include <memory>
48#include <stdexcept>
49#include <string_view>
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 explicit ParameterCache(Scalar maxOilSat = 1.0, unsigned regionIdx = 0)
177 : maxOilSat_(maxOilSat)
178 , regionIdx_(regionIdx)
179 {
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
295 static void setUseSaturatedTables(bool yesno)
296 { useSaturatedTables_ = yesno; }
297
301 static void setGasPvt(std::shared_ptr<GasPvt> pvtObj)
302 { gasPvt_ = pvtObj; }
303
307 static void setOilPvt(std::shared_ptr<OilPvt> pvtObj)
308 { oilPvt_ = pvtObj; }
309
313 static void setWaterPvt(std::shared_ptr<WaterPvt> pvtObj)
314 { waterPvt_ = pvtObj; }
315
316 static void setVapPars(const Scalar par1, const Scalar par2)
317 {
318 if (gasPvt_) {
319 gasPvt_->setVapPars(par1, par2);
320 }
321 if (oilPvt_) {
322 oilPvt_->setVapPars(par1, par2);
323 }
324 if (waterPvt_) {
325 waterPvt_->setVapPars(par1, par2);
326 }
327 }
328
336 static void setReferenceDensities(Scalar rhoOil,
337 Scalar rhoWater,
338 Scalar rhoGas,
339 unsigned regionIdx);
340
344 static void initEnd();
345
346 static bool isInitialized()
347 { return isInitialized_; }
348
349 /****************************************
350 * Generic phase properties
351 ****************************************/
352
354 static constexpr unsigned numPhases = 3;
355
357 static constexpr unsigned waterPhaseIdx = IndexTraits::waterPhaseIdx;
359 static constexpr unsigned oilPhaseIdx = IndexTraits::oilPhaseIdx;
361 static constexpr unsigned gasPhaseIdx = IndexTraits::gasPhaseIdx;
362
365
368
370 static std::string_view phaseName(unsigned phaseIdx);
371
373 static bool isLiquid(unsigned phaseIdx)
374 {
375 assert(phaseIdx < numPhases);
376 return phaseIdx != gasPhaseIdx;
377 }
378
379 /****************************************
380 * Generic component related properties
381 ****************************************/
382
384 static constexpr unsigned numComponents = 3;
385
387 static constexpr unsigned oilCompIdx = IndexTraits::oilCompIdx;
389 static constexpr unsigned waterCompIdx = IndexTraits::waterCompIdx;
391 static constexpr unsigned gasCompIdx = IndexTraits::gasCompIdx;
392
393protected:
394 static unsigned char numActivePhases_;
395 static std::array<bool,numPhases> phaseIsActive_;
396
397public:
399 static unsigned numActivePhases()
400 { return numActivePhases_; }
401
403 static bool phaseIsActive(unsigned phaseIdx)
404 {
405 assert(phaseIdx < numPhases);
406 return phaseIsActive_[phaseIdx];
407 }
408
410 static unsigned solventComponentIndex(unsigned phaseIdx);
411
413 static unsigned soluteComponentIndex(unsigned phaseIdx);
414
416 static std::string_view componentName(unsigned compIdx);
417
419 static Scalar molarMass(unsigned compIdx, unsigned regionIdx = 0)
420 { return molarMass_[regionIdx][compIdx]; }
421
423 static bool isIdealMixture(unsigned /*phaseIdx*/)
424 {
425 // fugacity coefficients are only pressure dependent -> we
426 // have an ideal mixture
427 return true;
428 }
429
431 static bool isCompressible(unsigned /*phaseIdx*/)
432 { return true; /* all phases are compressible */ }
433
435 static bool isIdealGas(unsigned /*phaseIdx*/)
436 { return false; }
437
438
439 /****************************************
440 * Black-oil specific properties
441 ****************************************/
447 static std::size_t numRegions()
448 { return molarMass_.size(); }
449
456 static bool enableDissolvedGas()
457 { return enableDissolvedGas_; }
458
459
467 { return enableDissolvedGasInWater_; }
468
475 static bool enableVaporizedOil()
476 { return enableVaporizedOil_; }
477
485 { return enableVaporizedWater_; }
486
492 static bool enableDiffusion()
493 { return enableDiffusion_; }
494
500 static bool useSaturatedTables()
501 { return useSaturatedTables_; }
502
508 static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
509 { return referenceDensity_[regionIdx][phaseIdx]; }
510
511 /****************************************
512 * thermodynamic quantities (generic version)
513 ****************************************/
515 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
516 static LhsEval density(const FluidState& fluidState,
517 const ParameterCache<ParamCacheEval>& paramCache,
518 unsigned phaseIdx)
519 { return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
520
522 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
523 static LhsEval fugacityCoefficient(const FluidState& fluidState,
524 const ParameterCache<ParamCacheEval>& paramCache,
525 unsigned phaseIdx,
526 unsigned compIdx)
527 {
528 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
529 phaseIdx,
530 compIdx,
531 paramCache.regionIndex());
532 }
533
535 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
536 static LhsEval viscosity(const FluidState& fluidState,
537 const ParameterCache<ParamCacheEval>& paramCache,
538 unsigned phaseIdx)
539 { return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
540
542 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
543 static LhsEval enthalpy(const FluidState& fluidState,
544 const ParameterCache<ParamCacheEval>& paramCache,
545 unsigned phaseIdx)
546 { return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
547
548 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
549 static LhsEval internalEnergy(const FluidState& fluidState,
550 const ParameterCache<ParamCacheEval>& paramCache,
551 unsigned phaseIdx)
552 { return internalEnergy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
553
554 /****************************************
555 * thermodynamic quantities (black-oil specific version: Note that the PVT region
556 * index is explicitly passed instead of a parameter cache object)
557 ****************************************/
559 template <class FluidState, class LhsEval = typename FluidState::Scalar>
560 static LhsEval density(const FluidState& fluidState,
561 unsigned phaseIdx,
562 unsigned regionIdx)
563 {
564 assert(phaseIdx <= numPhases);
565 assert(regionIdx <= numRegions());
566
567 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
568 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
569 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
570
571 switch (phaseIdx) {
572 case oilPhaseIdx: {
573 if (enableDissolvedGas()) {
574 // miscible oil
575 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
576 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
577
578 return
579 bo*referenceDensity(oilPhaseIdx, regionIdx)
580 + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
581 }
582
583 // immiscible oil
584 const LhsEval Rs(0.0);
585 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
586
587 return referenceDensity(phaseIdx, regionIdx)*bo;
588 }
589
590 case gasPhaseIdx: {
592 // gas containing vaporized oil and vaporized water
593 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
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 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
600 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
601 }
602 if (enableVaporizedOil()) {
603 // miscible gas
604 const LhsEval Rvw(0.0);
605 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
606 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
607
608 return
609 bg*referenceDensity(gasPhaseIdx, regionIdx)
610 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
611 }
612 if (enableVaporizedWater()) {
613 // gas containing vaporized water
614 const LhsEval Rv(0.0);
615 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
616 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
617
618 return
619 bg*referenceDensity(gasPhaseIdx, regionIdx)
620 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
621 }
622
623 // immiscible gas
624 const LhsEval Rv(0.0);
625 const LhsEval Rvw(0.0);
626 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
627 return bg*referenceDensity(phaseIdx, regionIdx);
628 }
629
630 case waterPhaseIdx:
632 // gas miscible in water
633 const LhsEval& Rsw =BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
634 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
635 return
636 bw*referenceDensity(waterPhaseIdx, regionIdx)
637 + Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
638 }
639 const LhsEval Rsw(0.0);
640 return
642 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
643 }
644
645 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
646 }
647
657 template <class FluidState, class LhsEval = typename FluidState::Scalar>
658 static LhsEval saturatedDensity(const FluidState& fluidState,
659 unsigned phaseIdx,
660 unsigned regionIdx)
661 {
662 assert(phaseIdx <= numPhases);
663 assert(regionIdx <= numRegions());
664
665 const auto& p = fluidState.pressure(phaseIdx);
666 const auto& T = fluidState.temperature(phaseIdx);
667
668 switch (phaseIdx) {
669 case oilPhaseIdx: {
670 if (enableDissolvedGas()) {
671 // miscible oil
672 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
673 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
674
675 return
676 bo*referenceDensity(oilPhaseIdx, regionIdx)
677 + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
678 }
679
680 // immiscible oil
681 const LhsEval Rs(0.0);
682 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
683 return referenceDensity(phaseIdx, regionIdx)*bo;
684 }
685
686 case gasPhaseIdx: {
688 // gas containing vaporized oil and vaporized water
689 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
690 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
691 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
692
693 return
694 bg*referenceDensity(gasPhaseIdx, regionIdx)
695 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
696 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx) ;
697 }
698
699 if (enableVaporizedOil()) {
700 // miscible gas
701 const LhsEval Rvw(0.0);
702 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
703 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
704
705 return
706 bg*referenceDensity(gasPhaseIdx, regionIdx)
707 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
708 }
709
710 if (enableVaporizedWater()) {
711 // gas containing vaporized water
712 const LhsEval Rv(0.0);
713 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
714 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
715
716 return
717 bg*referenceDensity(gasPhaseIdx, regionIdx)
718 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
719 }
720
721 // immiscible gas
722 const LhsEval Rv(0.0);
723 const LhsEval Rvw(0.0);
724 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
725
726 return referenceDensity(phaseIdx, regionIdx)*bg;
727
728 }
729
730 case waterPhaseIdx:
731 {
733 // miscible in water
734 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
735 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
736 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
737 return
738 bw*referenceDensity(waterPhaseIdx, regionIdx)
739 + Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
740 }
741 return
743 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
744 }
745 }
746
747 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
748 }
749
758 template <class FluidState, class LhsEval = typename FluidState::Scalar>
759 static LhsEval inverseFormationVolumeFactor(const FluidState& fluidState,
760 unsigned phaseIdx,
761 unsigned regionIdx)
762 {
763 OPM_TIMEBLOCK_LOCAL(inverseFormationVolumeFactor);
764 assert(phaseIdx <= numPhases);
765 assert(regionIdx <= numRegions());
766
767 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
768 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
769
770 switch (phaseIdx) {
771 case oilPhaseIdx: {
772 if (enableDissolvedGas()) {
773 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
774 if (useSaturatedTables() && fluidState.saturation(gasPhaseIdx) > 0.0
775 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
776 {
777 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
778 } else {
779 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
780 }
781 }
782
783 const LhsEval Rs(0.0);
784 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
785 }
786 case gasPhaseIdx: {
788 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
789 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
790 if (useSaturatedTables() && fluidState.saturation(waterPhaseIdx) > 0.0
791 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
792 && fluidState.saturation(oilPhaseIdx) > 0.0
793 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
794 {
795 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
796 } else {
797 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
798 }
799 }
800
801 if (enableVaporizedOil()) {
802 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
803 if (useSaturatedTables() && fluidState.saturation(oilPhaseIdx) > 0.0
804 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
805 {
806 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
807 } else {
808 const LhsEval Rvw(0.0);
809 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
810 }
811 }
812
813 if (enableVaporizedWater()) {
814 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
815 if (useSaturatedTables() && fluidState.saturation(waterPhaseIdx) > 0.0
816 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
817 {
818 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
819 } else {
820 const LhsEval Rv(0.0);
821 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
822 }
823 }
824
825 const LhsEval Rv(0.0);
826 const LhsEval Rvw(0.0);
827 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
828 }
829 case waterPhaseIdx:
830 {
831 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
833 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
834 if (useSaturatedTables() && fluidState.saturation(gasPhaseIdx) > 0.0
835 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
836 {
837 return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
838 } else {
839 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
840 }
841 }
842 const LhsEval Rsw(0.0);
843 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
844 }
845 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
846 }
847 }
848
858 template <class FluidState, class LhsEval = typename FluidState::Scalar>
859 static LhsEval saturatedInverseFormationVolumeFactor(const FluidState& fluidState,
860 unsigned phaseIdx,
861 unsigned regionIdx)
862 {
863 OPM_TIMEBLOCK_LOCAL(saturatedInverseFormationVolumeFactor);
864 assert(phaseIdx <= numPhases);
865 assert(regionIdx <= numRegions());
866
867 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
868 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
869 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
870
871 switch (phaseIdx) {
872 case oilPhaseIdx: return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
873 case gasPhaseIdx: return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
874 case waterPhaseIdx: return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
875 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
876 }
877 }
878
880 template <class FluidState, class LhsEval = typename FluidState::Scalar>
881 static LhsEval fugacityCoefficient(const FluidState& fluidState,
882 unsigned phaseIdx,
883 unsigned compIdx,
884 unsigned regionIdx)
885 {
886 assert(phaseIdx <= numPhases);
887 assert(compIdx <= numComponents);
888 assert(regionIdx <= numRegions());
889
890 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
891 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
892
893 // for the fugacity coefficient of the oil component in the oil phase, we use
894 // some pseudo-realistic value for the vapor pressure to ease physical
895 // interpretation of the results
896 const LhsEval phi_oO = 20e3/p;
897
898 // for the gas component in the gas phase, assume it to be an ideal gas
899 constexpr const Scalar phi_gG = 1.0;
900
901 // for the fugacity coefficient of the water component in the water phase, we use
902 // the same approach as for the oil component in the oil phase
903 const LhsEval phi_wW = 30e3/p;
904
905 switch (phaseIdx) {
906 case gasPhaseIdx: // fugacity coefficients for all components in the gas phase
907 switch (compIdx) {
908 case gasCompIdx:
909 return phi_gG;
910
911 // for the oil component, we calculate the Rv value for saturated gas and Rs
912 // for saturated oil, and compute the fugacity coefficients at the
913 // equilibrium. for this, we assume that in equilibrium the fugacities of the
914 // oil component is the same in both phases.
915 case oilCompIdx: {
916 if (!enableVaporizedOil())
917 // if there's no vaporized oil, the gas phase is assumed to be
918 // immiscible with the oil component
919 return phi_gG*1e6;
920
921 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
922 const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
923 const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
924
925 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
926 const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
927 const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
928 const auto& x_oOSat = 1.0 - x_oGSat;
929
930 const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
931 const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
932
933 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
934 }
935
936 case waterCompIdx:
937 // the water component is assumed to be never miscible with the gas phase
938 return phi_gG*1e6;
939
940 default:
941 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
942 }
943
944 case oilPhaseIdx: // fugacity coefficients for all components in the oil phase
945 switch (compIdx) {
946 case oilCompIdx:
947 return phi_oO;
948
949 // for the oil and water components, we proceed analogous to the gas and
950 // water components in the gas phase
951 case gasCompIdx: {
952 if (!enableDissolvedGas())
953 // if there's no dissolved gas, the oil phase is assumed to be
954 // immiscible with the gas component
955 return phi_oO*1e6;
956
957 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
958 const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
959 const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
960 const auto& x_gGSat = 1.0 - x_gOSat;
961
962 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
963 const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
964 const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
965
966 const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
967 const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
968
969 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
970 }
971
972 case waterCompIdx:
973 return phi_oO*1e6;
974
975 default:
976 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
977 }
978
979 case waterPhaseIdx: // fugacity coefficients for all components in the water phase
980 // the water phase fugacity coefficients are pretty simple: because the water
981 // phase is assumed to consist entirely from the water component, we just
982 // need to make sure that the fugacity coefficients for the other components
983 // are a few orders of magnitude larger than the one of the water
984 // component. (i.e., the affinity of the gas and oil components to the water
985 // phase is lower by a few orders of magnitude)
986 switch (compIdx) {
987 case waterCompIdx: return phi_wW;
988 case oilCompIdx: return 1.1e6*phi_wW;
989 case gasCompIdx: return 1e6*phi_wW;
990 default:
991 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
992 }
993
994 default:
995 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
996 }
997
998 throw std::logic_error("Unhandled phase or component index");
999 }
1000
1002 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1003 static LhsEval viscosity(const FluidState& fluidState,
1004 unsigned phaseIdx,
1005 unsigned regionIdx)
1006 {
1007 OPM_TIMEBLOCK_LOCAL(viscosity);
1008 assert(phaseIdx <= numPhases);
1009 assert(regionIdx <= numRegions());
1010
1011 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1012 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1013
1014 switch (phaseIdx) {
1015 case oilPhaseIdx: {
1016 if (enableDissolvedGas()) {
1017 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1018 if (useSaturatedTables() && fluidState.saturation(gasPhaseIdx) > 0.0
1019 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
1020 {
1021 return oilPvt_->saturatedViscosity(regionIdx, T, p);
1022 } else {
1023 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1024 }
1025 }
1026
1027 const LhsEval Rs(0.0);
1028 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1029 }
1030
1031 case gasPhaseIdx: {
1033 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1034 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1035 if (useSaturatedTables() && fluidState.saturation(waterPhaseIdx) > 0.0
1036 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1037 && fluidState.saturation(oilPhaseIdx) > 0.0
1038 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1039 {
1040 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1041 } else {
1042 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1043 }
1044 }
1045 if (enableVaporizedOil()) {
1046 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1047 if (useSaturatedTables() && fluidState.saturation(oilPhaseIdx) > 0.0
1048 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1049 {
1050 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1051 } else {
1052 const LhsEval Rvw(0.0);
1053 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1054 }
1055 }
1056 if (enableVaporizedWater()) {
1057 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1058 if (useSaturatedTables() && fluidState.saturation(waterPhaseIdx) > 0.0
1059 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1060 {
1061 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1062 } else {
1063 const LhsEval Rv(0.0);
1064 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1065 }
1066 }
1067
1068 const LhsEval Rv(0.0);
1069 const LhsEval Rvw(0.0);
1070 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1071 }
1072
1073 case waterPhaseIdx:
1074 {
1075 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1077 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1078 if (useSaturatedTables() && fluidState.saturation(gasPhaseIdx) > 0.0
1079 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
1080 {
1081 return waterPvt_->saturatedViscosity(regionIdx, T, p, saltConcentration);
1082 } else {
1083 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1084 }
1085 }
1086 const LhsEval Rsw(0.0);
1087 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1088 }
1089 }
1090
1091 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1092 }
1093
1094 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1095 static LhsEval internalEnergy(const FluidState& fluidState,
1096 unsigned phaseIdx,
1097 unsigned regionIdx){
1098 bool is_mixing = false;
1099 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1100 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1101
1102 switch (phaseIdx) {
1103 case oilPhaseIdx: {
1104 if(!oilPvt_->mixingEnergy()){
1105 const auto& oilEnergy =
1106 oilPvt_->internalEnergy(regionIdx, T, p,
1107 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1108
1109 return oilEnergy;
1110 }
1111 is_mixing = true;
1112 break;
1113 }
1114 case waterPhaseIdx: {
1115 if(!waterPvt_->mixingEnergy()){
1116 const auto waterEnergy =
1117 waterPvt_->internalEnergy(regionIdx, T, p,
1118 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1119 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1120
1121 return waterEnergy;
1122 }
1123 is_mixing = true;
1124 break;
1125 }
1126 case gasPhaseIdx: {
1127 if(!gasPvt_->mixingEnergy()){
1128 const auto gasEnergy =
1129 gasPvt_->internalEnergy(regionIdx, T, p,
1130 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1131 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1132
1133 return gasEnergy;
1134 }
1135 is_mixing = true;
1136 break;
1137 }
1138 }
1139 assert(is_mixing==true);
1140 const auto& energy = internalMixingTotalEnergy<FluidState,LhsEval>(fluidState, phaseIdx, regionIdx);
1141 const auto& phase_density = density<FluidState,LhsEval>(fluidState, phaseIdx, regionIdx);
1142 return energy/phase_density;
1143 }
1144
1145
1146 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1147 static LhsEval internalMixingTotalEnergy(const FluidState& fluidState,
1148 unsigned phaseIdx,
1149 unsigned regionIdx)
1150 {
1151 assert(phaseIdx <= numPhases);
1152 assert(regionIdx <= numRegions());
1153 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1154 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1155 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1156 // to avoid putting all thermal into the interface of the multiplexer
1157 switch (phaseIdx) {
1158 case oilPhaseIdx: {
1159 auto oilEnergy = oilPvt_->internalEnergy(regionIdx, T, p,
1160 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1161 assert(oilPvt_->mixingEnergy());
1162 //mixing energy adsed
1163 if (enableDissolvedGas()) {
1164 // miscible oil
1165 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1166 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
1167 const auto& gasEnergy =
1168 gasPvt_->internalEnergy(regionIdx, T, p,
1169 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1170 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1171 const auto hVapG = gasPvt_->hVap(regionIdx);// pressure correction ? assume equal to energy change
1172 return
1173 oilEnergy*bo*referenceDensity(oilPhaseIdx, regionIdx)
1174 + (gasEnergy-hVapG)*Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
1175 }
1176
1177 // immiscible oil
1178 const LhsEval Rs(0.0);
1179 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
1180
1181 return oilEnergy*referenceDensity(phaseIdx, regionIdx)*bo;
1182 }
1183
1184 case gasPhaseIdx: {
1185 const auto& gasEnergy =
1186 gasPvt_->internalEnergy(regionIdx, T, p,
1187 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1188 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1189 assert(gasPvt_->mixingEnergy());
1191 const auto& oilEnergy =
1192 oilPvt_->internalEnergy(regionIdx, T, p,
1193 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1194 const auto waterEnergy =
1195 waterPvt_->internalEnergy(regionIdx, T, p,
1196 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1197 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1198 // gas containing vaporized oil and vaporized water
1199 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1200 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1201 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1202 const auto hVapO = oilPvt_->hVap(regionIdx);
1203 const auto hVapW = waterPvt_->hVap(regionIdx);
1204 return
1205 gasEnergy*bg*referenceDensity(gasPhaseIdx, regionIdx)
1206 + (oilEnergy+hVapO)*Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
1207 + (waterEnergy+hVapW)*Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
1208 }
1209 if (enableVaporizedOil()) {
1210 const auto& oilEnergy =
1211 oilPvt_->internalEnergy(regionIdx, T, p,
1212 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1213 // miscible gas
1214 const LhsEval Rvw(0.0);
1215 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1216 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1217 const auto hVapO = oilPvt_->hVap(regionIdx);
1218 return
1219 gasEnergy*bg*referenceDensity(gasPhaseIdx, regionIdx)
1220 + (oilEnergy+hVapO)*Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
1221 }
1222 if (enableVaporizedWater()) {
1223 // gas containing vaporized water
1224 const LhsEval Rv(0.0);
1225 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1226 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1227 const auto waterEnergy =
1228 waterPvt_->internalEnergy(regionIdx, T, p,
1229 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1230 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1231 const auto hVapW = waterPvt_->hVap(regionIdx);
1232 return
1233 gasEnergy*bg*referenceDensity(gasPhaseIdx, regionIdx)
1234 + (waterEnergy+hVapW)*Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
1235 }
1236
1237 // immiscible gas
1238 const LhsEval Rv(0.0);
1239 const LhsEval Rvw(0.0);
1240 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1241 return gasEnergy*bg*referenceDensity(phaseIdx, regionIdx);
1242 }
1243
1244 case waterPhaseIdx:
1245 const auto waterEnergy =
1246 waterPvt_->internalEnergy(regionIdx, T, p,
1247 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1248 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1249 assert(waterPvt_->mixingEnergy());
1251 const auto& gasEnergy =
1252 gasPvt_->internalEnergy(regionIdx, T, p,
1253 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1254 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1255 // gas miscible in water
1256 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
1257 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
1258 return
1259 waterEnergy*bw*referenceDensity(waterPhaseIdx, regionIdx)
1260 + gasEnergy*Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
1261 }
1262 const LhsEval Rsw(0.0);
1263 return
1264 waterEnergy*referenceDensity(waterPhaseIdx, regionIdx)
1265 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
1266 }
1267 throw std::logic_error("Unhandled phase index " + std::to_string(phaseIdx));
1268 }
1269
1270
1271
1273 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1274 static LhsEval enthalpy(const FluidState& fluidState,
1275 unsigned phaseIdx,
1276 unsigned regionIdx)
1277 {
1278 // should preferably not be used values should be taken from intensive quantities fluid state.
1279 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1280 auto energy = internalEnergy<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1281 if(!enthalpy_eq_energy_){
1282 // used for simplified models
1283 energy += p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1284 }
1285 return energy;
1286 }
1287
1294 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1295 static LhsEval saturatedVaporizationFactor(const FluidState& fluidState,
1296 unsigned phaseIdx,
1297 unsigned regionIdx)
1298 {
1299 assert(phaseIdx <= numPhases);
1300 assert(regionIdx <= numRegions());
1301
1302 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1303 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1304 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
1305
1306 switch (phaseIdx) {
1307 case oilPhaseIdx: return 0.0;
1308 case gasPhaseIdx: return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p, saltConcentration);
1309 case waterPhaseIdx: return 0.0;
1310 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1311 }
1312 }
1313
1320 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1321 static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1322 unsigned phaseIdx,
1323 unsigned regionIdx,
1324 const LhsEval& maxOilSaturation)
1325 {
1326 OPM_TIMEBLOCK_LOCAL(saturatedDissolutionFactor);
1327 assert(phaseIdx <= numPhases);
1328 assert(regionIdx <= numRegions());
1329
1330 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1331 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1332 const auto& So = (phaseIdx == waterPhaseIdx) ? 0 : decay<LhsEval>(fluidState.saturation(oilPhaseIdx));
1333
1334 switch (phaseIdx) {
1335 case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1336 case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1337 case waterPhaseIdx: return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1338 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1339 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1340 }
1341 }
1342
1351 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1352 static LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1353 unsigned phaseIdx,
1354 unsigned regionIdx)
1355 {
1356 OPM_TIMEBLOCK_LOCAL(saturatedDissolutionFactor);
1357 assert(phaseIdx <= numPhases);
1358 assert(regionIdx <= numRegions());
1359
1360 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1361 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1362
1363 switch (phaseIdx) {
1364 case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1365 case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1366 case waterPhaseIdx: return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1367 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1368 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1369 }
1370 }
1371
1375 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1376 static LhsEval bubblePointPressure(const FluidState& fluidState,
1377 unsigned regionIdx)
1378 {
1379 return saturationPressure(fluidState, oilPhaseIdx, regionIdx);
1380 }
1381
1382
1386 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1387 static LhsEval dewPointPressure(const FluidState& fluidState,
1388 unsigned regionIdx)
1389 {
1390 return saturationPressure(fluidState, gasPhaseIdx, regionIdx);
1391 }
1392
1403 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1404 static LhsEval saturationPressure(const FluidState& fluidState,
1405 unsigned phaseIdx,
1406 unsigned regionIdx)
1407 {
1408 assert(phaseIdx <= numPhases);
1409 assert(regionIdx <= numRegions());
1410
1411 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1412
1413 switch (phaseIdx) {
1414 case oilPhaseIdx: return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1415 case gasPhaseIdx: return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1416 case waterPhaseIdx: return waterPvt_->saturationPressure(regionIdx, T,
1417 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1418 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1419 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1420 }
1421 }
1422
1423 /****************************************
1424 * Auxiliary and convenience methods for the black-oil model
1425 ****************************************/
1430 template <class LhsEval>
1431 static LhsEval convertXoGToRs(const LhsEval& XoG, unsigned regionIdx)
1432 {
1433 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1434 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1435
1436 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1437 }
1438
1443 template <class LhsEval>
1444 static LhsEval convertXwGToRsw(const LhsEval& XwG, unsigned regionIdx)
1445 {
1446 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1447 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1448
1449 return XwG/(1.0 - XwG)*(rho_wRef/rho_gRef);
1450 }
1451
1456 template <class LhsEval>
1457 static LhsEval convertXgOToRv(const LhsEval& XgO, unsigned regionIdx)
1458 {
1459 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1460 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1461
1462 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1463 }
1464
1469 template <class LhsEval>
1470 static LhsEval convertXgWToRvw(const LhsEval& XgW, unsigned regionIdx)
1471 {
1472 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1473 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1474
1475 return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1476 }
1477
1478
1483 template <class LhsEval>
1484 static LhsEval convertRsToXoG(const LhsEval& Rs, unsigned regionIdx)
1485 {
1486 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1487 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1488
1489 const LhsEval& rho_oG = Rs*rho_gRef;
1490 return rho_oG/(rho_oRef + rho_oG);
1491 }
1492
1497 template <class LhsEval>
1498 static LhsEval convertRswToXwG(const LhsEval& Rsw, unsigned regionIdx)
1499 {
1500 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1501 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1502
1503 const LhsEval& rho_wG = Rsw*rho_gRef;
1504 return rho_wG/(rho_wRef + rho_wG);
1505 }
1506
1511 template <class LhsEval>
1512 static LhsEval convertRvToXgO(const LhsEval& Rv, unsigned regionIdx)
1513 {
1514 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1515 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1516
1517 const LhsEval& rho_gO = Rv*rho_oRef;
1518 return rho_gO/(rho_gRef + rho_gO);
1519 }
1520
1525 template <class LhsEval>
1526 static LhsEval convertRvwToXgW(const LhsEval& Rvw, unsigned regionIdx)
1527 {
1528 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1529 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1530
1531 const LhsEval& rho_gW = Rvw*rho_wRef;
1532 return rho_gW/(rho_gRef + rho_gW);
1533 }
1534
1538 template <class LhsEval>
1539 static LhsEval convertXgWToxgW(const LhsEval& XgW, unsigned regionIdx)
1540 {
1541 Scalar MW = molarMass_[regionIdx][waterCompIdx];
1542 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1543
1544 return XgW*MG / (MW*(1 - XgW) + XgW*MG);
1545 }
1546
1550 template <class LhsEval>
1551 static LhsEval convertXwGToxwG(const LhsEval& XwG, unsigned regionIdx)
1552 {
1553 Scalar MW = molarMass_[regionIdx][waterCompIdx];
1554 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1555
1556 return XwG*MW / (MG*(1 - XwG) + XwG*MW);
1557 }
1558
1562 template <class LhsEval>
1563 static LhsEval convertXoGToxoG(const LhsEval& XoG, unsigned regionIdx)
1564 {
1565 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1566 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1567
1568 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1569 }
1570
1574 template <class LhsEval>
1575 static LhsEval convertxoGToXoG(const LhsEval& xoG, unsigned regionIdx)
1576 {
1577 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1578 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1579
1580 return xoG*MG / (xoG*(MG - MO) + MO);
1581 }
1582
1586 template <class LhsEval>
1587 static LhsEval convertXgOToxgO(const LhsEval& XgO, unsigned regionIdx)
1588 {
1589 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1590 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1591
1592 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1593 }
1594
1598 template <class LhsEval>
1599 static LhsEval convertxgOToXgO(const LhsEval& xgO, unsigned regionIdx)
1600 {
1601 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1602 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1603
1604 return xgO*MO / (xgO*(MO - MG) + MG);
1605 }
1606
1614 static const GasPvt& gasPvt()
1615 { return *gasPvt_; }
1616
1624 static const OilPvt& oilPvt()
1625 { return *oilPvt_; }
1626
1634 static const WaterPvt& waterPvt()
1635 { return *waterPvt_; }
1636
1642 static Scalar reservoirTemperature(unsigned = 0)
1643 { return reservoirTemperature_; }
1644
1651 { reservoirTemperature_ = value; }
1652
1653 static short activeToCanonicalPhaseIdx(unsigned activePhaseIdx);
1654
1655 static short canonicalToActivePhaseIdx(unsigned phaseIdx);
1656
1658 static Scalar diffusionCoefficient(unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1659 { return diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx]; }
1660
1662 static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1663 { diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx] = coefficient ; }
1664
1668 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
1669 static LhsEval diffusionCoefficient(const FluidState& fluidState,
1670 const ParameterCache<ParamCacheEval>& paramCache,
1671 unsigned phaseIdx,
1672 unsigned compIdx)
1673 {
1674 // diffusion is disabled by the user
1675 if(!enableDiffusion())
1676 return 0.0;
1677
1678 // diffusion coefficients are set, and we use them
1679 if(!diffusionCoefficients_.empty()) {
1680 return diffusionCoefficient(compIdx, phaseIdx, paramCache.regionIndex());
1681 }
1682
1683 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1684 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1685
1686 switch (phaseIdx) {
1687 case oilPhaseIdx: return oilPvt().diffusionCoefficient(T, p, compIdx);
1688 case gasPhaseIdx: return gasPvt().diffusionCoefficient(T, p, compIdx);
1689 case waterPhaseIdx: return waterPvt().diffusionCoefficient(T, p, compIdx);
1690 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1691 }
1692 }
1693 static void setEnergyEqualEnthalpy(bool enthalpy_eq_energy){
1694 enthalpy_eq_energy_ = enthalpy_eq_energy;
1695 }
1696
1697 static bool enthalpyEqualEnergy(){
1698 return enthalpy_eq_energy_;
1699 }
1700
1701private:
1702 static void resizeArrays_(std::size_t numRegions);
1703
1704 static Scalar reservoirTemperature_;
1705
1706 static std::shared_ptr<GasPvt> gasPvt_;
1707 static std::shared_ptr<OilPvt> oilPvt_;
1708 static std::shared_ptr<WaterPvt> waterPvt_;
1709
1710 static bool enableDissolvedGas_;
1711 static bool enableDissolvedGasInWater_;
1712 static bool enableVaporizedOil_;
1713 static bool enableVaporizedWater_;
1714 static bool enableDiffusion_;
1715
1716 // HACK for GCC 4.4: the array size has to be specified using the literal value '3'
1717 // here, because GCC 4.4 seems to be unable to determine the number of phases from
1718 // the BlackOil fluid system in the attribute declaration below...
1719 static std::vector<std::array<Scalar, /*numPhases=*/3> > referenceDensity_;
1720 static std::vector<std::array<Scalar, /*numComponents=*/3> > molarMass_;
1721 static std::vector<std::array<Scalar, /*numComponents=*/3 * /*numPhases=*/3> > diffusionCoefficients_;
1722
1723 static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1724 static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1725
1726 static bool isInitialized_;
1727 static bool useSaturatedTables_;
1728 inline static bool enthalpy_eq_energy_ = false;
1729};
1730
1731template <typename T> using BOFS = BlackOilFluidSystem<T, BlackOilDefaultIndexTraits>;
1732
1733#define DECLARE_INSTANCE(T) \
1734template<> unsigned char BOFS<T>::numActivePhases_; \
1735template<> std::array<bool, BOFS<T>::numPhases> BOFS<T>::phaseIsActive_; \
1736template<> std::array<short, BOFS<T>::numPhases> BOFS<T>::activeToCanonicalPhaseIdx_; \
1737template<> std::array<short, BOFS<T>::numPhases> BOFS<T>::canonicalToActivePhaseIdx_; \
1738template<> T BOFS<T>::surfaceTemperature; \
1739template<> T BOFS<T>::surfacePressure; \
1740template<> T BOFS<T>::reservoirTemperature_; \
1741template<> bool BOFS<T>::enableDissolvedGas_; \
1742template<> bool BOFS<T>::enableDissolvedGasInWater_; \
1743template<> bool BOFS<T>::enableVaporizedOil_; \
1744template<> bool BOFS<T>::enableVaporizedWater_; \
1745template<> bool BOFS<T>::enableDiffusion_; \
1746template<> std::shared_ptr<OilPvtMultiplexer<T>> BOFS<T>::oilPvt_; \
1747template<> std::shared_ptr<GasPvtMultiplexer<T>> BOFS<T>::gasPvt_; \
1748template<> std::shared_ptr<WaterPvtMultiplexer<T>> BOFS<T>::waterPvt_; \
1749template<> std::vector<std::array<T, 3>> BOFS<T>::referenceDensity_; \
1750template<> std::vector<std::array<T, 3>> BOFS<T>::molarMass_; \
1751template<> std::vector<std::array<T, 9>> BOFS<T>::diffusionCoefficients_; \
1752template<> bool BOFS<T>::isInitialized_; \
1753template<> bool BOFS<T>::useSaturatedTables_;
1754
1755DECLARE_INSTANCE(float)
1756DECLARE_INSTANCE(double)
1757
1758#undef DECLARE_INSTANCE
1759
1760} // namespace Opm
1761
1762#endif
The base class for all fluid systems.
The class which specifies the default phase and component indices for the black-oil fluid system.
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...
A parameter cache which does nothing.
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:43
ScalarT Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:48
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:389
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition BlackOilFluidSystem.hpp:399
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:1658
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:1431
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:1484
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition BlackOilFluidSystem.hpp:384
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:1526
static LhsEval convertXwGToxwG(const LhsEval &XwG, unsigned regionIdx)
Convert a gas mass fraction in the water phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1551
static void setUseSaturatedTables(bool yesno)
Specify whether the saturated tables should be used.
Definition BlackOilFluidSystem.hpp:295
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition BlackOilFluidSystem.hpp:859
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1563
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition BlackOilFluidSystem.hpp:1003
static void initEnd()
Finish initializing the black oil fluid system.
Definition BlackOilFluidSystem.cpp:267
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition BlackOilFluidSystem.hpp:1642
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:1498
static bool phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition BlackOilFluidSystem.hpp:403
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:1376
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition BlackOilFluidSystem.hpp:560
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition BlackOilFluidSystem.hpp:658
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:536
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition BlackOilFluidSystem.hpp:423
static bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition BlackOilFluidSystem.hpp:484
static constexpr unsigned numPhases
Number of fluid phases in the fluid system.
Definition BlackOilFluidSystem.hpp:354
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition BlackOilFluidSystem.hpp:313
static void initBegin(std::size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition BlackOilFluidSystem.cpp:229
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition BlackOilFluidSystem.hpp:1650
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:1321
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:1274
static constexpr unsigned oilCompIdx
Index of the oil component.
Definition BlackOilFluidSystem.hpp:387
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition BlackOilFluidSystem.hpp:373
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:543
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:881
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition BlackOilFluidSystem.hpp:1575
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition BlackOilFluidSystem.hpp:419
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition BlackOilFluidSystem.hpp:1352
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition BlackOilFluidSystem.hpp:1387
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:1404
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:1444
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition BlackOilFluidSystem.cpp:362
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1587
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:1669
static LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition BlackOilFluidSystem.hpp:1295
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition BlackOilFluidSystem.hpp:456
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition BlackOilFluidSystem.hpp:508
static Scalar surfaceTemperature
The temperature at the surface.
Definition BlackOilFluidSystem.hpp:367
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition BlackOilFluidSystem.hpp:307
static constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition BlackOilFluidSystem.hpp:357
static bool useSaturatedTables()
Returns whether the saturated tables should be used.
Definition BlackOilFluidSystem.hpp:500
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition BlackOilFluidSystem.hpp:1599
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition BlackOilFluidSystem.hpp:759
static LhsEval convertXgWToxgW(const LhsEval &XgW, unsigned regionIdx)
Convert a water mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1539
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition BlackOilFluidSystem.hpp:361
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:1457
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:523
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:516
static constexpr unsigned gasCompIdx
Index of the gas component.
Definition BlackOilFluidSystem.hpp:391
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:1512
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition BlackOilFluidSystem.hpp:1634
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition BlackOilFluidSystem.hpp:475
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition BlackOilFluidSystem.hpp:435
static unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition BlackOilFluidSystem.cpp:323
static Scalar surfacePressure
The pressure at the surface.
Definition BlackOilFluidSystem.hpp:364
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:359
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition BlackOilFluidSystem.hpp:1624
static bool enableDissolvedGasInWater()
Returns whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem.hpp:466
static unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition BlackOilFluidSystem.cpp:340
static bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem.hpp:492
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:1470
static void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition BlackOilFluidSystem.cpp:256
static std::size_t numRegions()
Returns the number of PVT regions which are considered.
Definition BlackOilFluidSystem.hpp:447
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 std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition BlackOilFluidSystem.cpp:306
static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition BlackOilFluidSystem.hpp:1662
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition BlackOilFluidSystem.hpp:1614
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition BlackOilFluidSystem.hpp:431
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition BlackOilFluidSystem.hpp:301
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition GasPvtMultiplexer.hpp:110
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:268
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:105
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:240
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition WaterPvtMultiplexer.hpp:90
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:231
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