My Project
DryHumidGasPvt.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_DRY_HUMID_GAS_PVT_HPP
28 #define OPM_DRY_HUMID_GAS_PVT_HPP
29 
34 
35 #if HAVE_ECL_INPUT
36 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
37 #include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
38 
39 #endif
40 
41 #if HAVE_OPM_COMMON
42 #include <opm/common/OpmLog/OpmLog.hpp>
43 #endif
44 
45 namespace Opm {
46 
51 template <class Scalar>
53 {
54  using SamplingPoints = std::vector<std::pair<Scalar, Scalar>>;
55 
56 public:
59 
61  {
62  vapPar1_ = 0.0;
63  }
64 
65  DryHumidGasPvt(const std::vector<Scalar>& gasReferenceDensity,
66  const std::vector<Scalar>& waterReferenceDensity,
67  const std::vector<TabulatedTwoDFunction>& inverseGasB,
68  const std::vector<TabulatedOneDFunction>& inverseSaturatedGasB,
69  const std::vector<TabulatedTwoDFunction>& gasMu,
70  const std::vector<TabulatedTwoDFunction>& inverseGasBMu,
71  const std::vector<TabulatedOneDFunction>& inverseSaturatedGasBMu,
72  const std::vector<TabulatedOneDFunction>& saturatedWaterVaporizationFactorTable,
73  const std::vector<TabulatedOneDFunction>& saturationPressure,
74  Scalar vapPar1)
75  : gasReferenceDensity_(gasReferenceDensity)
76  , waterReferenceDensity_(waterReferenceDensity)
77  , inverseGasB_(inverseGasB)
78  , inverseSaturatedGasB_(inverseSaturatedGasB)
79  , gasMu_(gasMu)
80  , inverseGasBMu_(inverseGasBMu)
81  , inverseSaturatedGasBMu_(inverseSaturatedGasBMu)
82  , saturatedWaterVaporizationFactorTable_(saturatedWaterVaporizationFactorTable)
83  , saturationPressure_(saturationPressure)
84  , vapPar1_(vapPar1)
85  {
86  }
87 
88 
89 #if HAVE_ECL_INPUT
95  void initFromState(const EclipseState& eclState, const Schedule&)
96  {
97  const auto& pvtgwTables = eclState.getTableManager().getPvtgwTables();
98  const auto& densityTable = eclState.getTableManager().getDensityTable();
99 
100  assert(pvtgwTables.size() == densityTable.size());
101 
102  size_t numRegions = pvtgwTables.size();
103  setNumRegions(numRegions);
104 
105  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
106  Scalar rhoRefO = densityTable[regionIdx].oil;
107  Scalar rhoRefG = densityTable[regionIdx].gas;
108  Scalar rhoRefW = densityTable[regionIdx].water;
109 
110  setReferenceDensities(regionIdx, rhoRefO, rhoRefG, rhoRefW);
111  }
112 
113  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
114  const auto& pvtgwTable = pvtgwTables[regionIdx];
115 
116  const auto& saturatedTable = pvtgwTable.getSaturatedTable();
117  assert(saturatedTable.numRows() > 1);
118 
119  auto& gasMu = gasMu_[regionIdx];
120  auto& invGasB = inverseGasB_[regionIdx];
121  auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
122  auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
123  auto& waterVaporizationFac = saturatedWaterVaporizationFactorTable_[regionIdx];
124 
125  waterVaporizationFac.setXYArrays(saturatedTable.numRows(),
126  saturatedTable.getColumn("PG"),
127  saturatedTable.getColumn("RW"));
128 
129  std::vector<Scalar> invSatGasBArray;
130  std::vector<Scalar> invSatGasBMuArray;
131 
132  // extract the table for the gas dissolution and the oil formation volume factors
133  for (unsigned outerIdx = 0; outerIdx < saturatedTable.numRows(); ++ outerIdx) {
134  Scalar pg = saturatedTable.get("PG" , outerIdx);
135  Scalar B = saturatedTable.get("BG" , outerIdx);
136  Scalar mu = saturatedTable.get("MUG" , outerIdx);
137 
138  invGasB.appendXPos(pg);
139  gasMu.appendXPos(pg);
140 
141  invSatGasBArray.push_back(1.0/B);
142  invSatGasBMuArray.push_back(1.0/(mu*B));
143 
144  assert(invGasB.numX() == outerIdx + 1);
145  assert(gasMu.numX() == outerIdx + 1);
146 
147  const auto& underSaturatedTable = pvtgwTable.getUnderSaturatedTable(outerIdx);
148  size_t numRows = underSaturatedTable.numRows();
149  for (size_t innerIdx = 0; innerIdx < numRows; ++ innerIdx) {
150  Scalar Rw = underSaturatedTable.get("RW" , innerIdx);
151  Scalar Bg = underSaturatedTable.get("BG" , innerIdx);
152  Scalar mug = underSaturatedTable.get("MUG" , innerIdx);
153 
154  invGasB.appendSamplePoint(outerIdx, Rw, 1.0/Bg);
155  gasMu.appendSamplePoint(outerIdx, Rw, mug);
156  }
157  }
158 
159  {
160  std::vector<double> tmpPressure = saturatedTable.getColumn("PG").vectorCopy( );
161 
162  invSatGasB.setXYContainers(tmpPressure, invSatGasBArray);
163  invSatGasBMu.setXYContainers(tmpPressure, invSatGasBMuArray);
164  }
165 
166  // make sure to have at least two sample points per gas pressure value
167  for (unsigned xIdx = 0; xIdx < invGasB.numX(); ++xIdx) {
168  // a single sample point is definitely needed
169  assert(invGasB.numY(xIdx) > 0);
170 
171  // everything is fine if the current table has two or more sampling points
172  // for a given mole fraction
173  if (invGasB.numY(xIdx) > 1)
174  continue;
175 
176  // find the master table which will be used as a template to extend the
177  // current line. We define master table as the first table which has values
178  // for undersaturated gas...
179  size_t masterTableIdx = xIdx + 1;
180  for (; masterTableIdx < saturatedTable.numRows(); ++masterTableIdx)
181  {
182  if (pvtgwTable.getUnderSaturatedTable(masterTableIdx).numRows() > 1)
183  break;
184  }
185 
186  if (masterTableIdx >= saturatedTable.numRows())
187  throw std::runtime_error("PVTGW tables are invalid: The last table must exhibit at least one "
188  "entry for undersaturated gas!");
189 
190 
191  // extend the current table using the master table.
192  extendPvtgwTable_(regionIdx,
193  xIdx,
194  pvtgwTable.getUnderSaturatedTable(xIdx),
195  pvtgwTable.getUnderSaturatedTable(masterTableIdx));
196  }
197  }
198 
199  vapPar1_ = 0.0;
200 
201  initEnd();
202  }
203 
204 private:
205  void extendPvtgwTable_(unsigned regionIdx,
206  unsigned xIdx,
207  const SimpleTable& curTable,
208  const SimpleTable& masterTable)
209  {
210  std::vector<double> RwArray = curTable.getColumn("RW").vectorCopy();
211  std::vector<double> gasBArray = curTable.getColumn("BG").vectorCopy();
212  std::vector<double> gasMuArray = curTable.getColumn("MUG").vectorCopy();
213 
214  auto& invGasB = inverseGasB_[regionIdx];
215  auto& gasMu = gasMu_[regionIdx];
216 
217  for (size_t newRowIdx = 1; newRowIdx < masterTable.numRows(); ++ newRowIdx) {
218  const auto& RWColumn = masterTable.getColumn("RW");
219  const auto& BGColumn = masterTable.getColumn("BG");
220  const auto& viscosityColumn = masterTable.getColumn("MUG");
221 
222  // compute the vaporized water factor Rw for the new entry
223  Scalar diffRw = RWColumn[newRowIdx] - RWColumn[newRowIdx - 1];
224  Scalar newRw = RwArray.back() + diffRw;
225 
226  // calculate the compressibility of the master table
227  Scalar B1 = BGColumn[newRowIdx];
228  Scalar B2 = BGColumn[newRowIdx - 1];
229  Scalar x = (B1 - B2)/( (B1 + B2)/2.0 );
230 
231  // calculate the gas formation volume factor which exhibits the same
232  // "compressibility" for the new value of Rw
233  Scalar newBg = gasBArray.back()*(1.0 + x/2.0)/(1.0 - x/2.0);
234 
235  // calculate the "viscosibility" of the master table
236  Scalar mu1 = viscosityColumn[newRowIdx];
237  Scalar mu2 = viscosityColumn[newRowIdx - 1];
238  Scalar xMu = (mu1 - mu2)/( (mu1 + mu2)/2.0 );
239 
240  // calculate the viscosity which exhibits the same
241  // "viscosibility" for the new Rw value
242  Scalar newMug = gasMuArray.back()*(1.0 + xMu/2)/(1.0 - xMu/2.0);
243 
244  // append the new values to the arrays which we use to compute the additional
245  // values ...
246  RwArray.push_back(newRw);
247  gasBArray.push_back(newBg);
248  gasMuArray.push_back(newMug);
249 
250  // ... and register them with the internal table objects
251  invGasB.appendSamplePoint(xIdx, newRw, 1.0/newBg);
252  gasMu.appendSamplePoint(xIdx, newRw, newMug);
253  }
254  }
255 
256 public:
257 #endif // HAVE_ECL_INPUT
258 
259  void setNumRegions(size_t numRegions)
260  {
261  waterReferenceDensity_.resize(numRegions);
262  gasReferenceDensity_.resize(numRegions);
263  inverseGasB_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::RightExtreme});
264  inverseGasBMu_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::RightExtreme});
265  inverseSaturatedGasB_.resize(numRegions);
266  inverseSaturatedGasBMu_.resize(numRegions);
267  gasMu_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::RightExtreme});
268  saturatedWaterVaporizationFactorTable_.resize(numRegions);
269  saturationPressure_.resize(numRegions);
270  }
271 
275  void setReferenceDensities(unsigned regionIdx,
276  Scalar /*rhoRefOil*/,
277  Scalar rhoRefGas,
278  Scalar rhoRefWater)
279  {
280  waterReferenceDensity_[regionIdx] = rhoRefWater;
281  gasReferenceDensity_[regionIdx] = rhoRefGas;
282  }
283 
289  void setSaturatedGasWaterVaporizationFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
290  { saturatedWaterVaporizationFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
291 
304  void setInverseGasFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction& invBg)
305  { inverseGasB_[regionIdx] = invBg; }
306 
312  void setGasViscosity(unsigned regionIdx, const TabulatedTwoDFunction& mug)
313  { gasMu_[regionIdx] = mug; }
314 
322  void setSaturatedGasViscosity(unsigned regionIdx, const SamplingPoints& samplePoints )
323  {
324  auto& waterVaporizationFac = saturatedWaterVaporizationFactorTable_[regionIdx];
325 
326  constexpr const Scalar RwMin = 0.0;
327  Scalar RwMax = waterVaporizationFac.eval(saturatedWaterVaporizationFactorTable_[regionIdx].xMax(), /*extrapolate=*/true);
328 
329  Scalar poMin = samplePoints.front().first;
330  Scalar poMax = samplePoints.back().first;
331 
332  constexpr const size_t nRw = 20;
333  size_t nP = samplePoints.size()*2;
334 
335  TabulatedOneDFunction mugTable;
336  mugTable.setContainerOfTuples(samplePoints);
337 
338  // calculate a table of estimated densities depending on pressure and gas mass
339  // fraction
340  for (size_t RwIdx = 0; RwIdx < nRw; ++RwIdx) {
341  Scalar Rw = RwMin + (RwMax - RwMin)*RwIdx/nRw;
342 
343  gasMu_[regionIdx].appendXPos(Rw);
344 
345  for (size_t pIdx = 0; pIdx < nP; ++pIdx) {
346  Scalar pg = poMin + (poMax - poMin)*pIdx/nP;
347  Scalar mug = mugTable.eval(pg, /*extrapolate=*/true);
348 
349  gasMu_[regionIdx].appendSamplePoint(RwIdx, pg, mug);
350  }
351  }
352  }
353 
357  void initEnd()
358  {
359  // calculate the final 2D functions which are used for interpolation.
360  size_t numRegions = gasMu_.size();
361  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
362  // calculate the table which stores the inverse of the product of the gas
363  // formation volume factor and the gas viscosity
364  const auto& gasMu = gasMu_[regionIdx];
365  const auto& invGasB = inverseGasB_[regionIdx];
366  assert(gasMu.numX() == invGasB.numX());
367 
368  auto& invGasBMu = inverseGasBMu_[regionIdx];
369  auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
370  auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
371 
372  std::vector<Scalar> satPressuresArray;
373  std::vector<Scalar> invSatGasBArray;
374  std::vector<Scalar> invSatGasBMuArray;
375  for (size_t pIdx = 0; pIdx < gasMu.numX(); ++pIdx) {
376  invGasBMu.appendXPos(gasMu.xAt(pIdx));
377 
378  assert(gasMu.numY(pIdx) == invGasB.numY(pIdx));
379 
380  size_t numRw = gasMu.numY(pIdx);
381  for (size_t RwIdx = 0; RwIdx < numRw; ++RwIdx)
382  invGasBMu.appendSamplePoint(pIdx,
383  gasMu.yAt(pIdx, RwIdx),
384  invGasB.valueAt(pIdx, RwIdx)
385  / gasMu.valueAt(pIdx, RwIdx));
386 
387  // the sampling points in UniformXTabulated2DFunction are always sorted
388  // in ascending order. Thus, the value for saturated gas is the last one
389  // (i.e., the one with the largest Rw value)
390  satPressuresArray.push_back(gasMu.xAt(pIdx));
391  invSatGasBArray.push_back(invGasB.valueAt(pIdx, numRw - 1));
392  invSatGasBMuArray.push_back(invGasBMu.valueAt(pIdx, numRw - 1));
393  }
394 
395  invSatGasB.setXYContainers(satPressuresArray, invSatGasBArray);
396  invSatGasBMu.setXYContainers(satPressuresArray, invSatGasBMuArray);
397 
398  updateSaturationPressure_(regionIdx);
399  }
400  }
401 
405  unsigned numRegions() const
406  { return gasReferenceDensity_.size(); }
407 
411  template <class Evaluation>
412  Evaluation internalEnergy(unsigned,
413  const Evaluation&,
414  const Evaluation&,
415  const Evaluation&) const
416  {
417  throw std::runtime_error("Requested the enthalpy of gas but the thermal option is not enabled");
418  }
419 
423  template <class Evaluation>
424  Evaluation viscosity(unsigned regionIdx,
425  const Evaluation& /*temperature*/,
426  const Evaluation& pressure,
427  const Evaluation& /*Rv*/,
428  const Evaluation& Rvw) const
429  {
430  const Evaluation& invBg = inverseGasB_[regionIdx].eval(pressure, Rvw, /*extrapolate=*/true);
431  const Evaluation& invMugBg = inverseGasBMu_[regionIdx].eval(pressure, Rvw, /*extrapolate=*/true);
432 
433  return invBg/invMugBg;
434  }
435 
439  template <class Evaluation>
440  Evaluation saturatedViscosity(unsigned regionIdx,
441  const Evaluation& /*temperature*/,
442  const Evaluation& pressure) const
443  {
444  const Evaluation& invBg = inverseSaturatedGasB_[regionIdx].eval(pressure, /*extrapolate=*/true);
445  const Evaluation& invMugBg = inverseSaturatedGasBMu_[regionIdx].eval(pressure, /*extrapolate=*/true);
446 
447  return invBg/invMugBg;
448  }
449 
453  template <class Evaluation>
454  Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
455  const Evaluation& /*temperature*/,
456  const Evaluation& pressure,
457  const Evaluation& /*Rv*/,
458  const Evaluation& Rvw) const
459  { return inverseGasB_[regionIdx].eval(pressure, Rvw, /*extrapolate=*/true); }
460 
464  template <class Evaluation>
465  Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
466  const Evaluation& /*temperature*/,
467  const Evaluation& pressure) const
468  { return inverseSaturatedGasB_[regionIdx].eval(pressure, /*extrapolate=*/true); }
469 
473  template <class Evaluation>
474  Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx,
475  const Evaluation& /*temperature*/,
476  const Evaluation& pressure) const
477  {
478  return saturatedWaterVaporizationFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
479  }
480 
484  template <class Evaluation>
485  Evaluation saturatedOilVaporizationFactor(unsigned /*regionIdx*/,
486  const Evaluation& /*temperature*/,
487  const Evaluation& /*pressure*/,
488  const Evaluation& /*oilSaturation*/,
489  const Evaluation& /*maxOilSaturation*/) const
490  { return 0.0; /* this is dry humid gas! */ }
491 
495  template <class Evaluation>
496  Evaluation saturatedOilVaporizationFactor(unsigned /*regionIdx*/,
497  const Evaluation& /*temperature*/,
498  const Evaluation& /*pressure*/) const
499  { return 0.0; /* this is dry humid gas! */ }
500 
508  template <class Evaluation>
509  Evaluation saturationPressure(unsigned regionIdx,
510  const Evaluation&,
511  const Evaluation& Rw) const
512  {
513  typedef MathToolbox<Evaluation> Toolbox;
514 
515  const auto& RwTable = saturatedWaterVaporizationFactorTable_[regionIdx];
516  const Scalar eps = std::numeric_limits<typename Toolbox::Scalar>::epsilon()*1e6;
517 
518  // use the tabulated saturation pressure function to get a pretty good initial value
519  Evaluation pSat = saturationPressure_[regionIdx].eval(Rw, /*extrapolate=*/true);
520 
521  // Newton method to do the remaining work. If the initial
522  // value is good, this should only take two to three
523  // iterations...
524  bool onProbation = false;
525  for (unsigned i = 0; i < 20; ++i) {
526  const Evaluation& f = RwTable.eval(pSat, /*extrapolate=*/true) - Rw;
527  const Evaluation& fPrime = RwTable.evalDerivative(pSat, /*extrapolate=*/true);
528 
529  // If the derivative is "zero" Newton will not converge,
530  // so simply return our initial guess.
531  if (std::abs(scalarValue(fPrime)) < 1.0e-30) {
532  return pSat;
533  }
534 
535  const Evaluation& delta = f/fPrime;
536 
537  pSat -= delta;
538 
539  if (pSat < 0.0) {
540  // if the pressure is lower than 0 Pascals, we set it back to 0. if this
541  // happens twice, we give up and just return 0 Pa...
542  if (onProbation)
543  return 0.0;
544 
545  onProbation = true;
546  pSat = 0.0;
547  }
548 
549  if (std::abs(scalarValue(delta)) < std::abs(scalarValue(pSat))*eps)
550  return pSat;
551  }
552 
553  std::stringstream errlog;
554  errlog << "Finding saturation pressure did not converge:"
555  << " pSat = " << pSat
556  << ", Rw = " << Rw;
557 #if HAVE_OPM_COMMON
558  OpmLog::debug("Wet gas saturation pressure", errlog.str());
559 #endif
560  throw NumericalIssue(errlog.str());
561  }
562 
563  template <class Evaluation>
564  Evaluation diffusionCoefficient(const Evaluation& /*temperature*/,
565  const Evaluation& /*pressure*/,
566  unsigned /*compIdx*/) const
567  {
568  throw std::runtime_error("Not implemented: The PVT model does not provide a diffusionCoefficient()");
569  }
570 
571  const Scalar gasReferenceDensity(unsigned regionIdx) const
572  { return gasReferenceDensity_[regionIdx]; }
573 
574  const Scalar waterReferenceDensity(unsigned regionIdx) const
575  { return waterReferenceDensity_[regionIdx]; }
576 
577  const std::vector<TabulatedTwoDFunction>& inverseGasB() const {
578  return inverseGasB_;
579  }
580 
581  const std::vector<TabulatedOneDFunction>& inverseSaturatedGasB() const {
582  return inverseSaturatedGasB_;
583  }
584 
585  const std::vector<TabulatedTwoDFunction>& gasMu() const {
586  return gasMu_;
587  }
588 
589  const std::vector<TabulatedTwoDFunction>& inverseGasBMu() const {
590  return inverseGasBMu_;
591  }
592 
593  const std::vector<TabulatedOneDFunction>& inverseSaturatedGasBMu() const {
594  return inverseSaturatedGasBMu_;
595  }
596 
597  const std::vector<TabulatedOneDFunction>& saturatedWaterVaporizationFactorTable() const {
598  return saturatedWaterVaporizationFactorTable_;
599  }
600 
601  const std::vector<TabulatedOneDFunction>& saturationPressure() const {
602  return saturationPressure_;
603  }
604 
605  Scalar vapPar1() const {
606  return vapPar1_;
607  }
608 
609  bool operator==(const DryHumidGasPvt<Scalar>& data) const
610  {
611  return this->gasReferenceDensity_ == data.gasReferenceDensity_ &&
612  this->waterReferenceDensity_ == data.waterReferenceDensity_ &&
613  this->inverseGasB() == data.inverseGasB() &&
614  this->inverseSaturatedGasB() == data.inverseSaturatedGasB() &&
615  this->gasMu() == data.gasMu() &&
616  this->inverseGasBMu() == data.inverseGasBMu() &&
617  this->inverseSaturatedGasBMu() == data.inverseSaturatedGasBMu() &&
618  this->saturatedWaterVaporizationFactorTable() == data.saturatedWaterVaporizationFactorTable() &&
619  this->saturationPressure() == data.saturationPressure() &&
620  this->vapPar1() == data.vapPar1();
621  }
622 
623 private:
624  void updateSaturationPressure_(unsigned regionIdx)
625  {
626  typedef std::pair<Scalar, Scalar> Pair;
627  const auto& waterVaporizationFac = saturatedWaterVaporizationFactorTable_[regionIdx];
628 
629  // create the taublated function representing saturation pressure depending of
630  // Rw
631  size_t n = waterVaporizationFac.numSamples();
632  Scalar delta = (waterVaporizationFac.xMax() - waterVaporizationFac.xMin())/Scalar(n + 1);
633 
634  SamplingPoints pSatSamplePoints;
635  Scalar Rw = 0;
636  for (size_t i = 0; i <= n; ++ i) {
637  Scalar pSat = waterVaporizationFac.xMin() + Scalar(i)*delta;
638  Rw = saturatedWaterVaporizationFactor(regionIdx, /*temperature=*/Scalar(1e30), pSat);
639 
640  Pair val(Rw, pSat);
641  pSatSamplePoints.push_back(val);
642  }
643 
644  //Prune duplicate Rv values (can occur, and will cause problems in further interpolation)
645  auto x_coord_comparator = [](const Pair& a, const Pair& b) { return a.first == b.first; };
646  auto last = std::unique(pSatSamplePoints.begin(), pSatSamplePoints.end(), x_coord_comparator);
647  if (std::distance(pSatSamplePoints.begin(), last) > 1) // only remove them if there are more than two points
648  pSatSamplePoints.erase(last, pSatSamplePoints.end());
649 
650  saturationPressure_[regionIdx].setContainerOfTuples(pSatSamplePoints);
651  }
652 
653  std::vector<Scalar> gasReferenceDensity_;
654  std::vector<Scalar> waterReferenceDensity_;
655  std::vector<TabulatedTwoDFunction> inverseGasB_;
656  std::vector<TabulatedOneDFunction> inverseSaturatedGasB_;
657  std::vector<TabulatedTwoDFunction> gasMu_;
658  std::vector<TabulatedTwoDFunction> inverseGasBMu_;
659  std::vector<TabulatedOneDFunction> inverseSaturatedGasBMu_;
660  std::vector<TabulatedOneDFunction> saturatedWaterVaporizationFactorTable_;
661  std::vector<TabulatedOneDFunction> saturationPressure_;
662 
663  Scalar vapPar1_;
664 };
665 
666 } // namespace Opm
667 
668 #endif
A central place for various physical constants occuring in some equations.
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
Implements a linearly interpolated scalar function that depends on one variable.
Implements a scalar function that depends on two variables and which is sampled uniformly in the X di...
This class represents the Pressure-Volume-Temperature relations of the gas phase with vaporized water...
Definition: DryHumidGasPvt.hpp:53
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &, const Evaluation &Rvw) const
Returns the formation volume factor [-] of the fluid phase.
Definition: DryHumidGasPvt.hpp:454
Evaluation saturatedOilVaporizationFactor(unsigned, const Evaluation &, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the oil vaporization factor [m^3/m^3] of the oil phase.
Definition: DryHumidGasPvt.hpp:485
void setGasViscosity(unsigned regionIdx, const TabulatedTwoDFunction &mug)
Initialize the viscosity of the gas phase.
Definition: DryHumidGasPvt.hpp:312
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &, const Evaluation &Rw) const
Returns the saturation pressure of the gas phase [Pa] depending on its mass fraction of the water com...
Definition: DryHumidGasPvt.hpp:509
void setInverseGasFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction &invBg)
Initialize the function for the gas formation volume factor.
Definition: DryHumidGasPvt.hpp:304
void initEnd()
Finish initializing the gas phase PVT properties.
Definition: DryHumidGasPvt.hpp:357
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: DryHumidGasPvt.hpp:405
void setReferenceDensities(unsigned regionIdx, Scalar, Scalar rhoRefGas, Scalar rhoRefWater)
Initialize the reference densities of all fluids for a given PVT region.
Definition: DryHumidGasPvt.hpp:275
Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the water vaporization factor [m^3/m^3] of the water phase.
Definition: DryHumidGasPvt.hpp:474
void setSaturatedGasWaterVaporizationFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the oil vaporization factor .
Definition: DryHumidGasPvt.hpp:289
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &, const Evaluation &Rvw) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: DryHumidGasPvt.hpp:424
void setSaturatedGasViscosity(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the phase viscosity for oil saturated gas.
Definition: DryHumidGasPvt.hpp:322
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of gas given a set of parameters.
Definition: DryHumidGasPvt.hpp:412
Evaluation saturatedOilVaporizationFactor(unsigned, const Evaluation &, const Evaluation &) const
Returns the oil vaporization factor [m^3/m^3] of the oil phase.
Definition: DryHumidGasPvt.hpp:496
Evaluation saturatedViscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of oil saturated gas at a given pressure.
Definition: DryHumidGasPvt.hpp:440
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of water saturated gas at a given pressure.
Definition: DryHumidGasPvt.hpp:465
Definition: Exceptions.hpp:46
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:48
Evaluation eval(const Evaluation &x, bool extrapolate=false) const
Evaluate the spline at a given position.
Definition: Tabulated1DFunction.hpp:259
void setContainerOfTuples(const XYContainer &points, bool sortInputs=true)
Set the sampling points of the piecewise linear function using a STL-compatible container of tuple-li...
Definition: Tabulated1DFunction.hpp:185
Implements a scalar function that depends on two variables and which is sampled uniformly in the X di...
Definition: UniformXTabulated2DFunction.hpp:54
Definition: MathToolbox.hpp:50