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_SATURATION_OVERLAY_FLUID_STATE_HPP
00028 #define OPM_SATURATION_OVERLAY_FLUID_STATE_HPP
00029
00030 #include <opm/common/Valgrind.hpp>
00031
00032 #include <array>
00033 #include <utility>
00034
00035 namespace Opm {
00036
00042 template <class FluidState>
00043 class SaturationOverlayFluidState
00044 {
00045 public:
00046 typedef typename FluidState::Scalar Scalar;
00047
00048 enum { numPhases = FluidState::numPhases };
00049 enum { numComponents = FluidState::numComponents };
00050
00058 SaturationOverlayFluidState(const FluidState& fs)
00059 : fs_(&fs)
00060 {
00061 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
00062 saturation_[phaseIdx] = fs.saturation(phaseIdx);
00063 }
00064
00065
00066 SaturationOverlayFluidState(const SaturationOverlayFluidState& fs)
00067 : fs_(fs.fs_)
00068 , saturation_(fs.saturation_)
00069 {}
00070
00071
00072 SaturationOverlayFluidState& operator=(const SaturationOverlayFluidState& fs)
00073 {
00074 fs_ = fs.fs_;
00075 saturation_ = fs.saturation_;
00076 return *this;
00077 }
00078
00079
00080
00081
00082
00086 auto saturation(unsigned phaseIdx) const
00087 -> decltype(std::declval<FluidState>().saturation(phaseIdx))
00088 { return saturation_[phaseIdx]; }
00089
00093 bool phaseIsPresent(unsigned phaseIdx) const
00094 { return saturation_[phaseIdx] > 0.0; }
00095
00099 auto moleFraction(unsigned phaseIdx, unsigned compIdx) const
00100 -> decltype(std::declval<FluidState>().moleFraction(phaseIdx, compIdx))
00101 { return fs_->moleFraction(phaseIdx, compIdx); }
00102
00106 auto massFraction(unsigned phaseIdx, unsigned compIdx) const
00107 -> decltype(std::declval<FluidState>().massFraction(phaseIdx, compIdx))
00108 { return fs_->massFraction(phaseIdx, compIdx); }
00109
00118 auto averageMolarMass(unsigned phaseIdx) const
00119 -> decltype(std::declval<FluidState>().averageMolarMass(phaseIdx))
00120 { return fs_->averageMolarMass(phaseIdx); }
00121
00131 auto molarity(unsigned phaseIdx, unsigned compIdx) const
00132 -> decltype(std::declval<FluidState>().molarity(phaseIdx, compIdx))
00133 { return fs_->molarity(phaseIdx, compIdx); }
00134
00138 auto fugacity(unsigned phaseIdx, unsigned compIdx) const
00139 -> decltype(std::declval<FluidState>().fugacity(phaseIdx, compIdx))
00140 { return fs_->fugacity(phaseIdx, compIdx); }
00141
00145 auto fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
00146 -> decltype(std::declval<FluidState>().fugacityCoefficient(phaseIdx, compIdx))
00147 { return fs_->fugacityCoefficient(phaseIdx, compIdx); }
00148
00152 auto molarVolume(unsigned phaseIdx) const
00153 -> decltype(std::declval<FluidState>().molarVolume(phaseIdx))
00154 { return fs_->molarVolume(phaseIdx); }
00155
00159 auto density(unsigned phaseIdx) const
00160 -> decltype(std::declval<FluidState>().density(phaseIdx))
00161 { return fs_->density(phaseIdx); }
00162
00166 auto molarDensity(unsigned phaseIdx) const
00167 -> decltype(std::declval<FluidState>().molarDensity(phaseIdx))
00168 { return fs_->molarDensity(phaseIdx); }
00169
00173 auto temperature(unsigned phaseIdx) const
00174 -> decltype(std::declval<FluidState>().temperature(phaseIdx))
00175 { return fs_->temperature(phaseIdx); }
00176
00180 auto pressure(unsigned phaseIdx) const
00181 -> decltype(std::declval<FluidState>().pressure(phaseIdx))
00182 { return fs_->pressure(phaseIdx); }
00183
00187 auto enthalpy(unsigned phaseIdx) const
00188 -> decltype(std::declval<FluidState>().enthalpy(phaseIdx))
00189 { return fs_->enthalpy(phaseIdx); }
00190
00194 auto internalEnergy(unsigned phaseIdx) const
00195 -> decltype(std::declval<FluidState>().internalEnergy(phaseIdx))
00196 { return fs_->internalEnergy(phaseIdx); }
00197
00201 auto viscosity(unsigned phaseIdx) const
00202 -> decltype(std::declval<FluidState>().viscosity(phaseIdx))
00203 { return fs_->viscosity(phaseIdx); }
00204
00205
00206
00207
00208
00209
00210
00214 void setSaturation(unsigned phaseIdx, const Scalar& value)
00215 { saturation_[phaseIdx] = value; }
00216
00225 void checkDefined() const
00226 {
00227 Valgrind::CheckDefined(saturation_);
00228 }
00229
00230 protected:
00231 const FluidState* fs_;
00232 std::array<Scalar, numPhases> saturation_;
00233 };
00234
00235 }
00236
00237 #endif