| Line | Branch | Exec | Source |
|---|---|---|---|
| 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 | // SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder | ||
| 5 | // SPDX-License-Identifier: GPL-3.0-or-later | ||
| 6 | // | ||
| 7 | /*! | ||
| 8 | * \file | ||
| 9 | * \ingroup EOS | ||
| 10 | * \brief The mixing rule for the oil and the gas phases of the SPE5 problem. | ||
| 11 | * | ||
| 12 | * This problem comprises \f$H_2O\f$, \f$C_1\f$, \f$C_3\f$, \f$C_6\f$, | ||
| 13 | * \f$C_10\f$, \f$C_15\f$ and \f$C_20\f$ as components. | ||
| 14 | * | ||
| 15 | * See: | ||
| 16 | * | ||
| 17 | * R. Reid, et al. (1987, pp. 43-44) \cite reid1987 <BR> | ||
| 18 | * | ||
| 19 | * and | ||
| 20 | * | ||
| 21 | * J.E. Killough, et al. (1987) \cite SPE5 | ||
| 22 | */ | ||
| 23 | #ifndef DUMUX_PENG_ROBINSON_PARAMS_MIXTURE_HH | ||
| 24 | #define DUMUX_PENG_ROBINSON_PARAMS_MIXTURE_HH | ||
| 25 | |||
| 26 | #include <algorithm> | ||
| 27 | #include <dumux/material/constants.hh> | ||
| 28 | |||
| 29 | #include "pengrobinsonparams.hh" | ||
| 30 | |||
| 31 | namespace Dumux { | ||
| 32 | |||
| 33 | /*! | ||
| 34 | * \ingroup EOS | ||
| 35 | * \brief The mixing rule for the oil and the gas phases of the SPE5 problem. | ||
| 36 | * | ||
| 37 | * This problem comprises \f$H_2O\f$, \f$C_1\f$, \f$C_3\f$, \f$C_6\f$, | ||
| 38 | * \f$C_10\f$, \f$C_15\f$ and \f$C_20\f$ as components. | ||
| 39 | * | ||
| 40 | * See: | ||
| 41 | * | ||
| 42 | * R. Reid, et al. (1987, pp. 43-44) \cite reid1987 <BR> | ||
| 43 | * | ||
| 44 | * and | ||
| 45 | * | ||
| 46 | * J.E. Killough, et al. (1987) \cite SPE5 | ||
| 47 | */ | ||
| 48 | template <class Scalar, class FluidSystem, int phaseIdx, bool useSpe5Relations=false> | ||
| 49 | class PengRobinsonParamsMixture | ||
| 50 | : public PengRobinsonParams<Scalar> | ||
| 51 | { | ||
| 52 | enum { numComponents = FluidSystem::numComponents }; | ||
| 53 | |||
| 54 | // Peng-Robinson parameters for pure substances | ||
| 55 | using PureParams = PengRobinsonParams<Scalar>; | ||
| 56 | |||
| 57 | public: | ||
| 58 | /*! | ||
| 59 | * \brief Update Peng-Robinson parameters for the pure components. | ||
| 60 | * \param fluidState Thermodynamic state of the fluids | ||
| 61 | * | ||
| 62 | */ | ||
| 63 | template <class FluidState> | ||
| 64 | void updatePure(const FluidState &fluidState) | ||
| 65 | { | ||
| 66 | updatePure(fluidState.temperature(phaseIdx), | ||
| 67 | fluidState.pressure(phaseIdx)); | ||
| 68 | } | ||
| 69 | |||
| 70 | /*! | ||
| 71 | * \brief Peng-Robinson parameters for the pure components. | ||
| 72 | * \param temperature Temperature in \f$\mathrm{[K]}\f$ | ||
| 73 | * \param pressure pressure in \f$\mathrm{[Pa]}\f$ | ||
| 74 | * This method is given by the SPE5 paper. | ||
| 75 | */ | ||
| 76 | 29490 | void updatePure(Scalar temperature, Scalar pressure) | |
| 77 | { | ||
| 78 | // Calculate the Peng-Robinson parameters of the pure | ||
| 79 | // components | ||
| 80 | // | ||
| 81 | // See: R. Reid, et al.: The Properties of Gases and Liquids, | ||
| 82 | // 4th edition, McGraw-Hill, 1987, p. 43 | ||
| 83 |
2/2✓ Branch 0 taken 103215 times.
✓ Branch 1 taken 14745 times.
|
235920 | for (int i = 0; i < numComponents; ++i) { |
| 84 |
2/2✓ Branch 0 taken 73725 times.
✓ Branch 1 taken 29490 times.
|
206430 | Scalar pc = FluidSystem::criticalPressure(i); |
| 85 |
2/2✓ Branch 0 taken 73725 times.
✓ Branch 1 taken 29490 times.
|
206430 | Scalar omega = FluidSystem::acentricFactor(i); |
| 86 |
2/2✓ Branch 0 taken 73725 times.
✓ Branch 1 taken 29490 times.
|
206430 | Scalar Tr = temperature/FluidSystem::criticalTemperature(i); |
| 87 | 206430 | Scalar RTc = Constants<Scalar>::R*FluidSystem::criticalTemperature(i); | |
| 88 | |||
| 89 | Scalar f_omega; | ||
| 90 | |||
| 91 | if (useSpe5Relations) { | ||
| 92 |
2/2✓ Branch 0 taken 73725 times.
✓ Branch 1 taken 29490 times.
|
206430 | if (omega < 0.49) f_omega = 0.37464 + omega*(1.54226 + omega*(-0.26992)); |
| 93 | 58980 | else f_omega = 0.379642 + omega*(1.48503 + omega*(-0.164423 + omega*0.016666)); | |
| 94 | } | ||
| 95 | else | ||
| 96 | f_omega = 0.37464 + omega*(1.54226 - omega*0.26992); | ||
| 97 | |||
| 98 | using std::sqrt; | ||
| 99 | 206430 | Scalar tmp = 1 + f_omega*(1 - sqrt(Tr)); | |
| 100 | 206430 | tmp = tmp*tmp; | |
| 101 | |||
| 102 | 206430 | Scalar a = 0.4572355*RTc*RTc/pc * tmp; | |
| 103 | 206430 | Scalar b = 0.0777961 * RTc / pc; | |
| 104 | using std::isfinite; | ||
| 105 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 103215 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 103215 times.
|
412860 | assert(isfinite(a)); |
| 106 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 103215 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 103215 times.
|
412860 | assert(isfinite(b)); |
| 107 | |||
| 108 | 206430 | this->pureParams_[i].setA(a); | |
| 109 | 412860 | this->pureParams_[i].setB(b); | |
| 110 | } | ||
| 111 | |||
| 112 | 29490 | updateACache_(); | |
| 113 | 29490 | } | |
| 114 | |||
| 115 | /*! | ||
| 116 | * \brief Calculates the "a" and "b" Peng-Robinson parameters for | ||
| 117 | * the mixture. | ||
| 118 | * | ||
| 119 | * The updatePure() method needs to be called _before_ calling | ||
| 120 | * this method! | ||
| 121 | * \param fs the thermodynamic state of the fluids | ||
| 122 | */ | ||
| 123 | template <class FluidState> | ||
| 124 | 76018 | void updateMix(const FluidState &fs) | |
| 125 | { | ||
| 126 | 76018 | Scalar sumx = 0.0; | |
| 127 |
2/2✓ Branch 0 taken 266063 times.
✓ Branch 1 taken 38009 times.
|
608144 | for (int compIdx = 0; compIdx < numComponents; ++compIdx) |
| 128 | 1063972 | sumx += fs.moleFraction(phaseIdx, compIdx); | |
| 129 | using std::max; | ||
| 130 |
1/2✓ Branch 0 taken 38009 times.
✗ Branch 1 not taken.
|
76018 | sumx = max(1e-10, sumx); |
| 131 | |||
| 132 | // Calculate the Peng-Robinson parameters of the mixture | ||
| 133 | // | ||
| 134 | // See: R. Reid, et al.: The Properties of Gases and Liquids, | ||
| 135 | // 4th edition, McGraw-Hill, 1987, p. 82 | ||
| 136 | 76018 | Scalar a = 0; | |
| 137 | 76018 | Scalar b = 0; | |
| 138 |
2/2✓ Branch 0 taken 266063 times.
✓ Branch 1 taken 38009 times.
|
608144 | for (int compIIdx = 0; compIIdx < numComponents; ++compIIdx) { |
| 139 | 532126 | Scalar xi = fs.moleFraction(phaseIdx, compIIdx) / sumx; | |
| 140 | using::std::isfinite; | ||
| 141 |
2/2✓ Branch 0 taken 1862441 times.
✓ Branch 1 taken 266063 times.
|
4257008 | for (int compJIdx = 0; compJIdx < numComponents; ++compJIdx) { |
| 142 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1861461 times.
|
3724882 | Scalar xj = fs.moleFraction(phaseIdx, compJIdx) / sumx; |
| 143 | |||
| 144 | // mixing rule from Reid, page 82 | ||
| 145 | 3724882 | a += xi * xj * aCache_[compIIdx][compJIdx]; | |
| 146 | |||
| 147 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 1862441 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1862441 times.
|
7449764 | assert(isfinite(a)); |
| 148 | } | ||
| 149 | |||
| 150 | // mixing rule from Reid, page 82 | ||
| 151 | 532126 | b += xi * this->pureParams_[compIIdx].b(); | |
| 152 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 266063 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 266063 times.
|
1064252 | assert(isfinite(b)); |
| 153 | } | ||
| 154 | |||
| 155 | 152036 | this->setA(a); | |
| 156 | 152036 | this->setB(b); | |
| 157 | 76018 | } | |
| 158 | |||
| 159 | /*! | ||
| 160 | * \brief Calculates the "a" and "b" Peng-Robinson parameters for | ||
| 161 | * the mixture provided that only a single mole fraction | ||
| 162 | * was changed. | ||
| 163 | * | ||
| 164 | * The updatePure() method needs to be called _before_ calling | ||
| 165 | * this method! | ||
| 166 | * \param fs the thermodynamic state of the fluids | ||
| 167 | * \param compIdx the component index | ||
| 168 | */ | ||
| 169 | template <class FluidState> | ||
| 170 | ✗ | void updateSingleMoleFraction(const FluidState &fs, | |
| 171 | int compIdx) | ||
| 172 | { | ||
| 173 | ✗ | updateMix(fs); | |
| 174 | ✗ | } | |
| 175 | |||
| 176 | /*! | ||
| 177 | * \brief Return the Peng-Robinson parameters of a pure substance, | ||
| 178 | * \param compIdx the component index | ||
| 179 | */ | ||
| 180 | const PureParams &pureParams(int compIdx) const | ||
| 181 |
7/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 12 times.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 12 times.
✓ Branch 7 taken 2 times.
|
2387000 | { return pureParams_[compIdx]; } |
| 182 | |||
| 183 | /*! | ||
| 184 | * \brief Returns the Peng-Robinson parameters for a pure component. | ||
| 185 | * \param compIdx the component index | ||
| 186 | */ | ||
| 187 | const PureParams &operator[](int compIdx) const | ||
| 188 | { | ||
| 189 | assert(0 <= compIdx && compIdx < numComponents); | ||
| 190 | return pureParams_[compIdx]; | ||
| 191 | } | ||
| 192 | |||
| 193 | protected: | ||
| 194 | PureParams pureParams_[numComponents]; | ||
| 195 | |||
| 196 | private: | ||
| 197 | 29490 | void updateACache_() | |
| 198 | { | ||
| 199 |
2/2✓ Branch 0 taken 103215 times.
✓ Branch 1 taken 14745 times.
|
235920 | for (int compIIdx = 0; compIIdx < numComponents; ++ compIIdx) { |
| 200 |
2/2✓ Branch 0 taken 722505 times.
✓ Branch 1 taken 103215 times.
|
1651440 | for (int compJIdx = 0; compJIdx < numComponents; ++ compJIdx) { |
| 201 | // interaction coefficient as given in SPE5 | ||
| 202 |
2/2✓ Branch 0 taken 309645 times.
✓ Branch 1 taken 412860 times.
|
1445010 | Scalar Psi = FluidSystem::interactionCoefficient(compIIdx, compJIdx); |
| 203 | using std::sqrt; | ||
| 204 | 1445010 | aCache_[compIIdx][compJIdx] = | |
| 205 | 1445010 | sqrt(this->pureParams_[compIIdx].a() | |
| 206 | 1445010 | * this->pureParams_[compJIdx].a()) | |
| 207 | 1445010 | * (1 - Psi); | |
| 208 | } | ||
| 209 | } | ||
| 210 | 29490 | } | |
| 211 | |||
| 212 | Scalar aCache_[numComponents][numComponents]; | ||
| 213 | }; | ||
| 214 | |||
| 215 | } // end namespace | ||
| 216 | |||
| 217 | #endif | ||
| 218 |