27 #ifndef OPM_WET_GAS_PVT_HPP 28 #define OPM_WET_GAS_PVT_HPP 36 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp> 37 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp> 46 template <
class Scalar>
51 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
65 void initFromDeck(
const Deck& deck,
const EclipseState& eclState)
67 const auto& pvtgTables = eclState.getTableManager().getPvtgTables();
68 const auto& densityKeyword = deck.getKeyword(
"DENSITY");
70 assert(pvtgTables.size() == densityKeyword.size());
75 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
76 Scalar rhoRefO = densityKeyword.getRecord(regionIdx).getItem(
"OIL").getSIDouble(0);
77 Scalar rhoRefG = densityKeyword.getRecord(regionIdx).getItem(
"GAS").getSIDouble(0);
78 Scalar rhoRefW = densityKeyword.getRecord(regionIdx).getItem(
"WATER").getSIDouble(0);
83 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
84 const auto& pvtgTable = pvtgTables[regionIdx];
86 const auto& saturatedTable = pvtgTable.getSaturatedTable();
87 assert(saturatedTable.numRows() > 1);
89 auto& gasMu = gasMu_[regionIdx];
90 auto& invGasB = inverseGasB_[regionIdx];
91 auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
92 auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
93 auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
95 oilVaporizationFac.setXYArrays(saturatedTable.numRows(),
96 saturatedTable.getColumn(
"PG"),
97 saturatedTable.getColumn(
"RV"));
99 std::vector<Scalar> invSatGasBArray;
100 std::vector<Scalar> invSatGasBMuArray;
103 for (
unsigned outerIdx = 0; outerIdx < saturatedTable.numRows(); ++ outerIdx) {
104 Scalar pg = saturatedTable.get(
"PG" , outerIdx);
105 Scalar B = saturatedTable.get(
"BG" , outerIdx);
106 Scalar mu = saturatedTable.get(
"MUG" , outerIdx);
108 invGasB.appendXPos(pg);
109 gasMu.appendXPos(pg);
111 invSatGasBArray.push_back(1.0/B);
112 invSatGasBMuArray.push_back(1.0/(mu*B));
114 assert(invGasB.numX() == outerIdx + 1);
115 assert(gasMu.numX() == outerIdx + 1);
117 const auto& underSaturatedTable = pvtgTable.getUnderSaturatedTable(outerIdx);
118 size_t numRows = underSaturatedTable.numRows();
119 for (
size_t innerIdx = 0; innerIdx < numRows; ++ innerIdx) {
120 Scalar Rv = underSaturatedTable.get(
"RV" , innerIdx);
121 Scalar Bg = underSaturatedTable.get(
"BG" , innerIdx);
122 Scalar mug = underSaturatedTable.get(
"MUG" , innerIdx);
124 invGasB.appendSamplePoint(outerIdx, Rv, 1.0/Bg);
125 gasMu.appendSamplePoint(outerIdx, Rv, mug);
130 std::vector<double> tmpPressure = saturatedTable.getColumn(
"PG").vectorCopy( );
132 invSatGasB.setXYContainers(tmpPressure, invSatGasBArray);
133 invSatGasBMu.setXYContainers(tmpPressure, invSatGasBMuArray);
137 for (
unsigned xIdx = 0; xIdx < invGasB.numX(); ++xIdx) {
139 assert(invGasB.numY(xIdx) > 0);
143 if (invGasB.numY(xIdx) > 1)
149 size_t masterTableIdx = xIdx + 1;
150 for (; masterTableIdx < saturatedTable.numRows(); ++masterTableIdx)
152 if (pvtgTable.getUnderSaturatedTable(masterTableIdx).numRows() > 1)
156 if (masterTableIdx >= saturatedTable.numRows())
157 OPM_THROW(std::runtime_error,
158 "PVTG tables are invalid: The last table must exhibit at least one " 159 "entry for undersaturated gas!");
163 extendPvtgTable_(regionIdx,
165 pvtgTable.getUnderSaturatedTable(xIdx),
166 pvtgTable.getUnderSaturatedTable(masterTableIdx));
171 if (deck.hasKeyword(
"VAPPARS")) {
172 const auto& vapParsKeyword = deck.getKeyword(
"VAPPARS");
173 vapPar1_ = vapParsKeyword.getRecord(0).getItem(
"OIL_VAP_PROPENSITY").template get<double>(0);
180 void extendPvtgTable_(
unsigned regionIdx,
182 const SimpleTable& curTable,
183 const SimpleTable& masterTable)
185 std::vector<double> RvArray = curTable.getColumn(
"RV").vectorCopy();
186 std::vector<double> gasBArray = curTable.getColumn(
"BG").vectorCopy();
187 std::vector<double> gasMuArray = curTable.getColumn(
"MUG").vectorCopy();
189 auto& invGasB = inverseGasB_[regionIdx];
190 auto& gasMu = gasMu_[regionIdx];
192 for (
size_t newRowIdx = 1; newRowIdx < masterTable.numRows(); ++ newRowIdx) {
193 const auto& RVColumn = masterTable.getColumn(
"RV");
194 const auto& BGColumn = masterTable.getColumn(
"BG");
195 const auto& viscosityColumn = masterTable.getColumn(
"MUG");
198 Scalar diffRv = RVColumn[newRowIdx] - RVColumn[newRowIdx - 1];
199 Scalar newRv = RvArray.back() + diffRv;
202 Scalar B1 = BGColumn[newRowIdx];
203 Scalar B2 = BGColumn[newRowIdx - 1];
204 Scalar x = (B1 - B2)/( (B1 + B2)/2.0 );
208 Scalar newBg = gasBArray.back()*(1.0 + x/2.0)/(1.0 - x/2.0);
211 Scalar mu1 = viscosityColumn[newRowIdx];
212 Scalar mu2 = viscosityColumn[newRowIdx - 1];
213 Scalar xMu = (mu1 - mu2)/( (mu1 + mu2)/2.0 );
217 Scalar newMug = gasMuArray.back()*(1.0 + xMu/2)/(1.0 - xMu/2.0);
221 RvArray.push_back(newRv);
222 gasBArray.push_back(newBg);
223 gasMuArray.push_back(newMug);
226 invGasB.appendSamplePoint(xIdx, newRv, 1.0/newBg);
227 gasMu.appendSamplePoint(xIdx, newRv, newMug);
232 #endif // HAVE_OPM_PARSER 243 saturatedOilVaporizationFactorTable_.resize(
numRegions);
255 oilReferenceDensity_[regionIdx] = rhoRefOil;
256 gasReferenceDensity_[regionIdx] = rhoRefGas;
265 { saturatedOilVaporizationFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
278 auto& invGasB = inverseGasB_[regionIdx];
280 const auto& RvTable = saturatedOilVaporizationFactorTable_[regionIdx];
282 Scalar T = 273.15 + 15.56;
285 Scalar RvMax = RvTable.eval(saturatedOilVaporizationFactorTable_[regionIdx].xMax(),
true);
287 Scalar poMin = samplePoints.front().first;
288 Scalar poMax = samplePoints.back().first;
291 size_t nP = samplePoints.size()*2;
293 Scalar rhogRef = gasReferenceDensity_[regionIdx];
294 Scalar rhooRef = oilReferenceDensity_[regionIdx];
299 updateSaturationPressure_(regionIdx);
305 for (
size_t RvIdx = 0; RvIdx < nRv; ++RvIdx) {
306 Scalar Rv = RvMin + (RvMax - RvMin)*RvIdx/nRv;
308 invGasB.appendXPos(Rv);
310 for (
size_t pIdx = 0; pIdx < nP; ++pIdx) {
311 Scalar pg = poMin + (poMax - poMin)*pIdx/nP;
314 Scalar BgSat = gasFormationVolumeFactor.
eval(poSat,
true);
315 Scalar drhoo_dp = (1.1200 - 1.1189)/((5000 - 4000)*6894.76);
316 Scalar rhoo = rhooRef/BgSat*(1 + drhoo_dp*(pg - poSat));
318 Scalar Bg = rhooRef/rhoo;
320 invGasB.appendSamplePoint(RvIdx, pg, 1.0/Bg);
338 { inverseGasB_[regionIdx] = invBg; }
346 { gasMu_[regionIdx] = mug; }
357 auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
360 Scalar RvMax = oilVaporizationFac.eval(saturatedOilVaporizationFactorTable_[regionIdx].xMax(),
true);
362 Scalar poMin = samplePoints.front().first;
363 Scalar poMax = samplePoints.back().first;
366 size_t nP = samplePoints.size()*2;
373 for (
size_t RvIdx = 0; RvIdx < nRv; ++RvIdx) {
374 Scalar Rv = RvMin + (RvMax - RvMin)*RvIdx/nRv;
376 gasMu_[regionIdx].appendXPos(Rv);
378 for (
size_t pIdx = 0; pIdx < nP; ++pIdx) {
379 Scalar pg = poMin + (poMax - poMin)*pIdx/nP;
380 Scalar mug = mugTable.
eval(pg,
true);
382 gasMu_[regionIdx].appendSamplePoint(RvIdx, pg, mug);
394 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
397 const auto& gasMu = gasMu_[regionIdx];
398 const auto& invGasB = inverseGasB_[regionIdx];
399 assert(gasMu.numX() == invGasB.numX());
401 auto& invGasBMu = inverseGasBMu_[regionIdx];
402 auto& invSatGasB = inverseSaturatedGasB_[regionIdx];
403 auto& invSatGasBMu = inverseSaturatedGasBMu_[regionIdx];
405 std::vector<Scalar> satPressuresArray;
406 std::vector<Scalar> invSatGasBArray;
407 std::vector<Scalar> invSatGasBMuArray;
408 for (
size_t pIdx = 0; pIdx < gasMu.numX(); ++pIdx) {
409 invGasBMu.appendXPos(gasMu.xAt(pIdx));
411 assert(gasMu.numY(pIdx) == invGasB.numY(pIdx));
413 size_t numRv = gasMu.numY(pIdx);
414 for (
size_t rvIdx = 0; rvIdx < numRv; ++rvIdx)
415 invGasBMu.appendSamplePoint(pIdx,
416 gasMu.yAt(pIdx, rvIdx),
417 invGasB.valueAt(pIdx, rvIdx)
418 / gasMu.valueAt(pIdx, rvIdx));
423 satPressuresArray.push_back(gasMu.xAt(pIdx));
424 invSatGasBArray.push_back(invGasB.valueAt(pIdx, numRv - 1));
425 invSatGasBMuArray.push_back(invGasBMu.valueAt(pIdx, numRv - 1));
428 invSatGasB.setXYContainers(satPressuresArray, invSatGasBArray);
429 invSatGasBMu.setXYContainers(satPressuresArray, invSatGasBMuArray);
431 updateSaturationPressure_(regionIdx);
439 {
return gasReferenceDensity_.size(); }
444 template <
class Evaluation>
447 const Evaluation& pressure,
448 const Evaluation& Rv)
const 450 const Evaluation& invBg = inverseGasB_[regionIdx].eval(pressure, Rv,
true);
451 const Evaluation& invMugBg = inverseGasBMu_[regionIdx].eval(pressure, Rv,
true);
453 return invBg/invMugBg;
459 template <
class Evaluation>
462 const Evaluation& pressure)
const 464 const Evaluation& invBg = inverseSaturatedGasB_[regionIdx].eval(pressure,
true);
465 const Evaluation& invMugBg = inverseSaturatedGasBMu_[regionIdx].eval(pressure,
true);
467 return invBg/invMugBg;
473 template <
class Evaluation>
476 const Evaluation& pressure,
477 const Evaluation& Rv)
const 478 {
return inverseGasB_[regionIdx].eval(pressure, Rv,
true); }
483 template <
class Evaluation>
486 const Evaluation& pressure)
const 487 {
return inverseSaturatedGasB_[regionIdx].eval(pressure,
true); }
492 template <
class Evaluation>
495 const Evaluation& pressure)
const 497 return saturatedOilVaporizationFactorTable_[regionIdx].eval(pressure,
true);
507 template <
class Evaluation>
510 const Evaluation& pressure,
511 const Evaluation& oilSaturation,
512 Scalar maxOilSaturation)
const 515 saturatedOilVaporizationFactorTable_[regionIdx].eval(pressure,
true);
519 maxOilSaturation = std::min(maxOilSaturation, Scalar(1.0));
520 if (vapPar1_ > 0.0 && maxOilSaturation > 0.01 && oilSaturation < maxOilSaturation) {
521 static const Scalar eps = 0.001;
522 const Evaluation& So = Opm::max(oilSaturation, eps);
523 tmp *= Opm::max(1e-3, Opm::pow(So/maxOilSaturation, vapPar1_));
539 template <
class Evaluation>
541 const Evaluation& temperature OPM_UNUSED,
542 const Evaluation& Rv)
const 546 const auto& RvTable = saturatedOilVaporizationFactorTable_[regionIdx];
547 const Scalar eps = std::numeric_limits<typename Toolbox::Scalar>::epsilon()*1e6;
550 Evaluation pSat = saturationPressure_[regionIdx].eval(Rv,
true);
555 bool onProbation =
false;
556 for (
unsigned i = 0; i < 20; ++i) {
557 const Evaluation& f = RvTable.eval(pSat,
true) - Rv;
558 const Evaluation& fPrime = RvTable.evalDerivative(pSat,
true);
562 if (std::abs(Opm::scalarValue(fPrime)) < 1.0e-30) {
566 const Evaluation& delta = f/fPrime;
580 if (std::abs(Opm::scalarValue(delta)) < std::abs(Opm::scalarValue(pSat))*eps)
584 std::stringstream errlog;
585 errlog <<
"Finding saturation pressure did not converge:" 586 <<
" pSat = " << pSat
588 OpmLog::debug(
"Wet gas saturation pressure", errlog.str());
589 OPM_THROW_NOLOG(NumericalProblem, errlog.str());
593 void updateSaturationPressure_(
unsigned regionIdx)
595 typedef std::pair<Scalar, Scalar> Pair;
596 const auto& oilVaporizationFac = saturatedOilVaporizationFactorTable_[regionIdx];
600 size_t n = oilVaporizationFac.numSamples();
601 Scalar delta = (oilVaporizationFac.xMax() - oilVaporizationFac.xMin())/Scalar(n + 1);
603 SamplingPoints pSatSamplePoints;
605 for (
size_t i = 0; i <= n; ++ i) {
606 Scalar pSat = oilVaporizationFac.xMin() + Scalar(i)*delta;
610 pSatSamplePoints.push_back(val);
614 auto x_coord_comparator = [](
const Pair& a,
const Pair& b) {
return a.first == b.first; };
615 auto last = std::unique(pSatSamplePoints.begin(), pSatSamplePoints.end(), x_coord_comparator);
616 pSatSamplePoints.erase(last, pSatSamplePoints.end());
618 saturationPressure_[regionIdx].setContainerOfTuples(pSatSamplePoints);
621 std::vector<Scalar> gasReferenceDensity_;
622 std::vector<Scalar> oilReferenceDensity_;
623 std::vector<TabulatedTwoDFunction> inverseGasB_;
624 std::vector<TabulatedOneDFunction> inverseSaturatedGasB_;
625 std::vector<TabulatedTwoDFunction> gasMu_;
626 std::vector<TabulatedTwoDFunction> inverseGasBMu_;
627 std::vector<TabulatedOneDFunction> inverseSaturatedGasBMu_;
628 std::vector<TabulatedOneDFunction> saturatedOilVaporizationFactorTable_;
629 std::vector<TabulatedOneDFunction> saturationPressure_;
void setSaturatedGasViscosity(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the phase viscosity for oil saturated gas.
Definition: WetGasPvt.hpp:355
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:484
void setInverseGasFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction &invBg)
Initialize the function for the gas formation volume factor.
Definition: WetGasPvt.hpp:337
Definition: Air_Mesitylene.hpp:33
Evaluation eval(const Evaluation &x, bool extrapolate=false) const
Evaluate the spline at a given position.
Definition: Tabulated1DFunction.hpp:252
void setGasViscosity(unsigned regionIdx, const TabulatedTwoDFunction &mug)
Initialize the viscosity of the gas phase.
Definition: WetGasPvt.hpp:345
void setSaturatedGasOilVaporizationFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the oil vaporization factor .
Definition: WetGasPvt.hpp:264
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar rhoRefGas, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: WetGasPvt.hpp:250
This class represents the Pressure-Volume-Temperature relations of the gas phas with vaporized oil...
Definition: WetGasPvt.hpp:47
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &temperature OPM_UNUSED, 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:540
Implements a linearly interpolated scalar function that depends on one variable.
A central place for various physical constants occuring in some equations.
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
void initEnd()
Finish initializing the gas phase PVT properties.
Definition: WetGasPvt.hpp:390
void setSaturatedGasFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the gas formation volume factor.
Definition: WetGasPvt.hpp:276
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:474
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: WetGasPvt.hpp:438
This file provides a wrapper around the "final" C++-2011 statement.
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:445
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47
Evaluation saturatedOilVaporizationFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &oilSaturation, Scalar maxOilSaturation) const
Returns the oil vaporization factor [m^3/m^3] of the gas phase.
Definition: WetGasPvt.hpp:508
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:493
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:460