]> https://gitweb.dealii.org/ - dealii.git/commitdiff
New ReadWriteVector::reinit(TrilinosWrappers::MPI::Vector) 5015/head
authorTimo Heister <timo.heister@gmail.com>
Sun, 3 Sep 2017 18:52:54 +0000 (14:52 -0400)
committerTimo Heister <timo.heister@gmail.com>
Tue, 5 Sep 2017 18:41:26 +0000 (14:41 -0400)
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
include/deal.II/lac/read_write_vector.templates.h
tests/trilinos/readwritevector_03.cc [new file with mode: 0644]
tests/trilinos/readwritevector_03.mpirun=2.output [new file with mode: 0644]

index d8a53ac8e53c6f52911f6faaf83f52660a7d7e2a..8a89f0229a4062e8cdcbe8fb34f7afdc7ddaa007 100644 (file)
@@ -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
index ffabed60e73a227dde01c1b4eb328bfa9e5e976b..e21a8ea3b1b45e16ab7ea1b21f98fb5601dfa47c 100644 (file)
@@ -137,6 +137,35 @@ namespace LinearAlgebra
 
 
 
+#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI)
+  template <typename Number>
+  void
+  ReadWriteVector<Number>::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 <typename Number>
   template <typename Functor>
   void
diff --git a/tests/trilinos/readwritevector_03.cc b/tests/trilinos/readwritevector_03.cc
new file mode 100644 (file)
index 0000000..86ec4ac
--- /dev/null
@@ -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 <deal.II/base/index_set.h>
+#include <deal.II/lac/read_write_vector.h>
+#include <deal.II/base/utilities.h>
+#include <deal.II/lac/trilinos_vector.h>
+
+#include <vector>
+
+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<double> readwrite;
+    readwrite.reinit(tril_vector);
+
+    deallog << "RWVector contents from tril_vector:" << std::endl;
+    readwrite.print(deallog.get_file_stream());
+  }
+  {
+    LinearAlgebra::ReadWriteVector<double> readwrite;
+    readwrite.reinit(tril_vector_ghosted);
+
+    deallog << "RWVector contents from tril_vector_ghosted:" << std::endl;
+    readwrite.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_03.mpirun=2.output b/tests/trilinos/readwritevector_03.mpirun=2.output
new file mode 100644 (file)
index 0000000..0fc5b0b
--- /dev/null
@@ -0,0 +1,73 @@
+
+DEAL:0::is: {[0,3]}
+DEAL:0::is_ghosted: {[2,5]}
+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 ghosted:
+size:8 local_size:6 :
+[0]: 0.000e+00
+[1]: 1.000e+00
+[2]: 2.000e+00
+[3]: 3.000e+00
+[4]: 4.000e+00
+[5]: 5.000e+00
+DEAL:0::ghosted IS: {[0,5]}
+DEAL:0::tril_vector_ghosted.owned_elements() {[0,3]}
+DEAL:0::RWVector contents from tril_vector:
+IndexSet: {[0,3]}
+
+0.000e+00
+1.000e+00
+2.000e+00
+3.000e+00
+DEAL:0::RWVector contents from tril_vector_ghosted:
+IndexSet: {[0,5]}
+
+0.000e+00
+1.000e+00
+2.000e+00
+3.000e+00
+4.000e+00
+5.000e+00
+DEAL:0::OK
+
+DEAL:1::is: {[4,7]}
+DEAL:1::is_ghosted: {[2,5]}
+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 ghosted:
+size:8 local_size:6 :
+[2]: 2.000e+00
+[3]: 3.000e+00
+[4]: 4.000e+00
+[5]: 5.000e+00
+[6]: 6.000e+00
+[7]: 7.000e+00
+DEAL:1::ghosted IS: {[2,7]}
+DEAL:1::tril_vector_ghosted.owned_elements() {[4,7]}
+DEAL:1::RWVector contents from tril_vector:
+IndexSet: {[4,7]}
+
+4.000e+00
+5.000e+00
+6.000e+00
+7.000e+00
+DEAL:1::RWVector contents from tril_vector_ghosted:
+IndexSet: {[2,7]}
+
+2.000e+00
+3.000e+00
+4.000e+00
+5.000e+00
+6.000e+00
+7.000e+00
+DEAL:1::OK
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.