From 169587a8e9ea30bd1a7f7e7216109f884bc6f4a5 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Fri, 28 Feb 2020 14:21:13 +0100 Subject: [PATCH] Add AnonymousConsensusAlgorithmProcess --- include/deal.II/base/mpi.h | 2 + .../deal.II/base/mpi_consensus_algorithm.h | 202 ++++++++++++++++++ tests/base/consensus_algorithm_01.cc | 71 ++++++ .../consensus_algorithm_01.mpirun=4.output | 15 ++ 4 files changed, 290 insertions(+) create mode 100644 include/deal.II/base/mpi_consensus_algorithm.h create mode 100644 tests/base/consensus_algorithm_01.cc create mode 100644 tests/base/consensus_algorithm_01.mpirun=4.output diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index 7e31c2f14e..ecb6ff11da 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -1069,6 +1069,8 @@ namespace Utilities const std::vector &recv_buffer); }; + + /** * A base class for algorithms that implement the task of coming up with * communication patterns to retrieve data from other processes in a diff --git a/include/deal.II/base/mpi_consensus_algorithm.h b/include/deal.II/base/mpi_consensus_algorithm.h new file mode 100644 index 0000000000..e9c13eb799 --- /dev/null +++ b/include/deal.II/base/mpi_consensus_algorithm.h @@ -0,0 +1,202 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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. +// +// --------------------------------------------------------------------- + +#ifndef dealii_mpi_consensus_algorithm_h +#define dealii_mpi_consensus_algorithm_h + +#include + +#include +#include + +DEAL_II_NAMESPACE_OPEN + + +namespace Utilities +{ + namespace MPI + { + /** + * This class implements Utilities::MPI::ConsensusAlgorithmProcess, + * using user-provided function wrappers. + * The advantage of this class is that users do not have to write their + * own implementation but can register lambda functions directly. + */ + template + class AnonymousConsensusAlgorithmProcess + : public ConsensusAlgorithmProcess + { + public: + /** + * Register functions that should be called for implementing the interface + * of ConsensusAlgorithmProcess. + * + * @param function_compute_targets called during `compute_targets`. + * @param function_create_request called during `create_request`. + * @param function_answer_request called during `answer_request`. + * @param function_prepare_buffer_for_answer called during + * `prepare_buffer_for_answer`. + * @param function_read_answer called during `read_answer`. + */ + AnonymousConsensusAlgorithmProcess( + const std::function()> + &function_compute_targets, + const std::function &)> + &function_create_request = + [](const unsigned int, std::vector &) {}, + const std::function &, + std::vector &)> &function_answer_request = + [](const unsigned int, const std::vector &, std::vector &) {}, + const std::function &)> + &function_prepare_buffer_for_answer = + [](const unsigned int, std::vector &) {}, + const std::function &)> + &function_read_answer = + [](const unsigned int, const std::vector &) {}); + + /** + * @copydoc ConsensusAlgorithmProcess::compute_targets() + */ + std::vector + compute_targets() override; + + /** + * @copydoc ConsensusAlgorithmProcess::create_request() + */ + void + create_request(const unsigned int other_rank, + std::vector & send_buffer) override; + + /** + * @copydoc ConsensusAlgorithmProcess::answer_request() + */ + void + answer_request(const unsigned int other_rank, + const std::vector &buffer_recv, + std::vector & request_buffer) override; + + /** + * @copydoc ConsensusAlgorithmProcess::prepare_buffer_for_answer() + */ + void + prepare_buffer_for_answer(const unsigned int other_rank, + std::vector & recv_buffer) override; + + /** + * @copydoc ConsensusAlgorithmProcess::read_answer() + */ + void + read_answer(const unsigned int other_rank, + const std::vector &recv_buffer) override; + + private: + const std::function()> function_compute_targets; + const std::function &)> + function_create_request; + const std::function< + void(const unsigned int, const std::vector &, std::vector &)> + function_answer_request; + const std::function &)> + function_prepare_buffer_for_answer; + const std::function &)> + function_read_answer; + }; + + + + template + AnonymousConsensusAlgorithmProcess:: + AnonymousConsensusAlgorithmProcess( + const std::function()> + &function_compute_targets, + const std::function &)> + & function_create_request, + const std::function &, + std::vector &)> &function_answer_request, + const std::function &)> + &function_prepare_buffer_for_answer, + const std::function &)> + &function_read_answer) + : function_compute_targets(function_compute_targets) + , function_create_request(function_create_request) + , function_answer_request(function_answer_request) + , function_prepare_buffer_for_answer(function_prepare_buffer_for_answer) + , function_read_answer(function_read_answer) + {} + + + + template + std::vector + AnonymousConsensusAlgorithmProcess::compute_targets() + { + return function_compute_targets(); + } + + + + template + void + AnonymousConsensusAlgorithmProcess::create_request( + const unsigned int other_rank, + std::vector & send_buffer) + { + function_create_request(other_rank, send_buffer); + } + + + + template + void + AnonymousConsensusAlgorithmProcess::answer_request( + const unsigned int other_rank, + const std::vector &buffer_recv, + std::vector & request_buffer) + { + function_answer_request(other_rank, buffer_recv, request_buffer); + } + + + + template + void + AnonymousConsensusAlgorithmProcess::prepare_buffer_for_answer( + const unsigned int other_rank, + std::vector & recv_buffer) + { + function_prepare_buffer_for_answer(other_rank, recv_buffer); + } + + + + template + void + AnonymousConsensusAlgorithmProcess::read_answer( + const unsigned int other_rank, + const std::vector &recv_buffer) + { + function_read_answer(other_rank, recv_buffer); + } + + + } // end of namespace MPI +} // end of namespace Utilities + + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/tests/base/consensus_algorithm_01.cc b/tests/base/consensus_algorithm_01.cc new file mode 100644 index 0000000000..562bbe72ba --- /dev/null +++ b/tests/base/consensus_algorithm_01.cc @@ -0,0 +1,71 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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 AnonymousConsensusAlgorithmProcess. + +#include + +#include "../tests.h" + +using namespace dealii; + +void +test(const MPI_Comm &comm) +{ + const unsigned int my_rank = dealii::Utilities::MPI::this_mpi_process(comm); + const unsigned int n_rank = dealii::Utilities::MPI::n_mpi_processes(comm); + + using T1 = unsigned int; + using T2 = unsigned int; + + dealii::Utilities::MPI::AnonymousConsensusAlgorithmProcess process( + [&]() { + std::vector result{(my_rank + 1) % n_rank}; + return result; + }, + [&](const unsigned int other_rank, std::vector &send_buffer) { + send_buffer.push_back(my_rank); + }, + [&](const unsigned int & other_rank, + const std::vector &buffer_recv, + std::vector & request_buffer) { + AssertDimension(other_rank, buffer_recv.front()); + deallog << "ConsensusAlgorithmProcess::answer_request() passed!" + << std::endl; + request_buffer.push_back(my_rank); + }, + [&](const unsigned int other_rank, std::vector &recv_buffer) { + recv_buffer.resize(1); + }, + [&](const unsigned int other_rank, const std::vector &recv_buffer) { + AssertDimension(other_rank, recv_buffer.front()); + deallog << "ConsensusAlgorithmProcess::function_read_answer() passed!" + << std::endl; + }); + dealii::Utilities::MPI::ConsensusAlgorithmSelector(process, comm) + .run(); +} + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll all; + + const MPI_Comm comm = MPI_COMM_WORLD; + + test(comm); +} diff --git a/tests/base/consensus_algorithm_01.mpirun=4.output b/tests/base/consensus_algorithm_01.mpirun=4.output new file mode 100644 index 0000000000..8faf9d944b --- /dev/null +++ b/tests/base/consensus_algorithm_01.mpirun=4.output @@ -0,0 +1,15 @@ + +DEAL:0::ConsensusAlgorithmProcess::answer_request() passed! +DEAL:0::ConsensusAlgorithmProcess::function_read_answer() passed! + +DEAL:1::ConsensusAlgorithmProcess::answer_request() passed! +DEAL:1::ConsensusAlgorithmProcess::function_read_answer() passed! + + +DEAL:2::ConsensusAlgorithmProcess::answer_request() passed! +DEAL:2::ConsensusAlgorithmProcess::function_read_answer() passed! + + +DEAL:3::ConsensusAlgorithmProcess::answer_request() passed! +DEAL:3::ConsensusAlgorithmProcess::function_read_answer() passed! + -- 2.39.5