search

SciML / ModelingToolkit.jl

An acausal modeling framework for automatically parallelized scientific machine learning (SciML) in Julia. A computer algebra system for integrated symbolics for physics-informed machine learning and automated transformations of differential equations
https://mtk.sciml.ai/dev/
Other
1.41k stars 204 forks source link

Unable to handle unit literals #2340

Open mtiller-jh opened 11 months ago

mtiller-jh commented 11 months ago

This is related to #2339. If I run the following MWE:

using ModelingToolkit
using DifferentialEquations
using Plots
using Unitful

@variables t x(t)
D = Differential(t)

@mtkmodel FirstOrder begin
    @structural_parameters begin
        f
    end
    @variables begin
        x(t), [unit = u"V"]
    end
    @equations begin
        D(x) ~ f(t) * 1 u"V"
    end
end

vf = function(t) sin(t) end
@mtkbuild root = FirstOrder(f=vf)
u0 = [
    root.x => 0
]
prob = ODEProblem(root, u0, (0, 2))
sol = solve(prob)
display(plot(sol, idxs=[root.x]))

Originally, I had:

@mtkmodel FirstOrder begin
    @structural_parameters begin
        f
    end
    @variables begin
        x(t), [unit = u"V"]
    end
    @equations begin
        D(x) ~ f(t) * 1
    end
end

But this leads to:

┌ Warning:  in eq. #1: units [V] for left and [] for right do not match.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/validation.jl:176
ERROR: ModelingToolkit.ValidationError("Some equations had invalid units. See warnings for details.")
Stacktrace:
 [1] check_units(eqs::Vector{Equation})
   @ ModelingToolkit ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/validation.jl:277
 [2] ODESystem(tag::UInt64, deqs::Vector{Equation}, iv::SymbolicUtils.BasicSymbolic{Real}, dvs::Vector{SymbolicUtils.BasicSymbolic{Real}}, ps::Vector{Any}, tspan::Nothing, var_to_name::Dict{Any, Any}, ctrls::Vector{Any}, observed::Vector{Equation}, tgrad::Base.RefValue{Vector{Num}}, jac::Base.RefValue{Any}, ctrl_jac::Base.RefValue{Any}, Wfact::Base.RefValue{Matrix{Num}}, Wfact_t::Base.RefValue{Matrix{Num}}, name::Symbol, systems::Vector{Any}, defaults::Dict{Any, Any}, torn_matching::Nothing, connector_type::Nothing, preface::Nothing, cevents::Vector{ModelingToolkit.SymbolicContinuousCallback}, devents::Vector{ModelingToolkit.SymbolicDiscreteCallback}, metadata::Nothing, gui_metadata::ModelingToolkit.GUIMetadata, tearing_state::Nothing, substitutions::Nothing, complete::Bool, discrete_subsystems::Nothing, unknown_states::Nothing, split_idxs::Nothing, parent::Nothing; checks::Bool)
   @ ModelingToolkit ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/diffeqs/odesystem.jl:166
 [3] ODESystem
   @ ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/diffeqs/odesystem.jl:151 [inlined]
 [4] ODESystem(deqs::Vector{Equation}, iv::Num, dvs::Vector{Num}, ps::Vector{Any}; controls::Vector{Num}, observed::Vector{Equation}, systems::Vector{Any}, tspan::Nothing, name::Symbol, default_u0::Dict{Any, Any}, default_p::Dict{Any, Any}, defaults::Dict{Any, Any}, connector_type::Nothing, preface::Nothing, continuous_events::Nothing, discrete_events::Nothing, checks::Bool, metadata::Nothing, gui_metadata::ModelingToolkit.GUIMetadata)
   @ ModelingToolkit ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/diffeqs/odesystem.jl:226
 [5] __FirstOrder__(; name::Symbol, f::var"#20#21", x::Nothing)
   @ Main ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/model_parsing.jl:82
 [6] __FirstOrder__
   @ ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/model_parsing.jl:82 [inlined]
 [7] #_#172
   @ ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/model_parsing.jl:25 [inlined]
 [8] top-level scope
   @ ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/abstractsystem.jl:1003

But when I try to correct this by assigning a unit to the constant involved (initial code), I get:

ERROR: syntax: "begin" at /Users/mtiller/Source/modeling-toolkit-introduction/srcs/unit_error.jl:16 expected "end", got "u"
Stacktrace:
 [1] top-level scope
   @ ~/Source/modeling-toolkit-introduction/srcs/unit_error.jl:17

In this case, it is clear what the origin of the issue is. But why should this be an error? If literal values appear in equations, we may need to assign them unit information. What is the alternative?

mtiller-jh commented 11 months ago

A further refactoring into:


@mtkmodel FirstOrder begin
    @structural_parameters begin
        f
    end
    @parameters begin 
        c = 1.0, [unit = u"V/s^-1"]
    end
    @variables begin
        x(t), [unit = u"V"]
    end
    @equations begin
        D(x) ~ f(t) * c
    end
end

@mtkbuild root = FirstOrder(f=sin)
u0 = [
    root.x => 0
]
prob = ODEProblem(root, u0, (0, 2))
sol = solve(prob)
display(plot(sol, idxs=[root.x]))

...which fixes an issue with t not having units and removed the extraneous x variable and introduced a constant to try and balance the units still ends with:

┌ Warning:  in eq. #1right: no method matching sin for arguments (Quantity{Int64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}},).
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/validation.jl:154
ERROR: ModelingToolkit.ValidationError("Some equations had invalid units. See warnings for details.")
Stacktrace:
 [1] check_units(eqs::Vector{Equation})
   @ ModelingToolkit ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/validation.jl:277
 [2] ODESystem(tag::UInt64, deqs::Vector{Equation}, iv::SymbolicUtils.BasicSymbolic{Real}, dvs::Vector{SymbolicUtils.BasicSymbolic{Real}}, ps::Vector{SymbolicUtils.BasicSymbolic{Real}}, tspan::Nothing, var_to_name::Dict{Any, Any}, ctrls::Vector{Any}, observed::Vector{Equation}, tgrad::Base.RefValue{Vector{Num}}, jac::Base.RefValue{Any}, ctrl_jac::Base.RefValue{Any}, Wfact::Base.RefValue{Matrix{Num}}, Wfact_t::Base.RefValue{Matrix{Num}}, name::Symbol, systems::Vector{Any}, defaults::Dict{Any, Any}, torn_matching::Nothing, connector_type::Nothing, preface::Nothing, cevents::Vector{ModelingToolkit.SymbolicContinuousCallback}, devents::Vector{ModelingToolkit.SymbolicDiscreteCallback}, metadata::Nothing, gui_metadata::ModelingToolkit.GUIMetadata, tearing_state::Nothing, substitutions::Nothing, complete::Bool, discrete_subsystems::Nothing, unknown_states::Nothing, split_idxs::Nothing, parent::Nothing; checks::Bool)
   @ ModelingToolkit ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/diffeqs/odesystem.jl:166
 [3] ODESystem
   @ ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/diffeqs/odesystem.jl:151 [inlined]
 [4] ODESystem(deqs::Vector{Equation}, iv::Num, dvs::Vector{Num}, ps::Vector{Num}; controls::Vector{Num}, observed::Vector{Equation}, systems::Vector{Any}, tspan::Nothing, name::Symbol, default_u0::Dict{Any, Any}, default_p::Dict{Any, Any}, defaults::Dict{Any, Any}, connector_type::Nothing, preface::Nothing, continuous_events::Nothing, discrete_events::Nothing, checks::Bool, metadata::Nothing, gui_metadata::ModelingToolkit.GUIMetadata)
   @ ModelingToolkit ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/diffeqs/odesystem.jl:226
 [5] __FirstOrder__(; name::Symbol, f::typeof(sin), c::Nothing, x::Nothing)
   @ Main ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/model_parsing.jl:82
 [6] __FirstOrder__
   @ ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/model_parsing.jl:82 [inlined]
 [7] #_#172
   @ ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/model_parsing.jl:25 [inlined]
 [8] top-level scope
   @ ~/.julia/packages/ModelingToolkit/7CI1l/src/systems/abstractsystem.jl:1003