My Project
ConstantCompressibilityBrinePvt.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_CONSTANT_COMPRESSIBILITY_BRINE_PVT_HPP
28 #define OPM_CONSTANT_COMPRESSIBILITY_BRINE_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/PvtwsaltTable.hpp>
36 #endif
37 
38 #include <vector>
39 
40 namespace Opm {
41 template <class Scalar, bool enableThermal, bool enableBrine>
42 class WaterPvtMultiplexer;
47 template <class Scalar>
49 {
50 public:
52 
54  ConstantCompressibilityBrinePvt(const std::vector<Scalar>& waterReferenceDensity,
55  const std::vector<Scalar>& referencePressure,
56  const std::vector<TabulatedFunction> formationVolumeTables,
57  const std::vector<TabulatedFunction> compressibilityTables,
58  const std::vector<TabulatedFunction> viscosityTables,
59  const std::vector<TabulatedFunction> viscosibilityTables)
60  : waterReferenceDensity_(waterReferenceDensity)
61  , referencePressure_(referencePressure)
62  , formationVolumeTables_(formationVolumeTables)
63  , compressibilityTables_(compressibilityTables)
64  , viscosityTables_(viscosityTables)
65  , viscosibilityTables_(viscosibilityTables)
66  { }
67 #if HAVE_ECL_INPUT
72  void initFromState(const EclipseState& eclState, const Schedule&)
73  {
74  const auto& tableManager = eclState.getTableManager();
75  size_t numRegions = tableManager.getTabdims().getNumPVTTables();
76  const auto& densityTable = tableManager.getDensityTable();
77 
78  formationVolumeTables_.resize(numRegions);
79  compressibilityTables_.resize(numRegions);
80  viscosityTables_.resize(numRegions);
81  viscosibilityTables_.resize(numRegions);
82  referencePressure_.resize(numRegions);
83 
84  const auto& pvtwsaltTables = tableManager.getPvtwSaltTables();
85  if(!pvtwsaltTables.empty()){
86  assert(numRegions == pvtwsaltTables.size());
87  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
88  const auto& pvtwsaltTable = pvtwsaltTables[regionIdx];
89  const auto& c = pvtwsaltTable.getSaltConcentrationColumn();
90 
91  const auto& B = pvtwsaltTable.getFormationVolumeFactorColumn();
92  formationVolumeTables_[regionIdx].setXYContainers(c, B);
93 
94  const auto& compressibility = pvtwsaltTable.getCompressibilityColumn();
95  compressibilityTables_[regionIdx].setXYContainers(c, compressibility);
96 
97  const auto& viscositytable = pvtwsaltTable.getViscosityColumn();
98  viscosityTables_[regionIdx].setXYContainers(c, viscositytable);
99 
100  const auto& viscosibility = pvtwsaltTable.getViscosibilityColumn();
101  viscosibilityTables_[regionIdx].setXYContainers(c, viscosibility);
102  referencePressure_[regionIdx] = pvtwsaltTable.getReferencePressureValue();
103  }
104  }
105  else {
106  throw std::runtime_error("PVTWSALT must be specified in BRINE runs\n");
107  }
108 
109 
110  size_t numPvtwRegions = numRegions;
111  setNumRegions(numPvtwRegions);
112 
113  for (unsigned regionIdx = 0; regionIdx < numPvtwRegions; ++ regionIdx) {
114 
115  waterReferenceDensity_[regionIdx] = densityTable[regionIdx].water;
116  }
117 
118  initEnd();
119  }
120 #endif
121 
122  void setNumRegions(size_t numRegions)
123  {
124  waterReferenceDensity_.resize(numRegions);
125 
126  for (unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
127  setReferenceDensities(regionIdx, 650.0, 1.0, 1000.0);
128  }
129  }
130 
134  void setReferenceDensities(unsigned regionIdx,
135  Scalar /*rhoRefOil*/,
136  Scalar /*rhoRefGas*/,
137  Scalar rhoRefWater)
138  { waterReferenceDensity_[regionIdx] = rhoRefWater; }
139 
143  void initEnd()
144  { }
145 
146 
150  unsigned numRegions() const
151  { return waterReferenceDensity_.size(); }
152 
156  template <class Evaluation>
157  Evaluation internalEnergy(unsigned,
158  const Evaluation&,
159  const Evaluation&,
160  const Evaluation&) const
161  {
162  throw std::runtime_error("Requested the enthalpy of water but the thermal option is not enabled");
163  }
164 
168  template <class Evaluation>
169  Evaluation viscosity(unsigned regionIdx,
170  const Evaluation& temperature,
171  const Evaluation& pressure,
172  const Evaluation& saltconcentration) const
173  {
174  // cf. ECLiPSE 2013.2 technical description, p. 114
175  Scalar pRef = referencePressure_[regionIdx];
176  const Evaluation C = compressibilityTables_[regionIdx].eval(saltconcentration, /*extrapolate=*/true);
177  const Evaluation Cv = viscosibilityTables_[regionIdx].eval(saltconcentration, /*extrapolate=*/true);
178  const Evaluation BwRef = formationVolumeTables_[regionIdx].eval(saltconcentration, /*extrapolate=*/true);
179  const Evaluation Y = (C-Cv)* (pressure - pRef);
180  Evaluation MuwRef = viscosityTables_[regionIdx].eval(saltconcentration, /*extrapolate=*/true);
181 
182  const Evaluation& bw = inverseFormationVolumeFactor(regionIdx, temperature, pressure, saltconcentration);
183 
184  return MuwRef*BwRef*bw/(1 + Y*(1 + Y/2));
185  }
186 
190  template <class Evaluation>
191  Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
192  const Evaluation& /*temperature*/,
193  const Evaluation& pressure,
194  const Evaluation& saltconcentration) const
195  {
196  Scalar pRef = referencePressure_[regionIdx];
197 
198  const Evaluation BwRef = formationVolumeTables_[regionIdx].eval(saltconcentration, /*extrapolate=*/true);
199  const Evaluation C = compressibilityTables_[regionIdx].eval(saltconcentration, /*extrapolate=*/true);
200  const Evaluation X = C * (pressure - pRef);
201 
202  return (1.0 + X*(1.0 + X/2.0))/BwRef;
203 
204  }
205 
206  const Scalar waterReferenceDensity(unsigned regionIdx) const
207  { return waterReferenceDensity_[regionIdx]; }
208 
209  const std::vector<Scalar>& referencePressure() const
210  { return referencePressure_; }
211 
212  const std::vector<TabulatedFunction>& formationVolumeTables() const
213  { return formationVolumeTables_; }
214 
215  const std::vector<TabulatedFunction>& compressibilityTables() const
216  { return compressibilityTables_; }
217 
218  const std::vector<TabulatedFunction>& viscosityTables() const
219  { return viscosityTables_; }
220 
221  const std::vector<TabulatedFunction>& viscosibilityTables() const
222  { return viscosibilityTables_; }
223 
224  bool operator==(const ConstantCompressibilityBrinePvt<Scalar>& data) const
225  {
226  return this->waterReferenceDensity_ == data.waterReferenceDensity_ &&
227  this->referencePressure() == data.referencePressure() &&
228  this->formationVolumeTables() == data.formationVolumeTables() &&
229  this->compressibilityTables() == data.compressibilityTables() &&
230  this->viscosityTables() == data.viscosityTables() &&
231  this->viscosibilityTables() == data.viscosibilityTables();
232  }
233 
234 private:
235  std::vector<Scalar> waterReferenceDensity_;
236  std::vector<Scalar> referencePressure_;
237  std::vector<TabulatedFunction> formationVolumeTables_;
238  std::vector<TabulatedFunction> compressibilityTables_;
239  std::vector<TabulatedFunction> viscosityTables_;
240  std::vector<TabulatedFunction> viscosibilityTables_;
241 };
242 
243 } // namespace Opm
244 
245 #endif
Implements a linearly interpolated scalar function that depends on one variable.
This class represents the Pressure-Volume-Temperature relations of the gas phase without vaporized oi...
Definition: ConstantCompressibilityBrinePvt.hpp:49
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: ConstantCompressibilityBrinePvt.hpp:150
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &saltconcentration) const
Returns the formation volume factor [-] of the fluid phase.
Definition: ConstantCompressibilityBrinePvt.hpp:191
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of water given a set of parameters.
Definition: ConstantCompressibilityBrinePvt.hpp:157
void setReferenceDensities(unsigned regionIdx, Scalar, Scalar, Scalar rhoRefWater)
Set the water reference density [kg / m^3].
Definition: ConstantCompressibilityBrinePvt.hpp:134
void initEnd()
Finish initializing the water phase PVT properties.
Definition: ConstantCompressibilityBrinePvt.hpp:143
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &saltconcentration) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: ConstantCompressibilityBrinePvt.hpp:169
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:48