StackOverFlowError with two usages of save_func

[MartinOtter]

Here is an ODE test problem:

Test_ode_save_func_callback.zip

with the essential statements:

...
prob1 = ODEProblem(simplePendulum1!, x1₀, tspan, model1)
sol1 = solve(prob1, Tsit5(), reltol = 1e-6, saveat=t_inc, callback=cb1)
sol2 = solve(prob1, Tsit5(), reltol = 1e-6, saveat=t_inc, callback=cb1)

When solve is called once and then at the second call without the callback, everything is fine. However, when the above two calls of solve are called with the callback cb1, then the following error occurs:

ERROR: LoadError: StackOverflowError:
Stacktrace:
 [1] ismutable(::Type) at HOME\.julia\packages\ArrayInterface\PQXWA\src\ArrayInterface.jl:15 (repeats 13032 times)
 [2] copyat_or_push! at HOME\.julia\packages\RecursiveArrayTools\qHvJz\src\utils.jl:59 [inlined]
 [3] (::DiffEqCallbacks.SavingAffect{typeof(Main.Test_ode_save_func_callback.simplePendulum2),Float64,Main.Test_ode_save_func_callback.Model,DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Array{Float64,1}})(::OrdinaryDiffEq.ODEIntegrator{OrdinaryDiffEq.Tsit5,true,Array{Float64,1},Float64,Main.Test_ode_save_func_callback.Model{Float64},Float64,Float64,Float64,Array{Array{Float64,1},1},DiffEqBase.ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},DiffEqBase.ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Main.Test_ode_save_func_callback.Model{Float64},DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},OrdinaryDiffEq.Tsit5,OrdinaryDiffEq.InterpolationData{DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}}},DiffEqBase.DEStats},DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),DiffEqBase.CallbackSet{Tuple{},Tuple{DiffEqBase.DiscreteCallback{DiffEqCallbacks.var"#30#31",DiffEqCallbacks.SavingAffect{typeof(Main.Test_ode_save_func_callback.simplePendulum2),Float64,Main.Test_ode_save_func_callback.Model,DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Array{Float64,1}},typeof(DiffEqCallbacks.saving_initialize)}}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Float64,Array{Float64,1}},Array{Float64,1},Float64,Nothing}, ::Bool) at HOME\.julia\dev\DiffEqCallbacks\src\saving.jl:65
 [4] SavingAffect at HOME\.julia\dev\DiffEqCallbacks\src\saving.jl:45 [inlined]
 [5] saving_initialize at HOME\.julia\dev\DiffEqCallbacks\src\saving.jl:85 [inlined]
 [6] initialize! at HOME\.julia\dev\DiffEqBase\src\callbacks.jl:273 [inlined]
 [7] initialize! at HOME\.julia\dev\DiffEqBase\src\callbacks.jl:263 [inlined]
 [8] initialize_callbacks!(::OrdinaryDiffEq.ODEIntegrator{OrdinaryDiffEq.Tsit5,true,Array{Float64,1},Float64,Main.Test_ode_save_func_callback.Model{Float64},Float64,Float64,Float64,Array{Array{Float64,1},1},DiffEqBase.ODESolution{Float64,2,Array{Array{Float64,1},1},Nothing,Nothing,Array{Float64,1},Array{Array{Array{Float64,1},1},1},DiffEqBase.ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Main.Test_ode_save_func_callback.Model{Float64},DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem},OrdinaryDiffEq.Tsit5,OrdinaryDiffEq.InterpolationData{DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Array{Float64,1},1},Array{Float64,1},Array{Array{Array{Float64,1},1},1},OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}}},DiffEqBase.DEStats},DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},OrdinaryDiffEq.Tsit5Cache{Array{Float64,1},Array{Float64,1},Array{Float64,1},OrdinaryDiffEq.Tsit5ConstantCache{Float64,Float64}},OrdinaryDiffEq.DEOptions{Float64,Float64,Float64,Float64,typeof(DiffEqBase.ODE_DEFAULT_NORM),typeof(LinearAlgebra.opnorm),DiffEqBase.CallbackSet{Tuple{},Tuple{DiffEqBase.DiscreteCallback{DiffEqCallbacks.var"#30#31",DiffEqCallbacks.SavingAffect{typeof(Main.Test_ode_save_func_callback.simplePendulum2),Float64,Main.Test_ode_save_func_callback.Model,DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Array{Float64,1}},typeof(DiffEqCallbacks.saving_initialize)}}},typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN),typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE),typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK),DataStructures.BinaryHeap{Float64,DataStructures.LessThan},DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Nothing,Nothing,Int64,Array{Float64,1},Float64,Array{Float64,1}},Array{Float64,1},Float64,Nothing}, ::Bool) at HOME\.julia\packages\OrdinaryDiffEq\6TJck\src\solve.jl:474
 [9] #__init#332(::Float64, ::Array{Float64,1}, ::Array{Float64,1}, ::Nothing, ::Bool, ::Bool, ::Bool, ::Bool, ::DiffEqBase.DiscreteCallback{DiffEqCallbacks.var"#30#31",DiffEqCallbacks.SavingAffect{typeof(Main.Test_ode_save_func_callback.simplePendulum2),Float64,Main.Test_ode_save_func_callback.Model,DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Array{Float64,1}},typeof(DiffEqCallbacks.saving_initialize)}, ::Bool, ::Bool, ::Float64, ::Float64, ::Float64, ::Bool, ::Bool, ::Rational{Int64}, ::Nothing, ::Float64, ::Rational{Int64}, ::Int64, ::Int64, ::Int64, ::Rational{Int64}, ::Bool, ::Int64, ::Nothing, ::Nothing, ::Int64, ::typeof(DiffEqBase.ODE_DEFAULT_NORM), ::typeof(LinearAlgebra.opnorm), ::typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), ::typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Int64, ::String, ::typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), ::Nothing, ::Bool, ::Bool, ::Bool, ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(DiffEqBase.__init), ::DiffEqBase.ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Main.Test_ode_save_func_callback.Model{Float64},DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::OrdinaryDiffEq.Tsit5, ::Array{Array{Float64,1},1}, ::Array{Float64,1}, ::Array{Any,1}, ::Type{Val{true}}) at HOME\.julia\packages\OrdinaryDiffEq\6TJck\src\solve.jl:356
 [10] (::DiffEqBase.var"#kw##__init")(::NamedTuple{(:reltol, :saveat, :callback),Tuple{Float64,Float64,DiffEqBase.DiscreteCallback{DiffEqCallbacks.var"#30#31",DiffEqCallbacks.SavingAffect{typeof(Main.Test_ode_save_func_callback.simplePendulum2),Float64,Main.Test_ode_save_func_callback.Model,DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Array{Float64,1}},typeof(DiffEqCallbacks.saving_initialize)}}}, ::typeof(DiffEqBase.__init), ::DiffEqBase.ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Main.Test_ode_save_func_callback.Model{Float64},DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::OrdinaryDiffEq.Tsit5, ::Array{Array{Float64,1},1}, ::Array{Float64,1}, ::Array{Any,1}, ::Type{Val{true}}) at .\none:0 (repeats 5 times)
 [11] #__solve#331 at HOME\.julia\packages\OrdinaryDiffEq\6TJck\src\solve.jl:4 [inlined]
 [12] #__solve at .\none:0 [inlined]
 [13] #solve_call#442(::Base.Iterators.Pairs{Symbol,Any,Tuple{Symbol,Symbol,Symbol},NamedTuple{(:reltol, :saveat, :callback),Tuple{Float64,Float64,DiffEqBase.DiscreteCallback{DiffEqCallbacks.var"#30#31",DiffEqCallbacks.SavingAffect{typeof(Main.Test_ode_save_func_callback.simplePendulum2),Float64,Main.Test_ode_save_func_callback.Model,DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Array{Float64,1}},typeof(DiffEqCallbacks.saving_initialize)}}}}, ::typeof(DiffEqBase.solve_call), ::DiffEqBase.ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Main.Test_ode_save_func_callback.Model{Float64},DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::OrdinaryDiffEq.Tsit5) at HOME\.julia\dev\DiffEqBase\src\solve.jl:40
 [14] #solve_call at .\none:0 [inlined]
 [15] #solve#443 at HOME\.julia\dev\DiffEqBase\src\solve.jl:57 [inlined]
 [16] (::DiffEqBase.var"#kw##solve")(::NamedTuple{(:reltol, :saveat, :callback),Tuple{Float64,Float64,DiffEqBase.DiscreteCallback{DiffEqCallbacks.var"#30#31",DiffEqCallbacks.SavingAffect{typeof(Main.Test_ode_save_func_callback.simplePendulum2),Float64,Main.Test_ode_save_func_callback.Model,DataStructures.BinaryHeap{Float64,DataStructures.LessThan},Array{Float64,1}},typeof(DiffEqCallbacks.saving_initialize)}}}, ::typeof(DiffEqBase.solve), ::DiffEqBase.ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Main.Test_ode_save_func_callback.Model{Float64},DiffEqBase.ODEFunction{true,typeof(Main.Test_ode_save_func_callback.simplePendulum1!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},DiffEqBase.StandardODEProblem}, ::OrdinaryDiffEq.Tsit5) at .\none:0
 [17] top-level scope at HOME\.julia\dev\Test_ode_save_func_callback.jl:63
 [18] include_string(::Module, ::String, ::String) at .\loading.jl:1075
 [19] (::Main._vscodeserver.var"#9#12"{String,Int64,Int64,String})() at HOME\.vscode\extensions\julialang.language-julia-0.13.1\scripts\terminalserver\terminalserver.jl:153
 [20] withpath(::Main._vscodeserver.var"#9#12"{String,Int64,Int64,String}, ::String) at HOME\.vscode\extensions\julialang.language-julia-0.13.1\scripts\terminalserver\repl.jl:62
 [21] (::Main._vscodeserver.var"#8#11"{String,Int64,Int64,String})() at HOME\.vscode\extensions\julialang.language-julia-0.13.1\scripts\terminalserver\terminalserver.jl:152
 [22] hideprompt(::Main._vscodeserver.var"#8#11"{String,Int64,Int64,String}) at HOME\.vscode\extensions\julialang.language-julia-0.13.1\scripts\terminalserver\repl.jl:28
 [23] macro expansion at HOME\.vscode\extensions\julialang.language-julia-0.13.1\scripts\terminalserver\terminalserver.jl:148 [inlined]
 [24] (::Main._vscodeserver.var"#7#10")() at .\task.jl:333
in expression starting at HOME\.julia\dev\Test_ode_save_func_callback.jl:63

This seems to be a bug related to the callback handling

[ChrisRackauckas]

You need a:

using ArrayInterface
ArrayInterface.ismutable(::Type{Model}) = true

This is a tricky one, since it doesn't know whether your Model is mutable since there's no subtyping information on it.

[MartinOtter]

You need a:

using ArrayInterface
ArrayInterface.ismutable(::Type{Model}) = true

This is a tricky one, since it doesn't know whether your Model is mutable since there's no subtyping information on it.

Thanks, this fixed the issue.

Can you give additional information in order that I understand in which situations the above approach needs to be used? Is this information necessary for the parameter "p" data structure but not for all the potential structs that might be used inside "p"? I am also a bit confused because "p" seems to be treated as an array, although it is actually a struct.

I wonder whether there is a simpler way? Is there an abstract mutable type in ArrayInterfaces and it would be fine to define something like:

mutable struct Model{R<:Real} <: XXXX
  ...
end

[ChrisRackauckas]

This is a very hard problem. For your u0, we specialize on whether it's mutable or not, which is separate from whether it's a mutable struct or a struct. For example, an Adjoint{Vector} is something that is a mutable vector even though it's held inside of a struct. Julia's Base Array interface is just missing this concept, but it can be pretty crucial for performance. It's only needed on values you're treating as u0.