Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #pragma once
11 :
12 : #include "SinglePhaseFluidProperties.h"
13 : #include "NewtonInversion.h"
14 : #include <array>
15 :
16 : #pragma GCC diagnostic push
17 : #pragma GCC diagnostic ignored "-Woverloaded-virtual"
18 :
19 : /**
20 : * Water (H2O) fluid properties as a function of pressure (Pa)
21 : * and temperature (K) from IAPWS-IF97:
22 : * Revised Release on the IAPWS Industrial Formulation 1997 for
23 : * the Thermodynamic Properties of Water and Steam.
24 : *
25 : * To avoid iteration in Region 3, the backwards equations from:
26 : * Revised Supplementary Release on Backward Equations for Specific Volume
27 : * as a Function of Pressure and Temperature v(p,T) for Region 3 of the
28 : * IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of Water
29 : * and Steam are implemented.
30 : *
31 : * Water viscosity as a function of density (kg/m^3) and temperature (K) from:
32 : * Release on the IAPWS Formulation 2008 for the Viscosity of Ordinary Water Substance.
33 : *
34 : * Thermal conductivity as a function of density (kg/m^3) and temperature (K) from:
35 : * Revised Release on the IAPS Formulation 1985 for the Thermal Conductivity of
36 : * Ordinary Water Substance
37 : */
38 : class Water97FluidProperties : public SinglePhaseFluidProperties
39 : {
40 : public:
41 : static InputParameters validParams();
42 :
43 : Water97FluidProperties(const InputParameters & parameters);
44 : virtual ~Water97FluidProperties();
45 :
46 : virtual std::string fluidName() const override;
47 :
48 : virtual Real molarMass() const override;
49 :
50 : virtual Real criticalPressure() const override;
51 :
52 : virtual Real criticalTemperature() const override;
53 :
54 : virtual Real criticalDensity() const override;
55 :
56 : virtual Real triplePointPressure() const override;
57 :
58 : virtual Real triplePointTemperature() const override;
59 :
60 : virtual Real rho_from_p_T(Real pressure, Real temperature) const override;
61 : virtual ADReal rho_from_p_T(const ADReal & pressure, const ADReal & temperature) const override;
62 :
63 : virtual Real v_from_p_T(Real pressure, Real temperature) const override;
64 : virtual ADReal v_from_p_T(const ADReal & pressure, const ADReal & temperature) const override;
65 : template <typename T>
66 : T v_from_p_T_template(const T & pressure, const T & temperature) const;
67 : virtual void
68 : v_from_p_T(Real pressure, Real temperature, Real & v, Real & dv_dp, Real & dv_dT) const override;
69 : virtual void v_from_p_T(const ADReal & pressure,
70 : const ADReal & temperature,
71 : ADReal & v,
72 : ADReal & dv_dp,
73 : ADReal & dv_dT) const override;
74 : template <typename T>
75 : void
76 : v_from_p_T_template(const T & pressure, const T & temperature, T & v, T & dv_dp, T & dv_dt) const;
77 :
78 : template <typename T>
79 : T rho_from_p_T_template(const T & pressure, const T & temperature) const;
80 : template <typename T>
81 : void rho_from_p_T_template(
82 : const T & pressure, const T & temperature, T & rho, T & drho_dp, T & drho_dt) const;
83 :
84 : virtual void rho_from_p_T(
85 : Real pressure, Real temperature, Real & rho, Real & drho_dp, Real & drho_dT) const override;
86 : virtual void rho_from_p_T(const ADReal & pressure,
87 : const ADReal & temperature,
88 : ADReal & rho,
89 : ADReal & drho_dp,
90 : ADReal & drho_dT) const override;
91 :
92 : Real p_from_v_e(Real v, Real e) const override;
93 : ADReal p_from_v_e(const ADReal & v, const ADReal & e) const override;
94 : template <typename T>
95 : T p_from_v_e_template(const T & v, const T & e) const;
96 :
97 : virtual Real e_from_p_rho(Real p, Real rho) const override;
98 : virtual ADReal e_from_p_rho(const ADReal & p, const ADReal & rho) const override;
99 : template <typename T>
100 : T e_from_p_rho_template(const T & p, const T & rho) const;
101 : void e_from_p_rho(Real p, Real rho, Real & e, Real & de_dp, Real & de_drho) const override;
102 : void e_from_p_rho(const ADReal & p,
103 : const ADReal & rho,
104 : ADReal & e,
105 : ADReal & de_dp,
106 : ADReal & de_drho) const override;
107 : template <typename T>
108 : void e_from_p_rho_template(const T & p, const T & rho, T & e, T & de_dp, T & de_drho) const;
109 :
110 : virtual Real e_from_p_T(Real pressure, Real temperature) const override;
111 : virtual ADReal e_from_p_T(const ADReal & pressure, const ADReal & temperature) const override;
112 : template <typename T>
113 : T e_from_p_T_template(const T & pressure, const T & temperature) const;
114 :
115 : virtual void
116 : e_from_p_T(Real pressure, Real temperature, Real & e, Real & de_dp, Real & de_dT) const override;
117 : virtual void e_from_p_T(const ADReal & pressure,
118 : const ADReal & temperature,
119 : ADReal & e,
120 : ADReal & de_dp,
121 : ADReal & de_dT) const override;
122 : template <typename T>
123 : void
124 : e_from_p_T_template(const T & pressure, const T & temperature, T & e, T & de_dp, T & de_dT) const;
125 :
126 : ADReal e_from_v_h(const ADReal & v, const ADReal & h) const override;
127 :
128 : Real T_from_v_e(Real v, Real e) const override;
129 : ADReal T_from_v_e(const ADReal & v, const ADReal & e) const override;
130 :
131 : virtual Real c_from_p_T(Real pressure, Real temperature) const override;
132 : virtual ADReal c_from_p_T(const ADReal & pressure, const ADReal & temperature) const override;
133 : template <typename T>
134 : T c_from_p_T_template(const T & pressure, const T & temperature) const;
135 : virtual ADReal c_from_v_e(const ADReal & v, const ADReal & e) const override;
136 :
137 : virtual Real cp_from_p_T(Real pressure, Real temperature) const override;
138 : virtual ADReal cp_from_p_T(const ADReal & pressure, const ADReal & temperature) const override;
139 : template <typename T>
140 : T cp_from_p_T_template(const T & pressure, const T & temperature) const;
141 :
142 : using SinglePhaseFluidProperties::cp_from_p_T;
143 :
144 : virtual ADReal cp_from_v_e(const ADReal & v, const ADReal & e) const override;
145 :
146 : virtual Real cv_from_p_T(Real pressure, Real temperature) const override;
147 : virtual ADReal cv_from_p_T(const ADReal & pressure, const ADReal & temperature) const override;
148 : template <typename T>
149 : T cv_from_p_T_template(const T & pressure, const T & temperature) const;
150 :
151 : virtual ADReal cv_from_v_e(const ADReal & v, const ADReal & e) const override;
152 :
153 : virtual Real mu_from_p_T(Real pressure, Real temperature) const override;
154 : virtual ADReal mu_from_p_T(const ADReal & pressure, const ADReal & temperature) const override;
155 : template <typename T>
156 : T mu_from_p_T_template(const T & pressure, const T & temperature) const;
157 :
158 : virtual void mu_from_p_T(
159 : Real pressure, Real temperature, Real & mu, Real & dmu_dp, Real & dmu_dT) const override;
160 :
161 : virtual Real mu_from_rho_T(Real density, Real temperature) const override;
162 : template <typename T>
163 : T mu_from_rho_T_template(const T & density, const T & temperature) const;
164 :
165 : void mu_from_rho_T(
166 : Real rho, Real temperature, Real drho_dT, Real & mu, Real & dmu_drho, Real & dmu_dT) const;
167 :
168 : ADReal mu_from_v_e(const ADReal & v, const ADReal & e) const override;
169 :
170 : virtual void
171 : rho_mu_from_p_T(Real pressure, Real temperature, Real & rho, Real & mu) const override;
172 :
173 : virtual void rho_mu_from_p_T(Real pressure,
174 : Real temperature,
175 : Real & rho,
176 : Real & drho_dp,
177 : Real & drho_dT,
178 : Real & mu,
179 : Real & dmu_dp,
180 : Real & dmu_dT) const override;
181 :
182 : virtual Real k_from_p_T(Real pressure, Real temperature) const override;
183 : virtual ADReal k_from_p_T(const ADReal & pressure, const ADReal & temperature) const override;
184 : template <typename T>
185 : T k_from_p_T_template(const T & pressure, const T & temperature) const;
186 :
187 : virtual void
188 : k_from_p_T(Real pressure, Real temperature, Real & k, Real & dk_dp, Real & dk_dT) const override;
189 :
190 : virtual Real k_from_rho_T(Real density, Real temperature) const override;
191 : template <typename T>
192 : T k_from_rho_T_template(const T & density, const T & temperature) const;
193 :
194 : Real k_from_v_e(Real v, Real e) const override;
195 : ADReal k_from_v_e(const ADReal & v, const ADReal & e) const override;
196 : template <typename T>
197 : T k_from_v_e_template(const T & v, const T & e) const;
198 :
199 : propfuncWithDefinitionOverride(s, p, T);
200 :
201 : virtual Real h_from_p_T(Real pressure, Real temperature) const override;
202 : virtual ADReal h_from_p_T(const ADReal & pressure, const ADReal & temperature) const override;
203 : template <typename T>
204 : T h_from_p_T_template(const T & pressure, const T & temperature) const;
205 :
206 : virtual void
207 : h_from_p_T(Real pressure, Real temperature, Real & h, Real & dh_dp, Real & dh_dT) const override;
208 : virtual void h_from_p_T(const ADReal & pressure,
209 : const ADReal & temperature,
210 : ADReal & h,
211 : ADReal & dh_dp,
212 : ADReal & dh_dT) const override;
213 : template <typename T>
214 : void
215 : h_from_p_T_template(const T & pressure, const T & temperature, T & h, T & dh_dp, T & dh_dT) const;
216 :
217 : virtual Real s_from_h_p(Real enthalpy, Real pressure) const override;
218 : virtual ADReal s_from_h_p(const ADReal & enthalpy, const ADReal & pressure) const override;
219 : void virtual s_from_h_p(
220 : Real enthalpy, Real pressure, Real & s, Real & ds_dh, Real & ds_dp) const override;
221 :
222 : virtual Real vaporPressure(Real temperature) const override;
223 :
224 : virtual void vaporPressure(Real temperature, Real & psat, Real & dpsat_dT) const override;
225 :
226 : template <typename T>
227 : void vaporPressureTemplate(const T & temperature, T & psat, T & dpsat_dT) const;
228 :
229 : /**
230 : * Saturation temperature as a function of pressure
231 : *
232 : * Eq. (31) from Revised Release on the IAPWS Industrial
233 : * Formulation 1997 for the Thermodynamic Properties of Water
234 : * and Steam
235 : *
236 : * Valid for 611.213 Pa <= p <= 22.064 MPa
237 : *
238 : * @param pressure water pressure (Pa)
239 : * @return saturation temperature (K)
240 : */
241 : Real vaporTemperature(Real pressure) const override;
242 : virtual void vaporTemperature(Real pressure, Real & Tsat, Real & dTsat_dp) const override;
243 :
244 : /**
245 : * Auxillary equation for the boundary between regions 2 and 3
246 : *
247 : * Eq. (5) from Revised Release on the IAPWS Industrial
248 : * Formulation 1997 for the Thermodynamic Properties of Water
249 : * and Steam
250 : *
251 : * @param temperature water temperature (K)
252 : * @return pressure (Pa) on the boundary between region 2 and 3
253 : */
254 : Real b23p(Real temperature) const;
255 :
256 : /**
257 : * Auxillary equation for the boundary between regions 2 and 3
258 : *
259 : * Eq. (6) from Revised Release on the IAPWS Industrial
260 : * Formulation 1997 for the Thermodynamic Properties of Water
261 : * and Steam
262 : *
263 : * @param pressure water pressure (Pa)
264 : * @return temperature (K) on the boundary between regions 2 and 3
265 : */
266 : Real b23T(Real pressure) const;
267 :
268 : /**
269 : * Determines the phase region that the given pressure and temperature values
270 : * lie in
271 : *
272 : * @param pressure water pressure (Pa)
273 : * @param temperature water temperature (K)
274 : * @return region phase region index
275 : */
276 : unsigned int inRegion(Real pressure, Real temperature) const;
277 :
278 : /**
279 : * Provides the correct subregion index for a (P,T) point in
280 : * region 3. From Revised Supplementary Release on Backward Equations for
281 : * Specific Volume as a Function of Pressure and Temperature v(p,T)
282 : * for Region 3 of the IAPWS Industrial Formulation 1997 for the
283 : * Thermodynamic Properties of Water and Steam
284 : *
285 : * @param pressure water pressure (Pa)
286 : * @param temperature water temperature (K)
287 : * @return subregion index
288 : */
289 : unsigned int subregion3(Real pressure, Real temperature) const;
290 :
291 : /**
292 : * Specific volume in all subregions of region 3 EXCEPT subregion n (13).
293 : *
294 : * @param pi scaled water pressure
295 : * @param theta scaled water temperature
296 : * @param a to e constants
297 : * @param sid subregion ID of the subregion
298 : * @return volume water specific volume (m^3/kg)
299 : */
300 : template <typename T>
301 : T subregionVolume(const T & pi,
302 : const T & theta,
303 : Real a,
304 : Real b,
305 : Real c,
306 : Real d,
307 : Real e,
308 : unsigned int sid) const;
309 :
310 : /**
311 : * Density function for Region 3 - supercritical water and steam
312 : *
313 : * To avoid iteration, use the backwards equations for region 3
314 : * from Revised Supplementary Release on Backward Equations for
315 : * Specific Volume as a Function of Pressure and Temperature v(p,T)
316 : * for Region 3 of the IAPWS Industrial Formulation 1997 for the
317 : * Thermodynamic Properties of Water and Steam.
318 : *
319 : * @param pressure water pressure (Pa)
320 : * @param temperature water temperature (K)
321 : * @return density (kg/m^3) in region 3
322 : */
323 : template <typename T>
324 : T densityRegion3(const T & pressure, const T & temperature) const;
325 :
326 : /**
327 : * Backwards equation T(p, h)
328 : * From Revised Release on the IAPWS Industrial Formulation 1997 for the
329 : * Thermodynamic Properties of Water and Steam
330 : *
331 : * @param pressure water pressure (Pa)
332 : * @param enthalpy water enthalpy (J/kg)
333 : * @return temperature water temperature (K)
334 : */
335 : virtual Real T_from_p_h(Real pressure, Real enthalpy) const override;
336 : virtual void T_from_p_h(Real p, Real h, Real & T, Real & dT_dp, Real & dT_dh) const override;
337 : virtual ADReal T_from_p_h(const ADReal & pressure, const ADReal & enthalpy) const override;
338 :
339 : /**
340 : * Boundary between subregions b and c in region 2.
341 : * Equation 21 from Revised Release on the IAPWS Industrial
342 : * Formulation 1997 for the Thermodynamic Properties of Water
343 : * and Steam.
344 : *
345 : * @param pressure water pressure (Pa)
346 : * @return enthalpy at boundary (J/kg)
347 : */
348 : Real b2bc(Real pressure) const;
349 :
350 : /**
351 : * Boundary between subregions a and b in region 3.
352 : * Equation 1 from Revised Supplementary Release on Backward Equations for
353 : * the Functions T(p,h), v(p,h) and T(p,s), v(p,s) for Region 3 of the IAPWS
354 : * Industrial Formulation 1997 for the Thermodynamic Properties of Water and
355 : * Steam
356 : *
357 : * @param pressure water pressure (Pa)
358 : * @return enthalpy at boundary (J/kg)
359 : */
360 : Real b3ab(Real pressure) const;
361 :
362 : /**
363 : * IAPWS formulation of Henry's law constant for dissolution in water
364 : * From Guidelines on the Henry's constant and vapour
365 : * liquid distribution constant for gases in H20 and D20 at high
366 : * temperatures, IAPWS (2004)
367 : * @param T fluid temperature (K)
368 : * @param coeffs Henry's constant coefficients of gas
369 : * @param[out] Kh Henry's constant
370 : * @param[out] dKh_dT derivative of Kh wrt temperature
371 : */
372 : Real henryConstant(Real temperature, const std::vector<Real> & coeffs) const;
373 : void
374 : henryConstant(Real temperature, const std::vector<Real> & coeffs, Real & Kh, Real & dKh_dT) const;
375 : ADReal henryConstant(const ADReal & temperature, const std::vector<Real> & coeffs) const;
376 :
377 : /**
378 : * Computes the pressure (first member of the pair) and temperature (second member of the pair) as
379 : * functions of specific volume and specific internal energy
380 : */
381 : template <typename T>
382 : std::pair<T, T> p_T_from_v_e(const T & v, const T & e) const;
383 :
384 : /**
385 : * Computes the density (first member of the pair) and temperature (second member of the pair) as
386 : * functions of specific volume and specific internal energy
387 : */
388 : template <typename T>
389 : std::pair<T, T> rho_T_from_v_e(const T & v, const T & e) const;
390 :
391 : /**
392 : * Computes the pressure (first member of the pair) and temperature (second member of the pair) as
393 : * functions of specific volume and specific enthalpy
394 : */
395 : template <typename T>
396 : std::pair<T, T> p_T_from_v_h(const T & v, const T & h) const;
397 : using SinglePhaseFluidProperties::p_T_from_v_h;
398 :
399 : protected:
400 : /**
401 : * Gibbs free energy in Region 1 - single phase liquid region
402 : *
403 : * From Eq. (7) From Revised Release on the IAPWS Industrial
404 : * Formulation 1997 for the Thermodynamic Properties of Water
405 : * and Steam, IAPWS 2007.
406 : *
407 : * @param pi reduced pressure (-)
408 : * @param tau reduced temperature (-)
409 : * @return Gibbs free energy (-)
410 : */
411 : template <typename T>
412 : T gamma1(const T & pi, const T & tau) const;
413 :
414 : /**
415 : * Derivative of Gibbs free energy in Region 1 wrt pi
416 : *
417 : * @param pi reduced pressure (-)
418 : * @param tau reduced temperature (-)
419 : * @return derivative of Gibbs free energy wrt pi (-)
420 : */
421 : template <typename T>
422 : T dgamma1_dpi(const T & pi, const T & tau) const;
423 :
424 : /**
425 : * Second derivative of Gibbs free energy in Region 1 wrt pi
426 : *
427 : * @param pi reduced pressure (-)
428 : * @param tau reduced temperature (-)
429 : * @return second derivative of Gibbs free energy wrt pi (-)
430 : */
431 : template <typename T>
432 : T d2gamma1_dpi2(const T & pi, const T & tau) const;
433 :
434 : /**
435 : * Derivative of Gibbs free energy in Region 1 wrt tau
436 : *
437 : * @param pi reduced pressure (-)
438 : * @param tau reduced temperature (-)
439 : * @return derivative of Gibbs free energy wrt tau (-)
440 : */
441 : template <typename T>
442 : T dgamma1_dtau(const T & pi, const T & tau) const;
443 :
444 : /**
445 : * Second derivative of Gibbs free energy in Region 1 wrt tau
446 : *
447 : * @param pi reduced pressure (-)
448 : * @param tau reduced temperature (-)
449 : * @return second derivative of Gibbs free energy wrt tau (-)
450 : */
451 : template <typename T>
452 : T d2gamma1_dtau2(const T & pi, const T & tau) const;
453 :
454 : /**
455 : * Second derivative of Gibbs free energy in Region 1 wrt pi and tau
456 : *
457 : * @param pi reduced pressure (-)
458 : * @param tau reduced temperature (-)
459 : * @return second derivative of Gibbs free energy wrt pi and tau (-)
460 : */
461 : template <typename T>
462 : T d2gamma1_dpitau(const T & pi, const T & tau) const;
463 :
464 : /**
465 : * Gibbs free energy in Region 2 - superheated steam
466 : *
467 : * From Eq. (15) From Revised Release on the IAPWS Industrial
468 : * Formulation 1997 for the Thermodynamic Properties of Water
469 : * and Steam, IAPWS 2007.
470 : *
471 : * @param pi reduced pressure (-)
472 : * @param tau reduced temperature (-)
473 : * @return Gibbs free energy (-)
474 : */
475 : template <typename T>
476 : T gamma2(const T & pi, const T & tau) const;
477 :
478 : /**
479 : * Derivative of Gibbs free energy in Region 2 wrt pi
480 : *
481 : * @param pi reduced pressure (-)
482 : * @param tau reduced temperature (-)
483 : * @return derivative of Gibbs free energy wrt pi (-)
484 : */
485 : template <typename T>
486 : T dgamma2_dpi(const T & pi, const T & tau) const;
487 :
488 : /**
489 : * Second derivative of Gibbs free energy in Region 2 wrt pi
490 : *
491 : * @param pi reduced pressure (-)
492 : * @param tau reduced temperature (-)
493 : * @return second derivative of Gibbs free energy wrt pi (-)
494 : */
495 : template <typename T>
496 : T d2gamma2_dpi2(const T & pi, const T & tau) const;
497 :
498 : /**
499 : * Derivative of Gibbs free energy in Region 2 wrt tau
500 : *
501 : * @param pi reduced pressure (-)
502 : * @param tau reduced temperature (-)
503 : * @return derivative of Gibbs free energy wrt tau (-)
504 : */
505 : template <typename T>
506 : T dgamma2_dtau(const T & pi, const T & tau) const;
507 :
508 : /**
509 : * Second derivative of Gibbs free energy in Region 2 wrt tau
510 : *
511 : * @param pi reduced pressure (-)
512 : * @param tau reduced temperature (-)
513 : * @return second derivative of Gibbs free energy wrt tau (-)
514 : */
515 : template <typename T>
516 : T d2gamma2_dtau2(const T & pi, const T & tau) const;
517 :
518 : /**
519 : * Second derivative of Gibbs free energy in Region 2 wrt pi and tau
520 : *
521 : * @param pi reduced pressure (-)
522 : * @param tau reduced temperature (-)
523 : * @return second derivative of Gibbs free energy wrt pi and tau (-)
524 : */
525 : template <typename T>
526 : T d2gamma2_dpitau(const T & pi, const T & tau) const;
527 :
528 : /**
529 : * Helmholtz free energy in Region 3
530 : *
531 : * From Eq. (28) From Revised Release on the IAPWS Industrial
532 : * Formulation 1997 for the Thermodynamic Properties of Water
533 : * and Steam, IAPWS 2007.
534 : *
535 : * @param delta reduced density (-)
536 : * @param tau reduced temperature (-)
537 : * @return Helmholtz free energy (-)
538 : */
539 : template <typename T>
540 : T phi3(const T & delta, const T & tau) const;
541 :
542 : /**
543 : * Derivative of Helmholtz free energy in Region 3 wrt delta
544 : *
545 : * @param delta reduced density (-)
546 : * @param tau reduced temperature (-)
547 : * @return derivative of Helmholtz free energy wrt delta (-)
548 : */
549 : template <typename T>
550 : T dphi3_ddelta(const T & delta, const T & tau) const;
551 :
552 : /**
553 : * Second derivative of Helmholtz free energy in Region 3 wrt delta
554 : *
555 : * @param delta reduced density (-)
556 : * @param tau reduced temperature (-)
557 : * @return second derivative of Helmholtz free energy wrt delta (-)
558 : */
559 : template <typename T>
560 : T d2phi3_ddelta2(const T & delta, const T & tau) const;
561 :
562 : /**
563 : * Derivative of Helmholtz free energy in Region 3 wrt tau
564 : *
565 : * @param delta reduced density (-)
566 : * @param tau reduced temperature (-)
567 : * @return derivative of Helmholtz free energy wrt tau (-)
568 : */
569 : template <typename T>
570 : T dphi3_dtau(const T & delta, const T & tau) const;
571 :
572 : /**
573 : * Second derivative of Helmholtz free energy in Region 3 wrt tau
574 : *
575 : * @param delta reduced density (-)
576 : * @param tau reduced temperature (-)
577 : * @return second derivative of Helmholtz free energy wrt tau (-)
578 : */
579 : template <typename T>
580 : T d2phi3_dtau2(const T & delta, const T & tau) const;
581 :
582 : /**
583 : * Second derivative of Helmholtz free energy in Region 3 wrt delta and tau
584 : *
585 : * @param delta reduced density (-)
586 : * @param tau reduced temperature (-)
587 : * @return second derivative of Helmholtz free energy wrt delta and tau (-)
588 : */
589 : template <typename T>
590 : T d2phi3_ddeltatau(const T & delta, const T & tau) const;
591 :
592 : /**
593 : * Gibbs free energy in Region 5
594 : *
595 : * From Eq. (32) From Revised Release on the IAPWS Industrial
596 : * Formulation 1997 for the Thermodynamic Properties of Water
597 : * and Steam, IAPWS 2007.
598 : *
599 : * @param pi reduced pressure (-)
600 : * @param tau reduced temperature (-)
601 : * @return Gibbs free energy (-)
602 : */
603 : template <typename T>
604 : T gamma5(const T & pi, const T & tau) const;
605 :
606 : /**
607 : * Derivative of Gibbs free energy in Region 5 wrt pi
608 : *
609 : * @param pi reduced pressure (-)
610 : * @param tau reduced temperature (-)
611 : * @return derivative of Gibbs free energy wrt pi (-)
612 : */
613 : template <typename T>
614 : T dgamma5_dpi(const T & pi, const T & tau) const;
615 :
616 : /**
617 : * Second derivative of Gibbs free energy in Region 5 wrt pi
618 : *
619 : * @param pi reduced pressure (-)
620 : * @param tau reduced temperature (-)
621 : * @return second derivative of Gibbs free energy wrt pi (-)
622 : */
623 : template <typename T>
624 : T d2gamma5_dpi2(const T & pi, const T & tau) const;
625 :
626 : /**
627 : * Derivative of Gibbs free energy in Region 5 wrt tau
628 : *
629 : * @param pi reduced pressure (-)
630 : * @param tau reduced temperature (-)
631 : * @return derivative of Gibbs free energy wrt tau (-)
632 : */
633 : template <typename T>
634 : T dgamma5_dtau(const T & pi, const T & tau) const;
635 :
636 : /**
637 : * Second derivative of Gibbs free energy in Region 5 wrt tau
638 : *
639 : * @param pi reduced pressure (-)
640 : * @param tau reduced temperature (-)
641 : * @return second derivative of Gibbs free energy wrt tau (-)
642 : */
643 : template <typename T>
644 : T d2gamma5_dtau2(const T & pi, const T & tau) const;
645 :
646 : /**
647 : * Second derivative of Gibbs free energy in Region 5 wrt pi and tau
648 : *
649 : * @param pi reduced pressure (-)
650 : * @param tau reduced temperature (-)
651 : * @return second derivative of Gibbs free energy wrt pi and tau (-)
652 : */
653 : template <typename T>
654 : T d2gamma5_dpitau(const T & pi, const T & tau) const;
655 :
656 : /// Enum of subregion ids for region 3
657 : enum subregionEnum
658 : {
659 : AB,
660 : CD,
661 : GH,
662 : IJ,
663 : JK,
664 : MN,
665 : OP,
666 : QU,
667 : RX,
668 : UV,
669 : WX,
670 : EF
671 : };
672 :
673 : /**
674 : * Boundaries between subregions in region 3.
675 : * From Revised Supplementary Release on Backward Equations for
676 : * Specific Volume as a Function of Pressure and Temperature v(p,T)
677 : * for Region 3 of the IAPWS Industrial Formulation 1997 for the
678 : * Thermodynamic Properties of Water and Steam
679 : *
680 : * @param pressure water pressure (Pa)
681 : * @param xy string to select the boundary between two subregions
682 : * @return temperature (K) along the boundary
683 : */
684 : template <typename T>
685 : T tempXY(const T & pressure, subregionEnum xy) const;
686 :
687 : /**
688 : * Determines the phase region that the given pressure and enthaply values
689 : * lie in.
690 : *
691 : * @param pressure water pressure (Pa)
692 : * @param enthalpy water enthalpy (J/kg)
693 : * @return region phase region index
694 : */
695 : unsigned int inRegionPH(Real pressure, Real enthalpy) const;
696 :
697 : /**
698 : * Provides the correct subregion index for a (P,h) point in
699 : * region 2. From Revised Release on the IAPWS Industrial
700 : * Formulation 1997 for the Thermodynamic Properties of Water
701 : * and Steam
702 : *
703 : * @param pressure water pressure (Pa)
704 : * @param enthalpy water enthalpy (J/kg)
705 : * @return subregion index
706 : */
707 : unsigned int subregion2ph(Real pressure, Real enthalpy) const;
708 :
709 : /**
710 : * Provides the correct subregion index for a (P,h) point in
711 : * region 3. From Revised Supplementary Release on Backward Equations for
712 : * the Functions T(p,h), v(p,h) and T(p,s), v(p,s) for Region 3 of the IAPWS
713 : * Industrial Formulation 1997 for the Thermodynamic Properties of Water and
714 : * Steam
715 : *
716 : * @param pressure water pressure (Pa)
717 : * @param enthalpy water enthalpy (J/kg)
718 : * @return subregion index
719 : */
720 : unsigned int subregion3ph(Real pressure, Real enthalpy) const;
721 :
722 : /**
723 : * AD version of backwards equation T(p, h) (used internally)
724 : * From Revised Release on the IAPWS Industrial Formulation 1997 for the
725 : * Thermodynamic Properties of Water and Steam
726 : *
727 : * @param pressure water pressure (Pa)
728 : * @param enthalpy water enthalpy (J/kg)
729 : * @return temperature water temperature (K)
730 : */
731 : ADReal T_from_p_h_ad(const ADReal & pressure, const ADReal & enthalpy) const;
732 :
733 : /**
734 : * Backwards equation T(p, h) in Region 1
735 : * Eq. (11) from Revised Release on the IAPWS Industrial
736 : * Formulation 1997 for the Thermodynamic Properties of Water
737 : * and Steam
738 : *
739 : * @param pressure water pressure (Pa)
740 : * @param enthalpy water enthalpy (J/kg)
741 : * @return temperature water temperature (K)
742 : */
743 : ADReal temperature_from_ph1(const ADReal & pressure, const ADReal & enthalpy) const;
744 :
745 : /**
746 : * Backwards equation T(p, h) in Region 2a
747 : * Eq. (22) from Revised Release on the IAPWS Industrial
748 : * Formulation 1997 for the Thermodynamic Properties of Water
749 : * and Steam
750 : *
751 : * @param pressure water pressure (Pa)
752 : * @param enthalpy water enthalpy (J/kg)
753 : * @return temperature water temperature (K)
754 : */
755 : ADReal temperature_from_ph2a(const ADReal & pressure, const ADReal & enthalpy) const;
756 :
757 : /**
758 : * Backwards equation T(p, h) in Region 2b
759 : * Eq. (23) from Revised Release on the IAPWS Industrial
760 : * Formulation 1997 for the Thermodynamic Properties of Water
761 : * and Steam
762 : *
763 : * @param pressure water pressure (Pa)
764 : * @param enthalpy water enthalpy (J/kg)
765 : * @return temperature water temperature (K)
766 : */
767 : ADReal temperature_from_ph2b(const ADReal & pressure, const ADReal & enthalpy) const;
768 :
769 : /**
770 : * Backwards equation T(p, h) in Region 2c
771 : * Eq. (24) from Revised Release on the IAPWS Industrial
772 : * Formulation 1997 for the Thermodynamic Properties of Water
773 : * and Steam
774 : *
775 : * @param pressure water pressure (Pa)
776 : * @param enthalpy water enthalpy (J/kg)
777 : * @return temperature water temperature (K)
778 : */
779 : ADReal temperature_from_ph2c(const ADReal & pressure, const ADReal & enthalpy) const;
780 :
781 : /**
782 : * Backwards equation T(p, h) in Region 3a
783 : * Eq. (2) from Revised Supplementary Release on Backward Equations for
784 : * the Functions T(p,h), v(p,h) and T(p,s), v(p,s) for Region 3 of the IAPWS
785 : * Industrial Formulation 1997 for the Thermodynamic Properties of Water and
786 : * Steam
787 : *
788 : * @param pressure water pressure (Pa)
789 : * @param enthalpy water enthalpy (J/kg)
790 : * @return temperature water temperature (K)
791 : */
792 : ADReal temperature_from_ph3a(const ADReal & pressure, const ADReal & enthalpy) const;
793 :
794 : /**
795 : * Backwards equation T(p, h) in Region 3b
796 : * Eq. (3) from Revised Supplementary Release on Backward Equations for
797 : * the Functions T(p,h), v(p,h) and T(p,s), v(p,s) for Region 3 of the IAPWS
798 : * Industrial Formulation 1997 for the Thermodynamic Properties of Water and
799 : * Steam
800 : *
801 : * @param pressure water pressure (Pa)
802 : * @param enthalpy water enthalpy (J/kg)
803 : * @return temperature water temperature (K)
804 : */
805 : ADReal temperature_from_ph3b(const ADReal & pressure, const ADReal & enthalpy) const;
806 :
807 : /**
808 : * AD version of saturation temperature as a function of pressure (used internally)
809 : *
810 : * Eq. (31) from Revised Release on the IAPWS Industrial
811 : * Formulation 1997 for the Thermodynamic Properties of Water
812 : * and Steam
813 : *
814 : * Valid for 611.213 Pa <= p <= 22.064 MPa
815 : *
816 : * @param pressure water pressure (Pa)
817 : * @return saturation temperature (K)
818 : */
819 : ADReal vaporTemperature_ad(const ADReal & pressure) const;
820 :
821 : /// Water molar mass (kg/mol)
822 : const Real _Mh2o;
823 : /// Specific gas constant for H2O (universal gas constant / molar mass of water - J/kg/K)
824 : const Real _Rw;
825 : /// Critical pressure (Pa)
826 : const Real _p_critical;
827 : /// Critical temperature (K)
828 : const Real _T_critical;
829 : /// Critical density (kg/m^3)
830 : const Real _rho_critical;
831 : /// Triple point pressure (Pa)
832 : const Real _p_triple;
833 : /// Triple point temperature (K)
834 : const Real _T_triple;
835 :
836 : /**
837 : * Reference constants used in to calculate thermophysical properties of water.
838 : * Taken from Revised Release on the IAPWS Industrial Formulation 1997 for the Thermodynamic
839 : * Properties of Water and Steam and from Revised Supplementary Release on Backward Equations for
840 : * Specific Volume as a Function of Pressure and Temperature v(p,T) for Region 3 of the IAPWS
841 : * Industrial Formulation 1997 for the Thermodynamic Properties of Water and Steam
842 : */
843 :
844 : /// Constants for region 1
845 : const std::array<Real, 34> _n1{
846 : {0.14632971213167e0, -0.84548187169114e0, -0.37563603672040e1, 0.33855169168385e1,
847 : -0.95791963387872e0, 0.15772038513228e0, -0.16616417199501e-1, 0.81214629983568e-3,
848 : 0.28319080123804e-3, -0.60706301565874e-3, -0.18990068218419e-1, -0.32529748770505e-1,
849 : -0.21841717175414e-1, -0.52838357969930e-4, -0.47184321073267e-3, -0.30001780793026e-3,
850 : 0.47661393906987e-4, -0.44141845330846e-5, -0.72694996297594e-15, -0.31679644845054e-4,
851 : -0.28270797985312e-5, -0.85205128120103e-9, -0.22425281908000e-5, -0.65171222895601e-6,
852 : -0.14341729937924e-12, -0.40516996860117e-6, -0.12734301741641e-8, -0.17424871230634e-9,
853 : -0.68762131295531e-18, 0.14478307828521e-19, 0.26335781662795e-22, -0.11947622640071e-22,
854 : 0.18228094581404e-23, -0.93537087292458e-25}};
855 :
856 : const std::array<int, 34> _I1{{0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 2, 2, 2,
857 : 2, 2, 3, 3, 3, 4, 4, 4, 5, 8, 8, 21, 23, 29, 30, 31, 32}};
858 :
859 : const std::array<int, 34> _J1{{-2, -1, 0, 1, 2, 3, 4, 5, -9, -7, -1, 0,
860 : 1, 3, -3, 0, 1, 3, 17, -4, 0, 6, -5, -2,
861 : 10, -8, -11, -6, -29, -31, -38, -39, -40, -41}};
862 :
863 : const std::array<Real, 20> _nph1{
864 : {-0.23872489924521e3, 0.40421188637945e3, 0.11349746881718e3, -0.58457616048039e1,
865 : -0.15285482413140e-3, -0.10866707695377e-5, -0.13391744872602e2, 0.43211039183559e2,
866 : -0.54010067170506e2, 0.30535892203916e2, -0.65964749423638e1, 0.93965400878363e-2,
867 : 0.11573647505340e-6, -0.25858641282073e-4, -0.40644363084799e-8, 0.66456186191635e-7,
868 : 0.80670734103027e-10, -0.93477771213947e-12, 0.58265442020601e-14, -0.15020185953503e-16}};
869 :
870 : const std::array<int, 20> _Iph1{{0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 2, 2, 3, 3, 4, 5, 6}};
871 :
872 : const std::array<int, 20> _Jph1{
873 : {0, 1, 2, 6, 22, 32, 0, 1, 2, 3, 4, 10, 32, 10, 32, 10, 32, 32, 32, 32}};
874 :
875 : /// Constants for region 2
876 : const std::array<Real, 9> _n02{{-0.96927686500217e1,
877 : 0.10086655968018e2,
878 : -0.56087911283020e-2,
879 : 0.71452738081455e-1,
880 : -0.40710498223928e0,
881 : 0.14240819171444e1,
882 : -0.43839511319450e1,
883 : -0.28408632460772e0,
884 : 0.21268463753307e-1}};
885 :
886 : const std::array<Real, 43> _n2{
887 : {-0.17731742473213e-2, -0.17834862292358e-1, -0.45996013696365e-1, -0.57581259083432e-1,
888 : -0.50325278727930e-1, -0.33032641670203e-4, -0.18948987516315e-3, -0.39392777243355e-2,
889 : -0.43797295650573e-1, -0.26674547914087e-4, 0.20481737692309e-7, 0.43870667284435e-6,
890 : -0.32277677238570e-4, -0.15033924542148e-2, -0.40668253562649e-1, -0.78847309559367e-9,
891 : 0.12790717852285e-7, 0.48225372718507e-6, 0.22922076337661e-5, -0.16714766451061e-10,
892 : -0.21171472321355e-2, -0.23895741934104e2, -0.59059564324270e-17, -0.12621808899101e-5,
893 : -0.38946842435739e-1, 0.11256211360459e-10, -0.82311340897998e1, 0.19809712802088e-7,
894 : 0.10406965210174e-18, -0.10234747095929e-12, -0.10018179379511e-8, -0.80882908646985e-10,
895 : 0.10693031879409e0, -0.33662250574171e0, 0.89185845355421e-24, 0.30629316876232e-12,
896 : -0.42002467698208e-5, -0.59056029685639e-25, 0.37826947613457e-5, -0.12768608934681e-14,
897 : 0.73087610595061e-28, 0.55414715350778e-16, -0.94369707241210e-6}};
898 :
899 : const std::array<int, 9> _J02{{0, 1, -5, -4, -3, -2, -1, 2, 3}};
900 :
901 : const std::array<int, 43> _I2{{1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3,
902 : 4, 4, 4, 5, 6, 6, 6, 7, 7, 7, 8, 8, 9, 10, 10,
903 : 10, 16, 16, 18, 20, 20, 20, 21, 22, 23, 24, 24, 24}};
904 :
905 : const std::array<int, 43> _J2{{0, 1, 2, 3, 6, 1, 2, 4, 7, 36, 0, 1, 3, 6, 35,
906 : 1, 2, 3, 7, 3, 16, 35, 0, 11, 25, 8, 36, 13, 4, 10,
907 : 14, 29, 50, 57, 20, 35, 48, 21, 53, 39, 26, 40, 58}};
908 :
909 : const std::array<Real, 36> _nph2a{
910 : {0.10898952318288e4, 0.84951654495535e3, -0.10781748091826e3, 0.33153654801263e2,
911 : -0.74232016790248e1, 0.11765048724356e2, 0.18445749355790e1, -0.41792700549624e1,
912 : 0.62478196935812e1, -0.17344563108114e2, -0.20058176862096e3, 0.27196065473796e3,
913 : -0.45511318285818e3, 0.30919688604755e4, 0.25226640357872e6, -0.61707422868339e-2,
914 : -0.31078046629583, 0.11670873077107e2, 0.12812798404046e9, -0.98554909623276e9,
915 : 0.28224546973002e10, -0.35948971410703e10, 0.17227349913197e10, -0.13551334240775e5,
916 : 0.12848734664650e8, 0.13865724283226e1, 0.23598832556514e6, -0.13105236545054e8,
917 : 0.73999835474766e4, -0.55196697030060e6, 0.37154085996233e7, 0.19127729239660e5,
918 : -0.41535164835634e6, -0.62459855192507e2}};
919 :
920 : const std::array<int, 36> _Iph2a{{0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2,
921 : 2, 2, 2, 2, 2, 2, 3, 3, 4, 4, 4, 5, 5, 5, 6, 6, 7}};
922 :
923 : const std::array<int, 36> _Jph2a{{0, 1, 2, 3, 7, 20, 0, 1, 2, 3, 7, 9,
924 : 11, 18, 44, 0, 2, 7, 36, 38, 40, 42, 44, 24,
925 : 44, 12, 32, 44, 32, 36, 42, 34, 44, 28}};
926 :
927 : const std::array<Real, 38> _nph2b{
928 : {0.14895041079516e4, 0.74307798314034e3, -0.97708318797837e2, 0.24742464705674e1,
929 : -0.63281320016026, 0.11385952129658e1, -0.47811863648625, 0.85208123431544e-2,
930 : 0.93747147377932, 0.33593118604916e1, 0.33809355601454e1, 0.16844539671904,
931 : 0.73875745236695, -0.47128737436186, 0.15020273139707, -0.21764114219750e-2,
932 : -0.21810755324761e-1, -0.10829784403677, -0.46333324635812e-1, 0.71280351959551e-4,
933 : 0.11032831789999e-3, 0.18955248387902e-3, 0.30891541160537e-2, 0.13555504554949e-2,
934 : 0.28640237477456e-6, -0.10779857357512e-4, -0.76462712454814e-4, 0.14052392818316e-4,
935 : -0.31083814331434e-4, -0.10302738212103e-5, 0.28217281635040e-6, 0.12704902271945e-5,
936 : 0.73803353468292e-7, -0.11030139238909e-7, -0.81456365207833e-13, -0.25180545682962e-10,
937 : -0.17565233969407e-17, 0.86934156344163e-14}};
938 :
939 : const std::array<int, 38> _Iph2b{{0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2,
940 : 2, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 5, 5, 5, 6, 7, 7, 9, 9}};
941 :
942 : const std::array<int, 38> _Jph2b{{0, 1, 2, 12, 18, 24, 28, 40, 0, 2, 6, 12, 18,
943 : 24, 28, 40, 2, 8, 18, 40, 1, 2, 12, 24, 2, 12,
944 : 18, 24, 28, 40, 18, 24, 40, 28, 2, 28, 1, 40}};
945 :
946 : const std::array<Real, 23> _nph2c{
947 : {-0.32368398555242e13, 0.73263350902181e13, 0.35825089945447e12, -0.58340131851590e12,
948 : -0.10783068217470e11, 0.20825544563171e11, 0.61074783564516e6, 0.85977722535580e6,
949 : -0.25745723604170e5, 0.31081088422714e5, 0.12082315865936e4, 0.48219755109255e3,
950 : 0.37966001272486e1, -0.10842984880077e2, -0.45364172676660e-1, 0.14559115658698e-12,
951 : 0.11261597407230e-11, -0.17804982240686e-10, 0.12324579690832e-6, -0.11606921130984e-5,
952 : 0.27846367088554e-4, -0.59270038474176e-3, 0.12918582991878e-2}};
953 :
954 : const std::array<int, 23> _Iph2c{
955 : {-7, -7, -6, -6, -5, -5, -2, -2, -1, -1, 0, 0, 1, 1, 2, 6, 6, 6, 6, 6, 6, 6, 6}};
956 :
957 : const std::array<int, 23> _Jph2c{
958 : {0, 4, 0, 2, 0, 2, 0, 1, 0, 2, 0, 1, 4, 8, 4, 0, 1, 4, 10, 12, 16, 20, 22}};
959 :
960 : /// Constants for the boundary between regions 2 and 3
961 : const std::array<Real, 5> _n23{{0.34805185628969e3,
962 : -0.11671859879975e1,
963 : 0.10192970039326e-2,
964 : 0.57254459862746e3,
965 : 0.13918839778870e2}};
966 :
967 : /// Constants for region 3
968 : const std::array<Real, 40> _n3{
969 : {0.10658070028513e1, -0.15732845290239e2, 0.20944396974307e2, -0.76867707878716e1,
970 : 0.26185947787954e1, -0.28080781148620e1, 0.12053369696517e1, -0.84566812812502e-2,
971 : -0.12654315477714e1, -0.11524407806681e1, 0.88521043984318e0, -0.64207765181607e0,
972 : 0.38493460186671e0, -0.85214708824206e0, 0.48972281541877e1, -0.30502617256965e1,
973 : 0.39420536879154e-1, 0.12558408424308e0, -0.27999329698710e0, 0.13899799569460e1,
974 : -0.20189915023570e1, -0.82147637173963e-2, -0.47596035734923e0, 0.43984074473500e-1,
975 : -0.44476435428739e0, 0.90572070719733e0, 0.70522450087967e0, 0.10770512626332e0,
976 : -0.32913623258954e0, -0.50871062041158e0, -0.22175400873096e-1, 0.94260751665092e-1,
977 : 0.16436278447961e0, -0.13503372241348e-1, -0.14834345352472e-1, 0.57922953628084e-3,
978 : 0.32308904703711e-2, 0.80964802996215e-4, -0.16557679795037e-3, -0.44923899061815e-4}};
979 :
980 : const std::array<int, 40> _I3{{0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 3, 3,
981 : 3, 3, 3, 4, 4, 4, 4, 5, 5, 5, 6, 6, 6, 7, 8, 9, 9, 10, 10, 11}};
982 :
983 : const std::array<int, 40> _J3{{0, 0, 1, 2, 7, 10, 12, 23, 2, 6, 15, 17, 0, 2,
984 : 6, 7, 22, 26, 0, 2, 4, 16, 26, 0, 2, 4, 26, 1,
985 : 3, 26, 0, 2, 26, 2, 26, 2, 26, 0, 1, 26}};
986 :
987 : const std::array<std::array<Real, 8>, 26> _par3{
988 : {{{0.0024, 100.0, 760.0, 0.085, 0.817, 1.0, 1.0, 1.0}},
989 : {{0.0041, 100.0, 860.0, 0.280, 0.779, 1.0, 1.0, 1.0}},
990 : {{0.0022, 40.0, 690.0, 0.259, 0.903, 1.0, 1.0, 1.0}},
991 : {{0.0029, 40.0, 690.0, 0.559, 0.939, 1.0, 1.0, 4.0}},
992 : {{0.0032, 40.0, 710.0, 0.587, 0.918, 1.0, 1.0, 1.0}},
993 : {{0.0064, 40.0, 730.0, 0.587, 0.891, 0.5, 1.0, 4.0}},
994 : {{0.0027, 25.0, 660.0, 0.872, 0.971, 1.0, 1.0, 4.0}},
995 : {{0.0032, 25.0, 660.0, 0.898, 0.983, 1.0, 1.0, 4.0}},
996 : {{0.0041, 25.0, 660.0, 0.910, 0.984, 0.5, 1.0, 4.0}},
997 : {{0.0054, 25.0, 670.0, 0.875, 0.964, 0.5, 1.0, 4.0}},
998 : {{0.0077, 25.0, 680.0, 0.802, 0.935, 1.0, 1.0, 1.0}},
999 : {{0.0026, 24.0, 650.0, 0.908, 0.989, 1.0, 1.0, 4.0}},
1000 : {{0.0028, 23.0, 650.0, 1.000, 0.997, 1.0, 0.25, 1.0}},
1001 : {{0.0031, 23.0, 650.0, 0.976, 0.997, 0.0, 0.0, 0.0}},
1002 : {{0.0034, 23.0, 650.0, 0.974, 0.996, 0.5, 1.0, 1.0}},
1003 : {{0.0041, 23.0, 650.0, 0.972, 0.997, 0.5, 1.0, 1.0}},
1004 : {{0.0022, 23.0, 650.0, 0.848, 0.983, 1.0, 1.0, 4.0}},
1005 : {{0.0054, 23.0, 650.0, 0.874, 0.982, 1.0, 1.0, 1.0}},
1006 : {{0.0022, 21.0, 640.0, 0.886, 0.990, 1.0, 1.0, 4.0}},
1007 : {{0.0088, 20.0, 650.0, 0.803, 1.020, 1.0, 1.0, 1.0}},
1008 : {{0.0026, 23.0, 650.0, 0.902, 0.988, 1.0, 1.0, 1.0}},
1009 : {{0.0031, 23.0, 650.0, 0.960, 0.995, 1.0, 1.0, 1.0}},
1010 : {{0.0039, 23.0, 650.0, 0.959, 0.995, 1.0, 1.0, 4.0}},
1011 : {{0.0049, 23.0, 650.0, 0.910, 0.988, 1.0, 1.0, 1.0}},
1012 : {{0.0031, 22.0, 650.0, 0.996, 0.994, 1.0, 1.0, 4.0}},
1013 : {{0.0038, 22.0, 650.0, 0.993, 0.994, 1.0, 1.0, 4.0}}}};
1014 :
1015 : const std::array<unsigned int, 26> _par3N{{30, 32, 35, 38, 29, 42, 38, 29, 42, 29, 34, 43, 40,
1016 : 39, 24, 27, 24, 27, 29, 33, 38, 39, 35, 36, 20, 23}};
1017 :
1018 : /// Constants for all 26 subregions in region 3
1019 : const std::vector<std::vector<Real>> _n3s{
1020 : {0.110879558823853e-2, 0.572616740810616e3, -0.767051948380852e5, -0.253321069529674e-1,
1021 : 0.628008049345689e4, 0.234105654131876e6, 0.216867826045856, -0.156237904341963e3,
1022 : -0.269893956176613e5, -0.180407100085505e-3, 0.116732227668261e-2, 0.266987040856040e2,
1023 : 0.282776617243286e5, -0.242431520029523e4, 0.435217323022733e-3, -0.122494831387441e-1,
1024 : 0.179357604019989e1, 0.442729521058314e2, -0.593223489018342e-2, 0.453186261685774,
1025 : 0.135825703129140e1, 0.408748415856745e-1, 0.474686397863312, 0.118646814997915e1,
1026 : 0.546987265727549, 0.195266770452643, -0.502268790869663e-1, -0.369645308193377,
1027 : 0.633828037528420e-2, 0.797441793901017e-1},
1028 : {-0.827670470003621e-1, 0.416887126010565e2, 0.483651982197059e-1, -0.291032084950276e5,
1029 : -0.111422582236948e3, -0.202300083904014e-1, 0.294002509338515e3, 0.140244997609658e3,
1030 : -0.344384158811459e3, 0.361182452612149e3, -0.140699677420738e4, -0.202023902676481e-2,
1031 : 0.171346792457471e3, -0.425597804058632e1, 0.691346085000334e-5, 0.151140509678925e-2,
1032 : -0.416375290166236e-1, -0.413754957011042e2, -0.506673295721637e2, -0.572212965569023e-3,
1033 : 0.608817368401785e1, 0.239600660256161e2, 0.122261479925384e-1, 0.216356057692938e1,
1034 : 0.398198903368642, -0.116892827834085, -0.102845919373532, -0.492676637589284,
1035 : 0.655540456406790e-1, -0.240462535078530, -0.269798180310075e-1, 0.128369435967012},
1036 : {0.311967788763030e1, 0.276713458847564e5, 0.322583103403269e8, -0.342416065095363e3,
1037 : -0.899732529907377e6, -0.793892049821251e8, 0.953193003217388e2, 0.229784742345072e4,
1038 : 0.175336675322499e6, 0.791214365222792e7, 0.319933345844209e-4, -0.659508863555767e2,
1039 : -0.833426563212851e6, 0.645734680583292e-1, -0.382031020570813e7, 0.406398848470079e-4,
1040 : 0.310327498492008e2, -0.892996718483724e-3, 0.234604891591616e3, 0.377515668966951e4,
1041 : 0.158646812591361e-1, 0.707906336241843, 0.126016225146570e2, 0.736143655772152,
1042 : 0.676544268999101, -0.178100588189137e2, -0.156531975531713, 0.117707430048158e2,
1043 : 0.840143653860447e-1, -0.186442467471949, -0.440170203949645e2, 0.123290423502494e7,
1044 : -0.240650039730845e-1, -0.107077716660869e7, 0.438319858566475e-1},
1045 : {-0.452484847171645e-9, 0.315210389538801e-4, -0.214991352047545e-2, 0.508058874808345e3,
1046 : -0.127123036845932e8, 0.115371133120497e13, -0.197805728776273e-15, 0.241554806033972e-10,
1047 : -0.156481703640525e-5, 0.277211346836625e-2, -0.203578994462286e2, 0.144369489909053e7,
1048 : -0.411254217946539e11, 0.623449786243773e-5, -0.221774281146038e2, -0.689315087933158e5,
1049 : -0.195419525060713e8, 0.316373510564015e4, 0.224040754426988e7, -0.436701347922356e-5,
1050 : -0.404213852833996e-3, -0.348153203414663e3, -0.385294213555289e6, 0.135203700099403e-6,
1051 : 0.134648383271089e-3, 0.125031835351736e6, 0.968123678455841e-1, 0.225660517512438e3,
1052 : -0.190102435341872e-3, -0.299628410819229e-1, 0.500833915372121e-2, 0.387842482998411,
1053 : -0.138535367777182e4, 0.870745245971773, 0.171946252068742e1, -0.326650121426383e-1,
1054 : 0.498044171727877e4, 0.551478022765087e-2},
1055 : {0.715815808404721e9, -0.114328360753449e12, 0.376531002015720e-11, -0.903983668691157e-4,
1056 : 0.665695908836252e6, 0.535364174960127e10, 0.794977402335603e11, 0.922230563421437e2,
1057 : -0.142586073991215e6, -0.111796381424162e7, 0.896121629640760e4, -0.669989239070491e4,
1058 : 0.451242538486834e-2, -0.339731325977713e2, -0.120523111552278e1, 0.475992667717124e5,
1059 : -0.266627750390341e6, -0.153314954386524e-3, 0.305638404828265, 0.123654999499486e3,
1060 : -0.104390794213011e4, -0.157496516174308e-1, 0.685331118940253, 0.178373462873903e1,
1061 : -0.544674124878910, 0.204529931318843e4, -0.228342359328752e5, 0.413197481515899,
1062 : -0.341931835910405e2},
1063 : {-0.251756547792325e-7, 0.601307193668763e-5, -0.100615977450049e-2,
1064 : 0.999969140252192, 0.214107759236486e1, -0.165175571959086e2,
1065 : -0.141987303638727e-2, 0.269251915156554e1, 0.349741815858722e2,
1066 : -0.300208695771783e2, -0.131546288252539e1, -0.839091277286169e1,
1067 : 0.181545608337015e-9, -0.591099206478909e-3, 0.152115067087106e1,
1068 : 0.252956470663225e-4, 0.100726265203786e-14, -0.14977453386065e1,
1069 : -0.793940970562969e-9, -0.150290891264717e-3, 0.151205531275133e1,
1070 : 0.470942606221652e-5, 0.195049710391712e-12, -0.911627886266077e-8,
1071 : 0.604374640201265e-3, -0.225132933900136e-15, 0.610916973582981e-11,
1072 : -0.303063908043404e-6, -0.137796070798409e-4, -0.919296736666106e-3,
1073 : 0.639288223132545e-9, 0.753259479898699e-6, -0.400321478682929e-12,
1074 : 0.756140294351614e-8, -0.912082054034891e-11, -0.237612381140539e-7,
1075 : 0.269586010591874e-4, -0.732828135157839e-10, 0.241995578306660e-9,
1076 : -0.405735532730322e-3, 0.189424143498011e-9, -0.486632965074563e-9},
1077 : {0.412209020652996e-4, -0.114987238280587e7, 0.948180885032080e10, -0.195788865718971e18,
1078 : 0.4962507048713e25, -0.105549884548496e29, -0.758642165988278e12, -0.922172769596101e23,
1079 : 0.725379072059348e30, -0.617718249205859e2, 0.107555033344858e5, -0.379545802336487e8,
1080 : 0.228646846221831e12, -0.499741093010619e7, -0.280214310054101e31, 0.104915406769586e7,
1081 : 0.613754229168619e28, 0.802056715528378e32, -0.298617819828065e8, -0.910782540134681e2,
1082 : 0.135033227281565e6, -0.712949383408211e19, -0.104578785289542e37, 0.304331584444093e2,
1083 : 0.593250797959445e10, -0.364174062110798e28, 0.921791403532461, -0.337693609657471,
1084 : -0.724644143758508e2, -0.110480239272601, 0.536516031875059e1, -0.291441872156205e4,
1085 : 0.616338176535305e40, -0.120889175861180e39, 0.818396024524612e23, 0.940781944835829e9,
1086 : -0.367279669545448e5, -0.837513931798655e16},
1087 : {0.561379678887577e-1, 0.774135421587083e10, 0.111482975877938e-8, -0.143987128208183e-2,
1088 : 0.193696558764920e4, -0.605971823585005e9, 0.171951568124337e14, -0.185461154985145e17,
1089 : 0.38785116807801e-16, -0.395464327846105e-13, -0.170875935679023e3, -0.21201062070122e4,
1090 : 0.177683337348191e8, 0.110177443629575e2, -0.234396091693313e6, -0.656174421999594e7,
1091 : 0.156362212977396e-4, -0.212946257021400e1, 0.135249306374858e2, 0.177189164145813,
1092 : 0.139499167345464e4, -0.703670932036388e-2, -0.152011044389648, 0.981916922991113e-4,
1093 : 0.147199658618076e-2, 0.202618487025578e2, 0.899345518944240, -0.211346402240858,
1094 : 0.249971752957491e2},
1095 : {0.106905684359136e1, -0.148620857922333e1, 0.259862256980408e15,
1096 : -0.446352055678749e-11, -0.566620757170032e-6, -0.235302885736849e-2,
1097 : -0.269226321968839, 0.922024992944392e1, 0.357633505503772e-11,
1098 : -0.173942565562222e2, 0.700681785556229e-5, -0.267050351075768e-3,
1099 : -0.231779669675624e1, -0.753533046979752e-12, 0.481337131452891e1,
1100 : -0.223286270422356e22, -0.118746004987383e-4, 0.646412934136496e-2,
1101 : -0.410588536330937e-9, 0.422739537057241e20, 0.313698180473812e-12,
1102 : 0.16439533434504e-23, -0.339823323754373e-5, -0.135268639905021e-1,
1103 : -0.723252514211625e-14, 0.184386437538366e-8, -0.463959533752385e-1,
1104 : -0.99226310037675e14, 0.688169154439335e-16, -0.222620998452197e-10,
1105 : -0.540843018624083e-7, 0.345570606200257e-2, 0.422275800304086e11,
1106 : -0.126974478770487e-14, 0.927237985153679e-9, 0.612670812016489e-13,
1107 : -0.722693924063497e-11, -0.383669502636822e-3, 0.374684572410204e-3,
1108 : -0.931976897511086e5, -0.247690616026922e-1, 0.658110546759474e2},
1109 : {-0.111371317395540e-3, 0.100342892423685e1, 0.530615581928979e1,
1110 : 0.179058760078792e-5, -0.728541958464774e-3, -0.187576133371704e2,
1111 : 0.199060874071849e-2, 0.243574755377290e2, -0.177040785499444e-3,
1112 : -0.25968038522713e-2, -0.198704578406823e3, 0.738627790224287e-4,
1113 : -0.236264692844138e-2, -0.161023121314333e1, 0.622322971786473e4,
1114 : -0.960754116701669e-8, -0.510572269720488e-10, 0.767373781404211e-2,
1115 : 0.663855469485254e-14, -0.717590735526745e-9, 0.146564542926508e-4,
1116 : 0.309029474277013e-11, -0.464216300971708e-15, -0.390499637961161e-13,
1117 : -0.236716126781431e-9, 0.454652854268717e-11, -0.422271787482497e-2,
1118 : 0.283911742354706e-10, 0.270929002720228e1},
1119 : {-0.401215699576099e9, 0.484501478318406e11, 0.394721471363678e-14,
1120 : 0.372629967374147e5, -0.369794374168666e-29, -0.380436407012452e-14,
1121 : 0.475361629970233e-6, -0.879148916140706e-3, 0.844317863844331,
1122 : 0.122433162656600e2, -0.104529634830279e3, 0.589702771277429e3,
1123 : -0.291026851164444e14, 0.170343072841850e-5, -0.277617606975748e-3,
1124 : -0.344709605486686e1, 0.221333862447095e2, -0.194646110037079e3,
1125 : 0.808354639772825e-15, -0.18084520914547e-10, -0.696664158132412e-5,
1126 : -0.181057560300994e-2, 0.255830298579027e1, 0.328913873658481e4,
1127 : -0.173270241249904e-18, -0.661876792558034e-6, -0.39568892342125e-2,
1128 : 0.604203299819132e-17, -0.400879935920517e-13, 0.160751107464958e-8,
1129 : 0.383719409025556e-4, -0.649565446702457e-14, -0.149095328506e-11,
1130 : 0.541449377329581e-8},
1131 : {0.260702058647537e10, -0.188277213604704e15, 0.554923870289667e19, -0.758966946387758e23,
1132 : 0.413865186848908e27, -0.81503800073806e12, -0.381458260489955e33, -0.123239564600519e-1,
1133 : 0.226095631437174e8, -0.49501780950672e12, 0.529482996422863e16, -0.444359478746295e23,
1134 : 0.521635864527315e35, -0.487095672740742e55, -0.714430209937547e6, 0.127868634615495,
1135 : -0.100752127917598e2, 0.777451437960990e7, -0.108105480796471e25, -0.357578581169659e-5,
1136 : -0.212857169423484e1, 0.270706111085238e30, -0.695953622348829e33, 0.110609027472280,
1137 : 0.721559163361354e2, -0.306367307532219e15, 0.265839618885530e-4, 0.253392392889754e-1,
1138 : -0.214443041836579e3, 0.937846601489667, 0.223184043101700e1, 0.338401222509191e2,
1139 : 0.494237237179718e21, -0.198068404154428, -0.141415349881140e31, -0.993862421613651e2,
1140 : 0.125070534142731e3, -0.996473529004439e3, 0.473137909872765e5, 0.116662121219322e33,
1141 : -0.315874976271533e16, -0.445703369196945e33, 0.642794932373694e33},
1142 : {0.811384363481847, -0.568199310990094e4, -0.178657198172556e11, 0.795537657613427e32,
1143 : -0.814568209346872e5, -0.659774567602874e8, -0.152861148659302e11, -0.560165667510446e12,
1144 : 0.458384828593949e6, -0.385754000383848e14, 0.453735800004273e8, 0.939454935735563e12,
1145 : 0.266572856432938e28, -0.547578313899097e10, 0.200725701112386e15, 0.185007245563239e13,
1146 : 0.185135446828337e9, -0.170451090076385e12, 0.157890366037614e15, -0.202530509748774e16,
1147 : 0.36819392618357e60, 0.170215539458936e18, 0.639234909918741e42, -0.821698160721956e15,
1148 : -0.795260241872306e24, 0.23341586947851e18, -0.600079934586803e23, 0.594584382273384e25,
1149 : 0.189461279349492e40, -0.810093428842645e46, 0.188813911076809e22, 0.111052244098768e36,
1150 : 0.291133958602503e46, -0.329421923951460e22, -0.137570282536696e26, 0.181508996303902e28,
1151 : -0.346865122768353e30, -0.21196114877426e38, -0.128617899887675e49, 0.479817895699239e65},
1152 : {0.280967799943151e-38, 0.614869006573609e-30, 0.582238667048942e-27,
1153 : 0.390628369238462e-22, 0.821445758255119e-20, 0.402137961842776e-14,
1154 : 0.651718171878301e-12, -0.211773355803058e-7, 0.264953354380072e-2,
1155 : -0.135031446451331e-31, -0.607246643970893e-23, -0.402352115234494e-18,
1156 : -0.744938506925544e-16, 0.189917206526237e-12, 0.364975183508473e-5,
1157 : 0.177274872361946e-25, -0.334952758812999e-18, -0.421537726098389e-8,
1158 : -0.391048167929649e-1, 0.541276911564176e-13, 0.705412100773699e-11,
1159 : 0.258585887897486e-8, -0.493111362030162e-10, -0.158649699894543e-5,
1160 : -0.525037427886100, 0.220019901729615e-2, -0.643064132636925e-2,
1161 : 0.629154149015048e2, 0.135147318617061e3, 0.240560808321713e-6,
1162 : -0.890763306701305e-3, -0.440209599407714e4, -0.302807107747776e3,
1163 : 0.159158748314599e4, 0.232534272709876e6, -0.792681207132600e6,
1164 : -0.869871364662769e11, 0.354542769185671e12, 0.400849240129329e15},
1165 : {0.128746023979718e-34, -0.735234770382342e-11, 0.28907869214915e-2,
1166 : 0.244482731907223, 0.141733492030985e-23, -0.354533853059476e-28,
1167 : -0.594539202901431e-17, -.585188401782779e-8, .201377325411803e-5,
1168 : 0.138647388209306e1, -0.173959365084772e-4, 0.137680878349369e-2,
1169 : 0.814897605805513e-14, 0.425596631351839e-25, -0.387449113787755e-17,
1170 : 0.13981474793024e-12, -0.171849638951521e-2, 0.641890529513296e-21,
1171 : 0.118960578072018e-10, -0.155282762571611e-17, 0.233907907347507e-7,
1172 : -0.174093247766213e-12, 0.377682649089149e-8, -0.516720236575302e-10},
1173 : {-0.982825342010366e-4, 0.105145700850612e1, 0.116033094095084e3, 0.324664750281543e4,
1174 : -0.123592348610137e4, -0.561403450013495e-1, 0.856677401640869e-7, 0.236313425393924e3,
1175 : 0.972503292350109e-2, -0.103001994531927e1, -0.149653706199162e-8, -0.215743778861592e-4,
1176 : -0.834452198291445e1, 0.586602660564988, 0.343480022104968e-25, 0.816256095947021e-5,
1177 : 0.294985697916798e-2, 0.711730466276584e-16, 0.400954763806941e-9, 0.107766027032853e2,
1178 : -0.409449599138182e-6, -0.729121307758902e-5, 0.677107970938909e-8, 0.602745973022975e-7,
1179 : -0.382323011855257e-10, 0.179946628317437e-2, -0.345042834640005e-3},
1180 : {-0.820433843259950e5, 0.473271518461586e11, -0.805950021005413e-1, 0.328600025435980e2,
1181 : -0.35661702998249e4, -0.172985781433335e10, 0.351769232729192e8, -0.775489259985144e6,
1182 : 0.710346691966018e-4, 0.993499883820274e5, -0.642094171904570, -0.612842816820083e4,
1183 : 0.232808472983776e3, -0.142808220416837e-4, -0.643596060678456e-2, -0.428577227475614e1,
1184 : 0.225689939161918e4, 0.100355651721510e-2, 0.333491455143516, 0.109697576888873e1,
1185 : 0.961917379376452, -0.838165632204598e-1, 0.247795908411492e1, -0.319114969006533e4},
1186 : {0.144165955660863e-2, -0.701438599628258e13, -0.830946716459219e-16, 0.261975135368109,
1187 : 0.393097214706245e3, -0.104334030654021e5, 0.490112654154211e9, -0.147104222772069e-3,
1188 : 0.103602748043408e1, 0.305308890065089e1, -0.399745276971264e7, 0.569233719593750e-11,
1189 : -0.464923504407778e-1, -0.535400396512906e-17, 0.399988795693162e-12, -0.536479560201811e-6,
1190 : 0.159536722411202e-1, 0.270303248860217e-14, 0.244247453858506e-7, -0.983430636716454e-5,
1191 : 0.663513144224454e-1, -0.993456957845006e1, 0.546491323528491e3, -0.143365406393758e5,
1192 : 0.150764974125511e6, -0.337209709340105e-9, 0.377501980025469e-8},
1193 : {-0.532466612140254e23, 0.100415480000824e32, -0.191540001821367e30, 0.105618377808847e17,
1194 : 0.202281884477061e59, 0.884585472596134e8, 0.166540181638363e23, -0.313563197669111e6,
1195 : -0.185662327545324e54, -0.624942093918942e-1, -0.50416072413259e10, 0.187514491833092e5,
1196 : 0.121399979993217e-2, 0.188317043049455e1, -0.167073503962060e4, 0.965961650599775,
1197 : 0.294885696802488e1, -0.653915627346115e5, 0.604012200163444e50, -0.198339358557937,
1198 : -0.175984090163501e58, 0.356314881403987e1, -0.575991255144384e3, 0.456213415338071e5,
1199 : -0.109174044987829e8, 0.437796099975134e34, -0.616552611135792e46, 0.193568768917797e10,
1200 : 0.950898170425042e54},
1201 : {0.155287249586268e1, 0.664235115009031e1, -0.289366236727210e4, -0.385923202309848e13,
1202 : -0.291002915783761e1, -0.829088246858083e12, 0.176814899675218e1, -0.534686695713469e9,
1203 : 0.160464608687834e18, 0.196435366560186e6, 0.156637427541729e13, -0.178154560260006e1,
1204 : -0.229746237623692e16, 0.385659001648006e8, 0.110554446790543e10, -0.677073830687349e14,
1205 : -0.327910592086523e31, -0.341552040860644e51, -0.527251339709047e21, 0.245375640937055e24,
1206 : -0.168776617209269e27, 0.358958955867578e29, -0.656475280339411e36, 0.355286045512301e39,
1207 : 0.569021454413270e58, -0.700584546433113e48, -0.705772623326374e65, 0.166861176200148e53,
1208 : -0.300475129680486e61, -0.668481295196808e51, 0.428432338620678e69, -0.444227367758304e72,
1209 : -0.281396013562745e77},
1210 : {0.122088349258355e18, 0.104216468608488e10, -.882666931564652e16, 0.259929510849499e20,
1211 : 0.222612779142211e15, -0.878473585050085e18, -0.314432577551552e22, -0.216934916996285e13,
1212 : 0.159079648196849e21, -0.339567617303423e3, 0.884387651337836e13, -0.843405926846418e21,
1213 : 0.114178193518022e2, -0.122708229235641e-3, -0.106201671767107e3, 0.903443213959313e25,
1214 : -0.693996270370852e28, 0.648916718965575e-8, 0.718957567127851e4, 0.105581745346187e-2,
1215 : -0.651903203602581e15, -0.160116813274676e25, -0.510254294237837e-8, -0.152355388953402,
1216 : 0.677143292290144e12, 0.276378438378930e15, 0.116862983141686e-1, -0.301426947980171e14,
1217 : 0.169719813884840e-7, 0.104674840020929e27, -0.10801690456014e5, -0.990623601934295e-12,
1218 : 0.536116483602738e7, 0.226145963747881e22, -0.488731565776210e-9, 0.151001548880670e-4,
1219 : -0.227700464643920e5, -0.781754507698846e28},
1220 : {-0.415652812061591e-54, 0.177441742924043e-60, -0.357078668203377e-54,
1221 : 0.359252213604114e-25, -0.259123736380269e2, 0.594619766193460e5,
1222 : -0.624184007103158e11, 0.313080299915944e17, 0.105006446192036e-8,
1223 : -0.192824336984852e-5, 0.654144373749937e6, 0.513117462865044e13,
1224 : -0.697595750347391e19, -0.103977184454767e29, 0.119563135540666e-47,
1225 : -0.436677034051655e-41, 0.926990036530639e-29, 0.587793105620748e21,
1226 : 0.280375725094731e-17, -0.192359972440634e23, 0.742705723302738e27,
1227 : -0.517429682450605e2, 0.820612048645469e7, -0.188214882341448e-8,
1228 : 0.184587261114837e-1, -0.135830407782663e-5, -0.723681885626348e17,
1229 : -0.223449194054124e27, -0.111526741826431e-34, 0.276032601145151e-28,
1230 : 0.134856491567853e15, 0.652440293345860e-9, 0.510655119774360e17,
1231 : -0.468138358908732e32, -0.760667491183279e16, -0.417247986986821e-18,
1232 : 0.312545677756104e14, -0.100375333864186e15, 0.247761392329058e27},
1233 : {-0.586219133817016e-7, -0.894460355005526e11, 0.531168037519774e-30,
1234 : 0.109892402329239, -0.575368389425212e-1, 0.228276853990249e5,
1235 : -0.158548609655002e19, 0.329865748576503e-27, -0.634987981190669e-24,
1236 : 0.615762068640611e-8, -0.961109240985747e8, -0.406274286652625e-44,
1237 : -0.471103725498077e-12, 0.725937724828145, 0.187768525763682e-38,
1238 : -0.103308436323771e4, -0.662552816342168e-1, 0.579514041765710e3,
1239 : 0.237416732616644e-26, 0.271700235739893e-14, -0.9078862134836e2,
1240 : -0.171242509570207e-36, 0.156792067854621e3, 0.923261357901470,
1241 : -0.597865988422577e1, 0.321988767636389e7, -0.399441390042203e-29,
1242 : 0.493429086046981e-7, 0.812036983370565e-19, -0.207610284654137e-11,
1243 : -0.340821291419719e-6, 0.542000573372233e-17, -0.856711586510214e-12,
1244 : 0.266170454405981e-13, 0.858133791857099e-5},
1245 : {0.377373741298151e19, -0.507100883722913e13, -0.10336322559886e16, 0.184790814320773e-5,
1246 : -0.924729378390945e-3, -0.425999562292738e24, -0.462307771873973e-12, 0.107319065855767e22,
1247 : 0.648662492280682e11, 0.244200600688281e1, -0.851535733484258e10, 0.169894481433592e22,
1248 : 0.215780222509020e-26, -0.320850551367334, -0.382642448458610e17, -0.275386077674421e-28,
1249 : -0.563199253391666e6, -0.326068646279314e21, 0.397949001553184e14, 0.100824008584757e-6,
1250 : 0.162234569738433e5, -0.432355225319745e11, -0.59287424559861e12, 0.133061647281106e1,
1251 : 0.157338197797544e7, 0.258189614270853e14, 0.262413209706358e25, -0.920011937431142e-1,
1252 : 0.220213765905426e-2, -0.110433759109547e2, 0.847004870612087e7, -0.592910695762536e9,
1253 : -0.183027173269660e-4, 0.181339603516302, -0.119228759669889e4, 0.430867658061468e7},
1254 : {-0.525597995024633e-9, 0.583441305228407e4, -0.134778968457925e17, 0.118973500934212e26,
1255 : -0.159096490904708e27, -0.315839902302021e-6, 0.496212197158239e3, 0.327777227273171e19,
1256 : -0.527114657850696e22, 0.210017506281863e-16, 0.705106224399834e21, -0.266713136106469e31,
1257 : -0.145370512554562e-7, 0.149333917053130e28, -0.149795620287641e8, -0.3818819062711e16,
1258 : 0.724660165585797e-4, -0.937808169550193e14, 0.514411468376383e10, -0.828198594040141e5},
1259 : {0.24400789229065e-10, -0.463057430331242e7, 0.728803274777712e10, 0.327776302858856e16,
1260 : -0.110598170118409e10, -0.323899915729957e13, 0.923814007023245e16, 0.842250080413712e-12,
1261 : 0.663221436245506e12, -0.167170186672139e15, 0.253749358701391e4, -0.819731559610523e-20,
1262 : 0.328380587890663e12, -0.625004791171543e8, 0.803197957462023e21, -0.204397011338353e-10,
1263 : -0.378391047055938e4, 0.97287654593862e-2, 0.154355721681459e2, -0.373962862928643e4,
1264 : -0.682859011374572e11, -0.248488015614543e-3, 0.394536049497068e7}};
1265 :
1266 : const std::vector<std::vector<int>> _I3s{
1267 : {-12, -12, -12, -10, -10, -10, -8, -8, -8, -6, -5, -5, -5, -4, -3,
1268 : -3, -3, -3, -2, -2, -2, -1, -1, -1, 0, 0, 1, 1, 2, 2},
1269 : {-12, -12, -10, -10, -8, -6, -6, -6, -5, -5, -5, -4, -4, -4, -3, -3,
1270 : -3, -3, -3, -2, -2, -2, -1, -1, 0, 0, 1, 1, 2, 3, 4, 4},
1271 : {-12, -12, -12, -10, -10, -10, -8, -8, -8, -6, -5, -5, -5, -4, -4, -3, -3, -2,
1272 : -2, -2, -1, -1, -1, 0, 0, 0, 1, 1, 2, 2, 2, 2, 3, 3, 8},
1273 : {-12, -12, -12, -12, -12, -12, -10, -10, -10, -10, -10, -10, -10, -8, -8, -8, -8, -6, -6,
1274 : -5, -5, -5, -5, -4, -4, -4, -3, -3, -2, -2, -1, -1, -1, 0, 0, 1, 1, 3},
1275 : {-12, -12, -10, -10, -10, -10, -10, -8, -8, -8, -6, -5, -4, -4, -3,
1276 : -3, -3, -2, -2, -2, -2, -1, 0, 0, 1, 1, 1, 2, 2},
1277 : {0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 3, 3, 3, 4, 5, 5, 6, 7, 7,
1278 : 10, 12, 12, 12, 14, 14, 14, 14, 14, 16, 16, 18, 18, 20, 20, 20, 22, 24, 24, 28, 32},
1279 : {-12, -12, -12, -12, -12, -12, -10, -10, -10, -8, -8, -8, -8, -6, -6, -5, -5, -4, -3,
1280 : -2, -2, -2, -2, -1, -1, -1, 0, 0, 0, 1, 1, 1, 3, 5, 6, 8, 10, 10},
1281 : {-12, -12, -10, -10, -10, -10, -10, -10, -8, -8, -8, -8, -8, -6, -6,
1282 : -6, -5, -5, -5, -4, -4, -3, -3, -2, -1, -1, 0, 1, 1},
1283 : {0, 0, 0, 1, 1, 1, 1, 2, 3, 3, 4, 4, 4, 5, 5, 5, 7, 7, 8, 8, 10,
1284 : 12, 12, 12, 14, 14, 14, 14, 18, 18, 18, 18, 18, 20, 20, 22, 24, 24, 32, 32, 36, 36},
1285 : {0, 0, 0, 1, 1, 1, 2, 2, 3, 4, 4, 5, 5, 5, 6,
1286 : 10, 12, 12, 14, 14, 14, 16, 18, 20, 20, 24, 24, 28, 28},
1287 : {-2, -2, -1, -1, 0, -0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1,
1288 : 1, 2, 2, 2, 2, 2, 2, 5, 5, 5, 6, 6, 6, 6, 8, 10, 12},
1289 : {-12, -12, -12, -12, -12, -10, -10, -8, -8, -8, -8, -8, -8, -8, -6,
1290 : -5, -5, -4, -4, -3, -3, -3, -3, -2, -2, -2, -1, -1, -1, 0,
1291 : 0, 0, 0, 1, 1, 2, 4, 5, 5, 6, 10, 10, 14},
1292 : {0, 3, 8, 20, 1, 3, 4, 5, 1, 6, 2, 4, 14, 2, 5, 3, 0, 1, 1, 1,
1293 : 28, 2, 16, 0, 5, 0, 3, 4, 12, 16, 1, 8, 14, 0, 2, 3, 4, 8, 14, 24},
1294 : {0, 3, 4, 6, 7, 10, 12, 14, 18, 0, 3, 5, 6, 8, 12, 0, 3, 7, 12, 2,
1295 : 3, 4, 2, 4, 7, 4, 3, 5, 6, 0, 0, 3, 1, 0, 1, 0, 1, 0, 1},
1296 : {0, 0, 0, 2, 3, 4, 4, 4, 4, 4, 5, 5, 6, 7, 8, 8, 8, 10, 10, 14, 14, 20, 20, 24},
1297 : {0, 0, 0, 0, 1, 2, 3, 3, 4, 6, 7, 7, 8, 10,
1298 : 12, 12, 12, 14, 14, 14, 16, 18, 20, 22, 24, 24, 36},
1299 : {-12, -12, -10, -10, -10, -10, -8, -6, -5, -5, -4, -4,
1300 : -3, -2, -2, -2, -2, -1, -1, -1, 0, 1, 1, 1},
1301 : {-8, -8, -3, -3, -3, -3, -3, 0, 0, 0, 0, 3, 3, 8,
1302 : 8, 8, 8, 10, 10, 10, 10, 10, 10, 10, 10, 12, 14},
1303 : {-12, -12, -10, -8, -6, -5, -5, -4, -4, -3, -3, -2, -1, -1, -1,
1304 : 0, 0, 0, 0, 1, 1, 3, 3, 3, 4, 4, 4, 5, 14},
1305 : {0, 0, 0, 0, 1, 1, 2, 2, 2, 3, 3, 4, 4, 7, 7, 7, 7,
1306 : 7, 10, 10, 10, 10, 10, 18, 20, 22, 22, 24, 28, 32, 32, 32, 36},
1307 : {-12, -10, -10, -10, -8, -8, -8, -6, -6, -5, -5, -5, -3, -1, -1, -1, -1, 0, 0,
1308 : 1, 2, 2, 3, 5, 5, 5, 6, 6, 8, 8, 10, 12, 12, 12, 14, 14, 14, 14},
1309 : {-10, -8, -6, -6, -6, -6, -6, -6, -5, -5, -5, -5, -5, -5, -4, -4, -4, -4, -3, -3,
1310 : -3, -2, -2, -1, -1, 0, 0, 0, 1, 1, 3, 4, 4, 4, 5, 8, 10, 12, 14},
1311 : {-12, -12, -10, -10, -8, -8, -8, -6, -6, -6, -6, -5, -4, -4, -3, -3, -2, -2,
1312 : -1, -1, -1, 0, 0, 1, 2, 2, 3, 3, 5, 5, 5, 8, 8, 10, 10},
1313 : {-8, -6, -5, -4, -4, -4, -3, -3, -1, 0, 0, 0, 1, 1, 2, 3, 3, 3,
1314 : 4, 5, 5, 5, 6, 8, 8, 8, 8, 10, 12, 12, 12, 12, 14, 14, 14, 14},
1315 : {0, 0, 0, 0, 1, 2, 2, 2, 2, 3, 3, 3, 4, 4, 5, 5, 8, 8, 10, 12},
1316 : {-8, -6, -5, -5, -4, -4, -4, -3, -3, -3, -2, -1, 0, 1, 2, 3, 3, 6, 6, 6, 6, 8, 8}};
1317 :
1318 : const std::vector<std::vector<int>> _J3s{
1319 : {5, 10, 12, 5, 10, 12, 5, 8, 10, 1, 1, 5, 10, 8, 0,
1320 : 1, 3, 6, 0, 2, 3, 0, 1, 2, 0, 1, 0, 2, 0, 2},
1321 : {10, 12, 8, 14, 8, 5, 6, 8, 5, 8, 10, 2, 4, 5, 0, 1,
1322 : 2, 3, 5, 0, 2, 5, 0, 2, 0, 1, 0, 2, 0, 2, 0, 1},
1323 : {6, 8, 10, 6, 8, 10, 5, 6, 7, 8, 1, 4, 7, 2, 8, 0, 3, 0,
1324 : 4, 5, 0, 1, 2, 0, 1, 2, 0, 2, 0, 1, 3, 7, 0, 7, 1},
1325 : {4, 6, 7, 10, 12, 16, 0, 2, 4, 6, 8, 10, 14, 3, 7, 8, 10, 6, 8,
1326 : 1, 2, 5, 7, 0, 1, 7, 2, 4, 0, 1, 0, 1, 5, 0, 2, 0, 6, 0},
1327 : {14, 16, 3, 6, 10, 14, 16, 7, 8, 10, 6, 6, 2, 4, 2, 6, 7, 0, 1, 3, 4, 0, 0, 1, 0, 4, 6, 0, 2},
1328 : {-3, -2, -1, 0, 1, 2, -1, 1, 2, 3, 0, 1, -5, -2,
1329 : 0, -3, -8, 1, -6, -4, 1, -6, -10, -8, -4, -12, -10, -8,
1330 : -6, -4, -10, -8, -12, -10, -12, -10, -6, -12, -12, -4, -12, -12},
1331 : {7, 12, 14, 18, 22, 24, 14, 20, 24, 7, 8, 10, 12, 8, 22, 7, 20, 22, 7,
1332 : 3, 5, 14, 24, 2, 8, 18, 0, 1, 2, 0, 1, 3, 24, 22, 12, 3, 0, 6},
1333 : {8, 12, 4, 6, 8, 10, 14, 16, 0, 1, 6, 7, 8, 4, 6, 8, 2, 3, 4, 2, 4, 1, 2, 0, 0, 2, 0, 0, 2},
1334 : {0, 1, 10, -4, -2, -1, 0, 0, -5, 0, -3, -2, -1, -6, -1, 12, -4, -3, -6, 10, -8,
1335 : -12, -6, -4, -10, -8, -4, 5, -12, -10, -8, -6, 2, -12, -10, -12, -12, -8, -10, -5, -10, -8},
1336 : {-1, 0, 1, -2, -1, 1, -1, 1, -2, -2, 2, -3, -2, 0, 3,
1337 : -6, -8, -3, -10, -8, -5, -10, -12, -12, -10, -12, -6, -12, -5},
1338 : {10, 12, -5, 6, -12, -6, -2, -1, 0, 1, 2, 3, 14, -3, -2, 0, 1,
1339 : 2, -8, -6, -3, -2, 0, 4, -12, -6, -3, -12, -10, -8, -5, -12, -12, -10},
1340 : {14, 16, 18, 20, 22, 14, 24, 6, 10, 12, 14, 18, 24, 36, 8, 4, 5, 7, 16, 1, 3, 18,
1341 : 20, 2, 3, 10, 0, 1, 3, 0, 1, 2, 12, 0, 16, 1, 0, 0, 1, 14, 4, 12, 10},
1342 : {0, 0, 0, 2, 5, 5, 5, 5, 6, 6, 7, 8, 8, 10, 10, 12, 14, 14, 18, 20,
1343 : 20, 22, 22, 24, 24, 28, 28, 28, 28, 28, 32, 32, 32, 36, 36, 36, 36, 36, 36, 36},
1344 : {-12, -12, -12, -12, -12, -12, -12, -12, -12, -10, -10, -10, -10,
1345 : -10, -10, -8, -8, -8, -8, -6, -6, -6, -5, -5, -5, -4,
1346 : -3, -3, -3, -2, -1, -1, 0, 1, 1, 2, 4, 5, 6},
1347 : {-12, -4, -1, -1, -10, -12, -8, -5, -4, -1, -4, -3,
1348 : -8, -12, -10, -8, -4, -12, -8, -12, -8, -12, -10, -12},
1349 : {-1, 0, 1, 2, 1, -1, -3, 0, -2, -2, -5, -4, -2, -3,
1350 : -12, -6, -5, -10, -8, -3, -8, -8, -10, -10, -12, -8, -12},
1351 : {10, 12, 6, 7, 8, 10, 8, 6, 2, 5, 3, 4, 3, 0, 1, 2, 4, 0, 1, 2, 0, 0, 1, 3},
1352 : {6, 14, -3, 3, 4, 5, 8, -1, 0, 1, 5, -6, -2, -12,
1353 : -10, -8, -5, -12, -10, -8, -6, -5, -4, -3, -2, -12, -12},
1354 : {20, 24, 22, 14, 36, 8, 16, 6, 32, 3, 8, 4, 1, 2, 3,
1355 : 0, 1, 4, 28, 0, 32, 0, 1, 2, 3, 18, 24, 4, 24},
1356 : {0, 1, 4, 12, 0, 10, 0, 6, 14, 3, 8, 0, 10, 3, 4, 7, 20,
1357 : 36, 10, 12, 14, 16, 22, 18, 32, 22, 36, 24, 28, 22, 32, 36, 36},
1358 : {14, 10, 12, 14, 10, 12, 14, 8, 12, 4, 8, 12, 2, -1, 1, 12, 14, -3, 1,
1359 : -2, 5, 10, -5, -4, 2, 3, -5, 2, -8, 8, -4, -12, -4, 4, -12, -10, -6, 6},
1360 : {-8, -12, -12, -3, 5, 6, 8, 10, 1, 2, 6, 8, 10, 14, -12, -10, -6, 10, -3, 10,
1361 : 12, 2, 4, -2, 0, -2, 6, 10, -12, -10, 3, -6, 3, 10, 2, -12, -2, -3, 1},
1362 : {8, 14, -1, 8, 6, 8, 14, -4, -3, 2, 8, -10, -1, 3, -10, 3, 1, 2,
1363 : -8, -4, 1, -12, 1, -1, -1, 2, -12, -5, -10, -8, -6, -12, -10, -12, -8},
1364 : {14, 10, 10, 1, 2, 14, -2, 12, 5, 0, 4, 10, -10, -1, 6, -12, 0, 8,
1365 : 3, -6, -2, 1, 1, -6, -3, 1, 8, -8, -10, -8, -5, -4, -12, -10, -8, -6},
1366 : {-3, 1, 5, 8, 8, -4, -1, 4, 5, -8, 4, 8, -6, 6, -2, 1, -8, -2, -5, -8},
1367 : {3, 6, 6, 8, 5, 6, 8, -2, 5, 6, 2, -6, 3, 1, 6, -6, -2, -6, -5, -4, -1, -8, -4}};
1368 :
1369 : const std::array<Real, 31> _nph3a{
1370 : {-0.133645667811215e-6, 0.455912656802978e-5, -0.146294640700979e-4, 0.639341312970080e-2,
1371 : 0.372783927268847e3, -0.718654377460447e4, 0.573494752103400e6, -0.267569329111439e7,
1372 : -0.334066283302614e-4, -0.245479214069597e-1, 0.478087847764996e2, 0.764664131818904e-5,
1373 : 0.128350627676972e-2, 0.171219081377331e-1, -0.851007304583213e1, -0.136513461629781e-1,
1374 : -0.384460997596657e-5, 0.337423807911655e-2, -0.551624873066791, 0.729202277107470,
1375 : -0.992522757376041e-2, -.119308831407288, .793929190615421, .454270731799386,
1376 : .20999859125991, -0.642109823904738e-2, -0.235155868604540e-1, 0.252233108341612e-2,
1377 : -0.764885133368119e-2, 0.136176427574291e-1, -0.133027883575669e-1}};
1378 :
1379 : const std::array<int, 31> _Iph3a{{-12, -12, -12, -12, -12, -12, -12, -12, -10, -10, -10,
1380 : -8, -8, -8, -8, -5, -3, -2, -2, -2, -1, -1,
1381 : 0, 0, 1, 3, 3, 4, 4, 10, 12}};
1382 :
1383 : const std::array<int, 31> _Jph3a{{0, 1, 2, 6, 14, 16, 20, 22, 1, 5, 12, 0, 2, 4, 10, 2,
1384 : 0, 1, 3, 4, 0, 2, 0, 1, 1, 0, 1, 0, 3, 4, 5}};
1385 :
1386 : const std::array<Real, 33> _nph3b{
1387 : {0.323254573644920e-4, -0.127575556587181e-3, -0.475851877356068e-3, 0.156183014181602e-2,
1388 : 0.105724860113781, -0.858514221132534e2, 0.724140095480911e3, 0.296475810273257e-2,
1389 : -0.592721983365988e-2, -0.126305422818666e-1, -0.115716196364853, 0.849000969739595e2,
1390 : -0.108602260086615e-1, 0.154304475328851e-1, 0.750455441524466e-1, 0.252520973612982e-1,
1391 : -0.602507901232996e-1, -0.307622221350501e1, -0.574011959864879e-1, 0.503471360939849e1,
1392 : -0.925081888584834, 0.391733882917546e1, -0.773146007130190e2, 0.949308762098587e4,
1393 : -0.141043719679409e7, 0.849166230819026e7, 0.861095729446704, 0.323346442811720,
1394 : 0.873281936020439, -0.436653048526683, 0.286596714529479, -0.131778331276228,
1395 : 0.676682064330275e-2}};
1396 :
1397 : const std::array<int, 33> _Iph3b{{-12, -12, -10, -10, -10, -10, -10, -8, -8, -8, -8,
1398 : -8, -6, -6, -6, -4, -4, -3, -2, -2, -1, -1,
1399 : -1, -1, -1, -1, 0, 0, 1, 3, 5, 6, 8}};
1400 :
1401 : const std::array<int, 33> _Jph3b{{0, 1, 0, 1, 5, 10, 12, 0, 1, 2, 4, 10, 0, 1, 2, 0, 1,
1402 : 5, 0, 4, 2, 4, 6, 10, 14, 16, 0, 2, 1, 1, 1, 1, 1}};
1403 :
1404 : /// Constants for region 4 (the saturation curve up to the critical point)
1405 : const std::array<Real, 10> _n4{{0.11670521452767e4,
1406 : -0.72421316703206e6,
1407 : -0.17073846940092e2,
1408 : 0.12020824702470e5,
1409 : -0.32325550322333e7,
1410 : 0.14915108613530e2,
1411 : -0.48232657361591e4,
1412 : 0.40511340542057e6,
1413 : -0.238555575678490,
1414 : 0.65017534844798e3}};
1415 :
1416 : /// Constants for region 5
1417 : const std::array<int, 6> _J05{{0, 1, -3, -2, -1, 2}};
1418 : const std::array<Real, 6> _n05{{-0.13179983674201e2,
1419 : 0.68540841634434e1,
1420 : -0.24805148933466e-1,
1421 : 0.36901534980333,
1422 : -0.31161318213925e1,
1423 : -0.32961626538917}};
1424 :
1425 : const std::array<int, 6> _I5{{1, 1, 1, 2, 2, 3}};
1426 : const std::array<int, 6> _J5{{1, 2, 3, 3, 9, 7}};
1427 : const std::array<Real, 6> _n5{{0.15736404855259e-2,
1428 : 0.90153761673944e-3,
1429 : -0.50270077677648e-2,
1430 : 0.22440037409485e-5,
1431 : -0.41163275453471e-5,
1432 : 0.37919454822955e-7}};
1433 :
1434 : /// Constnats for the tempXY() method
1435 : const std::vector<std::vector<int>> _tempXY_I{{0, 1, 2, -1, -2},
1436 : {0, 1, 2, 3},
1437 : {0, 1, 2, 3, 4},
1438 : {0, 1, 2, 3, 4},
1439 : {0, 1, 2, 3, 4},
1440 : {0, 1, 2, 3},
1441 : {0, 1, 2, -1, -2},
1442 : {0, 1, 2, 3},
1443 : {0, 1, 2, 3},
1444 : {0, 1, 2, 3, 4},
1445 : {0, 1, 2, -1, -2}};
1446 :
1447 : const std::vector<std::vector<Real>> _tempXY_n{
1448 : {0.154793642129415e4,
1449 : -0.187661219490113e3,
1450 : 0.213144632222113e2,
1451 : -0.191887498864292e4,
1452 : 0.918419702359447e3},
1453 : {0.585276966696349e3, 0.278233532206915e1, -0.127283549295878e-1, 0.159090746562729e-3},
1454 : {-0.249284240900418e5,
1455 : 0.428143584791546e4,
1456 : -0.269029173140130e3,
1457 : 0.751608051114157e1,
1458 : -0.787105249910383e-1},
1459 : {0.584814781649163e3,
1460 : -0.616179320924617,
1461 : 0.260763050899562,
1462 : -0.587071076864459e-2,
1463 : 0.515308185433082e-4},
1464 : {0.617229772068439e3,
1465 : -0.770600270141675e1,
1466 : 0.697072596851896,
1467 : -0.157391839848015e-1,
1468 : 0.137897492684194e-3},
1469 : {0.535339483742384e3, 0.761978122720128e1, -0.158365725441648, 0.192871054508108e-2},
1470 : {0.969461372400213e3,
1471 : -0.332500170441278e3,
1472 : 0.642859598466067e2,
1473 : 0.773845935768222e3,
1474 : -0.152313732937084e4},
1475 : {0.565603648239126e3, 0.529062258221222e1, -0.102020639611016, 0.122240301070145e-2},
1476 : {0.584561202520006e3, -0.102961025163669e1, 0.243293362700452, -0.294905044740799e-2},
1477 : {0.528199646263062e3, 0.890579602135307e1, -0.222814134903755, 0.286791682263697e-2},
1478 : {0.728052609145380e1,
1479 : 0.973505869861952e2,
1480 : 0.147370491183191e2,
1481 : 0.329196213998375e3,
1482 : 0.873371668682417e3}};
1483 :
1484 : /// Constants from Release on the IAPWS Formulation 2008 for the Viscosity of
1485 : /// Ordinary Water Substance
1486 : const std::array<Real, 4> _mu_H0{{1.67752, 2.20462, 0.6366564, -0.241605}};
1487 : const std::array<std::array<Real, 7>, 6> _mu_Hij{
1488 : {{{5.20094e-1, 2.22531e-1, -2.81378e-1, 1.61913e-1, -3.25372e-2, 0.0, 0.0}},
1489 : {{8.50895e-2, 9.99115e-1, -9.06851e-1, 2.57399e-1, 0.0, 0.0, 0.0}},
1490 : {{-1.08374, 1.88797, -7.72479e-1, 0.0, 0.0, 0.0, 0.0}},
1491 : {{-2.89555e-1, 1.26613, -4.89837e-1, 0.0, 6.98452e-2, 0.0, -4.35673e-3}},
1492 : {{0.0, 0.0, -2.57040e-1, 0.0, 0.0, 8.72102e-3, 0.0}},
1493 : {{0.0, 1.20573e-1, 0.0, 0.0, 0.0, 0.0, -5.93264e-4}}}};
1494 :
1495 : /// Constants for thermal conductivity
1496 : std::array<Real, 4> _k_a{{0.0102811, 0.0299621, 0.0156146, -0.00422464}};
1497 :
1498 : /// Temperature scale for each region
1499 : const std::array<Real, 5> _T_star{{1386.0, 540.0, _T_critical, 1.0, 1000.0}};
1500 : /// Pressure scale for each region
1501 : const std::array<Real, 5> _p_star{{16.53e6, 1.0e6, 1.0e6, 1.0e6, 1.0e6}};
1502 :
1503 : private:
1504 : /**
1505 : * Computes the temperature (first member of the pair) and the derivative of density (second
1506 : * member of the pair) with respect to temperature as a function of pressure and density
1507 : */
1508 : template <typename T>
1509 : std::pair<T, T> T_drhodT_from_p_rho(const T & p, const T & rho) const;
1510 : };
1511 :
1512 : #pragma GCC diagnostic pop
1513 :
1514 : template <typename T>
1515 : T
1516 653074 : Water97FluidProperties::dgamma1_dpi(const T & pi, const T & tau) const
1517 : {
1518 336 : T sum = 0.0;
1519 22857590 : for (std::size_t i = 0; i < _n1.size(); ++i)
1520 22227364 : sum += -_n1[i] * _I1[i] * MathUtils::pow(7.1 - pi, _I1[i] - 1) *
1521 22190508 : MathUtils::pow(tau - 1.222, _J1[i]);
1522 :
1523 653074 : return sum;
1524 : }
1525 :
1526 : template <typename T>
1527 : T
1528 9201 : Water97FluidProperties::dgamma2_dpi(const T & pi, const T & tau) const
1529 : {
1530 : // Ideal gas part of the Gibbs free energy
1531 9201 : T dg0 = 1.0 / pi;
1532 :
1533 : // Residual part of the Gibbs free energy
1534 0 : T dgr = 0.0;
1535 404844 : for (std::size_t i = 0; i < _n2.size(); ++i)
1536 395643 : dgr += _n2[i] * _I2[i] * MathUtils::pow(pi, _I2[i] - 1) * MathUtils::pow(tau - 0.5, _J2[i]);
1537 :
1538 9201 : return dg0 + dgr;
1539 : }
1540 :
1541 : template <typename T>
1542 : T
1543 102 : Water97FluidProperties::densityRegion3(const T & pressure, const T & temperature) const
1544 : {
1545 : using std::exp;
1546 :
1547 : // Region 3 is subdivided into 26 subregions, each with a given backwards equation
1548 : // to directly calculate density from pressure and temperature without the need for
1549 : // expensive iterations. Find the subregion id that the point is in:
1550 : unsigned int sid =
1551 102 : subregion3(MetaPhysicL::raw_value(pressure), MetaPhysicL::raw_value(temperature));
1552 :
1553 0 : T vstar, pi, theta;
1554 : Real a, b, c, d, e;
1555 : unsigned int N;
1556 :
1557 102 : vstar = _par3[sid][0];
1558 105 : pi = pressure / _par3[sid][1] / 1.0e6;
1559 102 : theta = temperature / _par3[sid][2];
1560 102 : a = _par3[sid][3];
1561 102 : b = _par3[sid][4];
1562 102 : c = _par3[sid][5];
1563 102 : d = _par3[sid][6];
1564 102 : e = _par3[sid][7];
1565 102 : N = _par3N[sid];
1566 :
1567 0 : T sum = 0.0;
1568 0 : T volume = 0.0;
1569 :
1570 : // Note that subregion 13 is the only different formulation
1571 102 : if (sid == 13)
1572 : {
1573 80 : for (std::size_t i = 0; i < N; ++i)
1574 78 : sum += _n3s[sid][i] * MathUtils::pow(pi - a, _I3s[sid][i]) *
1575 78 : MathUtils::pow(theta - b, _J3s[sid][i]);
1576 :
1577 2 : volume = vstar * exp(sum);
1578 : }
1579 : else
1580 103 : volume = vstar * subregionVolume(pi, theta, a, b, c, d, e, sid);
1581 :
1582 : // Density is the inverse of volume
1583 102 : return 1.0 / volume;
1584 : }
1585 :
1586 : template <typename T>
1587 : T
1588 242 : Water97FluidProperties::gamma1(const T & pi, const T & tau) const
1589 : {
1590 5 : T sum = 0.0;
1591 8470 : for (std::size_t i = 0; i < _n1.size(); ++i)
1592 8568 : sum += _n1[i] * MathUtils::pow(7.1 - pi, _I1[i]) * MathUtils::pow(tau - 1.222, _J1[i]);
1593 :
1594 242 : return sum;
1595 : }
1596 :
1597 : template <typename T>
1598 : T
1599 457 : Water97FluidProperties::d2gamma1_dpi2(const T & pi, const T & tau) const
1600 : {
1601 11 : T sum = 0.0;
1602 15995 : for (std::size_t i = 0; i < _n1.size(); ++i)
1603 16286 : sum += _n1[i] * _I1[i] * (_I1[i] - 1) * MathUtils::pow(7.1 - pi, _I1[i] - 2) *
1604 15164 : MathUtils::pow(tau - 1.222, _J1[i]);
1605 :
1606 457 : return sum;
1607 : }
1608 :
1609 : template <typename T>
1610 : T
1611 435224 : Water97FluidProperties::dgamma1_dtau(const T & pi, const T & tau) const
1612 : {
1613 61 : T g = 0.0;
1614 15232840 : for (std::size_t i = 0; i < _n1.size(); ++i)
1615 14801764 : g += _n1[i] * _J1[i] * MathUtils::pow(7.1 - pi, _I1[i]) *
1616 14794488 : MathUtils::pow(tau - 1.222, _J1[i] - 1);
1617 :
1618 435224 : return g;
1619 : }
1620 :
1621 : template <typename T>
1622 : T
1623 1119 : Water97FluidProperties::d2gamma1_dtau2(const T & pi, const T & tau) const
1624 : {
1625 16 : T dg = 0.0;
1626 39165 : for (std::size_t i = 0; i < _n1.size(); ++i)
1627 39134 : dg += _n1[i] * _J1[i] * (_J1[i] - 1) * MathUtils::pow(7.1 - pi, _I1[i]) *
1628 37502 : MathUtils::pow(tau - 1.222, _J1[i] - 2);
1629 :
1630 1119 : return dg;
1631 : }
1632 :
1633 : template <typename T>
1634 : T
1635 457 : Water97FluidProperties::d2gamma1_dpitau(const T & pi, const T & tau) const
1636 : {
1637 11 : T dg = 0.0;
1638 15995 : for (std::size_t i = 0; i < _n1.size(); ++i)
1639 16286 : dg += -_n1[i] * _I1[i] * _J1[i] * MathUtils::pow(7.1 - pi, _I1[i] - 1) *
1640 15164 : MathUtils::pow(tau - 1.222, _J1[i] - 1);
1641 :
1642 457 : return dg;
1643 : }
1644 :
1645 : template <typename T>
1646 : T
1647 12 : Water97FluidProperties::gamma2(const T & pi, const T & tau) const
1648 : {
1649 : using std::log;
1650 :
1651 : // Ideal gas part of the Gibbs free energy
1652 0 : T sum0 = 0.0;
1653 120 : for (std::size_t i = 0; i < _n02.size(); ++i)
1654 108 : sum0 += _n02[i] * MathUtils::pow(tau, _J02[i]);
1655 :
1656 13 : T g0 = log(pi) + sum0;
1657 :
1658 : // Residual part of the Gibbs free energy
1659 0 : T gr = 0.0;
1660 528 : for (std::size_t i = 0; i < _n2.size(); ++i)
1661 516 : gr += _n2[i] * MathUtils::pow(pi, _I2[i]) * MathUtils::pow(tau - 0.5, _J2[i]);
1662 :
1663 12 : return g0 + gr;
1664 : }
1665 :
1666 : template <typename T>
1667 : T
1668 5 : Water97FluidProperties::d2gamma2_dpi2(const T & pi, const T & tau) const
1669 : {
1670 : // Ideal gas part of the Gibbs free energy
1671 5 : T dg0 = -1.0 / pi / pi;
1672 :
1673 : // Residual part of the Gibbs free energy
1674 0 : T dgr = 0.0;
1675 220 : for (std::size_t i = 0; i < _n2.size(); ++i)
1676 215 : dgr += _n2[i] * _I2[i] * (_I2[i] - 1) * MathUtils::pow(pi, _I2[i] - 2) *
1677 215 : MathUtils::pow(tau - 0.5, _J2[i]);
1678 :
1679 5 : return dg0 + dgr;
1680 : }
1681 :
1682 : template <typename T>
1683 : T
1684 6190 : Water97FluidProperties::dgamma2_dtau(const T & pi, const T & tau) const
1685 : {
1686 : // Ideal gas part of the Gibbs free energy
1687 0 : T dg0 = 0.0;
1688 61900 : for (std::size_t i = 0; i < _n02.size(); ++i)
1689 55710 : dg0 += _n02[i] * _J02[i] * MathUtils::pow(tau, _J02[i] - 1);
1690 :
1691 : // Residual part of the Gibbs free energy
1692 0 : T dgr = 0.0;
1693 272360 : for (std::size_t i = 0; i < _n2.size(); ++i)
1694 266170 : dgr += _n2[i] * _J2[i] * MathUtils::pow(pi, _I2[i]) * MathUtils::pow(tau - 0.5, _J2[i] - 1);
1695 :
1696 6190 : return dg0 + dgr;
1697 : }
1698 :
1699 : template <typename T>
1700 : T
1701 8 : Water97FluidProperties::d2gamma2_dtau2(const T & pi, const T & tau) const
1702 : {
1703 : // Ideal gas part of the Gibbs free energy
1704 0 : T dg0 = 0.0;
1705 80 : for (std::size_t i = 0; i < _n02.size(); ++i)
1706 72 : dg0 += _n02[i] * _J02[i] * (_J02[i] - 1) * MathUtils::pow(tau, _J02[i] - 2);
1707 :
1708 : // Residual part of the Gibbs free energy
1709 0 : T dgr = 0.0;
1710 352 : for (std::size_t i = 0; i < _n2.size(); ++i)
1711 344 : dgr += _n2[i] * _J2[i] * (_J2[i] - 1) * MathUtils::pow(pi, _I2[i]) *
1712 344 : MathUtils::pow(tau - 0.5, _J2[i] - 2);
1713 :
1714 8 : return dg0 + dgr;
1715 : }
1716 :
1717 : template <typename T>
1718 : T
1719 5 : Water97FluidProperties::d2gamma2_dpitau(const T & pi, const T & tau) const
1720 : {
1721 : // Ideal gas part of the Gibbs free energy
1722 0 : T dg0 = 0.0;
1723 :
1724 : // Residual part of the Gibbs free energy
1725 0 : T dgr = 0.0;
1726 220 : for (std::size_t i = 0; i < _n2.size(); ++i)
1727 215 : dgr += _n2[i] * _I2[i] * _J2[i] * MathUtils::pow(pi, _I2[i] - 1) *
1728 215 : MathUtils::pow(tau - 0.5, _J2[i] - 1);
1729 :
1730 5 : return dg0 + dgr;
1731 : }
1732 :
1733 : template <typename T>
1734 : T
1735 12 : Water97FluidProperties::phi3(const T & delta, const T & tau) const
1736 : {
1737 : using std::log;
1738 :
1739 0 : T sum = 0.0;
1740 480 : for (std::size_t i = 1; i < _n3.size(); ++i)
1741 468 : sum += _n3[i] * MathUtils::pow(delta, _I3[i]) * MathUtils::pow(tau, _J3[i]);
1742 :
1743 12 : return _n3[0] * log(delta) + sum;
1744 : }
1745 :
1746 : template <typename T>
1747 : T
1748 19 : Water97FluidProperties::dphi3_ddelta(const T & delta, const T & tau) const
1749 : {
1750 0 : T sum = 0.0;
1751 760 : for (std::size_t i = 1; i < _n3.size(); ++i)
1752 741 : sum += _n3[i] * _I3[i] * MathUtils::pow(delta, _I3[i] - 1) * MathUtils::pow(tau, _J3[i]);
1753 :
1754 19 : return _n3[0] / delta + sum;
1755 : }
1756 :
1757 : template <typename T>
1758 : T
1759 7 : Water97FluidProperties::d2phi3_ddelta2(const T & delta, const T & tau) const
1760 : {
1761 0 : T sum = 0.0;
1762 280 : for (std::size_t i = 1; i < _n3.size(); ++i)
1763 273 : sum += _n3[i] * _I3[i] * (_I3[i] - 1) * MathUtils::pow(delta, _I3[i] - 2) *
1764 273 : MathUtils::pow(tau, _J3[i]);
1765 :
1766 7 : return -_n3[0] / delta / delta + sum;
1767 : }
1768 :
1769 : template <typename T>
1770 : T
1771 33 : Water97FluidProperties::dphi3_dtau(const T & delta, const T & tau) const
1772 : {
1773 0 : T sum = 0.0;
1774 1320 : for (std::size_t i = 1; i < _n3.size(); ++i)
1775 1287 : sum += _n3[i] * _J3[i] * MathUtils::pow(delta, _I3[i]) * MathUtils::pow(tau, _J3[i] - 1);
1776 :
1777 33 : return sum;
1778 : }
1779 :
1780 : template <typename T>
1781 : T
1782 8 : Water97FluidProperties::d2phi3_dtau2(const T & delta, const T & tau) const
1783 : {
1784 0 : T sum = 0.0;
1785 320 : for (std::size_t i = 1; i < _n3.size(); ++i)
1786 312 : sum += _n3[i] * _J3[i] * (_J3[i] - 1) * MathUtils::pow(delta, _I3[i]) *
1787 312 : MathUtils::pow(tau, _J3[i] - 2);
1788 :
1789 8 : return sum;
1790 : }
1791 :
1792 : template <typename T>
1793 : T
1794 7 : Water97FluidProperties::d2phi3_ddeltatau(const T & delta, const T & tau) const
1795 : {
1796 0 : T sum = 0.0;
1797 280 : for (std::size_t i = 1; i < _n3.size(); ++i)
1798 273 : sum += _n3[i] * _I3[i] * _J3[i] * MathUtils::pow(delta, _I3[i] - 1) *
1799 273 : MathUtils::pow(tau, _J3[i] - 1);
1800 :
1801 7 : return sum;
1802 : }
1803 :
1804 : template <typename T>
1805 : T
1806 9 : Water97FluidProperties::gamma5(const T & pi, const T & tau) const
1807 : {
1808 : using std::log;
1809 :
1810 : // Ideal gas part of the Gibbs free energy
1811 0 : T sum0 = 0.0;
1812 63 : for (std::size_t i = 0; i < _n05.size(); ++i)
1813 54 : sum0 += _n05[i] * MathUtils::pow(tau, _J05[i]);
1814 :
1815 10 : T g0 = log(pi) + sum0;
1816 :
1817 : // Residual part of the Gibbs free energy
1818 0 : T gr = 0.0;
1819 63 : for (std::size_t i = 0; i < _n5.size(); ++i)
1820 54 : gr += _n5[i] * MathUtils::pow(pi, _I5[i]) * MathUtils::pow(tau, _J5[i]);
1821 :
1822 9 : return g0 + gr;
1823 : }
1824 :
1825 : template <typename T>
1826 : T
1827 25 : Water97FluidProperties::dgamma5_dpi(const T & pi, const T & tau) const
1828 : {
1829 : // Ideal gas part of the Gibbs free energy
1830 25 : T dg0 = 1.0 / pi;
1831 :
1832 : // Residual part of the Gibbs free energy
1833 0 : T dgr = 0.0;
1834 175 : for (std::size_t i = 0; i < _n5.size(); ++i)
1835 150 : dgr += _n5[i] * _I5[i] * MathUtils::pow(pi, _I5[i] - 1) * MathUtils::pow(tau, _J5[i]);
1836 :
1837 25 : return dg0 + dgr;
1838 : }
1839 :
1840 : template <typename T>
1841 : T
1842 5 : Water97FluidProperties::d2gamma5_dpi2(const T & pi, const T & tau) const
1843 : {
1844 : // Ideal gas part of the Gibbs free energy
1845 5 : T dg0 = -1.0 / pi / pi;
1846 :
1847 : // Residual part of the Gibbs free energy
1848 0 : T dgr = 0.0;
1849 35 : for (std::size_t i = 0; i < _n5.size(); ++i)
1850 30 : dgr += _n5[i] * _I5[i] * (_I5[i] - 1) * MathUtils::pow(pi, _I5[i] - 2) *
1851 30 : MathUtils::pow(tau, _J5[i]);
1852 :
1853 5 : return dg0 + dgr;
1854 : }
1855 :
1856 : template <typename T>
1857 : T
1858 27 : Water97FluidProperties::dgamma5_dtau(const T & pi, const T & tau) const
1859 : {
1860 : // Ideal gas part of the Gibbs free energy
1861 0 : T dg0 = 0.0;
1862 189 : for (std::size_t i = 0; i < _n05.size(); ++i)
1863 162 : dg0 += _n05[i] * _J05[i] * MathUtils::pow(tau, _J05[i] - 1);
1864 :
1865 : // Residual part of the Gibbs free energy
1866 0 : T dgr = 0.0;
1867 189 : for (std::size_t i = 0; i < _n5.size(); ++i)
1868 162 : dgr += _n5[i] * _J5[i] * MathUtils::pow(pi, _I5[i]) * MathUtils::pow(tau, _J5[i] - 1);
1869 :
1870 27 : return dg0 + dgr;
1871 : }
1872 :
1873 : template <typename T>
1874 : T
1875 8 : Water97FluidProperties::d2gamma5_dtau2(const T & pi, const T & tau) const
1876 : {
1877 : // Ideal gas part of the Gibbs free energy
1878 0 : T dg0 = 0.0;
1879 56 : for (std::size_t i = 0; i < _n05.size(); ++i)
1880 48 : dg0 += _n05[i] * _J05[i] * (_J05[i] - 1) * MathUtils::pow(tau, _J05[i] - 2);
1881 :
1882 : // Residual part of the Gibbs free energy
1883 0 : T dgr = 0.0;
1884 56 : for (std::size_t i = 0; i < _n5.size(); ++i)
1885 48 : dgr += _n5[i] * _J5[i] * (_J5[i] - 1) * MathUtils::pow(pi, _I5[i]) *
1886 48 : MathUtils::pow(tau, _J5[i] - 2);
1887 :
1888 8 : return dg0 + dgr;
1889 : }
1890 :
1891 : template <typename T>
1892 : T
1893 5 : Water97FluidProperties::d2gamma5_dpitau(const T & pi, const T & tau) const
1894 : {
1895 : // Ideal gas part of the Gibbs free energy
1896 0 : T dg0 = 0.0;
1897 :
1898 : // Residual part of the Gibbs free energy
1899 0 : T dgr = 0.0;
1900 35 : for (std::size_t i = 0; i < _n5.size(); ++i)
1901 30 : dgr +=
1902 30 : _n5[i] * _I5[i] * _J5[i] * MathUtils::pow(pi, _I5[i] - 1) * MathUtils::pow(tau, _J5[i] - 1);
1903 :
1904 5 : return dg0 + dgr;
1905 : }
1906 :
1907 : template <typename T>
1908 : T
1909 682 : Water97FluidProperties::tempXY(const T & pressure, const subregionEnum xy) const
1910 : {
1911 682 : T pi = pressure / 1.0e6;
1912 :
1913 : // Choose the constants based on the string xy
1914 : unsigned int row;
1915 :
1916 : switch (xy)
1917 : {
1918 : case AB:
1919 : row = 0;
1920 : break;
1921 : case CD:
1922 : row = 1;
1923 : break;
1924 : case GH:
1925 : row = 2;
1926 : break;
1927 : case IJ:
1928 : row = 3;
1929 : break;
1930 : case JK:
1931 : row = 4;
1932 : break;
1933 : case MN:
1934 : row = 5;
1935 : break;
1936 : case OP:
1937 : row = 6;
1938 : break;
1939 : case QU:
1940 : row = 7;
1941 : break;
1942 : case RX:
1943 : row = 8;
1944 : break;
1945 : case UV:
1946 : row = 9;
1947 : break;
1948 : case WX:
1949 : row = 10;
1950 : break;
1951 : default:
1952 : row = 0;
1953 : }
1954 :
1955 : T sum = 0.0;
1956 :
1957 : using std::log;
1958 :
1959 682 : if (xy == AB || xy == OP || xy == WX)
1960 384 : for (std::size_t i = 0; i < _tempXY_n[row].size(); ++i)
1961 320 : sum += _tempXY_n[row][i] * MathUtils::pow(log(pi), _tempXY_I[row][i]);
1962 618 : else if (xy == EF)
1963 76 : sum += 3.727888004 * (pi - _p_critical / 1.0e6) + _T_critical;
1964 : else
1965 2826 : for (std::size_t i = 0; i < _tempXY_n[row].size(); ++i)
1966 2284 : sum += _tempXY_n[row][i] * MathUtils::pow(pi, _tempXY_I[row][i]);
1967 :
1968 682 : return sum;
1969 : }
1970 :
1971 : template <typename T>
1972 : void
1973 1 : Water97FluidProperties::vaporPressureTemplate(const T & temperature, T & psat, T & dpsat_dT) const
1974 : {
1975 : using std::sqrt, std::pow;
1976 :
1977 : // Check whether the input temperature is within the region of validity of this equation.
1978 : // Valid for 273.15 K <= t <= 647.096 K
1979 1 : if (temperature < 273.15 || temperature > _T_critical)
1980 0 : mooseException(name(),
1981 : ": vaporPressure(): Temperature is outside range 273.15 K <= T <= 647.096 K");
1982 :
1983 : T theta, dtheta_dT, theta2, a, b, c, da_dtheta, db_dtheta, dc_dtheta;
1984 1 : theta = temperature + _n4[8] / (temperature - _n4[9]);
1985 1 : dtheta_dT = 1.0 - _n4[8] / (temperature - _n4[9]) / (temperature - _n4[9]);
1986 1 : theta2 = theta * theta;
1987 :
1988 1 : a = theta2 + _n4[0] * theta + _n4[1];
1989 1 : b = _n4[2] * theta2 + _n4[3] * theta + _n4[4];
1990 1 : c = _n4[5] * theta2 + _n4[6] * theta + _n4[7];
1991 :
1992 1 : da_dtheta = 2.0 * theta + _n4[0];
1993 1 : db_dtheta = 2.0 * _n4[2] * theta + _n4[3];
1994 1 : dc_dtheta = 2.0 * _n4[5] * theta + _n4[6];
1995 :
1996 1 : T denominator = -b + sqrt(b * b - 4.0 * a * c);
1997 :
1998 1 : psat = Utility::pow<4>(2.0 * c / denominator) * 1.0e6;
1999 :
2000 : // The derivative wrt temperature is given by the chain rule
2001 1 : T dpsat = 4.0 * Utility::pow<3>(2.0 * c / denominator);
2002 1 : dpsat *= (2.0 * dc_dtheta / denominator -
2003 1 : 2.0 * c / denominator / denominator *
2004 1 : (-db_dtheta + pow(b * b - 4.0 * a * c, -0.5) *
2005 1 : (b * db_dtheta - 2.0 * da_dtheta * c - 2.0 * a * dc_dtheta)));
2006 1 : dpsat_dT = dpsat * dtheta_dT * 1.0e6;
2007 1 : }
2008 :
2009 : inline void
2010 1 : Water97FluidProperties::vaporPressure(const Real temperature, Real & psat, Real & dpsat_dT) const
2011 : {
2012 1 : vaporPressureTemplate(temperature, psat, dpsat_dT);
2013 1 : }
2014 :
2015 : template <typename T>
2016 : void
2017 256 : Water97FluidProperties::rho_from_p_T_template(
2018 : const T & pressure, const T & temperature, T & rho, T & drho_dp, T & drho_dT) const
2019 : {
2020 316 : auto functor = [this](const auto & pressure, const auto & temperature)
2021 316 : { return this->rho_from_p_T_template(pressure, temperature); };
2022 :
2023 276 : xyDerivatives(pressure, temperature, rho, drho_dp, drho_dT, functor);
2024 256 : }
2025 :
2026 : template <typename T>
2027 : T
2028 16 : Water97FluidProperties::v_from_p_T_template(const T & pressure, const T & temperature) const
2029 : {
2030 16 : return 1 / rho_from_p_T_template(pressure, temperature);
2031 : }
2032 :
2033 : template <typename T>
2034 : void
2035 12 : Water97FluidProperties::v_from_p_T_template(
2036 : const T & pressure, const T & temperature, T & v, T & dv_dp, T & dv_dT) const
2037 : {
2038 16 : auto functor = [this](const auto & pressure, const auto & temperature)
2039 16 : { return this->v_from_p_T_template(pressure, temperature); };
2040 :
2041 16 : xyDerivatives(pressure, temperature, v, dv_dp, dv_dT, functor);
2042 12 : }
2043 :
2044 : template <typename T>
2045 : T
2046 220316 : Water97FluidProperties::e_from_p_T_template(const T & pressure, const T & temperature) const
2047 : {
2048 : T internal_energy, pi, tau;
2049 :
2050 : // Determine which region the point is in
2051 : unsigned int region =
2052 220316 : inRegion(MetaPhysicL::raw_value(pressure), MetaPhysicL::raw_value(temperature));
2053 220316 : switch (region)
2054 : {
2055 217231 : case 1:
2056 217231 : pi = pressure / _p_star[0];
2057 217231 : tau = _T_star[0] / temperature;
2058 217311 : internal_energy =
2059 217231 : _Rw * temperature * (tau * dgamma1_dtau(pi, tau) - pi * dgamma1_dpi(pi, tau));
2060 217231 : break;
2061 :
2062 3067 : case 2:
2063 3067 : pi = pressure / _p_star[1];
2064 3067 : tau = _T_star[1] / temperature;
2065 3067 : internal_energy =
2066 3067 : _Rw * temperature * (tau * dgamma2_dtau(pi, tau) - pi * dgamma2_dpi(pi, tau));
2067 3067 : break;
2068 :
2069 9 : case 3:
2070 : {
2071 : // Calculate density first, then use that in Helmholtz free energy
2072 9 : T density3 = densityRegion3(pressure, temperature);
2073 9 : T delta = density3 / _rho_critical;
2074 9 : tau = _T_star[2] / temperature;
2075 9 : internal_energy = _Rw * temperature * tau * dphi3_dtau(delta, tau);
2076 : break;
2077 : }
2078 :
2079 9 : case 5:
2080 9 : pi = pressure / _p_star[4];
2081 9 : tau = _T_star[4] / temperature;
2082 9 : internal_energy =
2083 9 : _Rw * temperature * (tau * dgamma5_dtau(pi, tau) - pi * dgamma5_dpi(pi, tau));
2084 9 : break;
2085 :
2086 0 : default:
2087 0 : mooseError("inRegion() has given an incorrect region");
2088 : }
2089 : // Output in J/kg
2090 220316 : return internal_energy;
2091 : }
2092 :
2093 : template <typename T>
2094 : void
2095 8 : Water97FluidProperties::e_from_p_T_template(
2096 : const T & pressure, const T & temperature, T & e, T & de_dp, T & de_dT) const
2097 : {
2098 : T pi, tau, dinternal_energy_dp, dinternal_energy_dT;
2099 :
2100 : // Determine which region the point is in
2101 : unsigned int region =
2102 8 : inRegion(MetaPhysicL::raw_value(pressure), MetaPhysicL::raw_value(temperature));
2103 8 : switch (region)
2104 : {
2105 5 : case 1:
2106 : {
2107 5 : pi = pressure / _p_star[0];
2108 5 : tau = _T_star[0] / temperature;
2109 5 : T dgdp = dgamma1_dpi(pi, tau);
2110 5 : T d2gdpt = d2gamma1_dpitau(pi, tau);
2111 6 : dinternal_energy_dp =
2112 5 : _Rw * temperature * (tau * d2gdpt - dgdp - pi * d2gamma1_dpi2(pi, tau)) / _p_star[0];
2113 6 : dinternal_energy_dT =
2114 4 : _Rw * (pi * tau * d2gdpt - tau * tau * d2gamma1_dtau2(pi, tau) - pi * dgdp);
2115 4 : break;
2116 : }
2117 :
2118 1 : case 2:
2119 : {
2120 1 : pi = pressure / _p_star[1];
2121 1 : tau = _T_star[1] / temperature;
2122 1 : T dgdp = dgamma2_dpi(pi, tau);
2123 1 : T d2gdpt = d2gamma2_dpitau(pi, tau);
2124 1 : dinternal_energy_dp =
2125 1 : _Rw * temperature * (tau * d2gdpt - dgdp - pi * d2gamma2_dpi2(pi, tau)) / _p_star[1];
2126 1 : dinternal_energy_dT =
2127 1 : _Rw * (pi * tau * d2gdpt - tau * tau * d2gamma2_dtau2(pi, tau) - pi * dgdp);
2128 1 : break;
2129 : }
2130 :
2131 1 : case 3:
2132 : {
2133 : // Calculate density first, then use that in Helmholtz free energy
2134 1 : T density3 = densityRegion3(pressure, temperature);
2135 1 : T delta = density3 / _rho_critical;
2136 1 : tau = _T_star[2] / temperature;
2137 1 : T dpdd = dphi3_ddelta(delta, tau);
2138 1 : T d2pddt = d2phi3_ddeltatau(delta, tau);
2139 1 : T d2pdd2 = d2phi3_ddelta2(delta, tau);
2140 1 : dinternal_energy_dp =
2141 1 : _T_star[2] * d2pddt / _rho_critical /
2142 1 : (2.0 * temperature * delta * dpdd + temperature * delta * delta * d2pdd2);
2143 1 : dinternal_energy_dT =
2144 1 : -_Rw * (delta * tau * d2pddt * (dpdd - tau * d2pddt) / (2.0 * dpdd + delta * d2pdd2) +
2145 1 : tau * tau * d2phi3_dtau2(delta, tau));
2146 : break;
2147 : }
2148 :
2149 1 : case 5:
2150 : {
2151 1 : pi = pressure / _p_star[4];
2152 1 : tau = _T_star[4] / temperature;
2153 1 : T dgdp = dgamma5_dpi(pi, tau);
2154 1 : T d2gdpt = d2gamma5_dpitau(pi, tau);
2155 1 : dinternal_energy_dp =
2156 1 : _Rw * temperature * (tau * d2gdpt - dgdp - pi * d2gamma5_dpi2(pi, tau)) / _p_star[4];
2157 1 : dinternal_energy_dT =
2158 1 : _Rw * (pi * tau * d2gdpt - tau * tau * d2gamma5_dtau2(pi, tau) - pi * dgdp);
2159 1 : break;
2160 : }
2161 :
2162 0 : default:
2163 0 : mooseError("inRegion has given an incorrect region");
2164 : }
2165 :
2166 8 : e = this->e_from_p_T(pressure, temperature);
2167 8 : de_dp = dinternal_energy_dp;
2168 8 : de_dT = dinternal_energy_dT;
2169 8 : }
2170 :
2171 : template <typename T>
2172 : T
2173 100 : Water97FluidProperties::subregionVolume(
2174 : const T & pi, const T & theta, Real a, Real b, Real c, Real d, Real e, unsigned int sid) const
2175 : {
2176 : using std::pow;
2177 :
2178 0 : T sum = 0.0;
2179 :
2180 3394 : for (std::size_t i = 0; i < _n3s[sid].size(); ++i)
2181 3294 : sum += _n3s[sid][i] * MathUtils::pow(pow(pi - a, c), _I3s[sid][i]) *
2182 3189 : MathUtils::pow(pow(theta - b, d), _J3s[sid][i]);
2183 :
2184 100 : return pow(sum, e);
2185 : }
2186 :
2187 : template <typename T>
2188 : std::pair<T, T>
2189 74 : Water97FluidProperties::T_drhodT_from_p_rho(const T & p, const T & rho) const
2190 : {
2191 : // For NewtonSolve
2192 : // y = rho
2193 : // x = p
2194 : // z = T
2195 74 : auto lambda =
2196 256 : [this](const T & p, const T & temperature, T & rho, T & drho_dp, T & drho_dtemperature)
2197 296 : { rho_from_p_T_template(p, temperature, rho, drho_dp, drho_dtemperature); };
2198 10 : return FluidPropertiesUtils::NewtonSolve(
2199 148 : p, rho, _T_initial_guess, _tolerance, lambda, name() + "::T_drhodT_from_p_rho");
2200 : }
2201 :
2202 : template <typename T>
2203 : std::pair<T, T>
2204 18 : Water97FluidProperties::p_T_from_v_e(const T & v, const T & e) const
2205 : {
2206 18 : const auto rho = 1 / v;
2207 18 : const auto p = p_from_v_e(v, e);
2208 18 : return {p, T_drhodT_from_p_rho(p, rho).first};
2209 : }
2210 :
2211 : template <typename T>
2212 : std::pair<T, T>
2213 4 : Water97FluidProperties::rho_T_from_v_e(const T & v, const T & e) const
2214 : {
2215 4 : const auto rho = 1 / v;
2216 4 : const auto p = p_from_v_e(v, e);
2217 4 : return {rho, T_drhodT_from_p_rho(p, rho).first};
2218 : }
2219 :
2220 : template <typename T>
2221 : T
2222 236 : Water97FluidProperties::c_from_p_T_template(const T & pressure, const T & temperature) const
2223 : {
2224 : T speed2, pi, tau, delta;
2225 :
2226 : // Determine which region the point is in
2227 : unsigned int region =
2228 236 : inRegion(MetaPhysicL::raw_value(pressure), MetaPhysicL::raw_value(temperature));
2229 236 : switch (region)
2230 : {
2231 227 : case 1:
2232 227 : pi = pressure / _p_star[0];
2233 227 : tau = _T_star[0] / temperature;
2234 252 : speed2 = _Rw * temperature * Utility::pow<2>(dgamma1_dpi(pi, tau)) /
2235 222 : (Utility::pow<2>(dgamma1_dpi(pi, tau) - tau * d2gamma1_dpitau(pi, tau)) /
2236 222 : (tau * tau * d2gamma1_dtau2(pi, tau)) -
2237 222 : d2gamma1_dpi2(pi, tau));
2238 227 : break;
2239 :
2240 3 : case 2:
2241 3 : pi = pressure / _p_star[1];
2242 3 : tau = _T_star[1] / temperature;
2243 3 : speed2 = _Rw * temperature * Utility::pow<2>(pi * dgamma2_dpi(pi, tau)) /
2244 3 : ((-pi * pi * d2gamma2_dpi2(pi, tau)) +
2245 3 : Utility::pow<2>(pi * dgamma2_dpi(pi, tau) - tau * pi * d2gamma2_dpitau(pi, tau)) /
2246 3 : (tau * tau * d2gamma2_dtau2(pi, tau)));
2247 3 : break;
2248 :
2249 3 : case 3:
2250 : {
2251 : // Calculate density first, then use that in Helmholtz free energy
2252 3 : T density3 = densityRegion3(pressure, temperature);
2253 3 : delta = density3 / _rho_critical;
2254 3 : tau = _T_star[2] / temperature;
2255 3 : speed2 =
2256 3 : _Rw * temperature *
2257 3 : (2.0 * delta * dphi3_ddelta(delta, tau) + delta * delta * d2phi3_ddelta2(delta, tau) -
2258 3 : Utility::pow<2>(delta * dphi3_ddelta(delta, tau) -
2259 3 : delta * tau * d2phi3_ddeltatau(delta, tau)) /
2260 3 : (tau * tau * d2phi3_dtau2(delta, tau)));
2261 3 : break;
2262 : }
2263 :
2264 3 : case 5:
2265 3 : pi = pressure / _p_star[4];
2266 3 : tau = _T_star[4] / temperature;
2267 3 : speed2 = _Rw * temperature * Utility::pow<2>(pi * dgamma5_dpi(pi, tau)) /
2268 3 : ((-pi * pi * d2gamma5_dpi2(pi, tau)) +
2269 3 : Utility::pow<2>(pi * dgamma5_dpi(pi, tau) - tau * pi * d2gamma5_dpitau(pi, tau)) /
2270 3 : (tau * tau * d2gamma5_dtau2(pi, tau)));
2271 3 : break;
2272 :
2273 0 : default:
2274 0 : mooseError("inRegion() has given an incorrect region");
2275 : }
2276 :
2277 : using std::sqrt;
2278 236 : return sqrt(speed2);
2279 : }
2280 :
2281 : template <typename T>
2282 : T
2283 671 : Water97FluidProperties::cp_from_p_T_template(const T & pressure, const T & temperature) const
2284 : {
2285 : T specific_heat, pi, tau, delta;
2286 :
2287 : // Determine which region the point is in
2288 : unsigned int region =
2289 671 : inRegion(MetaPhysicL::raw_value(pressure), MetaPhysicL::raw_value(temperature));
2290 671 : switch (region)
2291 : {
2292 662 : case 1:
2293 662 : pi = pressure / _p_star[0];
2294 662 : tau = _T_star[0] / temperature;
2295 667 : specific_heat = -_Rw * tau * tau * d2gamma1_dtau2(pi, tau);
2296 662 : break;
2297 :
2298 3 : case 2:
2299 3 : pi = pressure / _p_star[1];
2300 3 : tau = _T_star[1] / temperature;
2301 3 : specific_heat = -_Rw * tau * tau * d2gamma2_dtau2(pi, tau);
2302 3 : break;
2303 :
2304 3 : case 3:
2305 : {
2306 : // Calculate density first, then use that in Helmholtz free energy
2307 3 : T density3 = densityRegion3(pressure, temperature);
2308 3 : delta = density3 / _rho_critical;
2309 3 : tau = _T_star[2] / temperature;
2310 3 : specific_heat =
2311 3 : _Rw *
2312 3 : (-tau * tau * d2phi3_dtau2(delta, tau) +
2313 3 : (delta * dphi3_ddelta(delta, tau) - delta * tau * d2phi3_ddeltatau(delta, tau)) *
2314 3 : (delta * dphi3_ddelta(delta, tau) - delta * tau * d2phi3_ddeltatau(delta, tau)) /
2315 3 : (2.0 * delta * dphi3_ddelta(delta, tau) +
2316 3 : delta * delta * d2phi3_ddelta2(delta, tau)));
2317 3 : break;
2318 : }
2319 :
2320 3 : case 5:
2321 3 : pi = pressure / _p_star[4];
2322 3 : tau = _T_star[4] / temperature;
2323 3 : specific_heat = -_Rw * tau * tau * d2gamma5_dtau2(pi, tau);
2324 3 : break;
2325 :
2326 0 : default:
2327 0 : mooseError("inRegion() has given an incorrect region");
2328 : }
2329 671 : return specific_heat;
2330 : }
2331 :
2332 : template <typename T>
2333 : T
2334 228 : Water97FluidProperties::cv_from_p_T_template(const T & pressure, const T & temperature) const
2335 : {
2336 : T specific_heat, pi, tau, delta;
2337 :
2338 : // Determine which region the point is in
2339 : unsigned int region =
2340 228 : inRegion(MetaPhysicL::raw_value(pressure), MetaPhysicL::raw_value(temperature));
2341 228 : switch (region)
2342 : {
2343 225 : case 1:
2344 225 : pi = pressure / _p_star[0];
2345 225 : tau = _T_star[0] / temperature;
2346 240 : specific_heat =
2347 225 : _Rw * (-tau * tau * d2gamma1_dtau2(pi, tau) +
2348 220 : Utility::pow<2>(dgamma1_dpi(pi, tau) - tau * d2gamma1_dpitau(pi, tau)) /
2349 220 : d2gamma1_dpi2(pi, tau));
2350 225 : break;
2351 :
2352 1 : case 2:
2353 1 : pi = pressure / _p_star[1];
2354 1 : tau = _T_star[1] / temperature;
2355 1 : specific_heat =
2356 1 : _Rw * (-tau * tau * d2gamma2_dtau2(pi, tau) +
2357 1 : Utility::pow<2>(dgamma2_dpi(pi, tau) - tau * d2gamma2_dpitau(pi, tau)) /
2358 1 : d2gamma2_dpi2(pi, tau));
2359 1 : break;
2360 :
2361 1 : case 3:
2362 : {
2363 : // Calculate density first, then use that in Helmholtz free energy
2364 1 : T density3 = densityRegion3(pressure, temperature);
2365 1 : delta = density3 / _rho_critical;
2366 1 : tau = _T_star[2] / temperature;
2367 1 : specific_heat = _Rw * (-tau * tau * d2phi3_dtau2(delta, tau));
2368 1 : break;
2369 : }
2370 :
2371 1 : case 5:
2372 1 : pi = pressure / _p_star[4];
2373 1 : tau = _T_star[4] / temperature;
2374 1 : specific_heat =
2375 1 : _Rw * (-tau * tau * d2gamma5_dtau2(pi, tau) +
2376 1 : Utility::pow<2>(dgamma5_dpi(pi, tau) - tau * d2gamma5_dpitau(pi, tau)) /
2377 1 : d2gamma5_dpi2(pi, tau));
2378 1 : break;
2379 :
2380 0 : default:
2381 0 : mooseError("inRegion() has given an incorrect region");
2382 : }
2383 228 : return specific_heat;
2384 : }
2385 :
2386 : template <typename T>
2387 : T
2388 241 : Water97FluidProperties::k_from_rho_T_template(const T & density, const T & temperature) const
2389 : {
2390 : using std::sqrt, std::exp, std::abs, std::pow;
2391 : // Scale the density and temperature. Note that the scales are slightly
2392 : // different to the critical values used in IAPWS-IF97
2393 241 : T Tbar = temperature / 647.26;
2394 241 : T rhobar = density / 317.7;
2395 :
2396 : // Ideal gas component
2397 8 : T sum0 = 0.0;
2398 :
2399 1205 : for (std::size_t i = 0; i < _k_a.size(); ++i)
2400 996 : sum0 += _k_a[i] * MathUtils::pow(Tbar, i);
2401 :
2402 249 : T lambda0 = sqrt(Tbar) * sum0;
2403 :
2404 : // The contribution due to finite density
2405 241 : T lambda1 =
2406 273 : -0.39707 + 0.400302 * rhobar + 1.06 * exp(-0.171587 * Utility::pow<2>(rhobar + 2.392190));
2407 :
2408 : // Critical enhancement
2409 241 : T DeltaT = abs(Tbar - 1.0) + 0.00308976;
2410 257 : T Q = 2.0 + 0.0822994 / pow(DeltaT, 0.6);
2411 249 : T S = (Tbar >= 1.0 ? 1.0 / DeltaT : 10.0932 / pow(DeltaT, 0.6));
2412 :
2413 265 : T lambda2 =
2414 249 : (0.0701309 / Utility::pow<10>(Tbar) + 0.011852) * pow(rhobar, 1.8) *
2415 249 : exp(0.642857 * (1.0 - pow(rhobar, 2.8))) +
2416 265 : 0.00169937 * S * pow(rhobar, Q) * exp((Q / (1.0 + Q)) * (1.0 - pow(rhobar, 1.0 + Q))) -
2417 265 : 1.02 * exp(-4.11717 * pow(Tbar, 1.5) - 6.17937 / Utility::pow<5>(rhobar));
2418 :
2419 241 : return lambda0 + lambda1 + lambda2;
2420 : }
2421 :
2422 : template <typename T>
2423 : T
2424 241 : Water97FluidProperties::k_from_p_T_template(const T & pressure, const T & temperature) const
2425 : {
2426 241 : T rho = this->rho_from_p_T(pressure, temperature);
2427 241 : return this->k_from_rho_T_template(rho, temperature);
2428 : }
2429 :
2430 : template <typename T>
2431 : T
2432 8 : Water97FluidProperties::k_from_v_e_template(const T & v, const T & e) const
2433 : {
2434 8 : const auto [p, temperature] = p_T_from_v_e(v, e);
2435 8 : return k_from_p_T_template(p, temperature);
2436 : }
2437 :
2438 : template <typename T>
2439 : T
2440 22 : Water97FluidProperties::p_from_v_e_template(const T & v, const T & e) const
2441 : {
2442 22 : const auto rho = 1 / v;
2443 : // For NewtonSolve
2444 : // y = e
2445 : // x = rho
2446 : // z = p
2447 66 : auto lambda = [this](const T & rho, const T & p, T & e, T & de_drho, T & de_dp)
2448 44 : { e_from_p_rho(p, rho, e, de_dp, de_drho); };
2449 7 : return FluidPropertiesUtils::NewtonSolve(
2450 22 : rho, e, _p_initial_guess, _tolerance, lambda, name() + "::p_from_v_e")
2451 37 : .first;
2452 : }
2453 :
2454 : template <typename T>
2455 : T
2456 52 : Water97FluidProperties::e_from_p_rho_template(const T & p, const T & rho) const
2457 : {
2458 52 : const auto temperature = T_drhodT_from_p_rho(p, rho).first;
2459 52 : return e_from_p_T_template(p, temperature);
2460 : }
2461 :
2462 : template <typename T>
2463 : void
2464 31 : Water97FluidProperties::e_from_p_rho_template(
2465 : const T & p, const T & rho, T & e, T & de_dp, T & de_drho) const
2466 : {
2467 46 : auto functor = [this](const auto & p, const auto & rho)
2468 46 : { return this->e_from_p_rho_template(p, rho); };
2469 :
2470 46 : xyDerivatives(p, rho, e, de_dp, de_drho, functor);
2471 31 : }
2472 :
2473 : template <typename T>
2474 : T
2475 220266 : Water97FluidProperties::mu_from_rho_T_template(const T & density, const T & temperature) const
2476 : {
2477 : using std::sqrt, std::exp;
2478 :
2479 5 : const T mu_star = 1.e-6;
2480 220266 : const T rhobar = density / _rho_critical;
2481 220266 : const T Tbar = temperature / _T_critical;
2482 :
2483 : // Viscosity in limit of zero density
2484 5 : T sum0 = 0.0;
2485 1101330 : for (std::size_t i = 0; i < _mu_H0.size(); ++i)
2486 881084 : sum0 += _mu_H0[i] / MathUtils::pow(Tbar, i);
2487 :
2488 220266 : const T mu0 = 100.0 * sqrt(Tbar) / sum0;
2489 :
2490 : // Residual component due to finite density
2491 5 : T sum1 = 0.0;
2492 1541862 : for (unsigned int i = 0; i < 6; ++i)
2493 : {
2494 1321626 : const T fact = MathUtils::pow(1.0 / Tbar - 1.0, i);
2495 10572768 : for (unsigned int j = 0; j < 7; ++j)
2496 9251382 : sum1 += fact * _mu_Hij[i][j] * MathUtils::pow(rhobar - 1.0, j);
2497 : }
2498 :
2499 220266 : const T mu1 = exp(rhobar * sum1);
2500 :
2501 : // The water viscosity (in Pa.s) is then given by
2502 220266 : return mu_star * mu0 * mu1;
2503 : }
2504 :
2505 : template <typename T>
2506 : std::pair<T, T>
2507 4 : Water97FluidProperties::p_T_from_v_h(const T & v, const T & h) const
2508 : {
2509 : T pressure, temperature;
2510 : bool conversion_succeeded;
2511 4 : p_T_from_v_h(
2512 4 : v, h, _p_initial_guess, _T_initial_guess, pressure, temperature, conversion_succeeded);
2513 4 : return {std::move(pressure), std::move(temperature)};
2514 : }
2515 :
2516 : template <typename T>
2517 : T
2518 220883 : Water97FluidProperties::h_from_p_T_template(const T & pressure, const T & temperature) const
2519 : {
2520 : T enthalpy, pi, tau, delta;
2521 :
2522 : // Determine which region the point is in
2523 : unsigned int region =
2524 220883 : inRegion(MetaPhysicL::raw_value(pressure), MetaPhysicL::raw_value(temperature));
2525 220883 : switch (region)
2526 : {
2527 217751 : case 1:
2528 217751 : pi = pressure / _p_star[0];
2529 217751 : tau = _T_star[0] / temperature;
2530 217767 : enthalpy = _Rw * _T_star[0] * dgamma1_dtau(pi, tau);
2531 217751 : break;
2532 :
2533 3111 : case 2:
2534 3111 : pi = pressure / _p_star[1];
2535 3111 : tau = _T_star[1] / temperature;
2536 3111 : enthalpy = _Rw * _T_star[1] * dgamma2_dtau(pi, tau);
2537 3111 : break;
2538 :
2539 12 : case 3:
2540 : {
2541 : // Calculate density first, then use that in Helmholtz free energy
2542 12 : T density3 = densityRegion3(pressure, temperature);
2543 12 : delta = density3 / _rho_critical;
2544 12 : tau = _T_star[2] / temperature;
2545 13 : enthalpy =
2546 12 : _Rw * temperature * (tau * dphi3_dtau(delta, tau) + delta * dphi3_ddelta(delta, tau));
2547 11 : break;
2548 : }
2549 :
2550 9 : case 5:
2551 9 : pi = pressure / _p_star[4];
2552 9 : tau = _T_star[4] / temperature;
2553 9 : enthalpy = _Rw * _T_star[4] * dgamma5_dtau(pi, tau);
2554 9 : break;
2555 :
2556 0 : default:
2557 0 : mooseError("Water97FluidProperties::inRegion has given an incorrect region");
2558 : }
2559 220883 : return enthalpy;
2560 : }
2561 :
2562 : template <typename T>
2563 : void
2564 15 : Water97FluidProperties::h_from_p_T_template(
2565 : const T & pressure, const T & temperature, T & h, T & dh_dp, T & dh_dT) const
2566 : {
2567 30 : auto functor = [this](const auto & pressure, const auto & temperature)
2568 30 : { return this->h_from_p_T_template(pressure, temperature); };
2569 :
2570 30 : xyDerivatives(pressure, temperature, h, dh_dp, dh_dT, functor);
2571 15 : }
2572 :
2573 : template <typename T>
2574 : T
2575 441530 : Water97FluidProperties::rho_from_p_T_template(const T & pressure, const T & temperature) const
2576 : {
2577 124 : T density, pi, tau;
2578 :
2579 : // Determine which region the point is in
2580 : unsigned int region =
2581 441530 : inRegion(MetaPhysicL::raw_value(pressure), MetaPhysicL::raw_value(temperature));
2582 :
2583 441530 : switch (region)
2584 : {
2585 435381 : case 1:
2586 435381 : pi = pressure / _p_star[0];
2587 435381 : tau = _T_star[0] / temperature;
2588 435661 : density = pressure / (pi * _Rw * temperature * dgamma1_dpi(pi, tau));
2589 435381 : break;
2590 :
2591 6129 : case 2:
2592 6129 : pi = pressure / _p_star[1];
2593 6129 : tau = _T_star[1] / temperature;
2594 6129 : density = pressure / (pi * _Rw * temperature * dgamma2_dpi(pi, tau));
2595 6129 : break;
2596 :
2597 9 : case 3:
2598 9 : density = densityRegion3(pressure, temperature);
2599 1 : break;
2600 :
2601 11 : case 5:
2602 11 : pi = pressure / _p_star[4];
2603 11 : tau = _T_star[4] / temperature;
2604 11 : density = pressure / (pi * _Rw * temperature * dgamma5_dpi(pi, tau));
2605 11 : break;
2606 :
2607 0 : default:
2608 0 : mooseError("inRegion() has given an incorrect region");
2609 : }
2610 441530 : return density;
2611 : }
2612 :
2613 : template <typename T>
2614 : T
2615 220247 : Water97FluidProperties::mu_from_p_T_template(const T & pressure, const T & temperature) const
2616 : {
2617 220247 : T rho = this->rho_from_p_T_template(pressure, temperature);
2618 220247 : return this->mu_from_rho_T_template(rho, temperature);
2619 : }
2620 :
2621 : template <typename T>
2622 : T
2623 275 : Water97FluidProperties::s_from_p_T_template(const T & pressure, const T & temperature) const
2624 : {
2625 : T entropy, pi, tau, density3, delta;
2626 :
2627 : // Determine which region the point is in
2628 : unsigned int region =
2629 275 : inRegion(MetaPhysicL::raw_value(pressure), MetaPhysicL::raw_value(temperature));
2630 275 : switch (region)
2631 : {
2632 242 : case 1:
2633 242 : pi = pressure / _p_star[0];
2634 242 : tau = _T_star[0] / temperature;
2635 256 : entropy = _Rw * (tau * dgamma1_dtau(pi, tau) - gamma1(pi, tau));
2636 242 : break;
2637 :
2638 12 : case 2:
2639 12 : pi = pressure / _p_star[1];
2640 12 : tau = _T_star[1] / temperature;
2641 13 : entropy = _Rw * (tau * dgamma2_dtau(pi, tau) - gamma2(pi, tau));
2642 12 : break;
2643 :
2644 12 : case 3:
2645 : // Calculate density first, then use that in Helmholtz free energy
2646 12 : density3 = densityRegion3(pressure, temperature);
2647 12 : delta = density3 / _rho_critical;
2648 12 : tau = _T_star[2] / temperature;
2649 13 : entropy = _Rw * (tau * dphi3_dtau(delta, tau) - phi3(delta, tau));
2650 12 : break;
2651 :
2652 9 : case 5:
2653 9 : pi = pressure / _p_star[4];
2654 9 : tau = _T_star[4] / temperature;
2655 10 : entropy = _Rw * (tau * dgamma5_dtau(pi, tau) - gamma5(pi, tau));
2656 9 : break;
2657 :
2658 0 : default:
2659 0 : mooseError("inRegion() has given an incorrect region");
2660 : }
2661 275 : return entropy;
2662 : }
2663 :
2664 : template <typename T>
2665 : void
2666 3 : Water97FluidProperties::s_from_p_T_template(
2667 : const T & pressure, const T & temperature, T & s, T & ds_dp, T & ds_dT) const
2668 : {
2669 10 : auto functor = [this](const auto & pressure, const auto & temperature)
2670 10 : { return this->s_from_p_T_template(pressure, temperature); };
2671 :
2672 10 : xyDerivatives(pressure, temperature, s, ds_dp, ds_dT, functor);
2673 3 : }
|
