]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Clean up compute_vertices_with_ghost_neighbors
authorPeter Munch <peterrmuench@gmail.com>
Wed, 18 Sep 2019 05:47:25 +0000 (07:47 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Wed, 18 Sep 2019 05:47:25 +0000 (07:47 +0200)
include/deal.II/distributed/fully_distributed_tria.h
include/deal.II/distributed/tria.h
source/distributed/fully_distributed_tria.cc
source/distributed/tria.cc
source/distributed/tria_base.cc

index eef0b55d2d5bc36c838552b15039734f750823c5..4c7843278c70778d0a6ed9f4a9a5653e4361b351 100644 (file)
@@ -327,13 +327,6 @@ namespace parallel
       virtual std::size_t
       memory_consumption() const override;
 
-      /**
-       * This function determines the neighboring subdomains that are adjacent
-       * to each vertex.
-       */
-      virtual std::map<unsigned int, std::set<dealii::types::subdomain_id>>
-      compute_vertices_with_ghost_neighbors() const override;
-
       virtual bool
       is_multilevel_hierarchy_constructed() const override;
 
index 467a21489cb798675ee13d3cd60458c9385bfd43..8376b466f60d37416236d9e2a184e524c8400748 100644 (file)
@@ -1193,16 +1193,6 @@ namespace parallel
       std::vector<unsigned int>
       get_cell_weights() const;
 
-      /**
-       * Override the implementation in parallel::TriangulationBase because
-       * we can ask p4est about ghost neighbors across periodic boundaries.
-       *
-       * Specifically, this function determines the neighboring subdomains that
-       * are adjacent to each vertex.
-       */
-      virtual std::map<unsigned int, std::set<dealii::types::subdomain_id>>
-      compute_vertices_with_ghost_neighbors() const override;
-
       /**
        * This method returns a bit vector of length tria.n_vertices()
        * indicating the locally active vertices on a level, i.e., the vertices
@@ -1368,13 +1358,6 @@ namespace parallel
        */
       Settings settings;
 
-      /**
-       * Like above, this method, which is only implemented for dim = 2 or 3,
-       * needs a stub because it is used in dof_handler_policy.cc
-       */
-      virtual std::map<unsigned int, std::set<dealii::types::subdomain_id>>
-      compute_vertices_with_ghost_neighbors() const override;
-
       /**
        * Like above, this method, which is only implemented for dim = 2 or 3,
        * needs a stub because it is used in dof_handler_policy.cc
index f4fe3720d48c23dce00e524e8206116b1ed0122e..b6625bc5419fb84e19429cb9d0fea20cf0580dde 100644 (file)
@@ -346,130 +346,6 @@ namespace parallel
 
 
 
-    template <int dim, int spacedim>
-    std::map<unsigned int, std::set<dealii::types::subdomain_id>>
-    Triangulation<dim, spacedim>::compute_vertices_with_ghost_neighbors() const
-    {
-      // 1) collect for each vertex on periodic faces all vertices it coincides
-      //    with
-      std::map<unsigned int, std::set<unsigned int>>
-        vertex_to_coinciding_vertices;
-      {
-        // 1a) collect nodes coinciding due to periodicity
-        std::map<unsigned int, unsigned int> vertex_to_coinciding_vertex;
-        for (auto &cell : this->active_cell_iterators())
-          if (cell->is_locally_owned() || cell->is_ghost())
-            for (auto i = 0u; i < GeometryInfo<dim>::faces_per_cell; ++i)
-              if (cell->has_periodic_neighbor(i) &&
-                  cell->periodic_neighbor(i)->active())
-                {
-                  auto face_t = cell->face(i);
-                  auto face_n = cell->periodic_neighbor(i)->face(
-                    cell->periodic_neighbor_face_no(i));
-                  for (auto j = 0u; j < GeometryInfo<dim>::vertices_per_face;
-                       ++j)
-                    {
-                      auto         v_t  = face_t->vertex_index(j);
-                      auto         v_n  = face_n->vertex_index(j);
-                      unsigned int temp = std::min(v_t, v_n);
-                      {
-                        auto it = vertex_to_coinciding_vertex.find(v_t);
-                        if (it != vertex_to_coinciding_vertex.end())
-                          temp = std::min(temp, it->second);
-                      }
-                      {
-                        auto it = vertex_to_coinciding_vertex.find(v_n);
-                        if (it != vertex_to_coinciding_vertex.end())
-                          temp = std::min(temp, it->second);
-                      }
-                      vertex_to_coinciding_vertex[v_t] = temp;
-                      vertex_to_coinciding_vertex[v_n] = temp;
-                    }
-                }
-
-        // 1b) compress map: let vertices point to the coinciding vertex with
-        //     the smallest index
-        for (auto &p : vertex_to_coinciding_vertex)
-          {
-            if (p.first == p.second)
-              continue;
-            unsigned int temp = p.second;
-            while (temp != vertex_to_coinciding_vertex[temp])
-              temp = vertex_to_coinciding_vertex[temp];
-            p.second = temp;
-          }
-
-#ifdef DEBUG
-        // check if map is actually compressed
-        for (auto p : vertex_to_coinciding_vertex)
-          {
-            if (p.first == p.second)
-              continue;
-            auto pp = vertex_to_coinciding_vertex.find(p.second);
-            if (pp->first == pp->second)
-              continue;
-            AssertThrow(false, ExcMessage("Map has to be compressed!"));
-          }
-#endif
-
-        // 1c) create a map: smallest index of coinciding index -> all
-        // coinciding indices
-        std::map<unsigned int, std::set<unsigned int>>
-          smallest_coinciding_vertex_to_coinciding_vertices;
-        for (auto p : vertex_to_coinciding_vertex)
-          smallest_coinciding_vertex_to_coinciding_vertices[p.second] =
-            std::set<unsigned int>();
-
-        for (auto p : vertex_to_coinciding_vertex)
-          smallest_coinciding_vertex_to_coinciding_vertices[p.second].insert(
-            p.first);
-
-        // 1d) create a map: vertex -> all coinciding indices
-        for (auto &s : smallest_coinciding_vertex_to_coinciding_vertices)
-          for (auto &ss : s.second)
-            vertex_to_coinciding_vertices[ss] = s.second;
-      }
-
-      // 2) collect vertices belonging to local cells
-      std::vector<bool> vertex_of_own_cell(this->n_vertices(), false);
-      for (const auto &cell : this->active_cell_iterators())
-        if (cell->is_locally_owned())
-          for (auto v = 0u; v < GeometryInfo<dim>::vertices_per_cell; ++v)
-            vertex_of_own_cell[cell->vertex_index(v)] = true;
-
-      // 3) for for each vertex belonging to a locally owned cell all ghost
-      //    neighbors (including the periodic own)
-      std::map<unsigned int, std::set<dealii::types::subdomain_id>> result;
-
-      // loop over all active ghost cells
-      for (const auto &cell : this->active_cell_iterators())
-        if (cell->is_ghost())
-          {
-            const auto owner = cell->subdomain_id();
-
-            // loop over all its vertices
-            for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell;
-                 ++v)
-              {
-                // set owner if vertex belongs to a local cell
-                if (vertex_of_own_cell[cell->vertex_index(v)])
-                  result[cell->vertex_index(v)].insert(owner);
-
-                // mark also nodes coinciding due to periodicity
-                auto coinciding_vertices =
-                  vertex_to_coinciding_vertices.find(cell->vertex_index(v));
-                if (coinciding_vertices != vertex_to_coinciding_vertices.end())
-                  for (auto coinciding_vertex : coinciding_vertices->second)
-                    if (vertex_of_own_cell[coinciding_vertex])
-                      result[coinciding_vertex].insert(owner);
-              }
-          }
-
-      return result;
-    }
-
-
-
     template <int dim, int spacedim>
     bool
     Triangulation<dim, spacedim>::is_multilevel_hierarchy_constructed() const
index aef7015732c32b8efca999cae9ccaee645ae7252..528e8234df1d333a9f5ed3e687355d7ea0016be6 100644 (file)
@@ -4194,7 +4194,8 @@ namespace parallel
       // Here, it is sufficient to collect all vertices that are located
       // at that boundary.
       const std::map<unsigned int, std::set<dealii::types::subdomain_id>>
-        vertices_with_ghost_neighbors = compute_vertices_with_ghost_neighbors();
+        vertices_with_ghost_neighbors =
+          this->compute_vertices_with_ghost_neighbors();
 
       // now collect cells and their vertices
       // for the interested neighbors
@@ -4470,19 +4471,6 @@ namespace parallel
 
 
 
-    template <int dim, int spacedim>
-    std::map<unsigned int, std::set<dealii::types::subdomain_id>>
-    Triangulation<dim, spacedim>::compute_vertices_with_ghost_neighbors() const
-    {
-      Assert(dim > 1, ExcNotImplemented());
-
-      return dealii::internal::p4est::compute_vertices_with_ghost_neighbors<
-        dim,
-        spacedim>(*this, this->parallel_forest, this->parallel_ghost);
-    }
-
-
-
     template <int dim, int spacedim>
     std::vector<bool>
     Triangulation<dim, spacedim>::mark_locally_active_vertices_on_level(
@@ -5035,16 +5023,6 @@ namespace parallel
 
 
 
-    template <int spacedim>
-    std::map<unsigned int, std::set<dealii::types::subdomain_id>>
-    Triangulation<1, spacedim>::compute_vertices_with_ghost_neighbors() const
-    {
-      Assert(false, ExcNotImplemented());
-      return std::map<unsigned int, std::set<dealii::types::subdomain_id>>();
-    }
-
-
-
     template <int spacedim>
     std::map<unsigned int, std::set<dealii::types::subdomain_id>>
     Triangulation<1, spacedim>::compute_level_vertices_with_ghost_neighbors(
index 44252fa12f287f8aefa03f70cf322b3525c84e37..283029a696f6dfde51a6b4ac663c98126de99905 100644 (file)
@@ -362,36 +362,125 @@ namespace parallel
   TriangulationBase<dim, spacedim>::compute_vertices_with_ghost_neighbors()
     const
   {
-    // TODO: we are not treating periodic neighbors correctly here. If we do
-    // we can remove the overriding implementation for p::d::Triangulation
-    // that is currently using a p4est callback to get correct ghost neighbors
-    // over periodic faces.
-    Assert(this->get_periodic_face_map().size() == 0, ExcNotImplemented());
+    // 1) collect for each vertex on periodic faces all vertices it coincides
+    //    with
+    std::map<unsigned int, std::set<unsigned int>>
+      vertex_to_coinciding_vertices;
+    {
+      // 1a) collect nodes coinciding due to periodicity
+      std::map<unsigned int, unsigned int> vertex_to_coinciding_vertex;
+      for (auto &cell : this->active_cell_iterators())
+        if (cell->is_locally_owned() || cell->is_ghost())
+          for (auto i = 0u; i < GeometryInfo<dim>::faces_per_cell; ++i)
+            if (cell->has_periodic_neighbor(i) &&
+                cell->periodic_neighbor(i)->active())
+              {
+                auto face_t = cell->face(i);
+                auto face_n = cell->periodic_neighbor(i)->face(
+                  cell->periodic_neighbor_face_no(i));
+                for (auto j = 0u; j < GeometryInfo<dim>::vertices_per_face; ++j)
+                  {
+                    auto         v_t  = face_t->vertex_index(j);
+                    auto         v_n  = face_n->vertex_index(j);
+                    unsigned int temp = std::min(v_t, v_n);
+                    {
+                      auto it = vertex_to_coinciding_vertex.find(v_t);
+                      if (it != vertex_to_coinciding_vertex.end())
+                        temp = std::min(temp, it->second);
+                    }
+                    {
+                      auto it = vertex_to_coinciding_vertex.find(v_n);
+                      if (it != vertex_to_coinciding_vertex.end())
+                        temp = std::min(temp, it->second);
+                    }
+                    vertex_to_coinciding_vertex[v_t] = temp;
+                    vertex_to_coinciding_vertex[v_n] = temp;
+                  }
+              }
+
+      // 1b) compress map: let vertices point to the coinciding vertex with
+      //     the smallest index
+      for (auto &p : vertex_to_coinciding_vertex)
+        {
+          if (p.first == p.second)
+            continue;
+          unsigned int temp = p.second;
+          while (temp != vertex_to_coinciding_vertex[temp])
+            temp = vertex_to_coinciding_vertex[temp];
+          p.second = temp;
+        }
 
+#ifdef DEBUG
+      // check if map is actually compressed
+      for (auto p : vertex_to_coinciding_vertex)
+        {
+          if (p.first == p.second)
+            continue;
+          auto pp = vertex_to_coinciding_vertex.find(p.second);
+          if (pp->first == pp->second)
+            continue;
+          AssertThrow(false, ExcMessage("Map has to be compressed!"));
+        }
+#endif
 
-    std::vector<bool> vertex_of_own_cell(this->n_vertices(), false);
+      // 1c) create a map: smallest index of coinciding index -> all
+      // coinciding indices
+      std::map<unsigned int, std::set<unsigned int>>
+        smallest_coinciding_vertex_to_coinciding_vertices;
+      for (auto p : vertex_to_coinciding_vertex)
+        smallest_coinciding_vertex_to_coinciding_vertices[p.second] =
+          std::set<unsigned int>();
+
+      for (auto p : vertex_to_coinciding_vertex)
+        smallest_coinciding_vertex_to_coinciding_vertices[p.second].insert(
+          p.first);
+
+      // 1d) create a map: vertex -> all coinciding indices
+      for (auto &s : smallest_coinciding_vertex_to_coinciding_vertices)
+        for (auto &ss : s.second)
+          vertex_to_coinciding_vertices[ss] = s.second;
+    }
 
+    // 2) collect vertices belonging to local cells
+    std::vector<bool> vertex_of_own_cell(this->n_vertices(), false);
     for (const auto &cell : this->active_cell_iterators())
       if (cell->is_locally_owned())
-        for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell; ++v)
+        for (auto v = 0u; v < GeometryInfo<dim>::vertices_per_cell; ++v)
           vertex_of_own_cell[cell->vertex_index(v)] = true;
 
+    // 3) for for each vertex belonging to a locally owned cell all ghost
+    //    neighbors (including the periodic own)
     std::map<unsigned int, std::set<dealii::types::subdomain_id>> result;
+
+    // loop over all active ghost cells
     for (const auto &cell : this->active_cell_iterators())
       if (cell->is_ghost())
         {
-          const types::subdomain_id owner = cell->subdomain_id();
+          const auto owner = cell->subdomain_id();
+
+          // loop over all its vertices
           for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell;
                ++v)
             {
+              // set owner if vertex belongs to a local cell
               if (vertex_of_own_cell[cell->vertex_index(v)])
                 result[cell->vertex_index(v)].insert(owner);
+
+              // mark also nodes coinciding due to periodicity
+              auto coinciding_vertices =
+                vertex_to_coinciding_vertices.find(cell->vertex_index(v));
+              if (coinciding_vertices != vertex_to_coinciding_vertices.end())
+                for (auto coinciding_vertex : coinciding_vertices->second)
+                  if (vertex_of_own_cell[coinciding_vertex])
+                    result[coinciding_vertex].insert(owner);
             }
         }
 
     return result;
   }
 
+
+
   template <int dim, int spacedim>
   DistributedTriangulationBase<dim, spacedim>::DistributedTriangulationBase(
     MPI_Comm mpi_communicator,

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.