From: Denis Davydov Date: Sat, 30 Sep 2017 21:34:17 +0000 (+0200) Subject: add Utilities::MPI::sum() for LAPACKFullMatrix objects X-Git-Tag: v9.0.0-rc1~1000^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=92ac8f4855c9528064e6e855203b1bb165e42ba9;p=dealii.git add Utilities::MPI::sum() for LAPACKFullMatrix objects --- diff --git a/doc/news/changes/minor/20170930DenisDavydov-2 b/doc/news/changes/minor/20170930DenisDavydov-2 new file mode 100644 index 0000000000..5e5ba2d94b --- /dev/null +++ b/doc/news/changes/minor/20170930DenisDavydov-2 @@ -0,0 +1,3 @@ +New: add Utilities::MPI::sum() for LAPACKFullMatrix objects. +
+(Denis Davydov, 2017/09/04) diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index dd8a512b5c..1c42e24a16 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -55,6 +55,8 @@ template class SymmetricTensor; template class Vector; //Forward type declaration to allow MPI sums over FullMatrix type template class FullMatrix; +//Forward type declaration to allow MPI sums over LAPACKFullMatrix type +template class LAPACKFullMatrix; namespace Utilities @@ -195,6 +197,14 @@ namespace Utilities const MPI_Comm &mpi_communicator, FullMatrix &sums); + /** + * Same as above but for LAPACKFullMatrix. + */ + template + void sum (const LAPACKFullMatrix &values, + const MPI_Comm &mpi_communicator, + LAPACKFullMatrix &sums); + /** * Perform an MPI sum of the entries of a symmetric tensor. * diff --git a/include/deal.II/base/mpi.templates.h b/include/deal.II/base/mpi.templates.h index ab961d1da4..f316f2de85 100644 --- a/include/deal.II/base/mpi.templates.h +++ b/include/deal.II/base/mpi.templates.h @@ -23,6 +23,7 @@ #include #include #include +#include #include @@ -227,6 +228,21 @@ namespace Utilities all_reduce(mpi_op, &values[0][0], mpi_communicator, &output[0][0], values.m() * values.n()); } + + + template + void all_reduce (const MPI_Op &mpi_op, + const LAPACKFullMatrix &values, + const MPI_Comm &mpi_communicator, + LAPACKFullMatrix &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 + void sum (const LAPACKFullMatrix &values, + const MPI_Comm &mpi_communicator, + LAPACKFullMatrix &sums) + { + internal::all_reduce(MPI_SUM, values, mpi_communicator, sums); + } + + + template Tensor sum (const Tensor &local, diff --git a/source/base/mpi.inst.in b/source/base/mpi.inst.in index d642d66412..504a385c58 100644 --- a/source/base/mpi.inst.in +++ b/source/base/mpi.inst.in @@ -13,6 +13,11 @@ // // --------------------------------------------------------------------- +for (S : REAL_SCALARS) +{ + template + void sum (const LAPACKFullMatrix &, const MPI_Comm &, LAPACKFullMatrix &); +} 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 index 0000000000..5394431f58 --- /dev/null +++ b/tests/mpi/collective_full_matrix_02.cc @@ -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 +#include + +template +void test(const unsigned int m = 13, const unsigned int n = 5) +{ + Assert( Utilities::MPI::job_supports_mpi(), ExcInternalError()); + + LAPACKFullMatrix 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 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(); + deallog.pop(); + deallog.push("double"); + test(); + deallog.pop(); + } + else + { + test(); + test(); + } + +} 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 index 0000000000..1d1ebcece8 --- /dev/null +++ b/tests/mpi/collective_full_matrix_02.mpirun=1.output @@ -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 index 0000000000..1d1ebcece8 --- /dev/null +++ b/tests/mpi/collective_full_matrix_02.mpirun=3.output @@ -0,0 +1,3 @@ + +DEAL:float::Ok +DEAL:double::Ok