Hi, I know someone may ask this question before.
I am trying to solve a SPD sparse system with sparse matrix M. However the rhs is also a very sparse matrix (<0.2% nnz). How can I solve this system efficiently? I use
Simplicial LLT to factorize M like:

SimplicialLLT<SparseMatrix<double>> solver;
solver.compute(M);
M=M.solve(B)

It turns out that, M. solve (B) takes an extremely long time comparing with X * B where X is some dense matrix with the same dimension of M.

Any way to improve the performance for M. solves (B)?

Thanks in advance!

The problem is that in general even if both M and B are very sparse, M^1 * B can be rather dense. Cannot you avoid computing M^-1 * B explicitly?

I am actually trying to get the matrix corresponding to B^T M^{-1} B, I guess I must explicitly compute M^{-1} B first or could you suggest any better walk-around?
Thank you for response.

and how do you use Q = B^T M^{-1} B ? If you use it for matrix-vector products like Q v or v^T Q v, then better keep Q "unpacked".

ggael wrote:and how do you use Q = B^T M^{-1} B ? If you use it for matrix-vector products like Q v or v^T Q v, then better keep Q "unpacked".

Hi ggael,
Thanks for your response and for developing Eigen - it is indeed helpful.

Q will be a square SPD dense, and it is used to solve a linear system... so I do need Q explicitly

This is what I am currently doing :

B is a sparse, so I manually extract columns in B that contain nonzeros in a std::vector<VectorXd>
Then I just use SparseLLT to solve each nonzero vectors in B to get another set of vectors corresponding to M^{-1} B
Finally I compute B^T (M^{-1}B)

This is much faster than directly solving B after factorizing M.
I am not sure if it is the best solution - because even M^{-1} corresponds to a dense matrix, a dense matrix times B (B is very sparse) is still much faster than what I did here. Right hand side vector in a linear system seems not contributing much performance improvement. Or do we have some tricks to handle sparse rhs?

Thanks

oh, when you mean that B is sparse, you actually mean that many of its columns are zeros? Then yes you should construct a "packed" version of B without empty column, and then compute M^{-1} B into a dense matrix.