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