https://mooseframework.inl.gov
Loading...
Searching...
No Matches
IdealRealGasMixtureFluidProperties.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
12#include "BrentsMethod.h"
13#include <numeric>
14
16
21#define define_from_v_e_using_T_v(prop) \
22 Real IdealRealGasMixtureFluidProperties::prop##_from_v_e( \
23 Real v, Real e, const std::vector<Real> & x) const \
24 { \
25 const Real T = T_from_v_e(v, e, x); \
26 return prop##_from_T_v(T, v, x); \
27 } \
28 \
29 void IdealRealGasMixtureFluidProperties::prop##_from_v_e(Real v, \
30 Real e, \
31 const std::vector<Real> & x, \
32 Real & y, \
33 Real & dy_dv, \
34 Real & dy_de, \
35 std::vector<Real> & dy_dx) const \
36 { \
37 Real T, dT_dv_e, dT_de_v; \
38 std::vector<Real> dT_dx_ve(_n_secondary_vapors); \
39 T_from_v_e(v, e, x, T, dT_dv_e, dT_de_v, dT_dx_ve); \
40 \
41 Real dy_dT_v, dy_dv_T; \
42 std::vector<Real> dy_dx_Tv; \
43 prop##_from_T_v(T, v, x, y, dy_dT_v, dy_dv_T, dy_dx_Tv); \
44 \
45 dy_dv = dy_dv_T + dy_dT_v * dT_dv_e; \
46 dy_de = dy_dT_v * dT_de_v; \
47 dy_dx.resize(_n_secondary_vapors); \
48 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
49 dy_dx[i] = dy_dx_Tv[i] + dy_dT_v * dT_dx_ve[i]; \
50 }
51
56#define define_from_p_T_using_T_v(prop) \
57 Real IdealRealGasMixtureFluidProperties::prop##_from_p_T( \
58 Real p, Real T, const std::vector<Real> & x) const \
59 { \
60 const Real v = v_from_p_T(p, T, x); \
61 return prop##_from_T_v(T, v, x); \
62 } \
63 \
64 void IdealRealGasMixtureFluidProperties::prop##_from_p_T(Real p, \
65 Real T, \
66 const std::vector<Real> & x, \
67 Real & y, \
68 Real & dy_dp, \
69 Real & dy_dT, \
70 std::vector<Real> & dy_dx) const \
71 { \
72 Real v, dv_dp_T, dv_dT_p; \
73 std::vector<Real> dv_dx_pT; \
74 v_from_p_T(p, T, x, v, dv_dp_T, dv_dT_p, dv_dx_pT); \
75 \
76 Real dy_dT_v, dy_dv_T; \
77 std::vector<Real> dy_dx_Tv; \
78 prop##_from_T_v(T, v, x, y, dy_dT_v, dy_dv_T, dy_dx_Tv); \
79 \
80 dy_dp = dy_dv_T * dv_dp_T; \
81 dy_dT = dy_dT_v + dy_dv_T * dv_dT_p; \
82 dy_dx.resize(_n_secondary_vapors); \
83 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
84 dy_dx[i] = dy_dx_Tv[i] + dy_dv_T * dv_dx_pT[i]; \
85 }
86
90#define define_mass_specific_prop_from_T_v(prop) \
91 Real IdealRealGasMixtureFluidProperties::prop##_from_T_v( \
92 Real T, Real v, const std::vector<Real> & x) const \
93 { \
94 const Real x_primary = primaryMassFraction(x); \
95 Real y = x_primary * _fp_primary->prop##_from_T_v(T, v / x_primary); \
96 \
97 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
98 y += x[i] * _fp_secondary[i]->prop##_from_T_v(T, v / x[i]); \
99 \
100 return y; \
101 } \
102 \
103 void IdealRealGasMixtureFluidProperties::prop##_from_T_v(Real T, \
104 Real v, \
105 const std::vector<Real> & x, \
106 Real & y, \
107 Real & dy_dT, \
108 Real & dy_dv, \
109 std::vector<Real> & dy_dx) const \
110 { \
111 const Real x_primary = primaryMassFraction(x); \
112 mooseAssert(!MooseUtils::absoluteFuzzyEqual(x_primary, 0.0), "Mass fraction may not be zero"); \
113 \
114 Real y_primary, dy_dT_primary, dy_dv_primary; \
115 _fp_primary->prop##_from_T_v(T, v / x_primary, y_primary, dy_dT_primary, dy_dv_primary); \
116 y = x_primary * y_primary; \
117 dy_dT = x_primary * dy_dT_primary; \
118 dy_dv = dy_dv_primary; \
119 \
120 Real dy_dT_sec; \
121 std::vector<Real> y_sec(_n_secondary_vapors), dy_dv_sec(_n_secondary_vapors); \
122 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
123 { \
124 mooseAssert(!MooseUtils::absoluteFuzzyEqual(x[i], 0.0), "Mass fraction may not be zero"); \
125 _fp_secondary[i]->prop##_from_T_v(T, v / x[i], y_sec[i], dy_dT_sec, dy_dv_sec[i]); \
126 y += x[i] * y_sec[i]; \
127 dy_dT += x[i] * dy_dT_sec; \
128 dy_dv += dy_dv_sec[i]; \
129 } \
130 \
131 dy_dx.resize(_n_secondary_vapors); \
132 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
133 { \
134 dy_dx[i] = y_sec[i] - x[i] * dy_dv_sec[i] * v / (x[i] * x[i]); \
135 for (unsigned int j = 0; j < _n_secondary_vapors; j++) \
136 { \
137 if (j == i) \
138 continue; \
139 if (_n_secondary_vapors > 1) \
140 imperfectJacobianMessage( \
141 "The mass fraction derivatives in the following " \
142 "function have not been tested for mixtures of 3 or more components:\n\n", \
143 __PRETTY_FUNCTION__); \
144 const Real dxj_dxi = 0; \
145 dy_dx[i] += dxj_dxi * y_sec[j] - x[j] * dy_dv_sec[j] * v / (x[j] * x[j]) * dxj_dxi; \
146 } \
147 const Real dx_primary_dxi = -x_primary / (1. - x[i]); \
148 dy_dx[i] += dx_primary_dxi * y_primary - \
149 x_primary * dy_dv_primary * v / (x_primary * x_primary) * dx_primary_dxi; \
150 } \
151 }
152
156#define define_transport_prop_from_T_v(prop) \
157 Real IdealRealGasMixtureFluidProperties::prop##_from_T_v( \
158 Real T, Real v, const std::vector<Real> & x) const \
159 { \
160 const Real x_primary = primaryMassFraction(x); \
161 Real M_primary = _fp_primary->molarMass(); \
162 \
163 Real sum = x_primary / M_primary; \
164 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
165 sum += x[i] / _fp_secondary[i]->molarMass(); \
166 const Real M_star = 1. / sum; \
167 \
168 const Real vp = v / x_primary; \
169 const Real ep = _fp_primary->e_from_T_v(T, vp); \
170 const Real yp = _fp_primary->prop##_from_v_e(vp, ep); \
171 Real y = x_primary * M_star / M_primary * yp; \
172 \
173 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
174 { \
175 const Real vi = v / x[i]; \
176 const Real ei = _fp_secondary[i]->e_from_T_v(T, vp); \
177 const Real Mi = _fp_secondary[i]->molarMass(); \
178 const Real yi = _fp_secondary[i]->prop##_from_v_e(vi, ei); \
179 y += x[i] * M_star / Mi * yi; \
180 } \
181 \
182 return y; \
183 } \
184 \
185 void IdealRealGasMixtureFluidProperties::prop##_from_T_v(Real T, \
186 Real v, \
187 const std::vector<Real> & x, \
188 Real & y, \
189 Real & dy_dT, \
190 Real & dy_dv, \
191 std::vector<Real> & dy_dx) const \
192 { \
193 const Real x_primary = primaryMassFraction(x); \
194 Real M_primary = _fp_primary->molarMass(); \
195 \
196 Real sum = x_primary / M_primary; \
197 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
198 sum += x[i] / _fp_secondary[i]->molarMass(); \
199 const Real M_star = 1. / sum; \
200 \
201 const Real vp = v / x_primary; \
202 const Real ep = _fp_primary->e_from_T_v(T, vp); \
203 const Real yp = _fp_primary->prop##_from_v_e(vp, ep); \
204 y = x_primary * M_star / M_primary * yp; \
205 \
206 imperfectJacobianMessage("The temperature and specific volume derivatives in the following " \
207 "function are currently neglected:\n\n", \
208 __PRETTY_FUNCTION__); \
209 dy_dT = 0; \
210 dy_dv = 0; \
211 \
212 Real sum_yj = 0; \
213 dy_dx.resize(_n_secondary_vapors); \
214 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
215 { \
216 const Real vi = v / x[i]; \
217 const Real ei = _fp_secondary[i]->e_from_T_v(T, vi); \
218 const Real Mi = _fp_secondary[i]->molarMass(); \
219 const Real yi = _fp_secondary[i]->prop##_from_v_e(vi, ei); \
220 y += x[i] * M_star / Mi * yi; \
221 dy_dx[i] = -M_star / M_primary * yp - \
222 x_primary * (1 / Mi - 1 / M_primary) * M_star * M_star / M_primary * yp + \
223 M_star / Mi * yi; \
224 sum_yj += x[i] * yi / Mi; \
225 } \
226 \
227 for (unsigned int i = 0; i < _n_secondary_vapors; i++) \
228 { \
229 const Real Mi = _fp_secondary[i]->molarMass(); \
230 dy_dx[i] += -M_star * M_star * (1 / Mi - 1 / M_primary) * sum_yj; \
231 } \
232 }
233
234// clang-format off
237
238define_transport_prop_from_T_v(mu)
239define_transport_prop_from_T_v(k)
240
241define_from_p_T_using_T_v(e)
242define_from_p_T_using_T_v(s)
243define_from_p_T_using_T_v(c)
244define_from_p_T_using_T_v(cp)
245define_from_p_T_using_T_v(cv)
246define_from_p_T_using_T_v(mu)
247define_from_p_T_using_T_v(k)
248
249define_from_v_e_using_T_v(p)
250define_from_v_e_using_T_v(c)
251
253// clang-format on
254{
256 params += NaNInterface::validParams();
257
258 params.addClassDescription("Class for fluid properties of an arbitrary vapor mixture");
259
260 params.addRequiredParam<UserObjectName>(
261 "fp_primary", "Name of fluid properties user object for primary vapor component");
262 params.addRequiredParam<std::vector<UserObjectName>>(
263 "fp_secondary", "Name of fluid properties user object(s) for secondary vapor component(s)");
264 params.addParam<Real>("_T_mix_max", 1300., "Maximum temperature of the mixture");
265
266 return params;
267}
268
270 const InputParameters & parameters)
271 : VaporMixtureFluidProperties(parameters),
272 NaNInterface(this),
273 _fp_primary(&getUserObject<SinglePhaseFluidProperties>("fp_primary")),
274 _fp_secondary_names(getParam<std::vector<UserObjectName>>("fp_secondary")),
275 _n_secondary_vapors(_fp_secondary_names.size()),
276 _T_mix_max(getParam<Real>("_T_mix_max"))
277{
279 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
280 _fp_secondary[i] = &getUserObjectByName<SinglePhaseFluidProperties>(_fp_secondary_names[i]);
281}
282
288
291{
292 mooseAssert(i < getNumberOfSecondaryVapors(), "Requested secondary index too high.");
293 return *_fp_secondary[i];
294}
295
296Real
297IdealRealGasMixtureFluidProperties::T_from_v_e(Real v, Real e, const std::vector<Real> & x) const
298{
299 Real v_primary = v / primaryMassFraction(x);
300 static const Real vc = 1. / _fp_primary->criticalDensity();
301 static const Real ec = _fp_primary->criticalInternalEnergy();
302
303 // Initial estimate of a bracketing interval for the temperature
304 Real lower_temperature, upper_temperature;
305 if (v_primary > vc)
306 {
307 Real e_sat_primary = _fp_primary->e_spndl_from_v(v_primary);
308 lower_temperature = _fp_primary->T_from_v_e(v_primary, e_sat_primary);
309 }
310 else
311 lower_temperature = _fp_primary->T_from_v_e(v_primary, ec);
312
313 upper_temperature = _T_mix_max;
314
315 // Use BrentsMethod to find temperature
316 auto energy_diff = [&v, &e, &x, this](Real T) { return this->e_from_T_v(T, v, x) - e; };
317
318 BrentsMethod::bracket(energy_diff, lower_temperature, upper_temperature);
319 return BrentsMethod::root(energy_diff, lower_temperature, upper_temperature);
320}
321
322void
323IdealRealGasMixtureFluidProperties::T_from_v_e(Real v,
324 Real e,
325 const std::vector<Real> & x,
326 Real & T,
327 Real & dT_dv,
328 Real & dT_de,
329 std::vector<Real> & dT_dx) const
330{
331 T = T_from_v_e(v, e, x);
332
333 Real e_unused, de_dT_v, de_dv_T;
334 std::vector<Real> de_dx_Tv;
335 e_from_T_v(T, v, x, e_unused, de_dT_v, de_dv_T, de_dx_Tv);
336
337 dT_dv = -de_dv_T / de_dT_v;
338 dT_de = 1.0 / de_dT_v;
339 dT_dx.resize(_n_secondary_vapors);
340 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
341 dT_dx[i] = -de_dx_Tv[i] / de_dT_v;
342}
343
344Real
345IdealRealGasMixtureFluidProperties::rho_from_p_T(Real p, Real T, const std::vector<Real> & x) const
346{
347 return 1.0 / v_from_p_T(p, T, x);
348}
349
350void
351IdealRealGasMixtureFluidProperties::rho_from_p_T(Real p,
352 Real T,
353 const std::vector<Real> & x,
354 Real & rho,
355 Real & drho_dp,
356 Real & drho_dT,
357 std::vector<Real> & drho_dx) const
358{
359 Real v, dv_dp, dv_dT;
360 std::vector<Real> dv_dx;
361 v_from_p_T(p, T, x, v, dv_dp, dv_dT, dv_dx);
362
363 rho = 1.0 / v;
364 const Real drho_dv = -1.0 / (v * v);
365 drho_dp = drho_dv * dv_dp;
366 drho_dT = drho_dv * dv_dT;
367 drho_dx.resize(_n_secondary_vapors);
368 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
369 drho_dx[i] = drho_dv * dv_dx[i];
370}
371
372Real
373IdealRealGasMixtureFluidProperties::v_from_p_T(Real p, Real T, const std::vector<Real> & x) const
374{
375 const Real x_primary = primaryMassFraction(x);
376 Real M_primary = _fp_primary->molarMass();
377
378 Real sum = x_primary / M_primary;
379 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
380 sum += x[i] / _fp_secondary[i]->molarMass();
381 Real M_star = 1. / sum;
382 Real v_ideal = R_molar * T / (M_star * p);
383
384 // check range of validity for primary (condensable) component
385 static const Real Tc = _fp_primary->criticalTemperature();
386 static const Real vc = 1. / _fp_primary->criticalDensity();
387 Real v_spndl, e_spndl;
388 if (T < Tc)
389 _fp_primary->v_e_spndl_from_T(T, v_spndl, e_spndl);
390 else
391 v_spndl = vc;
392
393 Real lower_spec_volume = v_spndl * x_primary;
394 Real upper_spec_volume = v_ideal; // p*v/(RT) <= 1
395
396 // Initial estimate of a bracketing interval for the temperature
397 Real p_max = p_from_T_v(T, lower_spec_volume, x);
398 if (p > p_max || upper_spec_volume < lower_spec_volume)
399 return getNaN();
400
401 // Use BrentsMethod to find temperature
402 auto pressure_diff = [&T, &p, &x, this](Real v) { return this->p_from_T_v(T, v, x) - p; };
403
404 BrentsMethod::bracket(pressure_diff, lower_spec_volume, upper_spec_volume);
405 return BrentsMethod::root(pressure_diff, lower_spec_volume, upper_spec_volume);
406}
407
408void
409IdealRealGasMixtureFluidProperties::v_from_p_T(Real p,
410 Real T,
411 const std::vector<Real> & x,
412 Real & v,
413 Real & dv_dp,
414 Real & dv_dT,
415 std::vector<Real> & dv_dx) const
416{
417 v = v_from_p_T(p, T, x);
418
419 Real p_unused, dp_dT, dp_dv;
420 std::vector<Real> dp_dx;
421 p_from_T_v(T, v, x, p_unused, dp_dT, dp_dv, dp_dx);
422
423 dv_dp = 1. / dp_dv;
424 dv_dT = -dp_dT / dp_dv;
425 dv_dx.resize(_n_secondary_vapors);
426 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
427 dv_dx[i] = -dp_dx[i] / dp_dv;
428}
429
430Real
431IdealRealGasMixtureFluidProperties::e_from_p_rho(Real p,
432 Real rho,
433 const std::vector<Real> & x) const
434{
435 Real v = 1. / rho;
436 Real T = T_from_p_v(p, v, x);
437
438 return e_from_T_v(T, v, x);
439}
440
441void
442IdealRealGasMixtureFluidProperties::e_from_p_rho(Real p,
443 Real rho,
444 const std::vector<Real> & x,
445 Real & e,
446 Real & de_dp,
447 Real & de_drho,
448 std::vector<Real> & de_dx) const
449{
450 Real v = 1. / rho;
451 Real T = T_from_p_v(p, v, x);
452 Real p_, dp_dT, dp_dv, de_dT, de_dv;
453 std::vector<Real> dp_dx_Tv;
454 std::vector<Real> de_dx_Tv;
455
456 p_from_T_v(T, v, x, p_, dp_dT, dp_dv, dp_dx_Tv);
457 e_from_T_v(T, v, x, e, de_dT, de_dv, de_dx_Tv);
458
459 de_dp = de_dT / dp_dT;
460 de_drho = (-v * v) * (de_dv - de_dT * dp_dv / dp_dT);
461
462 // Derivatives with respect to mass fractions:
463 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
464 {
465 de_dx[i] = de_dx_Tv[i] - de_dT * dp_dx_Tv[i] / dp_dT;
466 }
467}
468
469Real
470IdealRealGasMixtureFluidProperties::T_from_p_v(Real p, Real v, const std::vector<Real> & x) const
471{
472 Real v_primary = v / primaryMassFraction(x);
473 static const Real vc = 1. / _fp_primary->criticalDensity();
474 static const Real ec = _fp_primary->criticalInternalEnergy();
475
476 // Initial estimate of a bracketing interval for the temperature
477 Real lower_temperature, upper_temperature;
478 if (v_primary > vc)
479 {
480 Real e_sat_primary = _fp_primary->e_spndl_from_v(v_primary);
481 lower_temperature = _fp_primary->T_from_v_e(v_primary, e_sat_primary);
482 }
483 else
484 lower_temperature = _fp_primary->T_from_v_e(v_primary, ec);
485
486 upper_temperature = _T_mix_max;
487
488 // Use BrentsMethod to find temperature
489 auto pressure_diff = [&p, &v, &x, this](Real T) { return this->p_from_T_v(T, v, x) - p; };
490
491 BrentsMethod::bracket(pressure_diff, lower_temperature, upper_temperature);
492 return BrentsMethod::root(pressure_diff, lower_temperature, upper_temperature);
493}
494
495void
496IdealRealGasMixtureFluidProperties::T_from_p_v(Real p,
497 Real v,
498 const std::vector<Real> & x,
499 Real & T,
500 Real & dT_dp,
501 Real & dT_dv,
502 std::vector<Real> & dT_dx) const
503{
504 T = T_from_p_v(p, v, x);
505
506 // pressure and derivatives
507 Real p_unused, dp_dT_v, dp_dv_T;
508 std::vector<Real> dp_dx_Tv;
509 p_from_T_v(T, v, x, p_unused, dp_dT_v, dp_dv_T, dp_dx_Tv);
510
511 // Compute derivatives using the following rules:
512 dT_dp = 1. / dp_dT_v;
513 dT_dv = -dp_dv_T / dp_dT_v;
514
515 // Derivatives with respect to mass fractions:
516 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
517 dT_dx[i] = -dp_dx_Tv[i] / dp_dT_v;
518}
519
520Real
521IdealRealGasMixtureFluidProperties::p_from_T_v(Real T, Real v, const std::vector<Real> & x) const
522{
523 const Real x_primary = primaryMassFraction(x);
524 Real p = _fp_primary->p_from_T_v(T, v / x_primary);
525
526 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
527 p += _fp_secondary[i]->p_from_T_v(T, v / x[i]);
528
529 return p;
530}
531
532void
533IdealRealGasMixtureFluidProperties::p_from_T_v(Real T,
534 Real v,
535 const std::vector<Real> & x,
536 Real & p,
537 Real & dp_dT,
538 Real & dp_dv,
539 std::vector<Real> & dp_dx) const
540{
541 const Real x_primary = primaryMassFraction(x);
542
543 Real p_primary, dp_dT_primary, dp_dv_primary;
544 _fp_primary->p_from_T_v(T, v / x_primary, p_primary, dp_dT_primary, dp_dv_primary);
545 p = p_primary;
546 dp_dT = dp_dT_primary;
547 dp_dv = dp_dv_primary / x_primary;
548
549 // get the partial pressures and their derivatives first
550 Real p_sec, dp_dT_sec, dxj_dxi;
551 std::vector<Real> dp_dv_sec(_n_secondary_vapors);
552 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
553 {
554 _fp_secondary[i]->p_from_T_v(T, v / x[i], p_sec, dp_dT_sec, dp_dv_sec[i]);
555 p += p_sec;
556 dp_dT += dp_dT_sec;
557 dp_dv += dp_dv_sec[i] / x[i];
558 }
559
560 // get the composition dependent derivatives of the secondary vapors
561 dp_dx.resize(_n_secondary_vapors);
562 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
563 {
564 dp_dx[i] = -dp_dv_sec[i] * v / (x[i] * x[i]);
565 for (unsigned int j = 0; j < _n_secondary_vapors; j++)
566 {
567 if (j == i)
568 continue;
569 // Note: this was previously as follows, but we were unable to understand
570 // why and this is not currently tested (requires 3 or more components in
571 // the mixture):
572 // dxj_dxi = -x[j] / (1. - x[i]);
573 if (_n_secondary_vapors > 1)
575 "The mass fraction derivatives in the following "
576 "function have not been tested for mixtures of 3 or more components:\n\n",
577 __PRETTY_FUNCTION__);
578 dxj_dxi = 0;
579 dp_dx[i] += -dp_dv_sec[j] * v / (x[j] * x[j]) * dxj_dxi;
580 }
581 dp_dx[i] += -dp_dv_primary * v / (x_primary * x_primary) * (-x_primary / (1. - x[i]));
582 }
583}
584
585Real
586IdealRealGasMixtureFluidProperties::c_from_T_v(Real T, Real v, const std::vector<Real> & x) const
587{
588 Real p, dp_dT, dp_dv;
589 std::vector<Real> dp_dx;
590 p_from_T_v(T, v, x, p, dp_dT, dp_dv, dp_dx);
591
592 Real s, ds_dT, ds_dv;
593 std::vector<Real> ds_dx;
594 s_from_T_v(T, v, x, s, ds_dT, ds_dv, ds_dx);
595
596 const Real dp_dv_s = dp_dv - dp_dT * ds_dv / ds_dT;
597
598 if (dp_dv_s >= 0)
599 mooseWarning("c_from_T_v(): dp_dv_s = ", dp_dv_s, ". Should be negative.");
600 return v * std::sqrt(-dp_dv_s);
601}
602
603void
604IdealRealGasMixtureFluidProperties::c_from_T_v(Real T,
605 Real v,
606 const std::vector<Real> & x,
607 Real & c,
608 Real & dc_dT,
609 Real & dc_dv,
610 std::vector<Real> & dc_dx) const
611{
612 c = c_from_T_v(T, v, x);
613
614 // For derived properties, we would need higher order derivatives;
615 // therefore, numerical derivatives are used here.
616 const Real dT = 1.e-6;
617 const Real T_perturbed = T + dT;
618 Real c_perturbed = c_from_T_v(T_perturbed, v, x);
619 dc_dT = (c_perturbed - c) / dT;
620
621 const Real dv = v * 1.e-6;
622 const Real v_perturbed = v + dv;
623 c_perturbed = c_from_T_v(T, v_perturbed, x);
624 dc_dv = (c_perturbed - c) / dv;
625
626 dc_dx.resize(_n_secondary_vapors);
627 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
628 {
629 std::vector<Real> x_perturbed(x);
630 const Real dx_i = 1e-6;
631 for (unsigned int j = 0; j < _n_secondary_vapors; j++)
632 {
633 if (j != i)
634 x_perturbed[j] =
635 x[j] * (1.0 - (x[i] + dx_i)) / (1.0 - x[i]); // recalculate new mass fractions
636 }
637 x_perturbed[i] += dx_i;
638 c_perturbed = c_from_T_v(T, v, x_perturbed);
639 dc_dx[i] = ((c_perturbed - c) / dx_i);
640 }
641}
642
643Real
644IdealRealGasMixtureFluidProperties::cp_from_T_v(Real T, Real v, const std::vector<Real> & x) const
645{
646 const Real x_primary = primaryMassFraction(x);
647
648 Real p, dp_dT, dp_dv;
649 std::vector<Real> dp_dx;
650 p_from_T_v(T, v, x, p, dp_dT, dp_dv, dp_dx);
651
652 Real h, dh_dT, dh_dv;
653 _fp_primary->h_from_T_v(T, v / x_primary, h, dh_dT, dh_dv);
654 const Real cp_primary = dh_dT - dh_dv * dp_dT / dp_dv;
655 Real cp = x_primary * cp_primary;
656
657 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
658 {
659 _fp_secondary[i]->h_from_T_v(T, v / x[i], h, dh_dT, dh_dv);
660 const Real cp_sec = dh_dT - dh_dv * dp_dT / dp_dv;
661 cp += x[i] * cp_sec;
662 }
663
664 return cp;
665}
666
667void
668IdealRealGasMixtureFluidProperties::cp_from_T_v(Real T,
669 Real v,
670 const std::vector<Real> & x,
671 Real & cp,
672 Real & dcp_dT,
673 Real & dcp_dv,
674 std::vector<Real> & dcp_dx) const
675{
676 const Real x_primary = primaryMassFraction(x);
677
678 Real p, dp_dT, dp_dv;
679 std::vector<Real> dp_dx;
680 p_from_T_v(T, v, x, p, dp_dT, dp_dv, dp_dx);
681
682 Real h, dh_dT, dh_dv;
683 _fp_primary->h_from_T_v(T, v / x_primary, h, dh_dT, dh_dv);
684 const Real cp_primary = dh_dT - dh_dv * dp_dT / dp_dv;
685 cp = x_primary * cp_primary;
686
687 // Neglect these for now. These require higher-order derivatives, so a finite
688 // difference should probably be used, as for sound speed.
689 imperfectJacobianMessage("The temperature and specific volume derivatives in the following "
690 "function are currently neglected:\n\n",
691 __PRETTY_FUNCTION__);
692 dcp_dT = 0;
693 dcp_dv = 0;
694
695 dcp_dx.resize(_n_secondary_vapors);
696 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
697 {
698 _fp_secondary[i]->h_from_T_v(T, v / x[i], h, dh_dT, dh_dv);
699 const Real cp_sec = dh_dT - dh_dv * dp_dT / dp_dv;
700 cp += x[i] * cp_sec;
701 dcp_dx[i] = cp_sec - cp_primary;
702 }
703}
704
705Real
706IdealRealGasMixtureFluidProperties::cv_from_T_v(Real T, Real v, const std::vector<Real> & x) const
707{
708 const Real x_primary = primaryMassFraction(x);
709 Real cv = x_primary * _fp_primary->cv_from_T_v(T, v / x_primary);
710
711 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
712 cv += x[i] * _fp_secondary[i]->cv_from_T_v(T, v / x[i]);
713
714 return cv;
715}
716
717void
718IdealRealGasMixtureFluidProperties::cv_from_T_v(Real T,
719 Real v,
720 const std::vector<Real> & x,
721 Real & cv,
722 Real & dcv_dT,
723 Real & dcv_dv,
724 std::vector<Real> & dcv_dx) const
725{
726 const Real x_primary = primaryMassFraction(x);
727 const Real cv_primary = _fp_primary->cv_from_T_v(T, v / x_primary);
728 cv = x_primary * cv_primary;
729
730 // Neglect these for now. These require higher-order derivatives, so a finite
731 // difference should probably be used, as for sound speed.
732 imperfectJacobianMessage("The temperature and specific volume derivatives in the following "
733 "function are currently neglected:\n\n",
734 __PRETTY_FUNCTION__);
735 dcv_dT = 0;
736 dcv_dv = 0;
737
738 dcv_dx.resize(_n_secondary_vapors);
739 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
740 {
741 const Real cv_sec = _fp_secondary[i]->cv_from_T_v(T, v / x[i]);
742 cv += x[i] * cv_sec;
743 dcv_dx[i] = cv_sec - cv_primary;
744 }
745}
746
747Real
749 Real T,
750 const std::vector<Real> & x) const
751{
752 Real T_c = _fp_primary->criticalTemperature();
753 Real xs;
754 if (T > T_c)
755 // return 1. to indicate that the primary fluid will not condense for
756 // given (p,T)
757 xs = 1.;
758 else
759 {
760 Real pp_sat = _fp_primary->pp_sat_from_p_T(p, T);
761 if (pp_sat < 0.)
762 {
763 // return 1. to indicate that the primary fluid will not condense for
764 // given (p,T)
765 xs = 1.;
766 return xs;
767 }
768 Real v_primary = _fp_primary->v_from_p_T(pp_sat, T);
769 Real pp_sat_secondary = p - pp_sat;
770
771 Real v_secondary;
772 if (_n_secondary_vapors == 1)
773 v_secondary = _fp_secondary[0]->v_from_p_T(pp_sat_secondary, T);
774 else
775 {
776 std::vector<Real> x_sec(_n_secondary_vapors);
777 const Real x_primary = primaryMassFraction(x);
778 Real sum = 0.;
779 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
780 {
781 x_sec[i] = x[i] / (1. - x_primary);
782 sum += x_sec[i] / _fp_secondary[i]->molarMass();
783 }
784 Real M_star = 1. / sum;
785 v_secondary = R_molar * T / (M_star * pp_sat_secondary);
786 int it = 0;
787 double f = 1., df_dvs, pp_sec, p_sec, dp_dT_sec, dp_dv_sec, dp_dv;
788 double tol_p = 1.e-8;
789 while (std::fabs(f / pp_sat_secondary) > tol_p)
790 {
791 pp_sec = 0.;
792 dp_dv = 0.;
793 for (unsigned int i = 0; i < _n_secondary_vapors; i++)
794 {
795 _fp_secondary[i]->p_from_T_v(T, v_secondary / x_sec[i], p_sec, dp_dT_sec, dp_dv_sec);
796 pp_sec += p_sec;
797 dp_dv += dp_dv_sec / x_sec[i];
798 }
799 f = pp_sec - pp_sat_secondary;
800 df_dvs = dp_dv;
801 v_secondary -= f / df_dvs;
802 if (it++ > 15)
803 return getNaN();
804 }
805 }
806
807 xs = v_secondary / (v_primary + v_secondary);
808 }
809
810 return xs;
811}
Real f(Real x)
Test function for Brents method.
const std::vector< double > x
const double mu
const Real p
const double rho
const double T
const double v
registerMooseObject("FluidPropertiesApp", IdealRealGasMixtureFluidProperties)
Class for fluid properties of an arbitrary vapor mixture.
static constexpr const Real R_molar
molar (or universal) gas constant
virtual const SinglePhaseFluidProperties & getSecondaryFluidProperties(unsigned int i=0) const override
Gets a secondary component single-phase fluid properties.
const Real _T_mix_max
maximum temperature of all components
Real xs_prim_from_p_T(Real p, Real T, const std::vector< Real > &x) const
Mass fraction of primary (condensable) component at saturation from pressure and temperature.
const unsigned int _n_secondary_vapors
Number of secondary vapors.
IdealRealGasMixtureFluidProperties(const InputParameters &parameters)
std::vector< const SinglePhaseFluidProperties * > _fp_secondary
Secondary vapor fluid properties.
const SinglePhaseFluidProperties *const _fp_primary
Primary vapor fluid properties.
virtual const SinglePhaseFluidProperties & getPrimaryFluidProperties() const override
Gets the primary component single-phase fluid properties.
static InputParameters validParams()
const std::vector< UserObjectName > _fp_secondary_names
Names of secondary vapor fluid properties.
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 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()
Common class for single phase fluid properties.
virtual Real criticalInternalEnergy() const
Critical specific internal energy.
virtual void v_e_spndl_from_T(Real T, Real &v, Real &e) const
Specific internal energy from temperature and specific volume.
virtual Real e_spndl_from_v(Real v) const
Specific internal energy from temperature and specific volume.
virtual Real criticalDensity() const
Critical density.
virtual Real molarMass() const
Molar mass [kg/mol].
virtual Real criticalTemperature() const
Critical temperature.
Base class for fluid properties of vapor mixtures.
Real primaryMassFraction(const std::vector< Real > &x) const
Computes the mass fraction of the primary vapor given mass fractions of the secondary vapors.
void imperfectJacobianMessage(Args... args) const
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.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real