My Project
DeadOilPvt.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 */
27 #ifndef OPM_DEAD_OIL_PVT_HPP
28 #define OPM_DEAD_OIL_PVT_HPP
29 
33 
34 #if HAVE_ECL_INPUT
35 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
36 #include <opm/parser/eclipse/EclipseState/Tables/PvdoTable.hpp>
37 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
38 #endif
39 
40 namespace Opm {
45 template <class Scalar>
47 {
48  typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
49 
50 public:
52 
53  DeadOilPvt() = default;
54  DeadOilPvt(const std::vector<Scalar>& oilReferenceDensity,
55  const std::vector<TabulatedOneDFunction>& inverseOilB,
56  const std::vector<TabulatedOneDFunction>& oilMu,
57  const std::vector<TabulatedOneDFunction>& inverseOilBMu)
58  : oilReferenceDensity_(oilReferenceDensity)
59  , inverseOilB_(inverseOilB)
60  , oilMu_(oilMu)
61  , inverseOilBMu_(inverseOilBMu)
62  { }
63 
64 #if HAVE_ECL_INPUT
68  void initFromState(const EclipseState& eclState, const Schedule&)
69  {
70  const auto& pvdoTables = eclState.getTableManager().getPvdoTables();
71  const auto& densityTable = eclState.getTableManager().getDensityTable();
72 
73  assert(pvdoTables.size() == densityTable.size());
74 
75  size_t numRegions = pvdoTables.size();
76  setNumRegions(numRegions);
77 
78  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
79  Scalar rhoRefO = densityTable[regionIdx].oil;
80  Scalar rhoRefG = densityTable[regionIdx].gas;
81  Scalar rhoRefW = densityTable[regionIdx].water;
82 
83  setReferenceDensities(regionIdx, rhoRefO, rhoRefG, rhoRefW);
84 
85  const auto& pvdoTable = pvdoTables.getTable<PvdoTable>(regionIdx);
86 
87  const auto& BColumn(pvdoTable.getFormationFactorColumn());
88  std::vector<Scalar> invBColumn(BColumn.size());
89  for (unsigned i = 0; i < invBColumn.size(); ++i)
90  invBColumn[i] = 1/BColumn[i];
91 
92  inverseOilB_[regionIdx].setXYArrays(pvdoTable.numRows(),
93  pvdoTable.getPressureColumn(),
94  invBColumn);
95  oilMu_[regionIdx].setXYArrays(pvdoTable.numRows(),
96  pvdoTable.getPressureColumn(),
97  pvdoTable.getViscosityColumn());
98  }
99 
100  initEnd();
101  }
102 #endif // HAVE_ECL_INPUT
103 
104  void setNumRegions(size_t numRegions)
105  {
106  oilReferenceDensity_.resize(numRegions);
107  inverseOilB_.resize(numRegions);
108  inverseOilBMu_.resize(numRegions);
109  oilMu_.resize(numRegions);
110  }
111 
115  void setReferenceDensities(unsigned regionIdx,
116  Scalar rhoRefOil,
117  Scalar /*rhoRefGas*/,
118  Scalar /*rhoRefWater*/)
119  {
120  oilReferenceDensity_[regionIdx] = rhoRefOil;
121  }
122 
133  void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedOneDFunction& invBo)
134  { inverseOilB_[regionIdx] = invBo; }
135 
141  void setOilViscosity(unsigned regionIdx, const TabulatedOneDFunction& muo)
142  { oilMu_[regionIdx] = muo; }
143 
147  void initEnd()
148  {
149  // calculate the final 2D functions which are used for interpolation.
150  size_t numRegions = oilMu_.size();
151  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
152  // calculate the table which stores the inverse of the product of the oil
153  // formation volume factor and the oil viscosity
154  const auto& oilMu = oilMu_[regionIdx];
155  const auto& invOilB = inverseOilB_[regionIdx];
156  assert(oilMu.numSamples() == invOilB.numSamples());
157 
158  std::vector<Scalar> invBMuColumn;
159  std::vector<Scalar> pressureColumn;
160  invBMuColumn.resize(oilMu.numSamples());
161  pressureColumn.resize(oilMu.numSamples());
162 
163  for (unsigned pIdx = 0; pIdx < oilMu.numSamples(); ++pIdx) {
164  pressureColumn[pIdx] = invOilB.xAt(pIdx);
165  invBMuColumn[pIdx] = invOilB.valueAt(pIdx)*1/oilMu.valueAt(pIdx);
166  }
167 
168  inverseOilBMu_[regionIdx].setXYArrays(pressureColumn.size(),
169  pressureColumn,
170  invBMuColumn);
171  }
172  }
173 
177  unsigned numRegions() const
178  { return inverseOilBMu_.size(); }
179 
183  template <class Evaluation>
184  Evaluation internalEnergy(unsigned,
185  const Evaluation&,
186  const Evaluation&,
187  const Evaluation&) const
188  {
189  throw std::runtime_error("Requested the enthalpy of oil but the thermal option is not enabled");
190  }
191 
195  template <class Evaluation>
196  Evaluation viscosity(unsigned regionIdx,
197  const Evaluation& temperature,
198  const Evaluation& pressure,
199  const Evaluation& /*Rs*/) const
200  { return saturatedViscosity(regionIdx, temperature, pressure); }
201 
205  template <class Evaluation>
206  Evaluation saturatedViscosity(unsigned regionIdx,
207  const Evaluation& /*temperature*/,
208  const Evaluation& pressure) const
209  {
210  const Evaluation& invBo = inverseOilB_[regionIdx].eval(pressure, /*extrapolate=*/true);
211  const Evaluation& invMuoBo = inverseOilBMu_[regionIdx].eval(pressure, /*extrapolate=*/true);
212 
213  return invBo/invMuoBo;
214  }
215 
219  template <class Evaluation>
220  Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
221  const Evaluation& /*temperature*/,
222  const Evaluation& pressure,
223  const Evaluation& /*Rs*/) const
224  { return inverseOilB_[regionIdx].eval(pressure, /*extrapolate=*/true); }
225 
231  template <class Evaluation>
232  Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
233  const Evaluation& /*temperature*/,
234  const Evaluation& pressure) const
235  { return inverseOilB_[regionIdx].eval(pressure, /*extrapolate=*/true); }
236 
240  template <class Evaluation>
241  Evaluation saturatedGasDissolutionFactor(unsigned /*regionIdx*/,
242  const Evaluation& /*temperature*/,
243  const Evaluation& /*pressure*/) const
244  { return 0.0; /* this is dead oil! */ }
245 
249  template <class Evaluation>
250  Evaluation saturatedGasDissolutionFactor(unsigned /*regionIdx*/,
251  const Evaluation& /*temperature*/,
252  const Evaluation& /*pressure*/,
253  const Evaluation& /*oilSaturation*/,
254  const Evaluation& /*maxOilSaturation*/) const
255  { return 0.0; /* this is dead oil! */ }
256 
263  template <class Evaluation>
264  Evaluation saturationPressure(unsigned /*regionIdx*/,
265  const Evaluation& /*temperature*/,
266  const Evaluation& /*Rs*/) const
267  { return 0.0; /* this is dead oil, so there isn't any meaningful saturation pressure! */ }
268 
269  template <class Evaluation>
270  Evaluation saturatedGasMassFraction(unsigned /*regionIdx*/,
271  const Evaluation& /*temperature*/,
272  const Evaluation& /*pressure*/) const
273  { return 0.0; /* this is dead oil! */ }
274 
275  template <class Evaluation>
276  Evaluation saturatedGasMoleFraction(unsigned /*regionIdx*/,
277  const Evaluation& /*temperature*/,
278  const Evaluation& /*pressure*/) const
279  { return 0.0; /* this is dead oil! */ }
280 
281  template <class Evaluation>
282  Evaluation diffusionCoefficient(const Evaluation& /*temperature*/,
283  const Evaluation& /*pressure*/,
284  unsigned /*compIdx*/) const
285  {
286  throw std::runtime_error("Not implemented: The PVT model does not provide a diffusionCoefficient()");
287  }
288 
289  const Scalar oilReferenceDensity(unsigned regionIdx) const
290  { return oilReferenceDensity_[regionIdx]; }
291 
292  const std::vector<TabulatedOneDFunction>& inverseOilB() const
293  { return inverseOilB_; }
294 
295  const std::vector<TabulatedOneDFunction>& oilMu() const
296  { return oilMu_; }
297 
298  const std::vector<TabulatedOneDFunction>& inverseOilBMu() const
299  { return inverseOilBMu_; }
300 
301  bool operator==(const DeadOilPvt<Scalar>& data) const
302  {
303  return this->oilReferenceDensity_ == data.oilReferenceDensity_ &&
304  this->inverseOilB() == data.inverseOilB() &&
305  this->oilMu() == data.oilMu() &&
306  this->inverseOilBMu() == data.inverseOilBMu();
307  }
308 
309 private:
310  std::vector<Scalar> oilReferenceDensity_;
311  std::vector<TabulatedOneDFunction> inverseOilB_;
312  std::vector<TabulatedOneDFunction> oilMu_;
313  std::vector<TabulatedOneDFunction> inverseOilBMu_;
314 };
315 
316 } // namespace Opm
317 
318 #endif
Class implementing cubic splines.
Implements a linearly interpolated scalar function that depends on one variable.
Implements a scalar function that depends on two variables and which is sampled uniformly in the X di...
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
Definition: DeadOilPvt.hpp:47
void setOilViscosity(unsigned regionIdx, const TabulatedOneDFunction &muo)
Initialize the viscosity of the oil phase.
Definition: DeadOilPvt.hpp:141
Evaluation saturatedGasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: DeadOilPvt.hpp:241
void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedOneDFunction &invBo)
Initialize the function for the oil formation volume factor.
Definition: DeadOilPvt.hpp:133
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: DeadOilPvt.hpp:196
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: DeadOilPvt.hpp:115
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: DeadOilPvt.hpp:147
Evaluation saturatedViscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of gas saturated oil given a pressure.
Definition: DeadOilPvt.hpp:206
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of oil given a set of parameters.
Definition: DeadOilPvt.hpp:184
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of saturated oil.
Definition: DeadOilPvt.hpp:232
Evaluation saturationPressure(unsigned, const Evaluation &, const Evaluation &) const
Returns the saturation pressure of the oil phase [Pa] depending on its mass fraction of the gas compo...
Definition: DeadOilPvt.hpp:264
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &) const
Returns the formation volume factor [-] of the fluid phase.
Definition: DeadOilPvt.hpp:220
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: DeadOilPvt.hpp:177
Evaluation saturatedGasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: DeadOilPvt.hpp:250
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47