[Numpy-discussion] Zero row in SVD's unitary matrix on some Mac's
Pauli Virtanen
pav at iki.fi
Tue Apr 26 12:49:11 EDT 2011
Tue, 26 Apr 2011 11:36:19 -0500, Jason Grout wrote:
[clip]
> Okay, just one more data point. Our people that are seeing the problem
> with numpy returning a non-unitary V also see a non-unitary V being
> returned by the test C call to zgesdd. In other words, it really
> appears that zgesdd follows the Lapack docs, and if rows<columns, the
> returned V is not necessarily unitary, but may contain a zero row. This
> makes numpy's assumptions in using zgesdd false.
>
> You can see this report at
> http://trac.sagemath.org/sage_trac/ticket/11248#comment:25
What LAPACK promises is not very clear. Earlier on the the man page:
JOBZ (input) CHARACTER*1
Specifies options for computing all or part of the matrix U:
= 'A': all M columns of U and all N rows of V**H are returned
in the arrays U and VT; = 'S': the first min(M,N) columns of U
and the first min(M,N) rows of V**H are returned in the arrays
U and VT; = 'O': If M >= N, the first N columns of U are over‐
written in the array A and all rows of V**H are returned in the
array VT; otherwise, all columns of U are returned in the array
U and the first M rows of V**H are overwritten in the array A;
= 'N': no columns of U or rows of V**H are computed.
Looks to me the statement should be parsed as:
return_n_rows = (jobz == 'A') or (jobz == 'O' and m >= n)
So the current usage should be OK. (At least, as long as jobz == 'A',
in the other cases, I don't think the wrapper does the right thing.)
But apparently either there's a bug, or the LAPACK man page needs to
be understood as you say.
Pauli
More information about the NumPy-Discussion
mailing list