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::StorageIndex StorageIndex;
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 * \tparam _Scalar Numeric type, i.e. float, double, int
165 * \tparam _Rows Number of rows, or \b Dynamic
166 * \tparam _Cols Number of columns, or \b Dynamic
167 * \tparam _Supers Number of super diagonal
168 * \tparam _Subs Number of sub diagonal
169 * \tparam _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 Eigen::Index StorageIndex;
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>::StorageIndex StorageIndex;
205 typedef typename internal::traits<BandMatrix>::CoefficientsType CoefficientsType;
206
207 explicit 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::StorageIndex StorageIndex;
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>::StorageIndex StorageIndex;
268
269 explicit 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 * \tparam Scalar Numeric type, i.e. float, double, int
306 * \tparam Size Number of rows and cols, or \b Dynamic
307 * \tparam 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::StorageIndex StorageIndex;
316 public:
317 explicit 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
331struct BandShape {};
332
333template<typename _Scalar, int _Rows, int _Cols, int _Supers, int _Subs, int _Options>
334struct evaluator_traits<BandMatrix<_Scalar,_Rows,_Cols,_Supers,_Subs,_Options> >
335 : public evaluator_traits_base<BandMatrix<_Scalar,_Rows,_Cols,_Supers,_Subs,_Options> >
336{
337 typedef BandShape Shape;
338};
339
340template<typename _CoefficientsType,int _Rows, int _Cols, int _Supers, int _Subs,int _Options>
341struct evaluator_traits<BandMatrixWrapper<_CoefficientsType,_Rows,_Cols,_Supers,_Subs,_Options> >
342 : public evaluator_traits_base<BandMatrixWrapper<_CoefficientsType,_Rows,_Cols,_Supers,_Subs,_Options> >
343{
344 typedef BandShape Shape;
345};
346
347template<> struct AssignmentKind<DenseShape,BandShape> { typedef EigenBase2EigenBase Kind; };
348
349} // end namespace internal
350
351} // end namespace Eigen
352
353#endif // EIGEN_BANDMATRIX_H
354