]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix a problem with combining MPI::LargeCount and PETSc. 13885/head
authorDavid Wells <drwells@email.unc.edu>
Wed, 1 Jun 2022 15:58:53 +0000 (11:58 -0400)
committerDavid Wells <drwells@email.unc.edu>
Wed, 1 Jun 2022 17:13:34 +0000 (13:13 -0400)
PETSc, in petsclog.h, implements logging by redefining all MPI functions as
macros, e.g.,

    #define MPI_Send_c(buf,count,datatype,dest,tag,comm) \
    ((petsc_send_ct++,0) || PetscMPITypeSize((count),(datatype),(&petsc_send_len)) || MPI_Send_c((buf),(count),(datatype),(dest),(tag),(comm)))

This only works when all MPI declarations are made before these macro
definitions. Hence this library won't work with PETSc unless we include headers
in the right order.

Work around this by removing the MPI_ prefix from these function definitions.
Since we put these functions in our own namespace there is no need for them to
have the same names as the normal MPI functions.

include/deal.II/base/mpi_large_count.h
source/base/data_out_base.cc
source/base/mpi.cc
source/distributed/tria_base.cc

index 8577f366f05625c2c047e0a50043c3104c9773e1..8a6c4a90138416b9a8896e80561af35952b0a872 100644 (file)
@@ -64,12 +64,12 @@ namespace Utilities
        * See the MPI 4.x standard for details.
        */
       inline int
-      MPI_Type_contiguous_c(MPI_Count     count,
-                            MPI_Datatype  oldtype,
-                            MPI_Datatype *newtype)
+      Type_contiguous_c(MPI_Count     count,
+                        MPI_Datatype  oldtype,
+                        MPI_Datatype *newtype)
       {
 #  if MPI_VERSION >= 4
-        return ::MPI_Type_contiguous_c(count, oldtype, newtype);
+        return MPI_Type_contiguous_c(count, oldtype, newtype);
 #  else
         if (count <= LargeCount::mpi_max_int_count)
           return MPI_Type_contiguous(count, oldtype, newtype);
@@ -151,22 +151,22 @@ namespace Utilities
        * See the MPI 4.x standard for details.
        */
       inline int
-      MPI_Send_c(const void * buf,
-                 MPI_Count    count,
-                 MPI_Datatype datatype,
-                 int          dest,
-                 int          tag,
-                 MPI_Comm     comm)
+      Send_c(const void * buf,
+             MPI_Count    count,
+             MPI_Datatype datatype,
+             int          dest,
+             int          tag,
+             MPI_Comm     comm)
       {
 #  if MPI_VERSION >= 4
-        return ::MPI_Send_c(buf, count, datatype, dest, tag, comm);
+        return MPI_Send_c(buf, count, datatype, dest, tag, comm);
 #  else
         if (count <= LargeCount::mpi_max_int_count)
           return MPI_Send(buf, count, datatype, dest, tag, comm);
 
         MPI_Datatype bigtype;
         int          ierr;
-        ierr = MPI_Type_contiguous_c(count, datatype, &bigtype);
+        ierr = Type_contiguous_c(count, datatype, &bigtype);
         if (ierr != MPI_SUCCESS)
           return ierr;
         ierr = MPI_Type_commit(&bigtype);
@@ -190,23 +190,23 @@ namespace Utilities
        * See the MPI 4.x standard for details.
        */
       inline int
-      MPI_Recv_c(void *       buf,
-                 MPI_Count    count,
-                 MPI_Datatype datatype,
-                 int          source,
-                 int          tag,
-                 MPI_Comm     comm,
-                 MPI_Status * status)
+      Recv_c(void *       buf,
+             MPI_Count    count,
+             MPI_Datatype datatype,
+             int          source,
+             int          tag,
+             MPI_Comm     comm,
+             MPI_Status * status)
       {
 #  if MPI_VERSION >= 4
-        return ::MPI_Recv_c(buf, count, datatype, source, tag, comm, status);
+        return MPI_Recv_c(buf, count, datatype, source, tag, comm, status);
 #  else
         if (count <= LargeCount::mpi_max_int_count)
           return MPI_Recv(buf, count, datatype, source, tag, comm, status);
 
         MPI_Datatype bigtype;
         int          ierr;
-        ierr = MPI_Type_contiguous_c(count, datatype, &bigtype);
+        ierr = Type_contiguous_c(count, datatype, &bigtype);
         if (ierr != MPI_SUCCESS)
           return ierr;
 
@@ -232,21 +232,21 @@ namespace Utilities
        * See the MPI 4.x standard for details.
        */
       inline int
-      MPI_Bcast_c(void *       buf,
-                  MPI_Count    count,
-                  MPI_Datatype datatype,
-                  unsigned int root_mpi_rank,
-                  MPI_Comm     comm)
+      Bcast_c(void *       buf,
+              MPI_Count    count,
+              MPI_Datatype datatype,
+              unsigned int root_mpi_rank,
+              MPI_Comm     comm)
       {
 #  if MPI_VERSION >= 4
-        return ::MPI_Bcast_c(buf, count, datatype, root_mpi_rank, comm);
+        return MPI_Bcast_c(buf, count, datatype, root_mpi_rank, comm);
 #  else
         if (count <= LargeCount::mpi_max_int_count)
           return MPI_Bcast(buf, count, datatype, root_mpi_rank, comm);
 
         MPI_Datatype bigtype;
         int          ierr;
-        ierr = MPI_Type_contiguous_c(count, datatype, &bigtype);
+        ierr = Type_contiguous_c(count, datatype, &bigtype);
         if (ierr != MPI_SUCCESS)
           return ierr;
         ierr = MPI_Type_commit(&bigtype);
@@ -270,19 +270,19 @@ namespace Utilities
        * See the MPI 4.x standard for details.
        */
       inline int
-      MPI_File_write_at_c(MPI_File     fh,
-                          MPI_Offset   offset,
-                          const void * buf,
-                          MPI_Count    count,
-                          MPI_Datatype datatype,
-                          MPI_Status * status)
+      File_write_at_c(MPI_File     fh,
+                      MPI_Offset   offset,
+                      const void * buf,
+                      MPI_Count    count,
+                      MPI_Datatype datatype,
+                      MPI_Status * status)
       {
         if (count <= LargeCount::mpi_max_int_count)
           return MPI_File_write_at(fh, offset, buf, count, datatype, status);
 
         MPI_Datatype bigtype;
         int          ierr;
-        ierr = MPI_Type_contiguous_c(count, datatype, &bigtype);
+        ierr = Type_contiguous_c(count, datatype, &bigtype);
         if (ierr != MPI_SUCCESS)
           return ierr;
         ierr = MPI_Type_commit(&bigtype);
@@ -306,12 +306,12 @@ namespace Utilities
        * See the MPI 4.x standard for details.
        */
       inline int
-      MPI_File_write_at_all_c(MPI_File     fh,
-                              MPI_Offset   offset,
-                              const void * buf,
-                              MPI_Count    count,
-                              MPI_Datatype datatype,
-                              MPI_Status * status)
+      File_write_at_all_c(MPI_File     fh,
+                          MPI_Offset   offset,
+                          const void * buf,
+                          MPI_Count    count,
+                          MPI_Datatype datatype,
+                          MPI_Status * status)
       {
         if (count <= LargeCount::mpi_max_int_count)
           return MPI_File_write_at_all(
@@ -319,7 +319,7 @@ namespace Utilities
 
         MPI_Datatype bigtype;
         int          ierr;
-        ierr = MPI_Type_contiguous_c(count, datatype, &bigtype);
+        ierr = Type_contiguous_c(count, datatype, &bigtype);
         if (ierr != MPI_SUCCESS)
           return ierr;
         ierr = MPI_Type_commit(&bigtype);
@@ -342,18 +342,18 @@ namespace Utilities
        * See the MPI 4.x standard for details.
        */
       inline int
-      MPI_File_write_ordered_c(MPI_File     fh,
-                               const void * buf,
-                               MPI_Count    count,
-                               MPI_Datatype datatype,
-                               MPI_Status * status)
+      File_write_ordered_c(MPI_File     fh,
+                           const void * buf,
+                           MPI_Count    count,
+                           MPI_Datatype datatype,
+                           MPI_Status * status)
       {
         if (count <= LargeCount::mpi_max_int_count)
           return MPI_File_write_ordered(fh, buf, count, datatype, status);
 
         MPI_Datatype bigtype;
         int          ierr;
-        ierr = MPI_Type_contiguous_c(count, datatype, &bigtype);
+        ierr = Type_contiguous_c(count, datatype, &bigtype);
         if (ierr != MPI_SUCCESS)
           return ierr;
         ierr = MPI_Type_commit(&bigtype);
@@ -377,19 +377,19 @@ namespace Utilities
        * See the MPI 4.x standard for details.
        */
       inline int
-      MPI_File_read_at_c(MPI_File     fh,
-                         MPI_Offset   offset,
-                         void *       buf,
-                         MPI_Count    count,
-                         MPI_Datatype datatype,
-                         MPI_Status * status)
+      File_read_at_c(MPI_File     fh,
+                     MPI_Offset   offset,
+                     void *       buf,
+                     MPI_Count    count,
+                     MPI_Datatype datatype,
+                     MPI_Status * status)
       {
         if (count <= LargeCount::mpi_max_int_count)
           return MPI_File_read_at(fh, offset, buf, count, datatype, status);
 
         MPI_Datatype bigtype;
         int          ierr;
-        ierr = MPI_Type_contiguous_c(count, datatype, &bigtype);
+        ierr = Type_contiguous_c(count, datatype, &bigtype);
         if (ierr != MPI_SUCCESS)
           return ierr;
         ierr = MPI_Type_commit(&bigtype);
@@ -413,19 +413,19 @@ namespace Utilities
        * See the MPI 4.x standard for details.
        */
       inline int
-      MPI_File_read_at_all_c(MPI_File     fh,
-                             MPI_Offset   offset,
-                             void *       buf,
-                             MPI_Count    count,
-                             MPI_Datatype datatype,
-                             MPI_Status * status)
+      File_read_at_all_c(MPI_File     fh,
+                         MPI_Offset   offset,
+                         void *       buf,
+                         MPI_Count    count,
+                         MPI_Datatype datatype,
+                         MPI_Status * status)
       {
         if (count <= LargeCount::mpi_max_int_count)
           return MPI_File_read_at_all(fh, offset, buf, count, datatype, status);
 
         MPI_Datatype bigtype;
         int          ierr;
-        ierr = MPI_Type_contiguous_c(count, datatype, &bigtype);
+        ierr = Type_contiguous_c(count, datatype, &bigtype);
         if (ierr != MPI_SUCCESS)
           return ierr;
         ierr = MPI_Type_commit(&bigtype);
index e1712365b30bb1247c068f54120c254567ab1b59..42484a5ef632ffb3fdf1198a13b456aa6273be67 100644 (file)
@@ -7733,7 +7733,7 @@ DataOutInterface<dim, spacedim>::write_vtu_in_parallel(
       DataOutBase::write_vtu_header(ss, vtk_flags);
       header_size = ss.str().size();
       // Write the header on rank 0 at the start of a file, i.e., offset 0.
-      ierr = Utilities::MPI::LargeCount::MPI_File_write_at_c(
+      ierr = Utilities::MPI::LargeCount::File_write_at_c(
         fh, 0, ss.str().c_str(), header_size, MPI_CHAR, MPI_STATUS_IGNORE);
       AssertThrowMPI(ierr);
     }
@@ -7769,13 +7769,12 @@ DataOutInterface<dim, spacedim>::write_vtu_in_parallel(
     // Locate specific offset for each processor.
     const MPI_Offset offset = static_cast<MPI_Offset>(header_size) + prefix_sum;
 
-    ierr =
-      Utilities::MPI::LargeCount::MPI_File_write_at_all_c(fh,
-                                                          offset,
-                                                          ss.str().c_str(),
-                                                          ss.str().size(),
-                                                          MPI_CHAR,
-                                                          MPI_STATUS_IGNORE);
+    ierr = Utilities::MPI::LargeCount::File_write_at_all_c(fh,
+                                                           offset,
+                                                           ss.str().c_str(),
+                                                           ss.str().size(),
+                                                           MPI_CHAR,
+                                                           MPI_STATUS_IGNORE);
     AssertThrowMPI(ierr);
 
     if (myrank == n_ranks - 1)
@@ -7788,13 +7787,12 @@ DataOutInterface<dim, spacedim>::write_vtu_in_parallel(
         const unsigned int footer_size = ss.str().size();
 
         // Writing footer:
-        ierr =
-          Utilities::MPI::LargeCount::MPI_File_write_at_c(fh,
-                                                          footer_offset,
-                                                          ss.str().c_str(),
-                                                          footer_size,
-                                                          MPI_CHAR,
-                                                          MPI_STATUS_IGNORE);
+        ierr = Utilities::MPI::LargeCount::File_write_at_c(fh,
+                                                           footer_offset,
+                                                           ss.str().c_str(),
+                                                           footer_size,
+                                                           MPI_CHAR,
+                                                           MPI_STATUS_IGNORE);
         AssertThrowMPI(ierr);
       }
   }
index 6973cacf6baf370ee06a24452fe4f8fcd546e6a6..6678a8cc1059b8a30bc375698066cb35c4cf2600 100644 (file)
@@ -255,7 +255,7 @@ namespace Utilities
     create_mpi_data_type_n_bytes(const std::size_t n_bytes)
     {
       MPI_Datatype result;
-      int ierr = LargeCount::MPI_Type_contiguous_c(n_bytes, MPI_BYTE, &result);
+      int ierr = LargeCount::Type_contiguous_c(n_bytes, MPI_BYTE, &result);
       AssertThrowMPI(ierr);
       ierr = MPI_Type_commit(&result);
       AssertThrowMPI(ierr);
index 4f43fa4057585d2c1228fd0f32d7aec78921917b..5975375cccf96256451abd1f2364bd3365d530a7 100644 (file)
@@ -1443,7 +1443,7 @@ namespace parallel
       // this task.
       if (myrank == 0)
         {
-          ierr = Utilities::MPI::LargeCount::MPI_File_write_at_c(
+          ierr = Utilities::MPI::LargeCount::File_write_at_c(
             fh,
             0,
             sizes_fixed_cumulative.data(),
@@ -1464,12 +1464,12 @@ namespace parallel
         static_cast<MPI_Offset>(global_first_cell) * bytes_per_cell;
 
       ierr =
-        Utilities::MPI::LargeCount::MPI_File_write_at_c(fh,
-                                                        my_global_file_position,
-                                                        src_data_fixed.data(),
-                                                        src_data_fixed.size(),
-                                                        MPI_BYTE,
-                                                        MPI_STATUS_IGNORE);
+        Utilities::MPI::LargeCount::File_write_at_c(fh,
+                                                    my_global_file_position,
+                                                    src_data_fixed.data(),
+                                                    src_data_fixed.size(),
+                                                    MPI_BYTE,
+                                                    MPI_STATUS_IGNORE);
       AssertThrowMPI(ierr);
 
       ierr = MPI_File_close(&fh);
@@ -1519,7 +1519,7 @@ namespace parallel
                           std::numeric_limits<int>::max()),
                       ExcNotImplemented());
 
-          ierr = Utilities::MPI::LargeCount::MPI_File_write_at_c(
+          ierr = Utilities::MPI::LargeCount::File_write_at_c(
             fh,
             my_global_file_position,
             src_sizes_variable.data(),
@@ -1547,13 +1547,13 @@ namespace parallel
           prefix_sum;
 
         // Write data consecutively into file.
-        ierr = Utilities::MPI::LargeCount::MPI_File_write_at_c(
-          fh,
-          my_global_file_position,
-          src_data_variable.data(),
-          src_data_variable.size(),
-          MPI_BYTE,
-          MPI_STATUS_IGNORE);
+        ierr =
+          Utilities::MPI::LargeCount::File_write_at_c(fh,
+                                                      my_global_file_position,
+                                                      src_data_variable.data(),
+                                                      src_data_variable.size(),
+                                                      MPI_BYTE,
+                                                      MPI_STATUS_IGNORE);
         AssertThrowMPI(ierr);
 
 
@@ -1615,7 +1615,7 @@ namespace parallel
       // location in the file.
       sizes_fixed_cumulative.resize(1 + n_attached_deserialize_fixed +
                                     (variable_size_data_stored ? 1 : 0));
-      ierr = Utilities::MPI::LargeCount::MPI_File_read_at_c(
+      ierr = Utilities::MPI::LargeCount::File_read_at_c(
         fh,
         0,
         sizes_fixed_cumulative.data(),
@@ -1639,13 +1639,12 @@ namespace parallel
         size_header +
         static_cast<MPI_Offset>(global_first_cell) * bytes_per_cell;
 
-      ierr =
-        Utilities::MPI::LargeCount::MPI_File_read_at_c(fh,
-                                                       my_global_file_position,
-                                                       dest_data_fixed.data(),
-                                                       dest_data_fixed.size(),
-                                                       MPI_BYTE,
-                                                       MPI_STATUS_IGNORE);
+      ierr = Utilities::MPI::LargeCount::File_read_at_c(fh,
+                                                        my_global_file_position,
+                                                        dest_data_fixed.data(),
+                                                        dest_data_fixed.size(),
+                                                        MPI_BYTE,
+                                                        MPI_STATUS_IGNORE);
       AssertThrowMPI(ierr);
 
 
@@ -1679,7 +1678,7 @@ namespace parallel
         const MPI_Offset my_global_file_position_sizes =
           static_cast<MPI_Offset>(global_first_cell) * sizeof(unsigned int);
 
-        ierr = Utilities::MPI::LargeCount::MPI_File_read_at_c(
+        ierr = Utilities::MPI::LargeCount::File_read_at_c(
           fh,
           my_global_file_position_sizes,
           dest_sizes_variable.data(),
@@ -1711,13 +1710,13 @@ namespace parallel
 
         dest_data_variable.resize(size_on_proc);
 
-        ierr = Utilities::MPI::LargeCount::MPI_File_read_at_c(
-          fh,
-          my_global_file_position,
-          dest_data_variable.data(),
-          dest_data_variable.size(),
-          MPI_BYTE,
-          MPI_STATUS_IGNORE);
+        ierr =
+          Utilities::MPI::LargeCount::File_read_at_c(fh,
+                                                     my_global_file_position,
+                                                     dest_data_variable.data(),
+                                                     dest_data_variable.size(),
+                                                     MPI_BYTE,
+                                                     MPI_STATUS_IGNORE);
         AssertThrowMPI(ierr);
 
         ierr = MPI_File_close(&fh);

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.