]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add gather to Utilities::MPI 5784/head
authorBenjamin Brands <benjamin.brands@fau.de>
Tue, 23 Jan 2018 08:47:26 +0000 (09:47 +0100)
committerBenjamin Brands <benjamin.brands@fau.de>
Thu, 25 Jan 2018 21:23:02 +0000 (22:23 +0100)
index_set: replace MPI_Gather by Utilities::MPI::gather

doc/news/changes/minor/20180125BenjaminBrands [new file with mode: 0644]
include/deal.II/base/mpi.h
source/base/index_set.cc
tests/base/mpi_gather_01.cc [new file with mode: 0644]
tests/base/mpi_gather_01.mpirun=2.output [new file with mode: 0644]
tests/base/mpi_gather_01.mpirun=6.output [new file with mode: 0644]
tests/base/mpi_gather_01.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20180125BenjaminBrands b/doc/news/changes/minor/20180125BenjaminBrands
new file mode 100644 (file)
index 0000000..d8a151c
--- /dev/null
@@ -0,0 +1,3 @@
+New: Added Utilities::MPI::gather wrapper and tests to gather objects from all to one MPI process. 
+<br>
+(Benjamin Brands, 2018/01/25)
index 2902516e8fd8d77cb3a35b6c999d97ac188477b6..1851cd1c285fb2780400c24a611b5cb066e95b73 100644 (file)
@@ -534,6 +534,28 @@ namespace Utilities
     all_gather(const MPI_Comm &comm,
                const T        &object_to_send);
 
+    /**
+     * A generalization of the classic MPI_Gather function, that accepts
+     * arbitrary data types T, as long as boost::serialize accepts T as an
+     * argument.
+     *
+     * @param[in] comm MPI communicator.
+     * @param[in] object_to_send an object to send to the root process
+     * @param[in] root_process The process, which receives the objects from all processes.
+     *  By default the process with rank 0 is the root process.
+     *
+     * @return The @p root_process receives a vector of objects, with size equal to the number of
+     *  processes in the MPI communicator. Each entry contains the object
+     *  received from the processor with the corresponding rank within the
+     *  communicator. All other processes receive an empty vector.
+     *
+     * @author Benjamin Brands, 2017
+     */
+    template <typename T>
+    std::vector<T>
+    gather(const MPI_Comm &comm,
+           const T &object_to_send,
+           const unsigned int root_process=0);
 
 #ifndef DOXYGEN
     // declaration for an internal function that lives in mpi.templates.h
@@ -573,6 +595,76 @@ namespace Utilities
       internal::all_reduce(MPI_MIN, ArrayView<const T>(values, N),
                            mpi_communicator, ArrayView<T>(minima, N));
     }
+
+    template <typename T>
+    std::vector<T>
+    gather(const MPI_Comm &comm,
+           const T &object_to_send,
+           const unsigned int root_process)
+    {
+#ifndef DEAL_II_WITH_MPI
+      (void)comm;
+      (void)root_process;
+      std::vector<T> v(1, object_to_send);
+      return v;
+#else
+      const auto n_procs = dealii::Utilities::MPI::n_mpi_processes(comm);
+      const auto my_rank = dealii::Utilities::MPI::this_mpi_process(comm);
+
+      Assert(root_process < n_procs, ExcIndexRange(root_process,0,n_procs));
+
+      std::vector<char> buffer = Utilities::pack(object_to_send);
+      int n_local_data = buffer.size();
+
+      // Vector to store the size of loc_data_array for every process
+      // only the root process needs to allocate memory for that purpose
+      std::vector<int> size_all_data;
+      if (my_rank==root_process)
+        size_all_data.resize(n_procs,0);
+
+      // Exchanging the size of each buffer
+      int ierr = MPI_Gather(&n_local_data, 1, MPI_INT,
+                            size_all_data.data(), 1, MPI_INT,
+                            root_process, comm);
+      AssertThrowMPI(ierr);
+
+      // Now computing the displacement, relative to recvbuf,
+      // at which to store the incoming buffer; only for root
+      std::vector<int> rdispls;
+      if (my_rank==root_process)
+        {
+          rdispls.resize(n_procs,0);
+          for (unsigned int i=1; i<n_procs; ++i)
+            rdispls[i] = rdispls[i-1] + size_all_data[i-1];
+        }
+      // exchange the buffer:
+      std::vector<char> received_unrolled_buffer;
+      if (my_rank==root_process)
+        received_unrolled_buffer.resize(rdispls.back() + size_all_data.back());
+
+      ierr = MPI_Gatherv(buffer.data(), n_local_data, MPI_CHAR,
+                         received_unrolled_buffer.data(), size_all_data.data(),
+                         rdispls.data(), MPI_CHAR,
+                         root_process, comm);
+      AssertThrowMPI(ierr);
+
+      std::vector<T> received_objects;
+
+      if (my_rank==root_process)
+        {
+          received_objects.resize(n_procs);
+
+          for (unsigned int i=0; i<n_procs; ++i)
+            {
+              const std::vector<char> local_buffer(received_unrolled_buffer.begin()+rdispls[i],
+                                                   received_unrolled_buffer.begin()+rdispls[i]+size_all_data[i]);
+              received_objects[i] = Utilities::unpack<T>(local_buffer);
+            }
+        }
+      return received_objects;
+#endif
+    }
+
 #endif
   } // end of namespace MPI
 } // end of namespace Utilities
index 44d587a467da4b70c47d8164a14ea2caa5ede72a..adb7b16a44b45ff45b285c82d81be3e94c9aead1 100644 (file)
@@ -637,15 +637,7 @@ IndexSet::is_ascending_and_one_to_one (const MPI_Comm &communicator) const
       : numbers::invalid_dof_index ;
 
   const unsigned int my_rank = Utilities::MPI::this_mpi_process(communicator);
-  const unsigned int n_ranks = Utilities::MPI::n_mpi_processes(communicator);
-  // first gather all information on process 0
-  const unsigned int gather_size = (my_rank==0)?n_ranks:1;
-  std::vector<types::global_dof_index> global_dofs(gather_size);
-
-  int ierr = MPI_Gather(&first_local_dof, 1, DEAL_II_DOF_INDEX_MPI_TYPE,
-                        &(global_dofs[0]), 1, DEAL_II_DOF_INDEX_MPI_TYPE, 0,
-                        communicator);
-  AssertThrowMPI(ierr);
+  const std::vector<types::global_dof_index> global_dofs = Utilities::MPI::gather(communicator,first_local_dof,0);
 
   if (my_rank == 0)
     {
@@ -672,7 +664,7 @@ IndexSet::is_ascending_and_one_to_one (const MPI_Comm &communicator) const
 
   // now broadcast the result
   int is_ascending = is_globally_ascending ? 1 : 0;
-  ierr = MPI_Bcast(&is_ascending, 1, MPI_INT, 0, communicator);
+  int ierr = MPI_Bcast(&is_ascending, 1, MPI_INT, 0, communicator);
   AssertThrowMPI(ierr);
 
   return (is_ascending==1);
diff --git a/tests/base/mpi_gather_01.cc b/tests/base/mpi_gather_01.cc
new file mode 100644 (file)
index 0000000..f347068
--- /dev/null
@@ -0,0 +1,63 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// Test for Utilities::MPI::gather
+
+#include "../tests.h"
+#include <deal.II/base/point.h>
+#include <deal.II/base/mpi.h>
+#include <vector>
+
+int main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  initlog();
+
+  unsigned int root_process=0;
+  auto n_procs = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
+  auto my_proc = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+  // Creating the local array of points
+  std::vector<Point<3> > local_points(my_proc + 1);
+  for (unsigned int i=0; i<my_proc+1; ++i)
+    local_points[i] = Point<3>(my_proc, -my_proc, i);
+
+  auto gathered_points = Utilities::MPI::gather(MPI_COMM_WORLD, local_points, root_process);
+
+  if (my_proc==root_process)
+    {
+      bool test_passed = true;
+      for (unsigned int i=0; i<n_procs; ++i)
+        {
+          if (gathered_points[i].size() != i+1)
+            {
+              test_passed = false;
+              deallog << "Error: Points received from rank " << i << " have wrong size. " << std::endl;
+            }
+          for (unsigned int p=0; p< gathered_points[i].size(); ++p)
+            if ( gathered_points[i][p][0] != (double) i  ||
+                 gathered_points[i][p][1] != (double) -i ||
+                 gathered_points[i][p][2] != (double) p)
+              {
+                test_passed = false;
+                deallog << "Error with point " << p << " from rank " << i << std::endl;
+              }
+        }
+      if (test_passed)
+        deallog << "Test: ok" << std::endl;
+      else
+        deallog << "Test: FAILED" << std::endl;
+    }
+}
diff --git a/tests/base/mpi_gather_01.mpirun=2.output b/tests/base/mpi_gather_01.mpirun=2.output
new file mode 100644 (file)
index 0000000..186b780
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::Test: ok
diff --git a/tests/base/mpi_gather_01.mpirun=6.output b/tests/base/mpi_gather_01.mpirun=6.output
new file mode 100644 (file)
index 0000000..186b780
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::Test: ok
diff --git a/tests/base/mpi_gather_01.output b/tests/base/mpi_gather_01.output
new file mode 100644 (file)
index 0000000..186b780
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::Test: ok

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.