Skip to content

Commit 6c3b3aa

Browse files
Merge pull request #2517 from SciML/fb/continuousclock
feat: add `SolverStepClock`
2 parents a64aad8 + 7835e94 commit 6c3b3aa

File tree

3 files changed

+67
-0
lines changed

3 files changed

+67
-0
lines changed

src/clock.jl

Lines changed: 23 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -124,3 +124,26 @@ function Base.:(==)(c1::Clock, c2::Clock)
124124
end
125125

126126
is_concrete_time_domain(x) = x isa Union{AbstractClock, Continuous}
127+
128+
"""
129+
SolverStepClock <: AbstractClock
130+
SolverStepClock()
131+
SolverStepClock(t)
132+
133+
A clock that ticks at each solver step (sometimes referred to as "continuous sample time"). This clock **does generally not have equidistant tick intervals**, instead, the tick interval depends on the adaptive step-size slection of the continuous solver, as well as any continuous event handling. If adaptivity of the solver is turned off and there are no continuous events, the tick interval will be given by the fixed solver time step `dt`.
134+
135+
Due to possibly non-equidistant tick intervals, this clock should typically not be used with discrete-time systems that assume a fixed sample time, such as PID controllers and digital filters.
136+
"""
137+
struct SolverStepClock <: AbstractClock
138+
"Independent variable"
139+
t::Union{Nothing, Symbolic}
140+
"Period"
141+
SolverStepClock(t::Union{Num, Symbolic}) = new(value(t))
142+
end
143+
SolverStepClock() = SolverStepClock(nothing)
144+
145+
sampletime(c) = nothing
146+
Base.hash(c::SolverStepClock, seed::UInt) = seed 0x953d7b9a18874b91
147+
function Base.:(==)(c1::SolverStepClock, c2::SolverStepClock)
148+
((c1.t === nothing || c2.t === nothing) || isequal(c1.t, c2.t))
149+
end

src/systems/diffeqs/abstractodesystem.jl

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1053,6 +1053,10 @@ function DiffEqBase.ODEProblem{iip, specialize}(sys::AbstractODESystem, u0map =
10531053
discrete_cbs = map(affects, clocks, svs) do affect, clock, sv
10541054
if clock isa Clock
10551055
PeriodicCallback(DiscreteSaveAffect(affect, sv), clock.dt)
1056+
elseif clock isa SolverStepClock
1057+
affect = DiscreteSaveAffect(affect, sv)
1058+
DiscreteCallback(Returns(true), affect,
1059+
initialize = (c, u, t, integrator) -> affect(integrator))
10561060
else
10571061
error("$clock is not a supported clock type.")
10581062
end

test/clock.jl

Lines changed: 40 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -477,3 +477,43 @@ y = res.y[:]
477477
# y ~ x(k-1) + kp * u
478478
# Instead. This should be equivalent to the above, but gve me an error when I tried
479479
end
480+
481+
## Test continuous clock
482+
483+
c = ModelingToolkit.SolverStepClock(t)
484+
k = ShiftIndex(c)
485+
486+
@mtkmodel CounterSys begin
487+
@variables begin
488+
count(t) = 0
489+
u(t) = 0
490+
end
491+
@equations begin
492+
count(k + 1) ~ Sample(c)(u)
493+
end
494+
end
495+
496+
@mtkmodel FirstOrderSys begin
497+
@variables begin
498+
x(t) = 0
499+
end
500+
@equations begin
501+
D(x) ~ -x + sin(t)
502+
end
503+
end
504+
505+
@mtkmodel FirstOrderWithStepCounter begin
506+
@components begin
507+
counter = CounterSys()
508+
fo = FirstOrderSys()
509+
end
510+
@equations begin
511+
counter.u ~ fo.x
512+
end
513+
end
514+
515+
@mtkbuild model = FirstOrderWithStepCounter()
516+
prob = ODEProblem(model, [], (0.0, 10.0))
517+
sol = solve(prob, Tsit5(), kwargshandle = KeywordArgSilent)
518+
519+
@test sol.prob.kwargs[:disc_saved_values][1].t == sol.t[1:2:end] # Test that the discrete-tiem system executed at every step of the continuous solver. The solver saves each time step twice, one state value before discrete affect and one after.

0 commit comments

Comments
 (0)