From 4d4d66842306c45005e6e97e82412c69c4e3bc59 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 26 Sep 2011 18:26:56 +0000 Subject: [PATCH] Add array versions of the sum/max operators. git-svn-id: https://svn.dealii.org/trunk@24434 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/base/utilities.h | 68 +++++++++++++++++++ tests/mpi/collective_02_array.cc | 64 +++++++++++++++++ .../collective_02_array/ncpu_1/cmp/generic | 2 + .../collective_02_array/ncpu_10/cmp/generic | 2 + .../collective_02_array/ncpu_4/cmp/generic | 2 + tests/mpi/collective_03_array.cc | 64 +++++++++++++++++ .../collective_03_array/ncpu_1/cmp/generic | 2 + .../collective_03_array/ncpu_10/cmp/generic | 2 + .../collective_03_array/ncpu_4/cmp/generic | 2 + 9 files changed, 208 insertions(+) create mode 100644 tests/mpi/collective_02_array.cc create mode 100644 tests/mpi/collective_02_array/ncpu_1/cmp/generic create mode 100644 tests/mpi/collective_02_array/ncpu_10/cmp/generic create mode 100644 tests/mpi/collective_02_array/ncpu_4/cmp/generic create mode 100644 tests/mpi/collective_03_array.cc create mode 100644 tests/mpi/collective_03_array/ncpu_1/cmp/generic create mode 100644 tests/mpi/collective_03_array/ncpu_10/cmp/generic create mode 100644 tests/mpi/collective_03_array/ncpu_4/cmp/generic diff --git a/deal.II/include/deal.II/base/utilities.h b/deal.II/include/deal.II/base/utilities.h index 9b16ce4434..f6e48b48fa 100644 --- a/deal.II/include/deal.II/base/utilities.h +++ b/deal.II/include/deal.II/base/utilities.h @@ -383,6 +383,22 @@ namespace Utilities T sum (const T &t, const MPI_Comm &mpi_communicator); + /** + * 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 is the sum over + * the i-th entries of the input + * arrays from each processor. + */ + template + inline + void sum (const T (&values)[N], + const MPI_Comm &mpi_communicator, + T (&sums)[N]); + /** * Return the maximum over all processors of the value @p t. This function * is collective over all processors given in the communicator. If @@ -399,6 +415,22 @@ namespace Utilities T max (const T &t, 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. + */ + template + inline + void max (const T (&values)[N], + const MPI_Comm &mpi_communicator, + T (&maxima)[N]); + /** * Data structure to store the result of * min_max_avg(). @@ -1018,6 +1050,24 @@ namespace Utilities } + template + inline + void sum (const T (&values)[N], + const MPI_Comm &mpi_communicator, + T (&sums)[N]) + { +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + MPI_Allreduce (const_cast(static_cast(&values[0])), + &sums[0], N, internal::mpi_type_id(values), MPI_SUM, + mpi_communicator); +#else + (void)mpi_communicator; + for (unsigned int i=0; i inline T max (const T &t, @@ -1034,6 +1084,24 @@ namespace Utilities return t; #endif } + + + template + inline + void max (const T (&values)[N], + const MPI_Comm &mpi_communicator, + T (&maxima)[N]) + { +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + MPI_Allreduce (const_cast(static_cast(&values[0])), + &maxima[0], N, internal::mpi_type_id(values), MPI_MAX, + mpi_communicator); +#else + (void)mpi_communicator; + for (unsigned int i=0; i +#include +#include + +void test() +{ + unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + const unsigned int numprocs = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); + + unsigned int values[2] = { 1, 2 }; + unsigned int sums[2]; + Utilities::MPI::sum (values, + MPI_COMM_WORLD, + sums); + Assert (sums[0] == numprocs, ExcInternalError()); + Assert (sums[1] == 2*numprocs, ExcInternalError()); + + if (myid==0) + deallog << sums[0] << ' ' << sums[1] << std::endl; +} + + +int main(int argc, char *argv[]) +{ +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + Utilities::MPI::MPI_InitFinalize mpi (argc, argv); +#else + (void)argc; + (void)argv; + compile_time_error; + +#endif + + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + { + std::ofstream logfile(output_file_for_mpi("collective_02_array").c_str()); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + deallog.push("mpi"); + test(); + deallog.pop(); + } + else + test(); +} diff --git a/tests/mpi/collective_02_array/ncpu_1/cmp/generic b/tests/mpi/collective_02_array/ncpu_1/cmp/generic new file mode 100644 index 0000000000..6e10373d77 --- /dev/null +++ b/tests/mpi/collective_02_array/ncpu_1/cmp/generic @@ -0,0 +1,2 @@ + +DEAL:mpi::1 2 diff --git a/tests/mpi/collective_02_array/ncpu_10/cmp/generic b/tests/mpi/collective_02_array/ncpu_10/cmp/generic new file mode 100644 index 0000000000..97f11ae084 --- /dev/null +++ b/tests/mpi/collective_02_array/ncpu_10/cmp/generic @@ -0,0 +1,2 @@ + +DEAL:mpi::10 20 diff --git a/tests/mpi/collective_02_array/ncpu_4/cmp/generic b/tests/mpi/collective_02_array/ncpu_4/cmp/generic new file mode 100644 index 0000000000..ad500b96a3 --- /dev/null +++ b/tests/mpi/collective_02_array/ncpu_4/cmp/generic @@ -0,0 +1,2 @@ + +DEAL:mpi::4 8 diff --git a/tests/mpi/collective_03_array.cc b/tests/mpi/collective_03_array.cc new file mode 100644 index 0000000000..b41551427d --- /dev/null +++ b/tests/mpi/collective_03_array.cc @@ -0,0 +1,64 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2009, 2010, 2011 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//--------------------------------------------------------------------------- + + +// check Utilities::MPI::max() for arrays + +#include "../tests.h" +#include +#include +#include + +void test() +{ + unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + const unsigned int numprocs = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); + + unsigned int values[2] = { 1, 2 }; + unsigned int maxima[2]; + Utilities::MPI::max (values, + MPI_COMM_WORLD, + maxima); + Assert (maxima[0] == 1, ExcInternalError()); + Assert (maxima[1] == 2, ExcInternalError()); + + if (myid==0) + deallog << maxima[0] << ' ' << maxima[1] << std::endl; +} + + +int main(int argc, char *argv[]) +{ +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + Utilities::MPI::MPI_InitFinalize mpi (argc, argv); +#else + (void)argc; + (void)argv; + compile_time_error; + +#endif + + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + { + std::ofstream logfile(output_file_for_mpi("collective_03_array").c_str()); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + deallog.push("mpi"); + test(); + deallog.pop(); + } + else + test(); +} diff --git a/tests/mpi/collective_03_array/ncpu_1/cmp/generic b/tests/mpi/collective_03_array/ncpu_1/cmp/generic new file mode 100644 index 0000000000..6e10373d77 --- /dev/null +++ b/tests/mpi/collective_03_array/ncpu_1/cmp/generic @@ -0,0 +1,2 @@ + +DEAL:mpi::1 2 diff --git a/tests/mpi/collective_03_array/ncpu_10/cmp/generic b/tests/mpi/collective_03_array/ncpu_10/cmp/generic new file mode 100644 index 0000000000..6e10373d77 --- /dev/null +++ b/tests/mpi/collective_03_array/ncpu_10/cmp/generic @@ -0,0 +1,2 @@ + +DEAL:mpi::1 2 diff --git a/tests/mpi/collective_03_array/ncpu_4/cmp/generic b/tests/mpi/collective_03_array/ncpu_4/cmp/generic new file mode 100644 index 0000000000..6e10373d77 --- /dev/null +++ b/tests/mpi/collective_03_array/ncpu_4/cmp/generic @@ -0,0 +1,2 @@ + +DEAL:mpi::1 2 -- 2.39.5