JuliaDiff / JuliaDiff/ForwardDiff.jl
`InexactError` when trying to differentiate `ODEProblem` using `LabelledArrays`
Nobody has claimed this yet.
- Dominant language
- Julia
- Stars
- 1k
- Forks
- 160
- PR merge metrics
- No merged PRs in 30d
Description
It looks like the source of the problem is that recursive_unitless_bottom_eltype from RecursiveArrayTools is inferring an Int64 type rather than a Float64 type from the Dual of a LabelledArray. This is happening somewhere inside ODEProblem's __init.
Note that this works fine if you switch from LVector to Vector for the parameter/initial-condition arrays.
using LabelledArrays
using ModelingToolkit
using DifferentialEquations
using ForwardDiff
pars = @parameters β α
@variables t E(t)
D = Differential(t)
eqs = [D(E) ~ α*t + β]
sys = ODESystem(eqs)
u0 = LVector(E = 1.0)
p = LVector(β = 0.8,
α = 0.75)
tspan = (0.0,1.0)
sys = structural_simplify(sys)
prob = ODEProblem(sys,u0,tspan,p)
function loss(inputs)
N = length(prob.p)
_p = inputs[1:N]
_u0 = inputs[N:end]
tmp_sol = solve(prob,Tsit5(),p=_p,u0=_u0,saveat=0.1)
sum(abs2,tmp_sol[E])
end
#Substitute a subset of the parameters/initial-conditions into the full set
inputter(par_vec) = LVector(β = p[1], α = par_vec[1], E = par_vec[2])
p0 = LVector(α = 1.0, E = 0.5)
dummy_loss(par_vec) = sum(abs2,par_vec)
ForwardDiff.gradient(x->dummy_loss(inputter(x)),p0) #works
ForwardDiff.gradient(dummy_loss,inputter(p0)) #Works
ForwardDiff.gradient(loss,inputter(p0)) #Works
ForwardDiff.gradient(x->loss(inputter(x)),p0) # Error
Error message:
ERROR: InexactError: Int64(161//1000)
Stacktrace:
[1] Integer
@ ./rational.jl:109 [inlined]
[2] convert(#unused#::Type{Int64}, x::Rational{Int64})
@ Base /Applications/Julia-1.6.app/Contents/Resources/julia/share/julia/base/number.jl:7
[3] OrdinaryDiffEq.Tsit5ConstantCache(T::Type, T2::Type)
@ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/a5qvy/src/tableaus/low_order_rk_tableaus.jl:610
[4] alg_cache(alg::Tsit5, u::Vector{Real}, rate_prototype::Vector{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}}, uEltypeNoUnits::Type, uBottomEltypeNoUnits::Type, tTypeNoUnits::Type, uprev::Vector{Real}, uprev2::Vector{Real}, f::ODEFunction{true, ModelingToolkit.var"#f#165"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xbabc66d3, 0x633a919b, 0xe4683138, 0x0e254e2d, 0x36925fe8)}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x7eb6766a, 0x02a244e1, 0x54678655, 0x9f7f1370, 0x31962eb1)}}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Vector{Symbol}, Symbol, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, t::Float64, dt::Float64, reltol::Vector{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}}, p::Vector{Real}, calck::Bool, #unused#::Val{true})
@ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/a5qvy/src/caches/low_order_rk_caches.jl:348
[5] __init(prob::ODEProblem{Vector{Real}, Tuple{Float64, Float64}, true, Vector{Real}, ODEFunction{true, ModelingToolkit.var"#f#165"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xbabc66d3, 0x633a919b, 0xe4683138, 0x0e254e2d, 0x36925fe8)}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x7eb6766a, 0x02a244e1, 0x54678655, 0x9f7f1370, 0x31962eb1)}}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Vector{Symbol}, Symbol, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, alg::Tsit5, timeseries_init::Tuple{}, ts_init::Tuple{}, ks_init::Tuple{}, recompile::Type{Val{true}}; saveat::Float64, tstops::Tuple{}, d_discontinuities::Tuple{}, save_idxs::Nothing, save_everystep::Bool, save_on::Bool, save_start::Bool, save_end::Nothing, callback::Nothing, dense::Bool, calck::Bool, dt::Float64, dtmin::Nothing, dtmax::Float64, force_dtmin::Bool, adaptive::Bool, gamma::Rational{Int64}, abstol::Nothing, reltol::Nothing, qmin::Rational{Int64}, qmax::Int64, qsteady_min::Int64, qsteady_max::Int64, beta1::Nothing, beta2::Nothing, qoldinit::Rational{Int64}, controller::Nothing, fullnormalize::Bool, failfactor::Int64, maxiters::Int64, internalnorm::typeof(DiffEqBase.ODE_DEFAULT_NORM), internalopnorm::typeof(LinearAlgebra.opnorm), isoutofdomain::typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), unstable_check::typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), verbose::Bool, timeseries_errors::Bool, dense_errors::Bool, advance_to_tstop::Bool, stop_at_next_tstop::Bool, initialize_save::Bool, progress::Bool, progress_steps::Int64, progress_name::String, progress_message::typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), userdata::Nothing, allow_extrapolation::Bool, initialize_integrator::Bool, alias_u0::Bool, alias_du0::Bool, initializealg::OrdinaryDiffEq.DefaultInit, kwargs::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
@ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/a5qvy/src/solve.jl:295
[6] __solve(::ODEProblem{Vector{Real}, Tuple{Float64, Float64}, true, Vector{Real}, ODEFunction{true, ModelingToolkit.var"#f#165"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xbabc66d3, 0x633a919b, 0xe4683138, 0x0e254e2d, 0x36925fe8)}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x7eb6766a, 0x02a244e1, 0x54678655, 0x9f7f1370, 0x31962eb1)}}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Vector{Symbol}, Symbol, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, ::Tsit5; kwargs::Base.Iterators.Pairs{Symbol, Float64, Tuple{Symbol}, NamedTuple{(:saveat,), Tuple{Float64}}})
@ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/a5qvy/src/solve.jl:4
[7] solve_call(_prob::ODEProblem{Vector{Real}, Tuple{Float64, Float64}, true, Vector{Real}, ODEFunction{true, ModelingToolkit.var"#f#165"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xbabc66d3, 0x633a919b, 0xe4683138, 0x0e254e2d, 0x36925fe8)}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x7eb6766a, 0x02a244e1, 0x54678655, 0x9f7f1370, 0x31962eb1)}}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Vector{Symbol}, Symbol, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, args::Tsit5; merge_callbacks::Bool, kwargs::Base.Iterators.Pairs{Symbol, Float64, Tuple{Symbol}, NamedTuple{(:saveat,), Tuple{Float64}}})
@ DiffEqBase ~/.julia/packages/DiffEqBase/oe7VF/src/solve.jl:61
[8] solve_up(prob::ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, ModelingToolkit.var"#f#165"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xbabc66d3, 0x633a919b, 0xe4683138, 0x0e254e2d, 0x36925fe8)}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x7eb6766a, 0x02a244e1, 0x54678655, 0x9f7f1370, 0x31962eb1)}}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Vector{Symbol}, Symbol, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, sensealg::Nothing, u0::Vector{Real}, p::Vector{Real}, args::Tsit5; kwargs::Base.Iterators.Pairs{Symbol, Float64, Tuple{Symbol}, NamedTuple{(:saveat,), Tuple{Float64}}})
@ DiffEqBase ~/.julia/packages/DiffEqBase/oe7VF/src/solve.jl:82
[9] #solve#59
@ ~/.julia/packages/DiffEqBase/oe7VF/src/solve.jl:70 [inlined]
[10] loss(inputs::LArray{Real, 1, Vector{Real}, (:β, :α, :E)})
@ Main ./REPL[227]:5
[11] (::var"#65#66")(x::LArray{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}, 1, Vector{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}}, (:α, :E)})
@ Main ./REPL[234]:1
[12] vector_mode_dual_eval
@ ~/.julia/packages/ForwardDiff/QOqCN/src/apiutils.jl:37 [inlined]
[13] vector_mode_gradient(f::var"#65#66", x::LArray{Float64, 1, Vector{Float64}, (:α, :E)}, cfg::ForwardDiff.GradientConfig{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2, LArray{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}, 1, Vector{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}}, (:α, :E)}})
@ ForwardDiff ~/.julia/packages/ForwardDiff/QOqCN/src/gradient.jl:106
[14] gradient(f::Function, x::LArray{Float64, 1, Vector{Float64}, (:α, :E)}, cfg::ForwardDiff.GradientConfig{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2, LArray{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}, 1, Vector{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}}, (:α, :E)}}, ::Val{true})
@ ForwardDiff ~/.julia/packages/ForwardDiff/QOqCN/src/gradient.jl:19
[15] gradient(f::Function, x::LArray{Float64, 1, Vector{Float64}, (:α, :E)}, cfg::ForwardDiff.GradientConfig{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2, LArray{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}, 1, Vector{ForwardDiff.Dual{ForwardDiff.Tag{var"#65#66", Float64}, Float64, 2}}, (:α, :E)}}) (repeats 2 times)
@ ForwardDiff ~/.julia/packages/ForwardDiff/QOqCN/src/gradient.jl:17
[16] top-level scope
@ REPL[234]:1
Contributor guide
No contributing guide indexed for this repository
First steps
- Read the whole issue, then the project's contributing guide.
- Comment on the issue to say you are picking it up — it saves two people doing the same work.
- Fork the repository and make your change on a branch.
- Open a pull request that references the issue number.
Research direction
The report points to recursive_unitless_bottom_eltype in RecursiveArrayTools and ODEProblem's init; begin by reproducing the supplied LVector/Vector comparison and follow the stack trace into OrdinaryDiffEq's Tsit5 cache. Done means the final ForwardDiff.gradient call completes without InexactError while the Vector case remains working.
Written by the indexing model from the issue text.
Assessment
- Tech stack
- julia
- Domain
- devtools
- Issue type
- Bug
- Difficulty
- 4/5
- Estimated time
- 3-5 days
- Activity status
- Stale
- Clarity
- Mostly clear
- Newbie friendliness
- 35/100