www.mooseframework.org
DistributedGeneratedMesh.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 
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/partitioner.h"
18 #include "libmesh/metis_csr_graph.h"
19 #include "libmesh/edge_edge2.h"
20 #include "libmesh/edge_edge3.h"
21 #include "libmesh/edge_edge4.h"
22 #include "libmesh/mesh_communication.h"
23 #include "libmesh/remote_elem.h"
24 #include "libmesh/face_quad4.h"
25 #include "libmesh/cell_hex8.h"
26 
27 // C++ includes
28 #include <cmath> // provides round, not std::round (see http://www.cplusplus.com/reference/cmath/round/)
29 
30 #ifdef LIBMESH_HAVE_METIS
31 // MIPSPro 7.4.2 gets confused about these nested namespaces
32 #ifdef __sgi
33 #include <cstdarg>
34 #endif
35 namespace Metis
36 {
37 extern "C"
38 {
39 #include "libmesh/ignore_warnings.h"
40 #include "metis.h"
41 #include "libmesh/restore_warnings.h"
42 }
43 }
44 #else
45 #include "libmesh/sfc_partitioner.h"
46 #endif
47 
49 
50 template <>
53 {
55 
56  params.addParam<bool>("verbose", false, "Turn on verbose printing for the mesh generation");
57 
58  MooseEnum elem_types(
59  "EDGE EDGE2 EDGE3 EDGE4 QUAD QUAD4 QUAD8 QUAD9 TRI3 TRI6 HEX HEX8 HEX20 HEX27 TET4 TET10 "
60  "PRISM6 PRISM15 PRISM18 PYRAMID5 PYRAMID13 PYRAMID14"); // no default
61 
62  MooseEnum dims("1=1 2 3");
64  "dim", dims, "The dimension of the mesh to be generated"); // Make this parameter required
65 
66  params.addParam<dof_id_type>("nx", 1, "Number of elements in the X direction");
67  params.addParam<dof_id_type>("ny", 1, "Number of elements in the Y direction");
68  params.addParam<dof_id_type>("nz", 1, "Number of elements in the Z direction");
69  params.addParam<Real>("xmin", 0.0, "Lower X Coordinate of the generated mesh");
70  params.addParam<Real>("ymin", 0.0, "Lower Y Coordinate of the generated mesh");
71  params.addParam<Real>("zmin", 0.0, "Lower Z Coordinate of the generated mesh");
72  params.addParam<Real>("xmax", 1.0, "Upper X Coordinate of the generated mesh");
73  params.addParam<Real>("ymax", 1.0, "Upper Y Coordinate of the generated mesh");
74  params.addParam<Real>("zmax", 1.0, "Upper Z Coordinate of the generated mesh");
75  params.addParam<MooseEnum>("elem_type",
76  elem_types,
77  "The type of element from libMesh to "
78  "generate (default: linear element for "
79  "requested dimension)");
80  params.addRangeCheckedParam<Real>(
81  "bias_x",
82  1.,
83  "bias_x>=0.5 & bias_x<=2",
84  "The amount by which to grow (or shrink) the cells in the x-direction.");
85  params.addRangeCheckedParam<Real>(
86  "bias_y",
87  1.,
88  "bias_y>=0.5 & bias_y<=2",
89  "The amount by which to grow (or shrink) the cells in the y-direction.");
90  params.addRangeCheckedParam<Real>(
91  "bias_z",
92  1.,
93  "bias_z>=0.5 & bias_z<=2",
94  "The amount by which to grow (or shrink) the cells in the z-direction.");
95 
96  params.addParamNamesToGroup("dim", "Main");
97 
98  params.addClassDescription(
99  "Create a line, square, or cube mesh with uniformly spaced or biased elements.");
100 
101  // This mesh is always distributed
102  params.set<MooseEnum>("parallel_type") = "DISTRIBUTED";
103 
104  return params;
105 }
106 
108  : MooseMesh(parameters),
109  _verbose(getParam<bool>("verbose")),
110  _dim(getParam<MooseEnum>("dim")),
111  _nx(getParam<dof_id_type>("nx")),
112  _ny(getParam<dof_id_type>("ny")),
113  _nz(getParam<dof_id_type>("nz")),
114  _xmin(getParam<Real>("xmin")),
115  _xmax(getParam<Real>("xmax")),
116  _ymin(getParam<Real>("ymin")),
117  _ymax(getParam<Real>("ymax")),
118  _zmin(getParam<Real>("zmin")),
119  _zmax(getParam<Real>("zmax")),
120  _bias_x(getParam<Real>("bias_x")),
121  _bias_y(getParam<Real>("bias_y")),
122  _bias_z(getParam<Real>("bias_z")),
123  _dims_may_have_changed(false)
124 {
125  // All generated meshes are regular orthogonal meshes
127 }
128 
129 void
131 {
132  MooseMesh::prepared(state);
133 
134  // Fall back on scanning the mesh for coordinates instead of using input parameters for queries
135  if (!state)
136  _dims_may_have_changed = true;
137 }
138 
139 Real
140 DistributedGeneratedMesh::getMinInDimension(unsigned int component) const
141 {
143  return MooseMesh::getMinInDimension(component);
144 
145  switch (component)
146  {
147  case 0:
148  return _xmin;
149  case 1:
150  return _dim > 1 ? _ymin : 0;
151  case 2:
152  return _dim > 2 ? _zmin : 0;
153  default:
154  mooseError("Invalid component");
155  }
156 }
157 
158 Real
159 DistributedGeneratedMesh::getMaxInDimension(unsigned int component) const
160 {
162  return MooseMesh::getMaxInDimension(component);
163 
164  switch (component)
165  {
166  case 0:
167  return _xmax;
168  case 1:
169  return _dim > 1 ? _ymax : 0;
170  case 2:
171  return _dim > 2 ? _zmax : 0;
172  default:
173  mooseError("Invalid component");
174  }
175 }
176 
177 std::unique_ptr<MooseMesh>
179 {
180  return libmesh_make_unique<DistributedGeneratedMesh>(*this);
181 }
182 
183 namespace
184 {
185 
196 template <typename T>
197 inline dof_id_type
198 elem_id(const dof_id_type /*nx*/,
199  const dof_id_type /*ny*/,
200  const dof_id_type /*i*/,
201  const dof_id_type /*j*/,
202  const dof_id_type /*k*/)
203 {
204  mooseError("elem_id not implemented for this element type in DistributedGeneratedMesh");
205 }
206 
218 template <typename T>
219 inline dof_id_type
220 num_neighbors(const dof_id_type /*nx*/,
221  const dof_id_type /*ny*/,
222  const dof_id_type /*nz*/,
223  const dof_id_type /*i*/,
224  const dof_id_type /*j*/,
225  const dof_id_type /*k*/)
226 {
227  mooseError("num_neighbors not implemented for this element type in DistributedGeneratedMesh");
228 }
229 
242 template <typename T>
243 inline void
244 get_neighbors(const dof_id_type /*nx*/,
245  const dof_id_type /*ny*/,
246  const dof_id_type /*nz*/,
247  const dof_id_type /*i*/,
248  const dof_id_type /*j*/,
249  const dof_id_type /*k*/,
250  std::vector<dof_id_type> & /*neighbors*/)
251 {
252  mooseError("get_neighbors not implemented for this element type in DistributedGeneratedMesh");
253 }
254 
266 template <typename T>
267 inline dof_id_type
268 node_id(const ElemType /*type*/,
269  const dof_id_type /*nx*/,
270  const dof_id_type /*ny*/,
271  const dof_id_type /*i*/,
272  const dof_id_type /*j*/,
273  const dof_id_type /*k*/)
274 {
275  mooseError("node_id not implemented for this element type in DistributedGeneratedMesh");
276 }
277 
290 template <typename T>
291 Node *
292 add_point(const dof_id_type /*nx*/,
293  const dof_id_type /*ny*/,
294  const dof_id_type /*nz*/,
295  const dof_id_type /*i*/,
296  const dof_id_type /*j*/,
297  const dof_id_type /*k*/,
298  const ElemType /*type*/,
299  MeshBase & /*mesh*/)
300 {
301  mooseError("add_point not implemented for this element type in DistributedGeneratedMesh");
302 }
303 
319 template <typename T>
320 void
321 add_element(const dof_id_type /*nx*/,
322  const dof_id_type /*ny*/,
323  const dof_id_type /*nz*/,
324  const dof_id_type /*i*/,
325  const dof_id_type /*j*/,
326  const dof_id_type /*k*/,
327  const dof_id_type /*elem_id*/,
328  const processor_id_type /*pid*/,
329  const ElemType /*type*/,
330  MeshBase & /*mesh*/,
331  bool /*verbose*/)
332 {
333  mooseError("add_element not implemented for this element type in DistributedGeneratedMesh");
334 }
335 
346 template <typename T>
347 inline void
348 get_indices(const dof_id_type /*nx*/,
349  const dof_id_type /*ny*/,
350  const dof_id_type /*elem_id*/,
351  dof_id_type & /*i*/,
352  dof_id_type & /*j*/,
353  dof_id_type & /*k*/)
354 {
355  mooseError("get_indices not implemented for this element type in DistributedGeneratedMesh");
356 }
357 
366 template <typename T>
367 inline void
368 get_ghost_neighbors(const dof_id_type /*nx*/,
369  const dof_id_type /*ny*/,
370  const dof_id_type /*nz*/,
371  const MeshBase & /*mesh*/,
372  std::set<dof_id_type> & /*ghost_elems*/)
373 {
374  mooseError(
375  "get_ghost_neighbors not implemented for this element type in DistributedGeneratedMesh");
376 }
377 
383 template <typename T>
384 void
385 set_boundary_names(BoundaryInfo & /*boundary_info*/)
386 {
387  mooseError(
388  "set_boundary_names not implemented for this element type in DistributedGeneratedMesh");
389 }
390 
395 template <typename T>
396 void
397 scale_nodal_positions(dof_id_type /*nx*/,
398  dof_id_type /*ny*/,
399  dof_id_type /*nz*/,
400  Real /*xmin*/,
401  Real /*xmax*/,
402  Real /*ymin*/,
403  Real /*ymax*/,
404  Real /*zmin*/,
405  Real /*zmax*/,
406  MeshBase & /*mesh*/)
407 {
408  mooseError(
409  "scale_nodal_positions not implemented for this element type in DistributedGeneratedMesh");
410 }
411 
412 template <>
413 inline dof_id_type
414 num_neighbors<Edge2>(const dof_id_type nx,
415  const dof_id_type /*ny*/,
416  const dof_id_type /*nz*/,
417  const dof_id_type i,
418  const dof_id_type /*j*/,
419  const dof_id_type /*k*/)
420 {
421  // The ends only have one neighbor
422  if (i == 0 || i == nx - 1)
423  return 1;
424 
425  return 2;
426 }
427 
428 template <>
429 inline void
430 get_neighbors<Edge2>(const dof_id_type nx,
431  const dof_id_type /*ny*/,
432  const dof_id_type /*nz*/,
433  const dof_id_type i,
434  const dof_id_type /*j*/,
435  const dof_id_type /*k*/,
436  std::vector<dof_id_type> & neighbors)
437 
438 {
439  neighbors[0] = i - 1;
440  neighbors[1] = i + 1;
441 
442  // First element doesn't have a left neighbor
443  if (i == 0)
444  neighbors[0] = Elem::invalid_id;
445 
446  // Last element doesn't have a right neighbor
447  if (i == nx - 1)
448  neighbors[1] = Elem::invalid_id;
449 }
450 
451 template <>
452 inline void
453 get_indices<Edge2>(const dof_id_type /*nx*/,
454  const dof_id_type /*ny*/,
455  const dof_id_type elem_id,
456  dof_id_type & i,
457  dof_id_type & /*j*/,
458  dof_id_type & /*k*/)
459 {
460  i = elem_id;
461 }
462 
463 template <>
464 inline void
465 get_ghost_neighbors<Edge2>(const dof_id_type nx,
466  const dof_id_type /*ny*/,
467  const dof_id_type /*nz*/,
468  const MeshBase & mesh,
469  std::set<dof_id_type> & ghost_elems)
470 {
471  auto & boundary_info = mesh.get_boundary_info();
472 
473  std::vector<dof_id_type> neighbors(2);
474 
475  for (auto elem_ptr : mesh.element_ptr_range())
476  {
477  for (unsigned int s = 0; s < elem_ptr->n_sides(); s++)
478  {
479  // No current neighbor
480  if (!elem_ptr->neighbor_ptr(s))
481  {
482  // Not on a boundary
483  if (!boundary_info.n_boundary_ids(elem_ptr, s))
484  {
485  get_neighbors<Edge2>(nx, 0, 0, elem_ptr->id(), 0, 0, neighbors);
486 
487  ghost_elems.insert(neighbors[s]);
488  }
489  }
490  }
491  }
492 }
493 
494 template <>
495 inline dof_id_type
496 elem_id<Edge2>(const dof_id_type /*nx*/,
497  const dof_id_type /*ny*/,
498  const dof_id_type i,
499  const dof_id_type /*j*/,
500  const dof_id_type /*k*/)
501 {
502  return i;
503 }
504 
505 template <>
506 void
507 add_element<Edge2>(const dof_id_type nx,
508  const dof_id_type /*ny*/,
509  const dof_id_type /*nz*/,
510  const dof_id_type /*i*/,
511  const dof_id_type /*j*/,
512  const dof_id_type /*k*/,
513  const dof_id_type elem_id,
514  const processor_id_type pid,
515  const ElemType /*type*/,
516  MeshBase & mesh,
517  bool verbose)
518 {
519  BoundaryInfo & boundary_info = mesh.get_boundary_info();
520 
521  if (verbose)
522  Moose::out << "Adding elem: " << elem_id << " pid: " << pid << std::endl;
523 
524  auto node_offset = elem_id;
525 
526  auto node0_ptr = mesh.add_point(Point(static_cast<Real>(node_offset) / nx, 0, 0), node_offset);
527  node0_ptr->set_unique_id() = node_offset;
528 
529  auto node1_ptr =
530  mesh.add_point(Point(static_cast<Real>(node_offset + 1) / nx, 0, 0), node_offset + 1);
531  node1_ptr->set_unique_id() = node_offset + 1;
532 
533  if (verbose)
534  Moose::out << "Adding elem: " << elem_id << std::endl;
535 
536  Elem * elem = new Edge2;
537  elem->set_id(elem_id);
538  elem->processor_id() = pid;
539  elem->set_unique_id() = elem_id;
540  elem = mesh.add_elem(elem);
541  elem->set_node(0) = node0_ptr;
542  elem->set_node(1) = node1_ptr;
543 
544  if (elem_id == 0)
545  boundary_info.add_side(elem, 0, 0);
546 
547  if (elem_id == nx - 1)
548  boundary_info.add_side(elem, 1, 1);
549 }
550 
551 template <>
552 void
553 set_boundary_names<Edge2>(BoundaryInfo & boundary_info)
554 {
555  boundary_info.sideset_name(0) = "left";
556  boundary_info.sideset_name(1) = "right";
557 }
558 
559 template <>
560 void
561 scale_nodal_positions<Edge2>(dof_id_type /*nx*/,
562  dof_id_type /*ny*/,
563  dof_id_type /*nz*/,
564  Real xmin,
565  Real xmax,
566  Real /*ymin*/,
567  Real /*ymax*/,
568  Real /*zmin*/,
569  Real /*zmax*/,
570  MeshBase & mesh)
571 {
572  for (auto & node_ptr : mesh.node_ptr_range())
573  (*node_ptr)(0) = (*node_ptr)(0) * (xmax - xmin) + xmin;
574 }
575 
576 template <>
577 inline dof_id_type
578 num_neighbors<Quad4>(const dof_id_type nx,
579  const dof_id_type ny,
580  const dof_id_type /*nz*/,
581  const dof_id_type i,
582  const dof_id_type j,
583  const dof_id_type /*k*/)
584 {
585  dof_id_type n = 4;
586 
587  if (i == 0)
588  n--;
589 
590  if (i == nx - 1)
591  n--;
592 
593  if (j == 0)
594  n--;
595 
596  if (j == ny - 1)
597  n--;
598 
599  return n;
600 }
601 
602 template <>
603 inline dof_id_type
604 elem_id<Quad4>(const dof_id_type nx,
605  const dof_id_type /*nx*/,
606  const dof_id_type i,
607  const dof_id_type j,
608  const dof_id_type /*k*/)
609 {
610  return (j * nx) + i;
611 }
612 
613 template <>
614 inline void
615 get_neighbors<Quad4>(const dof_id_type nx,
616  const dof_id_type ny,
617  const dof_id_type /*nz*/,
618  const dof_id_type i,
619  const dof_id_type j,
620  const dof_id_type /*k*/,
621  std::vector<dof_id_type> & neighbors)
622 {
623  std::fill(neighbors.begin(), neighbors.end(), Elem::invalid_id);
624 
625  // Bottom
626  if (j != 0)
627  neighbors[0] = elem_id<Quad4>(nx, 0, i, j - 1, 0);
628 
629  // Right
630  if (i != nx - 1)
631  neighbors[1] = elem_id<Quad4>(nx, 0, i + 1, j, 0);
632 
633  // Top
634  if (j != ny - 1)
635  neighbors[2] = elem_id<Quad4>(nx, 0, i, j + 1, 0);
636 
637  // Left
638  if (i != 0)
639  neighbors[3] = elem_id<Quad4>(nx, 0, i - 1, j, 0);
640 }
641 
642 template <>
643 inline void
644 get_indices<Quad4>(const dof_id_type nx,
645  const dof_id_type /*ny*/,
646  const dof_id_type elem_id,
647  dof_id_type & i,
648  dof_id_type & j,
649  dof_id_type & /*k*/)
650 {
651  i = elem_id % nx;
652  j = (elem_id - i) / nx;
653 }
654 
655 template <>
656 inline void
657 get_ghost_neighbors<Quad4>(const dof_id_type nx,
658  const dof_id_type ny,
659  const dof_id_type /*nz*/,
660  const MeshBase & mesh,
661  std::set<dof_id_type> & ghost_elems)
662 {
663  auto & boundary_info = mesh.get_boundary_info();
664 
665  dof_id_type i, j, k;
666 
667  std::vector<dof_id_type> neighbors(4);
668 
669  for (auto elem_ptr : mesh.element_ptr_range())
670  {
671  for (unsigned int s = 0; s < elem_ptr->n_sides(); s++)
672  {
673  // No current neighbor
674  if (!elem_ptr->neighbor_ptr(s))
675  {
676  // Not on a boundary
677  if (!boundary_info.n_boundary_ids(elem_ptr, s))
678  {
679  auto elem_id = elem_ptr->id();
680 
681  get_indices<Quad4>(nx, 0, elem_id, i, j, k);
682 
683  get_neighbors<Quad4>(nx, ny, 0, i, j, 0, neighbors);
684 
685  ghost_elems.insert(neighbors[s]);
686  }
687  }
688  }
689  }
690 }
691 
692 template <>
693 inline dof_id_type
694 node_id<Quad4>(const ElemType /*type*/,
695  const dof_id_type nx,
696  const dof_id_type /*ny*/,
697  const dof_id_type i,
698  const dof_id_type j,
699  const dof_id_type /*k*/)
700 
701 {
702  return i + j * (nx + 1);
703 }
704 
705 template <>
706 void
707 add_element<Quad4>(const dof_id_type nx,
708  const dof_id_type ny,
709  const dof_id_type /*nz*/,
710  const dof_id_type i,
711  const dof_id_type j,
712  const dof_id_type /*k*/,
713  const dof_id_type elem_id,
714  const processor_id_type pid,
715  const ElemType type,
716  MeshBase & mesh,
717  bool verbose)
718 {
719  BoundaryInfo & boundary_info = mesh.get_boundary_info();
720 
721  if (verbose)
722  Moose::out << "Adding elem: " << elem_id << " pid: " << pid << std::endl;
723 
724  // Bottom Left
725  auto node0_ptr = mesh.add_point(Point(static_cast<Real>(i) / nx, static_cast<Real>(j) / ny, 0),
726  node_id<Quad4>(type, nx, 0, i, j, 0));
727  node0_ptr->set_unique_id() = node_id<Quad4>(type, nx, 0, i, j, 0);
728 
729  // Bottom Right
730  auto node1_ptr =
731  mesh.add_point(Point(static_cast<Real>(i + 1) / nx, static_cast<Real>(j) / ny, 0),
732  node_id<Quad4>(type, nx, 0, i + 1, j, 0));
733  node1_ptr->set_unique_id() = node_id<Quad4>(type, nx, 0, i + 1, j, 0);
734 
735  // Top Right
736  auto node2_ptr =
737  mesh.add_point(Point(static_cast<Real>(i + 1) / nx, static_cast<Real>(j + 1) / ny, 0),
738  node_id<Quad4>(type, nx, 0, i + 1, j + 1, 0));
739  node2_ptr->set_unique_id() = node_id<Quad4>(type, nx, 0, i + 1, j + 1, 0);
740 
741  // Top Left
742  auto node3_ptr =
743  mesh.add_point(Point(static_cast<Real>(i) / nx, static_cast<Real>(j + 1) / ny, 0),
744  node_id<Quad4>(type, nx, 0, i, j + 1, 0));
745  node3_ptr->set_unique_id() = node_id<Quad4>(type, nx, 0, i, j + 1, 0);
746 
747  Elem * elem = new Quad4;
748  elem->set_id(elem_id);
749  elem->processor_id() = pid;
750  elem->set_unique_id() = elem_id;
751  elem = mesh.add_elem(elem);
752  elem->set_node(0) = node0_ptr;
753  elem->set_node(1) = node1_ptr;
754  elem->set_node(2) = node2_ptr;
755  elem->set_node(3) = node3_ptr;
756 
757  // Bottom
758  if (j == 0)
759  boundary_info.add_side(elem, 0, 0);
760 
761  // Right
762  if (i == nx - 1)
763  boundary_info.add_side(elem, 1, 1);
764 
765  // Top
766  if (j == ny - 1)
767  boundary_info.add_side(elem, 2, 2);
768 
769  // Left
770  if (i == 0)
771  boundary_info.add_side(elem, 3, 3);
772 }
773 
774 template <>
775 void
776 set_boundary_names<Quad4>(BoundaryInfo & boundary_info)
777 {
778  boundary_info.sideset_name(0) = "bottom";
779  boundary_info.sideset_name(1) = "right";
780  boundary_info.sideset_name(2) = "top";
781  boundary_info.sideset_name(3) = "left";
782 }
783 
784 template <>
785 void
786 scale_nodal_positions<Quad4>(dof_id_type /*nx*/,
787  dof_id_type /*ny*/,
788  dof_id_type /*nz*/,
789  Real xmin,
790  Real xmax,
791  Real ymin,
792  Real ymax,
793  Real /*zmin*/,
794  Real /*zmax*/,
795  MeshBase & mesh)
796 {
797  for (auto & node_ptr : mesh.node_ptr_range())
798  {
799  (*node_ptr)(0) = (*node_ptr)(0) * (xmax - xmin) + xmin;
800  (*node_ptr)(1) = (*node_ptr)(1) * (ymax - ymin) + ymin;
801  }
802 }
803 
804 template <>
805 inline dof_id_type
806 elem_id<Hex8>(const dof_id_type nx,
807  const dof_id_type ny,
808  const dof_id_type i,
809  const dof_id_type j,
810  const dof_id_type k)
811 {
812  return i + (j * nx) + (k * nx * ny);
813 }
814 
815 template <>
816 inline dof_id_type
817 num_neighbors<Hex8>(const dof_id_type nx,
818  const dof_id_type ny,
819  const dof_id_type nz,
820  const dof_id_type i,
821  const dof_id_type j,
822  const dof_id_type k)
823 {
824  dof_id_type n = 6;
825 
826  if (i == 0)
827  n--;
828 
829  if (i == nx - 1)
830  n--;
831 
832  if (j == 0)
833  n--;
834 
835  if (j == ny - 1)
836  n--;
837 
838  if (k == 0)
839  n--;
840 
841  if (k == nz - 1)
842  n--;
843 
844  return n;
845 }
846 
847 template <>
848 inline void
849 get_neighbors<Hex8>(const dof_id_type nx,
850  const dof_id_type ny,
851  const dof_id_type nz,
852  const dof_id_type i,
853  const dof_id_type j,
854  const dof_id_type k,
855  std::vector<dof_id_type> & neighbors)
856 {
857  std::fill(neighbors.begin(), neighbors.end(), Elem::invalid_id);
858 
859  // Back
860  if (k != 0)
861  neighbors[0] = elem_id<Hex8>(nx, ny, i, j, k - 1);
862 
863  // Bottom
864  if (j != 0)
865  neighbors[1] = elem_id<Hex8>(nx, ny, i, j - 1, k);
866 
867  // Right
868  if (i != nx - 1)
869  neighbors[2] = elem_id<Hex8>(nx, ny, i + 1, j, k);
870 
871  // Top
872  if (j != ny - 1)
873  neighbors[3] = elem_id<Hex8>(nx, ny, i, j + 1, k);
874 
875  // Left
876  if (i != 0)
877  neighbors[4] = elem_id<Hex8>(nx, ny, i - 1, j, k);
878 
879  // Front
880  if (k != nz - 1)
881  neighbors[5] = elem_id<Hex8>(nx, ny, i, j, k + 1);
882 }
883 
884 template <>
885 inline dof_id_type
886 node_id<Hex8>(const ElemType /*type*/,
887  const dof_id_type nx,
888  const dof_id_type ny,
889  const dof_id_type i,
890  const dof_id_type j,
891  const dof_id_type k)
892 {
893  return i + (nx + 1) * (j + k * (ny + 1));
894 }
895 
896 template <>
897 Node *
898 add_point<Hex8>(const dof_id_type nx,
899  const dof_id_type ny,
900  const dof_id_type nz,
901  const dof_id_type i,
902  const dof_id_type j,
903  const dof_id_type k,
904  const ElemType type,
905  MeshBase & mesh)
906 {
907  auto id = node_id<Hex8>(type, nx, ny, i, j, k);
908  auto node_ptr = mesh.add_point(
909  Point(static_cast<Real>(i) / nx, static_cast<Real>(j) / ny, static_cast<Real>(k) / nz), id);
910  node_ptr->set_unique_id() = id;
911 
912  return node_ptr;
913 }
914 
915 template <>
916 void
917 add_element<Hex8>(const dof_id_type nx,
918  const dof_id_type ny,
919  const dof_id_type nz,
920  const dof_id_type i,
921  const dof_id_type j,
922  const dof_id_type k,
923  const dof_id_type elem_id,
924  const processor_id_type pid,
925  const ElemType type,
926  MeshBase & mesh,
927  bool verbose)
928 {
929  BoundaryInfo & boundary_info = mesh.get_boundary_info();
930 
931  if (verbose)
932  {
933  Moose::out << "Adding elem: " << elem_id << " pid: " << pid << std::endl;
934  Moose::out << "Type: " << type << " " << HEX8 << std::endl;
935  }
936 
937  // This ordering was picked to match the ordering in mesh_generation.C
938  auto node0_ptr = add_point<Hex8>(nx, ny, nz, i, j, k, type, mesh);
939  auto node1_ptr = add_point<Hex8>(nx, ny, nz, i + 1, j, k, type, mesh);
940  auto node2_ptr = add_point<Hex8>(nx, ny, nz, i + 1, j + 1, k, type, mesh);
941  auto node3_ptr = add_point<Hex8>(nx, ny, nz, i, j + 1, k, type, mesh);
942  auto node4_ptr = add_point<Hex8>(nx, ny, nz, i, j, k + 1, type, mesh);
943  auto node5_ptr = add_point<Hex8>(nx, ny, nz, i + 1, j, k + 1, type, mesh);
944  auto node6_ptr = add_point<Hex8>(nx, ny, nz, i + 1, j + 1, k + 1, type, mesh);
945  auto node7_ptr = add_point<Hex8>(nx, ny, nz, i, j + 1, k + 1, type, mesh);
946 
947  Elem * elem = new Hex8;
948  elem->set_id(elem_id);
949  elem->processor_id() = pid;
950  elem->set_unique_id() = elem_id;
951  elem = mesh.add_elem(elem);
952  elem->set_node(0) = node0_ptr;
953  elem->set_node(1) = node1_ptr;
954  elem->set_node(2) = node2_ptr;
955  elem->set_node(3) = node3_ptr;
956  elem->set_node(4) = node4_ptr;
957  elem->set_node(5) = node5_ptr;
958  elem->set_node(6) = node6_ptr;
959  elem->set_node(7) = node7_ptr;
960 
961  if (k == 0)
962  boundary_info.add_side(elem, 0, 0);
963 
964  if (k == (nz - 1))
965  boundary_info.add_side(elem, 5, 5);
966 
967  if (j == 0)
968  boundary_info.add_side(elem, 1, 1);
969 
970  if (j == (ny - 1))
971  boundary_info.add_side(elem, 3, 3);
972 
973  if (i == 0)
974  boundary_info.add_side(elem, 4, 4);
975 
976  if (i == (nx - 1))
977  boundary_info.add_side(elem, 2, 2);
978 }
979 
980 template <>
981 inline void
982 get_indices<Hex8>(const dof_id_type nx,
983  const dof_id_type ny,
984  const dof_id_type elem_id,
985  dof_id_type & i,
986  dof_id_type & j,
987  dof_id_type & k)
988 {
989  i = elem_id % nx;
990  j = (((elem_id - i) / nx) % ny);
991  k = ((elem_id - i) - (j * nx)) / (nx * ny);
992 }
993 
994 template <>
995 inline void
996 get_ghost_neighbors<Hex8>(const dof_id_type nx,
997  const dof_id_type ny,
998  const dof_id_type nz,
999  const MeshBase & mesh,
1000  std::set<dof_id_type> & ghost_elems)
1001 {
1002  auto & boundary_info = mesh.get_boundary_info();
1003 
1004  dof_id_type i, j, k;
1005 
1006  std::vector<dof_id_type> neighbors(6);
1007 
1008  for (auto elem_ptr : mesh.element_ptr_range())
1009  {
1010  for (unsigned int s = 0; s < elem_ptr->n_sides(); s++)
1011  {
1012  // No current neighbor
1013  if (!elem_ptr->neighbor_ptr(s))
1014  {
1015  // Not on a boundary
1016  if (!boundary_info.n_boundary_ids(elem_ptr, s))
1017  {
1018  auto elem_id = elem_ptr->id();
1019 
1020  get_indices<Hex8>(nx, ny, elem_id, i, j, k);
1021 
1022  get_neighbors<Hex8>(nx, ny, nz, i, j, k, neighbors);
1023 
1024  ghost_elems.insert(neighbors[s]);
1025  }
1026  }
1027  }
1028  }
1029 }
1030 
1031 template <>
1032 void
1033 set_boundary_names<Hex8>(BoundaryInfo & boundary_info)
1034 {
1035  boundary_info.sideset_name(0) = "back";
1036  boundary_info.sideset_name(1) = "bottom";
1037  boundary_info.sideset_name(2) = "right";
1038  boundary_info.sideset_name(3) = "top";
1039  boundary_info.sideset_name(4) = "left";
1040  boundary_info.sideset_name(5) = "front";
1041 }
1042 
1043 template <>
1044 void
1045 scale_nodal_positions<Hex8>(dof_id_type /*nx*/,
1046  dof_id_type /*ny*/,
1047  dof_id_type /*nz*/,
1048  Real xmin,
1049  Real xmax,
1050  Real ymin,
1051  Real ymax,
1052  Real zmin,
1053  Real zmax,
1054  MeshBase & mesh)
1055 {
1056  for (auto & node_ptr : mesh.node_ptr_range())
1057  {
1058  (*node_ptr)(0) = (*node_ptr)(0) * (xmax - xmin) + xmin;
1059  (*node_ptr)(1) = (*node_ptr)(1) * (ymax - ymin) + ymin;
1060  (*node_ptr)(2) = (*node_ptr)(2) * (zmax - zmin) + zmin;
1061  }
1062 }
1063 }
1064 
1065 template <typename T>
1066 void
1067 build_cube(UnstructuredMesh & mesh,
1068  const unsigned int nx,
1069  unsigned int ny,
1070  unsigned int nz,
1071  const Real xmin,
1072  const Real xmax,
1073  const Real ymin,
1074  const Real ymax,
1075  const Real zmin,
1076  const Real zmax,
1077  const ElemType type,
1078  bool verbose)
1079 {
1080  if (verbose)
1081  Moose::out << "nx: " << nx << "\n ny: " << ny << "\n nz: " << nz << std::endl;
1082 
1087 
1088  dof_id_type num_elems = nx * ny * nz;
1089  const auto n_pieces = mesh.comm().size();
1090  const auto pid = mesh.comm().rank();
1091 
1092  std::unique_ptr<Elem> canonical_elem = libmesh_make_unique<T>();
1093 
1094  // Will get used to find the neighbors of an element
1095  std::vector<dof_id_type> neighbors(canonical_elem->n_neighbors());
1096 
1097  // Data structure that Metis will fill up on processor 0 and broadcast.
1098  std::vector<Metis::idx_t> part(num_elems);
1099 
1100  if (mesh.processor_id() == 0)
1101  {
1102  // Data structures and parameters needed only on processor 0 by Metis.
1103  // std::vector<Metis::idx_t> options(5);
1104  // Weight by the number of nodes
1105  std::vector<Metis::idx_t> vwgt(num_elems, canonical_elem->n_nodes());
1106 
1107  Metis::idx_t n = num_elems, // number of "nodes" (elements) in the graph
1108  nparts = n_pieces, // number of subdomains to create
1109  edgecut = 0; // the numbers of edges cut by the resulting partition
1110 
1111  METIS_CSR_Graph<Metis::idx_t> csr_graph;
1112 
1113  csr_graph.offsets.resize(num_elems + 1, 0);
1114 
1115  for (dof_id_type k = 0; k < nz; k++)
1116  {
1117  for (dof_id_type j = 0; j < ny; j++)
1118  {
1119  for (dof_id_type i = 0; i < nx; i++)
1120  {
1121  auto n_neighbors = num_neighbors<T>(nx, ny, nz, i, j, k);
1122 
1123  auto e_id = elem_id<T>(nx, ny, i, j, k);
1124 
1125  if (verbose)
1126  Moose::out << e_id << " num_neighbors: " << n_neighbors << std::endl;
1127 
1128  csr_graph.prep_n_nonzeros(e_id, n_neighbors);
1129  }
1130  }
1131  }
1132 
1133  csr_graph.prepare_for_use();
1134 
1135  if (verbose)
1136  for (auto offset : csr_graph.offsets)
1137  Moose::out << "offset: " << offset << std::endl;
1138 
1139  for (dof_id_type k = 0; k < nz; k++)
1140  {
1141  for (dof_id_type j = 0; j < ny; j++)
1142  {
1143  for (dof_id_type i = 0; i < nx; i++)
1144  {
1145  auto e_id = elem_id<T>(nx, ny, i, j, k);
1146 
1147  dof_id_type connection = 0;
1148 
1149  get_neighbors<T>(nx, ny, nz, i, j, k, neighbors);
1150 
1151  for (auto neighbor : neighbors)
1152  {
1153  if (neighbor != Elem::invalid_id)
1154  {
1155  if (verbose)
1156  Moose::out << e_id << ": " << connection << " = " << neighbor << std::endl;
1157 
1158  csr_graph(e_id, connection++) = neighbor;
1159  }
1160  }
1161  }
1162  }
1163  }
1164 
1165  if (n_pieces == 1)
1166  {
1167  // Just assign them all to proc 0
1168  for (auto & elem_pid : part)
1169  elem_pid = 0;
1170  }
1171  else
1172  {
1173  Metis::idx_t ncon = 1;
1174 
1175  // Use recursive if the number of partitions is less than or equal to 8
1176  if (n_pieces <= 8)
1177  Metis::METIS_PartGraphRecursive(&n,
1178  &ncon,
1179  &csr_graph.offsets[0],
1180  &csr_graph.vals[0],
1181  &vwgt[0],
1182  libmesh_nullptr,
1183  libmesh_nullptr,
1184  &nparts,
1185  libmesh_nullptr,
1186  libmesh_nullptr,
1187  libmesh_nullptr,
1188  &edgecut,
1189  &part[0]);
1190 
1191  // Otherwise use kway
1192  else
1193  Metis::METIS_PartGraphKway(&n,
1194  &ncon,
1195  &csr_graph.offsets[0],
1196  &csr_graph.vals[0],
1197  &vwgt[0],
1198  libmesh_nullptr,
1199  libmesh_nullptr,
1200  &nparts,
1201  libmesh_nullptr,
1202  libmesh_nullptr,
1203  libmesh_nullptr,
1204  &edgecut,
1205  &part[0]);
1206  }
1207  } // end processor 0 part
1208 
1209  // Broadcast the resulting partition
1210  mesh.comm().broadcast(part);
1211 
1212  if (verbose)
1213  for (auto proc_id : part)
1214  Moose::out << "Part: " << proc_id << std::endl;
1215 
1216  BoundaryInfo & boundary_info = mesh.get_boundary_info();
1217 
1218  // Add elements this processor owns
1219  for (dof_id_type k = 0; k < nz; k++)
1220  {
1221  for (dof_id_type j = 0; j < ny; j++)
1222  {
1223  for (dof_id_type i = 0; i < nx; i++)
1224  {
1225  auto e_id = elem_id<Hex8>(nx, ny, i, j, k);
1226 
1227  if (static_cast<processor_id_type>(part[e_id]) == pid)
1228  add_element<T>(nx, ny, nz, i, j, k, e_id, pid, type, mesh, verbose);
1229  }
1230  }
1231  }
1232 
1233  if (verbose)
1234  for (auto & elem_ptr : mesh.element_ptr_range())
1235  for (unsigned int s = 0; s < elem_ptr->n_sides(); s++)
1236  Moose::out << "Elem neighbor: " << elem_ptr->neighbor_ptr(s) << " is remote "
1237  << (elem_ptr->neighbor_ptr(s) == remote_elem) << std::endl;
1238 
1239  // Need to link up the local elements before we can know what's missing
1240  mesh.find_neighbors();
1241 
1242  if (verbose)
1243  Moose::out << "After first find_neighbors" << std::endl;
1244 
1245  if (verbose)
1246  for (auto & elem_ptr : mesh.element_ptr_range())
1247  for (unsigned int s = 0; s < elem_ptr->n_sides(); s++)
1248  Moose::out << "Elem neighbor: " << elem_ptr->neighbor_ptr(s) << " is remote "
1249  << (elem_ptr->neighbor_ptr(s) == remote_elem) << std::endl;
1250 
1251  // Get the ghosts (missing face neighbors)
1252  std::set<dof_id_type> ghost_elems;
1253  get_ghost_neighbors<T>(nx, ny, nz, mesh, ghost_elems);
1254 
1255  // Add the ghosts to the mesh
1256  for (auto & ghost_id : ghost_elems)
1257  {
1258  dof_id_type i = 0;
1259  dof_id_type j = 0;
1260  dof_id_type k = 0;
1261 
1262  get_indices<T>(nx, ny, ghost_id, i, j, k);
1263 
1264  add_element<T>(nx, ny, nz, i, j, k, ghost_id, part[ghost_id], type, mesh, verbose);
1265  }
1266 
1267  if (verbose)
1268  Moose::out << "After adding ghosts" << std::endl;
1269 
1270  if (verbose)
1271  for (auto & elem_ptr : mesh.element_ptr_range())
1272  for (unsigned int s = 0; s < elem_ptr->n_sides(); s++)
1273  Moose::out << "Elem neighbor: " << elem_ptr->neighbor_ptr(s) << " is remote "
1274  << (elem_ptr->neighbor_ptr(s) == remote_elem) << std::endl;
1275 
1276  mesh.find_neighbors(true);
1277 
1278  if (verbose)
1279  Moose::out << "After second find neighbors " << std::endl;
1280 
1281  if (verbose)
1282  for (auto & elem_ptr : mesh.element_ptr_range())
1283  for (unsigned int s = 0; s < elem_ptr->n_sides(); s++)
1284  Moose::out << "Elem neighbor: " << elem_ptr->neighbor_ptr(s) << " is remote "
1285  << (elem_ptr->neighbor_ptr(s) == remote_elem) << std::endl;
1286 
1287  // Set RemoteElem neighbors
1288  for (auto & elem_ptr : mesh.element_ptr_range())
1289  for (unsigned int s = 0; s < elem_ptr->n_sides(); s++)
1290  if (!elem_ptr->neighbor_ptr(s) && !boundary_info.n_boundary_ids(elem_ptr, s))
1291  elem_ptr->set_neighbor(s, const_cast<RemoteElem *>(remote_elem));
1292 
1293  if (verbose)
1294  Moose::out << "After adding remote elements" << std::endl;
1295 
1296  if (verbose)
1297  for (auto & elem_ptr : mesh.element_ptr_range())
1298  for (unsigned int s = 0; s < elem_ptr->n_sides(); s++)
1299  Moose::out << "Elem neighbor: " << elem_ptr->neighbor_ptr(s) << " is remote "
1300  << (elem_ptr->neighbor_ptr(s) == remote_elem) << std::endl;
1301 
1302  set_boundary_names<T>(boundary_info);
1303 
1304  Partitioner::set_node_processor_ids(mesh);
1305 
1306  if (verbose)
1307  Moose::out << "mesh dim: " << mesh.mesh_dimension() << std::endl;
1308 
1309  if (verbose)
1310  for (auto & node_ptr : mesh.node_ptr_range())
1311  Moose::out << node_ptr->id() << ":" << node_ptr->processor_id() << std::endl;
1312 
1313  // Already partitioned!
1314  mesh.skip_partitioning(true);
1315 
1316  if (verbose)
1317  for (auto & elem_ptr : mesh.element_ptr_range())
1318  Moose::out << "Elem: " << elem_ptr->id() << " pid: " << elem_ptr->processor_id()
1319  << " uid: " << elem_ptr->unique_id() << std::endl;
1320 
1321  if (verbose)
1322  Moose::out << "Getting ready to prepare for use" << std::endl;
1323 
1324  // No need to renumber or find neighbors - done did it.
1325  // Avoid deprecation message/error by _also_ setting
1326  // allow_renumbering(false). This is a bit silly, but we want to
1327  // catch cases where people are purely using the old "skip"
1328  // interface and not the new flag setting one.
1329  mesh.allow_renumbering(false);
1330  mesh.prepare_for_use(/*skip_renumber (ignored!) = */ false,
1331  /*skip_find_neighbors = */ true);
1332 
1333  if (verbose)
1334  for (auto & elem_ptr : mesh.element_ptr_range())
1335  Moose::out << "Elem: " << elem_ptr->id() << " pid: " << elem_ptr->processor_id() << std::endl;
1336 
1337  if (verbose)
1338  for (auto & node_ptr : mesh.node_ptr_range())
1339  Moose::out << node_ptr->id() << ":" << node_ptr->processor_id() << std::endl;
1340 
1341  if (verbose)
1342  Moose::out << "mesh dim: " << mesh.mesh_dimension() << std::endl;
1343 
1344  // Scale the nodal positions
1345  scale_nodal_positions<T>(nx, ny, nz, xmin, xmax, ymin, ymax, zmin, zmax, mesh);
1346 
1347  if (verbose)
1348  mesh.print_info();
1349 }
1350 
1351 void
1353 {
1354  // Reference to the libmesh mesh
1355  MeshBase & mesh = getMesh();
1356 
1357  MooseEnum elem_type_enum = getParam<MooseEnum>("elem_type");
1358 
1359  if (!isParamValid("elem_type"))
1360  {
1361  // Switching on MooseEnum
1362  switch (_dim)
1363  {
1364  case 1:
1365  elem_type_enum = "EDGE2";
1366  break;
1367  case 2:
1368  elem_type_enum = "QUAD4";
1369  break;
1370  case 3:
1371  elem_type_enum = "HEX8";
1372  break;
1373  }
1374  }
1375 
1376  _elem_type = Utility::string_to_enum<ElemType>(elem_type_enum);
1377 
1378  mesh.set_mesh_dimension(_dim);
1379  mesh.set_spatial_dimension(_dim);
1380 
1381  // Switching on MooseEnum
1382  switch (_dim)
1383  {
1384  // The build_XYZ mesh generation functions take an
1385  // UnstructuredMesh& as the first argument, hence the dynamic_cast.
1386  case 1:
1387  build_cube<Edge2>(dynamic_cast<UnstructuredMesh &>(getMesh()),
1388  _nx,
1389  1,
1390  1,
1391  _xmin,
1392  _xmax,
1393  0,
1394  0,
1395  0,
1396  0,
1397  _elem_type,
1398  _verbose);
1399  break;
1400  case 2:
1401  build_cube<Quad4>(dynamic_cast<UnstructuredMesh &>(getMesh()),
1402  _nx,
1403  _ny,
1404  1,
1405  _xmin,
1406  _xmax,
1407  _ymin,
1408  _ymax,
1409  0,
1410  0,
1411  _elem_type,
1412  _verbose);
1413  break;
1414  case 3:
1415  build_cube<Hex8>(dynamic_cast<UnstructuredMesh &>(getMesh()),
1416  _nx,
1417  _ny,
1418  _nz,
1419  _xmin,
1420  _xmax,
1421  _ymin,
1422  _ymax,
1423  _zmin,
1424  _zmax,
1425  _elem_type,
1426  _verbose);
1427  break;
1428  default:
1429  mooseError(getParam<MooseEnum>("elem_type"),
1430  " is not a currently supported element type for DistributedGeneratedMesh");
1431  }
1432 
1433  // Apply the bias if any exists
1434  if (_bias_x != 1.0 || _bias_y != 1.0 || _bias_z != 1.0)
1435  {
1436  // Biases
1437  Real bias[3] = {_bias_x, _bias_y, _bias_z};
1438 
1439  // "width" of the mesh in each direction
1440  Real width[3] = {_xmax - _xmin, _ymax - _ymin, _zmax - _zmin};
1441 
1442  // Min mesh extent in each direction.
1443  Real mins[3] = {_xmin, _ymin, _zmin};
1444 
1445  // Number of elements in each direction.
1446  dof_id_type nelem[3] = {_nx, _ny, _nz};
1447 
1448  // We will need the biases raised to integer powers in each
1449  // direction, so let's pre-compute those...
1450  std::vector<std::vector<Real>> pows(LIBMESH_DIM);
1451  for (dof_id_type dir = 0; dir < LIBMESH_DIM; ++dir)
1452  {
1453  pows[dir].resize(nelem[dir] + 1);
1454  for (dof_id_type i = 0; i < pows[dir].size(); ++i)
1455  pows[dir][i] = std::pow(bias[dir], static_cast<int>(i));
1456  }
1457 
1458  // Loop over the nodes and move them to the desired location
1459  for (auto & node_ptr : mesh.node_ptr_range())
1460  {
1461  Node & node = *node_ptr;
1462 
1463  for (dof_id_type dir = 0; dir < LIBMESH_DIM; ++dir)
1464  {
1465  if (width[dir] != 0. && bias[dir] != 1.)
1466  {
1467  // Compute the scaled "index" of the current point. This
1468  // will either be close to a whole integer or a whole
1469  // integer+0.5 for quadratic nodes.
1470  Real float_index = (node(dir) - mins[dir]) * nelem[dir] / width[dir];
1471 
1472  Real integer_part = 0;
1473  Real fractional_part = std::modf(float_index, &integer_part);
1474 
1475  // Figure out where to move the node...
1476  if (std::abs(fractional_part) < TOLERANCE || std::abs(fractional_part - 1.0) < TOLERANCE)
1477  {
1478  // If the fractional_part ~ 0.0 or 1.0, this is a vertex node, so
1479  // we don't need an average.
1480  //
1481  // Compute the "index" we are at in the current direction. We
1482  // round to the nearest integral value to get this instead
1483  // of using "integer_part", since that could be off by a
1484  // lot (e.g. we want 3.9999 to map to 4.0 instead of 3.0).
1485  int index = round(float_index);
1486 
1487  // Move node to biased location.
1488  node(dir) =
1489  mins[dir] + width[dir] * (1. - pows[dir][index]) / (1. - pows[dir][nelem[dir]]);
1490  }
1491  else if (std::abs(fractional_part - 0.5) < TOLERANCE)
1492  {
1493  // If the fractional_part ~ 0.5, this is a midedge/face
1494  // (i.e. quadratic) node. We don't move those with the same
1495  // bias as the vertices, instead we put them midway between
1496  // their respective vertices.
1497  //
1498  // Also, since the fractional part is nearly 0.5, we know that
1499  // the integer_part will be the index of the vertex to the
1500  // left, and integer_part+1 will be the index of the
1501  // vertex to the right.
1502  node(dir) = mins[dir] +
1503  width[dir] *
1504  (1. - 0.5 * (pows[dir][integer_part] + pows[dir][integer_part + 1])) /
1505  (1. - pows[dir][nelem[dir]]);
1506  }
1507  else
1508  {
1509  // We don't yet handle anything higher order than quadratic...
1510  mooseError("Unable to bias node at node(", dir, ")=", node(dir));
1511  }
1512  }
1513  }
1514  }
1515  }
1516 }
virtual Real getMaxInDimension(unsigned int component) const
Definition: MooseMesh.C:1475
Real _bias_x
The amount by which to bias the cells in the x,y,z directions.
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)
DistributedGeneratedMesh(const InputParameters &parameters)
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:207
Mesh generated from parameters.
bool _verbose
Print output during generation.
Real round(Real x)
Definition: MathUtils.h:22
T & set(const std::string &name, bool quiet_mode=false)
Returns a writable reference to the named parameters.
virtual Real getMaxInDimension(unsigned int component) const override
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
bool _dims_may_have_changed
Boolean to indicate that dimensions may have changed.
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 _xmin
The min/max values for x,y,z component.
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
registerMooseObject("MooseApp", DistributedGeneratedMesh)
virtual void buildMesh() override
Must be overridden by child classes.
virtual Real getMinInDimension(unsigned int component) const override
Returns the min or max of the requested dimension respectively.
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
unsigned int _nx
Number of elements in x, y, z direction.
InputParameters validParams< DistributedGeneratedMesh >()
MatType type
PetscInt n
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...
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)
MooseEnum _dim
The dimension of the mesh.
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseObject.h:89
ElemType _elem_type
The type of element to build.
bool _regular_orthogonal_mesh
Boolean indicating whether this mesh was detected to be regular and orthogonal.
Definition: MooseMesh.h:1009
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...
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...