]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Simplify code. 6007/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 8 Mar 2018 14:37:19 +0000 (07:37 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 8 Mar 2018 14:37:19 +0000 (07:37 -0700)
The existing code seemed to suggest that something different had to happen if n_levels()>0.
But it was really that there was nothing to do if the triangulation is empty. Reflect this
in both the code and a comment.

source/distributed/tria_base.cc

index c75d95a00493fdbfd64940a0b25195959fc3a30e..ea1a5af0b2730c970ef0c1aef457116ed712d623 100644 (file)
@@ -224,65 +224,68 @@ namespace parallel
   {
     number_cache.level_ghost_owners.clear ();
 
-    if (this->n_levels() > 0)
-      {
-        // find level ghost owners
-        for (typename Triangulation<dim,spacedim>::cell_iterator
-             cell = this->begin();
-             cell != this->end();
-             ++cell)
-          if (cell->level_subdomain_id() != numbers::artificial_subdomain_id
-              && cell->level_subdomain_id() != this->locally_owned_subdomain())
-            this->number_cache.level_ghost_owners.insert(cell->level_subdomain_id());
+    // if there is nothing to do, then do nothing
+    if (this->n_levels() == 0)
+      return;
 
-#ifdef DEBUG
-        // Check that level_ghost_owners is symmetric by sending a message to everyone
-        {
-          int ierr = MPI_Barrier(this->mpi_communicator);
-          AssertThrowMPI(ierr);
+    // find level ghost owners
+    for (typename Triangulation<dim,spacedim>::cell_iterator
+         cell = this->begin();
+         cell != this->end();
+         ++cell)
+      if (cell->level_subdomain_id() != numbers::artificial_subdomain_id
+          && cell->level_subdomain_id() != this->locally_owned_subdomain())
+        this->number_cache.level_ghost_owners.insert(cell->level_subdomain_id());
 
-          // important: preallocate to avoid (re)allocation:
-          std::vector<MPI_Request> requests (this->number_cache.level_ghost_owners.size());
-          unsigned int dummy = 0;
-          unsigned int req_counter = 0;
+#ifdef DEBUG
+    // Check that level_ghost_owners is symmetric by sending a message to everyone
+    {
+      int ierr = MPI_Barrier(this->mpi_communicator);
+      AssertThrowMPI(ierr);
 
-          for (std::set<types::subdomain_id>::iterator it = this->number_cache.level_ghost_owners.begin();
-               it != this->number_cache.level_ghost_owners.end();
-               ++it, ++req_counter)
-            {
-              ierr = MPI_Isend(&dummy, 1, MPI_UNSIGNED,
-                               *it, 9001, this->mpi_communicator,
-                               &requests[req_counter]);
-              AssertThrowMPI(ierr);
-            }
+      // important: preallocate to avoid (re)allocation:
+      std::vector<MPI_Request> requests (this->number_cache.level_ghost_owners.size());
+      unsigned int dummy = 0;
+      unsigned int req_counter = 0;
 
-          for (std::set<types::subdomain_id>::iterator it = this->number_cache.level_ghost_owners.begin();
-               it != this->number_cache.level_ghost_owners.end();
-               ++it)
-            {
-              unsigned int dummy;
-              ierr = MPI_Recv(&dummy, 1, MPI_UNSIGNED,
-                              *it, 9001, this->mpi_communicator,
-                              MPI_STATUS_IGNORE);
-              AssertThrowMPI(ierr);
-            }
+      for (std::set<types::subdomain_id>::iterator it = this->number_cache.level_ghost_owners.begin();
+           it != this->number_cache.level_ghost_owners.end();
+           ++it, ++req_counter)
+        {
+          ierr = MPI_Isend(&dummy, 1, MPI_UNSIGNED,
+                           *it, 9001, this->mpi_communicator,
+                           &requests[req_counter]);
+          AssertThrowMPI(ierr);
+        }
 
-          if (requests.size() > 0)
-            {
-              ierr = MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
-              AssertThrowMPI(ierr);
-            }
+      for (std::set<types::subdomain_id>::iterator it = this->number_cache.level_ghost_owners.begin();
+           it != this->number_cache.level_ghost_owners.end();
+           ++it)
+        {
+          unsigned int dummy;
+          ierr = MPI_Recv(&dummy, 1, MPI_UNSIGNED,
+                          *it, 9001, this->mpi_communicator,
+                          MPI_STATUS_IGNORE);
+          AssertThrowMPI(ierr);
+        }
 
-          ierr = MPI_Barrier(this->mpi_communicator);
+      if (requests.size() > 0)
+        {
+          ierr = MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
           AssertThrowMPI(ierr);
         }
+
+      ierr = MPI_Barrier(this->mpi_communicator);
+      AssertThrowMPI(ierr);
+    }
 #endif
 
-        Assert(this->number_cache.level_ghost_owners.size() < Utilities::MPI::n_mpi_processes(this->mpi_communicator),
-               ExcInternalError());
-      }
+    Assert(this->number_cache.level_ghost_owners.size() < Utilities::MPI::n_mpi_processes(this->mpi_communicator),
+           ExcInternalError());
   }
+
 #else
+
   template <int dim, int spacedim>
   void
   Triangulation<dim,spacedim>::update_number_cache ()

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.