https://mooseframework.inl.gov
Loading...
Searching...
No Matches
CaloricallyImperfectGas.C
Go to the documentation of this file.
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
11#include "Conversion.h"
12#include "metaphysicl/raw_type.h"
13#include "Function.h"
14#include "BrentsMethod.h"
15
17
20{
22 params += NaNInterface::validParams();
23
24 params.addRequiredParam<Real>("molar_mass", "Constant molar mass of the fluid [kg/mol]");
25 params.addRequiredParam<FunctionName>("mu", "Dynamic viscosity, [Pa-s]");
26 params.addRequiredParam<FunctionName>("k", "Thermal conductivity, [W/(m-K)]");
27 params.addParam<Real>("T_c", 0, "Critical temperature, [K]");
28 params.addParam<Real>("rho_c", 0, "Critical density, [kg/m3]");
29 params.addParam<Real>("e_c", 0, "Specific internal energy at the critical point, [J/kg]");
30 params.addRequiredParam<FunctionName>(
31 "e", "Specific internal energy as a function of temperature [J/kg]");
32 params.addRequiredParam<Real>("min_temperature", "Smallest temperature for lookup tables [K]");
33 params.addRequiredParam<Real>("max_temperature", "Largest temperature for lookup tables [K]");
34 params.addParam<Real>(
35 "temperature_resolution", 1, "Step size in temperature for creating the inverse lookup T(e)");
36 params.addParam<bool>("out_of_bound_error", true, "Error if out of bounds");
37 params.addClassDescription("Fluid properties for an ideal gas with imperfect caloric behavior.");
38
39 return params;
40}
41
43 : SinglePhaseFluidProperties(parameters),
44 NaNInterface(this),
45 _molar_mass(getParam<Real>("molar_mass")),
46 _R_specific(_R / _molar_mass),
47 _T_c(getParam<Real>("T_c")),
48 _rho_c(getParam<Real>("rho_c")),
49 _e_c(getParam<Real>("e_c")),
50 _min_temperature(getParam<Real>("min_temperature")),
51 _max_temperature(getParam<Real>("max_temperature")),
52 _delta_T(getParam<Real>("temperature_resolution")),
53 _out_of_bound_error(getParam<bool>("out_of_bound_error")),
54 _tol(1e-4)
55{
56}
57
58void
68
69void
71{
72 // estimate number of points in inverse lookup & then adjust
73 // the _delta_T so (n - 1) * _delta_T = _max_temperature - _min_temperature
74 unsigned int n = std::floor((_max_temperature - _min_temperature) / _delta_T) + 1;
75 _delta_T = (_max_temperature - _min_temperature) / ((Real)n - 1.0);
76
77 // ensure e(T) is monotonic
78 for (unsigned int j = 0; j < n; ++j)
79 {
80 Real temperature = _min_temperature + j * _delta_T;
81
82 // Note that the function behavior at the end points
83 // can lead to de/dT = 0 & we want to allow that
84 Real cv = _e_T->timeDerivative(temperature, Point());
85 if (cv < 0 || (cv == 0 && j > 0 && j < n - 1))
86 mooseError("e(T) is not monotonically increasing with T");
87 }
88
89 // backward lookup tables
90 _min_e = _e_T->value(_min_temperature, Point());
91 _max_e = _e_T->value(_max_temperature, Point());
92 _delta_e = (_max_e - _min_e) / ((Real)n - 1.0);
95 _delta_h = (_max_h - _min_h) / ((Real)n - 1.0);
96
97 _T_e_lookup.resize(n);
98 _T_h_lookup.resize(n);
99 _Z_T_lookup.resize(n);
100 for (unsigned int j = 0; j < n; ++j)
101 {
102 // internal energy
103 {
104 Real internal_energy = _min_e + j * _delta_e;
105
106 // guarding against roundoff when summing
107 // at least once I saw roundoff cause an error here
108 internal_energy = std::min(internal_energy, _max_e);
109
110 // The temperature is found by e - e(T) = 0
111 Real min = _min_temperature;
112 Real max = _max_temperature;
113 auto e_diff = [&internal_energy, this](Real x)
114 { return this->e_from_p_T(0.0, x) - internal_energy; };
115 BrentsMethod::bracket(e_diff, min, max);
116 Real temperature = BrentsMethod::root(e_diff, min, max);
117 _T_e_lookup[j] = temperature;
118 }
119
120 // enthalpy
121 {
122 Real enthalpy = _min_h + j * _delta_h;
123
124 // guarding against roundoff when summing
125 // at least once I saw roundoff cause an error here
126 enthalpy = std::min(enthalpy, _max_h);
127
128 // The temperature is found by h - h(T) = 0
129 Real min = _min_temperature;
130 Real max = _max_temperature;
131 auto h_diff = [&enthalpy, this](Real x) { return this->h_from_p_T(0.0, x) - enthalpy; };
132 BrentsMethod::bracket(h_diff, min, max);
133 Real temperature = BrentsMethod::root(h_diff, min, max);
134 _T_h_lookup[j] = temperature;
135 }
136
137 // Z(T)
138 {
139 if (j == 0)
140 _Z_T_lookup[j] = 0;
141 else
142 {
143 Real temperature = _min_temperature + j * _delta_T;
144 Real temperature_prev = _min_temperature + (j - 1) * _delta_T;
145 Real f1 = cv_from_T(temperature) / temperature;
146 Real f2 = cv_from_T(temperature_prev) / temperature_prev;
147 _Z_T_lookup[j] = _Z_T_lookup[j - 1] + 0.5 * _delta_T * (f1 + f2);
148 }
149 }
150 }
151}
152
153std::string
155{
156 return "caloric_imperfect_gas";
157}
158
159Real
161{
162 if (T < _min_temperature || T > _max_temperature)
164
165 return _e_T->value(T, Point());
166}
167
168Real
170{
171 if (T < _min_temperature || T > _max_temperature)
173 return _e_T->value(T, Point()) + _R_specific * T;
174}
175
176Real
178{
179 if (T < _min_temperature || T > _max_temperature)
181
182 return _e_T->timeDerivative(T, Point()) + _R_specific;
183}
184
185void
186CaloricallyImperfectGas::cv_from_T(Real T, Real & cv, Real & dcv_dT) const
187{
188 if (T < _min_temperature || T > _max_temperature)
189 outOfBounds("cv_from_T (3 args)", T, _min_temperature, _max_temperature);
190 cv = cv_from_T(T);
191 const Real dT = temperatureStep(T);
192 dcv_dT = (cv_from_T(T + dT) - cv) / dT;
193}
194
195void
196CaloricallyImperfectGas::cp_from_T(Real T, Real & cp, Real & dcp_dT) const
197{
198 if (T < _min_temperature || T > _max_temperature)
199 outOfBounds("cp_from_T (3 args)", T, _min_temperature, _max_temperature);
200 cp = cp_from_T(T);
201 const Real dT = temperatureStep(T);
202 dcp_dT = (cp_from_T(T + dT) - cp) / dT;
203}
204
205Real
207{
208 if (e < _min_e || e > _max_e)
209 outOfBounds("T_from_e", e, _min_e, _max_e);
210
211 if (e < _min_e)
212 return _T_e_lookup[0];
213 else if (e > _max_e)
214 return _T_e_lookup.back();
215
216 // e == _max_e falls in the last interval, not past the end of the table
217 const unsigned int index =
218 std::min<unsigned int>(std::floor((e - _min_e) / _delta_e), _T_e_lookup.size() - 2);
219 Real x = (e - _min_e - index * _delta_e) / _delta_e;
220 return x * _T_e_lookup[index + 1] + (1 - x) * _T_e_lookup[index];
221}
222
223Real
225{
226 if (h < _min_h || h > _max_h)
227 outOfBounds("T_from_h", h, _min_h, _max_h);
228
229 if (h < _min_h)
230 return _T_h_lookup[0];
231 else if (h > _max_h)
232 return _T_h_lookup.back();
233
234 // h == _max_h falls in the last interval, not past the end of the table
235 const unsigned int index =
236 std::min<unsigned int>(std::floor((h - _min_h) / _delta_h), _T_h_lookup.size() - 2);
237 Real x = (h - _min_h - index * _delta_h) / _delta_h;
238 return x * _T_h_lookup[index + 1] + (1 - x) * _T_h_lookup[index];
239}
240
241Real
243{
244 if (T < _min_temperature || T > _max_temperature)
246
247 if (T <= _min_temperature)
248 return _Z_T_lookup[0];
249 else if (T >= _max_temperature)
250 return _Z_T_lookup.back();
251
252 const unsigned int index =
253 std::min<unsigned int>(std::floor((T - _min_temperature) / _delta_T), _Z_T_lookup.size() - 2);
254 Real x = (T - _min_temperature - index * _delta_T) / _delta_T;
255 return x * _Z_T_lookup[index + 1] + (1 - x) * _Z_T_lookup[index];
256}
257
258Real
260{
261 if (v == 0.0)
262 return getNaN("Invalid value of specific volume detected (v = " + Moose::stringify(v) + ").");
263
264 return _R_specific * T_from_e(e) / v;
265}
266
267ADReal
269{
270 if (v.value() == 0.0)
271 return getNaN("Invalid value of specific volume detected (v = " + Moose::stringify(v.value()) +
272 ").");
273
274 Real x = 0;
275 Real raw1 = v.value();
276 Real raw2 = e.value();
277 Real dxd1 = 0;
278 Real dxd2 = 0;
279 p_from_v_e(raw1, raw2, x, dxd1, dxd2);
280
281 ADReal result = x;
282 result.derivatives() = v.derivatives() * dxd1 + e.derivatives() * dxd2;
283 return result;
284}
285
286void
287CaloricallyImperfectGas::p_from_v_e(Real v, Real e, Real & p, Real & dp_dv, Real & dp_de) const
288{
289 p = p_from_v_e(v, e);
290 Real T = T_from_e(e);
291 dp_dv = -_R_specific * T / v / v;
292 dp_de = _R_specific / v / cv_from_T(T);
293}
294
295Real
296CaloricallyImperfectGas::T_from_v_e(Real /*v*/, Real e) const
297{
298 return T_from_e(e);
299}
300
301ADReal
303{
304 if (v.value() == 0.0)
305 return getNaN("Invalid value of specific volume detected (v = " + Moose::stringify(v.value()) +
306 ").");
307
308 Real x = 0;
309 Real raw1 = v.value();
310 Real raw2 = e.value();
311 Real dxd1 = 0;
312 Real dxd2 = 0;
313 T_from_v_e(raw1, raw2, x, dxd1, dxd2);
314
315 ADReal result = x;
316 result.derivatives() = v.derivatives() * dxd1 + e.derivatives() * dxd2;
317 return result;
318}
319
320void
321CaloricallyImperfectGas::T_from_v_e(Real v, Real e, Real & T, Real & dT_dv, Real & dT_de) const
322{
323 T = T_from_e(e);
324 dT_dv = 0.0;
325 dT_de = 1.0 / cv_from_v_e(v, e);
326}
327
328Real
330{
331 Real T = T_from_v_e(v, e);
332
333 const Real c2 = gamma_from_v_e(v, e) * _R_specific * T;
334 if (c2 < 0)
335 return getNaN("Sound speed squared (gamma * R * T) is negative: c2 = " + Moose::stringify(c2) +
336 ".");
337
338 return std::sqrt(c2);
339}
340
341ADReal
343{
344 using std::sqrt;
345 const auto T = T_from_v_e(v, e);
346
347 const auto c2 = gamma_from_v_e(v, e) * _R_specific * T;
348 if (MetaPhysicL::raw_value(c2) < 0)
349 return getNaN("Sound speed squared (gamma * R * T) is negative: c2 = " +
351
352 return sqrt(c2);
353}
354
355void
356CaloricallyImperfectGas::c_from_v_e(Real v, Real e, Real & c, Real & dc_dv, Real & dc_de) const
357{
358 Real T, dT_dv, dT_de;
359 T_from_v_e(v, e, T, dT_dv, dT_de);
360
361 Real gamma, dgamma_dv, dgamma_de;
362 gamma_from_v_e(v, e, gamma, dgamma_dv, dgamma_de);
363 c = std::sqrt(gamma * _R_specific * T);
364
365 dc_dv = 0.5 / c * _R_specific * (gamma * dT_dv + T * dgamma_dv);
366 dc_de = 0.5 / c * _R_specific * (gamma * dT_de + T * dgamma_de);
367}
368
369Real
371{
372 return std::sqrt(gamma_from_p_T(p, T) * _R_specific * T);
373}
374
375ADReal
377{
378 using std::sqrt;
379 return sqrt(gamma_from_p_T(p, T) * _R_specific * T);
380}
381
382void
384 const Real p, const Real T, Real & c, Real & dc_dp, Real & dc_dT) const
385{
386 Real gamma, dgamma_dp, dgamma_dT;
387 gamma_from_p_T(p, T, gamma, dgamma_dp, dgamma_dT);
388 c = std::sqrt(gamma * _R_specific * T);
389 dc_dp = 0;
390 dc_dT = 0.5 / c * _R_specific * (gamma + T * dgamma_dT);
391}
392
393Real
395{
396 Real T = T_from_v_e(v, e);
397 return cp_from_T(T);
398}
399
400void
401CaloricallyImperfectGas::cp_from_v_e(Real v, Real e, Real & cp, Real & dcp_dv, Real & dcp_de) const
402{
403 Real T = T_from_v_e(v, e);
404 Real dcp_dT;
405 cp_from_T(T, cp, dcp_dT);
406 dcp_dv = 0.0;
407 dcp_de = dcp_dT / cv_from_T(T);
408}
409
410Real
412{
413 Real T = T_from_v_e(v, e);
414 return cv_from_T(T);
415}
416
417ADReal
419{
420 Real x = 0;
421 Real raw1 = v.value();
422 Real raw2 = e.value();
423 Real dxd1 = 0;
424 Real dxd2 = 0;
425 cv_from_v_e(raw1, raw2, x, dxd1, dxd2);
426
427 ADReal result = x;
428 result.derivatives() = v.derivatives() * dxd1 + e.derivatives() * dxd2;
429 return result;
430}
431
432void
433CaloricallyImperfectGas::cv_from_v_e(Real v, Real e, Real & cv, Real & dcv_dv, Real & dcv_de) const
434{
435 Real T = T_from_v_e(v, e);
436 Real dcv_dT;
437 cv_from_T(T, cv, dcv_dT);
438 dcv_dv = 0.0;
439 dcv_de = dcv_dT / cv;
440}
441
442Real
444{
445 return cp_from_v_e(v, e) / cv_from_v_e(v, e);
446}
447
448ADReal
450{
451 Real x = 0;
452 Real raw1 = v.value();
453 Real raw2 = e.value();
454 Real dxd1 = 0;
455 Real dxd2 = 0;
456 gamma_from_v_e(raw1, raw2, x, dxd1, dxd2);
457
458 ADReal result = x;
459 result.derivatives() = v.derivatives() * dxd1 + e.derivatives() * dxd2;
460 return result;
461}
462
463void
465 Real v, Real e, Real & gamma, Real & dgamma_dv, Real & dgamma_de) const
466{
467 Real cp, dcp_dv, dcp_de;
468 cp_from_v_e(v, e, cp, dcp_dv, dcp_de);
469 Real cv, dcv_dv, dcv_de;
470 cv_from_v_e(v, e, cv, dcv_dv, dcv_de);
471 gamma = cp / cv;
472 dgamma_dv = 0.0;
473 dgamma_de = 1.0 / cv * dcp_de - gamma / cv * dcv_de;
474}
475
476void
478 Real p, Real T, Real & gamma, Real & dgamma_dp, Real & dgamma_dT) const
479{
480 Real cp, dcp_dp, dcp_dT;
481 cp_from_p_T(p, T, cp, dcp_dp, dcp_dT);
482 Real cv, dcv_dp, dcv_dT;
483 cv_from_p_T(p, T, cv, dcv_dp, dcv_dT);
484 gamma = cp / cv;
485 dgamma_dp = 0.0;
486 dgamma_dT = 1.0 / cv * dcp_dT - gamma / cv * dcv_dT;
487}
488
489Real
491{
492 return cp_from_p_T(p, T) / cv_from_p_T(p, T);
493}
494
495ADReal
497{
498 Real x = 0;
499 Real raw1 = p.value();
500 Real raw2 = T.value();
501 Real dxd1 = 0;
502 Real dxd2 = 0;
503 gamma_from_p_T(raw1, raw2, x, dxd1, dxd2);
504 ADReal result = x;
505 result.derivatives() = p.derivatives() * dxd1 + T.derivatives() * dxd2;
506 return result;
507}
508
509Real
511{
512 const Real T = T_from_v_e(v, e);
513 const Real p = p_from_v_e(v, e);
514 return mu_from_p_T(p, T);
515}
516
517void
518CaloricallyImperfectGas::mu_from_v_e(Real v, Real e, Real & mu, Real & dmu_dv, Real & dmu_de) const
519{
520 const Real T = T_from_v_e(v, e);
521 const Real p = p_from_v_e(v, e);
522 Real dmu_dp, dmu_dT;
523 mu_from_p_T(p, T, mu, dmu_dp, dmu_dT);
524 // only dmu_de = dmu/dT * dT/de = 1/cv * dmu/dT is nonzero
525 // (dmu/dv)_e = 0 because e only depends on T and e is constant
526 dmu_dv = 0.0;
527 dmu_de = dmu_dT / cv_from_p_T(p, T);
528}
529
530Real
532{
533 const Real T = T_from_v_e(v, e);
534 const Real p = p_from_v_e(v, e);
535 return k_from_p_T(p, T);
536}
537
538void
539CaloricallyImperfectGas::k_from_v_e(Real v, Real e, Real & k, Real & dk_dv, Real & dk_de) const
540{
541 const Real T = T_from_v_e(v, e);
542 const Real p = p_from_v_e(v, e);
543 Real dk_dp, dk_dT;
544 k_from_p_T(p, T, k, dk_dp, dk_dT);
545 // only dk_de = dk/dT * dT/de = 1/cv * dk/dT is nonzero
546 // (dk/dv)_e = 0 because e only depends on T and e is constant
547 dk_dv = 0.0;
548 dk_de = dk_dT / cv_from_p_T(p, T);
549}
550
551Real
553{
554 Real T = T_from_v_e(v, e);
555 Real Z = Z_from_T(T);
556 return Z + _R_specific * std::log(v);
557}
558
559void
560CaloricallyImperfectGas::s_from_v_e(Real v, Real e, Real & s, Real & ds_dv, Real & ds_de) const
561{
562 s = s_from_v_e(v, e);
563 // see documentation for derivation
564 ds_dv = _R_specific / v;
565 ds_de = 1.0 / T_from_v_e(v, e);
566}
567
568Real
570{
571 Real Z = Z_from_T(T);
572 Real v = 1.0 / rho_from_p_T(p, T);
573 return Z + _R_specific * std::log(v);
574}
575
576void
577CaloricallyImperfectGas::s_from_p_T(Real p, Real T, Real & s, Real & ds_dp, Real & ds_dT) const
578{
579 s = s_from_p_T(p, T);
580 ds_dp = -_R_specific / p;
581 ds_dT = cp_from_p_T(p, T) / T;
582}
583
584Real
586{
587 Real T = T_from_p_h(p, h);
588 Real v = 1.0 / rho_from_p_T(p, T);
589 Real Z = Z_from_T(T);
590 return Z + _R_specific * std::log(v);
591}
592
593void
594CaloricallyImperfectGas::s_from_h_p(Real h, Real p, Real & s, Real & ds_dh, Real & ds_dp) const
595{
596 s = s_from_h_p(h, p);
597 ds_dh = 1.0 / T_from_p_h(p, h);
598 ds_dp = -_R_specific / p;
599}
600
601Real
603{
604 auto s_diff = [&p, &s, this](Real x) { return this->s_from_p_T(p, x) - s; };
605 Real min = _min_temperature;
606 Real max = _max_temperature;
607 BrentsMethod::bracket(s_diff, min, max);
608 Real T = BrentsMethod::root(s_diff, min, max);
609 return rho_from_p_T(p, T);
610}
611
612Real
613CaloricallyImperfectGas::e_from_v_h(Real /*v*/, Real h) const
614{
615 Real T = T_from_h(h);
616 return e_from_T(T);
617}
618
619void
620CaloricallyImperfectGas::e_from_v_h(Real v, Real h, Real & e, Real & de_dv, Real & de_dh) const
621{
622 e = e_from_v_h(v, h);
623 Real cv = cv_from_v_e(v, e);
624 Real cp = cp_from_v_e(v, e);
625 de_dv = 0.0;
626 de_dh = cv / cp;
627}
628
629Real
631{
632 return p * _molar_mass / (_R * T);
633}
634
635ADReal
637{
638 return p * _molar_mass / (_R * T);
639}
640
641void
643 const ADReal & p, const ADReal & T, ADReal & rho, ADReal & drho_dp, ADReal & drho_dT) const
644{
645 rho = rho_from_p_T(p, T);
646 drho_dp = _molar_mass / (_R * T);
647 drho_dT = -p * _molar_mass / (_R * T * T);
648}
649
650void
652 Real p, Real T, Real & rho, Real & drho_dp, Real & drho_dT) const
653{
654 rho = rho_from_p_T(p, T);
655 drho_dp = _molar_mass / (_R * T);
656 drho_dT = -p * _molar_mass / (_R * T * T);
657}
658
659Real
660CaloricallyImperfectGas::e_from_p_rho(Real p, Real rho) const
661{
662 return e_from_p_rho_template(p, rho);
663}
664
665ADReal
666CaloricallyImperfectGas::e_from_p_rho(const ADReal & p, const ADReal & rho) const
667{
668 return e_from_p_rho_template(p, rho);
669}
670
671void
672CaloricallyImperfectGas::e_from_p_rho(
673 Real p, Real rho, Real & e, Real & de_dp, Real & de_drho) const
674{
675 e_from_p_rho_template(p, rho, e, de_dp, de_drho);
676}
677
678void
679CaloricallyImperfectGas::e_from_p_rho(
680 const ADReal & p, const ADReal & rho, ADReal & e, ADReal & de_dp, ADReal & de_drho) const
681{
682 e_from_p_rho_template(p, rho, e, de_dp, de_drho);
683}
684
685Real
687{
688 return e_from_T(T);
689}
690
691void
692CaloricallyImperfectGas::e_from_T_v(Real T, Real v, Real & e, Real & de_dT, Real & de_dv) const
693{
694 e = e_from_T_v(T, v);
695 de_dT = cv_from_T_v(T, v);
696 de_dv = 0.0;
697}
698
699ADReal
701{
702 Real x = 0;
703 Real raw1 = T.value();
704 Real raw2 = v.value();
705 Real dxd1 = 0;
706 Real dxd2 = 0;
707 e_from_T_v(raw1, raw2, x, dxd1, dxd2);
708
709 ADReal result = x;
710 result.derivatives() = T.derivatives() * dxd1 + v.derivatives() * dxd2;
711 return result;
712}
713
714Real
716{
717 return _R_specific * T / v;
718}
719
720void
721CaloricallyImperfectGas::p_from_T_v(Real T, Real v, Real & p, Real & dp_dT, Real & dp_dv) const
722{
723 p = _R_specific * T / v;
724 dp_dT = _R_specific / v;
725 dp_dv = -_R_specific * T / (v * v);
726}
727
728Real
730{
731 return h_from_T(T);
732}
733
734void
735CaloricallyImperfectGas::h_from_T_v(Real T, Real /*v*/, Real & h, Real & dh_dT, Real & dh_dv) const
736{
737 h = h_from_T(T);
738 dh_dT = cp_from_T(T);
739 dh_dv = 0.0;
740}
741
742Real
744{
745 Real p = p_from_T_v(T, v);
746 return s_from_p_T(p, T);
747}
748
749void
750CaloricallyImperfectGas::s_from_T_v(Real T, Real v, Real & s, Real & ds_dT, Real & ds_dv) const
751{
752 Real p, dp_dT_v, dp_dv_T;
753 Real ds_dp_T, ds_dT_p;
754 p_from_T_v(T, v, p, dp_dT_v, dp_dv_T);
755 s_from_p_T(p, T, s, ds_dp_T, ds_dT_p);
756 ds_dT = ds_dT_p + ds_dp_T * dp_dT_v;
757 ds_dv = ds_dp_T * dp_dv_T;
758}
759
760Real
762{
763 return cv_from_T(T);
764}
765
766Real
768{
769 return _e_c;
770}
771
772void
773CaloricallyImperfectGas::v_e_spndl_from_T(Real /*T*/, Real & v, Real & e) const
774{
775 v = 1. / _rho_c;
776 e = _e_c;
777}
778
779Real
781{
782 return h_from_T(T);
783}
784
785void
786CaloricallyImperfectGas::h_from_p_T(Real p, Real T, Real & h, Real & dh_dp, Real & dh_dT) const
787{
788 h = h_from_p_T(p, T);
789 dh_dp = 0.0;
790 dh_dT = cp_from_p_T(p, T);
791}
792
793Real
795{
796 return e_from_T(T);
797}
798
799ADReal
801{
802 Real x = 0;
803 Real raw1 = p.value();
804 Real raw2 = T.value();
805 Real dxd1 = 0;
806 Real dxd2 = 0;
807 e_from_p_T(raw1, raw2, x, dxd1, dxd2);
808
809 ADReal result = x;
810 result.derivatives() = p.derivatives() * dxd1 + T.derivatives() * dxd2;
811 return result;
812}
813
814void
815CaloricallyImperfectGas::e_from_p_T(Real p, Real T, Real & e, Real & de_dp, Real & de_dT) const
816{
817 e = e_from_p_T(p, T);
818 de_dp = 0.0;
819 de_dT = cv_from_p_T(p, T);
820}
821
822Real
827
828Real
833
834Real
839
840Real
845
846Real
847CaloricallyImperfectGas::T_from_p_h(Real /*p*/, Real h) const
848{
849 return T_from_h(h);
850}
851
852void
853CaloricallyImperfectGas::T_from_p_h(Real p, Real h, Real & T, Real & dT_dp, Real & dT_dh) const
854{
855 T = T_from_p_h(p, h);
856 dT_dp = 0;
857 dT_dh = 1.0 / cp_from_p_T(p, T);
858}
859
860Real
861CaloricallyImperfectGas::cv_from_p_T(Real /* pressure */, Real temperature) const
862{
863 return cv_from_T(temperature);
864}
865
866void
868 Real /*p*/, Real T, Real & cv, Real & dcv_dp, Real & dcv_dT) const
869{
870 cv_from_T(T, cv, dcv_dT);
871 dcv_dp = 0.0;
872}
873
874Real
875CaloricallyImperfectGas::cp_from_p_T(Real /* pressure */, Real temperature) const
876{
877 return cp_from_T(temperature);
878}
879
880void
882 Real /*p*/, Real T, Real & cp, Real & dcp_dp, Real & dcp_dT) const
883{
884 cp_from_T(T, cp, dcp_dT);
885 dcp_dp = 0.0;
886}
887
888Real
889CaloricallyImperfectGas::mu_from_p_T(Real /* pressure */, Real temperature) const
890{
891 return _mu_T->value(temperature, Point());
892}
893
894void
895CaloricallyImperfectGas::mu_from_p_T(Real p, Real T, Real & mu, Real & dmu_dp, Real & dmu_dT) const
896{
897 mu = this->mu_from_p_T(p, T);
898 dmu_dp = 0.0;
899 const Real dT = temperatureStep(T);
900 dmu_dT = (this->mu_from_p_T(p, T + dT) - mu) / dT;
901}
902
903Real
904CaloricallyImperfectGas::k_from_p_T(Real /* pressure */, Real temperature) const
905{
906 return _k_T->value(temperature, Point());
907}
908
909void
910CaloricallyImperfectGas::k_from_p_T(Real p, Real T, Real & k, Real & dk_dp, Real & dk_dT) const
911{
912 k = k_from_p_T(p, T);
913 dk_dp = 0.0;
914 const Real dT = temperatureStep(T);
915 dk_dT = (this->k_from_p_T(p, T + dT) - k) / dT;
916}
917
918Real
920{
921 // use the definition of the Gibbs free energy
922 Real T = T_from_v_e(v, e);
923 Real p = p_from_v_e(v, e);
924 return h_from_p_T(p, T) - T * s_from_p_T(p, T);
925}
926
927Real
929{
930 Real T = T_from_h(h);
931 Real e = e_from_T(T);
932 Real Z = Z_from_T(T);
933 Real v = std::exp((s - Z) / _R_specific);
934 return p_from_v_e(v, e);
935}
936
937void
938CaloricallyImperfectGas::p_from_h_s(Real h, Real s, Real & p, Real & dp_dh, Real & dp_ds) const
939{
940 p = p_from_h_s(h, s);
941 // From dh = T ds + v dp
942 const Real T = T_from_h(h);
943 const Real rho = rho_from_p_T(p, T);
944 dp_dh = rho;
945 dp_ds = -T * rho;
946}
947
948Real
950{
951 // Relative step for one-sided finite differences, taken backward near the upper bound of
952 // the tabulated temperature range so that the perturbed state stays inside it
953 const Real dT = 1.0e-7 * T;
954 return T + dT > _max_temperature ? -dT : dT;
955}
956
957void
958CaloricallyImperfectGas::outOfBounds(const std::string & function,
959 Real value,
960 Real min,
961 Real max) const
962{
963 std::stringstream ss;
964 ss << "Function " << function << " encountered argument value of " << value
965 << " which is outside of the bounds of " << min << " .. " << max;
967 mooseError(ss.str());
968 else
969 mooseWarning(ss.str());
970}
DualNumber< Real, DNDerivativeType, true > ADReal
registerMooseObject("FluidPropertiesApp", CaloricallyImperfectGas)
const std::vector< double > x
const double mu
const Real p
const double rho
const double T
const double v
A calorically imperfect gas fluid property class This fluid property assumes that internal energy is ...
Real _min_temperature
temperature limits when creating lookup tables
const Function * _k_T
Thermal conductivity.
virtual Real e_from_v_h(Real v, Real h) const override
virtual Real p_from_v_e(Real v, Real e) const override
virtual Real s_from_v_e(Real v, Real e) const override
virtual Real c_from_p_T(Real p, Real T) const override
virtual Real cv_from_v_e(Real v, Real e) const override
virtual Real gamma_from_p_T(Real p, Real T) const override
virtual Real rho_from_p_s(Real p, Real s) const override
const bool _out_of_bound_error
Flag to error if out of bounds.
virtual void initialSetup() override
virtual Real mu_from_p_T(Real p, Real T) const override
virtual Real criticalDensity() const override
Critical density.
virtual Real criticalInternalEnergy() const override
Critical specific internal energy.
virtual Real cv_from_T_v(Real T, Real v) const override
virtual void v_e_spndl_from_T(Real T, Real &v, Real &e) const override
Specific internal energy from temperature and specific volume.
virtual Real cp_from_v_e(Real v, Real e) const override
virtual Real e_spndl_from_v(Real v) const override
Specific internal energy from temperature and specific volume.
virtual Real k_from_p_T(Real pressure, Real temperature) const override
virtual Real e_from_p_T(Real p, Real T) const override
const Function * _mu_T
Dynamic viscosity.
Real e_from_T(Real T) const
helper functions for e(T) and h(T), T(e), T(h), cv(T), cp(T)
std::vector< Real > _T_h_lookup
inverse lookup table data
Real _delta_T
temperature interval in lookup tables
virtual Real s_from_T_v(Real T, Real v) const override
virtual Real h_from_p_T(Real p, Real T) const override
virtual Real e_from_T_v(Real T, Real v) const override
virtual Real h_from_T_v(Real T, Real v) const override
virtual Real T_from_v_e(Real v, Real e) const override
Real cv_from_T(const CppType &T) const
virtual Real g_from_v_e(Real v, Real e) const override
virtual Real s_from_h_p(Real h, Real p) const override
const Real & _molar_mass
molar mass
static InputParameters validParams()
Real _min_e
internal energy and enthalpy limits when creating lookup tables
CaloricallyImperfectGas(const InputParameters &parameters)
virtual Real cv_from_p_T(Real p, Real T) const override
virtual std::string fluidName() const override
Fluid name.
virtual Real T_from_p_h(Real p, Real h) const override
Real _delta_e
step size in internal energy and enthalpy
virtual Real c_from_v_e(Real v, Real e) const override
void outOfBounds(const std::string &function, Real value, Real min, Real max) const
function that handles exceeding parameter limits
virtual Real k_from_v_e(Real v, Real e) const override
virtual Real p_from_T_v(Real T, Real v) const override
virtual Real cp_from_p_T(Real p, Real T) const override
virtual Real criticalTemperature() const override
Critical temperature.
void setupLookupTables()
sets up the T(e) reverse lookup table
virtual Real rho_from_p_T(Real p, Real T) const override
virtual Real gamma_from_v_e(Real v, Real e) const override
const Function * _e_T
Internal energy as a function of temperature.
Real temperatureStep(Real T) const
Signed temperature step for one-sided finite differences that stays within the table.
virtual Real molarMass() const override
Molar mass [kg/mol].
virtual Real s_from_p_T(Real p, Real T) const override
std::vector< Real > _Z_T_lookup
Z(T) lookup table on uniform grid between _min_temperature and _max_temperature.
virtual Real mu_from_v_e(Real v, Real e) const override
virtual Real p_from_h_s(Real h, Real s) const override
const Real _R_specific
Specific gas constant (R / molar mass)
static const Real _R
Universal gas constant (J/mol/K)
const Function & getFunction(const std::string &name) const
virtual Real value(Real t, const Point &p) const
virtual Real timeDerivative(Real t, const Point &p) const
void addRequiredParam(const std::string &name, const std::string &doc_string)
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
void addClassDescription(const std::string &doc_string)
void mooseError(Args &&... args) const
void mooseWarning(Args &&... args) const
Interface class for producing errors, warnings, or just quiet NaNs.
Real getNaN() const
Throws an error or returns a NaN with or without a warning, with a default message.
static InputParameters validParams()
virtual void initialSetup()
Common class for single phase fluid properties.
static InputParameters validParams()
e e e e s T T T T T rho v v T e h
Real root(std::function< Real(Real)> const &f, Real x1, Real x2, Real tol=1.0e-12)
Finds the root of a function using Brent's method.
void bracket(std::function< Real(Real)> const &f, Real &x1, Real &x2)
Function to bracket a root of a given function.
auto raw_value(const Eigen::Map< T > &in)
std::string stringify(const T &t)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real