ModelingToolkit - DiscreteUpdate

jondavis847 · February 7, 2023, 3:11am

Hello. I’m trying to implement discrete control in MTK using the DiscreteUpdate operator. It failed to update on my more complex model, seems to fail for a very basic MWE as well.

using ModelingToolkit, DifferentialEquations

@variables t x(t)=0 y(t)=0
D = Differential(t)
U = DiscreteUpdate(t; dt = 1)

eqs = [D(x) ~ 1
       U(y) ~ x]

@named de = ODESystem(eqs,t)
sys = structural_simplify(de)
prob = ODEProblem(sys,[],(0,10))
sol = solve(prob)

Results in

julia> sol[x]
7-element Vector{Float64}:
 0.0
 9.999999999999996e-5
 0.0010999999999999996
 0.011099999999999994
 0.11109999999999993
 1.1110999999999993
 9.999999999999996

julia> sol[y]
7-element Vector{Int64}:
 0
 0
 0
 0
 0
 0
 0

Am I doing something wrong?

I found Coupling MTK with time discrete control · Issue #1180 · SciML/ModelingToolkit.jl · GitHub

Looks like a workaround using PeriodicCallback was found back in 2021 but is this still the best way to implement discrete control in MTK? What would that look like for models that compose multiple ODESystems together, some discrete and some continuous, since the callback is set at the ODEProblem or solve() level?

I have a similar model developed with only DifferentialEquations and a PeriodicCallback, but MTK with DiscreteUpdate seems like a much better long term solution.

Thanks!

ChrisRackauckas · February 7, 2023, 11:01am

Open an issue, the mixed discrete stuff still has a few things WIP so any non-working example is useful…

baggepinnen · February 7, 2023, 11:39am

We’ve been working on and off for the last couple of months adding support for hybrid continuous/discrete systems. The functionality is still not in a usable state and is as a consequence of that not yet documented. You can get a feeling for the current capabilities by looking at this test file, there are a number of increasingly more complicated hybrid systems

github.com

SciML/ModelingToolkit.jl/blob/master/test/clock.jl#L17


      
              ModelingToolkit.infer_clocks!(ci), Dict(ci.ts.fullvars .=> ci.var_domain)
          end
          
          
@info "Testing hybrid system"
          dt = 0.1
          @variables t x(t) y(t) u(t) yd(t) ud(t) r(t)
          @parameters kp
          D = Differential(t)
          # u(n + 1) := f(u(n))
          
          
eqs = [yd ~ Sample(t, dt)(y)
                 ud ~ kp * (r - yd)
                 r ~ 1.0
          
          
       # plant (time continuous part)
                 u ~ Hold(ud)
                 D(x) ~ -x + u
                 y ~ x]
          @named sys = ODESystem(eqs)
          # compute equation and variables' time domains
          #TODO: test linearize

Please consider this functionality experimental at this stage, and subject to breaking changes. Here’s a TODO list for things that needs to fall into place before this stuff is properly usable.

jondavis847 · February 7, 2023, 7:04pm

Thank you gentlemen for the response! Glad to hear it’s on the radar and being worked. Looks like Sample(t,dt)(f) is the functionality I’m looking for, though it looks like it still needs some work. Perhaps I’ll just consider DiscreteUpdate OBE and update issue #1180 with any findings I have from experimenting with Sample? Or maybe just have patience and let you guys work your magic.