I’m trying to save the derivative values from an ODE using SavingCallback
from DiffEqCallbacks.jl. Why is manual_du
not the same as save_du.saveval
in this code? Am I misunderstanding when these values are saved, or is this is an issue?
using OrdinaryDiffEq, DiffEqCallbacks
function f!(du, u, p, t)
du[1] = 0.3u[1]
du[2] = 0.2u[1] + u[2]
return nothing
end
save_function = (u, t, integrator) -> (_u = similar(u); get_du!(_u, integrator); return _u)
save_du = SavedValues(Float64, Vector{Float64})
t0, t1 = 0.0, 1.0
saveat = t0:0.1:t1
cb = SavingCallback(save_function, save_du; saveat=saveat)
prob = ODEProblem(f!, [1.0,1.0], (t0,t1), callback=cb)
sol = solve(prob, Tsit5(), saveat=saveat)
manual_du = [zeros(2) for _ in eachindex(saveat)]
for j in eachindex(sol)
u = sol.u[j]
f!(manual_du[j], u, nothing, nothing)
end
manual_du[1] - [0.3 * 1.0, 0.2 * 1.0 + 1.0]
save_du.saveval[1] - [0.3 * 1.0, 0.2 * 1.0 + 1.0]
manual_du - save_du.saveval
julia> manual_du[1] - [0.3 * 1.0, 0.2 * 1.0 + 1.0]
2-element Vector{Float64}:
0.0
0.0
julia> save_du.saveval[1] - [0.3 * 1.0, 0.2 * 1.0 + 1.0]
2-element Vector{Float64}:
0.009350154485183504
0.13582658198324649
julia> manual_du - save_du.saveval
11-element Vector{Vector{Float64}}:
[-0.009350154485183504, -0.13582658198324649]
[-0.0002137942991308761, -0.0032172187877890224]
[-0.015279228394231514, -0.26106071340989256]
[-0.005577907363575996, -0.09867494524327491]
[-0.0320137675756863, -0.6656221051911233]
[-0.021712550094837635, -0.46684142492757896]
[-0.01109761449724933, -0.2469357592801935]
[-0.00015940524436997006, -0.0036734640814164443]
[-0.023582896950197274, -0.6267856404645369]
[-0.011968307129786249, -0.32916801361346115]
[0.0, 0.0]