My Project
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 
37 namespace Opm {
38 
39 OPM_GENERATE_HAS_MEMBER(pvtRegionIndex, ) // Creates 'HasMember_pvtRegionIndex<T>'.
40 
41 template <class FluidState>
42 unsigned getPvtRegionIndex_(typename std::enable_if<HasMember_pvtRegionIndex<FluidState>::value,
43  const FluidState&>::type fluidState)
44 { return fluidState.pvtRegionIndex(); }
45 
46 template <class FluidState>
47 unsigned getPvtRegionIndex_(typename std::enable_if<!HasMember_pvtRegionIndex<FluidState>::value,
48  const FluidState&>::type)
49 { return 0; }
50 
51 OPM_GENERATE_HAS_MEMBER(invB, /*phaseIdx=*/0) // Creates 'HasMember_invB<T>'.
52 
53 template <class FluidSystem, class FluidState, class LhsEval>
54 auto 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 
61 template <class FluidSystem, class FluidState, class LhsEval>
62 LhsEval 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 
77 OPM_GENERATE_HAS_MEMBER(saltConcentration, ) // Creates 'HasMember_saltConcentration<T>'.
78 
79 template <class FluidState>
80 auto getSaltConcentration_(typename std::enable_if<HasMember_saltConcentration<FluidState>::value,
81  const FluidState&>::type fluidState)
82 { return fluidState.saltConcentration(); }
83 
84 template <class FluidState>
85 auto getSaltConcentration_(typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
86  const FluidState&>::type)
87 { return 0.0; }
88 
89 OPM_GENERATE_HAS_MEMBER(saltSaturation, ) // Creates 'HasMember_saltSaturation<T>'.
90 
91 template <class FluidState>
92 auto getSaltSaturation_(typename std::enable_if<HasMember_saltSaturation<FluidState>::value,
93  const FluidState&>::type fluidState)
94 { return fluidState.saltSaturation(); }
95 
96 
97 template <class FluidState>
98 auto getSaltSaturation_(typename std::enable_if<!HasMember_saltSaturation<FluidState>::value,
99  const FluidState&>::type)
100 { return 0.0; }
101 
109 template <class ScalarT,
110  class FluidSystem,
111  bool enableTemperature = false,
112  bool enableEnergy = false,
113  bool enableDissolution = true,
114  bool enableEvaporation = false,
115  bool enableBrine = false,
116  bool enableSaltPrecipitation = false,
117  unsigned numStoragePhases = FluidSystem::numPhases>
119 {
120  enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
121  enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
122  enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
123 
124  enum { waterCompIdx = FluidSystem::waterCompIdx };
125  enum { gasCompIdx = FluidSystem::gasCompIdx };
126  enum { oilCompIdx = FluidSystem::oilCompIdx };
127 
128 public:
129  using Scalar = ScalarT;
130  enum { numPhases = FluidSystem::numPhases };
131  enum { numComponents = FluidSystem::numComponents };
132 
141  void checkDefined() const
142  {
143 #ifndef NDEBUG
144  Valgrind::CheckDefined(pvtRegionIdx_);
145 
146  for (unsigned storagePhaseIdx = 0; storagePhaseIdx < numStoragePhases; ++ storagePhaseIdx) {
147  Valgrind::CheckDefined(saturation_[storagePhaseIdx]);
148  Valgrind::CheckDefined(pressure_[storagePhaseIdx]);
149  Valgrind::CheckDefined(density_[storagePhaseIdx]);
150  Valgrind::CheckDefined(invB_[storagePhaseIdx]);
151 
152  if constexpr (enableEnergy)
153  Valgrind::CheckDefined((*enthalpy_)[storagePhaseIdx]);
154  }
155 
156  if constexpr (enableDissolution) {
157  Valgrind::CheckDefined(*Rs_);
158  Valgrind::CheckDefined(*Rv_);
159  }
160 
161  if constexpr (enableEvaporation) {
162  Valgrind::CheckDefined(*Rvw_);
163  }
164 
165  if constexpr (enableBrine) {
166  Valgrind::CheckDefined(*saltConcentration_);
167  }
168 
169  if constexpr (enableSaltPrecipitation) {
170  Valgrind::CheckDefined(*saltSaturation_);
171  }
172 
173  if constexpr (enableTemperature || enableEnergy)
174  Valgrind::CheckDefined(*temperature_);
175 #endif // NDEBUG
176  }
177 
182  template <class FluidState>
183  void assign(const FluidState& fs)
184  {
185  if constexpr (enableTemperature || enableEnergy)
186  setTemperature(fs.temperature(/*phaseIdx=*/0));
187 
188  unsigned pvtRegionIdx = getPvtRegionIndex_<FluidState>(fs);
189  setPvtRegionIndex(pvtRegionIdx);
190 
191  if constexpr (enableDissolution) {
192  setRs(BlackOil::getRs_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
193  setRv(BlackOil::getRv_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
194  }
195  if constexpr (enableEvaporation) {
196  setRvw(BlackOil::getRvw_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
197  }
198  if constexpr (enableBrine){
199  setSaltConcentration(BlackOil::getSaltConcentration_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
200  }
201  if constexpr (enableSaltPrecipitation){
202  setSaltSaturation(BlackOil::getSaltSaturation_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
203  }
204  for (unsigned storagePhaseIdx = 0; storagePhaseIdx < numStoragePhases; ++storagePhaseIdx) {
205  unsigned phaseIdx = storageToCanonicalPhaseIndex_(storagePhaseIdx);
206  setSaturation(phaseIdx, fs.saturation(phaseIdx));
207  setPressure(phaseIdx, fs.pressure(phaseIdx));
208  setDensity(phaseIdx, fs.density(phaseIdx));
209 
210  if constexpr (enableEnergy)
211  setEnthalpy(phaseIdx, fs.enthalpy(phaseIdx));
212 
213  setInvB(phaseIdx, getInvB_<FluidSystem, FluidState, Scalar>(fs, phaseIdx, pvtRegionIdx));
214  }
215  }
216 
223  void setPvtRegionIndex(unsigned newPvtRegionIdx)
224  { pvtRegionIdx_ = static_cast<unsigned short>(newPvtRegionIdx); }
225 
229  void setPressure(unsigned phaseIdx, const Scalar& p)
230  { pressure_[canonicalToStoragePhaseIndex_(phaseIdx)] = p; }
231 
235  void setSaturation(unsigned phaseIdx, const Scalar& S)
236  { saturation_[canonicalToStoragePhaseIndex_(phaseIdx)] = S; }
237 
241  void setPc(unsigned phaseIdx, const Scalar& pc)
242  { pc_[canonicalToStoragePhaseIndex_(phaseIdx)] = pc; }
243 
247  void setTotalSaturation(const Scalar& value)
248  {
249  totalSaturation_ = value;
250  }
251 
258  void setTemperature(const Scalar& value)
259  {
260  assert(enableTemperature || enableEnergy);
261 
262  (*temperature_) = value;
263  }
264 
271  void setEnthalpy(unsigned phaseIdx, const Scalar& value)
272  {
273  assert(enableTemperature || enableEnergy);
274 
275  (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)] = value;
276  }
277 
281  void setInvB(unsigned phaseIdx, const Scalar& b)
282  { invB_[canonicalToStoragePhaseIndex_(phaseIdx)] = b; }
283 
287  void setDensity(unsigned phaseIdx, const Scalar& rho)
288  { density_[canonicalToStoragePhaseIndex_(phaseIdx)] = rho; }
289 
295  void setRs(const Scalar& newRs)
296  { *Rs_ = newRs; }
297 
303  void setRv(const Scalar& newRv)
304  { *Rv_ = newRv; }
305 
311  void setRvw(const Scalar& newRvw)
312  { *Rvw_ = newRvw; }
313 
317  void setSaltConcentration(const Scalar& newSaltConcentration)
318  { *saltConcentration_ = newSaltConcentration; }
319 
323  void setSaltSaturation(const Scalar& newSaltSaturation)
324  { *saltSaturation_ = newSaltSaturation; }
325 
329  const Scalar& pressure(unsigned phaseIdx) const
330  { return pressure_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
331 
335  const Scalar& saturation(unsigned phaseIdx) const
336  { return saturation_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
337 
341  const Scalar& pc(unsigned phaseIdx) const
342  { return pc_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
343 
347  const Scalar& totalSaturation() const
348  {
349  return totalSaturation_;
350  }
351 
355  const Scalar& temperature(unsigned) const
356  {
357  if constexpr (enableTemperature || enableEnergy) {
358  return *temperature_;
359  } else {
360  static Scalar tmp(FluidSystem::reservoirTemperature(pvtRegionIdx_));
361  return tmp;
362  }
363  }
364 
371  const Scalar& invB(unsigned phaseIdx) const
372  { return invB_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
373 
381  const Scalar& Rs() const
382  {
383  if constexpr (enableDissolution) {
384  return *Rs_;
385  } else {
386  static Scalar null = 0.0;
387  return null;
388  }
389  }
390 
398  const Scalar& Rv() const
399  {
400  if constexpr (!enableDissolution) {
401  static Scalar null = 0.0;
402  return null;
403  } else {
404  return *Rv_;
405  }
406  }
407 
415  const Scalar& Rvw() const
416  {
417  if constexpr (enableEvaporation) {
418  return *Rvw_;
419  } else {
420  static Scalar null = 0.0;
421  return null;
422  }
423  }
424 
428  const Scalar& saltConcentration() const
429  {
430  if constexpr (enableBrine) {
431  return *saltConcentration_;
432  } else {
433  static Scalar null = 0.0;
434  return null;
435  }
436  }
437 
441  const Scalar& saltSaturation() const
442  {
443  if constexpr (enableSaltPrecipitation) {
444  return *saltSaturation_;
445  } else {
446  static Scalar null = 0.0;
447  return null;
448  }
449  }
450 
459  unsigned short pvtRegionIndex() const
460  { return pvtRegionIdx_; }
461 
465  Scalar density(unsigned phaseIdx) const
466  { return density_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
467 
474  const Scalar& enthalpy(unsigned phaseIdx) const
475  { return (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)]; }
476 
483  Scalar internalEnergy(unsigned phaseIdx) const
484  { return (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)] - pressure(phaseIdx)/density(phaseIdx); }
485 
487  // slow methods
489 
493  Scalar molarDensity(unsigned phaseIdx) const
494  {
495  const auto& rho = density(phaseIdx);
496 
497  if (phaseIdx == waterPhaseIdx)
498  return rho/FluidSystem::molarMass(waterCompIdx, pvtRegionIdx_);
499 
500  return
501  rho*(moleFraction(phaseIdx, gasCompIdx)/FluidSystem::molarMass(gasCompIdx, pvtRegionIdx_)
502  + moleFraction(phaseIdx, oilCompIdx)/FluidSystem::molarMass(oilCompIdx, pvtRegionIdx_));
503 
504  }
505 
511  Scalar molarVolume(unsigned phaseIdx) const
512  { return 1.0/molarDensity(phaseIdx); }
513 
517  Scalar viscosity(unsigned phaseIdx) const
518  { return FluidSystem::viscosity(*this, phaseIdx, pvtRegionIdx_); }
519 
523  Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
524  {
525  switch (phaseIdx) {
526  case waterPhaseIdx:
527  if (compIdx == waterCompIdx)
528  return 1.0;
529  return 0.0;
530 
531  case oilPhaseIdx:
532  if (compIdx == waterCompIdx)
533  return 0.0;
534  else if (compIdx == oilCompIdx)
535  return 1.0 - FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_);
536  else {
537  assert(compIdx == gasCompIdx);
538  return FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_);
539  }
540  break;
541 
542  case gasPhaseIdx:
543  if (compIdx == waterCompIdx)
544  return 0.0;
545  else if (compIdx == oilCompIdx)
546  return FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_);
547  else {
548  assert(compIdx == gasCompIdx);
549  return 1.0 - FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_);
550  }
551  break;
552  }
553 
554  throw std::logic_error("Invalid phase or component index!");
555  }
556 
560  Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
561  {
562  switch (phaseIdx) {
563  case waterPhaseIdx:
564  if (compIdx == waterCompIdx)
565  return 1.0;
566  return 0.0;
567 
568  case oilPhaseIdx:
569  if (compIdx == waterCompIdx)
570  return 0.0;
571  else if (compIdx == oilCompIdx)
572  return 1.0 - FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_),
573  pvtRegionIdx_);
574  else {
575  assert(compIdx == gasCompIdx);
576  return FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_),
577  pvtRegionIdx_);
578  }
579  break;
580 
581  case gasPhaseIdx:
582  if (compIdx == waterCompIdx)
583  return 0.0;
584  else if (compIdx == oilCompIdx)
585  return FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_),
586  pvtRegionIdx_);
587  else {
588  assert(compIdx == gasCompIdx);
589  return 1.0 - FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_),
590  pvtRegionIdx_);
591  }
592  break;
593  }
594 
595  throw std::logic_error("Invalid phase or component index!");
596  }
597 
601  Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
602  { return moleFraction(phaseIdx, compIdx)*molarDensity(phaseIdx); }
603 
607  Scalar averageMolarMass(unsigned phaseIdx) const
608  {
609  Scalar result(0.0);
610  for (unsigned compIdx = 0; compIdx < numComponents; ++ compIdx)
611  result += FluidSystem::molarMass(compIdx, pvtRegionIdx_)*moleFraction(phaseIdx, compIdx);
612  return result;
613  }
614 
618  Scalar fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
619  { return FluidSystem::fugacityCoefficient(*this, phaseIdx, compIdx, pvtRegionIdx_); }
620 
624  Scalar fugacity(unsigned phaseIdx, unsigned compIdx) const
625  {
626  return
627  fugacityCoefficient(phaseIdx, compIdx)
628  *moleFraction(phaseIdx, compIdx)
629  *pressure(phaseIdx);
630  }
631 
632 private:
633  static unsigned storageToCanonicalPhaseIndex_(unsigned storagePhaseIdx)
634  {
635  if constexpr (numStoragePhases == 3)
636  return storagePhaseIdx;
637  else
638  return FluidSystem::activeToCanonicalPhaseIdx(storagePhaseIdx);
639  }
640 
641  static unsigned canonicalToStoragePhaseIndex_(unsigned canonicalPhaseIdx)
642  {
643  if constexpr (numStoragePhases == 3)
644  return canonicalPhaseIdx;
645  else
646  return FluidSystem::canonicalToActivePhaseIdx(canonicalPhaseIdx);
647  }
648 
649  ConditionalStorage<enableTemperature || enableEnergy, Scalar> temperature_;
650  ConditionalStorage<enableEnergy, std::array<Scalar, numStoragePhases> > enthalpy_;
651  Scalar totalSaturation_;
652  std::array<Scalar, numStoragePhases> pressure_;
653  std::array<Scalar, numStoragePhases> pc_;
654  std::array<Scalar, numStoragePhases> saturation_;
655  std::array<Scalar, numStoragePhases> invB_;
656  std::array<Scalar, numStoragePhases> density_;
657  ConditionalStorage<enableDissolution,Scalar> Rs_;
658  ConditionalStorage<enableDissolution, Scalar> Rv_;
659  ConditionalStorage<enableEvaporation,Scalar> Rvw_;
660  ConditionalStorage<enableBrine, Scalar> saltConcentration_;
661  ConditionalStorage<enableSaltPrecipitation, Scalar> saltSaturation_;
662  unsigned short pvtRegionIdx_;
663 };
664 
665 } // namespace Opm
666 
667 #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:119
void setTemperature(const Scalar &value)
Set the temperature [K].
Definition: BlackOilFluidState.hpp:258
const Scalar & Rv() const
Return the oil vaporization factor of gas [m^3/m^3].
Definition: BlackOilFluidState.hpp:398
void setSaltSaturation(const Scalar &newSaltSaturation)
Set the solid salt saturation.
Definition: BlackOilFluidState.hpp:323
void setPc(unsigned phaseIdx, const Scalar &pc)
Set the capillary pressure of a fluid phase [-].
Definition: BlackOilFluidState.hpp:241
const Scalar & Rvw() const
Return the water vaporization factor of gas [m^3/m^3].
Definition: BlackOilFluidState.hpp:415
unsigned short pvtRegionIndex() const
Return the PVT region where the current fluid state is assumed to be part of.
Definition: BlackOilFluidState.hpp:459
void setDensity(unsigned phaseIdx, const Scalar &rho)
\ brief Set the density of a fluid phase
Definition: BlackOilFluidState.hpp:287
const Scalar & invB(unsigned phaseIdx) const
Return the inverse formation volume factor of a fluid phase [-].
Definition: BlackOilFluidState.hpp:371
const Scalar & enthalpy(unsigned phaseIdx) const
Return the specific enthalpy [J/kg] of a given fluid phase.
Definition: BlackOilFluidState.hpp:474
void setRs(const Scalar &newRs)
Set the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: BlackOilFluidState.hpp:295
Scalar fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
Return the fugacity coefficient of a component in a fluid phase [-].
Definition: BlackOilFluidState.hpp:618
void setSaturation(unsigned phaseIdx, const Scalar &S)
Set the saturation of a fluid phase [-].
Definition: BlackOilFluidState.hpp:235
void setInvB(unsigned phaseIdx, const Scalar &b)
\ brief Set the inverse formation volume factor of a fluid phase
Definition: BlackOilFluidState.hpp:281
Scalar averageMolarMass(unsigned phaseIdx) const
Return the partial molar density of a fluid phase [kg / mol].
Definition: BlackOilFluidState.hpp:607
const Scalar & pc(unsigned phaseIdx) const
Return the capillary pressure of a fluid phase [-].
Definition: BlackOilFluidState.hpp:341
void setRv(const Scalar &newRv)
Set the oil vaporization factor [m^3/m^3] of the gas phase.
Definition: BlackOilFluidState.hpp:303
const Scalar & saltSaturation() const
Return the saturation of solid salt.
Definition: BlackOilFluidState.hpp:441
void setPressure(unsigned phaseIdx, const Scalar &p)
Set the pressure of a fluid phase [-].
Definition: BlackOilFluidState.hpp:229
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: BlackOilFluidState.hpp:183
Scalar density(unsigned phaseIdx) const
Return the density [kg/m^3] of a given fluid phase.
Definition: BlackOilFluidState.hpp:465
void checkDefined() const
Make sure that all attributes are defined.
Definition: BlackOilFluidState.hpp:141
Scalar molarVolume(unsigned phaseIdx) const
Return the molar volume of a fluid phase [m^3/mol].
Definition: BlackOilFluidState.hpp:511
const Scalar & pressure(unsigned phaseIdx) const
Return the pressure of a fluid phase [Pa].
Definition: BlackOilFluidState.hpp:329
const Scalar & temperature(unsigned) const
Return the temperature [K].
Definition: BlackOilFluidState.hpp:355
const Scalar & saturation(unsigned phaseIdx) const
Return the saturation of a fluid phase [-].
Definition: BlackOilFluidState.hpp:335
Scalar molarDensity(unsigned phaseIdx) const
Return the molar density of a fluid phase [mol/m^3].
Definition: BlackOilFluidState.hpp:493
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
Return the mole fraction of a component in a fluid phase [-].
Definition: BlackOilFluidState.hpp:560
void setSaltConcentration(const Scalar &newSaltConcentration)
Set the salt concentration.
Definition: BlackOilFluidState.hpp:317
Scalar internalEnergy(unsigned phaseIdx) const
Return the specific internal energy [J/kg] of a given fluid phase.
Definition: BlackOilFluidState.hpp:483
const Scalar & totalSaturation() const
Return the total saturation needed for sequential.
Definition: BlackOilFluidState.hpp:347
Scalar viscosity(unsigned phaseIdx) const
Return the dynamic viscosity of a fluid phase [Pa s].
Definition: BlackOilFluidState.hpp:517
void setPvtRegionIndex(unsigned newPvtRegionIdx)
Set the index of the fluid region.
Definition: BlackOilFluidState.hpp:223
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:601
void setRvw(const Scalar &newRvw)
Set the water vaporization factor [m^3/m^3] of the gas phase.
Definition: BlackOilFluidState.hpp:311
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
Return the mass fraction of a component in a fluid phase [-].
Definition: BlackOilFluidState.hpp:523
const Scalar & Rs() const
Return the gas dissulition factor of oil [m^3/m^3].
Definition: BlackOilFluidState.hpp:381
void setTotalSaturation(const Scalar &value)
Set the total saturation used for sequential methods.
Definition: BlackOilFluidState.hpp:247
const Scalar & saltConcentration() const
Return the concentration of salt in water.
Definition: BlackOilFluidState.hpp:428
Scalar fugacity(unsigned phaseIdx, unsigned compIdx) const
Return the fugacity of a component in a fluid phase [Pa].
Definition: BlackOilFluidState.hpp:624
void setEnthalpy(unsigned phaseIdx, const Scalar &value)
Set the specific enthalpy [J/kg] of a given fluid phase.
Definition: BlackOilFluidState.hpp:271
A simple class which only stores a given member attribute if a boolean condition is true.
Definition: ConditionalStorage.hpp:46