From 5eb5a70191ed346d34c04752851c4efa24f3229b Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Fri, 11 Feb 2011 02:49:29 +0000 Subject: [PATCH] Fix a bug and add a testcase for VectorTools::compute_mean_value. git-svn-id: https://svn.dealii.org/trunk@23325 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/numerics/vectors.templates.h | 6 +- tests/mpi/compute_mean_value.cc | 132 ++++++++++++++++++ .../compute_mean_value/ncpu_10/cmp/generic | 3 + .../mpi/compute_mean_value/ncpu_4/cmp/generic | 3 + 4 files changed, 141 insertions(+), 3 deletions(-) create mode 100644 tests/mpi/compute_mean_value.cc create mode 100644 tests/mpi/compute_mean_value/ncpu_10/cmp/generic create mode 100644 tests/mpi/compute_mean_value/ncpu_4/cmp/generic diff --git a/deal.II/include/deal.II/numerics/vectors.templates.h b/deal.II/include/deal.II/numerics/vectors.templates.h index cc2abca4bd..6074cc5147 100644 --- a/deal.II/include/deal.II/numerics/vectors.templates.h +++ b/deal.II/include/deal.II/numerics/vectors.templates.h @@ -4768,9 +4768,9 @@ VectorTools::compute_mean_value (const Mapping &mapping, double my_values[2] = { mean, area }; double global_values[2]; - MPI_Reduce (&my_values, &global_values, 2, MPI_DOUBLE, - MPI_SUM, 0, - p_d_triangulation->get_communicator()); + MPI_Allreduce (&my_values, &global_values, 2, MPI_DOUBLE, + MPI_SUM, + p_d_triangulation->get_communicator()); mean = global_values[0]; area = global_values[1]; diff --git a/tests/mpi/compute_mean_value.cc b/tests/mpi/compute_mean_value.cc new file mode 100644 index 0000000000..396a4de79c --- /dev/null +++ b/tests/mpi/compute_mean_value.cc @@ -0,0 +1,132 @@ +//--------------------------------------------------------------------------- +// $Id: 2d_coarse_grid_01.cc 17444 2008-10-31 19:35:14Z bangerth $ +// 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. +// +//--------------------------------------------------------------------------- + + +// Test VectorTools::compute_mean_value for parallel computations. we +// interpolate a linear function onto the grid with a symmetric mesh. the mean +// value of the interpolation must be the mean of the linear function + +#include "../tests.h" +#include "coarse_grid_common.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +#include + + +template +class LinearFunction : public Function +{ + public: + double value (const Point &p, + const unsigned int) const + { + return p[0] + 2; + } +}; + + +template +void test() +{ + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + GridGenerator::hyper_ball(tr); + tr.refine_global (2); + + const FE_Q fe(2); + DoFHandler dofh(tr); + dofh.distribute_dofs (fe); + + TrilinosWrappers::MPI::Vector interpolated(dofh.locally_owned_dofs(), + MPI_COMM_WORLD); + VectorTools::interpolate (dofh, + LinearFunction(), + interpolated); + + IndexSet relevant_set; + DoFTools::extract_locally_relevant_dofs (dofh, relevant_set); + TrilinosWrappers::MPI::Vector x_rel(relevant_set, MPI_COMM_WORLD); + x_rel = interpolated; + + const double mean + = VectorTools::compute_mean_value (dofh, QGauss(2), x_rel, 0); + + if (Utilities::System::get_this_mpi_process (MPI_COMM_WORLD) == 0) + deallog << "mean=" << mean + << std::endl; + + Assert (std::fabs(mean - 2) < 1e-3, ExcInternalError()); +} + + +int main(int argc, char *argv[]) +{ +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + MPI_Init (&argc,&argv); +#else + (void)argc; + (void)argv; +#endif + + unsigned int myid = Utilities::System::get_this_mpi_process (MPI_COMM_WORLD); + + + deallog.push(Utilities::int_to_string(myid)); + + if (myid == 0) + { + std::ofstream logfile(output_file_for_mpi("compute_mean_value").c_str()); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + deallog.push("2d"); + test<2>(); + deallog.pop(); + + deallog.push("3d"); + test<3>(); + deallog.pop(); + } + else + { + deallog.push("2d"); + test<2>(); + deallog.pop(); + + deallog.push("3d"); + test<3>(); + deallog.pop(); + } + +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + MPI_Finalize(); +#endif +} diff --git a/tests/mpi/compute_mean_value/ncpu_10/cmp/generic b/tests/mpi/compute_mean_value/ncpu_10/cmp/generic new file mode 100644 index 0000000000..75c95df4c9 --- /dev/null +++ b/tests/mpi/compute_mean_value/ncpu_10/cmp/generic @@ -0,0 +1,3 @@ + +DEAL:0:2d::mean=2.00000 +DEAL:0:3d::mean=2.00000 diff --git a/tests/mpi/compute_mean_value/ncpu_4/cmp/generic b/tests/mpi/compute_mean_value/ncpu_4/cmp/generic new file mode 100644 index 0000000000..75c95df4c9 --- /dev/null +++ b/tests/mpi/compute_mean_value/ncpu_4/cmp/generic @@ -0,0 +1,3 @@ + +DEAL:0:2d::mean=2.00000 +DEAL:0:3d::mean=2.00000 -- 2.39.5