New user self-registration is currently disabled. Please email eigen-core-team @ lists.tuxfamily.org if you need an account.
Bug 1053 - SparseLU factorization error: THE MATRIX IS STRUCTURALLY SINGULAR ... ZERO COLUMN
SparseLU factorization error: THE MATRIX IS STRUCTURALLY SINGULAR ... ZERO CO...
Status: RESOLVED FIXED
Product: Eigen
Classification: Unclassified
Component: Sparse
3.2
All All
: Normal Wrong Result
Assigned To: Nobody
:
Depends on:
Blocks:
  Show dependency treegraph
 
Reported: 2015-08-07 15:22 UTC by cDc
Modified: 2015-09-03 09:26 UTC (History)
2 users (show)



Attachments
test (4.81 KB, application/zip)
2015-08-07 15:22 UTC, cDc
no flags Details

Description cDc 2015-08-07 15:22:26 UTC
Created attachment 597 [details]
test

SparseLU factorization step fails always (except for very small and particular patrices) with an error message similar to:

THE MATRIX IS STRUCTURALLY SINGULAR ... ZERO COLUMN AT 428

The problem appears only when compiling with RowMajor as define:

#define EIGEN_DEFAULT_TO_ROW_MAJOR

Please find attached a small code that reproduces the problem with a test file (which is not important, all problems give a similar error). Please run as:

SparseLU.exe coeff.bin

If the EIGEN_DEFAULT_TO_ROW_MAJOR define is commented the factorization completes as expected.
Comment 1 Christoph Hertzberg 2015-08-07 21:15:17 UTC
EIGEN_DEFAULT_TO_ROW_MAJOR is documented as experimental/"should not be used by real-word code".

Nevertheless, I fixed the root of this bug in 3.2 and devel branch:
https://bitbucket.org/eigen/eigen/commits/890ac1744
https://bitbucket.org/eigen/eigen/commits/9b933a200

Btw: Your attachment was a 7z not a zip. Also, try to make more minimal examples, if possible.
Comment 2 cDc 2015-08-08 09:49:36 UTC
Thank you for the quick fix. The factorization works now as expected, but now the solve() is failing (crash).

PS: I can generate as many minimal examples as you want, but they are all above 1MB and the limit for attaching a file is 0.25MB
Comment 3 Christoph Hertzberg 2015-08-08 15:20:14 UTC
I meant "'more minimal' examples" not "more 'minimal examples'", i.e. try to reduce the problem as much as possible -- in many cases this already reveals the root of the problem (and makes resolving it easier for us).
Could you at least point out where the crash occurs?
Also, are the results correct in cases it does not crash?
Comment 4 cDc 2015-08-08 15:39:09 UTC
The error is at, where 'r' is a bad pointer:

Eigen::internal::triangular_solve_matrix<float,__int64,2,6,0,1,0>::run(__int64 size, __int64 otherSize, const float * _tri, __int64 triStride, float * _other, __int64 otherStride, Eigen::internal::level3_blocking<float,float> & blocking) Line 303




For smaller problems, it asserts at:

Eigen::DenseCoeffsBase<Eigen::Matrix<float,-1,1,0,-1,1>,1>::operator()(__int64 row, __int64 col) Line 336
Eigen::SparseLUMatrixUReturnType<Eigen::internal::MappedSuperNodalMatrix<float,int>,Eigen::MappedSparseMatrix<float,0,int> >::solveInPlace<Eigen::Matrix<float,-1,1,0,-1,1> >(Eigen::MatrixBase<Eigen::Matrix<float,-1,1,0,-1,1> > & X) Line 765
Eigen::SparseLU<Eigen::SparseMatrix<float,0,int>,Eigen::COLAMDOrdering<int> >::_solve<Eigen::Map<Eigen::Matrix<float,-1,1,0,-1,1>,0,Eigen::Stride<0,3> >,Eigen::Matrix<float,-1,1,0,-1,1> >(const Eigen::MatrixBase<Eigen::Map<Eigen::Matrix<float,-1,1,0,-1,1>,0,Eigen::Stride<0,3> > > & B, Eigen::MatrixBase<Eigen::Matrix<float,-1,1,0,-1,1> > & X_base) Line 242
Eigen::internal::solve_retval<Eigen::SparseLU<Eigen::SparseMatrix<float,0,int>,Eigen::COLAMDOrdering<int> >,Eigen::Map<Eigen::Matrix<float,-1,1,0,-1,1>,0,Eigen::Stride<0,3> > >::evalTo<Eigen::Matrix<float,-1,1,0,-1,1> >(Eigen::Matrix<float,-1,1,0,-1,1> & dst) Line 787
Eigen::internal::solve_retval_base<Eigen::SparseLU<Eigen::SparseMatrix<float,0,int>,Eigen::COLAMDOrdering<int> >,Eigen::Map<Eigen::Matrix<float,-1,1,0,-1,1>,0,Eigen::Stride<0,3> > >::evalTo<Eigen::Matrix<float,-1,1,0,-1,1> >(Eigen::Matrix<float,-1,1,0,-1,1> & dst) Line 52
Eigen::ReturnByValue<Eigen::internal::solve_retval_base<Eigen::SparseLU<Eigen::SparseMatrix<float,0,int>,Eigen::COLAMDOrdering<int> >,Eigen::Map<Eigen::Matrix<float,-1,1,0,-1,1>,0,Eigen::Stride<0,3> > > >::evalTo<Eigen::Matrix<float,-1,1,0,-1,1> >(Eigen::Matrix<float,-1,1,0,-1,1> & dst) Line 61
Eigen::Matrix<float,-1,1,0,-1,1>::Matrix<float,-1,1,0,-1,1><Eigen::internal::solve_retval_base<Eigen::SparseLU<Eigen::SparseMatrix<float,0,int>,Eigen::COLAMDOrdering<int> >,Eigen::Map<Eigen::Matrix<float,-1,1,0,-1,1>,0,Eigen::Stride<0,3> > > >(const Eigen::ReturnByValue<Eigen::internal::solve_retval_base<Eigen::SparseLU<Eigen::SparseMatrix<float,0,int>,Eigen::COLAMDOrdering<int> >,Eigen::Map<Eigen::Matrix<float,-1,1,0,-1,1>,0,Eigen::Stride<0,3> > > > & other) Line 296



For the sake of completness, there are these warnings as well:

SparseLU.h(733): warning C4244: 'initializing': conversion from '__int64' to 'int', possible loss of data
SparseLU.h(734): warning C4244: 'initializing': conversion from '__int64' to 'int', possible loss of data
SparseLU_SupernodalMatrix.h(236): warning C4244: 'initializing': conversion from '__int64' to 'int', possible loss of data
SparseLU_SupernodalMatrix.h(237): warning C4244: 'initializing': conversion from '__int64' to 'int', possible loss of data


I never get a complete solve to know if the results are correct or not.
Comment 5 Gael Guennebaud 2015-09-03 09:26:26 UTC
https://bitbucket.org/eigen/eigen/commits/2b2f1f084aa2/
Summary:     Bug 1053: fix SuplerLU::solve with EIGEN_DEFAULT_TO_ROW_MAJOR

Note You need to log in before you can comment on or make changes to this bug.