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 298 of file mesh_generation.C.

{
  START_LOG("build_cube()", "MeshTools::Generation");
 
  // Declare that we are using the indexing utility routine
  // in the "Private" part of our current namespace.  If this doesn't
  // work in GCC 2.95.3 we can either remove it or stop supporting
  // 2.95.3 altogether.
  // Changing this to import the whole namespace... just importing idx
  // causes an internal compiler error for Intel Compiler 11.0 on Linux
  // in debug mode.
  using namespace MeshTools::Generation::Private;
 
  // Clear the mesh and start from scratch
  mesh.clear();
 
  BoundaryInfo & boundary_info = mesh.get_boundary_info();
 
  if (nz != 0)
    {
      mesh.set_mesh_dimension(3);
      mesh.set_spatial_dimension(3);
    }
  else if (ny != 0)
    {
      mesh.set_mesh_dimension(2);
      mesh.set_spatial_dimension(2);
    }
  else if (nx != 0)
    {
      mesh.set_mesh_dimension(1);
      mesh.set_spatial_dimension(1);
    }
  else
    {
      // Will we get here?
      mesh.set_mesh_dimension(0);
      mesh.set_spatial_dimension(0);
    }
 
  switch (mesh.mesh_dimension())
    {
      //---------------------------------------------------------------------
      // Build a 0D point
    case 0:
      {
        libmesh_assert_equal_to (nx, 0);
        libmesh_assert_equal_to (ny, 0);
        libmesh_assert_equal_to (nz, 0);
 
        libmesh_assert (type == INVALID_ELEM || type == NODEELEM);
 
        // Build one nodal element for the mesh
        mesh.add_point (Point(0, 0, 0), 0);
        Elem * elem = mesh.add_elem (new NodeElem);
        elem->set_node(0) = mesh.node_ptr(0);
 
        break;
      }
 
 
 
      //---------------------------------------------------------------------
      // Build a 1D line
    case 1:
      {
        libmesh_assert_not_equal_to (nx, 0);
        libmesh_assert_equal_to (ny, 0);
        libmesh_assert_equal_to (nz, 0);
        libmesh_assert_less (xmin, xmax);
 
        // Reserve elements
        switch (type)
          {
          case INVALID_ELEM:
          case EDGE2:
          case EDGE3:
          case EDGE4:
            {
              mesh.reserve_elem (nx);
              break;
            }
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 1D element type.");
          }
 
        // Reserve nodes
        switch (type)
          {
          case INVALID_ELEM:
          case EDGE2:
            {
              mesh.reserve_nodes(nx+1);
              break;
            }
 
          case EDGE3:
            {
              mesh.reserve_nodes(2*nx+1);
              break;
            }
 
          case EDGE4:
            {
              mesh.reserve_nodes(3*nx+1);
              break;
            }
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 1D element type.");
          }
 
 
        // Build the nodes, depends on whether we're using linears,
        // quadratics or cubics and whether using uniform grid or Gauss-Lobatto
        unsigned int node_id = 0;
        switch(type)
          {
          case INVALID_ELEM:
          case EDGE2:
            {
              for (unsigned int i=0; i<=nx; i++)
                mesh.add_point (Point(static_cast<Real>(i)/nx, 0, 0), node_id++);
 
              break;
            }
 
          case EDGE3:
            {
              for (unsigned int i=0; i<=2*nx; i++)
                mesh.add_point (Point(static_cast<Real>(i)/(2*nx), 0, 0), node_id++);
              break;
            }
 
          case EDGE4:
            {
              for (unsigned int i=0; i<=3*nx; i++)
                mesh.add_point (Point(static_cast<Real>(i)/(3*nx), 0, 0), node_id++);
 
              break;
            }
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 1D element type.");
 
          }
 
        // Build the elements of the mesh
        switch(type)
          {
          case INVALID_ELEM:
          case EDGE2:
            {
              for (unsigned int i=0; i<nx; i++)
                {
                  Elem * elem = new Edge2;
                  elem->set_id(i);
                  elem = mesh.add_elem (elem);
                  elem->set_node(0) = mesh.node_ptr(i);
                  elem->set_node(1) = mesh.node_ptr(i+1);
 
                  if (i == 0)
                    boundary_info.add_side(elem, 0, 0);
 
                  if (i == (nx-1))
                    boundary_info.add_side(elem, 1, 1);
                }
              break;
            }
 
          case EDGE3:
            {
              for (unsigned int i=0; i<nx; i++)
                {
                  Elem * elem = new Edge3;
                  elem->set_id(i);
                  elem = mesh.add_elem (elem);
                  elem->set_node(0) = mesh.node_ptr(2*i);
                  elem->set_node(2) = mesh.node_ptr(2*i+1);
                  elem->set_node(1) = mesh.node_ptr(2*i+2);
 
                  if (i == 0)
                    boundary_info.add_side(elem, 0, 0);
 
                  if (i == (nx-1))
                    boundary_info.add_side(elem, 1, 1);
                }
              break;
            }
 
          case EDGE4:
            {
              for (unsigned int i=0; i<nx; i++)
                {
                  Elem * elem = new Edge4;
                  elem->set_id(i);
                  elem = mesh.add_elem (elem);
                  elem->set_node(0) = mesh.node_ptr(3*i);
                  elem->set_node(2) = mesh.node_ptr(3*i+1);
                  elem->set_node(3) = mesh.node_ptr(3*i+2);
                  elem->set_node(1) = mesh.node_ptr(3*i+3);
 
                  if (i == 0)
                    boundary_info.add_side(elem, 0, 0);
 
                  if (i == (nx-1))
                    boundary_info.add_side(elem, 1, 1);
                }
              break;
            }
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 1D element type.");
          }
 
        // Move the nodes to their final locations.
        if (gauss_lobatto_grid)
          {
            GaussLobattoRedistributionFunction func(nx, xmin, xmax);
            MeshTools::Modification::redistribute(mesh, func);
          }
        else // !gauss_lobatto_grid
          {
            for (Node * node : mesh.node_ptr_range())
              (*node)(0) = (*node)(0)*(xmax-xmin) + xmin;
          }
 
        // Add sideset names to boundary info
        boundary_info.sideset_name(0) = "left";
        boundary_info.sideset_name(1) = "right";
 
        // Add nodeset names to boundary info
        boundary_info.nodeset_name(0) = "left";
        boundary_info.nodeset_name(1) = "right";
 
        break;
      }
 
 
 
 
 
 
 
 
 
 
      //---------------------------------------------------------------------
      // Build a 2D quadrilateral
    case 2:
      {
        libmesh_assert_not_equal_to (nx, 0);
        libmesh_assert_not_equal_to (ny, 0);
        libmesh_assert_equal_to (nz, 0);
        libmesh_assert_less (xmin, xmax);
        libmesh_assert_less (ymin, ymax);
 
        // Reserve elements.  The TRI3 and TRI6 meshes
        // have twice as many elements...
        switch (type)
          {
          case INVALID_ELEM:
          case QUAD4:
          case QUAD8:
          case QUAD9:
            {
              mesh.reserve_elem (nx*ny);
              break;
            }
 
          case TRI3:
          case TRI6:
            {
              mesh.reserve_elem (2*nx*ny);
              break;
            }
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 2D element type.");
          }
 
 
 
        // Reserve nodes.  The quadratic element types
        // need to reserve more nodes than the linear types.
        switch (type)
          {
          case INVALID_ELEM:
          case QUAD4:
          case TRI3:
            {
              mesh.reserve_nodes( (nx+1)*(ny+1) );
              break;
            }
 
          case QUAD8:
          case QUAD9:
          case TRI6:
            {
              mesh.reserve_nodes( (2*nx+1)*(2*ny+1) );
              break;
            }
 
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 2D element type.");
          }
 
 
 
        // Build the nodes. Depends on whether you are using a linear
        // or quadratic element, and whether you are using a uniform
        // grid or the Gauss-Lobatto grid points.
        unsigned int node_id = 0;
        switch (type)
          {
          case INVALID_ELEM:
          case QUAD4:
          case TRI3:
            {
              for (unsigned int j=0; j<=ny; j++)
                for (unsigned int i=0; i<=nx; i++)
                  mesh.add_point (Point(static_cast<Real>(i)/static_cast<Real>(nx),
                                        static_cast<Real>(j)/static_cast<Real>(ny),
                                        0.), node_id++);
 
              break;
            }
 
          case QUAD8:
          case QUAD9:
          case TRI6:
            {
              for (unsigned int j=0; j<=(2*ny); j++)
                for (unsigned int i=0; i<=(2*nx); i++)
                  mesh.add_point (Point(static_cast<Real>(i)/static_cast<Real>(2*nx),
                                        static_cast<Real>(j)/static_cast<Real>(2*ny),
                                        0), node_id++);
 
              break;
            }
 
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 2D element type.");
          }
 
 
 
 
 
 
        // Build the elements.  Each one is a bit different.
        unsigned int elem_id = 0;
        switch (type)
          {
 
          case INVALID_ELEM:
          case QUAD4:
            {
              for (unsigned int j=0; j<ny; j++)
                for (unsigned int i=0; i<nx; i++)
                  {
                    Elem * elem = new Quad4;
                    elem->set_id(elem_id++);
                    elem = mesh.add_elem (elem);
 
                    elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j)    );
                    elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+1,j)  );
                    elem->set_node(2) = mesh.node_ptr(idx(type,nx,i+1,j+1));
                    elem->set_node(3) = mesh.node_ptr(idx(type,nx,i,j+1)  );
 
                    if (j == 0)
                      boundary_info.add_side(elem, 0, 0);
 
                    if (j == (ny-1))
                      boundary_info.add_side(elem, 2, 2);
 
                    if (i == 0)
                      boundary_info.add_side(elem, 3, 3);
 
                    if (i == (nx-1))
                      boundary_info.add_side(elem, 1, 1);
                  }
              break;
            }
 
 
          case TRI3:
            {
              for (unsigned int j=0; j<ny; j++)
                for (unsigned int i=0; i<nx; i++)
                  {
                    // Add first Tri3
                    Elem * elem = new Tri3;
                    elem->set_id(elem_id++);
                    elem = mesh.add_elem (elem);
 
                    elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j)    );
                    elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+1,j)  );
                    elem->set_node(2) = mesh.node_ptr(idx(type,nx,i+1,j+1));
 
                    if (j == 0)
                      boundary_info.add_side(elem, 0, 0);
 
                    if (i == (nx-1))
                      boundary_info.add_side(elem, 1, 1);
 
                    // Add second Tri3
                    elem = new Tri3;
                    elem->set_id(elem_id++);
                    elem = mesh.add_elem (elem);
 
                    elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j)    );
                    elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+1,j+1));
                    elem->set_node(2) = mesh.node_ptr(idx(type,nx,i,j+1)  );
 
                    if (j == (ny-1))
                      boundary_info.add_side(elem, 1, 2);
 
                    if (i == 0)
                      boundary_info.add_side(elem, 2, 3);
                  }
              break;
            }
 
 
 
          case QUAD8:
          case QUAD9:
            {
              for (unsigned int j=0; j<(2*ny); j += 2)
                for (unsigned int i=0; i<(2*nx); i += 2)
                  {
                    Elem * elem = (type == QUAD8) ?
                      static_cast<Elem *>(new Quad8) :
                      static_cast<Elem *>(new Quad9);
                    elem->set_id(elem_id++);
                    elem = mesh.add_elem (elem);
 
                    elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j)    );
                    elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+2,j)  );
                    elem->set_node(2) = mesh.node_ptr(idx(type,nx,i+2,j+2));
                    elem->set_node(3) = mesh.node_ptr(idx(type,nx,i,j+2)  );
                    elem->set_node(4) = mesh.node_ptr(idx(type,nx,i+1,j)  );
                    elem->set_node(5) = mesh.node_ptr(idx(type,nx,i+2,j+1));
                    elem->set_node(6) = mesh.node_ptr(idx(type,nx,i+1,j+2));
                    elem->set_node(7) = mesh.node_ptr(idx(type,nx,i,j+1)  );
                    if (type == QUAD9)
                      elem->set_node(8) = mesh.node_ptr(idx(type,nx,i+1,j+1));
 
 
                    if (j == 0)
                      boundary_info.add_side(elem, 0, 0);
 
                    if (j == 2*(ny-1))
                      boundary_info.add_side(elem, 2, 2);
 
                    if (i == 0)
                      boundary_info.add_side(elem, 3, 3);
 
                    if (i == 2*(nx-1))
                      boundary_info.add_side(elem, 1, 1);
                  }
              break;
            }
 
 
          case TRI6:
            {
              for (unsigned int j=0; j<(2*ny); j += 2)
                for (unsigned int i=0; i<(2*nx); i += 2)
                  {
                    // Add first Tri6
                    Elem * elem = new Tri6;
                    elem->set_id(elem_id++);
                    elem = mesh.add_elem (elem);
 
                    elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j)    );
                    elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+2,j)  );
                    elem->set_node(2) = mesh.node_ptr(idx(type,nx,i+2,j+2));
                    elem->set_node(3) = mesh.node_ptr(idx(type,nx,i+1,j)  );
                    elem->set_node(4) = mesh.node_ptr(idx(type,nx,i+2,j+1));
                    elem->set_node(5) = mesh.node_ptr(idx(type,nx,i+1,j+1));
 
                    if (j == 0)
                      boundary_info.add_side(elem, 0, 0);
 
                    if (i == 2*(nx-1))
                      boundary_info.add_side(elem, 1, 1);
 
                    // Add second Tri6
                    elem = new Tri6;
                    elem->set_id(elem_id++);
                    elem = mesh.add_elem (elem);
 
                    elem->set_node(0) = mesh.node_ptr(idx(type,nx,i,j)    );
                    elem->set_node(1) = mesh.node_ptr(idx(type,nx,i+2,j+2));
                    elem->set_node(2) = mesh.node_ptr(idx(type,nx,i,j+2)  );
                    elem->set_node(3) = mesh.node_ptr(idx(type,nx,i+1,j+1));
                    elem->set_node(4) = mesh.node_ptr(idx(type,nx,i+1,j+2));
                    elem->set_node(5) = mesh.node_ptr(idx(type,nx,i,j+1)  );
 
                    if (j == 2*(ny-1))
                      boundary_info.add_side(elem, 1, 2);
 
                    if (i == 0)
                      boundary_info.add_side(elem, 2, 3);
 
                  }
              break;
            };
 
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 2D element type.");
          }
 
 
 
 
        // Scale the nodal positions
        if (gauss_lobatto_grid)
          {
            GaussLobattoRedistributionFunction func(nx, xmin, xmax,
                                                    ny, ymin, ymax);
            MeshTools::Modification::redistribute(mesh, func);
          }
        else // !gauss_lobatto_grid
          {
            for (Node * node : mesh.node_ptr_range())
              {
                (*node)(0) = ((*node)(0))*(xmax-xmin) + xmin;
                (*node)(1) = ((*node)(1))*(ymax-ymin) + ymin;
              }
          }
 
        // Add sideset names to boundary info
        boundary_info.sideset_name(0) = "bottom";
        boundary_info.sideset_name(1) = "right";
        boundary_info.sideset_name(2) = "top";
        boundary_info.sideset_name(3) = "left";
 
        // Add nodeset names to boundary info
        boundary_info.nodeset_name(0) = "bottom";
        boundary_info.nodeset_name(1) = "right";
        boundary_info.nodeset_name(2) = "top";
        boundary_info.nodeset_name(3) = "left";
 
        break;
      }
 
 
 
 
 
 
 
 
 
 
 
      //---------------------------------------------------------------------
      // Build a 3D mesh using hexes, tets, prisms, or pyramids.
    case 3:
      {
        libmesh_assert_not_equal_to (nx, 0);
        libmesh_assert_not_equal_to (ny, 0);
        libmesh_assert_not_equal_to (nz, 0);
        libmesh_assert_less (xmin, xmax);
        libmesh_assert_less (ymin, ymax);
        libmesh_assert_less (zmin, zmax);
 
 
        // Reserve elements.  Meshes with prismatic elements require
        // twice as many elements.
        switch (type)
          {
          case INVALID_ELEM:
          case HEX8:
          case HEX20:
          case HEX27:
          case TET4:  // TET4's are created from an initial HEX27 discretization
          case TET10: // TET10's are created from an initial HEX27 discretization
          case PYRAMID5: // PYRAMIDs are created from an initial HEX27 discretization
          case PYRAMID13:
          case PYRAMID14:
            {
              mesh.reserve_elem(nx*ny*nz);
              break;
            }
 
          case PRISM6:
          case PRISM15:
          case PRISM18:
            {
              mesh.reserve_elem(2*nx*ny*nz);
              break;
            }
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 3D element type.");
          }
 
 
 
 
 
        // Reserve nodes.  Quadratic elements need twice as many nodes as linear elements.
        switch (type)
          {
          case INVALID_ELEM:
          case HEX8:
          case PRISM6:
            {
              mesh.reserve_nodes( (nx+1)*(ny+1)*(nz+1) );
              break;
            }
 
          case HEX20:
          case HEX27:
          case TET4: // TET4's are created from an initial HEX27 discretization
          case TET10: // TET10's are created from an initial HEX27 discretization
          case PYRAMID5: // PYRAMIDs are created from an initial HEX27 discretization
          case PYRAMID13:
          case PYRAMID14:
          case PRISM15:
          case PRISM18:
            {
              // FYI: The resulting TET4 mesh will have exactly
              // 5*(nx*ny*nz) + 2*(nx*ny + nx*nz + ny*nz) + (nx+ny+nz) + 1
              // nodes once the additional mid-edge nodes for the HEX27 discretization
              // have been deleted.
              mesh.reserve_nodes( (2*nx+1)*(2*ny+1)*(2*nz+1) );
              break;
            }
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 3D element type.");
          }
 
 
 
 
        // Build the nodes.
        unsigned int node_id = 0;
        switch (type)
          {
          case INVALID_ELEM:
          case HEX8:
          case PRISM6:
            {
              for (unsigned int k=0; k<=nz; k++)
                for (unsigned int j=0; j<=ny; j++)
                  for (unsigned int i=0; i<=nx; i++)
                    mesh.add_point(Point(static_cast<Real>(i)/static_cast<Real>(nx),
                                         static_cast<Real>(j)/static_cast<Real>(ny),
                                         static_cast<Real>(k)/static_cast<Real>(nz)), node_id++);
 
              break;
            }
 
          case HEX20:
          case HEX27:
          case TET4: // TET4's are created from an initial HEX27 discretization
          case TET10: // TET10's are created from an initial HEX27 discretization
          case PYRAMID5: // PYRAMIDs are created from an initial HEX27 discretization
          case PYRAMID13:
          case PYRAMID14:
          case PRISM15:
          case PRISM18:
            {
              for (unsigned int k=0; k<=(2*nz); k++)
                for (unsigned int j=0; j<=(2*ny); j++)
                  for (unsigned int i=0; i<=(2*nx); i++)
                    mesh.add_point(Point(static_cast<Real>(i)/static_cast<Real>(2*nx),
                                         static_cast<Real>(j)/static_cast<Real>(2*ny),
                                         static_cast<Real>(k)/static_cast<Real>(2*nz)), node_id++);
 
              break;
            }
 
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 3D element type.");
          }
 
 
 
 
        // Build the elements.
        unsigned int elem_id = 0;
        switch (type)
          {
          case INVALID_ELEM:
          case HEX8:
            {
              for (unsigned int k=0; k<nz; k++)
                for (unsigned int j=0; j<ny; j++)
                  for (unsigned int i=0; i<nx; i++)
                    {
                      Elem * elem = new Hex8;
                      elem->set_id(elem_id++);
                      elem = mesh.add_elem (elem);
 
                      elem->set_node(0) = mesh.node_ptr(idx(type,nx,ny,i,j,k)      );
                      elem->set_node(1) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k)    );
                      elem->set_node(2) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k)  );
                      elem->set_node(3) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k)    );
                      elem->set_node(4) = mesh.node_ptr(idx(type,nx,ny,i,j,k+1)    );
                      elem->set_node(5) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k+1)  );
                      elem->set_node(6) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
                      elem->set_node(7) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k+1)  );
 
                      if (k == 0)
                        boundary_info.add_side(elem, 0, 0);
 
                      if (k == (nz-1))
                        boundary_info.add_side(elem, 5, 5);
 
                      if (j == 0)
                        boundary_info.add_side(elem, 1, 1);
 
                      if (j == (ny-1))
                        boundary_info.add_side(elem, 3, 3);
 
                      if (i == 0)
                        boundary_info.add_side(elem, 4, 4);
 
                      if (i == (nx-1))
                        boundary_info.add_side(elem, 2, 2);
                    }
              break;
            }
 
 
 
 
          case PRISM6:
            {
              for (unsigned int k=0; k<nz; k++)
                for (unsigned int j=0; j<ny; j++)
                  for (unsigned int i=0; i<nx; i++)
                    {
                      // First Prism
                      Elem * elem = new Prism6;
                      elem->set_id(elem_id++);
                      elem = mesh.add_elem (elem);
 
                      elem->set_node(0) = mesh.node_ptr(idx(type,nx,ny,i,j,k)      );
                      elem->set_node(1) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k)    );
                      elem->set_node(2) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k)    );
                      elem->set_node(3) = mesh.node_ptr(idx(type,nx,ny,i,j,k+1)    );
                      elem->set_node(4) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k+1)  );
                      elem->set_node(5) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k+1)  );
 
                      // Add sides for first prism to boundary info object
                      if (i==0)
                        boundary_info.add_side(elem, 3, 4);
 
                      if (j==0)
                        boundary_info.add_side(elem, 1, 1);
 
                      if (k==0)
                        boundary_info.add_side(elem, 0, 0);
 
                      if (k == (nz-1))
                        boundary_info.add_side(elem, 4, 5);
 
                      // Second Prism
                      elem = new Prism6;
                      elem->set_id(elem_id++);
                      elem = mesh.add_elem (elem);
 
                      elem->set_node(0) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k)    );
                      elem->set_node(1) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k)  );
                      elem->set_node(2) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k)    );
                      elem->set_node(3) = mesh.node_ptr(idx(type,nx,ny,i+1,j,k+1)  );
                      elem->set_node(4) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
                      elem->set_node(5) = mesh.node_ptr(idx(type,nx,ny,i,j+1,k+1)  );
 
                      // Add sides for second prism to boundary info object
                      if (i == (nx-1))
                        boundary_info.add_side(elem, 1, 2);
 
                      if (j == (ny-1))
                        boundary_info.add_side(elem, 2, 3);
 
                      if (k==0)
                        boundary_info.add_side(elem, 0, 0);
 
                      if (k == (nz-1))
                        boundary_info.add_side(elem, 4, 5);
                    }
              break;
            }
 
 
 
 
 
 
          case HEX20:
          case HEX27:
          case TET4: // TET4's are created from an initial HEX27 discretization
          case TET10: // TET10's are created from an initial HEX27 discretization
          case PYRAMID5: // PYRAMIDs are created from an initial HEX27 discretization
          case PYRAMID13:
          case PYRAMID14:
            {
              for (unsigned int k=0; k<(2*nz); k += 2)
                for (unsigned int j=0; j<(2*ny); j += 2)
                  for (unsigned int i=0; i<(2*nx); i += 2)
                    {
                      Elem * elem = (type == HEX20) ?
                        static_cast<Elem *>(new Hex20) :
                        static_cast<Elem *>(new Hex27);
                      elem->set_id(elem_id++);
                      elem = mesh.add_elem (elem);
 
                      elem->set_node(0)  = mesh.node_ptr(idx(type,nx,ny,i,  j,  k)  );
                      elem->set_node(1)  = mesh.node_ptr(idx(type,nx,ny,i+2,j,  k)  );
                      elem->set_node(2)  = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k)  );
                      elem->set_node(3)  = mesh.node_ptr(idx(type,nx,ny,i,  j+2,k)  );
                      elem->set_node(4)  = mesh.node_ptr(idx(type,nx,ny,i,  j,  k+2));
                      elem->set_node(5)  = mesh.node_ptr(idx(type,nx,ny,i+2,j,  k+2));
                      elem->set_node(6)  = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k+2));
                      elem->set_node(7)  = mesh.node_ptr(idx(type,nx,ny,i,  j+2,k+2));
                      elem->set_node(8)  = mesh.node_ptr(idx(type,nx,ny,i+1,j,  k)  );
                      elem->set_node(9)  = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k)  );
                      elem->set_node(10) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k)  );
                      elem->set_node(11) = mesh.node_ptr(idx(type,nx,ny,i,  j+1,k)  );
                      elem->set_node(12) = mesh.node_ptr(idx(type,nx,ny,i,  j,  k+1));
                      elem->set_node(13) = mesh.node_ptr(idx(type,nx,ny,i+2,j,  k+1));
                      elem->set_node(14) = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k+1));
                      elem->set_node(15) = mesh.node_ptr(idx(type,nx,ny,i,  j+2,k+1));
                      elem->set_node(16) = mesh.node_ptr(idx(type,nx,ny,i+1,j,  k+2));
                      elem->set_node(17) = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k+2));
                      elem->set_node(18) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k+2));
                      elem->set_node(19) = mesh.node_ptr(idx(type,nx,ny,i,  j+1,k+2));
                      if ((type == HEX27) || (type == TET4) || (type == TET10) ||
                          (type == PYRAMID5) || (type == PYRAMID13) || (type == PYRAMID14))
                        {
                          elem->set_node(20) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k)  );
                          elem->set_node(21) = mesh.node_ptr(idx(type,nx,ny,i+1,j,  k+1));
                          elem->set_node(22) = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k+1));
                          elem->set_node(23) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k+1));
                          elem->set_node(24) = mesh.node_ptr(idx(type,nx,ny,i,  j+1,k+1));
                          elem->set_node(25) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+2));
                          elem->set_node(26) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
                        }
 
 
                      if (k == 0)
                        boundary_info.add_side(elem, 0, 0);
 
                      if (k == 2*(nz-1))
                        boundary_info.add_side(elem, 5, 5);
 
                      if (j == 0)
                        boundary_info.add_side(elem, 1, 1);
 
                      if (j == 2*(ny-1))
                        boundary_info.add_side(elem, 3, 3);
 
                      if (i == 0)
                        boundary_info.add_side(elem, 4, 4);
 
                      if (i == 2*(nx-1))
                        boundary_info.add_side(elem, 2, 2);
                    }
              break;
            }
 
 
 
 
          case PRISM15:
          case PRISM18:
            {
              for (unsigned int k=0; k<(2*nz); k += 2)
                for (unsigned int j=0; j<(2*ny); j += 2)
                  for (unsigned int i=0; i<(2*nx); i += 2)
                    {
                      // First Prism
                      Elem * elem = (type == PRISM15) ?
                        static_cast<Elem *>(new Prism15) :
                        static_cast<Elem *>(new Prism18);
                      elem->set_id(elem_id++);
                      elem = mesh.add_elem (elem);
 
                      elem->set_node(0)  = mesh.node_ptr(idx(type,nx,ny,i,  j,  k)  );
                      elem->set_node(1)  = mesh.node_ptr(idx(type,nx,ny,i+2,j,  k)  );
                      elem->set_node(2)  = mesh.node_ptr(idx(type,nx,ny,i,  j+2,k)  );
                      elem->set_node(3)  = mesh.node_ptr(idx(type,nx,ny,i,  j,  k+2));
                      elem->set_node(4)  = mesh.node_ptr(idx(type,nx,ny,i+2,j,  k+2));
                      elem->set_node(5)  = mesh.node_ptr(idx(type,nx,ny,i,  j+2,k+2));
                      elem->set_node(6)  = mesh.node_ptr(idx(type,nx,ny,i+1,j,  k)  );
                      elem->set_node(7)  = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k)  );
                      elem->set_node(8)  = mesh.node_ptr(idx(type,nx,ny,i,  j+1,k)  );
                      elem->set_node(9)  = mesh.node_ptr(idx(type,nx,ny,i,  j,  k+1));
                      elem->set_node(10) = mesh.node_ptr(idx(type,nx,ny,i+2,j,  k+1));
                      elem->set_node(11) = mesh.node_ptr(idx(type,nx,ny,i,  j+2,k+1));
                      elem->set_node(12) = mesh.node_ptr(idx(type,nx,ny,i+1,j,  k+2));
                      elem->set_node(13) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+2));
                      elem->set_node(14) = mesh.node_ptr(idx(type,nx,ny,i,  j+1,k+2));
                      if (type == PRISM18)
                        {
                          elem->set_node(15) = mesh.node_ptr(idx(type,nx,ny,i+1,j,  k+1));
                          elem->set_node(16) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
                          elem->set_node(17) = mesh.node_ptr(idx(type,nx,ny,i,  j+1,k+1));
                        }
 
                      // Add sides for first prism to boundary info object
                      if (i==0)
                        boundary_info.add_side(elem, 3, 4);
 
                      if (j==0)
                        boundary_info.add_side(elem, 1, 1);
 
                      if (k==0)
                        boundary_info.add_side(elem, 0, 0);
 
                      if (k == 2*(nz-1))
                        boundary_info.add_side(elem, 4, 5);
 
 
                      // Second Prism
                      elem = (type == PRISM15) ?
                        static_cast<Elem *>(new Prism15) :
                        static_cast<Elem *>(new Prism18);
                      elem->set_id(elem_id++);
                      elem = mesh.add_elem (elem);
 
                      elem->set_node(0)  = mesh.node_ptr(idx(type,nx,ny,i+2,j,k)     );
                      elem->set_node(1)  = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k)   );
                      elem->set_node(2)  = mesh.node_ptr(idx(type,nx,ny,i,j+2,k)     );
                      elem->set_node(3)  = mesh.node_ptr(idx(type,nx,ny,i+2,j,k+2)   );
                      elem->set_node(4)  = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k+2) );
                      elem->set_node(5)  = mesh.node_ptr(idx(type,nx,ny,i,j+2,k+2)   );
                      elem->set_node(6)  = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k)  );
                      elem->set_node(7)  = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k)  );
                      elem->set_node(8)  = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k)  );
                      elem->set_node(9)  = mesh.node_ptr(idx(type,nx,ny,i+2,j,k+1)  );
                      elem->set_node(10) = mesh.node_ptr(idx(type,nx,ny,i+2,j+2,k+1));
                      elem->set_node(11) = mesh.node_ptr(idx(type,nx,ny,i,j+2,k+1)  );
                      elem->set_node(12) = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k+2));
                      elem->set_node(13) = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k+2));
                      elem->set_node(14) = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+2));
                      if (type == PRISM18)
                        {
                          elem->set_node(15)  = mesh.node_ptr(idx(type,nx,ny,i+2,j+1,k+1));
                          elem->set_node(16)  = mesh.node_ptr(idx(type,nx,ny,i+1,j+2,k+1));
                          elem->set_node(17)  = mesh.node_ptr(idx(type,nx,ny,i+1,j+1,k+1));
                        }
 
                      // Add sides for second prism to boundary info object
                      if (i == 2*(nx-1))
                        boundary_info.add_side(elem, 1, 2);
 
                      if (j == 2*(ny-1))
                        boundary_info.add_side(elem, 2, 3);
 
                      if (k==0)
                        boundary_info.add_side(elem, 0, 0);
 
                      if (k == 2*(nz-1))
                        boundary_info.add_side(elem, 4, 5);
 
                    }
              break;
            }
 
 
 
 
 
          default:
            libmesh_error_msg("ERROR: Unrecognized 3D element type.");
          }
 
 
 
 
        //.......................................
        // Scale the nodal positions
        if (gauss_lobatto_grid)
          {
            GaussLobattoRedistributionFunction func(nx, xmin, xmax,
                                                    ny, ymin, ymax,
                                                    nz, zmin, zmax);
            MeshTools::Modification::redistribute(mesh, func);
          }
        else // !gauss_lobatto_grid
          {
            for (unsigned int p=0; p<mesh.n_nodes(); p++)
              {
                mesh.node_ref(p)(0) = (mesh.node_ref(p)(0))*(xmax-xmin) + xmin;
                mesh.node_ref(p)(1) = (mesh.node_ref(p)(1))*(ymax-ymin) + ymin;
                mesh.node_ref(p)(2) = (mesh.node_ref(p)(2))*(zmax-zmin) + zmin;
              }
          }
 
 
 
        // Additional work for tets and pyramids: we take the existing
        // HEX27 discretization and split each element into 24
        // sub-tets or 6 sub-pyramids.
        //
        // 24 isn't the minimum-possible number of tets, but it
        // obviates any concerns about the edge orientations between
        // the various elements.
        if ((type == TET4) ||
            (type == TET10) ||
            (type == PYRAMID5) ||
            (type == PYRAMID13) ||
            (type == PYRAMID14))
          {
            // Temporary storage for new elements. (24 tets per hex, 6 pyramids)
            std::vector<Elem *> new_elements;
 
            if ((type == TET4) || (type == TET10))
              new_elements.reserve(24*mesh.n_elem());
            else
              new_elements.reserve(6*mesh.n_elem());
 
            // Create tetrahedra or pyramids
            for (auto & base_hex : mesh.element_ptr_range())
              {
                // Get a pointer to the node located at the HEX27 centroid
                Node * apex_node = base_hex->node_ptr(26);
 
                // Container to catch ids handed back from BoundaryInfo
                std::vector<boundary_id_type> ids;
 
                for (auto s : base_hex->side_index_range())
                  {
                    // Get the boundary ID(s) for this side
                    boundary_info.boundary_ids(base_hex, s, ids);
 
                    // We're creating this Mesh, so there should be 0 or 1 boundary IDs.
                    libmesh_assert(ids.size() <= 1);
 
                    // A convenient name for the side's ID.
                    boundary_id_type b_id = ids.empty() ? BoundaryInfo::invalid_id : ids[0];
 
                    // Need to build the full-ordered side!
                    std::unique_ptr<Elem> side = base_hex->build_side_ptr(s);
 
                    if ((type == TET4) || (type == TET10))
                      {
                        // Build 4 sub-tets per side
                        for (unsigned int sub_tet=0; sub_tet<4; ++sub_tet)
                          {
                            new_elements.push_back( new Tet4 );
                            Elem * sub_elem = new_elements.back();
                            sub_elem->set_node(0) = side->node_ptr(sub_tet);
                            sub_elem->set_node(1) = side->node_ptr(8);                           // centroid of the face
                            sub_elem->set_node(2) = side->node_ptr(sub_tet==3 ? 0 : sub_tet+1 ); // wrap-around
                            sub_elem->set_node(3) = apex_node;                                   // apex node always used!
 
                            // If the original hex was a boundary hex, add the new sub_tet's side
                            // 0 with the same b_id.  Note: the tets are all aligned so that their
                            // side 0 is on the boundary.
                            if (b_id != BoundaryInfo::invalid_id)
                              boundary_info.add_side(sub_elem, 0, b_id);
                          }
                      } // end if ((type == TET4) || (type == TET10))
 
                    else // type==PYRAMID5 || type==PYRAMID13 || type==PYRAMID14
                      {
                        // Build 1 sub-pyramid per side.
                        new_elements.push_back(new Pyramid5);
                        Elem * sub_elem = new_elements.back();
 
                        // Set the base.  Note that since the apex is *inside* the base_hex,
                        // and the pyramid uses a counter-clockwise base numbering, we need to
                        // reverse the [1] and [3] node indices.
                        sub_elem->set_node(0) = side->node_ptr(0);
                        sub_elem->set_node(1) = side->node_ptr(3);
                        sub_elem->set_node(2) = side->node_ptr(2);
                        sub_elem->set_node(3) = side->node_ptr(1);
 
                        // Set the apex
                        sub_elem->set_node(4) = apex_node;
 
                        // If the original hex was a boundary hex, add the new sub_pyr's side
                        // 4 (the square base) with the same b_id.
                        if (b_id != BoundaryInfo::invalid_id)
                          boundary_info.add_side(sub_elem, 4, b_id);
                      } // end else type==PYRAMID5 || type==PYRAMID13 || type==PYRAMID14
                  }
              }
 
 
            // Delete the original HEX27 elements from the mesh, and the boundary info structure.
            for (auto & elem : mesh.element_ptr_range())
              {
                boundary_info.remove(elem); // Safe even if elem has no boundary info.
                mesh.delete_elem(elem);
              }
 
            // Add the new elements
            for (dof_id_type i=0,
                 n_new = cast_int<dof_id_type>(new_elements.size());
                 i != n_new; ++i)
              {
                new_elements[i]->set_id(i);
                mesh.add_elem(new_elements[i]);
              }
 
          } // end if (type == TET4,TET10,PYRAMID5,PYRAMID13,PYRAMID14
 
 
        // Use all_second_order to convert the TET4's to TET10's or PYRAMID5's to PYRAMID14's
        if ((type == TET10) || (type == PYRAMID14))
          mesh.all_second_order();
 
        else if (type == PYRAMID13)
          mesh.all_second_order(/*full_ordered=*/false);
 
 
        // Add sideset names to boundary info (Z axis out of the screen)
        boundary_info.sideset_name(0) = "back";
        boundary_info.sideset_name(1) = "bottom";
        boundary_info.sideset_name(2) = "right";
        boundary_info.sideset_name(3) = "top";
        boundary_info.sideset_name(4) = "left";
        boundary_info.sideset_name(5) = "front";
 
        // Add nodeset names to boundary info
        boundary_info.nodeset_name(0) = "back";
        boundary_info.nodeset_name(1) = "bottom";
        boundary_info.nodeset_name(2) = "right";
        boundary_info.nodeset_name(3) = "top";
        boundary_info.nodeset_name(4) = "left";
        boundary_info.nodeset_name(5) = "front";
 
        break;
      } // end case dim==3
 
    default:
      libmesh_error_msg("Unknown dimension " << mesh.mesh_dimension());
    }
 
  STOP_LOG("build_cube()", "MeshTools::Generation");
 
 
 
  // Done building the mesh.  Now prepare it for use.
  mesh.prepare_for_use (/*skip_renumber =*/ false);
}

References libMesh::BoundaryInfo::add_side(), libMesh::BoundaryInfo::boundary_ids(), libMesh::EDGE2, libMesh::EDGE3, libMesh::EDGE4, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::MeshTools::Generation::Private::idx(), libMesh::INVALID_ELEM, libMesh::BoundaryInfo::invalid_id, libMesh::libmesh_assert(), mesh, libMesh::NODEELEM, libMesh::BoundaryInfo::nodeset_name(), libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM6, libMesh::PYRAMID13, libMesh::PYRAMID14, libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::MeshTools::Modification::redistribute(), libMesh::BoundaryInfo::remove(), libMesh::DofObject::set_id(), libMesh::Elem::set_node(), libMesh::BoundaryInfo::sideset_name(), libMesh::TET10, libMesh::TET4, libMesh::TRI3, and libMesh::TRI6.

Referenced by add_cube_convex_hull_to_mesh(), Biharmonic::Biharmonic(), build_line(), ExtraIntegersTest::build_mesh(), build_point(), build_square(), main(), ElemTest< elem_type >::setUp(), ParsedFEMFunctionTest::setUp(), RationalMapTest< elem_type >::setUp(), setup(), FETest< order, family, elem_type >::setUp(), DistortTest::test_helper_3D(), AllTriTest::test_helper_3D(), MeshFunctionTest::test_p_level(), SystemsTest::testAssemblyWithDgFemContext(), SystemsTest::testBlockRestrictedVarNDofs(), SystemsTest::testBoundaryProjectCube(), DefaultCouplingTest::testCoupling(), PointNeighborCouplingTest::testCoupling(), SystemsTest::testDofCouplingWithVarGroups(), DofMapTest::testDofOwner(), BoundaryInfoTest::testEdgeBoundaryConditions(), MeshInputTest::testExodusWriteElementDataFromDiscontinuousNodalData(), PointLocatorTest::testLocator(), MeshStitchTest::testMeshStitch(), PartitionerTest< PartitionerSubclass, MeshClass >::testPartition(), SystemsTest::testProjectCube(), SystemsTest::testProjectCubeWithMeshFunction(), SystemsTest::testProjectMatrix3D(), InfFERadialTest::testSingleOrder(), and WriteEdgesetData::testWrite().

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 2321 of file mesh_generation.C.

{
  // Check for reasonable size
  libmesh_assert_greater_equal (nx, 1); // need at least 1 element in x-direction
  libmesh_assert_greater_equal (ny, 1); // need at least 1 element in y-direction
  libmesh_assert_less (xmin, xmax);
  libmesh_assert_less (ymin, ymax);
 
  // Clear out any data which may have been in the Mesh
  mesh.clear();
 
  BoundaryInfo & boundary_info = mesh.get_boundary_info();
 
  // Make sure the new Mesh will be 2D
  mesh.set_mesh_dimension(2);
 
  // The x and y spacing between boundary points
  const Real delta_x = (xmax-xmin) / static_cast<Real>(nx);
  const Real delta_y = (ymax-ymin) / static_cast<Real>(ny);
 
  // Bottom
  for (unsigned int p=0; p<=nx; ++p)
    mesh.add_point(Point(xmin + p*delta_x, ymin));
 
  // Right side
  for (unsigned int p=1; p<ny; ++p)
    mesh.add_point(Point(xmax, ymin + p*delta_y));
 
  // Top
  for (unsigned int p=0; p<=nx; ++p)
    mesh.add_point(Point(xmax - p*delta_x, ymax));
 
  // Left side
  for (unsigned int p=1; p<ny; ++p)
    mesh.add_point(Point(xmin,  ymax - p*delta_y));
 
  // Be sure we added as many points as we thought we did
  libmesh_assert_equal_to (mesh.n_nodes(), 2*(nx+ny));
 
  // Construct the Triangle Interface object
  TriangleInterface t(mesh);
 
  // Set custom variables for the triangulation
  t.desired_area()       = 0.5 * (xmax-xmin)*(ymax-ymin) / static_cast<Real>(nx*ny);
  t.triangulation_type() = TriangleInterface::PSLG;
  t.elem_type()          = type;
 
  if (holes != nullptr)
    t.attach_hole_list(holes);
 
  // Triangulate!
  t.triangulate();
 
  // The mesh is now generated, but we still need to mark the boundaries
  // to be consistent with the other build_square routines.  Note that all
  // hole boundary elements get the same ID, 4.
  for (auto & elem : mesh.element_ptr_range())
    for (auto s : elem->side_index_range())
      if (elem->neighbor_ptr(s) == nullptr)
        {
          std::unique_ptr<const Elem> side (elem->build_side_ptr(s));
 
          // Check the location of the side's midpoint.  Since
          // the square has straight sides, the midpoint is not
          // on the corner and thus it is uniquely on one of the
          // sides.
          Point side_midpoint= 0.5f*( side->point(0) + side->point(1) );
 
          // The boundary ids are set following the same convention as Quad4 sides
          // bottom = 0
          // right  = 1
          // top = 2
          // left = 3
          // hole = 4
          boundary_id_type bc_id=4;
 
          // bottom
          if      (std::fabs(side_midpoint(1) - ymin) < TOLERANCE)
            bc_id=0;
 
          // right
          else if (std::fabs(side_midpoint(0) - xmax) < TOLERANCE)
            bc_id=1;
 
          // top
          else if (std::fabs(side_midpoint(1) - ymax) < TOLERANCE)
            bc_id=2;
 
          // left
          else if (std::fabs(side_midpoint(0) - xmin) < TOLERANCE)
            bc_id=3;
 
          // If the point is not on any of the external boundaries, it
          // is on one of the holes....
 
          // Finally, add this element's information to the boundary info object.
          boundary_info.add_side(elem->id(), s, bc_id);
        }
 
} // end build_delaunay_square

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

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 1994 of file mesh_generation.C.

{
  if (!cross_section.n_elem())
    return;
 
  START_LOG("build_extrusion()", "MeshTools::Generation");
 
  dof_id_type orig_elem = cross_section.n_elem();
  dof_id_type orig_nodes = cross_section.n_nodes();
 
#ifdef LIBMESH_ENABLE_UNIQUE_ID
  unique_id_type orig_unique_ids = cross_section.parallel_max_unique_id();
#endif
 
  unsigned int order = 1;
 
  BoundaryInfo & boundary_info = mesh.get_boundary_info();
  const BoundaryInfo & cross_section_boundary_info = cross_section.get_boundary_info();
 
  // If cross_section is distributed, so is its extrusion
  if (!cross_section.is_serial())
    mesh.delete_remote_elements();
 
  // We know a priori how many elements we'll need
  mesh.reserve_elem(nz*orig_elem);
 
  // For straightforward meshes we need one or two additional layers per
  // element.
  if (cross_section.elements_begin() != cross_section.elements_end() &&
      (*cross_section.elements_begin())->default_order() == SECOND)
    order = 2;
  mesh.comm().max(order);
 
  mesh.reserve_nodes((order*nz+1)*orig_nodes);
 
  // Container to catch the boundary IDs handed back by the BoundaryInfo object
  std::vector<boundary_id_type> ids_to_copy;
 
  for (const auto & node : cross_section.node_ptr_range())
    {
      for (unsigned int k=0; k != order*nz+1; ++k)
        {
          Node * new_node =
            mesh.add_point(*node +
                           (extrusion_vector * k / nz / order),
                           node->id() + (k * orig_nodes),
                           node->processor_id());
 
#ifdef LIBMESH_ENABLE_UNIQUE_ID
          // Let's give the base of the extruded mesh the same
          // unique_ids as the source mesh, in case anyone finds that
          // a useful map to preserve.
          const unique_id_type uid = (k == 0) ?
            node->unique_id() :
            orig_unique_ids + (k-1)*(orig_nodes + orig_elem) + node->id();
 
          new_node->set_unique_id() = uid;
#endif
 
          cross_section_boundary_info.boundary_ids(node, ids_to_copy);
          boundary_info.add_node(new_node, ids_to_copy);
        }
    }
 
  const std::set<boundary_id_type> & side_ids =
    cross_section_boundary_info.get_side_boundary_ids();
 
  boundary_id_type next_side_id = side_ids.empty() ?
    0 : cast_int<boundary_id_type>(*side_ids.rbegin() + 1);
 
  // side_ids may not include ids from remote elements, in which case
  // some processors may have underestimated the next_side_id; let's
  // fix that.
  cross_section.comm().max(next_side_id);
 
  for (const auto & elem : cross_section.element_ptr_range())
    {
      const ElemType etype = elem->type();
 
      // build_extrusion currently only works on coarse meshes
      libmesh_assert (!elem->parent());
 
      for (unsigned int k=0; k != nz; ++k)
        {
          Elem * new_elem;
          switch (etype)
            {
            case EDGE2:
              {
                new_elem = new Quad4;
                new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (k * orig_nodes));
                new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (k * orig_nodes));
                new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(1)->id() + ((k+1) * orig_nodes));
                new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(0)->id() + ((k+1) * orig_nodes));
 
                if (elem->neighbor_ptr(0) == remote_elem)
                  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(1) == remote_elem)
                  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
 
                break;
              }
            case EDGE3:
              {
                new_elem = new Quad9;
                new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (2*k * orig_nodes));
                new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (2*k * orig_nodes));
                new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(2)->id() + (2*k * orig_nodes));
                new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(6) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(7) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(8) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+1) * orig_nodes));
 
                if (elem->neighbor_ptr(0) == remote_elem)
                  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(1) == remote_elem)
                  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
 
                break;
              }
            case TRI3:
              {
                new_elem = new Prism6;
                new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (k * orig_nodes));
                new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (k * orig_nodes));
                new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(2)->id() + (k * orig_nodes));
                new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(0)->id() + ((k+1) * orig_nodes));
                new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(1)->id() + ((k+1) * orig_nodes));
                new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(2)->id() + ((k+1) * orig_nodes));
 
                if (elem->neighbor_ptr(0) == remote_elem)
                  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(1) == remote_elem)
                  new_elem->set_neighbor(2, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(2) == remote_elem)
                  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
 
                break;
              }
            case TRI6:
              {
                new_elem = new Prism18;
                new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (2*k * orig_nodes));
                new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (2*k * orig_nodes));
                new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(2)->id() + (2*k * orig_nodes));
                new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(6) = mesh.node_ptr(elem->node_ptr(3)->id() + (2*k * orig_nodes));
                new_elem->set_node(7) = mesh.node_ptr(elem->node_ptr(4)->id() + (2*k * orig_nodes));
                new_elem->set_node(8) = mesh.node_ptr(elem->node_ptr(5)->id() + (2*k * orig_nodes));
                new_elem->set_node(9) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(10) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(11) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(12) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(13) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(14) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(15) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(16) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(17) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+1) * orig_nodes));
 
                if (elem->neighbor_ptr(0) == remote_elem)
                  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(1) == remote_elem)
                  new_elem->set_neighbor(2, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(2) == remote_elem)
                  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
 
                break;
              }
            case QUAD4:
              {
                new_elem = new Hex8;
                new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (k * orig_nodes));
                new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (k * orig_nodes));
                new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(2)->id() + (k * orig_nodes));
                new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(3)->id() + (k * orig_nodes));
                new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(0)->id() + ((k+1) * orig_nodes));
                new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(1)->id() + ((k+1) * orig_nodes));
                new_elem->set_node(6) = mesh.node_ptr(elem->node_ptr(2)->id() + ((k+1) * orig_nodes));
                new_elem->set_node(7) = mesh.node_ptr(elem->node_ptr(3)->id() + ((k+1) * orig_nodes));
 
                if (elem->neighbor_ptr(0) == remote_elem)
                  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(1) == remote_elem)
                  new_elem->set_neighbor(2, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(2) == remote_elem)
                  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(3) == remote_elem)
                  new_elem->set_neighbor(4, const_cast<RemoteElem *>(remote_elem));
 
                break;
              }
            case QUAD9:
              {
                new_elem = new Hex27;
                new_elem->set_node(0) = mesh.node_ptr(elem->node_ptr(0)->id() + (2*k * orig_nodes));
                new_elem->set_node(1) = mesh.node_ptr(elem->node_ptr(1)->id() + (2*k * orig_nodes));
                new_elem->set_node(2) = mesh.node_ptr(elem->node_ptr(2)->id() + (2*k * orig_nodes));
                new_elem->set_node(3) = mesh.node_ptr(elem->node_ptr(3)->id() + (2*k * orig_nodes));
                new_elem->set_node(4) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(5) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(6) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(7) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(8) = mesh.node_ptr(elem->node_ptr(4)->id() + (2*k * orig_nodes));
                new_elem->set_node(9) = mesh.node_ptr(elem->node_ptr(5)->id() + (2*k * orig_nodes));
                new_elem->set_node(10) = mesh.node_ptr(elem->node_ptr(6)->id() + (2*k * orig_nodes));
                new_elem->set_node(11) = mesh.node_ptr(elem->node_ptr(7)->id() + (2*k * orig_nodes));
                new_elem->set_node(12) = mesh.node_ptr(elem->node_ptr(0)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(13) = mesh.node_ptr(elem->node_ptr(1)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(14) = mesh.node_ptr(elem->node_ptr(2)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(15) = mesh.node_ptr(elem->node_ptr(3)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(16) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(17) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(18) = mesh.node_ptr(elem->node_ptr(6)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(19) = mesh.node_ptr(elem->node_ptr(7)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(20) = mesh.node_ptr(elem->node_ptr(8)->id() + (2*k * orig_nodes));
                new_elem->set_node(21) = mesh.node_ptr(elem->node_ptr(4)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(22) = mesh.node_ptr(elem->node_ptr(5)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(23) = mesh.node_ptr(elem->node_ptr(6)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(24) = mesh.node_ptr(elem->node_ptr(7)->id() + ((2*k+1) * orig_nodes));
                new_elem->set_node(25) = mesh.node_ptr(elem->node_ptr(8)->id() + ((2*k+2) * orig_nodes));
                new_elem->set_node(26) = mesh.node_ptr(elem->node_ptr(8)->id() + ((2*k+1) * orig_nodes));
 
                if (elem->neighbor_ptr(0) == remote_elem)
                  new_elem->set_neighbor(1, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(1) == remote_elem)
                  new_elem->set_neighbor(2, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(2) == remote_elem)
                  new_elem->set_neighbor(3, const_cast<RemoteElem *>(remote_elem));
                if (elem->neighbor_ptr(3) == remote_elem)
                  new_elem->set_neighbor(4, const_cast<RemoteElem *>(remote_elem));
 
                break;
              }
            default:
              {
                libmesh_not_implemented();
                break;
              }
            }
 
          new_elem->set_id(elem->id() + (k * orig_elem));
          new_elem->processor_id() = elem->processor_id();
 
#ifdef LIBMESH_ENABLE_UNIQUE_ID
          // Let's give the base of the extruded mesh the same
          // unique_ids as the source mesh, in case anyone finds that
          // a useful map to preserve.
          const unique_id_type uid = (k == 0) ?
            elem->unique_id() :
            orig_unique_ids + (k-1)*(orig_nodes + orig_elem) + orig_nodes + elem->id();
 
          new_elem->set_unique_id() = uid;
#endif
 
          if (!elem_subdomain)
            // maintain the subdomain_id
            new_elem->subdomain_id() = elem->subdomain_id();
          else
            // Allow the user to choose new subdomain_ids
            new_elem->subdomain_id() = elem_subdomain->get_subdomain_for_layer(elem, k);
 
          new_elem = mesh.add_elem(new_elem);
 
          // Copy any old boundary ids on all sides
          for (auto s : elem->side_index_range())
            {
              cross_section_boundary_info.boundary_ids(elem, s, ids_to_copy);
 
              if (new_elem->dim() == 3)
                {
                  // For 2D->3D extrusion, we give the boundary IDs
                  // for side s on the old element to side s+1 on the
                  // new element.  This is just a happy coincidence as
                  // far as I can tell...
                  boundary_info.add_side(new_elem,
                                         cast_int<unsigned short>(s+1),
                                         ids_to_copy);
                }
              else
                {
                  // For 1D->2D extrusion, the boundary IDs map as:
                  // Old elem -> New elem
                  // 0        -> 3
                  // 1        -> 1
                  libmesh_assert_less(s, 2);
                  const unsigned short sidemap[2] = {3, 1};
                  boundary_info.add_side(new_elem, sidemap[s], ids_to_copy);
                }
            }
 
          // Give new boundary ids to bottom and top
          if (k == 0)
            boundary_info.add_side(new_elem, 0, next_side_id);
          if (k == nz-1)
            {
              // For 2D->3D extrusion, the "top" ID is 1+the original
              // element's number of sides.  For 1D->2D extrusion, the
              // "top" ID is side 2.
              const unsigned short top_id = new_elem->dim() == 3 ?
                cast_int<unsigned short>(elem->n_sides()+1) : 2;
              boundary_info.add_side
                (new_elem, top_id,
                 cast_int<boundary_id_type>(next_side_id+1));
            }
        }
    }
 
  STOP_LOG("build_extrusion()", "MeshTools::Generation");
 
  // Done building the mesh.  Now prepare it for use.
  mesh.prepare_for_use(/*skip_renumber =*/ false);
}

References libMesh::MeshBase::add_elem(), libMesh::BoundaryInfo::add_node(), libMesh::MeshBase::add_point(), libMesh::BoundaryInfo::add_side(), libMesh::BoundaryInfo::boundary_ids(), libMesh::ParallelObject::comm(), libMesh::MeshBase::delete_remote_elements(), libMesh::Elem::dim(), libMesh::EDGE2, libMesh::EDGE3, libMesh::MeshBase::element_ptr_range(), libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), libMesh::MeshBase::get_boundary_info(), libMesh::BoundaryInfo::get_side_boundary_ids(), libMesh::MeshTools::Generation::QueryElemSubdomainIDBase::get_subdomain_for_layer(), libMesh::MeshBase::is_serial(), libMesh::libmesh_assert(), mesh, libMesh::MeshBase::n_elem(), libMesh::MeshBase::n_nodes(), libMesh::MeshBase::node_ptr(), libMesh::MeshBase::node_ptr_range(), libMesh::MeshBase::parallel_max_unique_id(), libMesh::MeshBase::prepare_for_use(), libMesh::DofObject::processor_id(), libMesh::QUAD4, libMesh::QUAD9, libMesh::remote_elem, libMesh::MeshBase::reserve_elem(), libMesh::MeshBase::reserve_nodes(), libMesh::SECOND, libMesh::DofObject::set_id(), libMesh::Elem::set_neighbor(), libMesh::Elem::set_node(), libMesh::DofObject::set_unique_id(), libMesh::Elem::subdomain_id(), libMesh::TRI3, libMesh::TRI6, and libMesh::DofObject::unique_id().

Referenced by MeshExtruderTest::testExtruder().

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 1480 of file mesh_generation.C.

{
  // This method only makes sense in 1D!
  // But we now just turn a non-1D mesh into a 1D mesh
  //libmesh_assert_equal_to (mesh.mesh_dimension(), 1);
 
  build_cube(mesh,
             nx, 0, 0,
             xmin, xmax,
             0., 0.,
             0., 0.,
             type,
             gauss_lobatto_grid);
}

References build_cube(), and mesh.

Referenced by Biharmonic::Biharmonic(), NodalNeighborsTest::do_test(), main(), MeshSpatialDimensionTest::test1D(), EquationSystemsTest::testPostInitAddElem(), SystemsTest::testProjectLine(), SystemsTest::testProjectMatrix1D(), and EquationSystemsTest::testReinitWithNodeElem().

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 1462 of file mesh_generation.C.

{
  // This method only makes sense in 0D!
  // But we now just turn a non-0D mesh into a 0D mesh
  //libmesh_assert_equal_to (mesh.mesh_dimension(), 1);
 
  build_cube(mesh,
             0, 0, 0,
             0., 0.,
             0., 0.,
             0., 0.,
             type,
             gauss_lobatto_grid);
}

References build_cube(), and mesh.

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

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 1532 of file mesh_generation.C.

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

Referenced by main().

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 1501 of file mesh_generation.C.

{
  // This method only makes sense in 2D!
  // But we now just turn a non-2D mesh into a 2D mesh
  //libmesh_assert_equal_to (mesh.mesh_dimension(), 2);
 
  // Call the build_cube() member to actually do the work for us.
  build_cube (mesh,
              nx, ny, 0,
              xmin, xmax,
              ymin, ymax,
              0., 0.,
              type,
              gauss_lobatto_grid);
}

References build_cube(), and mesh.

Referenced by AllSecondOrderTest::allSecondOrder(), Biharmonic::Biharmonic(), BoundaryMeshTest::build_mesh(), main(), MeshSpatialDimensionTest::test2D(), DistortTest::test_helper_2D(), AllTriTest::test_helper_2D(), DofMapTest::testConstraintLoopDetection(), MeshInputTest::testExodusCopyElementSolution(), MeshExtruderTest::testExtruder(), MappedSubdomainPartitionerTest::testMappedSubdomainPartitioner(), BoundaryInfoTest::testMesh(), MeshInputTest::testMeshMoveConstructor(), SystemsTest::testProjectMatrix2D(), SystemsTest::testProjectSquare(), EquationSystemsTest::testRefineThenReinitPreserveFlags(), EquationSystemsTest::testRepartitionThenReinit(), CheckpointIOTest::testSplitter(), WriteSidesetData::testWrite(), and WriteVecAndScalar::testWrite().