finite_difference routines temporarily mutate inputs

[BobPortmann]

tl;dr The finite_difference methods for the gradient and jacobian (but not hessian) temporarily mutate the input x vector in order to compute the derivatives. I think this is unwise since there is no way to know inside these functions if the input f function holds a reference to x and the user will not suspect this will happen. A copy of x should be made to make these routines safer.

More details for those interested:

This recently happened to me and it was difficult to track down. It occurred because I called the curve_fit routine from LsqFit.jl:

function curve_fit(model::Function, xpts, ydata, p0; kwargs...)

and inputted the same array for ydata and p0. I know this seems unusual but it makes sense for the problem I was solving (I can provide details if you want). This caused the fit to fail because curve_fit makes a function that includes ydata as a closure that is used in finite_difference_jacobian!. I don't see how a user of curve_fit would suspect that somewhere under the hood an input array will be mutated (even temporarily).

I see that the Calculus.jl routines are going to be superseded by the FiniteDiff.jl and this concern affects routines there as well (for computing the gradient and hessian; I don't see a jacobian routine there).

Thanks for your consideration.

[johnmyleswhite]

I think the way to do this is to require users to provide us with a buffer that we can use instead of x. I don't think it's acceptable for us to always allocate memory given that finite differencing is already a potential bottleneck in many use cases. But we can provide an interface that allocates memory automatically.

I'm sorry you had such trouble tracking down your bug. Your problem is a great example of the reason that Rust's borrow checker is so helpful: Julia doesn't really let us express the fact that we borrow x as a mutable container. So there was no way for us to communicate that essential fact to you programmatically.

[BobPortmann]

I'm not sure what you have in mind for an interface that allocates memory? But it sounds promising.

[johnmyleswhite]

We can just use the current interface you're used to in FiniteDiff while always allocating a buffer. I don't see any reason to hesitate doing that if we make the contract more explicit. I'm less happy doing that same thing for Calculus since it will make some code that's currently correct slower.

[BobPortmann]

Sounds good if FiniteDiff.jl is going to be the way forward.

[ChrisRackauckas]

We ended up implementing non-allocating versions of this over in DiffEqDiffTools.jl