BandMatrix.h revision c981c48f5bc9aefeffc0bcb0cc3934c2fae179dd 
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@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
13namespace Eigen {
14
15namespace internal {
16
17template<typename Derived>
18class 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
177template<typename _Scalar, int _Rows, int _Cols, int _Supers, int _Subs, int _Options>
178struct 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
198template<typename _Scalar, int Rows, int Cols, int Supers, int Subs, int Options>
199class 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
236template<typename _CoefficientsType,int _Rows, int _Cols, int _Supers, int _Subs,int _Options>
237class BandMatrixWrapper;
238
239template<typename _CoefficientsType,int _Rows, int _Cols, int _Supers, int _Subs,int _Options>
240struct 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
260template<typename _CoefficientsType,int _Rows, int _Cols, int _Supers, int _Subs,int _Options>
261class 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  */
311template<typename Scalar, int Size, int Options>
312class 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