Time Domain Simulation

Time domain simulation.

This file contains a collection of functions that calculate time responses for linear systems.

Convention for Time Series

This is a convention for function arguments and return values that represent time series: sequences of values that change over time. It is used throughout the library, for example in the functions forced_response(), step_response(), impulse_response(), and initial_response().

Note

This convention is different from the convention used in the library scipy.signal. In Scipy’s convention the meaning of rows and columns is interchanged. Thus, all 2D values must be transposed when they are used with functions from scipy.signal.

Types:

  • Arguments can be arrays, matrices, or nested lists.
  • Return values are arrays (not matrices).

The time vector is either 1D, or 2D with shape (1, n):

T = [[t1,     t2,     t3,     ..., tn    ]]

Input, state, and output all follow the same convention. Columns are different points in time, rows are different components. When there is only one row, a 1D object is accepted or returned, which adds convenience for SISO systems:

U = [[u1(t1), u1(t2), u1(t3), ..., u1(tn)]
     [u2(t1), u2(t2), u2(t3), ..., u2(tn)]
     ...
     ...
     [ui(t1), ui(t2), ui(t3), ..., ui(tn)]]

Same for X, Y

So, U[:,2] is the system’s input at the third point in time; and U[1] or U[1,:] is the sequence of values for the system’s second input.

The initial conditions are either 1D, or 2D with shape (j, 1):

X0 = [[x1]
      [x2]
      ...
      ...
      [xj]]

As all simulation functions return arrays, plotting is convenient:

t, y = step(sys)
plot(t, y)

The output of a MIMO system can be plotted like this:

t, y, x = lsim(sys, u, t)
plot(t, y[0], label='y_0')
plot(t, y[1], label='y_1')

The convention also works well with the state space form of linear systems. If D is the feedthrough matrix of a linear system, and U is its input (matrix or array), then the feedthrough part of the system’s response, can be computed like this:

ft = D * U

Time responses

control.forced_response(sys, T=None, U=0.0, X0=0.0, transpose=False, **keywords)

Simulate the output of a linear system.

As a convenience for parameters U, X0: Numbers (scalars) are converted to constant arrays with the correct shape. The correct shape is inferred from arguments sys and T.

For information on the shape of parameters U, T, X0 and return values T, yout, xout see: Convention for Time Series

Parameters:

sys: Lti (StateSpace, or TransferFunction) :

LTI system to simulate

T: array-like :

Time steps at which the input is defined; values must be evenly spaced.

U: array-like or number, optional :

Input array giving input at each time T (default = 0).

If U is None or 0, a special algorithm is used. This special algorithm is faster than the general algorithm, which is used otherwise.

X0: array-like or number, optional :

Initial condition (default = 0).

transpose: bool :

If True, transpose all input and output arrays (for backward compatibility with MATLAB and scipy.signal.lsim)

**keywords: :

Additional keyword arguments control the solution algorithm for the differential equations. These arguments are passed on to the function scipy.integrate.odeint(). See the documentation for scipy.integrate.odeint() for information about these arguments.

Returns:

T: array :

Time values of the output.

yout: array :

Response of the system.

xout: array :

Time evolution of the state vector.

See also

step_response, initial_response, impulse_response

Examples

>>> T, yout, xout = forced_response(sys, T, u, X0)
control.initial_response(sys, T=None, X0=0.0, input=0, output=None, transpose=False, **keywords)

Initial condition response of a linear system

If the system has multiple outputs (MIMO), optionally, one output may be selected. If no selection is made for the output, all outputs are given.

For information on the shape of parameters T, X0 and return values T, yout see: Convention for Time Series

Parameters:

sys: StateSpace, or TransferFunction :

LTI system to simulate

T: array-like object, optional :

Time vector (argument is autocomputed if not given)

X0: array-like object or number, optional :

Initial condition (default = 0)

Numbers are converted to constant arrays with the correct shape.

input: int :

Ignored, has no meaning in initial condition calculation. Parameter ensures compatibility with step_response and impulse_response

output: int :

Index of the output that will be used in this simulation. Set to None to not trim outputs

transpose: bool :

If True, transpose all input and output arrays (for backward compatibility with MATLAB and scipy.signal.lsim)

**keywords: :

Additional keyword arguments control the solution algorithm for the differential equations. These arguments are passed on to the function lsim(), which in turn passes them on to scipy.integrate.odeint(). See the documentation for scipy.integrate.odeint() for information about these arguments.

Returns:

T: array :

Time values of the output

yout: array :

Response of the system

See also

forced_response, impulse_response, step_response

Examples

>>> T, yout = initial_response(sys, T, X0)
control.step_response(sys, T=None, X0=0.0, input=None, output=None, transpose=False, **keywords)

Step response of a linear system

If the system has multiple inputs or outputs (MIMO), one input has to be selected for the simulation. Optionally, one output may be selected. The parameters input and output do this. All other inputs are set to 0, all other outputs are ignored.

For information on the shape of parameters T, X0 and return values T, yout see: Convention for Time Series

Parameters:

sys: StateSpace, or TransferFunction :

LTI system to simulate

T: array-like object, optional :

Time vector (argument is autocomputed if not given)

X0: array-like or number, optional :

Initial condition (default = 0)

Numbers are converted to constant arrays with the correct shape.

input: int :

Index of the input that will be used in this simulation.

output: int :

Index of the output that will be used in this simulation. Set to None to not trim outputs

transpose: bool :

If True, transpose all input and output arrays (for backward compatibility with MATLAB and scipy.signal.lsim)

**keywords: :

Additional keyword arguments control the solution algorithm for the differential equations. These arguments are passed on to the function lsim(), which in turn passes them on to scipy.integrate.odeint(). See the documentation for scipy.integrate.odeint() for information about these arguments.

Returns:

T: array :

Time values of the output

yout: array :

Response of the system

See also

forced_response, initial_response, impulse_response

Examples

>>> T, yout = step_response(sys, T, X0)

Phase portraits

control.phaseplot.phase_plot(odefun, X=None, Y=None, scale=1, X0=None, T=None, lingrid=None, lintime=None, logtime=None, timepts=None, parms=(), verbose=True)

Phase plot for 2D dynamical systems

Produces a vector field or stream line plot for a planar system.

Call signatures:
phase_plot(func, X, Y, ...) - display vector field on meshgrid phase_plot(func, X, Y, scale, ...) - scale arrows phase_plot(func. X0=(...), T=Tmax, ...) - display stream lines phase_plot(func, X, Y, X0=[...], T=Tmax, ...) - plot both phase_plot(func, X0=[...], T=Tmax, lingrid=N, ...) - plot both phase_plot(func, X0=[...], lintime=N, ...) - stream lines with arrows
Parameters:

func : callable(x, t, ...)

Computes the time derivative of y (compatible with odeint). The function should be the same for as used for scipy.integrate. Namely, it should be a function of the form dxdt = F(x, t) that accepts a state x of dimension 2 and returns a derivative dx/dt of dimension 2.

X, Y: ndarray, optional :

Two 1-D arrays representing x and y coordinates of a grid. These arguments are passed to meshgrid and generate the lists of points at which the vector field is plotted. If absent (or None), the vector field is not plotted.

scale: float, optional :

Scale size of arrows; default = 1

X0: ndarray of initial conditions, optional :

List of initial conditions from which streamlines are plotted. Each initial condition should be a pair of numbers.

T: array-like or number, optional :

Length of time to run simulations that generate streamlines. If a single number, the same simulation time is used for all initial conditions. Otherwise, should be a list of length len(X0) that gives the simulation time for each initial condition. Default value = 50.

lingrid = N or (N, M): integer or 2-tuple of integers, optional :

If X0 is given and X, Y are missing, a grid of arrows is produced using the limits of the initial conditions, with N grid points in each dimension or N grid points in x and M grid points in y.

lintime = N: integer, optional :

Draw N arrows using equally space time points

logtime = (N, lambda): (integer, float), optional :

Draw N arrows using exponential time constant lambda

timepts = [t1, t2, ...]: array-like, optional :

Draw arrows at the given list times

parms: tuple, optional :

List of parameters to pass to vector field: func(x, t, *parms)

See also

box_grid, Y

control.phaseplot.box_grid(xlimp, ylimp)

box_grid generate list of points on edge of box

list = box_grid([xmin xmax xnum], [ymin ymax ynum]) generates a list of points that correspond to a uniform grid at the end of the box defined by the corners [xmin ymin] and [xmax ymax].