My Project
Loading...
Searching...
No Matches
BlackOilFluidSystem_macrotemplate.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*/
32
33#include <opm/common/TimingMacros.hpp>
34#include <opm/common/utility/VectorWithDefaultAllocator.hpp>
35
38
44
45#include <array>
46#include <cstddef>
47#include <memory>
48#include <stdexcept>
49#include <string>
50#include <string_view>
51#include <vector>
52
53namespace Opm {
54
55#if HAVE_ECL_INPUT
56class EclipseState;
57class Schedule;
58#endif
59
60
61
62
69template <class Scalar, class IndexTraits = BlackOilDefaultIndexTraits,
70 template<typename> typename Storage = VectorWithDefaultAllocator,
71 template<typename> typename SmartPointer = std::shared_ptr>
72class FLUIDSYSTEM_CLASSNAME : public BaseFluidSystem<Scalar, FLUIDSYSTEM_CLASSNAME<Scalar, IndexTraits, Storage, SmartPointer> >
73{
75
76public:
80
81 #ifdef COMPILING_STATIC_FLUID_SYSTEM
83 template <class EvaluationT>
84 struct ParameterCache : public NullParameterCache<EvaluationT>
85 {
86 using Evaluation = EvaluationT;
87
88 public:
89 explicit ParameterCache(Scalar maxOilSat = 1.0, unsigned regionIdx = 0)
90 : maxOilSat_(maxOilSat)
91 , regionIdx_(regionIdx)
92 {
93 }
94
102 template <class OtherCache>
103 void assignPersistentData(const OtherCache& other)
104 {
105 regionIdx_ = other.regionIndex();
106 maxOilSat_ = other.maxOilSat();
107 }
108
116 unsigned regionIndex() const
117 { return regionIdx_; }
118
126 void setRegionIndex(unsigned val)
127 { regionIdx_ = val; }
128
129 const Evaluation& maxOilSat() const
130 { return maxOilSat_; }
131
132 void setMaxOilSat(const Evaluation& val)
133 { maxOilSat_ = val; }
134
135 private:
136 Evaluation maxOilSat_;
137 unsigned regionIdx_;
138 };
139
140 #else
141
142 // We want to use the exact same ParameterCache class for both the static and nonstatic versions
143 template<class EvaluationT>
144 using ParameterCache =
146 IndexTraits,
147 Storage,
148 SmartPointer>::template ParameterCache<EvaluationT>;
149 #endif
150
151 #ifndef COMPILING_STATIC_FLUID_SYSTEM
157 template<template<typename> typename StorageT, template<typename> typename SmartPointerT>
161 , numActivePhases_(other.numActivePhases_)
162 , phaseIsActive_(other.phaseIsActive_)
163 , reservoirTemperature_(other.reservoirTemperature_)
164 , gasPvt_(SmartPointerT<typename decltype(gasPvt_)::element_type>(other.gasPvt_))
165 , oilPvt_(SmartPointerT<typename decltype(oilPvt_)::element_type>(other.oilPvt_))
166 , waterPvt_(SmartPointerT<typename decltype(waterPvt_)::element_type>(other.waterPvt_))
167 , enableDissolvedGas_(other.enableDissolvedGas_)
168 , enableDissolvedGasInWater_(other.enableDissolvedGasInWater_)
169 , enableVaporizedOil_(other.enableVaporizedOil_)
170 , enableVaporizedWater_(other.enableVaporizedWater_)
171 , enableDiffusion_(other.enableDiffusion_)
172 , referenceDensity_(StorageT<typename decltype(referenceDensity_)::value_type>(other.referenceDensity_))
173 , molarMass_(StorageT<typename decltype(molarMass_)::value_type>(other.molarMass_))
174 , diffusionCoefficients_(StorageT<typename decltype(diffusionCoefficients_)::value_type>(other.diffusionCoefficients_))
175 , activeToCanonicalPhaseIdx_(other.activeToCanonicalPhaseIdx_)
176 , canonicalToActivePhaseIdx_(other.canonicalToActivePhaseIdx_)
177 , isInitialized_(other.isInitialized_)
178 , useSaturatedTables_(other.useSaturatedTables_)
179 , enthalpy_eq_energy_(other.enthalpy_eq_energy_)
180 {
181 }
182 #endif
183
184 #ifdef COMPILING_STATIC_FLUID_SYSTEM
192 template<template<typename> typename StorageT = VectorWithDefaultAllocator, template<typename> typename SmartPointerT = std::shared_ptr>
194 {
196 return instance;
197
198 }
199 #endif
200
201 /****************************************
202 * Initialization
203 ****************************************/
204#if HAVE_ECL_INPUT
208 STATIC_OR_DEVICE void initFromState(const EclipseState& eclState, const Schedule& schedule);
209#endif // HAVE_ECL_INPUT
210
219 STATIC_OR_DEVICE void initBegin(std::size_t numPvtRegions);
220
227 STATIC_OR_DEVICE void setEnableDissolvedGas(bool yesno)
228 { enableDissolvedGas_ = yesno; }
229
236 STATIC_OR_DEVICE void setEnableVaporizedOil(bool yesno)
237 { enableVaporizedOil_ = yesno; }
238
245 STATIC_OR_DEVICE void setEnableVaporizedWater(bool yesno)
246 { enableVaporizedWater_ = yesno; }
247
254 STATIC_OR_DEVICE void setEnableDissolvedGasInWater(bool yesno)
255 { enableDissolvedGasInWater_ = yesno; }
261 STATIC_OR_DEVICE void setEnableDiffusion(bool yesno)
262 { enableDiffusion_ = yesno; }
263
269 STATIC_OR_DEVICE void setUseSaturatedTables(bool yesno)
270 { useSaturatedTables_ = yesno; }
271
275 STATIC_OR_DEVICE void setGasPvt(SmartPointer<GasPvt> pvtObj)
276 { gasPvt_ = pvtObj; }
277
281 STATIC_OR_DEVICE void setOilPvt(SmartPointer<OilPvt> pvtObj)
282 { oilPvt_ = pvtObj; }
283
287 STATIC_OR_DEVICE void setWaterPvt(SmartPointer<WaterPvt> pvtObj)
288 { waterPvt_ = pvtObj; }
289
290 STATIC_OR_DEVICE void setVapPars(const Scalar par1, const Scalar par2)
291 {
292 if (gasPvt_) {
293 gasPvt_->setVapPars(par1, par2);
294 }
295 if (oilPvt_) {
296 oilPvt_->setVapPars(par1, par2);
297 }
298 if (waterPvt_) {
299 waterPvt_->setVapPars(par1, par2);
300 }
301 }
302
310 STATIC_OR_DEVICE void setReferenceDensities(Scalar rhoOil,
311 Scalar rhoWater,
312 Scalar rhoGas,
313 unsigned regionIdx);
314
318 STATIC_OR_DEVICE void initEnd();
319
320 STATIC_OR_DEVICE bool isInitialized()
321 { return isInitialized_; }
322
323 /****************************************
324 * Generic phase properties
325 ****************************************/
326
328 static constexpr unsigned numPhases = 3;
329
331 static constexpr unsigned waterPhaseIdx = IndexTraits::waterPhaseIdx;
333 static constexpr unsigned oilPhaseIdx = IndexTraits::oilPhaseIdx;
335 static constexpr unsigned gasPhaseIdx = IndexTraits::gasPhaseIdx;
336
338 STATIC_OR_NOTHING Scalar surfacePressure;
339
341 STATIC_OR_NOTHING Scalar surfaceTemperature;
342
344 STATIC_OR_DEVICE std::string_view phaseName(unsigned phaseIdx);
345
347 STATIC_OR_DEVICE bool isLiquid(unsigned phaseIdx)
348 {
349 assert(phaseIdx < numPhases);
350 return phaseIdx != gasPhaseIdx;
351 }
352
353 /****************************************
354 * Generic component related properties
355 ****************************************/
356
358 static constexpr unsigned numComponents = 3;
359
361 static constexpr int oilCompIdx = IndexTraits::oilCompIdx;
363 static constexpr int waterCompIdx = IndexTraits::waterCompIdx;
365 static constexpr int gasCompIdx = IndexTraits::gasCompIdx;
366
367protected:
368 STATIC_OR_NOTHING unsigned char numActivePhases_;
369 STATIC_OR_NOTHING std::array<bool,numPhases> phaseIsActive_;
370
371public:
373 STATIC_OR_DEVICE unsigned numActivePhases()
374 { return numActivePhases_; }
375
377 STATIC_OR_DEVICE bool phaseIsActive(unsigned phaseIdx)
378 {
379 assert(phaseIdx < numPhases);
380 return phaseIsActive_[phaseIdx];
381 }
382
384 STATIC_OR_DEVICE unsigned solventComponentIndex(unsigned phaseIdx);
385
387 STATIC_OR_DEVICE unsigned soluteComponentIndex(unsigned phaseIdx);
388
390 STATIC_OR_DEVICE std::string_view componentName(unsigned compIdx);
391
393 STATIC_OR_DEVICE Scalar molarMass(unsigned compIdx, unsigned regionIdx = 0)
394 { return molarMass_[regionIdx][compIdx]; }
395
397 STATIC_OR_DEVICE bool isIdealMixture(unsigned /*phaseIdx*/)
398 {
399 // fugacity coefficients are only pressure dependent -> we
400 // have an ideal mixture
401 return true;
402 }
403
405 STATIC_OR_DEVICE bool isCompressible(unsigned /*phaseIdx*/)
406 { return true; /* all phases are compressible */ }
407
409 STATIC_OR_DEVICE bool isIdealGas(unsigned /*phaseIdx*/)
410 { return false; }
411
412
413 /****************************************
414 * Black-oil specific properties
415 ****************************************/
421 STATIC_OR_DEVICE std::size_t numRegions()
422 { return molarMass_.size(); }
423
430 STATIC_OR_DEVICE bool enableDissolvedGas()
431 { return enableDissolvedGas_; }
432
433
440 STATIC_OR_DEVICE bool enableDissolvedGasInWater()
441 { return enableDissolvedGasInWater_; }
442
449 STATIC_OR_DEVICE bool enableVaporizedOil()
450 { return enableVaporizedOil_; }
451
458 STATIC_OR_DEVICE bool enableVaporizedWater()
459 { return enableVaporizedWater_; }
460
466 STATIC_OR_DEVICE bool enableDiffusion()
467 { return enableDiffusion_; }
468
474 STATIC_OR_DEVICE bool useSaturatedTables()
475 { return useSaturatedTables_; }
476
482 STATIC_OR_DEVICE Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
483 { return referenceDensity_[regionIdx][phaseIdx]; }
484
485 /****************************************
486 * thermodynamic quantities (generic version)
487 ****************************************/
489 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
490 STATIC_OR_DEVICE LhsEval density(const FluidState& fluidState,
491 const ParameterCache<ParamCacheEval>& paramCache,
492 unsigned phaseIdx)
493 { return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
494
496 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
497 STATIC_OR_DEVICE LhsEval fugacityCoefficient(const FluidState& fluidState,
498 const ParameterCache<ParamCacheEval>& paramCache,
499 unsigned phaseIdx,
500 unsigned compIdx)
501 {
502 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
503 phaseIdx,
504 compIdx,
505 paramCache.regionIndex());
506 }
507
509 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
510 STATIC_OR_DEVICE LhsEval viscosity(const FluidState& fluidState,
511 const ParameterCache<ParamCacheEval>& paramCache,
512 unsigned phaseIdx)
513 { return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
514
516 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
517 STATIC_OR_DEVICE LhsEval enthalpy(const FluidState& fluidState,
518 const ParameterCache<ParamCacheEval>& paramCache,
519 unsigned phaseIdx)
520 { return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
521
522 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
523 STATIC_OR_DEVICE LhsEval internalEnergy(const FluidState& fluidState,
524 const ParameterCache<ParamCacheEval>& paramCache,
525 unsigned phaseIdx)
526 { return internalEnergy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
527
528 /****************************************
529 * thermodynamic quantities (black-oil specific version: Note that the PVT region
530 * index is explicitly passed instead of a parameter cache object)
531 ****************************************/
533 template <class FluidState, class LhsEval = typename FluidState::Scalar>
534 STATIC_OR_DEVICE LhsEval density(const FluidState& fluidState,
535 unsigned phaseIdx,
536 unsigned regionIdx)
537 {
538 assert(phaseIdx <= numPhases);
539 assert(regionIdx <= numRegions());
540
541 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
542 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
543 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
544
545 switch (phaseIdx) {
546 case oilPhaseIdx: {
547 if (enableDissolvedGas()) {
548 // miscible oil
549 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
550 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
551
552 return
553 bo*referenceDensity(oilPhaseIdx, regionIdx)
554 + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
555 }
556
557 // immiscible oil
558 const LhsEval Rs(0.0);
559 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
560
561 return referenceDensity(phaseIdx, regionIdx)*bo;
562 }
563
564 case gasPhaseIdx: {
566 // gas containing vaporized oil and vaporized water
567 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
568 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
569 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
570
571 return
572 bg*referenceDensity(gasPhaseIdx, regionIdx)
573 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
574 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
575 }
576 if (enableVaporizedOil()) {
577 // miscible gas
578 const LhsEval Rvw(0.0);
579 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
580 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
581
582 return
583 bg*referenceDensity(gasPhaseIdx, regionIdx)
584 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
585 }
586 if (enableVaporizedWater()) {
587 // gas containing vaporized water
588 const LhsEval Rv(0.0);
589 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
590 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
591
592 return
593 bg*referenceDensity(gasPhaseIdx, regionIdx)
594 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
595 }
596
597 // immiscible gas
598 const LhsEval Rv(0.0);
599 const LhsEval Rvw(0.0);
600 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
601 return bg*referenceDensity(phaseIdx, regionIdx);
602 }
603
604 case waterPhaseIdx:
606 // gas miscible in water
607 const LhsEval& Rsw =BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
608 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
609 return
610 bw*referenceDensity(waterPhaseIdx, regionIdx)
611 + Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
612 }
613 const LhsEval Rsw(0.0);
614 return
616 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
617 }
618
619 throw std::logic_error("Unhandled phase index " + std::to_string(phaseIdx));
620 }
621
631 template <class FluidState, class LhsEval = typename FluidState::Scalar>
632 STATIC_OR_DEVICE LhsEval saturatedDensity(const FluidState& fluidState,
633 unsigned phaseIdx,
634 unsigned regionIdx)
635 {
636 assert(phaseIdx <= numPhases);
637 assert(regionIdx <= numRegions());
638
639 const auto& p = fluidState.pressure(phaseIdx);
640 const auto& T = fluidState.temperature(phaseIdx);
641
642 switch (phaseIdx) {
643 case oilPhaseIdx: {
644 if (enableDissolvedGas()) {
645 // miscible oil
646 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, oilPhaseIdx, regionIdx);
647 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
648
649 return
650 bo*referenceDensity(oilPhaseIdx, regionIdx)
651 + Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
652 }
653
654 // immiscible oil
655 const LhsEval Rs(0.0);
656 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
657 return referenceDensity(phaseIdx, regionIdx)*bo;
658 }
659
660 case gasPhaseIdx: {
662 // gas containing vaporized oil and vaporized water
663 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
664 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
665 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
666
667 return
668 bg*referenceDensity(gasPhaseIdx, regionIdx)
669 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
670 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx) ;
671 }
672
673 if (enableVaporizedOil()) {
674 // miscible gas
675 const LhsEval Rvw(0.0);
676 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
677 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
678
679 return
680 bg*referenceDensity(gasPhaseIdx, regionIdx)
681 + Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
682 }
683
684 if (enableVaporizedWater()) {
685 // gas containing vaporized water
686 const LhsEval Rv(0.0);
687 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState, gasPhaseIdx, regionIdx);
688 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
689
690 return
691 bg*referenceDensity(gasPhaseIdx, regionIdx)
692 + Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
693 }
694
695 // immiscible gas
696 const LhsEval Rv(0.0);
697 const LhsEval Rvw(0.0);
698 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
699
700 return referenceDensity(phaseIdx, regionIdx)*bg;
701
702 }
703
704 case waterPhaseIdx:
705 {
707 // miscible in water
708 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
709 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
710 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
711 return
712 bw*referenceDensity(waterPhaseIdx, regionIdx)
713 + Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
714 }
715 return
717 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
718 }
719 }
720
721 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
722 }
723
732 template <class FluidState, class LhsEval = typename FluidState::Scalar>
733 STATIC_OR_DEVICE LhsEval inverseFormationVolumeFactor(const FluidState& fluidState,
734 unsigned phaseIdx,
735 unsigned regionIdx)
736 {
737 OPM_TIMEBLOCK_LOCAL(inverseFormationVolumeFactor);
738 assert(phaseIdx <= numPhases);
739 assert(regionIdx <= numRegions());
740
741 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
742 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
743
744 switch (phaseIdx) {
745 case oilPhaseIdx: {
746 if (enableDissolvedGas()) {
747 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
748 if (useSaturatedTables() && fluidState.saturation(gasPhaseIdx) > 0.0
749 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
750 {
751 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
752 } else {
753 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
754 }
755 }
756
757 const LhsEval Rs(0.0);
758 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
759 }
760 case gasPhaseIdx: {
762 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
763 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
764 if (useSaturatedTables() && fluidState.saturation(waterPhaseIdx) > 0.0
765 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
766 && fluidState.saturation(oilPhaseIdx) > 0.0
767 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
768 {
769 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
770 } else {
771 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
772 }
773 }
774
775 if (enableVaporizedOil()) {
776 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
777 if (useSaturatedTables() && fluidState.saturation(oilPhaseIdx) > 0.0
778 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
779 {
780 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
781 } else {
782 const LhsEval Rvw(0.0);
783 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
784 }
785 }
786
787 if (enableVaporizedWater()) {
788 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
789 if (useSaturatedTables() && fluidState.saturation(waterPhaseIdx) > 0.0
790 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
791 {
792 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
793 } else {
794 const LhsEval Rv(0.0);
795 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
796 }
797 }
798
799 const LhsEval Rv(0.0);
800 const LhsEval Rvw(0.0);
801 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
802 }
803 case waterPhaseIdx:
804 {
805 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
807 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
808 if (useSaturatedTables() && fluidState.saturation(gasPhaseIdx) > 0.0
809 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
810 {
811 return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
812 } else {
813 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
814 }
815 }
816 const LhsEval Rsw(0.0);
817 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
818 }
819 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
820 }
821 }
822
832 template <class FluidState, class LhsEval = typename FluidState::Scalar>
833 STATIC_OR_DEVICE LhsEval saturatedInverseFormationVolumeFactor(const FluidState& fluidState,
834 unsigned phaseIdx,
835 unsigned regionIdx)
836 {
837 OPM_TIMEBLOCK_LOCAL(saturatedInverseFormationVolumeFactor);
838 assert(phaseIdx <= numPhases);
839 assert(regionIdx <= numRegions());
840
841 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
842 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
843 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
844
845 switch (phaseIdx) {
846 case oilPhaseIdx: return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
847 case gasPhaseIdx: return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
848 case waterPhaseIdx: return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
849 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
850 }
851 }
852
854 template <class FluidState, class LhsEval = typename FluidState::Scalar>
855 STATIC_OR_DEVICE LhsEval fugacityCoefficient(const FluidState& fluidState,
856 unsigned phaseIdx,
857 unsigned compIdx,
858 unsigned regionIdx)
859 {
860 assert(phaseIdx <= numPhases);
861 assert(compIdx <= numComponents);
862 assert(regionIdx <= numRegions());
863
864 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
865 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
866
867 // for the fugacity coefficient of the oil component in the oil phase, we use
868 // some pseudo-realistic value for the vapor pressure to ease physical
869 // interpretation of the results
870 const LhsEval phi_oO = 20e3/p;
871
872 // for the gas component in the gas phase, assume it to be an ideal gas
873 constexpr const Scalar phi_gG = 1.0;
874
875 // for the fugacity coefficient of the water component in the water phase, we use
876 // the same approach as for the oil component in the oil phase
877 const LhsEval phi_wW = 30e3/p;
878
879 switch (phaseIdx) {
880 case gasPhaseIdx: // fugacity coefficients for all components in the gas phase
881 switch (compIdx) {
882 case gasCompIdx:
883 return phi_gG;
884
885 // for the oil component, we calculate the Rv value for saturated gas and Rs
886 // for saturated oil, and compute the fugacity coefficients at the
887 // equilibrium. for this, we assume that in equilibrium the fugacities of the
888 // oil component is the same in both phases.
889 case oilCompIdx: {
890 if (!enableVaporizedOil())
891 // if there's no vaporized oil, the gas phase is assumed to be
892 // immiscible with the oil component
893 return phi_gG*1e6;
894
895 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
896 const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
897 const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
898
899 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
900 const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
901 const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
902 const auto& x_oOSat = 1.0 - x_oGSat;
903
904 const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
905 const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
906
907 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
908 }
909
910 case waterCompIdx:
911 // the water component is assumed to be never miscible with the gas phase
912 return phi_gG*1e6;
913
914 default:
915 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
916 }
917
918 case oilPhaseIdx: // fugacity coefficients for all components in the oil phase
919 switch (compIdx) {
920 case oilCompIdx:
921 return phi_oO;
922
923 // for the oil and water components, we proceed analogous to the gas and
924 // water components in the gas phase
925 case gasCompIdx: {
926 if (!enableDissolvedGas())
927 // if there's no dissolved gas, the oil phase is assumed to be
928 // immiscible with the gas component
929 return phi_oO*1e6;
930
931 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
932 const auto& X_gOSat = convertRvToXgO(R_vSat, regionIdx);
933 const auto& x_gOSat = convertXgOToxgO(X_gOSat, regionIdx);
934 const auto& x_gGSat = 1.0 - x_gOSat;
935
936 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
937 const auto& X_oGSat = convertRsToXoG(R_sSat, regionIdx);
938 const auto& x_oGSat = convertXoGToxoG(X_oGSat, regionIdx);
939
940 const auto& p_o = decay<LhsEval>(fluidState.pressure(oilPhaseIdx));
941 const auto& p_g = decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
942
943 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
944 }
945
946 case waterCompIdx:
947 return phi_oO*1e6;
948
949 default:
950 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
951 }
952
953 case waterPhaseIdx: // fugacity coefficients for all components in the water phase
954 // the water phase fugacity coefficients are pretty simple: because the water
955 // phase is assumed to consist entirely from the water component, we just
956 // need to make sure that the fugacity coefficients for the other components
957 // are a few orders of magnitude larger than the one of the water
958 // component. (i.e., the affinity of the gas and oil components to the water
959 // phase is lower by a few orders of magnitude)
960 switch (compIdx) {
961 case waterCompIdx: return phi_wW;
962 case oilCompIdx: return 1.1e6*phi_wW;
963 case gasCompIdx: return 1e6*phi_wW;
964 default:
965 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
966 }
967
968 default:
969 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
970 }
971
972 throw std::logic_error("Unhandled phase or component index");
973 }
974
976 template <class FluidState, class LhsEval = typename FluidState::Scalar>
977 STATIC_OR_DEVICE LhsEval viscosity(const FluidState& fluidState,
978 unsigned phaseIdx,
979 unsigned regionIdx)
980 {
981 OPM_TIMEBLOCK_LOCAL(viscosity);
982 assert(phaseIdx <= numPhases);
983 assert(regionIdx <= numRegions());
984
985 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
986 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
987
988 switch (phaseIdx) {
989 case oilPhaseIdx: {
990 if (enableDissolvedGas()) {
991 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
992 if (useSaturatedTables() && fluidState.saturation(gasPhaseIdx) > 0.0
993 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
994 {
995 return oilPvt_->saturatedViscosity(regionIdx, T, p);
996 } else {
997 return oilPvt_->viscosity(regionIdx, T, p, Rs);
998 }
999 }
1000
1001 const LhsEval Rs(0.0);
1002 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1003 }
1004
1005 case gasPhaseIdx: {
1007 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1008 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1009 if (useSaturatedTables() && fluidState.saturation(waterPhaseIdx) > 0.0
1010 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1011 && fluidState.saturation(oilPhaseIdx) > 0.0
1012 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1013 {
1014 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1015 } else {
1016 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1017 }
1018 }
1019 if (enableVaporizedOil()) {
1020 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1021 if (useSaturatedTables() && fluidState.saturation(oilPhaseIdx) > 0.0
1022 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1023 {
1024 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1025 } else {
1026 const LhsEval Rvw(0.0);
1027 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1028 }
1029 }
1030 if (enableVaporizedWater()) {
1031 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1032 if (useSaturatedTables() && fluidState.saturation(waterPhaseIdx) > 0.0
1033 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1034 {
1035 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1036 } else {
1037 const LhsEval Rv(0.0);
1038 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1039 }
1040 }
1041
1042 const LhsEval Rv(0.0);
1043 const LhsEval Rvw(0.0);
1044 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1045 }
1046
1047 case waterPhaseIdx:
1048 {
1049 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1051 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1052 if (useSaturatedTables() && fluidState.saturation(gasPhaseIdx) > 0.0
1053 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
1054 {
1055 return waterPvt_->saturatedViscosity(regionIdx, T, p, saltConcentration);
1056 } else {
1057 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1058 }
1059 }
1060 const LhsEval Rsw(0.0);
1061 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1062 }
1063 }
1064
1065 throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1066 }
1067
1068 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1069 STATIC_OR_DEVICE LhsEval internalEnergy(const FluidState& fluidState,
1070 const unsigned phaseIdx,
1071 const unsigned regionIdx)
1072 {
1073 const auto p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1074 const auto T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1075
1076 switch (phaseIdx) {
1077 case oilPhaseIdx:
1078 if (!oilPvt_->mixingEnergy()) {
1079 return oilPvt_->internalEnergy
1080 (regionIdx, T, p,
1081 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1082 }
1083 break;
1084
1085 case waterPhaseIdx:
1086 if (!waterPvt_->mixingEnergy()) {
1087 return waterPvt_->internalEnergy
1088 (regionIdx, T, p,
1089 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1090 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1091 }
1092 break;
1093
1094 case gasPhaseIdx:
1095 if (!gasPvt_->mixingEnergy()) {
1096 return gasPvt_->internalEnergy
1097 (regionIdx, T, p,
1098 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1099 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1100 }
1101 break;
1102
1103 default:
1104 throw std::logic_error {
1105 "Phase index " + std::to_string(phaseIdx) + " does not support internal energy"
1106 };
1107 }
1108
1109 return internalMixingTotalEnergy<FluidState,LhsEval>(fluidState, phaseIdx, regionIdx)
1110 / density<FluidState,LhsEval>(fluidState, phaseIdx, regionIdx);
1111 }
1112
1113
1114 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1115 STATIC_OR_DEVICE LhsEval internalMixingTotalEnergy(const FluidState& fluidState,
1116 unsigned phaseIdx,
1117 unsigned regionIdx)
1118 {
1119 assert(phaseIdx <= numPhases);
1120 assert(regionIdx <= numRegions());
1121 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1122 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1123 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1124 // to avoid putting all thermal into the interface of the multiplexer
1125 switch (phaseIdx) {
1126 case oilPhaseIdx: {
1127 auto oilEnergy = oilPvt_->internalEnergy(regionIdx, T, p,
1128 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1129 assert(oilPvt_->mixingEnergy());
1130 //mixing energy adsed
1131 if (enableDissolvedGas()) {
1132 // miscible oil
1133 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1134 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
1135 const auto& gasEnergy =
1136 gasPvt_->internalEnergy(regionIdx, T, p,
1137 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1138 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1139 const auto hVapG = gasPvt_->hVap(regionIdx);// pressure correction ? assume equal to energy change
1140 return
1141 oilEnergy*bo*referenceDensity(oilPhaseIdx, regionIdx)
1142 + (gasEnergy-hVapG)*Rs*bo*referenceDensity(gasPhaseIdx, regionIdx);
1143 }
1144
1145 // immiscible oil
1146 const LhsEval Rs(0.0);
1147 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
1148
1149 return oilEnergy*referenceDensity(phaseIdx, regionIdx)*bo;
1150 }
1151
1152 case gasPhaseIdx: {
1153 const auto& gasEnergy =
1154 gasPvt_->internalEnergy(regionIdx, T, p,
1155 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1156 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1157 assert(gasPvt_->mixingEnergy());
1159 const auto& oilEnergy =
1160 oilPvt_->internalEnergy(regionIdx, T, p,
1161 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1162 const auto waterEnergy =
1163 waterPvt_->internalEnergy(regionIdx, T, p,
1164 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1165 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1166 // gas containing vaporized oil and vaporized water
1167 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1168 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1169 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1170 const auto hVapO = oilPvt_->hVap(regionIdx);
1171 const auto hVapW = waterPvt_->hVap(regionIdx);
1172 return
1173 gasEnergy*bg*referenceDensity(gasPhaseIdx, regionIdx)
1174 + (oilEnergy+hVapO)*Rv*bg*referenceDensity(oilPhaseIdx, regionIdx)
1175 + (waterEnergy+hVapW)*Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
1176 }
1177 if (enableVaporizedOil()) {
1178 const auto& oilEnergy =
1179 oilPvt_->internalEnergy(regionIdx, T, p,
1180 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1181 // miscible gas
1182 const LhsEval Rvw(0.0);
1183 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1184 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1185 const auto hVapO = oilPvt_->hVap(regionIdx);
1186 return
1187 gasEnergy*bg*referenceDensity(gasPhaseIdx, regionIdx)
1188 + (oilEnergy+hVapO)*Rv*bg*referenceDensity(oilPhaseIdx, regionIdx);
1189 }
1190 if (enableVaporizedWater()) {
1191 // gas containing vaporized water
1192 const LhsEval Rv(0.0);
1193 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1194 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1195 const auto waterEnergy =
1196 waterPvt_->internalEnergy(regionIdx, T, p,
1197 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1198 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1199 const auto hVapW = waterPvt_->hVap(regionIdx);
1200 return
1201 gasEnergy*bg*referenceDensity(gasPhaseIdx, regionIdx)
1202 + (waterEnergy+hVapW)*Rvw*bg*referenceDensity(waterPhaseIdx, regionIdx);
1203 }
1204
1205 // immiscible gas
1206 const LhsEval Rv(0.0);
1207 const LhsEval Rvw(0.0);
1208 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1209 return gasEnergy*bg*referenceDensity(phaseIdx, regionIdx);
1210 }
1211
1212 case waterPhaseIdx:
1213 const auto waterEnergy =
1214 waterPvt_->internalEnergy(regionIdx, T, p,
1215 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1216 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1217 assert(waterPvt_->mixingEnergy());
1219 const auto& gasEnergy =
1220 gasPvt_->internalEnergy(regionIdx, T, p,
1221 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1222 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1223 // gas miscible in water
1224 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState, waterPhaseIdx, regionIdx);
1225 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
1226 return
1227 waterEnergy*bw*referenceDensity(waterPhaseIdx, regionIdx)
1228 + gasEnergy*Rsw*bw*referenceDensity(gasPhaseIdx, regionIdx);
1229 }
1230 const LhsEval Rsw(0.0);
1231 return
1232 waterEnergy*referenceDensity(waterPhaseIdx, regionIdx)
1233 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
1234 }
1235 throw std::logic_error("Unhandled phase index " + std::to_string(phaseIdx));
1236 }
1237
1238
1239
1241 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1242 STATIC_OR_DEVICE LhsEval enthalpy(const FluidState& fluidState,
1243 unsigned phaseIdx,
1244 unsigned regionIdx)
1245 {
1246 // should preferably not be used values should be taken from intensive quantities fluid state.
1247 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1248 auto energy = internalEnergy<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1249 if(!enthalpy_eq_energy_){
1250 // used for simplified models
1251 energy += p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1252 }
1253 return energy;
1254 }
1255
1262 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1263 STATIC_OR_DEVICE LhsEval saturatedVaporizationFactor(const FluidState& fluidState,
1264 unsigned phaseIdx,
1265 unsigned regionIdx)
1266 {
1267 assert(phaseIdx <= numPhases);
1268 assert(regionIdx <= numRegions());
1269
1270 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1271 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1272 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
1273
1274 switch (phaseIdx) {
1275 case oilPhaseIdx: return 0.0;
1276 case gasPhaseIdx: return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p, saltConcentration);
1277 case waterPhaseIdx: return 0.0;
1278 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1279 }
1280 }
1281
1288 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1289 STATIC_OR_DEVICE LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1290 unsigned phaseIdx,
1291 unsigned regionIdx,
1292 const LhsEval& maxOilSaturation)
1293 {
1294 OPM_TIMEBLOCK_LOCAL(saturatedDissolutionFactor);
1295 assert(phaseIdx <= numPhases);
1296 assert(regionIdx <= numRegions());
1297
1298 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1299 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1300 const auto& So = (phaseIdx == waterPhaseIdx) ? 0 : decay<LhsEval>(fluidState.saturation(oilPhaseIdx));
1301
1302 switch (phaseIdx) {
1303 case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1304 case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1305 case waterPhaseIdx: return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1306 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1307 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1308 }
1309 }
1310
1319 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1320 STATIC_OR_DEVICE LhsEval saturatedDissolutionFactor(const FluidState& fluidState,
1321 unsigned phaseIdx,
1322 unsigned regionIdx)
1323 {
1324 OPM_TIMEBLOCK_LOCAL(saturatedDissolutionFactor);
1325 assert(phaseIdx <= numPhases);
1326 assert(regionIdx <= numRegions());
1327
1328 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1329 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1330
1331 switch (phaseIdx) {
1332 case oilPhaseIdx: return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1333 case gasPhaseIdx: return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1334 case waterPhaseIdx: return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1335 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1336 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1337 }
1338 }
1339
1343 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1344 STATIC_OR_DEVICE LhsEval bubblePointPressure(const FluidState& fluidState,
1345 unsigned regionIdx)
1346 {
1347 return saturationPressure(fluidState, oilPhaseIdx, regionIdx);
1348 }
1349
1350
1354 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1355 STATIC_OR_DEVICE LhsEval dewPointPressure(const FluidState& fluidState,
1356 unsigned regionIdx)
1357 {
1358 return saturationPressure(fluidState, gasPhaseIdx, regionIdx);
1359 }
1360
1371 template <class FluidState, class LhsEval = typename FluidState::Scalar>
1372 STATIC_OR_DEVICE LhsEval saturationPressure(const FluidState& fluidState,
1373 unsigned phaseIdx,
1374 unsigned regionIdx)
1375 {
1376 assert(phaseIdx <= numPhases);
1377 assert(regionIdx <= numRegions());
1378
1379 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1380
1381 switch (phaseIdx) {
1382 case oilPhaseIdx: return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1383 case gasPhaseIdx: return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1384 case waterPhaseIdx: return waterPvt_->saturationPressure(regionIdx, T,
1385 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1386 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1387 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1388 }
1389 }
1390
1391 /****************************************
1392 * Auxiliary and convenience methods for the black-oil model
1393 ****************************************/
1398 template <class LhsEval>
1399 STATIC_OR_DEVICE LhsEval convertXoGToRs(const LhsEval& XoG, unsigned regionIdx)
1400 {
1401 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1402 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1403
1404 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1405 }
1406
1411 template <class LhsEval>
1412 STATIC_OR_DEVICE LhsEval convertXwGToRsw(const LhsEval& XwG, unsigned regionIdx)
1413 {
1414 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1415 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1416
1417 return XwG/(1.0 - XwG)*(rho_wRef/rho_gRef);
1418 }
1419
1424 template <class LhsEval>
1425 STATIC_OR_DEVICE LhsEval convertXgOToRv(const LhsEval& XgO, unsigned regionIdx)
1426 {
1427 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1428 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1429
1430 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1431 }
1432
1437 template <class LhsEval>
1438 STATIC_OR_DEVICE LhsEval convertXgWToRvw(const LhsEval& XgW, unsigned regionIdx)
1439 {
1440 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1441 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1442
1443 return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1444 }
1445
1446
1451 template <class LhsEval>
1452 STATIC_OR_DEVICE LhsEval convertRsToXoG(const LhsEval& Rs, unsigned regionIdx)
1453 {
1454 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1455 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1456
1457 const LhsEval& rho_oG = Rs*rho_gRef;
1458 return rho_oG/(rho_oRef + rho_oG);
1459 }
1460
1465 template <class LhsEval>
1466 STATIC_OR_DEVICE LhsEval convertRswToXwG(const LhsEval& Rsw, unsigned regionIdx)
1467 {
1468 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1469 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1470
1471 const LhsEval& rho_wG = Rsw*rho_gRef;
1472 return rho_wG/(rho_wRef + rho_wG);
1473 }
1474
1479 template <class LhsEval>
1480 STATIC_OR_DEVICE LhsEval convertRvToXgO(const LhsEval& Rv, unsigned regionIdx)
1481 {
1482 Scalar rho_oRef = referenceDensity_[regionIdx][oilPhaseIdx];
1483 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1484
1485 const LhsEval& rho_gO = Rv*rho_oRef;
1486 return rho_gO/(rho_gRef + rho_gO);
1487 }
1488
1493 template <class LhsEval>
1494 STATIC_OR_DEVICE LhsEval convertRvwToXgW(const LhsEval& Rvw, unsigned regionIdx)
1495 {
1496 Scalar rho_wRef = referenceDensity_[regionIdx][waterPhaseIdx];
1497 Scalar rho_gRef = referenceDensity_[regionIdx][gasPhaseIdx];
1498
1499 const LhsEval& rho_gW = Rvw*rho_wRef;
1500 return rho_gW/(rho_gRef + rho_gW);
1501 }
1502
1506 template <class LhsEval>
1507 STATIC_OR_DEVICE LhsEval convertXgWToxgW(const LhsEval& XgW, unsigned regionIdx)
1508 {
1509 Scalar MW = molarMass_[regionIdx][waterCompIdx];
1510 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1511
1512 return XgW*MG / (MW*(1 - XgW) + XgW*MG);
1513 }
1514
1518 template <class LhsEval>
1519 STATIC_OR_DEVICE LhsEval convertXwGToxwG(const LhsEval& XwG, unsigned regionIdx)
1520 {
1521 Scalar MW = molarMass_[regionIdx][waterCompIdx];
1522 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1523
1524 return XwG*MW / (MG*(1 - XwG) + XwG*MW);
1525 }
1526
1530 template <class LhsEval>
1531 STATIC_OR_DEVICE LhsEval convertXoGToxoG(const LhsEval& XoG, unsigned regionIdx)
1532 {
1533 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1534 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1535
1536 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1537 }
1538
1542 template <class LhsEval>
1543 STATIC_OR_DEVICE LhsEval convertxoGToXoG(const LhsEval& xoG, unsigned regionIdx)
1544 {
1545 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1546 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1547
1548 return xoG*MG / (xoG*(MG - MO) + MO);
1549 }
1550
1554 template <class LhsEval>
1555 STATIC_OR_DEVICE LhsEval convertXgOToxgO(const LhsEval& XgO, unsigned regionIdx)
1556 {
1557 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1558 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1559
1560 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1561 }
1562
1566 template <class LhsEval>
1567 STATIC_OR_DEVICE LhsEval convertxgOToXgO(const LhsEval& xgO, unsigned regionIdx)
1568 {
1569 Scalar MO = molarMass_[regionIdx][oilCompIdx];
1570 Scalar MG = molarMass_[regionIdx][gasCompIdx];
1571
1572 return xgO*MO / (xgO*(MO - MG) + MG);
1573 }
1574
1582 STATIC_OR_DEVICE const GasPvt& gasPvt()
1583 { return *gasPvt_; }
1584
1592 STATIC_OR_DEVICE const OilPvt& oilPvt()
1593 { return *oilPvt_; }
1594
1602 STATIC_OR_DEVICE const WaterPvt& waterPvt()
1603 { return *waterPvt_; }
1604
1610 STATIC_OR_DEVICE Scalar reservoirTemperature(unsigned = 0)
1611 { return reservoirTemperature_; }
1612
1618 STATIC_OR_DEVICE void setReservoirTemperature(Scalar value)
1619 { reservoirTemperature_ = value; }
1620
1621 STATIC_OR_DEVICE short activeToCanonicalPhaseIdx(unsigned activePhaseIdx);
1622
1623 STATIC_OR_DEVICE short canonicalToActivePhaseIdx(unsigned phaseIdx);
1624
1626 STATIC_OR_DEVICE Scalar diffusionCoefficient(unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1627 { return diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx]; }
1628
1630 STATIC_OR_DEVICE void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx = 0)
1631 { diffusionCoefficients_[regionIdx][numPhases*compIdx + phaseIdx] = coefficient ; }
1632
1636 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
1637 STATIC_OR_DEVICE LhsEval diffusionCoefficient(const FluidState& fluidState,
1638 const ParameterCache<ParamCacheEval>& paramCache,
1639 unsigned phaseIdx,
1640 unsigned compIdx)
1641 {
1642 // diffusion is disabled by the user
1643 if(!enableDiffusion())
1644 return 0.0;
1645
1646 // diffusion coefficients are set, and we use them
1647 if(!diffusionCoefficients_.empty()) {
1648 return diffusionCoefficient(compIdx, phaseIdx, paramCache.regionIndex());
1649 }
1650
1651 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1652 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1653
1654 switch (phaseIdx) {
1655 case oilPhaseIdx: return oilPvt().diffusionCoefficient(T, p, compIdx);
1656 case gasPhaseIdx: return gasPvt().diffusionCoefficient(T, p, compIdx);
1657 case waterPhaseIdx: return waterPvt().diffusionCoefficient(T, p, compIdx);
1658 default: throw std::logic_error("Unhandled phase index "+std::to_string(phaseIdx));
1659 }
1660 }
1661 STATIC_OR_DEVICE void setEnergyEqualEnthalpy(bool enthalpy_eq_energy){
1662 enthalpy_eq_energy_ = enthalpy_eq_energy;
1663 }
1664
1665 STATIC_OR_DEVICE bool enthalpyEqualEnergy(){
1666 return enthalpy_eq_energy_;
1667 }
1668
1669private:
1670 STATIC_OR_NOTHING void resizeArrays_(std::size_t numRegions);
1671
1672 STATIC_OR_NOTHING Scalar reservoirTemperature_;
1673
1674 STATIC_OR_NOTHING SmartPointer<GasPvt> gasPvt_;
1675 STATIC_OR_NOTHING SmartPointer<OilPvt> oilPvt_;
1676 STATIC_OR_NOTHING SmartPointer<WaterPvt> waterPvt_;
1677
1678 STATIC_OR_NOTHING bool enableDissolvedGas_;
1679 STATIC_OR_NOTHING bool enableDissolvedGasInWater_;
1680 STATIC_OR_NOTHING bool enableVaporizedOil_;
1681 STATIC_OR_NOTHING bool enableVaporizedWater_;
1682 STATIC_OR_NOTHING bool enableDiffusion_;
1683
1684 // HACK for GCC 4.4: the array size has to be specified using the literal value '3'
1685 // here, because GCC 4.4 seems to be unable to determine the number of phases from
1686 // the BlackOil fluid system in the attribute declaration below...
1687 STATIC_OR_NOTHING Storage<std::array<Scalar, /*numPhases=*/3> > referenceDensity_;
1688 STATIC_OR_NOTHING Storage<std::array<Scalar, /*numComponents=*/3> > molarMass_;
1689 STATIC_OR_NOTHING Storage<std::array<Scalar, /*numComponents=*/3 * /*numPhases=*/3> > diffusionCoefficients_;
1690
1691 STATIC_OR_NOTHING std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1692 STATIC_OR_NOTHING std::array<short, numPhases> canonicalToActivePhaseIdx_;
1693
1694 STATIC_OR_NOTHING bool isInitialized_;
1695 STATIC_OR_NOTHING bool useSaturatedTables_;
1696 STATIC_OR_NOTHING bool enthalpy_eq_energy_;
1697
1698 #ifndef COMPILING_STATIC_FLUID_SYSTEM
1699 template<template<typename> typename StorageT, template<typename> typename SmartPointerT>
1700 explicit FLUIDSYSTEM_CLASSNAME(const FLUIDSYSTEM_CLASSNAME_STATIC<Scalar, IndexTraits, StorageT, SmartPointerT>& other)
1703 , numActivePhases_(other.numActivePhases_)
1704 , phaseIsActive_(other.phaseIsActive_)
1705 , reservoirTemperature_(other.reservoirTemperature_)
1706 , gasPvt_(SmartPointerT<typename decltype(gasPvt_)::element_type>(other.gasPvt_))
1707 , oilPvt_(SmartPointerT<typename decltype(oilPvt_)::element_type>(other.oilPvt_))
1708 , waterPvt_(SmartPointerT<typename decltype(waterPvt_)::element_type>(other.waterPvt_))
1709 , enableDissolvedGas_(other.enableDissolvedGas_)
1710 , enableDissolvedGasInWater_(other.enableDissolvedGasInWater_)
1711 , enableVaporizedOil_(other.enableVaporizedOil_)
1712 , enableVaporizedWater_(other.enableVaporizedWater_)
1713 , enableDiffusion_(other.enableDiffusion_)
1714 , referenceDensity_(StorageT<typename decltype(referenceDensity_)::value_type>(other.referenceDensity_))
1715 , molarMass_(StorageT<typename decltype(molarMass_)::value_type>(other.molarMass_))
1716 , diffusionCoefficients_(StorageT<typename decltype(diffusionCoefficients_)::value_type>(other.diffusionCoefficients_))
1717 , activeToCanonicalPhaseIdx_(other.activeToCanonicalPhaseIdx_)
1718 , canonicalToActivePhaseIdx_(other.canonicalToActivePhaseIdx_)
1719 , isInitialized_(other.isInitialized_)
1720 , useSaturatedTables_(other.useSaturatedTables_)
1721 , enthalpy_eq_energy_(other.enthalpy_eq_energy_)
1722 {
1723 OPM_ERROR_IF(!other.isInitialized(), "The fluid system must be initialized before it can be copied.");
1724 }
1725
1726 template<class ScalarT, class IndexTraitsT, template<typename> typename StorageT, template<typename> typename SmartPointerT>
1727 friend class FLUIDSYSTEM_CLASSNAME_STATIC;
1728 #else
1729 template<class ScalarT, class IndexTraitsT, template<typename> typename StorageT, template<typename> typename SmartPointerT>
1730 friend class FLUIDSYSTEM_CLASSNAME_NONSTATIC;
1731 #endif
1732};
1733
1734template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1736initBegin(std::size_t numPvtRegions)
1737{
1738 isInitialized_ = false;
1739 useSaturatedTables_ = true;
1740
1741 enableDissolvedGas_ = true;
1742 enableDissolvedGasInWater_ = false;
1743 enableVaporizedOil_ = false;
1744 enableVaporizedWater_ = false;
1745 enableDiffusion_ = false;
1746
1747 oilPvt_ = nullptr;
1748 gasPvt_ = nullptr;
1749 waterPvt_ = nullptr;
1750
1751 surfaceTemperature = 273.15 + 15.56; // [K]
1752 surfacePressure = 1.01325e5; // [Pa]
1753 setReservoirTemperature(surfaceTemperature);
1754
1755 numActivePhases_ = numPhases;
1756 std::fill_n(&phaseIsActive_[0], numPhases, true);
1757
1758 resizeArrays_(numPvtRegions);
1759}
1760
1761template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1764 Scalar rhoWater,
1765 Scalar rhoGas,
1766 unsigned regionIdx)
1767{
1768 referenceDensity_[regionIdx][oilPhaseIdx] = rhoOil;
1769 referenceDensity_[regionIdx][waterPhaseIdx] = rhoWater;
1770 referenceDensity_[regionIdx][gasPhaseIdx] = rhoGas;
1771}
1772
1773template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1775{
1776 // calculate the final 2D functions which are used for interpolation.
1777 const std::size_t num_regions = molarMass_.size();
1778 for (unsigned regionIdx = 0; regionIdx < num_regions; ++regionIdx) {
1779 // calculate molar masses
1780
1781 // water is simple: 18 g/mol
1782 molarMass_[regionIdx][waterCompIdx] = 18e-3;
1783
1784 if (phaseIsActive(gasPhaseIdx)) {
1785 // for gas, we take the density at standard conditions and assume it to be ideal
1786 Scalar p = surfacePressure;
1787 Scalar T = surfaceTemperature;
1788 Scalar rho_g = referenceDensity_[/*regionIdx=*/0][gasPhaseIdx];
1789 molarMass_[regionIdx][gasCompIdx] = Constants<Scalar>::R*T*rho_g / p;
1790 }
1791 else
1792 // hydrogen gas. we just set this do avoid NaNs later
1793 molarMass_[regionIdx][gasCompIdx] = 2e-3;
1794
1795 // finally, for oil phase, we take the molar mass from the spe9 paper
1796 molarMass_[regionIdx][oilCompIdx] = 175e-3; // kg/mol
1797 }
1798
1799
1800 int activePhaseIdx = 0;
1801 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1802 if(phaseIsActive(phaseIdx)){
1803 canonicalToActivePhaseIdx_[phaseIdx] = activePhaseIdx;
1804 activeToCanonicalPhaseIdx_[activePhaseIdx] = phaseIdx;
1805 activePhaseIdx++;
1806 }
1807 }
1808 isInitialized_ = true;
1809}
1810
1811template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1813phaseName(unsigned phaseIdx)
1814{
1815 switch (phaseIdx) {
1816 case waterPhaseIdx:
1817 return "water";
1818 case oilPhaseIdx:
1819 return "oil";
1820 case gasPhaseIdx:
1821 return "gas";
1822
1823 default:
1824 throw std::logic_error(std::string("Phase index ") + std::to_string(phaseIdx) + " is unknown");
1825 }
1826}
1827
1828template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1830solventComponentIndex(unsigned phaseIdx)
1831{
1832 switch (phaseIdx) {
1833 case waterPhaseIdx:
1834 return waterCompIdx;
1835 case oilPhaseIdx:
1836 return oilCompIdx;
1837 case gasPhaseIdx:
1838 return gasCompIdx;
1839
1840 default:
1841 throw std::logic_error(std::string("Phase index ") + std::to_string(phaseIdx) + " is unknown");
1842 }
1843}
1844
1845template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1847soluteComponentIndex(unsigned phaseIdx)
1848{
1849 switch (phaseIdx) {
1850 case waterPhaseIdx:
1851 if (enableDissolvedGasInWater())
1852 return gasCompIdx;
1853 throw std::logic_error("The water phase does not have any solutes in the black oil model!");
1854 case oilPhaseIdx:
1855 return gasCompIdx;
1856 case gasPhaseIdx:
1857 if (enableVaporizedWater()) {
1858 return waterCompIdx;
1859 }
1860 return oilCompIdx;
1861
1862 default:
1863 throw std::logic_error(std::string("Phase index ") + std::to_string(phaseIdx) + " is unknown");
1864 }
1865}
1866
1867template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1869componentName(unsigned compIdx)
1870{
1871 switch (compIdx) {
1872 case waterCompIdx:
1873 return "Water";
1874 case oilCompIdx:
1875 return "Oil";
1876 case gasCompIdx:
1877 return "Gas";
1878
1879 default:
1880 throw std::logic_error(std::string("Component index ") + std::to_string(compIdx) + " is unknown");
1881 }
1882}
1883
1884template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1886activeToCanonicalPhaseIdx(unsigned activePhaseIdx)
1887{
1888 assert(activePhaseIdx<numActivePhases());
1889 return activeToCanonicalPhaseIdx_[activePhaseIdx];
1890}
1891
1892template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1893short FLUIDSYSTEM_CLASSNAME<Scalar,IndexTraits, Storage, SmartPointer>::
1894canonicalToActivePhaseIdx(unsigned phaseIdx)
1895{
1896 assert(phaseIdx<numPhases);
1897 assert(phaseIsActive(phaseIdx));
1898 return canonicalToActivePhaseIdx_[phaseIdx];
1899}
1900
1901template <class Scalar, class IndexTraits, template<typename> typename Storage, template<typename> typename SmartPointer>
1902void FLUIDSYSTEM_CLASSNAME<Scalar,IndexTraits,Storage, SmartPointer>::
1903resizeArrays_(std::size_t numRegions)
1904{
1905 molarMass_.resize(numRegions);
1906 referenceDensity_.resize(numRegions);
1907}
1908
1909#ifdef COMPILING_STATIC_FLUID_SYSTEM
1910template <typename T> using BOFS = FLUIDSYSTEM_CLASSNAME<T, BlackOilDefaultIndexTraits, VectorWithDefaultAllocator, std::shared_ptr>;
1911
1912#define DECLARE_INSTANCE(T) \
1913template<> unsigned char BOFS<T>::numActivePhases_; \
1914template<> std::array<bool, BOFS<T>::numPhases> BOFS<T>::phaseIsActive_; \
1915template<> std::array<short, BOFS<T>::numPhases> BOFS<T>::activeToCanonicalPhaseIdx_; \
1916template<> std::array<short, BOFS<T>::numPhases> BOFS<T>::canonicalToActivePhaseIdx_; \
1917template<> T BOFS<T>::surfaceTemperature; \
1918template<> T BOFS<T>::surfacePressure; \
1919template<> T BOFS<T>::reservoirTemperature_; \
1920template<> bool BOFS<T>::enableDissolvedGas_; \
1921template<> bool BOFS<T>::enableDissolvedGasInWater_; \
1922template<> bool BOFS<T>::enableVaporizedOil_; \
1923template<> bool BOFS<T>::enableVaporizedWater_; \
1924template<> bool BOFS<T>::enableDiffusion_; \
1925template<> std::shared_ptr<OilPvtMultiplexer<T>> BOFS<T>::oilPvt_; \
1926template<> std::shared_ptr<GasPvtMultiplexer<T>> BOFS<T>::gasPvt_; \
1927template<> std::shared_ptr<WaterPvtMultiplexer<T>> BOFS<T>::waterPvt_; \
1928template<> std::vector<std::array<T, 3>> BOFS<T>::referenceDensity_; \
1929template<> std::vector<std::array<T, 3>> BOFS<T>::molarMass_; \
1930template<> std::vector<std::array<T, 9>> BOFS<T>::diffusionCoefficients_; \
1931template<> bool BOFS<T>::isInitialized_; \
1932template<> bool BOFS<T>::useSaturatedTables_; \
1933template<> bool BOFS<T>::enthalpy_eq_energy_;
1934
1935DECLARE_INSTANCE(float)
1936DECLARE_INSTANCE(double)
1937
1938#undef DECLARE_INSTANCE
1939#endif
1940
1941} // namespace Opm
The base class for all fluid systems.
The class which specifies the default phase and component indices for the black-oil fluid system.
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
This macro generates a class HasMember_${MEMBER_NAME} which can be used for template specialization.
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
Definition Constants.hpp:40
Definition BlackOilFluidSystem.hpp:68
Definition BlackOilFluidSystemNonStatic.hpp:70
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
Definition BlackOilFluidSystem_macrotemplate.hpp:73
STATIC_OR_DEVICE 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_macrotemplate.hpp:1242
STATIC_OR_DEVICE void setEnableVaporizedOil(bool yesno)
Specify whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition BlackOilFluidSystem_macrotemplate.hpp:236
STATIC_OR_DEVICE LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1567
STATIC_OR_DEVICE bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition BlackOilFluidSystem_macrotemplate.hpp:405
STATIC_OR_DEVICE LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition BlackOilFluidSystem_macrotemplate.hpp:1344
STATIC_OR_DEVICE std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition BlackOilFluidSystem_macrotemplate.hpp:1869
STATIC_OR_DEVICE unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition BlackOilFluidSystem_macrotemplate.hpp:373
static constexpr int waterCompIdx
Index of the water component.
Definition BlackOilFluidSystem_macrotemplate.hpp:363
STATIC_OR_DEVICE 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_macrotemplate.hpp:1480
STATIC_OR_DEVICE unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition BlackOilFluidSystem_macrotemplate.hpp:1830
STATIC_OR_DEVICE 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_macrotemplate.hpp:517
STATIC_OR_DEVICE LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:534
STATIC_OR_DEVICE LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1531
STATIC_OR_DEVICE bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition BlackOilFluidSystem_macrotemplate.hpp:449
STATIC_OR_DEVICE bool phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition BlackOilFluidSystem_macrotemplate.hpp:377
STATIC_OR_DEVICE Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition BlackOilFluidSystem_macrotemplate.hpp:1610
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:335
STATIC_OR_DEVICE bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition BlackOilFluidSystem_macrotemplate.hpp:430
STATIC_OR_DEVICE LhsEval saturationPressure(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the saturation pressure of a given phase [Pa] depending on its composition.
Definition BlackOilFluidSystem_macrotemplate.hpp:1372
STATIC_OR_DEVICE 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_macrotemplate.hpp:1637
STATIC_OR_DEVICE void setEnableVaporizedWater(bool yesno)
Specify whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition BlackOilFluidSystem_macrotemplate.hpp:245
STATIC_OR_DEVICE LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1555
STATIC_OR_NOTHING Scalar surfacePressure
The pressure at the surface.
Definition BlackOilFluidSystem_macrotemplate.hpp:338
STATIC_OR_DEVICE bool useSaturatedTables()
Returns whether the saturated tables should be used.
Definition BlackOilFluidSystem_macrotemplate.hpp:474
STATIC_OR_DEVICE LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition BlackOilFluidSystem_macrotemplate.hpp:510
STATIC_OR_DEVICE bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition BlackOilFluidSystem_macrotemplate.hpp:347
STATIC_OR_DEVICE 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_macrotemplate.hpp:1452
STATIC_OR_DEVICE LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:833
STATIC_OR_DEVICE LhsEval convertXgWToxgW(const LhsEval &XgW, unsigned regionIdx)
Convert a water mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1507
STATIC_OR_DEVICE std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:1813
STATIC_OR_DEVICE 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_macrotemplate.hpp:855
STATIC_OR_DEVICE bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition BlackOilFluidSystem_macrotemplate.hpp:397
STATIC_OR_DEVICE LhsEval convertXwGToxwG(const LhsEval &XwG, unsigned regionIdx)
Convert a gas mass fraction in the water phase the corresponding mole fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1519
STATIC_OR_DEVICE void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition BlackOilFluidSystem_macrotemplate.hpp:1630
STATIC_OR_DEVICE 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_macrotemplate.hpp:1412
static constexpr int gasCompIdx
Index of the gas component.
Definition BlackOilFluidSystem_macrotemplate.hpp:365
STATIC_OR_DEVICE void setUseSaturatedTables(bool yesno)
Specify whether the saturated tables should be used.
Definition BlackOilFluidSystem_macrotemplate.hpp:269
STATIC_OR_DEVICE LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition BlackOilFluidSystem_macrotemplate.hpp:1355
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition BlackOilFluidSystem_macrotemplate.hpp:358
STATIC_OR_NOTHING Scalar surfaceTemperature
The temperature at the surface.
Definition BlackOilFluidSystem_macrotemplate.hpp:341
STATIC_OR_DEVICE 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_macrotemplate.hpp:1466
STATIC_OR_DEVICE unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition BlackOilFluidSystem_macrotemplate.hpp:1847
STATIC_OR_DEVICE 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_macrotemplate.hpp:1626
STATIC_OR_DEVICE bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition BlackOilFluidSystem_macrotemplate.hpp:458
STATIC_OR_DEVICE Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition BlackOilFluidSystem_macrotemplate.hpp:482
STATIC_OR_DEVICE void initEnd()
Finish initializing the black oil fluid system.
Definition BlackOilFluidSystem_macrotemplate.hpp:1774
STATIC_OR_DEVICE void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem_macrotemplate.hpp:261
STATIC_OR_DEVICE LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition BlackOilFluidSystem_macrotemplate.hpp:977
STATIC_OR_DEVICE LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1543
STATIC_OR_DEVICE bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition BlackOilFluidSystem_macrotemplate.hpp:409
STATIC_OR_DEVICE LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:1263
STATIC_OR_DEVICE std::size_t numRegions()
Returns the number of PVT regions which are considered.
Definition BlackOilFluidSystem_macrotemplate.hpp:421
static constexpr unsigned numPhases
Number of fluid phases in the fluid system.
Definition BlackOilFluidSystem_macrotemplate.hpp:328
STATIC_OR_DEVICE LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:490
static constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:331
STATIC_OR_DEVICE bool enableDissolvedGasInWater()
Returns whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem_macrotemplate.hpp:440
static constexpr unsigned oilPhaseIdx
Index of the oil phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:333
STATIC_OR_DEVICE 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_macrotemplate.hpp:1399
STATIC_OR_DEVICE LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx, const LhsEval &maxOilSaturation)
Returns the dissolution factor of a saturated fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:1289
FLUIDSYSTEM_CLASSNAME(const FLUIDSYSTEM_CLASSNAME< Scalar, IndexTraits, StorageT, SmartPointerT > &other)
Default copy constructor.
Definition BlackOilFluidSystem_macrotemplate.hpp:158
STATIC_OR_DEVICE void setEnableDissolvedGasInWater(bool yesno)
Specify whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem_macrotemplate.hpp:254
STATIC_OR_DEVICE 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_macrotemplate.hpp:1425
STATIC_OR_DEVICE void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition BlackOilFluidSystem_macrotemplate.hpp:1763
STATIC_OR_DEVICE LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:632
STATIC_OR_DEVICE const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition BlackOilFluidSystem_macrotemplate.hpp:1582
STATIC_OR_DEVICE LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:1320
STATIC_OR_DEVICE const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition BlackOilFluidSystem_macrotemplate.hpp:1592
STATIC_OR_DEVICE 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_macrotemplate.hpp:1494
STATIC_OR_DEVICE const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition BlackOilFluidSystem_macrotemplate.hpp:1602
STATIC_OR_DEVICE 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_macrotemplate.hpp:497
STATIC_OR_DEVICE void initBegin(std::size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition BlackOilFluidSystem_macrotemplate.hpp:1736
STATIC_OR_DEVICE void setWaterPvt(SmartPointer< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:287
STATIC_OR_DEVICE 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_macrotemplate.hpp:1438
static constexpr int oilCompIdx
Index of the oil component.
Definition BlackOilFluidSystem_macrotemplate.hpp:361
STATIC_OR_DEVICE bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem_macrotemplate.hpp:466
STATIC_OR_DEVICE void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition BlackOilFluidSystem_macrotemplate.hpp:1618
STATIC_OR_DEVICE void setEnableDissolvedGas(bool yesno)
Specify whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition BlackOilFluidSystem_macrotemplate.hpp:227
STATIC_OR_DEVICE void setGasPvt(SmartPointer< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:275
STATIC_OR_DEVICE LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:733
STATIC_OR_DEVICE Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition BlackOilFluidSystem_macrotemplate.hpp:393
STATIC_OR_DEVICE void setOilPvt(SmartPointer< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:281
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition GasPvtMultiplexer.hpp:111
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:269
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