Eigen  3.2.5
SelfAdjointEigenSolver.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008-2010 Gael Guennebaud <gael.guennebaud@inria.fr>
5 // Copyright (C) 2010 Jitse Niesen <jitse@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_SELFADJOINTEIGENSOLVER_H
12 #define EIGEN_SELFADJOINTEIGENSOLVER_H
13 
14 #include "./Tridiagonalization.h"
15 
16 namespace Eigen {
17 
18 template<typename _MatrixType>
19 class GeneralizedSelfAdjointEigenSolver;
20 
21 namespace internal {
22 template<typename SolverType,int Size,bool IsComplex> struct direct_selfadjoint_eigenvalues;
23 }
24 
68 template<typename _MatrixType> class SelfAdjointEigenSolver
69 {
70  public:
71 
72  typedef _MatrixType MatrixType;
73  enum {
74  Size = MatrixType::RowsAtCompileTime,
75  ColsAtCompileTime = MatrixType::ColsAtCompileTime,
76  Options = MatrixType::Options,
77  MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
78  };
79 
81  typedef typename MatrixType::Scalar Scalar;
82  typedef typename MatrixType::Index Index;
83 
91 
92  friend struct internal::direct_selfadjoint_eigenvalues<SelfAdjointEigenSolver,Size,NumTraits<Scalar>::IsComplex>;
93 
99  typedef typename internal::plain_col_type<MatrixType, RealScalar>::type RealVectorType;
100  typedef Tridiagonalization<MatrixType> TridiagonalizationType;
101 
113  : m_eivec(),
114  m_eivalues(),
115  m_subdiag(),
116  m_isInitialized(false)
117  { }
118 
132  : m_eivec(size, size),
133  m_eivalues(size),
134  m_subdiag(size > 1 ? size - 1 : 1),
135  m_isInitialized(false)
136  {}
137 
153  SelfAdjointEigenSolver(const MatrixType& matrix, int options = ComputeEigenvectors)
154  : m_eivec(matrix.rows(), matrix.cols()),
155  m_eivalues(matrix.cols()),
156  m_subdiag(matrix.rows() > 1 ? matrix.rows() - 1 : 1),
157  m_isInitialized(false)
158  {
159  compute(matrix, options);
160  }
161 
192  SelfAdjointEigenSolver& compute(const MatrixType& matrix, int options = ComputeEigenvectors);
193 
208  SelfAdjointEigenSolver& computeDirect(const MatrixType& matrix, int options = ComputeEigenvectors);
209 
228  const MatrixType& eigenvectors() const
229  {
230  eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
231  eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
232  return m_eivec;
233  }
234 
251  {
252  eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
253  return m_eivalues;
254  }
255 
274  MatrixType operatorSqrt() const
275  {
276  eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
277  eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
278  return m_eivec * m_eivalues.cwiseSqrt().asDiagonal() * m_eivec.adjoint();
279  }
280 
299  MatrixType operatorInverseSqrt() const
300  {
301  eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
302  eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
303  return m_eivec * m_eivalues.cwiseInverse().cwiseSqrt().asDiagonal() * m_eivec.adjoint();
304  }
305 
311  {
312  eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
313  return m_info;
314  }
315 
321  static const int m_maxIterations = 30;
322 
323  #ifdef EIGEN2_SUPPORT
324  SelfAdjointEigenSolver(const MatrixType& matrix, bool computeEigenvectors)
325  : m_eivec(matrix.rows(), matrix.cols()),
326  m_eivalues(matrix.cols()),
327  m_subdiag(matrix.rows() > 1 ? matrix.rows() - 1 : 1),
328  m_isInitialized(false)
329  {
330  compute(matrix, computeEigenvectors);
331  }
332 
333  SelfAdjointEigenSolver(const MatrixType& matA, const MatrixType& matB, bool computeEigenvectors = true)
334  : m_eivec(matA.cols(), matA.cols()),
335  m_eivalues(matA.cols()),
336  m_subdiag(matA.cols() > 1 ? matA.cols() - 1 : 1),
337  m_isInitialized(false)
338  {
339  static_cast<GeneralizedSelfAdjointEigenSolver<MatrixType>*>(this)->compute(matA, matB, computeEigenvectors ? ComputeEigenvectors : EigenvaluesOnly);
340  }
341 
342  void compute(const MatrixType& matrix, bool computeEigenvectors)
343  {
344  compute(matrix, computeEigenvectors ? ComputeEigenvectors : EigenvaluesOnly);
345  }
346 
347  void compute(const MatrixType& matA, const MatrixType& matB, bool computeEigenvectors = true)
348  {
349  compute(matA, matB, computeEigenvectors ? ComputeEigenvectors : EigenvaluesOnly);
350  }
351  #endif // EIGEN2_SUPPORT
352 
353  protected:
354  static void check_template_parameters()
355  {
356  EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar);
357  }
358 
359  MatrixType m_eivec;
360  RealVectorType m_eivalues;
361  typename TridiagonalizationType::SubDiagonalType m_subdiag;
362  ComputationInfo m_info;
363  bool m_isInitialized;
364  bool m_eigenvectorsOk;
365 };
366 
383 namespace internal {
384 template<int StorageOrder,typename RealScalar, typename Scalar, typename Index>
385 static void tridiagonal_qr_step(RealScalar* diag, RealScalar* subdiag, Index start, Index end, Scalar* matrixQ, Index n);
386 }
387 
388 template<typename MatrixType>
389 SelfAdjointEigenSolver<MatrixType>& SelfAdjointEigenSolver<MatrixType>
390 ::compute(const MatrixType& matrix, int options)
391 {
392  check_template_parameters();
393 
394  using std::abs;
395  eigen_assert(matrix.cols() == matrix.rows());
396  eigen_assert((options&~(EigVecMask|GenEigMask))==0
397  && (options&EigVecMask)!=EigVecMask
398  && "invalid option parameter");
399  bool computeEigenvectors = (options&ComputeEigenvectors)==ComputeEigenvectors;
400  Index n = matrix.cols();
401  m_eivalues.resize(n,1);
402 
403  if(n==1)
404  {
405  m_eivalues.coeffRef(0,0) = numext::real(matrix.coeff(0,0));
406  if(computeEigenvectors)
407  m_eivec.setOnes(n,n);
408  m_info = Success;
409  m_isInitialized = true;
410  m_eigenvectorsOk = computeEigenvectors;
411  return *this;
412  }
413 
414  // declare some aliases
415  RealVectorType& diag = m_eivalues;
416  MatrixType& mat = m_eivec;
417 
418  // map the matrix coefficients to [-1:1] to avoid over- and underflow.
419  mat = matrix.template triangularView<Lower>();
420  RealScalar scale = mat.cwiseAbs().maxCoeff();
421  if(scale==RealScalar(0)) scale = RealScalar(1);
422  mat.template triangularView<Lower>() /= scale;
423  m_subdiag.resize(n-1);
424  internal::tridiagonalization_inplace(mat, diag, m_subdiag, computeEigenvectors);
425 
426  Index end = n-1;
427  Index start = 0;
428  Index iter = 0; // total number of iterations
429 
430  while (end>0)
431  {
432  for (Index i = start; i<end; ++i)
433  if (internal::isMuchSmallerThan(abs(m_subdiag[i]),(abs(diag[i])+abs(diag[i+1]))))
434  m_subdiag[i] = 0;
435 
436  // find the largest unreduced block
437  while (end>0 && m_subdiag[end-1]==0)
438  {
439  end--;
440  }
441  if (end<=0)
442  break;
443 
444  // if we spent too many iterations, we give up
445  iter++;
446  if(iter > m_maxIterations * n) break;
447 
448  start = end - 1;
449  while (start>0 && m_subdiag[start-1]!=0)
450  start--;
451 
452  internal::tridiagonal_qr_step<MatrixType::Flags&RowMajorBit ? RowMajor : ColMajor>(diag.data(), m_subdiag.data(), start, end, computeEigenvectors ? m_eivec.data() : (Scalar*)0, n);
453  }
454 
455  if (iter <= m_maxIterations * n)
456  m_info = Success;
457  else
458  m_info = NoConvergence;
459 
460  // Sort eigenvalues and corresponding vectors.
461  // TODO make the sort optional ?
462  // TODO use a better sort algorithm !!
463  if (m_info == Success)
464  {
465  for (Index i = 0; i < n-1; ++i)
466  {
467  Index k;
468  m_eivalues.segment(i,n-i).minCoeff(&k);
469  if (k > 0)
470  {
471  std::swap(m_eivalues[i], m_eivalues[k+i]);
472  if(computeEigenvectors)
473  m_eivec.col(i).swap(m_eivec.col(k+i));
474  }
475  }
476  }
477 
478  // scale back the eigen values
479  m_eivalues *= scale;
480 
481  m_isInitialized = true;
482  m_eigenvectorsOk = computeEigenvectors;
483  return *this;
484 }
485 
486 
487 namespace internal {
488 
489 template<typename SolverType,int Size,bool IsComplex> struct direct_selfadjoint_eigenvalues
490 {
491  static inline void run(SolverType& eig, const typename SolverType::MatrixType& A, int options)
492  { eig.compute(A,options); }
493 };
494 
495 template<typename SolverType> struct direct_selfadjoint_eigenvalues<SolverType,3,false>
496 {
497  typedef typename SolverType::MatrixType MatrixType;
498  typedef typename SolverType::RealVectorType VectorType;
499  typedef typename SolverType::Scalar Scalar;
500  typedef typename MatrixType::Index Index;
501 
506  static inline void computeRoots(const MatrixType& m, VectorType& roots)
507  {
508  using std::sqrt;
509  using std::atan2;
510  using std::cos;
511  using std::sin;
512  const Scalar s_inv3 = Scalar(1.0)/Scalar(3.0);
513  const Scalar s_sqrt3 = sqrt(Scalar(3.0));
514 
515  // The characteristic equation is x^3 - c2*x^2 + c1*x - c0 = 0. The
516  // eigenvalues are the roots to this equation, all guaranteed to be
517  // real-valued, because the matrix is symmetric.
518  Scalar c0 = m(0,0)*m(1,1)*m(2,2) + Scalar(2)*m(1,0)*m(2,0)*m(2,1) - m(0,0)*m(2,1)*m(2,1) - m(1,1)*m(2,0)*m(2,0) - m(2,2)*m(1,0)*m(1,0);
519  Scalar c1 = m(0,0)*m(1,1) - m(1,0)*m(1,0) + m(0,0)*m(2,2) - m(2,0)*m(2,0) + m(1,1)*m(2,2) - m(2,1)*m(2,1);
520  Scalar c2 = m(0,0) + m(1,1) + m(2,2);
521 
522  // Construct the parameters used in classifying the roots of the equation
523  // and in solving the equation for the roots in closed form.
524  Scalar c2_over_3 = c2*s_inv3;
525  Scalar a_over_3 = (c2*c2_over_3 - c1)*s_inv3;
526  if(a_over_3<Scalar(0))
527  a_over_3 = Scalar(0);
528 
529  Scalar half_b = Scalar(0.5)*(c0 + c2_over_3*(Scalar(2)*c2_over_3*c2_over_3 - c1));
530 
531  Scalar q = a_over_3*a_over_3*a_over_3 - half_b*half_b;
532  if(q<Scalar(0))
533  q = Scalar(0);
534 
535  // Compute the eigenvalues by solving for the roots of the polynomial.
536  Scalar rho = sqrt(a_over_3);
537  Scalar theta = atan2(sqrt(q),half_b)*s_inv3; // since sqrt(q) > 0, atan2 is in [0, pi] and theta is in [0, pi/3]
538  Scalar cos_theta = cos(theta);
539  Scalar sin_theta = sin(theta);
540  // roots are already sorted, since cos is monotonically decreasing on [0, pi]
541  roots(0) = c2_over_3 - rho*(cos_theta + s_sqrt3*sin_theta); // == 2*rho*cos(theta+2pi/3)
542  roots(1) = c2_over_3 - rho*(cos_theta - s_sqrt3*sin_theta); // == 2*rho*cos(theta+ pi/3)
543  roots(2) = c2_over_3 + Scalar(2)*rho*cos_theta;
544  }
545 
546  static inline bool extract_kernel(MatrixType& mat, Ref<VectorType> res, Ref<VectorType> representative)
547  {
548  using std::abs;
549  Index i0;
550  // Find non-zero column i0 (by construction, there must exist a non zero coefficient on the diagonal):
551  mat.diagonal().cwiseAbs().maxCoeff(&i0);
552  // mat.col(i0) is a good candidate for an orthogonal vector to the current eigenvector,
553  // so let's save it:
554  representative = mat.col(i0);
555  Scalar n0, n1;
556  VectorType c0, c1;
557  n0 = (c0 = representative.cross(mat.col((i0+1)%3))).squaredNorm();
558  n1 = (c1 = representative.cross(mat.col((i0+2)%3))).squaredNorm();
559  if(n0>n1) res = c0/std::sqrt(n0);
560  else res = c1/std::sqrt(n1);
561 
562  return true;
563  }
564 
565  static inline void run(SolverType& solver, const MatrixType& mat, int options)
566  {
567  eigen_assert(mat.cols() == 3 && mat.cols() == mat.rows());
568  eigen_assert((options&~(EigVecMask|GenEigMask))==0
569  && (options&EigVecMask)!=EigVecMask
570  && "invalid option parameter");
571  bool computeEigenvectors = (options&ComputeEigenvectors)==ComputeEigenvectors;
572 
573  MatrixType& eivecs = solver.m_eivec;
574  VectorType& eivals = solver.m_eivalues;
575 
576  // Shift the matrix to the mean eigenvalue and map the matrix coefficients to [-1:1] to avoid over- and underflow.
577  Scalar shift = mat.trace() / Scalar(3);
578  // TODO Avoid this copy. Currently it is necessary to suppress bogus values when determining maxCoeff and for computing the eigenvectors later
579  MatrixType scaledMat = mat.template selfadjointView<Lower>();
580  scaledMat.diagonal().array() -= shift;
581  Scalar scale = scaledMat.cwiseAbs().maxCoeff();
582  if(scale > 0) scaledMat /= scale; // TODO for scale==0 we could save the remaining operations
583 
584  // compute the eigenvalues
585  computeRoots(scaledMat,eivals);
586 
587  // compute the eigenvectors
588  if(computeEigenvectors)
589  {
590  if((eivals(2)-eivals(0))<=Eigen::NumTraits<Scalar>::epsilon())
591  {
592  // All three eigenvalues are numerically the same
593  eivecs.setIdentity();
594  }
595  else
596  {
597  MatrixType tmp;
598  tmp = scaledMat;
599 
600  // Compute the eigenvector of the most distinct eigenvalue
601  Scalar d0 = eivals(2) - eivals(1);
602  Scalar d1 = eivals(1) - eivals(0);
603  Index k(0), l(2);
604  if(d0 > d1)
605  {
606  std::swap(k,l);
607  d0 = d1;
608  }
609 
610  // Compute the eigenvector of index k
611  {
612  tmp.diagonal().array () -= eivals(k);
613  // By construction, 'tmp' is of rank 2, and its kernel corresponds to the respective eigenvector.
614  extract_kernel(tmp, eivecs.col(k), eivecs.col(l));
615  }
616 
617  // Compute eigenvector of index l
618  if(d0<=2*Eigen::NumTraits<Scalar>::epsilon()*d1)
619  {
620  // If d0 is too small, then the two other eigenvalues are numerically the same,
621  // and thus we only have to ortho-normalize the near orthogonal vector we saved above.
622  eivecs.col(l) -= eivecs.col(k).dot(eivecs.col(l))*eivecs.col(l);
623  eivecs.col(l).normalize();
624  }
625  else
626  {
627  tmp = scaledMat;
628  tmp.diagonal().array () -= eivals(l);
629 
630  VectorType dummy;
631  extract_kernel(tmp, eivecs.col(l), dummy);
632  }
633 
634  // Compute last eigenvector from the other two
635  eivecs.col(1) = eivecs.col(2).cross(eivecs.col(0)).normalized();
636  }
637  }
638 
639  // Rescale back to the original size.
640  eivals *= scale;
641  eivals.array() += shift;
642 
643  solver.m_info = Success;
644  solver.m_isInitialized = true;
645  solver.m_eigenvectorsOk = computeEigenvectors;
646  }
647 };
648 
649 // 2x2 direct eigenvalues decomposition, code from Hauke Heibel
650 template<typename SolverType> struct direct_selfadjoint_eigenvalues<SolverType,2,false>
651 {
652  typedef typename SolverType::MatrixType MatrixType;
653  typedef typename SolverType::RealVectorType VectorType;
654  typedef typename SolverType::Scalar Scalar;
655 
656  static inline void computeRoots(const MatrixType& m, VectorType& roots)
657  {
658  using std::sqrt;
659  const Scalar t0 = Scalar(0.5) * sqrt( numext::abs2(m(0,0)-m(1,1)) + Scalar(4)*numext::abs2(m(1,0)));
660  const Scalar t1 = Scalar(0.5) * (m(0,0) + m(1,1));
661  roots(0) = t1 - t0;
662  roots(1) = t1 + t0;
663  }
664 
665  static inline void run(SolverType& solver, const MatrixType& mat, int options)
666  {
667  using std::sqrt;
668  using std::abs;
669 
670  eigen_assert(mat.cols() == 2 && mat.cols() == mat.rows());
671  eigen_assert((options&~(EigVecMask|GenEigMask))==0
672  && (options&EigVecMask)!=EigVecMask
673  && "invalid option parameter");
674  bool computeEigenvectors = (options&ComputeEigenvectors)==ComputeEigenvectors;
675 
676  MatrixType& eivecs = solver.m_eivec;
677  VectorType& eivals = solver.m_eivalues;
678 
679  // map the matrix coefficients to [-1:1] to avoid over- and underflow.
680  Scalar scale = mat.cwiseAbs().maxCoeff();
681  scale = (std::max)(scale,Scalar(1));
682  MatrixType scaledMat = mat / scale;
683 
684  // Compute the eigenvalues
685  computeRoots(scaledMat,eivals);
686 
687  // compute the eigen vectors
688  if(computeEigenvectors)
689  {
690  if((eivals(1)-eivals(0))<=abs(eivals(1))*Eigen::NumTraits<Scalar>::epsilon())
691  {
692  eivecs.setIdentity();
693  }
694  else
695  {
696  scaledMat.diagonal().array () -= eivals(1);
697  Scalar a2 = numext::abs2(scaledMat(0,0));
698  Scalar c2 = numext::abs2(scaledMat(1,1));
699  Scalar b2 = numext::abs2(scaledMat(1,0));
700  if(a2>c2)
701  {
702  eivecs.col(1) << -scaledMat(1,0), scaledMat(0,0);
703  eivecs.col(1) /= sqrt(a2+b2);
704  }
705  else
706  {
707  eivecs.col(1) << -scaledMat(1,1), scaledMat(1,0);
708  eivecs.col(1) /= sqrt(c2+b2);
709  }
710 
711  eivecs.col(0) << eivecs.col(1).unitOrthogonal();
712  }
713  }
714 
715  // Rescale back to the original size.
716  eivals *= scale;
717 
718  solver.m_info = Success;
719  solver.m_isInitialized = true;
720  solver.m_eigenvectorsOk = computeEigenvectors;
721  }
722 };
723 
724 }
725 
726 template<typename MatrixType>
727 SelfAdjointEigenSolver<MatrixType>& SelfAdjointEigenSolver<MatrixType>
728 ::computeDirect(const MatrixType& matrix, int options)
729 {
730  internal::direct_selfadjoint_eigenvalues<SelfAdjointEigenSolver,Size,NumTraits<Scalar>::IsComplex>::run(*this,matrix,options);
731  return *this;
732 }
733 
734 namespace internal {
735 template<int StorageOrder,typename RealScalar, typename Scalar, typename Index>
736 static void tridiagonal_qr_step(RealScalar* diag, RealScalar* subdiag, Index start, Index end, Scalar* matrixQ, Index n)
737 {
738  using std::abs;
739  RealScalar td = (diag[end-1] - diag[end])*RealScalar(0.5);
740  RealScalar e = subdiag[end-1];
741  // Note that thanks to scaling, e^2 or td^2 cannot overflow, however they can still
742  // underflow thus leading to inf/NaN values when using the following commented code:
743 // RealScalar e2 = numext::abs2(subdiag[end-1]);
744 // RealScalar mu = diag[end] - e2 / (td + (td>0 ? 1 : -1) * sqrt(td*td + e2));
745  // This explain the following, somewhat more complicated, version:
746  RealScalar mu = diag[end];
747  if(td==0)
748  mu -= abs(e);
749  else
750  {
751  RealScalar e2 = numext::abs2(subdiag[end-1]);
752  RealScalar h = numext::hypot(td,e);
753  if(e2==0) mu -= (e / (td + (td>0 ? 1 : -1))) * (e / h);
754  else mu -= e2 / (td + (td>0 ? h : -h));
755  }
756 
757  RealScalar x = diag[start] - mu;
758  RealScalar z = subdiag[start];
759  for (Index k = start; k < end; ++k)
760  {
761  JacobiRotation<RealScalar> rot;
762  rot.makeGivens(x, z);
763 
764  // do T = G' T G
765  RealScalar sdk = rot.s() * diag[k] + rot.c() * subdiag[k];
766  RealScalar dkp1 = rot.s() * subdiag[k] + rot.c() * diag[k+1];
767 
768  diag[k] = rot.c() * (rot.c() * diag[k] - rot.s() * subdiag[k]) - rot.s() * (rot.c() * subdiag[k] - rot.s() * diag[k+1]);
769  diag[k+1] = rot.s() * sdk + rot.c() * dkp1;
770  subdiag[k] = rot.c() * sdk - rot.s() * dkp1;
771 
772 
773  if (k > start)
774  subdiag[k - 1] = rot.c() * subdiag[k-1] - rot.s() * z;
775 
776  x = subdiag[k];
777 
778  if (k < end - 1)
779  {
780  z = -rot.s() * subdiag[k+1];
781  subdiag[k + 1] = rot.c() * subdiag[k+1];
782  }
783 
784  // apply the givens rotation to the unit matrix Q = Q * G
785  if (matrixQ)
786  {
787  // FIXME if StorageOrder == RowMajor this operation is not very efficient
788  Map<Matrix<Scalar,Dynamic,Dynamic,StorageOrder> > q(matrixQ,n,n);
789  q.applyOnTheRight(k,k+1,rot);
790  }
791  }
792 }
793 
794 } // end namespace internal
795 
796 } // end namespace Eigen
797 
798 #endif // EIGEN_SELFADJOINTEIGENSOLVER_H
ComputationInfo info() const
Reports whether previous computation was successful.
Definition: SelfAdjointEigenSolver.h:310
MatrixType::Scalar Scalar
Scalar type for matrices of type _MatrixType.
Definition: SelfAdjointEigenSolver.h:81
MatrixType operatorInverseSqrt() const
Computes the inverse square root of the matrix.
Definition: SelfAdjointEigenSolver.h:299
Definition: Constants.h:339
const RealVectorType & eigenvalues() const
Returns the eigenvalues of given matrix.
Definition: SelfAdjointEigenSolver.h:250
SelfAdjointEigenSolver(const MatrixType &matrix, int options=ComputeEigenvectors)
Constructor; computes eigendecomposition of given matrix.
Definition: SelfAdjointEigenSolver.h:153
SelfAdjointEigenSolver(Index size)
Constructor, pre-allocates memory for dynamic-size matrices.
Definition: SelfAdjointEigenSolver.h:131
Computes eigenvalues and eigenvectors of selfadjoint matrices.
Definition: SelfAdjointEigenSolver.h:68
Definition: LDLT.h:16
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:88
Tridiagonal decomposition of a selfadjoint matrix.
Definition: Tridiagonalization.h:61
SelfAdjointEigenSolver & computeDirect(const MatrixType &matrix, int options=ComputeEigenvectors)
Computes eigendecomposition of given matrix using a direct algorithm.
Definition: SelfAdjointEigenSolver.h:728
Computes eigenvalues and eigenvectors of the generalized selfadjoint eigen problem.
Definition: GeneralizedSelfAdjointEigenSolver.h:48
static const int m_maxIterations
Maximum number of iterations.
Definition: SelfAdjointEigenSolver.h:321
MatrixType operatorSqrt() const
Computes the positive-definite square root of the matrix.
Definition: SelfAdjointEigenSolver.h:274
NumTraits< Scalar >::Real RealScalar
Real scalar type for _MatrixType.
Definition: SelfAdjointEigenSolver.h:90
internal::plain_col_type< MatrixType, RealScalar >::type RealVectorType
Type for vector of eigenvalues as returned by eigenvalues().
Definition: SelfAdjointEigenSolver.h:99
Definition: Eigen_Colamd.h:54
const MatrixType & eigenvectors() const
Returns the eigenvectors of given matrix.
Definition: SelfAdjointEigenSolver.h:228
Definition: Constants.h:380
Definition: Constants.h:376
SelfAdjointEigenSolver & compute(const MatrixType &matrix, int options=ComputeEigenvectors)
Computes eigendecomposition of given matrix.
Definition: SelfAdjointEigenSolver.h:390
ComputationInfo
Definition: Constants.h:374
SelfAdjointEigenSolver()
Default constructor for fixed-size matrices.
Definition: SelfAdjointEigenSolver.h:112
Definition: Constants.h:336