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_CHECK_FLUIDSYSTEM_HPP
00028 #define OPM_CHECK_FLUIDSYSTEM_HPP
00029
00030
00031 #include <opm/material/fluidsystems/SinglePhaseFluidSystem.hpp>
00032 #include <opm/material/fluidsystems/TwoPhaseImmiscibleFluidSystem.hpp>
00033 #include <opm/material/fluidsystems/H2ON2FluidSystem.hpp>
00034 #include <opm/material/fluidsystems/H2ON2LiquidPhaseFluidSystem.hpp>
00035 #include <opm/material/fluidsystems/H2OAirFluidSystem.hpp>
00036 #include <opm/material/fluidsystems/H2OAirMesityleneFluidSystem.hpp>
00037 #include <opm/material/fluidsystems/H2OAirXyleneFluidSystem.hpp>
00038 #include <opm/material/fluidsystems/Spe5FluidSystem.hpp>
00039
00040
00041 #include <opm/material/fluidstates/PressureOverlayFluidState.hpp>
00042 #include <opm/material/fluidstates/SaturationOverlayFluidState.hpp>
00043 #include <opm/material/fluidstates/TemperatureOverlayFluidState.hpp>
00044 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
00045 #include <opm/material/fluidstates/NonEquilibriumFluidState.hpp>
00046 #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
00047
00048 #include <opm/common/Unused.hpp>
00049 #include <dune/common/classname.hh>
00050
00051 #include <iostream>
00052 #include <string>
00053
00058 template <class ScalarT,
00059 class FluidSystem,
00060 class BaseFluidState = Opm::CompositionalFluidState<ScalarT, FluidSystem> >
00061 class HairSplittingFluidState
00062 : protected BaseFluidState
00063 {
00064 public:
00065 enum { numPhases = FluidSystem::numPhases };
00066 enum { numComponents = FluidSystem::numComponents };
00067
00068 typedef ScalarT Scalar;
00069
00070 HairSplittingFluidState()
00071 {
00072
00073 allowTemperature(false);
00074 allowPressure(false);
00075 allowComposition(false);
00076 allowDensity(false);
00077
00078
00079 restrictToPhase(1000);
00080 }
00081
00082 void allowTemperature(bool yesno)
00083 { allowTemperature_ = yesno; }
00084
00085 void allowPressure(bool yesno)
00086 { allowPressure_ = yesno; }
00087
00088 void allowComposition(bool yesno)
00089 { allowComposition_ = yesno; }
00090
00091 void allowDensity(bool yesno)
00092 { allowDensity_ = yesno; }
00093
00094 void restrictToPhase(int phaseIdx)
00095 { restrictPhaseIdx_ = phaseIdx; }
00096
00097 BaseFluidState& base()
00098 { return *static_cast<BaseFluidState*>(this); }
00099
00100 const BaseFluidState& base() const
00101 { return *static_cast<const BaseFluidState*>(this); }
00102
00103 auto temperature(unsigned phaseIdx) const
00104 -> decltype(this->base().temperature(phaseIdx))
00105 {
00106 assert(allowTemperature_);
00107 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
00108 return this->base().temperature(phaseIdx);
00109 }
00110
00111 auto pressure(unsigned phaseIdx) const
00112 -> decltype(this->base().pressure(phaseIdx))
00113 {
00114 assert(allowPressure_);
00115 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
00116 return this->base().pressure(phaseIdx);
00117 }
00118
00119 auto moleFraction(unsigned phaseIdx, unsigned compIdx) const
00120 -> decltype(this->base().moleFraction(phaseIdx, compIdx))
00121 {
00122 assert(allowComposition_);
00123 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
00124 return this->base().moleFraction(phaseIdx, compIdx);
00125 }
00126
00127 auto massFraction(unsigned phaseIdx, unsigned compIdx) const
00128 -> decltype(this->base().massFraction(phaseIdx, compIdx))
00129 {
00130 assert(allowComposition_);
00131 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
00132 return this->base().massFraction(phaseIdx, compIdx);
00133 }
00134
00135 auto averageMolarMass(unsigned phaseIdx) const
00136 -> decltype(this->base().averageMolarMass(phaseIdx))
00137 {
00138 assert(allowComposition_);
00139 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
00140 return this->base().averageMolarMass(phaseIdx);
00141 }
00142
00143 auto molarity(unsigned phaseIdx, unsigned compIdx) const
00144 -> decltype(this->base().molarity(phaseIdx, compIdx))
00145 {
00146 assert(allowDensity_ && allowComposition_);
00147 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
00148 return this->base().molarity(phaseIdx, compIdx);
00149 }
00150
00151 auto molarDensity(unsigned phaseIdx) const
00152 -> decltype(this->base().molarDensity(phaseIdx))
00153 {
00154 assert(allowDensity_);
00155 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
00156 return this->base().molarDensity(phaseIdx);
00157 }
00158
00159 auto molarVolume(unsigned phaseIdx) const
00160 -> decltype(this->base().molarVolume(phaseIdx))
00161 {
00162 assert(allowDensity_);
00163 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
00164 return this->base().molarVolume(phaseIdx);
00165 }
00166
00167 auto density(unsigned phaseIdx) const
00168 -> decltype(this->base().density(phaseIdx))
00169 {
00170 assert(allowDensity_);
00171 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
00172 return this->base().density(phaseIdx);
00173 }
00174
00175 auto saturation(unsigned phaseIdx) const
00176 -> decltype(this->base().saturation(phaseIdx))
00177 {
00178 assert(false);
00179 return this->base().saturation(phaseIdx);
00180 }
00181
00182 auto phaseIsPresent(unsigned phaseIdx) const
00183 -> decltype(this->base().phaseIsPresent(phaseIdx))
00184 {
00185 assert(false);
00186 return this->base().phaseIsPresent(phaseIdx);
00187 }
00188
00189 auto fugacity(unsigned phaseIdx, unsigned compIdx) const
00190 -> decltype(this->base().fugacity(phaseIdx, compIdx))
00191 {
00192 assert(false);
00193 return this->base().fugacity(phaseIdx, compIdx);
00194 }
00195
00196 auto fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
00197 -> decltype(this->base().fugacityCoefficient(phaseIdx, compIdx))
00198 {
00199 assert(false);
00200 return this->base().fugacityCoefficient(phaseIdx, compIdx);
00201 }
00202
00203 auto enthalpy(unsigned phaseIdx) const
00204 -> decltype(this->base().enthalpy(phaseIdx))
00205 {
00206 assert(false);
00207 return this->base().enthalpy(phaseIdx);
00208 }
00209
00210 auto internalEnergy(unsigned phaseIdx) const
00211 -> decltype(this->base().internalEnergy(phaseIdx))
00212 {
00213 assert(false);
00214 return this->base().internalEnergy(phaseIdx);
00215 }
00216
00217 auto viscosity(unsigned phaseIdx) const
00218 -> decltype(this->base().viscosity(phaseIdx))
00219 {
00220 assert(false);
00221 return this->base().viscosity(phaseIdx);
00222 }
00223
00224 private:
00225 bool allowSaturation_;
00226 bool allowTemperature_;
00227 bool allowPressure_;
00228 bool allowComposition_;
00229 bool allowDensity_;
00230 int restrictPhaseIdx_;
00231 };
00232
00233 template <class Scalar, class BaseFluidState>
00234 void checkFluidState(const BaseFluidState& fs)
00235 {
00236
00237 BaseFluidState tmpFs(fs);
00238 tmpFs = fs;
00239
00240
00241 fs.checkDefined();
00242
00243
00244 typedef typename BaseFluidState::Scalar FsScalar;
00245 static_assert(std::is_same<FsScalar, Scalar>::value,
00246 "Fluid states must export the type they are given as scalar in an unmodified way");
00247
00248
00249 while (false) {
00250 Scalar val = 1.0;
00251
00252 val = 2*val;
00253
00254 val = fs.temperature(0);
00255 val = fs.pressure(0);
00256 val = fs.moleFraction(0, 0);
00257 val = fs.massFraction(0, 0);
00258 val = fs.averageMolarMass(0);
00259 val = fs.molarity(0, 0);
00260 val = fs.molarDensity(0);
00261 val = fs.molarVolume(0);
00262 val = fs.density(0);
00263 val = fs.saturation(0);
00264 bool b OPM_UNUSED = fs.phaseIsPresent(0);
00265 val = fs.fugacity(0, 0);
00266 val = fs.fugacityCoefficient(0, 0);
00267 val = fs.enthalpy(0);
00268 val = fs.internalEnergy(0);
00269 val = fs.viscosity(0);
00270 };
00271 }
00272
00276 template <class Scalar, class FluidSystem, class RhsEval, class LhsEval>
00277 void checkFluidSystem()
00278 {
00279 std::cout << "Testing fluid system '" << Dune::className<FluidSystem>() << "'\n";
00280
00281
00282
00283 enum { numPhases = FluidSystem::numPhases };
00284 enum { numComponents = FluidSystem::numComponents };
00285
00286 typedef HairSplittingFluidState<RhsEval, FluidSystem> FluidState;
00287 FluidState fs;
00288 fs.allowTemperature(true);
00289 fs.allowPressure(true);
00290 fs.allowComposition(true);
00291 fs.restrictToPhase(-1);
00292
00293
00294 fs.base().setTemperature(273.15 + 20.0);
00295 for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
00296 fs.base().setPressure(phaseIdx, 1e5);
00297 fs.base().setSaturation(phaseIdx, 1.0/numPhases);
00298 for (int compIdx = 0; compIdx < numComponents; ++ compIdx) {
00299 fs.base().setMoleFraction(phaseIdx, compIdx, 1.0/numComponents);
00300 }
00301 }
00302
00303 static_assert(std::is_same<typename FluidSystem::Scalar, Scalar>::value,
00304 "The type used for floating point used by the fluid system must be the same"
00305 " as the one passed to the checkFluidSystem() function");
00306
00307
00308 typedef typename FluidSystem::template ParameterCache<LhsEval> ParameterCache;
00309
00310 ParameterCache paramCache;
00311 try { paramCache.updateAll(fs); } catch (...) {};
00312 try { paramCache.updateAll(fs, ParameterCache::None); } catch (...) {};
00313 try { paramCache.updateAll(fs, ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); } catch (...) {};
00314 try { paramCache.updateAllPressures(fs); } catch (...) {};
00315
00316 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
00317 fs.restrictToPhase(static_cast<int>(phaseIdx));
00318 try { paramCache.updatePhase(fs, phaseIdx); } catch (...) {};
00319 try { paramCache.updatePhase(fs, phaseIdx, ParameterCache::None); } catch (...) {};
00320 try { paramCache.updatePhase(fs, phaseIdx, ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); } catch (...) {};
00321 try { paramCache.updateTemperature(fs, phaseIdx); } catch (...) {};
00322 try { paramCache.updatePressure(fs, phaseIdx); } catch (...) {};
00323 try { paramCache.updateComposition(fs, phaseIdx); } catch (...) {};
00324 try { paramCache.updateSingleMoleFraction(fs, phaseIdx, 0); } catch (...) {};
00325 }
00326
00327
00328
00329 LhsEval val = 0.0;
00330 Scalar scalarVal = 0.0;
00331
00332 scalarVal = 2*scalarVal;
00333 val = 2*val;
00334
00335
00336 try { FluidSystem::init(); } catch (...) {};
00337 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
00338 fs.restrictToPhase(static_cast<int>(phaseIdx));
00339 fs.allowPressure(FluidSystem::isCompressible(phaseIdx));
00340 fs.allowComposition(true);
00341 fs.allowDensity(false);
00342 try { auto tmpVal = FluidSystem::density(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value, "The default return value must be the scalar used by the fluid state!"); } catch (...) {};
00343 try { val = FluidSystem::template density<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
00344 try { scalarVal = FluidSystem::template density<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
00345
00346 fs.allowPressure(true);
00347 fs.allowDensity(true);
00348 try { auto tmpVal = FluidSystem::viscosity(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value, "The default return value must be the scalar used by the fluid state!"); } catch (...) {};
00349 try { auto tmpVal = FluidSystem::enthalpy(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value, "The default return value must be the scalar used by the fluid state!"); } catch (...) {};
00350 try { auto tmpVal = FluidSystem::heatCapacity(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value, "The default return value must be the scalar used by the fluid state!"); } catch (...) {};
00351 try { auto tmpVal = FluidSystem::thermalConductivity(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value, "The default return value must be the scalar used by the fluid state!"); } catch (...) {};
00352 try { val = FluidSystem::template viscosity<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
00353 try { val = FluidSystem::template enthalpy<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
00354 try { val = FluidSystem::template heatCapacity<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
00355 try { val = FluidSystem::template thermalConductivity<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
00356 try { scalarVal = FluidSystem::template viscosity<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
00357 try { scalarVal = FluidSystem::template enthalpy<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
00358 try { scalarVal = FluidSystem::template heatCapacity<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
00359 try { scalarVal = FluidSystem::template thermalConductivity<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
00360
00361 for (unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
00362 fs.allowComposition(!FluidSystem::isIdealMixture(phaseIdx));
00363 try { auto tmpVal = FluidSystem::fugacityCoefficient(fs, paramCache, phaseIdx, compIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value, "The default return value must be the scalar used by the fluid state!"); } catch (...) {};
00364 try { val = FluidSystem::template fugacityCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); } catch (...) {};
00365 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); } catch (...) {};
00366 fs.allowComposition(true);
00367 try { auto tmpVal = FluidSystem::diffusionCoefficient(fs, paramCache, phaseIdx, compIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value, "The default return value must be the scalar used by the fluid state!"); } catch (...) {};
00368 try { val = FluidSystem::template diffusionCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); } catch (...) {};
00369 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); } catch (...) {};
00370 }
00371 }
00372
00373
00374 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
00375 std::string name OPM_UNUSED = FluidSystem::phaseName(phaseIdx);
00376 bool bVal = FluidSystem::isLiquid(phaseIdx);
00377 bVal = FluidSystem::isIdealGas(phaseIdx);
00378 bVal = !bVal;
00379 }
00380
00381
00382 for (unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
00383 val = FluidSystem::molarMass(compIdx);
00384 std::string name = FluidSystem::componentName(compIdx);
00385 }
00386
00387 std::cout << "----------------------------------\n";
00388 }
00389
00390 #endif