From 888fcf659a419a4601492f6bdbf357e80fc13bda Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Mon, 3 Feb 2020 18:35:17 +0100 Subject: [PATCH] Add Utilities::MPI::min_max_avg() for vectors --- doc/news/changes/minor/20200203PeterMunch | 4 + include/deal.II/base/mpi.h | 32 ++++ source/base/mpi.cc | 173 +++++++++++------- tests/mpi/minmaxavg_vector_01.cc | 68 +++++++ tests/mpi/minmaxavg_vector_01.mpirun=4.output | 31 ++++ 5 files changed, 245 insertions(+), 63 deletions(-) create mode 100644 doc/news/changes/minor/20200203PeterMunch create mode 100644 tests/mpi/minmaxavg_vector_01.cc create mode 100644 tests/mpi/minmaxavg_vector_01.mpirun=4.output diff --git a/doc/news/changes/minor/20200203PeterMunch b/doc/news/changes/minor/20200203PeterMunch new file mode 100644 index 0000000000..75d66ea17d --- /dev/null +++ b/doc/news/changes/minor/20200203PeterMunch @@ -0,0 +1,4 @@ +New: The function Utilities::MPI::min_max_avg() can now compute the sum, average, minimum, maximum, and the +process id of the minimum and maximum for each entry of a ArrayView. +
+(Peter Munch, 2020/02/03) diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index fd954ec8cc..7e31c2f14e 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -719,6 +719,38 @@ namespace Utilities MinMaxAvg min_max_avg(const double my_value, const MPI_Comm &mpi_communicator); + /** + * Same as above but returning the sum, average, minimum, maximum, + * process id of minimum and maximum as a collective operation on the + * given MPI @ref GlossMPICommunicator "communicator" @p mpi_communicator + * for each entry of the vector. + * + * @note This function performs a single reduction sweep. + * + * @pre Size of the input vector has to be the same on all processes. + */ + std::vector + min_max_avg(const std::vector &my_value, + const MPI_Comm & mpi_communicator); + + + /** + * Same as above but returning the sum, average, minimum, maximum, + * process id of minimum and maximum as a collective operation on the + * given MPI @ref GlossMPICommunicator "communicator" @p mpi_communicator + * for each entry of the ArrayView. + * + * @note This function performs a single reduction sweep. + * + * @pre Size of the input ArrayView has to be the same on all processes + * and the input and output ArrayVew have to have the same size. + */ + void + min_max_avg(const ArrayView &my_values, + const ArrayView & result, + 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 diff --git a/source/base/mpi.cc b/source/base/mpi.cc index ae7df36c50..c7c46c4021 100644 --- a/source/base/mpi.cc +++ b/source/base/mpi.cc @@ -70,6 +70,31 @@ namespace Utilities { namespace MPI { + MinMaxAvg + min_max_avg(const double my_value, const MPI_Comm &mpi_communicator) + { + MinMaxAvg result; + min_max_avg(ArrayView(my_value), + ArrayView(result), + mpi_communicator); + + return result; + } + + + + std::vector + min_max_avg(const std::vector &my_values, + const MPI_Comm & mpi_communicator) + { + std::vector results(my_values.size()); + min_max_avg(my_values, results, mpi_communicator); + + return results; + } + + + #ifdef DEAL_II_WITH_MPI unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator) @@ -491,68 +516,80 @@ namespace Utilities int * len, MPI_Datatype *) { - (void)len; const MinMaxAvg *in_lhs = static_cast(in_lhs_); MinMaxAvg * inout_rhs = static_cast(inout_rhs_); - Assert(*len == 1, ExcInternalError()); - - inout_rhs->sum += in_lhs->sum; - if (inout_rhs->min > in_lhs->min) + for (int i = 0; i < *len; i++) { - inout_rhs->min = in_lhs->min; - inout_rhs->min_index = in_lhs->min_index; - } - else if (inout_rhs->min == in_lhs->min) - { - // choose lower cpu index when tied to make operator commutative - if (inout_rhs->min_index > in_lhs->min_index) - inout_rhs->min_index = in_lhs->min_index; - } + inout_rhs[i].sum += in_lhs[i].sum; + if (inout_rhs[i].min > in_lhs[i].min) + { + inout_rhs[i].min = in_lhs[i].min; + inout_rhs[i].min_index = in_lhs[i].min_index; + } + else if (inout_rhs[i].min == in_lhs[i].min) + { + // choose lower cpu index when tied to make operator commutative + if (inout_rhs[i].min_index > in_lhs[i].min_index) + inout_rhs[i].min_index = in_lhs[i].min_index; + } - if (inout_rhs->max < in_lhs->max) - { - inout_rhs->max = in_lhs->max; - inout_rhs->max_index = in_lhs->max_index; - } - else if (inout_rhs->max == in_lhs->max) - { - // choose lower cpu index when tied to make operator commutative - if (inout_rhs->max_index > in_lhs->max_index) - inout_rhs->max_index = in_lhs->max_index; + if (inout_rhs[i].max < in_lhs[i].max) + { + inout_rhs[i].max = in_lhs[i].max; + inout_rhs[i].max_index = in_lhs[i].max_index; + } + else if (inout_rhs[i].max == in_lhs[i].max) + { + // choose lower cpu index when tied to make operator commutative + if (inout_rhs[i].max_index > in_lhs[i].max_index) + inout_rhs[i].max_index = in_lhs[i].max_index; + } } } } // namespace - MinMaxAvg - min_max_avg(const double my_value, const MPI_Comm &mpi_communicator) + void + min_max_avg(const ArrayView &my_values, + const ArrayView & result, + const MPI_Comm & mpi_communicator) { // If MPI was not started, we have a serial computation and cannot run // the other MPI commands if (job_supports_mpi() == false || Utilities::MPI::n_mpi_processes(mpi_communicator) <= 1) { - MinMaxAvg result; - result.sum = my_value; - result.avg = my_value; - result.min = my_value; - result.max = my_value; - result.min_index = 0; - result.max_index = 0; - - return result; + for (unsigned int i = 0; i < my_values.size(); i++) + { + result[i].sum = my_values[i]; + result[i].avg = my_values[i]; + result[i].min = my_values[i]; + result[i].max = my_values[i]; + result[i].min_index = 0; + result[i].max_index = 0; + } } + AssertDimension(Utilities::MPI::min(my_values.size(), mpi_communicator), + Utilities::MPI::max(my_values.size(), mpi_communicator)); + + AssertDimension(my_values.size(), result.size()); + + + // To avoid uninitialized values on some MPI implementations, provide // result with a default value already... - MinMaxAvg result = {0., - std::numeric_limits::max(), - -std::numeric_limits::max(), - 0, - 0, - 0.}; + MinMaxAvg dummy = {0., + std::numeric_limits::max(), + -std::numeric_limits::max(), + 0, + 0, + 0.}; + + for (auto &i : result) + i = dummy; const unsigned int my_id = dealii::Utilities::MPI::this_mpi_process(mpi_communicator); @@ -566,21 +603,28 @@ namespace Utilities &op); AssertThrowMPI(ierr); - MinMaxAvg in; - in.sum = in.min = in.max = my_value; - in.min_index = in.max_index = my_id; + std::vector in(my_values.size()); + + for (unsigned int i = 0; i < my_values.size(); i++) + { + in[i].sum = in[i].min = in[i].max = my_values[i]; + in[i].min_index = in[i].max_index = my_id; + } MPI_Datatype type; - int lengths[] = {3, 2}; - MPI_Aint displacements[] = {0, offsetof(MinMaxAvg, min_index)}; - MPI_Datatype types[] = {MPI_DOUBLE, MPI_INT}; + int lengths[] = {3, 2, 1}; + MPI_Aint displacements[] = {0, + offsetof(MinMaxAvg, min_index), + offsetof(MinMaxAvg, avg)}; + MPI_Datatype types[] = {MPI_DOUBLE, MPI_INT, MPI_DOUBLE}; - ierr = MPI_Type_create_struct(2, lengths, displacements, types, &type); + ierr = MPI_Type_create_struct(3, lengths, displacements, types, &type); AssertThrowMPI(ierr); ierr = MPI_Type_commit(&type); AssertThrowMPI(ierr); - ierr = MPI_Allreduce(&in, &result, 1, type, op, mpi_communicator); + ierr = MPI_Allreduce( + in.data(), result.data(), my_values.size(), type, op, mpi_communicator); AssertThrowMPI(ierr); ierr = MPI_Type_free(&type); @@ -589,11 +633,11 @@ namespace Utilities ierr = MPI_Op_free(&op); AssertThrowMPI(ierr); - result.avg = result.sum / numproc; - - return result; + for (auto &r : result) + r.avg = r.sum / numproc; } + #else unsigned int @@ -635,19 +679,22 @@ namespace Utilities - MinMaxAvg - min_max_avg(const double my_value, const MPI_Comm &) + void + min_max_avg(const ArrayView &my_values, + const ArrayView & result, + const MPI_Comm &) { - MinMaxAvg result; - - result.sum = my_value; - result.avg = my_value; - result.min = my_value; - result.max = my_value; - result.min_index = 0; - result.max_index = 0; + AssertDimension(my_values.size(), result.size()); - return result; + for (unsigned int i = 0; i < my_values.size(); i++) + { + result[i].sum = my_values[i]; + result[i].avg = my_values[i]; + result[i].min = my_values[i]; + result[i].max = my_values[i]; + result[i].min_index = 0; + result[i].max_index = 0; + } } #endif diff --git a/tests/mpi/minmaxavg_vector_01.cc b/tests/mpi/minmaxavg_vector_01.cc new file mode 100644 index 0000000000..253561d00d --- /dev/null +++ b/tests/mpi/minmaxavg_vector_01.cc @@ -0,0 +1,68 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + + +// check Utilities::MPI::min_max_avg() for std::vector and ArrayView + +#include + +#include "../tests.h" + +void +print_it(const Utilities::MPI::MinMaxAvg &result) +{ + deallog << "sum: " << result.sum << " avg: " << result.avg + << " min: " << result.min << " @" << result.min_index + << " max: " << result.max << " @" << result.max_index << std::endl; +} + +void +test() +{ + unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + std::vector values; + + if (myid == 0) + values = std::vector{1.0, 1.0, 3.0}; + else if (myid == 1) + values = std::vector{2.0, 3.0, 2.0}; + else if (myid == 2) + values = std::vector{3.0, 2.0, 1.0}; + else if (myid == 3) + values = std::vector{3.0, 3.0, 0.0}; + + std::vector results(values.size()); + + // test ArrayView + Utilities::MPI::min_max_avg(values, results, MPI_COMM_WORLD); + for (const auto &result : results) + print_it(result); + + // test std::vector + for (const auto &result : Utilities::MPI::min_max_avg(values, MPI_COMM_WORLD)) + print_it(result); +} + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll all; + + test(); +} diff --git a/tests/mpi/minmaxavg_vector_01.mpirun=4.output b/tests/mpi/minmaxavg_vector_01.mpirun=4.output new file mode 100644 index 0000000000..40bc8d44b8 --- /dev/null +++ b/tests/mpi/minmaxavg_vector_01.mpirun=4.output @@ -0,0 +1,31 @@ + +DEAL:0::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @2 +DEAL:0::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @1 +DEAL:0::sum: 6.00000 avg: 1.50000 min: 0.00000 @3 max: 3.00000 @0 +DEAL:0::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @2 +DEAL:0::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @1 +DEAL:0::sum: 6.00000 avg: 1.50000 min: 0.00000 @3 max: 3.00000 @0 + +DEAL:1::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @2 +DEAL:1::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @1 +DEAL:1::sum: 6.00000 avg: 1.50000 min: 0.00000 @3 max: 3.00000 @0 +DEAL:1::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @2 +DEAL:1::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @1 +DEAL:1::sum: 6.00000 avg: 1.50000 min: 0.00000 @3 max: 3.00000 @0 + + +DEAL:2::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @2 +DEAL:2::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @1 +DEAL:2::sum: 6.00000 avg: 1.50000 min: 0.00000 @3 max: 3.00000 @0 +DEAL:2::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @2 +DEAL:2::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @1 +DEAL:2::sum: 6.00000 avg: 1.50000 min: 0.00000 @3 max: 3.00000 @0 + + +DEAL:3::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @2 +DEAL:3::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @1 +DEAL:3::sum: 6.00000 avg: 1.50000 min: 0.00000 @3 max: 3.00000 @0 +DEAL:3::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @2 +DEAL:3::sum: 9.00000 avg: 2.25000 min: 1.00000 @0 max: 3.00000 @1 +DEAL:3::sum: 6.00000 avg: 1.50000 min: 0.00000 @3 max: 3.00000 @0 + -- 2.39.5