Error on indixing into solution by Num

[bgctw]

After update to Julia 1.7.2 and associated update of DifferentialEquations and ModelingToolkit, I get errors when indexing into a solution by a symbolic:

    @variables t 
    D = Differential(t) 
    sts = @variables x(t) RHS(t)  # RHS is observed
    ps = @parameters τ=0.3      # parameters
    @named m = ODESystem([ RHS  ~ (1 - x)/τ, D(x) ~ RHS ], t, sts, ps)
    ms = structural_simplify(m)
prob = ODEProblem(ms, [x => 1.1], (0.0,1.0), [])
sol = solve(prob)

sol(0.3, idxs=[RHS])                                                                                                                                                                     

ERROR: ArgumentError: invalid index: RHS(t) of type Num                                                                                                                                         
Stacktrace:                                                                                                                                                                                     
  [1] to_index(i::Num)                                                                                                                                                                          
    @ Base ./indices.jl:300                                                                                                                                                                     
  [2] to_index(A::Vector{Float64}, i::Num)                                                                                                                                                      
    @ Base ./indices.jl:277                                                                                                                                                                     
  [3] to_indices                                                                                                                                                                                
    @ ./indices.jl:333 [inlined]                                                                                                                                                                
  [4] to_indices                                                                                                                                                                                
    @ ./indices.jl:325 [inlined]                                                                                                                                                                
  [5] getindex                                                                                                                                                                                  
    @ ./abstractarray.jl:1218 [inlined]                                                                                                                                                         
  [6] ode_interpolant                                                                                                                                                                           
    @ ~/scratch/twutz/julia_cluster_depots/packages/OrdinaryDiffEq/RbH1W/src/dense/generic_dense.jl:447 [inlined]                                                                               
  [7] ode_interpolation(tval::Float64, id::OrdinaryDiffEq.CompositeInterpolationData{ODEFunction{true, ModelingToolkit.var"#f#307"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :
ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xe4dfee93, 0x49199599, 0xc5b34c5d, 0x9de7865e, 0xd9dc5b00)}, RuntimeGeneratedFunctions.RuntimeGeneratedFunct
ion{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xd957969e, 0x589ebfe1, 0xc66b6fd7, 0xb7e92edf, 0x9d8be09e)}}, Matrix{Float64}, Nothi$
g, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Vector{Symbol}, Symbol, ModelingToolkit.var"#294#generated_observed#314"{Bool, ODESystem, Dict{Any, Any}}, N$
thing}, Vector{Vector{Float64}}, Vector{Float64}, Vector{Vector{Vector{Float64}}}, OrdinaryDiffEq.CompositeCache{Tuple{OrdinaryDiffEq.Tsit5Cache{Vector{Float64}, Vector{Float64}, Vector{Float6$
}, OrdinaryDiffEq.Tsit5ConstantCache{Float64, Float64}, typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, OrdinaryDiffEq.Rosenbrock23Cache{Vector$
Float64}, Vector{Float64}, Vector{Float64}, Matrix{Float64}, Matrix{Float64}, OrdinaryDiffEq.Rosenbrock23Tableau{Float64}, SciMLBase.TimeGradientWrapper{ODEFunction{true, ModelingToolkit.var"#$
#307"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xe4dfee93, 0x49199599, 0xc5b34c5d, 0x9d$
7865e, 0xd9dc5b00)}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xd957969e, 0x589$
bfe1, 0xc66b6fd7, 0xb7e92edf, 0x9d8be09e)}}, Matrix{Float64}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Vector{Symbol}, Symbol, ModelingToolkit.$
ar"#294#generated_observed#314"{Bool, ODESystem, Dict{Any, Any}}, Nothing}, Vector{Float64}, Vector{Float64}}, SciMLBase.UJacobianWrapper{ODEFunction{true, ModelingToolkit.var"#f#307"{RuntimeG$
neratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xe4dfee93, 0x49199599, 0xc5b34c5d, 0x9de7865e, 0xd9dc5$
00)}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xd957969e, 0x589ebfe1, 0xc66b6f$
7, 0xb7e92edf, 0x9d8be09e)}}, Matrix{Float64}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Vector{Symbol}, Symbol, ModelingToolkit.var"#294#genera$
ed_observed#314"{Bool, ODESystem, Dict{Any, Any}}, Nothing}, Float64, Vector{Float64}}, LinearSolve.LinearCache{Matrix{Float64}, Vector{Float64}, Vector{Float64}, SciMLBase.NullParameters, Gen$
ricFactorization{LinearSolve.RFWrapper{true, true}}, LinearAlgebra.LU{Float64, Matrix{Float64}}, LinearSolve.InvPreconditioner{LinearAlgebra.Diagonal{Float64, Vector{Float64}}}, LinearAlgebra.$
iagonal{Float64, Vector{Float64}}, Float64}, FiniteDiff.JacobianCache{Vector{Float64}, Vector{Float64}, Vector{Float64}, UnitRange{Int64}, Nothing, Val{:forward}(), Float64}, FiniteDiff.Gradie$
tCache{Nothing, Vector{Float64}, Vector{Float64}, Float64, Val{:forward}(), Float64, Val{true}()}, Float64, Rosenbrock23{1, false, GenericFactorization{LinearSolve.RFWrapper{true, true}}, type$
f(OrdinaryDiffEq.DEFAULT_PRECS), Val{:forward}, true, nothing}}}, OrdinaryDiffEq.AutoSwitchCache{Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.$
alse}, Rosenbrock23{0, false, Nothing, typeof(OrdinaryDiffEq.DEFAULT_PRECS), Val{:forward}, true, nothing}, Rational{Int64}, Int64}}}, idxs::Vector{Num}, deriv::Type{Val{0}}, p::Vector{Float64}
, continuity::Symbol)
    @ OrdinaryDiffEq ~/scratch/twutz/julia_cluster_depots/packages/OrdinaryDiffEq/RbH1W/src/dense/generic_dense.jl:353
  [8] CompositeInterpolationData
    @ ~/scratch/twutz/julia_cluster_depots/packages/OrdinaryDiffEq/RbH1W/src/interp_func.jl:82 [inlined]
  [9] #_#496
    @ ~/scratch/twutz/julia_cluster_depots/packages/OrdinaryDiffEq/RbH1W/src/composite_solution.jl:16 [inlined]
 [10] top-level scope
    @ REPL[64]:1

With Julia 1.7.2 and packages

  [2b5f629d] DiffEqBase v6.81.3
  [0c46a032] DifferentialEquations v7.1.0
  [2ee39098] LabelledArrays v1.7.0
  [961ee093] ModelingToolkit v8.3.2
  [ae029012] Requires v1.3.0
  [90137ffa] StaticArrays v1.3.3

[YingboMa]

I'll make a new release. This is fixed on master already.

[bgctw]

For me the error persists even after updating all packages and dev ModelingToolkit.

      Status `/Net/Groups/BGI/people/twutz/julia/dev/MTKHelpers/Project.toml`
  [2b5f629d] DiffEqBase v6.81.3
  [0c46a032] DifferentialEquations v7.1.0
  [2ee39098] LabelledArrays v1.7.1
  [961ee093] ModelingToolkit v8.4.0 `~/julia/dev/ModelingToolkit`
  [ae029012] Requires v1.3.0
  [90137ffa] StaticArrays v1.3.4

[YingboMa]

Ah, that's a different problem than I imagined. CompositeInterpolationData doesn't handle symbolic indices. @sharanry do you mind to handle that?

[bgctw]

Is there a workaround that I can use in the meantime?

I can interpolate states from the solution. How can I compute an observable given a state vector?

[YingboMa]

Just specify a solver.

[bgctw]

Could you, please, give an example?
I get the interpolated state, and how to go on then?

s = sol(0.3) 
rhs = ?

[YingboMa]

Specify an algorithm that's not a composite algorithm. Like sol = solve(prob, Tsit5()) or sol = solve(prob, Rodas4())

[sharanry]

Ah, that's a different problem than I imagined. CompositeInterpolationData doesn't handle symbolic indices. @sharanry do you mind to handle that?

I will take a look.

[bgctw]

Thanks. That works. Now I can follow and understand that the problem depends on the chosen solver.

[ChrisRackauckas]

SciML/SciMLBase.jl#141 + SciML/OrdinaryDiffEq.jl#1612 should fix this.

[ChrisRackauckas]

using ModelingToolkit, OrdinaryDiffEq

@variables t
D = Differential(t)
sts = @variables x(t) RHS(t)  # RHS is observed
ps = @parameters τ = 0.3      # parameters
@named m = ODESystem([RHS ~ (1 - x) / τ, D(x) ~ RHS], t, sts, ps)
ms = structural_simplify(m)
prob = ODEProblem(ms, [x => 1.1], (0.0, 1.0), [])
sol = solve(prob,AutoVern6(Rodas5()))
sol(0.3, idxs = [RHS]) 

Composite now works on latest releases.