From: Peter Munch Date: Mon, 6 Dec 2021 21:26:29 +0000 (+0100) Subject: p::f::T: allow repartitioning X-Git-Tag: v9.4.0-rc1~738^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cc43200461a4d0d682877a9177b981b97a84aeb2;p=dealii.git p::f::T: allow repartitioning --- diff --git a/include/deal.II/distributed/fully_distributed_tria.h b/include/deal.II/distributed/fully_distributed_tria.h index 2de0ba3d75..318e99f84e 100644 --- a/include/deal.II/distributed/fully_distributed_tria.h +++ b/include/deal.II/distributed/fully_distributed_tria.h @@ -39,6 +39,12 @@ namespace GridTools template struct PeriodicFacePair; } + +namespace RepartitioningPolicyTools +{ + template + class Base; +} #endif namespace parallel @@ -196,6 +202,20 @@ namespace parallel const unsigned int)> &partitioner, const TriangulationDescription::Settings & settings); + /** + * TODO + */ + void + set_partitioner( + const RepartitioningPolicyTools::Base &partitioner, + const TriangulationDescription::Settings & settings); + + /** + * TODO + */ + void + repartition(); + /** * Coarsen and refine the mesh according to refinement and coarsening * flags set. @@ -301,6 +321,12 @@ namespace parallel const unsigned int)> partitioner; + /** + * TODO + */ + SmartPointer> + partitioner_distributed; + /** * Sorted list of pairs of coarse-cell ids and their indices. */ diff --git a/include/deal.II/distributed/repartitioning_policy_tools.h b/include/deal.II/distributed/repartitioning_policy_tools.h index 81d412c45c..8240a72086 100644 --- a/include/deal.II/distributed/repartitioning_policy_tools.h +++ b/include/deal.II/distributed/repartitioning_policy_tools.h @@ -46,7 +46,7 @@ namespace RepartitioningPolicyTools * See the description of RepartitioningPolicyTools for more information. */ template - class Base + class Base : public Subscriptor { public: /** diff --git a/source/distributed/fully_distributed_tria.cc b/source/distributed/fully_distributed_tria.cc index 098f7951fc..2a7496d5ec 100644 --- a/source/distributed/fully_distributed_tria.cc +++ b/source/distributed/fully_distributed_tria.cc @@ -18,6 +18,7 @@ #include #include +#include #include @@ -310,6 +311,37 @@ namespace parallel + template + void + Triangulation::set_partitioner( + const RepartitioningPolicyTools::Base &partitioner, + const TriangulationDescription::Settings & settings) + { + this->partitioner_distributed = &partitioner; + this->settings = settings; + } + + + + template + void + Triangulation::repartition() + { + this->signals.pre_distributed_repartition(); + + const auto construction_data = TriangulationDescription::Utilities:: + create_description_from_triangulation( + *this, + this->partitioner_distributed->partition(*this), + this->settings); + + this->create_triangulation(construction_data); + + this->signals.post_distributed_repartition(); + } + + + template void Triangulation::execute_coarsening_and_refinement() diff --git a/source/grid/tria_description.cc b/source/grid/tria_description.cc index 3bb3e28a01..1606c0df64 100644 --- a/source/grid/tria_description.cc +++ b/source/grid/tria_description.cc @@ -17,6 +17,7 @@ #include #include +#include #include #include @@ -109,7 +110,8 @@ namespace TriangulationDescription collect( const std::vector & relevant_processes, const std::vector> &description_temp, - const MPI_Comm & comm) + const MPI_Comm & comm, + const bool vertices_have_unique_ids) { dealii::Utilities::MPI::ConsensusAlgorithms::AnonymousProcess @@ -135,7 +137,8 @@ namespace TriangulationDescription std::vector &) { this->merge( dealii::Utilities::unpack>( - recv_buffer, false)); + recv_buffer, false), + vertices_have_unique_ids); }); dealii::Utilities::MPI::ConsensusAlgorithms::Selector( @@ -149,17 +152,66 @@ namespace TriangulationDescription * This done by the reduce() function. */ void - merge(const DescriptionTemp &other) + merge(const DescriptionTemp &other, + const bool vertices_have_unique_ids) { this->cell_infos.resize( std::max(other.cell_infos.size(), this->cell_infos.size())); - this->coarse_cells.insert(this->coarse_cells.end(), - other.coarse_cells.begin(), - other.coarse_cells.end()); - this->coarse_cell_vertices.insert(this->coarse_cell_vertices.end(), - other.coarse_cell_vertices.begin(), - other.coarse_cell_vertices.end()); + if (vertices_have_unique_ids == false) + { + const auto comp = [](const auto &a, const auto &b) { + const double tolerance = 1e-10; + + if (a.distance(b) < tolerance) + return false; + + for (unsigned int i = 0; i < spacedim; ++i) + if (std::abs(a[i] - b[i]) > tolerance) + return a[i] < b[i]; + + return false; + }; + + std::map, unsigned int, decltype(comp)> map(comp); + std::map map_i; + + for (unsigned int i = 0; i < this->coarse_cell_vertices.size(); + ++i) + map[coarse_cell_vertices[i].second] = i; + + unsigned int counter = coarse_cell_vertices.size(); + + for (const auto p : other.coarse_cell_vertices) + if (map.find(p.second) == map.end()) + { + this->coarse_cell_vertices.emplace_back(counter, p.second); + map[p.second] = map_i[p.first] = counter++; + } + else + map_i[p.first] = map[p.second]; + + for (auto cell : other.coarse_cells) + { + for (auto &v : cell.vertices) + v = map_i[v]; + } + + this->coarse_cells.insert(this->coarse_cells.end(), + other.coarse_cells.begin(), + other.coarse_cells.end()); + } + else + { + this->coarse_cells.insert(this->coarse_cells.end(), + other.coarse_cells.begin(), + other.coarse_cells.end()); + this->coarse_cell_vertices.insert( + this->coarse_cell_vertices.end(), + other.coarse_cell_vertices.begin(), + other.coarse_cell_vertices.end()); + } + this->coarse_cell_index_to_coarse_cell_id.insert( this->coarse_cell_index_to_coarse_cell_id.end(), other.coarse_cell_index_to_coarse_cell_id.begin(), @@ -1025,9 +1077,13 @@ namespace TriangulationDescription // collect description from all processes that used to own locally-owned // active cells of this process in a single description DescriptionTemp description_merged; - description_merged.collect(relevant_processes, - descriptions_per_rank, - partition.get_mpi_communicator()); + description_merged.collect( + relevant_processes, + descriptions_per_rank, + partition.get_mpi_communicator(), + dynamic_cast< + const parallel::fullydistributed::Triangulation *>( + &tria) != nullptr); // remove redundant entries description_merged.reduce();