From 3ecc8b27d18f070d6d22f67642ecacf1877c3bf7 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Sat, 9 Oct 2021 19:40:44 +0200 Subject: [PATCH] Extend RepartitioningPolicyTools::DefaultPolicy --- .../distributed/repartitioning_policy_tools.h | 11 +++ .../repartitioning_policy_tools.cc | 53 +++++++++- .../repartitioning_07.cc | 96 +++++++++++++++++++ ...tioning_07.with_p4est=true.mpirun=4.output | 11 +++ 4 files changed, 169 insertions(+), 2 deletions(-) create mode 100644 tests/fullydistributed_grids/repartitioning_07.cc create mode 100644 tests/fullydistributed_grids/repartitioning_07.with_p4est=true.mpirun=4.output diff --git a/include/deal.II/distributed/repartitioning_policy_tools.h b/include/deal.II/distributed/repartitioning_policy_tools.h index 3bdcf77945..ff283c56fd 100644 --- a/include/deal.II/distributed/repartitioning_policy_tools.h +++ b/include/deal.II/distributed/repartitioning_policy_tools.h @@ -67,9 +67,20 @@ namespace RepartitioningPolicyTools class DefaultPolicy : public Base { public: + /** + * Constructor. + * + * @param tighten allows to renumber of subdomains so that empty ranks are + * positioned at the end. + */ + DefaultPolicy(const bool tighten = false); + virtual LinearAlgebra::distributed::Vector partition( const Triangulation &tria_coarse_in) const override; + + private: + const bool tighten; }; /** diff --git a/source/distributed/repartitioning_policy_tools.cc b/source/distributed/repartitioning_policy_tools.cc index 92d39e215c..9b203f93c2 100644 --- a/source/distributed/repartitioning_policy_tools.cc +++ b/source/distributed/repartitioning_policy_tools.cc @@ -46,12 +46,61 @@ namespace RepartitioningPolicyTools } } // namespace + + template + DefaultPolicy::DefaultPolicy(const bool tighten) + : tighten(tighten) + {} + template LinearAlgebra::distributed::Vector DefaultPolicy::partition( - const Triangulation &) const + const Triangulation &tria_in) const { - return {}; // nothing to do + if (tighten == false) + return {}; // nothing to do + + const auto tria = + dynamic_cast *>( + &tria_in); + + if (tria == nullptr) + return {}; // nothing to do, since serial triangulation + +#ifndef DEAL_II_WITH_MPI + Assert(false, ExcNeedsMPI()); + return {}; +#else + + const auto comm = tria->get_communicator(); + + const unsigned int process_has_active_locally_owned_cells = + tria->n_locally_owned_active_cells() > 0; + const unsigned int n_processes_with_active_locally_owned_cells = + Utilities::MPI::sum(process_has_active_locally_owned_cells, comm); + + if (n_processes_with_active_locally_owned_cells == + Utilities::MPI::n_mpi_processes(comm)) + return {}; // nothing to do, since all processes have cells + + unsigned int offset = 0; + + const int ierr = MPI_Exscan(&process_has_active_locally_owned_cells, + &offset, + 1, + Utilities::MPI::internal::mpi_type_id( + &process_has_active_locally_owned_cells), + MPI_SUM, + comm); + AssertThrowMPI(ierr); + + LinearAlgebra::distributed::Vector partition( + tria->global_active_cell_index_partitioner().lock()); + + partition = offset; + + return partition; +#endif } diff --git a/tests/fullydistributed_grids/repartitioning_07.cc b/tests/fullydistributed_grids/repartitioning_07.cc new file mode 100644 index 0000000000..091107ad7e --- /dev/null +++ b/tests/fullydistributed_grids/repartitioning_07.cc @@ -0,0 +1,96 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +// Test RepartitioningPolicyTools::DefaultPolicy's tightening functionality. + +#include + +#include +#include +#include + +#include + +#include +#include +#include +#include + +#include + +#include + +#include "./tests.h" + +using namespace dealii; + +template +void +test(const MPI_Comm comm) +{ + parallel::shared::Triangulation tria( + comm, + Triangulation::none, + true, + parallel::shared::Triangulation::partition_custom_signal); + + tria.signals.create.connect([&tria]() { + for (const auto &cell : tria.active_cell_iterators()) + { + switch (cell->active_cell_index()) + { + case 0: + cell->set_subdomain_id(0); + break; + case 1: + case 2: + cell->set_subdomain_id(2); + break; + case 3: + cell->set_subdomain_id(3); + break; + } + } + }); + + GridGenerator::subdivided_hyper_cube(tria, 2); + + const auto construction_data = + TriangulationDescription::Utilities::create_description_from_triangulation( + tria, + RepartitioningPolicyTools::DefaultPolicy(true).partition(tria)); + + parallel::fullydistributed::Triangulation tria_pft(comm); + tria_pft.create_triangulation(construction_data); + + deallog << tria.n_locally_owned_active_cells() << " " + << tria_pft.n_locally_owned_active_cells() << std::endl; +} + + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1); + MPILogInitAll all; + + MPI_Comm comm = MPI_COMM_WORLD; + + deallog.push("all"); + test<2>(comm); + deallog.pop(); +} diff --git a/tests/fullydistributed_grids/repartitioning_07.with_p4est=true.mpirun=4.output b/tests/fullydistributed_grids/repartitioning_07.with_p4est=true.mpirun=4.output new file mode 100644 index 0000000000..17d8de7aac --- /dev/null +++ b/tests/fullydistributed_grids/repartitioning_07.with_p4est=true.mpirun=4.output @@ -0,0 +1,11 @@ + +DEAL:0:all::1 1 + +DEAL:1:all::0 2 + + +DEAL:2:all::2 1 + + +DEAL:3:all::1 0 + -- 2.39.5