add unit tests for exceptions/warnings + cleanup

murrayrm · murrayrm · commit e180716898eb · 2021-01-28T12:23:48.000-08:00
diff --git a/control/descfcn.py b/control/descfcn.py
@@ -254,9 +254,15 @@ def describing_function_plot(
             if refine:
                 # Refine the answer to get more accuracy
                 def _cost(x):
+                    # If arguments are invalid, return a "large" value
+                    # Note: imposing bounds messed up the optimization (?)
+                    if x[0] < 0 or x[1] < 0:
+                        return 1
                     return abs(1 + H(1j * x[1]) *
                                describing_function(F, x[0]))**2
-                res = scipy.optimize.minimize(_cost, [a_guess, omega_guess])
+                res = scipy.optimize.minimize(
+                    _cost, [a_guess, omega_guess])
+                # bounds=[(A[i], A[i+1]), (H_omega[j], H_omega[j+1])])
 
                 if not res.success:
                     warn("not able to refine result; returning estimate")
@@ -322,6 +328,9 @@ class saturation_nonlinearity(DescribingFunctionNonlinearity):
 
     """
     def __init__(self, ub=1, lb=None):
+        # Create the describing function nonlinearity object
+        super(saturation_nonlinearity, self).__init__()
+
         # Process arguments
         if lb == None:
             # Only received one argument; assume symmetric around zero
@@ -341,6 +350,10 @@ def isstatic(self):
         return True
 
     def describing_function(self, A):
+        # Check to make sure the amplitude is positive
+        if A < 0:
+            raise ValueError("cannot evaluate describing function for A < 0")
+
         if self.lb <= A and A <= self.ub:
             return 1.
         else:
@@ -368,6 +381,9 @@ class relay_hysteresis_nonlinearity(DescribingFunctionNonlinearity):
 
     """
     def __init__(self, b, c):
+        # Create the describing function nonlinearity object
+        super(relay_hysteresis_nonlinearity, self).__init__()
+
         # Initialize the state to bottom branch
         self.branch = -1        # lower branch
         self.b = b              # relay output value
@@ -389,16 +405,16 @@ def __call__(self, x):
     def isstatic(self):
         return False
 
-    def describing_function(self, a):
-        def f(x):
-            return math.copysign(1, x) if abs(x) > 1 else \
-                (math.asin(x) + x * math.sqrt(1 - x**2)) * 2 / math.pi
+    def describing_function(self, A):
+        # Check to make sure the amplitude is positive
+        if A < 0:
+            raise ValueError("cannot evaluate describing function for A < 0")
 
-        if a < self.c:
+        if A < self.c:
             return np.nan
 
-        df_real = 4 * self.b * math.sqrt(1 - (self.c/a)**2) / (a * math.pi)
-        df_imag = -4 * self.b * self.c / (math.pi * a**2)
+        df_real = 4 * self.b * math.sqrt(1 - (self.c/A)**2) / (A * math.pi)
+        df_imag = -4 * self.b * self.c / (math.pi * A**2)
         return df_real + 1j * df_imag
 
 
@@ -421,6 +437,9 @@ class backlash_nonlinearity(DescribingFunctionNonlinearity):
     """
 
     def __init__(self, b):
+        # Create the describing function nonlinearity object
+        super(backlash_nonlinearity, self).__init__()
+
         self.b = b              # backlash distance
         self.center = 0         # current center position
 
@@ -442,6 +461,10 @@ def isstatic(self):
         return False
 
     def describing_function(self, A):
+        # Check to make sure the amplitude is positive
+        if A < 0:
+            raise ValueError("cannot evaluate describing function for A < 0")
+
         if A <= self.b/2:
             return 0
 
diff --git a/control/iosys.py b/control/iosys.py
@@ -821,7 +821,7 @@ def __call__(sys, u, squeeze=None, params=None):
         # If we received any parameters, update them before calling _out()
         if params is not None:
             sys._update_params(params)
-            
+
         # Evaluate the function on the argument
         out = sys._out(0, np.array((0,)), np.asarray(u))
         _, out = _process_time_response(sys, [], out, [], squeeze=squeeze)
@@ -2004,7 +2004,7 @@ def interconnect(syslist, connections=[], inplist=[], outlist=[],
     a warning is generated a copy of the system is created with the
     name of the new system determined by adding the prefix and suffix
     strings in config.defaults['iosys.linearized_system_name_prefix']
-    and config.defaults['iosys.linearized_system_name_suffix'], with the 
+    and config.defaults['iosys.linearized_system_name_suffix'], with the
     default being to add the suffix '$copy'$ to the system name.
 
     It is possible to replace lists in most of arguments with tuples instead,
diff --git a/control/tests/descfcn_test.py b/control/tests/descfcn_test.py
@@ -12,8 +12,12 @@
 import numpy as np
 import control as ct
 import math
+from control.descfcn import saturation_nonlinearity, backlash_nonlinearity, \
+    relay_hysteresis_nonlinearity
+
 
-class saturation():
+# Static function via a class
+class saturation_class():
     # Static nonlinear saturation function
     def __call__(self, x, lb=-1, ub=1):
         return np.maximum(lb, np.minimum(x, ub))
@@ -27,10 +31,15 @@ def describing_function(self, a):
             return 2/math.pi * (math.asin(b) + b * math.sqrt(1 - b**2))
 
 
+# Static function without a class
+def saturation(x):
+    return np.maximum(-1, np.minimum(x, 1))
+
+
 # Static nonlinear system implementing saturation
 @pytest.fixture
 def satsys():
-    satfcn = saturation()
+    satfcn = saturation_class()
     def _satfcn(t, x, u, params):
         return satfcn(u)
     return ct.NonlinearIOSystem(None, outfcn=_satfcn, input=1, output=1)
@@ -65,16 +74,16 @@ def _misofcn(t, x, u, params={}):
     np.testing.assert_array_equal(miso_sys([0, 0]), [0])
     np.testing.assert_array_equal(miso_sys([0, 0]), [0])
     np.testing.assert_array_equal(miso_sys([0, 0], squeeze=True), [0])
-    
+
 
 # Test saturation describing function in multiple ways
 def test_saturation_describing_function(satsys):
-    satfcn = saturation()
-    
+    satfcn = saturation_class()
+
     # Store the analytic describing function for comparison
     amprange = np.linspace(0, 10, 100)
     df_anal = [satfcn.describing_function(a) for a in amprange]
-    
+
     # Compute describing function for a static function
     df_fcn = [ct.describing_function(satfcn, a) for a in amprange]
     np.testing.assert_almost_equal(df_fcn, df_anal, decimal=3)
@@ -87,8 +96,9 @@ def test_saturation_describing_function(satsys):
     df_arr = ct.describing_function(satsys, amprange)
     np.testing.assert_almost_equal(df_arr, df_anal, decimal=3)
 
-from control.descfcn import saturation_nonlinearity, backlash_nonlinearity, \
-    relay_hysteresis_nonlinearity
+    # Evaluate static function at a negative amplitude
+    with pytest.raises(ValueError, match="cannot evaluate"):
+        ct.describing_function(saturation, -1)
 
 
 @pytest.mark.parametrize("fcn, amin, amax", [
@@ -100,7 +110,7 @@ def test_describing_function(fcn, amin, amax):
     # Store the analytic describing function for comparison
     amprange = np.linspace(amin, amax, 100)
     df_anal = [fcn.describing_function(a) for a in amprange]
-    
+
     # Compute describing function on an array of values
     df_arr = ct.describing_function(
         fcn, amprange, zero_check=False, try_method=False)
@@ -110,6 +120,11 @@ def test_describing_function(fcn, amin, amax):
     df_meth = ct.describing_function(fcn, amprange, zero_check=False)
     np.testing.assert_almost_equal(df_meth, df_anal, decimal=1)
 
+    # Make sure that evaluation at negative amplitude generates an exception
+    with pytest.raises(ValueError, match="cannot evaluate"):
+        ct.describing_function(fcn, -1)
+
+
 def test_describing_function_plot():
     # Simple linear system with at most 1 intersection
     H_simple = ct.tf([1], [1, 2, 2, 1])
@@ -141,3 +156,29 @@ def test_describing_function_plot():
         np.testing.assert_almost_equal(
             -1/ct.describing_function(F_backlash, a),
             H_multiple(1j*w), decimal=5)
+
+def test_describing_function_exceptions():
+    # Describing function with non-zero bias
+    with pytest.warns(UserWarning, match="asymmetric"):
+        saturation = ct.descfcn.saturation_nonlinearity(lb=-1, ub=2)
+        assert saturation(-3) == -1
+        assert saturation(3) == 2
+
+    # Turn off the bias check
+    bias = ct.describing_function(saturation, 0, zero_check=False)
+
+    # Function should evaluate to zero at zero amplitude
+    f = lambda x: x + 0.5
+    with pytest.raises(ValueError, match="must evaluate to zero"):
+        bias = ct.describing_function(f, 0, zero_check=True)
+
+    # Evaluate at a negative amplitude
+    with pytest.raises(ValueError, match="cannot evaluate"):
+        ct.describing_function(saturation, -1)
+
+    # Describing function with bad label
+    H_simple = ct.tf([8], [1, 2, 2, 1])
+    F_saturation = ct.descfcn.saturation_nonlinearity(1)
+    amp = np.linspace(1, 4, 10)
+    with pytest.raises(ValueError, match="formatting string"):
+        ct.describing_function_plot(H_simple, F_saturation, amp, label=1)
