My Project
Loading...
Searching...
No Matches
PengRobinsonMixture.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_MIXTURE_HPP
28#define OPM_PENG_ROBINSON_MIXTURE_HPP
29
30#include "PengRobinson.hpp"
31
33
34namespace Opm {
39template <class Scalar, class StaticParameters>
41{
42 enum { numComponents = StaticParameters::numComponents };
43
44 // this class cannot be instantiated!
45 PengRobinsonMixture() = default;
46
47 // the ideal gas constant
48 static const Scalar R;
49
50 // the u and w parameters as given by the Peng-Robinson EOS
51 static const Scalar u;
52 static const Scalar w;
53
54public:
55
73 template <class FluidState, class Params, class LhsEval = typename FluidState::Scalar>
74 static LhsEval computeFugacityCoefficient(const FluidState& fs,
75 const Params& params,
76 unsigned phaseIdx,
77 unsigned compIdx)
78 {
79 // note that we normalize the component mole fractions, so
80 // that their sum is 100%. This increases numerical stability
81 // considerably if the fluid state is not physical.
82 LhsEval Vm = params.molarVolume(phaseIdx);
83
84 // Calculate b_i / b
85 LhsEval bi_b = params.bPure(phaseIdx, compIdx) / params.b(phaseIdx);
86
87 // Calculate the compressibility factor
88 LhsEval RT = R*fs.temperature(phaseIdx);
89 LhsEval p = fs.pressure(phaseIdx); // molar volume in [bar]
90 LhsEval Z = p*Vm/RT; // compressibility factor
91
92 // Calculate A^* and B^* (see: Reid, p. 42)
93 LhsEval Astar = params.a(phaseIdx)*p/(RT*RT);
94 LhsEval Bstar = params.b(phaseIdx)*p/(RT);
95
96 LhsEval A_s = 0.0;
97 for (unsigned compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
98 A_s += params.aCache(phaseIdx, compIdx, compJIdx) * fs.moleFraction(phaseIdx, compJIdx) * p / (RT * RT);
99 }
100
101 LhsEval alpha;
102 LhsEval betta;
103 LhsEval gamma;
104 LhsEval ln_phi;
105 LhsEval fugCoeff;
106
107 Scalar m1;
108 Scalar m2;
109
110 m1 = 0.5*(u + std::sqrt(u*u - 4*w));
111 m2 = 0.5*(u - std::sqrt(u*u - 4*w));
112
113 alpha = -log(Z - Bstar) + bi_b * (Z - 1);
114 betta = log((Z + m2 * Bstar) / (Z + m1 * Bstar)) * Astar / ((m1 - m2) * Bstar);
115 gamma = (2 / Astar ) * A_s - bi_b;
116 ln_phi = alpha + (betta * gamma);
117
118 fugCoeff = exp(ln_phi);
119
121 // limit the fugacity coefficient to a reasonable range:
122 //
123 // on one side, we want the mole fraction to be at
124 // least 10^-3 if the fugacity is at the current pressure
125 //
126 fugCoeff = min(1e10, fugCoeff);
127 //
128 // on the other hand, if the mole fraction of the component is 100%, we want the
129 // fugacity to be at least 10^-3 Pa
130 //
131 fugCoeff = max(1e-10, fugCoeff);
133
134 return fugCoeff;
135 }
136
137};
138
139template <class Scalar, class StaticParameters>
140const Scalar PengRobinsonMixture<Scalar, StaticParameters>::R = Constants<Scalar>::R;
141template<class Scalar, class StaticParameters>
142const Scalar PengRobinsonMixture<Scalar, StaticParameters>::u = 2.0;
143template<class Scalar, class StaticParameters>
144const Scalar PengRobinsonMixture<Scalar, StaticParameters>::w = -1.0;
145
146} // namespace Opm
147
148#endif
A central place for various physical constants occuring in some equations.
Implements the Peng-Robinson equation of state for liquids and gases.
static const Scalar R
The ideal gas constant [J/(mol K)].
Definition Constants.hpp:45
Implements the Peng-Robinson equation of state for a mixture.
Definition PengRobinsonMixture.hpp:41
static LhsEval computeFugacityCoefficient(const FluidState &fs, const Params &params, unsigned phaseIdx, unsigned compIdx)
Returns the fugacity coefficient of an individual component in the phase.
Definition PengRobinsonMixture.hpp:74
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30