My Project
Loading...
Searching...
No Matches
checkFluidSystem.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_CHECK_FLUIDSYSTEM_HPP
28#define OPM_CHECK_FLUIDSYSTEM_HPP
29
30#include <opm/common/utility/Demangle.hpp>
31
32// include all fluid systems in opm-material
41
42// include all fluid states
49
50#include <iostream>
51#include <string>
52#include <string_view>
53#include <typeinfo>
54
59template <class ScalarT,
60 class FluidSystem,
63 : protected BaseFluidState
64{
65public:
66 enum { numPhases = FluidSystem::numPhases };
67 enum { numComponents = FluidSystem::numComponents };
68
69 typedef ScalarT Scalar;
70
72 {
73 // initially, do not allow anything
74 allowTemperature(false);
75 allowPressure(false);
76 allowComposition(false);
77 allowDensity(false);
78
79 // do not allow accessing any phase
80 restrictToPhase(1000);
81 }
82
83 void allowTemperature(bool yesno)
84 { allowTemperature_ = yesno; }
85
86 void allowPressure(bool yesno)
87 { allowPressure_ = yesno; }
88
89 void allowComposition(bool yesno)
90 { allowComposition_ = yesno; }
91
92 void allowDensity(bool yesno)
93 { allowDensity_ = yesno; }
94
95 void restrictToPhase(int phaseIdx)
96 { restrictPhaseIdx_ = phaseIdx; }
97
98 BaseFluidState& base()
99 { return *static_cast<BaseFluidState*>(this); }
100
101 const BaseFluidState& base() const
102 { return *static_cast<const BaseFluidState*>(this); }
103
104 auto temperature(unsigned phaseIdx) const
105 -> decltype(this->base().temperature(phaseIdx))
106 {
107 assert(allowTemperature_);
108 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
109 return this->base().temperature(phaseIdx);
110 }
111
112 auto pressure(unsigned phaseIdx) const
113 -> decltype(this->base().pressure(phaseIdx))
114 {
115 assert(allowPressure_);
116 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
117 return this->base().pressure(phaseIdx);
118 }
119
120 auto moleFraction(unsigned phaseIdx, unsigned compIdx) const
121 -> decltype(this->base().moleFraction(phaseIdx, compIdx))
122 {
123 assert(allowComposition_);
124 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
125 return this->base().moleFraction(phaseIdx, compIdx);
126 }
127
128 auto massFraction(unsigned phaseIdx, unsigned compIdx) const
129 -> decltype(this->base().massFraction(phaseIdx, compIdx))
130 {
131 assert(allowComposition_);
132 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
133 return this->base().massFraction(phaseIdx, compIdx);
134 }
135
136 auto averageMolarMass(unsigned phaseIdx) const
137 -> decltype(this->base().averageMolarMass(phaseIdx))
138 {
139 assert(allowComposition_);
140 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
141 return this->base().averageMolarMass(phaseIdx);
142 }
143
144 auto molarity(unsigned phaseIdx, unsigned compIdx) const
145 -> decltype(this->base().molarity(phaseIdx, compIdx))
146 {
147 assert(allowDensity_ && allowComposition_);
148 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
149 return this->base().molarity(phaseIdx, compIdx);
150 }
151
152 auto molarDensity(unsigned phaseIdx) const
153 -> decltype(this->base().molarDensity(phaseIdx))
154 {
155 assert(allowDensity_);
156 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
157 return this->base().molarDensity(phaseIdx);
158 }
159
160 auto molarVolume(unsigned phaseIdx) const
161 -> decltype(this->base().molarVolume(phaseIdx))
162 {
163 assert(allowDensity_);
164 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
165 return this->base().molarVolume(phaseIdx);
166 }
167
168 auto density(unsigned phaseIdx) const
169 -> decltype(this->base().density(phaseIdx))
170 {
171 assert(allowDensity_);
172 assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
173 return this->base().density(phaseIdx);
174 }
175
176 auto saturation(unsigned phaseIdx) const
177 -> decltype(this->base().saturation(phaseIdx))
178 {
179 assert(false);
180 return this->base().saturation(phaseIdx);
181 }
182
183 auto fugacity(unsigned phaseIdx, unsigned compIdx) const
184 -> decltype(this->base().fugacity(phaseIdx, compIdx))
185 {
186 assert(false);
187 return this->base().fugacity(phaseIdx, compIdx);
188 }
189
190 auto fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
191 -> decltype(this->base().fugacityCoefficient(phaseIdx, compIdx))
192 {
193 assert(false);
194 return this->base().fugacityCoefficient(phaseIdx, compIdx);
195 }
196
197 auto enthalpy(unsigned phaseIdx) const
198 -> decltype(this->base().enthalpy(phaseIdx))
199 {
200 assert(false);
201 return this->base().enthalpy(phaseIdx);
202 }
203
204 auto internalEnergy(unsigned phaseIdx) const
205 -> decltype(this->base().internalEnergy(phaseIdx))
206 {
207 assert(false);
208 return this->base().internalEnergy(phaseIdx);
209 }
210
211 auto viscosity(unsigned phaseIdx) const
212 -> decltype(this->base().viscosity(phaseIdx))
213 {
214 assert(false);
215 return this->base().viscosity(phaseIdx);
216 }
217
218 auto compressFactor(unsigned phaseIdx) const
219 -> decltype(this->base().compressFactor(phaseIdx))
220 {
221 return this->base().compressFactor(phaseIdx);
222 }
223
224private:
225 bool allowSaturation_;
226 bool allowTemperature_;
227 bool allowPressure_;
228 bool allowComposition_;
229 bool allowDensity_;
230 int restrictPhaseIdx_;
231};
232
233template <class Scalar, class BaseFluidState>
234void checkFluidState(const BaseFluidState& fs)
235{
236 // fluid states must be copy-able
237 [[maybe_unused]] BaseFluidState tmpFs(fs);
238 tmpFs = fs;
239
240 // a fluid state must provide a checkDefined() method
241 fs.checkDefined();
242
243 // fluid states must export the types which they use as Scalars
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");
247
248 // make sure the fluid state provides all mandatory methods
249 while (false) {
250 Scalar val = 1.0;
251
252 val = 2*val; // get rid of GCC warning (only occurs with paranoid warning flags)
253
254 val = fs.temperature(/*phaseIdx=*/0);
255 val = fs.pressure(/*phaseIdx=*/0);
256 val = fs.moleFraction(/*phaseIdx=*/0, /*compIdx=*/0);
257 val = fs.massFraction(/*phaseIdx=*/0, /*compIdx=*/0);
258 val = fs.averageMolarMass(/*phaseIdx=*/0);
259 val = fs.molarity(/*phaseIdx=*/0, /*compIdx=*/0);
260 val = fs.molarDensity(/*phaseIdx=*/0);
261 val = fs.molarVolume(/*phaseIdx=*/0);
262 val = fs.density(/*phaseIdx=*/0);
263 val = fs.saturation(/*phaseIdx=*/0);
264 val = fs.fugacity(/*phaseIdx=*/0, /*compIdx=*/0);
265 val = fs.fugacityCoefficient(/*phaseIdx=*/0, /*compIdx=*/0);
266 val = fs.enthalpy(/*phaseIdx=*/0);
267 val = fs.internalEnergy(/*phaseIdx=*/0);
268 val = fs.viscosity(/*phaseIdx=*/0);
269 };
270}
271
275template <class Scalar, class FluidSystem, class RhsEval, class LhsEval>
277{
278 std::cout << "Testing fluid system '"
279 << Opm::demangle(typeid(FluidSystem).name())
280 << ", RhsEval = " << Opm::demangle(typeid(RhsEval).name())
281 << ", LhsEval = " << Opm::demangle(typeid(LhsEval).name()) << "'\n";
282
283 // make sure the fluid system provides the number of phases and
284 // the number of components
285 enum { numPhases = FluidSystem::numPhases };
286 enum { numComponents = FluidSystem::numComponents };
287
289 FluidState fs;
290 fs.allowTemperature(true);
291 fs.allowPressure(true);
292 fs.allowComposition(true);
293 fs.restrictToPhase(-1);
294
295 // initialize memory the fluid state
296 fs.base().setTemperature(273.15 + 20.0);
297 for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
298 fs.base().setPressure(phaseIdx, 1e5);
299 fs.base().setSaturation(phaseIdx, 1.0/numPhases);
300 for (int compIdx = 0; compIdx < numComponents; ++ compIdx) {
301 fs.base().setMoleFraction(phaseIdx, compIdx, 1.0/numComponents);
302 }
303 }
304
305 static_assert(std::is_same<typename FluidSystem::Scalar, Scalar>::value,
306 "The type used for floating point used by the fluid system must be the same"
307 " as the one passed to the checkFluidSystem() function");
308
309 // check whether the parameter cache adheres to the API
310 typedef typename FluidSystem::template ParameterCache<LhsEval> ParameterCache;
311
312 ParameterCache paramCache;
313 try { paramCache.updateAll(fs); } catch (...) {};
314 try { paramCache.updateAll(fs, /*except=*/ParameterCache::None); } catch (...) {};
315 try { paramCache.updateAll(fs, /*except=*/ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); } catch (...) {};
316 try { paramCache.updateAllPressures(fs); } catch (...) {};
317
318 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
319 fs.restrictToPhase(static_cast<int>(phaseIdx));
320 try { paramCache.updatePhase(fs, phaseIdx); } catch (...) {};
321 try { paramCache.updatePhase(fs, phaseIdx, /*except=*/ParameterCache::None); } catch (...) {};
322 try { paramCache.updatePhase(fs, phaseIdx, /*except=*/ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); } catch (...) {};
323 try { paramCache.updateTemperature(fs, phaseIdx); } catch (...) {};
324 try { paramCache.updatePressure(fs, phaseIdx); } catch (...) {};
325 try { paramCache.updateComposition(fs, phaseIdx); } catch (...) {};
326 try { paramCache.updateSingleMoleFraction(fs, phaseIdx, /*compIdx=*/0); } catch (...) {};
327 }
328
329 // some value to make sure the return values of the fluid system
330 // are convertible to scalars
331 LhsEval val = 0.0;
332 Scalar scalarVal = 0.0;
333
334 scalarVal = 2*scalarVal; // get rid of GCC warning (only occurs with paranoid warning flags)
335 val = 2*val; // get rid of GCC warning (only occurs with paranoid warning flags)
336
337 // actually check the fluid system API
338 try { FluidSystem::init(); } catch (...) {};
339 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
340 fs.restrictToPhase(static_cast<int>(phaseIdx));
341 fs.allowPressure(FluidSystem::isCompressible(phaseIdx));
342 fs.allowComposition(true);
343 fs.allowDensity(false);
344 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 (...) {};
345 try { val = FluidSystem::template density<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
346 try { scalarVal = FluidSystem::template density<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
347
348 fs.allowPressure(true);
349 fs.allowDensity(true);
350 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 (...) {};
351 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 (...) {};
352 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 (...) {};
353 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 (...) {};
354 try { val = FluidSystem::template viscosity<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
355 try { val = FluidSystem::template enthalpy<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
356 try { val = FluidSystem::template heatCapacity<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
357 try { val = FluidSystem::template thermalConductivity<FluidState, LhsEval>(fs, paramCache, phaseIdx); } catch (...) {};
358 try { scalarVal = FluidSystem::template viscosity<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
359 try { scalarVal = FluidSystem::template enthalpy<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
360 try { scalarVal = FluidSystem::template heatCapacity<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
361 try { scalarVal = FluidSystem::template thermalConductivity<FluidState, Scalar>(fs, paramCache, phaseIdx); } catch (...) {};
362
363 for (unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
364 fs.allowComposition(!FluidSystem::isIdealMixture(phaseIdx));
365 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 (...) {};
366 try { val = FluidSystem::template fugacityCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); } catch (...) {};
367 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); } catch (...) {};
368 fs.allowComposition(true);
369 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 (...) {};
370 try { val = FluidSystem::template diffusionCoefficient<FluidState, LhsEval>(fs, paramCache, phaseIdx, compIdx); } catch (...) {};
371 try { scalarVal = FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fs, paramCache, phaseIdx, compIdx); } catch (...) {};
372 }
373 }
374
375 // test for phaseName(), isLiquid() and isIdealGas()
376 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
377 [[maybe_unused]] std::string name{FluidSystem::phaseName(phaseIdx)};
378 bool bVal = FluidSystem::isLiquid(phaseIdx);
379 bVal = FluidSystem::isIdealGas(phaseIdx);
380 bVal = !bVal; // get rid of GCC warning (only occurs with paranoid warning flags)
381 }
382
383 // test for molarMass() and componentName()
384 for (unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
385 val = FluidSystem::molarMass(compIdx);
386 std::string{FluidSystem::componentName(compIdx)};
387 }
388}
389
390#endif
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:276
This is a fluid state which makes sure that only the quantities allowed are accessed.
Definition checkFluidSystem.hpp:64
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system assum...
Definition CompositionalFluidState.hpp:46
std::string demangle(const char *mangled_symbol)
Returns demangled name of symbol.
Definition Demangle.cpp:32