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 "FunctorThermalResistanceBC.h"
11 : #include "HeatConductionNames.h"
12 :
13 : registerMooseObject("HeatTransferApp", FunctorThermalResistanceBC);
14 :
15 : InputParameters
16 76 : FunctorThermalResistanceBC::validParams()
17 : {
18 76 : auto params = FVFluxBC::validParams();
19 152 : params.addParam<MooseFunctorName>("temperature", "temperature variable");
20 76 : params.addRequiredParam<Real>(HeatConduction::T_ambient, "constant ambient temperature");
21 152 : params.addRequiredParam<MooseFunctorName>("htc", "heat transfer coefficient");
22 :
23 152 : params.addRequiredRangeCheckedParam<Real>(HeatConduction::emissivity,
24 76 : HeatConduction::emissivity + " >= 0.0 & " +
25 76 : HeatConduction::emissivity + " <= 1.0",
26 : "emissivity of the surface");
27 :
28 152 : params.addRequiredParam<std::vector<Real>>(
29 : "thermal_conductivities",
30 : "vector of thermal conductivity values used for the conduction layers");
31 152 : params.addRequiredParam<std::vector<Real>>("conduction_thicknesses",
32 : "vector of conduction layer thicknesses");
33 :
34 152 : MooseEnum geometry("cartesian cylindrical", "cartesian");
35 152 : params.addParam<MooseEnum>("geometry", geometry, "type of geometry");
36 152 : params.addRangeCheckedParam<Real>("inner_radius",
37 : "inner_radius > 0.0",
38 : "coordinate corresponding to the first resistance layer");
39 :
40 228 : params.addRangeCheckedParam<Real>(
41 152 : "step_size", 0.1, "step_size > 0.0", "underrelaxation step size");
42 :
43 228 : params.addRangeCheckedParam<unsigned int>(
44 152 : "max_iterations", 100, "max_iterations >= 0", "maximum iterations");
45 :
46 228 : params.addRangeCheckedParam<Real>(
47 152 : "tolerance", 1E-3, "tolerance > 0.0", "tolerance to converge iterations");
48 76 : params.addClassDescription("Thermal resistance heat flux boundary condition for the "
49 : "fluid and solid energy equations");
50 76 : return params;
51 76 : }
52 :
53 40 : FunctorThermalResistanceBC::FunctorThermalResistanceBC(const InputParameters & parameters)
54 : : FVFluxBC(parameters),
55 40 : _geometry(getParam<MooseEnum>("geometry").getEnum<Moose::CoordinateSystemType>()),
56 40 : _inner_radius(
57 50 : _geometry == Moose::CoordinateSystemType::COORD_RZ ? getParam<Real>("inner_radius") : 1.0),
58 120 : _T(isParamValid("temperature") ? getFunctor<ADReal>("temperature")
59 80 : : getFunctor<ADReal>("variable")),
60 40 : _T_ambient(getParam<Real>(HeatConduction::T_ambient)),
61 80 : _k(getParam<std::vector<Real>>("thermal_conductivities")),
62 80 : _dx(getParam<std::vector<Real>>("conduction_thicknesses")),
63 80 : _h(getFunctor<ADReal>(getParam<MooseFunctorName>("htc"))),
64 40 : _emissivity(getParam<Real>(HeatConduction::emissivity)),
65 80 : _max_iterations(getParam<unsigned int>("max_iterations")),
66 80 : _tolerance(getParam<Real>("tolerance")),
67 80 : _alpha(getParam<Real>("step_size")),
68 40 : _T_surface(0.0),
69 40 : _outer_radius(_inner_radius),
70 40 : _conduction_resistance(0.0),
71 80 : _parallel_resistance(0.0)
72 : {
73 40 : if (_k.size() != _dx.size())
74 0 : paramError("conduction_thicknesses",
75 : "Number of specified thermal conductivities must match "
76 : "the number of conduction layers!");
77 :
78 40 : if (_geometry == Moose::CoordinateSystemType::COORD_RZ)
79 40 : for (const auto & d : _dx)
80 : _outer_radius += d;
81 :
82 : // because the thermal conductivities are constant, we only need to compute
83 : // the conduction resistance one time
84 40 : computeConductionResistance();
85 40 : }
86 :
87 : void
88 40 : FunctorThermalResistanceBC::computeConductionResistance()
89 : {
90 40 : Real r = _inner_radius;
91 :
92 160 : for (const auto i : index_range(_k))
93 : {
94 120 : switch (_geometry)
95 : {
96 90 : case Moose::CoordinateSystemType::COORD_XYZ:
97 90 : _conduction_resistance += _dx[i] / _k[i];
98 90 : break;
99 30 : case Moose::CoordinateSystemType::COORD_RZ:
100 : {
101 30 : _conduction_resistance += std::log((_dx[i] + r) / r) / _k[i];
102 30 : r += _dx[i];
103 30 : break;
104 : }
105 0 : default:
106 0 : mooseError("Unhandled 'GeometryEnum' in 'FunctorThermalResistanceBC'!");
107 : }
108 : }
109 40 : }
110 :
111 : ADReal
112 600 : FunctorThermalResistanceBC::computeQpResidual()
113 : {
114 : using std::abs;
115 :
116 : // Evaluate material properties on the face
117 600 : const auto face_arg = singleSidedFaceArg();
118 :
119 : // radiation resistance has to be solved iteratively, since we don't know the
120 : // surface temperature. We do know that the heat flux in the conduction layers
121 : // must match the heat flux in the parallel convection-radiation segment. For a
122 : // first guess, take the surface temperature as the average of _T and T_ambient.
123 1800 : _T_surface = 0.5 * (_T(face_arg, determineState()) + _T_ambient);
124 :
125 : // total flux perpendicular to boundary
126 : ADReal flux;
127 :
128 : // resistance component representing the sum of the convection and radiation
129 : // resistances, in parallel
130 600 : computeParallelResistance();
131 :
132 : // other iteration requirements
133 : unsigned int iteration = 0;
134 600 : ADReal norm = 2 * _tolerance;
135 : ADReal T_surface_previous;
136 :
137 : // iterate to find the approximate surface temperature needed for evaluating the
138 : // radiation resistance. We only do this iteration if we have radiation transfer.
139 600 : if (_emissivity > 1e-8)
140 36300 : while (norm > (_tolerance * _alpha))
141 : {
142 36300 : T_surface_previous = _T_surface;
143 :
144 : // compute the flux based on the conduction part of the circuit
145 72600 : flux = (_T(face_arg, determineState()) - _T_surface) / _conduction_resistance;
146 :
147 36300 : computeParallelResistance();
148 :
149 : // use the flux computed from the conduction half to update T_surface
150 72600 : _T_surface = flux * _parallel_resistance + _T_ambient;
151 108900 : _T_surface = _alpha * _T_surface + (1 - _alpha) * T_surface_previous;
152 108900 : norm = abs(_T_surface - T_surface_previous) / abs(T_surface_previous);
153 :
154 36300 : if (iteration == _max_iterations)
155 : {
156 300 : mooseWarning("Maximum number of iterations reached in 'FunctorThermalResistanceBC'!");
157 : break;
158 : }
159 : else
160 36000 : iteration += 1;
161 : }
162 :
163 : // once we have determined T_surface, we can finally evaluate the complete
164 : // resistance to find the overall heat flux. For Cartesian, dividing by the
165 : // 'inner_radius' has no effect, but it is required for correct normalization
166 : // for cylindrical geometries.
167 600 : flux = (_T(face_arg, determineState()) - _T_ambient) /
168 1200 : (_conduction_resistance + _parallel_resistance) / _inner_radius;
169 600 : return flux;
170 : }
171 :
172 : void
173 36900 : FunctorThermalResistanceBC::computeParallelResistance()
174 : {
175 36900 : const auto face_arg = singleSidedFaceArg();
176 :
177 : // compute the parallel convection and radiation resistances, assuming they
178 : // act on the same surface area size
179 36900 : ADReal hr = _emissivity * HeatConduction::Constants::sigma *
180 36900 : (_T_surface * _T_surface + _T_ambient * _T_ambient) * (_T_surface + _T_ambient);
181 :
182 : // for Cartesian, dividing by the 'outer_radius' has no effect, but it is
183 : // required for correct normalization for cylindrical geometries
184 110700 : _parallel_resistance = 1.0 / (hr + _h(face_arg, determineState())) / _outer_radius;
185 36900 : }
|
