From: Matthias Maier <tamiko@43-1.org> Date: Thu, 14 Feb 2019 19:12:06 +0000 (-0600) Subject: lac/linear_operator.h: Prevent users from creating LinearOperators with temoporaries X-Git-Tag: v9.1.0-rc1~332^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=2c2d43580f7cfc06e8551b88de783eeec70fb6df;p=dealii.git lac/linear_operator.h: Prevent users from creating LinearOperators with temoporaries --- diff --git a/include/deal.II/lac/linear_operator.h b/include/deal.II/lac/linear_operator.h index f0115ff257..64dff9a099 100644 --- a/include/deal.II/lac/linear_operator.h +++ b/include/deal.II/lac/linear_operator.h @@ -1473,6 +1473,125 @@ linear_operator(const LinearOperator<Range, Domain, Payload> &operator_exemplar, //@} +#ifndef DOXYGEN + +// +// Ensure that we never capture a reference to a temporary by accident. +// This ensures that instead of silently allowing a "stack use after free", +// we at least bail out with a runtime error message that is slightly more +// understandable. +// + +template < + typename Range = Vector<double>, + typename Domain = Range, + typename Payload = internal::LinearOperatorImplementation::EmptyPayload, + typename OperatorExemplar, + typename Matrix, + typename = + typename std::enable_if<!std::is_lvalue_reference<Matrix>::value>::type> +LinearOperator<Range, Domain, Payload> +linear_operator(const OperatorExemplar &, Matrix &&) +{ + Assert(false, + ExcMessage( + "You are trying to construct a linear operator from a temporary " + "matrix object.")); +} + +template < + typename Range = Vector<double>, + typename Domain = Range, + typename Payload = internal::LinearOperatorImplementation::EmptyPayload, + typename OperatorExemplar, + typename Matrix, + typename = typename std::enable_if< + !std::is_lvalue_reference<OperatorExemplar>::value>::type> +LinearOperator<Range, Domain, Payload> +linear_operator(OperatorExemplar &&, const Matrix &) +{ + Assert(false, + ExcMessage( + "You are trying to construct a linear operator with a temporary " + "operator_exemplar object.")); +} + +template < + typename Range = Vector<double>, + typename Domain = Range, + typename Payload = internal::LinearOperatorImplementation::EmptyPayload, + typename OperatorExemplar, + typename Matrix, + typename = + typename std::enable_if<!std::is_lvalue_reference<Matrix>::value>::type, + typename = typename std::enable_if< + !std::is_lvalue_reference<OperatorExemplar>::value>::type> +LinearOperator<Range, Domain, Payload> +linear_operator(OperatorExemplar &&, Matrix &&) +{ + Assert(false, + ExcMessage( + "You are trying to construct a linear operator from a temporary " + "matrix object.")); +} + + +template < + typename Range = Vector<double>, + typename Domain = Range, + typename Payload = internal::LinearOperatorImplementation::EmptyPayload, + typename Matrix, + typename = + typename std::enable_if<!std::is_lvalue_reference<Matrix>::value>::type> +LinearOperator<Range, Domain, Payload> +linear_operator(const LinearOperator<Range, Domain, Payload> &, Matrix &&) +{ + Assert(false, + ExcMessage( + "You are trying to construct a linear operator from a temporary " + "matrix object.")); +} + +template < + typename Range = Vector<double>, + typename Domain = Range, + typename Payload = internal::LinearOperatorImplementation::EmptyPayload, + typename Matrix, + typename = + typename std::enable_if<!std::is_lvalue_reference<Matrix>::value>::type> +LinearOperator<Range, Domain, Payload> +linear_operator(Matrix &&) +{ + Assert(false, + ExcMessage( + "You are trying to construct a linear operator from a temporary " + "matrix object.")); +} + +template <typename Payload, + typename Solver, + typename Preconditioner, + typename Range = typename Solver::vector_type, + typename Domain = Range, + typename = typename std::enable_if< + !std::is_lvalue_reference<Preconditioner>::value>::type, + typename = typename std::enable_if< + !std::is_same<Preconditioner, PreconditionIdentity>::value>::type, + typename = typename std::enable_if< + !std::is_same<Preconditioner, + LinearOperator<Range, Domain, Payload>>::value>::type> +LinearOperator<Domain, Range, Payload> +inverse_operator(const LinearOperator<Range, Domain, Payload> &, + Solver &, + Preconditioner &&) +{ + Assert(false, + ExcMessage( + "You are trying to construct an inverse operator with a temporary " + "preconditioner object.")); +} + +#endif // DOXYGEN DEAL_II_NAMESPACE_CLOSE diff --git a/tests/lac/linear_operator_18.cc b/tests/lac/linear_operator_18.cc new file mode 100644 index 0000000000..08fbfb1801 --- /dev/null +++ b/tests/lac/linear_operator_18.cc @@ -0,0 +1,55 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 - 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. +// +// --------------------------------------------------------------------- + +// Test that we cannot accidentally create a linear operator from a +// temporary object + +#include <deal.II/lac/linear_operator.h> +#include <deal.II/lac/solver_cg.h> +#include <deal.II/lac/sparse_matrix.h> + +#include "../tests.h" + +#define CATCH(call) \ + try \ + { \ + call; \ + deallog << "Error: Something went wrong!" << std::endl; \ + } \ + catch (ExceptionBase & e) \ + { \ + deallog << e.get_exc_name() << std::endl; \ + } + +using namespace dealii; + +int +main() +{ + initlog(); + deal_II_exceptions::disable_abort_on_exception(); + + SparseMatrix<double> A; + const auto op_A = linear_operator(A); + + CATCH(linear_operator(SparseMatrix<double>())) + CATCH(linear_operator(SparseMatrix<double>(), A)) + CATCH(linear_operator(A, SparseMatrix<double>())) + CATCH(linear_operator(SparseMatrix<double>(), SparseMatrix<double>())) + + SolverControl solver_control; + SolverCG<Vector<double>> solver(solver_control); + CATCH(inverse_operator(op_A, solver, SparseMatrix<double>())) +} diff --git a/tests/lac/linear_operator_18.output b/tests/lac/linear_operator_18.output new file mode 100644 index 0000000000..8aad5753bc --- /dev/null +++ b/tests/lac/linear_operator_18.output @@ -0,0 +1,6 @@ + +DEAL::ExcMessage( "You are trying to construct a linear operator from a temporary " "matrix object.") +DEAL::ExcMessage( "You are trying to construct a linear operator with a temporary " "operator_exemplar object.") +DEAL::ExcMessage( "You are trying to construct a linear operator from a temporary " "matrix object.") +DEAL::ExcMessage( "You are trying to construct a linear operator from a temporary " "matrix object.") +DEAL::ExcMessage( "You are trying to construct an inverse operator with a temporary " "preconditioner object.")