My Project
BaseFluidSystem.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_BASE_FLUID_SYSTEM_HPP
28#define OPM_BASE_FLUID_SYSTEM_HPP
29
30#include <opm/common/utility/Demangle.hpp>
31
32#include <stdexcept>
33
35
36#include <typeinfo>
37
38namespace Opm {
39
44template <class ScalarT, class Implementation>
46{
47public:
51 typedef ScalarT Scalar;
52
60 template <class Evaluation>
62 ParameterCache() = delete; // derived fluid systems must specify this class!
63 };
64
66 static const int numComponents = -1000;
67
69 static const int numPhases = -2000;
70
76 static char* phaseName(unsigned /*phaseIdx*/)
77 {
78 throw std::runtime_error(not_implemented("phaseName"));
79 }
80
86 static bool isLiquid(unsigned /*phaseIdx*/)
87 {
88 throw std::runtime_error(not_implemented("isLiquid"));
89 }
90
105 static bool isIdealMixture(unsigned /*phaseIdx*/)
106 {
107 throw std::runtime_error(not_implemented("isIdealMixture"));
108 }
109
119 static bool isCompressible(unsigned /*phaseIdx*/)
120 {
121 throw std::runtime_error(not_implemented("isCompressible"));
122 }
123
130 static bool isIdealGas(unsigned /*phaseIdx*/)
131 {
132 throw std::runtime_error(not_implemented("isIdealGas"));
133 }
134
140 static const char* componentName(unsigned /*compIdx*/)
141 {
142 throw std::runtime_error(not_implemented("componentName"));
143 }
144
150 static Scalar molarMass(unsigned /*compIdx*/)
151 {
152 throw std::runtime_error(not_implemented("molarMass"));
153 }
154
160 static Scalar acentricFactor(unsigned /*compIdx*/)
161 {
162 throw std::runtime_error(not_implemented("acentricFactor"));
163 }
164
168 static void init()
169 { }
170
177 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCache>
178 static LhsEval density(const FluidState& /*fluidState*/,
179 const ParamCache& /*paramCache*/,
180 unsigned /*phaseIdx*/)
181 {
182 throw std::runtime_error(not_implemented("density"));
183 }
184
199 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCache>
200 static LhsEval fugacityCoefficient(const FluidState& /*fluidState*/,
201 ParamCache& /*paramCache*/,
202 unsigned /*phaseIdx*/,
203 unsigned /*compIdx*/)
204 {
205 throw std::runtime_error(not_implemented("fugacityCoefficient"));
206 }
207
214 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCache>
215 static LhsEval viscosity(const FluidState& /*fluidState*/,
216 ParamCache& /*paramCache*/,
217 unsigned /*phaseIdx*/)
218 {
219 throw std::runtime_error(not_implemented("viscosity"));
220 }
221
239 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCache>
240 static LhsEval diffusionCoefficient(const FluidState& /*fluidState*/,
241 ParamCache& /*paramCache*/,
242 unsigned /*phaseIdx*/,
243 unsigned /*compIdx*/)
244 {
245 throw std::runtime_error(not_implemented("diffusionCoefficient"));
246 }
247
255 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCache>
256 static LhsEval enthalpy(const FluidState& /*fluidState*/,
257 ParamCache& /*paramCache*/,
258 unsigned /*phaseIdx*/)
259 {
260 throw std::runtime_error(not_implemented("enthalpy"));
261 }
262
269 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCache>
270 static LhsEval thermalConductivity(const FluidState& /*fluidState*/,
271 ParamCache& /*paramCache*/,
272 unsigned /*phaseIdx*/)
273 {
274 throw std::runtime_error(not_implemented("thermalConductivity"));
275 }
276
283 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCache>
284 static LhsEval heatCapacity(const FluidState& /*fluidState*/,
285 ParamCache& /*paramCache*/,
286 unsigned /*phaseIdx*/)
287 {
288 throw std::runtime_error(not_implemented("heatCapacity"));
289 }
290
291
293 static unsigned phaseIsActive(unsigned /*phaseIdx*/)
294 {
295 return true;
296 }
297
298private:
299 static std::string not_implemented(const char* method)
300 {
301 return "Not implemented: The fluid system '" +
302 demangle(typeid(Implementation).name()) +
303 "' does not provide a " + method + "() method!";
304 }
305};
306
307} // namespace Opm
308
309#endif
The base class for all fluid systems.
Definition: BaseFluidSystem.hpp:46
static char * phaseName(unsigned)
Return the human readable name of a fluid phase.
Definition: BaseFluidSystem.hpp:76
static LhsEval heatCapacity(const FluidState &, ParamCache &, unsigned)
Specific isobaric heat capacity of a fluid phase [J/kg].
Definition: BaseFluidSystem.hpp:284
static unsigned phaseIsActive(unsigned)
Returns whether a fluid phase is active.
Definition: BaseFluidSystem.hpp:293
static LhsEval enthalpy(const FluidState &, ParamCache &, unsigned)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy [J/kg...
Definition: BaseFluidSystem.hpp:256
static const int numPhases
Number of fluid phases in the fluid system.
Definition: BaseFluidSystem.hpp:69
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: BaseFluidSystem.hpp:130
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: BaseFluidSystem.hpp:119
static Scalar molarMass(unsigned)
Return the molar mass of a component in [kg/mol].
Definition: BaseFluidSystem.hpp:150
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: BaseFluidSystem.hpp:105
static const char * componentName(unsigned)
Return the human readable name of a component.
Definition: BaseFluidSystem.hpp:140
static LhsEval diffusionCoefficient(const FluidState &, ParamCache &, unsigned, unsigned)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition: BaseFluidSystem.hpp:240
static bool isLiquid(unsigned)
Return whether a phase is liquid.
Definition: BaseFluidSystem.hpp:86
static void init()
Initialize the fluid system's static parameters.
Definition: BaseFluidSystem.hpp:168
ScalarT Scalar
The type used for scalar quantities.
Definition: BaseFluidSystem.hpp:51
static const int numComponents
Number of chemical species in the fluid system.
Definition: BaseFluidSystem.hpp:66
static LhsEval density(const FluidState &, const ParamCache &, unsigned)
Calculate the density [kg/m^3] of a fluid phase.
Definition: BaseFluidSystem.hpp:178
static Scalar acentricFactor(unsigned)
Return the acetntric factor of a component.
Definition: BaseFluidSystem.hpp:160
static LhsEval fugacityCoefficient(const FluidState &, ParamCache &, unsigned, unsigned)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition: BaseFluidSystem.hpp:200
static LhsEval viscosity(const FluidState &, ParamCache &, unsigned)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: BaseFluidSystem.hpp:215
static LhsEval thermalConductivity(const FluidState &, ParamCache &, unsigned)
Thermal conductivity of a fluid phase [W/(m K)].
Definition: BaseFluidSystem.hpp:270
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:30
std::string demangle(const char *mangled_symbol)
Returns demangled name of symbol.
The type of the fluid system's parameter cache.
Definition: BaseFluidSystem.hpp:61