libMesh
Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | Private Member Functions | Private Attributes | List of all members
libMesh::ParmetisPartitioner Class Reference

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

#include <parmetis_partitioner.h>

Inheritance diagram for libMesh::ParmetisPartitioner:
[legend]

Public Member Functions

 ParmetisPartitioner ()
 Constructor. More...
 
 ~ParmetisPartitioner ()
 Destructor. More...
 
virtual UniquePtr< Partitionerclone () const libmesh_override
 
virtual void partition (MeshBase &mesh, const unsigned int n)
 Partitions the MeshBase into n parts by setting processor_id() on Nodes and Elems. More...
 
virtual void partition (MeshBase &mesh)
 Partitions the MeshBase into mesh.n_processors() by setting processor_id() on Nodes and Elems. More...
 
virtual void partition_range (MeshBase &, MeshBase::element_iterator, MeshBase::element_iterator, const unsigned int)
 Partitions elements in the range (it, end) into n parts. More...
 
void repartition (MeshBase &mesh, const unsigned int n)
 Repartitions the MeshBase into n parts. More...
 
void repartition (MeshBase &mesh)
 Repartitions the MeshBase into mesh.n_processors() parts. More...
 
virtual void attach_weights (ErrorVector *)
 Attach weights that can be used for partitioning. More...
 

Static Public Member Functions

static void partition_unpartitioned_elements (MeshBase &mesh)
 These functions assign processor IDs to newly-created elements (in parallel) which are currently assigned to processor 0. More...
 
static void partition_unpartitioned_elements (MeshBase &mesh, const unsigned int n)
 
static void set_parent_processor_ids (MeshBase &mesh)
 This function is called after partitioning to set the processor IDs for the inactive parent elements. More...
 
static void set_node_processor_ids (MeshBase &mesh)
 This function is called after partitioning to set the processor IDs for the nodes. More...
 

Protected Member Functions

virtual void _do_repartition (MeshBase &mesh, const unsigned int n) libmesh_override
 Parmetis can handle dynamically repartitioning a mesh such that the redistribution costs are minimized. More...
 
virtual void _do_partition (MeshBase &mesh, const unsigned int n) libmesh_override
 Partition the MeshBase into n subdomains. More...
 
void single_partition (MeshBase &mesh)
 Trivially "partitions" the mesh for one processor. More...
 
void single_partition_range (MeshBase::element_iterator it, MeshBase::element_iterator end)
 Slightly generalized version of single_partition which acts on a range of elements defined by the pair of iterators (it, end). More...
 

Protected Attributes

ErrorVector_weights
 The weights that might be used for partitioning. More...
 

Static Protected Attributes

static const dof_id_type communication_blocksize = 1000000
 The blocksize to use when doing blocked parallel communication. More...
 

Private Member Functions

void initialize (const MeshBase &mesh, const unsigned int n_sbdmns)
 Initialize data structures. More...
 
void build_graph (const MeshBase &mesh)
 Build the graph. More...
 
void assign_partitioning (MeshBase &mesh)
 Assign the computed partitioning to the mesh. More...
 

Private Attributes

std::vector< dof_id_type_n_active_elem_on_proc
 The number of active elements on each processor. More...
 
vectormap< dof_id_type, dof_id_type_global_index_by_pid_map
 Maps active element ids into a contiguous range, as needed by ParMETIS. More...
 
UniquePtr< ParmetisHelper_pmetis
 Pointer to the Parmetis-specific data structures. More...
 

Detailed Description

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

Author
Benjamin S. Kirk
Date
2003 Partitioner which provides an interface to ParMETIS.

Definition at line 46 of file parmetis_partitioner.h.

Constructor & Destructor Documentation

libMesh::ParmetisPartitioner::ParmetisPartitioner ( )

Constructor.

Definition at line 62 of file parmetis_partitioner.C.

Referenced by clone().

64  : _pmetis(new ParmetisHelper)
65 #endif
66 {}
67 
68 
69 
71 {
72 }
UniquePtr< ParmetisHelper > _pmetis
Pointer to the Parmetis-specific data structures.
~ParmetisPartitioner()
Destructor.
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 76 of file parmetis_partitioner.C.

Referenced by clone().

78 {
79  this->_do_repartition (mesh, n_sbdmns);
80 }
virtual void _do_repartition(MeshBase &mesh, const unsigned int n) libmesh_override
Parmetis can handle dynamically repartitioning a mesh such that the redistribution costs are minimize...
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 84 of file parmetis_partitioner.C.

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

Referenced by clone().

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

References libMesh::MeshBase::active_element_ptr_range(), libMesh::ParallelObject::comm(), libMesh::libmesh_assert(), libMesh::MeshBase::n_active_local_elem(), libMesh::ParallelObject::n_processors(), libMesh::ParallelObject::processor_id(), and libMesh::Parallel::Communicator::send_receive().

Referenced by clone().

593 {
594  LOG_SCOPE("assign_partitioning()", "ParmetisPartitioner");
595 
596  // This function must be run on all processors at once
597  libmesh_parallel_only(mesh.comm());
598 
599  const dof_id_type
600  first_local_elem = _pmetis->vtxdist[mesh.processor_id()];
601 
602 #ifndef NDEBUG
603  const dof_id_type n_active_local_elem = mesh.n_active_local_elem();
604 #endif
605 
606  std::vector<std::vector<dof_id_type>>
607  requested_ids(mesh.n_processors()),
608  requests_to_fill(mesh.n_processors());
609 
610  for (auto & elem : mesh.active_element_ptr_range())
611  {
612  // we need to get the index from the owning processor
613  // (note we cannot assign it now -- we are iterating
614  // over elements again and this will be bad!)
615  libmesh_assert_less (elem->processor_id(), requested_ids.size());
616  requested_ids[elem->processor_id()].push_back(elem->id());
617  }
618 
619  // Trade with all processors (including self) to get their indices
620  for (processor_id_type pid=0; pid<mesh.n_processors(); pid++)
621  {
622  // Trade my requests with processor procup and procdown
623  const processor_id_type procup = (mesh.processor_id() + pid) % mesh.n_processors();
624  const processor_id_type procdown = (mesh.n_processors() +
625  mesh.processor_id() - pid) % mesh.n_processors();
626 
627  mesh.comm().send_receive (procup, requested_ids[procup],
628  procdown, requests_to_fill[procdown]);
629 
630  // we can overwrite these requested ids in-place.
631  for (std::size_t i=0; i<requests_to_fill[procdown].size(); i++)
632  {
633  const dof_id_type requested_elem_index =
634  requests_to_fill[procdown][i];
635 
636  libmesh_assert(_global_index_by_pid_map.count(requested_elem_index));
637 
638  const dof_id_type global_index_by_pid =
639  _global_index_by_pid_map[requested_elem_index];
640 
641  const dof_id_type local_index =
642  global_index_by_pid - first_local_elem;
643 
644  libmesh_assert_less (local_index, _pmetis->part.size());
645  libmesh_assert_less (local_index, n_active_local_elem);
646 
647  const unsigned int elem_procid =
648  static_cast<unsigned int>(_pmetis->part[local_index]);
649 
650  libmesh_assert_less (elem_procid, static_cast<unsigned int>(_pmetis->nparts));
651 
652  requests_to_fill[procdown][i] = elem_procid;
653  }
654 
655  // Trade back
656  mesh.comm().send_receive (procdown, requests_to_fill[procdown],
657  procup, requested_ids[procup]);
658  }
659 
660  // and finally assign the partitioning.
661  // note we are iterating in exactly the same order
662  // used to build up the request, so we can expect the
663  // required entries to be in the proper sequence.
664  std::vector<unsigned int> counters(mesh.n_processors(), 0);
665  for (auto & elem : mesh.active_element_ptr_range())
666  {
667  const processor_id_type current_pid = elem->processor_id();
668 
669  libmesh_assert_less (counters[current_pid], requested_ids[current_pid].size());
670 
671  const processor_id_type elem_procid =
672  requested_ids[current_pid][counters[current_pid]++];
673 
674  libmesh_assert_less (elem_procid, static_cast<unsigned int>(_pmetis->nparts));
675  elem->processor_id() = elem_procid;
676  }
677 }
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
Maps active element ids into a contiguous range, as needed by ParMETIS.
libmesh_assert(j)
UniquePtr< ParmetisHelper > _pmetis
Pointer to the Parmetis-specific data structures.
uint8_t dof_id_type
Definition: id_types.h:64
virtual void libMesh::Partitioner::attach_weights ( ErrorVector )
virtualinherited

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 403 of file parmetis_partitioner.C.

References libMesh::Elem::active(), libMesh::MeshBase::active_element_ptr_range(), libMesh::MeshBase::active_local_element_ptr_range(), libMesh::as_range(), libMesh::DofObject::id(), libMesh::libmesh_assert(), libmesh_nullptr, libMesh::MeshBase::n_active_local_elem(), libMesh::Elem::neighbor_ptr(), and libMesh::ParallelObject::processor_id().

Referenced by clone().

404 {
405  LOG_SCOPE("build_graph()", "ParmetisPartitioner");
406 
407  // build the graph in distributed CSR format. Note that
408  // the edges in the graph will correspond to
409  // face neighbors
410  const dof_id_type n_active_local_elem = mesh.n_active_local_elem();
411 
412  // If we have boundary elements in this mesh, we want to account for
413  // the connectivity between them and interior elements. We can find
414  // interior elements from boundary elements, but we need to build up
415  // a lookup map to do the reverse.
416 
417  typedef std::unordered_multimap<const Elem *, const Elem *> map_type;
418  map_type interior_to_boundary_map;
419 
420  for (const auto & elem : mesh.active_element_ptr_range())
421  {
422  // If we don't have an interior_parent then there's nothing to look us
423  // up.
424  if ((elem->dim() >= LIBMESH_DIM) ||
425  !elem->interior_parent())
426  continue;
427 
428  // get all relevant interior elements
429  std::set<const Elem *> neighbor_set;
430  elem->find_interior_neighbors(neighbor_set);
431 
432  std::set<const Elem *>::iterator n_it = neighbor_set.begin();
433  for (; n_it != neighbor_set.end(); ++n_it)
434  interior_to_boundary_map.insert(std::make_pair(*n_it, elem));
435  }
436 
437 #ifdef LIBMESH_ENABLE_AMR
438  std::vector<const Elem *> neighbors_offspring;
439 #endif
440 
441  std::vector<std::vector<dof_id_type>> graph(n_active_local_elem);
442  dof_id_type graph_size=0;
443 
444  const dof_id_type first_local_elem = _pmetis->vtxdist[mesh.processor_id()];
445 
446  for (const auto & elem : mesh.active_local_element_ptr_range())
447  {
449  const dof_id_type global_index_by_pid =
450  _global_index_by_pid_map[elem->id()];
451 
452  const dof_id_type local_index =
453  global_index_by_pid - first_local_elem;
454  libmesh_assert_less (local_index, n_active_local_elem);
455 
456  std::vector<dof_id_type> & graph_row = graph[local_index];
457 
458  // Loop over the element's neighbors. An element
459  // adjacency corresponds to a face neighbor
460  for (auto neighbor : elem->neighbor_ptr_range())
461  {
462  if (neighbor != libmesh_nullptr)
463  {
464  // If the neighbor is active treat it
465  // as a connection
466  if (neighbor->active())
467  {
469  const dof_id_type neighbor_global_index_by_pid =
470  _global_index_by_pid_map[neighbor->id()];
471 
472  graph_row.push_back(neighbor_global_index_by_pid);
473  graph_size++;
474  }
475 
476 #ifdef LIBMESH_ENABLE_AMR
477 
478  // Otherwise we need to find all of the
479  // neighbor's children that are connected to
480  // us and add them
481  else
482  {
483  // The side of the neighbor to which
484  // we are connected
485  const unsigned int ns =
486  neighbor->which_neighbor_am_i (elem);
487  libmesh_assert_less (ns, neighbor->n_neighbors());
488 
489  // Get all the active children (& grandchildren, etc...)
490  // of the neighbor
491 
492  // FIXME - this is the wrong thing, since we
493  // should be getting the active family tree on
494  // our side only. But adding too many graph
495  // links may cause hanging nodes to tend to be
496  // on partition interiors, which would reduce
497  // communication overhead for constraint
498  // equations, so we'll leave it.
499 
500  neighbor->active_family_tree (neighbors_offspring);
501 
502  // Get all the neighbor's children that
503  // live on that side and are thus connected
504  // to us
505  for (std::size_t nc=0; nc<neighbors_offspring.size(); nc++)
506  {
507  const Elem * child =
508  neighbors_offspring[nc];
509 
510  // This does not assume a level-1 mesh.
511  // Note that since children have sides numbered
512  // coincident with the parent then this is a sufficient test.
513  if (child->neighbor_ptr(ns) == elem)
514  {
515  libmesh_assert (child->active());
517  const dof_id_type child_global_index_by_pid =
518  _global_index_by_pid_map[child->id()];
519 
520  graph_row.push_back(child_global_index_by_pid);
521  graph_size++;
522  }
523  }
524  }
525 
526 #endif /* ifdef LIBMESH_ENABLE_AMR */
527 
528 
529  }
530  }
531 
532  if ((elem->dim() < LIBMESH_DIM) &&
533  elem->interior_parent())
534  {
535  // get all relevant interior elements
536  std::set<const Elem *> neighbor_set;
537  elem->find_interior_neighbors(neighbor_set);
538 
539  std::set<const Elem *>::iterator n_it = neighbor_set.begin();
540  for (; n_it != neighbor_set.end(); ++n_it)
541  {
542  // FIXME - non-const versions of the Elem set methods
543  // would be nice
544  Elem * neighbor = const_cast<Elem *>(*n_it);
545 
546  const dof_id_type neighbor_global_index_by_pid =
547  _global_index_by_pid_map[neighbor->id()];
548 
549  graph_row.push_back(neighbor_global_index_by_pid);
550  graph_size++;
551  }
552  }
553 
554  // Check for any boundary neighbors
555  for (const auto & pr : as_range(interior_to_boundary_map.equal_range(elem)))
556  {
557  const Elem * neighbor = pr.second;
558 
559  const dof_id_type neighbor_global_index_by_pid =
560  _global_index_by_pid_map[neighbor->id()];
561 
562  graph_row.push_back(neighbor_global_index_by_pid);
563  graph_size++;
564  }
565  }
566 
567  // Reserve space in the adjacency array
568  _pmetis->xadj.clear();
569  _pmetis->xadj.reserve (n_active_local_elem + 1);
570  _pmetis->adjncy.clear();
571  _pmetis->adjncy.reserve (graph_size);
572 
573  for (std::size_t r=0; r<graph.size(); r++)
574  {
575  _pmetis->xadj.push_back(_pmetis->adjncy.size());
576  std::vector<dof_id_type> graph_row; // build this empty
577  graph_row.swap(graph[r]); // this will deallocate at the end of scope
578  _pmetis->adjncy.insert(_pmetis->adjncy.end(),
579  graph_row.begin(),
580  graph_row.end());
581  }
582 
583  // The end of the adjacency array for the last elem
584  _pmetis->xadj.push_back(_pmetis->adjncy.size());
585 
586  libmesh_assert_equal_to (_pmetis->xadj.size(), n_active_local_elem+1);
587  libmesh_assert_equal_to (_pmetis->adjncy.size(), graph_size);
588 }
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
Maps active element ids into a contiguous range, as needed by ParMETIS.
libmesh_assert(j)
UniquePtr< ParmetisHelper > _pmetis
Pointer to the Parmetis-specific data structures.
SimpleRange< I > as_range(const std::pair< I, I > &p)
Helper function that allows us to treat a homogenous pair as a range.
Definition: simple_range.h:57
uint8_t dof_id_type
Definition: id_types.h:64
virtual UniquePtr<Partitioner> libMesh::ParmetisPartitioner::clone ( ) const
virtual
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 194 of file parmetis_partitioner.C.

References libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::MeshBase::active_local_element_ptr_range(), 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(), distance(), 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::ParallelObject::n_processors(), and libMesh::ParallelObject::processor_id().

Referenced by clone().

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

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)
These functions assign processor IDs to newly-created elements (in parallel) which are currently assi...
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(), distance(), 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)
Count up the number of elements of a specific type (as defined by an iterator range).
Definition: mesh_tools.C:656
libMesh::BoundingBox create_bounding_box(const MeshBase &mesh)
The same functionality as the deprecated MeshTools::bounding_box().
Definition: mesh_tools.C:329
MeshBase & mesh
uint8_t processor_id_type
Definition: id_types.h:99
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end...
Real distance(const Point &p)
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)
Trivially "partitions" the mesh for one processor.
Definition: partitioner.C:151
static void set_node_processor_ids(MeshBase &mesh)
This function is called after partitioning to set the processor IDs for the nodes.
Definition: partitioner.C:416
MeshBase & mesh
virtual void _do_repartition(MeshBase &mesh, const unsigned int n)
This is the actual re-partitioning method which can be overridden in derived classes.
Definition: partitioner.h:204
static void partition_unpartitioned_elements(MeshBase &mesh)
These functions assign processor IDs to newly-created elements (in parallel) which are currently assi...
Definition: partitioner.C:175
static void set_parent_processor_ids(MeshBase &mesh)
This function is called after partitioning to set the processor IDs for the inactive parent elements...
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)
Repartitions the MeshBase into n parts.
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 416 of file partitioner.C.

References libMesh::MeshBase::active_element_ptr_range(), libMesh::ParallelObject::comm(), libMesh::DofObject::invalid_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_ptr_range(), libMesh::MeshBase::node_ref(), 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().

417 {
418  LOG_SCOPE("set_node_processor_ids()","Partitioner");
419 
420  // This function must be run on all processors at once
421  libmesh_parallel_only(mesh.comm());
422 
423  // If we have any unpartitioned elements at this
424  // stage there is a problem
425  libmesh_assert (MeshTools::n_elem(mesh.unpartitioned_elements_begin(),
426  mesh.unpartitioned_elements_end()) == 0);
427 
428 
429  // const dof_id_type orig_n_local_nodes = mesh.n_local_nodes();
430 
431  // libMesh::err << "[" << mesh.processor_id() << "]: orig_n_local_nodes="
432  // << orig_n_local_nodes << std::endl;
433 
434  // Build up request sets. Each node is currently owned by a processor because
435  // it is connected to an element owned by that processor. However, during the
436  // repartitioning phase that element may have been assigned a new processor id, but
437  // it is still resident on the original processor. We need to know where to look
438  // for new ids before assigning new ids, otherwise we may be asking the wrong processors
439  // for the wrong information.
440  //
441  // The only remaining issue is what to do with unpartitioned nodes. Since they are required
442  // to live on all processors we can simply rely on ourselves to number them properly.
443  std::vector<std::vector<dof_id_type>>
444  requested_node_ids(mesh.n_processors());
445 
446  // Loop over all the nodes, count the ones on each processor. We can skip ourself
447  std::vector<dof_id_type> ghost_nodes_from_proc(mesh.n_processors(), 0);
448 
449  for (auto & node : mesh.node_ptr_range())
450  {
451  libmesh_assert(node);
452  const processor_id_type current_pid = node->processor_id();
453  if (current_pid != mesh.processor_id() &&
454  current_pid != DofObject::invalid_processor_id)
455  {
456  libmesh_assert_less (current_pid, ghost_nodes_from_proc.size());
457  ghost_nodes_from_proc[current_pid]++;
458  }
459  }
460 
461  // We know how many objects live on each processor, so reserve()
462  // space for each.
463  for (processor_id_type pid=0; pid != mesh.n_processors(); ++pid)
464  requested_node_ids[pid].reserve(ghost_nodes_from_proc[pid]);
465 
466  // We need to get the new pid for each node from the processor
467  // which *currently* owns the node. We can safely skip ourself
468  for (auto & node : mesh.node_ptr_range())
469  {
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, requested_node_ids.size());
476  libmesh_assert_less (requested_node_ids[current_pid].size(),
477  ghost_nodes_from_proc[current_pid]);
478  requested_node_ids[current_pid].push_back(node->id());
479  }
480 
481  // Unset any previously-set node processor ids
482  node->invalidate_processor_id();
483  }
484 
485  // Loop over all the active elements
486  for (auto & elem : mesh.active_element_ptr_range())
487  {
488  libmesh_assert(elem);
489 
490  libmesh_assert_not_equal_to (elem->processor_id(), DofObject::invalid_processor_id);
491 
492  // For each node, set the processor ID to the min of
493  // its current value and this Element's processor id.
494  //
495  // TODO: we would probably get better parallel partitioning if
496  // we did something like "min for even numbered nodes, max for
497  // odd numbered". We'd need to be careful about how that would
498  // affect solution ordering for I/O, though.
499  for (unsigned int n=0; n<elem->n_nodes(); ++n)
500  elem->node_ptr(n)->processor_id() = std::min(elem->node_ptr(n)->processor_id(),
501  elem->processor_id());
502  }
503 
504  // And loop over the subactive elements, but don't reassign
505  // nodes that are already active on another processor.
506  MeshBase::element_iterator sub_it = mesh.subactive_elements_begin();
507  const MeshBase::element_iterator sub_end = mesh.subactive_elements_end();
508 
509  for ( ; sub_it != sub_end; ++sub_it)
510  {
511  Elem * elem = *sub_it;
512  libmesh_assert(elem);
513 
514  libmesh_assert_not_equal_to (elem->processor_id(), DofObject::invalid_processor_id);
515 
516  for (unsigned int n=0; n<elem->n_nodes(); ++n)
517  if (elem->node_ptr(n)->processor_id() == DofObject::invalid_processor_id)
518  elem->node_ptr(n)->processor_id() = elem->processor_id();
519  }
520 
521  // Same for the inactive elements -- we will have already gotten most of these
522  // nodes, *except* for the case of a parent with a subset of children which are
523  // ghost elements. In that case some of the parent nodes will not have been
524  // properly handled yet
525  MeshBase::element_iterator not_it = mesh.not_active_elements_begin();
526  const MeshBase::element_iterator not_end = mesh.not_active_elements_end();
527 
528  for ( ; not_it != not_end; ++not_it)
529  {
530  Elem * elem = *not_it;
531  libmesh_assert(elem);
532 
533  libmesh_assert_not_equal_to (elem->processor_id(), DofObject::invalid_processor_id);
534 
535  for (unsigned int n=0; n<elem->n_nodes(); ++n)
536  if (elem->node_ptr(n)->processor_id() == DofObject::invalid_processor_id)
537  elem->node_ptr(n)->processor_id() = elem->processor_id();
538  }
539 
540  // We can't assert that all nodes are connected to elements, because
541  // a DistributedMesh with NodeConstraints might have pulled in some
542  // remote nodes solely for evaluating those constraints.
543  // MeshTools::libmesh_assert_connected_nodes(mesh);
544 
545  // For such nodes, we'll do a sanity check later when making sure
546  // that we successfully reset their processor ids to something
547  // valid.
548 
549  // Next set node ids from other processors, excluding self
550  for (processor_id_type p=1; p != mesh.n_processors(); ++p)
551  {
552  // Trade my requests with processor procup and procdown
553  processor_id_type procup = cast_int<processor_id_type>
554  ((mesh.processor_id() + p) % mesh.n_processors());
555  processor_id_type procdown = cast_int<processor_id_type>
556  ((mesh.n_processors() + mesh.processor_id() - p) %
557  mesh.n_processors());
558  std::vector<dof_id_type> request_to_fill;
559  mesh.comm().send_receive(procup, requested_node_ids[procup],
560  procdown, request_to_fill);
561 
562  // Fill those requests in-place
563  for (std::size_t i=0; i != request_to_fill.size(); ++i)
564  {
565  Node & node = mesh.node_ref(request_to_fill[i]);
566  const processor_id_type new_pid = node.processor_id();
567 
568  // We may have an invalid processor_id() on nodes that have been
569  // "detached" from coarsened-away elements but that have not yet
570  // themselves been removed.
571  // libmesh_assert_not_equal_to (new_pid, DofObject::invalid_processor_id);
572  // libmesh_assert_less (new_pid, mesh.n_partitions()); // this is the correct test --
573  request_to_fill[i] = new_pid; // the number of partitions may
574  } // not equal the number of processors
575 
576  // Trade back the results
577  std::vector<dof_id_type> filled_request;
578  mesh.comm().send_receive(procdown, request_to_fill,
579  procup, filled_request);
580  libmesh_assert_equal_to (filled_request.size(), requested_node_ids[procup].size());
581 
582  // And copy the id changes we've now been informed of
583  for (std::size_t i=0; i != filled_request.size(); ++i)
584  {
585  Node & node = mesh.node_ref(requested_node_ids[procup][i]);
586 
587  // this is the correct test -- the number of partitions may
588  // not equal the number of processors
589 
590  // But: we may have an invalid processor_id() on nodes that
591  // have been "detached" from coarsened-away elements but
592  // that have not yet themselves been removed.
593  // libmesh_assert_less (filled_request[i], mesh.n_partitions());
594 
595  node.processor_id(cast_int<processor_id_type>(filled_request[i]));
596  }
597  }
598 
599 #ifdef DEBUG
600  MeshTools::libmesh_assert_valid_procids<Node>(mesh);
601 #endif
602 }
dof_id_type n_elem(const MeshBase::const_element_iterator &begin, const MeshBase::const_element_iterator &end)
Count up the number of elements of a specific type (as defined by an iterator range).
Definition: mesh_tools.C:656
MeshBase & mesh
uint8_t processor_id_type
Definition: id_types.h:99
libmesh_assert(j)
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:335
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_element_ptr_range(), libMesh::Elem::active_family_tree(), libMesh::MeshBase::ancestor_elements_begin(), libMesh::MeshBase::ancestor_elements_end(), libMesh::Elem::child_ref_range(), libMesh::ParallelObject::comm(), libMesh::Partitioner::communication_blocksize, libMesh::DofObject::id(), libMesh::DofObject::invalid_processor_id, libMesh::DofObject::invalidate_processor_id(), libMesh::MeshBase::is_serial(), libMesh::libmesh_assert(), libMesh::libmesh_ignore(), libMesh::MeshBase::max_elem_id(), std::min(), libMesh::Parallel::Communicator::min(), 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  for (auto & child : mesh.active_element_ptr_range())
276  {
277  // First set descendents
278  std::vector<const Elem *> subactive_family;
279  child->total_family_tree(subactive_family);
280  for (std::size_t i = 0; i != subactive_family.size(); ++i)
281  const_cast<Elem *>(subactive_family[i])->processor_id() = child->processor_id();
282 
283  // Then set ancestors
284  Elem * parent = child->parent();
285 
286  while (parent)
287  {
288  // invalidate the parent id, otherwise the min below
289  // will not work if the current parent id is less
290  // than all the children!
291  parent->invalidate_processor_id();
292 
293  for (auto & child : parent->child_ref_range())
294  {
295  libmesh_assert(!child.is_remote());
296  libmesh_assert_not_equal_to (child.processor_id(), DofObject::invalid_processor_id);
297  parent->processor_id() = std::min(parent->processor_id(),
298  child.processor_id());
299  }
300  parent = parent->parent();
301  }
302  }
303  }
304 
305  // When the mesh is parallel we cannot guarantee that parents have access to
306  // all their children.
307  else
308  {
309  // Setting subactive processor ids is easy: we can guarantee
310  // that children have access to all their parents.
311 
312  // Loop over all the active elements in the mesh
313  for (auto & child : mesh.active_element_ptr_range())
314  {
315  std::vector<const Elem *> subactive_family;
316  child->total_family_tree(subactive_family);
317  for (std::size_t i = 0; i != subactive_family.size(); ++i)
318  const_cast<Elem *>(subactive_family[i])->processor_id() = child->processor_id();
319  }
320 
321  // When the mesh is parallel we cannot guarantee that parents have access to
322  // all their children.
323 
324  // We will use a brute-force approach here. Each processor finds its parent
325  // elements and sets the parent pid to the minimum of its
326  // semilocal descendants.
327  // A global reduction is then performed to make sure the true minimum is found.
328  // As noted, this is required because we cannot guarantee that a parent has
329  // access to all its children on any single processor.
330  libmesh_parallel_only(mesh.comm());
331  libmesh_assert(MeshTools::n_elem(mesh.unpartitioned_elements_begin(),
332  mesh.unpartitioned_elements_end()) == 0);
333 
334  const dof_id_type max_elem_id = mesh.max_elem_id();
335 
336  std::vector<processor_id_type>
337  parent_processor_ids (std::min(communication_blocksize,
338  max_elem_id));
339 
340  for (dof_id_type blk=0, last_elem_id=0; last_elem_id<max_elem_id; blk++)
341  {
342  last_elem_id =
343  std::min(static_cast<dof_id_type>((blk+1)*communication_blocksize),
344  max_elem_id);
345  const dof_id_type first_elem_id = blk*communication_blocksize;
346 
347  std::fill (parent_processor_ids.begin(),
348  parent_processor_ids.end(),
350 
351  // first build up local contributions to parent_processor_ids
352  MeshBase::element_iterator not_it = mesh.ancestor_elements_begin();
353  const MeshBase::element_iterator not_end = mesh.ancestor_elements_end();
354 
355  bool have_parent_in_block = false;
356 
357  for ( ; not_it != not_end; ++not_it)
358  {
359  Elem * parent = *not_it;
360 
361  const dof_id_type parent_idx = parent->id();
362  libmesh_assert_less (parent_idx, max_elem_id);
363 
364  if ((parent_idx >= first_elem_id) &&
365  (parent_idx < last_elem_id))
366  {
367  have_parent_in_block = true;
369 
370  std::vector<const Elem *> active_family;
371  parent->active_family_tree(active_family);
372  for (std::size_t i = 0; i != active_family.size(); ++i)
373  parent_pid = std::min (parent_pid, active_family[i]->processor_id());
374 
375  const dof_id_type packed_idx = parent_idx - first_elem_id;
376  libmesh_assert_less (packed_idx, parent_processor_ids.size());
377 
378  parent_processor_ids[packed_idx] = parent_pid;
379  }
380  }
381 
382  // then find the global minimum
383  mesh.comm().min (parent_processor_ids);
384 
385  // and assign the ids, if we have a parent in this block.
386  if (have_parent_in_block)
387  for (not_it = mesh.ancestor_elements_begin();
388  not_it != not_end; ++not_it)
389  {
390  Elem * parent = *not_it;
391 
392  const dof_id_type parent_idx = parent->id();
393 
394  if ((parent_idx >= first_elem_id) &&
395  (parent_idx < last_elem_id))
396  {
397  const dof_id_type packed_idx = parent_idx - first_elem_id;
398  libmesh_assert_less (packed_idx, parent_processor_ids.size());
399 
400  const processor_id_type parent_pid =
401  parent_processor_ids[packed_idx];
402 
403  libmesh_assert_not_equal_to (parent_pid, DofObject::invalid_processor_id);
404 
405  parent->processor_id() = parent_pid;
406  }
407  }
408  }
409  }
410 
411 #endif // LIBMESH_ENABLE_AMR
412 }
dof_id_type n_elem(const MeshBase::const_element_iterator &begin, const MeshBase::const_element_iterator &end)
Count up the number of elements of a specific type (as defined by an iterator range).
Definition: mesh_tools.C:656
MeshBase & mesh
uint8_t processor_id_type
Definition: id_types.h:99
libmesh_assert(j)
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:335
void libmesh_ignore(const T &)
static const dof_id_type communication_blocksize
The blocksize to use when doing blocked parallel communication.
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)
Slightly generalized version of single_partition which acts on a range of elements defined by the pai...
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
Also have a polymorphic pointer to the end object, this prevents iterating past the 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: