Just checked out the last SVN version of LAPACK.
The issue I have is related to the computation of workspace size in ?GESVD subroutines
Here below is a test case ZGESVD (file: test_zgesvd.f):
- Code: Select all
* .. Parameters ..
INTEGER NIN, NOUT
PARAMETER (NIN=5,NOUT=6)
INTEGER MMAX, NMAX
PARAMETER (MMAX=10,NMAX=24)
INTEGER LDA, LDVT
PARAMETER (LDA=MMAX,LDVT=NMAX)
* .. Local Scalars ..
INTEGER INFO, M, N
* .. Local Arrays ..
COMPLEX *16 A(LDA,NMAX), DUMMY(1,1), WORK(1,1)
DOUBLE PRECISION RWORK(5*MIN(MMAX,NMAX))
* .. Executable Statements ..
WRITE (NOUT,*) 'F08KPF Example Program Results'
WRITE (NOUT,*)
* Skip heading in data file
READ (NIN,*)
READ (NIN,*) M, N
WRITE (NOUT,*) 'READ>> M: ', M, ' - N: ', N
IF (M.LE.MMAX .AND. N.LE.NMAX) THEN
*
* Read the m by n matrix A from data file
*
READ (NIN,*) ((A(I,J),J=1,N),I=1,M)
*
* Compute the singular values and left and right singular vectors
* of A (A = U*S*(V**H), m.ge.n)
*
CALL ZGESVD('N','N',M,N,A,
+ LDA,S,DUMMY,1,DUMMY,1,WORK,-1,RWORK,INFO)
ELSE
WRITE (NOUT,*) 'MMAX and/or NMAX too small'
END IF
END
and the input matrix is (file: test_zgesvd.d):
- Code: Select all
6 12 :Values of M and N
(-1.377823724188512,0) (0,0) (0,0) (1,0) (0,0) (0,0) (0,0) (0,0) (0,0) (0,0) (0,0) (0,0)
(0,0) (-1.377823724188512,0) (0,0) (0,0) (1,0) (0,0) (0,0) (0,0) (0,0) (0,0) (0,0) (0,0)
(0,0) (0,0) (-1.377823724188512,0) (0,0) (0,0) (1,0) (0,0) (0,0) (0,0) (0,0) (0,0) (0,0)
(0.01588062702035706,0) (0,0) (0,0) (-1.382229846126757,0) (0,0) (0,0) (0.005286842551427753,0) (-0.003762200585945338,0) (-0.007151269102833426,0) (0.004018458066065118,0) (-0.002859598171060272,0) (-0.005435583664416818,0)
(0,0) (0.07630229180736924,0) (0,0) (0,0) (-1.28781829905337,0) (0,0) (0.01156945407605109,0) (-0.008233006086446525,0) (-0.01564946916158397,0) (0.04498579155453522,0) (-0.03201259914577823,0) (-0.0608502141081553,0)
(0,0) (0,0) (0.4070646177524311,0) (0,0) (0,0) (-0.2954402733863342,0) (-0.5243179610007159,0) (0.373112934782985,0) (0.7092208247345488,0) (0.02138338930221806,0) (-0.01521675725723726,0) (-0.02892432860319349,0) :End of matrix A
The output I get is the error:
- Code: Select all
** On entry to ZUNGLQ parameter number 5 had an illegal value
which is issued by the call to ZUNGLQ done at line 476 of ZGESVD.f (i.e., when calculating the space for ZUNGLQ)
I have not tried, but similar instructions are present in the other variant of ?GESVD.
I think the problem is related to the fact that since JOBVT='N' I pass LDVT as 1 which seems to be legal since the ?GESVD doc claims that LDVT>=1 and if JOBVT='N' the VT is not dereferenced.
A possible workaround that seems to work for me is to replace lines 476,477:
- Code: Select all
CALL ZUNGLQ( N, N, M, VT, LDVT, DUM(1), DUM(1), -1, IERR )
LWORK_ZUNGLQ_N=DUM(1)
with:
- Code: Select all
IF( .NOT.WNTVN ) THEN
CALL ZUNGLQ( N, N, M, VT, LDVT, DUM(1), DUM(1), -1,
$ IERR )
END IF
What do you think?
Thank you very much,
-- Marco
