Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Sparse eigenvalues using eigen3/sparse

Is there an distinct and effective way of finding eigenvalues and eigenvectors of a real, symmetrical, very large, let's say 10000x10000, sparse matrix in Eigen3? There is an eigenvalue solver for dense matrices but that doesn't make use of the property of the matrix e.g. it's symmetry. Furthermore I don't want to store the matrix in dense.

Or (alternative) is there a better (+better documented) library to do that?

like image 605
Philipp Avatar asked May 12 '15 10:05

Philipp


1 Answers

For Eigen, there's a library named Spectra. As is described on its web page, Spectra is a redesign of the ARPACK library using C++ language.

Unlike Armadillo, suggested in another answer, Spectra does support long double and any other real floating-point type (e.g. boost::multiprecision::float128).

Here's an example of usage (same as the version in documentation, but adapted for experiments with different floating-point types):

#include <Eigen/Core>
#include <SymEigsSolver.h>  // Also includes <MatOp/DenseSymMatProd.h>
#include <iostream>
#include <limits>

int main()
{
    using Real=long double;
    using Matrix=Eigen::Matrix<Real, Eigen::Dynamic, Eigen::Dynamic>;

    // We are going to calculate the eigenvalues of M
    const auto A = Matrix::Random(10, 10);
    const Matrix M = A + A.transpose();

    // Construct matrix operation object using the wrapper class DenseGenMatProd
    Spectra::DenseSymMatProd<Real> op(M);

    // Construct eigen solver object, requesting the largest three eigenvalues
    Spectra::SymEigsSolver<Real,
                           Spectra::LARGEST_ALGE,
                           Spectra::DenseSymMatProd<Real>> eigs(&op, 3, 6);

    // Initialize and compute
    eigs.init();
    const auto nconv = eigs.compute();
    std::cout << nconv << " eigenvalues converged.\n";

    // Retrieve results
    if(eigs.info() == Spectra::SUCCESSFUL)
    {
        const auto evalues = eigs.eigenvalues();
        std::cout.precision(std::numeric_limits<Real>::digits10);
        std::cout << "Eigenvalues found:\n" << evalues << '\n';
    }
}
like image 196
Ruslan Avatar answered Oct 09 '22 11:10

Ruslan