We need more helpful error messages in ModelingToolkit

dodoplus · June 17, 2022, 1:09pm

Hi,

Trying to solve a model, I’m getting a really unhelpful message:

@named x = MixingVolume()
sys = structural_simplify(x)

prob = ODAEProblem(sys, [], (0, 10.0))
sol = solve(prob, Tsit5())

MethodError: no method matching oneunit(::Type{Any})
Closest candidates are:
  oneunit(::Type{Union{Missing, T}}) where T at missing.jl:105
  oneunit(::Type{T}) where T at number.jl:358
  oneunit(::T) where T at number.jl:357
  ...

Stacktrace:
  [1] oneunit(#unused#::Type{Any})
    @ Base ./missing.jl:106
  [2] __init(prob::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, alg::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, timeseries_init::Tuple{}, ts_init::Tuple{}, ks_init::Tuple{}, recompile::Type{Val{true}}; saveat::Tuple{}, 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.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/ZgJ9s/src/solve.jl:159
  [3] __init(prob::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, alg::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, timeseries_init::Tuple{}, ts_init::Tuple{}, ks_init::Tuple{}, recompile::Type{Val{true}}) (repeats 5 times)
    @ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/ZgJ9s/src/solve.jl:9
  [4] __solve(::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, ::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}; kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/ZgJ9s/src/solve.jl:4
  [5] __solve(::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, ::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False})
    @ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/ZgJ9s/src/solve.jl:1
  [6] solve_call(_prob::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, args::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}; merge_callbacks::Bool, kwargshandle::KeywordArgError, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ DiffEqBase ~/.julia/packages/DiffEqBase/hZncn/src/solve.jl:428
  [7] solve_call(_prob::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, args::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False})
    @ DiffEqBase ~/.julia/packages/DiffEqBase/hZncn/src/solve.jl:409
  [8] solve_up(prob::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, sensealg::Nothing, u0::Nothing, p::Vector{Float64}, args::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}; kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ DiffEqBase ~/.julia/packages/DiffEqBase/hZncn/src/solve.jl:726
  [9] solve_up(prob::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, sensealg::Nothing, u0::Nothing, p::Vector{Float64}, args::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False})
    @ DiffEqBase ~/.julia/packages/DiffEqBase/hZncn/src/solve.jl:713
 [10] solve(prob::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, args::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}; sensealg::Nothing, u0::Nothing, p::Nothing, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ DiffEqBase ~/.julia/packages/DiffEqBase/hZncn/src/solve.jl:710
 [11] solve(prob::ODEProblem{Nothing, Tuple{Float64, Float64}, true, Vector{Float64}, ODEFunction{true, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(Symbol("##out#386"), Symbol("##arg#8441539597060760774"), Symbol("##arg#10374657382969776836"), :t), ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", ModelingToolkit.StructuralTransformations.var"#_RGF_ModTag", (0xcaccd448, 0x42cc35f6, 0xedb19974, 0xe85f4fda, 0x771b6431)}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, SparseArrays.SparseMatrixCSC{Bool, Int64}, Nothing, Nothing, Nothing, Vector{Symbol}, Nothing, ModelingToolkit.StructuralTransformations.var"#generated_observed#226"{Bool, ModelingToolkit.BipartiteGraphs.Matching{ModelingToolkit.BipartiteGraphs.Unassigned, Vector{Union{ModelingToolkit.BipartiteGraphs.Unassigned, Int64}}}, TearingState{ODESystem}, Dict{Any, Any}, BitVector, Vector{SymbolicUtils.Code.Assignment}, Tuple{Vector{Vector{Int64}}, Vector{BitSet}}, SymbolicUtils.Code.NameState, Dict{Any, Int64}}, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, args::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False})
    @ DiffEqBase ~/.julia/packages/DiffEqBase/hZncn/src/solve.jl:703

Yes, I get it that the model is not fully specified (at least, I think that’s what this is about).
However, this message doesn’t help in any way to understand the problem or help me resolve it…

DD

Sukera · June 17, 2022, 3:03pm

That array is untyped, so it has eltype of Any. I think that’s what it ends up complaining about.

dodoplus · June 17, 2022, 3:34pm

Right.

My point wasn’t the specific issue but the unhelpful error message.

ODAEProblem() should either force the input to comply with its requirements (e.g. ::Vector{Equation}, or check its parameters (like a Python library would).

ChrisRackauckas · June 17, 2022, 5:30pm

I agree. File an issue for every case that hits a bad error message and we will get those cleaned up.

dodoplus · June 18, 2022, 5:11am

https://github.com/SciML/ModelingToolkit.jl/issues/1654

dodoplus · June 18, 2022, 5:33am

This code produces the same error:

using ModelingToolkit, OrdinaryDiffEq

@variables t
@named model = ODESystem(Equation[], t, [], [])
sys = structural_simplify(model)
prob = ODEProblem(sys, [], (0, 1.0))
sol = solve(prob, Tsit5())

This doesn’t help either:

...
@variables x(t)
@parameters u
@named model = ODESystem(Equation[], t, [x], [u])
...

odow · June 18, 2022, 5:35am

We’re trying to develop the JuMP style guide to address MethodErrors like this one.

Here’s the in progress PR if anyone reading this has thoughts or suggestions:

https://github.com/jump-dev/JuMP.jl/pull/2995

ChrisRackauckas · June 18, 2022, 6:32am

All of that is already in the SciML Style Guide.

A guide isn’t the help here. Look at the actual issue. It has nothing to do with ModelingToolkit. The MWE is just:

using OrdinaryDiffEq
prob = ODEProblem((u,p,t)->u, [], (0, 1.0))
sol = solve(prob, Tsit5())

The problem is that we haven’t handled the “null ODE” case, which shouldn’t error but actually directly build a solution. No matter what style guide you tell someone to look at, people can miss that yes it’s possible to pass an ODE with no ODEs.

This only recently came to the forefront though because of ModelingToolkit’s ability to simplify equations, and so someone can end up with a null ODE from a real code simply because it ends up replacing it with observed variables for the analytical solution. The correct thing to do is not to error either, because then the user doesn’t get the analytical solution, so instead it needs to do a direct build of the ODESolution with the observed variables so the user can seamlessly interact with the solution as though it had actually had to solve something. Otherwise the user won’t get the answer, so any error is bad here.

Our style guide that promotes defense programming and has many principles that you explicitly reject that are for defensive programming. We have an entire system that catches high level issues before things ever hit solvers:

github.com

SciML/DiffEqBase.jl/blob/master/src/solve.jl#L113-L370


      
          function Base.showerror(io::IO, e::CommonKwargError)
              println(io, KWARGERROR_MESSAGE)
              notin = collect(map(x -> x ∉ allowedkeywords, keys(e.kwargs)))
              unrecognized = collect(keys(e.kwargs))[notin]
              print(io, "Unrecognized keyword arguments: ")
              printstyled(io, unrecognized; bold = true, color = :red)
              print(io, "\n\n")
              println(io, TruncatedStacktraces.VERBOSE_MSG)
          end
          
          @enum KeywordArgError KeywordArgWarn KeywordArgSilent
          
          const INCOMPATIBLE_U0_MESSAGE = """
                                          Initial condition incompatible with functional form.
                                          Detected an in-place function with an initial condition of type Number or SArray.
                                          This is incompatible because Numbers cannot be mutated, i.e.
                                          `x = 2.0; y = 2.0; x .= y` will error.
          
                                          If using a immutable initial condition type, please use the out-of-place form.
                                          I.e. define the function `du=f(u,p,t)` instead of attempting to "mutate" the immutable `du`.

This file has been truncated. show original

is just a taste of a piece of it. Sometimes, something can slip by, and no amount of telling someone to do defensive programming will capture every case. Sometimes, there’s a case that no one thought about, sorry.

FWIW, JuMP misses this case too, and silently just returns the null model is not optimized instead of erroring or saying it’s optimized.

using JuMP, Ipopt
model = Model(Ipopt.Optimizer)
optimize!(model)
termination_status(model) == LOCALLY_SOLVED # false, should be true or error

So anyways, @dodoplus I’m on this one right now but it’s a bit less trivial to fix because ODESolution needs to have its show printouts and everything robust to the case where sol.u = nothing or an empty array, a case that quite frankly we didn’t plan for before but now that we have symbolic simplification can show up in real-world codes.

P.S. and all of this talk about error messages, @dodoplus is using an old version. If your packages were up to date your error message would’ve been:

julia> sol = solve(prob, Tsit5())
ERROR: Non-concrete element type inside of an `Array` detected.
Arrays with non-concrete element types, such as
`Array{Union{Float32,Float64}}`, are not supported by the
differential equation solvers. Anyways, this is bad for
performance so you don't want to be doing this!

If this was a mistake, promote the element types to be
all the same. If this was intentional, for example,
using Unitful.jl with different unit values, then use
an array type which has fast broadcast support for
heterogeneous values such as the ArrayPartition
from RecursiveArrayTools.jl. For example:

using RecursiveArrayTools
x = ArrayPartition([1.0,2.0],[1f0,2f0])
y = ArrayPartition([3.0,4.0],[3f0,4f0])
x .+ y # fast, stable, and usable as u0 into DiffEq!

Element type:
Any

So clearly we have done a lot of defense programming errors (and have improved it a lot over the last 2 months), but we missed a case so even the defense warning was a bit off.

dodoplus · June 18, 2022, 7:10am

Thanks.

ccoffrin · June 19, 2022, 4:58am

ChrisRackauckas:

FWIW, JuMP misses this case too, and silently just returns the null model is not optimized instead of erroring or saying it’s optimized.
using JuMP, Ipopt
model = Model(Ipopt.Optimizer)
optimize!(model)
termination_status(model) == LOCALLY_SOLVED # false, should be true or error

In this example JuMP returns INVALID_MODEL::TerminationStatusCode = 21 I think that is a suitable status for this case, what am I missing?

ChrisRackauckas · June 19, 2022, 9:55am

If in the near future you build add a connection to a symbolic system that is able to eliminate variables and those are captured as observation variables, and then the user solves this model and checks a termination code to see if it solved correctly, this will silently return that it’s invalid and not a correct solution with just some hidden observed variables a user code check for validity would make them think it failed. That seems silly, but that’s exactly the problem we have here, where in this case we throw an error because we didn’t think about this case ever existing.

The point is that even the most defensive programming will miss cases that are unexpected because you can’t really be on a defense for unknown scope expansions, except for trying to throw errors, and ModelingToolkit was a massive recent scope expansion so it’s still shaking up some of the interfaces which made assumptions that are now not correct (like that solving a model with 0 equations is invalid or a non-interesting use case).

ChrisRackauckas · June 29, 2022, 4:13pm

Fix a typo in tearing by YingboMa · Pull Request #1676 · SciML/ModelingToolkit.jl · GitHub
Handle null differential equation and steady state cases by ChrisRackauckas · Pull Request #787 · SciML/DiffEqBase.jl · GitHub

That’s a full solution to this issue. So then this works:

using ModelingToolkit, OrdinaryDiffEq, SteadyStateDiffEq, Test

@variables t x(t) y(t)
eqs = [
    0 ~ x - y
    0 ~ y - x
 ]

 @named sys = ODESystem(eqs, t)
 sys = structural_simplify(sys)

 prob = ODEProblem(sys, Pair[], (0, 10.0))
 sol = solve(prob, Tsit5())

 @test sol[x] == [0.0, 0.0]
 @test sol[y] == [0.0, 0.0]

and it will just use its analytical solution, so plot(sol,vars = x) and everything will magically work. Just don’t tell anyone no ODEs were solved and we’re good.