My Project
Loading...
Searching...
No Matches
TwoPhaseImmiscibleFluidSystem.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_TWO_PHASE_IMMISCIBLE_FLUID_SYSTEM_HPP
28#define OPM_TWO_PHASE_IMMISCIBLE_FLUID_SYSTEM_HPP
29
30#include <cassert>
31#include <limits>
32#include <string_view>
33
37
38#include "BaseFluidSystem.hpp"
40
41namespace Opm {
42
56template <class Scalar, class WettingPhase, class NonwettingPhase>
58 : public BaseFluidSystem<Scalar, TwoPhaseImmiscibleFluidSystem<Scalar, WettingPhase, NonwettingPhase> >
59{
60 // do not try to instanciate this class, it has only static members!
62 {}
63
66
67public:
68 template <class Evaluation>
69 struct ParameterCache : public NullParameterCache<Evaluation>
70 {};
71
72 /****************************************
73 * Fluid phase related static parameters
74 ****************************************/
75
77 static const int numPhases = 2;
78
80 static const int wettingPhaseIdx = 0;
82 static const int nonWettingPhaseIdx = 1;
83
85 static std::string_view phaseName(unsigned phaseIdx)
86 {
87 assert(phaseIdx < numPhases);
88
89 static const std::string_view name[] = {
90 "wetting",
91 "nonwetting"
92 };
93 return name[phaseIdx];
94 }
95
97 static bool isLiquid(unsigned phaseIdx)
98 {
99 //assert(0 <= phaseIdx && phaseIdx < numPhases);
100 return
101 (phaseIdx == wettingPhaseIdx)
102 ? WettingPhase::isLiquid()
103 : NonwettingPhase::isLiquid();
104 }
105
107 static bool isCompressible(unsigned phaseIdx)
108 {
109 //assert(0 <= phaseIdx && phaseIdx < numPhases);
110
111 return
112 (phaseIdx == wettingPhaseIdx)
113 ? WettingPhase::isCompressible()
114 : NonwettingPhase::isCompressible();
115 }
116
118 static bool isIdealGas(unsigned phaseIdx)
119 {
120 //assert(0 <= phaseIdx && phaseIdx < numPhases);
121
122 // let the fluids decide
123 return
124 (phaseIdx == wettingPhaseIdx)
125 ? WettingPhase::isIdealGas()
126 : NonwettingPhase::isIdealGas();
127 }
128
130 static bool isIdealMixture(unsigned /*phaseIdx*/)
131 {
132 //assert(0 <= phaseIdx && phaseIdx < numPhases);
133
134 // we assume immisibility
135 return true;
136 }
137
138 /****************************************
139 * Component related static parameters
140 ****************************************/
141
143 static const int numComponents = 2;
144
146 static const int wettingCompIdx = 0;
148 static const int nonWettingCompIdx = 1;
149
151 static std::string_view componentName(unsigned compIdx)
152 {
153 assert(compIdx < numComponents);
154
155 if (compIdx == wettingCompIdx)
156 return WettingPhase::name();
157 return NonwettingPhase::name();
158 }
159
161 static Scalar molarMass(unsigned compIdx)
162 {
163 //assert(0 <= compIdx && compIdx < numComponents);
164
165 // let the fluids decide
166 return
167 (compIdx == wettingCompIdx)
168 ? WettingPhase::molarMass()
169 : NonwettingPhase::molarMass();
170 }
171
175 static Scalar criticalTemperature(unsigned compIdx)
176 {
177 //assert(0 <= compIdx && compIdx < numComponents);
178 // let the fluids decide
179 return
180 (compIdx == wettingCompIdx)
181 ? WettingPhase::criticalTemperature()
182 : NonwettingPhase::criticalTemperature();
183 }
184
188 static Scalar criticalPressure(unsigned compIdx)
189 {
190 //assert(0 <= compIdx && compIdx < numComponents);
191 // let the fluids decide
192 return
193 (compIdx == wettingCompIdx)
194 ? WettingPhase::criticalPressure()
195 : NonwettingPhase::criticalPressure();
196 }
197
201 static Scalar acentricFactor(unsigned compIdx)
202 {
203 //assert(0 <= compIdx && compIdx < numComponents);
204 // let the fluids decide
205 return
206 (compIdx == wettingCompIdx)
207 ? WettingPhase::acentricFactor()
208 : NonwettingPhase::acentricFactor();
209 }
210
211 /****************************************
212 * thermodynamic relations
213 ****************************************/
214
216 static void init()
217 {
218 // two gaseous phases at once do not make sense physically!
219 // (But two liquids are fine)
220 assert(WettingPhase::isLiquid() || NonwettingPhase::isLiquid());
221 }
222
224 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
225 static LhsEval density(const FluidState& fluidState,
226 const ParameterCache<ParamCacheEval>& /*paramCache*/,
227 unsigned phaseIdx)
228 {
229 assert(phaseIdx < numPhases);
230
231 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
232 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
233 if (phaseIdx == wettingPhaseIdx)
234 return WettingPhase::density(temperature, pressure);
235 return NonwettingPhase::density(temperature, pressure);
236 }
237
239 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
240 static LhsEval viscosity(const FluidState& fluidState,
241 const ParameterCache<ParamCacheEval>& /*paramCache*/,
242 unsigned phaseIdx)
243 {
244 assert(phaseIdx < numPhases);
245
246 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
247 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
248 if (phaseIdx == wettingPhaseIdx)
249 return WettingPhase::viscosity(temperature, pressure);
250 return NonwettingPhase::viscosity(temperature, pressure);
251 }
252
254 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
255 static LhsEval fugacityCoefficient(const FluidState& /*fluidState*/,
256 const ParameterCache<ParamCacheEval>& /*paramCache*/,
257 unsigned phaseIdx,
258 unsigned compIdx)
259 {
260 assert(phaseIdx < numPhases);
261 assert(compIdx < numComponents);
262
263 if (phaseIdx == compIdx)
264 // TODO (?): calculate the real fugacity coefficient of
265 // the component in the fluid. Probably that's not worth
266 // the effort, since the fugacity coefficient of the other
267 // component is infinite anyway...
268 return 1.0;
269 return std::numeric_limits<Scalar>::infinity();
270 }
271
273 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
274 static LhsEval enthalpy(const FluidState& fluidState,
275 const ParameterCache<ParamCacheEval>& /*paramCache*/,
276 unsigned phaseIdx)
277 {
278 assert(phaseIdx < numPhases);
279
280 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
281 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
282 if (phaseIdx == wettingPhaseIdx)
283 return WettingPhase::enthalpy(temperature, pressure);
284 return NonwettingPhase::enthalpy(temperature, pressure);
285 }
286
288 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
289 static LhsEval thermalConductivity(const FluidState& fluidState,
290 const ParameterCache<ParamCacheEval>& /*paramCache*/,
291 unsigned phaseIdx)
292 {
293 assert(phaseIdx < numPhases);
294
295 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
296 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
297 if (phaseIdx == wettingPhaseIdx)
298 return WettingPhase::thermalConductivity(temperature, pressure);
299 return NonwettingPhase::thermalConductivity(temperature, pressure);
300 }
301
303 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
304 static LhsEval heatCapacity(const FluidState& fluidState,
305 const ParameterCache<ParamCacheEval>& /*paramCache*/,
306 unsigned phaseIdx)
307 {
308 assert(phaseIdx < numPhases);
309
310 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
311 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
312 if (phaseIdx == wettingPhaseIdx)
313 return WettingPhase::heatCapacity(temperature, pressure);
314 return NonwettingPhase::heatCapacity(temperature, pressure);
315 }
316};
317
318} // namespace Opm
319
320#endif
The base class for all fluid systems.
Represents the gas phase of a single (pseudo-) component.
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system assum...
Represents the liquid phase of a single (pseudo-) component.
A parameter cache which does nothing.
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 two-phase models assuming immiscibility and thermodynamic equilibrium.
Definition TwoPhaseImmiscibleFluidSystem.hpp:59
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition TwoPhaseImmiscibleFluidSystem.hpp:130
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition TwoPhaseImmiscibleFluidSystem.hpp:201
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition TwoPhaseImmiscibleFluidSystem.hpp:289
static void init()
Initialize the fluid system's static parameters.
Definition TwoPhaseImmiscibleFluidSystem.hpp:216
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition TwoPhaseImmiscibleFluidSystem.hpp:97
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 TwoPhaseImmiscibleFluidSystem.hpp:255
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition TwoPhaseImmiscibleFluidSystem.hpp:161
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition TwoPhaseImmiscibleFluidSystem.hpp:151
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition TwoPhaseImmiscibleFluidSystem.hpp:118
static const int numPhases
Number of fluid phases in the fluid system.
Definition TwoPhaseImmiscibleFluidSystem.hpp:77
static const int wettingCompIdx
Index of the wetting phase's component.
Definition TwoPhaseImmiscibleFluidSystem.hpp:146
static const int nonWettingCompIdx
Index of the non-wetting phase's component.
Definition TwoPhaseImmiscibleFluidSystem.hpp:148
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 TwoPhaseImmiscibleFluidSystem.hpp:274
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition TwoPhaseImmiscibleFluidSystem.hpp:240
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition TwoPhaseImmiscibleFluidSystem.hpp:85
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition TwoPhaseImmiscibleFluidSystem.hpp:175
static LhsEval heatCapacity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Specific isobaric heat capacity of a fluid phase [J/kg].
Definition TwoPhaseImmiscibleFluidSystem.hpp:304
static const int nonWettingPhaseIdx
Index of the non-wetting phase.
Definition TwoPhaseImmiscibleFluidSystem.hpp:82
static const int numComponents
Number of chemical species in the fluid system.
Definition TwoPhaseImmiscibleFluidSystem.hpp:143
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition TwoPhaseImmiscibleFluidSystem.hpp:225
static const int wettingPhaseIdx
Index of the wetting phase.
Definition TwoPhaseImmiscibleFluidSystem.hpp:80
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition TwoPhaseImmiscibleFluidSystem.hpp:188
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition TwoPhaseImmiscibleFluidSystem.hpp:107
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition TwoPhaseImmiscibleFluidSystem.hpp:70