3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
region2.hh
Go to the documentation of this file.
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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
19
30#ifndef DUMUX_IAPWS_REGION2_HH
31#define DUMUX_IAPWS_REGION2_HH
32
33#include <cmath>
34#include <iostream>
36
37namespace Dumux {
38namespace IAPWS {
39
50template <class Scalar>
52{
53public:
62 static void checkValidityRange(Scalar temperature, Scalar pressure,
63 const std::string& propertyName = "This property")
64 {
65 // actually this is:
66 /* (273.15 <= temperature && temperature <= 623.15 && pressure <= vaporPressure(temperature)) ||
67 (623.15 < temperature && temperature <= 863.15 && pressure <= auxPressure(temperature)) ||
68 (863.15 < temperature && temperature <= 1073.15 && pressure < 100e6); */
69 if ((temperature <= 623.15 && pressure <= 100e6) ||
70 (temperature > 623.15 && temperature <= 1073.15 && pressure <= 16.532e6 ))
71 return;
72
73 DUNE_THROW(NumericalProblem,
74 propertyName << " of steam is only implemented for temperatures below 623.15K and "
75 "pressures below 100MPa. (T=" << temperature << ", p=" << pressure << ")");
76 }
77
83 static constexpr Scalar tau(Scalar temperature)
84 { return 540.0 / temperature; }
85
92 static constexpr Scalar dTau_dt(Scalar temperature)
93 { return - 540.0 / (temperature*temperature); }
94
100 static constexpr Scalar pi(Scalar pressure)
101 { return pressure / 1e6; }
102
109 static constexpr Scalar dPi_dp(Scalar pressure)
110 { return 1.0 / 1e6; }
111
118 static constexpr Scalar dp_dPi(Scalar pressure)
119 { return 1e6; }
120
132 static Scalar gamma(Scalar temperature, Scalar pressure)
133 {
134 Scalar tau = tau(temperature); /* reduced temperature */
135 Scalar pi = pi(pressure); /* reduced pressure */
136
137 Scalar result;
138
139 // ideal gas part
140 using std::pow;
141 result = ln(pi);
142 for (int i = 0; i < 9; ++i)
143 result += n_g(i)*pow(tau, J_g(i));
144
145 // residual part
146 for (int i = 0; i < 43; ++i)
147 result += n_r(i)*
148 pow(pi, I_r(i))*
149 pow(tau - 0.5, J_r(i));
150 return result;
151 }
152
165 static Scalar dGamma_dTau(Scalar temperature, Scalar pressure)
166 {
167 Scalar tau_ = tau(temperature); /* reduced temperature */
168 Scalar pi_ = pi(pressure); /* reduced pressure */
169
170 // ideal gas part
171 using std::pow;
172 Scalar result = 0;
173 for (int i = 0; i < 9; i++) {
174 result += n_g(i) *
175 J_g(i) *
176 pow(tau_, J_g(i) - 1);
177 }
178
179 // residual part
180 for (int i = 0; i < 43; i++) {
181 result += n_r(i) *
182 pow(pi_, I_r(i)) *
183 J_r(i) *
184 pow(tau_ - 0.5, J_r(i) - 1);
185 }
186
187 return result;
188 }
189
202 static Scalar dGamma_dPi(Scalar temperature, Scalar pressure)
203 {
204 Scalar tau_ = tau(temperature); /* reduced temperature */
205 Scalar pi_ = pi(pressure); /* reduced pressure */
206
207 // ideal gas part
208 Scalar result = 1/pi_;
209
210 // residual part
211 using std::pow;
212 for (int i = 0; i < 43; i++) {
213 result += n_r(i) *
214 I_r(i) *
215 pow(pi_, I_r(i) - 1) *
216 pow(tau_ - 0.5, J_r(i));
217 }
218
219 return result;
220 }
221
234 static Scalar ddGamma_dTaudPi(Scalar temperature, Scalar pressure)
235 {
236 Scalar tau_ = tau(temperature); /* reduced temperature */
237 Scalar pi_ = pi(pressure); /* reduced pressure */
238
239 // ideal gas part
240 Scalar result = 0;
241
242 // residual part
243 using std::pow;
244 for (int i = 0; i < 43; i++) {
245 result += n_r(i) *
246 I_r(i) *
247 J_r(i) *
248 pow(pi_, I_r(i) - 1) *
249 pow(tau_ - 0.5, J_r(i) - 1);
250 }
251
252 return result;
253 }
254
267 static Scalar ddGamma_ddPi(Scalar temperature, Scalar pressure)
268 {
269 Scalar tau_ = tau(temperature); /* reduced temperature */
270 Scalar pi_ = pi(pressure); /* reduced pressure */
271
272 // ideal gas part
273 Scalar result = -1/(pi_*pi_);
274
275 // residual part
276 using std::pow;
277 for (int i = 0; i < 43; i++) {
278 result += n_r(i) *
279 I_r(i) *
280 (I_r(i) - 1) *
281 pow(pi_, I_r(i) - 2) *
282 pow(tau_ - 0.5, J_r(i));
283 }
284
285 return result;
286 }
287
300 static Scalar ddGamma_ddTau(Scalar temperature, Scalar pressure)
301 {
302 Scalar tau_ = tau(temperature); /* reduced temperature */
303 Scalar pi_ = pi(pressure); /* reduced pressure */
304
305 // ideal gas part
306 using std::pow;
307 Scalar result = 0;
308 for (int i = 0; i < 9; i++) {
309 result += n_g(i) *
310 J_g(i) *
311 (J_g(i) - 1) *
312 pow(tau_, J_g(i) - 2);
313 }
314
315 // residual part
316 using std::pow;
317 for (int i = 0; i < 43; i++) {
318 result += n_r(i) *
319 pow(pi_, I_r(i)) *
320 J_r(i) *
321 (J_r(i) - 1.) *
322 pow(tau_ - 0.5, J_r(i) - 2.);
323 }
324
325 return result;
326 }
327
328private:
329 static Scalar n_g(int i)
330 {
331 constexpr const Scalar n[9] = {
332 -0.96927686500217e1, 0.10086655968018e2, -0.56087911283020e-2,
333 0.71452738081455e-1, -0.40710498223928, 0.14240819171444e1,
334 -0.43839511319450e1, -0.28408632460772, 0.21268463753307e-1
335 };
336 return n[i];
337 }
338
339 static Scalar n_r(int i)
340 {
341 constexpr const Scalar n[43] = {
342 -0.17731742473213e-2, -0.17834862292358e-1, -0.45996013696365e-1,
343 -0.57581259083432e-1, -0.50325278727930e-1, -0.33032641670203e-4,
344 -0.18948987516315e-3, -0.39392777243355e-2, -0.43797295650573e-1,
345 -0.26674547914087e-4, 0.20481737692309e-7, 0.43870667284435e-6,
346 -0.32277677238570e-4, -0.15033924542148e-2, -0.40668253562649e-1,
347 -0.78847309559367e-9, 0.12790717852285e-7, 0.48225372718507e-6,
348 0.22922076337661e-5, -0.16714766451061e-10, -0.21171472321355e-2,
349 -0.23895741934104e2, -0.59059564324270e-17, -0.12621808899101e-5,
350 -0.38946842435739e-1, 0.11256211360459e-10, -0.82311340897998e1,
351 0.19809712802088e-7, 0.10406965210174e-18, -0.10234747095929e-12,
352 -0.10018179379511e-8, -0.80882908646985e-10, 0.10693031879409,
353 -0.33662250574171, 0.89185845355421e-24, 0.30629316876232e-12,
354 -0.42002467698208e-5, -0.59056029685639e-25, 0.37826947613457e-5,
355 -0.12768608934681e-14, 0.73087610595061e-28, 0.55414715350778e-16,
356 -0.94369707241210e-6
357 };
358 return n[i];
359 }
360
361 static Scalar I_r(int i)
362 {
363 constexpr const short int I[43] = {
364 1, 1, 1,
365 1, 1, 2,
366 2, 2, 2,
367 2, 3, 3,
368 3, 3, 3,
369 4, 4, 4,
370 5, 6, 6,
371 6, 7, 7,
372 7, 8, 8,
373 9, 10, 10,
374 10, 16, 16,
375 18, 20, 20,
376 20, 21, 22,
377 23, 24, 24,
378 24
379 };
380 return I[i];
381 }
382
383 static Scalar J_g(int i)
384 {
385 constexpr const short int J[9] = {
386 0, 1, -5,
387 -4, -3, -2,
388 -1, 2, 3
389 };
390 return J[i];
391 }
392
393 static Scalar J_r(int i)
394 {
395 constexpr const short int J[43] = {
396 0, 1, 2,
397 3, 6, 1,
398 2, 4, 7,
399 36, 0, 1,
400 3, 6, 35,
401 1, 2, 3,
402 7, 3, 16,
403 35, 0, 11,
404 25, 8, 36,
405 13, 4, 10,
406 14, 29, 50,
407 57, 20, 35,
408 48, 21, 53,
409 39, 26, 40,
410 58
411 };
412 return J[i];
413 }
414
415};
416
417} // end namespace IAPWS
418} // end namespace Dumux
419
420#endif
Some exceptions thrown in DuMux
Definition: adapt.hh:29
std::string temperature() noexcept
I/O name of temperature for equilibrium models.
Definition: name.hh:51
std::string pressure(int phaseIdx) noexcept
I/O name of pressure for multiphase systems.
Definition: name.hh:34
Exception thrown if a fixable numerical problem occurs.
Definition: exceptions.hh:39
Implements the equations for region 2 of the IAPWS '97 formulation.
Definition: region2.hh:52
static Scalar ddGamma_ddTau(Scalar temperature, Scalar pressure)
The second partial derivative of the Gibbs free energy to the normalized temperature for IAPWS region...
Definition: region2.hh:300
static Scalar dGamma_dTau(Scalar temperature, Scalar pressure)
The partial derivative of the Gibbs free energy to the normalized temperature for IAPWS region 2 (i....
Definition: region2.hh:165
static Scalar ddGamma_dTaudPi(Scalar temperature, Scalar pressure)
The partial derivative of the Gibbs free energy to the normalized pressure and to the normalized temp...
Definition: region2.hh:234
static Scalar gamma(Scalar temperature, Scalar pressure)
The Gibbs free energy for IAPWS region 2 (i.e. sub-critical steam) (dimensionless).
Definition: region2.hh:132
static Scalar ddGamma_ddPi(Scalar temperature, Scalar pressure)
The second partial derivative of the Gibbs free energy to the normalized pressure for IAPWS region 2 ...
Definition: region2.hh:267
static void checkValidityRange(Scalar temperature, Scalar pressure, const std::string &propertyName="This property")
Returns true if IAPWS region 2 applies for a (temperature, pressure) pair.
Definition: region2.hh:62
static constexpr Scalar dTau_dt(Scalar temperature)
Returns the derivative of the reduced temperature to the temperature for IAPWS region 2.
Definition: region2.hh:92
static constexpr Scalar pi(Scalar pressure)
Returns the reduced pressure (dimensionless) for IAPWS region 2.
Definition: region2.hh:100
static constexpr Scalar dPi_dp(Scalar pressure)
Returns the derivative of the reduced pressure to the pressure for IAPWS region 2 in .
Definition: region2.hh:109
static constexpr Scalar tau(Scalar temperature)
Returns the reduced temperature (dimensionless) for IAPWS region 2.
Definition: region2.hh:83
static constexpr Scalar dp_dPi(Scalar pressure)
Returns the derivative of the pressure to the reduced pressure for IAPWS region 2 (dimensionless).
Definition: region2.hh:118
static Scalar dGamma_dPi(Scalar temperature, Scalar pressure)
The partial derivative of the Gibbs free energy to the normalized pressure for IAPWS region 2 (i....
Definition: region2.hh:202