Merge pull request #100 "hinfsyn and h2syn fixes"

slivingston · slivingston · commit d4e015598cc8 · 2016-12-26T05:24:26.000-05:00
#100 Changes are from branch `rory-robust-test` of https://github.com/roryyorke/python-control.git
diff --git a/control/robust.py b/control/robust.py
@@ -90,8 +90,8 @@ def h2syn(P,nmeas,ncon):
 
     n = np.size(P.A,0)
     m = np.size(P.B,1)
-    np = np.size(P.C,0)
-    out = sb10hd(n,m,np,ncon,nmeas,P.A,P.B,P.C,P.D)
+    np_ = np.size(P.C,0)
+    out = sb10hd(n,m,np_,ncon,nmeas,P.A,P.B,P.C,P.D)
     Ak = out[0]
     Bk = out[1]
     Ck = out[2]
@@ -115,7 +115,12 @@ def hinfsyn(P,nmeas,ncon):
     K: controller to stabilize P (State-space sys)
     CL: closed loop system (State-space sys)
     gam: infinity norm of closed loop system
-    info: info returned from siycot routine
+    rcond: 4-vector, reciprocal condition estimates of:
+           1: control transformation matrix
+           2: measurement transformation matrix
+           3: X-Ricatti equation
+           4: Y-Ricatti equation
+      TODO: document significance of rcond
 
     Raises
     ------
@@ -128,7 +133,7 @@ def hinfsyn(P,nmeas,ncon):
 
     Examples
     --------
-    >>> K, CL, gam, info = hinfsyn(P,nmeas,ncon)
+    >>> K, CL, gam, rcond = hinfsyn(P,nmeas,ncon)
 
     """
 
@@ -147,12 +152,11 @@ def hinfsyn(P,nmeas,ncon):
     except ImportError:
         raise ControlSlycot("can't find slycot subroutine sb10ad")
 
-    job = 3
     n = np.size(P.A,0)
     m = np.size(P.B,1)
-    np = np.size(P.C,0)
+    np_ = np.size(P.C,0)
     gamma = 1.e100
-    out = sb10ad(job,n,m,np,ncon,nmeas,gamma,P.A,P.B,P.C,P.D)
+    out = sb10ad(n,m,np_,ncon,nmeas,gamma,P.A,P.B,P.C,P.D)
     gam = out[0]
     Ak = out[1]
     Bk = out[2]
@@ -163,10 +167,8 @@ def hinfsyn(P,nmeas,ncon):
     Cc = out[7]
     Dc = out[8]
     rcond = out[9]
-    info = out[10]
 
     K = StateSpace(Ak, Bk, Ck, Dk)
     CL = StateSpace(Ac, Bc, Cc, Dc)
 
-    return K, CL, gam, info
-
+    return K, CL, gam, rcond
diff --git a/control/tests/robust_test.py b/control/tests/robust_test.py
@@ -0,0 +1,46 @@
+import unittest
+import numpy as np
+import control
+import control.robust
+from control.exception import slycot_check
+
+
+class TestHinf(unittest.TestCase):
+    @unittest.skipIf(not slycot_check(), "slycot not installed")
+    def testHinfsyn(self):
+        "Test hinfsyn"
+        p = control.ss(-1, [1, 1], [[1], [1]], [[0, 1], [1, 0]])
+        k, cl, gam, rcond = control.robust.hinfsyn(p, 1, 1)
+        # from Octave, which also uses SB10AD:
+        #   a= -1; b1= 1; b2= 1; c1= 1; c2= 1; d11= 0; d12= 1; d21= 1; d22= 0;
+        #   g = ss(a,[b1,b2],[c1;c2],[d11,d12;d21,d22]);
+        #   [k,cl] = hinfsyn(g,1,1);
+        np.testing.assert_array_almost_equal(k.A, [[-3]])
+        np.testing.assert_array_almost_equal(k.B, [[1]])
+        np.testing.assert_array_almost_equal(k.C, [[-1]])
+        np.testing.assert_array_almost_equal(k.D, [[0]])
+        np.testing.assert_array_almost_equal(cl.A, [[-1, -1], [1, -3]])
+        np.testing.assert_array_almost_equal(cl.B, [[1], [1]])
+        np.testing.assert_array_almost_equal(cl.C, [[1, -1]])
+        np.testing.assert_array_almost_equal(cl.D, [[0]])
+
+    # TODO: add more interesting examples
+
+class TestH2(unittest.TestCase):
+    @unittest.skipIf(not slycot_check(), "slycot not installed")
+    def testH2syn(self):
+        "Test h2syn"
+        p = control.ss(-1, [1, 1], [[1], [1]], [[0, 1], [1, 0]])
+        k = control.robust.h2syn(p, 1, 1)
+        # from Octave, which also uses SB10HD for H-2 synthesis:
+        #   a= -1; b1= 1; b2= 1; c1= 1; c2= 1; d11= 0; d12= 1; d21= 1; d22= 0;
+        #   g = ss(a,[b1,b2],[c1;c2],[d11,d12;d21,d22]);
+        #   k = h2syn(g,1,1);
+        # the solution is the same as for the hinfsyn test
+        np.testing.assert_array_almost_equal(k.A, [[-3]])
+        np.testing.assert_array_almost_equal(k.B, [[1]])
+        np.testing.assert_array_almost_equal(k.C, [[-1]])
+        np.testing.assert_array_almost_equal(k.D, [[0]])
+
+if __name__ == "__main__":
+    unittest.main()
