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 "GeneralFunctorFluidProps.h"
11 : #include "NS.h" // Variable Term Names
12 : #include "HeatTransferUtils.h"
13 : #include "NavierStokesMethods.h"
14 : #include "SinglePhaseFluidProperties.h"
15 :
16 : registerMooseObject("NavierStokesApp", GeneralFunctorFluidProps);
17 : registerMooseObject("NavierStokesApp", NonADGeneralFunctorFluidProps);
18 :
19 : template <bool is_ad>
20 : InputParameters
21 428 : GeneralFunctorFluidPropsTempl<is_ad>::validParams()
22 : {
23 428 : auto params = FunctorMaterial::validParams();
24 428 : params.addRequiredParam<UserObjectName>(NS::fluid, "Fluid properties object");
25 428 : params.addClassDescription("Creates functor fluid properties using a (P, T) formulation");
26 :
27 428 : params.addRequiredParam<MooseFunctorName>(NS::pressure, "Pressure");
28 428 : params.addRequiredParam<MooseFunctorName>(NS::T_fluid, "Fluid temperature");
29 428 : params.addRequiredParam<MooseFunctorName>(NS::speed, "Velocity norm");
30 428 : params.addParam<MooseFunctorName>(NS::density, "Density");
31 856 : params.addParam<bool>(
32 : "force_define_density",
33 856 : false,
34 : "Whether to force the definition of a density functor from the fluid properties");
35 856 : params.addParam<bool>("neglect_derivatives_of_density_time_derivative",
36 856 : true,
37 : "Whether to neglect the derivatives with regards to nonlinear variables "
38 : "of the density time derivatives");
39 :
40 856 : params.addParam<FunctionName>(
41 856 : "mu_rampdown", 1, "A function describing a ramp down of viscosity over time");
42 428 : params.addRequiredParam<MooseFunctorName>(NS::porosity, "porosity");
43 856 : params.addRequiredParam<MooseFunctorName>(
44 : "characteristic_length", "characteristic length for Reynolds number calculation");
45 :
46 : // Dynamic pressure
47 856 : params.addParam<bool>("solving_for_dynamic_pressure",
48 856 : false,
49 : "Whether to solve for the dynamic pressure instead of the total pressure");
50 428 : params.addParam<Point>("reference_pressure_point",
51 428 : Point(0, 0, 0),
52 : "Point at which the gravity term for the static pressure is zero");
53 856 : params.addParam<Real>("reference_pressure", 1e5, "Total pressure at the reference point");
54 856 : params.addParam<Point>("gravity", Point(0, 0, -9.81), "Gravity vector");
55 :
56 856 : params.addParamNamesToGroup(
57 : "solving_for_dynamic_pressure reference_pressure_point reference_pressure",
58 : "Dynamic pressure");
59 :
60 : // Property names
61 856 : params.addParam<MooseFunctorName>(
62 : "density_name", NS::density, "Name to give to the density functor");
63 856 : params.addParam<MooseFunctorName>(
64 : "dynamic_viscosity_name", NS::mu, "Name to give to the dynamic viscosity functor");
65 856 : params.addParam<MooseFunctorName>(
66 : "specific_heat_name", NS::cp, "Name to give to the specific heat (cp) functor");
67 856 : params.addParam<MooseFunctorName>(
68 : "thermal_conductivity_name", NS::k, "Name to give to the thermal conductivity functor");
69 856 : params.addParamNamesToGroup(
70 : "density_name dynamic_viscosity_name specific_heat_name thermal_conductivity_name",
71 : "Functor property names");
72 428 : return params;
73 0 : }
74 :
75 : template <bool is_ad>
76 225 : GeneralFunctorFluidPropsTempl<is_ad>::GeneralFunctorFluidPropsTempl(
77 : const InputParameters & parameters)
78 : : FunctorMaterial(parameters),
79 450 : _fluid(UserObjectInterface::getUserObject<SinglePhaseFluidProperties>(NS::fluid)),
80 225 : _eps(getFunctor<GenericReal<is_ad>>(NS::porosity)),
81 450 : _d(getFunctor<Real>("characteristic_length")),
82 :
83 450 : _pressure_is_dynamic(getParam<bool>("solving_for_dynamic_pressure")),
84 450 : _reference_pressure_point(getParam<Point>("reference_pressure_point")),
85 450 : _reference_pressure_value(getParam<Real>("reference_pressure")),
86 450 : _gravity_vec(getParam<Point>("gravity")),
87 :
88 225 : _pressure(getFunctor<GenericReal<is_ad>>(NS::pressure)),
89 225 : _T_fluid(getFunctor<GenericReal<is_ad>>(NS::T_fluid)),
90 225 : _speed(getFunctor<GenericReal<is_ad>>(NS::speed)),
91 450 : _force_define_density(getParam<bool>("force_define_density")),
92 225 : _rho(getFunctor<GenericReal<is_ad>>(NS::density)),
93 225 : _mu_rampdown(getFunction("mu_rampdown")),
94 225 : _neglect_derivatives_of_density_time_derivative(
95 225 : getParam<bool>("neglect_derivatives_of_density_time_derivative")),
96 :
97 450 : _density_name(getParam<MooseFunctorName>("density_name")),
98 450 : _dynamic_viscosity_name(getParam<MooseFunctorName>("dynamic_viscosity_name")),
99 450 : _specific_heat_name(getParam<MooseFunctorName>("specific_heat_name")),
100 675 : _thermal_conductivity_name(getParam<MooseFunctorName>("thermal_conductivity_name"))
101 : {
102 : // Check parameters
103 442 : if (!_pressure_is_dynamic &&
104 1093 : (isParamSetByUser("reference_pressure_point") || isParamSetByUser("reference_pressure")))
105 2 : paramError("solving_for_dynamic_pressure",
106 : "'reference_pressure_point' and 'reference_pressure' should not be set unless "
107 : "solving for the dynamic pressure");
108 :
109 : //
110 : // Set material properties functors
111 : //
112 :
113 : // If pressure is dynamic (else case), we must obtain the total pressure
114 : // We duplicate all this code. An alternative would be to redefine the pressure functor.
115 : // The issue with that is that you need the density functor to define the pressure,
116 : // and the pressure functor to define the density.
117 : // This could be solved by keeping a pointer to the pressure functor as an attribute and set
118 : // the pressure functor after the density functor has been defined.
119 223 : if (!_pressure_is_dynamic)
120 : {
121 215 : if (!isParamValid(NS::density) || _force_define_density)
122 615 : addFunctorProperty<GenericReal<is_ad>>(
123 205 : _density_name,
124 26837480 : [this](const auto & r, const auto & t) -> GenericReal<is_ad>
125 26837480 : { return _fluid.rho_from_p_T(_pressure(r, t), _T_fluid(r, t)); });
126 :
127 645 : addFunctorProperty<GenericReal<is_ad>>(
128 : NS::cv,
129 300 : [this](const auto & r, const auto & t) -> GenericReal<is_ad>
130 300 : { return _fluid.cv_from_p_T(_pressure(r, t), _T_fluid(r, t)); });
131 :
132 645 : const auto & cp = addFunctorProperty<GenericReal<is_ad>>(
133 215 : _specific_heat_name,
134 5734822 : [this](const auto & r, const auto & t) -> GenericReal<is_ad>
135 5734822 : { return _fluid.cp_from_p_T(_pressure(r, t), _T_fluid(r, t)); });
136 :
137 645 : const auto & mu = addFunctorProperty<GenericReal<is_ad>>(
138 215 : _dynamic_viscosity_name,
139 12047838 : [this](const auto & r, const auto & t) -> GenericReal<is_ad>
140 12047838 : { return _mu_rampdown(r, t) * _fluid.mu_from_p_T(_pressure(r, t), _T_fluid(r, t)); });
141 :
142 645 : const auto & k = addFunctorProperty<GenericReal<is_ad>>(
143 215 : _thermal_conductivity_name,
144 6122903 : [this](const auto & r, const auto & t) -> GenericReal<is_ad>
145 6122903 : { return _fluid.k_from_p_T(_pressure(r, t), _T_fluid(r, t)); });
146 :
147 : //
148 : // Time derivatives of fluid properties
149 : //
150 215 : if (_neglect_derivatives_of_density_time_derivative)
151 : {
152 820 : addFunctorProperty<GenericReal<is_ad>>(
153 205 : NS::time_deriv(_density_name),
154 1419772 : [this](const auto & r, const auto & t) -> GenericReal<is_ad>
155 : {
156 : Real rho, drho_dp, drho_dT;
157 1419772 : _fluid.rho_from_p_T(MetaPhysicL::raw_value(_pressure(r, t)),
158 1419772 : MetaPhysicL::raw_value(_T_fluid(r, t)),
159 : rho,
160 : drho_dp,
161 : drho_dT);
162 4259316 : return drho_dp * _pressure.dot(r, t) + drho_dT * _T_fluid.dot(r, t);
163 : });
164 : }
165 : else
166 : {
167 40 : addFunctorProperty<GenericReal<is_ad>>(
168 10 : NS::time_deriv(_density_name),
169 17080 : [this](const auto & r, const auto & t) -> GenericReal<is_ad>
170 : {
171 : GenericReal<is_ad> rho, drho_dp, drho_dT;
172 17080 : _fluid.rho_from_p_T(_pressure(r, t), _T_fluid(r, t), rho, drho_dp, drho_dT);
173 51240 : return drho_dp * _pressure.dot(r, t) + drho_dT * _T_fluid.dot(r, t);
174 : });
175 : }
176 :
177 860 : addFunctorProperty<GenericReal<is_ad>>(
178 215 : NS::time_deriv(_specific_heat_name),
179 300 : [this](const auto & r, const auto & t) -> GenericReal<is_ad>
180 : {
181 : Real dcp_dp, dcp_dT, dummy;
182 300 : const auto raw_pressure = MetaPhysicL::raw_value(_pressure(r, t));
183 300 : const auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
184 300 : _fluid.cp_from_p_T(raw_pressure, raw_T_fluid, dummy, dcp_dp, dcp_dT);
185 :
186 600 : return dcp_dp * _pressure.dot(r, t) + dcp_dT * _T_fluid.dot(r, t);
187 : });
188 :
189 : //
190 : // Temperature and pressure derivatives, to help with computing time derivatives
191 : //
192 :
193 860 : const auto & drho_dp = addFunctorProperty<Real>(
194 215 : derivativePropertyNameFirst(_density_name, NS::pressure),
195 600 : [this](const auto & r, const auto & t) -> Real
196 : {
197 : Real drho_dp, drho_dT, dummy;
198 600 : auto raw_pressure = MetaPhysicL::raw_value(_pressure(r, t));
199 600 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
200 :
201 600 : _fluid.rho_from_p_T(raw_pressure, raw_T_fluid, dummy, drho_dp, drho_dT);
202 600 : return drho_dp;
203 : });
204 :
205 860 : const auto & drho_dT = addFunctorProperty<Real>(
206 215 : derivativePropertyNameFirst(_density_name, NS::T_fluid),
207 600 : [this](const auto & r, const auto & t) -> Real
208 : {
209 : Real drho_dp, drho_dT, dummy;
210 600 : auto raw_pressure = MetaPhysicL::raw_value(_pressure(r, t));
211 600 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
212 :
213 600 : _fluid.rho_from_p_T(raw_pressure, raw_T_fluid, dummy, drho_dp, drho_dT);
214 600 : return drho_dT;
215 : });
216 :
217 860 : const auto & dcp_dp = addFunctorProperty<Real>(
218 215 : derivativePropertyNameFirst(_specific_heat_name, NS::pressure),
219 600 : [this](const auto & r, const auto & t) -> Real
220 : {
221 : Real dcp_dp, dcp_dT, dummy;
222 600 : auto raw_pressure = MetaPhysicL::raw_value(_pressure(r, t));
223 600 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
224 :
225 600 : _fluid.cp_from_p_T(raw_pressure, raw_T_fluid, dummy, dcp_dp, dcp_dT);
226 600 : return dcp_dp;
227 : });
228 :
229 860 : const auto & dcp_dT = addFunctorProperty<Real>(
230 215 : derivativePropertyNameFirst(_specific_heat_name, NS::T_fluid),
231 600 : [this](const auto & r, const auto & t) -> Real
232 : {
233 : Real dcp_dp, dcp_dT, dummy;
234 600 : auto raw_pressure = MetaPhysicL::raw_value(_pressure(r, t));
235 600 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
236 :
237 600 : _fluid.cp_from_p_T(raw_pressure, raw_T_fluid, dummy, dcp_dp, dcp_dT);
238 600 : return dcp_dT;
239 : });
240 :
241 860 : const auto & dmu_dp = addFunctorProperty<Real>(
242 215 : derivativePropertyNameFirst(_dynamic_viscosity_name, NS::pressure),
243 900 : [this](const auto & r, const auto & t) -> Real
244 : {
245 : Real dmu_dp, dmu_dT, dummy;
246 900 : auto raw_pressure = MetaPhysicL::raw_value(_pressure(r, t));
247 900 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
248 :
249 900 : _fluid.mu_from_p_T(raw_pressure, raw_T_fluid, dummy, dmu_dp, dmu_dT);
250 900 : return _mu_rampdown(r, t) * dmu_dp;
251 : });
252 :
253 860 : const auto & dmu_dT = addFunctorProperty<Real>(
254 215 : derivativePropertyNameFirst(_dynamic_viscosity_name, NS::T_fluid),
255 900 : [this](const auto & r, const auto & t) -> Real
256 : {
257 : Real dmu_dp, dmu_dT, dummy;
258 900 : auto raw_pressure = MetaPhysicL::raw_value(_pressure(r, t));
259 900 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
260 :
261 900 : _fluid.mu_from_p_T(raw_pressure, raw_T_fluid, dummy, dmu_dp, dmu_dT);
262 900 : return _mu_rampdown(r, t) * dmu_dT;
263 : });
264 :
265 860 : const auto & dk_dp = addFunctorProperty<Real>(
266 215 : derivativePropertyNameFirst(_thermal_conductivity_name, NS::pressure),
267 600 : [this](const auto & r, const auto & t) -> Real
268 : {
269 : Real dk_dp, dk_dT, dummy;
270 600 : auto raw_pressure = MetaPhysicL::raw_value(_pressure(r, t));
271 600 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
272 :
273 600 : _fluid.k_from_p_T(raw_pressure, raw_T_fluid, dummy, dk_dp, dk_dT);
274 600 : return dk_dp;
275 : });
276 :
277 860 : const auto & dk_dT = addFunctorProperty<Real>(
278 215 : derivativePropertyNameFirst(_thermal_conductivity_name, NS::T_fluid),
279 600 : [this](const auto & r, const auto & t) -> Real
280 : {
281 : Real dk_dp, dk_dT, dummy;
282 600 : auto raw_pressure = MetaPhysicL::raw_value(_pressure(r, t));
283 600 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
284 :
285 600 : _fluid.k_from_p_T(raw_pressure, raw_T_fluid, dummy, dk_dp, dk_dT);
286 600 : return dk_dT;
287 : });
288 :
289 : //
290 : // Fluid adimensional quantities, used in numerous correlations
291 : //
292 :
293 645 : addFunctorProperty<GenericReal<is_ad>>(
294 : NS::Prandtl,
295 3600 : [&cp, &mu, &k](const auto & r, const auto & t) -> GenericReal<is_ad>
296 : {
297 : static constexpr Real small_number = 1e-8;
298 :
299 : using std::max;
300 3600 : return HeatTransferUtils::prandtl(cp(r, t), mu(r, t), max(k(r, t), small_number));
301 : });
302 :
303 645 : addFunctorProperty<Real>(
304 : derivativePropertyNameFirst(NS::Prandtl, NS::pressure),
305 300 : [&mu, &cp, &k, &dmu_dp, &dcp_dp, &dk_dp](const auto & r, const auto & t) -> Real
306 : {
307 300 : return NS::prandtlPropertyDerivative(MetaPhysicL::raw_value(mu(r, t)),
308 300 : MetaPhysicL::raw_value(cp(r, t)),
309 300 : MetaPhysicL::raw_value(k(r, t)),
310 300 : dmu_dp(r, t),
311 300 : dcp_dp(r, t),
312 600 : dk_dp(r, t));
313 : });
314 :
315 645 : addFunctorProperty<Real>(
316 : derivativePropertyNameFirst(NS::Prandtl, NS::T_fluid),
317 300 : [&mu, &cp, &k, &dmu_dT, &dcp_dT, &dk_dT](const auto & r, const auto & t) -> Real
318 : {
319 300 : return NS::prandtlPropertyDerivative(MetaPhysicL::raw_value(mu(r, t)),
320 300 : MetaPhysicL::raw_value(cp(r, t)),
321 300 : MetaPhysicL::raw_value(k(r, t)),
322 300 : dmu_dT(r, t),
323 300 : dcp_dT(r, t),
324 600 : dk_dT(r, t));
325 : });
326 :
327 : //
328 : // (pore / particle) Reynolds number based on superficial velocity and
329 : // characteristic length. Only call Reynolds() one time to compute all three so that
330 : // we don't redundantly check that viscosity is not too close to zero.
331 : //
332 :
333 645 : const auto & Re = addFunctorProperty<GenericReal<is_ad>>(
334 : NS::Reynolds,
335 4800 : [this, &mu](const auto & r, const auto & t) -> GenericReal<is_ad>
336 : {
337 : static constexpr Real small_number = 1e-8;
338 : using std::max;
339 : return max(
340 3300 : HeatTransferUtils::reynolds(
341 4800 : _rho(r, t), _eps(r, t) * _speed(r, t), _d(r, t), max(mu(r, t), small_number)),
342 3000 : small_number);
343 : });
344 :
345 645 : addFunctorProperty<Real>(
346 : derivativePropertyNameFirst(NS::Reynolds, NS::pressure),
347 300 : [this, &Re, &mu, &drho_dp, &dmu_dp](const auto & r, const auto & t) -> Real
348 : {
349 300 : return NS::reynoldsPropertyDerivative(MetaPhysicL::raw_value(Re(r, t)),
350 300 : MetaPhysicL::raw_value(_rho(r, t)),
351 300 : MetaPhysicL::raw_value(mu(r, t)),
352 300 : drho_dp(r, t),
353 600 : dmu_dp(r, t));
354 : });
355 :
356 645 : addFunctorProperty<Real>(
357 : derivativePropertyNameFirst(NS::Reynolds, NS::T_fluid),
358 300 : [this, &Re, &mu, &drho_dT, &dmu_dT](const auto & r, const auto & t) -> Real
359 : {
360 300 : return NS::reynoldsPropertyDerivative(MetaPhysicL::raw_value(Re(r, t)),
361 300 : MetaPhysicL::raw_value(_rho(r, t)),
362 300 : MetaPhysicL::raw_value(mu(r, t)),
363 300 : drho_dT(r, t),
364 600 : dmu_dT(r, t));
365 : });
366 :
367 : // (hydraulic) Reynolds number
368 645 : addFunctorProperty<GenericReal<is_ad>>(
369 : NS::Reynolds_hydraulic,
370 300 : [this, &Re](const auto & r, const auto & t) -> GenericReal<is_ad>
371 : {
372 : static constexpr Real small_number = 1e-8;
373 : using std::max;
374 900 : return Re(r, t) / max(1 - _eps(r, t), small_number);
375 : });
376 :
377 : // (interstitial) Reynolds number
378 645 : addFunctorProperty<GenericReal<is_ad>>(
379 : NS::Reynolds_interstitial,
380 300 : [this, &Re](const auto & r, const auto & t) -> GenericReal<is_ad>
381 300 : { return Re(r, t) / _eps(r, t); });
382 : }
383 : else
384 : {
385 :
386 : const auto & rho =
387 0 : (!isParamValid(NS::density) || _force_define_density)
388 48 : ? addFunctorProperty<GenericReal<is_ad>>(
389 8 : _density_name,
390 7200 : [this](const auto & r, const auto & t) -> GenericReal<is_ad>
391 : {
392 7200 : auto total_pressure = _pressure(r, t) + _reference_pressure_value;
393 : // TODO: we should be integrating this term
394 7200 : const auto rho_approx = _fluid.rho_from_p_T(total_pressure, _T_fluid(r, t));
395 7200 : total_pressure +=
396 7200 : rho_approx * _gravity_vec * (r.getPoint() - _reference_pressure_point);
397 7200 : return _fluid.rho_from_p_T(total_pressure, _T_fluid(r, t));
398 : })
399 0 : : getFunctor<GenericReal<is_ad>>(_density_name);
400 :
401 24 : addFunctorProperty<GenericReal<is_ad>>(
402 : NS::cv,
403 150 : [this, &rho](const auto & r, const auto & t) -> GenericReal<is_ad>
404 : {
405 450 : const auto total_pressure =
406 150 : _pressure(r, t) + _reference_pressure_value +
407 150 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
408 150 : return _fluid.cv_from_p_T(total_pressure, _T_fluid(r, t));
409 : });
410 :
411 24 : const auto & cp = addFunctorProperty<GenericReal<is_ad>>(
412 8 : _specific_heat_name,
413 600 : [this, &rho](const auto & r, const auto & t) -> GenericReal<is_ad>
414 : {
415 1800 : const auto total_pressure =
416 600 : _pressure(r, t) + _reference_pressure_value +
417 600 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
418 600 : return _fluid.cp_from_p_T(total_pressure, _T_fluid(r, t));
419 : });
420 :
421 24 : const auto & mu = addFunctorProperty<GenericReal<is_ad>>(
422 8 : _dynamic_viscosity_name,
423 1650 : [this, &rho](const auto & r, const auto & t) -> GenericReal<is_ad>
424 : {
425 4950 : const auto total_pressure =
426 1650 : _pressure(r, t) + _reference_pressure_value +
427 1650 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
428 1650 : return _mu_rampdown(r, t) * _fluid.mu_from_p_T(total_pressure, _T_fluid(r, t));
429 : });
430 :
431 24 : const auto & k = addFunctorProperty<GenericReal<is_ad>>(
432 8 : _thermal_conductivity_name,
433 600 : [this, &rho](const auto & r, const auto & t) -> GenericReal<is_ad>
434 : {
435 1800 : const auto total_pressure =
436 600 : _pressure(r, t) + _reference_pressure_value +
437 600 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
438 600 : return _fluid.k_from_p_T(total_pressure, _T_fluid(r, t));
439 : });
440 :
441 : //
442 : // Time derivatives of fluid properties
443 : //
444 8 : if (_neglect_derivatives_of_density_time_derivative)
445 : {
446 32 : addFunctorProperty<GenericReal<is_ad>>(
447 8 : NS::time_deriv(_density_name),
448 150 : [this, &rho](const auto & r, const auto & t) -> GenericReal<is_ad>
449 : {
450 450 : const auto total_pressure =
451 150 : _pressure(r, t) + _reference_pressure_value +
452 150 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
453 : Real rho, drho_dp, drho_dT;
454 150 : _fluid.rho_from_p_T(MetaPhysicL::raw_value(total_pressure),
455 150 : MetaPhysicL::raw_value(_T_fluid(r, t)),
456 : rho,
457 : drho_dp,
458 : drho_dT);
459 300 : return drho_dp * _pressure.dot(r, t) + drho_dT * _T_fluid.dot(r, t);
460 : });
461 : }
462 : else
463 : {
464 0 : addFunctorProperty<GenericReal<is_ad>>(
465 0 : NS::time_deriv(_density_name),
466 0 : [this, &rho](const auto & r, const auto & t) -> GenericReal<is_ad>
467 : {
468 0 : const auto total_pressure =
469 0 : _pressure(r, t) + _reference_pressure_value +
470 0 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
471 : GenericReal<is_ad> rho, drho_dp, drho_dT;
472 0 : _fluid.rho_from_p_T(total_pressure, _T_fluid(r, t), rho, drho_dp, drho_dT);
473 0 : return drho_dp * _pressure.dot(r, t) + drho_dT * _T_fluid.dot(r, t);
474 : });
475 : }
476 :
477 32 : addFunctorProperty<GenericReal<is_ad>>(
478 8 : NS::time_deriv(_specific_heat_name),
479 150 : [this, &rho](const auto & r, const auto & t) -> GenericReal<is_ad>
480 : {
481 450 : const auto total_pressure =
482 150 : _pressure(r, t) + _reference_pressure_value +
483 300 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
484 : Real dcp_dp, dcp_dT, dummy;
485 : auto raw_pressure = MetaPhysicL::raw_value(total_pressure);
486 150 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
487 150 : _fluid.cp_from_p_T(raw_pressure, raw_T_fluid, dummy, dcp_dp, dcp_dT);
488 :
489 300 : return dcp_dp * _pressure.dot(r, t) + dcp_dT * _T_fluid.dot(r, t);
490 : });
491 :
492 : //
493 : // Temperature and pressure derivatives, to help with computing time derivatives
494 : //
495 :
496 32 : const auto & drho_dp = addFunctorProperty<Real>(
497 8 : derivativePropertyNameFirst(_density_name, NS::pressure),
498 300 : [this, &rho](const auto & r, const auto & t) -> Real
499 : {
500 900 : const auto total_pressure =
501 300 : _pressure(r, t) + _reference_pressure_value +
502 600 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
503 : Real drho_dp, drho_dT, dummy;
504 : auto raw_pressure = MetaPhysicL::raw_value(total_pressure);
505 300 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
506 :
507 300 : _fluid.rho_from_p_T(raw_pressure, raw_T_fluid, dummy, drho_dp, drho_dT);
508 300 : return drho_dp;
509 : });
510 :
511 32 : const auto & drho_dT = addFunctorProperty<Real>(
512 8 : derivativePropertyNameFirst(_density_name, NS::T_fluid),
513 300 : [this, &rho](const auto & r, const auto & t) -> Real
514 : {
515 900 : const auto total_pressure =
516 300 : _pressure(r, t) + _reference_pressure_value +
517 600 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
518 : Real drho_dp, drho_dT, dummy;
519 : auto raw_pressure = MetaPhysicL::raw_value(total_pressure);
520 300 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
521 :
522 300 : _fluid.rho_from_p_T(raw_pressure, raw_T_fluid, dummy, drho_dp, drho_dT);
523 300 : return drho_dT;
524 : });
525 :
526 32 : const auto & dcp_dp = addFunctorProperty<Real>(
527 8 : derivativePropertyNameFirst(_specific_heat_name, NS::pressure),
528 300 : [this, &rho](const auto & r, const auto & t) -> Real
529 : {
530 900 : const auto total_pressure =
531 300 : _pressure(r, t) + _reference_pressure_value +
532 600 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
533 : Real dcp_dp, dcp_dT, dummy;
534 : auto raw_pressure = MetaPhysicL::raw_value(total_pressure);
535 300 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
536 :
537 300 : _fluid.cp_from_p_T(raw_pressure, raw_T_fluid, dummy, dcp_dp, dcp_dT);
538 300 : return dcp_dp;
539 : });
540 :
541 32 : const auto & dcp_dT = addFunctorProperty<Real>(
542 8 : derivativePropertyNameFirst(_specific_heat_name, NS::T_fluid),
543 300 : [this, &rho](const auto & r, const auto & t) -> Real
544 : {
545 900 : const auto total_pressure =
546 300 : _pressure(r, t) + _reference_pressure_value +
547 600 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
548 : Real dcp_dp, dcp_dT, dummy;
549 : auto raw_pressure = MetaPhysicL::raw_value(total_pressure);
550 300 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
551 :
552 300 : _fluid.cp_from_p_T(raw_pressure, raw_T_fluid, dummy, dcp_dp, dcp_dT);
553 300 : return dcp_dT;
554 : });
555 :
556 32 : const auto & dmu_dp = addFunctorProperty<Real>(
557 8 : derivativePropertyNameFirst(_dynamic_viscosity_name, NS::pressure),
558 450 : [this, &rho](const auto & r, const auto & t) -> Real
559 : {
560 1350 : const auto total_pressure =
561 450 : _pressure(r, t) + _reference_pressure_value +
562 900 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
563 : Real dmu_dp, dmu_dT, dummy;
564 : auto raw_pressure = MetaPhysicL::raw_value(total_pressure);
565 450 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
566 :
567 450 : _fluid.mu_from_p_T(raw_pressure, raw_T_fluid, dummy, dmu_dp, dmu_dT);
568 450 : return _mu_rampdown(r, t) * dmu_dp;
569 : });
570 :
571 32 : const auto & dmu_dT = addFunctorProperty<Real>(
572 8 : derivativePropertyNameFirst(_dynamic_viscosity_name, NS::T_fluid),
573 450 : [this, &rho](const auto & r, const auto & t) -> Real
574 : {
575 1350 : const auto total_pressure =
576 450 : _pressure(r, t) + _reference_pressure_value +
577 900 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
578 : Real dmu_dp, dmu_dT, dummy;
579 : auto raw_pressure = MetaPhysicL::raw_value(total_pressure);
580 450 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
581 :
582 450 : _fluid.mu_from_p_T(raw_pressure, raw_T_fluid, dummy, dmu_dp, dmu_dT);
583 450 : return _mu_rampdown(r, t) * dmu_dT;
584 : });
585 :
586 32 : const auto & dk_dp = addFunctorProperty<Real>(
587 8 : derivativePropertyNameFirst(_thermal_conductivity_name, NS::pressure),
588 300 : [this, &rho](const auto & r, const auto & t) -> Real
589 : {
590 900 : const auto total_pressure =
591 300 : _pressure(r, t) + _reference_pressure_value +
592 600 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
593 : Real dk_dp, dk_dT, dummy;
594 : auto raw_pressure = MetaPhysicL::raw_value(total_pressure);
595 300 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
596 :
597 300 : _fluid.k_from_p_T(raw_pressure, raw_T_fluid, dummy, dk_dp, dk_dT);
598 300 : return dk_dp;
599 : });
600 :
601 32 : const auto & dk_dT = addFunctorProperty<Real>(
602 8 : derivativePropertyNameFirst(_thermal_conductivity_name, NS::T_fluid),
603 300 : [this, &rho](const auto & r, const auto & t) -> Real
604 : {
605 900 : const auto total_pressure =
606 300 : _pressure(r, t) + _reference_pressure_value +
607 600 : rho(r, t) * _gravity_vec * (r.getPoint() - _reference_pressure_point);
608 : Real dk_dp, dk_dT, dummy;
609 : auto raw_pressure = MetaPhysicL::raw_value(total_pressure);
610 300 : auto raw_T_fluid = MetaPhysicL::raw_value(_T_fluid(r, t));
611 :
612 300 : _fluid.k_from_p_T(raw_pressure, raw_T_fluid, dummy, dk_dp, dk_dT);
613 300 : return dk_dT;
614 : });
615 :
616 : //
617 : // Fluid adimensional quantities, used in numerous correlations
618 : //
619 :
620 24 : addFunctorProperty<GenericReal<is_ad>>(
621 : NS::Prandtl,
622 150 : [&cp, &mu, &k](const auto & r, const auto & t) -> GenericReal<is_ad>
623 : {
624 : static constexpr Real small_number = 1e-8;
625 : using std::max;
626 :
627 150 : return HeatTransferUtils::prandtl(cp(r, t), mu(r, t), max(k(r, t), small_number));
628 : });
629 :
630 24 : addFunctorProperty<Real>(
631 : derivativePropertyNameFirst(NS::Prandtl, NS::pressure),
632 150 : [&mu, &cp, &k, &dmu_dp, &dcp_dp, &dk_dp](const auto & r, const auto & t) -> Real
633 : {
634 150 : return NS::prandtlPropertyDerivative(MetaPhysicL::raw_value(mu(r, t)),
635 150 : MetaPhysicL::raw_value(cp(r, t)),
636 150 : MetaPhysicL::raw_value(k(r, t)),
637 150 : dmu_dp(r, t),
638 150 : dcp_dp(r, t),
639 300 : dk_dp(r, t));
640 : });
641 :
642 24 : addFunctorProperty<Real>(
643 : derivativePropertyNameFirst(NS::Prandtl, NS::T_fluid),
644 150 : [&mu, &cp, &k, &dmu_dT, &dcp_dT, &dk_dT](const auto & r, const auto & t) -> Real
645 : {
646 150 : return NS::prandtlPropertyDerivative(MetaPhysicL::raw_value(mu(r, t)),
647 150 : MetaPhysicL::raw_value(cp(r, t)),
648 150 : MetaPhysicL::raw_value(k(r, t)),
649 150 : dmu_dT(r, t),
650 150 : dcp_dT(r, t),
651 300 : dk_dT(r, t));
652 : });
653 :
654 : //
655 : // (pore / particle) Reynolds number based on superficial velocity and
656 : // characteristic length. Only call Reynolds() one time to compute all three so that
657 : // we don't redundantly check that viscosity is not too close to zero.
658 : //
659 :
660 24 : const auto & Re = addFunctorProperty<GenericReal<is_ad>>(
661 : NS::Reynolds,
662 750 : [this, &mu](const auto & r, const auto & t) -> GenericReal<is_ad>
663 : {
664 : static constexpr Real small_number = 1e-8;
665 : using std::max;
666 :
667 : return max(
668 0 : HeatTransferUtils::reynolds(
669 750 : _rho(r, t), _eps(r, t) * _speed(r, t), _d(r, t), max(mu(r, t), small_number)),
670 1500 : small_number);
671 : });
672 :
673 24 : addFunctorProperty<Real>(
674 : derivativePropertyNameFirst(NS::Reynolds, NS::pressure),
675 150 : [this, &Re, &mu, &drho_dp, &dmu_dp](const auto & r, const auto & t) -> Real
676 : {
677 150 : return NS::reynoldsPropertyDerivative(MetaPhysicL::raw_value(Re(r, t)),
678 150 : MetaPhysicL::raw_value(_rho(r, t)),
679 150 : MetaPhysicL::raw_value(mu(r, t)),
680 150 : drho_dp(r, t),
681 300 : dmu_dp(r, t));
682 : });
683 :
684 24 : addFunctorProperty<Real>(
685 : derivativePropertyNameFirst(NS::Reynolds, NS::T_fluid),
686 150 : [this, &Re, &mu, &drho_dT, &dmu_dT](const auto & r, const auto & t) -> Real
687 : {
688 150 : return NS::reynoldsPropertyDerivative(MetaPhysicL::raw_value(Re(r, t)),
689 150 : MetaPhysicL::raw_value(_rho(r, t)),
690 150 : MetaPhysicL::raw_value(mu(r, t)),
691 150 : drho_dT(r, t),
692 300 : dmu_dT(r, t));
693 : });
694 :
695 : // (hydraulic) Reynolds number
696 24 : addFunctorProperty<GenericReal<is_ad>>(
697 : NS::Reynolds_hydraulic,
698 150 : [this, &Re](const auto & r, const auto & t) -> GenericReal<is_ad>
699 : {
700 : static constexpr Real small_number = 1e-8;
701 : using std::max;
702 450 : return Re(r, t) / max(1 - _eps(r, t), small_number);
703 : });
704 :
705 : // (interstitial) Reynolds number
706 24 : addFunctorProperty<GenericReal<is_ad>>(
707 : NS::Reynolds_interstitial,
708 150 : [this, &Re](const auto & r, const auto & t) -> GenericReal<is_ad>
709 150 : { return Re(r, t) / _eps(r, t); });
710 : }
711 5342 : }
712 :
713 : template class GeneralFunctorFluidPropsTempl<false>;
714 : template class GeneralFunctorFluidPropsTempl<true>;
|
