JuliaControl / ControlSystems.jl

A Control Systems Toolbox for Julia
https://juliacontrol.github.io/ControlSystems.jl/stable/
Other
503 stars 85 forks source link

OutOfMemoryError on step without time step size #880

Closed mzaffalon closed 9 months ago

mzaffalon commented 9 months ago

I am not sure whether I should report it here, so apologies if this is the wrong place.

A =
[ -6538.4         -40.38          0.0                  0.0;
   409.6               0.0                 0.0                  0.0;
     0.0               0.9164            -28.28       -20.0;
     0.0               0.0                20.0                  0.0]
B =
[ 3846.15    0.0;
    0.0             -288.0;
    0.0                0.0;
    0.0                0.0]
C =[ 0.0  0.0  0.0  1.0]
D =[ 0.0  0.0]
step(ss(A,B,C,D))

throws an OutOfMemoryError because the temporal step size is too small. Indeed size(step(ss(A,B,C,D)).x) = (4, 153847, 2), a bit of an overkill for this stable slow system.

baggepinnen commented 9 months ago

The problem here is that this system has one very fast pole, and some very fast poles. When no time-vector is provided, the default time vector tries to be long enough to capture the slow dynamics, but dense enough to capture the fast dynamics.

With the change in

Your step response looks like this image

baggepinnen commented 9 months ago

This system can be reduced to third order effectively using balreal:

sys = ss(A,B,C,D)
sysr, G, T = balreal(sys)
plot(step.([sys, sysr], 10))

the step response is almost identical: image

mzaffalon commented 9 months ago

Oh, I understand what you mean: the fast pole is removed.

julia> poles(sys)
4-element Vector{ComplexF64}:
  -6535.930676480235 + 0.0im
 -2.5308619813034015 + 0.0im
 -14.140000000000017 + 14.144270925007058im
 -14.140000000000017 - 14.144270925007058im

julia> poles(sysr)
3-element Vector{ComplexF64}:
 -2.5308619557301677 + 0.0im
 -14.139998503418344 + 14.144270060245073im
 -14.139998503418344 - 14.144270060245073im
baggepinnen commented 9 months ago

Yeah, the fast pole wasn't observable

julia> gram(sys, :o)
4×4 Matrix{Float64}:
 1.48543e-6  2.37117e-5   6.53234e-5  4.98637e-5
 2.37117e-5  0.000378507  0.00104483  0.000799158
 6.53234e-5  0.00104483   0.0176803   0.025
 4.98637e-5  0.000799158  0.025       0.0530303

julia> gram(sys, :o) |> svdvals
4-element Vector{Float64}:
 0.06599408956924761
 0.004809087557537785
 0.0002874945673967246
 3.685244848996841e-20