SciML
Improve interactive system browsing and add system description
I am trying to improve system printing, interactive model inspection and discoverability in the REPL. My idea has been to show() a little more information about systems in the REPL:
- Model name (as before)
- Model description (new; a human-readable string that can provide a high-level description of what a system represents)
I think combined use of name + description opens the door to (but does not enforce) shorter variable names for more readable equations, while the description can "elaborate" when needed. For example, Ohm's law for the third resistor in a circuit could be written circuit.R3.V ~ circuit.R3.R * circuit.R3.I (instead of with R3 => resistor3), ideally with a description that circuit.R3 is a "Resistor".
- Subsystems and their descriptions (new; making it easier to read a model's hierarchy)
- Number of equations (as before, but now categorized by type, and including (optional) observed equations, too)
- Unknowns, parameters and their defaults and descriptions (as before)
- Everything is shown "lazily" such that if description/subsystems/unknowns/... is/are not present or empty, it is skipped
- When there are too many subsystems/unknowns/... to show, a hint is shown to the inspection method that gives them all
My intention is to make it easier to "browse" a system in the REPL, hopefully without being overwhelming.
Any thoughts on this? If you think something like this is useful, I think it would need some more work first.
Checklist
- [ ] Appropriate tests were added
- [x] Any code changes were done in a way that does not break public API
- [ ] All documentation related to code changes were updated
- [x] The new code follows the contributor guidelines, in particular the SciML Style Guide and COLPRAC.
- [ ] Any new documentation only uses public API
For example, consider a simplified example of the big RC circuit tutorial:
using ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as D
# Basic electrical components
@connector function Pin(; name)
@variables v(t) i(t) [connect = Flow]
description = "Pin"
ODESystem(Equation[], t, [v, i], []; name, description)
end
function Ground(; name)
@named g = Pin()
eqs = [g.v ~ 0]
description = "Ground connection"
compose(ODESystem(eqs, t, [], []; name, description), g)
end
function ConstantVoltage(; name)
@named p = Pin()
@named n = Pin()
pars = @parameters V
eqs = [V ~ p.v - n.v
0 ~ p.i + n.i]
description = "Constant voltage source"
compose(ODESystem(eqs, t, [], pars; name, description), p, n)
end
@connector function HeatPort(; name)
vars = @variables T(t) Qflow(t) [connect = Flow]
description = "Heat port"
ODESystem(Equation[], t, vars, []; name, description)
end
function HeatingResistor(; name)
@named p = Pin()
@named n = Pin()
@named h = HeatPort()
vars = @variables v(t) RTherm(t)
pars = @parameters R TAmbient α
eqs = [RTherm ~ R * (1 + α * (h.T - TAmbient))
v ~ p.i * RTherm
h.Qflow ~ -v * p.i
v ~ p.v - n.v
0 ~ p.i + n.i]
description = "Heating resistor"
compose(ODESystem(eqs, t, vars, pars; name, description), p, n, h)
end
function HeatCapacitor(; name)
pars = @parameters ρ V cp C = ρ * V * cp
@named h = HeatPort()
eqs = [
D(h.T) ~ h.Qflow / C
]
description = "Heat capacitor"
compose(ODESystem(eqs, t, [], pars; name, description), h)
end
function Capacitor(; name)
@named p = Pin()
@named n = Pin()
vars = @variables v(t) = 0.0
pars = @parameters C
eqs = [v ~ p.v - n.v
0 ~ p.i + n.i
D(v) ~ p.i / C]
description = "Capacitor"
compose(ODESystem(eqs, t, vars, pars; name, description), p, n)
end
function RCCircuit(S, G; name)
@named R = HeatingResistor()
@named C = Capacitor()
@named HC = HeatCapacitor()
eqs = [connect(S.p, R.p)
connect(R.n, C.p)
connect(C.n, S.n, G.g)
connect(R.h, HC.h)]
description = "A circuit with a resistor and capacitor"
compose(ODESystem(eqs, t; name, description), R, C, S, G, HC)
end
function BigRCCircuit(N; name)
@named S = ConstantVoltage()
@named G = Ground()
RCs = [RCCircuit(S, G; name = Symbol(:RC, i)) for i in 1:N]
vars = @variables E(t)
circuit = ODESystem(
[D(E) ~ sum(RC.R.h.Qflow for RC in RCs)], t, vars, [];
defaults = [
[RC.R.R for RC in RCs] .=> 10 .^ range(0, -4, N);
[RC.C.C for RC in RCs] .=> 10 .^ range(-3, 0, N);
S.V => 2.0
], initialization_eqs = [E ~ 0],
name, description = "A circuit consisting of parallell connected RC circuits"
)
compose(circuit, RCs)
end
@named bigRC = BigRCCircuit(50)
bigRC = structural_simplify(bigRC)
Now inspecting bigRC in the REPL shows
Model bigRC: A circuit consisting of parallell connected RC circuits
Subsystems (50):
RC1: A circuit with a resistor and capacitor
RC2: A circuit with a resistor and capacitor
RC3: A circuit with a resistor and capacitor
RC4: A circuit with a resistor and capacitor
RC5: A circuit with a resistor and capacitor
RC6: A circuit with a resistor and capacitor
RC7: A circuit with a resistor and capacitor
RC8: A circuit with a resistor and capacitor
RC9: A circuit with a resistor and capacitor
RC10: A circuit with a resistor and capacitor
⋮
see hierarchy(sys) for all
Equations (1051):
151 solvable … see equations(sys) for all
900 observed … see observed(sys) for all
Unknowns (151):
E(t) [defaults to 0.0]
RC1₊C₊v(t) [defaults to 0.0]
RC1₊HC₊h₊T(t)
RC2₊C₊v(t) [defaults to 0.0]
RC2₊HC₊h₊T(t)
RC3₊C₊v(t) [defaults to 0.0]
RC3₊HC₊h₊T(t)
RC4₊C₊v(t) [defaults to 0.0]
RC4₊HC₊h₊T(t)
RC5₊C₊v(t) [defaults to 0.0]
⋮
see unknowns(sys) for all
Parameters (450):
RC23₊HC₊C [defaults to RC23₊HC₊V*RC23₊HC₊cp*RC23₊HC₊ρ]
RC40₊C₊C [defaults to 0.244205]
RC9₊HC₊V
RC10₊R₊TAmbient
RC3₊HC₊C [defaults to RC3₊HC₊V*RC3₊HC₊cp*RC3₊HC₊ρ]
RC36₊HC₊C [defaults to RC36₊HC₊V*RC36₊HC₊cp*RC36₊HC₊ρ]
RC33₊R₊α
RC34₊HC₊cp
RC37₊HC₊C [defaults to RC37₊HC₊V*RC37₊HC₊cp*RC37₊HC₊ρ]
RC12₊HC₊C [defaults to RC12₊HC₊V*RC12₊HC₊cp*RC12₊HC₊ρ]
⋮
see parameters(sys) for all
This shows me, in English, which 5 components RC3 consists of.
@hersle @ChrisRackauckas This pull request is a good feature. I also wrote a module Explore using the API functions provided by MTK. Explore has two tool functions, Explore.sys() and Explore.model(). The former recursively checks all the details of ODESystem and its subsystems, while the latter checks the details of MTK.Model. These details will be printed to the specified file that can be overwritten. The file path variable is specified in Julia's startup.jl. The details of ODESystem are as follows:
Explore.sys: sys ...
sys, 5 subs [:fixedTemperature, :resistor, :capacitor, :source, :ground]
fixedTemperature, 1 subs [:port]
port, 0 subs
resistor, 3 subs [:p, :n, :conditionalHeatPort]
p, 0 subs
n, 0 subs
conditionalHeatPort, 1 subs [:heatPort]
heatPort, 0 subs
capacitor, 2 subs [:p, :n]
p, 0 subs
n, 0 subs
source, 2 subs [:p, :n]
p, 0 subs
n, 0 subs
ground, 1 subs [:p]
p, 0 subs
sys, 5 subs [:fixedTemperature, :resistor, :capacitor, :source, :ground]
Unknowns (2)
capacitor₊v(t) ▪ Unassigned ▪ Voltage drop of the two pins (= p.v - n.v)
resistor₊i(t) ▪ Unassigned ▪ Current flowing from pin p to pin n
Parameters (7)
fixedTemperature₊T ▪ 400.15 ▪ Fixed temperature at port
resistor₊R ▪ 1.0 ▪ Resistance at temperature T_ref
resistor₊T_ref ▪ 300.15 ▪ Reference temperature
resistor₊alpha ▪ 0 ▪ Temperature coefficient of resistance R_actual = R*(1 + alpha*(T_heatPort - T_ref))
resistor₊conditionalHeatPort₊T ▪ 293.15 ▪ Default environment temperature under which the model user is, if useHeatPort = false
capacitor₊C ▪ 1.0 ▪ Capacitance
source₊V ▪ 1.0 ▪ Value of constant voltage
Equations (2)
Differential(t)(capacitor₊v(t)) ~ capacitor₊i(t) / capacitor₊C
0 ~ -resistor₊v(t) + resistor₊R_actual(t)*resistor₊i(t)
......
capacitor, 2 subs [:p, :n]
Unknowns (2)
v(t) ▪ Unassigned ▪ Voltage drop of the two pins (= p.v - n.v)
i(t) ▪ Unassigned ▪ Current flowing from pin p to pin n
Parameters (1)
C ▪ 1.0 ▪ Capacitance
Equations (4)
v(t) ~ p₊v(t) - n₊v(t)
0 ~ n₊i(t) + p₊i(t)
i(t) ~ p₊i(t)
i(t) ~ C*Differential(t)(v(t))
p, 0 subs
Unknowns (2)
v(t) ▪ Potential at the pin
i(t) ▪ Current flowing into the pin
Parameters (0)
Equations (0)
.......
Equations_substitutions (25)
fixedTemperature₊port₊T(t) ~ fixedTemperature₊T
source₊v(t) ~ source₊V
resistor₊v(t) ~ -capacitor₊v(t) + source₊v(t)
resistor₊p₊i(t) ~ resistor₊i(t)
......
Equations_algebraic (15)
fixedTemperature (1)
port₊T(t) ~ T
resistor (6)
v(t) ~ p₊v(t) - n₊v(t)
0 ~ n₊i(t) + p₊i(t)
i(t) ~ p₊i(t)
R_actual(t) ~ R*(1 + (-T_ref + conditionalHeatPort₊T_heatPort(t))*alpha)
v(t) ~ R_actual(t)*i(t)
conditionalHeatPort₊LossPower(t) ~ v(t)*i(t)
.......
Equations_differential (1)
capacitor (1)
i(t) ~ C*Differential(t)(v(t))
Model bigRC: A circuit consisting of parallell connected RC circuits Subsystems (50): RC1: A circuit with a resistor and capacitor RC2: A circuit with a resistor and capacitor RC3: A circuit with a resistor and capacitor RC4: A circuit with a resistor and capacitor RC5: A circuit with a resistor and capacitor RC6: A circuit with a resistor and capacitor RC7: A circuit with a resistor and capacitor RC8: A circuit with a resistor and capacitor RC9: A circuit with a resistor and capacitor RC10: A circuit with a resistor and capacitor ⋮
@wang890 That's nice! Maybe too much for the default show() method, though. I noticed AbstractSystem already implements most of the AbstractTrees interface. I have tried to make this easier to discover by hinting to hierarchy(sys) in show(sys), so you can now do e.g. (continued from first example):
@named bigRC = BigRCCircuit(2)
hierarchy(bigRC; describe = true)
bigRC: A circuit consisting of parallell connected RC circuits
├─ RC1: A circuit with a resistor and capacitor
│ ├─ R: Heating resistor
│ │ ├─ p: Pin
│ │ ├─ n: Pin
│ │ └─ h: Heat port
│ ├─ C: Capacitor
│ │ ├─ p: Pin
│ │ └─ n: Pin
│ ├─ S: Constant voltage source
│ │ ├─ p: Pin
│ │ └─ n: Pin
│ ├─ G: Ground connection
│ │ └─ g: Pin
│ └─ HC: Heat capacitor
│ └─ h: Heat port
└─ RC2: A circuit with a resistor and capacitor
├─ R: Heating resistor
│ ├─ p: Pin
│ ├─ n: Pin
│ └─ h: Heat port
├─ C: Capacitor
│ ├─ p: Pin
│ └─ n: Pin
├─ S: Constant voltage source
│ ├─ p: Pin
│ └─ n: Pin
├─ G: Ground connection
│ └─ g: Pin
└─ HC: Heat capacitor
└─ h: Heat port
Question: what are "extra" equations? Before this PR (and still, I haven't changed this), show() counts
neqs = count(eq -> !(eq.lhs isa Connection), eqs)
next = n_extra_equations(sys)
and shows the equation count $neqs ($next). For example, continuing my 1st example with @named bigRC = BigRCCircuit(2) shows
Model bigRC with 25 (43) equations ...
Then I do equations(bigRC) and (to my surprise) get a
33-element Vector{Equation}: ...
I expected 25 or 43 elements (but exactly 8 of the 33 equations are of the form connect(..., ...), and 33-8=25).
Is there a more sensible way to do this? How about just
neqs = length(eqs)
or alternatively distinguish between
ncon = count(eq -> eq.lhs isa Connection, eqs)
neqs = length(eqs) - ncon
?
I am, yes. That just comes on top. I left them out of my question, since there are 0 (since it is before structural simplification). To be complete, I would like to print equation counts simply as
neqs = length(equations(sys))
nobs = length(observed(sys))
or alternatively
ncon = count(eq -> eq.lhs isa Connection, equations(sys))
neqs = length(equations(sys)) - ncon
nobs = length(observed(sys))
and get rid of n_extra_equations() altogether. Would that be ok?
My best guess is that n_extra_equations() counts how many equations all the connector equations expands to; i.e. equivalent to
function n_extra_equations(sys)
nexpand = length(equations(expand_connections(sys))) # number of equations after connection expansion
nnoconn = length(filter(eq -> !(eq.lhs isa Connection), equations(sys))) # number of non-connector equations before connection expansion
return nexpand - nnoconn
end
Why is the description cached?
Maybe that could be avoided. What would be a good and simple way to do it? Compute it on-demand with a function?
Compute it on-demand with a function?
Yeah make getdescription just a function that computes it, and call that during show.
I could foresee it being nice if getdescription returned a common SystemDescription that tabulates a bunch of information and has some nice properties, like show, DataFrame, writing information for serialization, etc. being very useful.
Hmm, in that case, wouldn't also the user need to pass the function? For example:
function Planet(; name)
@parameters m
description(sys) = "A planet with mass $m"
return ODESystem(Equation[], t, [], [m]; name, description)
end
I think that's less intuitive than passing just a string.
I think that's less intuitive than passing just a string.
Oh I thought you were saving the whole description 😅 . I see, if it's just a string that's fine.
This is the right direction. Let's go with it and see how it feels. It might need some adjustments but any minor changes I won't know until I've seen it a bit. Thanks @hersle!
