25 #ifndef OPM_WET_GAS_PVT_HPP
26 #define OPM_WET_GAS_PVT_HPP
36 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
41 template <
class Scalar>
42 class OilPvtMultiplexer;
48 template <
class Scalar>
56 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
65 void initFromDeck(DeckConstPtr deck, EclipseStateConstPtr eclState)
67 const auto& pvtgTables = eclState->getTableManager()->getPvtgTables();
68 DeckKeywordConstPtr densityKeyword = deck->getKeyword(
"DENSITY");
70 assert(pvtgTables.size() == densityKeyword->size());
72 size_t numRegions = pvtgTables.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);
84 Scalar T = 273.15 + 15.56;
92 const auto& pvtgTable = pvtgTables[regionIdx];
94 const auto saturatedTable = pvtgTable.getOuterTable();
95 assert(saturatedTable->numRows() > 1);
97 auto& gasMu = gasMu_[regionIdx];
98 auto& invGasB = inverseGasB_[regionIdx];
99 auto& oilVaporizationFac = oilVaporizationFactorTable_[regionIdx];
101 oilVaporizationFac.setXYArrays(saturatedTable->numRows(),
102 saturatedTable->getPressureColumn(),
103 saturatedTable->getOilSolubilityColumn());
106 for (
size_t outerIdx = 0; outerIdx < saturatedTable->numRows(); ++ outerIdx) {
107 Scalar pg = saturatedTable->getPressureColumn()[outerIdx];
109 invGasB.appendXPos(pg);
110 gasMu.appendXPos(pg);
112 assert(invGasB.numX() == outerIdx + 1);
113 assert(gasMu.numX() == outerIdx + 1);
115 const auto underSaturatedTable = pvtgTable.getInnerTable(outerIdx);
116 size_t numRows = underSaturatedTable->numRows();
117 for (
size_t innerIdx = 0; innerIdx < numRows; ++ innerIdx) {
118 Scalar Rv = underSaturatedTable->getOilSolubilityColumn()[innerIdx];
119 Scalar Bg = underSaturatedTable->getGasFormationFactorColumn()[innerIdx];
120 Scalar mug = underSaturatedTable->getGasViscosityColumn()[innerIdx];
122 invGasB.appendSamplePoint(outerIdx, Rv, 1.0/Bg);
123 gasMu.appendSamplePoint(outerIdx, Rv, mug);
128 for (
size_t xIdx = 0; xIdx < invGasB.numX(); ++xIdx) {
130 assert(invGasB.numY(xIdx) > 0);
134 if (invGasB.numY(xIdx) > 1)
140 size_t masterTableIdx = xIdx + 1;
141 for (; masterTableIdx < pvtgTable.getOuterTable()->numRows(); ++masterTableIdx)
143 if (pvtgTable.getInnerTable(masterTableIdx)->numRows() > 1)
147 if (masterTableIdx >= pvtgTable.getOuterTable()->numRows())
148 OPM_THROW(std::runtime_error,
149 "PVTG tables are invalid: The last table must exhibit at least one "
150 "entry for undersaturated gas!");
156 const auto masterTable = pvtgTable.getInnerTable(masterTableIdx);
157 const auto curTable = pvtgTable.getInnerTable(xIdx);
158 for (
size_t newRowIdx = 1; newRowIdx < masterTable->numRows(); ++ newRowIdx) {
160 masterTable->getOilSolubilityColumn()[newRowIdx]
161 / masterTable->getOilSolubilityColumn()[0];
164 masterTable->getGasFormationFactorColumn()[newRowIdx]
165 / masterTable->getGasFormationFactorColumn()[0];
168 masterTable->getGasViscosityColumn()[newRowIdx]
169 / masterTable->getGasViscosityColumn()[0];
171 Scalar newRv = curTable->getOilSolubilityColumn()[0]*alphaRv;
172 Scalar newBg = curTable->getGasFormationFactorColumn()[0]*alphaBg;
173 Scalar newMug = curTable->getGasViscosityColumn()[0]*alphaMug;
175 invGasB.appendSamplePoint(xIdx, newRv, 1.0/newBg);
176 gasMu.appendSamplePoint(xIdx, newRv, newMug);
181 #endif // HAVE_OPM_PARSER
185 oilMolarMass_.resize(numRegions);
186 gasMolarMass_.resize(numRegions);
187 oilReferenceDensity_.resize(numRegions);
188 gasReferenceDensity_.resize(numRegions);
189 inverseGasB_.resize(numRegions);
190 inverseGasBMu_.resize(numRegions);
191 gasMu_.resize(numRegions);
192 oilVaporizationFactorTable_.resize(numRegions);
193 saturationPressureSpline_.resize(numRegions);
204 oilReferenceDensity_[regionIdx] = rhoRefOil;
205 gasReferenceDensity_[regionIdx] = rhoRefGas;
216 oilMolarMass_[regionIdx] = MOil;
217 gasMolarMass_[regionIdx] = MGas;
226 { oilVaporizationFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
239 auto& invGasB = inverseGasB_[regionIdx];
241 auto &RvTable = oilVaporizationFactorTable_[regionIdx];
243 Scalar T = 273.15 + 15.56;
246 Scalar RvMax = RvTable.eval(oilVaporizationFactorTable_[regionIdx].xMax(),
true);
248 Scalar poMin = samplePoints.front().first;
249 Scalar poMax = samplePoints.back().first;
252 size_t nP = samplePoints.size()*2;
254 Scalar rhogRef = gasReferenceDensity_[regionIdx];
255 Scalar rhooRef = oilReferenceDensity_[regionIdx];
257 Spline gasFormationVolumeFactorSpline;
260 updateSaturationPressureSpline_(regionIdx);
266 for (
size_t RvIdx = 0; RvIdx < nRv; ++RvIdx) {
267 Scalar Rv = RvMin + (RvMax - RvMin)*RvIdx/nRv;
268 Scalar XgO = Rv/(rhooRef/rhogRef + Rv);
270 invGasB.appendXPos(Rv);
272 for (
size_t pIdx = 0; pIdx < nP; ++pIdx) {
273 Scalar pg = poMin + (poMax - poMin)*pIdx/nP;
276 Scalar BgSat = gasFormationVolumeFactorSpline.
eval(poSat,
true);
277 Scalar drhoo_dp = (1.1200 - 1.1189)/((5000 - 4000)*6894.76);
278 Scalar rhoo = rhooRef/BgSat*(1 + drhoo_dp*(pg - poSat));
280 Scalar Bg = rhooRef/rhoo;
282 invGasB.appendSamplePoint(RvIdx, pg, 1.0/Bg);
300 { inverseGasB_[regionIdx] = invBg; }
308 { gasMu_[regionIdx] = mug; }
319 auto& oilVaporizationFac = oilVaporizationFactorTable_[regionIdx];
322 Scalar RvMax = oilVaporizationFac.eval(oilVaporizationFactorTable_[regionIdx].xMax(),
true);
324 Scalar poMin = samplePoints.front().first;
325 Scalar poMax = samplePoints.back().first;
328 size_t nP = samplePoints.size()*2;
335 for (
size_t RvIdx = 0; RvIdx < nRv; ++RvIdx) {
336 Scalar Rv = RvMin + (RvMax - RvMin)*RvIdx/nRv;
338 gasMu_[regionIdx].appendXPos(Rv);
340 for (
size_t pIdx = 0; pIdx < nP; ++pIdx) {
341 Scalar pg = poMin + (poMax - poMin)*pIdx/nP;
342 Scalar mug = mugSpline.
eval(pg,
true);
344 gasMu_[regionIdx].appendSamplePoint(RvIdx, pg, mug);
357 size_t numRegions = gasMu_.size();
358 for (
unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
361 const auto& gasMu = gasMu_[regionIdx];
362 const auto& invGasB = inverseGasB_[regionIdx];
363 assert(gasMu.numX() == invGasB.numX());
365 auto& invGasBMu = inverseGasBMu_[regionIdx];
367 for (
size_t pIdx = 0; pIdx < gasMu.numX(); ++pIdx) {
368 invGasBMu.appendXPos(gasMu.xAt(pIdx));
370 assert(gasMu.numY(pIdx) == invGasB.numY(pIdx));
372 size_t numPressures = gasMu.numY(pIdx);
373 for (
size_t rvIdx = 0; rvIdx < numPressures; ++rvIdx)
374 invGasBMu.appendSamplePoint(pIdx,
375 gasMu.yAt(pIdx, rvIdx),
376 invGasB.valueAt(pIdx, rvIdx)*
377 1/gasMu.valueAt(pIdx, rvIdx));
380 updateSaturationPressureSpline_(regionIdx);
387 template <
class Evaluation>
390 const Evaluation& pressure,
391 const Evaluation& XgO)
const
393 Scalar rhooRef = oilReferenceDensity_[regionIdx];
394 Scalar rhogRef = gasReferenceDensity_[regionIdx];
395 const Evaluation& Rv = XgO/(1 - XgO)*(rhogRef/rhooRef);
397 const Evaluation& invBg = inverseGasB_[regionIdx].eval(pressure, Rv,
true);
398 const Evaluation& invMugBg = inverseGasBMu_[regionIdx].eval(pressure, Rv,
true);
400 return invBg/invMugBg;
406 template <
class Evaluation>
408 const Evaluation& temperature,
409 const Evaluation& pressure,
410 const Evaluation& XgO)
const
412 Scalar rhooRef = oilReferenceDensity_[regionIdx];
413 Scalar rhogRef = gasReferenceDensity_[regionIdx];
417 Evaluation rhog = rhogRef/Bg;
422 const Evaluation& Rv = XgO/(1 - XgO)*(rhogRef/rhooRef);
423 rhog += (rhogRef*Rv)/Bg;
431 template <
class Evaluation>
434 const Evaluation& pressure,
435 const Evaluation& XgO)
const
437 Scalar rhooRef = oilReferenceDensity_[regionIdx];
438 Scalar rhogRef = gasReferenceDensity_[regionIdx];
440 const Evaluation& Rv = XgO/(1-XgO)*(rhogRef/rhooRef);
442 return 1.0 / inverseGasB_[regionIdx].eval(pressure, Rv,
true);
449 template <
class Evaluation>
452 const Evaluation& )
const
459 template <
class Evaluation>
461 const Evaluation& temperature,
462 const Evaluation& pressure)
const
474 return phi_oO / x_gOSat;
477 template <
class Evaluation>
479 const Evaluation& temperature,
480 const Evaluation& pressure)
const
490 template <
class Evaluation>
493 const Evaluation& pressure)
const
494 {
return oilVaporizationFactorTable_[regionIdx].eval(pressure,
true); }
502 template <
class Evaluation>
504 const Evaluation& temperature,
505 const Evaluation& XgO)
const
510 Evaluation pSat = saturationPressureSpline_[regionIdx].eval(XgO,
true);
511 const Evaluation& eps = pSat*1e-11;
516 for (
unsigned i = 0; i < 20; ++i) {
520 const Evaluation& delta = f/fPrime;
523 if (
std::abs(Toolbox::value(delta)) <
std::abs(Toolbox::value(pSat)) * 1e-10)
527 OPM_THROW(NumericalProblem,
"Could find the gas saturation pressure for X_g^O = " << XgO);
530 template <
class Evaluation>
532 const Evaluation& temperature,
533 const Evaluation& pressure)
const
535 Scalar rho_gRef = gasReferenceDensity_[regionIdx];
536 Scalar rho_oRef = oilReferenceDensity_[regionIdx];
542 const Evaluation& rho_gO = Rv * rho_oRef;
546 return rho_gO/(rho_gRef + rho_gO);
549 template <
class Evaluation>
551 const Evaluation& temperature,
552 const Evaluation& pressure)
const
559 Scalar MG = gasMolarMass_[regionIdx];
560 Scalar MO = oilMolarMass_[regionIdx];
562 const Evaluation& avgMolarMass = MO/(1 + (1 - XgO)*(MO/MG - 1));
563 return XgO*avgMolarMass/MO;
567 void updateSaturationPressureSpline_(
unsigned regionIdx)
569 auto& oilVaporizationFac = oilVaporizationFactorTable_[regionIdx];
573 size_t n = oilVaporizationFac.numSamples()*5;
574 Scalar delta = (oilVaporizationFac.xMax() - oilVaporizationFac.xMin())/(n + 1);
576 SamplingPoints pSatSamplePoints;
578 for (
size_t i = 0; i <= n; ++ i) {
579 Scalar pSat = oilVaporizationFac.xMin() + i*delta;
582 std::pair<Scalar, Scalar> val(XgO, pSat);
583 pSatSamplePoints.push_back(val);
585 saturationPressureSpline_[regionIdx].setContainerOfTuples(pSatSamplePoints,
589 const OilPvtMultiplexer *oilPvt_;
591 std::vector<Scalar> gasMolarMass_;
592 std::vector<Scalar> oilMolarMass_;
593 std::vector<Scalar> gasReferenceDensity_;
594 std::vector<Scalar> oilReferenceDensity_;
595 std::vector<TabulatedTwoDFunction> inverseGasB_;
596 std::vector<TabulatedTwoDFunction> gasMu_;
597 std::vector<TabulatedTwoDFunction> inverseGasBMu_;
598 std::vector<TabulatedOneDFunction> oilVaporizationFactorTable_;
599 std::vector<Spline> saturationPressureSpline_;
Evaluation saturatedGasOilMassFraction(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Definition: WetGasPvt.hpp:531
Evaluation oilVaporizationFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: WetGasPvt.hpp:491
Evaluation formationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &XgO) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WetGasPvt.hpp:432
Evaluation fugacityCoefficientGas(unsigned, const Evaluation &, const Evaluation &) const
Returns the fugacity coefficient [Pa] of a component in the fluid phase given a set of parameters...
Definition: WetGasPvt.hpp:450
Evaluation fugacityCoefficientOil(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the fugacity coefficient [-] of the oil component in the oil phase given a pressure and a tem...
Definition: OilPvtMultiplexer.hpp:166
void initEnd(const OilPvtMultiplexer *oilPvt)
Finish initializing the gas phase PVT properties.
Definition: WetGasPvt.hpp:352
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
Definition: DryGasPvt.hpp:44
Definition: Spline.hpp:104
Definition: Air_Mesitylene.hpp:31
void setSaturatedGasViscosity(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the phase viscosity for oil saturated gas.
Definition: WetGasPvt.hpp:317
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &XgO) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WetGasPvt.hpp:388
Evaluation fugacityCoefficientOil(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Definition: WetGasPvt.hpp:460
void setContainerOfTuples(const XYContainer &points, Scalar m0, Scalar m1, bool sortInputs=false)
Set the sampling points and the boundary slopes of a full spline using a STL-compatible container of ...
Definition: Spline.hpp:474
Evaluation gasSaturationPressure(unsigned regionIdx, const Evaluation &temperature, const Evaluation &XgO) const
Returns the saturation pressure of the gas phase [Pa] depending on its mass fraction of the oil compo...
Definition: WetGasPvt.hpp:503
Evaluation fugacityCoefficientWater(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Definition: WetGasPvt.hpp:478
Class implementing cubic splines.
void setSaturatedGasOilVaporizationFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the oil vaporization factor .
Definition: WetGasPvt.hpp:225
Implements a linearly interpolated scalar function that depends on one variable.
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar rhoRefGas, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: WetGasPvt.hpp:199
A central place for various physical constants occuring in some equations.
This class represents the Pressure-Volume-Temperature relations of the gas phas with vaporized oil...
Definition: WetGasPvt.hpp:49
This file provides a wrapper around the "final" C++-2011 statement.
void setInverseGasFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction &invBg)
Initialize the function for the gas formation volume factor.
Definition: WetGasPvt.hpp:299
Evaluation< Scalar, VarSetTag, numVars > abs(const Evaluation< Scalar, VarSetTag, numVars > &)
Definition: Math.hpp:41
Class implementing cubic splines.
Definition: Spline.hpp:89
void setGasViscosity(unsigned regionIdx, const TabulatedTwoDFunction &mug)
Initialize the viscosity of the gas phase.
Definition: WetGasPvt.hpp:307
Evaluation saturatedGasOilMoleFraction(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Definition: WetGasPvt.hpp:550
void setNumRegions(size_t numRegions)
Definition: WetGasPvt.hpp:183
void setMolarMasses(unsigned regionIdx, Scalar MOil, Scalar MGas, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition: WetGasPvt.hpp:211
Scalar eval(Scalar x, bool extrapolate=false) const
Evaluate the spline at a given position.
Definition: Spline.hpp:809
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:44
A central place for various physical constants occuring in some equations.
Definition: Constants.hpp:39
Evaluation density(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &XgO) const
Returns the density [kg/m^3] of the fluid phase given a set of parameters.
Definition: WetGasPvt.hpp:407
void setSaturatedGasFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the gas formation volume factor.
Definition: WetGasPvt.hpp:237