]> https://gitweb.dealii.org/ - dealii.git/commitdiff
address comments
authorTimo Heister <timo.heister@gmail.com>
Thu, 18 Nov 2021 19:32:56 +0000 (14:32 -0500)
committerTimo Heister <timo.heister@gmail.com>
Thu, 18 Nov 2021 19:32:56 +0000 (14:32 -0500)
include/deal.II/base/mpi.h
source/base/mpi.cc
tests/mpi/create_mpi_datatype_01.cc

index 406dde250afa62fd57af384ca5ff4a5d537e6950..3449bebec3de70c0aa7689b456b458f3afbcf07c 100644 (file)
@@ -521,15 +521,15 @@ namespace Utilities
      *   MPI_Send(buffer.data(), buffer.size(), MPI_BYTE, dest, tag, comm);
      * else
      * {
-     *   MPI_Datatype bigtype;
-     *   Utilities::MPI::create_mpi_data_type_n_bytes(bigtype, buffer.size());
+     *   MPI_Datatype bigtype =
+     * Utilities::MPI::create_mpi_data_type_n_bytes(buffer.size());
      *   MPI_Send(buffer.data(), 1, bigtype, dest, tag, comm);
      *   MPI_Type_free(&bigtype);
      * }
      * </code>
      */
-    void
-    create_mpi_data_type_n_bytes(MPI_Datatype &result, std::size_t n_bytes);
+    MPI_Datatype
+    create_mpi_data_type_n_bytes(std::size_t n_bytes);
 
     /**
      * Return the sum over all processors of the value @p t. This function is
index b39d298fcee8c78a7299533bee3c76f5d89571c6..bcb7bcb0eeb653f46b32b1d895cf75884ba6276c 100644 (file)
@@ -297,18 +297,18 @@ namespace Utilities
 
 
 
-    void
-    create_mpi_data_type_n_bytes(MPI_Datatype &result, std::size_t n_bytes)
+    MPI_Datatype
+    create_mpi_data_type_n_bytes(std::size_t n_bytes)
     {
 #  ifdef DEAL_II_WITH_MPI
       // Simplified version from BigMPI repository, see
       // https://github.com/jeffhammond/BigMPI/blob/5300b18cc8ec1b2431bf269ee494054ee7bd9f72/src/type_contiguous_x.c#L74
       // (code is MIT licensed)
 
-      // We create an MPI datatype that has the layout A*nB where A is
+      // We create an MPI datatype that has the layout A*n+B where A is
       // max_signed_int bytes repeated n times and B is the remainder.
 
-      const MPI_Count max_signed_int = (1U << 31) - 1;
+      const MPI_Count max_signed_int = std::numeric_limits<int>::max();
 
       const MPI_Count n_chunks          = n_bytes / max_signed_int;
       const MPI_Count n_bytes_remainder = n_bytes % max_signed_int;
@@ -340,6 +340,8 @@ namespace Utilities
         ExcMessage(
           "ERROR in create_mpi_data_type_n_bytes(), size too big to support."));
 
+      MPI_Datatype result;
+
       const MPI_Datatype types[2] = {chunks, remainder};
       ierr =
         MPI_Type_create_struct(2, blocklengths, displacements, types, &result);
@@ -363,10 +365,12 @@ namespace Utilities
       Assert(size64 == static_cast<MPI_Count>(n_bytes), ExcInternalError());
 #      endif
 #    endif
-
+      return result;
 #  else
-      (void)result;
       (void)n_bytes;
+      Assert(false,
+             ExcMessage("This operation is not useful without MPI support."));
+      return MPI_Result(0);
 #  endif
     }
 
index f7564bda40488c34ee6591fdeaa730d235de96ed..48e14b426e5d71cdadbb95016b41ef3c91a1d27c 100644 (file)
@@ -27,10 +27,9 @@ const bool run_big = false;
 using namespace dealii;
 
 void
-test_data_type(long long n_bytes)
+test_data_type(const std::uint64_t n_bytes)
 {
-  MPI_Datatype bigtype;
-  Utilities::MPI::create_mpi_data_type_n_bytes(bigtype, n_bytes);
+  MPI_Datatype bigtype = Utilities::MPI::create_mpi_data_type_n_bytes(n_bytes);
 
   deallog << "checking size " << n_bytes << ":";
 
@@ -60,15 +59,15 @@ test_data_type(long long n_bytes)
 void
 test_send_recv(MPI_Comm comm)
 {
-  unsigned int    myid    = Utilities::MPI::this_mpi_process(comm);
-  const long long n_bytes = (run_big) ? ((1ULL << 31) + 37) : (37ULL);
+  unsigned int        myid    = Utilities::MPI::this_mpi_process(comm);
+  const std::uint64_t n_bytes = (run_big) ? ((1ULL << 31) + 37) : (37ULL);
 
   if (myid == 0)
     {
       std::vector<char> buffer(n_bytes, 'A');
       buffer[n_bytes - 1] = 'B';
-      MPI_Datatype bigtype;
-      Utilities::MPI::create_mpi_data_type_n_bytes(bigtype, buffer.size());
+      MPI_Datatype bigtype =
+        Utilities::MPI::create_mpi_data_type_n_bytes(buffer.size());
       int ierr =
         MPI_Send(buffer.data(), 1, bigtype, 1 /* dest */, 0 /* tag */, comm);
       AssertThrowMPI(ierr);
@@ -78,8 +77,8 @@ test_send_recv(MPI_Comm comm)
   else if (myid == 1)
     {
       std::vector<char> buffer(n_bytes, '?');
-      MPI_Datatype      bigtype;
-      Utilities::MPI::create_mpi_data_type_n_bytes(bigtype, buffer.size());
+      MPI_Datatype      bigtype =
+        Utilities::MPI::create_mpi_data_type_n_bytes(buffer.size());
       int ierr = MPI_Recv(buffer.data(),
                           1,
                           bigtype,

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.