From: daniel.arndt Date: Fri, 13 Sep 2013 10:53:37 +0000 (+0000) Subject: Splitted test for periodicity into two tests for periodicity and point_value X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=9b9ad17b22a0166395d205dcb790d6813ab783d6;p=dealii-svn.git Splitted test for periodicity into two tests for periodicity and point_value git-svn-id: https://svn.dealii.org/trunk@30684 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/mpi/periodicity_01.cc b/tests/mpi/periodicity_01.cc index f243ef3dab..413029deeb 100644 --- a/tests/mpi/periodicity_01.cc +++ b/tests/mpi/periodicity_01.cc @@ -298,16 +298,16 @@ namespace Step40 void LaplaceProblem::get_point_value (const Point point, const int proc, Vector &value) const { - std::vector tmp (value.size()); - try - { + typename DoFHandler::active_cell_iterator cell + = GridTools::find_active_cell_around_point (dof_handler, point); + + if (cell->is_locally_owned()) VectorTools::point_value (dof_handler, locally_relevant_solution, point, value); - for (unsigned int i=0; i tmp (value.size()); + for (unsigned int i=0; i tmp2 (value.size()); MPI_Reduce(&(tmp[0]), &(tmp2[0]), value.size(), MPI_DOUBLE, @@ -491,9 +491,9 @@ namespace Step40 (triangulation, /*b_id1*/ 2*i, /*b_id2*/2*i+1, /*direction*/ i, periodicity_vector); - + triangulation.add_periodicity(periodicity_vector); - + triangulation.refine_global (1); } else diff --git a/tests/mpi/point_value_01.cc b/tests/mpi/point_value_01.cc new file mode 100644 index 0000000000..2e1de3432d --- /dev/null +++ b/tests/mpi/point_value_01.cc @@ -0,0 +1,130 @@ +// --------------------------------------------------------------------- +// $Id: point_value_01.cc $ +// +// Copyright (C) 2001 - 2013 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. +// +// --------------------------------------------------------------------- + +// +// Checks whether VectorTools::point_value shows the expected behavior on +// distributed meshes. If the point is locally available everything should be +// fine, else the exception ExcPointNotAvailableHere should be thrown. +// + +#include "../tests.h" +#include +#include + +#include + +#include +#include +#include + +#include +#include +#include +#include +#include + +#include + +#include +#include + +#include +#include + +template +void test() +{ + parallel::distributed::Triangulation triangulation(MPI_COMM_WORLD); + + FE_Q fe(1); + + DoFHandler dof_handler (triangulation); + + GridGenerator::hyper_cube(triangulation); + triangulation.refine_global (1); + + dof_handler.distribute_dofs (fe); + + IndexSet locally_owned_dofs = dof_handler.locally_owned_dofs (); + parallel::distributed::Vector locally_owned_solution + (locally_owned_dofs, MPI_COMM_WORLD); + + locally_owned_solution=1; + + std::vector > points; + + for (int i=0; i<2; i++) + for (int j=0; j<2; j++) + if (dim==3) + for (int k=0; k<2; k++) + points.push_back(Point(.25+.5*i,.25+.5*j,.25+.5*k)); + else + points.push_back(Point(.25+.5*i,.25+.5*j)); + + typename std::vector >::iterator point_iterator, points_end; + point_iterator=points.begin(); + points_end=points.end(); + + Vector value(1); + + for(; point_iterator!=points_end; point_iterator++) + { + try + { + VectorTools::point_value (dof_handler, locally_owned_solution, + *point_iterator, value); + if (std::abs(value[0]-1.)>1e-8) + ExcInternalError(); + } + catch (const VectorTools::ExcPointNotAvailableHere &) + {} + + MPI_Barrier(MPI_COMM_WORLD); + + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + deallog << *point_iterator << " OK" << std::endl; + } +} + + +int main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + if (myid == 0) + { + std::ofstream logfile(output_file_for_mpi("point_value_01").c_str()); + deallog.attach(logfile, false); + 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 + { + test<2>(); + test<3>(); + } + + return 0; +} + diff --git a/tests/mpi/point_value_01/ncpu_2/cmp/generic b/tests/mpi/point_value_01/ncpu_2/cmp/generic new file mode 100644 index 0000000000..fbb1c377f2 --- /dev/null +++ b/tests/mpi/point_value_01/ncpu_2/cmp/generic @@ -0,0 +1,12 @@ +DEAL:2d::0.250000 0.250000 OK +DEAL:2d::0.250000 0.750000 OK +DEAL:2d::0.750000 0.250000 OK +DEAL:2d::0.750000 0.750000 OK +DEAL:3d::0.250000 0.250000 0.250000 OK +DEAL:3d::0.250000 0.250000 0.750000 OK +DEAL:3d::0.250000 0.750000 0.250000 OK +DEAL:3d::0.250000 0.750000 0.750000 OK +DEAL:3d::0.750000 0.250000 0.250000 OK +DEAL:3d::0.750000 0.250000 0.750000 OK +DEAL:3d::0.750000 0.750000 0.250000 OK +DEAL:3d::0.750000 0.750000 0.750000 OK diff --git a/tests/mpi/point_value_01/ncpu_4/cmp/generic b/tests/mpi/point_value_01/ncpu_4/cmp/generic new file mode 100644 index 0000000000..fbb1c377f2 --- /dev/null +++ b/tests/mpi/point_value_01/ncpu_4/cmp/generic @@ -0,0 +1,12 @@ +DEAL:2d::0.250000 0.250000 OK +DEAL:2d::0.250000 0.750000 OK +DEAL:2d::0.750000 0.250000 OK +DEAL:2d::0.750000 0.750000 OK +DEAL:3d::0.250000 0.250000 0.250000 OK +DEAL:3d::0.250000 0.250000 0.750000 OK +DEAL:3d::0.250000 0.750000 0.250000 OK +DEAL:3d::0.250000 0.750000 0.750000 OK +DEAL:3d::0.750000 0.250000 0.250000 OK +DEAL:3d::0.750000 0.250000 0.750000 OK +DEAL:3d::0.750000 0.750000 0.250000 OK +DEAL:3d::0.750000 0.750000 0.750000 OK diff --git a/tests/mpi/point_value_01/ncpu_6/cmp/generic b/tests/mpi/point_value_01/ncpu_6/cmp/generic new file mode 100644 index 0000000000..fbb1c377f2 --- /dev/null +++ b/tests/mpi/point_value_01/ncpu_6/cmp/generic @@ -0,0 +1,12 @@ +DEAL:2d::0.250000 0.250000 OK +DEAL:2d::0.250000 0.750000 OK +DEAL:2d::0.750000 0.250000 OK +DEAL:2d::0.750000 0.750000 OK +DEAL:3d::0.250000 0.250000 0.250000 OK +DEAL:3d::0.250000 0.250000 0.750000 OK +DEAL:3d::0.250000 0.750000 0.250000 OK +DEAL:3d::0.250000 0.750000 0.750000 OK +DEAL:3d::0.750000 0.250000 0.250000 OK +DEAL:3d::0.750000 0.250000 0.750000 OK +DEAL:3d::0.750000 0.750000 0.250000 OK +DEAL:3d::0.750000 0.750000 0.750000 OK