www.mooseframework.org
GeneratedMesh.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 "GeneratedMesh.h"
11 
12 #include "libmesh/mesh_generation.h"
13 #include "libmesh/string_to_enum.h"
14 #include "libmesh/periodic_boundaries.h"
15 #include "libmesh/periodic_boundary_base.h"
16 #include "libmesh/unstructured_mesh.h"
17 #include "libmesh/node.h"
18 
19 // C++ includes
20 #include <cmath> // provides round, not std::round (see http://www.cplusplus.com/reference/cmath/round/)
21 #include <array>
22 
24 
25 template <>
28 {
30 
31  MooseEnum elem_types(
32  "EDGE EDGE2 EDGE3 EDGE4 QUAD QUAD4 QUAD8 QUAD9 TRI3 TRI6 HEX HEX8 HEX20 HEX27 TET4 TET10 "
33  "PRISM6 PRISM15 PRISM18 PYRAMID5 PYRAMID13 PYRAMID14"); // no default
34 
35  MooseEnum dims("1=1 2 3");
37  "dim", dims, "The dimension of the mesh to be generated"); // Make this parameter required
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  return params;
78 }
79 
81  : MooseMesh(parameters),
82  _dim(getParam<MooseEnum>("dim")),
83  _nx(getParam<unsigned int>("nx")),
84  _ny(getParam<unsigned int>("ny")),
85  _nz(getParam<unsigned int>("nz")),
86  _xmin(getParam<Real>("xmin")),
87  _xmax(getParam<Real>("xmax")),
88  _ymin(getParam<Real>("ymin")),
89  _ymax(getParam<Real>("ymax")),
90  _zmin(getParam<Real>("zmin")),
91  _zmax(getParam<Real>("zmax")),
92  _gauss_lobatto_grid(getParam<bool>("gauss_lobatto_grid")),
93  _bias_x(getParam<Real>("bias_x")),
94  _bias_y(getParam<Real>("bias_y")),
95  _bias_z(getParam<Real>("bias_z")),
96  _dims_may_have_changed(false)
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 - until they get modified ;)
103 }
104 
105 void
107 {
108  MooseMesh::prepared(state);
109 
110  // Fall back on scanning the mesh for coordinates instead of using input parameters for queries
111  if (!state)
112  _dims_may_have_changed = true;
113 }
114 
115 Real
116 GeneratedMesh::getMinInDimension(unsigned int component) const
117 {
119  return MooseMesh::getMinInDimension(component);
120 
121  switch (component)
122  {
123  case 0:
124  return _xmin;
125  case 1:
126  return _dim > 1 ? _ymin : 0;
127  case 2:
128  return _dim > 2 ? _zmin : 0;
129  default:
130  mooseError("Invalid component");
131  }
132 }
133 
134 Real
135 GeneratedMesh::getMaxInDimension(unsigned int component) const
136 {
138  return MooseMesh::getMaxInDimension(component);
139 
140  switch (component)
141  {
142  case 0:
143  return _xmax;
144  case 1:
145  return _dim > 1 ? _ymax : 0;
146  case 2:
147  return _dim > 2 ? _zmax : 0;
148  default:
149  mooseError("Invalid component");
150  }
151 }
152 
153 std::unique_ptr<MooseMesh>
155 {
156  return libmesh_make_unique<GeneratedMesh>(*this);
157 }
158 
159 void
161 {
162  MooseEnum elem_type_enum = getParam<MooseEnum>("elem_type");
163 
164  if (!isParamValid("elem_type"))
165  {
166  // Switching on MooseEnum
167  switch (_dim)
168  {
169  case 1:
170  elem_type_enum = "EDGE2";
171  break;
172  case 2:
173  elem_type_enum = "QUAD4";
174  break;
175  case 3:
176  elem_type_enum = "HEX8";
177  break;
178  }
179  }
180 
181  ElemType elem_type = Utility::string_to_enum<ElemType>(elem_type_enum);
182 
183  // Switching on MooseEnum
184  switch (_dim)
185  {
186  // The build_XYZ mesh generation functions take an
187  // UnstructuredMesh& as the first argument, hence the dynamic_cast.
188  case 1:
189  MeshTools::Generation::build_line(dynamic_cast<UnstructuredMesh &>(getMesh()),
190  _nx,
191  _xmin,
192  _xmax,
193  elem_type,
195  break;
196  case 2:
197  MeshTools::Generation::build_square(dynamic_cast<UnstructuredMesh &>(getMesh()),
198  _nx,
199  _ny,
200  _xmin,
201  _xmax,
202  _ymin,
203  _ymax,
204  elem_type,
206  break;
207  case 3:
208  MeshTools::Generation::build_cube(dynamic_cast<UnstructuredMesh &>(getMesh()),
209  _nx,
210  _ny,
211  _nz,
212  _xmin,
213  _xmax,
214  _ymin,
215  _ymax,
216  _zmin,
217  _zmax,
218  elem_type,
220  break;
221  }
222 
223  // Apply the bias if any exists
224  if (_bias_x != 1.0 || _bias_y != 1.0 || _bias_z != 1.0)
225  {
226  // Reference to the libmesh mesh
227  MeshBase & mesh = getMesh();
228 
229  // Biases
230  std::array<Real, LIBMESH_DIM> bias = {
231  {_bias_x, _dim > 1 ? _bias_y : 1.0, _dim > 2 ? _bias_z : 1.0}};
232 
233  // "width" of the mesh in each direction
234  std::array<Real, LIBMESH_DIM> width = {
235  {_xmax - _xmin, _dim > 1 ? _ymax - _ymin : 0, _dim > 2 ? _zmax - _zmin : 0}};
236 
237  // Min mesh extent in each direction.
238  std::array<Real, LIBMESH_DIM> mins = {
240 
241  // Number of elements in each direction.
242  std::array<unsigned int, LIBMESH_DIM> nelem = {{_nx, _dim > 1 ? _ny : 1, _dim > 2 ? _nz : 1}};
243 
244  // We will need the biases raised to integer powers in each
245  // direction, so let's pre-compute those...
246  std::array<std::vector<Real>, LIBMESH_DIM> pows;
247  for (unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
248  {
249  pows[dir].resize(nelem[dir] + 1);
250  for (unsigned int i = 0; i < pows[dir].size(); ++i)
251  pows[dir][i] = std::pow(bias[dir], static_cast<int>(i));
252  }
253 
254  // Loop over the nodes and move them to the desired location
255  for (auto & node_ptr : mesh.node_ptr_range())
256  {
257  Node & node = *node_ptr;
258 
259  for (unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
260  {
261  if (width[dir] != 0. && bias[dir] != 1.)
262  {
263  // Compute the scaled "index" of the current point. This
264  // will either be close to a whole integer or a whole
265  // integer+0.5 for quadratic nodes.
266  Real float_index = (node(dir) - mins[dir]) * nelem[dir] / width[dir];
267 
268  Real integer_part = 0;
269  Real fractional_part = std::modf(float_index, &integer_part);
270 
271  // Figure out where to move the node...
272  if (std::abs(fractional_part) < TOLERANCE || std::abs(fractional_part - 1.0) < TOLERANCE)
273  {
274  // If the fractional_part ~ 0.0 or 1.0, this is a vertex node, so
275  // we don't need an average.
276  //
277  // Compute the "index" we are at in the current direction. We
278  // round to the nearest integral value to get this instead
279  // of using "integer_part", since that could be off by a
280  // lot (e.g. we want 3.9999 to map to 4.0 instead of 3.0).
281  int index = round(float_index);
282 
283  mooseAssert(index >= static_cast<int>(0) && index < static_cast<int>(pows[dir].size()),
284  "Scaled \"index\" out of range");
285 
286  // Move node to biased location.
287  node(dir) =
288  mins[dir] + width[dir] * (1. - pows[dir][index]) / (1. - pows[dir][nelem[dir]]);
289  }
290  else if (std::abs(fractional_part - 0.5) < TOLERANCE)
291  {
292  // If the fractional_part ~ 0.5, this is a midedge/face
293  // (i.e. quadratic) node. We don't move those with the same
294  // bias as the vertices, instead we put them midway between
295  // their respective vertices.
296  //
297  // Also, since the fractional part is nearly 0.5, we know that
298  // the integer_part will be the index of the vertex to the
299  // left, and integer_part+1 will be the index of the
300  // vertex to the right.
301  node(dir) = mins[dir] +
302  width[dir] *
303  (1. - 0.5 * (pows[dir][integer_part] + pows[dir][integer_part + 1])) /
304  (1. - pows[dir][nelem[dir]]);
305  }
306  else
307  {
308  // We don't yet handle anything higher order than quadratic...
309  mooseError("Unable to bias node at node(", dir, ")=", node(dir));
310  }
311  }
312  }
313  }
314  }
315 }
virtual Real getMaxInDimension(unsigned int component) const
Definition: MooseMesh.C:1475
virtual Real getMinInDimension(unsigned int component) const
Returns the min or max of the requested dimension respectively.
Definition: MooseMesh.C:1466
bool prepared() const
Setter/getter for the _is_prepared flag.
Definition: MooseMesh.C:2291
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 round(Real x)
Definition: MathUtils.h:22
virtual std::unique_ptr< MooseMesh > safeClone() const override
A safer version of the clone() method that hands back an allocated object wrapped in a smart pointer...
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:40
virtual Real getMaxInDimension(unsigned int component) const override
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
bool _dims_may_have_changed
Boolean to indicate that dimensions may have changed.
Definition: GeneratedMesh.h:63
virtual Real getMinInDimension(unsigned int component) const override
Returns the min or max of the requested dimension respectively.
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.
Mesh generated from parameters.
Definition: GeneratedMesh.h:22
Real pow(Real x, int e)
Definition: MathUtils.C:211
InputParameters validParams< MooseMesh >()
Definition: MooseMesh.C:63
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:2567
InputParameters validParams< GeneratedMesh >()
Definition: GeneratedMesh.C:27
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
unsigned int _ny
Definition: GeneratedMesh.h:43
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
Real _xmin
The min/max values for x,y,z component.
Definition: GeneratedMesh.h:46
GeneratedMesh(const InputParameters &parameters)
Definition: GeneratedMesh.C:80
registerMooseObject("MooseApp", GeneratedMesh)
unsigned int _nz
Definition: GeneratedMesh.h:43
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:53
virtual const Node & node(const dof_id_type i) const
Various accessors (pointers/references) for Node "i".
Definition: MooseMesh.C:422
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:43
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)
Real _bias_x
The amount by which to bias the cells in the x,y,z directions.
Definition: GeneratedMesh.h:60
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseObject.h:89
bool _regular_orthogonal_mesh
Boolean indicating whether this mesh was detected to be regular and orthogonal.
Definition: MooseMesh.h:1009
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...