27 #ifndef OPM_CHECK_FLUIDSYSTEM_HPP
28 #define OPM_CHECK_FLUIDSYSTEM_HPP
48 #include <opm/common/Unused.hpp>
49 #include <dune/common/classname.hh>
58 template <
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 phaseIsPresent(
unsigned phaseIdx)
const
183 -> decltype(this->base().phaseIsPresent(phaseIdx))
186 return this->base().phaseIsPresent(phaseIdx);
189 auto fugacity(
unsigned phaseIdx,
unsigned compIdx)
const
190 -> decltype(this->base().fugacity(phaseIdx, compIdx))
193 return this->base().fugacity(phaseIdx, compIdx);
196 auto fugacityCoefficient(
unsigned phaseIdx,
unsigned compIdx)
const
197 -> decltype(this->base().fugacityCoefficient(phaseIdx, compIdx))
200 return this->base().fugacityCoefficient(phaseIdx, compIdx);
203 auto enthalpy(
unsigned phaseIdx)
const
204 -> decltype(this->base().enthalpy(phaseIdx))
207 return this->base().enthalpy(phaseIdx);
210 auto internalEnergy(
unsigned phaseIdx)
const
211 -> decltype(this->base().internalEnergy(phaseIdx))
214 return this->base().internalEnergy(phaseIdx);
217 auto viscosity(
unsigned phaseIdx)
const
218 -> decltype(this->base().viscosity(phaseIdx))
221 return this->base().viscosity(phaseIdx);
225 bool allowSaturation_;
226 bool allowTemperature_;
228 bool allowComposition_;
230 int restrictPhaseIdx_;
233 template <
class Scalar,
class BaseFlu
idState>
234 void checkFluidState(
const BaseFluidState& fs)
237 BaseFluidState tmpFs(fs);
244 typedef typename BaseFluidState::Scalar FsScalar;
245 static_assert(std::is_same<FsScalar, Scalar>::value,
246 "Fluid states must export the type they are given as scalar in an unmodified way");
254 val = fs.temperature(0);
255 val = fs.pressure(0);
256 val = fs.moleFraction(0, 0);
257 val = fs.massFraction(0, 0);
258 val = fs.averageMolarMass(0);
259 val = fs.molarity(0, 0);
260 val = fs.molarDensity(0);
261 val = fs.molarVolume(0);
263 val = fs.saturation(0);
264 bool b OPM_UNUSED = fs.phaseIsPresent(0);
265 val = fs.fugacity(0, 0);
266 val = fs.fugacityCoefficient(0, 0);
267 val = fs.enthalpy(0);
268 val = fs.internalEnergy(0);
269 val = fs.viscosity(0);
276 template <
class Scalar,
class Flu
idSystem,
class RhsEval,
class LhsEval>
279 std::cout <<
"Testing fluid system '" << Dune::className<FluidSystem>() <<
"'\n";
283 enum { numPhases = FluidSystem::numPhases };
284 enum { numComponents = FluidSystem::numComponents };
288 fs.allowTemperature(
true);
289 fs.allowPressure(
true);
290 fs.allowComposition(
true);
291 fs.restrictToPhase(-1);
294 fs.base().setTemperature(273.15 + 20.0);
295 for (
int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
296 fs.base().setPressure(phaseIdx, 1e5);
297 fs.base().setSaturation(phaseIdx, 1.0/numPhases);
298 for (
int compIdx = 0; compIdx < numComponents; ++ compIdx) {
299 fs.base().setMoleFraction(phaseIdx, compIdx, 1.0/numComponents);
303 static_assert(std::is_same<typename FluidSystem::Scalar, Scalar>::value,
304 "The type used for floating point used by the fluid system must be the same"
305 " as the one passed to the checkFluidSystem() function");
308 typedef typename FluidSystem::template ParameterCache<LhsEval> ParameterCache;
310 ParameterCache paramCache;
311 try { paramCache.updateAll(fs); }
catch (...) {};
312 try { paramCache.updateAll(fs, ParameterCache::None); }
catch (...) {};
313 try { paramCache.updateAll(fs, ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); }
catch (...) {};
314 try { paramCache.updateAllPressures(fs); }
catch (...) {};
316 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
317 fs.restrictToPhase(static_cast<int>(phaseIdx));
318 try { paramCache.updatePhase(fs, phaseIdx); }
catch (...) {};
319 try { paramCache.updatePhase(fs, phaseIdx, ParameterCache::None); }
catch (...) {};
320 try { paramCache.updatePhase(fs, phaseIdx, ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); }
catch (...) {};
321 try { paramCache.updateTemperature(fs, phaseIdx); }
catch (...) {};
322 try { paramCache.updatePressure(fs, phaseIdx); }
catch (...) {};
323 try { paramCache.updateComposition(fs, phaseIdx); }
catch (...) {};
324 try { paramCache.updateSingleMoleFraction(fs, phaseIdx, 0); }
catch (...) {};
330 Scalar scalarVal = 0.0;
332 scalarVal = 2*scalarVal;
336 try { FluidSystem::init(); }
catch (...) {};
337 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
338 fs.restrictToPhase(static_cast<int>(phaseIdx));
339 fs.allowPressure(FluidSystem::isCompressible(phaseIdx));
340 fs.allowComposition(
true);
341 fs.allowDensity(
false);
342 try {
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 (...) {};
343 try { val = FluidSystem::template density<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
344 try { scalarVal = FluidSystem::template density<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
346 fs.allowPressure(
true);
347 fs.allowDensity(
true);
348 try {
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 (...) {};
349 try {
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 (...) {};
350 try {
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 (...) {};
351 try {
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 (...) {};
352 try { val = FluidSystem::template viscosity<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
353 try { val = FluidSystem::template enthalpy<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
354 try { val = FluidSystem::template heatCapacity<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
355 try { val = FluidSystem::template thermalConductivity<FluidState, LhsEval>(fs, paramCache, phaseIdx); }
catch (...) {};
356 try { scalarVal = FluidSystem::template viscosity<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
357 try { scalarVal = FluidSystem::template enthalpy<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
358 try { scalarVal = FluidSystem::template heatCapacity<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
359 try { scalarVal = FluidSystem::template thermalConductivity<FluidState, Scalar>(fs, paramCache, phaseIdx); }
catch (...) {};
361 for (
unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
362 fs.allowComposition(!FluidSystem::isIdealMixture(phaseIdx));
363 try {
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 (...) {};
364 try { val = FluidSystem::template fugacityCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
365 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
366 fs.allowComposition(
true);
367 try {
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 (...) {};
368 try { val = FluidSystem::template diffusionCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
369 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); }
catch (...) {};
374 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
375 std::string name OPM_UNUSED = FluidSystem::phaseName(phaseIdx);
376 bool bVal = FluidSystem::isLiquid(phaseIdx);
377 bVal = FluidSystem::isIdealGas(phaseIdx);
382 for (
unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
383 val = FluidSystem::molarMass(compIdx);
384 std::string name = FluidSystem::componentName(compIdx);
387 std::cout <<
"----------------------------------\n";
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system assum...
This is a fluid state which allows to set the fluid temperatures and takes all other quantities from ...
A fluid system with water, gas and NAPL as phases and water, air and NAPL (contaminant) as components...
A fluid system for single phase models.
This is a fluid state which makes sure that only the quantities allowed are accessed.
Definition: checkFluidSystem.hpp:61
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system assum...
A fluid system with water, gas and NAPL as phases and water, air and mesitylene (DNAPL) as components...
A two-phase fluid system with water and nitrogen as components.
This is a fluid state which allows to set the fluid saturations and takes all other quantities from a...
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 pressures and takes all other quantities from an ...
A liquid-phase-only fluid system with water and nitrogen as components.
A fluid system with a liquid and a gaseous phase and water and air as components. ...
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:277
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system not a...
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system assum...
Definition: CompositionalFluidState.hpp:46