From: Marc Fehling Date: Thu, 4 Mar 2021 05:54:42 +0000 (-0700) Subject: Added Utilities::MPI::logical_or(). X-Git-Tag: v9.3.0-rc1~366^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=832d3df0439d2c6a799c4dca16d4da82e5692288;p=dealii.git Added Utilities::MPI::logical_or(). --- diff --git a/doc/news/changes/minor/20210304Fehling b/doc/news/changes/minor/20210304Fehling new file mode 100644 index 0000000000..fdc0dd2d06 --- /dev/null +++ b/doc/news/changes/minor/20210304Fehling @@ -0,0 +1,3 @@ +New: Utilities::MPI::logical_or() for collective logical or operations. +
+(Marc Fehling, 2021/03/04) diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index e690e74ad5..cb63553c3b 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -53,6 +53,9 @@ using MPI_Op = int; # ifndef MPI_SUM # define MPI_SUM 0 # endif +# ifndef MPI_LOR +# define MPI_LOR 0 +# endif #endif @@ -678,6 +681,66 @@ namespace Utilities const MPI_Comm & mpi_communicator, const ArrayView & minima); + /** + * Performs a logical or operation over all processors of the value + * @p t. The logical or operator `||` returns the boolean value + * `true` if either or all operands are `true` and returns `false` + * otherwise. If the provided value @p t corresponds to `0` in its + * associated data type `T`, it will be interpreted as `false`, and `true` + * otherwise. Data type `T` must be of type `integral`, i.e., `bool`, + * `char`, `short`, `int`, `long`, or any of their variations. + * + * This function is collective over all processors given in the + * @ref GlossMPICommunicator "communicator". + * If deal.II is not configured for use of MPI, this function simply + * returns the value of @p value. This function corresponds to the + * MPI_Allreduce function, i.e., all processors receive the + * result of this operation. + * + * @note Sometimes, not all processors need a result and in that case one + * would call the MPI_Reduce function instead of the + * MPI_Allreduce function. The latter is at most twice as + * expensive, so if you are concerned about performance, it may be + * worthwhile investigating whether your algorithm indeed needs the result + * everywhere. + */ + template + T + logical_or(const T &t, const MPI_Comm &mpi_communicator); + + /** + * Like the previous function, but performs the logical or operation + * on each element of an array. In other words, the i-th element of the + * results array is the result of the logical or operation applied on + * the i-th entries of the input arrays from each processor. T and U must + * decay to the same type, e.g., they just differ by one of them having a + * const type qualifier and the other not. + * + * Input and output arrays may be the same. + * + * @note Depending on your standard library, this function may not work with + * specializations of `std::vector` for the data type `bool`. In that + * case, use a different container or data type. + */ + template + void + logical_or(const T &values, const MPI_Comm &mpi_communicator, U &results); + + /** + * Like the previous function, but performs the logical or operation + * on each element of an array as specified by the ArrayView arguments. + * In other words, the i-th element of the results array is the result of + * the logical or operation applied on the i-th entries of the input + * arrays from each processor. + * + * Input and output arrays may be the same. + */ + template + void + logical_or(const ArrayView &values, + const MPI_Comm & mpi_communicator, + const ArrayView & results); + /** * A data structure to store the result of the min_max_avg() function. * The structure stores the minimum, maximum, and average of one @@ -1131,6 +1194,21 @@ namespace Utilities ArrayView(minima, N)); } + template + void + logical_or(const T (&values)[N], + const MPI_Comm &mpi_communicator, + T (&results)[N]) + { + static_assert(std::is_integral::value, + "The MPI_LOR operation only allows integral data types."); + + internal::all_reduce(MPI_LOR, + ArrayView(values, N), + mpi_communicator, + ArrayView(results, N)); + } + template std::map some_to_some(const MPI_Comm & comm, diff --git a/include/deal.II/base/mpi.templates.h b/include/deal.II/base/mpi.templates.h index fa3cf3a31d..1e2c6d960c 100644 --- a/include/deal.II/base/mpi.templates.h +++ b/include/deal.II/base/mpi.templates.h @@ -43,6 +43,18 @@ namespace Utilities /** * Return the corresponding MPI data type id for the argument given. */ + inline MPI_Datatype + mpi_type_id(const bool *) + { +# if DEAL_II_MPI_VERSION_GTE(2, 2) + return MPI_CXX_BOOL; +# else + return MPI_C_BOOL; +# endif + } + + + inline MPI_Datatype mpi_type_id(const char *) { @@ -456,6 +468,66 @@ namespace Utilities + template + T + logical_or(const T &t, const MPI_Comm &mpi_communicator) + { + static_assert(std::is_integral::value, + "The MPI_LOR operation only allows integral data types."); + + T return_value{}; + internal::all_reduce(MPI_LOR, + ArrayView(&t, 1), + mpi_communicator, + ArrayView(&return_value, 1)); + return return_value; + } + + + + template + void + logical_or(const T &values, const MPI_Comm &mpi_communicator, U &results) + { + static_assert(std::is_same::type, + typename std::decay::type>::value, + "Input and output arguments must have the same type!"); + + static_assert(std::is_integral::value, + "The MPI_LOR operation only allows integral data types."); + + // Specializations of std containers for the data type bool do not + // necessarily store its elements as a contiguous array. Thus we will use + // the make_array_view() function with iterators here, which verifies + // this. + const auto array_view_values = + make_array_view(values.cbegin(), values.cend()); + const auto array_view_results = + make_array_view(results.begin(), results.end()); + + using const_type = + ArrayView; + logical_or(static_cast(array_view_values), + mpi_communicator, + array_view_results); + } + + + + template + void + logical_or(const ArrayView &values, + const MPI_Comm & mpi_communicator, + const ArrayView & results) + { + static_assert(std::is_integral::value, + "The MPI_LOR operation only allows integral data types."); + + internal::all_reduce(MPI_LOR, values, mpi_communicator, results); + } + + + template T all_reduce(const T & vec, diff --git a/source/base/mpi.cc b/source/base/mpi.cc index febffcd3a7..91b3ab81d4 100644 --- a/source/base/mpi.cc +++ b/source/base/mpi.cc @@ -1188,6 +1188,16 @@ namespace Utilities } + template bool + logical_or(const bool &, const MPI_Comm &); + + + template void + logical_or(const ArrayView &, + const MPI_Comm &, + const ArrayView &); + + template std::vector compute_set_union(const std::vector &vec, const MPI_Comm & comm); diff --git a/tests/mpi/collective_logical_or.cc b/tests/mpi/collective_logical_or.cc new file mode 100644 index 0000000000..fe4edba30b --- /dev/null +++ b/tests/mpi/collective_logical_or.cc @@ -0,0 +1,63 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 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::logical_or() + +#include + +#include "../tests.h" + +void +test() +{ + const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + bool collective_false = Utilities::MPI::logical_or(false, MPI_COMM_WORLD); + bool collective_alternating = + Utilities::MPI::logical_or(((myid % 2) == 0) ? true : false, + MPI_COMM_WORLD); + + if (myid == 0) + deallog << std::boolalpha << collective_false << ' ' + << collective_alternating << std::endl; +} + + +int +main(int argc, char *argv[]) +{ +#ifdef DEAL_II_WITH_MPI + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); +#else + (void)argc; + (void)argv; + compile_time_error; + +#endif + + if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0) + { + initlog(); + + deallog.push("mpi"); + test(); + deallog.pop(); + } + else + test(); +} diff --git a/tests/mpi/collective_logical_or.mpirun=1.output b/tests/mpi/collective_logical_or.mpirun=1.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or.mpirun=1.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or.mpirun=10.output b/tests/mpi/collective_logical_or.mpirun=10.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or.mpirun=10.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or.mpirun=4.output b/tests/mpi/collective_logical_or.mpirun=4.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or.mpirun=4.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_array.cc b/tests/mpi/collective_logical_or_array.cc new file mode 100644 index 0000000000..89a2627cc8 --- /dev/null +++ b/tests/mpi/collective_logical_or_array.cc @@ -0,0 +1,63 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 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::logical_or() for arrays + +#include + +#include "../tests.h" + +void +test() +{ + const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + bool values[2] = {false, (myid % 2) == 0}; + bool results[2]; + Utilities::MPI::logical_or(values, MPI_COMM_WORLD, results); + Assert(results[0] == false, ExcInternalError()); + Assert(results[1] == true, ExcInternalError()); + + if (myid == 0) + deallog << std::boolalpha << results[0] << ' ' << results[1] << std::endl; +} + + +int +main(int argc, char *argv[]) +{ +#ifdef DEAL_II_WITH_MPI + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); +#else + (void)argc; + (void)argv; + compile_time_error; + +#endif + + if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0) + { + initlog(); + + deallog.push("mpi"); + test(); + deallog.pop(); + } + else + test(); +} diff --git a/tests/mpi/collective_logical_or_array.mpirun=1.output b/tests/mpi/collective_logical_or_array.mpirun=1.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_array.mpirun=1.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_array.mpirun=10.output b/tests/mpi/collective_logical_or_array.mpirun=10.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_array.mpirun=10.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_array.mpirun=4.output b/tests/mpi/collective_logical_or_array.mpirun=4.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_array.mpirun=4.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_array_in_place.cc b/tests/mpi/collective_logical_or_array_in_place.cc new file mode 100644 index 0000000000..5d75d740cb --- /dev/null +++ b/tests/mpi/collective_logical_or_array_in_place.cc @@ -0,0 +1,62 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 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::logical_or() for arrays, but with input=output + +#include + +#include "../tests.h" + +void +test() +{ + const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + bool results[2] = {false, (myid % 2) == 0}; + Utilities::MPI::logical_or(results, MPI_COMM_WORLD, results); + Assert(results[0] == false, ExcInternalError()); + Assert(results[1] == true, ExcInternalError()); + + if (myid == 0) + deallog << std::boolalpha << results[0] << ' ' << results[1] << std::endl; +} + + +int +main(int argc, char *argv[]) +{ +#ifdef DEAL_II_WITH_MPI + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); +#else + (void)argc; + (void)argv; + compile_time_error; + +#endif + + if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0) + { + initlog(); + + deallog.push("mpi"); + test(); + deallog.pop(); + } + else + test(); +} diff --git a/tests/mpi/collective_logical_or_array_in_place.mpirun=1.output b/tests/mpi/collective_logical_or_array_in_place.mpirun=1.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_array_in_place.mpirun=1.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_array_in_place.mpirun=10.output b/tests/mpi/collective_logical_or_array_in_place.mpirun=10.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_array_in_place.mpirun=10.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_array_in_place.mpirun=4.output b/tests/mpi/collective_logical_or_array_in_place.mpirun=4.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_array_in_place.mpirun=4.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_vector.cc b/tests/mpi/collective_logical_or_vector.cc new file mode 100644 index 0000000000..62076cd263 --- /dev/null +++ b/tests/mpi/collective_logical_or_vector.cc @@ -0,0 +1,68 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 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::logical_or() for vectors +// +// std::vector does not necessarily store elements in contiguous array, +// use std::vector instead + +#include +#include + +#include "../tests.h" + +void +test() +{ + const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + std::vector values = {0, (myid % 2) == 0}; + std::vector results(2); + Utilities::MPI::logical_or(values, MPI_COMM_WORLD, results); + Assert(results[0] == false, ExcInternalError()); + Assert(results[1] == true, ExcInternalError()); + + if (myid == 0) + deallog << std::boolalpha << static_cast(results[0]) << ' ' + << static_cast(results[1]) << std::endl; +} + + +int +main(int argc, char *argv[]) +{ +#ifdef DEAL_II_WITH_MPI + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); +#else + (void)argc; + (void)argv; + compile_time_error; + +#endif + + if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0) + { + initlog(); + + deallog.push("mpi"); + test(); + deallog.pop(); + } + else + test(); +} diff --git a/tests/mpi/collective_logical_or_vector.mpirun=1.output b/tests/mpi/collective_logical_or_vector.mpirun=1.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_vector.mpirun=1.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_vector.mpirun=10.output b/tests/mpi/collective_logical_or_vector.mpirun=10.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_vector.mpirun=10.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_vector.mpirun=4.output b/tests/mpi/collective_logical_or_vector.mpirun=4.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_vector.mpirun=4.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_vector_in_place.cc b/tests/mpi/collective_logical_or_vector_in_place.cc new file mode 100644 index 0000000000..cb134887f4 --- /dev/null +++ b/tests/mpi/collective_logical_or_vector_in_place.cc @@ -0,0 +1,67 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 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::logical_or for vectors, but with input=output +// +// std::vector does not necessarily store elements in contiguous array, +// use std::vector instead + +#include +#include + +#include "../tests.h" + +void +test() +{ + const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + std::vector results = {0, (myid % 2) == 0}; + Utilities::MPI::logical_or(results, MPI_COMM_WORLD, results); + Assert(results[0] == false, ExcInternalError()); + Assert(results[1] == true, ExcInternalError()); + + if (myid == 0) + deallog << std::boolalpha << static_cast(results[0]) << ' ' + << static_cast(results[1]) << std::endl; +} + + +int +main(int argc, char *argv[]) +{ +#ifdef DEAL_II_WITH_MPI + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); +#else + (void)argc; + (void)argv; + compile_time_error; + +#endif + + if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0) + { + initlog(); + + deallog.push("mpi"); + test(); + deallog.pop(); + } + else + test(); +} diff --git a/tests/mpi/collective_logical_or_vector_in_place.mpirun=1.output b/tests/mpi/collective_logical_or_vector_in_place.mpirun=1.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_vector_in_place.mpirun=1.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_vector_in_place.mpirun=10.output b/tests/mpi/collective_logical_or_vector_in_place.mpirun=10.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_vector_in_place.mpirun=10.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true diff --git a/tests/mpi/collective_logical_or_vector_in_place.mpirun=4.output b/tests/mpi/collective_logical_or_vector_in_place.mpirun=4.output new file mode 100644 index 0000000000..314bcbafde --- /dev/null +++ b/tests/mpi/collective_logical_or_vector_in_place.mpirun=4.output @@ -0,0 +1,2 @@ + +DEAL:mpi::false true