diff --git a/src/callbacks.jl b/src/callbacks.jl index a84d11e4c..bd9c09978 100644 --- a/src/callbacks.jl +++ b/src/callbacks.jl @@ -140,7 +140,7 @@ end integrator, callbacks[$i], $i) - if event_occurred2 && (tmin2 < tmin || !event_occurred) + if event_occurred2 && (!event_occurred || integrator.tdir * tmin2 < integrator.tdir * tmin) tmin = tmin2 upcrossing = upcrossing2 event_occurred = true diff --git a/test/downstream/community_callback_tests.jl b/test/downstream/community_callback_tests.jl index fa68ffb7d..673e651e4 100644 --- a/test/downstream/community_callback_tests.jl +++ b/test/downstream/community_callback_tests.jl @@ -230,3 +230,34 @@ sol = solve(prob, DFBDF()) # test that the callback flipping p caused u[2] to get flipped. first_t = findfirst(isequal(0.5), sol.t) @test sol.u[first_t][2] == -sol.u[first_t + 1][2] + +# https://github.com/SciML/DiffEqBase.jl/issues/1231 +@testset "Successive callbacks in same integration step" begin + cb = ContinuousCallback( + (u, t, integrator) -> t - 0e-8, + (integrator) -> push!(record, 0) + ) + + vcb = VectorContinuousCallback( + (out, u, t, integrator) -> out .= (t - 1e-8, t - 2e-8), + (integrator, event_index) -> push!(record, event_index), + 2 + ) + + f(u, p, t) = 1.0 + u0 = 0.0 + + # Forward propagation with successive events + record = [] + tspan = (-1.0, 1.0) + prob = ODEProblem(f, u0, tspan) + sol = solve(prob, Tsit5(), dt = 2.0, callback = CallbackSet(cb, vcb)) + @test record == [0, 1, 2] + + # Backward propagation with successive events + record = [] + tspan = (1.0, -1.0) + prob = ODEProblem(f, u0, tspan) + sol = solve(prob, Tsit5(), dt = 2.0, callback = CallbackSet(cb, vcb)) + @test record == [2, 1, 0] +end