www.mooseframework.org
GeneratedMeshGenerator.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 "GeneratedMeshGenerator.h"
11 #include "CastUniquePointer.h"
12 
13 #include "libmesh/replicated_mesh.h"
14 #include "libmesh/mesh_generation.h"
15 #include "libmesh/mesh_modification.h"
16 #include "libmesh/string_to_enum.h"
17 #include "libmesh/periodic_boundaries.h"
18 #include "libmesh/periodic_boundary_base.h"
19 #include "libmesh/unstructured_mesh.h"
20 #include "libmesh/elem.h"
21 
22 // C++ includes
23 #include <cmath> // provides round, not std::round (see http://www.cplusplus.com/reference/cmath/round/)
24 
26 
29 {
31 
32  MooseEnum elem_types(LIST_GEOM_ELEM); // no default
33 
34  MooseEnum dims("1=1 2 3");
35  params.addRequiredParam<MooseEnum>("dim", dims, "The dimension of the mesh to be generated");
36 
37  params.addParam<unsigned int>("nx", 1, "Number of elements in the X direction");
38  params.addParam<unsigned int>("ny", 1, "Number of elements in the Y direction");
39  params.addParam<unsigned int>("nz", 1, "Number of elements in the Z direction");
40  params.addParam<Real>("xmin", 0.0, "Lower X Coordinate of the generated mesh");
41  params.addParam<Real>("ymin", 0.0, "Lower Y Coordinate of the generated mesh");
42  params.addParam<Real>("zmin", 0.0, "Lower Z Coordinate of the generated mesh");
43  params.addParam<Real>("xmax", 1.0, "Upper X Coordinate of the generated mesh");
44  params.addParam<Real>("ymax", 1.0, "Upper Y Coordinate of the generated mesh");
45  params.addParam<Real>("zmax", 1.0, "Upper Z Coordinate of the generated mesh");
46  params.addParam<MooseEnum>("elem_type",
47  elem_types,
48  "The type of element from libMesh to "
49  "generate (default: linear element for "
50  "requested dimension)");
51  params.addParam<std::vector<SubdomainID>>("subdomain_ids",
52  "Subdomain IDs for each element, default to all zero");
53 
54  params.addParam<bool>(
55  "gauss_lobatto_grid",
56  false,
57  "Grade mesh into boundaries according to Gauss-Lobatto quadrature spacing.");
58  params.addRangeCheckedParam<Real>(
59  "bias_x",
60  1.,
61  "bias_x>=0.5 & bias_x<=2",
62  "The amount by which to grow (or shrink) the cells in the x-direction.");
63  params.addRangeCheckedParam<Real>(
64  "bias_y",
65  1.,
66  "bias_y>=0.5 & bias_y<=2",
67  "The amount by which to grow (or shrink) the cells in the y-direction.");
68  params.addRangeCheckedParam<Real>(
69  "bias_z",
70  1.,
71  "bias_z>=0.5 & bias_z<=2",
72  "The amount by which to grow (or shrink) the cells in the z-direction.");
73 
74  params.addParam<std::string>("boundary_name_prefix",
75  "If provided, prefix the built in boundary names with this string");
76  params.addParam<boundary_id_type>(
77  "boundary_id_offset", 0, "This offset is added to the generated boundary IDs");
78 
79  params.addParam<std::vector<ExtraElementIDName>>("extra_element_integers",
80  "Names of extra element integers");
81 
82  params.addClassDescription(
83  "Create a line, square, or cube mesh with uniformly spaced or biased elements.");
84 
85  return params;
86 }
87 
89  : MeshGenerator(parameters),
90  _dim(getParam<MooseEnum>("dim")),
91  _nx(declareMeshProperty("num_elements_x", getParam<unsigned int>("nx"))),
92  _ny(declareMeshProperty("num_elements_y", getParam<unsigned int>("ny"))),
93  _nz(declareMeshProperty("num_elements_z", getParam<unsigned int>("nz"))),
94  _xmin(declareMeshProperty("xmin", getParam<Real>("xmin"))),
95  _xmax(declareMeshProperty("xmax", getParam<Real>("xmax"))),
96  _ymin(declareMeshProperty("ymin", getParam<Real>("ymin"))),
97  _ymax(declareMeshProperty("ymax", getParam<Real>("ymax"))),
98  _zmin(declareMeshProperty("zmin", getParam<Real>("zmin"))),
99  _zmax(declareMeshProperty("zmax", getParam<Real>("zmax"))),
100  _has_subdomain_ids(isParamValid("subdomain_ids")),
101  _gauss_lobatto_grid(getParam<bool>("gauss_lobatto_grid")),
102  _bias_x(getParam<Real>("bias_x")),
103  _bias_y(getParam<Real>("bias_y")),
104  _bias_z(getParam<Real>("bias_z")),
105  _boundary_name_prefix(isParamValid("boundary_name_prefix")
106  ? getParam<std::string>("boundary_name_prefix") + "_"
107  : ""),
108  _boundary_id_offset(getParam<boundary_id_type>("boundary_id_offset"))
109 {
110  if (_gauss_lobatto_grid && (_bias_x != 1.0 || _bias_y != 1.0 || _bias_z != 1.0))
111  mooseError("Cannot apply both Gauss-Lobatto mesh grading and biasing at the same time.");
112 }
113 
114 std::unique_ptr<MeshBase>
116 {
117  // Have MOOSE construct the correct libMesh::Mesh object using Mesh block and CLI parameters.
118  auto mesh = buildMeshBaseObject();
119 
120  if (isParamValid("extra_element_integers"))
121  {
122  for (auto & name : getParam<std::vector<ExtraElementIDName>>("extra_element_integers"))
123  mesh->add_elem_integer(name);
124  }
125 
126  MooseEnum elem_type_enum = getParam<MooseEnum>("elem_type");
127 
128  if (!isParamValid("elem_type"))
129  {
130  // Switching on MooseEnum
131  switch (_dim)
132  {
133  case 1:
134  elem_type_enum = "EDGE2";
135  break;
136  case 2:
137  elem_type_enum = "QUAD4";
138  break;
139  case 3:
140  elem_type_enum = "HEX8";
141  break;
142  }
143  }
144 
145  ElemType elem_type = Utility::string_to_enum<ElemType>(elem_type_enum);
146 
147  // Switching on MooseEnum
148  switch (_dim)
149  {
150  // The build_XYZ mesh generation functions take an
151  // UnstructuredMesh& as the first argument, hence the static_cast.
152  case 1:
153  MeshTools::Generation::build_line(static_cast<UnstructuredMesh &>(*mesh),
154  _nx,
155  _xmin,
156  _xmax,
157  elem_type,
159  break;
160  case 2:
161  MeshTools::Generation::build_square(static_cast<UnstructuredMesh &>(*mesh),
162  _nx,
163  _ny,
164  _xmin,
165  _xmax,
166  _ymin,
167  _ymax,
168  elem_type,
170  break;
171  case 3:
172  MeshTools::Generation::build_cube(static_cast<UnstructuredMesh &>(*mesh),
173  _nx,
174  _ny,
175  _nz,
176  _xmin,
177  _xmax,
178  _ymin,
179  _ymax,
180  _zmin,
181  _zmax,
182  elem_type,
184  break;
185  }
186 
187  if (_has_subdomain_ids)
188  {
189  auto & bids = getParam<std::vector<SubdomainID>>("subdomain_ids");
190  if (bids.size() != _nx * _ny * _nz)
191  paramError("subdomain_ids",
192  "Size must equal to the product of number of elements in all directions");
193  for (auto & elem : mesh->element_ptr_range())
194  {
195  const Point p = elem->vertex_average();
196  unsigned int ix = std::floor((p(0) - _xmin) / (_xmax - _xmin) * _nx);
197  unsigned int iy = std::floor((p(1) - _ymin) / (_ymax - _ymin) * _ny);
198  unsigned int iz = std::floor((p(2) - _zmin) / (_zmax - _zmin) * _nz);
199  unsigned int i = iz * _nx * _ny + iy * _nx + ix;
200  elem->subdomain_id() = bids[i];
201  }
202  }
203 
204  // rename and shift boundaries
205  BoundaryInfo & boundary_info = mesh->get_boundary_info();
206 
207  // Copy, since we're modifying the container mid-iteration
208  const auto mesh_boundary_ids = boundary_info.get_global_boundary_ids();
209  for (auto rit = mesh_boundary_ids.rbegin(); rit != mesh_boundary_ids.rend(); ++rit)
210  {
211  const std::string old_sideset_name = boundary_info.sideset_name(*rit);
212  const std::string old_nodeset_name = boundary_info.nodeset_name(*rit);
213 
214  MeshTools::Modification::change_boundary_id(*mesh, *rit, *rit + _boundary_id_offset);
215 
216  boundary_info.sideset_name(*rit + _boundary_id_offset) =
217  _boundary_name_prefix + old_sideset_name;
218  boundary_info.nodeset_name(*rit + _boundary_id_offset) =
219  _boundary_name_prefix + old_nodeset_name;
220  }
221 
222  // Apply the bias if any exists
223  if (_bias_x != 1.0 || _bias_y != 1.0 || _bias_z != 1.0)
224  {
225  const auto MIN = std::numeric_limits<Real>::max();
226 
227  // Biases
228  std::array<Real, LIBMESH_DIM> bias = {
229  {_bias_x, _dim > 1 ? _bias_y : 1.0, _dim > 2 ? _bias_z : 1.0}};
230 
231  // "width" of the mesh in each direction
232  std::array<Real, LIBMESH_DIM> width = {
233  {_xmax - _xmin, _dim > 1 ? _ymax - _ymin : 0, _dim > 2 ? _zmax - _zmin : 0}};
234 
235  // Min mesh extent in each direction.
236  std::array<Real, LIBMESH_DIM> mins = {{_xmin, _dim > 1 ? _ymin : MIN, _dim > 2 ? _zmin : MIN}};
237 
238  // Number of elements in each direction.
239  std::array<unsigned int, LIBMESH_DIM> nelem = {{_nx, _dim > 1 ? _ny : 1, _dim > 2 ? _nz : 1}};
240 
241  // We will need the biases raised to integer powers in each
242  // direction, so let's pre-compute those...
243  std::array<std::vector<Real>, LIBMESH_DIM> pows;
244  for (unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
245  {
246  pows[dir].resize(nelem[dir] + 1);
247  pows[dir][0] = 1.0;
248  for (unsigned int i = 1; i < pows[dir].size(); ++i)
249  pows[dir][i] = pows[dir][i - 1] * bias[dir];
250  }
251 
252  // Loop over the nodes and move them to the desired location
253  for (auto & node_ptr : mesh->node_ptr_range())
254  {
255  Node & node = *node_ptr;
256 
257  for (unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
258  {
259  if (width[dir] != 0. && bias[dir] != 1.)
260  {
261  // Compute the scaled "index" of the current point. This
262  // will either be close to a whole integer or a whole
263  // integer+0.5 for quadratic nodes.
264  Real float_index = (node(dir) - mins[dir]) * nelem[dir] / width[dir];
265 
266  Real integer_part = 0;
267  Real fractional_part = std::modf(float_index, &integer_part);
268 
269  // Figure out where to move the node...
270  if (std::abs(fractional_part) < TOLERANCE || std::abs(fractional_part - 1.0) < TOLERANCE)
271  {
272  // If the fractional_part ~ 0.0 or 1.0, this is a vertex node, so
273  // we don't need an average.
274  //
275  // Compute the "index" we are at in the current direction. We
276  // round to the nearest integral value to get this instead
277  // of using "integer_part", since that could be off by a
278  // lot (e.g. we want 3.9999 to map to 4.0 instead of 3.0).
279  int index = round(float_index);
280 
281  mooseAssert(index >= static_cast<int>(0) && index < static_cast<int>(pows[dir].size()),
282  "Scaled \"index\" out of range");
283 
284  // Move node to biased location.
285  node(dir) =
286  mins[dir] + width[dir] * (1. - pows[dir][index]) / (1. - pows[dir][nelem[dir]]);
287  }
288  else if (std::abs(fractional_part - 0.5) < TOLERANCE)
289  {
290  // If the fractional_part ~ 0.5, this is a midedge/face
291  // (i.e. quadratic) node. We don't move those with the same
292  // bias as the vertices, instead we put them midway between
293  // their respective vertices.
294  //
295  // Also, since the fractional part is nearly 0.5, we know that
296  // the integer_part will be the index of the vertex to the
297  // left, and integer_part+1 will be the index of the
298  // vertex to the right.
299  node(dir) = mins[dir] +
300  width[dir] *
301  (1. - 0.5 * (pows[dir][integer_part] + pows[dir][integer_part + 1])) /
302  (1. - pows[dir][nelem[dir]]);
303  }
304  else
305  {
306  // We don't yet handle anything higher order than quadratic...
307  mooseError("Unable to bias node at node(", dir, ")=", node(dir));
308  }
309  }
310  }
311  }
312  }
313 
314  mesh->set_isnt_prepared();
315  return dynamic_pointer_cast<MeshBase>(mesh);
316 }
const Real _bias_x
The amount by which to bias the cells in the x,y,z directions.
ElemType
void elem_types(const MeshBase &mesh, std::vector< ElemType > &et)
const std::string LIST_GEOM_ELEM
Definition: MooseMesh.h:58
bool _has_subdomain_ids
Whether or not subdomain_ids parameter is set.
GeneratedMeshGenerator(const InputParameters &parameters)
MeshBase & mesh
Generates a line, square, or cube mesh with uniformly spaced or biased elements.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
std::unique_ptr< T_DEST, T_DELETER > dynamic_pointer_cast(std::unique_ptr< T_SRC, T_DELETER > &src)
These are reworked from https://stackoverflow.com/a/11003103.
const boundary_id_type _boundary_id_offset
offset that is added to the boundary IDs
virtual const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:56
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...
auto max(const L &left, const R &right)
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
int8_t boundary_id_type
static InputParameters validParams()
T round(T x)
Definition: MathUtils.h:76
registerMooseObject("MooseApp", GeneratedMeshGenerator)
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
unsigned int & _nx
Number of elements in x, y, z direction.
Real & _xmin
The min/max values for x,y,z component.
void paramError(const std::string &param, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
static InputParameters validParams()
Definition: MeshGenerator.C:23
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const bool _gauss_lobatto_grid
All of the libmesh build_line/square/cube routines support an option to grade the mesh into the bound...
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
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...
std::unique_ptr< MeshBase > generate() override
Generate / modify the mesh.
void addRangeCheckedParam(const std::string &name, const T &value, const std::string &parsed_function, const std::string &doc_string)
MooseEnum _dim
The dimension of the mesh.
const std::string _boundary_name_prefix
prefix string for the boundary names
std::unique_ptr< MeshBase > buildMeshBaseObject(unsigned int dim=libMesh::invalid_uint)
Build a MeshBase object whose underlying type will be determined by the Mesh input file block...
MeshGenerators are objects that can modify or add to an existing mesh.
Definition: MeshGenerator.h:32
void ErrorVector unsigned int