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