My Project
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 
33 
34 #include <algorithm>
35 #include <array>
36 #include <cmath>
37 
38 namespace Opm {
39 
44 template <class Scalar,
45  class FluidSystem,
46  class Implementation>
48 {
49  enum { numPhases = FluidSystem::numPhases };
50  enum { numComponents = FluidSystem::numComponents };
51 
52 public:
54  {
55  for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
56  for (int compIdx = 0; compIdx < numComponents; ++compIdx)
57  moleFraction_[phaseIdx][compIdx] = 0.0;
58 
59  Valgrind::SetDefined(moleFraction_);
60  Valgrind::SetUndefined(averageMolarMass_);
61  Valgrind::SetUndefined(sumMoleFractions_);
62  }
63 
67  const Scalar& moleFraction(unsigned phaseIdx, unsigned compIdx) const
68  { return moleFraction_[phaseIdx][compIdx]; }
69 
73  Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
74  {
75  return
76  abs(sumMoleFractions_[phaseIdx])
77  *moleFraction_[phaseIdx][compIdx]
78  *FluidSystem::molarMass(compIdx)
79  / max(1e-40, abs(averageMolarMass_[phaseIdx]));
80  }
81 
90  const Scalar& averageMolarMass(unsigned phaseIdx) const
91  { return averageMolarMass_[phaseIdx]; }
92 
102  Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
103  { return asImp_().molarDensity(phaseIdx)*moleFraction(phaseIdx, compIdx); }
104 
110  void setMoleFraction(unsigned phaseIdx, unsigned compIdx, const Scalar& value)
111  {
112  Valgrind::CheckDefined(value);
113  Valgrind::SetUndefined(sumMoleFractions_[phaseIdx]);
114  Valgrind::SetUndefined(averageMolarMass_[phaseIdx]);
115  Valgrind::SetUndefined(moleFraction_[phaseIdx][compIdx]);
116 
117  moleFraction_[phaseIdx][compIdx] = value;
118 
119  // re-calculate the mean molar mass
120  sumMoleFractions_[phaseIdx] = 0.0;
121  averageMolarMass_[phaseIdx] = 0.0;
122  for (unsigned compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
123  sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
124  averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
125  }
126  }
127 
132  template <class FluidState>
133  void assign(const FluidState& fs)
134  {
135  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
136  averageMolarMass_[phaseIdx] = 0;
137  sumMoleFractions_[phaseIdx] = 0;
138  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
139  moleFraction_[phaseIdx][compIdx] =
140  decay<Scalar>(fs.moleFraction(phaseIdx, compIdx));
141 
142  averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
143  sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
144  }
145  }
146  }
147 
156  void checkDefined() const
157  {
158  Valgrind::CheckDefined(moleFraction_);
159  Valgrind::CheckDefined(averageMolarMass_);
160  Valgrind::CheckDefined(sumMoleFractions_);
161  }
162 
163 protected:
164  const Implementation& asImp_() const
165  { return *static_cast<const Implementation*>(this); }
166 
167  std::array<std::array<Scalar,numComponents>,numPhases> moleFraction_;
168  std::array<Scalar,numPhases> averageMolarMass_;
169  std::array<Scalar,numPhases> sumMoleFractions_;
170 };
171 
176 template <class Scalar,
177  class FluidSystem,
178  class Implementation>
180 {
181  enum { numPhases = FluidSystem::numPhases };
182 
183 public:
184  enum { numComponents = FluidSystem::numComponents };
185  static_assert(static_cast<int>(numPhases) == static_cast<int>(numComponents),
186  "The number of phases must be the same as the number of (pseudo-) components if you assume immiscibility");
187 
189  { }
190 
194  Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
195  { return (phaseIdx == compIdx)?1.0:0.0; }
196 
200  Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
201  { return (phaseIdx == compIdx)?1.0:0.0; }
202 
211  Scalar averageMolarMass(unsigned phaseIdx) const
212  { return FluidSystem::molarMass(/*compIdx=*/phaseIdx); }
213 
223  Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
224  { return asImp_().molarDensity(phaseIdx)*moleFraction(phaseIdx, compIdx); }
225 
230  template <class FluidState>
231  void assign(const FluidState& /* fs */)
232  { }
233 
242  void checkDefined() const
243  { }
244 
245 protected:
246  const Implementation& asImp_() const
247  { return *static_cast<const Implementation*>(this); }
248 };
249 
254 template <class Scalar>
256 {
257 public:
258  enum { numComponents = 0 };
259 
261  { }
262 
266  Scalar moleFraction(unsigned /* phaseIdx */, unsigned /* compIdx */) const
267  { throw std::logic_error("Mole fractions are not provided by this fluid state"); }
268 
272  Scalar massFraction(unsigned /* phaseIdx */, unsigned /* compIdx */) const
273  { throw std::logic_error("Mass fractions are not provided by this fluid state"); }
274 
283  Scalar averageMolarMass(unsigned /* phaseIdx */) const
284  { throw std::logic_error("Mean molar masses are not provided by this fluid state"); }
285 
295  Scalar molarity(unsigned /* phaseIdx */, unsigned /* compIdx */) const
296  { throw std::logic_error("Molarities are not provided by this fluid state"); }
297 
306  void checkDefined() const
307  { }
308 };
309 
310 
311 } // namespace Opm
312 
313 #endif
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
Some templates to wrap the valgrind client request macros.
Module for the modular fluid state which stores the phase compositions explicitly in terms of mole fr...
Definition: FluidStateCompositionModules.hpp:48
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:102
const Scalar & moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:67
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:133
const Scalar & averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:90
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:110
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:73
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:156
Module for the modular fluid state which provides the phase compositions assuming immiscibility.
Definition: FluidStateCompositionModules.hpp:180
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:200
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:223
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:242
void assign(const FluidState &)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:231
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:194
Scalar averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:211
Module for the modular fluid state which does not store the compositions but throws std::logic_error ...
Definition: FluidStateCompositionModules.hpp:256
Scalar moleFraction(unsigned, unsigned) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:266
Scalar averageMolarMass(unsigned) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:283
Scalar molarity(unsigned, unsigned) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:295
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:306
Scalar massFraction(unsigned, unsigned) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:272