My Project
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_ECL_INPUT
35 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
36 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
37 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
38 #include <opm/parser/eclipse/EclipseState/Tables/PvdsTable.hpp>
39 #endif
40 
41 #include <vector>
42 
43 namespace Opm {
48 template <class Scalar>
50 {
51  typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
52 
53 public:
55 
56  explicit SolventPvt() = default;
57  SolventPvt(const std::vector<Scalar>& solventReferenceDensity,
58  const std::vector<TabulatedOneDFunction>& inverseSolventB,
59  const std::vector<TabulatedOneDFunction>& solventMu,
60  const std::vector<TabulatedOneDFunction>& inverseSolventBMu)
61  : solventReferenceDensity_(solventReferenceDensity)
62  , inverseSolventB_(inverseSolventB)
63  , solventMu_(solventMu)
64  , inverseSolventBMu_(inverseSolventBMu)
65  {
66  }
67 
68 #if HAVE_ECL_INPUT
74  void initFromState(const EclipseState& eclState, const Schedule&)
75  {
76  const auto& pvdsTables = eclState.getTableManager().getPvdsTables();
77  const auto& sdensityTables = eclState.getTableManager().getSolventDensityTables();
78 
79  assert(pvdsTables.size() == sdensityTables.size());
80 
81  size_t numRegions = pvdsTables.size();
82  setNumRegions(numRegions);
83 
84  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
85  Scalar rhoRefS = sdensityTables[regionIdx].getSolventDensityColumn().front();
86 
87  setReferenceDensity(regionIdx, rhoRefS);
88 
89  const auto& pvdsTable = pvdsTables.getTable<PvdsTable>(regionIdx);
90 
91  // say 99.97% of all time: "premature optimization is the root of all
92  // evil". Eclipse does this "optimization" for apparently no good reason!
93  std::vector<Scalar> invB(pvdsTable.numRows());
94  const auto& Bg = pvdsTable.getFormationFactorColumn();
95  for (unsigned i = 0; i < Bg.size(); ++ i) {
96  invB[i] = 1.0/Bg[i];
97  }
98 
99  size_t numSamples = invB.size();
100  inverseSolventB_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), invB);
101  solventMu_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), pvdsTable.getViscosityColumn());
102  }
103 
104  initEnd();
105  }
106 #endif
107 
108  void setNumRegions(size_t numRegions)
109  {
110  solventReferenceDensity_.resize(numRegions);
111  inverseSolventB_.resize(numRegions);
112  inverseSolventBMu_.resize(numRegions);
113  solventMu_.resize(numRegions);
114  }
115 
119  void setReferenceDensity(unsigned regionIdx, Scalar rhoRefSolvent)
120  { solventReferenceDensity_[regionIdx] = rhoRefSolvent; }
121 
127  void setSolventViscosity(unsigned regionIdx, const TabulatedOneDFunction& mug)
128  { solventMu_[regionIdx] = mug; }
129 
135  void setSolventFormationVolumeFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
136  {
137  SamplingPoints tmp(samplePoints);
138  auto it = tmp.begin();
139  const auto& endIt = tmp.end();
140  for (; it != endIt; ++ it)
141  std::get<1>(*it) = 1.0/std::get<1>(*it);
142 
143  inverseSolventB_[regionIdx].setContainerOfTuples(tmp);
144  assert(inverseSolventB_[regionIdx].monotonic());
145  }
146 
150  void initEnd()
151  {
152  // calculate the final 2D functions which are used for interpolation.
153  size_t numRegions = solventMu_.size();
154  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
155  // calculate the table which stores the inverse of the product of the solvent
156  // formation volume factor and its viscosity
157  const auto& solventMu = solventMu_[regionIdx];
158  const auto& invSolventB = inverseSolventB_[regionIdx];
159  assert(solventMu.numSamples() == invSolventB.numSamples());
160 
161  std::vector<Scalar> pressureValues(solventMu.numSamples());
162  std::vector<Scalar> invSolventBMuValues(solventMu.numSamples());
163  for (unsigned pIdx = 0; pIdx < solventMu.numSamples(); ++pIdx) {
164  pressureValues[pIdx] = invSolventB.xAt(pIdx);
165  invSolventBMuValues[pIdx] = invSolventB.valueAt(pIdx) * (1.0/solventMu.valueAt(pIdx));
166  }
167 
168  inverseSolventBMu_[regionIdx].setXYContainers(pressureValues, invSolventBMuValues);
169  }
170  }
171 
175  unsigned numRegions() const
176  { return solventReferenceDensity_.size(); }
177 
181  template <class Evaluation>
182  Evaluation viscosity(unsigned regionIdx,
183  const Evaluation&,
184  const Evaluation& pressure) const
185  {
186  const Evaluation& invBg = inverseSolventB_[regionIdx].eval(pressure, /*extrapolate=*/true);
187  const Evaluation& invMugBg = inverseSolventBMu_[regionIdx].eval(pressure, /*extrapolate=*/true);
188 
189  return invBg/invMugBg;
190  }
191 
192  template <class Evaluation>
193  Evaluation diffusionCoefficient(const Evaluation& /*temperature*/,
194  const Evaluation& /*pressure*/,
195  unsigned /*compIdx*/) const
196  {
197  throw std::runtime_error("Not implemented: The PVT model does not provide a diffusionCoefficient()");
198  }
199 
203  Scalar referenceDensity(unsigned regionIdx) const
204  { return solventReferenceDensity_[regionIdx]; }
205 
209  template <class Evaluation>
210  Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
211  const Evaluation&,
212  const Evaluation& pressure) const
213  { return inverseSolventB_[regionIdx].eval(pressure, /*extrapolate=*/true); }
214 
215  const std::vector<Scalar>& solventReferenceDensity() const
216  { return solventReferenceDensity_; }
217 
218  const std::vector<TabulatedOneDFunction>& inverseSolventB() const
219  { return inverseSolventB_; }
220 
221  const std::vector<TabulatedOneDFunction>& solventMu() const
222  { return solventMu_; }
223 
224  const std::vector<TabulatedOneDFunction>& inverseSolventBMu() const
225  { return inverseSolventBMu_; }
226 
227  bool operator==(const SolventPvt<Scalar>& data) const
228  {
229  return solventReferenceDensity_ == data.solventReferenceDensity_ &&
230  inverseSolventB_ == data.inverseSolventB_ &&
231  solventMu_ == data.solventMu_ &&
232  inverseSolventBMu_ == data.inverseSolventBMu_;
233  }
234 
235 private:
236  std::vector<Scalar> solventReferenceDensity_;
237  std::vector<TabulatedOneDFunction> inverseSolventB_;
238  std::vector<TabulatedOneDFunction> solventMu_;
239  std::vector<TabulatedOneDFunction> inverseSolventBMu_;
240 };
241 
242 } // namespace Opm
243 
244 #endif
A central place for various physical constants occuring in some equations.
Implements a linearly interpolated scalar function that depends on one variable.
This class represents the Pressure-Volume-Temperature relations of the "second" gas phase in the of E...
Definition: SolventPvt.hpp:50
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: SolventPvt.hpp:182
void setReferenceDensity(unsigned regionIdx, Scalar rhoRefSolvent)
Initialize the reference density of the solvent gas for a given PVT region.
Definition: SolventPvt.hpp:119
void setSolventViscosity(unsigned regionIdx, const TabulatedOneDFunction &mug)
Initialize the viscosity of the solvent gas phase.
Definition: SolventPvt.hpp:127
Scalar referenceDensity(unsigned regionIdx) const
Return the reference density the solvent phase for a given PVT region.
Definition: SolventPvt.hpp:203
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: SolventPvt.hpp:150
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: SolventPvt.hpp:175
void setSolventFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the formation volume factor of solvent gas.
Definition: SolventPvt.hpp:135
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: SolventPvt.hpp:210
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47