C++ Duża prędkość rozkładania w eigendach

Question

Jako część mojego potoku muszę wykonać rozkład złożony z dużej matrycy rzędu 6000x6000. Matryca jest gęsta, więc chyba, że ​​uprościłem problem (w miarę możliwości sot), nie można zastosować żadnej rzadkiej metody.

W tej chwili gram z danymi o zabawkach. Używanie biblioteki Eigen dla macierzy 513x513 potrzebuje ~ 6,5 sekundy, podczas gdy dla matrycy 2049x2049 potrzebuję ~ 130 sekund, co brzmi wygórnie, ponieważ wzrost nie jest liniowy. Zostało to osiągnięte dzięki Eigen::SelfAdjointEigenSolver, podczas gdy innymi metodami, takimi jak Eigen::EigenSolver lub Eigen::ComplexEigenSolver, nie dostałem zauważalnej poprawy. To samo stało się, gdy wypróbowałem Armadillo z arma::eig_sym nawet z opcją "dc", która ma dać szybszy, ale przybliżony wynik. Armadillo ma kilka metod, które zwracają tylko pierwsze X wartości własnych dla przyspieszenia, ale jest to tylko dla rzadkich metod. W tej chwili prawdopodobnie mogę uciec z pierwszych 10-20 wartości własnych.

Czy istnieje sposób lub biblioteka/metoda, która może dać mi godne przyspieszenie?

Answer 1

Polecam wypróbować Arpack-Eigen. Wiem z Octave/Matlab, że może obliczyć największą wartość własną losowego 2049x2049 w ciągu sekundy, a największe 10 w ciągu 5-20 sekund, eigs(rand(2049), 10). Teraz jego dokumentacja help eigs wskazuje na ARPACK. Arpack-Eigen https://github.com/yixuan/arpack-eigen pozwala zażądać 10 wartości własnych z większej macierzy w następujący sposób: SymEigsSolver< double, LARGEST_ALGE, DenseGenMatProd<double> > eigs(&op, 10, 30);.

Answer 2

Spectra służy do pobierania kilku wartości własnych dużej matrycy.

Przykładowy kod do obliczania największej i najmniejszej wartości własne 10 może wyglądać następująco:

#include <Eigen/Core> 
#include <Eigen/Eigenvalues> 
#include <MatOp/DenseGenMatProd.h> 
#include <MatOp/DenseSymShiftSolve.h> 
#include <SymEigsSolver.h> 
#include <iostream> 

using namespace Spectra; 

int main() 
{ 
    srand(0); 
    // We are going to calculate the eigenvalues of M 
    Eigen::MatrixXd A = Eigen::MatrixXd::Random(1000, 1000); 
    Eigen::MatrixXd M = A.transpose() * A; 

    // Matrix operation objects 
    DenseGenMatProd<double> op_largest(M); 
    DenseSymShiftSolve<double> op_smallest(M); 

    // Construct solver object, requesting the largest 10 eigenvalues 
    SymEigsSolver< double, LARGEST_MAGN, DenseGenMatProd<double> > 
     eigs_largest(&op_largest, 10, 30); 

    // Initialize and compute 
    eigs_largest.init(); 
    eigs_largest.compute(); 

    std::cout << "Largest 10 Eigenvalues :\n" << 
     eigs_largest.eigenvalues() << std::endl; 

    // Construct solver object, requesting the smallest 10 eigenvalues 
    SymEigsShiftSolver< double, LARGEST_MAGN, DenseSymShiftSolve<double> > 
     eigs_smallest(&op_smallest, 10, 30, 0.0); 

    eigs_smallest.init(); 
    eigs_smallest.compute(); 
    std::cout << "Smallest 10 Eigenvalues :\n" << 
     eigs_smallest.eigenvalues() << std::endl; 

    return 0; 
}