I am trying to keep my solver autodiff compatible to be able to take advantage of the ForwardDiff capabilities in Julia. I have a function f(du, u, t, p) where computing du depends on dt. I want to be able to pass the Julia solver’s time step, dt, into my function to compute du to appropriately update du without needing to make some approximated work around. I need to use a native Julia solver to subsequently use the ForwardDiff tools. Is there an easy way to do this?
Could you explain the context here?
It doesn’t sound like you have an ODE, since by definition an ODE can be written in the form \frac{du}{dt} = \operatorname{somefunction}(u, t), and this right-hand-side is supposed to be what your ODEFunction computes in DifferentialEquations.jl
I am solving a population balance model where the change in the number density function is a function of dt. The update to the number density function is as follows where k is dt, h is dx for the internal coordinate, G is growth rate, and phi is an applied flux limiter:
Although I can factor out k, I am still left with a k within my step evaluation for f and need to pass k (which is dt) to compute it. This is a finite volume method.
That’s not an ODE, it’s a difference (not differential!) equation. You could try DifferenceEquations.jl, or you could just write a loop. Loops (properly written) are fast in Julia.
The complexity of DifferentialEquations.jl all comes from approximating the solution of continuous-time equations, but you’ve already eliminated the continuity here.
1 Like
Ive already made my own explicit euler integration, but I want to be able to use the automatic differentiation capabilities to later solve parameter estimation and optimization problems. Does DifferenceEquations.jl support autodiff or is it autodiff compatible? Is my problem autodiff compatible?
Autodiff compatibility isn’t about problems, it’s about code. You can certainly implement an Euler integration in a way that is compatible with AD, and ForwardDiff.jl is especially easy. It all depends on understanding the coding style that the AD system requires (with ForwardDiff.jl it mainly just means writing “type-generic” code that works regardless of the numeric type, which is a good idea anyway).
(I would advise learning a little about how AD works; if you have lots of parameters, e.g. \gg 10, you might want to switch to reverse mode. See our matrix calculus course notes on forward- vs. reverse-mode differentiation algorithms, and see also my notes on reverse/adjoint differentiation of recurrence/difference equations.)
Alternatively, you could take a step back and convert equations back into an ODE. In the \Delta t = k \to 0 limit, your equations turn back into the ODE:
\frac{d f_n}{dt} =- \frac{1}{2h} \left[ \left(G_{n+1/2} f_{n+1} - G_{n-1/2} f_n \right) \phi_n \right. \\ + \left. \left(G_{n-3/2} f_{n-1} - G_{n-1/2} f_n \right) (2 - \phi_{n-1}) \right]
(which seems like a discretized diffusion-like PDE); presumably you could also get this from your original PDE by just discretizing space and not time to begin with. Then you hand this to DifferentialEquations.jl and let it discretize time for you — this is called the method of lines (MOL).
The advantage of using MOL is that (1) you get DifferentialEquations.jl’s AD code, but more importantly (2) you can take advantage of more sophisticated time-discretization methods, including adaptive and implicit methods. The latter are especially important for diffusion-like equations, which are often stiff and require very small timesteps for Euler or other explicit schemes.
1 Like
Thank you very much and for taking the time to provide this detailed response. I have already made my code type-generic so I should be able to just use ForwardDiff directly. I thought I needed to use a native Julia solver for implementing ForwardDiff
If your code is written in Julia, then it is a “native Julia solver”. There’s no magic in libraries except that they are (hopefully) written by experts — this is unlike Python, where high-performance libraries often include large chunks of low-level code written in other languages.