PengRobinsonParamsMixture.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_PENG_ROBINSON_PARAMS_MIXTURE_HPP
28 #define OPM_PENG_ROBINSON_PARAMS_MIXTURE_HPP
29 
30 #include "PengRobinsonParams.hpp"
31 
34 
35 #include <algorithm>
36 
37 namespace Opm
38 {
39 
57 template <class Scalar, class FluidSystem, unsigned phaseIdx, bool useSpe5Relations=false>
59  : public PengRobinsonParams<Scalar>
60 {
61  enum { numComponents = FluidSystem::numComponents };
62 
63  // Peng-Robinson parameters for pure substances
65 
67 
68  // the ideal gas constant
69  static const Scalar R;
70 
71 public:
75  template <class FluidState>
76  void updatePure(const FluidState& fluidState)
77  {
78  updatePure(fluidState.temperature(phaseIdx),
79  fluidState.pressure(phaseIdx));
80  }
81 
87  void updatePure(Scalar temperature, Scalar pressure)
88  {
89  Valgrind::CheckDefined(temperature);
90  Valgrind::CheckDefined(pressure);
91 
92  // Calculate the Peng-Robinson parameters of the pure
93  // components
94  //
95  // See: R. Reid, et al.: The Properties of Gases and Liquids,
96  // 4th edition, McGraw-Hill, 1987, p. 43
97  for (unsigned i = 0; i < numComponents; ++i) {
98  Scalar pc = FluidSystem::criticalPressure(i);
99  Scalar omega = FluidSystem::acentricFactor(i);
100  Scalar Tr = temperature/FluidSystem::criticalTemperature(i);
101  Scalar RTc = R*FluidSystem::criticalTemperature(i);
102 
103  Scalar f_omega;
104 
105  if (useSpe5Relations) {
106  if (omega < 0.49) f_omega = 0.37464 + omega*(1.54226 + omega*(-0.26992));
107  else f_omega = 0.379642 + omega*(1.48503 + omega*(-0.164423 + omega*0.016666));
108  }
109  else
110  f_omega = 0.37464 + omega*(1.54226 - omega*0.26992);
111 
112  Valgrind::CheckDefined(f_omega);
113 
114  Scalar tmp = 1 + f_omega*(1 - Opm::sqrt(Tr));
115  tmp = tmp*tmp;
116 
117  Scalar newA = 0.4572355*RTc*RTc/pc * tmp;
118  Scalar newB = 0.0777961 * RTc / pc;
119  assert(std::isfinite(Opm::scalarValue(newA)));
120  assert(std::isfinite(Opm::scalarValue(newB)));
121 
122  this->pureParams_[i].setA(newA);
123  this->pureParams_[i].setB(newB);
124  Valgrind::CheckDefined(this->pureParams_[i].a());
125  Valgrind::CheckDefined(this->pureParams_[i].b());
126  }
127 
128  updateACache_();
129  }
130 
138  template <class FluidState>
139  void updateMix(const FluidState& fs)
140  {
141  Scalar sumx = 0.0;
142  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
143  sumx += fs.moleFraction(phaseIdx, compIdx);
144  sumx = std::max(Scalar(1e-10), sumx);
145 
146  // Calculate the Peng-Robinson parameters of the mixture
147  //
148  // See: R. Reid, et al.: The Properties of Gases and Liquids,
149  // 4th edition, McGraw-Hill, 1987, p. 82
150  Scalar newA = 0;
151  Scalar newB = 0;
152  for (unsigned compIIdx = 0; compIIdx < numComponents; ++compIIdx) {
153  const Scalar moleFracI = fs.moleFraction(phaseIdx, compIIdx);
154  Scalar xi = Opm::max(0.0, Opm::min(1.0, moleFracI));
155  Valgrind::CheckDefined(xi);
156 
157  for (unsigned compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
158  const Scalar moleFracJ = fs.moleFraction(phaseIdx, compJIdx );
159  Scalar xj = Opm::max(0.0, Opm::min(1.0, moleFracJ));
160  Valgrind::CheckDefined(xj);
161 
162  // mixing rule from Reid, page 82
163  newA += xi * xj * aCache_[compIIdx][compJIdx];
164 
165  assert(std::isfinite(Opm::scalarValue(newA)));
166  }
167 
168  // mixing rule from Reid, page 82
169  newB += Opm::max(0.0, xi) * this->pureParams_[compIIdx].b();
170  assert(std::isfinite(Opm::scalarValue(newB)));
171  }
172 
173  // assert(newB > 0);
174  this->setA(newA);
175  this->setB(newB);
176 
177  Valgrind::CheckDefined(this->a());
178  Valgrind::CheckDefined(this->b());
179 
180  }
181 
190  template <class FluidState>
191  void updateSingleMoleFraction(const FluidState& fs,
192  unsigned /*compIdx*/)
193  {
194  updateMix(fs);
195  }
196 
200  const PureParams& pureParams(unsigned compIdx) const
201  { return pureParams_[compIdx]; }
202 
206  const PureParams& operator[](unsigned compIdx) const
207  {
208  assert(0 <= compIdx && compIdx < numComponents);
209  return pureParams_[compIdx];
210  }
211 
216  void checkDefined() const
217  {
218 #ifndef NDEBUG
219  for (unsigned i = 0; i < numComponents; ++i)
220  pureParams_[i].checkDefined();
221 
222  Valgrind::CheckDefined(this->a());
223  Valgrind::CheckDefined(this->b());
224 #endif
225  }
226 
227 protected:
228  PureParams pureParams_[numComponents];
229 
230 private:
231  void updateACache_()
232  {
233  for (unsigned compIIdx = 0; compIIdx < numComponents; ++ compIIdx) {
234  for (unsigned compJIdx = 0; compJIdx < numComponents; ++ compJIdx) {
235  // interaction coefficient as given in SPE5
236  Scalar Psi = FluidSystem::interactionCoefficient(compIIdx, compJIdx);
237 
238  aCache_[compIIdx][compJIdx] =
239  Opm::sqrt(this->pureParams_[compIIdx].a()
240  * this->pureParams_[compJIdx].a())
241  * (1 - Psi);
242  }
243  }
244  }
245 
246  Scalar aCache_[numComponents][numComponents];
247 };
248 
249 template <class Scalar, class FluidSystem, unsigned phaseIdx, bool useSpe5Relations>
250 const Scalar PengRobinsonParamsMixture<Scalar, FluidSystem, phaseIdx, useSpe5Relations>::R = Opm::Constants<Scalar>::R;
251 
252 } // namespace Opm
253 
254 #endif
The mixing rule for the oil and the gas phases of the SPE5 problem.
Definition: PengRobinsonParamsMixture.hpp:58
void setA(Scalar value)
Set the attractive parameter &#39;a&#39; of the Peng-Robinson fluid.
Definition: PengRobinsonParams.hpp:76
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
const PureParams & operator[](unsigned compIdx) const
Returns the Peng-Robinson parameters for a pure component.
Definition: PengRobinsonParamsMixture.hpp:206
void setB(Scalar value)
Set the repulsive parameter &#39;b&#39; of the Peng-Robinson fluid.
Definition: PengRobinsonParams.hpp:83
Stores and provides access to the Peng-Robinson parameters.
void updateMix(const FluidState &fs)
Calculates the "a" and "b" Peng-Robinson parameters for the mixture.
Definition: PengRobinsonParamsMixture.hpp:139
Definition: Air_Mesitylene.hpp:33
A central place for various physical constants occuring in some equations.
Definition: Constants.hpp:40
void updatePure(const FluidState &fluidState)
Update Peng-Robinson parameters for the pure components.
Definition: PengRobinsonParamsMixture.hpp:76
void updatePure(Scalar temperature, Scalar pressure)
Peng-Robinson parameters for the pure components.
Definition: PengRobinsonParamsMixture.hpp:87
Definition: MathToolbox.hpp:48
void updateSingleMoleFraction(const FluidState &fs, unsigned)
Calculates the "a" and "b" Peng-Robinson parameters for the mixture provided that only a single mole ...
Definition: PengRobinsonParamsMixture.hpp:191
A central place for various physical constants occuring in some equations.
Scalar a() const
Returns the attractive parameter &#39;a&#39; of the Peng-Robinson fluid.
Definition: PengRobinsonParams.hpp:50
Stores and provides access to the Peng-Robinson parameters.
Definition: PengRobinsonParams.hpp:43
Scalar b() const
Returns the repulsive parameter &#39;b&#39; of the Peng-Robinson fluid.
Definition: PengRobinsonParams.hpp:57
void checkDefined() const
If run under valgrind, this method produces an warning if the parameters where not determined correct...
Definition: PengRobinsonParamsMixture.hpp:216
const PureParams & pureParams(unsigned compIdx) const
Return the Peng-Robinson parameters of a pure substance,.
Definition: PengRobinsonParamsMixture.hpp:200