# Time and Frequency response analysis

`ControlSystems.Simulator`

`ControlSystems.Simulator`

`ControlSystemsBase.BodemagWorkspace`

`ControlSystemsBase.LsimWorkspace`

`ControlSystemsBase.StepInfo`

`ControlSystemsBase.TransferFunction`

`Base.step`

`ControlSystemsBase.bode`

`ControlSystemsBase.bodemag!`

`ControlSystemsBase.evalfr`

`ControlSystemsBase.freqresp`

`ControlSystemsBase.freqresp!`

`ControlSystemsBase.impulse`

`ControlSystemsBase.lsim`

`ControlSystemsBase.lsim!`

`ControlSystemsBase.nyquist`

`ControlSystemsBase.sigma`

`ControlSystemsBase.stepinfo`

## Frequency response

Frequency responses are calculated using `freqresp`

, `bode`

, `sigma`

and `nyquist`

. Frequency-response plots are obtained using `bodeplot`

, `sigmaplot`

, `nyquistplot`

, `marginplot`

and `nicholsplot`

.

Any `TransferFunction`

can be evaluated at a point using `F(s)`

, `F(omega, true)`

, `F(z, false)`

`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(z,false)`

evaluates the discrete-time transfer function`F`

at`z`

A video demonstrating frequency-response analysis in ControlSystems.jl is available below.

## Time response (simulation)

Simulation with arbitrary inputs is primarily handled by the function `lsim`

, with `step`

and `impulse`

serving as convenience functions to simulate responses to particular inputs.

The function `lsim`

can take an input vector `u`

containing a sampled input trajectory, or an input function taking the state and time as arguments, `u(x,t)`

. This function can be used to easily simulate, e.g., ramp responses or saturated state-feedback control etc. See the docstring of `lsim`

for more details.

For more extensive nonlinear simulation capabilities, see the notes on ModelingToolkit and DifferentialEquations under The wider Julia ecosystem for control.

#### Example step response:

The following simulates a step response of a second-order system and plots the result.

```
using ControlSystemsBase, Plots
G = tf(1, [1, 1, 1])
res = step(G, 20) # Simulate 20 seconds step response
plot(res)
```

Using the function `stepinfo`

, we can compute characteristics of a step response:

`si = stepinfo(res)`

```
StepInfo:
Initial value: 0.000
Final value: 1.000
Step size: 1.000
Peak: 1.163
Peak time: 3.652 s
Overshoot: 16.30 %
Undershoot: 0.00 %
Settling time: 8.134 s
Rise time: 1.660 s
```

We can also plot the `StepInfo`

object

`plot(si)`

#### Example `lsim`

:

The function `lsim`

can take the control input as either

- An array of equidistantly sampled values, in this case the argument
`u`

is expected to have the shape`nu × n_time`

- A function of the state and time
`u(x,t)`

. This form allows simulation of state feedback, a step response at time $t_0$:`u(x, t) = amplitude * (t > t0)`

, or a ramp response:`u(x, t) = t`

etc.

The example below simulates state feedback with a step disturbance at $t=4$ by providing the function `u(x,t) = -L*x .+ (t > 4)`

to `lsim`

:

```
using ControlSystems
using LinearAlgebra: I
using Plots
A = [0 1; 0 0]
B = [0;1]
C = [1 0]
sys = ss(A,B,C,0)
Q = I
R = I
L = lqr(sys,Q,R)
u(x,t) = -L*x .+ (t > 4) # State feedback + step disturbance
t = 0:0.1:12
x0 = [1,0]
y, t, x, uout = lsim(sys,u,t,x0=x0)
plot(t,x', lab=["Position" "Velocity"], xlabel="Time [s]"); vline!([4], lab="Step disturbance", l=(:black, :dash, 0.5))
```

A video demonstrating time-domain simulation in ControlSystems.jl is available below.

## Docstrings

`ControlSystems.Simulator`

— Type`Simulator`

**Fields:**

```
P::StateSpace
f = (x,p,t) -> x
y = (x,t) -> y
```

`ControlSystems.Simulator`

— Method`Simulator(P::StateSpace, u = (x,t) -> 0)`

Used to simulate continuous-time systems. See function `?solve`

for additional info.

**Usage:**

```
using OrdinaryDiffEq, Plots
dt = 0.1
tfinal = 20
t = 0:dt:tfinal
P = ss(tf(1,[2,1])^2)
K = 5
reference(x,t) = [1.]
s = Simulator(P, reference)
x0 = [0.,0]
tspan = (0.0,tfinal)
sol = solve(s, x0, tspan, Tsit5())
plot(t, s.y(sol, t)[:], lab="Open loop step response")
```

`Base.step`

— Method```
y, t, x = step(sys[, tfinal])
y, t, x = step(sys[, t])
```

Calculate the response of the system `sys`

to a unit step at time `t = 0`

. If the final time `tfinal`

or time vector `t`

is not provided, one is calculated based on the system pole locations.

The return value is a structure of type `SimResult`

. A `SimResul`

can be plotted by `plot(result)`

, or destructured as `y, t, x = result`

.

`y`

has size `(ny, length(t), nu)`

, `x`

has size `(nx, length(t), nu)`

`ControlSystemsBase.impulse`

— Method```
y, t, x = impulse(sys[, tfinal])
y, t, x = impulse(sys[, t])
```

Calculate the response of the system `sys`

to an impulse at time `t = 0`

. For continous-time systems, the impulse is a unit Dirac impulse. For discrete-time systems, the impulse lasts one sample and has magnitude `1/Ts`

. If the final time `tfinal`

or time vector `t`

is not provided, one is calculated based on the system pole locations.

The return value is a structure of type `SimResult`

. A `SimResul`

can be plotted by `plot(result)`

, or destructured as `y, t, x = result`

.

`y`

has size `(ny, length(t), nu)`

, `x`

has size `(nx, length(t), nu)`

See also `lsim`

.

`ControlSystemsBase.lsim!`

— Method`res = lsim!(ws::LsimWorkspace, sys::AbstractStateSpace{<:Discrete}, u, [t]; x0)`

In-place version of `lsim`

that takes a workspace object created by calling `LsimWorkspace`

. *Notice*, if `u`

is a function, `res.u === ws.u`

. If `u`

is an array, `res.u === u`

.

`ControlSystemsBase.lsim`

— Method```
result = lsim(sys, u[, t]; x0, method])
result = lsim(sys, u::Function, t; x0, method)
```

Calculate the time response of system `sys`

to input `u`

. If `x0`

is omitted, a zero vector is used.

The result structure contains the fields `y, t, x, u`

and can be destructured automatically by iteration, e.g.,

`y, t, x, u = result`

`result::SimResult`

can also be plotted directly:

`plot(result, plotu=true, plotx=false)`

`y`

, `x`

, `u`

have time in the second dimension. Initial state `x0`

defaults to zero.

Continuous-time systems are simulated using an ODE solver if `u`

is a function (requires using ControlSystems). If `u`

is an array, the system is discretized (with `method=:zoh`

by default) before simulation. For a lower-level interface, see `?Simulator`

and `?solve`

. For continuous-time systems, keyword arguments are forwarded to the ODE solver. By default, the option `dtmax = t[2]-t[1]`

is used to prevent the solver from stepping over discontinuities in `u(x, t)`

. This prevents the solver from taking too large steps, but may also slow down the simulation when `u`

is smooth. To disable this behavior, set `dtmax = Inf`

.

`u`

can be a function or a *matrix* of precalculated control signals and must have dimensions `(nu, length(t))`

. If `u`

is a function, then `u(x,i)`

(for discrete systems) or `u(x,t)`

(for continuous ones) is called to calculate the control signal at every iteration (time instance used by solver). This can be used to provide a control law such as state feedback `u(x,t) = -L*x`

calculated by `lqr`

. To simulate a unit step at `t=t₀`

, use `(x,t)-> t ≥ t₀`

, for a ramp, use `(x,t)-> t`

, for a step at `t=5`

, use `(x,t)-> (t >= 5)`

etc.

*Note:* The function `u`

will be called once before simulating to verify that it returns an array of the correct dimensions. This can cause problems if `u`

is stateful. You can disable this check by passing `check_u = false`

.

For maximum performance, see function `lsim!`

, available for discrete-time systems only.

Usage example:

```
using ControlSystems
using LinearAlgebra: I
using Plots
A = [0 1; 0 0]
B = [0;1]
C = [1 0]
sys = ss(A,B,C,0)
Q = I
R = I
L = lqr(sys,Q,R)
u(x,t) = -L*x # Form control law
t = 0:0.1:5
x0 = [1,0]
y, t, x, uout = lsim(sys,u,t,x0=x0)
plot(t,x', lab=["Position" "Velocity"], xlabel="Time [s]")
# Alternative way of plotting
res = lsim(sys,u,t,x0=x0)
plot(res)
```

`ControlSystemsBase.stepinfo`

— Method`stepinfo(res::SimResult; y0 = nothing, yf = nothing, settling_th = 0.02, risetime_th = (0.1, 0.9))`

Compute the step response characteristics for a simulation result. The following information is computed and stored in a `StepInfo`

struct:

`y0`

: The initial value of the response`yf`

: The final value of the response`stepsize`

: The size of the step`peak`

: The peak value of the response`peaktime`

: The time at which the peak occurs`overshoot`

: The percentage overshoot of the response`undershoot`

: The percentage undershoot of the response. If the step response never reaches below the initial value, the undershoot is zero.`settlingtime`

: The time at which the response settles within`settling_th`

of the final value`settlingtimeind`

: The index at which the response settles within`settling_th`

of the final value`risetime`

: The time at which the response rises from`risetime_th[1]`

to`risetime_th[2]`

of the final value

**Arguments:**

`res`

: The result from a simulation using`step`

(or`lsim`

)`y0`

: The initial value, if not provided, the first value of the response is used.`yf`

: The final value, if not provided, the last value of the response is used. The simulation must have reached steady-state for an automatically computed value to make sense. If the simulation has not reached steady state, you may provide the final value manually.`settling_th`

: The threshold for computing the settling time. The settling time is the time at which the response settles within`settling_th`

of the final value.`risetime_th`

: The lower and upper threshold for computing the rise time. The rise time is the time at which the response rises from`risetime_th[1]`

to`risetime_th[2]`

of the final value.

**Example:**

```
G = tf([1], [1, 1, 1])
res = step(G, 15)
si = stepinfo(res)
plot(si)
```

`ControlSystemsBase.LsimWorkspace`

— Method```
LsimWorkspace(sys::AbstractStateSpace, N::Int)
LsimWorkspace(sys::AbstractStateSpace, u::AbstractMatrix)
LsimWorkspace{T}(ny, nu, nx, N)
```

Generate a workspace object for use with the in-place function `lsim!`

. `sys`

is the discrete-time system to be simulated and `N`

is the number of time steps, alternatively, the input `u`

can be provided instead of `N`

. Note: for threaded applications, create one workspace object per thread.

`ControlSystemsBase.StepInfo`

— Type`StepInfo`

Computed using `stepinfo`

**Fields:**

`y0`

: The initial value of the step response.`yf`

: The final value of the step response.`stepsize`

: The size of the step.`peak`

: The peak value of the step response.`peaktime`

: The time at which the peak occurs.`overshoot`

: The overshoot of the step response.`settlingtime`

: The time at which the step response has settled to within`settling_th`

of the final value.`settlingtimeind::Int`

: The index at which the step response has settled to within`settling_th`

of the final value.`risetime`

: The time at which the response rises from`risetime_th[1]`

to`risetime_th[2]`

of the final value`i10::Int`

: The index at which the response reaches`risetime_th[1]`

`i90::Int`

: The index at which the response reaches`risetime_th[2]`

`res::SimResult{SR}`

: The simulation result used to compute the step response characteristics.`settling_th`

: The threshold used to compute`settlingtime`

and`settlingtimeind`

.`risetime_th`

: The thresholds used to compute`risetime`

,`i10`

, and`i90`

.

`ControlSystemsBase.bode`

— Method`mag, phase, w = bode(sys[, w]; unwrap=true)`

Compute the magnitude and phase parts of the frequency response of system `sys`

at frequencies `w`

. See also `bodeplot`

`mag`

and `phase`

has size `(ny, nu, length(w))`

. If `unwrap`

is true (default), the function `unwrap!`

will be applied to the phase angles. This procedure is costly and can be avoided if the unwrapping is not required.

For higher performance, see the function `bodemag!`

that computes the magnitude only.

`ControlSystemsBase.bodemag!`

— Method`mag = bodemag!(ws::BodemagWorkspace, sys::LTISystem, w::AbstractVector)`

Compute the Bode magnitude operating in-place on an instance of `BodemagWorkspace`

. Note that the returned magnitude array is aliased with `ws.mag`

. The output array `mag`

is ∈ 𝐑(ny, nu, nω).

`ControlSystemsBase.evalfr`

— Method`evalfr(sys, x)`

Evaluate the transfer function of the LTI system sys at the complex number s=x (continuous-time) or z=x (discrete-time).

For many values of `x`

, use `freqresp`

instead.

`ControlSystemsBase.freqresp!`

— Method`freqresp!(R::Array{T, 3}, sys::LTISystem, w_vec::AbstractVector{<:Real})`

In-place version of `freqresp`

that takes a pre-allocated array `R`

of size (ny, nu, nw)`

`ControlSystemsBase.freqresp`

— Method`sys_fr = freqresp(sys, w)`

Evaluate the frequency response of a linear system

`w -> C*((iw*im*I - A)^-1)*B + D`

of system `sys`

over the frequency vector `w`

.

`ControlSystemsBase.nyquist`

— Method`re, im, w = nyquist(sys[, w])`

Compute the real and imaginary parts of the frequency response of system `sys`

at frequencies `w`

. See also `nyquistplot`

`re`

and `im`

has size `(ny, nu, length(w))`

`ControlSystemsBase.sigma`

— Method`sv, w = sigma(sys[, w])`

Compute the singular values `sv`

of the frequency response of system `sys`

at frequencies `w`

. See also `sigmaplot`

`sv`

has size `(min(ny, nu), length(w))`

`ControlSystemsBase.BodemagWorkspace`

— Method```
BodemagWorkspace(sys::LTISystem, N::Int)
BodemagWorkspace(sys::LTISystem, ω::AbstractVector)
BodemagWorkspace(R::Array{Complex{T}, 3}, mag::Array{T, 3})
BodemagWorkspace{T}(ny, nu, N)
```

Generate a workspace object for use with the in-place function `bodemag!`

. `N`

is the number of frequency points, alternatively, the input `ω`

can be provided instead of `N`

. Note: for threaded applications, create one workspace object per thread.

**Arguments:**

`mag`

: The output array ∈ 𝐑(ny, nu, nω)`R`

: Frequency-response array ∈ 𝐂(ny, nu, nω)

`ControlSystemsBase.TransferFunction`

— Method`F(s)`

, `F(omega, true)`

, `F(z, false)`

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(im
*Ts*omega) - F(z,false) evaluates the discrete-time transfer function F at z