]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Generalize DoFHandlerPolicy of p:d:t 8579/head
authorPeter Munch <peterrmuench@gmail.com>
Thu, 15 Aug 2019 04:00:03 +0000 (06:00 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Tue, 20 Aug 2019 06:03:30 +0000 (08:03 +0200)
doc/news/changes/minor/20190814PeterMunch-1 [new file with mode: 0644]
include/deal.II/distributed/tria.h
source/distributed/tria.cc
source/dofs/dof_handler_policy.cc
tests/mpi/mg_ghost_dofs_periodic_06.mpirun=7.output
tests/sharedtria/dof_05.cc

diff --git a/doc/news/changes/minor/20190814PeterMunch-1 b/doc/news/changes/minor/20190814PeterMunch-1
new file mode 100644 (file)
index 0000000..dd131aa
--- /dev/null
@@ -0,0 +1,4 @@
+New: Generalize internal::DoFHandlerImplementation::Policy::ParallelDistributed such that it uses 
+the new definition of CellId.
+<br>
+(Peter Munch, 2019/08/14)
index 4c4b9b0b64591e794fedf22897c449493c353391..cb70df8a589d08c57889c8f39fd6d0190509a0df 100644 (file)
@@ -1221,6 +1221,16 @@ namespace parallel
       std::vector<bool>
       mark_locally_active_vertices_on_level(const int level) const;
 
+      virtual unsigned int
+      coarse_cell_id_to_coarse_cell_index(
+        const types::coarse_cell_id coarse_cell_id) const override;
+
+      virtual types::coarse_cell_id
+      coarse_cell_index_to_coarse_cell_id(
+        const unsigned int coarse_cell_index) const override;
+
+
+
       template <typename>
       friend class dealii::internal::DoFHandlerImplementation::Policy::
         ParallelDistributed;
@@ -1388,6 +1398,14 @@ namespace parallel
        */
       virtual std::vector<bool>
       mark_locally_active_vertices_on_level(const unsigned int level) const;
+
+      virtual unsigned int
+      coarse_cell_id_to_coarse_cell_index(
+        const types::coarse_cell_id coarse_cell_id) const override;
+
+      virtual types::coarse_cell_id
+      coarse_cell_index_to_coarse_cell_id(
+        const unsigned int coarse_cell_index) const override;
     };
   } // namespace distributed
 } // namespace parallel
index 590ec84593b245e16032c502346c60f672586f82..a186866e3df7da7ccdc15c79bce319967b1e5888 100644 (file)
@@ -4551,6 +4551,26 @@ namespace parallel
 
 
 
+    template <int dim, int spacedim>
+    unsigned int
+    Triangulation<dim, spacedim>::coarse_cell_id_to_coarse_cell_index(
+      const types::coarse_cell_id coarse_cell_id) const
+    {
+      return p4est_tree_to_coarse_cell_permutation[coarse_cell_id];
+    }
+
+
+
+    template <int dim, int spacedim>
+    types::coarse_cell_id
+    Triangulation<dim, spacedim>::coarse_cell_index_to_coarse_cell_id(
+      const unsigned int coarse_cell_index) const
+    {
+      return coarse_cell_to_p4est_tree_permutation[coarse_cell_index];
+    }
+
+
+
     template <int dim, int spacedim>
     void
     Triangulation<dim, spacedim>::add_periodicity(
@@ -5047,6 +5067,27 @@ namespace parallel
 
 
 
+    template <int spacedim>
+    unsigned int
+    Triangulation<1, spacedim>::coarse_cell_id_to_coarse_cell_index(
+      const types::coarse_cell_id) const
+    {
+      Assert(false, ExcNotImplemented());
+      return 0;
+    }
+
+
+
+    template <int spacedim>
+    types::coarse_cell_id
+    Triangulation<1, spacedim>::coarse_cell_index_to_coarse_cell_id(
+      const unsigned int) const
+    {
+      Assert(false, ExcNotImplemented());
+      return 0;
+    }
+
+
     template <int spacedim>
     void
     Triangulation<1, spacedim>::load(const std::string &, const bool)
index f10dd33387604ae790b5ac9534fbb20719c72b46..4cb6d6228afe618155ab0cf29a4fa4428968720f 100644 (file)
@@ -4538,11 +4538,7 @@ namespace internal
         communicate_mg_ghost_cells(
           const typename parallel::distributed::Triangulation<dim, spacedim>
             &             tria,
-          DoFHandlerType &dof_handler,
-          const std::vector<dealii::types::global_dof_index>
-            &coarse_cell_to_p4est_tree_permutation,
-          const std::vector<dealii::types::global_dof_index>
-            &p4est_tree_to_coarse_cell_permutation)
+          DoFHandlerType &dof_handler)
         {
           // build list of cells to request for each neighbor
           std::set<dealii::types::subdomain_id> level_ghost_owners =
@@ -4565,7 +4561,7 @@ namespace internal
 
               find_marked_mg_ghost_cells_recursively<dim, spacedim>(
                 tria,
-                coarse_cell_to_p4est_tree_permutation[cell->index()],
+                cell->id().get_coarse_cell_id(),
                 cell,
                 p4est_coarse_cell,
                 neighbor_cell_list);
@@ -4634,11 +4630,14 @@ namespace internal
                    c < cell_data_transfer_buffer.tree_indices.size();
                    ++c)
                 {
+                  const auto temp =
+                    CellId(cell_data_transfer_buffer.tree_indices[c], 0, NULL)
+                      .to_cell(tria);
+
                   typename DoFHandlerType::level_cell_iterator cell(
                     &dof_handler.get_triangulation(),
                     0,
-                    p4est_tree_to_coarse_cell_permutation
-                      [cell_data_transfer_buffer.tree_indices[c]],
+                    temp->index(),
                     &dof_handler);
 
                   typename dealii::internal::p4est::types<dim>::quadrant
@@ -4703,12 +4702,12 @@ namespace internal
                    c < cell_data_transfer_buffer.tree_indices.size();
                    ++c, dofs += 1 + dofs[0])
                 {
+                  const auto temp =
+                    CellId(cell_data_transfer_buffer.tree_indices[c], 0, NULL)
+                      .to_cell(tria);
+
                   typename DoFHandlerType::level_cell_iterator cell(
-                    &tria,
-                    0,
-                    p4est_tree_to_coarse_cell_permutation
-                      [cell_data_transfer_buffer.tree_indices[c]],
-                    &dof_handler);
+                    &tria, 0, temp->index(), &dof_handler);
 
                   typename dealii::internal::p4est::types<dim>::quadrant
                     p4est_coarse_cell;
@@ -4750,9 +4749,7 @@ namespace internal
         void
         communicate_mg_ghost_cells(
           const typename parallel::distributed::Triangulation<1, spacedim> &,
-          DoFHandler<1, spacedim> &,
-          const std::vector<dealii::types::global_dof_index> &,
-          const std::vector<dealii::types::global_dof_index> &)
+          DoFHandler<1, spacedim> &)
         {
           Assert(false, ExcNotImplemented());
         }
@@ -4763,9 +4760,7 @@ namespace internal
         void
         communicate_mg_ghost_cells(
           const typename parallel::distributed::Triangulation<1, spacedim> &,
-          hp::DoFHandler<1, spacedim> &,
-          const std::vector<dealii::types::global_dof_index> &,
-          const std::vector<dealii::types::global_dof_index> &)
+          hp::DoFHandler<1, spacedim> &)
         {
           Assert(false, ExcNotImplemented());
         }
@@ -4794,9 +4789,7 @@ namespace internal
         void
         communicate_dof_indices_on_marked_cells(
           const DoFHandler<1, spacedim> &,
-          const std::map<unsigned int, std::set<dealii::types::subdomain_id>> &,
-          const std::vector<dealii::types::global_dof_index> &,
-          const std::vector<dealii::types::global_dof_index> &)
+          const std::map<unsigned int, std::set<dealii::types::subdomain_id>> &)
         {
           Assert(false, ExcNotImplemented());
         }
@@ -4807,9 +4800,7 @@ namespace internal
         void
         communicate_dof_indices_on_marked_cells(
           const hp::DoFHandler<1, spacedim> &,
-          const std::map<unsigned int, std::set<dealii::types::subdomain_id>> &,
-          const std::vector<dealii::types::global_dof_index> &,
-          const std::vector<dealii::types::global_dof_index> &)
+          const std::map<unsigned int, std::set<dealii::types::subdomain_id>> &)
         {
           Assert(false, ExcNotImplemented());
         }
@@ -4820,9 +4811,7 @@ namespace internal
         void
         communicate_dof_indices_on_marked_cells(
           const DoFHandlerType &dof_handler,
-          const std::map<unsigned int, std::set<dealii::types::subdomain_id>> &,
-          const std::vector<dealii::types::global_dof_index> &,
-          const std::vector<dealii::types::global_dof_index> &)
+          const std::map<unsigned int, std::set<dealii::types::subdomain_id>> &)
         {
 #  ifndef DEAL_II_WITH_MPI
           (void)vertices_with_ghost_neighbors;
@@ -5172,10 +5161,7 @@ namespace internal
           // as explained in the 'distributed' paper, this has to be
           // done twice
           communicate_dof_indices_on_marked_cells(
-            *dof_handler,
-            vertices_with_ghost_neighbors,
-            triangulation->coarse_cell_to_p4est_tree_permutation,
-            triangulation->p4est_tree_to_coarse_cell_permutation);
+            *dof_handler, vertices_with_ghost_neighbors);
 
           // in case of hp::DoFHandlers, we may have received valid
           // indices of degrees of freedom that are dominated by a fe
@@ -5190,10 +5176,7 @@ namespace internal
           //                    may still have invalid ones. thus, exchange
           //                    one more time.
           communicate_dof_indices_on_marked_cells(
-            *dof_handler,
-            vertices_with_ghost_neighbors,
-            triangulation->coarse_cell_to_p4est_tree_permutation,
-            triangulation->p4est_tree_to_coarse_cell_permutation);
+            *dof_handler, vertices_with_ghost_neighbors);
 
           // at this point, we must have taken care of the data transfer
           // on all cells we had previously marked. verify this
@@ -5400,11 +5383,7 @@ namespace internal
           // Phase 1. Request all marked cells from corresponding owners. If we
           // managed to get every DoF, remove the user_flag, otherwise we
           // will request them again in the step below.
-          communicate_mg_ghost_cells(
-            *triangulation,
-            *dof_handler,
-            triangulation->coarse_cell_to_p4est_tree_permutation,
-            triangulation->p4est_tree_to_coarse_cell_permutation);
+          communicate_mg_ghost_cells(*triangulation, *dof_handler);
 
           // have a barrier so that sends from above and below this
           // place are not mixed up.
@@ -5424,11 +5403,7 @@ namespace internal
 
           // Phase 2, only request the cells that were not completed
           // in Phase 1.
-          communicate_mg_ghost_cells(
-            *triangulation,
-            *dof_handler,
-            triangulation->coarse_cell_to_p4est_tree_permutation,
-            triangulation->p4est_tree_to_coarse_cell_permutation);
+          communicate_mg_ghost_cells(*triangulation, *dof_handler);
 
 #  ifdef DEBUG
           // make sure we have removed all flags:
@@ -5709,10 +5684,7 @@ namespace internal
               // as explained in the 'distributed' paper, this has to be
               // done twice
               communicate_dof_indices_on_marked_cells(
-                *dof_handler,
-                vertices_with_ghost_neighbors,
-                triangulation->coarse_cell_to_p4est_tree_permutation,
-                triangulation->p4est_tree_to_coarse_cell_permutation);
+                *dof_handler, vertices_with_ghost_neighbors);
 
               // in case of hp::DoFHandlers, we may have received valid
               // indices of degrees of freedom that are dominated by a fe
@@ -5723,10 +5695,7 @@ namespace internal
                 *dof_handler);
 
               communicate_dof_indices_on_marked_cells(
-                *dof_handler,
-                vertices_with_ghost_neighbors,
-                triangulation->coarse_cell_to_p4est_tree_permutation,
-                triangulation->p4est_tree_to_coarse_cell_permutation);
+                *dof_handler, vertices_with_ghost_neighbors);
 
               triangulation->load_user_flags(user_flags);
             }
index b3e28c9ef1b7a6078ecb9bfea30357ebd7727a1b..c04afa7b60e9f1cad6f5505122c09aa62defd9d9 100644 (file)
@@ -1,9 +1,9 @@
 
 DEAL:0::0_0: 0
 DEAL:0::1_0: 2
-DEAL:0::2_0: 3
+DEAL:0::4_0: 3
+DEAL:0::2_0: 2
 DEAL:0::3_0: 2
-DEAL:0::4_0: 2
 DEAL:0::5_0: 4
 DEAL:0::6_0: 4
 DEAL:0::7_0: 5
@@ -12,10 +12,10 @@ DEAL:0::0_1:0 0
 DEAL:0::0_1:1 1
 DEAL:0::0_1:2 1
 DEAL:0::0_1:3 2
-DEAL:0::2_1:0 4294967294
-DEAL:0::2_1:1 3
-DEAL:0::2_1:2 4294967294
-DEAL:0::2_1:3 3
+DEAL:0::4_1:0 4294967294
+DEAL:0::4_1:1 3
+DEAL:0::4_1:2 4294967294
+DEAL:0::4_1:3 3
 DEAL:0::6_1:0 4294967294
 DEAL:0::6_1:1 4294967294
 DEAL:0::6_1:2 4
@@ -31,9 +31,9 @@ DEAL:0::0_2:03 0
 
 DEAL:1::0_0: 0
 DEAL:1::1_0: 2
-DEAL:1::2_0: 4294967294
-DEAL:1::3_0: 2
 DEAL:1::4_0: 4294967294
+DEAL:1::2_0: 2
+DEAL:1::3_0: 4294967294
 DEAL:1::5_0: 4
 DEAL:1::6_0: 4294967294
 DEAL:1::7_0: 5
@@ -42,10 +42,10 @@ DEAL:1::0_1:0 0
 DEAL:1::0_1:1 1
 DEAL:1::0_1:2 1
 DEAL:1::0_1:3 2
-DEAL:1::2_1:0 4294967294
-DEAL:1::2_1:1 3
-DEAL:1::2_1:2 4294967294
-DEAL:1::2_1:3 3
+DEAL:1::4_1:0 4294967294
+DEAL:1::4_1:1 3
+DEAL:1::4_1:2 4294967294
+DEAL:1::4_1:3 3
 DEAL:1::6_1:0 4294967294
 DEAL:1::6_1:1 4294967294
 DEAL:1::6_1:2 4
@@ -62,9 +62,9 @@ DEAL:1::0_2:03 0
 
 DEAL:2::0_0: 0
 DEAL:2::1_0: 2
-DEAL:2::2_0: 3
+DEAL:2::4_0: 3
+DEAL:2::2_0: 2
 DEAL:2::3_0: 2
-DEAL:2::4_0: 2
 DEAL:2::5_0: 4
 DEAL:2::6_0: 4
 DEAL:2::7_0: 5
@@ -73,10 +73,10 @@ DEAL:2::0_1:0 0
 DEAL:2::0_1:1 1
 DEAL:2::0_1:2 1
 DEAL:2::0_1:3 2
-DEAL:2::2_1:0 3
-DEAL:2::2_1:1 4294967294
-DEAL:2::2_1:2 3
-DEAL:2::2_1:3 3
+DEAL:2::4_1:0 3
+DEAL:2::4_1:1 4294967294
+DEAL:2::4_1:2 3
+DEAL:2::4_1:3 3
 DEAL:2::6_1:0 4
 DEAL:2::6_1:1 4
 DEAL:2::6_1:2 4294967294
@@ -93,9 +93,9 @@ DEAL:2::0_2:03 0
 
 DEAL:3::0_0: 0
 DEAL:3::1_0: 2
-DEAL:3::2_0: 3
+DEAL:3::4_0: 3
+DEAL:3::2_0: 2
 DEAL:3::3_0: 2
-DEAL:3::4_0: 2
 DEAL:3::5_0: 4
 DEAL:3::6_0: 4
 DEAL:3::7_0: 5
@@ -104,10 +104,10 @@ DEAL:3::0_1:0 0
 DEAL:3::0_1:1 4294967294
 DEAL:3::0_1:2 1
 DEAL:3::0_1:3 4294967294
-DEAL:3::2_1:0 3
-DEAL:3::2_1:1 3
-DEAL:3::2_1:2 3
-DEAL:3::2_1:3 3
+DEAL:3::4_1:0 3
+DEAL:3::4_1:1 3
+DEAL:3::4_1:2 3
+DEAL:3::4_1:3 3
 DEAL:3::6_1:0 4294967294
 DEAL:3::6_1:1 4294967294
 DEAL:3::6_1:2 4
@@ -124,9 +124,9 @@ DEAL:3::0_2:03 4294967294
 
 DEAL:4::0_0: 0
 DEAL:4::1_0: 2
-DEAL:4::2_0: 3
+DEAL:4::4_0: 3
+DEAL:4::2_0: 2
 DEAL:4::3_0: 2
-DEAL:4::4_0: 2
 DEAL:4::5_0: 4
 DEAL:4::6_0: 4
 DEAL:4::7_0: 5
@@ -135,10 +135,10 @@ DEAL:4::0_1:0 0
 DEAL:4::0_1:1 1
 DEAL:4::0_1:2 1
 DEAL:4::0_1:3 4294967294
-DEAL:4::2_1:0 4294967294
-DEAL:4::2_1:1 3
-DEAL:4::2_1:2 3
-DEAL:4::2_1:3 3
+DEAL:4::4_1:0 4294967294
+DEAL:4::4_1:1 3
+DEAL:4::4_1:2 3
+DEAL:4::4_1:3 3
 DEAL:4::6_1:0 4
 DEAL:4::6_1:1 4
 DEAL:4::6_1:2 4
@@ -155,9 +155,9 @@ DEAL:4::0_2:03 4294967294
 
 DEAL:5::0_0: 0
 DEAL:5::1_0: 2
-DEAL:5::2_0: 3
+DEAL:5::4_0: 3
+DEAL:5::2_0: 2
 DEAL:5::3_0: 2
-DEAL:5::4_0: 2
 DEAL:5::5_0: 4
 DEAL:5::6_0: 4
 DEAL:5::7_0: 5
@@ -166,10 +166,10 @@ DEAL:5::0_1:0 4294967294
 DEAL:5::0_1:1 1
 DEAL:5::0_1:2 4294967294
 DEAL:5::0_1:3 4294967294
-DEAL:5::2_1:0 3
-DEAL:5::2_1:1 4294967294
-DEAL:5::2_1:2 4294967294
-DEAL:5::2_1:3 4294967294
+DEAL:5::4_1:0 3
+DEAL:5::4_1:1 4294967294
+DEAL:5::4_1:2 4294967294
+DEAL:5::4_1:3 4294967294
 DEAL:5::6_1:0 4294967294
 DEAL:5::6_1:1 4
 DEAL:5::6_1:2 4294967294
@@ -182,9 +182,9 @@ DEAL:5::8_1:3 4294967294
 
 DEAL:6::0_0: 0
 DEAL:6::1_0: 2
-DEAL:6::2_0: 3
+DEAL:6::4_0: 3
+DEAL:6::2_0: 2
 DEAL:6::3_0: 2
-DEAL:6::4_0: 2
 DEAL:6::5_0: 4
 DEAL:6::6_0: 4
 DEAL:6::7_0: 5
@@ -193,10 +193,10 @@ DEAL:6::0_1:0 0
 DEAL:6::0_1:1 4294967294
 DEAL:6::0_1:2 4294967294
 DEAL:6::0_1:3 4294967294
-DEAL:6::2_1:0 3
-DEAL:6::2_1:1 3
-DEAL:6::2_1:2 4294967294
-DEAL:6::2_1:3 4294967294
+DEAL:6::4_1:0 3
+DEAL:6::4_1:1 3
+DEAL:6::4_1:2 4294967294
+DEAL:6::4_1:3 4294967294
 DEAL:6::6_1:0 4
 DEAL:6::6_1:1 4294967294
 DEAL:6::6_1:2 4
index ddd7d7a5683ba14f2af34923d7b8d28560c05edd..8c92c472ecfd83fc1343da94d648998b00fd7ffc 100644 (file)
@@ -73,12 +73,10 @@ compare_meshes(DoFHandler<dim> &shared_dof_handler,
       if (cell->subdomain_id() == numbers::artificial_subdomain_id)
         continue;
 
-      typename Triangulation<dim>::active_cell_iterator tria_shared_cell =
-        cell->id().to_cell(shared_dof_handler.get_triangulation());
       typename DoFHandler<dim>::active_cell_iterator dof_shared_cell(
         &shared_dof_handler.get_triangulation(),
-        tria_shared_cell->level(),
-        tria_shared_cell->index(),
+        cell->level(),
+        cell->index(),
         &shared_dof_handler);
 
       std::vector<types::global_dof_index> distributed_cell_dofs(

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.