python-control · murrayrm · May 31, 2019 · Jul 14, 2018 · May 14, 2019 · May 18, 2019
diff --git a/control/bdalg.py b/control/bdalg.py
@@ -243,7 +243,7 @@ def feedback(sys1, sys2=1, sign=-1):
         elif isinstance(sys2, ss.StateSpace):
             sys1 = ss._convertToStateSpace(sys1)
         elif isinstance(sys2, frd.FRD):
-            sys1 = ss._convertToFRD(sys1)
+            sys1 = frd._convertToFRD(sys1, sys2.omega)
         else: # sys2 is a scalar.
             sys1 = tf._convert_to_transfer_function(sys1)
             sys2 = tf._convert_to_transfer_function(sys2)

diff --git a/control/canonical.py b/control/canonical.py
@@ -88,7 +88,9 @@ def reachable_form(xsys):
     # Transformation from one form to another
     Tzx = solve(Wrx.T, Wrz.T).T  # matrix right division, Tzx = Wrz * inv(Wrx)
 
-    if matrix_rank(Tzx) != xsys.states:
+    # Check to make sure inversion was OK.  Note that since we are inverting
+    # Wrx and we already checked its rank, this exception should never occur
+    if matrix_rank(Tzx) != xsys.states:         # pragma: no cover
         raise ValueError("Transformation matrix singular to working precision.")
 
     # Finally, compute the output matrix

diff --git a/control/config.py b/control/config.py
@@ -14,6 +14,16 @@
 bode_number_of_samples = None   # Bode plot number of samples
 bode_feature_periphery_decade = 1.0  # Bode plot feature periphery in decades
 
+
+def reset_defaults():
+    """Reset configuration values to their default values."""
+    global bode_dB; bode_dB = False
+    global bode_deg; bode_deg = True
+    global bode_Hz; bode_Hz = False
+    global bode_number_of_samples; bode_number_of_samples = None
+    global bode_feature_periphery_decade; bode_feature_periphery_decade = 1.0
+
+
 # Set defaults to match MATLAB
 def use_matlab_defaults():
     """
@@ -27,6 +37,7 @@ def use_matlab_defaults():
     global bode_deg; bode_deg = True
     global bode_Hz; bode_Hz = True
 
+
 # Set defaults to match FBS (Astrom and Murray)
 def use_fbs_defaults():
     """
@@ -39,4 +50,5 @@ def use_fbs_defaults():
     # Bode plot defaults
     global bode_dB; bode_dB = False
     global bode_deg; bode_deg = True
-    global bode_Hz; bode_Hz = True
+    global bode_Hz; bode_Hz = False
+
diff --git a/control/dtime.py b/control/dtime.py
@@ -112,6 +112,9 @@ def c2d(sysc, Ts, method='zoh'):
     '''
     #  Call the sample_system() function to do the work
     sysd = sample_system(sysc, Ts, method)
+
+    # TODO: is this check needed?  If sysc is  StateSpace, sysd is too?
     if isinstance(sysc, StateSpace) and not isinstance(sysd, StateSpace):
-        return _convertToStateSpace(sysd)
+        return _convertToStateSpace(sysd)       # pragma: no cover
+
     return sysd
diff --git a/control/exception.py b/control/exception.py
@@ -39,24 +39,24 @@
 #
 # $Id$
 
-class ControlSlycot(Exception):
+class ControlSlycot(ImportError):
     """Exception for Slycot import.  Used when we can't import a function
     from the slycot package"""
     pass
 
-class ControlDimension(Exception):
+class ControlDimension(ValueError):
     """Raised when dimensions of system objects are not correct"""
     pass
 
-class ControlArgument(Exception):
+class ControlArgument(TypeError):
     """Raised when arguments to a function are not correct"""
     pass
 
-class ControlMIMONotImplemented(Exception):
+class ControlMIMONotImplemented(NotImplementedError):
     """Function is not currently implemented for MIMO systems"""
     pass
 
-class ControlNotImplemented(Exception):
+class ControlNotImplemented(NotImplementedError):
     """Functionality is not yet implemented"""
     pass
 

diff --git a/control/frdata.py b/control/frdata.py
@@ -49,6 +49,7 @@
 """
 
 # External function declarations
+from warnings import warn
 import numpy as np
 from numpy import angle, array, empty, ones, \
     real, imag, matrix, absolute, eye, linalg, where, dot
@@ -187,9 +188,10 @@ def __add__(self, other):
 
         if isinstance(other, FRD):
             # verify that the frequencies match
-            if (other.omega != self.omega).any():
-                print("Warning: frequency points do not match; expect"
-                      " truncation and interpolation")
+            if len(other.omega) != len(self.omega) or \
+               (other.omega != self.omega).any():
+                warn("Frequency points do not match; expect"
+                      " truncation and interpolation.")
 
         # Convert the second argument to a frequency response function.
         # or re-base the frd to the current omega (if needed)
@@ -340,8 +342,9 @@ def evalfr(self, omega):
         intermediate values.
 
         """
-        warn("FRD.evalfr(omega) will be deprecated in a future release of python-control; use sys.eval(omega) instead", 
-             PendingDeprecationWarning)
+        warn("FRD.evalfr(omega) will be deprecated in a future release "
+             "of python-control; use sys.eval(omega) instead", 
+             PendingDeprecationWarning)         # pragma: no coverage
         return self._evalfr(omega)
 
     # Define the `eval` function to evaluate an FRD at a given (real)
@@ -352,7 +355,7 @@ def evalfr(self, omega):
     def eval(self, omega):
         """Evaluate a transfer function at a single angular frequency.
 
-        self._evalfr(omega) returns the value of the frequency response
+        self.evalfr(omega) returns the value of the frequency response
         at frequency omega.
 
         Note that a "normal" FRD only returns values for which there is an
@@ -462,7 +465,8 @@ def _convertToFRD(sys, omega, inputs=1, outputs=1):
 
     if isinstance(sys, FRD):
         omega.sort()
-        if (abs(omega - sys.omega) < FRD.epsw).all():
+        if len(omega) == len(sys.omega) and \
+           (abs(omega - sys.omega) < FRD.epsw).all():
             # frequencies match, and system was already frd; simply use
             return sys
 

diff --git a/control/freqplot.py b/control/freqplot.py
@@ -603,18 +603,24 @@ def default_frequency_range(syslist, Hz=None, number_of_samples=None,
     ----------
     syslist : list of LTI
         List of linear input/output systems (single system is OK)
-    Hz: boolean
+
+    Hz : bool
         If True, the limits (first and last value) of the frequencies
         are set to full decades in Hz so it fits plotting with logarithmic
         scale in Hz otherwise in rad/s. Omega is always returned in rad/sec.
-    number_of_samples: int
-        Number of samples to generate
-    feature_periphery_decade: float
+
+    number_of_samples : int, optional
+        Number of samples to generate.  Defaults to ``numpy.logspace`` default
+        value.
+
+    feature_periphery_decade : float, optional
         Defines how many decades shall be included in the frequency range on
         both sides of features (poles, zeros).
-        Example: If there is a feature, e.g. a pole, at 1Hz and feature_periphery_decade=1.
-        then the range of frequencies shall span 0.1 .. 10 Hz.
-        The default value is read from config.bode_feature_periphery_decade.
+
+        Example: If there is a feature, e.g. a pole, at 1 Hz and
+        feature_periphery_decade=1., then the range of frequencies shall span
+        0.1 .. 10 Hz.  The default value is read from
+        ``config.bode_feature_periphery_decade``.
 
     Returns
     -------
@@ -626,6 +632,7 @@ def default_frequency_range(syslist, Hz=None, number_of_samples=None,
     >>> from matlab import ss
     >>> sys = ss("1. -2; 3. -4", "5.; 7", "6. 8", "9.")
     >>> omega = default_frequency_range(sys)
+
     """
     # This code looks at the poles and zeros of all of the systems that
     # we are plotting and sets the frequency range to be one decade above

diff --git a/control/tests/bdalg_test.py b/control/tests/bdalg_test.py
@@ -253,9 +253,27 @@ def assert_equal(x,y):
         assert_equal(ref.C, tst.C)
         assert_equal(ref.D, tst.D)
 
+    def test_feedback_args(self):
+        # Added 25 May 2019 to cover missing exception handling in feedback()
+        # If first argument is not LTI or convertable, generate an exception
+        args = ([1], self.sys2)
+        self.assertRaises(TypeError, ctrl.feedback, *args)
+
+        # If second argument is not LTI or convertable, generate an exception
+        args = (self.sys1, np.array([1]))
+        self.assertRaises(TypeError, ctrl.feedback, *args)
+
+        # Convert first argument to FRD, if needed
+        h = TransferFunction([1], [1, 2, 2])
+        omega = np.logspace(-1, 2, 10)
+        frd = ctrl.FRD(h, omega)
+        sys = ctrl.feedback(1, frd)
+        self.assertTrue(isinstance(sys, ctrl.FRD))
+
 
 def suite():
     return unittest.TestLoader().loadTestsFromTestCase(TestFeedback)
 
+
 if __name__ == "__main__":
     unittest.main()
diff --git a/control/tests/canonical_test.py b/control/tests/canonical_test.py
@@ -2,9 +2,10 @@
 
 import unittest
 import numpy as np
-from control import ss
-from control.canonical import canonical_form
-
+from control import ss, tf
+from control.canonical import canonical_form, reachable_form, \
+    observable_form, modal_form
+from control.exception import ControlNotImplemented
 
 class TestCanonical(unittest.TestCase):
     """Tests for the canonical forms class"""
@@ -40,6 +41,11 @@ def test_reachable_form(self):
         np.testing.assert_array_almost_equal(sys_check.D, D_true)
         np.testing.assert_array_almost_equal(T_check, T_true)
 
+        # Reachable form only supports SISO
+        sys = tf([[ [1], [1] ]], [[ [1, 2, 1], [1, 2, 1] ]])
+        np.testing.assert_raises(ControlNotImplemented, reachable_form, sys)
+
+
     def test_unreachable_system(self):
         """Test reachable canonical form with an unreachable system"""
 
@@ -76,12 +82,84 @@ def test_modal_form(self):
         sys_check, T_check = canonical_form(ss(A, B, C, D), "modal")
 
         # Check against the true values
-        #TODO: Test in respect to ambiguous transformation (system characteristics?)
+        # TODO: Test in respect to ambiguous transformation (system characteristics?)
         np.testing.assert_array_almost_equal(sys_check.A, A_true)
         #np.testing.assert_array_almost_equal(sys_check.B, B_true)
         #np.testing.assert_array_almost_equal(sys_check.C, C_true)
         np.testing.assert_array_almost_equal(sys_check.D, D_true)
         #np.testing.assert_array_almost_equal(T_check, T_true)
+
+        # Check conversion when there are complex eigenvalues
+        A_true = np.array([[-1,  1,  0,  0],
+                           [-1, -1,  0,  0],
+                           [ 0,  0, -2,  0],
+                           [ 0,  0,  0, -3]])
+        B_true = np.array([[0], [1], [0], [1]])
+        C_true = np.array([[1, 0, 0, 1]])
+        D_true = np.array([[0]])
+
+        A = np.linalg.solve(T_true, A_true) * T_true
+        B = np.linalg.solve(T_true, B_true)
+        C = C_true * T_true
+        D = D_true
+
+        # Create state space system and convert to modal canonical form
+        sys_check, T_check = canonical_form(ss(A, B, C, D), 'modal')
+
+        # Check A and D matrix, which are uniquely defined
+        np.testing.assert_array_almost_equal(sys_check.A, A_true)
+        np.testing.assert_array_almost_equal(sys_check.D, D_true)
+
+        # B matrix should be all ones (or zero if not controllable)
+        # TODO: need to update modal_form() to implement this
+        if np.allclose(T_check, T_true):
+            np.testing.assert_array_almost_equal(sys_check.B, B_true)
+            np.testing.assert_array_almost_equal(sys_check.C, C_true)
+
+        # Make sure Hankel coefficients are OK
+        from numpy.linalg import matrix_power
+        for i in range(A.shape[0]):
+            np.testing.assert_almost_equal(
+                np.dot(np.dot(C_true, matrix_power(A_true, i)), B_true),
+                np.dot(np.dot(C, matrix_power(A, i)), B))
+
+        # Reorder rows to get complete coverage (real eigenvalue cxrtvfirst)
+        A_true = np.array([[-1,  0, 0,  0],
+                           [ 0, -2,  1,  0],
+                           [ 0, -1, -2,  0],
+                           [ 0,  0,  0, -3]])
+        B_true = np.array([[0], [0], [1], [1]])
+        C_true = np.array([[0, 1, 0, 1]])
+        D_true = np.array([[0]])
+
+        A = np.linalg.solve(T_true, A_true) * T_true
+        B = np.linalg.solve(T_true, B_true)
+        C = C_true * T_true
+        D = D_true
+
+        # Create state space system and convert to modal canonical form
+        sys_check, T_check = canonical_form(ss(A, B, C, D), 'modal')
+
+        # Check A and D matrix, which are uniquely defined
+        np.testing.assert_array_almost_equal(sys_check.A, A_true)
+        np.testing.assert_array_almost_equal(sys_check.D, D_true)
+
+        # B matrix should be all ones (or zero if not controllable)
+        # TODO: need to update modal_form() to implement this
+        if np.allclose(T_check, T_true):
+            np.testing.assert_array_almost_equal(sys_check.B, B_true)
+            np.testing.assert_array_almost_equal(sys_check.C, C_true)
+
+        # Make sure Hankel coefficients are OK
+        from numpy.linalg import matrix_power
+        for i in range(A.shape[0]):
+            np.testing.assert_almost_equal(
+                np.dot(np.dot(C_true, matrix_power(A_true, i)), B_true),
+                np.dot(np.dot(C, matrix_power(A, i)), B))
+
+        # Modal form only supports SISO
+        sys = tf([[ [1], [1] ]], [[ [1, 2, 1], [1, 2, 1] ]])
+        np.testing.assert_raises(ControlNotImplemented, modal_form, sys)
 
     def test_observable_form(self):
         """Test the observable canonical form"""
@@ -114,6 +192,11 @@ def test_observable_form(self):
         np.testing.assert_array_almost_equal(sys_check.D, D_true)
         np.testing.assert_array_almost_equal(T_check, T_true)
 
+        # Observable form only supports SISO
+        sys = tf([[ [1], [1] ]], [[ [1, 2, 1], [1, 2, 1] ]])
+        np.testing.assert_raises(ControlNotImplemented, observable_form, sys)
+
+
     def test_unobservable_system(self):
         """Test observable canonical form with an unobservable system"""
 
@@ -126,3 +209,18 @@ def test_unobservable_system(self):
 
         # Check if an exception is raised
         np.testing.assert_raises(ValueError, canonical_form, sys, "observable")
+
+    def test_arguments(self):
+        # Additional unit tests added on 25 May 2019 to increase coverage
+
+        # Unknown canonical forms should generate exception
+        sys = tf([1], [1, 2, 1])
+        np.testing.assert_raises(
+            ControlNotImplemented, canonical_form, sys, 'unknown')
+
+def suite():
+    return unittest.TestLoader().loadTestsFromTestCase(TestFeedback)
+
+
+if __name__ == "__main__":
+    unittest.main()