]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Bugfix: hp::Refinement::choose_p_over_h() in parallel 10309/head
authorMarc Fehling <marc.fehling@gmx.net>
Thu, 21 May 2020 20:56:56 +0000 (22:56 +0200)
committerMarc Fehling <marc.fehling@gmx.net>
Sat, 23 May 2020 22:59:24 +0000 (00:59 +0200)
doc/news/changes/minor/20200522Fehling [new file with mode: 0644]
source/hp/refinement.cc
tests/sharedtria/hp_choose_p_over_h.cc [new file with mode: 0644]
tests/sharedtria/hp_choose_p_over_h.with_mpi=true.mpirun=2.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20200522Fehling b/doc/news/changes/minor/20200522Fehling
new file mode 100644 (file)
index 0000000..b04788d
--- /dev/null
@@ -0,0 +1,3 @@
+Bugfix: hp::Refinement::choose_p_over_h() now works in parallel.
+<br>
+(Marc Fehling, 2020/05/22)
index a6191283d3dd3988027a1860fc6f36e610dde6dc..1b22490f9b472b90f91edf501a16e385f6a1dc9b 100644 (file)
@@ -23,6 +23,7 @@
 #include <deal.II/distributed/tria_base.h>
 
 #include <deal.II/grid/grid_refinement.h>
+#include <deal.II/grid/grid_tools.h>
 
 #include <deal.II/hp/dof_handler.h>
 #include <deal.II/hp/refinement.h>
@@ -635,6 +636,33 @@ namespace hp
     void
     choose_p_over_h(const hp::DoFHandler<dim, spacedim> &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<CellId, std::pair<bool, bool>> ghost_buffer;
+
+      if (dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
+            &dof_handler.get_triangulation()))
+        {
+          auto pack = [](
+                        const typename dealii::hp::DoFHandler<dim, spacedim>::
+                          active_cell_iterator &cell) -> std::pair<bool, bool> {
+            return {cell->coarsen_flag_set(), cell->future_fe_index_set()};
+          };
+
+          auto unpack = [&ghost_buffer](
+                          const typename dealii::hp::DoFHandler<dim, spacedim>::
+                            active_cell_iterator &    cell,
+                          const std::pair<bool, bool> pair) -> void {
+            ghost_buffer.emplace(cell->id(), pair);
+          };
+
+          GridTools::exchange_cell_data_to_ghosts<
+            std::pair<bool, bool>,
+            dealii::hp::DoFHandler<dim, spacedim>>(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<bool, bool> &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 (file)
index 0000000..0482b8a
--- /dev/null
@@ -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 <deal.II/distributed/shared_tria.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/hp/refinement.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test()
+{
+  // setup
+  const unsigned int n_cells = 2;
+
+  parallel::shared::Triangulation<dim> tr(
+    MPI_COMM_WORLD,
+    ::Triangulation<dim>::none,
+    false,
+    parallel::shared::Triangulation<dim>::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<unsigned int> rep(dim, 1);
+  rep[0] = n_cells;
+  Point<dim> 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<dim> fes;
+  for (unsigned int d = 1; d <= 2; ++d)
+    fes.push_back(FE_Q<dim>(d));
+
+  hp::DoFHandler<dim> 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 (file)
index 0000000..47a88e9
--- /dev/null
@@ -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

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.