Bug 485 - Diagonal Preconditioner's factorize() routine causes seg fault for (square) SparseMatrix of size 128 and greater
: Diagonal Preconditioner's factorize() routine causes seg fault for (square) ...
Status: RESOLVED FIXED
Product: Eigen
Classification: Unclassified
Component: Sparse
: 3.1
: x86 - 64-bit Linux
: Normal normal
Assigned To: Nobody
:
:
:
:
:
  Show dependency treegraph
 
Reported: 2012-06-30 03:50 UTC by Patrick
Modified: 2012-07-13 21:14 UTC (History)
2 users (show)



Attachments

Description Patrick 2012-06-30 03:50:12 UTC
I've been getting segfaults while attempting to use the ConjugateGradient
solver for a sparse system. gdb traces the fault back to the it.value() call
near line 83 of in the factorize() method in BasicPreconditioners.h , when
attempting to compute a new ConjugateGradient solver. Strangely this seems to
happen exactly when the problem size is 128 or greater. 

----This works:

int n = 127;
SparseMatrix<double> A(n,n);
// ...
Eigen::ConjugateGradient< SparseMatrix<double> > cg(A);

----This works:

int n = 128;
SparseMatrix<double> A(n,n);
// ...
Eigen::ConjugateGradient< SparseMatrix<double>, Eigen::Lower,               
Eigen::IdentityPreconditioner > cg(A);

-----This induces a seg fault:

int n = 128;
SparseMatrix<double> A(n,n);
// ...
Eigen::ConjugateGradient< SparseMatrix<double> > cg(A);
Comment 1 Desire NUENTSA 2012-07-12 17:23:15 UTC
Diagonal preconditioner takes the inverse of diagonal elements.
have you checked that all those elements are not null ?
Comment 2 Patrick 2012-07-12 20:50:38 UTC
I've checked again and I don't think that's the problem. Here's another test
program which shows the issue for the identity matrix. The particular matrix
size which produces the problem seems to depend on the scalar type for the
matrix. (Here I'm including Eigen 3.1, the latest stable release)

---------------------------------------------------------------------------

/*
g++ -I ../../../Desktop/eigen-eigen-ca142d0540d3/ -o showBug showBug.cpp
*/

#include<Eigen/Sparse>
using namespace Eigen;

typedef double myScalar;

int main(int argc, char* argv[]){
    int n = 128; // segfault if n is 128 or greater and myScalar is double, 
                 // segfault if n is 182 or greater and myScalar is float

    SparseMatrix<myScalar> A(n,n);
    for(int i=0;i<n;i++)
        A.insert(i,i) = (myScalar) 1.0 ;

    ConjugateGradient< SparseMatrix<myScalar> > cg(A);

    return 0;   
}
Comment 3 Desire NUENTSA 2012-07-13 10:21:01 UTC
The code is fine. I can't reproduce the error with Eigen 3.1 for any size of
the matrix.
Comment 4 Gael Guennebaud 2012-07-13 15:44:13 UTC
I cannot reproduce too. Make sure you use the 3.1 release, and maybe running
your problem under valgrind might help to discover the issue.
Comment 5 Patrick 2012-07-13 16:23:56 UTC
Good to know - I'll look for another cause - thanks for your time and sorry to
waste it!
Comment 6 Patrick 2012-07-13 17:48:20 UTC
I'm having no luck locating the problem. For reference, in case this comes up
again, here are valgrind and gdb outputs (path to library removed in gdb
output). "Invalid read of size 8" changes to "invalid read of size 4" when
using float instead of double. As far as I can tell, the error is only on the
attempt to read the first diagonal entry of the matrix, and the rest proceed
normally.

----------------------------------
valgrind  ./showBug
==21400== Memcheck, a memory error detector.
==21400== Copyright (C) 2002-2007, and GNU GPL'd, by Julian Seward et al.
==21400== Using LibVEX rev 1854, a library for dynamic binary translation.
==21400== Copyright (C) 2004-2007, and GNU GPL'd, by OpenWorks LLP.
==21400== Using valgrind-3.3.1, a dynamic binary instrumentation framework.
==21400== Copyright (C) 2000-2007, and GNU GPL'd, by Julian Seward et al.
==21400== For more details, rerun with: -v
==21400== Invalid read of size 8
==21400==    at 0x401797: Eigen::DenseBase<Eigen::Matrix<double, -1, -1, 0, -1,
-1> >::InnerIterator::value() const (CoreIterators.h:56)
==21400==    by 0x4046A4: Eigen::DiagonalPreconditioner<double>&
Eigen::DiagonalPreconditioner<double>::factorize<Eigen::SparseMatrix<double, 0,
int> >(Eigen::SparseMatrix<double, 0, int> const&) (BasicPreconditioners.h:83)
==21400==    by 0x404722: Eigen::DiagonalPreconditioner<double>&
Eigen::DiagonalPreconditioner<double>::compute<Eigen::SparseMatrix<double, 0,
int> >(Eigen::SparseMatrix<double, 0, int> const&) (BasicPreconditioners.h:94)
==21400==    by 0x404751:
Eigen::IterativeSolverBase<Eigen::ConjugateGradient<Eigen::SparseMatrix<double,
0, int>, 1, Eigen::DiagonalPreconditioner<double> >
>::compute(Eigen::SparseMatrix<double, 0, int> const&)
(IterativeSolverBase.h:121)
==21400==    by 0x4047BE:
Eigen::IterativeSolverBase<Eigen::ConjugateGradient<Eigen::SparseMatrix<double,
0, int>, 1, Eigen::DiagonalPreconditioner<double> >
>::IterativeSolverBase(Eigen::SparseMatrix<double, 0, int> const&)
(IterativeSolverBase.h:70)
==21400==    by 0x404830: Eigen::ConjugateGradient<Eigen::SparseMatrix<double,
0, int>, 1, Eigen::DiagonalPreconditioner<double>
>::ConjugateGradient(Eigen::SparseMatrix<double, 0, int> const&)
(ConjugateGradient.h:193)
==21400==    by 0x400FB6: main (showBug.cpp:28)
==21400==  Address 0x58fe528 is 0 bytes inside a block of size 131,072 free'd
==21400==    at 0x4C243AF: free (in
/usr/lib64/valgrind/amd64-linux/vgpreload_memcheck.so)
==21400==    by 0x403B86: Eigen::internal::aligned_free(void*) (Memory.h:240)
==21400==    by 0x403B9D: void
Eigen::internal::conditional_aligned_free<true>(void*) (Memory.h:314)
==21400==    by 0x403BB8: void
Eigen::internal::conditional_aligned_delete_auto<double, true>(double*,
unsigned long) (Memory.h:440)
==21400==    by 0x403BEA: Eigen::DenseStorage<double, -1, -1, -1,
0>::~DenseStorage() (DenseStorage.h:215)
==21400==    by 0x403C02: Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1,
0, -1, -1> >::~PlainObjectBase() (PlainObjectBase.h:89)
==21400==    by 0x403C1A: Eigen::Matrix<double, -1, -1, 0, -1, -1>::~Matrix()
(Matrix.h:144)
==21400==    by 0x40460B: Eigen::DiagonalPreconditioner<double>&
Eigen::DiagonalPreconditioner<double>::factorize<Eigen::SparseMatrix<double, 0,
int> >(Eigen::SparseMatrix<double, 0, int> const&) (BasicPreconditioners.h:80)
==21400==    by 0x404722: Eigen::DiagonalPreconditioner<double>&
Eigen::DiagonalPreconditioner<double>::compute<Eigen::SparseMatrix<double, 0,
int> >(Eigen::SparseMatrix<double, 0, int> const&) (BasicPreconditioners.h:94)
==21400==    by 0x404751:
Eigen::IterativeSolverBase<Eigen::ConjugateGradient<Eigen::SparseMatrix<double,
0, int>, 1, Eigen::DiagonalPreconditioner<double> >
>::compute(Eigen::SparseMatrix<double, 0, int> const&)
(IterativeSolverBase.h:121)
==21400==    by 0x4047BE:
Eigen::IterativeSolverBase<Eigen::ConjugateGradient<Eigen::SparseMatrix<double,
0, int>, 1, Eigen::DiagonalPreconditioner<double> >
>::IterativeSolverBase(Eigen::SparseMatrix<double, 0, int> const&)
(IterativeSolverBase.h:70)
==21400==    by 0x404830: Eigen::ConjugateGradient<Eigen::SparseMatrix<double,
0, int>, 1, Eigen::DiagonalPreconditioner<double>
>::ConjugateGradient(Eigen::SparseMatrix<double, 0, int> const&)
(ConjugateGradient.h:193)
==21400==
==21400== ERROR SUMMARY: 128 errors from 1 contexts (suppressed: 3 from 1)
==21400== malloc/free: in use at exit: 0 bytes in 0 blocks.
==21400== malloc/free: 133 allocs, 133 frees, 16,782,340 bytes allocated.
==21400== For counts of detected errors, rerun with: -v
==21400== All heap blocks were freed -- no leaks are possible.

-------------------------------------------------------------
> gdb ./showBug
.....................
(gdb) run
.....................
Program received signal SIGSEGV, Segmentation fault.
0x0000000000401797 in Eigen::DenseBase<Eigen::Matrix<double,
-0x00000000000000001, -0x00000000000000001, 0, -0x00000000000000001,
-0x00000000000000001> >::InnerIterator::value (this=0x7fffffffd640) at
PATH/TO/EIGEN/Eigen/src/SparseCore/CoreIterators.h:56
56                                : m_expression.coeff(m_inner, m_outer);
Comment 7 Patrick 2012-07-13 17:57:31 UTC
Also perhaps of interest is that valgrind reports an illegal read regardless of
the system size n, even though a segfault is thrown only for large enough
values.
Comment 8 Patrick 2012-07-13 19:46:29 UTC
Using gcc 4.6.1, the program above runs.

The crashes above were using gcc 4.3.2. Here's a dump of g++ -v :

> g++ -v
Using built-in specs.
Target: x86_64-suse-linux
Configured with: ../configure --prefix=/usr --infodir=/usr/share/info
--mandir=/usr/share/man --libdir=/usr/lib64 --libexecdir=/usr/lib64
--enable-languages=c,c++,objc,fortran,obj-c++,java,ada
--enable-checking=release --with-gxx-include-dir=/usr/include/c++/4.3
--enable-ssp --disable-libssp --with-bugurl=http://bugs.opensuse.org/
--with-pkgversion='SUSE Linux' --disable-libgcj --disable-libmudflap
--with-slibdir=/lib64 --with-system-zlib --enable-__cxa_atexit
--enable-libstdcxx-allocator=new --disable-libstdcxx-pch
--enable-version-specific-runtime-libs --program-suffix=-4.3
--enable-linux-futex --without-system-libunwind --with-cpu=generic
--build=x86_64-suse-linux
Thread model: posix
gcc version 4.3.2 [gcc-4_3-branch revision 141291] (SUSE Linux)
Comment 9 Gael Guennebaud 2012-07-13 20:08:42 UTC
thanbks for these report, I found and solved the problem. I'll commit it very
soon.
Comment 10 Gael Guennebaud 2012-07-13 21:14:15 UTC
https://bitbucket.org/eigen/eigen/changeset/28e7dfbc196f/
changeset:   28e7dfbc196f
user:        ggael
date:        2012-07-13 20:54:38
summary:     fix bug 485: conflict between a typedef and template type
parameter

https://bitbucket.org/eigen/eigen/changeset/60062b16e638/
changeset:   60062b16e638
branch:      3.1

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