handle non proper tf in _common_den() (#370)

bnavigator · web-flow · commit a09d059e57b0 · 2020-03-21T12:51:06.000-07:00
diff --git a/control/tests/xferfcn_test.py b/control/tests/xferfcn_test.py
@@ -546,6 +546,26 @@ def test_common_den(self):
                                              np.zeros((3, 5, 6)))
         np.testing.assert_array_almost_equal(den, denref)
 
+    def test_common_den_nonproper(self):
+        """ Test _common_den with order(num)>order(den) """
+
+        tf1 = TransferFunction(
+                [[[1., 2., 3.]], [[1., 2.]]],
+                [[[1., -2.]], [[1., -3.]]])
+        tf2 = TransferFunction(
+                [[[1., 2.]], [[1., 2., 3.]]],
+                [[[1., -2.]], [[1., -3.]]])
+
+        common_den_ref = np.array([[1., -5., 6.]])
+
+        np.testing.assert_raises(ValueError, tf1._common_den)
+        np.testing.assert_raises(ValueError, tf2._common_den)
+
+        _, den1, _ = tf1._common_den(allow_nonproper=True)
+        np.testing.assert_array_almost_equal(den1, common_den_ref)
+        _, den2, _ = tf2._common_den(allow_nonproper=True)
+        np.testing.assert_array_almost_equal(den2, common_den_ref)
+
     @unittest.skipIf(not slycot_check(), "slycot not installed")
     def test_pole_mimo(self):
         """Test for correct MIMO poles."""
@@ -557,6 +577,14 @@ def test_pole_mimo(self):
 
         np.testing.assert_array_almost_equal(p, [-2., -2., -7., -3., -2.])
 
+        # non proper transfer function
+        sys2 = TransferFunction(
+            [[[1., 2., 3., 4.], [1.]], [[1.], [1.]]],
+            [[[1., 2.], [1., 3.]], [[1., 4., 4.], [1., 9., 14.]]])
+        p2 = sys2.pole()
+
+        np.testing.assert_array_almost_equal(p2, [-2., -2., -7., -3., -2.])
+
     def test_double_cancelling_poles_siso(self):
 
         H = TransferFunction([1, 1], [1, 2, 1])
diff --git a/control/xferfcn.py b/control/xferfcn.py
@@ -679,7 +679,7 @@ def freqresp(self, omega):
 
     def pole(self):
         """Compute the poles of a transfer function."""
-        num, den, denorder = self._common_den()
+        _, den, denorder = self._common_den(allow_nonproper=True)
         rts = []
         for d, o in zip(den, denorder):
             rts.extend(roots(d[:o + 1]))
@@ -797,7 +797,7 @@ def returnScipySignalLTI(self):
 
         return out
 
-    def _common_den(self, imag_tol=None):
+    def _common_den(self, imag_tol=None, allow_nonproper=False):
         """
         Compute MIMO common denominators; return them and adjusted numerators.
 
@@ -813,6 +813,9 @@ def _common_den(self, imag_tol=None):
             Threshold for the imaginary part of a root to use in detecting
             complex poles
 
+        allow_nonproper : boolean
+            Do not enforce proper transfer functions
+
         Returns
         -------
         num: array
@@ -822,6 +825,8 @@ def _common_den(self, imag_tol=None):
             gives the numerator coefficient array for the ith output and jth
             input; padded for use in td04ad ('C' option); matches the
             denorder order; highest coefficient starts on the left.
+            If allow_nonproper=True and the order of a numerator exceeds the
+            order of the common denominator, num will be returned as None
 
         den: array
             sys.inputs by kd
@@ -906,6 +911,8 @@ def _common_den(self, imag_tol=None):
                     dtype=float)
         denorder = zeros((self.inputs,), dtype=int)
 
+        havenonproper = False
+
         for j in range(self.inputs):
             if not len(poles[j]):
                 # no poles matching this input; only one or more gains
@@ -930,11 +937,28 @@ def _common_den(self, imag_tol=None):
                         nwzeros.append(poles[j][ip])
 
                     numpoly = poleset[i][j][2] * np.atleast_1d(poly(nwzeros))
+
+                    # td04ad expects a proper transfer function. If the
+                    # numerater has a higher order than the denominator, the
+                    # padding will fail
+                    if len(numpoly) > maxindex + 1:
+                        if allow_nonproper:
+                            havenonproper = True
+                            break
+                        raise ValueError(
+                            self.__str__() +
+                            "is not a proper transfer function. "
+                            "The degree of the numerators must not exceed "
+                            "the degree of the denominators.")
+
                     # numerator polynomial should be padded on left and right
                     #   ending at maxindex to line up with what td04ad expects.
                     num[i, j, maxindex+1-len(numpoly):maxindex+1] = numpoly
                     # print(num[i, j])
 
+        if havenonproper:
+            num = None
+
         return num, den, denorder
 
     def sample(self, Ts, method='zoh', alpha=None):
