Home | History | Annotate | Download | only in Core 
      1 // This file is part of Eigen, a lightweight C++ template library
      2 // for linear algebra.
      3 //
      4 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud (at) inria.fr>
      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 EIGEN_BANDMATRIX_H
     11 #define EIGEN_BANDMATRIX_H
     12 
     13 namespace Eigen {
     14 
     15 namespace internal {
     16 
     17 template<typename Derived>
     18 class BandMatrixBase : public EigenBase<Derived>
     19 {
     20   public:
     21 
     22     enum {
     23       Flags = internal::traits<Derived>::Flags,
     24       CoeffReadCost = internal::traits<Derived>::CoeffReadCost,
     25       RowsAtCompileTime = internal::traits<Derived>::RowsAtCompileTime,
     26       ColsAtCompileTime = internal::traits<Derived>::ColsAtCompileTime,
     27       MaxRowsAtCompileTime = internal::traits<Derived>::MaxRowsAtCompileTime,
     28       MaxColsAtCompileTime = internal::traits<Derived>::MaxColsAtCompileTime,
     29       Supers = internal::traits<Derived>::Supers,
     30       Subs   = internal::traits<Derived>::Subs,
     31       Options = internal::traits<Derived>::Options
     32     };
     33     typedef typename internal::traits<Derived>::Scalar Scalar;
     34     typedef Matrix<Scalar,RowsAtCompileTime,ColsAtCompileTime> DenseMatrixType;
     35     typedef typename DenseMatrixType::Index Index;
     36     typedef typename internal::traits<Derived>::CoefficientsType CoefficientsType;
     37     typedef EigenBase<Derived> Base;
     38 
     39   protected:
     40     enum {
     41       DataRowsAtCompileTime = ((Supers!=Dynamic) && (Subs!=Dynamic))
     42                             ? 1 + Supers + Subs
     43                             : Dynamic,
     44       SizeAtCompileTime = EIGEN_SIZE_MIN_PREFER_DYNAMIC(RowsAtCompileTime,ColsAtCompileTime)
     45     };
     46 
     47   public:
     48 
     49     using Base::derived;
     50     using Base::rows;
     51     using Base::cols;
     52 
     53     /** \returns the number of super diagonals */
     54     inline Index supers() const { return derived().supers(); }
     55 
     56     /** \returns the number of sub diagonals */
     57     inline Index subs() const { return derived().subs(); }
     58 
     59     /** \returns an expression of the underlying coefficient matrix */
     60     inline const CoefficientsType& coeffs() const { return derived().coeffs(); }
     61 
     62     /** \returns an expression of the underlying coefficient matrix */
     63     inline CoefficientsType& coeffs() { return derived().coeffs(); }
     64 
     65     /** \returns a vector expression of the \a i -th column,
     66       * only the meaningful part is returned.
     67       * \warning the internal storage must be column major. */
     68     inline Block<CoefficientsType,Dynamic,1> col(Index i)
     69     {
     70       EIGEN_STATIC_ASSERT((Options&RowMajor)==0,THIS_METHOD_IS_ONLY_FOR_COLUMN_MAJOR_MATRICES);
     71       Index start = 0;
     72       Index len = coeffs().rows();
     73       if (i<=supers())
     74       {
     75         start = supers()-i;
     76         len = (std::min)(rows(),std::max<Index>(0,coeffs().rows() - (supers()-i)));
     77       }
     78       else if (i>=rows()-subs())
     79         len = std::max<Index>(0,coeffs().rows() - (i + 1 - rows() + subs()));
     80       return Block<CoefficientsType,Dynamic,1>(coeffs(), start, i, len, 1);
     81     }
     82 
     83     /** \returns a vector expression of the main diagonal */
     84     inline Block<CoefficientsType,1,SizeAtCompileTime> diagonal()
     85     { return Block<CoefficientsType,1,SizeAtCompileTime>(coeffs(),supers(),0,1,(std::min)(rows(),cols())); }
     86 
     87     /** \returns a vector expression of the main diagonal (const version) */
     88     inline const Block<const CoefficientsType,1,SizeAtCompileTime> diagonal() const
     89     { return Block<const CoefficientsType,1,SizeAtCompileTime>(coeffs(),supers(),0,1,(std::min)(rows(),cols())); }
     90 
     91     template<int Index> struct DiagonalIntReturnType {
     92       enum {
     93         ReturnOpposite = (Options&SelfAdjoint) && (((Index)>0 && Supers==0) || ((Index)<0 && Subs==0)),
     94         Conjugate = ReturnOpposite && NumTraits<Scalar>::IsComplex,
     95         ActualIndex = ReturnOpposite ? -Index : Index,
     96         DiagonalSize = (RowsAtCompileTime==Dynamic || ColsAtCompileTime==Dynamic)
     97                      ? Dynamic
     98                      : (ActualIndex<0
     99                      ? EIGEN_SIZE_MIN_PREFER_DYNAMIC(ColsAtCompileTime, RowsAtCompileTime + ActualIndex)
    100                      : EIGEN_SIZE_MIN_PREFER_DYNAMIC(RowsAtCompileTime, ColsAtCompileTime - ActualIndex))
    101       };
    102       typedef Block<CoefficientsType,1, DiagonalSize> BuildType;
    103       typedef typename internal::conditional<Conjugate,
    104                  CwiseUnaryOp<internal::scalar_conjugate_op<Scalar>,BuildType >,
    105                  BuildType>::type Type;
    106     };
    107 
    108     /** \returns a vector expression of the \a N -th sub or super diagonal */
    109     template<int N> inline typename DiagonalIntReturnType<N>::Type diagonal()
    110     {
    111       return typename DiagonalIntReturnType<N>::BuildType(coeffs(), supers()-N, (std::max)(0,N), 1, diagonalLength(N));
    112     }
    113 
    114     /** \returns a vector expression of the \a N -th sub or super diagonal */
    115     template<int N> inline const typename DiagonalIntReturnType<N>::Type diagonal() const
    116     {
    117       return typename DiagonalIntReturnType<N>::BuildType(coeffs(), supers()-N, (std::max)(0,N), 1, diagonalLength(N));
    118     }
    119 
    120     /** \returns a vector expression of the \a i -th sub or super diagonal */
    121     inline Block<CoefficientsType,1,Dynamic> diagonal(Index i)
    122     {
    123       eigen_assert((i<0 && -i<=subs()) || (i>=0 && i<=supers()));
    124       return Block<CoefficientsType,1,Dynamic>(coeffs(), supers()-i, std::max<Index>(0,i), 1, diagonalLength(i));
    125     }
    126 
    127     /** \returns a vector expression of the \a i -th sub or super diagonal */
    128     inline const Block<const CoefficientsType,1,Dynamic> diagonal(Index i) const
    129     {
    130       eigen_assert((i<0 && -i<=subs()) || (i>=0 && i<=supers()));
    131       return Block<const CoefficientsType,1,Dynamic>(coeffs(), supers()-i, std::max<Index>(0,i), 1, diagonalLength(i));
    132     }
    133 
    134     template<typename Dest> inline void evalTo(Dest& dst) const
    135     {
    136       dst.resize(rows(),cols());
    137       dst.setZero();
    138       dst.diagonal() = diagonal();
    139       for (Index i=1; i<=supers();++i)
    140         dst.diagonal(i) = diagonal(i);
    141       for (Index i=1; i<=subs();++i)
    142         dst.diagonal(-i) = diagonal(-i);
    143     }
    144 
    145     DenseMatrixType toDenseMatrix() const
    146     {
    147       DenseMatrixType res(rows(),cols());
    148       evalTo(res);
    149       return res;
    150     }
    151 
    152   protected:
    153 
    154     inline Index diagonalLength(Index i) const
    155     { return i<0 ? (std::min)(cols(),rows()+i) : (std::min)(rows(),cols()-i); }
    156 };
    157 
    158 /**
    159   * \class BandMatrix
    160   * \ingroup Core_Module
    161   *
    162   * \brief Represents a rectangular matrix with a banded storage
    163   *
    164   * \param _Scalar Numeric type, i.e. float, double, int
    165   * \param Rows Number of rows, or \b Dynamic
    166   * \param Cols Number of columns, or \b Dynamic
    167   * \param Supers Number of super diagonal
    168   * \param Subs Number of sub diagonal
    169   * \param _Options A combination of either \b #RowMajor or \b #ColMajor, and of \b #SelfAdjoint
    170   *                 The former controls \ref TopicStorageOrders "storage order", and defaults to
    171   *                 column-major. The latter controls whether the matrix represents a selfadjoint
    172   *                 matrix in which case either Supers of Subs have to be null.
    173   *
    174   * \sa class TridiagonalMatrix
    175   */
    176 
    177 template<typename _Scalar, int _Rows, int _Cols, int _Supers, int _Subs, int _Options>
    178 struct traits<BandMatrix<_Scalar,_Rows,_Cols,_Supers,_Subs,_Options> >
    179 {
    180   typedef _Scalar Scalar;
    181   typedef Dense StorageKind;
    182   typedef DenseIndex Index;
    183   enum {
    184     CoeffReadCost = NumTraits<Scalar>::ReadCost,
    185     RowsAtCompileTime = _Rows,
    186     ColsAtCompileTime = _Cols,
    187     MaxRowsAtCompileTime = _Rows,
    188     MaxColsAtCompileTime = _Cols,
    189     Flags = LvalueBit,
    190     Supers = _Supers,
    191     Subs = _Subs,
    192     Options = _Options,
    193     DataRowsAtCompileTime = ((Supers!=Dynamic) && (Subs!=Dynamic)) ? 1 + Supers + Subs : Dynamic
    194   };
    195   typedef Matrix<Scalar,DataRowsAtCompileTime,ColsAtCompileTime,Options&RowMajor?RowMajor:ColMajor> CoefficientsType;
    196 };
    197 
    198 template<typename _Scalar, int Rows, int Cols, int Supers, int Subs, int Options>
    199 class BandMatrix : public BandMatrixBase<BandMatrix<_Scalar,Rows,Cols,Supers,Subs,Options> >
    200 {
    201   public:
    202 
    203     typedef typename internal::traits<BandMatrix>::Scalar Scalar;
    204     typedef typename internal::traits<BandMatrix>::Index Index;
    205     typedef typename internal::traits<BandMatrix>::CoefficientsType CoefficientsType;
    206 
    207     inline BandMatrix(Index rows=Rows, Index cols=Cols, Index supers=Supers, Index subs=Subs)
    208       : m_coeffs(1+supers+subs,cols),
    209         m_rows(rows), m_supers(supers), m_subs(subs)
    210     {
    211     }
    212 
    213     /** \returns the number of columns */
    214     inline Index rows() const { return m_rows.value(); }
    215 
    216     /** \returns the number of rows */
    217     inline Index cols() const { return m_coeffs.cols(); }
    218 
    219     /** \returns the number of super diagonals */
    220     inline Index supers() const { return m_supers.value(); }
    221 
    222     /** \returns the number of sub diagonals */
    223     inline Index subs() const { return m_subs.value(); }
    224 
    225     inline const CoefficientsType& coeffs() const { return m_coeffs; }
    226     inline CoefficientsType& coeffs() { return m_coeffs; }
    227 
    228   protected:
    229 
    230     CoefficientsType m_coeffs;
    231     internal::variable_if_dynamic<Index, Rows>   m_rows;
    232     internal::variable_if_dynamic<Index, Supers> m_supers;
    233     internal::variable_if_dynamic<Index, Subs>   m_subs;
    234 };
    235 
    236 template<typename _CoefficientsType,int _Rows, int _Cols, int _Supers, int _Subs,int _Options>
    237 class BandMatrixWrapper;
    238 
    239 template<typename _CoefficientsType,int _Rows, int _Cols, int _Supers, int _Subs,int _Options>
    240 struct traits<BandMatrixWrapper<_CoefficientsType,_Rows,_Cols,_Supers,_Subs,_Options> >
    241 {
    242   typedef typename _CoefficientsType::Scalar Scalar;
    243   typedef typename _CoefficientsType::StorageKind StorageKind;
    244   typedef typename _CoefficientsType::Index Index;
    245   enum {
    246     CoeffReadCost = internal::traits<_CoefficientsType>::CoeffReadCost,
    247     RowsAtCompileTime = _Rows,
    248     ColsAtCompileTime = _Cols,
    249     MaxRowsAtCompileTime = _Rows,
    250     MaxColsAtCompileTime = _Cols,
    251     Flags = LvalueBit,
    252     Supers = _Supers,
    253     Subs = _Subs,
    254     Options = _Options,
    255     DataRowsAtCompileTime = ((Supers!=Dynamic) && (Subs!=Dynamic)) ? 1 + Supers + Subs : Dynamic
    256   };
    257   typedef _CoefficientsType CoefficientsType;
    258 };
    259 
    260 template<typename _CoefficientsType,int _Rows, int _Cols, int _Supers, int _Subs,int _Options>
    261 class BandMatrixWrapper : public BandMatrixBase<BandMatrixWrapper<_CoefficientsType,_Rows,_Cols,_Supers,_Subs,_Options> >
    262 {
    263   public:
    264 
    265     typedef typename internal::traits<BandMatrixWrapper>::Scalar Scalar;
    266     typedef typename internal::traits<BandMatrixWrapper>::CoefficientsType CoefficientsType;
    267     typedef typename internal::traits<BandMatrixWrapper>::Index Index;
    268 
    269     inline BandMatrixWrapper(const CoefficientsType& coeffs, Index rows=_Rows, Index cols=_Cols, Index supers=_Supers, Index subs=_Subs)
    270       : m_coeffs(coeffs),
    271         m_rows(rows), m_supers(supers), m_subs(subs)
    272     {
    273       EIGEN_UNUSED_VARIABLE(cols);
    274       //internal::assert(coeffs.cols()==cols() && (supers()+subs()+1)==coeffs.rows());
    275     }
    276 
    277     /** \returns the number of columns */
    278     inline Index rows() const { return m_rows.value(); }
    279 
    280     /** \returns the number of rows */
    281     inline Index cols() const { return m_coeffs.cols(); }
    282 
    283     /** \returns the number of super diagonals */
    284     inline Index supers() const { return m_supers.value(); }
    285 
    286     /** \returns the number of sub diagonals */
    287     inline Index subs() const { return m_subs.value(); }
    288 
    289     inline const CoefficientsType& coeffs() const { return m_coeffs; }
    290 
    291   protected:
    292 
    293     const CoefficientsType& m_coeffs;
    294     internal::variable_if_dynamic<Index, _Rows>   m_rows;
    295     internal::variable_if_dynamic<Index, _Supers> m_supers;
    296     internal::variable_if_dynamic<Index, _Subs>   m_subs;
    297 };
    298 
    299 /**
    300   * \class TridiagonalMatrix
    301   * \ingroup Core_Module
    302   *
    303   * \brief Represents a tridiagonal matrix with a compact banded storage
    304   *
    305   * \param _Scalar Numeric type, i.e. float, double, int
    306   * \param Size Number of rows and cols, or \b Dynamic
    307   * \param _Options Can be 0 or \b SelfAdjoint
    308   *
    309   * \sa class BandMatrix
    310   */
    311 template<typename Scalar, int Size, int Options>
    312 class TridiagonalMatrix : public BandMatrix<Scalar,Size,Size,Options&SelfAdjoint?0:1,1,Options|RowMajor>
    313 {
    314     typedef BandMatrix<Scalar,Size,Size,Options&SelfAdjoint?0:1,1,Options|RowMajor> Base;
    315     typedef typename Base::Index Index;
    316   public:
    317     TridiagonalMatrix(Index size = Size) : Base(size,size,Options&SelfAdjoint?0:1,1) {}
    318 
    319     inline typename Base::template DiagonalIntReturnType<1>::Type super()
    320     { return Base::template diagonal<1>(); }
    321     inline const typename Base::template DiagonalIntReturnType<1>::Type super() const
    322     { return Base::template diagonal<1>(); }
    323     inline typename Base::template DiagonalIntReturnType<-1>::Type sub()
    324     { return Base::template diagonal<-1>(); }
    325     inline const typename Base::template DiagonalIntReturnType<-1>::Type sub() const
    326     { return Base::template diagonal<-1>(); }
    327   protected:
    328 };
    329 
    330 } // end namespace internal
    331 
    332 } // end namespace Eigen
    333 
    334 #endif // EIGEN_BANDMATRIX_H
    335