| 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 IAPWS | ||
| 10 | * \brief Implements relations common for all regions of the IAPWS '97 | ||
| 11 | * formulation. | ||
| 12 | * See: | ||
| 13 | * | ||
| 14 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
| 15 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
| 16 | * http://www.iapws.org/relguide/IF97-Rev.pdf | ||
| 17 | */ | ||
| 18 | #ifndef DUMUX_IAPWS_COMMON_HH | ||
| 19 | #define DUMUX_IAPWS_COMMON_HH | ||
| 20 | |||
| 21 | #include <cmath> | ||
| 22 | #include <iostream> | ||
| 23 | |||
| 24 | #include <dune/common/math.hh> | ||
| 25 | |||
| 26 | #include <dumux/material/constants.hh> | ||
| 27 | |||
| 28 | namespace Dumux::IAPWS { | ||
| 29 | |||
| 30 | /*! | ||
| 31 | * \ingroup IAPWS | ||
| 32 | * \brief Implements relations which are common for all regions of the IAPWS '97 | ||
| 33 | * formulation. | ||
| 34 | * | ||
| 35 | * \tparam Scalar The type used for scalar values | ||
| 36 | * | ||
| 37 | * See: | ||
| 38 | * | ||
| 39 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
| 40 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
| 41 | * http://www.iapws.org/relguide/IF97-Rev.pdf | ||
| 42 | */ | ||
| 43 | template <class Scalar> | ||
| 44 | class Common | ||
| 45 | { | ||
| 46 | public: | ||
| 47 | //! The molar mass of water \f$\mathrm{[kg/mol]}\f$ | ||
| 48 | static constexpr Scalar molarMass = 18.01518e-3; | ||
| 49 | |||
| 50 | //! Specific gas constant of water \f$\mathrm{[J/(kg*K)]}\f$ | ||
| 51 | static constexpr Scalar Rs = Constants<Scalar>::R / molarMass; | ||
| 52 | |||
| 53 | //! Critical temperature of water \f$\mathrm{[K]}\f$ | ||
| 54 | static constexpr Scalar criticalTemperature = 647.096; | ||
| 55 | |||
| 56 | //! Critical pressure of water \f$\mathrm{[Pa]}\f$ | ||
| 57 | static constexpr Scalar criticalPressure = 22.064e6; | ||
| 58 | |||
| 59 | //! Critical molar volume of water \f$\mathrm{[m^3/mol]}\f$ | ||
| 60 | static constexpr Scalar criticalMolarVolume = molarMass / 322.0; | ||
| 61 | |||
| 62 | //! The acentric factor of water \f$\mathrm{[-]}\f$ | ||
| 63 | static constexpr Scalar acentricFactor = 0.344; | ||
| 64 | |||
| 65 | //! Density of water at the critical point \f$\mathrm{[kg/m^3]}\f$ | ||
| 66 | static constexpr Scalar criticalDensity = 322; | ||
| 67 | |||
| 68 | //! Triple temperature of water \f$\mathrm{[K]}\f$ | ||
| 69 | static constexpr Scalar tripleTemperature = 273.16; | ||
| 70 | |||
| 71 | //! Triple pressure of water \f$\mathrm{[Pa]}\f$ | ||
| 72 | static constexpr Scalar triplePressure = 611.657; | ||
| 73 | |||
| 74 | /*! | ||
| 75 | * \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of pure water. | ||
| 76 | * | ||
| 77 | * This relation is valid for all regions of the IAPWS '97 | ||
| 78 | * formulation. | ||
| 79 | * | ||
| 80 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
| 81 | * \param rho density of component in \f$\mathrm{[kg/m^3]}\f$ | ||
| 82 | * | ||
| 83 | * See: | ||
| 84 | * | ||
| 85 | * IAPWS: "Release on the IAPWS Formulation 2008 for the Viscosity | ||
| 86 | * of Ordinary Water Substance", http://www.iapws.org/relguide/visc.pdf | ||
| 87 | */ | ||
| 88 | 7859295 | static Scalar viscosity(Scalar temperature, Scalar rho) | |
| 89 | { | ||
| 90 | 7859295 | const Scalar rhoBar = rho/322.0; | |
| 91 | 7859295 | const Scalar TBar = temperature/criticalTemperature; | |
| 92 | |||
| 93 | // muBar = muBar_1 | ||
| 94 | 7859295 | constexpr Scalar Hij[6][7] = { | |
| 95 | { 5.20094e-1, 2.22531e-1,-2.81378e-1, 1.61913e-1,-3.25372e-2, 0, 0 }, | ||
| 96 | { 8.50895e-2, 9.99115e-1,-9.06851e-1, 2.57399e-1, 0, 0, 0 }, | ||
| 97 | {-1.08374 , 1.88797 ,-7.72479e-1, 0, 0, 0, 0 }, | ||
| 98 | {-2.89555e-1, 1.26613 ,-4.89837e-1, 0, 6.98452e-2, 0,-4.35673e-3 }, | ||
| 99 | { 0, 0,-2.57040e-1, 0, 0, 8.72102e-3, 0 }, | ||
| 100 | { 0, 1.20573e-1, 0, 0, 0, 0,-5.93264e-4 } | ||
| 101 | }; | ||
| 102 | |||
| 103 | 7859295 | Scalar tmp, tmp2, tmp3 = 1; | |
| 104 | 7859295 | Scalar muBar = 0; | |
| 105 |
2/2✓ Branch 0 taken 47155770 times.
✓ Branch 1 taken 7859295 times.
|
55015065 | for (int i = 0; i <= 5; ++i) { |
| 106 | tmp = 0; | ||
| 107 | tmp2 = 1; | ||
| 108 |
2/2✓ Branch 0 taken 330090390 times.
✓ Branch 1 taken 47155770 times.
|
377246160 | for (int j = 0; j <= 6; ++j) { |
| 109 | 330090390 | tmp += Hij[i][j]*tmp2; | |
| 110 | 330090390 | tmp2 *= (rhoBar - 1); | |
| 111 | } | ||
| 112 | 47155770 | muBar += tmp3 * tmp; | |
| 113 | 47155770 | tmp3 *= 1.0/TBar - 1; | |
| 114 | } | ||
| 115 | using std::exp; | ||
| 116 | 7859295 | muBar *= rhoBar; | |
| 117 | 7859295 | muBar = exp(muBar); | |
| 118 | |||
| 119 | // muBar *= muBar_0 | ||
| 120 | using std::sqrt; | ||
| 121 | 7859295 | muBar *= 100*sqrt(TBar); | |
| 122 | 7859295 | constexpr Scalar H[4] = { | |
| 123 | 1.67752, 2.20462, 0.6366564, -0.241605 | ||
| 124 | }; | ||
| 125 | |||
| 126 | 7859295 | tmp = 0, tmp2 = 1; | |
| 127 |
2/2✓ Branch 0 taken 31437180 times.
✓ Branch 1 taken 7859295 times.
|
39296475 | for (int i = 0; i < 4; ++i) { |
| 128 | 31437180 | tmp += H[i]/tmp2; | |
| 129 | 31437180 | tmp2 *= TBar; | |
| 130 | } | ||
| 131 | 7859295 | muBar /= tmp; | |
| 132 | |||
| 133 | 7859295 | return 1e-6*muBar; | |
| 134 | } | ||
| 135 | |||
| 136 | /*! | ||
| 137 | * \brief Thermal conductivity \f$\mathrm{[[W/(m*K)]}\f$ water (IAPWS) . | ||
| 138 | * | ||
| 139 | * Implementation taken from: | ||
| 140 | * freesteam - IAPWS-IF97 steam tables library | ||
| 141 | * copyright (C) 2004-2009 John Pye | ||
| 142 | * | ||
| 143 | * Appendix B: Recommended Interpolating equation for Industrial Use | ||
| 144 | * see http://www.iapws.org/relguide/thcond.pdf | ||
| 145 | * | ||
| 146 | * \param T absolute temperature in \f$\mathrm{[K]}\f$ | ||
| 147 | * \param rho density of water in \f$\mathrm{[kg/m^3]}\f$ | ||
| 148 | */ | ||
| 149 | 6604926 | static Scalar thermalConductivityIAPWS(const Scalar T, const Scalar rho) | |
| 150 | { | ||
| 151 | 6604926 | constexpr Scalar thcond_tstar = 647.26 ; | |
| 152 | 6604926 | constexpr Scalar thcond_rhostar = 317.7 ; | |
| 153 | /*static constexpr Scalar thcond_kstar = 1.0 ;*/ | ||
| 154 | |||
| 155 | 6604926 | constexpr Scalar thcond_b0 = -0.397070 ; | |
| 156 | 6604926 | constexpr Scalar thcond_b1 = 0.400302 ; | |
| 157 | 6604926 | constexpr Scalar thcond_b2 = 1.060000 ; | |
| 158 | 6604926 | constexpr Scalar thcond_B1 = -0.171587 ; | |
| 159 | 6604926 | constexpr Scalar thcond_B2 = 2.392190 ; | |
| 160 | |||
| 161 | 6604926 | constexpr Scalar thcond_c1 = 0.642857 ; | |
| 162 | 6604926 | constexpr Scalar thcond_c2 = -4.11717 ; | |
| 163 | 6604926 | constexpr Scalar thcond_c3 = -6.17937 ; | |
| 164 | 6604926 | constexpr Scalar thcond_c4 = 0.00308976 ; | |
| 165 | 6604926 | constexpr Scalar thcond_c5 = 0.0822994 ; | |
| 166 | 6604926 | constexpr Scalar thcond_c6 = 10.0932 ; | |
| 167 | |||
| 168 | 6604926 | constexpr Scalar thcond_d1 = 0.0701309 ; | |
| 169 | 6604926 | constexpr Scalar thcond_d2 = 0.0118520 ; | |
| 170 | 6604926 | constexpr Scalar thcond_d3 = 0.00169937 ; | |
| 171 | 6604926 | constexpr Scalar thcond_d4 = -1.0200 ; | |
| 172 | 6604926 | constexpr unsigned int thcond_a_count = 4; | |
| 173 | 6604926 | constexpr Scalar thcond_a[thcond_a_count] = { | |
| 174 | 0.0102811 | ||
| 175 | ,0.0299621 | ||
| 176 | ,0.0156146 | ||
| 177 | ,-0.00422464 | ||
| 178 | }; | ||
| 179 | |||
| 180 | 6604926 | const Scalar Tbar = T / thcond_tstar; | |
| 181 | 6604926 | const Scalar rhobar = rho / thcond_rhostar; | |
| 182 | |||
| 183 | /* fast implementation... minimised calls to 'pow' routine... */ | ||
| 184 | using std::sqrt; | ||
| 185 | 6604926 | const Scalar Troot = sqrt(Tbar); | |
| 186 | 6604926 | Scalar Tpow = Troot; | |
| 187 | 6604926 | Scalar lam = 0; | |
| 188 | |||
| 189 |
2/2✓ Branch 0 taken 26419704 times.
✓ Branch 1 taken 6604926 times.
|
33024630 | for(unsigned int k = 0; k < thcond_a_count; ++k) { |
| 190 | 26419704 | lam += thcond_a[k] * Tpow; | |
| 191 | 26419704 | Tpow *= Tbar; | |
| 192 | } | ||
| 193 | |||
| 194 | using std::exp; | ||
| 195 | 13209852 | lam += thcond_b0 + thcond_b1 | |
| 196 | 13209852 | * rhobar + thcond_b2 | |
| 197 | 6604926 | * exp(thcond_B1 * ((rhobar + thcond_B2)*(rhobar + thcond_B2))); | |
| 198 | |||
| 199 | using std::abs; | ||
| 200 | using std::pow; | ||
| 201 | using Dune::power; | ||
| 202 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 6604926 times.
|
6604926 | const Scalar DTbar = abs(Tbar - 1) + thcond_c4; |
| 203 | 6604926 | const Scalar DTbarpow = pow(DTbar, 3./5); | |
| 204 | 6604926 | const Scalar Q = 2. + thcond_c5 / DTbarpow; | |
| 205 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 6604926 times.
|
6604926 | const Scalar S = (Tbar >= 1) ? 1. / DTbar : thcond_c6 / DTbarpow; |
| 206 | |||
| 207 | 6604926 | const Scalar rhobar18 = pow(rhobar, 1.8); | |
| 208 | 6604926 | const Scalar rhobarQ = pow(rhobar, Q); | |
| 209 | |||
| 210 | 6604926 | lam += | |
| 211 | 6604926 | (thcond_d1 / power(Tbar,10) + thcond_d2) * rhobar18 * | |
| 212 | 6604926 | exp(thcond_c1 * (1 - rhobar * rhobar18)) | |
| 213 | 13209852 | + thcond_d3 * S * rhobarQ * | |
| 214 | 6604926 | exp((Q/(1+Q))*(1 - rhobar*rhobarQ)) | |
| 215 | 6604926 | + thcond_d4 * | |
| 216 | 13209852 | exp(thcond_c2 * power(Troot,3) + thcond_c3 / power(rhobar,5)); | |
| 217 | 6604926 | return /*thcond_kstar * */ lam; | |
| 218 | } | ||
| 219 | }; | ||
| 220 | |||
| 221 | } // end namespace Dumux::IAPWS | ||
| 222 | |||
| 223 | #endif | ||
| 224 |