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_PRESSURE_OVERLAY_FLUID_STATE_HPP
00028 #define OPM_PRESSURE_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 PressureOverlayFluidState
00044 {
00045 public:
00046 typedef typename FluidState::Scalar Scalar;
00047
00048 enum { numPhases = FluidState::numPhases };
00049 enum { numComponents = FluidState::numComponents };
00050
00058 PressureOverlayFluidState(const FluidState& fs)
00059 : fs_(&fs)
00060 {
00061 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
00062 pressure_[phaseIdx] = fs.pressure(phaseIdx);
00063 }
00064
00065
00066 PressureOverlayFluidState(const PressureOverlayFluidState& fs)
00067 : fs_(fs.fs_)
00068 , pressure_(fs.pressure_)
00069 {
00070 }
00071
00072
00073 PressureOverlayFluidState& operator=(const PressureOverlayFluidState& fs)
00074 {
00075 fs_ = fs.fs_;
00076 pressure_ = fs.pressure_;
00077 return *this;
00078 }
00079
00080
00081
00082
00083
00087 auto saturation(unsigned phaseIdx) const
00088 -> decltype(std::declval<FluidState>().saturation(phaseIdx))
00089 { return fs_->saturation(phaseIdx); }
00090
00094 bool phaseIsPresent(unsigned phaseIdx) const
00095 { return fs_->phaseIsPresent(phaseIdx); }
00096
00100 auto moleFraction(unsigned phaseIdx, unsigned compIdx) const
00101 -> decltype(std::declval<FluidState>().moleFraction(phaseIdx, compIdx))
00102 { return fs_->moleFraction(phaseIdx, compIdx); }
00103
00107 auto massFraction(unsigned phaseIdx, unsigned compIdx) const
00108 -> decltype(std::declval<FluidState>().massFraction(phaseIdx, compIdx))
00109 { return fs_->massFraction(phaseIdx, compIdx); }
00110
00119 auto averageMolarMass(unsigned phaseIdx) const
00120 -> decltype(std::declval<FluidState>().averageMolarMass(phaseIdx))
00121 { return fs_->averageMolarMass(phaseIdx); }
00122
00132 auto molarity(unsigned phaseIdx, unsigned compIdx) const
00133 -> decltype(std::declval<FluidState>().molarity(phaseIdx, compIdx))
00134 { return fs_->molarity(phaseIdx, compIdx); }
00135
00139 auto fugacity(unsigned phaseIdx, unsigned compIdx) const
00140 -> decltype(std::declval<FluidState>().fugacity(phaseIdx, compIdx))
00141 { return fs_->fugacity(phaseIdx, compIdx); }
00142
00146 auto fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
00147 -> decltype(std::declval<FluidState>().fugacityCoefficient(phaseIdx, compIdx))
00148 { return fs_->fugacityCoefficient(phaseIdx, compIdx); }
00149
00153 auto molarVolume(unsigned phaseIdx) const
00154 -> decltype(std::declval<FluidState>().molarVolume(phaseIdx))
00155 { return fs_->molarVolume(phaseIdx); }
00156
00160 auto density(unsigned phaseIdx) const
00161 -> decltype(std::declval<FluidState>().density(phaseIdx))
00162 { return fs_->density(phaseIdx); }
00163
00167 auto molarDensity(unsigned phaseIdx) const
00168 -> decltype(std::declval<FluidState>().molarDensity(phaseIdx))
00169 { return fs_->molarDensity(phaseIdx); }
00170
00174 auto temperature(unsigned phaseIdx) const
00175 -> decltype(std::declval<FluidState>().temperature(phaseIdx))
00176 { return fs_->temperature(phaseIdx); }
00177
00181 auto pressure(unsigned phaseIdx) const
00182 -> decltype(std::declval<FluidState>().pressure(phaseIdx))
00183 { return pressure_[phaseIdx]; }
00184
00188 auto enthalpy(unsigned phaseIdx) const
00189 -> decltype(std::declval<FluidState>().enthalpy(phaseIdx))
00190 { return fs_->enthalpy(phaseIdx); }
00191
00195 auto internalEnergy(unsigned phaseIdx) const
00196 -> decltype(std::declval<FluidState>().internalEnergy(phaseIdx))
00197 { return fs_->internalEnergy(phaseIdx); }
00198
00202 auto viscosity(unsigned phaseIdx) const
00203 -> decltype(std::declval<FluidState>().viscosity(phaseIdx))
00204 { return fs_->viscosity(phaseIdx); }
00205
00206
00207
00208
00209
00210
00211
00215 void setPressure(unsigned phaseIdx, const Scalar& value)
00216 { pressure_[phaseIdx] = value; }
00217
00226 void checkDefined() const
00227 {
00228 Valgrind::CheckDefined(pressure_);
00229 }
00230
00231 protected:
00232 const FluidState* fs_;
00233 std::array<Scalar, numPhases> pressure_;
00234 };
00235
00236 }
00237
00238 #endif