Skip to content

Commit 131c318

Browse files
murrayrmmurrayrm
committed
add documentation on predict keyword + input_output_response list processing
1 parent 86447f0 commit 131c318

3 files changed

Lines changed: 121 additions & 9 deletions

File tree

‎control/iosys.py‎

Lines changed: 64 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -1585,11 +1585,17 @@ def input_output_response(
15851585
T : array-like
15861586
Time steps at which the input is defined; values must be evenly spaced.
15871587
1588-
U : array-like or number, optional
1589-
Input array giving input at each time `T` (default = 0).
1590-
1591-
X0 : array-like or number, optional
1592-
Initial condition (default = 0).
1588+
U : array-like, list, or number, optional
1589+
Input array giving input at each time `T` (default = 0). If a list
1590+
is specified, each element in the list will be treated as a portion
1591+
of the input and broadcast (if necessary) to match the time vector.
1592+
1593+
X0 : array-like, list, or number, optional
1594+
Initial condition (default = 0). If a list is given, each element
1595+
in the list will be flattened and stacked into the initial
1596+
condition. If a smaller number of elements are given that the
1597+
number of states in the system, the initial condition will be padded
1598+
with zeros.
15931599
15941600
return_x : bool, optional
15951601
If True, return the state vector when assigning to a tuple (default =
@@ -1641,6 +1647,16 @@ def input_output_response(
16411647
ValueError
16421648
If time step does not match sampling time (for discrete time systems).
16431649
1650+
Notes
1651+
-----
1652+
1. If a smaller number of initial conditions are given than the number of
1653+
states in the system, the initial conditions will be padded with
1654+
zeros. This is often useful for interconnected control systems where
1655+
the process dynamics are the first system and all other components
1656+
start with zero initial condition since this can be specified as
1657+
[xsys_0, 0]. A warning is issued if the initial conditions are padded
1658+
and and the final listed initial state is not zero.
1659+
16441660
"""
16451661
#
16461662
# Process keyword arguments
@@ -1669,15 +1685,57 @@ def input_output_response(
16691685
n_steps = len(T)
16701686

16711687
# Check and convert the input, if needed
1672-
# TODO: improve MIMO ninputs check (choose from U)
16731688
if sys.ninputs is None or sys.ninputs == 1:
16741689
legal_shapes = [(n_steps,), (1, n_steps)]
16751690
else:
16761691
legal_shapes = [(sys.ninputs, n_steps)]
1692+
1693+
# If we were passed a list of input, concatenate them (w/ broadcast)
1694+
if isinstance(U, (tuple, list)):
1695+
U_elements = []
1696+
for i, u in enumerate(U):
1697+
u = np.array(u) # convert everyting to an array
1698+
# Process this input
1699+
if u.ndim == 0 or (u.ndim == 1 and u.shape[0] != T.shape[0]):
1700+
# Broadcast array to the length of the time input
1701+
u = np.outer(u, np.ones_like(T))
1702+
1703+
elif (u.ndim == 1 and u.shape[0] == T.shape[0]) or \
1704+
(u.ndim == 2 and u.shape[1] == T.shape[0]):
1705+
# No processing necessary; just stack
1706+
pass
1707+
1708+
else:
1709+
raise ValueError(f"Input element {i} has inconsistent shape")
1710+
1711+
# Append this input to our list
1712+
U_elements.append(u)
1713+
1714+
# Save the newly created input vector
1715+
U = np.vstack(U_elements)
1716+
1717+
# Make sure the input has the right shape
16771718
U = _check_convert_array(U, legal_shapes,
16781719
'Parameter ``U``: ', squeeze=False)
16791720
U = U.reshape(-1, n_steps)
16801721

1722+
# If we were passed a list of initial states, concatenate them
1723+
if isinstance(X0, (tuple, list)):
1724+
X0_list = []
1725+
for i, x0 in enumerate(X0):
1726+
x0 = np.array(x0).reshape(-1) # convert everyting to 1D array
1727+
X0_list += x0.tolist() # add elements to initial state
1728+
1729+
# Save the newly created input vector
1730+
X0 = np.array(X0_list)
1731+
1732+
# If the initial state is too short, make it longer (NB: sys.nstates
1733+
# could be None if nstates comes from size of initial condition)
1734+
if sys.nstates and isinstance(X0, np.ndarray) and X0.size < sys.nstates:
1735+
if X0[-1] != 0:
1736+
warn("initial state too short; padding with zeros")
1737+
X0 = np.hstack([X0, np.zeros(sys.nstates - X0.size)])
1738+
16811739
# Check to make sure this is not a static function
16821740
nstates = _find_size(sys.nstates, X0)
16831741
if nstates == 0:

‎control/stochsys.py‎

Lines changed: 15 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -323,9 +323,11 @@ def create_estimator_iosystem(
323323
324324
estim = ct.create_estimator_iosystem(sys, QN, RN)
325325
326-
where ``sys`` is the process dynamics and QN and RN are the covariance of
327-
the disturbance noise and sensor noise. The function returns the
328-
estimator ``estim`` as I/O systems.
326+
where ``sys`` is the process dynamics and QN and RN are the covariance
327+
of the disturbance noise and sensor noise. The function returns the
328+
estimator ``estim`` as I/O system with a parameter ``correct`` that can
329+
be used to turn off the correction term in the estimation (for forward
330+
predictions).
329331
330332
Parameters
331333
----------
@@ -368,6 +370,16 @@ def create_estimator_iosystem(
368370
est = ct.create_estimator_iosystem(sys, QN, RN, P0)
369371
ctrl, clsys = ct.create_statefbk_iosystem(sys, K, estimator=est)
370372
373+
The estimator can also be run on its own to process a noisy signal:
374+
375+
resp = ct.input_output_response(est, T, [Y, U], [X0, P0])
376+
377+
If desired, the ``correct`` parameter can be set to ``False`` to allow
378+
prediction with no additional sensor information:
379+
380+
resp = ct.input_output_response(
381+
est, T, 0, [X0, P0], param={'correct': False)
382+
371383
"""
372384

373385
# Make sure that we were passed an I/O system as an input

‎control/tests/iosys_test.py‎

Lines changed: 42 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1697,3 +1697,45 @@ def test_interconnect_unused_output():
16971697
inputs=['r'],
16981698
outputs=['y'],
16991699
ignore_outputs=['v'])
1700+
1701+
1702+
def test_input_output_broadcasting():
1703+
# Create a system, time vector, and noisy input
1704+
sys = ct.rss(6, 2, 3)
1705+
T = np.linspace(0, 10, 10)
1706+
U = np.zeros((sys.ninputs, T.size))
1707+
U[0, :] = np.sin(T)
1708+
U[1, :] = np.zeros_like(U[1, :])
1709+
U[2, :] = np.ones_like(U[2, :])
1710+
X0 = np.array([1, 2])
1711+
P0 = np.array([[3.11, 3.12], [3.21, 3.3]])
1712+
1713+
# Simulate the system with nominal input to establish baseline
1714+
resp_base = ct.input_output_response(
1715+
sys, T, U, np.hstack([X0, P0.reshape(-1)]))
1716+
1717+
# Split up the inputs into two pieces
1718+
resp_inp1 = ct.input_output_response(sys, T, [U[:1], U[1:]], [X0, P0])
1719+
np.testing.assert_equal(resp_base.states, resp_inp1.states)
1720+
1721+
# Specify two of the inputs as constants
1722+
resp_inp2 = ct.input_output_response(sys, T, [U[0], 0, 1], [X0, P0])
1723+
np.testing.assert_equal(resp_base.states, resp_inp2.states)
1724+
1725+
# Specify two of the inputs as constant vector
1726+
resp_inp3 = ct.input_output_response(sys, T, [U[0], [0, 1]], [X0, P0])
1727+
np.testing.assert_equal(resp_base.states, resp_inp3.states)
1728+
1729+
# Specify only some of the initial conditions
1730+
resp_init = ct.input_output_response(sys, T, [U[0], [0, 1]], [X0, 0])
1731+
resp_cov0 = ct.input_output_response(sys, T, U, [X0, P0 * 0])
1732+
np.testing.assert_equal(resp_cov0.states, resp_init.states)
1733+
1734+
# Specify only some of the initial conditions
1735+
with pytest.warns(UserWarning, match="initial state too short; padding"):
1736+
resp_short = ct.input_output_response(sys, T, [U[0], [0, 1]], [X0, 1])
1737+
1738+
# Make sure that inconsistent settings don't work
1739+
with pytest.raises(ValueError, match="inconsistent"):
1740+
resp_bad = ct.input_output_response(
1741+
sys, T, (U[0, :], U[:2, :-1]), [X0, P0])

0 commit comments

Comments
 (0)