mirror of
https://gitlab.com/libeigen/eigen.git
synced 2025-10-20 20:11:07 +08:00
301 lines
11 KiB
C++
301 lines
11 KiB
C++
// This file is part of Eigen, a lightweight C++ template library
|
|
// for linear algebra.
|
|
//
|
|
// Copyright (C) 2010 Gael Guennebaud <gael.guennebaud@inria.fr>
|
|
//
|
|
// 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/>.
|
|
|
|
#ifndef EIGEN_TRANSPOSITIONS_H
|
|
#define EIGEN_TRANSPOSITIONS_H
|
|
|
|
/** \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<typename TranspositionType, typename MatrixType, int Side, bool Transposed=false> struct transposition_matrix_product_retval;
|
|
}
|
|
|
|
template<int SizeAtCompileTime, int MaxSizeAtCompileTime>
|
|
class Transpositions
|
|
{
|
|
public:
|
|
|
|
typedef Matrix<DenseIndex, SizeAtCompileTime, 1, 0, MaxSizeAtCompileTime, 1> IndicesType;
|
|
typedef typename IndicesType::Index Index;
|
|
|
|
inline Transpositions() {}
|
|
|
|
/** Copy constructor. */
|
|
template<int OtherSize, int OtherMaxSize>
|
|
inline Transpositions(const Transpositions<OtherSize, OtherMaxSize>& 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<typename Other>
|
|
explicit inline Transpositions(const MatrixBase<Other>& indices) : m_indices(indices)
|
|
{}
|
|
|
|
/** Copies the \a other transpositions into \c *this */
|
|
template<int OtherSize, int OtherMaxSize>
|
|
Transpositions& operator=(const Transpositions<OtherSize, OtherMaxSize>& other)
|
|
{
|
|
m_indices = other.indices();
|
|
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=.
|
|
*/
|
|
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)
|
|
{}
|
|
|
|
/** \returns the number of transpositions */
|
|
inline Index size() const { return m_indices.size(); }
|
|
|
|
/** Direct access to the underlying index vector */
|
|
inline const Index& coeff(Index i) const { return m_indices.coeff(i); }
|
|
/** Direct access to the underlying index vector */
|
|
inline Index& coeffRef(Index i) { return m_indices.coeffRef(i); }
|
|
/** Direct access to the underlying index vector */
|
|
inline const Index& operator()(Index i) const { return m_indices(i); }
|
|
/** Direct access to the underlying index vector */
|
|
inline Index& operator()(Index i) { return m_indices(i); }
|
|
/** Direct access to the underlying index vector */
|
|
inline const Index& operator[](Index i) const { return m_indices(i); }
|
|
/** Direct access to the underlying index vector */
|
|
inline Index& operator[](Index i) { return m_indices(i); }
|
|
|
|
/** 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; }
|
|
|
|
/** Resizes to given size. */
|
|
inline void resize(int size)
|
|
{
|
|
m_indices.resize(size);
|
|
}
|
|
|
|
/** Sets \c *this to represents an identity transformation */
|
|
void setIdentity()
|
|
{
|
|
for(int i = 0; i < m_indices.size(); ++i)
|
|
m_indices.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<typename MatrixType>
|
|
void applyForwardToRows(MatrixType& mat) const
|
|
{
|
|
for(Index k=0 ; k<size() ; ++k)
|
|
if(m_indices(k)!=k)
|
|
mat.row(k).swap(mat.row(m_indices(k)));
|
|
}
|
|
|
|
template<typename MatrixType>
|
|
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<Transpositions> inverse() const
|
|
{ return *this; }
|
|
|
|
/** \returns the tranpose transformation */
|
|
inline Transpose<Transpositions> transpose() const
|
|
{ return *this; }
|
|
|
|
#ifndef EIGEN_PARSED_BY_DOXYGEN
|
|
template<int OtherSize, int OtherMaxSize>
|
|
Transpositions(const Transpose<Transpositions<OtherSize,OtherMaxSize> >& other)
|
|
: m_indices(other.size())
|
|
{
|
|
Index n = size();
|
|
Index j = size-1;
|
|
for(Index i=0; i<n;++i,--j)
|
|
m_indices.coeffRef(j) = other.nestedTranspositions().indices().coeff(i);
|
|
}
|
|
#endif
|
|
|
|
protected:
|
|
|
|
IndicesType m_indices;
|
|
};
|
|
|
|
/** \returns the \a matrix with the \a transpositions applied to the columns.
|
|
*/
|
|
template<typename Derived, int SizeAtCompileTime, int MaxSizeAtCompileTime>
|
|
inline const internal::transposition_matrix_product_retval<Transpositions<SizeAtCompileTime, MaxSizeAtCompileTime>, Derived, OnTheRight>
|
|
operator*(const MatrixBase<Derived>& matrix,
|
|
const Transpositions<SizeAtCompileTime, MaxSizeAtCompileTime> &transpositions)
|
|
{
|
|
return internal::transposition_matrix_product_retval
|
|
<Transpositions<SizeAtCompileTime, MaxSizeAtCompileTime>, Derived, OnTheRight>
|
|
(transpositions, matrix.derived());
|
|
}
|
|
|
|
/** \returns the \a matrix with the \a transpositions applied to the rows.
|
|
*/
|
|
template<typename Derived, int SizeAtCompileTime, int MaxSizeAtCompileTime>
|
|
inline const internal::transposition_matrix_product_retval
|
|
<Transpositions<SizeAtCompileTime, MaxSizeAtCompileTime>, Derived, OnTheLeft>
|
|
operator*(const Transpositions<SizeAtCompileTime, MaxSizeAtCompileTime> &transpositions,
|
|
const MatrixBase<Derived>& matrix)
|
|
{
|
|
return internal::transposition_matrix_product_retval
|
|
<Transpositions<SizeAtCompileTime, MaxSizeAtCompileTime>, Derived, OnTheLeft>
|
|
(transpositions, matrix.derived());
|
|
}
|
|
|
|
namespace internal {
|
|
|
|
template<typename TranspositionType, typename MatrixType, int Side, bool Transposed>
|
|
struct traits<transposition_matrix_product_retval<TranspositionType, MatrixType, Side, Transposed> >
|
|
{
|
|
typedef typename MatrixType::PlainObject ReturnType;
|
|
};
|
|
|
|
template<typename TranspositionType, typename MatrixType, int Side, bool Transposed>
|
|
struct transposition_matrix_product_retval
|
|
: public ReturnByValue<transposition_matrix_product_retval<TranspositionType, MatrixType, Side, Transposed> >
|
|
{
|
|
typedef typename remove_all<typename MatrixType::Nested>::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<typename Dest> inline void evalTo(Dest& dst) const
|
|
{
|
|
const int size = m_transpositions.size();
|
|
Index j = 0;
|
|
|
|
if(!(is_same<MatrixTypeNestedCleaned,Dest>::value && extract_data(dst) == extract_data(m_matrix)))
|
|
dst = m_matrix;
|
|
|
|
for(int k=(Transposed?size-1:0) ; Transposed?k>=0:k<size ; Transposed?--k:++k)
|
|
if((j=m_transpositions.coeff(k))!=k)
|
|
{
|
|
if(Side==OnTheLeft)
|
|
dst.row(k).swap(dst.row(j));
|
|
else if(Side==OnTheRight)
|
|
dst.col(k).swap(dst.col(j));
|
|
}
|
|
}
|
|
|
|
protected:
|
|
const TranspositionType& m_transpositions;
|
|
const typename MatrixType::Nested m_matrix;
|
|
};
|
|
|
|
} // end namespace internal
|
|
|
|
/* Template partial specialization for transposed/inverse transpositions */
|
|
|
|
template<int SizeAtCompileTime, int MaxSizeAtCompileTime>
|
|
class Transpose<Transpositions<SizeAtCompileTime, MaxSizeAtCompileTime> >
|
|
{
|
|
typedef Transpositions<SizeAtCompileTime, MaxSizeAtCompileTime> 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<typename Derived> friend
|
|
inline const internal::transposition_matrix_product_retval<TranspositionType, Derived, OnTheRight, true>
|
|
operator*(const MatrixBase<Derived>& matrix, const Transpose& trt)
|
|
{
|
|
return internal::transposition_matrix_product_retval<TranspositionType, Derived, OnTheRight, true>(trt.m_transpositions, matrix.derived());
|
|
}
|
|
|
|
/** \returns the \a matrix with the inverse transpositions applied to the rows.
|
|
*/
|
|
template<typename Derived>
|
|
inline const internal::transposition_matrix_product_retval<TranspositionType, Derived, OnTheLeft, true>
|
|
operator*(const MatrixBase<Derived>& matrix) const
|
|
{
|
|
return internal::transposition_matrix_product_retval<TranspositionType, Derived, OnTheLeft, true>(m_transpositions, matrix.derived());
|
|
}
|
|
|
|
const TranspositionType& nestedTranspositions() const { return m_transpositions; }
|
|
|
|
protected:
|
|
const TranspositionType& m_transpositions;
|
|
};
|
|
|
|
#endif // EIGEN_TRANSPOSITIONS_H
|