From dad245af56b4597f27d7c219e5f5f0b38cd62d24 Mon Sep 17 00:00:00 2001 From: Benoit Jacob Date: Tue, 25 Dec 2007 17:20:58 +0000 Subject: [PATCH] - eigen2 now fully enforces constness! found a way to achieve that with minimal code duplication. There now are only two (2) const_cast remaining in the whole source code. - eigen2 now fully allows copying a row-vector into a column-vector. added a unit-test for that. - split unit tests, improve docs, various improvements. --- Eigen/Core/Block.h | 17 ++++-- Eigen/Core/Coeffs.h | 8 +-- Eigen/Core/Column.h | 14 ++++- Eigen/Core/DiagonalCoeffs.h | 11 +++- Eigen/Core/DiagonalMatrix.h | 2 +- Eigen/Core/Dot.h | 2 +- Eigen/Core/DynBlock.h | 14 ++++- Eigen/Core/Fuzzy.h | 6 +- Eigen/Core/Map.h | 12 ++-- Eigen/Core/Matrix.h | 15 ++++- Eigen/Core/MatrixBase.h | 40 ++++++++++---- Eigen/Core/MatrixRef.h | 6 +- Eigen/Core/Minor.h | 10 +++- Eigen/Core/Ones.h | 2 +- Eigen/Core/OperatorEquals.h | 2 +- Eigen/Core/Random.h | 2 +- Eigen/Core/Row.h | 10 +++- Eigen/Core/Transpose.h | 10 +++- Eigen/Core/Zero.h | 2 +- doc/Doxyfile.in | 2 +- test/CMakeLists.txt | 3 + test/basicstuff.cpp | 96 +++----------------------------- test/linearstructure.cpp | 106 ++++++++++++++++++++++++++++++++++++ test/main.h | 3 + test/map.cpp | 60 ++++++++++++++++++++ test/product.cpp | 100 ++++++++++++++++++++++++++++++++++ 26 files changed, 414 insertions(+), 141 deletions(-) create mode 100644 test/linearstructure.cpp create mode 100644 test/map.cpp create mode 100644 test/product.cpp diff --git a/Eigen/Core/Block.h b/Eigen/Core/Block.h index 7dbd14e5c..eb80fc790 100644 --- a/Eigen/Core/Block.h +++ b/Eigen/Core/Block.h @@ -95,8 +95,9 @@ template class Block /** \returns a fixed-size expression of a block in *this. * - * \param blockRows the number of rows in the block - * \param blockCols the number of columns in the block + * The template parameters \a blockRows and \a blockCols are the number of + * rows and columns in the block + * * \param startRow the first row in the block * \param startCol the first column in the block * @@ -108,10 +109,18 @@ template class Block template template Block MatrixBase + ::block(int startRow, int startCol) +{ + return Block(ref(), startRow, startCol); +} + +/** This is the const version of block(). */ +template +template +const Block MatrixBase ::block(int startRow, int startCol) const { - return Block - (static_cast(const_cast(this))->ref(), startRow, startCol); + return Block(ref(), startRow, startCol); } #endif // EIGEN_BLOCK_H diff --git a/Eigen/Core/Coeffs.h b/Eigen/Core/Coeffs.h index fa8f44935..6c48a40c5 100644 --- a/Eigen/Core/Coeffs.h +++ b/Eigen/Core/Coeffs.h @@ -116,7 +116,7 @@ template Scalar MatrixBase ::coeff(int index) const { - eigen_internal_assert(IsVector); + eigen_internal_assert(IsVectorAtCompileTime); if(RowsAtCompileTime == 1) { eigen_internal_assert(index >= 0 && index < cols()); @@ -140,7 +140,7 @@ template Scalar MatrixBase ::operator[](int index) const { - assert(IsVector); + assert(IsVectorAtCompileTime); if(RowsAtCompileTime == 1) { assert(index >= 0 && index < cols()); @@ -171,7 +171,7 @@ template Scalar& MatrixBase ::coeffRef(int index) { - eigen_internal_assert(IsVector); + eigen_internal_assert(IsVectorAtCompileTime); if(RowsAtCompileTime == 1) { eigen_internal_assert(index >= 0 && index < cols()); @@ -194,7 +194,7 @@ template Scalar& MatrixBase ::operator[](int index) { - assert(IsVector); + assert(IsVectorAtCompileTime); if(RowsAtCompileTime == 1) { assert(index >= 0 && index < cols()); diff --git a/Eigen/Core/Column.h b/Eigen/Core/Column.h index a50be7114..e0da29638 100644 --- a/Eigen/Core/Column.h +++ b/Eigen/Core/Column.h @@ -89,15 +89,23 @@ template class Column /** \returns an expression of the \a i-th column of *this. Note that the numbering starts at 0. * - * Example: \include MatrixBase_col.cpp - * Output: \verbinclude MatrixBase_col.out + * Example: \include MatrixBase_column.cpp + * Output: \verbinclude MatrixBase_column.out * * \sa row(), class Column */ template Column +MatrixBase::col(int i) +{ + return Column(ref(), i); +} + +/** This is the const version of col(). */ +template +const Column MatrixBase::col(int i) const { - return Column(static_cast(const_cast(this))->ref(), i); + return Column(ref(), i); } #endif // EIGEN_COLUMN_H diff --git a/Eigen/Core/DiagonalCoeffs.h b/Eigen/Core/DiagonalCoeffs.h index b20e188b7..c8a0d43ed 100644 --- a/Eigen/Core/DiagonalCoeffs.h +++ b/Eigen/Core/DiagonalCoeffs.h @@ -64,10 +64,17 @@ template class DiagonalCoeffs template DiagonalCoeffs +MatrixBase::diagonal() +{ + return DiagonalCoeffs(ref()); +} + +/** This is the const version of diagonal(). */ +template +const DiagonalCoeffs MatrixBase::diagonal() const { - return DiagonalCoeffs - (static_cast(const_cast(this))->ref()); + return DiagonalCoeffs(ref()); } #endif // EIGEN_DIAGONALCOEFFS_H diff --git a/Eigen/Core/DiagonalMatrix.h b/Eigen/Core/DiagonalMatrix.h index 0ae9399e8..125fa58d9 100644 --- a/Eigen/Core/DiagonalMatrix.h +++ b/Eigen/Core/DiagonalMatrix.h @@ -38,7 +38,7 @@ class DiagonalMatrix : NoOperatorEquals, DiagonalMatrix(const CoeffsVecRef& coeffs) : m_coeffs(coeffs) { - assert(CoeffsVectorType::IsVector + assert(CoeffsVectorType::IsVectorAtCompileTime && _RowsAtCompileTime == _ColsAtCompileTime && _RowsAtCompileTime == CoeffsVectorType::SizeAtCompileTime && coeffs.size() > 0); diff --git a/Eigen/Core/Dot.h b/Eigen/Core/Dot.h index 1c0573a9a..ff6852af2 100644 --- a/Eigen/Core/Dot.h +++ b/Eigen/Core/Dot.h @@ -62,7 +62,7 @@ template template Scalar MatrixBase::dot(const OtherDerived& other) const { - assert(IsVector && OtherDerived::IsVector && size() == other.size()); + assert(IsVectorAtCompileTime && OtherDerived::IsVectorAtCompileTime && size() == other.size()); Scalar res; if(EIGEN_UNROLLED_LOOPS && SizeAtCompileTime != Dynamic && SizeAtCompileTime <= 16) DotUnroller diff --git a/Eigen/Core/DynBlock.h b/Eigen/Core/DynBlock.h index 2fd66746b..c06a33b07 100644 --- a/Eigen/Core/DynBlock.h +++ b/Eigen/Core/DynBlock.h @@ -109,10 +109,18 @@ template class DynBlock */ template DynBlock MatrixBase - ::dynBlock(int startRow, int startCol, int blockRows, int blockCols) const + ::dynBlock(int startRow, int startCol, int blockRows, int blockCols) { - return DynBlock(static_cast(const_cast(this))->ref(), - startRow, startCol, blockRows, blockCols); + return DynBlock(ref(), startRow, startCol, blockRows, blockCols); } +/** This is the const version of dynBlock(). */ +template +const DynBlock MatrixBase + ::dynBlock(int startRow, int startCol, int blockRows, int blockCols) const +{ + return DynBlock(ref(), startRow, startCol, blockRows, blockCols); +} + + #endif // EIGEN_DYNBLOCK_H diff --git a/Eigen/Core/Fuzzy.h b/Eigen/Core/Fuzzy.h index c1e4ddd74..bd3219568 100644 --- a/Eigen/Core/Fuzzy.h +++ b/Eigen/Core/Fuzzy.h @@ -34,7 +34,7 @@ bool MatrixBase::isApprox( ) const { assert(rows() == other.rows() && cols() == other.cols()); - if(IsVector) + if(IsVectorAtCompileTime) { return((*this - other).norm2() <= std::min(norm2(), other.norm2()) * prec * prec); } @@ -54,7 +54,7 @@ bool MatrixBase::isMuchSmallerThan( const typename NumTraits::Real& prec ) const { - if(IsVector) + if(IsVectorAtCompileTime) { return(norm2() <= abs2(other * prec)); } @@ -75,7 +75,7 @@ bool MatrixBase::isMuchSmallerThan( ) const { assert(rows() == other.rows() && cols() == other.cols()); - if(IsVector) + if(IsVectorAtCompileTime) { return(norm2() <= other.norm2() * prec * prec); } diff --git a/Eigen/Core/Map.h b/Eigen/Core/Map.h index 79951b226..c11ffcfce 100644 --- a/Eigen/Core/Map.h +++ b/Eigen/Core/Map.h @@ -66,23 +66,23 @@ template class Map template const Map MatrixBase::map(const Scalar* data, int rows, int cols) { - return Map(const_cast(data), rows, cols); + return Map(data, rows, cols); } template const Map MatrixBase::map(const Scalar* data, int size) { - assert(IsVector); + assert(IsVectorAtCompileTime); if(ColsAtCompileTime == 1) - return Map(const_cast(data), size, 1); + return Map(data, size, 1); else - return Map(const_cast(data), 1, size); + return Map(data, 1, size); } template const Map MatrixBase::map(const Scalar* data) { - return Map(const_cast(data), RowsAtCompileTime, ColsAtCompileTime); + return Map(data, RowsAtCompileTime, ColsAtCompileTime); } template @@ -94,7 +94,7 @@ Map MatrixBase::map(Scalar* data, int rows, int cols) template Map MatrixBase::map(Scalar* data, int size) { - assert(IsVector); + assert(IsVectorAtCompileTime); if(ColsAtCompileTime == 1) return Map(data, size, 1); else diff --git a/Eigen/Core/Matrix.h b/Eigen/Core/Matrix.h index 989b1a550..73d61a978 100644 --- a/Eigen/Core/Matrix.h +++ b/Eigen/Core/Matrix.h @@ -63,14 +63,23 @@ class Matrix : public MatrixBase<_Scalar, Matrix<_Scalar, _Rows, _Cols> >, template Matrix& operator=(const MatrixBase& other) { - resize(other.rows(), other.cols()); + if(_RowsAtCompileTime == 1) + { + assert(other.isVector()); + resize(1, other.size()); + } + else if(_ColsAtCompileTime == 1) + { + assert(other.isVector()); + resize(other.size(), 1); + } + else resize(other.rows(), other.cols()); return Base::operator=(other); } Matrix& operator=(const Matrix& other) { - resize(other.rows(), other.cols()); - return Base::operator=(other); + return operator=(other); } EIGEN_INHERIT_ASSIGNMENT_OPERATOR(Matrix, +=) diff --git a/Eigen/Core/MatrixBase.h b/Eigen/Core/MatrixBase.h index 85fca57af..4550deabf 100644 --- a/Eigen/Core/MatrixBase.h +++ b/Eigen/Core/MatrixBase.h @@ -77,7 +77,7 @@ template class MatrixBase * columns is known at compile-time to be equal to 1. Indeed, in that case, * we are dealing with a column-vector (if there is only one column) or with * a row-vector (if there is only one row). */ - static const bool IsVector = RowsAtCompileTime == 1 || ColsAtCompileTime == 1; + static const bool IsVectorAtCompileTime = RowsAtCompileTime == 1 || ColsAtCompileTime == 1; /** This is the "reference type" used to pass objects of type MatrixBase as arguments * to functions. If this MatrixBase type represents an expression, then \a Ref @@ -98,9 +98,11 @@ template class MatrixBase /** \returns the number of columns. \sa row(), ColsAtCompileTime*/ int cols() const { return static_cast(this)->_cols(); } /** \returns the number of coefficients, which is \a rows()*cols(). - * \sa rows(), cols(). */ + * \sa rows(), cols(), SizeAtCompileTime. */ int size() const { return rows() * cols(); } - + /** \returns true if either the number of rows or the number of columns is equal to 1. + * \sa rows(), cols(), IsVectorAtCompileTime. */ + bool isVector() const { return rows()==1 || cols()==1; } /** \returns a Ref to *this. \sa Ref */ Ref ref() const { return static_cast(this)->_ref(); } @@ -118,15 +120,27 @@ template class MatrixBase template const Cast cast() const; - Row row(int i) const; - Column col(int i) const; - Minor minor(int row, int col) const; - DynBlock dynBlock(int startRow, int startCol, - int blockRows, int blockCols) const; - template - Block block(int startRow, int startCol) const; + Row row(int i); + const Row row(int i) const; + + Column col(int i); + const Column col(int i) const; + + Minor minor(int row, int col); + const Minor minor(int row, int col) const; + + DynBlock dynBlock(int startRow, int startCol, int blockRows, int blockCols); + const DynBlock + dynBlock(int startRow, int startCol, int blockRows, int blockCols) const; + + template + Block block(int startRow, int startCol); + template + const Block block(int startRow, int startCol) const; + + Transpose transpose(); + const Transpose transpose() const; - Transpose transpose() const; const Conjugate conjugate() const; const Transpose > adjoint() const; Scalar trace() const; @@ -151,7 +165,9 @@ template class MatrixBase template static const DiagonalMatrix diagonal(const OtherDerived& coeffs); - DiagonalCoeffs diagonal() const; + + DiagonalCoeffs diagonal(); + const DiagonalCoeffs diagonal() const; static const Map map(const Scalar* array, int rows, int cols); static const Map map(const Scalar* array, int size); diff --git a/Eigen/Core/MatrixRef.h b/Eigen/Core/MatrixRef.h index ac4fe4bc0..ffdac9d8b 100644 --- a/Eigen/Core/MatrixRef.h +++ b/Eigen/Core/MatrixRef.h @@ -33,7 +33,7 @@ template class MatrixRef typedef typename MatrixType::Scalar Scalar; friend class MatrixBase; - MatrixRef(const MatrixType& matrix) : m_matrix(*const_cast(&matrix)) {} + MatrixRef(const MatrixType& matrix) : m_matrix(matrix) {} MatrixRef(const MatrixRef& other) : m_matrix(other.m_matrix) {} ~MatrixRef() {} @@ -53,11 +53,11 @@ template class MatrixRef Scalar& _coeffRef(int row, int col) { - return m_matrix.coeffRef(row, col); + return const_cast(&m_matrix)->_coeffRef(row, col); } protected: - MatrixType& m_matrix; + const MatrixType& m_matrix; }; #endif // EIGEN_MATRIXREF_H diff --git a/Eigen/Core/Minor.h b/Eigen/Core/Minor.h index 50422c5dd..1d2bf178c 100644 --- a/Eigen/Core/Minor.h +++ b/Eigen/Core/Minor.h @@ -78,9 +78,17 @@ template class Minor * row and column. */ template Minor +MatrixBase::minor(int row, int col) +{ + return Minor(ref(), row, col); +} + +/** This is the const version of minor(). */ +template +const Minor MatrixBase::minor(int row, int col) const { - return Minor(static_cast(const_cast(this))->ref(), row, col); + return Minor(ref(), row, col); } #endif // EIGEN_MINOR_H diff --git a/Eigen/Core/Ones.h b/Eigen/Core/Ones.h index 76a3665c6..36fbe58b0 100644 --- a/Eigen/Core/Ones.h +++ b/Eigen/Core/Ones.h @@ -64,7 +64,7 @@ const Ones MatrixBase::ones(int rows, int cols) template const Ones MatrixBase::ones(int size) { - assert(IsVector); + assert(IsVectorAtCompileTime); if(RowsAtCompileTime == 1) return Ones(1, size); else return Ones(size, 1); } diff --git a/Eigen/Core/OperatorEquals.h b/Eigen/Core/OperatorEquals.h index efda4c033..30d42b670 100644 --- a/Eigen/Core/OperatorEquals.h +++ b/Eigen/Core/OperatorEquals.h @@ -101,7 +101,7 @@ template Derived& MatrixBase ::operator=(const MatrixBase& other) { - if(IsVector && OtherDerived::IsVector) // copying a vector expression into a vector + if(IsVectorAtCompileTime && OtherDerived::IsVectorAtCompileTime) // copying a vector expression into a vector { assert(size() == other.size()); if(EIGEN_UNROLLED_LOOPS && SizeAtCompileTime != Dynamic && SizeAtCompileTime <= 25) diff --git a/Eigen/Core/Random.h b/Eigen/Core/Random.h index 769e737ee..24ef90655 100644 --- a/Eigen/Core/Random.h +++ b/Eigen/Core/Random.h @@ -64,7 +64,7 @@ Eval > MatrixBase::random(int rows, int cols) template Eval > MatrixBase::random(int size) { - assert(IsVector); + assert(IsVectorAtCompileTime); if(RowsAtCompileTime == 1) return Random(1, size).eval(); else return Random(size, 1).eval(); } diff --git a/Eigen/Core/Row.h b/Eigen/Core/Row.h index 84dc6f54c..34c6fe3b0 100644 --- a/Eigen/Core/Row.h +++ b/Eigen/Core/Row.h @@ -103,9 +103,17 @@ template class Row * \sa col(), class Row */ template Row +MatrixBase::row(int i) +{ + return Row(ref(), i); +} + +/** This is the const version of row(). */ +template +const Row MatrixBase::row(int i) const { - return Row(static_cast(const_cast(this))->ref(), i); + return Row(ref(), i); } #endif // EIGEN_ROW_H diff --git a/Eigen/Core/Transpose.h b/Eigen/Core/Transpose.h index 7cf989886..c3723ad82 100644 --- a/Eigen/Core/Transpose.h +++ b/Eigen/Core/Transpose.h @@ -65,9 +65,17 @@ template class Transpose template Transpose +MatrixBase::transpose() +{ + return Transpose(ref()); +} + +/** This is the const version of transpose(). */ +template +const Transpose MatrixBase::transpose() const { - return Transpose(static_cast(const_cast(this))->ref()); + return Transpose(ref()); } #endif // EIGEN_TRANSPOSE_H diff --git a/Eigen/Core/Zero.h b/Eigen/Core/Zero.h index 95685df48..2d3ec264e 100644 --- a/Eigen/Core/Zero.h +++ b/Eigen/Core/Zero.h @@ -64,7 +64,7 @@ const Zero MatrixBase::zero(int rows, int cols) template const Zero MatrixBase::zero(int size) { - assert(IsVector); + assert(IsVectorAtCompileTime); if(RowsAtCompileTime == 1) return Zero(1, size); else return Zero(size, 1); } diff --git a/doc/Doxyfile.in b/doc/Doxyfile.in index 36e19f6d0..721c3de81 100644 --- a/doc/Doxyfile.in +++ b/doc/Doxyfile.in @@ -37,7 +37,7 @@ SEPARATE_MEMBER_PAGES = NO TAB_SIZE = 8 ALIASES = \ "only_for_vectors=This is only for vectors (either row-vectors or column-vectors), \ -as determined by \link MatrixBase::IsVector \endlink." +as determined by \link MatrixBase::IsVectorAtCompileTime \endlink." OPTIMIZE_OUTPUT_FOR_C = NO OPTIMIZE_OUTPUT_JAVA = NO BUILTIN_STL_SUPPORT = NO diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 72a9f8fb7..6fff26b43 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -7,10 +7,13 @@ INCLUDE_DIRECTORIES( ${QT_INCLUDE_DIR} ) SET(test_SRCS main.cpp basicstuff.cpp + linearstructure.cpp + product.cpp adjoint.cpp submatrices.cpp miscmatrices.cpp smallvectors.cpp + map.cpp ) QT4_AUTOMOC(${test_SRCS}) diff --git a/test/basicstuff.cpp b/test/basicstuff.cpp index 04fefc25e..dcc909b4f 100644 --- a/test/basicstuff.cpp +++ b/test/basicstuff.cpp @@ -29,18 +29,9 @@ namespace Eigen { template void basicStuff(const MatrixType& m) { - /* this test covers the following files: - 1) Explicitly (see comments below): - Random.h Zero.h Identity.h Fuzzy.h Sum.h Difference.h - Opposite.h Product.h ScalarMultiple.h Map.h - - 2) Implicitly (the core stuff): - MatrixBase.h Matrix.h MatrixStorage.h CopyHelper.h MatrixRef.h - NumTraits.h Util.h MathFunctions.h OperatorEquals.h Coeffs.h - */ - typedef typename MatrixType::Scalar Scalar; typedef Matrix VectorType; + int rows = m.rows(); int cols = m.cols(); @@ -58,13 +49,9 @@ template void basicStuff(const MatrixType& m) v2 = VectorType::random(rows), vzero = VectorType::zero(rows); - Scalar s1 = random(), - s2 = random(); - int r = random(0, rows-1), c = random(0, cols-1); - // test Fuzzy.h and Zero.h. VERIFY_IS_APPROX( v1, v1); VERIFY_IS_NOT_APPROX( v1, 2*v1); VERIFY_IS_MUCH_SMALLER_THAN( vzero, v1); @@ -86,80 +73,13 @@ template void basicStuff(const MatrixType& m) // operator() that gets called, which in turn calls _read(). VERIFY_IS_MUCH_SMALLER_THAN(MatrixType::zero(rows,cols)(r,c), static_cast(1)); - // test the linear structure, i.e. the following files: - // Sum.h Difference.h Opposite.h ScalarMultiple.h - VERIFY_IS_APPROX(-(-m1), m1); - VERIFY_IS_APPROX(m1+m1, 2*m1); - VERIFY_IS_APPROX(m1+m2-m1, m2); - VERIFY_IS_APPROX(-m2+m1+m2, m1); - VERIFY_IS_APPROX(m1*s1, s1*m1); - VERIFY_IS_APPROX((m1+m2)*s1, s1*m1+s1*m2); - VERIFY_IS_APPROX((s1+s2)*m1, m1*s1+m1*s2); - VERIFY_IS_APPROX((m1-m2)*s1, s1*m1-s1*m2); - VERIFY_IS_APPROX((s1-s2)*m1, m1*s1-m1*s2); - VERIFY_IS_APPROX((-m1+m2)*s1, -s1*m1+s1*m2); - VERIFY_IS_APPROX((-s1+s2)*m1, -m1*s1+m1*s2); - m3 = m2; m3 += m1; - VERIFY_IS_APPROX(m3, m1+m2); - m3 = m2; m3 -= m1; - VERIFY_IS_APPROX(m3, m2-m1); - m3 = m2; m3 *= s1; - VERIFY_IS_APPROX(m3, s1*m2); - if(NumTraits::HasFloatingPoint) - { - m3 = m2; m3 /= s1; - VERIFY_IS_APPROX(m3, m2/s1); - } - - // again, test operator() to check const-qualification - VERIFY_IS_APPROX((-m1)(r,c), -(m1(r,c))); - VERIFY_IS_APPROX((m1-m2)(r,c), (m1(r,c))-(m2(r,c))); - VERIFY_IS_APPROX((m1+m2)(r,c), (m1(r,c))+(m2(r,c))); - VERIFY_IS_APPROX((s1*m1)(r,c), s1*(m1(r,c))); - VERIFY_IS_APPROX((m1*s1)(r,c), (m1(r,c))*s1); - if(NumTraits::HasFloatingPoint) - VERIFY_IS_APPROX((m1/s1)(r,c), (m1(r,c))/s1); - - // begin testing Product.h: only associativity for now - // (we use Transpose.h but this doesn't count as a test for it) - VERIFY_IS_APPROX((m1*m1.transpose())*m2, m1*(m1.transpose()*m2)); - m3 = m1; - m3 *= (m1.transpose() * m2); - VERIFY_IS_APPROX(m3, m1*(m1.transpose()*m2)); - VERIFY_IS_APPROX(m3, m1.lazyProduct(m1.transpose()*m2)); - - // continue testing Product.h: distributivity - VERIFY_IS_APPROX(square*(m1 + m2), square*m1+square*m2); - VERIFY_IS_APPROX(square*(m1 - m2), square*m1-square*m2); - - // continue testing Product.h: compatibility with ScalarMultiple.h - VERIFY_IS_APPROX(s1*(square*m1), (s1*square)*m1); - VERIFY_IS_APPROX(s1*(square*m1), square*(m1*s1)); - - // continue testing Product.h: lazyProduct - VERIFY_IS_APPROX(square.lazyProduct(m1), square*m1); - // again, test operator() to check const-qualification - s1 += square.lazyProduct(m1)(r,c); - - // test Product.h together with Identity.h - VERIFY_IS_APPROX(m1, identity*m1); - VERIFY_IS_APPROX(v1, identity*v1); - // again, test operator() to check const-qualification - VERIFY_IS_APPROX(MatrixType::identity(std::max(rows,cols))(r,c), static_cast(r==c)); - - // test Map.h - Scalar* array1 = new Scalar[rows]; - Scalar* array2 = new Scalar[rows]; - typedef Matrix VectorX; - VectorX::map(array1, rows) = VectorX::random(rows); - VectorX::map(array2, rows) = VectorX::map(array1, rows); - VectorX ma1 = VectorX::map(array1, rows); - VectorX ma2 = VectorX::map(array2, rows); - VERIFY_IS_APPROX(ma1, ma2); - VERIFY_IS_APPROX(ma1, VectorX(array2, rows)); - - delete[] array1; - delete[] array2; + // now test copying a row-vector into a (column-)vector and conversely. + square.col(r) = square.row(r).eval(); + Matrix rv(rows); + Matrix cv(rows); + rv = square.col(r); + cv = square.row(r); + VERIFY_IS_APPROX(rv, cv.transpose()); } void EigenTest::testBasicStuff() diff --git a/test/linearstructure.cpp b/test/linearstructure.cpp new file mode 100644 index 000000000..2bb7a1d81 --- /dev/null +++ b/test/linearstructure.cpp @@ -0,0 +1,106 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. Eigen itself is part of the KDE project. +// +// Copyright (C) 2006-2007 Benoit Jacob +// +// Eigen is free software; 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 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 General Public License for more +// details. +// +// You should have received a copy of the GNU General Public License along +// with Eigen; if not, write to the Free Software Foundation, Inc., 51 +// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. +// +// As a special exception, if other files instantiate templates or use macros +// or functions from this file, or you compile this file and link it +// with other works to produce a work based on this file, this file does not +// by itself cause the resulting work to be covered by the GNU General Public +// License. This exception does not invalidate any other reasons why a work +// based on this file might be covered by the GNU General Public License. + +#include "main.h" + +namespace Eigen { + +template void linearStructure(const MatrixType& m) +{ + /* this test covers the following files: + Sum.h Difference.h Opposite.h ScalarMultiple.h + */ + + typedef typename MatrixType::Scalar Scalar; + typedef Matrix VectorType; + + int rows = m.rows(); + int cols = m.cols(); + + // this test relies a lot on Random.h, and there's not much more that we can do + // to test it, hence I consider that we will have tested Random.h + MatrixType m1 = MatrixType::random(rows, cols), + m2 = MatrixType::random(rows, cols), + m3(rows, cols), + mzero = MatrixType::zero(rows, cols), + identity = Matrix + ::identity(rows), + square = Matrix + ::random(rows, rows); + VectorType v1 = VectorType::random(rows), + v2 = VectorType::random(rows), + vzero = VectorType::zero(rows); + + Scalar s1 = random(), + s2 = random(); + + int r = random(0, rows-1), + c = random(0, cols-1); + + VERIFY_IS_APPROX(-(-m1), m1); + VERIFY_IS_APPROX(m1+m1, 2*m1); + VERIFY_IS_APPROX(m1+m2-m1, m2); + VERIFY_IS_APPROX(-m2+m1+m2, m1); + VERIFY_IS_APPROX(m1*s1, s1*m1); + VERIFY_IS_APPROX((m1+m2)*s1, s1*m1+s1*m2); + VERIFY_IS_APPROX((s1+s2)*m1, m1*s1+m1*s2); + VERIFY_IS_APPROX((m1-m2)*s1, s1*m1-s1*m2); + VERIFY_IS_APPROX((s1-s2)*m1, m1*s1-m1*s2); + VERIFY_IS_APPROX((-m1+m2)*s1, -s1*m1+s1*m2); + VERIFY_IS_APPROX((-s1+s2)*m1, -m1*s1+m1*s2); + m3 = m2; m3 += m1; + VERIFY_IS_APPROX(m3, m1+m2); + m3 = m2; m3 -= m1; + VERIFY_IS_APPROX(m3, m2-m1); + m3 = m2; m3 *= s1; + VERIFY_IS_APPROX(m3, s1*m2); + if(NumTraits::HasFloatingPoint) + { + m3 = m2; m3 /= s1; + VERIFY_IS_APPROX(m3, m2/s1); + } + + // again, test operator() to check const-qualification + VERIFY_IS_APPROX((-m1)(r,c), -(m1(r,c))); + VERIFY_IS_APPROX((m1-m2)(r,c), (m1(r,c))-(m2(r,c))); + VERIFY_IS_APPROX((m1+m2)(r,c), (m1(r,c))+(m2(r,c))); + VERIFY_IS_APPROX((s1*m1)(r,c), s1*(m1(r,c))); + VERIFY_IS_APPROX((m1*s1)(r,c), (m1(r,c))*s1); + if(NumTraits::HasFloatingPoint) + VERIFY_IS_APPROX((m1/s1)(r,c), (m1(r,c))/s1); +} + +void EigenTest::testLinearStructure() +{ + for(int i = 0; i < m_repeat; i++) { + linearStructure(Matrix()); + linearStructure(Matrix4d()); + linearStructure(MatrixXcf(3, 3)); + linearStructure(MatrixXi(8, 12)); + linearStructure(MatrixXcd(20, 20)); + } +} + +} // namespace Eigen diff --git a/test/main.h b/test/main.h index 8cab47331..59c6b2de7 100644 --- a/test/main.h +++ b/test/main.h @@ -114,10 +114,13 @@ class EigenTest : public QObject private slots: void testBasicStuff(); + void testLinearStructure(); + void testProduct(); void testAdjoint(); void testSubmatrices(); void testMiscMatrices(); void testSmallVectors(); + void testMap(); protected: int m_repeat; }; diff --git a/test/map.cpp b/test/map.cpp new file mode 100644 index 000000000..fece15853 --- /dev/null +++ b/test/map.cpp @@ -0,0 +1,60 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. Eigen itself is part of the KDE project. +// +// Copyright (C) 2006-2007 Benoit Jacob +// +// Eigen is free software; 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 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 General Public License for more +// details. +// +// You should have received a copy of the GNU General Public License along +// with Eigen; if not, write to the Free Software Foundation, Inc., 51 +// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. +// +// As a special exception, if other files instantiate templates or use macros +// or functions from this file, or you compile this file and link it +// with other works to produce a work based on this file, this file does not +// by itself cause the resulting work to be covered by the GNU General Public +// License. This exception does not invalidate any other reasons why a work +// based on this file might be covered by the GNU General Public License. + +#include "main.h" + +namespace Eigen { + +template void tmap(const VectorType& m) +{ + typedef typename VectorType::Scalar Scalar; + + int size = m.size(); + + // test Map.h + Scalar* array1 = new Scalar[size]; + Scalar* array2 = new Scalar[size]; + VectorType::map(array1, size) = VectorType::random(size); + VectorType::map(array2, size) = VectorType::map(array1, size); + VectorType ma1 = VectorType::map(array1, size); + VectorType ma2 = VectorType::map(array2, size); + VERIFY_IS_APPROX(ma1, ma2); + VERIFY_IS_APPROX(ma1, VectorType(array2, size)); + delete[] array1; + delete[] array2; +} + +void EigenTest::testMap() +{ + for(int i = 0; i < m_repeat; i++) { + tmap(Matrix()); + tmap(Vector4d()); + tmap(RowVector4f()); + tmap(VectorXcf(8)); + tmap(VectorXi(12)); + } +} + +} // namespace Eigen diff --git a/test/product.cpp b/test/product.cpp new file mode 100644 index 000000000..81610d64c --- /dev/null +++ b/test/product.cpp @@ -0,0 +1,100 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. Eigen itself is part of the KDE project. +// +// Copyright (C) 2006-2007 Benoit Jacob +// +// Eigen is free software; 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 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 General Public License for more +// details. +// +// You should have received a copy of the GNU General Public License along +// with Eigen; if not, write to the Free Software Foundation, Inc., 51 +// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. +// +// As a special exception, if other files instantiate templates or use macros +// or functions from this file, or you compile this file and link it +// with other works to produce a work based on this file, this file does not +// by itself cause the resulting work to be covered by the GNU General Public +// License. This exception does not invalidate any other reasons why a work +// based on this file might be covered by the GNU General Public License. + +#include "main.h" + +namespace Eigen { + +template void product(const MatrixType& m) +{ + /* this test covers the following files: + Identity.h Product.h + */ + + typedef typename MatrixType::Scalar Scalar; + typedef Matrix VectorType; + + int rows = m.rows(); + int cols = m.cols(); + + // this test relies a lot on Random.h, and there's not much more that we can do + // to test it, hence I consider that we will have tested Random.h + MatrixType m1 = MatrixType::random(rows, cols), + m2 = MatrixType::random(rows, cols), + m3(rows, cols), + mzero = MatrixType::zero(rows, cols), + identity = Matrix + ::identity(rows), + square = Matrix + ::random(rows, rows); + VectorType v1 = VectorType::random(rows), + v2 = VectorType::random(rows), + vzero = VectorType::zero(rows); + + Scalar s1 = random(); + + int r = random(0, rows-1), + c = random(0, cols-1); + + // begin testing Product.h: only associativity for now + // (we use Transpose.h but this doesn't count as a test for it) + VERIFY_IS_APPROX((m1*m1.transpose())*m2, m1*(m1.transpose()*m2)); + m3 = m1; + m3 *= (m1.transpose() * m2); + VERIFY_IS_APPROX(m3, m1*(m1.transpose()*m2)); + VERIFY_IS_APPROX(m3, m1.lazyProduct(m1.transpose()*m2)); + + // continue testing Product.h: distributivity + VERIFY_IS_APPROX(square*(m1 + m2), square*m1+square*m2); + VERIFY_IS_APPROX(square*(m1 - m2), square*m1-square*m2); + + // continue testing Product.h: compatibility with ScalarMultiple.h + VERIFY_IS_APPROX(s1*(square*m1), (s1*square)*m1); + VERIFY_IS_APPROX(s1*(square*m1), square*(m1*s1)); + + // continue testing Product.h: lazyProduct + VERIFY_IS_APPROX(square.lazyProduct(m1), square*m1); + // again, test operator() to check const-qualification + s1 += square.lazyProduct(m1)(r,c); + + // test Product.h together with Identity.h + VERIFY_IS_APPROX(m1, identity*m1); + VERIFY_IS_APPROX(v1, identity*v1); + // again, test operator() to check const-qualification + VERIFY_IS_APPROX(MatrixType::identity(std::max(rows,cols))(r,c), static_cast(r==c)); +} + +void EigenTest::testProduct() +{ + for(int i = 0; i < m_repeat; i++) { + product(Matrix()); + product(Matrix4d()); + product(MatrixXcf(3, 3)); + product(MatrixXi(8, 12)); + product(MatrixXcd(20, 20)); + } +} + +} // namespace Eigen