When I use `DifferentialEquations.jl`

to solve a set of stiff coupled ODEs, my function is something like

```
function f!(du, u, par, t)
V = g(u, par) # calculate V based on u and some other parameters
du .= ??? # some calculations to update du based on u, par and V
return nothing
```

`du`

is updated based on the dependent variable `u`

, some variables and parameters stored inside `par`

, and also an intermediate variable `V`

which is calculated by calling some other functions that use `u`

and `par`

. The problem is that calculating `V=g(u,par)`

takes quite a long time and may require calling some external library. Since many stiff ODE solvers like the recommended `FBDF`

requires calculating the Jacobian, it is not quite feasible to calculate the Jacobian of `f!`

due to `V`

. However, for my specific problem I think it might be okay that within each time step `V`

is just treated as a variable independent of the current `u`

, and just approximates `V`

using the `u`

from the previous step. Since the system should evolve to a steady state, I suspect this method is good enough for the required accuracy.

So mu question is, is it possible to just compute `f!`

(and possibly the Jacobian) by ‘ignoring’ `V=g(u,par)`

? Or as I suggested, is it possible that at the current time step `k`

, calculate `V=g(u_{k-1},par)`

where `u_{k-1}`

is the solution at the previous time step? I’m wondering if I can have a callback to just compute `V`

after each time step and store it in `par`

, so that at the next time step the update of `du`

does not require `u`

in the current step, and so possibly the Jacobian calculation won’t take ages.