From ae9b1ba84fa5a86066a23aaf135fccbc5360e09d Mon Sep 17 00:00:00 2001 From: kronbichler Date: Thu, 27 Dec 2012 14:15:04 +0000 Subject: [PATCH] Avoid calls to MPI_Allreduce for vectors where only one processor participates. git-svn-id: https://svn.dealii.org/trunk@27859 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/lac/parallel_vector.h | 93 +++++++++++++------ tests/mpi/parallel_vector_06.cc | 6 +- .../mpi/parallel_vector_06/ncpu_1/cmp/generic | 24 +++++ 3 files changed, 92 insertions(+), 31 deletions(-) create mode 100644 tests/mpi/parallel_vector_06/ncpu_1/cmp/generic diff --git a/deal.II/include/deal.II/lac/parallel_vector.h b/deal.II/include/deal.II/lac/parallel_vector.h index d16f47870d..a5cc144f9c 100644 --- a/deal.II/include/deal.II/lac/parallel_vector.h +++ b/deal.II/include/deal.II/lac/parallel_vector.h @@ -1171,7 +1171,11 @@ namespace parallel // use int instead of bool int local_result = (partitioner->local_size()>0 ? -vector_view.all_zero () : -1); - return - Utilities::MPI::max(local_result,partitioner->get_communicator()); + if (partitioner->n_mpi_processes() > 1) + return -Utilities::MPI::max(local_result, + partitioner->get_communicator()); + else + return local_result; } @@ -1181,10 +1185,17 @@ namespace parallel bool Vector::is_non_negative () const { - // use int instead of bool + // use int instead of bool. in order to make this operation work also + // when MPI_Init was not called, only call MPI_Allreduce commands when + // there is more than one processor (note that reinit() functions handle + // this case correctly through the job_supports_mpi() query) int local_result = (partitioner->local_size()>0 ? -vector_view.is_non_negative () : -1); - return - Utilities::MPI::max(local_result,partitioner->get_communicator()); + if (partitioner->n_mpi_processes() > 1) + return -Utilities::MPI::max(local_result, + partitioner->get_communicator()); + else + return local_result; } @@ -1197,15 +1208,18 @@ namespace parallel { AssertDimension (local_size(), v.local_size()); - // MPI does not support bools, so use unsigned - // int instead. Two vectors are equal if the - // check for non-equal fails on all processors + // MPI does not support bools, so use unsigned int instead. Two vectors + // are equal if the check for non-equal fails on all processors unsigned int local_result = (partitioner->local_size()>0 ? vector_view.template operator != (v.vector_view) : 0 ); - unsigned int result = Utilities::MPI::max(local_result, - partitioner->get_communicator()); + unsigned int result = + partitioner->n_mpi_processes() > 1 + ? + Utilities::MPI::max(local_result, partitioner->get_communicator()) + : + local_result; return result==0; } @@ -1231,7 +1245,11 @@ namespace parallel Number local_result = (partitioner->local_size()>0 ? vector_view.operator* (V.vector_view) : 0); - return Utilities::MPI::sum (local_result, partitioner->get_communicator()); + if (partitioner->n_mpi_processes() > 1) + return Utilities::MPI::sum (local_result, + partitioner->get_communicator()); + else + return local_result; } @@ -1244,10 +1262,14 @@ namespace parallel // on some processors, the size might be zero, // which is not allowed by the base // class. Therefore, insert a check here - Number local_result = (partitioner->local_size()>0 ? - vector_view.norm_sqr() - : 0); - return Utilities::MPI::sum(local_result,partitioner->get_communicator()); + real_type local_result = (partitioner->local_size()>0 ? + vector_view.norm_sqr() + : 0); + if (partitioner->n_mpi_processes() > 1) + return Utilities::MPI::sum(local_result, + partitioner->get_communicator()); + else + return local_result; } @@ -1262,7 +1284,11 @@ namespace parallel vector_view.mean_value() : 0) *((real_type)partitioner->local_size()/(real_type)partitioner->size()); - return Utilities::MPI::sum (local_result, partitioner->get_communicator()); + if (partitioner->n_mpi_processes() > 1) + return Utilities::MPI::sum (local_result, + partitioner->get_communicator()); + else + return local_result; } @@ -1272,10 +1298,14 @@ namespace parallel typename Vector::real_type Vector::l1_norm () const { - Number local_result = (partitioner->local_size()>0 ? - vector_view.l1_norm() - : 0); - return Utilities::MPI::sum(local_result, partitioner->get_communicator()); + real_type local_result = (partitioner->local_size()>0 ? + vector_view.l1_norm() + : 0); + if (partitioner->n_mpi_processes() > 1) + return Utilities::MPI::sum(local_result, + partitioner->get_communicator()); + else + return local_result; } @@ -1295,12 +1325,15 @@ namespace parallel typename Vector::real_type Vector::lp_norm (const real_type p) const { - const Number local_result = (partitioner->local_size()>0 ? - std::pow(vector_view.lp_norm(p),p) - : 0); - return std::pow (Utilities::MPI::sum(local_result, - partitioner->get_communicator()), - static_cast(1.0/p)); + const real_type local_result = (partitioner->local_size()>0 ? + vector_view.lp_norm(p) + : 0); + if (partitioner->n_mpi_processes() > 1) + return std::pow (Utilities::MPI::sum(std::pow(local_result,p), + partitioner->get_communicator()), + static_cast(1.0/p)); + else + return local_result; } @@ -1310,10 +1343,14 @@ namespace parallel typename Vector::real_type Vector::linfty_norm () const { - const Number local_result = (partitioner->local_size()>0 ? - vector_view.linfty_norm() - : 0); - return Utilities::MPI::max (local_result, partitioner->get_communicator()); + const real_type local_result = (partitioner->local_size()>0 ? + vector_view.linfty_norm() + : 0); + if (partitioner->n_mpi_processes() > 1) + return Utilities::MPI::max (local_result, + partitioner->get_communicator()); + else + return local_result; } diff --git a/tests/mpi/parallel_vector_06.cc b/tests/mpi/parallel_vector_06.cc index 366b619c24..2e1380d482 100644 --- a/tests/mpi/parallel_vector_06.cc +++ b/tests/mpi/parallel_vector_06.cc @@ -90,8 +90,8 @@ void test () } // check mean value (should be equal to l1 - // norm here since we have no negative - // entries) + // norm divided by vector size here since we + // have no negative entries) { const double mean = v.mean_value(); if (myid == 0) @@ -203,7 +203,7 @@ void test () // only one processor has non-negative entry v3 = v; - if (myid == 1) + if (myid == 1 || numproc==1) v3.local_element(0) = -1; allnonneg = v3.is_non_negative(); if (myid == 0) diff --git a/tests/mpi/parallel_vector_06/ncpu_1/cmp/generic b/tests/mpi/parallel_vector_06/ncpu_1/cmp/generic new file mode 100644 index 0000000000..d16831617e --- /dev/null +++ b/tests/mpi/parallel_vector_06/ncpu_1/cmp/generic @@ -0,0 +1,24 @@ + +DEAL:0::numproc=4 +DEAL:0::l2 norm: 23.66 +DEAL:0::l1 norm: 56.00 +DEAL:0::linfty norm: 14.00 +DEAL:0::l2.2 norm: 22.04 +DEAL:0::Mean value: 7.000 +DEAL:0::Inner product: 312.0 +DEAL:0:: v==v2 ? 1 +DEAL:0:: v!=v2 ? 0 +DEAL:0:: v==v2 ? 0 +DEAL:0:: v!=v2 ? 1 +DEAL:0:: v==v2 ? 1 +DEAL:0:: v!=v2 ? 0 +DEAL:0:: v==v2 ? 0 +DEAL:0:: v!=v2 ? 1 +DEAL:0:: v==0 ? 0 +DEAL:0:: v2==0 ? 1 +DEAL:0:: v2==0 ? 0 +DEAL:0:: v>=0 ? 1 +DEAL:0:: v2>=0 ? 0 +DEAL:0:: v3>=0 ? 1 +DEAL:0:: v3>=0 ? 0 +DEAL:0::OK -- 2.39.5