https://mooseframework.inl.gov
CartesianMeshGenerator.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
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 "CartesianMeshGenerator.h"
11 #include "CastUniquePointer.h"
12 
13 // libMesh includes
14 #include "libmesh/mesh_generation.h"
15 #include "libmesh/unstructured_mesh.h"
16 #include "libmesh/replicated_mesh.h"
17 #include "libmesh/point.h"
18 #include "libmesh/elem.h"
19 #include "libmesh/node.h"
20 
22 
25 {
27  MooseEnum dims("1=1 2 3");
28  params.addRequiredParam<MooseEnum>("dim", dims, "The dimension of the mesh to be generated");
29 
30  params.addRequiredParam<std::vector<Real>>("dx", "Intervals in the X direction");
31  params.addParam<std::vector<unsigned int>>(
32  "ix", "Number of grids in all intervals in the X direction (default to all one)");
33  params.addParam<std::vector<Real>>(
34  "dy", "Intervals in the Y direction (required when dim>1 otherwise ignored)");
35  params.addParam<std::vector<unsigned int>>(
36  "iy", "Number of grids in all intervals in the Y direction (default to all one)");
37  params.addParam<std::vector<Real>>(
38  "dz", "Intervals in the Z direction (required when dim>2 otherwise ignored)");
39  params.addParam<std::vector<unsigned int>>(
40  "iz", "Number of grids in all intervals in the Z direction (default to all one)");
41  params.addParam<std::vector<unsigned int>>("subdomain_id", "Block IDs (default to all zero)");
42  params.addClassDescription("This CartesianMeshGenerator creates a non-uniform Cartesian mesh.");
43  return params;
44 }
45 
47  : MeshGenerator(parameters),
48  _dim(getParam<MooseEnum>("dim")),
49  _dx(getParam<std::vector<Real>>("dx"))
50 {
51  // get all other parameters if provided and check their sizes
52  if (isParamValid("ix"))
53  {
54  _ix = getParam<std::vector<unsigned int>>("ix");
55  if (_ix.size() != _dx.size())
56  mooseError("ix must be in the same size of dx");
57  for (unsigned int i = 0; i < _ix.size(); ++i)
58  if (_ix[i] == 0)
59  mooseError("ix cannot be zero");
60  }
61  else
62  _ix = std::vector<unsigned int>(_dx.size(), 1);
63 
64  for (unsigned int i = 0; i < _dx.size(); ++i)
65  if (_dx[i] <= 0)
66  mooseError("dx must be greater than zero. dx: ", Moose::stringify(_dx));
67 
68  if (isParamValid("dy"))
69  {
70  _dy = getParam<std::vector<Real>>("dy");
71  for (unsigned int i = 0; i < _dy.size(); ++i)
72  if (_dy[i] <= 0)
73  mooseError("dy must be greater than zero. dy: ", Moose::stringify(_dy));
74  }
75 
76  if (isParamValid("iy"))
77  {
78  _iy = getParam<std::vector<unsigned int>>("iy");
79  if (_iy.size() != _dy.size())
80  mooseError("iy must be in the same size of dy");
81  for (unsigned int i = 0; i < _iy.size(); ++i)
82  if (_iy[i] == 0)
83  mooseError("iy cannot be zero");
84  }
85  else
86  _iy = std::vector<unsigned int>(_dy.size(), 1);
87 
88  if (isParamValid("dz"))
89  {
90  _dz = getParam<std::vector<Real>>("dz");
91  for (unsigned int i = 0; i < _dz.size(); ++i)
92  if (_dz[i] <= 0)
93  mooseError("dz must be greater than zero dz: ", Moose::stringify(_dz));
94  }
95 
96  if (isParamValid("iz"))
97  {
98  _iz = getParam<std::vector<unsigned int>>("iz");
99  if (_iz.size() != _dz.size())
100  mooseError("iz must be in the same size of dz");
101  for (unsigned int i = 0; i < _iz.size(); ++i)
102  if (_iz[i] == 0)
103  mooseError("iz cannot be zero");
104  }
105  else
106  _iz = std::vector<unsigned int>(_dz.size(), 1);
107 
108  if (isParamValid("subdomain_id"))
109  {
110  _subdomain_id = getParam<std::vector<unsigned int>>("subdomain_id");
111  if (isParamValid("dz") && isParamValid("dy"))
112  {
113  if (_subdomain_id.size() != _dx.size() * _dy.size() * _dz.size())
114  mooseError(
115  "subdomain_id must be in the size of product of sizes of dx, dy and dz. Sizes are '" +
116  std::to_string(_subdomain_id.size()) + "' and '" +
117  std::to_string(_dx.size() * _dy.size() * _dz.size()) + "' respectively");
118  }
119  else if (isParamValid("dy"))
120  {
121  if (_subdomain_id.size() != _dx.size() * _dy.size())
122  mooseError(
123  "subdomain_id must be in the size of product of sizes of dx and dy. Sizes are '" +
124  std::to_string(_subdomain_id.size()) + "' and '" +
125  std::to_string(_dx.size() * _dy.size()) + "' respectively");
126  }
127  else
128  {
129  if (_subdomain_id.size() != _dx.size())
130  mooseError("subdomain_id must be in the size of product of sizes of dx. Sizes are '" +
131  std::to_string(_subdomain_id.size()) + "' and '" + std::to_string(_dx.size()) +
132  "' respectively");
133  }
134  }
135  else
136  {
137  if (isParamValid("dz"))
138  _subdomain_id = std::vector<unsigned int>(_dx.size() * _dy.size() * _dz.size(), 0);
139  else if (isParamValid("dy"))
140  _subdomain_id = std::vector<unsigned int>(_dx.size() * _dy.size(), 0);
141  else
142  _subdomain_id = std::vector<unsigned int>(_dx.size(), 0);
143  }
144 
145  // do dimension checks and expand block IDs for all sub-grids
146  switch (_dim)
147  {
148  case 1:
149  {
150  _nx = 0;
151  for (unsigned int i = 0; i < _dx.size(); ++i)
152  _nx += _ix[i];
153  _ny = 1;
154  _nz = 1;
155 
156  std::vector<unsigned int> new_id;
157  for (unsigned int i = 0; i < _dx.size(); ++i)
158  for (unsigned int ii = 0; ii < _ix[i]; ++ii)
159  new_id.push_back(_subdomain_id[i]);
160 
161  _subdomain_id = new_id;
162 
163  if (isParamValid("dy"))
164  mooseWarning("dy provided for 1D");
165  if (isParamValid("iy"))
166  mooseWarning("iy provided for 1D");
167  if (isParamValid("dz"))
168  mooseWarning("dz provided for 1D");
169  if (isParamValid("iz"))
170  mooseWarning("iz provided for 1D");
171  break;
172  }
173  case 2:
174  {
175  _nx = 0;
176  for (unsigned int i = 0; i < _dx.size(); ++i)
177  _nx += _ix[i];
178  _ny = 0;
179  for (unsigned int i = 0; i < _dy.size(); ++i)
180  _ny += _iy[i];
181  _nz = 1;
182 
183  std::vector<unsigned int> new_id;
184  for (unsigned int j = 0; j < _dy.size(); ++j)
185  for (unsigned int jj = 0; jj < _iy[j]; ++jj)
186  for (unsigned int i = 0; i < _dx.size(); ++i)
187  for (unsigned int ii = 0; ii < _ix[i]; ++ii)
188  new_id.push_back(_subdomain_id[j * _dx.size() + i]);
189 
190  _subdomain_id = new_id;
191 
192  if (!isParamValid("dy"))
193  mooseError("dy is not provided for 2D");
194  if (isParamValid("dz"))
195  mooseWarning("dz provided for 2D");
196  if (isParamValid("iz"))
197  mooseWarning("iz provided for 2D");
198  break;
199  }
200  case 3:
201  {
202  _nx = 0;
203  for (unsigned int i = 0; i < _dx.size(); ++i)
204  _nx += _ix[i];
205  _ny = 0;
206  for (unsigned int i = 0; i < _dy.size(); ++i)
207  _ny += _iy[i];
208  _nz = 0;
209  for (unsigned int i = 0; i < _dz.size(); ++i)
210  _nz += _iz[i];
211 
212  std::vector<unsigned int> new_id;
213  for (unsigned int k = 0; k < _dz.size(); ++k)
214  for (unsigned int kk = 0; kk < _iz[k]; ++kk)
215  for (unsigned int j = 0; j < _dy.size(); ++j)
216  for (unsigned int jj = 0; jj < _iy[j]; ++jj)
217  for (unsigned int i = 0; i < _dx.size(); ++i)
218  for (unsigned int ii = 0; ii < _ix[i]; ++ii)
219  new_id.push_back(_subdomain_id[k * _dx.size() * _dy.size() + j * _dx.size() + i]);
220 
221  _subdomain_id = new_id;
222 
223  if (!isParamValid("dy"))
224  mooseError("dy is not provided for 3D");
225  if (!isParamValid("dz"))
226  mooseError("dz is not provided for 3D");
227  break;
228  }
229  }
230 }
231 
232 std::unique_ptr<MeshBase>
234 {
235  auto mesh = buildMeshBaseObject();
236 
237  // switching on MooseEnum to generate the reference mesh
238  // Note: element type is fixed
239  switch (_dim)
240  {
241  case 1:
242  MeshTools::Generation::build_line(static_cast<UnstructuredMesh &>(*mesh), _nx, 0, _nx, EDGE2);
243  break;
244  case 2:
245  MeshTools::Generation::build_square(
246  static_cast<UnstructuredMesh &>(*mesh), _nx, _ny, 0, _nx, 0, _ny, QUAD4);
247  break;
248  case 3:
249  MeshTools::Generation::build_cube(
250  static_cast<UnstructuredMesh &>(*mesh), _nx, _ny, _nz, 0, _nx, 0, _ny, 0, _nz, HEX8);
251  break;
252  }
253 
254  // assign block IDs
255  MeshBase::element_iterator el = mesh->active_elements_begin();
256  MeshBase::element_iterator el_end = mesh->active_elements_end();
257  for (; el != el_end; ++el)
258  {
259  const Point p = (*el)->vertex_average();
260  unsigned int ix = std::floor(p(0));
261  unsigned int iy = std::floor(p(1));
262  unsigned int iz = std::floor(p(2));
263  unsigned int i = iz * _nx * _ny + iy * _nx + ix;
264  (*el)->subdomain_id() = _subdomain_id[i];
265  }
266 
267  // adjust node coordinates
268  switch (_dim)
269  {
270  case 1:
271  {
272  Real base;
273 
274  std::vector<Real> mapx;
275  // Note: the starting coordinate is zero
276  base = 0;
277  mapx.push_back(base);
278  for (unsigned int i = 0; i < _dx.size(); ++i)
279  {
280  for (unsigned int j = 1; j <= _ix[i]; ++j)
281  mapx.push_back(base + _dx[i] / _ix[i] * j);
282  base += _dx[i];
283  }
284 
285  MeshBase::node_iterator node = mesh->active_nodes_begin();
286  MeshBase::node_iterator node_end = mesh->active_nodes_end();
287  for (; node != node_end; ++node)
288  {
289  unsigned int i = (*(*node))(0) + 0.5;
290  (*(*node))(0) = mapx.at(i);
291  }
292  break;
293  }
294  case 2:
295  {
296  Real base;
297 
298  std::vector<Real> mapx;
299  base = 0;
300  mapx.push_back(base);
301  for (unsigned int i = 0; i < _dx.size(); ++i)
302  {
303  for (unsigned int j = 1; j <= _ix[i]; ++j)
304  mapx.push_back(base + _dx[i] / _ix[i] * j);
305  base += _dx[i];
306  }
307 
308  std::vector<Real> mapy;
309  base = 0;
310  mapy.push_back(base);
311  for (unsigned int i = 0; i < _dy.size(); ++i)
312  {
313  for (unsigned int j = 1; j <= _iy[i]; ++j)
314  mapy.push_back(base + _dy[i] / _iy[i] * j);
315  base += _dy[i];
316  }
317 
318  MeshBase::node_iterator node = mesh->active_nodes_begin();
319  MeshBase::node_iterator node_end = mesh->active_nodes_end();
320  for (; node != node_end; ++node)
321  {
322  unsigned int i = (*(*node))(0) + 0.5;
323  (*(*node))(0) = mapx.at(i);
324  unsigned int j = (*(*node))(1) + 0.5;
325  (*(*node))(1) = mapy.at(j);
326  }
327  break;
328  }
329  case 3:
330  {
331  Real base;
332 
333  std::vector<Real> mapx;
334  base = 0;
335  mapx.push_back(base);
336  for (unsigned int i = 0; i < _dx.size(); ++i)
337  {
338  for (unsigned int j = 1; j <= _ix[i]; ++j)
339  mapx.push_back(base + _dx[i] / _ix[i] * j);
340  base += _dx[i];
341  }
342 
343  std::vector<Real> mapy;
344  base = 0;
345  mapy.push_back(base);
346  for (unsigned int i = 0; i < _dy.size(); ++i)
347  {
348  for (unsigned int j = 1; j <= _iy[i]; ++j)
349  mapy.push_back(base + _dy[i] / _iy[i] * j);
350  base += _dy[i];
351  }
352 
353  std::vector<Real> mapz;
354  base = 0;
355  mapz.push_back(base);
356  for (unsigned int i = 0; i < _dz.size(); ++i)
357  {
358  for (unsigned int j = 1; j <= _iz[i]; ++j)
359  mapz.push_back(base + _dz[i] / _iz[i] * j);
360  base += _dz[i];
361  }
362 
363  MeshBase::node_iterator node = mesh->active_nodes_begin();
364  MeshBase::node_iterator node_end = mesh->active_nodes_end();
365  for (; node != node_end; ++node)
366  {
367  unsigned int i = (*(*node))(0) + 0.5;
368  (*(*node))(0) = mapx.at(i);
369  unsigned int j = (*(*node))(1) + 0.5;
370  (*(*node))(1) = mapy.at(j);
371  unsigned int k = (*(*node))(2) + 0.5;
372  (*(*node))(2) = mapz.at(k);
373  }
374  break;
375  }
376  }
377 
378  mesh->prepare_for_use();
379  return dynamic_pointer_cast<MeshBase>(mesh);
380 }
MooseEnum _dim
The dimension of the mesh.
HEX8
std::unique_ptr< MeshBase > generate() override
Generate / modify the mesh.
CartesianMeshGenerator(const InputParameters &parameters)
MeshBase & mesh
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.
int _nx
Number of elements in x, y, z direction.
void mooseWarning(Args &&... args) const
Emits a warning prefixed with object name and type.
std::vector< Real > _dy
Intervals in y direction.
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...
QUAD4
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
std::vector< unsigned int > _ix
Number of grids in all intervals in x direction.
std::vector< unsigned int > _subdomain_id
Block IDs.
registerMooseObject("MooseApp", CartesianMeshGenerator)
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:33
std::vector< unsigned int > _iy
Number of grids in all intervals in y direction.
static InputParameters validParams()
Definition: MeshGenerator.C:23
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:64
std::vector< Real > _dx
Intervals in x direction.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
EDGE2
std::vector< unsigned int > _iz
Number of grids in all intervals in z direction.
std::vector< Real > _dz
Intervals in z direction.
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 optional parameter and a documentation string to the InputParameters object...
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
static InputParameters validParams()