From: Wolfgang Bangerth Date: Mon, 29 Sep 2014 15:12:57 +0000 (-0500) Subject: MPI_Allreduce is apparently unable to accept input and output arguments at the same... X-Git-Tag: v8.2.0-rc1~119^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=9739e56818c47e01e6de32dcfb937cfc5f7c1167;p=dealii.git MPI_Allreduce is apparently unable to accept input and output arguments at the same location. Fix this by determining when &input=&output and in that case pass MPI_IN_PLACE. --- diff --git a/doc/news/changes.h b/doc/news/changes.h index 4bb5c9fc1b..62fbcd1efd 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -301,6 +301,13 @@ inconvenience this causes.

Specific improvements

    +
  1. Fixed: The vector and array versions of Utilities::MPI::sum() and + Utilities::MPI::max() produced segmentation faults with some MPI implementations + if the input and output arguments were the same. This is now fixed. +
    + (Wolfgang Bangerth, 2014/09/29) +
  2. +
  3. Fixed: Trying to have FE_Q(p) and FE_DGQ(r) elements next to each other in an hp::DoFHandler object led to assertions saying that these two elements don't know how to compute interface constraints where such diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index 5bbac50539..a6943e3691 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -138,6 +138,8 @@ namespace Utilities * 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. + * + * Input and output arrays may be the same. */ template inline @@ -149,6 +151,8 @@ namespace Utilities * Like the previous function, but take the sums over the elements of a * std::vector. 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. + * + * Input and output vectors may be the same. */ template inline @@ -184,6 +188,8 @@ namespace Utilities * 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. + * + * Input and output arrays may be the same. */ template inline @@ -196,6 +202,8 @@ namespace Utilities * std::vector. In other words, the i-th element of the results array is * the maximum over the i-th entries of the input arrays from each * processor. + * + * Input and output vectors may be the same. */ template inline @@ -415,7 +423,11 @@ namespace Utilities T (&sums)[N]) { #ifdef DEAL_II_WITH_MPI - MPI_Allreduce (const_cast(static_cast(&values[0])), + 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 @@ -434,7 +446,11 @@ namespace Utilities { #ifdef DEAL_II_WITH_MPI sums.resize (values.size()); - MPI_Allreduce (const_cast(static_cast(&values[0])), + 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 @@ -469,7 +485,11 @@ namespace Utilities T (&maxima)[N]) { #ifdef DEAL_II_WITH_MPI - MPI_Allreduce (const_cast(static_cast(&values[0])), + 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 @@ -488,7 +508,11 @@ namespace Utilities { #ifdef DEAL_II_WITH_MPI maxima.resize (values.size()); - MPI_Allreduce (const_cast(static_cast(&values[0])), + 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 diff --git a/tests/mpi/collective_02_array_in_place.cc b/tests/mpi/collective_02_array_in_place.cc new file mode 100644 index 0000000000..49f30056e7 --- /dev/null +++ b/tests/mpi/collective_02_array_in_place.cc @@ -0,0 +1,66 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2009 - 2014 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 arrays, but with input=output + +#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 sums[2] = { 1, 2 }; + Utilities::MPI::sum (sums, + 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_WITH_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"); + 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_in_place.mpirun=1.output b/tests/mpi/collective_02_array_in_place.mpirun=1.output new file mode 100644 index 0000000000..6e10373d77 --- /dev/null +++ b/tests/mpi/collective_02_array_in_place.mpirun=1.output @@ -0,0 +1,2 @@ + +DEAL:mpi::1 2 diff --git a/tests/mpi/collective_02_array_in_place.mpirun=10.output b/tests/mpi/collective_02_array_in_place.mpirun=10.output new file mode 100644 index 0000000000..97f11ae084 --- /dev/null +++ b/tests/mpi/collective_02_array_in_place.mpirun=10.output @@ -0,0 +1,2 @@ + +DEAL:mpi::10 20 diff --git a/tests/mpi/collective_02_array_in_place.mpirun=4.output b/tests/mpi/collective_02_array_in_place.mpirun=4.output new file mode 100644 index 0000000000..ad500b96a3 --- /dev/null +++ b/tests/mpi/collective_02_array_in_place.mpirun=4.output @@ -0,0 +1,2 @@ + +DEAL:mpi::4 8 diff --git a/tests/mpi/collective_02_vector.cc b/tests/mpi/collective_02_vector.cc new file mode 100644 index 0000000000..f9125af7ad --- /dev/null +++ b/tests/mpi/collective_02_vector.cc @@ -0,0 +1,68 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2009 - 2014 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 vectors + +#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 }; + std::vector values(&values_[0], &values_[2]); + std::vector 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_WITH_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"); + 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_vector.mpirun=1.output b/tests/mpi/collective_02_vector.mpirun=1.output new file mode 100644 index 0000000000..6e10373d77 --- /dev/null +++ b/tests/mpi/collective_02_vector.mpirun=1.output @@ -0,0 +1,2 @@ + +DEAL:mpi::1 2 diff --git a/tests/mpi/collective_02_vector.mpirun=10.output b/tests/mpi/collective_02_vector.mpirun=10.output new file mode 100644 index 0000000000..97f11ae084 --- /dev/null +++ b/tests/mpi/collective_02_vector.mpirun=10.output @@ -0,0 +1,2 @@ + +DEAL:mpi::10 20 diff --git a/tests/mpi/collective_02_vector.mpirun=4.output b/tests/mpi/collective_02_vector.mpirun=4.output new file mode 100644 index 0000000000..ad500b96a3 --- /dev/null +++ b/tests/mpi/collective_02_vector.mpirun=4.output @@ -0,0 +1,2 @@ + +DEAL:mpi::4 8 diff --git a/tests/mpi/collective_02_vector_in_place.cc b/tests/mpi/collective_02_vector_in_place.cc new file mode 100644 index 0000000000..d56c02978d --- /dev/null +++ b/tests/mpi/collective_02_vector_in_place.cc @@ -0,0 +1,67 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2009 - 2014 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 vectors, but with input=output + +#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 }; + std::vector sums(&values_[0], &values_[2]); + Utilities::MPI::sum (sums, + 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_WITH_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"); + 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_vector_in_place.mpirun=1.output b/tests/mpi/collective_02_vector_in_place.mpirun=1.output new file mode 100644 index 0000000000..6e10373d77 --- /dev/null +++ b/tests/mpi/collective_02_vector_in_place.mpirun=1.output @@ -0,0 +1,2 @@ + +DEAL:mpi::1 2 diff --git a/tests/mpi/collective_02_vector_in_place.mpirun=10.output b/tests/mpi/collective_02_vector_in_place.mpirun=10.output new file mode 100644 index 0000000000..97f11ae084 --- /dev/null +++ b/tests/mpi/collective_02_vector_in_place.mpirun=10.output @@ -0,0 +1,2 @@ + +DEAL:mpi::10 20 diff --git a/tests/mpi/collective_02_vector_in_place.mpirun=4.output b/tests/mpi/collective_02_vector_in_place.mpirun=4.output new file mode 100644 index 0000000000..ad500b96a3 --- /dev/null +++ b/tests/mpi/collective_02_vector_in_place.mpirun=4.output @@ -0,0 +1,2 @@ + +DEAL:mpi::4 8