libMesh::ParmetisPartitioner Class Reference

Partitioner which provides an interface to ParMETIS. More...

#include <parmetis_partitioner.h>

Inheritance diagram for libMesh::ParmetisPartitioner:

Public Member Functions

 ParmetisPartitioner ()
 
 ~ParmetisPartitioner ()
 
virtual UniquePtr< Partitionerclone () const libmesh_override
 
virtual void partition (MeshBase &mesh, const unsigned int n)
 
virtual void partition (MeshBase &mesh)
 
virtual void partition_range (MeshBase &, MeshBase::element_iterator, MeshBase::element_iterator, const unsigned int)
 
void repartition (MeshBase &mesh, const unsigned int n)
 
void repartition (MeshBase &mesh)
 
virtual void attach_weights (ErrorVector *)
 

Static Public Member Functions

static void partition_unpartitioned_elements (MeshBase &mesh)
 
static void partition_unpartitioned_elements (MeshBase &mesh, const unsigned int n)
 
static void set_parent_processor_ids (MeshBase &mesh)
 
static void set_node_processor_ids (MeshBase &mesh)
 

Protected Member Functions

virtual void _do_repartition (MeshBase &mesh, const unsigned int n) libmesh_override
 
virtual void _do_partition (MeshBase &mesh, const unsigned int n) libmesh_override
 
void single_partition (MeshBase &mesh)
 
void single_partition_range (MeshBase::element_iterator it, MeshBase::element_iterator end)
 

Protected Attributes

ErrorVector_weights
 

Static Protected Attributes

static const dof_id_type communication_blocksize = 1000000
 

Private Member Functions

void initialize (const MeshBase &mesh, const unsigned int n_sbdmns)
 
void build_graph (const MeshBase &mesh)
 
void assign_partitioning (MeshBase &mesh)
 

Private Attributes

std::vector< dof_id_type_n_active_elem_on_proc
 
vectormap< dof_id_type, dof_id_type_global_index_by_pid_map
 
UniquePtr< ParmetisHelper_pmetis
 

Detailed Description

Partitioner which provides an interface to ParMETIS.

The ParmetisPartitioner uses the Parmetis graph partitioner to partition the elements.

Author
Benjamin S. Kirk
Date
2003

Definition at line 46 of file parmetis_partitioner.h.

Constructor & Destructor Documentation

libMesh::ParmetisPartitioner::ParmetisPartitioner ( )

Constructor.

Definition at line 64 of file parmetis_partitioner.C.

Referenced by clone().

66  : _pmetis(new ParmetisHelper)
67 #endif
68 {}
69 
70 
71 
73 {
74 }
UniquePtr< ParmetisHelper > _pmetis
libMesh::ParmetisPartitioner::~ParmetisPartitioner ( )

Destructor.

Member Function Documentation

void libMesh::ParmetisPartitioner::_do_partition ( MeshBase mesh,
const unsigned int  n 
)
protectedvirtual

Partition the MeshBase into n subdomains.

Implements libMesh::Partitioner.

Definition at line 78 of file parmetis_partitioner.C.

Referenced by clone().

80 {
81  this->_do_repartition (mesh, n_sbdmns);
82 }
virtual void _do_repartition(MeshBase &mesh, const unsigned int n) libmesh_override
MeshBase & mesh
void libMesh::ParmetisPartitioner::_do_repartition ( MeshBase mesh,
const unsigned int  n 
)
protectedvirtual

Parmetis can handle dynamically repartitioning a mesh such that the redistribution costs are minimized. This method takes a previously partitioned mesh (which may have then been adaptively refined) and repartitions it.

Reimplemented from libMesh::Partitioner.

Definition at line 86 of file parmetis_partitioner.C.

References libMesh::MeshBase::allgather(), libMesh::ParallelObject::comm(), libMesh::err, libMesh::Parallel::Communicator::get(), libmesh_nullptr, libMesh::ParallelObject::n_processors(), and libMesh::Partitioner::partition().

Referenced by clone().

88 {
89  libmesh_assert_greater (n_sbdmns, 0);
90 
91  // Check for an easy return
92  if (n_sbdmns == 1)
93  {
94  this->single_partition(mesh);
95  return;
96  }
97 
98  // This function must be run on all processors at once
99  libmesh_parallel_only(mesh.comm());
100 
101  // What to do if the Parmetis library IS NOT present
102 #ifndef LIBMESH_HAVE_PARMETIS
103 
104  libmesh_here();
105  libMesh::err << "ERROR: The library has been built without" << std::endl
106  << "Parmetis support. Using a Metis" << std::endl
107  << "partitioner instead!" << std::endl;
108 
109  MetisPartitioner mp;
110 
111  mp.partition (mesh, n_sbdmns);
112 
113  // What to do if the Parmetis library IS present
114 #else
115 
116  // Revert to METIS on one processor.
117  if (mesh.n_processors() == 1)
118  {
119  // Make sure the mesh knows it's serial
120  mesh.allgather();
121 
122  MetisPartitioner mp;
123  mp.partition (mesh, n_sbdmns);
124  return;
125  }
126 
127  LOG_SCOPE("repartition()", "ParmetisPartitioner");
128 
129  // Initialize the data structures required by ParMETIS
130  this->initialize (mesh, n_sbdmns);
131 
132  // Make sure all processors have enough active local elements.
133  // Parmetis tends to crash when it's given only a couple elements
134  // per partition.
135  {
136  bool all_have_enough_elements = true;
137  for (std::size_t pid=0; pid<_n_active_elem_on_proc.size(); pid++)
139  all_have_enough_elements = false;
140 
141  // Parmetis will not work unless each processor has some
142  // elements. Specifically, it will abort when passed a NULL
143  // partition array on *any* of the processors.
144  if (!all_have_enough_elements)
145  {
146  // FIXME: revert to METIS, although this requires a serial mesh
147  MeshSerializer serialize(mesh);
148  MetisPartitioner mp;
149  mp.partition (mesh, n_sbdmns);
150  return;
151  }
152  }
153 
154  // build the graph corresponding to the mesh
155  this->build_graph (mesh);
156 
157 
158  // Partition the graph
159  std::vector<Parmetis::idx_t> vsize(_pmetis->vwgt.size(), 1);
160  Parmetis::real_t itr = 1000000.0;
161  MPI_Comm mpi_comm = mesh.comm().get();
162 
163  // Call the ParMETIS adaptive repartitioning method. This respects the
164  // original partitioning when computing the new partitioning so as to
165  // minimize the required data redistribution.
166  Parmetis::ParMETIS_V3_AdaptiveRepart(_pmetis->vtxdist.empty() ? libmesh_nullptr : &_pmetis->vtxdist[0],
167  _pmetis->xadj.empty() ? libmesh_nullptr : &_pmetis->xadj[0],
168  _pmetis->adjncy.empty() ? libmesh_nullptr : &_pmetis->adjncy[0],
169  _pmetis->vwgt.empty() ? libmesh_nullptr : &_pmetis->vwgt[0],
170  vsize.empty() ? libmesh_nullptr : &vsize[0],
172  &_pmetis->wgtflag,
173  &_pmetis->numflag,
174  &_pmetis->ncon,
175  &_pmetis->nparts,
176  _pmetis->tpwgts.empty() ? libmesh_nullptr : &_pmetis->tpwgts[0],
177  _pmetis->ubvec.empty() ? libmesh_nullptr : &_pmetis->ubvec[0],
178  &itr,
179  &_pmetis->options[0],
180  &_pmetis->edgecut,
181  _pmetis->part.empty() ? libmesh_nullptr : &_pmetis->part[0],
182  &mpi_comm);
183 
184  // Assign the returned processor ids
185  this->assign_partitioning (mesh);
186 
187 #endif // #ifndef LIBMESH_HAVE_PARMETIS ... else ...
188 
189 }
void single_partition(MeshBase &mesh)
Definition: partitioner.C:151
void initialize(const MeshBase &mesh, const unsigned int n_sbdmns)
MeshBase & mesh
const class libmesh_nullptr_t libmesh_nullptr
std::vector< dof_id_type > _n_active_elem_on_proc
UniquePtr< ParmetisHelper > _pmetis
void assign_partitioning(MeshBase &mesh)
void build_graph(const MeshBase &mesh)
OStreamProxy err(std::cerr)
const unsigned int MIN_ELEM_PER_PROC
void libMesh::ParmetisPartitioner::assign_partitioning ( MeshBase mesh)
private

Assign the computed partitioning to the mesh.

Definition at line 634 of file parmetis_partitioner.C.

References libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::ParallelObject::comm(), libMesh::DofObject::id(), libMesh::libmesh_assert(), libMesh::MeshBase::n_active_local_elem(), libMesh::ParallelObject::n_processors(), libMesh::ParallelObject::processor_id(), libMesh::DofObject::processor_id(), and libMesh::Parallel::Communicator::send_receive().

Referenced by clone().

635 {
636  LOG_SCOPE("assign_partitioning()", "ParmetisPartitioner");
637 
638  // This function must be run on all processors at once
639  libmesh_parallel_only(mesh.comm());
640 
641  const dof_id_type
642  first_local_elem = _pmetis->vtxdist[mesh.processor_id()];
643 
644 #ifndef NDEBUG
645  const dof_id_type n_active_local_elem = mesh.n_active_local_elem();
646 #endif
647 
648  std::vector<std::vector<dof_id_type> >
649  requested_ids(mesh.n_processors()),
650  requests_to_fill(mesh.n_processors());
651 
652  MeshBase::element_iterator elem_it = mesh.active_elements_begin();
653  MeshBase::element_iterator elem_end = mesh.active_elements_end();
654 
655  for (; elem_it != elem_end; ++elem_it)
656  {
657  Elem * elem = *elem_it;
658 
659  // we need to get the index from the owning processor
660  // (note we cannot assign it now -- we are iterating
661  // over elements again and this will be bad!)
662  libmesh_assert_less (elem->processor_id(), requested_ids.size());
663  requested_ids[elem->processor_id()].push_back(elem->id());
664  }
665 
666  // Trade with all processors (including self) to get their indices
667  for (processor_id_type pid=0; pid<mesh.n_processors(); pid++)
668  {
669  // Trade my requests with processor procup and procdown
670  const processor_id_type procup = (mesh.processor_id() + pid) % mesh.n_processors();
671  const processor_id_type procdown = (mesh.n_processors() +
672  mesh.processor_id() - pid) % mesh.n_processors();
673 
674  mesh.comm().send_receive (procup, requested_ids[procup],
675  procdown, requests_to_fill[procdown]);
676 
677  // we can overwrite these requested ids in-place.
678  for (std::size_t i=0; i<requests_to_fill[procdown].size(); i++)
679  {
680  const dof_id_type requested_elem_index =
681  requests_to_fill[procdown][i];
682 
683  libmesh_assert(_global_index_by_pid_map.count(requested_elem_index));
684 
685  const dof_id_type global_index_by_pid =
686  _global_index_by_pid_map[requested_elem_index];
687 
688  const dof_id_type local_index =
689  global_index_by_pid - first_local_elem;
690 
691  libmesh_assert_less (local_index, _pmetis->part.size());
692  libmesh_assert_less (local_index, n_active_local_elem);
693 
694  const unsigned int elem_procid =
695  static_cast<unsigned int>(_pmetis->part[local_index]);
696 
697  libmesh_assert_less (elem_procid, static_cast<unsigned int>(_pmetis->nparts));
698 
699  requests_to_fill[procdown][i] = elem_procid;
700  }
701 
702  // Trade back
703  mesh.comm().send_receive (procdown, requests_to_fill[procdown],
704  procup, requested_ids[procup]);
705  }
706 
707  // and finally assign the partitioning.
708  // note we are iterating in exactly the same order
709  // used to build up the request, so we can expect the
710  // required entries to be in the proper sequence.
711  elem_it = mesh.active_elements_begin();
712  elem_end = mesh.active_elements_end();
713 
714  for (std::vector<unsigned int> counters(mesh.n_processors(), 0);
715  elem_it != elem_end; ++elem_it)
716  {
717  Elem * elem = *elem_it;
718 
719  const processor_id_type current_pid = elem->processor_id();
720 
721  libmesh_assert_less (counters[current_pid], requested_ids[current_pid].size());
722 
723  const processor_id_type elem_procid =
724  requested_ids[current_pid][counters[current_pid]++];
725 
726  libmesh_assert_less (elem_procid, static_cast<unsigned int>(_pmetis->nparts));
727  elem->processor_id() = elem_procid;
728  }
729 }
difference_type count(const key_type &key) const
Definition: vectormap.h:210
MeshBase & mesh
uint8_t processor_id_type
Definition: id_types.h:99
vectormap< dof_id_type, dof_id_type > _global_index_by_pid_map
libmesh_assert(j)
UniquePtr< ParmetisHelper > _pmetis
uint8_t dof_id_type
Definition: id_types.h:64
virtual void libMesh::Partitioner::attach_weights ( ErrorVector )
inlinevirtualinherited

Attach weights that can be used for partitioning. This ErrorVector should be exactly the same on every processor and should have mesh->max_elem_id() entries.

Reimplemented in libMesh::MetisPartitioner.

Definition at line 170 of file partitioner.h.

References libMesh::Partitioner::_do_partition(), end, libMesh::Partitioner::single_partition(), and libMesh::Partitioner::single_partition_range().

170 { libmesh_not_implemented(); }
void libMesh::ParmetisPartitioner::build_graph ( const MeshBase mesh)
private

Build the graph.

Definition at line 410 of file parmetis_partitioner.C.

References libMesh::Elem::active(), libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::Elem::active_family_tree(), libMesh::MeshBase::active_local_elements_begin(), libMesh::MeshBase::active_local_elements_end(), libMesh::Elem::dim(), libMesh::Elem::find_interior_neighbors(), libMesh::DofObject::id(), libMesh::Elem::interior_parent(), libMesh::libmesh_assert(), libmesh_nullptr, libMesh::MeshBase::n_active_local_elem(), libMesh::Elem::n_neighbors(), libMesh::Elem::neighbor_ptr(), libMesh::ParallelObject::processor_id(), and libMesh::Elem::which_neighbor_am_i().

Referenced by clone().

411 {
412  LOG_SCOPE("build_graph()", "ParmetisPartitioner");
413 
414  // build the graph in distributed CSR format. Note that
415  // the edges in the graph will correspond to
416  // face neighbors
417  const dof_id_type n_active_local_elem = mesh.n_active_local_elem();
418 
419  // If we have boundary elements in this mesh, we want to account for
420  // the connectivity between them and interior elements. We can find
421  // interior elements from boundary elements, but we need to build up
422  // a lookup map to do the reverse.
423 
424  typedef LIBMESH_BEST_UNORDERED_MULTIMAP<const Elem *, const Elem *>
425  map_type;
426  map_type interior_to_boundary_map;
427 
428  {
429  MeshBase::const_element_iterator elem_it = mesh.active_elements_begin();
430  const MeshBase::const_element_iterator elem_end = mesh.active_elements_end();
431 
432  for (; elem_it != elem_end; ++elem_it)
433  {
434  const Elem * elem = *elem_it;
435 
436  // If we don't have an interior_parent then there's nothing to look us
437  // up.
438  if ((elem->dim() >= LIBMESH_DIM) ||
439  !elem->interior_parent())
440  continue;
441 
442  // get all relevant interior elements
443  std::set<const Elem *> neighbor_set;
444  elem->find_interior_neighbors(neighbor_set);
445 
446  std::set<const Elem *>::iterator n_it = neighbor_set.begin();
447  for (; n_it != neighbor_set.end(); ++n_it)
448  {
449  // FIXME - non-const versions of the Elem set methods
450  // would be nice
451  Elem * neighbor = const_cast<Elem *>(*n_it);
452 
453 #if defined(LIBMESH_HAVE_UNORDERED_MULTIMAP) || \
454  defined(LIBMESH_HAVE_TR1_UNORDERED_MAP) || \
455  defined(LIBMESH_HAVE_HASH_MAP) || \
456  defined(LIBMESH_HAVE_EXT_HASH_MAP)
457  interior_to_boundary_map.insert
458  (std::make_pair(neighbor, elem));
459 #else
460  interior_to_boundary_map.insert
461  (interior_to_boundary_map.begin(),
462  std::make_pair(neighbor, elem));
463 #endif
464  }
465  }
466  }
467 
468 #ifdef LIBMESH_ENABLE_AMR
469  std::vector<const Elem *> neighbors_offspring;
470 #endif
471 
472  std::vector<std::vector<dof_id_type> > graph(n_active_local_elem);
473  dof_id_type graph_size=0;
474 
475  const dof_id_type first_local_elem = _pmetis->vtxdist[mesh.processor_id()];
476 
477  MeshBase::const_element_iterator elem_it = mesh.active_local_elements_begin();
478  const MeshBase::const_element_iterator elem_end = mesh.active_local_elements_end();
479 
480  for (; elem_it != elem_end; ++elem_it)
481  {
482  const Elem * elem = *elem_it;
483 
485  const dof_id_type global_index_by_pid =
486  _global_index_by_pid_map[elem->id()];
487 
488  const dof_id_type local_index =
489  global_index_by_pid - first_local_elem;
490  libmesh_assert_less (local_index, n_active_local_elem);
491 
492  std::vector<dof_id_type> & graph_row = graph[local_index];
493 
494  // Loop over the element's neighbors. An element
495  // adjacency corresponds to a face neighbor
496  for (unsigned int ms=0; ms<elem->n_neighbors(); ms++)
497  {
498  const Elem * neighbor = elem->neighbor_ptr(ms);
499 
500  if (neighbor != libmesh_nullptr)
501  {
502  // If the neighbor is active treat it
503  // as a connection
504  if (neighbor->active())
505  {
507  const dof_id_type neighbor_global_index_by_pid =
508  _global_index_by_pid_map[neighbor->id()];
509 
510  graph_row.push_back(neighbor_global_index_by_pid);
511  graph_size++;
512  }
513 
514 #ifdef LIBMESH_ENABLE_AMR
515 
516  // Otherwise we need to find all of the
517  // neighbor's children that are connected to
518  // us and add them
519  else
520  {
521  // The side of the neighbor to which
522  // we are connected
523  const unsigned int ns =
524  neighbor->which_neighbor_am_i (elem);
525  libmesh_assert_less (ns, neighbor->n_neighbors());
526 
527  // Get all the active children (& grandchildren, etc...)
528  // of the neighbor
529 
530  // FIXME - this is the wrong thing, since we
531  // should be getting the active family tree on
532  // our side only. But adding too many graph
533  // links may cause hanging nodes to tend to be
534  // on partition interiors, which would reduce
535  // communication overhead for constraint
536  // equations, so we'll leave it.
537 
538  neighbor->active_family_tree (neighbors_offspring);
539 
540  // Get all the neighbor's children that
541  // live on that side and are thus connected
542  // to us
543  for (std::size_t nc=0; nc<neighbors_offspring.size(); nc++)
544  {
545  const Elem * child =
546  neighbors_offspring[nc];
547 
548  // This does not assume a level-1 mesh.
549  // Note that since children have sides numbered
550  // coincident with the parent then this is a sufficient test.
551  if (child->neighbor_ptr(ns) == elem)
552  {
553  libmesh_assert (child->active());
555  const dof_id_type child_global_index_by_pid =
556  _global_index_by_pid_map[child->id()];
557 
558  graph_row.push_back(child_global_index_by_pid);
559  graph_size++;
560  }
561  }
562  }
563 
564 #endif /* ifdef LIBMESH_ENABLE_AMR */
565 
566 
567  }
568  }
569 
570  if ((elem->dim() < LIBMESH_DIM) &&
571  elem->interior_parent())
572  {
573  // get all relevant interior elements
574  std::set<const Elem *> neighbor_set;
575  elem->find_interior_neighbors(neighbor_set);
576 
577  std::set<const Elem *>::iterator n_it = neighbor_set.begin();
578  for (; n_it != neighbor_set.end(); ++n_it)
579  {
580  // FIXME - non-const versions of the Elem set methods
581  // would be nice
582  Elem * neighbor = const_cast<Elem *>(*n_it);
583 
584  const dof_id_type neighbor_global_index_by_pid =
585  _global_index_by_pid_map[neighbor->id()];
586 
587  graph_row.push_back(neighbor_global_index_by_pid);
588  graph_size++;
589  }
590  }
591 
592  // Check for any boundary neighbors
593  typedef map_type::iterator map_it_type;
594  std::pair<map_it_type, map_it_type>
595  bounds = interior_to_boundary_map.equal_range(elem);
596 
597  for (map_it_type it = bounds.first; it != bounds.second; ++it)
598  {
599  const Elem * neighbor = it->second;
600 
601  const dof_id_type neighbor_global_index_by_pid =
602  _global_index_by_pid_map[neighbor->id()];
603 
604  graph_row.push_back(neighbor_global_index_by_pid);
605  graph_size++;
606  }
607  }
608 
609  // Reserve space in the adjacency array
610  _pmetis->xadj.clear();
611  _pmetis->xadj.reserve (n_active_local_elem + 1);
612  _pmetis->adjncy.clear();
613  _pmetis->adjncy.reserve (graph_size);
614 
615  for (std::size_t r=0; r<graph.size(); r++)
616  {
617  _pmetis->xadj.push_back(_pmetis->adjncy.size());
618  std::vector<dof_id_type> graph_row; // build this empty
619  graph_row.swap(graph[r]); // this will deallocate at the end of scope
620  _pmetis->adjncy.insert(_pmetis->adjncy.end(),
621  graph_row.begin(),
622  graph_row.end());
623  }
624 
625  // The end of the adjacency array for the last elem
626  _pmetis->xadj.push_back(_pmetis->adjncy.size());
627 
628  libmesh_assert_equal_to (_pmetis->xadj.size(), n_active_local_elem+1);
629  libmesh_assert_equal_to (_pmetis->adjncy.size(), graph_size);
630 }
difference_type count(const key_type &key) const
Definition: vectormap.h:210
MeshBase & mesh
const class libmesh_nullptr_t libmesh_nullptr
vectormap< dof_id_type, dof_id_type > _global_index_by_pid_map
libmesh_assert(j)
UniquePtr< ParmetisHelper > _pmetis
uint8_t dof_id_type
Definition: id_types.h:64
virtual UniquePtr<Partitioner> libMesh::ParmetisPartitioner::clone ( ) const
inlinevirtual
Returns
A copy of this partitioner wrapped in a smart pointer.

Implements libMesh::Partitioner.

Definition at line 63 of file parmetis_partitioner.h.

References _do_partition(), _do_repartition(), assign_partitioning(), build_graph(), initialize(), mesh, and ParmetisPartitioner().

64  {
65  return UniquePtr<Partitioner>(new ParmetisPartitioner());
66  }
void libMesh::ParmetisPartitioner::initialize ( const MeshBase mesh,
const unsigned int  n_sbdmns 
)
private

Initialize data structures.

Definition at line 196 of file parmetis_partitioner.C.

References libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::MeshBase::active_local_elements_begin(), libMesh::MeshBase::active_local_elements_end(), libMesh::MeshBase::active_pid_elements_begin(), libMesh::MeshBase::active_pid_elements_end(), libMesh::Parallel::Communicator::allgather(), libMesh::ParallelObject::comm(), libMesh::vectormap< Key, Tp >::count(), libMesh::MeshTools::create_bounding_box(), end, libMesh::MeshCommunication::find_global_indices(), libMesh::DofObject::id(), libMesh::vectormap< Key, Tp >::insert(), libMesh::libmesh_assert(), std::min(), libMesh::MeshBase::n_active_local_elem(), libMesh::Elem::n_nodes(), libMesh::ParallelObject::n_processors(), and libMesh::ParallelObject::processor_id().

Referenced by clone().

198 {
199  LOG_SCOPE("initialize()", "ParmetisPartitioner");
200 
201  const dof_id_type n_active_local_elem = mesh.n_active_local_elem();
202 
203  // Set parameters.
204  _pmetis->wgtflag = 2; // weights on vertices only
205  _pmetis->ncon = 1; // one weight per vertex
206  _pmetis->numflag = 0; // C-style 0-based numbering
207  _pmetis->nparts = static_cast<Parmetis::idx_t>(n_sbdmns); // number of subdomains to create
208  _pmetis->edgecut = 0; // the numbers of edges cut by the
209  // partition
210 
211  // Initialize data structures for ParMETIS
212  _pmetis->vtxdist.assign (mesh.n_processors()+1, 0);
213  _pmetis->tpwgts.assign (_pmetis->nparts, 1./_pmetis->nparts);
214  _pmetis->ubvec.assign (_pmetis->ncon, 1.05);
215  _pmetis->part.assign (n_active_local_elem, 0);
216  _pmetis->options.resize (5);
217  _pmetis->vwgt.resize (n_active_local_elem);
218 
219  // Set the options
220  _pmetis->options[0] = 1; // don't use default options
221  _pmetis->options[1] = 0; // default (level of timing)
222  _pmetis->options[2] = 15; // random seed (default)
223  _pmetis->options[3] = 2; // processor distribution and subdomain distribution are decoupled
224 
225  // Find the number of active elements on each processor. We cannot use
226  // mesh.n_active_elem_on_proc(pid) since that only returns the number of
227  // elements assigned to pid which are currently stored on the calling
228  // processor. This will not in general be correct for parallel meshes
229  // when (pid!=mesh.processor_id()).
230  _n_active_elem_on_proc.resize(mesh.n_processors());
231  mesh.comm().allgather(n_active_local_elem, _n_active_elem_on_proc);
232 
233  // count the total number of active elements in the mesh. Note we cannot
234  // use mesh.n_active_elem() in general since this only returns the number
235  // of active elements which are stored on the calling processor.
236  // We should not use n_active_elem for any allocation because that will
237  // be inherently unscalable, but it can be useful for libmesh_assertions.
238  dof_id_type n_active_elem=0;
239 
240  // Set up the vtxdist array. This will be the same on each processor.
241  // ***** Consult the Parmetis documentation. *****
242  libmesh_assert_equal_to (_pmetis->vtxdist.size(),
243  cast_int<std::size_t>(mesh.n_processors()+1));
244  libmesh_assert_equal_to (_pmetis->vtxdist[0], 0);
245 
246  for (processor_id_type pid=0; pid<mesh.n_processors(); pid++)
247  {
248  _pmetis->vtxdist[pid+1] = _pmetis->vtxdist[pid] + _n_active_elem_on_proc[pid];
249  n_active_elem += _n_active_elem_on_proc[pid];
250  }
251  libmesh_assert_equal_to (_pmetis->vtxdist.back(), static_cast<Parmetis::idx_t>(n_active_elem));
252 
253  // ParMetis expects the elements to be numbered in contiguous blocks
254  // by processor, i.e. [0, ne0), [ne0, ne0+ne1), ...
255  // Since we only partition active elements we should have no expectation
256  // that we currently have such a distribution. So we need to create it.
257  // Also, at the same time we are going to map all the active elements into a globally
258  // unique range [0,n_active_elem) which is *independent* of the current partitioning.
259  // This can be fed to ParMetis as the initial partitioning of the subdomains (decoupled
260  // from the partitioning of the objects themselves). This allows us to get the same
261  // resultant partitioning independent of the input partitioning.
262  libMesh::BoundingBox bbox =
264 
265  _global_index_by_pid_map.clear();
266 
267  // Maps active element ids into a contiguous range independent of partitioning.
268  // (only needs local scope)
269  vectormap<dof_id_type, dof_id_type> global_index_map;
270 
271  {
272  std::vector<dof_id_type> global_index;
273 
274  // create the mapping which is contiguous by processor
275  dof_id_type pid_offset=0;
276  for (processor_id_type pid=0; pid<mesh.n_processors(); pid++)
277  {
278  MeshBase::const_element_iterator it = mesh.active_pid_elements_begin(pid);
279  const MeshBase::const_element_iterator end = mesh.active_pid_elements_end(pid);
280 
281  // note that we may not have all (or any!) the active elements which belong on this processor,
282  // but by calling this on all processors a unique range in [0,_n_active_elem_on_proc[pid])
283  // is constructed. Only the indices for the elements we pass in are returned in the array.
284  MeshCommunication().find_global_indices (mesh.comm(),
285  bbox, it, end,
286  global_index);
287 
288  for (dof_id_type cnt=0; it != end; ++it)
289  {
290  const Elem * elem = *it;
291  // vectormap::count forces a sort, which is too expensive
292  // in a loop
293  // libmesh_assert (!_global_index_by_pid_map.count(elem->id()));
294  libmesh_assert_less (cnt, global_index.size());
295  libmesh_assert_less (global_index[cnt], _n_active_elem_on_proc[pid]);
296 
297  _global_index_by_pid_map.insert(std::make_pair(elem->id(), global_index[cnt++] + pid_offset));
298  }
299 
300  pid_offset += _n_active_elem_on_proc[pid];
301  }
302 
303  // create the unique mapping for all active elements independent of partitioning
304  {
305  MeshBase::const_element_iterator it = mesh.active_elements_begin();
306  const MeshBase::const_element_iterator end = mesh.active_elements_end();
307 
308  // Calling this on all processors a unique range in [0,n_active_elem) is constructed.
309  // Only the indices for the elements we pass in are returned in the array.
310  MeshCommunication().find_global_indices (mesh.comm(),
311  bbox, it, end,
312  global_index);
313 
314  for (dof_id_type cnt=0; it != end; ++it)
315  {
316  const Elem * elem = *it;
317  // vectormap::count forces a sort, which is too expensive
318  // in a loop
319  // libmesh_assert (!global_index_map.count(elem->id()));
320  libmesh_assert_less (cnt, global_index.size());
321  libmesh_assert_less (global_index[cnt], n_active_elem);
322 
323  global_index_map.insert(std::make_pair(elem->id(), global_index[cnt++]));
324  }
325  }
326  // really, shouldn't be close!
327  libmesh_assert_less_equal (global_index_map.size(), n_active_elem);
328  libmesh_assert_less_equal (_global_index_by_pid_map.size(), n_active_elem);
329 
330  // At this point the two maps should be the same size. If they are not
331  // then the number of active elements is not the same as the sum over all
332  // processors of the number of active elements per processor, which means
333  // there must be some unpartitioned objects out there.
334  if (global_index_map.size() != _global_index_by_pid_map.size())
335  libmesh_error_msg("ERROR: ParmetisPartitioner cannot handle unpartitioned objects!");
336  }
337 
338  // Finally, we need to initialize the vertex (partition) weights and the initial subdomain
339  // mapping. The subdomain mapping will be independent of the processor mapping, and is
340  // defined by a simple mapping of the global indices we just found.
341  {
342  std::vector<dof_id_type> subdomain_bounds(mesh.n_processors());
343 
344  const dof_id_type first_local_elem = _pmetis->vtxdist[mesh.processor_id()];
345 
346  for (processor_id_type pid=0; pid<mesh.n_processors(); pid++)
347  {
348  dof_id_type tgt_subdomain_size = 0;
349 
350  // watch out for the case that n_subdomains < n_processors
351  if (pid < static_cast<unsigned int>(_pmetis->nparts))
352  {
353  tgt_subdomain_size = n_active_elem/std::min
354  (cast_int<Parmetis::idx_t>(mesh.n_processors()), _pmetis->nparts);
355 
356  if (pid < n_active_elem%_pmetis->nparts)
357  tgt_subdomain_size++;
358  }
359  if (pid == 0)
360  subdomain_bounds[0] = tgt_subdomain_size;
361  else
362  subdomain_bounds[pid] = subdomain_bounds[pid-1] + tgt_subdomain_size;
363  }
364 
365  libmesh_assert_equal_to (subdomain_bounds.back(), n_active_elem);
366 
367  MeshBase::const_element_iterator elem_it = mesh.active_local_elements_begin();
368  const MeshBase::const_element_iterator elem_end = mesh.active_local_elements_end();
369 
370  for (; elem_it != elem_end; ++elem_it)
371  {
372  const Elem * elem = *elem_it;
373 
375  const dof_id_type global_index_by_pid =
376  _global_index_by_pid_map[elem->id()];
377  libmesh_assert_less (global_index_by_pid, n_active_elem);
378 
379  const dof_id_type local_index =
380  global_index_by_pid - first_local_elem;
381 
382  libmesh_assert_less (local_index, n_active_local_elem);
383  libmesh_assert_less (local_index, _pmetis->vwgt.size());
384 
385  // TODO:[BSK] maybe there is a better weight?
386  _pmetis->vwgt[local_index] = elem->n_nodes();
387 
388  // find the subdomain this element belongs in
389  libmesh_assert (global_index_map.count(elem->id()));
390  const dof_id_type global_index =
391  global_index_map[elem->id()];
392 
393  libmesh_assert_less (global_index, subdomain_bounds.back());
394 
395  const unsigned int subdomain_id =
396  std::distance(subdomain_bounds.begin(),
397  std::lower_bound(subdomain_bounds.begin(),
398  subdomain_bounds.end(),
399  global_index));
400  libmesh_assert_less (subdomain_id, static_cast<unsigned int>(_pmetis->nparts));
401  libmesh_assert_less (local_index, _pmetis->part.size());
402 
403  _pmetis->part[local_index] = subdomain_id;
404  }
405  }
406 }
difference_type count(const key_type &key) const
Definition: vectormap.h:210
libMesh::BoundingBox create_bounding_box(const MeshBase &mesh)
Definition: mesh_tools.C:341
MeshBase & mesh
uint8_t processor_id_type
Definition: id_types.h:99
std::vector< dof_id_type > _n_active_elem_on_proc
IterBase * end
vectormap< dof_id_type, dof_id_type > _global_index_by_pid_map
libmesh_assert(j)
UniquePtr< ParmetisHelper > _pmetis
void insert(const value_type &x)
Definition: vectormap.h:116
long double min(long double a, double b)
uint8_t dof_id_type
Definition: id_types.h:64
void libMesh::Partitioner::partition ( MeshBase mesh,
const unsigned int  n 
)
virtualinherited

Partitions the MeshBase into n parts by setting processor_id() on Nodes and Elems.

Note
If you are implementing a new type of Partitioner, you most likely do not want to override the partition() function, see instead the protected virtual _do_partition() method below. The partition() function is responsible for doing a lot of libmesh-internals-specific setup and finalization before and after the _do_partition() function is called. The only responsibility of the _do_partition() function, on the other hand, is to set the processor IDs of the elements according to a specific partitioning algorithm. See, e.g. MetisPartitioner for an example.

Definition at line 49 of file partitioner.C.

References libMesh::Partitioner::_do_partition(), libMesh::ParallelObject::comm(), libMesh::MeshTools::libmesh_assert_valid_remote_elems(), mesh, std::min(), libMesh::MeshBase::n_active_elem(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::MeshBase::redistribute(), libMesh::MeshBase::set_n_partitions(), libMesh::Partitioner::set_node_processor_ids(), libMesh::Partitioner::set_parent_processor_ids(), libMesh::Partitioner::single_partition(), and libMesh::MeshBase::update_post_partitioning().

Referenced by _do_repartition(), libMesh::Partitioner::partition(), and libMesh::Partitioner::~Partitioner().

51 {
52  libmesh_parallel_only(mesh.comm());
53 
54  // BSK - temporary fix while redistribution is integrated 6/26/2008
55  // Uncomment this to not repartition in parallel
56  // if (!mesh.is_serial())
57  // return;
58 
59  // we cannot partition into more pieces than we have
60  // active elements!
61  const unsigned int n_parts =
62  static_cast<unsigned int>
63  (std::min(mesh.n_active_elem(), static_cast<dof_id_type>(n)));
64 
65  // Set the number of partitions in the mesh
66  mesh.set_n_partitions()=n_parts;
67 
68  if (n_parts == 1)
69  {
70  this->single_partition (mesh);
71  return;
72  }
73 
74  // First assign a temporary partitioning to any unpartitioned elements
76 
77  // Call the partitioning function
78  this->_do_partition(mesh,n_parts);
79 
80  // Set the parent's processor ids
82 
83  // Redistribute elements if necessary, before setting node processor
84  // ids, to make sure those will be set consistently
85  mesh.redistribute();
86 
87 #ifdef DEBUG
89 
90  // Messed up elem processor_id()s can leave us without the child
91  // elements we need to restrict vectors on a distributed mesh
92  MeshTools::libmesh_assert_valid_procids<Elem>(mesh);
93 #endif
94 
95  // Set the node's processor ids
97 
98 #ifdef DEBUG
99  MeshTools::libmesh_assert_valid_procids<Elem>(mesh);
100 #endif
101 
102  // Give derived Mesh classes a chance to update any cached data to
103  // reflect the new partitioning
104  mesh.update_post_partitioning();
105 }
void single_partition(MeshBase &mesh)
Definition: partitioner.C:151
void libmesh_assert_valid_remote_elems(const MeshBase &mesh)
Definition: mesh_tools.C:1099
static void set_node_processor_ids(MeshBase &mesh)
Definition: partitioner.C:431
MeshBase & mesh
virtual void _do_partition(MeshBase &mesh, const unsigned int n)=0
static void partition_unpartitioned_elements(MeshBase &mesh)
Definition: partitioner.C:175
static void set_parent_processor_ids(MeshBase &mesh)
Definition: partitioner.C:256
long double min(long double a, double b)
uint8_t dof_id_type
Definition: id_types.h:64
void libMesh::Partitioner::partition ( MeshBase mesh)
virtualinherited

Partitions the MeshBase into mesh.n_processors() by setting processor_id() on Nodes and Elems.

Note
If you are implementing a new type of Partitioner, you most likely do not want to override the partition() function, see instead the protected virtual _do_partition() method below. The partition() function is responsible for doing a lot of libmesh-internals-specific setup and finalization before and after the _do_partition() function is called. The only responsibility of the _do_partition() function, on the other hand, is to set the processor IDs of the elements according to a specific partitioning algorithm. See, e.g. MetisPartitioner for an example.

Definition at line 42 of file partitioner.C.

References libMesh::ParallelObject::n_processors(), and libMesh::Partitioner::partition().

43 {
44  this->partition(mesh,mesh.n_processors());
45 }
MeshBase & mesh
virtual void partition(MeshBase &mesh, const unsigned int n)
Definition: partitioner.C:49
virtual void libMesh::Partitioner::partition_range ( MeshBase ,
MeshBase::element_iterator  ,
MeshBase::element_iterator  ,
const unsigned int   
)
inlinevirtualinherited

Partitions elements in the range (it, end) into n parts. The mesh from which the iterators are created must also be passed in, since it is a parallel object and has other useful information in it.

Although partition_range() is part of the public Partitioner interface, it should not generally be called by applications. Its main purpose is to support the SubdomainPartitioner, which uses it internally to individually partition ranges of elements before combining them into the final partitioning. Most of the time, the protected _do_partition() function is implemented in terms of partition_range() by passing a range which includes all the elements of the Mesh.

Reimplemented in libMesh::CentroidPartitioner, libMesh::MappedSubdomainPartitioner, libMesh::SFCPartitioner, libMesh::LinearPartitioner, and libMesh::MetisPartitioner.

Definition at line 120 of file partitioner.h.

References libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::Partitioner::repartition(), libMesh::Partitioner::set_node_processor_ids(), and libMesh::Partitioner::set_parent_processor_ids().

124  { libmesh_not_implemented(); }
void libMesh::Partitioner::partition_unpartitioned_elements ( MeshBase mesh)
staticinherited

These functions assign processor IDs to newly-created elements (in parallel) which are currently assigned to processor 0.

Definition at line 175 of file partitioner.C.

References libMesh::ParallelObject::n_processors().

Referenced by libMesh::Partitioner::partition(), libMesh::Partitioner::partition_range(), and libMesh::Partitioner::repartition().

176 {
178 }
MeshBase & mesh
static void partition_unpartitioned_elements(MeshBase &mesh)
Definition: partitioner.C:175
void libMesh::Partitioner::partition_unpartitioned_elements ( MeshBase mesh,
const unsigned int  n 
)
staticinherited

Definition at line 182 of file partitioner.C.

References libMesh::ParallelObject::comm(), libMesh::MeshTools::create_bounding_box(), end, libMesh::MeshCommunication::find_global_indices(), libMesh::MeshTools::n_elem(), libMesh::ParallelObject::n_processors(), libMesh::DofObject::processor_id(), libMesh::MeshBase::unpartitioned_elements_begin(), and libMesh::MeshBase::unpartitioned_elements_end().

184 {
185  MeshBase::element_iterator it = mesh.unpartitioned_elements_begin();
186  const MeshBase::element_iterator end = mesh.unpartitioned_elements_end();
187 
188  const dof_id_type n_unpartitioned_elements = MeshTools::n_elem (it, end);
189 
190  // the unpartitioned elements must exist on all processors. If the range is empty on one
191  // it is empty on all, and we can quit right here.
192  if (!n_unpartitioned_elements) return;
193 
194  // find the target subdomain sizes
195  std::vector<dof_id_type> subdomain_bounds(mesh.n_processors());
196 
197  for (processor_id_type pid=0; pid<mesh.n_processors(); pid++)
198  {
199  dof_id_type tgt_subdomain_size = 0;
200 
201  // watch out for the case that n_subdomains < n_processors
202  if (pid < n_subdomains)
203  {
204  tgt_subdomain_size = n_unpartitioned_elements/n_subdomains;
205 
206  if (pid < n_unpartitioned_elements%n_subdomains)
207  tgt_subdomain_size++;
208 
209  }
210 
211  //libMesh::out << "pid, #= " << pid << ", " << tgt_subdomain_size << std::endl;
212  if (pid == 0)
213  subdomain_bounds[0] = tgt_subdomain_size;
214  else
215  subdomain_bounds[pid] = subdomain_bounds[pid-1] + tgt_subdomain_size;
216  }
217 
218  libmesh_assert_equal_to (subdomain_bounds.back(), n_unpartitioned_elements);
219 
220  // create the unique mapping for all unpartitioned elements independent of partitioning
221  // determine the global indexing for all the unpartitioned elements
222  std::vector<dof_id_type> global_indices;
223 
224  // Calling this on all processors a unique range in [0,n_unpartitioned_elements) is constructed.
225  // Only the indices for the elements we pass in are returned in the array.
226  MeshCommunication().find_global_indices (mesh.comm(),
228  global_indices);
229 
230  for (dof_id_type cnt=0; it != end; ++it)
231  {
232  Elem * elem = *it;
233 
234  libmesh_assert_less (cnt, global_indices.size());
235  const dof_id_type global_index =
236  global_indices[cnt++];
237 
238  libmesh_assert_less (global_index, subdomain_bounds.back());
239  libmesh_assert_less (global_index, n_unpartitioned_elements);
240 
241  const processor_id_type subdomain_id =
242  cast_int<processor_id_type>
243  (std::distance(subdomain_bounds.begin(),
244  std::upper_bound(subdomain_bounds.begin(),
245  subdomain_bounds.end(),
246  global_index)));
247  libmesh_assert_less (subdomain_id, n_subdomains);
248 
249  elem->processor_id() = subdomain_id;
250  //libMesh::out << "assigning " << global_index << " to " << subdomain_id << std::endl;
251  }
252 }
dof_id_type n_elem(const MeshBase::const_element_iterator &begin, const MeshBase::const_element_iterator &end)
Definition: mesh_tools.C:684
libMesh::BoundingBox create_bounding_box(const MeshBase &mesh)
Definition: mesh_tools.C:341
MeshBase & mesh
uint8_t processor_id_type
Definition: id_types.h:99
IterBase * end
uint8_t dof_id_type
Definition: id_types.h:64
void libMesh::Partitioner::repartition ( MeshBase mesh,
const unsigned int  n 
)
inherited

Repartitions the MeshBase into n parts. (Some partitioning algorithms can repartition more efficiently than computing a new partitioning from scratch.) The default behavior is to simply call this->partition(mesh,n).

Definition at line 116 of file partitioner.C.

References libMesh::Partitioner::_do_repartition(), std::min(), libMesh::MeshBase::n_active_elem(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::MeshBase::set_n_partitions(), libMesh::Partitioner::set_node_processor_ids(), libMesh::Partitioner::set_parent_processor_ids(), and libMesh::Partitioner::single_partition().

Referenced by libMesh::Partitioner::partition_range(), and libMesh::Partitioner::repartition().

118 {
119  // we cannot partition into more pieces than we have
120  // active elements!
121  const unsigned int n_parts =
122  static_cast<unsigned int>
123  (std::min(mesh.n_active_elem(), static_cast<dof_id_type>(n)));
124 
125  // Set the number of partitions in the mesh
126  mesh.set_n_partitions()=n_parts;
127 
128  if (n_parts == 1)
129  {
130  this->single_partition (mesh);
131  return;
132  }
133 
134  // First assign a temporary partitioning to any unpartitioned elements
136 
137  // Call the partitioning function
138  this->_do_repartition(mesh,n_parts);
139 
140  // Set the parent's processor ids
142 
143  // Set the node's processor ids
145 }
void single_partition(MeshBase &mesh)
Definition: partitioner.C:151
static void set_node_processor_ids(MeshBase &mesh)
Definition: partitioner.C:431
MeshBase & mesh
virtual void _do_repartition(MeshBase &mesh, const unsigned int n)
Definition: partitioner.h:204
static void partition_unpartitioned_elements(MeshBase &mesh)
Definition: partitioner.C:175
static void set_parent_processor_ids(MeshBase &mesh)
Definition: partitioner.C:256
long double min(long double a, double b)
uint8_t dof_id_type
Definition: id_types.h:64
void libMesh::Partitioner::repartition ( MeshBase mesh)
inherited

Repartitions the MeshBase into mesh.n_processors() parts. This is required since some partitioning algorithms can repartition more efficiently than computing a new partitioning from scratch.

Definition at line 109 of file partitioner.C.

References libMesh::ParallelObject::n_processors(), and libMesh::Partitioner::repartition().

110 {
111  this->repartition(mesh,mesh.n_processors());
112 }
MeshBase & mesh
void repartition(MeshBase &mesh, const unsigned int n)
Definition: partitioner.C:116
void libMesh::Partitioner::set_node_processor_ids ( MeshBase mesh)
staticinherited

This function is called after partitioning to set the processor IDs for the nodes. By definition, a Node's processor ID is the minimum processor ID for all of the elements which share the node.

Definition at line 431 of file partitioner.C.

References libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::ParallelObject::comm(), libMesh::DofObject::id(), libMesh::DofObject::invalid_processor_id, libMesh::DofObject::invalidate_processor_id(), libMesh::libmesh_assert(), mesh, std::min(), libMesh::MeshTools::n_elem(), libMesh::Elem::n_nodes(), libMesh::ParallelObject::n_processors(), libMesh::Elem::node_ptr(), libMesh::MeshBase::node_ref(), libMesh::MeshBase::nodes_begin(), libMesh::MeshBase::nodes_end(), libMesh::MeshBase::not_active_elements_begin(), libMesh::MeshBase::not_active_elements_end(), libMesh::ParallelObject::processor_id(), libMesh::DofObject::processor_id(), libMesh::Parallel::Communicator::send_receive(), libMesh::MeshBase::subactive_elements_begin(), libMesh::MeshBase::subactive_elements_end(), libMesh::MeshBase::unpartitioned_elements_begin(), and libMesh::MeshBase::unpartitioned_elements_end().

Referenced by libMesh::UnstructuredMesh::all_first_order(), libMesh::Partitioner::partition(), libMesh::MeshBase::partition(), libMesh::Partitioner::partition_range(), libMesh::XdrIO::read(), libMesh::Partitioner::repartition(), and libMesh::BoundaryInfo::sync().

432 {
433  LOG_SCOPE("set_node_processor_ids()","Partitioner");
434 
435  // This function must be run on all processors at once
436  libmesh_parallel_only(mesh.comm());
437 
438  // If we have any unpartitioned elements at this
439  // stage there is a problem
440  libmesh_assert (MeshTools::n_elem(mesh.unpartitioned_elements_begin(),
441  mesh.unpartitioned_elements_end()) == 0);
442 
443 
444  // const dof_id_type orig_n_local_nodes = mesh.n_local_nodes();
445 
446  // libMesh::err << "[" << mesh.processor_id() << "]: orig_n_local_nodes="
447  // << orig_n_local_nodes << std::endl;
448 
449  // Build up request sets. Each node is currently owned by a processor because
450  // it is connected to an element owned by that processor. However, during the
451  // repartitioning phase that element may have been assigned a new processor id, but
452  // it is still resident on the original processor. We need to know where to look
453  // for new ids before assigning new ids, otherwise we may be asking the wrong processors
454  // for the wrong information.
455  //
456  // The only remaining issue is what to do with unpartitioned nodes. Since they are required
457  // to live on all processors we can simply rely on ourselves to number them properly.
458  std::vector<std::vector<dof_id_type> >
459  requested_node_ids(mesh.n_processors());
460 
461  // Loop over all the nodes, count the ones on each processor. We can skip ourself
462  std::vector<dof_id_type> ghost_nodes_from_proc(mesh.n_processors(), 0);
463 
464  MeshBase::node_iterator node_it = mesh.nodes_begin();
465  const MeshBase::node_iterator node_end = mesh.nodes_end();
466 
467  for (; node_it != node_end; ++node_it)
468  {
469  Node * node = *node_it;
470  libmesh_assert(node);
471  const processor_id_type current_pid = node->processor_id();
472  if (current_pid != mesh.processor_id() &&
473  current_pid != DofObject::invalid_processor_id)
474  {
475  libmesh_assert_less (current_pid, ghost_nodes_from_proc.size());
476  ghost_nodes_from_proc[current_pid]++;
477  }
478  }
479 
480  // We know how many objects live on each processor, so reserve()
481  // space for each.
482  for (processor_id_type pid=0; pid != mesh.n_processors(); ++pid)
483  requested_node_ids[pid].reserve(ghost_nodes_from_proc[pid]);
484 
485  // We need to get the new pid for each node from the processor
486  // which *currently* owns the node. We can safely skip ourself
487  for (node_it = mesh.nodes_begin(); node_it != node_end; ++node_it)
488  {
489  Node * node = *node_it;
490  libmesh_assert(node);
491  const processor_id_type current_pid = node->processor_id();
492  if (current_pid != mesh.processor_id() &&
493  current_pid != DofObject::invalid_processor_id)
494  {
495  libmesh_assert_less (current_pid, requested_node_ids.size());
496  libmesh_assert_less (requested_node_ids[current_pid].size(),
497  ghost_nodes_from_proc[current_pid]);
498  requested_node_ids[current_pid].push_back(node->id());
499  }
500 
501  // Unset any previously-set node processor ids
502  node->invalidate_processor_id();
503  }
504 
505  // Loop over all the active elements
506  MeshBase::element_iterator elem_it = mesh.active_elements_begin();
507  const MeshBase::element_iterator elem_end = mesh.active_elements_end();
508 
509  for ( ; elem_it != elem_end; ++elem_it)
510  {
511  Elem * elem = *elem_it;
512  libmesh_assert(elem);
513 
514  libmesh_assert_not_equal_to (elem->processor_id(), DofObject::invalid_processor_id);
515 
516  // For each node, set the processor ID to the min of
517  // its current value and this Element's processor id.
518  //
519  // TODO: we would probably get better parallel partitioning if
520  // we did something like "min for even numbered nodes, max for
521  // odd numbered". We'd need to be careful about how that would
522  // affect solution ordering for I/O, though.
523  for (unsigned int n=0; n<elem->n_nodes(); ++n)
524  elem->node_ptr(n)->processor_id() = std::min(elem->node_ptr(n)->processor_id(),
525  elem->processor_id());
526  }
527 
528  // And loop over the subactive elements, but don't reassign
529  // nodes that are already active on another processor.
530  MeshBase::element_iterator sub_it = mesh.subactive_elements_begin();
531  const MeshBase::element_iterator sub_end = mesh.subactive_elements_end();
532 
533  for ( ; sub_it != sub_end; ++sub_it)
534  {
535  Elem * elem = *sub_it;
536  libmesh_assert(elem);
537 
538  libmesh_assert_not_equal_to (elem->processor_id(), DofObject::invalid_processor_id);
539 
540  for (unsigned int n=0; n<elem->n_nodes(); ++n)
541  if (elem->node_ptr(n)->processor_id() == DofObject::invalid_processor_id)
542  elem->node_ptr(n)->processor_id() = elem->processor_id();
543  }
544 
545  // Same for the inactive elements -- we will have already gotten most of these
546  // nodes, *except* for the case of a parent with a subset of children which are
547  // ghost elements. In that case some of the parent nodes will not have been
548  // properly handled yet
549  MeshBase::element_iterator not_it = mesh.not_active_elements_begin();
550  const MeshBase::element_iterator not_end = mesh.not_active_elements_end();
551 
552  for ( ; not_it != not_end; ++not_it)
553  {
554  Elem * elem = *not_it;
555  libmesh_assert(elem);
556 
557  libmesh_assert_not_equal_to (elem->processor_id(), DofObject::invalid_processor_id);
558 
559  for (unsigned int n=0; n<elem->n_nodes(); ++n)
560  if (elem->node_ptr(n)->processor_id() == DofObject::invalid_processor_id)
561  elem->node_ptr(n)->processor_id() = elem->processor_id();
562  }
563 
564  // We can't assert that all nodes are connected to elements, because
565  // a DistributedMesh with NodeConstraints might have pulled in some
566  // remote nodes solely for evaluating those constraints.
567  // MeshTools::libmesh_assert_connected_nodes(mesh);
568 
569  // For such nodes, we'll do a sanity check later when making sure
570  // that we successfully reset their processor ids to something
571  // valid.
572 
573  // Next set node ids from other processors, excluding self
574  for (processor_id_type p=1; p != mesh.n_processors(); ++p)
575  {
576  // Trade my requests with processor procup and procdown
577  processor_id_type procup = cast_int<processor_id_type>
578  ((mesh.processor_id() + p) % mesh.n_processors());
579  processor_id_type procdown = cast_int<processor_id_type>
580  ((mesh.n_processors() + mesh.processor_id() - p) %
581  mesh.n_processors());
582  std::vector<dof_id_type> request_to_fill;
583  mesh.comm().send_receive(procup, requested_node_ids[procup],
584  procdown, request_to_fill);
585 
586  // Fill those requests in-place
587  for (std::size_t i=0; i != request_to_fill.size(); ++i)
588  {
589  Node & node = mesh.node_ref(request_to_fill[i]);
590  const processor_id_type new_pid = node.processor_id();
591 
592  // We may have an invalid processor_id() on nodes that have been
593  // "detached" from coarsened-away elements but that have not yet
594  // themselves been removed.
595  // libmesh_assert_not_equal_to (new_pid, DofObject::invalid_processor_id);
596  // libmesh_assert_less (new_pid, mesh.n_partitions()); // this is the correct test --
597  request_to_fill[i] = new_pid; // the number of partitions may
598  } // not equal the number of processors
599 
600  // Trade back the results
601  std::vector<dof_id_type> filled_request;
602  mesh.comm().send_receive(procdown, request_to_fill,
603  procup, filled_request);
604  libmesh_assert_equal_to (filled_request.size(), requested_node_ids[procup].size());
605 
606  // And copy the id changes we've now been informed of
607  for (std::size_t i=0; i != filled_request.size(); ++i)
608  {
609  Node & node = mesh.node_ref(requested_node_ids[procup][i]);
610 
611  // this is the correct test -- the number of partitions may
612  // not equal the number of processors
613 
614  // But: we may have an invalid processor_id() on nodes that
615  // have been "detached" from coarsened-away elements but
616  // that have not yet themselves been removed.
617  // libmesh_assert_less (filled_request[i], mesh.n_partitions());
618 
619  node.processor_id(cast_int<processor_id_type>(filled_request[i]));
620  }
621  }
622 
623 #ifdef DEBUG
624  MeshTools::libmesh_assert_valid_procids<Node>(mesh);
625 #endif
626 }
dof_id_type n_elem(const MeshBase::const_element_iterator &begin, const MeshBase::const_element_iterator &end)
Definition: mesh_tools.C:684
MeshBase & mesh
uint8_t processor_id_type
Definition: id_types.h:99
libmesh_assert(j)
static const processor_id_type invalid_processor_id
Definition: dof_object.h:345
long double min(long double a, double b)
void libMesh::Partitioner::set_parent_processor_ids ( MeshBase mesh)
staticinherited

This function is called after partitioning to set the processor IDs for the inactive parent elements. A parent's processor ID is the same as its first child.

Definition at line 256 of file partitioner.C.

References libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::Elem::active_family_tree(), libMesh::MeshBase::ancestor_elements_begin(), libMesh::MeshBase::ancestor_elements_end(), libMesh::Elem::child_ptr(), libMesh::ParallelObject::comm(), libMesh::Partitioner::communication_blocksize, end, libMesh::DofObject::id(), libMesh::DofObject::invalid_processor_id, libMesh::DofObject::invalidate_processor_id(), libMesh::Elem::is_remote(), libMesh::MeshBase::is_serial(), libMesh::libmesh_assert(), libMesh::libmesh_ignore(), libMesh::MeshBase::max_elem_id(), std::min(), libMesh::Parallel::Communicator::min(), libMesh::Elem::n_children(), libMesh::MeshTools::n_elem(), libMesh::Elem::parent(), libMesh::processor_id(), libMesh::DofObject::processor_id(), libMesh::Elem::total_family_tree(), libMesh::MeshBase::unpartitioned_elements_begin(), and libMesh::MeshBase::unpartitioned_elements_end().

Referenced by libMesh::Partitioner::partition(), libMesh::Partitioner::partition_range(), and libMesh::Partitioner::repartition().

257 {
258  // Ignore the parameter when !LIBMESH_ENABLE_AMR
260 
261  LOG_SCOPE("set_parent_processor_ids()", "Partitioner");
262 
263 #ifdef LIBMESH_ENABLE_AMR
264 
265  // If the mesh is serial we have access to all the elements,
266  // in particular all the active ones. We can therefore set
267  // the parent processor ids indirectly through their children, and
268  // set the subactive processor ids while examining their active
269  // ancestors.
270  // By convention a parent is assigned to the minimum processor
271  // of all its children, and a subactive is assigned to the processor
272  // of its active ancestor.
273  if (mesh.is_serial())
274  {
275  // Loop over all the active elements in the mesh
276  MeshBase::element_iterator it = mesh.active_elements_begin();
277  const MeshBase::element_iterator end = mesh.active_elements_end();
278 
279  for ( ; it!=end; ++it)
280  {
281  Elem * child = *it;
282 
283  // First set descendents
284 
285  std::vector<const Elem *> subactive_family;
286  child->total_family_tree(subactive_family);
287  for (std::size_t i = 0; i != subactive_family.size(); ++i)
288  const_cast<Elem *>(subactive_family[i])->processor_id() = child->processor_id();
289 
290  // Then set ancestors
291 
292  Elem * parent = child->parent();
293 
294  while (parent)
295  {
296  // invalidate the parent id, otherwise the min below
297  // will not work if the current parent id is less
298  // than all the children!
299  parent->invalidate_processor_id();
300 
301  for (unsigned int c=0; c<parent->n_children(); c++)
302  {
303  child = parent->child_ptr(c);
304  libmesh_assert(child);
305  libmesh_assert(!child->is_remote());
306  libmesh_assert_not_equal_to (child->processor_id(), DofObject::invalid_processor_id);
307  parent->processor_id() = std::min(parent->processor_id(),
308  child->processor_id());
309  }
310  parent = parent->parent();
311  }
312  }
313  }
314 
315  // When the mesh is parallel we cannot guarantee that parents have access to
316  // all their children.
317  else
318  {
319  // Setting subactive processor ids is easy: we can guarantee
320  // that children have access to all their parents.
321 
322  // Loop over all the active elements in the mesh
323  MeshBase::element_iterator it = mesh.active_elements_begin();
324  const MeshBase::element_iterator end = mesh.active_elements_end();
325 
326  for ( ; it!=end; ++it)
327  {
328  Elem * child = *it;
329 
330  std::vector<const Elem *> subactive_family;
331  child->total_family_tree(subactive_family);
332  for (std::size_t i = 0; i != subactive_family.size(); ++i)
333  const_cast<Elem *>(subactive_family[i])->processor_id() = child->processor_id();
334  }
335 
336  // When the mesh is parallel we cannot guarantee that parents have access to
337  // all their children.
338 
339  // We will use a brute-force approach here. Each processor finds its parent
340  // elements and sets the parent pid to the minimum of its
341  // semilocal descendants.
342  // A global reduction is then performed to make sure the true minimum is found.
343  // As noted, this is required because we cannot guarantee that a parent has
344  // access to all its children on any single processor.
345  libmesh_parallel_only(mesh.comm());
346  libmesh_assert(MeshTools::n_elem(mesh.unpartitioned_elements_begin(),
347  mesh.unpartitioned_elements_end()) == 0);
348 
349  const dof_id_type max_elem_id = mesh.max_elem_id();
350 
351  std::vector<processor_id_type>
352  parent_processor_ids (std::min(communication_blocksize,
353  max_elem_id));
354 
355  for (dof_id_type blk=0, last_elem_id=0; last_elem_id<max_elem_id; blk++)
356  {
357  last_elem_id =
358  std::min(static_cast<dof_id_type>((blk+1)*communication_blocksize),
359  max_elem_id);
360  const dof_id_type first_elem_id = blk*communication_blocksize;
361 
362  std::fill (parent_processor_ids.begin(),
363  parent_processor_ids.end(),
365 
366  // first build up local contributions to parent_processor_ids
367  MeshBase::element_iterator not_it = mesh.ancestor_elements_begin();
368  const MeshBase::element_iterator not_end = mesh.ancestor_elements_end();
369 
370  bool have_parent_in_block = false;
371 
372  for ( ; not_it != not_end; ++not_it)
373  {
374  Elem * parent = *not_it;
375 
376  const dof_id_type parent_idx = parent->id();
377  libmesh_assert_less (parent_idx, max_elem_id);
378 
379  if ((parent_idx >= first_elem_id) &&
380  (parent_idx < last_elem_id))
381  {
382  have_parent_in_block = true;
384 
385  std::vector<const Elem *> active_family;
386  parent->active_family_tree(active_family);
387  for (std::size_t i = 0; i != active_family.size(); ++i)
388  parent_pid = std::min (parent_pid, active_family[i]->processor_id());
389 
390  const dof_id_type packed_idx = parent_idx - first_elem_id;
391  libmesh_assert_less (packed_idx, parent_processor_ids.size());
392 
393  parent_processor_ids[packed_idx] = parent_pid;
394  }
395  }
396 
397  // then find the global minimum
398  mesh.comm().min (parent_processor_ids);
399 
400  // and assign the ids, if we have a parent in this block.
401  if (have_parent_in_block)
402  for (not_it = mesh.ancestor_elements_begin();
403  not_it != not_end; ++not_it)
404  {
405  Elem * parent = *not_it;
406 
407  const dof_id_type parent_idx = parent->id();
408 
409  if ((parent_idx >= first_elem_id) &&
410  (parent_idx < last_elem_id))
411  {
412  const dof_id_type packed_idx = parent_idx - first_elem_id;
413  libmesh_assert_less (packed_idx, parent_processor_ids.size());
414 
415  const processor_id_type parent_pid =
416  parent_processor_ids[packed_idx];
417 
418  libmesh_assert_not_equal_to (parent_pid, DofObject::invalid_processor_id);
419 
420  parent->processor_id() = parent_pid;
421  }
422  }
423  }
424  }
425 
426 #endif // LIBMESH_ENABLE_AMR
427 }
dof_id_type n_elem(const MeshBase::const_element_iterator &begin, const MeshBase::const_element_iterator &end)
Definition: mesh_tools.C:684
MeshBase & mesh
uint8_t processor_id_type
Definition: id_types.h:99
IterBase * end
libmesh_assert(j)
static const processor_id_type invalid_processor_id
Definition: dof_object.h:345
void libmesh_ignore(const T &)
static const dof_id_type communication_blocksize
Definition: partitioner.h:211
processor_id_type processor_id()
Definition: libmesh_base.h:96
long double min(long double a, double b)
uint8_t dof_id_type
Definition: id_types.h:64
void libMesh::Partitioner::single_partition ( MeshBase mesh)
protectedinherited

Trivially "partitions" the mesh for one processor. Simply loops through the elements and assigns all of them to processor 0. Is is provided as a separate function so that derived classes may use it without reimplementing it.

Definition at line 151 of file partitioner.C.

References libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), and libMesh::Partitioner::single_partition_range().

Referenced by libMesh::SubdomainPartitioner::_do_partition(), libMesh::Partitioner::attach_weights(), libMesh::Partitioner::partition(), and libMesh::Partitioner::repartition().

152 {
153  this->single_partition_range(mesh.elements_begin(),
154  mesh.elements_end());
155 }
MeshBase & mesh
void single_partition_range(MeshBase::element_iterator it, MeshBase::element_iterator end)
Definition: partitioner.C:159
void libMesh::Partitioner::single_partition_range ( MeshBase::element_iterator  it,
MeshBase::element_iterator  end 
)
protectedinherited

Slightly generalized version of single_partition which acts on a range of elements defined by the pair of iterators (it, end).

Definition at line 159 of file partitioner.C.

References end, libMesh::Elem::n_nodes(), libMesh::Elem::node_ptr(), and libMesh::DofObject::processor_id().

Referenced by libMesh::Partitioner::attach_weights(), libMesh::LinearPartitioner::partition_range(), libMesh::MappedSubdomainPartitioner::partition_range(), libMesh::CentroidPartitioner::partition_range(), and libMesh::Partitioner::single_partition().

161 {
162  LOG_SCOPE("single_partition_range()", "Partitioner");
163 
164  for ( ; it != end; ++it)
165  {
166  Elem * elem = *it;
167  elem->processor_id() = 0;
168 
169  // Assign all this element's nodes to processor 0 as well.
170  for (unsigned int n=0; n<elem->n_nodes(); ++n)
171  elem->node_ptr(n)->processor_id() = 0;
172  }
173 }
IterBase * end

Member Data Documentation

vectormap<dof_id_type, dof_id_type> libMesh::ParmetisPartitioner::_global_index_by_pid_map
private

Maps active element ids into a contiguous range, as needed by ParMETIS.

Definition at line 119 of file parmetis_partitioner.h.

std::vector<dof_id_type> libMesh::ParmetisPartitioner::_n_active_elem_on_proc
private

The number of active elements on each processor.

Note
ParMETIS requires that each processor have some active elements; it will abort if any processor passes a NULL _part array.

Definition at line 114 of file parmetis_partitioner.h.

UniquePtr<ParmetisHelper> libMesh::ParmetisPartitioner::_pmetis
private

Pointer to the Parmetis-specific data structures. Lets us avoid including parmetis.h here.

Definition at line 125 of file parmetis_partitioner.h.

ErrorVector* libMesh::Partitioner::_weights
protectedinherited

The weights that might be used for partitioning.

Definition at line 216 of file partitioner.h.

Referenced by libMesh::MetisPartitioner::attach_weights().

const dof_id_type libMesh::Partitioner::communication_blocksize = 1000000
staticprotectedinherited

The blocksize to use when doing blocked parallel communication. This limits the maximum vector size which can be used in a single communication step.

Definition at line 211 of file partitioner.h.

Referenced by libMesh::Partitioner::set_parent_processor_ids().


The documentation for this class was generated from the following files: