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