1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud (at) inria.fr> 5 // Copyright (C) 2010 Jitse Niesen <jitse (at) maths.leeds.ac.uk> 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_MATRIXBASEEIGENVALUES_H 12 #define EIGEN_MATRIXBASEEIGENVALUES_H 13 14 namespace Eigen { 15 16 namespace internal { 17 18 template<typename Derived, bool IsComplex> 19 struct eigenvalues_selector 20 { 21 // this is the implementation for the case IsComplex = true 22 static inline typename MatrixBase<Derived>::EigenvaluesReturnType const 23 run(const MatrixBase<Derived>& m) 24 { 25 typedef typename Derived::PlainObject PlainObject; 26 PlainObject m_eval(m); 27 return ComplexEigenSolver<PlainObject>(m_eval, false).eigenvalues(); 28 } 29 }; 30 31 template<typename Derived> 32 struct eigenvalues_selector<Derived, false> 33 { 34 static inline typename MatrixBase<Derived>::EigenvaluesReturnType const 35 run(const MatrixBase<Derived>& m) 36 { 37 typedef typename Derived::PlainObject PlainObject; 38 PlainObject m_eval(m); 39 return EigenSolver<PlainObject>(m_eval, false).eigenvalues(); 40 } 41 }; 42 43 } // end namespace internal 44 45 /** \brief Computes the eigenvalues of a matrix 46 * \returns Column vector containing the eigenvalues. 47 * 48 * \eigenvalues_module 49 * This function computes the eigenvalues with the help of the EigenSolver 50 * class (for real matrices) or the ComplexEigenSolver class (for complex 51 * matrices). 52 * 53 * The eigenvalues are repeated according to their algebraic multiplicity, 54 * so there are as many eigenvalues as rows in the matrix. 55 * 56 * The SelfAdjointView class provides a better algorithm for selfadjoint 57 * matrices. 58 * 59 * Example: \include MatrixBase_eigenvalues.cpp 60 * Output: \verbinclude MatrixBase_eigenvalues.out 61 * 62 * \sa EigenSolver::eigenvalues(), ComplexEigenSolver::eigenvalues(), 63 * SelfAdjointView::eigenvalues() 64 */ 65 template<typename Derived> 66 inline typename MatrixBase<Derived>::EigenvaluesReturnType 67 MatrixBase<Derived>::eigenvalues() const 68 { 69 typedef typename internal::traits<Derived>::Scalar Scalar; 70 return internal::eigenvalues_selector<Derived, NumTraits<Scalar>::IsComplex>::run(derived()); 71 } 72 73 /** \brief Computes the eigenvalues of a matrix 74 * \returns Column vector containing the eigenvalues. 75 * 76 * \eigenvalues_module 77 * This function computes the eigenvalues with the help of the 78 * SelfAdjointEigenSolver class. The eigenvalues are repeated according to 79 * their algebraic multiplicity, so there are as many eigenvalues as rows in 80 * the matrix. 81 * 82 * Example: \include SelfAdjointView_eigenvalues.cpp 83 * Output: \verbinclude SelfAdjointView_eigenvalues.out 84 * 85 * \sa SelfAdjointEigenSolver::eigenvalues(), MatrixBase::eigenvalues() 86 */ 87 template<typename MatrixType, unsigned int UpLo> 88 inline typename SelfAdjointView<MatrixType, UpLo>::EigenvaluesReturnType 89 SelfAdjointView<MatrixType, UpLo>::eigenvalues() const 90 { 91 typedef typename SelfAdjointView<MatrixType, UpLo>::PlainObject PlainObject; 92 PlainObject thisAsMatrix(*this); 93 return SelfAdjointEigenSolver<PlainObject>(thisAsMatrix, false).eigenvalues(); 94 } 95 96 97 98 /** \brief Computes the L2 operator norm 99 * \returns Operator norm of the matrix. 100 * 101 * \eigenvalues_module 102 * This function computes the L2 operator norm of a matrix, which is also 103 * known as the spectral norm. The norm of a matrix \f$ A \f$ is defined to be 104 * \f[ \|A\|_2 = \max_x \frac{\|Ax\|_2}{\|x\|_2} \f] 105 * where the maximum is over all vectors and the norm on the right is the 106 * Euclidean vector norm. The norm equals the largest singular value, which is 107 * the square root of the largest eigenvalue of the positive semi-definite 108 * matrix \f$ A^*A \f$. 109 * 110 * The current implementation uses the eigenvalues of \f$ A^*A \f$, as computed 111 * by SelfAdjointView::eigenvalues(), to compute the operator norm of a 112 * matrix. The SelfAdjointView class provides a better algorithm for 113 * selfadjoint matrices. 114 * 115 * Example: \include MatrixBase_operatorNorm.cpp 116 * Output: \verbinclude MatrixBase_operatorNorm.out 117 * 118 * \sa SelfAdjointView::eigenvalues(), SelfAdjointView::operatorNorm() 119 */ 120 template<typename Derived> 121 inline typename MatrixBase<Derived>::RealScalar 122 MatrixBase<Derived>::operatorNorm() const 123 { 124 using std::sqrt; 125 typename Derived::PlainObject m_eval(derived()); 126 // FIXME if it is really guaranteed that the eigenvalues are already sorted, 127 // then we don't need to compute a maxCoeff() here, comparing the 1st and last ones is enough. 128 return sqrt((m_eval*m_eval.adjoint()) 129 .eval() 130 .template selfadjointView<Lower>() 131 .eigenvalues() 132 .maxCoeff() 133 ); 134 } 135 136 /** \brief Computes the L2 operator norm 137 * \returns Operator norm of the matrix. 138 * 139 * \eigenvalues_module 140 * This function computes the L2 operator norm of a self-adjoint matrix. For a 141 * self-adjoint matrix, the operator norm is the largest eigenvalue. 142 * 143 * The current implementation uses the eigenvalues of the matrix, as computed 144 * by eigenvalues(), to compute the operator norm of the matrix. 145 * 146 * Example: \include SelfAdjointView_operatorNorm.cpp 147 * Output: \verbinclude SelfAdjointView_operatorNorm.out 148 * 149 * \sa eigenvalues(), MatrixBase::operatorNorm() 150 */ 151 template<typename MatrixType, unsigned int UpLo> 152 inline typename SelfAdjointView<MatrixType, UpLo>::RealScalar 153 SelfAdjointView<MatrixType, UpLo>::operatorNorm() const 154 { 155 return eigenvalues().cwiseAbs().maxCoeff(); 156 } 157 158 } // end namespace Eigen 159 160 #endif 161
