]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add Utilities::MPI::sum() for FullMatrix objects
authorDenis Davydov <davydden@gmail.com>
Mon, 4 Sep 2017 10:28:23 +0000 (12:28 +0200)
committerDenis Davydov <davydden@gmail.com>
Mon, 4 Sep 2017 10:28:48 +0000 (12:28 +0200)
doc/news/changes/minor/20170904DenisDavydov [new file with mode: 0644]
include/deal.II/base/mpi.h
include/deal.II/base/mpi.templates.h
source/base/mpi.inst.in
tests/mpi/collective_full_matrix.cc [new file with mode: 0644]
tests/mpi/collective_full_matrix.mpirun=1.output [new file with mode: 0644]
tests/mpi/collective_full_matrix.mpirun=10.output [new file with mode: 0644]
tests/mpi/collective_full_matrix.mpirun=3.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20170904DenisDavydov b/doc/news/changes/minor/20170904DenisDavydov
new file mode 100644 (file)
index 0000000..c70f198
--- /dev/null
@@ -0,0 +1,3 @@
+New: add Utilities::MPI::sum() for FullMatrix objects.
+<br>
+(Denis Davydov, 2017/09/04)
index 31c5f5d2522699a927b72d340eb4fde5af593d38..0dd3d11c5ae57c0df06f9f7cdbafb4327622868f 100644 (file)
@@ -43,6 +43,8 @@ template <int rank, int dim, typename Number> class SymmetricTensor;
 
 //Forward type declaration to allow MPI sums over Vector<number> type
 template <typename Number> class Vector;
+//Forward type declaration to allow MPI sums over FullMatrix<number> type
+template <typename Number> class FullMatrix;
 
 
 namespace Utilities
@@ -172,6 +174,16 @@ namespace Utilities
               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);
 
     /**
      * Perform an MPI sum of the entries of a symmetric tensor.
index 87dcfed76d456453a43da0e860ec61c0f7ef5625..ab961d1da4f781a478ee2d0f0fc68b09baf560f7 100644 (file)
@@ -22,6 +22,7 @@
 #include <deal.II/base/tensor.h>
 #include <deal.II/base/symmetric_tensor.h>
 #include <deal.II/lac/vector.h>
+#include <deal.II/lac/full_matrix.h>
 
 #include <vector>
 
@@ -43,48 +44,57 @@ namespace Utilities
       }
 
 
+
       inline MPI_Datatype mpi_type_id (const long int *)
       {
         return MPI_LONG;
       }
 
 
+
       inline MPI_Datatype mpi_type_id (const unsigned int *)
       {
         return MPI_UNSIGNED;
       }
 
 
+
       inline MPI_Datatype mpi_type_id (const unsigned long int *)
       {
         return MPI_UNSIGNED_LONG;
       }
 
 
+
       inline MPI_Datatype mpi_type_id (const unsigned long long int *)
       {
         return MPI_UNSIGNED_LONG_LONG;
       }
 
 
+
       inline MPI_Datatype mpi_type_id (const float *)
       {
         return MPI_FLOAT;
       }
 
 
+
       inline MPI_Datatype mpi_type_id (const double *)
       {
         return MPI_DOUBLE;
       }
 
 
+
       inline MPI_Datatype mpi_type_id (const long double *)
       {
         return MPI_LONG_DOUBLE;
       }
 #endif
 
+
+
       template <typename T>
       void all_reduce (const MPI_Op      &mpi_op,
                        const T *const    values,
@@ -123,6 +133,8 @@ namespace Utilities
           }
       }
 
+
+
       template <typename T>
       void all_reduce (const MPI_Op                   &mpi_op,
                        const std::complex<T> *const    values,
@@ -162,6 +174,8 @@ namespace Utilities
           }
       }
 
+
+
       template <typename T>
       T all_reduce (const MPI_Op &mpi_op,
                     const T &t,
@@ -172,6 +186,8 @@ namespace Utilities
         return output;
       }
 
+
+
       template <typename T>
       void all_reduce (const MPI_Op &mpi_op,
                        const std::vector<T> &values,
@@ -183,6 +199,8 @@ namespace Utilities
         all_reduce(mpi_op, &values[0], mpi_communicator, &output[0], values.size());
       }
 
+
+
       template <typename T>
       void all_reduce (const MPI_Op    &mpi_op,
                        const Vector<T> &values,
@@ -193,9 +211,26 @@ namespace Utilities
                ExcDimensionMismatch(values.size(), output.size()));
         all_reduce(mpi_op, values.begin(), mpi_communicator, output.begin(), values.size());
       }
+
+
+
+      template <typename T>
+      void all_reduce (const MPI_Op    &mpi_op,
+                       const FullMatrix<T> &values,
+                       const MPI_Comm  &mpi_communicator,
+                       FullMatrix<T>  &output)
+      {
+        Assert(values.m() == output.m(),
+               ExcDimensionMismatch(values.m(), output.m()));
+        Assert(values.n() == output.n(),
+               ExcDimensionMismatch(values.n(), output.n()));
+        all_reduce(mpi_op, &values[0][0], mpi_communicator, &output[0][0], values.m() * values.n());
+      }
+
     }
 
 
+
     template <typename T>
     T sum (const T &t,
            const MPI_Comm &mpi_communicator)
@@ -204,6 +239,7 @@ namespace Utilities
     }
 
 
+
     template <typename T>
     void sum (const std::vector<T> &values,
               const MPI_Comm       &mpi_communicator,
@@ -212,6 +248,8 @@ namespace Utilities
       internal::all_reduce(MPI_SUM, values, mpi_communicator, sums);
     }
 
+
+
     template <typename T>
     void sum (const Vector<T> &values,
               const MPI_Comm &mpi_communicator,
@@ -221,6 +259,17 @@ namespace Utilities
     }
 
 
+
+    template <typename T>
+    void sum (const FullMatrix<T> &values,
+              const MPI_Comm &mpi_communicator,
+              FullMatrix<T> &sums)
+    {
+      internal::all_reduce(MPI_SUM, values, mpi_communicator, sums);
+    }
+
+
+
     template <int rank, int dim, typename Number>
     Tensor<rank,dim,Number>
     sum (const Tensor<rank,dim,Number> &local,
@@ -242,6 +291,8 @@ namespace Utilities
       return global;
     }
 
+
+
     template <int rank, int dim, typename Number>
     SymmetricTensor<rank,dim,Number>
     sum (const SymmetricTensor<rank,dim,Number> &local,
@@ -263,6 +314,8 @@ namespace Utilities
       return global;
     }
 
+
+
     template <typename T>
     T max (const T &t,
            const MPI_Comm &mpi_communicator)
@@ -271,6 +324,7 @@ namespace Utilities
     }
 
 
+
     template <typename T>
     void max (const std::vector<T> &values,
               const MPI_Comm       &mpi_communicator,
@@ -280,6 +334,7 @@ namespace Utilities
     }
 
 
+
     template <typename T>
     T min (const T &t,
            const MPI_Comm &mpi_communicator)
@@ -288,6 +343,7 @@ namespace Utilities
     }
 
 
+
     template <typename T>
     void min (const std::vector<T> &values,
               const MPI_Comm       &mpi_communicator,
index 2420c84fea2bf279602e1ce6ca7c35d35ef59b7c..56cb0c2850db4cf7271cf21bc6bbba1c8c7a1d96 100644 (file)
@@ -19,6 +19,9 @@ for (S : MPI_SCALARS)
     template
     void sum<S> (const Vector<S> &, const MPI_Comm &, Vector<S> &);
 
+    template
+    void sum<S> (const FullMatrix<S> &, const MPI_Comm &, FullMatrix<S> &);
+
     template
     S sum<S> (const S &, const MPI_Comm &);
 
diff --git a/tests/mpi/collective_full_matrix.cc b/tests/mpi/collective_full_matrix.cc
new file mode 100644 (file)
index 0000000..3d8730d
--- /dev/null
@@ -0,0 +1,73 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+// check Utilities::MPI::sum() for FullMatrix objects
+
+#include "../tests.h"
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+
+template <typename NumberType>
+void test(const unsigned int m = 13, const unsigned int n = 5)
+{
+  Assert( Utilities::MPI::job_supports_mpi(), ExcInternalError());
+
+  FullMatrix<NumberType> full_matrix(m,n);
+  {
+    unsigned int index = 0;
+    for (unsigned int i = 0; i < full_matrix.m(); ++i)
+      for (unsigned int j = 0; j < full_matrix.n(); ++j)
+        full_matrix(i,j) = index++;
+  }
+
+  FullMatrix<NumberType> full_matrix_original(m,n);
+  full_matrix_original = full_matrix;
+
+  // inplace
+  Utilities::MPI::sum(full_matrix, MPI_COMM_WORLD, full_matrix);
+
+  const unsigned int numprocs = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD);
+  for (unsigned int i = 0; i < full_matrix.m(); ++i)
+    for (unsigned int j = 0; j < full_matrix.n(); ++j)
+      Assert (full_matrix(i,j) == full_matrix_original(i,j) * numprocs, ExcInternalError());
+
+  if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0)
+    deallog << "Ok" << std::endl;
+}
+
+
+int main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads());
+
+  if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0)
+    {
+      initlog();
+      deallog.push("float");
+      test<float>();
+      deallog.pop();
+      deallog.push("double");
+      test<double>();
+      deallog.pop();
+    }
+  else
+    {
+      test<float>();
+      test<double>();
+    }
+
+}
diff --git a/tests/mpi/collective_full_matrix.mpirun=1.output b/tests/mpi/collective_full_matrix.mpirun=1.output
new file mode 100644 (file)
index 0000000..1d1ebce
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:float::Ok
+DEAL:double::Ok
diff --git a/tests/mpi/collective_full_matrix.mpirun=10.output b/tests/mpi/collective_full_matrix.mpirun=10.output
new file mode 100644 (file)
index 0000000..1d1ebce
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:float::Ok
+DEAL:double::Ok
diff --git a/tests/mpi/collective_full_matrix.mpirun=3.output b/tests/mpi/collective_full_matrix.mpirun=3.output
new file mode 100644 (file)
index 0000000..1d1ebce
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:float::Ok
+DEAL:double::Ok

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.