Home | History | Annotate | Download | only in Eigen2Support 
      1 // This file is part of Eigen, a lightweight C++ template library
      2 // for linear algebra.
      3 //
      4 // Copyright (C) 2011 Benoit Jacob <jacob.benoit.1 (at) gmail.com>
      5 //
      6 // This Source Code Form is subject to the terms of the Mozilla
      7 // Public License v. 2.0. If a copy of the MPL was not distributed
      8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
      9 
     10 #ifndef EIGEN2_LU_H
     11 #define EIGEN2_LU_H
     12 
     13 namespace Eigen {
     14 
     15 template<typename MatrixType>
     16 class LU : public FullPivLU<MatrixType>
     17 {
     18   public:
     19 
     20     typedef typename MatrixType::Scalar Scalar;
     21     typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar;
     22     typedef Matrix<int, 1, MatrixType::ColsAtCompileTime, MatrixType::Options, 1, MatrixType::MaxColsAtCompileTime> IntRowVectorType;
     23     typedef Matrix<int, MatrixType::RowsAtCompileTime, 1, MatrixType::Options, MatrixType::MaxRowsAtCompileTime, 1> IntColVectorType;
     24     typedef Matrix<Scalar, 1, MatrixType::ColsAtCompileTime, MatrixType::Options, 1, MatrixType::MaxColsAtCompileTime> RowVectorType;
     25     typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1, MatrixType::Options, MatrixType::MaxRowsAtCompileTime, 1> ColVectorType;
     26 
     27     typedef Matrix<typename MatrixType::Scalar,
     28                   MatrixType::ColsAtCompileTime, // the number of rows in the "kernel matrix" is the number of cols of the original matrix
     29                                                  // so that the product "matrix * kernel = zero" makes sense
     30                   Dynamic,                       // we don't know at compile-time the dimension of the kernel
     31                   MatrixType::Options,
     32                   MatrixType::MaxColsAtCompileTime, // see explanation for 2nd template parameter
     33                   MatrixType::MaxColsAtCompileTime // the kernel is a subspace of the domain space, whose dimension is the number
     34                                                    // of columns of the original matrix
     35     > KernelResultType;
     36 
     37     typedef Matrix<typename MatrixType::Scalar,
     38                    MatrixType::RowsAtCompileTime, // the image is a subspace of the destination space, whose dimension is the number
     39                                                   // of rows of the original matrix
     40                    Dynamic,                       // we don't know at compile time the dimension of the image (the rank)
     41                    MatrixType::Options,
     42                    MatrixType::MaxRowsAtCompileTime, // the image matrix will consist of columns from the original matrix,
     43                    MatrixType::MaxColsAtCompileTime  // so it has the same number of rows and at most as many columns.
     44     > ImageResultType;
     45 
     46     typedef FullPivLU<MatrixType> Base;
     47 
     48     template<typename T>
     49     explicit LU(const T& t) : Base(t), m_originalMatrix(t) {}
     50 
     51     template<typename OtherDerived, typename ResultType>
     52     bool solve(const MatrixBase<OtherDerived>& b, ResultType *result) const
     53     {
     54       *result = static_cast<const Base*>(this)->solve(b);
     55       return true;
     56     }
     57 
     58     template<typename ResultType>
     59     inline void computeInverse(ResultType *result) const
     60     {
     61       solve(MatrixType::Identity(this->rows(), this->cols()), result);
     62     }
     63 
     64     template<typename KernelMatrixType>
     65     void computeKernel(KernelMatrixType *result) const
     66     {
     67       *result = static_cast<const Base*>(this)->kernel();
     68     }
     69 
     70     template<typename ImageMatrixType>
     71     void computeImage(ImageMatrixType *result) const
     72     {
     73       *result = static_cast<const Base*>(this)->image(m_originalMatrix);
     74     }
     75 
     76     const ImageResultType image() const
     77     {
     78       return static_cast<const Base*>(this)->image(m_originalMatrix);
     79     }
     80 
     81     const MatrixType& m_originalMatrix;
     82 };
     83 
     84 #if EIGEN2_SUPPORT_STAGE < STAGE20_RESOLVE_API_CONFLICTS
     85 /** \lu_module
     86   *
     87   * Synonym of partialPivLu().
     88   *
     89   * \return the partial-pivoting LU decomposition of \c *this.
     90   *
     91   * \sa class PartialPivLU
     92   */
     93 template<typename Derived>
     94 inline const LU<typename MatrixBase<Derived>::PlainObject>
     95 MatrixBase<Derived>::lu() const
     96 {
     97   return LU<PlainObject>(eval());
     98 }
     99 #endif
    100 
    101 #ifdef EIGEN2_SUPPORT
    102 /** \lu_module
    103   *
    104   * Synonym of partialPivLu().
    105   *
    106   * \return the partial-pivoting LU decomposition of \c *this.
    107   *
    108   * \sa class PartialPivLU
    109   */
    110 template<typename Derived>
    111 inline const LU<typename MatrixBase<Derived>::PlainObject>
    112 MatrixBase<Derived>::eigen2_lu() const
    113 {
    114   return LU<PlainObject>(eval());
    115 }
    116 #endif
    117 
    118 } // end namespace Eigen
    119 
    120 #endif // EIGEN2_LU_H
    121