fix up rounding errors and matrix/array

murrayrm · murrayrm · commit 1d1c2890966f · 2023-03-24T23:28:10.000-07:00
diff --git a/control/canonical.py b/control/canonical.py
@@ -17,6 +17,7 @@
 __all__ = ['canonical_form', 'reachable_form', 'observable_form', 'modal_form',
            'similarity_transform', 'bdschur']
 
+
 def canonical_form(xsys, form='reachable'):
     """Convert a system into canonical form
 
@@ -41,16 +42,16 @@ def canonical_form(xsys, form='reachable'):
     --------
     >>> Gs = ct.tf2ss([1], [1, 3, 2])
     >>> Gc, T = ct.canonical_form(Gs)  # default reachable
-    >>> Gc.B
+    >>> Gc.B.round()
     array([[1.],
            [0.]])
 
     >>> Gc, T = ct.canonical_form(Gs, 'observable')
-    >>> Gc.C
+    >>> Gc.C.round()
     array([[1., 0.]])
 
     >>> Gc, T = ct.canonical_form(Gs, 'modal')
-    >>> Gc.A
+    >>> Gc.A.round()                                            # doctest: +SKIP
     array([[-2.,  0.],
            [ 0., -1.]])
 
@@ -88,9 +89,9 @@ def reachable_form(xsys):
     --------
     >>> Gs = ct.tf2ss([1], [1, 3, 2])
     >>> Gc, T = ct.reachable_form(Gs)  # default reachable
-    >>> Gc.B                                                    # doctest: +SKIP
-    matrix([[1.],
-            [0.]])
+    >>> Gc.B.round()
+    array([[1.],
+           [0.]])
 
     """
     # Check to make sure we have a SISO system
@@ -151,8 +152,8 @@ def observable_form(xsys):
     --------
     >>> Gs = ct.tf2ss([1], [1, 3, 2])
     >>> Gc, T = ct.observable_form(Gs)
-    >>> Gc.C                                                    # doctest: +SKIP
-    matrix([[1., 0.]])
+    >>> Gc.C.round()
+    array([[1., 0.]])
 
     """
     # Check to make sure we have a SISO system
@@ -216,15 +217,15 @@ def similarity_transform(xsys, T, timescale=1, inverse=False):
     Examples
     --------
     >>> Gs = ct.tf2ss([1], [1, 3, 2])
-    >>> Gs.A                                                    # doctest: +SKIP
-    matrix([[-3., -2.],
-            [ 1.,  0.]])
+    >>> Gs.A.round()
+    array([[-3., -2.],
+           [ 1.,  0.]])
 
-    >>> T = np.array([[0,1],[1,0]])
+    >>> T = np.array([[0, 1], [1, 0]])
     >>> Gt = ct.similarity_transform(Gs, T)
-    >>> Gt.A                                                    # doctest: +SKIP
-    matrix([[ 0.,  1.],
-            [-2., -3.]])
+    >>> Gt.A.round()
+    array([[ 0.,  1.],
+           [-2., -3.]])
 
     """
     # Create a new system, starting with a copy of the old one
@@ -293,7 +294,8 @@ def _bdschur_condmax_search(aschur, tschur, condmax):
 
     Iterates mb03rd with different pmax values until:
       - result is non-defective;
-      - or condition number of similarity transform is unchanging despite large pmax;
+      - or condition number of similarity transform is unchanging
+        despite large pmax;
       - or condition number of similarity transform is close to condmax.
 
     Parameters
@@ -332,22 +334,25 @@ def _bdschur_condmax_search(aschur, tschur, condmax):
 
     # get lower bound; try condmax ** 0.5 first
     pmaxlower = condmax ** 0.5
-    amodal, tmodal, blksizes, eigvals = mb03rd(aschur.shape[0], aschur, tschur, pmax=pmaxlower)
+    amodal, tmodal, blksizes, eigvals = mb03rd(
+        aschur.shape[0], aschur, tschur, pmax=pmaxlower)
     if np.linalg.cond(tmodal) <= condmax:
         reslower = amodal, tmodal, blksizes, eigvals
     else:
         pmaxlower = 1.0
-        amodal, tmodal, blksizes, eigvals = mb03rd(aschur.shape[0], aschur, tschur, pmax=pmaxlower)
+        amodal, tmodal, blksizes, eigvals = mb03rd(
+            aschur.shape[0], aschur, tschur, pmax=pmaxlower)
         cond = np.linalg.cond(tmodal)
         if cond > condmax:
-            msg = 'minimum cond={} > condmax={}; try increasing condmax'.format(cond, condmax)
+            msg = f"minimum {cond=} > {condmax=}; try increasing condmax"
             raise RuntimeError(msg)
 
     pmax = pmaxlower
 
     # phase 1: search for upper bound on pmax
     for i in range(50):
-        amodal, tmodal, blksizes, eigvals = mb03rd(aschur.shape[0], aschur, tschur, pmax=pmax)
+        amodal, tmodal, blksizes, eigvals = mb03rd(
+            aschur.shape[0], aschur, tschur, pmax=pmax)
         cond = np.linalg.cond(tmodal)
         if cond < condmax:
             pmaxlower = pmax
@@ -368,7 +373,8 @@ def _bdschur_condmax_search(aschur, tschur, condmax):
     # phase 2: bisection search
     for i in range(50):
         pmax = (pmaxlower * pmaxupper) ** 0.5
-        amodal, tmodal, blksizes, eigvals = mb03rd(aschur.shape[0], aschur, tschur, pmax=pmax)
+        amodal, tmodal, blksizes, eigvals = mb03rd(
+            aschur.shape[0], aschur, tschur, pmax=pmax)
         cond = np.linalg.cond(tmodal)
 
         if cond < condmax:
@@ -440,7 +446,8 @@ def bdschur(a, condmax=None, sort=None):
         return a.copy(), np.eye(a.shape[1], a.shape[0]), np.array([])
 
     aschur, tschur = schur(a)
-    amodal, tmodal, blksizes, eigvals = _bdschur_condmax_search(aschur, tschur, condmax)
+    amodal, tmodal, blksizes, eigvals = _bdschur_condmax_search(
+        aschur, tschur, condmax)
 
     if sort in ('continuous', 'discrete'):
 
@@ -484,9 +491,11 @@ def modal_form(xsys, condmax=None, sort=False):
     xsys : StateSpace object
         System to be transformed, with state `x`
     condmax : None or float, optional
-        An upper bound on individual transformations.  If None, use `bdschur` default.
+        An upper bound on individual transformations.  If None, use
+        `bdschur` default.
     sort : bool, optional
-        If False (default), Schur blocks will not be sorted.  See `bdschur` for sort order.
+        If False (default), Schur blocks will not be sorted.  See `bdschur`
+        for sort order.
 
     Returns
     -------
@@ -499,9 +508,9 @@ def modal_form(xsys, condmax=None, sort=False):
     --------
     >>> Gs = ct.tf2ss([1], [1, 3, 2])
     >>> Gc, T = ct.modal_form(Gs)  # default reachable
-    >>> Gc.A                                                    # doctest: +SKIP
-    matrix([[-2.,  0.],
-            [ 0., -1.]])
+    >>> Gc.A.round()                                            # doctest: +SKIP
+    array([[-2.,  0.],
+           [ 0., -1.]])
 
     """
 
diff --git a/control/descfcn.py b/control/descfcn.py
@@ -249,7 +249,7 @@ def describing_function_plot(
     >>> H_simple = ct.tf([8], [1, 2, 2, 1])
     >>> F_saturation = ct.saturation_nonlinearity(1)
     >>> amp = np.linspace(1, 4, 10)
-    >>> ct.describing_function_plot(H_simple, F_saturation, amp)
+    >>> ct.describing_function_plot(H_simple, F_saturation, amp)  # doctest: +SKIP
     [(3.343844998258643, 1.4142293090899216)]
 
     """
