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" 34 "dim", dims,
"The dimension of the mesh to be generated");
36 params.
addParam<
unsigned int>(
"nx", 1,
"Number of elements in the X direction");
37 params.
addParam<
unsigned int>(
"ny", 1,
"Number of elements in the Y direction");
38 params.
addParam<
unsigned int>(
"nz", 1,
"Number of elements in the Z direction");
39 params.
addParam<
Real>(
"xmin", 0.0,
"Lower X Coordinate of the generated mesh");
40 params.
addParam<
Real>(
"ymin", 0.0,
"Lower Y Coordinate of the generated mesh");
41 params.
addParam<
Real>(
"zmin", 0.0,
"Lower Z Coordinate of the generated mesh");
42 params.
addParam<
Real>(
"xmax", 1.0,
"Upper X Coordinate of the generated mesh");
43 params.
addParam<
Real>(
"ymax", 1.0,
"Upper Y Coordinate of the generated mesh");
44 params.
addParam<
Real>(
"zmax", 1.0,
"Upper Z Coordinate of the generated mesh");
47 "The type of element from libMesh to " 48 "generate (default: linear element for " 49 "requested dimension)");
53 "Grade mesh into boundaries according to Gauss-Lobatto quadrature spacing.");
57 "bias_x>=0.5 & bias_x<=2",
58 "The amount by which to grow (or shrink) the cells in the x-direction.");
62 "bias_y>=0.5 & bias_y<=2",
63 "The amount by which to grow (or shrink) the cells in the y-direction.");
67 "bias_z>=0.5 & bias_z<=2",
68 "The amount by which to grow (or shrink) the cells in the z-direction.");
73 "Create a line, square, or cube mesh with uniformly spaced or biased elements.");
80 _nx(getParam<unsigned
int>(
"nx")),
81 _ny(getParam<unsigned
int>(
"ny")),
82 _nz(getParam<unsigned
int>(
"nz")),
83 _xmin(getParam<
Real>(
"xmin")),
84 _xmax(getParam<
Real>(
"xmax")),
85 _ymin(getParam<
Real>(
"ymin")),
86 _ymax(getParam<
Real>(
"ymax")),
87 _zmin(getParam<
Real>(
"zmin")),
88 _zmax(getParam<
Real>(
"zmax")),
89 _gauss_lobatto_grid(getParam<bool>(
"gauss_lobatto_grid")),
90 _bias_x(getParam<
Real>(
"bias_x")),
91 _bias_y(getParam<
Real>(
"bias_y")),
92 _bias_z(getParam<
Real>(
"bias_z")),
93 _dims_may_have_changed(false)
96 mooseError(
"Cannot apply both Gauss-Lobatto mesh grading and biasing at the same time.");
150 std::unique_ptr<MooseMesh>
153 return std::make_unique<GeneratedMesh>(*this);
159 MooseEnum elem_type_enum = getParam<MooseEnum>(
"elem_type");
167 elem_type_enum =
"EDGE2";
170 elem_type_enum =
"QUAD4";
173 elem_type_enum =
"HEX8";
178 ElemType elem_type = Utility::string_to_enum<ElemType>(elem_type_enum);
186 MeshTools::Generation::build_line(dynamic_cast<UnstructuredMesh &>(
getMesh()),
194 MeshTools::Generation::build_square(dynamic_cast<UnstructuredMesh &>(
getMesh()),
205 MeshTools::Generation::build_cube(dynamic_cast<UnstructuredMesh &>(
getMesh()),
227 std::array<Real, LIBMESH_DIM> bias = {
231 std::array<Real, LIBMESH_DIM> width = {
235 std::array<Real, LIBMESH_DIM> mins = {
239 std::array<unsigned int, LIBMESH_DIM> nelem = {{
_nx,
_dim > 1 ?
_ny : 1,
_dim > 2 ?
_nz : 1}};
243 std::array<std::vector<Real>, LIBMESH_DIM> pows;
244 for (
unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
246 pows[dir].resize(nelem[dir] + 1);
247 for (
unsigned int i = 0; i < pows[dir].size(); ++i)
248 pows[dir][i] =
std::pow(bias[dir], static_cast<int>(i));
252 for (
auto & node_ptr :
mesh.node_ptr_range())
254 Node &
node = *node_ptr;
256 for (
unsigned int dir = 0; dir < LIBMESH_DIM; ++dir)
258 if (width[dir] != 0. && bias[dir] != 1.)
263 Real float_index = (
node(dir) - mins[dir]) * nelem[dir] / width[dir];
265 Real integer_part = 0;
266 Real fractional_part = std::modf(float_index, &integer_part);
269 if (
std::abs(fractional_part) < TOLERANCE ||
std::abs(fractional_part - 1.0) < TOLERANCE)
278 int index =
round(float_index);
280 mooseAssert(index >= static_cast<int>(0) && index < static_cast<int>(pows[dir].size()),
281 "Scaled \"index\" out of range");
285 mins[dir] + width[dir] * (1. - pows[dir][index]) / (1. - pows[dir][nelem[dir]]);
287 else if (
std::abs(fractional_part - 0.5) < TOLERANCE)
298 node(dir) = mins[dir] +
300 (1. - 0.5 * (pows[dir][integer_part] + pows[dir][integer_part + 1])) /
301 (1. - pows[dir][nelem[dir]]);
static InputParameters validParams()
Typical "Moose-style" constructor and copy constructor.
virtual Real getMaxInDimension(unsigned int component) const
virtual Real getMinInDimension(unsigned int component) const
Returns the min or max of the requested dimension respectively.
bool prepared() const
Setter/getter for whether the mesh is prepared.
const std::string LIST_GEOM_ELEM
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...
MooseEnum _dim
The dimension of the mesh.
virtual Real getMaxInDimension(unsigned int component) const override
bool _dims_may_have_changed
Boolean to indicate that dimensions may have changed.
virtual Real getMinInDimension(unsigned int component) const override
Returns the min or max of the requested dimension respectively.
virtual void buildMesh() override
Must be overridden by child classes.
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
Mesh generated from parameters.
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Real _xmin
The min/max values for x,y,z component.
GeneratedMesh(const InputParameters ¶meters)
registerMooseObject("MooseApp", GeneratedMesh)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool _gauss_lobatto_grid
All of the libmesh build_line/square/cube routines support an option to grade the mesh into the bound...
virtual const Node & node(const dof_id_type i) const
Various accessors (pointers/references) for Node "i".
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
unsigned int _nx
Number of elements in x, y, z direction.
Real _bias_x
The amount by which to bias the cells in the x,y,z directions.
MooseUnits pow(const MooseUnits &, int)
void ErrorVector unsigned int
bool _regular_orthogonal_mesh
Boolean indicating whether this mesh was detected to be regular and orthogonal.
static InputParameters validParams()