27#ifndef OPM_CHECK_FLUIDSYSTEM_HPP
28#define OPM_CHECK_FLUIDSYSTEM_HPP
30#include <opm/common/utility/Demangle.hpp>
58template <
class ScalarT,
62 :
protected BaseFluidState
65 enum { numPhases = FluidSystem::numPhases };
66 enum { numComponents = FluidSystem::numComponents };
68 typedef ScalarT Scalar;
73 allowTemperature(
false);
75 allowComposition(
false);
79 restrictToPhase(1000);
82 void allowTemperature(
bool yesno)
83 { allowTemperature_ = yesno; }
85 void allowPressure(
bool yesno)
86 { allowPressure_ = yesno; }
88 void allowComposition(
bool yesno)
89 { allowComposition_ = yesno; }
91 void allowDensity(
bool yesno)
92 { allowDensity_ = yesno; }
94 void restrictToPhase(
int phaseIdx)
95 { restrictPhaseIdx_ = phaseIdx; }
97 BaseFluidState& base()
98 {
return *
static_cast<BaseFluidState*
>(
this); }
100 const BaseFluidState& base()
const
101 {
return *
static_cast<const BaseFluidState*
>(
this); }
103 auto temperature(
unsigned phaseIdx)
const
104 ->
decltype(this->base().temperature(phaseIdx))
106 assert(allowTemperature_);
107 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
108 return this->base().temperature(phaseIdx);
111 auto pressure(
unsigned phaseIdx)
const
112 ->
decltype(this->base().pressure(phaseIdx))
114 assert(allowPressure_);
115 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
116 return this->base().pressure(phaseIdx);
119 auto moleFraction(
unsigned phaseIdx,
unsigned compIdx)
const
120 ->
decltype(this->base().moleFraction(phaseIdx, compIdx))
122 assert(allowComposition_);
123 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
124 return this->base().moleFraction(phaseIdx, compIdx);
127 auto massFraction(
unsigned phaseIdx,
unsigned compIdx)
const
128 ->
decltype(this->base().massFraction(phaseIdx, compIdx))
130 assert(allowComposition_);
131 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
132 return this->base().massFraction(phaseIdx, compIdx);
135 auto averageMolarMass(
unsigned phaseIdx)
const
136 ->
decltype(this->base().averageMolarMass(phaseIdx))
138 assert(allowComposition_);
139 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
140 return this->base().averageMolarMass(phaseIdx);
143 auto molarity(
unsigned phaseIdx,
unsigned compIdx)
const
144 ->
decltype(this->base().molarity(phaseIdx, compIdx))
146 assert(allowDensity_ && allowComposition_);
147 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
148 return this->base().molarity(phaseIdx, compIdx);
151 auto molarDensity(
unsigned phaseIdx)
const
152 ->
decltype(this->base().molarDensity(phaseIdx))
154 assert(allowDensity_);
155 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
156 return this->base().molarDensity(phaseIdx);
159 auto molarVolume(
unsigned phaseIdx)
const
160 ->
decltype(this->base().molarVolume(phaseIdx))
162 assert(allowDensity_);
163 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
164 return this->base().molarVolume(phaseIdx);
167 auto density(
unsigned phaseIdx)
const
168 ->
decltype(this->base().density(phaseIdx))
170 assert(allowDensity_);
171 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ ==
static_cast<int>(phaseIdx));
172 return this->base().density(phaseIdx);
175 auto saturation(
unsigned phaseIdx)
const
176 ->
decltype(this->base().saturation(phaseIdx))
179 return this->base().saturation(phaseIdx);
182 auto fugacity(
unsigned phaseIdx,
unsigned compIdx)
const
183 ->
decltype(this->base().fugacity(phaseIdx, compIdx))
186 return this->base().fugacity(phaseIdx, compIdx);
189 auto fugacityCoefficient(
unsigned phaseIdx,
unsigned compIdx)
const
190 ->
decltype(this->base().fugacityCoefficient(phaseIdx, compIdx))
193 return this->base().fugacityCoefficient(phaseIdx, compIdx);
196 auto enthalpy(
unsigned phaseIdx)
const
197 ->
decltype(this->base().enthalpy(phaseIdx))
200 return this->base().enthalpy(phaseIdx);
203 auto internalEnergy(
unsigned phaseIdx)
const
204 ->
decltype(this->base().internalEnergy(phaseIdx))
207 return this->base().internalEnergy(phaseIdx);
210 auto viscosity(
unsigned phaseIdx)
const
211 ->
decltype(this->base().viscosity(phaseIdx))
214 return this->base().viscosity(phaseIdx);
218 bool allowSaturation_;
219 bool allowTemperature_;
221 bool allowComposition_;
223 int restrictPhaseIdx_;
226template <
class Scalar,
class BaseFlu
idState>
227void checkFluidState(
const BaseFluidState& fs)
230 [[maybe_unused]] BaseFluidState tmpFs(fs);
237 typedef typename BaseFluidState::Scalar FsScalar;
238 static_assert(std::is_same<FsScalar, Scalar>::value,
239 "Fluid states must export the type they are given as scalar in an unmodified way");
247 val = fs.temperature(0);
248 val = fs.pressure(0);
249 val = fs.moleFraction(0, 0);
250 val = fs.massFraction(0, 0);
251 val = fs.averageMolarMass(0);
252 val = fs.molarity(0, 0);
253 val = fs.molarDensity(0);
254 val = fs.molarVolume(0);
256 val = fs.saturation(0);
257 val = fs.fugacity(0, 0);
258 val = fs.fugacityCoefficient(0, 0);
259 val = fs.enthalpy(0);
260 val = fs.internalEnergy(0);
261 val = fs.viscosity(0);
268template <
class Scalar,
class Flu
idSystem,
class RhsEval,
class LhsEval>
271 std::cout <<
"Testing fluid system '" <<
Opm::demangle(
typeid(FluidSystem).name()) <<
"'\n";
275 enum { numPhases = FluidSystem::numPhases };
276 enum { numComponents = FluidSystem::numComponents };
280 fs.allowTemperature(
true);
281 fs.allowPressure(
true);
282 fs.allowComposition(
true);
283 fs.restrictToPhase(-1);
286 fs.base().setTemperature(273.15 + 20.0);
287 for (
int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
288 fs.base().setPressure(phaseIdx, 1e5);
289 fs.base().setSaturation(phaseIdx, 1.0/numPhases);
290 for (
int compIdx = 0; compIdx < numComponents; ++ compIdx) {
291 fs.base().setMoleFraction(phaseIdx, compIdx, 1.0/numComponents);
295 static_assert(std::is_same<typename FluidSystem::Scalar, Scalar>::value,
296 "The type used for floating point used by the fluid system must be the same"
297 " as the one passed to the checkFluidSystem() function");
300 typedef typename FluidSystem::template ParameterCache<LhsEval> ParameterCache;
302 ParameterCache paramCache;
303 try { paramCache.updateAll(fs); }
catch (...) {};
304 try { paramCache.updateAll(fs, ParameterCache::None); }
catch (...) {};
305 try { paramCache.updateAll(fs, ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); }
catch (...) {};
306 try { paramCache.updateAllPressures(fs); }
catch (...) {};
308 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
309 fs.restrictToPhase(
static_cast<int>(phaseIdx));
310 try { paramCache.updatePhase(fs, phaseIdx); }
catch (...) {};
311 try { paramCache.updatePhase(fs, phaseIdx, ParameterCache::None); }
catch (...) {};
312 try { paramCache.updatePhase(fs, phaseIdx, ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); }
catch (...) {};
313 try { paramCache.updateTemperature(fs, phaseIdx); }
catch (...) {};
314 try { paramCache.updatePressure(fs, phaseIdx); }
catch (...) {};
315 try { paramCache.updateComposition(fs, phaseIdx); }
catch (...) {};
316 try { paramCache.updateSingleMoleFraction(fs, phaseIdx, 0); }
catch (...) {};
322 Scalar scalarVal = 0.0;
324 scalarVal = 2*scalarVal;
328 try { FluidSystem::init(); }
catch (...) {};
329 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
330 fs.restrictToPhase(
static_cast<int>(phaseIdx));
331 fs.allowPressure(FluidSystem::isCompressible(phaseIdx));
332 fs.allowComposition(
true);
333 fs.allowDensity(
false);
334 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 (...) {};
335 try { val = FluidSystem::template density<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
336 try { scalarVal = FluidSystem::template density<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
338 fs.allowPressure(
true);
339 fs.allowDensity(
true);
340 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 (...) {};
341 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 (...) {};
342 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 (...) {};
343 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 (...) {};
344 try { val = FluidSystem::template viscosity<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
345 try { val = FluidSystem::template enthalpy<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
346 try { val = FluidSystem::template heatCapacity<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
347 try { val = FluidSystem::template thermalConductivity<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
348 try { scalarVal = FluidSystem::template viscosity<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
349 try { scalarVal = FluidSystem::template enthalpy<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
350 try { scalarVal = FluidSystem::template heatCapacity<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
351 try { scalarVal = FluidSystem::template thermalConductivity<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
353 for (
unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
354 fs.allowComposition(!FluidSystem::isIdealMixture(phaseIdx));
355 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 (...) {};
356 try { val = FluidSystem::template fugacityCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
357 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
358 fs.allowComposition(
true);
359 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 (...) {};
360 try { val = FluidSystem::template diffusionCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
361 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
366 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
367 [[maybe_unused]] std::string name = FluidSystem::phaseName(phaseIdx);
368 bool bVal = FluidSystem::isLiquid(phaseIdx);
369 bVal = FluidSystem::isIdealGas(phaseIdx);
374 for (
unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
375 val = FluidSystem::molarMass(compIdx);
376 std::string name = FluidSystem::componentName(compIdx);
379 std::cout <<
"----------------------------------\n";
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 assum...
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system not a...
void checkFluidSystem()
Checks whether a fluid system adheres to the specification.
Definition: checkFluidSystem.hpp:269
This is a fluid state which makes sure that only the quantities allowed are accessed.
Definition: checkFluidSystem.hpp:63
std::string demangle(const char *mangled_symbol)
Returns demangled name of symbol.