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