dgejsv bug?

[cottrell]

I'm trying out the new Jacobi solver, but having a hard time figuring out
the workspace requirements. As a general comment, it would be very
helpful if this routine had the apparatus offered by many others, where
you can give LWORK = -1 and get a response indicating a suitable workspace
size.

But anyway, take a 3x3 input matrix (M = N = 3), and say we want the left or
right singular vectors but not both. By my reading of dgejsv.f, the minimal
LWORK is supposed to be max(2*N+M,7) = 9. But if I give LWORK = 9, my
program aborts with
** On entry to DGEQP3 parameter number 8 had an illegal value

Param 8 is dgeqp3's LWORK, passed by dgejsv; and looking inside
dgejsv, the first such call is

CALL DGEQP3( M,N,A,LDA, IWORK,WORK, WORK(N+1),LWORK-N, IERR )

With N = 3 and LWORK = 9 on input, this sends dgeqp3 an LWORK of 6
= 2*N. Yet inside dgeqp3 we find the requirement that LWORK >= 3*N+1.
So the abort is not surprising: dgejsv is making an illegal call to
dgeqp3. It seems to me that dgejsv has to specify a bigger minimal
LWORK.

[cottrell]

To be more specific: I believe that on account of the xGEQP3 calls,
the minimal LWORK for dgejsv, for the case of one or other but
not both of the left- and right singular vectors wanted, should be
max(2*N+M,5*N+1,7), rather than max(2*N+M,7) as given.

[cottrell]

Hello again -- any confirmation/denial that this is a bug?

[Julien Langou]

Hello, thanks for your help on the forum and this bug report. Sorry for the delay. Yes there is a bug in DGEJSV, I confirm. When I perform the same experiment as you did, I have the same crash and your fix fixes it. So that's for the confirmation. Zlatko (original developer of the code) should have a look soon. Please use DGEJSV with your fix for the time being and we'll come back to you. Julien.

[cottrell]

Thanks for the confirmation! This function seems to work nicely when the
workspace issue is handled correctly.

[admin]

This has been corrected in our repository.
See http://www.netlib.org/lapack/Errata/index2.html
Thank you