From 4caf62bfcd997b2358bfacb9400f63823fac9f03 Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Sun, 3 Sep 2017 14:52:54 -0400 Subject: [PATCH] New ReadWriteVector::reinit(TrilinosWrappers::MPI::Vector) Let a ReadWriteVector be a view for a (potentially ghosted) TrilinosWrappers::MPI::Vector. This is the first step towards using RWV inside the library. --- include/deal.II/lac/read_write_vector.h | 18 ++++ .../deal.II/lac/read_write_vector.templates.h | 29 ++++++ tests/trilinos/readwritevector_03.cc | 92 +++++++++++++++++++ .../readwritevector_03.mpirun=2.output | 73 +++++++++++++++ 4 files changed, 212 insertions(+) create mode 100644 tests/trilinos/readwritevector_03.cc create mode 100644 tests/trilinos/readwritevector_03.mpirun=2.output diff --git a/include/deal.II/lac/read_write_vector.h b/include/deal.II/lac/read_write_vector.h index d8a53ac8e5..8a89f0229a 100644 --- a/include/deal.II/lac/read_write_vector.h +++ b/include/deal.II/lac/read_write_vector.h @@ -206,6 +206,24 @@ namespace LinearAlgebra virtual void reinit (const IndexSet &locally_stored_indices, const bool omit_zeroing_entries = false); + +#ifdef DEAL_II_WITH_TRILINOS +#ifdef DEAL_II_WITH_MPI + /** + * Initialize this ReadWriteVector by supplying access to all locally available + * entries in the given ghosted or non-ghosted vector. + * + * @note This function currently copies the values from the argument into + * the ReadWriteVector, so modifications here will not modify @p trilinos_vec. + * + * This function is mainly written for backwards-compatibility to get + * element access to a ghosted TrilinosWrappers::MPI::Vector inside the + * library. + */ + void reinit(const TrilinosWrappers::MPI::Vector &trilinos_vec); +#endif +#endif + /** * Apply the functor @p func to each element of the vector. The functor * should look like diff --git a/include/deal.II/lac/read_write_vector.templates.h b/include/deal.II/lac/read_write_vector.templates.h index ffabed60e7..e21a8ea3b1 100644 --- a/include/deal.II/lac/read_write_vector.templates.h +++ b/include/deal.II/lac/read_write_vector.templates.h @@ -137,6 +137,35 @@ namespace LinearAlgebra +#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI) + template + void + ReadWriteVector::reinit(const TrilinosWrappers::MPI::Vector &trilinos_vec) + { + // TODO: We could avoid copying the data by just using a view into the + // trilinos data but only if Number=double. Also update documentation that + // the argument's lifetime needs to be longer then. If we do this, we need + // to think about whether the view should be read/write. + + stored_elements = IndexSet(trilinos_vec.vector_partitioner()); + + resize_val(stored_elements.n_elements()); + + TrilinosScalar *start_ptr; + int leading_dimension; + int ierr = trilinos_vec.trilinos_vector().ExtractView (&start_ptr, &leading_dimension); + AssertThrow (ierr == 0, ExcTrilinosError(ierr)); + + std::copy(start_ptr, start_ptr + leading_dimension, val); + + // reset the communication pattern + source_stored_elements.clear(); + comm_pattern.reset(); + } +#endif + + + template template void diff --git a/tests/trilinos/readwritevector_03.cc b/tests/trilinos/readwritevector_03.cc new file mode 100644 index 0000000000..86ec4ac7ba --- /dev/null +++ b/tests/trilinos/readwritevector_03.cc @@ -0,0 +1,92 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 - 2016 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. +// +// --------------------------------------------------------------------- + +// test TrilinosWrappers::MPI::Vector -> RWV + +#include "../tests.h" +#include +#include +#include +#include + +#include + +void test() +{ + IndexSet is(8); + unsigned int rank = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + if (rank==0) + is.add_range(0,4); + if (rank==1) + is.add_range(4,8); + is.compress(); + + deallog << "is: "; + is.print(deallog); + + IndexSet is_ghosted(8); + is_ghosted.add_range(2,6); + deallog << "is_ghosted: "; + is_ghosted.print(deallog); + + TrilinosWrappers::MPI::Vector tril_vector(is); + TrilinosWrappers::MPI::Vector tril_vector_ghosted; + tril_vector_ghosted.reinit(is, is_ghosted, MPI_COMM_WORLD); + for (unsigned int i=0; i<8; ++i) + tril_vector[i] = i; + + tril_vector.compress(VectorOperation::insert); + deallog << "trilinos vec:" << std::endl; + tril_vector.print(deallog.get_file_stream()); + + tril_vector_ghosted = tril_vector; + + deallog << "trilinos vec ghosted:" << std::endl; + tril_vector_ghosted.print(deallog.get_file_stream()); + + + IndexSet readwrite_is (tril_vector_ghosted.vector_partitioner()); + deallog << "ghosted IS: "; + readwrite_is.print(deallog); + + deallog << "tril_vector_ghosted.owned_elements() "; + tril_vector_ghosted.locally_owned_elements().print(deallog); + + { + LinearAlgebra::ReadWriteVector readwrite; + readwrite.reinit(tril_vector); + + deallog << "RWVector contents from tril_vector:" << std::endl; + readwrite.print(deallog.get_file_stream()); + } + { + LinearAlgebra::ReadWriteVector readwrite; + readwrite.reinit(tril_vector_ghosted); + + deallog << "RWVector contents from tril_vector_ghosted:" << std::endl; + readwrite.print(deallog.get_file_stream()); + } + + deallog << "OK" <