From 65325c31942cd57815e7d53cd099b3301e359799 Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Sun, 16 Sep 2018 19:50:56 -0400 Subject: [PATCH] Add TrilinosWrappers::MPI::Vector::import(ReadWriteVector) --- doc/news/changes/minor/20180916Heister | 3 + include/deal.II/lac/trilinos_vector.h | 17 ++++ source/lac/trilinos_vector.cc | 31 +++++++ tests/trilinos/readwritevector_04.cc | 81 +++++++++++++++++++ .../readwritevector_04.mpirun=2.output | 57 +++++++++++++ 5 files changed, 189 insertions(+) create mode 100644 doc/news/changes/minor/20180916Heister create mode 100644 tests/trilinos/readwritevector_04.cc create mode 100644 tests/trilinos/readwritevector_04.mpirun=2.output diff --git a/doc/news/changes/minor/20180916Heister b/doc/news/changes/minor/20180916Heister new file mode 100644 index 0000000000..aeb1d3aa5b --- /dev/null +++ b/doc/news/changes/minor/20180916Heister @@ -0,0 +1,3 @@ +New: added the function TrilinosWrappers::MPI::Vector::import() to get data from a ReadWriteVector. +
+(Timo Heister, 2018/09/16) diff --git a/include/deal.II/lac/trilinos_vector.h b/include/deal.II/lac/trilinos_vector.h index 865159f3e0..29a9555829 100644 --- a/include/deal.II/lac/trilinos_vector.h +++ b/include/deal.II/lac/trilinos_vector.h @@ -47,6 +47,12 @@ DEAL_II_NAMESPACE_OPEN +namespace LinearAlgebra +{ + // Forward declaration + template + class ReadWriteVector; +} // namespace LinearAlgebra /** * @addtogroup TrilinosWrappers @@ -685,6 +691,17 @@ namespace TrilinosWrappers const dealii::TrilinosWrappers::SparseMatrix &matrix, const Vector & vector); + /** + * Imports all the elements present in the vector's IndexSet from the + * input vector @p rwv. VectorOperation::values @p operation is used to decide if + * the elements in @p rwv should be added to the current vector or replace the + * current elements. + */ + void + import(const LinearAlgebra::ReadWriteVector &rwv, + const VectorOperation::values operation); + + /** * Test for equality. This function assumes that the present vector and * the one to compare with have the same size already, since comparing diff --git a/source/lac/trilinos_vector.cc b/source/lac/trilinos_vector.cc index d0ac6a9b3d..7d084d7634 100644 --- a/source/lac/trilinos_vector.cc +++ b/source/lac/trilinos_vector.cc @@ -20,6 +20,7 @@ # include # include +# include # include # include # include @@ -544,6 +545,36 @@ namespace TrilinosWrappers } + void + Vector::import(const LinearAlgebra::ReadWriteVector &rwv, + const VectorOperation::values operation) + { + Assert( + this->size() == rwv.size(), + ExcMessage( + "Both vectors need to have the same size for import() to work!")); + // TODO: a generic import() function should handle any kind of data layout + // in ReadWriteVector, but this function is of limited use as this class + // will (hopefully) be retired eventually. + Assert(this->locally_owned_elements() == rwv.get_stored_elements(), + ExcNotImplemented()); + + if (operation == VectorOperation::insert) + { + for (const auto idx : this->locally_owned_elements()) + (*this)[idx] = rwv[idx]; + } + else if (operation == VectorOperation::add) + { + for (const auto idx : this->locally_owned_elements()) + (*this)[idx] += rwv[idx]; + } + else + AssertThrow(false, ExcNotImplemented()); + + this->compress(operation); + } + void Vector::compress(::dealii::VectorOperation::values given_last_action) diff --git a/tests/trilinos/readwritevector_04.cc b/tests/trilinos/readwritevector_04.cc new file mode 100644 index 0000000000..4b17c7b98a --- /dev/null +++ b/tests/trilinos/readwritevector_04.cc @@ -0,0 +1,81 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 - 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +// test: RWV::import() from TrilinosWrappers::MPI::Vector + +#include +#include + +#include +#include + +#include + +#include "../tests.h" + +void +test() +{ + IndexSet is(8); + unsigned int rank = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + if (rank == 0) + is.add_range(0, 4); + else if (rank == 1) + is.add_range(4, 8); + else + AssertThrow(false, ExcNotImplemented()); + + is.compress(); + + deallog << "IS: "; + is.print(deallog); + + TrilinosWrappers::MPI::Vector tril_vector(is); + LinearAlgebra::ReadWriteVector readwrite(is); + + for (auto idx : is) + readwrite[idx] = idx; + + deallog << "RWVector contents:" << std::endl; + readwrite.print(deallog.get_file_stream()); + + // import RWV->Trilinos + tril_vector.import(readwrite, VectorOperation::insert); + deallog << "trilinos vec:" << std::endl; + tril_vector.print(deallog.get_file_stream()); + + // test that ::add also works + tril_vector.import(readwrite, VectorOperation::add); + deallog << "trilinos vec (2x):" << std::endl; + tril_vector.print(deallog.get_file_stream()); + + // import again overwriting the contents + tril_vector.import(readwrite, VectorOperation::insert); + deallog << "trilinos vec (1x):" << std::endl; + tril_vector.print(deallog.get_file_stream()); + + deallog << "OK" << std::endl; +} + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + MPILogInitAll log; + + test(); +} diff --git a/tests/trilinos/readwritevector_04.mpirun=2.output b/tests/trilinos/readwritevector_04.mpirun=2.output new file mode 100644 index 0000000000..821a244c71 --- /dev/null +++ b/tests/trilinos/readwritevector_04.mpirun=2.output @@ -0,0 +1,57 @@ + +DEAL:0::IS: {[0,3]} +DEAL:0::RWVector contents: +IndexSet: {[0,3]} + +[0]: 0.000e+00 +[1]: 1.000e+00 +[2]: 2.000e+00 +[3]: 3.000e+00 +DEAL:0::trilinos vec: +size:8 local_size:4 : +[0]: 0.000e+00 +[1]: 1.000e+00 +[2]: 2.000e+00 +[3]: 3.000e+00 +DEAL:0::trilinos vec (2x): +size:8 local_size:4 : +[0]: 0.000e+00 +[1]: 2.000e+00 +[2]: 4.000e+00 +[3]: 6.000e+00 +DEAL:0::trilinos vec (1x): +size:8 local_size:4 : +[0]: 0.000e+00 +[1]: 1.000e+00 +[2]: 2.000e+00 +[3]: 3.000e+00 +DEAL:0::OK + +DEAL:1::IS: {[4,7]} +DEAL:1::RWVector contents: +IndexSet: {[4,7]} + +[4]: 4.000e+00 +[5]: 5.000e+00 +[6]: 6.000e+00 +[7]: 7.000e+00 +DEAL:1::trilinos vec: +size:8 local_size:4 : +[4]: 4.000e+00 +[5]: 5.000e+00 +[6]: 6.000e+00 +[7]: 7.000e+00 +DEAL:1::trilinos vec (2x): +size:8 local_size:4 : +[4]: 8.000e+00 +[5]: 1.000e+01 +[6]: 1.200e+01 +[7]: 1.400e+01 +DEAL:1::trilinos vec (1x): +size:8 local_size:4 : +[4]: 4.000e+00 +[5]: 5.000e+00 +[6]: 6.000e+00 +[7]: 7.000e+00 +DEAL:1::OK + -- 2.39.5