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_SPE5_PARAMETER_CACHE_HPP
00028 #define OPM_SPE5_PARAMETER_CACHE_HPP
00029
00030 #include <cassert>
00031
00032 #include <opm/material/components/H2O.hpp>
00033 #include <opm/material/fluidsystems/ParameterCacheBase.hpp>
00034
00035 #include <opm/material/eos/PengRobinson.hpp>
00036 #include <opm/material/eos/PengRobinsonParamsMixture.hpp>
00037
00038 namespace Opm {
00039
00044 template <class Scalar, class FluidSystem>
00045 class Spe5ParameterCache
00046 : public Opm::ParameterCacheBase<Spe5ParameterCache<Scalar, FluidSystem> >
00047 {
00048 typedef Spe5ParameterCache<Scalar, FluidSystem> ThisType;
00049 typedef Opm::ParameterCacheBase<ThisType> ParentType;
00050
00051 typedef Opm::PengRobinson<Scalar> PengRobinson;
00052
00053 enum { numPhases = FluidSystem::numPhases };
00054
00055 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
00056 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
00057 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
00058
00059 public:
00061 typedef Opm::PengRobinsonParamsMixture<Scalar, FluidSystem, oilPhaseIdx, true> OilPhaseParams;
00063 typedef Opm::PengRobinsonParamsMixture<Scalar, FluidSystem, gasPhaseIdx, true> GasPhaseParams;
00064
00065 Spe5ParameterCache()
00066 {
00067 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
00068 VmUpToDate_[phaseIdx] = false;
00069 Valgrind::SetUndefined(Vm_[phaseIdx]);
00070 }
00071 }
00072
00074 template <class FluidState>
00075 void updatePhase(const FluidState& fluidState,
00076 unsigned phaseIdx,
00077 int exceptQuantities = ParentType::None)
00078 {
00079 updateEosParams(fluidState, phaseIdx, exceptQuantities);
00080
00081
00082
00083 if (VmUpToDate_[phaseIdx])
00084 return;
00085
00086
00087 updateMolarVolume_(fluidState, phaseIdx);
00088 }
00089
00091 template <class FluidState>
00092 void updateSingleMoleFraction(const FluidState& fluidState,
00093 unsigned phaseIdx,
00094 unsigned compIdx)
00095 {
00096 if (phaseIdx == oilPhaseIdx)
00097 oilPhaseParams_.updateSingleMoleFraction(fluidState, compIdx);
00098 else if (phaseIdx == gasPhaseIdx)
00099 gasPhaseParams_.updateSingleMoleFraction(fluidState, compIdx);
00100
00101
00102 updateMolarVolume_(fluidState, phaseIdx);
00103 }
00104
00110 Scalar a(unsigned phaseIdx) const
00111 {
00112 switch (phaseIdx)
00113 {
00114 case oilPhaseIdx: return oilPhaseParams_.a();
00115 case gasPhaseIdx: return gasPhaseParams_.a();
00116 default:
00117 OPM_THROW(std::logic_error,
00118 "The a() parameter is only defined for "
00119 "oil and gas phases");
00120 };
00121 }
00122
00128 Scalar b(unsigned phaseIdx) const
00129 {
00130 switch (phaseIdx)
00131 {
00132 case oilPhaseIdx: return oilPhaseParams_.b();
00133 case gasPhaseIdx: return gasPhaseParams_.b();
00134 default:
00135 OPM_THROW(std::logic_error,
00136 "The b() parameter is only defined for "
00137 "oil and gas phases");
00138 };
00139 }
00140
00149 Scalar aPure(unsigned phaseIdx, unsigned compIdx) const
00150 {
00151 switch (phaseIdx)
00152 {
00153 case oilPhaseIdx: return oilPhaseParams_.pureParams(compIdx).a();
00154 case gasPhaseIdx: return gasPhaseParams_.pureParams(compIdx).a();
00155 default:
00156 OPM_THROW(std::logic_error,
00157 "The a() parameter is only defined for "
00158 "oil and gas phases");
00159 };
00160 }
00161
00169 Scalar bPure(unsigned phaseIdx, unsigned compIdx) const
00170 {
00171 switch (phaseIdx)
00172 {
00173 case oilPhaseIdx: return oilPhaseParams_.pureParams(compIdx).b();
00174 case gasPhaseIdx: return gasPhaseParams_.pureParams(compIdx).b();
00175 default:
00176 OPM_THROW(std::logic_error,
00177 "The b() parameter is only defined for "
00178 "oil and gas phases");
00179 };
00180 }
00181
00187 Scalar molarVolume(unsigned phaseIdx) const
00188 { assert(VmUpToDate_[phaseIdx]); return Vm_[phaseIdx]; }
00189
00190
00195 const OilPhaseParams& oilPhaseParams() const
00196 { return oilPhaseParams_; }
00197
00202 const GasPhaseParams& gasPhaseParams() const
00203 { return gasPhaseParams_; }
00204
00213 template <class FluidState>
00214 void updateEosParams(const FluidState& fluidState,
00215 unsigned phaseIdx,
00216 int exceptQuantities = ParentType::None)
00217 {
00218 if (!(exceptQuantities & ParentType::Temperature))
00219 {
00220 updatePure_(fluidState, phaseIdx);
00221 updateMix_(fluidState, phaseIdx);
00222 VmUpToDate_[phaseIdx] = false;
00223 }
00224 else if (!(exceptQuantities & ParentType::Composition))
00225 {
00226 updateMix_(fluidState, phaseIdx);
00227 VmUpToDate_[phaseIdx] = false;
00228 }
00229 else if (!(exceptQuantities & ParentType::Pressure)) {
00230 VmUpToDate_[phaseIdx] = false;
00231 }
00232 }
00233
00234 protected:
00241 template <class FluidState>
00242 void updatePure_(const FluidState& fluidState, unsigned phaseIdx)
00243 {
00244 Scalar T = fluidState.temperature(phaseIdx);
00245 Scalar p = fluidState.pressure(phaseIdx);
00246
00247 switch (phaseIdx)
00248 {
00249 case oilPhaseIdx: oilPhaseParams_.updatePure(T, p); break;
00250 case gasPhaseIdx: gasPhaseParams_.updatePure(T, p); break;
00251
00252 }
00253 }
00254
00262 template <class FluidState>
00263 void updateMix_(const FluidState& fluidState, unsigned phaseIdx)
00264 {
00265 Valgrind::CheckDefined(fluidState.averageMolarMass(phaseIdx));
00266 switch (phaseIdx)
00267 {
00268 case oilPhaseIdx:
00269 oilPhaseParams_.updateMix(fluidState);
00270 break;
00271 case gasPhaseIdx:
00272 gasPhaseParams_.updateMix(fluidState);
00273 break;
00274 case waterPhaseIdx:
00275 break;
00276 }
00277 }
00278
00279 template <class FluidState>
00280 void updateMolarVolume_(const FluidState& fluidState,
00281 unsigned phaseIdx)
00282 {
00283 VmUpToDate_[phaseIdx] = true;
00284
00285
00286
00287 switch (phaseIdx) {
00288 case gasPhaseIdx: {
00289
00290
00291
00292
00293
00294 Vm_[gasPhaseIdx] =
00295 PengRobinson::computeMolarVolume(fluidState,
00296 *this,
00297 phaseIdx,
00298 true);
00299 break;
00300 }
00301 case oilPhaseIdx: {
00302
00303
00304
00305
00306
00307 Vm_[oilPhaseIdx] =
00308 PengRobinson::computeMolarVolume(fluidState,
00309 *this,
00310 phaseIdx,
00311 false);
00312
00313 break;
00314 }
00315 case waterPhaseIdx: {
00316
00317
00318
00319 const Scalar stockTankWaterDensity = 62.4 * 0.45359237 / 0.028316847;
00320
00321
00322 Scalar overPressure = fluidState.pressure(waterPhaseIdx) - 1.013e5;
00323 Scalar waterDensity =
00324 stockTankWaterDensity * (1 + 3.3e-6*overPressure/6894.7573);
00325
00326
00327 Vm_[waterPhaseIdx] = fluidState.averageMolarMass(waterPhaseIdx)/waterDensity;
00328 break;
00329 };
00330 };
00331 }
00332
00333 bool VmUpToDate_[numPhases];
00334 Scalar Vm_[numPhases];
00335
00336 OilPhaseParams oilPhaseParams_;
00337 GasPhaseParams gasPhaseParams_;
00338 };
00339
00340 }
00341
00342 #endif