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 "AdvancedConcentricCircleGenerator.h"
11 : #include "PolygonalMeshGenerationUtils.h"
12 :
13 : // C++ includes
14 : #include <cmath>
15 :
16 : registerMooseObject("ReactorApp", AdvancedConcentricCircleGenerator);
17 :
18 : InputParameters
19 704 : AdvancedConcentricCircleGenerator::validParams()
20 : {
21 704 : InputParameters params = ConcentricCircleGeneratorBase::validParams();
22 :
23 704 : params.makeParamRequired<std::vector<Real>>("ring_radii");
24 704 : params.makeParamRequired<std::vector<unsigned int>>("ring_intervals");
25 :
26 1408 : params.addRangeCheckedParam<unsigned int>(
27 : "num_sectors",
28 : "num_sectors>2",
29 : "Number of azimuthal sectors of the circular mesh to be generated.");
30 1408 : params.addRangeCheckedParam<std::vector<Real>>(
31 : "customized_azimuthal_angles",
32 : "customized_azimuthal_angles>=0&customized_azimuthal_angles<360",
33 : "List of the user-specified azimuthal angles of the nodes.");
34 :
35 1408 : params.addParamNamesToGroup("num_sectors customized_azimuthal_angles", "Azimuthal Mesh Density");
36 :
37 704 : params.addClassDescription("This AdvancedConcentricCircleGenerator object is designed to mesh a "
38 : "concentric circular geometry.");
39 :
40 704 : return params;
41 0 : }
42 :
43 365 : AdvancedConcentricCircleGenerator::AdvancedConcentricCircleGenerator(
44 365 : const InputParameters & parameters)
45 : : ConcentricCircleGeneratorBase(parameters),
46 517 : _azimuthal_angles(isParamValid("customized_azimuthal_angles")
47 353 : ? getParam<std::vector<Real>>("customized_azimuthal_angles")
48 : : std::vector<Real>()),
49 1248 : _num_sectors(isParamValid("num_sectors") ? getParam<unsigned int>("num_sectors")
50 365 : : _azimuthal_angles.size())
51 : {
52 353 : const unsigned short tri_order = _tri_elem_type == TRI_ELEM_TYPE::TRI3 ? 1 : 2;
53 353 : const unsigned short quad_order = _quad_elem_type == QUAD_ELEM_TYPE::QUAD4 ? 1 : 2;
54 : // 1. If the generated mesh has only one ring layer of triangular elements, then no
55 : // quad elements are generated;
56 : // 2. Otherwise, both types of elements are generated.
57 353 : _order = tri_order;
58 110 : if (_ring_radii.size() == 1 && _ring_intervals.front() == 1 &&
59 463 : _ring_inner_boundary_layer_params.intervals.front() == 0 &&
60 110 : _ring_outer_boundary_layer_params.intervals.front() == 0)
61 : {
62 110 : if (tri_order != quad_order)
63 28 : _quad_elem_type = tri_order == 1 ? QUAD_ELEM_TYPE::QUAD4 : QUAD_ELEM_TYPE::QUAD9;
64 : }
65 243 : else if (tri_order != quad_order)
66 2 : paramError("tri_element_type",
67 : "the element types of triangular and quadrilateral elements must be compatible if "
68 : "both types of elements are generated.");
69 :
70 351 : if (_num_sectors == 0)
71 2 : paramError(
72 : "num_sectors",
73 : "this parameter must be specified if 'customized_azimuthal_angles' is not provided.");
74 349 : if (_azimuthal_angles.empty())
75 : {
76 2622 : for (unsigned int i = 0; i < _num_sectors; i++)
77 : {
78 2355 : _azimuthal_angles.push_back((Real)i * 360.0 / (Real)_num_sectors);
79 2355 : _virtual_nums_sectors.push_back((Real)_num_sectors);
80 : }
81 : }
82 : else
83 : {
84 82 : if (_num_sectors != _azimuthal_angles.size())
85 2 : paramError("num_sectors",
86 : "this parameter must be equal to the size of 'customized_azimuthal_angles' if "
87 : "both parameters are provided.");
88 7397 : for (unsigned int i = 0; i < _azimuthal_angles.size(); i++)
89 : {
90 7321 : const Real azi_angle_interval = i == _azimuthal_angles.size() - 1
91 7321 : ? 360.0 + _azimuthal_angles[0] - _azimuthal_angles[i]
92 7245 : : _azimuthal_angles[i + 1] - _azimuthal_angles[i];
93 7321 : if (azi_angle_interval <= 0.0)
94 2 : paramError("customized_azimuthal_angles",
95 : "the azimuthal angles provided must be strictly increasing.");
96 7319 : else if (azi_angle_interval >= 120.0)
97 2 : paramError("customized_azimuthal_angles",
98 : "please make sure the circle azimuthal discretization angles are less than "
99 : "120.0 to avert awkward polygonization.");
100 7317 : _virtual_nums_sectors.push_back(360.0 / azi_angle_interval);
101 : }
102 : }
103 : // Customized interface boundary id/name related error messages
104 343 : if (_inward_interface_boundary_names.size() > 0 &&
105 2 : _inward_interface_boundary_names.size() != _ring_radii.size() - 1)
106 2 : paramError("inward_interface_boundary_names",
107 : "If provided, the length of this parameter must be identical to the total number of "
108 : "interfaces.");
109 341 : if (_outward_interface_boundary_names.size() > 0 &&
110 2 : _outward_interface_boundary_names.size() != _ring_radii.size() - 1)
111 2 : paramError("outward_interface_boundary_names",
112 : "If provided, the length of this parameter must be identical to the total number of "
113 : "interfaces.");
114 :
115 : const unsigned int num_innermost_ring_layers =
116 339 : _ring_inner_boundary_layer_params.intervals.front() + _ring_intervals.front() +
117 339 : _ring_outer_boundary_layer_params.intervals.front();
118 339 : if (!_ring_block_ids.empty() &&
119 : _ring_block_ids.size() !=
120 302 : (_ring_intervals.size() + (unsigned int)(num_innermost_ring_layers != 1)))
121 0 : paramError("ring_block_ids",
122 : "This parameter must have the appropriate size if it is provided. The size should "
123 : "be the same as the size of 'ring_intervals' if the innermost ring interval "
124 : "(including boundary layers) is unity; otherwise the size should be greater than "
125 : "the size of 'ring_intervals' by one. If 'quad_center_elements' is true, it is "
126 : "optional to only provide this parameter with the same size as 'ring_intervals'");
127 339 : if (!_ring_block_names.empty() &&
128 : _ring_block_names.size() !=
129 302 : (_ring_intervals.size() + (unsigned int)(num_innermost_ring_layers != 1)))
130 0 : paramError("ring_block_names",
131 : "This parameter must have the appropriate size if it is set. The size should be the "
132 : "same as the size of 'ring_intervals' if the innermost ring interval (including "
133 : "boundary layers) is unity; otherwise the size should be greater than the size of "
134 : "'ring_intervals' by one. If 'quad_center_elements' is true, it is optional to only "
135 : "provide this parameter with the same size as 'ring_intervals'");
136 907 : for (unsigned int i = 0; i < _ring_radii.size(); i++)
137 : {
138 568 : const Real layer_width = _ring_radii[i] - (i == 0 ? 0.0 : _ring_radii[i - 1]);
139 568 : _ring_inner_boundary_layer_params.fractions.push_back(
140 568 : _ring_inner_boundary_layer_params.widths[i] / layer_width);
141 568 : _ring_outer_boundary_layer_params.fractions.push_back(
142 568 : _ring_outer_boundary_layer_params.widths[i] / layer_width);
143 : }
144 907 : for (unsigned int i = 0; i < _ring_inner_boundary_layer_params.fractions.size(); i++)
145 568 : if (MooseUtils::absoluteFuzzyEqual(_ring_inner_boundary_layer_params.fractions[i], 0.0) &&
146 568 : _ring_inner_boundary_layer_params.intervals[i] > 0)
147 0 : paramError("ring_inner_boundary_layer_intervals",
148 : "Ring inner boundary layer must have zero interval if its thickness is zero.");
149 : else if (MooseUtils::absoluteFuzzyGreaterThan(_ring_inner_boundary_layer_params.fractions[i],
150 568 : 0.0) &&
151 0 : _ring_inner_boundary_layer_params.intervals[i] == 0)
152 0 : paramError(
153 : "ring_inner_boundary_layer_intervals",
154 : "Ring inner boundary layer must have non-zero interval if its thickness is not zero.");
155 907 : for (unsigned int i = 0; i < _ring_outer_boundary_layer_params.fractions.size(); i++)
156 : {
157 568 : if (MooseUtils::absoluteFuzzyEqual(_ring_outer_boundary_layer_params.fractions[i], 0.0) &&
158 568 : _ring_outer_boundary_layer_params.intervals[i] > 0)
159 0 : paramError("ring_outer_boundary_layer_intervals",
160 : "Ring outer boundary layer must have zero interval if its thickness is zero.");
161 : else if (MooseUtils::absoluteFuzzyGreaterThan(_ring_outer_boundary_layer_params.fractions[i],
162 568 : 0.0) &&
163 0 : _ring_outer_boundary_layer_params.intervals[i] == 0)
164 0 : paramError(
165 : "ring_outer_boundary_layer_intervals",
166 : "Ring outer boundary layer must have non-zero interval if its thickness is not zero.");
167 568 : if (_ring_inner_boundary_layer_params.fractions[i] +
168 : _ring_outer_boundary_layer_params.fractions[i] >=
169 : 1.0)
170 0 : paramError("ring_inner_boundary_layer_widths",
171 : "Summation of ring_inner_boundary_layer_widths and "
172 : "ring_outer_boundary_layer_widths cannot exceeds the ring layer width.");
173 : }
174 :
175 907 : for (unsigned int i = 0; i < _ring_radii.size(); i++)
176 : {
177 568 : const Real layer_width = _ring_radii[i] - (i == 0 ? 0.0 : _ring_radii[i - 1]);
178 568 : _ring_inner_boundary_layer_params.fractions.push_back(
179 568 : _ring_inner_boundary_layer_params.widths[i] / layer_width);
180 568 : _ring_outer_boundary_layer_params.fractions.push_back(
181 568 : _ring_outer_boundary_layer_params.widths[i] / layer_width);
182 : }
183 339 : }
184 :
185 : std::unique_ptr<MeshBase>
186 263 : AdvancedConcentricCircleGenerator::generate()
187 : {
188 : std::vector<Real> ring_radii_corr;
189 : std::vector<Real> mod_azimuthal_angles;
190 :
191 8942 : for (unsigned int i = 1; i < _azimuthal_angles.size(); i++)
192 : {
193 8679 : mod_azimuthal_angles.push_back(_azimuthal_angles[i - 1]);
194 8679 : if (_order == 2)
195 231 : mod_azimuthal_angles.push_back((_azimuthal_angles[i - 1] + _azimuthal_angles[i]) / 2.0);
196 : }
197 263 : mod_azimuthal_angles.push_back(_azimuthal_angles.back());
198 263 : if (_order == 2)
199 28 : mod_azimuthal_angles.push_back((_azimuthal_angles.back() + _azimuthal_angles.front() + 360.0) /
200 : 2.0);
201 :
202 : const Real corr_factor =
203 263 : _preserve_volumes
204 263 : ? PolygonalMeshGenerationUtils::radiusCorrectionFactor(mod_azimuthal_angles, true, _order)
205 : : 1.0;
206 :
207 705 : for (const auto & ring_radius : _ring_radii)
208 442 : ring_radii_corr.push_back(ring_radius * corr_factor);
209 :
210 : const multiBdryLayerParams empty_params = {
211 : std::vector<Real>(), std::vector<Real>(), std::vector<unsigned int>(), std::vector<Real>()};
212 263 : const singleBdryLayerParams empty_param = {0.0, 0.0, 0, 1.0};
213 :
214 : // A dummy pitch number is needed for callind buildSlice()
215 : // Any value larger than twice of the largest ring radius will work
216 263 : const Real dummy_pitch = _ring_radii.back() * 3.0;
217 : const unsigned int num_sectors_per_side = 1;
218 : const unsigned int background_intervals = 0;
219 : const Real background_radial_bias = 1.0;
220 :
221 : dof_id_type node_id_background_meta;
222 :
223 : auto mesh = buildSlice(ring_radii_corr,
224 263 : _ring_intervals,
225 263 : _ring_radial_biases,
226 263 : _ring_inner_boundary_layer_params,
227 263 : _ring_outer_boundary_layer_params,
228 526 : std::vector<Real>(),
229 263 : std::vector<unsigned int>(),
230 263 : std::vector<Real>(),
231 : empty_params,
232 : empty_params,
233 : dummy_pitch,
234 : num_sectors_per_side,
235 : background_intervals,
236 : background_radial_bias,
237 : empty_param,
238 : empty_param,
239 : node_id_background_meta,
240 : _virtual_nums_sectors[0],
241 : /*side_index*/ 1,
242 263 : std::vector<Real>(),
243 263 : _block_id_shift,
244 : /* quad_center_elements */ false,
245 : /* center_quad_factor */ 0.0,
246 263 : _create_inward_interface_boundaries,
247 263 : _create_outward_interface_boundaries,
248 263 : _interface_boundary_id_shift,
249 : 1.0,
250 : true,
251 : _tri_elem_type,
252 263 : _quad_elem_type);
253 263 : MeshTools::Modification::rotate(*mesh, -_azimuthal_angles[0], 0, 0);
254 :
255 8942 : for (unsigned int i = 1; i < _num_sectors; i++)
256 : {
257 : auto mesh_tmp = buildSlice(ring_radii_corr,
258 : _ring_intervals,
259 : _ring_radial_biases,
260 : _ring_inner_boundary_layer_params,
261 : _ring_outer_boundary_layer_params,
262 17358 : std::vector<Real>(),
263 8679 : std::vector<unsigned int>(),
264 8679 : std::vector<Real>(),
265 : empty_params,
266 : empty_params,
267 : dummy_pitch,
268 : num_sectors_per_side,
269 : background_intervals,
270 : background_radial_bias,
271 : empty_param,
272 : empty_param,
273 : node_id_background_meta,
274 8679 : _virtual_nums_sectors[i],
275 : /*side_index*/ i + 1,
276 8679 : std::vector<Real>(),
277 8679 : _block_id_shift,
278 : /* quad_center_elements */ false,
279 : /* center_quad_factor */ 0.0,
280 8679 : _create_inward_interface_boundaries,
281 8679 : _create_outward_interface_boundaries,
282 8679 : _interface_boundary_id_shift,
283 : 1.0,
284 : true,
285 : _tri_elem_type,
286 8679 : _quad_elem_type);
287 :
288 8679 : ReplicatedMesh other_mesh(*mesh_tmp);
289 8679 : MeshTools::Modification::rotate(other_mesh, -_azimuthal_angles[i], 0, 0);
290 8679 : mesh->prepare_for_use();
291 8679 : other_mesh.prepare_for_use();
292 : // As we rotate the mesh in the negative direction (see "-" before _azimuthal_angles[i]),
293 : // the order of SLICE_END and SLICE_BEGIN should be reversed compared to the similar call in
294 : // PolygonConcentricCircleMeshGeneratorBase.C
295 8679 : mesh->stitch_meshes(other_mesh, SLICE_END, SLICE_BEGIN, TOLERANCE, true, false);
296 8679 : other_mesh.clear();
297 8679 : }
298 :
299 263 : if (!_generate_side_specific_boundaries)
300 9205 : for (unsigned int i = 0; i < _num_sectors; i++)
301 8942 : mesh->get_boundary_info().remove_id(i + 1 + OUTER_SIDESET_ID_ALT);
302 :
303 : // An extra step to stich the first and last slices together
304 263 : mesh->stitch_surfaces(SLICE_END, SLICE_BEGIN, TOLERANCE, true, false);
305 :
306 263 : mesh->prepare_for_use();
307 :
308 : // Set up customized Block Names and/or IDs
309 263 : unsigned int block_it = 0;
310 263 : unsigned int ring_block_num = 0;
311 : std::vector<subdomain_id_type> block_ids_old;
312 : std::vector<subdomain_id_type> block_ids_new;
313 : std::vector<SubdomainName> block_names;
314 :
315 263 : ringBlockIdsNamesPreparer(block_it, ring_block_num, block_ids_old, block_ids_new, block_names);
316 :
317 263 : assignBlockIdsNames(
318 : *mesh, block_ids_old, block_ids_new, block_names, "AdvancedConcentricCircleGenerator");
319 :
320 : // Customized boundary ids and names
321 263 : if (_external_boundary_id > 0)
322 124 : MooseMesh::changeBoundaryId(*mesh, OUTER_SIDESET_ID, _external_boundary_id, true);
323 : else
324 139 : MooseMesh::changeBoundaryId(*mesh, OUTER_SIDESET_ID, 0, true);
325 263 : if (!_external_boundary_name.empty())
326 : {
327 124 : mesh->get_boundary_info().sideset_name(_external_boundary_id > 0 ? _external_boundary_id : 0) =
328 124 : _external_boundary_name;
329 124 : mesh->get_boundary_info().nodeset_name(_external_boundary_id > 0 ? _external_boundary_id : 0) =
330 : _external_boundary_name;
331 : }
332 :
333 263 : assignInterfaceBoundaryNames(*mesh);
334 :
335 : mesh->unset_is_prepared();
336 526 : return dynamic_pointer_cast<MeshBase>(mesh);
337 263 : }
|
