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 1132 : PatternedCartesianMeshGenerator::validParams()
23 : {
24 1132 : InputParameters params = PolygonMeshGeneratorBase::validParams();
25 2264 : params.addRequiredParam<std::vector<MeshGeneratorName>>(
26 : "inputs", "The names of the meshes forming the pattern.");
27 2264 : 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 2264 : MooseEnum cartesian_pattern_boundary("none expanded", "expanded");
33 2264 : params.addParam<MooseEnum>(
34 : "pattern_boundary", cartesian_pattern_boundary, "The boundary shape of the patterned mesh.");
35 2264 : params.addParam<bool>(
36 : "generate_core_metadata",
37 2264 : 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 3396 : params.addRangeCheckedParam<unsigned int>("background_intervals",
41 2264 : 3,
42 : "background_intervals>0",
43 : "Radial intervals in the assembly peripheral region.");
44 2264 : 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 2264 : params.addRangeCheckedParam<std::vector<Real>>(
50 : "duct_sizes", "duct_sizes>0.0", "Distance(s) from center to duct(s) inner boundaries.");
51 2264 : MooseEnum duct_sizes_style("apothem radius", "apothem");
52 2264 : 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 2264 : params.addRangeCheckedParam<std::vector<unsigned int>>(
57 : "duct_intervals", "duct_intervals>0", "Number of meshing intervals in each enclosing duct.");
58 2264 : params.addParam<bool>("uniform_mesh_on_sides",
59 2264 : false,
60 : "Whether the side elements are reorganized to have a uniform size.");
61 2264 : params.addParam<bool>("generate_control_drum_positions_file",
62 2264 : false,
63 : "Whether a positions file is generated in the core mesh mode.");
64 2264 : params.addParam<bool>("assign_control_drum_id",
65 2264 : false,
66 : "Whether control drum id is assigned to the mesh as an extra integer.");
67 1132 : std::string position_file_default = "positions_meta.data";
68 2264 : 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 2264 : params.addParam<Real>(
73 : "rotate_angle",
74 2264 : 0.0,
75 : "Rotate the entire patterned mesh by a certain degrees that is defined here.");
76 2264 : 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 2264 : 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 2264 : 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 1132 : params.addParam<std::vector<SubdomainName>>(
89 : "duct_block_names",
90 1132 : 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 2264 : params.addRangeCheckedParam<boundary_id_type>("external_boundary_id",
94 : "external_boundary_id>0",
95 : "Optional customized external boundary id.");
96 2264 : params.addParam<bool>("create_inward_interface_boundaries",
97 2264 : false,
98 : "Whether the inward interface boundary sidesets are created.");
99 2264 : params.addParam<bool>("create_outward_interface_boundaries",
100 2264 : true,
101 : "Whether the outward interface boundary sidesets are created.");
102 3396 : params.addParam<BoundaryName>(
103 : "stitching_boundary_name",
104 2264 : std::to_string(OUTER_SIDESET_ID),
105 : "Name of the boundary used for stitching pins togethers and with the background region");
106 2264 : params.addParam<BoundaryName>(
107 1132 : "external_boundary_name", BoundaryName(), "Optional customized external boundary name.");
108 2264 : params.addParam<bool>("deform_non_circular_region",
109 2264 : true,
110 : "Whether the non-circular region (outside the rings) can be deformed.");
111 2264 : params.addParam<std::vector<std::string>>("id_name", "List of extra integer ID set names");
112 2264 : params.addParam<std::vector<MeshGeneratorName>>(
113 : "exclude_id", "Name of input meshes to be excluded in ID generation");
114 2264 : std::vector<MooseEnum> option = {MooseEnum("cell pattern manual", "cell")};
115 2264 : params.addParam<std::vector<MooseEnum>>(
116 : "assign_type", option, "List of integer ID assignment types");
117 2264 : 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 2264 : 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 2264 : MooseEnum quad_elem_type("QUAD4 QUAD8 QUAD9", "QUAD4");
126 2264 : 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 2264 : params.addParam<bool>(
131 : "allow_unused_inputs",
132 2264 : false,
133 : "Whether additional input assemblies can be part of inputs without being used in lattice");
134 2264 : params.addParam<bool>(
135 : "verbose_stitching",
136 2264 : false,
137 : "Whether to output the number of nodes stitched when stitching pin and background meshes");
138 :
139 2264 : 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 2264 : params.addParamNamesToGroup(
145 : "generate_control_drum_positions_file assign_control_drum_id position_file", "Control Drum");
146 2264 : params.addParamNamesToGroup(
147 : "background_intervals duct_intervals uniform_mesh_on_sides deform_non_circular_region",
148 : "Mesh Density");
149 2264 : params.addParamNamesToGroup("id_name exclude_id assign_type id_pattern", "Reporting ID");
150 1132 : 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 1132 : return params;
156 4528 : }
157 :
158 575 : PatternedCartesianMeshGenerator::PatternedCartesianMeshGenerator(const InputParameters & parameters)
159 : : PolygonMeshGeneratorBase(parameters),
160 575 : _mesh_ptrs(getMeshes("inputs")),
161 1150 : _input_names(getParam<std::vector<MeshGeneratorName>>("inputs")),
162 1150 : _pattern(getParam<std::vector<std::vector<unsigned int>>>("pattern")),
163 1150 : _pattern_boundary(getParam<MooseEnum>("pattern_boundary")),
164 1150 : _generate_core_metadata(getParam<bool>("generate_core_metadata")),
165 1150 : _background_intervals(getParam<unsigned int>("background_intervals")),
166 1150 : _has_assembly_duct(isParamValid("duct_sizes")),
167 1813 : _duct_sizes(isParamValid("duct_sizes") ? getParam<std::vector<Real>>("duct_sizes")
168 : : std::vector<Real>()),
169 1150 : _duct_sizes_style(getParam<MooseEnum>("duct_sizes_style").template getEnum<PolygonSizeStyle>()),
170 1326 : _duct_intervals(isParamValid("duct_intervals")
171 575 : ? getParam<std::vector<unsigned int>>("duct_intervals")
172 : : std::vector<unsigned int>()),
173 1150 : _uniform_mesh_on_sides(getParam<bool>("uniform_mesh_on_sides")),
174 1150 : _generate_control_drum_positions_file(getParam<bool>("generate_control_drum_positions_file")),
175 1150 : _assign_control_drum_id(getParam<bool>("assign_control_drum_id")),
176 1150 : _rotate_angle(getParam<Real>("rotate_angle")),
177 1318 : _duct_block_ids(isParamValid("duct_block_ids")
178 575 : ? getParam<std::vector<subdomain_id_type>>("duct_block_ids")
179 : : std::vector<subdomain_id_type>()),
180 1150 : _duct_block_names(getParam<std::vector<SubdomainName>>("duct_block_names")),
181 575 : _stitching_boundary_name(getParam<BoundaryName>("stitching_boundary_name")),
182 1335 : _external_boundary_id(isParamValid("external_boundary_id")
183 945 : ? getParam<boundary_id_type>("external_boundary_id")
184 : : 0),
185 575 : _external_boundary_name(getParam<BoundaryName>("external_boundary_name")),
186 1150 : _create_inward_interface_boundaries(getParam<bool>("create_inward_interface_boundaries")),
187 1150 : _create_outward_interface_boundaries(getParam<bool>("create_outward_interface_boundaries")),
188 1150 : _deform_non_circular_region(getParam<bool>("deform_non_circular_region")),
189 1150 : _use_reporting_id(isParamValid("id_name")),
190 1150 : _use_exclude_id(isParamValid("exclude_id")),
191 1150 : _use_interface_boundary_id_shift(isParamValid("interface_boundary_id_shift_pattern")),
192 575 : _boundary_quad_elem_type(
193 575 : getParam<MooseEnum>("boundary_region_element_type").template getEnum<QUAD_ELEM_TYPE>()),
194 1150 : _allow_unused_inputs(getParam<bool>("allow_unused_inputs")),
195 1725 : _verbose_stitching(getParam<bool>("verbose_stitching"))
196 : {
197 575 : declareMeshProperty("pattern_pitch_meta", 0.0);
198 575 : declareMeshProperty("input_pitch_meta", 0.0);
199 1150 : declareMeshProperty<bool>("is_control_drum_meta", false);
200 1150 : declareMeshProperty<std::vector<Point>>("control_drum_positions", std::vector<Point>());
201 1150 : declareMeshProperty<std::vector<Real>>("control_drum_angles", std::vector<Real>());
202 575 : declareMeshProperty<std::vector<std::vector<Real>>>("control_drums_azimuthal_meta",
203 575 : std::vector<std::vector<Real>>());
204 1725 : declareMeshProperty<std::string>("position_file_name", getParam<std::string>("position_file"));
205 575 : declareMeshProperty<bool>("square_peripheral_trimmability", !_generate_core_metadata);
206 575 : declareMeshProperty<bool>("square_center_trimmability", true);
207 575 : declareMeshProperty<bool>("peripheral_modifier_compatible", _pattern_boundary == "expanded");
208 :
209 575 : const unsigned int n_pattern_layers = _pattern.size();
210 575 : declareMeshProperty("pattern_size", n_pattern_layers);
211 575 : 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 2256 : for (const auto & pattern_elem : _pattern)
218 : {
219 1685 : if (pattern_elem.empty())
220 2 : paramError("pattern",
221 : "The element of the two-dimensional array parameter pattern must not be empty.");
222 1683 : pattern_elem_size.emplace(pattern_elem.size());
223 1683 : pattern_max_array.push_back(*std::max_element(pattern_elem.begin(), pattern_elem.end()));
224 1683 : pattern_1d.insert(pattern_1d.end(), pattern_elem.begin(), pattern_elem.end());
225 : }
226 571 : 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 569 : 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 1134 : if (std::set<unsigned int>(pattern_1d.begin(), pattern_1d.end()).size() < _input_names.size() &&
235 44 : !_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 1130 : if (isParamValid("background_block_id"))
241 : {
242 621 : _peripheral_block_ids.push_back(getParam<subdomain_id_type>("background_block_id"));
243 207 : _peripheral_block_ids.insert(
244 : _peripheral_block_ids.end(), _duct_block_ids.begin(), _duct_block_ids.end());
245 : }
246 358 : 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 1126 : if (isParamValid("background_block_name"))
251 : {
252 492 : _peripheral_block_names.push_back(getParam<SubdomainName>("background_block_name"));
253 164 : _peripheral_block_names.insert(
254 : _peripheral_block_names.end(), _duct_block_names.begin(), _duct_block_names.end());
255 : }
256 399 : 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 561 : if (_pattern_boundary == "expanded")
261 : {
262 317 : for (unsigned int i = 1; i < _duct_sizes.size(); i++)
263 60 : if (_duct_sizes[i] <= _duct_sizes[i - 1])
264 2 : paramError("duct_sizes", "This parameter must be strictly ascending.");
265 257 : 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 255 : if (!_peripheral_block_names.empty() &&
269 162 : _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 506 : 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 302 : 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 600 : if (isParamValid("square_size"))
283 2 : paramError("square_size",
284 : "This parameter must not be provided when pattern_boundary is none.");
285 : }
286 :
287 549 : if (_use_interface_boundary_id_shift)
288 : {
289 : // check "interface_boundary_id_shift_pattern" parameter
290 : _interface_boundary_id_shift_pattern =
291 21 : getParam<std::vector<std::vector<boundary_id_type>>>("interface_boundary_id_shift_pattern");
292 7 : 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 28 : for (const auto i : make_range(_pattern.size()))
305 21 : 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 1098 : declareMeshProperty<bool>("interface_boundaries", false);
321 1098 : declareMeshProperty<std::set<boundary_id_type>>("interface_boundary_ids", {});
322 :
323 549 : if (_use_reporting_id)
324 : {
325 : // get reporting id name input
326 765 : _reporting_id_names = getParam<std::vector<std::string>>("id_name");
327 255 : const unsigned int num_reporting_ids = _reporting_id_names.size();
328 : // get reporting id assign type input
329 765 : const auto input_assign_types = getParam<std::vector<MooseEnum>>("assign_type");
330 255 : 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 531 : for (const auto i : make_range(num_reporting_ids))
335 : {
336 276 : _assign_types.push_back(
337 276 : input_assign_types[i].getEnum<ReportingIDGeneratorUtils::AssignType>());
338 276 : if (_assign_types[i] == ReportingIDGeneratorUtils::AssignType::manual)
339 21 : manual_ids.push_back(_reporting_id_names[i]);
340 : }
341 : // processing "id_pattern" input parameter
342 297 : if (manual_ids.size() > 0 && !isParamValid("id_pattern"))
343 0 : paramError("id_pattern", "required when 'manual' is defined in \"assign_type\"");
344 510 : if (isParamValid("id_pattern"))
345 : {
346 : const auto input_id_patterns =
347 42 : getParam<std::vector<std::vector<std::vector<dof_id_type>>>>("id_pattern");
348 14 : 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 35 : for (unsigned int i = 0; i < manual_ids.size(); ++i)
353 21 : _id_patterns[manual_ids[i]] = input_id_patterns[i];
354 14 : }
355 : // processing exlude id
356 255 : _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 255 : if (_use_exclude_id)
360 : {
361 : std::vector<MeshGeneratorName> exclude_id_name =
362 118 : getParam<std::vector<MeshGeneratorName>>("exclude_id");
363 256 : for (unsigned int i = 0; i < _input_names.size(); ++i)
364 : {
365 : _exclude_ids[i] = false;
366 335 : for (auto input_name : exclude_id_name)
367 197 : if (_input_names[i] == input_name)
368 : {
369 : _exclude_ids[i] = true;
370 : break;
371 : }
372 : }
373 59 : }
374 : else
375 545 : for (unsigned int i = 0; i < _input_names.size(); ++i)
376 : _exclude_ids[i] = false;
377 255 : }
378 549 : }
379 :
380 : std::unique_ptr<MeshBase>
381 533 : PatternedCartesianMeshGenerator::generate()
382 : {
383 533 : std::vector<std::unique_ptr<ReplicatedMesh>> meshes(_input_names.size());
384 1533 : for (const auto i : index_range(_input_names))
385 : {
386 : mooseAssert(_mesh_ptrs[i] && (*_mesh_ptrs[i]).get(), "nullptr mesh");
387 2004 : meshes[i] = dynamic_pointer_cast<ReplicatedMesh>(std::move(*_mesh_ptrs[i]));
388 1002 : 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 2004 : if (hasMeshProperty<bool>("flat_side_up", _input_names[i]))
392 484 : 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 531 : if (_pattern_boundary == "none" && _generate_core_metadata)
411 : {
412 : // Extract & check pitch and drum metadata
413 482 : 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 708 : 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 708 : 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 354 : 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 352 : is_control_drum_array.push_back(
449 352 : getMeshProperty<bool>("is_control_drum_meta", _input_names[i]));
450 352 : control_drum_azimuthal_array.push_back(
451 352 : getMeshProperty<bool>("is_control_drum_meta", _input_names[i])
452 711 : ? getMeshProperty<std::vector<Real>>("azimuthal_angle_meta", _input_names[i])
453 : : std::vector<Real>());
454 : }
455 128 : 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 126 : _pattern_pitch = pattern_pitch_array.front();
472 252 : setMeshProperty("input_pitch_meta", _pattern_pitch);
473 : }
474 : }
475 : else
476 : {
477 401 : if (_pattern_boundary == "expanded")
478 486 : _pattern_pitch = getParam<Real>("square_size");
479 :
480 1043 : 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 1288 : 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 642 : pitch_array.push_back(getMeshProperty<Real>("pitch_meta", _input_names[i]));
491 :
492 : num_sectors_per_side_array_tmp =
493 1284 : getMeshProperty<std::vector<unsigned int>>("num_sectors_per_side_meta", _input_names[i]);
494 642 : if (*std::max_element(num_sectors_per_side_array_tmp.begin(),
495 : num_sectors_per_side_array_tmp.end()) !=
496 642 : *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 642 : num_sectors_per_side_array.push_back(*num_sectors_per_side_array_tmp.begin());
503 642 : background_intervals_array.push_back(
504 642 : getMeshProperty<unsigned int>("background_intervals_meta", _input_names[i]));
505 642 : node_id_background_array.push_back(
506 642 : getMeshProperty<dof_id_type>("node_id_background_meta", _input_names[i]));
507 1284 : max_radius_array.push_back(getMeshProperty<Real>("max_radius_meta", _input_names[i]));
508 : }
509 399 : max_radius_global = *max_element(max_radius_array.begin(), max_radius_array.end());
510 399 : 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 798 : setMeshProperty("input_pitch_meta", pitch_array.front());
518 399 : if (*std::max_element(num_sectors_per_side_array.begin(), num_sectors_per_side_array.end()) !=
519 399 : *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 399 : if (_pattern_boundary != "expanded")
526 156 : _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 525 : const Real extra_dist_tol = _pattern_boundary == "expanded" ? pitch_array.front() / 10.0 : 0.0;
535 525 : 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 525 : if (_pattern_boundary == "expanded")
538 : {
539 243 : if (_has_assembly_duct)
540 : {
541 210 : for (unsigned int i = 0; i < _duct_sizes.size(); i++)
542 : {
543 132 : 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 132 : extra_dist.push_back(0.5 * (_duct_sizes[i] * 2.0 - pitch_array.front() * _pattern.size()));
547 132 : peripheral_duct_intervals.push_back(_duct_intervals[i]);
548 : }
549 78 : 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 241 : extra_dist.push_back(0.5 * (_pattern_pitch - pitch_array.front() * _pattern.size()));
557 241 : 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 241 : if (extra_dist.front() <= extra_dist_tol)
565 : {
566 71 : extra_dist_shift = extra_dist_shift_0 - extra_dist.front();
567 258 : for (Real & d : extra_dist)
568 187 : d += extra_dist_shift;
569 71 : y_min = _deform_non_circular_region
570 71 : ? max_radius_global // Currently use this, ideally this should be the max of the
571 : // outer layer radii
572 9 : : (pitch_array.front() / 2.0);
573 71 : y_max_0 = pitch_array.front() / 2.0 + extra_dist.front();
574 71 : y_max_n = y_max_0 - extra_dist_shift;
575 71 : 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 1042 : 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 521 : input_interface_boundary_ids.resize(_input_names.size());
589 521 : if (_use_interface_boundary_id_shift)
590 : {
591 21 : for (const auto i : make_range(_input_names.size()))
592 : {
593 28 : 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 14 : 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 28 : 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 1503 : for (const auto i : make_range(_input_names.size()))
605 1964 : if (hasMeshProperty<std::set<boundary_id_type>>("interface_boundary_ids", _input_names[i]))
606 : input_interface_boundary_ids[i] =
607 1074 : getMeshProperty<std::set<boundary_id_type>>("interface_boundary_ids", _input_names[i]);
608 :
609 803 : const Real input_pitch((_pattern_boundary == "expanded" || !_generate_core_metadata)
610 803 : ? 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 521 : std::unique_ptr<ReplicatedMesh> out_mesh;
617 :
618 2036 : for (unsigned i = 0; i < _pattern.size(); i++)
619 : {
620 : const Real deltax = 0.0;
621 1517 : Real deltay = -(Real)(i)*input_pitch;
622 :
623 6368 : for (unsigned int j = 0; j < _pattern[i].size(); j++)
624 : {
625 4853 : const auto pattern = _pattern[i][j];
626 4853 : 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 4853 : if (_pattern_boundary == "none")
631 : {
632 1965 : if (_generate_core_metadata && is_control_drum_array[pattern] == 1)
633 : {
634 7 : control_drum_positions_x.push_back(deltax + j * input_pitch);
635 7 : control_drum_positions_y.push_back(deltay);
636 7 : control_drum_azimuthals.push_back(control_drum_azimuthal_array[pattern]);
637 : }
638 :
639 1965 : if (j == 0 && i == 0)
640 : {
641 564 : out_mesh = dynamic_pointer_cast<ReplicatedMesh>(pattern_mesh.clone());
642 282 : if (_assign_control_drum_id && _generate_core_metadata)
643 118 : 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 282 : 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 2550 : 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 2888 : if (j == 0 && i == 0)
664 : {
665 : rotation_angle = 90.;
666 : mesh_type = CORNER_MESH;
667 : orientation = 0.;
668 : }
669 2649 : else if (j == 0 && i == _pattern.size() - 1)
670 : {
671 : rotation_angle = 180.;
672 : mesh_type = CORNER_MESH;
673 : orientation = -90.;
674 : }
675 2410 : else if (j == _pattern[i].size() - 1 && i == 0)
676 : {
677 : rotation_angle = 0.;
678 : mesh_type = CORNER_MESH;
679 : orientation = 90.;
680 : }
681 2171 : 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 1932 : else if (i == 0)
688 : {
689 : rotation_angle = 0.;
690 : mesh_type = BOUNDARY_MESH;
691 : orientation = 0.;
692 : }
693 1604 : else if (i == _pattern.size() - 1)
694 : {
695 : rotation_angle = 180.;
696 : mesh_type = BOUNDARY_MESH;
697 : orientation = -180.;
698 : }
699 1276 : else if (j == 0)
700 : {
701 : rotation_angle = 90.;
702 : mesh_type = BOUNDARY_MESH;
703 : orientation = -90.;
704 : }
705 948 : 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 2888 : if (on_periphery)
715 : {
716 2268 : auto tmp_peripheral_mesh = dynamic_pointer_cast<ReplicatedMesh>(pattern_mesh.clone());
717 2268 : 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 2268 : if (extra_dist_shift != 0)
727 828 : 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 2268 : if (_use_interface_boundary_id_shift)
732 56 : reassignBoundaryIDs(*tmp_peripheral_mesh,
733 : _interface_boundary_id_shift_pattern[i][j],
734 : input_interface_boundary_ids[pattern]);
735 :
736 2268 : 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 2029 : 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 2029 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, *out_mesh);
748 : const auto stitching_boundary_id_periph =
749 2029 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, *tmp_peripheral_mesh);
750 :
751 2029 : MeshTools::Modification::translate(
752 2029 : *tmp_peripheral_mesh, deltax + j * input_pitch, deltay, 0);
753 2029 : 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 2029 : _verbose_stitching);
759 : }
760 :
761 : continue;
762 2268 : }
763 : }
764 :
765 2303 : if (_assign_control_drum_id && _pattern_boundary == "none" && _generate_core_metadata)
766 461 : 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 2303 : 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 2303 : 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 10711 : for (auto const & id_name_pair : main_subdomain_map)
783 8408 : main_subdomain_map_name_list.emplace(id_name_pair.second);
784 2303 : 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 2301 : if (_use_interface_boundary_id_shift)
788 7 : 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 2301 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, *out_mesh);
794 : const auto stitching_boundary_id_2 =
795 2301 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, pattern_mesh);
796 2301 : 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 2301 : _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 2301 : 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 2301 : if (_use_interface_boundary_id_shift)
808 7 : 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 519 : auto side_list = out_mesh->get_boundary_info().build_side_list();
817 : const auto stitching_boundary_id =
818 519 : MooseMeshUtils::getBoundaryID(_stitching_boundary_name, *out_mesh);
819 161855 : for (auto & sl : side_list)
820 : {
821 : // Remove it as an internal boundary
822 161336 : if (std::get<2>(sl) == stitching_boundary_id)
823 34972 : if (out_mesh->elem_ptr(std::get<0>(sl))->neighbor_ptr(std::get<1>(sl)) != nullptr)
824 17592 : out_mesh->get_boundary_info().remove_side(
825 17592 : out_mesh->elem_ptr(std::get<0>(sl)), std::get<1>(sl), std::get<2>(sl));
826 : }
827 :
828 519 : out_mesh->get_boundary_info().clear_boundary_node_ids();
829 :
830 1038 : out_mesh->get_boundary_info().build_node_list_from_side_list();
831 519 : 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 519 : const Point origin_pt = MooseMeshUtils::meshCentroidCalculator(*out_mesh);
837 519 : const Real origin_x = origin_pt(0);
838 519 : const Real origin_y = origin_pt(1);
839 :
840 519 : MeshTools::Modification::translate(*out_mesh, -origin_x, -origin_y, 0);
841 :
842 519 : if (_uniform_mesh_on_sides)
843 : {
844 69 : MeshTools::Modification::rotate(*out_mesh, -45.0, 0.0, 0.0);
845 : const Real azi_tol = 1E-8;
846 18336 : for (unsigned int i = 0; i < node_list.size(); ++i)
847 : {
848 18267 : if (std::get<1>(node_list[i]) == stitching_boundary_id)
849 : {
850 : node_azi_list.push_back(
851 4464 : std::make_pair(atan2(out_mesh->node_ref(std::get<0>(node_list[i]))(1),
852 2232 : out_mesh->node_ref(std::get<0>(node_list[i]))(0)) *
853 2232 : 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 2232 : if (node_azi_list.back().first + 180.0 <= azi_tol)
857 33 : node_azi_list.back().first = 180;
858 : }
859 : }
860 69 : std::sort(node_azi_list.begin(), node_azi_list.end());
861 69 : const unsigned int side_intervals = node_azi_list.size() / SQUARE_NUM_SIDES;
862 345 : for (unsigned int i = 0; i < SQUARE_NUM_SIDES; i++)
863 : {
864 2508 : for (unsigned int j = 1; j <= side_intervals; j++)
865 : {
866 2232 : Real azi_corr_tmp = atan2((Real)j * 2.0 / (Real)side_intervals - 1.0, 1.0);
867 2232 : Real x_tmp = _pattern_pitch / 2.0;
868 2232 : Real y_tmp = x_tmp * std::tan(azi_corr_tmp);
869 2232 : nodeCoordRotate(x_tmp, y_tmp, (Real)i * 90.0 - 135.0);
870 2232 : Point p_tmp = Point(x_tmp, y_tmp, 0.0);
871 2232 : 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 69 : if (_boundary_quad_elem_type != QUAD_ELEM_TYPE::QUAD4)
879 14 : adjustPeripheralQuadraticElements(*out_mesh, _boundary_quad_elem_type);
880 :
881 69 : MeshTools::Modification::rotate(*out_mesh, 45.0, 0.0, 0.0);
882 : }
883 :
884 519 : 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 519 : 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 133 : for (unsigned int i = 0; i < control_drum_positions_x.size(); ++i)
894 : {
895 7 : control_drum_positions_x[i] -= origin_x;
896 7 : control_drum_positions_y[i] -= origin_y;
897 7 : nodeCoordRotate(control_drum_positions_x[i], control_drum_positions_y[i], _rotate_angle);
898 7 : Real cd_angle = atan2(control_drum_positions_y[i], control_drum_positions_x[i]);
899 :
900 147 : for (unsigned int j = 0; j < control_drum_azimuthals[i].size(); j++)
901 : {
902 140 : control_drum_azimuthals[i][j] += _rotate_angle;
903 140 : control_drum_azimuthals[i][j] =
904 140 : atan2(std::sin(control_drum_azimuthals[i][j] / 180.0 * M_PI),
905 140 : std::cos(control_drum_azimuthals[i][j] / 180.0 * M_PI)) /
906 140 : M_PI * 180.0; // quick way to move from -M_PI to M_PI
907 : }
908 7 : std::sort(control_drum_azimuthals[i].begin(), control_drum_azimuthals[i].end());
909 :
910 7 : if (std::abs(cd_angle) < azi_tol)
911 : cd_angle = 0;
912 7 : else if (cd_angle < 0.0)
913 7 : cd_angle += 2 * M_PI;
914 7 : control_drum_tmp.push_back(
915 7 : std::make_tuple(cd_angle,
916 7 : Point(control_drum_positions_x[i], control_drum_positions_y[i], 0.0),
917 : control_drum_azimuthals[i],
918 7 : i + 1)); // control drum index to help sort control_drum_id
919 : }
920 126 : 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 133 : for (unsigned int i = 0; i < control_drum_tmp.size(); ++i)
925 : {
926 7 : control_drum_positions.push_back(std::get<1>(control_drum_tmp[i]));
927 7 : control_drum_angles.push_back(std::get<0>(control_drum_tmp[i]));
928 7 : control_drums_azimuthal_meta.push_back(std::get<2>(control_drum_tmp[i]));
929 7 : control_drum_id_sorted.push_back(std::get<3>(control_drum_tmp[i]));
930 : }
931 126 : setMeshProperty("control_drum_positions", control_drum_positions);
932 126 : setMeshProperty("control_drum_angles", control_drum_angles);
933 126 : setMeshProperty("control_drums_azimuthal_meta", control_drums_azimuthal_meta);
934 :
935 126 : if (_assign_control_drum_id)
936 : {
937 59 : unsigned int drum_integer_index = out_mesh->get_elem_integer_index("control_drum_id");
938 138586 : for (const auto & elem : out_mesh->element_ptr_range())
939 : {
940 69234 : dof_id_type unsorted_control_drum_id = elem->get_extra_integer(drum_integer_index);
941 69234 : if (unsorted_control_drum_id != 0)
942 : {
943 280 : auto sorted_iter = std::find(control_drum_id_sorted.begin(),
944 : control_drum_id_sorted.end(),
945 : unsorted_control_drum_id);
946 280 : elem->set_extra_integer(drum_integer_index,
947 280 : std::distance(control_drum_id_sorted.begin(), sorted_iter) + 1);
948 : }
949 59 : }
950 : }
951 :
952 126 : if (_generate_control_drum_positions_file)
953 : {
954 7 : std::string position_file_name = getMeshProperty<std::string>("position_file_name", name());
955 7 : std::ofstream pos_file(position_file_name);
956 14 : for (unsigned int i = 0; i < control_drum_positions.size(); ++i)
957 7 : pos_file << control_drum_positions[i](0) << " " << control_drum_positions[i](1) << " 0.0\n";
958 7 : pos_file.close();
959 7 : }
960 126 : }
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 519 : if (_use_reporting_id)
965 247 : addReportingIDs(*out_mesh, meshes);
966 :
967 : // Assign customized peripheral block ids and names
968 519 : if (!_peripheral_block_ids.empty())
969 129570 : for (const auto & elem : out_mesh->active_element_ptr_range())
970 154080 : for (subdomain_id_type i = PERIPHERAL_ID_SHIFT; i <= PERIPHERAL_ID_SHIFT + _duct_sizes.size();
971 : i++)
972 109024 : if (elem->subdomain_id() == i)
973 : {
974 19536 : elem->subdomain_id() = _peripheral_block_ids[i - PERIPHERAL_ID_SHIFT];
975 19536 : break;
976 193 : }
977 519 : if (!_peripheral_block_names.empty())
978 : {
979 346 : for (unsigned i = 0; i < _peripheral_block_names.size(); i++)
980 194 : out_mesh->subdomain_name(_peripheral_block_ids.empty() ? (PERIPHERAL_ID_SHIFT + i)
981 : : _peripheral_block_ids[i]) =
982 194 : _peripheral_block_names[i];
983 : }
984 : // Assign customized outer surface boundary id and name
985 519 : if (_external_boundary_id > 0)
986 177 : MooseMesh::changeBoundaryId(*out_mesh, stitching_boundary_id, _external_boundary_id, false);
987 519 : if (!_external_boundary_name.empty())
988 : {
989 177 : out_mesh->get_boundary_info().sideset_name(_external_boundary_id > 0
990 : ? _external_boundary_id
991 : : (boundary_id_type)stitching_boundary_id) =
992 177 : _external_boundary_name;
993 177 : 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 1497 : 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 978 : 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 978 : 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 519 : _interface_boundary_id_shift_pattern,
1015 : boundary_ids,
1016 : input_interface_boundary_ids,
1017 519 : _use_interface_boundary_id_shift,
1018 519 : _create_inward_interface_boundaries || _create_outward_interface_boundaries,
1019 1038 : extra_dist.size());
1020 519 : if (interface_boundary_ids.size() > 0)
1021 : {
1022 345 : setMeshProperty("interface_boundaries", true);
1023 690 : setMeshProperty("interface_boundary_ids", interface_boundary_ids);
1024 : }
1025 :
1026 : out_mesh->unset_is_prepared();
1027 519 : auto mesh = dynamic_pointer_cast<MeshBase>(out_mesh);
1028 519 : return mesh;
1029 1038 : }
1030 :
1031 : void
1032 2268 : 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 2268 : const auto stitching_boundary_id = MooseMeshUtils::getBoundaryID(_stitching_boundary_name, mesh);
1050 :
1051 2268 : 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 3824 : peripheral_point_index = {{0, 1, 2}, {2, 3, 4}};
1056 : else
1057 : // side mesh has one side that needs a peripheral mesh.
1058 3936 : peripheral_point_index = {{0, 1, 5}};
1059 :
1060 : // extra_dist includes background and ducts.
1061 : // Loop to calculate the positions of the boundaries.
1062 6012 : 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 3744 : positionSetup(
1067 1476 : 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 8964 : 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 20880 : for (unsigned int vector_index = 0; vector_index < 3; vector_index++)
1075 : {
1076 15660 : sub_positions_inner.push_back(
1077 15660 : positions_inner[peripheral_point_index[peripheral_index][vector_index]]);
1078 15660 : sub_d_positions_outer.push_back(
1079 15660 : d_positions_outer[peripheral_point_index[peripheral_index][vector_index]]);
1080 : }
1081 5220 : 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 5220 : _create_inward_interface_boundaries,
1088 5220 : (i != extra_dist.size() - 1) &&
1089 5220 : _create_outward_interface_boundaries);
1090 :
1091 : // The other_mesh must be prepared before stitching
1092 5220 : meshp0->prepare_for_use();
1093 :
1094 : // rotate the peripheral mesh to the desired side of the hexagon.
1095 5220 : MeshTools::Modification::rotate(*meshp0, rotation_angle, 0, 0);
1096 5220 : mesh.stitch_meshes(
1097 5220 : *meshp0, stitching_boundary_id, OUTER_SIDESET_ID, TOLERANCE, true, _verbose_stitching);
1098 5220 : sub_positions_inner.resize(0);
1099 5220 : sub_d_positions_outer.resize(0);
1100 5220 : }
1101 : }
1102 6804 : }
1103 :
1104 : void
1105 3744 : 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 3744 : positions_inner.resize(0);
1113 3744 : 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 3744 : positions_inner.push_back(std::make_pair(-pitch / 2.0, pitch / 2.0 + extra_dist_in));
1139 3744 : positions_inner.push_back(std::make_pair(0.0, pitch / 2.0 + extra_dist_in));
1140 : positions_inner.push_back(
1141 3744 : std::make_pair(pitch / 2.0 + extra_dist_in, pitch / 2.0 + extra_dist_in));
1142 3744 : positions_inner.push_back(std::make_pair(pitch / 2.0 + extra_dist_in, 0.0));
1143 3744 : positions_inner.push_back(std::make_pair(pitch / 2.0 + extra_dist_in, -pitch / 2.0));
1144 3744 : 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 3744 : d_positions_outer.push_back(std::make_pair(0.0, extra_dist_out - extra_dist_in));
1149 3744 : d_positions_outer.push_back(std::make_pair(0.0, extra_dist_out - extra_dist_in));
1150 : d_positions_outer.push_back(
1151 3744 : std::make_pair(extra_dist_out - extra_dist_in, extra_dist_out - extra_dist_in));
1152 3744 : d_positions_outer.push_back(std::make_pair(extra_dist_out - extra_dist_in, 0.0));
1153 3744 : d_positions_outer.push_back(std::make_pair(extra_dist_out - extra_dist_in, 0.0));
1154 3744 : d_positions_outer.push_back(std::make_pair(0.0, extra_dist_out - extra_dist_in));
1155 3744 : }
1156 :
1157 : void
1158 247 : PatternedCartesianMeshGenerator::addReportingIDs(
1159 : MeshBase & mesh, const std::vector<std::unique_ptr<ReplicatedMesh>> & from_meshes) const
1160 : {
1161 247 : const unsigned int num_reporting_ids = _reporting_id_names.size();
1162 515 : for (unsigned int i = 0; i < num_reporting_ids; ++i)
1163 : {
1164 268 : const std::string element_id_name = _reporting_id_names[i];
1165 : unsigned int extra_id_index;
1166 268 : if (!mesh.has_elem_integer(element_id_name))
1167 536 : 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 804 : (isParamValid("background_block_id")) ? std::set<subdomain_id_type>({PERIPHERAL_ID_SHIFT})
1180 268 : : std::set<subdomain_id_type>();
1181 :
1182 : const bool using_manual_id =
1183 268 : (_assign_types[i] == ReportingIDGeneratorUtils::AssignType::manual);
1184 268 : ReportingIDGeneratorUtils::assignReportingIDs(mesh,
1185 : extra_id_index,
1186 : _assign_types[i],
1187 268 : _use_exclude_id,
1188 268 : _exclude_ids,
1189 268 : _pattern_boundary == "expanded",
1190 : background_block_ids,
1191 : from_meshes,
1192 : _pattern,
1193 : (using_manual_id)
1194 536 : ? _id_patterns.at(element_id_name)
1195 : : std::vector<std::vector<dof_id_type>>());
1196 : }
1197 247 : }
|
