From 69960fa700aebb2f73183a56b55db8d856cb7b20 Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Mon, 2 Dec 2019 11:51:20 +0100 Subject: [PATCH] Fix some to some. --- doc/news/changes/minor/20191202LucaHeltai | 3 + include/deal.II/base/mpi.h | 50 ++++++++------ tests/base/mpi_some_to_some_03.cc | 65 +++++++++++++++++++ .../base/mpi_some_to_some_03.mpirun=1.output | 5 ++ .../base/mpi_some_to_some_03.mpirun=2.output | 11 ++++ .../base/mpi_some_to_some_03.mpirun=3.output | 17 +++++ .../base/mpi_some_to_some_03.mpirun=4.output | 23 +++++++ 7 files changed, 153 insertions(+), 21 deletions(-) create mode 100644 doc/news/changes/minor/20191202LucaHeltai create mode 100644 tests/base/mpi_some_to_some_03.cc create mode 100644 tests/base/mpi_some_to_some_03.mpirun=1.output create mode 100644 tests/base/mpi_some_to_some_03.mpirun=2.output create mode 100644 tests/base/mpi_some_to_some_03.mpirun=3.output create mode 100644 tests/base/mpi_some_to_some_03.mpirun=4.output diff --git a/doc/news/changes/minor/20191202LucaHeltai b/doc/news/changes/minor/20191202LucaHeltai new file mode 100644 index 0000000000..cabb3d0e4c --- /dev/null +++ b/doc/news/changes/minor/20191202LucaHeltai @@ -0,0 +1,3 @@ +Fixed: Utilities::MPI::some_to_some() now allows sending to our own mpi process. +
+(Luca Heltai, 2019/12/02) diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index f52898371c..3845b752c2 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -1481,21 +1481,24 @@ namespace Utilities { # ifndef DEAL_II_WITH_MPI (void)comm; - Assert(objects_to_send.size() == 0, + Assert(objects_to_send.size() < 2, ExcMessage("Cannot send to more than one processor.")); Assert(objects_to_send.find(0) != objects_to_send.end() || objects_to_send.size() == 0, ExcMessage("Can only send to myself or to nobody.")); return objects_to_send; # else + const auto my_proc = 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()); + std::vector send_to; + send_to.reserve(objects_to_send.size()); + for (const auto &m : objects_to_send) + if (m.first != my_proc) + send_to.emplace_back(m.first); + + // Shortcut, when running in serial + if (send_to.size() == 0) + return objects_to_send; const auto receive_from = Utilities::MPI::compute_point_to_point_communication_pattern(comm, @@ -1514,19 +1517,20 @@ namespace Utilities { 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, - mpi_tag, - comm, - &buffer_send_requests[i]); - AssertThrowMPI(ierr); - ++i; - } + if (rank_obj.first != my_proc) + { + 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, + mpi_tag, + comm, + &buffer_send_requests[i]); + AssertThrowMPI(ierr); + ++i; + } } // Receiving buffers @@ -1571,6 +1575,10 @@ namespace Utilities buffer_send_requests.data(), MPI_STATUSES_IGNORE); + // If necessary, insert my own objects into the received_objects + if (objects_to_send.find(my_proc) != objects_to_send.end()) + received_objects[my_proc] = objects_to_send.at(my_proc); + return received_objects; # endif // deal.II with MPI } diff --git a/tests/base/mpi_some_to_some_03.cc b/tests/base/mpi_some_to_some_03.cc new file mode 100644 index 0000000000..1978ef6604 --- /dev/null +++ b/tests/base/mpi_some_to_some_03.cc @@ -0,0 +1,65 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2019 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +// Test Utilities::MPI::some_to_some when sending to myself as well. + +#include +#include +#include + +#include +#include + +#include "../tests.h" + + +void +test() +{ + const auto n_procs = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + const auto my_proc = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + std::map> to_send; + + // send to myself and to my next one + to_send[my_proc].emplace_back(my_proc + 1); + to_send[(my_proc + 1) % n_procs].emplace_back((my_proc + 1) * 10); + + const auto received = Utilities::MPI::some_to_some(MPI_COMM_WORLD, to_send); + + const auto original = Utilities::MPI::some_to_some(MPI_COMM_WORLD, received); + + deallog << "Sent : " + << Patterns::Tools::to_string(to_send) << std::endl; + deallog << "Received : " + << Patterns::Tools::to_string(received) << std::endl; + deallog << "Received(Received) == Sent: " + << Patterns::Tools::to_string(original) << std::endl; + + // now check that to_send and original are the same + if (original == to_send) + deallog << "OK" << std::endl; + else + deallog << "Not OK" << std::endl; +} + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; + + test(); +} diff --git a/tests/base/mpi_some_to_some_03.mpirun=1.output b/tests/base/mpi_some_to_some_03.mpirun=1.output new file mode 100644 index 0000000000..81a581e751 --- /dev/null +++ b/tests/base/mpi_some_to_some_03.mpirun=1.output @@ -0,0 +1,5 @@ + +DEAL:0::Sent : 0:1, 10 +DEAL:0::Received : 0:1, 10 +DEAL:0::Received(Received) == Sent: 0:1, 10 +DEAL:0::OK diff --git a/tests/base/mpi_some_to_some_03.mpirun=2.output b/tests/base/mpi_some_to_some_03.mpirun=2.output new file mode 100644 index 0000000000..b33812ead6 --- /dev/null +++ b/tests/base/mpi_some_to_some_03.mpirun=2.output @@ -0,0 +1,11 @@ + +DEAL:0::Sent : 0:1; 1:10 +DEAL:0::Received : 0:1; 1:20 +DEAL:0::Received(Received) == Sent: 0:1; 1:10 +DEAL:0::OK + +DEAL:1::Sent : 0:20; 1:2 +DEAL:1::Received : 0:10; 1:2 +DEAL:1::Received(Received) == Sent: 0:20; 1:2 +DEAL:1::OK + diff --git a/tests/base/mpi_some_to_some_03.mpirun=3.output b/tests/base/mpi_some_to_some_03.mpirun=3.output new file mode 100644 index 0000000000..09ae742d78 --- /dev/null +++ b/tests/base/mpi_some_to_some_03.mpirun=3.output @@ -0,0 +1,17 @@ + +DEAL:0::Sent : 0:1; 1:10 +DEAL:0::Received : 0:1; 2:30 +DEAL:0::Received(Received) == Sent: 0:1; 1:10 +DEAL:0::OK + +DEAL:1::Sent : 1:2; 2:20 +DEAL:1::Received : 0:10; 1:2 +DEAL:1::Received(Received) == Sent: 1:2; 2:20 +DEAL:1::OK + + +DEAL:2::Sent : 0:30; 2:3 +DEAL:2::Received : 1:20; 2:3 +DEAL:2::Received(Received) == Sent: 0:30; 2:3 +DEAL:2::OK + diff --git a/tests/base/mpi_some_to_some_03.mpirun=4.output b/tests/base/mpi_some_to_some_03.mpirun=4.output new file mode 100644 index 0000000000..3b280f8b21 --- /dev/null +++ b/tests/base/mpi_some_to_some_03.mpirun=4.output @@ -0,0 +1,23 @@ + +DEAL:0::Sent : 0:1; 1:10 +DEAL:0::Received : 0:1; 3:40 +DEAL:0::Received(Received) == Sent: 0:1; 1:10 +DEAL:0::OK + +DEAL:1::Sent : 1:2; 2:20 +DEAL:1::Received : 0:10; 1:2 +DEAL:1::Received(Received) == Sent: 1:2; 2:20 +DEAL:1::OK + + +DEAL:2::Sent : 2:3; 3:30 +DEAL:2::Received : 1:20; 2:3 +DEAL:2::Received(Received) == Sent: 2:3; 3:30 +DEAL:2::OK + + +DEAL:3::Sent : 0:40; 3:4 +DEAL:3::Received : 2:30; 3:4 +DEAL:3::Received(Received) == Sent: 0:40; 3:4 +DEAL:3::OK + -- 2.39.5