| 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 FluidSystems | ||
| 10 | * \brief @copybrief Dumux::FluidSystems::H2OHeavyOil | ||
| 11 | */ | ||
| 12 | #ifndef DUMUX_H2O_HEAVYOIL_FLUID_SYSTEM_HH | ||
| 13 | #define DUMUX_H2O_HEAVYOIL_FLUID_SYSTEM_HH | ||
| 14 | |||
| 15 | #include <dumux/material/idealgas.hh> | ||
| 16 | #include <dumux/material/components/h2o.hh> | ||
| 17 | #include <dumux/material/components/tabulatedcomponent.hh> | ||
| 18 | #include <dumux/material/components/heavyoil.hh> | ||
| 19 | |||
| 20 | #include <dumux/material/binarycoefficients/h2o_heavyoil.hh> | ||
| 21 | |||
| 22 | #include <dumux/material/fluidsystems/base.hh> | ||
| 23 | |||
| 24 | #include <dumux/io/name.hh> | ||
| 25 | |||
| 26 | namespace Dumux { | ||
| 27 | namespace FluidSystems { | ||
| 28 | |||
| 29 | /*! | ||
| 30 | * \ingroup FluidSystems | ||
| 31 | * \brief A compositional fluid system with water and heavy oil | ||
| 32 | * components in both the liquid and the gas phase. | ||
| 33 | */ | ||
| 34 | template <class Scalar, | ||
| 35 | class H2OType = Dumux::Components::TabulatedComponent<Dumux::Components::H2O<Scalar> > > | ||
| 36 | class H2OHeavyOil | ||
| 37 | : public Base<Scalar, H2OHeavyOil<Scalar, H2OType> > | ||
| 38 | { | ||
| 39 | using ThisType = H2OHeavyOil<Scalar, H2OType>; | ||
| 40 | |||
| 41 | public: | ||
| 42 | using HeavyOil = Dumux::Components::HeavyOil<Scalar>; | ||
| 43 | using H2O = H2OType; | ||
| 44 | |||
| 45 | |||
| 46 | static const int numPhases = 3; | ||
| 47 | static const int numComponents = 2; | ||
| 48 | |||
| 49 | static const int wPhaseIdx = 0; // index of the water phase | ||
| 50 | static const int nPhaseIdx = 1; // index of the NAPL phase | ||
| 51 | static const int gPhaseIdx = 2; // index of the gas phase | ||
| 52 | |||
| 53 | static const int H2OIdx = 0; | ||
| 54 | static const int NAPLIdx = 1; | ||
| 55 | |||
| 56 | // export component indices to indicate the main component | ||
| 57 | // of the corresponding phase at atmospheric pressure 1 bar | ||
| 58 | // and room temperature 20°C: | ||
| 59 | static const int wCompIdx = H2OIdx; | ||
| 60 | static const int nCompIdx = NAPLIdx; | ||
| 61 | |||
| 62 | /*! | ||
| 63 | * \brief Initialize the fluid system's static parameters generically | ||
| 64 | * | ||
| 65 | * If a tabulated H2O component is used, we do our best to create | ||
| 66 | * tables that always work. | ||
| 67 | */ | ||
| 68 | static void init() | ||
| 69 | { | ||
| 70 |
1/3✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
|
2 | init(/*tempMin=*/273.15, |
| 71 | /*tempMax=*/623.15, | ||
| 72 | /*numTemp=*/100, | ||
| 73 | /*pMin=*/0.0, | ||
| 74 | /*pMax=*/20e6, | ||
| 75 | /*numP=*/200); | ||
| 76 | } | ||
| 77 | |||
| 78 | /*! | ||
| 79 | * \brief Initialize the fluid system's static parameters using | ||
| 80 | * problem specific temperature and pressure ranges | ||
| 81 | * | ||
| 82 | * \param tempMin The minimum temperature used for tabulation of water [K] | ||
| 83 | * \param tempMax The maximum temperature used for tabulation of water [K] | ||
| 84 | * \param nTemp The number of ticks on the temperature axis of the table of water | ||
| 85 | * \param pressMin The minimum pressure used for tabulation of water [Pa] | ||
| 86 | * \param pressMax The maximum pressure used for tabulation of water [Pa] | ||
| 87 | * \param nPress The number of ticks on the pressure axis of the table of water | ||
| 88 | */ | ||
| 89 | static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, | ||
| 90 | Scalar pressMin, Scalar pressMax, unsigned nPress) | ||
| 91 | { | ||
| 92 | if (H2O::isTabulated) | ||
| 93 | { | ||
| 94 |
1/2✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
|
2 | H2O::init(tempMin, tempMax, nTemp, |
| 95 | pressMin, pressMax, nPress); | ||
| 96 | } | ||
| 97 | } | ||
| 98 | |||
| 99 | /*! | ||
| 100 | * \brief Get the main component of a given phase | ||
| 101 | * \param phaseIdx The index of the fluid phase to consider | ||
| 102 | */ | ||
| 103 | static constexpr int getMainComponent(int phaseIdx) | ||
| 104 | { | ||
| 105 | // For the gas phase, choosing a main component appears to be | ||
| 106 | // rather arbitrary. Motivated by the fact that the thermal conductivity | ||
| 107 | // of the gas phase is set to the thermal conductivity of pure water, | ||
| 108 | // water is chosen for now. | ||
| 109 |
6/6✓ Branch 0 taken 15974400 times.
✓ Branch 1 taken 7987200 times.
✓ Branch 2 taken 31948800 times.
✓ Branch 3 taken 15974400 times.
✓ Branch 4 taken 15974400 times.
✓ Branch 5 taken 7987200 times.
|
95846400 | if (phaseIdx == nPhaseIdx) |
| 110 | return nCompIdx; | ||
| 111 | else | ||
| 112 | return wCompIdx; | ||
| 113 | } | ||
| 114 | |||
| 115 | /*! | ||
| 116 | * \brief Returns whether the fluids are miscible | ||
| 117 | */ | ||
| 118 | static constexpr bool isMiscible() | ||
| 119 | { return true; } | ||
| 120 | |||
| 121 | /*! | ||
| 122 | * \brief Return whether a phase is gaseous | ||
| 123 | * | ||
| 124 | * \param phaseIdx The index of the fluid phase to consider | ||
| 125 | */ | ||
| 126 | static constexpr bool isGas(int phaseIdx) | ||
| 127 | { | ||
| 128 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
|
2 | assert(0 <= phaseIdx && phaseIdx < numPhases); |
| 129 | return phaseIdx == gPhaseIdx; | ||
| 130 | } | ||
| 131 | |||
| 132 | /*! | ||
| 133 | * \brief Returns true if and only if a fluid phase is assumed to | ||
| 134 | * be an ideal gas. | ||
| 135 | * | ||
| 136 | * \param phaseIdx The index of the fluid phase to consider | ||
| 137 | */ | ||
| 138 | ✗ | static bool isIdealGas(int phaseIdx) | |
| 139 | ✗ | { return phaseIdx == gPhaseIdx && H2O::gasIsIdeal() && HeavyOil::gasIsIdeal(); } | |
| 140 | |||
| 141 | /*! | ||
| 142 | * \brief Returns true if and only if a fluid phase is assumed to | ||
| 143 | * be an ideal mixture. | ||
| 144 | * | ||
| 145 | * We define an ideal mixture as a fluid phase where the fugacity | ||
| 146 | * coefficients of all components times the pressure of the phase | ||
| 147 | * are independent on the fluid composition. This assumption is true | ||
| 148 | * if Henry's law and Raoult's law apply. If you are unsure what | ||
| 149 | * this function should return, it is safe to return false. The | ||
| 150 | * only damage done will be (slightly) increased computation times | ||
| 151 | * in some cases. | ||
| 152 | * | ||
| 153 | * \param phaseIdx The index of the fluid phase to consider | ||
| 154 | */ | ||
| 155 | static bool isIdealMixture(int phaseIdx) | ||
| 156 | { | ||
| 157 | assert(0 <= phaseIdx && phaseIdx < numPhases); | ||
| 158 | // we assume Henry's and Raoult's laws for the water phase and | ||
| 159 | // and no interaction between gas molecules of different | ||
| 160 | // components, so all phases are ideal mixtures! | ||
| 161 | return true; | ||
| 162 | } | ||
| 163 | |||
| 164 | /*! | ||
| 165 | * \brief Returns true if and only if a fluid phase is assumed to | ||
| 166 | * be compressible. | ||
| 167 | * | ||
| 168 | * Compressible means that the partial derivative of the density | ||
| 169 | * to the fluid pressure is always larger than zero. | ||
| 170 | * | ||
| 171 | * \param phaseIdx The index of the fluid phase to consider | ||
| 172 | */ | ||
| 173 | static constexpr bool isCompressible(int phaseIdx) | ||
| 174 | { | ||
| 175 | assert(0 <= phaseIdx && phaseIdx < numPhases); | ||
| 176 | // gases are always compressible | ||
| 177 | if (phaseIdx == gPhaseIdx) | ||
| 178 | return true; | ||
| 179 | else if (phaseIdx == wPhaseIdx) | ||
| 180 | // the water component decides for the water phase... | ||
| 181 | return H2O::liquidIsCompressible(); | ||
| 182 | |||
| 183 | // the NAPL component decides for the napl phase... | ||
| 184 | return HeavyOil::liquidIsCompressible(); | ||
| 185 | } | ||
| 186 | |||
| 187 | /*! | ||
| 188 | * \brief Return the human readable name of a phase (used in indices) | ||
| 189 | */ | ||
| 190 | 24 | static std::string phaseName(int phaseIdx) | |
| 191 | { | ||
| 192 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 24 times.
|
24 | assert(0 <= phaseIdx && phaseIdx < numPhases); |
| 193 |
3/4✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 8 times.
✗ Branch 3 not taken.
|
24 | switch (phaseIdx) |
| 194 | { | ||
| 195 | 8 | case wPhaseIdx: return IOName::aqueousPhase(); | |
| 196 | 8 | case nPhaseIdx: return IOName::naplPhase(); | |
| 197 | 8 | case gPhaseIdx: return IOName::gaseousPhase(); | |
| 198 | } | ||
| 199 | ✗ | DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx); | |
| 200 | } | ||
| 201 | |||
| 202 | /*! | ||
| 203 | * \brief Return the human readable name of a component (used in indices) | ||
| 204 | */ | ||
| 205 | 8 | static std::string componentName(int compIdx) | |
| 206 | { | ||
| 207 |
2/3✓ Branch 0 taken 4 times.
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
|
8 | switch (compIdx) { |
| 208 | 4 | case H2OIdx: return H2O::name(); | |
| 209 | 4 | case NAPLIdx: return HeavyOil::name(); | |
| 210 | }; | ||
| 211 | ✗ | DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx); | |
| 212 | } | ||
| 213 | |||
| 214 | /*! | ||
| 215 | * \brief Return the molar mass of a component in [kg/mol]. | ||
| 216 | */ | ||
| 217 | 197738884 | static Scalar molarMass(int compIdx) | |
| 218 | { | ||
| 219 |
2/3✓ Branch 0 taken 114843842 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 82895042 times.
|
197738884 | switch (compIdx) { |
| 220 | case H2OIdx: return H2O::molarMass(); | ||
| 221 | 114843842 | case NAPLIdx: return HeavyOil::molarMass(); | |
| 222 | }; | ||
| 223 | ✗ | DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx); | |
| 224 | } | ||
| 225 | |||
| 226 | using Base<Scalar, ThisType>::density; | ||
| 227 | //! \copydoc Base<Scalar,ThisType>::density(const FluidState&,int) | ||
| 228 | template <class FluidState> | ||
| 229 | 11564843 | static Scalar density(const FluidState &fluidState, int phaseIdx) | |
| 230 | { | ||
| 231 |
2/2✓ Branch 0 taken 3854949 times.
✓ Branch 1 taken 7709894 times.
|
11564843 | if (phaseIdx == wPhaseIdx) { |
| 232 | // See: doctoral thesis of Steffen Ochs 2007 | ||
| 233 | // Steam injection into saturated porous media : process analysis including experimental and numerical investigations | ||
| 234 | // http://elib.uni-stuttgart.de/bitstream/11682/271/1/Diss_Ochs_OPUS.pdf | ||
| 235 | |||
| 236 | // This assumes each gas molecule displaces exactly one | ||
| 237 | // molecule in the liquid. | ||
| 238 | |||
| 239 | 11564845 | return H2O::liquidMolarDensity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)) | |
| 240 | 3854947 | * (H2O::molarMass()*fluidState.moleFraction(wPhaseIdx, H2OIdx) | |
| 241 | 7709893 | + HeavyOil::molarMass()*fluidState.moleFraction(wPhaseIdx, NAPLIdx)); | |
| 242 | } | ||
| 243 |
2/2✓ Branch 0 taken 3854947 times.
✓ Branch 1 taken 3854947 times.
|
7709894 | else if (phaseIdx == nPhaseIdx) { |
| 244 | // assume pure NAPL for the NAPL phase | ||
| 245 | 3854947 | Scalar pressure = HeavyOil::liquidIsCompressible()?fluidState.pressure(phaseIdx):1e100; | |
| 246 | 11564839 | return HeavyOil::liquidDensity(fluidState.temperature(phaseIdx), pressure); | |
| 247 | } | ||
| 248 | |||
| 249 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 3854947 times.
|
3854947 | assert (phaseIdx == gPhaseIdx); |
| 250 | 3854947 | Scalar pH2O = | |
| 251 | 7709893 | fluidState.moleFraction(gPhaseIdx, H2OIdx) * | |
| 252 | fluidState.pressure(gPhaseIdx); | ||
| 253 | 3854947 | Scalar pNAPL = | |
| 254 | 7709893 | fluidState.moleFraction(gPhaseIdx, NAPLIdx) * | |
| 255 | fluidState.pressure(gPhaseIdx); | ||
| 256 | 7709893 | return H2O::gasDensity(fluidState.temperature(phaseIdx), pH2O) | |
| 257 | 11564839 | + HeavyOil::gasDensity(fluidState.temperature(phaseIdx), pNAPL); | |
| 258 | } | ||
| 259 | |||
| 260 | using Base<Scalar, ThisType>::molarDensity; | ||
| 261 | //! \copydoc Base<Scalar,ThisType>::molarDensity(const FluidState&,int) | ||
| 262 | template <class FluidState> | ||
| 263 | 7709895 | static Scalar molarDensity(const FluidState &fluidState, int phaseIdx) | |
| 264 | { | ||
| 265 |
2/2✓ Branch 0 taken 3854946 times.
✓ Branch 1 taken 3854946 times.
|
11564841 | Scalar temperature = fluidState.temperature(phaseIdx); |
| 266 |
2/2✓ Branch 0 taken 3854946 times.
✓ Branch 1 taken 3854946 times.
|
11564841 | Scalar pressure = fluidState.pressure(phaseIdx); |
| 267 |
2/2✓ Branch 0 taken 3854947 times.
✓ Branch 1 taken 3854948 times.
|
7709895 | if (phaseIdx == nPhaseIdx) |
| 268 | { | ||
| 269 | 7709894 | return HeavyOil::liquidMolarDensity(temperature, pressure); | |
| 270 | } | ||
| 271 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 3854947 times.
|
3854948 | else if (phaseIdx == wPhaseIdx) |
| 272 | { // This assumes each gas molecule displaces exactly one | ||
| 273 | // molecule in the liquid. | ||
| 274 | 3854947 | return H2O::liquidMolarDensity(temperature, pressure); | |
| 275 | } | ||
| 276 | else | ||
| 277 | { | ||
| 278 | 7709893 | return H2O::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, H2OIdx)) | |
| 279 | 11564839 | + HeavyOil::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, NAPLIdx)); | |
| 280 | } | ||
| 281 | } | ||
| 282 | |||
| 283 | using Base<Scalar, ThisType>::viscosity; | ||
| 284 | //! \copydoc Base<Scalar,ThisType>::viscosity(const FluidState&,int) | ||
| 285 | template <class FluidState> | ||
| 286 | 11564841 | static Scalar viscosity(const FluidState &fluidState, | |
| 287 | int phaseIdx) | ||
| 288 | { | ||
| 289 |
2/2✓ Branch 0 taken 3854947 times.
✓ Branch 1 taken 7709894 times.
|
11564841 | if (phaseIdx == wPhaseIdx) { |
| 290 | // assume pure water viscosity | ||
| 291 | 11564839 | return H2O::liquidViscosity(fluidState.temperature(phaseIdx), | |
| 292 | fluidState.pressure(phaseIdx)); | ||
| 293 | } | ||
| 294 |
2/2✓ Branch 0 taken 3854947 times.
✓ Branch 1 taken 3854947 times.
|
7709894 | else if (phaseIdx == nPhaseIdx) { |
| 295 | // assume pure NAPL viscosity | ||
| 296 | 11564839 | return HeavyOil::liquidViscosity(fluidState.temperature(phaseIdx), | |
| 297 | 3854947 | fluidState.pressure(phaseIdx)); | |
| 298 | } | ||
| 299 | |||
| 300 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 3854947 times.
|
3854947 | assert (phaseIdx == gPhaseIdx); |
| 301 | |||
| 302 | /* Wilke method. See: | ||
| 303 | * | ||
| 304 | * See: R. Reid, et al.: The Properties of Gases and Liquids, | ||
| 305 | * 4th edition, McGraw-Hill, 1987, 407-410 | ||
| 306 | * 5th edition, McGraw-Hill, 20001, p. 9.21/22 | ||
| 307 | * | ||
| 308 | * in this case, we use a simplified version in order to avoid | ||
| 309 | * computationally costly evaluation of sqrt and pow functions and | ||
| 310 | * divisions | ||
| 311 | * -- compare e.g. with Promo Class p. 32/33 | ||
| 312 | */ | ||
| 313 | 7709894 | const Scalar mu[numComponents] = { | |
| 314 | 11564839 | h2oGasViscosityInMixture(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)), | |
| 315 | 15419785 | HeavyOil::gasViscosity(fluidState.temperature(phaseIdx), HeavyOil::vaporPressure(fluidState.temperature(phaseIdx))) | |
| 316 | }; | ||
| 317 | |||
| 318 | 3854947 | return mu[H2OIdx]*fluidState.moleFraction(gPhaseIdx, H2OIdx) | |
| 319 | 7709893 | + mu[NAPLIdx]*fluidState.moleFraction(gPhaseIdx, NAPLIdx); | |
| 320 | } | ||
| 321 | |||
| 322 | using Base<Scalar, ThisType>::binaryDiffusionCoefficient; | ||
| 323 | /*! | ||
| 324 | * \brief Given a phase's composition, temperature and pressure, | ||
| 325 | * return the binary diffusion coefficient \f$\mathrm{[m^2/s]}\f$ for components | ||
| 326 | * \f$\mathrm{i}\f$ and \f$\mathrm{j}\f$ in this phase. | ||
| 327 | * | ||
| 328 | * Be aware that in this case there are only two components in three phases. | ||
| 329 | * Therefore, we assume the diffusion to simply be the binary diffusion coefficients. | ||
| 330 | * This was previously implemented in diffusionCoefficient(), but is now moved. | ||
| 331 | * | ||
| 332 | * \param fluidState The fluid state | ||
| 333 | * \param phaseIdx Index of the fluid phase | ||
| 334 | * \param compIIdx Index of the component i | ||
| 335 | * \param compJIdx Index of the component j | ||
| 336 | */ | ||
| 337 | template <class FluidState> | ||
| 338 | 7709904 | static Scalar binaryDiffusionCoefficient(const FluidState &fluidState, | |
| 339 | int phaseIdx, | ||
| 340 | int compIIdx, | ||
| 341 | int compJIdx) | ||
| 342 | |||
| 343 | { | ||
| 344 |
2/2✓ Branch 0 taken 3854946 times.
✓ Branch 1 taken 3854946 times.
|
7709904 | const Scalar T = fluidState.temperature(phaseIdx); |
| 345 |
2/2✓ Branch 0 taken 3854946 times.
✓ Branch 1 taken 3854946 times.
|
7709904 | const Scalar p = fluidState.pressure(phaseIdx); |
| 346 | |||
| 347 | // liquid phase | ||
| 348 |
2/2✓ Branch 0 taken 3854954 times.
✓ Branch 1 taken 3854950 times.
|
7709904 | if (phaseIdx == wPhaseIdx) |
| 349 | return BinaryCoeff::H2O_HeavyOil::liquidDiffCoeff(T, p); | ||
| 350 | |||
| 351 | // gas phase | ||
| 352 |
2/2✓ Branch 0 taken 4 times.
✓ Branch 1 taken 3854950 times.
|
3854954 | else if (phaseIdx == gPhaseIdx) |
| 353 | return BinaryCoeff::H2O_HeavyOil::gasDiffCoeff(T, p); | ||
| 354 | else | ||
| 355 |
7/16✓ Branch 2 taken 4 times.
✗ Branch 3 not taken.
✓ Branch 11 taken 4 times.
✗ Branch 12 not taken.
✓ Branch 16 taken 4 times.
✗ Branch 17 not taken.
✓ Branch 19 taken 4 times.
✗ Branch 20 not taken.
✓ Branch 22 taken 4 times.
✗ Branch 23 not taken.
✓ Branch 24 taken 4 times.
✗ Branch 25 not taken.
✗ Branch 27 not taken.
✓ Branch 28 taken 4 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
|
48 | DUNE_THROW(Dune::InvalidStateException, |
| 356 | "Non-existent binary diffusion coefficient for phase index " | ||
| 357 | << phaseIdx); | ||
| 358 | } | ||
| 359 | |||
| 360 | using Base<Scalar, ThisType>::diffusionCoefficient; | ||
| 361 | /*! | ||
| 362 | * \brief Calculate the molecular diffusion coefficient for | ||
| 363 | * a component in a fluid phase \f$\mathrm{[mol^2 * s / (kg*m^3)]}\f$ | ||
| 364 | * \param fluidState The fluid state | ||
| 365 | * \param phaseIdx Index of the fluid phase | ||
| 366 | * \todo This signature does not match the base fluid system's signature, | ||
| 367 | * which takes an additional component index. | ||
| 368 | */ | ||
| 369 | template <class FluidState> | ||
| 370 | 7709892 | static Scalar diffusionCoefficient(const FluidState &fluidState, int phaseIdx) | |
| 371 | { | ||
| 372 | // liquid phase | ||
| 373 |
2/2✓ Branch 0 taken 3854946 times.
✓ Branch 1 taken 3854946 times.
|
7709892 | if (phaseIdx == wPhaseIdx) |
| 374 | 3854946 | return binaryDiffusionCoefficient(fluidState, phaseIdx, H2OIdx, NAPLIdx); | |
| 375 | // gas phase | ||
| 376 |
1/2✓ Branch 0 taken 3854946 times.
✗ Branch 1 not taken.
|
3854946 | else if (phaseIdx == gPhaseIdx) |
| 377 | 3854946 | return binaryDiffusionCoefficient(fluidState, phaseIdx, NAPLIdx, H2OIdx); | |
| 378 | else | ||
| 379 | ✗ | DUNE_THROW(Dune::InvalidStateException, | |
| 380 | "Non-existent diffusion coefficient for phase index "<< phaseIdx); | ||
| 381 | } | ||
| 382 | |||
| 383 | /*! | ||
| 384 | * \brief Henry coefficients \f$[N/m^2]\f$ of a component in a phase. | ||
| 385 | * \param fluidState An arbitrary fluid state | ||
| 386 | * \param phaseIdx The index of the phase to consider | ||
| 387 | * \param compIdx The index of the component to consider | ||
| 388 | */ | ||
| 389 | template <class FluidState> | ||
| 390 | static Scalar henryCoefficient(const FluidState &fluidState, | ||
| 391 | int phaseIdx, | ||
| 392 | int compIdx) | ||
| 393 | { | ||
| 394 | assert(0 <= phaseIdx && phaseIdx < numPhases); | ||
| 395 | assert(0 <= compIdx && compIdx < numComponents); | ||
| 396 | |||
| 397 | const Scalar T = fluidState.temperature(phaseIdx); | ||
| 398 | |||
| 399 | if (compIdx == NAPLIdx && phaseIdx == wPhaseIdx) | ||
| 400 | return Dumux::BinaryCoeff::H2O_HeavyOil::henryOilInWater(T); | ||
| 401 | |||
| 402 | else if (phaseIdx == nPhaseIdx && compIdx == H2OIdx) | ||
| 403 | return Dumux::BinaryCoeff::H2O_HeavyOil::henryWaterInOil(T); | ||
| 404 | |||
| 405 | else | ||
| 406 | DUNE_THROW(Dune::InvalidStateException, "non-existent henry coefficient for phase index " << phaseIdx | ||
| 407 | << " and component index " << compIdx); | ||
| 408 | } | ||
| 409 | |||
| 410 | /*! | ||
| 411 | * \brief Partial pressures in the gas phase, taken from saturation vapor pressures. | ||
| 412 | * \param fluidState An arbitrary fluid state | ||
| 413 | * \param phaseIdx The index of the phase to consider | ||
| 414 | * \param compIdx The index of the component to consider | ||
| 415 | */ | ||
| 416 | template <class FluidState> | ||
| 417 | 7712470 | static Scalar partialPressureGas(const FluidState &fluidState, int phaseIdx, | |
| 418 | int compIdx) | ||
| 419 | { | ||
| 420 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 7712470 times.
|
7712470 | assert(0 <= compIdx && compIdx < numComponents); |
| 421 | |||
| 422 |
2/2✓ Branch 0 taken 3856050 times.
✓ Branch 1 taken 3856420 times.
|
7712470 | const Scalar T = fluidState.temperature(phaseIdx); |
| 423 |
2/2✓ Branch 0 taken 3856050 times.
✓ Branch 1 taken 3856420 times.
|
7712470 | if (compIdx == NAPLIdx) |
| 424 | 7712100 | return HeavyOil::vaporPressure(T); | |
| 425 | else if (compIdx == H2OIdx) | ||
| 426 | 3856420 | return H2O::vaporPressure(T); | |
| 427 | else | ||
| 428 | DUNE_THROW(Dune::InvalidStateException, "non-existent component index " << compIdx); | ||
| 429 | } | ||
| 430 | |||
| 431 | /*! | ||
| 432 | * \brief Inverse vapor pressures, taken from inverse saturation vapor pressures | ||
| 433 | * \param fluidState An arbitrary fluid state | ||
| 434 | * \param phaseIdx The index of the phase to consider | ||
| 435 | * \param compIdx The index of the component to consider | ||
| 436 | */ | ||
| 437 | template <class FluidState> | ||
| 438 | 750 | static Scalar inverseVaporPressureCurve(const FluidState &fluidState, | |
| 439 | int phaseIdx, | ||
| 440 | int compIdx) | ||
| 441 | { | ||
| 442 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 750 times.
|
750 | assert(0 <= compIdx && compIdx < numComponents); |
| 443 | |||
| 444 |
2/2✓ Branch 0 taken 370 times.
✓ Branch 1 taken 380 times.
|
750 | const Scalar pressure = fluidState.pressure(phaseIdx); |
| 445 |
2/2✓ Branch 0 taken 370 times.
✓ Branch 1 taken 380 times.
|
750 | if (compIdx == NAPLIdx) |
| 446 | 740 | return HeavyOil::vaporTemperature(pressure); | |
| 447 | else if (compIdx == H2OIdx) | ||
| 448 | 760 | return H2O::vaporTemperature(pressure); | |
| 449 | else | ||
| 450 | DUNE_THROW(Dune::InvalidStateException, "non-existent component index " << compIdx); | ||
| 451 | } | ||
| 452 | |||
| 453 | using Base<Scalar, ThisType>::enthalpy; | ||
| 454 | /*! | ||
| 455 | * \brief Given all mole fractions in a phase, return the specific | ||
| 456 | * phase enthalpy\f$\mathrm{[J/kg]}\f$. | ||
| 457 | * \param fluidState An arbitrary fluid state | ||
| 458 | * \param phaseIdx The index of the phase to consider | ||
| 459 | * \todo This system neglects the contribution of gas-molecules in the liquid phase. | ||
| 460 | * This contribution is probably not big. | ||
| 461 | * Somebody would have to find out the enthalpy of solution for this system. ... | ||
| 462 | */ | ||
| 463 | template <class FluidState> | ||
| 464 | 11564841 | static Scalar enthalpy(const FluidState &fluidState, | |
| 465 | int phaseIdx) | ||
| 466 | { | ||
| 467 |
2/2✓ Branch 0 taken 3854947 times.
✓ Branch 1 taken 7709894 times.
|
11564841 | if (phaseIdx == wPhaseIdx) { |
| 468 | 11564839 | return H2O::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)); | |
| 469 | } | ||
| 470 |
2/2✓ Branch 0 taken 3854947 times.
✓ Branch 1 taken 3854947 times.
|
7709894 | else if (phaseIdx == nPhaseIdx) { |
| 471 | 15419785 | return HeavyOil::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)); | |
| 472 | } | ||
| 473 |
1/2✓ Branch 0 taken 3854947 times.
✗ Branch 1 not taken.
|
3854947 | else if (phaseIdx == gPhaseIdx) { // gas phase enthalpy depends strongly on composition |
| 474 | 11564839 | Scalar hgc = HeavyOil::gasEnthalpy(fluidState.temperature(phaseIdx), | |
| 475 | fluidState.pressure(phaseIdx)); | ||
| 476 | 11564839 | Scalar hgw = H2O::gasEnthalpy(fluidState.temperature(phaseIdx), | |
| 477 | fluidState.pressure(phaseIdx)); | ||
| 478 | |||
| 479 | 3854947 | Scalar result = 0; | |
| 480 | 3854947 | result += hgw * fluidState.massFraction(gPhaseIdx, H2OIdx); | |
| 481 | 3854947 | result += hgc * fluidState.massFraction(gPhaseIdx, NAPLIdx); | |
| 482 | |||
| 483 | 3854947 | return result; | |
| 484 | } | ||
| 485 | ✗ | DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx); | |
| 486 | } | ||
| 487 | |||
| 488 | /*! | ||
| 489 | * \brief Returns the specific enthalpy \f$\mathrm{[J/kg]}\f$ of a component in a specific phase | ||
| 490 | * \param fluidState The fluid state | ||
| 491 | * \param phaseIdx The index of the phase | ||
| 492 | * \param componentIdx The index of the component | ||
| 493 | */ | ||
| 494 | template <class FluidState> | ||
| 495 | static Scalar componentEnthalpy(const FluidState& fluidState, int phaseIdx, int componentIdx) | ||
| 496 | { | ||
| 497 | const Scalar T = fluidState.temperature(phaseIdx); | ||
| 498 | const Scalar p = fluidState.pressure(phaseIdx); | ||
| 499 | |||
| 500 | if (phaseIdx == wPhaseIdx) | ||
| 501 | { | ||
| 502 | if (componentIdx == H2OIdx) | ||
| 503 | return H2O::liquidEnthalpy(T, p); | ||
| 504 | else if (componentIdx == NAPLIdx) | ||
| 505 | DUNE_THROW(Dune::NotImplemented, "The component enthalpy for NAPL in water is not implemented."); | ||
| 506 | DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx); | ||
| 507 | } | ||
| 508 | else if (phaseIdx == nPhaseIdx) | ||
| 509 | { | ||
| 510 | if (componentIdx == H2OIdx) | ||
| 511 | DUNE_THROW(Dune::NotImplemented, "The component enthalpy for water in NAPL is not implemented."); | ||
| 512 | else if (componentIdx == NAPLIdx) | ||
| 513 | return HeavyOil::liquidEnthalpy(T, p); | ||
| 514 | DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx); | ||
| 515 | } | ||
| 516 | else if (phaseIdx == gPhaseIdx) | ||
| 517 | { | ||
| 518 | if (componentIdx == H2OIdx) | ||
| 519 | return H2O::gasEnthalpy(T, p); | ||
| 520 | else if (componentIdx == NAPLIdx) | ||
| 521 | return HeavyOil::gasEnthalpy(T, p); | ||
| 522 | DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx); | ||
| 523 | } | ||
| 524 | DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx); | ||
| 525 | } | ||
| 526 | |||
| 527 | using Base<Scalar, ThisType>::heatCapacity; | ||
| 528 | //! \copydoc Base<Scalar,ThisType>::heatCapacity(const FluidState&,int) | ||
| 529 | template <class FluidState> | ||
| 530 | 3 | static Scalar heatCapacity(const FluidState &fluidState, | |
| 531 | int phaseIdx) | ||
| 532 | { | ||
| 533 |
7/16✓ Branch 2 taken 3 times.
✗ Branch 3 not taken.
✓ Branch 11 taken 3 times.
✗ Branch 12 not taken.
✓ Branch 15 taken 3 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 3 times.
✗ Branch 19 not taken.
✓ Branch 21 taken 3 times.
✗ Branch 22 not taken.
✓ Branch 23 taken 3 times.
✗ Branch 24 not taken.
✗ Branch 26 not taken.
✓ Branch 27 taken 3 times.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
|
33 | DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2ONAPL::heatCapacity()"); |
| 534 | } | ||
| 535 | |||
| 536 | using Base<Scalar, ThisType>::thermalConductivity; | ||
| 537 | /*! | ||
| 538 | * \brief Thermal conductivity of a fluid phase \f$\mathrm{[W/(m K)]}\f$. | ||
| 539 | * | ||
| 540 | * Use the conductivity of water (wPhase and gPhase) and oil (nPhase) as a first approximation. | ||
| 541 | * | ||
| 542 | * \param fluidState An arbitrary fluid state | ||
| 543 | * \param phaseIdx The index of the fluid phase to consider | ||
| 544 | */ | ||
| 545 | template <class FluidState> | ||
| 546 | 11564841 | static Scalar thermalConductivity(const FluidState &fluidState, | |
| 547 | int phaseIdx) | ||
| 548 | { | ||
| 549 |
2/2✓ Branch 0 taken 3854946 times.
✓ Branch 1 taken 7709892 times.
|
11564841 | const Scalar temperature = fluidState.temperature(phaseIdx) ; |
| 550 |
2/2✓ Branch 0 taken 3854946 times.
✓ Branch 1 taken 7709892 times.
|
11564841 | const Scalar pressure = fluidState.pressure(phaseIdx); |
| 551 |
2/2✓ Branch 0 taken 3854947 times.
✓ Branch 1 taken 7709894 times.
|
11564841 | if (phaseIdx == wPhaseIdx) |
| 552 | { | ||
| 553 | 3854947 | return H2O::liquidThermalConductivity(temperature, pressure); | |
| 554 | } | ||
| 555 |
2/2✓ Branch 0 taken 3854947 times.
✓ Branch 1 taken 3854947 times.
|
7709894 | else if (phaseIdx == nPhaseIdx) |
| 556 | { | ||
| 557 | return HeavyOil::liquidThermalConductivity(temperature, pressure); | ||
| 558 | } | ||
| 559 |
1/2✓ Branch 0 taken 3854947 times.
✗ Branch 1 not taken.
|
3854947 | else if (phaseIdx == gPhaseIdx) |
| 560 | { | ||
| 561 | 3854947 | return H2O::gasThermalConductivity(temperature, pressure); | |
| 562 | } | ||
| 563 | ✗ | DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx); | |
| 564 | } | ||
| 565 | |||
| 566 | }; | ||
| 567 | } // end namespace FluidSystems | ||
| 568 | } // end namespace Dumux | ||
| 569 | |||
| 570 | #endif | ||
| 571 |