SolventPvt.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_SOLVENT_PVT_HPP
28 #define OPM_SOLVENT_PVT_HPP
29 
31 
33 
34 #if HAVE_OPM_PARSER
35 #include <opm/parser/eclipse/Deck/Deck.hpp>
36 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
37 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
38 #include <opm/parser/eclipse/EclipseState/Tables/PvdsTable.hpp>
39 #include <opm/parser/eclipse/Deck/DeckKeyword.hpp>
40 #include <opm/parser/eclipse/Deck/DeckRecord.hpp>
41 #endif
42 
43 #include <vector>
44 
45 namespace Opm {
50 template <class Scalar>
52 {
54  typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
55 
56 public:
57 #if HAVE_OPM_PARSER
58 
63  void initFromDeck(const Deck& deck, const EclipseState& eclState)
64  {
65  const auto& pvdsTables = eclState.getTableManager().getPvdsTables();
66  const auto& sdensityKeyword = deck.getKeyword("SDENSITY");
67 
68  assert(pvdsTables.size() == sdensityKeyword.size());
69 
70  size_t numRegions = pvdsTables.size();
71  setNumRegions(numRegions);
72 
73  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
74  Scalar rhoRefS = sdensityKeyword.getRecord(regionIdx).getItem("SOLVENT_DENSITY").getSIDouble(0);
75 
76  setReferenceDensity(regionIdx, rhoRefS);
77 
78  const auto& pvdsTable = pvdsTables.getTable<PvdsTable>(regionIdx);
79 
80  // say 99.97% of all time: "premature optimization is the root of all
81  // evil". Eclipse does this "optimization" for apparently no good reason!
82  std::vector<Scalar> invB(pvdsTable.numRows());
83  const auto& Bg = pvdsTable.getFormationFactorColumn();
84  for (unsigned i = 0; i < Bg.size(); ++ i) {
85  invB[i] = 1.0/Bg[i];
86  }
87 
88  size_t numSamples = invB.size();
89  inverseSolventB_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), invB);
90  solventMu_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), pvdsTable.getViscosityColumn());
91  }
92 
93  initEnd();
94  }
95 #endif
96 
97  void setNumRegions(size_t numRegions)
98  {
99  solventReferenceDensity_.resize(numRegions);
100  inverseSolventB_.resize(numRegions);
101  inverseSolventBMu_.resize(numRegions);
102  solventMu_.resize(numRegions);
103  }
104 
108  void setReferenceDensity(unsigned regionIdx, Scalar rhoRefSolvent)
109  { solventReferenceDensity_[regionIdx] = rhoRefSolvent; }
110 
116  void setSolventViscosity(unsigned regionIdx, const TabulatedOneDFunction& mug)
117  { solventMu_[regionIdx] = mug; }
118 
124  void setSolventFormationVolumeFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
125  {
126  SamplingPoints tmp(samplePoints);
127  auto it = tmp.begin();
128  const auto& endIt = tmp.end();
129  for (; it != endIt; ++ it)
130  std::get<1>(*it) = 1.0/std::get<1>(*it);
131 
132  inverseSolventB_[regionIdx].setContainerOfTuples(tmp);
133  assert(inverseSolventB_[regionIdx].monotonic());
134  }
135 
139  void initEnd()
140  {
141  // calculate the final 2D functions which are used for interpolation.
142  size_t numRegions = solventMu_.size();
143  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
144  // calculate the table which stores the inverse of the product of the solvent
145  // formation volume factor and its viscosity
146  const auto& solventMu = solventMu_[regionIdx];
147  const auto& invSolventB = inverseSolventB_[regionIdx];
148  assert(solventMu.numSamples() == invSolventB.numSamples());
149 
150  std::vector<Scalar> pressureValues(solventMu.numSamples());
151  std::vector<Scalar> invSolventBMuValues(solventMu.numSamples());
152  for (unsigned pIdx = 0; pIdx < solventMu.numSamples(); ++pIdx) {
153  pressureValues[pIdx] = invSolventB.xAt(pIdx);
154  invSolventBMuValues[pIdx] = invSolventB.valueAt(pIdx) * (1.0/solventMu.valueAt(pIdx));
155  }
156 
157  inverseSolventBMu_[regionIdx].setXYContainers(pressureValues, invSolventBMuValues);
158  }
159  }
160 
164  unsigned numRegions() const
165  { return solventReferenceDensity_.size(); }
166 
170  template <class Evaluation>
171  Evaluation viscosity(unsigned regionIdx,
172  const Evaluation& temperature OPM_UNUSED,
173  const Evaluation& pressure) const
174  {
175  const Evaluation& invBg = inverseSolventB_[regionIdx].eval(pressure, /*extrapolate=*/true);
176  const Evaluation& invMugBg = inverseSolventBMu_[regionIdx].eval(pressure, /*extrapolate=*/true);
177 
178  return invBg/invMugBg;
179  }
180 
184  Scalar referenceDensity(unsigned regionIdx) const
185  { return solventReferenceDensity_[regionIdx]; }
186 
190  template <class Evaluation>
191  Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
192  const Evaluation& temperature OPM_UNUSED,
193  const Evaluation& pressure) const
194  { return inverseSolventB_[regionIdx].eval(pressure, /*extrapolate=*/true); }
195 
196 private:
197  std::vector<Scalar> solventReferenceDensity_;
198  std::vector<TabulatedOneDFunction> inverseSolventB_;
199  std::vector<TabulatedOneDFunction> solventMu_;
200  std::vector<TabulatedOneDFunction> inverseSolventBMu_;
201 };
202 
203 } // namespace Opm
204 
205 #endif
Scalar referenceDensity(unsigned regionIdx) const
Return the reference density the solvent phase for a given PVT region.
Definition: SolventPvt.hpp:184
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature OPM_UNUSED, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: SolventPvt.hpp:191
This class represents the Pressure-Volume-Temperature relations of the "second" gas phase in the of E...
Definition: SolventPvt.hpp:51
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: SolventPvt.hpp:164
void setReferenceDensity(unsigned regionIdx, Scalar rhoRefSolvent)
Initialize the reference density of the solvent gas for a given PVT region.
Definition: SolventPvt.hpp:108
Definition: Air_Mesitylene.hpp:33
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: SolventPvt.hpp:139
void setSolventViscosity(unsigned regionIdx, const TabulatedOneDFunction &mug)
Initialize the viscosity of the solvent gas phase.
Definition: SolventPvt.hpp:116
Implements a linearly interpolated scalar function that depends on one variable.
A central place for various physical constants occuring in some equations.
void setSolventFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the formation volume factor of solvent gas.
Definition: SolventPvt.hpp:124
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature OPM_UNUSED, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: SolventPvt.hpp:171