OrdinaryDiffEq: unexpected results in DEDataArray

You’re relying on behavior that cannot and should not be guaranteed. Basically, you’re relying on the assumption that the final f call of the Runge-Kutta method is at t + dt. That is not the case in many methods (for example, it’s not required in explicit or implicit Runge-Kutta methods). So it’s returning the t that corresponds to that final f call, which then is the “bug”.

Why is RK4 of all algorithms doing that? Because, if you look at the documentation, this adaptive RK4 is:

  • RK4 - The canonical Runge-Kutta Order 4 method. Uses a defect control for adaptive stepping using maximum error over the whole interval.

Defect control is done via:

Reference: http://faculty.smu.edu/shampine/residuals.pdf

You can see two f calls in there. This algorithm has many nice properties of course since it can bound the interval error, and adaptive=false then returns you to the simple RK4 that most people know. But that right there is the reason for why the t is not what you think it is.

In general, this kind of thing is the reason why DEDataArray is being deprecated, because relying on hidden values implicit in the array type also relies on other assumptions about the algorithms which are not necessarily correct. Instead, for all of these use cases with observables and such, ModelingToolkit.jl is a much better solution that we are now focusing on.

https://mtk.sciml.ai/dev/

1 Like