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_DEAD_OIL_PVT_HPP
00028 #define OPM_DEAD_OIL_PVT_HPP
00029
00030 #include <opm/material/common/UniformXTabulated2DFunction.hpp>
00031 #include <opm/material/common/Tabulated1DFunction.hpp>
00032 #include <opm/material/common/Spline.hpp>
00033
00034 #if HAVE_OPM_PARSER
00035 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
00036 #include <opm/parser/eclipse/EclipseState/Tables/PvdoTable.hpp>
00037 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
00038 #endif
00039
00040 namespace Opm {
00045 template <class Scalar>
00046 class DeadOilPvt
00047 {
00048 typedef Opm::Tabulated1DFunction<Scalar> TabulatedOneDFunction;
00049 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
00050
00051 public:
00052 #if HAVE_OPM_PARSER
00053
00056 void initFromDeck(const Deck& deck, const EclipseState& eclState)
00057 {
00058 const auto& pvdoTables = eclState.getTableManager().getPvdoTables();
00059 const auto& densityKeyword = deck.getKeyword("DENSITY");
00060
00061 assert(pvdoTables.size() == densityKeyword.size());
00062
00063 size_t numRegions = pvdoTables.size();
00064 setNumRegions(numRegions);
00065
00066 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
00067 Scalar rhoRefO = densityKeyword.getRecord(regionIdx).getItem("OIL").getSIDouble(0);
00068 Scalar rhoRefG = densityKeyword.getRecord(regionIdx).getItem("GAS").getSIDouble(0);
00069 Scalar rhoRefW = densityKeyword.getRecord(regionIdx).getItem("WATER").getSIDouble(0);
00070
00071 setReferenceDensities(regionIdx, rhoRefO, rhoRefG, rhoRefW);
00072
00073 const auto& pvdoTable = pvdoTables.getTable<PvdoTable>(regionIdx);
00074
00075 const auto& BColumn(pvdoTable.getFormationFactorColumn());
00076 std::vector<Scalar> invBColumn(BColumn.size());
00077 for (unsigned i = 0; i < invBColumn.size(); ++i)
00078 invBColumn[i] = 1/BColumn[i];
00079
00080 inverseOilB_[regionIdx].setXYArrays(pvdoTable.numRows(),
00081 pvdoTable.getPressureColumn(),
00082 invBColumn);
00083 oilMu_[regionIdx].setXYArrays(pvdoTable.numRows(),
00084 pvdoTable.getPressureColumn(),
00085 pvdoTable.getViscosityColumn());
00086 }
00087
00088 initEnd();
00089 }
00090 #endif // HAVE_OPM_PARSER
00091
00092 void setNumRegions(size_t numRegions)
00093 {
00094 oilReferenceDensity_.resize(numRegions);
00095 inverseOilB_.resize(numRegions);
00096 inverseOilBMu_.resize(numRegions);
00097 oilMu_.resize(numRegions);
00098 }
00099
00103 void setReferenceDensities(unsigned regionIdx,
00104 Scalar rhoRefOil,
00105 Scalar ,
00106 Scalar )
00107 {
00108 oilReferenceDensity_[regionIdx] = rhoRefOil;
00109 }
00110
00121 void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedOneDFunction& invBo)
00122 { inverseOilB_[regionIdx] = invBo; }
00123
00129 void setOilViscosity(unsigned regionIdx, const TabulatedOneDFunction& muo)
00130 { oilMu_[regionIdx] = muo; }
00131
00135 void initEnd()
00136 {
00137
00138 size_t numRegions = oilMu_.size();
00139 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
00140
00141
00142 const auto& oilMu = oilMu_[regionIdx];
00143 const auto& invOilB = inverseOilB_[regionIdx];
00144 assert(oilMu.numSamples() == invOilB.numSamples());
00145
00146 std::vector<Scalar> invBMuColumn;
00147 std::vector<Scalar> pressureColumn;
00148 invBMuColumn.resize(oilMu.numSamples());
00149 pressureColumn.resize(oilMu.numSamples());
00150
00151 for (unsigned pIdx = 0; pIdx < oilMu.numSamples(); ++pIdx) {
00152 pressureColumn[pIdx] = invOilB.xAt(pIdx);
00153 invBMuColumn[pIdx] = invOilB.valueAt(pIdx)*1/oilMu.valueAt(pIdx);
00154 }
00155
00156 inverseOilBMu_[regionIdx].setXYArrays(pressureColumn.size(),
00157 pressureColumn,
00158 invBMuColumn);
00159 }
00160 }
00161
00165 unsigned numRegions() const
00166 { return inverseOilBMu_.size(); }
00167
00171 template <class Evaluation>
00172 Evaluation viscosity(unsigned regionIdx,
00173 const Evaluation& temperature,
00174 const Evaluation& pressure,
00175 const Evaluation& ) const
00176 { return saturatedViscosity(regionIdx, temperature, pressure); }
00177
00181 template <class Evaluation>
00182 Evaluation saturatedViscosity(unsigned regionIdx,
00183 const Evaluation& ,
00184 const Evaluation& pressure) const
00185 {
00186 const Evaluation& invBo = inverseOilB_[regionIdx].eval(pressure, true);
00187 const Evaluation& invMuoBo = inverseOilBMu_[regionIdx].eval(pressure, true);
00188
00189 return invBo/invMuoBo;
00190 }
00191
00195 template <class Evaluation>
00196 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
00197 const Evaluation& ,
00198 const Evaluation& pressure,
00199 const Evaluation& ) const
00200 { return inverseOilB_[regionIdx].eval(pressure, true); }
00201
00207 template <class Evaluation>
00208 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
00209 const Evaluation& ,
00210 const Evaluation& pressure) const
00211 { return inverseOilB_[regionIdx].eval(pressure, true); }
00212
00216 template <class Evaluation>
00217 Evaluation saturatedGasDissolutionFactor(unsigned ,
00218 const Evaluation& ,
00219 const Evaluation& ) const
00220 { return 0.0; }
00221
00225 template <class Evaluation>
00226 Evaluation saturatedGasDissolutionFactor(unsigned ,
00227 const Evaluation& ,
00228 const Evaluation& ,
00229 const Evaluation& ,
00230 Scalar ) const
00231 { return 0.0; }
00232
00239 template <class Evaluation>
00240 Evaluation saturationPressure(unsigned ,
00241 const Evaluation& ,
00242 const Evaluation& ) const
00243 { return 0.0; }
00244
00245 template <class Evaluation>
00246 Evaluation saturatedGasMassFraction(unsigned ,
00247 const Evaluation& ,
00248 const Evaluation& ) const
00249 { return 0.0; }
00250
00251 template <class Evaluation>
00252 Evaluation saturatedGasMoleFraction(unsigned ,
00253 const Evaluation& ,
00254 const Evaluation& ) const
00255 { return 0.0; }
00256
00257 private:
00258 std::vector<Scalar> oilReferenceDensity_;
00259 std::vector<TabulatedOneDFunction> inverseOilB_;
00260 std::vector<TabulatedOneDFunction> oilMu_;
00261 std::vector<TabulatedOneDFunction> inverseOilBMu_;
00262 };
00263
00264 }
00265
00266 #endif