From bc934b3278477835dcd86414eeac481ec4781b5f Mon Sep 17 00:00:00 2001
From: Timo Heister <timo.heister@gmail.com>
Date: Sat, 12 Aug 2017 18:20:22 -0600
Subject: [PATCH] clean up p4est::iterate

---
 include/deal.II/distributed/p4est_wrappers.h | 21 +++---
 source/distributed/p4est_wrappers.cc         | 69 +++++++++++++-------
 source/distributed/p4est_wrappers.inst.in    | 10 +--
 source/distributed/tria.cc                   | 15 +----
 4 files changed, 63 insertions(+), 52 deletions(-)

diff --git a/include/deal.II/distributed/p4est_wrappers.h b/include/deal.II/distributed/p4est_wrappers.h
index 42d9e0f27e..a6d687e2ef 100644
--- a/include/deal.II/distributed/p4est_wrappers.h
+++ b/include/deal.II/distributed/p4est_wrappers.h
@@ -472,16 +472,14 @@ namespace internal
       static
       size_t (&connectivity_memory_used) (types<3>::connectivity *p4est);
 
-      template <int spacedim>
-      static void iterate(dealii::internal::p4est::types<3>::forest *parallel_forest,
-                          dealii::internal::p4est::types<3>::ghost *parallel_ghost,
-                          void *user_data);
+
 
       static const unsigned int max_level = P8EST_MAXLEVEL;
     };
 
 
 
+
     /**
      * This struct templatizes the p4est iterate structs and function
      * prototypes, which are used to execute callback functions for faces,
@@ -565,18 +563,15 @@ namespace internal
                          const typename types<dim>::topidx coarse_grid_cell);
 
 
+
     /**
-     * This is the callback data structure used to fill
-     * vertices_with_ghost_neighbors via the p4est_iterate tool
+     * Compute the ghost neighbors surrounding each vertex by querying p4est
      */
     template <int dim, int spacedim>
-    struct FindGhosts
-    {
-      const typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation;
-      sc_array_t                                                                *subids;
-      std::map<unsigned int, std::set<dealii::types::subdomain_id> >            *vertices_with_ghost_neighbors;
-    };
-
+    std::map<unsigned int, std::set<dealii::types::subdomain_id> >
+    compute_vertices_with_ghost_neighbors(const dealii::parallel::distributed::Triangulation<dim,spacedim> &tria,
+                                          typename dealii::internal::p4est::types<dim>::forest *parallel_forest,
+                                          typename dealii::internal::p4est::types<dim>::ghost *parallel_ghost);
 
   }
 }
diff --git a/source/distributed/p4est_wrappers.cc b/source/distributed/p4est_wrappers.cc
index 41cb8266aa..155d0e0729 100644
--- a/source/distributed/p4est_wrappers.cc
+++ b/source/distributed/p4est_wrappers.cc
@@ -35,6 +35,18 @@ namespace internal
         return out_cell;
       }
 
+      /**
+       * This is the callback data structure used to fill
+       * vertices_with_ghost_neighbors via the p4est_iterate tool
+       */
+      template <int dim, int spacedim>
+      struct FindGhosts
+      {
+        const typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation;
+        sc_array_t                                                                *subids;
+        std::map<unsigned int, std::set<dealii::types::subdomain_id> >            *vertices_with_ghost_neighbors;
+      };
+
 
       /** At a corner (vertex), determine if any of the neighboring cells are
        * ghosts.  If there are, find out their subdomain ids, and if this is a
@@ -436,22 +448,48 @@ namespace internal
     size_t (&functions<2>::connectivity_memory_used) (types<2>::connectivity *p4est)
       = p4est_connectivity_memory_used;
 
-    template <int spacedim>
-    void functions<2>::iterate(dealii::internal::p4est::types<2>::forest *parallel_forest,
-                               dealii::internal::p4est::types<2>::ghost *parallel_ghost,
-                               void *user_data)
+    template <int dim, int spacedim>
+    std::map<unsigned int, std::set<dealii::types::subdomain_id> >
+    compute_vertices_with_ghost_neighbors(const typename dealii::parallel::distributed::Triangulation<dim,spacedim> &tria,
+                                          typename dealii::internal::p4est::types<dim>::forest *parallel_forest,
+                                          typename dealii::internal::p4est::types<dim>::ghost *parallel_ghost)
     {
-      p4est_iterate (parallel_forest, parallel_ghost, user_data,
-                     nullptr,
-                     find_ghosts_face<2,spacedim>,
-                     find_ghosts_corner<2,spacedim>);
+      std::map<unsigned int, std::set<dealii::types::subdomain_id> > vertices_with_ghost_neighbors;
+
+      dealii::internal::p4est::FindGhosts<dim,spacedim> fg;
+      fg.subids = sc_array_new (sizeof (dealii::types::subdomain_id));
+      fg.triangulation = &tria;
+      fg.vertices_with_ghost_neighbors = &vertices_with_ghost_neighbors;
+
+      switch (dim)
+        {
+        case 2:
+          p4est_iterate (reinterpret_cast<dealii::internal::p4est::types<2>::forest *>(parallel_forest),
+                         reinterpret_cast<dealii::internal::p4est::types<2>::ghost *>(parallel_ghost),
+                         static_cast<void *>(&fg),
+                         NULL, find_ghosts_face<2,spacedim>, find_ghosts_corner<2,spacedim>);
+          break;
+
+        case 3:
+          p8est_iterate (reinterpret_cast<dealii::internal::p4est::types<3>::forest *>(parallel_forest),
+                         reinterpret_cast<dealii::internal::p4est::types<3>::ghost *>(parallel_ghost),
+                         static_cast<void *>(&fg),
+                         NULL, find_ghosts_face<3,3>, find_ghosts_edge<3,3>, find_ghosts_corner<3,3>);
+          break;
+
+        default:
+          Assert (false, ExcNotImplemented());
+        }
+
+      sc_array_destroy (fg.subids);
+
+      return vertices_with_ghost_neighbors;
     }
 
     const unsigned int functions<2>::max_level;
 
 
 
-
     int (&functions<3>::quadrant_compare) (const void *v1, const void *v2)
       = p8est_quadrant_compare;
 
@@ -638,23 +676,10 @@ namespace internal
     size_t (&functions<3>::connectivity_memory_used) (types<3>::connectivity *p4est)
       = p8est_connectivity_memory_used;
 
-    template <int spacedim>
-    void functions<3>::iterate(dealii::internal::p4est::types<3>::forest *parallel_forest,
-                               dealii::internal::p4est::types<3>::ghost *parallel_ghost,
-                               void *user_data)
-    {
-      p8est_iterate (parallel_forest, parallel_ghost, user_data,
-                     nullptr,
-                     find_ghosts_face<3,spacedim>,
-                     find_ghosts_edge<3,spacedim>,
-                     find_ghosts_corner<3,spacedim>);
-    }
-
     const unsigned int functions<3>::max_level;
 
 
 
-
     template <int dim>
     void
     init_quadrant_children
diff --git a/source/distributed/p4est_wrappers.inst.in b/source/distributed/p4est_wrappers.inst.in
index 464352a06c..69df43a82b 100644
--- a/source/distributed/p4est_wrappers.inst.in
+++ b/source/distributed/p4est_wrappers.inst.in
@@ -71,11 +71,11 @@ for (deal_II_dimension, deal_II_space_dimension : DIMENSIONS)
 #if     deal_II_dimension > 1
 #if       deal_II_space_dimension >= deal_II_dimension
     template
-    void
-    functions<deal_II_dimension>::
-    iterate<deal_II_space_dimension>(dealii::internal::p4est::types<deal_II_dimension>::forest *parallel_forest,
-                                     dealii::internal::p4est::types<deal_II_dimension>::ghost *parallel_ghost,
-                                     void *user_data);
+    std::map<unsigned int, std::set<dealii::types::subdomain_id> >
+    compute_vertices_with_ghost_neighbors(const dealii::parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &tria,
+                                          typename dealii::internal::p4est::types<deal_II_dimension>::forest *parallel_forest,
+                                          typename dealii::internal::p4est::types<deal_II_dimension>::ghost *parallel_ghost);
+
 #endif
 #endif
     \}
diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc
index 0b46eaf08e..3333a29ac6 100644
--- a/source/distributed/tria.cc
+++ b/source/distributed/tria.cc
@@ -3387,20 +3387,11 @@ namespace parallel
     {
       Assert (dim>1, ExcNotImplemented());
 
-      std::map<unsigned int, std::set<dealii::types::subdomain_id> > vertices_with_ghost_neighbors;
+      return dealii::internal::p4est::compute_vertices_with_ghost_neighbors<dim, spacedim> (*this,
+             this->parallel_forest,
+             this->parallel_ghost);
 
-      dealii::internal::p4est::FindGhosts<dim,spacedim> fg;
-      fg.subids = sc_array_new (sizeof (dealii::types::subdomain_id));
-      fg.triangulation = this;
-      fg.vertices_with_ghost_neighbors = &vertices_with_ghost_neighbors;
 
-      dealii::internal::p4est::functions<dim>::template iterate<spacedim> (this->parallel_forest,
-          this->parallel_ghost,
-          static_cast<void *>(&fg));
-
-      sc_array_destroy (fg.subids);
-
-      return vertices_with_ghost_neighbors;
     }
 
 
-- 
2.39.5