Skip to content

Commit 60d6c2c

Browse files
author
arnold
committed
Removed dependence on python-control and scipy. Made dependence of lzwork and ldwork on n explicit to address travis failing to build with python 3.3.5.
1 parent d6f3f61 commit 60d6c2c

File tree

3 files changed

+66
-87
lines changed

3 files changed

+66
-87
lines changed

‎slycot/src/transform.pyf‎

Lines changed: 12 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -159,10 +159,10 @@ subroutine tb05ad_ag(baleig,inita,n,m,p,freq,a,lda,b,ldb,c,ldc,rcond,g,ldg,evre,
159159
integer intent(hide),depend(n) :: ldhinv = n
160160
! cache variables
161161
integer intent(hide,cache),dimension(n) :: iwork
162-
double precision intent(hide,cache),dimension(2*n),depend(n) :: dwork
163-
integer optional :: ldwork = 2*n
164-
complex*16 intent(hide,cache),dimension(lzwork),depend(n) :: zwork
165-
integer optional :: lzwork = n*n+2*n
162+
double precision intent(hide,cache),dimension(ldwork),depend(ldwork) :: dwork
163+
integer optional,depend(n) :: ldwork = 2*n
164+
complex*16 intent(hide,cache),dimension(lzwork),depend(lzwork) :: zwork
165+
integer optional,depend(n) :: lzwork = n*n+2*n
166166
integer intent(out) :: info
167167
end subroutine tb05ad_ag
168168
subroutine tb05ad_ng(baleig,inita,n,m,p,freq,a,lda,b,ldb,c,ldc,rcond,g,ldg,evre,evim,hinvb,ldhinv,iwork,dwork,ldwork,zwork,lzwork,info)
@@ -188,10 +188,10 @@ subroutine tb05ad_ng(baleig,inita,n,m,p,freq,a,lda,b,ldb,c,ldc,rcond,g,ldg,evre,
188188
integer intent(hide),depend(n) :: ldhinv = n
189189
! cache variables
190190
integer intent(hide,cache),dimension(n) :: iwork
191-
double precision intent(hide,cache),dimension(2*n),depend(n) :: dwork
192-
integer optional :: ldwork = 2*n
193-
complex*16 intent(hide,cache),dimension(lzwork),depend(n) :: zwork
194-
integer optional :: lzwork = n*n+2*n
191+
double precision intent(hide,cache),dimension(ldwork),depend(ldwork) :: dwork
192+
integer optional,depend(n) :: ldwork = 2*n
193+
complex*16 intent(hide,cache),dimension(lzwork),depend(lzwork) :: zwork
194+
integer optional,depend(n) :: lzwork = n*n+2*n
195195
integer intent(out) :: info
196196
end subroutine tb05ad_ng
197197

@@ -219,10 +219,10 @@ subroutine tb05ad_nh(baleig,inita,n,m,p,freq,a,lda,b,ldb,c,ldc,rcond,g,ldg,evre,
219219
integer intent(hide),depend(n) :: ldhinv = n
220220
! cache variables
221221
integer intent(hide,cache),dimension(n) :: iwork
222-
double precision intent(hide,cache),dimension(2*n),depend(n) :: dwork
223-
integer optional :: ldwork = 2*n
224-
complex*16 intent(hide,cache),dimension(lzwork),depend(n) :: zwork
225-
integer optional :: lzwork = n*n+2*n
222+
double precision intent(hide,cache),dimension(ldwork),depend(ldwork) :: dwork
223+
integer optional,depend(n) :: ldwork = 2*n
224+
complex*16 intent(hide,cache),dimension(lzwork),depend(lzwork) :: zwork
225+
integer optional,depend(n) :: lzwork = n*n+2*n
226226
integer intent(out) :: info
227227
end subroutine tb05ad_h
228228

‎slycot/tests/test.py‎

Lines changed: 42 additions & 66 deletions
Original file line numberDiff line numberDiff line change
@@ -67,129 +67,105 @@ def test_td04ad_static(self):
6767
# ===================================================
6868
# Begin tb05ad tests
6969
import numpy as np
70-
from control import ss
7170

72-
from scipy import linalg
7371
from numpy.testing import assert_raises, assert_almost_equal
7472
from numpy.testing import assert_array_almost_equal, dec
7573

7674
CASES = {}
7775

78-
CASES['fail1'] = ss(np.array([[-0.5, 0., 0., 0. ],
76+
CASES['fail1'] = {'A': np.array([[-0.5, 0., 0., 0. ],
7977
[ 0., -1., 0. , 0. ],
8078
[ 1., 0., -0.5, 0. ],
8179
[ 0., 1., 0., -1. ]]),
82-
np.array([[ 1., 0.],
80+
'B': np.array([[ 1., 0.],
8381
[ 0., 1.],
8482
[ 0., 0.],
8583
[ 0., 0.]]),
86-
np.array([[ 0., 1., 1., 0.],
84+
'C': np.array([[ 0., 1., 1., 0.],
8785
[ 0., 1., 0., 1.],
88-
[ 0., 1., 1., 1.]]),
89-
np.zeros([3,2]) )
86+
[ 0., 1., 1., 1.]])}
9087

9188
n = 20
9289
p = 10
9390
m = 14
9491
np.random.seed(40)
95-
CASES['pass1'] = ss(np.random.rand(n, n),
96-
np.random.rand(n, m),
97-
np.random.rand(p, n),
98-
np.zeros([p,m]))
92+
93+
CASES['pass1'] = {'A': np.random.rand(n, n),
94+
'B': np.random.rand(n, m),
95+
'C': np.random.rand(p, n)}
96+
9997

10098
class test_tb05ad(unittest.TestCase):
10199

102100
def setUp(self):
103101
pass
104102

105-
CASES = {}
106-
107-
CASES['fail1'] = ss(np.array([[-0.5, 0., 0., 0. ],
108-
[ 0., -1., 0. , 0. ],
109-
[ 1., 0., -0.5, 0. ],
110-
[ 0., 1., 0., -1. ]]),
111-
np.array([[ 1., 0.],
112-
[ 0., 1.],
113-
[ 0., 0.],
114-
[ 0., 0.]]),
115-
np.array([[ 0., 1., 1., 0.],
116-
[ 0., 1., 0., 1.],
117-
[ 0., 1., 1., 1.]]),
118-
np.zeros([3,2]) )
119-
120-
n = 20
121-
p = 10
122-
m = 14
123-
np.random.seed(40)
124-
CASES['pass1'] = ss(np.random.rand(n, n),
125-
np.random.rand(n, m),
126-
np.random.rand(p, n),
127-
np.zeros([p,m]))
128-
129-
130103
def test_tb05ad_ng(self):
131104
for key in CASES:
132105
sys = CASES[key]
133106
self.check_tb05ad_AG_NG(sys, 10*1j, 'NG')
134107

135-
136-
@dec.knownfailureif
108+
@dec.knownfailureif(True)
137109
def test_tb05ad_ag_failure(self):
138110
# Test the failure when we do balancing on certain A matrices.
139-
self.check_tb05ad_AG_NG(CASES['fail1'], 'AG')
140-
111+
self.check_tb05ad_AG_NG(CASES['fail1'], 10*1j, 'AG')
141112

142113
def test_tb05ad_nh(self):
143114
# Test the conversion to Hessenberg form and subsequent solution.
144115
jomega = 10*1j
145116
for key in CASES:
146117
sys = CASES[key]
147-
sys2 = self.check_tb05ad_AG_NG(sys, jomega, 'NG')
148-
self.check_tb05ad_NH(sys2, sys, jomega)
149-
118+
sys_transformed = self.check_tb05ad_AG_NG(sys, jomega, 'NG')
119+
self.check_tb05ad_NH(sys_transformed, sys, jomega)
150120

151121
# Check error handling.
152122
def test_tb05ad_errors(self):
153123
self.check_tb05ad_errors(CASES['pass1'])
154124

155-
156125
def check_tb05ad_AG_NG(self, sys, jomega, job):
157-
result = transform.tb05ad(sys.states, sys.inputs, sys.outputs, jomega,
158-
sys.A, sys.B, sys.C, job=job)
126+
n, m = sys['B'].shape
127+
p = sys['C'].shape[0]
128+
result = transform.tb05ad(n, m, p, jomega,
129+
sys['A'], sys['B'], sys['C'], job=job)
159130
g_i = result[3]
160-
hinvb = linalg.solve(np.eye(sys.states) * jomega - sys.A, sys.B)
161-
g_i_solve = sys.C.dot(hinvb)
131+
hinvb = np.linalg.solve(np.eye(n) * jomega - sys['A'], sys['B'])
132+
g_i_solve = sys['C'].dot(hinvb)
162133
np.testing.assert_almost_equal(g_i_solve, g_i)
163-
return ss(result[0], result[1], result[2], sys.D)
164-
134+
sys_transformed = {'A': result[0], 'B': result[1], 'C': result[2]}
135+
return sys_transformed
165136

166137
def check_tb05ad_NH(self, sys_transformed, sys, jomega):
167-
# When input matrices are already Hessenberg, output format changes.
168-
result = transform.tb05ad(sys_transformed.states, sys_transformed.inputs,
169-
sys_transformed.outputs, jomega,
170-
sys_transformed.A, sys_transformed.B,
171-
sys_transformed.C, job='NH')
138+
# Check that computing frequency response on the transformed
139+
# system gives the same result as computing C(sI - A)^-1B with
140+
# the original system.
141+
142+
n, m = sys_transformed['B'].shape
143+
p = sys_transformed['C'].shape[0]
144+
result = transform.tb05ad(n, m, p, jomega, sys_transformed['A'],
145+
sys_transformed['B'], sys_transformed['C'],
146+
job='NH')
172147
g_i = result[0]
173-
hinvb = linalg.solve(np.eye(sys.states) * jomega - sys.A, sys.B)
174-
g_i_solve = sys.C.dot(hinvb)
148+
hinvb = np.linalg.solve(np.eye(n) * jomega - sys['A'], sys['B'])
149+
g_i_solve = sys['C'].dot(hinvb)
175150
np.testing.assert_almost_equal(g_i_solve, g_i)
176151

177-
178152
def check_tb05ad_errors(self, sys):
153+
n, m = sys['B'].shape
154+
p = sys['C'].shape[0]
179155
jomega = 10*1j
180156
# test error handling
181157
# wrong size A
182-
assert_raises(ValueError, transform.tb05ad, sys.states+1,
183-
sys.inputs, sys.outputs, jomega, sys.A, sys.B, sys.C, job='NH')
158+
assert_raises(ValueError, transform.tb05ad, n+1, m, p,
159+
jomega, sys['A'], sys['B'], sys['C'], job='NH')
184160
# wrong size B
185-
assert_raises(ValueError, transform.tb05ad, sys.states,
186-
sys.inputs+1, sys.outputs, jomega, sys.A, sys.B, sys.C, job='NH')
161+
assert_raises(ValueError, transform.tb05ad, n, m+1, p,
162+
jomega, sys['A'], sys['B'], sys['C'], job='NH')
187163
# wrong size C
188-
assert_raises(ValueError, transform.tb05ad,sys.states,
189-
sys.inputs, sys.outputs+1, jomega, sys.A, sys.B, sys.C, job='NH')
164+
assert_raises(ValueError, transform.tb05ad, n, m, p+1,
165+
jomega, sys['A'], sys['B'], sys['C'], job='NH')
190166
# unrecognized job
191-
assert_raises(ValueError, transform.tb05ad,sys.states,
192-
sys.inputs, sys.outputs, jomega, sys.A, sys.B, sys.C, job='a')
167+
assert_raises(ValueError, transform.tb05ad, n, m, p, jomega,
168+
sys['A'], sys['B'], sys['C'], job='a')
193169

194170

195171

‎slycot/transform.py‎

Lines changed: 12 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -334,7 +334,7 @@ def tb04ad(n,m,p,A,B,C,D,tol1=0.0,tol2=0.0,ldwork=None):
334334
return A[:Nr,:Nr],B[:Nr,:m],C[:p,:Nr],Nr,index,dcoeff[:porm,:kdcoef],ucoeff[:porm,:porp,:kdcoef]
335335

336336

337-
def tb05ad(n, m, p, jomega, A, B, C, job='all'):
337+
def tb05ad(n, m, p, jomega, A, B, C, job='NG'):
338338
"""At, Bt, Ct, g_jw, rcond, ev, hinvb = tb05ad_a(n, m, p, jomega, A, B, C):
339339
340340
To find the complex frequency response matrix (transfer matrix)
@@ -483,19 +483,22 @@ def tb05ad(n, m, p, jomega, A, B, C, job='all'):
483483
Example
484484
-------
485485
>>> A = np.array([[0.0, 1.0],
486-
[-100.0, -20.1]),
486+
[-100.0, -20.1]])
487487
>>> B = np.array([[0.],[100]])
488488
>>> C = np.array([[1., 0.]])
489489
>>> n = np.shape(A)[0]
490490
>>> m = np.shape(B)[1]
491491
>>> p = np.shape(C)[0]
492492
>>> jw_s = [1j*10, 1j*15]
493-
>>> at, bt, ct, g_1, rcond, ev, hinvb = tb05ad(n, m, p, jw_s[0], A, B, C, job='NG')
494-
>>> g_2, hinv2 = tb05ad(n, m, p, jw_s[1], at, bt, ct, job='NH')
493+
>>> at, bt, ct, g_1, hinvb, info = slycot.tb05ad(n, m, p, jw_s[0],
494+
A, B, C, job='NG')
495+
>>> g_2, hinv2,info = slycot.tb05ad(n, m, p, jw_s[1], at, bt, ct, job='NH')
495496
496497
"""
497498
def error_handler(out, arg_list):
498499
if out[-1] < 0:
500+
# Conform fortran 1-based argument indexing to
501+
# to python zero indexing.
499502
error_text = ("The following argument had an illegal value: "
500503
+ arg_list[-out[-1]-1])
501504
e = ValueError(error_text)
@@ -525,25 +528,25 @@ def error_handler(out, arg_list):
525528
# TB05AD(baleig,inita,n,m,p,freq,a,lda,b,ldb,c,ldc,rcond,g,ldg,evre,evim,hinvb,ldhinv,
526529
# iwork,dwork,ldwork,zwork,lzwork,info)
527530

528-
# Sanity check on matrix dimenations
531+
# Sanity check on matrix dimensions
529532
if A.shape != (n, n):
530533
e = ValueError("The shape of A is (" + str(A.shape[0]) + "," +
531534
str(A.shape[1]) + "), but expected (" + str(n) +
532-
"," + str(m) + ")")
533-
e.info = -6
535+
"," + str(n) + ")")
536+
e.info = -6
534537
raise e
535538

536539
if B.shape != (n, m):
537540
e = ValueError("The shape of B is (" + str(B.shape[0]) + "," +
538541
str(B.shape[1]) + "), but expected (" + str(n) +
539542
"," + str(m) + ")")
540-
e.info = -6
543+
e.info = -8
541544
raise e
542545
if C.shape != (p, n):
543546
e = ValueError("The shape of C is (" + str(C.shape[0]) + "," +
544547
str(C.shape[1]) + "), but expected (" + str(p) +
545548
"," + str(n) + ")")
546-
e.info = -8
549+
e.info = -10
547550
raise e
548551

549552
# ----------------------------------------------------

0 commit comments

Comments
 (0)