Is this code from linalg.lstsq for the complex case correct?
lapack_routine = lapack_lite.zgelsd
lwork = 1
rwork = zeros((lwork,), real_t)
work = zeros((lwork,),t)
results = lapack_routine(m, n, n_rhs, a, m, bstar, ldb, s, rcond,
0, work, -1, rwork, iwork, 0)
lwork = int(abs(work[0]))
rwork = zeros((lwork,),real_t)
a_real = zeros((m,n),real_t)
bstar_real = zeros((ldb,n_rhs,),real_t)
results = lapack_lite.dgelsd(m, n, n_rhs, a_real, m,
bstar_real, ldb, s, rcond,
0, rwork, -1, iwork, 0)
lrwork = int(rwork[0])
work = zeros((lwork,), t)
rwork = zeros((lrwork,), real_t)
results = lapack_routine(m, n, n_rhs, a, m, bstar, ldb, s, rcond,
The middle call to dgelsd looks unnecessary to me. At the very least,
allocating astar_real and bstar_real shouldn't be necessary since they
aren't referenced anywhere else in the lstsq function. The lapack
documentation for zgelsd also doesn't mention any need to call dgelsd
to compute the size of the work array.
--bb