Gradient computations

[astamm]

The issue to brainstorm on numerical derivatives vs automatic differentiation and whether or not we should do something about it in nloptr.

Just attended a talk on pnd which might be worth considering.

[aadler]

While I haven't been involved in many years, I was a bit more active in the Stan community early on, and AD adds a level of complexity which we may not want to consider. Personally, my workflow when things get really hairy is to use Maxima to calculate the derivatives in closed form when possible, instead of numeric differentiation. I think that the primary focus should be to implement as much as nlopt as we can accurately, perhaps exposing more of its internals, but that's just my opinion.

[aadler]

For that matter, we may wish to consider the five-point stencil, as I suggested a while ago, which increases the accuracy to fourth order from second order. It will be a breaking change since the gard/jacobian functions will return different—ostensibly more accurate—values.

If you want even more precision, we can implement Fornberg's algorithm. That would allow calculation of the derivative to any level of precision we want, so long as we select enough points. Although, to be fair, the Fornberg algorithm returns the weights of the five (and three) point stencils for five (and three) points symmetrically around the point in question—as it should, as those are the best linear estimates for those number of points.

I've implemented the algorithm in my personal R package (not on CRAN/private Github repo), so it's not hard, but it is a triply nested loop, for which I use a three dimensional array and just read off the last calculated value. I think it'soverkill, but it's something to consider and significantly simpler than trying to implement automatic differentiation, to the best of my understanding.