C++ MeshBase::add_elem方法代码示例

void add_cube_convex_hull_to_mesh(MeshBase& mesh, Point lower_limit, Point upper_limit)
{
#ifdef LIBMESH_HAVE_TETGEN
  SerialMesh cube_mesh(mesh.comm(),3);

  unsigned n_elem = 1;

  MeshTools::Generation::build_cube(cube_mesh,
                                    n_elem,n_elem,n_elem, // n. elements in each direction
                                    lower_limit(0), upper_limit(0),
                                    lower_limit(1), upper_limit(1),
                                    lower_limit(2), upper_limit(2),
                                    HEX8);

  // The pointset_convexhull() algorithm will ignore the Hex8s
  // in the Mesh, and just construct the triangulation
  // of the convex hull.
  TetGenMeshInterface t(cube_mesh);
  t.pointset_convexhull();

  // Now add all nodes from the boundary of the cube_mesh to the input mesh.

  // Map from "node id in cube_mesh" -> "node id in mesh".  Initially inserted
  // with a dummy value, later to be assigned a value by the input mesh.
  std::map<unsigned,unsigned> node_id_map;
  typedef std::map<unsigned,unsigned>::iterator iterator;

  {
    MeshBase::element_iterator it = cube_mesh.elements_begin();
    const MeshBase::element_iterator end = cube_mesh.elements_end();
    for ( ; it != end; ++it)
      {
        Elem* elem = *it;

        for (unsigned s=0; s<elem->n_sides(); ++s)
          if (elem->neighbor(s) == NULL)
            {
              // Add the node IDs of this side to the set
              AutoPtr<Elem> side = elem->side(s);

              for (unsigned n=0; n<side->n_nodes(); ++n)
                node_id_map.insert( std::make_pair(side->node(n), /*dummy_value=*/0) );
            }
      }
  }

  // For each node in the map, insert it into the input mesh and keep
  // track of the ID assigned.
  for (iterator it=node_id_map.begin(); it != node_id_map.end(); ++it)
    {
      // Id of the node in the cube mesh
      unsigned id = (*it).first;

      // Pointer to node in the cube mesh
      Node* old_node = cube_mesh.node_ptr(id);

      // Add geometric point to input mesh
      Node* new_node = mesh.add_point ( *old_node );

      // Track ID value of new_node in map
      (*it).second = new_node->id();
    }

  // With the points added and the map data structure in place, we are
  // ready to add each TRI3 element of the cube_mesh to the input Mesh
  // with proper node assignments
  {
    MeshBase::element_iterator       el     = cube_mesh.elements_begin();
    const MeshBase::element_iterator end_el = cube_mesh.elements_end();

    for (; el != end_el; ++el)
      {
        Elem* old_elem = *el;

        if (old_elem->type() == TRI3)
          {
            Elem* new_elem = mesh.add_elem(new Tri3);

            // Assign nodes in new elements.  Since this is an example,
            // we'll do it in several steps.
            for (unsigned i=0; i<old_elem->n_nodes(); ++i)
              {
                // Locate old node ID in the map
                iterator it = node_id_map.find(old_elem->node(i));

                // Check for not found
                if (it == node_id_map.end())
                  {
                    libMesh::err << "Node id " << old_elem->node(i) << " not found in map!" << std::endl;
                    libmesh_error();
                  }

                // Mapping to node ID in input mesh
                unsigned new_node_id = (*it).second;

                // Node pointer assigned from input mesh
                new_elem->set_node(i) = mesh.node_ptr(new_node_id);
              }
          }
      }
//.........这里部分代码省略.........

void MeshTools::Subdivision::add_boundary_ghosts(MeshBase & mesh)
{
  static const Real tol = 1e-5;

  // add the mirrored ghost elements (without using iterators, because the mesh is modified in the course)
  std::vector<Tri3Subdivision *> ghost_elems;
  std::vector<Node *> ghost_nodes;
  const unsigned int n_elem = mesh.n_elem();
  for (unsigned int eid = 0; eid < n_elem; ++eid)
    {
      Elem * elem = mesh.elem(eid);
      libmesh_assert_equal_to(elem->type(), TRI3SUBDIVISION);

      // If the triangle happens to be in a corner (two boundary
      // edges), we perform a counter-clockwise loop by mirroring the
      // previous triangle until we come back to the original
      // triangle.  This prevents degenerated triangles in the mesh
      // corners and guarantees that the node in the middle of the
      // loop is of valence=6.
      for (unsigned int i = 0; i < elem->n_sides(); ++i)
        {
          libmesh_assert_not_equal_to(elem->neighbor(i), elem);

          if (elem->neighbor(i) == libmesh_nullptr &&
              elem->neighbor(next[i]) == libmesh_nullptr)
            {
              Elem * nelem = elem;
              unsigned int k = i;
              for (unsigned int l=0;l<4;l++)
                {
                  // this is the vertex to be mirrored
                  Point point = nelem->point(k) + nelem->point(next[k]) - nelem->point(prev[k]);

                  // Check if the proposed vertex doesn't coincide
                  // with one of the existing vertices.  This is
                  // necessary because for some triangulations, it can
                  // happen that two mirrored ghost vertices coincide,
                  // which would then lead to a zero size ghost
                  // element below.
                  Node * node = libmesh_nullptr;
                  for (unsigned int j = 0; j < ghost_nodes.size(); ++j)
                    {
                      if ((*ghost_nodes[j] - point).size() < tol * (elem->point(k) - point).size())
                        {
                          node = ghost_nodes[j];
                          break;
                        }
                    }

                  // add the new vertex only if no other is nearby
                  if (node == libmesh_nullptr)
                    {
                      node = mesh.add_point(point);
                      ghost_nodes.push_back(node);
                    }

                  Tri3Subdivision * newelem = new Tri3Subdivision();

                  // add the first new ghost element to the list just as in the non-corner case
                  if (l == 0)
                    ghost_elems.push_back(newelem);

                  newelem->set_node(0) = nelem->get_node(next[k]);
                  newelem->set_node(1) = nelem->get_node(k);
                  newelem->set_node(2) = node;
                  newelem->set_neighbor(0, nelem);
                  newelem->set_ghost(true);
                  if (l>0)
                    newelem->set_neighbor(2, libmesh_nullptr);
                  nelem->set_neighbor(k, newelem);

                  mesh.add_elem(newelem);
                  mesh.get_boundary_info().add_node(nelem->get_node(k), 1);
                  mesh.get_boundary_info().add_node(nelem->get_node(next[k]), 1);
                  mesh.get_boundary_info().add_node(nelem->get_node(prev[k]), 1);
                  mesh.get_boundary_info().add_node(node, 1);

                  nelem = newelem;
                  k = 2 ;
                }

              Tri3Subdivision * newelem = new Tri3Subdivision();

              newelem->set_node(0) = elem->get_node(next[i]);
              newelem->set_node(1) = nelem->get_node(2);
              newelem->set_node(2) = elem->get_node(prev[i]);
              newelem->set_neighbor(0, nelem);
              nelem->set_neighbor(2, newelem);
              newelem->set_ghost(true);
              newelem->set_neighbor(2, elem);
              elem->set_neighbor(next[i],newelem);

              mesh.add_elem(newelem);

              break;
            }
        }

      for (unsigned int i = 0; i < elem->n_sides(); ++i)
        {
//.........这里部分代码省略.........