diff options
author | Benoit Jacob <jacob.benoit.1@gmail.com> | 2011-11-18 11:10:27 -0500 |
---|---|---|
committer | Benoit Jacob <jacob.benoit.1@gmail.com> | 2011-11-18 11:10:27 -0500 |
commit | bc6d78982fb2b8d07827246e682acdf47d0e8944 (patch) | |
tree | be84292862f0c1fc54743a273c7c6e139a8a44c9 | |
parent | de22ad117cb5324e3e1d0700dde7466921fdf9ca (diff) |
Bugs 157 and 377 - General tightening/testing of vectorwise ops:
* add lots of static assertions making it very explicit when all these ops
are supposed to work:
** all ops require the rhs vector to go in the right direction
** all ops already require that the lhs and rhs are of the same kind
(matrix vs vector) otherwise we'd have to do complex work
** multiplicative ops (introduced Kibeom's patch) are restricted to arrays, if only because for matrices they could be ambiguous.
* add a new test, vectorwiseop.cpp.
* these compound-assign operators used to be implemented with for loops:
for(Index j=0; j<subVectors(); ++j)
subVector(j).array() += other.derived().array();
This didn't seem to be needed; replaced by using expressions like operator+ and operator- did.
-rw-r--r-- | Eigen/src/Core/VectorwiseOp.h | 54 | ||||
-rw-r--r-- | Eigen/src/Core/util/StaticAssert.h | 16 | ||||
-rw-r--r-- | test/CMakeLists.txt | 2 | ||||
-rw-r--r-- | test/vectorwiseop.cpp | 187 |
4 files changed, 235 insertions, 24 deletions
diff --git a/Eigen/src/Core/VectorwiseOp.h b/Eigen/src/Core/VectorwiseOp.h index b61234b57..7819cb489 100644 --- a/Eigen/src/Core/VectorwiseOp.h +++ b/Eigen/src/Core/VectorwiseOp.h @@ -237,7 +237,10 @@ template<typename ExpressionType, int Direction> class VectorwiseOp typename ExtendedType<OtherDerived>::Type extendedTo(const DenseBase<OtherDerived>& other) const { - EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived); + EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(Direction==Vertical, OtherDerived::MaxColsAtCompileTime==1), + YOU_PASSED_A_ROW_VECTOR_BUT_A_COLUMN_VECTOR_WAS_EXPECTED) + EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(Direction==Horizontal, OtherDerived::MaxRowsAtCompileTime==1), + YOU_PASSED_A_COLUMN_VECTOR_BUT_A_ROW_VECTOR_WAS_EXPECTED) return typename ExtendedType<OtherDerived>::Type (other.derived(), Direction==Vertical ? 1 : m_matrix.rows(), @@ -418,10 +421,9 @@ template<typename ExpressionType, int Direction> class VectorwiseOp ExpressionType& operator=(const DenseBase<OtherDerived>& other) { EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) //eigen_assert((m_matrix.isNull()) == (other.isNull())); FIXME - for(Index j=0; j<subVectors(); ++j) - subVector(j) = other; - return const_cast<ExpressionType&>(m_matrix); + return const_cast<ExpressionType&>(m_matrix = extendedTo(other.derived())); } /** Adds the vector \a other to each subvector of \c *this */ @@ -429,9 +431,8 @@ template<typename ExpressionType, int Direction> class VectorwiseOp ExpressionType& operator+=(const DenseBase<OtherDerived>& other) { EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) - for(Index j=0; j<subVectors(); ++j) - subVector(j) += other.derived(); - return const_cast<ExpressionType&>(m_matrix); + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) + return const_cast<ExpressionType&>(m_matrix += extendedTo(other.derived())); } /** Substracts the vector \a other to each subvector of \c *this */ @@ -439,28 +440,29 @@ template<typename ExpressionType, int Direction> class VectorwiseOp ExpressionType& operator-=(const DenseBase<OtherDerived>& other) { EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) - for(Index j=0; j<subVectors(); ++j) - subVector(j) -= other.derived(); - return const_cast<ExpressionType&>(m_matrix); + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) + return const_cast<ExpressionType&>(m_matrix -= extendedTo(other.derived())); } - /** Multiplies the vector \a other to each subvector of \c *this */ + /** Multiples each subvector of \c *this by the vector \a other */ template<typename OtherDerived> ExpressionType& operator*=(const DenseBase<OtherDerived>& other) { EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) - for(Index j=0; j<subVectors(); ++j) - subVector(j).array() *= other.derived().array(); + EIGEN_STATIC_ASSERT_ARRAYXPR(ExpressionType) + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) + m_matrix *= extendedTo(other.derived()); return const_cast<ExpressionType&>(m_matrix); } - /** Divides the vector \a other to each subvector of \c *this */ + /** Divides each subvector of \c *this by the vector \a other */ template<typename OtherDerived> ExpressionType& operator/=(const DenseBase<OtherDerived>& other) { EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) - for(Index j=0; j<subVectors(); ++j) - subVector(j).array() /= other.derived(); + EIGEN_STATIC_ASSERT_ARRAYXPR(ExpressionType) + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) + m_matrix /= extendedTo(other.derived()); return const_cast<ExpressionType&>(m_matrix); } @@ -471,7 +473,8 @@ template<typename ExpressionType, int Direction> class VectorwiseOp const typename ExtendedType<OtherDerived>::Type> operator+(const DenseBase<OtherDerived>& other) const { - EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived); + EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) return m_matrix + extendedTo(other.derived()); } @@ -482,29 +485,36 @@ template<typename ExpressionType, int Direction> class VectorwiseOp const typename ExtendedType<OtherDerived>::Type> operator-(const DenseBase<OtherDerived>& other) const { - EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived); + EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) return m_matrix - extendedTo(other.derived()); } - /** Returns the expression of the multiplication of the vector \a other to each subvector of \c *this */ + /** Returns the expression where each subvector is the product of the vector \a other + * by the corresponding subvector of \c *this */ template<typename OtherDerived> EIGEN_STRONG_INLINE CwiseBinaryOp<internal::scalar_product_op<Scalar>, const ExpressionTypeNestedCleaned, const typename ExtendedType<OtherDerived>::Type> operator*(const DenseBase<OtherDerived>& other) const { - EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived); + EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) + EIGEN_STATIC_ASSERT_ARRAYXPR(ExpressionType) + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) return m_matrix * extendedTo(other.derived()); } - /** Returns the expression of the division between each subvector of \c *this and the vector \a other */ + /** Returns the expression where each subvector is the quotient of the corresponding + * subvector of \c *this by the vector \a other */ template<typename OtherDerived> CwiseBinaryOp<internal::scalar_quotient_op<Scalar>, const ExpressionTypeNestedCleaned, const typename ExtendedType<OtherDerived>::Type> operator/(const DenseBase<OtherDerived>& other) const { - EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived); + EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived) + EIGEN_STATIC_ASSERT_ARRAYXPR(ExpressionType) + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(ExpressionType, OtherDerived) return m_matrix / extendedTo(other.derived()); } diff --git a/Eigen/src/Core/util/StaticAssert.h b/Eigen/src/Core/util/StaticAssert.h index d14c87c7c..672c6ae96 100644 --- a/Eigen/src/Core/util/StaticAssert.h +++ b/Eigen/src/Core/util/StaticAssert.h @@ -97,7 +97,10 @@ THIS_EXPRESSION_IS_NOT_A_LVALUE__IT_IS_READ_ONLY, YOU_ARE_TRYING_TO_USE_AN_INDEX_BASED_ACCESSOR_ON_AN_EXPRESSION_THAT_DOES_NOT_SUPPORT_THAT, THIS_METHOD_IS_ONLY_FOR_1x1_EXPRESSIONS, - THIS_METHOD_IS_ONLY_FOR_EXPRESSIONS_OF_BOOL + THIS_METHOD_IS_ONLY_FOR_EXPRESSIONS_OF_BOOL, + THIS_METHOD_IS_ONLY_FOR_ARRAYS_NOT_MATRICES, + YOU_PASSED_A_ROW_VECTOR_BUT_A_COLUMN_VECTOR_WAS_EXPECTED, + YOU_PASSED_A_COLUMN_VECTOR_BUT_A_ROW_VECTOR_WAS_EXPECTED }; }; @@ -197,4 +200,15 @@ EIGEN_STATIC_ASSERT(internal::is_lvalue<Derived>::value, \ THIS_EXPRESSION_IS_NOT_A_LVALUE__IT_IS_READ_ONLY) +#define EIGEN_STATIC_ASSERT_ARRAYXPR(Derived) \ + EIGEN_STATIC_ASSERT((internal::is_same<typename internal::traits<Derived>::XprKind, ArrayXpr>::value), \ + THIS_METHOD_IS_ONLY_FOR_ARRAYS_NOT_MATRICES) + +#define EIGEN_STATIC_ASSERT_SAME_XPR_KIND(Derived1, Derived2) \ + EIGEN_STATIC_ASSERT((internal::is_same<typename internal::traits<Derived1>::XprKind, \ + typename internal::traits<Derived2>::XprKind \ + >::value), \ + YOU_CANNOT_MIX_ARRAYS_AND_MATRICES) + + #endif // EIGEN_STATIC_ASSERT_H diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 7086f5251..e1379d78b 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -178,7 +178,7 @@ ei_add_test(dontalign) ei_add_test(evaluators) ei_add_test(sizeoverflow) ei_add_test(prec_inverse_4x4) - +ei_add_test(vectorwiseop) ei_add_test(simplicial_cholesky) ei_add_test(conjugate_gradient) diff --git a/test/vectorwiseop.cpp b/test/vectorwiseop.cpp new file mode 100644 index 000000000..d3518b7ec --- /dev/null +++ b/test/vectorwiseop.cpp @@ -0,0 +1,187 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2011 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/>. + +#define EIGEN_NO_STATIC_ASSERT + +#include "main.h" + +template<typename ArrayType> void vectorwiseop_array(const ArrayType& m) +{ + typedef typename ArrayType::Index Index; + typedef typename ArrayType::Scalar Scalar; + typedef typename NumTraits<Scalar>::Real RealScalar; + typedef Array<Scalar, ArrayType::RowsAtCompileTime, 1> ColVectorType; + typedef Array<Scalar, 1, ArrayType::ColsAtCompileTime> RowVectorType; + + Index rows = m.rows(); + Index cols = m.cols(); + Index r = internal::random<Index>(0, rows-1), + c = internal::random<Index>(0, cols-1); + + ArrayType m1 = ArrayType::Random(rows, cols), + m2(rows, cols), + m3(rows, cols); + + ColVectorType colvec = ColVectorType::Random(rows); + RowVectorType rowvec = RowVectorType::Random(cols); + + // test addition + + m2 = m1; + m2.colwise() += colvec; + VERIFY_IS_APPROX(m2, m1.colwise() + colvec); + VERIFY_IS_APPROX(m2.col(c), m1.col(c) + colvec); + + VERIFY_RAISES_ASSERT(m2.colwise() += colvec.transpose()); + VERIFY_RAISES_ASSERT(m1.colwise() + colvec.transpose()); + + m2 = m1; + m2.rowwise() += rowvec; + VERIFY_IS_APPROX(m2, m1.rowwise() + rowvec); + VERIFY_IS_APPROX(m2.row(r), m1.row(r) + rowvec); + + VERIFY_RAISES_ASSERT(m2.rowwise() += rowvec.transpose()); + VERIFY_RAISES_ASSERT(m1.rowwise() + rowvec.transpose()); + + // test substraction + + m2 = m1; + m2.colwise() -= colvec; + VERIFY_IS_APPROX(m2, m1.colwise() - colvec); + VERIFY_IS_APPROX(m2.col(c), m1.col(c) - colvec); + + VERIFY_RAISES_ASSERT(m2.colwise() -= colvec.transpose()); + VERIFY_RAISES_ASSERT(m1.colwise() - colvec.transpose()); + + m2 = m1; + m2.rowwise() -= rowvec; + VERIFY_IS_APPROX(m2, m1.rowwise() - rowvec); + VERIFY_IS_APPROX(m2.row(r), m1.row(r) - rowvec); + + VERIFY_RAISES_ASSERT(m2.rowwise() -= rowvec.transpose()); + VERIFY_RAISES_ASSERT(m1.rowwise() - rowvec.transpose()); + + // test multiplication + + m2 = m1; + m2.colwise() *= colvec; + VERIFY_IS_APPROX(m2, m1.colwise() * colvec); + VERIFY_IS_APPROX(m2.col(c), m1.col(c) * colvec); + + VERIFY_RAISES_ASSERT(m2.colwise() *= colvec.transpose()); + VERIFY_RAISES_ASSERT(m1.colwise() * colvec.transpose()); + + m2 = m1; + m2.rowwise() *= rowvec; + VERIFY_IS_APPROX(m2, m1.rowwise() * rowvec); + VERIFY_IS_APPROX(m2.row(r), m1.row(r) * rowvec); + + VERIFY_RAISES_ASSERT(m2.rowwise() *= rowvec.transpose()); + VERIFY_RAISES_ASSERT(m1.rowwise() * rowvec.transpose()); + + // test quotient + + m2 = m1; + m2.colwise() /= colvec; + VERIFY_IS_APPROX(m2, m1.colwise() / colvec); + VERIFY_IS_APPROX(m2.col(c), m1.col(c) / colvec); + + VERIFY_RAISES_ASSERT(m2.colwise() /= colvec.transpose()); + VERIFY_RAISES_ASSERT(m1.colwise() / colvec.transpose()); + + m2 = m1; + m2.rowwise() /= rowvec; + VERIFY_IS_APPROX(m2, m1.rowwise() / rowvec); + VERIFY_IS_APPROX(m2.row(r), m1.row(r) / rowvec); + + VERIFY_RAISES_ASSERT(m2.rowwise() /= rowvec.transpose()); + VERIFY_RAISES_ASSERT(m1.rowwise() / rowvec.transpose()); +} + +template<typename MatrixType> void vectorwiseop_matrix(const MatrixType& m) +{ + typedef typename MatrixType::Index Index; + typedef typename MatrixType::Scalar Scalar; + typedef typename NumTraits<Scalar>::Real RealScalar; + typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> ColVectorType; + typedef Matrix<Scalar, 1, MatrixType::ColsAtCompileTime> RowVectorType; + + Index rows = m.rows(); + Index cols = m.cols(); + Index r = internal::random<Index>(0, rows-1), + c = internal::random<Index>(0, cols-1); + + MatrixType m1 = MatrixType::Random(rows, cols), + m2(rows, cols), + m3(rows, cols); + + ColVectorType colvec = ColVectorType::Random(rows); + RowVectorType rowvec = RowVectorType::Random(cols); + + // test addition + + m2 = m1; + m2.colwise() += colvec; + VERIFY_IS_APPROX(m2, m1.colwise() + colvec); + VERIFY_IS_APPROX(m2.col(c), m1.col(c) + colvec); + + VERIFY_RAISES_ASSERT(m2.colwise() += colvec.transpose()); + VERIFY_RAISES_ASSERT(m1.colwise() + colvec.transpose()); + + m2 = m1; + m2.rowwise() += rowvec; + VERIFY_IS_APPROX(m2, m1.rowwise() + rowvec); + VERIFY_IS_APPROX(m2.row(r), m1.row(r) + rowvec); + + VERIFY_RAISES_ASSERT(m2.rowwise() += rowvec.transpose()); + VERIFY_RAISES_ASSERT(m1.rowwise() + rowvec.transpose()); + + // test substraction + + m2 = m1; + m2.colwise() -= colvec; + VERIFY_IS_APPROX(m2, m1.colwise() - colvec); + VERIFY_IS_APPROX(m2.col(c), m1.col(c) - colvec); + + VERIFY_RAISES_ASSERT(m2.colwise() -= colvec.transpose()); + VERIFY_RAISES_ASSERT(m1.colwise() - colvec.transpose()); + + m2 = m1; + m2.rowwise() -= rowvec; + VERIFY_IS_APPROX(m2, m1.rowwise() - rowvec); + VERIFY_IS_APPROX(m2.row(r), m1.row(r) - rowvec); + + VERIFY_RAISES_ASSERT(m2.rowwise() -= rowvec.transpose()); + VERIFY_RAISES_ASSERT(m1.rowwise() - rowvec.transpose()); +} + +void test_vectorwiseop() +{ + CALL_SUBTEST_1(vectorwiseop_array(Array22cd())); + CALL_SUBTEST_2(vectorwiseop_array(Array<double, 3, 2>())); + CALL_SUBTEST_3(vectorwiseop_array(ArrayXXf(3, 4))); + CALL_SUBTEST_4(vectorwiseop_matrix(Matrix4cf())); + CALL_SUBTEST_5(vectorwiseop_matrix(Matrix<float,4,5>())); + CALL_SUBTEST_6(vectorwiseop_matrix(MatrixXd(7,2))); +} |