]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make sure MPI BigData objects are always destroyed. 12987/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Wed, 29 Dec 2021 00:16:00 +0000 (17:16 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 29 Dec 2021 00:24:16 +0000 (17:24 -0700)
include/deal.II/base/mpi.h
source/base/mpi.cc
source/distributed/tria_base.cc
tests/mpi/create_mpi_datatype_01.cc

index ea98a5c31754492f78b98432a29d119b257ccad2..43177039579dfbc5e43b013e9b8973b7e3ba78ee 100644 (file)
@@ -500,35 +500,53 @@ namespace Utilities
 
 
     /**
-     * Create an MPI_Datatype that consists of @p n_bytes bytes.
+     * Create a object that contains an `MPI_Datatype` that represents @p n_bytes bytes.
      *
-     * The resulting data type can be used in MPI send/recv or MPI IO to process
-     * messages of sizes larger than 2 GB with MPI_Byte as the underlying data
-     * type. This helper is required for MPI versions before 4.0 because
-     * routines like MPI_Send
+     * The resulting data type can be used in MPI send/receive or MPI IO to
+     * process messages of sizes larger than 2 GB with MPI_Byte as the
+     * underlying data type. This helper is required for MPI versions before 4.0
+     * because routines like `MPI_Send`
      * use a signed interger for the @p count variable. Instead, you can use this
      * data type with the appropriate size set to the size of your message and
      * by passing
      * 1 as the @p count.
      *
-     * @note You need to free this data type after you are done using it by calling
-     * <code>MPI_Type_free(&result)</code>.
+     * @note The function does not just return an object of type `MPI_Datatype`
+     *   because such objects need to be destroyed by a call to `MPI_Type_free`
+     *   and it is easy to forget to do so (thereby creating a resource leak).
+     *   Rather, the function returns an object that *points* to such an
+     *   `MPI_Datatype` object, but also has a "deleter" function that ensures
+     *   that `MPI_Type_free` is called whenever the object returned by this
+     *   function goes out of scope.
      *
      * Usage example:
      * <code>
      * std::vector<char> buffer;
+     * [...]
      * if (buffer.size()<(1U<<31))
+     * {                               // less than 2GB of data
      *   MPI_Send(buffer.data(), buffer.size(), MPI_BYTE, dest, tag, comm);
+     * }
      * else
-     * {
-     *   MPI_Datatype bigtype =
+     * {                               // more than 2GB of data
+     *   const auto bigtype =
      *     Utilities::MPI::create_mpi_data_type_n_bytes(buffer.size());
-     *   MPI_Send(buffer.data(), 1, bigtype, dest, tag, comm);
-     *   MPI_Type_free(&bigtype);
+     *   MPI_Send(buffer.data(), 1, *bigtype, dest, tag, comm);
+     * }
+     * </code>
+     * Alternatively, the code in the `else` branch can be simplified to
+     * the following:
+     * <code>
+     * [...]
+     * else
+     * {                               // more than 2GB of data
+     *   MPI_Send(buffer.data(), 1,
+     *            *Utilities::MPI::create_mpi_data_type_n_bytes(buffer.size()),
+     *            dest, tag, comm);
      * }
      * </code>
      */
-    MPI_Datatype
+    std::unique_ptr<MPI_Datatype, void (*)(MPI_Datatype *)>
     create_mpi_data_type_n_bytes(const std::size_t n_bytes);
 
     /**
index 5cb26fc5ee5deaf5b170e1578fe4c92c4ce877c1..a6ba44f83ddcfe07f53dc993df2ec69009bb9773 100644 (file)
@@ -297,7 +297,7 @@ namespace Utilities
 
 
 
-    MPI_Datatype
+    std::unique_ptr<MPI_Datatype, void (*)(MPI_Datatype *)>
     create_mpi_data_type_n_bytes(const std::size_t n_bytes)
     {
       // Simplified version from BigMPI repository, see
@@ -337,7 +337,7 @@ namespace Utilities
       AssertThrow(
         displacements[1] == n_chunks * max_signed_int,
         ExcMessage(
-          "ERROR in create_mpi_data_type_n_bytes(), size too big to support."));
+          "Error in create_mpi_data_type_n_bytes(): the size is too big to support."));
 
       MPI_Datatype result;
 
@@ -357,14 +357,33 @@ namespace Utilities
 #  ifdef DEBUG
 #    if DEAL_II_MPI_VERSION_GTE(3, 0)
       MPI_Count size64;
-      // this function is only available starting MPI 3.0:
+      // this function is only available starting with MPI 3.0:
       ierr = MPI_Type_size_x(result, &size64);
       AssertThrowMPI(ierr);
 
       Assert(size64 == static_cast<MPI_Count>(n_bytes), ExcInternalError());
 #    endif
 #  endif
-      return result;
+
+      // Now put the new data type into a std::unique_ptr with a custom
+      // deleter. We call the std::unique_ptr constructor that as first
+      // argument takes a pointer (here, a pointer to a copy of the `result`
+      // object, and as second argument a pointer-to-function, for which
+      // we here use a lambda function without captures that acts as the
+      // 'deleter' object: it calls `MPI_Type_free` and then deletes the
+      // pointer.
+      return {// The copy of the object:
+              new MPI_Datatype(result),
+              // The deleter:
+              [](MPI_Datatype *p) {
+                if (p != nullptr)
+                  {
+                    const int ierr = MPI_Type_free(p);
+                    AssertThrowMPI(ierr);
+
+                    delete p;
+                  }
+              }};
     }
 
 
index 7c622f5cc02b61f1b1ab8cda65e026c8d6ce28fa..66b75f803dd5ddc25790ac22f03bd8d4409af745 100644 (file)
@@ -1466,24 +1466,23 @@ namespace parallel
       else
         {
           // Writes bigger than 2GB require some extra care:
-          MPI_Datatype bigtype =
-            Utilities::MPI::create_mpi_data_type_n_bytes(src_data_fixed.size());
           ierr =
             MPI_File_write_at(fh,
                               my_global_file_position,
                               DEAL_II_MPI_CONST_CAST(src_data_fixed.data()),
                               1,
-                              bigtype,
+                              *Utilities::MPI::create_mpi_data_type_n_bytes(
+                                src_data_fixed.size()),
                               MPI_STATUS_IGNORE);
           AssertThrowMPI(ierr);
-          ierr = MPI_Type_free(&bigtype);
-          AssertThrowMPI(ierr);
         }
 
       ierr = MPI_File_close(&fh);
       AssertThrowMPI(ierr);
     }
 
+
+
     //
     // ---------- Variable size data ----------
     //
@@ -1568,17 +1567,15 @@ namespace parallel
         else
           {
             // Writes bigger than 2GB require some extra care:
-            MPI_Datatype bigtype = Utilities::MPI::create_mpi_data_type_n_bytes(
-              src_data_variable.size());
-            ierr = MPI_File_write_at(fh,
-                                     my_global_file_position,
-                                     DEAL_II_MPI_CONST_CAST(
-                                       src_data_variable.data()),
-                                     1,
-                                     bigtype,
-                                     MPI_STATUS_IGNORE);
-            AssertThrowMPI(ierr);
-            ierr = MPI_Type_free(&bigtype);
+            ierr =
+              MPI_File_write_at(fh,
+                                my_global_file_position,
+                                DEAL_II_MPI_CONST_CAST(
+                                  src_data_variable.data()),
+                                1,
+                                *Utilities::MPI::create_mpi_data_type_n_bytes(
+                                  src_data_variable.size()),
+                                MPI_STATUS_IGNORE);
             AssertThrowMPI(ierr);
           }
 
@@ -1681,17 +1678,14 @@ namespace parallel
       else
         {
           // Reads bigger than 2GB require some extra care:
-          MPI_Datatype bigtype = Utilities::MPI::create_mpi_data_type_n_bytes(
-            dest_data_fixed.size());
           ierr = MPI_File_read_at(fh,
                                   my_global_file_position,
                                   dest_data_fixed.data(),
                                   1,
-                                  bigtype,
+                                  *Utilities::MPI::create_mpi_data_type_n_bytes(
+                                    dest_data_fixed.size()),
                                   MPI_STATUS_IGNORE);
           AssertThrowMPI(ierr);
-          ierr = MPI_Type_free(&bigtype);
-          AssertThrowMPI(ierr);
         }
 
       ierr = MPI_File_close(&fh);
@@ -1772,16 +1766,14 @@ namespace parallel
         else
           {
             // Reads bigger than 2GB require some extra care:
-            MPI_Datatype bigtype = Utilities::MPI::create_mpi_data_type_n_bytes(
-              src_data_fixed.size());
-            ierr = MPI_File_read_at(fh,
-                                    my_global_file_position,
-                                    dest_data_variable.data(),
-                                    1,
-                                    bigtype,
-                                    MPI_STATUS_IGNORE);
-            AssertThrowMPI(ierr);
-            ierr = MPI_Type_free(&bigtype);
+            ierr =
+              MPI_File_read_at(fh,
+                               my_global_file_position,
+                               dest_data_variable.data(),
+                               1,
+                               *Utilities::MPI::create_mpi_data_type_n_bytes(
+                                 src_data_fixed.size()),
+                               MPI_STATUS_IGNORE);
             AssertThrowMPI(ierr);
           }
 
index 48e14b426e5d71cdadbb95016b41ef3c91a1d27c..57c9340d1f38177e548962861a44f00d6de6242c 100644 (file)
@@ -29,12 +29,12 @@ using namespace dealii;
 void
 test_data_type(const std::uint64_t n_bytes)
 {
-  MPI_Datatype bigtype = Utilities::MPI::create_mpi_data_type_n_bytes(n_bytes);
+  const auto bigtype = Utilities::MPI::create_mpi_data_type_n_bytes(n_bytes);
 
   deallog << "checking size " << n_bytes << ":";
 
   int size32;
-  int ierr = MPI_Type_size(bigtype, &size32);
+  int ierr = MPI_Type_size(*bigtype, &size32);
   AssertThrowMPI(ierr);
 
   if (size32 == MPI_UNDEFINED)
@@ -45,17 +45,17 @@ test_data_type(const std::uint64_t n_bytes)
 
 #if DEAL_II_MPI_VERSION_GTE(3, 0)
   MPI_Count size64;
-  ierr = MPI_Type_size_x(bigtype, &size64);
+  ierr = MPI_Type_size_x(*bigtype, &size64);
   AssertThrowMPI(ierr);
 
   deallog << " size64=" << size64;
 #endif
 
   deallog << std::endl;
-
-  MPI_Type_free(&bigtype);
 }
 
+
+
 void
 test_send_recv(MPI_Comm comm)
 {
@@ -66,29 +66,25 @@ test_send_recv(MPI_Comm comm)
     {
       std::vector<char> buffer(n_bytes, 'A');
       buffer[n_bytes - 1] = 'B';
-      MPI_Datatype bigtype =
+      const auto 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);
-      ierr = MPI_Type_free(&bigtype);
+        MPI_Send(buffer.data(), 1, *bigtype, 1 /* dest */, 0 /* tag */, comm);
       AssertThrowMPI(ierr);
     }
   else if (myid == 1)
     {
       std::vector<char> buffer(n_bytes, '?');
-      MPI_Datatype      bigtype =
+      const auto        bigtype =
         Utilities::MPI::create_mpi_data_type_n_bytes(buffer.size());
       int ierr = MPI_Recv(buffer.data(),
                           1,
-                          bigtype,
+                          *bigtype,
                           0 /* src */,
                           0 /* tag */,
                           comm,
                           MPI_STATUS_IGNORE);
       AssertThrowMPI(ierr);
-      ierr = MPI_Type_free(&bigtype);
-      AssertThrowMPI(ierr);
 
       AssertThrow(buffer[0] == 'A', ExcInternalError());
       AssertThrow(buffer[n_bytes - 1] == 'B', ExcInternalError());

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.