WaterPvtThermal.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_WATER_PVT_THERMAL_HPP
28 #define OPM_WATER_PVT_THERMAL_HPP
29 
31 
36 
37 #if HAVE_OPM_PARSER
38 #include <opm/parser/eclipse/Deck/Deck.hpp>
39 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
40 #include <opm/parser/eclipse/EclipseState/Tables/SimpleTable.hpp>
41 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
42 #endif
43 
44 namespace Opm {
45 template <class Scalar, bool enableThermal>
46 class WaterPvtMultiplexer;
47 
54 template <class Scalar>
56 {
58  typedef WaterPvtMultiplexer<Scalar, /*enableThermal=*/false> IsothermalPvt;
59 
60 public:
62  { delete isothermalPvt_; }
63 
64 #if HAVE_OPM_PARSER
65 
68  void initFromDeck(const Deck& deck,
69  const EclipseState& eclState)
70  {
72  // initialize the isothermal part
74  isothermalPvt_ = new IsothermalPvt;
75  isothermalPvt_->initFromDeck(deck, eclState);
76 
78  // initialize the thermal part
80  const auto& tables = eclState.getTableManager();
81 
82  enableThermalDensity_ = deck.hasKeyword("WATDENT");
83  enableThermalViscosity_ = deck.hasKeyword("VISCREF");
84 
85  unsigned numRegions = isothermalPvt_->numRegions();
86  setNumRegions(numRegions);
87 
88  if (enableThermalDensity_) {
89  const auto& watdentKeyword = deck.getKeyword("WATDENT");
90 
91  assert(watdentKeyword.size() == numRegions);
92  for (unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
93  const auto& watdentRecord = watdentKeyword.getRecord(regionIdx);
94 
95  watdentRefTemp_[regionIdx] = watdentRecord.getItem("REFERENCE_TEMPERATURE").getSIDouble(0);
96  watdentCT1_[regionIdx] = watdentRecord.getItem("EXPANSION_COEFF_LINEAR").getSIDouble(0);
97  watdentCT2_[regionIdx] = watdentRecord.getItem("EXPANSION_COEFF_QUADRATIC").getSIDouble(0);
98  }
99  }
100 
101  if (enableThermalViscosity_) {
102  const auto& viscrefKeyword = deck.getKeyword("VISCREF");
103 
104  const auto& watvisctTables = tables.getWatvisctTables();
105 
106  assert(watvisctTables.size() == numRegions);
107  assert(viscrefKeyword.size() == numRegions);
108 
109  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
110  const auto& T = watvisctTables[regionIdx].getColumn("Temperature").vectorCopy();
111  const auto& mu = watvisctTables[regionIdx].getColumn("Viscosity").vectorCopy();
112  watvisctCurves_[regionIdx].setXYContainers(T, mu);
113 
114  const auto& viscrefRecord = viscrefKeyword.getRecord(regionIdx);
115  viscrefPress_[regionIdx] = viscrefRecord.getItem("REFERENCE_PRESSURE").getSIDouble(0);
116  }
117  }
118  }
119 #endif // HAVE_OPM_PARSER
120 
124  void setNumRegions(size_t numRegions)
125  {
126  pvtwRefPress_.resize(numRegions);
127  pvtwRefB_.resize(numRegions);
128  pvtwCompressibility_.resize(numRegions);
129  pvtwViscosity_.resize(numRegions);
130  pvtwViscosibility_.resize(numRegions);
131  viscrefPress_.resize(numRegions);
132  watvisctCurves_.resize(numRegions);
133  watdentRefTemp_.resize(numRegions);
134  watdentCT1_.resize(numRegions);
135  watdentCT2_.resize(numRegions);
136  }
137 
141  void initEnd()
142  { }
143 
147  bool enableThermalDensity() const
148  { return enableThermalDensity_; }
149 
154  { return enableThermalViscosity_; }
155 
156  size_t numRegions() const
157  { return pvtwRefPress_.size(); }
158 
162  template <class Evaluation>
163  Evaluation viscosity(unsigned regionIdx,
164  const Evaluation& temperature,
165  const Evaluation& pressure) const
166  {
167  const auto& isothermalMu = isothermalPvt_->viscosity(regionIdx, temperature, pressure);
168  if (!enableThermalViscosity())
169  return isothermalMu;
170 
171  Scalar x = -pvtwViscosibility_[regionIdx]*(viscrefPress_[regionIdx] - pvtwRefPress_[regionIdx]);
172  Scalar muRef = pvtwViscosity_[regionIdx]/(1.0 + x + 0.5*x*x);
173 
174  // compute the viscosity deviation due to temperature
175  const auto& muWatvisct = watvisctCurves_[regionIdx].eval(temperature);
176  return isothermalMu * muWatvisct/muRef;
177  }
178 
182  template <class Evaluation>
183  Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
184  const Evaluation& temperature,
185  const Evaluation& pressure) const
186  {
187  if (!enableThermalDensity())
188  return isothermalPvt_->inverseFormationVolumeFactor(regionIdx, temperature, pressure);
189 
190  Scalar BwRef = pvtwRefB_[regionIdx];
191  Scalar TRef = watdentRefTemp_[regionIdx];
192  const Evaluation& X = pvtwCompressibility_[regionIdx]*(pressure - pvtwRefPress_[regionIdx]);
193  Scalar cT1 = watdentCT1_[regionIdx];
194  Scalar cT2 = watdentCT2_[regionIdx];
195  const Evaluation& Y = temperature - TRef;
196 
197  return ((1 - X)*(1 + cT1*Y + cT2*Y*Y))/BwRef;
198  }
199 
200 private:
201  IsothermalPvt* isothermalPvt_;
202 
203  // The PVT properties needed for temperature dependence. We need to store one
204  // value per PVT region.
205  std::vector<Scalar> viscrefPress_;
206 
207  std::vector<Scalar> watdentRefTemp_;
208  std::vector<Scalar> watdentCT1_;
209  std::vector<Scalar> watdentCT2_;
210 
211  std::vector<Scalar> pvtwRefPress_;
212  std::vector<Scalar> pvtwRefB_;
213  std::vector<Scalar> pvtwCompressibility_;
214  std::vector<Scalar> pvtwViscosity_;
215  std::vector<Scalar> pvtwViscosibility_;
216 
217  std::vector<TabulatedOneDFunction> watvisctCurves_;
218 
219  bool enableThermalDensity_;
220  bool enableThermalViscosity_;
221 };
222 
223 } // namespace Opm
224 
225 #endif
This class implements temperature dependence of the PVT properties of water.
Definition: WaterPvtThermal.hpp:55
void initEnd()
Finish initializing the thermal part of the water phase PVT properties.
Definition: WaterPvtThermal.hpp:141
Implements a scalar function that depends on two variables and which is sampled uniformly in the X di...
void setNumRegions(size_t numRegions)
Set the number of PVT-regions considered by this object.
Definition: WaterPvtThermal.hpp:124
Definition: Air_Mesitylene.hpp:33
Class implementing cubic splines.
bool enableThermalViscosity() const
Returns true iff the viscosity of the water phase is temperature dependent.
Definition: WaterPvtThermal.hpp:153
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WaterPvtThermal.hpp:163
bool enableThermalDensity() const
Returns true iff the density of the water phase is temperature dependent.
Definition: WaterPvtThermal.hpp:147
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WaterPvtMultiplexer.hpp:123
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition: WaterPvtMultiplexer.hpp:55
Implements a linearly interpolated scalar function that depends on one variable.
A central place for various physical constants occuring in some equations.
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WaterPvtThermal.hpp:183
This file provides a wrapper around the "final" C++-2011 statement.
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: WaterPvtMultiplexer.hpp:116
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WaterPvtMultiplexer.hpp:132