8000 Merge pull request #889 from SCLiao47/feature_bandwidth · python-control/python-control@15ed92c · GitHub
[go: up one dir, main page]
Skip to content

Commit 15ed92c

Browse files
murrayrmmurrayrm
authored
Merge pull request #889 from SCLiao47/feature_bandwidth
bandwidth feature
2 parents c06a205 + ab68562 commit 15ed92c

File tree

4 files changed

+144
-3
lines changed

4 files changed

+144
-3
lines changed

‎control/lti.py

Lines changed: 108 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -20,7 +20,7 @@
2020
from .namedio import NamedIOSystem, isdtime
2121

2222
__all__ = ['poles', 'zeros', 'damp', 'evalfr', 'frequency_response',
23-
'freqresp', 'dcgain', 'pole', 'zero']
23+
'freqresp', 'dcgain', 'bandwidth', 'pole', 'zero']
2424

2525

2626
class LTI(NamedIOSystem):
@@ -202,6 +202,68 @@ def _dcgain(self, warn_infinite):
202202
else:
203203
return zeroresp
204204

205+
def bandwidth(self, dbdrop=-3):
206+
"""Evaluate the bandwidth of the LTI system for a given dB drop.
207+
208+
Evaluate the first frequency that the response magnitude is lower than
209+
DC gain by dbdrop dB.
210+
211+
Parameters
212+
----------
213+
dpdrop : float, optional
214+
A strictly negative scalar in dB (default = -3) defines the
215+
amount of gain drop for deciding bandwidth.
216+
217+
Returns
218+
-------
219+
bandwidth : ndarray
220+
The first frequency (rad/time-unit) where the gain drops below
221+
dbdrop of the dc gain of the system, or nan if the system has
222+
infinite dc gain, inf if the gain does not drop for all frequency
223+
224+
Raises
225+
------
226+
TypeError
227+
if 'sys' is not an SISO LTI instance
228+
ValueError
229+
if 'dbdrop' is not a negative scalar
230+
"""
231+
# check if system is SISO and dbdrop is a negative scalar
232+
if not self.issiso():
233+
raise TypeError("system should be a SISO system")
234+
235+
if (not np.isscalar(dbdrop)) or dbdrop >= 0:
236+
raise ValueError("expecting dbdrop be a negative scalar in dB")
237+
238+
dcgain = self.dcgain()
239+
if np.isinf(dcgain):
240+
# infinite dcgain, return np.nan
241+
return np.nan
242+
243+
# use frequency range to identify the 0-crossing (dbdrop) bracket
244+
from control.freqplot import _default_frequency_range
245+
omega = _default_frequency_range(self)
246+
mag, phase, omega = self.frequency_response(omega)
247+
idx_dropped = np.nonzero(mag - dcgain*10**(dbdrop/20) < 0)[0]
248+
249+
if idx_dropped.shape[0] == 0:
250+
# no frequency response is dbdrop below the dc gain, return np.inf
251+
return np.inf
252+
else:
253+
# solve for the bandwidth, use scipy.optimize.root_scalar() to
254+
# solve using bisection
255+
import scipy
256+
result = scipy.optimize.root_scalar(
257+
lambda w: np.abs(self(w*1j)) - np.abs(dcgain)*10**(dbdrop/20),
258+
bracket=[omega[idx_dropped[0] - 1], omega[idx_dropped[0]]],
259+
method='bisect')
260+
261+
# check solution
262+
if result.converged:
263+
return np.abs(result.root)
264+
else:
265+
raise Exception(result.message)
266+
205267
def ispassive(self):
206268
# importing here prevents circular dependancy
207269
from control.passivity import ispassive
@@ -499,6 +561,51 @@ def dcgain(sys):
499561
return sys.dcgain()
500562

501563

564+
def bandwidth(sys, dbdrop=-3):
565+
"""Return the first freqency where the gain drop by dbdrop of the system.
566+
567+
Parameters
568+
----------
569+
sys: StateSpace or TransferFunction
570+
Linear system
571+
dbdrop : float, optional
572+
By how much the gain drop in dB (default = -3) that defines the
573+
bandwidth. Should be a negative scalar
574+
575+
Returns
576+
-------
577+
bandwidth : ndarray
578+
The first frequency (rad/time-unit) where the gain drops below dbdrop
579+
of the dc gain of the system, or nan if the system has infinite dc
580+
gain, inf if the gain does not drop for all frequency
581+
582+
Raises
583+
------
584+
TypeError
585+
if 'sys' is not an SISO LTI instance
586+
ValueError
587+
if 'dbdrop' is not a negative scalar
588+
589+
Example
590+
-------
591+
>>> G = ct.tf([1], [1, 1])
592+
>>> ct.bandwidth(G)
593+
0.9976
594+
595+
>>> G1 = ct.tf(0.1, [1, 0.1])
596+
>>> wn2 = 1
597+
>>> zeta2 = 0.001
598+
>>> G2 = ct.tf(wn2**2, [1, 2*zeta2*wn2, wn2**2])
599+
>>> ct.bandwidth(G1*G2)
600+
0.1018
601+
602+
"""
603+
if not isinstance(sys, LTI):
604+
raise TypeError("sys must be a LTI instance.")
605+
606+
return sys.bandwidth(dbdrop)
607+
608+
502609
# Process frequency responses in a uniform way
503610
def _process_frequency_response(sys, omega, out, squeeze=None):
504611
# Set value of squeeze argument if not set

‎control/matlab/__init__.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -189,7 +189,7 @@
189189
190190
== ========================== ============================================
191191
\* :func:`dcgain` steady-state (D.C.) gain
192-
\ lti/bandwidth system bandwidth
192+
\* :func:`bandwidth` system bandwidth
193193
\ lti/norm h2 and Hinfinity norms of LTI models
194194
\* :func:`pole` system poles
195195
\* :func:`zero` system (transmission) zeros

‎control/tests/lti_test.py

Lines changed: 33 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -6,7 +6,7 @@
66

77
import control as ct
88
from control import c2d, tf, ss, tf2ss, NonlinearIOSystem
9-
from control.lti import LTI, evalfr, damp, dcgain, zeros, poles
9+
from control.lti import LTI, evalfr, damp, dcgain, zeros, poles, bandwidth
1010
from control import common_timebase, isctime, isdtime, issiso, timebaseEqual
1111
from control.tests.conftest import slycotonly
1212
from control.exception import slycot_check
@@ -104,6 +104,38 @@ def test_dcgain(self):
104104
np.testing.assert_allclose(sys.dcgain(), 42)
105105
np.testing.assert_allclose(dcgain(sys), 42)
106106

107+
def test_bandwidth(self):
108+
# test a first-order system, compared with matlab
109+
sys1 = tf(0.1, [1, 0.1])
110+
np.testing.assert_allclose(sys1.bandwidth(), 0.099762834511098)
111+
np.testing.assert_allclose(bandwidth(sys1), 0.099762834511098)
112+
113+
# test a second-order system, compared with matlab
114+
wn2 = 1
115+
zeta2 = 0.001
116+
sys2 = sys1 * tf(wn2**2, [1, 2*zeta2*wn2, wn2**2])
117+
np.testing.assert_allclose(sys2.bandwidth(), 0.101848388240241)
118+
np.testing.assert_allclose(bandwidth(sys2), 0.101848388240241)
119+
120+
# test constant gain, bandwidth should be infinity
121+
sysAP = tf(1,1)
122+
np.testing.assert_allclose(bandwidth(sysAP), np.inf)
123+
124+
# test integrator, bandwidth should return np.nan
125+
sysInt = tf(1, [1, 0])
126+
np.testing.assert_allclose(bandwidth(sysInt), np.nan)
127+
128+
# test exception for system other than LTI
129+
np.testing.assert_raises(TypeError, bandwidth, 1)
130+
131+
# test exception for system other than SISO system
132+
sysMIMO = tf([[[-1, 41], [1]], [[1, 2], [3, 4]]],
133+
[[[1, 10], [1, 20]], [[1, 30], [1, 40]]])
134+
np.testing.assert_raises(TypeError, bandwidth, sysMIMO)
135+
136+
# test if raise exception if dbdrop is positive scalar
137+
np.testing.assert_raises(ValueError, bandwidth, sys1, 3)
138+
107139
@pytest.mark.parametrize("dt1, dt2, expected",
108140
[(None, None, True),
109141
(None, 0, True),

‎control/xferfcn.py

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -74,6 +74,7 @@
7474
'xferfcn.floating_point_format': '.4g'
7575
}
7676

77+
7778
def _float2str(value):
7879
_num_format = config.defaults.get('xferfcn.floating_point_format', ':.4g')
7980
return f"{value:{_num_format}}"
@@ -1407,6 +1408,7 @@ def _tf_factorized_polynomial_to_string(roots, gain=1, var='s'):
14071408

14081409
return multiplier + " ".join(factors)
14091410

1411+
14101412
def _tf_string_to_latex(thestr, var='s'):
14111413
""" make sure to superscript all digits in a polynomial string
14121414
and convert float coefficients in scientific notation

0 commit comments

Comments
 (0)
0