ARPACK support

Submitted by Gael Guennebaud @ggael

Assigned to Nobody

Link to original bugzilla bug (#230)
Version: 3.0

Description

One way to get eigen decomposition on sparse matrices is to wrap ARPACK++, here is a start:
http://forum.kde.org/viewtopic.php?f=74&t=84238&p=192268#p192198

and here is a copy of schizofrenikh's message from the forum:

#include <arpack++/arssym.h>  
#include &lt;vector&gt;  
#include <eigen2/Eigen/Eigen>  

class TMatrixForArpackpp {  
private:  
    Eigen::SparseMatrix&lt; double &gt; &M;  
public:  
    TMatrixForArpackpp( Eigen::SparseMatrix&lt; double &gt; &_M ) :  
        M(_M) {  
        ;  
    };  
    void MultMv( double *v, double *w ) {  
        Eigen::Map< Eigen::VectorXd > v_map = Eigen::VectorXd::Map( v, M.rows() );  
        Eigen::Map< Eigen::VectorXd > w_map = Eigen::VectorXd::Map( w, M.cols() );  
        w_map = M.marked<Eigen::SelfAdjoint|Eigen::LowerTriangular>()*v_map;  
    }  
};  


void Lanczos_symm( Eigen::SparseMatrix&lt; double &gt; &M, std::vector&lt; double &gt; &eigenValues, std::vector< Eigen::VectorXd > &eigenVectors ) {  
    int n = M.cols();  
    TMatrixForArpackpp M_arp( M );  
    ARSymStdEig< double, TMatrixForArpackpp > prob( n, n-1, &M_arp, &TMatrixForArpackpp::MultMv, "SM" ); //"SM" to get smallest eigenvalues  
    prob.FindEigenvectors();  
    //prob.FindEigenvalues();  
    int m = prob.ConvergedEigenvalues();  
    eigenValues.resize(m);  
    eigenVectors.resize(m);  
    eigenValues.assign( prob.RawEigenvalues(), prob.RawEigenvalues()+m );  
    for( int i=0; i<m; i++ )  
        eigenVectors[i] = Eigen::VectorXd::Map( prob.RawEigenvector(i), n);  
    /*eigenValues = *( prob.StlEigenvalues() );  
    eigenVectors = *( prob.StlEigenvectors() );*/  
}  

To compile with arpack on my debian system I had to add -lblas -lgfortran -lgfortranbegin -lnsl -larpack -larpack++