]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add Utilities::MPI::sum() for LAPACKFullMatrix objects
authorDenis Davydov <davydden@gmail.com>
Sat, 30 Sep 2017 21:34:17 +0000 (23:34 +0200)
committerDenis Davydov <davydden@gmail.com>
Tue, 3 Oct 2017 05:58:07 +0000 (07:58 +0200)
doc/news/changes/minor/20170930DenisDavydov-2 [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_02.cc [new file with mode: 0644]
tests/mpi/collective_full_matrix_02.mpirun=1.output [new file with mode: 0644]
tests/mpi/collective_full_matrix_02.mpirun=3.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20170930DenisDavydov-2 b/doc/news/changes/minor/20170930DenisDavydov-2
new file mode 100644 (file)
index 0000000..5e5ba2d
--- /dev/null
@@ -0,0 +1,3 @@
+New: add Utilities::MPI::sum() for LAPACKFullMatrix objects.
+<br>
+(Denis Davydov, 2017/09/04)
index dd8a512b5c2e1294639d9121b91dea50648b6b2a..1c42e24a16c7625ae09b436c299608d76b6529bc 100644 (file)
@@ -55,6 +55,8 @@ template <int rank, int dim, typename Number> class SymmetricTensor;
 template <typename Number> class Vector;
 //Forward type declaration to allow MPI sums over FullMatrix<number> type
 template <typename Number> class FullMatrix;
+//Forward type declaration to allow MPI sums over LAPACKFullMatrix<number> type
+template <typename Number> class LAPACKFullMatrix;
 
 
 namespace Utilities
@@ -195,6 +197,14 @@ namespace Utilities
               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.
      *
index ab961d1da4f781a478ee2d0f0fc68b09baf560f7..f316f2de851492d2f0136c58c32e68eb6deb8033 100644 (file)
@@ -23,6 +23,7 @@
 #include <deal.II/base/symmetric_tensor.h>
 #include <deal.II/lac/vector.h>
 #include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/lapack_full_matrix.h>
 
 #include <vector>
 
@@ -227,6 +228,21 @@ namespace Utilities
         all_reduce(mpi_op, &values[0][0], mpi_communicator, &output[0][0], values.m() * values.n());
       }
 
+
+
+      template <typename T>
+      void all_reduce (const MPI_Op    &mpi_op,
+                       const LAPACKFullMatrix<T> &values,
+                       const MPI_Comm  &mpi_communicator,
+                       LAPACKFullMatrix<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());
+      }
+
     }
 
 
@@ -270,6 +286,16 @@ namespace Utilities
 
 
 
+    template <typename T>
+    void sum (const LAPACKFullMatrix<T> &values,
+              const MPI_Comm &mpi_communicator,
+              LAPACKFullMatrix<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,
index d642d664127f83783d082e7ade05e4d92a887c50..504a385c589f0ad353bce9041ea78ef49beb12d5 100644 (file)
 //
 // ---------------------------------------------------------------------
 
+for (S : REAL_SCALARS)
+{
+    template
+    void sum<S> (const LAPACKFullMatrix<S> &, const MPI_Comm &, LAPACKFullMatrix<S> &);
+}
 
 for (S : MPI_SCALARS)
 {
diff --git a/tests/mpi/collective_full_matrix_02.cc b/tests/mpi/collective_full_matrix_02.cc
new file mode 100644 (file)
index 0000000..5394431
--- /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 LAPACKFullMatrix objects
+
+#include "../tests.h"
+#include <deal.II/lac/lapack_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());
+
+  LAPACKFullMatrix<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++;
+  }
+
+  LAPACKFullMatrix<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) * double(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_02.mpirun=1.output b/tests/mpi/collective_full_matrix_02.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_02.mpirun=3.output b/tests/mpi/collective_full_matrix_02.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.