00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00027 #ifndef OPM_DRY_GAS_PVT_HPP
00028 #define OPM_DRY_GAS_PVT_HPP
00029
00030 #include <opm/material/Constants.hpp>
00031
00032 #include <opm/material/common/Tabulated1DFunction.hpp>
00033
00034 #if HAVE_OPM_PARSER
00035 #include <opm/parser/eclipse/Deck/Deck.hpp>
00036 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
00037 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
00038 #include <opm/parser/eclipse/EclipseState/Tables/PvdgTable.hpp>
00039 #include <opm/parser/eclipse/Deck/DeckKeyword.hpp>
00040 #include <opm/parser/eclipse/Deck/DeckRecord.hpp>
00041 #endif
00042
00043 #include <vector>
00044
00045 namespace Opm {
00050 template <class Scalar>
00051 class DryGasPvt
00052 {
00053 typedef Opm::Tabulated1DFunction<Scalar> TabulatedOneDFunction;
00054 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
00055
00056 public:
00057 #if HAVE_OPM_PARSER
00058
00063 void initFromDeck(const Deck& deck, const EclipseState& eclState)
00064 {
00065 const auto& pvdgTables = eclState.getTableManager().getPvdgTables();
00066 const auto& densityKeyword = deck.getKeyword("DENSITY");
00067
00068 assert(pvdgTables.size() == densityKeyword.size());
00069
00070 size_t numRegions = pvdgTables.size();
00071 setNumRegions(numRegions);
00072
00073 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
00074 Scalar rhoRefO = densityKeyword.getRecord(regionIdx).getItem("OIL").getSIDouble(0);
00075 Scalar rhoRefG = densityKeyword.getRecord(regionIdx).getItem("GAS").getSIDouble(0);
00076 Scalar rhoRefW = densityKeyword.getRecord(regionIdx).getItem("WATER").getSIDouble(0);
00077
00078 setReferenceDensities(regionIdx, rhoRefO, rhoRefG, rhoRefW);
00079
00080
00081 Scalar p = 1.01325e5;
00082 Scalar T = 273.15 + 15.56;
00083 Scalar MO = 175e-3;
00084 Scalar MG = Opm::Constants<Scalar>::R*T*rhoRefG / p;
00085 Scalar MW = 18.0e-3;
00086
00087
00088 setMolarMasses(regionIdx, MO, MG, MW);
00089
00090 const auto& pvdgTable = pvdgTables.getTable<PvdgTable>(regionIdx);
00091
00092
00093
00094 std::vector<Scalar> invB(pvdgTable.numRows());
00095 const auto& Bg = pvdgTable.getFormationFactorColumn();
00096 for (unsigned i = 0; i < Bg.size(); ++ i) {
00097 invB[i] = 1.0/Bg[i];
00098 }
00099
00100 size_t numSamples = invB.size();
00101 inverseGasB_[regionIdx].setXYArrays(numSamples, pvdgTable.getPressureColumn(), invB);
00102 gasMu_[regionIdx].setXYArrays(numSamples, pvdgTable.getPressureColumn(), pvdgTable.getViscosityColumn());
00103 }
00104
00105 initEnd();
00106 }
00107 #endif
00108
00109 void setNumRegions(size_t numRegions)
00110 {
00111 gasReferenceDensity_.resize(numRegions);
00112 inverseGasB_.resize(numRegions);
00113 inverseGasBMu_.resize(numRegions);
00114 gasMu_.resize(numRegions);
00115 }
00116
00117
00121 void setReferenceDensities(unsigned regionIdx,
00122 Scalar ,
00123 Scalar rhoRefGas,
00124 Scalar )
00125 {
00126 gasReferenceDensity_[regionIdx] = rhoRefGas;
00127 }
00128
00132 void setMolarMasses(unsigned ,
00133 Scalar ,
00134 Scalar ,
00135 Scalar )
00136 { }
00137
00143 void setGasViscosity(unsigned regionIdx, const TabulatedOneDFunction& mug)
00144 { gasMu_[regionIdx] = mug; }
00145
00151 void setGasFormationVolumeFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
00152 {
00153 SamplingPoints tmp(samplePoints);
00154 auto it = tmp.begin();
00155 const auto& endIt = tmp.end();
00156 for (; it != endIt; ++ it)
00157 std::get<1>(*it) = 1.0/std::get<1>(*it);
00158
00159 inverseGasB_[regionIdx].setContainerOfTuples(tmp);
00160 assert(inverseGasB_[regionIdx].monotonic());
00161 }
00162
00166 void initEnd()
00167 {
00168
00169 size_t numRegions = gasMu_.size();
00170 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
00171
00172
00173 const auto& gasMu = gasMu_[regionIdx];
00174 const auto& invGasB = inverseGasB_[regionIdx];
00175 assert(gasMu.numSamples() == invGasB.numSamples());
00176
00177 std::vector<Scalar> pressureValues(gasMu.numSamples());
00178 std::vector<Scalar> invGasBMuValues(gasMu.numSamples());
00179 for (unsigned pIdx = 0; pIdx < gasMu.numSamples(); ++pIdx) {
00180 pressureValues[pIdx] = invGasB.xAt(pIdx);
00181 invGasBMuValues[pIdx] = invGasB.valueAt(pIdx) * (1.0/gasMu.valueAt(pIdx));
00182 }
00183
00184 inverseGasBMu_[regionIdx].setXYContainers(pressureValues, invGasBMuValues);
00185 }
00186 }
00187
00191 unsigned numRegions() const
00192 { return gasReferenceDensity_.size(); }
00193
00197 template <class Evaluation>
00198 Evaluation viscosity(unsigned regionIdx,
00199 const Evaluation& temperature,
00200 const Evaluation& pressure,
00201 const Evaluation& ) const
00202 { return saturatedViscosity(regionIdx, temperature, pressure); }
00203
00207 template <class Evaluation>
00208 Evaluation saturatedViscosity(unsigned regionIdx,
00209 const Evaluation& ,
00210 const Evaluation& pressure) const
00211 {
00212 const Evaluation& invBg = inverseGasB_[regionIdx].eval(pressure, true);
00213 const Evaluation& invMugBg = inverseGasBMu_[regionIdx].eval(pressure, true);
00214
00215 return invBg/invMugBg;
00216 }
00217
00221 template <class Evaluation>
00222 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
00223 const Evaluation& temperature,
00224 const Evaluation& pressure,
00225 const Evaluation& ) const
00226 { return saturatedInverseFormationVolumeFactor(regionIdx, temperature, pressure); }
00227
00231 template <class Evaluation>
00232 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
00233 const Evaluation& ,
00234 const Evaluation& pressure) const
00235 { return inverseGasB_[regionIdx].eval(pressure, true); }
00236
00243 template <class Evaluation>
00244 Evaluation saturationPressure(unsigned ,
00245 const Evaluation& ,
00246 const Evaluation& ) const
00247 { return 0.0; }
00248
00252 template <class Evaluation>
00253 Evaluation saturatedOilVaporizationFactor(unsigned ,
00254 const Evaluation& ,
00255 const Evaluation& ,
00256 const Evaluation& ,
00257 Scalar ) const
00258 { return 0.0; }
00259
00263 template <class Evaluation>
00264 Evaluation saturatedOilVaporizationFactor(unsigned ,
00265 const Evaluation& ,
00266 const Evaluation& ) const
00267 { return 0.0; }
00268
00269 private:
00270 std::vector<Scalar> gasReferenceDensity_;
00271 std::vector<TabulatedOneDFunction> inverseGasB_;
00272 std::vector<TabulatedOneDFunction> gasMu_;
00273 std::vector<TabulatedOneDFunction> inverseGasBMu_;
00274 };
00275
00276 }
00277
00278 #endif