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 
97 template <class ScalarT,
98  class FluidSystem,
99  bool enableTemperature = false,
100  bool enableEnergy = false,
101  bool enableDissolution = true,
102  bool enableBrine = false,
103  unsigned numStoragePhases = FluidSystem::numPhases>
105 {
106  enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
107  enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
108  enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
109 
110  enum { waterCompIdx = FluidSystem::waterCompIdx };
111  enum { gasCompIdx = FluidSystem::gasCompIdx };
112  enum { oilCompIdx = FluidSystem::oilCompIdx };
113 
114 public:
115  typedef ScalarT Scalar;
116  enum { numPhases = FluidSystem::numPhases };
117  enum { numComponents = FluidSystem::numComponents };
118 
127  void checkDefined() const
128  {
129 #ifndef NDEBUG
130  Valgrind::CheckDefined(pvtRegionIdx_);
131 
132  for (unsigned storagePhaseIdx = 0; storagePhaseIdx < numStoragePhases; ++ storagePhaseIdx) {
133  Valgrind::CheckDefined(saturation_[storagePhaseIdx]);
134  Valgrind::CheckDefined(pressure_[storagePhaseIdx]);
135  Valgrind::CheckDefined(density_[storagePhaseIdx]);
136  Valgrind::CheckDefined(invB_[storagePhaseIdx]);
137 
138  if (enableEnergy)
139  Valgrind::CheckDefined((*enthalpy_)[storagePhaseIdx]);
140  }
141 
142  if (enableDissolution) {
143  Valgrind::CheckDefined(*Rs_);
144  Valgrind::CheckDefined(*Rv_);
145  }
146 
147  if (enableBrine) {
148  Valgrind::CheckDefined(*saltConcentration_);
149  }
150 
151  if (enableTemperature || enableEnergy)
152  Valgrind::CheckDefined(*temperature_);
153 #endif // NDEBUG
154  }
155 
160  template <class FluidState>
161  void assign(const FluidState& fs)
162  {
163  if (enableTemperature || enableEnergy)
164  setTemperature(fs.temperature(/*phaseIdx=*/0));
165 
166  unsigned pvtRegionIdx = getPvtRegionIndex_<FluidState>(fs);
167  setPvtRegionIndex(pvtRegionIdx);
168 
169  if (enableDissolution) {
170  setRs(BlackOil::getRs_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
171  setRv(BlackOil::getRv_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
172  }
173 
174  if (enableBrine){
175  setSaltConcentration(BlackOil::getSaltConcentration_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
176  }
177  for (unsigned storagePhaseIdx = 0; storagePhaseIdx < numStoragePhases; ++storagePhaseIdx) {
178  unsigned phaseIdx = storageToCanonicalPhaseIndex_(storagePhaseIdx);
179  setSaturation(phaseIdx, fs.saturation(phaseIdx));
180  setPressure(phaseIdx, fs.pressure(phaseIdx));
181  setDensity(phaseIdx, fs.density(phaseIdx));
182 
183  if (enableEnergy)
184  setEnthalpy(phaseIdx, fs.enthalpy(phaseIdx));
185 
186  setInvB(phaseIdx, getInvB_<FluidSystem, FluidState, Scalar>(fs, phaseIdx, pvtRegionIdx));
187  }
188  }
189 
196  void setPvtRegionIndex(unsigned newPvtRegionIdx)
197  { pvtRegionIdx_ = static_cast<unsigned short>(newPvtRegionIdx); }
198 
202  void setPressure(unsigned phaseIdx, const Scalar& p)
203  { pressure_[canonicalToStoragePhaseIndex_(phaseIdx)] = p; }
204 
208  void setSaturation(unsigned phaseIdx, const Scalar& S)
209  { saturation_[canonicalToStoragePhaseIndex_(phaseIdx)] = S; }
210 
214  void setPc(unsigned phaseIdx, const Scalar& pc)
215  { pc_[canonicalToStoragePhaseIndex_(phaseIdx)] = pc; }
216 
220  void setTotalSaturation(const Scalar& value)
221  {
222  totalSaturation_ = value;
223  }
224 
231  void setTemperature(const Scalar& value)
232  {
233  assert(enableTemperature || enableEnergy);
234 
235  (*temperature_) = value;
236  }
237 
244  void setEnthalpy(unsigned phaseIdx, const Scalar& value)
245  {
246  assert(enableTemperature || enableEnergy);
247 
248  (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)] = value;
249  }
250 
254  void setInvB(unsigned phaseIdx, const Scalar& b)
255  { invB_[canonicalToStoragePhaseIndex_(phaseIdx)] = b; }
256 
260  void setDensity(unsigned phaseIdx, const Scalar& rho)
261  { density_[canonicalToStoragePhaseIndex_(phaseIdx)] = rho; }
262 
268  void setRs(const Scalar& newRs)
269  { *Rs_ = newRs; }
270 
276  void setRv(const Scalar& newRv)
277  { *Rv_ = newRv; }
278 
282  void setSaltConcentration(const Scalar& newSaltConcentration)
283  { *saltConcentration_ = newSaltConcentration; }
284 
288  const Scalar& pressure(unsigned phaseIdx) const
289  { return pressure_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
290 
294  const Scalar& saturation(unsigned phaseIdx) const
295  { return saturation_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
296 
300  const Scalar& pc(unsigned phaseIdx) const
301  { return pc_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
302 
306  const Scalar& totalSaturation() const
307  {
308  return totalSaturation_;
309  }
310 
314  const Scalar& temperature(unsigned) const
315  {
316  if (!enableTemperature && !enableEnergy) {
317  static Scalar tmp(FluidSystem::reservoirTemperature(pvtRegionIdx_));
318  return tmp;
319  }
320 
321  return *temperature_;
322  }
323 
330  const Scalar& invB(unsigned phaseIdx) const
331  { return invB_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
332 
340  const Scalar& Rs() const
341  {
342  if (!enableDissolution) {
343  static Scalar null = 0.0;
344  return null;
345  }
346 
347  return *Rs_;
348  }
349 
357  const Scalar& Rv() const
358  {
359  if (!enableDissolution) {
360  static Scalar null = 0.0;
361  return null;
362  }
363 
364  return *Rv_;
365  }
366 
370  const Scalar& saltConcentration() const
371  {
372  if (!enableBrine) {
373  static Scalar null = 0.0;
374  return null;
375  }
376 
377  return *saltConcentration_;
378  }
379 
388  unsigned short pvtRegionIndex() const
389  { return pvtRegionIdx_; }
390 
394  Scalar density(unsigned phaseIdx) const
395  { return density_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
396 
403  const Scalar& enthalpy(unsigned phaseIdx) const
404  { return (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)]; }
405 
412  Scalar internalEnergy(unsigned phaseIdx) const
413  { return (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)] - pressure(phaseIdx)/density(phaseIdx); }
414 
416  // slow methods
418 
422  Scalar molarDensity(unsigned phaseIdx) const
423  {
424  const auto& rho = density(phaseIdx);
425 
426  if (phaseIdx == waterPhaseIdx)
427  return rho/FluidSystem::molarMass(waterCompIdx, pvtRegionIdx_);
428 
429  return
430  rho*(moleFraction(phaseIdx, gasCompIdx)/FluidSystem::molarMass(gasCompIdx, pvtRegionIdx_)
431  + moleFraction(phaseIdx, oilCompIdx)/FluidSystem::molarMass(oilCompIdx, pvtRegionIdx_));
432 
433  }
434 
440  Scalar molarVolume(unsigned phaseIdx) const
441  { return 1.0/molarDensity(phaseIdx); }
442 
446  Scalar viscosity(unsigned phaseIdx) const
447  { return FluidSystem::viscosity(*this, phaseIdx, pvtRegionIdx_); }
448 
452  Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
453  {
454  switch (phaseIdx) {
455  case waterPhaseIdx:
456  if (compIdx == waterCompIdx)
457  return 1.0;
458  return 0.0;
459 
460  case oilPhaseIdx:
461  if (compIdx == waterCompIdx)
462  return 0.0;
463  else if (compIdx == oilCompIdx)
464  return 1.0 - FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_);
465  else {
466  assert(compIdx == gasCompIdx);
467  return FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_);
468  }
469  break;
470 
471  case gasPhaseIdx:
472  if (compIdx == waterCompIdx)
473  return 0.0;
474  else if (compIdx == oilCompIdx)
475  return FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_);
476  else {
477  assert(compIdx == gasCompIdx);
478  return 1.0 - FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_);
479  }
480  break;
481  }
482 
483  throw std::logic_error("Invalid phase or component index!");
484  }
485 
489  Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
490  {
491  switch (phaseIdx) {
492  case waterPhaseIdx:
493  if (compIdx == waterCompIdx)
494  return 1.0;
495  return 0.0;
496 
497  case oilPhaseIdx:
498  if (compIdx == waterCompIdx)
499  return 0.0;
500  else if (compIdx == oilCompIdx)
501  return 1.0 - FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_),
502  pvtRegionIdx_);
503  else {
504  assert(compIdx == gasCompIdx);
505  return FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs(), pvtRegionIdx_),
506  pvtRegionIdx_);
507  }
508  break;
509 
510  case gasPhaseIdx:
511  if (compIdx == waterCompIdx)
512  return 0.0;
513  else if (compIdx == oilCompIdx)
514  return FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_),
515  pvtRegionIdx_);
516  else {
517  assert(compIdx == gasCompIdx);
518  return 1.0 - FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv(), pvtRegionIdx_),
519  pvtRegionIdx_);
520  }
521  break;
522  }
523 
524  throw std::logic_error("Invalid phase or component index!");
525  }
526 
530  Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
531  { return moleFraction(phaseIdx, compIdx)*molarDensity(phaseIdx); }
532 
536  Scalar averageMolarMass(unsigned phaseIdx) const
537  {
538  Scalar result(0.0);
539  for (unsigned compIdx = 0; compIdx < numComponents; ++ compIdx)
540  result += FluidSystem::molarMass(compIdx, pvtRegionIdx_)*moleFraction(phaseIdx, compIdx);
541  return result;
542  }
543 
547  Scalar fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
548  { return FluidSystem::fugacityCoefficient(*this, phaseIdx, compIdx, pvtRegionIdx_); }
549 
553  Scalar fugacity(unsigned phaseIdx, unsigned compIdx) const
554  {
555  return
556  fugacityCoefficient(phaseIdx, compIdx)
557  *moleFraction(phaseIdx, compIdx)
558  *pressure(phaseIdx);
559  }
560 
561 private:
562  static unsigned storageToCanonicalPhaseIndex_(unsigned storagePhaseIdx)
563  {
564  if (numStoragePhases == 3)
565  return storagePhaseIdx;
566 
567  return FluidSystem::activeToCanonicalPhaseIdx(storagePhaseIdx);
568  }
569 
570  static unsigned canonicalToStoragePhaseIndex_(unsigned canonicalPhaseIdx)
571  {
572  if (numStoragePhases == 3)
573  return canonicalPhaseIdx;
574 
575  return FluidSystem::canonicalToActivePhaseIdx(canonicalPhaseIdx);
576  }
577 
578  ConditionalStorage<enableTemperature || enableEnergy, Scalar> temperature_;
579  ConditionalStorage<enableEnergy, std::array<Scalar, numStoragePhases> > enthalpy_;
580  Scalar totalSaturation_;
581  std::array<Scalar, numStoragePhases> pressure_;
582  std::array<Scalar, numStoragePhases> pc_;
583  std::array<Scalar, numStoragePhases> saturation_;
584  std::array<Scalar, numStoragePhases> invB_;
585  std::array<Scalar, numStoragePhases> density_;
586  ConditionalStorage<enableDissolution,Scalar> Rs_;
587  ConditionalStorage<enableDissolution, Scalar> Rv_;
588  ConditionalStorage<enableBrine, Scalar> saltConcentration_;
589  unsigned short pvtRegionIdx_;
590 };
591 
592 } // namespace Opm
593 
594 #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:105
void setDensity(unsigned phaseIdx, const Scalar &rho)
\ brief Set the density of a fluid phase
Definition: BlackOilFluidState.hpp:260
void setInvB(unsigned phaseIdx, const Scalar &b)
\ brief Set the inverse formation volume factor of a fluid phase
Definition: BlackOilFluidState.hpp:254
const Scalar & temperature(unsigned) const
Return the temperature [K].
Definition: BlackOilFluidState.hpp:314
void checkDefined() const
Make sure that all attributes are defined.
Definition: BlackOilFluidState.hpp:127
const Scalar & Rs() const
Return the gas dissulition factor of oil [m^3/m^3].
Definition: BlackOilFluidState.hpp:340
Scalar fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
Return the fugacity coefficient of a component in a fluid phase [-].
Definition: BlackOilFluidState.hpp:547
Scalar molarDensity(unsigned phaseIdx) const
Return the molar density of a fluid phase [mol/m^3].
Definition: BlackOilFluidState.hpp:422
Scalar molarVolume(unsigned phaseIdx) const
Return the molar volume of a fluid phase [m^3/mol].
Definition: BlackOilFluidState.hpp:440
const Scalar & totalSaturation() const
Return the total saturation needed for sequential.
Definition: BlackOilFluidState.hpp:306
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: BlackOilFluidState.hpp:161
Scalar internalEnergy(unsigned phaseIdx) const
Return the specific internal energy [J/kg] of a given fluid phase.
Definition: BlackOilFluidState.hpp:412
void setPc(unsigned phaseIdx, const Scalar &pc)
Set the capillary pressure of a fluid phase [-].
Definition: BlackOilFluidState.hpp:214
void setRs(const Scalar &newRs)
Set the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: BlackOilFluidState.hpp:268
void setEnthalpy(unsigned phaseIdx, const Scalar &value)
Set the specific enthalpy [J/kg] of a given fluid phase.
Definition: BlackOilFluidState.hpp:244
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
Return the mole fraction of a component in a fluid phase [-].
Definition: BlackOilFluidState.hpp:489
void setRv(const Scalar &newRv)
Set the oil vaporization factor [m^3/m^3] of the gas phase.
Definition: BlackOilFluidState.hpp:276
const Scalar & Rv() const
Return the oil vaporization factor of gas [m^3/m^3].
Definition: BlackOilFluidState.hpp:357
const Scalar & saturation(unsigned phaseIdx) const
Return the saturation of a fluid phase [-].
Definition: BlackOilFluidState.hpp:294
Scalar fugacity(unsigned phaseIdx, unsigned compIdx) const
Return the fugacity of a component in a fluid phase [Pa].
Definition: BlackOilFluidState.hpp:553
void setSaltConcentration(const Scalar &newSaltConcentration)
Set the salt concentration.
Definition: BlackOilFluidState.hpp:282
unsigned short pvtRegionIndex() const
Return the PVT region where the current fluid state is assumed to be part of.
Definition: BlackOilFluidState.hpp:388
void setPressure(unsigned phaseIdx, const Scalar &p)
Set the pressure of a fluid phase [-].
Definition: BlackOilFluidState.hpp:202
void setSaturation(unsigned phaseIdx, const Scalar &S)
Set the saturation of a fluid phase [-].
Definition: BlackOilFluidState.hpp:208
const Scalar & enthalpy(unsigned phaseIdx) const
Return the specific enthalpy [J/kg] of a given fluid phase.
Definition: BlackOilFluidState.hpp:403
const Scalar & pressure(unsigned phaseIdx) const
Return the pressure of a fluid phase [Pa].
Definition: BlackOilFluidState.hpp:288
Scalar viscosity(unsigned phaseIdx) const
Return the dynamic viscosity of a fluid phase [Pa s].
Definition: BlackOilFluidState.hpp:446
const Scalar & pc(unsigned phaseIdx) const
Return the capillary pressure of a fluid phase [-].
Definition: BlackOilFluidState.hpp:300
void setPvtRegionIndex(unsigned newPvtRegionIdx)
Set the index of the fluid region.
Definition: BlackOilFluidState.hpp:196
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:530
const Scalar & invB(unsigned phaseIdx) const
Return the inverse formation volume factor of a fluid phase [-].
Definition: BlackOilFluidState.hpp:330
Scalar density(unsigned phaseIdx) const
Return the density [kg/m^3] of a given fluid phase.
Definition: BlackOilFluidState.hpp:394
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
Return the mass fraction of a component in a fluid phase [-].
Definition: BlackOilFluidState.hpp:452
Scalar averageMolarMass(unsigned phaseIdx) const
Return the partial molar density of a fluid phase [kg / mol].
Definition: BlackOilFluidState.hpp:536
void setTemperature(const Scalar &value)
Set the temperature [K].
Definition: BlackOilFluidState.hpp:231
void setTotalSaturation(const Scalar &value)
Set the total saturation used for sequential methods.
Definition: BlackOilFluidState.hpp:220
const Scalar & saltConcentration() const
Return the concentration of salt in water.
Definition: BlackOilFluidState.hpp:370