Skip to content

Run Runic on the disco-optimizations files #4064

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
ChrisRackauckas-Claude wants to merge 1 commit into
base: master
Choose a base branch
from
+38 −34
Draft

Run Runic on the disco-optimizations files #4064

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 8 additions & 4 deletions lib/OrdinaryDiffEqCore/src/disco.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -54,8 +54,10 @@ function find_discontinuity(integrator)
if (disco_zero.out_low[j] * disco_zero.out_high[j] < zero(disco_zero.out_low[j]))
if (!addsteps_called)
addsteps_called = true
_ode_addsteps!(disco_zero.k, disco_zero.tprev, disco_zero.uprev, disco_zero.u,
disco_zero.dt, disco_zero.f, disco_zero.p, disco_zero.cache, false, true, false)
_ode_addsteps!(
disco_zero.k, disco_zero.tprev, disco_zero.uprev, disco_zero.u,
disco_zero.dt, disco_zero.f, disco_zero.p, disco_zero.cache, false, true, false
)
end
disco_zero.ind = j
disco_prob.tspan[2] = breakpointθ
Expand All @@ -74,8 +76,10 @@ function find_discontinuity(integrator)
if (disco_zero.out_low[1] * disco_zero.out_high[1] < zero(disco_zero.out_low[1]))
if (!addsteps_called)
addsteps_called = true
_ode_addsteps!(disco_zero.k, disco_zero.tprev, disco_zero.uprev, disco_zero.u,
disco_zero.dt, disco_zero.f, disco_zero.p, disco_zero.cache, false, true, false)
_ode_addsteps!(
disco_zero.k, disco_zero.tprev, disco_zero.uprev, disco_zero.u,
disco_zero.dt, disco_zero.f, disco_zero.p, disco_zero.cache, false, true, false
)
end
disco_prob.tspan[2] = breakpointθ
sol = solve(disco_prob)
Expand Down
2 changes: 1 addition & 1 deletion lib/OrdinaryDiffEqCore/src/integrators/controllers.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -571,7 +571,7 @@ for (accessor, default) in (
end
end

# Discontinuity-detection helpers
# Discontinuity-detection helpers
# Shared logic in every step_accept/reject_controller! below.

function handle_disco_accept!(integrator, controller_basic, t, nominal_new_dt)
Expand Down
58 changes: 29 additions & 29 deletions test/Integrators_I/disco_tests.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,23 +19,23 @@ prob = ODEProblem(f, u0, tspan)
condition(u, t, integrator) = u[1] - 1
cb = ContinuousCallback(condition, default_affect!; maybe_discontinuity = true)

sol_disco_radau = solve(prob, RadauIIA5(); callback = cb, reltol = 1e-6, controller = predictive_disco_controller(RadauIIA5()))
sol_disco_radau = solve(prob, RadauIIA5(); callback = cb, reltol = 1.0e-6, controller = predictive_disco_controller(RadauIIA5()))
# 314.500 μs (8909 allocations: 278.73 KiB)
sol_no_disco_radau = solve(prob, RadauIIA5(); callback = cb, reltol = 1e-6)
sol_no_disco_radau = solve(prob, RadauIIA5(); callback = cb, reltol = 1.0e-6)
# 376.125 μs (10628 allocations: 331.61 KiB)
@test sol_disco_radau.retcode == ReturnCode.Success
@test sol_disco_radau.stats.nreject <= sol_no_disco_radau.stats.nreject

sol_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1e-6, controller = PI_disco_controller(Rodas5P()))
sol_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1.0e-6, controller = PI_disco_controller(Rodas5P()))
# 444.500 μs (15596 allocations: 547.62 KiB)
sol_no_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1e-6)
sol_no_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1.0e-6)
# 620.709 μs (21830 allocations: 765.77 KiB)
@test sol_disco_rosenbrock.retcode == ReturnCode.Success
@test sol_disco_rosenbrock.stats.nreject <= sol_no_disco_rosenbrock.stats.nreject

sol_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1e-6, controller = PI_disco_controller(Tsit5()))
sol_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1.0e-6, controller = PI_disco_controller(Tsit5()))
# 55.625 μs (7576 allocations: 240.67 KiB)
sol_no_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1e-6)
sol_no_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1.0e-6)
# 55.083 μs (7570 allocations: 240.58 KiB)
@test sol_disco_tsit5.retcode == ReturnCode.Success
@test sol_disco_tsit5.stats.nreject <= sol_no_disco_tsit5.stats.nreject
Expand Down Expand Up @@ -64,16 +64,16 @@ condition2(u, t, integrator) = u[1] - 2
cb2 = ContinuousCallback(condition2, default_affect!; maybe_discontinuity = true)
cb = CallbackSet(cb1, cb2)

sol_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1e-6, controller = PI_disco_controller(Rodas5P()))
sol_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1.0e-6, controller = PI_disco_controller(Rodas5P()))
# 1.224 ms (42698 allocations: 1.46 MiB)
sol_no_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1e-6)
sol_no_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1.0e-6)
# 1.551 ms (52905 allocations: 1.80 MiB)
@test sol_disco_rosenbrock.retcode == ReturnCode.Success
@test sol_disco_rosenbrock.stats.nreject <= sol_no_disco_rosenbrock.stats.nreject

sol_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1e-6, controller = PI_disco_controller(Tsit5()))
sol_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1.0e-6, controller = PI_disco_controller(Tsit5()))
# 267.167 μs (34737 allocations: 1.07 MiB)
sol_no_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1e-6)
sol_no_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1.0e-6)
# 279.416 μs (39645 allocations: 1.22 MiB)
@test sol_disco_tsit5.retcode == ReturnCode.Success
@test sol_disco_tsit5.stats.nreject <= sol_no_disco_tsit5.stats.nreject
Expand Down Expand Up @@ -103,31 +103,31 @@ cb_multi_2 = ContinuousCallback(cond_multi_2, default_affect!; maybe_discontinui

cb_multi = CallbackSet(cb_multi_1, cb_multi_2)

sol_disco_radau = solve(prob_multi, RadauIIA5(); callback=cb_multi, reltol=1e-7, abstol=1e-9, controller = PI_disco_controller(RadauIIA5()))
sol_disco_radau = solve(prob_multi, RadauIIA5(); callback = cb_multi, reltol = 1.0e-7, abstol = 1.0e-9, controller = PI_disco_controller(RadauIIA5()))
# 148.875 μs (1564 allocations: 64.56 KiB)
sol_no_disco_radau = solve(prob_multi, RadauIIA5(); callback=cb_multi, reltol=1e-7, abstol=1e-9)
sol_no_disco_radau = solve(prob_multi, RadauIIA5(); callback = cb_multi, reltol = 1.0e-7, abstol = 1.0e-9)
# 129.875 μs (1073 allocations: 48.44 KiB)
@test sol_disco_radau.retcode == ReturnCode.Success
@test sol_disco_radau.stats.nreject <= sol_no_disco_radau.stats.nreject

sol_disco_rosenbrock = solve(prob_multi, Rodas5P(); callback=cb_multi, reltol=1e-7, abstol=1e-9, controller = PI_disco_controller(Rodas5P()))
sol_disco_rosenbrock = solve(prob_multi, Rodas5P(); callback = cb_multi, reltol = 1.0e-7, abstol = 1.0e-9, controller = PI_disco_controller(Rodas5P()))
# 209.000 μs (1474 allocations: 60.94 KiB)
sol_no_disco_rosenbrock = solve(prob_multi, Rodas5P(); callback=cb_multi, reltol=1e-7, abstol=1e-9)
sol_no_disco_rosenbrock = solve(prob_multi, Rodas5P(); callback = cb_multi, reltol = 1.0e-7, abstol = 1.0e-9)
# 237.250 μs (969 allocations: 45.44 KiB)
@test sol_disco_rosenbrock.retcode == ReturnCode.Success
@test sol_disco_rosenbrock.stats.nreject <= sol_no_disco_rosenbrock.stats.nreject

sol_disco_tsit5 = solve(prob_multi, Tsit5(); callback=cb_multi, reltol=1e-7, abstol=1e-9, controller = PI_disco_controller(Tsit5()))
sol_disco_tsit5 = solve(prob_multi, Tsit5(); callback = cb_multi, reltol = 1.0e-7, abstol = 1.0e-9, controller = PI_disco_controller(Tsit5()))
# 105.250 μs (1592 allocations: 60.34 KiB)
sol_no_disco_tsit5 = solve(prob_multi, Tsit5(); callback=cb_multi, reltol = 1e-7, abstol = 1e-9)
sol_no_disco_tsit5 = solve(prob_multi, Tsit5(); callback = cb_multi, reltol = 1.0e-7, abstol = 1.0e-9)
# 88.000 μs (1133 allocations: 46.84 KiB)
@test sol_disco_tsit5.retcode == ReturnCode.Success
@test sol_disco_tsit5.stats.nreject <= sol_no_disco_tsit5.stats.nreject

#TEST 4: STIFF MULTI-COMPONENT DISCONTINUITY
function f_stiff_disc!(du, u, p, t)
λ = p[1]
if u[1] < 0.5
return if u[1] < 0.5
du[1] = -λ * u[1] + λ * exp(-t)
du[2] = -λ * (u[2] - u[1])
du[3] = -(u[3] - u[1]^2)
Expand All @@ -145,9 +145,9 @@ prob_stiff = ODEProblem(f_stiff_disc!, u0_stiff, tspan_stiff, [500.0])
cond_stiff(u, t, integrator) = u[1] - 0.5
cb_stiff = ContinuousCallback(cond_stiff, default_affect!; maybe_discontinuity = true)

sol_disco_rosenbrock = solve(prob_stiff, Rodas5P(); callback=cb_stiff, reltol=1e-9, abstol=1e-11, controller = PI_disco_controller(Rodas5P()))
sol_disco_rosenbrock = solve(prob_stiff, Rodas5P(); callback = cb_stiff, reltol = 1.0e-9, abstol = 1.0e-11, controller = PI_disco_controller(Rodas5P()))
# 364.833 μs (2030 allocations: 85.34 KiB)
sol_no_disco_rosenbrock = solve(prob_stiff, Rodas5P(); callback=cb_stiff, reltol=1e-9, abstol=1e-11)
sol_no_disco_rosenbrock = solve(prob_stiff, Rodas5P(); callback = cb_stiff, reltol = 1.0e-9, abstol = 1.0e-11)
# 332.750 μs (1765 allocations: 79.78 KiB)
@test sol_disco_rosenbrock.retcode == ReturnCode.Success
@test sol_disco_rosenbrock.stats.nreject <= sol_no_disco_rosenbrock.stats.nreject
Expand Down Expand Up @@ -178,9 +178,9 @@ prob_dae = ODEProblem(f_dae_func, u0_dae, tspan_dae)
cond_dae(u, t, integrator) = u[1] - 0.5
cb_dae = ContinuousCallback(cond_dae, default_affect!; maybe_discontinuity = true)

sol_disco_rosenbrock = solve(prob_dae, Rodas5P(); callback=cb_dae, reltol=1e-8, abstol=1e-10, controller = PI_disco_controller(Rodas5P()))
sol_disco_rosenbrock = solve(prob_dae, Rodas5P(); callback = cb_dae, reltol = 1.0e-8, abstol = 1.0e-10, controller = PI_disco_controller(Rodas5P()))
# 127.083 μs (831 allocations: 35.53 KiB)
sol_no_disco_rosenbrock = solve(prob_dae, Rodas5P(); callback=cb_dae, reltol=1e-8, abstol=1e-10)
sol_no_disco_rosenbrock = solve(prob_dae, Rodas5P(); callback = cb_dae, reltol = 1.0e-8, abstol = 1.0e-10)
# 121.625 μs (625 allocations: 31.69 KiB)
@test sol_disco_rosenbrock.retcode == ReturnCode.Success
@test sol_disco_rosenbrock.stats.nreject <= sol_no_disco_rosenbrock.stats.nreject
Expand All @@ -203,16 +203,16 @@ end
default_affect_vec!(integrator, idx) = nothing
cb_vec = VectorContinuousCallback(condition_vec!, default_affect_vec!, 2; maybe_discontinuity = true)

sol_disco_rosenbrock = solve(prob_vec, Rodas5P(); callback = cb_vec, reltol=1e-7, abstol=1e-9, controller = PI_disco_controller(Rodas5P()))
sol_disco_rosenbrock = solve(prob_vec, Rodas5P(); callback = cb_vec, reltol = 1.0e-7, abstol = 1.0e-9, controller = PI_disco_controller(Rodas5P()))
# 191.250 μs (1212 allocations: 53.62 KiB)
sol_no_disco_rosenbrock = solve(prob_vec, Rodas5P(); callback = cb_vec, reltol=1e-7, abstol=1e-9)
sol_no_disco_rosenbrock = solve(prob_vec, Rodas5P(); callback = cb_vec, reltol = 1.0e-7, abstol = 1.0e-9)
# 202.000 μs (939 allocations: 49.47 KiB)
@test sol_disco_rosenbrock.retcode == ReturnCode.Success
@test sol_disco_rosenbrock.stats.nreject <= sol_no_disco_rosenbrock.stats.nreject

sol_disco_tsit5 = solve(prob_vec, Tsit5(); callback = cb_vec, reltol=1e-7, abstol=1e-9, controller = PI_disco_controller(Tsit5()))
sol_disco_tsit5 = solve(prob_vec, Tsit5(); callback = cb_vec, reltol = 1.0e-7, abstol = 1.0e-9, controller = PI_disco_controller(Tsit5()))
# 94.459 μs (1626 allocations: 67.73 KiB)
sol_no_disco_tsit5 = solve(prob_vec, Tsit5(); callback = cb_vec, reltol=1e-7, abstol=1e-9)
sol_no_disco_tsit5 = solve(prob_vec, Tsit5(); callback = cb_vec, reltol = 1.0e-7, abstol = 1.0e-9)
# 93.084 μs (1356 allocations: 63.04 KiB)
@test sol_disco_tsit5.retcode == ReturnCode.Success
@test sol_disco_tsit5.stats.nreject <= sol_no_disco_tsit5.stats.nreject
Expand All @@ -235,16 +235,16 @@ prob = ODEProblem(f!, u, tspan)
cond(u, t, integrator) = u[2]
cb = ContinuousCallback(cond, default_affect!; maybe_discontinuity = true)

sol_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1e-8, abstol = 1e-10, controller = PI_disco_controller(Rodas5P()))
sol_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1.0e-8, abstol = 1.0e-10, controller = PI_disco_controller(Rodas5P()))
# 169.208 μs (1190 allocations: 51.78 KiB)
sol_no_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1e-8, abstol = 1e-10)
sol_no_disco_rosenbrock = solve(prob, Rodas5P(); callback = cb, reltol = 1.0e-8, abstol = 1.0e-10)
# 184.792 μs (877 allocations: 45.72 KiB)
@test sol_disco_rosenbrock.retcode == ReturnCode.Success
@test sol_disco_rosenbrock.stats.nreject <= sol_no_disco_rosenbrock.stats.nreject

sol_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1e-8, abstol = 1e-10, controller = PI_disco_controller(Tsit5()))
sol_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1.0e-8, abstol = 1.0e-10, controller = PI_disco_controller(Tsit5()))
# 59.625 μs (1248 allocations: 53.79 KiB)
sol_no_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1e-8, abstol = 1e-10)
sol_no_disco_tsit5 = solve(prob, Tsit5(); callback = cb, reltol = 1.0e-8, abstol = 1.0e-10)
# 54.041 μs (1155 allocations: 54.55 KiB)
@test sol_disco_tsit5.retcode == ReturnCode.Success
@test sol_disco_tsit5.stats.nreject <= sol_no_disco_tsit5.stats.nreject
Loading