Skip to content

Commit

Permalink
Merge 8379266 into a57d602
Browse files Browse the repository at this point in the history
  • Loading branch information
olof3 committed Sep 11, 2020
2 parents a57d602 + 8379266 commit 454f858
Show file tree
Hide file tree
Showing 26 changed files with 220 additions and 254 deletions.
2 changes: 1 addition & 1 deletion Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -27,5 +27,5 @@ IterTools = "1.0"
LaTeXStrings = "1.0"
OrdinaryDiffEq = "5.2"
Plots = "0.24, 0.25, 0.26, 0.27, 0.28, 0.29, 1.0"
Polynomials = "0.6.0, 0.7"
Polynomials = "0.6.0"
julia = "1.0"
6 changes: 3 additions & 3 deletions docs/src/man/creating_systems.md
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@ The syntax for creating a transfer function is
```julia
tf(num, den, Ts=0)
```
where `num` and `den` are the polinomial coefficients of the numerator and denominator of the polynomial and `Ts` is the sample time.
where `num` and `den` are the polynomial coefficients of the numerator and denominator of the polynomial and `Ts` is the sample time.
### Example:
```julia
tf([1.0],[1,2,1])
Expand All @@ -29,7 +29,7 @@ Continuous-time transfer function model
The transfer functions created using this method will be of type `TransferFunction{SisoRational}`.

## zpk - Pole-Zero-Gain Representation
Sometimes it's better to represent the transferfunction by its poles, zeros and gain, this can be done using
Sometimes it's better to represent the transfer function by its poles, zeros and gain, this can be done using
```julia
zpk(zeros, poles, gain, Ts=0)
```
Expand Down Expand Up @@ -72,7 +72,7 @@ A state-space system is created using
```julia
ss(A,B,C,D,Ts=0)
```
and they behave similarily to transfer functions. State-space systems with heterogeneous matrix types are also available, which can be used to create systems with static or sized matrices, e.g.,
and they behave similarly to transfer functions. State-space systems with heterogeneous matrix types are also available, which can be used to create systems with static or sized matrices, e.g.,
```jldoctest HSS
using StaticArrays
import ControlSystems.HeteroStateSpace
Expand Down
23 changes: 3 additions & 20 deletions src/ControlSystems.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,6 @@ export LTISystem,
DelayLtiSystem,
Continuous,
Discrete,
Static,
ss,
tf,
zpk,
Expand Down Expand Up @@ -88,10 +87,8 @@ export LTISystem,
denvec,
numpoly,
denpoly,
sampletime,
iscontinuous,
isdiscrete,
isstatic
isdiscrete


# QUESTION: are these used? LaTeXStrings, Requires, IterTools
Expand All @@ -112,8 +109,8 @@ abstract type AbstractSystem end

include("types/TimeType.jl")
## Added interface:
# sampletime(Lti) -> Number
# sampletype(Lti) -> TimeType
# time(Lti) -> TimeType (not exported)


include("types/Lti.jl")

Expand Down Expand Up @@ -174,20 +171,6 @@ function covar(D::Union{AbstractMatrix,UniformScaling}, R)
D*R*D'
end

function Base.getproperty(sys::Union{StateSpace,HeteroStateSpace,TransferFunction}, s::Symbol)
if s === :Ts
# if !isdiscrete(sys) # NOTE this line seems to be breaking inference of isdiscrete
if !isdiscrete(sys.time) # We use this instead
@warn "Getting sampletime 0.0 for non-discrete systems is deprecated. Check `isdiscrete` before trying to access sampletime."
return 0.0
else
return sampletime(sys)
end
else
return getfield(sys, s)
end
end

# The path has to be evaluated upon initial import
const __CONTROLSYSTEMS_SOURCE_DIR__ = dirname(Base.source_path())

Expand Down
9 changes: 4 additions & 5 deletions src/analysis.jl
Original file line number Diff line number Diff line change
Expand Up @@ -118,7 +118,7 @@ Compute the zeros, poles, and gains of system `sys`.
`k` : Matrix{Float64}, (ny x nu)"""
function zpkdata(sys::LTISystem)
G = convert(TransferFunction{sampletype(sys),SisoZpk}, sys)
G = convert(TransferFunction{typeof(time(sys)),SisoZpk}, sys)

zs = map(x -> x.z, G.matrix)
ps = map(x -> x.p, G.matrix)
Expand All @@ -133,9 +133,8 @@ Compute the natural frequencies, `Wn`, and damping ratios, `zeta`, of the
poles, `ps`, of `sys`"""
function damp(sys::LTISystem)
ps = pole(sys)
if !iscontinuous(sys)
Ts = isstatic(sys) ? 1 : sampletime(sys)
ps = log(ps)/Ts
if isdiscrete(sys)
ps = log(ps)/sys.Ts
end
Wn = abs.(ps)
order = sortperm(Wn; by=z->(abs(z), real(z), imag(z)))
Expand Down Expand Up @@ -447,7 +446,7 @@ function delaymargin(G::LTISystem)
ωϕₘ = m[3][i]
dₘ = ϕₘ/ωϕₘ
if isdiscrete(G)
dₘ /= sampletime(G) # Give delay margin in number of sample times, as matlab does
dₘ /= G.Ts # Give delay margin in number of sample times, as matlab does
end
dₘ
end
Expand Down
32 changes: 16 additions & 16 deletions src/connections.jl
Original file line number Diff line number Diff line change
Expand Up @@ -22,18 +22,18 @@ Append systems in block diagonal form
"""
function append(systems::(ST where ST<:AbstractStateSpace)...)
ST = promote_type(typeof.(systems)...)
Ts = common_sample_time(s.time for s in systems)
time = common_time(systems...)
A = blockdiag([s.A for s in systems]...)
B = blockdiag([s.B for s in systems]...)
C = blockdiag([s.C for s in systems]...)
D = blockdiag([s.D for s in systems]...)
return ST(A, B, C, D, Ts)
return ST(A, B, C, D, time)
end

function append(systems::TransferFunction...)
Ts = common_sample_time(s.time for s in systems)
time = common_time(systems...)
mat = blockdiag([s.matrix for s in systems]...)
return TransferFunction(mat, Ts)
return TransferFunction(mat, time)
end

append(systems::LTISystem...) = append(promote(systems...)...)
Expand Down Expand Up @@ -62,8 +62,8 @@ function Base.vcat(systems::ST...) where ST <: AbstractStateSpace
B = vcat([s.B for s in systems]...)
C = blockdiag([s.C for s in systems]...)
D = vcat([s.D for s in systems]...)
Ts = common_sample_time(s.time for s in systems)
return ST(A, B, C, D, Ts)
time = common_time(systems...)
return ST(A, B, C, D, time)
end

function Base.vcat(systems::TransferFunction...)
Expand All @@ -72,10 +72,10 @@ function Base.vcat(systems::TransferFunction...)
if !all(s.nu == nu for s in systems)
error("All systems must have same input dimension")
end
Ts = common_sample_time(s.time for s in systems)
time = common_time(systems...)
mat = vcat([s.matrix for s in systems]...)

return TransferFunction(mat, Ts)
return TransferFunction(mat, time)
end

Base.vcat(systems::LTISystem...) = vcat(promote(systems...)...)
Expand All @@ -86,13 +86,13 @@ function Base.hcat(systems::ST...) where ST <: AbstractStateSpace
if !all(s.ny == ny for s in systems)
error("All systems must have same output dimension")
end
Ts = common_sample_time(s.time for s in systems)
time = common_time(systems...)
A = blockdiag([s.A for s in systems]...)
B = blockdiag([s.B for s in systems]...)
C = hcat([s.C for s in systems]...)
D = hcat([s.D for s in systems]...)

return ST(A, B, C, D, Ts)
return ST(A, B, C, D, time)
end

function Base.hcat(systems::TransferFunction...)
Expand All @@ -101,10 +101,10 @@ function Base.hcat(systems::TransferFunction...)
if !all(s.ny == ny for s in systems)
error("All systems must have same output dimension")
end
Ts = common_sample_time(s.time for s in systems)
time = common_time(systems...)
mat = hcat([s.matrix for s in systems]...)

return TransferFunction(mat, Ts)
return TransferFunction(mat, time)
end

Base.hcat(systems::LTISystem...) = hcat(promote(systems...)...)
Expand Down Expand Up @@ -201,13 +201,13 @@ function feedback(sys::Union{StateSpace, DelayLtiSystem})
end

function feedback(sys1::StateSpace,sys2::StateSpace)
Ts = common_sample_time(sys1.time,sys2.time)
time = common_time(sys1,sys2)
!(iszero(sys1.D) || iszero(sys2.D)) && error("There cannot be a direct term (D) in both sys1 and sys2")
A = [sys1.A+sys1.B*(-sys2.D)*sys1.C sys1.B*(-sys2.C);
sys2.B*sys1.C sys2.A+sys2.B*sys1.D*(-sys2.C)]
B = [sys1.B; sys2.B*sys1.D]
C = [sys1.C sys1.D*(-sys2.C)]
ss(A, B, C, sys1.D, Ts)
ss(A, B, C, sys1.D, time)
end


Expand All @@ -230,7 +230,7 @@ See Zhou etc. for similar (somewhat less symmetric) formulas.
U1=:, Y1=:, U2=:, Y2=:, W1=:, Z1=:, W2=Int[], Z2=Int[],
Wperm=:, Zperm=:, pos_feedback::Bool=false)

Ts = common_sample_time(sys1.time,sys2.time)
time = common_time(sys1,sys2)

if !(isa(Y1, Colon) || allunique(Y1)); @warn "Connecting single output to multiple inputs Y1=$Y1"; end
if !(isa(Y2, Colon) || allunique(Y2)); @warn "Connecting single output to multiple inputs Y2=$Y2"; end
Expand Down Expand Up @@ -299,7 +299,7 @@ See Zhou etc. for similar (somewhat less symmetric) formulas.
s2_D12*R2*s1_D21 s2_D11 + α*s2_D12*R2*s1_D22*s2_D21]
end

return StateSpace(A, B[:, Wperm], C[Zperm,:], D[Zperm, Wperm], Ts)
return StateSpace(A, B[:, Wperm], C[Zperm,:], D[Zperm, Wperm], time)
end


Expand Down
2 changes: 1 addition & 1 deletion src/delay_systems.jl
Original file line number Diff line number Diff line change
Expand Up @@ -161,7 +161,7 @@ function _bounds_and_features(sys::DelayLtiSystem, plot::Val)
end

# Againm we have to do something for default vectors, more or less a copy from timeresp.jl
function _default_Ts(sys::DelayLtiSystem)
function _default_dt(sys::DelayLtiSystem)
if !isstable(sys.P.P)
return 0.05 # Something small
else
Expand Down
2 changes: 1 addition & 1 deletion src/discrete.jl
Original file line number Diff line number Diff line change
Expand Up @@ -222,7 +222,7 @@ function lsima(sys::StateSpace, t::AbstractVector, r::AbstractVector, control_si
if iscontinuous(sys)
dsys = c2d(sys, dt, :zoh)[1]
else
if sampletime(sys) != dt
if sys.Ts != dt
error("Time vector must match sample time for discrete system")
end
dsys = sys
Expand Down
15 changes: 7 additions & 8 deletions src/freqresp.jl
Original file line number Diff line number Diff line change
Expand Up @@ -7,11 +7,10 @@ Evaluate the frequency response of a linear system
of system `sys` over the frequency vector `w`."""
function freqresp(sys::LTISystem, w_vec::AbstractVector{<:Real})
# Create imaginary freq vector s
if !iscontinuous(sys)
Ts = isstatic(sys) ? 1.0 : sampletime(sys)
s_vec = exp.(w_vec*(im*Ts))
else
if iscontinuous(sys)
s_vec = im*w_vec
else
s_vec = exp.(w_vec*(im*sys.Ts))
end
if isa(sys, StateSpace)
sys = _preprocess_for_freqresp(sys)
Expand Down Expand Up @@ -65,7 +64,7 @@ end
Notation for frequency response evaluation.
- F(s) evaluates the continuous-time transfer function F at s.
- F(omega,true) evaluates the discrete-time transfer function F at exp(i*Ts*omega)
- F(omega,true) evaluates the discrete-time transfer function F at exp(im*Ts*omega)
- F(z,false) evaluates the discrete-time transfer function F at z
"""
function (sys::TransferFunction)(s)
Expand All @@ -75,7 +74,7 @@ end
function (sys::TransferFunction)(z_or_omega::Number, map_to_unit_circle::Bool)
@assert isdiscrete(sys) "It only makes no sense to call this function with discrete systems"
if map_to_unit_circle
isreal(z_or_omega) ? evalfr(sys,exp(im*z_or_omega.*sampletime(sys))) : error("To map to the unit circle, omega should be real")
isreal(z_or_omega) ? evalfr(sys,exp(im*z_or_omega.*sys.Ts)) : error("To map to the unit circle, omega should be real")
else
evalfr(sys,z_or_omega)
end
Expand Down Expand Up @@ -169,8 +168,8 @@ function _bounds_and_features(sys::LTISystem, plot::Val)
w1 = 0.0
w2 = 2.0
end
if !iscontinuous(sys) && !isstatic(sys) # Do not draw above Nyquist freq for disc. systems
w2 = min(w2, log10/sampletime(sys)))
if isdiscrete(sys) # Do not draw above Nyquist freq for disc. systems
w2 = min(w2, log10/sys.Ts))
end
return [w1, w2], zp
end
17 changes: 6 additions & 11 deletions src/matrix_comps.jl
Original file line number Diff line number Diff line change
Expand Up @@ -241,13 +241,8 @@ state space systems in continuous and discrete time', American Control Conferenc
See also [`linfnorm`](@ref).
"""
function hinfnorm(sys::AbstractStateSpace; tol=1e-6)
if !isdiscrete(sys) # Continuous or Static
return _infnorm_two_steps_ct(sys, :hinf, tol)
else
return _infnorm_two_steps_dt(sys, :hinf, tol)
end
end
hinfnorm(sys::AbstractStateSpace{<:Continuous}; tol=1e-6) = _infnorm_two_steps_ct(sys, :hinf, tol)
hinfnorm(sys::AbstractStateSpace{<:Discrete}; tol=1e-6) = _infnorm_two_steps_dt(sys, :hinf, tol)
hinfnorm(sys::TransferFunction; tol=1e-6) = hinfnorm(ss(sys); tol=tol)

"""
Expand Down Expand Up @@ -374,8 +369,8 @@ function _infnorm_two_steps_dt(sys::AbstractStateSpace, normtype::Symbol, tol=1e

on_unit_circle = z -> abs(abs(z) - 1) < approxcirc # Helper fcn for readability

T = promote_type(real(numeric_type(sys)), Float64, typeof(true/sampletime(sys)))
Tw = typeof(one(T)/sampletime(sys))
T = promote_type(real(numeric_type(sys)), Float64, typeof(true/sys.Ts))
Tw = typeof(one(T)/sys.Ts)

if sys.nx == 0 # static gain
return (T(opnorm(sys.D)), Tw(0))
Expand All @@ -386,7 +381,7 @@ function _infnorm_two_steps_dt(sys::AbstractStateSpace, normtype::Symbol, tol=1e
# Check if there is a pole on the unit circle
pidx = findfirst(on_unit_circle, pole_vec)
if !(pidx isa Nothing)
return T(Inf), Tw(angle(pole_vec[pidx])/sampletime(sys))
return T(Inf), Tw(angle(pole_vec[pidx])/sys.Ts)
end

if normtype == :hinf && any(z -> abs(z) > 1, pole_vec)
Expand Down Expand Up @@ -439,7 +434,7 @@ function _infnorm_two_steps_dt(sys::AbstractStateSpace, normtype::Symbol, tol=1e
sort!(θ_vec)

if isempty(θ_vec)
return T((1+tol)*lb), Tw(θ_peak/sampletime(sys))
return T((1+tol)*lb), Tw(θ_peak/sys.Ts)
end

# Improve the lower bound
Expand Down

0 comments on commit 454f858

Please sign in to comment.