8000 Merge pull request #982 from murrayrm/save_params-30Mar2024 · python-control/python-control@c8d9c8f · GitHub
[go: up one dir, main page]
Skip to content

Commit c8d9c8f

Browse files
murrayrmmurrayrm
authored
Merge pull request #982 from murrayrm/save_params-30Mar2024
Allow passing and saving of params in time responses
2 parents ea3e6a5 + 85568a9 commit c8d9c8f

File tree

4 files changed

+84
-14
lines changed

4 files changed

+84
-14
lines changed

‎control/nlsys.py

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1367,6 +1367,8 @@ def input_output_response(
13671367
13681368
* inputs (array): Input(s) to the system, indexed by input and time.
13691369
1370+
* params (dict): Parameters values used for the simulation.
1371+
13701372
The return value of the system can also be accessed by assigning the
13711373
function to a tuple of length 2 (time, output) or of length 3 (time,
13721374
output, state) if ``return_x`` is ``True``. If the input/output
@@ -1643,7 +1645,7 @@ def ivp_rhs(t, x):
16431645
raise TypeError("Can't determine system type")
16441646

16451647
return TimeResponseData(
1646-
soln.t, y, soln.y, u, issiso=sys.issiso(),
1648+
soln.t, y, soln.y, u, params=params, issiso=sys.issiso(),
16471649
output_labels=sys.output_labels, input_labels=sys.input_labels,
16481650
state_labels=sys.state_labels, sysname=sys.name,
16491651
title="Input/output response for " + sys.name,

‎control/tests/nlsys_test.py

Lines changed: 17 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -75,3 +75,20 @@ def test_lti_nlsys_response(nin, nout, input, output):
7575
resp_nl = ct.forced_response(sys_nl, timepts, U, X0=X0)
7676
np.testing.assert_equal(resp_ss.time, resp_nl.time)
7777
np.testing.assert_allclose(resp_ss.states, resp_nl.states, atol=0.05)
78+
79+
80+
# Test to make sure that impulse responses are not allowed
81+
def test_nlsys_impulse():
82+
sys_ss = ct.rss(4, 1, 1, strictly_proper=True)
83+
sys_nl = ct.nlsys(
84+
lambda t, x, u, params: sys_ss.A @ x + sys_ss.B @ u,
85+
lambda t, x, u, params: sys_ss.C @ x + sys_ss.D @ u,
86+
inputs=1, outputs=1, states=4)
87+
88+
# Figure out the time to use from the linear impulse response
89+
resp_ss = ct.impulse_response(sys_ss)
90+
timepts = np.linspace(0, resp_ss.time[-1]/10, 100)
91+
92+
# Impulse_response (not implemented)
93+
with pytest.raises(ValueError, match="system must be LTI"):
94+
resp_nl = ct.impulse_response(sys_nl, timepts)

‎control/tests/trdata_test.py

Lines changed: 33 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -236,7 +236,8 @@ def test_trdata_labels():
236236
np.testing.assert_equal(
237237
response.input_labels, ["u[%d]" % i for i in range(sys.ninputs)])
238238

239-
# Make sure the selected input and output are both correctly transferred to the response
239+
# Make sure the selected input and output are both correctly
240+
# transferred to the response
240241
for nu in range(sys.ninputs):
241242
for ny in range(sys.noutputs):
242243
step_response = ct.step_response(sys, T, input=nu, output=ny)
@@ -339,6 +340,37 @@ def test_trdata_multitrace():
339340
np.zeros(5), np.ones(5), np.zeros((1, 5)), np.zeros(6))
340341

341342

343+
@pytest.mark.parametrize("func, args", [
344+
(ct.step_response, ()),
345+
(ct.initial_response, (1, )),
346+
(ct.forced_response, (0, 1)),
347+
(ct.input_output_response, (0, 1)),
348+
])
349+
@pytest.mark.parametrize("dt", [0, 1])
350+
def test_trdata_params(func, args, dt):
351+
# Create a nonlinear system with parameters, neutrally stable
352+
nlsys = ct.nlsys(
353+
lambda t, x, u, params: params['a'] * x[0] + u[0],
354+
states = 1, inputs = 1, outputs = 1, params={'a': 0}, dt=dt)
355+
lnsys = ct.ss([[-0.5]], [[1]], [[1]], 0, dt=dt)
356+
357+
# Compute the response, setting parameters to make things stable
358+
timevec = np.linspace(0, 1) if dt == 0 else np.arange(0, 10, 1)
359+
nlresp = func(nlsys, timevec, *args, params={'a': -0.5})
360+
lnresp = func(lnsys, timevec, *args)
361+
362+
# Make sure the modified system was stable
363+
np.testing.assert_allclose(
364+
nlresp.states, lnresp.states, rtol=1e-3, atol=1e-5)
365+
assert lnresp.params == None
366+
assert nlresp.params['a'] == -0.5
367+
368+
# Make sure the match was not accidental
369+
bdresp = func(nlsys, timevec, *args)
370+
assert not np.allclose(
371+
bdresp.states, nlresp.states, rtol=1e-3, atol=1e-5)
372+
373+
342374
def test_trdata_exceptions():
343375
# Incorrect dimension for time vector
344376
with pytest.raises(ValueError, match="Time vector must be 1D"):

‎control/timeresp.py

Lines changed: 31 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -173,6 +173,9 @@ class TimeResponseData:
173173
sysname : str, optional
174174
Name of the system that created the data.
175175
176+
params : dict, optional
177+
If system is a nonlinear I/O system, set parameter values.
178+
176179
plot_inputs : bool, optional
177180
Whether or not to plot the inputs by default (can be overridden in
178181
the plot() method)
@@ -227,7 +230,7 @@ def __init__(
227230
output_labels=None, state_labels=None, input_labels=None,
228231
title=None, transpose=False, return_x=False, squeeze=None,
229232
multi_trace=False, trace_labels=None, trace_types=None,
230-
plot_inputs=True, sysname=None
233+
plot_inputs=True, sysname=None, params=None
231234
):
232235
"""Create an input/output time response object.
233236
@@ -315,6 +318,7 @@ def __init__(
315318
raise ValueError("Time vector must be 1D array")
316319
self.title = title
317320
self.sysname = sysname
321+
self.params = params
318322

319323
#
320324
# Output vector (and number of traces)
@@ -856,7 +860,7 @@ def shape_matches(s_legal, s_actual):
856860

857861

858862
# Forced response of a linear system
859-
def forced_response(sys, T=None, U=0., X0=0., transpose=False,
863+
def forced_response(sys, T=None, U=0., X0=0., transpose=False, params=None,
860864
interpolate=False, return_x=None, squeeze=None):
861865
"""Compute the output of a linear system given the input.
862866
@@ -889,6 +893,9 @@ def forced_response(sys, T=None, U=0., X0=0., transpose=False,
889893
X0 : array_like or float, default=0.
890894
Initial condition.
891895
896+
params : dict, optional
897+
If system is a nonlinear I/O system, set parameter values.
898+
892899
transpose : bool, default=False
893900
If True, transpose all input and output arrays (for backward
894901
compatibility with MATLAB and :func:`scipy.signal.lsim`).
@@ -981,7 +988,7 @@ def forced_response(sys, T=None, U=0., X0=0., transpose=False,
981988
warnings.warn(
982989
"interpolation not supported for nonlinear I/O systems")
983990
return input_output_response(
984-
sys, T, U, X0, transpose=transpose,
991+
sys, T, U, X0, params=params, transpose=transpose,
985992
return_x=return_x, squeeze=squeeze)
986993
else:
987994
raise TypeError('Parameter ``sys``: must be a ``StateSpace`` or'
@@ -1169,7 +1176,7 @@ def forced_response(sys, T=None, U=0., X0=0., transpose=False,
11691176
yout = np.transpose(yout)
11701177

11711178
return TimeResponseData(
1172-
tout, yout, xout, U, issiso=sys.issiso(),
1179+
tout, yout, xout, U, params=params, issiso=sys.issiso(),
11731180
output_labels=sys.output_labels, input_labels=sys.input_labels,
11741181
state_labels=sys.state_labels, sysname=sys.name, plot_inputs=True,
11751182
title="Forced response for " + sys.name, trace_types=['forced'],
@@ -1256,7 +1263,7 @@ def _process_time_response(
12561263

12571264

12581265
def step_response(sys, T=None, X0=0, input=None, output=None, T_num=None,
1259-
transpose=False, return_x=False, squeeze=None):
1266+
transpose=False, return_x=False, squeeze=None, params=None):
12601267
# pylint: disable=W0622
12611268
"""Compute the step response for a linear system.
12621269
@@ -1298,6 +1305,9 @@ def step_response(sys, T=None, X0=0, input=None, output=None, T_num=None,
12981305
Only report the step response for the listed output. If not
12991306
specified, all outputs are reported.
13001307
1308+
params : dict, optional
1309+
If system is a nonlinear I/O system, set parameter values.
1310+
13011311
T_num : int, optional
13021312
Number of time steps to use in simulation if T is not provided as an
13031313
array (autocomputed if not given); ignored if sys is discrete-time.
@@ -1403,7 +1413,7 @@ def step_response(sys, T=None, X0=0, input=None, output=None, T_num=None,
14031413
U = np.zeros((sys.ninputs, T.size))
14041414
U[i, :] = np.ones_like(T)
14051415

1406-
response = forced_response(sys, T, U, X0, squeeze=True)
1416+
response = forced_response(sys, T, U, X0, squeeze=True, params=params)
14071417
inpidx = i if input is None else 0
14081418
yout[:, inpidx, :] = response.y if output is None \
14091419
else response.y[output]
@@ -1424,11 +1434,11 @@ def step_response(sys, T=None, X0=0, input=None, output=None, T_num=None,
14241434
output_labels=output_labels, input_labels=input_labels,
14251435
state_labels=sys.state_labels, title="Step response for " + sys.name,
14261436
transpose=transpose, return_x=return_x, squeeze=squeeze,
1427-
sysname=sys.name, trace_labels=trace_labels,
1437+
sysname=sys.name, params=params, trace_labels=trace_labels,
14281438
trace_types=trace_types, plot_inputs=False)
14291439

14301440

1431-
def step_info(sysdata, T=None, T_num=None, yfinal=None,
1441+
def step_info(sysdata, T=None, T_num=None, yfinal=None, params=None,
14321442
SettlingTimeThreshold=0.02, RiseTimeLimits=(0.1, 0.9)):
14331443
"""
14341444
Step response characteristics (Rise time, Settling Time, Peak and others).
@@ -1451,6 +1461,8 @@ def step_info(sysdata, T=None, T_num=None, yfinal=None,
14511461
systems to simulate and the last value of the the response data is
14521462
used for a given time series of response data. Scalar for SISO,
14531463
(noutputs, ninputs) array_like for MIMO systems.
1464+
params : dict, optional
1465+
If system is a nonlinear I/O system, set parameter values.
14541466
SettlingTimeThreshold : float, optional
14551467
Defines the error to compute settling time (default = 0.02)
14561468
RiseTimeLimits : tuple (lower_threshold, upper_theshold)
@@ -1536,7 +1548,7 @@ def step_info(sysdata, T=None, T_num=None, yfinal=None,
15361548
from .nlsys import NonlinearIOSystem
15371549

15381550
if isinstance(sysdata, (StateSpace, TransferFunction, NonlinearIOSystem)):
1539-
T, Yout = step_response(sysdata, T, squeeze=False)
1551+
T, Yout = step_response(sysdata, T, squeeze=False, params=params)
15401552
if yfinal:
15411553
InfValues = np.atleast_2d(yfinal)
15421554
else:
@@ -1651,7 +1663,7 @@ def step_info(sysdata, T=None, T_num=None, yfinal=None,
16511663
return ret[0][0] if retsiso else ret
16521664

16531665

1654-
def initial_response(sys, T=None, X0=0, output=None, T_num=None,
1666+
def initial_response(sys, T=None, X0=0, output=None, T_num=None, params=None,
16551667
transpose=False, return_x=False, squeeze=None):
16561668
# pylint: disable=W0622
16571669
"""Compute the initial condition response for a linear system.
@@ -1684,6 +1696,9 @@ def initial_response(sys, T=None, X0=0, output=None, T_num=None,
16841696
Number of time steps to use in simulation if T is not provided as an
16851697
array (autocomputed if not given); ignored if sys is discrete-time.
16861698
1699+
params : dict, optional
1700+
If system is a nonlinear I/O system, set parameter values.
1701+
16871702
transpose : bool, optional
16881703
If True, transpose all input and output arrays (for backward
16891704
compatibility with MATLAB and :func:`scipy.signal.lsim`). Default
@@ -1748,7 +1763,7 @@ def initial_response(sys, T=None, X0=0, output=None, T_num=None,
17481763
raise ValueError("invalid value of T for this type of system")
17491764

17501765
# Compute the forced response
1751-
response = forced_response(sys, T, 0, X0)
1766+
response = forced_response(sys, T, 0, X0, params=params)
17521767

17531768
# Figure out if the system is SISO or not
17541769
issiso = sys.issiso() or output is not None
@@ -1760,7 +1775,7 @@ def initial_response(sys, T=None, X0=0, output=None, T_num=None,
17601775

17611776
# Store the response without an input
17621777
return TimeResponseData(
1763-
response.t, yout, response.x, None, issiso=issiso,
1778+
response.t, yout, response.x, None, params=params, issiso=issiso,
17641779
output_labels=output_labels, input_labels=None,
17651780
state_labels=sys.state_labels, sysname=sys.name,
17661781
title="Initial response for " + sys.name, trace_types=['initial'],
@@ -1863,6 +1878,10 @@ def impulse_response(sys, T=None, input=None, output=None, T_num=None,
18631878
from .statesp import _convert_to_statespace
18641879
from .lti import LTI
18651880

1881+
# Make sure we have an LTI system
1882+
if not isinstance(sys, LTI):
1883+
raise ValueError("system must be LTI system for impulse response")
1884+
18661885
# Create the time and input vectors
18671886
if T is None or np.asarray(T).size == 1:
18681887
T = _default_time_vector(sys, N=T_num, tfinal=T, is_step=False)

0 commit comments

Comments
 (0)
0