My Project
WetGasPvt.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_WET_GAS_PVT_HPP
28 #define OPM_WET_GAS_PVT_HPP
29 
35 
36 #if HAVE_ECL_INPUT
37 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
38 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
39 #include <opm/parser/eclipse/EclipseState/Schedule/OilVaporizationProperties.hpp>
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>
53 class WetGasPvt
54 {
55  typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
56 
57 public:
60 
61  WetGasPvt()
62  {
63  vapPar1_ = 0.0;
64  }
65 
66  WetGasPvt(const std::vector<Scalar>& gasReferenceDensity,
67  const std::vector<Scalar>& oilReferenceDensity,
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>& saturatedOilVaporizationFactorTable,
74  const std::vector<TabulatedOneDFunction>& saturationPressure,
75  Scalar vapPar1)
76  : gasReferenceDensity_(gasReferenceDensity)
77  , oilReferenceDensity_(oilReferenceDensity)
78  , inverseGasB_(inverseGasB)
79  , inverseSaturatedGasB_(inverseSaturatedGasB)
80  , gasMu_(gasMu)
81  , inverseGasBMu_(inverseGasBMu)
82  , inverseSaturatedGasBMu_(inverseSaturatedGasBMu)
83  , saturatedOilVaporizationFactorTable_(saturatedOilVaporizationFactorTable)
84  , saturationPressure_(saturationPressure)
85  , vapPar1_(vapPar1)
86  {
87  }
88 
89 
90 #if HAVE_ECL_INPUT
96  void initFromState(const EclipseState& eclState, const Schedule& schedule)
97  {
98  const auto& pvtgTables = eclState.getTableManager().getPvtgTables();
99  const auto& densityTable = eclState.getTableManager().getDensityTable();
100 
101  assert(pvtgTables.size() == densityTable.size());
102 
103  size_t numRegions = pvtgTables.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& pvtgTable = pvtgTables[regionIdx];
116 
117  const auto& saturatedTable = pvtgTable.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& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
125 
126  oilVaporizationFac.setXYArrays(saturatedTable.numRows(),
127  saturatedTable.getColumn("PG"),
128  saturatedTable.getColumn("RV"));
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 = pvtgTable.getUnderSaturatedTable(outerIdx);
149  size_t numRows = underSaturatedTable.numRows();
150  for (size_t innerIdx = 0; innerIdx < numRows; ++ innerIdx) {
151  Scalar Rv = underSaturatedTable.get("RV" , innerIdx);
152  Scalar Bg = underSaturatedTable.get("BG" , innerIdx);
153  Scalar mug = underSaturatedTable.get("MUG" , innerIdx);
154 
155  invGasB.appendSamplePoint(outerIdx, Rv, 1.0/Bg);
156  gasMu.appendSamplePoint(outerIdx, Rv, 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 (pvtgTable.getUnderSaturatedTable(masterTableIdx).numRows() > 1)
184  break;
185  }
186 
187  if (masterTableIdx >= saturatedTable.numRows())
188  throw std::runtime_error("PVTG 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  extendPvtgTable_(regionIdx,
194  xIdx,
195  pvtgTable.getUnderSaturatedTable(xIdx),
196  pvtgTable.getUnderSaturatedTable(masterTableIdx));
197  }
198  }
199 
200  vapPar1_ = 0.0;
201  const auto& oilVap = schedule[0].oilvap();
202  if (oilVap.getType() == OilVaporizationProperties::OilVaporization::VAPPARS) {
203  vapPar1_ = oilVap.vap1();
204  }
205 
206  initEnd();
207  }
208 
209 private:
210  void extendPvtgTable_(unsigned regionIdx,
211  unsigned xIdx,
212  const SimpleTable& curTable,
213  const SimpleTable& masterTable)
214  {
215  std::vector<double> RvArray = curTable.getColumn("RV").vectorCopy();
216  std::vector<double> gasBArray = curTable.getColumn("BG").vectorCopy();
217  std::vector<double> gasMuArray = curTable.getColumn("MUG").vectorCopy();
218 
219  auto& invGasB = inverseGasB_[regionIdx];
220  auto& gasMu = gasMu_[regionIdx];
221 
222  for (size_t newRowIdx = 1; newRowIdx < masterTable.numRows(); ++ newRowIdx) {
223  const auto& RVColumn = masterTable.getColumn("RV");
224  const auto& BGColumn = masterTable.getColumn("BG");
225  const auto& viscosityColumn = masterTable.getColumn("MUG");
226 
227  // compute the gas pressure for the new entry
228  Scalar diffRv = RVColumn[newRowIdx] - RVColumn[newRowIdx - 1];
229  Scalar newRv = RvArray.back() + diffRv;
230 
231  // calculate the compressibility of the master table
232  Scalar B1 = BGColumn[newRowIdx];
233  Scalar B2 = BGColumn[newRowIdx - 1];
234  Scalar x = (B1 - B2)/( (B1 + B2)/2.0 );
235 
236  // calculate the gas formation volume factor which exhibits the same
237  // "compressibility" for the new value of Rv
238  Scalar newBg = gasBArray.back()*(1.0 + x/2.0)/(1.0 - x/2.0);
239 
240  // calculate the "viscosibility" of the master table
241  Scalar mu1 = viscosityColumn[newRowIdx];
242  Scalar mu2 = viscosityColumn[newRowIdx - 1];
243  Scalar xMu = (mu1 - mu2)/( (mu1 + mu2)/2.0 );
244 
245  // calculate the gas formation volume factor which exhibits the same
246  // compressibility for the new pressure
247  Scalar newMug = gasMuArray.back()*(1.0 + xMu/2)/(1.0 - xMu/2.0);
248 
249  // append the new values to the arrays which we use to compute the additional
250  // values ...
251  RvArray.push_back(newRv);
252  gasBArray.push_back(newBg);
253  gasMuArray.push_back(newMug);
254 
255  // ... and register them with the internal table objects
256  invGasB.appendSamplePoint(xIdx, newRv, 1.0/newBg);
257  gasMu.appendSamplePoint(xIdx, newRv, newMug);
258  }
259  }
260 
261 public:
262 #endif // HAVE_ECL_INPUT
263 
264  void setNumRegions(size_t numRegions)
265  {
266  oilReferenceDensity_.resize(numRegions);
267  gasReferenceDensity_.resize(numRegions);
268  inverseGasB_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::RightExtreme});
269  inverseGasBMu_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::RightExtreme});
270  inverseSaturatedGasB_.resize(numRegions);
271  inverseSaturatedGasBMu_.resize(numRegions);
272  gasMu_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::RightExtreme});
273  saturatedOilVaporizationFactorTable_.resize(numRegions);
274  saturationPressure_.resize(numRegions);
275  }
276 
280  void setReferenceDensities(unsigned regionIdx,
281  Scalar rhoRefOil,
282  Scalar rhoRefGas,
283  Scalar /*rhoRefWater*/)
284  {
285  oilReferenceDensity_[regionIdx] = rhoRefOil;
286  gasReferenceDensity_[regionIdx] = rhoRefGas;
287  }
288 
294  void setSaturatedGasOilVaporizationFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
295  { saturatedOilVaporizationFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
296 
306  void setSaturatedGasFormationVolumeFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
307  {
308  auto& invGasB = inverseGasB_[regionIdx];
309 
310  const auto& RvTable = saturatedOilVaporizationFactorTable_[regionIdx];
311 
312  Scalar T = 273.15 + 15.56; // [K]
313 
314  Scalar RvMin = 0.0;
315  Scalar RvMax = RvTable.eval(saturatedOilVaporizationFactorTable_[regionIdx].xMax(), /*extrapolate=*/true);
316 
317  Scalar poMin = samplePoints.front().first;
318  Scalar poMax = samplePoints.back().first;
319 
320  size_t nRv = 20;
321  size_t nP = samplePoints.size()*2;
322 
323  Scalar rhooRef = oilReferenceDensity_[regionIdx];
324 
325  TabulatedOneDFunction gasFormationVolumeFactor;
326  gasFormationVolumeFactor.setContainerOfTuples(samplePoints);
327 
328  updateSaturationPressure_(regionIdx);
329 
330  // calculate a table of estimated densities depending on pressure and gas mass
331  // fraction. note that this assumes oil of constant compressibility. (having said
332  // that, if only the saturated gas densities are available, there's not much
333  // choice.)
334  for (size_t RvIdx = 0; RvIdx < nRv; ++RvIdx) {
335  Scalar Rv = RvMin + (RvMax - RvMin)*RvIdx/nRv;
336 
337  invGasB.appendXPos(Rv);
338 
339  for (size_t pIdx = 0; pIdx < nP; ++pIdx) {
340  Scalar pg = poMin + (poMax - poMin)*pIdx/nP;
341 
342  Scalar poSat = saturationPressure(regionIdx, T, Rv);
343  Scalar BgSat = gasFormationVolumeFactor.eval(poSat, /*extrapolate=*/true);
344  Scalar drhoo_dp = (1.1200 - 1.1189)/((5000 - 4000)*6894.76);
345  Scalar rhoo = rhooRef/BgSat*(1 + drhoo_dp*(pg - poSat));
346 
347  Scalar Bg = rhooRef/rhoo;
348 
349  invGasB.appendSamplePoint(RvIdx, pg, 1.0/Bg);
350  }
351  }
352  }
353 
366  void setInverseGasFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction& invBg)
367  { inverseGasB_[regionIdx] = invBg; }
368 
374  void setGasViscosity(unsigned regionIdx, const TabulatedTwoDFunction& mug)
375  { gasMu_[regionIdx] = mug; }
376 
384  void setSaturatedGasViscosity(unsigned regionIdx, const SamplingPoints& samplePoints )
385  {
386  auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
387 
388  Scalar RvMin = 0.0;
389  Scalar RvMax = oilVaporizationFac.eval(saturatedOilVaporizationFactorTable_[regionIdx].xMax(), /*extrapolate=*/true);
390 
391  Scalar poMin = samplePoints.front().first;
392  Scalar poMax = samplePoints.back().first;
393 
394  size_t nRv = 20;
395  size_t nP = samplePoints.size()*2;
396 
397  TabulatedOneDFunction mugTable;
398  mugTable.setContainerOfTuples(samplePoints);
399 
400  // calculate a table of estimated densities depending on pressure and gas mass
401  // fraction
402  for (size_t RvIdx = 0; RvIdx < nRv; ++RvIdx) {
403  Scalar Rv = RvMin + (RvMax - RvMin)*RvIdx/nRv;
404 
405  gasMu_[regionIdx].appendXPos(Rv);
406 
407  for (size_t pIdx = 0; pIdx < nP; ++pIdx) {
408  Scalar pg = poMin + (poMax - poMin)*pIdx/nP;
409  Scalar mug = mugTable.eval(pg, /*extrapolate=*/true);
410 
411  gasMu_[regionIdx].appendSamplePoint(RvIdx, pg, mug);
412  }
413  }
414  }
415 
419  void initEnd()
420  {
421  // calculate the final 2D functions which are used for interpolation.
422  size_t numRegions = gasMu_.size();
423  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
424  // calculate the table which stores the inverse of the product of the gas
425  // formation volume factor and the gas viscosity
426  const auto& gasMu = gasMu_[regionIdx];
427  const auto& invGasB = inverseGasB_[regionIdx];
428  assert(gasMu.numX() == invGasB.numX());
429 
430  auto& invGasBMu = inverseGasBMu_[regionIdx];
431  auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
432  auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
433 
434  std::vector<Scalar> satPressuresArray;
435  std::vector<Scalar> invSatGasBArray;
436  std::vector<Scalar> invSatGasBMuArray;
437  for (size_t pIdx = 0; pIdx < gasMu.numX(); ++pIdx) {
438  invGasBMu.appendXPos(gasMu.xAt(pIdx));
439 
440  assert(gasMu.numY(pIdx) == invGasB.numY(pIdx));
441 
442  size_t numRv = gasMu.numY(pIdx);
443  for (size_t rvIdx = 0; rvIdx < numRv; ++rvIdx)
444  invGasBMu.appendSamplePoint(pIdx,
445  gasMu.yAt(pIdx, rvIdx),
446  invGasB.valueAt(pIdx, rvIdx)
447  / gasMu.valueAt(pIdx, rvIdx));
448 
449  // the sampling points in UniformXTabulated2DFunction are always sorted
450  // in ascending order. Thus, the value for saturated gas is the last one
451  // (i.e., the one with the largest Rv value)
452  satPressuresArray.push_back(gasMu.xAt(pIdx));
453  invSatGasBArray.push_back(invGasB.valueAt(pIdx, numRv - 1));
454  invSatGasBMuArray.push_back(invGasBMu.valueAt(pIdx, numRv - 1));
455  }
456 
457  invSatGasB.setXYContainers(satPressuresArray, invSatGasBArray);
458  invSatGasBMu.setXYContainers(satPressuresArray, invSatGasBMuArray);
459 
460  updateSaturationPressure_(regionIdx);
461  }
462  }
463 
467  unsigned numRegions() const
468  { return gasReferenceDensity_.size(); }
469 
473  template <class Evaluation>
474  Evaluation internalEnergy(unsigned,
475  const Evaluation&,
476  const Evaluation&,
477  const Evaluation&) const
478  {
479  throw std::runtime_error("Requested the enthalpy of gas but the thermal option is not enabled");
480  }
481 
485  template <class Evaluation>
486  Evaluation viscosity(unsigned regionIdx,
487  const Evaluation& /*temperature*/,
488  const Evaluation& pressure,
489  const Evaluation& Rv) const
490  {
491  const Evaluation& invBg = inverseGasB_[regionIdx].eval(pressure, Rv, /*extrapolate=*/true);
492  const Evaluation& invMugBg = inverseGasBMu_[regionIdx].eval(pressure, Rv, /*extrapolate=*/true);
493 
494  return invBg/invMugBg;
495  }
496 
500  template <class Evaluation>
501  Evaluation saturatedViscosity(unsigned regionIdx,
502  const Evaluation& /*temperature*/,
503  const Evaluation& pressure) const
504  {
505  const Evaluation& invBg = inverseSaturatedGasB_[regionIdx].eval(pressure, /*extrapolate=*/true);
506  const Evaluation& invMugBg = inverseSaturatedGasBMu_[regionIdx].eval(pressure, /*extrapolate=*/true);
507 
508  return invBg/invMugBg;
509  }
510 
514  template <class Evaluation>
515  Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
516  const Evaluation& /*temperature*/,
517  const Evaluation& pressure,
518  const Evaluation& Rv) const
519  { return inverseGasB_[regionIdx].eval(pressure, Rv, /*extrapolate=*/true); }
520 
524  template <class Evaluation>
525  Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
526  const Evaluation& /*temperature*/,
527  const Evaluation& pressure) const
528  { return inverseSaturatedGasB_[regionIdx].eval(pressure, /*extrapolate=*/true); }
529 
533  template <class Evaluation>
534  Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
535  const Evaluation& /*temperature*/,
536  const Evaluation& pressure) const
537  {
538  return saturatedOilVaporizationFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
539  }
540 
548  template <class Evaluation>
549  Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
550  const Evaluation& /*temperature*/,
551  const Evaluation& pressure,
552  const Evaluation& oilSaturation,
553  Evaluation maxOilSaturation) const
554  {
555  Evaluation tmp =
556  saturatedOilVaporizationFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
557 
558  // apply the vaporization parameters for the gas phase (cf. the Eclipse VAPPARS
559  // keyword)
560  maxOilSaturation = min(maxOilSaturation, Scalar(1.0));
561  if (vapPar1_ > 0.0 && maxOilSaturation > 0.01 && oilSaturation < maxOilSaturation) {
562  static const Scalar eps = 0.001;
563  const Evaluation& So = max(oilSaturation, eps);
564  tmp *= max(1e-3, pow(So/maxOilSaturation, vapPar1_));
565  }
566 
567  return tmp;
568  }
569 
580  template <class Evaluation>
581  Evaluation saturationPressure(unsigned regionIdx,
582  const Evaluation&,
583  const Evaluation& Rv) const
584  {
585  typedef MathToolbox<Evaluation> Toolbox;
586 
587  const auto& RvTable = saturatedOilVaporizationFactorTable_[regionIdx];
588  const Scalar eps = std::numeric_limits<typename Toolbox::Scalar>::epsilon()*1e6;
589 
590  // use the tabulated saturation pressure function to get a pretty good initial value
591  Evaluation pSat = saturationPressure_[regionIdx].eval(Rv, /*extrapolate=*/true);
592 
593  // Newton method to do the remaining work. If the initial
594  // value is good, this should only take two to three
595  // iterations...
596  bool onProbation = false;
597  for (unsigned i = 0; i < 20; ++i) {
598  const Evaluation& f = RvTable.eval(pSat, /*extrapolate=*/true) - Rv;
599  const Evaluation& fPrime = RvTable.evalDerivative(pSat, /*extrapolate=*/true);
600 
601  // If the derivative is "zero" Newton will not converge,
602  // so simply return our initial guess.
603  if (std::abs(scalarValue(fPrime)) < 1.0e-30) {
604  return pSat;
605  }
606 
607  const Evaluation& delta = f/fPrime;
608 
609  pSat -= delta;
610 
611  if (pSat < 0.0) {
612  // if the pressure is lower than 0 Pascals, we set it back to 0. if this
613  // happens twice, we give up and just return 0 Pa...
614  if (onProbation)
615  return 0.0;
616 
617  onProbation = true;
618  pSat = 0.0;
619  }
620 
621  if (std::abs(scalarValue(delta)) < std::abs(scalarValue(pSat))*eps)
622  return pSat;
623  }
624 
625  std::stringstream errlog;
626  errlog << "Finding saturation pressure did not converge:"
627  << " pSat = " << pSat
628  << ", Rv = " << Rv;
629 #if HAVE_OPM_COMMON
630  OpmLog::debug("Wet gas saturation pressure", errlog.str());
631 #endif
632  throw NumericalIssue(errlog.str());
633  }
634 
635  template <class Evaluation>
636  Evaluation diffusionCoefficient(const Evaluation& /*temperature*/,
637  const Evaluation& /*pressure*/,
638  unsigned /*compIdx*/) const
639  {
640  throw std::runtime_error("Not implemented: The PVT model does not provide a diffusionCoefficient()");
641  }
642 
643  const Scalar gasReferenceDensity(unsigned regionIdx) const
644  { return gasReferenceDensity_[regionIdx]; }
645 
646  const Scalar oilReferenceDensity(unsigned regionIdx) const
647  { return oilReferenceDensity_[regionIdx]; }
648 
649  const std::vector<TabulatedTwoDFunction>& inverseGasB() const {
650  return inverseGasB_;
651  }
652 
653  const std::vector<TabulatedOneDFunction>& inverseSaturatedGasB() const {
654  return inverseSaturatedGasB_;
655  }
656 
657  const std::vector<TabulatedTwoDFunction>& gasMu() const {
658  return gasMu_;
659  }
660 
661  const std::vector<TabulatedTwoDFunction>& inverseGasBMu() const {
662  return inverseGasBMu_;
663  }
664 
665  const std::vector<TabulatedOneDFunction>& inverseSaturatedGasBMu() const {
666  return inverseSaturatedGasBMu_;
667  }
668 
669  const std::vector<TabulatedOneDFunction>& saturatedOilVaporizationFactorTable() const {
670  return saturatedOilVaporizationFactorTable_;
671  }
672 
673  const std::vector<TabulatedOneDFunction>& saturationPressure() const {
674  return saturationPressure_;
675  }
676 
677  Scalar vapPar1() const {
678  return vapPar1_;
679  }
680 
681  bool operator==(const WetGasPvt<Scalar>& data) const
682  {
683  return this->gasReferenceDensity_ == data.gasReferenceDensity_ &&
684  this->oilReferenceDensity_ == data.oilReferenceDensity_ &&
685  this->inverseGasB() == data.inverseGasB() &&
686  this->inverseSaturatedGasB() == data.inverseSaturatedGasB() &&
687  this->gasMu() == data.gasMu() &&
688  this->inverseGasBMu() == data.inverseGasBMu() &&
689  this->inverseSaturatedGasBMu() == data.inverseSaturatedGasBMu() &&
690  this->saturatedOilVaporizationFactorTable() == data.saturatedOilVaporizationFactorTable() &&
691  this->saturationPressure() == data.saturationPressure() &&
692  this->vapPar1() == data.vapPar1();
693  }
694 
695 private:
696  void updateSaturationPressure_(unsigned regionIdx)
697  {
698  typedef std::pair<Scalar, Scalar> Pair;
699  const auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
700 
701  // create the taublated function representing saturation pressure depending of
702  // Rv
703  size_t n = oilVaporizationFac.numSamples();
704  Scalar delta = (oilVaporizationFac.xMax() - oilVaporizationFac.xMin())/Scalar(n + 1);
705 
706  SamplingPoints pSatSamplePoints;
707  Scalar Rv = 0;
708  for (size_t i = 0; i <= n; ++ i) {
709  Scalar pSat = oilVaporizationFac.xMin() + Scalar(i)*delta;
710  Rv = saturatedOilVaporizationFactor(regionIdx, /*temperature=*/Scalar(1e30), pSat);
711 
712  Pair val(Rv, pSat);
713  pSatSamplePoints.push_back(val);
714  }
715 
716  //Prune duplicate Rv values (can occur, and will cause problems in further interpolation)
717  auto x_coord_comparator = [](const Pair& a, const Pair& b) { return a.first == b.first; };
718  auto last = std::unique(pSatSamplePoints.begin(), pSatSamplePoints.end(), x_coord_comparator);
719  pSatSamplePoints.erase(last, pSatSamplePoints.end());
720 
721  saturationPressure_[regionIdx].setContainerOfTuples(pSatSamplePoints);
722  }
723 
724  std::vector<Scalar> gasReferenceDensity_;
725  std::vector<Scalar> oilReferenceDensity_;
726  std::vector<TabulatedTwoDFunction> inverseGasB_;
727  std::vector<TabulatedOneDFunction> inverseSaturatedGasB_;
728  std::vector<TabulatedTwoDFunction> gasMu_;
729  std::vector<TabulatedTwoDFunction> inverseGasBMu_;
730  std::vector<TabulatedOneDFunction> inverseSaturatedGasBMu_;
731  std::vector<TabulatedOneDFunction> saturatedOilVaporizationFactorTable_;
732  std::vector<TabulatedOneDFunction> saturationPressure_;
733 
734  Scalar vapPar1_;
735 };
736 
737 } // namespace Opm
738 
739 #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...
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
This class represents the Pressure-Volume-Temperature relations of the gas phas with vaporized oil.
Definition: WetGasPvt.hpp:54
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rv) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WetGasPvt.hpp:515
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of oil saturated gas at a given pressure.
Definition: WetGasPvt.hpp:525
Evaluation saturatedOilVaporizationFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the oil vaporization factor [m^3/m^3] of the gas phase.
Definition: WetGasPvt.hpp:534
void initEnd()
Finish initializing the gas phase PVT properties.
Definition: WetGasPvt.hpp:419
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &, const Evaluation &Rv) const
Returns the saturation pressure of the gas phase [Pa] depending on its mass fraction of the oil compo...
Definition: WetGasPvt.hpp:581
Evaluation saturatedOilVaporizationFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &oilSaturation, Evaluation maxOilSaturation) const
Returns the oil vaporization factor [m^3/m^3] of the gas phase.
Definition: WetGasPvt.hpp:549
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rv) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WetGasPvt.hpp:486
void setGasViscosity(unsigned regionIdx, const TabulatedTwoDFunction &mug)
Initialize the viscosity of the gas phase.
Definition: WetGasPvt.hpp:374
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of gas given a set of parameters.
Definition: WetGasPvt.hpp:474
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: WetGasPvt.hpp:467
void setSaturatedGasViscosity(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the phase viscosity for oil saturated gas.
Definition: WetGasPvt.hpp:384
void setSaturatedGasFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the gas formation volume factor.
Definition: WetGasPvt.hpp:306
void setSaturatedGasOilVaporizationFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the oil vaporization factor .
Definition: WetGasPvt.hpp:294
void setInverseGasFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction &invBg)
Initialize the function for the gas formation volume factor.
Definition: WetGasPvt.hpp:366
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: WetGasPvt.hpp:501
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar rhoRefGas, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: WetGasPvt.hpp:280
Definition: MathToolbox.hpp:52