From 10c2ebdfc5c60f495c5c34387d9a3c3fbadee843 Mon Sep 17 00:00:00 2001 From: Werner Saar Date: Mon, 7 Mar 2016 10:34:04 +0100 Subject: [PATCH] BUGFIX: removed fixes for bugs #148 and #149, because info for xerbla is wrong --- lapack-netlib/SRC/cgesvdx.f | 85 +++++++++++----------------- lapack-netlib/SRC/dgesvdx.f | 74 ++++++------------------- lapack-netlib/SRC/sgesvdx.f | 74 ++++++------------------- lapack-netlib/SRC/zgesvdx.f | 132 ++++++++++++++++++-------------------------- 4 files changed, 124 insertions(+), 241 deletions(-) diff --git a/lapack-netlib/SRC/cgesvdx.f b/lapack-netlib/SRC/cgesvdx.f index 87ea986..235426a 100644 --- a/lapack-netlib/SRC/cgesvdx.f +++ b/lapack-netlib/SRC/cgesvdx.f @@ -170,7 +170,7 @@ *> vectors, stored columnwise) as specified by RANGE; if *> JOBU = 'N', U is not referenced. *> Note: The user must ensure that UCOL >= NS; if RANGE = 'V', -*> the exact value of NS is not known in advance and an upper +*> the exact value of NS is not known ILQFin advance and an upper *> bound must be used. *> \endverbatim *> @@ -294,8 +294,8 @@ CHARACTER JOBZ, RNGTGK LOGICAL ALLS, INDS, LQUERY, VALS, WANTU, WANTVT INTEGER I, ID, IE, IERR, ILQF, ILTGK, IQRF, ISCL, - $ ITAU, ITAUP, ITAUQ, ITEMP, ITEMPR, ITGKZ, - $ IUTGK, J, K, MAXWRK, MINMN, MINWRK, MNTHR + $ ITAU, ITAUP, ITAUQ, ITEMP, ITGKZ, IUTGK, + $ J, K, MAXWRK, MINMN, MINWRK, MNTHR REAL ABSTOL, ANRM, BIGNUM, EPS, SMLNUM * .. * .. Local Arrays .. @@ -367,14 +367,8 @@ IF( INFO.EQ.0 ) THEN IF( WANTU .AND. LDU.LT.M ) THEN INFO = -15 - ELSE IF( WANTVT ) THEN - IF( INDS ) THEN - IF( LDVT.LT.IU-IL+1 ) THEN - INFO = -17 - END IF - ELSE IF( LDVT.LT.MINMN ) THEN - INFO = -17 - END IF + ELSE IF( WANTVT .AND. LDVT.LT.MINMN ) THEN + INFO = -16 END IF END IF END IF @@ -396,24 +390,18 @@ * * Path 1 (M much larger than N) * - MINWRK = N*(N+5) - MAXWRK = N + N*ILAENV(1,'CGEQRF',' ',M,N,-1,-1) - MAXWRK = MAX(MAXWRK, - $ N*N+2*N+2*N*ILAENV(1,'CGEBRD',' ',N,N,-1,-1)) - IF (WANTU .OR. WANTVT) THEN - MAXWRK = MAX(MAXWRK, - $ N*N+2*N+N*ILAENV(1,'CUNMQR','LN',N,N,N,-1)) - END IF + MAXWRK = N + N* + $ ILAENV( 1, 'SGEQRF', ' ', M, N, -1, -1 ) + MAXWRK = MAX( MAXWRK, N*N + N + 2*N* + $ ILAENV( 1, 'SGEBRD', ' ', N, N, -1, -1 ) ) + MINWRK = N*(N+4) ELSE * * Path 2 (M at least N, but not much larger) * - MINWRK = 3*N + M - MAXWRK = 2*N + (M+N)*ILAENV(1,'CGEBRD',' ',M,N,-1,-1) - IF (WANTU .OR. WANTVT) THEN - MAXWRK = MAX(MAXWRK, - $ 2*N+N*ILAENV(1,'CUNMQR','LN',N,N,N,-1)) - END IF + MAXWRK = 2*N + ( M+N )* + $ ILAENV( 1, 'CGEBRD', ' ', M, N, -1, -1 ) + MINWRK = 2*N + M END IF ELSE MNTHR = ILAENV( 6, 'CGESVD', JOBU // JOBVT, M, N, 0, 0 ) @@ -421,25 +409,18 @@ * * Path 1t (N much larger than M) * - MINWRK = M*(M+5) - MAXWRK = M + M*ILAENV(1,'CGELQF',' ',M,N,-1,-1) - MAXWRK = MAX(MAXWRK, - $ M*M+2*M+2*M*ILAENV(1,'CGEBRD',' ',M,M,-1,-1)) - IF (WANTU .OR. WANTVT) THEN - MAXWRK = MAX(MAXWRK, - $ M*M+2*M+M*ILAENV(1,'CUNMQR','LN',M,M,M,-1)) - END IF + MAXWRK = M + M* + $ ILAENV( 1, 'CGELQF', ' ', M, N, -1, -1 ) + MAXWRK = MAX( MAXWRK, M*M + M + 2*M* + $ ILAENV( 1, 'CGEBRD', ' ', M, M, -1, -1 ) ) + MINWRK = M*(M+4) ELSE * * Path 2t (N greater than M, but not much larger) * -* - MINWRK = 3*M + N - MAXWRK = 2*M + (M+N)*ILAENV(1,'CGEBRD',' ',M,N,-1,-1) - IF (WANTU .OR. WANTVT) THEN - MAXWRK = MAX(MAXWRK, - $ 2*M+M*ILAENV(1,'CUNMQR','LN',M,M,M,-1)) - END IF + MAXWRK = M*(M*2+19) + ( M+N )* + $ ILAENV( 1, 'CGEBRD', ' ', M, N, -1, -1 ) + MINWRK = 2*M + N END IF END IF END IF @@ -537,14 +518,14 @@ CALL CGEBRD( N, N, WORK( IQRF ), N, RWORK( ID ), $ RWORK( IE ), WORK( ITAUQ ), WORK( ITAUP ), $ WORK( ITEMP ), LWORK-ITEMP+1, INFO ) - ITEMPR = ITGKZ + N*(N*2+1) + ITEMP = ITGKZ + N*(N*2+1) * * Solve eigenvalue problem TGK*Z=Z*S. * (Workspace: need 2*N*N+14*N) * CALL SBDSVDX( 'U', JOBZ, RNGTGK, N, RWORK( ID ), $ RWORK( IE ), VL, VU, ILTGK, IUTGK, NS, S, - $ RWORK( ITGKZ ), N*2, RWORK( ITEMPR ), + $ RWORK( ITGKZ ), N*2, RWORK( ITEMP ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -558,7 +539,7 @@ END DO K = K + N END DO - CALL CLASET( 'A', M-N, NS, CZERO, CZERO, U( N+1,1 ), LDU) + CALL CLASET( 'A', M-N, N, CZERO, CZERO, U( N+1,1 ), LDU ) * * Call CUNMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -613,14 +594,14 @@ CALL CGEBRD( M, N, A, LDA, RWORK( ID ), RWORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), WORK( ITEMP ), $ LWORK-ITEMP+1, INFO ) - ITEMPR = ITGKZ + N*(N*2+1) + ITEMP = ITGKZ + N*(N*2+1) * * Solve eigenvalue problem TGK*Z=Z*S. * (Workspace: need 2*N*N+14*N) * CALL SBDSVDX( 'U', JOBZ, RNGTGK, N, RWORK( ID ), $ RWORK( IE ), VL, VU, ILTGK, IUTGK, NS, S, - $ RWORK( ITGKZ ), N*2, RWORK( ITEMPR ), + $ RWORK( ITGKZ ), N*2, RWORK( ITEMP ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -634,7 +615,7 @@ END DO K = K + N END DO - CALL CLASET( 'A', M-N, NS, CZERO, CZERO, U( N+1,1 ), LDU) + CALL CLASET( 'A', M-N, N, CZERO, CZERO, U( N+1,1 ), LDU ) * * Call CUNMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -700,14 +681,14 @@ CALL CGEBRD( M, M, WORK( ILQF ), M, RWORK( ID ), $ RWORK( IE ), WORK( ITAUQ ), WORK( ITAUP ), $ WORK( ITEMP ), LWORK-ITEMP+1, INFO ) - ITEMPR = ITGKZ + M*(M*2+1) + ITEMP = ITGKZ + M*(M*2+1) * * Solve eigenvalue problem TGK*Z=Z*S. * (Workspace: need 2*M*M+14*M) * CALL SBDSVDX( 'U', JOBZ, RNGTGK, M, RWORK( ID ), $ RWORK( IE ), VL, VU, ILTGK, IUTGK, NS, S, - $ RWORK( ITGKZ ), M*2, RWORK( ITEMPR ), + $ RWORK( ITGKZ ), M*2, RWORK( ITEMP ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -741,7 +722,7 @@ END DO K = K + M END DO - CALL CLASET( 'A', NS, N-M, CZERO, CZERO, + CALL CLASET( 'A', M, N-M, CZERO, CZERO, $ VT( 1,M+1 ), LDVT ) * * Call CUNMBR to compute (VB**T)*(PB**T) @@ -777,14 +758,14 @@ CALL CGEBRD( M, N, A, LDA, RWORK( ID ), RWORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), WORK( ITEMP ), $ LWORK-ITEMP+1, INFO ) - ITEMPR = ITGKZ + M*(M*2+1) + ITEMP = ITGKZ + M*(M*2+1) * * Solve eigenvalue problem TGK*Z=Z*S. * (Workspace: need 2*M*M+14*M) * CALL SBDSVDX( 'L', JOBZ, RNGTGK, M, RWORK( ID ), $ RWORK( IE ), VL, VU, ILTGK, IUTGK, NS, S, - $ RWORK( ITGKZ ), M*2, RWORK( ITEMPR ), + $ RWORK( ITGKZ ), M*2, RWORK( ITEMP ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -818,7 +799,7 @@ END DO K = K + M END DO - CALL CLASET( 'A', NS, N-M, CZERO, CZERO, + CALL CLASET( 'A', M, N-M, CZERO, CZERO, $ VT( 1,M+1 ), LDVT ) * * Call CUNMBR to compute VB**T * PB**T diff --git a/lapack-netlib/SRC/dgesvdx.f b/lapack-netlib/SRC/dgesvdx.f index 4588083..cfa2ff0 100644 --- a/lapack-netlib/SRC/dgesvdx.f +++ b/lapack-netlib/SRC/dgesvdx.f @@ -169,7 +169,7 @@ *> vectors, stored columnwise) as specified by RANGE; if *> JOBU = 'N', U is not referenced. *> Note: The user must ensure that UCOL >= NS; if RANGE = 'V', -*> the exact value of NS is not known in advance and an upper +*> the exact value of NS is not known ILQFin advance and an upper *> bound must be used. *> \endverbatim *> @@ -357,14 +357,8 @@ IF( INFO.EQ.0 ) THEN IF( WANTU .AND. LDU.LT.M ) THEN INFO = -15 - ELSE IF( WANTVT ) THEN - IF( INDS ) THEN - IF( LDVT.LT.IU-IL+1 ) THEN - INFO = -17 - END IF - ELSE IF( LDVT.LT.MINMN ) THEN - INFO = -17 - END IF + ELSE IF( WANTVT .AND. LDVT.LT.MINMN ) THEN + INFO = -16 END IF END IF END IF @@ -386,34 +380,18 @@ * * Path 1 (M much larger than N) * - MAXWRK = N + + MAXWRK = N*(N*2+16) + $ N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - MAXWRK = MAX( MAXWRK, N*(N+5) + 2*N* + MAXWRK = MAX( MAXWRK, N*(N*2+20) + 2*N* $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - IF (WANTU) THEN - MAXWRK = MAX(MAXWRK,N*(N*3+6)+N* - $ ILAENV( 1, 'DORMQR', ' ', N, N, -1, -1 ) ) - END IF - IF (WANTVT) THEN - MAXWRK = MAX(MAXWRK,N*(N*3+6)+N* - $ ILAENV( 1, 'DORMLQ', ' ', N, N, -1, -1 ) ) - END IF - MINWRK = N*(N*3+20) + MINWRK = N*(N*2+21) ELSE * * Path 2 (M at least N, but not much larger) * - MAXWRK = 4*N + ( M+N )* + MAXWRK = N*(N*2+19) + ( M+N )* $ ILAENV( 1, 'DGEBRD', ' ', M, N, -1, -1 ) - IF (WANTU) THEN - MAXWRK = MAX(MAXWRK,N*(N*2+5)+N* - $ ILAENV( 1, 'DORMQR', ' ', N, N, -1, -1 ) ) - END IF - IF (WANTVT) THEN - MAXWRK = MAX(MAXWRK,N*(N*2+5)+N* - $ ILAENV( 1, 'DORMLQ', ' ', N, N, -1, -1 ) ) - END IF - MINWRK = MAX(N*(N*2+19),4*N+M) + MINWRK = N*(N*2+20) + M END IF ELSE MNTHR = ILAENV( 6, 'DGESVD', JOBU // JOBVT, M, N, 0, 0 ) @@ -421,34 +399,18 @@ * * Path 1t (N much larger than M) * - MAXWRK = M + + MAXWRK = M*(M*2+16) + $ M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - MAXWRK = MAX( MAXWRK, M*(M+5) + 2*M* + MAXWRK = MAX( MAXWRK, M*(M*2+20) + 2*M* $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - IF (WANTU) THEN - MAXWRK = MAX(MAXWRK,M*(M*3+6)+M* - $ ILAENV( 1, 'DORMQR', ' ', M, M, -1, -1 ) ) - END IF - IF (WANTVT) THEN - MAXWRK = MAX(MAXWRK,M*(M*3+6)+M* - $ ILAENV( 1, 'DORMLQ', ' ', M, M, -1, -1 ) ) - END IF - MINWRK = M*(M*3+20) + MINWRK = M*(M*2+21) ELSE * -* Path 2t (N at least M, but not much larger) +* Path 2t (N greater than M, but not much larger) * - MAXWRK = 4*M + ( M+N )* + MAXWRK = M*(M*2+19) + ( M+N )* $ ILAENV( 1, 'DGEBRD', ' ', M, N, -1, -1 ) - IF (WANTU) THEN - MAXWRK = MAX(MAXWRK,M*(M*2+5)+M* - $ ILAENV( 1, 'DORMQR', ' ', M, M, -1, -1 ) ) - END IF - IF (WANTVT) THEN - MAXWRK = MAX(MAXWRK,M*(M*2+5)+M* - $ ILAENV( 1, 'DORMLQ', ' ', M, M, -1, -1 ) ) - END IF - MINWRK = MAX(M*(M*2+19),4*M+N) + MINWRK = M*(M*2+20) + N END IF END IF END IF @@ -560,7 +522,7 @@ CALL DCOPY( N, WORK( J ), 1, U( 1,I ), 1 ) J = J + N*2 END DO - CALL DLASET( 'A', M-N, NS, ZERO, ZERO, U( N+1,1 ), LDU ) + CALL DLASET( 'A', M-N, N, ZERO, ZERO, U( N+1,1 ), LDU ) * * Call DORMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -629,7 +591,7 @@ CALL DCOPY( N, WORK( J ), 1, U( 1,I ), 1 ) J = J + N*2 END DO - CALL DLASET( 'A', M-N, NS, ZERO, ZERO, U( N+1,1 ), LDU ) + CALL DLASET( 'A', M-N, N, ZERO, ZERO, U( N+1,1 ), LDU ) * * Call DORMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -725,7 +687,7 @@ CALL DCOPY( M, WORK( J ), 1, VT( I,1 ), LDVT ) J = J + M*2 END DO - CALL DLASET( 'A', NS, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT) + CALL DLASET( 'A', M, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT ) * * Call DORMBR to compute (VB**T)*(PB**T) * (Workspace in WORK( ITEMP ): need M, prefer M*NB) @@ -794,7 +756,7 @@ CALL DCOPY( M, WORK( J ), 1, VT( I,1 ), LDVT ) J = J + M*2 END DO - CALL DLASET( 'A', NS, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT) + CALL DLASET( 'A', M, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT ) * * Call DORMBR to compute VB**T * PB**T * (Workspace in WORK( ITEMP ): need M, prefer M*NB) diff --git a/lapack-netlib/SRC/sgesvdx.f b/lapack-netlib/SRC/sgesvdx.f index 9128a7c..aae8b07 100644 --- a/lapack-netlib/SRC/sgesvdx.f +++ b/lapack-netlib/SRC/sgesvdx.f @@ -169,7 +169,7 @@ *> vectors, stored columnwise) as specified by RANGE; if *> JOBU = 'N', U is not referenced. *> Note: The user must ensure that UCOL >= NS; if RANGE = 'V', -*> the exact value of NS is not known in advance and an upper +*> the exact value of NS is not known ILQFin advance and an upper *> bound must be used. *> \endverbatim *> @@ -357,14 +357,8 @@ IF( INFO.EQ.0 ) THEN IF( WANTU .AND. LDU.LT.M ) THEN INFO = -15 - ELSE IF( WANTVT ) THEN - IF( INDS ) THEN - IF( LDVT.LT.IU-IL+1 ) THEN - INFO = -17 - END IF - ELSE IF( LDVT.LT.MINMN ) THEN - INFO = -17 - END IF + ELSE IF( WANTVT .AND. LDVT.LT.MINMN ) THEN + INFO = -16 END IF END IF END IF @@ -386,34 +380,18 @@ * * Path 1 (M much larger than N) * - MAXWRK = N + + MAXWRK = N*(N*2+16) + $ N*ILAENV( 1, 'SGEQRF', ' ', M, N, -1, -1 ) - MAXWRK = MAX( MAXWRK, N*(N+5) + 2*N* + MAXWRK = MAX( MAXWRK, N*(N*2+20) + 2*N* $ ILAENV( 1, 'SGEBRD', ' ', N, N, -1, -1 ) ) - IF (WANTU) THEN - MAXWRK = MAX(MAXWRK,N*(N*3+6)+N* - $ ILAENV( 1, 'SORMQR', ' ', N, N, -1, -1 ) ) - END IF - IF (WANTVT) THEN - MAXWRK = MAX(MAXWRK,N*(N*3+6)+N* - $ ILAENV( 1, 'SORMLQ', ' ', N, N, -1, -1 ) ) - END IF - MINWRK = N*(N*3+20) + MINWRK = N*(N*2+21) ELSE * * Path 2 (M at least N, but not much larger) * - MAXWRK = 4*N + ( M+N )* + MAXWRK = N*(N*2+19) + ( M+N )* $ ILAENV( 1, 'SGEBRD', ' ', M, N, -1, -1 ) - IF (WANTU) THEN - MAXWRK = MAX(MAXWRK,N*(N*2+5)+N* - $ ILAENV( 1, 'SORMQR', ' ', N, N, -1, -1 ) ) - END IF - IF (WANTVT) THEN - MAXWRK = MAX(MAXWRK,N*(N*2+5)+N* - $ ILAENV( 1, 'SORMLQ', ' ', N, N, -1, -1 ) ) - END IF - MINWRK = MAX(N*(N*2+19),4*N+M) + MINWRK = N*(N*2+20) + M END IF ELSE MNTHR = ILAENV( 6, 'SGESVD', JOBU // JOBVT, M, N, 0, 0 ) @@ -421,34 +399,18 @@ * * Path 1t (N much larger than M) * - MAXWRK = M + + MAXWRK = M*(M*2+16) + $ M*ILAENV( 1, 'SGELQF', ' ', M, N, -1, -1 ) - MAXWRK = MAX( MAXWRK, M*(M+5) + 2*M* + MAXWRK = MAX( MAXWRK, M*(M*2+20) + 2*M* $ ILAENV( 1, 'SGEBRD', ' ', M, M, -1, -1 ) ) - IF (WANTU) THEN - MAXWRK = MAX(MAXWRK,M*(M*3+6)+M* - $ ILAENV( 1, 'SORMQR', ' ', M, M, -1, -1 ) ) - END IF - IF (WANTVT) THEN - MAXWRK = MAX(MAXWRK,M*(M*3+6)+M* - $ ILAENV( 1, 'SORMLQ', ' ', M, M, -1, -1 ) ) - END IF - MINWRK = M*(M*3+20) + MINWRK = M*(M*2+21) ELSE * -* Path 2t (N at least M, but not much larger) +* Path 2t (N greater than M, but not much larger) * - MAXWRK = 4*M + ( M+N )* + MAXWRK = M*(M*2+19) + ( M+N )* $ ILAENV( 1, 'SGEBRD', ' ', M, N, -1, -1 ) - IF (WANTU) THEN - MAXWRK = MAX(MAXWRK,M*(M*2+5)+M* - $ ILAENV( 1, 'SORMQR', ' ', M, M, -1, -1 ) ) - END IF - IF (WANTVT) THEN - MAXWRK = MAX(MAXWRK,M*(M*2+5)+M* - $ ILAENV( 1, 'SORMLQ', ' ', M, M, -1, -1 ) ) - END IF - MINWRK = MAX(M*(M*2+19),4*M+N) + MINWRK = M*(M*2+20) + N END IF END IF END IF @@ -560,7 +522,7 @@ CALL SCOPY( N, WORK( J ), 1, U( 1,I ), 1 ) J = J + N*2 END DO - CALL SLASET( 'A', M-N, NS, ZERO, ZERO, U( N+1,1 ), LDU ) + CALL SLASET( 'A', M-N, N, ZERO, ZERO, U( N+1,1 ), LDU ) * * Call SORMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -629,7 +591,7 @@ CALL SCOPY( N, WORK( J ), 1, U( 1,I ), 1 ) J = J + N*2 END DO - CALL SLASET( 'A', M-N, NS, ZERO, ZERO, U( N+1,1 ), LDU ) + CALL SLASET( 'A', M-N, N, ZERO, ZERO, U( N+1,1 ), LDU ) * * Call SORMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -725,7 +687,7 @@ CALL SCOPY( M, WORK( J ), 1, VT( I,1 ), LDVT ) J = J + M*2 END DO - CALL SLASET( 'A', NS, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT) + CALL SLASET( 'A', M, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT ) * * Call SORMBR to compute (VB**T)*(PB**T) * (Workspace in WORK( ITEMP ): need M, prefer M*NB) @@ -794,7 +756,7 @@ CALL SCOPY( M, WORK( J ), 1, VT( I,1 ), LDVT ) J = J + M*2 END DO - CALL SLASET( 'A', NS, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT) + CALL SLASET( 'A', M, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT ) * * Call SORMBR to compute VB**T * PB**T * (Workspace in WORK( ITEMP ): need M, prefer M*NB) diff --git a/lapack-netlib/SRC/zgesvdx.f b/lapack-netlib/SRC/zgesvdx.f index c9509e4..6f7d5ba 100644 --- a/lapack-netlib/SRC/zgesvdx.f +++ b/lapack-netlib/SRC/zgesvdx.f @@ -36,30 +36,27 @@ * .. * * -*> \par Purpose: -* ============= -*> -*> \verbatim -*> -*> ZGESVDX computes the singular value decomposition (SVD) of a complex -*> M-by-N matrix A, optionally computing the left and/or right singular -*> vectors. The SVD is written -*> -*> A = U * SIGMA * transpose(V) -*> -*> where SIGMA is an M-by-N matrix which is zero except for its -*> min(m,n) diagonal elements, U is an M-by-M unitary matrix, and -*> V is an N-by-N unitary matrix. The diagonal elements of SIGMA -*> are the singular values of A; they are real and non-negative, and -*> are returned in descending order. The first min(m,n) columns of -*> U and V are the left and right singular vectors of A. -*> -*> ZGESVDX uses an eigenvalue problem for obtaining the SVD, which -*> allows for the computation of a subset of singular values and -*> vectors. See DBDSVDX for details. -*> -*> Note that the routine returns V**T, not V. -*> \endverbatim +* Purpose +* ======= +* +* ZGESVDX computes the singular value decomposition (SVD) of a complex +* M-by-N matrix A, optionally computing the left and/or right singular +* vectors. The SVD is written +* +* A = U * SIGMA * transpose(V) +* +* where SIGMA is an M-by-N matrix which is zero except for its +* min(m,n) diagonal elements, U is an M-by-M unitary matrix, and +* V is an N-by-N unitary matrix. The diagonal elements of SIGMA +* are the singular values of A; they are real and non-negative, and +* are returned in descending order. The first min(m,n) columns of +* U and V are the left and right singular vectors of A. +* +* ZGESVDX uses an eigenvalue problem for obtaining the SVD, which +* allows for the computation of a subset of singular values and +* vectors. See DBDSVDX for details. +* +* Note that the routine returns V**T, not V. * * Arguments: * ========== @@ -110,7 +107,7 @@ *> *> \param[in,out] A *> \verbatim -*> A is COMPLEX*16 array, dimension (LDA,N) +*> A is COMPLEX array, dimension (LDA,N) *> On entry, the M-by-N matrix A. *> On exit, the contents of A are destroyed. *> \endverbatim @@ -170,7 +167,7 @@ *> vectors, stored columnwise) as specified by RANGE; if *> JOBU = 'N', U is not referenced. *> Note: The user must ensure that UCOL >= NS; if RANGE = 'V', -*> the exact value of NS is not known in advance and an upper +*> the exact value of NS is not known ILQFin advance and an upper *> bound must be used. *> \endverbatim *> @@ -294,8 +291,8 @@ CHARACTER JOBZ, RNGTGK LOGICAL ALLS, INDS, LQUERY, VALS, WANTU, WANTVT INTEGER I, ID, IE, IERR, ILQF, ILTGK, IQRF, ISCL, - $ ITAU, ITAUP, ITAUQ, ITEMP, ITEMPR, ITGKZ, - $ IUTGK, J, K, MAXWRK, MINMN, MINWRK, MNTHR + $ ITAU, ITAUP, ITAUQ, ITEMP, ITGKZ, IUTGK, + $ J, K, MAXWRK, MINMN, MINWRK, MNTHR DOUBLE PRECISION ABSTOL, ANRM, BIGNUM, EPS, SMLNUM * .. * .. Local Arrays .. @@ -367,14 +364,8 @@ IF( INFO.EQ.0 ) THEN IF( WANTU .AND. LDU.LT.M ) THEN INFO = -15 - ELSE IF( WANTVT ) THEN - IF( INDS ) THEN - IF( LDVT.LT.IU-IL+1 ) THEN - INFO = -17 - END IF - ELSE IF( LDVT.LT.MINMN ) THEN - INFO = -17 - END IF + ELSE IF( WANTVT .AND. LDVT.LT.MINMN ) THEN + INFO = -16 END IF END IF END IF @@ -396,24 +387,18 @@ * * Path 1 (M much larger than N) * - MINWRK = N*(N+5) - MAXWRK = N + N*ILAENV(1,'ZGEQRF',' ',M,N,-1,-1) - MAXWRK = MAX(MAXWRK, - $ N*N+2*N+2*N*ILAENV(1,'ZGEBRD',' ',N,N,-1,-1)) - IF (WANTU .OR. WANTVT) THEN - MAXWRK = MAX(MAXWRK, - $ N*N+2*N+N*ILAENV(1,'ZUNMQR','LN',N,N,N,-1)) - END IF + MAXWRK = N + N* + $ ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) + MAXWRK = MAX( MAXWRK, N*N + N + 2*N* + $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) + MINWRK = N*(N+4) ELSE * * Path 2 (M at least N, but not much larger) * - MINWRK = 3*N + M - MAXWRK = 2*N + (M+N)*ILAENV(1,'ZGEBRD',' ',M,N,-1,-1) - IF (WANTU .OR. WANTVT) THEN - MAXWRK = MAX(MAXWRK, - $ 2*N+N*ILAENV(1,'ZUNMQR','LN',N,N,N,-1)) - END IF + MAXWRK = 2*N + ( M+N )* + $ ILAENV( 1, 'ZGEBRD', ' ', M, N, -1, -1 ) + MINWRK = 2*N + M END IF ELSE MNTHR = ILAENV( 6, 'ZGESVD', JOBU // JOBVT, M, N, 0, 0 ) @@ -421,25 +406,18 @@ * * Path 1t (N much larger than M) * - MINWRK = M*(M+5) - MAXWRK = M + M*ILAENV(1,'ZGELQF',' ',M,N,-1,-1) - MAXWRK = MAX(MAXWRK, - $ M*M+2*M+2*M*ILAENV(1,'ZGEBRD',' ',M,M,-1,-1)) - IF (WANTU .OR. WANTVT) THEN - MAXWRK = MAX(MAXWRK, - $ M*M+2*M+M*ILAENV(1,'ZUNMQR','LN',M,M,M,-1)) - END IF + MAXWRK = M + M* + $ ILAENV( 1, 'ZGELQF', ' ', M, N, -1, -1 ) + MAXWRK = MAX( MAXWRK, M*M + M + 2*M* + $ ILAENV( 1, 'ZGEBRD', ' ', M, M, -1, -1 ) ) + MINWRK = M*(M+4) ELSE * * Path 2t (N greater than M, but not much larger) * -* - MINWRK = 3*M + N - MAXWRK = 2*M + (M+N)*ILAENV(1,'ZGEBRD',' ',M,N,-1,-1) - IF (WANTU .OR. WANTVT) THEN - MAXWRK = MAX(MAXWRK, - $ 2*M+M*ILAENV(1,'ZUNMQR','LN',M,M,M,-1)) - END IF + MAXWRK = M*(M*2+19) + ( M+N )* + $ ILAENV( 1, 'ZGEBRD', ' ', M, N, -1, -1 ) + MINWRK = 2*M + N END IF END IF END IF @@ -537,14 +515,14 @@ CALL ZGEBRD( N, N, WORK( IQRF ), N, RWORK( ID ), $ RWORK( IE ), WORK( ITAUQ ), WORK( ITAUP ), $ WORK( ITEMP ), LWORK-ITEMP+1, INFO ) - ITEMPR = ITGKZ + N*(N*2+1) + ITEMP = ITGKZ + N*(N*2+1) * * Solve eigenvalue problem TGK*Z=Z*S. * (Workspace: need 2*N*N+14*N) * CALL DBDSVDX( 'U', JOBZ, RNGTGK, N, RWORK( ID ), $ RWORK( IE ), VL, VU, ILTGK, IUTGK, NS, S, - $ RWORK( ITGKZ ), N*2, RWORK( ITEMPR ), + $ RWORK( ITGKZ ), N*2, RWORK( ITEMP ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -558,7 +536,7 @@ END DO K = K + N END DO - CALL ZLASET( 'A', M-N, NS, CZERO, CZERO, U( N+1,1 ), LDU) + CALL ZLASET( 'A', M-N, N, CZERO, CZERO, U( N+1,1 ), LDU ) * * Call ZUNMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -613,14 +591,14 @@ CALL ZGEBRD( M, N, A, LDA, RWORK( ID ), RWORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), WORK( ITEMP ), $ LWORK-ITEMP+1, INFO ) - ITEMPR = ITGKZ + N*(N*2+1) + ITEMP = ITGKZ + N*(N*2+1) * * Solve eigenvalue problem TGK*Z=Z*S. * (Workspace: need 2*N*N+14*N) * CALL DBDSVDX( 'U', JOBZ, RNGTGK, N, RWORK( ID ), $ RWORK( IE ), VL, VU, ILTGK, IUTGK, NS, S, - $ RWORK( ITGKZ ), N*2, RWORK( ITEMPR ), + $ RWORK( ITGKZ ), N*2, RWORK( ITEMP ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -634,7 +612,7 @@ END DO K = K + N END DO - CALL ZLASET( 'A', M-N, NS, CZERO, CZERO, U( N+1,1 ), LDU) + CALL ZLASET( 'A', M-N, N, CZERO, CZERO, U( N+1,1 ), LDU ) * * Call ZUNMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -700,14 +678,14 @@ CALL ZGEBRD( M, M, WORK( ILQF ), M, RWORK( ID ), $ RWORK( IE ), WORK( ITAUQ ), WORK( ITAUP ), $ WORK( ITEMP ), LWORK-ITEMP+1, INFO ) - ITEMPR = ITGKZ + M*(M*2+1) + ITEMP = ITGKZ + M*(M*2+1) * * Solve eigenvalue problem TGK*Z=Z*S. * (Workspace: need 2*M*M+14*M) * CALL DBDSVDX( 'U', JOBZ, RNGTGK, M, RWORK( ID ), $ RWORK( IE ), VL, VU, ILTGK, IUTGK, NS, S, - $ RWORK( ITGKZ ), M*2, RWORK( ITEMPR ), + $ RWORK( ITGKZ ), M*2, RWORK( ITEMP ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -741,7 +719,7 @@ END DO K = K + M END DO - CALL ZLASET( 'A', NS, N-M, CZERO, CZERO, + CALL ZLASET( 'A', M, N-M, CZERO, CZERO, $ VT( 1,M+1 ), LDVT ) * * Call ZUNMBR to compute (VB**T)*(PB**T) @@ -777,14 +755,14 @@ CALL ZGEBRD( M, N, A, LDA, RWORK( ID ), RWORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), WORK( ITEMP ), $ LWORK-ITEMP+1, INFO ) - ITEMPR = ITGKZ + M*(M*2+1) + ITEMP = ITGKZ + M*(M*2+1) * * Solve eigenvalue problem TGK*Z=Z*S. * (Workspace: need 2*M*M+14*M) * CALL DBDSVDX( 'L', JOBZ, RNGTGK, M, RWORK( ID ), $ RWORK( IE ), VL, VU, ILTGK, IUTGK, NS, S, - $ RWORK( ITGKZ ), M*2, RWORK( ITEMPR ), + $ RWORK( ITGKZ ), M*2, RWORK( ITEMP ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -818,7 +796,7 @@ END DO K = K + M END DO - CALL ZLASET( 'A', NS, N-M, CZERO, CZERO, + CALL ZLASET( 'A', M, N-M, CZERO, CZERO, $ VT( 1,M+1 ), LDVT ) * * Call ZUNMBR to compute VB**T * PB**T -- 2.7.4