Reference-LAPACK / lapack

LAPACK development repository
Other
1.49k stars 434 forks source link

dgges() not working as expected on scipy #689

Closed marcdelabarrera closed 2 years ago

marcdelabarrera commented 2 years ago

Describe your issue.

I'm a python user using scipy to perform a QZ decomposition. I'm testing the scipy.linalg.qz function, with complex results, which calls the LATPACK function dgges().

If I perform the operation VSL@S@VSR.T, I would expect to recover A, but it's not the case. If instead, I do VSL.conjugate()@S.conjugate()@VSR.T, then I do recover A. The same happens with B. Moreover, I would expect VSL to be orthonormal so VSL@VSL.T=I, but it's not the case either.

Is this behavior expected? Here you have a reproducible code (in python)

Reproducing Code Example

A = np.array([[1,2],[3,4]]).astype('D')
B = np.array([[3,4],[4,4]]).astype('D')
gges, = get_lapack_funcs(('gges',), (A, B))
S, T,_,_,_,VSL,VSR,_,info = gges(lambda x: None, A, B)
print(VSL@S@VSR.T)
print(VSL.conjugate()@S.conjugate()@VSR.T)
ilayn commented 2 years ago

You are asking a for complex QZ decomposition hence the result will be a triangular complex array which necessitates a conjugation.

There is no check if the data has 0. imaginary parts and treats it as a complex array.

marcdelabarrera commented 2 years ago

Extremely helpful, thank you. Now VSL@S@VSR.conjugate().T=A and VSL@VSL.conjugate().T=I