Poster of Linux kernelThe best gift for a Linux geek
MeshCommunication

MeshCommunication

Section: C Library Functions (3) Updated: Thu Apr 7 2011
Local index Up
 

NAME

MeshCommunication -  

SYNOPSIS


#include <mesh_communication.h>  

Public Member Functions


MeshCommunication ()

~MeshCommunication ()

void clear ()

void broadcast (MeshBase &) const

void redistribute (ParallelMesh &) const

void gather_neighboring_elements (ParallelMesh &) const

void allgather (ParallelMesh &) const

void delete_remote_elements (ParallelMesh &) const

void assign_global_indices (MeshBase &) const

template<typename ForwardIterator > void find_global_indices (const MeshTools::BoundingBox &, const ForwardIterator &, const ForwardIterator &, std::vector< unsigned int > &) const

void make_elems_parallel_consistent (MeshBase &)

void make_node_ids_parallel_consistent (MeshBase &, LocationMap< Node > &)

void make_node_proc_ids_parallel_consistent (MeshBase &, LocationMap< Node > &)

void make_nodes_parallel_consistent (MeshBase &, LocationMap< Node > &)
 

Private Member Functions


void broadcast_mesh (MeshBase &) const

void broadcast_bcs (const MeshBase &, BoundaryInfo &) const

void allgather_mesh (ParallelMesh &) const

void allgather_bcs (const ParallelMesh &, BoundaryInfo &) const
 

Detailed Description

This is the MeshCommunication class. It handles all the details of communicating mesh information from one processor to another. All parallelization of the Mesh data structures is done via this class.

Author:

Benjamin S. Kirk, 2003

Definition at line 52 of file mesh_communication.h.  

Constructor & Destructor Documentation

 

MeshCommunication::MeshCommunication () [inline]Constructor.

Definition at line 59 of file mesh_communication.h.

{}
 

MeshCommunication::~MeshCommunication () [inline]Destructor.

Definition at line 64 of file mesh_communication.h.

{}
 

Member Function Documentation

 

void MeshCommunication::allgather (ParallelMesh &mesh) constThis method takes an input ParallelMesh which may be distributed among all the processors. Each processor then sends its local nodes and elements to the other processors. The end result is that a previously distributed ParallelMesh will be serialized on each processor. Since this method is collective it must be called by all processors.

Definition at line 1729 of file mesh_communication.C.

References allgather_bcs(), allgather_mesh(), MeshBase::boundary_info, UnstructuredMesh::find_neighbors(), and ParallelMesh::is_serial().

Referenced by ParallelMesh::allgather(), assign_global_indices(), and find_global_indices().

{
  START_LOG('allgather()','MeshCommunication');

  // The mesh should know it's about to be serialized
  libmesh_assert (mesh.is_serial());

  this->allgather_mesh (mesh);
  this->allgather_bcs  (mesh, *(mesh.boundary_info));

  // Inform new elements of their neighbors,
  // while resetting all remote_elem links
  mesh.find_neighbors(true);

  STOP_LOG('allgather()','MeshCommunication');
}
 

void MeshCommunication::allgather_bcs (const ParallelMesh &mesh, BoundaryInfo &boundary_info) const [private]

Definition at line 1756 of file mesh_communication.C.

Referenced by allgather().

{
  // NO MPI == one processor, no need for this method
  return;
}
 

void MeshCommunication::allgather_mesh (ParallelMesh &mesh) const [private]

Definition at line 1749 of file mesh_communication.C.

Referenced by allgather().

{
  // NO MPI == one processor, no need for this method
  return;
}
 

void MeshCommunication::assign_global_indices (MeshBase &mesh) constThis method assigns globally unique, partition-agnostic indices to the nodes and elements in the mesh. The approach is to compute the Hilbert space-filling curve key and use its value to assign an index in [0,N_global). Since the Hilbert key is unique for each spatial location, two objects occupying the same location will be assigned the same global id. Thus, this method can also be useful for identifying duplicate nodes which may occur during parallel refinement.

Definition at line 172 of file mesh_communication_global_indices.C.

References allgather(), Parallel::Sort< KeyType >::bin(), MeshTools::bounding_box(), Elem::centroid(), libMesh::COMM_WORLD, MeshBase::elements_begin(), MeshBase::elements_end(), MeshTools::Generation::Private::idx(), MeshBase::local_elements_begin(), MeshBase::local_elements_end(), MeshBase::local_nodes_begin(), MeshBase::local_nodes_end(), MeshBase::n_elem(), MeshBase::n_nodes(), libMesh::n_processors(), MeshBase::nodes_begin(), MeshBase::nodes_end(), Threads::parallel_for(), libMesh::processor_id(), DofObject::set_id(), and Parallel::Sort< KeyType >::sort().

Referenced by MeshTools::Private::globally_renumber_nodes_and_elements().

{
  START_LOG ('assign_global_indices()', 'MeshCommunication');

  // This method determines partition-agnostic global indices
  // for nodes and elements.

  // Algorithm:
  // (1) compute the Hilbert key for each local node/element
  // (2) perform a parallel sort of the Hilbert key
  // (3) get the min/max value on each processor
  // (4) determine the position in the global ranking for
  //     each local object

  // Global bounding box
  MeshTools::BoundingBox bbox =
    MeshTools::bounding_box (mesh);

  // Set up a derived MPI datatype to handle communication of HilbertIndices
  MPI_Datatype hilbert_type;
  MPI_Type_contiguous (3, MPI_UNSIGNED, &hilbert_type);
  MPI_Type_commit     (&hilbert_type);


  //-------------------------------------------------------------
  // (1) compute Hilbert keys
  std::vector<Hilbert::HilbertIndices>
    node_keys, elem_keys;
  
  {
    // Nodes first
    {
      ConstNodeRange nr (mesh.local_nodes_begin(),
                         mesh.local_nodes_end());
      node_keys.resize (nr.size()); 
      Threads::parallel_for (nr, ComputeHilbertKeys (bbox, node_keys));
    }
    
    // Elements next
    {
      ConstElemRange er (mesh.local_elements_begin(),
                         mesh.local_elements_end());
      elem_keys.resize (er.size());
      Threads::parallel_for (er, ComputeHilbertKeys (bbox, elem_keys));
    }    
  } // done computing Hilbert keys

  
 
  //-------------------------------------------------------------
  // (2) parallel sort the Hilbert keys
  Parallel::Sort<Hilbert::HilbertIndices> node_sorter (node_keys);
  node_sorter.sort(); /* done with node_keys */ //node_keys.clear();
    
  const std::vector<Hilbert::HilbertIndices> &my_node_bin = 
    node_sorter.bin();

  Parallel::Sort<Hilbert::HilbertIndices> elem_sorter (elem_keys);
  elem_sorter.sort(); /* done with elem_keys */ //elem_keys.clear();
    
  const std::vector<Hilbert::HilbertIndices> &my_elem_bin = 
    elem_sorter.bin();



  //-------------------------------------------------------------
  // (3) get the max value on each processor
  std::vector<Hilbert::HilbertIndices>    
    node_upper_bounds(libMesh::n_processors()),
    elem_upper_bounds(libMesh::n_processors());

  { // limit scope of temporaries
    std::vector<Hilbert::HilbertIndices> recvbuf(2*libMesh::n_processors());
    std::vector<unsigned short int> /* do not use a vector of bools here since it is not always so! */
      empty_nodes (libMesh::n_processors()),
      empty_elem  (libMesh::n_processors());
    Hilbert::HilbertIndices my_max[2];
    
    Parallel::allgather (static_cast<unsigned short int>(my_node_bin.empty()), empty_nodes);
    Parallel::allgather (static_cast<unsigned short int>(my_elem_bin.empty()),  empty_elem);
               
    if (!my_node_bin.empty()) my_max[0] = my_node_bin.back();
    if (!my_elem_bin.empty()) my_max[1] = my_elem_bin.back();

    MPI_Allgather (my_max,      2, hilbert_type,
                   &recvbuf[0], 2, hilbert_type,
                   libMesh::COMM_WORLD);

    // Be cereful here.  The *_upper_bounds will be used to find the processor
    // a given object belongs to.  So, if a processor contains no objects (possible!)
    // then copy the bound from the lower processor id.
    for (unsigned int p=0; p<libMesh::n_processors(); p++)
      {
        node_upper_bounds[p] = recvbuf[2*p+0];
        elem_upper_bounds[p] = recvbuf[2*p+1];

        if (p > 0) // default hilbert index value is the OK upper bound for processor 0.
          {
            if (empty_nodes[p]) node_upper_bounds[p] = node_upper_bounds[p-1];
            if (empty_elem[p])  elem_upper_bounds[p] = elem_upper_bounds[p-1];
          }
      }
  }



  //-------------------------------------------------------------
  // (4) determine the position in the global ranking for
  //     each local object
  {
    //----------------------------------------------
    // Nodes first -- all nodes, not just local ones
    {
      // Request sets to send to each processor
      std::vector<std::vector<Hilbert::HilbertIndices> > 
        requested_ids (libMesh::n_processors());
      // Results to gather from each processor
      std::vector<std::vector<unsigned int> >
        filled_request (libMesh::n_processors());

      MeshBase::const_node_iterator       it  = mesh.nodes_begin();
      const MeshBase::const_node_iterator end = mesh.nodes_end();

      // build up list of requests
      for (; it != end; ++it)
        {
          const Node* node = (*it);
          libmesh_assert (node != NULL);
          const Hilbert::HilbertIndices hi = 
            get_hilbert_index (*node, bbox);
          const unsigned int pid = 
            std::distance (node_upper_bounds.begin(), 
                           std::lower_bound(node_upper_bounds.begin(), 
                                            node_upper_bounds.end(),
                                            hi));

          libmesh_assert (pid < libMesh::n_processors());

          requested_ids[pid].push_back(hi);
        }

      // The number of objects in my_node_bin on each processor
      std::vector<unsigned int> node_bin_sizes(libMesh::n_processors());
      Parallel::allgather (static_cast<unsigned int>(my_node_bin.size()), node_bin_sizes);

      // The offset of my first global index
      unsigned int my_offset = 0;
      for (unsigned int pid=0; pid<libMesh::processor_id(); pid++)
        my_offset += node_bin_sizes[pid];

      // start with pid=0, so that we will trade with ourself
      for (unsigned int pid=0; pid<libMesh::n_processors(); pid++)
        {
          // Trade my requests with processor procup and procdown
          const unsigned int procup = (libMesh::processor_id() + pid) %
                                       libMesh::n_processors();
          const unsigned int procdown = (libMesh::n_processors() +
                                         libMesh::processor_id() - pid) %
                                         libMesh::n_processors();

          std::vector<Hilbert::HilbertIndices> request_to_fill;
          Parallel::send_receive(procup, requested_ids[procup],
                                 procdown, request_to_fill,
                                 hilbert_type);   

          // Fill the requests
          std::vector<unsigned int> global_ids;  global_ids.reserve(request_to_fill.size());
          for (unsigned int idx=0; idx<request_to_fill.size(); idx++)
            {
              const Hilbert::HilbertIndices &hi = request_to_fill[idx];
              libmesh_assert (hi <= node_upper_bounds[libMesh::processor_id()]);
              
              // find the requested index in my node bin
              std::vector<Hilbert::HilbertIndices>::const_iterator pos =
                 std::lower_bound (my_node_bin.begin(), my_node_bin.end(), hi);
              libmesh_assert (pos != my_node_bin.end());
              libmesh_assert (*pos == hi);
              
              // Finally, assign the global index based off the position of the index
              // in my array, properly offset.
              global_ids.push_back (std::distance(my_node_bin.begin(), pos) + my_offset);
            }

          // and trade back
          Parallel::send_receive (procdown, global_ids,
                                  procup,   filled_request[procup]);
        }

      // We now have all the filled requests, so we can loop through our
      // nodes once and assign the global index to each one.
      {
        std::vector<std::vector<unsigned int>::const_iterator>
          next_obj_on_proc; next_obj_on_proc.reserve(libMesh::n_processors());
        for (unsigned int pid=0; pid<libMesh::n_processors(); pid++)
          next_obj_on_proc.push_back(filled_request[pid].begin());

        MeshBase::node_iterator       it  = mesh.nodes_begin();
        const MeshBase::node_iterator end = mesh.nodes_end();

        for (; it != end; ++it)
          {
            Node* node = (*it);
            libmesh_assert (node != NULL);
            const Hilbert::HilbertIndices hi = 
              get_hilbert_index (*node, bbox);
            const unsigned int pid = 
              std::distance (node_upper_bounds.begin(), 
                             std::lower_bound(node_upper_bounds.begin(), 
                                              node_upper_bounds.end(),
                                              hi));

            libmesh_assert (pid < libMesh::n_processors());
            libmesh_assert (next_obj_on_proc[pid] != filled_request[pid].end());

            const unsigned int global_index = *next_obj_on_proc[pid];
            libmesh_assert (global_index < mesh.n_nodes());
            node->set_id() = global_index;
            
            ++next_obj_on_proc[pid];
          }
      }
    }

    //---------------------------------------------------
    // elements next -- all elements, not just local ones
    {
      // Request sets to send to each processor
      std::vector<std::vector<Hilbert::HilbertIndices> > 
        requested_ids (libMesh::n_processors());
      // Results to gather from each processor
      std::vector<std::vector<unsigned int> >
        filled_request (libMesh::n_processors());
      
      MeshBase::const_element_iterator       it  = mesh.elements_begin();
      const MeshBase::const_element_iterator end = mesh.elements_end();

      for (; it != end; ++it)
        {
          const Elem* elem = (*it);
          libmesh_assert (elem != NULL);
          const Hilbert::HilbertIndices hi = 
            get_hilbert_index (elem->centroid(), bbox);
          const unsigned int pid = 
            std::distance (elem_upper_bounds.begin(), 
                           std::lower_bound(elem_upper_bounds.begin(), 
                                            elem_upper_bounds.end(),
                                            hi));

          libmesh_assert (pid < libMesh::n_processors());

          requested_ids[pid].push_back(hi);
        }

      // The number of objects in my_elem_bin on each processor
      std::vector<unsigned int> elem_bin_sizes(libMesh::n_processors());
      Parallel::allgather (static_cast<unsigned int>(my_elem_bin.size()), elem_bin_sizes);

      // The offset of my first global index
      unsigned int my_offset = 0;
      for (unsigned int pid=0; pid<libMesh::processor_id(); pid++)
        my_offset += elem_bin_sizes[pid];

      // start with pid=0, so that we will trade with ourself
      for (unsigned int pid=0; pid<libMesh::n_processors(); pid++)
        {
          // Trade my requests with processor procup and procdown
          const unsigned int procup = (libMesh::processor_id() + pid) %
                                       libMesh::n_processors();
          const unsigned int procdown = (libMesh::n_processors() +
                                         libMesh::processor_id() - pid) %
                                         libMesh::n_processors();

          std::vector<Hilbert::HilbertIndices> request_to_fill;
          Parallel::send_receive(procup, requested_ids[procup],
                                 procdown, request_to_fill,
                                 hilbert_type);   

          // Fill the requests
          std::vector<unsigned int> global_ids;  global_ids.reserve(request_to_fill.size());
          for (unsigned int idx=0; idx<request_to_fill.size(); idx++)
            {
              const Hilbert::HilbertIndices &hi = request_to_fill[idx];
              libmesh_assert (hi <= elem_upper_bounds[libMesh::processor_id()]);
              
              // find the requested index in my elem bin
              std::vector<Hilbert::HilbertIndices>::const_iterator pos =
                std::lower_bound (my_elem_bin.begin(), my_elem_bin.end(), hi);
              libmesh_assert (pos != my_elem_bin.end());
              libmesh_assert (*pos == hi);
              
              // Finally, assign the global index based off the position of the index
              // in my array, properly offset.
              global_ids.push_back (std::distance(my_elem_bin.begin(), pos) + my_offset);
            }

          // and trade back
          Parallel::send_receive (procdown, global_ids,
                                  procup,   filled_request[procup]);
        }

      // We now have all the filled requests, so we can loop through our
      // elements once and assign the global index to each one.
      {
        std::vector<std::vector<unsigned int>::const_iterator>
          next_obj_on_proc; next_obj_on_proc.reserve(libMesh::n_processors());
        for (unsigned int pid=0; pid<libMesh::n_processors(); pid++)
          next_obj_on_proc.push_back(filled_request[pid].begin());

        MeshBase::element_iterator       it  = mesh.elements_begin();
        const MeshBase::element_iterator end = mesh.elements_end();

        for (; it != end; ++it)
          {
            Elem* elem = (*it);
            libmesh_assert (elem != NULL);
            const Hilbert::HilbertIndices hi = 
              get_hilbert_index (elem->centroid(), bbox);
            const unsigned int pid = 
              std::distance (elem_upper_bounds.begin(), 
                             std::lower_bound(elem_upper_bounds.begin(), 
                                              elem_upper_bounds.end(),
                                              hi));

            libmesh_assert (pid < libMesh::n_processors());
            libmesh_assert (next_obj_on_proc[pid] != filled_request[pid].end());

            const unsigned int global_index = *next_obj_on_proc[pid];
            libmesh_assert (global_index < mesh.n_elem());
            elem->set_id() = global_index;
            
            ++next_obj_on_proc[pid];
          }
      }
    }        
  }


  // Clean up
  MPI_Type_free (&hilbert_type);

  STOP_LOG ('assign_global_indices()', 'MeshCommunication');
}
 

void MeshCommunication::broadcast (MeshBase &mesh) constFinds all the processors that may contain elements that neighbor my elements. This list is guaranteed to include all processors that border any of my elements, but may include additional ones as well. This method computes bounding boxes for the elements on each processor and checks for overlaps. This method takes a mesh (which is assumed to reside on processor 0) and broadcasts it to all the other processors. It also broadcasts any boundary information the mesh has associated with it.

Definition at line 1333 of file mesh_communication.C.

References MeshBase::boundary_info, broadcast_bcs(), broadcast_mesh(), and libMesh::n_processors().

Referenced by UnstructuredMesh::read().

{
  // Don't need to do anything if there is
  // only one processor.
  if (libMesh::n_processors() == 1)
    return;
  
  this->broadcast_mesh (mesh);
  this->broadcast_bcs  (mesh, *(mesh.boundary_info));
}
 

void MeshCommunication::broadcast_bcs (const MeshBase &mesh, BoundaryInfo &boundary_info) const [private]

Definition at line 1596 of file mesh_communication.C.

Referenced by broadcast().

{
}
 

void MeshCommunication::broadcast_mesh (MeshBase &mesh) const [private]

Definition at line 1348 of file mesh_communication.C.

Referenced by broadcast().

{
  // no MPI, no need for this method...
  return;
}
 

void MeshCommunication::clear ()Clears all data structures and returns to a pristine state.

Definition at line 83 of file mesh_communication.C.

{
  //  _neighboring_processors.clear();
}
 

void MeshCommunication::delete_remote_elements (ParallelMesh &mesh) constThis method takes an input ParallelMesh which may be distributed among all the processors. Each processor deletes all elements which are neither local elements nor 'ghost' elements which touch local elements, and deletes all nodes which are not contained in local or ghost elements. The end result is that a previously serial ParallelMesh will be distributed between processors. Since this method is collective it must be called by all processors.

Definition at line 2412 of file mesh_communication.C.

References ParallelMesh::delete_elem(), DofObject::id(), ParallelMesh::is_serial(), ParallelMesh::level_elements_begin(), ParallelMesh::level_elements_end(), ParallelMesh::local_elements_begin(), ParallelMesh::local_elements_end(), Elem::make_links_to_me_remote(), ParallelMesh::max_elem_id(), ParallelMesh::max_node_id(), MeshTools::n_levels(), Elem::n_nodes(), Elem::node(), ParallelMesh::not_local_elements_begin(), ParallelMesh::not_local_elements_end(), Elem::parent(), ParallelMesh::renumber_nodes_and_elements(), ParallelMesh::unpartitioned_elements_begin(), and ParallelMesh::unpartitioned_elements_end().

Referenced by ParallelMesh::delete_remote_elements().

{
  // The mesh should know it's about to be parallelized
  libmesh_assert (!mesh.is_serial());

  START_LOG('delete_remote_elements()', 'MeshCommunication');

  // FIXME - should these be 'unsorted_set's?  O(N) is O(N)...
  std::vector<bool> local_nodes(mesh.max_node_id(), false);
  std::vector<bool> semilocal_elems(mesh.max_elem_id(), false);

  // We don't want to delete any element that shares a node
  // with or is an ancestor of a local element.
  MeshBase::const_element_iterator l_elem_it = mesh.local_elements_begin(),
                                   l_end     = mesh.local_elements_end();
  for (; l_elem_it != l_end; ++l_elem_it)
    {
      const Elem *elem = *l_elem_it;
      for (unsigned int n=0; n != elem->n_nodes(); ++n)
        local_nodes[elem->node(n)] = true;
      while (elem)
        {
          semilocal_elems[elem->id()] = true;
          elem = elem->parent();
        }
    }

  // We don't want to delete any element that shares a node
  // with or is an ancestor of an unpartitioned element either.
  MeshBase::const_element_iterator
    u_elem_it = mesh.unpartitioned_elements_begin(),
    u_end     = mesh.unpartitioned_elements_end();
  
  for (; u_elem_it != u_end; ++u_elem_it)
    {
      const Elem *elem = *u_elem_it;
      for (unsigned int n=0; n != elem->n_nodes(); ++n)
        local_nodes[elem->node(n)] = true;
      while (elem)
        {
          semilocal_elems[elem->id()] = true;
          elem = elem->parent();
        }
    }

  // Flag all the elements that share nodes with
  // local and unpartitioned elements, along with their ancestors
  MeshBase::element_iterator nl_elem_it = mesh.not_local_elements_begin(),
                             nl_end     = mesh.not_local_elements_end();
  for (; nl_elem_it != nl_end; ++nl_elem_it)
    {
      const Elem *elem = *nl_elem_it;
      for (unsigned int n=0; n != elem->n_nodes(); ++n)
        if (local_nodes[elem->node(n)])
          {
            while (elem)
              {
                semilocal_elems[elem->id()] = true;
                elem = elem->parent();
              }
            break;
          }
    }

  // Delete all the elements we have no reason to save,
  // starting with the most refined so that the mesh
  // is valid at all intermediate steps
  unsigned int n_levels = MeshTools::n_levels(mesh);

  for (int l = n_levels - 1; l >= 0; --l)
    {
      MeshBase::element_iterator lev_elem_it = mesh.level_elements_begin(l),
                                 lev_end     = mesh.level_elements_end(l);
      for (; lev_elem_it != lev_end; ++lev_elem_it)
        {
          Elem *elem = *lev_elem_it;
          libmesh_assert (elem);
          if (!semilocal_elems[elem->id()])
            {
              // Make sure we don't leave any invalid pointers
              elem->make_links_to_me_remote();

              // delete_elem doesn't currently invalidate element
              // iterators... that had better not change
              mesh.delete_elem(elem);
            }
        }
    }

  STOP_LOG('delete_remote_elements()', 'MeshCommunication');

  // Now make sure the containers actually shrink - strip
  // any newly-created NULL voids out of the element array
  mesh.renumber_nodes_and_elements();
}
 

template<typename ForwardIterator > void MeshCommunication::find_global_indices (const MeshTools::BoundingBox &bbox, const ForwardIterator &begin, const ForwardIterator &end, std::vector< unsigned int > &index_map) constThis method determines a globally unique, partition-agnostic index for each object in the input range.

Definition at line 524 of file mesh_communication_global_indices.C.

References allgather(), Parallel::Sort< KeyType >::bin(), libMesh::COMM_WORLD, MeshTools::Generation::Private::idx(), DofObject::invalid_processor_id, libMesh::n_processors(), libMesh::processor_id(), and Parallel::Sort< KeyType >::sort().

Referenced by MetisPartitioner::_do_partition(), ParmetisPartitioner::initialize(), and Partitioner::partition_unpartitioned_elements().

{
  START_LOG ('find_global_indices()', 'MeshCommunication');

  // This method determines partition-agnostic global indices
  // for nodes and elements.

  // Algorithm:
  // (1) compute the Hilbert key for each local node/element
  // (2) perform a parallel sort of the Hilbert key
  // (3) get the min/max value on each processor
  // (4) determine the position in the global ranking for
  //     each local object
  index_map.clear();
  index_map.reserve(std::distance (begin, end));

  // Set up a derived MPI datatype to handle communication of HilbertIndices
  MPI_Datatype hilbert_type;
  MPI_Type_contiguous (3, MPI_UNSIGNED, &hilbert_type);
  MPI_Type_commit     (&hilbert_type);


  //-------------------------------------------------------------
  // (1) compute Hilbert keys
  // These aren't trivial to compute, and we will need them again.
  // But the binsort will sort the input vector, trashing the order
  // that we'd like to rely on.  So, two vectors...
  std::vector<Hilbert::HilbertIndices> 
    sorted_hilbert_keys,
    hilbert_keys;
  sorted_hilbert_keys.reserve(index_map.capacity());
  hilbert_keys.reserve(index_map.capacity());
  {
    START_LOG('compute_hilbert_indices()', 'MeshCommunication');
    for (ForwardIterator it=begin; it!=end; ++it)
      {
        const Hilbert::HilbertIndices hi(get_hilbert_index (*it, bbox));
        hilbert_keys.push_back(hi);
        
        if ((*it)->processor_id() == libMesh::processor_id())
          sorted_hilbert_keys.push_back(hi);

        // someone needs to take care of unpartitioned objects!
        if ((libMesh::processor_id() == 0) &&
            ((*it)->processor_id() == DofObject::invalid_processor_id))
          sorted_hilbert_keys.push_back(hi);
      }    
    STOP_LOG('compute_hilbert_indices()', 'MeshCommunication');
  }
  
  //-------------------------------------------------------------
  // (2) parallel sort the Hilbert keys
  START_LOG ('parallel_sort()', 'MeshCommunication');  
  Parallel::Sort<Hilbert::HilbertIndices> sorter (sorted_hilbert_keys);
  sorter.sort(); 
  STOP_LOG ('parallel_sort()', 'MeshCommunication');
  const std::vector<Hilbert::HilbertIndices> &my_bin = sorter.bin();

  // The number of objects in my_bin on each processor
  std::vector<unsigned int> bin_sizes(libMesh::n_processors());
  Parallel::allgather (static_cast<unsigned int>(my_bin.size()), bin_sizes);
    
  // The offset of my first global index
  unsigned int my_offset = 0;
  for (unsigned int pid=0; pid<libMesh::processor_id(); pid++)
    my_offset += bin_sizes[pid];
  
  //-------------------------------------------------------------
  // (3) get the max value on each processor
  std::vector<Hilbert::HilbertIndices>    
    upper_bounds(libMesh::n_processors());
    
  // limit scope of temporaries
  {
    Hilbert::HilbertIndices my_max;
    
    if (!my_bin.empty()) my_max = my_bin.back();
    
    MPI_Allgather (&my_max,          1, hilbert_type,
                   &upper_bounds[0], 1, hilbert_type,
                   libMesh::COMM_WORLD);

    // Be cereful here.  The *_upper_bounds will be used to find the processor
    // a given object belongs to.  So, if a processor contains no objects (possible!)
    // then copy the bound from the lower processor id.
    for (unsigned int p=1; p<libMesh::n_processors(); p++)
      if (!bin_sizes[p]) upper_bounds[p] = upper_bounds[p-1];
  }


  //-------------------------------------------------------------
  // (4) determine the position in the global ranking for
  //     each local object
  {
    //----------------------------------------------
    // all objects, not just local ones
    
    // Request sets to send to each processor
    std::vector<std::vector<Hilbert::HilbertIndices> > 
      requested_ids (libMesh::n_processors());
    // Results to gather from each processor
    std::vector<std::vector<unsigned int> >
      filled_request (libMesh::n_processors());

    // build up list of requests
    std::vector<Hilbert::HilbertIndices>::const_iterator hi = 
      hilbert_keys.begin();

    for (ForwardIterator it = begin; it != end; ++it)
      {
        libmesh_assert (hi != hilbert_keys.end());
        const unsigned int pid = 
          std::distance (upper_bounds.begin(), 
                         std::lower_bound(upper_bounds.begin(), 
                                          upper_bounds.end(),
                                          *hi));

        libmesh_assert (pid < libMesh::n_processors());

        requested_ids[pid].push_back(*hi);

        hi++;
        // go ahead and put pid in index_map, that way we 
        // don't have to repeat the std::lower_bound()
        index_map.push_back(pid);
      }

    // start with pid=0, so that we will trade with ourself
    std::vector<Hilbert::HilbertIndices> request_to_fill;
    std::vector<unsigned int> global_ids;
    for (unsigned int pid=0; pid<libMesh::n_processors(); pid++)
      {
        // Trade my requests with processor procup and procdown
        const unsigned int procup = (libMesh::processor_id() + pid) %
                                     libMesh::n_processors();
        const unsigned int procdown = (libMesh::n_processors() +
                                       libMesh::processor_id() - pid) %
                                       libMesh::n_processors();

        Parallel::send_receive(procup, requested_ids[procup],
                               procdown, request_to_fill,
                               hilbert_type);     

        // Fill the requests
        global_ids.clear();  global_ids.reserve(request_to_fill.size());
        for (unsigned int idx=0; idx<request_to_fill.size(); idx++)
          {
            const Hilbert::HilbertIndices &hi = request_to_fill[idx];
            libmesh_assert (hi <= upper_bounds[libMesh::processor_id()]);
            
            // find the requested index in my node bin
            std::vector<Hilbert::HilbertIndices>::const_iterator pos =
              std::lower_bound (my_bin.begin(), my_bin.end(), hi);
            libmesh_assert (pos != my_bin.end());
            libmesh_assert (*pos == hi);
            
            // Finally, assign the global index based off the position of the index
            // in my array, properly offset.
            global_ids.push_back (std::distance(my_bin.begin(), pos) + my_offset);
          }
        
        // and trade back
        Parallel::send_receive (procdown, global_ids,
                                procup,   filled_request[procup]);
      }

    // We now have all the filled requests, so we can loop through our
    // nodes once and assign the global index to each one.
    {
      std::vector<std::vector<unsigned int>::const_iterator>
        next_obj_on_proc; next_obj_on_proc.reserve(libMesh::n_processors());
      for (unsigned int pid=0; pid<libMesh::n_processors(); pid++)
        next_obj_on_proc.push_back(filled_request[pid].begin());
      
      unsigned int cnt=0;
      for (ForwardIterator it = begin; it != end; ++it, cnt++)
        {
          const unsigned int pid = index_map[cnt];
          
          libmesh_assert (pid < libMesh::n_processors());
          libmesh_assert (next_obj_on_proc[pid] != filled_request[pid].end());
          
          const unsigned int global_index = *next_obj_on_proc[pid];
          index_map[cnt] = global_index;

          ++next_obj_on_proc[pid];
        }
    }
  }


  // Clean up
  MPI_Type_free (&hilbert_type);

  STOP_LOG ('find_global_indices()', 'MeshCommunication');
}
 

void MeshCommunication::gather_neighboring_elements (ParallelMesh &mesh) const

Definition at line 673 of file mesh_communication.C.

{
  // no MPI == one processor, no need for this method...
  return;
}
 

void MeshCommunication::make_elems_parallel_consistent (MeshBase &mesh)Copy ids of ghost elements from their local processors.

Definition at line 2266 of file mesh_communication.C.

References MeshBase::active_elements_begin(), MeshBase::active_elements_end(), MeshBase::renumber_elem(), and Parallel::sync_element_data_by_parent_id().

Referenced by MeshRefinement::_refine_elements().

{
  // This function must be run on all processors at once
  parallel_only();

  START_LOG ('make_elems_parallel_consistent()', 'MeshCommunication');

  SyncIds syncids(mesh, &MeshBase::renumber_elem);
  Parallel::sync_element_data_by_parent_id
    (mesh, mesh.active_elements_begin(),
     mesh.active_elements_end(), syncids);

  STOP_LOG ('make_elems_parallel_consistent()', 'MeshCommunication');
}
 

void MeshCommunication::make_node_ids_parallel_consistent (MeshBase &mesh, LocationMap< Node > &loc_map)Assuming all ids on local nodes are globally unique, and assuming all processor ids are parallel consistent, this function makes all other ids parallel consistent.

Definition at line 2247 of file mesh_communication.C.

References MeshBase::nodes_begin(), MeshBase::nodes_end(), MeshBase::renumber_node(), and Parallel::sync_dofobject_data_by_xyz().

{
  // This function must be run on all processors at once
  parallel_only();

  START_LOG ('make_node_ids_parallel_consistent()', 'MeshCommunication');

  SyncIds syncids(mesh, &MeshBase::renumber_node);
  Parallel::sync_dofobject_data_by_xyz
    (mesh.nodes_begin(), mesh.nodes_end(),
     loc_map, syncids);

  STOP_LOG ('make_node_ids_parallel_consistent()', 'MeshCommunication');
}
 

void MeshCommunication::make_node_proc_ids_parallel_consistent (MeshBase &mesh, LocationMap< Node > &loc_map)Assuming all processor ids on nodes touching local elements are parallel consistent, this function makes all other processor ids parallel consistent as well.

Definition at line 2332 of file mesh_communication.C.

References MeshBase::nodes_begin(), MeshBase::nodes_end(), and Parallel::sync_dofobject_data_by_xyz().

Referenced by MeshTools::correct_node_proc_ids().

{
  START_LOG ('make_node_proc_ids_parallel_consistent()', 'MeshCommunication');

  // This function must be run on all processors at once
  parallel_only();

  // When this function is called, each section of a parallelized mesh
  // should be in the following state:
  //
  // All nodes should have the exact same physical location on every
  // processor where they exist.
  //
  // Local nodes should have unique authoritative ids,
  // and processor ids consistent with all processors which own
  // an element touching them.
  //
  // Ghost nodes touching local elements should have processor ids
  // consistent with all processors which own an element touching
  // them.
  
  SyncProcIds sync(mesh);
  Parallel::sync_dofobject_data_by_xyz
    (mesh.nodes_begin(), mesh.nodes_end(), loc_map, sync);

  STOP_LOG ('make_node_proc_ids_parallel_consistent()', 'MeshCommunication');
}
 

void MeshCommunication::make_nodes_parallel_consistent (MeshBase &mesh, LocationMap< Node > &loc_map)Copy processor_ids and ids on ghost nodes from their local processors. This is an internal function of MeshRefinement which turns out to be useful for other code which wants to add nodes to a distributed mesh.

Definition at line 2365 of file mesh_communication.C.

References MeshTools::correct_node_proc_ids(), LocationMap< T >::empty(), LocationMap< T >::init(), std::max(), MeshBase::nodes_begin(), and MeshBase::nodes_end().

Referenced by MeshRefinement::_coarsen_elements(), MeshRefinement::_refine_elements(), and UnstructuredMesh::all_second_order().

{
  // This function must be run on all processors at once
  parallel_only();

  // Create the loc_map if it hasn't been done already
  bool need_map_update = (mesh.nodes_begin() != mesh.nodes_end() && 
                          loc_map.empty());
  Parallel::max(need_map_update);

  if (need_map_update)
    loc_map.init(mesh);

  // When this function is called, each section of a parallelized mesh
  // should be in the following state:
  //
  // All nodes should have the exact same physical location on every
  // processor where they exist.
  //
  // Local nodes should have unique authoritative ids,
  // and processor ids consistent with all processors which own
  // an element touching them.
  //
  // Ghost nodes touching local elements should have processor ids
  // consistent with all processors which own an element touching
  // them.
  //
  // Ghost nodes should have ids which are either already correct
  // or which are in the 'unpartitioned' id space.

  // First, let's sync up processor ids.  Some of these processor ids
  // may be 'wrong' from coarsening, but they're right in the sense
  // that they'll tell us who has the authoritative dofobject ids for
  // each node.
  this->make_node_proc_ids_parallel_consistent(mesh, loc_map);

  // Second, sync up dofobject ids.
  this->make_node_ids_parallel_consistent(mesh, loc_map);

  // Finally, correct the processor ids to make DofMap happy
  MeshTools::correct_node_proc_ids(mesh, loc_map);
}
 

void MeshCommunication::redistribute (ParallelMesh &mesh) constThis method takes a parallel distributed mesh and redistributes the elements. Specifically, any elements stored on a given processor are sent to the processor which 'owns' them. Similarly, any elements assigned to the current processor but stored on another are received. Once this step is completed any required ghost elements are updated. The final result is that each processor stores only the elements it actually owns and any ghost elements required to satisfy data dependencies. This method can be invoked after a partitioning step to affect the new partitioning.

Definition at line 154 of file mesh_communication.C.

{
  // no MPI == one processor, no redistribution
  return;
}

 

Author

Generated automatically by Doxygen for libMesh from the source code.


 

Index

NAME
SYNOPSIS
Public Member Functions
Private Member Functions
Detailed Description
Constructor & Destructor Documentation
MeshCommunication::MeshCommunication () [inline]Constructor.
MeshCommunication::~MeshCommunication () [inline]Destructor.
Member Function Documentation
void MeshCommunication::allgather (ParallelMesh &mesh) constThis method takes an input ParallelMesh which may be distributed among all the processors. Each processor then sends its local nodes and elements to the other processors. The end result is that a previously distributed ParallelMesh will be serialized on each processor. Since this method is collective it must be called by all processors.
void MeshCommunication::allgather_bcs (const ParallelMesh &mesh, BoundaryInfo &boundary_info) const [private]
void MeshCommunication::allgather_mesh (ParallelMesh &mesh) const [private]
void MeshCommunication::assign_global_indices (MeshBase &mesh) constThis method assigns globally unique, partition-agnostic indices to the nodes and elements in the mesh. The approach is to compute the Hilbert space-filling curve key and use its value to assign an index in [0,N_global). Since the Hilbert key is unique for each spatial location, two objects occupying the same location will be assigned the same global id. Thus, this method can also be useful for identifying duplicate nodes which may occur during parallel refinement.
void MeshCommunication::broadcast (MeshBase &mesh) constFinds all the processors that may contain elements that neighbor my elements. This list is guaranteed to include all processors that border any of my elements, but may include additional ones as well. This method computes bounding boxes for the elements on each processor and checks for overlaps. This method takes a mesh (which is assumed to reside on processor 0) and broadcasts it to all the other processors. It also broadcasts any boundary information the mesh has associated with it.
void MeshCommunication::broadcast_bcs (const MeshBase &mesh, BoundaryInfo &boundary_info) const [private]
void MeshCommunication::broadcast_mesh (MeshBase &mesh) const [private]
void MeshCommunication::clear ()Clears all data structures and returns to a pristine state.
void MeshCommunication::delete_remote_elements (ParallelMesh &mesh) constThis method takes an input ParallelMesh which may be distributed among all the processors. Each processor deletes all elements which are neither local elements nor 'ghost' elements which touch local elements, and deletes all nodes which are not contained in local or ghost elements. The end result is that a previously serial ParallelMesh will be distributed between processors. Since this method is collective it must be called by all processors.
template<typename ForwardIterator > void MeshCommunication::find_global_indices (const MeshTools::BoundingBox &bbox, const ForwardIterator &begin, const ForwardIterator &end, std::vector< unsigned int > &index_map) constThis method determines a globally unique, partition-agnostic index for each object in the input range.
void MeshCommunication::gather_neighboring_elements (ParallelMesh &mesh) const
void MeshCommunication::make_elems_parallel_consistent (MeshBase &mesh)Copy ids of ghost elements from their local processors.
void MeshCommunication::make_node_ids_parallel_consistent (MeshBase &mesh, LocationMap< Node > &loc_map)Assuming all ids on local nodes are globally unique, and assuming all processor ids are parallel consistent, this function makes all other ids parallel consistent.
void MeshCommunication::make_node_proc_ids_parallel_consistent (MeshBase &mesh, LocationMap< Node > &loc_map)Assuming all processor ids on nodes touching local elements are parallel consistent, this function makes all other processor ids parallel consistent as well.
void MeshCommunication::make_nodes_parallel_consistent (MeshBase &mesh, LocationMap< Node > &loc_map)Copy processor_ids and ids on ghost nodes from their local processors. This is an internal function of MeshRefinement which turns out to be useful for other code which wants to add nodes to a distributed mesh.
void MeshCommunication::redistribute (ParallelMesh &mesh) constThis method takes a parallel distributed mesh and redistributes the elements. Specifically, any elements stored on a given processor are sent to the processor which 'owns' them. Similarly, any elements assigned to the current processor but stored on another are received. Once this step is completed any required ghost elements are updated. The final result is that each processor stores only the elements it actually owns and any ghost elements required to satisfy data dependencies. This method can be invoked after a partitioning step to affect the new partitioning.
Author

This document was created by man2html, using the manual pages.
Time: 21:50:10 GMT, April 16, 2011