search

kaskr / adcomp

AD computation with Template Model Builder (TMB)
Other
178 stars 80 forks source link

Misuse of `M_chm_sparse_to_SEXP` by `tmb_invQ_tril_halfdiag` #379

Closed jaganmn closed 1 year ago

jaganmn commented 1 year ago

Forthcoming improvements to subscript methods in Matrix have uncovered that you are misusing M_chm_sparse_to_SEXP in tmb_invQ_tril_halfdiag:

https://github.com/kaskr/adcomp/blob/bb2e1432638a8af60aee9b6f1868cc38cb5731dc/TMB/src/solve_subset.c#L318

Here, the fifth argument is used to set the diag slot of the resulting dtCMatrix and so must be "N" (non-unit triangular) or "U" (unit triangular) but never "". Please correct, then test your patch against Matrix-devel:

$ svn checkout svn://svn.r-forge.r-project.org/svnroot/matrix/pkg/Matrix
$ R CMD build Matrix
$ R CMD INSTALL Matrix_1.5-5.tar.gz

Trouble with the invalid ihessian (having diag = "") seems to start inside of this lookup call and result in a gradient with NA in the "simple" example (and for several packages depending on TMB):

https://github.com/kaskr/adcomp/blob/bb2e1432638a8af60aee9b6f1868cc38cb5731dc/TMB/R/TMB.R#L767

jaganmn commented 1 year ago

BTW, going forward the subscript method for CsparseMatrix will not construct an intermediate dgCMatrix for <dsCMatrix>[i, i] when i is strictly increasing or strictly decreasing. So lookup can in principle be made faster if you know that r is monotone.

jaganmn commented 1 year ago

Somehow, I didn't notice that the bug appears twice, i.e., also here:

https://github.com/kaskr/adcomp/blob/5de4200fd8ca0c74dfb297e228c31dd6ce1e8cc2/TMB/src/solve_subset.c#L297