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/string_to_enum.h"
16 #include "libmesh/periodic_boundaries.h"
17 #include "libmesh/periodic_boundary_base.h"
18 #include "libmesh/unstructured_mesh.h"
19 #include "libmesh/elem.h"
20 
21 // C++ includes
22 #include <cmath> // provides round, not std::round (see http://www.cplusplus.com/reference/cmath/round/)
23 
25 
26 template <>
29 {
31 
32  MooseEnum elem_types(
33  "EDGE EDGE2 EDGE3 EDGE4 QUAD QUAD4 QUAD8 QUAD9 TRI3 TRI6 HEX HEX8 HEX20 HEX27 TET4 TET10 "
34  "PRISM6 PRISM15 PRISM18 PYRAMID5 PYRAMID13 PYRAMID14"); // no default
35 
36  MooseEnum dims("1=1 2 3");
37  params.addRequiredParam<MooseEnum>("dim", dims, "The dimension of the mesh to be generated");
38 
39  params.addParam<unsigned int>("nx", 1, "Number of elements in the X direction");
40  params.addParam<unsigned int>("ny", 1, "Number of elements in the Y direction");
41  params.addParam<unsigned int>("nz", 1, "Number of elements in the Z direction");
42  params.addParam<Real>("xmin", 0.0, "Lower X Coordinate of the generated mesh");
43  params.addParam<Real>("ymin", 0.0, "Lower Y Coordinate of the generated mesh");
44  params.addParam<Real>("zmin", 0.0, "Lower Z Coordinate of the generated mesh");
45  params.addParam<Real>("xmax", 1.0, "Upper X Coordinate of the generated mesh");
46  params.addParam<Real>("ymax", 1.0, "Upper Y Coordinate of the generated mesh");
47  params.addParam<Real>("zmax", 1.0, "Upper Z Coordinate of the generated mesh");
48  params.addParam<MooseEnum>("elem_type",
49  elem_types,
50  "The type of element from libMesh to "
51  "generate (default: linear element for "
52  "requested dimension)");
53  params.addParam<bool>(
54  "gauss_lobatto_grid",
55  false,
56  "Grade mesh into boundaries according to Gauss-Lobatto quadrature spacing.");
57  params.addRangeCheckedParam<Real>(
58  "bias_x",
59  1.,
60  "bias_x>=0.5 & bias_x<=2",
61  "The amount by which to grow (or shrink) the cells in the x-direction.");
62  params.addRangeCheckedParam<Real>(
63  "bias_y",
64  1.,
65  "bias_y>=0.5 & bias_y<=2",
66  "The amount by which to grow (or shrink) the cells in the y-direction.");
67  params.addRangeCheckedParam<Real>(
68  "bias_z",
69  1.,
70  "bias_z>=0.5 & bias_z<=2",
71  "The amount by which to grow (or shrink) the cells in the z-direction.");
72 
73  params.addParamNamesToGroup("dim", "Main");
74 
75  params.addClassDescription(
76  "Create a line, square, or cube mesh with uniformly spaced or biased elements.");
77 
78  return params;
79 }
80 
82  : MeshGenerator(parameters),
83  _dim(getParam<MooseEnum>("dim")),
84  _nx(getParam<unsigned int>("nx")),
85  _ny(getParam<unsigned int>("ny")),
86  _nz(getParam<unsigned int>("nz")),
87  _xmin(getParam<Real>("xmin")),
88  _xmax(getParam<Real>("xmax")),
89  _ymin(getParam<Real>("ymin")),
90  _ymax(getParam<Real>("ymax")),
91  _zmin(getParam<Real>("zmin")),
92  _zmax(getParam<Real>("zmax")),
93  _gauss_lobatto_grid(getParam<bool>("gauss_lobatto_grid")),
94  _bias_x(getParam<Real>("bias_x")),
95  _bias_y(getParam<Real>("bias_y")),
96  _bias_z(getParam<Real>("bias_z"))
97 {
98  if (_gauss_lobatto_grid && (_bias_x != 1.0 || _bias_y != 1.0 || _bias_z != 1.0))
99  mooseError("Cannot apply both Gauss-Lobatto mesh grading and biasing at the same time.");
100 }
101 
102 std::unique_ptr<MeshBase>
104 {
105  // Have MOOSE construct the correct libMesh::Mesh object using Mesh block and CLI parameters.
106  auto mesh = _mesh->buildMeshBaseObject();
107 
108  MooseEnum elem_type_enum = getParam<MooseEnum>("elem_type");
109 
110  if (!isParamValid("elem_type"))
111  {
112  // Switching on MooseEnum
113  switch (_dim)
114  {
115  case 1:
116  elem_type_enum = "EDGE2";
117  break;
118  case 2:
119  elem_type_enum = "QUAD4";
120  break;
121  case 3:
122  elem_type_enum = "HEX8";
123  break;
124  }
125  }
126 
127  ElemType elem_type = Utility::string_to_enum<ElemType>(elem_type_enum);
128 
129  // Switching on MooseEnum
130  switch (_dim)
131  {
132  // The build_XYZ mesh generation functions take an
133  // UnstructuredMesh& as the first argument, hence the dynamic_cast.
134  case 1:
135  MeshTools::Generation::build_line(static_cast<UnstructuredMesh &>(*mesh),
136  _nx,
137  _xmin,
138  _xmax,
139  elem_type,
141  break;
142  case 2:
143  MeshTools::Generation::build_square(static_cast<UnstructuredMesh &>(*mesh),
144  _nx,
145  _ny,
146  _xmin,
147  _xmax,
148  _ymin,
149  _ymax,
150  elem_type,
152  break;
153  case 3:
154  MeshTools::Generation::build_cube(static_cast<UnstructuredMesh &>(*mesh),
155  _nx,
156  _ny,
157  _nz,
158  _xmin,
159  _xmax,
160  _ymin,
161  _ymax,
162  _zmin,
163  _zmax,
164  elem_type,
166  break;
167  }
168 
169  // Apply the bias if any exists
170  if (_bias_x != 1.0 || _bias_y != 1.0 || _bias_z != 1.0)
171  {
172  const auto MIN = std::numeric_limits<Real>::max();
173 
174  // Biases
175  std::array<Real, LIBMESH_DIM> bias = {
176  {_bias_x, _dim > 1 ? _bias_y : 1.0, _dim > 2 ? _bias_z : 1.0}};
177 
178  // "width" of the mesh in each direction
179  std::array<Real, LIBMESH_DIM> width = {
180  {_xmax - _xmin, _dim > 1 ? _ymax - _ymin : 0, _dim > 2 ? _zmax - _zmin : 0}};
181 
182  // Min mesh extent in each direction.
183  std::array<Real, LIBMESH_DIM> mins = {{_xmin, _dim > 1 ? _ymin : MIN, _dim > 2 ? _zmin : MIN}};
184 
185  // Number of elements in each direction.
186  std::array<unsigned int, LIBMESH_DIM> nelem = {{_nx, _dim > 1 ? _ny : 1, _dim > 2 ? _nz : 1}};
187 
188  // We will need the biases raised to integer powers in each
189  // direction, so let's pre-compute those...
190  std::array<std::vector<Real>, LIBMESH_DIM> pows;
191  for (unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
192  {
193  pows[dir].resize(nelem[dir] + 1);
194  for (unsigned int i = 0; i < pows[dir].size(); ++i)
195  pows[dir][i] = std::pow(bias[dir], static_cast<int>(i));
196  }
197 
198  // Loop over the nodes and move them to the desired location
199  for (auto & node_ptr : mesh->node_ptr_range())
200  {
201  Node & node = *node_ptr;
202 
203  for (unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
204  {
205  if (width[dir] != 0. && bias[dir] != 1.)
206  {
207  // Compute the scaled "index" of the current point. This
208  // will either be close to a whole integer or a whole
209  // integer+0.5 for quadratic nodes.
210  Real float_index = (node(dir) - mins[dir]) * nelem[dir] / width[dir];
211 
212  Real integer_part = 0;
213  Real fractional_part = std::modf(float_index, &integer_part);
214 
215  // Figure out where to move the node...
216  if (std::abs(fractional_part) < TOLERANCE || std::abs(fractional_part - 1.0) < TOLERANCE)
217  {
218  // If the fractional_part ~ 0.0 or 1.0, this is a vertex node, so
219  // we don't need an average.
220  //
221  // Compute the "index" we are at in the current direction. We
222  // round to the nearest integral value to get this instead
223  // of using "integer_part", since that could be off by a
224  // lot (e.g. we want 3.9999 to map to 4.0 instead of 3.0).
225  int index = round(float_index);
226 
227  mooseAssert(index >= static_cast<int>(0) && index < static_cast<int>(pows[dir].size()),
228  "Scaled \"index\" out of range");
229 
230  // Move node to biased location.
231  node(dir) =
232  mins[dir] + width[dir] * (1. - pows[dir][index]) / (1. - pows[dir][nelem[dir]]);
233  }
234  else if (std::abs(fractional_part - 0.5) < TOLERANCE)
235  {
236  // If the fractional_part ~ 0.5, this is a midedge/face
237  // (i.e. quadratic) node. We don't move those with the same
238  // bias as the vertices, instead we put them midway between
239  // their respective vertices.
240  //
241  // Also, since the fractional part is nearly 0.5, we know that
242  // the integer_part will be the index of the vertex to the
243  // left, and integer_part+1 will be the index of the
244  // vertex to the right.
245  node(dir) = mins[dir] +
246  width[dir] *
247  (1. - 0.5 * (pows[dir][integer_part] + pows[dir][integer_part + 1])) /
248  (1. - pows[dir][nelem[dir]]);
249  }
250  else
251  {
252  // We don't yet handle anything higher order than quadratic...
253  mooseError("Unable to bias node at node(", dir, ")=", node(dir));
254  }
255  }
256  }
257  }
258  }
259 
260  return dynamic_pointer_cast<MeshBase>(mesh);
261 }
unsigned int _nx
Number of elements in x, y, z direction.
MetaPhysicL::DualNumber< T, D > abs(const MetaPhysicL::DualNumber< T, D > &in)
void build_cube(UnstructuredMesh &mesh, const unsigned int nx, unsigned int ny, unsigned int nz, const Real xmin, const Real xmax, const Real ymin, const Real ymax, const Real zmin, const Real zmax, const ElemType type, bool verbose)
Real _xmin
The min/max values for x,y,z component.
Real round(Real x)
Definition: MathUtils.h:22
GeneratedMeshGenerator(const InputParameters &parameters)
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...
Real _bias_x
The amount by which to bias the cells in the x,y,z directions.
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.
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
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...
Real pow(Real x, int e)
Definition: MathUtils.C:211
registerMooseObject("MooseApp", GeneratedMeshGenerator)
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
bool _gauss_lobatto_grid
All of the libmesh build_line/square/cube routines support an option to grade the mesh into the bound...
InputParameters validParams< GeneratedMeshGenerator >()
std::shared_ptr< MooseMesh > & _mesh
References to the mesh and displaced mesh (currently in the ActionWarehouse)
Definition: MeshGenerator.h:72
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.
InputParameters validParams< MeshGenerator >()
Definition: MeshGenerator.C:16
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseObject.h:89
MeshGenerators are objects that can modify or add to an existing mesh.
Definition: MeshGenerator.h:30
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
This method takes a space delimited list of parameter names and adds them to the specified group name...