]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Clean up some p4est-related code. 3631/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 24 Nov 2016 16:55:41 +0000 (09:55 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 24 Nov 2016 18:41:52 +0000 (11:41 -0700)
Specifically, commit daf3146 (via #3625) introduced a dimension dependent dispatch
mechanism that allowed writing some code in a more generic way. The idea was right,
but the approach duplicated the dispatching because we already have a dispatching
mechanism via the internal::p4est::functions classes.

This patch simply merges the two approaches. It also allows to move some
code out of distributed/tria.cc into distributed/p4est_wrappers.cc, at the
cost of some code churn. The patch does not actually change any kind of
functionality -- it just moves things.

include/deal.II/distributed/p4est_wrappers.h
source/distributed/p4est_wrappers.cc
source/distributed/p4est_wrappers.inst.in
source/distributed/tria.cc

index 33deee0b0afec2157f490cc04657a95209972ded..795a5619479d1159594270c6c26b060d39633d99 100644 (file)
 #  include <p8est_communication.h>
 #  include <p8est_iterate.h>
 
+#include <map>
+#include <set>
+
+
 DEAL_II_NAMESPACE_OPEN
 
+namespace parallel
+{
+  namespace distributed
+  {
+    template <int dim, int spacedim> class Triangulation;
+  }
+}
+
+
 namespace internal
 {
   namespace p4est
@@ -272,6 +285,11 @@ namespace internal
       static
       size_t (&connectivity_memory_used) (types<2>::connectivity *p4est);
 
+      template <int spacedim>
+      static void iterate(dealii::internal::p4est::types<2>::forest *parallel_forest,
+                          dealii::internal::p4est::types<2>::ghost *parallel_ghost,
+                          void *user_data);
+
       static const unsigned int max_level = P4EST_MAXLEVEL;
     };
 
@@ -454,6 +472,11 @@ 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;
     };
 
@@ -539,6 +562,21 @@ namespace internal
     bool
     tree_exists_locally (const typename types<dim>::forest *parallel_forest,
                          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
+     */
+    template <int dim, int spacedim>
+    struct FindGhosts
+    {
+      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;
+    };
+
+
   }
 }
 
index d2e6d4c614cb0971d56ccab43f4b0c37e8884db2..c88da22d20eacd6a67d6bb5593433632e812b3fd 100644 (file)
@@ -1,4 +1,5 @@
 #include <deal.II/distributed/p4est_wrappers.h>
+#include <deal.II/distributed/tria.h>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -8,6 +9,251 @@ namespace internal
 {
   namespace p4est
   {
+    namespace
+    {
+      template <int dim, int spacedim>
+      typename dealii::Triangulation<dim,spacedim>::cell_iterator
+      cell_from_quad
+      (typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation,
+       typename dealii::internal::p4est::types<dim>::topidx treeidx,
+       typename dealii::internal::p4est::types<dim>::quadrant &quad)
+      {
+        int i, l = quad.level;
+        int child_id;
+        dealii::types::global_dof_index dealii_index =
+          triangulation->get_p4est_tree_to_coarse_cell_permutation()[treeidx];
+
+        for (i = 0; i < l; i++)
+          {
+            typename dealii::Triangulation<dim,spacedim>::cell_iterator cell (triangulation, i, dealii_index);
+            child_id = dealii::internal::p4est::functions<dim>::quadrant_ancestor_id (&quad, i + 1);
+            Assert (cell->has_children (), ExcMessage ("p4est quadrant does not correspond to a cell!"));
+            dealii_index = cell->child_index(child_id);
+          }
+
+        typename dealii::Triangulation<dim,spacedim>::cell_iterator out_cell (triangulation, l, dealii_index);
+
+        return out_cell;
+      }
+
+
+      /** 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
+       * local vertex, then add these subdomain ids to the map
+       * vertices_with_ghost_neighbors of that index
+       */
+      template <int dim, int spacedim>
+      void
+      find_ghosts_corner
+      (typename dealii::internal::p4est::iter<dim>::corner_info *info,
+       void *user_data)
+      {
+        int i, j;
+        int nsides = info->sides.elem_count;
+        typename dealii::internal::p4est::iter<dim>::corner_side *sides =
+          (typename dealii::internal::p4est::iter<dim>::corner_side *)
+          (info->sides.array);
+        FindGhosts<dim,spacedim> *fg = static_cast<FindGhosts<dim,spacedim> *>(user_data);
+        sc_array_t *subids = fg->subids;
+        typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation = fg->triangulation;
+        int nsubs;
+        dealii::types::subdomain_id *subdomain_ids;
+        std::map<unsigned int, std::set<dealii::types::subdomain_id> >
+        *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
+
+        subids->elem_count = 0;
+        for (i = 0; i < nsides; i++)
+          {
+            if (sides[i].is_ghost)
+              {
+                typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].quad));
+                Assert (cell->is_ghost(), ExcMessage ("ghost quad did not find ghost cell"));
+                dealii::types::subdomain_id *subid =
+                  static_cast<dealii::types::subdomain_id *>(sc_array_push (subids));
+                *subid = cell->subdomain_id();
+              }
+          }
+
+        if (!subids->elem_count)
+          {
+            return;
+          }
+
+        nsubs = (int) subids->elem_count;
+        subdomain_ids = (dealii::types::subdomain_id *) (subids->array);
+
+        for (i = 0; i < nsides; i++)
+          {
+            if (!sides[i].is_ghost)
+              {
+                typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].quad));
+
+                Assert (!cell->is_ghost(), ExcMessage ("local quad found ghost cell"));
+
+                for (j = 0; j < nsubs; j++)
+                  {
+                    (*vertices_with_ghost_neighbors)[cell->vertex_index(sides[i].corner)]
+                    .insert (subdomain_ids[j]);
+                  }
+              }
+          }
+
+        subids->elem_count = 0;
+      }
+
+      /** Similar to find_ghosts_corner, but for the hanging vertex in the
+       * middle of an edge
+       */
+      template <int dim, int spacedim>
+      void
+      find_ghosts_edge
+      (typename dealii::internal::p4est::iter<dim>::edge_info *info,
+       void *user_data)
+      {
+        int i, j, k;
+        int nsides = info->sides.elem_count;
+        typename dealii::internal::p4est::iter<dim>::edge_side *sides =
+          (typename dealii::internal::p4est::iter<dim>::edge_side *)
+          (info->sides.array);
+        FindGhosts<dim,spacedim> *fg = static_cast<FindGhosts<dim,spacedim> *>(user_data);
+        sc_array_t *subids = fg->subids;
+        typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation = fg->triangulation;
+        int nsubs;
+        dealii::types::subdomain_id *subdomain_ids;
+        std::map<unsigned int, std::set<dealii::types::subdomain_id> >
+        *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
+
+        subids->elem_count = 0;
+        for (i = 0; i < nsides; i++)
+          {
+            if (sides[i].is_hanging)
+              {
+                for (j = 0; j < 2; j++)
+                  {
+                    if (sides[i].is.hanging.is_ghost[j])
+                      {
+                        typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
+                        dealii::types::subdomain_id *subid =
+                          static_cast<dealii::types::subdomain_id *>(sc_array_push (subids));
+                        *subid = cell->subdomain_id();
+                      }
+                  }
+              }
+          }
+
+        if (!subids->elem_count)
+          {
+            return;
+          }
+
+        nsubs = (int) subids->elem_count;
+        subdomain_ids = (dealii::types::subdomain_id *) (subids->array);
+
+        for (i = 0; i < nsides; i++)
+          {
+            if (sides[i].is_hanging)
+              {
+                for (j = 0; j < 2; j++)
+                  {
+                    if (!sides[i].is.hanging.is_ghost[j])
+                      {
+                        typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
+
+                        for (k = 0; k < nsubs; k++)
+                          {
+                            (*vertices_with_ghost_neighbors)[cell->vertex_index(p8est_edge_corners[sides[i].edge][1^j])]
+                            .insert (subdomain_ids[k]);
+                          }
+                      }
+                  }
+              }
+          }
+
+        subids->elem_count = 0;
+      }
+
+      /** Similar to find_ghosts_corner, but for the hanging vertex in the
+       * middle of a face
+       */
+      template <int dim, int spacedim>
+      void
+      find_ghosts_face
+      (typename dealii::internal::p4est::iter<dim>::face_info *info,
+       void *user_data)
+      {
+        int i, j, k;
+        int nsides = info->sides.elem_count;
+        typename dealii::internal::p4est::iter<dim>::face_side *sides =
+          (typename dealii::internal::p4est::iter<dim>::face_side *)
+          (info->sides.array);
+        FindGhosts<dim,spacedim> *fg = static_cast<FindGhosts<dim,spacedim> *>(user_data);
+        sc_array_t *subids = fg->subids;
+        typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation = fg->triangulation;
+        int nsubs;
+        dealii::types::subdomain_id *subdomain_ids;
+        std::map<unsigned int, std::set<dealii::types::subdomain_id> >
+        *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
+        int limit = (dim == 2) ? 2 : 4;
+
+        subids->elem_count = 0;
+        for (i = 0; i < nsides; i++)
+          {
+            if (sides[i].is_hanging)
+              {
+                for (j = 0; j < limit; j++)
+                  {
+                    if (sides[i].is.hanging.is_ghost[j])
+                      {
+                        typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
+                        dealii::types::subdomain_id *subid =
+                          static_cast<dealii::types::subdomain_id *>(sc_array_push (subids));
+                        *subid = cell->subdomain_id();
+                      }
+                  }
+              }
+          }
+
+        if (!subids->elem_count)
+          {
+            return;
+          }
+
+        nsubs = (int) subids->elem_count;
+        subdomain_ids = (dealii::types::subdomain_id *) (subids->array);
+
+        for (i = 0; i < nsides; i++)
+          {
+            if (sides[i].is_hanging)
+              {
+                for (j = 0; j < limit; j++)
+                  {
+                    if (!sides[i].is.hanging.is_ghost[j])
+                      {
+                        typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
+
+                        for (k = 0; k < nsubs; k++)
+                          {
+                            if (dim == 2)
+                              {
+                                (*vertices_with_ghost_neighbors)[cell->vertex_index(p4est_face_corners[sides[i].face][(limit - 1)^j])]
+                                .insert (subdomain_ids[k]);
+                              }
+                            else
+                              {
+                                (*vertices_with_ghost_neighbors)[cell->vertex_index(p8est_face_corners[sides[i].face][(limit - 1)^j])]
+                                .insert (subdomain_ids[k]);
+                              }
+                          }
+                      }
+                  }
+              }
+          }
+
+        subids->elem_count = 0;
+      }
+    }
+
+
     int (&functions<2>::quadrant_compare) (const void *v1, const void *v2)
       = p4est_quadrant_compare;
 
@@ -191,8 +437,22 @@ 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)
+    {
+      p4est_iterate (parallel_forest, parallel_ghost, user_data,
+                     NULL,
+                     find_ghosts_face<2,spacedim>,
+                     find_ghosts_corner<2,spacedim>);
+    }
+
     const unsigned int functions<2>::max_level;
 
+
+
+
     int (&functions<3>::quadrant_compare) (const void *v1, const void *v2)
       = p8est_quadrant_compare;
 
@@ -379,8 +639,23 @@ 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,
+                     NULL,
+                     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
index 1806dc820163183f3917c2a97f0da8679f4d88f3..464352a06cc2d8d66665b0daa88067ec1a0d8200 100644 (file)
@@ -58,3 +58,28 @@ for (deal_II_dimension : DIMENSIONS)
     \}
 #   endif // DEAL_II_WITH_P4EST
 }
+
+
+for (deal_II_dimension, deal_II_space_dimension : DIMENSIONS)
+{
+#   ifdef DEAL_II_WITH_P4EST
+
+    namespace internal
+    \{
+    namespace p4est
+    \{
+#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);
+#endif
+#endif
+    \}
+    \}
+#   endif // DEAL_II_WITH_P4EST
+}
+
index 646375d75cbfaa8f3b68efda3df0566e5488ffae..3e175ffdb432b7c04290fe29c032cca85bd1fa6f 100644 (file)
@@ -3380,280 +3380,6 @@ namespace parallel
     }
 
 
-
-    namespace
-    {
-      /**
-       * 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
-      {
-        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
-       * local vertex, then add these subdomain ids to the map
-       * vertices_with_ghost_neighbors of that index
-       */
-      template <int dim, int spacedim>
-      void
-      find_ghosts_corner
-      (typename dealii::internal::p4est::iter<dim>::corner_info *info,
-       void *user_data)
-      {
-        int i, j;
-        int nsides = info->sides.elem_count;
-        typename dealii::internal::p4est::iter<dim>::corner_side *sides =
-          (typename dealii::internal::p4est::iter<dim>::corner_side *)
-          (info->sides.array);
-        FindGhosts<dim,spacedim> *fg = static_cast<FindGhosts<dim,spacedim> *>(user_data);
-        sc_array_t *subids = fg->subids;
-        typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation = fg->triangulation;
-        int nsubs;
-        dealii::types::subdomain_id *subdomain_ids;
-        std::map<unsigned int, std::set<dealii::types::subdomain_id> >
-        *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
-
-        subids->elem_count = 0;
-        for (i = 0; i < nsides; i++)
-          {
-            if (sides[i].is_ghost)
-              {
-                typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].quad));
-                Assert (cell->is_ghost(), ExcMessage ("ghost quad did not find ghost cell"));
-                dealii::types::subdomain_id *subid =
-                  static_cast<dealii::types::subdomain_id *>(sc_array_push (subids));
-                *subid = cell->subdomain_id();
-              }
-          }
-
-        if (!subids->elem_count)
-          {
-            return;
-          }
-
-        nsubs = (int) subids->elem_count;
-        subdomain_ids = (dealii::types::subdomain_id *) (subids->array);
-
-        for (i = 0; i < nsides; i++)
-          {
-            if (!sides[i].is_ghost)
-              {
-                typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].quad));
-
-                Assert (!cell->is_ghost(), ExcMessage ("local quad found ghost cell"));
-
-                for (j = 0; j < nsubs; j++)
-                  {
-                    (*vertices_with_ghost_neighbors)[cell->vertex_index(sides[i].corner)]
-                    .insert (subdomain_ids[j]);
-                  }
-              }
-          }
-
-        subids->elem_count = 0;
-      }
-
-      /** Similar to find_ghosts_corner, but for the hanging vertex in the
-       * middle of an edge
-       */
-      template <int dim, int spacedim>
-      void
-      find_ghosts_edge
-      (typename dealii::internal::p4est::iter<dim>::edge_info *info,
-       void *user_data)
-      {
-        int i, j, k;
-        int nsides = info->sides.elem_count;
-        typename dealii::internal::p4est::iter<dim>::edge_side *sides =
-          (typename dealii::internal::p4est::iter<dim>::edge_side *)
-          (info->sides.array);
-        FindGhosts<dim,spacedim> *fg = static_cast<FindGhosts<dim,spacedim> *>(user_data);
-        sc_array_t *subids = fg->subids;
-        typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation = fg->triangulation;
-        int nsubs;
-        dealii::types::subdomain_id *subdomain_ids;
-        std::map<unsigned int, std::set<dealii::types::subdomain_id> >
-        *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
-
-        subids->elem_count = 0;
-        for (i = 0; i < nsides; i++)
-          {
-            if (sides[i].is_hanging)
-              {
-                for (j = 0; j < 2; j++)
-                  {
-                    if (sides[i].is.hanging.is_ghost[j])
-                      {
-                        typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
-                        dealii::types::subdomain_id *subid =
-                          static_cast<dealii::types::subdomain_id *>(sc_array_push (subids));
-                        *subid = cell->subdomain_id();
-                      }
-                  }
-              }
-          }
-
-        if (!subids->elem_count)
-          {
-            return;
-          }
-
-        nsubs = (int) subids->elem_count;
-        subdomain_ids = (dealii::types::subdomain_id *) (subids->array);
-
-        for (i = 0; i < nsides; i++)
-          {
-            if (sides[i].is_hanging)
-              {
-                for (j = 0; j < 2; j++)
-                  {
-                    if (!sides[i].is.hanging.is_ghost[j])
-                      {
-                        typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
-
-                        for (k = 0; k < nsubs; k++)
-                          {
-                            (*vertices_with_ghost_neighbors)[cell->vertex_index(p8est_edge_corners[sides[i].edge][1^j])]
-                            .insert (subdomain_ids[k]);
-                          }
-                      }
-                  }
-              }
-          }
-
-        subids->elem_count = 0;
-      }
-
-      /** Similar to find_ghosts_corner, but for the hanging vertex in the
-       * middle of a face
-       */
-      template <int dim, int spacedim>
-      void
-      find_ghosts_face
-      (typename dealii::internal::p4est::iter<dim>::face_info *info,
-       void *user_data)
-      {
-        int i, j, k;
-        int nsides = info->sides.elem_count;
-        typename dealii::internal::p4est::iter<dim>::face_side *sides =
-          (typename dealii::internal::p4est::iter<dim>::face_side *)
-          (info->sides.array);
-        FindGhosts<dim,spacedim> *fg = static_cast<FindGhosts<dim,spacedim> *>(user_data);
-        sc_array_t *subids = fg->subids;
-        typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation = fg->triangulation;
-        int nsubs;
-        dealii::types::subdomain_id *subdomain_ids;
-        std::map<unsigned int, std::set<dealii::types::subdomain_id> >
-        *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
-        int limit = (dim == 2) ? 2 : 4;
-
-        subids->elem_count = 0;
-        for (i = 0; i < nsides; i++)
-          {
-            if (sides[i].is_hanging)
-              {
-                for (j = 0; j < limit; j++)
-                  {
-                    if (sides[i].is.hanging.is_ghost[j])
-                      {
-                        typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
-                        dealii::types::subdomain_id *subid =
-                          static_cast<dealii::types::subdomain_id *>(sc_array_push (subids));
-                        *subid = cell->subdomain_id();
-                      }
-                  }
-              }
-          }
-
-        if (!subids->elem_count)
-          {
-            return;
-          }
-
-        nsubs = (int) subids->elem_count;
-        subdomain_ids = (dealii::types::subdomain_id *) (subids->array);
-
-        for (i = 0; i < nsides; i++)
-          {
-            if (sides[i].is_hanging)
-              {
-                for (j = 0; j < limit; j++)
-                  {
-                    if (!sides[i].is.hanging.is_ghost[j])
-                      {
-                        typename dealii::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
-
-                        for (k = 0; k < nsubs; k++)
-                          {
-                            if (dim == 2)
-                              {
-                                (*vertices_with_ghost_neighbors)[cell->vertex_index(p4est_face_corners[sides[i].face][(limit - 1)^j])]
-                                .insert (subdomain_ids[k]);
-                              }
-                            else
-                              {
-                                (*vertices_with_ghost_neighbors)[cell->vertex_index(p8est_face_corners[sides[i].face][(limit - 1)^j])]
-                                .insert (subdomain_ids[k]);
-                              }
-                          }
-                      }
-                  }
-              }
-          }
-
-        subids->elem_count = 0;
-      }
-
-
-      /**
-       * This struct wraps the p4est_iterate / p8est_iterate in a templated
-       * fashion so that we can avoid the compiler complaining about invalid
-       * types being referred, which would be the case inside switch or if
-       * statements.
-       */
-      template <int dim, int spacedim>
-      struct pXest_iterate_wrapper
-      {
-        static void iterate(typename dealii::internal::p4est::types<dim>::forest *parallel_forest,
-                            typename dealii::internal::p4est::types<dim>::ghost *parallel_ghost,
-                            void *user_data);
-      };
-
-      template <int spacedim>
-      struct pXest_iterate_wrapper<2,spacedim>
-      {
-        static void iterate(dealii::internal::p4est::types<2>::forest *parallel_forest,
-                            dealii::internal::p4est::types<2>::ghost *parallel_ghost,
-                            void *user_data)
-        {
-          p4est_iterate (parallel_forest, parallel_ghost, user_data,
-                         NULL, find_ghosts_face<2,spacedim>,
-                         find_ghosts_corner<2,spacedim>);
-        }
-      };
-
-      template <int spacedim>
-      struct pXest_iterate_wrapper<3,spacedim>
-      {
-        static void 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,
-                         NULL, find_ghosts_face<3,spacedim>, find_ghosts_edge<3,spacedim>,
-                         find_ghosts_corner<3,spacedim>);
-        }
-      };
-    }
-
-
     /**
      * Determine the neighboring subdomains that are adjacent to each vertex.
      * This is achieved via the p4est_iterate/p8est_iterate tool
@@ -3667,14 +3393,14 @@ namespace parallel
     {
       Assert (dim>1, ExcNotImplemented());
 
-      FindGhosts<dim,spacedim> fg;
+      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;
 
-      pXest_iterate_wrapper<dim,spacedim>::iterate (this->parallel_forest,
-                                                    this->parallel_ghost,
-                                                    static_cast<void *>(&fg));
+      dealii::internal::p4est::functions<dim>::template iterate<spacedim> (this->parallel_forest,
+          this->parallel_ghost,
+          static_cast<void *>(&fg));
 
       sc_array_destroy (fg.subids);
     }
@@ -4138,31 +3864,6 @@ namespace parallel
       return weights;
     }
 
-    template <int dim, int spacedim>
-    typename dealii::Triangulation<dim,spacedim>::cell_iterator
-    cell_from_quad
-    (typename dealii::parallel::distributed::Triangulation<dim,spacedim> *triangulation,
-     typename dealii::internal::p4est::types<dim>::topidx treeidx,
-     typename dealii::internal::p4est::types<dim>::quadrant &quad)
-    {
-      int i, l = quad.level;
-      int child_id;
-      types::global_dof_index dealii_index =
-        triangulation->get_p4est_tree_to_coarse_cell_permutation()[treeidx];
-
-      for (i = 0; i < l; i++)
-        {
-          typename dealii::Triangulation<dim,spacedim>::cell_iterator cell (triangulation, i, dealii_index);
-          child_id = dealii::internal::p4est::functions<dim>::quadrant_ancestor_id (&quad, i + 1);
-          Assert (cell->has_children (), ExcMessage ("p4est quadrant does not correspond to a cell!"));
-          dealii_index = cell->child_index(child_id);
-        }
-
-      typename dealii::Triangulation<dim,spacedim>::cell_iterator out_cell (triangulation, l, dealii_index);
-
-      return out_cell;
-    }
-
 
 
     template <int spacedim>

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.