From ad36e1bf78e76c0bc4ae9a0b1c10e13be7d324ea Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Wed, 8 Mar 2017 10:36:29 +0100 Subject: [PATCH] Add test for communication patterns supported now. --- tests/mpi/communicator_validity.cc | 76 +++++++++++++++++++ .../mpi/communicator_validity.mpirun=2.output | 5 ++ 2 files changed, 81 insertions(+) create mode 100644 tests/mpi/communicator_validity.cc create mode 100644 tests/mpi/communicator_validity.mpirun=2.output diff --git a/tests/mpi/communicator_validity.cc b/tests/mpi/communicator_validity.cc new file mode 100644 index 0000000000..6a4b5c6160 --- /dev/null +++ b/tests/mpi/communicator_validity.cc @@ -0,0 +1,76 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 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 that LinearAlgebra::distributed::Vector::operator= does not carry +// over any state from a duplicated MPI communicator in conjunction with a +// distributed mesh + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +template +void do_test(MPI_Comm communicator) +{ + const int dim = 2; + parallel::distributed::Triangulation tria(communicator); + GridGenerator::hyper_cube(tria); + tria.refine_global(2); + FE_DGQ fe(0); + DoFHandler dof(tria); + dof.distribute_dofs(fe); + + VectorType v1; + v1.reinit(dof.locally_owned_dofs(), communicator); + if (dof.locally_owned_dofs().n_elements() > 0) + v1.local_element(0) = 1; + + GrowingVectorMemory memory; + typename VectorMemory::Pointer v3(memory); + *v3 = v1; + + deallog << v1.l2_norm() << " " << v3->l2_norm() << std::endl; +} + + + +int main (int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads()); + + mpi_initlog(); + + do_test >(MPI_COMM_WORLD); + + { + MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD); + do_test >(communicator); + MPI_Comm_free (&communicator); + } + { + MPI_Comm communicator = Utilities::MPI::duplicate_communicator(MPI_COMM_WORLD); + do_test >(communicator); + MPI_Comm_free (&communicator); + } + + do_test >(MPI_COMM_WORLD); +} diff --git a/tests/mpi/communicator_validity.mpirun=2.output b/tests/mpi/communicator_validity.mpirun=2.output new file mode 100644 index 0000000000..b5fba2e3b9 --- /dev/null +++ b/tests/mpi/communicator_validity.mpirun=2.output @@ -0,0 +1,5 @@ + +DEAL::1.41421 1.41421 +DEAL::1.41421 1.41421 +DEAL::1.41421 1.41421 +DEAL::1.41421 1.41421 -- 2.39.5