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_SOLVENT_PVT_HPP
00028 #define OPM_SOLVENT_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/PvdsTable.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 SolventPvt
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& pvdsTables = eclState.getTableManager().getPvdsTables();
00066 const auto& sdensityKeyword = deck.getKeyword("SDENSITY");
00067
00068 assert(pvdsTables.size() == sdensityKeyword.size());
00069
00070 size_t numRegions = pvdsTables.size();
00071 setNumRegions(numRegions);
00072
00073 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
00074 Scalar rhoRefS = sdensityKeyword.getRecord(regionIdx).getItem("SOLVENT_DENSITY").getSIDouble(0);
00075
00076 setReferenceDensity(regionIdx, rhoRefS);
00077
00078 const auto& pvdsTable = pvdsTables.getTable<PvdsTable>(regionIdx);
00079
00080
00081
00082 std::vector<Scalar> invB(pvdsTable.numRows());
00083 const auto& Bg = pvdsTable.getFormationFactorColumn();
00084 for (unsigned i = 0; i < Bg.size(); ++ i) {
00085 invB[i] = 1.0/Bg[i];
00086 }
00087
00088 size_t numSamples = invB.size();
00089 inverseSolventB_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), invB);
00090 solventMu_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), pvdsTable.getViscosityColumn());
00091 }
00092
00093 initEnd();
00094 }
00095 #endif
00096
00097 void setNumRegions(size_t numRegions)
00098 {
00099 solventReferenceDensity_.resize(numRegions);
00100 inverseSolventB_.resize(numRegions);
00101 inverseSolventBMu_.resize(numRegions);
00102 solventMu_.resize(numRegions);
00103 }
00104
00108 void setReferenceDensity(unsigned regionIdx, Scalar rhoRefSolvent)
00109 { solventReferenceDensity_[regionIdx] = rhoRefSolvent; }
00110
00116 void setSolventViscosity(unsigned regionIdx, const TabulatedOneDFunction& mug)
00117 { solventMu_[regionIdx] = mug; }
00118
00124 void setSolventFormationVolumeFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
00125 {
00126 SamplingPoints tmp(samplePoints);
00127 auto it = tmp.begin();
00128 const auto& endIt = tmp.end();
00129 for (; it != endIt; ++ it)
00130 std::get<1>(*it) = 1.0/std::get<1>(*it);
00131
00132 inverseSolventB_[regionIdx].setContainerOfTuples(tmp);
00133 assert(inverseSolventB_[regionIdx].monotonic());
00134 }
00135
00139 void initEnd()
00140 {
00141
00142 size_t numRegions = solventMu_.size();
00143 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
00144
00145
00146 const auto& solventMu = solventMu_[regionIdx];
00147 const auto& invSolventB = inverseSolventB_[regionIdx];
00148 assert(solventMu.numSamples() == invSolventB.numSamples());
00149
00150 std::vector<Scalar> pressureValues(solventMu.numSamples());
00151 std::vector<Scalar> invSolventBMuValues(solventMu.numSamples());
00152 for (unsigned pIdx = 0; pIdx < solventMu.numSamples(); ++pIdx) {
00153 pressureValues[pIdx] = invSolventB.xAt(pIdx);
00154 invSolventBMuValues[pIdx] = invSolventB.valueAt(pIdx) * (1.0/solventMu.valueAt(pIdx));
00155 }
00156
00157 inverseSolventBMu_[regionIdx].setXYContainers(pressureValues, invSolventBMuValues);
00158 }
00159 }
00160
00164 unsigned numRegions() const
00165 { return solventReferenceDensity_.size(); }
00166
00170 template <class Evaluation>
00171 Evaluation viscosity(unsigned regionIdx,
00172 const Evaluation& temperature OPM_UNUSED,
00173 const Evaluation& pressure) const
00174 {
00175 const Evaluation& invBg = inverseSolventB_[regionIdx].eval(pressure, true);
00176 const Evaluation& invMugBg = inverseSolventBMu_[regionIdx].eval(pressure, true);
00177
00178 return invBg/invMugBg;
00179 }
00180
00184 Scalar referenceDensity(unsigned regionIdx) const
00185 { return solventReferenceDensity_[regionIdx]; }
00186
00190 template <class Evaluation>
00191 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
00192 const Evaluation& temperature OPM_UNUSED,
00193 const Evaluation& pressure) const
00194 { return inverseSolventB_[regionIdx].eval(pressure, true); }
00195
00196 private:
00197 std::vector<Scalar> solventReferenceDensity_;
00198 std::vector<TabulatedOneDFunction> inverseSolventB_;
00199 std::vector<TabulatedOneDFunction> solventMu_;
00200 std::vector<TabulatedOneDFunction> inverseSolventBMu_;
00201 };
00202
00203 }
00204
00205 #endif