ArpackSelfAdjointEigenSolver.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2012 David Harmon <dharmon@gmail.com>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #ifndef EIGEN_ARPACKGENERALIZEDSELFADJOINTEIGENSOLVER_H
11 #define EIGEN_ARPACKGENERALIZEDSELFADJOINTEIGENSOLVER_H
12 
13 #include "../../../../Eigen/Dense"
14 
15 namespace Eigen {
16 
17 namespace internal {
18  template<typename Scalar, typename RealScalar> struct arpack_wrapper;
19  template<typename MatrixSolver, typename MatrixType, typename Scalar, bool BisSPD> struct OP;
20 }
21 
22 
23 
24 template<typename MatrixType, typename MatrixSolver=SimplicialLLT<MatrixType>, bool BisSPD=false>
26 {
27 public:
28  //typedef typename MatrixSolver::MatrixType MatrixType;
29 
31  typedef typename MatrixType::Scalar Scalar;
32  typedef typename MatrixType::Index Index;
33 
41 
48 
56  : m_eivec(),
57  m_eivalues(),
58  m_isInitialized(false),
59  m_eigenvectorsOk(false),
60  m_nbrConverged(0),
62  { }
63 
87  Index nbrEigenvalues, std::string eigs_sigma="LM",
89  : m_eivec(),
90  m_eivalues(),
91  m_isInitialized(false),
92  m_eigenvectorsOk(false),
93  m_nbrConverged(0),
95  {
96  compute(A, B, nbrEigenvalues, eigs_sigma, options, tol);
97  }
98 
122  Index nbrEigenvalues, std::string eigs_sigma="LM",
124  : m_eivec(),
125  m_eivalues(),
126  m_isInitialized(false),
127  m_eigenvectorsOk(false),
128  m_nbrConverged(0),
129  m_nbrIterations(0)
130  {
131  compute(A, nbrEigenvalues, eigs_sigma, options, tol);
132  }
133 
134 
159  Index nbrEigenvalues, std::string eigs_sigma="LM",
161 
185  Index nbrEigenvalues, std::string eigs_sigma="LM",
187 
188 
209  {
210  eigen_assert(m_isInitialized && "ArpackGeneralizedSelfAdjointEigenSolver is not initialized.");
211  eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
212  return m_eivec;
213  }
214 
231  {
232  eigen_assert(m_isInitialized && "ArpackGeneralizedSelfAdjointEigenSolver is not initialized.");
233  return m_eivalues;
234  }
235 
255  {
256  eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
257  eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
258  return m_eivec * m_eivalues.cwiseSqrt().asDiagonal() * m_eivec.adjoint();
259  }
260 
280  {
281  eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
282  eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
283  return m_eivec * m_eivalues.cwiseInverse().cwiseSqrt().asDiagonal() * m_eivec.adjoint();
284  }
285 
291  {
292  eigen_assert(m_isInitialized && "ArpackGeneralizedSelfAdjointEigenSolver is not initialized.");
293  return m_info;
294  }
295 
297  { return m_nbrConverged; }
298 
299  size_t getNbrIterations() const
300  { return m_nbrIterations; }
301 
302 protected:
308 
311 };
312 
313 
314 
315 
316 
317 template<typename MatrixType, typename MatrixSolver, bool BisSPD>
320 ::compute(const MatrixType& A, Index nbrEigenvalues,
321  std::string eigs_sigma, int options, RealScalar tol)
322 {
323  MatrixType B(0,0);
324  compute(A, B, nbrEigenvalues, eigs_sigma, options, tol);
325 
326  return *this;
327 }
328 
329 
330 template<typename MatrixType, typename MatrixSolver, bool BisSPD>
333 ::compute(const MatrixType& A, const MatrixType& B, Index nbrEigenvalues,
334  std::string eigs_sigma, int options, RealScalar tol)
335 {
336  eigen_assert(A.cols() == A.rows());
337  eigen_assert(B.cols() == B.rows());
338  eigen_assert(B.rows() == 0 || A.cols() == B.rows());
340  && (options & EigVecMask) != EigVecMask
341  && "invalid option parameter");
342 
343  bool isBempty = (B.rows() == 0) || (B.cols() == 0);
344 
345  // For clarity, all parameters match their ARPACK name
346  //
347  // Always 0 on the first call
348  //
349  int ido = 0;
350 
351  int n = (int)A.cols();
352 
353  // User options: "LA", "SA", "SM", "LM", "BE"
354  //
355  char whch[3] = "LM";
356 
357  // Specifies the shift if iparam[6] = { 3, 4, 5 }, not used if iparam[6] = { 1, 2 }
358  //
359  RealScalar sigma = 0.0;
360 
361  if (eigs_sigma.length() >= 2 && isalpha(eigs_sigma[0]) && isalpha(eigs_sigma[1]))
362  {
363  eigs_sigma[0] = toupper(eigs_sigma[0]);
364  eigs_sigma[1] = toupper(eigs_sigma[1]);
365 
366  // In the following special case we're going to invert the problem, since solving
367  // for larger magnitude is much much faster
368  // i.e., if 'SM' is specified, we're going to really use 'LM', the default
369  //
370  if (eigs_sigma.substr(0,2) != "SM")
371  {
372  whch[0] = eigs_sigma[0];
373  whch[1] = eigs_sigma[1];
374  }
375  }
376  else
377  {
378  eigen_assert(false && "Specifying clustered eigenvalues is not yet supported!");
379 
380  // If it's not scalar values, then the user may be explicitly
381  // specifying the sigma value to cluster the evs around
382  //
383  sigma = atof(eigs_sigma.c_str());
384 
385  // If atof fails, it returns 0.0, which is a fine default
386  //
387  }
388 
389  // "I" means normal eigenvalue problem, "G" means generalized
390  //
391  char bmat[2] = "I";
392  if (eigs_sigma.substr(0,2) == "SM" || !(isalpha(eigs_sigma[0]) && isalpha(eigs_sigma[1])) || (!isBempty && !BisSPD))
393  bmat[0] = 'G';
394 
395  // Now we determine the mode to use
396  //
397  int mode = (bmat[0] == 'G') + 1;
398  if (eigs_sigma.substr(0,2) == "SM" || !(isalpha(eigs_sigma[0]) && isalpha(eigs_sigma[1])))
399  {
400  // We're going to use shift-and-invert mode, and basically find
401  // the largest eigenvalues of the inverse operator
402  //
403  mode = 3;
404  }
405 
406  // The user-specified number of eigenvalues/vectors to compute
407  //
408  int nev = (int)nbrEigenvalues;
409 
410  // Allocate space for ARPACK to store the residual
411  //
412  Scalar *resid = new Scalar[n];
413 
414  // Number of Lanczos vectors, must satisfy nev < ncv <= n
415  // Note that this indicates that nev != n, and we cannot compute
416  // all eigenvalues of a mtrix
417  //
418  int ncv = std::min(std::max(2*nev, 20), n);
419 
420  // The working n x ncv matrix, also store the final eigenvectors (if computed)
421  //
422  Scalar *v = new Scalar[n*ncv];
423  int ldv = n;
424 
425  // Working space
426  //
427  Scalar *workd = new Scalar[3*n];
428  int lworkl = ncv*ncv+8*ncv; // Must be at least this length
429  Scalar *workl = new Scalar[lworkl];
430 
431  int *iparam= new int[11];
432  iparam[0] = 1; // 1 means we let ARPACK perform the shifts, 0 means we'd have to do it
433  iparam[2] = std::max(300, (int)std::ceil(2*n/std::max(ncv,1)));
434  iparam[6] = mode; // The mode, 1 is standard ev problem, 2 for generalized ev, 3 for shift-and-invert
435 
436  // Used during reverse communicate to notify where arrays start
437  //
438  int *ipntr = new int[11];
439 
440  // Error codes are returned in here, initial value of 0 indicates a random initial
441  // residual vector is used, any other values means resid contains the initial residual
442  // vector, possibly from a previous run
443  //
444  int info = 0;
445 
446  Scalar scale = 1.0;
447  //if (!isBempty)
448  //{
449  //Scalar scale = B.norm() / std::sqrt(n);
450  //scale = std::pow(2, std::floor(std::log(scale+1)));
452  //for (size_t i=0; i<(size_t)B.outerSize(); i++)
453  // for (typename MatrixType::InnerIterator it(B, i); it; ++it)
454  // it.valueRef() /= scale;
455  //}
456 
457  MatrixSolver OP;
458  if (mode == 1 || mode == 2)
459  {
460  if (!isBempty)
461  OP.compute(B);
462  }
463  else if (mode == 3)
464  {
465  if (sigma == 0.0)
466  {
467  OP.compute(A);
468  }
469  else
470  {
471  // Note: We will never enter here because sigma must be 0.0
472  //
473  if (isBempty)
474  {
475  MatrixType AminusSigmaB(A);
476  for (Index i=0; i<A.rows(); ++i)
477  AminusSigmaB.coeffRef(i,i) -= sigma;
478 
479  OP.compute(AminusSigmaB);
480  }
481  else
482  {
483  MatrixType AminusSigmaB = A - sigma * B;
484  OP.compute(AminusSigmaB);
485  }
486  }
487  }
488 
489  if (!(mode == 1 && isBempty) && !(mode == 2 && isBempty) && OP.info() != Success)
490  std::cout << "Error factoring matrix" << std::endl;
491 
492  do
493  {
494  internal::arpack_wrapper<Scalar, RealScalar>::saupd(&ido, bmat, &n, whch, &nev, &tol, resid,
495  &ncv, v, &ldv, iparam, ipntr, workd, workl,
496  &lworkl, &info);
497 
498  if (ido == -1 || ido == 1)
499  {
500  Scalar *in = workd + ipntr[0] - 1;
501  Scalar *out = workd + ipntr[1] - 1;
502 
503  if (ido == 1 && mode != 2)
504  {
505  Scalar *out2 = workd + ipntr[2] - 1;
506  if (isBempty || mode == 1)
508  else
510 
511  in = workd + ipntr[2] - 1;
512  }
513 
514  if (mode == 1)
515  {
516  if (isBempty)
517  {
518  // OP = A
519  //
521  }
522  else
523  {
524  // OP = L^{-1}AL^{-T}
525  //
527  }
528  }
529  else if (mode == 2)
530  {
531  if (ido == 1)
533 
534  // OP = B^{-1} A
535  //
537  }
538  else if (mode == 3)
539  {
540  // OP = (A-\sigmaB)B (\sigma could be 0, and B could be I)
541  // The B * in is already computed and stored at in if ido == 1
542  //
543  if (ido == 1 || isBempty)
545  else
547  }
548  }
549  else if (ido == 2)
550  {
551  Scalar *in = workd + ipntr[0] - 1;
552  Scalar *out = workd + ipntr[1] - 1;
553 
554  if (isBempty || mode == 1)
556  else
558  }
559  } while (ido != 99);
560 
561  if (info == 1)
562  m_info = NoConvergence;
563  else if (info == 3)
564  m_info = NumericalIssue;
565  else if (info < 0)
566  m_info = InvalidInput;
567  else if (info != 0)
568  eigen_assert(false && "Unknown ARPACK return value!");
569  else
570  {
571  // Do we compute eigenvectors or not?
572  //
574 
575  // "A" means "All", use "S" to choose specific eigenvalues (not yet supported in ARPACK))
576  //
577  char howmny[2] = "A";
578 
579  // if howmny == "S", specifies the eigenvalues to compute (not implemented in ARPACK)
580  //
581  int *select = new int[ncv];
582 
583  // Final eigenvalues
584  //
585  m_eivalues.resize(nev, 1);
586 
587  internal::arpack_wrapper<Scalar, RealScalar>::seupd(&rvec, howmny, select, m_eivalues.data(), v, &ldv,
588  &sigma, bmat, &n, whch, &nev, &tol, resid, &ncv,
589  v, &ldv, iparam, ipntr, workd, workl, &lworkl, &info);
590 
591  if (info == -14)
592  m_info = NoConvergence;
593  else if (info != 0)
594  m_info = InvalidInput;
595  else
596  {
597  if (rvec)
598  {
599  m_eivec.resize(A.rows(), nev);
600  for (int i=0; i<nev; i++)
601  for (int j=0; j<n; j++)
602  m_eivec(j,i) = v[i*n+j] / scale;
603 
604  if (mode == 1 && !isBempty && BisSPD)
606 
607  m_eigenvectorsOk = true;
608  }
609 
610  m_nbrIterations = iparam[2];
611  m_nbrConverged = iparam[4];
612 
613  m_info = Success;
614  }
615 
616  delete[] select;
617  }
618 
619  delete[] v;
620  delete[] iparam;
621  delete[] ipntr;
622  delete[] workd;
623  delete[] workl;
624  delete[] resid;
625 
626  m_isInitialized = true;
627 
628  return *this;
629 }
630 
631 
632 // Single precision
633 //
634 extern "C" void ssaupd_(int *ido, char *bmat, int *n, char *which,
635  int *nev, float *tol, float *resid, int *ncv,
636  float *v, int *ldv, int *iparam, int *ipntr,
637  float *workd, float *workl, int *lworkl,
638  int *info);
639 
640 extern "C" void sseupd_(int *rvec, char *All, int *select, float *d,
641  float *z, int *ldz, float *sigma,
642  char *bmat, int *n, char *which, int *nev,
643  float *tol, float *resid, int *ncv, float *v,
644  int *ldv, int *iparam, int *ipntr, float *workd,
645  float *workl, int *lworkl, int *ierr);
646 
647 // Double precision
648 //
649 extern "C" void dsaupd_(int *ido, char *bmat, int *n, char *which,
650  int *nev, double *tol, double *resid, int *ncv,
651  double *v, int *ldv, int *iparam, int *ipntr,
652  double *workd, double *workl, int *lworkl,
653  int *info);
654 
655 extern "C" void dseupd_(int *rvec, char *All, int *select, double *d,
656  double *z, int *ldz, double *sigma,
657  char *bmat, int *n, char *which, int *nev,
658  double *tol, double *resid, int *ncv, double *v,
659  int *ldv, int *iparam, int *ipntr, double *workd,
660  double *workl, int *lworkl, int *ierr);
661 
662 
663 namespace internal {
664 
665 template<typename Scalar, typename RealScalar> struct arpack_wrapper
666 {
667  static inline void saupd(int *ido, char *bmat, int *n, char *which,
668  int *nev, RealScalar *tol, Scalar *resid, int *ncv,
669  Scalar *v, int *ldv, int *iparam, int *ipntr,
670  Scalar *workd, Scalar *workl, int *lworkl, int *info)
671  {
672  EIGEN_STATIC_ASSERT(!NumTraits<Scalar>::IsComplex, NUMERIC_TYPE_MUST_BE_REAL)
673  }
674 
675  static inline void seupd(int *rvec, char *All, int *select, Scalar *d,
676  Scalar *z, int *ldz, RealScalar *sigma,
677  char *bmat, int *n, char *which, int *nev,
678  RealScalar *tol, Scalar *resid, int *ncv, Scalar *v,
679  int *ldv, int *iparam, int *ipntr, Scalar *workd,
680  Scalar *workl, int *lworkl, int *ierr)
681  {
682  EIGEN_STATIC_ASSERT(!NumTraits<Scalar>::IsComplex, NUMERIC_TYPE_MUST_BE_REAL)
683  }
684 };
685 
686 template <> struct arpack_wrapper<float, float>
687 {
688  static inline void saupd(int *ido, char *bmat, int *n, char *which,
689  int *nev, float *tol, float *resid, int *ncv,
690  float *v, int *ldv, int *iparam, int *ipntr,
691  float *workd, float *workl, int *lworkl, int *info)
692  {
693  ssaupd_(ido, bmat, n, which, nev, tol, resid, ncv, v, ldv, iparam, ipntr, workd, workl, lworkl, info);
694  }
695 
696  static inline void seupd(int *rvec, char *All, int *select, float *d,
697  float *z, int *ldz, float *sigma,
698  char *bmat, int *n, char *which, int *nev,
699  float *tol, float *resid, int *ncv, float *v,
700  int *ldv, int *iparam, int *ipntr, float *workd,
701  float *workl, int *lworkl, int *ierr)
702  {
703  sseupd_(rvec, All, select, d, z, ldz, sigma, bmat, n, which, nev, tol, resid, ncv, v, ldv, iparam, ipntr,
704  workd, workl, lworkl, ierr);
705  }
706 };
707 
708 template <> struct arpack_wrapper<double, double>
709 {
710  static inline void saupd(int *ido, char *bmat, int *n, char *which,
711  int *nev, double *tol, double *resid, int *ncv,
712  double *v, int *ldv, int *iparam, int *ipntr,
713  double *workd, double *workl, int *lworkl, int *info)
714  {
715  dsaupd_(ido, bmat, n, which, nev, tol, resid, ncv, v, ldv, iparam, ipntr, workd, workl, lworkl, info);
716  }
717 
718  static inline void seupd(int *rvec, char *All, int *select, double *d,
719  double *z, int *ldz, double *sigma,
720  char *bmat, int *n, char *which, int *nev,
721  double *tol, double *resid, int *ncv, double *v,
722  int *ldv, int *iparam, int *ipntr, double *workd,
723  double *workl, int *lworkl, int *ierr)
724  {
725  dseupd_(rvec, All, select, d, v, ldv, sigma, bmat, n, which, nev, tol, resid, ncv, v, ldv, iparam, ipntr,
726  workd, workl, lworkl, ierr);
727  }
728 };
729 
730 
731 template<typename MatrixSolver, typename MatrixType, typename Scalar, bool BisSPD>
732 struct OP
733 {
734  static inline void applyOP(MatrixSolver &OP, const MatrixType &A, int n, Scalar *in, Scalar *out);
735  static inline void project(MatrixSolver &OP, int n, int k, Scalar *vecs);
736 };
737 
738 template<typename MatrixSolver, typename MatrixType, typename Scalar>
739 struct OP<MatrixSolver, MatrixType, Scalar, true>
740 {
741  static inline void applyOP(MatrixSolver &OP, const MatrixType &A, int n, Scalar *in, Scalar *out)
742 {
743  // OP = L^{-1} A L^{-T} (B = LL^T)
744  //
745  // First solve L^T out = in
746  //
749 
750  // Then compute out = A out
751  //
753 
754  // Then solve L out = out
755  //
758 }
759 
760  static inline void project(MatrixSolver &OP, int n, int k, Scalar *vecs)
761 {
762  // Solve L^T out = in
763  //
765  Matrix<Scalar, Dynamic, Dynamic>::Map(vecs, n, k) = OP.permutationPinv() * Matrix<Scalar, Dynamic, Dynamic>::Map(vecs, n, k);
766 }
767 
768 };
769 
770 template<typename MatrixSolver, typename MatrixType, typename Scalar>
771 struct OP<MatrixSolver, MatrixType, Scalar, false>
772 {
773  static inline void applyOP(MatrixSolver &OP, const MatrixType &A, int n, Scalar *in, Scalar *out)
774 {
775  eigen_assert(false && "Should never be in here...");
776 }
777 
778  static inline void project(MatrixSolver &OP, int n, int k, Scalar *vecs)
779 {
780  eigen_assert(false && "Should never be in here...");
781 }
782 
783 };
784 
785 } // end namespace internal
786 
787 } // end namespace Eigen
788 
789 #endif // EIGEN_ARPACKSELFADJOINTEIGENSOLVER_H
790 
gtsam.examples.DogLegOptimizerExample.int
int
Definition: DogLegOptimizerExample.py:111
Eigen::NumericalIssue
@ NumericalIssue
Definition: Constants.h:444
Eigen
Namespace containing all symbols from the Eigen library.
Definition: jet.h:637
Eigen::internal::arpack_wrapper::seupd
static void seupd(int *rvec, char *All, int *select, Scalar *d, Scalar *z, int *ldz, RealScalar *sigma, char *bmat, int *n, char *which, int *nev, RealScalar *tol, Scalar *resid, int *ncv, Scalar *v, int *ldv, int *iparam, int *ipntr, Scalar *workd, Scalar *workl, int *lworkl, int *ierr)
Definition: ArpackSelfAdjointEigenSolver.h:675
Eigen::internal::arpack_wrapper
Definition: ArpackSelfAdjointEigenSolver.h:18
d
static const double d[K][N]
Definition: igam.h:11
MatrixType
MatrixXf MatrixType
Definition: benchmark-blocking-sizes.cpp:52
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::info
ComputationInfo info() const
Reports whether previous computation was successful.
Definition: ArpackSelfAdjointEigenSolver.h:290
eigen_assert
#define eigen_assert(x)
Definition: Macros.h:1037
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::ArpackGeneralizedSelfAdjointEigenSolver
ArpackGeneralizedSelfAdjointEigenSolver()
Default constructor.
Definition: ArpackSelfAdjointEigenSolver.h:55
B
Definition: test_numpy_dtypes.cpp:299
Eigen::internal::OP::applyOP
static void applyOP(MatrixSolver &OP, const MatrixType &A, int n, Scalar *in, Scalar *out)
Eigen::Success
@ Success
Definition: Constants.h:442
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::getNbrConvergedEigenValues
size_t getNbrConvergedEigenValues() const
Definition: ArpackSelfAdjointEigenSolver.h:296
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::m_isInitialized
bool m_isInitialized
Definition: ArpackSelfAdjointEigenSolver.h:306
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::m_nbrConverged
size_t m_nbrConverged
Definition: ArpackSelfAdjointEigenSolver.h:309
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::m_eigenvectorsOk
bool m_eigenvectorsOk
Definition: ArpackSelfAdjointEigenSolver.h:307
Eigen::ComputeEigenvectors
@ ComputeEigenvectors
Definition: Constants.h:405
sampling::sigma
static const double sigma
Definition: testGaussianBayesNet.cpp:170
Eigen::internal::OP< MatrixSolver, MatrixType, Scalar, false >::project
static void project(MatrixSolver &OP, int n, int k, Scalar *vecs)
Definition: ArpackSelfAdjointEigenSolver.h:778
n
int n
Definition: BiCGSTAB_simple.cpp:1
Eigen::internal::arpack_wrapper< float, float >::seupd
static void seupd(int *rvec, char *All, int *select, float *d, float *z, int *ldz, float *sigma, char *bmat, int *n, char *which, int *nev, float *tol, float *resid, int *ncv, float *v, int *ldv, int *iparam, int *ipntr, float *workd, float *workl, int *lworkl, int *ierr)
Definition: ArpackSelfAdjointEigenSolver.h:696
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::m_info
ComputationInfo m_info
Definition: ArpackSelfAdjointEigenSolver.h:305
Eigen::internal::true_type
Definition: Meta.h:96
Eigen::internal::arpack_wrapper< float, float >::saupd
static void saupd(int *ido, char *bmat, int *n, char *which, int *nev, float *tol, float *resid, int *ncv, float *v, int *ldv, int *iparam, int *ipntr, float *workd, float *workl, int *lworkl, int *info)
Definition: ArpackSelfAdjointEigenSolver.h:688
Eigen::internal::arpack_wrapper< double, double >::seupd
static void seupd(int *rvec, char *All, int *select, double *d, double *z, int *ldz, double *sigma, char *bmat, int *n, char *which, int *nev, double *tol, double *resid, int *ncv, double *v, int *ldv, int *iparam, int *ipntr, double *workd, double *workl, int *lworkl, int *ierr)
Definition: ArpackSelfAdjointEigenSolver.h:718
Eigen::NoConvergence
@ NoConvergence
Definition: Constants.h:446
A
Definition: test_numpy_dtypes.cpp:298
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::eigenvalues
const Matrix< Scalar, Dynamic, 1 > & eigenvalues() const
Returns the eigenvalues of given matrix.
Definition: ArpackSelfAdjointEigenSolver.h:230
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::compute
ArpackGeneralizedSelfAdjointEigenSolver & compute(const MatrixType &A, const MatrixType &B, Index nbrEigenvalues, std::string eigs_sigma="LM", int options=ComputeEigenvectors, RealScalar tol=0.0)
Computes generalized eigenvalues / eigenvectors of given matrix using the external ARPACK library.
Definition: ArpackSelfAdjointEigenSolver.h:333
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::m_eivalues
Matrix< Scalar, Dynamic, 1 > m_eivalues
Definition: ArpackSelfAdjointEigenSolver.h:304
scale
set noclip points set clip one set noclip two set bar set border lt lw set xdata set ydata set zdata set x2data set y2data set boxwidth set dummy y set format x g set format y g set format x2 g set format y2 g set format z g set angles radians set nogrid set key title set key left top Right noreverse box linetype linewidth samplen spacing width set nolabel set noarrow set nologscale set logscale x set set pointsize set encoding default set nopolar set noparametric set set set set surface set nocontour set clabel set mapping cartesian set nohidden3d set cntrparam order set cntrparam linear set cntrparam levels auto set cntrparam points set size set set xzeroaxis lt lw set x2zeroaxis lt lw set yzeroaxis lt lw set y2zeroaxis lt lw set tics in set ticslevel set tics scale
Definition: gnuplot_common_settings.hh:54
j
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2
Eigen::internal::OP< MatrixSolver, MatrixType, Scalar, true >::project
static void project(MatrixSolver &OP, int n, int k, Scalar *vecs)
Definition: ArpackSelfAdjointEigenSolver.h:760
mode
static const DiscreteKey mode(modeKey, 2)
Eigen::ssaupd_
void ssaupd_(int *ido, char *bmat, int *n, char *which, int *nev, float *tol, float *resid, int *ncv, float *v, int *ldv, int *iparam, int *ipntr, float *workd, float *workl, int *lworkl, int *info)
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::m_nbrIterations
size_t m_nbrIterations
Definition: ArpackSelfAdjointEigenSolver.h:310
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::Index
MatrixType::Index Index
Definition: ArpackSelfAdjointEigenSolver.h:32
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::Scalar
MatrixType::Scalar Scalar
Scalar type for matrices of type MatrixType.
Definition: ArpackSelfAdjointEigenSolver.h:31
Eigen::GenEigMask
@ GenEigMask
Definition: Constants.h:418
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::RealVectorType
internal::plain_col_type< MatrixType, RealScalar >::type RealVectorType
Type for vector of eigenvalues as returned by eigenvalues().
Definition: ArpackSelfAdjointEigenSolver.h:47
Eigen::internal::OP< MatrixSolver, MatrixType, Scalar, false >::applyOP
static void applyOP(MatrixSolver &OP, const MatrixType &A, int n, Scalar *in, Scalar *out)
Definition: ArpackSelfAdjointEigenSolver.h:773
Eigen::sseupd_
void sseupd_(int *rvec, char *All, int *select, float *d, float *z, int *ldz, float *sigma, char *bmat, int *n, char *which, int *nev, float *tol, float *resid, int *ncv, float *v, int *ldv, int *iparam, int *ipntr, float *workd, float *workl, int *lworkl, int *ierr)
info
else if n * info
Definition: 3rdparty/Eigen/lapack/cholesky.cpp:18
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::ArpackGeneralizedSelfAdjointEigenSolver
ArpackGeneralizedSelfAdjointEigenSolver(const MatrixType &A, Index nbrEigenvalues, std::string eigs_sigma="LM", int options=ComputeEigenvectors, RealScalar tol=0.0)
Constructor; computes eigenvalues of given matrix.
Definition: ArpackSelfAdjointEigenSolver.h:121
pybind_wrapper_test_script.z
z
Definition: pybind_wrapper_test_script.py:61
Eigen::internal::arpack_wrapper::saupd
static void saupd(int *ido, char *bmat, int *n, char *which, int *nev, RealScalar *tol, Scalar *resid, int *ncv, Scalar *v, int *ldv, int *iparam, int *ipntr, Scalar *workd, Scalar *workl, int *lworkl, int *info)
Definition: ArpackSelfAdjointEigenSolver.h:667
compute
EIGEN_DONT_INLINE void compute(Solver &solver, const MatrixType &A)
Definition: dense_solvers.cpp:25
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::operatorSqrt
Matrix< Scalar, Dynamic, Dynamic > operatorSqrt() const
Computes the positive-definite square root of the matrix.
Definition: ArpackSelfAdjointEigenSolver.h:254
out
std::ofstream out("Result.txt")
Eigen::PlainObjectBase< Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > >::Map
static ConstMapType Map(const Scalar *data)
Definition: PlainObjectBase.h:644
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::operatorInverseSqrt
Matrix< Scalar, Dynamic, Dynamic > operatorInverseSqrt() const
Computes the inverse square root of the matrix.
Definition: ArpackSelfAdjointEigenSolver.h:279
RealScalar
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:47
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::RealScalar
NumTraits< Scalar >::Real RealScalar
Real scalar type for MatrixType.
Definition: ArpackSelfAdjointEigenSolver.h:40
Eigen::dsaupd_
void dsaupd_(int *ido, char *bmat, int *n, char *which, int *nev, double *tol, double *resid, int *ncv, double *v, int *ldv, int *iparam, int *ipntr, double *workd, double *workl, int *lworkl, int *info)
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::m_eivec
Matrix< Scalar, Dynamic, Dynamic > m_eivec
Definition: ArpackSelfAdjointEigenSolver.h:303
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::eigenvectors
const Matrix< Scalar, Dynamic, Dynamic > & eigenvectors() const
Returns the eigenvectors of given matrix.
Definition: ArpackSelfAdjointEigenSolver.h:208
EIGEN_STATIC_ASSERT
#define EIGEN_STATIC_ASSERT(CONDITION, MSG)
Definition: StaticAssert.h:127
OP
#define OP(X)
Definition: gtsam/3rdparty/Eigen/blas/common.h:47
v
Array< int, Dynamic, 1 > v
Definition: Array_initializer_list_vector_cxx11.cpp:1
Eigen::ArpackGeneralizedSelfAdjointEigenSolver
Definition: ArpackSelfAdjointEigenSolver.h:25
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::getNbrIterations
size_t getNbrIterations() const
Definition: ArpackSelfAdjointEigenSolver.h:299
Eigen::internal::OP
Definition: ArpackSelfAdjointEigenSolver.h:19
min
#define min(a, b)
Definition: datatypes.h:19
gtsam::tol
const G double tol
Definition: Group.h:79
ceil
const EIGEN_DEVICE_FUNC CeilReturnType ceil() const
Definition: ArrayCwiseUnaryOps.h:495
gtsam.examples.DogLegOptimizerExample.float
float
Definition: DogLegOptimizerExample.py:113
Eigen::EigVecMask
@ EigVecMask
Definition: Constants.h:407
Eigen::Matrix< Scalar, Dynamic, Dynamic >
Eigen::InvalidInput
@ InvalidInput
Definition: Constants.h:449
internal
Definition: BandTriangularSolver.h:13
Eigen::ArpackGeneralizedSelfAdjointEigenSolver::ArpackGeneralizedSelfAdjointEigenSolver
ArpackGeneralizedSelfAdjointEigenSolver(const MatrixType &A, const MatrixType &B, Index nbrEigenvalues, std::string eigs_sigma="LM", int options=ComputeEigenvectors, RealScalar tol=0.0)
Constructor; computes generalized eigenvalues of given matrix with respect to another matrix.
Definition: ArpackSelfAdjointEigenSolver.h:86
Eigen::ComputationInfo
ComputationInfo
Definition: Constants.h:440
Eigen::internal::OP< MatrixSolver, MatrixType, Scalar, true >::applyOP
static void applyOP(MatrixSolver &OP, const MatrixType &A, int n, Scalar *in, Scalar *out)
Definition: ArpackSelfAdjointEigenSolver.h:741
Eigen::internal::OP::project
static void project(MatrixSolver &OP, int n, int k, Scalar *vecs)
max
#define max(a, b)
Definition: datatypes.h:20
Eigen::dseupd_
void dseupd_(int *rvec, char *All, int *select, double *d, double *z, int *ldz, double *sigma, char *bmat, int *n, char *which, int *nev, double *tol, double *resid, int *ncv, double *v, int *ldv, int *iparam, int *ipntr, double *workd, double *workl, int *lworkl, int *ierr)
options
Definition: options.h:16
Eigen::NumTraits
Holds information about the various numeric (i.e. scalar) types allowed by Eigen.
Definition: NumTraits.h:232
Eigen::internal::arpack_wrapper< double, double >::saupd
static void saupd(int *ido, char *bmat, int *n, char *which, int *nev, double *tol, double *resid, int *ncv, double *v, int *ldv, int *iparam, int *ipntr, double *workd, double *workl, int *lworkl, int *info)
Definition: ArpackSelfAdjointEigenSolver.h:710
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9
Scalar
SCALAR Scalar
Definition: bench_gemm.cpp:46
Eigen::Index
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74


gtsam
Author(s):
autogenerated on Fri Nov 1 2024 03:31:52