]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add some more AssertThrowMPIs. 12309/head
authorDavid Wells <drwells@email.unc.edu>
Tue, 25 May 2021 22:14:41 +0000 (18:14 -0400)
committerDavid Wells <drwells@email.unc.edu>
Tue, 25 May 2021 22:16:44 +0000 (18:16 -0400)
include/deal.II/base/aligned_vector.h
include/deal.II/base/mpi.h
include/deal.II/base/mpi_compute_index_owner_internal.h
include/deal.II/base/mpi_consensus_algorithms.templates.h
include/deal.II/base/mpi_noncontiguous_partitioner.templates.h
include/deal.II/base/mpi_remote_point_evaluation.h
include/deal.II/lac/la_parallel_vector.templates.h
include/deal.II/matrix_free/face_setup_internal.h
include/deal.II/matrix_free/matrix_free.h
include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h

index 0d9385108b5e750fa841d36d1f52a121d28c9193..9f216af047e77d3495585dbd2bf34149017cf441 100644 (file)
@@ -1536,7 +1536,8 @@ AlignedVector<T>::replicate_across_communicator(const MPI_Comm &   communicator,
 
   // Make sure that the shared memory host has copied the data before we try to
   // access it.
-  MPI_Barrier(shmem_group_communicator);
+  const int ierr = MPI_Barrier(shmem_group_communicator);
+  AssertThrowMPI(ierr);
 
   // **** Step 7 ****
   // Finally, we need to set the pointers of this object to what we just
index 7d99cddb9baa3c296e9b4334e2f307b07aef2fec..d0f416a5a31011eff0a7d5d8cdeb659045bf48d1 100644 (file)
@@ -1398,8 +1398,9 @@ namespace Utilities
       std::vector<int> size_all_data(n_procs, 0);
 
       // Exchanging the size of each buffer
-      MPI_Allgather(
+      int ierr = MPI_Allgather(
         &n_local_data, 1, MPI_INT, size_all_data.data(), 1, MPI_INT, comm);
+      AssertThrowMPI(ierr);
 
       // Now computing the displacement, relative to recvbuf,
       // at which to store the incoming buffer
@@ -1412,14 +1413,15 @@ namespace Utilities
       std::vector<char> received_unrolled_buffer(rdispls.back() +
                                                  size_all_data.back());
 
-      MPI_Allgatherv(buffer.data(),
-                     n_local_data,
-                     MPI_CHAR,
-                     received_unrolled_buffer.data(),
-                     size_all_data.data(),
-                     rdispls.data(),
-                     MPI_CHAR,
-                     comm);
+      ierr = MPI_Allgatherv(buffer.data(),
+                            n_local_data,
+                            MPI_CHAR,
+                            received_unrolled_buffer.data(),
+                            size_all_data.data(),
+                            rdispls.data(),
+                            MPI_CHAR,
+                            comm);
+      AssertThrowMPI(ierr);
 
       std::vector<T> received_objects(n_procs);
       for (unsigned int i = 0; i < n_procs; ++i)
index 7b00cfea4f81073014859c29bb5cc91d9561246e..71cdc86d2144709fd4395aca4dad55d79595b367 100644 (file)
@@ -359,7 +359,7 @@ namespace Utilities
                   {
                     // wait for an incoming message
                     MPI_Status status;
-                    auto       ierr =
+                    int        ierr =
                       MPI_Probe(MPI_ANY_SOURCE, mpi_tag, comm, &status);
                     AssertThrowMPI(ierr);
 
@@ -876,9 +876,8 @@ namespace Utilities
                  ++c)
               {
                 // wait for an incoming message
-                MPI_Status   status;
-                unsigned int ierr =
-                  MPI_Probe(MPI_ANY_SOURCE, mpi_tag, comm, &status);
+                MPI_Status status;
+                int ierr = MPI_Probe(MPI_ANY_SOURCE, mpi_tag, comm, &status);
                 AssertThrowMPI(ierr);
 
                 // retrieve size of incoming message
index eb6945da09331dcb00e2c7749f7feeb3c18e6a49..2fd8e28279d42acf8735e6f8ca64c5c35db39784 100644 (file)
@@ -346,19 +346,20 @@ namespace Utilities
             }
 
 
-          const int ierr = MPI_Wait(&barrier_request, MPI_STATUS_IGNORE);
+          int ierr = MPI_Wait(&barrier_request, MPI_STATUS_IGNORE);
           AssertThrowMPI(ierr);
 
           for (auto &i : request_requests)
             {
-              const auto ierr = MPI_Wait(i.get(), MPI_STATUS_IGNORE);
+              ierr = MPI_Wait(i.get(), MPI_STATUS_IGNORE);
               AssertThrowMPI(ierr);
             }
 
 #  ifdef DEBUG
           // note: IBarrier seems to make problem during testing, this
           // additional Barrier seems to help
-          MPI_Barrier(this->comm);
+          ierr = MPI_Barrier(this->comm);
+          AssertThrowMPI(ierr);
 #  endif
         }
 
index e3a1447d703492d59cc447f4352e7d18c8eff87e..38a15b4151fcef2ea0f0edd5e60595f9d50e045e 100644 (file)
@@ -106,7 +106,7 @@ namespace Utilities
       AssertIndexRange(recv_ranks.size(), recv_ptr.size());
       for (types::global_dof_index i = 0; i < recv_ranks.size(); i++)
         {
-          const auto ierr =
+          const int ierr =
             MPI_Irecv(buffers.data() + recv_ptr[i],
                       recv_ptr[i + 1] - recv_ptr[i],
                       Utilities::MPI::internal::mpi_type_id(buffers.data()),
@@ -135,7 +135,7 @@ namespace Utilities
                    (send_ptr[i] == buffers.size() &&
                     send_ptr[i + 1] == send_ptr[i]),
                  ExcMessage("The input buffer doesn't contain enough entries"));
-          const auto ierr =
+          const int ierr =
             MPI_Isend(buffers.data() + send_ptr[i],
                       send_ptr[i + 1] - send_ptr[i],
                       Utilities::MPI::internal::mpi_type_id(buffers.data()),
@@ -182,7 +182,7 @@ namespace Utilities
         }
 
       // wait that all data packages have been sent
-      const auto ierr =
+      const int ierr =
         MPI_Waitall(send_ranks.size(), requests.data(), MPI_STATUSES_IGNORE);
       AssertThrowMPI(ierr);
 #endif
index ebd3f3a5b8a513e27c7cc6020e62900570ca869a..74dc477e78072b980d0e4d0c93d5ed4fd4032c93 100644 (file)
@@ -310,13 +310,14 @@ namespace Utilities
 
           requests.push_back(MPI_Request());
 
-          MPI_Isend(buffer.data(),
-                    buffer.size(),
-                    MPI_CHAR,
-                    send_ranks[i],
-                    internal::Tags::remote_point_evaluation,
-                    tria->get_communicator(),
-                    &requests.back());
+          const int ierr = MPI_Isend(buffer.data(),
+                                     buffer.size(),
+                                     MPI_CHAR,
+                                     send_ranks[i],
+                                     internal::Tags::remote_point_evaluation,
+                                     tria->get_communicator(),
+                                     &requests.back());
+          AssertThrowMPI(ierr);
         }
 
       for (const auto recv_rank : recv_ranks)
@@ -325,30 +326,35 @@ namespace Utilities
             continue;
 
           MPI_Status status;
-          MPI_Probe(MPI_ANY_SOURCE,
-                    internal::Tags::remote_point_evaluation,
-                    tria->get_communicator(),
-                    &status);
+          int        ierr = MPI_Probe(MPI_ANY_SOURCE,
+                               internal::Tags::remote_point_evaluation,
+                               tria->get_communicator(),
+                               &status);
+          AssertThrowMPI(ierr);
 
           int message_length;
-          MPI_Get_count(&status, MPI_CHAR, &message_length);
+          ierr = MPI_Get_count(&status, MPI_CHAR, &message_length);
+          AssertThrowMPI(ierr);
 
           std::vector<char> buffer(message_length);
 
-          MPI_Recv(buffer.data(),
-                   buffer.size(),
-                   MPI_CHAR,
-                   status.MPI_SOURCE,
-                   internal::Tags::remote_point_evaluation,
-                   tria->get_communicator(),
-                   MPI_STATUS_IGNORE);
+          ierr = MPI_Recv(buffer.data(),
+                          buffer.size(),
+                          MPI_CHAR,
+                          status.MPI_SOURCE,
+                          internal::Tags::remote_point_evaluation,
+                          tria->get_communicator(),
+                          MPI_STATUS_IGNORE);
+          AssertThrowMPI(ierr);
 
           temp_recv_map[status.MPI_SOURCE] =
             Utilities::unpack<std::vector<T>>(buffer, false);
         }
 
       // make sure all messages have been sent
-      MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
+      const int ierr =
+        MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
+      AssertThrowMPI(ierr);
 
       // copy received data into output vector
       auto it = recv_permutation.begin();
@@ -443,13 +449,14 @@ namespace Utilities
 
           requests.push_back(MPI_Request());
 
-          MPI_Isend(buffer_send.data(),
-                    buffer_send.size(),
-                    MPI_CHAR,
-                    recv_rank,
-                    internal::Tags::remote_point_evaluation,
-                    tria->get_communicator(),
-                    &requests.back());
+          const int ierr = MPI_Isend(buffer_send.data(),
+                                     buffer_send.size(),
+                                     MPI_CHAR,
+                                     recv_rank,
+                                     internal::Tags::remote_point_evaluation,
+                                     tria->get_communicator(),
+                                     &requests.back());
+          AssertThrowMPI(ierr);
         }
 
       for (unsigned int i = 0; i < send_ranks.size(); ++i)
@@ -474,24 +481,26 @@ namespace Utilities
             }
 
           MPI_Status status;
-          MPI_Probe(MPI_ANY_SOURCE,
-                    internal::Tags::remote_point_evaluation,
-                    tria->get_communicator(),
-                    &status);
+          int        ierr = MPI_Probe(MPI_ANY_SOURCE,
+                               internal::Tags::remote_point_evaluation,
+                               tria->get_communicator(),
+                               &status);
+          AssertThrowMPI(ierr);
 
           int message_length;
-          MPI_Get_count(&status, MPI_CHAR, &message_length);
+          ierr = MPI_Get_count(&status, MPI_CHAR, &message_length);
+          AssertThrowMPI(ierr);
 
           std::vector<char> recv_buffer(message_length);
 
-          MPI_Recv(recv_buffer.data(),
-                   recv_buffer.size(),
-                   MPI_CHAR,
-                   status.MPI_SOURCE,
-                   internal::Tags::remote_point_evaluation,
-                   tria->get_communicator(),
-                   MPI_STATUS_IGNORE);
-
+          ierr = MPI_Recv(recv_buffer.data(),
+                          recv_buffer.size(),
+                          MPI_CHAR,
+                          status.MPI_SOURCE,
+                          internal::Tags::remote_point_evaluation,
+                          tria->get_communicator(),
+                          MPI_STATUS_IGNORE);
+          AssertThrowMPI(ierr);
 
           const auto recv_buffer_unpacked =
             Utilities::unpack<std::vector<T>>(recv_buffer, false);
@@ -515,7 +524,9 @@ namespace Utilities
             }
         }
 
-      MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
+      const int ierr =
+        MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
+      AssertThrowMPI(ierr);
 
       // sort for easy access during function call
       for (unsigned int i = 0; i < send_permutation.size(); ++i)
index bf55d336bdc87e93d2cba8b7a1b66f763c5ddd9d..4073de08297ce491ebcb741e954a294527826ed7 100644 (file)
@@ -164,9 +164,11 @@ namespace LinearAlgebra
               std::vector<Number *> others(size_sm);
 
               MPI_Info info;
-              MPI_Info_create(&info);
+              int      ierr = MPI_Info_create(&info);
+              AssertThrowMPI(ierr);
 
-              MPI_Info_set(info, "alloc_shared_noncontig", "true");
+              ierr = MPI_Info_set(info, "alloc_shared_noncontig", "true");
+              AssertThrowMPI(ierr);
 
               const std::size_t align_by = 64;
 
@@ -175,7 +177,7 @@ namespace LinearAlgebra
                  sizeof(Number)) *
                 sizeof(Number);
 
-              auto ierr = MPI_Win_allocate_shared(
+              ierr = MPI_Win_allocate_shared(
                 s, sizeof(Number), info, comm_shared, &data_this, &mpi_window);
               AssertThrowMPI(ierr);
 
index 89935c6f4caba28cf0e9c1d6179d15472cac052a..504bd6dc634a6f5856f5b6af9d198de7e9572f32 100644 (file)
@@ -299,46 +299,50 @@ namespace internal
               MPI_Status   status;
               unsigned int mysize    = inner_face.second.shared_faces.size();
               unsigned int othersize = numbers::invalid_unsigned_int;
-              MPI_Sendrecv(&mysize,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           600 + my_domain,
-                           &othersize,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           600 + inner_face.first,
-                           comm,
-                           &status);
+
+              int ierr = MPI_Sendrecv(&mysize,
+                                      1,
+                                      MPI_UNSIGNED,
+                                      inner_face.first,
+                                      600 + my_domain,
+                                      &othersize,
+                                      1,
+                                      MPI_UNSIGNED,
+                                      inner_face.first,
+                                      600 + inner_face.first,
+                                      comm,
+                                      &status);
+              AssertThrowMPI(ierr);
               AssertDimension(mysize, othersize);
               mysize = inner_face.second.n_hanging_faces_smaller_subdomain;
-              MPI_Sendrecv(&mysize,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           700 + my_domain,
-                           &othersize,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           700 + inner_face.first,
-                           comm,
-                           &status);
+              ierr   = MPI_Sendrecv(&mysize,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  700 + my_domain,
+                                  &othersize,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  700 + inner_face.first,
+                                  comm,
+                                  &status);
+              AssertThrowMPI(ierr);
               AssertDimension(mysize, othersize);
               mysize = inner_face.second.n_hanging_faces_larger_subdomain;
-              MPI_Sendrecv(&mysize,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           800 + my_domain,
-                           &othersize,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           800 + inner_face.first,
-                           comm,
-                           &status);
+              ierr   = MPI_Sendrecv(&mysize,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  800 + my_domain,
+                                  &othersize,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  800 + inner_face.first,
+                                  comm,
+                                  &status);
+              AssertThrowMPI(ierr);
               AssertDimension(mysize, othersize);
 #  endif
 
@@ -463,44 +467,47 @@ namespace internal
 
                 // make sure the splitting is consistent between both sides
 #  if defined(DEAL_II_WITH_MPI) && defined(DEBUG)
-              MPI_Sendrecv(&split_index,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           900 + my_domain,
-                           &othersize,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           900 + inner_face.first,
-                           comm,
-                           &status);
+              ierr = MPI_Sendrecv(&split_index,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  900 + my_domain,
+                                  &othersize,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  900 + inner_face.first,
+                                  comm,
+                                  &status);
+              AssertThrowMPI(ierr);
               AssertDimension(split_index, othersize);
-              MPI_Sendrecv(&n_faces_lower_proc,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           1000 + my_domain,
-                           &othersize,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           1000 + inner_face.first,
-                           comm,
-                           &status);
+              ierr = MPI_Sendrecv(&n_faces_lower_proc,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  1000 + my_domain,
+                                  &othersize,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  1000 + inner_face.first,
+                                  comm,
+                                  &status);
+              AssertThrowMPI(ierr);
               AssertDimension(n_faces_lower_proc, othersize);
-              MPI_Sendrecv(&n_faces_higher_proc,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           1100 + my_domain,
-                           &othersize,
-                           1,
-                           MPI_UNSIGNED,
-                           inner_face.first,
-                           1100 + inner_face.first,
-                           comm,
-                           &status);
+              ierr = MPI_Sendrecv(&n_faces_higher_proc,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  1100 + my_domain,
+                                  &othersize,
+                                  1,
+                                  MPI_UNSIGNED,
+                                  inner_face.first,
+                                  1100 + inner_face.first,
+                                  comm,
+                                  &status);
+              AssertThrowMPI(ierr);
               AssertDimension(n_faces_higher_proc, othersize);
 #  endif
 
index 53fe0202befc1c1d777588206f28156964906f44..54ed7d273df85cf1f0606778751233c88736d9f5 100644 (file)
@@ -3741,7 +3741,11 @@ namespace internal
             }
 
           if (Utilities::MPI::job_supports_mpi())
-            MPI_Barrier(matrix_free.get_task_info().communicator_sm);
+            {
+              const int ierr =
+                MPI_Barrier(matrix_free.get_task_info().communicator_sm);
+              AssertThrowMPI(ierr);
+            }
 #  endif
         }
     }
index d38b2eb6e35a602914ef148b18807a5da3cbb894..ae60aa915480edc2b6d332a395daa1ff44707737 100644 (file)
@@ -723,7 +723,7 @@ namespace internal
         for (unsigned int i = 0; i < ranks.size(); ++i)
           {
             MPI_Status status;
-            const auto ierr_1 = MPI_Probe(
+            const int  ierr_1 = MPI_Probe(
               MPI_ANY_SOURCE,
               Utilities::MPI::internal::Tags::fine_dof_handler_view_reinit,
               communicator,
@@ -732,8 +732,8 @@ namespace internal
 
             std::vector<types::global_dof_index> buffer;
 
-            int        message_length;
-            const auto ierr_2 =
+            int       message_length;
+            const int ierr_2 =
               MPI_Get_count(&status,
                             Utilities::MPI::internal::mpi_type_id(
                               buffer.data()),
@@ -742,7 +742,7 @@ namespace internal
 
             buffer.resize(message_length);
 
-            const auto ierr_3 = MPI_Recv(
+            const int ierr_3 = MPI_Recv(
               buffer.data(),
               buffer.size(),
               Utilities::MPI::internal::mpi_type_id(buffer.data()),
@@ -779,7 +779,7 @@ namespace internal
               }
           }
 
-        const auto ierr_1 =
+        const int ierr_1 =
           MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
         AssertThrowMPI(ierr_1);
       }

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.