Open
Description
A bug was reported in python-control issue #935 in which tf2ss returns the incorrect poles.
It looks like the underlying problem is in td04ad. The following code illustrates the bug:
from numpy import array, roots, sort
from numpy.linalg import eig
from numpy.testing import assert_allclose
from slycot import td04ad
from scipy.signal import tf2ss
# Test case #1 (works)
num = array([[[1, 1, 1, 1, 1]]])
den = array([[1, 2, 1, 2, 1]])
denorder = np.array([den.size - 1])
ssout = td04ad('C', 1, 1, denorder, den, num, tol=0)
tf_poles = sort(roots(den[0]))
ss_poles = sort(eig(ssout[1])[0])
assert_allclose(tf_poles, ss_poles)
# Test case #2 (doesn't work)
num = array([[[1.05351324e-02, 2.83623596e-01, 2.44891698e+00, 1.22011448e+01,
5.34846570e+01, 1.42092310e+02, 2.24627039e+02, 2.54305034e+02,
2.18224432e+02, 1.24102950e+02, 4.99448809e+01, 1.26651596e+01,
2.02239050e+00, 1.06681479e-01, 1.56212741e-03, 8.90622718e-06,
1.77205330e-08]]])
den = array([[1.00000000e+00, 6.41551095e+02, 8.40213016e+03, 5.55121354e+04,
1.93428590e+05, 1.31800818e+05, 9.20183802e+04, 2.90211403e+04,
4.54824418e+03, 4.08636727e+02, 1.72367976e+01, 8.91452881e-02,
1.04426479e-04, 4.32051569e-08, 5.88752872e-12, 3.24913608e-16,
6.34764375e-21]])
denorder = np.array([den.size - 1])
# Try scipy
A, B, C, D = tf2ss(num[0, 0], den[0])
tf_poles = sort(roots(den[0]))
ss_poles = sort(eig(A)[0])
assert_allclose(tf_poles, ss_poles)
# Now try slycot
ssout = td04ad('C', 1, 1, denorder, den, num, tol=0)
ss_poles = sort(eig(ssout[1])[0])
assert_allclose(tf_poles, ss_poles)
For the second system, the scipy-generated poles match, but some of the the slycot-generated poles are not the same (and one of them is unstable).
It's possible this is a numerical issue (the coefficients span 21 orders of magnitude!), but odd the SciPy gets it right and slycot doesn't. (SciPy uses an observable canonical form-based realization.)
Metadata
Metadata
Assignees
Labels
No labels