]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Query the size of the reply message from MPI.
authorWolfgang Bangerth <bangerth@colostate.edu>
Tue, 11 Jan 2022 19:09:20 +0000 (12:09 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Tue, 18 Jan 2022 20:47:51 +0000 (13:47 -0700)
This avoids requiring the user to do so by hand.

include/deal.II/base/mpi_consensus_algorithms.templates.h

index 54aa6e8daf7c480c487f5933d5da429ea9527762..c6d21aa18936c8133946f5286b6bf3540790475c 100644 (file)
@@ -424,36 +424,45 @@ namespace Utilities
               const int tag_deliver = Utilities::MPI::internal::Tags::
                 consensus_algorithm_nbx_process_deliver;
 
-#  ifdef DEBUG
               // Just to be sure, double check that the message really
               // is already here -- in other words, that the following
-              // MPI_Recv is going to return immediately
+              // MPI_Recv is going to return immediately. But as part of
+              // this, we also use the status object to query the size
+              // of the message so that we can resize the receive buffer.
+              int        request_is_pending;
+              MPI_Status status;
               {
-                int        request_is_pending;
-                const auto ierr = MPI_Iprobe(target,
-                                             tag_deliver,
-                                             this->comm,
-                                             &request_is_pending,
-                                             MPI_STATUS_IGNORE);
+                const int ierr = MPI_Iprobe(target,
+                                            tag_deliver,
+                                            this->comm,
+                                            &request_is_pending,
+                                            &status);
                 AssertThrowMPI(ierr);
-
-                (void)request_is_pending;
-                Assert(request_is_pending, ExcInternalError());
               }
-#  endif
+
+              (void)request_is_pending;
+              Assert(request_is_pending, ExcInternalError());
 
               // OK, so yes, a message is here. Receive it.
-              std::vector<T2> recv_buffer;
-              this->process.prepare_buffer_for_answer(target, recv_buffer);
-
-              const int ierr = MPI_Recv(recv_buffer.data(),
-                                        recv_buffer.size() * sizeof(T2),
-                                        MPI_BYTE,
-                                        target,
-                                        tag_deliver,
-                                        this->comm,
-                                        MPI_STATUS_IGNORE);
-              AssertThrowMPI(ierr);
+              int message_size;
+              {
+                const int ierr =
+                  MPI_Get_count(&status, MPI_BYTE, &message_size);
+                AssertThrowMPI(ierr);
+              }
+              Assert(message_size % sizeof(T2) == 0, ExcInternalError());
+              std::vector<T2> recv_buffer(message_size / sizeof(T2));
+
+              {
+                const int ierr = MPI_Recv(recv_buffer.data(),
+                                          recv_buffer.size() * sizeof(T2),
+                                          MPI_BYTE,
+                                          target,
+                                          tag_deliver,
+                                          this->comm,
+                                          MPI_STATUS_IGNORE);
+                AssertThrowMPI(ierr);
+              }
 
               this->process.read_answer(target, recv_buffer);
             }

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.