From: danshapero Date: Thu, 24 Mar 2016 00:47:10 +0000 (-0700) Subject: Added move constructor to SparseMatrix X-Git-Tag: v8.5.0-rc1~1167^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ccc51e04b48085a6fe0dd21500617ca5c4f1297b;p=dealii.git Added move constructor to SparseMatrix --- diff --git a/include/deal.II/lac/sparse_matrix.h b/include/deal.II/lac/sparse_matrix.h index 1cde275e27..50ea51f48a 100644 --- a/include/deal.II/lac/sparse_matrix.h +++ b/include/deal.II/lac/sparse_matrix.h @@ -543,6 +543,20 @@ public: */ SparseMatrix (const SparseMatrix &); +#ifdef DEAL_II_WITH_CXX11 + /** + * Move constructor. Construct a new sparse matrix by transferring the + * internal data of the matrix @p m into a new object. + * + * Move construction allows an object to be returned from a function or + * packed into a tuple even when the class cannot be copy-constructed. + * + * @note this constructor is only available if deal.II is configured with + * C++11 support. + */ + SparseMatrix (SparseMatrix &&m); +#endif + /** * Constructor. Takes the given matrix sparsity structure to represent the * sparsity pattern of this matrix. You can change the sparsity pattern @@ -584,6 +598,17 @@ public: */ SparseMatrix &operator = (const SparseMatrix &); +#ifdef DEAL_II_WITH_CXX11 + /** + * Move assignment operator. This operator replaces the present matrix with + * @p m by transferring the internal data of @p m. + * + * @note This operator is only available if deal.II is configured with C++11 + * support. + */ + SparseMatrix &operator = (SparseMatrix &&m); +#endif + /** * Copy operator: initialize the matrix with the identity matrix. This * operator will throw an exception if the sizes of the sparsity pattern and diff --git a/include/deal.II/lac/sparse_matrix.templates.h b/include/deal.II/lac/sparse_matrix.templates.h index 4f145c555c..0cda654f7a 100644 --- a/include/deal.II/lac/sparse_matrix.templates.h +++ b/include/deal.II/lac/sparse_matrix.templates.h @@ -69,6 +69,23 @@ SparseMatrix::SparseMatrix (const SparseMatrix &m) +#ifdef DEAL_II_WITH_CXX11 +template +SparseMatrix::SparseMatrix (SparseMatrix &&m) + : + Subscriptor (m), + cols(m.cols), + val(m.val), + max_len(m.max_len) +{ + m.cols = 0; + m.val = 0; + m.max_len = 0; +} +#endif + + + template SparseMatrix & SparseMatrix::operator = (const SparseMatrix &m) @@ -83,6 +100,25 @@ SparseMatrix::operator = (const SparseMatrix &m) +#ifdef DEAL_II_WITH_CXX11 +template +SparseMatrix & +SparseMatrix::operator = (SparseMatrix &&m) +{ + cols = m.cols; + val = m.val; + max_len = m.max_len; + + m.cols = 0; + m.val = 0; + m.max_len = 0; + + return *this; +} +#endif + + + template SparseMatrix::SparseMatrix (const SparsityPattern &c) : diff --git a/tests/lac/sparse_matrix_move.cc b/tests/lac/sparse_matrix_move.cc new file mode 100644 index 0000000000..6d67eb3708 --- /dev/null +++ b/tests/lac/sparse_matrix_move.cc @@ -0,0 +1,100 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 1998 - 2015 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. +// +// --------------------------------------------------------------------- + + + +#include "../tests.h" +#include "testmatrix.h" +#include +#include + +//#include + +void graph_laplacian(const SparsityPattern& sparsity, + SparseMatrix& matrix) +{ + matrix = 0.0; + + for (SparsityPattern::const_iterator it = sparsity.begin(); + it != sparsity.end(); ++it) + { + const auto i = (*it).row(); + const auto j = (*it).column(); + matrix.add(i, j, -1); + matrix.add(i, i, 1); + } +} + + +SparseMatrix graph_laplacian(const SparsityPattern& sparsity) +{ + SparseMatrix A(sparsity); + graph_laplacian(sparsity, A); + + return A; +} + + +SparseMatrix +graph_laplacian_move_return(const SparsityPattern& sparsity) +{ + SparseMatrix A(sparsity); + graph_laplacian(sparsity, A); + + return std::move(A); +} + + +int main() +{ + initlog(); + deallog << std::setprecision(3); + + const unsigned int size = 5; + + FDMatrix testproblem (size, size); + unsigned int dim = (size-1) * (size-1); + + SparsityPattern sparsity(dim, dim, size); + testproblem.five_point_structure(sparsity); + sparsity.compress(); + + // Return a sparse matrix, possibly using RVO or move constructor + SparseMatrix A = graph_laplacian(sparsity); + deallog << A.n_nonzero_elements() << std::endl; + + Vector x(dim), y(dim); + x = 1.0; + y = 1.0; + A.vmult(y, x); + + deallog << y.l2_norm() << std::endl; + + // Return a sparse matrix using the move constructor + SparseMatrix B = graph_laplacian_move_return(sparsity); + deallog << B.n_nonzero_elements() << std::endl; + y = 1.0; + B.vmult(y, x); + + deallog << y.l2_norm() << std::endl; + + // Explicitly move a sparse matrix + SparseMatrix C; + C = std::move(B); + deallog << C.m() << std::endl; + deallog << B.empty() << std::endl; + + return 0; +} diff --git a/tests/lac/sparse_matrix_move.with_cxx11=on.output b/tests/lac/sparse_matrix_move.with_cxx11=on.output new file mode 100644 index 0000000000..7239c9cda8 --- /dev/null +++ b/tests/lac/sparse_matrix_move.with_cxx11=on.output @@ -0,0 +1,7 @@ + +DEAL::64 +DEAL::0.00 +DEAL::64 +DEAL::0.00 +DEAL::16 +DEAL::1