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_WATER_PVT_THERMAL_HPP
00028 #define OPM_WATER_PVT_THERMAL_HPP
00029
00030 #include <opm/material/Constants.hpp>
00031
00032 #include <opm/material/common/OpmFinal.hpp>
00033 #include <opm/material/common/UniformXTabulated2DFunction.hpp>
00034 #include <opm/material/common/Tabulated1DFunction.hpp>
00035 #include <opm/material/common/Spline.hpp>
00036
00037 #if HAVE_OPM_PARSER
00038 #include <opm/parser/eclipse/Deck/Deck.hpp>
00039 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
00040 #include <opm/parser/eclipse/EclipseState/Tables/SimpleTable.hpp>
00041 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
00042 #endif
00043
00044 namespace Opm {
00045 template <class Scalar, bool enableThermal>
00046 class WaterPvtMultiplexer;
00047
00054 template <class Scalar>
00055 class WaterPvtThermal
00056 {
00057 typedef Opm::Tabulated1DFunction<Scalar> TabulatedOneDFunction;
00058 typedef WaterPvtMultiplexer<Scalar, false> IsothermalPvt;
00059
00060 public:
00061 ~WaterPvtThermal()
00062 { delete isothermalPvt_; }
00063
00064 #if HAVE_OPM_PARSER
00065
00068 void initFromDeck(const Deck& deck,
00069 const EclipseState& eclState)
00070 {
00072
00074 isothermalPvt_ = new IsothermalPvt;
00075 isothermalPvt_->initFromDeck(deck, eclState);
00076
00078
00080 const auto& tables = eclState.getTableManager();
00081
00082 enableThermalDensity_ = deck.hasKeyword("WATDENT");
00083 enableThermalViscosity_ = deck.hasKeyword("VISCREF");
00084
00085 unsigned numRegions = isothermalPvt_->numRegions();
00086 setNumRegions(numRegions);
00087
00088 if (enableThermalDensity_) {
00089 const auto& watdentKeyword = deck.getKeyword("WATDENT");
00090
00091 assert(watdentKeyword.size() == numRegions);
00092 for (unsigned regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
00093 const auto& watdentRecord = watdentKeyword.getRecord(regionIdx);
00094
00095 watdentRefTemp_[regionIdx] = watdentRecord.getItem("REFERENCE_TEMPERATURE").getSIDouble(0);
00096 watdentCT1_[regionIdx] = watdentRecord.getItem("EXPANSION_COEFF_LINEAR").getSIDouble(0);
00097 watdentCT2_[regionIdx] = watdentRecord.getItem("EXPANSION_COEFF_QUADRATIC").getSIDouble(0);
00098 }
00099 }
00100
00101 if (enableThermalViscosity_) {
00102 const auto& viscrefKeyword = deck.getKeyword("VISCREF");
00103
00104 const auto& watvisctTables = tables.getWatvisctTables();
00105
00106 assert(watvisctTables.size() == numRegions);
00107 assert(viscrefKeyword.size() == numRegions);
00108
00109 for (unsigned regionIdx = 0; regionIdx < numRegions; ++ regionIdx) {
00110 const auto& T = watvisctTables[regionIdx].getColumn("Temperature").vectorCopy();
00111 const auto& mu = watvisctTables[regionIdx].getColumn("Viscosity").vectorCopy();
00112 watvisctCurves_[regionIdx].setXYContainers(T, mu);
00113
00114 const auto& viscrefRecord = viscrefKeyword.getRecord(regionIdx);
00115 viscrefPress_[regionIdx] = viscrefRecord.getItem("REFERENCE_PRESSURE").getSIDouble(0);
00116 }
00117 }
00118 }
00119 #endif // HAVE_OPM_PARSER
00120
00124 void setNumRegions(size_t numRegions)
00125 {
00126 pvtwRefPress_.resize(numRegions);
00127 pvtwRefB_.resize(numRegions);
00128 pvtwCompressibility_.resize(numRegions);
00129 pvtwViscosity_.resize(numRegions);
00130 pvtwViscosibility_.resize(numRegions);
00131 viscrefPress_.resize(numRegions);
00132 watvisctCurves_.resize(numRegions);
00133 watdentRefTemp_.resize(numRegions);
00134 watdentCT1_.resize(numRegions);
00135 watdentCT2_.resize(numRegions);
00136 }
00137
00141 void initEnd()
00142 { }
00143
00147 bool enableThermalDensity() const
00148 { return enableThermalDensity_; }
00149
00153 bool enableThermalViscosity() const
00154 { return enableThermalViscosity_; }
00155
00156 size_t numRegions() const
00157 { return pvtwRefPress_.size(); }
00158
00162 template <class Evaluation>
00163 Evaluation viscosity(unsigned regionIdx,
00164 const Evaluation& temperature,
00165 const Evaluation& pressure) const
00166 {
00167 const auto& isothermalMu = isothermalPvt_->viscosity(regionIdx, temperature, pressure);
00168 if (!enableThermalViscosity())
00169 return isothermalMu;
00170
00171 Scalar x = -pvtwViscosibility_[regionIdx]*(viscrefPress_[regionIdx] - pvtwRefPress_[regionIdx]);
00172 Scalar muRef = pvtwViscosity_[regionIdx]/(1.0 + x + 0.5*x*x);
00173
00174
00175 const auto& muWatvisct = watvisctCurves_[regionIdx].eval(temperature);
00176 return isothermalMu * muWatvisct/muRef;
00177 }
00178
00182 template <class Evaluation>
00183 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
00184 const Evaluation& temperature,
00185 const Evaluation& pressure) const
00186 {
00187 if (!enableThermalDensity())
00188 return isothermalPvt_->inverseFormationVolumeFactor(regionIdx, temperature, pressure);
00189
00190 Scalar BwRef = pvtwRefB_[regionIdx];
00191 Scalar TRef = watdentRefTemp_[regionIdx];
00192 const Evaluation& X = pvtwCompressibility_[regionIdx]*(pressure - pvtwRefPress_[regionIdx]);
00193 Scalar cT1 = watdentCT1_[regionIdx];
00194 Scalar cT2 = watdentCT2_[regionIdx];
00195 const Evaluation& Y = temperature - TRef;
00196
00197 return ((1 - X)*(1 + cT1*Y + cT2*Y*Y))/BwRef;
00198 }
00199
00200 private:
00201 IsothermalPvt* isothermalPvt_;
00202
00203
00204
00205 std::vector<Scalar> viscrefPress_;
00206
00207 std::vector<Scalar> watdentRefTemp_;
00208 std::vector<Scalar> watdentCT1_;
00209 std::vector<Scalar> watdentCT2_;
00210
00211 std::vector<Scalar> pvtwRefPress_;
00212 std::vector<Scalar> pvtwRefB_;
00213 std::vector<Scalar> pvtwCompressibility_;
00214 std::vector<Scalar> pvtwViscosity_;
00215 std::vector<Scalar> pvtwViscosibility_;
00216
00217 std::vector<TabulatedOneDFunction> watvisctCurves_;
00218
00219 bool enableThermalDensity_;
00220 bool enableThermalViscosity_;
00221 };
00222
00223 }
00224
00225 #endif