SciML / SciMLBase.jl

The Base interface of the SciML ecosystem
https://docs.sciml.ai/SciMLBase/stable
MIT License
131 stars 97 forks source link

Symbolic Indexing mostly broken for EnsembleProblems (remake works though) #661

Closed TorkelE closed 5 months ago

TorkelE commented 5 months ago

As title, here's a list of cases to check.Technically there are Ensemble solutions as well, however, I am not sure of symbolic indexing should be a thing here or not (one could make it work, but would have to return the values across all solutions).

using OrdinaryDiffEq, ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as D
using ModelingToolkit: getu, getp, setu, setp

@parameters p d
@variables X(t) X2(t)
eqs = [
    D(X) ~ p - d*X,
    X2 ~ 2*X
]
@mtkbuild osys = ODESystem(eqs, t)

u0 = [X => 0.1]
tspan = (0.0, 10.0)
p = [p => 1.0, d => 0.2]
oprob = ODEProblem(osys, u0, tspan, p)
eprob = EnsembleProblem(oprob)

# Get stuff.
eprob[X] # ERROR: type EnsembleProblem has no field f
eprob[:X] # ERROR: type EnsembleProblem has no field f
eprob[osys.X] # ERROR: type EnsembleProblem has no field f
eprob.ps[p] # ERROR: type EnsembleProblem has no field f
eprob.ps[:p] # ERROR: type EnsembleProblem has no field f
eprob.ps[osys.p] # ERROR: type EnsembleProblem has no field f
getu(eprob, X)(eprob) # ERROR: type EnsembleProblem has no field f
getu(eprob, :X)(eprob) # ERROR: type EnsembleProblem has no field f
getu(eprob, osys.X)(eprob) # ERROR: type EnsembleProblem has no field f
getp(eprob, p)(eprob) # ERROR: type EnsembleProblem has no field f
getp(eprob, :p)(eprob) # ERROR: type EnsembleProblem has no field f
getp(eprob, osys.p)(eprob) # ERROR: type EnsembleProblem has no field f

# Set stuff.
eprob[X] = 0.0 # ERROR: type EnsembleProblem has no field f
eprob[:X] = 0.0 # ERROR: type EnsembleProblem has no field f
eprob[osys.X] = 0.0 # ERROR: type EnsembleProblem has no field f
eprob.ps[p] = 0.0 # ERROR: type EnsembleProblem has no field f
eprob.ps[:p] = 0.0 # ERROR: type EnsembleProblem has no field f
eprob.ps[osys.p] = 0.0 # ERROR: type EnsembleProblem has no field f
setu(eprob, X)(eprob, 0.0) # ERROR: type EnsembleProblem has no field f
setu(eprob, :X)(eprob, 0.0) # ERROR: type EnsembleProblem has no field f
setu(eprob, osys.X)(eprob, 0.0) # ERROR: type EnsembleProblem has no field f
setp(eprob, p)(eprob, 0.0) # ERROR: type EnsembleProblem has no field f
setp(eprob, :p)(eprob, 0.0) # ERROR: type EnsembleProblem has no field f
setp(eprob, osys.p)(eprob, 0.0) # ERROR: type EnsembleProblem has no field f

# Remake works.
remake(eprob; u0 = [X => 0.1]).prob.u0 # Gives [0.1]
remake(eprob; u0 = [:X => 0.2]).prob.u0 # Gives [0.2]
remake(eprob; u0 = [osys.X => 0.3]).prob.u0 # Gives [0.3]

remake(eprob; p = [d => 0.4]) # Probably works (no errors, but cannot access parameter directly).
remake(eprob; p = [:d => 0.5]) # Probably works (no errors, but cannot access parameter directly).
remake(eprob; p = [osys.d => 0.6]) # Probably works (no errors, but cannot access parameter directly).  
ChrisRackauckas commented 5 months ago

Yeah, I hadn't considered this case. That's definitely an oversight we need to fix.

TorkelE commented 5 months ago

Quite a few of these issues I'm raising are not that critical (e.g. https://github.com/SciML/SciMLBase.jl/issues/660). However, I think I have gotten quite good at finding these now, so I am just running large scans of cases that I think we might have missed. Then we can just tick of the list of corner cases as there is time to sort them out.