www.mooseframework.org
ConcentricCircleMesh.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
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 "ConcentricCircleMesh.h"
11 #include "libmesh/face_quad4.h"
12 #include "MooseMesh.h"
13 #include "libmesh/mesh_modification.h"
14 #include "libmesh/serial_mesh.h"
15 #include "libmesh/boundary_info.h"
16 #include "libmesh/utility.h"
17 // C++ includes
18 #include <cmath> // provides round, not std::round (see http://www.cplusplus.com/reference/cmath/round/)
19 
21 
24 {
26  MooseEnum portion(
27  "full top_right top_left bottom_left bottom_right right_half left_half top_half bottom_half",
28  "full");
29  params.addRequiredParam<unsigned int>("num_sectors",
30  "num_sectors % 2 = 0, num_sectors > 0"
31  "Number of azimuthal sectors in each quadrant"
32  "'num_sectors' must be an even number.");
33  params.addRequiredParam<std::vector<Real>>("radii", "Radii of major concentric circles");
34  params.addRequiredParam<std::vector<unsigned int>>(
35  "rings", "Number of rings in each circle or in the moderator");
36  params.addRequiredParam<Real>("inner_mesh_fraction",
37  "Length of inner square / radius of the innermost circle");
38  params.addRequiredParam<bool>(
39  "has_outer_square",
40  "It determines if meshes for a outer square are added to concentric circle meshes.");
41  params.addRangeCheckedParam<Real>(
42  "pitch",
43  0.0,
44  "pitch>=0.0",
45  "The moderator can be added to complete meshes for one unit cell of fuel assembly."
46  "Elements are quad meshes.");
47  params.addParam<MooseEnum>("portion", portion, "Control of which part of mesh is created");
48  params.addRequiredParam<bool>(
49  "preserve_volumes", "Volume of concentric circles can be preserved using this function.");
50  params.addClassDescription("This ConcentricCircleMesh source code is to generate concentric "
51  "circle meshes.");
52  return params;
53 }
54 
56  : MooseMesh(parameters),
57  _num_sectors(getParam<unsigned int>("num_sectors")),
58  _radii(getParam<std::vector<Real>>("radii")),
59  _rings(getParam<std::vector<unsigned int>>("rings")),
60  _inner_mesh_fraction(getParam<Real>("inner_mesh_fraction")),
61  _has_outer_square(getParam<bool>("has_outer_square")),
62  _pitch(getParam<Real>("pitch")),
63  _preserve_volumes(getParam<bool>("preserve_volumes")),
64  _portion(getParam<MooseEnum>("portion"))
65 {
66 
67  if (_num_sectors % 2 != 0)
68  mooseError("ConcentricCircleMesh: num_sectors must be an even number.");
69 
70  // radii data check
71  for (unsigned i = 0; i < _radii.size() - 1; ++i)
72  if (_radii[i] > _radii[i + 1])
73  mooseError("Radii must be provided in order by starting with the smallest radius and "
74  "providing the following gradual radii.");
75 
76  // condition for setting the size of inner squares.
77  if (_inner_mesh_fraction > std::cos(M_PI / 4))
78  mooseError("The aspect ratio can not be larger than cos(PI/4).");
79 
80  // size of 'rings' check
82  {
83  if (_rings.size() != _radii.size() + 1)
84  mooseError("The size of 'rings' must be equal to the size of 'radii' plus 1.");
85  }
86  else
87  {
88  if (_rings.size() != _radii.size())
89  mooseError("The size of 'rings' must be equal to the size of 'radii'.");
90  }
91  // pitch / 2 must be bigger than any raddi.
93  for (unsigned i = 0; i < _radii.size(); ++i)
94  if (_pitch / 2 < _radii[i])
95  mooseError("The pitch / 2 must be larger than any radii.");
96 }
97 
98 void
100 {
101  // Get the actual libMesh mesh
102  ReplicatedMesh & mesh = cast_ref<ReplicatedMesh &>(getMesh());
103  // Set dimension of mesh
104  mesh.set_mesh_dimension(2);
105  mesh.set_spatial_dimension(2);
106  BoundaryInfo & boundary_info = mesh.get_boundary_info();
107 
108  // Creating real mesh concentric circles
109  // i: index for _rings, j: index for _radii
110  std::vector<Real> total_concentric_circles;
111  unsigned int j = 0;
112  while (j < _radii.size())
113  {
114  unsigned int i = 0;
115  if (j == 0)
116  while (i < _rings[j])
117  {
118  total_concentric_circles.push_back(_inner_mesh_fraction * _radii[j] +
119  (_radii[j] - _inner_mesh_fraction * _radii[j]) /
120  _rings[j] * (i + 1));
121  ++i;
122  }
123  else
124  while (i < _rings[j])
125  {
126  total_concentric_circles.push_back(_radii[j - 1] +
127  (_radii[j] - _radii[j - 1]) / _rings[j] * (i + 1));
128  ++i;
129  }
130  ++j;
131  }
132 
133  // volume preserving function is used to conserve volume.
134  const Real d_angle = M_PI / 2 / _num_sectors;
135 
136  if (_preserve_volumes)
137  {
138  Real original_radius = 0.0;
139  for (unsigned i = 0; i < total_concentric_circles.size(); ++i)
140  {
141  // volume preserving function for the center circle
142  if (i == 0)
143  {
144  const Real target_area = M_PI * Utility::pow<2>(total_concentric_circles[i]);
145  Real modified_radius = std::sqrt(2 * target_area / std::sin(d_angle) / _num_sectors / 4);
146  original_radius = total_concentric_circles[i];
147  total_concentric_circles[i] = modified_radius;
148  }
149  else
150  {
151  // volume preserving functions for outer circles
152  const Real target_area = M_PI * (Utility::pow<2>(total_concentric_circles[i]) -
153  Utility::pow<2>(original_radius));
154  Real modified_radius = std::sqrt(target_area / std::sin(d_angle) / _num_sectors / 2 +
155  Utility::pow<2>(total_concentric_circles[i - 1]));
156  original_radius = total_concentric_circles[i];
157  total_concentric_circles[i] = modified_radius;
158  }
159  }
160  }
161 
162  // number of total nodes
163  unsigned num_total_nodes = 0;
164  if (_has_outer_square)
165  num_total_nodes = Utility::pow<2>(_num_sectors / 2 + 1) +
166  (_num_sectors + 1) * (total_concentric_circles.size() + _rings.back()) +
167  (_num_sectors + 1);
168  else
169  num_total_nodes = Utility::pow<2>(_num_sectors / 2 + 1) +
170  (_num_sectors + 1) * total_concentric_circles.size();
171 
172  std::vector<Node *> nodes(num_total_nodes);
173  unsigned node_id = 0;
174 
175  // for adding nodes for the square at the center of the circle
176  for (unsigned i = 0; i <= _num_sectors / 2; ++i)
177  {
178  const Real x = i * _inner_mesh_fraction * total_concentric_circles[0] / (_num_sectors / 2);
179  for (unsigned j = 0; j <= _num_sectors / 2; ++j)
180  {
181  const Real y = j * _inner_mesh_fraction * total_concentric_circles[0] / (_num_sectors / 2);
182  nodes[node_id] = mesh.add_point(Point(x, y, 0.0), node_id);
183  ++node_id;
184  }
185  }
186 
187  // for adding the outer nodes of the square
188  Real current_radius = 0.0;
189 
190  for (unsigned layers = 0; layers < total_concentric_circles.size(); ++layers)
191  {
192  current_radius = total_concentric_circles[layers];
193  for (unsigned num_outer_nodes = 0; num_outer_nodes <= _num_sectors; ++num_outer_nodes)
194  {
195  const Real x = current_radius * std::cos(num_outer_nodes * d_angle);
196  const Real y = current_radius * std::sin(num_outer_nodes * d_angle);
197  nodes[node_id] = mesh.add_point(Point(x, y, 0.0), node_id);
198  ++node_id;
199  }
200  }
201 
202  // adding nodes for the unit cell of fuel assembly.
203  if (_has_outer_square)
204  {
205  Real current_radius_moderator = 0.0;
206  for (unsigned i = 1; i <= _rings.back(); ++i)
207  {
208  current_radius_moderator =
209  _radii.back() + i * (_pitch / 2 - _radii.back()) / (_rings.back() + 1);
210  total_concentric_circles.push_back(current_radius_moderator);
211  for (unsigned num_outer_nodes = 0; num_outer_nodes <= _num_sectors; ++num_outer_nodes)
212  {
213  const Real x = current_radius_moderator * std::cos(num_outer_nodes * d_angle);
214  const Real y = current_radius_moderator * std::sin(num_outer_nodes * d_angle);
215  nodes[node_id] = mesh.add_point(Point(x, y, 0.0), node_id);
216  ++node_id;
217  }
218  }
219 
220  for (unsigned j = 0; j < _num_sectors / 2 + 1; ++j)
221  {
222  const Real x = _pitch / 2;
223  const Real y = _pitch / 2 * std::tan(j * d_angle);
224  nodes[node_id] = mesh.add_point(Point(x, y, 0.0), node_id);
225  ++node_id;
226  }
227 
228  for (unsigned i = 0; i < _num_sectors / 2; ++i)
229  {
230  const Real x = _pitch / 2 * std::cos((i + _num_sectors / 2 + 1) * d_angle) /
231  std::sin((i + _num_sectors / 2 + 1) * d_angle);
232  const Real y = _pitch / 2;
233  nodes[node_id] = mesh.add_point(Point(x, y, 0.0), node_id);
234  ++node_id;
235  }
236  }
237  // Currently, index, limit, counter variables use the int type because of the 'modulo' function.
238  // adding elements
239  int index = 0;
240  int limit = 0;
241  int standard = static_cast<int>(_num_sectors);
242 
243  // This is to set the limit for the index
244  if (standard > 4)
245  {
246  int additional_term = 0;
247  int counter = standard;
248  while (counter > 4)
249  {
250  counter = counter - 2;
251  additional_term = additional_term + counter;
252  }
253  limit = standard + additional_term;
254  }
255  else if (standard == 4)
256  limit = standard;
257 
258  // SubdomainIDs set up
259  std::vector<unsigned int> subdomainIDs;
260  for (unsigned int i = 0; i < _rings.size(); ++i)
261  for (unsigned int j = 0; j < _rings[i]; ++j)
262  subdomainIDs.push_back(i + 1);
263 
264  if (_has_outer_square)
265  subdomainIDs.push_back(subdomainIDs.back());
266  // adding elements in the square
267  while (index <= limit)
268  {
269  Elem * elem = mesh.add_elem(new Quad4);
270  elem->set_node(0) = nodes[index];
271  elem->set_node(1) = nodes[index + _num_sectors / 2 + 1];
272  elem->set_node(2) = nodes[index + _num_sectors / 2 + 2];
273  elem->set_node(3) = nodes[index + 1];
274  elem->subdomain_id() = subdomainIDs[0];
275 
276  if (index < standard / 2)
277  boundary_info.add_side(elem, 3, 1);
278  if (index % (standard / 2 + 1) == 0)
279  boundary_info.add_side(elem, 0, 2);
280 
281  ++index;
282  if ((index - standard / 2) % (standard / 2 + 1) == 0)
283  ++index;
284  }
285 
286  index = (_num_sectors / 2 + 1) * (_num_sectors / 2);
287  limit = (_num_sectors / 2) * (_num_sectors / 2 + 2);
288 
289  // adding elements in one outer layer of the square (right side)
290  while (index < limit)
291  {
292  Elem * elem = mesh.add_elem(new Quad4);
293  elem->set_node(0) = nodes[index];
294  elem->set_node(1) = nodes[index + _num_sectors / 2 + 1];
295  elem->set_node(2) = nodes[index + _num_sectors / 2 + 2];
296  elem->set_node(3) = nodes[index + 1];
297  elem->subdomain_id() = subdomainIDs[0];
298 
299  if (index == (standard / 2 + 1) * (standard / 2))
300  boundary_info.add_side(elem, 0, 2);
301 
302  ++index;
303  }
304 
305  // adding elements in one outer layer of the square (left side)
306  int counter = 0;
307  while (index != standard / 2)
308  {
309  Elem * elem = mesh.add_elem(new Quad4);
310  elem->set_node(0) = nodes[index];
311  elem->set_node(1) = nodes[index + (_num_sectors / 2 + 1) + counter * (_num_sectors / 2 + 2)];
312  elem->set_node(2) =
313  nodes[index + (_num_sectors / 2 + 1) + counter * (_num_sectors / 2 + 2) + 1];
314  elem->set_node(3) = nodes[index - _num_sectors / 2 - 1];
315  elem->subdomain_id() = subdomainIDs[0];
316 
317  if (index == standard + 1)
318  boundary_info.add_side(elem, 2, 1);
319 
320  index = index - _num_sectors / 2 - 1;
321  ++counter;
322  }
323 
324  // adding elements for other concentric circles
325  index = Utility::pow<2>(_num_sectors / 2 + 1);
326  limit = static_cast<int>(num_total_nodes) - standard - 2;
327  int num_nodes_boundary = Utility::pow<2>(_num_sectors / 2 + 1) + _num_sectors + 1;
328 
329  counter = 0;
330  while (index < limit)
331  {
332 
333  Elem * elem = mesh.add_elem(new Quad4);
334  elem->set_node(0) = nodes[index];
335  elem->set_node(1) = nodes[index + _num_sectors + 1];
336  elem->set_node(2) = nodes[index + _num_sectors + 2];
337  elem->set_node(3) = nodes[index + 1];
338 
339  for (int i = 0; i < static_cast<int>(subdomainIDs.size()) - 1; ++i)
340  if (index < limit - (standard + 1) * i && index >= limit - (standard + 1) * (i + 1))
341  elem->subdomain_id() = subdomainIDs[subdomainIDs.size() - 1 - i];
342 
343  int const initial = Utility::pow<2>(standard / 2 + 1);
344  int const final = Utility::pow<2>(standard / 2 + 1) + _num_sectors - 1;
345 
346  if ((index - initial) % (standard + 1) == 0)
347  boundary_info.add_side(elem, 0, 2);
348  if ((index - final) % (standard + 1) == 0)
349  boundary_info.add_side(elem, 2, 1);
350  if (index > limit - (standard + 1))
351  {
352  if (_has_outer_square)
353  {
354  if (index < limit - standard + standard / 2)
355  boundary_info.add_side(elem, 1, 3);
356  else
357  boundary_info.add_side(elem, 1, 4);
358  }
359  else
360  {
361  boundary_info.add_side(elem, 1, 3);
362  }
363  }
364  ++index;
365  if (index == (num_nodes_boundary + counter * (standard + 1)) - 1)
366  {
367  ++index;
368  ++counter;
369  }
370  }
371 
372  // This is to set boundary names.
373  boundary_info.sideset_name(1) = "left";
374  boundary_info.sideset_name(2) = "bottom";
375 
376  if (!_has_outer_square)
377  boundary_info.sideset_name(3) = "outer";
378  else
379  {
380  boundary_info.sideset_name(3) = "right";
381  boundary_info.sideset_name(4) = "top";
382  }
383 
384  if (_portion == "top_left")
385  {
386  MeshTools::Modification::rotate(mesh, 90, 0, 0);
387  boundary_info.sideset_name(1) = "bottom";
388  boundary_info.sideset_name(2) = "right";
389 
390  if (!_has_outer_square)
391  boundary_info.sideset_name(3) = "outer";
392  else
393  {
394  boundary_info.sideset_name(3) = "top";
395  boundary_info.sideset_name(4) = "left";
396  }
397  }
398  else if (_portion == "bottom_left")
399  {
400  MeshTools::Modification::rotate(mesh, 180, 0, 0);
401  boundary_info.sideset_name(1) = "right";
402  boundary_info.sideset_name(2) = "top";
403 
404  if (!_has_outer_square)
405  boundary_info.sideset_name(3) = "outer";
406  else
407  {
408  boundary_info.sideset_name(3) = "left";
409  boundary_info.sideset_name(4) = "bottom";
410  }
411  }
412  else if (_portion == "bottom_right")
413  {
414  MeshTools::Modification::rotate(mesh, 270, 0, 0);
415  boundary_info.sideset_name(1) = "top";
416  boundary_info.sideset_name(2) = "left";
417 
418  if (!_has_outer_square)
419  boundary_info.sideset_name(3) = "outer";
420  else
421  {
422  boundary_info.sideset_name(3) = "bottom";
423  boundary_info.sideset_name(4) = "right";
424  }
425  }
426 
427  else if (_portion == "top_half")
428  {
429  ReplicatedMesh other_mesh(mesh);
430  // This is to rotate the mesh and also to reset boundary IDs.
431  MeshTools::Modification::rotate(other_mesh, 90, 0, 0);
432  if (_has_outer_square)
433  {
434  MeshTools::Modification::change_boundary_id(other_mesh, 1, 5);
435  MeshTools::Modification::change_boundary_id(other_mesh, 2, 6);
436  MeshTools::Modification::change_boundary_id(other_mesh, 3, 7);
437  MeshTools::Modification::change_boundary_id(other_mesh, 4, 1);
438  MeshTools::Modification::change_boundary_id(other_mesh, 5, 2);
439  MeshTools::Modification::change_boundary_id(other_mesh, 6, 3);
440  MeshTools::Modification::change_boundary_id(other_mesh, 7, 4);
441  mesh.prepare_for_use();
442  other_mesh.prepare_for_use();
443  mesh.stitch_meshes(other_mesh, 1, 3, TOLERANCE, true);
444  mesh.get_boundary_info().sideset_name(1) = "left";
445  mesh.get_boundary_info().sideset_name(2) = "bottom";
446  mesh.get_boundary_info().sideset_name(3) = "right";
447  mesh.get_boundary_info().sideset_name(4) = "top";
448  }
449  else
450  {
451  MeshTools::Modification::change_boundary_id(other_mesh, 1, 5);
452  MeshTools::Modification::change_boundary_id(other_mesh, 2, 1);
453  MeshTools::Modification::change_boundary_id(other_mesh, 5, 2);
454  mesh.prepare_for_use();
455  other_mesh.prepare_for_use();
456  mesh.stitch_meshes(other_mesh, 1, 1, TOLERANCE, true);
457 
458  MeshTools::Modification::change_boundary_id(mesh, 2, 1);
459  MeshTools::Modification::change_boundary_id(mesh, 3, 2);
460  mesh.get_boundary_info().sideset_name(1) = "bottom";
461  mesh.get_boundary_info().sideset_name(2) = "outer";
462  }
463  other_mesh.clear();
464  }
465 
466  else if (_portion == "right_half")
467  {
468  ReplicatedMesh other_mesh(mesh);
469  // This is to rotate the mesh and also to reset boundary IDs.
470  MeshTools::Modification::rotate(other_mesh, 270, 0, 0);
471  if (_has_outer_square)
472  {
473  MeshTools::Modification::change_boundary_id(other_mesh, 1, 5);
474  MeshTools::Modification::change_boundary_id(other_mesh, 2, 6);
475  MeshTools::Modification::change_boundary_id(other_mesh, 3, 7);
476  MeshTools::Modification::change_boundary_id(other_mesh, 4, 3);
477  MeshTools::Modification::change_boundary_id(other_mesh, 5, 4);
478  MeshTools::Modification::change_boundary_id(other_mesh, 6, 1);
479  MeshTools::Modification::change_boundary_id(other_mesh, 7, 2);
480  mesh.prepare_for_use();
481  other_mesh.prepare_for_use();
482  mesh.stitch_meshes(other_mesh, 2, 4, TOLERANCE, true);
483  mesh.get_boundary_info().sideset_name(1) = "left";
484  mesh.get_boundary_info().sideset_name(2) = "bottom";
485  mesh.get_boundary_info().sideset_name(3) = "right";
486  mesh.get_boundary_info().sideset_name(4) = "top";
487  }
488  else
489  {
490  MeshTools::Modification::change_boundary_id(other_mesh, 1, 5);
491  MeshTools::Modification::change_boundary_id(other_mesh, 2, 1);
492  MeshTools::Modification::change_boundary_id(other_mesh, 5, 2);
493  mesh.prepare_for_use();
494  other_mesh.prepare_for_use();
495  mesh.stitch_meshes(other_mesh, 2, 2, TOLERANCE, true);
496 
497  MeshTools::Modification::change_boundary_id(mesh, 3, 2);
498  mesh.get_boundary_info().sideset_name(1) = "left";
499  mesh.get_boundary_info().sideset_name(2) = "outer";
500  }
501  other_mesh.clear();
502  }
503  else if (_portion == "left_half")
504  {
505  ReplicatedMesh other_mesh(mesh);
506 
507  // This is to rotate the mesh and to reset boundary IDs.
508  MeshTools::Modification::rotate(other_mesh, 90, 0, 0);
509  MeshTools::Modification::rotate(mesh, 180, 0, 0);
510  if (_has_outer_square)
511  {
512  // The other mesh is created by rotating the original mesh about 90 degrees.
513  MeshTools::Modification::change_boundary_id(other_mesh, 1, 5);
514  MeshTools::Modification::change_boundary_id(other_mesh, 2, 6);
515  MeshTools::Modification::change_boundary_id(other_mesh, 3, 7);
516  MeshTools::Modification::change_boundary_id(other_mesh, 4, 1);
517  MeshTools::Modification::change_boundary_id(other_mesh, 5, 2);
518  MeshTools::Modification::change_boundary_id(other_mesh, 6, 3);
519  MeshTools::Modification::change_boundary_id(other_mesh, 7, 4);
520  // The original mesh is then rotated about 180 degrees.
521  MeshTools::Modification::change_boundary_id(mesh, 1, 5);
522  MeshTools::Modification::change_boundary_id(mesh, 2, 6);
523  MeshTools::Modification::change_boundary_id(mesh, 3, 7);
524  MeshTools::Modification::change_boundary_id(mesh, 4, 2);
525  MeshTools::Modification::change_boundary_id(mesh, 5, 3);
526  MeshTools::Modification::change_boundary_id(mesh, 6, 4);
527  MeshTools::Modification::change_boundary_id(mesh, 7, 1);
528  mesh.prepare_for_use();
529  other_mesh.prepare_for_use();
530  mesh.stitch_meshes(other_mesh, 4, 2, TOLERANCE, true);
531  mesh.get_boundary_info().sideset_name(1) = "left";
532  mesh.get_boundary_info().sideset_name(2) = "bottom";
533  mesh.get_boundary_info().sideset_name(3) = "right";
534  mesh.get_boundary_info().sideset_name(4) = "top";
535  }
536  else
537  {
538  MeshTools::Modification::change_boundary_id(mesh, 1, 5);
539  MeshTools::Modification::change_boundary_id(mesh, 2, 1);
540  MeshTools::Modification::change_boundary_id(mesh, 5, 2);
541  mesh.prepare_for_use();
542  other_mesh.prepare_for_use();
543  mesh.stitch_meshes(other_mesh, 1, 1, TOLERANCE, true);
544 
545  MeshTools::Modification::change_boundary_id(mesh, 2, 1);
546  MeshTools::Modification::change_boundary_id(mesh, 3, 2);
547  mesh.get_boundary_info().sideset_name(1) = "right";
548  mesh.get_boundary_info().sideset_name(2) = "outer";
549  }
550  other_mesh.clear();
551  }
552  else if (_portion == "bottom_half")
553  {
554  ReplicatedMesh other_mesh(mesh);
555  // This is to rotate the mesh and also to reset boundary IDs.
556  MeshTools::Modification::rotate(other_mesh, 180, 0, 0);
557  MeshTools::Modification::rotate(mesh, 270, 0, 0);
558  if (_has_outer_square)
559  {
560  // The other mesh is created by rotating the original mesh about 180 degrees.
561  MeshTools::Modification::change_boundary_id(other_mesh, 1, 5);
562  MeshTools::Modification::change_boundary_id(other_mesh, 2, 6);
563  MeshTools::Modification::change_boundary_id(other_mesh, 3, 7);
564  MeshTools::Modification::change_boundary_id(other_mesh, 4, 2);
565  MeshTools::Modification::change_boundary_id(other_mesh, 5, 3);
566  MeshTools::Modification::change_boundary_id(other_mesh, 6, 4);
567  MeshTools::Modification::change_boundary_id(other_mesh, 7, 1);
568  // The original mesh is rotated about 270 degrees.
569  MeshTools::Modification::change_boundary_id(mesh, 1, 5);
570  MeshTools::Modification::change_boundary_id(mesh, 2, 6);
571  MeshTools::Modification::change_boundary_id(mesh, 3, 7);
572  MeshTools::Modification::change_boundary_id(mesh, 4, 3);
573  MeshTools::Modification::change_boundary_id(mesh, 5, 4);
574  MeshTools::Modification::change_boundary_id(mesh, 6, 1);
575  MeshTools::Modification::change_boundary_id(mesh, 7, 2);
576  mesh.prepare_for_use();
577  other_mesh.prepare_for_use();
578  mesh.stitch_meshes(other_mesh, 1, 3, TOLERANCE, true);
579  mesh.get_boundary_info().sideset_name(1) = "left";
580  mesh.get_boundary_info().sideset_name(2) = "bottom";
581  mesh.get_boundary_info().sideset_name(3) = "right";
582  mesh.get_boundary_info().sideset_name(4) = "top";
583  }
584  else
585  {
586  MeshTools::Modification::change_boundary_id(mesh, 1, 5);
587  MeshTools::Modification::change_boundary_id(mesh, 2, 1);
588  MeshTools::Modification::change_boundary_id(mesh, 5, 2);
589  mesh.prepare_for_use();
590  other_mesh.prepare_for_use();
591  mesh.stitch_meshes(other_mesh, 1, 1, TOLERANCE, true);
592 
593  MeshTools::Modification::change_boundary_id(mesh, 2, 1);
594  MeshTools::Modification::change_boundary_id(mesh, 3, 2);
595  mesh.get_boundary_info().sideset_name(1) = "top";
596  mesh.get_boundary_info().sideset_name(2) = "outer";
597  }
598  other_mesh.clear();
599  }
600  else if (_portion == "full")
601  {
602  ReplicatedMesh portion_two(mesh);
603 
604  // This is to rotate the mesh and also to reset boundary IDs.
605  MeshTools::Modification::rotate(portion_two, 90, 0, 0);
606 
607  if (_has_outer_square)
608  {
609  // Portion 2: 2nd quadrant
610  MeshTools::Modification::change_boundary_id(portion_two, 1, 5);
611  MeshTools::Modification::change_boundary_id(portion_two, 2, 6);
612  MeshTools::Modification::change_boundary_id(portion_two, 3, 7);
613  MeshTools::Modification::change_boundary_id(portion_two, 4, 1);
614  MeshTools::Modification::change_boundary_id(portion_two, 5, 2);
615  MeshTools::Modification::change_boundary_id(portion_two, 6, 3);
616  MeshTools::Modification::change_boundary_id(portion_two, 7, 4);
617  mesh.prepare_for_use();
618  portion_two.prepare_for_use();
619  // 'top_half'
620  mesh.stitch_meshes(portion_two, 1, 3, TOLERANCE, true);
621 
622  // 'bottom_half'
623  ReplicatedMesh portion_bottom(mesh);
624  MeshTools::Modification::rotate(portion_bottom, 180, 0, 0);
625  MeshTools::Modification::change_boundary_id(portion_bottom, 1, 5);
626  MeshTools::Modification::change_boundary_id(portion_bottom, 2, 6);
627  MeshTools::Modification::change_boundary_id(portion_bottom, 3, 7);
628  MeshTools::Modification::change_boundary_id(portion_bottom, 4, 2);
629  MeshTools::Modification::change_boundary_id(portion_bottom, 5, 3);
630  MeshTools::Modification::change_boundary_id(portion_bottom, 6, 4);
631  MeshTools::Modification::change_boundary_id(portion_bottom, 7, 1);
632  mesh.prepare_for_use();
633  portion_bottom.prepare_for_use();
634  // 'full'
635  mesh.stitch_meshes(portion_bottom, 2, 4, TOLERANCE, true);
636 
637  mesh.get_boundary_info().sideset_name(1) = "left";
638  mesh.get_boundary_info().sideset_name(2) = "bottom";
639  mesh.get_boundary_info().sideset_name(3) = "right";
640  mesh.get_boundary_info().sideset_name(4) = "top";
641  portion_bottom.clear();
642  }
643  else
644  {
645  MeshTools::Modification::change_boundary_id(portion_two, 1, 5);
646  MeshTools::Modification::change_boundary_id(portion_two, 2, 1);
647  MeshTools::Modification::change_boundary_id(portion_two, 5, 2);
648  // 'top half'
649  mesh.prepare_for_use();
650  portion_two.prepare_for_use();
651  mesh.stitch_meshes(portion_two, 1, 1, TOLERANCE, true);
652  // 'bottom half'
653  ReplicatedMesh portion_bottom(mesh);
654  MeshTools::Modification::rotate(portion_bottom, 180, 0, 0);
655  // 'full'
656  mesh.prepare_for_use();
657  portion_bottom.prepare_for_use();
658  mesh.stitch_meshes(portion_bottom, 2, 2, TOLERANCE, true);
659  MeshTools::Modification::change_boundary_id(mesh, 3, 1);
660  mesh.get_boundary_info().sideset_name(1) = "outer";
661  portion_bottom.clear();
662  }
663  portion_two.clear();
664  }
665  if (_portion != "top_half" && _portion != "right_half" && _portion != "left_half" &&
666  _portion != "bottom_half" && _portion != "full")
667  mesh.prepare_for_use();
668 }
static InputParameters validParams()
Typical "Moose-style" constructor and copy constructor.
Definition: MooseMesh.C:78
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template * sin(_arg) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(tan
Real _inner_mesh_fraction
Size of inner square in relation to radius of the innermost concentric circle.
virtual void buildMesh() override
Must be overridden by child classes.
static InputParameters validParams()
unsigned int _num_sectors
Number of sectors in one quadrant.
MeshBase & mesh
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
std::vector< Real > _radii
Radii of concentric circles.
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template cos(_arg) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(cos
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:3198
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:88
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
ConcentricCircleMesh(const InputParameters &parameters)
bool _has_outer_square
Adding the moderator is optional.
std::vector< unsigned int > _rings
Number of rings in each circle or in the moderator.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
registerMooseObject("MooseApp", ConcentricCircleMesh)
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
Mesh generated from parameters.
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
void addRangeCheckedParam(const std::string &name, const T &value, const std::string &parsed_function, const std::string &doc_string)
virtual Elem * elem(const dof_id_type i)
Various accessors (pointers/references) for Elem "i".
Definition: MooseMesh.C:2849
bool _preserve_volumes
Volume preserving function is optional.
void ErrorVector unsigned int
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template pow< 2 >(tan(_arg))+1.0) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(sqrt
MooseEnum _portion
Control of which portion of mesh will be developed.