// This file is part of Eigen, a lightweight C++ template library // for linear algebra. // // Copyright (C) 2010-2011 Gael Guennebaud // // This Source Code Form is subject to the terms of the Mozilla // Public License v. 2.0. If a copy of the MPL was not distributed // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. #ifndef EIGEN_TRANSPOSITIONS_H #define EIGEN_TRANSPOSITIONS_H namespace Eigen { /** \class Transpositions * \ingroup Core_Module * * \brief Represents a sequence of transpositions (row/column interchange) * * \param SizeAtCompileTime the number of transpositions, or Dynamic * \param MaxSizeAtCompileTime the maximum number of transpositions, or Dynamic. This optional parameter defaults to SizeAtCompileTime. Most of the time, you should not have to specify it. * * This class represents a permutation transformation as a sequence of \em n transpositions * \f$[T_{n-1} \ldots T_{i} \ldots T_{0}]\f$. It is internally stored as a vector of integers \c indices. * Each transposition \f$T_{i} \f$ applied on the left of a matrix (\f$T_{i} M\f$) interchanges * the rows \c i and \c indices[i] of the matrix \c M. * A transposition applied on the right (e.g., \f$M T_{i}\f$) yields a column interchange. * * Compared to the class PermutationMatrix, such a sequence of transpositions is what is * computed during a decomposition with pivoting, and it is faster when applying the permutation in-place. * * To apply a sequence of transpositions to a matrix, simply use the operator * as in the following example: * \code * Transpositions tr; * MatrixXf mat; * mat = tr * mat; * \endcode * In this example, we detect that the matrix appears on both side, and so the transpositions * are applied in-place without any temporary or extra copy. * * \sa class PermutationMatrix */ namespace internal { template struct transposition_matrix_product_retval; } template class TranspositionsBase { typedef internal::traits Traits; public: typedef typename Traits::IndicesType IndicesType; typedef typename IndicesType::Scalar Index; Derived& derived() { return *static_cast(this); } const Derived& derived() const { return *static_cast(this); } /** Copies the \a other transpositions into \c *this */ template Derived& operator=(const TranspositionsBase& other) { indices() = other.indices(); return derived(); } #ifndef EIGEN_PARSED_BY_DOXYGEN /** This is a special case of the templated operator=. Its purpose is to * prevent a default operator= from hiding the templated operator=. */ Derived& operator=(const TranspositionsBase& other) { indices() = other.indices(); return derived(); } #endif /** \returns the number of transpositions */ inline Index size() const { return indices().size(); } /** Direct access to the underlying index vector */ inline const Index& coeff(Index i) const { return indices().coeff(i); } /** Direct access to the underlying index vector */ inline Index& coeffRef(Index i) { return indices().coeffRef(i); } /** Direct access to the underlying index vector */ inline const Index& operator()(Index i) const { return indices()(i); } /** Direct access to the underlying index vector */ inline Index& operator()(Index i) { return indices()(i); } /** Direct access to the underlying index vector */ inline const Index& operator[](Index i) const { return indices()(i); } /** Direct access to the underlying index vector */ inline Index& operator[](Index i) { return indices()(i); } /** const version of indices(). */ const IndicesType& indices() const { return derived().indices(); } /** \returns a reference to the stored array representing the transpositions. */ IndicesType& indices() { return derived().indices(); } /** Resizes to given size. */ inline void resize(int newSize) { indices().resize(newSize); } /** Sets \c *this to represents an identity transformation */ void setIdentity() { for(int i = 0; i < indices().size(); ++i) coeffRef(i) = i; } // FIXME: do we want such methods ? // might be usefull when the target matrix expression is complex, e.g.: // object.matrix().block(..,..,..,..) = trans * object.matrix().block(..,..,..,..); /* template void applyForwardToRows(MatrixType& mat) const { for(Index k=0 ; k void applyBackwardToRows(MatrixType& mat) const { for(Index k=size()-1 ; k>=0 ; --k) if(m_indices(k)!=k) mat.row(k).swap(mat.row(m_indices(k))); } */ /** \returns the inverse transformation */ inline Transpose inverse() const { return Transpose(derived()); } /** \returns the tranpose transformation */ inline Transpose transpose() const { return Transpose(derived()); } protected: }; namespace internal { template struct traits > { typedef IndexType Index; typedef Matrix IndicesType; }; } template class Transpositions : public TranspositionsBase > { typedef internal::traits Traits; public: typedef TranspositionsBase Base; typedef typename Traits::IndicesType IndicesType; typedef typename IndicesType::Scalar Index; inline Transpositions() {} /** Copy constructor. */ template inline Transpositions(const TranspositionsBase& other) : m_indices(other.indices()) {} #ifndef EIGEN_PARSED_BY_DOXYGEN /** Standard copy constructor. Defined only to prevent a default copy constructor * from hiding the other templated constructor */ inline Transpositions(const Transpositions& other) : m_indices(other.indices()) {} #endif /** Generic constructor from expression of the transposition indices. */ template explicit inline Transpositions(const MatrixBase& a_indices) : m_indices(a_indices) {} /** Copies the \a other transpositions into \c *this */ template Transpositions& operator=(const TranspositionsBase& other) { return Base::operator=(other); } #ifndef EIGEN_PARSED_BY_DOXYGEN /** This is a special case of the templated operator=. Its purpose is to * prevent a default operator= from hiding the templated operator=. */ Transpositions& operator=(const Transpositions& other) { m_indices = other.m_indices; return *this; } #endif /** Constructs an uninitialized permutation matrix of given size. */ inline Transpositions(Index size) : m_indices(size) {} /** const version of indices(). */ const IndicesType& indices() const { return m_indices; } /** \returns a reference to the stored array representing the transpositions. */ IndicesType& indices() { return m_indices; } protected: IndicesType m_indices; }; namespace internal { template struct traits,_PacketAccess> > { typedef IndexType Index; typedef Map, _PacketAccess> IndicesType; }; } template class Map,PacketAccess> : public TranspositionsBase,PacketAccess> > { typedef internal::traits Traits; public: typedef TranspositionsBase Base; typedef typename Traits::IndicesType IndicesType; typedef typename IndicesType::Scalar Index; inline Map(const Index* indicesPtr) : m_indices(indicesPtr) {} inline Map(const Index* indicesPtr, Index size) : m_indices(indicesPtr,size) {} /** Copies the \a other transpositions into \c *this */ template Map& operator=(const TranspositionsBase& other) { return Base::operator=(other); } #ifndef EIGEN_PARSED_BY_DOXYGEN /** This is a special case of the templated operator=. Its purpose is to * prevent a default operator= from hiding the templated operator=. */ Map& operator=(const Map& other) { m_indices = other.m_indices; return *this; } #endif /** const version of indices(). */ const IndicesType& indices() const { return m_indices; } /** \returns a reference to the stored array representing the transpositions. */ IndicesType& indices() { return m_indices; } protected: IndicesType m_indices; }; namespace internal { template struct traits > { typedef typename _IndicesType::Scalar Index; typedef _IndicesType IndicesType; }; } template class TranspositionsWrapper : public TranspositionsBase > { typedef internal::traits Traits; public: typedef TranspositionsBase Base; typedef typename Traits::IndicesType IndicesType; typedef typename IndicesType::Scalar Index; inline TranspositionsWrapper(IndicesType& a_indices) : m_indices(a_indices) {} /** Copies the \a other transpositions into \c *this */ template TranspositionsWrapper& operator=(const TranspositionsBase& other) { return Base::operator=(other); } #ifndef EIGEN_PARSED_BY_DOXYGEN /** This is a special case of the templated operator=. Its purpose is to * prevent a default operator= from hiding the templated operator=. */ TranspositionsWrapper& operator=(const TranspositionsWrapper& other) { m_indices = other.m_indices; return *this; } #endif /** const version of indices(). */ const IndicesType& indices() const { return m_indices; } /** \returns a reference to the stored array representing the transpositions. */ IndicesType& indices() { return m_indices; } protected: const typename IndicesType::Nested m_indices; }; /** \returns the \a matrix with the \a transpositions applied to the columns. */ template inline const internal::transposition_matrix_product_retval operator*(const MatrixBase& matrix, const TranspositionsBase &transpositions) { return internal::transposition_matrix_product_retval (transpositions.derived(), matrix.derived()); } /** \returns the \a matrix with the \a transpositions applied to the rows. */ template inline const internal::transposition_matrix_product_retval operator*(const TranspositionsBase &transpositions, const MatrixBase& matrix) { return internal::transposition_matrix_product_retval (transpositions.derived(), matrix.derived()); } namespace internal { template struct traits > { typedef typename MatrixType::PlainObject ReturnType; }; template struct transposition_matrix_product_retval : public ReturnByValue > { typedef typename remove_all::type MatrixTypeNestedCleaned; typedef typename TranspositionType::Index Index; transposition_matrix_product_retval(const TranspositionType& tr, const MatrixType& matrix) : m_transpositions(tr), m_matrix(matrix) {} inline int rows() const { return m_matrix.rows(); } inline int cols() const { return m_matrix.cols(); } template inline void evalTo(Dest& dst) const { const int size = m_transpositions.size(); Index j = 0; if(!(is_same::value && extract_data(dst) == extract_data(m_matrix))) dst = m_matrix; for(int k=(Transposed?size-1:0) ; Transposed?k>=0:k class Transpose > { typedef TranspositionsDerived TranspositionType; typedef typename TranspositionType::IndicesType IndicesType; public: Transpose(const TranspositionType& t) : m_transpositions(t) {} inline int size() const { return m_transpositions.size(); } /** \returns the \a matrix with the inverse transpositions applied to the columns. */ template friend inline const internal::transposition_matrix_product_retval operator*(const MatrixBase& matrix, const Transpose& trt) { return internal::transposition_matrix_product_retval(trt.m_transpositions, matrix.derived()); } /** \returns the \a matrix with the inverse transpositions applied to the rows. */ template inline const internal::transposition_matrix_product_retval operator*(const MatrixBase& matrix) const { return internal::transposition_matrix_product_retval(m_transpositions, matrix.derived()); } protected: const TranspositionType& m_transpositions; }; } // end namespace Eigen #endif // EIGEN_TRANSPOSITIONS_H