ODE solvers and Autodiff

Hi,

I’ve seen multiple discussions in the forum but they are all 3-5 years old (at least the ones I spotted). Do you guys have any recommendations for ODE solvers (anything close to DifferentialEquations.jl in performance and broadness). Similarly, any autodiff packages that might be compatible?

I don’t really code Julia; I come from Python, but scipy’s solvers can be notably slow (mostly due to the Python side generating the equations). On the other hand, using something like jax requires that I rewrite most of the code for my use case, which I’m currently doing. I wanted to try Fortran and make a comparison, in terms of speed, and also how long it takes me to get this simple code working. But since Julia seems to have the upper hand in differential equations, I was thinking perhaps I should go there directly, even though I’ve been avoiding Julia for a few years now

I saw some recommendations here like FLINT and FOODIE, but it seems like those are unmaintained.

I’ve developed an automatic differentiation package in Fortran:

It can be used up to the 3d order derivatives and there is no limit in terms of independent variables. You can get gradient, jacobian…

I checked Beliavsky’s list of Fortran code on GitHub and I could only see one ODE solver library that mentioned anything which sounded like AD: FATODE. I have never used this library and it appears to be unmaintained.

Many of these libraries were developed well before automatic differentiation was popularized, so it’s not likely many are compatible with Fortran derived type based AD.

My suggestion would be to pick an AD library you like and write your own minimal ODE solver that works with it. The ODE solver part will be a lot easier than the AD part. This is roughly what I did in my Nerf blaster simulator, which uses a forward-mode AD library I wrote for another project, a new RK4 integrator (which is pretty deeply embedded and can not really be extracted into a standalone library), and a form of event detection.

Developing a Fortran ODE solver library that works with a particular AD library would be a good idea in my view. The various AD libraries I’ve seen seem to have broadly compatible operators, but incompatible ways of initializing variables. There would be benefit to standardizing around a particular way to initialize AD variables in Fortran to encourage interoperability of AD libraries.