00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 #ifndef EIGEN_DIAGONALMATRIX_H
00027 #define EIGEN_DIAGONALMATRIX_H
00028
00029 namespace Eigen {
00030
00031 #ifndef EIGEN_PARSED_BY_DOXYGEN
00032 template<typename Derived>
00033 class DiagonalBase : public EigenBase<Derived>
00034 {
00035 public:
00036 typedef typename internal::traits<Derived>::DiagonalVectorType DiagonalVectorType;
00037 typedef typename DiagonalVectorType::Scalar Scalar;
00038 typedef typename internal::traits<Derived>::StorageKind StorageKind;
00039 typedef typename internal::traits<Derived>::Index Index;
00040
00041 enum {
00042 RowsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
00043 ColsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
00044 MaxRowsAtCompileTime = DiagonalVectorType::MaxSizeAtCompileTime,
00045 MaxColsAtCompileTime = DiagonalVectorType::MaxSizeAtCompileTime,
00046 IsVectorAtCompileTime = 0,
00047 Flags = 0
00048 };
00049
00050 typedef Matrix<Scalar, RowsAtCompileTime, ColsAtCompileTime, 0, MaxRowsAtCompileTime, MaxColsAtCompileTime> DenseMatrixType;
00051 typedef DenseMatrixType DenseType;
00052 typedef DiagonalMatrix<Scalar,DiagonalVectorType::SizeAtCompileTime,DiagonalVectorType::MaxSizeAtCompileTime> PlainObject;
00053
00054 inline const Derived& derived() const { return *static_cast<const Derived*>(this); }
00055 inline Derived& derived() { return *static_cast<Derived*>(this); }
00056
00057 DenseMatrixType toDenseMatrix() const { return derived(); }
00058 template<typename DenseDerived>
00059 void evalTo(MatrixBase<DenseDerived> &other) const;
00060 template<typename DenseDerived>
00061 void addTo(MatrixBase<DenseDerived> &other) const
00062 { other.diagonal() += diagonal(); }
00063 template<typename DenseDerived>
00064 void subTo(MatrixBase<DenseDerived> &other) const
00065 { other.diagonal() -= diagonal(); }
00066
00067 inline const DiagonalVectorType& diagonal() const { return derived().diagonal(); }
00068 inline DiagonalVectorType& diagonal() { return derived().diagonal(); }
00069
00070 inline Index rows() const { return diagonal().size(); }
00071 inline Index cols() const { return diagonal().size(); }
00072
00073 template<typename MatrixDerived>
00074 const DiagonalProduct<MatrixDerived, Derived, OnTheLeft>
00075 operator*(const MatrixBase<MatrixDerived> &matrix) const;
00076
00077 inline const DiagonalWrapper<const CwiseUnaryOp<internal::scalar_inverse_op<Scalar>, const DiagonalVectorType> >
00078 inverse() const
00079 {
00080 return diagonal().cwiseInverse();
00081 }
00082
00083 #ifdef EIGEN2_SUPPORT
00084 template<typename OtherDerived>
00085 bool isApprox(const DiagonalBase<OtherDerived>& other, typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision()) const
00086 {
00087 return diagonal().isApprox(other.diagonal(), precision);
00088 }
00089 template<typename OtherDerived>
00090 bool isApprox(const MatrixBase<OtherDerived>& other, typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision()) const
00091 {
00092 return toDenseMatrix().isApprox(other, precision);
00093 }
00094 #endif
00095 };
00096
00097 template<typename Derived>
00098 template<typename DenseDerived>
00099 void DiagonalBase<Derived>::evalTo(MatrixBase<DenseDerived> &other) const
00100 {
00101 other.setZero();
00102 other.diagonal() = diagonal();
00103 }
00104 #endif
00105
00119 namespace internal {
00120 template<typename _Scalar, int SizeAtCompileTime, int MaxSizeAtCompileTime>
00121 struct traits<DiagonalMatrix<_Scalar,SizeAtCompileTime,MaxSizeAtCompileTime> >
00122 : traits<Matrix<_Scalar,SizeAtCompileTime,SizeAtCompileTime,0,MaxSizeAtCompileTime,MaxSizeAtCompileTime> >
00123 {
00124 typedef Matrix<_Scalar,SizeAtCompileTime,1,0,MaxSizeAtCompileTime,1> DiagonalVectorType;
00125 typedef Dense StorageKind;
00126 typedef DenseIndex Index;
00127 enum {
00128 Flags = LvalueBit
00129 };
00130 };
00131 }
00132 template<typename _Scalar, int SizeAtCompileTime, int MaxSizeAtCompileTime>
00133 class DiagonalMatrix
00134 : public DiagonalBase<DiagonalMatrix<_Scalar,SizeAtCompileTime,MaxSizeAtCompileTime> >
00135 {
00136 public:
00137 #ifndef EIGEN_PARSED_BY_DOXYGEN
00138 typedef typename internal::traits<DiagonalMatrix>::DiagonalVectorType DiagonalVectorType;
00139 typedef const DiagonalMatrix& Nested;
00140 typedef _Scalar Scalar;
00141 typedef typename internal::traits<DiagonalMatrix>::StorageKind StorageKind;
00142 typedef typename internal::traits<DiagonalMatrix>::Index Index;
00143 #endif
00144
00145 protected:
00146
00147 DiagonalVectorType m_diagonal;
00148
00149 public:
00150
00152 inline const DiagonalVectorType& diagonal() const { return m_diagonal; }
00154 inline DiagonalVectorType& diagonal() { return m_diagonal; }
00155
00157 inline DiagonalMatrix() {}
00158
00160 inline DiagonalMatrix(Index dim) : m_diagonal(dim) {}
00161
00163 inline DiagonalMatrix(const Scalar& x, const Scalar& y) : m_diagonal(x,y) {}
00164
00166 inline DiagonalMatrix(const Scalar& x, const Scalar& y, const Scalar& z) : m_diagonal(x,y,z) {}
00167
00169 template<typename OtherDerived>
00170 inline DiagonalMatrix(const DiagonalBase<OtherDerived>& other) : m_diagonal(other.diagonal()) {}
00171
00172 #ifndef EIGEN_PARSED_BY_DOXYGEN
00173
00174 inline DiagonalMatrix(const DiagonalMatrix& other) : m_diagonal(other.diagonal()) {}
00175 #endif
00176
00178 template<typename OtherDerived>
00179 explicit inline DiagonalMatrix(const MatrixBase<OtherDerived>& other) : m_diagonal(other)
00180 {}
00181
00183 template<typename OtherDerived>
00184 DiagonalMatrix& operator=(const DiagonalBase<OtherDerived>& other)
00185 {
00186 m_diagonal = other.diagonal();
00187 return *this;
00188 }
00189
00190 #ifndef EIGEN_PARSED_BY_DOXYGEN
00191
00194 DiagonalMatrix& operator=(const DiagonalMatrix& other)
00195 {
00196 m_diagonal = other.diagonal();
00197 return *this;
00198 }
00199 #endif
00200
00202 inline void resize(Index size) { m_diagonal.resize(size); }
00204 inline void setZero() { m_diagonal.setZero(); }
00206 inline void setZero(Index size) { m_diagonal.setZero(size); }
00208 inline void setIdentity() { m_diagonal.setOnes(); }
00210 inline void setIdentity(Index size) { m_diagonal.setOnes(size); }
00211 };
00212
00227 namespace internal {
00228 template<typename _DiagonalVectorType>
00229 struct traits<DiagonalWrapper<_DiagonalVectorType> >
00230 {
00231 typedef _DiagonalVectorType DiagonalVectorType;
00232 typedef typename DiagonalVectorType::Scalar Scalar;
00233 typedef typename DiagonalVectorType::Index Index;
00234 typedef typename DiagonalVectorType::StorageKind StorageKind;
00235 enum {
00236 RowsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
00237 ColsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
00238 MaxRowsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
00239 MaxColsAtCompileTime = DiagonalVectorType::SizeAtCompileTime,
00240 Flags = traits<DiagonalVectorType>::Flags & LvalueBit
00241 };
00242 };
00243 }
00244
00245 template<typename _DiagonalVectorType>
00246 class DiagonalWrapper
00247 : public DiagonalBase<DiagonalWrapper<_DiagonalVectorType> >, internal::no_assignment_operator
00248 {
00249 public:
00250 #ifndef EIGEN_PARSED_BY_DOXYGEN
00251 typedef _DiagonalVectorType DiagonalVectorType;
00252 typedef DiagonalWrapper Nested;
00253 #endif
00254
00256 inline DiagonalWrapper(DiagonalVectorType& diagonal) : m_diagonal(diagonal) {}
00257
00259 const DiagonalVectorType& diagonal() const { return m_diagonal; }
00260
00261 protected:
00262 typename DiagonalVectorType::Nested m_diagonal;
00263 };
00264
00274 template<typename Derived>
00275 inline const DiagonalWrapper<const Derived>
00276 MatrixBase<Derived>::asDiagonal() const
00277 {
00278 return derived();
00279 }
00280
00289 template<typename Derived>
00290 bool MatrixBase<Derived>::isDiagonal(RealScalar prec) const
00291 {
00292 if(cols() != rows()) return false;
00293 RealScalar maxAbsOnDiagonal = static_cast<RealScalar>(-1);
00294 for(Index j = 0; j < cols(); ++j)
00295 {
00296 RealScalar absOnDiagonal = internal::abs(coeff(j,j));
00297 if(absOnDiagonal > maxAbsOnDiagonal) maxAbsOnDiagonal = absOnDiagonal;
00298 }
00299 for(Index j = 0; j < cols(); ++j)
00300 for(Index i = 0; i < j; ++i)
00301 {
00302 if(!internal::isMuchSmallerThan(coeff(i, j), maxAbsOnDiagonal, prec)) return false;
00303 if(!internal::isMuchSmallerThan(coeff(j, i), maxAbsOnDiagonal, prec)) return false;
00304 }
00305 return true;
00306 }
00307
00308 }
00309
00310 #endif // EIGEN_DIAGONALMATRIX_H