From: Denis Davydov Date: Mon, 4 Sep 2017 10:28:23 +0000 (+0200) Subject: add Utilities::MPI::sum() for FullMatrix objects X-Git-Tag: v9.0.0-rc1~1128^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1f9a055f4701220200d1a1dba16f7da4ef37e106;p=dealii.git add Utilities::MPI::sum() for FullMatrix objects --- diff --git a/doc/news/changes/minor/20170904DenisDavydov b/doc/news/changes/minor/20170904DenisDavydov new file mode 100644 index 0000000000..c70f1987b1 --- /dev/null +++ b/doc/news/changes/minor/20170904DenisDavydov @@ -0,0 +1,3 @@ +New: add Utilities::MPI::sum() for FullMatrix objects. +
+(Denis Davydov, 2017/09/04) diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index 31c5f5d252..0dd3d11c5a 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -43,6 +43,8 @@ template class SymmetricTensor; //Forward type declaration to allow MPI sums over Vector type template class Vector; +//Forward type declaration to allow MPI sums over FullMatrix type +template class FullMatrix; namespace Utilities @@ -172,6 +174,16 @@ namespace Utilities const MPI_Comm &mpi_communicator, Vector &sums); + /** + * Like the previous function, but take the sums over the elements of a + * FullMatrix. + * + * Input and output matrices may be the same. + */ + template + void sum (const FullMatrix &values, + const MPI_Comm &mpi_communicator, + FullMatrix &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 87dcfed76d..ab961d1da4 100644 --- a/include/deal.II/base/mpi.templates.h +++ b/include/deal.II/base/mpi.templates.h @@ -22,6 +22,7 @@ #include #include #include +#include #include @@ -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 void all_reduce (const MPI_Op &mpi_op, const T *const values, @@ -123,6 +133,8 @@ namespace Utilities } } + + template void all_reduce (const MPI_Op &mpi_op, const std::complex *const values, @@ -162,6 +174,8 @@ namespace Utilities } } + + template T all_reduce (const MPI_Op &mpi_op, const T &t, @@ -172,6 +186,8 @@ namespace Utilities return output; } + + template void all_reduce (const MPI_Op &mpi_op, const std::vector &values, @@ -183,6 +199,8 @@ namespace Utilities all_reduce(mpi_op, &values[0], mpi_communicator, &output[0], values.size()); } + + template void all_reduce (const MPI_Op &mpi_op, const Vector &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 + void all_reduce (const MPI_Op &mpi_op, + const FullMatrix &values, + const MPI_Comm &mpi_communicator, + FullMatrix &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 T sum (const T &t, const MPI_Comm &mpi_communicator) @@ -204,6 +239,7 @@ namespace Utilities } + template void sum (const std::vector &values, const MPI_Comm &mpi_communicator, @@ -212,6 +248,8 @@ namespace Utilities internal::all_reduce(MPI_SUM, values, mpi_communicator, sums); } + + template void sum (const Vector &values, const MPI_Comm &mpi_communicator, @@ -221,6 +259,17 @@ namespace Utilities } + + template + void sum (const FullMatrix &values, + const MPI_Comm &mpi_communicator, + FullMatrix &sums) + { + internal::all_reduce(MPI_SUM, values, mpi_communicator, sums); + } + + + template Tensor sum (const Tensor &local, @@ -242,6 +291,8 @@ namespace Utilities return global; } + + template SymmetricTensor sum (const SymmetricTensor &local, @@ -263,6 +314,8 @@ namespace Utilities return global; } + + template T max (const T &t, const MPI_Comm &mpi_communicator) @@ -271,6 +324,7 @@ namespace Utilities } + template void max (const std::vector &values, const MPI_Comm &mpi_communicator, @@ -280,6 +334,7 @@ namespace Utilities } + template T min (const T &t, const MPI_Comm &mpi_communicator) @@ -288,6 +343,7 @@ namespace Utilities } + template void min (const std::vector &values, const MPI_Comm &mpi_communicator, diff --git a/source/base/mpi.inst.in b/source/base/mpi.inst.in index 2420c84fea..56cb0c2850 100644 --- a/source/base/mpi.inst.in +++ b/source/base/mpi.inst.in @@ -19,6 +19,9 @@ for (S : MPI_SCALARS) template void sum (const Vector &, const MPI_Comm &, Vector &); + template + void sum (const FullMatrix &, const MPI_Comm &, FullMatrix &); + template S sum (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 index 0000000000..3d8730d27a --- /dev/null +++ b/tests/mpi/collective_full_matrix.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 FullMatrix 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()); + + FullMatrix 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 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(); + deallog.pop(); + deallog.push("double"); + test(); + deallog.pop(); + } + else + { + test(); + test(); + } + +} diff --git a/tests/mpi/collective_full_matrix.mpirun=1.output b/tests/mpi/collective_full_matrix.mpirun=1.output new file mode 100644 index 0000000000..1d1ebcece8 --- /dev/null +++ b/tests/mpi/collective_full_matrix.mpirun=1.output @@ -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 index 0000000000..1d1ebcece8 --- /dev/null +++ b/tests/mpi/collective_full_matrix.mpirun=10.output @@ -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 index 0000000000..1d1ebcece8 --- /dev/null +++ b/tests/mpi/collective_full_matrix.mpirun=3.output @@ -0,0 +1,3 @@ + +DEAL:float::Ok +DEAL:double::Ok