27 #ifndef OPM_WET_HUMID_GAS_PVT_HPP
28 #define OPM_WET_HUMID_GAS_PVT_HPP
35 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
36 #include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
41 #include <opm/common/OpmLog/OpmLog.hpp>
50 template <
class Scalar>
53 using SamplingPoints = std::vector<std::pair<Scalar, Scalar>>;
65 const std::vector<Scalar>& oilReferenceDensity,
66 const std::vector<Scalar>& waterReferenceDensity,
67 const std::vector<TabulatedTwoDFunction>& inverseGasBRvwSat,
68 const std::vector<TabulatedTwoDFunction>& inverseGasBRvSat,
69 const std::vector<TabulatedOneDFunction>& inverseSaturatedGasB,
70 const std::vector<TabulatedTwoDFunction>& gasMuRvwSat,
71 const std::vector<TabulatedTwoDFunction>& gasMuRvSat,
72 const std::vector<TabulatedTwoDFunction>& inverseGasBMuRvwSat,
73 const std::vector<TabulatedTwoDFunction>& inverseGasBMuRvSat,
74 const std::vector<TabulatedOneDFunction>& inverseSaturatedGasBMu,
75 const std::vector<TabulatedOneDFunction>& saturatedWaterVaporizationFactorTable,
76 const std::vector<TabulatedOneDFunction>& saturatedOilVaporizationFactorTable,
79 : gasReferenceDensity_(gasReferenceDensity)
80 , oilReferenceDensity_(oilReferenceDensity)
81 , waterReferenceDensity_(waterReferenceDensity)
82 , inverseGasBRvwSat_(inverseGasBRvwSat)
83 , inverseGasBRvSat_(inverseGasBRvSat)
84 , inverseSaturatedGasB_(inverseSaturatedGasB)
85 , gasMuRvwSat_(gasMuRvwSat)
86 , gasMuRvSat_(gasMuRvSat)
87 , inverseGasBMuRvwSat_(inverseGasBMuRvwSat)
88 , inverseGasBMuRvSat_(inverseGasBMuRvSat)
89 , inverseSaturatedGasBMu_(inverseSaturatedGasBMu)
90 , saturatedWaterVaporizationFactorTable_(saturatedWaterVaporizationFactorTable)
91 , saturatedOilVaporizationFactorTable_(saturatedOilVaporizationFactorTable)
92 , saturationPressure_(saturationPressure)
104 void initFromState(
const EclipseState& eclState,
const Schedule& schedule)
106 const auto& pvtgwTables = eclState.getTableManager().getPvtgwTables();
107 const auto& pvtgTables = eclState.getTableManager().getPvtgTables();
108 const auto& densityTable = eclState.getTableManager().getDensityTable();
110 assert(pvtgwTables.size() == densityTable.size());
111 assert(pvtgTables.size() == densityTable.size());
116 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
117 Scalar rhoRefO = densityTable[regionIdx].oil;
118 Scalar rhoRefG = densityTable[regionIdx].gas;
119 Scalar rhoRefW = densityTable[regionIdx].water;
124 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
125 const auto& pvtgwTable = pvtgwTables[regionIdx];
127 const auto& saturatedTable = pvtgwTable.getSaturatedTable();
128 assert(saturatedTable.numRows() > 1);
131 auto& gasMuRvSat = gasMuRvSat_[regionIdx];
132 auto& invGasBRvSat = inverseGasBRvSat_[regionIdx];
133 auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
134 auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
135 auto& waterVaporizationFac = saturatedWaterVaporizationFactorTable_[regionIdx];
137 waterVaporizationFac.setXYArrays(saturatedTable.numRows(),
138 saturatedTable.getColumn(
"PG"),
139 saturatedTable.getColumn(
"RW"));
141 std::vector<Scalar> invSatGasBArray;
142 std::vector<Scalar> invSatGasBMuArray;
145 for (
unsigned outerIdx = 0; outerIdx < saturatedTable.numRows(); ++ outerIdx) {
146 Scalar pg = saturatedTable.get(
"PG" , outerIdx);
147 Scalar B = saturatedTable.get(
"BG" , outerIdx);
148 Scalar mu = saturatedTable.get(
"MUG" , outerIdx);
150 invGasBRvSat.appendXPos(pg);
151 gasMuRvSat.appendXPos(pg);
153 invSatGasBArray.push_back(1.0/B);
154 invSatGasBMuArray.push_back(1.0/(mu*B));
156 assert(invGasBRvSat.numX() == outerIdx + 1);
157 assert(gasMuRvSat.numX() == outerIdx + 1);
159 const auto& underSaturatedTable = pvtgwTable.getUnderSaturatedTable(outerIdx);
160 size_t numRows = underSaturatedTable.numRows();
161 for (
size_t innerIdx = 0; innerIdx < numRows; ++ innerIdx) {
162 Scalar Rw = underSaturatedTable.get(
"RW" , innerIdx);
163 Scalar Bg = underSaturatedTable.get(
"BG" , innerIdx);
164 Scalar mug = underSaturatedTable.get(
"MUG" , innerIdx);
166 invGasBRvSat.appendSamplePoint(outerIdx, Rw, 1.0/Bg);
167 gasMuRvSat.appendSamplePoint(outerIdx, Rw, mug);
172 std::vector<double> tmpPressure = saturatedTable.getColumn(
"PG").vectorCopy( );
174 invSatGasB.setXYContainers(tmpPressure, invSatGasBArray);
175 invSatGasBMu.setXYContainers(tmpPressure, invSatGasBMuArray);
179 for (
unsigned xIdx = 0; xIdx < invGasBRvSat.numX(); ++xIdx) {
181 assert(invGasBRvSat.numY(xIdx) > 0);
185 if (invGasBRvSat.numY(xIdx) > 1)
191 size_t masterTableIdx = xIdx + 1;
192 for (; masterTableIdx < saturatedTable.numRows(); ++masterTableIdx)
194 if (pvtgwTable.getUnderSaturatedTable(masterTableIdx).numRows() > 1)
198 if (masterTableIdx >= saturatedTable.numRows())
199 throw std::runtime_error(
"PVTGW tables are invalid: The last table must exhibit at least one "
200 "entry for undersaturated gas!");
204 extendPvtgwTable_(regionIdx,
206 pvtgwTable.getUnderSaturatedTable(xIdx),
207 pvtgwTable.getUnderSaturatedTable(masterTableIdx));
212 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
213 const auto& pvtgTable = pvtgTables[regionIdx];
215 const auto& saturatedTable = pvtgTable.getSaturatedTable();
216 assert(saturatedTable.numRows() > 1);
218 auto& gasMuRvwSat = gasMuRvwSat_[regionIdx];
219 auto& invGasBRvwSat = inverseGasBRvwSat_[regionIdx];
220 auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
221 auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
222 auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
224 oilVaporizationFac.setXYArrays(saturatedTable.numRows(),
225 saturatedTable.getColumn(
"PG"),
226 saturatedTable.getColumn(
"RV"));
228 std::vector<Scalar> invSatGasBArray;
229 std::vector<Scalar> invSatGasBMuArray;
232 for (
unsigned outerIdx = 0; outerIdx < saturatedTable.numRows(); ++ outerIdx) {
233 Scalar pg = saturatedTable.get(
"PG" , outerIdx);
234 Scalar B = saturatedTable.get(
"BG" , outerIdx);
235 Scalar mu = saturatedTable.get(
"MUG" , outerIdx);
237 invGasBRvwSat.appendXPos(pg);
238 gasMuRvwSat.appendXPos(pg);
240 invSatGasBArray.push_back(1.0/B);
241 invSatGasBMuArray.push_back(1.0/(mu*B));
243 assert(invGasBRvwSat.numX() == outerIdx + 1);
244 assert(gasMuRvwSat.numX() == outerIdx + 1);
246 const auto& underSaturatedTable = pvtgTable.getUnderSaturatedTable(outerIdx);
247 size_t numRows = underSaturatedTable.numRows();
248 for (
size_t innerIdx = 0; innerIdx < numRows; ++ innerIdx) {
249 Scalar Rv = underSaturatedTable.get(
"RV" , innerIdx);
250 Scalar Bg = underSaturatedTable.get(
"BG" , innerIdx);
251 Scalar mug = underSaturatedTable.get(
"MUG" , innerIdx);
253 invGasBRvwSat.appendSamplePoint(outerIdx, Rv, 1.0/Bg);
254 gasMuRvwSat.appendSamplePoint(outerIdx, Rv, mug);
259 std::vector<double> tmpPressure = saturatedTable.getColumn(
"PG").vectorCopy( );
261 invSatGasB.setXYContainers(tmpPressure, invSatGasBArray);
262 invSatGasBMu.setXYContainers(tmpPressure, invSatGasBMuArray);
266 for (
unsigned xIdx = 0; xIdx < invGasBRvwSat.numX(); ++xIdx) {
268 assert(invGasBRvwSat.numY(xIdx) > 0);
272 if (invGasBRvwSat.numY(xIdx) > 1)
278 size_t masterTableIdx = xIdx + 1;
279 for (; masterTableIdx < saturatedTable.numRows(); ++masterTableIdx)
281 if (pvtgTable.getUnderSaturatedTable(masterTableIdx).numRows() > 1)
285 if (masterTableIdx >= saturatedTable.numRows())
286 throw std::runtime_error(
"PVTG tables are invalid: The last table must exhibit at least one "
287 "entry for undersaturated gas!");
291 extendPvtgTable_(regionIdx,
293 pvtgTable.getUnderSaturatedTable(xIdx),
294 pvtgTable.getUnderSaturatedTable(masterTableIdx));
298 const auto& oilVap = schedule[0].oilvap();
299 if (oilVap.getType() == OilVaporizationProperties::OilVaporization::VAPPARS) {
300 vapPar1_ = oilVap.vap1();
307 void extendPvtgwTable_(
unsigned regionIdx,
309 const SimpleTable& curTable,
310 const SimpleTable& masterTable)
312 std::vector<double> RvArray = curTable.getColumn(
"RW").vectorCopy();
313 std::vector<double> gasBArray = curTable.getColumn(
"BG").vectorCopy();
314 std::vector<double> gasMuArray = curTable.getColumn(
"MUG").vectorCopy();
316 auto& invGasBRvSat = inverseGasBRvSat_[regionIdx];
317 auto& gasMuRvSat = gasMuRvSat_[regionIdx];
319 for (
size_t newRowIdx = 1; newRowIdx < masterTable.numRows(); ++ newRowIdx) {
320 const auto& RVColumn = masterTable.getColumn(
"RW");
321 const auto& BGColumn = masterTable.getColumn(
"BG");
322 const auto& viscosityColumn = masterTable.getColumn(
"MUG");
325 Scalar diffRv = RVColumn[newRowIdx] - RVColumn[newRowIdx - 1];
326 Scalar newRv = RvArray.back() + diffRv;
329 Scalar B1 = BGColumn[newRowIdx];
330 Scalar B2 = BGColumn[newRowIdx - 1];
331 Scalar x = (B1 - B2)/( (B1 + B2)/2.0 );
335 Scalar newBg = gasBArray.back()*(1.0 + x/2.0)/(1.0 - x/2.0);
338 Scalar mu1 = viscosityColumn[newRowIdx];
339 Scalar mu2 = viscosityColumn[newRowIdx - 1];
340 Scalar xMu = (mu1 - mu2)/( (mu1 + mu2)/2.0 );
344 Scalar newMug = gasMuArray.back()*(1.0 + xMu/2)/(1.0 - xMu/2.0);
348 RvArray.push_back(newRv);
349 gasBArray.push_back(newBg);
350 gasMuArray.push_back(newMug);
353 invGasBRvSat.appendSamplePoint(xIdx, newRv, 1.0/newBg);
354 gasMuRvSat.appendSamplePoint(xIdx, newRv, newMug);
357 void extendPvtgTable_(
unsigned regionIdx,
359 const SimpleTable& curTable,
360 const SimpleTable& masterTable)
362 std::vector<double> RvArray = curTable.getColumn(
"RV").vectorCopy();
363 std::vector<double> gasBArray = curTable.getColumn(
"BG").vectorCopy();
364 std::vector<double> gasMuArray = curTable.getColumn(
"MUG").vectorCopy();
366 auto& invGasBRvwSat= inverseGasBRvwSat_[regionIdx];
367 auto& gasMuRvwSat = gasMuRvwSat_[regionIdx];
369 for (
size_t newRowIdx = 1; newRowIdx < masterTable.numRows(); ++ newRowIdx) {
370 const auto& RVColumn = masterTable.getColumn(
"RV");
371 const auto& BGColumn = masterTable.getColumn(
"BG");
372 const auto& viscosityColumn = masterTable.getColumn(
"MUG");
375 Scalar diffRv = RVColumn[newRowIdx] - RVColumn[newRowIdx - 1];
376 Scalar newRv = RvArray.back() + diffRv;
379 Scalar B1 = BGColumn[newRowIdx];
380 Scalar B2 = BGColumn[newRowIdx - 1];
381 Scalar x = (B1 - B2)/( (B1 + B2)/2.0 );
385 Scalar newBg = gasBArray.back()*(1.0 + x/2.0)/(1.0 - x/2.0);
388 Scalar mu1 = viscosityColumn[newRowIdx];
389 Scalar mu2 = viscosityColumn[newRowIdx - 1];
390 Scalar xMu = (mu1 - mu2)/( (mu1 + mu2)/2.0 );
394 Scalar newMug = gasMuArray.back()*(1.0 + xMu/2)/(1.0 - xMu/2.0);
398 RvArray.push_back(newRv);
399 gasBArray.push_back(newBg);
400 gasMuArray.push_back(newMug);
403 invGasBRvwSat.appendSamplePoint(xIdx, newRv, 1.0/newBg);
404 gasMuRvwSat.appendSamplePoint(xIdx, newRv, newMug);
424 saturatedWaterVaporizationFactorTable_.resize(
numRegions);
425 saturatedOilVaporizationFactorTable_.resize(
numRegions);
437 waterReferenceDensity_[regionIdx] = rhoRefWater;
438 oilReferenceDensity_[regionIdx] = rhoRefOil;
439 gasReferenceDensity_[regionIdx] = rhoRefGas;
448 { saturatedWaterVaporizationFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
456 { saturatedOilVaporizationFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
468 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
471 const auto& gasMuRvSat = gasMuRvSat_[regionIdx];
472 const auto& invGasBRvSat = inverseGasBRvSat_[regionIdx];
473 assert(gasMuRvSat.numX() == invGasBRvSat.numX());
475 auto& invGasBMuRvSat = inverseGasBMuRvSat_[regionIdx];
476 auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
477 auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
479 std::vector<Scalar> satPressuresArray;
480 std::vector<Scalar> invSatGasBArray;
481 std::vector<Scalar> invSatGasBMuArray;
482 for (
size_t pIdx = 0; pIdx < gasMuRvSat.numX(); ++pIdx) {
483 invGasBMuRvSat.appendXPos(gasMuRvSat.xAt(pIdx));
485 assert(gasMuRvSat.numY(pIdx) == invGasBRvSat.numY(pIdx));
487 size_t numRw = gasMuRvSat.numY(pIdx);
488 for (
size_t RwIdx = 0; RwIdx < numRw; ++RwIdx)
489 invGasBMuRvSat.appendSamplePoint(pIdx,
490 gasMuRvSat.yAt(pIdx, RwIdx),
491 invGasBRvSat.valueAt(pIdx, RwIdx)
492 / gasMuRvSat.valueAt(pIdx, RwIdx));
497 satPressuresArray.push_back(gasMuRvSat.xAt(pIdx));
498 invSatGasBArray.push_back(invGasBRvSat.valueAt(pIdx, numRw - 1));
499 invSatGasBMuArray.push_back(invGasBMuRvSat.valueAt(pIdx, numRw - 1));
502 invSatGasB.setXYContainers(satPressuresArray, invSatGasBArray);
503 invSatGasBMu.setXYContainers(satPressuresArray, invSatGasBMuArray);
509 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
512 const auto& gasMuRvwSat = gasMuRvwSat_[regionIdx];
513 const auto& invGasBRvwSat = inverseGasBRvwSat_[regionIdx];
514 assert(gasMuRvwSat.numX() == invGasBRvwSat.numX());
516 auto& invGasBMuRvwSat = inverseGasBMuRvwSat_[regionIdx];
517 auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
518 auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
520 std::vector<Scalar> satPressuresArray;
521 std::vector<Scalar> invSatGasBArray;
522 std::vector<Scalar> invSatGasBMuArray;
523 for (
size_t pIdx = 0; pIdx < gasMuRvwSat.numX(); ++pIdx) {
524 invGasBMuRvwSat.appendXPos(gasMuRvwSat.xAt(pIdx));
526 assert(gasMuRvwSat.numY(pIdx) == invGasBRvwSat.numY(pIdx));
528 size_t numRw = gasMuRvwSat.numY(pIdx);
529 for (
size_t RwIdx = 0; RwIdx < numRw; ++RwIdx)
530 invGasBMuRvwSat.appendSamplePoint(pIdx,
531 gasMuRvwSat.yAt(pIdx, RwIdx),
532 invGasBRvwSat.valueAt(pIdx, RwIdx)
533 / gasMuRvwSat.valueAt(pIdx, RwIdx));
538 satPressuresArray.push_back(gasMuRvwSat.xAt(pIdx));
539 invSatGasBArray.push_back(invGasBRvwSat.valueAt(pIdx, numRw - 1));
540 invSatGasBMuArray.push_back(invGasBMuRvwSat.valueAt(pIdx, numRw - 1));
543 invSatGasB.setXYContainers(satPressuresArray, invSatGasBArray);
544 invSatGasBMu.setXYContainers(satPressuresArray, invSatGasBMuArray);
546 updateSaturationPressure_(regionIdx);
554 {
return gasReferenceDensity_.size(); }
559 template <
class Evaluation>
563 const Evaluation&)
const
565 throw std::runtime_error(
"Requested the enthalpy of gas but the thermal option is not enabled");
571 template <
class Evaluation>
574 const Evaluation& pressure,
575 const Evaluation& Rv,
576 const Evaluation& Rvw)
const
578 const Evaluation& temperature = 1E30;
580 if (Rv >= (1.0 - 1e-10)*saturatedOilVaporizationFactor(regionIdx, temperature, pressure)) {
581 const Evaluation& invBg = inverseGasBRvSat_[regionIdx].eval(pressure, Rvw,
true);
582 const Evaluation& invMugBg = inverseGasBMuRvSat_[regionIdx].eval(pressure, Rvw,
true);
583 return invBg/invMugBg;
587 const Evaluation& invBg = inverseGasBRvwSat_[regionIdx].eval(pressure, Rv,
true);
588 const Evaluation& invMugBg = inverseGasBMuRvwSat_[regionIdx].eval(pressure, Rv,
true);
589 return invBg/invMugBg;
596 template <
class Evaluation>
599 const Evaluation& pressure)
const
601 const Evaluation& invBg = inverseSaturatedGasB_[regionIdx].eval(pressure,
true);
602 const Evaluation& invMugBg = inverseSaturatedGasBMu_[regionIdx].eval(pressure,
true);
604 return invBg/invMugBg;
617 template <
class Evaluation>
620 const Evaluation& pressure,
621 const Evaluation& Rv,
622 const Evaluation& Rvw)
const
624 const Evaluation& temperature = 1E30;
626 if (Rv >= (1.0 - 1e-10)*saturatedOilVaporizationFactor(regionIdx, temperature, pressure)) {
627 return inverseGasBRvSat_[regionIdx].eval(pressure, Rvw,
true);
631 return inverseGasBRvwSat_[regionIdx].eval(pressure, Rv,
true);
640 template <
class Evaluation>
643 const Evaluation& pressure)
const
644 {
return inverseSaturatedGasB_[regionIdx].eval(pressure,
true); }
649 template <
class Evaluation>
652 const Evaluation& pressure)
const
654 return saturatedWaterVaporizationFactorTable_[regionIdx].eval(pressure,
true);
657 template <
class Evaluation>
658 Evaluation saturatedOilVaporizationFactor(
unsigned regionIdx,
660 const Evaluation& pressure)
const
662 return saturatedOilVaporizationFactorTable_[regionIdx].eval(pressure,
true);
672 template <
class Evaluation>
675 const Evaluation& pressure,
676 const Evaluation& oilSaturation,
677 Evaluation maxOilSaturation)
const
680 saturatedOilVaporizationFactorTable_[regionIdx].eval(pressure,
true);
684 maxOilSaturation = min(maxOilSaturation, Scalar(1.0));
685 if (vapPar1_ > 0.0 && maxOilSaturation > 0.01 && oilSaturation < maxOilSaturation) {
686 constexpr
const Scalar eps = 0.001;
687 const Evaluation& So = max(oilSaturation, eps);
688 tmp *= max(1e-3, pow(So/maxOilSaturation, vapPar1_));
702 template <
class Evaluation>
705 const Evaluation& Rw)
const
709 const auto& RwTable = saturatedWaterVaporizationFactorTable_[regionIdx];
710 constexpr
const Scalar eps = std::numeric_limits<typename Toolbox::Scalar>::epsilon()*1e6;
713 Evaluation pSat = saturationPressure_[regionIdx].eval(Rw,
true);
718 bool onProbation =
false;
719 for (
unsigned i = 0; i < 20; ++i) {
720 const Evaluation& f = RwTable.eval(pSat,
true) - Rw;
721 const Evaluation& fPrime = RwTable.evalDerivative(pSat,
true);
725 if (std::abs(scalarValue(fPrime)) < 1.0e-30) {
729 const Evaluation& delta = f/fPrime;
743 if (std::abs(scalarValue(delta)) < std::abs(scalarValue(pSat))*eps)
747 std::stringstream errlog;
748 errlog <<
"Finding saturation pressure did not converge:"
749 <<
" pSat = " << pSat
752 OpmLog::debug(
"Wet gas saturation pressure", errlog.str());
757 template <
class Evaluation>
758 Evaluation diffusionCoefficient(
const Evaluation& ,
762 throw std::runtime_error(
"Not implemented: The PVT model does not provide a diffusionCoefficient()");
765 const Scalar gasReferenceDensity(
unsigned regionIdx)
const
766 {
return gasReferenceDensity_[regionIdx]; }
768 const Scalar oilReferenceDensity(
unsigned regionIdx)
const
769 {
return oilReferenceDensity_[regionIdx]; }
771 const Scalar waterReferenceDensity(
unsigned regionIdx)
const
772 {
return waterReferenceDensity_[regionIdx]; }
774 const std::vector<TabulatedTwoDFunction>& inverseGasB()
const {
775 return inverseGasBRvSat_;
778 const std::vector<TabulatedOneDFunction>& inverseSaturatedGasB()
const {
779 return inverseSaturatedGasB_;
782 const std::vector<TabulatedTwoDFunction>& gasMu()
const {
786 const std::vector<TabulatedTwoDFunction>& inverseGasBMu()
const {
787 return inverseGasBMuRvSat_;
790 const std::vector<TabulatedOneDFunction>& inverseSaturatedGasBMu()
const {
791 return inverseSaturatedGasBMu_;
794 const std::vector<TabulatedOneDFunction>& saturatedWaterVaporizationFactorTable()
const {
795 return saturatedWaterVaporizationFactorTable_;
798 const std::vector<TabulatedOneDFunction>& saturatedOilVaporizationFactorTable()
const {
799 return saturatedOilVaporizationFactorTable_;
803 const std::vector<TabulatedOneDFunction>& saturationPressure()
const {
804 return saturationPressure_;
807 Scalar vapPar1()
const {
811 bool operator==(
const WetHumidGasPvt<Scalar>& data)
const
813 return this->gasReferenceDensity_ == data.gasReferenceDensity_ &&
814 this->oilReferenceDensity_ == data.oilReferenceDensity_ &&
815 this->waterReferenceDensity_ == data.waterReferenceDensity_ &&
816 this->inverseGasB() == data.inverseGasB() &&
817 this->inverseSaturatedGasB() == data.inverseSaturatedGasB() &&
818 this->gasMu() == data.gasMu() &&
819 this->inverseGasBMu() == data.inverseGasBMu() &&
820 this->inverseSaturatedGasBMu() == data.inverseSaturatedGasBMu() &&
821 this->saturatedWaterVaporizationFactorTable() == data.saturatedWaterVaporizationFactorTable() &&
822 this->saturationPressure() == data.saturationPressure() &&
823 this->vapPar1() == data.vapPar1();
827 void updateSaturationPressure_(
unsigned regionIdx)
829 const auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
833 size_t n = oilVaporizationFac.numSamples();
834 Scalar delta = (oilVaporizationFac.xMax() - oilVaporizationFac.xMin())/Scalar(n + 1);
836 SamplingPoints pSatSamplePoints;
838 for (
size_t i = 0; i <= n; ++ i) {
839 Scalar pSat = oilVaporizationFac.xMin() + Scalar(i)*delta;
840 Rv = saturatedOilVaporizationFactor(regionIdx, Scalar(1e30), pSat);
842 pSatSamplePoints.emplace_back(Rv, pSat);
846 auto x_coord_comparator = [](
const auto& a,
const auto& b) {
return a.first == b.first; };
847 auto last = std::unique(pSatSamplePoints.begin(), pSatSamplePoints.end(), x_coord_comparator);
848 pSatSamplePoints.erase(last, pSatSamplePoints.end());
850 saturationPressure_[regionIdx].setContainerOfTuples(pSatSamplePoints);
853 std::vector<Scalar> gasReferenceDensity_;
854 std::vector<Scalar> waterReferenceDensity_;
855 std::vector<Scalar> oilReferenceDensity_;
856 std::vector<TabulatedTwoDFunction> inverseGasBRvwSat_;
857 std::vector<TabulatedTwoDFunction> inverseGasBRvSat_;
858 std::vector<TabulatedOneDFunction> inverseSaturatedGasB_;
859 std::vector<TabulatedTwoDFunction> gasMuRvwSat_;
860 std::vector<TabulatedTwoDFunction> gasMuRvSat_;
861 std::vector<TabulatedTwoDFunction> inverseGasBMuRvwSat_;
862 std::vector<TabulatedTwoDFunction> inverseGasBMuRvSat_;
863 std::vector<TabulatedOneDFunction> inverseSaturatedGasBMu_;
864 std::vector<TabulatedOneDFunction> saturatedWaterVaporizationFactorTable_;
865 std::vector<TabulatedOneDFunction> saturatedOilVaporizationFactorTable_;
866 std::vector<TabulatedOneDFunction> saturationPressure_;
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Exceptions.hpp:46
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:48
This class represents the Pressure-Volume-Temperature relations of the gas phase with vaporized oil a...
Definition: WetHumidGasPvt.hpp:52
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: WetHumidGasPvt.hpp:703
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: WetHumidGasPvt.hpp:553
void setSaturatedGasOilVaporizationFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the oil vaporization factor .
Definition: WetHumidGasPvt.hpp:455
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rv, const Evaluation &Rvw) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WetHumidGasPvt.hpp:572
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of water saturated gas at a given pressure.
Definition: WetHumidGasPvt.hpp:641
Evaluation saturatedWaterVaporizationFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the water vaporization factor [m^3/m^3] of the water phase.
Definition: WetHumidGasPvt.hpp:650
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: WetHumidGasPvt.hpp:597
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of gas given a set of parameters.
Definition: WetHumidGasPvt.hpp:560
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: WetHumidGasPvt.hpp:673
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar rhoRefGas, Scalar rhoRefWater)
Initialize the reference densities of all fluids for a given PVT region.
Definition: WetHumidGasPvt.hpp:432
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rv, const Evaluation &Rvw) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WetHumidGasPvt.hpp:618
void initEnd()
Finish initializing the gas phase PVT properties.
Definition: WetHumidGasPvt.hpp:462
void setSaturatedGasWaterVaporizationFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the water vaporization factor .
Definition: WetHumidGasPvt.hpp:447