Help debugging differential equation

[ChrisRackauckas]

@Ge0rges

I would also appreciate general debugging tips, as stepping through 1000 iterations of a model call is impractical and so is manual logging. What are some strategies used?

How to debug Julia simulation codes (ODEs, optimization, etc.!)

That has a lot of tactics. Note that the easiest way to step through 1000 iterations is to use the integrator interface.

https://diffeq.sciml.ai/stable/basics/integrator/

So it’s just a few lines of code

integ = init(prob, alg, ...)
for i in 1:1000
  step!(integ)
end

@show integ.u
@show integ.t
@show integ.dt

and you can start to play with it to see what’s diverging first and all. Check integ.k[1] to look at the derivatives: see if those make sense.

@Ge0rges

One thing I am aware of is that the model never seems to settle on zero values. Often it gives approximations such as 1e-300.

Yeah this is a funny thing. LU factorizations do not guarantee pure zeros, so you can get phantom values that revive. I’m going to write a blog post on this “soon”, but for now you can refer to Dead species coming back alive · Issue #1651 · SciML/OrdinaryDiffEq.jl · GitHub which describes techniques to use in such a case.