My Project
LiveOilPvt.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_LIVE_OIL_PVT_HPP
28 #define OPM_LIVE_OIL_PVT_HPP
29 
35 
36 #if HAVE_ECL_INPUT
37 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
38 #include <opm/parser/eclipse/EclipseState/Schedule/OilVaporizationProperties.hpp>
39 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
40 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
41 #endif
42 
43 #if HAVE_OPM_COMMON
44 #include <opm/common/OpmLog/OpmLog.hpp>
45 #endif
46 
47 namespace Opm {
52 template <class Scalar>
54 {
55  typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
56 
57 public:
60 
61  LiveOilPvt()
62  {
63  vapPar2_ = 0.0;
64  }
65 
66  LiveOilPvt(const std::vector<Scalar>& gasReferenceDensity,
67  const std::vector<Scalar>& oilReferenceDensity,
68  const std::vector<TabulatedTwoDFunction>& inverseOilBTable,
69  const std::vector<TabulatedTwoDFunction>& oilMuTable,
70  const std::vector<TabulatedTwoDFunction>& inverseOilBMuTable,
71  const std::vector<TabulatedOneDFunction>& saturatedOilMuTable,
72  const std::vector<TabulatedOneDFunction>& inverseSaturatedOilBTable,
73  const std::vector<TabulatedOneDFunction>& inverseSaturatedOilBMuTable,
74  const std::vector<TabulatedOneDFunction>& saturatedGasDissolutionFactorTable,
75  const std::vector<TabulatedOneDFunction>& saturationPressure,
76  Scalar vapPar2)
77  : gasReferenceDensity_(gasReferenceDensity)
78  , oilReferenceDensity_(oilReferenceDensity)
79  , inverseOilBTable_(inverseOilBTable)
80  , oilMuTable_(oilMuTable)
81  , inverseOilBMuTable_(inverseOilBMuTable)
82  , saturatedOilMuTable_(saturatedOilMuTable)
83  , inverseSaturatedOilBTable_(inverseSaturatedOilBTable)
84  , inverseSaturatedOilBMuTable_(inverseSaturatedOilBMuTable)
85  , saturatedGasDissolutionFactorTable_(saturatedGasDissolutionFactorTable)
86  , saturationPressure_(saturationPressure)
87  , vapPar2_(vapPar2)
88  { }
89 
90 #if HAVE_ECL_INPUT
94  void initFromState(const EclipseState& eclState, const Schedule& schedule)
95  {
96  const auto& pvtoTables = eclState.getTableManager().getPvtoTables();
97  const auto& densityTable = eclState.getTableManager().getDensityTable();
98 
99  assert(pvtoTables.size() == densityTable.size());
100 
101  size_t numRegions = pvtoTables.size();
102  setNumRegions(numRegions);
103 
104  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
105  Scalar rhoRefO = densityTable[regionIdx].oil;
106  Scalar rhoRefG = densityTable[regionIdx].gas;
107  Scalar rhoRefW = densityTable[regionIdx].water;
108 
109  setReferenceDensities(regionIdx, rhoRefO, rhoRefG, rhoRefW);
110  }
111 
112  // initialize the internal table objects
113  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
114  const auto& pvtoTable = pvtoTables[regionIdx];
115 
116  const auto& saturatedTable = pvtoTable.getSaturatedTable();
117  assert(saturatedTable.numRows() > 1);
118 
119  auto& oilMu = oilMuTable_[regionIdx];
120  auto& satOilMu = saturatedOilMuTable_[regionIdx];
121  auto& invOilB = inverseOilBTable_[regionIdx];
122  auto& invSatOilB = inverseSaturatedOilBTable_[regionIdx];
123  auto& gasDissolutionFac = saturatedGasDissolutionFactorTable_[regionIdx];
124  std::vector<Scalar> invSatOilBArray;
125  std::vector<Scalar> satOilMuArray;
126 
127  // extract the table for the gas dissolution and the oil formation volume factors
128  for (unsigned outerIdx = 0; outerIdx < saturatedTable.numRows(); ++ outerIdx) {
129  Scalar Rs = saturatedTable.get("RS", outerIdx);
130  Scalar BoSat = saturatedTable.get("BO", outerIdx);
131  Scalar muoSat = saturatedTable.get("MU", outerIdx);
132 
133  satOilMuArray.push_back(muoSat);
134  invSatOilBArray.push_back(1.0/BoSat);
135 
136  invOilB.appendXPos(Rs);
137  oilMu.appendXPos(Rs);
138 
139  assert(invOilB.numX() == outerIdx + 1);
140  assert(oilMu.numX() == outerIdx + 1);
141 
142  const auto& underSaturatedTable = pvtoTable.getUnderSaturatedTable(outerIdx);
143  size_t numRows = underSaturatedTable.numRows();
144  for (unsigned innerIdx = 0; innerIdx < numRows; ++ innerIdx) {
145  Scalar po = underSaturatedTable.get("P", innerIdx);
146  Scalar Bo = underSaturatedTable.get("BO", innerIdx);
147  Scalar muo = underSaturatedTable.get("MU", innerIdx);
148 
149  invOilB.appendSamplePoint(outerIdx, po, 1.0/Bo);
150  oilMu.appendSamplePoint(outerIdx, po, muo);
151  }
152  }
153 
154  // update the tables for the formation volume factor and for the gas
155  // dissolution factor of saturated oil
156  {
157  const auto& tmpPressureColumn = saturatedTable.getColumn("P");
158  const auto& tmpGasSolubilityColumn = saturatedTable.getColumn("RS");
159 
160  invSatOilB.setXYContainers(tmpPressureColumn, invSatOilBArray);
161  satOilMu.setXYContainers(tmpPressureColumn, satOilMuArray);
162  gasDissolutionFac.setXYContainers(tmpPressureColumn, tmpGasSolubilityColumn);
163  }
164 
165  updateSaturationPressure_(regionIdx);
166  // make sure to have at least two sample points per Rs value
167  for (unsigned xIdx = 0; xIdx < invOilB.numX(); ++xIdx) {
168  // a single sample point is definitely needed
169  assert(invOilB.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 (invOilB.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 oil...
179  size_t masterTableIdx = xIdx + 1;
180  for (; masterTableIdx < saturatedTable.numRows(); ++masterTableIdx)
181  {
182  if (pvtoTable.getUnderSaturatedTable(masterTableIdx).numRows() > 1)
183  break;
184  }
185 
186  if (masterTableIdx >= saturatedTable.numRows())
187  throw std::runtime_error("PVTO tables are invalid: The last table must exhibit at least one "
188  "entry for undersaturated oil!");
189 
190  // extend the current table using the master table.
191  extendPvtoTable_(regionIdx,
192  xIdx,
193  pvtoTable.getUnderSaturatedTable(xIdx),
194  pvtoTable.getUnderSaturatedTable(masterTableIdx));
195  }
196  }
197 
198  vapPar2_ = 0.0;
199  const auto& oilVap = schedule[0].oilvap();
200  if (oilVap.getType() == OilVaporizationProperties::OilVaporization::VAPPARS) {
201  vapPar2_ = oilVap.vap2();
202  }
203 
204  initEnd();
205  }
206 
207 private:
208  void extendPvtoTable_(unsigned regionIdx,
209  unsigned xIdx,
210  const SimpleTable& curTable,
211  const SimpleTable& masterTable)
212  {
213  std::vector<double> pressuresArray = curTable.getColumn("P").vectorCopy();
214  std::vector<double> oilBArray = curTable.getColumn("BO").vectorCopy();
215  std::vector<double> oilMuArray = curTable.getColumn("MU").vectorCopy();
216 
217  auto& invOilB = inverseOilBTable_[regionIdx];
218  auto& oilMu = oilMuTable_[regionIdx];
219 
220  for (unsigned newRowIdx = 1; newRowIdx < masterTable.numRows(); ++ newRowIdx) {
221  const auto& pressureColumn = masterTable.getColumn("P");
222  const auto& BOColumn = masterTable.getColumn("BO");
223  const auto& viscosityColumn = masterTable.getColumn("MU");
224 
225  // compute the oil pressure for the new entry
226  Scalar diffPo = pressureColumn[newRowIdx] - pressureColumn[newRowIdx - 1];
227  Scalar newPo = pressuresArray.back() + diffPo;
228 
229  // calculate the compressibility of the master table
230  Scalar B1 = BOColumn[newRowIdx];
231  Scalar B2 = BOColumn[newRowIdx - 1];
232  Scalar x = (B1 - B2)/( (B1 + B2)/2.0 );
233 
234  // calculate the oil formation volume factor which exhibits the same
235  // compressibility for the new pressure
236  Scalar newBo = oilBArray.back()*(1.0 + x/2.0)/(1.0 - x/2.0);
237 
238  // calculate the "viscosibility" of the master table
239  Scalar mu1 = viscosityColumn[newRowIdx];
240  Scalar mu2 = viscosityColumn[newRowIdx - 1];
241  Scalar xMu = (mu1 - mu2)/( (mu1 + mu2)/2.0 );
242 
243  // calculate the oil formation volume factor which exhibits the same
244  // compressibility for the new pressure
245  Scalar newMuo = oilMuArray.back()*(1.0 + xMu/2)/(1.0 - xMu/2.0);
246 
247  // append the new values to the arrays which we use to compute the additional
248  // values ...
249  pressuresArray.push_back(newPo);
250  oilBArray.push_back(newBo);
251  oilMuArray.push_back(newMuo);
252 
253  // ... and register them with the internal table objects
254  invOilB.appendSamplePoint(xIdx, newPo, 1.0/newBo);
255  oilMu.appendSamplePoint(xIdx, newPo, newMuo);
256  }
257  }
258 
259 public:
260 #endif // HAVE_ECL_INPUT
261 
262  void setNumRegions(size_t numRegions)
263  {
264  oilReferenceDensity_.resize(numRegions);
265  gasReferenceDensity_.resize(numRegions);
266  inverseOilBTable_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::LeftExtreme});
267  inverseOilBMuTable_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::LeftExtreme});
268  inverseSaturatedOilBTable_.resize(numRegions);
269  inverseSaturatedOilBMuTable_.resize(numRegions);
270  oilMuTable_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::LeftExtreme});
271  saturatedOilMuTable_.resize(numRegions);
272  saturatedGasDissolutionFactorTable_.resize(numRegions);
273  saturationPressure_.resize(numRegions);
274  }
275 
279  void setReferenceDensities(unsigned regionIdx,
280  Scalar rhoRefOil,
281  Scalar rhoRefGas,
282  Scalar /*rhoRefWater*/)
283  {
284  oilReferenceDensity_[regionIdx] = rhoRefOil;
285  gasReferenceDensity_[regionIdx] = rhoRefGas;
286  }
287 
293  void setSaturatedOilGasDissolutionFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
294  { saturatedGasDissolutionFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
295 
305  void setSaturatedOilFormationVolumeFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
306  {
307  Scalar T = 273.15 + 15.56; // [K]
308  auto& invOilB = inverseOilBTable_[regionIdx];
309 
310  updateSaturationPressure_(regionIdx);
311 
312  // calculate a table of estimated densities of undersatured gas
313  for (size_t pIdx = 0; pIdx < samplePoints.size(); ++pIdx) {
314  Scalar p1 = std::get<0>(samplePoints[pIdx]);
315  Scalar p2 = p1 * 2.0;
316 
317  Scalar Bo1 = std::get<1>(samplePoints[pIdx]);
318  Scalar drhoo_dp = (1.1200 - 1.1189)/((5000 - 4000)*6894.76);
319  Scalar Bo2 = Bo1/(1.0 + (p2 - p1)*drhoo_dp);
320 
321  Scalar Rs = saturatedGasDissolutionFactor(regionIdx, T, p1);
322 
323  invOilB.appendXPos(Rs);
324  invOilB.appendSamplePoint(pIdx, p1, 1.0/Bo1);
325  invOilB.appendSamplePoint(pIdx, p2, 1.0/Bo2);
326  }
327  }
328 
341  void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction& invBo)
342  { inverseOilBTable_[regionIdx] = invBo; }
343 
349  void setOilViscosity(unsigned regionIdx, const TabulatedTwoDFunction& muo)
350  { oilMuTable_[regionIdx] = muo; }
351 
359  void setSaturatedOilViscosity(unsigned regionIdx, const SamplingPoints& samplePoints)
360  {
361  Scalar T = 273.15 + 15.56; // [K]
362 
363  // update the table for the saturated oil
364  saturatedOilMuTable_[regionIdx].setContainerOfTuples(samplePoints);
365 
366  // calculate a table of estimated viscosities depending on pressure and gas mass
367  // fraction for untersaturated oil to make the other code happy
368  for (size_t pIdx = 0; pIdx < samplePoints.size(); ++pIdx) {
369  Scalar p1 = std::get<0>(samplePoints[pIdx]);
370  Scalar p2 = p1 * 2.0;
371 
372  // no pressure dependence of the viscosity
373  Scalar mu1 = std::get<1>(samplePoints[pIdx]);
374  Scalar mu2 = mu1;
375 
376  Scalar Rs = saturatedGasDissolutionFactor(regionIdx, T, p1);
377 
378  oilMuTable_[regionIdx].appendXPos(Rs);
379  oilMuTable_[regionIdx].appendSamplePoint(pIdx, p1, mu1);
380  oilMuTable_[regionIdx].appendSamplePoint(pIdx, p2, mu2);
381  }
382  }
383 
387  void initEnd()
388  {
389  // calculate the final 2D functions which are used for interpolation.
390  size_t numRegions = oilMuTable_.size();
391  for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
392  // calculate the table which stores the inverse of the product of the oil
393  // formation volume factor and the oil viscosity
394  const auto& oilMu = oilMuTable_[regionIdx];
395  const auto& satOilMu = saturatedOilMuTable_[regionIdx];
396  const auto& invOilB = inverseOilBTable_[regionIdx];
397  assert(oilMu.numX() == invOilB.numX());
398 
399  auto& invOilBMu = inverseOilBMuTable_[regionIdx];
400  auto& invSatOilB = inverseSaturatedOilBTable_[regionIdx];
401  auto& invSatOilBMu = inverseSaturatedOilBMuTable_[regionIdx];
402 
403  std::vector<Scalar> satPressuresArray;
404  std::vector<Scalar> invSatOilBArray;
405  std::vector<Scalar> invSatOilBMuArray;
406  for (unsigned rsIdx = 0; rsIdx < oilMu.numX(); ++rsIdx) {
407  invOilBMu.appendXPos(oilMu.xAt(rsIdx));
408 
409  assert(oilMu.numY(rsIdx) == invOilB.numY(rsIdx));
410 
411  size_t numPressures = oilMu.numY(rsIdx);
412  for (unsigned pIdx = 0; pIdx < numPressures; ++pIdx)
413  invOilBMu.appendSamplePoint(rsIdx,
414  oilMu.yAt(rsIdx, pIdx),
415  invOilB.valueAt(rsIdx, pIdx)
416  / oilMu.valueAt(rsIdx, pIdx));
417 
418  // the sampling points in UniformXTabulated2DFunction are always sorted
419  // in ascending order. Thus, the value for saturated oil is the first one
420  // (i.e., the one for the lowest pressure value)
421  satPressuresArray.push_back(oilMu.yAt(rsIdx, 0));
422  invSatOilBArray.push_back(invOilB.valueAt(rsIdx, 0));
423  invSatOilBMuArray.push_back(invSatOilBArray.back()/satOilMu.valueAt(rsIdx));
424  }
425 
426  invSatOilB.setXYContainers(satPressuresArray, invSatOilBArray);
427  invSatOilBMu.setXYContainers(satPressuresArray, invSatOilBMuArray);
428 
429  updateSaturationPressure_(regionIdx);
430  }
431  }
432 
436  unsigned numRegions() const
437  { return inverseOilBMuTable_.size(); }
438 
442  template <class Evaluation>
443  Evaluation internalEnergy(unsigned,
444  const Evaluation&,
445  const Evaluation&,
446  const Evaluation&) const
447  {
448  throw std::runtime_error("Requested the enthalpy of oil but the thermal option is not enabled");
449  }
450 
454  template <class Evaluation>
455  Evaluation viscosity(unsigned regionIdx,
456  const Evaluation& /*temperature*/,
457  const Evaluation& pressure,
458  const Evaluation& Rs) const
459  {
460  // ATTENTION: Rs is the first axis!
461  const Evaluation& invBo = inverseOilBTable_[regionIdx].eval(Rs, pressure, /*extrapolate=*/true);
462  const Evaluation& invMuoBo = inverseOilBMuTable_[regionIdx].eval(Rs, pressure, /*extrapolate=*/true);
463 
464  return invBo/invMuoBo;
465  }
466 
470  template <class Evaluation>
471  Evaluation saturatedViscosity(unsigned regionIdx,
472  const Evaluation& /*temperature*/,
473  const Evaluation& pressure) const
474  {
475  // ATTENTION: Rs is the first axis!
476  const Evaluation& invBo = inverseSaturatedOilBTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
477  const Evaluation& invMuoBo = inverseSaturatedOilBMuTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
478 
479  return invBo/invMuoBo;
480  }
481 
485  template <class Evaluation>
486  Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
487  const Evaluation& /*temperature*/,
488  const Evaluation& pressure,
489  const Evaluation& Rs) const
490  {
491  // ATTENTION: Rs is represented by the _first_ axis!
492  return inverseOilBTable_[regionIdx].eval(Rs, pressure, /*extrapolate=*/true);
493  }
494 
498  template <class Evaluation>
499  Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
500  const Evaluation& /*temperature*/,
501  const Evaluation& pressure) const
502  {
503  // ATTENTION: Rs is represented by the _first_ axis!
504  return inverseSaturatedOilBTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
505  }
506 
510  template <class Evaluation>
511  Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
512  const Evaluation& /*temperature*/,
513  const Evaluation& pressure) const
514  { return saturatedGasDissolutionFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true); }
515 
523  template <class Evaluation>
524  Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
525  const Evaluation& /*temperature*/,
526  const Evaluation& pressure,
527  const Evaluation& oilSaturation,
528  Evaluation maxOilSaturation) const
529  {
530  Evaluation tmp =
531  saturatedGasDissolutionFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
532 
533  // apply the vaporization parameters for the gas phase (cf. the Eclipse VAPPARS
534  // keyword)
535  maxOilSaturation = min(maxOilSaturation, Scalar(1.0));
536  if (vapPar2_ > 0.0 && maxOilSaturation > 0.01 && oilSaturation < maxOilSaturation) {
537  static const Scalar eps = 0.001;
538  const Evaluation& So = max(oilSaturation, eps);
539  tmp *= max(1e-3, pow(So/maxOilSaturation, vapPar2_));
540  }
541 
542  return tmp;
543  }
544 
551  template <class Evaluation>
552  Evaluation saturationPressure(unsigned regionIdx,
553  const Evaluation&,
554  const Evaluation& Rs) const
555  {
556  typedef MathToolbox<Evaluation> Toolbox;
557 
558  const auto& RsTable = saturatedGasDissolutionFactorTable_[regionIdx];
559  const Scalar eps = std::numeric_limits<typename Toolbox::Scalar>::epsilon()*1e6;
560 
561  // use the saturation pressure function to get a pretty good initial value
562  Evaluation pSat = saturationPressure_[regionIdx].eval(Rs, /*extrapolate=*/true);
563 
564  // Newton method to do the remaining work. If the initial
565  // value is good, this should only take two to three
566  // iterations...
567  bool onProbation = false;
568  for (int i = 0; i < 20; ++i) {
569  const Evaluation& f = RsTable.eval(pSat, /*extrapolate=*/true) - Rs;
570  const Evaluation& fPrime = RsTable.evalDerivative(pSat, /*extrapolate=*/true);
571 
572  // If the derivative is "zero" Newton will not converge,
573  // so simply return our initial guess.
574  if (std::abs(scalarValue(fPrime)) < 1.0e-30) {
575  return pSat;
576  }
577 
578  const Evaluation& delta = f/fPrime;
579 
580  pSat -= delta;
581 
582  if (pSat < 0.0) {
583  // if the pressure is lower than 0 Pascals, we set it back to 0. if this
584  // happens twice, we give up and just return 0 Pa...
585  if (onProbation)
586  return 0.0;
587 
588  onProbation = true;
589  pSat = 0.0;
590  }
591 
592  if (std::abs(scalarValue(delta)) < std::abs(scalarValue(pSat))*eps)
593  return pSat;
594  }
595 
596  std::stringstream errlog;
597  errlog << "Finding saturation pressure did not converge:"
598  << " pSat = " << pSat
599  << ", Rs = " << Rs;
600 #if HAVE_OPM_COMMON
601  OpmLog::debug("Live oil saturation pressure", errlog.str());
602 #endif
603  throw NumericalIssue(errlog.str());
604  }
605 
606  template <class Evaluation>
607  Evaluation diffusionCoefficient(const Evaluation& /*temperature*/,
608  const Evaluation& /*pressure*/,
609  unsigned /*compIdx*/) const
610  {
611  throw std::runtime_error("Not implemented: The PVT model does not provide a diffusionCoefficient()");
612  }
613  const Scalar gasReferenceDensity(unsigned regionIdx) const
614  { return gasReferenceDensity_[regionIdx]; }
615 
616  const Scalar oilReferenceDensity(unsigned regionIdx) const
617  { return oilReferenceDensity_[regionIdx]; }
618 
619  const std::vector<TabulatedTwoDFunction>& inverseOilBTable() const
620  { return inverseOilBTable_; }
621 
622  const std::vector<TabulatedTwoDFunction>& oilMuTable() const
623  { return oilMuTable_; }
624 
625  const std::vector<TabulatedTwoDFunction>& inverseOilBMuTable() const
626  { return inverseOilBMuTable_; }
627 
628  const std::vector<TabulatedOneDFunction>& saturatedOilMuTable() const
629  { return saturatedOilMuTable_; }
630 
631  const std::vector<TabulatedOneDFunction>& inverseSaturatedOilBTable() const
632  { return inverseSaturatedOilBTable_; }
633 
634  const std::vector<TabulatedOneDFunction>& inverseSaturatedOilBMuTable() const
635  { return inverseSaturatedOilBMuTable_; }
636 
637  const std::vector<TabulatedOneDFunction>& saturatedGasDissolutionFactorTable() const
638  { return saturatedGasDissolutionFactorTable_; }
639 
640  const std::vector<TabulatedOneDFunction>& saturationPressure() const
641  { return saturationPressure_; }
642 
643  Scalar vapPar2() const
644  { return vapPar2_; }
645 
646  bool operator==(const LiveOilPvt<Scalar>& data) const
647  {
648  return this->gasReferenceDensity_ == data.gasReferenceDensity_ &&
649  this->oilReferenceDensity_ == data.oilReferenceDensity_ &&
650  this->inverseOilBTable() == data.inverseOilBTable() &&
651  this->oilMuTable() == data.oilMuTable() &&
652  this->inverseOilBMuTable() == data.inverseOilBMuTable() &&
653  this->saturatedOilMuTable() == data.saturatedOilMuTable() &&
654  this->inverseSaturatedOilBTable() == data.inverseSaturatedOilBTable() &&
655  this->inverseSaturatedOilBMuTable() == data.inverseSaturatedOilBMuTable() &&
656  this->saturatedGasDissolutionFactorTable() == data.saturatedGasDissolutionFactorTable() &&
657  this->vapPar2() == data.vapPar2();
658  }
659 
660 private:
661  void updateSaturationPressure_(unsigned regionIdx)
662  {
663  typedef std::pair<Scalar, Scalar> Pair;
664  const auto& gasDissolutionFac = saturatedGasDissolutionFactorTable_[regionIdx];
665 
666  // create the function representing saturation pressure depending of the mass
667  // fraction in gas
668  size_t n = gasDissolutionFac.numSamples();
669  Scalar delta = (gasDissolutionFac.xMax() - gasDissolutionFac.xMin())/Scalar(n + 1);
670 
671  SamplingPoints pSatSamplePoints;
672  Scalar Rs = 0;
673  for (size_t i=0; i <= n; ++ i) {
674  Scalar pSat = gasDissolutionFac.xMin() + Scalar(i)*delta;
675  Rs = saturatedGasDissolutionFactor(regionIdx,
676  /*temperature=*/Scalar(1e30),
677  pSat);
678 
679  Pair val(Rs, pSat);
680  pSatSamplePoints.push_back(val);
681  }
682 
683  //Prune duplicate Rs values (can occur, and will cause problems in further interpolation)
684  auto x_coord_comparator = [](const Pair& a, const Pair& b) { return a.first == b.first; };
685  auto last = std::unique(pSatSamplePoints.begin(), pSatSamplePoints.end(), x_coord_comparator);
686  pSatSamplePoints.erase(last, pSatSamplePoints.end());
687 
688  saturationPressure_[regionIdx].setContainerOfTuples(pSatSamplePoints);
689  }
690 
691  std::vector<Scalar> gasReferenceDensity_;
692  std::vector<Scalar> oilReferenceDensity_;
693  std::vector<TabulatedTwoDFunction> inverseOilBTable_;
694  std::vector<TabulatedTwoDFunction> oilMuTable_;
695  std::vector<TabulatedTwoDFunction> inverseOilBMuTable_;
696  std::vector<TabulatedOneDFunction> saturatedOilMuTable_;
697  std::vector<TabulatedOneDFunction> inverseSaturatedOilBTable_;
698  std::vector<TabulatedOneDFunction> inverseSaturatedOilBMuTable_;
699  std::vector<TabulatedOneDFunction> saturatedGasDissolutionFactorTable_;
700  std::vector<TabulatedOneDFunction> saturationPressure_;
701 
702  Scalar vapPar2_;
703 };
704 
705 } // namespace Opm
706 
707 #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 oil phas with dissolved gas.
Definition: LiveOilPvt.hpp:54
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of oil given a set of parameters.
Definition: LiveOilPvt.hpp:443
void setSaturatedOilFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the oil formation volume factor.
Definition: LiveOilPvt.hpp:305
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rs) const
Returns the formation volume factor [-] of the fluid phase.
Definition: LiveOilPvt.hpp:486
Evaluation saturatedGasDissolutionFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &oilSaturation, Evaluation maxOilSaturation) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: LiveOilPvt.hpp:524
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: LiveOilPvt.hpp:499
void setSaturatedOilViscosity(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the phase viscosity for gas saturated oil.
Definition: LiveOilPvt.hpp:359
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &, const Evaluation &Rs) const
Returns the saturation pressure of the oil phase [Pa] depending on its mass fraction of the gas compo...
Definition: LiveOilPvt.hpp:552
void setSaturatedOilGasDissolutionFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the gas dissolution factor .
Definition: LiveOilPvt.hpp:293
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rs) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: LiveOilPvt.hpp:455
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar rhoRefGas, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: LiveOilPvt.hpp:279
void setOilViscosity(unsigned regionIdx, const TabulatedTwoDFunction &muo)
Initialize the viscosity of the oil phase.
Definition: LiveOilPvt.hpp:349
Evaluation saturatedGasDissolutionFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: LiveOilPvt.hpp:511
Evaluation saturatedViscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: LiveOilPvt.hpp:471
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: LiveOilPvt.hpp:436
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: LiveOilPvt.hpp:387
void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction &invBo)
Initialize the function for the oil formation volume factor.
Definition: LiveOilPvt.hpp:341
Definition: Exceptions.hpp:46
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47
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