[CppAD] Sparse Hessians and Jacobians

[bradbell at seanet.com]

> 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)

Yes, this would be more efficient, but it would be a marked change from
the other cases where the value is returned. Perhaps one should have the
option to pass in the allocated memory as you suggest above.

Perhaps all the CppAD routines should have an interface that uses input
and output iterators (for arbitrary containers), instead just simple
vectors. Harder for the beginning user, but more versatile, more
efficient, and expected by the C++ expert.

>
> 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.

The sparse Jacobians and Hessians only allocate space corresponding to one
row (forward mode) or one column (reverse mode). Then (using the sparsity)
multiple rows or columns are computed as once and then separated.

On the other hand, the one row or one column is for every intermediate
variable in the operation sequence; see size_var
      http://www.coin-or.org/CppAD/Doc/seq_property.htm#size_var
not just the independent and dependent variables.

>
>
>>
>>
>>   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/
>
> _______________________________________________
> CppAD mailing list
> CppAD at list.coin-or.org
> http://list.coin-or.org/mailman/listinfo/cppad
>