Attachment #714 for bug #707

This bugzilla service is closed. All entries have been migrated to https://gitlab.com/libeigen/eigen

View | Details | Raw Unified | Return to bug 707
Collapse All | Expand All




  */
template<typename _MatrixType, int _UpLo> class LDLT
{
  public:
    typedef _MatrixType MatrixType;
    enum {
      RowsAtCompileTime = MatrixType::RowsAtCompileTime,
      ColsAtCompileTime = MatrixType::ColsAtCompileTime,
      Options = MatrixType::Options & ~RowMajorBit, // these are the options for the TmpMatrixType, we need a ColMajor matrix here!
      MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
      MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
      UpLo = _UpLo
    };
    typedef typename MatrixType::Scalar Scalar;
    typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar;
    typedef Eigen::Index Index; ///< \deprecated since Eigen 3.3
    typedef typename MatrixType::StorageIndex StorageIndex;
    typedef Matrix<Scalar, RowsAtCompileTime, 1, 0, MaxRowsAtCompileTime, 1> TmpMatrixType;

    typedef Transpositions<RowsAtCompileTime, MaxRowsAtCompileTime> TranspositionType;
    typedef PermutationMatrix<RowsAtCompileTime, MaxRowsAtCompileTime> PermutationType;

    typedef internal::LDLT_Traits<MatrixType,UpLo> Traits;

    /** \brief Default Constructor.
      *

        m_temporary(size),
        m_sign(internal::ZeroSign),
        m_isInitialized(false)
    {}

    /** \brief Constructor with decomposition
      *
      * This calculates the decomposition for the input \a matrix.
      *
      * \sa LDLT(Index size)
      */
    template<typename InputType>
    explicit LDLT(const EigenBase<InputType>& matrix)
      : m_matrix(matrix.rows(), matrix.cols()),
        m_transpositions(matrix.rows()),
        m_temporary(matrix.rows()),
        m_sign(internal::ZeroSign),
        m_isInitialized(false)
    {
      compute(matrix.derived());
    }

    /** \brief Constructs a LDLT factorization from a given matrix
      *
      * This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref.
      *
      * \sa LDLT(const EigenBase&)
      */
    template<typename InputType>
    explicit LDLT(EigenBase<InputType>& matrix)
      : m_matrix(matrix.derived()),
        m_transpositions(matrix.rows()),
        m_temporary(matrix.rows()),
        m_sign(internal::ZeroSign),
        m_isInitialized(false)
    {
      compute(matrix.derived());
    }

    /** Clear any existing decomposition
     * \sa rankUpdate(w,sigma)
     */
    void setZero()
    {
      m_isInitialized = false;
    }





  */
template<typename _MatrixType, int _UpLo> class LLT
{
  public:
    typedef _MatrixType MatrixType;
    enum {
      RowsAtCompileTime = MatrixType::RowsAtCompileTime,
      ColsAtCompileTime = MatrixType::ColsAtCompileTime,
      Options = MatrixType::Options,
      MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
    };
    typedef typename MatrixType::Scalar Scalar;
    typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar;
    typedef Eigen::Index Index; ///< \deprecated since Eigen 3.3
    typedef typename MatrixType::StorageIndex StorageIndex;

    enum {

    template<typename InputType>
    explicit LLT(const EigenBase<InputType>& matrix)
      : m_matrix(matrix.rows(), matrix.cols()),
        m_isInitialized(false)
    {
      compute(matrix.derived());
    }

    /** \brief Constructs a LDLT factorization from a given matrix
      *
      * This overloaded constructor is provided for inplace solving when
      * \c MatrixType is a Eigen::Ref.
      *
      * \sa LLT(const EigenBase&)
      */
    template<typename InputType>
    explicit LLT(EigenBase<InputType>& matrix)
      : m_matrix(matrix.derived()),
        m_isInitialized(false)
    {
      compute(matrix.derived());
    }

    /** \returns a view of the upper triangular matrix U */
    inline typename Traits::MatrixU matrixU() const
    {
      eigen_assert(m_isInitialized && "LLT is not initialized.");
      return Traits::getU(m_matrix);
    }

    /** \returns a view of the lower triangular matrix L */




    /** Constructor.
      *
      * \param matrix the matrix of which to compute the LU decomposition.
      *               It is required to be nonzero.
      */
    template<typename InputType>
    explicit FullPivLU(const EigenBase<InputType>& matrix);

    /** \brief Constructs a LU factorization from a given matrix
      *
      * This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref.
      *
      * \sa FullPivLU(const EigenBase&)
      */
    template<typename InputType>
    explicit FullPivLU(EigenBase<InputType>& matrix);

    /** Computes the LU decomposition of the given matrix.
      *
      * \param matrix the matrix of which to compute the LU decomposition.
      *               It is required to be nonzero.
      *
      * \returns a reference to *this
      */
    template<typename InputType>

    m_isInitialized(false),
    m_usePrescribedThreshold(false)
{
  compute(matrix.derived());
}

template<typename MatrixType>
template<typename InputType>
FullPivLU<MatrixType>::FullPivLU(EigenBase<InputType>& matrix)
  : m_lu(matrix.derived()),
    m_p(matrix.rows()),
    m_q(matrix.cols()),
    m_rowsTranspositions(matrix.rows()),
    m_colsTranspositions(matrix.cols()),
    m_isInitialized(false),
    m_usePrescribedThreshold(false)
{
  compute(matrix.derived());
}

template<typename MatrixType>
template<typename InputType>
FullPivLU<MatrixType>& FullPivLU<MatrixType>::compute(const EigenBase<InputType>& matrix)
{
  check_template_parameters();

  // the permutations are stored as int indices, so just to be sure:
  eigen_assert(matrix.rows()<=NumTraits<int>::highest() && matrix.cols()<=NumTraits<int>::highest());

  m_lu = matrix.derived();




      * \param matrix the matrix of which to compute the LU decomposition.
      *
      * \warning The matrix should have full rank (e.g. if it's square, it should be invertible).
      * If you need to deal with non-full rank, use class FullPivLU instead.
      */
    template<typename InputType>
    explicit PartialPivLU(const EigenBase<InputType>& matrix);

    /** Constructor for inplace decomposition
      *
      * \param matrix the matrix of which to compute the LU decomposition.
      *
      * If \c MatrixType is an Eigen::Ref, then the storage of \a matrix will be shared
      * between \a matrix and \c *this and the decomposition will take place in-place.
      * The memory of \a matrix will be used througrough the lifetime of \c *this. In
      * particular, further calls to \c this->compute(A) will still operate on the memory
      * of \a matrix meaning. This also implies that the sizes of \c A must match the
      * ones of \a matrix.
      *
      * \warning The matrix should have full rank (e.g. if it's square, it should be invertible).
      * If you need to deal with non-full rank, use class FullPivLU instead.
      */
    template<typename InputType>
    explicit PartialPivLU(EigenBase<InputType>& matrix);

    template<typename InputType>
    PartialPivLU& compute(const EigenBase<InputType>& matrix);

    /** \returns the LU decomposition matrix: the upper-triangular part is U, the
      * unit-lower-triangular part is L (at least for square matrices; in the non-square
      * case, special care is needed, see the documentation of class FullPivLU).
      *
      * \sa matrixL(), matrixU()

    m_det_p(0),
    m_isInitialized(false)
{
}

template<typename MatrixType>
template<typename InputType>
PartialPivLU<MatrixType>::PartialPivLU(const EigenBase<InputType>& matrix)
  : m_lu(matrix.rows(),matrix.cols()),
    m_p(matrix.rows()),
    m_rowsTranspositions(matrix.rows()),
    m_l1_norm(0),
    m_det_p(0),
    m_isInitialized(false)
{
  compute(matrix.derived());
}

template<typename MatrixType>
template<typename InputType>
PartialPivLU<MatrixType>::PartialPivLU(EigenBase<InputType>& matrix)
  : m_lu(matrix.derived()),
    m_p(matrix.rows()),
    m_rowsTranspositions(matrix.rows()),
    m_l1_norm(0),
    m_det_p(0),
    m_isInitialized(false)
{
  compute(matrix.derived());
}




template<typename _MatrixType> class ColPivHouseholderQR
{
  public:

    typedef _MatrixType MatrixType;
    enum {
      RowsAtCompileTime = MatrixType::RowsAtCompileTime,
      ColsAtCompileTime = MatrixType::ColsAtCompileTime,
      Options = MatrixType::Options,
      MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
      MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
    };
    typedef typename MatrixType::Scalar Scalar;
    typedef typename MatrixType::RealScalar RealScalar;
    // FIXME should be int
    typedef typename MatrixType::StorageIndex StorageIndex;
    typedef Matrix<Scalar, RowsAtCompileTime, RowsAtCompileTime, Options, MaxRowsAtCompileTime, MaxRowsAtCompileTime> MatrixQType;
    typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType;
    typedef PermutationMatrix<ColsAtCompileTime, MaxColsAtCompileTime> PermutationType;
    typedef typename internal::plain_row_type<MatrixType, Index>::type IntRowVectorType;
    typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;
    typedef typename internal::plain_row_type<MatrixType, RealScalar>::type RealRowVectorType;
    typedef HouseholderSequence<MatrixType,typename internal::remove_all<typename HCoeffsType::ConjugateReturnType>::type> HouseholderSequenceType;
    typedef typename MatrixType::PlainObject PlainObject;


        m_colNormsUpdated(matrix.cols()),
        m_colNormsDirect(matrix.cols()),
        m_isInitialized(false),
        m_usePrescribedThreshold(false)
    {
      compute(matrix.derived());
    }

    /** \brief Constructs a QR factorization from a given matrix
      *
      * This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref.
      *
      * \sa ColPivHouseholderQR(const EigenBase&)
      */
    template<typename InputType>
    explicit ColPivHouseholderQR(EigenBase<InputType>& matrix)
      : m_qr(matrix.derived()),
        m_hCoeffs((std::min)(matrix.rows(),matrix.cols())),
        m_colsPermutation(PermIndexType(matrix.cols())),
        m_colsTranspositions(matrix.cols()),
        m_temp(matrix.cols()),
        m_colNormsUpdated(matrix.cols()),
        m_colNormsDirect(matrix.cols()),
        m_isInitialized(false),
        m_usePrescribedThreshold(false)
    {
      compute(matrix.derived());
    }

    /** This method finds a solution x to the equation Ax=b, where A is the matrix of which
      * *this is the QR decomposition, if any exists.
      *
      * \param b the right-hand-side of the equation to solve.
      *
      * \returns a solution.
      *
      * \note The case where b is a matrix is not yet implemented. Also, this




template<typename _MatrixType> class FullPivHouseholderQR
{
  public:

    typedef _MatrixType MatrixType;
    enum {
      RowsAtCompileTime = MatrixType::RowsAtCompileTime,
      ColsAtCompileTime = MatrixType::ColsAtCompileTime,
      Options = MatrixType::Options,
      MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
      MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
    };
    typedef typename MatrixType::Scalar Scalar;
    typedef typename MatrixType::RealScalar RealScalar;
    // FIXME should be int
    typedef typename MatrixType::StorageIndex StorageIndex;
    typedef internal::FullPivHouseholderQRMatrixQReturnType<MatrixType> MatrixQReturnType;

        m_cols_permutation(matrix.cols()),
        m_temp(matrix.cols()),
        m_isInitialized(false),
        m_usePrescribedThreshold(false)
    {
      compute(matrix.derived());
    }

    /** \brief Constructs a QR factorization from a given matrix
      *
      * This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref.
      *
      * \sa FullPivHouseholderQR(const EigenBase&)
      */
    template<typename InputType>
    explicit FullPivHouseholderQR(EigenBase<InputType>& matrix)
      : m_qr(matrix.derived()),
        m_hCoeffs((std::min)(matrix.rows(), matrix.cols())),
        m_rows_transpositions((std::min)(matrix.rows(), matrix.cols())),
        m_cols_transpositions((std::min)(matrix.rows(), matrix.cols())),
        m_cols_permutation(matrix.cols()),
        m_temp(matrix.cols()),
        m_isInitialized(false),
        m_usePrescribedThreshold(false)
    {
      compute(matrix.derived());
    }

    /** This method finds a solution x to the equation Ax=b, where A is the matrix of which
      * \c *this is the QR decomposition.
      *
      * \param b the right-hand-side of the equation to solve.
      *
      * \returns the exact or least-square solution if the rank is greater or equal to the number of columns of A,
      * and an arbitrary solution otherwise.
      *




template<typename _MatrixType> class HouseholderQR
{
  public:

    typedef _MatrixType MatrixType;
    enum {
      RowsAtCompileTime = MatrixType::RowsAtCompileTime,
      ColsAtCompileTime = MatrixType::ColsAtCompileTime,
      Options = MatrixType::Options,
      MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
      MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
    };
    typedef typename MatrixType::Scalar Scalar;
    typedef typename MatrixType::RealScalar RealScalar;
    // FIXME should be int
    typedef typename MatrixType::StorageIndex StorageIndex;
    typedef Matrix<Scalar, RowsAtCompileTime, RowsAtCompileTime, (MatrixType::Flags&RowMajorBit) ? RowMajor : ColMajor, MaxRowsAtCompileTime, MaxRowsAtCompileTime> MatrixQType;

      : m_qr(matrix.rows(), matrix.cols()),
        m_hCoeffs((std::min)(matrix.rows(),matrix.cols())),
        m_temp(matrix.cols()),
        m_isInitialized(false)
    {
      compute(matrix.derived());
    }


    /** \brief Constructs a QR factorization from a given matrix
      *
      * This overloaded constructor is provided for inplace solving when
      * \c MatrixType is a Eigen::Ref.
      *
      * \sa HouseholderQR(const EigenBase&)
      */
    template<typename InputType>
    explicit HouseholderQR(EigenBase<InputType>& matrix)
      : m_qr(matrix.derived()),
        m_hCoeffs((std::min)(matrix.rows(),matrix.cols())),
        m_temp(matrix.cols()),
        m_isInitialized(false)
    {
      compute(matrix.derived());
    }

    /** This method finds a solution x to the equation Ax=b, where A is the matrix of which
      * *this is the QR decomposition, if any exists.
      *
      * \param b the right-hand-side of the equation to solve.
      *
      * \returns a solution.
      *
      * \note The case where b is a matrix is not yet implemented. Also, this

Return to bug 707

Lines 47-72 namespace internal { Link Here

(-)a/Eigen/src/Cholesky/LDLT.h (-2 / +19 lines)
47	*/	47	*/
48	template<typename _MatrixType, int _UpLo> class LDLT	48	template<typename _MatrixType, int _UpLo> class LDLT
49	{	49	{
50	public:	50	public:
51	typedef _MatrixType MatrixType;	51	typedef _MatrixType MatrixType;
52	enum {	52	enum {
53	RowsAtCompileTime = MatrixType::RowsAtCompileTime,	53	RowsAtCompileTime = MatrixType::RowsAtCompileTime,
54	ColsAtCompileTime = MatrixType::ColsAtCompileTime,	54	ColsAtCompileTime = MatrixType::ColsAtCompileTime,
55	Options = MatrixType::Options & ~RowMajorBit, // these are the options for the TmpMatrixType, we need a ColMajor matrix here!
56	MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,	55	MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
57	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,	56	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
58	UpLo = _UpLo	57	UpLo = _UpLo
59	};	58	};
60	typedef typename MatrixType::Scalar Scalar;	59	typedef typename MatrixType::Scalar Scalar;
61	typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar;	60	typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar;
62	typedef Eigen::Index Index; ///< \deprecated since Eigen 3.3	61	typedef Eigen::Index Index; ///< \deprecated since Eigen 3.3
63	typedef typename MatrixType::StorageIndex StorageIndex;	62	typedef typename MatrixType::StorageIndex StorageIndex;
64	typedef Matrix<Scalar, RowsAtCompileTime, 1, Options, MaxRowsAtCompileTime, 1> TmpMatrixType;	63	typedef Matrix<Scalar, RowsAtCompileTime, 1, 0, MaxRowsAtCompileTime, 1> TmpMatrixType;
65		64
66	typedef Transpositions<RowsAtCompileTime, MaxRowsAtCompileTime> TranspositionType;	65	typedef Transpositions<RowsAtCompileTime, MaxRowsAtCompileTime> TranspositionType;
67	typedef PermutationMatrix<RowsAtCompileTime, MaxRowsAtCompileTime> PermutationType;	66	typedef PermutationMatrix<RowsAtCompileTime, MaxRowsAtCompileTime> PermutationType;
68		67
69	typedef internal::LDLT_Traits<MatrixType,UpLo> Traits;	68	typedef internal::LDLT_Traits<MatrixType,UpLo> Traits;
70		69
71	/** \brief Default Constructor.	70	/** \brief Default Constructor.
72	*	71	*
Lines 92-120 template<typename _MatrixType, int _UpLo Link Here
92	m_temporary(size),	91	m_temporary(size),
93	m_sign(internal::ZeroSign),	92	m_sign(internal::ZeroSign),
94	m_isInitialized(false)	93	m_isInitialized(false)
95	{}	94	{}
96		95
97	/** \brief Constructor with decomposition	96	/** \brief Constructor with decomposition
98	*	97	*
99	* This calculates the decomposition for the input \a matrix.	98	* This calculates the decomposition for the input \a matrix.
		99	*
100	* \sa LDLT(Index size)	100	* \sa LDLT(Index size)
101	*/	101	*/
102	template<typename InputType>	102	template<typename InputType>
103	explicit LDLT(const EigenBase<InputType>& matrix)	103	explicit LDLT(const EigenBase<InputType>& matrix)
104	: m_matrix(matrix.rows(), matrix.cols()),	104	: m_matrix(matrix.rows(), matrix.cols()),
105	m_transpositions(matrix.rows()),	105	m_transpositions(matrix.rows()),
106	m_temporary(matrix.rows()),	106	m_temporary(matrix.rows()),
107	m_sign(internal::ZeroSign),	107	m_sign(internal::ZeroSign),
108	m_isInitialized(false)	108	m_isInitialized(false)
109	{	109	{
110	compute(matrix.derived());	110	compute(matrix.derived());
111	}	111	}
112		112
		113	/** \brief Constructs a LDLT factorization from a given matrix
		114	*
		115	* This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref.
		116	*
		117	* \sa LDLT(const EigenBase&)
		118	*/
		119	template<typename InputType>
		120	explicit LDLT(EigenBase<InputType>& matrix)
		121	: m_matrix(matrix.derived()),
		122	m_transpositions(matrix.rows()),
		123	m_temporary(matrix.rows()),
		124	m_sign(internal::ZeroSign),
		125	m_isInitialized(false)
		126	{
		127	compute(matrix.derived());
		128	}
		129
113	/** Clear any existing decomposition	130	/** Clear any existing decomposition
114	* \sa rankUpdate(w,sigma)	131	* \sa rankUpdate(w,sigma)
115	*/	132	*/
116	void setZero()	133	void setZero()
117	{	134	{
118	m_isInitialized = false;	135	m_isInitialized = false;
119	}	136	}
120		137

Lines 49-65 template<typename MatrixType, int UpLo> Link Here

(-)a/Eigen/src/Cholesky/LLT.h (-1 / +15 lines)
49	*/	49	*/
50	template<typename _MatrixType, int _UpLo> class LLT	50	template<typename _MatrixType, int _UpLo> class LLT
51	{	51	{
52	public:	52	public:
53	typedef _MatrixType MatrixType;	53	typedef _MatrixType MatrixType;
54	enum {	54	enum {
55	RowsAtCompileTime = MatrixType::RowsAtCompileTime,	55	RowsAtCompileTime = MatrixType::RowsAtCompileTime,
56	ColsAtCompileTime = MatrixType::ColsAtCompileTime,	56	ColsAtCompileTime = MatrixType::ColsAtCompileTime,
57	Options = MatrixType::Options,
58	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime	57	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
59	};	58	};
60	typedef typename MatrixType::Scalar Scalar;	59	typedef typename MatrixType::Scalar Scalar;
61	typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar;	60	typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar;
62	typedef Eigen::Index Index; ///< \deprecated since Eigen 3.3	61	typedef Eigen::Index Index; ///< \deprecated since Eigen 3.3
63	typedef typename MatrixType::StorageIndex StorageIndex;	62	typedef typename MatrixType::StorageIndex StorageIndex;
64		63
65	enum {	64	enum {
Lines 90-105 template<typename _MatrixType, int _UpLo Link Here
90	template<typename InputType>	89	template<typename InputType>
91	explicit LLT(const EigenBase<InputType>& matrix)	90	explicit LLT(const EigenBase<InputType>& matrix)
92	: m_matrix(matrix.rows(), matrix.cols()),	91	: m_matrix(matrix.rows(), matrix.cols()),
93	m_isInitialized(false)	92	m_isInitialized(false)
94	{	93	{
95	compute(matrix.derived());	94	compute(matrix.derived());
96	}	95	}
97		96
		97	/** \brief Constructs a LDLT factorization from a given matrix
		98	*
		99	* This overloaded constructor is provided for inplace solving when
		100	* \c MatrixType is a Eigen::Ref.
		101	*
		102	* \sa LLT(const EigenBase&)
		103	*/
		104	template<typename InputType>
		105	explicit LLT(EigenBase<InputType>& matrix)
		106	: m_matrix(matrix.derived()),
		107	m_isInitialized(false)
		108	{
		109	compute(matrix.derived());
		110	}
		111
98	/** \returns a view of the upper triangular matrix U */	112	/** \returns a view of the upper triangular matrix U */
99	inline typename Traits::MatrixU matrixU() const	113	inline typename Traits::MatrixU matrixU() const
100	{	114	{
101	eigen_assert(m_isInitialized && "LLT is not initialized.");	115	eigen_assert(m_isInitialized && "LLT is not initialized.");
102	return Traits::getU(m_matrix);	116	return Traits::getU(m_matrix);
103	}	117	}
104		118
105	/** \returns a view of the lower triangular matrix L */	119	/** \returns a view of the lower triangular matrix L */

Lines 92-107 template<typename _MatrixType> class Ful Link Here

(-)a/Eigen/src/LU/FullPivLU.h (+23 lines)
92	/** Constructor.	92	/** Constructor.
93	*	93	*
94	* \param matrix the matrix of which to compute the LU decomposition.	94	* \param matrix the matrix of which to compute the LU decomposition.
95	* It is required to be nonzero.	95	* It is required to be nonzero.
96	*/	96	*/
97	template<typename InputType>	97	template<typename InputType>
98	explicit FullPivLU(const EigenBase<InputType>& matrix);	98	explicit FullPivLU(const EigenBase<InputType>& matrix);
99		99
		100	/** \brief Constructs a LU factorization from a given matrix
		101	*
		102	* This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref.
		103	*
		104	* \sa FullPivLU(const EigenBase&)
		105	*/
		106	template<typename InputType>
		107	explicit FullPivLU(EigenBase<InputType>& matrix);
		108
100	/** Computes the LU decomposition of the given matrix.	109	/** Computes the LU decomposition of the given matrix.
101	*	110	*
102	* \param matrix the matrix of which to compute the LU decomposition.	111	* \param matrix the matrix of which to compute the LU decomposition.
103	* It is required to be nonzero.	112	* It is required to be nonzero.
104	*	113	*
105	* \returns a reference to *this	114	* \returns a reference to *this
106	*/	115	*/
107	template<typename InputType>	116	template<typename InputType>
Lines 454-469 FullPivLU<MatrixType>::FullPivLU(const E Link Here
454	m_isInitialized(false),	463	m_isInitialized(false),
455	m_usePrescribedThreshold(false)	464	m_usePrescribedThreshold(false)
456	{	465	{
457	compute(matrix.derived());	466	compute(matrix.derived());
458	}	467	}
459		468
460	template<typename MatrixType>	469	template<typename MatrixType>
461	template<typename InputType>	470	template<typename InputType>
		471	FullPivLU<MatrixType>::FullPivLU(EigenBase<InputType>& matrix)
		472	: m_lu(matrix.derived()),
		473	m_p(matrix.rows()),
		474	m_q(matrix.cols()),
		475	m_rowsTranspositions(matrix.rows()),
		476	m_colsTranspositions(matrix.cols()),
		477	m_isInitialized(false),
		478	m_usePrescribedThreshold(false)
		479	{
		480	compute(matrix.derived());
		481	}
		482
		483	template<typename MatrixType>
		484	template<typename InputType>
462	FullPivLU<MatrixType>& FullPivLU<MatrixType>::compute(const EigenBase<InputType>& matrix)	485	FullPivLU<MatrixType>& FullPivLU<MatrixType>::compute(const EigenBase<InputType>& matrix)
463	{	486	{
464	check_template_parameters();	487	check_template_parameters();
465		488
466	// the permutations are stored as int indices, so just to be sure:	489	// the permutations are stored as int indices, so just to be sure:
467	eigen_assert(matrix.rows()<=NumTraits<int>::highest() && matrix.cols()<=NumTraits<int>::highest());	490	eigen_assert(matrix.rows()<=NumTraits<int>::highest() && matrix.cols()<=NumTraits<int>::highest());
468		491
469	m_lu = matrix.derived();	492	m_lu = matrix.derived();

Lines 97-112 template<typename _MatrixType> class Par Link Here

(-)a/Eigen/src/LU/PartialPivLU.h (-1 / +31 lines)
97	* \param matrix the matrix of which to compute the LU decomposition.	97	* \param matrix the matrix of which to compute the LU decomposition.
98	*	98	*
99	* \warning The matrix should have full rank (e.g. if it's square, it should be invertible).	99	* \warning The matrix should have full rank (e.g. if it's square, it should be invertible).
100	* If you need to deal with non-full rank, use class FullPivLU instead.	100	* If you need to deal with non-full rank, use class FullPivLU instead.
101	*/	101	*/
102	template<typename InputType>	102	template<typename InputType>
103	explicit PartialPivLU(const EigenBase<InputType>& matrix);	103	explicit PartialPivLU(const EigenBase<InputType>& matrix);
104		104
		105	/** Constructor for inplace decomposition
		106	*
		107	* \param matrix the matrix of which to compute the LU decomposition.
		108	*
		109	* If \c MatrixType is an Eigen::Ref, then the storage of \a matrix will be shared
		110	* between \a matrix and \c *this and the decomposition will take place in-place.
		111	* The memory of \a matrix will be used througrough the lifetime of \c *this. In
		112	* particular, further calls to \c this->compute(A) will still operate on the memory
		113	* of \a matrix meaning. This also implies that the sizes of \c A must match the
		114	* ones of \a matrix.
		115	*
		116	* \warning The matrix should have full rank (e.g. if it's square, it should be invertible).
		117	* If you need to deal with non-full rank, use class FullPivLU instead.
		118	*/
		119	template<typename InputType>
		120	explicit PartialPivLU(EigenBase<InputType>& matrix);
		121
105	template<typename InputType>	122	template<typename InputType>
106	PartialPivLU& compute(const EigenBase<InputType>& matrix);	123	PartialPivLU& compute(const EigenBase<InputType>& matrix);
107		124
108	/** \returns the LU decomposition matrix: the upper-triangular part is U, the	125	/** \returns the LU decomposition matrix: the upper-triangular part is U, the
109	* unit-lower-triangular part is L (at least for square matrices; in the non-square	126	* unit-lower-triangular part is L (at least for square matrices; in the non-square
110	* case, special care is needed, see the documentation of class FullPivLU).	127	* case, special care is needed, see the documentation of class FullPivLU).
111	*	128	*
112	* \sa matrixL(), matrixU()	129	* \sa matrixL(), matrixU()
Lines 279-295 PartialPivLU<MatrixType>::PartialPivLU(I Link Here
279	m_det_p(0),	296	m_det_p(0),
280	m_isInitialized(false)	297	m_isInitialized(false)
281	{	298	{
282	}	299	}
283		300
284	template<typename MatrixType>	301	template<typename MatrixType>
285	template<typename InputType>	302	template<typename InputType>
286	PartialPivLU<MatrixType>::PartialPivLU(const EigenBase<InputType>& matrix)	303	PartialPivLU<MatrixType>::PartialPivLU(const EigenBase<InputType>& matrix)
287	: m_lu(matrix.rows(), matrix.rows()),	304	: m_lu(matrix.rows(),matrix.cols()),
		305	m_p(matrix.rows()),
		306	m_rowsTranspositions(matrix.rows()),
		307	m_l1_norm(0),
		308	m_det_p(0),
		309	m_isInitialized(false)
		310	{
		311	compute(matrix.derived());
		312	}
		313
		314	template<typename MatrixType>
		315	template<typename InputType>
		316	PartialPivLU<MatrixType>::PartialPivLU(EigenBase<InputType>& matrix)
		317	: m_lu(matrix.derived()),
288	m_p(matrix.rows()),	318	m_p(matrix.rows()),
289	m_rowsTranspositions(matrix.rows()),	319	m_rowsTranspositions(matrix.rows()),
290	m_l1_norm(0),	320	m_l1_norm(0),
291	m_det_p(0),	321	m_det_p(0),
292	m_isInitialized(false)	322	m_isInitialized(false)
293	{	323	{
294	compute(matrix.derived());	324	compute(matrix.derived());
295	}	325	}

Lines 46-70 template<typename _MatrixType> struct tr Link Here

(-)a/Eigen/src/QR/ColPivHouseholderQR.h (-2 / +21 lines)
46	template<typename _MatrixType> class ColPivHouseholderQR	46	template<typename _MatrixType> class ColPivHouseholderQR
47	{	47	{
48	public:	48	public:
49		49
50	typedef _MatrixType MatrixType;	50	typedef _MatrixType MatrixType;
51	enum {	51	enum {
52	RowsAtCompileTime = MatrixType::RowsAtCompileTime,	52	RowsAtCompileTime = MatrixType::RowsAtCompileTime,
53	ColsAtCompileTime = MatrixType::ColsAtCompileTime,	53	ColsAtCompileTime = MatrixType::ColsAtCompileTime,
54	Options = MatrixType::Options,
55	MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,	54	MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
56	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime	55	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
57	};	56	};
58	typedef typename MatrixType::Scalar Scalar;	57	typedef typename MatrixType::Scalar Scalar;
59	typedef typename MatrixType::RealScalar RealScalar;	58	typedef typename MatrixType::RealScalar RealScalar;
60	// FIXME should be int	59	// FIXME should be int
61	typedef typename MatrixType::StorageIndex StorageIndex;	60	typedef typename MatrixType::StorageIndex StorageIndex;
62	typedef Matrix<Scalar, RowsAtCompileTime, RowsAtCompileTime, Options, MaxRowsAtCompileTime, MaxRowsAtCompileTime> MatrixQType;
63	typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType;	61	typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType;
64	typedef PermutationMatrix<ColsAtCompileTime, MaxColsAtCompileTime> PermutationType;	62	typedef PermutationMatrix<ColsAtCompileTime, MaxColsAtCompileTime> PermutationType;
65	typedef typename internal::plain_row_type<MatrixType, Index>::type IntRowVectorType;	63	typedef typename internal::plain_row_type<MatrixType, Index>::type IntRowVectorType;
66	typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;	64	typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;
67	typedef typename internal::plain_row_type<MatrixType, RealScalar>::type RealRowVectorType;	65	typedef typename internal::plain_row_type<MatrixType, RealScalar>::type RealRowVectorType;
68	typedef HouseholderSequence<MatrixType,typename internal::remove_all<typename HCoeffsType::ConjugateReturnType>::type> HouseholderSequenceType;	66	typedef HouseholderSequence<MatrixType,typename internal::remove_all<typename HCoeffsType::ConjugateReturnType>::type> HouseholderSequenceType;
69	typedef typename MatrixType::PlainObject PlainObject;	67	typedef typename MatrixType::PlainObject PlainObject;
70		68
Lines 130-145 template<typename _MatrixType> class Col Link Here
130	m_colNormsUpdated(matrix.cols()),	128	m_colNormsUpdated(matrix.cols()),
131	m_colNormsDirect(matrix.cols()),	129	m_colNormsDirect(matrix.cols()),
132	m_isInitialized(false),	130	m_isInitialized(false),
133	m_usePrescribedThreshold(false)	131	m_usePrescribedThreshold(false)
134	{	132	{
135	compute(matrix.derived());	133	compute(matrix.derived());
136	}	134	}
137		135
		136	/** \brief Constructs a QR factorization from a given matrix
		137	*
		138	* This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref.
		139	*
		140	* \sa ColPivHouseholderQR(const EigenBase&)
		141	*/
		142	template<typename InputType>
		143	explicit ColPivHouseholderQR(EigenBase<InputType>& matrix)
		144	: m_qr(matrix.derived()),
		145	m_hCoeffs((std::min)(matrix.rows(),matrix.cols())),
		146	m_colsPermutation(PermIndexType(matrix.cols())),
		147	m_colsTranspositions(matrix.cols()),
		148	m_temp(matrix.cols()),
		149	m_colNormsUpdated(matrix.cols()),
		150	m_colNormsDirect(matrix.cols()),
		151	m_isInitialized(false),
		152	m_usePrescribedThreshold(false)
		153	{
		154	compute(matrix.derived());
		155	}
		156
138	/** This method finds a solution x to the equation Ax=b, where A is the matrix of which	157	/** This method finds a solution x to the equation Ax=b, where A is the matrix of which
139	* *this is the QR decomposition, if any exists.	158	* *this is the QR decomposition, if any exists.
140	*	159	*
141	* \param b the right-hand-side of the equation to solve.	160	* \param b the right-hand-side of the equation to solve.
142	*	161	*
143	* \returns a solution.	162	* \returns a solution.
144	*	163	*
145	* \note The case where b is a matrix is not yet implemented. Also, this	164	* \note The case where b is a matrix is not yet implemented. Also, this

Lines 55-71 struct traits<FullPivHouseholderQRMatrix Link Here

(-)a/Eigen/src/QR/FullPivHouseholderQR.h (-1 / +20 lines)
55	template<typename _MatrixType> class FullPivHouseholderQR	55	template<typename _MatrixType> class FullPivHouseholderQR
56	{	56	{
57	public:	57	public:
58		58
59	typedef _MatrixType MatrixType;	59	typedef _MatrixType MatrixType;
60	enum {	60	enum {
61	RowsAtCompileTime = MatrixType::RowsAtCompileTime,	61	RowsAtCompileTime = MatrixType::RowsAtCompileTime,
62	ColsAtCompileTime = MatrixType::ColsAtCompileTime,	62	ColsAtCompileTime = MatrixType::ColsAtCompileTime,
63	Options = MatrixType::Options,
64	MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,	63	MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
65	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime	64	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
66	};	65	};
67	typedef typename MatrixType::Scalar Scalar;	66	typedef typename MatrixType::Scalar Scalar;
68	typedef typename MatrixType::RealScalar RealScalar;	67	typedef typename MatrixType::RealScalar RealScalar;
69	// FIXME should be int	68	// FIXME should be int
70	typedef typename MatrixType::StorageIndex StorageIndex;	69	typedef typename MatrixType::StorageIndex StorageIndex;
71	typedef internal::FullPivHouseholderQRMatrixQReturnType<MatrixType> MatrixQReturnType;	70	typedef internal::FullPivHouseholderQRMatrixQReturnType<MatrixType> MatrixQReturnType;
Lines 130-145 template<typename _MatrixType> class Ful Link Here
130	m_cols_permutation(matrix.cols()),	129	m_cols_permutation(matrix.cols()),
131	m_temp(matrix.cols()),	130	m_temp(matrix.cols()),
132	m_isInitialized(false),	131	m_isInitialized(false),
133	m_usePrescribedThreshold(false)	132	m_usePrescribedThreshold(false)
134	{	133	{
135	compute(matrix.derived());	134	compute(matrix.derived());
136	}	135	}
137		136
		137	/** \brief Constructs a QR factorization from a given matrix
		138	*
		139	* This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref.
		140	*
		141	* \sa FullPivHouseholderQR(const EigenBase&)
		142	*/
		143	template<typename InputType>
		144	explicit FullPivHouseholderQR(EigenBase<InputType>& matrix)
		145	: m_qr(matrix.derived()),
		146	m_hCoeffs((std::min)(matrix.rows(), matrix.cols())),
		147	m_rows_transpositions((std::min)(matrix.rows(), matrix.cols())),
		148	m_cols_transpositions((std::min)(matrix.rows(), matrix.cols())),
		149	m_cols_permutation(matrix.cols()),
		150	m_temp(matrix.cols()),
		151	m_isInitialized(false),
		152	m_usePrescribedThreshold(false)
		153	{
		154	compute(matrix.derived());
		155	}
		156
138	/** This method finds a solution x to the equation Ax=b, where A is the matrix of which	157	/** This method finds a solution x to the equation Ax=b, where A is the matrix of which
139	* \c *this is the QR decomposition.	158	* \c *this is the QR decomposition.
140	*	159	*
141	* \param b the right-hand-side of the equation to solve.	160	* \param b the right-hand-side of the equation to solve.
142	*	161	*
143	* \returns the exact or least-square solution if the rank is greater or equal to the number of columns of A,	162	* \returns the exact or least-square solution if the rank is greater or equal to the number of columns of A,
144	* and an arbitrary solution otherwise.	163	* and an arbitrary solution otherwise.
145	*	164	*

Lines 42-58 namespace Eigen { Link Here

(-)a/Eigen/src/QR/HouseholderQR.h (-1 / +18 lines)
42	template<typename _MatrixType> class HouseholderQR	42	template<typename _MatrixType> class HouseholderQR
43	{	43	{
44	public:	44	public:
45		45
46	typedef _MatrixType MatrixType;	46	typedef _MatrixType MatrixType;
47	enum {	47	enum {
48	RowsAtCompileTime = MatrixType::RowsAtCompileTime,	48	RowsAtCompileTime = MatrixType::RowsAtCompileTime,
49	ColsAtCompileTime = MatrixType::ColsAtCompileTime,	49	ColsAtCompileTime = MatrixType::ColsAtCompileTime,
50	Options = MatrixType::Options,
51	MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,	50	MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
52	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime	51	MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
53	};	52	};
54	typedef typename MatrixType::Scalar Scalar;	53	typedef typename MatrixType::Scalar Scalar;
55	typedef typename MatrixType::RealScalar RealScalar;	54	typedef typename MatrixType::RealScalar RealScalar;
56	// FIXME should be int	55	// FIXME should be int
57	typedef typename MatrixType::StorageIndex StorageIndex;	56	typedef typename MatrixType::StorageIndex StorageIndex;
58	typedef Matrix<Scalar, RowsAtCompileTime, RowsAtCompileTime, (MatrixType::Flags&RowMajorBit) ? RowMajor : ColMajor, MaxRowsAtCompileTime, MaxRowsAtCompileTime> MatrixQType;	57	typedef Matrix<Scalar, RowsAtCompileTime, RowsAtCompileTime, (MatrixType::Flags&RowMajorBit) ? RowMajor : ColMajor, MaxRowsAtCompileTime, MaxRowsAtCompileTime> MatrixQType;
Lines 97-112 template<typename _MatrixType> class Hou Link Here
97	: m_qr(matrix.rows(), matrix.cols()),	96	: m_qr(matrix.rows(), matrix.cols()),
98	m_hCoeffs((std::min)(matrix.rows(),matrix.cols())),	97	m_hCoeffs((std::min)(matrix.rows(),matrix.cols())),
99	m_temp(matrix.cols()),	98	m_temp(matrix.cols()),
100	m_isInitialized(false)	99	m_isInitialized(false)
101	{	100	{
102	compute(matrix.derived());	101	compute(matrix.derived());
103	}	102	}
104		103
		104
		105	/** \brief Constructs a QR factorization from a given matrix
		106	*
		107	* This overloaded constructor is provided for inplace solving when
		108	* \c MatrixType is a Eigen::Ref.
		109	*
		110	* \sa HouseholderQR(const EigenBase&)
		111	*/
		112	template<typename InputType>
		113	explicit HouseholderQR(EigenBase<InputType>& matrix)
		114	: m_qr(matrix.derived()),
		115	m_hCoeffs((std::min)(matrix.rows(),matrix.cols())),
		116	m_temp(matrix.cols()),
		117	m_isInitialized(false)
		118	{
		119	compute(matrix.derived());
		120	}
		121
105	/** This method finds a solution x to the equation Ax=b, where A is the matrix of which	122	/** This method finds a solution x to the equation Ax=b, where A is the matrix of which
106	* *this is the QR decomposition, if any exists.	123	* *this is the QR decomposition, if any exists.
107	*	124	*
108	* \param b the right-hand-side of the equation to solve.	125	* \param b the right-hand-side of the equation to solve.
109	*	126	*
110	* \returns a solution.	127	* \returns a solution.
111	*	128	*
112	* \note The case where b is a matrix is not yet implemented. Also, this	129	* \note The case where b is a matrix is not yet implemented. Also, this