]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a move assignment operator to LA::d::V. 14969/head
authorDavid Wells <drwells@email.unc.edu>
Fri, 24 Mar 2023 12:42:44 +0000 (08:42 -0400)
committerDavid Wells <drwells@email.unc.edu>
Wed, 19 Apr 2023 11:52:28 +0000 (07:52 -0400)
include/deal.II/lac/la_parallel_vector.h
include/deal.II/lac/la_parallel_vector.templates.h
tests/lac/parallel_vector_08.cc

index afd9a24bcdf34fffad98cd3376670bcfc3f7c5d2..79fbe965b91e26bbb7fd02b28bf68b50f4b363e8 100644 (file)
@@ -443,6 +443,15 @@ namespace LinearAlgebra
       void
       swap(Vector<Number, MemorySpace> &v);
 
+      /**
+       * Move assignment operator.
+       *
+       * @note This method may throw an exception (should an MPI check fail) and
+       * is consequently not `noexcept`.
+       */
+      Vector<Number, MemorySpace> &
+      operator=(Vector<Number, MemorySpace> &&in_vector); // NOLINT
+
       /**
        * Assigns the vector to the parallel partitioning of the input vector
        * @p in_vector, and copies all the data.
index 20962b80db2fb2dc0ded41b831a7b081a18f179d..a0de1a320e7cbe7d3ea1fda012f4162e9fc563c2 100644 (file)
@@ -1375,6 +1375,18 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number, typename MemorySpaceType>
+    Vector<Number, MemorySpaceType> &
+    Vector<Number, MemorySpaceType>::operator=( // NOLINT
+      Vector<Number, MemorySpaceType> &&v)
+    {
+      static_cast<Subscriptor &>(*this) = static_cast<Subscriptor &&>(v);
+      this->swap(v);
+      return *this;
+    }
+
+
+
     template <typename Number, typename MemorySpaceType>
     Vector<Number, MemorySpaceType> &
     Vector<Number, MemorySpaceType>::operator=(const Number s)
index b9c4b6cd97d5ef537807daab05d54627da99e002..cc7d46af41bfb827af8f4e91accb90d2fb2facf7 100644 (file)
@@ -14,8 +14,8 @@
 // ---------------------------------------------------------------------
 
 
-// check parallel_vector::copy_from to update ghost values. Same vector layout
-// as in parallel_vector_07.cc
+// check parallel_vector assignment and ghost values. Same vector layout as in
+// parallel_vector_07.cc
 
 #include <deal.II/base/cuda.h>
 #include <deal.II/base/index_set.h>
@@ -58,10 +58,12 @@ test()
   IndexSet local_owned(global_size);
   local_owned.add_range(my_start, my_start + local_size);
   IndexSet local_relevant(global_size);
-  local_relevant                 = local_owned;
-  unsigned int ghost_indices[10] = {
+  local_relevant                            = local_owned;
+  const std::size_t n_ghosts                = 10;
+  unsigned int      ghost_indices[n_ghosts] = {
     1, 2, 13, set - 3, set - 2, set - 1, set, set + 1, set + 2, set + 3};
-  local_relevant.add_indices(&ghost_indices[0], &ghost_indices[0] + 10);
+  local_relevant.add_indices(std::begin(ghost_indices),
+                             std::end(ghost_indices));
 
   // v has ghosts, w has none. set some entries
   // on w, copy into v and check if they are
@@ -80,6 +82,29 @@ test()
   v = w;
   v.update_ghost_values();
 
+  // Check local and ghost values with move assignment:
+  {
+    decltype(v) v_move, v_copy = v;
+    v_move = std::move(v_copy);
+
+    LinearAlgebra::ReadWriteVector<double> v_rw(local_relevant);
+    v_rw.import(v, VectorOperation::insert);
+    LinearAlgebra::ReadWriteVector<double> v_move_rw(local_relevant);
+    v_move_rw.import(v_move, VectorOperation::insert);
+
+    for (unsigned int i = 0; i < v_rw.locally_owned_size(); ++i)
+      AssertThrow(v_move_rw.local_element(i) == v_rw.local_element(i),
+                  ExcInternalError());
+
+    // v_copy should now be empty
+    AssertThrow(v_copy.locally_owned_size() == 0, ExcInternalError());
+    AssertThrow(v_copy.get_partitioner()->size() == 0, ExcInternalError());
+    AssertThrow(v_copy.get_partitioner()->n_ghost_indices() == 0,
+                ExcInternalError());
+    AssertThrow(v_copy.get_mpi_communicator() == MPI_COMM_SELF,
+                ExcInternalError());
+  }
+
   // check local values for correctness
   rw_vector.import(v, VectorOperation::insert);
   for (unsigned int i = 0; i < local_size; ++i)

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.