diff options
-rw-r--r-- | Eigen/src/Core/Matrix.h | 34 | ||||
-rw-r--r-- | Eigen/src/Core/SelfAdjointView.h | 3 | ||||
-rw-r--r-- | Eigen/src/Core/products/TriangularMatrixVector.h | 3 | ||||
-rw-r--r-- | test/CMakeLists.txt | 1 | ||||
-rw-r--r-- | test/product_notemporary.cpp | 119 | ||||
-rw-r--r-- | test/product_symm.cpp | 2 | ||||
-rw-r--r-- | test/product_trmm.cpp | 2 | ||||
-rw-r--r-- | test/triangular.cpp | 2 |
8 files changed, 147 insertions, 19 deletions
diff --git a/Eigen/src/Core/Matrix.h b/Eigen/src/Core/Matrix.h index 2b4c4634a..d5c508128 100644 --- a/Eigen/src/Core/Matrix.h +++ b/Eigen/src/Core/Matrix.h @@ -360,7 +360,7 @@ class Matrix /** \internal */ Matrix(ei_constructor_without_unaligned_array_assert) : m_storage(ei_constructor_without_unaligned_array_assert()) - {} + { _check_template_params(); } #endif /** Constructs a vector or row-vector with given dimension. \only_for_vectors @@ -436,7 +436,10 @@ class Matrix /** Copy constructor with in-place evaluation */ template<typename OtherDerived,typename OtherEvalType> EIGEN_STRONG_INLINE Matrix(const ReturnByValue<OtherDerived,OtherEvalType>& other) - { other.evalTo(*this); } + { + _check_template_params(); + other.evalTo(*this); + } /** Destructor */ inline ~Matrix() {} @@ -454,6 +457,7 @@ class Matrix EIGEN_STRONG_INLINE Matrix(const AnyMatrixBase<OtherDerived> &other) : m_storage(other.derived().rows() * other.derived().cols(), other.derived().rows(), other.derived().cols()) { + _check_template_params(); *this = other; } @@ -587,17 +591,21 @@ class Matrix static EIGEN_STRONG_INLINE void _check_template_params() { - EIGEN_STATIC_ASSERT(((_Rows >= _MaxRows) - && (_Cols >= _MaxCols) - && (_MaxRows >= 0) - && (_MaxCols >= 0) - && (_Rows <= Dynamic) - && (_Cols <= Dynamic) - && (_MaxRows == _Rows || _Rows==Dynamic) - && (_MaxCols == _Cols || _Cols==Dynamic) - && ((_MaxRows==Dynamic?1:_MaxRows)*(_MaxCols==Dynamic?1:_MaxCols)<Dynamic) - && (_Options & (DontAlign|RowMajor)) == _Options), - INVALID_MATRIX_TEMPLATE_PARAMETERS) + #ifdef EIGEN_DEBUG_MATRIX_CTOR + EIGEN_DEBUG_MATRIX_CTOR(Matrix); + #endif + + EIGEN_STATIC_ASSERT(((_Rows >= _MaxRows) + && (_Cols >= _MaxCols) + && (_MaxRows >= 0) + && (_MaxCols >= 0) + && (_Rows <= Dynamic) + && (_Cols <= Dynamic) + && (_MaxRows == _Rows || _Rows==Dynamic) + && (_MaxCols == _Cols || _Cols==Dynamic) + && ((_MaxRows==Dynamic?1:_MaxRows)*(_MaxCols==Dynamic?1:_MaxCols)<Dynamic) + && (_Options & (DontAlign|RowMajor)) == _Options), + INVALID_MATRIX_TEMPLATE_PARAMETERS) } diff --git a/Eigen/src/Core/SelfAdjointView.h b/Eigen/src/Core/SelfAdjointView.h index 3faebfc5d..5644e2fab 100644 --- a/Eigen/src/Core/SelfAdjointView.h +++ b/Eigen/src/Core/SelfAdjointView.h @@ -306,13 +306,14 @@ struct ei_selfadjoint_product_returntype<Lhs,LhsMode,false,Rhs,RhsMode,false> template<typename Dest> void evalTo(Dest& dst) const { - dst.resize(m_lhs.rows(), m_rhs.cols()); dst.setZero(); evalTo(dst,1); } template<typename Dest> void evalTo(Dest& dst, Scalar alpha) const { + ei_assert(dst.rows()==m_lhs.rows() && dst.cols()==m_rhs.cols()); + const ActualLhsType lhs = LhsBlasTraits::extract(m_lhs); const ActualRhsType rhs = RhsBlasTraits::extract(m_rhs); diff --git a/Eigen/src/Core/products/TriangularMatrixVector.h b/Eigen/src/Core/products/TriangularMatrixVector.h index c557e5237..de65b1ece 100644 --- a/Eigen/src/Core/products/TriangularMatrixVector.h +++ b/Eigen/src/Core/products/TriangularMatrixVector.h @@ -148,13 +148,14 @@ struct ei_triangular_product_returntype<Mode,true,Lhs,false,Rhs,true> template<typename Dest> void evalTo(Dest& dst) const { - dst.resize(m_lhs.rows(), m_rhs.cols()); dst.setZero(); evalTo(dst,1); } template<typename Dest> void evalTo(Dest& dst, Scalar alpha) const { + ei_assert(dst.rows()==m_lhs.rows() && dst.cols()==m_rhs.cols()); + const ActualLhsType lhs = LhsBlasTraits::extract(m_lhs); const ActualRhsType rhs = RhsBlasTraits::extract(m_rhs); diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 99224ff60..2f1aa21a3 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -115,6 +115,7 @@ ei_add_test(product_syrk ${EI_OFLAG}) ei_add_test(product_trmv ${EI_OFLAG}) ei_add_test(product_trmm ${EI_OFLAG}) ei_add_test(product_trsm ${EI_OFLAG}) +ei_add_test(product_notemporary ${EI_OFLAG}) ei_add_test(bandmatrix) ei_add_test(cholesky " " "${GSL_LIBRARIES}") ei_add_test(lu ${EI_OFLAG}) diff --git a/test/product_notemporary.cpp b/test/product_notemporary.cpp new file mode 100644 index 000000000..478bc2521 --- /dev/null +++ b/test/product_notemporary.cpp @@ -0,0 +1,119 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1@gmail.com> +// +// Eigen is free software; you can 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 3 of the License, or (at your option) any later version. +// +// Alternatively, you can redistribute it and/or +// modify it under the terms of the GNU General Public License as +// published by the Free Software Foundation; either version 2 of +// the License, or (at your option) any later version. +// +// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY +// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS +// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the +// GNU General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License and a copy of the GNU General Public License along with +// Eigen. If not, see <http://www.gnu.org/licenses/>. + +static int nb_temporaries; + +#define EIGEN_DEBUG_MATRIX_CTOR(MTYPE) { \ + if(MTYPE::SizeAtCompileTime==Dynamic) \ + nb_temporaries++; \ + } + +#include "main.h" +#include <Eigen/Array> + +#define VERIFY_EVALUATION_COUNT(XPR,N) {\ + nb_temporaries = 0; \ + XPR; \ + if(nb_temporaries!=N) std::cerr << "nb_temporaries == " << nb_temporaries << "\n"; \ + VERIFY( (#XPR) && nb_temporaries==N ); \ + } + +template<typename MatrixType> void product_notemporary(const MatrixType& m) +{ + /* This test checks the number of tempories created + * during the evaluation of a complex expression */ + + typedef typename MatrixType::Scalar Scalar; + typedef Matrix<Scalar, 1, Dynamic> RowVectorType; + typedef Matrix<Scalar, Dynamic, 1> ColVectorType; + typedef Matrix<Scalar, Dynamic, Dynamic, RowMajor> RowMajorMatrixType; + + int rows = m.rows(); + int cols = m.cols(); + + MatrixType m1 = MatrixType::Random(rows, cols), + m2 = MatrixType::Random(rows, cols), + m3(rows, cols); + RowVectorType rv1 = RowVectorType::Random(rows), rvres(rows); + ColVectorType vc2 = ColVectorType::Random(cols), cvres(cols); + RowMajorMatrixType rm3(rows, cols); + + Scalar s1 = ei_random<Scalar>(), + s2 = ei_random<Scalar>(), + s3 = ei_random<Scalar>(); + + int c0 = ei_random<int>(4,cols-8), + c1 = ei_random<int>(8,cols-c0), + r0 = ei_random<int>(4,cols-8), + r1 = ei_random<int>(8,rows-r0); + + VERIFY_EVALUATION_COUNT( m3 = (m1 * m2.adjoint()), 1); + VERIFY_EVALUATION_COUNT( m3 = (m1 * m2.adjoint()).lazy(), 0); + // NOTE in this case the slow product is used: + VERIFY_EVALUATION_COUNT( m3 = s1 * (m1 * m2.transpose()).lazy(), 0); + VERIFY_EVALUATION_COUNT( m3 = (s1 * m1 * s2 * m2.adjoint()).lazy(), 0); + VERIFY_EVALUATION_COUNT( m3 = (s1 * m1 * s2 * (m1*s3+m2*s2).adjoint()).lazy(), 1); + VERIFY_EVALUATION_COUNT( m3 = ((s1 * m1).adjoint() * s2 * m2).lazy(), 0); + VERIFY_EVALUATION_COUNT( m3 -= (s1 * (-m1*s3).adjoint() * (s2 * m2 * s3)).lazy(), 0); + VERIFY_EVALUATION_COUNT( m3 -= (s1 * (m1.transpose() * m2)).lazy(), 1); + + VERIFY_EVALUATION_COUNT(( m3.block(r0,r0,r1,r1) += (-m1.block(r0,c0,r1,c1) * (s2*m2.block(r0,c0,r1,c1)).adjoint()).lazy() ), 0); + VERIFY_EVALUATION_COUNT(( m3.block(r0,r0,r1,r1) -= (s1 * m1.block(r0,c0,r1,c1) * m2.block(c0,r0,c1,r1)).lazy() ), 0); + // NOTE this is because the Block expression is not handled yet by our expression analyser + VERIFY_EVALUATION_COUNT(( m3.block(r0,r0,r1,r1) = (s1 * m1.block(r0,c0,r1,c1) * (s1*m2).block(c0,r0,c1,r1)).lazy() ), 1); + + VERIFY_EVALUATION_COUNT( m3 -= (s1 * m1).template triangularView<LowerTriangular>() * m2, 0); + VERIFY_EVALUATION_COUNT( rm3 = (s1 * m1.adjoint()).template triangularView<UpperTriangular>() * (m2+m2), 1); + VERIFY_EVALUATION_COUNT( rm3 = (s1 * m1.adjoint()).template triangularView<UnitUpperTriangular>() * m2.adjoint(), 0); + + VERIFY_EVALUATION_COUNT( rm3.col(c0) = (s1 * m1.adjoint()).template triangularView<UnitUpperTriangular>() * (s2*m2.row(c0)).adjoint(), 0); + + VERIFY_EVALUATION_COUNT( m1.template triangularView<LowerTriangular>().solveInPlace(m3), 0); + // FIXME this is because the rhs/result must be column major: + VERIFY_EVALUATION_COUNT( m1.adjoint().template triangularView<LowerTriangular>().solveInPlace(m3.transpose()), 1); + + VERIFY_EVALUATION_COUNT( m3 -= (s1 * m1).adjoint().template selfadjointView<LowerTriangular>() * (-m2*s3).adjoint(), 0); + VERIFY_EVALUATION_COUNT( m3 = s2 * m2.adjoint() * (s1 * m1.adjoint()).template selfadjointView<UpperTriangular>(), 0); + VERIFY_EVALUATION_COUNT( rm3 = (s1 * m1.adjoint()).template selfadjointView<LowerTriangular>() * m2.adjoint(), 0); + + VERIFY_EVALUATION_COUNT( m3.col(c0) = (s1 * m1).adjoint().template selfadjointView<LowerTriangular>() * (-m2.row(c0)*s3).adjoint(), 0); + VERIFY_EVALUATION_COUNT( m3.col(c0) -= (s1 * m1).adjoint().template selfadjointView<UpperTriangular>() * (-m2.row(c0)*s3).adjoint(), 0); + + VERIFY_EVALUATION_COUNT(( m3.block(r0,r0,r1,r1) += m1.block(r0,r0,r1,r1).template selfadjointView<UpperTriangular>() * (s1*m2.block(c0,r0,c1,r1)) ), 0); + VERIFY_EVALUATION_COUNT(( m3.block(r0,r0,r1,r1) = m1.block(r0,r0,r1,r1).template selfadjointView<UpperTriangular>() * m2.block(c0,r0,c1,r1) ), 0); + VERIFY_EVALUATION_COUNT(( m3 = m1.block(r0,r0,r1,r1).template selfadjointView<LowerTriangular>() * m2.block(c0,r0,c1,r1) ), 0); + + VERIFY_EVALUATION_COUNT( m3.template selfadjointView<LowerTriangular>().rankUpdate(m2.adjoint()), 0); +} + +void test_product_notemporary() +{ + int s; + for(int i = 0; i < g_repeat; i++) { + s = ei_random<int>(16,320); + CALL_SUBTEST( product_notemporary(MatrixXf(s, s)) ); + s = ei_random<int>(16,120); + CALL_SUBTEST( product_notemporary(MatrixXcd(s,s)) ); + } +} diff --git a/test/product_symm.cpp b/test/product_symm.cpp index 54bf91fb9..88bac878b 100644 --- a/test/product_symm.cpp +++ b/test/product_symm.cpp @@ -38,7 +38,7 @@ template<int OtherSize> struct symm_extra { template<> struct symm_extra<1> { template<typename M1, typename M2, typename Scalar> - static void run(M1& m1, M1& m2, M2& rhs2, M2& rhs22, M2& rhs23, Scalar s1, Scalar s2) {} + static void run(M1&, M1&, M2&, M2&, M2&, Scalar, Scalar) {} }; template<typename Scalar, int Size, int OtherSize> void symm(int size = Size, int othersize = OtherSize) diff --git a/test/product_trmm.cpp b/test/product_trmm.cpp index 47ffb4af3..734d8c970 100644 --- a/test/product_trmm.cpp +++ b/test/product_trmm.cpp @@ -53,7 +53,7 @@ template<typename Scalar> void trmm(int size,int othersize) VERIFY_IS_APPROX(rge3 = tri.template triangularView<LowerTriangular>() * ge2.adjoint(), loTri * ge2.adjoint()); VERIFY_IS_APPROX( ge3 = tri.template triangularView<UpperTriangular>() * ge2.adjoint(), upTri * ge2.adjoint()); VERIFY_IS_APPROX(rge3 = tri.template triangularView<UpperTriangular>() * ge2.adjoint(), upTri * ge2.adjoint()); - VERIFY_IS_APPROX( ge3 = tri.adjoint().template triangularView<UpperTriangular>() * ge2.adjoint(), loTri.adjoint() * ge2.adjoint()); + VERIFY_IS_APPROX( ge3 = (s1*tri).adjoint().template triangularView<UpperTriangular>() * ge2.adjoint(), ei_conj(s1) * loTri.adjoint() * ge2.adjoint()); VERIFY_IS_APPROX(rge3 = tri.adjoint().template triangularView<UpperTriangular>() * ge2.adjoint(), loTri.adjoint() * ge2.adjoint()); VERIFY_IS_APPROX( ge3 = tri.adjoint().template triangularView<LowerTriangular>() * ge2.adjoint(), upTri.adjoint() * ge2.adjoint()); VERIFY_IS_APPROX(rge3 = tri.adjoint().template triangularView<LowerTriangular>() * ge2.adjoint(), upTri.adjoint() * ge2.adjoint()); diff --git a/test/triangular.cpp b/test/triangular.cpp index 6385bffd1..39430034d 100644 --- a/test/triangular.cpp +++ b/test/triangular.cpp @@ -51,8 +51,6 @@ template<typename MatrixType> void triangular(const MatrixType& m) v2 = VectorType::Random(rows), vzero = VectorType::Zero(rows); - Scalar s1 = ei_random<Scalar>(); - MatrixType m1up = m1.template triangularView<Eigen::UpperTriangular>(); MatrixType m2up = m2.template triangularView<Eigen::UpperTriangular>(); |