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 
128  void setCompressFactor(unsigned phaseIdx, const Scalar& value) {
129  Valgrind::CheckDefined(value);
130  Z_[phaseIdx] = value;
131  }
132 
133  Scalar compressFactor(unsigned phaseIdx) const {
134  return Z_[phaseIdx];
135  }
136 
141  template <class FluidState>
142  void assign(const FluidState& fs)
143  {
144  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
145  averageMolarMass_[phaseIdx] = 0;
146  sumMoleFractions_[phaseIdx] = 0;
147  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
148  moleFraction_[phaseIdx][compIdx] =
149  decay<Scalar>(fs.moleFraction(phaseIdx, compIdx));
150 
151  averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
152  sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
153  }
154  }
155  }
156 
165  void checkDefined() const
166  {
167  Valgrind::CheckDefined(moleFraction_);
168  Valgrind::CheckDefined(averageMolarMass_);
169  Valgrind::CheckDefined(sumMoleFractions_);
170  Valgrind::CheckDefined(K_);
171  Valgrind::CheckDefined(L_);
172  }
173 
174  const Scalar& K(unsigned compIdx) const
175  {
176  return K_[compIdx];
177  }
178 
182  void setKvalue(unsigned compIdx, const Scalar& value)
183  {
184  K_[compIdx] = value;
185  }
186 
190  const Scalar& L() const
191  {
192  return L_;
193  }
194 
198  void setLvalue(const Scalar& value)
199  {
200  L_ = value;
201  }
202 
207  Scalar wilsonK_(unsigned compIdx) const
208  {
209  const auto& acf = FluidSystem::acentricFactor(compIdx);
210  const auto& T_crit = FluidSystem::criticalTemperature(compIdx);
211  const auto& T = asImp_().temperature(0);
212  const auto& p_crit = FluidSystem::criticalPressure(compIdx);
213  const auto& p = asImp_().pressure(0); //for now assume no capillary pressure
214 
215  const auto tmp = exp(5.37 * (1+acf) * (1-T_crit/T)) * (p_crit/p);
216  return tmp;
217  }
218 
219 protected:
220  const Implementation& asImp_() const
221  {
222  return *static_cast<const Implementation*>(this);
223  }
224 
225  std::array<std::array<Scalar,numComponents>,numPhases> moleFraction_;
226  std::array<Scalar,numPhases> averageMolarMass_;
227  std::array<Scalar,numPhases> sumMoleFractions_;
228  std::array<Scalar,numPhases> Z_;
229  std::array<Scalar,numComponents> K_;
230  Scalar L_;
231 };
232 
237 template <class Scalar,
238  class FluidSystem,
239  class Implementation>
241 {
242  enum { numPhases = FluidSystem::numPhases };
243 
244 public:
245  enum { numComponents = FluidSystem::numComponents };
246  static_assert(static_cast<int>(numPhases) == static_cast<int>(numComponents),
247  "The number of phases must be the same as the number of (pseudo-) components if you assume immiscibility");
248 
250 
254  Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
255  { return (phaseIdx == compIdx)?1.0:0.0; }
256 
260  Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
261  { return (phaseIdx == compIdx)?1.0:0.0; }
262 
271  Scalar averageMolarMass(unsigned phaseIdx) const
272  { return FluidSystem::molarMass(/*compIdx=*/phaseIdx); }
273 
283  Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
284  { return asImp_().molarDensity(phaseIdx)*moleFraction(phaseIdx, compIdx); }
285 
290  template <class FluidState>
291  void assign(const FluidState& /* fs */)
292  { }
293 
302  void checkDefined() const
303  { }
304 
305 protected:
306  const Implementation& asImp_() const
307  { return *static_cast<const Implementation*>(this); }
308 };
309 
314 template <class Scalar>
316 {
317 public:
318  enum { numComponents = 0 };
319 
321 
325  Scalar moleFraction(unsigned /* phaseIdx */, unsigned /* compIdx */) const
326  { throw std::logic_error("Mole fractions are not provided by this fluid state"); }
327 
331  Scalar massFraction(unsigned /* phaseIdx */, unsigned /* compIdx */) const
332  { throw std::logic_error("Mass fractions are not provided by this fluid state"); }
333 
342  Scalar averageMolarMass(unsigned /* phaseIdx */) const
343  { throw std::logic_error("Mean molar masses are not provided by this fluid state"); }
344 
354  Scalar molarity(unsigned /* phaseIdx */, unsigned /* compIdx */) const
355  { throw std::logic_error("Molarities are not provided by this fluid state"); }
356 
365  void checkDefined() const
366  { }
367 };
368 
369 
370 } // namespace Opm
371 
372 #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
void setLvalue(const Scalar &value)
Set the L value [-].
Definition: FluidStateCompositionModules.hpp:198
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
Scalar wilsonK_(unsigned compIdx) const
Wilson formula to calculate K.
Definition: FluidStateCompositionModules.hpp:207
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:142
const Scalar & L() const
The L value of a composition [-].
Definition: FluidStateCompositionModules.hpp:190
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
void setKvalue(unsigned compIdx, const Scalar &value)
Set the K value of a component [-].
Definition: FluidStateCompositionModules.hpp:182
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:165
Module for the modular fluid state which provides the phase compositions assuming immiscibility.
Definition: FluidStateCompositionModules.hpp:241
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:260
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:283
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:302
void assign(const FluidState &)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:291
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:254
Scalar averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:271
Module for the modular fluid state which does not store the compositions but throws std::logic_error ...
Definition: FluidStateCompositionModules.hpp:316
Scalar moleFraction(unsigned, unsigned) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:325
Scalar averageMolarMass(unsigned) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:342
Scalar molarity(unsigned, unsigned) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:354
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:365
Scalar massFraction(unsigned, unsigned) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:331