21 #include "libmesh/parmetis_partitioner.h" 24 #include "libmesh/libmesh_config.h" 25 #include "libmesh/elem.h" 26 #include "libmesh/enum_partitioner_type.h" 27 #include "libmesh/libmesh_logging.h" 28 #include "libmesh/mesh_base.h" 29 #include "libmesh/mesh_serializer.h" 30 #include "libmesh/mesh_tools.h" 31 #include "libmesh/mesh_communication.h" 32 #include "libmesh/metis_partitioner.h" 33 #include "libmesh/parallel_only.h" 34 #include "libmesh/parmetis_helper.h" 41 #ifdef LIBMESH_HAVE_PARMETIS 49 # include "libmesh/ignore_warnings.h" 50 # include "parmetis.h" 51 # include "libmesh/restore_warnings.h" 59 #include <unordered_map> 67 #ifdef LIBMESH_HAVE_PARMETIS 74 #ifdef LIBMESH_HAVE_PARMETIS 75 : _pmetis(std::make_unique<ParmetisHelper>())
83 #ifdef LIBMESH_HAVE_PARMETIS
84 , _pmetis(
std::make_unique<ParmetisHelper>(*(other._pmetis)))
91 ParmetisPartitioner::~ParmetisPartitioner() =
default;
103 const unsigned int n_sbdmns)
105 this->_do_repartition (
mesh, n_sbdmns);
111 const unsigned int n_sbdmns)
114 libmesh_parallel_only(
mesh.comm());
122 this->single_partition(
mesh);
126 libmesh_assert_greater (n_sbdmns, 0);
129 #ifndef LIBMESH_HAVE_PARMETIS 132 libMesh::out <<
"ERROR: The library has been built without" << std::endl
133 <<
"Parmetis support. Using a Metis" << std::endl
134 <<
"partitioner instead!" << std::endl;);
145 mp.partition_range (
mesh,
mesh.active_elements_begin(),
146 mesh.active_elements_end(), n_sbdmns);
152 if (
mesh.n_processors() == 1)
162 mesh.active_elements_end(), n_sbdmns);
166 LOG_SCOPE(
"repartition()",
"ParmetisPartitioner");
172 this->build_graph (
mesh);
179 bool ready_for_parmetis =
true;
180 for (
const auto & nelem : _n_active_elem_on_proc)
182 ready_for_parmetis =
false;
184 std::size_t my_adjacency = _pmetis->adjncy.size();
185 mesh.comm().min(my_adjacency);
187 ready_for_parmetis =
false;
192 if (!ready_for_parmetis)
204 std::vector<Parmetis::idx_t> vsize(_pmetis->vwgt.size(), 1);
205 Parmetis::real_t itr = 1000000.0;
206 MPI_Comm mpi_comm =
mesh.comm().get();
211 Parmetis::ParMETIS_V3_AdaptiveRepart(_pmetis->vtxdist.empty() ? nullptr : _pmetis->vtxdist.data(),
212 _pmetis->xadj.empty() ? nullptr : _pmetis->xadj.data(),
213 _pmetis->adjncy.empty() ? nullptr : _pmetis->adjncy.data(),
214 _pmetis->vwgt.empty() ? nullptr : _pmetis->vwgt.data(),
215 vsize.empty() ? nullptr : vsize.data(),
221 _pmetis->tpwgts.empty() ? nullptr : _pmetis->tpwgts.data(),
222 _pmetis->ubvec.empty() ? nullptr : _pmetis->ubvec.data(),
224 _pmetis->options.data(),
226 _pmetis->part.empty() ? nullptr :
reinterpret_cast<Parmetis::idx_t *
>(_pmetis->part.data()),
230 this->assign_partitioning (
mesh, _pmetis->part);
232 #endif // #ifndef LIBMESH_HAVE_PARMETIS ... else ... 239 #ifdef LIBMESH_HAVE_PARMETIS 242 const unsigned int n_sbdmns)
244 LOG_SCOPE(
"initialize()",
"ParmetisPartitioner");
248 _pmetis->wgtflag = 2;
250 _pmetis->numflag = 0;
251 _pmetis->nparts =
static_cast<Parmetis::idx_t
>(n_sbdmns);
252 _pmetis->edgecut = 0;
256 _pmetis->vtxdist.assign (
mesh.n_processors()+1, 0);
257 _pmetis->tpwgts.assign (_pmetis->nparts, 1./_pmetis->nparts);
258 _pmetis->ubvec.assign (_pmetis->ncon, 1.05);
259 _pmetis->part.assign (n_active_local_elem, 0);
260 _pmetis->options.resize (5);
261 _pmetis->vwgt.resize (n_active_local_elem);
264 _pmetis->options[0] = 1;
265 _pmetis->options[1] = 0;
266 _pmetis->options[2] = 15;
267 _pmetis->options[3] = 2;
281 _find_global_index_by_pid_map(
mesh);
293 libmesh_assert_equal_to (_pmetis->vtxdist.size(),
294 cast_int<std::size_t>(
mesh.n_processors()+1));
295 libmesh_assert_equal_to (_pmetis->vtxdist[0], 0);
299 _pmetis->vtxdist[pid+1] = _pmetis->vtxdist[pid] + _n_active_elem_on_proc[pid];
300 n_active_elem += _n_active_elem_on_proc[pid];
302 libmesh_assert_equal_to (_pmetis->vtxdist.back(),
static_cast<Parmetis::idx_t
>(n_active_elem));
307 std::unordered_map<dof_id_type, dof_id_type> global_index_map;
310 std::vector<dof_id_type> global_index;
325 const Elem * elem = *it;
329 libmesh_assert_less (cnt, global_index.size());
330 libmesh_assert_less (global_index[cnt], n_active_elem);
332 global_index_map.emplace(elem->
id(), global_index[cnt++]);
336 libmesh_assert_less_equal (global_index_map.size(), n_active_elem);
337 libmesh_assert_less_equal (_global_index_by_pid_map.size(), n_active_elem);
343 libmesh_error_msg_if(global_index_map.size() != _global_index_by_pid_map.size(),
344 "ERROR: ParmetisPartitioner cannot handle unpartitioned objects!");
351 std::vector<dof_id_type> subdomain_bounds(
mesh.n_processors());
353 const dof_id_type first_local_elem = _pmetis->vtxdist[
mesh.processor_id()];
360 if (pid < static_cast<unsigned int>(_pmetis->nparts))
362 tgt_subdomain_size = n_active_elem/std::min
363 (cast_int<Parmetis::idx_t>(
mesh.n_processors()), _pmetis->nparts);
365 if (pid < n_active_elem%_pmetis->nparts)
366 tgt_subdomain_size++;
369 subdomain_bounds[0] = tgt_subdomain_size;
371 subdomain_bounds[pid] = subdomain_bounds[pid-1] + tgt_subdomain_size;
374 libmesh_assert_equal_to (subdomain_bounds.back(), n_active_elem);
376 for (
const auto & elem :
mesh.active_local_element_ptr_range())
380 _global_index_by_pid_map[elem->id()];
381 libmesh_assert_less (global_index_by_pid, n_active_elem);
384 global_index_by_pid - first_local_elem;
386 libmesh_assert_less (local_index, n_active_local_elem);
387 libmesh_assert_less (local_index, _pmetis->vwgt.size());
396 _pmetis->vwgt[local_index] = 50;
398 _pmetis->vwgt[local_index] = elem->n_nodes();
403 global_index_map[elem->id()];
405 libmesh_assert_less (global_index, subdomain_bounds.back());
407 const unsigned int subdomain_id =
408 cast_int<unsigned int>
410 std::lower_bound(subdomain_bounds.begin(),
411 subdomain_bounds.end(),
413 libmesh_assert_less (subdomain_id, _pmetis->nparts);
414 libmesh_assert_less (local_index, _pmetis->part.size());
416 _pmetis->part[local_index] = subdomain_id;
425 LOG_SCOPE(
"build_graph()",
"ParmetisPartitioner");
432 Partitioner::build_graph(
mesh);
436 for (
auto & row: _dual_graph)
437 graph_size += cast_int<dof_id_type>(row.size());
440 _pmetis->xadj.clear();
441 _pmetis->xadj.reserve (n_active_local_elem + 1);
442 _pmetis->adjncy.clear();
443 _pmetis->adjncy.reserve (graph_size);
445 for (
auto & graph_row : _dual_graph)
447 _pmetis->xadj.push_back(cast_int<int>(_pmetis->adjncy.size()));
448 _pmetis->adjncy.insert(_pmetis->adjncy.end(),
454 _pmetis->xadj.push_back(cast_int<int>(_pmetis->adjncy.size()));
456 libmesh_assert_equal_to (_pmetis->xadj.size(), n_active_local_elem+1);
457 libmesh_assert_equal_to (_pmetis->adjncy.size(), graph_size);
460 #endif // #ifdef LIBMESH_HAVE_PARMETIS void find_global_indices(const Parallel::Communicator &communicator, const libMesh::BoundingBox &, const ForwardIterator &, const ForwardIterator &, std::vector< dof_id_type > &) const
This method determines a globally unique, partition-agnostic index for each object in the input range...
The definition of the const_element_iterator struct.
virtual void partition_range(MeshBase &mesh, MeshBase::element_iterator it, MeshBase::element_iterator end, const unsigned int n) override
Called by the SubdomainPartitioner to partition elements in the range (it, end).
This is the base class from which all geometric element types are derived.
ParmetisPartitioner()
Default and copy ctors.
The MetisPartitioner uses the Metis graph partitioner to partition the elements.
The libMesh namespace provides an interface to certain functionality in the library.
void initialize(EquationSystems &es, const std::string &system_name)
Real distance(const Point &p)
This is the MeshBase class.
PartitionerType
Defines an enum for mesh partitioner types.
This is the MeshCommunication class.
Defines a Cartesian bounding box by the two corner extremum.
const unsigned int MIN_ELEM_PER_PROC
Temporarily serialize a DistributedMesh for non-distributed-mesh capable code paths.
virtual void partition(MeshBase &mesh, const unsigned int n)
Partitions the MeshBase into n parts by setting processor_id() on Nodes and Elems.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...