My Project
Loading...
Searching...
No Matches
SinglePhaseFluidSystem.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_SINGLE_PHASE_FLUIDSYSTEM_HPP
28#define OPM_SINGLE_PHASE_FLUIDSYSTEM_HPP
29
30#include "BaseFluidSystem.hpp"
32
39
40#include <limits>
41#include <cassert>
42#include <string_view>
43
44namespace Opm {
45
55template <class Scalar, class Fluid>
57 : public BaseFluidSystem<Scalar, SinglePhaseFluidSystem<Scalar, Fluid> >
58{
61
62public:
64 template <class Evaluation>
65 struct ParameterCache : public NullParameterCache<Evaluation>
66 {};
67
68 /****************************************
69 * Fluid phase related static parameters
70 ****************************************/
71
73 static const int numPhases = 1;
74
76 static std::string_view phaseName([[maybe_unused]] unsigned phaseIdx)
77 {
78 assert(phaseIdx < numPhases);
79
80 return Fluid::name();
81 }
82
84 static bool isLiquid(unsigned /*phaseIdx*/)
85 {
86 //assert(0 <= phaseIdx && phaseIdx < numPhases);
87
88 return Fluid::isLiquid();
89 }
90
92 static bool isCompressible(unsigned /*phaseIdx*/)
93 {
94 //assert(0 <= phaseIdx && phaseIdx < numPhases);
95
96 // let the fluid decide
97 return Fluid::isCompressible();
98 }
99
101 static bool isIdealMixture(unsigned /*phaseIdx*/)
102 {
103 //assert(0 <= phaseIdx && phaseIdx < numPhases);
104
105 // we assume immisibility
106 return true;
107 }
108
110 static bool isIdealGas(unsigned /*phaseIdx*/)
111 {
112 //assert(0 <= phaseIdx && phaseIdx < numPhases);
113
114 // let the fluid decide
115 return Fluid::isIdealGas();
116 }
117
118 /****************************************
119 * Component related static parameters
120 ****************************************/
121
123 static const int numComponents = 1;
124
126 static std::string_view componentName([[maybe_unused]] unsigned compIdx)
127 {
128 assert(compIdx < numComponents);
129
130 return Fluid::name();
131 }
132
134 static Scalar molarMass(unsigned /*compIdx*/)
135 {
136 //assert(0 <= compIdx && compIdx < numComponents);
137
138 return Fluid::molarMass();
139 }
140
146 static Scalar criticalTemperature(unsigned /*compIdx*/)
147 {
148 //assert(0 <= compIdx && compIdx < numComponents);
149
150 return Fluid::criticalTemperature();
151 }
152
158 static Scalar criticalPressure(unsigned /*compIdx*/)
159 {
160 //assert(0 <= compIdx && compIdx < numComponents);
161
162 return Fluid::criticalPressure();
163 }
164
170 static Scalar acentricFactor(unsigned /*compIdx*/)
171 {
172 //assert(0 <= compIdx && compIdx < numComponents);
173
174 return Fluid::acentricFactor();
175 }
176
177 /****************************************
178 * thermodynamic relations
179 ****************************************/
180
182 static void init()
183 { }
184
186 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
187 static LhsEval density(const FluidState& fluidState,
188 const ParameterCache<ParamCacheEval>& /*paramCache*/,
189 unsigned phaseIdx)
190 {
191 assert(phaseIdx < numPhases);
192
193 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
194 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
195 return Fluid::density(T, p);
196 }
197
199 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
200 static LhsEval viscosity(const FluidState& fluidState,
201 const ParameterCache<ParamCacheEval>& /*paramCache*/,
202 unsigned phaseIdx)
203 {
204 assert(phaseIdx < numPhases);
205
206 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
207 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
208 return Fluid::viscosity(T, p);
209 }
210
212 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
213 static LhsEval fugacityCoefficient(const FluidState& /*fluidState*/,
214 const ParameterCache<ParamCacheEval>& /*paramCache*/,
215 unsigned phaseIdx,
216 unsigned compIdx)
217 {
218 assert(phaseIdx < numPhases);
219 assert(compIdx < numComponents);
220
221 if (phaseIdx == compIdx)
222 // TODO (?): calculate the real fugacity coefficient of
223 // the component in the fluid. Probably that's not worth
224 // the effort, since the fugacity coefficient of the other
225 // component is infinite anyway...
226 return 1.0;
227 return std::numeric_limits<Scalar>::infinity();
228 }
229
231 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
232 static LhsEval enthalpy(const FluidState& fluidState,
233 const ParameterCache<ParamCacheEval>& /*paramCache*/,
234 unsigned phaseIdx)
235 {
236 assert(phaseIdx < numPhases);
237
238 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
239 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
240 return Fluid::enthalpy(T, p);
241 }
242
244 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
245 static LhsEval thermalConductivity(const FluidState& fluidState,
246 const ParameterCache<ParamCacheEval>& /*paramCache*/,
247 unsigned phaseIdx)
248 {
249 assert(phaseIdx < numPhases);
250
251 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
252 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
253 return Fluid::thermalConductivity(T, p);
254 }
255
257 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
258 static LhsEval heatCapacity(const FluidState& fluidState,
259 const ParameterCache<ParamCacheEval>& /*paramCache*/,
260 unsigned phaseIdx)
261 {
262 assert(phaseIdx < numPhases);
263
264 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
265 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
266 return Fluid::heatCapacity(T, p);
267 }
268};
269
270} // namespace Opm
271
272#endif
The base class for all fluid systems.
Represents the gas phase of a single (pseudo-) component.
Material properties of pure water .
Represents the liquid phase of a single (pseudo-) component.
Properties of pure molecular nitrogen .
A parameter cache which does nothing.
A simple version of pure water.
A generic class which tabulates all thermodynamic properties of a given component.
The base class for all fluid systems.
Definition BaseFluidSystem.hpp:47
Scalar Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:52
A parameter cache which does nothing.
Definition NullParameterCache.hpp:40
A fluid system for single phase models.
Definition SinglePhaseFluidSystem.hpp:58
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition SinglePhaseFluidSystem.hpp:245
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition SinglePhaseFluidSystem.hpp:200
static const int numComponents
Number of chemical species in the fluid system.
Definition SinglePhaseFluidSystem.hpp:123
static Scalar criticalTemperature(unsigned)
Critical temperature of a component [K].
Definition SinglePhaseFluidSystem.hpp:146
static Scalar acentricFactor(unsigned)
The acentric factor of a component [].
Definition SinglePhaseFluidSystem.hpp:170
static const int numPhases
Number of fluid phases in the fluid system.
Definition SinglePhaseFluidSystem.hpp:73
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition SinglePhaseFluidSystem.hpp:76
static LhsEval fugacityCoefficient(const FluidState &, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition SinglePhaseFluidSystem.hpp:213
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition SinglePhaseFluidSystem.hpp:110
static void init()
Initialize the fluid system's static parameters.
Definition SinglePhaseFluidSystem.hpp:182
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition SinglePhaseFluidSystem.hpp:187
static Scalar molarMass(unsigned)
Return the molar mass of a component in [kg/mol].
Definition SinglePhaseFluidSystem.hpp:134
static Scalar criticalPressure(unsigned)
Critical pressure of a component [Pa].
Definition SinglePhaseFluidSystem.hpp:158
static LhsEval heatCapacity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Specific isobaric heat capacity of a fluid phase [J/kg].
Definition SinglePhaseFluidSystem.hpp:258
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition SinglePhaseFluidSystem.hpp:92
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition SinglePhaseFluidSystem.hpp:101
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition SinglePhaseFluidSystem.hpp:126
static LhsEval enthalpy(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy [J/kg...
Definition SinglePhaseFluidSystem.hpp:232
static bool isLiquid(unsigned)
Return whether a phase is liquid.
Definition SinglePhaseFluidSystem.hpp:84
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
The type of the fluid system's parameter cache.
Definition SinglePhaseFluidSystem.hpp:66