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 : // Moose includes
11 : #include "CHTHandler.h"
12 : #include "LinearFVFluxKernel.h"
13 : #include "LinearFVDiffusion.h"
14 : #include "LinearFVAnisotropicDiffusion.h"
15 : #include "LinearFVBoundaryCondition.h"
16 : #include "LinearFVCHTBCInterface.h"
17 : #include "FEProblemBase.h"
18 :
19 : namespace NS
20 : {
21 : namespace FV
22 : {
23 :
24 : InputParameters
25 1344 : CHTHandler::validParams()
26 : {
27 1344 : auto params = emptyInputParameters();
28 2688 : params.addParam<std::vector<BoundaryName>>(
29 : "cht_interfaces",
30 : {},
31 : "The interfaces where we would like to add conjugate heat transfer handling.");
32 :
33 4032 : params.addRangeCheckedParam<unsigned int>(
34 : "max_cht_fpi",
35 2688 : 1,
36 : "max_cht_fpi >= 1",
37 : "Number of maximum fixed point iterations (FPI). Currently only applied to"
38 : " conjugate heat transfer simulations. The default value of 1 essentially keeps"
39 : " the FPI feature turned off. CHT iteration ends after this number of iteration even if the "
40 : "tolerance is not met.");
41 :
42 4032 : params.addRangeCheckedParam<Real>(
43 : "cht_heat_flux_tolerance",
44 2688 : 1e-5,
45 : "cht_heat_flux_tolerance > 0 & cht_heat_flux_tolerance <= 1.0",
46 : "The relative tolerance for terminating conjugate heat transfer iteration before the maximum "
47 : "number of CHT iterations. Relative tolerance is ignore if the maximum number of CHT "
48 : "iterations is reached.");
49 :
50 2688 : params.addParam<std::vector<Real>>(
51 : "cht_fluid_temperature_relaxation",
52 : {},
53 : "The relaxation factors for the boundary temperature when being updated on the fluid side.");
54 2688 : params.addParam<std::vector<Real>>(
55 : "cht_solid_temperature_relaxation",
56 : {},
57 : "The relaxation factors for the boundary temperature when being updated on the solid side.");
58 2688 : params.addParam<std::vector<Real>>(
59 : "cht_fluid_flux_relaxation",
60 : {},
61 : "The relaxation factors for the boundary flux when being updated on the fluid side.");
62 2688 : params.addParam<std::vector<Real>>(
63 : "cht_solid_flux_relaxation",
64 : {},
65 : "The relaxation factors for the boundary flux when being updated on the solid side.");
66 :
67 2688 : params.addParamNamesToGroup(
68 : "cht_interfaces max_cht_fpi cht_heat_flux_tolerance cht_fluid_temperature_relaxation "
69 : "cht_solid_temperature_relaxation cht_fluid_flux_relaxation "
70 : "cht_solid_flux_relaxation",
71 : "Conjugate Heat Transfer");
72 :
73 1344 : return params;
74 0 : }
75 :
76 672 : CHTHandler::CHTHandler(const InputParameters & params)
77 : : MooseObject(params),
78 1344 : _problem(*getCheckedPointerParam<FEProblemBase *>(
79 : "_fe_problem_base", "This might happen if you don't have a mesh")),
80 672 : _mesh(_problem.mesh()),
81 1344 : _cht_boundary_names(getParam<std::vector<BoundaryName>>("cht_interfaces")),
82 672 : _cht_boundary_ids(_mesh.getBoundaryIDs(_cht_boundary_names)),
83 1344 : _max_cht_fpi(getParam<unsigned int>("max_cht_fpi")),
84 2016 : _cht_heat_flux_tolerance(getParam<Real>("cht_heat_flux_tolerance"))
85 : {
86 2016 : if (isParamSetByUser("cht_interfaces") && !_cht_boundary_names.size())
87 0 : paramError("cht_interfaces", "You must declare at least one interface!");
88 672 : }
89 :
90 : void
91 40 : CHTHandler::linkEnergySystems(SystemBase * solid_energy_system,
92 : SystemBase * fluid_energy_system,
93 : std::vector<SystemBase *> pm_radiation_systems)
94 : {
95 40 : _energy_system = fluid_energy_system;
96 40 : _solid_energy_system = solid_energy_system;
97 40 : _pm_radiation_systems = pm_radiation_systems;
98 :
99 40 : if (!_energy_system || !_solid_energy_system)
100 0 : paramError("cht_interfaces",
101 : "You selected to do conjugate heat transfer treatment, but it needs two energy "
102 : "systems: a solid and a fluid. One of these systems is missing.");
103 40 : }
104 :
105 : void
106 40 : CHTHandler::deduceCHTBoundaryCoupling()
107 : {
108 40 : if (_solid_energy_system->nVariables() != 1)
109 0 : mooseError("We should have only one variable in the solid energy system: ",
110 0 : _energy_system->name(),
111 : "! Right now we have: ",
112 0 : Moose::stringify(_solid_energy_system->getVariableNames()));
113 40 : if (_energy_system->nVariables() != 1)
114 0 : mooseError("We should have only one variable in the fluid energy system: ",
115 0 : _energy_system->name(),
116 : "! Right now we have: ",
117 0 : Moose::stringify(_energy_system->getVariableNames()));
118 120 : const std::vector<std::string> solid_fluid({"solid", "fluid"});
119 :
120 : // We do some setup at the beginning to make sure the container sizes are good
121 : _cht_system_numbers =
122 40 : std::vector<unsigned int>({_solid_energy_system->number(), _energy_system->number()});
123 40 : _cht_conduction_kernels = std::vector<LinearFVFluxKernel *>({nullptr, nullptr});
124 40 : _cht_boundary_conditions.clear();
125 40 : _cht_boundary_conditions.resize(_cht_boundary_names.size(), {nullptr, nullptr});
126 :
127 : // Populate the PM radiation system numbers
128 40 : if (!_pm_radiation_systems.empty())
129 : {
130 32 : for (const auto sys_i : index_range(_pm_radiation_systems))
131 16 : _cht_pm_radiation_system_numbers.push_back(_pm_radiation_systems[sys_i]->number());
132 :
133 : // Reserve space for _cht_pm_radiation_kernels based on the size of
134 : // _cht_pm_radiation_system_numbers
135 16 : _cht_pm_radiation_kernels.reserve(_cht_pm_radiation_system_numbers.size());
136 : // Reserve space for pm radiation boundary conditions
137 16 : _cht_pm_radiation_boundary_conditions.clear();
138 16 : _cht_pm_radiation_boundary_conditions.resize(
139 : _cht_boundary_names.size(),
140 32 : std::vector<LinearFVBoundaryCondition *>(_cht_pm_radiation_system_numbers.size(), nullptr));
141 : }
142 :
143 : const auto flux_relaxation_param_names =
144 120 : std::vector<std::string>({"cht_solid_flux_relaxation", "cht_fluid_flux_relaxation"});
145 : const auto temperature_relaxation_param_names = std::vector<std::string>(
146 120 : {"cht_solid_temperature_relaxation", "cht_fluid_temperature_relaxation"});
147 40 : _cht_flux_relaxation_factor.clear();
148 40 : _cht_flux_relaxation_factor.resize(2, std::vector<Real>(_cht_boundary_names.size(), 1.0));
149 40 : _cht_temperature_relaxation_factor.clear();
150 40 : _cht_temperature_relaxation_factor.resize(2, std::vector<Real>(_cht_boundary_names.size(), 1.0));
151 :
152 120 : for (const auto region_index : index_range(solid_fluid))
153 : {
154 : // First thing, we fetch the relaxation parameter values
155 : const auto & flux_param_value =
156 : getParam<std::vector<Real>>(flux_relaxation_param_names[region_index]);
157 80 : if (flux_param_value.empty() || (flux_param_value.size() != _cht_boundary_names.size()))
158 0 : paramError(flux_relaxation_param_names[region_index],
159 : "The number of relaxation factors is not the same as the number of interfaces!");
160 :
161 80 : _cht_flux_relaxation_factor[region_index] = flux_param_value;
162 : // We have to do the range check here because the intput parameter check errors if the vector is
163 : // empty
164 160 : for (const auto param : _cht_flux_relaxation_factor[region_index])
165 80 : if (param <= 0 || param > 1.0)
166 0 : paramError(flux_relaxation_param_names[region_index],
167 : "The relaxation parameter should be between 0 and 1!");
168 :
169 : const auto & temperature_param_value =
170 : getParam<std::vector<Real>>(temperature_relaxation_param_names[region_index]);
171 80 : if (temperature_param_value.empty() ||
172 : (temperature_param_value.size() != _cht_boundary_names.size()))
173 0 : paramError(temperature_relaxation_param_names[region_index],
174 : "The number of relaxation factors is not the same as the number of interfaces!");
175 :
176 80 : _cht_temperature_relaxation_factor[region_index] = temperature_param_value;
177 : // We have to do the range check here because the intput parameter check errors if the vector is
178 : // empty
179 160 : for (const auto param : _cht_temperature_relaxation_factor[region_index])
180 80 : if (param <= 0 || param > 1.0)
181 0 : paramError(temperature_relaxation_param_names[region_index],
182 : "The relaxation parameter should be between 0 and 1!");
183 :
184 : // We then fetch the conduction kernels
185 : std::vector<LinearFVFluxKernel *> flux_kernels;
186 80 : _app.theWarehouse()
187 80 : .query()
188 80 : .template condition<AttribSystem>("LinearFVFluxKernel")
189 160 : .template condition<AttribVar>(0)
190 80 : .template condition<AttribSysNum>(_cht_system_numbers[region_index])
191 : .queryInto(flux_kernels);
192 :
193 : // We then fetch the radiation conduction kernels in the fluid region
194 80 : if (!_pm_radiation_systems.empty() && region_index == 1)
195 32 : for (const auto sys_i : index_range(_pm_radiation_systems))
196 : {
197 : // We then fetch the radiation conduction kernels
198 : std::vector<LinearFVFluxKernel *> radiation_kernels;
199 32 : _app.theWarehouse()
200 16 : .query()
201 16 : .template condition<AttribSystem>("LinearFVFluxKernel")
202 32 : .template condition<AttribVar>(0)
203 16 : .template condition<AttribSysNum>(_cht_pm_radiation_system_numbers[sys_i])
204 : .queryInto(radiation_kernels);
205 :
206 16 : if (radiation_kernels.size() > 1)
207 0 : mooseError(
208 : "We already have a kernel that describes the participating media radiation diffusion "
209 : "with the name: ",
210 0 : radiation_kernels[0]->name(),
211 : ". Make sure that you have only one conduction kernel.");
212 16 : else if (radiation_kernels.empty())
213 0 : mooseError("We did not find a diffusion kernel for the participating media radiation "
214 : "diffusion to compute the "
215 : "radiative heat flux. Please add a diffusion kernel.");
216 : else
217 16 : _cht_pm_radiation_kernels.push_back(radiation_kernels[0]);
218 16 : }
219 :
220 184 : for (auto kernel : flux_kernels)
221 : {
222 104 : auto check_diff = dynamic_cast<LinearFVDiffusion *>(kernel);
223 104 : auto check_aniso_diff = dynamic_cast<LinearFVAnisotropicDiffusion *>(kernel);
224 104 : if (_cht_conduction_kernels[region_index] && (check_diff || check_aniso_diff))
225 0 : mooseError("We already have a kernel that describes the heat conduction for the ",
226 : solid_fluid[region_index],
227 : " domain: ",
228 : _cht_conduction_kernels[region_index]->name(),
229 : " We found another one with the name: ",
230 : (check_diff ? check_diff->name() : check_aniso_diff->name()),
231 : " Make sure that you have only one conduction kernel on the ",
232 : solid_fluid[region_index],
233 : " side!");
234 :
235 104 : if (check_diff || check_aniso_diff)
236 80 : _cht_conduction_kernels[region_index] = kernel;
237 : }
238 :
239 : // Then we check the boundary conditions, to make sure at least there is something defined
240 : // from both sides
241 160 : for (const auto bd_index : index_range(_cht_boundary_names))
242 : {
243 : const auto & boundary_name = _cht_boundary_names[bd_index];
244 80 : const auto boundary_id = _cht_boundary_ids[bd_index];
245 :
246 : std::vector<LinearFVBoundaryCondition *> bcs;
247 80 : _problem.getMooseApp()
248 : .theWarehouse()
249 80 : .query()
250 80 : .template condition<AttribSystem>("LinearFVBoundaryCondition")
251 160 : .template condition<AttribVar>(0)
252 80 : .template condition<AttribSysNum>(_cht_system_numbers[region_index])
253 80 : .template condition<AttribBoundaries>(boundary_id)
254 : .queryInto(bcs);
255 :
256 : // We then fetch the radiation conduction bcs in the fluid region (i.e MarshakBC in P1)
257 80 : if (!_pm_radiation_systems.empty() && region_index == 1)
258 32 : for (const auto sys_i : index_range(_pm_radiation_systems))
259 : {
260 : std::vector<LinearFVBoundaryCondition *> rad_bcs;
261 16 : _app.theWarehouse()
262 16 : .query()
263 16 : .template condition<AttribSystem>("LinearFVBoundaryCondition")
264 32 : .template condition<AttribVar>(0)
265 16 : .template condition<AttribSysNum>(_cht_pm_radiation_system_numbers[sys_i])
266 16 : .template condition<AttribBoundaries>(boundary_id)
267 : .queryInto(rad_bcs);
268 :
269 16 : if (!rad_bcs.empty())
270 16 : _cht_pm_radiation_boundary_conditions[bd_index][sys_i] = rad_bcs[0];
271 : else
272 0 : mooseError("No LinearFVBoundaryCondition found for the given boundary or system.");
273 16 : }
274 :
275 80 : if (bcs.size() != 1)
276 0 : mooseError("We found multiple or no boundary conditions for solid energy on boundary ",
277 : boundary_name,
278 : " (ID: ",
279 : boundary_id,
280 : "). Make sure you define exactly one for conjugate heat transfer applications!");
281 80 : _cht_boundary_conditions[bd_index][region_index] = bcs[0];
282 :
283 80 : if (!dynamic_cast<LinearFVCHTBCInterface *>(_cht_boundary_conditions[bd_index][region_index]))
284 0 : mooseError("The selected boundary condition cannot be used with CHT problems! Make sure it "
285 : "inherits from LinearFVCHTBCInterface!");
286 80 : }
287 80 : }
288 240 : }
289 :
290 : void
291 40 : CHTHandler::setupConjugateHeatTransferContainers()
292 : {
293 : // We already error in initialSetup if we have more variables
294 : const auto * fluid_variable =
295 40 : dynamic_cast<const MooseLinearVariableFVReal *>(&_energy_system->getVariable(0, 0));
296 : const auto * solid_variable =
297 40 : dynamic_cast<const MooseLinearVariableFVReal *>(&_solid_energy_system->getVariable(0, 0));
298 :
299 40 : _cht_face_info.clear();
300 40 : _cht_face_info.resize(_cht_boundary_ids.size());
301 40 : _boundary_heat_flux.clear();
302 40 : _boundary_temperature.clear();
303 40 : _integrated_boundary_heat_flux.clear();
304 :
305 80 : for (const auto bd_index : index_range(_cht_boundary_ids))
306 : {
307 40 : const auto bd_id = _cht_boundary_ids[bd_index];
308 : const auto & bd_name = _cht_boundary_names[bd_index];
309 :
310 : // We populate the face infos for every interface
311 : auto & bd_fi_container = _cht_face_info[bd_index];
312 49734 : for (auto & fi : _problem.mesh().faceInfo())
313 49694 : if (fi->boundaryIDs().count(bd_id))
314 400 : bd_fi_container.push_back(fi);
315 :
316 : // We do this because the coupling functors should be evaluated on both sides
317 : // of the interface and there are rigorous checks if the functors don't support a subdomain
318 : std::set<SubdomainID> combined_set;
319 40 : std::set_union(solid_variable->blockIDs().begin(),
320 40 : solid_variable->blockIDs().end(),
321 40 : fluid_variable->blockIDs().begin(),
322 40 : fluid_variable->blockIDs().end(),
323 : std::inserter(combined_set, combined_set.begin()));
324 :
325 : // We instantiate the coupling fuctors for heat flux and temperature
326 : FaceCenteredMapFunctor<Real, std::unordered_map<dof_id_type, Real>> solid_bd_flux(
327 40 : _problem.mesh(), combined_set, "heat_flux_to_solid_" + bd_name);
328 : FaceCenteredMapFunctor<Real, std::unordered_map<dof_id_type, Real>> fluid_bd_flux(
329 40 : _problem.mesh(), combined_set, "heat_flux_to_fluid_" + bd_name);
330 :
331 : _boundary_heat_flux.push_back(
332 160 : std::vector<FaceCenteredMapFunctor<Real, std::unordered_map<dof_id_type, Real>>>(
333 80 : {std::move(solid_bd_flux), std::move(fluid_bd_flux)}));
334 : auto & flux_container = _boundary_heat_flux.back();
335 :
336 80 : _integrated_boundary_heat_flux.push_back(std::vector<Real>({0.0, 0.0}));
337 :
338 : FaceCenteredMapFunctor<Real, std::unordered_map<dof_id_type, Real>> solid_bd_temperature(
339 40 : _problem.mesh(), combined_set, "interface_temperature_solid_" + bd_name);
340 : FaceCenteredMapFunctor<Real, std::unordered_map<dof_id_type, Real>> fluid_bd_temperature(
341 40 : _problem.mesh(), combined_set, "interface_temperature_fluid_" + bd_name);
342 :
343 : _boundary_temperature.push_back(
344 160 : std::vector<FaceCenteredMapFunctor<Real, std::unordered_map<dof_id_type, Real>>>(
345 80 : {std::move(solid_bd_temperature), std::move(fluid_bd_temperature)}));
346 : auto & temperature_container = _boundary_temperature.back();
347 :
348 : // Time to register the functors on all of the threads
349 80 : for (const auto tid : make_range(libMesh::n_threads()))
350 : {
351 40 : _problem.addFunctor("heat_flux_to_solid_" + bd_name, flux_container[NS::CHTSide::SOLID], tid);
352 40 : _problem.addFunctor("heat_flux_to_fluid_" + bd_name, flux_container[NS::CHTSide::FLUID], tid);
353 40 : _problem.addFunctor(
354 40 : "interface_temperature_solid_" + bd_name, temperature_container[NS::CHTSide::SOLID], tid);
355 40 : _problem.addFunctor(
356 80 : "interface_temperature_fluid_" + bd_name, temperature_container[NS::CHTSide::FLUID], tid);
357 : }
358 :
359 : // Initialize the containers, they will be filled with correct values soon.
360 : // Before any solve happens.
361 120 : for (const auto region_index : make_range(2))
362 880 : for (auto & fi : bd_fi_container)
363 : {
364 800 : flux_container[region_index][fi->id()] = 0.0;
365 800 : temperature_container[region_index][fi->id()] = 0.0;
366 : }
367 40 : }
368 120 : }
369 :
370 : void
371 38 : CHTHandler::initializeCHTCouplingFields()
372 : {
373 76 : for (const auto bd_index : index_range(_cht_boundary_ids))
374 : {
375 : const auto & bd_fi_container = _cht_face_info[bd_index];
376 : auto & temperature_container = _boundary_temperature[bd_index];
377 :
378 114 : for (const auto region_index : make_range(2))
379 : {
380 : // Can't be const considering we will update members from here
381 76 : auto bc = _cht_boundary_conditions[bd_index][region_index];
382 872 : for (const auto & fi : bd_fi_container)
383 : {
384 796 : bc->setupFaceData(fi, fi->faceType(std::make_pair(0, _cht_system_numbers[region_index])));
385 796 : temperature_container[1 - region_index][fi->id()] = bc->computeBoundaryValue();
386 : }
387 : }
388 : }
389 38 : }
390 :
391 : void
392 14552 : CHTHandler::updateCHTBoundaryCouplingFields(const NS::CHTSide side)
393 : {
394 : // Well we can just use the face that this enum casts into int very nicely
395 : // we can use it to get the index of the other side
396 14552 : const NS::CHTSide other_side = static_cast<NS::CHTSide>(1 - side);
397 :
398 29104 : for (const auto bd_index : index_range(_cht_boundary_ids))
399 : {
400 14552 : auto & other_bc = _cht_boundary_conditions[bd_index][other_side];
401 : auto & other_kernel = _cht_conduction_kernels[other_side];
402 :
403 : // We get the relaxation from the other side, so if we are fluid side we get the solid
404 : // relaxation
405 14552 : const auto temperature_relaxation = _cht_flux_relaxation_factor[other_side][bd_index];
406 14552 : const auto flux_relaxation = _cht_temperature_relaxation_factor[other_side][bd_index];
407 :
408 : // Fetching the right container here, if side is fluid we fetch "heat_flux_to_fluid"
409 14552 : auto & flux_container = _boundary_heat_flux[bd_index][side];
410 : // Fetching the other side's contaienr here, if side is fluid we fetch the solid temperature
411 : auto & temperature_container = _boundary_temperature[bd_index][other_side];
412 : // We will also update the integrated flux for output info
413 : auto & integrated_flux = _integrated_boundary_heat_flux[bd_index][side];
414 : // We are recomputing this so, time to zero this out
415 14552 : integrated_flux = 0.0;
416 :
417 : const auto & bd_fi_container = _cht_face_info[bd_index];
418 :
419 : // We enter the face loop to update the coupling fields
420 181024 : for (const auto & fi : bd_fi_container)
421 : {
422 166472 : other_kernel->setupFaceData(fi);
423 : // We will want the flux in W/m2 for the coupling so no face integral for now,
424 : // this can cause issues if we start using face area in the kernels
425 : // for more than just face integral multipliers.
426 : // Also, if we decide to not require overlapping meshes on the boundary
427 : // this will probably have to change.
428 166472 : other_kernel->setCurrentFaceArea(1.0);
429 166472 : other_bc->setupFaceData(fi, fi->faceType(std::make_pair(0, _cht_system_numbers[other_side])));
430 :
431 : // T_new = relaxation * T_boundary + (1-relaxation) * T_old
432 332944 : temperature_container[fi->id()] =
433 166472 : temperature_relaxation * other_bc->computeBoundaryValue() +
434 166472 : (1 - temperature_relaxation) * temperature_container[fi->id()];
435 :
436 : // Flux_new = relaxation * Flux_boundary + (1-relaxation) * Flux_old,
437 : // minus sign is due to the normal differences
438 :
439 : // Conductive flux
440 166472 : auto flux = other_kernel->computeBoundaryFlux(*other_bc);
441 :
442 : // If participating media radiation system exists we add the heat flux from the fluid
443 : // to the solid region.
444 166472 : if (!_pm_radiation_systems.empty() && side == NS::CHTSide::SOLID)
445 4424 : for (const auto sys_i : index_range(_pm_radiation_systems))
446 : {
447 2212 : _cht_pm_radiation_kernels[sys_i]->setupFaceData(fi);
448 2212 : _cht_pm_radiation_kernels[sys_i]->setCurrentFaceArea(1.0);
449 2212 : _cht_pm_radiation_boundary_conditions[bd_index][sys_i]->setupFaceData(
450 : fi, fi->faceType(std::make_pair(0, _cht_pm_radiation_system_numbers[sys_i])));
451 2212 : flux += _cht_pm_radiation_kernels[sys_i]->computeBoundaryFlux(
452 : *_cht_pm_radiation_boundary_conditions[bd_index][sys_i]);
453 : }
454 :
455 332944 : flux_container[fi->id()] =
456 166472 : flux_relaxation * flux + (1 - flux_relaxation) * flux_container[fi->id()];
457 :
458 : // We do the integral here
459 166472 : integrated_flux += flux * fi->faceArea() * fi->faceCoord();
460 : }
461 : }
462 14552 : }
463 :
464 : void
465 7276 : CHTHandler::sumIntegratedFluxes()
466 : {
467 14552 : for (const auto i : index_range(_integrated_boundary_heat_flux))
468 : {
469 : auto & integrated_fluxes = _integrated_boundary_heat_flux[i];
470 7276 : _problem.comm().sum(integrated_fluxes[NS::CHTSide::SOLID]);
471 7276 : _problem.comm().sum(integrated_fluxes[NS::CHTSide::FLUID]);
472 : }
473 7276 : }
474 :
475 : void
476 7276 : CHTHandler::printIntegratedFluxes() const
477 : {
478 14552 : for (const auto i : index_range(_integrated_boundary_heat_flux))
479 : {
480 : auto & integrated_fluxes = _integrated_boundary_heat_flux[i];
481 7276 : _console << " Iteration " << _fpi_it << " Boundary " << _cht_boundary_names[i]
482 7276 : << " flux on solid side " << integrated_fluxes[NS::CHTSide::SOLID]
483 7276 : << " flux on fluid side: " << integrated_fluxes[NS::CHTSide::FLUID] << std::endl;
484 : }
485 7276 : }
486 :
487 : void
488 4572 : CHTHandler::resetIntegratedFluxes()
489 : {
490 9144 : for (const auto i : index_range(_integrated_boundary_heat_flux))
491 9144 : _integrated_boundary_heat_flux[i] = std::vector<Real>({0.0, 0.0});
492 4572 : }
493 :
494 : bool
495 56908 : CHTHandler::converged() const
496 : {
497 56908 : if (_fpi_it >= _max_cht_fpi)
498 : return true;
499 :
500 40582 : for (const auto & boundary_flux : _integrated_boundary_heat_flux)
501 : {
502 10378 : const Real f1 = boundary_flux[0];
503 10378 : const Real f2 = boundary_flux[1];
504 :
505 : // Special case: both are zero at startup not converged yet
506 10378 : if (_fpi_it != 0 && (f1 == 0.0 && f2 == 0.0))
507 : return true;
508 :
509 : // These fluxes should be of opposite sign
510 10378 : const Real diff = std::abs(f1 + f2);
511 10378 : const Real denom = std::max({std::fabs(f1), std::fabs(f2), Real(1e-14)});
512 10378 : const Real rel_diff = diff / denom;
513 :
514 10378 : if (rel_diff >= _cht_heat_flux_tolerance)
515 : return false;
516 : }
517 :
518 30204 : return _fpi_it;
519 : }
520 :
521 : } // End FV namespace
522 : } // End Moose namespace
|