]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add TrilinosWrappers::MPI::Vector::import(ReadWriteVector) 7202/head
authorTimo Heister <timo.heister@gmail.com>
Sun, 16 Sep 2018 23:50:56 +0000 (19:50 -0400)
committerTimo Heister <timo.heister@gmail.com>
Sun, 16 Sep 2018 23:52:23 +0000 (19:52 -0400)
doc/news/changes/minor/20180916Heister [new file with mode: 0644]
include/deal.II/lac/trilinos_vector.h
source/lac/trilinos_vector.cc
tests/trilinos/readwritevector_04.cc [new file with mode: 0644]
tests/trilinos/readwritevector_04.mpirun=2.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20180916Heister b/doc/news/changes/minor/20180916Heister
new file mode 100644 (file)
index 0000000..aeb1d3a
--- /dev/null
@@ -0,0 +1,3 @@
+New: added the function TrilinosWrappers::MPI::Vector::import() to get data from a ReadWriteVector.
+<br>
+(Timo Heister, 2018/09/16)
index 865159f3e0ac93bce5e96cc388e9470a2baaa168..29a95558294d46335d518102c0d5c71b7232f61f 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
+namespace LinearAlgebra
+{
+  // Forward declaration
+  template <typename Number>
+  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<double> &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
index d0ac6a9b3d15a3a21574394e2b39aadcaa4fc861..7d084d7634e69c8ffc03747884e135aba17a885a 100644 (file)
@@ -20,6 +20,7 @@
 #  include <deal.II/base/mpi.h>
 #  include <deal.II/base/std_cxx14/memory.h>
 
+#  include <deal.II/lac/read_write_vector.h>
 #  include <deal.II/lac/trilinos_index_access.h>
 #  include <deal.II/lac/trilinos_parallel_block_vector.h>
 #  include <deal.II/lac/trilinos_sparse_matrix.h>
@@ -544,6 +545,36 @@ namespace TrilinosWrappers
     }
 
 
+    void
+    Vector::import(const LinearAlgebra::ReadWriteVector<double> &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 (file)
index 0000000..4b17c7b
--- /dev/null
@@ -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 <deal.II/base/index_set.h>
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/lac/read_write_vector.h>
+#include <deal.II/lac/trilinos_vector.h>
+
+#include <vector>
+
+#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<double> 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 (file)
index 0000000..821a244
--- /dev/null
@@ -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
+

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.