From 89826ea9b1676e38873e5f517b5eb2614f2960c9 Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Thu, 5 Sep 2024 13:53:42 +0200 Subject: [PATCH 1/3] add slew-rate limiter --- src/ModelingToolkitSampledData.jl | 2 +- src/discrete_blocks.jl | 64 +++++++++++++++++++++++++++++++ test/test_discrete_blocks.jl | 28 ++++++++++++++ 3 files changed, 93 insertions(+), 1 deletion(-) diff --git a/src/ModelingToolkitSampledData.jl b/src/ModelingToolkitSampledData.jl index 00aa621..47e3377 100644 --- a/src/ModelingToolkitSampledData.jl +++ b/src/ModelingToolkitSampledData.jl @@ -8,7 +8,7 @@ export DiscreteIntegrator, DiscreteDerivative, Delay, Difference, ZeroOrderHold, ClockChanger, DiscretePIDParallel, DiscretePIDStandard, DiscreteStateSpace, DiscreteTransferFunction, NormalNoise, UniformNoise, Quantization, - ExponentialFilter + DiscreteSlewRateLimiter, ExponentialFilter export DiscreteOnOffController include("discrete_blocks.jl") diff --git a/src/discrete_blocks.jl b/src/discrete_blocks.jl index aa57a5a..035ae8e 100644 --- a/src/discrete_blocks.jl +++ b/src/discrete_blocks.jl @@ -878,6 +878,70 @@ See also [ControlSystemsMTK.jl](https://juliacontrol.github.io/ControlSystemsMTK # compose(ODESystem(eqs, t, sts, pars; name = name), input, output) # end + +""" + DiscreteSlewRateLimiter(rate = 1.0, rate_negative = rate) + +A discrete-time slew rate limiter that limits the rate of change of the input signal. + +Note, the sample interval is not taken into account when computing the rate of change, the difference between two consequetive samples is saturated. + +# Parameters: +- `rate`: Slew rate limit (in positive/increasing direction). Must be a positive number. +- `rate_negative`: Negative slew rate limit, defaults to `rate`. Must be a positive number. + +# Variables +- `d`: Unsaturated rate of change of the input signal +- `u`: Input signal +- `y`: Output signal (saturated slew rate) + +# Connectors: +- `input` +- `output` + +# Example +``` +cl = Clock(0.1) +z = ShiftIndex(cl) +@mtkmodel SlweRateLimiterModel begin + @components begin + input = Sine(amplitude=1, frequency=0.8) + limiter = DiscreteSlewRateLimiter(; z, rate=0.4, rate_negative = 0.3) + end + @variables begin + x(t) = 0 # Dummy variable to workaround JSCompiler bug + end + @equations begin + connect(input.output, limiter.input) + D(x) ~ 0.1x + Hold(limiter.y) + end +end +@named m = SlweRateLimiterModel() +m = complete(m) +ssys = structural_simplify(IRSystem(m)) +prob = ODEProblem(ssys, [m.limiter.y(z-1) => 0], (0.0, 2.0)) +sol = solve(prob, Tsit5(), dtmax=0.01) +plot(sol, idxs=[m.input.output.u, m.limiter.y], title="Slew rate limited sine wave") +``` +""" +@mtkmodel DiscreteSlewRateLimiter begin + @extend u, y = siso = SISO() + @structural_parameters begin + z = ShiftIndex() + end + @parameters begin + rate = 1.0, [description = "Slew rate limit"] + rate_negative = rate, [description = "Negative slew rate limit"] + end + @variables begin + d(t) + end + @equations begin + d(z) ~ u(z) - y(z-1) + y(z) ~ y(z-1) + clamp(d(z), -rate_negative, rate) + end +end + """ Quantization diff --git a/test/test_discrete_blocks.jl b/test/test_discrete_blocks.jl index c6c0dbd..5b4e53c 100644 --- a/test/test_discrete_blocks.jl +++ b/test/test_discrete_blocks.jl @@ -385,7 +385,35 @@ end # @test reduce(vcat, sol((0:10) .+ 1e-2))[:]≈[zeros(2); 1; zeros(8)] atol=1e-2 # end* +## +@testset "SlweRateLimiter" begin + @info "Testing SlweRateLimiter" + cl = Clock(0.1) + z = ShiftIndex(cl) + @mtkmodel SlweRateLimiterModel begin + @components begin + input = Sine(amplitude=1, frequency=0.8) + limiter = DiscreteSlewRateLimiter(; z, rate=0.4, rate_negative = 0.3) + end + @variables begin + x(t) = 0 # Dummy variable to workaround JSCompiler bug + end + @equations begin + connect(input.output, limiter.input) + D(x) ~ 0.1x + Hold(limiter.y) + end + end + @named m = SlweRateLimiterModel() + m = complete(m) + ssys = structural_simplify(IRSystem(m)) + prob = ODEProblem(ssys, [m.limiter.y(z-1) => 0], (0.0, 2.0)) + sol = solve(prob, Tsit5(), dtmax=0.01) + plot(sol, idxs=[m.input.output.u, m.limiter.y], title="Slew rate limited sine wave") + @test maximum(diff(sol[m.limiter.y])) ≈ 0.4 + @test minimum(diff(sol[m.limiter.y])) ≈ -0.3 +end +## @testset "quantization" begin @info "Testing quantization" From 6d1d134915524901a2ff452ed000a2a7b57110b3 Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Thu, 5 Sep 2024 13:54:42 +0200 Subject: [PATCH 2/3] add to docs --- docs/src/blocks.md | 1 + 1 file changed, 1 insertion(+) diff --git a/docs/src/blocks.md b/docs/src/blocks.md index abe8f44..d6cd23c 100644 --- a/docs/src/blocks.md +++ b/docs/src/blocks.md @@ -9,6 +9,7 @@ - [`Difference`](@ref) - [`DiscreteDerivative`](@ref) - [`DiscreteIntegrator`](@ref) +- [`DiscreteSlewRateLimiter`](@ref) - [`Sampler`](@ref) - [`ZeroOrderHold`](@ref) From d80330f0e34a6c9936f16c6de11f069045126e8f Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Thu, 5 Sep 2024 14:07:56 +0200 Subject: [PATCH 3/3] rm plot --- test/test_discrete_blocks.jl | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/test/test_discrete_blocks.jl b/test/test_discrete_blocks.jl index 5b4e53c..c7536c6 100644 --- a/test/test_discrete_blocks.jl +++ b/test/test_discrete_blocks.jl @@ -386,11 +386,11 @@ end # end* ## -@testset "SlweRateLimiter" begin - @info "Testing SlweRateLimiter" +@testset "SlewRateLimiter" begin + @info "Testing SlewRateLimiter" cl = Clock(0.1) z = ShiftIndex(cl) - @mtkmodel SlweRateLimiterModel begin + @mtkmodel SlewRateLimiterModel begin @components begin input = Sine(amplitude=1, frequency=0.8) limiter = DiscreteSlewRateLimiter(; z, rate=0.4, rate_negative = 0.3) @@ -403,12 +403,12 @@ end D(x) ~ 0.1x + Hold(limiter.y) end end - @named m = SlweRateLimiterModel() + @named m = SlewRateLimiterModel() m = complete(m) ssys = structural_simplify(IRSystem(m)) prob = ODEProblem(ssys, [m.limiter.y(z-1) => 0], (0.0, 2.0)) sol = solve(prob, Tsit5(), dtmax=0.01) - plot(sol, idxs=[m.input.output.u, m.limiter.y], title="Slew rate limited sine wave") + # plot(sol, idxs=[m.input.output.u, m.limiter.y], title="Slew rate limited sine wave") @test maximum(diff(sol[m.limiter.y])) ≈ 0.4 @test minimum(diff(sol[m.limiter.y])) ≈ -0.3 end