diff --git a/src/Blocks/Blocks.jl b/src/Blocks/Blocks.jl index 4c5d75b56..249d52f0a 100644 --- a/src/Blocks/Blocks.jl +++ b/src/Blocks/Blocks.jl @@ -11,7 +11,7 @@ D = Differential(t) export RealInput, RealOutput, SISO include("utils.jl") -export Gain, Sum, MatrixGain, Feedback, Add, Product, Division +export Gain, Sum, MatrixGain, Feedback, Add, Add3, Product, Division export Abs, Sign, Sqrt, Sin, Cos, Tan, Asin, Acos, Atan, Atan2, Sinh, Cosh, Tanh, Exp export Log, Log10 include("math.jl") diff --git a/src/Blocks/continuous.jl b/src/Blocks/continuous.jl index efda074a8..d700508c2 100644 --- a/src/Blocks/continuous.jl +++ b/src/Blocks/continuous.jl @@ -192,14 +192,12 @@ function PID(;name, k=1, Ti=false, Td=false, Nd=10, xi_start=0, xd_start=0) push!(eqs, connect(err_input, int.input)) push!(eqs, connect(int.output, addPID.input2)) else - push!(eqs, connect(err_input, Izero.input)) push!(eqs, connect(Izero.output, addPID.input2)) end if with_D push!(eqs, connect(err_input, der.input)) push!(eqs, connect(der.output, addPID.input3)) else - push!(eqs, connect(err_input, Dzero.input)) push!(eqs, connect(Dzero.output, addPID.input3)) end ODESystem(eqs, t, [], []; name=name, systems=sys) @@ -285,6 +283,7 @@ function LimPID(; name, k=1, Ti=false, Td=false, wp=1, wd=1, ) with_I = !isequal(Ti, false) with_D = !isequal(Td, false) + with_AWM = Ni != Inf if gains Ti = k / Ti Td = Td / k @@ -302,12 +301,16 @@ function LimPID(; name, k=1, Ti=false, Td=false, wp=1, wd=1, @named addP = Add(k1=wp, k2=-1) @named gainPID = Gain(k) @named addPID = Add3() + @named limiter = Limiter(y_max=u_max, y_min=u_min) if with_I - @named addI = Add3(k1=1, k2=-1, k3=1) + if with_AWM + @named addI = Add3(k1=1, k2=-1, k3=1) + @named addSat = Add(k1=1, k2=-1) + @named gainTrack = Gain(1/(k * Ni)) + else + @named addI = Add(k1=1, k2=-1) + end @named int = Integrator(k=1/Ti, x_start=xi_start) - @named limiter = Limiter(y_max=u_max, y_min=u_min) - @named addSat = Add(k1=1, k2=-1) - @named gainTrack = Gain(1/(k * Ni)) else @named Izero = Constant(k=0) end @@ -318,9 +321,12 @@ function LimPID(; name, k=1, Ti=false, Td=false, wp=1, wd=1, @named Dzero = Constant(k=0) end - sys = [reference, measurement, ctr_output, addP, gainPID, addPID] + sys = [reference, measurement, ctr_output, addP, gainPID, addPID, limiter] if with_I - push!(sys, [addI, int, limiter, addSat, gainTrack]...) + if with_AWM + push!(sys, [addSat, gainTrack]...) + end + push!(sys, [addI, int]...) else push!(sys, Izero) end @@ -335,16 +341,18 @@ function LimPID(; name, k=1, Ti=false, Td=false, wp=1, wd=1, connect(measurement, addP.input2), connect(addP.output, addPID.input1), connect(addPID.output, gainPID.input), + connect(gainPID.output, limiter.input), + connect(limiter.output, ctr_output), ] if with_I push!(eqs, connect(reference, addI.input1)) push!(eqs, connect(measurement, addI.input2)) - push!(eqs, connect(gainPID.output, limiter.input)) - push!(eqs, connect(limiter.output, ctr_output)) - push!(eqs, connect(limiter.input, addSat.input2)) - push!(eqs, connect(limiter.output, addSat.input1)) - push!(eqs, connect(addSat.output, gainTrack.input)) - push!(eqs, connect(gainTrack.output, addI.input3)) + if with_AWM + push!(eqs, connect(limiter.input, addSat.input2)) + push!(eqs, connect(limiter.output, addSat.input1)) + push!(eqs, connect(addSat.output, gainTrack.input)) + push!(eqs, connect(gainTrack.output, addI.input3)) + end push!(eqs, connect(addI.output, int.input)) push!(eqs, connect(int.output, addPID.input3)) else diff --git a/test/Blocks/continuous.jl b/test/Blocks/continuous.jl index cadc5fb08..1f865aa2c 100644 --- a/test/Blocks/continuous.jl +++ b/test/Blocks/continuous.jl @@ -17,7 +17,8 @@ an integrator with a constant input is often used together with the system under sys = structural_simplify(iosys) prob = ODEProblem(sys, Pair[int.x=>1.0], (0.0, 1.0)) sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 2 + @test all(sol[c.output.u] .≈ 1) + @test sol[int.output.u][end] .≈ 2 # expected solution end @testset "Derivative" begin @@ -34,7 +35,7 @@ end sys = structural_simplify(iosys) prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 10.0)) sol = solve(prob, Rodas4()) - @test isapprox(sol[source.output.u], sol[int.output.u], atol=1e-1) + @test all(isapprox.(sol[source.output.u], sol[int.output.u], atol=1e-1)) end @testset "PT1" begin @@ -88,6 +89,9 @@ end sys = structural_simplify(model) prob = ODEProblem(sys, Pair[], (0.0, 100.0)) sol = solve(prob, Rodas4()) + # initial condition + @test sol[ss.x[1]][1] ≈ 0 atol=1e-3 + @test sol[ss.x[2]][1] ≈ 0 atol=1e-3 # equilibrium point is at [1, 0] @test sol[ss.x[1]][end] ≈ 1 atol=1e-3 @test sol[ss.x[2]][end] ≈ 0 atol=1e-3 @@ -108,7 +112,8 @@ function Plant(;name, x_start=zeros(2)) end @testset "PI" begin - @named ref = Constant(; k=2) + re_val = 2 + @named ref = Constant(; k=re_val) @named pi_controller = PI(k=1, T=1) @named plant = Plant() @named fb = Feedback() @@ -125,11 +130,13 @@ end sys = structural_simplify(model) prob = ODEProblem(sys, Pair[], (0.0, 100.0)) sol = solve(prob, Rodas4()) - @test sol[ref.output.u - plant.output.u][end] ≈ 0 atol=1e-3 # zero control error after 100s + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s end @testset "PID" begin - @named ref = Constant(; k=2) + re_val = 2 + @named ref = Constant(; k=re_val) @named pid_controller = PID(k=3, Ti=0.5, Td=100) @named plant = Plant() @named fb = Feedback() @@ -146,7 +153,46 @@ end sys = structural_simplify(model) prob = ODEProblem(sys, Pair[], (0.0, 100.0)) sol = solve(prob, Rodas4()) - @test sol[ref.output.u - plant.output.u][end] ≈ 0 atol=1e-3 # zero control error after 100s + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s + + @testset "PI" begin + @named pid_controller = PID(k=3, Ti=0.5, Td=false) + @named model = ODESystem( + [ + connect(ref.output, fb.input1), + connect(plant.output, fb.input2), + connect(fb.output, pid_controller.err_input), + connect(pid_controller.ctr_output, plant.input), + ], + t, + systems=[pid_controller, plant, ref, fb] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[], (0.0, 100.0)) + sol = solve(prob, Rodas4()) + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s + end + + @testset "PD" begin + @named pid_controller = PID(k=10, Ti=false, Td=1) + @named model = ODESystem( + [ + connect(ref.output, fb.input1), + connect(plant.output, fb.input2), + connect(fb.output, pid_controller.err_input), + connect(pid_controller.ctr_output, plant.input), + ], + t, + systems=[pid_controller, plant, ref, fb] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[], (0.0, 100.0)) + sol = solve(prob, Rodas4()) + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] > 1 # without I there will be a steady-state error + end end @test_skip begin @@ -236,7 +282,8 @@ end end @testset "LimPI" begin - @named ref = Constant(; k=1) + re_val = 1 + @named ref = Constant(; k=re_val) @named pi_controller_lim = LimPI(k=3, T=0.5, u_max=1.5, u_min=-1.5, Ta=0.1) @named pi_controller = PI(k=3, T=0.5) @named sat = Limiter(y_max=1.5, y_min=-1.5) @@ -279,15 +326,19 @@ end sol = solve(prob, Rodas4()) end - @test sol[ref.output.u - plant.output.u][end] ≈ 0 atol=1e-3 # zero control error after 100s - @test sol_lim[ref.output.u - plant.output.u][end] ≈ 0 atol=1e-3 # zero control error after 100s + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test all(isapprox.(sol_lim[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s + @test sol_lim[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s + @test all(-1.5 .<= sol_lim[pi_controller_lim.ctr_output.u] .<= 1.5) # test limit # Plots.plot(sol; vars=[plant.output.u]) # without anti-windup measure # Plots.plot!(sol_lim; vars=[plant.output.u]) # with anti-windup measure end @testset "LimPID" begin - @named ref = Constant(; k=1) + re_val = 1 + @named ref = Constant(; k=re_val) @named pid_controller = LimPID(k=3, Ti=0.5, Td=100, u_max=1.5, u_min=-1.5, Ni=0.1/0.5) @named plant = Plant() @named model = ODESystem( @@ -304,5 +355,112 @@ end sol = solve(prob, Rodas4()) # Plots.plot(sol, vars=[plant.output.u, plant.input.u]) - @test sol[ref.output.u - plant.output.u][end] ≈ 0 atol=1e-3 # zero control error after 100s + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s + @test all(-1.5 .<= sol[pid_controller.ctr_output.u] .<= 1.5) # test limit + + @testset "PI" begin + @named pid_controller = LimPID(k=3, Ti=0.5, Td=false, u_max=1.5, u_min=-1.5, Ni=0.1/0.5) + @named model = ODESystem( + [ + connect(ref.output, pid_controller.reference), + connect(plant.output, pid_controller.measurement), + connect(pid_controller.ctr_output, plant.input), + ], + t, + systems=[pid_controller, plant, ref] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[], (0.0, 100.0)) + sol = solve(prob, Rodas4()) + + # Plots.plot(sol, vars=[plant.output.u, plant.input.u]) + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s + @test all(-1.5 .<= sol[pid_controller.ctr_output.u] .<= 1.5) # test limit + end + @testset "PD" begin + @named pid_controller = LimPID(k=10, Ti=false, Td=1, u_max=1.5, u_min=-1.5) + @named model = ODESystem( + [ + connect(ref.output, pid_controller.reference), + connect(plant.output, pid_controller.measurement), + connect(pid_controller.ctr_output, plant.input), + ], + t, + systems=[pid_controller, plant, ref] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[], (0.0, 100.0)) + sol = solve(prob, Rodas4()) + + # Plots.plot(sol, vars=[plant.output.u, plant.input.u]) + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] > 0.5 # without I there will be a steady-state error + @test all(-1.5 .<= sol[pid_controller.ctr_output.u] .<= 1.5) # test limit + end + @testset "set-point weights" begin + @testset "wp" begin + @named pid_controller = LimPID(k=3, Ti=0.5, Td=100, u_max=1.5, u_min=-1.5, Ni=0.1/0.5, wp=0, wd=1) + @named model = ODESystem( + [ + connect(ref.output, pid_controller.reference), + connect(plant.output, pid_controller.measurement), + connect(pid_controller.ctr_output, plant.input), + ], + t, + systems=[pid_controller, plant, ref] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[], (0.0, 100.0)) + sol = solve(prob, Rodas4()) + + # Plots.plot(sol, vars=[plant.output.u, plant.input.u]) + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + sol[pid_controller.addP.output.u] == -sol[pid_controller.measurement.u] + @test sol[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s + @test all(-1.5 .<= sol[pid_controller.ctr_output.u] .<= 1.5) # test limit + end + @testset "wd" begin + @named pid_controller = LimPID(k=3, Ti=0.5, Td=100, u_max=1.5, u_min=-1.5, Ni=0.1/0.5, wp=1, wd=0) + @named model = ODESystem( + [ + connect(ref.output, pid_controller.reference), + connect(plant.output, pid_controller.measurement), + connect(pid_controller.ctr_output, plant.input), + ], + t, + systems=[pid_controller, plant, ref] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[], (0.0, 100.0)) + sol = solve(prob, Rodas4()) + + # Plots.plot(sol, vars=[plant.output.u, plant.input.u]) + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s + sol[pid_controller.addD.output.u] == -sol[pid_controller.measurement.u] + @test all(-1.5 .<= sol[pid_controller.ctr_output.u] .<= 1.5) # test limit + end + end + @testset "PI without AWM" begin + @named pid_controller = LimPID(k=3, Ti=0.5, Td=false, u_max=1.5, u_min=-1.5, Ni=Inf) + @named model = ODESystem( + [ + connect(ref.output, pid_controller.reference), + connect(plant.output, pid_controller.measurement), + connect(pid_controller.ctr_output, plant.input), + ], + t, + systems=[pid_controller, plant, ref] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[], (0.0, 100.0)) + sol = solve(prob, Rodas4()) + + # Plots.plot(sol, vars=[plant.output.u, plant.input.u]) + @test all(isapprox.(sol[ref.output.u], re_val, atol=1e-3)) # check reference + @test sol[plant.output.u][end] ≈ re_val atol=1e-3 # zero control error after 100s + @test all(-1.5 .<= sol[pid_controller.ctr_output.u] .<= 1.5) # test limit + end end \ No newline at end of file diff --git a/test/Blocks/math.jl b/test/Blocks/math.jl index c6bf82fe7..f36ea9e1f 100644 --- a/test/Blocks/math.jl +++ b/test/Blocks/math.jl @@ -15,7 +15,8 @@ using ModelingToolkit, OrdinaryDiffEq sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 2 + @test all(sol[c.output.u] .≈ 1) + @test sol[int.output.u][end] ≈ 2 # expected solution after 1s end @testset "Feedback loop" begin @@ -38,12 +39,12 @@ end prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 100.0)) sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 2 + @test sol[int.output.u][end] ≈ 2 # expected solution after 1s end @testset "Add" begin @named c1 = Constant(; k=1) - @named c2 = Constant(; k=2) + @named c2 = Sine(; frequency=1) @named add = Add(;) @named int = Integrator(; k=1) @named model = ODESystem( @@ -56,16 +57,77 @@ end systems=[int, add, c1, c2] ) sys = structural_simplify(model) - prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) + sol = solve(prob, Rodas4()) + @test sol[add.output.u] ≈ 1 .+ sin.(2*pi*sol.t) + + @testset "weights" begin + k1 = -1 + k2 = 2 + @named add = Add(;k1=k1, k2=k2) + @named model = ODESystem( + [ + connect(c1.output, add.input1), + connect(c2.output, add.input2), + connect(add.output, int.input), + ], + t, + systems=[int, add, c1, c2] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) + sol = solve(prob, Rodas4()) + @test sol[add.output.u] ≈ k1 .* 1 .+ k2 .* sin.(2*pi*sol.t) + end +end +@testset "Add3" begin + @named c1 = Constant(; k=1) + @named c2 = Sine(; frequency=1) + @named c3 = Sine(; frequency=2) + @named add = Add3(;) + @named int = Integrator(; k=1) + @named model = ODESystem( + [ + connect(c1.output, add.input1), + connect(c2.output, add.input2), + connect(c3.output, add.input3), + connect(add.output, int.input), + ], + t, + systems=[int, add, c1, c2, c3] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 3 + @test sol[add.output.u] ≈ 1 .+ sin.(2*pi*sol.t) .+ sin.(2*pi*2*sol.t) + + @testset "weights" begin + k1 = -1 + k2 = 2 + k3 = -pi + @named add = Add3(;k1=k1, k2=k2, k3=k3) + @named model = ODESystem( + [ + connect(c1.output, add.input1), + connect(c2.output, add.input2), + connect(c3.output, add.input3), + connect(add.output, int.input), + ], + t, + systems=[int, add, c1, c2, c3] + ) + sys = structural_simplify(model) + prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) + sol = solve(prob, Rodas4()) + @test sol[add.output.u] ≈ k1 .* 1 .+ k2 .* sin.(2*pi*sol.t) .+ k3 .* sin.(2*pi*2*sol.t) + + end end @testset "Product" begin - @named c1 = Constant(; k=1) - @named c2 = Constant(; k=2) + @named c1 = Constant(; k=2) + @named c2 = Sine(; frequency=1) @named prod = Product(;) @named int = Integrator(; k=1) @named model = ODESystem( @@ -78,15 +140,13 @@ end systems=[int, prod, c1, c2] ) sys = structural_simplify(model) - prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) - sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 2 + @test sol[prod.output.u] ≈ 2 * sin.(2*pi*sol.t) end @testset "Division" begin - @named c1 = Constant(; k=1) + @named c1 = Sine(; frequency=1) @named c2 = Constant(; k=2) @named div = Division(;) @named int = Integrator(; k=1) @@ -100,71 +160,27 @@ end systems=[int, div, c1, c2] ) sys = structural_simplify(model) - prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) - sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 1/2 + @test sol[div.output.u] ≈ sin.(2*pi*sol.t) ./ 2 end @testset "Abs" begin - @named c = Constant(; k=-1) - @named abs = Abs(;) - @named int = Integrator(; k=1) - @named model = ODESystem( - [ - connect(c.output, abs.input), - connect(abs.output, int.input), - ], - t, - systems=[int, abs, c] - ) - sys = structural_simplify(model) - - prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) - - sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 1 -end - -@testset "Sqrt" begin - @named c = Constant(; k=4) - @named sqr = Sqrt(;) + @named c = Sine(; frequency=1) + @named absb = Abs(;) @named int = Integrator(; k=1) @named model = ODESystem( [ - connect(c.output, sqr.input), - connect(sqr.output, int.input), + connect(c.output, absb.input), + connect(absb.output, int.input), ], t, - systems=[int, sqr, c] + systems=[int, absb, c] ) sys = structural_simplify(model) - prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) - sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 2 -end - -@testset "Sign" begin - @named c = Constant(; k=3) - @named sig = Sign(;) - @named int = Integrator(; k=1) - @named model = ODESystem( - [ - connect(c.output, sig.input), - connect(sig.output, int.input), - ], - t, - systems=[int, sig, c] - ) - sys = structural_simplify(model) - - prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) - - sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 1 + @test sol[absb.output.u] ≈ abs.(sin.(2*pi*sol.t)) end @testset "MatrixGain" begin @@ -179,15 +195,13 @@ end end @testset "Math" begin - for (block, func) in [(Abs, abs), (Sin, sin), (Cos, cos), (Tan, tan), (Asin, asin), (Acos, acos), (Atan, atan), (Sinh, sinh), (Cosh, cosh), (Tanh, tanh), (Exp, exp)] + for (block, func) in [(Abs, abs), (Sign, sign), (Sin, sin), (Cos, cos), (Tan, tan), (Asin, asin), (Acos, acos), (Atan, atan), (Sinh, sinh), (Cosh, cosh), (Tanh, tanh), (Exp, exp)] @named source = Sine(frequency=1) @named b = block() @named int = Integrator() @named model = ODESystem([connect(source.output, b.input), connect(b.output, int.input)], t, systems=[int, b, source]) sys = structural_simplify(model) - prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) - sol = solve(prob, Rodas4()) @test sol[b.output.u] ≈ func.(sol[source.output.u]) end @@ -199,9 +213,7 @@ end @named int = Integrator() @named model = ODESystem([connect(source.output, b.input), connect(b.output, int.input)], t, systems=[int, b, source]) sys = structural_simplify(model) - prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) - sol = solve(prob, Rodas4()) @test sol[b.output.u] ≈ func.(sol[source.output.u]) end @@ -222,9 +234,7 @@ end systems=[int, b, c1, c2] ) sys = structural_simplify(model) - prob = ODEProblem(sys, Pair[int.x=>0.0], (0.0, 1.0)) - sol = solve(prob, Rodas4()) @test sol[int.output.u][end] ≈ atan(1, 2) end diff --git a/test/Blocks/nonlinear.jl b/test/Blocks/nonlinear.jl index 2c792930e..8957b97b6 100644 --- a/test/Blocks/nonlinear.jl +++ b/test/Blocks/nonlinear.jl @@ -66,7 +66,7 @@ end prob = ODEProblem(sys, [int.x=>1.0], (0.0, 1.0)) sol = solve(prob, Rodas4()) - @test sol[int.output.u][end] ≈ 2 + @test all(sol[int.output.u][end] .≈ 2) end @testset "Sine" begin @@ -108,6 +108,8 @@ end prob = ODEProblem(sys, Pair[], (0.0, 10.0)) - sol = solve(prob, Rodas4(), saveat=0.01, abstol=1e-10, reltol=1e-10) + tS = 0.01 + sol = solve(prob, Rodas4(), saveat=tS, abstol=1e-10, reltol=1e-10) @test all(abs.(sol[rl.output.u]) .<= 0.51) + @test all(-1 - 1e-5 .<= diff(sol[rl.output.u]) ./ tS .<= 1 + 1e-5) # just an approximation end \ No newline at end of file