]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added Utilities::MPI::logical_or(). 11850/head
authorMarc Fehling <mafehling.git@gmail.com>
Thu, 4 Mar 2021 05:54:42 +0000 (22:54 -0700)
committerMarc Fehling <mafehling.git@gmail.com>
Sun, 7 Mar 2021 03:57:33 +0000 (20:57 -0700)
24 files changed:
doc/news/changes/minor/20210304Fehling [new file with mode: 0644]
include/deal.II/base/mpi.h
include/deal.II/base/mpi.templates.h
source/base/mpi.cc
tests/mpi/collective_logical_or.cc [new file with mode: 0644]
tests/mpi/collective_logical_or.mpirun=1.output [new file with mode: 0644]
tests/mpi/collective_logical_or.mpirun=10.output [new file with mode: 0644]
tests/mpi/collective_logical_or.mpirun=4.output [new file with mode: 0644]
tests/mpi/collective_logical_or_array.cc [new file with mode: 0644]
tests/mpi/collective_logical_or_array.mpirun=1.output [new file with mode: 0644]
tests/mpi/collective_logical_or_array.mpirun=10.output [new file with mode: 0644]
tests/mpi/collective_logical_or_array.mpirun=4.output [new file with mode: 0644]
tests/mpi/collective_logical_or_array_in_place.cc [new file with mode: 0644]
tests/mpi/collective_logical_or_array_in_place.mpirun=1.output [new file with mode: 0644]
tests/mpi/collective_logical_or_array_in_place.mpirun=10.output [new file with mode: 0644]
tests/mpi/collective_logical_or_array_in_place.mpirun=4.output [new file with mode: 0644]
tests/mpi/collective_logical_or_vector.cc [new file with mode: 0644]
tests/mpi/collective_logical_or_vector.mpirun=1.output [new file with mode: 0644]
tests/mpi/collective_logical_or_vector.mpirun=10.output [new file with mode: 0644]
tests/mpi/collective_logical_or_vector.mpirun=4.output [new file with mode: 0644]
tests/mpi/collective_logical_or_vector_in_place.cc [new file with mode: 0644]
tests/mpi/collective_logical_or_vector_in_place.mpirun=1.output [new file with mode: 0644]
tests/mpi/collective_logical_or_vector_in_place.mpirun=10.output [new file with mode: 0644]
tests/mpi/collective_logical_or_vector_in_place.mpirun=4.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20210304Fehling b/doc/news/changes/minor/20210304Fehling
new file mode 100644 (file)
index 0000000..fdc0dd2
--- /dev/null
@@ -0,0 +1,3 @@
+New: Utilities::MPI::logical_or() for collective <i>logical or</i> operations.
+<br>
+(Marc Fehling, 2021/03/04)
index e690e74ad55cec5dd9429c51931afadb6f503c47..cb63553c3b339b7f5b905d7de33d4b851d3c94b3 100644 (file)
@@ -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<T> &      minima);
 
+    /**
+     * Performs a <i>logical or</i> operation over all processors of the value
+     * @p t. The <i>logical or</i> 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
+     * <code>MPI_Allreduce</code> 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 <code>MPI_Reduce</code> function instead of the
+     * <code>MPI_Allreduce</code> 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 <typename T>
+    T
+    logical_or(const T &t, const MPI_Comm &mpi_communicator);
+
+    /**
+     * Like the previous function, but performs the <i>logical or</i> operation
+     * on each element of an array. In other words, the i-th element of the
+     * results array is the result of the <i>logical or</i> 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 <typename T, typename U>
+    void
+    logical_or(const T &values, const MPI_Comm &mpi_communicator, U &results);
+
+    /**
+     * Like the previous function, but performs the <i>logical or</i> 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 <i>logical or</i> operation applied on the i-th entries of the input
+     * arrays from each processor.
+     *
+     * Input and output arrays may be the same.
+     */
+    template <typename T>
+    void
+    logical_or(const ArrayView<const T> &values,
+               const MPI_Comm &          mpi_communicator,
+               const ArrayView<T> &      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<T>(minima, N));
     }
 
+    template <typename T, unsigned int N>
+    void
+    logical_or(const T (&values)[N],
+               const MPI_Comm &mpi_communicator,
+               T (&results)[N])
+    {
+      static_assert(std::is_integral<T>::value,
+                    "The MPI_LOR operation only allows integral data types.");
+
+      internal::all_reduce(MPI_LOR,
+                           ArrayView<const T>(values, N),
+                           mpi_communicator,
+                           ArrayView<T>(results, N));
+    }
+
     template <typename T>
     std::map<unsigned int, T>
     some_to_some(const MPI_Comm &                 comm,
index fa3cf3a31d433cc163a30699fd127c6fcd6368d2..1e2c6d960cb687d62eb12e47e112c11df4a701ed 100644 (file)
@@ -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 <typename T>
+    T
+    logical_or(const T &t, const MPI_Comm &mpi_communicator)
+    {
+      static_assert(std::is_integral<T>::value,
+                    "The MPI_LOR operation only allows integral data types.");
+
+      T return_value{};
+      internal::all_reduce(MPI_LOR,
+                           ArrayView<const T>(&t, 1),
+                           mpi_communicator,
+                           ArrayView<T>(&return_value, 1));
+      return return_value;
+    }
+
+
+
+    template <typename T, typename U>
+    void
+    logical_or(const T &values, const MPI_Comm &mpi_communicator, U &results)
+    {
+      static_assert(std::is_same<typename std::decay<T>::type,
+                                 typename std::decay<U>::type>::value,
+                    "Input and output arguments must have the same type!");
+
+      static_assert(std::is_integral<typename T::value_type>::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<const typename decltype(array_view_values)::value_type>;
+      logical_or(static_cast<const_type>(array_view_values),
+                 mpi_communicator,
+                 array_view_results);
+    }
+
+
+
+    template <typename T>
+    void
+    logical_or(const ArrayView<const T> &values,
+               const MPI_Comm &          mpi_communicator,
+               const ArrayView<T> &      results)
+    {
+      static_assert(std::is_integral<T>::value,
+                    "The MPI_LOR operation only allows integral data types.");
+
+      internal::all_reduce(MPI_LOR, values, mpi_communicator, results);
+    }
+
+
+
     template <typename T>
     T
     all_reduce(const T &                                     vec,
index febffcd3a7f9364ed0b347e260d992f10b90999d..91b3ab81d4b817ecbbcb24a62480c8074c11ec75 100644 (file)
@@ -1188,6 +1188,16 @@ namespace Utilities
     }
 
 
+    template bool
+    logical_or<bool>(const bool &, const MPI_Comm &);
+
+
+    template void
+    logical_or<bool>(const ArrayView<const bool> &,
+                     const MPI_Comm &,
+                     const ArrayView<bool> &);
+
+
     template std::vector<unsigned int>
     compute_set_union(const std::vector<unsigned int> &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 (file)
index 0000000..fe4edba
--- /dev/null
@@ -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 <deal.II/base/utilities.h>
+
+#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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..89a2627
--- /dev/null
@@ -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 <deal.II/base/utilities.h>
+
+#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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..5d75d74
--- /dev/null
@@ -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 <deal.II/base/utilities.h>
+
+#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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..62076cd
--- /dev/null
@@ -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<bool> does not necessarily store elements in contiguous array,
+// use std::vector<char> instead
+
+#include <deal.II/base/mpi.templates.h>
+#include <deal.II/base/utilities.h>
+
+#include "../tests.h"
+
+void
+test()
+{
+  const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+
+  std::vector<char> values = {0, (myid % 2) == 0};
+  std::vector<char> 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<bool>(results[0]) << ' '
+            << static_cast<bool>(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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..cb13488
--- /dev/null
@@ -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<bool> does not necessarily store elements in contiguous array,
+// use std::vector<char> instead
+
+#include <deal.II/base/mpi.templates.h>
+#include <deal.II/base/utilities.h>
+
+#include "../tests.h"
+
+void
+test()
+{
+  const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+
+  std::vector<char> 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<bool>(results[0]) << ' '
+            << static_cast<bool>(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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -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 (file)
index 0000000..314bcba
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL:mpi::false true

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.