]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Reduce the volume of communication somewhat. 14138/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 15 Jul 2022 03:24:51 +0000 (21:24 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 15 Jul 2022 18:46:40 +0000 (12:46 -0600)
source/distributed/shared_tria.cc

index f9b999d28f73860016d3a0f0da2ecb09506105bf..f111aa38fe7c9eeda11b6f9debf7738c27f5f007 100644 (file)
@@ -14,6 +14,7 @@
 // ---------------------------------------------------------------------
 
 #include <deal.II/base/mpi.h>
+#include <deal.II/base/mpi.templates.h>
 #include <deal.II/base/utilities.h>
 
 #include <deal.II/distributed/shared_tria.h>
@@ -27,6 +28,8 @@
 
 #include <deal.II/lac/sparsity_tools.h>
 
+#include <type_traits>
+
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -354,15 +357,30 @@ namespace parallel
       // make sure that all refinement/coarsening flags are the same on all
       // processes
       {
-        std::vector<unsigned int> refinement_configurations(
-          this->n_active_cells() * 2, 0u);
+        // Obtain the type used to store the different possibilities
+        // a cell can be refined. This is a bit awkward because
+        // what `cell->refine_flag_set()` returns is a struct
+        // type, RefinementCase, which internally stores a
+        // std::uint8_t, which actually holds integers of
+        // enum type RefinementPossibilities<dim>::Possibilities.
+        // In the following, use the actual name of the enum, but
+        // make sure that it is in fact a `std::uint8_t` or
+        // equally sized type.
+        using int_type = std::underlying_type_t<
+          typename RefinementPossibilities<dim>::Possibilities>;
+        static_assert(sizeof(int_type) == sizeof(std::uint8_t),
+                      "Internal type mismatch.");
+
+        std::vector<int_type> refinement_configurations(this->n_active_cells() *
+                                                          2,
+                                                        int_type(0));
         for (const auto &cell : this->active_cell_iterators())
           if (cell->is_locally_owned())
             {
               refinement_configurations[cell->active_cell_index() * 2 + 0] =
-                cell->refine_flag_set();
+                static_cast<int_type>(cell->refine_flag_set());
               refinement_configurations[cell->active_cell_index() * 2 + 1] =
-                cell->coarsen_flag_set();
+                static_cast<int_type>(cell->coarsen_flag_set() ? 1 : 0);
             }
 
         Utilities::MPI::max(refinement_configurations,
@@ -385,7 +403,7 @@ namespace parallel
               ExcMessage(
                 "Refinement/coarsening flags of cells are not consistent in parallel!"));
 
-            if (refinement_configurations[cell->active_cell_index() * 2 + 0] >
+            if (refinement_configurations[cell->active_cell_index() * 2 + 0] !=
                 0)
               cell->set_refine_flag(RefinementCase<dim>(
                 refinement_configurations[cell->active_cell_index() * 2 + 0]));

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.