Skip to content

Commit

Permalink
MAINT: disable sort keyword of linalg.qz, due to Windows issues. See #…
Browse files Browse the repository at this point in the history
  • Loading branch information
rgommers committed Sep 23, 2012
1 parent 5fd5866 commit e1c61f1
Show file tree
Hide file tree
Showing 2 changed files with 57 additions and 47 deletions.
15 changes: 11 additions & 4 deletions scipy/linalg/_decomp_qz.py
Original file line number Original file line Diff line number Diff line change
Expand Up @@ -83,6 +83,8 @@ def qz(A, B, output='real', lwork=None, sort=None, overwrite_a=False,
lwork : int, optional lwork : int, optional
Work array size. If None or -1, it is automatically computed. Work array size. If None or -1, it is automatically computed.
sort : {None, callable, 'lhp', 'rhp', 'iuc', 'ouc'}, optional sort : {None, callable, 'lhp', 'rhp', 'iuc', 'ouc'}, optional
NOTE: THIS INPUT IS DISABLED FOR NOW, IT DOESN'T WORK WELL ON WINDOWS.
Specifies whether the upper eigenvalues should be sorted. A callable Specifies whether the upper eigenvalues should be sorted. A callable
may be passed that, given a eigenvalue, returns a boolean denoting may be passed that, given a eigenvalue, returns a boolean denoting
whether the eigenvalue should be sorted to the top-left (True). For whether the eigenvalue should be sorted to the top-left (True). For
Expand Down Expand Up @@ -146,6 +148,11 @@ def qz(A, B, output='real', lwork=None, sort=None, overwrite_a=False,
[-0.54861681, -0.6210585 , -0.55973739]]) [-0.54861681, -0.6210585 , -0.55973739]])
""" """
if sort is not None:
# Disabled due to segfaults on win32, see ticket 1717.
raise ValueError("The 'sort' input of qz() has to be None (will "
" change when this functionality is made more robust).")

if not output in ['real','complex','r','c']: if not output in ['real','complex','r','c']:
raise ValueError("argument must be 'real', or 'complex'") raise ValueError("argument must be 'real', or 'complex'")


Expand Down Expand Up @@ -201,14 +208,14 @@ def qz(A, B, output='real', lwork=None, sort=None, overwrite_a=False,
raise ValueError("Illegal value in argument %d of gges" % -info) raise ValueError("Illegal value in argument %d of gges" % -info)
elif info > 0 and info <= a_n: elif info > 0 and info <= a_n:
warnings.warn("The QZ iteration failed. (a,b) are not in Schur " warnings.warn("The QZ iteration failed. (a,b) are not in Schur "
"form, but ALPHAR(j), ALPHAI(j), and BETA(j) should be correct" "form, but ALPHAR(j), ALPHAI(j), and BETA(j) should be correct "
"for J=%d,...,N" % info-1, UserWarning) "for J=%d,...,N" % info-1, UserWarning)
elif info == a_n+1: elif info == a_n+1:
raise LinAlgError("Something other than QZ iteration failed") raise LinAlgError("Something other than QZ iteration failed")
elif info == a_n+2: elif info == a_n+2:
raise LinAlgError("After reordering, roundoff changed values of some" raise LinAlgError("After reordering, roundoff changed values of some "
"complex eigenvalues so that leading eigenvalues in the" "complex eigenvalues so that leading eigenvalues in the "
"Generalized Schur form no longer satisfy sort=True." "Generalized Schur form no longer satisfy sort=True. "
"This could also be caused due to scaling.") "This could also be caused due to scaling.")
elif info == a_n+3: elif info == a_n+3:
raise LinAlgError("Reordering failed in <s,d,c,z>tgsen") raise LinAlgError("Reordering failed in <s,d,c,z>tgsen")
Expand Down
89 changes: 46 additions & 43 deletions scipy/linalg/tests/test_decomp.py
Original file line number Original file line Diff line number Diff line change
Expand Up @@ -1715,34 +1715,36 @@ def test_qz_double_sort(self):


sort = lambda ar,ai,beta : ai == 0 sort = lambda ar,ai,beta : ai == 0


AA,BB,Q,Z,sdim = qz(A,B,sort=sort) assert_raises(ValueError, qz, A, B, sort=sort)
#assert_(sdim == 2) if False:
assert_(sdim == 4) AA,BB,Q,Z,sdim = qz(A,B,sort=sort)
assert_array_almost_equal(dot(dot(Q,AA),Z.T), A) #assert_(sdim == 2)
assert_array_almost_equal(dot(dot(Q,BB),Z.T), B) assert_(sdim == 4)

assert_array_almost_equal(dot(dot(Q,AA),Z.T), A)
# test absolute values bc the sign is ambiguous and might be platform assert_array_almost_equal(dot(dot(Q,BB),Z.T), B)
# dependent
assert_array_almost_equal(np.abs(AA), np.abs(np.array( # test absolute values bc the sign is ambiguous and might be platform
[[ 35.7864, -80.9061, -12.0629, -9.498 ], # dependent
[ 0. , 2.7638, -2.3505, 7.3256], assert_array_almost_equal(np.abs(AA), np.abs(np.array(
[ 0. , 0. , 0.6258, -0.0398], [[ 35.7864, -80.9061, -12.0629, -9.498 ],
[ 0. , 0. , 0. , -12.8217]])), 4) [ 0. , 2.7638, -2.3505, 7.3256],
assert_array_almost_equal(np.abs(BB), np.abs(np.array( [ 0. , 0. , 0.6258, -0.0398],
[[ 4.5324, -8.7878, 3.2357, -3.5526], [ 0. , 0. , 0. , -12.8217]])), 4)
[ 0. , 1.4314, -2.1894, 0.9709], assert_array_almost_equal(np.abs(BB), np.abs(np.array(
[ 0. , 0. , 1.3126, -0.3468], [[ 4.5324, -8.7878, 3.2357, -3.5526],
[ 0. , 0. , 0. , 0.559 ]])), 4) [ 0. , 1.4314, -2.1894, 0.9709],
assert_array_almost_equal(np.abs(Q), np.abs(np.array( [ 0. , 0. , 1.3126, -0.3468],
[[-0.4193, -0.605 , -0.1894, -0.6498], [ 0. , 0. , 0. , 0.559 ]])), 4)
[-0.5495, 0.6987, 0.2654, -0.3734], assert_array_almost_equal(np.abs(Q), np.abs(np.array(
[-0.4973, -0.3682, 0.6194, 0.4832], [[-0.4193, -0.605 , -0.1894, -0.6498],
[-0.5243, 0.1008, -0.7142, 0.4526]])), 4) [-0.5495, 0.6987, 0.2654, -0.3734],
assert_array_almost_equal(np.abs(Z), np.abs(np.array( [-0.4973, -0.3682, 0.6194, 0.4832],
[[-0.9471, -0.2971, -0.1217, 0.0055], [-0.5243, 0.1008, -0.7142, 0.4526]])), 4)
[-0.0367, 0.1209, 0.0358, 0.9913], assert_array_almost_equal(np.abs(Z), np.abs(np.array(
[ 0.3171, -0.9041, -0.2547, 0.1312], [[-0.9471, -0.2971, -0.1217, 0.0055],
[ 0.0346, 0.2824, -0.9587, 0.0014]])), 4) [-0.0367, 0.1209, 0.0358, 0.9913],
[ 0.3171, -0.9041, -0.2547, 0.1312],
[ 0.0346, 0.2824, -0.9587, 0.0014]])), 4)


# test absolute values bc the sign is ambiguous and might be platform # test absolute values bc the sign is ambiguous and might be platform
# dependent # dependent
Expand All @@ -1769,24 +1771,25 @@ def test_qz_double_sort(self):
# [0.0626, -0.6934, -0.7114, 0.0956]])), 4) # [0.0626, -0.6934, -0.7114, 0.0956]])), 4)
#assert_array_almost_equal(dot(Z,Z.T), eye(4)) #assert_array_almost_equal(dot(Z,Z.T), eye(4))


def test_qz_complex_sort(self): #def test_qz_complex_sort(self):
cA = np.array([ # cA = np.array([
[-21.10+22.50*1j, 53.50+-50.50*1j, -34.50+127.50*1j, 7.50+ 0.50*1j], # [-21.10+22.50*1j, 53.50+-50.50*1j, -34.50+127.50*1j, 7.50+ 0.50*1j],
[-0.46+ -7.78*1j, -3.50+-37.50*1j, -15.50+ 58.50*1j,-10.50+ -1.50*1j], # [-0.46+ -7.78*1j, -3.50+-37.50*1j, -15.50+ 58.50*1j,-10.50+ -1.50*1j],
[ 4.30+ -5.50*1j, 39.70+-17.10*1j, -68.50+ 12.50*1j, -7.50+ -3.50*1j], # [ 4.30+ -5.50*1j, 39.70+-17.10*1j, -68.50+ 12.50*1j, -7.50+ -3.50*1j],
[ 5.50+ 4.40*1j, 14.40+ 43.30*1j, -32.50+-46.00*1j,-19.00+-32.50*1j]]) # [ 5.50+ 4.40*1j, 14.40+ 43.30*1j, -32.50+-46.00*1j,-19.00+-32.50*1j]])

# cB = np.array([
# [1.00+ -5.00*1j, 1.60+ 1.20*1j,-3.00+ 0.00*1j, 0.00+ -1.00*1j],
# [0.80+ -0.60*1j, 3.00+ -5.00*1j,-4.00+ 3.00*1j,-2.40+ -3.20*1j],
# [1.00+ 0.00*1j, 2.40+ 1.80*1j,-4.00+ -5.00*1j, 0.00+ -3.00*1j],
# [0.00+ 1.00*1j,-1.80+ 2.40*1j, 0.00+ -4.00*1j, 4.00+ -5.00*1j]])


cB = np.array([ # AAS,BBS,QS,ZS,sdim = qz(cA,cB,sort='lhp')
[1.00+ -5.00*1j, 1.60+ 1.20*1j,-3.00+ 0.00*1j, 0.00+ -1.00*1j],
[0.80+ -0.60*1j, 3.00+ -5.00*1j,-4.00+ 3.00*1j,-2.40+ -3.20*1j],
[1.00+ 0.00*1j, 2.40+ 1.80*1j,-4.00+ -5.00*1j, 0.00+ -3.00*1j],
[0.00+ 1.00*1j,-1.80+ 2.40*1j, 0.00+ -4.00*1j, 4.00+ -5.00*1j]])


AAS,BBS,QS,ZS,sdim = qz(cA,cB,sort='lhp') # eigenvalues = diag(AAS)/diag(BBS)
# assert_(all(np.real(eigenvalues[:sdim] < 0)))
# assert_(all(np.real(eigenvalues[sdim:] > 0)))


eigenvalues = diag(AAS)/diag(BBS)
assert_(all(np.real(eigenvalues[:sdim] < 0)))
assert_(all(np.real(eigenvalues[sdim:] > 0)))


class TestDatacopied(TestCase): class TestDatacopied(TestCase):


Expand Down

0 comments on commit e1c61f1

Please sign in to comment.