www.mooseframework.org
GeneratedMesh.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "GeneratedMesh.h"
16 
17 #include "libmesh/mesh_generation.h"
18 #include "libmesh/string_to_enum.h"
19 #include "libmesh/periodic_boundaries.h"
20 #include "libmesh/periodic_boundary_base.h"
21 #include "libmesh/unstructured_mesh.h"
22 
23 // C++ includes
24 #include <cmath> // provides round, not std::round (see http://www.cplusplus.com/reference/cmath/round/)
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");
38  "dim", dims, "The dimension of the mesh to be generated"); // Make this parameter required
39 
40  params.addParam<unsigned int>("nx", 1, "Number of elements in the X direction");
41  params.addParam<unsigned int>("ny", 1, "Number of elements in the Y direction");
42  params.addParam<unsigned int>("nz", 1, "Number of elements in the Z direction");
43  params.addParam<Real>("xmin", 0.0, "Lower X Coordinate of the generated mesh");
44  params.addParam<Real>("ymin", 0.0, "Lower Y Coordinate of the generated mesh");
45  params.addParam<Real>("zmin", 0.0, "Lower Z Coordinate of the generated mesh");
46  params.addParam<Real>("xmax", 1.0, "Upper X Coordinate of the generated mesh");
47  params.addParam<Real>("ymax", 1.0, "Upper Y Coordinate of the generated mesh");
48  params.addParam<Real>("zmax", 1.0, "Upper Z Coordinate of the generated mesh");
49  params.addParam<MooseEnum>("elem_type",
50  elem_types,
51  "The type of element from libMesh to "
52  "generate (default: linear element for "
53  "requested dimension)");
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.addParamNamesToGroup("dim", "Main");
75 
76  params.addClassDescription(
77  "Create a line, square, or cube mesh with uniformly spaced or biased elements.");
78  return params;
79 }
80 
82  : MooseMesh(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  // All generated meshes are regular orthogonal meshes
103 }
104 
105 Real
106 GeneratedMesh::getMinInDimension(unsigned int component) const
107 {
108  switch (component)
109  {
110  case 0:
111  return _xmin;
112  case 1:
113  return _dim > 1 ? _ymin : 0;
114  case 2:
115  return _dim > 2 ? _zmin : 0;
116  default:
117  mooseError("Invalid component");
118  }
119 }
120 
121 Real
122 GeneratedMesh::getMaxInDimension(unsigned int component) const
123 {
124  switch (component)
125  {
126  case 0:
127  return _xmax;
128  case 1:
129  return _dim > 1 ? _ymax : 0;
130  case 2:
131  return _dim > 2 ? _zmax : 0;
132  default:
133  mooseError("Invalid component");
134  }
135 }
136 
137 MooseMesh &
139 {
140  return *(new GeneratedMesh(*this));
141 }
142 
143 void
145 {
146  MooseEnum elem_type_enum = getParam<MooseEnum>("elem_type");
147 
148  if (!isParamValid("elem_type"))
149  {
150  // Switching on MooseEnum
151  switch (_dim)
152  {
153  case 1:
154  elem_type_enum = "EDGE2";
155  break;
156  case 2:
157  elem_type_enum = "QUAD4";
158  break;
159  case 3:
160  elem_type_enum = "HEX8";
161  break;
162  }
163  }
164 
165  ElemType elem_type = Utility::string_to_enum<ElemType>(elem_type_enum);
166 
167  // Switching on MooseEnum
168  switch (_dim)
169  {
170  // The build_XYZ mesh generation functions take an
171  // UnstructuredMesh& as the first argument, hence the dynamic_cast.
172  case 1:
173  MeshTools::Generation::build_line(dynamic_cast<UnstructuredMesh &>(getMesh()),
174  _nx,
175  _xmin,
176  _xmax,
177  elem_type,
179  break;
180  case 2:
181  MeshTools::Generation::build_square(dynamic_cast<UnstructuredMesh &>(getMesh()),
182  _nx,
183  _ny,
184  _xmin,
185  _xmax,
186  _ymin,
187  _ymax,
188  elem_type,
190  break;
191  case 3:
192  MeshTools::Generation::build_cube(dynamic_cast<UnstructuredMesh &>(getMesh()),
193  _nx,
194  _ny,
195  _nz,
196  _xmin,
197  _xmax,
198  _ymin,
199  _ymax,
200  _zmin,
201  _zmax,
202  elem_type,
204  break;
205  }
206 
207  // Apply the bias if any exists
208  if (_bias_x != 1.0 || _bias_y != 1.0 || _bias_z != 1.0)
209  {
210  // Reference to the libmesh mesh
211  MeshBase & mesh = getMesh();
212 
213  // Biases
214  Real bias[3] = {_bias_x, _bias_y, _bias_z};
215 
216  // "width" of the mesh in each direction
217  Real width[3] = {_xmax - _xmin, _ymax - _ymin, _zmax - _zmin};
218 
219  // Min mesh extent in each direction.
220  Real mins[3] = {_xmin, _ymin, _zmin};
221 
222  // Number of elements in each direction.
223  unsigned int nelem[3] = {_nx, _ny, _nz};
224 
225  // We will need the biases raised to integer powers in each
226  // direction, so let's pre-compute those...
227  std::vector<std::vector<Real>> pows(LIBMESH_DIM);
228  for (unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
229  {
230  pows[dir].resize(nelem[dir] + 1);
231  for (unsigned int i = 0; i < pows[dir].size(); ++i)
232  pows[dir][i] = std::pow(bias[dir], static_cast<int>(i));
233  }
234 
235  // Loop over the nodes and move them to the desired location
236  for (auto & node_ptr : mesh.node_ptr_range())
237  {
238  Node & node = *node_ptr;
239 
240  for (unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
241  {
242  if (width[dir] != 0. && bias[dir] != 1.)
243  {
244  // Compute the scaled "index" of the current point. This
245  // will either be close to a whole integer or a whole
246  // integer+0.5 for quadratic nodes.
247  Real float_index = (node(dir) - mins[dir]) * nelem[dir] / width[dir];
248 
249  Real integer_part = 0;
250  Real fractional_part = std::modf(float_index, &integer_part);
251 
252  // Figure out where to move the node...
253  if (std::abs(fractional_part) < TOLERANCE || std::abs(fractional_part - 1.0) < TOLERANCE)
254  {
255  // If the fractional_part ~ 0.0 or 1.0, this is a vertex node, so
256  // we don't need an average.
257  //
258  // Compute the "index" we are at in the current direction. We
259  // round to the nearest integral value to get this instead
260  // of using "integer_part", since that could be off by a
261  // lot (e.g. we want 3.9999 to map to 4.0 instead of 3.0).
262  int index = round(float_index);
263 
264  // Move node to biased location.
265  node(dir) =
266  mins[dir] + width[dir] * (1. - pows[dir][index]) / (1. - pows[dir][nelem[dir]]);
267  }
268  else if (std::abs(fractional_part - 0.5) < TOLERANCE)
269  {
270  // If the fractional_part ~ 0.5, this is a midedge/face
271  // (i.e. quadratic) node. We don't move those with the same
272  // bias as the vertices, instead we put them midway between
273  // their respective vertices.
274  //
275  // Also, since the fractional part is nearly 0.5, we know that
276  // the integer_part will be the index of the vertex to the
277  // left, and integer_part+1 will be the index of the
278  // vertex to the right.
279  node(dir) = mins[dir] +
280  width[dir] *
281  (1. - 0.5 * (pows[dir][integer_part] + pows[dir][integer_part + 1])) /
282  (1. - pows[dir][nelem[dir]]);
283  }
284  else
285  {
286  // We don't yet handle anything higher order than quadratic...
287  mooseError("Unable to bias node at node(", dir, ")=", node(dir));
288  }
289  }
290  }
291  }
292  }
293 }
virtual MooseMesh & clone() const override
Clone method.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
MooseEnum _dim
The dimension of the mesh.
Definition: GeneratedMesh.h:44
virtual const Node & node(const dof_id_type i) const
Various accessors (pointers/references) for Node "i".
Definition: MooseMesh.C:384
virtual Real getMaxInDimension(unsigned int component) const override
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...
virtual void buildMesh() override
Must be overridden by child classes.
InputParameters validParams< MooseMesh >()
Definition: MooseMesh.C:65
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:2355
InputParameters validParams< GeneratedMesh >()
Definition: GeneratedMesh.C:28
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
virtual Real getMinInDimension(unsigned int component) const override
Returns the min or max of the requested dimension respectively.
unsigned int _ny
Definition: GeneratedMesh.h:47
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:37
Real _xmin
The min/max values for x,y,z component.
Definition: GeneratedMesh.h:50
GeneratedMesh(const InputParameters &parameters)
Definition: GeneratedMesh.C:81
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseObject.h:67
unsigned int _nz
Definition: GeneratedMesh.h:47
bool _gauss_lobatto_grid
All of the libmesh build_line/square/cube routines support an option to grade the mesh into the bound...
Definition: GeneratedMesh.h:57
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...
unsigned int _nx
Number of elements in x, y, z direction.
Definition: GeneratedMesh.h:47
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 mooseError(Args &&...args) const
Definition: MooseObject.h:80
void addRangeCheckedParam(const std::string &name, const T &value, const std::string &parsed_function, const std::string &doc_string)
Real _bias_x
The amount by which to bias the cells in the x,y,z directions.
Definition: GeneratedMesh.h:64
bool _regular_orthogonal_mesh
Boolean indicating whether this mesh was detected to be regular and orthogonal.
Definition: MooseMesh.h:961
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...