1*> \brief \b ZGESDD 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8*> \htmlonly 9*> Download ZGESDD + dependencies 10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgesdd.f"> 11*> [TGZ]</a> 12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgesdd.f"> 13*> [ZIP]</a> 14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgesdd.f"> 15*> [TXT]</a> 16*> \endhtmlonly 17* 18* Definition: 19* =========== 20* 21* SUBROUTINE ZGESDD( JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, 22* LWORK, RWORK, IWORK, INFO ) 23* 24* .. Scalar Arguments .. 25* CHARACTER JOBZ 26* INTEGER INFO, LDA, LDU, LDVT, LWORK, M, N 27* .. 28* .. Array Arguments .. 29* INTEGER IWORK( * ) 30* DOUBLE PRECISION RWORK( * ), S( * ) 31* COMPLEX*16 A( LDA, * ), U( LDU, * ), VT( LDVT, * ), 32* $ WORK( * ) 33* .. 34* 35* 36*> \par Purpose: 37* ============= 38*> 39*> \verbatim 40*> 41*> ZGESDD computes the singular value decomposition (SVD) of a complex 42*> M-by-N matrix A, optionally computing the left and/or right singular 43*> vectors, by using divide-and-conquer method. The SVD is written 44*> 45*> A = U * SIGMA * conjugate-transpose(V) 46*> 47*> where SIGMA is an M-by-N matrix which is zero except for its 48*> min(m,n) diagonal elements, U is an M-by-M unitary matrix, and 49*> V is an N-by-N unitary matrix. The diagonal elements of SIGMA 50*> are the singular values of A; they are real and non-negative, and 51*> are returned in descending order. The first min(m,n) columns of 52*> U and V are the left and right singular vectors of A. 53*> 54*> Note that the routine returns VT = V**H, not V. 55*> 56*> The divide and conquer algorithm makes very mild assumptions about 57*> floating point arithmetic. It will work on machines with a guard 58*> digit in add/subtract, or on those binary machines without guard 59*> digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or 60*> Cray-2. It could conceivably fail on hexadecimal or decimal machines 61*> without guard digits, but we know of none. 62*> \endverbatim 63* 64* Arguments: 65* ========== 66* 67*> \param[in] JOBZ 68*> \verbatim 69*> JOBZ is CHARACTER*1 70*> Specifies options for computing all or part of the matrix U: 71*> = 'A': all M columns of U and all N rows of V**H are 72*> returned in the arrays U and VT; 73*> = 'S': the first min(M,N) columns of U and the first 74*> min(M,N) rows of V**H are returned in the arrays U 75*> and VT; 76*> = 'O': If M >= N, the first N columns of U are overwritten 77*> in the array A and all rows of V**H are returned in 78*> the array VT; 79*> otherwise, all columns of U are returned in the 80*> array U and the first M rows of V**H are overwritten 81*> in the array A; 82*> = 'N': no columns of U or rows of V**H are computed. 83*> \endverbatim 84*> 85*> \param[in] M 86*> \verbatim 87*> M is INTEGER 88*> The number of rows of the input matrix A. M >= 0. 89*> \endverbatim 90*> 91*> \param[in] N 92*> \verbatim 93*> N is INTEGER 94*> The number of columns of the input matrix A. N >= 0. 95*> \endverbatim 96*> 97*> \param[in,out] A 98*> \verbatim 99*> A is COMPLEX*16 array, dimension (LDA,N) 100*> On entry, the M-by-N matrix A. 101*> On exit, 102*> if JOBZ = 'O', A is overwritten with the first N columns 103*> of U (the left singular vectors, stored 104*> columnwise) if M >= N; 105*> A is overwritten with the first M rows 106*> of V**H (the right singular vectors, stored 107*> rowwise) otherwise. 108*> if JOBZ .ne. 'O', the contents of A are destroyed. 109*> \endverbatim 110*> 111*> \param[in] LDA 112*> \verbatim 113*> LDA is INTEGER 114*> The leading dimension of the array A. LDA >= max(1,M). 115*> \endverbatim 116*> 117*> \param[out] S 118*> \verbatim 119*> S is DOUBLE PRECISION array, dimension (min(M,N)) 120*> The singular values of A, sorted so that S(i) >= S(i+1). 121*> \endverbatim 122*> 123*> \param[out] U 124*> \verbatim 125*> U is COMPLEX*16 array, dimension (LDU,UCOL) 126*> UCOL = M if JOBZ = 'A' or JOBZ = 'O' and M < N; 127*> UCOL = min(M,N) if JOBZ = 'S'. 128*> If JOBZ = 'A' or JOBZ = 'O' and M < N, U contains the M-by-M 129*> unitary matrix U; 130*> if JOBZ = 'S', U contains the first min(M,N) columns of U 131*> (the left singular vectors, stored columnwise); 132*> if JOBZ = 'O' and M >= N, or JOBZ = 'N', U is not referenced. 133*> \endverbatim 134*> 135*> \param[in] LDU 136*> \verbatim 137*> LDU is INTEGER 138*> The leading dimension of the array U. LDU >= 1; if 139*> JOBZ = 'S' or 'A' or JOBZ = 'O' and M < N, LDU >= M. 140*> \endverbatim 141*> 142*> \param[out] VT 143*> \verbatim 144*> VT is COMPLEX*16 array, dimension (LDVT,N) 145*> If JOBZ = 'A' or JOBZ = 'O' and M >= N, VT contains the 146*> N-by-N unitary matrix V**H; 147*> if JOBZ = 'S', VT contains the first min(M,N) rows of 148*> V**H (the right singular vectors, stored rowwise); 149*> if JOBZ = 'O' and M < N, or JOBZ = 'N', VT is not referenced. 150*> \endverbatim 151*> 152*> \param[in] LDVT 153*> \verbatim 154*> LDVT is INTEGER 155*> The leading dimension of the array VT. LDVT >= 1; if 156*> JOBZ = 'A' or JOBZ = 'O' and M >= N, LDVT >= N; 157*> if JOBZ = 'S', LDVT >= min(M,N). 158*> \endverbatim 159*> 160*> \param[out] WORK 161*> \verbatim 162*> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK)) 163*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK. 164*> \endverbatim 165*> 166*> \param[in] LWORK 167*> \verbatim 168*> LWORK is INTEGER 169*> The dimension of the array WORK. LWORK >= 1. 170*> if JOBZ = 'N', LWORK >= 2*min(M,N)+max(M,N). 171*> if JOBZ = 'O', 172*> LWORK >= 2*min(M,N)*min(M,N)+2*min(M,N)+max(M,N). 173*> if JOBZ = 'S' or 'A', 174*> LWORK >= min(M,N)*min(M,N)+2*min(M,N)+max(M,N). 175*> For good performance, LWORK should generally be larger. 176*> 177*> If LWORK = -1, a workspace query is assumed. The optimal 178*> size for the WORK array is calculated and stored in WORK(1), 179*> and no other work except argument checking is performed. 180*> \endverbatim 181*> 182*> \param[out] RWORK 183*> \verbatim 184*> RWORK is DOUBLE PRECISION array, dimension (MAX(1,LRWORK)) 185*> If JOBZ = 'N', LRWORK >= 7*min(M,N). 186*> Otherwise, 187*> LRWORK >= min(M,N)*max(5*min(M,N)+7,2*max(M,N)+2*min(M,N)+1) 188*> \endverbatim 189*> 190*> \param[out] IWORK 191*> \verbatim 192*> IWORK is INTEGER array, dimension (8*min(M,N)) 193*> \endverbatim 194*> 195*> \param[out] INFO 196*> \verbatim 197*> INFO is INTEGER 198*> = 0: successful exit. 199*> < 0: if INFO = -i, the i-th argument had an illegal value. 200*> > 0: The updating process of DBDSDC did not converge. 201*> \endverbatim 202* 203* Authors: 204* ======== 205* 206*> \author Univ. of Tennessee 207*> \author Univ. of California Berkeley 208*> \author Univ. of Colorado Denver 209*> \author NAG Ltd. 210* 211*> \date November 2015 212* 213*> \ingroup complex16GEsing 214* 215*> \par Contributors: 216* ================== 217*> 218*> Ming Gu and Huan Ren, Computer Science Division, University of 219*> California at Berkeley, USA 220*> 221* ===================================================================== 222 SUBROUTINE ZGESDD( JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, 223 $ LWORK, RWORK, IWORK, INFO ) 224* 225* -- LAPACK driver routine (version 3.6.0) -- 226* -- LAPACK is a software package provided by Univ. of Tennessee, -- 227* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 228* November 2015 229* 230* .. Scalar Arguments .. 231 CHARACTER JOBZ 232 INTEGER INFO, LDA, LDU, LDVT, LWORK, M, N 233* .. 234* .. Array Arguments .. 235 INTEGER IWORK( * ) 236 DOUBLE PRECISION RWORK( * ), S( * ) 237 COMPLEX*16 A( LDA, * ), U( LDU, * ), VT( LDVT, * ), 238 $ WORK( * ) 239* .. 240* 241* ===================================================================== 242* 243* .. Parameters .. 244 INTEGER LQUERV 245 PARAMETER ( LQUERV = -1 ) 246 COMPLEX*16 CZERO, CONE 247 PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ), 248 $ CONE = ( 1.0D+0, 0.0D+0 ) ) 249 DOUBLE PRECISION ZERO, ONE 250 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 251* .. 252* .. Local Scalars .. 253 LOGICAL WNTQA, WNTQAS, WNTQN, WNTQO, WNTQS 254 INTEGER BLK, CHUNK, I, IE, IERR, IL, IR, IRU, IRVT, 255 $ ISCL, ITAU, ITAUP, ITAUQ, IU, IVT, LDWKVT, 256 $ LDWRKL, LDWRKR, LDWRKU, MAXWRK, MINMN, MINWRK, 257 $ MNTHR1, MNTHR2, NRWORK, NWORK, WRKBL 258 DOUBLE PRECISION ANRM, BIGNUM, EPS, SMLNUM 259* .. 260* .. Local Arrays .. 261 INTEGER IDUM( 1 ) 262 DOUBLE PRECISION DUM( 1 ) 263* .. 264* .. External Subroutines .. 265 EXTERNAL DBDSDC, DLASCL, XERBLA, ZGEBRD, ZGELQF, ZGEMM, 266 $ ZGEQRF, ZLACP2, ZLACPY, ZLACRM, ZLARCM, ZLASCL, 267 $ ZLASET, ZUNGBR, ZUNGLQ, ZUNGQR, ZUNMBR 268* .. 269* .. External Functions .. 270 LOGICAL LSAME 271 INTEGER ILAENV 272 DOUBLE PRECISION DLAMCH, ZLANGE 273 EXTERNAL LSAME, ILAENV, DLAMCH, ZLANGE 274* .. 275* .. Intrinsic Functions .. 276 INTRINSIC INT, MAX, MIN, SQRT 277* .. 278* .. Executable Statements .. 279* 280* Test the input arguments 281* 282 INFO = 0 283 MINMN = MIN( M, N ) 284 MNTHR1 = INT( MINMN*17.0D0 / 9.0D0 ) 285 MNTHR2 = INT( MINMN*5.0D0 / 3.0D0 ) 286 WNTQA = LSAME( JOBZ, 'A' ) 287 WNTQS = LSAME( JOBZ, 'S' ) 288 WNTQAS = WNTQA .OR. WNTQS 289 WNTQO = LSAME( JOBZ, 'O' ) 290 WNTQN = LSAME( JOBZ, 'N' ) 291 MINWRK = 1 292 MAXWRK = 1 293* 294 IF( .NOT.( WNTQA .OR. WNTQS .OR. WNTQO .OR. WNTQN ) ) THEN 295 INFO = -1 296 ELSE IF( M.LT.0 ) THEN 297 INFO = -2 298 ELSE IF( N.LT.0 ) THEN 299 INFO = -3 300 ELSE IF( LDA.LT.MAX( 1, M ) ) THEN 301 INFO = -5 302 ELSE IF( LDU.LT.1 .OR. ( WNTQAS .AND. LDU.LT.M ) .OR. 303 $ ( WNTQO .AND. M.LT.N .AND. LDU.LT.M ) ) THEN 304 INFO = -8 305 ELSE IF( LDVT.LT.1 .OR. ( WNTQA .AND. LDVT.LT.N ) .OR. 306 $ ( WNTQS .AND. LDVT.LT.MINMN ) .OR. 307 $ ( WNTQO .AND. M.GE.N .AND. LDVT.LT.N ) ) THEN 308 INFO = -10 309 END IF 310* 311* Compute workspace 312* (Note: Comments in the code beginning "Workspace:" describe the 313* minimal amount of workspace needed at that point in the code, 314* as well as the preferred amount for good performance. 315* CWorkspace refers to complex workspace, and RWorkspace to 316* real workspace. NB refers to the optimal block size for the 317* immediately following subroutine, as returned by ILAENV.) 318* 319 IF( INFO.EQ.0 .AND. M.GT.0 .AND. N.GT.0 ) THEN 320 IF( M.GE.N ) THEN 321* 322* There is no complex work space needed for bidiagonal SVD 323* The real work space needed for bidiagonal SVD is BDSPAC 324* for computing singular values and singular vectors; BDSPAN 325* for computing singular values only. 326* BDSPAC = 5*N*N + 7*N 327* BDSPAN = MAX(7*N+4, 3*N+2+SMLSIZ*(SMLSIZ+8)) 328* 329 IF( M.GE.MNTHR1 ) THEN 330 IF( WNTQN ) THEN 331* 332* Path 1 (M much larger than N, JOBZ='N') 333* 334 MAXWRK = N + N*ILAENV( 1, 'ZGEQRF', ' ', M, N, -1, 335 $ -1 ) 336 MAXWRK = MAX( MAXWRK, 2*N+2*N* 337 $ ILAENV( 1, 'ZGEBRD', ' ', N, N, -1, -1 ) ) 338 MINWRK = 3*N 339 ELSE IF( WNTQO ) THEN 340* 341* Path 2 (M much larger than N, JOBZ='O') 342* 343 WRKBL = N + N*ILAENV( 1, 'ZGEQRF', ' ', M, N, -1, -1 ) 344 WRKBL = MAX( WRKBL, N+N*ILAENV( 1, 'ZUNGQR', ' ', M, 345 $ N, N, -1 ) ) 346 WRKBL = MAX( WRKBL, 2*N+2*N* 347 $ ILAENV( 1, 'ZGEBRD', ' ', N, N, -1, -1 ) ) 348 WRKBL = MAX( WRKBL, 2*N+N* 349 $ ILAENV( 1, 'ZUNMBR', 'QLN', N, N, N, -1 ) ) 350 WRKBL = MAX( WRKBL, 2*N+N* 351 $ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) ) 352 MAXWRK = M*N + N*N + WRKBL 353 MINWRK = 2*N*N + 3*N 354 ELSE IF( WNTQS ) THEN 355* 356* Path 3 (M much larger than N, JOBZ='S') 357* 358 WRKBL = N + N*ILAENV( 1, 'ZGEQRF', ' ', M, N, -1, -1 ) 359 WRKBL = MAX( WRKBL, N+N*ILAENV( 1, 'ZUNGQR', ' ', M, 360 $ N, N, -1 ) ) 361 WRKBL = MAX( WRKBL, 2*N+2*N* 362 $ ILAENV( 1, 'ZGEBRD', ' ', N, N, -1, -1 ) ) 363 WRKBL = MAX( WRKBL, 2*N+N* 364 $ ILAENV( 1, 'ZUNMBR', 'QLN', N, N, N, -1 ) ) 365 WRKBL = MAX( WRKBL, 2*N+N* 366 $ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) ) 367 MAXWRK = N*N + WRKBL 368 MINWRK = N*N + 3*N 369 ELSE IF( WNTQA ) THEN 370* 371* Path 4 (M much larger than N, JOBZ='A') 372* 373 WRKBL = N + N*ILAENV( 1, 'ZGEQRF', ' ', M, N, -1, -1 ) 374 WRKBL = MAX( WRKBL, N+M*ILAENV( 1, 'ZUNGQR', ' ', M, 375 $ M, N, -1 ) ) 376 WRKBL = MAX( WRKBL, 2*N+2*N* 377 $ ILAENV( 1, 'ZGEBRD', ' ', N, N, -1, -1 ) ) 378 WRKBL = MAX( WRKBL, 2*N+N* 379 $ ILAENV( 1, 'ZUNMBR', 'QLN', N, N, N, -1 ) ) 380 WRKBL = MAX( WRKBL, 2*N+N* 381 $ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) ) 382 MAXWRK = N*N + WRKBL 383 MINWRK = N*N + 2*N + M 384 END IF 385 ELSE IF( M.GE.MNTHR2 ) THEN 386* 387* Path 5 (M much larger than N, but not as much as MNTHR1) 388* 389 MAXWRK = 2*N + ( M+N )*ILAENV( 1, 'ZGEBRD', ' ', M, N, 390 $ -1, -1 ) 391 MINWRK = 2*N + M 392 IF( WNTQO ) THEN 393 MAXWRK = MAX( MAXWRK, 2*N+N* 394 $ ILAENV( 1, 'ZUNGBR', 'P', N, N, N, -1 ) ) 395 MAXWRK = MAX( MAXWRK, 2*N+N* 396 $ ILAENV( 1, 'ZUNGBR', 'Q', M, N, N, -1 ) ) 397 MAXWRK = MAXWRK + M*N 398 MINWRK = MINWRK + N*N 399 ELSE IF( WNTQS ) THEN 400 MAXWRK = MAX( MAXWRK, 2*N+N* 401 $ ILAENV( 1, 'ZUNGBR', 'P', N, N, N, -1 ) ) 402 MAXWRK = MAX( MAXWRK, 2*N+N* 403 $ ILAENV( 1, 'ZUNGBR', 'Q', M, N, N, -1 ) ) 404 ELSE IF( WNTQA ) THEN 405 MAXWRK = MAX( MAXWRK, 2*N+N* 406 $ ILAENV( 1, 'ZUNGBR', 'P', N, N, N, -1 ) ) 407 MAXWRK = MAX( MAXWRK, 2*N+M* 408 $ ILAENV( 1, 'ZUNGBR', 'Q', M, M, N, -1 ) ) 409 END IF 410 ELSE 411* 412* Path 6 (M at least N, but not much larger) 413* 414 MAXWRK = 2*N + ( M+N )*ILAENV( 1, 'ZGEBRD', ' ', M, N, 415 $ -1, -1 ) 416 MINWRK = 2*N + M 417 IF( WNTQO ) THEN 418 MAXWRK = MAX( MAXWRK, 2*N+N* 419 $ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) ) 420 MAXWRK = MAX( MAXWRK, 2*N+N* 421 $ ILAENV( 1, 'ZUNMBR', 'QLN', M, N, N, -1 ) ) 422 MAXWRK = MAXWRK + M*N 423 MINWRK = MINWRK + N*N 424 ELSE IF( WNTQS ) THEN 425 MAXWRK = MAX( MAXWRK, 2*N+N* 426 $ ILAENV( 1, 'ZUNMBR', 'PRC', N, N, N, -1 ) ) 427 MAXWRK = MAX( MAXWRK, 2*N+N* 428 $ ILAENV( 1, 'ZUNMBR', 'QLN', M, N, N, -1 ) ) 429 ELSE IF( WNTQA ) THEN 430 MAXWRK = MAX( MAXWRK, 2*N+N* 431 $ ILAENV( 1, 'ZUNGBR', 'PRC', N, N, N, -1 ) ) 432 MAXWRK = MAX( MAXWRK, 2*N+M* 433 $ ILAENV( 1, 'ZUNGBR', 'QLN', M, M, N, -1 ) ) 434 END IF 435 END IF 436 ELSE 437* 438* There is no complex work space needed for bidiagonal SVD 439* The real work space needed for bidiagonal SVD is BDSPAC 440* for computing singular values and singular vectors; BDSPAN 441* for computing singular values only. 442* BDSPAC = 5*M*M + 7*M 443* BDSPAN = MAX(7*M+4, 3*M+2+SMLSIZ*(SMLSIZ+8)) 444* 445 IF( N.GE.MNTHR1 ) THEN 446 IF( WNTQN ) THEN 447* 448* Path 1t (N much larger than M, JOBZ='N') 449* 450 MAXWRK = M + M*ILAENV( 1, 'ZGELQF', ' ', M, N, -1, 451 $ -1 ) 452 MAXWRK = MAX( MAXWRK, 2*M+2*M* 453 $ ILAENV( 1, 'ZGEBRD', ' ', M, M, -1, -1 ) ) 454 MINWRK = 3*M 455 ELSE IF( WNTQO ) THEN 456* 457* Path 2t (N much larger than M, JOBZ='O') 458* 459 WRKBL = M + M*ILAENV( 1, 'ZGELQF', ' ', M, N, -1, -1 ) 460 WRKBL = MAX( WRKBL, M+M*ILAENV( 1, 'ZUNGLQ', ' ', M, 461 $ N, M, -1 ) ) 462 WRKBL = MAX( WRKBL, 2*M+2*M* 463 $ ILAENV( 1, 'ZGEBRD', ' ', M, M, -1, -1 ) ) 464 WRKBL = MAX( WRKBL, 2*M+M* 465 $ ILAENV( 1, 'ZUNMBR', 'PRC', M, M, M, -1 ) ) 466 WRKBL = MAX( WRKBL, 2*M+M* 467 $ ILAENV( 1, 'ZUNMBR', 'QLN', M, M, M, -1 ) ) 468 MAXWRK = M*N + M*M + WRKBL 469 MINWRK = 2*M*M + 3*M 470 ELSE IF( WNTQS ) THEN 471* 472* Path 3t (N much larger than M, JOBZ='S') 473* 474 WRKBL = M + M*ILAENV( 1, 'ZGELQF', ' ', M, N, -1, -1 ) 475 WRKBL = MAX( WRKBL, M+M*ILAENV( 1, 'ZUNGLQ', ' ', M, 476 $ N, M, -1 ) ) 477 WRKBL = MAX( WRKBL, 2*M+2*M* 478 $ ILAENV( 1, 'ZGEBRD', ' ', M, M, -1, -1 ) ) 479 WRKBL = MAX( WRKBL, 2*M+M* 480 $ ILAENV( 1, 'ZUNMBR', 'PRC', M, M, M, -1 ) ) 481 WRKBL = MAX( WRKBL, 2*M+M* 482 $ ILAENV( 1, 'ZUNMBR', 'QLN', M, M, M, -1 ) ) 483 MAXWRK = M*M + WRKBL 484 MINWRK = M*M + 3*M 485 ELSE IF( WNTQA ) THEN 486* 487* Path 4t (N much larger than M, JOBZ='A') 488* 489 WRKBL = M + M*ILAENV( 1, 'ZGELQF', ' ', M, N, -1, -1 ) 490 WRKBL = MAX( WRKBL, M+N*ILAENV( 1, 'ZUNGLQ', ' ', N, 491 $ N, M, -1 ) ) 492 WRKBL = MAX( WRKBL, 2*M+2*M* 493 $ ILAENV( 1, 'ZGEBRD', ' ', M, M, -1, -1 ) ) 494 WRKBL = MAX( WRKBL, 2*M+M* 495 $ ILAENV( 1, 'ZUNMBR', 'PRC', M, M, M, -1 ) ) 496 WRKBL = MAX( WRKBL, 2*M+M* 497 $ ILAENV( 1, 'ZUNMBR', 'QLN', M, M, M, -1 ) ) 498 MAXWRK = M*M + WRKBL 499 MINWRK = M*M + 2*M + N 500 END IF 501 ELSE IF( N.GE.MNTHR2 ) THEN 502* 503* Path 5t (N much larger than M, but not as much as MNTHR1) 504* 505 MAXWRK = 2*M + ( M+N )*ILAENV( 1, 'ZGEBRD', ' ', M, N, 506 $ -1, -1 ) 507 MINWRK = 2*M + N 508 IF( WNTQO ) THEN 509 MAXWRK = MAX( MAXWRK, 2*M+M* 510 $ ILAENV( 1, 'ZUNGBR', 'P', M, N, M, -1 ) ) 511 MAXWRK = MAX( MAXWRK, 2*M+M* 512 $ ILAENV( 1, 'ZUNGBR', 'Q', M, M, N, -1 ) ) 513 MAXWRK = MAXWRK + M*N 514 MINWRK = MINWRK + M*M 515 ELSE IF( WNTQS ) THEN 516 MAXWRK = MAX( MAXWRK, 2*M+M* 517 $ ILAENV( 1, 'ZUNGBR', 'P', M, N, M, -1 ) ) 518 MAXWRK = MAX( MAXWRK, 2*M+M* 519 $ ILAENV( 1, 'ZUNGBR', 'Q', M, M, N, -1 ) ) 520 ELSE IF( WNTQA ) THEN 521 MAXWRK = MAX( MAXWRK, 2*M+N* 522 $ ILAENV( 1, 'ZUNGBR', 'P', N, N, M, -1 ) ) 523 MAXWRK = MAX( MAXWRK, 2*M+M* 524 $ ILAENV( 1, 'ZUNGBR', 'Q', M, M, N, -1 ) ) 525 END IF 526 ELSE 527* 528* Path 6t (N greater than M, but not much larger) 529* 530 MAXWRK = 2*M + ( M+N )*ILAENV( 1, 'ZGEBRD', ' ', M, N, 531 $ -1, -1 ) 532 MINWRK = 2*M + N 533 IF( WNTQO ) THEN 534 MAXWRK = MAX( MAXWRK, 2*M+M* 535 $ ILAENV( 1, 'ZUNMBR', 'PRC', M, N, M, -1 ) ) 536 MAXWRK = MAX( MAXWRK, 2*M+M* 537 $ ILAENV( 1, 'ZUNMBR', 'QLN', M, M, N, -1 ) ) 538 MAXWRK = MAXWRK + M*N 539 MINWRK = MINWRK + M*M 540 ELSE IF( WNTQS ) THEN 541 MAXWRK = MAX( MAXWRK, 2*M+M* 542 $ ILAENV( 1, 'ZUNGBR', 'PRC', M, N, M, -1 ) ) 543 MAXWRK = MAX( MAXWRK, 2*M+M* 544 $ ILAENV( 1, 'ZUNGBR', 'QLN', M, M, N, -1 ) ) 545 ELSE IF( WNTQA ) THEN 546 MAXWRK = MAX( MAXWRK, 2*M+N* 547 $ ILAENV( 1, 'ZUNGBR', 'PRC', N, N, M, -1 ) ) 548 MAXWRK = MAX( MAXWRK, 2*M+M* 549 $ ILAENV( 1, 'ZUNGBR', 'QLN', M, M, N, -1 ) ) 550 END IF 551 END IF 552 END IF 553 MAXWRK = MAX( MAXWRK, MINWRK ) 554 END IF 555 IF( INFO.EQ.0 ) THEN 556 WORK( 1 ) = MAXWRK 557 IF( LWORK.LT.MINWRK .AND. LWORK.NE.LQUERV ) 558 $ INFO = -13 559 END IF 560* 561* Quick returns 562* 563 IF( INFO.NE.0 ) THEN 564 CALL XERBLA( 'ZGESDD', -INFO ) 565 RETURN 566 END IF 567 IF( LWORK.EQ.LQUERV ) 568 $ RETURN 569 IF( M.EQ.0 .OR. N.EQ.0 ) THEN 570 RETURN 571 END IF 572* 573* Get machine constants 574* 575 EPS = DLAMCH( 'P' ) 576 SMLNUM = SQRT( DLAMCH( 'S' ) ) / EPS 577 BIGNUM = ONE / SMLNUM 578* 579* Scale A if max element outside range [SMLNUM,BIGNUM] 580* 581 ANRM = ZLANGE( 'M', M, N, A, LDA, DUM ) 582 ISCL = 0 583 IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN 584 ISCL = 1 585 CALL ZLASCL( 'G', 0, 0, ANRM, SMLNUM, M, N, A, LDA, IERR ) 586 ELSE IF( ANRM.GT.BIGNUM ) THEN 587 ISCL = 1 588 CALL ZLASCL( 'G', 0, 0, ANRM, BIGNUM, M, N, A, LDA, IERR ) 589 END IF 590* 591 IF( M.GE.N ) THEN 592* 593* A has at least as many rows as columns. If A has sufficiently 594* more rows than columns, first reduce using the QR 595* decomposition (if sufficient workspace available) 596* 597 IF( M.GE.MNTHR1 ) THEN 598* 599 IF( WNTQN ) THEN 600* 601* Path 1 (M much larger than N, JOBZ='N') 602* No singular vectors to be computed 603* 604 ITAU = 1 605 NWORK = ITAU + N 606* 607* Compute A=Q*R 608* (CWorkspace: need 2*N, prefer N+N*NB) 609* (RWorkspace: need 0) 610* 611 CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ), 612 $ LWORK-NWORK+1, IERR ) 613* 614* Zero out below R 615* 616 CALL ZLASET( 'L', N-1, N-1, CZERO, CZERO, A( 2, 1 ), 617 $ LDA ) 618 IE = 1 619 ITAUQ = 1 620 ITAUP = ITAUQ + N 621 NWORK = ITAUP + N 622* 623* Bidiagonalize R in A 624* (CWorkspace: need 3*N, prefer 2*N+2*N*NB) 625* (RWorkspace: need N) 626* 627 CALL ZGEBRD( N, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ), 628 $ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1, 629 $ IERR ) 630 NRWORK = IE + N 631* 632* Perform bidiagonal SVD, compute singular values only 633* (CWorkspace: 0) 634* (RWorkspace: need BDSPAN) 635* 636 CALL DBDSDC( 'U', 'N', N, S, RWORK( IE ), DUM, 1, DUM, 1, 637 $ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO ) 638* 639 ELSE IF( WNTQO ) THEN 640* 641* Path 2 (M much larger than N, JOBZ='O') 642* N left singular vectors to be overwritten on A and 643* N right singular vectors to be computed in VT 644* 645 IU = 1 646* 647* WORK(IU) is N by N 648* 649 LDWRKU = N 650 IR = IU + LDWRKU*N 651 IF( LWORK.GE.M*N+N*N+3*N ) THEN 652* 653* WORK(IR) is M by N 654* 655 LDWRKR = M 656 ELSE 657 LDWRKR = ( LWORK-N*N-3*N ) / N 658 END IF 659 ITAU = IR + LDWRKR*N 660 NWORK = ITAU + N 661* 662* Compute A=Q*R 663* (CWorkspace: need N*N+2*N, prefer M*N+N+N*NB) 664* (RWorkspace: 0) 665* 666 CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ), 667 $ LWORK-NWORK+1, IERR ) 668* 669* Copy R to WORK( IR ), zeroing out below it 670* 671 CALL ZLACPY( 'U', N, N, A, LDA, WORK( IR ), LDWRKR ) 672 CALL ZLASET( 'L', N-1, N-1, CZERO, CZERO, WORK( IR+1 ), 673 $ LDWRKR ) 674* 675* Generate Q in A 676* (CWorkspace: need 2*N, prefer N+N*NB) 677* (RWorkspace: 0) 678* 679 CALL ZUNGQR( M, N, N, A, LDA, WORK( ITAU ), 680 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 681 IE = 1 682 ITAUQ = ITAU 683 ITAUP = ITAUQ + N 684 NWORK = ITAUP + N 685* 686* Bidiagonalize R in WORK(IR) 687* (CWorkspace: need N*N+3*N, prefer M*N+2*N+2*N*NB) 688* (RWorkspace: need N) 689* 690 CALL ZGEBRD( N, N, WORK( IR ), LDWRKR, S, RWORK( IE ), 691 $ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ), 692 $ LWORK-NWORK+1, IERR ) 693* 694* Perform bidiagonal SVD, computing left singular vectors 695* of R in WORK(IRU) and computing right singular vectors 696* of R in WORK(IRVT) 697* (CWorkspace: need 0) 698* (RWorkspace: need BDSPAC) 699* 700 IRU = IE + N 701 IRVT = IRU + N*N 702 NRWORK = IRVT + N*N 703 CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ), 704 $ N, RWORK( IRVT ), N, DUM, IDUM, 705 $ RWORK( NRWORK ), IWORK, INFO ) 706* 707* Copy real matrix RWORK(IRU) to complex matrix WORK(IU) 708* Overwrite WORK(IU) by the left singular vectors of R 709* (CWorkspace: need 2*N*N+3*N, prefer M*N+N*N+2*N+N*NB) 710* (RWorkspace: 0) 711* 712 CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, WORK( IU ), 713 $ LDWRKU ) 714 CALL ZUNMBR( 'Q', 'L', 'N', N, N, N, WORK( IR ), LDWRKR, 715 $ WORK( ITAUQ ), WORK( IU ), LDWRKU, 716 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 717* 718* Copy real matrix RWORK(IRVT) to complex matrix VT 719* Overwrite VT by the right singular vectors of R 720* (CWorkspace: need N*N+3*N, prefer M*N+2*N+N*NB) 721* (RWorkspace: 0) 722* 723 CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT ) 724 CALL ZUNMBR( 'P', 'R', 'C', N, N, N, WORK( IR ), LDWRKR, 725 $ WORK( ITAUP ), VT, LDVT, WORK( NWORK ), 726 $ LWORK-NWORK+1, IERR ) 727* 728* Multiply Q in A by left singular vectors of R in 729* WORK(IU), storing result in WORK(IR) and copying to A 730* (CWorkspace: need 2*N*N, prefer N*N+M*N) 731* (RWorkspace: 0) 732* 733 DO 10 I = 1, M, LDWRKR 734 CHUNK = MIN( M-I+1, LDWRKR ) 735 CALL ZGEMM( 'N', 'N', CHUNK, N, N, CONE, A( I, 1 ), 736 $ LDA, WORK( IU ), LDWRKU, CZERO, 737 $ WORK( IR ), LDWRKR ) 738 CALL ZLACPY( 'F', CHUNK, N, WORK( IR ), LDWRKR, 739 $ A( I, 1 ), LDA ) 740 10 CONTINUE 741* 742 ELSE IF( WNTQS ) THEN 743* 744* Path 3 (M much larger than N, JOBZ='S') 745* N left singular vectors to be computed in U and 746* N right singular vectors to be computed in VT 747* 748 IR = 1 749* 750* WORK(IR) is N by N 751* 752 LDWRKR = N 753 ITAU = IR + LDWRKR*N 754 NWORK = ITAU + N 755* 756* Compute A=Q*R 757* (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) 758* (RWorkspace: 0) 759* 760 CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ), 761 $ LWORK-NWORK+1, IERR ) 762* 763* Copy R to WORK(IR), zeroing out below it 764* 765 CALL ZLACPY( 'U', N, N, A, LDA, WORK( IR ), LDWRKR ) 766 CALL ZLASET( 'L', N-1, N-1, CZERO, CZERO, WORK( IR+1 ), 767 $ LDWRKR ) 768* 769* Generate Q in A 770* (CWorkspace: need 2*N, prefer N+N*NB) 771* (RWorkspace: 0) 772* 773 CALL ZUNGQR( M, N, N, A, LDA, WORK( ITAU ), 774 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 775 IE = 1 776 ITAUQ = ITAU 777 ITAUP = ITAUQ + N 778 NWORK = ITAUP + N 779* 780* Bidiagonalize R in WORK(IR) 781* (CWorkspace: need N*N+3*N, prefer N*N+2*N+2*N*NB) 782* (RWorkspace: need N) 783* 784 CALL ZGEBRD( N, N, WORK( IR ), LDWRKR, S, RWORK( IE ), 785 $ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ), 786 $ LWORK-NWORK+1, IERR ) 787* 788* Perform bidiagonal SVD, computing left singular vectors 789* of bidiagonal matrix in RWORK(IRU) and computing right 790* singular vectors of bidiagonal matrix in RWORK(IRVT) 791* (CWorkspace: need 0) 792* (RWorkspace: need BDSPAC) 793* 794 IRU = IE + N 795 IRVT = IRU + N*N 796 NRWORK = IRVT + N*N 797 CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ), 798 $ N, RWORK( IRVT ), N, DUM, IDUM, 799 $ RWORK( NRWORK ), IWORK, INFO ) 800* 801* Copy real matrix RWORK(IRU) to complex matrix U 802* Overwrite U by left singular vectors of R 803* (CWorkspace: need N*N+3*N, prefer N*N+2*N+N*NB) 804* (RWorkspace: 0) 805* 806 CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, U, LDU ) 807 CALL ZUNMBR( 'Q', 'L', 'N', N, N, N, WORK( IR ), LDWRKR, 808 $ WORK( ITAUQ ), U, LDU, WORK( NWORK ), 809 $ LWORK-NWORK+1, IERR ) 810* 811* Copy real matrix RWORK(IRVT) to complex matrix VT 812* Overwrite VT by right singular vectors of R 813* (CWorkspace: need N*N+3*N, prefer N*N+2*N+N*NB) 814* (RWorkspace: 0) 815* 816 CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT ) 817 CALL ZUNMBR( 'P', 'R', 'C', N, N, N, WORK( IR ), LDWRKR, 818 $ WORK( ITAUP ), VT, LDVT, WORK( NWORK ), 819 $ LWORK-NWORK+1, IERR ) 820* 821* Multiply Q in A by left singular vectors of R in 822* WORK(IR), storing result in U 823* (CWorkspace: need N*N) 824* (RWorkspace: 0) 825* 826 CALL ZLACPY( 'F', N, N, U, LDU, WORK( IR ), LDWRKR ) 827 CALL ZGEMM( 'N', 'N', M, N, N, CONE, A, LDA, WORK( IR ), 828 $ LDWRKR, CZERO, U, LDU ) 829* 830 ELSE IF( WNTQA ) THEN 831* 832* Path 4 (M much larger than N, JOBZ='A') 833* M left singular vectors to be computed in U and 834* N right singular vectors to be computed in VT 835* 836 IU = 1 837* 838* WORK(IU) is N by N 839* 840 LDWRKU = N 841 ITAU = IU + LDWRKU*N 842 NWORK = ITAU + N 843* 844* Compute A=Q*R, copying result to U 845* (CWorkspace: need 2*N, prefer N+N*NB) 846* (RWorkspace: 0) 847* 848 CALL ZGEQRF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ), 849 $ LWORK-NWORK+1, IERR ) 850 CALL ZLACPY( 'L', M, N, A, LDA, U, LDU ) 851* 852* Generate Q in U 853* (CWorkspace: need N+M, prefer N+M*NB) 854* (RWorkspace: 0) 855* 856 CALL ZUNGQR( M, M, N, U, LDU, WORK( ITAU ), 857 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 858* 859* Produce R in A, zeroing out below it 860* 861 CALL ZLASET( 'L', N-1, N-1, CZERO, CZERO, A( 2, 1 ), 862 $ LDA ) 863 IE = 1 864 ITAUQ = ITAU 865 ITAUP = ITAUQ + N 866 NWORK = ITAUP + N 867* 868* Bidiagonalize R in A 869* (CWorkspace: need 3*N, prefer 2*N+2*N*NB) 870* (RWorkspace: need N) 871* 872 CALL ZGEBRD( N, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ), 873 $ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1, 874 $ IERR ) 875 IRU = IE + N 876 IRVT = IRU + N*N 877 NRWORK = IRVT + N*N 878* 879* Perform bidiagonal SVD, computing left singular vectors 880* of bidiagonal matrix in RWORK(IRU) and computing right 881* singular vectors of bidiagonal matrix in RWORK(IRVT) 882* (CWorkspace: need 0) 883* (RWorkspace: need BDSPAC) 884* 885 CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ), 886 $ N, RWORK( IRVT ), N, DUM, IDUM, 887 $ RWORK( NRWORK ), IWORK, INFO ) 888* 889* Copy real matrix RWORK(IRU) to complex matrix WORK(IU) 890* Overwrite WORK(IU) by left singular vectors of R 891* (CWorkspace: need N*N+3*N, prefer N*N+2*N+N*NB) 892* (RWorkspace: 0) 893* 894 CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, WORK( IU ), 895 $ LDWRKU ) 896 CALL ZUNMBR( 'Q', 'L', 'N', N, N, N, A, LDA, 897 $ WORK( ITAUQ ), WORK( IU ), LDWRKU, 898 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 899* 900* Copy real matrix RWORK(IRVT) to complex matrix VT 901* Overwrite VT by right singular vectors of R 902* (CWorkspace: need 3*N, prefer 2*N+N*NB) 903* (RWorkspace: 0) 904* 905 CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT ) 906 CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA, 907 $ WORK( ITAUP ), VT, LDVT, WORK( NWORK ), 908 $ LWORK-NWORK+1, IERR ) 909* 910* Multiply Q in U by left singular vectors of R in 911* WORK(IU), storing result in A 912* (CWorkspace: need N*N) 913* (RWorkspace: 0) 914* 915 CALL ZGEMM( 'N', 'N', M, N, N, CONE, U, LDU, WORK( IU ), 916 $ LDWRKU, CZERO, A, LDA ) 917* 918* Copy left singular vectors of A from A to U 919* 920 CALL ZLACPY( 'F', M, N, A, LDA, U, LDU ) 921* 922 END IF 923* 924 ELSE IF( M.GE.MNTHR2 ) THEN 925* 926* MNTHR2 <= M < MNTHR1 927* 928* Path 5 (M much larger than N, but not as much as MNTHR1) 929* Reduce to bidiagonal form without QR decomposition, use 930* ZUNGBR and matrix multiplication to compute singular vectors 931* 932 IE = 1 933 NRWORK = IE + N 934 ITAUQ = 1 935 ITAUP = ITAUQ + N 936 NWORK = ITAUP + N 937* 938* Bidiagonalize A 939* (CWorkspace: need 2*N+M, prefer 2*N+(M+N)*NB) 940* (RWorkspace: need N) 941* 942 CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ), 943 $ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1, 944 $ IERR ) 945 IF( WNTQN ) THEN 946* 947* Compute singular values only 948* (Cworkspace: 0) 949* (Rworkspace: need BDSPAN) 950* 951 CALL DBDSDC( 'U', 'N', N, S, RWORK( IE ), DUM, 1, DUM, 1, 952 $ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO ) 953 ELSE IF( WNTQO ) THEN 954 IU = NWORK 955 IRU = NRWORK 956 IRVT = IRU + N*N 957 NRWORK = IRVT + N*N 958* 959* Copy A to VT, generate P**H 960* (Cworkspace: need 2*N, prefer N+N*NB) 961* (Rworkspace: 0) 962* 963 CALL ZLACPY( 'U', N, N, A, LDA, VT, LDVT ) 964 CALL ZUNGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), 965 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 966* 967* Generate Q in A 968* (CWorkspace: need 2*N, prefer N+N*NB) 969* (RWorkspace: 0) 970* 971 CALL ZUNGBR( 'Q', M, N, N, A, LDA, WORK( ITAUQ ), 972 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 973* 974 IF( LWORK.GE.M*N+3*N ) THEN 975* 976* WORK( IU ) is M by N 977* 978 LDWRKU = M 979 ELSE 980* 981* WORK(IU) is LDWRKU by N 982* 983 LDWRKU = ( LWORK-3*N ) / N 984 END IF 985 NWORK = IU + LDWRKU*N 986* 987* Perform bidiagonal SVD, computing left singular vectors 988* of bidiagonal matrix in RWORK(IRU) and computing right 989* singular vectors of bidiagonal matrix in RWORK(IRVT) 990* (CWorkspace: need 0) 991* (RWorkspace: need BDSPAC) 992* 993 CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ), 994 $ N, RWORK( IRVT ), N, DUM, IDUM, 995 $ RWORK( NRWORK ), IWORK, INFO ) 996* 997* Multiply real matrix RWORK(IRVT) by P**H in VT, 998* storing the result in WORK(IU), copying to VT 999* (Cworkspace: need 0) 1000* (Rworkspace: need 3*N*N) 1001* 1002 CALL ZLARCM( N, N, RWORK( IRVT ), N, VT, LDVT, 1003 $ WORK( IU ), LDWRKU, RWORK( NRWORK ) ) 1004 CALL ZLACPY( 'F', N, N, WORK( IU ), LDWRKU, VT, LDVT ) 1005* 1006* Multiply Q in A by real matrix RWORK(IRU), storing the 1007* result in WORK(IU), copying to A 1008* (CWorkspace: need N*N, prefer M*N) 1009* (Rworkspace: need 3*N*N, prefer N*N+2*M*N) 1010* 1011 NRWORK = IRVT 1012 DO 20 I = 1, M, LDWRKU 1013 CHUNK = MIN( M-I+1, LDWRKU ) 1014 CALL ZLACRM( CHUNK, N, A( I, 1 ), LDA, RWORK( IRU ), 1015 $ N, WORK( IU ), LDWRKU, RWORK( NRWORK ) ) 1016 CALL ZLACPY( 'F', CHUNK, N, WORK( IU ), LDWRKU, 1017 $ A( I, 1 ), LDA ) 1018 20 CONTINUE 1019* 1020 ELSE IF( WNTQS ) THEN 1021* 1022* Copy A to VT, generate P**H 1023* (Cworkspace: need 2*N, prefer N+N*NB) 1024* (Rworkspace: 0) 1025* 1026 CALL ZLACPY( 'U', N, N, A, LDA, VT, LDVT ) 1027 CALL ZUNGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), 1028 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1029* 1030* Copy A to U, generate Q 1031* (Cworkspace: need 2*N, prefer N+N*NB) 1032* (Rworkspace: 0) 1033* 1034 CALL ZLACPY( 'L', M, N, A, LDA, U, LDU ) 1035 CALL ZUNGBR( 'Q', M, N, N, U, LDU, WORK( ITAUQ ), 1036 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1037* 1038* Perform bidiagonal SVD, computing left singular vectors 1039* of bidiagonal matrix in RWORK(IRU) and computing right 1040* singular vectors of bidiagonal matrix in RWORK(IRVT) 1041* (CWorkspace: need 0) 1042* (RWorkspace: need BDSPAC) 1043* 1044 IRU = NRWORK 1045 IRVT = IRU + N*N 1046 NRWORK = IRVT + N*N 1047 CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ), 1048 $ N, RWORK( IRVT ), N, DUM, IDUM, 1049 $ RWORK( NRWORK ), IWORK, INFO ) 1050* 1051* Multiply real matrix RWORK(IRVT) by P**H in VT, 1052* storing the result in A, copying to VT 1053* (Cworkspace: need 0) 1054* (Rworkspace: need 3*N*N) 1055* 1056 CALL ZLARCM( N, N, RWORK( IRVT ), N, VT, LDVT, A, LDA, 1057 $ RWORK( NRWORK ) ) 1058 CALL ZLACPY( 'F', N, N, A, LDA, VT, LDVT ) 1059* 1060* Multiply Q in U by real matrix RWORK(IRU), storing the 1061* result in A, copying to U 1062* (CWorkspace: need 0) 1063* (Rworkspace: need N*N+2*M*N) 1064* 1065 NRWORK = IRVT 1066 CALL ZLACRM( M, N, U, LDU, RWORK( IRU ), N, A, LDA, 1067 $ RWORK( NRWORK ) ) 1068 CALL ZLACPY( 'F', M, N, A, LDA, U, LDU ) 1069 ELSE 1070* 1071* Copy A to VT, generate P**H 1072* (Cworkspace: need 2*N, prefer N+N*NB) 1073* (Rworkspace: 0) 1074* 1075 CALL ZLACPY( 'U', N, N, A, LDA, VT, LDVT ) 1076 CALL ZUNGBR( 'P', N, N, N, VT, LDVT, WORK( ITAUP ), 1077 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1078* 1079* Copy A to U, generate Q 1080* (Cworkspace: need 2*N, prefer N+N*NB) 1081* (Rworkspace: 0) 1082* 1083 CALL ZLACPY( 'L', M, N, A, LDA, U, LDU ) 1084 CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ), 1085 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1086* 1087* Perform bidiagonal SVD, computing left singular vectors 1088* of bidiagonal matrix in RWORK(IRU) and computing right 1089* singular vectors of bidiagonal matrix in RWORK(IRVT) 1090* (CWorkspace: need 0) 1091* (RWorkspace: need BDSPAC) 1092* 1093 IRU = NRWORK 1094 IRVT = IRU + N*N 1095 NRWORK = IRVT + N*N 1096 CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ), 1097 $ N, RWORK( IRVT ), N, DUM, IDUM, 1098 $ RWORK( NRWORK ), IWORK, INFO ) 1099* 1100* Multiply real matrix RWORK(IRVT) by P**H in VT, 1101* storing the result in A, copying to VT 1102* (Cworkspace: need 0) 1103* (Rworkspace: need 3*N*N) 1104* 1105 CALL ZLARCM( N, N, RWORK( IRVT ), N, VT, LDVT, A, LDA, 1106 $ RWORK( NRWORK ) ) 1107 CALL ZLACPY( 'F', N, N, A, LDA, VT, LDVT ) 1108* 1109* Multiply Q in U by real matrix RWORK(IRU), storing the 1110* result in A, copying to U 1111* (CWorkspace: 0) 1112* (Rworkspace: need 3*N*N) 1113* 1114 NRWORK = IRVT 1115 CALL ZLACRM( M, N, U, LDU, RWORK( IRU ), N, A, LDA, 1116 $ RWORK( NRWORK ) ) 1117 CALL ZLACPY( 'F', M, N, A, LDA, U, LDU ) 1118 END IF 1119* 1120 ELSE 1121* 1122* M .LT. MNTHR2 1123* 1124* Path 6 (M at least N, but not much larger) 1125* Reduce to bidiagonal form without QR decomposition 1126* Use ZUNMBR to compute singular vectors 1127* 1128 IE = 1 1129 NRWORK = IE + N 1130 ITAUQ = 1 1131 ITAUP = ITAUQ + N 1132 NWORK = ITAUP + N 1133* 1134* Bidiagonalize A 1135* (CWorkspace: need 2*N+M, prefer 2*N+(M+N)*NB) 1136* (RWorkspace: need N) 1137* 1138 CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ), 1139 $ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1, 1140 $ IERR ) 1141 IF( WNTQN ) THEN 1142* 1143* Compute singular values only 1144* (Cworkspace: 0) 1145* (Rworkspace: need BDSPAN) 1146* 1147 CALL DBDSDC( 'U', 'N', N, S, RWORK( IE ), DUM, 1, DUM, 1, 1148 $ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO ) 1149 ELSE IF( WNTQO ) THEN 1150 IU = NWORK 1151 IRU = NRWORK 1152 IRVT = IRU + N*N 1153 NRWORK = IRVT + N*N 1154 IF( LWORK.GE.M*N+3*N ) THEN 1155* 1156* WORK( IU ) is M by N 1157* 1158 LDWRKU = M 1159 ELSE 1160* 1161* WORK( IU ) is LDWRKU by N 1162* 1163 LDWRKU = ( LWORK-3*N ) / N 1164 END IF 1165 NWORK = IU + LDWRKU*N 1166* 1167* Perform bidiagonal SVD, computing left singular vectors 1168* of bidiagonal matrix in RWORK(IRU) and computing right 1169* singular vectors of bidiagonal matrix in RWORK(IRVT) 1170* (CWorkspace: need 0) 1171* (RWorkspace: need BDSPAC) 1172* 1173 CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ), 1174 $ N, RWORK( IRVT ), N, DUM, IDUM, 1175 $ RWORK( NRWORK ), IWORK, INFO ) 1176* 1177* Copy real matrix RWORK(IRVT) to complex matrix VT 1178* Overwrite VT by right singular vectors of A 1179* (Cworkspace: need 2*N, prefer N+N*NB) 1180* (Rworkspace: need 0) 1181* 1182 CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT ) 1183 CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA, 1184 $ WORK( ITAUP ), VT, LDVT, WORK( NWORK ), 1185 $ LWORK-NWORK+1, IERR ) 1186* 1187 IF( LWORK.GE.M*N+3*N ) THEN 1188* 1189* Copy real matrix RWORK(IRU) to complex matrix WORK(IU) 1190* Overwrite WORK(IU) by left singular vectors of A, copying 1191* to A 1192* (Cworkspace: need M*N+2*N, prefer M*N+N+N*NB) 1193* (Rworkspace: need 0) 1194* 1195 CALL ZLASET( 'F', M, N, CZERO, CZERO, WORK( IU ), 1196 $ LDWRKU ) 1197 CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, WORK( IU ), 1198 $ LDWRKU ) 1199 CALL ZUNMBR( 'Q', 'L', 'N', M, N, N, A, LDA, 1200 $ WORK( ITAUQ ), WORK( IU ), LDWRKU, 1201 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1202 CALL ZLACPY( 'F', M, N, WORK( IU ), LDWRKU, A, LDA ) 1203 ELSE 1204* 1205* Generate Q in A 1206* (Cworkspace: need 2*N, prefer N+N*NB) 1207* (Rworkspace: need 0) 1208* 1209 CALL ZUNGBR( 'Q', M, N, N, A, LDA, WORK( ITAUQ ), 1210 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1211* 1212* Multiply Q in A by real matrix RWORK(IRU), storing the 1213* result in WORK(IU), copying to A 1214* (CWorkspace: need N*N, prefer M*N) 1215* (Rworkspace: need 3*N*N, prefer N*N+2*M*N) 1216* 1217 NRWORK = IRVT 1218 DO 30 I = 1, M, LDWRKU 1219 CHUNK = MIN( M-I+1, LDWRKU ) 1220 CALL ZLACRM( CHUNK, N, A( I, 1 ), LDA, 1221 $ RWORK( IRU ), N, WORK( IU ), LDWRKU, 1222 $ RWORK( NRWORK ) ) 1223 CALL ZLACPY( 'F', CHUNK, N, WORK( IU ), LDWRKU, 1224 $ A( I, 1 ), LDA ) 1225 30 CONTINUE 1226 END IF 1227* 1228 ELSE IF( WNTQS ) THEN 1229* 1230* Perform bidiagonal SVD, computing left singular vectors 1231* of bidiagonal matrix in RWORK(IRU) and computing right 1232* singular vectors of bidiagonal matrix in RWORK(IRVT) 1233* (CWorkspace: need 0) 1234* (RWorkspace: need BDSPAC) 1235* 1236 IRU = NRWORK 1237 IRVT = IRU + N*N 1238 NRWORK = IRVT + N*N 1239 CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ), 1240 $ N, RWORK( IRVT ), N, DUM, IDUM, 1241 $ RWORK( NRWORK ), IWORK, INFO ) 1242* 1243* Copy real matrix RWORK(IRU) to complex matrix U 1244* Overwrite U by left singular vectors of A 1245* (CWorkspace: need 3*N, prefer 2*N+N*NB) 1246* (RWorkspace: 0) 1247* 1248 CALL ZLASET( 'F', M, N, CZERO, CZERO, U, LDU ) 1249 CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, U, LDU ) 1250 CALL ZUNMBR( 'Q', 'L', 'N', M, N, N, A, LDA, 1251 $ WORK( ITAUQ ), U, LDU, WORK( NWORK ), 1252 $ LWORK-NWORK+1, IERR ) 1253* 1254* Copy real matrix RWORK(IRVT) to complex matrix VT 1255* Overwrite VT by right singular vectors of A 1256* (CWorkspace: need 3*N, prefer 2*N+N*NB) 1257* (RWorkspace: 0) 1258* 1259 CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT ) 1260 CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA, 1261 $ WORK( ITAUP ), VT, LDVT, WORK( NWORK ), 1262 $ LWORK-NWORK+1, IERR ) 1263 ELSE 1264* 1265* Perform bidiagonal SVD, computing left singular vectors 1266* of bidiagonal matrix in RWORK(IRU) and computing right 1267* singular vectors of bidiagonal matrix in RWORK(IRVT) 1268* (CWorkspace: need 0) 1269* (RWorkspace: need BDSPAC) 1270* 1271 IRU = NRWORK 1272 IRVT = IRU + N*N 1273 NRWORK = IRVT + N*N 1274 CALL DBDSDC( 'U', 'I', N, S, RWORK( IE ), RWORK( IRU ), 1275 $ N, RWORK( IRVT ), N, DUM, IDUM, 1276 $ RWORK( NRWORK ), IWORK, INFO ) 1277* 1278* Set the right corner of U to identity matrix 1279* 1280 CALL ZLASET( 'F', M, M, CZERO, CZERO, U, LDU ) 1281 IF( M.GT.N ) THEN 1282 CALL ZLASET( 'F', M-N, M-N, CZERO, CONE, 1283 $ U( N+1, N+1 ), LDU ) 1284 END IF 1285* 1286* Copy real matrix RWORK(IRU) to complex matrix U 1287* Overwrite U by left singular vectors of A 1288* (CWorkspace: need 2*N+M, prefer 2*N+M*NB) 1289* (RWorkspace: 0) 1290* 1291 CALL ZLACP2( 'F', N, N, RWORK( IRU ), N, U, LDU ) 1292 CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA, 1293 $ WORK( ITAUQ ), U, LDU, WORK( NWORK ), 1294 $ LWORK-NWORK+1, IERR ) 1295* 1296* Copy real matrix RWORK(IRVT) to complex matrix VT 1297* Overwrite VT by right singular vectors of A 1298* (CWorkspace: need 3*N, prefer 2*N+N*NB) 1299* (RWorkspace: 0) 1300* 1301 CALL ZLACP2( 'F', N, N, RWORK( IRVT ), N, VT, LDVT ) 1302 CALL ZUNMBR( 'P', 'R', 'C', N, N, N, A, LDA, 1303 $ WORK( ITAUP ), VT, LDVT, WORK( NWORK ), 1304 $ LWORK-NWORK+1, IERR ) 1305 END IF 1306* 1307 END IF 1308* 1309 ELSE 1310* 1311* A has more columns than rows. If A has sufficiently more 1312* columns than rows, first reduce using the LQ decomposition (if 1313* sufficient workspace available) 1314* 1315 IF( N.GE.MNTHR1 ) THEN 1316* 1317 IF( WNTQN ) THEN 1318* 1319* Path 1t (N much larger than M, JOBZ='N') 1320* No singular vectors to be computed 1321* 1322 ITAU = 1 1323 NWORK = ITAU + M 1324* 1325* Compute A=L*Q 1326* (CWorkspace: need 2*M, prefer M+M*NB) 1327* (RWorkspace: 0) 1328* 1329 CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ), 1330 $ LWORK-NWORK+1, IERR ) 1331* 1332* Zero out above L 1333* 1334 CALL ZLASET( 'U', M-1, M-1, CZERO, CZERO, A( 1, 2 ), 1335 $ LDA ) 1336 IE = 1 1337 ITAUQ = 1 1338 ITAUP = ITAUQ + M 1339 NWORK = ITAUP + M 1340* 1341* Bidiagonalize L in A 1342* (CWorkspace: need 3*M, prefer 2*M+2*M*NB) 1343* (RWorkspace: need M) 1344* 1345 CALL ZGEBRD( M, M, A, LDA, S, RWORK( IE ), WORK( ITAUQ ), 1346 $ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1, 1347 $ IERR ) 1348 NRWORK = IE + M 1349* 1350* Perform bidiagonal SVD, compute singular values only 1351* (CWorkspace: 0) 1352* (RWorkspace: need BDSPAN) 1353* 1354 CALL DBDSDC( 'U', 'N', M, S, RWORK( IE ), DUM, 1, DUM, 1, 1355 $ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO ) 1356* 1357 ELSE IF( WNTQO ) THEN 1358* 1359* Path 2t (N much larger than M, JOBZ='O') 1360* M right singular vectors to be overwritten on A and 1361* M left singular vectors to be computed in U 1362* 1363 IVT = 1 1364 LDWKVT = M 1365* 1366* WORK(IVT) is M by M 1367* 1368 IL = IVT + LDWKVT*M 1369 IF( LWORK.GE.M*N+M*M+3*M ) THEN 1370* 1371* WORK(IL) M by N 1372* 1373 LDWRKL = M 1374 CHUNK = N 1375 ELSE 1376* 1377* WORK(IL) is M by CHUNK 1378* 1379 LDWRKL = M 1380 CHUNK = ( LWORK-M*M-3*M ) / M 1381 END IF 1382 ITAU = IL + LDWRKL*CHUNK 1383 NWORK = ITAU + M 1384* 1385* Compute A=L*Q 1386* (CWorkspace: need 2*M, prefer M+M*NB) 1387* (RWorkspace: 0) 1388* 1389 CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ), 1390 $ LWORK-NWORK+1, IERR ) 1391* 1392* Copy L to WORK(IL), zeroing about above it 1393* 1394 CALL ZLACPY( 'L', M, M, A, LDA, WORK( IL ), LDWRKL ) 1395 CALL ZLASET( 'U', M-1, M-1, CZERO, CZERO, 1396 $ WORK( IL+LDWRKL ), LDWRKL ) 1397* 1398* Generate Q in A 1399* (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) 1400* (RWorkspace: 0) 1401* 1402 CALL ZUNGLQ( M, N, M, A, LDA, WORK( ITAU ), 1403 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1404 IE = 1 1405 ITAUQ = ITAU 1406 ITAUP = ITAUQ + M 1407 NWORK = ITAUP + M 1408* 1409* Bidiagonalize L in WORK(IL) 1410* (CWorkspace: need M*M+3*M, prefer M*M+2*M+2*M*NB) 1411* (RWorkspace: need M) 1412* 1413 CALL ZGEBRD( M, M, WORK( IL ), LDWRKL, S, RWORK( IE ), 1414 $ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ), 1415 $ LWORK-NWORK+1, IERR ) 1416* 1417* Perform bidiagonal SVD, computing left singular vectors 1418* of bidiagonal matrix in RWORK(IRU) and computing right 1419* singular vectors of bidiagonal matrix in RWORK(IRVT) 1420* (CWorkspace: need 0) 1421* (RWorkspace: need BDSPAC) 1422* 1423 IRU = IE + M 1424 IRVT = IRU + M*M 1425 NRWORK = IRVT + M*M 1426 CALL DBDSDC( 'U', 'I', M, S, RWORK( IE ), RWORK( IRU ), 1427 $ M, RWORK( IRVT ), M, DUM, IDUM, 1428 $ RWORK( NRWORK ), IWORK, INFO ) 1429* 1430* Copy real matrix RWORK(IRU) to complex matrix WORK(IU) 1431* Overwrite WORK(IU) by the left singular vectors of L 1432* (CWorkspace: need N*N+3*N, prefer M*N+2*N+N*NB) 1433* (RWorkspace: 0) 1434* 1435 CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU ) 1436 CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, WORK( IL ), LDWRKL, 1437 $ WORK( ITAUQ ), U, LDU, WORK( NWORK ), 1438 $ LWORK-NWORK+1, IERR ) 1439* 1440* Copy real matrix RWORK(IRVT) to complex matrix WORK(IVT) 1441* Overwrite WORK(IVT) by the right singular vectors of L 1442* (CWorkspace: need N*N+3*N, prefer M*N+2*N+N*NB) 1443* (RWorkspace: 0) 1444* 1445 CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, WORK( IVT ), 1446 $ LDWKVT ) 1447 CALL ZUNMBR( 'P', 'R', 'C', M, M, M, WORK( IL ), LDWRKL, 1448 $ WORK( ITAUP ), WORK( IVT ), LDWKVT, 1449 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1450* 1451* Multiply right singular vectors of L in WORK(IL) by Q 1452* in A, storing result in WORK(IL) and copying to A 1453* (CWorkspace: need 2*M*M, prefer M*M+M*N)) 1454* (RWorkspace: 0) 1455* 1456 DO 40 I = 1, N, CHUNK 1457 BLK = MIN( N-I+1, CHUNK ) 1458 CALL ZGEMM( 'N', 'N', M, BLK, M, CONE, WORK( IVT ), M, 1459 $ A( 1, I ), LDA, CZERO, WORK( IL ), 1460 $ LDWRKL ) 1461 CALL ZLACPY( 'F', M, BLK, WORK( IL ), LDWRKL, 1462 $ A( 1, I ), LDA ) 1463 40 CONTINUE 1464* 1465 ELSE IF( WNTQS ) THEN 1466* 1467* Path 3t (N much larger than M, JOBZ='S') 1468* M right singular vectors to be computed in VT and 1469* M left singular vectors to be computed in U 1470* 1471 IL = 1 1472* 1473* WORK(IL) is M by M 1474* 1475 LDWRKL = M 1476 ITAU = IL + LDWRKL*M 1477 NWORK = ITAU + M 1478* 1479* Compute A=L*Q 1480* (CWorkspace: need 2*M, prefer M+M*NB) 1481* (RWorkspace: 0) 1482* 1483 CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ), 1484 $ LWORK-NWORK+1, IERR ) 1485* 1486* Copy L to WORK(IL), zeroing out above it 1487* 1488 CALL ZLACPY( 'L', M, M, A, LDA, WORK( IL ), LDWRKL ) 1489 CALL ZLASET( 'U', M-1, M-1, CZERO, CZERO, 1490 $ WORK( IL+LDWRKL ), LDWRKL ) 1491* 1492* Generate Q in A 1493* (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) 1494* (RWorkspace: 0) 1495* 1496 CALL ZUNGLQ( M, N, M, A, LDA, WORK( ITAU ), 1497 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1498 IE = 1 1499 ITAUQ = ITAU 1500 ITAUP = ITAUQ + M 1501 NWORK = ITAUP + M 1502* 1503* Bidiagonalize L in WORK(IL) 1504* (CWorkspace: need M*M+3*M, prefer M*M+2*M+2*M*NB) 1505* (RWorkspace: need M) 1506* 1507 CALL ZGEBRD( M, M, WORK( IL ), LDWRKL, S, RWORK( IE ), 1508 $ WORK( ITAUQ ), WORK( ITAUP ), WORK( NWORK ), 1509 $ LWORK-NWORK+1, IERR ) 1510* 1511* Perform bidiagonal SVD, computing left singular vectors 1512* of bidiagonal matrix in RWORK(IRU) and computing right 1513* singular vectors of bidiagonal matrix in RWORK(IRVT) 1514* (CWorkspace: need 0) 1515* (RWorkspace: need BDSPAC) 1516* 1517 IRU = IE + M 1518 IRVT = IRU + M*M 1519 NRWORK = IRVT + M*M 1520 CALL DBDSDC( 'U', 'I', M, S, RWORK( IE ), RWORK( IRU ), 1521 $ M, RWORK( IRVT ), M, DUM, IDUM, 1522 $ RWORK( NRWORK ), IWORK, INFO ) 1523* 1524* Copy real matrix RWORK(IRU) to complex matrix U 1525* Overwrite U by left singular vectors of L 1526* (CWorkspace: need M*M+3*M, prefer M*M+2*M+M*NB) 1527* (RWorkspace: 0) 1528* 1529 CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU ) 1530 CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, WORK( IL ), LDWRKL, 1531 $ WORK( ITAUQ ), U, LDU, WORK( NWORK ), 1532 $ LWORK-NWORK+1, IERR ) 1533* 1534* Copy real matrix RWORK(IRVT) to complex matrix VT 1535* Overwrite VT by left singular vectors of L 1536* (CWorkspace: need M*M+3*M, prefer M*M+2*M+M*NB) 1537* (RWorkspace: 0) 1538* 1539 CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, VT, LDVT ) 1540 CALL ZUNMBR( 'P', 'R', 'C', M, M, M, WORK( IL ), LDWRKL, 1541 $ WORK( ITAUP ), VT, LDVT, WORK( NWORK ), 1542 $ LWORK-NWORK+1, IERR ) 1543* 1544* Copy VT to WORK(IL), multiply right singular vectors of L 1545* in WORK(IL) by Q in A, storing result in VT 1546* (CWorkspace: need M*M) 1547* (RWorkspace: 0) 1548* 1549 CALL ZLACPY( 'F', M, M, VT, LDVT, WORK( IL ), LDWRKL ) 1550 CALL ZGEMM( 'N', 'N', M, N, M, CONE, WORK( IL ), LDWRKL, 1551 $ A, LDA, CZERO, VT, LDVT ) 1552* 1553 ELSE IF( WNTQA ) THEN 1554* 1555* Path 9t (N much larger than M, JOBZ='A') 1556* N right singular vectors to be computed in VT and 1557* M left singular vectors to be computed in U 1558* 1559 IVT = 1 1560* 1561* WORK(IVT) is M by M 1562* 1563 LDWKVT = M 1564 ITAU = IVT + LDWKVT*M 1565 NWORK = ITAU + M 1566* 1567* Compute A=L*Q, copying result to VT 1568* (CWorkspace: need 2*M, prefer M+M*NB) 1569* (RWorkspace: 0) 1570* 1571 CALL ZGELQF( M, N, A, LDA, WORK( ITAU ), WORK( NWORK ), 1572 $ LWORK-NWORK+1, IERR ) 1573 CALL ZLACPY( 'U', M, N, A, LDA, VT, LDVT ) 1574* 1575* Generate Q in VT 1576* (CWorkspace: need M+N, prefer M+N*NB) 1577* (RWorkspace: 0) 1578* 1579 CALL ZUNGLQ( N, N, M, VT, LDVT, WORK( ITAU ), 1580 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1581* 1582* Produce L in A, zeroing out above it 1583* 1584 CALL ZLASET( 'U', M-1, M-1, CZERO, CZERO, A( 1, 2 ), 1585 $ LDA ) 1586 IE = 1 1587 ITAUQ = ITAU 1588 ITAUP = ITAUQ + M 1589 NWORK = ITAUP + M 1590* 1591* Bidiagonalize L in A 1592* (CWorkspace: need M*M+3*M, prefer M*M+2*M+2*M*NB) 1593* (RWorkspace: need M) 1594* 1595 CALL ZGEBRD( M, M, A, LDA, S, RWORK( IE ), WORK( ITAUQ ), 1596 $ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1, 1597 $ IERR ) 1598* 1599* Perform bidiagonal SVD, computing left singular vectors 1600* of bidiagonal matrix in RWORK(IRU) and computing right 1601* singular vectors of bidiagonal matrix in RWORK(IRVT) 1602* (CWorkspace: need 0) 1603* (RWorkspace: need BDSPAC) 1604* 1605 IRU = IE + M 1606 IRVT = IRU + M*M 1607 NRWORK = IRVT + M*M 1608 CALL DBDSDC( 'U', 'I', M, S, RWORK( IE ), RWORK( IRU ), 1609 $ M, RWORK( IRVT ), M, DUM, IDUM, 1610 $ RWORK( NRWORK ), IWORK, INFO ) 1611* 1612* Copy real matrix RWORK(IRU) to complex matrix U 1613* Overwrite U by left singular vectors of L 1614* (CWorkspace: need 3*M, prefer 2*M+M*NB) 1615* (RWorkspace: 0) 1616* 1617 CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU ) 1618 CALL ZUNMBR( 'Q', 'L', 'N', M, M, M, A, LDA, 1619 $ WORK( ITAUQ ), U, LDU, WORK( NWORK ), 1620 $ LWORK-NWORK+1, IERR ) 1621* 1622* Copy real matrix RWORK(IRVT) to complex matrix WORK(IVT) 1623* Overwrite WORK(IVT) by right singular vectors of L 1624* (CWorkspace: need M*M+3*M, prefer M*M+2*M+M*NB) 1625* (RWorkspace: 0) 1626* 1627 CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, WORK( IVT ), 1628 $ LDWKVT ) 1629 CALL ZUNMBR( 'P', 'R', 'C', M, M, M, A, LDA, 1630 $ WORK( ITAUP ), WORK( IVT ), LDWKVT, 1631 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1632* 1633* Multiply right singular vectors of L in WORK(IVT) by 1634* Q in VT, storing result in A 1635* (CWorkspace: need M*M) 1636* (RWorkspace: 0) 1637* 1638 CALL ZGEMM( 'N', 'N', M, N, M, CONE, WORK( IVT ), LDWKVT, 1639 $ VT, LDVT, CZERO, A, LDA ) 1640* 1641* Copy right singular vectors of A from A to VT 1642* 1643 CALL ZLACPY( 'F', M, N, A, LDA, VT, LDVT ) 1644* 1645 END IF 1646* 1647 ELSE IF( N.GE.MNTHR2 ) THEN 1648* 1649* MNTHR2 <= N < MNTHR1 1650* 1651* Path 5t (N much larger than M, but not as much as MNTHR1) 1652* Reduce to bidiagonal form without QR decomposition, use 1653* ZUNGBR and matrix multiplication to compute singular vectors 1654* 1655* 1656 IE = 1 1657 NRWORK = IE + M 1658 ITAUQ = 1 1659 ITAUP = ITAUQ + M 1660 NWORK = ITAUP + M 1661* 1662* Bidiagonalize A 1663* (CWorkspace: need 2*M+N, prefer 2*M+(M+N)*NB) 1664* (RWorkspace: M) 1665* 1666 CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ), 1667 $ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1, 1668 $ IERR ) 1669* 1670 IF( WNTQN ) THEN 1671* 1672* Compute singular values only 1673* (Cworkspace: 0) 1674* (Rworkspace: need BDSPAN) 1675* 1676 CALL DBDSDC( 'L', 'N', M, S, RWORK( IE ), DUM, 1, DUM, 1, 1677 $ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO ) 1678 ELSE IF( WNTQO ) THEN 1679 IRVT = NRWORK 1680 IRU = IRVT + M*M 1681 NRWORK = IRU + M*M 1682 IVT = NWORK 1683* 1684* Copy A to U, generate Q 1685* (Cworkspace: need 2*M, prefer M+M*NB) 1686* (Rworkspace: 0) 1687* 1688 CALL ZLACPY( 'L', M, M, A, LDA, U, LDU ) 1689 CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ), 1690 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1691* 1692* Generate P**H in A 1693* (Cworkspace: need 2*M, prefer M+M*NB) 1694* (Rworkspace: 0) 1695* 1696 CALL ZUNGBR( 'P', M, N, M, A, LDA, WORK( ITAUP ), 1697 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1698* 1699 LDWKVT = M 1700 IF( LWORK.GE.M*N+3*M ) THEN 1701* 1702* WORK( IVT ) is M by N 1703* 1704 NWORK = IVT + LDWKVT*N 1705 CHUNK = N 1706 ELSE 1707* 1708* WORK( IVT ) is M by CHUNK 1709* 1710 CHUNK = ( LWORK-3*M ) / M 1711 NWORK = IVT + LDWKVT*CHUNK 1712 END IF 1713* 1714* Perform bidiagonal SVD, computing left singular vectors 1715* of bidiagonal matrix in RWORK(IRU) and computing right 1716* singular vectors of bidiagonal matrix in RWORK(IRVT) 1717* (CWorkspace: need 0) 1718* (RWorkspace: need BDSPAC) 1719* 1720 CALL DBDSDC( 'L', 'I', M, S, RWORK( IE ), RWORK( IRU ), 1721 $ M, RWORK( IRVT ), M, DUM, IDUM, 1722 $ RWORK( NRWORK ), IWORK, INFO ) 1723* 1724* Multiply Q in U by real matrix RWORK(IRVT) 1725* storing the result in WORK(IVT), copying to U 1726* (Cworkspace: need 0) 1727* (Rworkspace: need 2*M*M) 1728* 1729 CALL ZLACRM( M, M, U, LDU, RWORK( IRU ), M, WORK( IVT ), 1730 $ LDWKVT, RWORK( NRWORK ) ) 1731 CALL ZLACPY( 'F', M, M, WORK( IVT ), LDWKVT, U, LDU ) 1732* 1733* Multiply RWORK(IRVT) by P**H in A, storing the 1734* result in WORK(IVT), copying to A 1735* (CWorkspace: need M*M, prefer M*N) 1736* (Rworkspace: need 2*M*M, prefer 2*M*N) 1737* 1738 NRWORK = IRU 1739 DO 50 I = 1, N, CHUNK 1740 BLK = MIN( N-I+1, CHUNK ) 1741 CALL ZLARCM( M, BLK, RWORK( IRVT ), M, A( 1, I ), LDA, 1742 $ WORK( IVT ), LDWKVT, RWORK( NRWORK ) ) 1743 CALL ZLACPY( 'F', M, BLK, WORK( IVT ), LDWKVT, 1744 $ A( 1, I ), LDA ) 1745 50 CONTINUE 1746 ELSE IF( WNTQS ) THEN 1747* 1748* Copy A to U, generate Q 1749* (Cworkspace: need 2*M, prefer M+M*NB) 1750* (Rworkspace: 0) 1751* 1752 CALL ZLACPY( 'L', M, M, A, LDA, U, LDU ) 1753 CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ), 1754 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1755* 1756* Copy A to VT, generate P**H 1757* (Cworkspace: need 2*M, prefer M+M*NB) 1758* (Rworkspace: 0) 1759* 1760 CALL ZLACPY( 'U', M, N, A, LDA, VT, LDVT ) 1761 CALL ZUNGBR( 'P', M, N, M, VT, LDVT, WORK( ITAUP ), 1762 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1763* 1764* Perform bidiagonal SVD, computing left singular vectors 1765* of bidiagonal matrix in RWORK(IRU) and computing right 1766* singular vectors of bidiagonal matrix in RWORK(IRVT) 1767* (CWorkspace: need 0) 1768* (RWorkspace: need BDSPAC) 1769* 1770 IRVT = NRWORK 1771 IRU = IRVT + M*M 1772 NRWORK = IRU + M*M 1773 CALL DBDSDC( 'L', 'I', M, S, RWORK( IE ), RWORK( IRU ), 1774 $ M, RWORK( IRVT ), M, DUM, IDUM, 1775 $ RWORK( NRWORK ), IWORK, INFO ) 1776* 1777* Multiply Q in U by real matrix RWORK(IRU), storing the 1778* result in A, copying to U 1779* (CWorkspace: need 0) 1780* (Rworkspace: need 3*M*M) 1781* 1782 CALL ZLACRM( M, M, U, LDU, RWORK( IRU ), M, A, LDA, 1783 $ RWORK( NRWORK ) ) 1784 CALL ZLACPY( 'F', M, M, A, LDA, U, LDU ) 1785* 1786* Multiply real matrix RWORK(IRVT) by P**H in VT, 1787* storing the result in A, copying to VT 1788* (Cworkspace: need 0) 1789* (Rworkspace: need M*M+2*M*N) 1790* 1791 NRWORK = IRU 1792 CALL ZLARCM( M, N, RWORK( IRVT ), M, VT, LDVT, A, LDA, 1793 $ RWORK( NRWORK ) ) 1794 CALL ZLACPY( 'F', M, N, A, LDA, VT, LDVT ) 1795 ELSE 1796* 1797* Copy A to U, generate Q 1798* (Cworkspace: need 2*M, prefer M+M*NB) 1799* (Rworkspace: 0) 1800* 1801 CALL ZLACPY( 'L', M, M, A, LDA, U, LDU ) 1802 CALL ZUNGBR( 'Q', M, M, N, U, LDU, WORK( ITAUQ ), 1803 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1804* 1805* Copy A to VT, generate P**H 1806* (Cworkspace: need 2*M, prefer M+M*NB) 1807* (Rworkspace: 0) 1808* 1809 CALL ZLACPY( 'U', M, N, A, LDA, VT, LDVT ) 1810 CALL ZUNGBR( 'P', N, N, M, VT, LDVT, WORK( ITAUP ), 1811 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1812* 1813* Perform bidiagonal SVD, computing left singular vectors 1814* of bidiagonal matrix in RWORK(IRU) and computing right 1815* singular vectors of bidiagonal matrix in RWORK(IRVT) 1816* (CWorkspace: need 0) 1817* (RWorkspace: need BDSPAC) 1818* 1819 IRVT = NRWORK 1820 IRU = IRVT + M*M 1821 NRWORK = IRU + M*M 1822 CALL DBDSDC( 'L', 'I', M, S, RWORK( IE ), RWORK( IRU ), 1823 $ M, RWORK( IRVT ), M, DUM, IDUM, 1824 $ RWORK( NRWORK ), IWORK, INFO ) 1825* 1826* Multiply Q in U by real matrix RWORK(IRU), storing the 1827* result in A, copying to U 1828* (CWorkspace: need 0) 1829* (Rworkspace: need 3*M*M) 1830* 1831 CALL ZLACRM( M, M, U, LDU, RWORK( IRU ), M, A, LDA, 1832 $ RWORK( NRWORK ) ) 1833 CALL ZLACPY( 'F', M, M, A, LDA, U, LDU ) 1834* 1835* Multiply real matrix RWORK(IRVT) by P**H in VT, 1836* storing the result in A, copying to VT 1837* (Cworkspace: need 0) 1838* (Rworkspace: need M*M+2*M*N) 1839* 1840 CALL ZLARCM( M, N, RWORK( IRVT ), M, VT, LDVT, A, LDA, 1841 $ RWORK( NRWORK ) ) 1842 CALL ZLACPY( 'F', M, N, A, LDA, VT, LDVT ) 1843 END IF 1844* 1845 ELSE 1846* 1847* N .LT. MNTHR2 1848* 1849* Path 6t (N greater than M, but not much larger) 1850* Reduce to bidiagonal form without LQ decomposition 1851* Use ZUNMBR to compute singular vectors 1852* 1853 IE = 1 1854 NRWORK = IE + M 1855 ITAUQ = 1 1856 ITAUP = ITAUQ + M 1857 NWORK = ITAUP + M 1858* 1859* Bidiagonalize A 1860* (CWorkspace: need 2*M+N, prefer 2*M+(M+N)*NB) 1861* (RWorkspace: M) 1862* 1863 CALL ZGEBRD( M, N, A, LDA, S, RWORK( IE ), WORK( ITAUQ ), 1864 $ WORK( ITAUP ), WORK( NWORK ), LWORK-NWORK+1, 1865 $ IERR ) 1866 IF( WNTQN ) THEN 1867* 1868* Compute singular values only 1869* (Cworkspace: 0) 1870* (Rworkspace: need BDSPAN) 1871* 1872 CALL DBDSDC( 'L', 'N', M, S, RWORK( IE ), DUM, 1, DUM, 1, 1873 $ DUM, IDUM, RWORK( NRWORK ), IWORK, INFO ) 1874 ELSE IF( WNTQO ) THEN 1875 LDWKVT = M 1876 IVT = NWORK 1877 IF( LWORK.GE.M*N+3*M ) THEN 1878* 1879* WORK( IVT ) is M by N 1880* 1881 CALL ZLASET( 'F', M, N, CZERO, CZERO, WORK( IVT ), 1882 $ LDWKVT ) 1883 NWORK = IVT + LDWKVT*N 1884 ELSE 1885* 1886* WORK( IVT ) is M by CHUNK 1887* 1888 CHUNK = ( LWORK-3*M ) / M 1889 NWORK = IVT + LDWKVT*CHUNK 1890 END IF 1891* 1892* Perform bidiagonal SVD, computing left singular vectors 1893* of bidiagonal matrix in RWORK(IRU) and computing right 1894* singular vectors of bidiagonal matrix in RWORK(IRVT) 1895* (CWorkspace: need 0) 1896* (RWorkspace: need BDSPAC) 1897* 1898 IRVT = NRWORK 1899 IRU = IRVT + M*M 1900 NRWORK = IRU + M*M 1901 CALL DBDSDC( 'L', 'I', M, S, RWORK( IE ), RWORK( IRU ), 1902 $ M, RWORK( IRVT ), M, DUM, IDUM, 1903 $ RWORK( NRWORK ), IWORK, INFO ) 1904* 1905* Copy real matrix RWORK(IRU) to complex matrix U 1906* Overwrite U by left singular vectors of A 1907* (Cworkspace: need 2*M, prefer M+M*NB) 1908* (Rworkspace: need 0) 1909* 1910 CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU ) 1911 CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA, 1912 $ WORK( ITAUQ ), U, LDU, WORK( NWORK ), 1913 $ LWORK-NWORK+1, IERR ) 1914* 1915 IF( LWORK.GE.M*N+3*M ) THEN 1916* 1917* Copy real matrix RWORK(IRVT) to complex matrix WORK(IVT) 1918* Overwrite WORK(IVT) by right singular vectors of A, 1919* copying to A 1920* (Cworkspace: need M*N+2*M, prefer M*N+M+M*NB) 1921* (Rworkspace: need 0) 1922* 1923 CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, WORK( IVT ), 1924 $ LDWKVT ) 1925 CALL ZUNMBR( 'P', 'R', 'C', M, N, M, A, LDA, 1926 $ WORK( ITAUP ), WORK( IVT ), LDWKVT, 1927 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1928 CALL ZLACPY( 'F', M, N, WORK( IVT ), LDWKVT, A, LDA ) 1929 ELSE 1930* 1931* Generate P**H in A 1932* (Cworkspace: need 2*M, prefer M+M*NB) 1933* (Rworkspace: need 0) 1934* 1935 CALL ZUNGBR( 'P', M, N, M, A, LDA, WORK( ITAUP ), 1936 $ WORK( NWORK ), LWORK-NWORK+1, IERR ) 1937* 1938* Multiply Q in A by real matrix RWORK(IRU), storing the 1939* result in WORK(IU), copying to A 1940* (CWorkspace: need M*M, prefer M*N) 1941* (Rworkspace: need 3*M*M, prefer M*M+2*M*N) 1942* 1943 NRWORK = IRU 1944 DO 60 I = 1, N, CHUNK 1945 BLK = MIN( N-I+1, CHUNK ) 1946 CALL ZLARCM( M, BLK, RWORK( IRVT ), M, A( 1, I ), 1947 $ LDA, WORK( IVT ), LDWKVT, 1948 $ RWORK( NRWORK ) ) 1949 CALL ZLACPY( 'F', M, BLK, WORK( IVT ), LDWKVT, 1950 $ A( 1, I ), LDA ) 1951 60 CONTINUE 1952 END IF 1953 ELSE IF( WNTQS ) THEN 1954* 1955* Perform bidiagonal SVD, computing left singular vectors 1956* of bidiagonal matrix in RWORK(IRU) and computing right 1957* singular vectors of bidiagonal matrix in RWORK(IRVT) 1958* (CWorkspace: need 0) 1959* (RWorkspace: need BDSPAC) 1960* 1961 IRVT = NRWORK 1962 IRU = IRVT + M*M 1963 NRWORK = IRU + M*M 1964 CALL DBDSDC( 'L', 'I', M, S, RWORK( IE ), RWORK( IRU ), 1965 $ M, RWORK( IRVT ), M, DUM, IDUM, 1966 $ RWORK( NRWORK ), IWORK, INFO ) 1967* 1968* Copy real matrix RWORK(IRU) to complex matrix U 1969* Overwrite U by left singular vectors of A 1970* (CWorkspace: need 3*M, prefer 2*M+M*NB) 1971* (RWorkspace: M*M) 1972* 1973 CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU ) 1974 CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA, 1975 $ WORK( ITAUQ ), U, LDU, WORK( NWORK ), 1976 $ LWORK-NWORK+1, IERR ) 1977* 1978* Copy real matrix RWORK(IRVT) to complex matrix VT 1979* Overwrite VT by right singular vectors of A 1980* (CWorkspace: need 3*M, prefer 2*M+M*NB) 1981* (RWorkspace: M*M) 1982* 1983 CALL ZLASET( 'F', M, N, CZERO, CZERO, VT, LDVT ) 1984 CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, VT, LDVT ) 1985 CALL ZUNMBR( 'P', 'R', 'C', M, N, M, A, LDA, 1986 $ WORK( ITAUP ), VT, LDVT, WORK( NWORK ), 1987 $ LWORK-NWORK+1, IERR ) 1988 ELSE 1989* 1990* Perform bidiagonal SVD, computing left singular vectors 1991* of bidiagonal matrix in RWORK(IRU) and computing right 1992* singular vectors of bidiagonal matrix in RWORK(IRVT) 1993* (CWorkspace: need 0) 1994* (RWorkspace: need BDSPAC) 1995* 1996 IRVT = NRWORK 1997 IRU = IRVT + M*M 1998 NRWORK = IRU + M*M 1999* 2000 CALL DBDSDC( 'L', 'I', M, S, RWORK( IE ), RWORK( IRU ), 2001 $ M, RWORK( IRVT ), M, DUM, IDUM, 2002 $ RWORK( NRWORK ), IWORK, INFO ) 2003* 2004* Copy real matrix RWORK(IRU) to complex matrix U 2005* Overwrite U by left singular vectors of A 2006* (CWorkspace: need 3*M, prefer 2*M+M*NB) 2007* (RWorkspace: M*M) 2008* 2009 CALL ZLACP2( 'F', M, M, RWORK( IRU ), M, U, LDU ) 2010 CALL ZUNMBR( 'Q', 'L', 'N', M, M, N, A, LDA, 2011 $ WORK( ITAUQ ), U, LDU, WORK( NWORK ), 2012 $ LWORK-NWORK+1, IERR ) 2013* 2014* Set all of VT to identity matrix 2015* 2016 CALL ZLASET( 'F', N, N, CZERO, CONE, VT, LDVT ) 2017* 2018* Copy real matrix RWORK(IRVT) to complex matrix VT 2019* Overwrite VT by right singular vectors of A 2020* (CWorkspace: need 2*M+N, prefer 2*M+N*NB) 2021* (RWorkspace: M*M) 2022* 2023 CALL ZLACP2( 'F', M, M, RWORK( IRVT ), M, VT, LDVT ) 2024 CALL ZUNMBR( 'P', 'R', 'C', N, N, M, A, LDA, 2025 $ WORK( ITAUP ), VT, LDVT, WORK( NWORK ), 2026 $ LWORK-NWORK+1, IERR ) 2027 END IF 2028* 2029 END IF 2030* 2031 END IF 2032* 2033* Undo scaling if necessary 2034* 2035 IF( ISCL.EQ.1 ) THEN 2036 IF( ANRM.GT.BIGNUM ) 2037 $ CALL DLASCL( 'G', 0, 0, BIGNUM, ANRM, MINMN, 1, S, MINMN, 2038 $ IERR ) 2039 IF( INFO.NE.0 .AND. ANRM.GT.BIGNUM ) 2040 $ CALL DLASCL( 'G', 0, 0, BIGNUM, ANRM, MINMN-1, 1, 2041 $ RWORK( IE ), MINMN, IERR ) 2042 IF( ANRM.LT.SMLNUM ) 2043 $ CALL DLASCL( 'G', 0, 0, SMLNUM, ANRM, MINMN, 1, S, MINMN, 2044 $ IERR ) 2045 IF( INFO.NE.0 .AND. ANRM.LT.SMLNUM ) 2046 $ CALL DLASCL( 'G', 0, 0, SMLNUM, ANRM, MINMN-1, 1, 2047 $ RWORK( IE ), MINMN, IERR ) 2048 END IF 2049* 2050* Return optimal workspace in WORK(1) 2051* 2052 WORK( 1 ) = MAXWRK 2053* 2054 RETURN 2055* 2056* End of ZGESDD 2057* 2058 END 2059