From 2fcbc7d9c3e5930ff81adc544ca007885e2eca50 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Thu, 21 May 2020 22:56:56 +0200 Subject: [PATCH] Bugfix: hp::Refinement::choose_p_over_h() in parallel --- doc/news/changes/minor/20200522Fehling | 3 + source/hp/refinement.cc | 91 +++++++++-- tests/sharedtria/hp_choose_p_over_h.cc | 152 ++++++++++++++++++ ...ose_p_over_h.with_mpi=true.mpirun=2.output | 7 + 4 files changed, 238 insertions(+), 15 deletions(-) create mode 100644 doc/news/changes/minor/20200522Fehling create mode 100644 tests/sharedtria/hp_choose_p_over_h.cc create mode 100644 tests/sharedtria/hp_choose_p_over_h.with_mpi=true.mpirun=2.output diff --git a/doc/news/changes/minor/20200522Fehling b/doc/news/changes/minor/20200522Fehling new file mode 100644 index 0000000000..b04788de08 --- /dev/null +++ b/doc/news/changes/minor/20200522Fehling @@ -0,0 +1,3 @@ +Bugfix: hp::Refinement::choose_p_over_h() now works in parallel. +
+(Marc Fehling, 2020/05/22) diff --git a/source/hp/refinement.cc b/source/hp/refinement.cc index a6191283d3..1b22490f9b 100644 --- a/source/hp/refinement.cc +++ b/source/hp/refinement.cc @@ -23,6 +23,7 @@ #include #include +#include #include #include @@ -635,6 +636,33 @@ namespace hp void choose_p_over_h(const hp::DoFHandler &dof_handler) { + // Siblings of cells to be coarsened may not be owned by the same + // processor. We will exchange coarsening flags on ghost cells and + // temporarily store them. + std::map> ghost_buffer; + + if (dynamic_cast *>( + &dof_handler.get_triangulation())) + { + auto pack = []( + const typename dealii::hp::DoFHandler:: + active_cell_iterator &cell) -> std::pair { + return {cell->coarsen_flag_set(), cell->future_fe_index_set()}; + }; + + auto unpack = [&ghost_buffer]( + const typename dealii::hp::DoFHandler:: + active_cell_iterator & cell, + const std::pair pair) -> void { + ghost_buffer.emplace(cell->id(), pair); + }; + + GridTools::exchange_cell_data_to_ghosts< + std::pair, + dealii::hp::DoFHandler>(dof_handler, pack, unpack); + } + + for (const auto &cell : dof_handler.active_cell_iterators()) if (cell->is_locally_owned() && cell->future_fe_index_set()) { @@ -656,27 +684,60 @@ namespace hp const auto &child = parent->child(child_index); if (child->is_active()) { - if (child->coarsen_flag_set()) - ++h_flagged_children; - if (child->future_fe_index_set()) - ++p_flagged_children; + if (child->is_locally_owned()) + { + if (child->coarsen_flag_set()) + ++h_flagged_children; + if (child->future_fe_index_set()) + ++p_flagged_children; + } + else if (child->is_ghost()) + { + const std::pair &flags = + ghost_buffer[child->id()]; + + if (flags.first) + ++h_flagged_children; + if (flags.second) + ++p_flagged_children; + } + else + { + // Siblings of locally owned cells are all + // either also locally owned or ghost cells. + Assert(false, ExcInternalError()); + } } } if (h_flagged_children == n_children && p_flagged_children != n_children) - // Perform pure h coarsening and - // drop all p adaptation flags. - for (unsigned int child_index = 0; child_index < n_children; - ++child_index) - parent->child(child_index)->clear_future_fe_index(); + { + // Perform pure h coarsening and + // drop all p adaptation flags. + for (unsigned int child_index = 0; child_index < n_children; + ++child_index) + { + const auto &child = parent->child(child_index); + // h_flagged_children == n_children implies + // that all children are active + Assert(child->is_active(), ExcInternalError()); + if (child->is_locally_owned()) + child->clear_future_fe_index(); + } + } else - // Perform p adaptation on all children and - // drop all h coarsening flags. - for (unsigned int child_index = 0; child_index < n_children; - ++child_index) - if (parent->child(child_index)->is_active()) - parent->child(child_index)->clear_coarsen_flag(); + { + // Perform p adaptation on all children and + // drop all h coarsening flags. + for (unsigned int child_index = 0; child_index < n_children; + ++child_index) + { + const auto &child = parent->child(child_index); + if (child->is_active() && child->is_locally_owned()) + child->clear_coarsen_flag(); + } + } } } } diff --git a/tests/sharedtria/hp_choose_p_over_h.cc b/tests/sharedtria/hp_choose_p_over_h.cc new file mode 100644 index 0000000000..0482b8a9c6 --- /dev/null +++ b/tests/sharedtria/hp_choose_p_over_h.cc @@ -0,0 +1,152 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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. +// +// --------------------------------------------------------------------- + + + +// hp::Refinement::choose_p_over_h called future_fe_index_set() on cells that +// are not locally owned and triggered an assertion at some point. + + +#include + +#include + +#include + +#include +#include +#include + +#include "../tests.h" + + +template +void +test() +{ + // setup + const unsigned int n_cells = 2; + + parallel::shared::Triangulation tr( + MPI_COMM_WORLD, + ::Triangulation::none, + false, + parallel::shared::Triangulation::partition_custom_signal); + tr.signals.post_refinement.connect([&tr]() { + // partition the triangulation by hand + for (const auto &cell : tr.active_cell_iterators()) + cell->set_subdomain_id(cell->active_cell_index() % + Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD)); + }); + + std::vector rep(dim, 1); + rep[0] = n_cells; + Point p1, p2; + for (unsigned int d = 0; d < dim; ++d) + { + p1[d] = 0; + p2[d] = (d == 0) ? n_cells : 1; + } + GridGenerator::subdivided_hyper_rectangle(tr, rep, p1, p2); + tr.refine_global(1); + + hp::FECollection fes; + for (unsigned int d = 1; d <= 2; ++d) + fes.push_back(FE_Q(d)); + + hp::DoFHandler dh(tr); + dh.set_fe(fes); + + // set flags + for (auto cell = dh.begin(0); cell != dh.end(0); ++cell) + { + if (cell->id().to_string() == "0_0:") + { + // all children will be flagged for both h- and p-adaptation. + // choose_p_over_h() will decide in favor of p-adaptation. + for (unsigned int i = 0; i < cell->n_children(); ++i) + { + const auto &child = cell->child(i); + if (child->is_locally_owned()) + { + child->set_future_fe_index(1); + child->set_coarsen_flag(); + } + } + } + else if (cell->id().to_string() == "1_0:") + { + // all children will be flagged for both h-adaptation + // and only one of them for p-adaptation. + // choose_p_over_h() will decide in favor of h-adaptation. + for (unsigned int i = 0; i < cell->n_children(); ++i) + { + const auto &child = cell->child(i); + if (child->is_locally_owned()) + { + if (i == 0) + child->set_future_fe_index(1); + child->set_coarsen_flag(); + } + } + } + } + + // decide between p and h flags + hp::Refinement::choose_p_over_h(dh); + + // verify + unsigned int h_flagged_cells = 0, p_flagged_cells = 0; + for (const auto &cell : dh.active_cell_iterators()) + if (cell->is_locally_owned()) + { + if (cell->coarsen_flag_set()) + ++h_flagged_cells; + if (cell->future_fe_index_set()) + ++p_flagged_cells; + } + const unsigned int global_h_flagged_cells = + Utilities::MPI::sum(h_flagged_cells, MPI_COMM_WORLD), + global_p_flagged_cells = + Utilities::MPI::sum(p_flagged_cells, MPI_COMM_WORLD); + + deallog << "h-flags:" << global_h_flagged_cells << std::endl + << "p-flags:" << global_p_flagged_cells << std::endl + << "OK" << std::endl; +} + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0) + { + initlog(); + + deallog.push("2d"); + test<2>(); + deallog.pop(); + deallog.push("3d"); + test<3>(); + deallog.pop(); + } + else + { + test<2>(); + test<3>(); + } +} diff --git a/tests/sharedtria/hp_choose_p_over_h.with_mpi=true.mpirun=2.output b/tests/sharedtria/hp_choose_p_over_h.with_mpi=true.mpirun=2.output new file mode 100644 index 0000000000..47a88e9794 --- /dev/null +++ b/tests/sharedtria/hp_choose_p_over_h.with_mpi=true.mpirun=2.output @@ -0,0 +1,7 @@ + +DEAL:2d::h-flags:4 +DEAL:2d::p-flags:4 +DEAL:2d::OK +DEAL:3d::h-flags:8 +DEAL:3d::p-flags:8 +DEAL:3d::OK -- 2.39.5