27 #ifndef OPM_PENG_ROBINSON_PARAMS_MIXTURE_HPP 28 #define OPM_PENG_ROBINSON_PARAMS_MIXTURE_HPP 57 template <
class Scalar,
class Flu
idSystem,
unsigned phaseIdx,
bool useSpe5Relations=false>
61 enum { numComponents = FluidSystem::numComponents };
69 static const Scalar R;
75 template <
class Flu
idState>
79 fluidState.pressure(phaseIdx));
89 Valgrind::CheckDefined(temperature);
90 Valgrind::CheckDefined(pressure);
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);
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));
110 f_omega = 0.37464 + omega*(1.54226 - omega*0.26992);
112 Valgrind::CheckDefined(f_omega);
114 Scalar tmp = 1 + f_omega*(1 - Opm::sqrt(Tr));
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)));
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());
138 template <
class Flu
idState>
142 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
143 sumx += fs.moleFraction(phaseIdx, compIdx);
144 sumx = std::max(Scalar(1e-10), sumx);
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);
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);
163 newA += xi * xj * aCache_[compIIdx][compJIdx];
165 assert(std::isfinite(Opm::scalarValue(newA)));
169 newB += Opm::max(0.0, xi) * this->pureParams_[compIIdx].
b();
170 assert(std::isfinite(Opm::scalarValue(newB)));
177 Valgrind::CheckDefined(this->
a());
178 Valgrind::CheckDefined(this->
b());
190 template <
class Flu
idState>
201 {
return pureParams_[compIdx]; }
208 assert(0 <= compIdx && compIdx < numComponents);
209 return pureParams_[compIdx];
219 for (
unsigned i = 0; i < numComponents; ++i)
222 Valgrind::CheckDefined(this->
a());
223 Valgrind::CheckDefined(this->
b());
228 PureParams pureParams_[numComponents];
233 for (
unsigned compIIdx = 0; compIIdx < numComponents; ++ compIIdx) {
234 for (
unsigned compJIdx = 0; compJIdx < numComponents; ++ compJIdx) {
236 Scalar Psi = FluidSystem::interactionCoefficient(compIIdx, compJIdx);
238 aCache_[compIIdx][compJIdx] =
239 Opm::sqrt(this->pureParams_[compIIdx].
a()
240 * this->pureParams_[compJIdx].
a())
246 Scalar aCache_[numComponents][numComponents];
249 template <
class Scalar,
class Flu
idSystem,
unsigned phaseIdx,
bool useSpe5Relations>
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 'a' of the Peng-Robinson fluid.
Definition: PengRobinsonParams.hpp:76
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 'b' 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
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 'a' 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 'b' 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