diff options
Diffstat (limited to 'Eigen/src/Core/PermutationMatrix.h')
| -rw-r--r-- | Eigen/src/Core/PermutationMatrix.h | 133 |
1 files changed, 89 insertions, 44 deletions
diff --git a/Eigen/src/Core/PermutationMatrix.h b/Eigen/src/Core/PermutationMatrix.h index 8aa4c8bc5..4846f2ae1 100644 --- a/Eigen/src/Core/PermutationMatrix.h +++ b/Eigen/src/Core/PermutationMatrix.h @@ -13,7 +13,8 @@ namespace Eigen { -template<int RowCol,typename IndicesType,typename MatrixType, typename StorageKind> class PermutedImpl; +// TODO: this does not seems to be needed at all: +// template<int RowCol,typename IndicesType,typename MatrixType, typename StorageKind> class PermutedImpl; /** \class PermutationBase * \ingroup Core_Module @@ -60,7 +61,6 @@ class PermutationBase : public EigenBase<Derived> typedef typename Traits::IndicesType IndicesType; enum { Flags = Traits::Flags, - CoeffReadCost = Traits::CoeffReadCost, RowsAtCompileTime = Traits::RowsAtCompileTime, ColsAtCompileTime = Traits::ColsAtCompileTime, MaxRowsAtCompileTime = Traits::MaxRowsAtCompileTime, @@ -73,6 +73,7 @@ class PermutationBase : public EigenBase<Derived> typedef PermutationMatrix<IndicesType::SizeAtCompileTime,IndicesType::MaxSizeAtCompileTime,StorageIndexType> PlainPermutationType; using Base::derived; + typedef Transpose<PermutationBase> TransposeReturnType; #endif /** Copies the other permutation into *this */ @@ -118,7 +119,7 @@ class PermutationBase : public EigenBase<Derived> void evalTo(MatrixBase<DenseDerived>& other) const { other.setZero(); - for (int i=0; i<rows();++i) + for (Index i=0; i<rows(); ++i) other.coeffRef(indices().coeff(i),i) = typename DenseDerived::Scalar(1); } #endif @@ -163,10 +164,10 @@ class PermutationBase : public EigenBase<Derived> * * \returns a reference to *this. * - * \warning This is much slower than applyTranspositionOnTheRight(int,int): + * \warning This is much slower than applyTranspositionOnTheRight(Index,Index): * this has linear complexity and requires a lot of branching. * - * \sa applyTranspositionOnTheRight(int,int) + * \sa applyTranspositionOnTheRight(Index,Index) */ Derived& applyTranspositionOnTheLeft(Index i, Index j) { @@ -185,7 +186,7 @@ class PermutationBase : public EigenBase<Derived> * * This is a fast operation, it only consists in swapping two indices. * - * \sa applyTranspositionOnTheLeft(int,int) + * \sa applyTranspositionOnTheLeft(Index,Index) */ Derived& applyTranspositionOnTheRight(Index i, Index j) { @@ -198,14 +199,14 @@ class PermutationBase : public EigenBase<Derived> * * \note \note_try_to_help_rvo */ - inline Transpose<PermutationBase> inverse() const - { return derived(); } + inline TransposeReturnType inverse() const + { return TransposeReturnType(derived()); } /** \returns the tranpose permutation matrix. * * \note \note_try_to_help_rvo */ - inline Transpose<PermutationBase> transpose() const - { return derived(); } + inline TransposeReturnType transpose() const + { return TransposeReturnType(derived()); } /**** multiplication helpers to hopefully get RVO ****/ @@ -215,13 +216,13 @@ class PermutationBase : public EigenBase<Derived> template<typename OtherDerived> void assignTranspose(const PermutationBase<OtherDerived>& other) { - for (int i=0; i<rows();++i) indices().coeffRef(other.indices().coeff(i)) = i; + for (Index i=0; i<rows();++i) indices().coeffRef(other.indices().coeff(i)) = i; } template<typename Lhs,typename Rhs> void assignProduct(const Lhs& lhs, const Rhs& rhs) { eigen_assert(lhs.cols() == rhs.rows()); - for (int i=0; i<rows();++i) indices().coeffRef(i) = lhs.indices().coeff(rhs.indices().coeff(i)); + for (Index i=0; i<rows();++i) indices().coeffRef(i) = lhs.indices().coeff(rhs.indices().coeff(i)); } #endif @@ -274,6 +275,7 @@ template<int SizeAtCompileTime, int MaxSizeAtCompileTime, typename _StorageIndex struct traits<PermutationMatrix<SizeAtCompileTime, MaxSizeAtCompileTime, _StorageIndexType> > : traits<Matrix<_StorageIndexType,SizeAtCompileTime,SizeAtCompileTime,0,MaxSizeAtCompileTime,MaxSizeAtCompileTime> > { + typedef PermutationStorage StorageKind; typedef Matrix<_StorageIndexType, SizeAtCompileTime, 1, 0, MaxSizeAtCompileTime, 1> IndicesType; typedef typename IndicesType::Index Index; typedef _StorageIndexType StorageIndexType; @@ -287,6 +289,8 @@ class PermutationMatrix : public PermutationBase<PermutationMatrix<SizeAtCompile typedef internal::traits<PermutationMatrix> Traits; public: + typedef const PermutationMatrix& Nested; + #ifndef EIGEN_PARSED_BY_DOXYGEN typedef typename Traits::IndicesType IndicesType; typedef typename Traits::StorageIndexType StorageIndexType; @@ -298,8 +302,10 @@ class PermutationMatrix : public PermutationBase<PermutationMatrix<SizeAtCompile /** Constructs an uninitialized permutation matrix of given size. */ - inline PermutationMatrix(Index size) : m_indices(size) - {} + explicit inline PermutationMatrix(Index size) : m_indices(size) + { + eigen_internal_assert(size <= NumTraits<StorageIndexType>::highest()); + } /** Copy constructor. */ template<typename OtherDerived> @@ -370,7 +376,10 @@ class PermutationMatrix : public PermutationBase<PermutationMatrix<SizeAtCompile PermutationMatrix(const Transpose<PermutationBase<Other> >& other) : m_indices(other.nestedPermutation().size()) { - for (int i=0; i<m_indices.size();++i) m_indices.coeffRef(other.nestedPermutation().indices().coeff(i)) = i; + eigen_internal_assert(m_indices.size() <= NumTraits<StorageIndexType>::highest()); + StorageIndexType end = StorageIndexType(m_indices.size()); + for (StorageIndexType i=0; i<end;++i) + m_indices.coeffRef(other.nestedPermutation().indices().coeff(i)) = i; } template<typename Lhs,typename Rhs> PermutationMatrix(internal::PermPermProduct_t, const Lhs& lhs, const Rhs& rhs) @@ -391,6 +400,7 @@ template<int SizeAtCompileTime, int MaxSizeAtCompileTime, typename _StorageIndex struct traits<Map<PermutationMatrix<SizeAtCompileTime, MaxSizeAtCompileTime, _StorageIndexType>,_PacketAccess> > : traits<Matrix<_StorageIndexType,SizeAtCompileTime,SizeAtCompileTime,0,MaxSizeAtCompileTime,MaxSizeAtCompileTime> > { + typedef PermutationStorage StorageKind; typedef Map<const Matrix<_StorageIndexType, SizeAtCompileTime, 1, 0, MaxSizeAtCompileTime, 1>, _PacketAccess> IndicesType; typedef typename IndicesType::Index Index; typedef _StorageIndexType StorageIndexType; @@ -462,8 +472,6 @@ class Map<PermutationMatrix<SizeAtCompileTime, MaxSizeAtCompileTime, _StorageInd * \sa class PermutationBase, class PermutationMatrix */ -struct PermutationStorage {}; - template<typename _IndicesType> class TranspositionsWrapper; namespace internal { template<typename _IndicesType> @@ -477,10 +485,9 @@ struct traits<PermutationWrapper<_IndicesType> > enum { RowsAtCompileTime = _IndicesType::SizeAtCompileTime, ColsAtCompileTime = _IndicesType::SizeAtCompileTime, - MaxRowsAtCompileTime = IndicesType::MaxRowsAtCompileTime, - MaxColsAtCompileTime = IndicesType::MaxColsAtCompileTime, - Flags = 0, - CoeffReadCost = _IndicesType::CoeffReadCost + MaxRowsAtCompileTime = IndicesType::MaxSizeAtCompileTime, + MaxColsAtCompileTime = IndicesType::MaxSizeAtCompileTime, + Flags = 0 }; }; } @@ -509,35 +516,39 @@ class PermutationWrapper : public PermutationBase<PermutationWrapper<_IndicesTyp typename IndicesType::Nested m_indices; }; + +// TODO: Do we need to define these operator* functions? Would it be better to have them inherited +// from MatrixBase? + /** \returns the matrix with the permutation applied to the columns. */ -template<typename Derived, typename PermutationDerived> -inline const internal::permut_matrix_product_retval<PermutationDerived, Derived, OnTheRight> -operator*(const MatrixBase<Derived>& matrix, - const PermutationBase<PermutationDerived> &permutation) +template<typename MatrixDerived, typename PermutationDerived> +EIGEN_DEVICE_FUNC +const Product<MatrixDerived, PermutationDerived, DefaultProduct> +operator*(const MatrixBase<MatrixDerived> &matrix, + const PermutationBase<PermutationDerived>& permutation) { - return internal::permut_matrix_product_retval - <PermutationDerived, Derived, OnTheRight> - (permutation.derived(), matrix.derived()); + return Product<MatrixDerived, PermutationDerived, DefaultProduct> + (matrix.derived(), permutation.derived()); } /** \returns the matrix with the permutation applied to the rows. */ -template<typename Derived, typename PermutationDerived> -inline const internal::permut_matrix_product_retval - <PermutationDerived, Derived, OnTheLeft> +template<typename PermutationDerived, typename MatrixDerived> +EIGEN_DEVICE_FUNC +const Product<PermutationDerived, MatrixDerived, DefaultProduct> operator*(const PermutationBase<PermutationDerived> &permutation, - const MatrixBase<Derived>& matrix) + const MatrixBase<MatrixDerived>& matrix) { - return internal::permut_matrix_product_retval - <PermutationDerived, Derived, OnTheLeft> - (permutation.derived(), matrix.derived()); + return Product<PermutationDerived, MatrixDerived, DefaultProduct> + (permutation.derived(), matrix.derived()); } namespace internal { template<typename PermutationType, typename MatrixType, int Side, bool Transposed> struct traits<permut_matrix_product_retval<PermutationType, MatrixType, Side, Transposed> > + : traits<typename MatrixType::PlainObject> { typedef typename MatrixType::PlainObject ReturnType; }; @@ -590,7 +601,7 @@ struct permut_matrix_product_retval } else { - for(int i = 0; i < n; ++i) + for(Index i = 0; i < n; ++i) { Block<Dest, Side==OnTheLeft ? 1 : Dest::RowsAtCompileTime, Side==OnTheRight ? 1 : Dest::ColsAtCompileTime> (dst, ((Side==OnTheLeft) ^ Transposed) ? m_permutation.indices().coeff(i) : i) @@ -617,6 +628,8 @@ struct traits<Transpose<PermutationBase<Derived> > > } // end namespace internal +// TODO: the specificties should be handled by the evaluator, +// at the very least we should only specialize TransposeImpl template<typename Derived> class Transpose<PermutationBase<Derived> > : public EigenBase<Transpose<PermutationBase<Derived> > > @@ -631,26 +644,26 @@ class Transpose<PermutationBase<Derived> > typedef typename Derived::DenseMatrixType DenseMatrixType; enum { Flags = Traits::Flags, - CoeffReadCost = Traits::CoeffReadCost, RowsAtCompileTime = Traits::RowsAtCompileTime, ColsAtCompileTime = Traits::ColsAtCompileTime, MaxRowsAtCompileTime = Traits::MaxRowsAtCompileTime, MaxColsAtCompileTime = Traits::MaxColsAtCompileTime }; typedef typename Traits::Scalar Scalar; + typedef typename Traits::Index Index; #endif Transpose(const PermutationType& p) : m_permutation(p) {} - inline int rows() const { return m_permutation.rows(); } - inline int cols() const { return m_permutation.cols(); } + inline Index rows() const { return m_permutation.rows(); } + inline Index cols() const { return m_permutation.cols(); } #ifndef EIGEN_PARSED_BY_DOXYGEN template<typename DenseDerived> void evalTo(MatrixBase<DenseDerived>& other) const { other.setZero(); - for (int i=0; i<rows();++i) + for (Index i=0; i<rows();++i) other.coeffRef(i, m_permutation.indices().coeff(i)) = typename DenseDerived::Scalar(1); } #endif @@ -663,19 +676,19 @@ class Transpose<PermutationBase<Derived> > /** \returns the matrix with the inverse permutation applied to the columns. */ template<typename OtherDerived> friend - inline const internal::permut_matrix_product_retval<PermutationType, OtherDerived, OnTheRight, true> + const Product<OtherDerived, Transpose, DefaultProduct> operator*(const MatrixBase<OtherDerived>& matrix, const Transpose& trPerm) { - return internal::permut_matrix_product_retval<PermutationType, OtherDerived, OnTheRight, true>(trPerm.m_permutation, matrix.derived()); + return Product<OtherDerived, Transpose, DefaultProduct>(matrix.derived(), trPerm.derived()); } /** \returns the matrix with the inverse permutation applied to the rows. */ template<typename OtherDerived> - inline const internal::permut_matrix_product_retval<PermutationType, OtherDerived, OnTheLeft, true> + const Product<Transpose, OtherDerived, DefaultProduct> operator*(const MatrixBase<OtherDerived>& matrix) const { - return internal::permut_matrix_product_retval<PermutationType, OtherDerived, OnTheLeft, true>(m_permutation, matrix.derived()); + return Product<Transpose, OtherDerived, DefaultProduct>(*this, matrix.derived()); } const PermutationType& nestedPermutation() const { return m_permutation; } @@ -690,6 +703,38 @@ const PermutationWrapper<const Derived> MatrixBase<Derived>::asPermutation() con return derived(); } +namespace internal { + +// TODO currently a permutation matrix expression has the form PermutationMatrix or PermutationWrapper +// or their transpose; in the future shape should be defined by the expression traits +template<int SizeAtCompileTime, int MaxSizeAtCompileTime, typename IndexType> +struct evaluator_traits<PermutationMatrix<SizeAtCompileTime, MaxSizeAtCompileTime, IndexType> > +{ + typedef typename storage_kind_to_evaluator_kind<Dense>::Kind Kind; + typedef PermutationShape Shape; + static const int AssumeAliasing = 0; +}; + +template<typename IndicesType> +struct evaluator_traits<PermutationWrapper<IndicesType> > +{ + typedef typename storage_kind_to_evaluator_kind<Dense>::Kind Kind; + typedef PermutationShape Shape; + static const int AssumeAliasing = 0; +}; + +template<typename Derived> +struct evaluator_traits<Transpose<PermutationBase<Derived> > > +{ + typedef typename storage_kind_to_evaluator_kind<Dense>::Kind Kind; + typedef PermutationShape Shape; + static const int AssumeAliasing = 0; +}; + +template<> struct AssignmentKind<DenseShape,PermutationShape> { typedef EigenBase2EigenBase Kind; }; + +} // end namespace internal + } // end namespace Eigen #endif // EIGEN_PERMUTATIONMATRIX_H |