Merge pull request #490 from SciML/kg/cbfix3

Callback hotfix

Author: ChrisRackauckas

src/callbacks.jl

@@ -765,12 +765,12 @@ end
 
 function apply_callback!(integrator,callback::Union{ContinuousCallback,VectorContinuousCallback},cb_time,prev_sign,event_idx)
 
-  change_t_via_interpolation!(integrator,integrator.tprev+cb_time)
-
   if integrator.opts.adaptive
-    set_proposed_dt!(integrator, max(integrator.opts.dtmin+eps(integrator.dt), callback.dtrelax * integrator.dt))
+    set_proposed_dt!(integrator, max(nextfloat(integrator.opts.dtmin), callback.dtrelax * integrator.dt))
   end
 
+  change_t_via_interpolation!(integrator,integrator.tprev+cb_time)
+
   # handle saveat
   _, savedexactly = savevalues!(integrator)
   saved_in_cb = true

test/downstream/Project.toml

@@ -1,7 +1,9 @@
 [deps]
+Calculus = "49dc2e85-a5d0-5ad3-a950-438e2897f1b9"
 DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
 DiffEqCallbacks = "459566f4-90b8-5000-8ac3-15dfb0a30def"
 DiffEqProblemLibrary = "a077e3f3-b75c-5d7f-a0c6-6bc4c8ec64a9"
+ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 SteadyStateDiffEq = "9672c7b4-1e72-59bd-8a11-6ac3964bc41f"
 StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0"

test/downstream/ad_tests.jl

@@ -0,0 +1,154 @@
+using Test
+using OrdinaryDiffEq, Calculus, ForwardDiff
+
+function f(du,u,p,t)
+  du[1] = -p[1]
+  du[2] = p[2]
+end
+
+for x in 0:0.001:5
+  called = false
+  function test_f(p)
+    cb = ContinuousCallback((u,t,i) -> u[1], (integrator)->(called=true;integrator.p[2]=zero(integrator.p[2])))
+    prob = ODEProblem(f,eltype(p).([1.0,0.0]),eltype(p).((0.0,1.0)),copy(p))
+    integrator = init(prob,Tsit5(),abstol=1e-14,reltol=1e-14,callback=cb)
+    step!(integrator)
+    solve!(integrator).u[end]
+  end
+  p = [2.0, x]
+  called = false
+  findiff = Calculus.finite_difference_jacobian(test_f,p)
+  @test called
+  called = false
+  fordiff = ForwardDiff.jacobian(test_f,p)
+  @test called
+  @test findiff ≈ fordiff
+end
+
+function f2(du,u,p,t)
+  du[1] = -u[2]
+  du[2] = p[2]
+end
+
+for x in 2.1:0.001:5
+  called = false
+  function test_f2(p)
+    cb = ContinuousCallback((u,t,i) -> u[1], (integrator)->(called=true;integrator.p[2]=zero(integrator.p[2])))
+    prob = ODEProblem(f2,eltype(p).([1.0,0.0]),eltype(p).((0.0,1.0)),copy(p))
+    integrator = init(prob,Tsit5(),abstol=1e-12,reltol=1e-12,callback=cb)
+    step!(integrator)
+    solve!(integrator).u[end]
+  end
+  p = [2.0, x]
+  findiff = Calculus.finite_difference_jacobian(test_f2,p)
+  @test called
+  called = false
+  fordiff = ForwardDiff.jacobian(test_f2,p)
+  @test called
+  @test findiff ≈ fordiff
+end
+
+#=
+#x = 2.0 is an interesting case
+
+x = 2.0
+
+function test_f2(p)
+  cb = ContinuousCallback((u,t,i) -> u[1], (integrator)->(@show(x,integrator.t);called=true;integrator.p[2]=zero(integrator.p[2])))
+  prob = ODEProblem(f2,eltype(p).([1.0,0.0]),eltype(p).((0.0,1.0)),copy(p))
+  integrator = init(prob,Tsit5(),abstol=1e-12,reltol=1e-12,callback=cb)
+  step!(integrator)
+  solve!(integrator).u[end]
+end
+
+p = [2.0, x]
+findiff = Calculus.finite_difference_jacobian(test_f2,p)
+@test called
+called = false
+fordiff = ForwardDiff.jacobian(test_f2,p)
+@test called
+
+# At that value, it shouldn't be called, but a small perturbation will make it called, so finite difference is wrong!
+=#
+
+for x in 1.0:0.001:2.5
+  function lotka_volterra(du,u,p,t)
+    x, y = u
+    α, β, δ, γ = p
+    du[1] = dx = α*x - β*x*y
+    du[2] = dy = -δ*y + γ*x*y
+  end
+  u0 = [1.0,1.0]
+  tspan = (0.0,10.0)
+  p = [x,1.0,3.0,1.0]
+  prob = ODEProblem(lotka_volterra,u0,tspan,p)
+  sol = solve(prob,Tsit5())
+
+  called=false
+  function test_lotka(p)
+    cb = ContinuousCallback((u,t,i) -> u[1]-2.5, (integrator)->(called=true;integrator.p[4]=1.5))
+    prob = ODEProblem(lotka_volterra,eltype(p).([1.0,1.0]),eltype(p).((0.0,10.0)),copy(p))
+    integrator = init(prob,Tsit5(),abstol=1e-12,reltol=1e-12,callback=cb)
+    step!(integrator)
+    solve!(integrator).u[end]
+  end
+
+  findiff = Calculus.finite_difference_jacobian(test_lotka,p)
+  @test called
+  called = false
+  fordiff = ForwardDiff.jacobian(test_lotka,p)
+  @test called
+  @test findiff ≈ fordiff
+end
+
+# Gradients and Hessians
+
+function myobj(θ)
+  f(u,p,t) = -θ[1]*u
+  u0, _ = promote(10.0, θ[1])
+  prob = ODEProblem(f, u0, (0.0, 1.0))
+  sol = solve(prob, Tsit5())
+  diff = sol.u - 10*exp.(-sol.t)
+  return diff'diff
+end
+
+ForwardDiff.gradient(myobj, [1.0])
+ForwardDiff.hessian(myobj, [1.0])
+
+function myobj2(θ)
+  f(du,u,p,t) = (du[1]=-θ[1]*u[1])
+  u0, _ = promote(10.0, θ[1])
+  prob = ODEProblem(f, [u0], (0.0, 1.0))
+  sol = solve(prob, Tsit5())
+  diff = sol[:,1] .- 10 .*exp.(-sol.t)
+  return diff'diff
+end
+
+ForwardDiff.gradient(myobj2, [1.0])
+ForwardDiff.hessian(myobj2, [1.0])
+
+function myobj3(θ)
+  f(u,p,t) = -θ[1]*u
+  u0, _ = promote(10.0, θ[1])
+  tspan_end, _ = promote(1.0, θ[1])
+  prob = ODEProblem(f, u0, (0.0, tspan_end))
+  sol = solve(prob, Tsit5())
+  diff = sol.u - 10*exp.(-sol.t)
+  return diff'diff
+end
+
+ForwardDiff.gradient(myobj3, [1.0])
+ForwardDiff.hessian(myobj3, [1.0])
+
+function myobj4(θ)
+  f(du,u,p,t) = (du[1] = -θ[1]*u[1])
+  u0, _ = promote(10.0, θ[1])
+  tspan_end, _ = promote(1.0, θ[1])
+  prob = ODEProblem(f, [u0], (0.0, tspan_end))
+  sol = solve(prob, Tsit5())
+  diff = sol[:,1] .- 10 .* exp.(-sol.t)
+  return diff'diff
+end
+
+ForwardDiff.gradient(myobj4, [1.0])
+ForwardDiff.hessian(myobj4, [1.0])

test/runtests.jl

@@ -47,6 +47,7 @@ if !is_APPVEYOR && GROUP == "Downstream"
     @time @safetestset "DE stats" begin include("downstream/destats_tests.jl") end
     @time @safetestset "DEDataArray" begin include("downstream/data_array_regression_tests.jl") end
     @time @safetestset "Concrete_solve Tests" begin include("downstream/concrete_solve_tests.jl") end
+    @time @safetestset "AD Tests" begin include("downstream/ad_tests.jl") end
 end
 
 if !is_APPVEYOR && GROUP == "GPU"