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
33
34#if HAVE_ECL_INPUT
35#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
36#include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
37#include <opm/input/eclipse/Schedule/OilVaporizationProperties.hpp>
38#endif
39
40#if HAVE_OPM_COMMON
41#include <opm/common/OpmLog/OpmLog.hpp>
42#endif
43
44namespace Opm {
45
50template <class Scalar>
52{
53 using SamplingPoints = std::vector<std::pair<Scalar, Scalar>>;
54
55public:
58
59 WetGasPvt()
60 {
61 vapPar1_ = 0.0;
62 }
63
64 WetGasPvt(const std::vector<Scalar>& gasReferenceDensity,
65 const std::vector<Scalar>& oilReferenceDensity,
66 const std::vector<TabulatedTwoDFunction>& inverseGasB,
67 const std::vector<TabulatedOneDFunction>& inverseSaturatedGasB,
68 const std::vector<TabulatedTwoDFunction>& gasMu,
69 const std::vector<TabulatedTwoDFunction>& inverseGasBMu,
70 const std::vector<TabulatedOneDFunction>& inverseSaturatedGasBMu,
71 const std::vector<TabulatedOneDFunction>& saturatedOilVaporizationFactorTable,
72 const std::vector<TabulatedOneDFunction>& saturationPressure,
73 Scalar vapPar1)
74 : gasReferenceDensity_(gasReferenceDensity)
75 , oilReferenceDensity_(oilReferenceDensity)
76 , inverseGasB_(inverseGasB)
77 , inverseSaturatedGasB_(inverseSaturatedGasB)
78 , gasMu_(gasMu)
79 , inverseGasBMu_(inverseGasBMu)
80 , inverseSaturatedGasBMu_(inverseSaturatedGasBMu)
81 , saturatedOilVaporizationFactorTable_(saturatedOilVaporizationFactorTable)
82 , saturationPressure_(saturationPressure)
83 , vapPar1_(vapPar1)
84 {
85 }
86
87
88#if HAVE_ECL_INPUT
94 void initFromState(const EclipseState& eclState, const Schedule& schedule)
95 {
96 const auto& pvtgTables = eclState.getTableManager().getPvtgTables();
97 const auto& densityTable = eclState.getTableManager().getDensityTable();
98
99 assert(pvtgTables.size() == densityTable.size());
100
101 size_t numRegions = pvtgTables.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 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
113 const auto& pvtgTable = pvtgTables[regionIdx];
114
115 const auto& saturatedTable = pvtgTable.getSaturatedTable();
116 assert(saturatedTable.numRows() > 1);
117
118 auto& gasMu = gasMu_[regionIdx];
119 auto& invGasB = inverseGasB_[regionIdx];
120 auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
121 auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
122 auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
123
124 oilVaporizationFac.setXYArrays(saturatedTable.numRows(),
125 saturatedTable.getColumn("PG"),
126 saturatedTable.getColumn("RV"));
127
128 std::vector<Scalar> invSatGasBArray;
129 std::vector<Scalar> invSatGasBMuArray;
130
131 // extract the table for the gas dissolution and the oil formation volume factors
132 for (unsigned outerIdx = 0; outerIdx < saturatedTable.numRows(); ++ outerIdx) {
133 Scalar pg = saturatedTable.get("PG" , outerIdx);
134 Scalar B = saturatedTable.get("BG" , outerIdx);
135 Scalar mu = saturatedTable.get("MUG" , outerIdx);
136
137 invGasB.appendXPos(pg);
138 gasMu.appendXPos(pg);
139
140 invSatGasBArray.push_back(1.0/B);
141 invSatGasBMuArray.push_back(1.0/(mu*B));
142
143 assert(invGasB.numX() == outerIdx + 1);
144 assert(gasMu.numX() == outerIdx + 1);
145
146 const auto& underSaturatedTable = pvtgTable.getUnderSaturatedTable(outerIdx);
147 size_t numRows = underSaturatedTable.numRows();
148 for (size_t innerIdx = 0; innerIdx < numRows; ++ innerIdx) {
149 Scalar Rv = underSaturatedTable.get("RV" , innerIdx);
150 Scalar Bg = underSaturatedTable.get("BG" , innerIdx);
151 Scalar mug = underSaturatedTable.get("MUG" , innerIdx);
152
153 invGasB.appendSamplePoint(outerIdx, Rv, 1.0/Bg);
154 gasMu.appendSamplePoint(outerIdx, Rv, mug);
155 }
156 }
157
158 {
159 std::vector<double> tmpPressure = saturatedTable.getColumn("PG").vectorCopy( );
160
161 invSatGasB.setXYContainers(tmpPressure, invSatGasBArray);
162 invSatGasBMu.setXYContainers(tmpPressure, invSatGasBMuArray);
163 }
164
165 // make sure to have at least two sample points per gas pressure value
166 for (unsigned xIdx = 0; xIdx < invGasB.numX(); ++xIdx) {
167 // a single sample point is definitely needed
168 assert(invGasB.numY(xIdx) > 0);
169
170 // everything is fine if the current table has two or more sampling points
171 // for a given mole fraction
172 if (invGasB.numY(xIdx) > 1)
173 continue;
174
175 // find the master table which will be used as a template to extend the
176 // current line. We define master table as the first table which has values
177 // for undersaturated gas...
178 size_t masterTableIdx = xIdx + 1;
179 for (; masterTableIdx < saturatedTable.numRows(); ++masterTableIdx)
180 {
181 if (pvtgTable.getUnderSaturatedTable(masterTableIdx).numRows() > 1)
182 break;
183 }
184
185 if (masterTableIdx >= saturatedTable.numRows())
186 throw std::runtime_error("PVTG tables are invalid: The last table must exhibit at least one "
187 "entry for undersaturated gas!");
188
189
190 // extend the current table using the master table.
191 extendPvtgTable_(regionIdx,
192 xIdx,
193 pvtgTable.getUnderSaturatedTable(xIdx),
194 pvtgTable.getUnderSaturatedTable(masterTableIdx));
195 }
196 }
197
198 vapPar1_ = 0.0;
199 const auto& oilVap = schedule[0].oilvap();
200 if (oilVap.getType() == OilVaporizationProperties::OilVaporization::VAPPARS) {
201 vapPar1_ = oilVap.vap1();
202 }
203
204 initEnd();
205 }
206
207private:
208 void extendPvtgTable_(unsigned regionIdx,
209 unsigned xIdx,
210 const SimpleTable& curTable,
211 const SimpleTable& masterTable)
212 {
213 std::vector<double> RvArray = curTable.getColumn("RV").vectorCopy();
214 std::vector<double> gasBArray = curTable.getColumn("BG").vectorCopy();
215 std::vector<double> gasMuArray = curTable.getColumn("MUG").vectorCopy();
216
217 auto& invGasB = inverseGasB_[regionIdx];
218 auto& gasMu = gasMu_[regionIdx];
219
220 for (size_t newRowIdx = 1; newRowIdx < masterTable.numRows(); ++ newRowIdx) {
221 const auto& RVColumn = masterTable.getColumn("RV");
222 const auto& BGColumn = masterTable.getColumn("BG");
223 const auto& viscosityColumn = masterTable.getColumn("MUG");
224
225 // compute the gas pressure for the new entry
226 Scalar diffRv = RVColumn[newRowIdx] - RVColumn[newRowIdx - 1];
227 Scalar newRv = RvArray.back() + diffRv;
228
229 // calculate the compressibility of the master table
230 Scalar B1 = BGColumn[newRowIdx];
231 Scalar B2 = BGColumn[newRowIdx - 1];
232 Scalar x = (B1 - B2)/( (B1 + B2)/2.0 );
233
234 // calculate the gas formation volume factor which exhibits the same
235 // "compressibility" for the new value of Rv
236 Scalar newBg = gasBArray.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 gas formation volume factor which exhibits the same
244 // compressibility for the new pressure
245 Scalar newMug = gasMuArray.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 RvArray.push_back(newRv);
250 gasBArray.push_back(newBg);
251 gasMuArray.push_back(newMug);
252
253 // ... and register them with the internal table objects
254 invGasB.appendSamplePoint(xIdx, newRv, 1.0/newBg);
255 gasMu.appendSamplePoint(xIdx, newRv, newMug);
256 }
257 }
258
259public:
260#endif // HAVE_ECL_INPUT
261
262 void setNumRegions(size_t numRegions)
263 {
264 oilReferenceDensity_.resize(numRegions);
265 gasReferenceDensity_.resize(numRegions);
266 inverseGasB_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::RightExtreme});
267 inverseGasBMu_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::RightExtreme});
268 inverseSaturatedGasB_.resize(numRegions);
269 inverseSaturatedGasBMu_.resize(numRegions);
270 gasMu_.resize(numRegions, TabulatedTwoDFunction{TabulatedTwoDFunction::InterpolationPolicy::RightExtreme});
271 saturatedOilVaporizationFactorTable_.resize(numRegions);
272 saturationPressure_.resize(numRegions);
273 }
274
278 void setReferenceDensities(unsigned regionIdx,
279 Scalar rhoRefOil,
280 Scalar rhoRefGas,
281 Scalar /*rhoRefWater*/)
282 {
283 oilReferenceDensity_[regionIdx] = rhoRefOil;
284 gasReferenceDensity_[regionIdx] = rhoRefGas;
285 }
286
292 void setSaturatedGasOilVaporizationFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
293 { saturatedOilVaporizationFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
294
304 void setSaturatedGasFormationVolumeFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
305 {
306 auto& invGasB = inverseGasB_[regionIdx];
307
308 const auto& RvTable = saturatedOilVaporizationFactorTable_[regionIdx];
309
310 constexpr const Scalar T = 273.15 + 15.56; // [K]
311
312 constexpr const Scalar RvMin = 0.0;
313 Scalar RvMax = RvTable.eval(saturatedOilVaporizationFactorTable_[regionIdx].xMax(), /*extrapolate=*/true);
314
315 Scalar poMin = samplePoints.front().first;
316 Scalar poMax = samplePoints.back().first;
317
318 constexpr const size_t nRv = 20;
319 size_t nP = samplePoints.size()*2;
320
321 Scalar rhooRef = oilReferenceDensity_[regionIdx];
322
323 TabulatedOneDFunction gasFormationVolumeFactor;
324 gasFormationVolumeFactor.setContainerOfTuples(samplePoints);
325
326 updateSaturationPressure_(regionIdx);
327
328 // calculate a table of estimated densities depending on pressure and gas mass
329 // fraction. note that this assumes oil of constant compressibility. (having said
330 // that, if only the saturated gas densities are available, there's not much
331 // choice.)
332 for (size_t RvIdx = 0; RvIdx < nRv; ++RvIdx) {
333 Scalar Rv = RvMin + (RvMax - RvMin)*RvIdx/nRv;
334
335 invGasB.appendXPos(Rv);
336
337 for (size_t pIdx = 0; pIdx < nP; ++pIdx) {
338 Scalar pg = poMin + (poMax - poMin)*pIdx/nP;
339
340 Scalar poSat = saturationPressure(regionIdx, T, Rv);
341 Scalar BgSat = gasFormationVolumeFactor.eval(poSat, /*extrapolate=*/true);
342 Scalar drhoo_dp = (1.1200 - 1.1189)/((5000 - 4000)*6894.76);
343 Scalar rhoo = rhooRef/BgSat*(1 + drhoo_dp*(pg - poSat));
344
345 Scalar Bg = rhooRef/rhoo;
346
347 invGasB.appendSamplePoint(RvIdx, pg, 1.0/Bg);
348 }
349 }
350 }
351
364 void setInverseGasFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction& invBg)
365 { inverseGasB_[regionIdx] = invBg; }
366
372 void setGasViscosity(unsigned regionIdx, const TabulatedTwoDFunction& mug)
373 { gasMu_[regionIdx] = mug; }
374
382 void setSaturatedGasViscosity(unsigned regionIdx, const SamplingPoints& samplePoints )
383 {
384 auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
385
386 constexpr const Scalar RvMin = 0.0;
387 Scalar RvMax = oilVaporizationFac.eval(saturatedOilVaporizationFactorTable_[regionIdx].xMax(), /*extrapolate=*/true);
388
389 Scalar poMin = samplePoints.front().first;
390 Scalar poMax = samplePoints.back().first;
391
392 constexpr const size_t nRv = 20;
393 size_t nP = samplePoints.size()*2;
394
395 TabulatedOneDFunction mugTable;
396 mugTable.setContainerOfTuples(samplePoints);
397
398 // calculate a table of estimated densities depending on pressure and gas mass
399 // fraction
400 for (size_t RvIdx = 0; RvIdx < nRv; ++RvIdx) {
401 Scalar Rv = RvMin + (RvMax - RvMin)*RvIdx/nRv;
402
403 gasMu_[regionIdx].appendXPos(Rv);
404
405 for (size_t pIdx = 0; pIdx < nP; ++pIdx) {
406 Scalar pg = poMin + (poMax - poMin)*pIdx/nP;
407 Scalar mug = mugTable.eval(pg, /*extrapolate=*/true);
408
409 gasMu_[regionIdx].appendSamplePoint(RvIdx, pg, mug);
410 }
411 }
412 }
413
417 void initEnd()
418 {
419 // calculate the final 2D functions which are used for interpolation.
420 size_t numRegions = gasMu_.size();
421 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
422 // calculate the table which stores the inverse of the product of the gas
423 // formation volume factor and the gas viscosity
424 const auto& gasMu = gasMu_[regionIdx];
425 const auto& invGasB = inverseGasB_[regionIdx];
426 assert(gasMu.numX() == invGasB.numX());
427
428 auto& invGasBMu = inverseGasBMu_[regionIdx];
429 auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
430 auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
431
432 std::vector<Scalar> satPressuresArray;
433 std::vector<Scalar> invSatGasBArray;
434 std::vector<Scalar> invSatGasBMuArray;
435 for (size_t pIdx = 0; pIdx < gasMu.numX(); ++pIdx) {
436 invGasBMu.appendXPos(gasMu.xAt(pIdx));
437
438 assert(gasMu.numY(pIdx) == invGasB.numY(pIdx));
439
440 size_t numRv = gasMu.numY(pIdx);
441 for (size_t rvIdx = 0; rvIdx < numRv; ++rvIdx)
442 invGasBMu.appendSamplePoint(pIdx,
443 gasMu.yAt(pIdx, rvIdx),
444 invGasB.valueAt(pIdx, rvIdx)
445 / gasMu.valueAt(pIdx, rvIdx));
446
447 // the sampling points in UniformXTabulated2DFunction are always sorted
448 // in ascending order. Thus, the value for saturated gas is the last one
449 // (i.e., the one with the largest Rv value)
450 satPressuresArray.push_back(gasMu.xAt(pIdx));
451 invSatGasBArray.push_back(invGasB.valueAt(pIdx, numRv - 1));
452 invSatGasBMuArray.push_back(invGasBMu.valueAt(pIdx, numRv - 1));
453 }
454
455 invSatGasB.setXYContainers(satPressuresArray, invSatGasBArray);
456 invSatGasBMu.setXYContainers(satPressuresArray, invSatGasBMuArray);
457
458 updateSaturationPressure_(regionIdx);
459 }
460 }
461
465 unsigned numRegions() const
466 { return gasReferenceDensity_.size(); }
467
471 template <class Evaluation>
472 Evaluation internalEnergy(unsigned,
473 const Evaluation&,
474 const Evaluation&,
475 const Evaluation&) const
476 {
477 throw std::runtime_error("Requested the enthalpy of gas but the thermal option is not enabled");
478 }
479
483 template <class Evaluation>
484 Evaluation viscosity(unsigned regionIdx,
485 const Evaluation& /*temperature*/,
486 const Evaluation& pressure,
487 const Evaluation& Rv,
488 const Evaluation& /*Rvw*/) const
489 {
490 const Evaluation& invBg = inverseGasB_[regionIdx].eval(pressure, Rv, /*extrapolate=*/true);
491 const Evaluation& invMugBg = inverseGasBMu_[regionIdx].eval(pressure, Rv, /*extrapolate=*/true);
492
493 return invBg/invMugBg;
494 }
495
499 template <class Evaluation>
500 Evaluation saturatedViscosity(unsigned regionIdx,
501 const Evaluation& /*temperature*/,
502 const Evaluation& pressure) const
503 {
504 const Evaluation& invBg = inverseSaturatedGasB_[regionIdx].eval(pressure, /*extrapolate=*/true);
505 const Evaluation& invMugBg = inverseSaturatedGasBMu_[regionIdx].eval(pressure, /*extrapolate=*/true);
506
507 return invBg/invMugBg;
508 }
509
513 template <class Evaluation>
514 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
515 const Evaluation& /*temperature*/,
516 const Evaluation& pressure,
517 const Evaluation& Rv,
518 const Evaluation& /*Rvw*/) 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 saturatedWaterVaporizationFactor(unsigned /*regionIdx*/,
535 const Evaluation& /*temperature*/,
536 const Evaluation& /*pressure*/) const
537 { return 0.0; /* this is non-humid gas! */ }
538
542 template <class Evaluation>
543 Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
544 const Evaluation& /*temperature*/,
545 const Evaluation& pressure) const
546 {
547 return saturatedOilVaporizationFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
548 }
549
557 template <class Evaluation>
558 Evaluation saturatedOilVaporizationFactor(unsigned regionIdx,
559 const Evaluation& /*temperature*/,
560 const Evaluation& pressure,
561 const Evaluation& oilSaturation,
562 Evaluation maxOilSaturation) const
563 {
564 Evaluation tmp =
565 saturatedOilVaporizationFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
566
567 // apply the vaporization parameters for the gas phase (cf. the Eclipse VAPPARS
568 // keyword)
569 maxOilSaturation = min(maxOilSaturation, Scalar(1.0));
570 if (vapPar1_ > 0.0 && maxOilSaturation > 0.01 && oilSaturation < maxOilSaturation) {
571 constexpr const Scalar eps = 0.001;
572 const Evaluation& So = max(oilSaturation, eps);
573 tmp *= max(1e-3, pow(So/maxOilSaturation, vapPar1_));
574 }
575
576 return tmp;
577 }
578
589 template <class Evaluation>
590 Evaluation saturationPressure(unsigned regionIdx,
591 const Evaluation&,
592 const Evaluation& Rv) const
593 {
594 typedef MathToolbox<Evaluation> Toolbox;
595
596 const auto& RvTable = saturatedOilVaporizationFactorTable_[regionIdx];
597 constexpr const Scalar eps = std::numeric_limits<typename Toolbox::Scalar>::epsilon()*1e6;
598
599 // use the tabulated saturation pressure function to get a pretty good initial value
600 Evaluation pSat = saturationPressure_[regionIdx].eval(Rv, /*extrapolate=*/true);
601
602 // Newton method to do the remaining work. If the initial
603 // value is good, this should only take two to three
604 // iterations...
605 bool onProbation = false;
606 for (unsigned i = 0; i < 20; ++i) {
607 const Evaluation& f = RvTable.eval(pSat, /*extrapolate=*/true) - Rv;
608 const Evaluation& fPrime = RvTable.evalDerivative(pSat, /*extrapolate=*/true);
609
610 // If the derivative is "zero" Newton will not converge,
611 // so simply return our initial guess.
612 if (std::abs(scalarValue(fPrime)) < 1.0e-30) {
613 return pSat;
614 }
615
616 const Evaluation& delta = f/fPrime;
617
618 pSat -= delta;
619
620 if (pSat < 0.0) {
621 // if the pressure is lower than 0 Pascals, we set it back to 0. if this
622 // happens twice, we give up and just return 0 Pa...
623 if (onProbation)
624 return 0.0;
625
626 onProbation = true;
627 pSat = 0.0;
628 }
629
630 if (std::abs(scalarValue(delta)) < std::abs(scalarValue(pSat))*eps)
631 return pSat;
632 }
633
634 std::stringstream errlog;
635 errlog << "Finding saturation pressure did not converge:"
636 << " pSat = " << pSat
637 << ", Rv = " << Rv;
638#if HAVE_OPM_COMMON
639 OpmLog::debug("Wet gas saturation pressure", errlog.str());
640#endif
641 throw NumericalIssue(errlog.str());
642 }
643
644 template <class Evaluation>
645 Evaluation diffusionCoefficient(const Evaluation& /*temperature*/,
646 const Evaluation& /*pressure*/,
647 unsigned /*compIdx*/) const
648 {
649 throw std::runtime_error("Not implemented: The PVT model does not provide a diffusionCoefficient()");
650 }
651
652 const Scalar gasReferenceDensity(unsigned regionIdx) const
653 { return gasReferenceDensity_[regionIdx]; }
654
655 const Scalar oilReferenceDensity(unsigned regionIdx) const
656 { return oilReferenceDensity_[regionIdx]; }
657
658 const std::vector<TabulatedTwoDFunction>& inverseGasB() const {
659 return inverseGasB_;
660 }
661
662 const std::vector<TabulatedOneDFunction>& inverseSaturatedGasB() const {
663 return inverseSaturatedGasB_;
664 }
665
666 const std::vector<TabulatedTwoDFunction>& gasMu() const {
667 return gasMu_;
668 }
669
670 const std::vector<TabulatedTwoDFunction>& inverseGasBMu() const {
671 return inverseGasBMu_;
672 }
673
674 const std::vector<TabulatedOneDFunction>& inverseSaturatedGasBMu() const {
675 return inverseSaturatedGasBMu_;
676 }
677
678 const std::vector<TabulatedOneDFunction>& saturatedOilVaporizationFactorTable() const {
679 return saturatedOilVaporizationFactorTable_;
680 }
681
682 const std::vector<TabulatedOneDFunction>& saturationPressure() const {
683 return saturationPressure_;
684 }
685
686 Scalar vapPar1() const {
687 return vapPar1_;
688 }
689
690 bool operator==(const WetGasPvt<Scalar>& data) const
691 {
692 return this->gasReferenceDensity_ == data.gasReferenceDensity_ &&
693 this->oilReferenceDensity_ == data.oilReferenceDensity_ &&
694 this->inverseGasB() == data.inverseGasB() &&
695 this->inverseSaturatedGasB() == data.inverseSaturatedGasB() &&
696 this->gasMu() == data.gasMu() &&
697 this->inverseGasBMu() == data.inverseGasBMu() &&
698 this->inverseSaturatedGasBMu() == data.inverseSaturatedGasBMu() &&
699 this->saturatedOilVaporizationFactorTable() == data.saturatedOilVaporizationFactorTable() &&
700 this->saturationPressure() == data.saturationPressure() &&
701 this->vapPar1() == data.vapPar1();
702 }
703
704private:
705 void updateSaturationPressure_(unsigned regionIdx)
706 {
707 const auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
708
709 // create the taublated function representing saturation pressure depending of
710 // Rv
711 size_t n = oilVaporizationFac.numSamples();
712 Scalar delta = (oilVaporizationFac.xMax() - oilVaporizationFac.xMin())/Scalar(n + 1);
713
714 SamplingPoints pSatSamplePoints;
715 Scalar Rv = 0;
716 for (size_t i = 0; i <= n; ++ i) {
717 Scalar pSat = oilVaporizationFac.xMin() + Scalar(i)*delta;
718 Rv = saturatedOilVaporizationFactor(regionIdx, /*temperature=*/Scalar(1e30), pSat);
719
720 pSatSamplePoints.emplace_back(Rv, pSat);
721 }
722
723 //Prune duplicate Rv values (can occur, and will cause problems in further interpolation)
724 auto x_coord_comparator = [](const auto& a, const auto& b) { return a.first == b.first; };
725 auto last = std::unique(pSatSamplePoints.begin(), pSatSamplePoints.end(), x_coord_comparator);
726 if (std::distance(pSatSamplePoints.begin(), last) > 1) // only remove them if there are more than two points
727 pSatSamplePoints.erase(last, pSatSamplePoints.end());
728
729 saturationPressure_[regionIdx].setContainerOfTuples(pSatSamplePoints);
730 }
731
732 std::vector<Scalar> gasReferenceDensity_;
733 std::vector<Scalar> oilReferenceDensity_;
734 std::vector<TabulatedTwoDFunction> inverseGasB_;
735 std::vector<TabulatedOneDFunction> inverseSaturatedGasB_;
736 std::vector<TabulatedTwoDFunction> gasMu_;
737 std::vector<TabulatedTwoDFunction> inverseGasBMu_;
738 std::vector<TabulatedOneDFunction> inverseSaturatedGasBMu_;
739 std::vector<TabulatedOneDFunction> saturatedOilVaporizationFactorTable_;
740 std::vector<TabulatedOneDFunction> saturationPressure_;
741
742 Scalar vapPar1_;
743};
744
745} // namespace Opm
746
747#endif
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
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: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
This class represents the Pressure-Volume-Temperature relations of the gas phas with vaporized oil.
Definition: WetGasPvt.hpp:52
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:543
void initEnd()
Finish initializing the gas phase PVT properties.
Definition: WetGasPvt.hpp:417
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:590
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:558
void setGasViscosity(unsigned regionIdx, const TabulatedTwoDFunction &mug)
Initialize the viscosity of the gas phase.
Definition: WetGasPvt.hpp:372
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:472
Evaluation saturatedWaterVaporizationFactor(unsigned, const Evaluation &, const Evaluation &) const
Returns the water vaporization factor [m^3/m^3] of the gasphase.
Definition: WetGasPvt.hpp:534
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: WetGasPvt.hpp:465
void setSaturatedGasViscosity(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the phase viscosity for oil saturated gas.
Definition: WetGasPvt.hpp:382
void setSaturatedGasFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the gas formation volume factor.
Definition: WetGasPvt.hpp:304
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rv, const Evaluation &) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WetGasPvt.hpp:514
void setSaturatedGasOilVaporizationFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the oil vaporization factor .
Definition: WetGasPvt.hpp:292
void setInverseGasFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction &invBg)
Initialize the function for the gas formation volume factor.
Definition: WetGasPvt.hpp:364
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:500
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rv, const Evaluation &) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WetGasPvt.hpp:484
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar rhoRefGas, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: WetGasPvt.hpp:278
Definition: MathToolbox.hpp:50