My Project
Loading...
Searching...
No Matches
BlackOilFluidState.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*/
28#ifndef OPM_BLACK_OIL_FLUID_STATE_HH
29#define OPM_BLACK_OIL_FLUID_STATE_HH
30
33
36
37namespace Opm {
38
39OPM_GENERATE_HAS_MEMBER(pvtRegionIndex, ) // Creates 'HasMember_pvtRegionIndex<T>'.
40
41template <class FluidState>
42unsigned getPvtRegionIndex_(typename std::enable_if<HasMember_pvtRegionIndex<FluidState>::value,
43 const FluidState&>::type fluidState)
44{ return fluidState.pvtRegionIndex(); }
45
46template <class FluidState>
47unsigned getPvtRegionIndex_(typename std::enable_if<!HasMember_pvtRegionIndex<FluidState>::value,
48 const FluidState&>::type)
49{ return 0; }
50
51OPM_GENERATE_HAS_MEMBER(invB, /*phaseIdx=*/0) // Creates 'HasMember_invB<T>'.
52
53template <class FluidSystem, class FluidState, class LhsEval>
54auto getInvB_(typename std::enable_if<HasMember_invB<FluidState>::value,
55 const FluidState&>::type fluidState,
56 unsigned phaseIdx,
57 unsigned)
58 -> decltype(decay<LhsEval>(fluidState.invB(phaseIdx)))
59{ return decay<LhsEval>(fluidState.invB(phaseIdx)); }
60
61template <class FluidSystem, class FluidState, class LhsEval>
62LhsEval getInvB_(typename std::enable_if<!HasMember_invB<FluidState>::value,
63 const FluidState&>::type fluidState,
64 unsigned phaseIdx,
65 unsigned pvtRegionIdx)
66{
67 const auto& rho = fluidState.density(phaseIdx);
68 const auto& Xsolvent =
69 fluidState.massFraction(phaseIdx, FluidSystem::solventComponentIndex(phaseIdx));
70
71 return
72 decay<LhsEval>(rho)
73 *decay<LhsEval>(Xsolvent)
74 /FluidSystem::referenceDensity(phaseIdx, pvtRegionIdx);
75}
76
77OPM_GENERATE_HAS_MEMBER(saltConcentration, ) // Creates 'HasMember_saltConcentration<T>'.
78
79template <class FluidState>
80auto getSaltConcentration_(typename std::enable_if<HasMember_saltConcentration<FluidState>::value,
81 const FluidState&>::type fluidState)
82{ return fluidState.saltConcentration(); }
83
84template <class FluidState>
85auto getSaltConcentration_(typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
86 const FluidState&>::type)
87{ return 0.0; }
88
89OPM_GENERATE_HAS_MEMBER(saltSaturation, ) // Creates 'HasMember_saltSaturation<T>'.
90
91template <class FluidState>
92auto getSaltSaturation_(typename std::enable_if<HasMember_saltSaturation<FluidState>::value,
93 const FluidState&>::type fluidState)
94{ return fluidState.saltSaturation(); }
95
96
97template <class FluidState>
98auto getSaltSaturation_(typename std::enable_if<!HasMember_saltSaturation<FluidState>::value,
99 const FluidState&>::type)
100{ return 0.0; }
101
109template <class ScalarT,
110 class FluidSystem,
111 bool enableTemperature = false,
112 bool enableEnergy = false,
113 bool enableDissolution = true,
114 bool enableVapwat = false,
115 bool enableBrine = false,
116 bool enableSaltPrecipitation = false,
117 bool enableDissolutionInWater = false,
118 unsigned numStoragePhases = FluidSystem::numPhases>
120{
121 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
122 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
123 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
124
125 enum { waterCompIdx = FluidSystem::waterCompIdx };
126 enum { gasCompIdx = FluidSystem::gasCompIdx };
127 enum { oilCompIdx = FluidSystem::oilCompIdx };
128
129public:
130 using Scalar = ScalarT;
131 enum { numPhases = FluidSystem::numPhases };
132 enum { numComponents = FluidSystem::numComponents };
133
142 void checkDefined() const
143 {
144#ifndef NDEBUG
145 Valgrind::CheckDefined(pvtRegionIdx_);
146
147 for (unsigned storagePhaseIdx = 0; storagePhaseIdx < numStoragePhases; ++ storagePhaseIdx) {
148 Valgrind::CheckDefined(saturation_[storagePhaseIdx]);
149 Valgrind::CheckDefined(pressure_[storagePhaseIdx]);
150 Valgrind::CheckDefined(density_[storagePhaseIdx]);
151 Valgrind::CheckDefined(invB_[storagePhaseIdx]);
152
153 if constexpr (enableEnergy)
154 Valgrind::CheckDefined((*enthalpy_)[storagePhaseIdx]);
155 }
156
157 if constexpr (enableDissolution) {
158 Valgrind::CheckDefined(*Rs_);
159 Valgrind::CheckDefined(*Rv_);
160 }
161
162 if constexpr (enableVapwat) {
163 Valgrind::CheckDefined(*Rvw_);
164 }
165
166 if constexpr (enableDissolutionInWater) {
167 Valgrind::CheckDefined(*Rsw_);
168 }
169
170 if constexpr (enableBrine) {
171 Valgrind::CheckDefined(*saltConcentration_);
172 }
173
174 if constexpr (enableSaltPrecipitation) {
175 Valgrind::CheckDefined(*saltSaturation_);
176 }
177
178 if constexpr (enableTemperature || enableEnergy)
179 Valgrind::CheckDefined(*temperature_);
180#endif // NDEBUG
181 }
182
187 template <class FluidState>
188 void assign(const FluidState& fs)
189 {
190 if constexpr (enableTemperature || enableEnergy)
191 setTemperature(fs.temperature(/*phaseIdx=*/0));
192
193 unsigned pvtRegionIdx = getPvtRegionIndex_<FluidState>(fs);
194 setPvtRegionIndex(pvtRegionIdx);
195
196 if constexpr (enableDissolution) {
197 setRs(BlackOil::getRs_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
198 setRv(BlackOil::getRv_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
199 }
200 if constexpr (enableVapwat) {
201 setRvw(BlackOil::getRvw_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
202 }
203 if constexpr (enableDissolutionInWater) {
204 setRsw(BlackOil::getRsw_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
205 }
206 if constexpr (enableBrine){
207 setSaltConcentration(BlackOil::getSaltConcentration_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
208 }
209 if constexpr (enableSaltPrecipitation){
210 setSaltSaturation(BlackOil::getSaltSaturation_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
211 }
212 for (unsigned storagePhaseIdx = 0; storagePhaseIdx < numStoragePhases; ++storagePhaseIdx) {
213 unsigned phaseIdx = storageToCanonicalPhaseIndex_(storagePhaseIdx);
214 setSaturation(phaseIdx, fs.saturation(phaseIdx));
215 setPressure(phaseIdx, fs.pressure(phaseIdx));
216 setDensity(phaseIdx, fs.density(phaseIdx));
217
218 if constexpr (enableEnergy)
219 setEnthalpy(phaseIdx, fs.enthalpy(phaseIdx));
220
221 setInvB(phaseIdx, getInvB_<FluidSystem, FluidState, Scalar>(fs, phaseIdx, pvtRegionIdx));
222 }
223 }
224
231 void setPvtRegionIndex(unsigned newPvtRegionIdx)
232 { pvtRegionIdx_ = static_cast<unsigned short>(newPvtRegionIdx); }
233
237 void setPressure(unsigned phaseIdx, const Scalar& p)
238 { pressure_[canonicalToStoragePhaseIndex_(phaseIdx)] = p; }
239
243 void setSaturation(unsigned phaseIdx, const Scalar& S)
244 { saturation_[canonicalToStoragePhaseIndex_(phaseIdx)] = S; }
245
249 void setPc(unsigned phaseIdx, const Scalar& pc)
250 { pc_[canonicalToStoragePhaseIndex_(phaseIdx)] = pc; }
251
255 void setTotalSaturation(const Scalar& value)
256 {
257 totalSaturation_ = value;
258 }
259
266 void setTemperature(const Scalar& value)
267 {
268 assert(enableTemperature || enableEnergy);
269
270 (*temperature_) = value;
271 }
272
279 void setEnthalpy(unsigned phaseIdx, const Scalar& value)
280 {
281 assert(enableTemperature || enableEnergy);
282
283 (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)] = value;
284 }
285
289 void setInvB(unsigned phaseIdx, const Scalar& b)
290 { invB_[canonicalToStoragePhaseIndex_(phaseIdx)] = b; }
291
295 void setDensity(unsigned phaseIdx, const Scalar& rho)
296 { density_[canonicalToStoragePhaseIndex_(phaseIdx)] = rho; }
297
303 void setRs(const Scalar& newRs)
304 { *Rs_ = newRs; }
305
311 void setRv(const Scalar& newRv)
312 { *Rv_ = newRv; }
313
319 void setRvw(const Scalar& newRvw)
320 { *Rvw_ = newRvw; }
321
327 void setRsw(const Scalar& newRsw)
328 { *Rsw_ = newRsw; }
329
333 void setSaltConcentration(const Scalar& newSaltConcentration)
334 { *saltConcentration_ = newSaltConcentration; }
335
339 void setSaltSaturation(const Scalar& newSaltSaturation)
340 { *saltSaturation_ = newSaltSaturation; }
341
345 const Scalar& pressure(unsigned phaseIdx) const
346 { return pressure_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
347
351 const Scalar& saturation(unsigned phaseIdx) const
352 { return saturation_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
353
357 const Scalar& pc(unsigned phaseIdx) const
358 { return pc_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
359
363 const Scalar& totalSaturation() const
364 {
365 return totalSaturation_;
366 }
367
371 const Scalar& temperature(unsigned) const
372 {
373 if constexpr (enableTemperature || enableEnergy) {
374 return *temperature_;
375 } else {
376 static Scalar tmp(FluidSystem::reservoirTemperature(pvtRegionIdx_));
377 return tmp;
378 }
379 }
380
387 const Scalar& invB(unsigned phaseIdx) const
388 { return invB_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
389
397 const Scalar& Rs() const
398 {
399 if constexpr (enableDissolution) {
400 return *Rs_;
401 } else {
402 static Scalar null = 0.0;
403 return null;
404 }
405 }
406
414 const Scalar& Rv() const
415 {
416 if constexpr (!enableDissolution) {
417 static Scalar null = 0.0;
418 return null;
419 } else {
420 return *Rv_;
421 }
422 }
423
431 const Scalar& Rvw() const
432 {
433 if constexpr (enableVapwat) {
434 return *Rvw_;
435 } else {
436 static Scalar null = 0.0;
437 return null;
438 }
439 }
440
448 const Scalar& Rsw() const
449 {
450 if constexpr (enableDissolutionInWater) {
451 return *Rsw_;
452 } else {
453 static Scalar null = 0.0;
454 return null;
455 }
456 }
457
461 const Scalar& saltConcentration() const
462 {
463 if constexpr (enableBrine) {
464 return *saltConcentration_;
465 } else {
466 static Scalar null = 0.0;
467 return null;
468 }
469 }
470
474 const Scalar& saltSaturation() const
475 {
476 if constexpr (enableSaltPrecipitation) {
477 return *saltSaturation_;
478 } else {
479 static Scalar null = 0.0;
480 return null;
481 }
482 }
483
488 unsigned short pvtRegionIndex() const
489 { return pvtRegionIdx_; }
490
494 Scalar density(unsigned phaseIdx) const
495 { return density_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
496
503 const Scalar& enthalpy(unsigned phaseIdx) const
504 { return (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)]; }
505
512 Scalar internalEnergy(unsigned phaseIdx) const
513 { return (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)] - pressure(phaseIdx)/density(phaseIdx); }
514
516 // slow methods
518
522 Scalar molarDensity(unsigned phaseIdx) const
523 {
524 const auto& rho = density(phaseIdx);
525
526 if (phaseIdx == waterPhaseIdx)
527 return rho/FluidSystem::molarMass(waterCompIdx, pvtRegionIdx_);
528
529 return
530 rho*(moleFraction(phaseIdx, gasCompIdx)/FluidSystem::molarMass(gasCompIdx, pvtRegionIdx_)
531 + moleFraction(phaseIdx, oilCompIdx)/FluidSystem::molarMass(oilCompIdx, pvtRegionIdx_));
532
533 }
534
540 Scalar molarVolume(unsigned phaseIdx) const
541 { return 1.0/molarDensity(phaseIdx); }
542
546 Scalar viscosity(unsigned phaseIdx) const
547 { return FluidSystem::viscosity(*this, phaseIdx, pvtRegionIdx_); }
548
552 Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
553 {
554 switch (phaseIdx) {
555 case waterPhaseIdx:
556 if (compIdx == waterCompIdx)
557 return 1.0;
558 return 0.0;
559
560 case oilPhaseIdx:
561 if (compIdx == waterCompIdx)
562 return 0.0;
563 else if (compIdx == oilCompIdx)
564 return 1.0 - FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_);
565 else {
566 assert(compIdx == gasCompIdx);
567 return FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_);
568 }
569 break;
570
571 case gasPhaseIdx:
572 if (compIdx == waterCompIdx)
573 return 0.0;
574 else if (compIdx == oilCompIdx)
575 return FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_);
576 else {
577 assert(compIdx == gasCompIdx);
578 return 1.0 - FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_);
579 }
580 break;
581 }
582
583 throw std::logic_error("Invalid phase or component index!");
584 }
585
589 Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
590 {
591 switch (phaseIdx) {
592 case waterPhaseIdx:
593 if (compIdx == waterCompIdx)
594 return 1.0;
595 return 0.0;
596
597 case oilPhaseIdx:
598 if (compIdx == waterCompIdx)
599 return 0.0;
600 else if (compIdx == oilCompIdx)
601 return 1.0 - FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_),
602 pvtRegionIdx_);
603 else {
604 assert(compIdx == gasCompIdx);
605 return FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_),
606 pvtRegionIdx_);
607 }
608 break;
609
610 case gasPhaseIdx:
611 if (compIdx == waterCompIdx)
612 return 0.0;
613 else if (compIdx == oilCompIdx)
614 return FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_),
615 pvtRegionIdx_);
616 else {
617 assert(compIdx == gasCompIdx);
618 return 1.0 - FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_),
619 pvtRegionIdx_);
620 }
621 break;
622 }
623
624 throw std::logic_error("Invalid phase or component index!");
625 }
626
630 Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
631 { return moleFraction(phaseIdx, compIdx)*molarDensity(phaseIdx); }
632
636 Scalar averageMolarMass(unsigned phaseIdx) const
637 {
638 Scalar result(0.0);
639 for (unsigned compIdx = 0; compIdx < numComponents; ++ compIdx)
640 result += FluidSystem::molarMass(compIdx, pvtRegionIdx_)*moleFraction(phaseIdx, compIdx);
641 return result;
642 }
643
647 Scalar fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
648 { return FluidSystem::fugacityCoefficient(*this, phaseIdx, compIdx, pvtRegionIdx_); }
649
653 Scalar fugacity(unsigned phaseIdx, unsigned compIdx) const
654 {
655 return
656 fugacityCoefficient(phaseIdx, compIdx)
657 *moleFraction(phaseIdx, compIdx)
658 *pressure(phaseIdx);
659 }
660
661private:
662 static unsigned storageToCanonicalPhaseIndex_(unsigned storagePhaseIdx)
663 {
664 if constexpr (numStoragePhases == 3)
665 return storagePhaseIdx;
666 else
667 return FluidSystem::activeToCanonicalPhaseIdx(storagePhaseIdx);
668 }
669
670 static unsigned canonicalToStoragePhaseIndex_(unsigned canonicalPhaseIdx)
671 {
672 if constexpr (numStoragePhases == 3)
673 return canonicalPhaseIdx;
674 else
675 return FluidSystem::canonicalToActivePhaseIdx(canonicalPhaseIdx);
676 }
677
678 ConditionalStorage<enableTemperature || enableEnergy, Scalar> temperature_;
679 ConditionalStorage<enableEnergy, std::array<Scalar, numStoragePhases> > enthalpy_;
680 Scalar totalSaturation_;
681 std::array<Scalar, numStoragePhases> pressure_;
682 std::array<Scalar, numStoragePhases> pc_;
683 std::array<Scalar, numStoragePhases> saturation_;
684 std::array<Scalar, numStoragePhases> invB_;
685 std::array<Scalar, numStoragePhases> density_;
686 ConditionalStorage<enableDissolution,Scalar> Rs_;
687 ConditionalStorage<enableDissolution, Scalar> Rv_;
688 ConditionalStorage<enableVapwat,Scalar> Rvw_;
689 ConditionalStorage<enableDissolutionInWater,Scalar> Rsw_;
690 ConditionalStorage<enableBrine, Scalar> saltConcentration_;
691 ConditionalStorage<enableSaltPrecipitation, Scalar> saltSaturation_;
692 unsigned short pvtRegionIdx_;
693};
694
695} // namespace Opm
696
697#endif
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
A simple class which only stores a given member attribute if a boolean condition is true.
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
Some templates to wrap the valgrind client request macros.
Implements a "tailor-made" fluid state class for the black-oil model.
Definition BlackOilFluidState.hpp:120
void setEnthalpy(unsigned phaseIdx, const Scalar &value)
Set the specific enthalpy [J/kg] of a given fluid phase.
Definition BlackOilFluidState.hpp:279
const Scalar & pressure(unsigned phaseIdx) const
Return the pressure of a fluid phase [Pa].
Definition BlackOilFluidState.hpp:345
void setTemperature(const Scalar &value)
Set the temperature [K].
Definition BlackOilFluidState.hpp:266
const Scalar & invB(unsigned phaseIdx) const
Return the inverse formation volume factor of a fluid phase [-].
Definition BlackOilFluidState.hpp:387
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
Return the mole fraction of a component in a fluid phase [-].
Definition BlackOilFluidState.hpp:589
void setSaltSaturation(const Scalar &newSaltSaturation)
Set the solid salt saturation.
Definition BlackOilFluidState.hpp:339
void setSaltConcentration(const Scalar &newSaltConcentration)
Set the salt concentration.
Definition BlackOilFluidState.hpp:333
void setRs(const Scalar &newRs)
Set the gas dissolution factor [m^3/m^3] of the oil phase.
Definition BlackOilFluidState.hpp:303
const Scalar & Rv() const
Return the oil vaporization factor of gas [m^3/m^3].
Definition BlackOilFluidState.hpp:414
void setDensity(unsigned phaseIdx, const Scalar &rho)
\ brief Set the density of a fluid phase
Definition BlackOilFluidState.hpp:295
void checkDefined() const
Make sure that all attributes are defined.
Definition BlackOilFluidState.hpp:142
Scalar averageMolarMass(unsigned phaseIdx) const
Return the partial molar density of a fluid phase [kg / mol].
Definition BlackOilFluidState.hpp:636
void setInvB(unsigned phaseIdx, const Scalar &b)
\ brief Set the inverse formation volume factor of a fluid phase
Definition BlackOilFluidState.hpp:289
const Scalar & Rsw() const
Return the gas dissolution factor of water [m^3/m^3].
Definition BlackOilFluidState.hpp:448
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition BlackOilFluidState.hpp:188
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
Return the partial molar density of a component in a fluid phase [mol / m^3].
Definition BlackOilFluidState.hpp:630
Scalar internalEnergy(unsigned phaseIdx) const
Return the specific internal energy [J/kg] of a given fluid phase.
Definition BlackOilFluidState.hpp:512
void setSaturation(unsigned phaseIdx, const Scalar &S)
Set the saturation of a fluid phase [-].
Definition BlackOilFluidState.hpp:243
void setRv(const Scalar &newRv)
Set the oil vaporization factor [m^3/m^3] of the gas phase.
Definition BlackOilFluidState.hpp:311
void setPvtRegionIndex(unsigned newPvtRegionIdx)
Set the index of the fluid region.
Definition BlackOilFluidState.hpp:231
const Scalar & pc(unsigned phaseIdx) const
Return the capillary pressure of a fluid phase [-].
Definition BlackOilFluidState.hpp:357
Scalar fugacity(unsigned phaseIdx, unsigned compIdx) const
Return the fugacity of a component in a fluid phase [Pa].
Definition BlackOilFluidState.hpp:653
void setTotalSaturation(const Scalar &value)
Set the total saturation used for sequential methods.
Definition BlackOilFluidState.hpp:255
const Scalar & saltConcentration() const
Return the concentration of salt in water.
Definition BlackOilFluidState.hpp:461
const Scalar & saltSaturation() const
Return the saturation of solid salt.
Definition BlackOilFluidState.hpp:474
const Scalar & temperature(unsigned) const
Return the temperature [K].
Definition BlackOilFluidState.hpp:371
Scalar molarDensity(unsigned phaseIdx) const
Return the molar density of a fluid phase [mol/m^3].
Definition BlackOilFluidState.hpp:522
void setPressure(unsigned phaseIdx, const Scalar &p)
Set the pressure of a fluid phase [-].
Definition BlackOilFluidState.hpp:237
unsigned short pvtRegionIndex() const
Return the PVT region where the current fluid state is assumed to be part of.
Definition BlackOilFluidState.hpp:488
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
Return the mass fraction of a component in a fluid phase [-].
Definition BlackOilFluidState.hpp:552
void setPc(unsigned phaseIdx, const Scalar &pc)
Set the capillary pressure of a fluid phase [-].
Definition BlackOilFluidState.hpp:249
const Scalar & totalSaturation() const
Return the total saturation needed for sequential.
Definition BlackOilFluidState.hpp:363
const Scalar & saturation(unsigned phaseIdx) const
Return the saturation of a fluid phase [-].
Definition BlackOilFluidState.hpp:351
Scalar fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
Return the fugacity coefficient of a component in a fluid phase [-].
Definition BlackOilFluidState.hpp:647
const Scalar & enthalpy(unsigned phaseIdx) const
Return the specific enthalpy [J/kg] of a given fluid phase.
Definition BlackOilFluidState.hpp:503
Scalar molarVolume(unsigned phaseIdx) const
Return the molar volume of a fluid phase [m^3/mol].
Definition BlackOilFluidState.hpp:540
const Scalar & Rvw() const
Return the water vaporization factor of gas [m^3/m^3].
Definition BlackOilFluidState.hpp:431
const Scalar & Rs() const
Return the gas dissolution factor of oil [m^3/m^3].
Definition BlackOilFluidState.hpp:397
Scalar density(unsigned phaseIdx) const
Return the density [kg/m^3] of a given fluid phase.
Definition BlackOilFluidState.hpp:494
void setRsw(const Scalar &newRsw)
Set the gas dissolution factor [m^3/m^3] of the water phase.
Definition BlackOilFluidState.hpp:327
Scalar viscosity(unsigned phaseIdx) const
Return the dynamic viscosity of a fluid phase [Pa s].
Definition BlackOilFluidState.hpp:546
void setRvw(const Scalar &newRvw)
Set the water vaporization factor [m^3/m^3] of the gas phase.
Definition BlackOilFluidState.hpp:319
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30