From: Luca Heltai Date: Sat, 18 Nov 2017 20:49:54 +0000 (+0100) Subject: some_to_some and all_gather. X-Git-Tag: v9.0.0-rc1~748^2~5 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=edc7571aa2339e9cde96782c7546efc64d448fd0;p=dealii.git some_to_some and all_gather. --- diff --git a/doc/news/changes/minor/20171118GiovanniAlzetta-LucaHeltai b/doc/news/changes/minor/20171118GiovanniAlzetta-LucaHeltai new file mode 100644 index 0000000000..b54cfef10c --- /dev/null +++ b/doc/news/changes/minor/20171118GiovanniAlzetta-LucaHeltai @@ -0,0 +1,4 @@ +New: Utilities::MPI::all_gather and Utilities::MPI::some_to_some functions have been +added to perform general collective communications between processors. +
+(Giovanni Alzetta, Luca Heltai, 2017/11/18) diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index f931f23d63..f50191b62b 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -20,6 +20,7 @@ #include #include +#include #if !defined(DEAL_II_WITH_MPI) && !defined(DEAL_II_WITH_PETSC) // without MPI, we would still like to use @@ -395,8 +396,6 @@ namespace Utilities min_max_avg (const double my_value, const MPI_Comm &mpi_communicator); - - /** * A class that is used to initialize the MPI system at the beginning of a * program and to shut it down again at the end. It also allows you to @@ -491,6 +490,49 @@ namespace Utilities */ bool job_supports_mpi (); + /** + * Initiate a some-to-some communication, and exchange arbitrary objects + * (the class T should be serializable using boost::serialize) between + * processors. + * + * @param[in] comm MPI communicator. + * + * @param[in] objects_to_send A map from the rank (unsigned int) of the + * process meant to receive the data and the object to send (the type T + * must be serializable for this function to work properly). + * + * @param[out] A map from the rank (unsigned int) of the process + * which sent the data and object received. + * + * @author Giovanni Alzetta, Luca Heltai, 2017 + */ + template + std::map + some_to_some(const MPI_Comm &comm, + const std::map &objects_to_send); + + /** + * A generalization of the classic MPI_Allgather function, that accepts + * arbitrary data types T, as long as boost::serialize accepts T as an + * argument. + * + * @param[in] comm MPI communicator. + * @param[in] objects_to_send A map of `rank:object`, where the key `rank` + * is the rank of the processor we want to send to, and `object` is the + * object to send + * + * @param[out] 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. + * + * @author Giovanni Alzetta, Luca Heltai, 2017 + */ + template + std::vector + all_gather(const MPI_Comm &comm, + const T &object_to_send); + #ifndef DOXYGEN // declaration for an internal function that lives in mpi.templates.h namespace internal diff --git a/include/deal.II/base/mpi.templates.h b/include/deal.II/base/mpi.templates.h index 4ad8353035..4fab76fcb2 100644 --- a/include/deal.II/base/mpi.templates.h +++ b/include/deal.II/base/mpi.templates.h @@ -190,6 +190,142 @@ namespace Utilities + template + std::map + some_to_some(const MPI_Comm &comm, + const std::map &objects_to_send) + { +#ifndef DEAL_II_WITH_MPI + (void)comm; + (void)objects_to_send; + Assert (false, ExcMessage ("The function some_to_some doesn't make" + "any sense if you do not have MPI enabled. ")); +#else + auto n_procs = dealii::Utilities::MPI::n_mpi_processes(comm); + auto my_proc = dealii::Utilities::MPI::this_mpi_process(comm); + + std::vector send_to(objects_to_send.size()); + { + unsigned int i=0; + for (const auto &m: objects_to_send) + send_to[i++] = m.first; + } + AssertDimension(send_to.size(), objects_to_send.size()); + + const auto receive_from = + Utilities::MPI::compute_point_to_point_communication_pattern(comm, send_to); + + // Sending buffers + std::vector > buffers_to_send(send_to.size()); + std::vector buffer_send_requests(send_to.size()); + { + unsigned int i = 0; + for (const auto &rank_obj : objects_to_send) + { + const auto &rank = rank_obj.first; + buffers_to_send[i] = Utilities::pack(rank_obj.second); + const int ierr = MPI_Isend(buffers_to_send[i].data(), + buffers_to_send[i].size(), MPI_CHAR, + rank, 21, comm, &buffer_send_requests[i]); + AssertThrowMPI(ierr); + ++i; + } + } + + // Receiving buffers + std::map received_objects; + { + unsigned int i = 0; + std::vector buffer; + + // We do this on a first come/first served basis + while (i(buffer); + ++i; + } + } + + // Wait to have sent all objects. + MPI_Waitall(send_to.size(), buffer_send_requests.data(),MPI_STATUSES_IGNORE); + + return received_objects; +#endif // deal.II with MPI + } + + + + template + std::vector all_gather(const MPI_Comm &comm, + const T &object) + { +#ifndef DEAL_II_WITH_MPI + (void)comm; + (void)objects_to_send; + Assert (false, ExcMessage ("The function all_gather doesn't make" + "any sense if you do not have MPI enabled. ")); +#else + auto n_procs = dealii::Utilities::MPI::n_mpi_processes(comm); + + std::vector buffer = Utilities::pack(object); + + int n_local_data = buffer.size(); + + // Vector to store the size of loc_data_array for every process + std::vector size_all_data(n_procs,0); + + // Exchanging the size of each buffer + MPI_Allgather(&n_local_data, 1, MPI_INT, + &(size_all_data[0]), 1, MPI_INT, + comm); + + // Now computing the the displacement, relative to recvbuf, + // at which to store the incoming buffer + std::vector rdispls(n_procs); + rdispls[0] = 0; + for (unsigned int i=1; i < n_procs; ++i) + rdispls[i] = rdispls[i-1] + size_all_data[i-1]; + + // Step 3: exchange the buffer: + std::vector received_unrolled_buffer(rdispls.back() + size_all_data.back()); + + MPI_Allgatherv(buffer.data(), n_local_data, MPI_CHAR, + received_unrolled_buffer.data(), size_all_data.data(), + rdispls.data(), MPI_CHAR, comm); + + std::vector received_objects(n_procs); + for (unsigned int i= 0; i < n_procs; ++i) + { + std::vector local_buffer(received_unrolled_buffer.begin()+rdispls[i], + received_unrolled_buffer.begin()+rdispls[i]+size_all_data[i]); + received_objects[i] = Utilities::unpack(local_buffer); + } + + return received_objects; +#endif + } + + + template T sum (const T &t, const MPI_Comm &mpi_communicator) diff --git a/tests/base/mpi_send_and_receive_01.cc b/tests/base/mpi_send_and_receive_01.cc new file mode 100644 index 0000000000..3e8cff82ea --- /dev/null +++ b/tests/base/mpi_send_and_receive_01.cc @@ -0,0 +1,64 @@ +// --------------------------------------------------------------------- +// +// 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 the internal::send_and_receive function in the case of +// collective communication + +#include "../tests.h" +#include "../../include/deal.II/base/mpi.templates.h" + +#include +#include + +#include +#include + +int main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + initlog(); + + 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 > local_points(my_proc + 1); + for (unsigned int i=0; i(my_proc, -my_proc, i); + + auto gathered_points = Utilities::MPI::all_gather(MPI_COMM_WORLD, local_points); + 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_send_and_receive_01.mpirun=2.output b/tests/base/mpi_send_and_receive_01.mpirun=2.output new file mode 100644 index 0000000000..186b780267 --- /dev/null +++ b/tests/base/mpi_send_and_receive_01.mpirun=2.output @@ -0,0 +1,2 @@ + +DEAL::Test: ok diff --git a/tests/base/mpi_send_and_receive_01.mpirun=6.output b/tests/base/mpi_send_and_receive_01.mpirun=6.output new file mode 100644 index 0000000000..186b780267 --- /dev/null +++ b/tests/base/mpi_send_and_receive_01.mpirun=6.output @@ -0,0 +1,2 @@ + +DEAL::Test: ok diff --git a/tests/base/mpi_send_and_receive_02.cc b/tests/base/mpi_send_and_receive_02.cc new file mode 100644 index 0000000000..602c5f2248 --- /dev/null +++ b/tests/base/mpi_send_and_receive_02.cc @@ -0,0 +1,58 @@ +// --------------------------------------------------------------------- +// +// 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 the internal::send_and_receive function in the case of +// collective communication + +#include "../tests.h" +#include "../../include/deal.II/base/mpi.templates.h" + +#include +#include + +#include +#include + + + +int main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + initlog(); + + auto n_procs = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + auto my_proc = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + // Creating the map array of points to be sent + std::map > > m; + for (unsigned int i=0; i(my_proc, -2.0*my_proc)); + + auto received_pts = Utilities::MPI::some_to_some(MPI_COMM_WORLD, m); + + bool test_passed = true; + for (auto const &pt: received_pts) + if ( std::abs(pt.first - pt.second[0][0]) > 1e-12 || + std::abs(2.0*pt.first + pt.second[0][1]) > 1e-12 ) + { + test_passed = false; + deallog << "Error with point " << pt.second[0] << " received from rank " << pt.first << std::endl; + } + if (test_passed) + deallog << "Test: ok" << std::endl; + else + deallog << "Test: FAILED" << std::endl; +} diff --git a/tests/base/mpi_send_and_receive_02.mpirun=2.output b/tests/base/mpi_send_and_receive_02.mpirun=2.output new file mode 100644 index 0000000000..186b780267 --- /dev/null +++ b/tests/base/mpi_send_and_receive_02.mpirun=2.output @@ -0,0 +1,2 @@ + +DEAL::Test: ok diff --git a/tests/base/mpi_send_and_receive_02.mpirun=6.output b/tests/base/mpi_send_and_receive_02.mpirun=6.output new file mode 100644 index 0000000000..186b780267 --- /dev/null +++ b/tests/base/mpi_send_and_receive_02.mpirun=6.output @@ -0,0 +1,2 @@ + +DEAL::Test: ok