Skip to content

Commit 0ef2941

Browse files
murrayrmmurrayrm
committed
update time response functions to allow nonlinear systems
1 parent e9e3a8e commit 0ef2941

7 files changed

Lines changed: 158 additions & 154 deletions

File tree

‎control/matlab/timeresp.py‎

Lines changed: 4 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -6,7 +6,7 @@
66

77
__all__ = ['step', 'stepinfo', 'impulse', 'initial', 'lsim']
88

9-
def step(sys, T=None, X0=0., input=0, output=None, return_x=False):
9+
def step(sys, T=None, input=0, output=None, return_x=False):
1010
'''Step response of a linear system
1111
1212
If the system has multiple inputs or outputs (MIMO), one input has
@@ -22,9 +22,6 @@ def step(sys, T=None, X0=0., input=0, output=None, return_x=False):
2222
T: array-like or number, optional
2323
Time vector, or simulation time duration if a number (time vector is
2424
autocomputed if not given)
25-
X0: array-like or number, optional
26-
Initial condition (default = 0)
27-
Numbers are converted to constant arrays with the correct shape.
2825
input: int
2926
Index of the input that will be used in this simulation.
3027
output: int
@@ -55,7 +52,7 @@ def step(sys, T=None, X0=0., input=0, output=None, return_x=False):
5552
from ..timeresp import step_response
5653

5754
# Switch output argument order and transpose outputs
58-
out = step_response(sys, T, X0, input, output,
55+
out = step_response(sys, T, input=input, output=output,
5956
transpose=True, return_x=return_x)
6057
return (out[1], out[0], out[2]) if return_x else (out[1], out[0])
6158

@@ -134,7 +131,7 @@ def stepinfo(sysdata, T=None, yfinal=None, SettlingTimeThreshold=0.02,
134131

135132
return S
136133

137-
def impulse(sys, T=None, X0=0., input=0, output=None, return_x=False):
134+
def impulse(sys, T=None, input=0, output=None, return_x=False):
138135
'''Impulse response of a linear system
139136
140137
If the system has multiple inputs or outputs (MIMO), one input has
@@ -150,10 +147,6 @@ def impulse(sys, T=None, X0=0., input=0, output=None, return_x=False):
150147
T: array-like or number, optional
151148
Time vector, or simulation time duration if a number (time vector is
152149
autocomputed if not given)
153-
X0: array-like or number, optional
154-
Initial condition (default = 0)
155-
156-
Numbers are converted to constant arrays with the correct shape.
157150
input: int
158151
Index of the input that will be used in this simulation.
159152
output: int
@@ -183,7 +176,7 @@ def impulse(sys, T=None, X0=0., input=0, output=None, return_x=False):
183176
from ..timeresp import impulse_response
184177

185178
# Switch output argument order and transpose outputs
186-
out = impulse_response(sys, T, X0, input, output,
179+
out = impulse_response(sys, T, input, output,
187180
transpose = True, return_x=return_x)
188181
return (out[1], out[0], out[2]) if return_x else (out[1], out[0])
189182

@@ -203,8 +196,6 @@ def initial(sys, T=None, X0=0., input=None, output=None, return_x=False):
203196
autocomputed if not given)
204197
X0: array-like object or number, optional
205198
Initial condition (default = 0)
206-
207-
Numbers are converted to constant arrays with the correct shape.
208199
input: int
209200
This input is ignored, but present for compatibility with step
210201
and impulse.

‎control/tests/matlab2_test.py‎

Lines changed: 3 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -126,8 +126,7 @@ def test_step(self, SISO_mats, MIMO_mats, mplcleanup):
126126

127127
subplot2grid(plot_shape, (0, 1))
128128
T = linspace(0, 2, 100)
129-
X0 = array([1, 1])
130-
y, t = step(sys, T, X0)
129+
y, t = step(sys, T)
131130
plot(t, y)
132131

133132
# Test output of state vector
@@ -153,9 +152,8 @@ def test_impulse(self, SISO_mats, mplcleanup):
153152

154153
#supply time and X0
155154
T = linspace(0, 2, 100)
156-
X0 = [0.2, 0.2]
157-
t, y = impulse(sys, T, X0)
158-
plot(t, y, label='t=0..2, X0=[0.2, 0.2]')
155+
t, y = impulse(sys, T)
156+
plot(t, y, label='t=0..2')
159157

160158
#Test system with direct feed-though, the function should print a warning.
161159
D = [[0.5]]

‎control/tests/matlab_test.py‎

Lines changed: 4 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -195,16 +195,7 @@ def testStep(self, siso):
195195
np.testing.assert_array_almost_equal(tout, t)
196196

197197
# Play with arguments
198-
yout, tout = step(sys, T=t, X0=0)
199-
np.testing.assert_array_almost_equal(yout, youttrue, decimal=4)
200-
np.testing.assert_array_almost_equal(tout, t)
201-
202-
X0 = np.array([0, 0])
203-
yout, tout = step(sys, T=t, X0=X0)
204-
np.testing.assert_array_almost_equal(yout, youttrue, decimal=4)
205-
np.testing.assert_array_almost_equal(tout, t)
206-
207-
yout, tout, xout = step(sys, T=t, X0=0, return_x=True)
198+
yout, tout, xout = step(sys, T=t, return_x=True)
208199
np.testing.assert_array_almost_equal(yout, youttrue, decimal=4)
209200
np.testing.assert_array_almost_equal(tout, t)
210201

@@ -249,20 +240,19 @@ def testImpulse(self, siso):
249240
# produce a warning for a system with direct feedthrough
250241
with pytest.warns(UserWarning, match="System has direct feedthrough"):
251242
# Play with arguments
252-
yout, tout = impulse(sys, T=t, X0=0)
243+
yout, tout = impulse(sys, T=t)
253244
np.testing.assert_array_almost_equal(yout, youttrue, decimal=4)
254245
np.testing.assert_array_almost_equal(tout, t)
255246

256247
# produce a warning for a system with direct feedthrough
257248
with pytest.warns(UserWarning, match="System has direct feedthrough"):
258-
X0 = np.array([0, 0])
259-
yout, tout = impulse(sys, T=t, X0=X0)
249+
yout, tout = impulse(sys, T=t)
260250
np.testing.assert_array_almost_equal(yout, youttrue, decimal=4)
261251
np.testing.assert_array_almost_equal(tout, t)
262252

263253
# produce a warning for a system with direct feedthrough
264254
with pytest.warns(UserWarning, match="System has direct feedthrough"):
265-
yout, tout, xout = impulse(sys, T=t, X0=0, return_x=True)
255+
yout, tout, xout = impulse(sys, T=t, return_x=True)
266256
np.testing.assert_array_almost_equal(yout, youttrue, decimal=4)
267257
np.testing.assert_array_almost_equal(tout, t)
268258

‎control/tests/nlsys_test.py‎

Lines changed: 44 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -11,6 +11,7 @@
1111
import numpy as np
1212
import control as ct
1313

14+
# Basic test of nlsys()
1415
def test_nlsys_basic():
1516
def kincar_update(t, x, u, params):
1617
l = params.get('l', 1) # wheelbase
@@ -31,3 +32,46 @@ def kincar_output(t, x, u, params):
3132
assert kincar.input_labels == ['U[0]', 'U[1]']
3233
assert kincar.output_labels == ['y[0]', 'y[1]']
3334
assert kincar.state_labels == ['x', 'y', 'theta']
35+
36+
37+
# Test nonlinear initial, step, and forced response
38+
@pytest.mark.parametrize(
39+
"nin, nout, input, output", [
40+
( 1, 1, None, None),
41+
( 2, 2, None, None),
42+
( 2, 2, 0, None),
43+
( 2, 2, None, 1),
44+
( 2, 2, 1, 0),
45+
])
46+
def test_lti_nlsys_response(nin, nout, input, output):
47+
sys_ss = ct.rss(4, nin, nout, strictly_proper=True)
48+
sys_nl = ct.nlsys(
49+
lambda t, x, u, params: sys_ss.A @ x + sys_ss.B @ u,
50+
lambda t, x, u, params: sys_ss.C @ x + sys_ss.D @ u,
51+
inputs=nin, outputs=nout, states=4)
52+
53+
# Figure out the time to use from the linear impulse response
54+
resp_ss = ct.impulse_response(sys_ss)
55+
timepts = np.linspace(0, resp_ss.time[-1]/10, 100)
56+
57+
# Initial response
58+
resp_ss = ct.initial_response(sys_ss, timepts, output=output)
59+
resp_nl = ct.initial_response(sys_nl, timepts, output=output)
60+
np.testing.assert_equal(resp_ss.time, resp_nl.time)
61+
np.testing.assert_allclose(resp_ss.states, resp_nl.states, atol=0.01)
62+
63+
# Step response
64+
resp_ss = ct.step_response(sys_ss, timepts, input=input, output=output)
65+
resp_nl = ct.step_response(sys_nl, timepts, input=input, output=output)
66+
np.testing.assert_equal(resp_ss.time, resp_nl.time)
67+
np.testing.assert_allclose(resp_ss.states, resp_nl.states, atol=0.01)
68+
69+
# Forced response
70+
X0 = np.linspace(0, 1, sys_ss.nstates)
71+
U = np.zeros((nin, timepts.size))
72+
for i in range(nin):
73+
U[i] = 0.01 * np.sin(timepts + i)
74+
resp_ss = ct.forced_response(sys_ss, timepts, U, X0=X0)
75+
resp_nl = ct.forced_response(sys_nl, timepts, U, X0=X0)
76+
np.testing.assert_equal(resp_ss.time, resp_nl.time)
77+
np.testing.assert_allclose(resp_ss.states, resp_nl.states, atol=0.05)

‎control/tests/timeresp_test.py‎

Lines changed: 13 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -486,9 +486,7 @@ def test_step_pole_cancellation(self, tsystem):
486486
@pytest.mark.parametrize(
487487
"tsystem, kwargs",
488488
[("siso_ss2", {}),
489-
("siso_ss2", {'X0': 0}),
490-
("siso_ss2", {'X0': np.array([0, 0])}),
491-
("siso_ss2", {'X0': 0, 'return_x': True}),
489+
("siso_ss2", {'return_x': True}),
492490
("siso_dtf0", {})],
493491
indirect=["tsystem"])
494492
def test_impulse_response_siso(self, tsystem, kwargs):
@@ -567,9 +565,9 @@ def test_initial_response_mimo(self, tsystem):
567565
yref = tsystem.yinitial
568566
yref_notrim = np.broadcast_to(yref, (2, len(t)))
569567

570-
_t, y_00 = initial_response(sys, T=t, X0=x0, input=0, output=0)
568+
_t, y_00 = initial_response(sys, T=t, X0=x0, output=0)
571569
np.testing.assert_array_almost_equal(y_00, yref, decimal=4)
572-
_t, y_11 = initial_response(sys, T=t, X0=x0, input=0, output=1)
570+
_t, y_11 = initial_response(sys, T=t, X0=x0, output=1)
573571
np.testing.assert_array_almost_equal(y_11, yref, decimal=4)
574572
_t, yy = initial_response(sys, T=t, X0=x0)
575573
np.testing.assert_array_almost_equal(yy, yref_notrim, decimal=4)
@@ -874,7 +872,8 @@ def test_time_vector(self, tsystem, fun, squeeze):
874872
kw['U'] = np.vstack([np.sin(t) for i in range(sys.ninputs)])
875873
elif fun == forced_response and isctime(sys, strict=True):
876874
pytest.skip("No continuous forced_response without time vector.")
877-
if hasattr(sys, "nstates") and sys.nstates is not None:
875+
if hasattr(sys, "nstates") and sys.nstates is not None and \
876+
fun != impulse_response:
878877
kw['X0'] = np.arange(sys.nstates) + 1
879878
if sys.ninputs > 1 and fun in [step_response, impulse_response]:
880879
kw['input'] = 1
@@ -1047,8 +1046,9 @@ def test_squeeze(self, fcn, nstate, nout, ninp, squeeze, shape1, shape2):
10471046
_, yvec, xvec = ct.initial_response(
10481047
sys, tvec, 1, squeeze=squeeze, return_x=True)
10491048
assert xvec.shape == (sys.nstates, 8)
1050-
else:
1051-
_, yvec = ct.initial_response(sys, tvec, 1, squeeze=squeeze)
1049+
elif isinstance(sys, TransferFunction):
1050+
with pytest.warns(UserWarning, match="may not be consistent"):
1051+
_, yvec = ct.initial_response(sys, tvec, 1, squeeze=squeeze)
10521052
assert yvec.shape == shape2
10531053

10541054
# Forced response (only indexed by output)
@@ -1090,7 +1090,11 @@ def test_squeeze(self, fcn, nstate, nout, ninp, squeeze, shape1, shape2):
10901090
if squeeze is not True or sys.ninputs > 1 or sys.noutputs > 1:
10911091
assert yvec.shape == (sys.noutputs, sys.ninputs, 8)
10921092

1093-
_, yvec = ct.initial_response(sys, tvec, 1)
1093+
if isinstance(sys, TransferFunction):
1094+
with pytest.warns(UserWarning, match="may not be consistent"):
1095+
_, yvec = ct.initial_response(sys, tvec, 1)
1096+
else:
1097+
_, yvec = ct.initial_response(sys, tvec, 1)
10941098
if squeeze is not True or sys.noutputs > 1:
10951099
assert yvec.shape == (sys.noutputs, 8)
10961100

‎control/tests/trdata_test.py‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -243,7 +243,7 @@ def test_trdata_labels():
243243
assert step_response.input_labels == [sys.input_labels[nu]]
244244
assert step_response.output_labels == [sys.output_labels[ny]]
245245

246-
init_response = ct.initial_response(sys, T, input=nu, output=ny)
246+
init_response = ct.initial_response(sys, T, output=ny)
247247
assert init_response.input_labels == None
248248
assert init_response.output_labels == [sys.output_labels[ny]]
249249

0 commit comments

Comments
 (0)