[CppAD] Sparse Hessians and Jacobians

[Michael Braun]

If we only want to support one format, this triplet format can work.  It is intuitive, and there are routines available to convert to and from CSC and CSR formats.  

I do have a question about how the Jacobian is returned.  If I understand the code correctly, CppAD allocates new memory for a new Jacobian/Hessian every time the function is called.  If I am repeatedly computing gradients and Hessians, this is an awful lot of dynamic memory allocation.  Wouldn't it be more efficient for the user to allocate memory for the output, and then pass in a pointer or reference for that array?

In that case, the call would be something like:

f.SparseJacobian(x, p, r, c, *out)

Also, will the SparseJacobian / SparseMatrix functions still need to allocate a full m-by-n matrix internally?  The value of a sparse matrix algorithm is not just that the matrix is return in a sparse format, but the algorithm itself exploits the sparseness in terms of its memory usage.


> 
> 
>   1. Sparse Hessians and Jacobians (Brad Bell)
> 
> 
> ----------------------------------------------------------------------
> 
> Message: 1
> Date: Tue, 08 May 2012 07:45:05 -0700
> From: Brad Bell <bradbell at seanet.com>
> To: CppAD at list.coin-or.org
> Subject: [CppAD] Sparse Hessians and Jacobians
> Message-ID: <4FA93171.5040207 at seanet.com>
> Content-Type: text/plain; charset=ISO-8859-1; format=flowed
> 
> I am considering adding sparse Jacobian and Hessians returns to CppAD as 
> per the sparsity discussion under
>     http://list.coin-or.org/pipermail/cppad/2012q1/thread.html
> 
> Reading over the discussion, it seems to me that the format used by 
> ADOL-C is reasonable; i.e. three vectors rind[], cind[], and values[] 
> all of the same length could be used to represent Sparse Hessians and 
> Jacobians (where the corresponding element types are size_t, size_t, and 
> Base).
> 
> For example, the sparse Jacobian driver
>     http://www.coin-or.org/CppAD/Doc/sparse_jacobian.xml
> could be extended to have the syntax
>     jac = f.SparseJacobian(x, p, r, c)
> where
>     x is the point where the Jacobian is being evaluated
>     p is the CppAD sparsity pattern (for the entire Jacobian)
>     r is the row indices
>     c is the column indices
>     jac is the Jacobian values at the corresponding row and column indices.
> 
> Note that the sparsity pattern p can be computed by either
>     http://www.coin-or.org/CppAD/Doc/forsparsejac.xml
>     http://www.coin-or.org/CppAD/Doc/revsparsejac.xml
> Then the users can decide which non-zero indices they want back and in 
> what order. For example, if the Jacobian is know to be symmetric, one 
> only needs to ask for the upper or lower triangle.
> 
> A similar approach be used for the sparse Hessians.
>     http://www.coin-or.org/CppAD/Doc/sparse_hessian.xml
> 
> 
> 
> ------------------------------
> 
> _______________________________________________
> CppAD mailing list
> CppAD at list.coin-or.org
> http://list.coin-or.org/mailman/listinfo/cppad
> 
> 
> End of CppAD Digest, Vol 56, Issue 1
> ************************************

Michael Braun 
MIT Sloan School of Management
braunm at mit.edu
-------------------------------------------------- 
View my research at
http://braunm.scripts.mit.edu/

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://list.coin-or.org/pipermail/cppad/attachments/20120509/018f7c96/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: smime.p7s
Type: application/pkcs7-signature
Size: 1844 bytes
Desc: not available
URL: <http://list.coin-or.org/pipermail/cppad/attachments/20120509/018f7c96/attachment.p7s>