libMesh
Namespaces | Classes | Functions
libMesh::MeshTools::Generation Namespace Reference

Tools for Mesh generation. More...

Namespaces

 Private
 

Classes

class  QueryElemSubdomainIDBase
 Class for receiving the callback during extrusion generation and providing user-defined subdomains based on the old (existing) element id and the current layer. More...
 

Functions

void build_cube (UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
 Builds a \( nx \times ny \times nz \) (elements) cube. More...
 
void build_point (UnstructuredMesh &mesh, const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
 A specialized build_cube() for 0D meshes. More...
 
void build_line (UnstructuredMesh &mesh, const unsigned int nx, const Real xmin=0., const Real xmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
 A specialized build_cube() for 1D meshes. More...
 
void build_square (UnstructuredMesh &mesh, const unsigned int nx, const unsigned int ny, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
 A specialized build_cube() for 2D meshes. More...
 
void build_sphere (UnstructuredMesh &mesh, const Real rad=1, const unsigned int nr=2, const ElemType type=INVALID_ELEM, const unsigned int n_smooth=2, const bool flat=true)
 Meshes a spherical or mapped-spherical domain. More...
 
void build_extrusion (UnstructuredMesh &mesh, const MeshBase &cross_section, const unsigned int nz, RealVectorValue extrusion_vector, QueryElemSubdomainIDBase *elem_subdomain=nullptr)
 Meshes the tensor product of a 1D and a 1D-or-2D domain. More...
 
void build_delaunay_square (UnstructuredMesh &mesh, const unsigned int nx, const unsigned int ny, const Real xmin, const Real xmax, const Real ymin, const Real ymax, const ElemType type, const std::vector< TriangleInterface::Hole *> *holes=nullptr)
 Meshes a rectangular (2D) region (with or without holes) with a Delaunay triangulation. More...
 

Detailed Description

Tools for Mesh generation.

Author
Benjamin S. Kirk
Date
2004

Function Documentation

◆ build_cube()

void libMesh::MeshTools::Generation::build_cube ( UnstructuredMesh &  mesh,
const unsigned int  nx = 0,
const unsigned int  ny = 0,
const unsigned int  nz = 0,
const Real  xmin = 0.,
const Real  xmax = 1.,
const Real  ymin = 0.,
const Real  ymax = 1.,
const Real  zmin = 0.,
const Real  zmax = 1.,
const ElemType  type = INVALID_ELEM,
const bool  gauss_lobatto_grid = false 
)

Builds a \( nx \times ny \times nz \) (elements) cube.

Defaults to a unit cube (or line in 1D, square in 2D), but the dimensions can be specified through the optional arguments.

Boundary ids are set to be equal to the side indexing on a master hex

Definition at line 311 of file mesh_generation.C.

References libMesh::BoundaryInfo::add_node(), libMesh::BoundaryInfo::add_side(), libMesh::BoundaryInfo::boundary_ids(), libMesh::Elem::build(), libMesh::Elem::build_with_id(), libMesh::EDGE2, libMesh::EDGE3, libMesh::EDGE4, libMesh::Utility::enum_to_string(), libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::DofObject::id(), libMesh::MeshTools::Generation::Private::idx(), libMesh::index_range(), libMesh::INVALID_ELEM, libMesh::BoundaryInfo::invalid_id, libMesh::libmesh_assert(), mesh, libMesh::NODEELEM, libMesh::BoundaryInfo::nodeset_name(), libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM20, libMesh::PRISM21, libMesh::PRISM6, libMesh::PYRAMID13, libMesh::PYRAMID14, libMesh::PYRAMID18, libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::QUADSHELL4, libMesh::QUADSHELL8, libMesh::Real, libMesh::MeshTools::Modification::redistribute(), libMesh::BoundaryInfo::remove(), libMesh::Elem::set_node(), libMesh::BoundaryInfo::sideset_name(), libMesh::TET10, libMesh::TET14, libMesh::TET4, libMesh::TRI3, libMesh::TRI6, and libMesh::TRI7.

Referenced by add_cube_convex_hull_to_mesh(), AllSecondOrderTest::allCompleteOrderMixedFixing3D(), AllSecondOrderTest::allSecondOrderMixedFixing3D(), Biharmonic::Biharmonic(), build_line(), ExtraIntegersTest::build_mesh(), build_point(), build_square(), CopyNodesAndElementsTest::collectMeshes(), main(), PerElemTest< elem_type >::setUp(), ParsedFEMFunctionTest::setUp(), RationalMapTest< elem_type >::setUp(), setup(), ParallelGhostSyncTest::setUp(), FETestBase< order, family, elem_type, 1 >::setUp(), SystemsTest::test3DProjectVectorFE(), DistortTest::test_helper_3D(), AllTriTest::test_helper_3D(), MeshFunctionTest::test_p_level(), VolumeTest::test_true_centroid_and_volume(), MeshStitchTest::testAmbiguousRemappingStitch(), SystemsTest::testAssemblyWithDgFemContext(), SystemsTest::testBlockRestrictedVarNDofs(), MeshStitchTest::testBoundaryInfo(), SystemsTest::testBoundaryProjectCube(), MeshGenerationTest::testBuildCube(), DefaultCouplingTest::testCoupling(), PointNeighborCouplingTest::testCoupling(), SystemsTest::testDofCouplingWithVarGroups(), DofMapTest::testDofOwner(), BoundaryInfoTest::testEdgeBoundaryConditions(), MeshInputTest::testExodusWriteElementDataFromDiscontinuousNodalData(), VolumeTest::testHex20PLevelTrueCentroid(), InfFERadialTest::testInfQuants(), PointLocatorTest::testLocator(), MeshStitchTest::testMeshStitch(), MeshStitchTest::testMeshStitchElemsets(), PartitionerTest< PartitionerSubclass, MeshClass >::testPartition(), SystemsTest::testProjectCube(), SystemsTest::testProjectCubeWithMeshFunction(), SystemsTest::testProjectMatrix3D(), MeshStitchTest::testRemappingStitch(), InfFERadialTest::testSides(), InfFERadialTest::testSingleOrder(), VolumeTest::testTwistedVolume(), and WriteEdgesetData::testWriteImpl().

320 {
321  LOG_SCOPE("build_cube()", "MeshTools::Generation");
322 
323  // Declare that we are using the indexing utility routine
324  // in the "Private" part of our current namespace. If this doesn't
325  // work in GCC 2.95.3 we can either remove it or stop supporting
326  // 2.95.3 altogether.
327  // Changing this to import the whole namespace... just importing idx
328  // causes an internal compiler error for Intel Compiler 11.0 on Linux
329  // in debug mode.
330  using namespace MeshTools::Generation::Private;
331 
332  // Clear the mesh and start from scratch
333  mesh.clear();
334 
335  BoundaryInfo & boundary_info = mesh.get_boundary_info();
336 
337  if (nz != 0)
338  {
339  mesh.set_mesh_dimension(3);
340  mesh.set_spatial_dimension(3);
341  }
342  else if (ny != 0)
343  {
344  mesh.set_mesh_dimension(2);
345  mesh.set_spatial_dimension(2);
346  }
347  else if (nx != 0)
348  {
349  mesh.set_mesh_dimension(1);
350  mesh.set_spatial_dimension(1);
351  }
352  else
353  {
354  // Will we get here?
355  mesh.set_mesh_dimension(0);
356  mesh.set_spatial_dimension(0);
357  }
358 
359  switch (mesh.mesh_dimension())
360  {
361  //---------------------------------------------------------------------
362  // Build a 0D point
363  case 0:
364  {
365  libmesh_assert_equal_to (nx, 0);
366  libmesh_assert_equal_to (ny, 0);
367  libmesh_assert_equal_to (nz, 0);
368 
369  libmesh_assert (type == INVALID_ELEM || type == NODEELEM);
370 
371  // Build one nodal element for the mesh
372  mesh.add_point (Point(0, 0, 0), 0);
373  Elem * elem = mesh.add_elem(Elem::build(NODEELEM));
374  elem->set_node(0) = mesh.node_ptr(0);
375 
376  break;
377  }
378 
379 
380 
381  //---------------------------------------------------------------------
382  // Build a 1D line
383  case 1:
384  {
385  libmesh_assert_not_equal_to (nx, 0);
386  libmesh_assert_equal_to (ny, 0);
387  libmesh_assert_equal_to (nz, 0);
388  libmesh_assert_less (xmin, xmax);
389 
390  // Reserve elements
391  switch (type)
392  {
393  case INVALID_ELEM:
394  case EDGE2:
395  case EDGE3:
396  case EDGE4:
397  {
398  mesh.reserve_elem (nx);
399  break;
400  }
401 
402  default:
403  libmesh_error_msg("ERROR: Unrecognized 1D element type == " << Utility::enum_to_string(type));
404  }
405 
406  // Reserve nodes
407  switch (type)
408  {
409  case INVALID_ELEM:
410  case EDGE2:
411  {
412  mesh.reserve_nodes(nx+1);
413  break;
414  }
415 
416  case EDGE3:
417  {
418  mesh.reserve_nodes(2*nx+1);
419  break;
420  }
421 
422  case EDGE4:
423  {
424  mesh.reserve_nodes(3*nx+1);
425  break;
426  }
427 
428  default:
429  libmesh_error_msg("ERROR: Unrecognized 1D element type == " << Utility::enum_to_string(type));
430  }
431 
432 
433  // Build the nodes, depends on whether we're using linears,
434  // quadratics or cubics and whether using uniform grid or Gauss-Lobatto
435  unsigned int node_id = 0;
436  switch(type)
437  {
438  case INVALID_ELEM:
439  case EDGE2:
440  {
441  for (unsigned int i=0; i<=nx; i++)
442  {
443  const Node * const node = mesh.add_point (Point(static_cast<Real>(i)/nx, 0, 0), node_id++);
444  if (i == 0)
445  boundary_info.add_node(node, 0);
446  if (i == nx)
447  boundary_info.add_node(node, 1);
448  }
449 
450  break;
451  }
452 
453  case EDGE3:
454  {
455  for (unsigned int i=0; i<=2*nx; i++)
456  {
457  const Node * const node = mesh.add_point (Point(static_cast<Real>(i)/(2*nx), 0, 0), node_id++);
458  if (i == 0)
459  boundary_info.add_node(node, 0);
460  if (i == 2*nx)
461  boundary_info.add_node(node, 1);
462  }
463  break;
464  }
465 
466  case EDGE4:
467  {
468  for (unsigned int i=0; i<=3*nx; i++)
469  {
470  const Node * const node = mesh.add_point (Point(static_cast<Real>(i)/(3*nx), 0, 0), node_id++);
471  if (i == 0)
472  boundary_info.add_node(node, 0);
473  if (i == 3*nx)
474  boundary_info.add_node(node, 1);
475  }
476 
477  break;
478  }
479 
480  default:
481  libmesh_error_msg("ERROR: Unrecognized 1D element type == " << Utility::enum_to_string(type));
482 
483  }
484 
485  // Build the elements of the mesh
486  switch(type)
487  {
488  case INVALID_ELEM:
489  case EDGE2:
490  {
491  for (unsigned int i=0; i<nx; i++)
492  {
493  Elem * elem = mesh.add_elem(Elem::build_with_id(EDGE2, i));
494  elem->set_node(0) = mesh.node_ptr(i);
495  elem->set_node(1) = mesh.node_ptr(i+1);
496 
497  if (i == 0)
498  boundary_info.add_side(elem, 0, 0);
499 
500  if (i == (nx-1))
501  boundary_info.add_side(elem, 1, 1);
502 
503  }
504  break;
505  }
506 
507  case EDGE3:
508  {
509  for (unsigned int i=0; i<nx; i++)
510  {
511  Elem * elem = mesh.add_elem(Elem::build_with_id(EDGE3, i));
512  elem->set_node(0) = mesh.node_ptr(2*i);
513  elem->set_node(2) = mesh.node_ptr(2*i+1);
514  elem->set_node(1) = mesh.node_ptr(2*i+2);
515 
516  if (i == 0)
517  boundary_info.add_side(elem, 0, 0);
518 
519  if (i == (nx-1))
520  boundary_info.add_side(elem, 1, 1);
521  }
522  break;
523  }
524 
525  case EDGE4:
526  {
527  for (unsigned int i=0; i<nx; i++)
528  {
529  Elem * elem = mesh.add_elem(Elem::build_with_id(EDGE4, i));
530  elem->set_node(0) = mesh.node_ptr(3*i);
531  elem->set_node(2) = mesh.node_ptr(3*i+1);
532  elem->set_node(3) = mesh.node_ptr(3*i+2);
533  elem->set_node(1) = mesh.node_ptr(3*i+3);
534 
535  if (i == 0)
536  boundary_info.add_side(elem, 0, 0);
537 
538  if (i == (nx-1))
539  boundary_info.add_side(elem, 1, 1);
540  }
541  break;
542  }
543 
544  default:
545  libmesh_error_msg("ERROR: Unrecognized 1D element type == " << Utility::enum_to_string(type));
546  }
547 
548  // Move the nodes to their final locations.
549  if (gauss_lobatto_grid)
550  {
551  GaussLobattoRedistributionFunction func(nx, xmin, xmax);
552  MeshTools::Modification::redistribute(mesh, func);
553  }
554  else // !gauss_lobatto_grid
555  {
556  for (Node * node : mesh.node_ptr_range())
557  (*node)(0) = (*node)(0)*(xmax-xmin) + xmin;
558  }
559 
560  // Add sideset names to boundary info
561  boundary_info.sideset_name(0) = "left";
562  boundary_info.sideset_name(1) = "right";
563 
564  // Add nodeset names to boundary info
565  boundary_info.nodeset_name(0) = "left";
566  boundary_info.nodeset_name(1) = "right";
567 
568  break;
569  }
570 
571 
572 
573 
574 
575 
576 
577 
578 
579 
580  //---------------------------------------------------------------------
581  // Build a 2D quadrilateral
582  case 2:
583  {
584  libmesh_assert_not_equal_to (nx, 0);
585  libmesh_assert_not_equal_to (ny, 0);
586  libmesh_assert_equal_to (nz, 0);
587  libmesh_assert_less (xmin, xmax);
588  libmesh_assert_less (ymin, ymax);
589 
590  // Reserve elements. The TRI3 and TRI6 meshes
591  // have twice as many elements...
592  switch (type)
593  {
594  case INVALID_ELEM:
595  case QUAD4:
596  case QUADSHELL4:
597  case QUAD8:
598  case QUADSHELL8:
599  case QUAD9:
600  {
601  mesh.reserve_elem (nx*ny);
602  break;
603  }
604 
605  case TRI3:
606  case TRI6:
607  case TRI7:
608  {
609  mesh.reserve_elem (2*nx*ny);
610  break;
611  }
612 
613  default:
614  libmesh_error_msg("ERROR: Unrecognized 2D element type == " << Utility::enum_to_string(type));
615  }
616 
617 
618 
619  // Reserve nodes. The quadratic element types
620  // need to reserve more nodes than the linear types.
621  switch (type)
622  {
623  case INVALID_ELEM:
624  case QUAD4:
625  case QUADSHELL4:
626  case TRI3:
627  {
628  mesh.reserve_nodes( (nx+1)*(ny+1) );
629  break;
630  }
631 
632  case QUAD8:
633  case QUADSHELL8:
634  case QUAD9:
635  case TRI6:
636  {
637  mesh.reserve_nodes( (2*nx+1)*(2*ny+1) );
638  break;
639  }
640 
641  case TRI7:
642  {
643  mesh.reserve_nodes( (2*nx+1)*(2*ny+1) + 2*nx*ny );
644  break;
645  }
646 
647  default:
648  libmesh_error_msg("ERROR: Unrecognized 2D element type == " << Utility::enum_to_string(type));
649  }
650 
651 
652 
653  // Build the nodes. Depends on whether you are using a linear
654  // or quadratic element, and whether you are using a uniform
655  // grid or the Gauss-Lobatto grid points.
656  unsigned int node_id = 0;
657  switch (type)
658  {
659  case INVALID_ELEM:
660  case QUAD4:
661  case QUADSHELL4:
662  case TRI3:
663  {
664  for (unsigned int j=0; j<=ny; j++)
665  for (unsigned int i=0; i<=nx; i++)
666  {
667  const Node * const node =
668  mesh.add_point(Point(static_cast<Real>(i) / static_cast<Real>(nx),
669  static_cast<Real>(j) / static_cast<Real>(ny),
670  0.),
671  node_id++);
672  if (j == 0)
673  boundary_info.add_node(node, 0);
674  if (j == ny)
675  boundary_info.add_node(node, 2);
676  if (i == 0)
677  boundary_info.add_node(node, 3);
678  if (i == nx)
679  boundary_info.add_node(node, 1);
680  }
681 
682  break;
683  }
684 
685  case QUAD8:
686  case QUADSHELL8:
687  case QUAD9:
688  case TRI6:
689  case TRI7:
690  {
691  for (unsigned int j=0; j<=(2*ny); j++)
692  for (unsigned int i=0; i<=(2*nx); i++)
693  {
694  const Node * const node =
695  mesh.add_point(Point(static_cast<Real>(i) / static_cast<Real>(2 * nx),
696  static_cast<Real>(j) / static_cast<Real>(2 * ny),
697  0),
698  node_id++);
699  if (j == 0)
700  boundary_info.add_node(node, 0);
701  if (j == 2*ny)
702  boundary_info.add_node(node, 2);
703  if (i == 0)
704  boundary_info.add_node(node, 3);
705  if (i == 2*nx)
706  boundary_info.add_node(node, 1);
707  }
708 
709  // We'll add any interior Tri7 nodes last, to keep from
710  // messing with our idx function
711  if (type == TRI7)
712  for (unsigned int j=0; j<(3*ny); j += 3)
713  for (unsigned int i=0; i<(3*nx); i += 3)
714  {
715  // The bottom-right triangle's center node
716  mesh.add_point(Point(static_cast<Real>(i+2) / static_cast<Real>(3 * nx),
717  static_cast<Real>(j+1) / static_cast<Real>(3 * ny),
718  0),
719  node_id++);
720  // The top-left triangle's center node
721  mesh.add_point(Point(static_cast<Real>(i+1) / static_cast<Real>(3 * nx),
722  static_cast<Real>(j+2) / static_cast<Real>(3 * ny),
723  0),
724  node_id++);
725  }
726 
727  break;
728  }
729 
730 
731  default:
732  libmesh_error_msg("ERROR: Unrecognized 2D element type == " << Utility::enum_to_string(type));
733  }
734 
735 
736 
737 
738 
739 
740  // Build the elements. Each one is a bit different.
741  unsigned int elem_id = 0;
742  switch (type)
743  {
744 
745  case INVALID_ELEM:
746  case QUAD4:
747  case QUADSHELL4:
748  {
749  for (unsigned int j=0; j<ny; j++)
750  for (unsigned int i=0; i<nx; i++)
751  {
752  Elem * elem = mesh.add_elem(Elem::build_with_id(type == INVALID_ELEM ? QUAD4 : type, elem_id++));
753  elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j) );
754  elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+1,j) );
755  elem->set_node(2) = mesh.node_ptr(idx(type,nx,i+1,j+1));
756  elem->set_node(3) = mesh.node_ptr(idx(type,nx,i,j+1) );
757 
758  if (j == 0)
759  boundary_info.add_side(elem, 0, 0);
760 
761  if (j == (ny-1))
762  boundary_info.add_side(elem, 2, 2);
763 
764  if (i == 0)
765  boundary_info.add_side(elem, 3, 3);
766 
767  if (i == (nx-1))
768  boundary_info.add_side(elem, 1, 1);
769  }
770  break;
771  }
772 
773 
774  case TRI3:
775  {
776  for (unsigned int j=0; j<ny; j++)
777  for (unsigned int i=0; i<nx; i++)
778  {
779  // Add first Tri3
780  Elem * elem = mesh.add_elem(Elem::build_with_id(TRI3, elem_id++));
781  elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j) );
782  elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+1,j) );
783  elem->set_node(2) = mesh.node_ptr(idx(type,nx,i+1,j+1));
784 
785  if (j == 0)
786  boundary_info.add_side(elem, 0, 0);
787 
788  if (i == (nx-1))
789  boundary_info.add_side(elem, 1, 1);
790 
791  // Add second Tri3
792  elem = mesh.add_elem(Elem::build_with_id(TRI3, elem_id++));
793  elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j) );
794  elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+1,j+1));
795  elem->set_node(2) = mesh.node_ptr(idx(type,nx,i,j+1) );
796 
797  if (j == (ny-1))
798  boundary_info.add_side(elem, 1, 2);
799 
800  if (i == 0)
801  boundary_info.add_side(elem, 2, 3);
802  }
803  break;
804  }
805 
806 
807 
808  case QUAD8:
809  case QUADSHELL8:
810  case QUAD9:
811  {
812  for (unsigned int j=0; j<(2*ny); j += 2)
813  for (unsigned int i=0; i<(2*nx); i += 2)
814  {
815  Elem * elem = mesh.add_elem(Elem::build_with_id(type, elem_id++));
816  elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j) );
817  elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+2,j) );
818  elem->set_node(2) = mesh.node_ptr(idx(type,nx,i+2,j+2));
819  elem->set_node(3) = mesh.node_ptr(idx(type,nx,i,j+2) );
820  elem->set_node(4) = mesh.node_ptr(idx(type,nx,i+1,j) );
821  elem->set_node(5) = mesh.node_ptr(idx(type,nx,i+2,j+1));
822  elem->set_node(6) = mesh.node_ptr(idx(type,nx,i+1,j+2));
823  elem->set_node(7) = mesh.node_ptr(idx(type,nx,i,j+1) );
824 
825  if (type == QUAD9)
826  elem->set_node(8) = mesh.node_ptr(idx(type,nx,i+1,j+1));
827 
828  if (j == 0)
829  boundary_info.add_side(elem, 0, 0);
830 
831  if (j == 2*(ny-1))
832  boundary_info.add_side(elem, 2, 2);
833 
834  if (i == 0)
835  boundary_info.add_side(elem, 3, 3);
836 
837  if (i == 2*(nx-1))
838  boundary_info.add_side(elem, 1, 1);
839  }
840  break;
841  }
842 
843 
844  case TRI6:
845  case TRI7:
846  {
847  for (unsigned int j=0; j<(2*ny); j += 2)
848  for (unsigned int i=0; i<(2*nx); i += 2)
849  {
850  // Add first Tri in the bottom-right of its quad
851  Elem * elem = mesh.add_elem(Elem::build_with_id(type, elem_id++));
852  elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j) );
853  elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+2,j) );
854  elem->set_node(2) = mesh.node_ptr(idx(type,nx,i+2,j+2));
855  elem->set_node(3) = mesh.node_ptr(idx(type,nx,i+1,j) );
856  elem->set_node(4) = mesh.node_ptr(idx(type,nx,i+2,j+1));
857  elem->set_node(5) = mesh.node_ptr(idx(type,nx,i+1,j+1));
858 
859  if (type == TRI7)
860  elem->set_node(6) = mesh.node_ptr(elem->id()+(2*nx+1)*(2*ny+1));
861 
862  if (j == 0)
863  boundary_info.add_side(elem, 0, 0);
864 
865  if (i == 2*(nx-1))
866  boundary_info.add_side(elem, 1, 1);
867 
868  // Add second Tri in the top left of its quad
869  elem = mesh.add_elem(Elem::build_with_id(type, elem_id++));
870  elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j) );
871  elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+2,j+2));
872  elem->set_node(2) = mesh.node_ptr(idx(type,nx,i,j+2) );
873  elem->set_node(3) = mesh.node_ptr(idx(type,nx,i+1,j+1));
874  elem->set_node(4) = mesh.node_ptr(idx(type,nx,i+1,j+2));
875  elem->set_node(5) = mesh.node_ptr(idx(type,nx,i,j+1) );
876 
877  if (type == TRI7)
878  elem->set_node(6) = mesh.node_ptr(elem->id()+(2*nx+1)*(2*ny+1));
879 
880  if (j == 2*(ny-1))
881  boundary_info.add_side(elem, 1, 2);
882 
883  if (i == 0)
884  boundary_info.add_side(elem, 2, 3);
885  }
886  break;
887  };
888 
889 
890  default:
891  libmesh_error_msg("ERROR: Unrecognized 2D element type == " << Utility::enum_to_string(type));
892  }
893 
894 
895 
896 
897  // Scale the nodal positions
898  if (gauss_lobatto_grid)
899  {
900  GaussLobattoRedistributionFunction func(nx, xmin, xmax,
901  ny, ymin, ymax);
902  MeshTools::Modification::redistribute(mesh, func);
903  }
904  else // !gauss_lobatto_grid
905  {
906  for (Node * node : mesh.node_ptr_range())
907  {
908  (*node)(0) = ((*node)(0))*(xmax-xmin) + xmin;
909  (*node)(1) = ((*node)(1))*(ymax-ymin) + ymin;
910  }
911  }
912 
913  // Add sideset names to boundary info
914  boundary_info.sideset_name(0) = "bottom";
915  boundary_info.sideset_name(1) = "right";
916  boundary_info.sideset_name(2) = "top";
917  boundary_info.sideset_name(3) = "left";
918 
919  // Add nodeset names to boundary info
920  boundary_info.nodeset_name(0) = "bottom";
921  boundary_info.nodeset_name(1) = "right";
922  boundary_info.nodeset_name(2) = "top";
923  boundary_info.nodeset_name(3) = "left";
924 
925  break;
926  }
927 
928 
929 
930 
931 
932 
933 
934 
935 
936 
937 
938  //---------------------------------------------------------------------
939  // Build a 3D mesh using hexes, tets, prisms, or pyramids.
940  case 3:
941  {
942  libmesh_assert_not_equal_to (nx, 0);
943  libmesh_assert_not_equal_to (ny, 0);
944  libmesh_assert_not_equal_to (nz, 0);
945  libmesh_assert_less (xmin, xmax);
946  libmesh_assert_less (ymin, ymax);
947  libmesh_assert_less (zmin, zmax);
948 
949 
950  // Reserve elements. Meshes with prismatic elements require
951  // twice as many elements.
952  switch (type)
953  {
954  case INVALID_ELEM:
955  case HEX8:
956  case HEX20:
957  case HEX27:
958  case TET4: // TET4's are created from an initial HEX27 discretization
959  case TET10: // TET10's are created from an initial HEX27 discretization
960  case TET14: // TET14's are created from an initial HEX27 discretization
961  case PYRAMID5: // PYRAMIDs are created from an initial HEX27 discretization
962  case PYRAMID13:
963  case PYRAMID14:
964  case PYRAMID18:
965  {
966  mesh.reserve_elem(nx*ny*nz);
967  break;
968  }
969 
970  case PRISM6:
971  case PRISM15:
972  case PRISM18:
973  case PRISM20:
974  case PRISM21:
975  {
976  mesh.reserve_elem(2*nx*ny*nz);
977  break;
978  }
979 
980  default:
981  libmesh_error_msg("ERROR: Unrecognized 3D element type == " << Utility::enum_to_string(type));
982  }
983 
984 
985 
986 
987 
988  // Reserve nodes. Quadratic elements need twice as many nodes as linear elements.
989  switch (type)
990  {
991  case INVALID_ELEM:
992  case HEX8:
993  case PRISM6:
994  {
995  mesh.reserve_nodes( (nx+1)*(ny+1)*(nz+1) );
996  break;
997  }
998 
999  case HEX20:
1000  case HEX27:
1001  case TET4: // TET4's are created from an initial HEX27 discretization
1002  case TET10: // TET10's are created from an initial HEX27 discretization
1003  case PYRAMID5: // PYRAMIDs are created from an initial HEX27 discretization
1004  case PYRAMID13:
1005  case PYRAMID14:
1006  case PYRAMID18:
1007  case PRISM15:
1008  case PRISM18:
1009  {
1010  // FYI: The resulting TET4 mesh will have exactly
1011  // 5*(nx*ny*nz) + 2*(nx*ny + nx*nz + ny*nz) + (nx+ny+nz) + 1
1012  // nodes once the additional mid-edge nodes for the HEX27 discretization
1013  // have been deleted.
1014  mesh.reserve_nodes( (2*nx+1)*(2*ny+1)*(2*nz+1) );
1015  break;
1016  }
1017 
1018  case TET14:
1019  {
1020  mesh.reserve_nodes( (2*nx+1)*(2*ny+1)*(2*nz+1) +
1021  24*nx*ny*nz +
1022  4*(nx*ny + ny*nz + nx*nz) );
1023  break;
1024  }
1025 
1026  case PRISM20:
1027  {
1028  mesh.reserve_nodes( (2*nx+1)*(2*ny+1)*(2*nz+1) +
1029  2*nx*ny*(nz+1) );
1030  break;
1031  }
1032 
1033  case PRISM21:
1034  {
1035  mesh.reserve_nodes( (2*nx+1)*(2*ny+1)*(2*nz+1) +
1036  2*nx*ny*(2*nz+1) );
1037  break;
1038  }
1039 
1040  default:
1041  libmesh_error_msg("ERROR: Unrecognized 3D element type == " << Utility::enum_to_string(type));
1042  }
1043 
1044 
1045 
1046 
1047  // Build the nodes.
1048  unsigned int node_id = 0;
1049  switch (type)
1050  {
1051  case INVALID_ELEM:
1052  case HEX8:
1053  case PRISM6:
1054  {
1055  for (unsigned int k=0; k<=nz; k++)
1056  for (unsigned int j=0; j<=ny; j++)
1057  for (unsigned int i=0; i<=nx; i++)
1058  {
1059  const Node * const node =
1060  mesh.add_point(Point(static_cast<Real>(i) / static_cast<Real>(nx),
1061  static_cast<Real>(j) / static_cast<Real>(ny),
1062  static_cast<Real>(k) / static_cast<Real>(nz)),
1063  node_id++);
1064  if (k == 0)
1065  boundary_info.add_node(node, 0);
1066  if (k == nz)
1067  boundary_info.add_node(node, 5);
1068  if (j == 0)
1069  boundary_info.add_node(node, 1);
1070  if (j == ny)
1071  boundary_info.add_node(node, 3);
1072  if (i == 0)
1073  boundary_info.add_node(node, 4);
1074  if (i == nx)
1075  boundary_info.add_node(node, 2);
1076  }
1077 
1078  break;
1079  }
1080 
1081  case HEX20:
1082  case HEX27:
1083  case TET4: // TET4's are created from an initial HEX27 discretization
1084  case TET10: // TET10's are created from an initial HEX27 discretization
1085  case TET14: // TET14's are created from an initial HEX27 discretization
1086  case PYRAMID5: // PYRAMIDs are created from an initial HEX27 discretization
1087  case PYRAMID13:
1088  case PYRAMID14:
1089  case PYRAMID18:
1090  case PRISM15:
1091  case PRISM18:
1092  case PRISM20:
1093  case PRISM21:
1094  {
1095  for (unsigned int k=0; k<=(2*nz); k++)
1096  for (unsigned int j=0; j<=(2*ny); j++)
1097  for (unsigned int i=0; i<=(2*nx); i++)
1098  {
1099  const Node * const node =
1100  mesh.add_point(Point(static_cast<Real>(i) / static_cast<Real>(2 * nx),
1101  static_cast<Real>(j) / static_cast<Real>(2 * ny),
1102  static_cast<Real>(k) / static_cast<Real>(2 * nz)),
1103  node_id++);
1104  if (k == 0)
1105  boundary_info.add_node(node, 0);
1106  if (k == 2*nz)
1107  boundary_info.add_node(node, 5);
1108  if (j == 0)
1109  boundary_info.add_node(node, 1);
1110  if (j == 2*ny)
1111  boundary_info.add_node(node, 3);
1112  if (i == 0)
1113  boundary_info.add_node(node, 4);
1114  if (i == 2*nx)
1115  boundary_info.add_node(node, 2);
1116  }
1117 
1118  if (type == PRISM20 ||
1119  type == PRISM21)
1120  {
1121  const unsigned int kmax = (type == PRISM20) ? nz : 2*nz;
1122  for (unsigned int k=0; k<=kmax; k++)
1123  for (unsigned int j=0; j<ny; j++)
1124  for (unsigned int i=0; i<nx; i++)
1125  {
1126  const Node * const node1 =
1127  mesh.add_point(Point((static_cast<Real>(i)+1/Real(3)) / static_cast<Real>(nx),
1128  (static_cast<Real>(j)+1/Real(3)) / static_cast<Real>(ny),
1129  static_cast<Real>(k) / static_cast<Real>(kmax)),
1130  node_id++);
1131  if (k == 0)
1132  boundary_info.add_node(node1, 0);
1133  if (k == kmax)
1134  boundary_info.add_node(node1, 5);
1135 
1136  const Node * const node2 =
1137  mesh.add_point(Point((static_cast<Real>(i)+2/Real(3)) / static_cast<Real>(nx),
1138  (static_cast<Real>(j)+2/Real(3)) / static_cast<Real>(ny),
1139  static_cast<Real>(k) / static_cast<Real>(kmax)),
1140  node_id++);
1141  if (k == 0)
1142  boundary_info.add_node(node2, 0);
1143  if (k == kmax)
1144  boundary_info.add_node(node2, 5);
1145  }
1146  }
1147 
1148  break;
1149  }
1150 
1151 
1152  default:
1153  libmesh_error_msg("ERROR: Unrecognized 3D element type == " << Utility::enum_to_string(type));
1154  }
1155 
1156 
1157 
1158 
1159  // Build the elements.
1160  unsigned int elem_id = 0;
1161  switch (type)
1162  {
1163  case INVALID_ELEM:
1164  case HEX8:
1165  {
1166  for (unsigned int k=0; k<nz; k++)
1167  for (unsigned int j=0; j<ny; j++)
1168  for (unsigned int i=0; i<nx; i++)
1169  {
1170  Elem * elem = mesh.add_elem(Elem::build_with_id(HEX8, elem_id++));
1171  elem->set_node(0) = mesh.node_ptr(idx(type,nx,ny,i,j,k) );
1172  elem->set_node(1) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k) );
1173  elem->set_node(2) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k) );
1174  elem->set_node(3) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k) );
1175  elem->set_node(4) = mesh.node_ptr(idx(type,nx,ny,i,j,k+1) );
1176  elem->set_node(5) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k+1) );
1177  elem->set_node(6) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
1178  elem->set_node(7) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k+1) );
1179 
1180  if (k == 0)
1181  boundary_info.add_side(elem, 0, 0);
1182 
1183  if (k == (nz-1))
1184  boundary_info.add_side(elem, 5, 5);
1185 
1186  if (j == 0)
1187  boundary_info.add_side(elem, 1, 1);
1188 
1189  if (j == (ny-1))
1190  boundary_info.add_side(elem, 3, 3);
1191 
1192  if (i == 0)
1193  boundary_info.add_side(elem, 4, 4);
1194 
1195  if (i == (nx-1))
1196  boundary_info.add_side(elem, 2, 2);
1197  }
1198  break;
1199  }
1200 
1201 
1202 
1203 
1204  case PRISM6:
1205  {
1206  for (unsigned int k=0; k<nz; k++)
1207  for (unsigned int j=0; j<ny; j++)
1208  for (unsigned int i=0; i<nx; i++)
1209  {
1210  // First Prism
1211  Elem * elem = mesh.add_elem(Elem::build_with_id(PRISM6, elem_id++));
1212  elem->set_node(0) = mesh.node_ptr(idx(type,nx,ny,i,j,k) );
1213  elem->set_node(1) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k) );
1214  elem->set_node(2) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k) );
1215  elem->set_node(3) = mesh.node_ptr(idx(type,nx,ny,i,j,k+1) );
1216  elem->set_node(4) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k+1) );
1217  elem->set_node(5) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k+1) );
1218 
1219  // Add sides for first prism to boundary info object
1220  if (i==0)
1221  boundary_info.add_side(elem, 3, 4);
1222 
1223  if (j==0)
1224  boundary_info.add_side(elem, 1, 1);
1225 
1226  if (k==0)
1227  boundary_info.add_side(elem, 0, 0);
1228 
1229  if (k == (nz-1))
1230  boundary_info.add_side(elem, 4, 5);
1231 
1232  // Second Prism
1233  elem = mesh.add_elem(Elem::build_with_id(PRISM6, elem_id++));
1234  elem->set_node(0) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k) );
1235  elem->set_node(1) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k) );
1236  elem->set_node(2) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k) );
1237  elem->set_node(3) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k+1) );
1238  elem->set_node(4) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
1239  elem->set_node(5) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k+1) );
1240 
1241  // Add sides for second prism to boundary info object
1242  if (i == (nx-1))
1243  boundary_info.add_side(elem, 1, 2);
1244 
1245  if (j == (ny-1))
1246  boundary_info.add_side(elem, 2, 3);
1247 
1248  if (k==0)
1249  boundary_info.add_side(elem, 0, 0);
1250 
1251  if (k == (nz-1))
1252  boundary_info.add_side(elem, 4, 5);
1253  }
1254  break;
1255  }
1256 
1257 
1258 
1259 
1260 
1261 
1262  case HEX20:
1263  case HEX27:
1264  case TET4: // TET4's are created from an initial HEX27 discretization
1265  case TET10: // TET10's are created from an initial HEX27 discretization
1266  case TET14: // TET14's are created from an initial HEX27 discretization
1267  case PYRAMID5: // PYRAMIDs are created from an initial HEX27 discretization
1268  case PYRAMID13:
1269  case PYRAMID14:
1270  case PYRAMID18:
1271  {
1272  for (unsigned int k=0; k<(2*nz); k += 2)
1273  for (unsigned int j=0; j<(2*ny); j += 2)
1274  for (unsigned int i=0; i<(2*nx); i += 2)
1275  {
1276  ElemType build_type = (type == HEX20) ? HEX20 : HEX27;
1277  Elem * elem = mesh.add_elem(Elem::build_with_id(build_type, elem_id++));
1278 
1279  elem->set_node(0) = mesh.node_ptr(idx(type,nx,ny,i, j, k) );
1280  elem->set_node(1) = mesh.node_ptr(idx(type,nx,ny,i+2,j, k) );
1281  elem->set_node(2) = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k) );
1282  elem->set_node(3) = mesh.node_ptr(idx(type,nx,ny,i, j+2,k) );
1283  elem->set_node(4) = mesh.node_ptr(idx(type,nx,ny,i, j, k+2));
1284  elem->set_node(5) = mesh.node_ptr(idx(type,nx,ny,i+2,j, k+2));
1285  elem->set_node(6) = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k+2));
1286  elem->set_node(7) = mesh.node_ptr(idx(type,nx,ny,i, j+2,k+2));
1287  elem->set_node(8) = mesh.node_ptr(idx(type,nx,ny,i+1,j, k) );
1288  elem->set_node(9) = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k) );
1289  elem->set_node(10) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k) );
1290  elem->set_node(11) = mesh.node_ptr(idx(type,nx,ny,i, j+1,k) );
1291  elem->set_node(12) = mesh.node_ptr(idx(type,nx,ny,i, j, k+1));
1292  elem->set_node(13) = mesh.node_ptr(idx(type,nx,ny,i+2,j, k+1));
1293  elem->set_node(14) = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k+1));
1294  elem->set_node(15) = mesh.node_ptr(idx(type,nx,ny,i, j+2,k+1));
1295  elem->set_node(16) = mesh.node_ptr(idx(type,nx,ny,i+1,j, k+2));
1296  elem->set_node(17) = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k+2));
1297  elem->set_node(18) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k+2));
1298  elem->set_node(19) = mesh.node_ptr(idx(type,nx,ny,i, j+1,k+2));
1299 
1300  if ((type == HEX27) || (type == TET4) || (type == TET10) || (type == TET14) ||
1301  (type == PYRAMID5) || (type == PYRAMID13) || (type == PYRAMID14) ||
1302  (type == PYRAMID18))
1303  {
1304  elem->set_node(20) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k) );
1305  elem->set_node(21) = mesh.node_ptr(idx(type,nx,ny,i+1,j, k+1));
1306  elem->set_node(22) = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k+1));
1307  elem->set_node(23) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k+1));
1308  elem->set_node(24) = mesh.node_ptr(idx(type,nx,ny,i, j+1,k+1));
1309  elem->set_node(25) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+2));
1310  elem->set_node(26) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
1311  }
1312 
1313  if (k == 0)
1314  boundary_info.add_side(elem, 0, 0);
1315 
1316  if (k == 2*(nz-1))
1317  boundary_info.add_side(elem, 5, 5);
1318 
1319  if (j == 0)
1320  boundary_info.add_side(elem, 1, 1);
1321 
1322  if (j == 2*(ny-1))
1323  boundary_info.add_side(elem, 3, 3);
1324 
1325  if (i == 0)
1326  boundary_info.add_side(elem, 4, 4);
1327 
1328  if (i == 2*(nx-1))
1329  boundary_info.add_side(elem, 2, 2);
1330  }
1331  break;
1332  }
1333 
1334 
1335 
1336 
1337  case PRISM15:
1338  case PRISM18:
1339  case PRISM20:
1340  case PRISM21:
1341  {
1342  for (unsigned int k=0; k<(2*nz); k += 2)
1343  for (unsigned int j=0; j<(2*ny); j += 2)
1344  for (unsigned int i=0; i<(2*nx); i += 2)
1345  {
1346  // First Prism
1347  Elem * elem = mesh.add_elem(Elem::build_with_id(type, elem_id++));
1348  elem->set_node(0) = mesh.node_ptr(idx(type,nx,ny,i, j, k) );
1349  elem->set_node(1) = mesh.node_ptr(idx(type,nx,ny,i+2,j, k) );
1350  elem->set_node(2) = mesh.node_ptr(idx(type,nx,ny,i, j+2,k) );
1351  elem->set_node(3) = mesh.node_ptr(idx(type,nx,ny,i, j, k+2));
1352  elem->set_node(4) = mesh.node_ptr(idx(type,nx,ny,i+2,j, k+2));
1353  elem->set_node(5) = mesh.node_ptr(idx(type,nx,ny,i, j+2,k+2));
1354  elem->set_node(6) = mesh.node_ptr(idx(type,nx,ny,i+1,j, k) );
1355  elem->set_node(7) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k) );
1356  elem->set_node(8) = mesh.node_ptr(idx(type,nx,ny,i, j+1,k) );
1357  elem->set_node(9) = mesh.node_ptr(idx(type,nx,ny,i, j, k+1));
1358  elem->set_node(10) = mesh.node_ptr(idx(type,nx,ny,i+2,j, k+1));
1359  elem->set_node(11) = mesh.node_ptr(idx(type,nx,ny,i, j+2,k+1));
1360  elem->set_node(12) = mesh.node_ptr(idx(type,nx,ny,i+1,j, k+2));
1361  elem->set_node(13) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+2));
1362  elem->set_node(14) = mesh.node_ptr(idx(type,nx,ny,i, j+1,k+2));
1363 
1364  if (type == PRISM18 ||
1365  type == PRISM20 ||
1366  type == PRISM21)
1367  {
1368  elem->set_node(15) = mesh.node_ptr(idx(type,nx,ny,i+1,j, k+1));
1369  elem->set_node(16) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
1370  elem->set_node(17) = mesh.node_ptr(idx(type,nx,ny,i, j+1,k+1));
1371  }
1372 
1373  if (type == PRISM20)
1374  {
1375  const dof_id_type base_idx = (2*nx+1)*(2*ny+1)*(2*nz+1);
1376  elem->set_node(18) = mesh.node_ptr(base_idx+((k/2)*(nx*ny)+j/2*nx+i/2)*2);
1377  elem->set_node(19) = mesh.node_ptr(base_idx+(((k/2)+1)*(nx*ny)+j/2*nx+i/2)*2);
1378  }
1379 
1380  if (type == PRISM21)
1381  {
1382  const dof_id_type base_idx = (2*nx+1)*(2*ny+1)*(2*nz+1);
1383  elem->set_node(18) = mesh.node_ptr(base_idx+(k*(nx*ny)+j/2*nx+i/2)*2);
1384  elem->set_node(19) = mesh.node_ptr(base_idx+((k+2)*(nx*ny)+j/2*nx+i/2)*2);
1385  elem->set_node(20) = mesh.node_ptr(base_idx+((k+1)*(nx*ny)+j/2*nx+i/2)*2);
1386  }
1387 
1388  // Add sides for first prism to boundary info object
1389  if (i==0)
1390  boundary_info.add_side(elem, 3, 4);
1391 
1392  if (j==0)
1393  boundary_info.add_side(elem, 1, 1);
1394 
1395  if (k==0)
1396  boundary_info.add_side(elem, 0, 0);
1397 
1398  if (k == 2*(nz-1))
1399  boundary_info.add_side(elem, 4, 5);
1400 
1401 
1402  // Second Prism
1403  elem = mesh.add_elem(Elem::build_with_id(type, elem_id++));
1404  elem->set_node(0) = mesh.node_ptr(idx(type,nx,ny,i+2,j,k) );
1405  elem->set_node(1) = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k) );
1406  elem->set_node(2) = mesh.node_ptr(idx(type,nx,ny,i,j+2,k) );
1407  elem->set_node(3) = mesh.node_ptr(idx(type,nx,ny,i+2,j,k+2) );
1408  elem->set_node(4) = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k+2) );
1409  elem->set_node(5) = mesh.node_ptr(idx(type,nx,ny,i,j+2,k+2) );
1410  elem->set_node(6) = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k) );
1411  elem->set_node(7) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k) );
1412  elem->set_node(8) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k) );
1413  elem->set_node(9) = mesh.node_ptr(idx(type,nx,ny,i+2,j,k+1) );
1414  elem->set_node(10) = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k+1));
1415  elem->set_node(11) = mesh.node_ptr(idx(type,nx,ny,i,j+2,k+1) );
1416  elem->set_node(12) = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k+2));
1417  elem->set_node(13) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k+2));
1418  elem->set_node(14) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+2));
1419 
1420  if (type == PRISM18 ||
1421  type == PRISM20 ||
1422  type == PRISM21)
1423  {
1424  elem->set_node(15) = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k+1));
1425  elem->set_node(16) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k+1));
1426  elem->set_node(17) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
1427  }
1428 
1429  if (type == PRISM20)
1430  {
1431  const dof_id_type base_idx = (2*nx+1)*(2*ny+1)*(2*nz+1);
1432  elem->set_node(18) = mesh.node_ptr(base_idx+((k/2)*(nx*ny)+j/2*nx+i/2)*2+1);
1433  elem->set_node(19) = mesh.node_ptr(base_idx+(((k/2)+1)*(nx*ny)+j/2*nx+i/2)*2+1);
1434  }
1435 
1436  if (type == PRISM21)
1437  {
1438  const dof_id_type base_idx = (2*nx+1)*(2*ny+1)*(2*nz+1);
1439  elem->set_node(18) = mesh.node_ptr(base_idx+(k*(nx*ny)+j/2*nx+i/2)*2+1);
1440  elem->set_node(19) = mesh.node_ptr(base_idx+((k+2)*(nx*ny)+j/2*nx+i/2)*2+1);
1441  elem->set_node(20) = mesh.node_ptr(base_idx+((k+1)*(nx*ny)+j/2*nx+i/2)*2+1);
1442  }
1443 
1444  // Add sides for second prism to boundary info object
1445  if (i == 2*(nx-1))
1446  boundary_info.add_side(elem, 1, 2);
1447 
1448  if (j == 2*(ny-1))
1449  boundary_info.add_side(elem, 2, 3);
1450 
1451  if (k==0)
1452  boundary_info.add_side(elem, 0, 0);
1453 
1454  if (k == 2*(nz-1))
1455  boundary_info.add_side(elem, 4, 5);
1456 
1457  }
1458  break;
1459  }
1460 
1461 
1462 
1463 
1464 
1465  default:
1466  libmesh_error_msg("ERROR: Unrecognized 3D element type == " << Utility::enum_to_string(type));
1467  }
1468 
1469 
1470 
1471 
1472  //.......................................
1473  // Scale the nodal positions
1474  if (gauss_lobatto_grid)
1475  {
1476  GaussLobattoRedistributionFunction func(nx, xmin, xmax,
1477  ny, ymin, ymax,
1478  nz, zmin, zmax);
1479  MeshTools::Modification::redistribute(mesh, func);
1480  }
1481  else // !gauss_lobatto_grid
1482  {
1483  for (unsigned int p=0; p<mesh.n_nodes(); p++)
1484  {
1485  mesh.node_ref(p)(0) = (mesh.node_ref(p)(0))*(xmax-xmin) + xmin;
1486  mesh.node_ref(p)(1) = (mesh.node_ref(p)(1))*(ymax-ymin) + ymin;
1487  mesh.node_ref(p)(2) = (mesh.node_ref(p)(2))*(zmax-zmin) + zmin;
1488  }
1489  }
1490 
1491 
1492 
1493  // Additional work for tets and pyramids: we take the existing
1494  // HEX27 discretization and split each element into 24
1495  // sub-tets or 6 sub-pyramids.
1496  //
1497  // 24 isn't the minimum-possible number of tets, but it
1498  // obviates any concerns about the edge orientations between
1499  // the various elements.
1500  if ((type == TET4) ||
1501  (type == TET10) ||
1502  (type == TET14) ||
1503  (type == PYRAMID5) ||
1504  (type == PYRAMID13) ||
1505  (type == PYRAMID14) ||
1506  (type == PYRAMID18))
1507  {
1508  // Temporary storage for new elements. (24 tets per hex, 6 pyramids)
1509  std::vector<std::unique_ptr<Elem>> new_elements;
1510 
1511  // For avoiding extraneous construction of element sides
1512  std::unique_ptr<Elem> side;
1513 
1514  if ((type == TET4) || (type == TET10) || (type == TET14))
1515  new_elements.reserve(24*mesh.n_elem());
1516  else
1517  new_elements.reserve(6*mesh.n_elem());
1518 
1519  // Create tetrahedra or pyramids
1520  for (auto & base_hex : mesh.element_ptr_range())
1521  {
1522  // Get a pointer to the node located at the HEX27 center
1523  Node * apex_node = base_hex->node_ptr(26);
1524 
1525  // Container to catch ids handed back from BoundaryInfo
1526  std::vector<boundary_id_type> ids;
1527 
1528  for (auto s : base_hex->side_index_range())
1529  {
1530  // Get the boundary ID(s) for this side
1531  boundary_info.boundary_ids(base_hex, s, ids);
1532 
1533  // We're creating this Mesh, so there should be 0 or 1 boundary IDs.
1534  libmesh_assert(ids.size() <= 1);
1535 
1536  // A convenient name for the side's ID.
1537  boundary_id_type b_id = ids.empty() ? BoundaryInfo::invalid_id : ids[0];
1538 
1539  // Need to build the full-ordered side!
1540  base_hex->build_side_ptr(side, s);
1541 
1542  if ((type == TET4) || (type == TET10) || (type == TET14))
1543  {
1544  // Build 4 sub-tets per side
1545  for (unsigned int sub_tet=0; sub_tet<4; ++sub_tet)
1546  {
1547  new_elements.push_back( Elem::build(TET4) );
1548  auto & sub_elem = new_elements.back();
1549  sub_elem->set_node(0) = side->node_ptr(sub_tet);
1550  sub_elem->set_node(1) = side->node_ptr(8); // center of the face
1551  sub_elem->set_node(2) = side->node_ptr(sub_tet==3 ? 0 : sub_tet+1 ); // wrap-around
1552  sub_elem->set_node(3) = apex_node; // apex node always used!
1553 
1554  // If the original hex was a boundary hex, add the new sub_tet's side
1555  // 0 with the same b_id. Note: the tets are all aligned so that their
1556  // side 0 is on the boundary.
1557  if (b_id != BoundaryInfo::invalid_id)
1558  boundary_info.add_side(sub_elem.get(), 0, b_id);
1559  }
1560  } // end if ((type == TET4) || (type == TET10) || (type == TET14))
1561 
1562  else // type==PYRAMID*
1563  {
1564  // Build 1 sub-pyramid per side.
1565  new_elements.push_back( Elem::build(PYRAMID5) );
1566  auto & sub_elem = new_elements.back();
1567 
1568  // Set the base. Note that since the apex is *inside* the base_hex,
1569  // and the pyramid uses a counter-clockwise base numbering, we need to
1570  // reverse the [1] and [3] node indices.
1571  sub_elem->set_node(0) = side->node_ptr(0);
1572  sub_elem->set_node(1) = side->node_ptr(3);
1573  sub_elem->set_node(2) = side->node_ptr(2);
1574  sub_elem->set_node(3) = side->node_ptr(1);
1575 
1576  // Set the apex
1577  sub_elem->set_node(4) = apex_node;
1578 
1579  // If the original hex was a boundary hex, add the new sub_pyr's side
1580  // 4 (the square base) with the same b_id.
1581  if (b_id != BoundaryInfo::invalid_id)
1582  boundary_info.add_side(sub_elem.get(), 4, b_id);
1583  } // end else type==PYRAMID*
1584  }
1585  }
1586 
1587 
1588  // Delete the original HEX27 elements from the mesh, and the boundary info structure.
1589  for (auto & elem : mesh.element_ptr_range())
1590  {
1591  boundary_info.remove(elem); // Safe even if elem has no boundary info.
1592  mesh.delete_elem(elem);
1593  }
1594 
1595  // Add the new elements
1596  for (auto i : index_range(new_elements))
1597  {
1598  new_elements[i]->set_id(i);
1599  mesh.add_elem( std::move(new_elements[i]) );
1600  }
1601 
1602  } // end if (type == TET*,PYRAMID*)
1603 
1604 
1605  // Use all_second_order to convert the TET4's to TET10's or PYRAMID5's to PYRAMID14's
1606  if ((type == TET10) || (type == PYRAMID14))
1607  mesh.all_second_order();
1608 
1609  else if (type == PYRAMID13)
1610  mesh.all_second_order(/*full_ordered=*/false);
1611 
1612  else if ((type == TET14) || (type == PYRAMID18))
1613  mesh.all_complete_order();
1614 
1615 
1616  // Add sideset names to boundary info (Z axis out of the screen)
1617  boundary_info.sideset_name(0) = "back";
1618  boundary_info.sideset_name(1) = "bottom";
1619  boundary_info.sideset_name(2) = "right";
1620  boundary_info.sideset_name(3) = "top";
1621  boundary_info.sideset_name(4) = "left";
1622  boundary_info.sideset_name(5) = "front";
1623 
1624  // Add nodeset names to boundary info
1625  boundary_info.nodeset_name(0) = "back";
1626  boundary_info.nodeset_name(1) = "bottom";
1627  boundary_info.nodeset_name(2) = "right";
1628  boundary_info.nodeset_name(3) = "top";
1629  boundary_info.nodeset_name(4) = "left";
1630  boundary_info.nodeset_name(5) = "front";
1631 
1632  break;
1633  } // end case dim==3
1634 
1635  default:
1636  libmesh_error_msg("Unknown dimension " << mesh.mesh_dimension());
1637  }
1638 
1639  // Done building the mesh. Now prepare it for use.
1640  mesh.prepare_for_use ();
1641 }
libMesh::MeshTools::Generation::Private::idx
unsigned int idx(const ElemType type, const unsigned int nx, const unsigned int ny, const unsigned int i, const unsigned int j, const unsigned int k)
Definition: mesh_generation.C:113
libMesh::ElemType
ElemType
Defines an enum for geometric element types.
Definition: enum_elem_type.h:33
mesh
MeshBase & mesh
Definition: mesh_communication.C:1308
libMesh::boundary_id_type
int8_t boundary_id_type
Definition: id_types.h:51
libMesh::libmesh_assert
libmesh_assert(ctx)
libMesh::Utility::enum_to_string
std::string enum_to_string(const T e)
libMesh::Real
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Definition: libmesh_common.h:143
libMesh::MeshTools::Modification::redistribute
void redistribute(MeshBase &mesh, const FunctionBase< Real > &mapfunc)
Deterministically perturb the nodal locations.
Definition: mesh_modification.C:207
libMesh::index_range
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
libMesh::dof_id_type
uint8_t dof_id_type
Definition: id_types.h:67

◆ build_delaunay_square()

void libMesh::MeshTools::Generation::build_delaunay_square ( UnstructuredMesh &  mesh,
const unsigned int  nx,
const unsigned int  ny,
const Real  xmin,
const Real  xmax,
const Real  ymin,
const Real  ymax,
const ElemType  type,
const std::vector< TriangleInterface::Hole *> *  holes = nullptr 
)

Meshes a rectangular (2D) region (with or without holes) with a Delaunay triangulation.

This function internally calls the triangle library written by J.R. Shewchuk.

Definition at line 2625 of file mesh_generation.C.

References libMesh::MeshBase::add_point(), libMesh::BoundaryInfo::add_side(), libMesh::TriangulatorInterface::attach_hole_list(), libMesh::MeshBase::clear(), libMesh::TriangulatorInterface::desired_area(), libMesh::TriangulatorInterface::elem_type(), libMesh::MeshBase::get_boundary_info(), mesh, libMesh::MeshBase::n_nodes(), libMesh::TriangulatorInterface::PSLG, libMesh::Real, libMesh::MeshBase::set_mesh_dimension(), libMesh::TOLERANCE, libMesh::TriangleInterface::triangulate(), and libMesh::TriangulatorInterface::triangulation_type().

2632 {
2633  // Check for reasonable size
2634  libmesh_assert_greater_equal (nx, 1); // need at least 1 element in x-direction
2635  libmesh_assert_greater_equal (ny, 1); // need at least 1 element in y-direction
2636  libmesh_assert_less (xmin, xmax);
2637  libmesh_assert_less (ymin, ymax);
2638 
2639  // Clear out any data which may have been in the Mesh
2640  mesh.clear();
2641 
2642  BoundaryInfo & boundary_info = mesh.get_boundary_info();
2643 
2644  // Make sure the new Mesh will be 2D
2645  mesh.set_mesh_dimension(2);
2646 
2647  // The x and y spacing between boundary points
2648  const Real delta_x = (xmax-xmin) / static_cast<Real>(nx);
2649  const Real delta_y = (ymax-ymin) / static_cast<Real>(ny);
2650 
2651  // Bottom
2652  for (unsigned int p=0; p<=nx; ++p)
2653  mesh.add_point(Point(xmin + p*delta_x, ymin));
2654 
2655  // Right side
2656  for (unsigned int p=1; p<ny; ++p)
2657  mesh.add_point(Point(xmax, ymin + p*delta_y));
2658 
2659  // Top
2660  for (unsigned int p=0; p<=nx; ++p)
2661  mesh.add_point(Point(xmax - p*delta_x, ymax));
2662 
2663  // Left side
2664  for (unsigned int p=1; p<ny; ++p)
2665  mesh.add_point(Point(xmin, ymax - p*delta_y));
2666 
2667  // Be sure we added as many points as we thought we did
2668  libmesh_assert_equal_to (mesh.n_nodes(), 2*(nx+ny));
2669 
2670  // Construct the Triangle Interface object
2671  TriangleInterface t(mesh);
2672 
2673  // Set custom variables for the triangulation
2674  t.desired_area() = 0.5 * (xmax-xmin)*(ymax-ymin) / static_cast<Real>(nx*ny);
2675  t.triangulation_type() = TriangleInterface::PSLG;
2676  t.elem_type() = type;
2677 
2678  if (holes != nullptr)
2679  t.attach_hole_list(holes);
2680 
2681  // Triangulate!
2682  t.triangulate();
2683 
2684  // For avoiding extraneous side element construction
2685  std::unique_ptr<const Elem> side;
2686 
2687  // The mesh is now generated, but we still need to mark the boundaries
2688  // to be consistent with the other build_square routines. Note that all
2689  // hole boundary elements get the same ID, 4.
2690  for (auto & elem : mesh.element_ptr_range())
2691  for (auto s : elem->side_index_range())
2692  if (elem->neighbor_ptr(s) == nullptr)
2693  {
2694  elem->build_side_ptr(side, s);
2695 
2696  // Check the location of the side's midpoint. Since
2697  // the square has straight sides, the midpoint is not
2698  // on the corner and thus it is uniquely on one of the
2699  // sides.
2700  Point side_midpoint= 0.5f*( side->point(0) + side->point(1) );
2701 
2702  // The boundary ids are set following the same convention as Quad4 sides
2703  // bottom = 0
2704  // right = 1
2705  // top = 2
2706  // left = 3
2707  // hole = 4
2708  boundary_id_type bc_id=4;
2709 
2710  // bottom
2711  if (std::fabs(side_midpoint(1) - ymin) < TOLERANCE)
2712  bc_id=0;
2713 
2714  // right
2715  else if (std::fabs(side_midpoint(0) - xmax) < TOLERANCE)
2716  bc_id=1;
2717 
2718  // top
2719  else if (std::fabs(side_midpoint(1) - ymax) < TOLERANCE)
2720  bc_id=2;
2721 
2722  // left
2723  else if (std::fabs(side_midpoint(0) - xmin) < TOLERANCE)
2724  bc_id=3;
2725 
2726  // If the point is not on any of the external boundaries, it
2727  // is on one of the holes....
2728 
2729  // Finally, add this element's information to the boundary info object.
2730  boundary_info.add_side(elem->id(), s, bc_id);
2731  }
2732 
2733 } // end build_delaunay_square
libMesh::TOLERANCE
static constexpr Real TOLERANCE
Definition: libmesh_common.h:150
mesh
MeshBase & mesh
Definition: mesh_communication.C:1308
libMesh::boundary_id_type
int8_t boundary_id_type
Definition: id_types.h:51
libMesh::Real
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Definition: libmesh_common.h:143

◆ build_extrusion()

void libMesh::MeshTools::Generation::build_extrusion ( UnstructuredMesh &  mesh,
const MeshBase &  cross_section,
const unsigned int  nz,
RealVectorValue  extrusion_vector,
QueryElemSubdomainIDBase *  elem_subdomain = nullptr 
)

Meshes the tensor product of a 1D and a 1D-or-2D domain.

Definition at line 2252 of file mesh_generation.C.

References libMesh::MeshBase::add_elem(), libMesh::BoundaryInfo::add_node(), libMesh::MeshBase::add_node(), libMesh::BoundaryInfo::add_side(), libMesh::BoundaryInfo::boundary_ids(), libMesh::Node::build(), libMesh::Elem::build(), libMesh::ParallelObject::comm(), libMesh::MeshBase::delete_remote_elements(), libMesh::Elem::dim(), libMesh::EDGE2, libMesh::EDGE3, libMesh::MeshBase::get_boundary_info(), libMesh::BoundaryInfo::get_edgeset_name_map(), libMesh::BoundaryInfo::get_nodeset_name_map(), libMesh::BoundaryInfo::get_side_boundary_ids(), libMesh::BoundaryInfo::get_sideset_name_map(), libMesh::MeshTools::Generation::QueryElemSubdomainIDBase::get_subdomain_for_layer(), libMesh::HEX27, libMesh::HEX8, libMesh::MeshBase::is_serial(), libMesh::libmesh_assert(), TIMPI::Communicator::max(), mesh, libMesh::MeshBase::n_elem(), libMesh::MeshBase::n_nodes(), libMesh::MeshBase::node_ptr(), libMesh::MeshBase::parallel_max_unique_id(), libMesh::MeshBase::prepare_for_use(), libMesh::PRISM18, libMesh::PRISM21, libMesh::PRISM6, libMesh::QUAD4, libMesh::QUAD9, libMesh::MeshBase::query_node_ptr(), libMesh::remote_elem, libMesh::MeshBase::reserve_elem(), libMesh::MeshBase::reserve_nodes(), libMesh::SECOND, libMesh::BoundaryInfo::set_edgeset_name_map(), libMesh::BoundaryInfo::set_nodeset_name_map(), libMesh::BoundaryInfo::set_sideset_name_map(), top_id, libMesh::TRI3, libMesh::TRI6, and libMesh::TRI7.

Referenced by main(), and MeshExtruderTest::testExtruder().

2257 {
2258  LOG_SCOPE("build_extrusion()", "MeshTools::Generation");
2259 
2260  if (!cross_section.n_elem())
2261  return;
2262 
2263  dof_id_type orig_elem = cross_section.n_elem();
2264  dof_id_type orig_nodes = cross_section.n_nodes();
2265 
2266 #ifdef LIBMESH_ENABLE_UNIQUE_ID
2267  unique_id_type orig_unique_ids = cross_section.parallel_max_unique_id();
2268 #endif
2269 
2270  unsigned int order = 1;
2271 
2272  BoundaryInfo & boundary_info = mesh.get_boundary_info();
2273  const BoundaryInfo & cross_section_boundary_info = cross_section.get_boundary_info();
2274 
2275  // Copy name maps from old to new boundary. We won't copy the whole
2276  // BoundaryInfo because that copies bc ids too, and we need to set
2277  // those more carefully.
2278  boundary_info.set_sideset_name_map() = cross_section_boundary_info.get_sideset_name_map();
2279  boundary_info.set_nodeset_name_map() = cross_section_boundary_info.get_nodeset_name_map();
2280  boundary_info.set_edgeset_name_map() = cross_section_boundary_info.get_edgeset_name_map();
2281 
2282  // If cross_section is distributed, so is its extrusion
2283  if (!cross_section.is_serial())
2284  mesh.delete_remote_elements();
2285 
2286  // We know a priori how many elements we'll need
2287  mesh.reserve_elem(nz*orig_elem);
2288 
2289  // For straightforward meshes we need one or two additional layers per
2290  // element.
2291  if (cross_section.elements_begin() != cross_section.elements_end() &&
2292  (*cross_section.elements_begin())->default_order() == SECOND)
2293  order = 2;
2294  mesh.comm().max(order);
2295 
2296  mesh.reserve_nodes((order*nz+1)*orig_nodes);
2297 
2298  // Container to catch the boundary IDs handed back by the BoundaryInfo object
2299  std::vector<boundary_id_type> ids_to_copy;
2300 
2301  for (const auto & node : cross_section.node_ptr_range())
2302  {
2303  for (unsigned int k=0; k != order*nz+1; ++k)
2304  {
2305  const dof_id_type new_node_id = node->id() + k * orig_nodes;
2306  Node * my_node = mesh.query_node_ptr(new_node_id);
2307  if (!my_node)
2308  {
2309  std::unique_ptr<Node> new_node = Node::build
2310  (*node + (extrusion_vector * k / nz / order),
2311  new_node_id);
2312  new_node->processor_id() = node->processor_id();
2313 
2314 #ifdef LIBMESH_ENABLE_UNIQUE_ID
2315  // Let's give the base of the extruded mesh the same
2316  // unique_ids as the source mesh, in case anyone finds that
2317  // a useful map to preserve.
2318  const unique_id_type uid = (k == 0) ?
2319  node->unique_id() :
2320  orig_unique_ids + (k-1)*(orig_nodes + orig_elem) + node->id();
2321 
2322  new_node->set_unique_id(uid);
2323 #endif
2324 
2325  cross_section_boundary_info.boundary_ids(node, ids_to_copy);
2326  boundary_info.add_node(new_node.get(), ids_to_copy);
2327 
2328  mesh.add_node(std::move(new_node));
2329  }
2330  }
2331  }
2332 
2333  const std::set<boundary_id_type> & side_ids =
2334  cross_section_boundary_info.get_side_boundary_ids();
2335 
2336  boundary_id_type next_side_id = side_ids.empty() ?
2337  0 : cast_int<boundary_id_type>(*side_ids.rbegin() + 1);
2338 
2339  // side_ids may not include ids from remote elements, in which case
2340  // some processors may have underestimated the next_side_id; let's
2341  // fix that.
2342  cross_section.comm().max(next_side_id);
2343 
2344  for (const auto & elem : cross_section.element_ptr_range())
2345  {
2346  const ElemType etype = elem->type();
2347 
2348  // build_extrusion currently only works on coarse meshes
2349  libmesh_assert (!elem->parent());
2350 
2351  for (unsigned int k=0; k != nz; ++k)
2352  {
2353  std::unique_ptr<Elem> new_elem;
2354  switch (etype)
2355  {
2356  case EDGE2:
2357  {
2358  new_elem = Elem::build(QUAD4);
2359  new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (k * orig_nodes));
2360  new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (k * orig_nodes));
2361  new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(1)->id() + ((k+1) * orig_nodes));
2362  new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(0)->id() + ((k+1) * orig_nodes));
2363 
2364  if (elem->neighbor_ptr(0) == remote_elem)
2365  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
2366  if (elem->neighbor_ptr(1) == remote_elem)
2367  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
2368 
2369  break;
2370  }
2371  case EDGE3:
2372  {
2373  new_elem = Elem::build(QUAD9);
2374  new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (2*k * orig_nodes));
2375  new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (2*k * orig_nodes));
2376  new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+2) * orig_nodes));
2377  new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+2) * orig_nodes));
2378  new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(2)->id() + (2*k * orig_nodes));
2379  new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+1) * orig_nodes));
2380  new_elem->set_node(6) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+2) * orig_nodes));
2381  new_elem->set_node(7) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+1) * orig_nodes));
2382  new_elem->set_node(8) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+1) * orig_nodes));
2383 
2384  if (elem->neighbor_ptr(0) == remote_elem)
2385  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
2386  if (elem->neighbor_ptr(1) == remote_elem)
2387  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
2388 
2389  break;
2390  }
2391  case TRI3:
2392  {
2393  new_elem = Elem::build(PRISM6);
2394  new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (k * orig_nodes));
2395  new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (k * orig_nodes));
2396  new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(2)->id() + (k * orig_nodes));
2397  new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(0)->id() + ((k+1) * orig_nodes));
2398  new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(1)->id() + ((k+1) * orig_nodes));
2399  new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(2)->id() + ((k+1) * orig_nodes));
2400 
2401  if (elem->neighbor_ptr(0) == remote_elem)
2402  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
2403  if (elem->neighbor_ptr(1) == remote_elem)
2404  new_elem->set_neighbor(2, const_cast<RemoteElem *>(remote_elem));
2405  if (elem->neighbor_ptr(2) == remote_elem)
2406  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
2407 
2408  break;
2409  }
2410  case TRI6:
2411  {
2412  new_elem = Elem::build(PRISM18);
2413  new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (2*k * orig_nodes));
2414  new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (2*k * orig_nodes));
2415  new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(2)->id() + (2*k * orig_nodes));
2416  new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+2) * orig_nodes));
2417  new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+2) * orig_nodes));
2418  new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+2) * orig_nodes));
2419  new_elem->set_node(6) = mesh.node_ptr(elem->node_ptr(3)->id() + (2*k * orig_nodes));
2420  new_elem->set_node(7) = mesh.node_ptr(elem->node_ptr(4)->id() + (2*k * orig_nodes));
2421  new_elem->set_node(8) = mesh.node_ptr(elem->node_ptr(5)->id() + (2*k * orig_nodes));
2422  new_elem->set_node(9) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+1) * orig_nodes));
2423  new_elem->set_node(10) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+1) * orig_nodes));
2424  new_elem->set_node(11) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+1) * orig_nodes));
2425  new_elem->set_node(12) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+2) * orig_nodes));
2426  new_elem->set_node(13) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+2) * orig_nodes));
2427  new_elem->set_node(14) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+2) * orig_nodes));
2428  new_elem->set_node(15) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+1) * orig_nodes));
2429  new_elem->set_node(16) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+1) * orig_nodes));
2430  new_elem->set_node(17) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+1) * orig_nodes));
2431 
2432  if (elem->neighbor_ptr(0) == remote_elem)
2433  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
2434  if (elem->neighbor_ptr(1) == remote_elem)
2435  new_elem->set_neighbor(2, const_cast<RemoteElem *>(remote_elem));
2436  if (elem->neighbor_ptr(2) == remote_elem)
2437  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
2438 
2439  break;
2440  }
2441  case TRI7:
2442  {
2443  new_elem = Elem::build(PRISM21);
2444  new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (2*k * orig_nodes));
2445  new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (2*k * orig_nodes));
2446  new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(2)->id() + (2*k * orig_nodes));
2447  new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+2) * orig_nodes));
2448  new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+2) * orig_nodes));
2449  new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+2) * orig_nodes));
2450  new_elem->set_node(6) = mesh.node_ptr(elem->node_ptr(3)->id() + (2*k * orig_nodes));
2451  new_elem->set_node(7) = mesh.node_ptr(elem->node_ptr(4)->id() + (2*k * orig_nodes));
2452  new_elem->set_node(8) = mesh.node_ptr(elem->node_ptr(5)->id() + (2*k * orig_nodes));
2453  new_elem->set_node(9) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+1) * orig_nodes));
2454  new_elem->set_node(10) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+1) * orig_nodes));
2455  new_elem->set_node(11) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+1) * orig_nodes));
2456  new_elem->set_node(12) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+2) * orig_nodes));
2457  new_elem->set_node(13) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+2) * orig_nodes));
2458  new_elem->set_node(14) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+2) * orig_nodes));
2459  new_elem->set_node(15) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+1) * orig_nodes));
2460  new_elem->set_node(16) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+1) * orig_nodes));
2461  new_elem->set_node(17) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+1) * orig_nodes));
2462 
2463  new_elem->set_node(18) = mesh.node_ptr(elem->node_ptr(6)->id() + (2*k * orig_nodes));
2464  new_elem->set_node(19) = mesh.node_ptr(elem->node_ptr(6)->id() + ((2*k+2) * orig_nodes));
2465  new_elem->set_node(20) = mesh.node_ptr(elem->node_ptr(6)->id() + ((2*k+1) * orig_nodes));
2466 
2467  if (elem->neighbor_ptr(0) == remote_elem)
2468  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
2469  if (elem->neighbor_ptr(1) == remote_elem)
2470  new_elem->set_neighbor(2, const_cast<RemoteElem *>(remote_elem));
2471  if (elem->neighbor_ptr(2) == remote_elem)
2472  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
2473 
2474  break;
2475  }
2476  case QUAD4:
2477  {
2478  new_elem = Elem::build(HEX8);
2479  new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (k * orig_nodes));
2480  new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (k * orig_nodes));
2481  new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(2)->id() + (k * orig_nodes));
2482  new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(3)->id() + (k * orig_nodes));
2483  new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(0)->id() + ((k+1) * orig_nodes));
2484  new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(1)->id() + ((k+1) * orig_nodes));
2485  new_elem->set_node(6) = mesh.node_ptr(elem->node_ptr(2)->id() + ((k+1) * orig_nodes));
2486  new_elem->set_node(7) = mesh.node_ptr(elem->node_ptr(3)->id() + ((k+1) * orig_nodes));
2487 
2488  if (elem->neighbor_ptr(0) == remote_elem)
2489  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
2490  if (elem->neighbor_ptr(1) == remote_elem)
2491  new_elem->set_neighbor(2, const_cast<RemoteElem *>(remote_elem));
2492  if (elem->neighbor_ptr(2) == remote_elem)
2493  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
2494  if (elem->neighbor_ptr(3) == remote_elem)
2495  new_elem->set_neighbor(4, const_cast<RemoteElem *>(remote_elem));
2496 
2497  break;
2498  }
2499  case QUAD9:
2500  {
2501  new_elem = Elem::build(HEX27);
2502  new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (2*k * orig_nodes));
2503  new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (2*k * orig_nodes));
2504  new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(2)->id() + (2*k * orig_nodes));
2505  new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(3)->id() + (2*k * orig_nodes));
2506  new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+2) * orig_nodes));
2507  new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+2) * orig_nodes));
2508  new_elem->set_node(6) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+2) * orig_nodes));
2509  new_elem->set_node(7) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+2) * orig_nodes));
2510  new_elem->set_node(8) = mesh.node_ptr(elem->node_ptr(4)->id() + (2*k * orig_nodes));
2511  new_elem->set_node(9) = mesh.node_ptr(elem->node_ptr(5)->id() + (2*k * orig_nodes));
2512  new_elem->set_node(10) = mesh.node_ptr(elem->node_ptr(6)->id() + (2*k * orig_nodes));
2513  new_elem->set_node(11) = mesh.node_ptr(elem->node_ptr(7)->id() + (2*k * orig_nodes));
2514  new_elem->set_node(12) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+1) * orig_nodes));
2515  new_elem->set_node(13) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+1) * orig_nodes));
2516  new_elem->set_node(14) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+1) * orig_nodes));
2517  new_elem->set_node(15) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+1) * orig_nodes));
2518  new_elem->set_node(16) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+2) * orig_nodes));
2519  new_elem->set_node(17) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+2) * orig_nodes));
2520  new_elem->set_node(18) = mesh.node_ptr(elem->node_ptr(6)->id() + ((2*k+2) * orig_nodes));
2521  new_elem->set_node(19) = mesh.node_ptr(elem->node_ptr(7)->id() + ((2*k+2) * orig_nodes));
2522  new_elem->set_node(20) = mesh.node_ptr(elem->node_ptr(8)->id() + (2*k * orig_nodes));
2523  new_elem->set_node(21) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+1) * orig_nodes));
2524  new_elem->set_node(22) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+1) * orig_nodes));
2525  new_elem->set_node(23) = mesh.node_ptr(elem->node_ptr(6)->id() + ((2*k+1) * orig_nodes));
2526  new_elem->set_node(24) = mesh.node_ptr(elem->node_ptr(7)->id() + ((2*k+1) * orig_nodes));
2527  new_elem->set_node(25) = mesh.node_ptr(elem->node_ptr(8)->id() + ((2*k+2) * orig_nodes));
2528  new_elem->set_node(26) = mesh.node_ptr(elem->node_ptr(8)->id() + ((2*k+1) * orig_nodes));
2529 
2530  if (elem->neighbor_ptr(0) == remote_elem)
2531  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
2532  if (elem->neighbor_ptr(1) == remote_elem)
2533  new_elem->set_neighbor(2, const_cast<RemoteElem *>(remote_elem));
2534  if (elem->neighbor_ptr(2) == remote_elem)
2535  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
2536  if (elem->neighbor_ptr(3) == remote_elem)
2537  new_elem->set_neighbor(4, const_cast<RemoteElem *>(remote_elem));
2538 
2539  break;
2540  }
2541  default:
2542  {
2543  libmesh_not_implemented();
2544  break;
2545  }
2546  }
2547 
2548  new_elem->set_id(elem->id() + (k * orig_elem));
2549  new_elem->processor_id() = elem->processor_id();
2550 
2551 #ifdef LIBMESH_ENABLE_UNIQUE_ID
2552  // Let's give the base of the extruded mesh the same
2553  // unique_ids as the source mesh, in case anyone finds that
2554  // a useful map to preserve.
2555  const unique_id_type uid = (k == 0) ?
2556  elem->unique_id() :
2557  orig_unique_ids + (k-1)*(orig_nodes + orig_elem) + orig_nodes + elem->id();
2558 
2559  new_elem->set_unique_id(uid);
2560 #endif
2561 
2562  if (!elem_subdomain)
2563  // maintain the subdomain_id
2564  new_elem->subdomain_id() = elem->subdomain_id();
2565  else
2566  // Allow the user to choose new subdomain_ids
2567  new_elem->subdomain_id() = elem_subdomain->get_subdomain_for_layer(elem, k);
2568 
2569  Elem * added_elem = mesh.add_elem(std::move(new_elem));
2570 
2571  // Copy any old boundary ids on all sides
2572  for (auto s : elem->side_index_range())
2573  {
2574  cross_section_boundary_info.boundary_ids(elem, s, ids_to_copy);
2575 
2576  if (added_elem->dim() == 3)
2577  {
2578  // For 2D->3D extrusion, we give the boundary IDs
2579  // for side s on the old element to side s+1 on the
2580  // new element. This is just a happy coincidence as
2581  // far as I can tell...
2582  boundary_info.add_side(added_elem,
2583  cast_int<unsigned short>(s+1),
2584  ids_to_copy);
2585  }
2586  else
2587  {
2588  // For 1D->2D extrusion, the boundary IDs map as:
2589  // Old elem -> New elem
2590  // 0 -> 3
2591  // 1 -> 1
2592  libmesh_assert_less(s, 2);
2593  const unsigned short sidemap[2] = {3, 1};
2594  boundary_info.add_side(added_elem, sidemap[s], ids_to_copy);
2595  }
2596  }
2597 
2598  // Give new boundary ids to bottom and top
2599  if (k == 0)
2600  boundary_info.add_side(added_elem, 0, next_side_id);
2601  if (k == nz-1)
2602  {
2603  // For 2D->3D extrusion, the "top" ID is 1+the original
2604  // element's number of sides. For 1D->2D extrusion, the
2605  // "top" ID is side 2.
2606  const unsigned short top_id = added_elem->dim() == 3 ?
2607  cast_int<unsigned short>(elem->n_sides()+1) : 2;
2608  boundary_info.add_side
2609  (added_elem, top_id,
2610  cast_int<boundary_id_type>(next_side_id+1));
2611  }
2612  }
2613  }
2614 
2615  // Done building the mesh. Now prepare it for use.
2616  mesh.prepare_for_use();
2617 }
libMesh::ElemType
ElemType
Defines an enum for geometric element types.
Definition: enum_elem_type.h:33
mesh
MeshBase & mesh
Definition: mesh_communication.C:1308
libMesh::boundary_id_type
int8_t boundary_id_type
Definition: id_types.h:51
libMesh::libmesh_assert
libmesh_assert(ctx)
top_id
const boundary_id_type top_id
Definition: periodic_bc_test.C:156
libMesh::unique_id_type
uint8_t unique_id_type
Definition: id_types.h:86
libMesh::dof_id_type
uint8_t dof_id_type
Definition: id_types.h:67
libMesh::remote_elem
const RemoteElem * remote_elem
Definition: remote_elem.C:54

◆ build_line()

void libMesh::MeshTools::Generation::build_line ( UnstructuredMesh &  mesh,
const unsigned int  nx,
const Real  xmin = 0.,
const Real  xmax = 1.,
const ElemType  type = INVALID_ELEM,
const bool  gauss_lobatto_grid = false 
)

A specialized build_cube() for 1D meshes.

Boundary ids are set to be equal to the side indexing on a master edge

Definition at line 1663 of file mesh_generation.C.

References build_cube(), and mesh.

Referenced by Biharmonic::Biharmonic(), NodalNeighborsTest::do_test(), main(), DualShapeTest::setUp(), SystemsTest::simpleSetup(), SystemsTest::test100KVariables(), MeshSpatialDimensionTest::test1D(), ConstraintOperatorTest::test1DCoarseningNewNodes(), ConstraintOperatorTest::test1DCoarseningOperator(), MeshGenerationTest::testBuildLine(), VolumeTest::testEdge3Volume(), MeshSubdomainIDTest::testMultiple(), BoundaryInfoTest::testNameCopying(), EquationSystemsTest::testPostInitAddElem(), SystemsTest::testProjectLine(), SystemsTest::testProjectMatrix1D(), EquationSystemsTest::testReinitWithNodeElem(), DistributedMeshTest::testRemoteElemError(), and EquationSystemsTest::testSelectivePRefine().

1668 {
1669  // This method only makes sense in 1D!
1670  // But we now just turn a non-1D mesh into a 1D mesh
1671  //libmesh_assert_equal_to (mesh.mesh_dimension(), 1);
1672 
1673  build_cube(mesh,
1674  nx, 0, 0,
1675  xmin, xmax,
1676  0., 0.,
1677  0., 0.,
1678  type,
1679  gauss_lobatto_grid);
1680 }
mesh
MeshBase & mesh
Definition: mesh_communication.C:1308
libMesh::MeshTools::Generation::build_cube
void build_cube(UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
Builds a (elements) cube.
Definition: mesh_generation.C:311

◆ build_point()

void libMesh::MeshTools::Generation::build_point ( UnstructuredMesh &  mesh,
const ElemType  type = INVALID_ELEM,
const bool  gauss_lobatto_grid = false 
)

A specialized build_cube() for 0D meshes.

The resulting mesh is a single NodeElem suitable for ODE tests

Definition at line 1645 of file mesh_generation.C.

References build_cube(), and mesh.

Referenced by TimeSolverTestImplementation< NewmarkSolver >::run_test_with_exact_soln(), EquationSystemsTest::testInit(), EquationSystemsTest::testPostInitAddRealSystem(), and EquationSystemsTest::testPostInitAddSystem().

1648 {
1649  // This method only makes sense in 0D!
1650  // But we now just turn a non-0D mesh into a 0D mesh
1651  //libmesh_assert_equal_to (mesh.mesh_dimension(), 1);
1652 
1653  build_cube(mesh,
1654  0, 0, 0,
1655  0., 0.,
1656  0., 0.,
1657  0., 0.,
1658  type,
1659  gauss_lobatto_grid);
1660 }
mesh
MeshBase & mesh
Definition: mesh_communication.C:1308
libMesh::MeshTools::Generation::build_cube
void build_cube(UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
Builds a (elements) cube.
Definition: mesh_generation.C:311

◆ build_sphere()

void libMesh::MeshTools::Generation::build_sphere ( UnstructuredMesh &  mesh,
const Real  rad = 1,
const unsigned int  nr = 2,
const ElemType  type = INVALID_ELEM,
const unsigned int  n_smooth = 2,
const bool  flat = true 
)

Meshes a spherical or mapped-spherical domain.

Definition at line 1715 of file mesh_generation.C.

Referenced by main(), MeshGenerationTest::testBuildSphere(), InfFERadialTest::testInfQuants_numericDeriv(), and MeshTriangulationTest::testTriangulatorRoundHole().

1721 {
1722  libmesh_error_msg("Building a circle/sphere only works with AMR.");
1723 }

◆ build_square()

void libMesh::MeshTools::Generation::build_square ( UnstructuredMesh &  mesh,
const unsigned int  nx,
const unsigned int  ny,
const Real  xmin = 0.,
const Real  xmax = 1.,
const Real  ymin = 0.,
const Real  ymax = 1.,
const ElemType  type = INVALID_ELEM,
const bool  gauss_lobatto_grid = false 
)

A specialized build_cube() for 2D meshes.

Boundary ids are set to be equal to the side indexing on a master quad

Definition at line 1684 of file mesh_generation.C.

References build_cube(), and mesh.

Referenced by AllSecondOrderTest::allCompleteOrder(), AllSecondOrderTest::allCompleteOrderDoNothing(), AllSecondOrderTest::allCompleteOrderMixed(), AllSecondOrderTest::allCompleteOrderMixedFixing(), AllSecondOrderTest::allCompleteOrderRange(), AllSecondOrderTest::allSecondOrder(), AllSecondOrderTest::allSecondOrderDoNothing(), AllSecondOrderTest::allSecondOrderMixed(), AllSecondOrderTest::allSecondOrderMixedFixing(), AllSecondOrderTest::allSecondOrderRange(), Biharmonic::Biharmonic(), BoundaryMeshTest::build_mesh(), main(), WriteVecAndScalar::setupTests(), MultiEvaluablePredTest::test(), MeshSpatialDimensionTest::test2D(), SystemsTest::test2DProjectVectorFE(), DistortTest::test_helper_2D(), AllTriTest::test_helper_2D(), MeshFunctionTest::test_subdomain_id_sets(), DofMapTest::testBadElemFECombo(), EquationSystemsTest::testBadVarNames(), BoundaryInfoTest::testBoundaryOnChildrenBoundaryIDs(), BoundaryInfoTest::testBoundaryOnChildrenBoundarySides(), BoundaryInfoTest::testBoundaryOnChildrenElementsRefineCoarsen(), BoundaryInfoTest::testBoundaryOnChildrenErrors(), MeshGenerationTest::testBuildSquare(), DofMapTest::testConstraintLoopDetection(), MeshInputTest::testCopyElementSolutionImpl(), MeshInputTest::testCopyElementVectorImpl(), MeshInputTest::testCopyNodalSolutionImpl(), MeshDeletionsTest::testDeleteElem(), MeshInputTest::testExodusReadHeader(), MeshExtruderTest::testExtruder(), MixedOrderTest::testFindNeighbors(), MappedSubdomainPartitionerTest::testMappedSubdomainPartitioner(), BoundaryInfoTest::testMesh(), MeshBaseTest::testMeshBaseVerifyIsPrepared(), MeshAssignTest::testMeshMoveAssign(), MeshInputTest::testNemesisReadImpl(), PointLocatorTest::testPlanar(), SystemsTest::testProjectMatrix2D(), SystemsTest::testProjectSquare(), VolumeTest::testQuad4TrueCentroid(), EquationSystemsTest::testRefineThenReinitPreserveFlags(), BoundaryInfoTest::testRenumber(), EquationSystemsTest::testRepartitionThenReinit(), MeshInputTest::testSingleElementImpl(), CheckpointIOTest::testSplitter(), MixedOrderTest::testStitch(), MeshTriangulationTest::testTriangulatorMeshedHoles(), MeshInputTest::testVTKPreserveElemIds(), MeshInputTest::testVTKPreserveSubdomainIds(), WriteNodesetData::testWriteImpl(), WriteSidesetData::testWriteImpl(), and WriteElemsetData::testWriteImpl().

1691 {
1692  // This method only makes sense in 2D!
1693  // But we now just turn a non-2D mesh into a 2D mesh
1694  //libmesh_assert_equal_to (mesh.mesh_dimension(), 2);
1695 
1696  // Call the build_cube() member to actually do the work for us.
1697  build_cube (mesh,
1698  nx, ny, 0,
1699  xmin, xmax,
1700  ymin, ymax,
1701  0., 0.,
1702  type,
1703  gauss_lobatto_grid);
1704 }
mesh
MeshBase & mesh
Definition: mesh_communication.C:1308
libMesh::MeshTools::Generation::build_cube
void build_cube(UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
Builds a (elements) cube.
Definition: mesh_generation.C:311
Generated on Thu Apr 18 2024 11:05:28 for libMesh by   doxygen 1.8.14