TemperatureOverlayFluidState.hpp
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
27 #ifndef OPM_TEMPERATURE_OVERLAY_FLUID_STATE_HPP
28 #define OPM_TEMPERATURE_OVERLAY_FLUID_STATE_HPP
29 
30 #include <opm/common/Valgrind.hpp>
31 
32 #include <utility>
33 
34 namespace Opm {
35 
41 template <class FluidState>
43 {
44 public:
45  typedef typename FluidState::Scalar Scalar;
46 
47  enum { numPhases = FluidState::numPhases };
48  enum { numComponents = FluidState::numComponents };
49 
58  TemperatureOverlayFluidState(const FluidState& fs)
59  : fs_(&fs)
60  {
61  temperature_ = fs.temperature(/*phaseIdx=*/0);
62  }
63 
64  TemperatureOverlayFluidState(Scalar T, const FluidState& fs)
65  : temperature_(T), fs_(&fs)
66  { }
67 
68  // copy constructor
70  : fs_(fs.fs_)
71  , temperature_(fs.temperature_)
72  {}
73 
74  // assignment operator
76  {
77  fs_ = fs.fs_;
78  temperature_ = fs.temperature_;
79  return *this;
80  }
81 
82  /*****************************************************
83  * Generic access to fluid properties (No assumptions
84  * on thermodynamic equilibrium required)
85  *****************************************************/
89  auto saturation(unsigned phaseIdx) const
90  -> decltype(std::declval<FluidState>().saturation(phaseIdx))
91  { return fs_->saturation(phaseIdx); }
92 
96  bool phaseIsPresent(unsigned phaseIdx) const
97  { return fs_->phaseIsPresent(phaseIdx); }
98 
102  auto moleFraction(unsigned phaseIdx, unsigned compIdx) const
103  -> decltype(std::declval<FluidState>().moleFraction(phaseIdx, compIdx))
104  { return fs_->moleFraction(phaseIdx, compIdx); }
105 
109  auto massFraction(unsigned phaseIdx, unsigned compIdx) const
110  -> decltype(std::declval<FluidState>().massFraction(phaseIdx, compIdx))
111  { return fs_->massFraction(phaseIdx, compIdx); }
112 
121  auto averageMolarMass(unsigned phaseIdx) const
122  -> decltype(std::declval<FluidState>().averageMolarMass(phaseIdx))
123  { return fs_->averageMolarMass(phaseIdx); }
124 
134  auto molarity(unsigned phaseIdx, unsigned compIdx) const
135  -> decltype(std::declval<FluidState>().molarity(phaseIdx, compIdx))
136  { return fs_->molarity(phaseIdx, compIdx); }
137 
141  auto fugacity(unsigned phaseIdx, unsigned compIdx) const
142  -> decltype(std::declval<FluidState>().fugacity(phaseIdx, compIdx))
143  { return fs_->fugacity(phaseIdx, compIdx); }
144 
148  auto fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
149  -> decltype(std::declval<FluidState>().fugacityCoefficient(phaseIdx, compIdx))
150  { return fs_->fugacityCoefficient(phaseIdx, compIdx); }
151 
155  auto molarVolume(unsigned phaseIdx) const
156  -> decltype(std::declval<FluidState>().molarVolume(phaseIdx))
157  { return fs_->molarVolume(phaseIdx); }
158 
162  auto density(unsigned phaseIdx) const
163  -> decltype(std::declval<FluidState>().density(phaseIdx))
164  { return fs_->density(phaseIdx); }
165 
169  auto molarDensity(unsigned phaseIdx) const
170  -> decltype(std::declval<FluidState>().molarDensity(phaseIdx))
171  { return fs_->molarDensity(phaseIdx); }
172 
176  const Scalar& temperature(unsigned /*phaseIdx*/) const
177  { return temperature_; }
178 
182  auto pressure(unsigned phaseIdx) const
183  -> decltype(std::declval<FluidState>().pressure(phaseIdx))
184  { return fs_->pressure(phaseIdx); }
185 
189  auto enthalpy(unsigned phaseIdx) const
190  -> decltype(std::declval<FluidState>().enthalpy(phaseIdx))
191  { return fs_->enthalpy(phaseIdx); }
192 
196  auto internalEnergy(unsigned phaseIdx) const
197  -> decltype(std::declval<FluidState>().internalEnergy(phaseIdx))
198  { return fs_->internalEnergy(phaseIdx); }
199 
203  auto viscosity(unsigned phaseIdx) const
204  -> decltype(std::declval<FluidState>().viscosity(phaseIdx))
205  { return fs_->viscosity(phaseIdx); }
206 
207 
208  /*****************************************************
209  * Setter methods. Note that these are not part of the
210  * generic FluidState interface but specific for each
211  * implementation...
212  *****************************************************/
216  void setTemperature(const Scalar& value)
217  { temperature_ = value; }
218 
227  void checkDefined() const
228  {
229  Valgrind::CheckDefined(temperature_);
230  }
231 
232 protected:
233  const FluidState* fs_;
234  Scalar temperature_;
235 };
236 
237 } // namespace Opm
238 
239 #endif
auto fugacity(unsigned phaseIdx, unsigned compIdx) const -> decltype(std::declval< FluidState >().fugacity(phaseIdx, compIdx))
The fugacity of a component in a phase [Pa].
Definition: TemperatureOverlayFluidState.hpp:141
auto density(unsigned phaseIdx) const -> decltype(std::declval< FluidState >().density(phaseIdx))
The mass density of a fluid phase [kg/m^3].
Definition: TemperatureOverlayFluidState.hpp:162
void checkDefined() const
Make sure that all attributes are defined.
Definition: TemperatureOverlayFluidState.hpp:227
auto molarDensity(unsigned phaseIdx) const -> decltype(std::declval< FluidState >().molarDensity(phaseIdx))
The molar density of a fluid phase [mol/m^3].
Definition: TemperatureOverlayFluidState.hpp:169
This is a fluid state which allows to set the fluid temperatures and takes all other quantities from ...
Definition: TemperatureOverlayFluidState.hpp:42
TemperatureOverlayFluidState(const FluidState &fs)
Constructor.
Definition: TemperatureOverlayFluidState.hpp:58
auto internalEnergy(unsigned phaseIdx) const -> decltype(std::declval< FluidState >().internalEnergy(phaseIdx))
The specific internal energy of a fluid phase [J/kg].
Definition: TemperatureOverlayFluidState.hpp:196
Definition: Air_Mesitylene.hpp:33
bool phaseIsPresent(unsigned phaseIdx) const
Returns true iff a fluid phase shall be assumed to be present.
Definition: TemperatureOverlayFluidState.hpp:96
auto viscosity(unsigned phaseIdx) const -> decltype(std::declval< FluidState >().viscosity(phaseIdx))
The dynamic viscosity of a fluid phase [Pa s].
Definition: TemperatureOverlayFluidState.hpp:203
void setTemperature(const Scalar &value)
Set the temperature [K] of a fluid phase.
Definition: TemperatureOverlayFluidState.hpp:216
auto fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const -> decltype(std::declval< FluidState >().fugacityCoefficient(phaseIdx, compIdx))
The fugacity coefficient of a component in a phase [].
Definition: TemperatureOverlayFluidState.hpp:148
auto molarVolume(unsigned phaseIdx) const -> decltype(std::declval< FluidState >().molarVolume(phaseIdx))
The molar volume of a fluid phase [m^3/mol].
Definition: TemperatureOverlayFluidState.hpp:155
auto enthalpy(unsigned phaseIdx) const -> decltype(std::declval< FluidState >().enthalpy(phaseIdx))
The specific enthalpy of a fluid phase [J/kg].
Definition: TemperatureOverlayFluidState.hpp:189
auto averageMolarMass(unsigned phaseIdx) const -> decltype(std::declval< FluidState >().averageMolarMass(phaseIdx))
The average molar mass of a fluid phase [kg/mol].
Definition: TemperatureOverlayFluidState.hpp:121
auto moleFraction(unsigned phaseIdx, unsigned compIdx) const -> decltype(std::declval< FluidState >().moleFraction(phaseIdx, compIdx))
The mole fraction of a component in a phase [].
Definition: TemperatureOverlayFluidState.hpp:102
auto massFraction(unsigned phaseIdx, unsigned compIdx) const -> decltype(std::declval< FluidState >().massFraction(phaseIdx, compIdx))
The mass fraction of a component in a phase [].
Definition: TemperatureOverlayFluidState.hpp:109
auto saturation(unsigned phaseIdx) const -> decltype(std::declval< FluidState >().saturation(phaseIdx))
Returns the saturation of a phase [].
Definition: TemperatureOverlayFluidState.hpp:89
auto molarity(unsigned phaseIdx, unsigned compIdx) const -> decltype(std::declval< FluidState >().molarity(phaseIdx, compIdx))
The molar concentration of a component in a phase [mol/m^3].
Definition: TemperatureOverlayFluidState.hpp:134
auto pressure(unsigned phaseIdx) const -> decltype(std::declval< FluidState >().pressure(phaseIdx))
The pressure of a fluid phase [Pa].
Definition: TemperatureOverlayFluidState.hpp:182
const Scalar & temperature(unsigned) const
The temperature of a fluid phase [K].
Definition: TemperatureOverlayFluidState.hpp:176