blackoilfluidstate.hh
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 EWOMS_BLACK_OIL_FLUID_STATE_HH
29 #define EWOMS_BLACK_OIL_FLUID_STATE_HH
30 
31 #include "blackoilproperties.hh"
32 
33 #include <opm/common/Valgrind.hpp>
34 #include <opm/common/Unused.hpp>
35 #include <opm/common/ErrorMacros.hpp>
36 #include <opm/common/Exceptions.hpp>
37 
38 namespace Ewoms {
48 template <class TypeTag>
50 {
51  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
52  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
53 
54  enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
55  enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
56  enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
57 
58  enum { waterCompIdx = FluidSystem::waterCompIdx };
59  enum { gasCompIdx = FluidSystem::gasCompIdx };
60  enum { oilCompIdx = FluidSystem::oilCompIdx };
61 
62 public:
63  typedef Evaluation Scalar;
64  enum { numPhases = FluidSystem::numPhases };
65  enum { numComponents = FluidSystem::numComponents };
66 
75  void checkDefined() const
76  {
77 #ifndef NDEBUG
78  Opm::Valgrind::CheckDefined(pvtRegionIdx_);
79 
80  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
81  Opm::Valgrind::CheckDefined(saturation_[phaseIdx]);
82  Opm::Valgrind::CheckDefined(pressure_[phaseIdx]);
83  Opm::Valgrind::CheckDefined(invB_[phaseIdx]);
84  }
85 
86  Opm::Valgrind::CheckDefined(Rs_);
87  Opm::Valgrind::CheckDefined(Rv_);
88 #endif // NDEBUG
89  }
90 
95  template <class FluidState>
96  void assign(const FluidState& fs OPM_UNUSED)
97  {
98  assert(false); // not yet implemented
99  }
100 
101  void setPvtRegionIndex(unsigned newPvtRegionIdx)
102  { pvtRegionIdx_ = static_cast<unsigned short>(newPvtRegionIdx); }
103 
104  void setPressure(unsigned phaseIdx, const Evaluation& p)
105  { pressure_[phaseIdx] = p; }
106 
107  void setSaturation(unsigned phaseIdx, const Evaluation& S)
108  { saturation_[phaseIdx] = S; }
109 
110  void setInvB(unsigned phaseIdx, const Evaluation& b)
111  { invB_[phaseIdx] = b; }
112 
113  void setDensity(unsigned phaseIdx, const Evaluation& rho)
114  { density_[phaseIdx] = rho; }
115 
116  void setRs(const Evaluation& newRs)
117  { Rs_ = newRs; }
118 
119  void setRv(const Evaluation& newRv)
120  { Rv_ = newRv; }
121 
122  const Evaluation& pressure(unsigned phaseIdx) const
123  { return pressure_[phaseIdx]; }
124 
125  const Evaluation& saturation(unsigned phaseIdx) const
126  { return saturation_[phaseIdx]; }
127 
128  const Evaluation& temperature(unsigned phaseIdx OPM_UNUSED) const
129  { return temperature_; }
130 
131  const Evaluation& invB(unsigned phaseIdx) const
132  { return invB_[phaseIdx]; }
133 
134  const Evaluation& Rs() const
135  { return Rs_; }
136 
137  const Evaluation& Rv() const
138  { return Rv_; }
139 
140  unsigned short pvtRegionIndex() const
141  { return pvtRegionIdx_; }
142 
143  bool phaseIsPresent(unsigned phaseIdx) const
144  { return saturation_[phaseIdx] > 0.0; }
145 
147  // slow methods
149  Evaluation density(unsigned phaseIdx) const
150  { return density_[phaseIdx]; }
151 
152  Evaluation molarDensity(unsigned phaseIdx) const
153  {
154  const auto& rho = density(phaseIdx);
155 
156  if (phaseIdx == waterPhaseIdx)
157  return rho/FluidSystem::molarMass(waterCompIdx, pvtRegionIdx_);
158 
159  return
160  rho*(moleFraction(phaseIdx, gasCompIdx)/FluidSystem::molarMass(gasCompIdx, pvtRegionIdx_)
161  + moleFraction(phaseIdx, oilCompIdx)/FluidSystem::molarMass(oilCompIdx, pvtRegionIdx_));
162 
163  }
164 
165  Evaluation molarVolume(unsigned phaseIdx) const
166  { return 1.0/molarDensity(phaseIdx); }
167 
168  Evaluation viscosity(unsigned phaseIdx) const
169  { return FluidSystem::viscosity(*this, phaseIdx, pvtRegionIdx_); }
170 
171  Evaluation enthalpy(unsigned phaseIdx OPM_UNUSED) const
172  {
173  OPM_THROW(Opm::NotImplemented,
174  "The black-oil model does not support energy conservation yet.");
175  }
176 
177  Evaluation internalEnergy(unsigned phaseIdx OPM_UNUSED) const
178  {
179  OPM_THROW(Opm::NotImplemented,
180  "The black-oil model does not support energy conservation yet.");
181  }
182 
183  Evaluation massFraction(unsigned phaseIdx, unsigned compIdx) const
184  {
185  switch (phaseIdx) {
186  case waterPhaseIdx:
187  if (compIdx == waterCompIdx)
188  return 1.0;
189  return 0.0;
190 
191  case oilPhaseIdx:
192  if (compIdx == waterCompIdx)
193  return 0.0;
194  else if (compIdx == oilCompIdx)
195  return 1.0 - FluidSystem::convertRsToXoG(Rs_, pvtRegionIdx_);
196  else {
197  assert(compIdx == gasCompIdx);
198  return FluidSystem::convertRsToXoG(Rs_, pvtRegionIdx_);
199  }
200  break;
201 
202  case gasPhaseIdx:
203  if (compIdx == waterCompIdx)
204  return 0.0;
205  else if (compIdx == oilCompIdx)
206  return FluidSystem::convertRvToXgO(Rv_, pvtRegionIdx_);
207  else {
208  assert(compIdx == gasCompIdx);
209  return 1.0 - FluidSystem::convertRvToXgO(Rv_, pvtRegionIdx_);
210  }
211  break;
212  }
213 
214  OPM_THROW(std::logic_error,
215  "Invalid phase or component index!");
216  }
217 
218  Evaluation moleFraction(unsigned phaseIdx, unsigned compIdx) const
219  {
220  switch (phaseIdx) {
221  case waterPhaseIdx:
222  if (compIdx == waterCompIdx)
223  return 1.0;
224  return 0.0;
225 
226  case oilPhaseIdx:
227  if (compIdx == waterCompIdx)
228  return 0.0;
229  else if (compIdx == oilCompIdx)
230  return 1.0 - FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs_, pvtRegionIdx_),
231  pvtRegionIdx_);
232  else {
233  assert(compIdx == gasCompIdx);
234  return FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs_, pvtRegionIdx_),
235  pvtRegionIdx_);
236  }
237  break;
238 
239  case gasPhaseIdx:
240  if (compIdx == waterCompIdx)
241  return 0.0;
242  else if (compIdx == oilCompIdx)
243  return FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv_, pvtRegionIdx_),
244  pvtRegionIdx_);
245  else {
246  assert(compIdx == gasCompIdx);
247  return 1.0 - FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv_, pvtRegionIdx_),
248  pvtRegionIdx_);
249  }
250  break;
251  }
252 
253  OPM_THROW(std::logic_error,
254  "Invalid phase or component index!");
255  }
256 
257  Evaluation molarity(unsigned phaseIdx, unsigned compIdx) const
258  { return moleFraction(phaseIdx, compIdx)*molarDensity(phaseIdx); }
259 
260  Evaluation averageMolarMass(unsigned phaseIdx) const
261  {
262  Evaluation result(0.0);
263  for (unsigned compIdx = 0; compIdx < numComponents; ++ compIdx)
264  result += FluidSystem::molarMass(compIdx, pvtRegionIdx_)*moleFraction(phaseIdx, compIdx);
265  return result;
266  }
267 
268  Evaluation fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
269  { return FluidSystem::fugacityCoefficient(*this, phaseIdx, compIdx, pvtRegionIdx_); }
270 
271  Evaluation fugacity(unsigned phaseIdx, unsigned compIdx) const
272  {
273  return
274  fugacityCoefficient(phaseIdx, compIdx)
275  *moleFraction(phaseIdx, compIdx)
276  *pressure(phaseIdx);
277  }
278 
279 private:
280  static const Evaluation temperature_;
281  std::array<Evaluation, numPhases> pressure_;
282  std::array<Evaluation, numPhases> saturation_;
283  std::array<Evaluation, numPhases> invB_;
284  std::array<Evaluation, numPhases> density_;
285  Evaluation Rs_;
286  Evaluation Rv_;
287  unsigned short pvtRegionIdx_;
288 };
289 
290 template <class TypeTag>
291 const typename BlackOilFluidState<TypeTag>::Evaluation BlackOilFluidState<TypeTag>::temperature_ =
292  BlackOilFluidState<TypeTag>::FluidSystem::surfaceTemperature;
293 
294 } // namespace Ewoms
295 
296 #endif
Definition: baseauxiliarymodule.hh:37
void assign(const FluidState &fs OPM_UNUSED)
Retrieve all parameters from an arbitrary fluid state.
Definition: blackoilfluidstate.hh:96
void checkDefined() const
Make sure that all attributes are defined.
Definition: blackoilfluidstate.hh:75
Declares the properties required by the black oil model.
Implements a "taylor-made" fluid state class for the black-oil model.
Definition: blackoilfluidstate.hh:49