From 54c472d59df48aaccda8c54536efdc7d1170e819 Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Fri, 27 Jan 2017 09:50:32 -0500 Subject: [PATCH] Fix a bug in import function of distributed::Vector. --- .../lac/la_parallel_vector.templates.h | 21 +++-- .../deal.II/lac/read_write_vector.templates.h | 2 +- tests/mpi/parallel_vector_20.cc | 84 +++++++++++++++++++ tests/mpi/parallel_vector_20.mpirun=2.output | 2 + 4 files changed, 103 insertions(+), 6 deletions(-) create mode 100644 tests/mpi/parallel_vector_20.cc create mode 100644 tests/mpi/parallel_vector_20.mpirun=2.output diff --git a/include/deal.II/lac/la_parallel_vector.templates.h b/include/deal.II/lac/la_parallel_vector.templates.h index 4b8fdd1b51..2033688dea 100644 --- a/include/deal.II/lac/la_parallel_vector.templates.h +++ b/include/deal.II/lac/la_parallel_vector.templates.h @@ -827,13 +827,20 @@ namespace LinearAlgebra VectorOperation::values operation, std_cxx11::shared_ptr communication_pattern) { + IndexSet locally_owned_elem = locally_owned_elements(); // If no communication pattern is given, create one. Otherwise, use the // given one. std_cxx11::shared_ptr comm_pattern; if (communication_pattern.get() == NULL) { - comm_pattern.reset(new Utilities::MPI::Partitioner(locally_owned_elements(), - V.get_stored_elements(), + // Split the IndexSet of V in locally owned elements and ghost indices + // then create the communication pattern + IndexSet ghost_indices(V.get_stored_elements()); + ghost_indices.subtract_set(locally_owned_elem); + IndexSet local_indices(V.get_stored_elements()); + local_indices.subtract_set(ghost_indices); + comm_pattern.reset(new Utilities::MPI::Partitioner(local_indices, + ghost_indices, get_mpi_communicator())); } else @@ -856,9 +863,13 @@ namespace LinearAlgebra tmp_vector.compress(operation); - dealii::internal::VectorOperations::Vector_copy copier(tmp_vector.val, val); - internal::VectorOperations::parallel_for(copier, partitioner->local_size(), - thread_loop_partitioner); + // Copy the local elements of tmp_vector to the right place in val + IndexSet tmp_index_set = tmp_vector.locally_owned_elements(); + for (size_type i=0; ioperator= (Number()); - // reset the communication patter + // reset the communication pattern source_stored_elements.clear(); comm_pattern.reset(); } diff --git a/tests/mpi/parallel_vector_20.cc b/tests/mpi/parallel_vector_20.cc new file mode 100644 index 0000000000..fb9ed93be5 --- /dev/null +++ b/tests/mpi/parallel_vector_20.cc @@ -0,0 +1,84 @@ +// --------------------------------------------------------------------- +// +// 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 import function + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +void test() +{ + unsigned int my_id = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + unsigned int n_procs = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); + + IndexSet locally_owned(n_procs*2); + locally_owned.add_range(my_id*2, my_id*2+2); + IndexSet read_write_owned(4); + read_write_owned.add_range(my_id*2, my_id*2+2); + + LinearAlgebra::distributed::Vector v(locally_owned, MPI_COMM_WORLD); + LinearAlgebra::ReadWriteVector read_write_vector(read_write_owned); + read_write_vector.local_element(0) = 1.; + read_write_vector.local_element(1) = 2.; + + v.import(read_write_vector, VectorOperation::values::insert); + + AssertThrow(v.local_element(0) == 1., ExcInternalError()); + AssertThrow(v.local_element(1) == 2., ExcInternalError()); + + read_write_owned.clear(); + read_write_owned.add_index(1); + read_write_owned.add_index(2); + read_write_vector.reinit(read_write_owned); + read_write_vector.local_element(0) = 1.; + read_write_vector.local_element(1) = 2.; + + v.import(read_write_vector, VectorOperation::values::insert); + + AssertThrow(v.local_element(0) == my_id+1, ExcInternalError()); + AssertThrow(v.local_element(1) == my_id+1, ExcInternalError()); + + if (my_id == 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); + deallog.threshold_double(1.e-10); + + test(); + } + else + test(); + +} diff --git a/tests/mpi/parallel_vector_20.mpirun=2.output b/tests/mpi/parallel_vector_20.mpirun=2.output new file mode 100644 index 0000000000..be8d055f86 --- /dev/null +++ b/tests/mpi/parallel_vector_20.mpirun=2.output @@ -0,0 +1,2 @@ + +DEAL:0::OK -- 2.39.5