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