From addb8df58df4199dc054f675465d9495e3324908 Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Mon, 28 Sep 2020 15:14:52 -0400 Subject: [PATCH] reduce number of prepare_coarsening_and_refinement() calls I noticed that we call prepare_coarsening_and_refinement() many times when using a p::d::Triangulation. Many of these are certainly unnecessary: 1. distributed tria for mesh_reconstruction_after_repartitioning. This is used to remove all refinement and smoothing will be called inside the execute...() call below. No need to call our more expensive version. 2. serial execute_coarsening_and_refinement() should not call the derived version, which in turn calls the serial version. For the test tests/distributed_grids/refine_periodic_01: Before: 154 sequential, 154 parallel After: 126 sequential, 53 parallel fixup --- source/distributed/tria.cc | 16 ++++++++++------ source/grid/tria.cc | 7 ++++++- 2 files changed, 16 insertions(+), 7 deletions(-) diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc index 09c21bf910..e5a70ac970 100644 --- a/source/distributed/tria.cc +++ b/source/distributed/tria.cc @@ -2244,7 +2244,7 @@ namespace parallel void Triangulation::copy_local_forest_to_triangulation() { - // disable mesh smoothing for recreating the deal.II triangulation, + // Disable mesh smoothing for recreating the deal.II triangulation, // otherwise we might not be able to reproduce the p4est mesh // exactly. We restore the original smoothing at the end of this // function. Note that the smoothing flag is used in the normal @@ -2270,21 +2270,25 @@ namespace parallel bool mesh_changed = false; - // remove all deal.II refinements. Note that we could skip this and + // Remove all deal.II refinements. Note that we could skip this and // start from our current state, because the algorithm later coarsens as // necessary. This has the advantage of being faster when large parts // of the local partition changes (likely) and gives a deterministic // ordering of the cells (useful for snapshot/resume). // TODO: is there a more efficient way to do this? if (settings & mesh_reconstruction_after_repartitioning) - while (this->begin_active()->level() > 0) + while (this->n_levels() > 1) { - for (const auto &cell : this->active_cell_iterators()) + // Instead of marking all active cells, we slice off the finest + // level, one level at a time. This takes the same number of + // iterations but solves an issue where not all cells on a + // periodic boundary are indeed coarsened and we run into an + // irrelevant Assert() in update_periodic_face_map(). + for (const auto &cell : + this->active_cell_iterators_on_level(this->n_levels() - 1)) { cell->set_coarsen_flag(); } - - this->prepare_coarsening_and_refinement(); try { dealii::Triangulation:: diff --git a/source/grid/tria.cc b/source/grid/tria.cc index 10f106a297..33cfcd5252 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -13378,7 +13378,12 @@ template void Triangulation::execute_coarsening_and_refinement() { - prepare_coarsening_and_refinement(); + // Call our version of prepare_coarsening_and_refinement() even if a derived + // class like parallel::distributed::Triangulation overrides it. Their + // function will be called in their execute_coarsening_and_refinement() + // function. Even in a distributed computation our job here is to reconstruct + // the local part of the mesh and as such checking our flags is enough. + Triangulation::prepare_coarsening_and_refinement(); // verify a case with which we have had // some difficulty in the past (see the -- 2.39.5