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