27 #ifndef OPM_CHECK_FLUIDSYSTEM_HPP
28 #define OPM_CHECK_FLUIDSYSTEM_HPP
48 #include <dune/common/classname.hh>
57 template <
class ScalarT,
61 :
protected BaseFluidState
64 enum { numPhases = FluidSystem::numPhases };
65 enum { numComponents = FluidSystem::numComponents };
67 typedef ScalarT Scalar;
72 allowTemperature(
false);
74 allowComposition(
false);
78 restrictToPhase(1000);
81 void allowTemperature(
bool yesno)
82 { allowTemperature_ = yesno; }
84 void allowPressure(
bool yesno)
85 { allowPressure_ = yesno; }
87 void allowComposition(
bool yesno)
88 { allowComposition_ = yesno; }
90 void allowDensity(
bool yesno)
91 { allowDensity_ = yesno; }
93 void restrictToPhase(
int phaseIdx)
94 { restrictPhaseIdx_ = phaseIdx; }
96 BaseFluidState& base()
97 {
return *
static_cast<BaseFluidState*
>(
this); }
99 const BaseFluidState& base()
const
100 {
return *
static_cast<const BaseFluidState*
>(
this); }
102 auto temperature(
unsigned phaseIdx)
const
103 -> decltype(this->base().temperature(phaseIdx))
105 assert(allowTemperature_);
106 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
107 return this->base().temperature(phaseIdx);
110 auto pressure(
unsigned phaseIdx)
const
111 -> decltype(this->base().pressure(phaseIdx))
113 assert(allowPressure_);
114 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
115 return this->base().pressure(phaseIdx);
118 auto moleFraction(
unsigned phaseIdx,
unsigned compIdx)
const
119 -> decltype(this->base().moleFraction(phaseIdx, compIdx))
121 assert(allowComposition_);
122 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
123 return this->base().moleFraction(phaseIdx, compIdx);
126 auto massFraction(
unsigned phaseIdx,
unsigned compIdx)
const
127 -> decltype(this->base().massFraction(phaseIdx, compIdx))
129 assert(allowComposition_);
130 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
131 return this->base().massFraction(phaseIdx, compIdx);
134 auto averageMolarMass(
unsigned phaseIdx)
const
135 -> decltype(this->base().averageMolarMass(phaseIdx))
137 assert(allowComposition_);
138 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
139 return this->base().averageMolarMass(phaseIdx);
142 auto molarity(
unsigned phaseIdx,
unsigned compIdx)
const
143 -> decltype(this->base().molarity(phaseIdx, compIdx))
145 assert(allowDensity_ && allowComposition_);
146 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
147 return this->base().molarity(phaseIdx, compIdx);
150 auto molarDensity(
unsigned phaseIdx)
const
151 -> decltype(this->base().molarDensity(phaseIdx))
153 assert(allowDensity_);
154 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
155 return this->base().molarDensity(phaseIdx);
158 auto molarVolume(
unsigned phaseIdx)
const
159 -> decltype(this->base().molarVolume(phaseIdx))
161 assert(allowDensity_);
162 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
163 return this->base().molarVolume(phaseIdx);
166 auto density(
unsigned phaseIdx)
const
167 -> decltype(this->base().density(phaseIdx))
169 assert(allowDensity_);
170 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
171 return this->base().density(phaseIdx);
174 auto saturation(
unsigned phaseIdx)
const
175 -> decltype(this->base().saturation(phaseIdx))
178 return this->base().saturation(phaseIdx);
181 auto fugacity(
unsigned phaseIdx,
unsigned compIdx)
const
182 -> decltype(this->base().fugacity(phaseIdx, compIdx))
185 return this->base().fugacity(phaseIdx, compIdx);
188 auto fugacityCoefficient(
unsigned phaseIdx,
unsigned compIdx)
const
189 -> decltype(this->base().fugacityCoefficient(phaseIdx, compIdx))
192 return this->base().fugacityCoefficient(phaseIdx, compIdx);
195 auto enthalpy(
unsigned phaseIdx)
const
196 -> decltype(this->base().enthalpy(phaseIdx))
199 return this->base().enthalpy(phaseIdx);
202 auto internalEnergy(
unsigned phaseIdx)
const
203 -> decltype(this->base().internalEnergy(phaseIdx))
206 return this->base().internalEnergy(phaseIdx);
209 auto viscosity(
unsigned phaseIdx)
const
210 -> decltype(this->base().viscosity(phaseIdx))
213 return this->base().viscosity(phaseIdx);
217 bool allowSaturation_;
218 bool allowTemperature_;
220 bool allowComposition_;
222 int restrictPhaseIdx_;
225 template <
class Scalar,
class BaseFlu
idState>
226 void checkFluidState(
const BaseFluidState& fs)
229 [[maybe_unused]] BaseFluidState tmpFs(fs);
236 typedef typename BaseFluidState::Scalar FsScalar;
237 static_assert(std::is_same<FsScalar, Scalar>::value,
238 "Fluid states must export the type they are given as scalar in an unmodified way");
246 val = fs.temperature(0);
247 val = fs.pressure(0);
248 val = fs.moleFraction(0, 0);
249 val = fs.massFraction(0, 0);
250 val = fs.averageMolarMass(0);
251 val = fs.molarity(0, 0);
252 val = fs.molarDensity(0);
253 val = fs.molarVolume(0);
255 val = fs.saturation(0);
256 val = fs.fugacity(0, 0);
257 val = fs.fugacityCoefficient(0, 0);
258 val = fs.enthalpy(0);
259 val = fs.internalEnergy(0);
260 val = fs.viscosity(0);
267 template <
class Scalar,
class Flu
idSystem,
class RhsEval,
class LhsEval>
270 std::cout <<
"Testing fluid system '" << Dune::className<FluidSystem>() <<
"'\n";
274 enum { numPhases = FluidSystem::numPhases };
275 enum { numComponents = FluidSystem::numComponents };
279 fs.allowTemperature(
true);
280 fs.allowPressure(
true);
281 fs.allowComposition(
true);
282 fs.restrictToPhase(-1);
285 fs.base().setTemperature(273.15 + 20.0);
286 for (
int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
287 fs.base().setPressure(phaseIdx, 1e5);
288 fs.base().setSaturation(phaseIdx, 1.0/numPhases);
289 for (
int compIdx = 0; compIdx < numComponents; ++ compIdx) {
290 fs.base().setMoleFraction(phaseIdx, compIdx, 1.0/numComponents);
294 static_assert(std::is_same<typename FluidSystem::Scalar, Scalar>::value,
295 "The type used for floating point used by the fluid system must be the same"
296 " as the one passed to the checkFluidSystem() function");
299 typedef typename FluidSystem::template ParameterCache<LhsEval> ParameterCache;
301 ParameterCache paramCache;
302 try { paramCache.updateAll(fs); }
catch (...) {};
303 try { paramCache.updateAll(fs, ParameterCache::None); }
catch (...) {};
304 try { paramCache.updateAll(fs, ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); }
catch (...) {};
305 try { paramCache.updateAllPressures(fs); }
catch (...) {};
307 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
308 fs.restrictToPhase(
static_cast<int>(phaseIdx));
309 try { paramCache.updatePhase(fs, phaseIdx); }
catch (...) {};
310 try { paramCache.updatePhase(fs, phaseIdx, ParameterCache::None); }
catch (...) {};
311 try { paramCache.updatePhase(fs, phaseIdx, ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); }
catch (...) {};
312 try { paramCache.updateTemperature(fs, phaseIdx); }
catch (...) {};
313 try { paramCache.updatePressure(fs, phaseIdx); }
catch (...) {};
314 try { paramCache.updateComposition(fs, phaseIdx); }
catch (...) {};
315 try { paramCache.updateSingleMoleFraction(fs, phaseIdx, 0); }
catch (...) {};
321 Scalar scalarVal = 0.0;
323 scalarVal = 2*scalarVal;
327 try { FluidSystem::init(); }
catch (...) {};
328 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
329 fs.restrictToPhase(
static_cast<int>(phaseIdx));
330 fs.allowPressure(FluidSystem::isCompressible(phaseIdx));
331 fs.allowComposition(
true);
332 fs.allowDensity(
false);
333 try { [[maybe_unused]]
auto tmpVal = FluidSystem::density(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value,
"The default return value must be the scalar used by the fluid state!"); }
catch (...) {};
334 try { val = FluidSystem::template density<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
335 try { scalarVal = FluidSystem::template density<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
337 fs.allowPressure(
true);
338 fs.allowDensity(
true);
339 try { [[maybe_unused]]
auto tmpVal = FluidSystem::viscosity(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value,
"The default return value must be the scalar used by the fluid state!"); }
catch (...) {};
340 try { [[maybe_unused]]
auto tmpVal = FluidSystem::enthalpy(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value,
"The default return value must be the scalar used by the fluid state!"); }
catch (...) {};
341 try { [[maybe_unused]]
auto tmpVal = FluidSystem::heatCapacity(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value,
"The default return value must be the scalar used by the fluid state!"); }
catch (...) {};
342 try { [[maybe_unused]]
auto tmpVal = FluidSystem::thermalConductivity(fs, paramCache, phaseIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value,
"The default return value must be the scalar used by the fluid state!"); }
catch (...) {};
343 try { val = FluidSystem::template viscosity<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
344 try { val = FluidSystem::template enthalpy<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
345 try { val = FluidSystem::template heatCapacity<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
346 try { val = FluidSystem::template thermalConductivity<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
347 try { scalarVal = FluidSystem::template viscosity<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
348 try { scalarVal = FluidSystem::template enthalpy<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
349 try { scalarVal = FluidSystem::template heatCapacity<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
350 try { scalarVal = FluidSystem::template thermalConductivity<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
352 for (
unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
353 fs.allowComposition(!FluidSystem::isIdealMixture(phaseIdx));
354 try { [[maybe_unused]]
auto tmpVal = FluidSystem::fugacityCoefficient(fs, paramCache, phaseIdx, compIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value,
"The default return value must be the scalar used by the fluid state!"); }
catch (...) {};
355 try { val = FluidSystem::template fugacityCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
356 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
357 fs.allowComposition(
true);
358 try { [[maybe_unused]]
auto tmpVal = FluidSystem::diffusionCoefficient(fs, paramCache, phaseIdx, compIdx); static_assert(std::is_same<decltype(tmpVal), RhsEval>::value,
"The default return value must be the scalar used by the fluid state!"); }
catch (...) {};
359 try { val = FluidSystem::template diffusionCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
360 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
365 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
366 [[maybe_unused]] std::string name = FluidSystem::phaseName(phaseIdx);
367 bool bVal = FluidSystem::isLiquid(phaseIdx);
368 bVal = FluidSystem::isIdealGas(phaseIdx);
373 for (
unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
374 val = FluidSystem::molarMass(compIdx);
375 std::string name = FluidSystem::componentName(compIdx);
378 std::cout <<
"----------------------------------\n";
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system assum...
A fluid system with a liquid and a gaseous phase and water and air as components.
A fluid system with water, gas and NAPL as phases and water, air and mesitylene (DNAPL) as components...
A fluid system with water, gas and NAPL as phases and water, air and NAPL (contaminant) as components...
A two-phase fluid system with water and nitrogen as components.
A liquid-phase-only fluid system with water and nitrogen as components.
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system assum...
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system not a...
This is a fluid state which allows to set the fluid pressures and takes all other quantities from an ...
This is a fluid state which allows to set the fluid saturations and takes all other quantities from a...
A fluid system for single phase models.
The fluid system for the oil, gas and water phases of the SPE5 problem.
This is a fluid state which allows to set the fluid temperatures and takes all other quantities from ...
A fluid system for two-phase models assuming immiscibility and thermodynamic equilibrium.
void checkFluidSystem()
Checks whether a fluid system adheres to the specification.
Definition: checkFluidSystem.hpp:268
This is a fluid state which makes sure that only the quantities allowed are accessed.
Definition: checkFluidSystem.hpp:62