add sampletime operator

Author: baggepinnen

src/ModelingToolkit.jl

@@ -269,7 +269,7 @@ export debug_system
 #export Continuous, Discrete, sampletime, input_timedomain, output_timedomain
 #export has_discrete_domain, has_continuous_domain
 #export is_discrete_domain, is_continuous_domain, is_hybrid_domain
-export Sample, Hold, Shift, ShiftIndex
+export Sample, Hold, Shift, ShiftIndex, sampletime
 export Clock #, InferredDiscrete,
 
 end # module

src/clock.jl

@@ -117,7 +117,8 @@ end
 Clock(dt::Real) = Clock(nothing, dt)
 Clock() = Clock(nothing, nothing)
 
-sampletime(c) = isdefined(c, :dt) ? c.dt : nothing
+sampletime() = InferredSampleTime()
+sampletime(c) = something(isdefined(c, :dt) ? c.dt : nothing, InferredSampleTime())
 Base.hash(c::Clock, seed::UInt) = hash(c.dt, seed ⊻ 0x953d7a9a18874b90)
 function Base.:(==)(c1::Clock, c2::Clock)
     ((c1.t === nothing || c2.t === nothing) || isequal(c1.t, c2.t)) && c1.dt == c2.dt

src/discretedomain.jl

@@ -1,5 +1,14 @@
 using Symbolics: Operator, Num, Term, value, recursive_hasoperator
 
+struct InferredSampleTime <: Operator end
+function SymbolicUtils.promote_symtype(::Type{InferredSampleTime}, t...)
+    Real
+end
+function InferredSampleTime()
+    # Term{Real}(InferredSampleTime, Any[])
+    SymbolicUtils.term(InferredSampleTime, type = Real)
+end
+
 # Shift
 
 """
@@ -15,8 +24,6 @@ $(FIELDS)
 ```jldoctest
 julia> using Symbolics
 
-julia> @variables t;
-
 julia> Δ = Shift(t)
 (::Shift) (generic function with 2 methods)
 ```
@@ -176,16 +183,18 @@ end
 function (xn::Num)(k::ShiftIndex)
     @unpack clock, steps = k
     x = value(xn)
-    t = clock.t
     # Verify that the independent variables of k and x match and that the expression doesn't have multiple variables
     vars = Symbolics.get_variables(x)
     length(vars) == 1 ||
         error("Cannot shift a multivariate expression $x. Either create a new unknown and shift this, or shift the individual variables in the expression.")
     args = Symbolics.arguments(vars[]) # args should be one element vector with the t in x(t)
     length(args) == 1 ||
         error("Cannot shift an expression with multiple independent variables $x.")
-    isequal(args[], t) ||
-        error("Independent variable of $xn is not the same as that of the ShiftIndex $(k.t)")
+    t = args[]
+    if hasfield(typeof(clock), :t)
+        isequal(t, clock.t) ||
+            error("Independent variable of $xn is not the same as that of the ShiftIndex $(k.t)")
+    end
 
     # d, _ = propagate_time_domain(xn)
     # if d != clock # this is only required if the variable has another clock

src/systems/systemstructure.jl

@@ -652,6 +652,18 @@ function structural_simplify!(state::TearingState, io = nothing; simplify = fals
                 append!(appended_parameters, inputs[i], unknowns(ss))
                 discrete_subsystems[i] = ss
             end
+            for i in eachindex(discrete_subsystems)
+                discsys = discrete_subsystems[i]
+                eqs = collect(discsys.eqs)
+                for eqi in eachindex(eqs)
+                    clock = id_to_clock[i]
+                    clock isa AbstractDiscrete || continue
+                    Ts = sampletime(clock)
+                    eqs[eqi] = substitute(eqs[eqi], InferredSampleTime() => Ts)
+                end
+                @set discsys.eqs = eqs
+                discrete_subsystems[i] = discsys
+            end
             @set! sys.discrete_subsystems = discrete_subsystems, inputs, continuous_id,
             id_to_clock
             @set! sys.ps = appended_parameters

test/clock.jl

@@ -330,7 +330,7 @@ using ModelingToolkitStandardLibrary.Blocks
 
 dt = 0.05
 d = Clock(t, dt)
-k = ShiftIndex(d)
+k = ShiftIndex()
 
 @mtkmodel DiscretePI begin
     @components begin
@@ -347,7 +347,7 @@ k = ShiftIndex(d)
         y(t)
     end
     @equations begin
-        x(k) ~ x(k - 1) + ki * u(k)
+        x(k) ~ x(k - 1) + ki * u(k) * sampletime() / dt
         output.u(k) ~ y(k)
         input.u(k) ~ u(k)
         y(k) ~ x(k - 1) + kp * u(k)
@@ -364,21 +364,18 @@ end
     end
 end
 
-@mtkmodel Holder begin
-    @components begin
-        input = RealInput()
-        output = RealOutput()
-    end
+@mtkmodel ZeroOrderHold begin
+    @extend u, y = siso = Blocks.SISO()
     @equations begin
-        output.u ~ Hold(input.u)
+        y ~ Hold(u)
     end
 end
 
 @mtkmodel ClosedLoop begin
     @components begin
         plant = FirstOrder(k = 0.3, T = 1)
-        sampler = Sampler()
-        holder = Holder()
+        sampler = Blocks.Sampler(; clock = d)
+        holder = ZeroOrderHold()
         controller = DiscretePI(kp = 2, ki = 2)
         feedback = Feedback()
         ref = Constant(k = 0.5)