diff options
author | Benoit Jacob <jacob.benoit.1@gmail.com> | 2009-11-16 15:07:33 -0500 |
---|---|---|
committer | Benoit Jacob <jacob.benoit.1@gmail.com> | 2009-11-16 15:07:33 -0500 |
commit | e8d0dbf82e224fa57f0741bef58c3f57ec9f7895 (patch) | |
tree | 9000372e2e14205dc433aabc52068e30059fd2e7 /Eigen | |
parent | 8a1bada43d64820f5e973d1692cc51eaf865bf76 (diff) |
PermutationMatrix:
* make multiplication order not be reversed
* release-quality documentation
Diffstat (limited to 'Eigen')
-rw-r--r-- | Eigen/src/Core/PermutationMatrix.h | 72 |
1 files changed, 57 insertions, 15 deletions
diff --git a/Eigen/src/Core/PermutationMatrix.h b/Eigen/src/Core/PermutationMatrix.h index 1c66cde8e..f2bde0e71 100644 --- a/Eigen/src/Core/PermutationMatrix.h +++ b/Eigen/src/Core/PermutationMatrix.h @@ -25,18 +25,24 @@ #ifndef EIGEN_PERMUTATIONMATRIX_H #define EIGEN_PERMUTATIONMATRIX_H -/** \nonstableyet - * \class PermutationMatrix +/** \class PermutationMatrix * * \brief Permutation matrix * * \param SizeAtCompileTime the number of rows/cols, or Dynamic - * \param MaxSizeAtCompileTime the maximum number of rows/cols, or Dynamic. This optional parameter defaults to SizeAtCompileTime. + * \param MaxSizeAtCompileTime the maximum number of rows/cols, or Dynamic. This optional parameter defaults to SizeAtCompileTime. Most of the time, you should not have to specify it. * * This class represents a permutation matrix, internally stored as a vector of integers. - * The convention followed here is the same as on <a href="http://en.wikipedia.org/wiki/Permutation_matrix">Wikipedia</a>, - * namely: the matrix of permutation \a p is the matrix such that on each row \a i, the only nonzero coefficient is - * in column p(i). + * The convention followed here is that if \f$ \sigma \f$ is a permutation, the corresponding permutation matrix + * \f$ P_\sigma \f$ is such that if \f$ (e_1,\ldots,e_p) \f$ is the canonical basis, we have: + * \f[ P_\sigma(e_i) = e_{\sigma(i)}. \f] + * This convention ensures that for any two permutations \f$ \sigma, \tau \f$, we have: + * \f[ P_{\sigma\circ\tau} = P_\sigma P_\tau. \f] + * + * Permutation matrices are square and invertible. + * + * Notice that in addition to the member functions and operators listed here, there also are non-member + * operator* to multiply a PermutationMatrix with any kind of matrix expression (MatrixBase) on either side. * * \sa class DiagonalMatrix */ @@ -53,6 +59,7 @@ class PermutationMatrix : public AnyMatrixBase<PermutationMatrix<SizeAtCompileTi { public: + #ifndef EIGEN_PARSED_BY_DOXYGEN typedef ei_traits<PermutationMatrix> Traits; typedef Matrix<int,SizeAtCompileTime,SizeAtCompileTime,0,MaxSizeAtCompileTime,MaxSizeAtCompileTime> DenseMatrixType; @@ -65,25 +72,37 @@ class PermutationMatrix : public AnyMatrixBase<PermutationMatrix<SizeAtCompileTi MaxColsAtCompileTime = Traits::MaxColsAtCompileTime }; typedef typename Traits::Scalar Scalar; + #endif - typedef Matrix<int, RowsAtCompileTime, 1, 0, MaxRowsAtCompileTime, 1> IndicesType; + typedef Matrix<int, SizeAtCompileTime, 1, 0, MaxSizeAtCompileTime, 1> IndicesType; inline PermutationMatrix() { } + /** Copy constructor. */ template<int OtherSize, int OtherMaxSize> inline PermutationMatrix(const PermutationMatrix<OtherSize, OtherMaxSize>& other) : m_indices(other.indices()) {} - /** copy constructor. prevent a default copy constructor from hiding the other templated constructor */ + #ifndef EIGEN_PARSED_BY_DOXYGEN + /** Standard copy constructor. Defined only to prevent a default copy constructor + * from hiding the other templated constructor */ inline PermutationMatrix(const PermutationMatrix& other) : m_indices(other.indices()) {} - - /** generic constructor from expression of the indices */ + #endif + + /** Generic constructor from expression of the indices. The indices + * array has the meaning that the permutations sends each integer i to indices[i]. + * + * \warning It is your responsibility to check that the indices array that you passes actually + * describes a permutation, \ie each value between 0 and n-1 occurs exactly once, where n is the + * array's size. + */ template<typename Other> - explicit inline PermutationMatrix(const MatrixBase<Other>& other) : m_indices(other) + explicit inline PermutationMatrix(const MatrixBase<Other>& indices) : m_indices(indices) {} + /** Copies the other permutation into *this */ template<int OtherSize, int OtherMaxSize> PermutationMatrix& operator=(const PermutationMatrix<OtherSize, OtherMaxSize>& other) { @@ -91,6 +110,7 @@ class PermutationMatrix : public AnyMatrixBase<PermutationMatrix<SizeAtCompileTi return *this; } + #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=. */ @@ -99,8 +119,12 @@ class PermutationMatrix : public AnyMatrixBase<PermutationMatrix<SizeAtCompileTi m_indices = other.m_indices(); return *this; } + #endif - inline PermutationMatrix(int rows, int cols) : m_indices(rows) + /** Constructs an uninitialized permutation matrix of given size. Note that it is required + * that rows == cols, since permutation matrices are square. The \a cols parameter may be omitted. + */ + inline PermutationMatrix(int rows, int cols = rows) : m_indices(rows) { ei_assert(rows == cols); } @@ -111,24 +135,33 @@ class PermutationMatrix : public AnyMatrixBase<PermutationMatrix<SizeAtCompileTi /** \returns the number of columns */ inline int cols() const { return m_indices.size(); } + #ifndef EIGEN_PARSED_BY_DOXYGEN template<typename DenseDerived> void evalTo(MatrixBase<DenseDerived>& other) const { other.setZero(); for (int i=0; i<rows();++i) - other.coeffRef(i,m_indices.coeff(i)) = typename DenseDerived::Scalar(1); + other.coeffRef(m_indices.coeff(i),i) = typename DenseDerived::Scalar(1); } + #endif + /** \returns a Matrix object initialized from this permutation matrix. Notice that it + * is inefficient to return this Matrix object by value. For efficiency, favor using + * the Matrix constructor taking AnyMatrixBase objects. + */ DenseMatrixType toDenseMatrix() const { return *this; } + /** const version of indices(). */ const IndicesType& indices() const { return m_indices; } + /** \returns a reference to the stored array representing the permutation. */ IndicesType& indices() { return m_indices; } /**** inversion and multiplication helpers to hopefully get RVO ****/ +#ifndef EIGEN_PARSED_BY_DOXYGEN protected: enum Inverse_t {Inverse}; PermutationMatrix(Inverse_t, const PermutationMatrix& other) @@ -143,10 +176,19 @@ class PermutationMatrix : public AnyMatrixBase<PermutationMatrix<SizeAtCompileTi ei_assert(lhs.cols() == rhs.rows()); for (int i=0; i<rows();++i) m_indices.coeffRef(i) = lhs.m_indices.coeff(rhs.m_indices.coeff(i)); } +#endif public: + /** \returns the inverse permutation matrix. + * + * \note \note_try_to_help_rvo + */ inline PermutationMatrix inverse() const { return PermutationMatrix(Inverse, *this); } + /** \returns the product permutation matrix. + * + * \note \note_try_to_help_rvo + */ template<int OtherSize, int OtherMaxSize> inline PermutationMatrix operator*(const PermutationMatrix<OtherSize, OtherMaxSize>& other) const { return PermutationMatrix(Product, *this, other); } @@ -209,7 +251,7 @@ struct ei_permut_matrix_product_retval Dest, Side==OnTheLeft ? 1 : Dest::RowsAtCompileTime, Side==OnTheRight ? 1 : Dest::ColsAtCompileTime - >(dst, Side==OnTheRight ? m_permutation.indices().coeff(i) : i) + >(dst, Side==OnTheLeft ? m_permutation.indices().coeff(i) : i) = @@ -217,7 +259,7 @@ struct ei_permut_matrix_product_retval MatrixTypeNestedCleaned, Side==OnTheLeft ? 1 : MatrixType::RowsAtCompileTime, Side==OnTheRight ? 1 : MatrixType::ColsAtCompileTime - >(m_matrix, Side==OnTheLeft ? m_permutation.indices().coeff(i) : i); + >(m_matrix, Side==OnTheRight ? m_permutation.indices().coeff(i) : i); } } |