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 };
73 Scalar
getaCache(
unsigned compIIdx,
unsigned compJIdx )
const
75 return aCache_[compIIdx][compJIdx];
81 template <
class Flu
idState>
85 fluidState.pressure(phaseIdx));
95 Valgrind::CheckDefined(temperature);
96 Valgrind::CheckDefined(pressure);
103 for (
unsigned i = 0; i < numComponents; ++i) {
104 Scalar pc = FluidSystem::criticalPressure(i);
105 Scalar omega = FluidSystem::acentricFactor(i);
106 Scalar Tr = temperature/FluidSystem::criticalTemperature(i);
112 f_omega = 0.37464 + omega*(1.54226 + omega*(-0.26992));
114 f_omega = 0.379642 + omega*(1.48503 + omega*(-0.164423 + omega*0.016666));
116 Valgrind::CheckDefined(f_omega);
118 Scalar tmp = 1 + f_omega*(1 - sqrt(Tr));
121 Scalar newA = 0.4572355*RTc*RTc/pc * tmp;
122 Scalar newB = 0.0777961 * RTc / pc;
123 assert(std::isfinite(scalarValue(newA)));
124 assert(std::isfinite(scalarValue(newB)));
126 this->pureParams_[i].
setA(newA);
127 this->pureParams_[i].
setB(newB);
128 Valgrind::CheckDefined(this->pureParams_[i].
a());
129 Valgrind::CheckDefined(this->pureParams_[i].
b());
142 template <
class Flu
idState>
145 using FlashEval =
typename FluidState::Scalar;
147 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
148 sumx += fs.moleFraction(phaseIdx, compIdx);
149 sumx = std::max(Scalar(1e-10), sumx);
157 for (
unsigned compIIdx = 0; compIIdx < numComponents; ++compIIdx) {
158 const FlashEval moleFracI = fs.moleFraction(phaseIdx, compIIdx);
159 FlashEval xi = max(0.0, min(1.0, moleFracI));
160 Valgrind::CheckDefined(xi);
162 for (
unsigned compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
163 const FlashEval moleFracJ = fs.moleFraction(phaseIdx, compJIdx );
164 FlashEval xj = max(0.0, min(1.0, moleFracJ));
165 Valgrind::CheckDefined(xj);
168 newA += xi * xj * aCache_[compIIdx][compJIdx];
170 assert(std::isfinite(scalarValue(newA)));
174 newB += max(0.0, xi) * this->pureParams_[compIIdx].
b();
175 assert(std::isfinite(scalarValue(newB)));
182 Valgrind::CheckDefined(this->
a());
183 Valgrind::CheckDefined(this->
b());
195 template <
class Flu
idState>
206 {
return pureParams_[compIdx]; }
213 assert(0 <= compIdx && compIdx < numComponents);
214 return pureParams_[compIdx];
224 for (
unsigned i = 0; i < numComponents; ++i)
227 Valgrind::CheckDefined(this->
a());
228 Valgrind::CheckDefined(this->
b());
233 PureParams pureParams_[numComponents];
238 for (
unsigned compIIdx = 0; compIIdx < numComponents; ++ compIIdx) {
239 for (
unsigned compJIdx = 0; compJIdx < numComponents; ++ compJIdx) {
241 Scalar Psi = FluidSystem::interactionCoefficient(compIIdx, compJIdx);
243 aCache_[compIIdx][compJIdx] =
244 sqrt(this->pureParams_[compIIdx].
a()
245 * this->pureParams_[compJIdx].
a())
251 Scalar aCache_[numComponents][numComponents];
A central place for various physical constants occuring in some equations.
Stores and provides access to the Peng-Robinson parameters.
A central place for various physical constants occuring in some equations.
Definition: Constants.hpp:41
The mixing rule for the oil and the gas phases of the SPE5 problem.
Definition: PengRobinsonParamsMixture.hpp:60
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:196
void updatePure(Scalar temperature, Scalar pressure)
Peng-Robinson parameters for the pure components.
Definition: PengRobinsonParamsMixture.hpp:93
void updatePure(const FluidState &fluidState)
Update Peng-Robinson parameters for the pure components.
Definition: PengRobinsonParamsMixture.hpp:82
const PureParams & operator[](unsigned compIdx) const
Returns the Peng-Robinson parameters for a pure component.
Definition: PengRobinsonParamsMixture.hpp:211
Scalar getaCache(unsigned compIIdx, unsigned compJIdx) const
TODO.
Definition: PengRobinsonParamsMixture.hpp:73
void checkDefined() const
If run under valgrind, this method produces an warning if the parameters where not determined correct...
Definition: PengRobinsonParamsMixture.hpp:221
void updateMix(const FluidState &fs)
Calculates the "a" and "b" Peng-Robinson parameters for the mixture.
Definition: PengRobinsonParamsMixture.hpp:143
const PureParams & pureParams(unsigned compIdx) const
Return the Peng-Robinson parameters of a pure substance,.
Definition: PengRobinsonParamsMixture.hpp:205
Stores and provides access to the Peng-Robinson parameters.
Definition: PengRobinsonParams.hpp:44
Scalar a() const
Returns the attractive parameter 'a' of the Peng-Robinson fluid.
Definition: PengRobinsonParams.hpp:50
void setB(Scalar value)
Set the repulsive parameter 'b' of the Peng-Robinson fluid.
Definition: PengRobinsonParams.hpp:83
Scalar b() const
Returns the repulsive parameter 'b' of the Peng-Robinson fluid.
Definition: PengRobinsonParams.hpp:57
void setA(Scalar value)
Set the attractive parameter 'a' of the Peng-Robinson fluid.
Definition: PengRobinsonParams.hpp:76