Hi,

I'm trying to solve the system Ax = B taking advantage that A is sparse, but it didn't work so far.

Code: Select all

int mrow = 1000;
int mcol = 3000;
Eigen::DynamicSparseMatrix<double> data(mrow, mcol);

//**** fills data with coeffRef() +=
...;
Eigen::SparseMatrix<double> fdata = data;
//****

Eigen::SparseMatrix<double> mul = fdata*fdata.transpose();
Eigen::SparseLU<Eigen::SparseMatrix<double>, Eigen::SuperLU > lu(mul);


I've tried with superlu2, superlu3, superlu4, eigen2, eigen3, SparseLU, SparseLLT, SparseLDLT, UmfPack, Taucs, Cholmod and the default backend of SparseLDLT. Most of the time solve() or solveInPlace() return false without any error messages. Otherwise, I get these errors when it calls compute():

Eigen::SparseLLT<MatrixType, Backend>::compute(const MatrixType&) [with MatrixType = Eigen::SparseMatrix<double, 0, int>, int Backend = 0]: Assertion `it.index()==j && "matrix must has non zero diagonal entries and only the lower triangular part must be stored"' failed.

Choldmod : warning: matrix not positive definite

Using superlu4:

** On entry to dgssvx, parameter number 13 had an illegal value

Using superlu3:

free: Invalid pointer

I didn't check whether the matrix is positive definite or not, but it's symmetric with no zero values on the diagonal. A*A.transpose should be positive definite anyway.

It's interesting to note that the dense solvers, Eigen::LU, Eigen::LLT and Eigen::LDLT work fine but they are too slow for my problem.

Any ideas?

thanks for your help

Although your sparse matrix is symmetric with no zero values on the diagonal, maybe you need to delete the upper triangular part of your matrix and only left the lower part, i.e. your input matrix needs to be a lower triangular sparse matrix.

I use cholmod and it works well. Please let me know if SuperLU, UmfPack and Taucs work also for you as I got compile error when using these backends.

Thank you.