From a73f66338f841aa6584355f375b60b88a1c9729a Mon Sep 17 00:00:00 2001 From: David Wells Date: Fri, 25 Jun 2021 11:57:10 -0400 Subject: [PATCH] Add a move ctor for p::d::Vector. --- include/deal.II/lac/la_parallel_vector.h | 5 + .../lac/la_parallel_vector.templates.h | 11 ++ tests/mpi/parallel_vector_12a.cc | 138 ++++++++++++++++++ tests/mpi/parallel_vector_12a.mpirun=1.output | 6 + .../mpi/parallel_vector_12a.mpirun=10.output | 6 + tests/mpi/parallel_vector_12a.mpirun=4.output | 6 + 6 files changed, 172 insertions(+) create mode 100644 tests/mpi/parallel_vector_12a.cc create mode 100644 tests/mpi/parallel_vector_12a.mpirun=1.output create mode 100644 tests/mpi/parallel_vector_12a.mpirun=10.output create mode 100644 tests/mpi/parallel_vector_12a.mpirun=4.output diff --git a/include/deal.II/lac/la_parallel_vector.h b/include/deal.II/lac/la_parallel_vector.h index bf2fcff752..842dfa0b83 100644 --- a/include/deal.II/lac/la_parallel_vector.h +++ b/include/deal.II/lac/la_parallel_vector.h @@ -282,6 +282,11 @@ namespace LinearAlgebra */ Vector(const Vector &in_vector); + /** + * Move constructor. Uses the swap method. + */ + Vector(Vector &&in_vector); + /** * Construct a parallel vector of the given global size without any * actual parallel distribution. diff --git a/include/deal.II/lac/la_parallel_vector.templates.h b/include/deal.II/lac/la_parallel_vector.templates.h index 4073de0829..7050732328 100644 --- a/include/deal.II/lac/la_parallel_vector.templates.h +++ b/include/deal.II/lac/la_parallel_vector.templates.h @@ -711,6 +711,16 @@ namespace LinearAlgebra + template + Vector::Vector(Vector &&v) + : Vector() + { + static_cast(*this) = static_cast(v); + this->swap(v); + } + + + template Vector::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 index 0000000000..aa4c2ea4e6 --- /dev/null +++ b/tests/mpi/parallel_vector_12a.cc @@ -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 +#include + +#include + +#include +#include + +#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(my_start0), + static_cast(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 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 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 index 0000000000..0eba92b76c --- /dev/null +++ b/tests/mpi/parallel_vector_12a.mpirun=1.output @@ -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 index 0000000000..687a75f512 --- /dev/null +++ b/tests/mpi/parallel_vector_12a.mpirun=10.output @@ -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 index 0000000000..9fd8b026b3 --- /dev/null +++ b/tests/mpi/parallel_vector_12a.mpirun=4.output @@ -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 -- 2.39.5