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 "PatternedCartesianMeshGenerator.h"
11 : #include "ReportingIDGeneratorUtils.h"
12 : #include "MooseUtils.h"
13 : #include "MooseMeshUtils.h"
14 :
15 : // C++ includes
16 : #include <cmath> // provides round, not std::round (see http://www.cplusplus.com/reference/cmath/round/)
17 : #include <fstream> // used to generate the optional control drum position file
18 :
19 : registerMooseObject("ReactorApp", PatternedCartesianMeshGenerator);
20 :
21 : InputParameters
22 1400 : PatternedCartesianMeshGenerator::validParams()
23 : {
24 1400 : InputParameters params = PolygonMeshGeneratorBase::validParams();
25 2800 : params.addRequiredParam<std::vector<MeshGeneratorName>>(
26 : "inputs", "The names of the meshes forming the pattern.");
27 2800 : params.addRequiredRangeCheckedParam<std::vector<std::vector<unsigned int>>>(
28 : "pattern",
29 : "pattern>=0",
30 : "A two-dimensional cartesian (square-shaped) array starting with the upper-left corner."
31 : "It is composed of indexes into the inputs vector");
32 2800 : MooseEnum cartesian_pattern_boundary("none expanded", "expanded");
33 2800 : params.addParam<MooseEnum>(
34 : "pattern_boundary", cartesian_pattern_boundary, "The boundary shape of the patterned mesh.");
35 2800 : params.addParam<bool>(
36 : "generate_core_metadata",
37 2800 : false,
38 : "A Boolean parameter that controls whether the core related metadata "
39 : "is generated for other MOOSE objects such as 'MultiControlDrumFunction' or not.");
40 4200 : params.addRangeCheckedParam<unsigned int>("background_intervals",
41 2800 : 3,
42 : "background_intervals>0",
43 : "Radial intervals in the assembly peripheral region.");
44 2800 : params.addRangeCheckedParam<Real>(
45 : "square_size",
46 : "square_size>0.0",
47 : "Size (side length) of the outmost square boundary to be generated; this is "
48 : "required only when pattern type is 'expanded'.");
49 2800 : params.addRangeCheckedParam<std::vector<Real>>(
50 : "duct_sizes", "duct_sizes>0.0", "Distance(s) from center to duct(s) inner boundaries.");
51 2800 : MooseEnum duct_sizes_style("apothem radius", "apothem");
52 2800 : params.addParam<MooseEnum>("duct_sizes_style",
53 : duct_sizes_style,
54 : "Style in which square center to duct distance(s) is given (apothem "
55 : "= center-to-face, radius = center-to-vertex).");
56 2800 : params.addRangeCheckedParam<std::vector<unsigned int>>(
57 : "duct_intervals", "duct_intervals>0", "Number of meshing intervals in each enclosing duct.");
58 2800 : params.addParam<bool>("uniform_mesh_on_sides",
59 2800 : false,
60 : "Whether the side elements are reorganized to have a uniform size.");
61 2800 : params.addParam<bool>("generate_control_drum_positions_file",
62 2800 : false,
63 : "Whether a positions file is generated in the core mesh mode.");
64 2800 : params.addParam<bool>("assign_control_drum_id",
65 2800 : false,
66 : "Whether control drum id is assigned to the mesh as an extra integer.");
67 1400 : std::string position_file_default = "positions_meta.data";
68 2800 : params.addParam<std::string>(
69 : "position_file", position_file_default, "Data file name to store control drum positions.");
70 : // A square pattern_boundary mesh can be used in "inputs" of `PatternedCartesianMeshGenerator`
71 : // without rotation or with rotation of 90, 180, or 270 degrees.
72 2800 : params.addParam<Real>(
73 : "rotate_angle",
74 2800 : 0.0,
75 : "Rotate the entire patterned mesh by a certain degrees that is defined here.");
76 2800 : params.addParam<subdomain_id_type>(
77 : "background_block_id",
78 : "Optional customized block id for the background block in 'assembly' mode; must be provided "
79 : "along with 'duct_block_ids' if 'duct_sizes' is provided.");
80 2800 : params.addParam<SubdomainName>(
81 : "background_block_name",
82 : "Optional customized block name for the background block in 'assembly' mode; must be "
83 : "provided along with 'duct_block_names' if 'duct_sizes' is provided.");
84 2800 : params.addParam<std::vector<subdomain_id_type>>(
85 : "duct_block_ids",
86 : "Optional customized block ids for each duct geometry block in 'assembly' mode; must be "
87 : "provided along with 'background_block_id'.");
88 1400 : params.addParam<std::vector<SubdomainName>>(
89 : "duct_block_names",
90 1400 : std::vector<SubdomainName>(),
91 : "Optional customized block names for each duct geometry block in 'assembly' mode; must be "
92 : "provided along with 'background_block_name'.");
93 2800 : params.addRangeCheckedParam<boundary_id_type>("external_boundary_id",
94 : "external_boundary_id>0",
95 : "Optional customized external boundary id.");
96 2800 : params.addParam<bool>("create_inward_interface_boundaries",
97 2800 : false,
98 : "Whether the inward interface boundary sidesets are created.");
99 2800 : params.addParam<bool>("create_outward_interface_boundaries",
100 2800 : true,
101 : "Whether the outward interface boundary sidesets are created.");
102 4200 : params.addParam<BoundaryName>(
103 : "stitching_boundary_name",
104 2800 : std::to_string(OUTER_SIDESET_ID),
105 : "Name of the boundary used for stitching pins togethers and with the background region");
106 2800 : params.addParam<BoundaryName>(
107 1400 : "external_boundary_name", BoundaryName(), "Optional customized external boundary name.");
108 2800 : params.addParam<bool>("deform_non_circular_region",
109 2800 : true,
110 : "Whether the non-circular region (outside the rings) can be deformed.");
111 2800 : params.addParam<std::vector<std::string>>("id_name", "List of extra integer ID set names");
112 2800 : params.addParam<std::vector<MeshGeneratorName>>(
113 : "exclude_id", "Name of input meshes to be excluded in ID generation");
114 2800 : std::vector<MooseEnum> option = {MooseEnum("cell pattern manual", "cell")};
115 2800 : params.addParam<std::vector<MooseEnum>>(
116 : "assign_type", option, "List of integer ID assignment types");
117 2800 : params.addParam<std::vector<std::vector<std::vector<dof_id_type>>>>(
118 : "id_pattern",
119 : "User-defined element IDs. A double-indexed array starting with the upper-left corner. When "
120 : "providing multiple patterns, each pattern should be separated using '|'");
121 2800 : params.addParam<std::vector<std::vector<boundary_id_type>>>(
122 : "interface_boundary_id_shift_pattern",
123 : "User-defined shift values for each pattern cell. A double-indexed array starting with the "
124 : "upper-left corner.");
125 2800 : MooseEnum quad_elem_type("QUAD4 QUAD8 QUAD9", "QUAD4");
126 2800 : params.addParam<MooseEnum>(
127 : "boundary_region_element_type",
128 : quad_elem_type,
129 : "Type of the quadrilateral elements to be generated in the boundary region.");
130 2800 : params.addParam<bool>(
131 : "allow_unused_inputs",
132 2800 : false,
133 : "Whether additional input assemblies can be part of inputs without being used in lattice");
134 2800 : params.addParam<bool>(
135 : "verbose_stitching",
136 2800 : false,
137 : "Whether to output the number of nodes stitched when stitching pin and background meshes");
138 :
139 2800 : params.addParamNamesToGroup(
140 : "pattern_boundary background_block_id background_block_name duct_block_ids duct_block_names "
141 : "external_boundary_id external_boundary_name create_inward_interface_boundaries "
142 : "create_outward_interface_boundaries boundary_region_element_type",
143 : "Customized Subdomain/Boundary");
144 2800 : params.addParamNamesToGroup(
145 : "generate_control_drum_positions_file assign_control_drum_id position_file", "Control Drum");
146 2800 : params.addParamNamesToGroup(
147 : "background_intervals duct_intervals uniform_mesh_on_sides deform_non_circular_region",
148 : "Mesh Density");
149 2800 : params.addParamNamesToGroup("id_name exclude_id assign_type id_pattern", "Reporting ID");
150 1400 : params.addClassDescription(
151 : "This PatternedCartesianMeshGenerator source code assembles square meshes into a square "
152 : "grid "
153 : "and optionally forces the outer boundary to be square and/or adds a duct.");
154 :
155 1400 : return params;
156 5600 : }
157 :
158 709 : PatternedCartesianMeshGenerator::PatternedCartesianMeshGenerator(const InputParameters & parameters)
159 : : PolygonMeshGeneratorBase(parameters),
160 709 : _mesh_ptrs(getMeshes("inputs")),
161 1418 : _input_names(getParam<std::vector<MeshGeneratorName>>("inputs")),
162 1418 : _pattern(getParam<std::vector<std::vector<unsigned int>>>("pattern")),
163 1418 : _pattern_boundary(getParam<MooseEnum>("pattern_boundary")),
164 1418 : _generate_core_metadata(getParam<bool>("generate_core_metadata")),
165 1418 : _background_intervals(getParam<unsigned int>("background_intervals")),
166 1418 : _has_assembly_duct(isParamValid("duct_sizes")),
167 2235 : _duct_sizes(isParamValid("duct_sizes") ? getParam<std::vector<Real>>("duct_sizes")
168 : : std::vector<Real>()),
169 1418 : _duct_sizes_style(getParam<MooseEnum>("duct_sizes_style").template getEnum<PolygonSizeStyle>()),
170 1634 : _duct_intervals(isParamValid("duct_intervals")
171 709 : ? getParam<std::vector<unsigned int>>("duct_intervals")
172 : : std::vector<unsigned int>()),
173 1418 : _uniform_mesh_on_sides(getParam<bool>("uniform_mesh_on_sides")),
174 1418 : _generate_control_drum_positions_file(getParam<bool>("generate_control_drum_positions_file")),
175 1418 : _assign_control_drum_id(getParam<bool>("assign_control_drum_id")),
176 1418 : _rotate_angle(getParam<Real>("rotate_angle")),
177 1626 : _duct_block_ids(isParamValid("duct_block_ids")
178 709 : ? getParam<std::vector<subdomain_id_type>>("duct_block_ids")
179 : : std::vector<subdomain_id_type>()),
180 1418 : _duct_block_names(getParam<std::vector<SubdomainName>>("duct_block_names")),
181 709 : _stitching_boundary_name(getParam<BoundaryName>("stitching_boundary_name")),
182 1645 : _external_boundary_id(isParamValid("external_boundary_id")
183 1163 : ? getParam<boundary_id_type>("external_boundary_id")
184 : : 0),
185 709 : _external_boundary_name(getParam<BoundaryName>("external_boundary_name")),
186 1418 : _create_inward_interface_boundaries(getParam<bool>("create_inward_interface_boundaries")),
187 1418 : _create_outward_interface_boundaries(getParam<bool>("create_outward_interface_boundaries")),
188 1418 : _deform_non_circular_region(getParam<bool>("deform_non_circular_region")),
189 1418 : _use_reporting_id(isParamValid("id_name")),
190 1418 : _use_exclude_id(isParamValid("exclude_id")),
191 1418 : _use_interface_boundary_id_shift(isParamValid("interface_boundary_id_shift_pattern")),
192 709 : _boundary_quad_elem_type(
193 709 : getParam<MooseEnum>("boundary_region_element_type").template getEnum<QUAD_ELEM_TYPE>()),
194 1418 : _allow_unused_inputs(getParam<bool>("allow_unused_inputs")),
195 2127 : _verbose_stitching(getParam<bool>("verbose_stitching"))
196 : {
197 709 : declareMeshProperty("pattern_pitch_meta", 0.0);
198 709 : declareMeshProperty("input_pitch_meta", 0.0);
199 1418 : declareMeshProperty<bool>("is_control_drum_meta", false);
200 1418 : declareMeshProperty<std::vector<Point>>("control_drum_positions", std::vector<Point>());
201 1418 : declareMeshProperty<std::vector<Real>>("control_drum_angles", std::vector<Real>());
202 709 : declareMeshProperty<std::vector<std::vector<Real>>>("control_drums_azimuthal_meta",
203 709 : std::vector<std::vector<Real>>());
204 2127 : declareMeshProperty<std::string>("position_file_name", getParam<std::string>("position_file"));
205 709 : declareMeshProperty<bool>("square_peripheral_trimmability", !_generate_core_metadata);
206 709 : declareMeshProperty<bool>("square_center_trimmability", true);
207 709 : declareMeshProperty<bool>("peripheral_modifier_compatible", _pattern_boundary == "expanded");
208 :
209 709 : const unsigned int n_pattern_layers = _pattern.size();
210 709 : declareMeshProperty("pattern_size", n_pattern_layers);
211 709 : if (n_pattern_layers == 1)
212 2 : paramError("pattern", "The length (layer number) of this parameter must be larger than unity.");
213 : // Examine pattern for consistency
214 : std::vector<unsigned int> pattern_max_array;
215 : std::vector<unsigned int> pattern_1d;
216 : std::set<unsigned int> pattern_elem_size;
217 2782 : for (const auto & pattern_elem : _pattern)
218 : {
219 2077 : if (pattern_elem.empty())
220 2 : paramError("pattern",
221 : "The element of the two-dimensional array parameter pattern must not be empty.");
222 2075 : pattern_elem_size.emplace(pattern_elem.size());
223 2075 : pattern_max_array.push_back(*std::max_element(pattern_elem.begin(), pattern_elem.end()));
224 2075 : pattern_1d.insert(pattern_1d.end(), pattern_elem.begin(), pattern_elem.end());
225 : }
226 705 : if (pattern_elem_size.size() > 1 || *pattern_elem_size.begin() != _pattern.size())
227 2 : paramError("pattern",
228 : "The two-dimensional array parameter pattern must have a correct square shape.");
229 :
230 703 : if (*std::max_element(pattern_max_array.begin(), pattern_max_array.end()) >= _input_names.size())
231 2 : paramError(
232 : "pattern",
233 : "Elements of this parameter must be smaller than the length of inputs (0-indexing).");
234 1402 : if (std::set<unsigned int>(pattern_1d.begin(), pattern_1d.end()).size() < _input_names.size() &&
235 56 : !_allow_unused_inputs)
236 2 : paramError("pattern",
237 : "All the meshes provided in inputs must be used in the lattice pattern. To bypass "
238 : "this requirement, set 'allow_unused_inputs = true'");
239 :
240 1398 : if (isParamValid("background_block_id"))
241 : {
242 759 : _peripheral_block_ids.push_back(getParam<subdomain_id_type>("background_block_id"));
243 253 : _peripheral_block_ids.insert(
244 : _peripheral_block_ids.end(), _duct_block_ids.begin(), _duct_block_ids.end());
245 : }
246 446 : else if (!_duct_block_ids.empty())
247 2 : paramError("background_block_id",
248 : "This parameter and duct_block_ids must be "
249 : "provided simultaneously.");
250 1394 : if (isParamValid("background_block_name"))
251 : {
252 600 : _peripheral_block_names.push_back(getParam<SubdomainName>("background_block_name"));
253 200 : _peripheral_block_names.insert(
254 : _peripheral_block_names.end(), _duct_block_names.begin(), _duct_block_names.end());
255 : }
256 497 : else if (!_duct_block_names.empty())
257 2 : paramError("background_block_name",
258 : "This parameter and duct_block_names must be provided simultaneously.");
259 :
260 695 : if (_pattern_boundary == "expanded")
261 : {
262 387 : for (unsigned int i = 1; i < _duct_sizes.size(); i++)
263 72 : if (_duct_sizes[i] <= _duct_sizes[i - 1])
264 2 : paramError("duct_sizes", "This parameter must be strictly ascending.");
265 315 : if (!_peripheral_block_ids.empty() && _peripheral_block_ids.size() != _duct_sizes.size() + 1)
266 2 : paramError("duct_block_ids",
267 : "This parameter, if provided, must have a length equal to length of duct_sizes.");
268 313 : if (!_peripheral_block_names.empty() &&
269 198 : _peripheral_block_names.size() != _duct_sizes.size() + 1)
270 2 : paramError("duct_block_names",
271 : "This parameter, if provided, must have a length equal to length of duct_sizes.");
272 622 : if (!isParamValid("square_size"))
273 2 : paramError("square_size",
274 : "This parameter must be provided when pattern_boundary is expanded.");
275 : }
276 : else
277 : {
278 378 : if (!_peripheral_block_ids.empty() || !_peripheral_block_names.empty())
279 2 : paramError("background_block_id",
280 : "This parameter and background_block_name must not be set when the "
281 : "pattern_boundary is none.");
282 752 : if (isParamValid("square_size"))
283 2 : paramError("square_size",
284 : "This parameter must not be provided when pattern_boundary is none.");
285 : }
286 :
287 683 : if (_use_interface_boundary_id_shift)
288 : {
289 : // check "interface_boundary_id_shift_pattern" parameter
290 : _interface_boundary_id_shift_pattern =
291 27 : getParam<std::vector<std::vector<boundary_id_type>>>("interface_boundary_id_shift_pattern");
292 9 : if (_interface_boundary_id_shift_pattern.size() != _pattern.size())
293 : {
294 : std::string shape_pattern =
295 0 : "(" + std::to_string(_pattern.size()) + ", " + std::to_string(_pattern[0].size()) + ") ";
296 0 : paramError("interface_boundary_id_shift_pattern",
297 0 : "This parameter, if provided, should have the same two-dimensional array shape " +
298 0 : shape_pattern +
299 : "as "
300 0 : "the 'pattern' parameter. First dimension '" +
301 0 : std::to_string(_interface_boundary_id_shift_pattern.size()) +
302 : "' does not match.");
303 : }
304 36 : for (const auto i : make_range(_pattern.size()))
305 27 : if (_interface_boundary_id_shift_pattern[i].size() != _pattern[i].size())
306 : {
307 0 : std::string shape_pattern = "(" + std::to_string(_pattern.size()) + ", " +
308 0 : std::to_string(_pattern[0].size()) + ") ";
309 0 : paramError(
310 : "interface_boundary_id_shift_pattern",
311 0 : "This parameter, if provided, should have the same two-dimensional array shape " +
312 0 : shape_pattern +
313 : "as "
314 0 : "the 'pattern' parameter. Second dimension '" +
315 0 : std::to_string(_interface_boundary_id_shift_pattern[i].size()) +
316 : "' does not match.");
317 : }
318 : }
319 : // declare metadata for internal interface boundaries
320 1366 : declareMeshProperty<bool>("interface_boundaries", false);
321 1366 : declareMeshProperty<std::set<boundary_id_type>>("interface_boundary_ids", {});
322 :
323 683 : if (_use_reporting_id)
324 : {
325 : // get reporting id name input
326 951 : _reporting_id_names = getParam<std::vector<std::string>>("id_name");
327 317 : const unsigned int num_reporting_ids = _reporting_id_names.size();
328 : // get reporting id assign type input
329 951 : const auto input_assign_types = getParam<std::vector<MooseEnum>>("assign_type");
330 317 : if (input_assign_types.size() != num_reporting_ids)
331 0 : paramError("assign_type", "This parameter must have a length equal to length of id_name.");
332 : // list of reporting id names using manual id patterns;
333 : std::vector<std::string> manual_ids;
334 661 : for (const auto i : make_range(num_reporting_ids))
335 : {
336 344 : _assign_types.push_back(
337 344 : input_assign_types[i].getEnum<ReportingIDGeneratorUtils::AssignType>());
338 344 : if (_assign_types[i] == ReportingIDGeneratorUtils::AssignType::manual)
339 27 : manual_ids.push_back(_reporting_id_names[i]);
340 : }
341 : // processing "id_pattern" input parameter
342 371 : if (manual_ids.size() > 0 && !isParamValid("id_pattern"))
343 0 : paramError("id_pattern", "required when 'manual' is defined in \"assign_type\"");
344 634 : if (isParamValid("id_pattern"))
345 : {
346 : const auto input_id_patterns =
347 54 : getParam<std::vector<std::vector<std::vector<dof_id_type>>>>("id_pattern");
348 18 : if (input_id_patterns.size() != manual_ids.size())
349 0 : paramError("id_pattern",
350 : "The number of patterns must be equal to the number of 'manual' types defined "
351 : "in \"assign_type\".");
352 45 : for (unsigned int i = 0; i < manual_ids.size(); ++i)
353 27 : _id_patterns[manual_ids[i]] = input_id_patterns[i];
354 18 : }
355 : // processing exlude id
356 317 : _exclude_ids.resize(_input_names.size());
357 : // in case of using 'exclude_id', create a vector containg flag for each input tile to indicate
358 : // whether it is excluded from reporting id assignment
359 317 : if (_use_exclude_id)
360 : {
361 : std::vector<MeshGeneratorName> exclude_id_name =
362 146 : getParam<std::vector<MeshGeneratorName>>("exclude_id");
363 312 : for (unsigned int i = 0; i < _input_names.size(); ++i)
364 : {
365 : _exclude_ids[i] = false;
366 405 : for (auto input_name : exclude_id_name)
367 239 : if (_input_names[i] == input_name)
368 : {
369 : _exclude_ids[i] = true;
370 : break;
371 : }
372 : }
373 73 : }
374 : else
375 683 : for (unsigned int i = 0; i < _input_names.size(); ++i)
376 : _exclude_ids[i] = false;
377 317 : }
378 683 : }
379 :
380 : std::unique_ptr<MeshBase>
381 667 : PatternedCartesianMeshGenerator::generate()
382 : {
383 667 : std::vector<std::unique_ptr<ReplicatedMesh>> meshes(_input_names.size());
384 1919 : for (const auto i : index_range(_input_names))
385 : {
386 : mooseAssert(_mesh_ptrs[i] && (*_mesh_ptrs[i]).get(), "nullptr mesh");
387 2508 : meshes[i] = dynamic_pointer_cast<ReplicatedMesh>(std::move(*_mesh_ptrs[i]));
388 1254 : if (!meshes[i])
389 0 : paramError("inputs", "Mesh '", _input_names[i], "' is not a replicated mesh but it must be");
390 : // throw an error message if the input mesh does not have a flat side up
391 2508 : if (hasMeshProperty<bool>("flat_side_up", _input_names[i]))
392 608 : if (!getMeshProperty<bool>("flat_side_up", _input_names[i]))
393 2 : paramError("inputs",
394 : "Mesh '",
395 2 : _input_names[i],
396 : "' does not have a flat side facing up, which is not supported.");
397 : }
398 :
399 : std::vector<Real> pitch_array;
400 : std::vector<unsigned int> num_sectors_per_side_array;
401 : std::vector<unsigned int> num_sectors_per_side_array_tmp;
402 : std::vector<std::vector<Real>> control_drum_azimuthal_array;
403 : std::vector<unsigned int> background_intervals_array;
404 : std::vector<dof_id_type> node_id_background_array;
405 : std::vector<Real> max_radius_array;
406 : std::vector<bool> is_control_drum_array;
407 : Real max_radius_global(0.0);
408 : std::vector<Real> pattern_pitch_array;
409 :
410 665 : if (_pattern_boundary == "none" && _generate_core_metadata)
411 : {
412 : // Extract & check pitch and drum metadata
413 598 : for (MooseIndex(_input_names) i = 0; i < _input_names.size(); ++i)
414 : {
415 : // throw an error message if the input mesh does not contain the required meta data
416 876 : if (!hasMeshProperty<Real>("pattern_pitch_meta", _input_names[i]))
417 0 : mooseError(
418 : "In PatternedCartesianMeshGenerator ",
419 0 : _name,
420 : ": the unit square input mesh does not contain appropriate meta data "
421 : "required for generating a core mesh. Involved input mesh: ",
422 0 : _input_names[i],
423 : "; metadata issue: 'pattern_pitch_meta' is missing. Note that "
424 : "'generate_core_metadata' is set to true, which"
425 : "means that the mesh generator is producing a core mesh by stitching the input "
426 : "assembly meshes together. Therefore,"
427 : "the input meshes must contain the metadata of assembly meshes, which can "
428 : "usually be either automatically assigned "
429 : "by using another PatternedCartesianMeshGenerator with 'generate_core_metadata' set as "
430 : "false or manually assigned by AddMetaDataGenerator.");
431 876 : pattern_pitch_array.push_back(getMeshProperty<Real>("pattern_pitch_meta", _input_names[i]));
432 : // throw an error message if the input mesh contains non-sense meta data
433 438 : if (pattern_pitch_array.back() == 0.0)
434 2 : mooseError(
435 : "In PatternedCartesianMeshGenerator ",
436 2 : _name,
437 : ": the unit square input mesh does not contain appropriate meta data "
438 : "required for generating a core mesh. Involved input mesh: ",
439 2 : _input_names[i],
440 : "; metadata issue: 'pattern_pitch_meta' is zero. Note that "
441 : "'generate_core_metadata' is set to true, which"
442 : "means that the mesh generator is producing a core mesh by stitching the input "
443 : "assembly meshes together. Therefore,"
444 : "the input meshes must contain the metadata of assembly meshes, which can "
445 : "usually be either automatically assigned "
446 : "by using another PatternedCartesianMeshGenerator with 'generate_core_metadata' set as "
447 : "false or manually assigned by AddMetaDataGenerator.");
448 436 : is_control_drum_array.push_back(
449 436 : getMeshProperty<bool>("is_control_drum_meta", _input_names[i]));
450 436 : control_drum_azimuthal_array.push_back(
451 436 : getMeshProperty<bool>("is_control_drum_meta", _input_names[i])
452 881 : ? getMeshProperty<std::vector<Real>>("azimuthal_angle_meta", _input_names[i])
453 : : std::vector<Real>());
454 : }
455 160 : if (!MooseUtils::absoluteFuzzyEqual(
456 : *std::max_element(pattern_pitch_array.begin(), pattern_pitch_array.end()),
457 : *std::min_element(pattern_pitch_array.begin(), pattern_pitch_array.end())))
458 4 : mooseError(
459 : "In PatternedCartesianMeshGenerator ",
460 2 : _name,
461 : ": pattern_pitch metadata values of all input mesh generators must be identical when "
462 : "pattern_boundary is 'none' and generate_core_metadata is true. Please check the "
463 : "parameters of the mesh generators that produce the input meshes. "
464 : "Note that some of these mesh generator, such as "
465 : "CartesianConcentricCircleAdaptiveBoundaryMeshGenerator and FlexiblePatternGenerator,"
466 : "may have different definitions of square size in their input parameters. Please refer "
467 : "to the documentation of these mesh generators.",
468 2 : pitchMetaDataErrorGenerator(_input_names, pattern_pitch_array, "pattern_pitch_meta"));
469 : else
470 : {
471 158 : _pattern_pitch = pattern_pitch_array.front();
472 316 : setMeshProperty("input_pitch_meta", _pattern_pitch);
473 : }
474 : }
475 : else
476 : {
477 503 : if (_pattern_boundary == "expanded")
478 602 : _pattern_pitch = getParam<Real>("square_size");
479 :
480 1313 : for (MooseIndex(_input_names) i = 0; i < _input_names.size(); ++i)
481 : {
482 : // throw an error message if the input mesh does not contain the required meta data
483 1624 : if (!hasMeshProperty<Real>("pitch_meta", _input_names[i]))
484 2 : mooseError("In PatternedCartesianMeshGenerator ",
485 2 : _name,
486 : ": the unit square input mesh does not contain appropriate meta data "
487 : "required for generating an assembly. Involved input mesh: ",
488 2 : _input_names[i],
489 : "; metadata issue: 'pitch_meta' is missing");
490 810 : pitch_array.push_back(getMeshProperty<Real>("pitch_meta", _input_names[i]));
491 :
492 : num_sectors_per_side_array_tmp =
493 1620 : getMeshProperty<std::vector<unsigned int>>("num_sectors_per_side_meta", _input_names[i]);
494 810 : if (*std::max_element(num_sectors_per_side_array_tmp.begin(),
495 : num_sectors_per_side_array_tmp.end()) !=
496 810 : *std::min_element(num_sectors_per_side_array_tmp.begin(),
497 : num_sectors_per_side_array_tmp.end()))
498 0 : mooseError("In PatternedCartesianMeshGenerator ",
499 0 : _name,
500 : ": num_sectors_per_side metadata values of all four sides of each input mesh "
501 : "generator must be identical.");
502 810 : num_sectors_per_side_array.push_back(*num_sectors_per_side_array_tmp.begin());
503 810 : background_intervals_array.push_back(
504 810 : getMeshProperty<unsigned int>("background_intervals_meta", _input_names[i]));
505 810 : node_id_background_array.push_back(
506 810 : getMeshProperty<dof_id_type>("node_id_background_meta", _input_names[i]));
507 1620 : max_radius_array.push_back(getMeshProperty<Real>("max_radius_meta", _input_names[i]));
508 : }
509 501 : max_radius_global = *max_element(max_radius_array.begin(), max_radius_array.end());
510 501 : if (!MooseUtils::absoluteFuzzyEqual(*std::max_element(pitch_array.begin(), pitch_array.end()),
511 : *std::min_element(pitch_array.begin(), pitch_array.end())))
512 0 : mooseError("In PatternedCartesianMeshGenerator ",
513 0 : _name,
514 : ": pitch metadata values of all input mesh generators must be identical. Please "
515 : "check the parameters of the mesh generators that produce the input meshes.",
516 0 : pitchMetaDataErrorGenerator(_input_names, pitch_array, "pitch_meta"));
517 1002 : setMeshProperty("input_pitch_meta", pitch_array.front());
518 501 : if (*std::max_element(num_sectors_per_side_array.begin(), num_sectors_per_side_array.end()) !=
519 501 : *std::min_element(num_sectors_per_side_array.begin(), num_sectors_per_side_array.end()))
520 0 : mooseError(
521 : "In PatternedCartesianMeshGenerator ",
522 0 : _name,
523 : ": num_sectors_per_side metadata values of all input mesh generators must be identical.");
524 :
525 501 : if (_pattern_boundary != "expanded")
526 200 : _pattern_pitch = pitch_array.front() * (Real)_pattern.size();
527 : }
528 :
529 : std::vector<Real> extra_dist;
530 : Real extra_dist_shift(0.0);
531 : Real y_min(0.0);
532 : Real y_max_0(0.0);
533 : Real y_max_n(0.0);
534 659 : const Real extra_dist_tol = _pattern_boundary == "expanded" ? pitch_array.front() / 10.0 : 0.0;
535 659 : const Real extra_dist_shift_0 = _pattern_boundary == "expanded" ? pitch_array.front() / 5.0 : 0.0;
536 : std::vector<unsigned int> peripheral_duct_intervals;
537 659 : if (_pattern_boundary == "expanded")
538 : {
539 301 : if (_has_assembly_duct)
540 : {
541 262 : for (unsigned int i = 0; i < _duct_sizes.size(); i++)
542 : {
543 164 : if (_duct_sizes_style == PolygonSizeStyle::radius)
544 0 : _duct_sizes[i] *= std::cos(M_PI / 4.0);
545 : /// As square geometry is used here, size of patterned mesh can be straightforwardly calculated.
546 164 : extra_dist.push_back(0.5 * (_duct_sizes[i] * 2.0 - pitch_array.front() * _pattern.size()));
547 164 : peripheral_duct_intervals.push_back(_duct_intervals[i]);
548 : }
549 98 : if (_duct_sizes.back() >= _pattern_pitch / 2.0)
550 2 : paramError("duct_sizes",
551 : "The duct sizes should not exceed the size of the square boundary.");
552 : }
553 : // calculate the distance between the larger square boundary and the boundary of stitched unit
554 : // squares this is used to decide whether deformation is needed when cut-off happens or when
555 : // the distance is small.
556 299 : extra_dist.push_back(0.5 * (_pattern_pitch - pitch_array.front() * _pattern.size()));
557 299 : peripheral_duct_intervals.insert(peripheral_duct_intervals.begin(), _background_intervals);
558 :
559 : // In some cases, when the external square size is small enough, the external square boundary
560 : // may either be very close to the input square meshes that are near the boundary or even cut
561 : // off by these squares. As long as the ring regions are not cut off, the input squares can
562 : // be deformed to accomodate the external square shape. This block sets up the range of mesh
563 : // region that needs to be deformed.
564 299 : if (extra_dist.front() <= extra_dist_tol)
565 : {
566 89 : extra_dist_shift = extra_dist_shift_0 - extra_dist.front();
567 322 : for (Real & d : extra_dist)
568 233 : d += extra_dist_shift;
569 89 : y_min = _deform_non_circular_region
570 89 : ? max_radius_global // Currently use this, ideally this should be the max of the
571 : // outer layer radii
572 11 : : (pitch_array.front() / 2.0);
573 89 : y_max_0 = pitch_array.front() / 2.0 + extra_dist.front();
574 89 : y_max_n = y_max_0 - extra_dist_shift;
575 89 : if (y_max_n <= y_min)
576 2 : mooseError("In PatternedCartesianMeshGenerator ",
577 2 : _name,
578 : ": the assembly is cut off so much that the internal structure that should not "
579 : "be altered is compromised.");
580 : }
581 : }
582 :
583 1310 : setMeshProperty("pattern_pitch_meta", _pattern_pitch);
584 :
585 : // create a list of interface boundary ids for each input mesh
586 : // NOTE: list of interface boundary ids is stored in mesh metadata
587 : std::vector<std::set<boundary_id_type>> input_interface_boundary_ids;
588 655 : input_interface_boundary_ids.resize(_input_names.size());
589 655 : if (_use_interface_boundary_id_shift)
590 : {
591 27 : for (const auto i : make_range(_input_names.size()))
592 : {
593 36 : if (!hasMeshProperty<bool>("interface_boundaries", _input_names[i]))
594 0 : mooseError("Metadata 'interface_boundaries' could not be found on the input mesh: ",
595 0 : _input_names[i]);
596 18 : if (!getMeshProperty<bool>("interface_boundaries", _input_names[i]))
597 0 : mooseError("Interface boundary ids were not constructed in the input mesh",
598 0 : _input_names[i]);
599 36 : if (!hasMeshProperty<std::set<boundary_id_type>>("interface_boundary_ids", _input_names[i]))
600 0 : mooseError("Metadata 'interface_boundary_ids' could not be found on the input mesh: ",
601 0 : _input_names[i]);
602 : }
603 : }
604 1889 : for (const auto i : make_range(_input_names.size()))
605 2468 : if (hasMeshProperty<std::set<boundary_id_type>>("interface_boundary_ids", _input_names[i]))
606 : input_interface_boundary_ids[i] =
607 1358 : getMeshProperty<std::set<boundary_id_type>>("interface_boundary_ids", _input_names[i]);
608 :
609 1013 : const Real input_pitch((_pattern_boundary == "expanded" || !_generate_core_metadata)
610 1013 : ? pitch_array.front()
611 : : _pattern_pitch);
612 : std::vector<Real> control_drum_positions_x;
613 : std::vector<Real> control_drum_positions_y;
614 : std::vector<std::vector<Real>> control_drum_azimuthals;
615 :
616 655 : std::unique_ptr<ReplicatedMesh> out_mesh;
617 :
618 2562 : for (unsigned i = 0; i < _pattern.size(); i++)
619 : {
620 : const Real deltax = 0.0;
621 1909 : Real deltay = -(Real)(i)*input_pitch;
622 :
623 8024 : for (unsigned int j = 0; j < _pattern[i].size(); j++)
624 : {
625 6117 : const auto pattern = _pattern[i][j];
626 6117 : ReplicatedMesh & pattern_mesh = *meshes[pattern];
627 :
628 : // No pattern boundary, so no need to wrap with a peripheral mesh
629 : // Use the inputs as they stand and translate accordingly later
630 6117 : if (_pattern_boundary == "none")
631 : {
632 2487 : if (_generate_core_metadata && is_control_drum_array[pattern] == 1)
633 : {
634 9 : control_drum_positions_x.push_back(deltax + j * input_pitch);
635 9 : control_drum_positions_y.push_back(deltay);
636 9 : control_drum_azimuthals.push_back(control_drum_azimuthal_array[pattern]);
637 : }
638 :
639 2487 : if (j == 0 && i == 0)
640 : {
641 716 : out_mesh = dynamic_pointer_cast<ReplicatedMesh>(pattern_mesh.clone());
642 358 : if (_assign_control_drum_id && _generate_core_metadata)
643 146 : out_mesh->add_elem_integer(
644 : "control_drum_id",
645 : true,
646 : is_control_drum_array[pattern] ? control_drum_azimuthals.size() : 0);
647 : // shift interface boundary ids
648 358 : if (_use_interface_boundary_id_shift)
649 0 : reassignBoundaryIDs(*out_mesh,
650 : _interface_boundary_id_shift_pattern[i][j],
651 : input_interface_boundary_ids[pattern]);
652 3202 : continue;
653 : }
654 : }
655 : // Has a pattern boundary so we need to wrap with a peripheral mesh
656 : else // has a pattern boundary
657 : {
658 : Real rotation_angle = std::numeric_limits<Real>::max();
659 : Real orientation = std::numeric_limits<Real>::max();
660 : unsigned int mesh_type = std::numeric_limits<unsigned int>::max();
661 : bool on_periphery = true;
662 :
663 3630 : if (j == 0 && i == 0)
664 : {
665 : rotation_angle = 90.;
666 : mesh_type = CORNER_MESH;
667 : orientation = 0.;
668 : }
669 3333 : else if (j == 0 && i == _pattern.size() - 1)
670 : {
671 : rotation_angle = 180.;
672 : mesh_type = CORNER_MESH;
673 : orientation = -90.;
674 : }
675 3036 : else if (j == _pattern[i].size() - 1 && i == 0)
676 : {
677 : rotation_angle = 0.;
678 : mesh_type = CORNER_MESH;
679 : orientation = 90.;
680 : }
681 2739 : else if (j == _pattern[i].size() - 1 && i == _pattern.size() - 1)
682 : {
683 : rotation_angle = -90.;
684 : mesh_type = CORNER_MESH;
685 : orientation = 180.;
686 : }
687 2442 : else if (i == 0)
688 : {
689 : rotation_angle = 0.;
690 : mesh_type = BOUNDARY_MESH;
691 : orientation = 0.;
692 : }
693 2028 : else if (i == _pattern.size() - 1)
694 : {
695 : rotation_angle = 180.;
696 : mesh_type = BOUNDARY_MESH;
697 : orientation = -180.;
698 : }
699 1614 : else if (j == 0)
700 : {
701 : rotation_angle = 90.;
702 : mesh_type = BOUNDARY_MESH;
703 : orientation = -90.;
704 : }
705 1200 : else if (j == _pattern[i].size() - 1)
706 : {
707 : rotation_angle = -90.;
708 : mesh_type = BOUNDARY_MESH;
709 : orientation = 90;
710 : }
711 : else
712 : on_periphery = false;
713 :
714 3630 : if (on_periphery)
715 : {
716 2844 : auto tmp_peripheral_mesh = dynamic_pointer_cast<ReplicatedMesh>(pattern_mesh.clone());
717 2844 : addPeripheralMesh(*tmp_peripheral_mesh,
718 : pattern,
719 : pitch_array.front(),
720 : extra_dist,
721 : num_sectors_per_side_array,
722 : peripheral_duct_intervals,
723 : rotation_angle,
724 : mesh_type);
725 :
726 2844 : if (extra_dist_shift != 0)
727 1044 : cutOffPolyDeform(
728 : *tmp_peripheral_mesh, orientation, y_max_0, y_max_n, y_min, mesh_type, 90.0);
729 :
730 : // Reassign interface boundary ids
731 2844 : if (_use_interface_boundary_id_shift)
732 72 : reassignBoundaryIDs(*tmp_peripheral_mesh,
733 : _interface_boundary_id_shift_pattern[i][j],
734 : input_interface_boundary_ids[pattern]);
735 :
736 2844 : if (i == 0 && j == 0)
737 : out_mesh = std::move(tmp_peripheral_mesh);
738 : else
739 : {
740 : // Retrieve subdomain name map from the mesh to be stitched and insert it to the main
741 : // subdomain map
742 : const auto & increment_subdomain_map = tmp_peripheral_mesh->get_subdomain_name_map();
743 2547 : out_mesh->set_subdomain_name_map().insert(increment_subdomain_map.begin(),
744 : increment_subdomain_map.end());
745 :
746 : const auto stitching_boundary_id_base =
747 2547 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, *out_mesh);
748 : const auto stitching_boundary_id_periph =
749 2547 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, *tmp_peripheral_mesh);
750 :
751 2547 : MeshTools::Modification::translate(
752 2547 : *tmp_peripheral_mesh, deltax + j * input_pitch, deltay, 0);
753 2547 : out_mesh->stitch_meshes(*tmp_peripheral_mesh,
754 : stitching_boundary_id_base,
755 : stitching_boundary_id_periph,
756 : TOLERANCE,
757 : /*clear_stitched_boundary_ids=*/true,
758 2547 : _verbose_stitching);
759 : }
760 :
761 : continue;
762 2844 : }
763 : }
764 :
765 2915 : if (_assign_control_drum_id && _pattern_boundary == "none" && _generate_core_metadata)
766 547 : pattern_mesh.add_elem_integer(
767 : "control_drum_id",
768 : true,
769 : is_control_drum_array[pattern] ? control_drum_azimuthals.size() : 0);
770 :
771 : // Translate to correct position
772 2915 : MeshTools::Modification::translate(pattern_mesh, deltax + j * input_pitch, deltay, 0);
773 :
774 : // Define a reference map variable for subdomain map
775 : auto & main_subdomain_map = out_mesh->set_subdomain_name_map();
776 : // Retrieve subdomain name map from the mesh to be stitched and insert it to the main
777 : // subdomain map
778 : const auto & increment_subdomain_map = pattern_mesh.get_subdomain_name_map();
779 2915 : main_subdomain_map.insert(increment_subdomain_map.begin(), increment_subdomain_map.end());
780 : // Check if one SubdomainName is shared by more than one subdomain ids
781 : std::set<SubdomainName> main_subdomain_map_name_list;
782 13405 : for (auto const & id_name_pair : main_subdomain_map)
783 10490 : main_subdomain_map_name_list.emplace(id_name_pair.second);
784 2915 : if (main_subdomain_map.size() != main_subdomain_map_name_list.size())
785 2 : paramError("inputs", "The input meshes contain subdomain name maps with conflicts.");
786 : // Reassign interface boundary ids
787 2913 : if (_use_interface_boundary_id_shift)
788 9 : reassignBoundaryIDs(pattern_mesh,
789 : _interface_boundary_id_shift_pattern[i][j],
790 : input_interface_boundary_ids[pattern]);
791 :
792 : const auto stitching_boundary_id_1 =
793 2913 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, *out_mesh);
794 : const auto stitching_boundary_id_2 =
795 2913 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, pattern_mesh);
796 2913 : out_mesh->stitch_meshes(pattern_mesh,
797 : stitching_boundary_id_1,
798 : stitching_boundary_id_2,
799 : TOLERANCE,
800 : /*clear_stitched_boundary_ids=*/false,
801 2913 : _verbose_stitching);
802 :
803 : // Translate back now that we've stitched so that anyone else that uses this mesh has it at
804 : // the origin
805 2913 : MeshTools::Modification::translate(pattern_mesh, -(deltax + j * input_pitch), -deltay, 0);
806 : // Roll back the changes in interface boundary ids for the same reason
807 2913 : if (_use_interface_boundary_id_shift)
808 9 : reassignBoundaryIDs(pattern_mesh,
809 : _interface_boundary_id_shift_pattern[i][j],
810 : input_interface_boundary_ids[pattern],
811 : true);
812 : }
813 : }
814 :
815 : // Check if stitching_boundary_id is really the external boundary. Correct if needed.
816 653 : auto side_list = out_mesh->get_boundary_info().build_side_list();
817 : const auto stitching_boundary_id =
818 653 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, *out_mesh);
819 205013 : for (auto & sl : side_list)
820 : {
821 : // Remove it as an internal boundary
822 204360 : if (std::get<2>(sl) == stitching_boundary_id)
823 43988 : if (out_mesh->elem_ptr(std::get<0>(sl))->neighbor_ptr(std::get<1>(sl)) != nullptr)
824 22248 : out_mesh->get_boundary_info().remove_side(
825 22248 : out_mesh->elem_ptr(std::get<0>(sl)), std::get<1>(sl), std::get<2>(sl));
826 : }
827 :
828 653 : out_mesh->get_boundary_info().clear_boundary_node_ids();
829 :
830 653 : out_mesh->get_boundary_info().build_node_list_from_side_list();
831 653 : const auto node_list = out_mesh->get_boundary_info().build_node_list();
832 :
833 : std::vector<Real> bd_x_list;
834 : std::vector<Real> bd_y_list;
835 : std::vector<std::pair<Real, dof_id_type>> node_azi_list;
836 653 : const Point origin_pt = MooseMeshUtils::meshCentroidCalculator(*out_mesh);
837 653 : const Real origin_x = origin_pt(0);
838 653 : const Real origin_y = origin_pt(1);
839 :
840 653 : MeshTools::Modification::translate(*out_mesh, -origin_x, -origin_y, 0);
841 :
842 653 : if (_uniform_mesh_on_sides)
843 : {
844 83 : MeshTools::Modification::rotate(*out_mesh, -45.0, 0.0, 0.0);
845 : const Real azi_tol = 1E-8;
846 23262 : for (unsigned int i = 0; i < node_list.size(); ++i)
847 : {
848 23179 : if (std::get<1>(node_list[i]) == stitching_boundary_id)
849 : {
850 : node_azi_list.push_back(
851 5488 : std::make_pair(atan2(out_mesh->node_ref(std::get<0>(node_list[i]))(1),
852 2744 : out_mesh->node_ref(std::get<0>(node_list[i]))(0)) *
853 2744 : 180.0 / M_PI,
854 : std::get<0>(node_list[i])));
855 : // correct the last node's angle value (180.0) if it becomes a negative value.
856 2744 : if (node_azi_list.back().first + 180.0 <= azi_tol)
857 45 : node_azi_list.back().first = 180;
858 : }
859 : }
860 83 : std::sort(node_azi_list.begin(), node_azi_list.end());
861 83 : const unsigned int side_intervals = node_azi_list.size() / SQUARE_NUM_SIDES;
862 415 : for (unsigned int i = 0; i < SQUARE_NUM_SIDES; i++)
863 : {
864 3076 : for (unsigned int j = 1; j <= side_intervals; j++)
865 : {
866 2744 : Real azi_corr_tmp = atan2((Real)j * 2.0 / (Real)side_intervals - 1.0, 1.0);
867 2744 : Real x_tmp = _pattern_pitch / 2.0;
868 2744 : Real y_tmp = x_tmp * std::tan(azi_corr_tmp);
869 2744 : nodeCoordRotate(x_tmp, y_tmp, (Real)i * 90.0 - 135.0);
870 2744 : Point p_tmp = Point(x_tmp, y_tmp, 0.0);
871 2744 : out_mesh->add_point(p_tmp, node_azi_list[i * side_intervals + j - 1].second);
872 : }
873 : }
874 :
875 : // if quadratic elements are used, additional nodes need to be adjusted based on the new
876 : // boundary node locations. adjust side mid-edge nodes to the midpoints of the corner
877 : // points, and if QUAD9, adjust center point to new centroid.
878 83 : if (_boundary_quad_elem_type != QUAD_ELEM_TYPE::QUAD4)
879 18 : adjustPeripheralQuadraticElements(*out_mesh, _boundary_quad_elem_type);
880 :
881 83 : MeshTools::Modification::rotate(*out_mesh, 45.0, 0.0, 0.0);
882 : }
883 :
884 653 : MeshTools::Modification::rotate(*out_mesh, _rotate_angle, 0.0, 0.0);
885 :
886 : // This combination of input parameters is usually used for core mesh generation by stitching
887 : // assembly meshes together.
888 653 : if (_pattern_boundary == "none" && _generate_core_metadata)
889 : {
890 : const Real azi_tol = 1E-8;
891 : std::vector<std::tuple<Real, Point, std::vector<Real>, dof_id_type>> control_drum_tmp;
892 : std::vector<dof_id_type> control_drum_id_sorted;
893 167 : for (unsigned int i = 0; i < control_drum_positions_x.size(); ++i)
894 : {
895 9 : control_drum_positions_x[i] -= origin_x;
896 9 : control_drum_positions_y[i] -= origin_y;
897 9 : nodeCoordRotate(control_drum_positions_x[i], control_drum_positions_y[i], _rotate_angle);
898 9 : Real cd_angle = atan2(control_drum_positions_y[i], control_drum_positions_x[i]);
899 :
900 189 : for (unsigned int j = 0; j < control_drum_azimuthals[i].size(); j++)
901 : {
902 180 : control_drum_azimuthals[i][j] += _rotate_angle;
903 180 : control_drum_azimuthals[i][j] =
904 180 : atan2(std::sin(control_drum_azimuthals[i][j] / 180.0 * M_PI),
905 180 : std::cos(control_drum_azimuthals[i][j] / 180.0 * M_PI)) /
906 180 : M_PI * 180.0; // quick way to move from -M_PI to M_PI
907 : }
908 9 : std::sort(control_drum_azimuthals[i].begin(), control_drum_azimuthals[i].end());
909 :
910 9 : if (std::abs(cd_angle) < azi_tol)
911 : cd_angle = 0;
912 9 : else if (cd_angle < 0.0)
913 9 : cd_angle += 2 * M_PI;
914 9 : control_drum_tmp.push_back(
915 9 : std::make_tuple(cd_angle,
916 9 : Point(control_drum_positions_x[i], control_drum_positions_y[i], 0.0),
917 : control_drum_azimuthals[i],
918 9 : i + 1)); // control drum index to help sort control_drum_id
919 : }
920 158 : std::sort(control_drum_tmp.begin(), control_drum_tmp.end());
921 : std::vector<Point> control_drum_positions;
922 : std::vector<Real> control_drum_angles;
923 : std::vector<std::vector<Real>> control_drums_azimuthal_meta;
924 167 : for (unsigned int i = 0; i < control_drum_tmp.size(); ++i)
925 : {
926 9 : control_drum_positions.push_back(std::get<1>(control_drum_tmp[i]));
927 9 : control_drum_angles.push_back(std::get<0>(control_drum_tmp[i]));
928 9 : control_drums_azimuthal_meta.push_back(std::get<2>(control_drum_tmp[i]));
929 9 : control_drum_id_sorted.push_back(std::get<3>(control_drum_tmp[i]));
930 : }
931 158 : setMeshProperty("control_drum_positions", control_drum_positions);
932 158 : setMeshProperty("control_drum_angles", control_drum_angles);
933 158 : setMeshProperty("control_drums_azimuthal_meta", control_drums_azimuthal_meta);
934 :
935 158 : if (_assign_control_drum_id)
936 : {
937 73 : unsigned int drum_integer_index = out_mesh->get_elem_integer_index("control_drum_id");
938 150246 : for (const auto & elem : out_mesh->element_ptr_range())
939 : {
940 75050 : dof_id_type unsorted_control_drum_id = elem->get_extra_integer(drum_integer_index);
941 75050 : if (unsorted_control_drum_id != 0)
942 : {
943 360 : auto sorted_iter = std::find(control_drum_id_sorted.begin(),
944 : control_drum_id_sorted.end(),
945 : unsorted_control_drum_id);
946 360 : elem->set_extra_integer(drum_integer_index,
947 360 : std::distance(control_drum_id_sorted.begin(), sorted_iter) + 1);
948 : }
949 73 : }
950 : }
951 :
952 158 : if (_generate_control_drum_positions_file)
953 : {
954 9 : std::string position_file_name = getMeshProperty<std::string>("position_file_name", name());
955 9 : std::ofstream pos_file(position_file_name);
956 18 : for (unsigned int i = 0; i < control_drum_positions.size(); ++i)
957 9 : pos_file << control_drum_positions[i](0) << " " << control_drum_positions[i](1) << " 0.0\n";
958 9 : pos_file.close();
959 9 : }
960 158 : }
961 :
962 : // add reporting IDs if _use_reporting_id is set true
963 : // NOTE: addReportingIDs should be called before applying customized peripheral block ids
964 653 : if (_use_reporting_id)
965 309 : addReportingIDs(*out_mesh, meshes);
966 :
967 : // Assign customized peripheral block ids and names
968 653 : if (!_peripheral_block_ids.empty())
969 162078 : for (const auto & elem : out_mesh->active_element_ptr_range())
970 192336 : for (subdomain_id_type i = PERIPHERAL_ID_SHIFT; i <= PERIPHERAL_ID_SHIFT + _duct_sizes.size();
971 : i++)
972 135952 : if (elem->subdomain_id() == i)
973 : {
974 24416 : elem->subdomain_id() = _peripheral_block_ids[i - PERIPHERAL_ID_SHIFT];
975 24416 : break;
976 239 : }
977 653 : if (!_peripheral_block_names.empty())
978 : {
979 430 : for (unsigned i = 0; i < _peripheral_block_names.size(); i++)
980 242 : out_mesh->subdomain_name(_peripheral_block_ids.empty() ? (PERIPHERAL_ID_SHIFT + i)
981 : : _peripheral_block_ids[i]) =
982 242 : _peripheral_block_names[i];
983 : }
984 : // Assign customized outer surface boundary id and name
985 653 : if (_external_boundary_id > 0)
986 219 : MooseMesh::changeBoundaryId(*out_mesh, stitching_boundary_id, _external_boundary_id, false);
987 653 : if (!_external_boundary_name.empty())
988 : {
989 219 : out_mesh->get_boundary_info().sideset_name(_external_boundary_id > 0
990 : ? _external_boundary_id
991 : : (boundary_id_type)stitching_boundary_id) =
992 219 : _external_boundary_name;
993 219 : out_mesh->get_boundary_info().nodeset_name(_external_boundary_id > 0
994 : ? _external_boundary_id
995 : : (boundary_id_type)stitching_boundary_id) =
996 : _external_boundary_name;
997 : }
998 : // Merge the boundary name maps of all the input meshed to generate the output mesh's boundary
999 : // name maps
1000 : auto & new_sideset_map = out_mesh->get_boundary_info().set_sideset_name_map();
1001 : auto & new_nodeset_map = out_mesh->get_boundary_info().set_nodeset_name_map();
1002 1883 : for (unsigned int i = 0; i < meshes.size(); i++)
1003 : {
1004 : const auto input_sideset_map = meshes[i]->get_boundary_info().get_sideset_name_map();
1005 1230 : new_sideset_map.insert(input_sideset_map.begin(), input_sideset_map.end());
1006 : const auto input_nodeset_map = meshes[i]->get_boundary_info().get_nodeset_name_map();
1007 1230 : new_nodeset_map.insert(input_nodeset_map.begin(), input_nodeset_map.end());
1008 : }
1009 :
1010 : // set mesh metadata related with interface boundary ids
1011 : const std::set<boundary_id_type> boundary_ids = out_mesh->get_boundary_info().get_boundary_ids();
1012 : const std::set<boundary_id_type> interface_boundary_ids = getInterfaceBoundaryIDs(
1013 : _pattern,
1014 653 : _interface_boundary_id_shift_pattern,
1015 : boundary_ids,
1016 : input_interface_boundary_ids,
1017 653 : _use_interface_boundary_id_shift,
1018 653 : _create_inward_interface_boundaries || _create_outward_interface_boundaries,
1019 1306 : extra_dist.size());
1020 653 : if (interface_boundary_ids.size() > 0)
1021 : {
1022 435 : setMeshProperty("interface_boundaries", true);
1023 870 : setMeshProperty("interface_boundary_ids", interface_boundary_ids);
1024 : }
1025 :
1026 : out_mesh->set_isnt_prepared();
1027 653 : auto mesh = dynamic_pointer_cast<MeshBase>(out_mesh);
1028 653 : return mesh;
1029 1306 : }
1030 :
1031 : void
1032 2844 : PatternedCartesianMeshGenerator::addPeripheralMesh(
1033 : ReplicatedMesh & mesh,
1034 : const unsigned int pattern,
1035 : const Real pitch,
1036 : const std::vector<Real> & extra_dist,
1037 : const std::vector<unsigned int> & num_sectors_per_side_array,
1038 : const std::vector<unsigned int> & peripheral_duct_intervals,
1039 : const Real rotation_angle,
1040 : const unsigned int mesh_type)
1041 : {
1042 : std::vector<std::pair<Real, Real>> positions_inner;
1043 : std::vector<std::pair<Real, Real>> d_positions_outer;
1044 :
1045 : std::vector<std::vector<unsigned int>> peripheral_point_index;
1046 : std::vector<std::pair<Real, Real>> sub_positions_inner;
1047 : std::vector<std::pair<Real, Real>> sub_d_positions_outer;
1048 :
1049 2844 : const auto stitching_boundary_id = MooseMeshUtils::getBoundaryID(_stitching_boundary_name, mesh);
1050 :
1051 2844 : if (mesh_type == CORNER_MESH)
1052 : // corner mesh has two sides that need peripheral meshes.
1053 : // each element has three sub-elements, representing beginning, middle, and ending azimuthal
1054 : // points
1055 4752 : peripheral_point_index = {{0, 1, 2}, {2, 3, 4}};
1056 : else
1057 : // side mesh has one side that needs a peripheral mesh.
1058 4968 : peripheral_point_index = {{0, 1, 5}};
1059 :
1060 : // extra_dist includes background and ducts.
1061 : // Loop to calculate the positions of the boundaries.
1062 7524 : for (unsigned int i = 0; i < extra_dist.size(); i++)
1063 : {
1064 : // Generate the node positions for the peripheral meshes
1065 : // The node positions are generated for the two possible cases: the edge and the corner
1066 4680 : positionSetup(
1067 1836 : positions_inner, d_positions_outer, i == 0 ? 0.0 : extra_dist[i - 1], extra_dist[i], pitch);
1068 :
1069 : // Loop for all applicable sides that need peripheral mesh (2 for corner and 1 for edge)
1070 11196 : for (unsigned int peripheral_index = 0; peripheral_index < peripheral_point_index.size();
1071 : peripheral_index++)
1072 : {
1073 : // Loop for beginning, middle and ending positions of a side
1074 26064 : for (unsigned int vector_index = 0; vector_index < 3; vector_index++)
1075 : {
1076 19548 : sub_positions_inner.push_back(
1077 19548 : positions_inner[peripheral_point_index[peripheral_index][vector_index]]);
1078 19548 : sub_d_positions_outer.push_back(
1079 19548 : d_positions_outer[peripheral_point_index[peripheral_index][vector_index]]);
1080 : }
1081 6516 : auto meshp0 = buildSimplePeripheral(num_sectors_per_side_array[pattern],
1082 : peripheral_duct_intervals[i],
1083 : sub_positions_inner,
1084 : sub_d_positions_outer,
1085 : i,
1086 : _boundary_quad_elem_type,
1087 6516 : _create_inward_interface_boundaries,
1088 6516 : (i != extra_dist.size() - 1) &&
1089 6516 : _create_outward_interface_boundaries);
1090 :
1091 : // The other_mesh must be prepared before stitching
1092 6516 : meshp0->prepare_for_use();
1093 :
1094 : // rotate the peripheral mesh to the desired side of the hexagon.
1095 6516 : MeshTools::Modification::rotate(*meshp0, rotation_angle, 0, 0);
1096 6516 : mesh.stitch_meshes(
1097 6516 : *meshp0, stitching_boundary_id, OUTER_SIDESET_ID, TOLERANCE, true, _verbose_stitching);
1098 6516 : sub_positions_inner.resize(0);
1099 6516 : sub_d_positions_outer.resize(0);
1100 6516 : }
1101 : }
1102 8532 : }
1103 :
1104 : void
1105 4680 : PatternedCartesianMeshGenerator::positionSetup(
1106 : std::vector<std::pair<Real, Real>> & positions_inner,
1107 : std::vector<std::pair<Real, Real>> & d_positions_outer,
1108 : const Real extra_dist_in,
1109 : const Real extra_dist_out,
1110 : const Real pitch) const
1111 : {
1112 4680 : positions_inner.resize(0);
1113 4680 : d_positions_outer.resize(0);
1114 : // Nine sets of positions are generated here, as shown below.
1115 : // CORNER MESH Peripheral {0 1 2} and {2 3 4}
1116 : //
1117 : // 0 1 2
1118 : // | | /
1119 : // | | /
1120 : // |____|____/
1121 : // | |
1122 : // | |____ 3
1123 : // | |
1124 : // |_________|____ 4
1125 : //
1126 : // EDGE MESH Peripheral {0 1 5}
1127 : //
1128 : // 0 1 5
1129 : // | | |
1130 : // | | |
1131 : // |____|____|
1132 : // | |
1133 : // | |
1134 : // | |
1135 : // |_________|
1136 :
1137 : // Inner positions defined from index 0 through 5 as shown in the above cartoon
1138 4680 : positions_inner.push_back(std::make_pair(-pitch / 2.0, pitch / 2.0 + extra_dist_in));
1139 4680 : positions_inner.push_back(std::make_pair(0.0, pitch / 2.0 + extra_dist_in));
1140 : positions_inner.push_back(
1141 4680 : std::make_pair(pitch / 2.0 + extra_dist_in, pitch / 2.0 + extra_dist_in));
1142 4680 : positions_inner.push_back(std::make_pair(pitch / 2.0 + extra_dist_in, 0.0));
1143 4680 : positions_inner.push_back(std::make_pair(pitch / 2.0 + extra_dist_in, -pitch / 2.0));
1144 4680 : positions_inner.push_back(std::make_pair(pitch / 2.0, pitch / 2.0 + extra_dist_in));
1145 :
1146 : // Outer positions (relative displacement from inner ones) defined from index 0 through 5 as shown
1147 : // in the above cartoon
1148 4680 : d_positions_outer.push_back(std::make_pair(0.0, extra_dist_out - extra_dist_in));
1149 4680 : d_positions_outer.push_back(std::make_pair(0.0, extra_dist_out - extra_dist_in));
1150 : d_positions_outer.push_back(
1151 4680 : std::make_pair(extra_dist_out - extra_dist_in, extra_dist_out - extra_dist_in));
1152 4680 : d_positions_outer.push_back(std::make_pair(extra_dist_out - extra_dist_in, 0.0));
1153 4680 : d_positions_outer.push_back(std::make_pair(extra_dist_out - extra_dist_in, 0.0));
1154 4680 : d_positions_outer.push_back(std::make_pair(0.0, extra_dist_out - extra_dist_in));
1155 4680 : }
1156 :
1157 : void
1158 309 : PatternedCartesianMeshGenerator::addReportingIDs(
1159 : MeshBase & mesh, const std::vector<std::unique_ptr<ReplicatedMesh>> & from_meshes) const
1160 : {
1161 309 : const unsigned int num_reporting_ids = _reporting_id_names.size();
1162 645 : for (unsigned int i = 0; i < num_reporting_ids; ++i)
1163 : {
1164 336 : const std::string element_id_name = _reporting_id_names[i];
1165 : unsigned int extra_id_index;
1166 336 : if (!mesh.has_elem_integer(element_id_name))
1167 672 : extra_id_index = mesh.add_elem_integer(element_id_name);
1168 : else
1169 : {
1170 0 : extra_id_index = mesh.get_elem_integer_index(element_id_name);
1171 0 : paramWarning(
1172 : "id_name", "An element integer with the name '", element_id_name, "' already exists");
1173 : }
1174 :
1175 : // assign reporting IDs to individual elements
1176 : // NOTE: background block id should be set "PERIPHERAL_ID_SHIFT" because this function is called
1177 : // before assigning the user-defined background block id
1178 : std::set<subdomain_id_type> background_block_ids =
1179 1008 : (isParamValid("background_block_id")) ? std::set<subdomain_id_type>({PERIPHERAL_ID_SHIFT})
1180 336 : : std::set<subdomain_id_type>();
1181 :
1182 : const bool using_manual_id =
1183 336 : (_assign_types[i] == ReportingIDGeneratorUtils::AssignType::manual);
1184 336 : ReportingIDGeneratorUtils::assignReportingIDs(mesh,
1185 : extra_id_index,
1186 : _assign_types[i],
1187 336 : _use_exclude_id,
1188 336 : _exclude_ids,
1189 336 : _pattern_boundary == "expanded",
1190 : background_block_ids,
1191 : from_meshes,
1192 : _pattern,
1193 : (using_manual_id)
1194 672 : ? _id_patterns.at(element_id_name)
1195 : : std::vector<std::vector<dof_id_type>>());
1196 : }
1197 309 : }
|
