diff --git a/lapack-netlib/LAPACKE/src/lapacke_clantr.c b/lapack-netlib/LAPACKE/src/lapacke_clantr.c index 4ba436753..33e6e57ff 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_clantr.c +++ b/lapack-netlib/LAPACKE/src/lapacke_clantr.c @@ -51,8 +51,7 @@ float LAPACKE_clantr( int matrix_layout, char norm, char uplo, char diag, } #endif /* Allocate memory for working array(s) */ - if( LAPACKE_lsame( norm, 'i' ) || LAPACKE_lsame( norm, '1' ) || - LAPACKE_lsame( norm, 'O' ) ) { + if( LAPACKE_lsame( norm, 'i' ) ) { work = (float*)LAPACKE_malloc( sizeof(float) * MAX(1,MAX(m,n)) ); if( work == NULL ) { info = LAPACK_WORK_MEMORY_ERROR; @@ -63,8 +62,7 @@ float LAPACKE_clantr( int matrix_layout, char norm, char uplo, char diag, res = LAPACKE_clantr_work( matrix_layout, norm, uplo, diag, m, n, a, lda, work ); /* Release memory and exit */ - if( LAPACKE_lsame( norm, 'i' ) || LAPACKE_lsame( norm, '1' ) || - LAPACKE_lsame( norm, 'O' ) ) { + if( LAPACKE_lsame( norm, 'i' ) ) { LAPACKE_free( work ); } exit_level_0: diff --git a/lapack-netlib/LAPACKE/src/lapacke_dlantr.c b/lapack-netlib/LAPACKE/src/lapacke_dlantr.c index 65802f3e3..8fd112084 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_dlantr.c +++ b/lapack-netlib/LAPACKE/src/lapacke_dlantr.c @@ -51,8 +51,7 @@ double LAPACKE_dlantr( int matrix_layout, char norm, char uplo, char diag, } #endif /* Allocate memory for working array(s) */ - if( LAPACKE_lsame( norm, 'i' ) || LAPACKE_lsame( norm, '1' ) || - LAPACKE_lsame( norm, 'O' ) ) { + if( LAPACKE_lsame( norm, 'i' ) ) { work = (double*)LAPACKE_malloc( sizeof(double) * MAX(1,MAX(m,n)) ); if( work == NULL ) { info = LAPACK_WORK_MEMORY_ERROR; @@ -63,8 +62,7 @@ double LAPACKE_dlantr( int matrix_layout, char norm, char uplo, char diag, res = LAPACKE_dlantr_work( matrix_layout, norm, uplo, diag, m, n, a, lda, work ); /* Release memory and exit */ - if( LAPACKE_lsame( norm, 'i' ) || LAPACKE_lsame( norm, '1' ) || - LAPACKE_lsame( norm, 'O' ) ) { + if( LAPACKE_lsame( norm, 'i' ) ) { LAPACKE_free( work ); } exit_level_0: diff --git a/lapack-netlib/LAPACKE/src/lapacke_dlantr_work.c b/lapack-netlib/LAPACKE/src/lapacke_dlantr_work.c index 59eef3801..ee5e665ee 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_dlantr_work.c +++ b/lapack-netlib/LAPACKE/src/lapacke_dlantr_work.c @@ -38,10 +38,10 @@ double LAPACKE_dlantr_work( int matrix_layout, char norm, char uplo, const double* a, lapack_int lda, double* work ) { lapack_int info = 0; - double res = 0.; + double res = 0.; if( matrix_layout == LAPACK_COL_MAJOR ) { /* Call LAPACK function and adjust info */ - LAPACK_dlantr( &norm, &uplo, &diag, &m, &n, a, &lda, work ); + res = LAPACK_dlantr( &norm, &uplo, &diag, &m, &n, a, &lda, work ); if( info < 0 ) { info = info - 1; } diff --git a/lapack-netlib/LAPACKE/src/lapacke_dormbr_work.c b/lapack-netlib/LAPACKE/src/lapacke_dormbr_work.c index 9db92ce98..dcd8842fa 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_dormbr_work.c +++ b/lapack-netlib/LAPACKE/src/lapacke_dormbr_work.c @@ -74,11 +74,10 @@ lapack_int LAPACKE_dormbr_work( int matrix_layout, char vect, char side, } /* Allocate memory for temporary array(s) */ if( LAPACKE_lsame( vect, 'q' ) ) { - a_t = (double*)LAPACKE_malloc( sizeof(double) * lda_t * k ); + a_t = (double*)LAPACKE_malloc( sizeof(double) * lda_t * MAX(1,k) ); } else { - a_t = (double*)LAPACKE_malloc( sizeof(double) * lda_t * nq ); + a_t = (double*)LAPACKE_malloc( sizeof(double) * lda_t * MAX(1,nq) ); } - if( a_t == NULL ) { info = LAPACK_TRANSPOSE_MEMORY_ERROR; goto exit_level_0; @@ -89,11 +88,7 @@ lapack_int LAPACKE_dormbr_work( int matrix_layout, char vect, char side, goto exit_level_1; } /* Transpose input matrices */ - if( LAPACKE_lsame( vect, 'q' ) ) { - LAPACKE_dge_trans( matrix_layout, nq, k, a, lda, a_t, lda_t ); - } else { - LAPACKE_dge_trans( matrix_layout, k, nq, a, lda, a_t, lda_t ); - } + LAPACKE_dge_trans( matrix_layout, r, MIN(nq,k), a, lda, a_t, lda_t ); LAPACKE_dge_trans( matrix_layout, m, n, c, ldc, c_t, ldc_t ); /* Call LAPACK function and adjust info */ LAPACK_dormbr( &vect, &side, &trans, &m, &n, &k, a_t, &lda_t, tau, c_t, diff --git a/lapack-netlib/LAPACKE/src/lapacke_dormlq_work.c b/lapack-netlib/LAPACKE/src/lapacke_dormlq_work.c index 2a59cd56a..f46c6d3b1 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_dormlq_work.c +++ b/lapack-netlib/LAPACKE/src/lapacke_dormlq_work.c @@ -87,12 +87,7 @@ lapack_int LAPACKE_dormlq_work( int matrix_layout, char side, char trans, goto exit_level_1; } /* Transpose input matrices */ - if( LAPACKE_lsame( side, 'l' ) ){ - LAPACKE_dge_trans( matrix_layout, k, m, a, lda, a_t, lda_t ); - } else { - LAPACKE_dge_trans( matrix_layout, k, n, a, lda, a_t, lda_t ); - } - + LAPACKE_dge_trans( matrix_layout, k, m, a, lda, a_t, lda_t ); LAPACKE_dge_trans( matrix_layout, m, n, c, ldc, c_t, ldc_t ); /* Call LAPACK function and adjust info */ LAPACK_dormlq( &side, &trans, &m, &n, &k, a_t, &lda_t, tau, c_t, &ldc_t, diff --git a/lapack-netlib/LAPACKE/src/lapacke_slantr.c b/lapack-netlib/LAPACKE/src/lapacke_slantr.c index 99942b352..23b19b15d 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_slantr.c +++ b/lapack-netlib/LAPACKE/src/lapacke_slantr.c @@ -51,8 +51,7 @@ float LAPACKE_slantr( int matrix_layout, char norm, char uplo, char diag, } #endif /* Allocate memory for working array(s) */ - if( LAPACKE_lsame( norm, 'i' ) || LAPACKE_lsame( norm, '1' ) || - LAPACKE_lsame( norm, 'O' ) ) { + if( LAPACKE_lsame( norm, 'i' ) ) { work = (float*)LAPACKE_malloc( sizeof(float) * MAX(1,MAX(m,n)) ); if( work == NULL ) { info = LAPACK_WORK_MEMORY_ERROR; @@ -63,8 +62,7 @@ float LAPACKE_slantr( int matrix_layout, char norm, char uplo, char diag, res = LAPACKE_slantr_work( matrix_layout, norm, uplo, diag, m, n, a, lda, work ); /* Release memory and exit */ - if( LAPACKE_lsame( norm, 'i' ) || LAPACKE_lsame( norm, '1' ) || - LAPACKE_lsame( norm, 'O' ) ) { + if( LAPACKE_lsame( norm, 'i' ) ) { LAPACKE_free( work ); } exit_level_0: diff --git a/lapack-netlib/LAPACKE/src/lapacke_slantr_work.c b/lapack-netlib/LAPACKE/src/lapacke_slantr_work.c index 79c71a00f..92a0e4017 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_slantr_work.c +++ b/lapack-netlib/LAPACKE/src/lapacke_slantr_work.c @@ -41,7 +41,7 @@ float LAPACKE_slantr_work( int matrix_layout, char norm, char uplo, float res = 0.; if( matrix_layout == LAPACK_COL_MAJOR ) { /* Call LAPACK function and adjust info */ - LAPACK_slantr( &norm, &uplo, &diag, &m, &n, a, &lda, work ); + res = LAPACK_slantr( &norm, &uplo, &diag, &m, &n, a, &lda, work ); if( info < 0 ) { info = info - 1; } diff --git a/lapack-netlib/LAPACKE/src/lapacke_sormbr_work.c b/lapack-netlib/LAPACKE/src/lapacke_sormbr_work.c index f7a2ff3d0..34bc41f30 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_sormbr_work.c +++ b/lapack-netlib/LAPACKE/src/lapacke_sormbr_work.c @@ -73,8 +73,11 @@ lapack_int LAPACKE_sormbr_work( int matrix_layout, char vect, char side, return (info < 0) ? (info - 1) : info; } /* Allocate memory for temporary array(s) */ - a_t = (float*) - LAPACKE_malloc( sizeof(float) * lda_t * MAX(1,MIN(nq,k)) ); + if( LAPACKE_lsame( vect, 'q' ) ) { + a_t = (float*)LAPACKE_malloc( sizeof(float) * lda_t * MAX(1,k) ); + } else { + a_t = (float*)LAPACKE_malloc( sizeof(float) * lda_t * MAX(1,nq) ); + } if( a_t == NULL ) { info = LAPACK_TRANSPOSE_MEMORY_ERROR; goto exit_level_0; diff --git a/lapack-netlib/LAPACKE/src/lapacke_sormlq_work.c b/lapack-netlib/LAPACKE/src/lapacke_sormlq_work.c index a277436cd..b02a2d100 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_sormlq_work.c +++ b/lapack-netlib/LAPACKE/src/lapacke_sormlq_work.c @@ -72,7 +72,11 @@ lapack_int LAPACKE_sormlq_work( int matrix_layout, char side, char trans, return (info < 0) ? (info - 1) : info; } /* Allocate memory for temporary array(s) */ - a_t = (float*)LAPACKE_malloc( sizeof(float) * lda_t * MAX(1,m) ); + if( LAPACKE_lsame( side, 'l' ) ) { + a_t = (float*)LAPACKE_malloc( sizeof(float) * lda_t * MAX(1,m) ); + } else { + a_t = (float*)LAPACKE_malloc( sizeof(float) * lda_t * MAX(1,n) ); + } if( a_t == NULL ) { info = LAPACK_TRANSPOSE_MEMORY_ERROR; goto exit_level_0; diff --git a/lapack-netlib/LAPACKE/src/lapacke_zlantr.c b/lapack-netlib/LAPACKE/src/lapacke_zlantr.c index b2c637101..2b645e750 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_zlantr.c +++ b/lapack-netlib/LAPACKE/src/lapacke_zlantr.c @@ -51,8 +51,7 @@ double LAPACKE_zlantr( int matrix_layout, char norm, char uplo, char diag, } #endif /* Allocate memory for working array(s) */ - if( LAPACKE_lsame( norm, 'i' ) || LAPACKE_lsame( norm, '1' ) || - LAPACKE_lsame( norm, 'O' ) ) { + if( LAPACKE_lsame( norm, 'i' ) ) { work = (double*)LAPACKE_malloc( sizeof(double) * MAX(1,MAX(m,n)) ); if( work == NULL ) { info = LAPACK_WORK_MEMORY_ERROR; @@ -63,8 +62,7 @@ double LAPACKE_zlantr( int matrix_layout, char norm, char uplo, char diag, res = LAPACKE_zlantr_work( matrix_layout, norm, uplo, diag, m, n, a, lda, work ); /* Release memory and exit */ - if( LAPACKE_lsame( norm, 'i' ) || LAPACKE_lsame( norm, '1' ) || - LAPACKE_lsame( norm, 'O' ) ) { + if( LAPACKE_lsame( norm, 'i' ) ) { LAPACKE_free( work ); } exit_level_0: diff --git a/lapack-netlib/LAPACKE/src/lapacke_zlantr_work.c b/lapack-netlib/LAPACKE/src/lapacke_zlantr_work.c index 5fd9f8442..0988bf6e8 100644 --- a/lapack-netlib/LAPACKE/src/lapacke_zlantr_work.c +++ b/lapack-netlib/LAPACKE/src/lapacke_zlantr_work.c @@ -39,7 +39,7 @@ double LAPACKE_zlantr_work( int matrix_layout, char norm, char uplo, double* work ) { lapack_int info = 0; - double res = 0.; + double res = 0.; if( matrix_layout == LAPACK_COL_MAJOR ) { /* Call LAPACK function and adjust info */ res = LAPACK_zlantr( &norm, &uplo, &diag, &m, &n, a, &lda, work ); diff --git a/lapack-netlib/SRC/cgeev.f b/lapack-netlib/SRC/cgeev.f index b79b64544..0f48322a8 100644 --- a/lapack-netlib/SRC/cgeev.f +++ b/lapack-netlib/SRC/cgeev.f @@ -405,9 +405,9 @@ $ WORK( IWRK ), LWORK-IWRK+1, INFO ) END IF * -* If INFO > 0 from CHSEQR, then quit +* If INFO .NE. 0 from CHSEQR, then quit * - IF( INFO.GT.0 ) + IF( INFO.NE.0 ) $ GO TO 50 * IF( WANTVL .OR. WANTVR ) THEN diff --git a/lapack-netlib/SRC/cgesvdx.f b/lapack-netlib/SRC/cgesvdx.f index 235426ad4..87ea9861d 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 ILQFin advance and an upper +*> the exact value of NS is not known in 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, ITGKZ, IUTGK, - $ J, K, MAXWRK, MINMN, MINWRK, MNTHR + $ ITAU, ITAUP, ITAUQ, ITEMP, ITEMPR, ITGKZ, + $ IUTGK, J, K, MAXWRK, MINMN, MINWRK, MNTHR REAL ABSTOL, ANRM, BIGNUM, EPS, SMLNUM * .. * .. Local Arrays .. @@ -367,8 +367,14 @@ IF( INFO.EQ.0 ) THEN IF( WANTU .AND. LDU.LT.M ) THEN INFO = -15 - ELSE IF( WANTVT .AND. LDVT.LT.MINMN ) THEN - INFO = -16 + 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 END IF END IF END IF @@ -390,18 +396,24 @@ * * Path 1 (M much larger than N) * - 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) + 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 ELSE * * Path 2 (M at least N, but not much larger) * - MAXWRK = 2*N + ( M+N )* - $ ILAENV( 1, 'CGEBRD', ' ', M, N, -1, -1 ) - MINWRK = 2*N + M + 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 END IF ELSE MNTHR = ILAENV( 6, 'CGESVD', JOBU // JOBVT, M, N, 0, 0 ) @@ -409,18 +421,25 @@ * * Path 1t (N much larger than M) * - 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) + 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 ELSE * * Path 2t (N greater than M, but not much larger) * - MAXWRK = M*(M*2+19) + ( M+N )* - $ ILAENV( 1, 'CGEBRD', ' ', M, N, -1, -1 ) - MINWRK = 2*M + N +* + 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 END IF END IF END IF @@ -518,14 +537,14 @@ CALL CGEBRD( N, N, WORK( IQRF ), N, RWORK( ID ), $ RWORK( IE ), WORK( ITAUQ ), WORK( ITAUP ), $ WORK( ITEMP ), LWORK-ITEMP+1, INFO ) - ITEMP = ITGKZ + N*(N*2+1) + ITEMPR = 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( ITEMP ), + $ RWORK( ITGKZ ), N*2, RWORK( ITEMPR ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -539,7 +558,7 @@ END DO K = K + N END DO - CALL CLASET( 'A', M-N, N, CZERO, CZERO, U( N+1,1 ), LDU ) + CALL CLASET( 'A', M-N, NS, CZERO, CZERO, U( N+1,1 ), LDU) * * Call CUNMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -594,14 +613,14 @@ CALL CGEBRD( M, N, A, LDA, RWORK( ID ), RWORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), WORK( ITEMP ), $ LWORK-ITEMP+1, INFO ) - ITEMP = ITGKZ + N*(N*2+1) + ITEMPR = 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( ITEMP ), + $ RWORK( ITGKZ ), N*2, RWORK( ITEMPR ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -615,7 +634,7 @@ END DO K = K + N END DO - CALL CLASET( 'A', M-N, N, CZERO, CZERO, U( N+1,1 ), LDU ) + CALL CLASET( 'A', M-N, NS, CZERO, CZERO, U( N+1,1 ), LDU) * * Call CUNMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -681,14 +700,14 @@ CALL CGEBRD( M, M, WORK( ILQF ), M, RWORK( ID ), $ RWORK( IE ), WORK( ITAUQ ), WORK( ITAUP ), $ WORK( ITEMP ), LWORK-ITEMP+1, INFO ) - ITEMP = ITGKZ + M*(M*2+1) + ITEMPR = 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( ITEMP ), + $ RWORK( ITGKZ ), M*2, RWORK( ITEMPR ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -722,7 +741,7 @@ END DO K = K + M END DO - CALL CLASET( 'A', M, N-M, CZERO, CZERO, + CALL CLASET( 'A', NS, N-M, CZERO, CZERO, $ VT( 1,M+1 ), LDVT ) * * Call CUNMBR to compute (VB**T)*(PB**T) @@ -758,14 +777,14 @@ CALL CGEBRD( M, N, A, LDA, RWORK( ID ), RWORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), WORK( ITEMP ), $ LWORK-ITEMP+1, INFO ) - ITEMP = ITGKZ + M*(M*2+1) + ITEMPR = 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( ITEMP ), + $ RWORK( ITGKZ ), M*2, RWORK( ITEMPR ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -799,7 +818,7 @@ END DO K = K + M END DO - CALL CLASET( 'A', M, N-M, CZERO, CZERO, + CALL CLASET( 'A', NS, N-M, CZERO, CZERO, $ VT( 1,M+1 ), LDVT ) * * Call CUNMBR to compute VB**T * PB**T diff --git a/lapack-netlib/SRC/cgetc2.f b/lapack-netlib/SRC/cgetc2.f index fac6b56820..99eb69d92 100644 --- a/lapack-netlib/SRC/cgetc2.f +++ b/lapack-netlib/SRC/cgetc2.f @@ -145,15 +145,33 @@ INTRINSIC ABS, CMPLX, MAX * .. * .. Executable Statements .. +* + INFO = 0 +* +* Quick return if possible +* + IF( N.EQ.0 ) + $ RETURN * * Set constants to control overflow * - INFO = 0 EPS = SLAMCH( 'P' ) SMLNUM = SLAMCH( 'S' ) / EPS BIGNUM = ONE / SMLNUM CALL SLABAD( SMLNUM, BIGNUM ) * +* Handle the case N=1 by itself +* + IF( N.EQ.1 ) THEN + IPIV( 1 ) = 1 + JPIV( 1 ) = 1 + IF( ABS( A( 1, 1 ) ).LT.SMLNUM ) THEN + INFO = 1 + A( 1, 1 ) = CMPLX( SMLNUM, ZERO ) + END IF + RETURN + END IF +* * Factorize A using complete pivoting. * Set pivots less than SMIN to SMIN * diff --git a/lapack-netlib/SRC/cggev3.f b/lapack-netlib/SRC/cggev3.f index 4a000fe10..decdae509 100644 --- a/lapack-netlib/SRC/cggev3.f +++ b/lapack-netlib/SRC/cggev3.f @@ -339,16 +339,16 @@ $ LDVL, VR, LDVR, WORK, -1, IERR ) LWKOPT = MAX( LWKOPT, N+INT( WORK( 1 ) ) ) CALL CHGEQZ( 'S', JOBVL, JOBVR, N, 1, N, A, LDA, B, LDB, - $ ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, - $ -1, WORK, IERR ) + $ ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, -1, + $ RWORK, IERR ) LWKOPT = MAX( LWKOPT, N+INT( WORK( 1 ) ) ) ELSE CALL CGGHD3( 'N', 'N', N, 1, N, A, LDA, B, LDB, VL, LDVL, $ VR, LDVR, WORK, -1, IERR ) LWKOPT = MAX( LWKOPT, N+INT( WORK( 1 ) ) ) CALL CHGEQZ( 'E', JOBVL, JOBVR, N, 1, N, A, LDA, B, LDB, - $ ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, - $ -1, WORK, IERR ) + $ ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, -1, + $ RWORK, IERR ) LWKOPT = MAX( LWKOPT, N+INT( WORK( 1 ) ) ) END IF WORK( 1 ) = CMPLX( LWKOPT ) diff --git a/lapack-netlib/SRC/dgeev.f b/lapack-netlib/SRC/dgeev.f index dd60db69e..328eaa39c 100644 --- a/lapack-netlib/SRC/dgeev.f +++ b/lapack-netlib/SRC/dgeev.f @@ -418,9 +418,9 @@ $ WORK( IWRK ), LWORK-IWRK+1, INFO ) END IF * -* If INFO > 0 from DHSEQR, then quit +* If INFO .NE. 0 from DHSEQR, then quit * - IF( INFO.GT.0 ) + IF( INFO.NE.0 ) $ GO TO 50 * IF( WANTVL .OR. WANTVR ) THEN diff --git a/lapack-netlib/SRC/dgesvdx.f b/lapack-netlib/SRC/dgesvdx.f index cfa2ff05d..4588083f8 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 ILQFin advance and an upper +*> the exact value of NS is not known in advance and an upper *> bound must be used. *> \endverbatim *> @@ -357,8 +357,14 @@ IF( INFO.EQ.0 ) THEN IF( WANTU .AND. LDU.LT.M ) THEN INFO = -15 - ELSE IF( WANTVT .AND. LDVT.LT.MINMN ) THEN - INFO = -16 + 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 END IF END IF END IF @@ -380,18 +386,34 @@ * * Path 1 (M much larger than N) * - MAXWRK = N*(N*2+16) + + MAXWRK = N + $ N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - MAXWRK = MAX( MAXWRK, N*(N*2+20) + 2*N* + MAXWRK = MAX( MAXWRK, N*(N+5) + 2*N* $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - MINWRK = N*(N*2+21) + 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) ELSE * * Path 2 (M at least N, but not much larger) * - MAXWRK = N*(N*2+19) + ( M+N )* + MAXWRK = 4*N + ( M+N )* $ ILAENV( 1, 'DGEBRD', ' ', M, N, -1, -1 ) - MINWRK = N*(N*2+20) + M + 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) END IF ELSE MNTHR = ILAENV( 6, 'DGESVD', JOBU // JOBVT, M, N, 0, 0 ) @@ -399,18 +421,34 @@ * * Path 1t (N much larger than M) * - MAXWRK = M*(M*2+16) + + MAXWRK = M + $ M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - MAXWRK = MAX( MAXWRK, M*(M*2+20) + 2*M* + MAXWRK = MAX( MAXWRK, M*(M+5) + 2*M* $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - MINWRK = M*(M*2+21) + 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) ELSE * -* Path 2t (N greater than M, but not much larger) +* Path 2t (N at least M, but not much larger) * - MAXWRK = M*(M*2+19) + ( M+N )* + MAXWRK = 4*M + ( M+N )* $ ILAENV( 1, 'DGEBRD', ' ', M, N, -1, -1 ) - MINWRK = M*(M*2+20) + N + 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) END IF END IF END IF @@ -522,7 +560,7 @@ CALL DCOPY( N, WORK( J ), 1, U( 1,I ), 1 ) J = J + N*2 END DO - CALL DLASET( 'A', M-N, N, ZERO, ZERO, U( N+1,1 ), LDU ) + CALL DLASET( 'A', M-N, NS, ZERO, ZERO, U( N+1,1 ), LDU ) * * Call DORMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -591,7 +629,7 @@ CALL DCOPY( N, WORK( J ), 1, U( 1,I ), 1 ) J = J + N*2 END DO - CALL DLASET( 'A', M-N, N, ZERO, ZERO, U( N+1,1 ), LDU ) + CALL DLASET( 'A', M-N, NS, ZERO, ZERO, U( N+1,1 ), LDU ) * * Call DORMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -687,7 +725,7 @@ CALL DCOPY( M, WORK( J ), 1, VT( I,1 ), LDVT ) J = J + M*2 END DO - CALL DLASET( 'A', M, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT ) + CALL DLASET( 'A', NS, 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) @@ -756,7 +794,7 @@ CALL DCOPY( M, WORK( J ), 1, VT( I,1 ), LDVT ) J = J + M*2 END DO - CALL DLASET( 'A', M, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT ) + CALL DLASET( 'A', NS, 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/dgetc2.f b/lapack-netlib/SRC/dgetc2.f index 7e43a0236..3cd7eeb2b 100644 --- a/lapack-netlib/SRC/dgetc2.f +++ b/lapack-netlib/SRC/dgetc2.f @@ -145,15 +145,33 @@ INTRINSIC ABS, MAX * .. * .. Executable Statements .. +* + INFO = 0 +* +* Quick return if possible +* + IF( N.EQ.0 ) + $ RETURN * * Set constants to control overflow * - INFO = 0 EPS = DLAMCH( 'P' ) SMLNUM = DLAMCH( 'S' ) / EPS BIGNUM = ONE / SMLNUM CALL DLABAD( SMLNUM, BIGNUM ) * +* Handle the case N=1 by itself +* + IF( N.EQ.1 ) THEN + IPIV( 1 ) = 1 + JPIV( 1 ) = 1 + IF( ABS( A( 1, 1 ) ).LT.SMLNUM ) THEN + INFO = 1 + A( 1, 1 ) = SMLNUM + END IF + RETURN + END IF +* * Factorize A using complete pivoting. * Set pivots less than SMIN to SMIN. * diff --git a/lapack-netlib/SRC/sgeev.f b/lapack-netlib/SRC/sgeev.f index 89dbe08c8..667de0afe 100644 --- a/lapack-netlib/SRC/sgeev.f +++ b/lapack-netlib/SRC/sgeev.f @@ -418,9 +418,9 @@ $ WORK( IWRK ), LWORK-IWRK+1, INFO ) END IF * -* If INFO > 0 from SHSEQR, then quit +* If INFO .NE. 0 from SHSEQR, then quit * - IF( INFO.GT.0 ) + IF( INFO.NE.0 ) $ GO TO 50 * IF( WANTVL .OR. WANTVR ) THEN diff --git a/lapack-netlib/SRC/sgesvdx.f b/lapack-netlib/SRC/sgesvdx.f index aae8b0764..9128a7c0a 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 ILQFin advance and an upper +*> the exact value of NS is not known in advance and an upper *> bound must be used. *> \endverbatim *> @@ -357,8 +357,14 @@ IF( INFO.EQ.0 ) THEN IF( WANTU .AND. LDU.LT.M ) THEN INFO = -15 - ELSE IF( WANTVT .AND. LDVT.LT.MINMN ) THEN - INFO = -16 + 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 END IF END IF END IF @@ -380,18 +386,34 @@ * * Path 1 (M much larger than N) * - MAXWRK = N*(N*2+16) + + MAXWRK = N + $ N*ILAENV( 1, 'SGEQRF', ' ', M, N, -1, -1 ) - MAXWRK = MAX( MAXWRK, N*(N*2+20) + 2*N* + MAXWRK = MAX( MAXWRK, N*(N+5) + 2*N* $ ILAENV( 1, 'SGEBRD', ' ', N, N, -1, -1 ) ) - MINWRK = N*(N*2+21) + 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) ELSE * * Path 2 (M at least N, but not much larger) * - MAXWRK = N*(N*2+19) + ( M+N )* + MAXWRK = 4*N + ( M+N )* $ ILAENV( 1, 'SGEBRD', ' ', M, N, -1, -1 ) - MINWRK = N*(N*2+20) + M + 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) END IF ELSE MNTHR = ILAENV( 6, 'SGESVD', JOBU // JOBVT, M, N, 0, 0 ) @@ -399,18 +421,34 @@ * * Path 1t (N much larger than M) * - MAXWRK = M*(M*2+16) + + MAXWRK = M + $ M*ILAENV( 1, 'SGELQF', ' ', M, N, -1, -1 ) - MAXWRK = MAX( MAXWRK, M*(M*2+20) + 2*M* + MAXWRK = MAX( MAXWRK, M*(M+5) + 2*M* $ ILAENV( 1, 'SGEBRD', ' ', M, M, -1, -1 ) ) - MINWRK = M*(M*2+21) + 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) ELSE * -* Path 2t (N greater than M, but not much larger) +* Path 2t (N at least M, but not much larger) * - MAXWRK = M*(M*2+19) + ( M+N )* + MAXWRK = 4*M + ( M+N )* $ ILAENV( 1, 'SGEBRD', ' ', M, N, -1, -1 ) - MINWRK = M*(M*2+20) + N + 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) END IF END IF END IF @@ -522,7 +560,7 @@ CALL SCOPY( N, WORK( J ), 1, U( 1,I ), 1 ) J = J + N*2 END DO - CALL SLASET( 'A', M-N, N, ZERO, ZERO, U( N+1,1 ), LDU ) + CALL SLASET( 'A', M-N, NS, ZERO, ZERO, U( N+1,1 ), LDU ) * * Call SORMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -591,7 +629,7 @@ CALL SCOPY( N, WORK( J ), 1, U( 1,I ), 1 ) J = J + N*2 END DO - CALL SLASET( 'A', M-N, N, ZERO, ZERO, U( N+1,1 ), LDU ) + CALL SLASET( 'A', M-N, NS, ZERO, ZERO, U( N+1,1 ), LDU ) * * Call SORMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -687,7 +725,7 @@ CALL SCOPY( M, WORK( J ), 1, VT( I,1 ), LDVT ) J = J + M*2 END DO - CALL SLASET( 'A', M, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT ) + CALL SLASET( 'A', NS, 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) @@ -756,7 +794,7 @@ CALL SCOPY( M, WORK( J ), 1, VT( I,1 ), LDVT ) J = J + M*2 END DO - CALL SLASET( 'A', M, N-M, ZERO, ZERO, VT( 1,M+1 ), LDVT ) + CALL SLASET( 'A', NS, 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/sgetc2.f b/lapack-netlib/SRC/sgetc2.f index 3c3880d4e..598446519 100644 --- a/lapack-netlib/SRC/sgetc2.f +++ b/lapack-netlib/SRC/sgetc2.f @@ -145,15 +145,33 @@ INTRINSIC ABS, MAX * .. * .. Executable Statements .. +* + INFO = 0 +* +* Quick return if possible +* + IF( N.EQ.0 ) + $ RETURN * * Set constants to control overflow * - INFO = 0 EPS = SLAMCH( 'P' ) SMLNUM = SLAMCH( 'S' ) / EPS BIGNUM = ONE / SMLNUM CALL SLABAD( SMLNUM, BIGNUM ) * +* Handle the case N=1 by itself +* + IF( N.EQ.1 ) THEN + IPIV( 1 ) = 1 + JPIV( 1 ) = 1 + IF( ABS( A( 1, 1 ) ).LT.SMLNUM ) THEN + INFO = 1 + A( 1, 1 ) = SMLNUM + END IF + RETURN + END IF +* * Factorize A using complete pivoting. * Set pivots less than SMIN to SMIN. * diff --git a/lapack-netlib/SRC/zgeev.f b/lapack-netlib/SRC/zgeev.f index d4520805f..a518b4cd9 100644 --- a/lapack-netlib/SRC/zgeev.f +++ b/lapack-netlib/SRC/zgeev.f @@ -404,9 +404,9 @@ $ WORK( IWRK ), LWORK-IWRK+1, INFO ) END IF * -* If INFO > 0 from ZHSEQR, then quit +* If INFO .NE. 0 from ZHSEQR, then quit * - IF( INFO.GT.0 ) + IF( INFO.NE.0 ) $ GO TO 50 * IF( WANTVL .OR. WANTVR ) THEN diff --git a/lapack-netlib/SRC/zgesvdx.f b/lapack-netlib/SRC/zgesvdx.f index 6f7d5ba04..c9509e458 100644 --- a/lapack-netlib/SRC/zgesvdx.f +++ b/lapack-netlib/SRC/zgesvdx.f @@ -36,27 +36,30 @@ * .. * * -* 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. +*> \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 * * Arguments: * ========== @@ -107,7 +110,7 @@ *> *> \param[in,out] A *> \verbatim -*> A is COMPLEX array, dimension (LDA,N) +*> A is COMPLEX*16 array, dimension (LDA,N) *> On entry, the M-by-N matrix A. *> On exit, the contents of A are destroyed. *> \endverbatim @@ -167,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 ILQFin advance and an upper +*> the exact value of NS is not known in advance and an upper *> bound must be used. *> \endverbatim *> @@ -291,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, ITGKZ, IUTGK, - $ J, K, MAXWRK, MINMN, MINWRK, MNTHR + $ ITAU, ITAUP, ITAUQ, ITEMP, ITEMPR, ITGKZ, + $ IUTGK, J, K, MAXWRK, MINMN, MINWRK, MNTHR DOUBLE PRECISION ABSTOL, ANRM, BIGNUM, EPS, SMLNUM * .. * .. Local Arrays .. @@ -364,8 +367,14 @@ IF( INFO.EQ.0 ) THEN IF( WANTU .AND. LDU.LT.M ) THEN INFO = -15 - ELSE IF( WANTVT .AND. LDVT.LT.MINMN ) THEN - INFO = -16 + 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 END IF END IF END IF @@ -387,18 +396,24 @@ * * Path 1 (M much larger than N) * - 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) + 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 ELSE * * Path 2 (M at least N, but not much larger) * - MAXWRK = 2*N + ( M+N )* - $ ILAENV( 1, 'ZGEBRD', ' ', M, N, -1, -1 ) - MINWRK = 2*N + M + 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 END IF ELSE MNTHR = ILAENV( 6, 'ZGESVD', JOBU // JOBVT, M, N, 0, 0 ) @@ -406,18 +421,25 @@ * * Path 1t (N much larger than M) * - 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) + 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 ELSE * * Path 2t (N greater than M, but not much larger) * - MAXWRK = M*(M*2+19) + ( M+N )* - $ ILAENV( 1, 'ZGEBRD', ' ', M, N, -1, -1 ) - MINWRK = 2*M + N +* + 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 END IF END IF END IF @@ -515,14 +537,14 @@ CALL ZGEBRD( N, N, WORK( IQRF ), N, RWORK( ID ), $ RWORK( IE ), WORK( ITAUQ ), WORK( ITAUP ), $ WORK( ITEMP ), LWORK-ITEMP+1, INFO ) - ITEMP = ITGKZ + N*(N*2+1) + ITEMPR = 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( ITEMP ), + $ RWORK( ITGKZ ), N*2, RWORK( ITEMPR ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -536,7 +558,7 @@ END DO K = K + N END DO - CALL ZLASET( 'A', M-N, N, CZERO, CZERO, U( N+1,1 ), LDU ) + CALL ZLASET( 'A', M-N, NS, CZERO, CZERO, U( N+1,1 ), LDU) * * Call ZUNMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -591,14 +613,14 @@ CALL ZGEBRD( M, N, A, LDA, RWORK( ID ), RWORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), WORK( ITEMP ), $ LWORK-ITEMP+1, INFO ) - ITEMP = ITGKZ + N*(N*2+1) + ITEMPR = 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( ITEMP ), + $ RWORK( ITGKZ ), N*2, RWORK( ITEMPR ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -612,7 +634,7 @@ END DO K = K + N END DO - CALL ZLASET( 'A', M-N, N, CZERO, CZERO, U( N+1,1 ), LDU ) + CALL ZLASET( 'A', M-N, NS, CZERO, CZERO, U( N+1,1 ), LDU) * * Call ZUNMBR to compute QB*UB. * (Workspace in WORK( ITEMP ): need N, prefer N*NB) @@ -678,14 +700,14 @@ CALL ZGEBRD( M, M, WORK( ILQF ), M, RWORK( ID ), $ RWORK( IE ), WORK( ITAUQ ), WORK( ITAUP ), $ WORK( ITEMP ), LWORK-ITEMP+1, INFO ) - ITEMP = ITGKZ + M*(M*2+1) + ITEMPR = 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( ITEMP ), + $ RWORK( ITGKZ ), M*2, RWORK( ITEMPR ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -719,7 +741,7 @@ END DO K = K + M END DO - CALL ZLASET( 'A', M, N-M, CZERO, CZERO, + CALL ZLASET( 'A', NS, N-M, CZERO, CZERO, $ VT( 1,M+1 ), LDVT ) * * Call ZUNMBR to compute (VB**T)*(PB**T) @@ -755,14 +777,14 @@ CALL ZGEBRD( M, N, A, LDA, RWORK( ID ), RWORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), WORK( ITEMP ), $ LWORK-ITEMP+1, INFO ) - ITEMP = ITGKZ + M*(M*2+1) + ITEMPR = 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( ITEMP ), + $ RWORK( ITGKZ ), M*2, RWORK( ITEMPR ), $ IWORK, INFO) * * If needed, compute left singular vectors. @@ -796,7 +818,7 @@ END DO K = K + M END DO - CALL ZLASET( 'A', M, N-M, CZERO, CZERO, + CALL ZLASET( 'A', NS, N-M, CZERO, CZERO, $ VT( 1,M+1 ), LDVT ) * * Call ZUNMBR to compute VB**T * PB**T diff --git a/lapack-netlib/SRC/zgetc2.f b/lapack-netlib/SRC/zgetc2.f index 3179612f5..bf59415b5 100644 --- a/lapack-netlib/SRC/zgetc2.f +++ b/lapack-netlib/SRC/zgetc2.f @@ -145,15 +145,33 @@ INTRINSIC ABS, DCMPLX, MAX * .. * .. Executable Statements .. +* + INFO = 0 +* +* Quick return if possible +* + IF( N.EQ.0 ) + $ RETURN * * Set constants to control overflow * - INFO = 0 EPS = DLAMCH( 'P' ) SMLNUM = DLAMCH( 'S' ) / EPS BIGNUM = ONE / SMLNUM CALL DLABAD( SMLNUM, BIGNUM ) * +* Handle the case N=1 by itself +* + IF( N.EQ.1 ) THEN + IPIV( 1 ) = 1 + JPIV( 1 ) = 1 + IF( ABS( A( 1, 1 ) ).LT.SMLNUM ) THEN + INFO = 1 + A( 1, 1 ) = DCMPLX( SMLNUM, ZERO ) + END IF + RETURN + END IF +* * Factorize A using complete pivoting. * Set pivots less than SMIN to SMIN * diff --git a/lapack-netlib/SRC/zggev3.f b/lapack-netlib/SRC/zggev3.f index 1c4e832af..78337fd07 100644 --- a/lapack-netlib/SRC/zggev3.f +++ b/lapack-netlib/SRC/zggev3.f @@ -340,7 +340,7 @@ LWKOPT = MAX( LWKOPT, N+INT( WORK( 1 ) ) ) CALL ZHGEQZ( 'S', JOBVL, JOBVR, N, 1, N, A, LDA, B, LDB, $ ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, -1, - $ WORK, IERR ) + $ RWORK, IERR ) LWKOPT = MAX( LWKOPT, N+INT( WORK( 1 ) ) ) ELSE CALL ZGGHD3( JOBVL, JOBVR, N, 1, N, A, LDA, B, LDB, VL, @@ -348,7 +348,7 @@ LWKOPT = MAX( LWKOPT, N+INT( WORK( 1 ) ) ) CALL ZHGEQZ( 'E', JOBVL, JOBVR, N, 1, N, A, LDA, B, LDB, $ ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, -1, - $ WORK, IERR ) + $ RWORK, IERR ) LWKOPT = MAX( LWKOPT, N+INT( WORK( 1 ) ) ) END IF WORK( 1 ) = DCMPLX( LWKOPT )