]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added `p:d:TemporarilyMatchRefineFlags`. 11981/head
authorMarc Fehling <mafehling.git@gmail.com>
Wed, 24 Mar 2021 22:22:12 +0000 (16:22 -0600)
committerMarc Fehling <mafehling.git@gmail.com>
Thu, 1 Apr 2021 18:54:01 +0000 (12:54 -0600)
include/deal.II/distributed/tria.h
source/distributed/tria.cc
source/distributed/tria.inst.in
tests/mpi/limit_p_level_difference_01.cc
tests/mpi/limit_p_level_difference_02.cc
tests/mpi/limit_p_level_difference_03.cc
tests/mpi/limit_p_level_difference_04.cc

index 6deb092a97e2cab02f3b890709123be2c305c008..3727241811a13eea53a54a6c060ab24c159a2951 100644 (file)
@@ -73,6 +73,18 @@ namespace GridTools
   template <typename CellIterator>
   struct PeriodicFacePair;
 }
+
+namespace internal
+{
+  namespace parallel
+  {
+    namespace distributed
+    {
+      template <int, int>
+      class TemporarilyMatchRefineFlags;
+    }
+  } // namespace parallel
+} // namespace internal
 #  endif
 
 namespace parallel
@@ -803,6 +815,10 @@ namespace parallel
 
       template <int, int, class>
       friend class dealii::FETools::internal::ExtrapolateImplementation;
+
+      template <int, int>
+      friend class dealii::internal::parallel::distributed::
+        TemporarilyMatchRefineFlags;
     };
 
 
@@ -910,6 +926,10 @@ namespace parallel
       virtual types::coarse_cell_id
       coarse_cell_index_to_coarse_cell_id(
         const unsigned int coarse_cell_index) const override;
+
+      template <int, int>
+      friend class dealii::internal::parallel::distributed::
+        TemporarilyMatchRefineFlags;
     };
   } // namespace distributed
 } // namespace parallel
@@ -934,7 +954,7 @@ namespace parallel
      */
     template <int dim, int spacedim = dim>
     class Triangulation
-      : public dealii::parallel::TriangulationBase<dim, spacedim>
+      : public dealii::parallel::DistributedTriangulationBase<dim, spacedim>
     {
     public:
       /**
@@ -950,6 +970,75 @@ namespace parallel
 #endif
 
 
+namespace internal
+{
+  namespace parallel
+  {
+    namespace distributed
+    {
+      /**
+       * This class temporarily modifies the refine and coarsen flags of all
+       * active cells to match the p4est oracle.
+       *
+       * The modification only happens on parallel::distributed::Triangulation
+       * objects, and persists for the lifetime of an instantiation of this
+       * class.
+       *
+       * The TemporarilyMatchRefineFlags class should only be used in
+       * combination with the Triangulation::Signals::post_p4est_refinement
+       * signal. At this stage, the p4est orcale already has been refined, but
+       * the triangulation is still unchanged. After the modification, all
+       * refine and coarsen flags describe how the traingulation will acutally
+       * be refined.
+       */
+      template <int dim, int spacedim = dim>
+      class TemporarilyMatchRefineFlags : public Subscriptor
+      {
+      public:
+        /**
+         * Constructor.
+         *
+         * Stores the refine and coarsen flags of all active cells if the
+         * provided Triangulation is of type
+         * parallel::distributed::Triangulation.
+         *
+         * Adjusts them to be consistent with the p4est oracle.
+         */
+        TemporarilyMatchRefineFlags(Triangulation<dim, spacedim> &tria);
+
+        /**
+         * Destructor.
+         *
+         * Returns the refine and coarsen flags of all active cells on the
+         * parallel::distributed::Triangulation into their previous state.
+         */
+        ~TemporarilyMatchRefineFlags();
+
+      private:
+        /**
+         * The modified parallel::distributed::Triangulation.
+         */
+        const SmartPointer<
+          dealii::parallel::distributed::Triangulation<dim, spacedim>>
+          distributed_tria;
+
+        /**
+         * A vector that temporarily stores the refine flags before they have
+         * been modified on the parallel::distributed::Triangulation.
+         */
+        std::vector<bool> saved_refine_flags;
+
+        /**
+         * A vector that temporarily stores the coarsen flags before they have
+         * been modified on the parallel::distributed::Triangulation.
+         */
+        std::vector<bool> saved_coarsen_flags;
+      };
+    } // namespace distributed
+  }   // namespace parallel
+} // namespace internal
+
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 65b34777b0181d1c5a938738fb63899875560ac9..09c21bf910be1aa9f84a5fb54a79c6c37c08d12d 100644 (file)
@@ -3559,6 +3559,92 @@ namespace parallel
 
 
 
+namespace internal
+{
+  namespace parallel
+  {
+    namespace distributed
+    {
+      template <int dim, int spacedim>
+      TemporarilyMatchRefineFlags<dim, spacedim>::TemporarilyMatchRefineFlags(
+        Triangulation<dim, spacedim> &tria)
+        : distributed_tria(
+            dynamic_cast<
+              dealii::parallel::distributed::Triangulation<dim, spacedim> *>(
+              &tria))
+      {
+#ifdef DEAL_II_WITH_P4EST
+        if (distributed_tria != nullptr)
+          {
+            // Save the current set of refinement flags, and adjust the
+            // refinement flags to be consistent with the p4est oracle.
+            distributed_tria->save_coarsen_flags(saved_coarsen_flags);
+            distributed_tria->save_refine_flags(saved_refine_flags);
+
+            for (const auto &pair : distributed_tria->local_cell_relations)
+              {
+                const auto &cell   = pair.first;
+                const auto &status = pair.second;
+
+                switch (status)
+                  {
+                    case dealii::Triangulation<dim, spacedim>::CELL_PERSIST:
+                      // cell remains unchanged
+                      cell->clear_refine_flag();
+                      cell->clear_coarsen_flag();
+                      break;
+
+                    case dealii::Triangulation<dim, spacedim>::CELL_REFINE:
+                      // cell will be refined
+                      cell->clear_coarsen_flag();
+                      cell->set_refine_flag();
+                      break;
+
+                    case dealii::Triangulation<dim, spacedim>::CELL_COARSEN:
+                      // children of this cell will be coarsened
+                      for (const auto &child : cell->child_iterators())
+                        {
+                          child->clear_refine_flag();
+                          child->set_coarsen_flag();
+                        }
+                      break;
+
+                    case dealii::Triangulation<dim, spacedim>::CELL_INVALID:
+                      // do nothing as cell does not exist yet
+                      break;
+
+                    default:
+                      Assert(false, ExcInternalError());
+                      break;
+                  }
+              }
+          }
+#endif
+      }
+
+
+
+      template <int dim, int spacedim>
+      TemporarilyMatchRefineFlags<dim, spacedim>::~TemporarilyMatchRefineFlags()
+      {
+#ifdef DEAL_II_WITH_P4EST
+        if (distributed_tria)
+          {
+            // Undo the refinement flags modification.
+            distributed_tria->load_coarsen_flags(saved_coarsen_flags);
+            distributed_tria->load_refine_flags(saved_refine_flags);
+          }
+#else
+        // pretend that this destructor does something to silence clang-tidy
+        (void)distributed_tria;
+#endif
+      }
+    } // namespace distributed
+  }   // namespace parallel
+} // namespace internal
+
+
+
 /*-------------- Explicit Instantiations -------------------------------*/
 #include "tria.inst"
 
index 141e63fc17c97165eca86ab9eb1b582f4743e164..c0e43213581830c77a083dd936182798992b13c6 100644 (file)
 
 
 
-for (deal_II_dimension : DIMENSIONS)
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
   {
     namespace parallel
     \{
       namespace distributed
       \{
-        template class Triangulation<deal_II_dimension>;
-#if deal_II_dimension < 3
-        template class Triangulation<deal_II_dimension, deal_II_dimension + 1>;
+#if deal_II_dimension <= deal_II_space_dimension
+        template class Triangulation<deal_II_dimension,
+                                     deal_II_space_dimension>;
 #endif
-#if deal_II_dimension < 2
-        template class Triangulation<deal_II_dimension, deal_II_dimension + 2>;
+      \}
+    \}
+
+    namespace internal
+    \{
+      namespace parallel
+      \{
+        namespace distributed
+        \{
+#if deal_II_dimension <= deal_II_space_dimension
+          template class TemporarilyMatchRefineFlags<deal_II_dimension,
+                                                     deal_II_space_dimension>;
 #endif
+        \}
       \}
     \}
   }
index 9ca79b1e5f7992b01bcd8c06cb36db7d71fe7166..a21913a9bd3b67079499d7b06e1e5a6bb30a6cd6 100644 (file)
@@ -75,13 +75,20 @@ test(const unsigned int fes_size, const unsigned int max_difference)
     if (cell->is_locally_owned() && cell->center() == Point<dim>())
       cell->set_active_fe_index(sequence.back());
 
-  const bool fe_indices_changed =
-    hp::Refinement::limit_p_level_difference(dofh,
-                                             max_difference,
-                                             contains_fe_index);
+  bool fe_indices_changed = false;
+  tria.signals.post_p4est_refinement.connect(
+    [&]() {
+      const internal::parallel::distributed::TemporarilyMatchRefineFlags<dim>
+        refine_modifier(tria);
+      fe_indices_changed =
+        hp::Refinement::limit_p_level_difference(dofh,
+                                                 max_difference,
+                                                 contains_fe_index);
+    },
+    boost::signals2::at_front);
+
   tria.execute_coarsening_and_refinement();
 
-  (void)fe_indices_changed;
   Assert(fe_indices_changed, ExcInternalError());
 
   // display number of cells for each FE index
index f87bad819050eef666f99ede8b57a8a91dc4a36f..5c7cf6c1af91b5c9331195843f5711a761a9e336 100644 (file)
@@ -70,6 +70,18 @@ test(const unsigned int fes_size, const unsigned int max_difference)
   DoFHandler<dim> dofh(tria);
   dofh.distribute_dofs(fes);
 
+  bool fe_indices_changed;
+  tria.signals.post_p4est_refinement.connect(
+    [&]() {
+      const internal::parallel::distributed::TemporarilyMatchRefineFlags<dim>
+        refine_modifier(tria);
+      fe_indices_changed =
+        hp::Refinement::limit_p_level_difference(dofh,
+                                                 max_difference,
+                                                 contains_fe_index);
+    },
+    boost::signals2::at_front);
+
   // increase p-level in center subsequently
   for (unsigned int i = 0; i < sequence.size() - 1; ++i)
     {
@@ -79,13 +91,9 @@ test(const unsigned int fes_size, const unsigned int max_difference)
           cell->set_future_fe_index(
             fes.next_in_hierarchy(cell->active_fe_index()));
 
-      const bool fe_indices_changed =
-        hp::Refinement::limit_p_level_difference(dofh,
-                                                 max_difference,
-                                                 contains_fe_index);
+      fe_indices_changed = false;
       tria.execute_coarsening_and_refinement();
 
-      (void)fe_indices_changed;
       if (i >= max_difference)
         Assert(fe_indices_changed, ExcInternalError());
 
index 127e339b2a564c379f58ebdb55bea100b54d8949..8ac8ed9308b6070f2286ac648bf415f90353e5bb 100644 (file)
@@ -90,13 +90,20 @@ test(const unsigned int fes_size, const unsigned int max_difference)
       cell->set_active_fe_index(sequence.back());
   dofh.distribute_dofs(fes);
 
-  const bool fe_indices_changed =
-    hp::Refinement::limit_p_level_difference(dofh,
-                                             max_difference,
-                                             contains_fe_index);
+  bool fe_indices_changed = false;
+  tria.signals.post_p4est_refinement.connect(
+    [&]() {
+      const internal::parallel::distributed::TemporarilyMatchRefineFlags<dim>
+        refine_modifier(tria);
+      fe_indices_changed =
+        hp::Refinement::limit_p_level_difference(dofh,
+                                                 max_difference,
+                                                 contains_fe_index);
+    },
+    boost::signals2::at_front);
+
   tria.execute_coarsening_and_refinement();
 
-  (void)fe_indices_changed;
   Assert(fe_indices_changed, ExcInternalError());
 
 #ifdef DEBUG
index f3e2d1c5f488ce14c45cbe7f7cda52dde450864a..07546f3b6234f7ebeb1bd0b79f5b8df44095019f 100644 (file)
@@ -82,15 +82,22 @@ test(const unsigned int max_difference)
       cell->set_active_fe_index(sequence.back());
   dofh.distribute_dofs(fes);
 
-  const bool fe_indices_changed =
-    hp::Refinement::limit_p_level_difference(dofh,
-                                             max_difference,
-                                             contains_fe_index);
-  (void)fe_indices_changed;
-  Assert(fe_indices_changed, ExcInternalError());
+  bool fe_indices_changed = false;
+  tria.signals.post_p4est_refinement.connect(
+    [&]() {
+      const internal::parallel::distributed::TemporarilyMatchRefineFlags<dim>
+        refine_modifier(tria);
+      fe_indices_changed =
+        hp::Refinement::limit_p_level_difference(dofh,
+                                                 max_difference,
+                                                 contains_fe_index);
+    },
+    boost::signals2::at_front);
 
   tria.execute_coarsening_and_refinement();
 
+  Assert(fe_indices_changed, ExcInternalError());
+
   deallog << "active FE indices after adaptation:" << std::endl;
   for (const auto &cell : dofh.active_cell_iterators())
     if (cell->is_locally_owned())

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.