From: bangerth Date: Mon, 26 Sep 2011 17:54:20 +0000 (+0000) Subject: Add Utilities::MPI::max. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=aa1d857adb8c7ebfe404370741de01d0a5832780;p=dealii-svn.git Add Utilities::MPI::max. git-svn-id: https://svn.dealii.org/trunk@24433 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index e030ab9e50..5ea9804271 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -336,10 +336,8 @@ its "local" index or numbers::invalid_unsigned_int.
(Guido Kanschat, 2011/09/26) -
  • New: There is now a function GridTools::volume() computing the volume - -
  • New: The function Utilities::MPI::sum -function can be used to compute the sum of values over a number of +
  • New: The functions Utilities::MPI::sum and Utilities::MPI::max +function can be used to compute the sum or maximum of values over a number of MPI processes without having to deal with the underlying MPI functions.
    (Wolfgang Bangerth, 2011/09/25) @@ -366,7 +364,7 @@ threads kept pouncing on it. This is now fixed.
    (Patrick Sodré, 2011/09/07) -
  • New: There is now a function GridTools::volume computing the volume +
  • New: There is now a function GridTools::volume() computing the volume of a triangulation.
    (Wolfgang Bangerth, 2011/09/02) diff --git a/deal.II/include/deal.II/base/utilities.h b/deal.II/include/deal.II/base/utilities.h index 803d62e78f..9b16ce4434 100644 --- a/deal.II/include/deal.II/base/utilities.h +++ b/deal.II/include/deal.II/base/utilities.h @@ -381,6 +381,22 @@ namespace Utilities */ template T sum (const T &t, + const MPI_Comm &mpi_communicator); + + /** + * Return the maximum over all processors of the value @p t. This function + * is collective over all processors given in the communicator. If + * deal.II is not configured for use of MPI, this function simply + * returns the value of @p t. This function corresponds to the + * MPI_Allreduce function, i.e. all processors receive + * the result of this operation. + * + * @note This function is only implemented for certain template + * arguments T, namely float, double, int, + * unsigned int. + */ + template + T max (const T &t, const MPI_Comm &mpi_communicator); /** @@ -983,6 +999,7 @@ namespace Utilities #endif } + template inline T sum (const T &t, @@ -999,6 +1016,24 @@ namespace Utilities return t; #endif } + + + template + inline + T max (const T &t, + const MPI_Comm &mpi_communicator) + { +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + T sum; + MPI_Allreduce (const_cast(static_cast(&t)), + &sum, 1, internal::mpi_type_id(&t), MPI_MAX, + mpi_communicator); + return sum; +#else + (void)mpi_communicator; + return t; +#endif + } } } diff --git a/tests/mpi/collective_02.cc b/tests/mpi/collective_02.cc index d053862530..4f1cc1ea93 100644 --- a/tests/mpi/collective_02.cc +++ b/tests/mpi/collective_02.cc @@ -12,7 +12,7 @@ //--------------------------------------------------------------------------- -// check Utilities::System::calculate_collective_mpi_sum() +// check Utilities::MPI::sum() #include "../tests.h" #include diff --git a/tests/mpi/collective_03.cc b/tests/mpi/collective_03.cc new file mode 100644 index 0000000000..89a6bee8a7 --- /dev/null +++ b/tests/mpi/collective_03.cc @@ -0,0 +1,74 @@ +//--------------------------------------------------------------------------- +// $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() + +#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); + + int int_sum, uint_sum, double_sum, float_sum; + int_sum + = + Utilities::MPI::max(myid+1, + MPI_COMM_WORLD); + uint_sum + = + Utilities::MPI::max(myid+1, + MPI_COMM_WORLD); + float_sum + = + Utilities::MPI::max(myid+1, + MPI_COMM_WORLD); + double_sum + = + Utilities::MPI::max(myid+1, + MPI_COMM_WORLD); + + if (myid==0) + deallog << int_sum << ' ' << uint_sum << ' ' << double_sum << ' ' << float_sum << 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").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/ncpu_1/cmp/generic b/tests/mpi/collective_03/ncpu_1/cmp/generic new file mode 100644 index 0000000000..4723d36fe0 --- /dev/null +++ b/tests/mpi/collective_03/ncpu_1/cmp/generic @@ -0,0 +1,2 @@ + +DEAL:mpi::1 1 1 1 diff --git a/tests/mpi/collective_03/ncpu_10/cmp/generic b/tests/mpi/collective_03/ncpu_10/cmp/generic new file mode 100644 index 0000000000..c72ae8e67c --- /dev/null +++ b/tests/mpi/collective_03/ncpu_10/cmp/generic @@ -0,0 +1,2 @@ + +DEAL:mpi::10 10 10 10 diff --git a/tests/mpi/collective_03/ncpu_4/cmp/generic b/tests/mpi/collective_03/ncpu_4/cmp/generic new file mode 100644 index 0000000000..b0dbea0c21 --- /dev/null +++ b/tests/mpi/collective_03/ncpu_4/cmp/generic @@ -0,0 +1,2 @@ + +DEAL:mpi::4 4 4 4