Eigne Matrix Class Library
Dependents: Eigen_test Odometry_test AttitudeEstimation_usingTicker MPU9250_Quaternion_Binary_Serial ... more
SelfAdjointView.h
00001 // This file is part of Eigen, a lightweight C++ template library 00002 // for linear algebra. 00003 // 00004 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr> 00005 // 00006 // This Source Code Form is subject to the terms of the Mozilla 00007 // Public License v. 2.0. If a copy of the MPL was not distributed 00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00009 00010 #ifndef EIGEN_SELFADJOINTMATRIX_H 00011 #define EIGEN_SELFADJOINTMATRIX_H 00012 00013 namespace Eigen { 00014 00015 /** \class SelfAdjointView 00016 * \ingroup Core_Module 00017 * 00018 * 00019 * \brief Expression of a selfadjoint matrix from a triangular part of a dense matrix 00020 * 00021 * \param MatrixType the type of the dense matrix storing the coefficients 00022 * \param TriangularPart can be either \c #Lower or \c #Upper 00023 * 00024 * This class is an expression of a sefladjoint matrix from a triangular part of a matrix 00025 * with given dense storage of the coefficients. It is the return type of MatrixBase::selfadjointView() 00026 * and most of the time this is the only way that it is used. 00027 * 00028 * \sa class TriangularBase, MatrixBase::selfadjointView() 00029 */ 00030 00031 namespace internal { 00032 template<typename MatrixType, unsigned int UpLo> 00033 struct traits<SelfAdjointView<MatrixType, UpLo> > : traits<MatrixType> 00034 { 00035 typedef typename nested<MatrixType>::type MatrixTypeNested; 00036 typedef typename remove_all<MatrixTypeNested>::type MatrixTypeNestedCleaned; 00037 typedef MatrixType ExpressionType; 00038 typedef typename MatrixType::PlainObject DenseMatrixType; 00039 enum { 00040 Mode = UpLo | SelfAdjoint, 00041 Flags = MatrixTypeNestedCleaned::Flags & (HereditaryBits) 00042 & (~(PacketAccessBit | DirectAccessBit | LinearAccessBit)), // FIXME these flags should be preserved 00043 CoeffReadCost = MatrixTypeNestedCleaned::CoeffReadCost 00044 }; 00045 }; 00046 } 00047 00048 template <typename Lhs, int LhsMode, bool LhsIsVector, 00049 typename Rhs, int RhsMode, bool RhsIsVector> 00050 struct SelfadjointProductMatrix; 00051 00052 // FIXME could also be called SelfAdjointWrapper to be consistent with DiagonalWrapper ?? 00053 template<typename MatrixType, unsigned int UpLo> class SelfAdjointView 00054 : public TriangularBase<SelfAdjointView<MatrixType, UpLo> > 00055 { 00056 public: 00057 00058 typedef TriangularBase<SelfAdjointView> Base; 00059 typedef typename internal::traits<SelfAdjointView>::MatrixTypeNested MatrixTypeNested; 00060 typedef typename internal::traits<SelfAdjointView>::MatrixTypeNestedCleaned MatrixTypeNestedCleaned; 00061 00062 /** \brief The type of coefficients in this matrix */ 00063 typedef typename internal::traits<SelfAdjointView>::Scalar Scalar; 00064 00065 typedef typename MatrixType::Index Index; 00066 00067 enum { 00068 Mode = internal::traits<SelfAdjointView>::Mode 00069 }; 00070 typedef typename MatrixType::PlainObject PlainObject; 00071 00072 inline SelfAdjointView(MatrixType& matrix) : m_matrix(matrix) 00073 {} 00074 00075 inline Index rows () const { return m_matrix.rows(); } 00076 inline Index cols () const { return m_matrix.cols(); } 00077 inline Index outerStride() const { return m_matrix.outerStride(); } 00078 inline Index innerStride() const { return m_matrix.innerStride(); } 00079 00080 /** \sa MatrixBase::coeff() 00081 * \warning the coordinates must fit into the referenced triangular part 00082 */ 00083 inline Scalar coeff(Index row, Index col) const 00084 { 00085 Base::check_coordinates_internal(row, col); 00086 return m_matrix.coeff(row, col); 00087 } 00088 00089 /** \sa MatrixBase::coeffRef() 00090 * \warning the coordinates must fit into the referenced triangular part 00091 */ 00092 inline Scalar& coeffRef(Index row, Index col) 00093 { 00094 Base::check_coordinates_internal(row, col); 00095 return m_matrix.const_cast_derived().coeffRef(row, col); 00096 } 00097 00098 /** \internal */ 00099 const MatrixTypeNestedCleaned& _expression() const { return m_matrix; } 00100 00101 const MatrixTypeNestedCleaned& nestedExpression() const { return m_matrix; } 00102 MatrixTypeNestedCleaned& nestedExpression() { return *const_cast<MatrixTypeNestedCleaned*>(&m_matrix); } 00103 00104 /** Efficient self-adjoint matrix times vector/matrix product */ 00105 template<typename OtherDerived> 00106 SelfadjointProductMatrix<MatrixType,Mode,false,OtherDerived,0,OtherDerived::IsVectorAtCompileTime> 00107 operator* (const MatrixBase<OtherDerived>& rhs) const 00108 { 00109 return SelfadjointProductMatrix 00110 <MatrixType,Mode,false,OtherDerived,0,OtherDerived::IsVectorAtCompileTime> 00111 (m_matrix, rhs.derived()); 00112 } 00113 00114 /** Efficient vector/matrix times self-adjoint matrix product */ 00115 template<typename OtherDerived> friend 00116 SelfadjointProductMatrix<OtherDerived,0,OtherDerived::IsVectorAtCompileTime,MatrixType,Mode,false> 00117 operator* (const MatrixBase<OtherDerived>& lhs, const SelfAdjointView& rhs) 00118 { 00119 return SelfadjointProductMatrix 00120 <OtherDerived,0,OtherDerived::IsVectorAtCompileTime,MatrixType,Mode,false> 00121 (lhs.derived(),rhs.m_matrix); 00122 } 00123 00124 /** Perform a symmetric rank 2 update of the selfadjoint matrix \c *this: 00125 * \f$ this = this + \alpha u v^* + conj(\alpha) v u^* \f$ 00126 * \returns a reference to \c *this 00127 * 00128 * The vectors \a u and \c v \b must be column vectors, however they can be 00129 * a adjoint expression without any overhead. Only the meaningful triangular 00130 * part of the matrix is updated, the rest is left unchanged. 00131 * 00132 * \sa rankUpdate(const MatrixBase<DerivedU>&, Scalar) 00133 */ 00134 template<typename DerivedU, typename DerivedV> 00135 SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const MatrixBase<DerivedV>& v, const Scalar& alpha = Scalar(1)); 00136 00137 /** Perform a symmetric rank K update of the selfadjoint matrix \c *this: 00138 * \f$ this = this + \alpha ( u u^* ) \f$ where \a u is a vector or matrix. 00139 * 00140 * \returns a reference to \c *this 00141 * 00142 * Note that to perform \f$ this = this + \alpha ( u^* u ) \f$ you can simply 00143 * call this function with u.adjoint(). 00144 * 00145 * \sa rankUpdate(const MatrixBase<DerivedU>&, const MatrixBase<DerivedV>&, Scalar) 00146 */ 00147 template<typename DerivedU> 00148 SelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, const Scalar& alpha = Scalar(1)); 00149 00150 /////////// Cholesky module /////////// 00151 00152 const LLT<PlainObject, UpLo> llt() const; 00153 const LDLT<PlainObject, UpLo> ldlt() const; 00154 00155 /////////// Eigenvalue module /////////// 00156 00157 /** Real part of #Scalar */ 00158 typedef typename NumTraits<Scalar>::Real RealScalar; 00159 /** Return type of eigenvalues() */ 00160 typedef Matrix<RealScalar, internal::traits<MatrixType>::ColsAtCompileTime, 1> EigenvaluesReturnType; 00161 00162 EigenvaluesReturnType eigenvalues() const; 00163 RealScalar operatorNorm() const; 00164 00165 #ifdef EIGEN2_SUPPORT 00166 template<typename OtherDerived> 00167 SelfAdjointView& operator=(const MatrixBase<OtherDerived>& other) 00168 { 00169 enum { 00170 OtherPart = UpLo == Upper ? StrictlyLower : StrictlyUpper 00171 }; 00172 m_matrix.const_cast_derived().template triangularView<UpLo>() = other; 00173 m_matrix.const_cast_derived().template triangularView<OtherPart>() = other.adjoint (); 00174 return *this; 00175 } 00176 template<typename OtherMatrixType, unsigned int OtherMode> 00177 SelfAdjointView& operator=(const TriangularView<OtherMatrixType, OtherMode>& other) 00178 { 00179 enum { 00180 OtherPart = UpLo == Upper ? StrictlyLower : StrictlyUpper 00181 }; 00182 m_matrix.const_cast_derived().template triangularView<UpLo>() = other.toDenseMatrix(); 00183 m_matrix.const_cast_derived().template triangularView<OtherPart>() = other.toDenseMatrix().adjoint(); 00184 return *this; 00185 } 00186 #endif 00187 00188 protected: 00189 MatrixTypeNested m_matrix; 00190 }; 00191 00192 00193 // template<typename OtherDerived, typename MatrixType, unsigned int UpLo> 00194 // internal::selfadjoint_matrix_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> > 00195 // operator*(const MatrixBase<OtherDerived>& lhs, const SelfAdjointView<MatrixType,UpLo>& rhs) 00196 // { 00197 // return internal::matrix_selfadjoint_product_returntype<OtherDerived,SelfAdjointView<MatrixType,UpLo> >(lhs.derived(),rhs); 00198 // } 00199 00200 // selfadjoint to dense matrix 00201 00202 namespace internal { 00203 00204 template<typename Derived1, typename Derived2, int UnrollCount, bool ClearOpposite> 00205 struct triangular_assignment_selector<Derived1, Derived2, (SelfAdjoint|Upper), UnrollCount, ClearOpposite> 00206 { 00207 enum { 00208 col = (UnrollCount-1) / Derived1::RowsAtCompileTime, 00209 row = (UnrollCount-1) % Derived1::RowsAtCompileTime 00210 }; 00211 00212 static inline void run(Derived1 &dst, const Derived2 &src) 00213 { 00214 triangular_assignment_selector<Derived1, Derived2, (SelfAdjoint|Upper), UnrollCount-1, ClearOpposite>::run(dst, src); 00215 00216 if(row == col) 00217 dst.coeffRef(row, col) = numext::real(src.coeff(row, col)); 00218 else if(row < col) 00219 dst.coeffRef(col, row) = numext::conj(dst.coeffRef(row, col) = src.coeff(row, col)); 00220 } 00221 }; 00222 00223 template<typename Derived1, typename Derived2, bool ClearOpposite> 00224 struct triangular_assignment_selector<Derived1, Derived2, SelfAdjoint|Upper, 0, ClearOpposite> 00225 { 00226 static inline void run(Derived1 &, const Derived2 &) {} 00227 }; 00228 00229 template<typename Derived1, typename Derived2, int UnrollCount, bool ClearOpposite> 00230 struct triangular_assignment_selector<Derived1, Derived2, (SelfAdjoint|Lower), UnrollCount, ClearOpposite> 00231 { 00232 enum { 00233 col = (UnrollCount-1) / Derived1::RowsAtCompileTime, 00234 row = (UnrollCount-1) % Derived1::RowsAtCompileTime 00235 }; 00236 00237 static inline void run(Derived1 &dst, const Derived2 &src) 00238 { 00239 triangular_assignment_selector<Derived1, Derived2, (SelfAdjoint|Lower), UnrollCount-1, ClearOpposite>::run(dst, src); 00240 00241 if(row == col) 00242 dst.coeffRef(row, col) = numext::real(src.coeff(row, col)); 00243 else if(row > col) 00244 dst.coeffRef(col, row) = numext::conj(dst.coeffRef(row, col) = src.coeff(row, col)); 00245 } 00246 }; 00247 00248 template<typename Derived1, typename Derived2, bool ClearOpposite> 00249 struct triangular_assignment_selector<Derived1, Derived2, SelfAdjoint|Lower, 0, ClearOpposite> 00250 { 00251 static inline void run(Derived1 &, const Derived2 &) {} 00252 }; 00253 00254 template<typename Derived1, typename Derived2, bool ClearOpposite> 00255 struct triangular_assignment_selector<Derived1, Derived2, SelfAdjoint|Upper, Dynamic, ClearOpposite> 00256 { 00257 typedef typename Derived1::Index Index; 00258 static inline void run(Derived1 &dst, const Derived2 &src) 00259 { 00260 for(Index j = 0; j < dst.cols(); ++j) 00261 { 00262 for(Index i = 0; i < j; ++i) 00263 { 00264 dst.copyCoeff(i, j, src); 00265 dst.coeffRef(j,i) = numext::conj(dst.coeff(i,j)); 00266 } 00267 dst.copyCoeff(j, j, src); 00268 } 00269 } 00270 }; 00271 00272 template<typename Derived1, typename Derived2, bool ClearOpposite> 00273 struct triangular_assignment_selector<Derived1, Derived2, SelfAdjoint|Lower, Dynamic, ClearOpposite> 00274 { 00275 static inline void run(Derived1 &dst, const Derived2 &src) 00276 { 00277 typedef typename Derived1::Index Index; 00278 for(Index i = 0; i < dst.rows(); ++i) 00279 { 00280 for(Index j = 0; j < i; ++j) 00281 { 00282 dst.copyCoeff(i, j, src); 00283 dst.coeffRef(j,i) = numext::conj(dst.coeff(i,j)); 00284 } 00285 dst.copyCoeff(i, i, src); 00286 } 00287 } 00288 }; 00289 00290 } // end namespace internal 00291 00292 /*************************************************************************** 00293 * Implementation of MatrixBase methods 00294 ***************************************************************************/ 00295 00296 template<typename Derived> 00297 template<unsigned int UpLo> 00298 typename MatrixBase<Derived>::template ConstSelfAdjointViewReturnType<UpLo>::Type 00299 MatrixBase<Derived>::selfadjointView() const 00300 { 00301 return derived(); 00302 } 00303 00304 template<typename Derived> 00305 template<unsigned int UpLo> 00306 typename MatrixBase<Derived>::template SelfAdjointViewReturnType<UpLo>::Type 00307 MatrixBase<Derived>::selfadjointView() 00308 { 00309 return derived(); 00310 } 00311 00312 } // end namespace Eigen 00313 00314 #endif // EIGEN_SELFADJOINTMATRIX_H
Generated on Tue Jul 12 2022 17:47:00 by
![doxygen](doxygen.png)