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_CONSTANT_COMPRESSIBILITY_WATER_PVT_HPP
00028 #define OPM_CONSTANT_COMPRESSIBILITY_WATER_PVT_HPP
00029
00030 #include <opm/material/common/Tabulated1DFunction.hpp>
00031
00032 #if HAVE_OPM_PARSER
00033 #include <opm/parser/eclipse/Deck/Deck.hpp>
00034 #include <opm/parser/eclipse/Deck/DeckKeyword.hpp>
00035 #include <opm/parser/eclipse/Deck/DeckRecord.hpp>
00036 #include <opm/parser/eclipse/Deck/DeckItem.hpp>
00037 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
00038 #endif
00039
00040 #include <vector>
00041
00042 namespace Opm {
00047 template <class Scalar>
00048 class ConstantCompressibilityWaterPvt
00049 {
00050 typedef Opm::Tabulated1DFunction<Scalar> TabulatedOneDFunction;
00051 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
00052
00053 public:
00054 #if HAVE_OPM_PARSER
00055
00059 void initFromDeck(const Deck& deck, const EclipseState& )
00060 {
00061 const auto& pvtwKeyword = deck.getKeyword("PVTW");
00062 const auto& densityKeyword = deck.getKeyword("DENSITY");
00063
00064 assert(pvtwKeyword.size() == densityKeyword.size());
00065
00066 size_t numRegions = pvtwKeyword.size();
00067 setNumRegions(numRegions);
00068
00069 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
00070 auto pvtwRecord = pvtwKeyword.getRecord(regionIdx);
00071 auto densityRecord = densityKeyword.getRecord(regionIdx);
00072
00073 waterReferenceDensity_[regionIdx] =
00074 densityRecord.getItem("WATER").getSIDouble(0);
00075
00076 waterReferencePressure_[regionIdx] =
00077 pvtwRecord.getItem("P_REF").getSIDouble(0);
00078 waterReferenceFormationVolumeFactor_[regionIdx] =
00079 pvtwRecord.getItem("WATER_VOL_FACTOR").getSIDouble(0);
00080 waterCompressibility_[regionIdx] =
00081 pvtwRecord.getItem("WATER_COMPRESSIBILITY").getSIDouble(0);
00082 waterViscosity_[regionIdx] =
00083 pvtwRecord.getItem("WATER_VISCOSITY").getSIDouble(0);
00084 waterViscosibility_[regionIdx] =
00085 pvtwRecord.getItem("WATER_VISCOSIBILITY").getSIDouble(0);
00086 }
00087
00088 initEnd();
00089 }
00090 #endif
00091
00092 void setNumRegions(size_t numRegions)
00093 {
00094 waterReferenceDensity_.resize(numRegions);
00095 waterReferencePressure_.resize(numRegions);
00096 waterReferenceFormationVolumeFactor_.resize(numRegions);
00097 waterCompressibility_.resize(numRegions);
00098 waterViscosity_.resize(numRegions);
00099 waterViscosibility_.resize(numRegions);
00100
00101 for (unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
00102 setReferenceDensities(regionIdx, 650.0, 1.0, 1000.0);
00103 setReferenceFormationVolumeFactor(regionIdx, 1.0);
00104 setReferencePressure(regionIdx, 1e5);
00105 }
00106 }
00107
00111 void setReferenceDensities(unsigned regionIdx,
00112 Scalar ,
00113 Scalar ,
00114 Scalar rhoRefWater)
00115 { waterReferenceDensity_[regionIdx] = rhoRefWater; }
00116
00120 void setReferencePressure(unsigned regionIdx, Scalar p)
00121 { waterReferencePressure_[regionIdx] = p; }
00122
00126 void setViscosity(unsigned regionIdx, Scalar muw, Scalar waterViscosibility = 0.0)
00127 {
00128 waterViscosity_[regionIdx] = muw;
00129 waterViscosibility_[regionIdx] = waterViscosibility;
00130 }
00131
00135 void setCompressibility(unsigned regionIdx, Scalar waterCompressibility)
00136 { waterCompressibility_[regionIdx] = waterCompressibility; }
00137
00141 void setReferenceFormationVolumeFactor(unsigned regionIdx, Scalar BwRef)
00142 { waterReferenceFormationVolumeFactor_[regionIdx] = BwRef; }
00143
00147 void setViscosibility(unsigned regionIdx, Scalar muComp)
00148 { waterViscosibility_[regionIdx] = muComp; }
00149
00153 void initEnd()
00154 { }
00155
00159 unsigned numRegions() const
00160 { return waterReferenceDensity_.size(); }
00161
00165 template <class Evaluation>
00166 Evaluation viscosity(unsigned regionIdx,
00167 const Evaluation& temperature,
00168 const Evaluation& pressure) const
00169 {
00170 Scalar BwMuwRef = waterViscosity_[regionIdx]*waterReferenceFormationVolumeFactor_[regionIdx];
00171 const Evaluation& bw = inverseFormationVolumeFactor(regionIdx, temperature, pressure);
00172
00173 Scalar pRef = waterReferencePressure_[regionIdx];
00174 const Evaluation& Y =
00175 (waterCompressibility_[regionIdx] - waterViscosibility_[regionIdx])
00176 * (pressure - pRef);
00177 return BwMuwRef*bw/(1 + Y*(1 + Y/2));
00178 }
00179
00183 template <class Evaluation>
00184 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
00185 const Evaluation& ,
00186 const Evaluation& pressure) const
00187 {
00188
00189 Scalar pRef = waterReferencePressure_[regionIdx];
00190 const Evaluation& X = waterCompressibility_[regionIdx]*(pressure - pRef);
00191
00192 Scalar BwRef = waterReferenceFormationVolumeFactor_[regionIdx];
00193
00194
00195 return (1.0 + X*(1.0 + X/2.0))/BwRef;
00196 }
00197
00198 private:
00199 std::vector<Scalar> waterReferenceDensity_;
00200 std::vector<Scalar> waterReferencePressure_;
00201 std::vector<Scalar> waterReferenceFormationVolumeFactor_;
00202 std::vector<Scalar> waterCompressibility_;
00203 std::vector<Scalar> waterViscosity_;
00204 std::vector<Scalar> waterViscosibility_;
00205 };
00206
00207 }
00208
00209 #endif