Dear Eigen Folks,
Keir and I are quite interested in Eigen having world class sparse matrix support. We spent some time talking about it and here are some thoughts. Some of this stuff is from my previous email on the subject.
1. Construction.
Eigen should support a (i,j,s) triplet based constructors of the form
SparseMatrix* M = BuildSparseMatrix( vector<int> rows, vector<int> cols, vector<double> values, int num_rows, int num_cols);
SparseMatrix* M = BuildSparseMatrix( int* rows, int* cols, double* values, int nnz, int num_rows, int num_cols);
I think this a much more useful way of building a sparse matrix rather than the DynamicSparseMatrix route.
2. BLAS
We should first focus on building a sparse BLAS, but we don't actually need the whole thing. The most common operations that we should focus on are
vector * scalar
vector / scalar
vector + scalar (maybe)
vector - scalar (maybe)
vector.norm() (all p norms)
vector +/- vector
dot(vector, vector)
Matrix * vector
Matrix^T * vector
If we have well optimized robust implementations of these operations, a majority of sparse matrix usage will be satisfied.. In particular it gives us enough operators to build iterative solvers with. Which brings up the issue of efficient implementation of SpMV - sparse matrix dense vector product. This is the source of much activity in the high performance computing lately, and a naive implementation of spmv has extremely poor performance. It is hard to hard code a strategy that works well for all matrices.
One direction we could look into is wrapping the OSKI library
3. Iterative Methods
A library of sparse iterative solvers built on top of the sparse matrix support. This may or may not be part of Eigen proper. A very useful thing to have here are implementations of a variety of basic preconditioning schemes. Not just for square/symmetric systems but also for rectangular systems and non-symmetric systems.
4. Direct methods
Here we need better support for exposing the interface to symbolic factorization reuse, and the various ordering options that solvers provide.
Both Keir and I are interested in helping with this.
Sameer