DiagonalMatrix.h
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 // Copyright (C) 2007-2009 Benoit Jacob <jacob.benoit.1@gmail.com>
6 //
7 // This Source Code Form is subject to the terms of the Mozilla
8 // Public License v. 2.0. If a copy of the MPL was not distributed
9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
10 
11 #ifndef EIGEN_DIAGONALMATRIX_H
12 #define EIGEN_DIAGONALMATRIX_H
13 
14 namespace Eigen {
15 
16 #ifndef EIGEN_PARSED_BY_DOXYGEN
17 template<typename Derived>
18 class DiagonalBase : public EigenBase<Derived>
19 {
20  public:
21  typedef typename internal::traits<Derived>::DiagonalVectorType DiagonalVectorType;
22  typedef typename DiagonalVectorType::Scalar Scalar;
23  typedef typename DiagonalVectorType::RealScalar RealScalar;
24  typedef typename internal::traits<Derived>::StorageKind StorageKind;
25  typedef typename internal::traits<Derived>::Index Index;
26 
27  enum {
28  RowsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
29  ColsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
30  MaxRowsAtCompileTime = DiagonalVectorType::MaxSizeAtCompileTime,
31  MaxColsAtCompileTime = DiagonalVectorType::MaxSizeAtCompileTime,
32  IsVectorAtCompileTime = 0,
33  Flags = 0
34  };
35 
36  typedef Matrix<Scalar, RowsAtCompileTime, ColsAtCompileTime, 0, MaxRowsAtCompileTime, MaxColsAtCompileTime> DenseMatrixType;
37  typedef DenseMatrixType DenseType;
38  typedef DiagonalMatrix<Scalar,DiagonalVectorType::SizeAtCompileTime,DiagonalVectorType::MaxSizeAtCompileTime> PlainObject;
39 
40  inline const Derived& derived() const { return *static_cast<const Derived*>(this); }
41  inline Derived& derived() { return *static_cast<Derived*>(this); }
42 
43  DenseMatrixType toDenseMatrix() const { return derived(); }
44  template<typename DenseDerived>
45  void evalTo(MatrixBase<DenseDerived> &other) const;
46  template<typename DenseDerived>
47  void addTo(MatrixBase<DenseDerived> &other) const
48  { other.diagonal() += diagonal(); }
49  template<typename DenseDerived>
50  void subTo(MatrixBase<DenseDerived> &other) const
51  { other.diagonal() -= diagonal(); }
52 
53  inline const DiagonalVectorType& diagonal() const { return derived().diagonal(); }
54  inline DiagonalVectorType& diagonal() { return derived().diagonal(); }
55 
56  inline Index rows() const { return diagonal().size(); }
57  inline Index cols() const { return diagonal().size(); }
58 
59  template<typename MatrixDerived>
60  const DiagonalProduct<MatrixDerived, Derived, OnTheLeft>
61  operator*(const MatrixBase<MatrixDerived> &matrix) const;
62 
63  inline const DiagonalWrapper<const CwiseUnaryOp<internal::scalar_inverse_op<Scalar>, const DiagonalVectorType> >
64  inverse() const
65  {
66  return diagonal().cwiseInverse();
67  }
68 
69  inline const DiagonalWrapper<const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DiagonalVectorType> >
70  operator*(const Scalar& scalar) const
71  {
72  return diagonal() * scalar;
73  }
74  friend inline const DiagonalWrapper<const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DiagonalVectorType> >
75  operator*(const Scalar& scalar, const DiagonalBase& other)
76  {
77  return other.diagonal() * scalar;
78  }
79 
80  #ifdef EIGEN2_SUPPORT
81  template<typename OtherDerived>
82  bool isApprox(const DiagonalBase<OtherDerived>& other, typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision()) const
83  {
84  return diagonal().isApprox(other.diagonal(), precision);
85  }
86  template<typename OtherDerived>
87  bool isApprox(const MatrixBase<OtherDerived>& other, typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision()) const
88  {
89  return toDenseMatrix().isApprox(other, precision);
90  }
91  #endif
92 };
93 
94 template<typename Derived>
95 template<typename DenseDerived>
96 void DiagonalBase<Derived>::evalTo(MatrixBase<DenseDerived> &other) const
97 {
98  other.setZero();
99  other.diagonal() = diagonal();
100 }
101 #endif
102 
116 namespace internal {
117 template<typename _Scalar, int SizeAtCompileTime, int MaxSizeAtCompileTime>
118 struct traits<DiagonalMatrix<_Scalar,SizeAtCompileTime,MaxSizeAtCompileTime> >
119  : traits<Matrix<_Scalar,SizeAtCompileTime,SizeAtCompileTime,0,MaxSizeAtCompileTime,MaxSizeAtCompileTime> >
120 {
121  typedef Matrix<_Scalar,SizeAtCompileTime,1,0,MaxSizeAtCompileTime,1> DiagonalVectorType;
122  typedef Dense StorageKind;
123  typedef DenseIndex Index;
124  enum {
125  Flags = LvalueBit
126  };
127 };
128 }
129 template<typename _Scalar, int SizeAtCompileTime, int MaxSizeAtCompileTime>
131  : public DiagonalBase<DiagonalMatrix<_Scalar,SizeAtCompileTime,MaxSizeAtCompileTime> >
132 {
133  public:
134  #ifndef EIGEN_PARSED_BY_DOXYGEN
135  typedef typename internal::traits<DiagonalMatrix>::DiagonalVectorType DiagonalVectorType;
136  typedef const DiagonalMatrix& Nested;
137  typedef _Scalar Scalar;
138  typedef typename internal::traits<DiagonalMatrix>::StorageKind StorageKind;
139  typedef typename internal::traits<DiagonalMatrix>::Index Index;
140  #endif
141 
142  protected:
143 
144  DiagonalVectorType m_diagonal;
145 
146  public:
147 
149  inline const DiagonalVectorType& diagonal() const { return m_diagonal; }
151  inline DiagonalVectorType& diagonal() { return m_diagonal; }
152 
154  inline DiagonalMatrix() {}
155 
157  inline DiagonalMatrix(Index dim) : m_diagonal(dim) {}
158 
160  inline DiagonalMatrix(const Scalar& x, const Scalar& y) : m_diagonal(x,y) {}
161 
163  inline DiagonalMatrix(const Scalar& x, const Scalar& y, const Scalar& z) : m_diagonal(x,y,z) {}
164 
166  template<typename OtherDerived>
167  inline DiagonalMatrix(const DiagonalBase<OtherDerived>& other) : m_diagonal(other.diagonal()) {}
168 
169  #ifndef EIGEN_PARSED_BY_DOXYGEN
170 
171  inline DiagonalMatrix(const DiagonalMatrix& other) : m_diagonal(other.diagonal()) {}
172  #endif
173 
175  template<typename OtherDerived>
176  explicit inline DiagonalMatrix(const MatrixBase<OtherDerived>& other) : m_diagonal(other)
177  {}
178 
180  template<typename OtherDerived>
181  DiagonalMatrix& operator=(const DiagonalBase<OtherDerived>& other)
182  {
183  m_diagonal = other.diagonal();
184  return *this;
185  }
186 
187  #ifndef EIGEN_PARSED_BY_DOXYGEN
188 
192  {
193  m_diagonal = other.diagonal();
194  return *this;
195  }
196  #endif
197 
199  inline void resize(Index size) { m_diagonal.resize(size); }
201  inline void setZero() { m_diagonal.setZero(); }
203  inline void setZero(Index size) { m_diagonal.setZero(size); }
205  inline void setIdentity() { m_diagonal.setOnes(); }
207  inline void setIdentity(Index size) { m_diagonal.setOnes(size); }
208 };
209 
224 namespace internal {
225 template<typename _DiagonalVectorType>
226 struct traits<DiagonalWrapper<_DiagonalVectorType> >
227 {
228  typedef _DiagonalVectorType DiagonalVectorType;
229  typedef typename DiagonalVectorType::Scalar Scalar;
230  typedef typename DiagonalVectorType::Index Index;
231  typedef typename DiagonalVectorType::StorageKind StorageKind;
232  enum {
233  RowsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
234  ColsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
235  MaxRowsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
236  MaxColsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
237  Flags = traits<DiagonalVectorType>::Flags & LvalueBit
238  };
239 };
240 }
241 
242 template<typename _DiagonalVectorType>
244  : public DiagonalBase<DiagonalWrapper<_DiagonalVectorType> >, internal::no_assignment_operator
245 {
246  public:
247  #ifndef EIGEN_PARSED_BY_DOXYGEN
248  typedef _DiagonalVectorType DiagonalVectorType;
249  typedef DiagonalWrapper Nested;
250  #endif
251 
253  inline DiagonalWrapper(DiagonalVectorType& diagonal) : m_diagonal(diagonal) {}
254 
256  const DiagonalVectorType& diagonal() const { return m_diagonal; }
257 
258  protected:
259  typename DiagonalVectorType::Nested m_diagonal;
260 };
261 
271 template<typename Derived>
272 inline const DiagonalWrapper<const Derived>
274 {
275  return derived();
276 }
277 
286 template<typename Derived>
287 bool MatrixBase<Derived>::isDiagonal(RealScalar prec) const
288 {
289  if(cols() != rows()) return false;
290  RealScalar maxAbsOnDiagonal = static_cast<RealScalar>(-1);
291  for(Index j = 0; j < cols(); ++j)
292  {
293  RealScalar absOnDiagonal = internal::abs(coeff(j,j));
294  if(absOnDiagonal > maxAbsOnDiagonal) maxAbsOnDiagonal = absOnDiagonal;
295  }
296  for(Index j = 0; j < cols(); ++j)
297  for(Index i = 0; i < j; ++i)
298  {
299  if(!internal::isMuchSmallerThan(coeff(i, j), maxAbsOnDiagonal, prec)) return false;
300  if(!internal::isMuchSmallerThan(coeff(j, i), maxAbsOnDiagonal, prec)) return false;
301  }
302  return true;
303 }
304 
305 } // end namespace Eigen
306 
307 #endif // EIGEN_DIAGONALMATRIX_H