All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
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_OPM_PARSER
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 {
49  typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
50 
51 public:
52 #if HAVE_OPM_PARSER
53 
56  void initFromDeck(const Deck& deck, const EclipseState& eclState)
57  {
58  const auto& pvdoTables = eclState.getTableManager().getPvdoTables();
59  const auto& densityKeyword = deck.getKeyword("DENSITY");
60 
61  assert(pvdoTables.size() == densityKeyword.size());
62 
63  size_t numRegions = pvdoTables.size();
64  setNumRegions(numRegions);
65 
66  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
67  Scalar rhoRefO = densityKeyword.getRecord(regionIdx).getItem("OIL").getSIDouble(0);
68  Scalar rhoRefG = densityKeyword.getRecord(regionIdx).getItem("GAS").getSIDouble(0);
69  Scalar rhoRefW = densityKeyword.getRecord(regionIdx).getItem("WATER").getSIDouble(0);
70 
71  setReferenceDensities(regionIdx, rhoRefO, rhoRefG, rhoRefW);
72 
73  const auto& pvdoTable = pvdoTables.getTable<PvdoTable>(regionIdx);
74 
75  const auto& BColumn(pvdoTable.getFormationFactorColumn());
76  std::vector<Scalar> invBColumn(BColumn.size());
77  for (unsigned i = 0; i < invBColumn.size(); ++i)
78  invBColumn[i] = 1/BColumn[i];
79 
80  inverseOilB_[regionIdx].setXYArrays(pvdoTable.numRows(),
81  pvdoTable.getPressureColumn(),
82  invBColumn);
83  oilMu_[regionIdx].setXYArrays(pvdoTable.numRows(),
84  pvdoTable.getPressureColumn(),
85  pvdoTable.getViscosityColumn());
86  }
87 
88  initEnd();
89  }
90 #endif // HAVE_OPM_PARSER
91 
92  void setNumRegions(size_t numRegions)
93  {
94  oilReferenceDensity_.resize(numRegions);
95  inverseOilB_.resize(numRegions);
96  inverseOilBMu_.resize(numRegions);
97  oilMu_.resize(numRegions);
98  }
99 
103  void setReferenceDensities(unsigned regionIdx,
104  Scalar rhoRefOil,
105  Scalar /*rhoRefGas*/,
106  Scalar /*rhoRefWater*/)
107  {
108  oilReferenceDensity_[regionIdx] = rhoRefOil;
109  }
110 
121  void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedOneDFunction& invBo)
122  { inverseOilB_[regionIdx] = invBo; }
123 
129  void setOilViscosity(unsigned regionIdx, const TabulatedOneDFunction& muo)
130  { oilMu_[regionIdx] = muo; }
131 
135  void initEnd()
136  {
137  // calculate the final 2D functions which are used for interpolation.
138  size_t numRegions = oilMu_.size();
139  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
140  // calculate the table which stores the inverse of the product of the oil
141  // formation volume factor and the oil viscosity
142  const auto& oilMu = oilMu_[regionIdx];
143  const auto& invOilB = inverseOilB_[regionIdx];
144  assert(oilMu.numSamples() == invOilB.numSamples());
145 
146  std::vector<Scalar> invBMuColumn;
147  std::vector<Scalar> pressureColumn;
148  invBMuColumn.resize(oilMu.numSamples());
149  pressureColumn.resize(oilMu.numSamples());
150 
151  for (unsigned pIdx = 0; pIdx < oilMu.numSamples(); ++pIdx) {
152  pressureColumn[pIdx] = invOilB.xAt(pIdx);
153  invBMuColumn[pIdx] = invOilB.valueAt(pIdx)*1/oilMu.valueAt(pIdx);
154  }
155 
156  inverseOilBMu_[regionIdx].setXYArrays(pressureColumn.size(),
157  pressureColumn,
158  invBMuColumn);
159  }
160  }
161 
165  unsigned numRegions() const
166  { return inverseOilBMu_.size(); }
167 
171  template <class Evaluation>
172  Evaluation viscosity(unsigned regionIdx,
173  const Evaluation& temperature,
174  const Evaluation& pressure,
175  const Evaluation& /*Rs*/) const
176  { return saturatedViscosity(regionIdx, temperature, pressure); }
177 
181  template <class Evaluation>
182  Evaluation saturatedViscosity(unsigned regionIdx,
183  const Evaluation& /*temperature*/,
184  const Evaluation& pressure) const
185  {
186  const Evaluation& invBo = inverseOilB_[regionIdx].eval(pressure, /*extrapolate=*/true);
187  const Evaluation& invMuoBo = inverseOilBMu_[regionIdx].eval(pressure, /*extrapolate=*/true);
188 
189  return invBo/invMuoBo;
190  }
191 
195  template <class Evaluation>
196  Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
197  const Evaluation& /*temperature*/,
198  const Evaluation& pressure,
199  const Evaluation& /*Rs*/) const
200  { return inverseOilB_[regionIdx].eval(pressure, /*extrapolate=*/true); }
201 
207  template <class Evaluation>
208  Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
209  const Evaluation& /*temperature*/,
210  const Evaluation& pressure) const
211  { return inverseOilB_[regionIdx].eval(pressure, /*extrapolate=*/true); }
212 
216  template <class Evaluation>
217  Evaluation saturatedGasDissolutionFactor(unsigned /*regionIdx*/,
218  const Evaluation& /*temperature*/,
219  const Evaluation& /*pressure*/) const
220  { return 0.0; /* this is dead oil! */ }
221 
225  template <class Evaluation>
226  Evaluation saturatedGasDissolutionFactor(unsigned /*regionIdx*/,
227  const Evaluation& /*temperature*/,
228  const Evaluation& /*pressure*/,
229  const Evaluation& /*oilSaturation*/,
230  Scalar /*maxOilSaturation*/) const
231  { return 0.0; /* this is dead oil! */ }
232 
239  template <class Evaluation>
240  Evaluation saturationPressure(unsigned /*regionIdx*/,
241  const Evaluation& /*temperature*/,
242  const Evaluation& /*Rs*/) const
243  { return 0.0; /* this is dead oil, so there isn't any meaningful saturation pressure! */ }
244 
245  template <class Evaluation>
246  Evaluation saturatedGasMassFraction(unsigned /*regionIdx*/,
247  const Evaluation& /*temperature*/,
248  const Evaluation& /*pressure*/) const
249  { return 0.0; /* this is dead oil! */ }
250 
251  template <class Evaluation>
252  Evaluation saturatedGasMoleFraction(unsigned /*regionIdx*/,
253  const Evaluation& /*temperature*/,
254  const Evaluation& /*pressure*/) const
255  { return 0.0; /* this is dead oil! */ }
256 
257 private:
258  std::vector<Scalar> oilReferenceDensity_;
259  std::vector<TabulatedOneDFunction> inverseOilB_;
260  std::vector<TabulatedOneDFunction> oilMu_;
261  std::vector<TabulatedOneDFunction> inverseOilBMu_;
262 };
263 
264 } // namespace Opm
265 
266 #endif
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
Definition: DeadOilPvt.hpp:46
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: DeadOilPvt.hpp:165
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &) const
Returns the formation volume factor [-] of the fluid phase.
Definition: DeadOilPvt.hpp:196
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:172
Implements a scalar function that depends on two variables and which is sampled uniformly in the X di...
void setOilViscosity(unsigned regionIdx, const TabulatedOneDFunction &muo)
Initialize the viscosity of the oil phase.
Definition: DeadOilPvt.hpp:129
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:240
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: DeadOilPvt.hpp:135
Class implementing cubic splines.
Evaluation saturatedGasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: DeadOilPvt.hpp:217
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of saturated oil.
Definition: DeadOilPvt.hpp:208
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: DeadOilPvt.hpp:103
Evaluation saturatedGasDissolutionFactor(unsigned, const Evaluation &, const Evaluation &, const Evaluation &, Scalar) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: DeadOilPvt.hpp:226
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:182
Implements a linearly interpolated scalar function that depends on one variable.
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47
void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedOneDFunction &invBo)
Initialize the function for the oil formation volume factor.
Definition: DeadOilPvt.hpp:121