From d915ad1a5f209a38d86a702af4ea3fc59d4131fa Mon Sep 17 00:00:00 2001 From: Denis Davydov Date: Mon, 26 Mar 2018 15:32:20 +0200 Subject: [PATCH] bugfix: la::p::Vector::compress(insert) now correctly zeroes ghosts --- doc/news/changes/minor/20180326DenisDavydov | 4 + .../lac/la_parallel_vector.templates.h | 5 +- tests/mpi/parallel_vector_03a.cc | 123 ++++++++++++++++++ tests/mpi/parallel_vector_03a.mpirun=2.output | 2 + 4 files changed, 132 insertions(+), 2 deletions(-) create mode 100644 doc/news/changes/minor/20180326DenisDavydov create mode 100644 tests/mpi/parallel_vector_03a.cc create mode 100644 tests/mpi/parallel_vector_03a.mpirun=2.output diff --git a/doc/news/changes/minor/20180326DenisDavydov b/doc/news/changes/minor/20180326DenisDavydov new file mode 100644 index 0000000000..1846744fb5 --- /dev/null +++ b/doc/news/changes/minor/20180326DenisDavydov @@ -0,0 +1,4 @@ +Fixed: LinearAlgebra::distributed::Vector::compress(VectorOperation::insert) now flashes +ghost part of the vector in Release mode. +
+(Denis Davydov, 2018/03/17) diff --git a/include/deal.II/lac/la_parallel_vector.templates.h b/include/deal.II/lac/la_parallel_vector.templates.h index 04e95f50f6..6a1a0c53f9 100644 --- a/include/deal.II/lac/la_parallel_vector.templates.h +++ b/include/deal.II/lac/la_parallel_vector.templates.h @@ -560,8 +560,9 @@ namespace LinearAlgebra #ifdef DEAL_II_WITH_MPI vector_is_ghosted = false; - if (compress_requests.size() == 0) - return; + // in order to zero ghost part of the vector, we need to call + // import_from_ghosted_array_finish() regardless of + // compress_requests.size() == 0 // make this function thread safe Threads::Mutex::ScopedLock lock (mutex); diff --git a/tests/mpi/parallel_vector_03a.cc b/tests/mpi/parallel_vector_03a.cc new file mode 100644 index 0000000000..6d913b82cd --- /dev/null +++ b/tests/mpi/parallel_vector_03a.cc @@ -0,0 +1,123 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2011 - 2015 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. +// +// --------------------------------------------------------------------- + + +// similar to parallel_sparse_vector_03.cc, but make sure +// compress(insert) zeroes out ghosts in Release mode + +#include "../tests.h" +#include +#include +#include +#include +#include +#include + +void check(const unsigned int myid, + const LinearAlgebra::distributed::Vector &v) +{ + if (myid==0) + { + AssertThrow(v(10) == 10.0, ExcInternalError()); + AssertThrow(v(11) == 0., ExcInternalError()); + AssertThrow(v(12) == 0., ExcInternalError()); + AssertThrow(v(14) == 14., ExcInternalError()); + + AssertThrow(v(5) == 55., ExcInternalError()); + } + else + { + AssertThrow(v(4) == 0., ExcInternalError()); + AssertThrow(v(5) == 55., ExcInternalError()); + AssertThrow(v(6) == 66., ExcInternalError()); + } +} + + +void test () +{ + unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + unsigned int numproc = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); + + Assert (numproc==2, ExcNotImplemented()); + + const unsigned int size = 20; + IndexSet local_owned(size); + IndexSet local_nonzero(size); + IndexSet local_relevant(size); + if (myid==0) + { + local_owned.add_range(0,10); + local_nonzero.add_range(5,10); + local_relevant = local_owned; + local_relevant.add_range(10,13); + local_relevant.add_range(14,15); + } + else + { + local_owned.add_range(10,size); + local_nonzero.add_range(10,11); + local_nonzero.add_range(13,15); + local_relevant = local_owned; + local_relevant.add_range(4,7); + } + + LinearAlgebra::distributed::Vector v(local_owned, local_relevant, MPI_COMM_WORLD); + v = 0.; + + // set local values + for (unsigned int i = 0; i < local_nonzero.n_elements(); i++) + v(local_nonzero.nth_index_in_set(i)) = local_nonzero.nth_index_in_set(i); + + // set value from processor which does not own it: + v(5) = 55.; + v.compress(VectorOperation::insert); + + // add to value from processor which has it as a ghost + if (myid == 1) + v(6) = 60; + v.compress(VectorOperation::add); // 60 + 6 + // compress(insert) used to leave ghosts un-touched which resulted in + // the wrong 55+55 for this compress(add) operation. + + v.update_ghost_values(); + + check(myid,v); + + if (myid == 0) + deallog << "OK" << std::endl; +} + + + +int main (int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads()); + + unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + deallog.push(Utilities::int_to_string(myid)); + + if (myid == 0) + { + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog << std::setprecision(4); + + test(); + } + else + test(); + +} diff --git a/tests/mpi/parallel_vector_03a.mpirun=2.output b/tests/mpi/parallel_vector_03a.mpirun=2.output new file mode 100644 index 0000000000..be8d055f86 --- /dev/null +++ b/tests/mpi/parallel_vector_03a.mpirun=2.output @@ -0,0 +1,2 @@ + +DEAL:0::OK -- 2.39.5