FluidStateCompositionModules.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 */
28 #ifndef OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
29 #define OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
30 
31 #include <opm/common/Valgrind.hpp>
33 
34 #include <opm/common/ErrorMacros.hpp>
35 #include <opm/common/Exceptions.hpp>
36 
37 #include <algorithm>
38 #include <cmath>
39 
40 namespace Opm {
41 
46 template <class Scalar,
47  class FluidSystem,
48  class Implementation>
50 {
51  enum { numPhases = FluidSystem::numPhases };
52  enum { numComponents = FluidSystem::numComponents };
53 
54 public:
56  {
57  for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
58  for (int compIdx = 0; compIdx < numComponents; ++compIdx)
59  moleFraction_[phaseIdx][compIdx] = 0.0;
60 
61  Valgrind::SetDefined(moleFraction_);
62  Valgrind::SetUndefined(averageMolarMass_);
63  Valgrind::SetUndefined(sumMoleFractions_);
64  }
65 
69  const Scalar& moleFraction(unsigned phaseIdx, unsigned compIdx) const
70  { return moleFraction_[phaseIdx][compIdx]; }
71 
75  Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
76  {
77  return
78  Opm::abs(sumMoleFractions_[phaseIdx])
79  *moleFraction_[phaseIdx][compIdx]
80  *FluidSystem::molarMass(compIdx)
81  / Opm::max(1e-40, Opm::abs(averageMolarMass_[phaseIdx]));
82  }
83 
92  const Scalar& averageMolarMass(unsigned phaseIdx) const
93  { return averageMolarMass_[phaseIdx]; }
94 
104  Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
105  { return asImp_().molarDensity(phaseIdx)*moleFraction(phaseIdx, compIdx); }
106 
112  void setMoleFraction(unsigned phaseIdx, unsigned compIdx, const Scalar& value)
113  {
114  Valgrind::CheckDefined(value);
115  Valgrind::SetUndefined(sumMoleFractions_[phaseIdx]);
116  Valgrind::SetUndefined(averageMolarMass_[phaseIdx]);
117  Valgrind::SetUndefined(moleFraction_[phaseIdx][compIdx]);
118 
119  moleFraction_[phaseIdx][compIdx] = value;
120 
121  // re-calculate the mean molar mass
122  sumMoleFractions_[phaseIdx] = 0.0;
123  averageMolarMass_[phaseIdx] = 0.0;
124  for (unsigned compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
125  sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
126  averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
127  }
128  }
129 
134  template <class FluidState>
135  void assign(const FluidState& fs)
136  {
137  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
138  averageMolarMass_[phaseIdx] = 0;
139  sumMoleFractions_[phaseIdx] = 0;
140  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
141  moleFraction_[phaseIdx][compIdx] =
142  Opm::decay<Scalar>(fs.moleFraction(phaseIdx, compIdx));
143 
144  averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
145  sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
146  }
147  }
148  }
149 
158  void checkDefined() const
159  {
160  Valgrind::CheckDefined(moleFraction_);
161  Valgrind::CheckDefined(averageMolarMass_);
162  Valgrind::CheckDefined(sumMoleFractions_);
163  }
164 
165 protected:
166  const Implementation& asImp_() const
167  { return *static_cast<const Implementation*>(this); }
168 
169  Scalar moleFraction_[numPhases][numComponents];
170  Scalar averageMolarMass_[numPhases];
171  Scalar sumMoleFractions_[numPhases];
172 };
173 
178 template <class Scalar,
179  class FluidSystem,
180  class Implementation>
182 {
183  enum { numPhases = FluidSystem::numPhases };
184 
185 public:
186  enum { numComponents = FluidSystem::numComponents };
187  static_assert(static_cast<int>(numPhases) == static_cast<int>(numComponents),
188  "The number of phases must be the same as the number of (pseudo-) components if you assume immiscibility");
189 
191  { }
192 
196  Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
197  { return (phaseIdx == compIdx)?1.0:0.0; }
198 
202  Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
203  { return (phaseIdx == compIdx)?1.0:0.0; }
204 
213  Scalar averageMolarMass(unsigned phaseIdx) const
214  { return FluidSystem::molarMass(/*compIdx=*/phaseIdx); }
215 
225  Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
226  { return asImp_().molarDensity(phaseIdx)*moleFraction(phaseIdx, compIdx); }
227 
232  template <class FluidState>
233  void assign(const FluidState& /* fs */)
234  { }
235 
244  void checkDefined() const
245  { }
246 
247 protected:
248  const Implementation& asImp_() const
249  { return *static_cast<const Implementation*>(this); }
250 };
251 
256 template <class Scalar>
258 {
259 public:
260  enum { numComponents = 0 };
261 
263  { }
264 
268  Scalar moleFraction(unsigned /* phaseIdx */, unsigned /* compIdx */) const
269  { OPM_THROW(std::logic_error, "Mole fractions are not provided by this fluid state"); }
270 
274  Scalar massFraction(unsigned /* phaseIdx */, unsigned /* compIdx */) const
275  { OPM_THROW(std::logic_error, "Mass fractions are not provided by this fluid state"); }
276 
285  Scalar averageMolarMass(unsigned /* phaseIdx */) const
286  { OPM_THROW(std::logic_error, "Mean molar masses are not provided by this fluid state"); }
287 
297  Scalar molarity(unsigned /* phaseIdx */, unsigned /* compIdx */) const
298  { OPM_THROW(std::logic_error, "Molarities are not provided by this fluid state"); }
299 
308  void checkDefined() const
309  { }
310 };
311 
312 
313 } // namespace Opm
314 
315 #endif
Scalar molarity(unsigned, unsigned) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:297
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:135
Module for the modular fluid state which does not store the compositions but throws std::logic_error ...
Definition: FluidStateCompositionModules.hpp:257
const Scalar & moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:69
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
Module for the modular fluid state which stores the phase compositions explicitly in terms of mole fr...
Definition: FluidStateCompositionModules.hpp:49
Scalar moleFraction(unsigned, unsigned) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:268
Scalar averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:213
const Scalar & averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:92
Definition: Air_Mesitylene.hpp:33
void assign(const FluidState &)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:233
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:196
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:308
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:158
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:225
Module for the modular fluid state which provides the phase compositions assuming immiscibility...
Definition: FluidStateCompositionModules.hpp:181
Scalar massFraction(unsigned, unsigned) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:274
Scalar averageMolarMass(unsigned) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:285
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:104
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:244
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:75
void setMoleFraction(unsigned phaseIdx, unsigned compIdx, const Scalar &value)
Set the mole fraction of a component in a phase [] and update the average molar mass [kg/mol] accordi...
Definition: FluidStateCompositionModules.hpp:112
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:202