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 Components | ||
10 | * \brief Material properties of pure water \f$H_2O\f$. | ||
11 | */ | ||
12 | #ifndef DUMUX_H2O_HH | ||
13 | #define DUMUX_H2O_HH | ||
14 | |||
15 | #include <cmath> | ||
16 | #include <cassert> | ||
17 | |||
18 | #include <dumux/nonlinear/findscalarroot.hh> | ||
19 | #include <dumux/material/idealgas.hh> | ||
20 | #include <dumux/common/exceptions.hh> | ||
21 | |||
22 | #include "iapws/common.hh" | ||
23 | #include "iapws/region1.hh" | ||
24 | #include "iapws/region2.hh" | ||
25 | #include "iapws/region4.hh" | ||
26 | |||
27 | #include <dumux/material/components/base.hh> | ||
28 | #include <dumux/material/components/liquid.hh> | ||
29 | #include <dumux/material/components/gas.hh> | ||
30 | |||
31 | namespace Dumux { | ||
32 | namespace Components { | ||
33 | |||
34 | /*! | ||
35 | * \ingroup Components | ||
36 | * \brief Material properties of pure water \f$H_2O\f$. | ||
37 | * \tparam Scalar The type used for scalar values | ||
38 | * | ||
39 | * See: | ||
40 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
41 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
42 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
43 | */ | ||
44 | template <class Scalar> | ||
45 | class H2O | ||
46 | : public Components::Base<Scalar, H2O<Scalar> > | ||
47 | , public Components::Liquid<Scalar, H2O<Scalar> > | ||
48 | , public Components::Gas<Scalar, H2O<Scalar> > | ||
49 | { | ||
50 | |||
51 | using Common = IAPWS::Common<Scalar>; | ||
52 | using Region1 = IAPWS::Region1<Scalar>; | ||
53 | using Region2 = IAPWS::Region2<Scalar>; | ||
54 | using Region4 = IAPWS::Region4<Scalar>; | ||
55 | |||
56 | // specific gas constant of water | ||
57 | static constexpr Scalar Rs = Common::Rs; | ||
58 | |||
59 | public: | ||
60 | /*! | ||
61 | * \brief A human readable name for the water. | ||
62 | */ | ||
63 | static std::string name() | ||
64 |
54/110✓ Branch 3 taken 95 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 95 times.
✗ Branch 7 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✓ Branch 19 taken 34 times.
✗ Branch 20 not taken.
✓ Branch 22 taken 34 times.
✗ Branch 23 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
✓ Branch 33 taken 1 times.
✗ Branch 34 not taken.
✓ Branch 35 taken 57 times.
✓ Branch 36 taken 1 times.
✓ Branch 37 taken 1 times.
✓ Branch 38 taken 57 times.
✓ Branch 39 taken 1 times.
✓ Branch 40 taken 1 times.
✗ Branch 41 not taken.
✓ Branch 42 taken 1 times.
✗ Branch 43 not taken.
✗ Branch 44 not taken.
✓ Branch 45 taken 2 times.
✗ Branch 46 not taken.
✗ Branch 47 not taken.
✓ Branch 48 taken 2 times.
✗ Branch 49 not taken.
✓ Branch 51 taken 2 times.
✗ Branch 52 not taken.
✓ Branch 54 taken 2 times.
✗ Branch 55 not taken.
✓ Branch 57 taken 1 times.
✗ Branch 58 not taken.
✗ Branch 59 not taken.
✓ Branch 60 taken 1 times.
✗ Branch 61 not taken.
✗ Branch 62 not taken.
✓ Branch 63 taken 1 times.
✗ Branch 64 not taken.
✓ Branch 66 taken 1 times.
✓ Branch 67 taken 2 times.
✗ Branch 68 not taken.
✓ Branch 69 taken 1 times.
✓ Branch 70 taken 2 times.
✗ Branch 71 not taken.
✓ Branch 72 taken 1 times.
✗ Branch 73 not taken.
✓ Branch 75 taken 1 times.
✗ Branch 76 not taken.
✓ Branch 78 taken 1 times.
✗ Branch 79 not taken.
✓ Branch 81 taken 1 times.
✗ Branch 82 not taken.
✓ Branch 84 taken 1 times.
✗ Branch 85 not taken.
✓ Branch 87 taken 1 times.
✗ Branch 88 not taken.
✓ Branch 90 taken 1 times.
✗ Branch 91 not taken.
✓ Branch 93 taken 1 times.
✗ Branch 94 not taken.
✓ Branch 96 taken 1 times.
✗ Branch 97 not taken.
✓ Branch 99 taken 2 times.
✗ Branch 100 not taken.
✓ Branch 102 taken 2 times.
✗ Branch 103 not taken.
✗ Branch 104 not taken.
✓ Branch 105 taken 2 times.
✗ Branch 106 not taken.
✓ Branch 107 taken 1 times.
✓ Branch 108 taken 1 times.
✓ Branch 109 taken 1 times.
✗ Branch 110 not taken.
✓ Branch 111 taken 1 times.
✓ Branch 112 taken 1 times.
✗ Branch 113 not taken.
✓ Branch 114 taken 1 times.
✓ Branch 115 taken 1 times.
✗ Branch 116 not taken.
✓ Branch 117 taken 2 times.
✗ Branch 118 not taken.
✓ Branch 120 taken 1 times.
✗ Branch 121 not taken.
✓ Branch 123 taken 1 times.
✗ Branch 124 not taken.
✓ Branch 126 taken 1 times.
✗ Branch 127 not taken.
✓ Branch 129 taken 1 times.
✗ Branch 130 not taken.
✓ Branch 132 taken 1 times.
✗ Branch 133 not taken.
✓ Branch 135 taken 1 times.
✗ Branch 136 not taken.
✓ Branch 138 taken 1 times.
✗ Branch 139 not taken.
✓ Branch 141 taken 1 times.
✗ Branch 142 not taken.
✓ Branch 144 taken 1 times.
✗ Branch 145 not taken.
✓ Branch 147 taken 1 times.
✗ Branch 148 not taken.
✓ Branch 150 taken 1 times.
✗ Branch 151 not taken.
|
798938 | { return "H2O"; } |
65 | |||
66 | /*! | ||
67 | * \brief The molar mass in \f$\mathrm{[kg/mol]}\f$ of water. | ||
68 | */ | ||
69 | static constexpr Scalar molarMass() | ||
70 | { return Common::molarMass; } | ||
71 | |||
72 | /*! | ||
73 | * \brief The acentric factor \f$\mathrm{[-]}\f$ of water. | ||
74 | */ | ||
75 | static constexpr Scalar acentricFactor() | ||
76 | { return Common::acentricFactor; } | ||
77 | |||
78 | /*! | ||
79 | * \brief Returns the critical temperature \f$\mathrm{[K]}\f$ of water | ||
80 | */ | ||
81 | static constexpr Scalar criticalTemperature() | ||
82 | { return Common::criticalTemperature; } | ||
83 | |||
84 | /*! | ||
85 | * \brief Returns the critical pressure \f$\mathrm{[Pa]}\f$ of water. | ||
86 | */ | ||
87 | static constexpr Scalar criticalPressure() | ||
88 | { return Common::criticalPressure; } | ||
89 | |||
90 | /*! | ||
91 | * \brief Returns the molar volume \f$\mathrm{[m^3/mol]}\f$ of water at the critical point | ||
92 | */ | ||
93 | static constexpr Scalar criticalMolarVolume() | ||
94 | { return Common::criticalMolarVolume; } | ||
95 | |||
96 | /*! | ||
97 | * \brief Returns the temperature \f$\mathrm{[K]}\f$ at water's triple point. | ||
98 | */ | ||
99 | static constexpr Scalar tripleTemperature() | ||
100 | { return Common::tripleTemperature; } | ||
101 | |||
102 | /*! | ||
103 | * \brief Returns the pressure \f$\mathrm{[Pa]}\f$ at water's triple point. | ||
104 | */ | ||
105 | static constexpr Scalar triplePressure() | ||
106 | { return Common::triplePressure; } | ||
107 | |||
108 | /*! | ||
109 | * \brief The vapor pressure in \f$\mathrm{[Pa]}\f$ of pure water | ||
110 | * at a given temperature. | ||
111 | * | ||
112 | *\param T temperature of component in \f$\mathrm{[K]}\f$ | ||
113 | * | ||
114 | * See: | ||
115 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
116 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
117 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
118 | */ | ||
119 | 230553591 | static Scalar vaporPressure(Scalar T) | |
120 | { | ||
121 | using std::min; | ||
122 | using std::max; | ||
123 |
2/2✓ Branch 0 taken 2 times.
✓ Branch 1 taken 230553589 times.
|
230553591 | T = min(T, criticalTemperature()); |
124 |
2/2✓ Branch 0 taken 1061883 times.
✓ Branch 1 taken 229491708 times.
|
230553591 | T = max(T,tripleTemperature()); |
125 | |||
126 | 230553592 | return Region4::saturationPressure(T); | |
127 | } | ||
128 | |||
129 | /*! | ||
130 | * \brief The vapor temperature in \f$\mathrm{[K]}\f$ of pure water | ||
131 | * at a given pressure. | ||
132 | * | ||
133 | *\param pressure pressure in \f$\mathrm{[Pa]}\f$ | ||
134 | * | ||
135 | * See: | ||
136 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
137 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
138 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
139 | */ | ||
140 | 101966 | static Scalar vaporTemperature(Scalar pressure) | |
141 | { | ||
142 | using std::min; | ||
143 | using std::max; | ||
144 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 101966 times.
|
101966 | pressure = min(pressure, criticalPressure()); |
145 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 101966 times.
|
101966 | pressure = max(pressure, triplePressure()); |
146 | |||
147 | 101966 | return Region4::vaporTemperature(pressure); | |
148 | } | ||
149 | |||
150 | /*! | ||
151 | * \brief Specific enthalpy of water steam \f$\mathrm{[J/kg]}\f$. | ||
152 | * | ||
153 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
154 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
155 | * | ||
156 | * See: | ||
157 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
158 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
159 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
160 | */ | ||
161 | 3278701 | static const Scalar gasEnthalpy(Scalar temperature, | |
162 | Scalar pressure) | ||
163 | { | ||
164 |
6/12✓ Branch 1 taken 3278701 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3278701 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 3278701 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 3278701 times.
✓ Branch 11 taken 97768 times.
✓ Branch 12 taken 3180933 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
|
6557402 | Region2::checkValidityRange(temperature, pressure, "Enthalpy"); |
165 | |||
166 | // regularization | ||
167 |
2/2✓ Branch 0 taken 97768 times.
✓ Branch 1 taken 3180933 times.
|
3278701 | if (pressure < triplePressure() - 100) { |
168 | // We assume an ideal gas for low pressures to avoid the | ||
169 | // 0/0 for the gas enthalpy at very low pressures. The | ||
170 | // enthalpy of an ideal gas does not exhibit any | ||
171 | // dependence on pressure, so we can just return the | ||
172 | // specific enthalpy at the point of regularization, i.e. | ||
173 | // the triple pressure - 100Pa | ||
174 | 195536 | return enthalpyRegion2_(temperature, triplePressure() - 100); | |
175 | } | ||
176 | 3180933 | Scalar pv = vaporPressure(temperature); | |
177 |
2/2✓ Branch 0 taken 673373 times.
✓ Branch 1 taken 2507560 times.
|
3180933 | if (pressure > pv) { |
178 | // the pressure is too high, in this case we use the slope | ||
179 | // of the enthalpy at the vapor pressure to regularize | ||
180 | 673373 | Scalar dh_dp = | |
181 | 1346746 | Rs*temperature* | |
182 | 673373 | Region2::tau(temperature)* | |
183 | 673373 | Region2::dPi_dp(pv)* | |
184 | 673373 | Region2::ddGamma_dTaudPi(temperature, pv); | |
185 | |||
186 | return | ||
187 | 1346746 | enthalpyRegion2_(temperature, pv) + | |
188 | 673373 | (pressure - pv)*dh_dp; | |
189 | } | ||
190 | |||
191 | 5015120 | return enthalpyRegion2_(temperature, pressure); | |
192 | } | ||
193 | |||
194 | /*! | ||
195 | * \brief Specific enthalpy of liquid water \f$\mathrm{[J/kg]}\f$. | ||
196 | * | ||
197 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
198 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
199 | * | ||
200 | * See: | ||
201 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
202 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
203 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
204 | */ | ||
205 | 4858064 | static const Scalar liquidEnthalpy(Scalar temperature, | |
206 | Scalar pressure) | ||
207 | { | ||
208 |
4/10✓ Branch 1 taken 4858064 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4858064 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 4858064 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 4858064 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
9716128 | Region1::checkValidityRange(temperature, pressure, "Enthalpy"); |
209 | |||
210 | // regularization | ||
211 | 4858064 | Scalar pv = vaporPressure(temperature); | |
212 |
2/2✓ Branch 0 taken 81496 times.
✓ Branch 1 taken 4776568 times.
|
4858064 | if (pressure < pv) { |
213 | // the pressure is too low, in this case we use the slope | ||
214 | // of the enthalpy at the vapor pressure to regularize | ||
215 | 81496 | Scalar dh_dp = | |
216 | 162992 | Rs * temperature* | |
217 | 81496 | Region1::tau(temperature)* | |
218 | 81496 | Region1::dPi_dp(pv)* | |
219 | 81496 | Region1::ddGamma_dTaudPi(temperature, pv); | |
220 | |||
221 | return | ||
222 | 162992 | enthalpyRegion1_(temperature, pv) + | |
223 | 81496 | (pressure - pv)*dh_dp; | |
224 | } | ||
225 | |||
226 | 9553136 | return enthalpyRegion1_(temperature, pressure); | |
227 | } | ||
228 | |||
229 | /*! | ||
230 | * \brief Specific isobaric heat capacity of water steam \f$\mathrm{[J/(kg*K)}\f$. | ||
231 | * | ||
232 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
233 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
234 | * | ||
235 | * See: | ||
236 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
237 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
238 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
239 | */ | ||
240 | 2764815 | static const Scalar gasHeatCapacity(Scalar temperature, | |
241 | Scalar pressure) | ||
242 | { | ||
243 |
6/12✓ Branch 1 taken 2764815 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2764815 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 2764815 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 2764815 times.
✓ Branch 11 taken 97741 times.
✓ Branch 12 taken 2667074 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
|
5529630 | Region2::checkValidityRange(temperature, pressure, "Heat capacity"); |
244 | |||
245 | // regularization | ||
246 |
2/2✓ Branch 0 taken 97741 times.
✓ Branch 1 taken 2667074 times.
|
2764815 | if (pressure < triplePressure() - 100) { |
247 | 195482 | return heatCap_p_Region2_(temperature, triplePressure() - 100); | |
248 | } | ||
249 | 2667074 | Scalar pv = vaporPressure(temperature); | |
250 |
2/2✓ Branch 0 taken 364355 times.
✓ Branch 1 taken 2302719 times.
|
2667074 | if (pressure > pv) { |
251 | // the pressure is too high, in this case we use the heat | ||
252 | // cap at the vapor pressure to regularize | ||
253 | 728710 | return heatCap_p_Region2_(temperature, pv); | |
254 | } | ||
255 | 4605438 | return heatCap_p_Region2_(temperature, pressure); | |
256 | } | ||
257 | |||
258 | /*! | ||
259 | * \brief Specific isobaric heat capacity of liquid water \f$\mathrm{[J/(kg*K)]}\f$. | ||
260 | * | ||
261 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
262 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
263 | * | ||
264 | * See: | ||
265 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
266 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
267 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
268 | */ | ||
269 | 2765353 | static const Scalar liquidHeatCapacity(Scalar temperature, | |
270 | Scalar pressure) | ||
271 | { | ||
272 |
4/10✓ Branch 1 taken 2765353 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2765353 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 2765353 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 2765353 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
5530706 | Region1::checkValidityRange(temperature, pressure, "Heat capacity"); |
273 | |||
274 | // regularization | ||
275 | 2765353 | Scalar pv = vaporPressure(temperature); | |
276 |
2/2✓ Branch 0 taken 81252 times.
✓ Branch 1 taken 2684101 times.
|
2765353 | if (pressure < pv) { |
277 | // the pressure is too low, in this case we use the heat cap at the vapor pressure to regularize | ||
278 | 162504 | return heatCap_p_Region1_(temperature, pv); | |
279 | } | ||
280 | |||
281 | 5368202 | return heatCap_p_Region1_(temperature, pressure); | |
282 | } | ||
283 | |||
284 | /*! | ||
285 | * \brief Specific internal energy of liquid water \f$\mathrm{[J/kg]}\f$. | ||
286 | * | ||
287 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
288 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
289 | * | ||
290 | * See: | ||
291 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
292 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
293 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
294 | */ | ||
295 | 1235383 | static const Scalar liquidInternalEnergy(Scalar temperature, | |
296 | Scalar pressure) | ||
297 | { | ||
298 |
4/10✓ Branch 1 taken 1235383 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1235383 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1235383 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 1235383 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
2470766 | Region1::checkValidityRange(temperature, pressure, "Internal energy"); |
299 | |||
300 | // regularization | ||
301 | 1235383 | Scalar pv = vaporPressure(temperature); | |
302 |
2/2✓ Branch 0 taken 217 times.
✓ Branch 1 taken 1235166 times.
|
1235383 | if (pressure < pv) { |
303 | // the pressure is too low, in this case we use the slope | ||
304 | // of the internal energy at the vapor pressure to | ||
305 | // regularize | ||
306 | |||
307 | /* | ||
308 | // calculate the partial derivative of the internal energy | ||
309 | // to the pressure at the vapor pressure. | ||
310 | Scalar tau = Region1::tau(temperature); | ||
311 | Scalar dGamma_dPi = Region1::dGamma_dPi(temperature, pv); | ||
312 | Scalar ddGamma_dTaudPi = Region1::ddGamma_dTaudPi(temperature, pv); | ||
313 | Scalar ddGamma_ddPi = Region1::ddGamma_ddPi(temperature, pv); | ||
314 | Scalar pi = Region1::pi(pv); | ||
315 | Scalar dPi_dp = Region1::dPi_dp(pv); | ||
316 | Scalar du_dp = | ||
317 | Rs*temperature* | ||
318 | (tau*dPi_dp*ddGamma_dTaudPi + dPi_dp*dPi_dp*dGamma_dPi + pi*dPi_dp*ddGamma_ddPi); | ||
319 | */ | ||
320 | |||
321 | // use a straight line for extrapolation. use forward | ||
322 | // differences to calculate the partial derivative to the | ||
323 | // pressure at the vapor pressure | ||
324 | static const Scalar eps = 1e-7; | ||
325 | 217 | Scalar uv = internalEnergyRegion1_(temperature, pv); | |
326 | 217 | Scalar uvPEps = internalEnergyRegion1_(temperature, pv + eps); | |
327 | 217 | Scalar du_dp = (uvPEps - uv)/eps; | |
328 | 217 | return uv + du_dp*(pressure - pv); | |
329 | } | ||
330 | |||
331 | 1235166 | return internalEnergyRegion1_(temperature, pressure); | |
332 | } | ||
333 | |||
334 | /*! | ||
335 | * \brief Specific internal energy of steam and water vapor \f$\mathrm{[J/kg]}\f$. | ||
336 | * | ||
337 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
338 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
339 | * | ||
340 | * See: | ||
341 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
342 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
343 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
344 | */ | ||
345 | 205023 | static Scalar gasInternalEnergy(Scalar temperature, Scalar pressure) | |
346 | { | ||
347 |
5/12✓ Branch 1 taken 205023 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 205023 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 205023 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 205023 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 205023 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
|
410046 | Region2::checkValidityRange(temperature, pressure, "Internal energy"); |
348 | |||
349 | // regularization | ||
350 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 205023 times.
|
205023 | if (pressure < triplePressure() - 100) { |
351 | // We assume an ideal gas for low pressures to avoid the | ||
352 | // 0/0 for the internal energy of gas at very low | ||
353 | // pressures. The enthalpy of an ideal gas does not | ||
354 | // exhibit any dependence on pressure, so we can just | ||
355 | // return the specific enthalpy at the point of | ||
356 | // regularization, i.e. the triple pressure - 100Pa, and | ||
357 | // subtract the work required to change the volume for an | ||
358 | // ideal gas. | ||
359 | return | ||
360 | ✗ | enthalpyRegion2_(temperature, triplePressure() - 100) | |
361 | - | ||
362 | ✗ | Rs*temperature; // = p*v for an ideal gas! | |
363 | } | ||
364 | 205023 | Scalar pv = vaporPressure(temperature); | |
365 |
2/2✓ Branch 0 taken 189 times.
✓ Branch 1 taken 204834 times.
|
205023 | if (pressure > pv) { |
366 | // the pressure is too high, in this case we use the slope | ||
367 | // of the internal energy at the vapor pressure to | ||
368 | // regularize | ||
369 | |||
370 | /* | ||
371 | // calculate the partial derivative of the internal energy | ||
372 | // to the pressure at the vapor pressure. | ||
373 | Scalar tau = Region2::tau(temperature); | ||
374 | Scalar dGamma_dPi = Region2::dGamma_dPi(temperature, pv); | ||
375 | Scalar ddGamma_dTaudPi = Region2::ddGamma_dTaudPi(temperature, pv); | ||
376 | Scalar ddGamma_ddPi = Region2::ddGamma_ddPi(temperature, pv); | ||
377 | Scalar pi = Region2::pi(pv); | ||
378 | Scalar dPi_dp = Region2::dPi_dp(pv); | ||
379 | Scalar du_dp = | ||
380 | Rs*temperature* | ||
381 | (tau*dPi_dp*ddGamma_dTaudPi + dPi_dp*dPi_dp*dGamma_dPi + pi*dPi_dp*ddGamma_ddPi); | ||
382 | |||
383 | // use a straight line for extrapolation | ||
384 | Scalar uv = internalEnergyRegion2_(temperature, pv); | ||
385 | return uv + du_dp*(pressure - pv); | ||
386 | */ | ||
387 | |||
388 | // use a straight line for extrapolation. use backward | ||
389 | // differences to calculate the partial derivative to the | ||
390 | // pressure at the vapor pressure | ||
391 | static const Scalar eps = 1e-7; | ||
392 | 189 | Scalar uv = internalEnergyRegion2_(temperature, pv); | |
393 | 189 | Scalar uvMEps = internalEnergyRegion2_(temperature, pv - eps); | |
394 | 189 | Scalar du_dp = (uv - uvMEps)/eps; | |
395 | 189 | return uv + du_dp*(pressure - pv); | |
396 | } | ||
397 | |||
398 | 204834 | return internalEnergyRegion2_(temperature, pressure); | |
399 | } | ||
400 | |||
401 | /*! | ||
402 | * \brief Specific isochoric heat capacity of liquid water \f$\mathrm{[J/(m^3*K)]}\f$. | ||
403 | * | ||
404 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
405 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
406 | * | ||
407 | * See: | ||
408 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
409 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
410 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
411 | */ | ||
412 | static const Scalar liquidHeatCapacityConstVolume(Scalar temperature, | ||
413 | Scalar pressure) | ||
414 | { | ||
415 | Region1::checkValidityRange(temperature, pressure, "Heat capacity for a constant volume"); | ||
416 | |||
417 | // regularization | ||
418 | Scalar pv = vaporPressure(temperature); | ||
419 | if (pressure < pv) { | ||
420 | // the pressure is too low, in this case we use the heat cap at the vapor pressure to regularize | ||
421 | |||
422 | return heatCap_v_Region1_(temperature, pv); | ||
423 | } | ||
424 | |||
425 | return heatCap_v_Region1_(temperature, pressure); | ||
426 | } | ||
427 | |||
428 | /*! | ||
429 | * \brief Specific isochoric heat capacity of steam and water vapor \f$\mathrm{[J/(kg*K)}\f$. | ||
430 | * | ||
431 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
432 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
433 | * | ||
434 | * See: | ||
435 | * | ||
436 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
437 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
438 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
439 | */ | ||
440 | static Scalar gasHeatCapacityConstVolume(Scalar temperature, Scalar pressure) | ||
441 | { | ||
442 | Region2::checkValidityRange(temperature, pressure, "Heat capacity for a constant volume"); | ||
443 | |||
444 | // regularization | ||
445 | if (pressure < triplePressure() - 100) { | ||
446 | return | ||
447 | heatCap_v_Region2_(temperature, triplePressure() - 100); | ||
448 | } | ||
449 | Scalar pv = vaporPressure(temperature); | ||
450 | if (pressure > pv) { | ||
451 | return heatCap_v_Region2_(temperature, pv); | ||
452 | } | ||
453 | |||
454 | return heatCap_v_Region2_(temperature, pressure); | ||
455 | } | ||
456 | |||
457 | /*! | ||
458 | * \brief Returns true if the gas phase is assumed to be compressible | ||
459 | */ | ||
460 | static constexpr bool gasIsCompressible() | ||
461 | { return true; } | ||
462 | |||
463 | /*! | ||
464 | * \brief Returns true if the liquid phase is assumed to be compressible | ||
465 | */ | ||
466 | static constexpr bool liquidIsCompressible() | ||
467 | { return true; } | ||
468 | |||
469 | /*! | ||
470 | * \brief The density of steam in \f$\mathrm{[kg/m^3]}\f$ at a given pressure and temperature. | ||
471 | * | ||
472 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
473 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
474 | * | ||
475 | * See: | ||
476 | * | ||
477 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
478 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
479 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
480 | */ | ||
481 | 41481744 | static Scalar gasDensity(Scalar temperature, Scalar pressure) | |
482 | { | ||
483 |
6/12✓ Branch 1 taken 41481744 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 41481744 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 41481744 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 41481744 times.
✓ Branch 11 taken 738824 times.
✓ Branch 12 taken 40742920 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
|
82963488 | Region2::checkValidityRange(temperature, pressure, "Density"); |
484 | |||
485 | // regularization | ||
486 |
2/2✓ Branch 0 taken 738824 times.
✓ Branch 1 taken 40742920 times.
|
41481744 | if (pressure < triplePressure() - 100) { |
487 | // We assume an ideal gas for low pressures to avoid the | ||
488 | // 0/0 for the internal energy and enthalpy. | ||
489 | 1477648 | Scalar rho0IAPWS = 1.0/volumeRegion2_(temperature, | |
490 | triplePressure() - 100); | ||
491 | 1477648 | Scalar rho0Id = IdealGas<Scalar>::density(molarMass(), | |
492 | temperature, | ||
493 | triplePressure() - 100); | ||
494 | return | ||
495 | 738824 | rho0IAPWS/rho0Id * | |
496 | 738824 | IdealGas<Scalar>::density(molarMass(), | |
497 | temperature, | ||
498 | 738824 | pressure); | |
499 | } | ||
500 | 40742920 | Scalar pv = vaporPressure(temperature); | |
501 |
2/2✓ Branch 0 taken 6097038 times.
✓ Branch 1 taken 34645882 times.
|
40742920 | if (pressure > pv) { |
502 | // the pressure is too high, in this case we use the slope | ||
503 | // of the density energy at the vapor pressure to | ||
504 | // regularize | ||
505 | |||
506 | // calculate the partial derivative of the specific volume | ||
507 | // to the pressure at the vapor pressure. | ||
508 | 6097038 | const Scalar eps = pv*1e-8; | |
509 | 12194076 | Scalar v0 = volumeRegion2_(temperature, pv); | |
510 | 12194076 | Scalar v1 = volumeRegion2_(temperature, pv + eps); | |
511 | 6097038 | Scalar dv_dp = (v1 - v0)/eps; | |
512 | /* | ||
513 | Scalar pi = Region2::pi(pv); | ||
514 | Scalar dp_dPi = Region2::dp_dPi(pv); | ||
515 | Scalar dGamma_dPi = Region2::dGamma_dPi(temperature, pv); | ||
516 | Scalar ddGamma_ddPi = Region2::ddGamma_ddPi(temperature, pv); | ||
517 | |||
518 | Scalar RT = Rs*temperature; | ||
519 | Scalar dv_dp = | ||
520 | RT/(dp_dPi*pv) | ||
521 | * | ||
522 | (dGamma_dPi + pi*ddGamma_ddPi - v0*dp_dPi/RT); | ||
523 | */ | ||
524 | |||
525 | // calculate the partial derivative of the density to the | ||
526 | // pressure at vapor pressure | ||
527 | 6097038 | Scalar drho_dp = - 1/(v0*v0)*dv_dp; | |
528 | |||
529 | // use a straight line for extrapolation | ||
530 | 6097038 | return 1.0/v0 + (pressure - pv)*drho_dp; | |
531 | } | ||
532 | |||
533 | 69291764 | return 1.0/volumeRegion2_(temperature, pressure); | |
534 | } | ||
535 | |||
536 | /*! | ||
537 | * \brief The molar density of steam in \f$\mathrm{[mol/m^3]}\f$ at a given pressure and temperature. | ||
538 | * | ||
539 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
540 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
541 | * | ||
542 | */ | ||
543 | static Scalar gasMolarDensity(Scalar temperature, Scalar pressure) | ||
544 | 2303502 | { return gasDensity(temperature, pressure)/molarMass(); } | |
545 | |||
546 | /*! | ||
547 | * \brief Returns true if the gas phase is assumed to be ideal | ||
548 | */ | ||
549 | static constexpr bool gasIsIdeal() | ||
550 | { return false; } | ||
551 | |||
552 | /*! | ||
553 | * \brief The pressure of steam in \f$\mathrm{[Pa]}\f$ at a given density and temperature. | ||
554 | * | ||
555 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
556 | * \param density of component in \f$\mathrm{[kg/m^3]}\f$ | ||
557 | * | ||
558 | * See: | ||
559 | * | ||
560 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
561 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
562 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
563 | */ | ||
564 | 2969632 | static Scalar gasPressure(Scalar temperature, Scalar density) | |
565 | { | ||
566 | // We use the newton method for this. For the initial value we | ||
567 | // assume steam to be an ideal gas | ||
568 | 2969632 | Scalar pressure = IdealGas<Scalar>::pressure(temperature, density/molarMass()); | |
569 | 2969632 | Scalar eps = pressure*1e-7; | |
570 | |||
571 | 2969632 | Scalar deltaP = pressure*2; | |
572 | using std::abs; | ||
573 |
8/8✓ Branch 0 taken 11823500 times.
✓ Branch 1 taken 217058 times.
✓ Branch 2 taken 9070926 times.
✓ Branch 3 taken 2752574 times.
✓ Branch 4 taken 9070926 times.
✓ Branch 5 taken 2752574 times.
✓ Branch 6 taken 9070926 times.
✓ Branch 7 taken 2752574 times.
|
12040558 | for (int i = 0; i < 5 && abs(pressure*1e-9) < abs(deltaP); ++i) |
574 | { | ||
575 | 9070926 | const Scalar f = gasDensity(temperature, pressure) - density; | |
576 | |||
577 | Scalar df_dp; | ||
578 | 9070926 | df_dp = gasDensity(temperature, pressure + eps); | |
579 | 9070926 | df_dp -= gasDensity(temperature, pressure - eps); | |
580 | 9070926 | df_dp /= 2*eps; | |
581 | |||
582 | 9070926 | deltaP = - f/df_dp; | |
583 | |||
584 | 9070926 | pressure += deltaP; | |
585 | } | ||
586 | |||
587 | 2969632 | return pressure; | |
588 | } | ||
589 | |||
590 | /*! | ||
591 | * \brief The density of pure water in \f$\mathrm{[kg/m^3]}\f$ at a given pressure and temperature. | ||
592 | * | ||
593 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
594 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
595 | * | ||
596 | * See: | ||
597 | * | ||
598 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
599 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
600 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
601 | */ | ||
602 | 144681088 | static Scalar liquidDensity(Scalar temperature, Scalar pressure) | |
603 | { | ||
604 |
6/10✓ Branch 1 taken 144681088 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 144681088 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 144681082 times.
✓ Branch 8 taken 6 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 144681082 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 6 times.
|
289362176 | Region1::checkValidityRange(temperature, pressure, "Density"); |
605 | |||
606 | // regularization | ||
607 | 144681082 | Scalar pv = vaporPressure(temperature); | |
608 |
2/2✓ Branch 0 taken 7727607 times.
✓ Branch 1 taken 136953475 times.
|
144681082 | if (pressure < pv) { |
609 | // the pressure is too low, in this case we use the slope | ||
610 | // of the density at the vapor pressure to regularize | ||
611 | |||
612 | // calculate the partial derivative of the specific volume | ||
613 | // to the pressure at the vapor pressure. | ||
614 | 7727607 | const Scalar eps = pv*1e-8; | |
615 | 15455214 | Scalar v0 = volumeRegion1_(temperature, pv); | |
616 | 15455214 | Scalar v1 = volumeRegion1_(temperature, pv + eps); | |
617 | 7727607 | Scalar dv_dp = (v1 - v0)/eps; | |
618 | |||
619 | /* | ||
620 | Scalar v0 = volumeRegion1_(temperature, pv); | ||
621 | Scalar pi = Region1::pi(pv); | ||
622 | Scalar dp_dPi = Region1::dp_dPi(pv); | ||
623 | Scalar dGamma_dPi = Region1::dGamma_dPi(temperature, pv); | ||
624 | Scalar ddGamma_ddPi = Region1::ddGamma_ddPi(temperature, pv); | ||
625 | |||
626 | Scalar RT = Rs*temperature; | ||
627 | Scalar dv_dp = | ||
628 | RT/(dp_dPi*pv) | ||
629 | * | ||
630 | (dGamma_dPi + pi*ddGamma_ddPi - v0*dp_dPi/RT); | ||
631 | */ | ||
632 | |||
633 | // calculate the partial derivative of the density to the | ||
634 | // pressure at vapor pressure | ||
635 | 7727607 | Scalar drho_dp = - 1/(v0*v0)*dv_dp; | |
636 | |||
637 | // use a straight line for extrapolation | ||
638 | 7727607 | return 1.0/v0 + (pressure - pv)*drho_dp; | |
639 | } | ||
640 | |||
641 | 273906950 | return 1/volumeRegion1_(temperature, pressure); | |
642 | } | ||
643 | |||
644 | /*! | ||
645 | * \brief The molar density of water in \f$\mathrm{[mol/m^3]}\f$ at a given pressure and temperature. | ||
646 | * | ||
647 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
648 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
649 | * | ||
650 | */ | ||
651 | static Scalar liquidMolarDensity(Scalar temperature, Scalar pressure) | ||
652 | 14 | { return liquidDensity(temperature, pressure)/molarMass(); } | |
653 | |||
654 | /*! | ||
655 | * \brief The pressure of liquid water in \f$\mathrm{[Pa]}\f$ at a given density and | ||
656 | * temperature. | ||
657 | * | ||
658 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
659 | * \param density density of component in \f$\mathrm{[kg/m^3]}\f$ | ||
660 | * | ||
661 | * See: | ||
662 | * | ||
663 | * IAPWS: "Revised Release on the IAPWS Industrial Formulation | ||
664 | * 1997 for the Thermodynamic Properties of Water and Steam", | ||
665 | * http://www.iapws.org/relguide/IF97-Rev.pdf \cite IAPWS1997 | ||
666 | */ | ||
667 | 3999992 | static Scalar liquidPressure(Scalar temperature, Scalar density) | |
668 | { | ||
669 | // We use Brent's method for this, with a pressure range including the surroundings of the | ||
670 | // vapor pressure line | ||
671 | 3999992 | Scalar minPressure = vaporPressure(temperature)/1.11; | |
672 | 3999992 | Scalar maxPressure = 100e6; | |
673 | 3999992 | const auto residualFunction = [&] (const Scalar pressure) { | |
674 |
3/4✗ Branch 3 not taken.
✓ Branch 4 taken 3999992 times.
✓ Branch 6 taken 110416689 times.
✓ Branch 7 taken 7737416 times.
|
122154097 | return liquidDensity(temperature, pressure) - density; |
675 | }; | ||
676 | try | ||
677 | { | ||
678 |
1/2✓ Branch 1 taken 3999992 times.
✗ Branch 2 not taken.
|
3999992 | return findScalarRootBrent(minPressure, maxPressure, residualFunction); |
679 | } | ||
680 | ✗ | catch (const NumericalProblem& e) | |
681 | { | ||
682 | ✗ | DUNE_THROW(NumericalProblem, | |
683 | "searched for pressure(T=" << temperature << ",rho=" << density | ||
684 | <<") in [" << minPressure << ", " << maxPressure << "]: " | ||
685 | << e.what()); | ||
686 | } | ||
687 | } | ||
688 | |||
689 | /*! | ||
690 | * \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of steam. | ||
691 | * | ||
692 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
693 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
694 | * | ||
695 | * We assume pure water vapor here. For water in a mixture of other gaseous | ||
696 | * components, consider the free function h2oGasViscosityInMixture. | ||
697 | * | ||
698 | * We use the IAPWS Formulation, see: | ||
699 | * IAPWS: "Release on the IAPWS Formulation 2008 for the Viscosity | ||
700 | * of Ordinary Water Substance", http://www.iapws.org/relguide/visc.pdf \cite cooper2008 | ||
701 | * This method is only valid if pressure is below or at the vapor | ||
702 | * pressure of water. | ||
703 | */ | ||
704 | 2969840 | static Scalar gasViscosity(Scalar temperature, Scalar pressure) | |
705 | { | ||
706 |
4/10✓ Branch 1 taken 2969840 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2969840 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 2969840 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 2969840 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
5939680 | Region2::checkValidityRange(temperature, pressure, "Viscosity"); |
707 | |||
708 | 2969840 | Scalar rho = gasDensity(temperature, pressure); | |
709 | 2969840 | return Common::viscosity(temperature, rho); | |
710 | } | ||
711 | |||
712 | |||
713 | /*! | ||
714 | * \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of pure water. | ||
715 | * | ||
716 | * \param temperature temperature of component in \f$\mathrm{[K]}\f$ | ||
717 | * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$ | ||
718 | * | ||
719 | * See: | ||
720 | * IAPWS: "Release on the IAPWS Formulation 2008 for the Viscosity | ||
721 | * of Ordinary Water Substance", http://www.iapws.org/relguide/visc.pdf \cite cooper2008 | ||
722 | */ | ||
723 | 4889455 | static Scalar liquidViscosity(Scalar temperature, Scalar pressure) | |
724 | { | ||
725 |
4/10✓ Branch 1 taken 4889455 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4889455 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 4889455 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 4889455 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
9778910 | Region1::checkValidityRange(temperature, pressure, "Viscosity"); |
726 | |||
727 | 4889455 | Scalar rho = liquidDensity(temperature, pressure); | |
728 | 4889455 | return Common::viscosity(temperature, rho); | |
729 | } | ||
730 | |||
731 | /*! | ||
732 | * \brief Thermal conductivity \f$\mathrm{[[W/(m*K)]}\f$ of water (IAPWS) . | ||
733 | * | ||
734 | * Implementation taken from: | ||
735 | * freesteam - IAPWS-IF97 steam tables library | ||
736 | * copyright (C) 2004-2009 John Pye | ||
737 | * | ||
738 | * See: | ||
739 | * IAPWS: "Release on the IAPWS Formulation 2011 for the Thermal Conductivity | ||
740 | * of Ordinary Water Substance", http://www.iapws.org/relguide/ThCond.pdf | ||
741 | * \cite IAPWS_ThCond | ||
742 | * | ||
743 | * \param temperature absolute temperature in \f$\mathrm{[K]}\f$ | ||
744 | * \param pressure of the phase in \f$\mathrm{[Pa]}\f$ | ||
745 | */ | ||
746 | 3597563 | static Scalar liquidThermalConductivity( Scalar temperature, Scalar pressure) | |
747 | { | ||
748 | // Thermal conductivity of water is empirically fit. | ||
749 | // Evaluating that fitting-function outside the area of validity does not make sense. | ||
750 |
6/8✓ Branch 0 taken 3597563 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 3597546 times.
✓ Branch 3 taken 17 times.
✓ Branch 4 taken 1619600 times.
✓ Branch 5 taken 1977946 times.
✓ Branch 6 taken 17 times.
✗ Branch 7 not taken.
|
3597580 | if ( !( (pressure <= 400e6 && (273.15 <= temperature) && (temperature <= 398.15)) |
751 |
5/6✓ Branch 0 taken 1619617 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1619600 times.
✓ Branch 3 taken 17 times.
✓ Branch 4 taken 727400 times.
✓ Branch 5 taken 892200 times.
|
1619617 | || (pressure <= 200e6 && (398.15 < temperature) && (temperature <= 523.15)) |
752 |
4/6✓ Branch 0 taken 727417 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 727400 times.
✓ Branch 3 taken 17 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 727400 times.
|
727417 | || (pressure <= 150e6 && (523.15 < temperature) && (temperature <= 673.15)) |
753 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 17 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
17 | || (pressure <= 100e6 && (673.15 < temperature) && (temperature <= 1073.15)) )) |
754 | { | ||
755 |
9/20✓ Branch 2 taken 17 times.
✗ Branch 3 not taken.
✓ Branch 11 taken 17 times.
✗ Branch 12 not taken.
✓ Branch 16 taken 17 times.
✗ Branch 17 not taken.
✓ Branch 20 taken 17 times.
✗ Branch 21 not taken.
✓ Branch 23 taken 17 times.
✗ Branch 24 not taken.
✓ Branch 26 taken 17 times.
✗ Branch 27 not taken.
✓ Branch 29 taken 17 times.
✗ Branch 30 not taken.
✓ Branch 31 taken 17 times.
✗ Branch 32 not taken.
✓ Branch 35 taken 17 times.
✗ Branch 36 not taken.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
|
204 | DUNE_THROW(NumericalProblem, |
756 | "Evaluating the IAPWS fit function for thermal conductivity outside range of applicability." | ||
757 | "(T=" << temperature << ", p=" << pressure << ")"); | ||
758 | } | ||
759 | |||
760 | 3597546 | Scalar rho = liquidDensity(temperature, pressure); | |
761 | 3597546 | return Common::thermalConductivityIAPWS(temperature, rho); | |
762 | } | ||
763 | |||
764 | /*! | ||
765 | * \brief Thermal conductivity \f$\mathrm{[[W/(m*K)]}\f$ of steam (IAPWS) . | ||
766 | * | ||
767 | * Implementation taken from: | ||
768 | * freesteam - IAPWS-IF97 steam tables library | ||
769 | * copyright (C) 2004-2009 John Pye | ||
770 | * | ||
771 | * See: | ||
772 | * IAPWS: "Release on the IAPWS Formulation 2011 for the Thermal Conductivity | ||
773 | * of Ordinary Water Substance", http://www.iapws.org/relguide/ThCond.pdf | ||
774 | * \cite IAPWS_ThCond | ||
775 | * | ||
776 | * \param temperature absolute temperature in \f$\mathrm{[K]}\f$ | ||
777 | * \param pressure of the phase in \f$\mathrm{[Pa]}\f$ | ||
778 | */ | ||
779 | 3007380 | static Scalar gasThermalConductivity(const Scalar temperature, const Scalar pressure) | |
780 | { | ||
781 | // Thermal conductivity of water is empirically fit. | ||
782 | // Evaluating that fitting-function outside the area of validity does not make sense. | ||
783 |
4/8✓ Branch 0 taken 3007380 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 3007380 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 1619600 times.
✓ Branch 5 taken 1387780 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
3007380 | if ( !( (pressure <= 400e6 && (273.15 <= temperature) && (temperature <= 398.15)) |
784 |
4/6✓ Branch 0 taken 1619600 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1619600 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 727400 times.
✓ Branch 5 taken 892200 times.
|
1619600 | || (pressure <= 200e6 && (398.15 < temperature) && (temperature <= 523.15)) |
785 |
3/6✓ Branch 0 taken 727400 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 727400 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 727400 times.
|
727400 | || (pressure <= 150e6 && (523.15 < temperature) && (temperature <= 673.15)) |
786 | ✗ | || (pressure <= 100e6 && (673.15 < temperature) && (temperature <= 1073.15)) )) | |
787 | { | ||
788 | ✗ | DUNE_THROW(NumericalProblem, | |
789 | "Evaluating the IAPWS fit function for thermal conductivity outside range of applicability." | ||
790 | "(T=" << temperature << ", p=" << pressure << ")"); | ||
791 | } | ||
792 | |||
793 | 3007380 | Scalar rho = gasDensity(temperature, pressure); | |
794 | 3007380 | return Common::thermalConductivityIAPWS(temperature, rho); | |
795 | } | ||
796 | |||
797 | private: | ||
798 | // the unregularized specific enthalpy for liquid water | ||
799 | static constexpr Scalar enthalpyRegion1_(Scalar temperature, Scalar pressure) | ||
800 | { | ||
801 | return | ||
802 | 4776568 | Region1::tau(temperature) * | |
803 | 4858064 | Region1::dGamma_dTau(temperature, pressure) * | |
804 | 4858064 | Rs*temperature; | |
805 | } | ||
806 | |||
807 | // the unregularized specific isobaric heat capacity | ||
808 | static constexpr Scalar heatCap_p_Region1_(Scalar temperature, Scalar pressure) | ||
809 | { | ||
810 | return | ||
811 | 2765353 | - Region1::tau(temperature) * Region1::tau(temperature) * | |
812 | 2765353 | Region1::ddGamma_ddTau(temperature, pressure) * | |
813 | 2765353 | Rs; | |
814 | } | ||
815 | |||
816 | // the unregularized specific isochoric heat capacity | ||
817 | static Scalar heatCap_v_Region1_(Scalar temperature, Scalar pressure) | ||
818 | { | ||
819 | Scalar tau = Region1::tau(temperature); | ||
820 | Scalar num = Region1::dGamma_dPi(temperature, pressure) - tau * Region1::ddGamma_dTaudPi(temperature, pressure); | ||
821 | Scalar diff = num * num / Region1::ddGamma_ddPi(temperature, pressure); | ||
822 | |||
823 | return | ||
824 | - tau * tau * | ||
825 | Region1::ddGamma_ddTau(temperature, pressure) * Rs + | ||
826 | diff; | ||
827 | } | ||
828 | |||
829 | // the unregularized specific internal energy for liquid water | ||
830 | 1235600 | static constexpr Scalar internalEnergyRegion1_(Scalar temperature, Scalar pressure) | |
831 | { | ||
832 | return | ||
833 | 1235600 | Rs * temperature * | |
834 | 2471200 | ( Region1::tau(temperature)*Region1::dGamma_dTau(temperature, pressure) - | |
835 | 2471200 | Region1::pi(pressure)*Region1::dGamma_dPi(temperature, pressure)); | |
836 | } | ||
837 | |||
838 | // the unregularized specific volume for liquid water | ||
839 | static constexpr Scalar volumeRegion1_(Scalar temperature, Scalar pressure) | ||
840 | { | ||
841 | return | ||
842 | 144681082 | Region1::pi(pressure)* | |
843 | 144681082 | Region1::dGamma_dPi(temperature, pressure) * | |
844 | 144681082 | Rs * temperature / pressure; | |
845 | } | ||
846 | |||
847 | // the unregularized specific enthalpy for steam | ||
848 | static constexpr Scalar enthalpyRegion2_(Scalar temperature, Scalar pressure) | ||
849 | { | ||
850 | return | ||
851 | 2605328 | Region2::tau(temperature) * | |
852 | 3278701 | Region2::dGamma_dTau(temperature, pressure) * | |
853 | 3278701 | Rs*temperature; | |
854 | } | ||
855 | |||
856 | // the unregularized specific internal energy for steam | ||
857 | 205212 | static constexpr Scalar internalEnergyRegion2_(Scalar temperature, Scalar pressure) | |
858 | { | ||
859 | return | ||
860 | 205212 | Rs * temperature * | |
861 | 410424 | ( Region2::tau(temperature)*Region2::dGamma_dTau(temperature, pressure) - | |
862 | 410424 | Region2::pi(pressure)*Region2::dGamma_dPi(temperature, pressure)); | |
863 | } | ||
864 | |||
865 | // the unregularized specific isobaric heat capacity | ||
866 | static constexpr Scalar heatCap_p_Region2_(Scalar temperature, Scalar pressure) | ||
867 | { | ||
868 | return | ||
869 | 2764815 | - Region2::tau(temperature) * Region2::tau(temperature) * | |
870 | 2764815 | Region2::ddGamma_ddTau(temperature, pressure) * | |
871 | 2764815 | Rs; | |
872 | } | ||
873 | |||
874 | // the unregularized specific isochoric heat capacity | ||
875 | static Scalar heatCap_v_Region2_(Scalar temperature, Scalar pressure) | ||
876 | { | ||
877 | Scalar tau = Region2::tau(temperature); | ||
878 | Scalar pi = Region2::pi(pressure); | ||
879 | Scalar num = 1 + pi * Region2::dGamma_dPi(temperature, pressure) + tau * pi * Region2::ddGamma_dTaudPi(temperature, pressure); | ||
880 | Scalar diff = num * num / (1 - pi * pi * Region2::ddGamma_ddPi(temperature, pressure)); | ||
881 | return | ||
882 | - tau * tau * | ||
883 | Region2::ddGamma_ddTau(temperature, pressure) * Rs | ||
884 | - diff; | ||
885 | } | ||
886 | |||
887 | // the unregularized specific volume for steam | ||
888 | static constexpr Scalar volumeRegion2_(Scalar temperature, Scalar pressure) | ||
889 | { | ||
890 | return | ||
891 | 40742920 | Region2::pi(pressure)* | |
892 | 41481744 | Region2::dGamma_dPi(temperature, pressure) * | |
893 | 41481744 | Rs * temperature / pressure; | |
894 | } | ||
895 | }; | ||
896 | |||
897 | template <class Scalar> | ||
898 | struct IsAqueous<H2O<Scalar>> : public std::true_type {}; | ||
899 | |||
900 | } // end namespace Components | ||
901 | |||
902 | namespace FluidSystems::Detail { | ||
903 | // viscosity according to Reid, R.C. | ||
904 | template <class Scalar> | ||
905 | 59043789 | Scalar viscosityReid_(Scalar temperature) | |
906 | { | ||
907 | 59043789 | constexpr Scalar tc = 647.3; | |
908 | using std::max; | ||
909 |
2/2✓ Branch 0 taken 14 times.
✓ Branch 1 taken 59043775 times.
|
59043789 | const Scalar tr = max(temperature/tc, 1e-8); |
910 | |||
911 | 59043789 | const Scalar fp0 = 1.0 + 0.221*(0.96 + 0.1*(tr - 0.7)); | |
912 | 59043789 | constexpr Scalar xi = 3.334e-3; | |
913 | using std::pow; | ||
914 | using std::exp; | ||
915 | 118087578 | const Scalar eta_xi = (0.807*pow(tr, 0.618) - 0.357*exp((-0.449)*tr) | |
916 | 59043789 | + 0.34*exp((-4.058)*tr) + 0.018)*fp0; | |
917 | |||
918 | 59043789 | return 1.0e-7*eta_xi/xi; | |
919 | } | ||
920 | |||
921 | // viscosity according to Nagel, T. et al. | ||
922 | template <class Scalar> | ||
923 | Scalar viscosityNagel_(Scalar temperature) | ||
924 | { | ||
925 | 25428 | constexpr Scalar a1 = -4.4189440e-6; | |
926 | 25428 | constexpr Scalar a2 = 4.6876380e-8; | |
927 | 25428 | constexpr Scalar a3 = -5.3894310e-12; | |
928 | 25428 | constexpr Scalar a4 = 3.2028560e-16; | |
929 | 25428 | constexpr Scalar a5 = 4.9191790e-22; | |
930 | |||
931 | 25428 | return a1 + a2*temperature + a3*temperature*temperature | |
932 | 25428 | + a4*temperature*temperature*temperature | |
933 | 25428 | + a5*temperature*temperature*temperature*temperature; | |
934 | } | ||
935 | } // end FluidSystems::Detail | ||
936 | namespace FluidSystems { | ||
937 | /*! | ||
938 | * \brief The dynamic viscosity \f$\mathrm{[Pa*s]}\f$ of steam in a gas mixture. | ||
939 | * | ||
940 | * \param temperature temperature in \f$\mathrm{[K]}\f$ | ||
941 | * \param pressure pressure | ||
942 | * | ||
943 | * We assume here that water is in mixture with other gaseous components. | ||
944 | * For pure water, use the gasViscosity function of Components::H2O. | ||
945 | * | ||
946 | * We apply two different laws depending on the gas temperature. | ||
947 | * | ||
948 | * For temperatures below 480 K see: | ||
949 | * "Reid, R.C., Prausnitz, J.M., Poling, B.E.: The Properties of | ||
950 | * Gases and Liquids (1987)" | ||
951 | * Lucas corresponding states method | ||
952 | * https://www.osti.gov/scitech/biblio/6504847 \cite reid1987 | ||
953 | * | ||
954 | * For temperatures above 500 K see: | ||
955 | * Nagel, T. et al.: THC-Processes (2018) | ||
956 | * https://doi.org/10.1007/978-3-319-68225-9_12 | ||
957 | * | ||
958 | * In the range 480 - 500 K, we interpolate between the two laws. | ||
959 | */ | ||
960 | template <class Scalar> | ||
961 | 59069112 | Scalar h2oGasViscosityInMixture(Scalar temperature, Scalar pressure) | |
962 | { | ||
963 |
2/2✓ Branch 0 taken 59043684 times.
✓ Branch 1 taken 25428 times.
|
59069112 | if (temperature < 480.0) |
964 | 59043684 | return Detail::viscosityReid_(temperature); | |
965 | |||
966 |
2/2✓ Branch 0 taken 25323 times.
✓ Branch 1 taken 105 times.
|
25428 | else if (temperature > 500.0) |
967 | 25323 | return Detail::viscosityNagel_(temperature); | |
968 | |||
969 | else // interpolate | ||
970 | { | ||
971 | 105 | const Scalar op = (500.0 - temperature)/20.0; | |
972 | |||
973 | 105 | return op*Detail::viscosityReid_(temperature) | |
974 | 105 | + (1.0 - op)*Detail::viscosityNagel_(temperature); | |
975 | } | ||
976 | } | ||
977 | } // end namespace FluidSystems | ||
978 | } // end namespace Dumux | ||
979 | |||
980 | #endif | ||
981 |