From fa5c0a0acf45deeb26c6d3f5fd2c03510af87d9d Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Mon, 8 Mar 2021 20:01:35 -0700 Subject: [PATCH] Support subfaces for DoFHandler::prepare_coarsening_and_refinement(). --- source/dofs/dof_handler.cc | 62 +++++-- .../prepare_coarsening_and_refinement_03.cc | 124 +++++++++++++ ...repare_coarsening_and_refinement_03.output | 4 + .../prepare_coarsening_and_refinement_03.cc | 166 ++++++++++++++++++ ...inement_03.with_p4est=true.mpirun=1.output | 13 ++ ...inement_03.with_p4est=true.mpirun=4.output | 13 ++ 6 files changed, 363 insertions(+), 19 deletions(-) create mode 100644 tests/hp/prepare_coarsening_and_refinement_03.cc create mode 100644 tests/hp/prepare_coarsening_and_refinement_03.output create mode 100644 tests/mpi/prepare_coarsening_and_refinement_03.cc create mode 100644 tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=1.output create mode 100644 tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=4.output diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc index 2e36ac73e2..9f07be99e7 100644 --- a/source/dofs/dof_handler.cc +++ b/source/dofs/dof_handler.cc @@ -2741,6 +2741,35 @@ DoFHandler::prepare_coarsening_and_refinement( // - always raise levels to match criterion, never lower them // - exchange level indices on ghost cells + // Function that updates the level of neighbor to fulfill difference + // criterion, and returns whether it was changed. + const auto update_neighbor_level = + [&future_levels, max_difference](const active_cell_iterator &neighbor, + const level_type cell_level) -> bool { + // We only care about locally owned neighbors. If neighbor is a ghost cell, + // its future FE index will be updated on the owning process and + // communicated at the next loop iteration. + if (neighbor->is_locally_owned()) + { + const level_type neighbor_level = static_cast( + future_levels[neighbor->global_active_cell_index()]); + + // ignore neighbors that are not part of the hierarchy + if (neighbor_level == invalid_level) + return false; + + if ((cell_level - max_difference) > neighbor_level) + { + future_levels[neighbor->global_active_cell_index()] = + cell_level - max_difference; + + return true; + } + } + + return false; + }; + bool levels_changed = false; bool levels_changed_in_cycle; do @@ -2767,29 +2796,24 @@ DoFHandler::prepare_coarsening_and_refinement( for (unsigned int f = 0; f < cell->n_faces(); ++f) if (cell->face(f)->at_boundary() == false) { - const auto neighbor = cell->neighbor(f); - - // We only care about locally owned neighbors. If neighbor is - // a ghost cell, its future FE index will be updated on the - // owning process and communicated at the next loop iteration. - if (neighbor->is_locally_owned()) + if (cell->face(f)->has_children()) { - const level_type neighbor_level = static_cast( - future_levels[neighbor->global_active_cell_index()]); - - // ignore neighbors that are not part of the hierarchy - if (neighbor_level == invalid_level) - continue; - - if ((cell_level - max_difference) > neighbor_level) + for (unsigned int sf = 0; + sf < cell->face(f)->n_children(); + ++sf) { - // update future level - future_levels[neighbor->global_active_cell_index()] = - cell_level - max_difference; - - levels_changed_in_cycle = true; + const auto neighbor = + cell->neighbor_child_on_subface(f, sf); + levels_changed_in_cycle |= + update_neighbor_level(neighbor, cell_level); } } + else + { + const auto neighbor = cell->neighbor(f); + levels_changed_in_cycle |= + update_neighbor_level(neighbor, cell_level); + } } } diff --git a/tests/hp/prepare_coarsening_and_refinement_03.cc b/tests/hp/prepare_coarsening_and_refinement_03.cc new file mode 100644 index 0000000000..806e656307 --- /dev/null +++ b/tests/hp/prepare_coarsening_and_refinement_03.cc @@ -0,0 +1,124 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + +// verify restrictions on level differences imposed by +// DoFHandler::prepare_coarsening_and_refinement() +// on h-refined grids + + +#include + +#include + +#include + +#include + +#include + +#include "../tests.h" + +#include "../test_grids.h" + + +unsigned int +divide_and_ceil(unsigned int x, unsigned int y) +{ + return x / y + (x % y > 0); +} + + +template +void +test(const unsigned int fes_size, const unsigned int max_difference) +{ + Assert(max_difference > 0, ExcInternalError()); + Assert(fes_size > 1, ExcInternalError()); + + // setup FE collection + hp::FECollection fes; + while (fes.size() < fes_size) + fes.push_back(FE_Q(1)); + + const unsigned int contains_fe_index = 0; + const auto sequence = fes.get_hierarchy_sequence(contains_fe_index); + + // set up line grid + // - refine leftmost column of cells consecutively + // - assign highest p-level to rightmost cell + // + // +++-+---+------+ + // +++ | | | + // +++-+---+ | + // +++ | | | + // +++-+---+------+ + // + // after prepare_coarsening_and_refinement(), each p-level will correspond to + // a unique h-level + + Triangulation tria; + TestGrids::hyper_line(tria, 2); + const unsigned int n_refinements = + divide_and_ceil(sequence.size() - 1, max_difference); + for (unsigned int i = 0; i < n_refinements; ++i) + { + for (const auto &cell : tria.active_cell_iterators()) + for (unsigned int v = 0; v < cell->n_vertices(); ++v) + if (cell->vertex(v)[0] == 0.) + cell->set_refine_flag(); + + tria.execute_coarsening_and_refinement(); + } + + DoFHandler dofh(tria); + for (const auto &cell : dofh.cell_iterators_on_level(0)) + if (cell->is_active()) + cell->set_active_fe_index(sequence.back()); + dofh.distribute_dofs(fes); + + const bool fe_indices_changed = + dofh.prepare_coarsening_and_refinement(max_difference, contains_fe_index); + tria.execute_coarsening_and_refinement(); + + (void)fe_indices_changed; + Assert(fe_indices_changed, ExcInternalError()); + +#ifdef DEBUG + // check each cell's active FE by its h-level + for (unsigned int l = 0; l < tria.n_levels(); ++l) + for (const auto &cell : dofh.cell_iterators_on_level(l)) + if (cell->is_active()) + { + const unsigned int expected_level = std::max( + 0, static_cast(sequence.size() - 1 - l * max_difference)); + Assert(cell->active_fe_index() == sequence[expected_level], + ExcInternalError()); + } +#endif + + deallog << "OK" << std::endl; +} + + +int +main() +{ + initlog(); + + test<2>(5, 1); + test<2>(10, 2); + test<2>(15, 3); +} diff --git a/tests/hp/prepare_coarsening_and_refinement_03.output b/tests/hp/prepare_coarsening_and_refinement_03.output new file mode 100644 index 0000000000..fb71de2867 --- /dev/null +++ b/tests/hp/prepare_coarsening_and_refinement_03.output @@ -0,0 +1,4 @@ + +DEAL::OK +DEAL::OK +DEAL::OK diff --git a/tests/mpi/prepare_coarsening_and_refinement_03.cc b/tests/mpi/prepare_coarsening_and_refinement_03.cc new file mode 100644 index 0000000000..1c84fa5f89 --- /dev/null +++ b/tests/mpi/prepare_coarsening_and_refinement_03.cc @@ -0,0 +1,166 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + +// verify restrictions on level differences imposed by +// DoFHandler::prepare_coarsening_and_refinement() +// on h-refined grids + + +#include + +#include +#include +#include + +#include + +#include + +#include + +#include "../tests.h" + +#include "../test_grids.h" + + +unsigned int +divide_and_ceil(unsigned int x, unsigned int y) +{ + return x / y + (x % y > 0); +} + + +template +void +test(parallel::TriangulationBase &tria, + const unsigned int fes_size, + const unsigned int max_difference) +{ + Assert(tria.n_levels() == 0, ExcInternalError()); + Assert(max_difference > 0, ExcInternalError()); + Assert(fes_size > 1, ExcInternalError()); + + // setup FE collection + hp::FECollection fes; + while (fes.size() < fes_size) + fes.push_back(FE_Q(1)); + + const unsigned int contains_fe_index = 0; + const auto sequence = fes.get_hierarchy_sequence(contains_fe_index); + + // set up line grid + // - refine leftmost column of cells consecutively + // - assign highest p-level to rightmost cell + // + // +++-+---+------+ + // +++ | | | + // +++-+---+ | + // +++ | | | + // +++-+---+------+ + // + // after prepare_coarsening_and_refinement(), each p-level will correspond to + // a unique column of cells and thus h-level + + TestGrids::hyper_line(tria, 2); + const unsigned int n_refinements = + divide_and_ceil(sequence.size() - 1, max_difference); + for (unsigned int i = 0; i < n_refinements; ++i) + { + for (const auto &cell : tria.active_cell_iterators()) + if (cell->is_locally_owned()) + for (unsigned int v = 0; v < cell->n_vertices(); ++v) + if (cell->vertex(v)[0] == 0.) + cell->set_refine_flag(); + + tria.execute_coarsening_and_refinement(); + } + + DoFHandler dofh(tria); + for (const auto &cell : dofh.cell_iterators_on_level(0)) + if (cell->is_active() && cell->is_locally_owned()) + cell->set_active_fe_index(sequence.back()); + dofh.distribute_dofs(fes); + + const bool fe_indices_changed = + dofh.prepare_coarsening_and_refinement(max_difference, contains_fe_index); + tria.execute_coarsening_and_refinement(); + + (void)fe_indices_changed; + Assert(fe_indices_changed, ExcInternalError()); + +#ifdef DEBUG + // check each cell's active FE by its h-level + for (unsigned int l = 0; l < tria.n_global_levels(); ++l) + for (const auto &cell : dofh.cell_iterators_on_level(l)) + if (cell->is_active() && cell->is_locally_owned()) + { + const unsigned int expected_level = std::max( + 0, static_cast(sequence.size() - 1 - l * max_difference)); + Assert(cell->active_fe_index() == sequence[expected_level], + ExcInternalError()); + } +#endif + + deallog << "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(); + + constexpr const unsigned int dim = 2; + + deallog << "parallel::shared::Triangulation" << std::endl; + { + parallel::shared::Triangulation tria(MPI_COMM_WORLD); + + test(tria, 5, 1); + tria.clear(); + test(tria, 10, 2); + tria.clear(); + test(tria, 15, 3); + } + + deallog << "parallel::shared::Triangulation with artificial cells" + << std::endl; + { + parallel::shared::Triangulation tria(MPI_COMM_WORLD, + Triangulation::none, + /*allow_artificial_cells=*/true); + + test(tria, 5, 1); + tria.clear(); + test(tria, 10, 2); + tria.clear(); + test(tria, 15, 3); + } + + deallog << "parallel::distributed::Triangulation" << std::endl; + { + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + + test(tria, 5, 1); + tria.clear(); + test(tria, 10, 2); + tria.clear(); + test(tria, 15, 3); + } +} diff --git a/tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=1.output b/tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=1.output new file mode 100644 index 0000000000..7715c474e9 --- /dev/null +++ b/tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=1.output @@ -0,0 +1,13 @@ + +DEAL::parallel::shared::Triangulation +DEAL::OK +DEAL::OK +DEAL::OK +DEAL::parallel::shared::Triangulation with artificial cells +DEAL::OK +DEAL::OK +DEAL::OK +DEAL::parallel::distributed::Triangulation +DEAL::OK +DEAL::OK +DEAL::OK diff --git a/tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=4.output b/tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=4.output new file mode 100644 index 0000000000..7715c474e9 --- /dev/null +++ b/tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=4.output @@ -0,0 +1,13 @@ + +DEAL::parallel::shared::Triangulation +DEAL::OK +DEAL::OK +DEAL::OK +DEAL::parallel::shared::Triangulation with artificial cells +DEAL::OK +DEAL::OK +DEAL::OK +DEAL::parallel::distributed::Triangulation +DEAL::OK +DEAL::OK +DEAL::OK -- 2.39.5