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.