Re: [eigen] Re: Map and pointer freeing in unsupported/Eigen/src/SparseExtra/CholmodSupport.h

[Douglas Bates <bates@xxxxxxxxxxxxx>]

On Thu, Aug 18, 2011 at 10:08 AM, Gael Guennebaud
<gael.guennebaud@xxxxxxxxx> wrote:
>> I would say there won't be an easy and save fix, because one can't be sure
>> what kind of allocator CholMod uses on different platforms or in future
>> versions.
>
> Yes that's why swapping is probably not the right solution. Als,o
> swapping would only allow to get rid of the extra copy, not of the
> extra allocation.

I apologize for being late on rejoining this discussion.  I have been
travelling with often unreliable internet access.

Perhaps I didn't phrase my original request well.  I find that running
make check on the development branch of eigen ends with


99% tests passed, 2 tests failed out of 586

Total Test time (real) = 192.22 sec

The following tests FAILED:
	557 - sparse_llt_1 (Failed)
	559 - sparse_ldlt_1 (Failed)
Errors while running CTest

on an Ubuntu 11.04 amd_64 system where I have the CHOLMOD libraries
and headers installed.  Running these specific tests produces

$ ./sparse_llt_1
Initializing random number generator with seed 1314030613
Repeating each test 10 times
CHOLMOD warning: not positive definite
Test (sparse_llt<double,long int>(s,s)) failed in
/home/bates/devel/eigen/unsupported/test/sparse_llt.cpp (72)
    (m3*x).isApprox(b,test_precision<Scalar>()) && "LLT legacy:
cholmod solveInPlace"

Aborted
$ ./sparse_ldlt_1
Initializing random number generator with seed 1314030835
Repeating each test 10 times
CHOLMOD error: all methods failed
CHOLMOD error: argument missing
Segmentation fault

Even if these errors were fixed (and I will try to do so) these tests
are based on the code in CholmodSupportLegacy.h file and not on
CholmodSupport.h

The solve methods in both of those header files have the problem that
they map memory allocated by cholmod functions that is then freed by a
call to cholmod_free_dense or cholmod_free_sparse.

I think that trying to swap pointers will be unsuitable in the general
case because cholmod can be configured to use whatever malloc/free
functions are available, such as Matlab's alloc/free combination when
used in Matlab.  Without tweaking the CHOLMOD code at a low level I
think it will be necessary to create an eigen object, including
allocating the memory, and copy the contents of the cholmod result
then free the cholmod result.

My particular problem is that I still get errors showing up through
valgrind even when I do that.

Regarding whether this is only an issue for the supernodal
decomposition, I needed to switch to the cholmod interface instead of
using the SimplicialCholesky.h code for two reasons.  First, I need to
determine the Cholesky factor of many sparse, positive definite
symmetric matrices that have the same structure but different values
of the non-zeros so I want to keep the same symbolic factorization and
simply update it.  However, the matrix is generated as a product of
sparse matrices in which some of the "non-zeros" can end up with the
value of zero.  Because I am unable to control pruning of zeros in a
product of sparse matrices I need to detect this situation and perform
the symbolic factorization on the reduced matrix.

The second issue is performance.  Cholmod provides the colamd as well
as the amd methods for determining a fill-reducing permutation and
also provides ways of factoring A A' + k I directly from A.  Both of
these are important to me as they result in a non-trivial change in
execution time (as much as an order of magnitude).

I am convinced that the current solve methods in both CholmodSupport.h
and CholmodSupportLegacy.h are in error.   If you would like me to
document just that part, I will.   I am willing to work on this if I
can have some guidance about where my efforts would be most helpful.

>> A save way would be to do the triangular solve operations with own code (or
>> is there a CholMod function which works on views instead of allocating own
>> memory?).
>
> This would indeed be better. I quickly checked the user manual, and
> there are low levels solvers routines for the supernodal mode only,
> but that's all we need because otherwise me can directly use Eigen's
> own routines:
>
> cholmod_super_lsolve(...)
> cholmod_super_ltsolve(...)
>
> this is however a clear low priority for I because the benefits will
> be very marginal.
>
> cheers,
> gael
>
>
>>> viewAsEigen will return a MappedSparseMatrix with pointers to contents
>>> that will be freed in the subsequent line, unless I am missing
>>> something subtle.
>>>
>>> The comments about swapping lead me to believe that the dest argument
>>> has storage allocated on entry to the function, in which case the
>>> values could be copied from the structs&x_cd or&x_cs.  I must admit
>>> I get a little lost in the details of the sparse_solve_retval struct
>>> so if anyone could tell me the status of the Dest argument when the
>>> _solve method is called, I would appreciate it.  In particular, for a
>>> dense return value, has the storage for the data already been
>>> allocated?
>>
>> cholmod_solve allocates its own result vector (stored in x_cs). dest is a
>> parameter of _solve(). Currently the values _are_ copied (which technically
>> is not necessary).
>> If one assigns a solve_retval object to a vector x, then solve_retval's
>> evalTo function is called with parameter x.
>>
>>
>> Christoph
>>
>> --
>> ----------------------------------------------
>> Dipl.-Inf. Christoph Hertzberg
>> Cartesium 0.051
>> Universität Bremen
>> Enrique-Schmidt-Straße 5
>> 28359 Bremen
>>
>> Tel: (+49) 421-218-64252
>> ----------------------------------------------
>>
>>
>>
>
>
>