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 all fluid systems in opm-material
39 
40 // include all fluid states
47 
48 #include <opm/common/Unused.hpp>
49 #include <dune/common/classname.hh>
50 
51 #include <iostream>
52 #include <string>
53 
58 template <class ScalarT,
59  class FluidSystem,
62  : protected BaseFluidState
63 {
64 public:
65  enum { numPhases = FluidSystem::numPhases };
66  enum { numComponents = FluidSystem::numComponents };
67 
68  typedef ScalarT Scalar;
69 
71  {
72  // initially, do not allow anything
73  allowTemperature(false);
74  allowPressure(false);
75  allowComposition(false);
76  allowDensity(false);
77 
78  // do not allow accessing any phase
79  restrictToPhase(1000);
80  }
81 
82  void allowTemperature(bool yesno)
83  { allowTemperature_ = yesno; }
84 
85  void allowPressure(bool yesno)
86  { allowPressure_ = yesno; }
87 
88  void allowComposition(bool yesno)
89  { allowComposition_ = yesno; }
90 
91  void allowDensity(bool yesno)
92  { allowDensity_ = yesno; }
93 
94  void restrictToPhase(int phaseIdx)
95  { restrictPhaseIdx_ = phaseIdx; }
96 
97  BaseFluidState& base()
98  { return *static_cast<BaseFluidState*>(this); }
99 
100  const BaseFluidState& base() const
101  { return *static_cast<const BaseFluidState*>(this); }
102 
103  auto temperature(unsigned phaseIdx) const
104  -> decltype(this->base().temperature(phaseIdx))
105  {
106  assert(allowTemperature_);
107  assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
108  return this->base().temperature(phaseIdx);
109  }
110 
111  auto pressure(unsigned phaseIdx) const
112  -> decltype(this->base().pressure(phaseIdx))
113  {
114  assert(allowPressure_);
115  assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
116  return this->base().pressure(phaseIdx);
117  }
118 
119  auto moleFraction(unsigned phaseIdx, unsigned compIdx) const
120  -> decltype(this->base().moleFraction(phaseIdx, compIdx))
121  {
122  assert(allowComposition_);
123  assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
124  return this->base().moleFraction(phaseIdx, compIdx);
125  }
126 
127  auto massFraction(unsigned phaseIdx, unsigned compIdx) const
128  -> decltype(this->base().massFraction(phaseIdx, compIdx))
129  {
130  assert(allowComposition_);
131  assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
132  return this->base().massFraction(phaseIdx, compIdx);
133  }
134 
135  auto averageMolarMass(unsigned phaseIdx) const
136  -> decltype(this->base().averageMolarMass(phaseIdx))
137  {
138  assert(allowComposition_);
139  assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
140  return this->base().averageMolarMass(phaseIdx);
141  }
142 
143  auto molarity(unsigned phaseIdx, unsigned compIdx) const
144  -> decltype(this->base().molarity(phaseIdx, compIdx))
145  {
146  assert(allowDensity_ && allowComposition_);
147  assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
148  return this->base().molarity(phaseIdx, compIdx);
149  }
150 
151  auto molarDensity(unsigned phaseIdx) const
152  -> decltype(this->base().molarDensity(phaseIdx))
153  {
154  assert(allowDensity_);
155  assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
156  return this->base().molarDensity(phaseIdx);
157  }
158 
159  auto molarVolume(unsigned phaseIdx) const
160  -> decltype(this->base().molarVolume(phaseIdx))
161  {
162  assert(allowDensity_);
163  assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
164  return this->base().molarVolume(phaseIdx);
165  }
166 
167  auto density(unsigned phaseIdx) const
168  -> decltype(this->base().density(phaseIdx))
169  {
170  assert(allowDensity_);
171  assert(restrictPhaseIdx_ < 0 || restrictPhaseIdx_ == static_cast<int>(phaseIdx));
172  return this->base().density(phaseIdx);
173  }
174 
175  auto saturation(unsigned phaseIdx) const
176  -> decltype(this->base().saturation(phaseIdx))
177  {
178  assert(false);
179  return this->base().saturation(phaseIdx);
180  }
181 
182  auto phaseIsPresent(unsigned phaseIdx) const
183  -> decltype(this->base().phaseIsPresent(phaseIdx))
184  {
185  assert(false);
186  return this->base().phaseIsPresent(phaseIdx);
187  }
188 
189  auto fugacity(unsigned phaseIdx, unsigned compIdx) const
190  -> decltype(this->base().fugacity(phaseIdx, compIdx))
191  {
192  assert(false);
193  return this->base().fugacity(phaseIdx, compIdx);
194  }
195 
196  auto fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
197  -> decltype(this->base().fugacityCoefficient(phaseIdx, compIdx))
198  {
199  assert(false);
200  return this->base().fugacityCoefficient(phaseIdx, compIdx);
201  }
202 
203  auto enthalpy(unsigned phaseIdx) const
204  -> decltype(this->base().enthalpy(phaseIdx))
205  {
206  assert(false);
207  return this->base().enthalpy(phaseIdx);
208  }
209 
210  auto internalEnergy(unsigned phaseIdx) const
211  -> decltype(this->base().internalEnergy(phaseIdx))
212  {
213  assert(false);
214  return this->base().internalEnergy(phaseIdx);
215  }
216 
217  auto viscosity(unsigned phaseIdx) const
218  -> decltype(this->base().viscosity(phaseIdx))
219  {
220  assert(false);
221  return this->base().viscosity(phaseIdx);
222  }
223 
224 private:
225  bool allowSaturation_;
226  bool allowTemperature_;
227  bool allowPressure_;
228  bool allowComposition_;
229  bool allowDensity_;
230  int restrictPhaseIdx_;
231 };
232 
233 template <class Scalar, class BaseFluidState>
234 void checkFluidState(const BaseFluidState& fs)
235 {
236  // fluid states must be copy-able
237  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  bool b OPM_UNUSED = fs.phaseIsPresent(/*phaseIdx=*/0);
265  val = fs.fugacity(/*phaseIdx=*/0, /*compIdx=*/0);
266  val = fs.fugacityCoefficient(/*phaseIdx=*/0, /*compIdx=*/0);
267  val = fs.enthalpy(/*phaseIdx=*/0);
268  val = fs.internalEnergy(/*phaseIdx=*/0);
269  val = fs.viscosity(/*phaseIdx=*/0);
270  };
271 }
272 
276 template <class Scalar, class FluidSystem, class RhsEval, class LhsEval>
278 {
279  std::cout << "Testing fluid system '" << Dune::className<FluidSystem>() << "'\n";
280 
281  // make sure the fluid system provides the number of phases and
282  // the number of components
283  enum { numPhases = FluidSystem::numPhases };
284  enum { numComponents = FluidSystem::numComponents };
285 
287  FluidState fs;
288  fs.allowTemperature(true);
289  fs.allowPressure(true);
290  fs.allowComposition(true);
291  fs.restrictToPhase(-1);
292 
293  // initialize memory the fluid state
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);
300  }
301  }
302 
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");
306 
307  // check whether the parameter cache adheres to the API
308  typedef typename FluidSystem::template ParameterCache<LhsEval> ParameterCache;
309 
310  ParameterCache paramCache;
311  try { paramCache.updateAll(fs); } catch (...) {};
312  try { paramCache.updateAll(fs, /*except=*/ParameterCache::None); } catch (...) {};
313  try { paramCache.updateAll(fs, /*except=*/ParameterCache::Temperature | ParameterCache::Pressure | ParameterCache::Composition); } catch (...) {};
314  try { paramCache.updateAllPressures(fs); } catch (...) {};
315 
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, /*except=*/ParameterCache::None); } catch (...) {};
320  try { paramCache.updatePhase(fs, phaseIdx, /*except=*/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, /*compIdx=*/0); } catch (...) {};
325  }
326 
327  // some value to make sure the return values of the fluid system
328  // are convertible to scalars
329  LhsEval val = 0.0;
330  Scalar scalarVal = 0.0;
331 
332  scalarVal = 2*scalarVal; // get rid of GCC warning (only occurs with paranoid warning flags)
333  val = 2*val; // get rid of GCC warning (only occurs with paranoid warning flags)
334 
335  // actually check the fluid system API
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 (...) {};
345 
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 (...) {};
360 
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 (...) {};
370  }
371  }
372 
373  // test for phaseName(), isLiquid() and isIdealGas()
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);
378  bVal = !bVal; // get rid of GCC warning (only occurs with paranoid warning flags)
379  }
380 
381  // test for molarMass() and componentName()
382  for (unsigned compIdx = 0; compIdx < numComponents; ++ compIdx) {
383  val = FluidSystem::molarMass(compIdx);
384  std::string name = FluidSystem::componentName(compIdx);
385  }
386 
387  std::cout << "----------------------------------\n";
388 }
389 
390 #endif
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