]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a move ctor for p::d::Vector.
authorDavid Wells <drwells@email.unc.edu>
Fri, 25 Jun 2021 15:57:10 +0000 (11:57 -0400)
committerDavid Wells <drwells@email.unc.edu>
Sat, 26 Jun 2021 13:04:46 +0000 (09:04 -0400)
include/deal.II/lac/la_parallel_vector.h
include/deal.II/lac/la_parallel_vector.templates.h
tests/mpi/parallel_vector_12a.cc [new file with mode: 0644]
tests/mpi/parallel_vector_12a.mpirun=1.output [new file with mode: 0644]
tests/mpi/parallel_vector_12a.mpirun=10.output [new file with mode: 0644]
tests/mpi/parallel_vector_12a.mpirun=4.output [new file with mode: 0644]

index bf2fcff752d39867d14a83e1b1941408d511c156..842dfa0b839a12ac6e461d4f99f5db2ac035f3b6 100644 (file)
@@ -282,6 +282,11 @@ namespace LinearAlgebra
        */
       Vector(const Vector<Number, MemorySpace> &in_vector);
 
+      /**
+       * Move constructor. Uses the swap method.
+       */
+      Vector(Vector<Number, MemorySpace> &&in_vector);
+
       /**
        * Construct a parallel vector of the given global size without any
        * actual parallel distribution.
index 4073de08297ce491ebcb741e954a294527826ed7..7050732328d479c2b0e0d1fb7610ade5a7d3eb44 100644 (file)
@@ -711,6 +711,16 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number, typename MemorySpaceType>
+    Vector<Number, MemorySpaceType>::Vector(Vector<Number, MemorySpaceType> &&v)
+      : Vector()
+    {
+      static_cast<Subscriptor &>(*this) = static_cast<Subscriptor &&>(v);
+      this->swap(v);
+    }
+
+
+
     template <typename Number, typename MemorySpaceType>
     Vector<Number, MemorySpaceType>::Vector(const IndexSet &local_range,
                                             const IndexSet &ghost_indices,
@@ -1411,6 +1421,7 @@ namespace LinearAlgebra
 
       std::swap(compress_requests, v.compress_requests);
       std::swap(update_ghost_values_requests, v.update_ghost_values_requests);
+      std::swap(comm_sm, v.comm_sm);
 #endif
 
       std::swap(partitioner, v.partitioner);
diff --git a/tests/mpi/parallel_vector_12a.cc b/tests/mpi/parallel_vector_12a.cc
new file mode 100644 (file)
index 0000000..aa4c2ea
--- /dev/null
@@ -0,0 +1,138 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2012 - 2018 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.
+//
+// ---------------------------------------------------------------------
+
+
+// check LinearAlgebra::distributed::Vector's move constructor. Should work like
+// swap.
+
+#include <deal.II/base/index_set.h>
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/lac/la_parallel_vector.h>
+
+#include <iostream>
+#include <vector>
+
+#include "../tests.h"
+
+DeclException2(ExcNonEqual,
+               double,
+               double,
+               << "Left compare: " << arg1 << ", right compare: " << arg2);
+
+void
+test()
+{
+  unsigned int myid    = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+  unsigned int numproc = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
+
+  if (myid == 0)
+    deallog << "numproc=" << numproc << std::endl;
+
+  // vector 0:
+  // global size: 20, locally_owned_size: 3 as long as
+  // less than 20
+  const unsigned int locally_owned_size0 = 3;
+  const unsigned int global_size0 =
+    std::min(20U, locally_owned_size0 * numproc);
+  const unsigned int my_start0 =
+    std::min(locally_owned_size0 * myid, global_size0);
+  const unsigned int my_end0 =
+    std::min(locally_owned_size0 * (myid + 1), global_size0);
+  const unsigned int actual_locally_owned_size0 = my_end0 - my_start0;
+
+  IndexSet local_owned0(global_size0);
+  if (my_end0 > my_start0)
+    local_owned0.add_range(static_cast<unsigned int>(my_start0),
+                           static_cast<unsigned int>(my_end0));
+  IndexSet local_relevant0(global_size0);
+  local_relevant0 = local_owned0;
+  local_relevant0.add_index(2);
+  if (numproc > 2)
+    local_relevant0.add_index(8);
+
+  LinearAlgebra::distributed::Vector<double> v0(local_owned0,
+                                                local_relevant0,
+                                                MPI_COMM_WORLD);
+  v0 = 1;
+  // check assignment in initial state
+  for (unsigned int i = 0; i < v0.locally_owned_size(); ++i)
+    AssertThrow(v0.local_element(i) == 1.,
+                ExcNonEqual(v0.local_element(i), 1.));
+
+  // check ghost elements in initial state
+  v0.update_ghost_values();
+  AssertThrow(v0(2) == 1., ExcNonEqual(v0(2), 1.));
+  if (numproc > 2)
+    AssertThrow(v0(8) == 1., ExcNonEqual(v0(8), 1.));
+  MPI_Barrier(MPI_COMM_WORLD);
+  if (myid == 0)
+    deallog << "Initial set and ghost update OK" << std::endl;
+
+  // now move.
+  LinearAlgebra::distributed::Vector<double> v1 = std::move(v0);
+  // Make sure we actually moved and not copied
+  AssertDimension(v0.locally_owned_size(), 0);
+  AssertDimension(v1.locally_owned_size(), actual_locally_owned_size0);
+  AssertDimension(v0.size(), 0);
+  AssertDimension(v1.size(), global_size0);
+  MPI_Barrier(MPI_COMM_WORLD);
+  if (myid == 0)
+    deallog << "First move: dimensions OK" << std::endl;
+  for (unsigned int i = 0; i < actual_locally_owned_size0; ++i)
+    AssertThrow(v1.local_element(i) == 1.,
+                ExcNonEqual(v1.local_element(i), 1.));
+  // Since we moved the ghost values should be present
+  for (const auto &ghost_index : v1.get_partitioner()->ghost_indices())
+    {
+      AssertThrow(v1(ghost_index) == 1., ExcNonEqual(v1(ghost_index), 1.));
+    }
+
+  MPI_Barrier(MPI_COMM_WORLD);
+  if (myid == 0)
+    deallog << "First move: local values OK" << std::endl;
+  v0.update_ghost_values();
+  v1.update_ghost_values();
+  AssertThrow(v1(2) == 1., ExcNonEqual(v1(2), 1.));
+  if (numproc > 2)
+    AssertThrow(v1(8) == 1., ExcNonEqual(v1(8), 1.));
+  if (numproc > 2)
+    AssertThrow(v1(8) == 1., ExcNonEqual(v1(8), 1.));
+  MPI_Barrier(MPI_COMM_WORLD);
+  if (myid == 0)
+    deallog << "Ghost values after first move OK" << std::endl;
+}
+
+
+
+int
+main(int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(
+    argc, argv, testing_max_num_threads());
+
+  unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+  deallog.push(Utilities::int_to_string(myid));
+
+  if (myid == 0)
+    {
+      initlog();
+      deallog << std::setprecision(4);
+
+      test();
+    }
+  else
+    test();
+}
diff --git a/tests/mpi/parallel_vector_12a.mpirun=1.output b/tests/mpi/parallel_vector_12a.mpirun=1.output
new file mode 100644 (file)
index 0000000..0eba92b
--- /dev/null
@@ -0,0 +1,6 @@
+
+DEAL:0::numproc=1
+DEAL:0::Initial set and ghost update OK
+DEAL:0::First move: dimensions OK
+DEAL:0::First move: local values OK
+DEAL:0::Ghost values after first move OK
diff --git a/tests/mpi/parallel_vector_12a.mpirun=10.output b/tests/mpi/parallel_vector_12a.mpirun=10.output
new file mode 100644 (file)
index 0000000..687a75f
--- /dev/null
@@ -0,0 +1,6 @@
+
+DEAL:0::numproc=10
+DEAL:0::Initial set and ghost update OK
+DEAL:0::First move: dimensions OK
+DEAL:0::First move: local values OK
+DEAL:0::Ghost values after first move OK
diff --git a/tests/mpi/parallel_vector_12a.mpirun=4.output b/tests/mpi/parallel_vector_12a.mpirun=4.output
new file mode 100644 (file)
index 0000000..9fd8b02
--- /dev/null
@@ -0,0 +1,6 @@
+
+DEAL:0::numproc=4
+DEAL:0::Initial set and ghost update OK
+DEAL:0::First move: dimensions OK
+DEAL:0::First move: local values OK
+DEAL:0::Ghost values after first move 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.