3.6-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
fluidsystems/base.hh
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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program 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 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program 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 this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_BASE_FLUID_SYSTEM_HH
25#define DUMUX_BASE_FLUID_SYSTEM_HH
26
27#include <string>
28
29#include <dune/common/exceptions.hh>
31#include "nullparametercache.hh"
32
33namespace Dumux {
34namespace FluidSystems {
35
43template <class ScalarType, class Implementation>
44class Base
45{
46public:
48 using Scalar = ScalarType;
49
52
56 // \{
57
59 static constexpr bool isTracerFluidSystem()
60 { return false; }
61
71 template<class I = Implementation, std::enable_if_t<!I::isTracerFluidSystem(), int> = 0>
72 static constexpr int getMainComponent(int phaseIdx)
73 { return phaseIdx; }
74
84 template<class T = Implementation>
85 static constexpr bool isCompressible(int phaseIdx)
86 {
87 static_assert(AlwaysFalse<T>::value, "Mandatory function not implemented: isCompressible(phaseIdx)");
88 return true;
89 }
90
94 template<class T = Implementation>
95 static constexpr bool isMiscible()
96 {
97 static_assert(AlwaysFalse<T>::value, "Mandatory function not implemented: isMiscible()");
98 return true;
99 }
100
107 static constexpr bool viscosityIsConstant(int phaseIdx)
108 { return false; }
109
115 static std::string phaseName(int phaseIdx)
116 { return "DefaultPhaseName"; }
117
123 static std::string componentName(int phaseIdx)
124 { return "DefaultComponentName"; }
125
126 // \}
127
133 template <class FluidState>
134 static Scalar density(const FluidState &fluidState,
135 int phaseIdx)
136 {
137 DUNE_THROW(Dune::NotImplemented, "density() method not implemented by the fluid system.");
138 }
139
146 template <class FluidState>
147 static Scalar density(const FluidState &fluidState,
148 const ParameterCache &paramCache,
149 int phaseIdx)
150 {
151 return Implementation::density(fluidState, phaseIdx);
152 }
153
159 template <class FluidState>
160 static Scalar molarDensity(const FluidState &fluidState,
161 int phaseIdx)
162 {
163 DUNE_THROW(Dune::NotImplemented, "molarDensity() method not implemented by the fluid system.");
164 }
165
172 template <class FluidState>
173 static Scalar molarDensity(const FluidState &fluidState,
174 const ParameterCache &paramCache,
175 int phaseIdx)
176 {
177 return Implementation::molarDensity(fluidState, phaseIdx);
178 }
179
196 template <class FluidState>
197 static Scalar fugacityCoefficient(const FluidState &fluidState,
198 int phaseIdx,
199 int compIdx)
200 {
201 DUNE_THROW(Dune::NotImplemented, "fugacityCoefficient() method not implemented by the fluid system.");
202 }
203
221 template <class FluidState>
222 static Scalar fugacityCoefficient(const FluidState &fluidState,
223 const ParameterCache &paramCache,
224 int phaseIdx,
225 int compIdx)
226 {
227 return Implementation::fugacityCoefficient(fluidState, phaseIdx, compIdx);
228 }
229
235 template <class FluidState>
236 static Scalar viscosity(const FluidState &fluidState,
237 int phaseIdx)
238 {
239 DUNE_THROW(Dune::NotImplemented, "viscosity() method not implemented by the fluid system.");
240 }
241
248 template <class FluidState>
249 static Scalar viscosity(const FluidState &fluidState,
250 const ParameterCache &paramCache,
251 int phaseIdx)
252 {
253 return Implementation::viscosity(fluidState, phaseIdx);
254 }
255
277 template <class FluidState>
278 static Scalar diffusionCoefficient(const FluidState &fluidState,
279 int phaseIdx,
280 int compIdx)
281 {
282 DUNE_THROW(Dune::NotImplemented, "diffusionCoefficient() method not implemented by the fluid system.");
283 }
284
307 template <class FluidState>
308 static Scalar diffusionCoefficient(const FluidState &fluidState,
309 const ParameterCache &paramCache,
310 int phaseIdx,
311 int compIdx)
312 {
313 return Implementation::diffusionCoefficient(fluidState, phaseIdx, compIdx);
314 }
315
325 template <class FluidState>
326 static Scalar binaryDiffusionCoefficient(const FluidState &fluidState,
327 int phaseIdx,
328 int compIIdx,
329 int compJIdx)
330
331 {
332 DUNE_THROW(Dune::NotImplemented, "binaryDiffusionCoefficient() method not implemented by the fluid system.");
333 }
334
345 template <class FluidState>
346 static Scalar binaryDiffusionCoefficient(const FluidState &fluidState,
347 const ParameterCache &paramCache,
348 int phaseIdx,
349 int compIIdx,
350 int compJIdx)
351
352 {
353 return Implementation::binaryDiffusionCoefficient(fluidState, phaseIdx, compIIdx, compJIdx);
355
362 template <class FluidState>
363 static Scalar enthalpy(const FluidState &fluidState,
364 int phaseIdx)
365 {
366 DUNE_THROW(Dune::NotImplemented, "enthalpy() method not implemented by the fluid system.");
367 }
368
376 template <class FluidState>
377 static Scalar enthalpy(const FluidState &fluidState,
378 const ParameterCache &paramCache,
379 int phaseIdx)
380 {
381 return Implementation::enthalpy(fluidState, phaseIdx);
382 }
383
389 template <class FluidState>
390 static Scalar thermalConductivity(const FluidState &fluidState,
391 int phaseIdx)
392 {
393 DUNE_THROW(Dune::NotImplemented, "thermalConductivity() method not implemented by the fluid system.");
394 }
395
402 template <class FluidState>
403 static Scalar thermalConductivity(const FluidState &fluidState,
404 const ParameterCache &paramCache,
405 int phaseIdx)
406 {
407 return Implementation::thermalConductivity(fluidState, phaseIdx);
408 }
409
423 template <class FluidState>
424 static Scalar heatCapacity(const FluidState &fluidState,
425 int phaseIdx)
426 {
427 DUNE_THROW(Dune::NotImplemented, "heatCapacity() method not implemented by the fluid system.");
428 }
429
444 template <class FluidState>
445 static Scalar heatCapacity(const FluidState &fluidState,
446 const ParameterCache &paramCache,
447 int phaseIdx)
448 {
449 return Implementation::heatCapacity(fluidState, phaseIdx);
451};
452
453} // end namespace FluidSystems
454
455} // end namespace Dumux
456
457#endif
Type traits.
The a parameter cache which does nothing.
Adaption of the non-isothermal two-phase two-component flow model to problems with CO2.
Definition: adapt.hh:29
std::string viscosity(int phaseIdx) noexcept
I/O name of viscosity for multiphase systems.
Definition: name.hh:74
std::string molarDensity(int phaseIdx) noexcept
I/O name of molar density for multiphase systems.
Definition: name.hh:83
std::string density(int phaseIdx) noexcept
I/O name of density for multiphase systems.
Definition: name.hh:65
Template which always yields a false value.
Definition: typetraits.hh:35
Fluid system base class.
Definition: fluidsystems/base.hh:45
static constexpr bool viscosityIsConstant(int phaseIdx)
Returns true if and only if a fluid phase is assumed to have a constant viscosity.
Definition: fluidsystems/base.hh:107
ScalarType Scalar
export the scalar type
Definition: fluidsystems/base.hh:48
static constexpr bool isTracerFluidSystem()
Some properties of the fluid system.
Definition: fluidsystems/base.hh:59
static Scalar binaryDiffusionCoefficient(const FluidState &fluidState, const ParameterCache &paramCache, int phaseIdx, int compIIdx, int compJIdx)
Given a phase's composition, temperature and pressure, return the binary diffusion coefficient for c...
Definition: fluidsystems/base.hh:346
static Scalar density(const FluidState &fluidState, const ParameterCache &paramCache, int phaseIdx)
Calculate the density of a fluid phase.
Definition: fluidsystems/base.hh:147
static Scalar fugacityCoefficient(const FluidState &fluidState, const ParameterCache &paramCache, int phaseIdx, int compIdx)
Calculate the fugacity coefficient of an individual component in a fluid phase.
Definition: fluidsystems/base.hh:222
static constexpr bool isMiscible()
Returns whether the fluids are miscible.
Definition: fluidsystems/base.hh:95
static Scalar density(const FluidState &fluidState, int phaseIdx)
Calculate the density of a fluid phase.
Definition: fluidsystems/base.hh:134
static Scalar heatCapacity(const FluidState &fluidState, const ParameterCache &paramCache, int phaseIdx)
Specific isobaric heat capacity of a fluid phase .
Definition: fluidsystems/base.hh:445
static Scalar thermalConductivity(const FluidState &fluidState, int phaseIdx)
Thermal conductivity of a fluid phase .
Definition: fluidsystems/base.hh:390
static std::string phaseName(int phaseIdx)
Return the human readable name of a fluid phase.
Definition: fluidsystems/base.hh:115
static Scalar molarDensity(const FluidState &fluidState, const ParameterCache &paramCache, int phaseIdx)
Calculate the molar density of a fluid phase.
Definition: fluidsystems/base.hh:173
static constexpr bool isCompressible(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: fluidsystems/base.hh:85
static Scalar fugacityCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Calculate the fugacity coefficient of an individual component in a fluid phase.
Definition: fluidsystems/base.hh:197
static Scalar diffusionCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase .
Definition: fluidsystems/base.hh:278
static Scalar binaryDiffusionCoefficient(const FluidState &fluidState, int phaseIdx, int compIIdx, int compJIdx)
Given a phase's composition, temperature and pressure, return the binary diffusion coefficient for c...
Definition: fluidsystems/base.hh:326
static std::string componentName(int phaseIdx)
Return the human readable name of a fluid phase.
Definition: fluidsystems/base.hh:123
static Scalar enthalpy(const FluidState &fluidState, int phaseIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy .
Definition: fluidsystems/base.hh:363
static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
Calculate the molar density of a fluid phase.
Definition: fluidsystems/base.hh:160
static constexpr int getMainComponent(int phaseIdx)
Get the main component of a given phase if possible.
Definition: fluidsystems/base.hh:72
static Scalar enthalpy(const FluidState &fluidState, const ParameterCache &paramCache, int phaseIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy .
Definition: fluidsystems/base.hh:377
static Scalar viscosity(const FluidState &fluidState, int phaseIdx)
Calculate the dynamic viscosity of a fluid phase .
Definition: fluidsystems/base.hh:236
static Scalar viscosity(const FluidState &fluidState, const ParameterCache &paramCache, int phaseIdx)
Calculate the dynamic viscosity of a fluid phase .
Definition: fluidsystems/base.hh:249
static Scalar heatCapacity(const FluidState &fluidState, int phaseIdx)
Specific isobaric heat capacity of a fluid phase .
Definition: fluidsystems/base.hh:424
static Scalar diffusionCoefficient(const FluidState &fluidState, const ParameterCache &paramCache, int phaseIdx, int compIdx)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase .
Definition: fluidsystems/base.hh:308
static Scalar thermalConductivity(const FluidState &fluidState, const ParameterCache &paramCache, int phaseIdx)
Thermal conductivity of a fluid phase .
Definition: fluidsystems/base.hh:403
The a parameter cache which does nothing.
Definition: nullparametercache.hh:34