From: Daniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Date: Sat, 20 Jan 2018 14:51:00 +0000 (+0100)
Subject: Provide general template interface for Utilities::MPI::sum/min/max
X-Git-Tag: v9.0.0-rc1~532^2~1
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=af0d3a75b5edd997e221d422dedc8d7af151ba41;p=dealii.git

Provide general template interface for Utilities::MPI::sum/min/max
---

diff --git a/include/deal.II/base/array_view.h b/include/deal.II/base/array_view.h
index 013e86a034..f003e5bad4 100644
--- a/include/deal.II/base/array_view.h
+++ b/include/deal.II/base/array_view.h
@@ -19,6 +19,7 @@
 #include <deal.II/base/config.h>
 #include <deal.II/base/exceptions.h>
 #include <deal.II/base/table.h>
+#include <deal.II/lac/lapack_full_matrix.h>
 
 #include <type_traits>
 #include <vector>
@@ -420,6 +421,30 @@ ArrayView<ElementType>::operator[](const std::size_t i) const
 }
 
 
+
+/**
+ * Create a view to an entire C-style array. This is equivalent to
+ * initializing an ArrayView object with a pointer to the first element and
+ * the size of the given argument.
+ *
+ * Whether the resulting ArrayView is writable or not depends on the
+ * ElementType being a const type or not.
+ *
+ * @param[in] array The C-style array for which we want to have an ArrayView
+ * object. The ArrayView corresponds to the <em>entire</em> vector.
+ *
+ * @relates ArrayView
+ */
+template <typename ElementType, int N>
+inline
+ArrayView<ElementType>
+make_array_view (ElementType (&array)[N])
+{
+  return ArrayView<ElementType> (array, N);
+}
+
+
+
 /**
  * Create a view to an entire std::vector object. This is equivalent to
  * initializing an ArrayView object with a pointer to the first element and
@@ -444,6 +469,53 @@ make_array_view (std::vector<ElementType> &vector)
 
 
 
+/**
+ * Create a view to an entire Vector object. This is equivalent to
+ * initializing an ArrayView object with a pointer to the first element and
+ * the size of the given argument.
+ *
+ * This function is used for @p const references to objects of Vector type
+ * because they contain immutable elements. Consequently, the return type of
+ * this function is a view to a set of @p const objects.
+ *
+ * @param[in] vector The Vector for which we want to have an array view
+ * object. The array view corresponds to the <em>entire</em> Vector.
+ *
+ * @relates ArrayView
+ */
+template <typename ElementType>
+inline
+ArrayView<ElementType>
+make_array_view (Vector<ElementType> &vector)
+{
+  return ArrayView<ElementType> (vector.begin(), vector.size());
+}
+
+
+/**
+ * Create a view to an entire Vector object. This is equivalent to
+ * initializing an ArrayView object with a pointer to the first element and
+ * the size of the given argument.
+ *
+ * This function is used for @p const references to objects of Vector type
+ * because they contain immutable elements. Consequently, the return type of
+ * this function is a view to a set of @p const objects.
+ *
+ * @param[in] vector The Vector for which we want to have an array view
+ * object. The array view corresponds to the <em>entire</em> Vector.
+ *
+ * @relates ArrayView
+ */
+template <typename ElementType>
+inline
+ArrayView<const ElementType>
+make_array_view (const Vector<ElementType> &vector)
+{
+  return ArrayView<const ElementType> (vector.begin(), vector.size());
+}
+
+
+
 /**
  * Create a view to an entire std::vector object. This is equivalent to
  * initializing an ArrayView object with a pointer to the first element and
@@ -643,9 +715,10 @@ make_array_view (Table<2,ElementType> &table)
  * initializing an ArrayView object with a pointer to the first element of the
  * given table, and the number of table entries as the length of the view.
  *
- * This function is used for @p const references to objects of Table type
- * because they contain immutable elements. Consequently, the return type of
- * this function is a view to a set of @p const objects.
+ * This function is used for non-@p const references to objects of
+ * LAPACKFullMatrix type. Such objects contain elements that can be
+ * written to. Consequently, the return type of this function is a
+ * view to a set of writable objects.
  *
  * @param[in] table The Table for which we want to have an array view object.
  * The array view corresponds to the <em>entire</em> table.
@@ -661,6 +734,51 @@ make_array_view (const Table<2,ElementType> &table)
 }
 
 
+/**
+ * Create a view to an entire LAPACKFullMatrix object. This is equivalent to
+ * initializing an ArrayView object with a pointer to the first element of the
+ * given object, and the number entries as the length of the view.
+ *
+ * This function is used for @p non-const references to objects of Table type
+ * because they contain immutable elements. Consequently, the return type of
+ * this function is a view to a set of @p non-const objects.
+ *
+ * @param[in] table The LAPACKFullMatrix for which we want to have an array view
+ * object. The array view corresponds to the <em>entire</em> object.
+ *
+ * @relates ArrayView
+ */
+template <typename ElementType>
+inline
+ArrayView<ElementType>
+make_array_view (LAPACKFullMatrix<ElementType> &matrix)
+{
+  return ArrayView<ElementType> (&matrix(0,0), matrix.n_elements());
+}
+
+
+/**
+ * Create a view to an entire LAPACKFullMatrix object. This is equivalent to
+ * initializing an ArrayView object with a pointer to the first element of the
+ * given object, and the number of entries as the length of the view.
+ *
+ * This function is used for @p const references to objects of LAPACKFullMatrix
+ * type because they contain immutable elements. Consequently, the return type
+ * of this function is a view to a set of @p const objects.
+ *
+ * @param[in] table The LAPACKFullMatrix for which we want to have an array view
+ * object. The array view corresponds to the <em>entire</em> object.
+ *
+ * @relates ArrayView
+ */
+template <typename ElementType>
+inline
+ArrayView<const ElementType>
+make_array_view (const LAPACKFullMatrix<ElementType> &matrix)
+{
+  return ArrayView<const ElementType> (&matrix(0,0), matrix.n_elements());
+}
+
 
 /**
  * Create a view to an entire row of a Table<2> object. This is equivalent to
diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h
index 2790c9e8b0..6854cef75b 100644
--- a/include/deal.II/base/mpi.h
+++ b/include/deal.II/base/mpi.h
@@ -154,16 +154,16 @@ namespace Utilities
 
     /**
      * Like the previous function, but take the sums over the elements of an
-     * array of length N. In other words, the i-th element of the results
+     * array of type T. In other words, the i-th element of the results
      * array is the sum over the i-th entries of the input arrays from each
-     * processor.
+     * processor. T and U must decay to the same type.
      *
      * Input and output arrays may be the same.
      */
-    template <typename T, unsigned int N>
-    void sum (const T (&values)[N],
+    template <typename T, typename U>
+    void sum (const T        &values,
               const MPI_Comm &mpi_communicator,
-              T (&sums)[N]);
+              U              &sums);
 
     /**
      * Like the previous function, but take the sums over the elements of an
@@ -179,48 +179,6 @@ namespace Utilities
               const MPI_Comm           &mpi_communicator,
               const ArrayView<T>       &sums);
 
-    /**
-     * Like the previous function, but take the sums over the elements of a
-     * std::vector. In other words, the i-th element of the results array is
-     * the sum over the i-th entries of the input arrays from each processor.
-     *
-     * Input and output vectors may be the same.
-     */
-    template <typename T>
-    void sum (const std::vector<T> &values,
-              const MPI_Comm &mpi_communicator,
-              std::vector<T> &sums);
-
-    /**
-     * Like the previous function, but take the sums over the elements of a
-     * Vector<T>.
-     *
-     * Input and output vectors may be the same.
-     */
-    template <typename T>
-    void sum (const Vector<T> &values,
-              const MPI_Comm &mpi_communicator,
-              Vector<T> &sums);
-
-    /**
-     * Like the previous function, but take the sums over the elements of a
-     * FullMatrix<T>.
-     *
-     * Input and output matrices may be the same.
-     */
-    template <typename T>
-    void sum (const FullMatrix<T> &values,
-              const MPI_Comm &mpi_communicator,
-              FullMatrix<T> &sums);
-
-    /**
-     * Same as above but for LAPACKFullMatrix.
-     */
-    template <typename T>
-    void sum (const LAPACKFullMatrix<T> &values,
-              const MPI_Comm &mpi_communicator,
-              LAPACKFullMatrix<T> &sums);
-
     /**
      * Perform an MPI sum of the entries of a symmetric tensor.
      *
@@ -265,30 +223,31 @@ namespace Utilities
            const MPI_Comm &mpi_communicator);
 
     /**
-     * Like the previous function, but take the maxima over the elements of an
-     * array of length N. In other words, the i-th element of the results
-     * array is the maximum of the i-th entries of the input arrays from each
-     * processor.
+     * Like the previous function, but take the maximum over the elements of an
+     * array of type T. In other words, the i-th element of the results array is
+     * the maximum over the i-th entries of the input arrays from each
+     * processor. T and U must decay to the same type.
      *
-     * Input and output arrays may be the same.
+     * Input and output vectors may be the same.
      */
-    template <typename T, unsigned int N>
-    void max (const T (&values)[N],
+    template <typename T, typename U>
+    void max (const T        &values,
               const MPI_Comm &mpi_communicator,
-              T (&maxima)[N]);
+              U              &maxima);
 
     /**
-     * Like the previous function, but take the maximum over the elements of a
-     * std::vector. In other words, the i-th element of the results array is
-     * the maximum over the i-th entries of the input arrays from each
+     * Like the previous function, but take the maximum over the elements of an
+     * array as specified by the ArrayView arguments.
+     * In other words, the i-th element of the results
+     * array is the maximum over the i-th entries of the input arrays from each
      * processor.
      *
-     * Input and output vectors may be the same.
+     * Input and output arrays may be the same.
      */
     template <typename T>
-    void max (const std::vector<T> &values,
-              const MPI_Comm &mpi_communicator,
-              std::vector<T> &maxima);
+    void max (const ArrayView<const T> &values,
+              const MPI_Comm           &mpi_communicator,
+              const ArrayView<T>       &maxima);
 
     /**
      * Return the minimum over all processors of the value @p t. This function
@@ -315,30 +274,30 @@ namespace Utilities
 
     /**
      * Like the previous function, but take the minima over the elements of an
-     * array of length N. In other words, the i-th element of the results
+     * array of type T. In other words, the i-th element of the results
      * array is the minimum of the i-th entries of the input arrays from each
-     * processor.
+     * processor. T and U must decay to the same type.
      *
      * Input and output arrays may be the same.
      */
-    template <typename T, unsigned int N>
-    void min (const T (&values)[N],
+    template <typename T, typename U>
+    void min (const T        &values,
               const MPI_Comm &mpi_communicator,
-              T (&minima)[N]);
+              U              &minima);
 
     /**
-     * Like the previous function, but take the minimum over the elements of a
-     * std::vector. In other words, the i-th element of the results array is
-     * the minimum over the i-th entries of the input arrays from each
+     * Like the previous function, but take the minimum over the elements of an
+     * array as specified by the ArrayView arguments.
+     * In other words, the i-th element of the results
+     * array is the minimum over the i-th entries of the input arrays from each
      * processor.
      *
-     * Input and output vectors may be the same.
+     * Input and output arrays may be the same.
      */
     template <typename T>
-    void min (const std::vector<T> &values,
-              const MPI_Comm &mpi_communicator,
-              std::vector<T> &minima);
-
+    void min (const ArrayView<const T> &values,
+              const MPI_Comm           &mpi_communicator,
+              const ArrayView<T>       &minima);
 
     /**
      * A data structure to store the result of the min_max_avg() function.
diff --git a/include/deal.II/base/mpi.templates.h b/include/deal.II/base/mpi.templates.h
index db4b33f35c..02ee73f80b 100644
--- a/include/deal.II/base/mpi.templates.h
+++ b/include/deal.II/base/mpi.templates.h
@@ -336,13 +336,19 @@ namespace Utilities
 
 
 
-    template <typename T>
-    void sum (const std::vector<T> &values,
-              const MPI_Comm       &mpi_communicator,
-              std::vector<T>       &sums)
+    template <typename T, typename U>
+    void sum (const T        &values,
+              const MPI_Comm &mpi_communicator,
+              U              &sums)
     {
-      internal::all_reduce(MPI_SUM, ArrayView<const T>(values),
-                           mpi_communicator, ArrayView<T>(sums));
+      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!");
+      const auto array_view_values = make_array_view(values);
+      using const_type = ArrayView<const typename decltype(array_view_values)::value_type>;
+      sum(static_cast<const_type> (array_view_values),
+          mpi_communicator,
+          make_array_view(sums));
     }
 
 
@@ -357,44 +363,6 @@ namespace Utilities
 
 
 
-    template <typename T>
-    void sum (const Vector<T> &values,
-              const MPI_Comm &mpi_communicator,
-              Vector<T> &sums)
-    {
-      const auto &size = values.size();
-      internal::all_reduce(MPI_SUM, ArrayView<const T>(values.begin(), size),
-                           mpi_communicator, ArrayView<T>(sums.begin(), size));
-    }
-
-
-
-    template <typename T>
-    void sum (const FullMatrix<T> &values,
-              const MPI_Comm &mpi_communicator,
-              FullMatrix<T> &sums)
-    {
-      const auto size_values = values.n()*values.m();
-      const auto size_sums = sums.n()*sums.m();
-      internal::all_reduce(MPI_SUM, ArrayView<const T>(&values[0][0], size_values),
-                           mpi_communicator, ArrayView<T>(&sums[0][0], size_sums));
-    }
-
-
-
-    template <typename T>
-    void sum (const LAPACKFullMatrix<T> &values,
-              const MPI_Comm &mpi_communicator,
-              LAPACKFullMatrix<T> &sums)
-    {
-      const auto size_values = values.n()*values.m();
-      const auto size_sums = sums.n()*sums.m();
-      internal::all_reduce(MPI_SUM, ArrayView<const T>(&values(0,0), size_values),
-                           mpi_communicator, ArrayView<T>(&sums(0,0), size_sums));
-    }
-
-
-
     template <int rank, int dim, typename Number>
     Tensor<rank,dim,Number>
     sum (const Tensor<rank,dim,Number> &local,
@@ -453,13 +421,29 @@ namespace Utilities
 
 
 
+    template <typename T, typename U>
+    void max (const T         &values,
+              const MPI_Comm  &mpi_communicator,
+              U               &maxima)
+    {
+      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!");
+      const auto array_view_values = make_array_view(values);
+      using const_type = ArrayView<const typename decltype(array_view_values)::value_type>;
+      max(static_cast<const_type> (array_view_values),
+          mpi_communicator,
+          make_array_view(maxima));
+    }
+
+
+
     template <typename T>
-    void max (const std::vector<T> &values,
-              const MPI_Comm       &mpi_communicator,
-              std::vector<T>       &maxima)
+    void max (const ArrayView<const T> &values,
+              const MPI_Comm           &mpi_communicator,
+              const ArrayView<T>       &maxima)
     {
-      internal::all_reduce(MPI_MAX, ArrayView<const T>(values),
-                           mpi_communicator, ArrayView<T> (maxima));
+      internal::all_reduce(MPI_MAX, values, mpi_communicator, maxima);
     }
 
 
@@ -476,13 +460,29 @@ namespace Utilities
 
 
 
+    template <typename T, typename U>
+    void min (const T        &values,
+              const MPI_Comm &mpi_communicator,
+              U              &minima)
+    {
+      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!");
+      const auto array_view_values = make_array_view(values);
+      using const_type = ArrayView<const typename decltype(array_view_values)::value_type>;
+      min(static_cast<const_type> (array_view_values),
+          mpi_communicator,
+          make_array_view(minima));
+    }
+
+
+
     template <typename T>
-    void min (const std::vector<T> &values,
-              const MPI_Comm       &mpi_communicator,
-              std::vector<T>       &minima)
+    void min (const ArrayView<const T> &values,
+              const MPI_Comm           &mpi_communicator,
+              const ArrayView<T>       &minima)
     {
-      internal::all_reduce(MPI_MIN, ArrayView<const T>(values),
-                           mpi_communicator, ArrayView<T> (minima));
+      internal::all_reduce(MPI_MIN, values, mpi_communicator, minima);
     }
   } // end of namespace MPI
 } // end of namespace Utilities
diff --git a/source/base/mpi.inst.in b/source/base/mpi.inst.in
index c96e796efd..a6acf7774e 100644
--- a/source/base/mpi.inst.in
+++ b/source/base/mpi.inst.in
@@ -17,13 +17,13 @@
 for (S : MPI_SCALARS)
 {
     template
-    void sum<S> (const LAPACKFullMatrix<S> &, const MPI_Comm &, LAPACKFullMatrix<S> &);
+    void sum<LAPACKFullMatrix<S> > (const LAPACKFullMatrix<S> &, const MPI_Comm &, LAPACKFullMatrix<S> &);
 
     template
-    void sum<S> (const Vector<S> &, const MPI_Comm &, Vector<S> &);
+    void sum<Vector<S> > (const Vector<S> &, const MPI_Comm &, Vector<S> &);
 
     template
-    void sum<S> (const FullMatrix<S> &, const MPI_Comm &, FullMatrix<S> &);
+    void sum<FullMatrix<S> > (const FullMatrix<S> &, const MPI_Comm &, FullMatrix<S> &);
 
     template
     void sum<S> (const ArrayView<const S> &, const MPI_Comm &, const ArrayView<S> &);
@@ -32,19 +32,19 @@ for (S : MPI_SCALARS)
     S sum<S> (const S &, const MPI_Comm &);
 
     template
-    void sum<S> (const std::vector<S> &, const MPI_Comm &, std::vector<S> &);
+    void sum<std::vector<S> > (const std::vector<S> &, const MPI_Comm &, std::vector<S> &);
 
     template
     S max<S> (const S &, const MPI_Comm &);
 
     template
-    void max<S> (const std::vector<S> &, const MPI_Comm &, std::vector<S> &);
+    void max<std::vector<S> > (const std::vector<S> &, const MPI_Comm &, std::vector<S> &);
 
     template
     S min<S> (const S &, const MPI_Comm &);
 
     template
-    void min<S> (const std::vector<S> &, const MPI_Comm &, std::vector<S> &);
+    void min<std::vector<S> > (const std::vector<S> &, const MPI_Comm &, std::vector<S> &);
 
     // The fixed-length array (i.e., things declared like T(&values)[N])
     // versions of the functions above live in the header file mpi.h since the