From 6452dc3b3b0ba35225aaef25e2568b89c614162b Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Tue, 13 Jan 2015 22:39:18 +0100 Subject: [PATCH] Check for job_supports_mpi() inside MPI Allgather functions. MPI calls are only allowed if MPI_Init has been called. To simplify calls to MPI::sum, MPI::max etc. and avoid putting job_supports_mpi() around them, place this hopefully very cheap test inside the functions. Also moved job_supports_mpi() from the namespace Utilities::System to Utilities::MPI. --- doc/news/changes.h | 15 ++ include/deal.II/base/mpi.h | 177 ++++++++++++------ include/deal.II/base/utilities.h | 14 +- .../deal.II/lac/parallel_vector.templates.h | 2 +- .../matrix_free/matrix_free.templates.h | 25 +-- source/base/mpi.cc | 22 ++- source/base/partitioner.cc | 1 - source/base/timer.cc | 50 ++--- source/base/utilities.cc | 11 +- tests/mpi/collective_01.cc | 2 +- tests/mpi/distribute_sp_01.cc | 2 +- tests/mpi/distribute_sp_02.cc | 2 +- tests/mpi/point_to_point_pattern_01.cc | 2 +- tests/mpi/simple_mpi_01.cc | 2 +- 14 files changed, 179 insertions(+), 148 deletions(-) diff --git a/doc/news/changes.h b/doc/news/changes.h index 2feac5d682..bceba1927b 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -202,6 +202,21 @@ inconvenience this causes.

Specific improvements

    +
  1. Improved: MPI collective operations such as MPI::sum, MPI::max now + check for job_supports_mpi() internally, which allows running them also + without a call to MPI_Init. +
    + (Martin Kronbichler, 2015/01/13) +
  2. + +
  3. Changed: The method job_supports_mpi() now resides in the namespace + Utilities::MPI instead of Utilities::System for consistency with other MPI + methods. The old method has been marked deprecated and will be removed in + a future version. +
    + (Martin Kronbichler, 2015/01/13) +
  4. +
  5. Fixed: The update of ghost values in parallel::distributed::Vector when calling the assignment operator is now active when one of the two vector had its ghost values updated before or when the layout of the right hand side diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index 24a3ad6cfe..25934a1dcb 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2011 - 2014 by the deal.II authors +// Copyright (C) 2011 - 2015 by the deal.II authors // // This file is part of the deal.II library. // @@ -316,6 +316,19 @@ namespace Utilities ~MPI_InitFinalize(); }; + /** + * Return whether (i) deal.II has been compiled to support MPI (for + * example by compiling with CXX=mpiCC) and if so whether + * (ii) MPI_Init() has been called (for example using the + * Utilities::MPI::MPI_InitFinalize class). In other words, the result + * indicates whether the current job is running under MPI. + * + * @note The function does not take into account whether an MPI job + * actually runs on more than one processor or is, in fact, a single-node + * job that happens to run under MPI. + */ + bool job_supports_mpi (); + namespace internal { #ifdef DEAL_II_WITH_MPI @@ -378,15 +391,20 @@ namespace Utilities const MPI_Comm &mpi_communicator) { #ifdef DEAL_II_WITH_MPI - T sum; - MPI_Allreduce (const_cast(static_cast(&t)), - &sum, 1, internal::mpi_type_id(&t), MPI_SUM, - mpi_communicator); - return sum; -#else - (void)mpi_communicator; - return t; + if (job_supports_mpi()) + { + T sum; + MPI_Allreduce (const_cast(static_cast(&t)), + &sum, 1, internal::mpi_type_id(&t), MPI_SUM, + mpi_communicator); + return sum; + } + else #endif + { + (void)mpi_communicator; + return t; + } } @@ -397,18 +415,23 @@ namespace Utilities T (&sums)[N]) { #ifdef DEAL_II_WITH_MPI - MPI_Allreduce ((&values[0] != &sums[0] - ? - const_cast(static_cast(&values[0])) - : - MPI_IN_PLACE), - &sums[0], N, internal::mpi_type_id(values), MPI_SUM, - mpi_communicator); -#else - (void)mpi_communicator; - for (unsigned int i=0; i(static_cast(&values[0])) + : + MPI_IN_PLACE), + &sums[0], N, internal::mpi_type_id(values), MPI_SUM, + mpi_communicator); + } + else #endif + { + (void)mpi_communicator; + for (unsigned int i=0; i &sums) { #ifdef DEAL_II_WITH_MPI - sums.resize (values.size()); - MPI_Allreduce ((&values[0] != &sums[0] - ? - const_cast(static_cast(&values[0])) - : - MPI_IN_PLACE), - &sums[0], values.size(), internal::mpi_type_id((T *)0), MPI_SUM, - mpi_communicator); -#else - (void)mpi_communicator; - sums = values; + if (job_supports_mpi()) + { + sums.resize (values.size()); + MPI_Allreduce ((&values[0] != &sums[0] + ? + const_cast(static_cast(&values[0])) + : + MPI_IN_PLACE), + &sums[0], values.size(), internal::mpi_type_id((T *)0), MPI_SUM, + mpi_communicator); + } + else #endif + { + (void)mpi_communicator; + sums = values; + } } @@ -440,15 +468,20 @@ namespace Utilities const MPI_Comm &mpi_communicator) { #ifdef DEAL_II_WITH_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; + if (job_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 #endif + { + (void)mpi_communicator; + return t; + } } @@ -459,18 +492,23 @@ namespace Utilities T (&maxima)[N]) { #ifdef DEAL_II_WITH_MPI - MPI_Allreduce ((&values[0] != &maxima[0] - ? - const_cast(static_cast(&values[0])) - : - MPI_IN_PLACE), - &maxima[0], N, internal::mpi_type_id(values), MPI_MAX, - mpi_communicator); -#else - (void)mpi_communicator; - for (unsigned int i=0; i(static_cast(&values[0])) + : + MPI_IN_PLACE), + &maxima[0], N, internal::mpi_type_id(values), MPI_MAX, + mpi_communicator); + } + else #endif + { + (void)mpi_communicator; + for (unsigned int i=0; i &maxima) { #ifdef DEAL_II_WITH_MPI - maxima.resize (values.size()); - MPI_Allreduce ((&values[0] != &maxima[0] - ? - const_cast(static_cast(&values[0])) - : - MPI_IN_PLACE), - &maxima[0], values.size(), internal::mpi_type_id((T *)0), MPI_MAX, - mpi_communicator); + if (job_supports_mpi()) + { + maxima.resize (values.size()); + MPI_Allreduce ((&values[0] != &maxima[0] + ? + const_cast(static_cast(&values[0])) + : + MPI_IN_PLACE), + &maxima[0], values.size(), internal::mpi_type_id((T *)0), MPI_MAX, + mpi_communicator); + } + else +#endif + { + (void)mpi_communicator; + maxima = values; + } + } + + + inline + bool job_supports_mpi () + { +#ifdef DEAL_II_WITH_MPI + int MPI_has_been_started = 0; + MPI_Initialized(&MPI_has_been_started); + + return (MPI_has_been_started > 0); #else - (void)mpi_communicator; - maxima = values; + return false; #endif } } // end of namespace MPI diff --git a/include/deal.II/base/utilities.h b/include/deal.II/base/utilities.h index a005889de2..52f0f25b62 100644 --- a/include/deal.II/base/utilities.h +++ b/include/deal.II/base/utilities.h @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2005 - 2013 by the deal.II authors +// Copyright (C) 2005 - 2015 by the deal.II authors // // This file is part of the deal.II library. // @@ -342,17 +342,9 @@ namespace Utilities std::string get_time (); /** - * Return whether (i) deal.II has been compiled to support MPI (for - * example by compiling with CXX=mpiCC) and if so whether - * (ii) MPI_Init() has been called (for example using the - * Utilities::System::MPI_InitFinalize class). In other words, the result - * indicates whether the current job is running under MPI. - * - * @note The function does not take into account whether an MPI job - * actually runs on more than one processor or is, in fact, a single-node - * job that happens to run under MPI. + * @deprecated Use Utilities::MPI::job_supports_mpi() instead. */ - bool job_supports_mpi (); + bool job_supports_mpi () DEAL_II_DEPRECATED; } diff --git a/include/deal.II/lac/parallel_vector.templates.h b/include/deal.II/lac/parallel_vector.templates.h index 015d3833ae..2a7f3d1b73 100644 --- a/include/deal.II/lac/parallel_vector.templates.h +++ b/include/deal.II/lac/parallel_vector.templates.h @@ -600,7 +600,7 @@ namespace parallel #ifdef DEAL_II_WITH_MPI #ifdef DEBUG - if (Utilities::System::job_supports_mpi()) + if (Utilities::MPI::job_supports_mpi()) { // make sure that there are not outstanding requests from updating // ghost values or compress diff --git a/include/deal.II/matrix_free/matrix_free.templates.h b/include/deal.II/matrix_free/matrix_free.templates.h index eafe233f0f..67f8264a09 100644 --- a/include/deal.II/matrix_free/matrix_free.templates.h +++ b/include/deal.II/matrix_free/matrix_free.templates.h @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2011 - 2014 by the deal.II authors +// Copyright (C) 2011 - 2015 by the deal.II authors // // This file is part of the deal.II library. // @@ -81,7 +81,7 @@ namespace internal dynamic_cast*>(&tria); if (dist_tria != 0) { - if (Utilities::System::job_supports_mpi()) + if (Utilities::MPI::job_supports_mpi()) { int communicators_same = 0; MPI_Comm_compare (dist_tria->get_communicator(), comm_mf, @@ -136,7 +136,7 @@ internal_reinit(const Mapping &mapping, internal::assert_communicator_equality (dof_handler[0]->get_tria(), additional_data.mpi_communicator); size_info.communicator = additional_data.mpi_communicator; - if (Utilities::System::job_supports_mpi() == true) + if (Utilities::MPI::job_supports_mpi() == true) { size_info.my_pid = Utilities::MPI::this_mpi_process(size_info.communicator); @@ -266,7 +266,7 @@ internal_reinit(const Mapping &mapping, internal::assert_communicator_equality (dof_handler[0]->get_tria(), additional_data.mpi_communicator); size_info.communicator = additional_data.mpi_communicator; - if (Utilities::System::job_supports_mpi() == true) + if (Utilities::MPI::job_supports_mpi() == true) { size_info.my_pid = Utilities::MPI::this_mpi_process(size_info.communicator); @@ -890,21 +890,8 @@ namespace internal void SizeInfo::print_memory_statistics (STREAM &out, std::size_t data_length) const { - Utilities::MPI::MinMaxAvg memory_c; - if (Utilities::System::job_supports_mpi() == true) - { - memory_c = Utilities::MPI::min_max_avg (1e-6*data_length, - communicator); - } - else - { - memory_c.sum = 1e-6*data_length; - memory_c.min = memory_c.sum; - memory_c.max = memory_c.sum; - memory_c.avg = memory_c.sum; - memory_c.min_index = 0; - memory_c.max_index = 0; - } + Utilities::MPI::MinMaxAvg memory_c + = Utilities::MPI::min_max_avg (1e-6*data_length, communicator); if (n_procs < 2) out << memory_c.min; else diff --git a/source/base/mpi.cc b/source/base/mpi.cc index e0d65cd1a1..b9a8a0ce76 100644 --- a/source/base/mpi.cc +++ b/source/base/mpi.cc @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2005 - 2014 by the deal.II authors +// Copyright (C) 2005 - 2015 by the deal.II authors // // This file is part of the deal.II library. // @@ -230,6 +230,21 @@ namespace Utilities min_max_avg(const double my_value, const MPI_Comm &mpi_communicator) { + // If MPI was not started, we have a serial computation and cannot run + // the other MPI commands + if (job_supports_mpi() == false) + { + MinMaxAvg result; + result.sum = my_value; + result.avg = my_value; + result.min = my_value; + result.max = my_value; + result.min_index = 0; + result.max_index = 0; + + return result; + } + // To avoid uninitialized values on some MPI implementations, provide // result with a default value already... MinMaxAvg result = { 0., std::numeric_limits::max(), @@ -501,10 +516,7 @@ namespace Utilities // when running PETSc, because we initialize MPI ourselves before calling // PetscInitialize #ifdef DEAL_II_WITH_MPI - int MPI_has_been_started = 0; - MPI_Initialized(&MPI_has_been_started); - if (Utilities::System::job_supports_mpi() == true && - MPI_has_been_started != 0) + if (job_supports_mpi() == true) { if (std::uncaught_exception()) { diff --git a/source/base/partitioner.cc b/source/base/partitioner.cc index a75dad66b5..39240e1e4e 100644 --- a/source/base/partitioner.cc +++ b/source/base/partitioner.cc @@ -141,7 +141,6 @@ namespace Utilities n_ghost_indices_data = ghost_indices_data.n_elements(); have_ghost_indices = - Utilities::System::job_supports_mpi() && Utilities::MPI::sum(n_ghost_indices_data, communicator) > 0; // In the rest of this function, we determine the point-to-point diff --git a/source/base/timer.cc b/source/base/timer.cc index a4ef54a23d..911358bc52 100644 --- a/source/base/timer.cc +++ b/source/base/timer.cc @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 1998 - 2014 by the deal.II authors +// Copyright (C) 1998 - 2015 by the deal.II authors // // This file is part of the deal.II library. // @@ -170,7 +170,7 @@ double Timer::stop () #endif #ifdef DEAL_II_WITH_MPI - if (sync_wall_time && Utilities::System::job_supports_mpi()) + if (sync_wall_time && Utilities::MPI::job_supports_mpi()) { this->mpi_data = Utilities::MPI::min_max_avg (last_lap_time, mpi_communicator); @@ -188,16 +188,8 @@ double Timer::stop () double Timer::get_lap_time() const { - // time already has the difference - // between the last start()/stop() - // cycle. -#ifdef DEAL_II_WITH_MPI - if (Utilities::System::job_supports_mpi()) - return Utilities::MPI::max (last_lap_time, mpi_communicator); - else -#endif - return last_lap_time; - + // time already has the difference between the last start()/stop() cycle. + return Utilities::MPI::max (last_lap_time, mpi_communicator); } @@ -219,16 +211,12 @@ double Timer::operator() () const const double running_time = dtime - start_time + dtime_children - start_time_children + cumulative_time; - if (Utilities::System::job_supports_mpi()) - // in case of MPI, need to get the time - // passed by summing the time over all - // processes in the network. works also - // in case we just want to have the time - // of a single thread, since then the - // communicator is MPI_COMM_SELF - return Utilities::MPI::sum (running_time, mpi_communicator); - else - return running_time; + // in case of MPI, need to get the time passed by summing the time over + // all processes in the network. works also in case we just want to have + // the time of a single thread, since then the communicator is + // MPI_COMM_SELF + return Utilities::MPI::sum (running_time, mpi_communicator); + #elif defined(DEAL_II_MSVC) const double running_time = windows::cpu_clock() - start_time + cumulative_time; return running_time; @@ -238,10 +226,7 @@ double Timer::operator() () const } else { - if (Utilities::System::job_supports_mpi()) - return Utilities::MPI::sum (cumulative_time, mpi_communicator); - else - return cumulative_time; + return Utilities::MPI::sum (cumulative_time, mpi_communicator); } } @@ -426,11 +411,7 @@ TimerOutput::leave_subsection (const std::string §ion_name) // initialize the Timers here according to that double cpu_time = sections[actual_section_name].timer(); sections[actual_section_name].total_cpu_time - += (Utilities::System::job_supports_mpi() - ? - Utilities::MPI::sum (cpu_time, mpi_communicator) - : - cpu_time); + += Utilities::MPI::sum (cpu_time, mpi_communicator); // in case we have to print out something, do that here... if ((output_frequency == every_call || output_frequency == every_call_and_summary) @@ -474,12 +455,7 @@ TimerOutput::print_summary () const // in case we want to write CPU times if (output_type != wall_times) { - double total_cpu_time = (Utilities::System::job_supports_mpi() - ? - Utilities::MPI::sum (timer_all(), - mpi_communicator) - : - timer_all()); + double total_cpu_time = Utilities::MPI::sum(timer_all(), mpi_communicator); // check that the sum of all times is // less or equal than the total diff --git a/source/base/utilities.cc b/source/base/utilities.cc index 2de85686a2..8af548922d 100644 --- a/source/base/utilities.cc +++ b/source/base/utilities.cc @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2005 - 2013 by the deal.II authors +// Copyright (C) 2005 - 2015 by the deal.II authors // // This file is part of the deal.II library. // @@ -619,14 +619,7 @@ namespace Utilities bool job_supports_mpi () { -#ifdef DEAL_II_WITH_MPI - int MPI_has_been_started = 0; - MPI_Initialized(&MPI_has_been_started); - - return true && (MPI_has_been_started > 0); -#else - return false; -#endif + return Utilities::MPI::job_supports_mpi(); } } diff --git a/tests/mpi/collective_01.cc b/tests/mpi/collective_01.cc index b6c06f8058..08fb965042 100644 --- a/tests/mpi/collective_01.cc +++ b/tests/mpi/collective_01.cc @@ -33,7 +33,7 @@ void print_it(Utilities::MPI::MinMaxAvg &result) void test() { - Assert( Utilities::System::job_supports_mpi(), ExcInternalError()); + Assert( Utilities::MPI::job_supports_mpi(), ExcInternalError()); unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); const unsigned int numprocs = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); diff --git a/tests/mpi/distribute_sp_01.cc b/tests/mpi/distribute_sp_01.cc index b2be59d991..216ea26b8f 100644 --- a/tests/mpi/distribute_sp_01.cc +++ b/tests/mpi/distribute_sp_01.cc @@ -30,7 +30,7 @@ void test_mpi() { - Assert( Utilities::System::job_supports_mpi(), ExcInternalError()); + Assert( Utilities::MPI::job_supports_mpi(), ExcInternalError()); unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); diff --git a/tests/mpi/distribute_sp_02.cc b/tests/mpi/distribute_sp_02.cc index db80ed35dd..931ea40064 100644 --- a/tests/mpi/distribute_sp_02.cc +++ b/tests/mpi/distribute_sp_02.cc @@ -31,7 +31,7 @@ void test_mpi() { - Assert( Utilities::System::job_supports_mpi(), ExcInternalError()); + Assert( Utilities::MPI::job_supports_mpi(), ExcInternalError()); unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); const unsigned int numprocs = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); diff --git a/tests/mpi/point_to_point_pattern_01.cc b/tests/mpi/point_to_point_pattern_01.cc index 19fbbe6cfc..aa94ae6763 100644 --- a/tests/mpi/point_to_point_pattern_01.cc +++ b/tests/mpi/point_to_point_pattern_01.cc @@ -27,7 +27,7 @@ void test_mpi() { - Assert( Utilities::System::job_supports_mpi(), ExcInternalError()); + Assert( Utilities::MPI::job_supports_mpi(), ExcInternalError()); unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); const unsigned int numprocs = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); diff --git a/tests/mpi/simple_mpi_01.cc b/tests/mpi/simple_mpi_01.cc index 69a9c35a5d..0ed26678cd 100644 --- a/tests/mpi/simple_mpi_01.cc +++ b/tests/mpi/simple_mpi_01.cc @@ -27,7 +27,7 @@ void test_mpi() { - Assert( Utilities::System::job_supports_mpi(), ExcInternalError()); + Assert( Utilities::MPI::job_supports_mpi(), ExcInternalError()); unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); -- 2.39.5