--- rpl/lapack/lapack/dgesvd.f 2010/01/26 15:22:45 1.1.1.1 +++ rpl/lapack/lapack/dgesvd.f 2017/06/17 10:53:49 1.17 @@ -1,10 +1,220 @@ - SUBROUTINE DGESVD( JOBU, JOBVT, M, N, A, LDA, S, U, LDU, VT, LDVT, - $ WORK, LWORK, INFO ) +*> \brief DGESVD computes the singular value decomposition (SVD) for GE matrices * -* -- LAPACK driver routine (version 3.2) -- +* =========== DOCUMENTATION =========== +* +* Online html documentation available at +* http://www.netlib.org/lapack/explore-html/ +* +*> \htmlonly +*> Download DGESVD + dependencies +*> +*> [TGZ] +*> +*> [ZIP] +*> +*> [TXT] +*> \endhtmlonly +* +* Definition: +* =========== +* +* SUBROUTINE DGESVD( JOBU, JOBVT, M, N, A, LDA, S, U, LDU, VT, LDVT, +* WORK, LWORK, INFO ) +* +* .. Scalar Arguments .. +* CHARACTER JOBU, JOBVT +* INTEGER INFO, LDA, LDU, LDVT, LWORK, M, N +* .. +* .. Array Arguments .. +* DOUBLE PRECISION A( LDA, * ), S( * ), U( LDU, * ), +* $ VT( LDVT, * ), WORK( * ) +* .. +* +* +*> \par Purpose: +* ============= +*> +*> \verbatim +*> +*> DGESVD computes the singular value decomposition (SVD) of a real +*> 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 orthogonal matrix, and +*> V is an N-by-N orthogonal 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. +*> +*> Note that the routine returns V**T, not V. +*> \endverbatim +* +* Arguments: +* ========== +* +*> \param[in] JOBU +*> \verbatim +*> JOBU is CHARACTER*1 +*> Specifies options for computing all or part of the matrix U: +*> = 'A': all M columns of U are returned in array U: +*> = 'S': the first min(m,n) columns of U (the left singular +*> vectors) are returned in the array U; +*> = 'O': the first min(m,n) columns of U (the left singular +*> vectors) are overwritten on the array A; +*> = 'N': no columns of U (no left singular vectors) are +*> computed. +*> \endverbatim +*> +*> \param[in] JOBVT +*> \verbatim +*> JOBVT is CHARACTER*1 +*> Specifies options for computing all or part of the matrix +*> V**T: +*> = 'A': all N rows of V**T are returned in the array VT; +*> = 'S': the first min(m,n) rows of V**T (the right singular +*> vectors) are returned in the array VT; +*> = 'O': the first min(m,n) rows of V**T (the right singular +*> vectors) are overwritten on the array A; +*> = 'N': no rows of V**T (no right singular vectors) are +*> computed. +*> +*> JOBVT and JOBU cannot both be 'O'. +*> \endverbatim +*> +*> \param[in] M +*> \verbatim +*> M is INTEGER +*> The number of rows of the input matrix A. M >= 0. +*> \endverbatim +*> +*> \param[in] N +*> \verbatim +*> N is INTEGER +*> The number of columns of the input matrix A. N >= 0. +*> \endverbatim +*> +*> \param[in,out] A +*> \verbatim +*> A is DOUBLE PRECISION array, dimension (LDA,N) +*> On entry, the M-by-N matrix A. +*> On exit, +*> if JOBU = 'O', A is overwritten with the first min(m,n) +*> columns of U (the left singular vectors, +*> stored columnwise); +*> if JOBVT = 'O', A is overwritten with the first min(m,n) +*> rows of V**T (the right singular vectors, +*> stored rowwise); +*> if JOBU .ne. 'O' and JOBVT .ne. 'O', the contents of A +*> are destroyed. +*> \endverbatim +*> +*> \param[in] LDA +*> \verbatim +*> LDA is INTEGER +*> The leading dimension of the array A. LDA >= max(1,M). +*> \endverbatim +*> +*> \param[out] S +*> \verbatim +*> S is DOUBLE PRECISION array, dimension (min(M,N)) +*> The singular values of A, sorted so that S(i) >= S(i+1). +*> \endverbatim +*> +*> \param[out] U +*> \verbatim +*> U is DOUBLE PRECISION array, dimension (LDU,UCOL) +*> (LDU,M) if JOBU = 'A' or (LDU,min(M,N)) if JOBU = 'S'. +*> If JOBU = 'A', U contains the M-by-M orthogonal matrix U; +*> if JOBU = 'S', U contains the first min(m,n) columns of U +*> (the left singular vectors, stored columnwise); +*> if JOBU = 'N' or 'O', U is not referenced. +*> \endverbatim +*> +*> \param[in] LDU +*> \verbatim +*> LDU is INTEGER +*> The leading dimension of the array U. LDU >= 1; if +*> JOBU = 'S' or 'A', LDU >= M. +*> \endverbatim +*> +*> \param[out] VT +*> \verbatim +*> VT is DOUBLE PRECISION array, dimension (LDVT,N) +*> If JOBVT = 'A', VT contains the N-by-N orthogonal matrix +*> V**T; +*> if JOBVT = 'S', VT contains the first min(m,n) rows of +*> V**T (the right singular vectors, stored rowwise); +*> if JOBVT = 'N' or 'O', VT is not referenced. +*> \endverbatim +*> +*> \param[in] LDVT +*> \verbatim +*> LDVT is INTEGER +*> The leading dimension of the array VT. LDVT >= 1; if +*> JOBVT = 'A', LDVT >= N; if JOBVT = 'S', LDVT >= min(M,N). +*> \endverbatim +*> +*> \param[out] WORK +*> \verbatim +*> WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK)) +*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK; +*> if INFO > 0, WORK(2:MIN(M,N)) contains the unconverged +*> superdiagonal elements of an upper bidiagonal matrix B +*> whose diagonal is in S (not necessarily sorted). B +*> satisfies A = U * B * VT, so it has the same singular values +*> as A, and singular vectors related by U and VT. +*> \endverbatim +*> +*> \param[in] LWORK +*> \verbatim +*> LWORK is INTEGER +*> The dimension of the array WORK. +*> LWORK >= MAX(1,5*MIN(M,N)) for the paths (see comments inside code): +*> - PATH 1 (M much larger than N, JOBU='N') +*> - PATH 1t (N much larger than M, JOBVT='N') +*> LWORK >= MAX(1,3*MIN(M,N) + MAX(M,N),5*MIN(M,N)) for the other paths +*> For good performance, LWORK should generally be larger. +*> +*> If LWORK = -1, then a workspace query is assumed; the routine +*> only calculates the optimal size of the WORK array, returns +*> this value as the first entry of the WORK array, and no error +*> message related to LWORK is issued by XERBLA. +*> \endverbatim +*> +*> \param[out] INFO +*> \verbatim +*> INFO is INTEGER +*> = 0: successful exit. +*> < 0: if INFO = -i, the i-th argument had an illegal value. +*> > 0: if DBDSQR did not converge, INFO specifies how many +*> superdiagonals of an intermediate bidiagonal form B +*> did not converge to zero. See the description of WORK +*> above for details. +*> \endverbatim +* +* Authors: +* ======== +* +*> \author Univ. of Tennessee +*> \author Univ. of California Berkeley +*> \author Univ. of Colorado Denver +*> \author NAG Ltd. +* +*> \date April 2012 +* +*> \ingroup doubleGEsing +* +* ===================================================================== + SUBROUTINE DGESVD( JOBU, JOBVT, M, N, A, LDA, S, U, LDU, + $ VT, LDVT, WORK, LWORK, INFO ) +* +* -- LAPACK driver routine (version 3.7.0) -- * -- LAPACK is a software package provided by Univ. of Tennessee, -- * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- -* November 2006 +* April 2012 * * .. Scalar Arguments .. CHARACTER JOBU, JOBVT @@ -15,122 +225,6 @@ $ VT( LDVT, * ), WORK( * ) * .. * -* Purpose -* ======= -* -* DGESVD computes the singular value decomposition (SVD) of a real -* 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 orthogonal matrix, and -* V is an N-by-N orthogonal 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. -* -* Note that the routine returns V**T, not V. -* -* Arguments -* ========= -* -* JOBU (input) CHARACTER*1 -* Specifies options for computing all or part of the matrix U: -* = 'A': all M columns of U are returned in array U: -* = 'S': the first min(m,n) columns of U (the left singular -* vectors) are returned in the array U; -* = 'O': the first min(m,n) columns of U (the left singular -* vectors) are overwritten on the array A; -* = 'N': no columns of U (no left singular vectors) are -* computed. -* -* JOBVT (input) CHARACTER*1 -* Specifies options for computing all or part of the matrix -* V**T: -* = 'A': all N rows of V**T are returned in the array VT; -* = 'S': the first min(m,n) rows of V**T (the right singular -* vectors) are returned in the array VT; -* = 'O': the first min(m,n) rows of V**T (the right singular -* vectors) are overwritten on the array A; -* = 'N': no rows of V**T (no right singular vectors) are -* computed. -* -* JOBVT and JOBU cannot both be 'O'. -* -* M (input) INTEGER -* The number of rows of the input matrix A. M >= 0. -* -* N (input) INTEGER -* The number of columns of the input matrix A. N >= 0. -* -* A (input/output) DOUBLE PRECISION array, dimension (LDA,N) -* On entry, the M-by-N matrix A. -* On exit, -* if JOBU = 'O', A is overwritten with the first min(m,n) -* columns of U (the left singular vectors, -* stored columnwise); -* if JOBVT = 'O', A is overwritten with the first min(m,n) -* rows of V**T (the right singular vectors, -* stored rowwise); -* if JOBU .ne. 'O' and JOBVT .ne. 'O', the contents of A -* are destroyed. -* -* LDA (input) INTEGER -* The leading dimension of the array A. LDA >= max(1,M). -* -* S (output) DOUBLE PRECISION array, dimension (min(M,N)) -* The singular values of A, sorted so that S(i) >= S(i+1). -* -* U (output) DOUBLE PRECISION array, dimension (LDU,UCOL) -* (LDU,M) if JOBU = 'A' or (LDU,min(M,N)) if JOBU = 'S'. -* If JOBU = 'A', U contains the M-by-M orthogonal matrix U; -* if JOBU = 'S', U contains the first min(m,n) columns of U -* (the left singular vectors, stored columnwise); -* if JOBU = 'N' or 'O', U is not referenced. -* -* LDU (input) INTEGER -* The leading dimension of the array U. LDU >= 1; if -* JOBU = 'S' or 'A', LDU >= M. -* -* VT (output) DOUBLE PRECISION array, dimension (LDVT,N) -* If JOBVT = 'A', VT contains the N-by-N orthogonal matrix -* V**T; -* if JOBVT = 'S', VT contains the first min(m,n) rows of -* V**T (the right singular vectors, stored rowwise); -* if JOBVT = 'N' or 'O', VT is not referenced. -* -* LDVT (input) INTEGER -* The leading dimension of the array VT. LDVT >= 1; if -* JOBVT = 'A', LDVT >= N; if JOBVT = 'S', LDVT >= min(M,N). -* -* WORK (workspace/output) DOUBLE PRECISION array, dimension (MAX(1,LWORK)) -* On exit, if INFO = 0, WORK(1) returns the optimal LWORK; -* if INFO > 0, WORK(2:MIN(M,N)) contains the unconverged -* superdiagonal elements of an upper bidiagonal matrix B -* whose diagonal is in S (not necessarily sorted). B -* satisfies A = U * B * VT, so it has the same singular values -* as A, and singular vectors related by U and VT. -* -* LWORK (input) INTEGER -* The dimension of the array WORK. -* LWORK >= MAX(1,3*MIN(M,N)+MAX(M,N),5*MIN(M,N)). -* For good performance, LWORK should generally be larger. -* -* If LWORK = -1, then a workspace query is assumed; the routine -* only calculates the optimal size of the WORK array, returns -* this value as the first entry of the WORK array, and no error -* message related to LWORK is issued by XERBLA. -* -* INFO (output) INTEGER -* = 0: successful exit. -* < 0: if INFO = -i, the i-th argument had an illegal value. -* > 0: if DBDSQR did not converge, INFO specifies how many -* superdiagonals of an intermediate bidiagonal form B -* did not converge to zero. See the description of WORK -* above for details. -* * ===================================================================== * * .. Parameters .. @@ -144,6 +238,9 @@ $ ITAU, ITAUP, ITAUQ, IU, IWORK, LDWRKR, LDWRKU, $ MAXWRK, MINMN, MINWRK, MNTHR, NCU, NCVT, NRU, $ NRVT, WRKBL + INTEGER LWORK_DGEQRF, LWORK_DORGQR_N, LWORK_DORGQR_M, + $ LWORK_DGEBRD, LWORK_DORGBR_P, LWORK_DORGBR_Q, + $ LWORK_DGELQF, LWORK_DORGLQ_N, LWORK_DORGLQ_M DOUBLE PRECISION ANRM, BIGNUM, EPS, SMLNUM * .. * .. Local Arrays .. @@ -215,163 +312,160 @@ * MNTHR = ILAENV( 6, 'DGESVD', JOBU // JOBVT, M, N, 0, 0 ) BDSPAC = 5*N +* Compute space needed for DGEQRF + CALL DGEQRF( M, N, A, LDA, DUM(1), DUM(1), -1, IERR ) + LWORK_DGEQRF = INT( DUM(1) ) +* Compute space needed for DORGQR + CALL DORGQR( M, N, N, A, LDA, DUM(1), DUM(1), -1, IERR ) + LWORK_DORGQR_N = INT( DUM(1) ) + CALL DORGQR( M, M, N, A, LDA, DUM(1), DUM(1), -1, IERR ) + LWORK_DORGQR_M = INT( DUM(1) ) +* Compute space needed for DGEBRD + CALL DGEBRD( N, N, A, LDA, S, DUM(1), DUM(1), + $ DUM(1), DUM(1), -1, IERR ) + LWORK_DGEBRD = INT( DUM(1) ) +* Compute space needed for DORGBR P + CALL DORGBR( 'P', N, N, N, A, LDA, DUM(1), + $ DUM(1), -1, IERR ) + LWORK_DORGBR_P = INT( DUM(1) ) +* Compute space needed for DORGBR Q + CALL DORGBR( 'Q', N, N, N, A, LDA, DUM(1), + $ DUM(1), -1, IERR ) + LWORK_DORGBR_Q = INT( DUM(1) ) +* IF( M.GE.MNTHR ) THEN IF( WNTUN ) THEN * * Path 1 (M much larger than N, JOBU='N') * - MAXWRK = N + N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, - $ -1 ) - MAXWRK = MAX( MAXWRK, 3*N+2*N* - $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) + MAXWRK = N + LWORK_DGEQRF + MAXWRK = MAX( MAXWRK, 3*N + LWORK_DGEBRD ) IF( WNTVO .OR. WNTVAS ) - $ MAXWRK = MAX( MAXWRK, 3*N+( N-1 )* - $ ILAENV( 1, 'DORGBR', 'P', N, N, N, -1 ) ) + $ MAXWRK = MAX( MAXWRK, 3*N + LWORK_DORGBR_P ) MAXWRK = MAX( MAXWRK, BDSPAC ) MINWRK = MAX( 4*N, BDSPAC ) ELSE IF( WNTUO .AND. WNTVN ) THEN * * Path 2 (M much larger than N, JOBU='O', JOBVT='N') * - WRKBL = N + N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, N+N*ILAENV( 1, 'DORGQR', ' ', M, - $ N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+2*N* - $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+N* - $ ILAENV( 1, 'DORGBR', 'Q', N, N, N, -1 ) ) + WRKBL = N + LWORK_DGEQRF + WRKBL = MAX( WRKBL, N + LWORK_DORGQR_N ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_Q ) WRKBL = MAX( WRKBL, BDSPAC ) - MAXWRK = MAX( N*N+WRKBL, N*N+M*N+N ) - MINWRK = MAX( 3*N+M, BDSPAC ) + MAXWRK = MAX( N*N + WRKBL, N*N + M*N + N ) + MINWRK = MAX( 3*N + M, BDSPAC ) ELSE IF( WNTUO .AND. WNTVAS ) THEN * * Path 3 (M much larger than N, JOBU='O', JOBVT='S' or * 'A') * - WRKBL = N + N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, N+N*ILAENV( 1, 'DORGQR', ' ', M, - $ N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+2*N* - $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+N* - $ ILAENV( 1, 'DORGBR', 'Q', N, N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+( N-1 )* - $ ILAENV( 1, 'DORGBR', 'P', N, N, N, -1 ) ) + WRKBL = N + LWORK_DGEQRF + WRKBL = MAX( WRKBL, N + LWORK_DORGQR_N ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_Q ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_P ) WRKBL = MAX( WRKBL, BDSPAC ) - MAXWRK = MAX( N*N+WRKBL, N*N+M*N+N ) - MINWRK = MAX( 3*N+M, BDSPAC ) + MAXWRK = MAX( N*N + WRKBL, N*N + M*N + N ) + MINWRK = MAX( 3*N + M, BDSPAC ) ELSE IF( WNTUS .AND. WNTVN ) THEN * * Path 4 (M much larger than N, JOBU='S', JOBVT='N') * - WRKBL = N + N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, N+N*ILAENV( 1, 'DORGQR', ' ', M, - $ N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+2*N* - $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+N* - $ ILAENV( 1, 'DORGBR', 'Q', N, N, N, -1 ) ) + WRKBL = N + LWORK_DGEQRF + WRKBL = MAX( WRKBL, N + LWORK_DORGQR_N ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_Q ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = N*N + WRKBL - MINWRK = MAX( 3*N+M, BDSPAC ) + MINWRK = MAX( 3*N + M, BDSPAC ) ELSE IF( WNTUS .AND. WNTVO ) THEN * * Path 5 (M much larger than N, JOBU='S', JOBVT='O') * - WRKBL = N + N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, N+N*ILAENV( 1, 'DORGQR', ' ', M, - $ N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+2*N* - $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+N* - $ ILAENV( 1, 'DORGBR', 'Q', N, N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+( N-1 )* - $ ILAENV( 1, 'DORGBR', 'P', N, N, N, -1 ) ) + WRKBL = N + LWORK_DGEQRF + WRKBL = MAX( WRKBL, N + LWORK_DORGQR_N ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_Q ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_P ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = 2*N*N + WRKBL - MINWRK = MAX( 3*N+M, BDSPAC ) + MINWRK = MAX( 3*N + M, BDSPAC ) ELSE IF( WNTUS .AND. WNTVAS ) THEN * * Path 6 (M much larger than N, JOBU='S', JOBVT='S' or * 'A') * - WRKBL = N + N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, N+N*ILAENV( 1, 'DORGQR', ' ', M, - $ N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+2*N* - $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+N* - $ ILAENV( 1, 'DORGBR', 'Q', N, N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+( N-1 )* - $ ILAENV( 1, 'DORGBR', 'P', N, N, N, -1 ) ) + WRKBL = N + LWORK_DGEQRF + WRKBL = MAX( WRKBL, N + LWORK_DORGQR_N ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_Q ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_P ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = N*N + WRKBL - MINWRK = MAX( 3*N+M, BDSPAC ) + MINWRK = MAX( 3*N + M, BDSPAC ) ELSE IF( WNTUA .AND. WNTVN ) THEN * * Path 7 (M much larger than N, JOBU='A', JOBVT='N') * - WRKBL = N + N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, N+M*ILAENV( 1, 'DORGQR', ' ', M, - $ M, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+2*N* - $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+N* - $ ILAENV( 1, 'DORGBR', 'Q', N, N, N, -1 ) ) + WRKBL = N + LWORK_DGEQRF + WRKBL = MAX( WRKBL, N + LWORK_DORGQR_M ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_Q ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = N*N + WRKBL - MINWRK = MAX( 3*N+M, BDSPAC ) + MINWRK = MAX( 3*N + M, BDSPAC ) ELSE IF( WNTUA .AND. WNTVO ) THEN * * Path 8 (M much larger than N, JOBU='A', JOBVT='O') * - WRKBL = N + N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, N+M*ILAENV( 1, 'DORGQR', ' ', M, - $ M, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+2*N* - $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+N* - $ ILAENV( 1, 'DORGBR', 'Q', N, N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+( N-1 )* - $ ILAENV( 1, 'DORGBR', 'P', N, N, N, -1 ) ) + WRKBL = N + LWORK_DGEQRF + WRKBL = MAX( WRKBL, N + LWORK_DORGQR_M ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_Q ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_P ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = 2*N*N + WRKBL - MINWRK = MAX( 3*N+M, BDSPAC ) + MINWRK = MAX( 3*N + M, BDSPAC ) ELSE IF( WNTUA .AND. WNTVAS ) THEN * * Path 9 (M much larger than N, JOBU='A', JOBVT='S' or * 'A') * - WRKBL = N + N*ILAENV( 1, 'DGEQRF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, N+M*ILAENV( 1, 'DORGQR', ' ', M, - $ M, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+2*N* - $ ILAENV( 1, 'DGEBRD', ' ', N, N, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+N* - $ ILAENV( 1, 'DORGBR', 'Q', N, N, N, -1 ) ) - WRKBL = MAX( WRKBL, 3*N+( N-1 )* - $ ILAENV( 1, 'DORGBR', 'P', N, N, N, -1 ) ) + WRKBL = N + LWORK_DGEQRF + WRKBL = MAX( WRKBL, N + LWORK_DORGQR_M ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_Q ) + WRKBL = MAX( WRKBL, 3*N + LWORK_DORGBR_P ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = N*N + WRKBL - MINWRK = MAX( 3*N+M, BDSPAC ) + MINWRK = MAX( 3*N + M, BDSPAC ) END IF ELSE * * Path 10 (M at least N, but not much larger) * - MAXWRK = 3*N + ( M+N )*ILAENV( 1, 'DGEBRD', ' ', M, N, - $ -1, -1 ) - IF( WNTUS .OR. WNTUO ) - $ MAXWRK = MAX( MAXWRK, 3*N+N* - $ ILAENV( 1, 'DORGBR', 'Q', M, N, N, -1 ) ) - IF( WNTUA ) - $ MAXWRK = MAX( MAXWRK, 3*N+M* - $ ILAENV( 1, 'DORGBR', 'Q', M, M, N, -1 ) ) - IF( .NOT.WNTVN ) - $ MAXWRK = MAX( MAXWRK, 3*N+( N-1 )* - $ ILAENV( 1, 'DORGBR', 'P', N, N, N, -1 ) ) + CALL DGEBRD( M, N, A, LDA, S, DUM(1), DUM(1), + $ DUM(1), DUM(1), -1, IERR ) + LWORK_DGEBRD = INT( DUM(1) ) + MAXWRK = 3*N + LWORK_DGEBRD + IF( WNTUS .OR. WNTUO ) THEN + CALL DORGBR( 'Q', M, N, N, A, LDA, DUM(1), + $ DUM(1), -1, IERR ) + LWORK_DORGBR_Q = INT( DUM(1) ) + MAXWRK = MAX( MAXWRK, 3*N + LWORK_DORGBR_Q ) + END IF + IF( WNTUA ) THEN + CALL DORGBR( 'Q', M, M, N, A, LDA, DUM(1), + $ DUM(1), -1, IERR ) + LWORK_DORGBR_Q = INT( DUM(1) ) + MAXWRK = MAX( MAXWRK, 3*N + LWORK_DORGBR_Q ) + END IF + IF( .NOT.WNTVN ) THEN + MAXWRK = MAX( MAXWRK, 3*N + LWORK_DORGBR_P ) + END IF MAXWRK = MAX( MAXWRK, BDSPAC ) - MINWRK = MAX( 3*N+M, BDSPAC ) + MINWRK = MAX( 3*N + M, BDSPAC ) END IF ELSE IF( MINMN.GT.0 ) THEN * @@ -379,163 +473,160 @@ * MNTHR = ILAENV( 6, 'DGESVD', JOBU // JOBVT, M, N, 0, 0 ) BDSPAC = 5*M +* Compute space needed for DGELQF + CALL DGELQF( M, N, A, LDA, DUM(1), DUM(1), -1, IERR ) + LWORK_DGELQF = INT( DUM(1) ) +* Compute space needed for DORGLQ + CALL DORGLQ( N, N, M, DUM(1), N, DUM(1), DUM(1), -1, IERR ) + LWORK_DORGLQ_N = INT( DUM(1) ) + CALL DORGLQ( M, N, M, A, LDA, DUM(1), DUM(1), -1, IERR ) + LWORK_DORGLQ_M = INT( DUM(1) ) +* Compute space needed for DGEBRD + CALL DGEBRD( M, M, A, LDA, S, DUM(1), DUM(1), + $ DUM(1), DUM(1), -1, IERR ) + LWORK_DGEBRD = INT( DUM(1) ) +* Compute space needed for DORGBR P + CALL DORGBR( 'P', M, M, M, A, N, DUM(1), + $ DUM(1), -1, IERR ) + LWORK_DORGBR_P = INT( DUM(1) ) +* Compute space needed for DORGBR Q + CALL DORGBR( 'Q', M, M, M, A, N, DUM(1), + $ DUM(1), -1, IERR ) + LWORK_DORGBR_Q = INT( DUM(1) ) IF( N.GE.MNTHR ) THEN IF( WNTVN ) THEN * * Path 1t(N much larger than M, JOBVT='N') * - MAXWRK = M + M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, - $ -1 ) - MAXWRK = MAX( MAXWRK, 3*M+2*M* - $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) + MAXWRK = M + LWORK_DGELQF + MAXWRK = MAX( MAXWRK, 3*M + LWORK_DGEBRD ) IF( WNTUO .OR. WNTUAS ) - $ MAXWRK = MAX( MAXWRK, 3*M+M* - $ ILAENV( 1, 'DORGBR', 'Q', M, M, M, -1 ) ) + $ MAXWRK = MAX( MAXWRK, 3*M + LWORK_DORGBR_Q ) MAXWRK = MAX( MAXWRK, BDSPAC ) MINWRK = MAX( 4*M, BDSPAC ) ELSE IF( WNTVO .AND. WNTUN ) THEN * * Path 2t(N much larger than M, JOBU='N', JOBVT='O') * - WRKBL = M + M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, M+M*ILAENV( 1, 'DORGLQ', ' ', M, - $ N, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+2*M* - $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+( M-1 )* - $ ILAENV( 1, 'DORGBR', 'P', M, M, M, -1 ) ) + WRKBL = M + LWORK_DGELQF + WRKBL = MAX( WRKBL, M + LWORK_DORGLQ_M ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_P ) WRKBL = MAX( WRKBL, BDSPAC ) - MAXWRK = MAX( M*M+WRKBL, M*M+M*N+M ) - MINWRK = MAX( 3*M+N, BDSPAC ) + MAXWRK = MAX( M*M + WRKBL, M*M + M*N + M ) + MINWRK = MAX( 3*M + N, BDSPAC ) ELSE IF( WNTVO .AND. WNTUAS ) THEN * * Path 3t(N much larger than M, JOBU='S' or 'A', * JOBVT='O') * - WRKBL = M + M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, M+M*ILAENV( 1, 'DORGLQ', ' ', M, - $ N, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+2*M* - $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+( M-1 )* - $ ILAENV( 1, 'DORGBR', 'P', M, M, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+M* - $ ILAENV( 1, 'DORGBR', 'Q', M, M, M, -1 ) ) + WRKBL = M + LWORK_DGELQF + WRKBL = MAX( WRKBL, M + LWORK_DORGLQ_M ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_P ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_Q ) WRKBL = MAX( WRKBL, BDSPAC ) - MAXWRK = MAX( M*M+WRKBL, M*M+M*N+M ) - MINWRK = MAX( 3*M+N, BDSPAC ) + MAXWRK = MAX( M*M + WRKBL, M*M + M*N + M ) + MINWRK = MAX( 3*M + N, BDSPAC ) ELSE IF( WNTVS .AND. WNTUN ) THEN * * Path 4t(N much larger than M, JOBU='N', JOBVT='S') * - WRKBL = M + M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, M+M*ILAENV( 1, 'DORGLQ', ' ', M, - $ N, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+2*M* - $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+( M-1 )* - $ ILAENV( 1, 'DORGBR', 'P', M, M, M, -1 ) ) + WRKBL = M + LWORK_DGELQF + WRKBL = MAX( WRKBL, M + LWORK_DORGLQ_M ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_P ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = M*M + WRKBL - MINWRK = MAX( 3*M+N, BDSPAC ) + MINWRK = MAX( 3*M + N, BDSPAC ) ELSE IF( WNTVS .AND. WNTUO ) THEN * * Path 5t(N much larger than M, JOBU='O', JOBVT='S') * - WRKBL = M + M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, M+M*ILAENV( 1, 'DORGLQ', ' ', M, - $ N, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+2*M* - $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+( M-1 )* - $ ILAENV( 1, 'DORGBR', 'P', M, M, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+M* - $ ILAENV( 1, 'DORGBR', 'Q', M, M, M, -1 ) ) + WRKBL = M + LWORK_DGELQF + WRKBL = MAX( WRKBL, M + LWORK_DORGLQ_M ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_P ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_Q ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = 2*M*M + WRKBL - MINWRK = MAX( 3*M+N, BDSPAC ) + MINWRK = MAX( 3*M + N, BDSPAC ) ELSE IF( WNTVS .AND. WNTUAS ) THEN * * Path 6t(N much larger than M, JOBU='S' or 'A', * JOBVT='S') * - WRKBL = M + M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, M+M*ILAENV( 1, 'DORGLQ', ' ', M, - $ N, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+2*M* - $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+( M-1 )* - $ ILAENV( 1, 'DORGBR', 'P', M, M, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+M* - $ ILAENV( 1, 'DORGBR', 'Q', M, M, M, -1 ) ) + WRKBL = M + LWORK_DGELQF + WRKBL = MAX( WRKBL, M + LWORK_DORGLQ_M ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_P ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_Q ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = M*M + WRKBL - MINWRK = MAX( 3*M+N, BDSPAC ) + MINWRK = MAX( 3*M + N, BDSPAC ) ELSE IF( WNTVA .AND. WNTUN ) THEN * * Path 7t(N much larger than M, JOBU='N', JOBVT='A') * - WRKBL = M + M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, M+N*ILAENV( 1, 'DORGLQ', ' ', N, - $ N, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+2*M* - $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+( M-1 )* - $ ILAENV( 1, 'DORGBR', 'P', M, M, M, -1 ) ) + WRKBL = M + LWORK_DGELQF + WRKBL = MAX( WRKBL, M + LWORK_DORGLQ_N ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_P ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = M*M + WRKBL - MINWRK = MAX( 3*M+N, BDSPAC ) + MINWRK = MAX( 3*M + N, BDSPAC ) ELSE IF( WNTVA .AND. WNTUO ) THEN * * Path 8t(N much larger than M, JOBU='O', JOBVT='A') * - WRKBL = M + M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, M+N*ILAENV( 1, 'DORGLQ', ' ', N, - $ N, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+2*M* - $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+( M-1 )* - $ ILAENV( 1, 'DORGBR', 'P', M, M, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+M* - $ ILAENV( 1, 'DORGBR', 'Q', M, M, M, -1 ) ) + WRKBL = M + LWORK_DGELQF + WRKBL = MAX( WRKBL, M + LWORK_DORGLQ_N ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_P ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_Q ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = 2*M*M + WRKBL - MINWRK = MAX( 3*M+N, BDSPAC ) + MINWRK = MAX( 3*M + N, BDSPAC ) ELSE IF( WNTVA .AND. WNTUAS ) THEN * * Path 9t(N much larger than M, JOBU='S' or 'A', * JOBVT='A') * - WRKBL = M + M*ILAENV( 1, 'DGELQF', ' ', M, N, -1, -1 ) - WRKBL = MAX( WRKBL, M+N*ILAENV( 1, 'DORGLQ', ' ', N, - $ N, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+2*M* - $ ILAENV( 1, 'DGEBRD', ' ', M, M, -1, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+( M-1 )* - $ ILAENV( 1, 'DORGBR', 'P', M, M, M, -1 ) ) - WRKBL = MAX( WRKBL, 3*M+M* - $ ILAENV( 1, 'DORGBR', 'Q', M, M, M, -1 ) ) + WRKBL = M + LWORK_DGELQF + WRKBL = MAX( WRKBL, M + LWORK_DORGLQ_N ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DGEBRD ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_P ) + WRKBL = MAX( WRKBL, 3*M + LWORK_DORGBR_Q ) WRKBL = MAX( WRKBL, BDSPAC ) MAXWRK = M*M + WRKBL - MINWRK = MAX( 3*M+N, BDSPAC ) + MINWRK = MAX( 3*M + N, BDSPAC ) END IF ELSE * * Path 10t(N greater than M, but not much larger) * - MAXWRK = 3*M + ( M+N )*ILAENV( 1, 'DGEBRD', ' ', M, N, - $ -1, -1 ) - IF( WNTVS .OR. WNTVO ) - $ MAXWRK = MAX( MAXWRK, 3*M+M* - $ ILAENV( 1, 'DORGBR', 'P', M, N, M, -1 ) ) - IF( WNTVA ) - $ MAXWRK = MAX( MAXWRK, 3*M+N* - $ ILAENV( 1, 'DORGBR', 'P', N, N, M, -1 ) ) - IF( .NOT.WNTUN ) - $ MAXWRK = MAX( MAXWRK, 3*M+( M-1 )* - $ ILAENV( 1, 'DORGBR', 'Q', M, M, M, -1 ) ) + CALL DGEBRD( M, N, A, LDA, S, DUM(1), DUM(1), + $ DUM(1), DUM(1), -1, IERR ) + LWORK_DGEBRD = INT( DUM(1) ) + MAXWRK = 3*M + LWORK_DGEBRD + IF( WNTVS .OR. WNTVO ) THEN +* Compute space needed for DORGBR P + CALL DORGBR( 'P', M, N, M, A, N, DUM(1), + $ DUM(1), -1, IERR ) + LWORK_DORGBR_P = INT( DUM(1) ) + MAXWRK = MAX( MAXWRK, 3*M + LWORK_DORGBR_P ) + END IF + IF( WNTVA ) THEN + CALL DORGBR( 'P', N, N, M, A, N, DUM(1), + $ DUM(1), -1, IERR ) + LWORK_DORGBR_P = INT( DUM(1) ) + MAXWRK = MAX( MAXWRK, 3*M + LWORK_DORGBR_P ) + END IF + IF( .NOT.WNTUN ) THEN + MAXWRK = MAX( MAXWRK, 3*M + LWORK_DORGBR_Q ) + END IF MAXWRK = MAX( MAXWRK, BDSPAC ) - MINWRK = MAX( 3*M+N, BDSPAC ) + MINWRK = MAX( 3*M + N, BDSPAC ) END IF END IF MAXWRK = MAX( MAXWRK, MINWRK ) @@ -594,21 +685,24 @@ IWORK = ITAU + N * * Compute A=Q*R -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Zero out below R * - CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, A( 2, 1 ), LDA ) + IF( N .GT. 1 ) THEN + CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, A( 2, 1 ), + $ LDA ) + END IF IE = 1 ITAUQ = IE + N ITAUP = ITAUQ + N IWORK = ITAUP + N * * Bidiagonalize R in A -* (Workspace: need 4*N, prefer 3*N+2*N*NB) +* (Workspace: need 4*N, prefer 3*N + 2*N*NB) * CALL DGEBRD( N, N, A, LDA, S, WORK( IE ), WORK( ITAUQ ), $ WORK( ITAUP ), WORK( IWORK ), LWORK-IWORK+1, @@ -617,7 +711,7 @@ IF( WNTVO .OR. WNTVAS ) THEN * * If right singular vectors desired, generate P'. -* (Workspace: need 4*N-1, prefer 3*N+(N-1)*NB) +* (Workspace: need 4*N-1, prefer 3*N + (N-1)*NB) * CALL DORGBR( 'P', N, N, N, A, LDA, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -648,13 +742,13 @@ * Sufficient workspace for a fast algorithm * IR = 1 - IF( LWORK.GE.MAX( WRKBL, LDA*N+N )+LDA*N ) THEN + IF( LWORK.GE.MAX( WRKBL, LDA*N + N ) + LDA*N ) THEN * * WORK(IU) is LDA by N, WORK(IR) is LDA by N * LDWRKU = LDA LDWRKR = LDA - ELSE IF( LWORK.GE.MAX( WRKBL, LDA*N+N )+N*N ) THEN + ELSE IF( LWORK.GE.MAX( WRKBL, LDA*N + N ) + N*N ) THEN * * WORK(IU) is LDA by N, WORK(IR) is N by N * @@ -671,7 +765,7 @@ IWORK = ITAU + N * * Compute A=Q*R -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -683,7 +777,7 @@ $ LDWRKR ) * * Generate Q in A -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DORGQR( M, N, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -693,14 +787,14 @@ IWORK = ITAUP + N * * Bidiagonalize R in WORK(IR) -* (Workspace: need N*N+4*N, prefer N*N+3*N+2*N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + 2*N*NB) * CALL DGEBRD( N, N, WORK( IR ), LDWRKR, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Generate left vectors bidiagonalizing R -* (Workspace: need N*N+4*N, prefer N*N+3*N+N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + N*NB) * CALL DORGBR( 'Q', N, N, N, WORK( IR ), LDWRKR, $ WORK( ITAUQ ), WORK( IWORK ), @@ -709,7 +803,7 @@ * * Perform bidiagonal QR iteration, computing left * singular vectors of R in WORK(IR) -* (Workspace: need N*N+BDSPAC) +* (Workspace: need N*N + BDSPAC) * CALL DBDSQR( 'U', N, 0, N, 0, S, WORK( IE ), DUM, 1, $ WORK( IR ), LDWRKR, DUM, 1, @@ -718,7 +812,7 @@ * * Multiply Q in A by left singular vectors of R in * WORK(IR), storing result in WORK(IU) and copying to A -* (Workspace: need N*N+2*N, prefer N*N+M*N+N) +* (Workspace: need N*N + 2*N, prefer N*N + M*N + N) * DO 10 I = 1, M, LDWRKU CHUNK = MIN( M-I+1, LDWRKU ) @@ -739,14 +833,14 @@ IWORK = ITAUP + N * * Bidiagonalize A -* (Workspace: need 3*N+M, prefer 3*N+(M+N)*NB) +* (Workspace: need 3*N + M, prefer 3*N + (M + N)*NB) * CALL DGEBRD( M, N, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Generate left vectors bidiagonalizing A -* (Workspace: need 4*N, prefer 3*N+N*NB) +* (Workspace: need 4*N, prefer 3*N + N*NB) * CALL DORGBR( 'Q', M, N, N, A, LDA, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -772,13 +866,13 @@ * Sufficient workspace for a fast algorithm * IR = 1 - IF( LWORK.GE.MAX( WRKBL, LDA*N+N )+LDA*N ) THEN + IF( LWORK.GE.MAX( WRKBL, LDA*N + N ) + LDA*N ) THEN * * WORK(IU) is LDA by N and WORK(IR) is LDA by N * LDWRKU = LDA LDWRKR = LDA - ELSE IF( LWORK.GE.MAX( WRKBL, LDA*N+N )+N*N ) THEN + ELSE IF( LWORK.GE.MAX( WRKBL, LDA*N + N ) + N*N ) THEN * * WORK(IU) is LDA by N and WORK(IR) is N by N * @@ -795,7 +889,7 @@ IWORK = ITAU + N * * Compute A=Q*R -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -808,7 +902,7 @@ $ VT( 2, 1 ), LDVT ) * * Generate Q in A -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DORGQR( M, N, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -818,7 +912,7 @@ IWORK = ITAUP + N * * Bidiagonalize R in VT, copying result to WORK(IR) -* (Workspace: need N*N+4*N, prefer N*N+3*N+2*N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + 2*N*NB) * CALL DGEBRD( N, N, VT, LDVT, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -826,14 +920,14 @@ CALL DLACPY( 'L', N, N, VT, LDVT, WORK( IR ), LDWRKR ) * * Generate left vectors bidiagonalizing R in WORK(IR) -* (Workspace: need N*N+4*N, prefer N*N+3*N+N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + N*NB) * CALL DORGBR( 'Q', N, N, N, WORK( IR ), LDWRKR, $ WORK( ITAUQ ), WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right vectors bidiagonalizing R in VT -* (Workspace: need N*N+4*N-1, prefer N*N+3*N+(N-1)*NB) +* (Workspace: need N*N + 4*N-1, prefer N*N + 3*N + (N-1)*NB) * CALL DORGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -842,7 +936,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of R in WORK(IR) and computing right * singular vectors of R in VT -* (Workspace: need N*N+BDSPAC) +* (Workspace: need N*N + BDSPAC) * CALL DBDSQR( 'U', N, N, N, 0, S, WORK( IE ), VT, LDVT, $ WORK( IR ), LDWRKR, DUM, 1, @@ -851,7 +945,7 @@ * * Multiply Q in A by left singular vectors of R in * WORK(IR), storing result in WORK(IU) and copying to A -* (Workspace: need N*N+2*N, prefer N*N+M*N+N) +* (Workspace: need N*N + 2*N, prefer N*N + M*N + N) * DO 20 I = 1, M, LDWRKU CHUNK = MIN( M-I+1, LDWRKU ) @@ -870,7 +964,7 @@ IWORK = ITAU + N * * Compute A=Q*R -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -883,7 +977,7 @@ $ VT( 2, 1 ), LDVT ) * * Generate Q in A -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DORGQR( M, N, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -893,21 +987,21 @@ IWORK = ITAUP + N * * Bidiagonalize R in VT -* (Workspace: need 4*N, prefer 3*N+2*N*NB) +* (Workspace: need 4*N, prefer 3*N + 2*N*NB) * CALL DGEBRD( N, N, VT, LDVT, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Multiply Q in A by left vectors bidiagonalizing R -* (Workspace: need 3*N+M, prefer 3*N+M*NB) +* (Workspace: need 3*N + M, prefer 3*N + M*NB) * CALL DORMBR( 'Q', 'R', 'N', M, N, N, VT, LDVT, $ WORK( ITAUQ ), A, LDA, WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right vectors bidiagonalizing R in VT -* (Workspace: need 4*N-1, prefer 3*N+(N-1)*NB) +* (Workspace: need 4*N-1, prefer 3*N + (N-1)*NB) * CALL DORGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -951,7 +1045,7 @@ IWORK = ITAU + N * * Compute A=Q*R -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -964,7 +1058,7 @@ $ WORK( IR+1 ), LDWRKR ) * * Generate Q in A -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DORGQR( M, N, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -974,7 +1068,7 @@ IWORK = ITAUP + N * * Bidiagonalize R in WORK(IR) -* (Workspace: need N*N+4*N, prefer N*N+3*N+2*N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + 2*N*NB) * CALL DGEBRD( N, N, WORK( IR ), LDWRKR, S, $ WORK( IE ), WORK( ITAUQ ), @@ -982,7 +1076,7 @@ $ LWORK-IWORK+1, IERR ) * * Generate left vectors bidiagonalizing R in WORK(IR) -* (Workspace: need N*N+4*N, prefer N*N+3*N+N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + N*NB) * CALL DORGBR( 'Q', N, N, N, WORK( IR ), LDWRKR, $ WORK( ITAUQ ), WORK( IWORK ), @@ -991,7 +1085,7 @@ * * Perform bidiagonal QR iteration, computing left * singular vectors of R in WORK(IR) -* (Workspace: need N*N+BDSPAC) +* (Workspace: need N*N + BDSPAC) * CALL DBDSQR( 'U', N, 0, N, 0, S, WORK( IE ), DUM, $ 1, WORK( IR ), LDWRKR, DUM, 1, @@ -1012,14 +1106,14 @@ IWORK = ITAU + N * * Compute A=Q*R, copying result to U -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'L', M, N, A, LDA, U, LDU ) * * Generate Q in U -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DORGQR( M, N, N, U, LDU, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1030,18 +1124,20 @@ * * Zero out below R in A * - CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, A( 2, 1 ), - $ LDA ) + IF( N .GT. 1 ) THEN + CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, + $ A( 2, 1 ), LDA ) + END IF * * Bidiagonalize R in A -* (Workspace: need 4*N, prefer 3*N+2*N*NB) +* (Workspace: need 4*N, prefer 3*N + 2*N*NB) * CALL DGEBRD( N, N, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Multiply Q in U by left vectors bidiagonalizing R -* (Workspace: need 3*N+M, prefer 3*N+M*NB) +* (Workspace: need 3*N + M, prefer 3*N + M*NB) * CALL DORMBR( 'Q', 'R', 'N', M, N, N, A, LDA, $ WORK( ITAUQ ), U, LDU, WORK( IWORK ), @@ -1076,7 +1172,7 @@ LDWRKU = LDA IR = IU + LDWRKU*N LDWRKR = LDA - ELSE IF( LWORK.GE.WRKBL+( LDA+N )*N ) THEN + ELSE IF( LWORK.GE.WRKBL+( LDA + N )*N ) THEN * * WORK(IU) is LDA by N and WORK(IR) is N by N * @@ -1095,7 +1191,7 @@ IWORK = ITAU + N * * Compute A=Q*R -* (Workspace: need 2*N*N+2*N, prefer 2*N*N+N+N*NB) +* (Workspace: need 2*N*N + 2*N, prefer 2*N*N + N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1108,7 +1204,7 @@ $ WORK( IU+1 ), LDWRKU ) * * Generate Q in A -* (Workspace: need 2*N*N+2*N, prefer 2*N*N+N+N*NB) +* (Workspace: need 2*N*N + 2*N, prefer 2*N*N + N + N*NB) * CALL DORGQR( M, N, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1119,7 +1215,7 @@ * * Bidiagonalize R in WORK(IU), copying result to * WORK(IR) -* (Workspace: need 2*N*N+4*N, +* (Workspace: need 2*N*N + 4*N, * prefer 2*N*N+3*N+2*N*NB) * CALL DGEBRD( N, N, WORK( IU ), LDWRKU, S, @@ -1130,14 +1226,14 @@ $ WORK( IR ), LDWRKR ) * * Generate left bidiagonalizing vectors in WORK(IU) -* (Workspace: need 2*N*N+4*N, prefer 2*N*N+3*N+N*NB) +* (Workspace: need 2*N*N + 4*N, prefer 2*N*N + 3*N + N*NB) * CALL DORGBR( 'Q', N, N, N, WORK( IU ), LDWRKU, $ WORK( ITAUQ ), WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right bidiagonalizing vectors in WORK(IR) -* (Workspace: need 2*N*N+4*N-1, +* (Workspace: need 2*N*N + 4*N-1, * prefer 2*N*N+3*N+(N-1)*NB) * CALL DORGBR( 'P', N, N, N, WORK( IR ), LDWRKR, @@ -1148,7 +1244,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of R in WORK(IU) and computing * right singular vectors of R in WORK(IR) -* (Workspace: need 2*N*N+BDSPAC) +* (Workspace: need 2*N*N + BDSPAC) * CALL DBDSQR( 'U', N, N, N, 0, S, WORK( IE ), $ WORK( IR ), LDWRKR, WORK( IU ), @@ -1175,14 +1271,14 @@ IWORK = ITAU + N * * Compute A=Q*R, copying result to U -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'L', M, N, A, LDA, U, LDU ) * * Generate Q in U -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DORGQR( M, N, N, U, LDU, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1193,25 +1289,27 @@ * * Zero out below R in A * - CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, A( 2, 1 ), - $ LDA ) + IF( N .GT. 1 ) THEN + CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, + $ A( 2, 1 ), LDA ) + END IF * * Bidiagonalize R in A -* (Workspace: need 4*N, prefer 3*N+2*N*NB) +* (Workspace: need 4*N, prefer 3*N + 2*N*NB) * CALL DGEBRD( N, N, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Multiply Q in U by left vectors bidiagonalizing R -* (Workspace: need 3*N+M, prefer 3*N+M*NB) +* (Workspace: need 3*N + M, prefer 3*N + M*NB) * CALL DORMBR( 'Q', 'R', 'N', M, N, N, A, LDA, $ WORK( ITAUQ ), U, LDU, WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right vectors bidiagonalizing R in A -* (Workspace: need 4*N-1, prefer 3*N+(N-1)*NB) +* (Workspace: need 4*N-1, prefer 3*N + (N-1)*NB) * CALL DORGBR( 'P', N, N, N, A, LDA, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1255,7 +1353,7 @@ IWORK = ITAU + N * * Compute A=Q*R -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1268,7 +1366,7 @@ $ WORK( IU+1 ), LDWRKU ) * * Generate Q in A -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DORGQR( M, N, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1278,7 +1376,7 @@ IWORK = ITAUP + N * * Bidiagonalize R in WORK(IU), copying result to VT -* (Workspace: need N*N+4*N, prefer N*N+3*N+2*N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + 2*N*NB) * CALL DGEBRD( N, N, WORK( IU ), LDWRKU, S, $ WORK( IE ), WORK( ITAUQ ), @@ -1288,14 +1386,14 @@ $ LDVT ) * * Generate left bidiagonalizing vectors in WORK(IU) -* (Workspace: need N*N+4*N, prefer N*N+3*N+N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + N*NB) * CALL DORGBR( 'Q', N, N, N, WORK( IU ), LDWRKU, $ WORK( ITAUQ ), WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right bidiagonalizing vectors in VT -* (Workspace: need N*N+4*N-1, +* (Workspace: need N*N + 4*N-1, * prefer N*N+3*N+(N-1)*NB) * CALL DORGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), @@ -1305,7 +1403,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of R in WORK(IU) and computing * right singular vectors of R in VT -* (Workspace: need N*N+BDSPAC) +* (Workspace: need N*N + BDSPAC) * CALL DBDSQR( 'U', N, N, N, 0, S, WORK( IE ), VT, $ LDVT, WORK( IU ), LDWRKU, DUM, 1, @@ -1326,14 +1424,14 @@ IWORK = ITAU + N * * Compute A=Q*R, copying result to U -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'L', M, N, A, LDA, U, LDU ) * * Generate Q in U -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DORGQR( M, N, N, U, LDU, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1350,7 +1448,7 @@ IWORK = ITAUP + N * * Bidiagonalize R in VT -* (Workspace: need 4*N, prefer 3*N+2*N*NB) +* (Workspace: need 4*N, prefer 3*N + 2*N*NB) * CALL DGEBRD( N, N, VT, LDVT, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -1358,14 +1456,14 @@ * * Multiply Q in U by left bidiagonalizing vectors * in VT -* (Workspace: need 3*N+M, prefer 3*N+M*NB) +* (Workspace: need 3*N + M, prefer 3*N + M*NB) * CALL DORMBR( 'Q', 'R', 'N', M, N, N, VT, LDVT, $ WORK( ITAUQ ), U, LDU, WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right bidiagonalizing vectors in VT -* (Workspace: need 4*N-1, prefer 3*N+(N-1)*NB) +* (Workspace: need 4*N-1, prefer 3*N + (N-1)*NB) * CALL DORGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1412,7 +1510,7 @@ IWORK = ITAU + N * * Compute A=Q*R, copying result to U -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1426,7 +1524,7 @@ $ WORK( IR+1 ), LDWRKR ) * * Generate Q in U -* (Workspace: need N*N+N+M, prefer N*N+N+M*NB) +* (Workspace: need N*N + N + M, prefer N*N + N + M*NB) * CALL DORGQR( M, M, N, U, LDU, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1436,7 +1534,7 @@ IWORK = ITAUP + N * * Bidiagonalize R in WORK(IR) -* (Workspace: need N*N+4*N, prefer N*N+3*N+2*N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + 2*N*NB) * CALL DGEBRD( N, N, WORK( IR ), LDWRKR, S, $ WORK( IE ), WORK( ITAUQ ), @@ -1444,7 +1542,7 @@ $ LWORK-IWORK+1, IERR ) * * Generate left bidiagonalizing vectors in WORK(IR) -* (Workspace: need N*N+4*N, prefer N*N+3*N+N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + N*NB) * CALL DORGBR( 'Q', N, N, N, WORK( IR ), LDWRKR, $ WORK( ITAUQ ), WORK( IWORK ), @@ -1453,7 +1551,7 @@ * * Perform bidiagonal QR iteration, computing left * singular vectors of R in WORK(IR) -* (Workspace: need N*N+BDSPAC) +* (Workspace: need N*N + BDSPAC) * CALL DBDSQR( 'U', N, 0, N, 0, S, WORK( IE ), DUM, $ 1, WORK( IR ), LDWRKR, DUM, 1, @@ -1478,14 +1576,14 @@ IWORK = ITAU + N * * Compute A=Q*R, copying result to U -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'L', M, N, A, LDA, U, LDU ) * * Generate Q in U -* (Workspace: need N+M, prefer N+M*NB) +* (Workspace: need N + M, prefer N + M*NB) * CALL DORGQR( M, M, N, U, LDU, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1496,11 +1594,13 @@ * * Zero out below R in A * - CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, A( 2, 1 ), - $ LDA ) + IF( N .GT. 1 ) THEN + CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, + $ A( 2, 1 ), LDA ) + END IF * * Bidiagonalize R in A -* (Workspace: need 4*N, prefer 3*N+2*N*NB) +* (Workspace: need 4*N, prefer 3*N + 2*N*NB) * CALL DGEBRD( N, N, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -1508,7 +1608,7 @@ * * Multiply Q in U by left bidiagonalizing vectors * in A -* (Workspace: need 3*N+M, prefer 3*N+M*NB) +* (Workspace: need 3*N + M, prefer 3*N + M*NB) * CALL DORMBR( 'Q', 'R', 'N', M, N, N, A, LDA, $ WORK( ITAUQ ), U, LDU, WORK( IWORK ), @@ -1543,7 +1643,7 @@ LDWRKU = LDA IR = IU + LDWRKU*N LDWRKR = LDA - ELSE IF( LWORK.GE.WRKBL+( LDA+N )*N ) THEN + ELSE IF( LWORK.GE.WRKBL+( LDA + N )*N ) THEN * * WORK(IU) is LDA by N and WORK(IR) is N by N * @@ -1562,14 +1662,14 @@ IWORK = ITAU + N * * Compute A=Q*R, copying result to U -* (Workspace: need 2*N*N+2*N, prefer 2*N*N+N+N*NB) +* (Workspace: need 2*N*N + 2*N, prefer 2*N*N + N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'L', M, N, A, LDA, U, LDU ) * * Generate Q in U -* (Workspace: need 2*N*N+N+M, prefer 2*N*N+N+M*NB) +* (Workspace: need 2*N*N + N + M, prefer 2*N*N + N + M*NB) * CALL DORGQR( M, M, N, U, LDU, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1587,7 +1687,7 @@ * * Bidiagonalize R in WORK(IU), copying result to * WORK(IR) -* (Workspace: need 2*N*N+4*N, +* (Workspace: need 2*N*N + 4*N, * prefer 2*N*N+3*N+2*N*NB) * CALL DGEBRD( N, N, WORK( IU ), LDWRKU, S, @@ -1598,14 +1698,14 @@ $ WORK( IR ), LDWRKR ) * * Generate left bidiagonalizing vectors in WORK(IU) -* (Workspace: need 2*N*N+4*N, prefer 2*N*N+3*N+N*NB) +* (Workspace: need 2*N*N + 4*N, prefer 2*N*N + 3*N + N*NB) * CALL DORGBR( 'Q', N, N, N, WORK( IU ), LDWRKU, $ WORK( ITAUQ ), WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right bidiagonalizing vectors in WORK(IR) -* (Workspace: need 2*N*N+4*N-1, +* (Workspace: need 2*N*N + 4*N-1, * prefer 2*N*N+3*N+(N-1)*NB) * CALL DORGBR( 'P', N, N, N, WORK( IR ), LDWRKR, @@ -1616,7 +1716,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of R in WORK(IU) and computing * right singular vectors of R in WORK(IR) -* (Workspace: need 2*N*N+BDSPAC) +* (Workspace: need 2*N*N + BDSPAC) * CALL DBDSQR( 'U', N, N, N, 0, S, WORK( IE ), $ WORK( IR ), LDWRKR, WORK( IU ), @@ -1646,14 +1746,14 @@ IWORK = ITAU + N * * Compute A=Q*R, copying result to U -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'L', M, N, A, LDA, U, LDU ) * * Generate Q in U -* (Workspace: need N+M, prefer N+M*NB) +* (Workspace: need N + M, prefer N + M*NB) * CALL DORGQR( M, M, N, U, LDU, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1664,11 +1764,13 @@ * * Zero out below R in A * - CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, A( 2, 1 ), - $ LDA ) + IF( N .GT. 1 ) THEN + CALL DLASET( 'L', N-1, N-1, ZERO, ZERO, + $ A( 2, 1 ), LDA ) + END IF * * Bidiagonalize R in A -* (Workspace: need 4*N, prefer 3*N+2*N*NB) +* (Workspace: need 4*N, prefer 3*N + 2*N*NB) * CALL DGEBRD( N, N, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -1676,14 +1778,14 @@ * * Multiply Q in U by left bidiagonalizing vectors * in A -* (Workspace: need 3*N+M, prefer 3*N+M*NB) +* (Workspace: need 3*N + M, prefer 3*N + M*NB) * CALL DORMBR( 'Q', 'R', 'N', M, N, N, A, LDA, $ WORK( ITAUQ ), U, LDU, WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right bidiagonalizing vectors in A -* (Workspace: need 4*N-1, prefer 3*N+(N-1)*NB) +* (Workspace: need 4*N-1, prefer 3*N + (N-1)*NB) * CALL DORGBR( 'P', N, N, N, A, LDA, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1727,14 +1829,14 @@ IWORK = ITAU + N * * Compute A=Q*R, copying result to U -* (Workspace: need N*N+2*N, prefer N*N+N+N*NB) +* (Workspace: need N*N + 2*N, prefer N*N + N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'L', M, N, A, LDA, U, LDU ) * * Generate Q in U -* (Workspace: need N*N+N+M, prefer N*N+N+M*NB) +* (Workspace: need N*N + N + M, prefer N*N + N + M*NB) * CALL DORGQR( M, M, N, U, LDU, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1751,7 +1853,7 @@ IWORK = ITAUP + N * * Bidiagonalize R in WORK(IU), copying result to VT -* (Workspace: need N*N+4*N, prefer N*N+3*N+2*N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + 2*N*NB) * CALL DGEBRD( N, N, WORK( IU ), LDWRKU, S, $ WORK( IE ), WORK( ITAUQ ), @@ -1761,14 +1863,14 @@ $ LDVT ) * * Generate left bidiagonalizing vectors in WORK(IU) -* (Workspace: need N*N+4*N, prefer N*N+3*N+N*NB) +* (Workspace: need N*N + 4*N, prefer N*N + 3*N + N*NB) * CALL DORGBR( 'Q', N, N, N, WORK( IU ), LDWRKU, $ WORK( ITAUQ ), WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right bidiagonalizing vectors in VT -* (Workspace: need N*N+4*N-1, +* (Workspace: need N*N + 4*N-1, * prefer N*N+3*N+(N-1)*NB) * CALL DORGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), @@ -1778,7 +1880,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of R in WORK(IU) and computing * right singular vectors of R in VT -* (Workspace: need N*N+BDSPAC) +* (Workspace: need N*N + BDSPAC) * CALL DBDSQR( 'U', N, N, N, 0, S, WORK( IE ), VT, $ LDVT, WORK( IU ), LDWRKU, DUM, 1, @@ -1803,14 +1905,14 @@ IWORK = ITAU + N * * Compute A=Q*R, copying result to U -* (Workspace: need 2*N, prefer N+N*NB) +* (Workspace: need 2*N, prefer N + N*NB) * CALL DGEQRF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'L', M, N, A, LDA, U, LDU ) * * Generate Q in U -* (Workspace: need N+M, prefer N+M*NB) +* (Workspace: need N + M, prefer N + M*NB) * CALL DORGQR( M, M, N, U, LDU, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1827,7 +1929,7 @@ IWORK = ITAUP + N * * Bidiagonalize R in VT -* (Workspace: need 4*N, prefer 3*N+2*N*NB) +* (Workspace: need 4*N, prefer 3*N + 2*N*NB) * CALL DGEBRD( N, N, VT, LDVT, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -1835,14 +1937,14 @@ * * Multiply Q in U by left bidiagonalizing vectors * in VT -* (Workspace: need 3*N+M, prefer 3*N+M*NB) +* (Workspace: need 3*N + M, prefer 3*N + M*NB) * CALL DORMBR( 'Q', 'R', 'N', M, N, N, VT, LDVT, $ WORK( ITAUQ ), U, LDU, WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate right bidiagonalizing vectors in VT -* (Workspace: need 4*N-1, prefer 3*N+(N-1)*NB) +* (Workspace: need 4*N-1, prefer 3*N + (N-1)*NB) * CALL DORGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1876,7 +1978,7 @@ IWORK = ITAUP + N * * Bidiagonalize A -* (Workspace: need 3*N+M, prefer 3*N+(M+N)*NB) +* (Workspace: need 3*N + M, prefer 3*N + (M + N)*NB) * CALL DGEBRD( M, N, A, LDA, S, WORK( IE ), WORK( ITAUQ ), $ WORK( ITAUP ), WORK( IWORK ), LWORK-IWORK+1, @@ -1885,7 +1987,7 @@ * * If left singular vectors desired in U, copy result to U * and generate left bidiagonalizing vectors in U -* (Workspace: need 3*N+NCU, prefer 3*N+NCU*NB) +* (Workspace: need 3*N + NCU, prefer 3*N + NCU*NB) * CALL DLACPY( 'L', M, N, A, LDA, U, LDU ) IF( WNTUS ) @@ -1899,7 +2001,7 @@ * * If right singular vectors desired in VT, copy result to * VT and generate right bidiagonalizing vectors in VT -* (Workspace: need 4*N-1, prefer 3*N+(N-1)*NB) +* (Workspace: need 4*N-1, prefer 3*N + (N-1)*NB) * CALL DLACPY( 'U', N, N, A, LDA, VT, LDVT ) CALL DORGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), @@ -1909,7 +2011,7 @@ * * If left singular vectors desired in A, generate left * bidiagonalizing vectors in A -* (Workspace: need 4*N, prefer 3*N+N*NB) +* (Workspace: need 4*N, prefer 3*N + N*NB) * CALL DORGBR( 'Q', M, N, N, A, LDA, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1918,7 +2020,7 @@ * * If right singular vectors desired in A, generate right * bidiagonalizing vectors in A -* (Workspace: need 4*N-1, prefer 3*N+(N-1)*NB) +* (Workspace: need 4*N-1, prefer 3*N + (N-1)*NB) * CALL DORGBR( 'P', N, N, N, A, LDA, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -1980,7 +2082,7 @@ IWORK = ITAU + M * * Compute A=L*Q -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), WORK( IWORK ), $ LWORK-IWORK+1, IERR ) @@ -1994,7 +2096,7 @@ IWORK = ITAUP + M * * Bidiagonalize L in A -* (Workspace: need 4*M, prefer 3*M+2*M*NB) +* (Workspace: need 4*M, prefer 3*M + 2*M*NB) * CALL DGEBRD( M, M, A, LDA, S, WORK( IE ), WORK( ITAUQ ), $ WORK( ITAUP ), WORK( IWORK ), LWORK-IWORK+1, @@ -2002,7 +2104,7 @@ IF( WNTUO .OR. WNTUAS ) THEN * * If left singular vectors desired, generate Q -* (Workspace: need 4*M, prefer 3*M+M*NB) +* (Workspace: need 4*M, prefer 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, A, LDA, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2035,14 +2137,14 @@ * Sufficient workspace for a fast algorithm * IR = 1 - IF( LWORK.GE.MAX( WRKBL, LDA*N+M )+LDA*M ) THEN + IF( LWORK.GE.MAX( WRKBL, LDA*N + M ) + LDA*M ) THEN * * WORK(IU) is LDA by N and WORK(IR) is LDA by M * LDWRKU = LDA CHUNK = N LDWRKR = LDA - ELSE IF( LWORK.GE.MAX( WRKBL, LDA*N+M )+M*M ) THEN + ELSE IF( LWORK.GE.MAX( WRKBL, LDA*N + M ) + M*M ) THEN * * WORK(IU) is LDA by N and WORK(IR) is M by M * @@ -2061,7 +2163,7 @@ IWORK = ITAU + M * * Compute A=L*Q -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2073,7 +2175,7 @@ $ WORK( IR+LDWRKR ), LDWRKR ) * * Generate Q in A -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DORGLQ( M, N, M, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2083,14 +2185,14 @@ IWORK = ITAUP + M * * Bidiagonalize L in WORK(IR) -* (Workspace: need M*M+4*M, prefer M*M+3*M+2*M*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + 2*M*NB) * CALL DGEBRD( M, M, WORK( IR ), LDWRKR, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Generate right vectors bidiagonalizing L -* (Workspace: need M*M+4*M-1, prefer M*M+3*M+(M-1)*NB) +* (Workspace: need M*M + 4*M-1, prefer M*M + 3*M + (M-1)*NB) * CALL DORGBR( 'P', M, M, M, WORK( IR ), LDWRKR, $ WORK( ITAUP ), WORK( IWORK ), @@ -2099,7 +2201,7 @@ * * Perform bidiagonal QR iteration, computing right * singular vectors of L in WORK(IR) -* (Workspace: need M*M+BDSPAC) +* (Workspace: need M*M + BDSPAC) * CALL DBDSQR( 'U', M, M, 0, 0, S, WORK( IE ), $ WORK( IR ), LDWRKR, DUM, 1, DUM, 1, @@ -2108,7 +2210,7 @@ * * Multiply right singular vectors of L in WORK(IR) by Q * in A, storing result in WORK(IU) and copying to A -* (Workspace: need M*M+2*M, prefer M*M+M*N+M) +* (Workspace: need M*M + 2*M, prefer M*M + M*N + M) * DO 30 I = 1, N, CHUNK BLK = MIN( N-I+1, CHUNK ) @@ -2129,14 +2231,14 @@ IWORK = ITAUP + M * * Bidiagonalize A -* (Workspace: need 3*M+N, prefer 3*M+(M+N)*NB) +* (Workspace: need 3*M + N, prefer 3*M + (M + N)*NB) * CALL DGEBRD( M, N, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Generate right vectors bidiagonalizing A -* (Workspace: need 4*M, prefer 3*M+M*NB) +* (Workspace: need 4*M, prefer 3*M + M*NB) * CALL DORGBR( 'P', M, N, M, A, LDA, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2162,14 +2264,14 @@ * Sufficient workspace for a fast algorithm * IR = 1 - IF( LWORK.GE.MAX( WRKBL, LDA*N+M )+LDA*M ) THEN + IF( LWORK.GE.MAX( WRKBL, LDA*N + M ) + LDA*M ) THEN * * WORK(IU) is LDA by N and WORK(IR) is LDA by M * LDWRKU = LDA CHUNK = N LDWRKR = LDA - ELSE IF( LWORK.GE.MAX( WRKBL, LDA*N+M )+M*M ) THEN + ELSE IF( LWORK.GE.MAX( WRKBL, LDA*N + M ) + M*M ) THEN * * WORK(IU) is LDA by N and WORK(IR) is M by M * @@ -2188,7 +2290,7 @@ IWORK = ITAU + M * * Compute A=L*Q -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2200,7 +2302,7 @@ $ LDU ) * * Generate Q in A -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DORGLQ( M, N, M, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2210,7 +2312,7 @@ IWORK = ITAUP + M * * Bidiagonalize L in U, copying result to WORK(IR) -* (Workspace: need M*M+4*M, prefer M*M+3*M+2*M*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + 2*M*NB) * CALL DGEBRD( M, M, U, LDU, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -2218,14 +2320,14 @@ CALL DLACPY( 'U', M, M, U, LDU, WORK( IR ), LDWRKR ) * * Generate right vectors bidiagonalizing L in WORK(IR) -* (Workspace: need M*M+4*M-1, prefer M*M+3*M+(M-1)*NB) +* (Workspace: need M*M + 4*M-1, prefer M*M + 3*M + (M-1)*NB) * CALL DORGBR( 'P', M, M, M, WORK( IR ), LDWRKR, $ WORK( ITAUP ), WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate left vectors bidiagonalizing L in U -* (Workspace: need M*M+4*M, prefer M*M+3*M+M*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, U, LDU, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2234,7 +2336,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of L in U, and computing right * singular vectors of L in WORK(IR) -* (Workspace: need M*M+BDSPAC) +* (Workspace: need M*M + BDSPAC) * CALL DBDSQR( 'U', M, M, M, 0, S, WORK( IE ), $ WORK( IR ), LDWRKR, U, LDU, DUM, 1, @@ -2243,7 +2345,7 @@ * * Multiply right singular vectors of L in WORK(IR) by Q * in A, storing result in WORK(IU) and copying to A -* (Workspace: need M*M+2*M, prefer M*M+M*N+M)) +* (Workspace: need M*M + 2*M, prefer M*M + M*N + M)) * DO 40 I = 1, N, CHUNK BLK = MIN( N-I+1, CHUNK ) @@ -2262,7 +2364,7 @@ IWORK = ITAU + M * * Compute A=L*Q -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2274,7 +2376,7 @@ $ LDU ) * * Generate Q in A -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DORGLQ( M, N, M, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2284,21 +2386,21 @@ IWORK = ITAUP + M * * Bidiagonalize L in U -* (Workspace: need 4*M, prefer 3*M+2*M*NB) +* (Workspace: need 4*M, prefer 3*M + 2*M*NB) * CALL DGEBRD( M, M, U, LDU, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Multiply right vectors bidiagonalizing L by Q in A -* (Workspace: need 3*M+N, prefer 3*M+N*NB) +* (Workspace: need 3*M + N, prefer 3*M + N*NB) * CALL DORMBR( 'P', 'L', 'T', M, N, M, U, LDU, $ WORK( ITAUP ), A, LDA, WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate left vectors bidiagonalizing L in U -* (Workspace: need 4*M, prefer 3*M+M*NB) +* (Workspace: need 4*M, prefer 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, U, LDU, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2342,7 +2444,7 @@ IWORK = ITAU + M * * Compute A=L*Q -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2355,7 +2457,7 @@ $ WORK( IR+LDWRKR ), LDWRKR ) * * Generate Q in A -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DORGLQ( M, N, M, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2365,7 +2467,7 @@ IWORK = ITAUP + M * * Bidiagonalize L in WORK(IR) -* (Workspace: need M*M+4*M, prefer M*M+3*M+2*M*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + 2*M*NB) * CALL DGEBRD( M, M, WORK( IR ), LDWRKR, S, $ WORK( IE ), WORK( ITAUQ ), @@ -2374,7 +2476,7 @@ * * Generate right vectors bidiagonalizing L in * WORK(IR) -* (Workspace: need M*M+4*M, prefer M*M+3*M+(M-1)*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + (M-1)*NB) * CALL DORGBR( 'P', M, M, M, WORK( IR ), LDWRKR, $ WORK( ITAUP ), WORK( IWORK ), @@ -2383,7 +2485,7 @@ * * Perform bidiagonal QR iteration, computing right * singular vectors of L in WORK(IR) -* (Workspace: need M*M+BDSPAC) +* (Workspace: need M*M + BDSPAC) * CALL DBDSQR( 'U', M, M, 0, 0, S, WORK( IE ), $ WORK( IR ), LDWRKR, DUM, 1, DUM, 1, @@ -2404,7 +2506,7 @@ IWORK = ITAU + M * * Compute A=L*Q -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2414,7 +2516,7 @@ CALL DLACPY( 'U', M, N, A, LDA, VT, LDVT ) * * Generate Q in VT -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DORGLQ( M, N, M, VT, LDVT, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2429,14 +2531,14 @@ $ LDA ) * * Bidiagonalize L in A -* (Workspace: need 4*M, prefer 3*M+2*M*NB) +* (Workspace: need 4*M, prefer 3*M + 2*M*NB) * CALL DGEBRD( M, M, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Multiply right vectors bidiagonalizing L by Q in VT -* (Workspace: need 3*M+N, prefer 3*M+N*NB) +* (Workspace: need 3*M + N, prefer 3*M + N*NB) * CALL DORMBR( 'P', 'L', 'T', M, N, M, A, LDA, $ WORK( ITAUP ), VT, LDVT, @@ -2471,7 +2573,7 @@ LDWRKU = LDA IR = IU + LDWRKU*M LDWRKR = LDA - ELSE IF( LWORK.GE.WRKBL+( LDA+M )*M ) THEN + ELSE IF( LWORK.GE.WRKBL+( LDA + M )*M ) THEN * * WORK(IU) is LDA by M and WORK(IR) is M by M * @@ -2490,7 +2592,7 @@ IWORK = ITAU + M * * Compute A=L*Q -* (Workspace: need 2*M*M+2*M, prefer 2*M*M+M+M*NB) +* (Workspace: need 2*M*M + 2*M, prefer 2*M*M + M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2503,7 +2605,7 @@ $ WORK( IU+LDWRKU ), LDWRKU ) * * Generate Q in A -* (Workspace: need 2*M*M+2*M, prefer 2*M*M+M+M*NB) +* (Workspace: need 2*M*M + 2*M, prefer 2*M*M + M + M*NB) * CALL DORGLQ( M, N, M, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2514,7 +2616,7 @@ * * Bidiagonalize L in WORK(IU), copying result to * WORK(IR) -* (Workspace: need 2*M*M+4*M, +* (Workspace: need 2*M*M + 4*M, * prefer 2*M*M+3*M+2*M*NB) * CALL DGEBRD( M, M, WORK( IU ), LDWRKU, S, @@ -2525,7 +2627,7 @@ $ WORK( IR ), LDWRKR ) * * Generate right bidiagonalizing vectors in WORK(IU) -* (Workspace: need 2*M*M+4*M-1, +* (Workspace: need 2*M*M + 4*M-1, * prefer 2*M*M+3*M+(M-1)*NB) * CALL DORGBR( 'P', M, M, M, WORK( IU ), LDWRKU, @@ -2533,7 +2635,7 @@ $ LWORK-IWORK+1, IERR ) * * Generate left bidiagonalizing vectors in WORK(IR) -* (Workspace: need 2*M*M+4*M, prefer 2*M*M+3*M+M*NB) +* (Workspace: need 2*M*M + 4*M, prefer 2*M*M + 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, WORK( IR ), LDWRKR, $ WORK( ITAUQ ), WORK( IWORK ), @@ -2543,7 +2645,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of L in WORK(IR) and computing * right singular vectors of L in WORK(IU) -* (Workspace: need 2*M*M+BDSPAC) +* (Workspace: need 2*M*M + BDSPAC) * CALL DBDSQR( 'U', M, M, M, 0, S, WORK( IE ), $ WORK( IU ), LDWRKU, WORK( IR ), @@ -2570,14 +2672,14 @@ IWORK = ITAU + M * * Compute A=L*Q, copying result to VT -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'U', M, N, A, LDA, VT, LDVT ) * * Generate Q in VT -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DORGLQ( M, N, M, VT, LDVT, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2592,21 +2694,21 @@ $ LDA ) * * Bidiagonalize L in A -* (Workspace: need 4*M, prefer 3*M+2*M*NB) +* (Workspace: need 4*M, prefer 3*M + 2*M*NB) * CALL DGEBRD( M, M, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Multiply right vectors bidiagonalizing L by Q in VT -* (Workspace: need 3*M+N, prefer 3*M+N*NB) +* (Workspace: need 3*M + N, prefer 3*M + N*NB) * CALL DORMBR( 'P', 'L', 'T', M, N, M, A, LDA, $ WORK( ITAUP ), VT, LDVT, $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Generate left bidiagonalizing vectors of L in A -* (Workspace: need 4*M, prefer 3*M+M*NB) +* (Workspace: need 4*M, prefer 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, A, LDA, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2650,7 +2752,7 @@ IWORK = ITAU + M * * Compute A=L*Q -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2663,7 +2765,7 @@ $ WORK( IU+LDWRKU ), LDWRKU ) * * Generate Q in A -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DORGLQ( M, N, M, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2673,7 +2775,7 @@ IWORK = ITAUP + M * * Bidiagonalize L in WORK(IU), copying result to U -* (Workspace: need M*M+4*M, prefer M*M+3*M+2*M*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + 2*M*NB) * CALL DGEBRD( M, M, WORK( IU ), LDWRKU, S, $ WORK( IE ), WORK( ITAUQ ), @@ -2683,7 +2785,7 @@ $ LDU ) * * Generate right bidiagonalizing vectors in WORK(IU) -* (Workspace: need M*M+4*M-1, +* (Workspace: need M*M + 4*M-1, * prefer M*M+3*M+(M-1)*NB) * CALL DORGBR( 'P', M, M, M, WORK( IU ), LDWRKU, @@ -2691,7 +2793,7 @@ $ LWORK-IWORK+1, IERR ) * * Generate left bidiagonalizing vectors in U -* (Workspace: need M*M+4*M, prefer M*M+3*M+M*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, U, LDU, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2700,7 +2802,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of L in U and computing right * singular vectors of L in WORK(IU) -* (Workspace: need M*M+BDSPAC) +* (Workspace: need M*M + BDSPAC) * CALL DBDSQR( 'U', M, M, M, 0, S, WORK( IE ), $ WORK( IU ), LDWRKU, U, LDU, DUM, 1, @@ -2721,14 +2823,14 @@ IWORK = ITAU + M * * Compute A=L*Q, copying result to VT -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'U', M, N, A, LDA, VT, LDVT ) * * Generate Q in VT -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DORGLQ( M, N, M, VT, LDVT, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2744,7 +2846,7 @@ IWORK = ITAUP + M * * Bidiagonalize L in U -* (Workspace: need 4*M, prefer 3*M+2*M*NB) +* (Workspace: need 4*M, prefer 3*M + 2*M*NB) * CALL DGEBRD( M, M, U, LDU, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -2752,14 +2854,14 @@ * * Multiply right bidiagonalizing vectors in U by Q * in VT -* (Workspace: need 3*M+N, prefer 3*M+N*NB) +* (Workspace: need 3*M + N, prefer 3*M + N*NB) * CALL DORMBR( 'P', 'L', 'T', M, N, M, U, LDU, $ WORK( ITAUP ), VT, LDVT, $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Generate left bidiagonalizing vectors in U -* (Workspace: need 4*M, prefer 3*M+M*NB) +* (Workspace: need 4*M, prefer 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, U, LDU, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2786,7 +2888,7 @@ * N right singular vectors to be computed in VT and * no left singular vectors to be computed * - IF( LWORK.GE.M*M+MAX( N+M, 4*M, BDSPAC ) ) THEN + IF( LWORK.GE.M*M+MAX( N + M, 4*M, BDSPAC ) ) THEN * * Sufficient workspace for a fast algorithm * @@ -2806,7 +2908,7 @@ IWORK = ITAU + M * * Compute A=L*Q, copying result to VT -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2820,7 +2922,7 @@ $ WORK( IR+LDWRKR ), LDWRKR ) * * Generate Q in VT -* (Workspace: need M*M+M+N, prefer M*M+M+N*NB) +* (Workspace: need M*M + M + N, prefer M*M + M + N*NB) * CALL DORGLQ( N, N, M, VT, LDVT, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2830,7 +2932,7 @@ IWORK = ITAUP + M * * Bidiagonalize L in WORK(IR) -* (Workspace: need M*M+4*M, prefer M*M+3*M+2*M*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + 2*M*NB) * CALL DGEBRD( M, M, WORK( IR ), LDWRKR, S, $ WORK( IE ), WORK( ITAUQ ), @@ -2838,7 +2940,7 @@ $ LWORK-IWORK+1, IERR ) * * Generate right bidiagonalizing vectors in WORK(IR) -* (Workspace: need M*M+4*M-1, +* (Workspace: need M*M + 4*M-1, * prefer M*M+3*M+(M-1)*NB) * CALL DORGBR( 'P', M, M, M, WORK( IR ), LDWRKR, @@ -2848,7 +2950,7 @@ * * Perform bidiagonal QR iteration, computing right * singular vectors of L in WORK(IR) -* (Workspace: need M*M+BDSPAC) +* (Workspace: need M*M + BDSPAC) * CALL DBDSQR( 'U', M, M, 0, 0, S, WORK( IE ), $ WORK( IR ), LDWRKR, DUM, 1, DUM, 1, @@ -2873,14 +2975,14 @@ IWORK = ITAU + M * * Compute A=L*Q, copying result to VT -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'U', M, N, A, LDA, VT, LDVT ) * * Generate Q in VT -* (Workspace: need M+N, prefer M+N*NB) +* (Workspace: need M + N, prefer M + N*NB) * CALL DORGLQ( N, N, M, VT, LDVT, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2895,7 +2997,7 @@ $ LDA ) * * Bidiagonalize L in A -* (Workspace: need 4*M, prefer 3*M+2*M*NB) +* (Workspace: need 4*M, prefer 3*M + 2*M*NB) * CALL DGEBRD( M, M, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -2903,7 +3005,7 @@ * * Multiply right bidiagonalizing vectors in A by Q * in VT -* (Workspace: need 3*M+N, prefer 3*M+N*NB) +* (Workspace: need 3*M + N, prefer 3*M + N*NB) * CALL DORMBR( 'P', 'L', 'T', M, N, M, A, LDA, $ WORK( ITAUP ), VT, LDVT, @@ -2926,7 +3028,7 @@ * N right singular vectors to be computed in VT and * M left singular vectors to be overwritten on A * - IF( LWORK.GE.2*M*M+MAX( N+M, 4*M, BDSPAC ) ) THEN + IF( LWORK.GE.2*M*M+MAX( N + M, 4*M, BDSPAC ) ) THEN * * Sufficient workspace for a fast algorithm * @@ -2938,7 +3040,7 @@ LDWRKU = LDA IR = IU + LDWRKU*M LDWRKR = LDA - ELSE IF( LWORK.GE.WRKBL+( LDA+M )*M ) THEN + ELSE IF( LWORK.GE.WRKBL+( LDA + M )*M ) THEN * * WORK(IU) is LDA by M and WORK(IR) is M by M * @@ -2957,14 +3059,14 @@ IWORK = ITAU + M * * Compute A=L*Q, copying result to VT -* (Workspace: need 2*M*M+2*M, prefer 2*M*M+M+M*NB) +* (Workspace: need 2*M*M + 2*M, prefer 2*M*M + M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'U', M, N, A, LDA, VT, LDVT ) * * Generate Q in VT -* (Workspace: need 2*M*M+M+N, prefer 2*M*M+M+N*NB) +* (Workspace: need 2*M*M + M + N, prefer 2*M*M + M + N*NB) * CALL DORGLQ( N, N, M, VT, LDVT, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -2982,7 +3084,7 @@ * * Bidiagonalize L in WORK(IU), copying result to * WORK(IR) -* (Workspace: need 2*M*M+4*M, +* (Workspace: need 2*M*M + 4*M, * prefer 2*M*M+3*M+2*M*NB) * CALL DGEBRD( M, M, WORK( IU ), LDWRKU, S, @@ -2993,7 +3095,7 @@ $ WORK( IR ), LDWRKR ) * * Generate right bidiagonalizing vectors in WORK(IU) -* (Workspace: need 2*M*M+4*M-1, +* (Workspace: need 2*M*M + 4*M-1, * prefer 2*M*M+3*M+(M-1)*NB) * CALL DORGBR( 'P', M, M, M, WORK( IU ), LDWRKU, @@ -3001,7 +3103,7 @@ $ LWORK-IWORK+1, IERR ) * * Generate left bidiagonalizing vectors in WORK(IR) -* (Workspace: need 2*M*M+4*M, prefer 2*M*M+3*M+M*NB) +* (Workspace: need 2*M*M + 4*M, prefer 2*M*M + 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, WORK( IR ), LDWRKR, $ WORK( ITAUQ ), WORK( IWORK ), @@ -3011,7 +3113,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of L in WORK(IR) and computing * right singular vectors of L in WORK(IU) -* (Workspace: need 2*M*M+BDSPAC) +* (Workspace: need 2*M*M + BDSPAC) * CALL DBDSQR( 'U', M, M, M, 0, S, WORK( IE ), $ WORK( IU ), LDWRKU, WORK( IR ), @@ -3041,14 +3143,14 @@ IWORK = ITAU + M * * Compute A=L*Q, copying result to VT -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'U', M, N, A, LDA, VT, LDVT ) * * Generate Q in VT -* (Workspace: need M+N, prefer M+N*NB) +* (Workspace: need M + N, prefer M + N*NB) * CALL DORGLQ( N, N, M, VT, LDVT, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -3063,7 +3165,7 @@ $ LDA ) * * Bidiagonalize L in A -* (Workspace: need 4*M, prefer 3*M+2*M*NB) +* (Workspace: need 4*M, prefer 3*M + 2*M*NB) * CALL DGEBRD( M, M, A, LDA, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -3071,14 +3173,14 @@ * * Multiply right bidiagonalizing vectors in A by Q * in VT -* (Workspace: need 3*M+N, prefer 3*M+N*NB) +* (Workspace: need 3*M + N, prefer 3*M + N*NB) * CALL DORMBR( 'P', 'L', 'T', M, N, M, A, LDA, $ WORK( ITAUP ), VT, LDVT, $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Generate left bidiagonalizing vectors in A -* (Workspace: need 4*M, prefer 3*M+M*NB) +* (Workspace: need 4*M, prefer 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, A, LDA, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -3102,7 +3204,7 @@ * N right singular vectors to be computed in VT and * M left singular vectors to be computed in U * - IF( LWORK.GE.M*M+MAX( N+M, 4*M, BDSPAC ) ) THEN + IF( LWORK.GE.M*M+MAX( N + M, 4*M, BDSPAC ) ) THEN * * Sufficient workspace for a fast algorithm * @@ -3122,14 +3224,14 @@ IWORK = ITAU + M * * Compute A=L*Q, copying result to VT -* (Workspace: need M*M+2*M, prefer M*M+M+M*NB) +* (Workspace: need M*M + 2*M, prefer M*M + M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'U', M, N, A, LDA, VT, LDVT ) * * Generate Q in VT -* (Workspace: need M*M+M+N, prefer M*M+M+N*NB) +* (Workspace: need M*M + M + N, prefer M*M + M + N*NB) * CALL DORGLQ( N, N, M, VT, LDVT, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -3146,7 +3248,7 @@ IWORK = ITAUP + M * * Bidiagonalize L in WORK(IU), copying result to U -* (Workspace: need M*M+4*M, prefer M*M+3*M+2*M*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + 2*M*NB) * CALL DGEBRD( M, M, WORK( IU ), LDWRKU, S, $ WORK( IE ), WORK( ITAUQ ), @@ -3156,14 +3258,14 @@ $ LDU ) * * Generate right bidiagonalizing vectors in WORK(IU) -* (Workspace: need M*M+4*M, prefer M*M+3*M+(M-1)*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + (M-1)*NB) * CALL DORGBR( 'P', M, M, M, WORK( IU ), LDWRKU, $ WORK( ITAUP ), WORK( IWORK ), $ LWORK-IWORK+1, IERR ) * * Generate left bidiagonalizing vectors in U -* (Workspace: need M*M+4*M, prefer M*M+3*M+M*NB) +* (Workspace: need M*M + 4*M, prefer M*M + 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, U, LDU, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -3172,7 +3274,7 @@ * Perform bidiagonal QR iteration, computing left * singular vectors of L in U and computing right * singular vectors of L in WORK(IU) -* (Workspace: need M*M+BDSPAC) +* (Workspace: need M*M + BDSPAC) * CALL DBDSQR( 'U', M, M, M, 0, S, WORK( IE ), $ WORK( IU ), LDWRKU, U, LDU, DUM, 1, @@ -3197,14 +3299,14 @@ IWORK = ITAU + M * * Compute A=L*Q, copying result to VT -* (Workspace: need 2*M, prefer M+M*NB) +* (Workspace: need 2*M, prefer M + M*NB) * CALL DGELQF( M, N, A, LDA, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) CALL DLACPY( 'U', M, N, A, LDA, VT, LDVT ) * * Generate Q in VT -* (Workspace: need M+N, prefer M+N*NB) +* (Workspace: need M + N, prefer M + N*NB) * CALL DORGLQ( N, N, M, VT, LDVT, WORK( ITAU ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -3220,7 +3322,7 @@ IWORK = ITAUP + M * * Bidiagonalize L in U -* (Workspace: need 4*M, prefer 3*M+2*M*NB) +* (Workspace: need 4*M, prefer 3*M + 2*M*NB) * CALL DGEBRD( M, M, U, LDU, S, WORK( IE ), $ WORK( ITAUQ ), WORK( ITAUP ), @@ -3228,14 +3330,14 @@ * * Multiply right bidiagonalizing vectors in U by Q * in VT -* (Workspace: need 3*M+N, prefer 3*M+N*NB) +* (Workspace: need 3*M + N, prefer 3*M + N*NB) * CALL DORMBR( 'P', 'L', 'T', M, N, M, U, LDU, $ WORK( ITAUP ), VT, LDVT, $ WORK( IWORK ), LWORK-IWORK+1, IERR ) * * Generate left bidiagonalizing vectors in U -* (Workspace: need 4*M, prefer 3*M+M*NB) +* (Workspace: need 4*M, prefer 3*M + M*NB) * CALL DORGBR( 'Q', M, M, M, U, LDU, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -3269,7 +3371,7 @@ IWORK = ITAUP + M * * Bidiagonalize A -* (Workspace: need 3*M+N, prefer 3*M+(M+N)*NB) +* (Workspace: need 3*M + N, prefer 3*M + (M + N)*NB) * CALL DGEBRD( M, N, A, LDA, S, WORK( IE ), WORK( ITAUQ ), $ WORK( ITAUP ), WORK( IWORK ), LWORK-IWORK+1, @@ -3278,7 +3380,7 @@ * * If left singular vectors desired in U, copy result to U * and generate left bidiagonalizing vectors in U -* (Workspace: need 4*M-1, prefer 3*M+(M-1)*NB) +* (Workspace: need 4*M-1, prefer 3*M + (M-1)*NB) * CALL DLACPY( 'L', M, M, A, LDA, U, LDU ) CALL DORGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ), @@ -3288,7 +3390,7 @@ * * If right singular vectors desired in VT, copy result to * VT and generate right bidiagonalizing vectors in VT -* (Workspace: need 3*M+NRVT, prefer 3*M+NRVT*NB) +* (Workspace: need 3*M + NRVT, prefer 3*M + NRVT*NB) * CALL DLACPY( 'U', M, N, A, LDA, VT, LDVT ) IF( WNTVA ) @@ -3302,7 +3404,7 @@ * * If left singular vectors desired in A, generate left * bidiagonalizing vectors in A -* (Workspace: need 4*M-1, prefer 3*M+(M-1)*NB) +* (Workspace: need 4*M-1, prefer 3*M + (M-1)*NB) * CALL DORGBR( 'Q', M, M, N, A, LDA, WORK( ITAUQ ), $ WORK( IWORK ), LWORK-IWORK+1, IERR ) @@ -3311,7 +3413,7 @@ * * If right singular vectors desired in A, generate right * bidiagonalizing vectors in A -* (Workspace: need 4*M, prefer 3*M+M*NB) +* (Workspace: need 4*M, prefer 3*M + M*NB) * CALL DORGBR( 'P', M, N, M, A, LDA, WORK( ITAUP ), $ WORK( IWORK ), LWORK-IWORK+1, IERR )