OrdinaryDiffEq.jl icon indicating copy to clipboard operation
OrdinaryDiffEq.jl copied to clipboard

Step size after tstop can be too small

Open andreasnoack opened this issue 3 years ago • 3 comments

This can cause unnecessary abortion. E.g.

julia> f = (u,p,t) -> -p.a*u
#3 (generic function with 1 method)

julia> prob = ODEProblem(f, [10.0], (0.0, 1.0), (a=-100.0,))
ODEProblem with uType Vector{Float64} and tType Float64. In-place: false
timespan: (0.0, 1.0)
u0: 1-element Vector{Float64}:
 10.0

julia> sol1 = solve(prob, Vern7());

julia> solve(prob, Vern7(), tstops=[sol1.t[2] + eps()])
┌ Warning: dt <= dtmin. Aborting. There is either an error in your model specification or the true solution is unstable.
└ @ SciMLBase ~/.julia/packages/SciMLBase/BoNUy/src/integrator_interface.jl:345
retcode: DtLessThanMin
Interpolation: specialized 7th order lazy interpolation
t: 2-element Vector{Float64}:
 0.0
 0.01
u: 2-element Vector{Vector{Float64}}:
 [10.0]
 [27.18282103271114]

I suppose that the step proposal after a tstop should always be at least the size of the initial step size proposal. Even in the cases where the current behavior doesn't cause an error, the step size might make the solving slower than needed.

It's not super easy to hit a case like this for Float64 inputs. I think this is because of the snapping-to-tstop in https://github.com/SciML/OrdinaryDiffEq.jl/blob/e06e744b2e290170154890dd3ecb075e77b2fe5f/src/integrators/integrator_utils.jl#L197-L206. However, before https://github.com/SciML/DiffEqBase.jl/pull/735 it could happen relatively frequently when passing ForwardDiff.Duals.

andreasnoack avatar Mar 14 '22 07:03 andreasnoack

abs(ttmp - tstop) < 100eps(float(max(integrator.t,tstop)/oneunit(integrator.t)))*oneunit(integrator.t)

Now from this issue, it seems to make more sense to do 100dtmin.

ChrisRackauckas avatar Mar 14 '22 14:03 ChrisRackauckas

Yes. That would leave some room when starting again after a tstop. However, isn't there still a risk that an artificially small step before a tstop makes the first step after the tstop (even if it's larger than 100dtmin) such that it would be better to choose the step size right after a tstop as if it was the initial step?

andreasnoack avatar Mar 15 '22 19:03 andreasnoack

such that it would be better to choose the step size right after a tstop as if it was the initial step?

No, first of all it's not like the initial stepsize algorithm is a robust thing. It's useful mostly just to minimize the number of times it needs to reject at the beginning. But secondly, that wouldn't solve this problem. It will still need to be the case that the maximum stepsize you can take is the distance between two tstops, because by definition you cannot step over the tstop. So if a user provides two tstops which are closer than dtmin, then it will error. We could do something like skip doing a real step and just update t in this case, but for the most part this is a user input error since it's directly asking for a dt<dtmin.

ChrisRackauckas avatar Mar 27 '22 21:03 ChrisRackauckas

So if a user provides two tstops which are closer than dtmin, then it will error.

Notice that this is not what happens here. The user only provides a single tstop value. The other one is chosen by the solver, so it's impossible to reason about for the user. I'm again trying to help a user who is hitting this issue.

andreasnoack avatar Oct 04 '23 06:10 andreasnoack