Редкие собственные значения с использованием собственных/разреженных

Question

Есть ли отличный и эффективный способ нахождения собственных значений и собственных векторов реальной, симметричной, очень большой, скажем, 10000x10000, разреженной матрицы в Eigen3? Для плотных матриц существует решатель на собственные значения, но он не использует свойство матрицы, например. это симметрия. Кроме того, я не хочу, чтобы матрица была плотной.

Или (альтернативно) есть ли библиотека (лучше документирована) для этого?

Answer 1

Armadillo будет делать это с помощью eigs_sym

Обратите внимание, что вычислительное все собственные очень дорогая операция, что бы вы ни делали, как правило, то, что делается, чтобы найти только K наибольшее или наименьшее собственные значения (что это будет делать).

Answer 2

Для Eigen существует библиотека с именем Spectra. Как описано на его веб-странице, Spectra представляет собой редизайн библиотеки ARPACK с использованием языка C++.

В отличие от Armadillo, предложенного в another answer, Spectra поддерживает long double и любые другие реальные типы с плавающей точкой (например, boost::multiprecision::float128).

Вот пример использования (такой же, как версия в документации, но адаптирована для экспериментов с различными типами с плавающей точкой):

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