Nested Labelled Arrays with Modelingtoolkitize

[jbiffl]

I'm having the same issues as in this previous issue but with nested labelled vectors with the latest version. The code below gives the same kind of error as the original example by fkrauer, but only when using nested LVectors now that single layer LVectors work.

using DifferentialEquations
using LabelledArrays
using ModelingToolkit

# ODE model: nonsense based on lorenz
function ODEFun!(du,u,p,t)

    # params in two layers does not work

    a = p.firstset.a
    b = p.firstset.b
    c = p.firstset.c

    d = p.secondset.d
    e = p.secondset.e
    f = p.secondset.f

    #params in one layer does work

    #=
    a = p.a
    b = p.b
    c = p.c
    d = p.d
    e = p.e
    f = p.f
    =#

    # change in states
    du[1] = a*(u[2] - u[1])+d
    du[2] = u[1]*(b-u[3]-u[2])+e
    du[3] = u[1]*u[2] - c*u[3]+f

end

# Solver settings
tmin = 0.0
tmax = 10.0
tspan = (tmin, tmax)
solvsettings = (abstol = 1.0e-6, 
                reltol = 1.0e-6, 
                saveat = 1.0,
                solver = Rodas4())

# ODE Parameters as nested LVector does not work
p = LVector(firstset = LVector(a = 10, b = 29, c = 8/3), secondset = LVector(d = 0, e = 0, f = 0)) #value of variables has no impact on error
# ODE Parameters as single layer LVector works
#p = LVector(a = 10, b = 28, c = 8/3, d = 0, e = 0, f = 0)
u0 = [1.0 0.0 0.0]

# Initiate ODE problem
problem = ODEProblem(ODEFun!,u0,tspan,p)
problem_acc = ODEProblem(modelingtoolkitize(problem), u0, tspan, p, jac=true, sparse=true)

sol = solve(problem_acc, 
            solvsettings.solver, 
            abstol=solvsettings.abstol, 
            reltol=solvsettings.reltol, 
            isoutofdomain=(u,p,t)->any(x->x<0.0,u), 
            saveat=solvsettings.saveat)

That code throws this error:

ERROR: LoadError: type Num has no field a
Stacktrace:
 [1] getproperty(x::Num, f::Symbol)
   @ Base .\Base.jl:33
 [2] ODEFun!(du::Matrix{Num}, u::Matrix{Num}, p::LArray{Num, 1, Vector{Num}, (:firstset, :secondset)}, t::Num)

When you inspect p while it is running, you obtain

2-element LArray{Num, 1, Vector{Num}, (:firstset, :secondset)}:
:firstset => firstset
:secondset => secondset

[ChrisRackauckas]

Thanks. This is definitely not a case I considered, so I'm pushing a reminder to the weekend.

[jonniedie]

BTW, ComponentArrays.jl is basically a version of of LabelledArrays that can be nested like this. However, it doesn't actually solve the modelingtoolkitize problem because that's a separate thing. I've had it on my TODO list to make it easier to modelingtoolkitize ComponentArrays problems and "componentarrayize" ModelingToolkit problems. The reason the latter is useful is it would let you more easily write callbacks and such for ModelingToolkit models. To be honest, it shouldn't be too hard to convert between the two, I just haven't gotten around to it yet. And I'm not sure where I want it to live (because I have a few other things I'd like to do to make it easier to convert between the two and gets to the point where it seems silly to do that much with Requires).

[ChrisRackauckas]

https://github.com/SciML/ModelingToolkit.jl/blob/master/src/systems/diffeqs/modelingtoolkitize.jl#L82-L93

Some of LabelledArrays is already specialized. I think we could just add a hook to ComponentArrays here? I think more and more it would be good to just integrate it fully into the SciML ecosystem. define_vars would just need to be made recursive for ComponentArrays to work well there though.