Allow the use of a combined function + Jacobian

[astrozot]

I need to solve a relatively complicated function (which depends on changing parameters) several thousands or millions of times, and I am trying to obtain a code that is as fast as possible. I am able to provide the Jacobian, but I would like to be able to compute the gradient and the function itself within the same function to avoid repeating, costly computations.

However, I do not see how one can use this package (and the sub package SimpleNonlinearSolve) with a setting where f and df share common computations. Ideally, I would like to be able to provide three separate functions: f, df, and fdf, with the last one used to compute at the same time the function and its Jacobian.

See also

https://discourse.julialang.org/t/specifying-analytical-jacobian-in-nonlinearproblem-a-la-optim-only-fg/115991

for a similar discussion.

P.S. For performance reasons, since my problem is in low dimensions I am using non-allocating functions and StaticArrays, so memoization is not an option.