FIX: transfer_to_state tripping over numerical noise

harold/_classes.py

@@ -1,6 +1,6 @@
 import numpy as np
 import warnings
-from numpy import zeros_like, kron, ndarray, zeros, exp, convolve
+from numpy import zeros_like, kron, ndarray, zeros, exp, convolve, spacing
 from numpy.random import rand, choice
 from numpy.linalg.linalg import _assertNdSquareness
 from scipy.linalg import (eigvals, svdvals, block_diag, qz, norm, solve, expm,
@@ -2898,6 +2898,9 @@ def transfer_to_state(G, output='system'):
         else:
             # Watch out for full cancellation !!
             NumOrEmpty, datanum = haroldpolydiv(num, den)
+            # Clean up the tiny entries
+            datanum[np.abs(datanum) < spacing(100.)] = 0.
+
             # If all cancelled datanum is returned empty
             if datanum.size == 0:
                 A = None
@@ -2934,6 +2937,9 @@ def transfer_to_state(G, output='system'):
                     pass  # D[x,y] is already 0.
                 else:
                     NumOrEmpty, datanum = haroldpolydiv(datanum, dataden)
+                    # Clean up the tiny entries
+                    datanum[np.abs(datanum) < spacing(100.)] = 0.
+
                     # Case 3: If all cancelled datanum is returned empty
                     if np.count_nonzero(datanum) == 0:
                         D[x, y] = NumOrEmpty

harold/tests/test_classes.py

@@ -8,9 +8,9 @@
 from numpy.testing import (assert_,
                            assert_equal,
                            assert_array_equal,
-                           assert_raises,
                            assert_almost_equal,
                            assert_array_almost_equal)
+from pytest import raises as assert_raises
 
 
 def test_concatenate_state_matrices():
@@ -759,8 +759,8 @@ def test_random_state_model():
     assert not G._isgain
     assert not G._isSISO
 
-    G = random_state_model(5, 2, 4, stable=False)
-    assert np.any(G.poles.real > 0)
+    G = random_state_model(5, 2, 4, stable=True)
+    assert not (G.poles.real > 0).any()
     G = random_state_model(11, stable=False, prob_dist=[0, 0, 0.5, 0.5])
     assert_array_almost_equal(np.abs(G.poles.real), np.zeros(11))
     assert np.any(G.poles.imag)
@@ -834,6 +834,20 @@ def test_transfer_to_state():
     assert_array_almost_equal(Gss.c, des_c)
     assert_array_almost_equal(Gss.d, np.zeros([2, 2]))
 
+    # reported in gh-#50
+    num = [[[61.79732492202783, 36.24988430260625, 0.7301196233698941],
+            [0.0377840674057878, 0.9974993795127982, 21.763622825733773]]]
+    den = [[[84.64, 18.4, 1.0], [1.0, 7.2, 144.0]]]
+
+    TF = transfer_to_state((num, den))
+    assert_array_almost_equal([-3.6-1.14472704e+01j,
+                               -3.6+1.14472704e+01j,
+                               -0.10869565-1.74405237e-07j,
+                               -0.10869565+1.74405237e-07j,
+                               ],
+                              np.sort(TF.poles))
+    assert TF.zeros.size == 0
+
 
 def test_state_to_transfer():
     G = State(-2*np.eye(2), np.eye(2), [[1, -3], [0, 0]], [[0, 1], [-1, 0]])