rpl/lapack/lapack/zgesdd.f - diff

Return to zgesdd.f CVS log

Up to [local] / rpl / lapack / lapack

Diff for /rpl/lapack/lapack/zgesdd.f between versions 1.8 and 1.20

version 1.8, 2010/12/21 13:53:44	version 1.20, 2023/08/07 08:39:19
Line 1	Line 1
SUBROUTINE ZGESDD( JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK,	*> \brief \b ZGESDD
$ LWORK, RWORK, IWORK, INFO )
*	*
* -- LAPACK driver routine (version 3.2.2) --	* =========== DOCUMENTATION ===========
	*
	* Online html documentation available at
	* http://www.netlib.org/lapack/explore-html/
	*
	*> \htmlonly
	*> Download ZGESDD + dependencies
	*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgesdd.f">
	*> [TGZ]</a>
	*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgesdd.f">
	*> [ZIP]</a>
	*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgesdd.f">
	*> [TXT]</a>
	*> \endhtmlonly
	*
	* Definition:
	* ===========
	*
	* SUBROUTINE ZGESDD( JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT,
	* WORK, LWORK, RWORK, IWORK, INFO )
	*
	* .. Scalar Arguments ..
	* CHARACTER JOBZ
	* INTEGER INFO, LDA, LDU, LDVT, LWORK, M, N
	* ..
	* .. Array Arguments ..
	* INTEGER IWORK( * )
	* DOUBLE PRECISION RWORK( * ), S( * )
	* COMPLEX16 A( LDA, ), U( LDU, * ), VT( LDVT, * ),
	* $ WORK( * )
	* ..
	*
	*
	*> \par Purpose:
	* =============
	*>
	*> \verbatim
	*>
	*> ZGESDD computes the singular value decomposition (SVD) of a complex
	*> M-by-N matrix A, optionally computing the left and/or right singular
	*> vectors, by using divide-and-conquer method. The SVD is written
	*>
	> A = U SIGMA * conjugate-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.
	*>
	> Note that the routine returns VT = V*H, not V.
	*>
	*> The divide and conquer algorithm makes very mild assumptions about
	*> floating point arithmetic. It will work on machines with a guard
	*> digit in add/subtract, or on those binary machines without guard
	*> digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
	*> Cray-2. It could conceivably fail on hexadecimal or decimal machines
	*> without guard digits, but we know of none.
	*> \endverbatim
	*
	* Arguments:
	* ==========
	*
	*> \param[in] JOBZ
	*> \verbatim
	> JOBZ is CHARACTER1
	*> Specifies options for computing all or part of the matrix U:
	> = 'A': all M columns of U and all N rows of V*H are
	*> returned in the arrays U and VT;
	*> = 'S': the first min(M,N) columns of U and the first
	> min(M,N) rows of V*H are returned in the arrays U
	*> and VT;
	*> = 'O': If M >= N, the first N columns of U are overwritten
	> in the array A and all rows of V*H are returned in
	*> the array VT;
	*> otherwise, all columns of U are returned in the
	> array U and the first M rows of V*H are overwritten
	*> in the array A;
	> = 'N': no columns of U or rows of V*H are computed.
	*> \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 COMPLEX16 array, dimension (LDA,N)
	*> On entry, the M-by-N matrix A.
	*> On exit,
	*> if JOBZ = 'O', A is overwritten with the first N columns
	*> of U (the left singular vectors, stored
	*> columnwise) if M >= N;
	*> A is overwritten with the first M rows
	> of V*H (the right singular vectors, stored
	*> rowwise) otherwise.
	*> if JOBZ .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 COMPLEX16 array, dimension (LDU,UCOL)
	*> UCOL = M if JOBZ = 'A' or JOBZ = 'O' and M < N;
	*> UCOL = min(M,N) if JOBZ = 'S'.
	*> If JOBZ = 'A' or JOBZ = 'O' and M < N, U contains the M-by-M
	*> unitary matrix U;
	*> if JOBZ = 'S', U contains the first min(M,N) columns of U
	*> (the left singular vectors, stored columnwise);
	*> if JOBZ = 'O' and M >= N, or JOBZ = 'N', U is not referenced.
	*> \endverbatim
	*>
	*> \param[in] LDU
	*> \verbatim
	*> LDU is INTEGER
	*> The leading dimension of the array U. LDU >= 1;
	*> if JOBZ = 'S' or 'A' or JOBZ = 'O' and M < N, LDU >= M.
	*> \endverbatim
	*>
	*> \param[out] VT
	*> \verbatim
	> VT is COMPLEX16 array, dimension (LDVT,N)
	*> If JOBZ = 'A' or JOBZ = 'O' and M >= N, VT contains the
	> N-by-N unitary matrix V*H;
	*> if JOBZ = 'S', VT contains the first min(M,N) rows of
	> V*H (the right singular vectors, stored rowwise);
	*> if JOBZ = 'O' and M < N, or JOBZ = 'N', VT is not referenced.
	*> \endverbatim
	*>
	*> \param[in] LDVT
	*> \verbatim
	*> LDVT is INTEGER
	*> The leading dimension of the array VT. LDVT >= 1;
	*> if JOBZ = 'A' or JOBZ = 'O' and M >= N, LDVT >= N;
	*> if JOBZ = 'S', LDVT >= min(M,N).
	*> \endverbatim
	*>
	*> \param[out] WORK
	*> \verbatim
	> WORK is COMPLEX16 array, dimension (MAX(1,LWORK))
	*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
	*> \endverbatim
	*>
	*> \param[in] LWORK
	*> \verbatim
	*> LWORK is INTEGER
	*> The dimension of the array WORK. LWORK >= 1.
	*> If LWORK = -1, a workspace query is assumed. The optimal
	*> size for the WORK array is calculated and stored in WORK(1),
	*> and no other work except argument checking is performed.
	*>
	*> Let mx = max(M,N) and mn = min(M,N).
	> If JOBZ = 'N', LWORK >= 2mn + mx.
	> If JOBZ = 'O', LWORK >= 2mnmn + 2mn + mx.
	> If JOBZ = 'S', LWORK >= mnmn + 3*mn.
	> If JOBZ = 'A', LWORK >= mnmn + 2*mn + mx.
	*> These are not tight minimums in all cases; see comments inside code.
	*> For good performance, LWORK should generally be larger;
	*> a query is recommended.
	*> \endverbatim
	*>
	*> \param[out] RWORK
	*> \verbatim
	*> RWORK is DOUBLE PRECISION array, dimension (MAX(1,LRWORK))
	*> Let mx = max(M,N) and mn = min(M,N).
	> If JOBZ = 'N', LRWORK >= 5mn (LAPACK <= 3.6 needs 7*mn);
	> else if mx >> mn, LRWORK >= 5mnmn + 5mn;
	> else LRWORK >= max( 5mnmn + 5mn,
	> 2mxmn + 2mn*mn + mn ).
	*> \endverbatim
	*>
	*> \param[out] IWORK
	*> \verbatim
	> IWORK is INTEGER array, dimension (8min(M,N))
	*> \endverbatim
	*>
	*> \param[out] INFO
	*> \verbatim
	*> INFO is INTEGER
	*> < 0: if INFO = -i, the i-th argument had an illegal value.
	*> = -4: if A had a NAN entry.
	*> > 0: The updating process of DBDSDC did not converge.
	*> = 0: successful exit.
	*> \endverbatim
	*
	* Authors:
	* ========
	*
	*> \author Univ. of Tennessee
	*> \author Univ. of California Berkeley
	*> \author Univ. of Colorado Denver
	*> \author NAG Ltd.
	*
	*> \ingroup complex16GEsing
	*
	*> \par Contributors:
	* ==================
	*>
	*> Ming Gu and Huan Ren, Computer Science Division, University of
	*> California at Berkeley, USA
	*>
	* =====================================================================
	SUBROUTINE ZGESDD( JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT,
	$ WORK, LWORK, RWORK, IWORK, INFO )
	implicit none
	*
	* -- LAPACK driver routine --
* -- LAPACK is a software package provided by Univ. of Tennessee, --	* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--	* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* June 2010
* 8-15-00: Improve consistency of WS calculations (eca)
*	*
* .. Scalar Arguments ..	* .. Scalar Arguments ..
CHARACTER JOBZ	CHARACTER JOBZ
Line 18	Line 241
$ WORK( * )	$ WORK( * )
* ..	* ..
*	*
* Purpose
* =======
*
* ZGESDD computes the singular value decomposition (SVD) of a complex
* M-by-N matrix A, optionally computing the left and/or right singular
* vectors, by using divide-and-conquer method. The SVD is written
*
* A = U * SIGMA * conjugate-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.
*
* Note that the routine returns VT = V**H, not V.
*
* The divide and conquer algorithm makes very mild assumptions about
* floating point arithmetic. It will work on machines with a guard
* digit in add/subtract, or on those binary machines without guard
* digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
* Cray-2. It could conceivably fail on hexadecimal or decimal machines
* without guard digits, but we know of none.
*
* Arguments
* =========
*
* JOBZ (input) CHARACTER*1
* Specifies options for computing all or part of the matrix U:
* = 'A': all M columns of U and all N rows of V**H are
* returned in the arrays U and VT;
* = 'S': the first min(M,N) columns of U and the first
* min(M,N) rows of V**H are returned in the arrays U
* and VT;
* = 'O': If M >= N, the first N columns of U are overwritten
* in the array A and all rows of V**H are returned in
* the array VT;
* otherwise, all columns of U are returned in the
* array U and the first M rows of V**H are overwritten
* in the array A;
* = 'N': no columns of U or rows of V**H are computed.
*
* 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) COMPLEX*16 array, dimension (LDA,N)
* On entry, the M-by-N matrix A.
* On exit,
* if JOBZ = 'O', A is overwritten with the first N columns
* of U (the left singular vectors, stored
* columnwise) if M >= N;
* A is overwritten with the first M rows
* of V**H (the right singular vectors, stored
* rowwise) otherwise.
* if JOBZ .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) COMPLEX*16 array, dimension (LDU,UCOL)
* UCOL = M if JOBZ = 'A' or JOBZ = 'O' and M < N;
* UCOL = min(M,N) if JOBZ = 'S'.
* If JOBZ = 'A' or JOBZ = 'O' and M < N, U contains the M-by-M
* unitary matrix U;
* if JOBZ = 'S', U contains the first min(M,N) columns of U
* (the left singular vectors, stored columnwise);
* if JOBZ = 'O' and M >= N, or JOBZ = 'N', U is not referenced.
*
* LDU (input) INTEGER
* The leading dimension of the array U. LDU >= 1; if
* JOBZ = 'S' or 'A' or JOBZ = 'O' and M < N, LDU >= M.
*
* VT (output) COMPLEX*16 array, dimension (LDVT,N)
* If JOBZ = 'A' or JOBZ = 'O' and M >= N, VT contains the
* N-by-N unitary matrix V**H;
* if JOBZ = 'S', VT contains the first min(M,N) rows of
* V**H (the right singular vectors, stored rowwise);
* if JOBZ = 'O' and M < N, or JOBZ = 'N', VT is not referenced.
*
* LDVT (input) INTEGER
* The leading dimension of the array VT. LDVT >= 1; if
* JOBZ = 'A' or JOBZ = 'O' and M >= N, LDVT >= N;
* if JOBZ = 'S', LDVT >= min(M,N).
*
* WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
* On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
*
* LWORK (input) INTEGER
* The dimension of the array WORK. LWORK >= 1.
* if JOBZ = 'N', LWORK >= 2*min(M,N)+max(M,N).
* if JOBZ = 'O',
* LWORK >= 2min(M,N)min(M,N)+2*min(M,N)+max(M,N).
* if JOBZ = 'S' or 'A',
* LWORK >= min(M,N)min(M,N)+2min(M,N)+max(M,N).
* For good performance, LWORK should generally be larger.
*
* If LWORK = -1, a workspace query is assumed. The optimal
* size for the WORK array is calculated and stored in WORK(1),
* and no other work except argument checking is performed.
*
* RWORK (workspace) DOUBLE PRECISION array, dimension (MAX(1,LRWORK))
* If JOBZ = 'N', LRWORK >= 5*min(M,N).
* Otherwise,
* LRWORK >= min(M,N)max(5min(M,N)+7,2max(M,N)+2min(M,N)+1)
*
* IWORK (workspace) INTEGER array, dimension (8*min(M,N))
*
* INFO (output) INTEGER
* = 0: successful exit.
* < 0: if INFO = -i, the i-th argument had an illegal value.
* > 0: The updating process of DBDSDC did not converge.
*
* Further Details
* ===============
*
* Based on contributions by
* Ming Gu and Huan Ren, Computer Science Division, University of
* California at Berkeley, USA
*
* =====================================================================	* =====================================================================
*	*
* .. Parameters ..	* .. Parameters ..
INTEGER LQUERV
PARAMETER ( LQUERV = -1 )
COMPLEX*16 CZERO, CONE	COMPLEX*16 CZERO, CONE
PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ),	PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ),
$ CONE = ( 1.0D+0, 0.0D+0 ) )	$ CONE = ( 1.0D+0, 0.0D+0 ) )
Line 156	Line 251
PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )	PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
* ..	* ..
* .. Local Scalars ..	* .. Local Scalars ..
LOGICAL WNTQA, WNTQAS, WNTQN, WNTQO, WNTQS	LOGICAL LQUERY, WNTQA, WNTQAS, WNTQN, WNTQO, WNTQS
INTEGER BLK, CHUNK, I, IE, IERR, IL, IR, IRU, IRVT,	INTEGER BLK, CHUNK, I, IE, IERR, IL, IR, IRU, IRVT,
$ ISCL, ITAU, ITAUP, ITAUQ, IU, IVT, LDWKVT,	$ ISCL, ITAU, ITAUP, ITAUQ, IU, IVT, LDWKVT,
$ LDWRKL, LDWRKR, LDWRKU, MAXWRK, MINMN, MINWRK,	$ LDWRKL, LDWRKR, LDWRKU, MAXWRK, MINMN, MINWRK,
$ MNTHR1, MNTHR2, NRWORK, NWORK, WRKBL	$ MNTHR1, MNTHR2, NRWORK, NWORK, WRKBL
	INTEGER LWORK_ZGEBRD_MN, LWORK_ZGEBRD_MM,
	$ LWORK_ZGEBRD_NN, LWORK_ZGELQF_MN,
	$ LWORK_ZGEQRF_MN,
	$ LWORK_ZUNGBR_P_MN, LWORK_ZUNGBR_P_NN,
	$ LWORK_ZUNGBR_Q_MN, LWORK_ZUNGBR_Q_MM,
	$ LWORK_ZUNGLQ_MN, LWORK_ZUNGLQ_NN,
	$ LWORK_ZUNGQR_MM, LWORK_ZUNGQR_MN,
	$ LWORK_ZUNMBR_PRC_MM, LWORK_ZUNMBR_QLN_MM,
	$ LWORK_ZUNMBR_PRC_MN, LWORK_ZUNMBR_QLN_MN,
	$ LWORK_ZUNMBR_PRC_NN, LWORK_ZUNMBR_QLN_NN
DOUBLE PRECISION ANRM, BIGNUM, EPS, SMLNUM	DOUBLE PRECISION ANRM, BIGNUM, EPS, SMLNUM
* ..	* ..
* .. Local Arrays ..	* .. Local Arrays ..
INTEGER IDUM( 1 )	INTEGER IDUM( 1 )
DOUBLE PRECISION DUM( 1 )	DOUBLE PRECISION DUM( 1 )
	COMPLEX*16 CDUM( 1 )
* ..	* ..
* .. External Subroutines ..	* .. External Subroutines ..
EXTERNAL DBDSDC, DLASCL, XERBLA, ZGEBRD, ZGELQF, ZGEMM,	EXTERNAL DBDSDC, DLASCL, XERBLA, ZGEBRD, ZGELQF, ZGEMM,
Line 173	Line 279
$ ZLASET, ZUNGBR, ZUNGLQ, ZUNGQR, ZUNMBR	$ ZLASET, ZUNGBR, ZUNGLQ, ZUNGQR, ZUNMBR
* ..	* ..
* .. External Functions ..	* .. External Functions ..
LOGICAL LSAME	LOGICAL LSAME, DISNAN
INTEGER ILAENV	DOUBLE PRECISION DLAMCH, ZLANGE, DROUNDUP_LWORK
DOUBLE PRECISION DLAMCH, ZLANGE	EXTERNAL LSAME, DLAMCH, ZLANGE, DISNAN,
EXTERNAL LSAME, ILAENV, DLAMCH, ZLANGE	$ DROUNDUP_LWORK
* ..	* ..
* .. Intrinsic Functions ..	* .. Intrinsic Functions ..
INTRINSIC INT, MAX, MIN, SQRT	INTRINSIC INT, MAX, MIN, SQRT
Line 185	Line 291
*	*
* Test the input arguments	* Test the input arguments
*	*
INFO = 0	INFO = 0
MINMN = MIN( M, N )	MINMN = MIN( M, N )
MNTHR1 = INT( MINMN*17.0D0 / 9.0D0 )	MNTHR1 = INT( MINMN*17.0D0 / 9.0D0 )
MNTHR2 = INT( MINMN*5.0D0 / 3.0D0 )	MNTHR2 = INT( MINMN*5.0D0 / 3.0D0 )
WNTQA = LSAME( JOBZ, 'A' )	WNTQA = LSAME( JOBZ, 'A' )
WNTQS = LSAME( JOBZ, 'S' )	WNTQS = LSAME( JOBZ, 'S' )
WNTQAS = WNTQA .OR. WNTQS	WNTQAS = WNTQA .OR. WNTQS
WNTQO = LSAME( JOBZ, 'O' )	WNTQO = LSAME( JOBZ, 'O' )
WNTQN = LSAME( JOBZ, 'N' )	WNTQN = LSAME( JOBZ, 'N' )
	LQUERY = ( LWORK.EQ.-1 )
MINWRK = 1	MINWRK = 1
MAXWRK = 1	MAXWRK = 1
*	*
Line 215	Line 322
END IF	END IF
*	*
* Compute workspace	* Compute workspace
* (Note: Comments in the code beginning "Workspace:" describe the	* Note: Comments in the code beginning "Workspace:" describe the
* minimal amount of workspace needed at that point in the code,	* minimal amount of workspace allocated at that point in the code,
* as well as the preferred amount for good performance.	* as well as the preferred amount for good performance.
* CWorkspace refers to complex workspace, and RWorkspace to	* CWorkspace refers to complex workspace, and RWorkspace to
* real workspace. NB refers to the optimal block size for the	* real workspace. NB refers to the optimal block size for the
* immediately following subroutine, as returned by ILAENV.)	* immediately following subroutine, as returned by ILAENV.)
*	*
IF( INFO.EQ.0 .AND. M.GT.0 .AND. N.GT.0 ) THEN	IF( INFO.EQ.0 ) THEN
IF( M.GE.N ) THEN	MINWRK = 1
	MAXWRK = 1
	IF( M.GE.N .AND. MINMN.GT.0 ) THEN
*	*
* There is no complex work space needed for bidiagonal SVD	* There is no complex work space needed for bidiagonal SVD
* The real work space needed for bidiagonal SVD is BDSPAC	* The real work space needed for bidiagonal SVD (dbdsdc) is
* for computing singular values and singular vectors; BDSPAN	* BDSPAC = 3NN + 4*N for singular values and vectors;
* for computing singular values only.	* BDSPAC = 4*N for singular values only;
* BDSPAC = 5NN + 7*N	* not including e, RU, and RVT matrices.
* BDSPAN = MAX(7N+4, 3N+2+SMLSIZ*(SMLSIZ+8))	*
	* Compute space preferred for each routine
	CALL ZGEBRD( M, N, CDUM(1), M, DUM(1), DUM(1), CDUM(1),
	$ CDUM(1), CDUM(1), -1, IERR )
	LWORK_ZGEBRD_MN = INT( CDUM(1) )
	*
	CALL ZGEBRD( N, N, CDUM(1), N, DUM(1), DUM(1), CDUM(1),
	$ CDUM(1), CDUM(1), -1, IERR )
	LWORK_ZGEBRD_NN = INT( CDUM(1) )
	*
	CALL ZGEQRF( M, N, CDUM(1), M, CDUM(1), CDUM(1), -1, IERR )
	LWORK_ZGEQRF_MN = INT( CDUM(1) )
	*
	CALL ZUNGBR( 'P', N, N, N, CDUM(1), N, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGBR_P_NN = INT( CDUM(1) )
	*
	CALL ZUNGBR( 'Q', M, M, N, CDUM(1), M, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGBR_Q_MM = INT( CDUM(1) )
	*
	CALL ZUNGBR( 'Q', M, N, N, CDUM(1), M, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGBR_Q_MN = INT( CDUM(1) )
	*
	CALL ZUNGQR( M, M, N, CDUM(1), M, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGQR_MM = INT( CDUM(1) )
	*
	CALL ZUNGQR( M, N, N, CDUM(1), M, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGQR_MN = INT( CDUM(1) )
	*
	CALL ZUNMBR( 'P', 'R', 'C', N, N, N, CDUM(1), N, CDUM(1),
	$ CDUM(1), N, CDUM(1), -1, IERR )
	LWORK_ZUNMBR_PRC_NN = INT( CDUM(1) )
	*
	CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, CDUM(1), M, CDUM(1),
	$ CDUM(1), M, CDUM(1), -1, IERR )
	LWORK_ZUNMBR_QLN_MM = INT( CDUM(1) )
	*
	CALL ZUNMBR( 'Q', 'L', 'N', M, N, N, CDUM(1), M, CDUM(1),
	$ CDUM(1), M, CDUM(1), -1, IERR )
	LWORK_ZUNMBR_QLN_MN = INT( CDUM(1) )
	*
	CALL ZUNMBR( 'Q', 'L', 'N', N, N, N, CDUM(1), N, CDUM(1),
	$ CDUM(1), N, CDUM(1), -1, IERR )
	LWORK_ZUNMBR_QLN_NN = INT( CDUM(1) )
*	*
IF( M.GE.MNTHR1 ) THEN	IF( M.GE.MNTHR1 ) THEN
IF( WNTQN ) THEN	IF( WNTQN ) THEN
*	*
* Path 1 (M much larger than N, JOBZ='N')	* Path 1 (M >> N, JOBZ='N')
*	*
MAXWRK = N + N*ILAENV( 1, 'ZGEQRF', ' ', M, N, -1,	MAXWRK = N + LWORK_ZGEQRF_MN
$ -1 )	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZGEBRD_NN )
MAXWRK = MAX( MAXWRK, 2N+2N*
$ ILAENV( 1, 'ZGEBRD', ' ', N, N, -1, -1 ) )
MINWRK = 3*N	MINWRK = 3*N
ELSE IF( WNTQO ) THEN	ELSE IF( WNTQO ) THEN
*	*
* Path 2 (M much larger than N, JOBZ='O')	* Path 2 (M >> N, JOBZ='O')
*	*
WRKBL = N + N*ILAENV( 1, 'ZGEQRF', ' ', M, N, -1, -1 )	WRKBL = N + LWORK_ZGEQRF_MN
WRKBL = MAX( WRKBL, N+N*ILAENV( 1, 'ZUNGQR', ' ', M,	WRKBL = MAX( WRKBL, N + LWORK_ZUNGQR_MN )
$ N, N, -1 ) )	WRKBL = MAX( WRKBL, 2*N + LWORK_ZGEBRD_NN )
WRKBL = MAX( WRKBL, 2N+2N*	WRKBL = MAX( WRKBL, 2*N + LWORK_ZUNMBR_QLN_NN )
$ ILAENV( 1, 'ZGEBRD', ' ', N, N, -1, -1 ) )	WRKBL = MAX( WRKBL, 2*N + LWORK_ZUNMBR_PRC_NN )
WRKBL = MAX( WRKBL, 2N+N
$ ILAENV( 1, 'ZUNMBR', 'QLN', N, N, N, -1 ) )
WRKBL = MAX( WRKBL, 2N+N
$ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) )
MAXWRK = MN + NN + WRKBL	MAXWRK = MN + NN + WRKBL
MINWRK = 2NN + 3*N	MINWRK = 2NN + 3*N
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
*	*
* Path 3 (M much larger than N, JOBZ='S')	* Path 3 (M >> N, JOBZ='S')
*	*
WRKBL = N + N*ILAENV( 1, 'ZGEQRF', ' ', M, N, -1, -1 )	WRKBL = N + LWORK_ZGEQRF_MN
WRKBL = MAX( WRKBL, N+N*ILAENV( 1, 'ZUNGQR', ' ', M,	WRKBL = MAX( WRKBL, N + LWORK_ZUNGQR_MN )
$ N, N, -1 ) )	WRKBL = MAX( WRKBL, 2*N + LWORK_ZGEBRD_NN )
WRKBL = MAX( WRKBL, 2N+2N*	WRKBL = MAX( WRKBL, 2*N + LWORK_ZUNMBR_QLN_NN )
$ ILAENV( 1, 'ZGEBRD', ' ', N, N, -1, -1 ) )	WRKBL = MAX( WRKBL, 2*N + LWORK_ZUNMBR_PRC_NN )
WRKBL = MAX( WRKBL, 2N+N
$ ILAENV( 1, 'ZUNMBR', 'QLN', N, N, N, -1 ) )
WRKBL = MAX( WRKBL, 2N+N
$ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) )
MAXWRK = N*N + WRKBL	MAXWRK = N*N + WRKBL
MINWRK = NN + 3N	MINWRK = NN + 3N
ELSE IF( WNTQA ) THEN	ELSE IF( WNTQA ) THEN
*	*
* Path 4 (M much larger than N, JOBZ='A')	* Path 4 (M >> N, JOBZ='A')
*	*
WRKBL = N + N*ILAENV( 1, 'ZGEQRF', ' ', M, N, -1, -1 )	WRKBL = N + LWORK_ZGEQRF_MN
WRKBL = MAX( WRKBL, N+M*ILAENV( 1, 'ZUNGQR', ' ', M,	WRKBL = MAX( WRKBL, N + LWORK_ZUNGQR_MM )
$ M, N, -1 ) )	WRKBL = MAX( WRKBL, 2*N + LWORK_ZGEBRD_NN )
WRKBL = MAX( WRKBL, 2N+2N*	WRKBL = MAX( WRKBL, 2*N + LWORK_ZUNMBR_QLN_NN )
$ ILAENV( 1, 'ZGEBRD', ' ', N, N, -1, -1 ) )	WRKBL = MAX( WRKBL, 2*N + LWORK_ZUNMBR_PRC_NN )
WRKBL = MAX( WRKBL, 2N+N
$ ILAENV( 1, 'ZUNMBR', 'QLN', N, N, N, -1 ) )
WRKBL = MAX( WRKBL, 2N+N
$ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) )
MAXWRK = N*N + WRKBL	MAXWRK = N*N + WRKBL
MINWRK = NN + 2N + M	MINWRK = NN + MAX( 3N, N + M )
END IF	END IF
ELSE IF( M.GE.MNTHR2 ) THEN	ELSE IF( M.GE.MNTHR2 ) THEN
*	*
* Path 5 (M much larger than N, but not as much as MNTHR1)	* Path 5 (M >> N, but not as much as MNTHR1)
*	*
MAXWRK = 2N + ( M+N )ILAENV( 1, 'ZGEBRD', ' ', M, N,	MAXWRK = 2*N + LWORK_ZGEBRD_MN
$ -1, -1 )
MINWRK = 2*N + M	MINWRK = 2*N + M
IF( WNTQO ) THEN	IF( WNTQO ) THEN
MAXWRK = MAX( MAXWRK, 2N+N	* Path 5o (M >> N, JOBZ='O')
$ ILAENV( 1, 'ZUNGBR', 'P', N, N, N, -1 ) )	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNGBR_P_NN )
MAXWRK = MAX( MAXWRK, 2N+N	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNGBR_Q_MN )
$ ILAENV( 1, 'ZUNGBR', 'Q', M, N, N, -1 ) )
MAXWRK = MAXWRK + M*N	MAXWRK = MAXWRK + M*N
MINWRK = MINWRK + N*N	MINWRK = MINWRK + N*N
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
MAXWRK = MAX( MAXWRK, 2N+N	* Path 5s (M >> N, JOBZ='S')
$ ILAENV( 1, 'ZUNGBR', 'P', N, N, N, -1 ) )	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNGBR_P_NN )
MAXWRK = MAX( MAXWRK, 2N+N	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNGBR_Q_MN )
$ ILAENV( 1, 'ZUNGBR', 'Q', M, N, N, -1 ) )
ELSE IF( WNTQA ) THEN	ELSE IF( WNTQA ) THEN
MAXWRK = MAX( MAXWRK, 2N+N	* Path 5a (M >> N, JOBZ='A')
$ ILAENV( 1, 'ZUNGBR', 'P', N, N, N, -1 ) )	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNGBR_P_NN )
MAXWRK = MAX( MAXWRK, 2N+M	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNGBR_Q_MM )
$ ILAENV( 1, 'ZUNGBR', 'Q', M, M, N, -1 ) )
END IF	END IF
ELSE	ELSE
*	*
* Path 6 (M at least N, but not much larger)	* Path 6 (M >= N, but not much larger)
*	*
MAXWRK = 2N + ( M+N )ILAENV( 1, 'ZGEBRD', ' ', M, N,	MAXWRK = 2*N + LWORK_ZGEBRD_MN
$ -1, -1 )
MINWRK = 2*N + M	MINWRK = 2*N + M
IF( WNTQO ) THEN	IF( WNTQO ) THEN
MAXWRK = MAX( MAXWRK, 2N+N	* Path 6o (M >= N, JOBZ='O')
$ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) )	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNMBR_PRC_NN )
MAXWRK = MAX( MAXWRK, 2N+N	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNMBR_QLN_MN )
$ ILAENV( 1, 'ZUNMBR', 'QLN', M, N, N, -1 ) )
MAXWRK = MAXWRK + M*N	MAXWRK = MAXWRK + M*N
MINWRK = MINWRK + N*N	MINWRK = MINWRK + N*N
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
MAXWRK = MAX( MAXWRK, 2N+N	* Path 6s (M >= N, JOBZ='S')
$ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) )	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNMBR_QLN_MN )
MAXWRK = MAX( MAXWRK, 2N+N	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNMBR_PRC_NN )
$ ILAENV( 1, 'ZUNMBR', 'QLN', M, N, N, -1 ) )
ELSE IF( WNTQA ) THEN	ELSE IF( WNTQA ) THEN
MAXWRK = MAX( MAXWRK, 2N+N	* Path 6a (M >= N, JOBZ='A')
$ ILAENV( 1, 'ZUNGBR', 'PRC', N, N, N, -1 ) )	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNMBR_QLN_MM )
MAXWRK = MAX( MAXWRK, 2N+M	MAXWRK = MAX( MAXWRK, 2*N + LWORK_ZUNMBR_PRC_NN )
$ ILAENV( 1, 'ZUNGBR', 'QLN', M, M, N, -1 ) )
END IF	END IF
END IF	END IF
ELSE	ELSE IF( MINMN.GT.0 ) THEN
*	*
* There is no complex work space needed for bidiagonal SVD	* There is no complex work space needed for bidiagonal SVD
* The real work space needed for bidiagonal SVD is BDSPAC	* The real work space needed for bidiagonal SVD (dbdsdc) is
* for computing singular values and singular vectors; BDSPAN	* BDSPAC = 3MM + 4*M for singular values and vectors;
* for computing singular values only.	* BDSPAC = 4*M for singular values only;
* BDSPAC = 5MM + 7*M	* not including e, RU, and RVT matrices.
* BDSPAN = MAX(7M+4, 3M+2+SMLSIZ*(SMLSIZ+8))	*
	* Compute space preferred for each routine
	CALL ZGEBRD( M, N, CDUM(1), M, DUM(1), DUM(1), CDUM(1),
	$ CDUM(1), CDUM(1), -1, IERR )
	LWORK_ZGEBRD_MN = INT( CDUM(1) )
	*
	CALL ZGEBRD( M, M, CDUM(1), M, DUM(1), DUM(1), CDUM(1),
	$ CDUM(1), CDUM(1), -1, IERR )
	LWORK_ZGEBRD_MM = INT( CDUM(1) )
	*
	CALL ZGELQF( M, N, CDUM(1), M, CDUM(1), CDUM(1), -1, IERR )
	LWORK_ZGELQF_MN = INT( CDUM(1) )
	*
	CALL ZUNGBR( 'P', M, N, M, CDUM(1), M, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGBR_P_MN = INT( CDUM(1) )
	*
	CALL ZUNGBR( 'P', N, N, M, CDUM(1), N, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGBR_P_NN = INT( CDUM(1) )
	*
	CALL ZUNGBR( 'Q', M, M, N, CDUM(1), M, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGBR_Q_MM = INT( CDUM(1) )
	*
	CALL ZUNGLQ( M, N, M, CDUM(1), M, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGLQ_MN = INT( CDUM(1) )
	*
	CALL ZUNGLQ( N, N, M, CDUM(1), N, CDUM(1), CDUM(1),
	$ -1, IERR )
	LWORK_ZUNGLQ_NN = INT( CDUM(1) )
	*
	CALL ZUNMBR( 'P', 'R', 'C', M, M, M, CDUM(1), M, CDUM(1),
	$ CDUM(1), M, CDUM(1), -1, IERR )
	LWORK_ZUNMBR_PRC_MM = INT( CDUM(1) )
	*
	CALL ZUNMBR( 'P', 'R', 'C', M, N, M, CDUM(1), M, CDUM(1),
	$ CDUM(1), M, CDUM(1), -1, IERR )
	LWORK_ZUNMBR_PRC_MN = INT( CDUM(1) )
	*
	CALL ZUNMBR( 'P', 'R', 'C', N, N, M, CDUM(1), N, CDUM(1),
	$ CDUM(1), N, CDUM(1), -1, IERR )
	LWORK_ZUNMBR_PRC_NN = INT( CDUM(1) )
	*
	CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, CDUM(1), M, CDUM(1),
	$ CDUM(1), M, CDUM(1), -1, IERR )
	LWORK_ZUNMBR_QLN_MM = INT( CDUM(1) )
*	*
IF( N.GE.MNTHR1 ) THEN	IF( N.GE.MNTHR1 ) THEN
IF( WNTQN ) THEN	IF( WNTQN ) THEN
*	*
* Path 1t (N much larger than M, JOBZ='N')	* Path 1t (N >> M, JOBZ='N')
*	*
MAXWRK = M + M*ILAENV( 1, 'ZGELQF', ' ', M, N, -1,	MAXWRK = M + LWORK_ZGELQF_MN
$ -1 )	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZGEBRD_MM )
MAXWRK = MAX( MAXWRK, 2M+2M*
$ ILAENV( 1, 'ZGEBRD', ' ', M, M, -1, -1 ) )
MINWRK = 3*M	MINWRK = 3*M
ELSE IF( WNTQO ) THEN	ELSE IF( WNTQO ) THEN
*	*
* Path 2t (N much larger than M, JOBZ='O')	* Path 2t (N >> M, JOBZ='O')
*	*
WRKBL = M + M*ILAENV( 1, 'ZGELQF', ' ', M, N, -1, -1 )	WRKBL = M + LWORK_ZGELQF_MN
WRKBL = MAX( WRKBL, M+M*ILAENV( 1, 'ZUNGLQ', ' ', M,	WRKBL = MAX( WRKBL, M + LWORK_ZUNGLQ_MN )
$ N, M, -1 ) )	WRKBL = MAX( WRKBL, 2*M + LWORK_ZGEBRD_MM )
WRKBL = MAX( WRKBL, 2M+2M*	WRKBL = MAX( WRKBL, 2*M + LWORK_ZUNMBR_QLN_MM )
$ ILAENV( 1, 'ZGEBRD', ' ', M, M, -1, -1 ) )	WRKBL = MAX( WRKBL, 2*M + LWORK_ZUNMBR_PRC_MM )
WRKBL = MAX( WRKBL, 2M+M
$ ILAENV( 1, 'ZUNMBR', 'PRC', M, M, M, -1 ) )
WRKBL = MAX( WRKBL, 2M+M
$ ILAENV( 1, 'ZUNMBR', 'QLN', M, M, M, -1 ) )
MAXWRK = MN + MM + WRKBL	MAXWRK = MN + MM + WRKBL
MINWRK = 2MM + 3*M	MINWRK = 2MM + 3*M
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
*	*
* Path 3t (N much larger than M, JOBZ='S')	* Path 3t (N >> M, JOBZ='S')
*	*
WRKBL = M + M*ILAENV( 1, 'ZGELQF', ' ', M, N, -1, -1 )	WRKBL = M + LWORK_ZGELQF_MN
WRKBL = MAX( WRKBL, M+M*ILAENV( 1, 'ZUNGLQ', ' ', M,	WRKBL = MAX( WRKBL, M + LWORK_ZUNGLQ_MN )
$ N, M, -1 ) )	WRKBL = MAX( WRKBL, 2*M + LWORK_ZGEBRD_MM )
WRKBL = MAX( WRKBL, 2M+2M*	WRKBL = MAX( WRKBL, 2*M + LWORK_ZUNMBR_QLN_MM )
$ ILAENV( 1, 'ZGEBRD', ' ', M, M, -1, -1 ) )	WRKBL = MAX( WRKBL, 2*M + LWORK_ZUNMBR_PRC_MM )
WRKBL = MAX( WRKBL, 2M+M
$ ILAENV( 1, 'ZUNMBR', 'PRC', M, M, M, -1 ) )
WRKBL = MAX( WRKBL, 2M+M
$ ILAENV( 1, 'ZUNMBR', 'QLN', M, M, M, -1 ) )
MAXWRK = M*M + WRKBL	MAXWRK = M*M + WRKBL
MINWRK = MM + 3M	MINWRK = MM + 3M
ELSE IF( WNTQA ) THEN	ELSE IF( WNTQA ) THEN
*	*
* Path 4t (N much larger than M, JOBZ='A')	* Path 4t (N >> M, JOBZ='A')
*	*
WRKBL = M + M*ILAENV( 1, 'ZGELQF', ' ', M, N, -1, -1 )	WRKBL = M + LWORK_ZGELQF_MN
WRKBL = MAX( WRKBL, M+N*ILAENV( 1, 'ZUNGLQ', ' ', N,	WRKBL = MAX( WRKBL, M + LWORK_ZUNGLQ_NN )
$ N, M, -1 ) )	WRKBL = MAX( WRKBL, 2*M + LWORK_ZGEBRD_MM )
WRKBL = MAX( WRKBL, 2M+2M*	WRKBL = MAX( WRKBL, 2*M + LWORK_ZUNMBR_QLN_MM )
$ ILAENV( 1, 'ZGEBRD', ' ', M, M, -1, -1 ) )	WRKBL = MAX( WRKBL, 2*M + LWORK_ZUNMBR_PRC_MM )
WRKBL = MAX( WRKBL, 2M+M
$ ILAENV( 1, 'ZUNMBR', 'PRC', M, M, M, -1 ) )
WRKBL = MAX( WRKBL, 2M+M
$ ILAENV( 1, 'ZUNMBR', 'QLN', M, M, M, -1 ) )
MAXWRK = M*M + WRKBL	MAXWRK = M*M + WRKBL
MINWRK = MM + 2M + N	MINWRK = MM + MAX( 3M, M + N )
END IF	END IF
ELSE IF( N.GE.MNTHR2 ) THEN	ELSE IF( N.GE.MNTHR2 ) THEN
*	*
* Path 5t (N much larger than M, but not as much as MNTHR1)	* Path 5t (N >> M, but not as much as MNTHR1)
*	*
MAXWRK = 2M + ( M+N )ILAENV( 1, 'ZGEBRD', ' ', M, N,	MAXWRK = 2*M + LWORK_ZGEBRD_MN
$ -1, -1 )
MINWRK = 2*M + N	MINWRK = 2*M + N
IF( WNTQO ) THEN	IF( WNTQO ) THEN
MAXWRK = MAX( MAXWRK, 2M+M	* Path 5to (N >> M, JOBZ='O')
$ ILAENV( 1, 'ZUNGBR', 'P', M, N, M, -1 ) )	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNGBR_Q_MM )
MAXWRK = MAX( MAXWRK, 2M+M	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNGBR_P_MN )
$ ILAENV( 1, 'ZUNGBR', 'Q', M, M, N, -1 ) )
MAXWRK = MAXWRK + M*N	MAXWRK = MAXWRK + M*N
MINWRK = MINWRK + M*M	MINWRK = MINWRK + M*M
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
MAXWRK = MAX( MAXWRK, 2M+M	* Path 5ts (N >> M, JOBZ='S')
$ ILAENV( 1, 'ZUNGBR', 'P', M, N, M, -1 ) )	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNGBR_Q_MM )
MAXWRK = MAX( MAXWRK, 2M+M	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNGBR_P_MN )
$ ILAENV( 1, 'ZUNGBR', 'Q', M, M, N, -1 ) )
ELSE IF( WNTQA ) THEN	ELSE IF( WNTQA ) THEN
MAXWRK = MAX( MAXWRK, 2M+N	* Path 5ta (N >> M, JOBZ='A')
$ ILAENV( 1, 'ZUNGBR', 'P', N, N, M, -1 ) )	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNGBR_Q_MM )
MAXWRK = MAX( MAXWRK, 2M+M	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNGBR_P_NN )
$ ILAENV( 1, 'ZUNGBR', 'Q', M, M, N, -1 ) )
END IF	END IF
ELSE	ELSE
*	*
* Path 6t (N greater than M, but not much larger)	* Path 6t (N > M, but not much larger)
*	*
MAXWRK = 2M + ( M+N )ILAENV( 1, 'ZGEBRD', ' ', M, N,	MAXWRK = 2*M + LWORK_ZGEBRD_MN
$ -1, -1 )
MINWRK = 2*M + N	MINWRK = 2*M + N
IF( WNTQO ) THEN	IF( WNTQO ) THEN
MAXWRK = MAX( MAXWRK, 2M+M	* Path 6to (N > M, JOBZ='O')
$ ILAENV( 1, 'ZUNMBR', 'PRC', M, N, M, -1 ) )	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNMBR_QLN_MM )
MAXWRK = MAX( MAXWRK, 2M+M	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNMBR_PRC_MN )
$ ILAENV( 1, 'ZUNMBR', 'QLN', M, M, N, -1 ) )
MAXWRK = MAXWRK + M*N	MAXWRK = MAXWRK + M*N
MINWRK = MINWRK + M*M	MINWRK = MINWRK + M*M
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
MAXWRK = MAX( MAXWRK, 2M+M	* Path 6ts (N > M, JOBZ='S')
$ ILAENV( 1, 'ZUNGBR', 'PRC', M, N, M, -1 ) )	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNMBR_QLN_MM )
MAXWRK = MAX( MAXWRK, 2M+M	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNMBR_PRC_MN )
$ ILAENV( 1, 'ZUNGBR', 'QLN', M, M, N, -1 ) )
ELSE IF( WNTQA ) THEN	ELSE IF( WNTQA ) THEN
MAXWRK = MAX( MAXWRK, 2M+N	* Path 6ta (N > M, JOBZ='A')
$ ILAENV( 1, 'ZUNGBR', 'PRC', N, N, M, -1 ) )	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNMBR_QLN_MM )
MAXWRK = MAX( MAXWRK, 2M+M	MAXWRK = MAX( MAXWRK, 2*M + LWORK_ZUNMBR_PRC_NN )
$ ILAENV( 1, 'ZUNGBR', 'QLN', M, M, N, -1 ) )
END IF	END IF
END IF	END IF
END IF	END IF
MAXWRK = MAX( MAXWRK, MINWRK )	MAXWRK = MAX( MAXWRK, MINWRK )
END IF	END IF
IF( INFO.EQ.0 ) THEN	IF( INFO.EQ.0 ) THEN
WORK( 1 ) = MAXWRK	WORK( 1 ) = DROUNDUP_LWORK( MAXWRK )
IF( LWORK.LT.MINWRK .AND. LWORK.NE.LQUERV )	IF( LWORK.LT.MINWRK .AND. .NOT. LQUERY ) THEN
$ INFO = -13	INFO = -12
	END IF
END IF	END IF
*	*
* Quick returns
*
IF( INFO.NE.0 ) THEN	IF( INFO.NE.0 ) THEN
CALL XERBLA( 'ZGESDD', -INFO )	CALL XERBLA( 'ZGESDD', -INFO )
RETURN	RETURN
	ELSE IF( LQUERY ) THEN
	RETURN
END IF	END IF
IF( LWORK.EQ.LQUERV )	*
$ RETURN	* Quick return if possible
	*
IF( M.EQ.0 .OR. N.EQ.0 ) THEN	IF( M.EQ.0 .OR. N.EQ.0 ) THEN
RETURN	RETURN
END IF	END IF
Line 485	Line 646
* Scale A if max element outside range [SMLNUM,BIGNUM]	* Scale A if max element outside range [SMLNUM,BIGNUM]
*	*
ANRM = ZLANGE( 'M', M, N, A, LDA, DUM )	ANRM = ZLANGE( 'M', M, N, A, LDA, DUM )
	IF( DISNAN( ANRM ) ) THEN
	INFO = -4
	RETURN
	END IF
ISCL = 0	ISCL = 0
IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN	IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN
ISCL = 1	ISCL = 1
Line 504	Line 669
*	*
IF( WNTQN ) THEN	IF( WNTQN ) THEN
*	*
* Path 1 (M much larger than N, JOBZ='N')	* Path 1 (M >> N, JOBZ='N')
* No singular vectors to be computed	* No singular vectors to be computed
*	*
ITAU = 1	ITAU = 1
NWORK = ITAU + N	NWORK = ITAU + N
*	*
* Compute A=Q*R	* Compute A=Q*R
* (CWorkspace: need 2N, prefer N+NNB)	* CWorkspace: need N [tau] + N [work]
* (RWorkspace: need 0)	* CWorkspace: prefer N [tau] + N*NB [work]
	* RWorkspace: need 0
*	*
CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),	CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
Line 527	Line 693
NWORK = ITAUP + N	NWORK = ITAUP + N
*	*
* Bidiagonalize R in A	* Bidiagonalize R in A
* (CWorkspace: need 3N, prefer 2N+2NNB)	* CWorkspace: need 2*N [tauq, taup] + N [work]
* (RWorkspace: need N)	* CWorkspace: prefer 2N [tauq, taup] + 2N*NB [work]
	* RWorkspace: need N [e]
*	*
CALL ZGEBRD( N, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),	CALL ZGEBRD( N, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),
$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,	$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,
Line 536	Line 703
NRWORK = IE + N	NRWORK = IE + N
*	*
* Perform bidiagonal SVD, compute singular values only	* Perform bidiagonal SVD, compute singular values only
* (CWorkspace: 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAN)	* RWorkspace: need N [e] + BDSPAC
*	*
CALL DBDSDC( 'U', 'N', N, S, RWORK( IE ), DUM, 1, DUM, 1,	CALL DBDSDC( 'U', 'N', N, S, RWORK( IE ), DUM,1,DUM,1,
$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )	$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )
*	*
ELSE IF( WNTQO ) THEN	ELSE IF( WNTQO ) THEN
*	*
* Path 2 (M much larger than N, JOBZ='O')	* Path 2 (M >> N, JOBZ='O')
* N left singular vectors to be overwritten on A and	* N left singular vectors to be overwritten on A and
* N right singular vectors to be computed in VT	* N right singular vectors to be computed in VT
*	*
Line 554	Line 721
*	*
LDWRKU = N	LDWRKU = N
IR = IU + LDWRKU*N	IR = IU + LDWRKU*N
IF( LWORK.GE.MN+NN+3*N ) THEN	IF( LWORK .GE. MN + NN + 3*N ) THEN
*	*
* WORK(IR) is M by N	* WORK(IR) is M by N
*	*
LDWRKR = M	LDWRKR = M
ELSE	ELSE
LDWRKR = ( LWORK-NN-3N ) / N	LDWRKR = ( LWORK - NN - 3N ) / N
END IF	END IF
ITAU = IR + LDWRKR*N	ITAU = IR + LDWRKR*N
NWORK = ITAU + N	NWORK = ITAU + N
*	*
* Compute A=Q*R	* Compute A=Q*R
* (CWorkspace: need NN+2N, prefer MN+N+NNB)	* CWorkspace: need NN [U] + NN [R] + N [tau] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [U] + NN [R] + N [tau] + N*NB [work]
	* RWorkspace: need 0
*	*
CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),	CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
Line 579	Line 747
$ LDWRKR )	$ LDWRKR )
*	*
* Generate Q in A	* Generate Q in A
* (CWorkspace: need 2N, prefer N+NNB)	* CWorkspace: need NN [U] + NN [R] + N [tau] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [U] + NN [R] + N [tau] + N*NB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGQR( M, N, N, A, LDA, WORK( ITAU ),	CALL ZUNGQR( M, N, N, A, LDA, WORK( ITAU ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
Line 590	Line 759
NWORK = ITAUP + N	NWORK = ITAUP + N
*	*
* Bidiagonalize R in WORK(IR)	* Bidiagonalize R in WORK(IR)
* (CWorkspace: need NN+3N, prefer MN+2N+2NNB)	* CWorkspace: need NN [U] + NN [R] + 2*N [tauq, taup] + N [work]
* (RWorkspace: need N)	* CWorkspace: prefer NN [U] + NN [R] + 2N [tauq, taup] + 2N*NB [work]
	* RWorkspace: need N [e]
*	*
CALL ZGEBRD( N, N, WORK( IR ), LDWRKR, S, RWORK( IE ),	CALL ZGEBRD( N, N, WORK( IR ), LDWRKR, S, RWORK( IE ),
$ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ),	$ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ),
Line 600	Line 770
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of R in WORK(IRU) and computing right singular vectors	* of R in WORK(IRU) and computing right singular vectors
* of R in WORK(IRVT)	* of R in WORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + BDSPAC
*	*
IRU = IE + N	IRU = IE + N
IRVT = IRU + N*N	IRVT = IRU + N*N
Line 612	Line 782
*	*
* Copy real matrix RWORK(IRU) to complex matrix WORK(IU)	* Copy real matrix RWORK(IRU) to complex matrix WORK(IU)
* Overwrite WORK(IU) by the left singular vectors of R	* Overwrite WORK(IU) by the left singular vectors of R
* (CWorkspace: need 2NN+3N, prefer MN+NN+2N+N*NB)	* CWorkspace: need NN [U] + NN [R] + 2*N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [U] + NN [R] + 2N [tauq, taup] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, WORK( IU ),	CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, WORK( IU ),
$ LDWRKU )	$ LDWRKU )
Line 623	Line 794
*	*
* Copy real matrix RWORK(IRVT) to complex matrix VT	* Copy real matrix RWORK(IRVT) to complex matrix VT
* Overwrite VT by the right singular vectors of R	* Overwrite VT by the right singular vectors of R
* (CWorkspace: need NN+3N, prefer MN+2N+N*NB)	* CWorkspace: need NN [U] + NN [R] + 2*N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [U] + NN [R] + 2N [tauq, taup] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )	CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )
CALL ZUNMBR( 'P', 'R', 'C', N, N, N, WORK( IR ), LDWRKR,	CALL ZUNMBR( 'P', 'R', 'C', N, N, N, WORK( IR ), LDWRKR,
Line 633	Line 805
*	*
* Multiply Q in A by left singular vectors of R in	* Multiply Q in A by left singular vectors of R in
* WORK(IU), storing result in WORK(IR) and copying to A	* WORK(IU), storing result in WORK(IR) and copying to A
* (CWorkspace: need 2NN, prefer NN+MN)	* CWorkspace: need NN [U] + NN [R]
* (RWorkspace: 0)	* CWorkspace: prefer NN [U] + MN [R]
	* RWorkspace: need 0
*	*
DO 10 I = 1, M, LDWRKR	DO 10 I = 1, M, LDWRKR
CHUNK = MIN( M-I+1, LDWRKR )	CHUNK = MIN( M-I+1, LDWRKR )
Line 647	Line 820
*	*
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
*	*
* Path 3 (M much larger than N, JOBZ='S')	* Path 3 (M >> N, JOBZ='S')
* N left singular vectors to be computed in U and	* N left singular vectors to be computed in U and
* N right singular vectors to be computed in VT	* N right singular vectors to be computed in VT
*	*
Line 660	Line 833
NWORK = ITAU + N	NWORK = ITAU + N
*	*
* Compute A=Q*R	* Compute A=Q*R
* (CWorkspace: need NN+2N, prefer NN+N+NNB)	* CWorkspace: need N*N [R] + N [tau] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [R] + N [tau] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),	CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
Line 673	Line 847
$ LDWRKR )	$ LDWRKR )
*	*
* Generate Q in A	* Generate Q in A
* (CWorkspace: need 2N, prefer N+NNB)	* CWorkspace: need N*N [R] + N [tau] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [R] + N [tau] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGQR( M, N, N, A, LDA, WORK( ITAU ),	CALL ZUNGQR( M, N, N, A, LDA, WORK( ITAU ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
Line 684	Line 859
NWORK = ITAUP + N	NWORK = ITAUP + N
*	*
* Bidiagonalize R in WORK(IR)	* Bidiagonalize R in WORK(IR)
* (CWorkspace: need NN+3N, prefer NN+2N+2NNB)	* CWorkspace: need NN [R] + 2N [tauq, taup] + N [work]
* (RWorkspace: need N)	* CWorkspace: prefer NN [R] + 2N [tauq, taup] + 2NNB [work]
	* RWorkspace: need N [e]
*	*
CALL ZGEBRD( N, N, WORK( IR ), LDWRKR, S, RWORK( IE ),	CALL ZGEBRD( N, N, WORK( IR ), LDWRKR, S, RWORK( IE ),
$ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ),	$ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ),
Line 694	Line 870
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + BDSPAC
*	*
IRU = IE + N	IRU = IE + N
IRVT = IRU + N*N	IRVT = IRU + N*N
Line 706	Line 882
*	*
* Copy real matrix RWORK(IRU) to complex matrix U	* Copy real matrix RWORK(IRU) to complex matrix U
* Overwrite U by left singular vectors of R	* Overwrite U by left singular vectors of R
* (CWorkspace: need NN+3N, prefer NN+2N+N*NB)	* CWorkspace: need NN [R] + 2N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [R] + 2N [tauq, taup] + N*NB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, U, LDU )	CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, U, LDU )
CALL ZUNMBR( 'Q', 'L', 'N', N, N, N, WORK( IR ), LDWRKR,	CALL ZUNMBR( 'Q', 'L', 'N', N, N, N, WORK( IR ), LDWRKR,
Line 716	Line 893
*	*
* Copy real matrix RWORK(IRVT) to complex matrix VT	* Copy real matrix RWORK(IRVT) to complex matrix VT
* Overwrite VT by right singular vectors of R	* Overwrite VT by right singular vectors of R
* (CWorkspace: need NN+3N, prefer NN+2N+N*NB)	* CWorkspace: need NN [R] + 2N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [R] + 2N [tauq, taup] + N*NB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )	CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )
CALL ZUNMBR( 'P', 'R', 'C', N, N, N, WORK( IR ), LDWRKR,	CALL ZUNMBR( 'P', 'R', 'C', N, N, N, WORK( IR ), LDWRKR,
Line 726	Line 904
*	*
* Multiply Q in A by left singular vectors of R in	* Multiply Q in A by left singular vectors of R in
* WORK(IR), storing result in U	* WORK(IR), storing result in U
* (CWorkspace: need N*N)	* CWorkspace: need N*N [R]
* (RWorkspace: 0)	* RWorkspace: need 0
*	*
CALL ZLACPY( 'F', N, N, U, LDU, WORK( IR ), LDWRKR )	CALL ZLACPY( 'F', N, N, U, LDU, WORK( IR ), LDWRKR )
CALL ZGEMM( 'N', 'N', M, N, N, CONE, A, LDA, WORK( IR ),	CALL ZGEMM( 'N', 'N', M, N, N, CONE, A, LDA, WORK( IR ),
Line 735	Line 913
*	*
ELSE IF( WNTQA ) THEN	ELSE IF( WNTQA ) THEN
*	*
* Path 4 (M much larger than N, JOBZ='A')	* Path 4 (M >> N, JOBZ='A')
* M left singular vectors to be computed in U and	* M left singular vectors to be computed in U and
* N right singular vectors to be computed in VT	* N right singular vectors to be computed in VT
*	*
Line 748	Line 926
NWORK = ITAU + N	NWORK = ITAU + N
*	*
* Compute A=Q*R, copying result to U	* Compute A=Q*R, copying result to U
* (CWorkspace: need 2N, prefer N+NNB)	* CWorkspace: need N*N [U] + N [tau] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [U] + N [tau] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),	CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
CALL ZLACPY( 'L', M, N, A, LDA, U, LDU )	CALL ZLACPY( 'L', M, N, A, LDA, U, LDU )
*	*
* Generate Q in U	* Generate Q in U
* (CWorkspace: need N+M, prefer N+M*NB)	* CWorkspace: need N*N [U] + N [tau] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [U] + N [tau] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGQR( M, M, N, U, LDU, WORK( ITAU ),	CALL ZUNGQR( M, M, N, U, LDU, WORK( ITAU ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
Line 772	Line 952
NWORK = ITAUP + N	NWORK = ITAUP + N
*	*
* Bidiagonalize R in A	* Bidiagonalize R in A
* (CWorkspace: need 3N, prefer 2N+2NNB)	* CWorkspace: need NN [U] + 2N [tauq, taup] + N [work]
* (RWorkspace: need N)	* CWorkspace: prefer NN [U] + 2N [tauq, taup] + 2NNB [work]
	* RWorkspace: need N [e]
*	*
CALL ZGEBRD( N, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),	CALL ZGEBRD( N, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),
$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,	$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,
Line 785	Line 966
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + BDSPAC
*	*
CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ),	CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ),
$ N, RWORK( IRVT ), N, DUM, IDUM,	$ N, RWORK( IRVT ), N, DUM, IDUM,
Line 794	Line 975
*	*
* Copy real matrix RWORK(IRU) to complex matrix WORK(IU)	* Copy real matrix RWORK(IRU) to complex matrix WORK(IU)
* Overwrite WORK(IU) by left singular vectors of R	* Overwrite WORK(IU) by left singular vectors of R
* (CWorkspace: need NN+3N, prefer NN+2N+N*NB)	* CWorkspace: need NN [U] + 2N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [U] + 2N [tauq, taup] + N*NB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, WORK( IU ),	CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, WORK( IU ),
$ LDWRKU )	$ LDWRKU )
Line 805	Line 987
*	*
* Copy real matrix RWORK(IRVT) to complex matrix VT	* Copy real matrix RWORK(IRVT) to complex matrix VT
* Overwrite VT by right singular vectors of R	* Overwrite VT by right singular vectors of R
* (CWorkspace: need 3N, prefer 2N+N*NB)	* CWorkspace: need NN [U] + 2N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer NN [U] + 2N [tauq, taup] + N*NB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )	CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )
CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA,	CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA,
Line 815	Line 998
*	*
* Multiply Q in U by left singular vectors of R in	* Multiply Q in U by left singular vectors of R in
* WORK(IU), storing result in A	* WORK(IU), storing result in A
* (CWorkspace: need N*N)	* CWorkspace: need N*N [U]
* (RWorkspace: 0)	* RWorkspace: need 0
*	*
CALL ZGEMM( 'N', 'N', M, N, N, CONE, U, LDU, WORK( IU ),	CALL ZGEMM( 'N', 'N', M, N, N, CONE, U, LDU, WORK( IU ),
$ LDWRKU, CZERO, A, LDA )	$ LDWRKU, CZERO, A, LDA )
Line 831	Line 1014
*	*
* MNTHR2 <= M < MNTHR1	* MNTHR2 <= M < MNTHR1
*	*
* Path 5 (M much larger than N, but not as much as MNTHR1)	* Path 5 (M >> N, but not as much as MNTHR1)
* Reduce to bidiagonal form without QR decomposition, use	* Reduce to bidiagonal form without QR decomposition, use
* ZUNGBR and matrix multiplication to compute singular vectors	* ZUNGBR and matrix multiplication to compute singular vectors
*	*
Line 842	Line 1025
NWORK = ITAUP + N	NWORK = ITAUP + N
*	*
* Bidiagonalize A	* Bidiagonalize A
* (CWorkspace: need 2N+M, prefer 2N+(M+N)*NB)	* CWorkspace: need 2*N [tauq, taup] + M [work]
* (RWorkspace: need N)	* CWorkspace: prefer 2N [tauq, taup] + (M+N)NB [work]
	* RWorkspace: need N [e]
*	*
CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),	CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),
$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,	$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,
$ IERR )	$ IERR )
IF( WNTQN ) THEN	IF( WNTQN ) THEN
*	*
	* Path 5n (M >> N, JOBZ='N')
* Compute singular values only	* Compute singular values only
* (Cworkspace: 0)	* CWorkspace: need 0
* (Rworkspace: need BDSPAN)	* RWorkspace: need N [e] + BDSPAC
*	*
CALL DBDSDC( 'U', 'N', N, S, RWORK( IE ), DUM, 1, DUM, 1,	CALL DBDSDC( 'U', 'N', N, S, RWORK( IE ), DUM, 1,DUM,1,
$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )	$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )
ELSE IF( WNTQO ) THEN	ELSE IF( WNTQO ) THEN
IU = NWORK	IU = NWORK
Line 862	Line 1047
IRVT = IRU + N*N	IRVT = IRU + N*N
NRWORK = IRVT + N*N	NRWORK = IRVT + N*N
*	*
	* Path 5o (M >> N, JOBZ='O')
* Copy A to VT, generate P**H	* Copy A to VT, generate P**H
* (Cworkspace: need 2N, prefer N+NNB)	* CWorkspace: need 2*N [tauq, taup] + N [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'U', N, N, A, LDA, VT, LDVT )	CALL ZLACPY( 'U', N, N, A, LDA, VT, LDVT )
CALL ZUNGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ),	CALL ZUNGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
* Generate Q in A	* Generate Q in A
* (CWorkspace: need 2N, prefer N+NNB)	* CWorkspace: need 2*N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGBR( 'Q', M, N, N, A, LDA, WORK( ITAUQ ),	CALL ZUNGBR( 'Q', M, N, N, A, LDA, WORK( ITAUQ ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
IF( LWORK.GE.MN+3N ) THEN	IF( LWORK .GE. MN + 3N ) THEN
*	*
* WORK( IU ) is M by N	* WORK( IU ) is M by N
*	*
Line 886	Line 1074
*	*
* WORK(IU) is LDWRKU by N	* WORK(IU) is LDWRKU by N
*	*
LDWRKU = ( LWORK-3*N ) / N	LDWRKU = ( LWORK - 3*N ) / N
END IF	END IF
NWORK = IU + LDWRKU*N	NWORK = IU + LDWRKU*N
*	*
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + BDSPAC
*	*
CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ),	CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ),
$ N, RWORK( IRVT ), N, DUM, IDUM,	$ N, RWORK( IRVT ), N, DUM, IDUM,
Line 902	Line 1090
*	*
* Multiply real matrix RWORK(IRVT) by P**H in VT,	* Multiply real matrix RWORK(IRVT) by P**H in VT,
* storing the result in WORK(IU), copying to VT	* storing the result in WORK(IU), copying to VT
* (Cworkspace: need 0)	* CWorkspace: need 2N [tauq, taup] + NN [U]
* (Rworkspace: need 3NN)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + 2NN [rwork]
*	*
CALL ZLARCM( N, N, RWORK( IRVT ), N, VT, LDVT,	CALL ZLARCM( N, N, RWORK( IRVT ), N, VT, LDVT,
$ WORK( IU ), LDWRKU, RWORK( NRWORK ) )	$ WORK( IU ), LDWRKU, RWORK( NRWORK ) )
Line 911	Line 1099
*	*
* Multiply Q in A by real matrix RWORK(IRU), storing the	* Multiply Q in A by real matrix RWORK(IRU), storing the
* result in WORK(IU), copying to A	* result in WORK(IU), copying to A
* (CWorkspace: need NN, prefer MN)	* CWorkspace: need 2N [tauq, taup] + NN [U]
* (Rworkspace: need 3NN, prefer NN+2M*N)	* CWorkspace: prefer 2N [tauq, taup] + MN [U]
	* RWorkspace: need N [e] + NN [RU] + 2N*N [rwork]
	* RWorkspace: prefer N [e] + NN [RU] + 2MN [rwork] < N + 5NN since M < 2N here
*	*
NRWORK = IRVT	NRWORK = IRVT
DO 20 I = 1, M, LDWRKU	DO 20 I = 1, M, LDWRKU
Line 925	Line 1115
*	*
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
*	*
	* Path 5s (M >> N, JOBZ='S')
* Copy A to VT, generate P**H	* Copy A to VT, generate P**H
* (Cworkspace: need 2N, prefer N+NNB)	* CWorkspace: need 2*N [tauq, taup] + N [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'U', N, N, A, LDA, VT, LDVT )	CALL ZLACPY( 'U', N, N, A, LDA, VT, LDVT )
CALL ZUNGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ),	CALL ZUNGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
* Copy A to U, generate Q	* Copy A to U, generate Q
* (Cworkspace: need 2N, prefer N+NNB)	* CWorkspace: need 2*N [tauq, taup] + N [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'L', M, N, A, LDA, U, LDU )	CALL ZLACPY( 'L', M, N, A, LDA, U, LDU )
CALL ZUNGBR( 'Q', M, N, N, U, LDU, WORK( ITAUQ ),	CALL ZUNGBR( 'Q', M, N, N, U, LDU, WORK( ITAUQ ),
Line 944	Line 1137
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + BDSPAC
*	*
IRU = NRWORK	IRU = NRWORK
IRVT = IRU + N*N	IRVT = IRU + N*N
Line 956	Line 1149
*	*
* Multiply real matrix RWORK(IRVT) by P**H in VT,	* Multiply real matrix RWORK(IRVT) by P**H in VT,
* storing the result in A, copying to VT	* storing the result in A, copying to VT
* (Cworkspace: need 0)	* CWorkspace: need 0
* (Rworkspace: need 3NN)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + 2NN [rwork]
*	*
CALL ZLARCM( N, N, RWORK( IRVT ), N, VT, LDVT, A, LDA,	CALL ZLARCM( N, N, RWORK( IRVT ), N, VT, LDVT, A, LDA,
$ RWORK( NRWORK ) )	$ RWORK( NRWORK ) )
Line 965	Line 1158
*	*
* Multiply Q in U by real matrix RWORK(IRU), storing the	* Multiply Q in U by real matrix RWORK(IRU), storing the
* result in A, copying to U	* result in A, copying to U
* (CWorkspace: need 0)	* CWorkspace: need 0
* (Rworkspace: need NN+2M*N)	* RWorkspace: need N [e] + NN [RU] + 2MN [rwork] < N + 5NN since M < 2N here
*	*
NRWORK = IRVT	NRWORK = IRVT
CALL ZLACRM( M, N, U, LDU, RWORK( IRU ), N, A, LDA,	CALL ZLACRM( M, N, U, LDU, RWORK( IRU ), N, A, LDA,
Line 974	Line 1167
CALL ZLACPY( 'F', M, N, A, LDA, U, LDU )	CALL ZLACPY( 'F', M, N, A, LDA, U, LDU )
ELSE	ELSE
*	*
	* Path 5a (M >> N, JOBZ='A')
* Copy A to VT, generate P**H	* Copy A to VT, generate P**H
* (Cworkspace: need 2N, prefer N+NNB)	* CWorkspace: need 2*N [tauq, taup] + N [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'U', N, N, A, LDA, VT, LDVT )	CALL ZLACPY( 'U', N, N, A, LDA, VT, LDVT )
CALL ZUNGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ),	CALL ZUNGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
* Copy A to U, generate Q	* Copy A to U, generate Q
* (Cworkspace: need 2N, prefer N+NNB)	* CWorkspace: need 2*N [tauq, taup] + M [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'L', M, N, A, LDA, U, LDU )	CALL ZLACPY( 'L', M, N, A, LDA, U, LDU )
CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ),	CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ),
Line 993	Line 1189
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + BDSPAC
*	*
IRU = NRWORK	IRU = NRWORK
IRVT = IRU + N*N	IRVT = IRU + N*N
Line 1005	Line 1201
*	*
* Multiply real matrix RWORK(IRVT) by P**H in VT,	* Multiply real matrix RWORK(IRVT) by P**H in VT,
* storing the result in A, copying to VT	* storing the result in A, copying to VT
* (Cworkspace: need 0)	* CWorkspace: need 0
* (Rworkspace: need 3NN)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + 2NN [rwork]
*	*
CALL ZLARCM( N, N, RWORK( IRVT ), N, VT, LDVT, A, LDA,	CALL ZLARCM( N, N, RWORK( IRVT ), N, VT, LDVT, A, LDA,
$ RWORK( NRWORK ) )	$ RWORK( NRWORK ) )
Line 1014	Line 1210
*	*
* Multiply Q in U by real matrix RWORK(IRU), storing the	* Multiply Q in U by real matrix RWORK(IRU), storing the
* result in A, copying to U	* result in A, copying to U
* (CWorkspace: 0)	* CWorkspace: need 0
* (Rworkspace: need 3NN)	* RWorkspace: need N [e] + NN [RU] + 2MN [rwork] < N + 5NN since M < 2N here
*	*
NRWORK = IRVT	NRWORK = IRVT
CALL ZLACRM( M, N, U, LDU, RWORK( IRU ), N, A, LDA,	CALL ZLACRM( M, N, U, LDU, RWORK( IRU ), N, A, LDA,
Line 1027	Line 1223
*	*
* M .LT. MNTHR2	* M .LT. MNTHR2
*	*
* Path 6 (M at least N, but not much larger)	* Path 6 (M >= N, but not much larger)
* Reduce to bidiagonal form without QR decomposition	* Reduce to bidiagonal form without QR decomposition
* Use ZUNMBR to compute singular vectors	* Use ZUNMBR to compute singular vectors
*	*
Line 1038	Line 1234
NWORK = ITAUP + N	NWORK = ITAUP + N
*	*
* Bidiagonalize A	* Bidiagonalize A
* (CWorkspace: need 2N+M, prefer 2N+(M+N)*NB)	* CWorkspace: need 2*N [tauq, taup] + M [work]
* (RWorkspace: need N)	* CWorkspace: prefer 2N [tauq, taup] + (M+N)NB [work]
	* RWorkspace: need N [e]
*	*
CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),	CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),
$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,	$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,
$ IERR )	$ IERR )
IF( WNTQN ) THEN	IF( WNTQN ) THEN
*	*
	* Path 6n (M >= N, JOBZ='N')
* Compute singular values only	* Compute singular values only
* (Cworkspace: 0)	* CWorkspace: need 0
* (Rworkspace: need BDSPAN)	* RWorkspace: need N [e] + BDSPAC
*	*
CALL DBDSDC( 'U', 'N', N, S, RWORK( IE ), DUM, 1, DUM, 1,	CALL DBDSDC( 'U', 'N', N, S, RWORK( IE ), DUM,1,DUM,1,
$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )	$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )
ELSE IF( WNTQO ) THEN	ELSE IF( WNTQO ) THEN
IU = NWORK	IU = NWORK
IRU = NRWORK	IRU = NRWORK
IRVT = IRU + N*N	IRVT = IRU + N*N
NRWORK = IRVT + N*N	NRWORK = IRVT + N*N
IF( LWORK.GE.MN+3N ) THEN	IF( LWORK .GE. MN + 3N ) THEN
*	*
* WORK( IU ) is M by N	* WORK( IU ) is M by N
*	*
Line 1066	Line 1264
*	*
* WORK( IU ) is LDWRKU by N	* WORK( IU ) is LDWRKU by N
*	*
LDWRKU = ( LWORK-3*N ) / N	LDWRKU = ( LWORK - 3*N ) / N
END IF	END IF
NWORK = IU + LDWRKU*N	NWORK = IU + LDWRKU*N
*	*
	* Path 6o (M >= N, JOBZ='O')
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + BDSPAC
*	*
CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ),	CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ),
$ N, RWORK( IRVT ), N, DUM, IDUM,	$ N, RWORK( IRVT ), N, DUM, IDUM,
Line 1082	Line 1281
*	*
* Copy real matrix RWORK(IRVT) to complex matrix VT	* Copy real matrix RWORK(IRVT) to complex matrix VT
* Overwrite VT by right singular vectors of A	* Overwrite VT by right singular vectors of A
* (Cworkspace: need 2N, prefer N+NNB)	* CWorkspace: need 2N [tauq, taup] + NN [U] + N [work]
* (Rworkspace: need 0)	* CWorkspace: prefer 2N [tauq, taup] + NN [U] + N*NB [work]
	* RWorkspace: need N [e] + NN [RU] + NN [RVT]
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )	CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )
CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA,	CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA,
$ WORK( ITAUP ), VT, LDVT, WORK( NWORK ),	$ WORK( ITAUP ), VT, LDVT, WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
*	*
IF( LWORK.GE.MN+3N ) THEN	IF( LWORK .GE. MN + 3N ) THEN
*	*
* Copy real matrix RWORK(IRU) to complex matrix WORK(IU)	* Path 6o-fast
* Overwrite WORK(IU) by left singular vectors of A, copying	* Copy real matrix RWORK(IRU) to complex matrix WORK(IU)
* to A	* Overwrite WORK(IU) by left singular vectors of A, copying
* (Cworkspace: need MN+2N, prefer MN+N+NNB)	* to A
* (Rworkspace: need 0)	* CWorkspace: need 2N [tauq, taup] + MN [U] + N [work]
	* CWorkspace: prefer 2N [tauq, taup] + MN [U] + N*NB [work]
	* RWorkspace: need N [e] + N*N [RU]
*	*
CALL ZLASET( 'F', M, N, CZERO, CZERO, WORK( IU ),	CALL ZLASET( 'F', M, N, CZERO, CZERO, WORK( IU ),
$ LDWRKU )	$ LDWRKU )
Line 1108	Line 1310
CALL ZLACPY( 'F', M, N, WORK( IU ), LDWRKU, A, LDA )	CALL ZLACPY( 'F', M, N, WORK( IU ), LDWRKU, A, LDA )
ELSE	ELSE
*	*
	* Path 6o-slow
* Generate Q in A	* Generate Q in A
* (Cworkspace: need 2N, prefer N+NNB)	* CWorkspace: need 2N [tauq, taup] + NN [U] + N [work]
* (Rworkspace: need 0)	* CWorkspace: prefer 2N [tauq, taup] + NN [U] + N*NB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGBR( 'Q', M, N, N, A, LDA, WORK( ITAUQ ),	CALL ZUNGBR( 'Q', M, N, N, A, LDA, WORK( ITAUQ ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
* Multiply Q in A by real matrix RWORK(IRU), storing the	* Multiply Q in A by real matrix RWORK(IRU), storing the
* result in WORK(IU), copying to A	* result in WORK(IU), copying to A
* (CWorkspace: need NN, prefer MN)	* CWorkspace: need 2N [tauq, taup] + NN [U]
* (Rworkspace: need 3NN, prefer NN+2M*N)	* CWorkspace: prefer 2N [tauq, taup] + MN [U]
	* RWorkspace: need N [e] + NN [RU] + 2N*N [rwork]
	* RWorkspace: prefer N [e] + NN [RU] + 2MN [rwork] < N + 5NN since M < 2N here
*	*
NRWORK = IRVT	NRWORK = IRVT
DO 30 I = 1, M, LDWRKU	DO 30 I = 1, M, LDWRKU
Line 1133	Line 1339
*	*
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
*	*
	* Path 6s (M >= N, JOBZ='S')
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + BDSPAC
*	*
IRU = NRWORK	IRU = NRWORK
IRVT = IRU + N*N	IRVT = IRU + N*N
Line 1148	Line 1355
*	*
* Copy real matrix RWORK(IRU) to complex matrix U	* Copy real matrix RWORK(IRU) to complex matrix U
* Overwrite U by left singular vectors of A	* Overwrite U by left singular vectors of A
* (CWorkspace: need 3N, prefer 2N+N*NB)	* CWorkspace: need 2*N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + NNB [work]
	* RWorkspace: need N [e] + NN [RU] + NN [RVT]
*	*
CALL ZLASET( 'F', M, N, CZERO, CZERO, U, LDU )	CALL ZLASET( 'F', M, N, CZERO, CZERO, U, LDU )
CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, U, LDU )	CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, U, LDU )
Line 1159	Line 1367
*	*
* Copy real matrix RWORK(IRVT) to complex matrix VT	* Copy real matrix RWORK(IRVT) to complex matrix VT
* Overwrite VT by right singular vectors of A	* Overwrite VT by right singular vectors of A
* (CWorkspace: need 3N, prefer 2N+N*NB)	* CWorkspace: need 2*N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + NNB [work]
	* RWorkspace: need N [e] + NN [RU] + NN [RVT]
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )	CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )
CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA,	CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA,
Line 1168	Line 1377
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
ELSE	ELSE
*	*
	* Path 6a (M >= N, JOBZ='A')
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need N [e] + NN [RU] + NN [RVT] + BDSPAC
*	*
IRU = NRWORK	IRU = NRWORK
IRVT = IRU + N*N	IRVT = IRU + N*N
Line 1191	Line 1401
*	*
* Copy real matrix RWORK(IRU) to complex matrix U	* Copy real matrix RWORK(IRU) to complex matrix U
* Overwrite U by left singular vectors of A	* Overwrite U by left singular vectors of A
* (CWorkspace: need 2N+M, prefer 2N+M*NB)	* CWorkspace: need 2*N [tauq, taup] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + MNB [work]
	* RWorkspace: need N [e] + NN [RU] + NN [RVT]
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, U, LDU )	CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, U, LDU )
CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA,	CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA,
Line 1201	Line 1412
*	*
* Copy real matrix RWORK(IRVT) to complex matrix VT	* Copy real matrix RWORK(IRVT) to complex matrix VT
* Overwrite VT by right singular vectors of A	* Overwrite VT by right singular vectors of A
* (CWorkspace: need 3N, prefer 2N+N*NB)	* CWorkspace: need 2*N [tauq, taup] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer 2N [tauq, taup] + NNB [work]
	* RWorkspace: need N [e] + NN [RU] + NN [RVT]
*	*
CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )	CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT )
CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA,	CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA,
Line 1222	Line 1434
*	*
IF( WNTQN ) THEN	IF( WNTQN ) THEN
*	*
* Path 1t (N much larger than M, JOBZ='N')	* Path 1t (N >> M, JOBZ='N')
* No singular vectors to be computed	* No singular vectors to be computed
*	*
ITAU = 1	ITAU = 1
NWORK = ITAU + M	NWORK = ITAU + M
*	*
* Compute A=L*Q	* Compute A=L*Q
* (CWorkspace: need 2M, prefer M+MNB)	* CWorkspace: need M [tau] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer M [tau] + M*NB [work]
	* RWorkspace: need 0
*	*
CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),	CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
Line 1245	Line 1458
NWORK = ITAUP + M	NWORK = ITAUP + M
*	*
* Bidiagonalize L in A	* Bidiagonalize L in A
* (CWorkspace: need 3M, prefer 2M+2MNB)	* CWorkspace: need 2*M [tauq, taup] + M [work]
* (RWorkspace: need M)	* CWorkspace: prefer 2M [tauq, taup] + 2M*NB [work]
	* RWorkspace: need M [e]
*	*
CALL ZGEBRD( M, M, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),	CALL ZGEBRD( M, M, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),
$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,	$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,
Line 1254	Line 1468
NRWORK = IE + M	NRWORK = IE + M
*	*
* Perform bidiagonal SVD, compute singular values only	* Perform bidiagonal SVD, compute singular values only
* (CWorkspace: 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAN)	* RWorkspace: need M [e] + BDSPAC
*	*
CALL DBDSDC( 'U', 'N', M, S, RWORK( IE ), DUM, 1, DUM, 1,	CALL DBDSDC( 'U', 'N', M, S, RWORK( IE ), DUM,1,DUM,1,
$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )	$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )
*	*
ELSE IF( WNTQO ) THEN	ELSE IF( WNTQO ) THEN
*	*
* Path 2t (N much larger than M, JOBZ='O')	* Path 2t (N >> M, JOBZ='O')
* M right singular vectors to be overwritten on A and	* M right singular vectors to be overwritten on A and
* M left singular vectors to be computed in U	* M left singular vectors to be computed in U
*	*
Line 1272	Line 1486
* WORK(IVT) is M by M	* WORK(IVT) is M by M
*	*
IL = IVT + LDWKVT*M	IL = IVT + LDWKVT*M
IF( LWORK.GE.MN+MM+3*M ) THEN	IF( LWORK .GE. MN + MM + 3*M ) THEN
*	*
* WORK(IL) M by N	* WORK(IL) M by N
*	*
Line 1283	Line 1497
* WORK(IL) is M by CHUNK	* WORK(IL) is M by CHUNK
*	*
LDWRKL = M	LDWRKL = M
CHUNK = ( LWORK-MM-3M ) / M	CHUNK = ( LWORK - MM - 3M ) / M
END IF	END IF
ITAU = IL + LDWRKL*CHUNK	ITAU = IL + LDWRKL*CHUNK
NWORK = ITAU + M	NWORK = ITAU + M
*	*
* Compute A=L*Q	* Compute A=L*Q
* (CWorkspace: need 2M, prefer M+MNB)	* CWorkspace: need MM [VT] + MM [L] + M [tau] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [VT] + MM [L] + M [tau] + M*NB [work]
	* RWorkspace: need 0
*	*
CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),	CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
Line 1302	Line 1517
$ WORK( IL+LDWRKL ), LDWRKL )	$ WORK( IL+LDWRKL ), LDWRKL )
*	*
* Generate Q in A	* Generate Q in A
* (CWorkspace: need MM+2M, prefer MM+M+MNB)	* CWorkspace: need MM [VT] + MM [L] + M [tau] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [VT] + MM [L] + M [tau] + M*NB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGLQ( M, N, M, A, LDA, WORK( ITAU ),	CALL ZUNGLQ( M, N, M, A, LDA, WORK( ITAU ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
Line 1313	Line 1529
NWORK = ITAUP + M	NWORK = ITAUP + M
*	*
* Bidiagonalize L in WORK(IL)	* Bidiagonalize L in WORK(IL)
* (CWorkspace: need MM+3M, prefer MM+2M+2MNB)	* CWorkspace: need MM [VT] + MM [L] + 2*M [tauq, taup] + M [work]
* (RWorkspace: need M)	* CWorkspace: prefer MM [VT] + MM [L] + 2M [tauq, taup] + 2M*NB [work]
	* RWorkspace: need M [e]
*	*
CALL ZGEBRD( M, M, WORK( IL ), LDWRKL, S, RWORK( IE ),	CALL ZGEBRD( M, M, WORK( IL ), LDWRKL, S, RWORK( IE ),
$ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ),	$ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ),
Line 1323	Line 1540
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need M [e] + MM [RU] + MM [RVT] + BDSPAC
*	*
IRU = IE + M	IRU = IE + M
IRVT = IRU + M*M	IRVT = IRU + M*M
Line 1335	Line 1552
*	*
* Copy real matrix RWORK(IRU) to complex matrix WORK(IU)	* Copy real matrix RWORK(IRU) to complex matrix WORK(IU)
* Overwrite WORK(IU) by the left singular vectors of L	* Overwrite WORK(IU) by the left singular vectors of L
* (CWorkspace: need NN+3N, prefer MN+2N+N*NB)	* CWorkspace: need MM [VT] + MM [L] + 2*M [tauq, taup] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [VT] + MM [L] + 2M [tauq, taup] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )	CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )
CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, WORK( IL ), LDWRKL,	CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, WORK( IL ), LDWRKL,
Line 1345	Line 1563
*	*
* Copy real matrix RWORK(IRVT) to complex matrix WORK(IVT)	* Copy real matrix RWORK(IRVT) to complex matrix WORK(IVT)
* Overwrite WORK(IVT) by the right singular vectors of L	* Overwrite WORK(IVT) by the right singular vectors of L
* (CWorkspace: need NN+3N, prefer MN+2N+N*NB)	* CWorkspace: need MM [VT] + MM [L] + 2*M [tauq, taup] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [VT] + MM [L] + 2M [tauq, taup] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, WORK( IVT ),	CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, WORK( IVT ),
$ LDWKVT )	$ LDWKVT )
Line 1356	Line 1575
*	*
* Multiply right singular vectors of L in WORK(IL) by Q	* Multiply right singular vectors of L in WORK(IL) by Q
* in A, storing result in WORK(IL) and copying to A	* in A, storing result in WORK(IL) and copying to A
* (CWorkspace: need 2MM, prefer MM+MN))	* CWorkspace: need MM [VT] + MM [L]
* (RWorkspace: 0)	* CWorkspace: prefer MM [VT] + MN [L]
	* RWorkspace: need 0
*	*
DO 40 I = 1, N, CHUNK	DO 40 I = 1, N, CHUNK
BLK = MIN( N-I+1, CHUNK )	BLK = MIN( N-I+1, CHUNK )
Line 1370	Line 1590
*	*
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
*	*
* Path 3t (N much larger than M, JOBZ='S')	* Path 3t (N >> M, JOBZ='S')
* M right singular vectors to be computed in VT and	* M right singular vectors to be computed in VT and
* M left singular vectors to be computed in U	* M left singular vectors to be computed in U
*	*
IL = 1	IL = 1
*	*
Line 1383	Line 1603
NWORK = ITAU + M	NWORK = ITAU + M
*	*
* Compute A=L*Q	* Compute A=L*Q
* (CWorkspace: need 2M, prefer M+MNB)	* CWorkspace: need M*M [L] + M [tau] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [L] + M [tau] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),	CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
Line 1396	Line 1617
$ WORK( IL+LDWRKL ), LDWRKL )	$ WORK( IL+LDWRKL ), LDWRKL )
*	*
* Generate Q in A	* Generate Q in A
* (CWorkspace: need MM+2M, prefer MM+M+MNB)	* CWorkspace: need M*M [L] + M [tau] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [L] + M [tau] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGLQ( M, N, M, A, LDA, WORK( ITAU ),	CALL ZUNGLQ( M, N, M, A, LDA, WORK( ITAU ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
Line 1407	Line 1629
NWORK = ITAUP + M	NWORK = ITAUP + M
*	*
* Bidiagonalize L in WORK(IL)	* Bidiagonalize L in WORK(IL)
* (CWorkspace: need MM+3M, prefer MM+2M+2MNB)	* CWorkspace: need MM [L] + 2M [tauq, taup] + M [work]
* (RWorkspace: need M)	* CWorkspace: prefer MM [L] + 2M [tauq, taup] + 2MNB [work]
	* RWorkspace: need M [e]
*	*
CALL ZGEBRD( M, M, WORK( IL ), LDWRKL, S, RWORK( IE ),	CALL ZGEBRD( M, M, WORK( IL ), LDWRKL, S, RWORK( IE ),
$ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ),	$ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ),
Line 1417	Line 1640
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need M [e] + MM [RU] + MM [RVT] + BDSPAC
*	*
IRU = IE + M	IRU = IE + M
IRVT = IRU + M*M	IRVT = IRU + M*M
Line 1429	Line 1652
*	*
* Copy real matrix RWORK(IRU) to complex matrix U	* Copy real matrix RWORK(IRU) to complex matrix U
* Overwrite U by left singular vectors of L	* Overwrite U by left singular vectors of L
* (CWorkspace: need MM+3M, prefer MM+2M+M*NB)	* CWorkspace: need MM [L] + 2M [tauq, taup] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [L] + 2M [tauq, taup] + M*NB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )	CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )
CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, WORK( IL ), LDWRKL,	CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, WORK( IL ), LDWRKL,
Line 1439	Line 1663
*	*
* Copy real matrix RWORK(IRVT) to complex matrix VT	* Copy real matrix RWORK(IRVT) to complex matrix VT
* Overwrite VT by left singular vectors of L	* Overwrite VT by left singular vectors of L
* (CWorkspace: need MM+3M, prefer MM+2M+M*NB)	* CWorkspace: need MM [L] + 2M [tauq, taup] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [L] + 2M [tauq, taup] + M*NB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, VT, LDVT )	CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, VT, LDVT )
CALL ZUNMBR( 'P', 'R', 'C', M, M, M, WORK( IL ), LDWRKL,	CALL ZUNMBR( 'P', 'R', 'C', M, M, M, WORK( IL ), LDWRKL,
Line 1449	Line 1674
*	*
* Copy VT to WORK(IL), multiply right singular vectors of L	* Copy VT to WORK(IL), multiply right singular vectors of L
* in WORK(IL) by Q in A, storing result in VT	* in WORK(IL) by Q in A, storing result in VT
* (CWorkspace: need M*M)	* CWorkspace: need M*M [L]
* (RWorkspace: 0)	* RWorkspace: need 0
*	*
CALL ZLACPY( 'F', M, M, VT, LDVT, WORK( IL ), LDWRKL )	CALL ZLACPY( 'F', M, M, VT, LDVT, WORK( IL ), LDWRKL )
CALL ZGEMM( 'N', 'N', M, N, M, CONE, WORK( IL ), LDWRKL,	CALL ZGEMM( 'N', 'N', M, N, M, CONE, WORK( IL ), LDWRKL,
Line 1458	Line 1683
*	*
ELSE IF( WNTQA ) THEN	ELSE IF( WNTQA ) THEN
*	*
* Path 9t (N much larger than M, JOBZ='A')	* Path 4t (N >> M, JOBZ='A')
* N right singular vectors to be computed in VT and	* N right singular vectors to be computed in VT and
* M left singular vectors to be computed in U	* M left singular vectors to be computed in U
*	*
Line 1471	Line 1696
NWORK = ITAU + M	NWORK = ITAU + M
*	*
* Compute A=L*Q, copying result to VT	* Compute A=L*Q, copying result to VT
* (CWorkspace: need 2M, prefer M+MNB)	* CWorkspace: need M*M [VT] + M [tau] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [VT] + M [tau] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),	CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
CALL ZLACPY( 'U', M, N, A, LDA, VT, LDVT )	CALL ZLACPY( 'U', M, N, A, LDA, VT, LDVT )
*	*
* Generate Q in VT	* Generate Q in VT
* (CWorkspace: need M+N, prefer M+N*NB)	* CWorkspace: need M*M [VT] + M [tau] + N [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [VT] + M [tau] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGLQ( N, N, M, VT, LDVT, WORK( ITAU ),	CALL ZUNGLQ( N, N, M, VT, LDVT, WORK( ITAU ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
Line 1495	Line 1722
NWORK = ITAUP + M	NWORK = ITAUP + M
*	*
* Bidiagonalize L in A	* Bidiagonalize L in A
* (CWorkspace: need MM+3M, prefer MM+2M+2MNB)	* CWorkspace: need MM [VT] + 2M [tauq, taup] + M [work]
* (RWorkspace: need M)	* CWorkspace: prefer MM [VT] + 2M [tauq, taup] + 2MNB [work]
	* RWorkspace: need M [e]
*	*
CALL ZGEBRD( M, M, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),	CALL ZGEBRD( M, M, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),
$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,	$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,
Line 1505	Line 1733
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need M [e] + MM [RU] + MM [RVT] + BDSPAC
*	*
IRU = IE + M	IRU = IE + M
IRVT = IRU + M*M	IRVT = IRU + M*M
Line 1517	Line 1745
*	*
* Copy real matrix RWORK(IRU) to complex matrix U	* Copy real matrix RWORK(IRU) to complex matrix U
* Overwrite U by left singular vectors of L	* Overwrite U by left singular vectors of L
* (CWorkspace: need 3M, prefer 2M+M*NB)	* CWorkspace: need MM [VT] + 2M [tauq, taup] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [VT] + 2M [tauq, taup] + M*NB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )	CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )
CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, A, LDA,	CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, A, LDA,
Line 1527	Line 1756
*	*
* Copy real matrix RWORK(IRVT) to complex matrix WORK(IVT)	* Copy real matrix RWORK(IRVT) to complex matrix WORK(IVT)
* Overwrite WORK(IVT) by right singular vectors of L	* Overwrite WORK(IVT) by right singular vectors of L
* (CWorkspace: need MM+3M, prefer MM+2M+M*NB)	* CWorkspace: need MM [VT] + 2M [tauq, taup] + M [work]
* (RWorkspace: 0)	* CWorkspace: prefer MM [VT] + 2M [tauq, taup] + M*NB [work]
	* RWorkspace: need 0
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, WORK( IVT ),	CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, WORK( IVT ),
$ LDWKVT )	$ LDWKVT )
Line 1538	Line 1768
*	*
* Multiply right singular vectors of L in WORK(IVT) by	* Multiply right singular vectors of L in WORK(IVT) by
* Q in VT, storing result in A	* Q in VT, storing result in A
* (CWorkspace: need M*M)	* CWorkspace: need M*M [VT]
* (RWorkspace: 0)	* RWorkspace: need 0
*	*
CALL ZGEMM( 'N', 'N', M, N, M, CONE, WORK( IVT ), LDWKVT,	CALL ZGEMM( 'N', 'N', M, N, M, CONE, WORK( IVT ), LDWKVT,
$ VT, LDVT, CZERO, A, LDA )	$ VT, LDVT, CZERO, A, LDA )
Line 1554	Line 1784
*	*
* MNTHR2 <= N < MNTHR1	* MNTHR2 <= N < MNTHR1
*	*
* Path 5t (N much larger than M, but not as much as MNTHR1)	* Path 5t (N >> M, but not as much as MNTHR1)
* Reduce to bidiagonal form without QR decomposition, use	* Reduce to bidiagonal form without QR decomposition, use
* ZUNGBR and matrix multiplication to compute singular vectors	* ZUNGBR and matrix multiplication to compute singular vectors
*	*
*
IE = 1	IE = 1
NRWORK = IE + M	NRWORK = IE + M
ITAUQ = 1	ITAUQ = 1
Line 1566	Line 1795
NWORK = ITAUP + M	NWORK = ITAUP + M
*	*
* Bidiagonalize A	* Bidiagonalize A
* (CWorkspace: need 2M+N, prefer 2M+(M+N)*NB)	* CWorkspace: need 2*M [tauq, taup] + N [work]
* (RWorkspace: M)	* CWorkspace: prefer 2M [tauq, taup] + (M+N)NB [work]
	* RWorkspace: need M [e]
*	*
CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),	CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),
$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,	$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,
Line 1575	Line 1805
*	*
IF( WNTQN ) THEN	IF( WNTQN ) THEN
*	*
	* Path 5tn (N >> M, JOBZ='N')
* Compute singular values only	* Compute singular values only
* (Cworkspace: 0)	* CWorkspace: need 0
* (Rworkspace: need BDSPAN)	* RWorkspace: need M [e] + BDSPAC
*	*
CALL DBDSDC( 'L', 'N', M, S, RWORK( IE ), DUM, 1, DUM, 1,	CALL DBDSDC( 'L', 'N', M, S, RWORK( IE ), DUM,1,DUM,1,
$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )	$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )
ELSE IF( WNTQO ) THEN	ELSE IF( WNTQO ) THEN
IRVT = NRWORK	IRVT = NRWORK
Line 1587	Line 1818
NRWORK = IRU + M*M	NRWORK = IRU + M*M
IVT = NWORK	IVT = NWORK
*	*
	* Path 5to (N >> M, JOBZ='O')
* Copy A to U, generate Q	* Copy A to U, generate Q
* (Cworkspace: need 2M, prefer M+MNB)	* CWorkspace: need 2*M [tauq, taup] + M [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2M [tauq, taup] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'L', M, M, A, LDA, U, LDU )	CALL ZLACPY( 'L', M, M, A, LDA, U, LDU )
CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ),	CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
* Generate P**H in A	* Generate P**H in A
* (Cworkspace: need 2M, prefer M+MNB)	* CWorkspace: need 2*M [tauq, taup] + M [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2M [tauq, taup] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGBR( 'P', M, N, M, A, LDA, WORK( ITAUP ),	CALL ZUNGBR( 'P', M, N, M, A, LDA, WORK( ITAUP ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
LDWKVT = M	LDWKVT = M
IF( LWORK.GE.MN+3M ) THEN	IF( LWORK .GE. MN + 3M ) THEN
*	*
* WORK( IVT ) is M by N	* WORK( IVT ) is M by N
*	*
Line 1613	Line 1847
*	*
* WORK( IVT ) is M by CHUNK	* WORK( IVT ) is M by CHUNK
*	*
CHUNK = ( LWORK-3*M ) / M	CHUNK = ( LWORK - 3*M ) / M
NWORK = IVT + LDWKVT*CHUNK	NWORK = IVT + LDWKVT*CHUNK
END IF	END IF
*	*
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need M [e] + MM [RVT] + MM [RU] + BDSPAC
*	*
CALL DBDSDC( 'L', 'I', M, S, RWORK( IE ), RWORK( IRU ),	CALL DBDSDC( 'L', 'I', M, S, RWORK( IE ), RWORK( IRU ),
$ M, RWORK( IRVT ), M, DUM, IDUM,	$ M, RWORK( IRVT ), M, DUM, IDUM,
Line 1629	Line 1863
*	*
* Multiply Q in U by real matrix RWORK(IRVT)	* Multiply Q in U by real matrix RWORK(IRVT)
* storing the result in WORK(IVT), copying to U	* storing the result in WORK(IVT), copying to U
* (Cworkspace: need 0)	* CWorkspace: need 2M [tauq, taup] + MM [VT]
* (Rworkspace: need 2MM)	* RWorkspace: need M [e] + MM [RVT] + MM [RU] + 2MM [rwork]
*	*
CALL ZLACRM( M, M, U, LDU, RWORK( IRU ), M, WORK( IVT ),	CALL ZLACRM( M, M, U, LDU, RWORK( IRU ), M, WORK( IVT ),
$ LDWKVT, RWORK( NRWORK ) )	$ LDWKVT, RWORK( NRWORK ) )
Line 1638	Line 1872
*	*
* Multiply RWORK(IRVT) by P**H in A, storing the	* Multiply RWORK(IRVT) by P**H in A, storing the
* result in WORK(IVT), copying to A	* result in WORK(IVT), copying to A
* (CWorkspace: need MM, prefer MN)	* CWorkspace: need 2M [tauq, taup] + MM [VT]
* (Rworkspace: need 2MM, prefer 2MN)	* CWorkspace: prefer 2M [tauq, taup] + MN [VT]
	* RWorkspace: need M [e] + MM [RVT] + 2M*M [rwork]
	* RWorkspace: prefer M [e] + MM [RVT] + 2MN [rwork] < M + 5MM since N < 2M here
*	*
NRWORK = IRU	NRWORK = IRU
DO 50 I = 1, N, CHUNK	DO 50 I = 1, N, CHUNK
Line 1651	Line 1887
50 CONTINUE	50 CONTINUE
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
*	*
	* Path 5ts (N >> M, JOBZ='S')
* Copy A to U, generate Q	* Copy A to U, generate Q
* (Cworkspace: need 2M, prefer M+MNB)	* CWorkspace: need 2*M [tauq, taup] + M [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2M [tauq, taup] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'L', M, M, A, LDA, U, LDU )	CALL ZLACPY( 'L', M, M, A, LDA, U, LDU )
CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ),	CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
* Copy A to VT, generate P**H	* Copy A to VT, generate P**H
* (Cworkspace: need 2M, prefer M+MNB)	* CWorkspace: need 2*M [tauq, taup] + M [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2M [tauq, taup] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'U', M, N, A, LDA, VT, LDVT )	CALL ZLACPY( 'U', M, N, A, LDA, VT, LDVT )
CALL ZUNGBR( 'P', M, N, M, VT, LDVT, WORK( ITAUP ),	CALL ZUNGBR( 'P', M, N, M, VT, LDVT, WORK( ITAUP ),
Line 1670	Line 1909
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need M [e] + MM [RVT] + MM [RU] + BDSPAC
*	*
IRVT = NRWORK	IRVT = NRWORK
IRU = IRVT + M*M	IRU = IRVT + M*M
Line 1682	Line 1921
*	*
* Multiply Q in U by real matrix RWORK(IRU), storing the	* Multiply Q in U by real matrix RWORK(IRU), storing the
* result in A, copying to U	* result in A, copying to U
* (CWorkspace: need 0)	* CWorkspace: need 0
* (Rworkspace: need 3MM)	* RWorkspace: need M [e] + MM [RVT] + MM [RU] + 2MM [rwork]
*	*
CALL ZLACRM( M, M, U, LDU, RWORK( IRU ), M, A, LDA,	CALL ZLACRM( M, M, U, LDU, RWORK( IRU ), M, A, LDA,
$ RWORK( NRWORK ) )	$ RWORK( NRWORK ) )
Line 1691	Line 1930
*	*
* Multiply real matrix RWORK(IRVT) by P**H in VT,	* Multiply real matrix RWORK(IRVT) by P**H in VT,
* storing the result in A, copying to VT	* storing the result in A, copying to VT
* (Cworkspace: need 0)	* CWorkspace: need 0
* (Rworkspace: need MM+2M*N)	* RWorkspace: need M [e] + MM [RVT] + 2MN [rwork] < M + 5MM since N < 2M here
*	*
NRWORK = IRU	NRWORK = IRU
CALL ZLARCM( M, N, RWORK( IRVT ), M, VT, LDVT, A, LDA,	CALL ZLARCM( M, N, RWORK( IRVT ), M, VT, LDVT, A, LDA,
Line 1700	Line 1939
CALL ZLACPY( 'F', M, N, A, LDA, VT, LDVT )	CALL ZLACPY( 'F', M, N, A, LDA, VT, LDVT )
ELSE	ELSE
*	*
	* Path 5ta (N >> M, JOBZ='A')
* Copy A to U, generate Q	* Copy A to U, generate Q
* (Cworkspace: need 2M, prefer M+MNB)	* CWorkspace: need 2*M [tauq, taup] + M [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2M [tauq, taup] + MNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'L', M, M, A, LDA, U, LDU )	CALL ZLACPY( 'L', M, M, A, LDA, U, LDU )
CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ),	CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
* Copy A to VT, generate P**H	* Copy A to VT, generate P**H
* (Cworkspace: need 2M, prefer M+MNB)	* CWorkspace: need 2*M [tauq, taup] + N [work]
* (Rworkspace: 0)	* CWorkspace: prefer 2M [tauq, taup] + NNB [work]
	* RWorkspace: need 0
*	*
CALL ZLACPY( 'U', M, N, A, LDA, VT, LDVT )	CALL ZLACPY( 'U', M, N, A, LDA, VT, LDVT )
CALL ZUNGBR( 'P', N, N, M, VT, LDVT, WORK( ITAUP ),	CALL ZUNGBR( 'P', N, N, M, VT, LDVT, WORK( ITAUP ),
Line 1719	Line 1961
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need M [e] + MM [RVT] + MM [RU] + BDSPAC
*	*
IRVT = NRWORK	IRVT = NRWORK
IRU = IRVT + M*M	IRU = IRVT + M*M
Line 1731	Line 1973
*	*
* Multiply Q in U by real matrix RWORK(IRU), storing the	* Multiply Q in U by real matrix RWORK(IRU), storing the
* result in A, copying to U	* result in A, copying to U
* (CWorkspace: need 0)	* CWorkspace: need 0
* (Rworkspace: need 3MM)	* RWorkspace: need M [e] + MM [RVT] + MM [RU] + 2MM [rwork]
*	*
CALL ZLACRM( M, M, U, LDU, RWORK( IRU ), M, A, LDA,	CALL ZLACRM( M, M, U, LDU, RWORK( IRU ), M, A, LDA,
$ RWORK( NRWORK ) )	$ RWORK( NRWORK ) )
Line 1740	Line 1982
*	*
* Multiply real matrix RWORK(IRVT) by P**H in VT,	* Multiply real matrix RWORK(IRVT) by P**H in VT,
* storing the result in A, copying to VT	* storing the result in A, copying to VT
* (Cworkspace: need 0)	* CWorkspace: need 0
* (Rworkspace: need MM+2M*N)	* RWorkspace: need M [e] + MM [RVT] + 2MN [rwork] < M + 5MM since N < 2M here
*	*
	NRWORK = IRU
CALL ZLARCM( M, N, RWORK( IRVT ), M, VT, LDVT, A, LDA,	CALL ZLARCM( M, N, RWORK( IRVT ), M, VT, LDVT, A, LDA,
$ RWORK( NRWORK ) )	$ RWORK( NRWORK ) )
CALL ZLACPY( 'F', M, N, A, LDA, VT, LDVT )	CALL ZLACPY( 'F', M, N, A, LDA, VT, LDVT )
Line 1752	Line 1995
*	*
* N .LT. MNTHR2	* N .LT. MNTHR2
*	*
* Path 6t (N greater than M, but not much larger)	* Path 6t (N > M, but not much larger)
* Reduce to bidiagonal form without LQ decomposition	* Reduce to bidiagonal form without LQ decomposition
* Use ZUNMBR to compute singular vectors	* Use ZUNMBR to compute singular vectors
*	*
Line 1763	Line 2006
NWORK = ITAUP + M	NWORK = ITAUP + M
*	*
* Bidiagonalize A	* Bidiagonalize A
* (CWorkspace: need 2M+N, prefer 2M+(M+N)*NB)	* CWorkspace: need 2*M [tauq, taup] + N [work]
* (RWorkspace: M)	* CWorkspace: prefer 2M [tauq, taup] + (M+N)NB [work]
	* RWorkspace: need M [e]
*	*
CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),	CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ),
$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,	$ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1,
$ IERR )	$ IERR )
IF( WNTQN ) THEN	IF( WNTQN ) THEN
*	*
	* Path 6tn (N > M, JOBZ='N')
* Compute singular values only	* Compute singular values only
* (Cworkspace: 0)	* CWorkspace: need 0
* (Rworkspace: need BDSPAN)	* RWorkspace: need M [e] + BDSPAC
*	*
CALL DBDSDC( 'L', 'N', M, S, RWORK( IE ), DUM, 1, DUM, 1,	CALL DBDSDC( 'L', 'N', M, S, RWORK( IE ), DUM,1,DUM,1,
$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )	$ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO )
ELSE IF( WNTQO ) THEN	ELSE IF( WNTQO ) THEN
	* Path 6to (N > M, JOBZ='O')
LDWKVT = M	LDWKVT = M
IVT = NWORK	IVT = NWORK
IF( LWORK.GE.MN+3M ) THEN	IF( LWORK .GE. MN + 3M ) THEN
*	*
* WORK( IVT ) is M by N	* WORK( IVT ) is M by N
*	*
Line 1791	Line 2037
*	*
* WORK( IVT ) is M by CHUNK	* WORK( IVT ) is M by CHUNK
*	*
CHUNK = ( LWORK-3*M ) / M	CHUNK = ( LWORK - 3*M ) / M
NWORK = IVT + LDWKVT*CHUNK	NWORK = IVT + LDWKVT*CHUNK
END IF	END IF
*	*
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need M [e] + MM [RVT] + MM [RU] + BDSPAC
*	*
IRVT = NRWORK	IRVT = NRWORK
IRU = IRVT + M*M	IRU = IRVT + M*M
Line 1810	Line 2056
*	*
* Copy real matrix RWORK(IRU) to complex matrix U	* Copy real matrix RWORK(IRU) to complex matrix U
* Overwrite U by left singular vectors of A	* Overwrite U by left singular vectors of A
* (Cworkspace: need 2M, prefer M+MNB)	* CWorkspace: need 2M [tauq, taup] + MM [VT] + M [work]
* (Rworkspace: need 0)	* CWorkspace: prefer 2M [tauq, taup] + MM [VT] + M*NB [work]
	* RWorkspace: need M [e] + MM [RVT] + MM [RU]
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )	CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )
CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA,	CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA,
$ WORK( ITAUQ ), U, LDU, WORK( NWORK ),	$ WORK( ITAUQ ), U, LDU, WORK( NWORK ),
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
*	*
IF( LWORK.GE.MN+3M ) THEN	IF( LWORK .GE. MN + 3M ) THEN
*	*
* Copy real matrix RWORK(IRVT) to complex matrix WORK(IVT)	* Path 6to-fast
* Overwrite WORK(IVT) by right singular vectors of A,	* Copy real matrix RWORK(IRVT) to complex matrix WORK(IVT)
* copying to A	* Overwrite WORK(IVT) by right singular vectors of A,
* (Cworkspace: need MN+2M, prefer MN+M+MNB)	* copying to A
* (Rworkspace: need 0)	* CWorkspace: need 2M [tauq, taup] + MN [VT] + M [work]
	* CWorkspace: prefer 2M [tauq, taup] + MN [VT] + M*NB [work]
	* RWorkspace: need M [e] + M*M [RVT]
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, WORK( IVT ),	CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, WORK( IVT ),
$ LDWKVT )	$ LDWKVT )
Line 1834	Line 2083
CALL ZLACPY( 'F', M, N, WORK( IVT ), LDWKVT, A, LDA )	CALL ZLACPY( 'F', M, N, WORK( IVT ), LDWKVT, A, LDA )
ELSE	ELSE
*	*
	* Path 6to-slow
* Generate P**H in A	* Generate P**H in A
* (Cworkspace: need 2M, prefer M+MNB)	* CWorkspace: need 2M [tauq, taup] + MM [VT] + M [work]
* (Rworkspace: need 0)	* CWorkspace: prefer 2M [tauq, taup] + MM [VT] + M*NB [work]
	* RWorkspace: need 0
*	*
CALL ZUNGBR( 'P', M, N, M, A, LDA, WORK( ITAUP ),	CALL ZUNGBR( 'P', M, N, M, A, LDA, WORK( ITAUP ),
$ WORK( NWORK ), LWORK-NWORK+1, IERR )	$ WORK( NWORK ), LWORK-NWORK+1, IERR )
*	*
* Multiply Q in A by real matrix RWORK(IRU), storing the	* Multiply Q in A by real matrix RWORK(IRU), storing the
* result in WORK(IU), copying to A	* result in WORK(IU), copying to A
* (CWorkspace: need MM, prefer MN)	* CWorkspace: need 2M [tauq, taup] + MM [VT]
* (Rworkspace: need 3MM, prefer MM+2M*N)	* CWorkspace: prefer 2M [tauq, taup] + MN [VT]
	* RWorkspace: need M [e] + MM [RVT] + 2M*M [rwork]
	* RWorkspace: prefer M [e] + MM [RVT] + 2MN [rwork] < M + 5MM since N < 2M here
*	*
NRWORK = IRU	NRWORK = IRU
DO 60 I = 1, N, CHUNK	DO 60 I = 1, N, CHUNK
Line 1858	Line 2111
END IF	END IF
ELSE IF( WNTQS ) THEN	ELSE IF( WNTQS ) THEN
*	*
	* Path 6ts (N > M, JOBZ='S')
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need M [e] + MM [RVT] + MM [RU] + BDSPAC
*	*
IRVT = NRWORK	IRVT = NRWORK
IRU = IRVT + M*M	IRU = IRVT + M*M
Line 1873	Line 2127
*	*
* Copy real matrix RWORK(IRU) to complex matrix U	* Copy real matrix RWORK(IRU) to complex matrix U
* Overwrite U by left singular vectors of A	* Overwrite U by left singular vectors of A
* (CWorkspace: need 3M, prefer 2M+M*NB)	* CWorkspace: need 2*M [tauq, taup] + M [work]
* (RWorkspace: M*M)	* CWorkspace: prefer 2M [tauq, taup] + MNB [work]
	* RWorkspace: need M [e] + MM [RVT] + MM [RU]
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )	CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )
CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA,	CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA,
Line 1883	Line 2138
*	*
* Copy real matrix RWORK(IRVT) to complex matrix VT	* Copy real matrix RWORK(IRVT) to complex matrix VT
* Overwrite VT by right singular vectors of A	* Overwrite VT by right singular vectors of A
* (CWorkspace: need 3M, prefer 2M+M*NB)	* CWorkspace: need 2*M [tauq, taup] + M [work]
* (RWorkspace: M*M)	* CWorkspace: prefer 2M [tauq, taup] + MNB [work]
	* RWorkspace: need M [e] + M*M [RVT]
*	*
CALL ZLASET( 'F', M, N, CZERO, CZERO, VT, LDVT )	CALL ZLASET( 'F', M, N, CZERO, CZERO, VT, LDVT )
CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, VT, LDVT )	CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, VT, LDVT )
Line 1893	Line 2149
$ LWORK-NWORK+1, IERR )	$ LWORK-NWORK+1, IERR )
ELSE	ELSE
*	*
	* Path 6ta (N > M, JOBZ='A')
* Perform bidiagonal SVD, computing left singular vectors	* Perform bidiagonal SVD, computing left singular vectors
* of bidiagonal matrix in RWORK(IRU) and computing right	* of bidiagonal matrix in RWORK(IRU) and computing right
* singular vectors of bidiagonal matrix in RWORK(IRVT)	* singular vectors of bidiagonal matrix in RWORK(IRVT)
* (CWorkspace: need 0)	* CWorkspace: need 0
* (RWorkspace: need BDSPAC)	* RWorkspace: need M [e] + MM [RVT] + MM [RU] + BDSPAC
*	*
IRVT = NRWORK	IRVT = NRWORK
IRU = IRVT + M*M	IRU = IRVT + M*M
Line 1909	Line 2166
*	*
* Copy real matrix RWORK(IRU) to complex matrix U	* Copy real matrix RWORK(IRU) to complex matrix U
* Overwrite U by left singular vectors of A	* Overwrite U by left singular vectors of A
* (CWorkspace: need 3M, prefer 2M+M*NB)	* CWorkspace: need 2*M [tauq, taup] + M [work]
* (RWorkspace: M*M)	* CWorkspace: prefer 2M [tauq, taup] + MNB [work]
	* RWorkspace: need M [e] + MM [RVT] + MM [RU]
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )	CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU )
CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA,	CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA,
Line 1923	Line 2181
*	*
* Copy real matrix RWORK(IRVT) to complex matrix VT	* Copy real matrix RWORK(IRVT) to complex matrix VT
* Overwrite VT by right singular vectors of A	* Overwrite VT by right singular vectors of A
* (CWorkspace: need 2M+N, prefer 2M+N*NB)	* CWorkspace: need 2*M [tauq, taup] + N [work]
* (RWorkspace: M*M)	* CWorkspace: prefer 2M [tauq, taup] + NNB [work]
	* RWorkspace: need M [e] + M*M [RVT]
*	*
CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, VT, LDVT )	CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, VT, LDVT )
CALL ZUNMBR( 'P', 'R', 'C', N, N, M, A, LDA,	CALL ZUNMBR( 'P', 'R', 'C', N, N, M, A, LDA,
Line 1955	Line 2214
*	*
* Return optimal workspace in WORK(1)	* Return optimal workspace in WORK(1)
*	*
WORK( 1 ) = MAXWRK	WORK( 1 ) = DROUNDUP_LWORK( MAXWRK )
*	*
RETURN	RETURN
*	*

CVSweb interface <joel.bertrand@systella.fr>

Removed from v.1.8
changed lines
	Added in v.1.20