Michael Ernst Peter / Eigen

src/Eigenvalues/MatrixBaseEigenvalues.h@1:3b8049da21b8, 2019-09-24 (annotated)

Committer:
jsoh91
Date:
Tue Sep 24 00:18:23 2019 +0000
Revision:
1:3b8049da21b8
Parent:
0:13a5d365ba16 
ignore and revise some of error parts

Who changed what in which revision?

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