]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use get_communicator() 11849/head
authorPeter Munch <peterrmuench@gmail.com>
Thu, 4 Mar 2021 20:34:35 +0000 (21:34 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Thu, 4 Mar 2021 20:34:35 +0000 (21:34 +0100)
include/deal.II/dofs/dof_handler.h
include/deal.II/multigrid/mg_transfer.h
include/deal.II/multigrid/mg_transfer.templates.h
include/deal.II/multigrid/mg_transfer_matrix_free.h
include/deal.II/numerics/vector_tools_integrate_difference.templates.h
source/multigrid/mg_level_global_transfer.cc
source/multigrid/mg_transfer_internal.cc
source/multigrid/mg_transfer_prebuilt.cc
source/particles/particle_handler.cc

index df5a1ec44bf9d3b891129510067909c934fff4b7..21687df85ee1a00fff3687ce37fa9121186ffc47 100644 (file)
@@ -1945,20 +1945,10 @@ DoFHandler<dim, spacedim>::n_locally_owned_dofs_per_processor() const
   if (number_cache.n_locally_owned_dofs_per_processor.empty() &&
       number_cache.n_global_dofs > 0)
     {
-      MPI_Comm comm;
-
-      const parallel::TriangulationBase<dim, spacedim> *tr =
-        (dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
-          &this->get_triangulation()));
-      if (tr != nullptr)
-        comm = tr->get_communicator();
-      else
-        comm = MPI_COMM_SELF;
-
       const_cast<dealii::internal::DoFHandlerImplementation::NumberCache &>(
         number_cache)
         .n_locally_owned_dofs_per_processor =
-        number_cache.get_n_locally_owned_dofs_per_processor(comm);
+        number_cache.get_n_locally_owned_dofs_per_processor(get_communicator());
     }
   return number_cache.n_locally_owned_dofs_per_processor;
 }
@@ -1972,20 +1962,10 @@ DoFHandler<dim, spacedim>::locally_owned_dofs_per_processor() const
   if (number_cache.locally_owned_dofs_per_processor.empty() &&
       number_cache.n_global_dofs > 0)
     {
-      MPI_Comm comm;
-
-      const parallel::TriangulationBase<dim, spacedim> *tr =
-        (dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
-          &this->get_triangulation()));
-      if (tr != nullptr)
-        comm = tr->get_communicator();
-      else
-        comm = MPI_COMM_SELF;
-
       const_cast<dealii::internal::DoFHandlerImplementation::NumberCache &>(
         number_cache)
         .locally_owned_dofs_per_processor =
-        number_cache.get_locally_owned_dofs_per_processor(comm);
+        number_cache.get_locally_owned_dofs_per_processor(get_communicator());
     }
   return number_cache.locally_owned_dofs_per_processor;
 }
@@ -2007,20 +1987,11 @@ DoFHandler<dim, spacedim>::locally_owned_mg_dofs_per_processor(
   if (mg_number_cache[level].locally_owned_dofs_per_processor.empty() &&
       mg_number_cache[level].n_global_dofs > 0)
     {
-      MPI_Comm comm;
-
-      const parallel::TriangulationBase<dim, spacedim> *tr =
-        (dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
-          &this->get_triangulation()));
-      if (tr != nullptr)
-        comm = tr->get_communicator();
-      else
-        comm = MPI_COMM_SELF;
-
       const_cast<dealii::internal::DoFHandlerImplementation::NumberCache &>(
         mg_number_cache[level])
         .locally_owned_dofs_per_processor =
-        mg_number_cache[level].get_locally_owned_dofs_per_processor(comm);
+        mg_number_cache[level].get_locally_owned_dofs_per_processor(
+          get_communicator());
     }
   return mg_number_cache[level].locally_owned_dofs_per_processor;
 }
index 370a3c156e8af9512c701c5ef54dd8861e936841..ab979647609ee276dbd3b4e0b8f591bc448eeeec 100644 (file)
@@ -84,12 +84,7 @@ namespace internal
            const SparsityPatternType &      sp,
            const DoFHandler<dim, spacedim> &dh)
     {
-      const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria =
-        dynamic_cast<
-          const dealii::parallel::TriangulationBase<dim, spacedim> *>(
-          &(dh.get_triangulation()));
-      MPI_Comm communicator =
-        dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF;
+      const MPI_Comm communicator = dh.get_communicator();
 
       matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
                     dh.locally_owned_mg_dofs(level),
@@ -115,12 +110,8 @@ namespace internal
            const SparsityPatternType &      sp,
            const DoFHandler<dim, spacedim> &dh)
     {
-      const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria =
-        dynamic_cast<
-          const dealii::parallel::TriangulationBase<dim, spacedim> *>(
-          &(dh.get_triangulation()));
-      MPI_Comm communicator =
-        dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF;
+      const MPI_Comm communicator = dh.get_communicator();
+
       matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
                     dh.locally_owned_mg_dofs(level),
                     sp,
@@ -147,12 +138,8 @@ namespace internal
            const SparsityPatternType &      sp,
            const DoFHandler<dim, spacedim> &dh)
     {
-      const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria =
-        dynamic_cast<
-          const dealii::parallel::TriangulationBase<dim, spacedim> *>(
-          &(dh.get_triangulation()));
-      MPI_Comm communicator =
-        dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF;
+      const MPI_Comm communicator = dh.get_communicator();
+
       matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
                     dh.locally_owned_mg_dofs(level),
                     sp,
@@ -178,12 +165,8 @@ namespace internal
            const SparsityPatternType &      sp,
            const DoFHandler<dim, spacedim> &dh)
     {
-      const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria =
-        dynamic_cast<
-          const dealii::parallel::TriangulationBase<dim, spacedim> *>(
-          &(dh.get_triangulation()));
-      MPI_Comm communicator =
-        dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF;
+      const MPI_Comm communicator = dh.get_communicator();
+
       matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
                     dh.locally_owned_mg_dofs(level),
                     sp,
@@ -238,12 +221,8 @@ namespace internal
            const SparsityPatternType &      sp,
            const DoFHandler<dim, spacedim> &dh)
     {
-      const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria =
-        dynamic_cast<
-          const dealii::parallel::TriangulationBase<dim, spacedim> *>(
-          &(dh.get_triangulation()));
-      MPI_Comm communicator =
-        dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF;
+      const MPI_Comm communicator = dh.get_communicator();
+
       // Reinit PETSc matrix
       matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
                     dh.locally_owned_mg_dofs(level),
index 3554957ffec0b387a9b54de9b4316e5f9dd005ab..9153d0c915888453829a42ed1bc5fc52a79a5b6b 100644 (file)
@@ -443,13 +443,8 @@ MGLevelGlobalTransfer<LinearAlgebra::distributed::Vector<Number>>::copy_to_mg(
           dst[level].reinit(ghosted_level_vector[level], false);
         else
           {
-            const dealii::parallel::TriangulationBase<dim, spacedim> *tria =
-              (dynamic_cast<
-                const dealii::parallel::TriangulationBase<dim, spacedim> *>(
-                &dof_handler.get_triangulation()));
             dst[level].reinit(dof_handler.locally_owned_mg_dofs(level),
-                              tria != nullptr ? tria->get_communicator() :
-                                                MPI_COMM_SELF);
+                              dof_handler.get_communicator());
           }
       }
     else if ((perform_plain_copy == false &&
index 210aa36a1bc089da8b49dbd61a9fe662d5d1ce3e..be7d526ea57d147eebd0757ee588c64e797fe72f 100644 (file)
@@ -674,8 +674,7 @@ MGTransferBlockMatrixFree<dim, Number>::copy_to_mg(
                 LinearAlgebra::distributed::Vector<Number> &v =
                   dst[level].block(b);
                 v.reinit(dof_handler[b]->locally_owned_mg_dofs(level),
-                         tria != nullptr ? tria->get_communicator() :
-                                           MPI_COMM_SELF);
+                         dof_handler[b]->get_communicator());
               }
             dst[level].collect_sizes();
           }
index 4e77bd2fced5fe0a580ed75b124e5eef94f1c973..29c7a84b9df4a98a1e36ad066cee03ccd8d9c75f 100644 (file)
@@ -1260,13 +1260,7 @@ namespace VectorTools
     }
 #endif
 
-    MPI_Comm comm = MPI_COMM_SELF;
-#ifdef DEAL_II_WITH_MPI
-    if (const parallel::TriangulationBase<dim, spacedim> *ptria =
-          dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
-            &tria))
-      comm = ptria->get_communicator();
-#endif
+    const MPI_Comm comm = tria.get_communicator();
 
     switch (norm)
       {
index 51ffc58d02efb0dba7b70f9613de7ee2693d97cc..abcc0d3d580682f30d32475f4fe2b1da5b7ae7b7 100644 (file)
@@ -284,11 +284,7 @@ void
 MGLevelGlobalTransfer<LinearAlgebra::distributed::Vector<Number>>::
   fill_and_communicate_copy_indices(const DoFHandler<dim, spacedim> &mg_dof)
 {
-  const parallel::TriangulationBase<dim, spacedim> *ptria =
-    dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
-      &mg_dof.get_triangulation());
-  const MPI_Comm mpi_communicator =
-    ptria != nullptr ? ptria->get_communicator() : MPI_COMM_SELF;
+  const MPI_Comm mpi_communicator = mg_dof.get_communicator();
 
   fill_internal(mg_dof,
                 mg_constrained_dofs,
index 99106b8e3ccdd8603f5405a8502b6d8700c9edab..d6405b17fc69c5e0fd0904737de2552e8bf9bf86 100644 (file)
@@ -946,18 +946,12 @@ namespace internal
           // the base class for keeping ghosted transfer indices. To avoid
           // keeping two very similar vectors, we keep one single ghosted
           // vector that is augmented/filled here.
-          const dealii::parallel::TriangulationBase<dim, dim> *ptria =
-            (dynamic_cast<
-              const dealii::parallel::TriangulationBase<dim, dim> *>(&tria));
-          const MPI_Comm communicator =
-            ptria != nullptr ? ptria->get_communicator() : MPI_COMM_SELF;
-
           reinit_level_partitioner(dof_handler.locally_owned_mg_dofs(level),
                                    ghosted_level_dofs,
                                    external_partitioners.empty() ?
                                      nullptr :
                                      external_partitioners[level],
-                                   communicator,
+                                   tria.get_communicator(),
                                    target_partitioners[level],
                                    copy_indices_global_mine[level]);
 
@@ -976,7 +970,7 @@ namespace internal
                                        external_partitioners.empty() ?
                                          nullptr :
                                          external_partitioners[0],
-                                       communicator,
+                                       tria.get_communicator(),
                                        target_partitioners[0],
                                        copy_indices_global_mine[0]);
 
index 472942fa46ef4780f791a70856f016e6ad464904..5c576eca6521c6cca60a76f8be23a7276420715e 100644 (file)
@@ -266,20 +266,11 @@ MGTransferPrebuilt<VectorType>::build(
           // complete sparsity patterns on their own, the sparsity pattern must
           // be manually distributed.
 
-          // Retrieve communicator from triangulation if it is parallel
-          const parallel::TriangulationBase<dim, spacedim> *dist_tria =
-            dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
-              &(dof_handler.get_triangulation()));
-
-          MPI_Comm communicator = dist_tria != nullptr ?
-                                    dist_tria->get_communicator() :
-                                    MPI_COMM_SELF;
-
           // Distribute sparsity pattern
           ::dealii::SparsityTools::distribute_sparsity_pattern(
             dsp,
             dof_handler.locally_owned_mg_dofs(level + 1),
-            communicator,
+            dof_handler.get_communicator(),
             dsp.row_index_set());
         }
 #endif
index c69838d84a846ec7f19a5d94ee8e8cecff27efc7..514ba68a181c2541eab11c2305f2317170973064 100644 (file)
@@ -491,11 +491,7 @@ namespace Particles
     if (!ids.empty())
       AssertDimension(ids.size(), positions.size());
 
-    const auto tria =
-      dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
-        &(*triangulation));
-    const auto comm =
-      (tria != nullptr ? tria->get_communicator() : MPI_COMM_WORLD);
+    const auto comm = triangulation->get_communicator();
 
     const auto n_mpi_processes = Utilities::MPI::n_mpi_processes(comm);
 

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.