engineeringkit.control_systems#

Domain hub · History · Examples · API reference

Classical stability, PID responses, state feedback, LQR and estimation.

engineeringkit.control_systems: feedback, stability, frequency response, state space and optimal control.

Black’s feedback amplifier and Minorsky’s PID controller with Ziegler-Nichols tuning; Maxwell’s characteristic roots, the Routh and Hurwitz tests and Lyapunov’s direct method; Nyquist’s encirclement criterion, Bode plots and stability margins; Evans’s root locus; Hazen’s second-order servo specifications; zero-order-hold discretisation; Kalman’s controllability and observability, pole placement and Luenberger observers; and the LQR, discrete LQR and Kalman-Bucy Riccati designs. Models are continuous-time and SISO unless stated; SciPy supplies the LTI objects, Riccati/Lyapunov solvers, pole placement and step responses.

class engineeringkit.control_systems.PID(kp, ki=0.0, kd=0.0)[source]#

Bases: object

Ideal parallel PID controller \(C(s) = K_p + K_i/s + K_d s\).

Derivative action is ideal (no roll-off filter), and there is no actuator saturation or anti-windup. With ki == 0 no integrator is included.

Examples

>>> from engineeringkit.control_systems import transfer_function
>>> r = PID(1, 1).simulate(transfer_function([1], [1, 1]), np.linspace(0, 5, 51))
>>> round(r.summary()["final"], 4)
0.9933
Parameters:

  • kp (float)

  • ki (float)

  • kd (float)

closed_loop(plant)[source]#

Unity-feedback closed loop \(CG/(1 + CG)\) from reference to output.

kd: float = 0.0#

Derivative gain [input s/error].

Type:

float

ki: float = 0.0#

Integral gain [input/(error s)].

Type:

float

kp: float#

Proportional gain [input/error].

Type:

float

simulate(plant, T=None)[source]#

Unit-step reference response of the closed loop, as a ResponseResult.

transfer_function()[source]#

Controller transfer function \(C(s)\); all-zero gains have none.

class engineeringkit.control_systems.ResponseResult(t, y)[source]#

Bases: object

Time response of a simulated control loop.

Parameters:

  • t (ndarray)

  • y (ndarray)

plot(ax=None)[source]#

Plot the response against time; return (fig, ax).

summary()[source]#

Final value and peak of the response.

Return type:

dict

t: ndarray#

Time samples [s].

Type:

ndarray

y: ndarray#

Output at each sample, in output units per unit reference.

Type:

ndarray

engineeringkit.control_systems.as_polynomial(coeffs, name='coeffs', min_degree=0)[source]#

Validate a descending-power polynomial and strip leading zeros.

Parameters:

  • coeffs (array_like) – Real coefficients, highest power first.

  • name (str, optional) – Argument name used in error messages.

  • min_degree (int, optional) – Smallest admissible degree after stripping leading zeros.

Returns:

Coefficients with a nonzero leading term.

Return type:

ndarray

Raises:

ValueError – If the coefficients are not a finite real vector, are all zero, or the degree is below min_degree.

Examples

>>> as_polynomial([0, 1, 3, 2])
array([1., 3., 2.])
engineeringkit.control_systems.bode(sys, w=None)[source]#

Bode magnitude and phase of a continuous-time model.

Parameters:

  • sys (scipy.signal.lti) – Continuous-time model.

  • w (array_like, optional) – Angular frequencies [rad/s]; chosen automatically if omitted.

Returns:

(w, magnitude, phase): angular frequency [rad/s], gain [dB] and phase [deg].

Return type:

tuple of ndarray

Notes

Wraps scipy.signal.bode(). See Control systems: matrix contracts and numerical stability.

References

H. W. Bode, Network Analysis and Feedback Amplifier Design (Van Nostrand, 1945).

Examples

>>> from engineeringkit.control_systems import transfer_function
>>> w, mag, phase = bode(transfer_function([1], [1, 1]), [1.0])
>>> round(float(mag[0]), 3), round(float(phase[0]), 1)
(-3.01, -45.0)
engineeringkit.control_systems.controllability_matrix(A, B)[source]#

Kalman controllability matrix \([B,\ AB,\ \dots,\ A^{n-1}B]\).

Parameters:

  • A (array_like) – State matrix, shape (n, n).

  • B (array_like) – Input matrix, shape (n, m).

Returns:

Shape (n, n m).

Return type:

ndarray

Notes

See Control systems: matrix contracts and numerical stability.

References

R. E. Kalman, “On the general theory of control systems,” Proc. 1st IFAC Congress (Moscow, 1960) 481-492.

Examples

>>> controllability_matrix([[0, 1], [0, 0]], [[0], [1]]).tolist()
[[0.0, 1.0], [1.0, 0.0]]
engineeringkit.control_systems.dc_gain(sys)[source]#

Steady-state gain \(G(0)\) for a constant input.

Parameters:

sys (scipy.signal.lti) – SISO continuous-time model.

Returns:

Output per unit input at zero frequency.

Return type:

float

Raises:

ValueError – If the model has a pole at \(s = 0\) (an integrator), so the gain is unbounded.

Notes

The steady state is reached only if the system is stable. See Control systems: matrix contracts and numerical stability.

Examples

>>> dc_gain(transfer_function([2, 2], [1, 3, 2]))
1.0
engineeringkit.control_systems.discretize_zoh(A, B, dt)[source]#

Exact zero-order-hold discretisation \((A_d, B_d)\) of \(\dot x = Ax + Bu\).

Parameters:

  • A (array_like) – Continuous state matrix, shape (n, n) [1/s].

  • B (array_like) – Continuous input matrix, shape (n, m).

  • dt (float) – Sampling period [s], positive.

Returns:

(Ad, Bd) for \(x_{k+1} = A_d x_k + B_d u_k\).

Return type:

tuple of ndarray

Notes

Both matrices come from one matrix exponential of the augmented matrix \(\begin{bmatrix}A & B\\ 0 & 0\end{bmatrix} T\), which is valid even when \(A\) is singular. See Control systems: matrix contracts and numerical stability.

References

J. R. Ragazzini and L. A. Zadeh, “The analysis of sampled-data systems,” Trans. AIEE 71 (1952) 225-234.

Examples

>>> Ad, Bd = discretize_zoh([[0.0]], [[1.0]], 0.5)
>>> Ad.tolist(), Bd.tolist()
([[1.0]], [[0.5]])
engineeringkit.control_systems.dlqr(A, B, Q, R)[source]#

Discrete-time infinite-horizon linear-quadratic regulator.

Minimises \(\sum_{k \ge 0} (x_k^\mathsf{T} Q x_k + u_k^\mathsf{T} R u_k)\) over \(u_k = -K x_k\) for \(x_{k+1} = A x_k + B u_k\); the steady state of Bellman’s backward dynamic-programming recursion.

Parameters:

  • A (array_like) – Discrete state matrix, shape (n, n).

  • B (array_like) – Discrete input matrix, shape (n, m).

  • Q (array_like) – Symmetric positive semidefinite state weight, shape (n, n).

  • R (array_like) – Symmetric positive definite input weight, shape (m, m).

Returns:

(K, P, poles): gain \(K = (R + B^\mathsf{T}PB)^{-1}B^\mathsf{T}PA\), the discrete Riccati solution \(P\) and the closed-loop eigenvalues (inside the unit circle when stabilising).

Return type:

tuple

Notes

See Control systems: matrix contracts and numerical stability.

References

  1. Bellman, Dynamic Programming (Princeton University Press, 1957).

Examples

>>> K, P, poles = dlqr([[1.0]], [[1.0]], [[1.0]], [[1.0]])
>>> round(float(P[0, 0]), 6), round(float(K[0, 0]), 6)
(1.618034, 0.618034)
engineeringkit.control_systems.feedback(G, H=1.0, sign=-1)[source]#

Closed-loop transfer function of a feedback loop, \(G/(1 - \sigma GH)\).

Parameters:

  • G (scipy.signal.lti) – Forward-path model.

  • H (scipy.signal.lti or float, optional) – Feedback-path model or constant gain (unity feedback by default).

  • sign ({-1, 1}, optional) – Feedback sign \(\sigma\); -1 is negative feedback.

Returns:

Closed-loop model from reference to output.

Return type:

scipy.signal.TransferFunction

Notes

Pole-zero cancellations are not removed. See Control systems: matrix contracts and numerical stability.

References

H. S. Black, “Stabilized feedback amplifiers,” Bell Syst. Tech. J. 13 (1934) 1-18.

Examples

>>> T = feedback(transfer_function([1000], [1]), 0.01)
>>> round(dc_gain(T), 3)
90.909
engineeringkit.control_systems.feedforward_gain(A, B, C, K)[source]#

Reference precompensator \(N\) giving unit DC gain under SISO state feedback.

With \(u = -Kx + Nr\), the closed loop \(\dot x = (A - BK)x + BNr\) tracks a constant reference when \(N = -1/\left(C (A - BK)^{-1} B\right)\).

Parameters:

  • A (array_like) – State matrix, shape (n, n).

  • B (array_like) – Input matrix, shape (n, 1).

  • C (array_like) – Output matrix, shape (1, n).

  • K (array_like) – State-feedback gain, shape (1, n).

Returns:

Precompensator gain [input/reference units].

Return type:

float

Raises:

ValueError – If the system is not SISO, or the closed loop has a zero at the origin so no finite \(N\) exists.

Notes

Exact only for a perfect model; integral action is needed for robust tracking. See Control systems: matrix contracts and numerical stability.

Examples

>>> feedforward_gain([[-1]], [[1]], [[1]], [[2]])
3.0
engineeringkit.control_systems.hurwitz_determinants(coeffs)[source]#

Leading principal minors \(\Delta_1, \dots, \Delta_n\) of the Hurwitz matrix.

For \(a_0 s^n + a_1 s^{n-1} + \dots + a_n\) with \(a_0 > 0\), all roots lie in the open left half-plane if and only if every \(\Delta_k > 0\). The Hurwitz matrix has entries \(H_{ij} = a_{2j-i}\) (1-based, \(a_k = 0\) outside \(0 \le k \le n\)).

Parameters:

coeffs (array_like) – Real coefficients in descending powers of \(s\), degree at least one. They are normalised by the leading coefficient.

Returns:

The \(n\) determinants.

Return type:

ndarray

Notes

On the stability boundary \(\Delta_{n-1} = 0\) (a pair of imaginary roots) or \(\Delta_n = 0\) (a root at the origin). See Control systems: matrix contracts and numerical stability.

References

A. Hurwitz, “Ueber die Bedingungen, unter welchen eine Gleichung nur Wurzeln mit negativen reellen Theilen besitzt,” Math. Ann. 46 (1895) 273-284.

Examples

>>> hurwitz_determinants([1, 3, 3, 1]).round(6).tolist()
[3.0, 8.0, 8.0]
engineeringkit.control_systems.impulse_response(sys, T=None)[source]#

Response to a unit impulse input from rest.

Parameters:

  • sys (scipy.signal.lti) – Continuous-time model.

  • T (array_like, optional) – Time samples [s]; chosen automatically if omitted.

Returns:

(t, h): time [s] and impulse response [output/(input s)].

Return type:

tuple of ndarray

Notes

Wraps scipy.signal.impulse(). See Control systems: matrix contracts and numerical stability.

Examples

>>> t, h = impulse_response(transfer_function([1], [1, 1]), [0.0, 1.0])
>>> round(float(h[-1]), 4)
0.3679
engineeringkit.control_systems.is_controllable(A, B)[source]#

Whether the controllability matrix has full state rank.

Parameters:

  • A (array_like) – State matrix, shape (n, n).

  • B (array_like) – Input matrix, shape (n, m).

Returns:

True if every state can be steered by the inputs.

Return type:

bool

Notes

Rank is computed numerically with SVD tolerances. See Control systems: matrix contracts and numerical stability.

Examples

>>> is_controllable([[-1, 0], [0, -1]], [[1], [1]])
False
engineeringkit.control_systems.is_observable(A, C)[source]#

Whether the observability matrix has full state rank.

Parameters:

  • A (array_like) – State matrix, shape (n, n).

  • C (array_like) – Output matrix, shape (p, n).

Returns:

True if every state can be reconstructed from the outputs.

Return type:

bool

Notes

See Control systems: matrix contracts and numerical stability.

Examples

>>> is_observable([[0, 1], [0, 0]], [[1, 0]])
True
engineeringkit.control_systems.kalman_gain(A, C, Q, R)[source]#

Steady-state Kalman-Bucy estimator gain for continuous-time measurements.

For \(\dot x = Ax + w\), \(y = Cx + v\) with white process and measurement noise of spectral densities \(Q\) and \(R\), the estimator \(\dot{\hat x} = A\hat x + L(y - C\hat x)\) with \(L = PC^\mathsf{T}R^{-1}\) minimises the steady error covariance \(P\).

Parameters:

  • A (array_like) – State matrix, shape (n, n) [1/s].

  • C (array_like) – Output matrix, shape (p, n).

  • Q (array_like) – Process-noise spectral density, symmetric positive semidefinite, shape (n, n).

  • R (array_like) – Measurement-noise spectral density, symmetric positive definite, shape (p, p).

Returns:

(L, P, poles): observer gain, error covariance and estimator poles [1/s].

Return type:

tuple

Notes

\(P\) solves \(AP + PA^\mathsf{T} - PC^\mathsf{T}R^{-1}CP + Q = 0\). \(Q\) and \(R\) are spectral densities, not sampled covariances. See Control systems: matrix contracts and numerical stability.

References

R. E. Kalman and R. S. Bucy, “New results in linear filtering and prediction theory,” J. Basic Eng. 83 (1961) 95-108.

Examples

>>> L, P, poles = kalman_gain([[0]], [[1]], [[4]], [[1]])
>>> L.round(6).tolist(), P.round(6).tolist()
([[2.0]], [[2.0]])
engineeringkit.control_systems.lqr(A, B, Q, R)[source]#

Continuous-time linear-quadratic regulator.

Minimises \(J = \int_0^\infty (x^\mathsf{T} Q x + u^\mathsf{T} R u)\,dt\) over \(u = -Kx\) for \(\dot x = Ax + Bu\).

Parameters:

  • A (array_like) – State matrix, shape (n, n) [1/s].

  • B (array_like) – Input matrix, shape (n, m).

  • Q (array_like) – Symmetric positive semidefinite state weight, shape (n, n).

  • R (array_like) – Symmetric positive definite input weight, shape (m, m).

Returns:

(K, P, poles): gain \(K = R^{-1} B^\mathsf{T} P\), the Riccati solution \(P\) (so the optimal cost from \(x_0\) is \(x_0^\mathsf{T} P x_0\)) and the closed-loop poles [1/s].

Return type:

tuple

Notes

\(P\) solves \(A^\mathsf{T}P + PA - PBR^{-1}B^\mathsf{T}P + Q = 0\). A stabilising solution needs \((A, B)\) stabilisable and \((Q^{1/2}, A)\) detectable; SciPy errors are not masked. See Control systems: matrix contracts and numerical stability.

References

R. E. Kalman, “Contributions to the theory of optimal control,” Bol. Soc. Mat. Mexicana 5 (1960) 102-119.

Examples

>>> K, P, poles = lqr([[0]], [[1]], [[4]], [[1]])
>>> K.round(6).tolist(), poles.real.round(6).tolist()
([[2.0]], [-2.0])
engineeringkit.control_systems.lyapunov_equation(A, Q=None)[source]#

Solve the continuous Lyapunov equation \(A^\mathsf{T} P + P A = -Q\).

If \(A\) is stable, \(V(x) = x^\mathsf{T} P x\) is a Lyapunov function: it is positive and decreases at the rate \(\dot V = -x^\mathsf{T} Q x\) along every trajectory of \(\dot x = Ax\).

Parameters:

  • A (array_like) – State matrix, shape (n, n) [1/s].

  • Q (array_like, optional) – Symmetric positive semidefinite decay-rate matrix, shape (n, n); the identity by default.

Returns:

Symmetric solution \(P\).

Return type:

ndarray

Raises:

ValueError – If two eigenvalues of \(A\) sum to zero, so the equation has no unique solution.

Notes

See Control systems: matrix contracts and numerical stability.

References

A. M. Lyapunov, The General Problem of the Stability of Motion (Kharkov, 1892; English translation, Taylor & Francis, 1992).

Examples

>>> lyapunov_equation([[-1.0]]).tolist()
[[0.5]]
engineeringkit.control_systems.lyapunov_stable(A)[source]#

Asymptotic stability of \(\dot x = Ax\) by Lyapunov’s direct method.

Parameters:

A (array_like) – State matrix, shape (n, n) [1/s].

Returns:

True if \(A^\mathsf{T} P + PA = -I\) has a positive definite solution \(P\).

Return type:

bool

Notes

A singular equation (eigenvalues summing to zero, e.g. an undamped oscillator) means the origin is not asymptotically stable. See Control systems: matrix contracts and numerical stability.

Examples

>>> lyapunov_stable([[0, 1], [-2, -3]]), lyapunov_stable([[0, 1], [-1, 0]])
(True, False)
engineeringkit.control_systems.num_den(sys)[source]#

Numerator and denominator polynomials of a SISO continuous-time model.

Parameters:

sys (scipy.signal.lti) – A TransferFunction, or any SciPy LTI model that converts to one (state space, zeros-poles-gain).

Returns:

(num, den) in descending powers of \(s\).

Return type:

tuple of ndarray

Raises:

ValueError – If the model is not single-input single-output.

Examples

>>> num_den(signal.TransferFunction([2, 2], [1, 3, 2]))
(array([2., 2.]), array([1., 3., 2.]))
engineeringkit.control_systems.nyquist(sys, w)[source]#

Complex frequency response \(G(j\omega)\).

Parameters:

  • sys (scipy.signal.lti) – SISO continuous-time model.

  • w (float or array_like) – Angular frequencies [rad/s]; negative values give the mirror image.

Returns:

Output/input gain at each frequency.

Return type:

complex or ndarray of complex

Notes

See Control systems: matrix contracts and numerical stability.

Examples

>>> from engineeringkit.control_systems import transfer_function
>>> complex(np.round(nyquist(transfer_function([2, 2], [1, 3, 2]), 2.0), 12))
(0.5-0.5j)
engineeringkit.control_systems.nyquist_encirclements(sys, w=None)[source]#

Net clockwise encirclements of \(-1\) by the Nyquist plot of a loop.

Nyquist’s criterion: the closed loop \(L/(1+L)\) has \(Z = N + P\) right-half-plane poles, where \(P\) is the number of open-loop right-half-plane poles and \(N\) the value returned here.

Parameters:

  • sys (scipy.signal.lti) – Proper open-loop model \(L(s)\) with no poles on the imaginary axis.

  • w (array_like, optional) – Positive angular frequencies [rad/s] spanning all features of \(L\); defaults to 40001 points from \(10^{-4}\) to \(10^4\).

Returns:

Clockwise encirclements of \(-1\) (negative if counter-clockwise).

Return type:

int

Raises:

ValueError – If \(L\) is improper, has imaginary-axis poles, or passes through \(-1\) on the grid.

Notes

The plot is traced for \(\omega\) from \(-\infty\) to \(\infty\); for a proper loop the infinite semicircle maps to a single point. A grid that misses a fast phase swing can miscount. See Control systems: matrix contracts and numerical stability.

References

H. Nyquist, “Regeneration theory,” Bell Syst. Tech. J. 11 (1932) 126-147.

Examples

>>> from engineeringkit.control_systems import transfer_function
>>> nyquist_encirclements(transfer_function([12], [1, 3, 3, 1]))
2
engineeringkit.control_systems.observability_matrix(A, C)[source]#

Kalman observability matrix \([C;\ CA;\ \dots;\ CA^{n-1}]\).

Parameters:

  • A (array_like) – State matrix, shape (n, n).

  • C (array_like) – Output matrix, shape (p, n).

Returns:

Shape (n p, n).

Return type:

ndarray

Notes

See Control systems: matrix contracts and numerical stability.

Examples

>>> observability_matrix([[0, 1], [0, 0]], [[1, 0]]).tolist()
[[1.0, 0.0], [0.0, 1.0]]
engineeringkit.control_systems.observer_gain(A, C, desired_poles)[source]#

Luenberger observer gain \(L\) placing the eigenvalues of \(A - LC\).

The observer \(\dot{\hat x} = A\hat x + Bu + L(y - C\hat x)\) makes the estimation error obey \(\dot e = (A - LC)e\).

Parameters:

  • A (array_like) – State matrix, shape (n, n).

  • C (array_like) – Output matrix, shape (p, n); (A, C) must be observable.

  • desired_poles (array_like) – n error-dynamics poles [1/s], conjugate pairs for complex values.

Returns:

Gain \(L\), shape (n, p).

Return type:

ndarray

Notes

Computed by duality, placing the poles of \(A^\mathsf{T} - C^\mathsf{T} L^\mathsf{T}\). See Control systems: matrix contracts and numerical stability.

References

D. G. Luenberger, “Observing the state of a linear system,” IEEE Trans. Mil. Electron. 8 (1964) 74-80.

Examples

>>> observer_gain([[0, 1], [0, 0]], [[1, 0]], [-1, -2]).round(6).tolist()
[[3.0], [2.0]]
engineeringkit.control_systems.poles(sys)[source]#

Poles of a SISO model, the roots of its denominator.

Parameters:

sys (scipy.signal.lti) – SISO continuous-time model.

Returns:

Poles [1/s].

Return type:

ndarray of complex

Notes

The system is asymptotically stable when every pole has a negative real part. See Control systems: matrix contracts and numerical stability.

Examples

>>> np.sort(poles(transfer_function([1], [1, 3, 2]))).real.tolist()
[-2.0, -1.0]
engineeringkit.control_systems.root_locus(num, den, gains)[source]#

Closed-loop characteristic roots of \(D(s) + K N(s)\) for each gain.

Parameters:

  • num (array_like) – Open-loop numerator and denominator, descending powers of \(s\). The open loop must be strictly proper so the number of roots is constant.

  • den (array_like) – Open-loop numerator and denominator, descending powers of \(s\). The open loop must be strictly proper so the number of roots is constant.

  • gains (array_like) – Nonempty vector of loop gains \(K\) [-].

Returns:

Shape (len(gains), len(den) - 1); row i holds the roots for gains[i]. Roots within a row are not tracked along branches.

Return type:

ndarray of complex

Notes

See Control systems: matrix contracts and numerical stability.

References

W. R. Evans, “Graphical analysis of control systems,” Trans. AIEE 67 (1948) 547-551.

Examples

>>> np.sort(root_locus([1], [1, 3, 2], [0.25])[0].real).tolist()
[-1.5, -1.5]
engineeringkit.control_systems.routh_first_column(coeffs)[source]#

First column of the Routh array of a real polynomial.

The number of sign changes down the column equals the number of roots in the open right half-plane.

Parameters:

coeffs (array_like) – Real coefficients in descending powers of \(s\), with a nonzero leading coefficient and degree at least one.

Returns:

First-column entries for rows \(s^n\) down to \(s^0\), normalised so the leading coefficient is one.

Return type:

ndarray

Notes

A row of zeros is replaced by the derivative of the auxiliary polynomial formed from the row above. A zero pivot is replaced by a small positive \(\varepsilon\), so this is a numerical teaching method rather than a symbolic proof for roots on or near the imaginary axis. See Control systems: matrix contracts and numerical stability.

References

E. J. Routh, A Treatise on the Stability of a Given State of Motion (Macmillan, 1877).

Examples

>>> routh_first_column([1, 3, 3, 1]).round(4).tolist()
[1.0, 3.0, 2.6667, 1.0]
engineeringkit.control_systems.routh_stable(coeffs)[source]#

Whether every root of a polynomial lies strictly in the left half-plane.

Parameters:

coeffs (array_like) – Real coefficients in descending powers of \(s\).

Returns:

True if the Routh column has no sign change and no root lies on the imaginary axis; marginal stability counts as unstable.

Return type:

bool

Notes

Routh signs alone cannot distinguish the auxiliary-row marginal case, so root locations confirm it with a relative tolerance. See Control systems: matrix contracts and numerical stability.

Examples

>>> routh_stable([1, 3, 3, 1]), routh_stable([1, 0, 1])
(True, False)
engineeringkit.control_systems.second_order_metrics(zeta, omega_n)[source]#

Peak time, percent overshoot and 2 % settling time of an underdamped second-order step response.

Parameters:

  • zeta (float or array_like) – Damping ratio [-], strictly between 0 and 1.

  • omega_n (float or array_like) – Undamped natural frequency [rad/s], positive.

Returns:

peak_time \(\pi/\omega_d\) [s], percent_overshoot \(100\,e^{-\pi\zeta/\sqrt{1-\zeta^2}}\) [%] and settling_time_2pct \(4/(\zeta\omega_n)\) [s], with \(\omega_d = \omega_n\sqrt{1-\zeta^2}\).

Return type:

dict

Notes

Peak time and overshoot are exact; the settling time is the usual envelope estimate, not the exact 2 % crossing. See Control systems: matrix contracts and numerical stability.

References

H. L. Hazen, “Theory of servo-mechanisms,” J. Franklin Inst. 218 (1934) 279-331.

Examples

>>> m = second_order_metrics(0.5, 2.0)
>>> round(float(m["percent_overshoot"]), 2), float(m["settling_time_2pct"])
(16.3, 4.0)
engineeringkit.control_systems.stability_margins(sys, w=None)[source]#

Gain and phase margins from the first crossovers of a loop’s frequency response.

Parameters:

  • sys (scipy.signal.lti) – SISO open-loop model \(L(s)\).

  • w (array_like, optional) – Positive angular frequencies [rad/s]; defaults to 20000 points from \(10^{-4}\) to \(10^4\).

Returns:

gain_margin [-] (inf without a phase crossover), phase_margin_deg [deg] (nan without a gain crossover), gain_cross_frequency and phase_cross_frequency [rad/s].

Return type:

dict

Notes

Crossings are interpolated linearly in log-frequency between grid points. Only the first crossover of each kind is reported, so loops with multiple crossovers, time delays or negative gain need separate analysis. See Control systems: matrix contracts and numerical stability.

References

H. W. Bode, “Relations between attenuation and phase in feedback amplifier design,” Bell Syst. Tech. J. 19 (1940) 421-454.

Examples

>>> from engineeringkit.control_systems import transfer_function
>>> m = stability_margins(transfer_function([2], [1, 3, 2, 0]))
>>> round(m["gain_margin"], 2), round(m["phase_cross_frequency"], 3)
(3.0, 1.414)
engineeringkit.control_systems.state_feedback_poles(A, B, desired_poles)[source]#

State-feedback gain \(K\) placing the eigenvalues of \(A - BK\).

Parameters:

  • A (array_like) – State matrix, shape (n, n).

  • B (array_like) – Input matrix, shape (n, m); (A, B) must be controllable.

  • desired_poles (array_like) – n closed-loop poles [1/s]; complex poles must come in conjugate pairs.

Returns:

Gain \(K\), shape (m, n), for the control law \(u = -Kx\).

Return type:

ndarray

Notes

Wraps scipy.signal.place_poles(). See Control systems: matrix contracts and numerical stability.

References

W. M. Wonham, “On pole assignment in multi-input controllable linear systems,” IEEE Trans. Autom. Control 12 (1967) 660-665.

Examples

>>> state_feedback_poles([[0, 1], [0, 0]], [[0], [1]], [-1, -2]).round(6).tolist()
[[2.0, 3.0]]
engineeringkit.control_systems.state_matrices(A, B=None, C=None)[source]#

Validate and convert state-space matrices to float arrays.

Parameters:

  • A (array_like) – Square state matrix, shape (n, n) [1/s].

  • B (array_like, optional) – Input matrix, shape (n, m).

  • C (array_like, optional) – Output matrix, shape (p, n).

Returns:

(A, B, C), with None for matrices not supplied.

Return type:

tuple of ndarray

Raises:

ValueError – If a matrix is nonfinite, empty or has an inconsistent shape.

Examples

>>> A, B, C = state_matrices([[0, 1], [-2, -3]], B=[[0], [1]])
>>> A.shape, B.shape, C
((2, 2), (2, 1), None)
engineeringkit.control_systems.state_space(A, B, C, D)[source]#

Continuous-time state-space model \(\dot x = Ax + Bu,\ y = Cx + Du\).

Parameters:

  • A (array_like) – State matrix, shape (n, n) [1/s].

  • B (array_like) – Input matrix, shape (n, m).

  • C (array_like) – Output matrix, shape (p, n).

  • D (array_like) – Feedthrough matrix, shape (p, m).

Returns:

Continuous-time model.

Return type:

scipy.signal.StateSpace

Notes

Matrix units follow the chosen state, input and output units; see Control systems: matrix contracts and numerical stability.

Examples

>>> sys = state_space([[-2]], [[2]], [[1]], [[0]])
>>> float((-sys.C @ np.linalg.solve(sys.A, sys.B) + sys.D).item())
1.0
engineeringkit.control_systems.step_response(sys, T=None)[source]#

Response to a unit step input from rest.

Parameters:

  • sys (scipy.signal.lti) – Continuous-time model.

  • T (array_like, optional) – Time samples [s]; chosen automatically if omitted.

Returns:

(t, y): time [s] and output per unit step.

Return type:

tuple of ndarray

Notes

Wraps scipy.signal.step(). See Control systems: matrix contracts and numerical stability.

Examples

>>> t, y = step_response(transfer_function([1], [1, 1]), [0.0, 1.0])
>>> round(float(y[-1]), 4)
0.6321
engineeringkit.control_systems.transfer_function(num, den)[source]#

Continuous-time transfer function \(N(s)/D(s)\).

Parameters:

  • num (array_like) – Numerator and denominator coefficients in descending powers of \(s\). The denominator must have degree at least zero and must not be identically zero.

  • den (array_like) – Numerator and denominator coefficients in descending powers of \(s\). The denominator must have degree at least zero and must not be identically zero.

Returns:

Continuous-time SISO model.

Return type:

scipy.signal.TransferFunction

Notes

Units of \(G\) are output units per input unit. See Control systems: matrix contracts and numerical stability.

Examples

>>> G = transfer_function([1], [1, 1])
>>> G.den.tolist()
[1.0, 1.0]
engineeringkit.control_systems.weight_matrix(matrix, size, name, strict=False)[source]#

Validate a symmetric positive (semi)definite weighting or covariance matrix.

Parameters:

  • matrix (array_like) – Candidate matrix, shape (size, size).

  • size (int) – Required dimension.

  • name (str) – Argument name used in error messages.

  • strict (bool, optional) – Require positive definite (True) rather than semidefinite.

Returns:

The matrix as a float array.

Return type:

ndarray

Raises:

ValueError – If the matrix has the wrong shape, is not symmetric, or is indefinite.

Examples

>>> weight_matrix([[2.0]], 1, "R", strict=True)
array([[2.]])
engineeringkit.control_systems.zeros(sys)[source]#

Zeros of a SISO model, the roots of its numerator.

Parameters:

sys (scipy.signal.lti) – SISO continuous-time model.

Returns:

Zeros [1/s].

Return type:

ndarray of complex

Notes

See Control systems: matrix contracts and numerical stability.

Examples

>>> zeros(transfer_function([2, 2], [1, 3, 2])).real.tolist()
[-1.0]
engineeringkit.control_systems.ziegler_nichols(ku, tu, controller='PID')[source]#

Ziegler-Nichols ultimate-sensitivity tuning.

Raise a proportional gain until the loop oscillates steadily: that gain is the ultimate gain \(K_u\) and the oscillation period is \(T_u\). The classic rules then set

type

\(K_p\)

\(T_i\)

\(T_d\)

P

\(0.5K_u\)

PI

\(0.45K_u\)

\(T_u/1.2\)

PID

\(0.6K_u\)

\(T_u/2\)

\(T_u/8\)

with \(K_i = K_p/T_i\) and \(K_d = K_p T_d\).

Parameters:

  • ku (float) – Ultimate proportional gain [input/error], positive.

  • tu (float) – Ultimate oscillation period [s], positive.

  • controller ({"P", "PI", "PID"}, optional) – Controller structure.

Returns:

Tuned controller.

Return type:

PID

Notes

The rules aim for roughly quarter-amplitude decay and are only a starting point; they give large overshoot on many plants. See Control systems: matrix contracts and numerical stability.

References

J. G. Ziegler and N. B. Nichols, “Optimum settings for automatic controllers,” Trans. ASME 64 (1942) 759-768.

Examples

>>> ziegler_nichols(8.0, 4.0)
PID(kp=4.8, ki=2.4, kd=2.4)