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 "SimpleFluidProperties.h"
11 : #include "NewtonInversion.h"
12 :
13 : registerMooseObject("FluidPropertiesApp", SimpleFluidProperties);
14 :
15 : InputParameters
16 20 : SimpleFluidProperties::validParams()
17 : {
18 20 : InputParameters params = SinglePhaseFluidProperties::validParams();
19 40 : params.addParam<Real>("molar_mass", 1.8E-2, "Constant molar mass of the fluid (kg/mol)");
20 40 : params.addParam<Real>(
21 40 : "thermal_expansion", 2.14E-4, "Constant coefficient of thermal expansion (1/K)");
22 40 : params.addParam<Real>(
23 40 : "cv", 4186.0, "Constant specific heat capacity at constant volume (J/kg/K)");
24 40 : params.addParam<Real>(
25 40 : "cp", 4194.0, "Constant specific heat capacity at constant pressure (J/kg/K)");
26 60 : params.addRangeCheckedParam<Real>(
27 40 : "bulk_modulus", 2.0E9, "bulk_modulus>0", "Constant bulk modulus (Pa)");
28 40 : params.addParam<Real>("thermal_conductivity", 0.6, "Constant thermal conductivity (W/m/K)");
29 40 : params.addParam<Real>("specific_entropy", 300.0, "Constant specific entropy (J/kg/K)");
30 40 : params.addParam<Real>("viscosity", 1.0E-3, "Constant dynamic viscosity (Pa.s)");
31 40 : params.addParam<Real>("density0", 1000.0, "Density at zero pressure and zero temperature");
32 40 : params.addParam<Real>("porepressure_coefficient",
33 40 : 1.0,
34 : "The enthalpy is internal_energy + P / density * "
35 : "porepressure_coefficient. Physically this should be 1.0, "
36 : "but analytic solutions are simplified when it is zero");
37 20 : params.addClassDescription("Fluid properties for a simple fluid with a constant bulk density");
38 20 : return params;
39 0 : }
40 :
41 10 : SimpleFluidProperties::SimpleFluidProperties(const InputParameters & parameters)
42 : : SinglePhaseFluidProperties(parameters),
43 10 : _molar_mass(getParam<Real>("molar_mass")),
44 20 : _thermal_expansion(getParam<Real>("thermal_expansion")),
45 20 : _cv(getParam<Real>("cv")),
46 20 : _cp(getParam<Real>("cp")),
47 20 : _bulk_modulus(getParam<Real>("bulk_modulus")),
48 20 : _thermal_conductivity(getParam<Real>("thermal_conductivity")),
49 20 : _specific_entropy(getParam<Real>("specific_entropy")),
50 20 : _viscosity(getParam<Real>("viscosity")),
51 20 : _density0(getParam<Real>("density0")),
52 30 : _pp_coeff(getParam<Real>("porepressure_coefficient"))
53 : {
54 10 : }
55 :
56 20 : SimpleFluidProperties::~SimpleFluidProperties() {}
57 :
58 : std::string
59 1 : SimpleFluidProperties::fluidName() const
60 : {
61 1 : return "simple_fluid";
62 : }
63 :
64 : Real
65 1 : SimpleFluidProperties::molarMass() const
66 : {
67 1 : return _molar_mass;
68 : }
69 :
70 : Real
71 2 : SimpleFluidProperties::beta_from_p_T(Real /*pressure*/, Real /*temperature*/) const
72 : {
73 2 : return _thermal_expansion;
74 : }
75 :
76 : void
77 0 : SimpleFluidProperties::beta_from_p_T(
78 : Real pressure, Real temperature, Real & beta, Real & dbeta_dp, Real & dbeta_dT) const
79 : {
80 0 : beta = beta_from_p_T(pressure, temperature);
81 0 : dbeta_dp = 0.0;
82 0 : dbeta_dT = 0.0;
83 0 : }
84 :
85 : Real
86 14 : SimpleFluidProperties::cp_from_p_T(Real /*pressure*/, Real /*temperature*/) const
87 : {
88 14 : return _cp;
89 : }
90 :
91 : void
92 2 : SimpleFluidProperties::cp_from_p_T(
93 : Real pressure, Real temperature, Real & cp, Real & dcp_dp, Real & dcp_dT) const
94 : {
95 2 : cp = cp_from_p_T(pressure, temperature);
96 2 : dcp_dp = 0.0;
97 2 : dcp_dT = 0.0;
98 2 : }
99 :
100 : Real
101 14 : SimpleFluidProperties::cp_from_v_e(Real /*v*/, Real /*e*/) const
102 : {
103 14 : return _cp;
104 : }
105 :
106 : void
107 2 : SimpleFluidProperties::cp_from_v_e(Real v, Real e, Real & cp, Real & dcp_dv, Real & dcp_de) const
108 : {
109 2 : cp = cp_from_v_e(v, e);
110 2 : dcp_dv = 0.0;
111 2 : dcp_de = 0.0;
112 2 : }
113 :
114 : Real
115 14 : SimpleFluidProperties::cv_from_p_T(Real /*pressure*/, Real /*temperature*/) const
116 : {
117 14 : return _cv;
118 : }
119 :
120 : void
121 2 : SimpleFluidProperties::cv_from_p_T(
122 : Real pressure, Real temperature, Real & cv, Real & dcv_dp, Real & dcv_dT) const
123 : {
124 2 : cv = cv_from_p_T(pressure, temperature);
125 2 : dcv_dp = 0.0;
126 2 : dcv_dT = 0.0;
127 2 : }
128 :
129 : Real
130 14 : SimpleFluidProperties::cv_from_v_e(Real /*v*/, Real /*e*/) const
131 : {
132 14 : return _cv;
133 : }
134 :
135 : void
136 2 : SimpleFluidProperties::cv_from_v_e(Real v, Real e, Real & cv, Real & dcv_dv, Real & dcv_de) const
137 : {
138 2 : cv = cv_from_v_e(v, e);
139 2 : dcv_dv = 0.0;
140 2 : dcv_de = 0.0;
141 2 : }
142 :
143 : Real
144 14 : SimpleFluidProperties::c_from_p_T(Real pressure, Real temperature) const
145 : {
146 14 : return std::sqrt(_bulk_modulus / rho_from_p_T(pressure, temperature));
147 : }
148 :
149 : void
150 2 : SimpleFluidProperties::c_from_p_T(
151 : Real pressure, Real temperature, Real & c, Real & dc_dp, Real & dc_dT) const
152 : {
153 2 : c = c_from_p_T(pressure, temperature);
154 2 : dc_dp =
155 2 : -std::exp(_thermal_expansion * temperature - pressure / _bulk_modulus) / 2 / _density0 /
156 2 : std::sqrt(_bulk_modulus *
157 2 : std::exp(_thermal_expansion * temperature - pressure / _bulk_modulus) / _density0);
158 2 : dc_dT =
159 2 : _thermal_expansion * _bulk_modulus *
160 2 : std::exp(_thermal_expansion * temperature - pressure / _bulk_modulus) / 2 / _density0 /
161 2 : std::sqrt(_bulk_modulus *
162 2 : std::exp(_thermal_expansion * temperature - pressure / _bulk_modulus) / _density0);
163 2 : }
164 :
165 : Real
166 12 : SimpleFluidProperties::c_from_v_e(Real v, Real /*e*/) const
167 : {
168 12 : return std::sqrt(_bulk_modulus * v);
169 : }
170 :
171 : void
172 2 : SimpleFluidProperties::c_from_v_e(Real v, Real /*e*/, Real & c, Real & dc_dv, Real & dc_de) const
173 : {
174 2 : c = std::sqrt(_bulk_modulus * v);
175 :
176 2 : dc_dv = 0.5 * std::sqrt(_bulk_modulus / v);
177 2 : dc_de = 0.0;
178 2 : }
179 :
180 : Real
181 14 : SimpleFluidProperties::k_from_p_T(Real /*pressure*/, Real /*temperature*/) const
182 : {
183 14 : return _thermal_conductivity;
184 : }
185 :
186 : void
187 2 : SimpleFluidProperties::k_from_p_T(
188 : Real /*pressure*/, Real /*temperature*/, Real & k, Real & dk_dp, Real & dk_dT) const
189 : {
190 2 : k = _thermal_conductivity;
191 2 : dk_dp = 0;
192 2 : dk_dT = 0;
193 2 : }
194 :
195 : Real
196 12 : SimpleFluidProperties::k_from_v_e(Real /*v*/, Real /*e*/) const
197 : {
198 12 : return _thermal_conductivity;
199 : }
200 :
201 : void
202 2 : SimpleFluidProperties::k_from_v_e(
203 : Real /*v*/, Real /*e*/, Real & k, Real & dk_dv, Real & dk_de) const
204 : {
205 2 : k = _thermal_conductivity;
206 2 : dk_dv = 0;
207 2 : dk_de = 0;
208 2 : }
209 :
210 : Real
211 12 : SimpleFluidProperties::s_from_p_T(Real /*pressure*/, Real /*temperature*/) const
212 : {
213 12 : return _specific_entropy;
214 : }
215 :
216 : void
217 2 : SimpleFluidProperties::s_from_p_T(
218 : Real /*p*/, Real /*T*/, Real & s, Real & ds_dp, Real & ds_dT) const
219 : {
220 2 : s = _specific_entropy;
221 2 : ds_dp = 0;
222 2 : ds_dT = 0;
223 2 : }
224 :
225 : Real
226 2 : SimpleFluidProperties::s_from_h_p(Real /*enthalpy*/, Real /*pressure*/) const
227 : {
228 2 : return _specific_entropy;
229 : }
230 :
231 : Real
232 12 : SimpleFluidProperties::s_from_v_e(Real /*v*/, Real /*e*/) const
233 : {
234 12 : return _specific_entropy;
235 : }
236 :
237 : void
238 2 : SimpleFluidProperties::s_from_v_e(
239 : Real /*v*/, Real /*e*/, Real & s, Real & ds_dv, Real & ds_de) const
240 : {
241 2 : s = _specific_entropy;
242 2 : ds_dv = 0;
243 2 : ds_de = 0;
244 2 : }
245 :
246 : Real
247 74 : SimpleFluidProperties::rho_from_p_T(Real pressure, Real temperature) const
248 : {
249 74 : return _density0 * std::exp(pressure / _bulk_modulus - _thermal_expansion * temperature);
250 : }
251 :
252 : void
253 13 : SimpleFluidProperties::rho_from_p_T(
254 : Real pressure, Real temperature, Real & rho, Real & drho_dp, Real & drho_dT) const
255 : {
256 13 : rho = this->rho_from_p_T(pressure, temperature);
257 13 : drho_dp = rho / _bulk_modulus;
258 13 : drho_dT = -_thermal_expansion * rho;
259 13 : }
260 :
261 : void
262 0 : SimpleFluidProperties::rho_from_p_T(const ADReal & pressure,
263 : const ADReal & temperature,
264 : ADReal & rho,
265 : ADReal & drho_dp,
266 : ADReal & drho_dT) const
267 : {
268 0 : rho = SinglePhaseFluidProperties::rho_from_p_T(pressure, temperature);
269 0 : drho_dp = rho / _bulk_modulus;
270 0 : drho_dT = -_thermal_expansion * rho;
271 0 : }
272 :
273 : Real
274 34 : SimpleFluidProperties::T_from_v_e(Real /*v*/, Real e) const
275 : {
276 : // NOTE: while e = _cv * T, h is not equal to _cp * T
277 34 : return e / _cv;
278 : }
279 :
280 : Real
281 28 : SimpleFluidProperties::T_from_v_h(Real v, Real h) const
282 : {
283 28 : return (std::log(1. / _density0 / v) - h / v / _bulk_modulus) / -_thermal_expansion /
284 28 : (1 + _cv / v / _thermal_expansion / _bulk_modulus);
285 : }
286 :
287 : void
288 4 : SimpleFluidProperties::T_from_v_h(Real v, Real h, Real & T, Real & dT_dv, Real & dT_dh) const
289 : {
290 4 : T = T_from_v_h(v, h);
291 4 : dT_dv =
292 4 : (1 / -_thermal_expansion) *
293 4 : ((-1. / v + h / v / v / _bulk_modulus) * (1 + _cv / v / _thermal_expansion / _bulk_modulus) -
294 4 : (std::log(1. / _density0 / v) - h / v / _bulk_modulus) *
295 4 : (-_cv / v / v / _thermal_expansion / _bulk_modulus)) /
296 : Utility::pow<2>(1 + _cv / v / _thermal_expansion / _bulk_modulus);
297 4 : dT_dh = (-1 / v / _bulk_modulus) / -_thermal_expansion /
298 4 : (1 + _cv / v / _thermal_expansion / _bulk_modulus);
299 4 : }
300 :
301 : void
302 2 : SimpleFluidProperties::T_from_v_e(Real v, Real e, Real & T, Real & dT_dv, Real & dT_de) const
303 : {
304 2 : T = T_from_v_e(v, e);
305 2 : dT_dv = 0.0;
306 2 : dT_de = 1.0 / _cv;
307 2 : }
308 :
309 : void
310 0 : SimpleFluidProperties::T_from_v_e(
311 : const ADReal & v, const ADReal & e, ADReal & T, ADReal & dT_dv, ADReal & dT_de) const
312 : {
313 0 : T = SinglePhaseFluidProperties::T_from_v_e(v, e);
314 0 : dT_dv = 0.0;
315 0 : dT_de = 1.0 / _cv;
316 0 : }
317 :
318 : Real
319 28 : SimpleFluidProperties::T_from_p_rho(Real p, Real rho) const
320 : {
321 : mooseAssert(rho > 0, "Density should be positive");
322 28 : return (std::log(rho / _density0) - p / _bulk_modulus) / -_thermal_expansion;
323 : }
324 :
325 : void
326 4 : SimpleFluidProperties::T_from_p_rho(Real p, Real rho, Real & T, Real & dT_dp, Real & dT_drho) const
327 : {
328 4 : T = T_from_p_rho(p, rho);
329 4 : dT_dp = 1 / (_thermal_expansion * _bulk_modulus);
330 4 : dT_drho = 1 / (-_thermal_expansion * rho);
331 4 : }
332 :
333 : Real
334 2 : SimpleFluidProperties::T_from_p_h(Real p, Real h) const
335 : {
336 : // Likely a better guess than user-selected
337 2 : Real T_initial = h / _cp;
338 :
339 : // exponential dependence in rho and linear dependence in e makes it challenging
340 : auto lambda = [&](Real p, Real current_T, Real & new_rho, Real & dh_dp, Real & dh_dT)
341 6 : { h_from_p_T(p, current_T, new_rho, dh_dp, dh_dT); };
342 2 : return FluidPropertiesUtils::NewtonSolve(
343 2 : p, h, T_initial, _tolerance, lambda, name() + "::T_from_p_h", _max_newton_its)
344 2 : .first;
345 : }
346 :
347 : Real
348 16 : SimpleFluidProperties::p_from_v_e(Real v, Real e) const
349 : {
350 16 : Real temperature = T_from_v_e(v, e);
351 16 : return _bulk_modulus * (_thermal_expansion * temperature + std::log(1 / (v * _density0)));
352 : }
353 :
354 : void
355 2 : SimpleFluidProperties::p_from_v_e(Real v, Real e, Real & p, Real & dp_dv, Real & dp_de) const
356 : {
357 2 : p = p_from_v_e(v, e);
358 2 : dp_dv = -_bulk_modulus / v;
359 2 : dp_de = _bulk_modulus * _thermal_expansion / _cv;
360 2 : }
361 :
362 : void
363 0 : SimpleFluidProperties::p_from_v_e(
364 : const ADReal & v, const ADReal & e, ADReal & p, ADReal & dp_dv, ADReal & dp_de) const
365 : {
366 0 : p = SinglePhaseFluidProperties::p_from_v_e(v, e);
367 0 : dp_dv = -_bulk_modulus / v;
368 0 : dp_de = _bulk_modulus * _thermal_expansion / _cv;
369 0 : }
370 :
371 : Real
372 12 : SimpleFluidProperties::p_from_v_h(Real v, Real h) const
373 : {
374 12 : Real T = T_from_v_h(v, h);
375 12 : return _bulk_modulus * (_thermal_expansion * T + std::log(1 / (v * _density0)));
376 : }
377 :
378 : void
379 2 : SimpleFluidProperties::p_from_v_h(Real v, Real h, Real & p, Real & dp_dv, Real & dp_dh) const
380 : {
381 : Real T, dT_dv, dT_dh;
382 2 : T_from_v_h(v, h, T, dT_dv, dT_dh);
383 2 : p = _bulk_modulus * (_thermal_expansion * T + std::log(1 / (v * _density0)));
384 2 : dp_dv = _bulk_modulus * (_thermal_expansion * dT_dv - 1. / v);
385 2 : dp_dh = _bulk_modulus * (_thermal_expansion * dT_dh);
386 2 : }
387 :
388 : Real
389 70 : SimpleFluidProperties::e_from_p_T(Real /*pressure*/, Real temperature) const
390 : {
391 70 : return _cv * temperature;
392 : }
393 :
394 : void
395 8 : SimpleFluidProperties::e_from_p_T(
396 : Real pressure, Real temperature, Real & e, Real & de_dp, Real & de_dT) const
397 : {
398 8 : e = this->e_from_p_T(pressure, temperature);
399 8 : de_dp = 0.0;
400 8 : de_dT = _cv;
401 8 : }
402 :
403 : Real
404 12 : SimpleFluidProperties::e_from_p_rho(Real p, Real rho) const
405 : {
406 12 : Real T = T_from_p_rho(p, rho);
407 12 : return e_from_p_T(p, T);
408 : }
409 :
410 : void
411 2 : SimpleFluidProperties::e_from_p_rho(Real p, Real rho, Real & e, Real & de_dp, Real & de_drho) const
412 : {
413 : // get temperature and derivatives
414 : Real T, dT_dp, dT_drho;
415 2 : T_from_p_rho(p, rho, T, dT_dp, dT_drho);
416 :
417 : // get energy and derivatives
418 : Real de_dT;
419 2 : e_from_p_T(p, T, e, de_dp, de_dT);
420 2 : de_dp = de_dT * dT_dp + de_dp;
421 2 : de_drho = de_dT * dT_drho + de_dp * dT_dp;
422 2 : }
423 :
424 : Real
425 12 : SimpleFluidProperties::e_from_v_h(Real v, Real h) const
426 : {
427 12 : Real T = T_from_v_h(v, h);
428 12 : Real p = p_from_v_h(v, h);
429 12 : return e_from_p_T(p, T);
430 : }
431 :
432 : void
433 2 : SimpleFluidProperties::e_from_v_h(Real v, Real h, Real & e, Real & de_dv, Real & de_dh) const
434 : {
435 : Real T, dT_dv, dT_dh;
436 : Real p, dp_dv, dp_dh;
437 2 : T_from_v_h(v, h, T, dT_dv, dT_dh);
438 2 : p_from_v_h(v, h, p, dp_dv, dp_dh);
439 :
440 : Real de_dp, de_dT;
441 2 : e_from_p_T(p, T, e, de_dp, de_dT);
442 2 : de_dv = de_dp * dp_dv + de_dT * dT_dv;
443 2 : de_dh = de_dp * dp_dh + de_dT * dT_dh;
444 2 : }
445 :
446 : Real
447 17 : SimpleFluidProperties::mu_from_p_T(Real /*pressure*/, Real /*temperature*/) const
448 : {
449 17 : return _viscosity;
450 : }
451 :
452 : void
453 4 : SimpleFluidProperties::mu_from_p_T(
454 : Real pressure, Real temperature, Real & mu, Real & dmu_dp, Real & dmu_dT) const
455 : {
456 4 : mu = this->mu_from_p_T(pressure, temperature);
457 4 : dmu_dp = 0.0;
458 4 : dmu_dT = 0.0;
459 4 : }
460 :
461 : Real
462 14 : SimpleFluidProperties::mu_from_v_e(Real /*v*/, Real /*e*/) const
463 : {
464 14 : return _viscosity;
465 : }
466 :
467 : void
468 2 : SimpleFluidProperties::mu_from_v_e(Real v, Real e, Real & mu, Real & dmu_dv, Real & dmu_de) const
469 : {
470 2 : mu = this->mu_from_v_e(v, e);
471 2 : dmu_dv = 0.0;
472 2 : dmu_de = 0.0;
473 2 : }
474 :
475 : Real
476 26 : SimpleFluidProperties::h_from_p_T(Real pressure, Real temperature) const
477 : {
478 26 : return e_from_p_T(pressure, temperature) +
479 26 : _pp_coeff * pressure / rho_from_p_T(pressure, temperature);
480 : }
481 :
482 : void
483 8 : SimpleFluidProperties::h_from_p_T(
484 : Real pressure, Real temperature, Real & h, Real & dh_dp, Real & dh_dT) const
485 : {
486 8 : h = this->h_from_p_T(pressure, temperature);
487 :
488 : Real density, ddensity_dp, ddensity_dT;
489 8 : rho_from_p_T(pressure, temperature, density, ddensity_dp, ddensity_dT);
490 :
491 8 : dh_dp = _pp_coeff / density - _pp_coeff * pressure * ddensity_dp / density / density;
492 8 : dh_dT = _cv - _pp_coeff * pressure * ddensity_dT / density / density;
493 8 : }
|
