Welcome to Skogestad Python’s documentation!¶
Indices and tables¶
Created on Jan 27, 2012
@author: Carl Sandrock
-
utils.
BoundKS
(G, poles, up, e=1e-05)¶ The functions uses equaption 6.24 (p229) to calculate the peak value for KS transfer function using the stable version of the plant.
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
- poles : numpy array (number of poles)
List of right-half plane poles.
- up : numpy array (number of poles)
List of input pole directions.
- e : float
Avoid division by zero. Let epsilon be very small (optional).
Returns: - KS_max : float
Minimum peak value.
-
utils.
BoundST
(G, poles, zeros, deadtime=None)¶ This function will calculate the minimum peak values of S and T if the system has zeros and poles in the input or output. For standard conditions Equation 6.8 (p224) is applied. Equation 6.16 (p226) is used when deadtime is included.
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
- poles : numpy array (number of zeros)
List of poles.
- zeros : numpy array (number of zeros)
List of zeros.
- deadtime : numpy matrix (n, n)
Deadtime or time delay for G.
Returns: - Ms_min : real
Minimum peak value.
-
utils.
Closed_loop
(Kz, Kp, Gz, Gp)¶ Return zero and pole polynomial for a closed loop function.
Parameters: - Kz & Gz : list
Polynomial constants in the numerator.
- Kz & Gz : list
Polynomial constants in the denominator.
Returns: - Zeros_poly : list
List of zero polynomial for closed loop function.
- Poles_poly : list
List of pole polynomial for closed loop function.
-
utils.
ControllerTuning
(G, method='ZN')¶ Calculates either the Ziegler-Nichols or Tyreus-Luyben tuning parameters for a PI controller based on the continuous cycling method.
Parameters: - G : tf
plant model
- method : Use ‘ZN’ for Ziegler-Nichols tuning parameters and
‘TT’ for Tyreus-Luyben parameters. The default is to return Ziegler-Nichols tuning parameters.
Returns: - Kc : array containing a real number
proportional gain
- Taui : array containing a real number
integral gain
- Ku : array containing a real number
ultimate P controller gain
- Pu : array containing a real number
corresponding period of oscillations
-
utils.
IterRGA
(A, n)¶ Computes the n’th iteration of the RGA.
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
Returns: - n’th iteration of RGA matrix : matrix
iterated RGA matrix of complex numbers.
-
utils.
RGA
(G)¶ Computes the RGA (Relative Gain Array) of a matrix.
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
Returns: - RGA matrix : matrix
RGA matrix of complex numbers.
-
utils.
RGAnumber
(G, I)¶ Computes the RGA (Relative Gain Array) number of a matrix.
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
- I : numpy matrix
Pairing matrix.
Returns: - RGA number : float
RGA number.
-
utils.
SVD
(G)¶ Returns the singular values (Sv) as well as the input and output singular vectors (V and U respectively).
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
Returns: - U : matrix of complex numbers
Unitary matrix of output singular vectors.
- Sv : array
Singular values of Gin arranged in decending order.
- V : matrix of complex numbers
Unitary matrix of input singular vectors.
-
utils.
Wp
(wB, M, A, s)¶ Computes the magnitude of the performance weighting function. Based on Equation 2.105 (p62).
Parameters: - wB : float
Approximate bandwidth requirement. Asymptote crosses 1 at this frequency.
- M : float
Maximum frequency.
- A : float
Maximum steady state tracking error. Typically 0.
- s : complex
Typically w*1j.
Returns: - `|Wp(s)|` : float
The magnitude of the performance weighting fucntion at a specific frequency (s).
-
utils.
arrayfun
(f, A)¶ Recurses down to scalar elements in A, then applies f, returning lists containing the result.
Parameters: - A : array
- f : function
Returns: - arrayfun : list
- >>> def f(x):
- … return 1.
- >>> arrayfun(f, numpy.array([1, 2, 3]))
- [1.0, 1.0, 1.0]
- >>> arrayfun(f, numpy.array([[1, 2, 3], [1, 2, 3]]))
- [[1.0, 1.0, 1.0], [1.0, 1.0, 1.0]]
- >>> arrayfun(f, 1)
- 1.0
-
utils.
astf
(maybetf)¶ Parameters: maybetf – something which could be a tf Returns: a transfer function object >>> G = tf( 1, [ 1, 1]) >>> astf(G) tf([1.], [1. 1.])
>>> astf(1) tf([1.], [1.])
>>> astf(numpy.matrix([[G, 1.], [0., G]])) matrix([[tf([1.], [1. 1.]), tf([1.], [1.])], [tf([0.], [1]), tf([1.], [1. 1.])]], dtype=object)
-
utils.
circle
(cx, cy, r)¶ Return the coordinates of a circle
Parameters: - cx : float
Center x coordinate.
- cy : float
Center y coordinate.
- r : float
Radius.
Returns: - x, y : float
Circle coordinates.
-
utils.
det
(A)¶ Calculate determinant via elementary operations
Parameters: A – Array-like object Returns: determinant >>> det(2.) 2.0
>>> A = [[1., 2.], ... [1., 2.]] >>> det(A) 0.0
>>> B = [[1., 2.], ... [3., 4.]] >>> det(B) -2.0
>>> C = [[1., 2., 3.], ... [1., 3., 2.], ... [3., 2., 1.]] >>> det(C) -12.0
# Can handle matrices of tf objects # TODO: make this a little more natural (without the .matrix) >>> G11 = tf([1], [1, 2]) >>> G = mimotf([[G11, G11], [G11, G11]]) >>> det(G.matrix) tf([0.], [1])
>>> G = mimotf([[G11, 2*G11], [G11**2, 3*G11]]) >>> det(G.matrix) tf([ 3. 16. 28. 16.], [ 1. 10. 40. 80. 80. 32.])
-
utils.
distRHPZ
(G, Gd, RHP_Z)¶ Applies equation 6.48 (p239) For performance requirements imposed by disturbances. Calculate the system’s zeros alignment with the disturbacne matrix.
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
- gd : numpy matrix (n x 1)
The transfer function Gd(s) of the distrurbances.
- RHP_Z : complex
Right-half plane zero
Returns: - Dist_RHPZ : float
Minimum peak value.
-
utils.
distRej
(G, gd)¶ Convenience wrapper for calculation of ||gd||2 (equation 6.42, p238) and the disturbance condition number (equation 6.43) for each disturbance.
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
- gd : numpy matrix (m x n)
The transfer function Gd(s) of the distrurbances.
Returns: - 1/||gd|| :math:`_2` : float
The inverse of the 2-norm of a single disturbance gd.
- distCondNum : float
The disturbance condition number \(\sigma\) (G) \(\sigma\) (G \(^{-1}\) yd)
- yd : numpy matrix
Disturbance direction.
-
utils.
feedback
(forward, backward=None, positive=False)¶ Defined for use in connect function Calculates a feedback loop This version is for trasnfer function objects Negative feedback is assumed, use positive=True for positive feedback Forward refers to the function that goes out of the comparator Backward refers to the function that goes into the comparator
-
utils.
feedback_mimo
(forward, backward=None, positive=False)¶ Calculates a feedback loop This version is for matrices Negative feedback is assumed, use positive=True for positive feedback Forward refers to the function that goes out of the comparator Backward refers to the function that goes into the comparator
-
utils.
findst
(G, K)¶ Find S and T given a value for G and K.
Parameters: - G : numpy array
Matrix of transfer functions.
- K : numpy array
Matrix of controller functions.
Returns: - S : numpy array
Matrix of sensitivities.
- T : numpy array
Matrix of complementary sensitivities.
-
utils.
freq
(G)¶ Calculate the frequency response for an optimisation problem
Parameters: - G : tf
plant model
Returns: - Gw : frequency response function
-
utils.
gaintf
(K)¶ Transform a gain value into a transfer function.
Parameters: - K : float
Gain.
Returns: - gaintf : tf
Transfer function.
-
utils.
kalman_controllable
(A, B, C, P=None, RP=None)¶ Computes the Kalman Controllable Canonical Form of the inout system A, B, C, making use of QR Decomposition. Can be used in sequentially with kalman_observable to obtain a minimal realisation.
Parameters: - A : numpy matrix
The system state matrix.
- B : numpy matrix
The system input matrix.
- C : numpy matrix
The system output matrix.
- P : (optional) numpy matrix
The controllability matrix
- RP : (optional int)
Returns: - Ac : numpy matrix
The state matrix of the controllable system
- Bc : nump matrix
The input matrix of the controllable system
- Cc : numpy matrix
The output matrix of the controllable system
-
utils.
kalman_observable
(A, B, C, Q=None, RQ=None)¶ Computes the Kalman Observable Canonical Form of the inout system A, B, C, making use of QR Decomposition. Can be used in sequentially with kalman_controllable to obtain a minimal realisation.
Parameters: - A : numpy matrix
The system state matrix.
- B : numpy matrix
The system input matrix.
- C : numpy matrix
The system output matrix.
- Q : (optional) numpy matrix
Observability matrix
- RQ : (optional) int
Rank of observability matrxi
Returns: - Ao : numpy matrix
The state matrix of the observable system
- Bo : nump matrix
The input matrix of the observable system
- Co : numpy matrix
The output matrix of the observable system
-
utils.
lcm_of_all_minors
(G)¶ Returns the lowest common multiple of all minors of G
-
utils.
listify
(A)¶ Transform a gain value into a transfer function.
Parameters: - K : float
Gain.
Returns: - gaintf : tf
Transfer function.
-
utils.
margins
(G)¶ Calculates the gain and phase margins, together with the gain and phase crossover frequency for a plant model
Parameters: - G : tf
plant model
Returns: - GM : array containing a real number
gain margin
- PM : array containing a real number
phase margin
- wc : array containing a real number
gain crossover frequency where |G(jwc)| = 1
- w_180 : array containing a real number
phase crossover frequency where angle[G(jw_180] = -180 deg
-
utils.
marginsclosedloop
(L)¶ Calculates the gain and phase margins, together with the gain and phase crossover frequency for a control model
Parameters: - L : tf
loop transfer function
Returns: - GM : real
gain margin
- PM : real
phase margin
- wc : real
gain crossover frequency for L
- wb : real
closed loop bandwidth for S
- wbt : real
closed loop bandwidth for T
-
utils.
matrix_as_scalar
(M)¶ Return a scalar from a 1x1 matrix
Parameters: M – matrix Returns: scalar part of matrix if it is 1x1 else just a matrix
-
utils.
maxpeak
(G, w_start=-2, w_end=2, points=1000)¶ Computes the maximum bode magnitude peak of a transfer function
-
class
utils.
mimotf
(matrix)¶ Represents MIMO transfer function matrix
This is a pretty basic wrapper around the numpy.matrix class which deals with most of the heavy lifting.
You can construct the object from siso tf objects similarly to calling numpy.matrix:
>>> G11 = G12 = G21 = G22 = tf(1, [1, 1]) >>> G = mimotf([[G11, G12], [G21, G22]]) >>> G mimotf([[tf([1.], [1. 1.]) tf([1.], [1. 1.])] [tf([1.], [1. 1.]) tf([1.], [1. 1.])]])
Some coersion will take place on the elements: >>> mimotf([[1]]) mimotf([[tf([1.], [1.])]])
The object knows how to do:
addition
>>> G + G mimotf([[tf([2.], [1. 1.]) tf([2.], [1. 1.])] [tf([2.], [1. 1.]) tf([2.], [1. 1.])]])
>>> 0 + G mimotf([[tf([1.], [1. 1.]) tf([1.], [1. 1.])] [tf([1.], [1. 1.]) tf([1.], [1. 1.])]])
>>> G + 0 mimotf([[tf([1.], [1. 1.]) tf([1.], [1. 1.])] [tf([1.], [1. 1.]) tf([1.], [1. 1.])]])
multiplication >>> G * G mimotf([[tf([2.], [1. 2. 1.]) tf([2.], [1. 2. 1.])]
[tf([2.], [1. 2. 1.]) tf([2.], [1. 2. 1.])]])>>> 1*G mimotf([[tf([1.], [1. 1.]) tf([1.], [1. 1.])] [tf([1.], [1. 1.]) tf([1.], [1. 1.])]])
>>> G*1 mimotf([[tf([1.], [1. 1.]) tf([1.], [1. 1.])] [tf([1.], [1. 1.]) tf([1.], [1. 1.])]])
>>> G*tf(1) mimotf([[tf([1.], [1. 1.]) tf([1.], [1. 1.])] [tf([1.], [1. 1.]) tf([1.], [1. 1.])]])
>>> tf(1)*G mimotf([[tf([1.], [1. 1.]) tf([1.], [1. 1.])] [tf([1.], [1. 1.]) tf([1.], [1. 1.])]])
exponentiation with positive integer constants
>>> G**2 mimotf([[tf([2.], [1. 2. 1.]) tf([2.], [1. 2. 1.])] [tf([2.], [1. 2. 1.]) tf([2.], [1. 2. 1.])]])
Methods
__call__
(s)>>> G = mimotf([[1]]) >>> G(0) matrix([[1.]])
inverse
()Calculate inverse of mimotf object poles
()Calculate poles >>> s = tf([1, 0], [1]) >>> G = mimotf([[(s - 1) / (s + 2), 4 / (s + 2)], … cofactor_mat det mimotf_slice zeros -
inverse
()¶ Calculate inverse of mimotf object
>>> s = tf([1, 0], 1) >>> G = mimotf([[(s - 1) / (s + 2), 4 / (s + 2)], ... [4.5 / (s + 2), 2 * (s - 1) / (s + 2)]]) >>> G.inverse() matrix([[tf([ 1. -1.], [ 1. -4.]), tf([-2.], [ 1. -4.])], [tf([-2.25], [ 1. -4.]), tf([ 0.5 -0.5], [ 1. -4.])]], dtype=object)
>>> G.inverse()*G.matrix matrix([[tf([1.], [1.]), tf([0.], [1])], [tf([0.], [1]), tf([1.], [1.])]], dtype=object)
-
poles
()¶ Calculate poles >>> s = tf([1, 0], [1]) >>> G = mimotf([[(s - 1) / (s + 2), 4 / (s + 2)], … [4.5 / (s + 2), 2 * (s - 1) / (s + 2)]]) >>> G.poles() array([-2.])
-
-
utils.
mimotf2sym
(G, deadtime=False)¶ Converts a mimotf object making use of individual tf objects to a transfer function system in sympy.Matrix form.
Parameters: - G : mimotf matrix
The mimotf system matrix.
- deadtime: boolean
Should deadtime be added to sympy matrix or not.
Returns: - Gs : sympy matrix
The sympy system matrix
- s : sympy symbol
Sympy symbol generated
-
utils.
minimal_realisation
(a, b, c)¶ “This function will obtain a minimal realisation for a state space model in the form given in Skogestad second edition p 119 equations 4.3 and 4.4
Parameters: - a – numpy matrix the A matrix in the state space model
- b – numpy matrix the B matrix in the state space model
- c – numpy matrix the C matrix in the state space model
Examples
Example 1:
>>> A = numpy.matrix([[0, 0, 0, 0], ... [0, -2, 0, 0], ... [2.5, 2.5, -1, 0], ... [2.5, 2.5, 0, -3]])
>>> B = numpy.matrix([[1], ... [1], ... [0], ... [0]])
>>> C = numpy.matrix([0, 0, 1, 1])
>>> Aco, Bco, Cco = minimal_realisation(A, B, C)
Add null to eliminate negatives null elements (-0.)
>>> Aco.round(decimals=3) + 0. array([[-2.038, 5.192], [ 0.377, -0.962]])
>>> Bco.round(decimals=3) + 0. array([[ 0. ], [-1.388]])
>>> Cco.round(decimals=3) + 0. array([[-1.388, 0. ]])
Example 2:
>>> A = numpy.matrix([[1, 1, 0], ... [0, 1, 0], ... [0, 1, 1]])
>>> B = numpy.matrix([[0, 1], ... [1, 0], ... [0, 1]])
>>> C = numpy.matrix([1, 1, 1])
>>> Aco, Bco, Cco = minimal_realisation(A, B, C)
Add null to eliminate negatives null elements (-0.)
>>> Aco.round(decimals=3) + 0. array([[ 1. , 0. ], [-1.414, 1. ]])
>>> Bco.round(decimals=3) + 0. array([[-1. , 0. ], [ 0. , 1.414]])
>>> Cco.round(decimals=3) + 0. array([[-1. , 1.414]])
-
utils.
minors
(G, order)¶ Returns the order minors of a MIMO tf G.
-
utils.
omega
(w_start, w_end)¶ Convenience wrapper Defines the frequency range for calculation of frequency response Frequency in rad/time where time is the time unit used in the model.
-
utils.
phase
(G, deg=False)¶ Return the phase angle in degrees or radians
Parameters: - G : tf
Plant of transfer functions.
- deg : booleans
True if radians result is required, otherwise degree is default (optional).
Returns: - phase : float
Phase angle.
-
utils.
pole_zero_directions
(G, vec, dir_type, display_type='a', e=1e-08, min_tol=1e-05, max_tol=10000.0)¶ Crude method to calculate the input and output direction of a pole or zero, from the SVD.
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
- vec : array
A vector containing all the transmission poles or zeros of a system.
- dir_type : string
Type of direction to calculate.
dir_type Choose ‘p’ Poles ‘z’ Zeros - display_type : string
Choose the type of directional data to return (optional).
display_type Directional data to return ‘a’ All data (default) ‘u’ Only input direction ‘y’ Only output direction - e : float
Avoid division by zero. Let epsilon be very small (optional).
- min_tol : float
Acceptable tolerance for zero validation. Let min_tol be very small (optional).
- max_tol : float
Acceptable tolerance for pole validation. Let max_tol be very small (optional).
Returns: - pz_dir : array
Pole or zero direction in the form: (pole/zero, input direction, output direction, valid)
- valid : integer array
If 1 the directions are valid, else if 0 the directions are not valid.
-
utils.
poles
(G=None, A=None)¶ If G is passed then return the poles of a multivariable transfer function system. Applies Theorem 4.4 (p135). If G is NOT specified but A is, returns the poles from the state space description as per section 4.4.2.
Parameters: - G : sympy or mimotf matrix (n x n)
The transfer function G(s) of the system.
- A : State Space A matrix
Returns: - pole : array
List of poles.
-
utils.
poles_and_zeros_of_square_tf_matrix
(G)¶ Determine poles and zeros of a square mimotf matrix, making use of the determinant. This method may fail in special cases. If terms cancel out during calculation of the determinant, not all poles and zeros will be determined.
Parameters: - G : mimotf matrix (n x n)
The transfer function of the system.
Returns: - z : array
List of zeros.
- p : array
List of poles.
- possible_cancel : boolean
Test whether terms were possibly cancelled out in determinant calculation.
-
utils.
polygcd
(a, b)¶ Find the approximate Greatest Common Divisor of two polynomials
>>> a = numpy.poly1d([1, 1]) * numpy.poly1d([1, 2]) >>> b = numpy.poly1d([1, 1]) * numpy.poly1d([1, 3]) >>> polygcd(a, b) poly1d([1., 1.])
>>> polygcd(numpy.poly1d([1, 1]), numpy.poly1d([1])) poly1d([1.])
-
utils.
scaling
(G_hat, e, u, input_type='symbolic', Gd_hat=None, d=None)¶ Receives symbolic matrix of plant and disturbance transfer functions as well as array of maximum deviations, scales plant variables according to eq () and ()
Parameters: - G_hat : matrix of plant WITHOUT deadtime
- e : array of maximum plant output variable deviations
in same order as G matrix plant outputs
- u : array of maximum plant input variable deviations
in same order as G matrix plant inputs
- input_type : specifies whether input is symbolic matrix or utils mimotf
- Gd_hat : optional
matrix of plant disturbance model WITHOUT deadtime
- d : optional
array of maximum plant disturbance variable deviations in same order as Gd matrix plant disturbances
- Returns
- ———-
- G_scaled : scaled plant function
- Gd_scaled : scaled plant disturbance function
-
utils.
sigmas
(A, position=None)¶ Returns the singular values of A
Parameters: - A : array
Transfer function matrix.
- position : string
Type of sigmas to return (optional).
position Type of sigmas to return max Maximum singular value min Minimal singular value
Returns: - :math:`sigma` (A) : array
Singular values of A arranged in decending order.
- This is a convenience wrapper to enable easy calculation of
- singular values over frequency
-
utils.
ssr_solve
(A, B, C, D)¶ Solves the zeros and poles of a state-space representation of a system.
Parameters: - A – System state matrix
- B – matrix
- C – matrix
- D – matrix
For information on the meanings of A, B, C, and D consult Skogestad 4.1.1
- Returns:
- zeros: The system’s zeros poles: the system’s poles
TODO: Add any other relevant values to solve for, for example, if coprime factorisations are useful somewhere add them to this function’s return dict rather than writing another function.
-
utils.
state_controllability
(A, B)¶ This method checks if the state space description of the system is state controllable according to Definition 4.1 (p127).
Parameters: - A : numpy matrix
Matrix A of state-space representation.
- B : numpy matrix
Matrix B of state-space representation.
Returns: - state_control : boolean
True if state controllable
- u_p : array
Input pole vectors for the states u_p_i
- control_matrix : numpy matrix
State Controllability Matrix
-
utils.
state_observability_matrix
(a, c)¶ calculate the observability matrix
Parameters: - a – numpy matrix the A matrix in the state space model
- c – numpy matrix the C matrix in the state space model
-
utils.
sv_dir
(G, table=False)¶ Returns the input and output singular vectors associated with the minimum and maximum singular values.
Parameters: - G : numpy matrix (n x n)
The transfer function G(s) of the system.
- table : True of False boolean
Default set to False.
Returns: - u : list of arrays containing complex numbers
Output vector associated with the maximum and minium singular values. The maximum singular output vector is the first entry u[0] and the minimum is the second u[1].
- v : list of arrays containing complex numbers
Input vector associated with the maximum and minium singular values. The maximum singular intput vector is the first entry u[0] and the minimum is the second u[1].
- table : If table is True then the output and input vectors are summarised
and returned as a table in the command window. Values are reported to five significant figures.
-
utils.
sym2mimotf
(Gmat, deadtime=None)¶ Converts a MIMO transfer function system in sympy.Matrix form to a mimotf object making use of individual tf objects.
Parameters: - Gmat : sympy matrix
The system transfer function matrix.
- deadtime: numpy matrix of same dimensions as Gmat
The dead times of Gmat with corresponding indexes.
Returns: - Gmimotf : sympy matrix
The mimotf system matrix
-
class
utils.
tf
(numerator, denominator=1, deadtime=0, name='', u='', y='', prec=3)¶ Very basic transfer function object
Construct with a numerator and denominator:
>>> G = tf(1, [1, 1]) >>> G tf([1.], [1. 1.])
>>> G2 = tf(1, [2, 1])
The object knows how to do:
addition
>>> G + G2 tf([1.5 1. ], [1. 1.5 0.5]) >>> G + G # check for simplification tf([2.], [1. 1.])
multiplication
>>> G * G2 tf([0.5], [1. 1.5 0.5])
division
>>> G / G2 tf([2. 1.], [1. 1.])
Deadtime is supported:
>>> G3 = tf(1, [1, 1], deadtime=2) >>> G3 tf([1.], [1. 1.], deadtime=2)
Note we can’t add transfer functions with different deadtime:
>>> G2 + G3 Traceback (most recent call last): ... ValueError: Transfer functions can only be added if their deadtimes are the same. self=tf([0.5], [1. 0.5]), other=tf([1.], [1. 1.], deadtime=2)
Although we can add a zero-gain tf to anything
>>> G2 + 0*G3 tf([0.5], [1. 0.5])
>>> 0*G2 + G3 tf([1.], [1. 1.], deadtime=2)
It is sometimes useful to define
>>> s = tf([1, 0]) >>> 1 + s tf([1. 1.], [1.])
>>> 1/(s + 1) tf([1.], [1. 1.])
Methods
__call__
(s)This allows the transfer function to be evaluated at particular values of s. exp
()If this is basically “D*s” defined as tf([D, 0], 1), return dead time inverse
()Inverse of the transfer function lsim
(*args)Negative step response step
(*args)Step response poles simplify zeros -
exp
()¶ - If this is basically “D*s” defined as tf([D, 0], 1),
- return dead time
>>> s = tf([1, 0], 1) >>> numpy.exp(-2*s) tf([1.], [1.], deadtime=2.0)
-
inverse
()¶ Inverse of the transfer function
-
lsim
(*args)¶ Negative step response
-
step
(*args)¶ Step response
-
-
utils.
tf2ss
(H)¶ Converts a mimotf object to the controllable canonical form state space representation. This method and the examples were obtained from course work notes available at http://www.egr.msu.edu/classes/me851/jchoi/lecture/Lect_20.pdf which appears to derive the method from “A Linear Systems Primer” by Antsaklis and Birkhauser.
Parameters: - H : mimotf
The mimotf object transfer function form
Returns: - Ac : numpy matrix
The state matrix of the observable system
- Bc : nump matrix
The input matrix of the observable system
- Cc : numpy matrix
The output matrix of the observable system
- Dc : numpy matrix
The output matrix of the observable system
-
utils.
tf_step
(G, t_end=10, initial_val=0, points=1000, constraint=None, Y=None, method='numeric')¶ Validate the step response data of a transfer function by considering dead time and constraints. A unit step response is generated.
Parameters: - G : tf
Transfer function (input[u] or output[y]) to evauate step response.
- Y : tf
Transfer function output[y] to evaluate constrain step response (optional) (required if constraint is specified).
- t_end : integer
length of time to evaluate step response (optional).
- initial_val : integer
starting value to evalaute step response (optional).
- points : integer
number of iteration that will be calculated (optional).
- constraint : real
The upper limit the step response cannot exceed. Is only calculated if a value is specified (optional).
- method : [‘numeric’,’analytic’]
The method that is used to calculate a constrainted response. A constraint value is required (optional).
Returns: - timedata : array
Array of floating time values.
- process : array (1 or 2 dim)
1 or 2 dimensional array of floating process values.
-
utils.
zero_directions_ss
(A, B, C, D)¶ This function calculates the zeros with input and output directions from a state space representation using the method outlined on pg. 140
Parameters: - A : numpy matrix
A matrix of state space representation
- B : numpy matrix
B matrix of state space representation
- C : numpy matrix
C matrix of state space representation
- D : numpy matrix
D matrix of state space representation
Returns: - zeros_in_out : list
zeros_in_out[i] contains a zero, input direction vector and output direction vector
-
utils.
zeros
(G=None, A=None, B=None, C=None, D=None)¶ Return the zeros of a multivariable transfer function system for with transfer functions or state-space. For transfer functions, Theorem 4.5 (p139) is used. For state-space, the method from Equations 4.66 and 4.67 (p138) is applied.
Parameters: - G : sympy or mimotf matrix (n x n)
The transfer function G(s) of the system.
- A, B, C, D : numpy matrix
State space parameters
Returns: - zero : array
List of zeros.
Common features to plotting functions in this script
Default parameters¶
- axlim : list [xmin, xmax, ymin, ymax]
- A list containing the minimum and maximum limits for the x and y-axis. To autoscale a limit enter ‘None’ in its placeholder. The default is to allow autoscaling of the axes.
- w_start : float
- The x-axis valve at which to start the plot.
- w_end : float
- The x-axis value at which to stop plotting.
- points : float
- The number of data points to be used in generating the plot.
Example¶
- def G(s):
- return numpy.matrix([[s/(s+1), 1],
- [s**2 + 1, 1/(s+1)]])
plt.figure(‘Example 1’) your_utilsplot_functionA(G, w_start=-5, w_end=2, axlim=[None, None, 0, 1], more_parameters) plt.show()
plt.figure(‘Example 2’) plt.subplot(2, 1, 1) your_utilsplot_functionB(G) plt.subplot(2, 1, 2) your_utilsplot_functionC(G) plt.show()
-
utilsplot.
adjust_spine
(xlabel, ylabel, x0=0, y0=0, width=1, height=1)¶ General function to adjust the margins for subplots.
Parameters: - xlabel : string
Label on the main x-axis.
- ylabel : string
Label on the main x-axis.
- x0 : integer
Horizontal offset of xlabel.
- y0 : integer
Verticle offset of ylabel.
- width : float
Scaling factor of width of subplots.
- height : float
Scaling factor of height of subplots.
Returns: - fig : matplotlib subplot area
-
utilsplot.
bode
(G, w_start=-2, w_end=2, axlim=None, points=1000, margin=False)¶ Shows the bode plot for a plant model
Parameters: - G : tf
Plant transfer function.
- margin : boolean
Show the cross over frequencies on the plot (optional).
Returns: - GM : array containing a real number
Gain margin.
- PM : array containing a real number
Phase margin.
- Plot : matplotlib figure
-
utilsplot.
bodeclosedloop
(G, K, w_start=-2, w_end=2, axlim=None, points=1000, margin=False)¶ Shows the bode plot for a controller model
Parameters: - G : tf
Plant transfer function.
- K : tf
Controller transfer function.
- margin : boolean
Show the cross over frequencies on the plot (optional).
-
utilsplot.
complexplane
(args, color=True, marker='o', msize=5)¶ Plot up to 8 arguments on a complex plane (limited by the colors) Useful when you wish to compare sets of complex numbers graphically or plot your poles and zeros
Parameters: - args : A list of the list of numbers to plot
- color : True if every tuple of info must be a different color
False if all must be the same color
- marker : Type of marker to use
- msize : Size of the marker
- Example:
A = [1+2j, 1-2j, 1+1j, 2-1j] B = [1+2j, 3+2j, 1, 1+2j] complexplane([A, B, [1+3j, 2+5j]], marker=’+’, msize=8)
-
utilsplot.
condtn_nm_plot
(G, w_start=-2, w_end=2, axlim=None, points=1000)¶ Plot of the condition number, the maximum over the minimum singular value
Parameters: - G : numpy matrix
Plant model.
Returns: - Plot : matplotlib figure
-
utilsplot.
dis_rejctn_plot
(G, Gd, S=None, w_start=-2, w_end=2, axlim=None, points=1000)¶ A subplot of disturbance condition number to check for input saturation (equation 6.43, p238). Two more subplots indicate if the disturbances fall withing the bounds of S, applying equations 6.45 and 6.46 (p239).
Parameters: - G : numpy matrix
Plant model.
- Gd : numpy matrix
Plant disturbance model.
- S : numpy matrix
Sensitivity function (optional, if available).
- # TODO test S condition
- Returns
- ——-
- Plot : matplotlib figure
-
utilsplot.
freq_step_response_plot
(G, K, Kc, t_end=50, freqtype='S', w_start=-2, w_end=2, axlim=None, points=1000)¶ A subplot function for both the frequency response and step response for a controlled plant
Parameters: - G : tf
Plant transfer function.
- K : tf
Controller transfer function.
- Kc : integer
Controller constant.
- t_end : integer
Time period which the step response should occur.
- freqtype : string (optional)
Type of function to plot:
freqtype Type of function to plot S Sensitivity function T Complementary sensitivity function L Loop function
Returns: - Plot : matplotlib figure
-
utilsplot.
input_acceptable_const_plot
(G, Gd, w_start=-2, w_end=2, axlim=None, points=1000, modified=False)¶ Subbplots for input constraints for accepatable control. Applies equation 6.55 (p241).
Parameters: - G : numpy matrix
Plant model.
- Gd : numpy matrix
Plant disturbance model.
- modified : boolean
If true, the arguments in the equation are change to \(\sigma_1 (G) + 1 \geq |u_i^H g_d|\). This is to avoid a negative log scale.
Returns: - Plot : matplotlib figure
-
utilsplot.
input_perfect_const_plot
(G, Gd, w_start=-2, w_end=2, axlim=None, points=1000, simultaneous=False)¶ Plot for input constraints for perfect control. Applies equation 6.50 (p240).
Parameters: - G : numpy matrix
Plant model.
- Gd : numpy matrix
Plant disturbance model.
- simultaneous : boolean.
If true, the induced max-norm is calculated for simultaneous disturbances (optional).
Returns: - Plot : matplotlib figure
-
utilsplot.
mimo_bode
(G, w_start=-2, w_end=2, axlim=None, points=1000, Kin=None, text=False, sv_all=False)¶ Plots the max and min singular values of G and computes the crossover frequency.
If a controller is specified, the max and min singular values of S are also plotted and the bandwidth frequency computed (p81).
Parameters: - G : numpy matrix
Matrix of plant transfer functions.
- Kin : numpy matrix
Controller matrix (optional).
- text : boolean
If true, the crossover and bandwidth frequencies are plotted (optional).
- sv_all : boolean
If true, plot all the singular values of the plant (optional).
Returns: - wC : real
Crossover frequency.
- wB : real
Bandwidth frequency.
- Plot : matplotlib figure
-
utilsplot.
mimo_nyquist_plot
(L, w_start=-2, w_end=2, axlim=None, points=1000)¶ Nyquist stability plot for MIMO system.
Parameters: - L : numpy matrix
Closed loop transfer function matrix as a function of s, i.e. def L(s).
Returns: - Plot : matplotlib figure
-
utilsplot.
perf_Wp_plot
(S, wB_req, maxSSerror, w_start, w_end, axlim=None, points=1000)¶ MIMO sensitivity S and performance weight Wp plotting funtion.
Parameters: - S : numpy array
Sensitivity transfer function matrix as function of s => S(s)
- wB_req : float
The design or require bandwidth of the plant in rad/time. 1/time eg: wB_req = 1/20sec = 0.05rad/s
- maxSSerror : float
The maximum stead state tracking error required of the plant.
- wStart : float
Minimum power of w for the frequency range in rad/time. eg: for w starting at 10e-3, wStart = -3.
- wEnd : float
Maximum value of w for the frequency range in rad/time. eg: for w ending at 10e3, wStart = 3.
Returns: - wB : float
The actually plant bandwidth in rad/time given the specified controller used to generate the sensitivity matrix S(s).
- Plot : matplotlib figure
-
utilsplot.
ref_perfect_const_plot
(G, R, wr, w_start=-2, w_end=2, axlim=None, points=1000, plot_type='all')¶ Use these plots to determine the constraints for perfect control in terms of combined reference changes. Equation 6.52 (p241) calculates the minimal requirement for input saturation to check in terms of set point tracking. A more tighter bounds is calculated with equation 6.53 (p241).
Parameters: - G : tf
Plant transfer function.
- R : numpy matrix (n x n)
Reference changes (usually diagonal with all elements larger than 1)
- wr : float
Frequency up to witch reference tracking is required
- type_eq : string
Type of plot:
plot_type Type of plot minimal Minimal requirement, equation 6.52 tighter Tighter requirement, equation 6.53 allo All requirements
Returns: - Plot : matplotlib figure
-
utilsplot.
rga_nm_plot
(G, pairing_list=None, pairing_names=None, w_start=-2, w_end=2, axlim=None, points=1000, plot_type='all')¶ Plots the RGA number for a specified pairing
Parameters: - G : numpy matrix
Plant model.
- pairing_list : List of sparse numpy matrices of the same shape as G.
An array of zeros with a 1. at each required output-input pairing. The default is a diagonal pairing with 1.’s on the diagonal.
- plot_type : string
Type of plot:
plot_type Type of plot all All the pairings on one plot element Each pairing has its own plot
Returns: - Plot : matplotlib figure
-
utilsplot.
rga_plot
(G, w_start=-2, w_end=2, axlim=None, points=1000, fig=0, plot_type='elements', input_label=None, output_label=None)¶ Plots the relative gain interaction between each output and input pairing
Parameters: - G : numpy matrix
Plant model.
- plot_type : string
Type of plot.
plot_type Type of plot all All the RGAs on one plot output Plots grouped by output input Plots grouped by input element Each element has its own plot
Returns: - Plot : matplotlib figure
-
utilsplot.
step
(G, t_end=100, initial_val=0, input_label=None, output_label=None, points=1000)¶ This function is similar to the MatLab step function.
Parameters: - G : tf
Plant transfer function.
- t_end : integer
Time period which the step response should occur (optional).
- initial_val : integer
Starting value to evaluate step response (optional).
- input_label : array
List of input variable labels.
- output_label : array
List of output variable labels.
Returns: - Plot : matplotlib figure
-
utilsplot.
step_response_plot
(Y, U, t_end=50, initial_val=0, timedim='sec', axlim=None, points=1000, constraint=None, method='numeric')¶ A plot of the step response of a transfer function
Parameters: - Y : tf
Output transfer function.
- U : tf
Input transfer function.
- t_end : integer
Time period which the step response should occur (optional).
- initial_val : integer
Starting value to evaluate step response (optional).
- constraint : float
The upper limit the step response cannot exceed. is only calculated if a value is specified (optional).
- method : [‘numeric’,’analytic’]
The method that is used to calculate a constrained response. A constraint value is required (optional).
Returns: - Plot : matplotlib figure
-
utilsplot.
sv_dir_plot
(G, plot_type, w_start=-2, w_end=2, axlim=None, points=1000)¶ Plot the input and output singular vectors associated with the minimum and maximum singular values.
Parameters: - G : matrix
Plant model or sensitivity function.
- plot_type : string
Type of plot.
plot_type Type of plot input Plots input vectors output Plots output vectors
Returns: - Plot : matplotlib figure