Autovalori sparsi utilizzando eigen3/sparse

Question

Esiste un modo distinto ed efficace per trovare autovalori e autovettori di una matrice reale, simmetrica, molto grande, ad esempio 10000x10000, sparse in Eigen3? Esiste un risolutore di autovalori per matrici dense ma che non fa uso della proprietà della matrice, ad es. è simmetria. Inoltre non voglio memorizzare la matrice in modo denso.

Oppure (alternativa) c'è una libreria migliore (+ meglio documentata) per farlo?

Answer 1

Armadillo farà questo usando eigs_sym

Si noti che il calcolo tutti gli autovalori è un'operazione molto costosa qualunque cosa tu faccia, di solito quello che è fatto è quello di trovare solo il k più grande, o autovalori più piccoli (che è ciò che questo lo farà).

Answer 2

Per Eigen, esiste una libreria denominata Spectra. Come è descritto nella sua pagina web, Spectra è una riprogettazione della libreria ARPACK usando il linguaggio C++.

A differenza di Armadillo, suggerito in another answer, Spectra supporta long double e qualsiasi altro tipo a virgola mobile reale (ad esempio boost::multiprecision::float128).

Ecco un esempio di utilizzo (stessa della versione nella documentazione, ma adattato per gli esperimenti con diversi tipi in virgola mobile):

#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'; 
    } 
}