swiftnav: cgesvd.c Source File

Go to the documentation of this file.
00001 /* cgesvd.f -- translated by f2c (version 20061008).
00002    You must link the resulting object file with libf2c:
00003         on Microsoft Windows system, link with libf2c.lib;
00004         on Linux or Unix systems, link with .../path/to/libf2c.a -lm
00005         or, if you install libf2c.a in a standard place, with -lf2c -lm
00006         -- in that order, at the end of the command line, as in
00007                 cc *.o -lf2c -lm
00008         Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
00009 
00010                 http://www.netlib.org/f2c/libf2c.zip
00011 */
00012 
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015 
00016 /* Table of constant values */
00017 
00018 static complex c_b1 = {0.f,0.f};
00019 static complex c_b2 = {1.f,0.f};
00020 static integer c__6 = 6;
00021 static integer c__0 = 0;
00022 static integer c__2 = 2;
00023 static integer c__1 = 1;
00024 static integer c_n1 = -1;
00025 
00026 /* Subroutine */ int cgesvd_(char *jobu, char *jobvt, integer *m, integer *n, 
00027         complex *a, integer *lda, real *s, complex *u, integer *ldu, complex *
00028         vt, integer *ldvt, complex *work, integer *lwork, real *rwork, 
00029         integer *info)
00030 {
00031     /* System generated locals */
00032     address a__1[2];
00033     integer a_dim1, a_offset, u_dim1, u_offset, vt_dim1, vt_offset, i__1[2], 
00034             i__2, i__3, i__4;
00035     char ch__1[2];
00036 
00037     /* Builtin functions */
00038     /* Subroutine */ int s_cat(char *, char **, integer *, integer *, ftnlen);
00039     double sqrt(doublereal);
00040 
00041     /* Local variables */
00042     integer i__, ie, ir, iu, blk, ncu;
00043     real dum[1], eps;
00044     integer nru;
00045     complex cdum[1];
00046     integer iscl;
00047     real anrm;
00048     integer ierr, itau, ncvt, nrvt;
00049     extern /* Subroutine */ int cgemm_(char *, char *, integer *, integer *, 
00050             integer *, complex *, complex *, integer *, complex *, integer *, 
00051             complex *, complex *, integer *);
00052     extern logical lsame_(char *, char *);
00053     integer chunk, minmn, wrkbl, itaup, itauq, mnthr, iwork;
00054     logical wntua, wntva, wntun, wntuo, wntvn, wntvo, wntus, wntvs;
00055     extern /* Subroutine */ int cgebrd_(integer *, integer *, complex *, 
00056             integer *, real *, real *, complex *, complex *, complex *, 
00057             integer *, integer *);
00058     extern doublereal clange_(char *, integer *, integer *, complex *, 
00059             integer *, real *);
00060     extern /* Subroutine */ int cgelqf_(integer *, integer *, complex *, 
00061             integer *, complex *, complex *, integer *, integer *), clascl_(
00062             char *, integer *, integer *, real *, real *, integer *, integer *
00063 , complex *, integer *, integer *), cgeqrf_(integer *, 
00064             integer *, complex *, integer *, complex *, complex *, integer *, 
00065             integer *);
00066     extern doublereal slamch_(char *);
00067     extern /* Subroutine */ int clacpy_(char *, integer *, integer *, complex 
00068             *, integer *, complex *, integer *), claset_(char *, 
00069             integer *, integer *, complex *, complex *, complex *, integer *), cbdsqr_(char *, integer *, integer *, integer *, integer 
00070             *, real *, real *, complex *, integer *, complex *, integer *, 
00071             complex *, integer *, real *, integer *), xerbla_(char *, 
00072             integer *), cungbr_(char *, integer *, integer *, integer 
00073             *, complex *, integer *, complex *, complex *, integer *, integer 
00074             *);
00075     real bignum;
00076     extern /* Subroutine */ int slascl_(char *, integer *, integer *, real *, 
00077             real *, integer *, integer *, real *, integer *, integer *);
00078     extern integer ilaenv_(integer *, char *, char *, integer *, integer *, 
00079             integer *, integer *);
00080     extern /* Subroutine */ int cunmbr_(char *, char *, char *, integer *, 
00081             integer *, integer *, complex *, integer *, complex *, complex *, 
00082             integer *, complex *, integer *, integer *), cunglq_(integer *, integer *, integer *, complex *, 
00083             integer *, complex *, complex *, integer *, integer *), cungqr_(
00084             integer *, integer *, integer *, complex *, integer *, complex *, 
00085             complex *, integer *, integer *);
00086     integer ldwrkr, minwrk, ldwrku, maxwrk;
00087     real smlnum;
00088     integer irwork;
00089     logical lquery, wntuas, wntvas;
00090 
00091 
00092 /*  -- LAPACK driver routine (version 3.2) -- */
00093 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00094 /*     November 2006 */
00095 
00096 /*     .. Scalar Arguments .. */
00097 /*     .. */
00098 /*     .. Array Arguments .. */
00099 /*     .. */
00100 
00101 /*  Purpose */
00102 /*  ======= */
00103 
00104 /*  CGESVD computes the singular value decomposition (SVD) of a complex */
00105 /*  M-by-N matrix A, optionally computing the left and/or right singular */
00106 /*  vectors. The SVD is written */
00107 
00108 /*       A = U * SIGMA * conjugate-transpose(V) */
00109 
00110 /*  where SIGMA is an M-by-N matrix which is zero except for its */
00111 /*  min(m,n) diagonal elements, U is an M-by-M unitary matrix, and */
00112 /*  V is an N-by-N unitary matrix.  The diagonal elements of SIGMA */
00113 /*  are the singular values of A; they are real and non-negative, and */
00114 /*  are returned in descending order.  The first min(m,n) columns of */
00115 /*  U and V are the left and right singular vectors of A. */
00116 
00117 /*  Note that the routine returns V**H, not V. */
00118 
00119 /*  Arguments */
00120 /*  ========= */
00121 
00122 /*  JOBU    (input) CHARACTER*1 */
00123 /*          Specifies options for computing all or part of the matrix U: */
00124 /*          = 'A':  all M columns of U are returned in array U: */
00125 /*          = 'S':  the first min(m,n) columns of U (the left singular */
00126 /*                  vectors) are returned in the array U; */
00127 /*          = 'O':  the first min(m,n) columns of U (the left singular */
00128 /*                  vectors) are overwritten on the array A; */
00129 /*          = 'N':  no columns of U (no left singular vectors) are */
00130 /*                  computed. */
00131 
00132 /*  JOBVT   (input) CHARACTER*1 */
00133 /*          Specifies options for computing all or part of the matrix */
00134 /*          V**H: */
00135 /*          = 'A':  all N rows of V**H are returned in the array VT; */
00136 /*          = 'S':  the first min(m,n) rows of V**H (the right singular */
00137 /*                  vectors) are returned in the array VT; */
00138 /*          = 'O':  the first min(m,n) rows of V**H (the right singular */
00139 /*                  vectors) are overwritten on the array A; */
00140 /*          = 'N':  no rows of V**H (no right singular vectors) are */
00141 /*                  computed. */
00142 
00143 /*          JOBVT and JOBU cannot both be 'O'. */
00144 
00145 /*  M       (input) INTEGER */
00146 /*          The number of rows of the input matrix A.  M >= 0. */
00147 
00148 /*  N       (input) INTEGER */
00149 /*          The number of columns of the input matrix A.  N >= 0. */
00150 
00151 /*  A       (input/output) COMPLEX array, dimension (LDA,N) */
00152 /*          On entry, the M-by-N matrix A. */
00153 /*          On exit, */
00154 /*          if JOBU = 'O',  A is overwritten with the first min(m,n) */
00155 /*                          columns of U (the left singular vectors, */
00156 /*                          stored columnwise); */
00157 /*          if JOBVT = 'O', A is overwritten with the first min(m,n) */
00158 /*                          rows of V**H (the right singular vectors, */
00159 /*                          stored rowwise); */
00160 /*          if JOBU .ne. 'O' and JOBVT .ne. 'O', the contents of A */
00161 /*                          are destroyed. */
00162 
00163 /*  LDA     (input) INTEGER */
00164 /*          The leading dimension of the array A.  LDA >= max(1,M). */
00165 
00166 /*  S       (output) REAL array, dimension (min(M,N)) */
00167 /*          The singular values of A, sorted so that S(i) >= S(i+1). */
00168 
00169 /*  U       (output) COMPLEX array, dimension (LDU,UCOL) */
00170 /*          (LDU,M) if JOBU = 'A' or (LDU,min(M,N)) if JOBU = 'S'. */
00171 /*          If JOBU = 'A', U contains the M-by-M unitary matrix U; */
00172 /*          if JOBU = 'S', U contains the first min(m,n) columns of U */
00173 /*          (the left singular vectors, stored columnwise); */
00174 /*          if JOBU = 'N' or 'O', U is not referenced. */
00175 
00176 /*  LDU     (input) INTEGER */
00177 /*          The leading dimension of the array U.  LDU >= 1; if */
00178 /*          JOBU = 'S' or 'A', LDU >= M. */
00179 
00180 /*  VT      (output) COMPLEX array, dimension (LDVT,N) */
00181 /*          If JOBVT = 'A', VT contains the N-by-N unitary matrix */
00182 /*          V**H; */
00183 /*          if JOBVT = 'S', VT contains the first min(m,n) rows of */
00184 /*          V**H (the right singular vectors, stored rowwise); */
00185 /*          if JOBVT = 'N' or 'O', VT is not referenced. */
00186 
00187 /*  LDVT    (input) INTEGER */
00188 /*          The leading dimension of the array VT.  LDVT >= 1; if */
00189 /*          JOBVT = 'A', LDVT >= N; if JOBVT = 'S', LDVT >= min(M,N). */
00190 
00191 /*  WORK    (workspace/output) COMPLEX array, dimension (MAX(1,LWORK)) */
00192 /*          On exit, if INFO = 0, WORK(1) returns the optimal LWORK. */
00193 
00194 /*  LWORK   (input) INTEGER */
00195 /*          The dimension of the array WORK. */
00196 /*          LWORK >=  MAX(1,2*MIN(M,N)+MAX(M,N)). */
00197 /*          For good performance, LWORK should generally be larger. */
00198 
00199 /*          If LWORK = -1, then a workspace query is assumed; the routine */
00200 /*          only calculates the optimal size of the WORK array, returns */
00201 /*          this value as the first entry of the WORK array, and no error */
00202 /*          message related to LWORK is issued by XERBLA. */
00203 
00204 /*  RWORK   (workspace) REAL array, dimension (5*min(M,N)) */
00205 /*          On exit, if INFO > 0, RWORK(1:MIN(M,N)-1) contains the */
00206 /*          unconverged superdiagonal elements of an upper bidiagonal */
00207 /*          matrix B whose diagonal is in S (not necessarily sorted). */
00208 /*          B satisfies A = U * B * VT, so it has the same singular */
00209 /*          values as A, and singular vectors related by U and VT. */
00210 
00211 /*  INFO    (output) INTEGER */
00212 /*          = 0:  successful exit. */
00213 /*          < 0:  if INFO = -i, the i-th argument had an illegal value. */
00214 /*          > 0:  if CBDSQR did not converge, INFO specifies how many */
00215 /*                superdiagonals of an intermediate bidiagonal form B */
00216 /*                did not converge to zero. See the description of RWORK */
00217 /*                above for details. */
00218 
00219 /*  ===================================================================== */
00220 
00221 /*     .. Parameters .. */
00222 /*     .. */
00223 /*     .. Local Scalars .. */
00224 /*     .. */
00225 /*     .. Local Arrays .. */
00226 /*     .. */
00227 /*     .. External Subroutines .. */
00228 /*     .. */
00229 /*     .. External Functions .. */
00230 /*     .. */
00231 /*     .. Intrinsic Functions .. */
00232 /*     .. */
00233 /*     .. Executable Statements .. */
00234 
00235 /*     Test the input arguments */
00236 
00237     /* Parameter adjustments */
00238     a_dim1 = *lda;
00239     a_offset = 1 + a_dim1;
00240     a -= a_offset;
00241     --s;
00242     u_dim1 = *ldu;
00243     u_offset = 1 + u_dim1;
00244     u -= u_offset;
00245     vt_dim1 = *ldvt;
00246     vt_offset = 1 + vt_dim1;
00247     vt -= vt_offset;
00248     --work;
00249     --rwork;
00250 
00251     /* Function Body */
00252     *info = 0;
00253     minmn = min(*m,*n);
00254     wntua = lsame_(jobu, "A");
00255     wntus = lsame_(jobu, "S");
00256     wntuas = wntua || wntus;
00257     wntuo = lsame_(jobu, "O");
00258     wntun = lsame_(jobu, "N");
00259     wntva = lsame_(jobvt, "A");
00260     wntvs = lsame_(jobvt, "S");
00261     wntvas = wntva || wntvs;
00262     wntvo = lsame_(jobvt, "O");
00263     wntvn = lsame_(jobvt, "N");
00264     lquery = *lwork == -1;
00265 
00266     if (! (wntua || wntus || wntuo || wntun)) {
00267         *info = -1;
00268     } else if (! (wntva || wntvs || wntvo || wntvn) || wntvo && wntuo) {
00269         *info = -2;
00270     } else if (*m < 0) {
00271         *info = -3;
00272     } else if (*n < 0) {
00273         *info = -4;
00274     } else if (*lda < max(1,*m)) {
00275         *info = -6;
00276     } else if (*ldu < 1 || wntuas && *ldu < *m) {
00277         *info = -9;
00278     } else if (*ldvt < 1 || wntva && *ldvt < *n || wntvs && *ldvt < minmn) {
00279         *info = -11;
00280     }
00281 
00282 /*     Compute workspace */
00283 /*      (Note: Comments in the code beginning "Workspace:" describe the */
00284 /*       minimal amount of workspace needed at that point in the code, */
00285 /*       as well as the preferred amount for good performance. */
00286 /*       CWorkspace refers to complex workspace, and RWorkspace to */
00287 /*       real workspace. NB refers to the optimal block size for the */
00288 /*       immediately following subroutine, as returned by ILAENV.) */
00289 
00290     if (*info == 0) {
00291         minwrk = 1;
00292         maxwrk = 1;
00293         if (*m >= *n && minmn > 0) {
00294 
00295 /*           Space needed for CBDSQR is BDSPAC = 5*N */
00296 
00297 /* Writing concatenation */
00298             i__1[0] = 1, a__1[0] = jobu;
00299             i__1[1] = 1, a__1[1] = jobvt;
00300             s_cat(ch__1, a__1, i__1, &c__2, (ftnlen)2);
00301             mnthr = ilaenv_(&c__6, "CGESVD", ch__1, m, n, &c__0, &c__0);
00302             if (*m >= mnthr) {
00303                 if (wntun) {
00304 
00305 /*                 Path 1 (M much larger than N, JOBU='N') */
00306 
00307                     maxwrk = *n + *n * ilaenv_(&c__1, "CGEQRF", " ", m, n, &
00308                             c_n1, &c_n1);
00309 /* Computing MAX */
00310                     i__2 = maxwrk, i__3 = (*n << 1) + (*n << 1) * ilaenv_(&
00311                             c__1, "CGEBRD", " ", n, n, &c_n1, &c_n1);
00312                     maxwrk = max(i__2,i__3);
00313                     if (wntvo || wntvas) {
00314 /* Computing MAX */
00315                         i__2 = maxwrk, i__3 = (*n << 1) + (*n - 1) * ilaenv_(&
00316                                 c__1, "CUNGBR", "P", n, n, n, &c_n1);
00317                         maxwrk = max(i__2,i__3);
00318                     }
00319                     minwrk = *n * 3;
00320                 } else if (wntuo && wntvn) {
00321 
00322 /*                 Path 2 (M much larger than N, JOBU='O', JOBVT='N') */
00323 
00324                     wrkbl = *n + *n * ilaenv_(&c__1, "CGEQRF", " ", m, n, &
00325                             c_n1, &c_n1);
00326 /* Computing MAX */
00327                     i__2 = wrkbl, i__3 = *n + *n * ilaenv_(&c__1, "CUNGQR", 
00328                             " ", m, n, n, &c_n1);
00329                     wrkbl = max(i__2,i__3);
00330 /* Computing MAX */
00331                     i__2 = wrkbl, i__3 = (*n << 1) + (*n << 1) * ilaenv_(&
00332                             c__1, "CGEBRD", " ", n, n, &c_n1, &c_n1);
00333                     wrkbl = max(i__2,i__3);
00334 /* Computing MAX */
00335                     i__2 = wrkbl, i__3 = (*n << 1) + *n * ilaenv_(&c__1, 
00336                             "CUNGBR", "Q", n, n, n, &c_n1);
00337                     wrkbl = max(i__2,i__3);
00338 /* Computing MAX */
00339                     i__2 = *n * *n + wrkbl, i__3 = *n * *n + *m * *n;
00340                     maxwrk = max(i__2,i__3);
00341                     minwrk = (*n << 1) + *m;
00342                 } else if (wntuo && wntvas) {
00343 
00344 /*                 Path 3 (M much larger than N, JOBU='O', JOBVT='S' or */
00345 /*                 'A') */
00346 
00347                     wrkbl = *n + *n * ilaenv_(&c__1, "CGEQRF", " ", m, n, &
00348                             c_n1, &c_n1);
00349 /* Computing MAX */
00350                     i__2 = wrkbl, i__3 = *n + *n * ilaenv_(&c__1, "CUNGQR", 
00351                             " ", m, n, n, &c_n1);
00352                     wrkbl = max(i__2,i__3);
00353 /* Computing MAX */
00354                     i__2 = wrkbl, i__3 = (*n << 1) + (*n << 1) * ilaenv_(&
00355                             c__1, "CGEBRD", " ", n, n, &c_n1, &c_n1);
00356                     wrkbl = max(i__2,i__3);
00357 /* Computing MAX */
00358                     i__2 = wrkbl, i__3 = (*n << 1) + *n * ilaenv_(&c__1, 
00359                             "CUNGBR", "Q", n, n, n, &c_n1);
00360                     wrkbl = max(i__2,i__3);
00361 /* Computing MAX */
00362                     i__2 = wrkbl, i__3 = (*n << 1) + (*n - 1) * ilaenv_(&c__1, 
00363                              "CUNGBR", "P", n, n, n, &c_n1);
00364                     wrkbl = max(i__2,i__3);
00365 /* Computing MAX */
00366                     i__2 = *n * *n + wrkbl, i__3 = *n * *n + *m * *n;
00367                     maxwrk = max(i__2,i__3);
00368                     minwrk = (*n << 1) + *m;
00369                 } else if (wntus && wntvn) {
00370 
00371 /*                 Path 4 (M much larger than N, JOBU='S', JOBVT='N') */
00372 
00373                     wrkbl = *n + *n * ilaenv_(&c__1, "CGEQRF", " ", m, n, &
00374                             c_n1, &c_n1);
00375 /* Computing MAX */
00376                     i__2 = wrkbl, i__3 = *n + *n * ilaenv_(&c__1, "CUNGQR", 
00377                             " ", m, n, n, &c_n1);
00378                     wrkbl = max(i__2,i__3);
00379 /* Computing MAX */
00380                     i__2 = wrkbl, i__3 = (*n << 1) + (*n << 1) * ilaenv_(&
00381                             c__1, "CGEBRD", " ", n, n, &c_n1, &c_n1);
00382                     wrkbl = max(i__2,i__3);
00383 /* Computing MAX */
00384                     i__2 = wrkbl, i__3 = (*n << 1) + *n * ilaenv_(&c__1, 
00385                             "CUNGBR", "Q", n, n, n, &c_n1);
00386                     wrkbl = max(i__2,i__3);
00387                     maxwrk = *n * *n + wrkbl;
00388                     minwrk = (*n << 1) + *m;
00389                 } else if (wntus && wntvo) {
00390 
00391 /*                 Path 5 (M much larger than N, JOBU='S', JOBVT='O') */
00392 
00393                     wrkbl = *n + *n * ilaenv_(&c__1, "CGEQRF", " ", m, n, &
00394                             c_n1, &c_n1);
00395 /* Computing MAX */
00396                     i__2 = wrkbl, i__3 = *n + *n * ilaenv_(&c__1, "CUNGQR", 
00397                             " ", m, n, n, &c_n1);
00398                     wrkbl = max(i__2,i__3);
00399 /* Computing MAX */
00400                     i__2 = wrkbl, i__3 = (*n << 1) + (*n << 1) * ilaenv_(&
00401                             c__1, "CGEBRD", " ", n, n, &c_n1, &c_n1);
00402                     wrkbl = max(i__2,i__3);
00403 /* Computing MAX */
00404                     i__2 = wrkbl, i__3 = (*n << 1) + *n * ilaenv_(&c__1, 
00405                             "CUNGBR", "Q", n, n, n, &c_n1);
00406                     wrkbl = max(i__2,i__3);
00407 /* Computing MAX */
00408                     i__2 = wrkbl, i__3 = (*n << 1) + (*n - 1) * ilaenv_(&c__1, 
00409                              "CUNGBR", "P", n, n, n, &c_n1);
00410                     wrkbl = max(i__2,i__3);
00411                     maxwrk = (*n << 1) * *n + wrkbl;
00412                     minwrk = (*n << 1) + *m;
00413                 } else if (wntus && wntvas) {
00414 
00415 /*                 Path 6 (M much larger than N, JOBU='S', JOBVT='S' or */
00416 /*                 'A') */
00417 
00418                     wrkbl = *n + *n * ilaenv_(&c__1, "CGEQRF", " ", m, n, &
00419                             c_n1, &c_n1);
00420 /* Computing MAX */
00421                     i__2 = wrkbl, i__3 = *n + *n * ilaenv_(&c__1, "CUNGQR", 
00422                             " ", m, n, n, &c_n1);
00423                     wrkbl = max(i__2,i__3);
00424 /* Computing MAX */
00425                     i__2 = wrkbl, i__3 = (*n << 1) + (*n << 1) * ilaenv_(&
00426                             c__1, "CGEBRD", " ", n, n, &c_n1, &c_n1);
00427                     wrkbl = max(i__2,i__3);
00428 /* Computing MAX */
00429                     i__2 = wrkbl, i__3 = (*n << 1) + *n * ilaenv_(&c__1, 
00430                             "CUNGBR", "Q", n, n, n, &c_n1);
00431                     wrkbl = max(i__2,i__3);
00432 /* Computing MAX */
00433                     i__2 = wrkbl, i__3 = (*n << 1) + (*n - 1) * ilaenv_(&c__1, 
00434                              "CUNGBR", "P", n, n, n, &c_n1);
00435                     wrkbl = max(i__2,i__3);
00436                     maxwrk = *n * *n + wrkbl;
00437                     minwrk = (*n << 1) + *m;
00438                 } else if (wntua && wntvn) {
00439 
00440 /*                 Path 7 (M much larger than N, JOBU='A', JOBVT='N') */
00441 
00442                     wrkbl = *n + *n * ilaenv_(&c__1, "CGEQRF", " ", m, n, &
00443                             c_n1, &c_n1);
00444 /* Computing MAX */
00445                     i__2 = wrkbl, i__3 = *n + *m * ilaenv_(&c__1, "CUNGQR", 
00446                             " ", m, m, n, &c_n1);
00447                     wrkbl = max(i__2,i__3);
00448 /* Computing MAX */
00449                     i__2 = wrkbl, i__3 = (*n << 1) + (*n << 1) * ilaenv_(&
00450                             c__1, "CGEBRD", " ", n, n, &c_n1, &c_n1);
00451                     wrkbl = max(i__2,i__3);
00452 /* Computing MAX */
00453                     i__2 = wrkbl, i__3 = (*n << 1) + *n * ilaenv_(&c__1, 
00454                             "CUNGBR", "Q", n, n, n, &c_n1);
00455                     wrkbl = max(i__2,i__3);
00456                     maxwrk = *n * *n + wrkbl;
00457                     minwrk = (*n << 1) + *m;
00458                 } else if (wntua && wntvo) {
00459 
00460 /*                 Path 8 (M much larger than N, JOBU='A', JOBVT='O') */
00461 
00462                     wrkbl = *n + *n * ilaenv_(&c__1, "CGEQRF", " ", m, n, &
00463                             c_n1, &c_n1);
00464 /* Computing MAX */
00465                     i__2 = wrkbl, i__3 = *n + *m * ilaenv_(&c__1, "CUNGQR", 
00466                             " ", m, m, n, &c_n1);
00467                     wrkbl = max(i__2,i__3);
00468 /* Computing MAX */
00469                     i__2 = wrkbl, i__3 = (*n << 1) + (*n << 1) * ilaenv_(&
00470                             c__1, "CGEBRD", " ", n, n, &c_n1, &c_n1);
00471                     wrkbl = max(i__2,i__3);
00472 /* Computing MAX */
00473                     i__2 = wrkbl, i__3 = (*n << 1) + *n * ilaenv_(&c__1, 
00474                             "CUNGBR", "Q", n, n, n, &c_n1);
00475                     wrkbl = max(i__2,i__3);
00476 /* Computing MAX */
00477                     i__2 = wrkbl, i__3 = (*n << 1) + (*n - 1) * ilaenv_(&c__1, 
00478                              "CUNGBR", "P", n, n, n, &c_n1);
00479                     wrkbl = max(i__2,i__3);
00480                     maxwrk = (*n << 1) * *n + wrkbl;
00481                     minwrk = (*n << 1) + *m;
00482                 } else if (wntua && wntvas) {
00483 
00484 /*                 Path 9 (M much larger than N, JOBU='A', JOBVT='S' or */
00485 /*                 'A') */
00486 
00487                     wrkbl = *n + *n * ilaenv_(&c__1, "CGEQRF", " ", m, n, &
00488                             c_n1, &c_n1);
00489 /* Computing MAX */
00490                     i__2 = wrkbl, i__3 = *n + *m * ilaenv_(&c__1, "CUNGQR", 
00491                             " ", m, m, n, &c_n1);
00492                     wrkbl = max(i__2,i__3);
00493 /* Computing MAX */
00494                     i__2 = wrkbl, i__3 = (*n << 1) + (*n << 1) * ilaenv_(&
00495                             c__1, "CGEBRD", " ", n, n, &c_n1, &c_n1);
00496                     wrkbl = max(i__2,i__3);
00497 /* Computing MAX */
00498                     i__2 = wrkbl, i__3 = (*n << 1) + *n * ilaenv_(&c__1, 
00499                             "CUNGBR", "Q", n, n, n, &c_n1);
00500                     wrkbl = max(i__2,i__3);
00501 /* Computing MAX */
00502                     i__2 = wrkbl, i__3 = (*n << 1) + (*n - 1) * ilaenv_(&c__1, 
00503                              "CUNGBR", "P", n, n, n, &c_n1);
00504                     wrkbl = max(i__2,i__3);
00505                     maxwrk = *n * *n + wrkbl;
00506                     minwrk = (*n << 1) + *m;
00507                 }
00508             } else {
00509 
00510 /*              Path 10 (M at least N, but not much larger) */
00511 
00512                 maxwrk = (*n << 1) + (*m + *n) * ilaenv_(&c__1, "CGEBRD", 
00513                         " ", m, n, &c_n1, &c_n1);
00514                 if (wntus || wntuo) {
00515 /* Computing MAX */
00516                     i__2 = maxwrk, i__3 = (*n << 1) + *n * ilaenv_(&c__1, 
00517                             "CUNGBR", "Q", m, n, n, &c_n1);
00518                     maxwrk = max(i__2,i__3);
00519                 }
00520                 if (wntua) {
00521 /* Computing MAX */
00522                     i__2 = maxwrk, i__3 = (*n << 1) + *m * ilaenv_(&c__1, 
00523                             "CUNGBR", "Q", m, m, n, &c_n1);
00524                     maxwrk = max(i__2,i__3);
00525                 }
00526                 if (! wntvn) {
00527 /* Computing MAX */
00528                     i__2 = maxwrk, i__3 = (*n << 1) + (*n - 1) * ilaenv_(&
00529                             c__1, "CUNGBR", "P", n, n, n, &c_n1);
00530                     maxwrk = max(i__2,i__3);
00531                 }
00532                 minwrk = (*n << 1) + *m;
00533             }
00534         } else if (minmn > 0) {
00535 
00536 /*           Space needed for CBDSQR is BDSPAC = 5*M */
00537 
00538 /* Writing concatenation */
00539             i__1[0] = 1, a__1[0] = jobu;
00540             i__1[1] = 1, a__1[1] = jobvt;
00541             s_cat(ch__1, a__1, i__1, &c__2, (ftnlen)2);
00542             mnthr = ilaenv_(&c__6, "CGESVD", ch__1, m, n, &c__0, &c__0);
00543             if (*n >= mnthr) {
00544                 if (wntvn) {
00545 
00546 /*                 Path 1t(N much larger than M, JOBVT='N') */
00547 
00548                     maxwrk = *m + *m * ilaenv_(&c__1, "CGELQF", " ", m, n, &
00549                             c_n1, &c_n1);
00550 /* Computing MAX */
00551                     i__2 = maxwrk, i__3 = (*m << 1) + (*m << 1) * ilaenv_(&
00552                             c__1, "CGEBRD", " ", m, m, &c_n1, &c_n1);
00553                     maxwrk = max(i__2,i__3);
00554                     if (wntuo || wntuas) {
00555 /* Computing MAX */
00556                         i__2 = maxwrk, i__3 = (*m << 1) + *m * ilaenv_(&c__1, 
00557                                 "CUNGBR", "Q", m, m, m, &c_n1);
00558                         maxwrk = max(i__2,i__3);
00559                     }
00560                     minwrk = *m * 3;
00561                 } else if (wntvo && wntun) {
00562 
00563 /*                 Path 2t(N much larger than M, JOBU='N', JOBVT='O') */
00564 
00565                     wrkbl = *m + *m * ilaenv_(&c__1, "CGELQF", " ", m, n, &
00566                             c_n1, &c_n1);
00567 /* Computing MAX */
00568                     i__2 = wrkbl, i__3 = *m + *m * ilaenv_(&c__1, "CUNGLQ", 
00569                             " ", m, n, m, &c_n1);
00570                     wrkbl = max(i__2,i__3);
00571 /* Computing MAX */
00572                     i__2 = wrkbl, i__3 = (*m << 1) + (*m << 1) * ilaenv_(&
00573                             c__1, "CGEBRD", " ", m, m, &c_n1, &c_n1);
00574                     wrkbl = max(i__2,i__3);
00575 /* Computing MAX */
00576                     i__2 = wrkbl, i__3 = (*m << 1) + (*m - 1) * ilaenv_(&c__1, 
00577                              "CUNGBR", "P", m, m, m, &c_n1);
00578                     wrkbl = max(i__2,i__3);
00579 /* Computing MAX */
00580                     i__2 = *m * *m + wrkbl, i__3 = *m * *m + *m * *n;
00581                     maxwrk = max(i__2,i__3);
00582                     minwrk = (*m << 1) + *n;
00583                 } else if (wntvo && wntuas) {
00584 
00585 /*                 Path 3t(N much larger than M, JOBU='S' or 'A', */
00586 /*                 JOBVT='O') */
00587 
00588                     wrkbl = *m + *m * ilaenv_(&c__1, "CGELQF", " ", m, n, &
00589                             c_n1, &c_n1);
00590 /* Computing MAX */
00591                     i__2 = wrkbl, i__3 = *m + *m * ilaenv_(&c__1, "CUNGLQ", 
00592                             " ", m, n, m, &c_n1);
00593                     wrkbl = max(i__2,i__3);
00594 /* Computing MAX */
00595                     i__2 = wrkbl, i__3 = (*m << 1) + (*m << 1) * ilaenv_(&
00596                             c__1, "CGEBRD", " ", m, m, &c_n1, &c_n1);
00597                     wrkbl = max(i__2,i__3);
00598 /* Computing MAX */
00599                     i__2 = wrkbl, i__3 = (*m << 1) + (*m - 1) * ilaenv_(&c__1, 
00600                              "CUNGBR", "P", m, m, m, &c_n1);
00601                     wrkbl = max(i__2,i__3);
00602 /* Computing MAX */
00603                     i__2 = wrkbl, i__3 = (*m << 1) + *m * ilaenv_(&c__1, 
00604                             "CUNGBR", "Q", m, m, m, &c_n1);
00605                     wrkbl = max(i__2,i__3);
00606 /* Computing MAX */
00607                     i__2 = *m * *m + wrkbl, i__3 = *m * *m + *m * *n;
00608                     maxwrk = max(i__2,i__3);
00609                     minwrk = (*m << 1) + *n;
00610                 } else if (wntvs && wntun) {
00611 
00612 /*                 Path 4t(N much larger than M, JOBU='N', JOBVT='S') */
00613 
00614                     wrkbl = *m + *m * ilaenv_(&c__1, "CGELQF", " ", m, n, &
00615                             c_n1, &c_n1);
00616 /* Computing MAX */
00617                     i__2 = wrkbl, i__3 = *m + *m * ilaenv_(&c__1, "CUNGLQ", 
00618                             " ", m, n, m, &c_n1);
00619                     wrkbl = max(i__2,i__3);
00620 /* Computing MAX */
00621                     i__2 = wrkbl, i__3 = (*m << 1) + (*m << 1) * ilaenv_(&
00622                             c__1, "CGEBRD", " ", m, m, &c_n1, &c_n1);
00623                     wrkbl = max(i__2,i__3);
00624 /* Computing MAX */
00625                     i__2 = wrkbl, i__3 = (*m << 1) + (*m - 1) * ilaenv_(&c__1, 
00626                              "CUNGBR", "P", m, m, m, &c_n1);
00627                     wrkbl = max(i__2,i__3);
00628                     maxwrk = *m * *m + wrkbl;
00629                     minwrk = (*m << 1) + *n;
00630                 } else if (wntvs && wntuo) {
00631 
00632 /*                 Path 5t(N much larger than M, JOBU='O', JOBVT='S') */
00633 
00634                     wrkbl = *m + *m * ilaenv_(&c__1, "CGELQF", " ", m, n, &
00635                             c_n1, &c_n1);
00636 /* Computing MAX */
00637                     i__2 = wrkbl, i__3 = *m + *m * ilaenv_(&c__1, "CUNGLQ", 
00638                             " ", m, n, m, &c_n1);
00639                     wrkbl = max(i__2,i__3);
00640 /* Computing MAX */
00641                     i__2 = wrkbl, i__3 = (*m << 1) + (*m << 1) * ilaenv_(&
00642                             c__1, "CGEBRD", " ", m, m, &c_n1, &c_n1);
00643                     wrkbl = max(i__2,i__3);
00644 /* Computing MAX */
00645                     i__2 = wrkbl, i__3 = (*m << 1) + (*m - 1) * ilaenv_(&c__1, 
00646                              "CUNGBR", "P", m, m, m, &c_n1);
00647                     wrkbl = max(i__2,i__3);
00648 /* Computing MAX */
00649                     i__2 = wrkbl, i__3 = (*m << 1) + *m * ilaenv_(&c__1, 
00650                             "CUNGBR", "Q", m, m, m, &c_n1);
00651                     wrkbl = max(i__2,i__3);
00652                     maxwrk = (*m << 1) * *m + wrkbl;
00653                     minwrk = (*m << 1) + *n;
00654                 } else if (wntvs && wntuas) {
00655 
00656 /*                 Path 6t(N much larger than M, JOBU='S' or 'A', */
00657 /*                 JOBVT='S') */
00658 
00659                     wrkbl = *m + *m * ilaenv_(&c__1, "CGELQF", " ", m, n, &
00660                             c_n1, &c_n1);
00661 /* Computing MAX */
00662                     i__2 = wrkbl, i__3 = *m + *m * ilaenv_(&c__1, "CUNGLQ", 
00663                             " ", m, n, m, &c_n1);
00664                     wrkbl = max(i__2,i__3);
00665 /* Computing MAX */
00666                     i__2 = wrkbl, i__3 = (*m << 1) + (*m << 1) * ilaenv_(&
00667                             c__1, "CGEBRD", " ", m, m, &c_n1, &c_n1);
00668                     wrkbl = max(i__2,i__3);
00669 /* Computing MAX */
00670                     i__2 = wrkbl, i__3 = (*m << 1) + (*m - 1) * ilaenv_(&c__1, 
00671                              "CUNGBR", "P", m, m, m, &c_n1);
00672                     wrkbl = max(i__2,i__3);
00673 /* Computing MAX */
00674                     i__2 = wrkbl, i__3 = (*m << 1) + *m * ilaenv_(&c__1, 
00675                             "CUNGBR", "Q", m, m, m, &c_n1);
00676                     wrkbl = max(i__2,i__3);
00677                     maxwrk = *m * *m + wrkbl;
00678                     minwrk = (*m << 1) + *n;
00679                 } else if (wntva && wntun) {
00680 
00681 /*                 Path 7t(N much larger than M, JOBU='N', JOBVT='A') */
00682 
00683                     wrkbl = *m + *m * ilaenv_(&c__1, "CGELQF", " ", m, n, &
00684                             c_n1, &c_n1);
00685 /* Computing MAX */
00686                     i__2 = wrkbl, i__3 = *m + *n * ilaenv_(&c__1, "CUNGLQ", 
00687                             " ", n, n, m, &c_n1);
00688                     wrkbl = max(i__2,i__3);
00689 /* Computing MAX */
00690                     i__2 = wrkbl, i__3 = (*m << 1) + (*m << 1) * ilaenv_(&
00691                             c__1, "CGEBRD", " ", m, m, &c_n1, &c_n1);
00692                     wrkbl = max(i__2,i__3);
00693 /* Computing MAX */
00694                     i__2 = wrkbl, i__3 = (*m << 1) + (*m - 1) * ilaenv_(&c__1, 
00695                              "CUNGBR", "P", m, m, m, &c_n1);
00696                     wrkbl = max(i__2,i__3);
00697                     maxwrk = *m * *m + wrkbl;
00698                     minwrk = (*m << 1) + *n;
00699                 } else if (wntva && wntuo) {
00700 
00701 /*                 Path 8t(N much larger than M, JOBU='O', JOBVT='A') */
00702 
00703                     wrkbl = *m + *m * ilaenv_(&c__1, "CGELQF", " ", m, n, &
00704                             c_n1, &c_n1);
00705 /* Computing MAX */
00706                     i__2 = wrkbl, i__3 = *m + *n * ilaenv_(&c__1, "CUNGLQ", 
00707                             " ", n, n, m, &c_n1);
00708                     wrkbl = max(i__2,i__3);
00709 /* Computing MAX */
00710                     i__2 = wrkbl, i__3 = (*m << 1) + (*m << 1) * ilaenv_(&
00711                             c__1, "CGEBRD", " ", m, m, &c_n1, &c_n1);
00712                     wrkbl = max(i__2,i__3);
00713 /* Computing MAX */
00714                     i__2 = wrkbl, i__3 = (*m << 1) + (*m - 1) * ilaenv_(&c__1, 
00715                              "CUNGBR", "P", m, m, m, &c_n1);
00716                     wrkbl = max(i__2,i__3);
00717 /* Computing MAX */
00718                     i__2 = wrkbl, i__3 = (*m << 1) + *m * ilaenv_(&c__1, 
00719                             "CUNGBR", "Q", m, m, m, &c_n1);
00720                     wrkbl = max(i__2,i__3);
00721                     maxwrk = (*m << 1) * *m + wrkbl;
00722                     minwrk = (*m << 1) + *n;
00723                 } else if (wntva && wntuas) {
00724 
00725 /*                 Path 9t(N much larger than M, JOBU='S' or 'A', */
00726 /*                 JOBVT='A') */
00727 
00728                     wrkbl = *m + *m * ilaenv_(&c__1, "CGELQF", " ", m, n, &
00729                             c_n1, &c_n1);
00730 /* Computing MAX */
00731                     i__2 = wrkbl, i__3 = *m + *n * ilaenv_(&c__1, "CUNGLQ", 
00732                             " ", n, n, m, &c_n1);
00733                     wrkbl = max(i__2,i__3);
00734 /* Computing MAX */
00735                     i__2 = wrkbl, i__3 = (*m << 1) + (*m << 1) * ilaenv_(&
00736                             c__1, "CGEBRD", " ", m, m, &c_n1, &c_n1);
00737                     wrkbl = max(i__2,i__3);
00738 /* Computing MAX */
00739                     i__2 = wrkbl, i__3 = (*m << 1) + (*m - 1) * ilaenv_(&c__1, 
00740                              "CUNGBR", "P", m, m, m, &c_n1);
00741                     wrkbl = max(i__2,i__3);
00742 /* Computing MAX */
00743                     i__2 = wrkbl, i__3 = (*m << 1) + *m * ilaenv_(&c__1, 
00744                             "CUNGBR", "Q", m, m, m, &c_n1);
00745                     wrkbl = max(i__2,i__3);
00746                     maxwrk = *m * *m + wrkbl;
00747                     minwrk = (*m << 1) + *n;
00748                 }
00749             } else {
00750 
00751 /*              Path 10t(N greater than M, but not much larger) */
00752 
00753                 maxwrk = (*m << 1) + (*m + *n) * ilaenv_(&c__1, "CGEBRD", 
00754                         " ", m, n, &c_n1, &c_n1);
00755                 if (wntvs || wntvo) {
00756 /* Computing MAX */
00757                     i__2 = maxwrk, i__3 = (*m << 1) + *m * ilaenv_(&c__1, 
00758                             "CUNGBR", "P", m, n, m, &c_n1);
00759                     maxwrk = max(i__2,i__3);
00760                 }
00761                 if (wntva) {
00762 /* Computing MAX */
00763                     i__2 = maxwrk, i__3 = (*m << 1) + *n * ilaenv_(&c__1, 
00764                             "CUNGBR", "P", n, n, m, &c_n1);
00765                     maxwrk = max(i__2,i__3);
00766                 }
00767                 if (! wntun) {
00768 /* Computing MAX */
00769                     i__2 = maxwrk, i__3 = (*m << 1) + (*m - 1) * ilaenv_(&
00770                             c__1, "CUNGBR", "Q", m, m, m, &c_n1);
00771                     maxwrk = max(i__2,i__3);
00772                 }
00773                 minwrk = (*m << 1) + *n;
00774             }
00775         }
00776         maxwrk = max(minwrk,maxwrk);
00777         work[1].r = (real) maxwrk, work[1].i = 0.f;
00778 
00779         if (*lwork < minwrk && ! lquery) {
00780             *info = -13;
00781         }
00782     }
00783 
00784     if (*info != 0) {
00785         i__2 = -(*info);
00786         xerbla_("CGESVD", &i__2);
00787         return 0;
00788     } else if (lquery) {
00789         return 0;
00790     }
00791 
00792 /*     Quick return if possible */
00793 
00794     if (*m == 0 || *n == 0) {
00795         return 0;
00796     }
00797 
00798 /*     Get machine constants */
00799 
00800     eps = slamch_("P");
00801     smlnum = sqrt(slamch_("S")) / eps;
00802     bignum = 1.f / smlnum;
00803 
00804 /*     Scale A if max element outside range [SMLNUM,BIGNUM] */
00805 
00806     anrm = clange_("M", m, n, &a[a_offset], lda, dum);
00807     iscl = 0;
00808     if (anrm > 0.f && anrm < smlnum) {
00809         iscl = 1;
00810         clascl_("G", &c__0, &c__0, &anrm, &smlnum, m, n, &a[a_offset], lda, &
00811                 ierr);
00812     } else if (anrm > bignum) {
00813         iscl = 1;
00814         clascl_("G", &c__0, &c__0, &anrm, &bignum, m, n, &a[a_offset], lda, &
00815                 ierr);
00816     }
00817 
00818     if (*m >= *n) {
00819 
00820 /*        A has at least as many rows as columns. If A has sufficiently */
00821 /*        more rows than columns, first reduce using the QR */
00822 /*        decomposition (if sufficient workspace available) */
00823 
00824         if (*m >= mnthr) {
00825 
00826             if (wntun) {
00827 
00828 /*              Path 1 (M much larger than N, JOBU='N') */
00829 /*              No left singular vectors to be computed */
00830 
00831                 itau = 1;
00832                 iwork = itau + *n;
00833 
00834 /*              Compute A=Q*R */
00835 /*              (CWorkspace: need 2*N, prefer N+N*NB) */
00836 /*              (RWorkspace: need 0) */
00837 
00838                 i__2 = *lwork - iwork + 1;
00839                 cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[iwork], &
00840                         i__2, &ierr);
00841 
00842 /*              Zero out below R */
00843 
00844                 i__2 = *n - 1;
00845                 i__3 = *n - 1;
00846                 claset_("L", &i__2, &i__3, &c_b1, &c_b1, &a[a_dim1 + 2], lda);
00847                 ie = 1;
00848                 itauq = 1;
00849                 itaup = itauq + *n;
00850                 iwork = itaup + *n;
00851 
00852 /*              Bidiagonalize R in A */
00853 /*              (CWorkspace: need 3*N, prefer 2*N+2*N*NB) */
00854 /*              (RWorkspace: need N) */
00855 
00856                 i__2 = *lwork - iwork + 1;
00857                 cgebrd_(n, n, &a[a_offset], lda, &s[1], &rwork[ie], &work[
00858                         itauq], &work[itaup], &work[iwork], &i__2, &ierr);
00859                 ncvt = 0;
00860                 if (wntvo || wntvas) {
00861 
00862 /*                 If right singular vectors desired, generate P'. */
00863 /*                 (CWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) */
00864 /*                 (RWorkspace: 0) */
00865 
00866                     i__2 = *lwork - iwork + 1;
00867                     cungbr_("P", n, n, n, &a[a_offset], lda, &work[itaup], &
00868                             work[iwork], &i__2, &ierr);
00869                     ncvt = *n;
00870                 }
00871                 irwork = ie + *n;
00872 
00873 /*              Perform bidiagonal QR iteration, computing right */
00874 /*              singular vectors of A in A if desired */
00875 /*              (CWorkspace: 0) */
00876 /*              (RWorkspace: need BDSPAC) */
00877 
00878                 cbdsqr_("U", n, &ncvt, &c__0, &c__0, &s[1], &rwork[ie], &a[
00879                         a_offset], lda, cdum, &c__1, cdum, &c__1, &rwork[
00880                         irwork], info);
00881 
00882 /*              If right singular vectors desired in VT, copy them there */
00883 
00884                 if (wntvas) {
00885                     clacpy_("F", n, n, &a[a_offset], lda, &vt[vt_offset], 
00886                             ldvt);
00887                 }
00888 
00889             } else if (wntuo && wntvn) {
00890 
00891 /*              Path 2 (M much larger than N, JOBU='O', JOBVT='N') */
00892 /*              N left singular vectors to be overwritten on A and */
00893 /*              no right singular vectors to be computed */
00894 
00895                 if (*lwork >= *n * *n + *n * 3) {
00896 
00897 /*                 Sufficient workspace for a fast algorithm */
00898 
00899                     ir = 1;
00900 /* Computing MAX */
00901                     i__2 = wrkbl, i__3 = *lda * *n;
00902                     if (*lwork >= max(i__2,i__3) + *lda * *n) {
00903 
00904 /*                    WORK(IU) is LDA by N, WORK(IR) is LDA by N */
00905 
00906                         ldwrku = *lda;
00907                         ldwrkr = *lda;
00908                     } else /* if(complicated condition) */ {
00909 /* Computing MAX */
00910                         i__2 = wrkbl, i__3 = *lda * *n;
00911                         if (*lwork >= max(i__2,i__3) + *n * *n) {
00912 
00913 /*                    WORK(IU) is LDA by N, WORK(IR) is N by N */
00914 
00915                             ldwrku = *lda;
00916                             ldwrkr = *n;
00917                         } else {
00918 
00919 /*                    WORK(IU) is LDWRKU by N, WORK(IR) is N by N */
00920 
00921                             ldwrku = (*lwork - *n * *n) / *n;
00922                             ldwrkr = *n;
00923                         }
00924                     }
00925                     itau = ir + ldwrkr * *n;
00926                     iwork = itau + *n;
00927 
00928 /*                 Compute A=Q*R */
00929 /*                 (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
00930 /*                 (RWorkspace: 0) */
00931 
00932                     i__2 = *lwork - iwork + 1;
00933                     cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[iwork]
00934 , &i__2, &ierr);
00935 
00936 /*                 Copy R to WORK(IR) and zero out below it */
00937 
00938                     clacpy_("U", n, n, &a[a_offset], lda, &work[ir], &ldwrkr);
00939                     i__2 = *n - 1;
00940                     i__3 = *n - 1;
00941                     claset_("L", &i__2, &i__3, &c_b1, &c_b1, &work[ir + 1], &
00942                             ldwrkr);
00943 
00944 /*                 Generate Q in A */
00945 /*                 (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
00946 /*                 (RWorkspace: 0) */
00947 
00948                     i__2 = *lwork - iwork + 1;
00949                     cungqr_(m, n, n, &a[a_offset], lda, &work[itau], &work[
00950                             iwork], &i__2, &ierr);
00951                     ie = 1;
00952                     itauq = itau;
00953                     itaup = itauq + *n;
00954                     iwork = itaup + *n;
00955 
00956 /*                 Bidiagonalize R in WORK(IR) */
00957 /*                 (CWorkspace: need N*N+3*N, prefer N*N+2*N+2*N*NB) */
00958 /*                 (RWorkspace: need N) */
00959 
00960                     i__2 = *lwork - iwork + 1;
00961                     cgebrd_(n, n, &work[ir], &ldwrkr, &s[1], &rwork[ie], &
00962                             work[itauq], &work[itaup], &work[iwork], &i__2, &
00963                             ierr);
00964 
00965 /*                 Generate left vectors bidiagonalizing R */
00966 /*                 (CWorkspace: need N*N+3*N, prefer N*N+2*N+N*NB) */
00967 /*                 (RWorkspace: need 0) */
00968 
00969                     i__2 = *lwork - iwork + 1;
00970                     cungbr_("Q", n, n, n, &work[ir], &ldwrkr, &work[itauq], &
00971                             work[iwork], &i__2, &ierr);
00972                     irwork = ie + *n;
00973 
00974 /*                 Perform bidiagonal QR iteration, computing left */
00975 /*                 singular vectors of R in WORK(IR) */
00976 /*                 (CWorkspace: need N*N) */
00977 /*                 (RWorkspace: need BDSPAC) */
00978 
00979                     cbdsqr_("U", n, &c__0, n, &c__0, &s[1], &rwork[ie], cdum, 
00980                             &c__1, &work[ir], &ldwrkr, cdum, &c__1, &rwork[
00981                             irwork], info);
00982                     iu = itauq;
00983 
00984 /*                 Multiply Q in A by left singular vectors of R in */
00985 /*                 WORK(IR), storing result in WORK(IU) and copying to A */
00986 /*                 (CWorkspace: need N*N+N, prefer N*N+M*N) */
00987 /*                 (RWorkspace: 0) */
00988 
00989                     i__2 = *m;
00990                     i__3 = ldwrku;
00991                     for (i__ = 1; i__3 < 0 ? i__ >= i__2 : i__ <= i__2; i__ +=
00992                              i__3) {
00993 /* Computing MIN */
00994                         i__4 = *m - i__ + 1;
00995                         chunk = min(i__4,ldwrku);
00996                         cgemm_("N", "N", &chunk, n, n, &c_b2, &a[i__ + a_dim1]
00997 , lda, &work[ir], &ldwrkr, &c_b1, &work[iu], &
00998                                 ldwrku);
00999                         clacpy_("F", &chunk, n, &work[iu], &ldwrku, &a[i__ + 
01000                                 a_dim1], lda);
01001 /* L10: */
01002                     }
01003 
01004                 } else {
01005 
01006 /*                 Insufficient workspace for a fast algorithm */
01007 
01008                     ie = 1;
01009                     itauq = 1;
01010                     itaup = itauq + *n;
01011                     iwork = itaup + *n;
01012 
01013 /*                 Bidiagonalize A */
01014 /*                 (CWorkspace: need 2*N+M, prefer 2*N+(M+N)*NB) */
01015 /*                 (RWorkspace: N) */
01016 
01017                     i__3 = *lwork - iwork + 1;
01018                     cgebrd_(m, n, &a[a_offset], lda, &s[1], &rwork[ie], &work[
01019                             itauq], &work[itaup], &work[iwork], &i__3, &ierr);
01020 
01021 /*                 Generate left vectors bidiagonalizing A */
01022 /*                 (CWorkspace: need 3*N, prefer 2*N+N*NB) */
01023 /*                 (RWorkspace: 0) */
01024 
01025                     i__3 = *lwork - iwork + 1;
01026                     cungbr_("Q", m, n, n, &a[a_offset], lda, &work[itauq], &
01027                             work[iwork], &i__3, &ierr);
01028                     irwork = ie + *n;
01029 
01030 /*                 Perform bidiagonal QR iteration, computing left */
01031 /*                 singular vectors of A in A */
01032 /*                 (CWorkspace: need 0) */
01033 /*                 (RWorkspace: need BDSPAC) */
01034 
01035                     cbdsqr_("U", n, &c__0, m, &c__0, &s[1], &rwork[ie], cdum, 
01036                             &c__1, &a[a_offset], lda, cdum, &c__1, &rwork[
01037                             irwork], info);
01038 
01039                 }
01040 
01041             } else if (wntuo && wntvas) {
01042 
01043 /*              Path 3 (M much larger than N, JOBU='O', JOBVT='S' or 'A') */
01044 /*              N left singular vectors to be overwritten on A and */
01045 /*              N right singular vectors to be computed in VT */
01046 
01047                 if (*lwork >= *n * *n + *n * 3) {
01048 
01049 /*                 Sufficient workspace for a fast algorithm */
01050 
01051                     ir = 1;
01052 /* Computing MAX */
01053                     i__3 = wrkbl, i__2 = *lda * *n;
01054                     if (*lwork >= max(i__3,i__2) + *lda * *n) {
01055 
01056 /*                    WORK(IU) is LDA by N and WORK(IR) is LDA by N */
01057 
01058                         ldwrku = *lda;
01059                         ldwrkr = *lda;
01060                     } else /* if(complicated condition) */ {
01061 /* Computing MAX */
01062                         i__3 = wrkbl, i__2 = *lda * *n;
01063                         if (*lwork >= max(i__3,i__2) + *n * *n) {
01064 
01065 /*                    WORK(IU) is LDA by N and WORK(IR) is N by N */
01066 
01067                             ldwrku = *lda;
01068                             ldwrkr = *n;
01069                         } else {
01070 
01071 /*                    WORK(IU) is LDWRKU by N and WORK(IR) is N by N */
01072 
01073                             ldwrku = (*lwork - *n * *n) / *n;
01074                             ldwrkr = *n;
01075                         }
01076                     }
01077                     itau = ir + ldwrkr * *n;
01078                     iwork = itau + *n;
01079 
01080 /*                 Compute A=Q*R */
01081 /*                 (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
01082 /*                 (RWorkspace: 0) */
01083 
01084                     i__3 = *lwork - iwork + 1;
01085                     cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[iwork]
01086 , &i__3, &ierr);
01087 
01088 /*                 Copy R to VT, zeroing out below it */
01089 
01090                     clacpy_("U", n, n, &a[a_offset], lda, &vt[vt_offset], 
01091                             ldvt);
01092                     if (*n > 1) {
01093                         i__3 = *n - 1;
01094                         i__2 = *n - 1;
01095                         claset_("L", &i__3, &i__2, &c_b1, &c_b1, &vt[vt_dim1 
01096                                 + 2], ldvt);
01097                     }
01098 
01099 /*                 Generate Q in A */
01100 /*                 (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
01101 /*                 (RWorkspace: 0) */
01102 
01103                     i__3 = *lwork - iwork + 1;
01104                     cungqr_(m, n, n, &a[a_offset], lda, &work[itau], &work[
01105                             iwork], &i__3, &ierr);
01106                     ie = 1;
01107                     itauq = itau;
01108                     itaup = itauq + *n;
01109                     iwork = itaup + *n;
01110 
01111 /*                 Bidiagonalize R in VT, copying result to WORK(IR) */
01112 /*                 (CWorkspace: need N*N+3*N, prefer N*N+2*N+2*N*NB) */
01113 /*                 (RWorkspace: need N) */
01114 
01115                     i__3 = *lwork - iwork + 1;
01116                     cgebrd_(n, n, &vt[vt_offset], ldvt, &s[1], &rwork[ie], &
01117                             work[itauq], &work[itaup], &work[iwork], &i__3, &
01118                             ierr);
01119                     clacpy_("L", n, n, &vt[vt_offset], ldvt, &work[ir], &
01120                             ldwrkr);
01121 
01122 /*                 Generate left vectors bidiagonalizing R in WORK(IR) */
01123 /*                 (CWorkspace: need N*N+3*N, prefer N*N+2*N+N*NB) */
01124 /*                 (RWorkspace: 0) */
01125 
01126                     i__3 = *lwork - iwork + 1;
01127                     cungbr_("Q", n, n, n, &work[ir], &ldwrkr, &work[itauq], &
01128                             work[iwork], &i__3, &ierr);
01129 
01130 /*                 Generate right vectors bidiagonalizing R in VT */
01131 /*                 (CWorkspace: need N*N+3*N-1, prefer N*N+2*N+(N-1)*NB) */
01132 /*                 (RWorkspace: 0) */
01133 
01134                     i__3 = *lwork - iwork + 1;
01135                     cungbr_("P", n, n, n, &vt[vt_offset], ldvt, &work[itaup], 
01136                             &work[iwork], &i__3, &ierr);
01137                     irwork = ie + *n;
01138 
01139 /*                 Perform bidiagonal QR iteration, computing left */
01140 /*                 singular vectors of R in WORK(IR) and computing right */
01141 /*                 singular vectors of R in VT */
01142 /*                 (CWorkspace: need N*N) */
01143 /*                 (RWorkspace: need BDSPAC) */
01144 
01145                     cbdsqr_("U", n, n, n, &c__0, &s[1], &rwork[ie], &vt[
01146                             vt_offset], ldvt, &work[ir], &ldwrkr, cdum, &c__1, 
01147                              &rwork[irwork], info);
01148                     iu = itauq;
01149 
01150 /*                 Multiply Q in A by left singular vectors of R in */
01151 /*                 WORK(IR), storing result in WORK(IU) and copying to A */
01152 /*                 (CWorkspace: need N*N+N, prefer N*N+M*N) */
01153 /*                 (RWorkspace: 0) */
01154 
01155                     i__3 = *m;
01156                     i__2 = ldwrku;
01157                     for (i__ = 1; i__2 < 0 ? i__ >= i__3 : i__ <= i__3; i__ +=
01158                              i__2) {
01159 /* Computing MIN */
01160                         i__4 = *m - i__ + 1;
01161                         chunk = min(i__4,ldwrku);
01162                         cgemm_("N", "N", &chunk, n, n, &c_b2, &a[i__ + a_dim1]
01163 , lda, &work[ir], &ldwrkr, &c_b1, &work[iu], &
01164                                 ldwrku);
01165                         clacpy_("F", &chunk, n, &work[iu], &ldwrku, &a[i__ + 
01166                                 a_dim1], lda);
01167 /* L20: */
01168                     }
01169 
01170                 } else {
01171 
01172 /*                 Insufficient workspace for a fast algorithm */
01173 
01174                     itau = 1;
01175                     iwork = itau + *n;
01176 
01177 /*                 Compute A=Q*R */
01178 /*                 (CWorkspace: need 2*N, prefer N+N*NB) */
01179 /*                 (RWorkspace: 0) */
01180 
01181                     i__2 = *lwork - iwork + 1;
01182                     cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[iwork]
01183 , &i__2, &ierr);
01184 
01185 /*                 Copy R to VT, zeroing out below it */
01186 
01187                     clacpy_("U", n, n, &a[a_offset], lda, &vt[vt_offset], 
01188                             ldvt);
01189                     if (*n > 1) {
01190                         i__2 = *n - 1;
01191                         i__3 = *n - 1;
01192                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &vt[vt_dim1 
01193                                 + 2], ldvt);
01194                     }
01195 
01196 /*                 Generate Q in A */
01197 /*                 (CWorkspace: need 2*N, prefer N+N*NB) */
01198 /*                 (RWorkspace: 0) */
01199 
01200                     i__2 = *lwork - iwork + 1;
01201                     cungqr_(m, n, n, &a[a_offset], lda, &work[itau], &work[
01202                             iwork], &i__2, &ierr);
01203                     ie = 1;
01204                     itauq = itau;
01205                     itaup = itauq + *n;
01206                     iwork = itaup + *n;
01207 
01208 /*                 Bidiagonalize R in VT */
01209 /*                 (CWorkspace: need 3*N, prefer 2*N+2*N*NB) */
01210 /*                 (RWorkspace: N) */
01211 
01212                     i__2 = *lwork - iwork + 1;
01213                     cgebrd_(n, n, &vt[vt_offset], ldvt, &s[1], &rwork[ie], &
01214                             work[itauq], &work[itaup], &work[iwork], &i__2, &
01215                             ierr);
01216 
01217 /*                 Multiply Q in A by left vectors bidiagonalizing R */
01218 /*                 (CWorkspace: need 2*N+M, prefer 2*N+M*NB) */
01219 /*                 (RWorkspace: 0) */
01220 
01221                     i__2 = *lwork - iwork + 1;
01222                     cunmbr_("Q", "R", "N", m, n, n, &vt[vt_offset], ldvt, &
01223                             work[itauq], &a[a_offset], lda, &work[iwork], &
01224                             i__2, &ierr);
01225 
01226 /*                 Generate right vectors bidiagonalizing R in VT */
01227 /*                 (CWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) */
01228 /*                 (RWorkspace: 0) */
01229 
01230                     i__2 = *lwork - iwork + 1;
01231                     cungbr_("P", n, n, n, &vt[vt_offset], ldvt, &work[itaup], 
01232                             &work[iwork], &i__2, &ierr);
01233                     irwork = ie + *n;
01234 
01235 /*                 Perform bidiagonal QR iteration, computing left */
01236 /*                 singular vectors of A in A and computing right */
01237 /*                 singular vectors of A in VT */
01238 /*                 (CWorkspace: 0) */
01239 /*                 (RWorkspace: need BDSPAC) */
01240 
01241                     cbdsqr_("U", n, n, m, &c__0, &s[1], &rwork[ie], &vt[
01242                             vt_offset], ldvt, &a[a_offset], lda, cdum, &c__1, 
01243                             &rwork[irwork], info);
01244 
01245                 }
01246 
01247             } else if (wntus) {
01248 
01249                 if (wntvn) {
01250 
01251 /*                 Path 4 (M much larger than N, JOBU='S', JOBVT='N') */
01252 /*                 N left singular vectors to be computed in U and */
01253 /*                 no right singular vectors to be computed */
01254 
01255                     if (*lwork >= *n * *n + *n * 3) {
01256 
01257 /*                    Sufficient workspace for a fast algorithm */
01258 
01259                         ir = 1;
01260                         if (*lwork >= wrkbl + *lda * *n) {
01261 
01262 /*                       WORK(IR) is LDA by N */
01263 
01264                             ldwrkr = *lda;
01265                         } else {
01266 
01267 /*                       WORK(IR) is N by N */
01268 
01269                             ldwrkr = *n;
01270                         }
01271                         itau = ir + ldwrkr * *n;
01272                         iwork = itau + *n;
01273 
01274 /*                    Compute A=Q*R */
01275 /*                    (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
01276 /*                    (RWorkspace: 0) */
01277 
01278                         i__2 = *lwork - iwork + 1;
01279                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
01280                                 iwork], &i__2, &ierr);
01281 
01282 /*                    Copy R to WORK(IR), zeroing out below it */
01283 
01284                         clacpy_("U", n, n, &a[a_offset], lda, &work[ir], &
01285                                 ldwrkr);
01286                         i__2 = *n - 1;
01287                         i__3 = *n - 1;
01288                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &work[ir + 1]
01289 , &ldwrkr);
01290 
01291 /*                    Generate Q in A */
01292 /*                    (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
01293 /*                    (RWorkspace: 0) */
01294 
01295                         i__2 = *lwork - iwork + 1;
01296                         cungqr_(m, n, n, &a[a_offset], lda, &work[itau], &
01297                                 work[iwork], &i__2, &ierr);
01298                         ie = 1;
01299                         itauq = itau;
01300                         itaup = itauq + *n;
01301                         iwork = itaup + *n;
01302 
01303 /*                    Bidiagonalize R in WORK(IR) */
01304 /*                    (CWorkspace: need N*N+3*N, prefer N*N+2*N+2*N*NB) */
01305 /*                    (RWorkspace: need N) */
01306 
01307                         i__2 = *lwork - iwork + 1;
01308                         cgebrd_(n, n, &work[ir], &ldwrkr, &s[1], &rwork[ie], &
01309                                 work[itauq], &work[itaup], &work[iwork], &
01310                                 i__2, &ierr);
01311 
01312 /*                    Generate left vectors bidiagonalizing R in WORK(IR) */
01313 /*                    (CWorkspace: need N*N+3*N, prefer N*N+2*N+N*NB) */
01314 /*                    (RWorkspace: 0) */
01315 
01316                         i__2 = *lwork - iwork + 1;
01317                         cungbr_("Q", n, n, n, &work[ir], &ldwrkr, &work[itauq]
01318 , &work[iwork], &i__2, &ierr);
01319                         irwork = ie + *n;
01320 
01321 /*                    Perform bidiagonal QR iteration, computing left */
01322 /*                    singular vectors of R in WORK(IR) */
01323 /*                    (CWorkspace: need N*N) */
01324 /*                    (RWorkspace: need BDSPAC) */
01325 
01326                         cbdsqr_("U", n, &c__0, n, &c__0, &s[1], &rwork[ie], 
01327                                 cdum, &c__1, &work[ir], &ldwrkr, cdum, &c__1, 
01328                                 &rwork[irwork], info);
01329 
01330 /*                    Multiply Q in A by left singular vectors of R in */
01331 /*                    WORK(IR), storing result in U */
01332 /*                    (CWorkspace: need N*N) */
01333 /*                    (RWorkspace: 0) */
01334 
01335                         cgemm_("N", "N", m, n, n, &c_b2, &a[a_offset], lda, &
01336                                 work[ir], &ldwrkr, &c_b1, &u[u_offset], ldu);
01337 
01338                     } else {
01339 
01340 /*                    Insufficient workspace for a fast algorithm */
01341 
01342                         itau = 1;
01343                         iwork = itau + *n;
01344 
01345 /*                    Compute A=Q*R, copying result to U */
01346 /*                    (CWorkspace: need 2*N, prefer N+N*NB) */
01347 /*                    (RWorkspace: 0) */
01348 
01349                         i__2 = *lwork - iwork + 1;
01350                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
01351                                 iwork], &i__2, &ierr);
01352                         clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], 
01353                                 ldu);
01354 
01355 /*                    Generate Q in U */
01356 /*                    (CWorkspace: need 2*N, prefer N+N*NB) */
01357 /*                    (RWorkspace: 0) */
01358 
01359                         i__2 = *lwork - iwork + 1;
01360                         cungqr_(m, n, n, &u[u_offset], ldu, &work[itau], &
01361                                 work[iwork], &i__2, &ierr);
01362                         ie = 1;
01363                         itauq = itau;
01364                         itaup = itauq + *n;
01365                         iwork = itaup + *n;
01366 
01367 /*                    Zero out below R in A */
01368 
01369                         i__2 = *n - 1;
01370                         i__3 = *n - 1;
01371                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &a[a_dim1 + 
01372                                 2], lda);
01373 
01374 /*                    Bidiagonalize R in A */
01375 /*                    (CWorkspace: need 3*N, prefer 2*N+2*N*NB) */
01376 /*                    (RWorkspace: need N) */
01377 
01378                         i__2 = *lwork - iwork + 1;
01379                         cgebrd_(n, n, &a[a_offset], lda, &s[1], &rwork[ie], &
01380                                 work[itauq], &work[itaup], &work[iwork], &
01381                                 i__2, &ierr);
01382 
01383 /*                    Multiply Q in U by left vectors bidiagonalizing R */
01384 /*                    (CWorkspace: need 2*N+M, prefer 2*N+M*NB) */
01385 /*                    (RWorkspace: 0) */
01386 
01387                         i__2 = *lwork - iwork + 1;
01388                         cunmbr_("Q", "R", "N", m, n, n, &a[a_offset], lda, &
01389                                 work[itauq], &u[u_offset], ldu, &work[iwork], 
01390                                 &i__2, &ierr)
01391                                 ;
01392                         irwork = ie + *n;
01393 
01394 /*                    Perform bidiagonal QR iteration, computing left */
01395 /*                    singular vectors of A in U */
01396 /*                    (CWorkspace: 0) */
01397 /*                    (RWorkspace: need BDSPAC) */
01398 
01399                         cbdsqr_("U", n, &c__0, m, &c__0, &s[1], &rwork[ie], 
01400                                 cdum, &c__1, &u[u_offset], ldu, cdum, &c__1, &
01401                                 rwork[irwork], info);
01402 
01403                     }
01404 
01405                 } else if (wntvo) {
01406 
01407 /*                 Path 5 (M much larger than N, JOBU='S', JOBVT='O') */
01408 /*                 N left singular vectors to be computed in U and */
01409 /*                 N right singular vectors to be overwritten on A */
01410 
01411                     if (*lwork >= (*n << 1) * *n + *n * 3) {
01412 
01413 /*                    Sufficient workspace for a fast algorithm */
01414 
01415                         iu = 1;
01416                         if (*lwork >= wrkbl + (*lda << 1) * *n) {
01417 
01418 /*                       WORK(IU) is LDA by N and WORK(IR) is LDA by N */
01419 
01420                             ldwrku = *lda;
01421                             ir = iu + ldwrku * *n;
01422                             ldwrkr = *lda;
01423                         } else if (*lwork >= wrkbl + (*lda + *n) * *n) {
01424 
01425 /*                       WORK(IU) is LDA by N and WORK(IR) is N by N */
01426 
01427                             ldwrku = *lda;
01428                             ir = iu + ldwrku * *n;
01429                             ldwrkr = *n;
01430                         } else {
01431 
01432 /*                       WORK(IU) is N by N and WORK(IR) is N by N */
01433 
01434                             ldwrku = *n;
01435                             ir = iu + ldwrku * *n;
01436                             ldwrkr = *n;
01437                         }
01438                         itau = ir + ldwrkr * *n;
01439                         iwork = itau + *n;
01440 
01441 /*                    Compute A=Q*R */
01442 /*                    (CWorkspace: need 2*N*N+2*N, prefer 2*N*N+N+N*NB) */
01443 /*                    (RWorkspace: 0) */
01444 
01445                         i__2 = *lwork - iwork + 1;
01446                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
01447                                 iwork], &i__2, &ierr);
01448 
01449 /*                    Copy R to WORK(IU), zeroing out below it */
01450 
01451                         clacpy_("U", n, n, &a[a_offset], lda, &work[iu], &
01452                                 ldwrku);
01453                         i__2 = *n - 1;
01454                         i__3 = *n - 1;
01455                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &work[iu + 1]
01456 , &ldwrku);
01457 
01458 /*                    Generate Q in A */
01459 /*                    (CWorkspace: need 2*N*N+2*N, prefer 2*N*N+N+N*NB) */
01460 /*                    (RWorkspace: 0) */
01461 
01462                         i__2 = *lwork - iwork + 1;
01463                         cungqr_(m, n, n, &a[a_offset], lda, &work[itau], &
01464                                 work[iwork], &i__2, &ierr);
01465                         ie = 1;
01466                         itauq = itau;
01467                         itaup = itauq + *n;
01468                         iwork = itaup + *n;
01469 
01470 /*                    Bidiagonalize R in WORK(IU), copying result to */
01471 /*                    WORK(IR) */
01472 /*                    (CWorkspace: need   2*N*N+3*N, */
01473 /*                                 prefer 2*N*N+2*N+2*N*NB) */
01474 /*                    (RWorkspace: need   N) */
01475 
01476                         i__2 = *lwork - iwork + 1;
01477                         cgebrd_(n, n, &work[iu], &ldwrku, &s[1], &rwork[ie], &
01478                                 work[itauq], &work[itaup], &work[iwork], &
01479                                 i__2, &ierr);
01480                         clacpy_("U", n, n, &work[iu], &ldwrku, &work[ir], &
01481                                 ldwrkr);
01482 
01483 /*                    Generate left bidiagonalizing vectors in WORK(IU) */
01484 /*                    (CWorkspace: need 2*N*N+3*N, prefer 2*N*N+2*N+N*NB) */
01485 /*                    (RWorkspace: 0) */
01486 
01487                         i__2 = *lwork - iwork + 1;
01488                         cungbr_("Q", n, n, n, &work[iu], &ldwrku, &work[itauq]
01489 , &work[iwork], &i__2, &ierr);
01490 
01491 /*                    Generate right bidiagonalizing vectors in WORK(IR) */
01492 /*                    (CWorkspace: need   2*N*N+3*N-1, */
01493 /*                                 prefer 2*N*N+2*N+(N-1)*NB) */
01494 /*                    (RWorkspace: 0) */
01495 
01496                         i__2 = *lwork - iwork + 1;
01497                         cungbr_("P", n, n, n, &work[ir], &ldwrkr, &work[itaup]
01498 , &work[iwork], &i__2, &ierr);
01499                         irwork = ie + *n;
01500 
01501 /*                    Perform bidiagonal QR iteration, computing left */
01502 /*                    singular vectors of R in WORK(IU) and computing */
01503 /*                    right singular vectors of R in WORK(IR) */
01504 /*                    (CWorkspace: need 2*N*N) */
01505 /*                    (RWorkspace: need BDSPAC) */
01506 
01507                         cbdsqr_("U", n, n, n, &c__0, &s[1], &rwork[ie], &work[
01508                                 ir], &ldwrkr, &work[iu], &ldwrku, cdum, &c__1, 
01509                                  &rwork[irwork], info);
01510 
01511 /*                    Multiply Q in A by left singular vectors of R in */
01512 /*                    WORK(IU), storing result in U */
01513 /*                    (CWorkspace: need N*N) */
01514 /*                    (RWorkspace: 0) */
01515 
01516                         cgemm_("N", "N", m, n, n, &c_b2, &a[a_offset], lda, &
01517                                 work[iu], &ldwrku, &c_b1, &u[u_offset], ldu);
01518 
01519 /*                    Copy right singular vectors of R to A */
01520 /*                    (CWorkspace: need N*N) */
01521 /*                    (RWorkspace: 0) */
01522 
01523                         clacpy_("F", n, n, &work[ir], &ldwrkr, &a[a_offset], 
01524                                 lda);
01525 
01526                     } else {
01527 
01528 /*                    Insufficient workspace for a fast algorithm */
01529 
01530                         itau = 1;
01531                         iwork = itau + *n;
01532 
01533 /*                    Compute A=Q*R, copying result to U */
01534 /*                    (CWorkspace: need 2*N, prefer N+N*NB) */
01535 /*                    (RWorkspace: 0) */
01536 
01537                         i__2 = *lwork - iwork + 1;
01538                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
01539                                 iwork], &i__2, &ierr);
01540                         clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], 
01541                                 ldu);
01542 
01543 /*                    Generate Q in U */
01544 /*                    (CWorkspace: need 2*N, prefer N+N*NB) */
01545 /*                    (RWorkspace: 0) */
01546 
01547                         i__2 = *lwork - iwork + 1;
01548                         cungqr_(m, n, n, &u[u_offset], ldu, &work[itau], &
01549                                 work[iwork], &i__2, &ierr);
01550                         ie = 1;
01551                         itauq = itau;
01552                         itaup = itauq + *n;
01553                         iwork = itaup + *n;
01554 
01555 /*                    Zero out below R in A */
01556 
01557                         i__2 = *n - 1;
01558                         i__3 = *n - 1;
01559                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &a[a_dim1 + 
01560                                 2], lda);
01561 
01562 /*                    Bidiagonalize R in A */
01563 /*                    (CWorkspace: need 3*N, prefer 2*N+2*N*NB) */
01564 /*                    (RWorkspace: need N) */
01565 
01566                         i__2 = *lwork - iwork + 1;
01567                         cgebrd_(n, n, &a[a_offset], lda, &s[1], &rwork[ie], &
01568                                 work[itauq], &work[itaup], &work[iwork], &
01569                                 i__2, &ierr);
01570 
01571 /*                    Multiply Q in U by left vectors bidiagonalizing R */
01572 /*                    (CWorkspace: need 2*N+M, prefer 2*N+M*NB) */
01573 /*                    (RWorkspace: 0) */
01574 
01575                         i__2 = *lwork - iwork + 1;
01576                         cunmbr_("Q", "R", "N", m, n, n, &a[a_offset], lda, &
01577                                 work[itauq], &u[u_offset], ldu, &work[iwork], 
01578                                 &i__2, &ierr)
01579                                 ;
01580 
01581 /*                    Generate right vectors bidiagonalizing R in A */
01582 /*                    (CWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) */
01583 /*                    (RWorkspace: 0) */
01584 
01585                         i__2 = *lwork - iwork + 1;
01586                         cungbr_("P", n, n, n, &a[a_offset], lda, &work[itaup], 
01587                                  &work[iwork], &i__2, &ierr);
01588                         irwork = ie + *n;
01589 
01590 /*                    Perform bidiagonal QR iteration, computing left */
01591 /*                    singular vectors of A in U and computing right */
01592 /*                    singular vectors of A in A */
01593 /*                    (CWorkspace: 0) */
01594 /*                    (RWorkspace: need BDSPAC) */
01595 
01596                         cbdsqr_("U", n, n, m, &c__0, &s[1], &rwork[ie], &a[
01597                                 a_offset], lda, &u[u_offset], ldu, cdum, &
01598                                 c__1, &rwork[irwork], info);
01599 
01600                     }
01601 
01602                 } else if (wntvas) {
01603 
01604 /*                 Path 6 (M much larger than N, JOBU='S', JOBVT='S' */
01605 /*                         or 'A') */
01606 /*                 N left singular vectors to be computed in U and */
01607 /*                 N right singular vectors to be computed in VT */
01608 
01609                     if (*lwork >= *n * *n + *n * 3) {
01610 
01611 /*                    Sufficient workspace for a fast algorithm */
01612 
01613                         iu = 1;
01614                         if (*lwork >= wrkbl + *lda * *n) {
01615 
01616 /*                       WORK(IU) is LDA by N */
01617 
01618                             ldwrku = *lda;
01619                         } else {
01620 
01621 /*                       WORK(IU) is N by N */
01622 
01623                             ldwrku = *n;
01624                         }
01625                         itau = iu + ldwrku * *n;
01626                         iwork = itau + *n;
01627 
01628 /*                    Compute A=Q*R */
01629 /*                    (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
01630 /*                    (RWorkspace: 0) */
01631 
01632                         i__2 = *lwork - iwork + 1;
01633                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
01634                                 iwork], &i__2, &ierr);
01635 
01636 /*                    Copy R to WORK(IU), zeroing out below it */
01637 
01638                         clacpy_("U", n, n, &a[a_offset], lda, &work[iu], &
01639                                 ldwrku);
01640                         i__2 = *n - 1;
01641                         i__3 = *n - 1;
01642                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &work[iu + 1]
01643 , &ldwrku);
01644 
01645 /*                    Generate Q in A */
01646 /*                    (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
01647 /*                    (RWorkspace: 0) */
01648 
01649                         i__2 = *lwork - iwork + 1;
01650                         cungqr_(m, n, n, &a[a_offset], lda, &work[itau], &
01651                                 work[iwork], &i__2, &ierr);
01652                         ie = 1;
01653                         itauq = itau;
01654                         itaup = itauq + *n;
01655                         iwork = itaup + *n;
01656 
01657 /*                    Bidiagonalize R in WORK(IU), copying result to VT */
01658 /*                    (CWorkspace: need N*N+3*N, prefer N*N+2*N+2*N*NB) */
01659 /*                    (RWorkspace: need N) */
01660 
01661                         i__2 = *lwork - iwork + 1;
01662                         cgebrd_(n, n, &work[iu], &ldwrku, &s[1], &rwork[ie], &
01663                                 work[itauq], &work[itaup], &work[iwork], &
01664                                 i__2, &ierr);
01665                         clacpy_("U", n, n, &work[iu], &ldwrku, &vt[vt_offset], 
01666                                  ldvt);
01667 
01668 /*                    Generate left bidiagonalizing vectors in WORK(IU) */
01669 /*                    (CWorkspace: need N*N+3*N, prefer N*N+2*N+N*NB) */
01670 /*                    (RWorkspace: 0) */
01671 
01672                         i__2 = *lwork - iwork + 1;
01673                         cungbr_("Q", n, n, n, &work[iu], &ldwrku, &work[itauq]
01674 , &work[iwork], &i__2, &ierr);
01675 
01676 /*                    Generate right bidiagonalizing vectors in VT */
01677 /*                    (CWorkspace: need   N*N+3*N-1, */
01678 /*                                 prefer N*N+2*N+(N-1)*NB) */
01679 /*                    (RWorkspace: 0) */
01680 
01681                         i__2 = *lwork - iwork + 1;
01682                         cungbr_("P", n, n, n, &vt[vt_offset], ldvt, &work[
01683                                 itaup], &work[iwork], &i__2, &ierr)
01684                                 ;
01685                         irwork = ie + *n;
01686 
01687 /*                    Perform bidiagonal QR iteration, computing left */
01688 /*                    singular vectors of R in WORK(IU) and computing */
01689 /*                    right singular vectors of R in VT */
01690 /*                    (CWorkspace: need N*N) */
01691 /*                    (RWorkspace: need BDSPAC) */
01692 
01693                         cbdsqr_("U", n, n, n, &c__0, &s[1], &rwork[ie], &vt[
01694                                 vt_offset], ldvt, &work[iu], &ldwrku, cdum, &
01695                                 c__1, &rwork[irwork], info);
01696 
01697 /*                    Multiply Q in A by left singular vectors of R in */
01698 /*                    WORK(IU), storing result in U */
01699 /*                    (CWorkspace: need N*N) */
01700 /*                    (RWorkspace: 0) */
01701 
01702                         cgemm_("N", "N", m, n, n, &c_b2, &a[a_offset], lda, &
01703                                 work[iu], &ldwrku, &c_b1, &u[u_offset], ldu);
01704 
01705                     } else {
01706 
01707 /*                    Insufficient workspace for a fast algorithm */
01708 
01709                         itau = 1;
01710                         iwork = itau + *n;
01711 
01712 /*                    Compute A=Q*R, copying result to U */
01713 /*                    (CWorkspace: need 2*N, prefer N+N*NB) */
01714 /*                    (RWorkspace: 0) */
01715 
01716                         i__2 = *lwork - iwork + 1;
01717                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
01718                                 iwork], &i__2, &ierr);
01719                         clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], 
01720                                 ldu);
01721 
01722 /*                    Generate Q in U */
01723 /*                    (CWorkspace: need 2*N, prefer N+N*NB) */
01724 /*                    (RWorkspace: 0) */
01725 
01726                         i__2 = *lwork - iwork + 1;
01727                         cungqr_(m, n, n, &u[u_offset], ldu, &work[itau], &
01728                                 work[iwork], &i__2, &ierr);
01729 
01730 /*                    Copy R to VT, zeroing out below it */
01731 
01732                         clacpy_("U", n, n, &a[a_offset], lda, &vt[vt_offset], 
01733                                 ldvt);
01734                         if (*n > 1) {
01735                             i__2 = *n - 1;
01736                             i__3 = *n - 1;
01737                             claset_("L", &i__2, &i__3, &c_b1, &c_b1, &vt[
01738                                     vt_dim1 + 2], ldvt);
01739                         }
01740                         ie = 1;
01741                         itauq = itau;
01742                         itaup = itauq + *n;
01743                         iwork = itaup + *n;
01744 
01745 /*                    Bidiagonalize R in VT */
01746 /*                    (CWorkspace: need 3*N, prefer 2*N+2*N*NB) */
01747 /*                    (RWorkspace: need N) */
01748 
01749                         i__2 = *lwork - iwork + 1;
01750                         cgebrd_(n, n, &vt[vt_offset], ldvt, &s[1], &rwork[ie], 
01751                                  &work[itauq], &work[itaup], &work[iwork], &
01752                                 i__2, &ierr);
01753 
01754 /*                    Multiply Q in U by left bidiagonalizing vectors */
01755 /*                    in VT */
01756 /*                    (CWorkspace: need 2*N+M, prefer 2*N+M*NB) */
01757 /*                    (RWorkspace: 0) */
01758 
01759                         i__2 = *lwork - iwork + 1;
01760                         cunmbr_("Q", "R", "N", m, n, n, &vt[vt_offset], ldvt, 
01761                                 &work[itauq], &u[u_offset], ldu, &work[iwork], 
01762                                  &i__2, &ierr);
01763 
01764 /*                    Generate right bidiagonalizing vectors in VT */
01765 /*                    (CWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) */
01766 /*                    (RWorkspace: 0) */
01767 
01768                         i__2 = *lwork - iwork + 1;
01769                         cungbr_("P", n, n, n, &vt[vt_offset], ldvt, &work[
01770                                 itaup], &work[iwork], &i__2, &ierr)
01771                                 ;
01772                         irwork = ie + *n;
01773 
01774 /*                    Perform bidiagonal QR iteration, computing left */
01775 /*                    singular vectors of A in U and computing right */
01776 /*                    singular vectors of A in VT */
01777 /*                    (CWorkspace: 0) */
01778 /*                    (RWorkspace: need BDSPAC) */
01779 
01780                         cbdsqr_("U", n, n, m, &c__0, &s[1], &rwork[ie], &vt[
01781                                 vt_offset], ldvt, &u[u_offset], ldu, cdum, &
01782                                 c__1, &rwork[irwork], info);
01783 
01784                     }
01785 
01786                 }
01787 
01788             } else if (wntua) {
01789 
01790                 if (wntvn) {
01791 
01792 /*                 Path 7 (M much larger than N, JOBU='A', JOBVT='N') */
01793 /*                 M left singular vectors to be computed in U and */
01794 /*                 no right singular vectors to be computed */
01795 
01796 /* Computing MAX */
01797                     i__2 = *n + *m, i__3 = *n * 3;
01798                     if (*lwork >= *n * *n + max(i__2,i__3)) {
01799 
01800 /*                    Sufficient workspace for a fast algorithm */
01801 
01802                         ir = 1;
01803                         if (*lwork >= wrkbl + *lda * *n) {
01804 
01805 /*                       WORK(IR) is LDA by N */
01806 
01807                             ldwrkr = *lda;
01808                         } else {
01809 
01810 /*                       WORK(IR) is N by N */
01811 
01812                             ldwrkr = *n;
01813                         }
01814                         itau = ir + ldwrkr * *n;
01815                         iwork = itau + *n;
01816 
01817 /*                    Compute A=Q*R, copying result to U */
01818 /*                    (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
01819 /*                    (RWorkspace: 0) */
01820 
01821                         i__2 = *lwork - iwork + 1;
01822                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
01823                                 iwork], &i__2, &ierr);
01824                         clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], 
01825                                 ldu);
01826 
01827 /*                    Copy R to WORK(IR), zeroing out below it */
01828 
01829                         clacpy_("U", n, n, &a[a_offset], lda, &work[ir], &
01830                                 ldwrkr);
01831                         i__2 = *n - 1;
01832                         i__3 = *n - 1;
01833                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &work[ir + 1]
01834 , &ldwrkr);
01835 
01836 /*                    Generate Q in U */
01837 /*                    (CWorkspace: need N*N+N+M, prefer N*N+N+M*NB) */
01838 /*                    (RWorkspace: 0) */
01839 
01840                         i__2 = *lwork - iwork + 1;
01841                         cungqr_(m, m, n, &u[u_offset], ldu, &work[itau], &
01842                                 work[iwork], &i__2, &ierr);
01843                         ie = 1;
01844                         itauq = itau;
01845                         itaup = itauq + *n;
01846                         iwork = itaup + *n;
01847 
01848 /*                    Bidiagonalize R in WORK(IR) */
01849 /*                    (CWorkspace: need N*N+3*N, prefer N*N+2*N+2*N*NB) */
01850 /*                    (RWorkspace: need N) */
01851 
01852                         i__2 = *lwork - iwork + 1;
01853                         cgebrd_(n, n, &work[ir], &ldwrkr, &s[1], &rwork[ie], &
01854                                 work[itauq], &work[itaup], &work[iwork], &
01855                                 i__2, &ierr);
01856 
01857 /*                    Generate left bidiagonalizing vectors in WORK(IR) */
01858 /*                    (CWorkspace: need N*N+3*N, prefer N*N+2*N+N*NB) */
01859 /*                    (RWorkspace: 0) */
01860 
01861                         i__2 = *lwork - iwork + 1;
01862                         cungbr_("Q", n, n, n, &work[ir], &ldwrkr, &work[itauq]
01863 , &work[iwork], &i__2, &ierr);
01864                         irwork = ie + *n;
01865 
01866 /*                    Perform bidiagonal QR iteration, computing left */
01867 /*                    singular vectors of R in WORK(IR) */
01868 /*                    (CWorkspace: need N*N) */
01869 /*                    (RWorkspace: need BDSPAC) */
01870 
01871                         cbdsqr_("U", n, &c__0, n, &c__0, &s[1], &rwork[ie], 
01872                                 cdum, &c__1, &work[ir], &ldwrkr, cdum, &c__1, 
01873                                 &rwork[irwork], info);
01874 
01875 /*                    Multiply Q in U by left singular vectors of R in */
01876 /*                    WORK(IR), storing result in A */
01877 /*                    (CWorkspace: need N*N) */
01878 /*                    (RWorkspace: 0) */
01879 
01880                         cgemm_("N", "N", m, n, n, &c_b2, &u[u_offset], ldu, &
01881                                 work[ir], &ldwrkr, &c_b1, &a[a_offset], lda);
01882 
01883 /*                    Copy left singular vectors of A from A to U */
01884 
01885                         clacpy_("F", m, n, &a[a_offset], lda, &u[u_offset], 
01886                                 ldu);
01887 
01888                     } else {
01889 
01890 /*                    Insufficient workspace for a fast algorithm */
01891 
01892                         itau = 1;
01893                         iwork = itau + *n;
01894 
01895 /*                    Compute A=Q*R, copying result to U */
01896 /*                    (CWorkspace: need 2*N, prefer N+N*NB) */
01897 /*                    (RWorkspace: 0) */
01898 
01899                         i__2 = *lwork - iwork + 1;
01900                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
01901                                 iwork], &i__2, &ierr);
01902                         clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], 
01903                                 ldu);
01904 
01905 /*                    Generate Q in U */
01906 /*                    (CWorkspace: need N+M, prefer N+M*NB) */
01907 /*                    (RWorkspace: 0) */
01908 
01909                         i__2 = *lwork - iwork + 1;
01910                         cungqr_(m, m, n, &u[u_offset], ldu, &work[itau], &
01911                                 work[iwork], &i__2, &ierr);
01912                         ie = 1;
01913                         itauq = itau;
01914                         itaup = itauq + *n;
01915                         iwork = itaup + *n;
01916 
01917 /*                    Zero out below R in A */
01918 
01919                         i__2 = *n - 1;
01920                         i__3 = *n - 1;
01921                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &a[a_dim1 + 
01922                                 2], lda);
01923 
01924 /*                    Bidiagonalize R in A */
01925 /*                    (CWorkspace: need 3*N, prefer 2*N+2*N*NB) */
01926 /*                    (RWorkspace: need N) */
01927 
01928                         i__2 = *lwork - iwork + 1;
01929                         cgebrd_(n, n, &a[a_offset], lda, &s[1], &rwork[ie], &
01930                                 work[itauq], &work[itaup], &work[iwork], &
01931                                 i__2, &ierr);
01932 
01933 /*                    Multiply Q in U by left bidiagonalizing vectors */
01934 /*                    in A */
01935 /*                    (CWorkspace: need 2*N+M, prefer 2*N+M*NB) */
01936 /*                    (RWorkspace: 0) */
01937 
01938                         i__2 = *lwork - iwork + 1;
01939                         cunmbr_("Q", "R", "N", m, n, n, &a[a_offset], lda, &
01940                                 work[itauq], &u[u_offset], ldu, &work[iwork], 
01941                                 &i__2, &ierr)
01942                                 ;
01943                         irwork = ie + *n;
01944 
01945 /*                    Perform bidiagonal QR iteration, computing left */
01946 /*                    singular vectors of A in U */
01947 /*                    (CWorkspace: 0) */
01948 /*                    (RWorkspace: need BDSPAC) */
01949 
01950                         cbdsqr_("U", n, &c__0, m, &c__0, &s[1], &rwork[ie], 
01951                                 cdum, &c__1, &u[u_offset], ldu, cdum, &c__1, &
01952                                 rwork[irwork], info);
01953 
01954                     }
01955 
01956                 } else if (wntvo) {
01957 
01958 /*                 Path 8 (M much larger than N, JOBU='A', JOBVT='O') */
01959 /*                 M left singular vectors to be computed in U and */
01960 /*                 N right singular vectors to be overwritten on A */
01961 
01962 /* Computing MAX */
01963                     i__2 = *n + *m, i__3 = *n * 3;
01964                     if (*lwork >= (*n << 1) * *n + max(i__2,i__3)) {
01965 
01966 /*                    Sufficient workspace for a fast algorithm */
01967 
01968                         iu = 1;
01969                         if (*lwork >= wrkbl + (*lda << 1) * *n) {
01970 
01971 /*                       WORK(IU) is LDA by N and WORK(IR) is LDA by N */
01972 
01973                             ldwrku = *lda;
01974                             ir = iu + ldwrku * *n;
01975                             ldwrkr = *lda;
01976                         } else if (*lwork >= wrkbl + (*lda + *n) * *n) {
01977 
01978 /*                       WORK(IU) is LDA by N and WORK(IR) is N by N */
01979 
01980                             ldwrku = *lda;
01981                             ir = iu + ldwrku * *n;
01982                             ldwrkr = *n;
01983                         } else {
01984 
01985 /*                       WORK(IU) is N by N and WORK(IR) is N by N */
01986 
01987                             ldwrku = *n;
01988                             ir = iu + ldwrku * *n;
01989                             ldwrkr = *n;
01990                         }
01991                         itau = ir + ldwrkr * *n;
01992                         iwork = itau + *n;
01993 
01994 /*                    Compute A=Q*R, copying result to U */
01995 /*                    (CWorkspace: need 2*N*N+2*N, prefer 2*N*N+N+N*NB) */
01996 /*                    (RWorkspace: 0) */
01997 
01998                         i__2 = *lwork - iwork + 1;
01999                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
02000                                 iwork], &i__2, &ierr);
02001                         clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], 
02002                                 ldu);
02003 
02004 /*                    Generate Q in U */
02005 /*                    (CWorkspace: need 2*N*N+N+M, prefer 2*N*N+N+M*NB) */
02006 /*                    (RWorkspace: 0) */
02007 
02008                         i__2 = *lwork - iwork + 1;
02009                         cungqr_(m, m, n, &u[u_offset], ldu, &work[itau], &
02010                                 work[iwork], &i__2, &ierr);
02011 
02012 /*                    Copy R to WORK(IU), zeroing out below it */
02013 
02014                         clacpy_("U", n, n, &a[a_offset], lda, &work[iu], &
02015                                 ldwrku);
02016                         i__2 = *n - 1;
02017                         i__3 = *n - 1;
02018                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &work[iu + 1]
02019 , &ldwrku);
02020                         ie = 1;
02021                         itauq = itau;
02022                         itaup = itauq + *n;
02023                         iwork = itaup + *n;
02024 
02025 /*                    Bidiagonalize R in WORK(IU), copying result to */
02026 /*                    WORK(IR) */
02027 /*                    (CWorkspace: need   2*N*N+3*N, */
02028 /*                                 prefer 2*N*N+2*N+2*N*NB) */
02029 /*                    (RWorkspace: need   N) */
02030 
02031                         i__2 = *lwork - iwork + 1;
02032                         cgebrd_(n, n, &work[iu], &ldwrku, &s[1], &rwork[ie], &
02033                                 work[itauq], &work[itaup], &work[iwork], &
02034                                 i__2, &ierr);
02035                         clacpy_("U", n, n, &work[iu], &ldwrku, &work[ir], &
02036                                 ldwrkr);
02037 
02038 /*                    Generate left bidiagonalizing vectors in WORK(IU) */
02039 /*                    (CWorkspace: need 2*N*N+3*N, prefer 2*N*N+2*N+N*NB) */
02040 /*                    (RWorkspace: 0) */
02041 
02042                         i__2 = *lwork - iwork + 1;
02043                         cungbr_("Q", n, n, n, &work[iu], &ldwrku, &work[itauq]
02044 , &work[iwork], &i__2, &ierr);
02045 
02046 /*                    Generate right bidiagonalizing vectors in WORK(IR) */
02047 /*                    (CWorkspace: need   2*N*N+3*N-1, */
02048 /*                                 prefer 2*N*N+2*N+(N-1)*NB) */
02049 /*                    (RWorkspace: 0) */
02050 
02051                         i__2 = *lwork - iwork + 1;
02052                         cungbr_("P", n, n, n, &work[ir], &ldwrkr, &work[itaup]
02053 , &work[iwork], &i__2, &ierr);
02054                         irwork = ie + *n;
02055 
02056 /*                    Perform bidiagonal QR iteration, computing left */
02057 /*                    singular vectors of R in WORK(IU) and computing */
02058 /*                    right singular vectors of R in WORK(IR) */
02059 /*                    (CWorkspace: need 2*N*N) */
02060 /*                    (RWorkspace: need BDSPAC) */
02061 
02062                         cbdsqr_("U", n, n, n, &c__0, &s[1], &rwork[ie], &work[
02063                                 ir], &ldwrkr, &work[iu], &ldwrku, cdum, &c__1, 
02064                                  &rwork[irwork], info);
02065 
02066 /*                    Multiply Q in U by left singular vectors of R in */
02067 /*                    WORK(IU), storing result in A */
02068 /*                    (CWorkspace: need N*N) */
02069 /*                    (RWorkspace: 0) */
02070 
02071                         cgemm_("N", "N", m, n, n, &c_b2, &u[u_offset], ldu, &
02072                                 work[iu], &ldwrku, &c_b1, &a[a_offset], lda);
02073 
02074 /*                    Copy left singular vectors of A from A to U */
02075 
02076                         clacpy_("F", m, n, &a[a_offset], lda, &u[u_offset], 
02077                                 ldu);
02078 
02079 /*                    Copy right singular vectors of R from WORK(IR) to A */
02080 
02081                         clacpy_("F", n, n, &work[ir], &ldwrkr, &a[a_offset], 
02082                                 lda);
02083 
02084                     } else {
02085 
02086 /*                    Insufficient workspace for a fast algorithm */
02087 
02088                         itau = 1;
02089                         iwork = itau + *n;
02090 
02091 /*                    Compute A=Q*R, copying result to U */
02092 /*                    (CWorkspace: need 2*N, prefer N+N*NB) */
02093 /*                    (RWorkspace: 0) */
02094 
02095                         i__2 = *lwork - iwork + 1;
02096                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
02097                                 iwork], &i__2, &ierr);
02098                         clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], 
02099                                 ldu);
02100 
02101 /*                    Generate Q in U */
02102 /*                    (CWorkspace: need N+M, prefer N+M*NB) */
02103 /*                    (RWorkspace: 0) */
02104 
02105                         i__2 = *lwork - iwork + 1;
02106                         cungqr_(m, m, n, &u[u_offset], ldu, &work[itau], &
02107                                 work[iwork], &i__2, &ierr);
02108                         ie = 1;
02109                         itauq = itau;
02110                         itaup = itauq + *n;
02111                         iwork = itaup + *n;
02112 
02113 /*                    Zero out below R in A */
02114 
02115                         i__2 = *n - 1;
02116                         i__3 = *n - 1;
02117                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &a[a_dim1 + 
02118                                 2], lda);
02119 
02120 /*                    Bidiagonalize R in A */
02121 /*                    (CWorkspace: need 3*N, prefer 2*N+2*N*NB) */
02122 /*                    (RWorkspace: need N) */
02123 
02124                         i__2 = *lwork - iwork + 1;
02125                         cgebrd_(n, n, &a[a_offset], lda, &s[1], &rwork[ie], &
02126                                 work[itauq], &work[itaup], &work[iwork], &
02127                                 i__2, &ierr);
02128 
02129 /*                    Multiply Q in U by left bidiagonalizing vectors */
02130 /*                    in A */
02131 /*                    (CWorkspace: need 2*N+M, prefer 2*N+M*NB) */
02132 /*                    (RWorkspace: 0) */
02133 
02134                         i__2 = *lwork - iwork + 1;
02135                         cunmbr_("Q", "R", "N", m, n, n, &a[a_offset], lda, &
02136                                 work[itauq], &u[u_offset], ldu, &work[iwork], 
02137                                 &i__2, &ierr)
02138                                 ;
02139 
02140 /*                    Generate right bidiagonalizing vectors in A */
02141 /*                    (CWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) */
02142 /*                    (RWorkspace: 0) */
02143 
02144                         i__2 = *lwork - iwork + 1;
02145                         cungbr_("P", n, n, n, &a[a_offset], lda, &work[itaup], 
02146                                  &work[iwork], &i__2, &ierr);
02147                         irwork = ie + *n;
02148 
02149 /*                    Perform bidiagonal QR iteration, computing left */
02150 /*                    singular vectors of A in U and computing right */
02151 /*                    singular vectors of A in A */
02152 /*                    (CWorkspace: 0) */
02153 /*                    (RWorkspace: need BDSPAC) */
02154 
02155                         cbdsqr_("U", n, n, m, &c__0, &s[1], &rwork[ie], &a[
02156                                 a_offset], lda, &u[u_offset], ldu, cdum, &
02157                                 c__1, &rwork[irwork], info);
02158 
02159                     }
02160 
02161                 } else if (wntvas) {
02162 
02163 /*                 Path 9 (M much larger than N, JOBU='A', JOBVT='S' */
02164 /*                         or 'A') */
02165 /*                 M left singular vectors to be computed in U and */
02166 /*                 N right singular vectors to be computed in VT */
02167 
02168 /* Computing MAX */
02169                     i__2 = *n + *m, i__3 = *n * 3;
02170                     if (*lwork >= *n * *n + max(i__2,i__3)) {
02171 
02172 /*                    Sufficient workspace for a fast algorithm */
02173 
02174                         iu = 1;
02175                         if (*lwork >= wrkbl + *lda * *n) {
02176 
02177 /*                       WORK(IU) is LDA by N */
02178 
02179                             ldwrku = *lda;
02180                         } else {
02181 
02182 /*                       WORK(IU) is N by N */
02183 
02184                             ldwrku = *n;
02185                         }
02186                         itau = iu + ldwrku * *n;
02187                         iwork = itau + *n;
02188 
02189 /*                    Compute A=Q*R, copying result to U */
02190 /*                    (CWorkspace: need N*N+2*N, prefer N*N+N+N*NB) */
02191 /*                    (RWorkspace: 0) */
02192 
02193                         i__2 = *lwork - iwork + 1;
02194                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
02195                                 iwork], &i__2, &ierr);
02196                         clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], 
02197                                 ldu);
02198 
02199 /*                    Generate Q in U */
02200 /*                    (CWorkspace: need N*N+N+M, prefer N*N+N+M*NB) */
02201 /*                    (RWorkspace: 0) */
02202 
02203                         i__2 = *lwork - iwork + 1;
02204                         cungqr_(m, m, n, &u[u_offset], ldu, &work[itau], &
02205                                 work[iwork], &i__2, &ierr);
02206 
02207 /*                    Copy R to WORK(IU), zeroing out below it */
02208 
02209                         clacpy_("U", n, n, &a[a_offset], lda, &work[iu], &
02210                                 ldwrku);
02211                         i__2 = *n - 1;
02212                         i__3 = *n - 1;
02213                         claset_("L", &i__2, &i__3, &c_b1, &c_b1, &work[iu + 1]
02214 , &ldwrku);
02215                         ie = 1;
02216                         itauq = itau;
02217                         itaup = itauq + *n;
02218                         iwork = itaup + *n;
02219 
02220 /*                    Bidiagonalize R in WORK(IU), copying result to VT */
02221 /*                    (CWorkspace: need N*N+3*N, prefer N*N+2*N+2*N*NB) */
02222 /*                    (RWorkspace: need N) */
02223 
02224                         i__2 = *lwork - iwork + 1;
02225                         cgebrd_(n, n, &work[iu], &ldwrku, &s[1], &rwork[ie], &
02226                                 work[itauq], &work[itaup], &work[iwork], &
02227                                 i__2, &ierr);
02228                         clacpy_("U", n, n, &work[iu], &ldwrku, &vt[vt_offset], 
02229                                  ldvt);
02230 
02231 /*                    Generate left bidiagonalizing vectors in WORK(IU) */
02232 /*                    (CWorkspace: need N*N+3*N, prefer N*N+2*N+N*NB) */
02233 /*                    (RWorkspace: 0) */
02234 
02235                         i__2 = *lwork - iwork + 1;
02236                         cungbr_("Q", n, n, n, &work[iu], &ldwrku, &work[itauq]
02237 , &work[iwork], &i__2, &ierr);
02238 
02239 /*                    Generate right bidiagonalizing vectors in VT */
02240 /*                    (CWorkspace: need   N*N+3*N-1, */
02241 /*                                 prefer N*N+2*N+(N-1)*NB) */
02242 /*                    (RWorkspace: need   0) */
02243 
02244                         i__2 = *lwork - iwork + 1;
02245                         cungbr_("P", n, n, n, &vt[vt_offset], ldvt, &work[
02246                                 itaup], &work[iwork], &i__2, &ierr)
02247                                 ;
02248                         irwork = ie + *n;
02249 
02250 /*                    Perform bidiagonal QR iteration, computing left */
02251 /*                    singular vectors of R in WORK(IU) and computing */
02252 /*                    right singular vectors of R in VT */
02253 /*                    (CWorkspace: need N*N) */
02254 /*                    (RWorkspace: need BDSPAC) */
02255 
02256                         cbdsqr_("U", n, n, n, &c__0, &s[1], &rwork[ie], &vt[
02257                                 vt_offset], ldvt, &work[iu], &ldwrku, cdum, &
02258                                 c__1, &rwork[irwork], info);
02259 
02260 /*                    Multiply Q in U by left singular vectors of R in */
02261 /*                    WORK(IU), storing result in A */
02262 /*                    (CWorkspace: need N*N) */
02263 /*                    (RWorkspace: 0) */
02264 
02265                         cgemm_("N", "N", m, n, n, &c_b2, &u[u_offset], ldu, &
02266                                 work[iu], &ldwrku, &c_b1, &a[a_offset], lda);
02267 
02268 /*                    Copy left singular vectors of A from A to U */
02269 
02270                         clacpy_("F", m, n, &a[a_offset], lda, &u[u_offset], 
02271                                 ldu);
02272 
02273                     } else {
02274 
02275 /*                    Insufficient workspace for a fast algorithm */
02276 
02277                         itau = 1;
02278                         iwork = itau + *n;
02279 
02280 /*                    Compute A=Q*R, copying result to U */
02281 /*                    (CWorkspace: need 2*N, prefer N+N*NB) */
02282 /*                    (RWorkspace: 0) */
02283 
02284                         i__2 = *lwork - iwork + 1;
02285                         cgeqrf_(m, n, &a[a_offset], lda, &work[itau], &work[
02286                                 iwork], &i__2, &ierr);
02287                         clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], 
02288                                 ldu);
02289 
02290 /*                    Generate Q in U */
02291 /*                    (CWorkspace: need N+M, prefer N+M*NB) */
02292 /*                    (RWorkspace: 0) */
02293 
02294                         i__2 = *lwork - iwork + 1;
02295                         cungqr_(m, m, n, &u[u_offset], ldu, &work[itau], &
02296                                 work[iwork], &i__2, &ierr);
02297 
02298 /*                    Copy R from A to VT, zeroing out below it */
02299 
02300                         clacpy_("U", n, n, &a[a_offset], lda, &vt[vt_offset], 
02301                                 ldvt);
02302                         if (*n > 1) {
02303                             i__2 = *n - 1;
02304                             i__3 = *n - 1;
02305                             claset_("L", &i__2, &i__3, &c_b1, &c_b1, &vt[
02306                                     vt_dim1 + 2], ldvt);
02307                         }
02308                         ie = 1;
02309                         itauq = itau;
02310                         itaup = itauq + *n;
02311                         iwork = itaup + *n;
02312 
02313 /*                    Bidiagonalize R in VT */
02314 /*                    (CWorkspace: need 3*N, prefer 2*N+2*N*NB) */
02315 /*                    (RWorkspace: need N) */
02316 
02317                         i__2 = *lwork - iwork + 1;
02318                         cgebrd_(n, n, &vt[vt_offset], ldvt, &s[1], &rwork[ie], 
02319                                  &work[itauq], &work[itaup], &work[iwork], &
02320                                 i__2, &ierr);
02321 
02322 /*                    Multiply Q in U by left bidiagonalizing vectors */
02323 /*                    in VT */
02324 /*                    (CWorkspace: need 2*N+M, prefer 2*N+M*NB) */
02325 /*                    (RWorkspace: 0) */
02326 
02327                         i__2 = *lwork - iwork + 1;
02328                         cunmbr_("Q", "R", "N", m, n, n, &vt[vt_offset], ldvt, 
02329                                 &work[itauq], &u[u_offset], ldu, &work[iwork], 
02330                                  &i__2, &ierr);
02331 
02332 /*                    Generate right bidiagonalizing vectors in VT */
02333 /*                    (CWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) */
02334 /*                    (RWorkspace: 0) */
02335 
02336                         i__2 = *lwork - iwork + 1;
02337                         cungbr_("P", n, n, n, &vt[vt_offset], ldvt, &work[
02338                                 itaup], &work[iwork], &i__2, &ierr)
02339                                 ;
02340                         irwork = ie + *n;
02341 
02342 /*                    Perform bidiagonal QR iteration, computing left */
02343 /*                    singular vectors of A in U and computing right */
02344 /*                    singular vectors of A in VT */
02345 /*                    (CWorkspace: 0) */
02346 /*                    (RWorkspace: need BDSPAC) */
02347 
02348                         cbdsqr_("U", n, n, m, &c__0, &s[1], &rwork[ie], &vt[
02349                                 vt_offset], ldvt, &u[u_offset], ldu, cdum, &
02350                                 c__1, &rwork[irwork], info);
02351 
02352                     }
02353 
02354                 }
02355 
02356             }
02357 
02358         } else {
02359 
02360 /*           M .LT. MNTHR */
02361 
02362 /*           Path 10 (M at least N, but not much larger) */
02363 /*           Reduce to bidiagonal form without QR decomposition */
02364 
02365             ie = 1;
02366             itauq = 1;
02367             itaup = itauq + *n;
02368             iwork = itaup + *n;
02369 
02370 /*           Bidiagonalize A */
02371 /*           (CWorkspace: need 2*N+M, prefer 2*N+(M+N)*NB) */
02372 /*           (RWorkspace: need N) */
02373 
02374             i__2 = *lwork - iwork + 1;
02375             cgebrd_(m, n, &a[a_offset], lda, &s[1], &rwork[ie], &work[itauq], 
02376                     &work[itaup], &work[iwork], &i__2, &ierr);
02377             if (wntuas) {
02378 
02379 /*              If left singular vectors desired in U, copy result to U */
02380 /*              and generate left bidiagonalizing vectors in U */
02381 /*              (CWorkspace: need 2*N+NCU, prefer 2*N+NCU*NB) */
02382 /*              (RWorkspace: 0) */
02383 
02384                 clacpy_("L", m, n, &a[a_offset], lda, &u[u_offset], ldu);
02385                 if (wntus) {
02386                     ncu = *n;
02387                 }
02388                 if (wntua) {
02389                     ncu = *m;
02390                 }
02391                 i__2 = *lwork - iwork + 1;
02392                 cungbr_("Q", m, &ncu, n, &u[u_offset], ldu, &work[itauq], &
02393                         work[iwork], &i__2, &ierr);
02394             }
02395             if (wntvas) {
02396 
02397 /*              If right singular vectors desired in VT, copy result to */
02398 /*              VT and generate right bidiagonalizing vectors in VT */
02399 /*              (CWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) */
02400 /*              (RWorkspace: 0) */
02401 
02402                 clacpy_("U", n, n, &a[a_offset], lda, &vt[vt_offset], ldvt);
02403                 i__2 = *lwork - iwork + 1;
02404                 cungbr_("P", n, n, n, &vt[vt_offset], ldvt, &work[itaup], &
02405                         work[iwork], &i__2, &ierr);
02406             }
02407             if (wntuo) {
02408 
02409 /*              If left singular vectors desired in A, generate left */
02410 /*              bidiagonalizing vectors in A */
02411 /*              (CWorkspace: need 3*N, prefer 2*N+N*NB) */
02412 /*              (RWorkspace: 0) */
02413 
02414                 i__2 = *lwork - iwork + 1;
02415                 cungbr_("Q", m, n, n, &a[a_offset], lda, &work[itauq], &work[
02416                         iwork], &i__2, &ierr);
02417             }
02418             if (wntvo) {
02419 
02420 /*              If right singular vectors desired in A, generate right */
02421 /*              bidiagonalizing vectors in A */
02422 /*              (CWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) */
02423 /*              (RWorkspace: 0) */
02424 
02425                 i__2 = *lwork - iwork + 1;
02426                 cungbr_("P", n, n, n, &a[a_offset], lda, &work[itaup], &work[
02427                         iwork], &i__2, &ierr);
02428             }
02429             irwork = ie + *n;
02430             if (wntuas || wntuo) {
02431                 nru = *m;
02432             }
02433             if (wntun) {
02434                 nru = 0;
02435             }
02436             if (wntvas || wntvo) {
02437                 ncvt = *n;
02438             }
02439             if (wntvn) {
02440                 ncvt = 0;
02441             }
02442             if (! wntuo && ! wntvo) {
02443 
02444 /*              Perform bidiagonal QR iteration, if desired, computing */
02445 /*              left singular vectors in U and computing right singular */
02446 /*              vectors in VT */
02447 /*              (CWorkspace: 0) */
02448 /*              (RWorkspace: need BDSPAC) */
02449 
02450                 cbdsqr_("U", n, &ncvt, &nru, &c__0, &s[1], &rwork[ie], &vt[
02451                         vt_offset], ldvt, &u[u_offset], ldu, cdum, &c__1, &
02452                         rwork[irwork], info);
02453             } else if (! wntuo && wntvo) {
02454 
02455 /*              Perform bidiagonal QR iteration, if desired, computing */
02456 /*              left singular vectors in U and computing right singular */
02457 /*              vectors in A */
02458 /*              (CWorkspace: 0) */
02459 /*              (RWorkspace: need BDSPAC) */
02460 
02461                 cbdsqr_("U", n, &ncvt, &nru, &c__0, &s[1], &rwork[ie], &a[
02462                         a_offset], lda, &u[u_offset], ldu, cdum, &c__1, &
02463                         rwork[irwork], info);
02464             } else {
02465 
02466 /*              Perform bidiagonal QR iteration, if desired, computing */
02467 /*              left singular vectors in A and computing right singular */
02468 /*              vectors in VT */
02469 /*              (CWorkspace: 0) */
02470 /*              (RWorkspace: need BDSPAC) */
02471 
02472                 cbdsqr_("U", n, &ncvt, &nru, &c__0, &s[1], &rwork[ie], &vt[
02473                         vt_offset], ldvt, &a[a_offset], lda, cdum, &c__1, &
02474                         rwork[irwork], info);
02475             }
02476 
02477         }
02478 
02479     } else {
02480 
02481 /*        A has more columns than rows. If A has sufficiently more */
02482 /*        columns than rows, first reduce using the LQ decomposition (if */
02483 /*        sufficient workspace available) */
02484 
02485         if (*n >= mnthr) {
02486 
02487             if (wntvn) {
02488 
02489 /*              Path 1t(N much larger than M, JOBVT='N') */
02490 /*              No right singular vectors to be computed */
02491 
02492                 itau = 1;
02493                 iwork = itau + *m;
02494 
02495 /*              Compute A=L*Q */
02496 /*              (CWorkspace: need 2*M, prefer M+M*NB) */
02497 /*              (RWorkspace: 0) */
02498 
02499                 i__2 = *lwork - iwork + 1;
02500                 cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[iwork], &
02501                         i__2, &ierr);
02502 
02503 /*              Zero out above L */
02504 
02505                 i__2 = *m - 1;
02506                 i__3 = *m - 1;
02507                 claset_("U", &i__2, &i__3, &c_b1, &c_b1, &a[(a_dim1 << 1) + 1]
02508 , lda);
02509                 ie = 1;
02510                 itauq = 1;
02511                 itaup = itauq + *m;
02512                 iwork = itaup + *m;
02513 
02514 /*              Bidiagonalize L in A */
02515 /*              (CWorkspace: need 3*M, prefer 2*M+2*M*NB) */
02516 /*              (RWorkspace: need M) */
02517 
02518                 i__2 = *lwork - iwork + 1;
02519                 cgebrd_(m, m, &a[a_offset], lda, &s[1], &rwork[ie], &work[
02520                         itauq], &work[itaup], &work[iwork], &i__2, &ierr);
02521                 if (wntuo || wntuas) {
02522 
02523 /*                 If left singular vectors desired, generate Q */
02524 /*                 (CWorkspace: need 3*M, prefer 2*M+M*NB) */
02525 /*                 (RWorkspace: 0) */
02526 
02527                     i__2 = *lwork - iwork + 1;
02528                     cungbr_("Q", m, m, m, &a[a_offset], lda, &work[itauq], &
02529                             work[iwork], &i__2, &ierr);
02530                 }
02531                 irwork = ie + *m;
02532                 nru = 0;
02533                 if (wntuo || wntuas) {
02534                     nru = *m;
02535                 }
02536 
02537 /*              Perform bidiagonal QR iteration, computing left singular */
02538 /*              vectors of A in A if desired */
02539 /*              (CWorkspace: 0) */
02540 /*              (RWorkspace: need BDSPAC) */
02541 
02542                 cbdsqr_("U", m, &c__0, &nru, &c__0, &s[1], &rwork[ie], cdum, &
02543                         c__1, &a[a_offset], lda, cdum, &c__1, &rwork[irwork], 
02544                         info);
02545 
02546 /*              If left singular vectors desired in U, copy them there */
02547 
02548                 if (wntuas) {
02549                     clacpy_("F", m, m, &a[a_offset], lda, &u[u_offset], ldu);
02550                 }
02551 
02552             } else if (wntvo && wntun) {
02553 
02554 /*              Path 2t(N much larger than M, JOBU='N', JOBVT='O') */
02555 /*              M right singular vectors to be overwritten on A and */
02556 /*              no left singular vectors to be computed */
02557 
02558                 if (*lwork >= *m * *m + *m * 3) {
02559 
02560 /*                 Sufficient workspace for a fast algorithm */
02561 
02562                     ir = 1;
02563 /* Computing MAX */
02564                     i__2 = wrkbl, i__3 = *lda * *n;
02565                     if (*lwork >= max(i__2,i__3) + *lda * *m) {
02566 
02567 /*                    WORK(IU) is LDA by N and WORK(IR) is LDA by M */
02568 
02569                         ldwrku = *lda;
02570                         chunk = *n;
02571                         ldwrkr = *lda;
02572                     } else /* if(complicated condition) */ {
02573 /* Computing MAX */
02574                         i__2 = wrkbl, i__3 = *lda * *n;
02575                         if (*lwork >= max(i__2,i__3) + *m * *m) {
02576 
02577 /*                    WORK(IU) is LDA by N and WORK(IR) is M by M */
02578 
02579                             ldwrku = *lda;
02580                             chunk = *n;
02581                             ldwrkr = *m;
02582                         } else {
02583 
02584 /*                    WORK(IU) is M by CHUNK and WORK(IR) is M by M */
02585 
02586                             ldwrku = *m;
02587                             chunk = (*lwork - *m * *m) / *m;
02588                             ldwrkr = *m;
02589                         }
02590                     }
02591                     itau = ir + ldwrkr * *m;
02592                     iwork = itau + *m;
02593 
02594 /*                 Compute A=L*Q */
02595 /*                 (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
02596 /*                 (RWorkspace: 0) */
02597 
02598                     i__2 = *lwork - iwork + 1;
02599                     cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[iwork]
02600 , &i__2, &ierr);
02601 
02602 /*                 Copy L to WORK(IR) and zero out above it */
02603 
02604                     clacpy_("L", m, m, &a[a_offset], lda, &work[ir], &ldwrkr);
02605                     i__2 = *m - 1;
02606                     i__3 = *m - 1;
02607                     claset_("U", &i__2, &i__3, &c_b1, &c_b1, &work[ir + 
02608                             ldwrkr], &ldwrkr);
02609 
02610 /*                 Generate Q in A */
02611 /*                 (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
02612 /*                 (RWorkspace: 0) */
02613 
02614                     i__2 = *lwork - iwork + 1;
02615                     cunglq_(m, n, m, &a[a_offset], lda, &work[itau], &work[
02616                             iwork], &i__2, &ierr);
02617                     ie = 1;
02618                     itauq = itau;
02619                     itaup = itauq + *m;
02620                     iwork = itaup + *m;
02621 
02622 /*                 Bidiagonalize L in WORK(IR) */
02623 /*                 (CWorkspace: need M*M+3*M, prefer M*M+2*M+2*M*NB) */
02624 /*                 (RWorkspace: need M) */
02625 
02626                     i__2 = *lwork - iwork + 1;
02627                     cgebrd_(m, m, &work[ir], &ldwrkr, &s[1], &rwork[ie], &
02628                             work[itauq], &work[itaup], &work[iwork], &i__2, &
02629                             ierr);
02630 
02631 /*                 Generate right vectors bidiagonalizing L */
02632 /*                 (CWorkspace: need M*M+3*M-1, prefer M*M+2*M+(M-1)*NB) */
02633 /*                 (RWorkspace: 0) */
02634 
02635                     i__2 = *lwork - iwork + 1;
02636                     cungbr_("P", m, m, m, &work[ir], &ldwrkr, &work[itaup], &
02637                             work[iwork], &i__2, &ierr);
02638                     irwork = ie + *m;
02639 
02640 /*                 Perform bidiagonal QR iteration, computing right */
02641 /*                 singular vectors of L in WORK(IR) */
02642 /*                 (CWorkspace: need M*M) */
02643 /*                 (RWorkspace: need BDSPAC) */
02644 
02645                     cbdsqr_("U", m, m, &c__0, &c__0, &s[1], &rwork[ie], &work[
02646                             ir], &ldwrkr, cdum, &c__1, cdum, &c__1, &rwork[
02647                             irwork], info);
02648                     iu = itauq;
02649 
02650 /*                 Multiply right singular vectors of L in WORK(IR) by Q */
02651 /*                 in A, storing result in WORK(IU) and copying to A */
02652 /*                 (CWorkspace: need M*M+M, prefer M*M+M*N) */
02653 /*                 (RWorkspace: 0) */
02654 
02655                     i__2 = *n;
02656                     i__3 = chunk;
02657                     for (i__ = 1; i__3 < 0 ? i__ >= i__2 : i__ <= i__2; i__ +=
02658                              i__3) {
02659 /* Computing MIN */
02660                         i__4 = *n - i__ + 1;
02661                         blk = min(i__4,chunk);
02662                         cgemm_("N", "N", m, &blk, m, &c_b2, &work[ir], &
02663                                 ldwrkr, &a[i__ * a_dim1 + 1], lda, &c_b1, &
02664                                 work[iu], &ldwrku);
02665                         clacpy_("F", m, &blk, &work[iu], &ldwrku, &a[i__ * 
02666                                 a_dim1 + 1], lda);
02667 /* L30: */
02668                     }
02669 
02670                 } else {
02671 
02672 /*                 Insufficient workspace for a fast algorithm */
02673 
02674                     ie = 1;
02675                     itauq = 1;
02676                     itaup = itauq + *m;
02677                     iwork = itaup + *m;
02678 
02679 /*                 Bidiagonalize A */
02680 /*                 (CWorkspace: need 2*M+N, prefer 2*M+(M+N)*NB) */
02681 /*                 (RWorkspace: need M) */
02682 
02683                     i__3 = *lwork - iwork + 1;
02684                     cgebrd_(m, n, &a[a_offset], lda, &s[1], &rwork[ie], &work[
02685                             itauq], &work[itaup], &work[iwork], &i__3, &ierr);
02686 
02687 /*                 Generate right vectors bidiagonalizing A */
02688 /*                 (CWorkspace: need 3*M, prefer 2*M+M*NB) */
02689 /*                 (RWorkspace: 0) */
02690 
02691                     i__3 = *lwork - iwork + 1;
02692                     cungbr_("P", m, n, m, &a[a_offset], lda, &work[itaup], &
02693                             work[iwork], &i__3, &ierr);
02694                     irwork = ie + *m;
02695 
02696 /*                 Perform bidiagonal QR iteration, computing right */
02697 /*                 singular vectors of A in A */
02698 /*                 (CWorkspace: 0) */
02699 /*                 (RWorkspace: need BDSPAC) */
02700 
02701                     cbdsqr_("L", m, n, &c__0, &c__0, &s[1], &rwork[ie], &a[
02702                             a_offset], lda, cdum, &c__1, cdum, &c__1, &rwork[
02703                             irwork], info);
02704 
02705                 }
02706 
02707             } else if (wntvo && wntuas) {
02708 
02709 /*              Path 3t(N much larger than M, JOBU='S' or 'A', JOBVT='O') */
02710 /*              M right singular vectors to be overwritten on A and */
02711 /*              M left singular vectors to be computed in U */
02712 
02713                 if (*lwork >= *m * *m + *m * 3) {
02714 
02715 /*                 Sufficient workspace for a fast algorithm */
02716 
02717                     ir = 1;
02718 /* Computing MAX */
02719                     i__3 = wrkbl, i__2 = *lda * *n;
02720                     if (*lwork >= max(i__3,i__2) + *lda * *m) {
02721 
02722 /*                    WORK(IU) is LDA by N and WORK(IR) is LDA by M */
02723 
02724                         ldwrku = *lda;
02725                         chunk = *n;
02726                         ldwrkr = *lda;
02727                     } else /* if(complicated condition) */ {
02728 /* Computing MAX */
02729                         i__3 = wrkbl, i__2 = *lda * *n;
02730                         if (*lwork >= max(i__3,i__2) + *m * *m) {
02731 
02732 /*                    WORK(IU) is LDA by N and WORK(IR) is M by M */
02733 
02734                             ldwrku = *lda;
02735                             chunk = *n;
02736                             ldwrkr = *m;
02737                         } else {
02738 
02739 /*                    WORK(IU) is M by CHUNK and WORK(IR) is M by M */
02740 
02741                             ldwrku = *m;
02742                             chunk = (*lwork - *m * *m) / *m;
02743                             ldwrkr = *m;
02744                         }
02745                     }
02746                     itau = ir + ldwrkr * *m;
02747                     iwork = itau + *m;
02748 
02749 /*                 Compute A=L*Q */
02750 /*                 (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
02751 /*                 (RWorkspace: 0) */
02752 
02753                     i__3 = *lwork - iwork + 1;
02754                     cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[iwork]
02755 , &i__3, &ierr);
02756 
02757 /*                 Copy L to U, zeroing about above it */
02758 
02759                     clacpy_("L", m, m, &a[a_offset], lda, &u[u_offset], ldu);
02760                     i__3 = *m - 1;
02761                     i__2 = *m - 1;
02762                     claset_("U", &i__3, &i__2, &c_b1, &c_b1, &u[(u_dim1 << 1) 
02763                             + 1], ldu);
02764 
02765 /*                 Generate Q in A */
02766 /*                 (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
02767 /*                 (RWorkspace: 0) */
02768 
02769                     i__3 = *lwork - iwork + 1;
02770                     cunglq_(m, n, m, &a[a_offset], lda, &work[itau], &work[
02771                             iwork], &i__3, &ierr);
02772                     ie = 1;
02773                     itauq = itau;
02774                     itaup = itauq + *m;
02775                     iwork = itaup + *m;
02776 
02777 /*                 Bidiagonalize L in U, copying result to WORK(IR) */
02778 /*                 (CWorkspace: need M*M+3*M, prefer M*M+2*M+2*M*NB) */
02779 /*                 (RWorkspace: need M) */
02780 
02781                     i__3 = *lwork - iwork + 1;
02782                     cgebrd_(m, m, &u[u_offset], ldu, &s[1], &rwork[ie], &work[
02783                             itauq], &work[itaup], &work[iwork], &i__3, &ierr);
02784                     clacpy_("U", m, m, &u[u_offset], ldu, &work[ir], &ldwrkr);
02785 
02786 /*                 Generate right vectors bidiagonalizing L in WORK(IR) */
02787 /*                 (CWorkspace: need M*M+3*M-1, prefer M*M+2*M+(M-1)*NB) */
02788 /*                 (RWorkspace: 0) */
02789 
02790                     i__3 = *lwork - iwork + 1;
02791                     cungbr_("P", m, m, m, &work[ir], &ldwrkr, &work[itaup], &
02792                             work[iwork], &i__3, &ierr);
02793 
02794 /*                 Generate left vectors bidiagonalizing L in U */
02795 /*                 (CWorkspace: need M*M+3*M, prefer M*M+2*M+M*NB) */
02796 /*                 (RWorkspace: 0) */
02797 
02798                     i__3 = *lwork - iwork + 1;
02799                     cungbr_("Q", m, m, m, &u[u_offset], ldu, &work[itauq], &
02800                             work[iwork], &i__3, &ierr);
02801                     irwork = ie + *m;
02802 
02803 /*                 Perform bidiagonal QR iteration, computing left */
02804 /*                 singular vectors of L in U, and computing right */
02805 /*                 singular vectors of L in WORK(IR) */
02806 /*                 (CWorkspace: need M*M) */
02807 /*                 (RWorkspace: need BDSPAC) */
02808 
02809                     cbdsqr_("U", m, m, m, &c__0, &s[1], &rwork[ie], &work[ir], 
02810                              &ldwrkr, &u[u_offset], ldu, cdum, &c__1, &rwork[
02811                             irwork], info);
02812                     iu = itauq;
02813 
02814 /*                 Multiply right singular vectors of L in WORK(IR) by Q */
02815 /*                 in A, storing result in WORK(IU) and copying to A */
02816 /*                 (CWorkspace: need M*M+M, prefer M*M+M*N)) */
02817 /*                 (RWorkspace: 0) */
02818 
02819                     i__3 = *n;
02820                     i__2 = chunk;
02821                     for (i__ = 1; i__2 < 0 ? i__ >= i__3 : i__ <= i__3; i__ +=
02822                              i__2) {
02823 /* Computing MIN */
02824                         i__4 = *n - i__ + 1;
02825                         blk = min(i__4,chunk);
02826                         cgemm_("N", "N", m, &blk, m, &c_b2, &work[ir], &
02827                                 ldwrkr, &a[i__ * a_dim1 + 1], lda, &c_b1, &
02828                                 work[iu], &ldwrku);
02829                         clacpy_("F", m, &blk, &work[iu], &ldwrku, &a[i__ * 
02830                                 a_dim1 + 1], lda);
02831 /* L40: */
02832                     }
02833 
02834                 } else {
02835 
02836 /*                 Insufficient workspace for a fast algorithm */
02837 
02838                     itau = 1;
02839                     iwork = itau + *m;
02840 
02841 /*                 Compute A=L*Q */
02842 /*                 (CWorkspace: need 2*M, prefer M+M*NB) */
02843 /*                 (RWorkspace: 0) */
02844 
02845                     i__2 = *lwork - iwork + 1;
02846                     cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[iwork]
02847 , &i__2, &ierr);
02848 
02849 /*                 Copy L to U, zeroing out above it */
02850 
02851                     clacpy_("L", m, m, &a[a_offset], lda, &u[u_offset], ldu);
02852                     i__2 = *m - 1;
02853                     i__3 = *m - 1;
02854                     claset_("U", &i__2, &i__3, &c_b1, &c_b1, &u[(u_dim1 << 1) 
02855                             + 1], ldu);
02856 
02857 /*                 Generate Q in A */
02858 /*                 (CWorkspace: need 2*M, prefer M+M*NB) */
02859 /*                 (RWorkspace: 0) */
02860 
02861                     i__2 = *lwork - iwork + 1;
02862                     cunglq_(m, n, m, &a[a_offset], lda, &work[itau], &work[
02863                             iwork], &i__2, &ierr);
02864                     ie = 1;
02865                     itauq = itau;
02866                     itaup = itauq + *m;
02867                     iwork = itaup + *m;
02868 
02869 /*                 Bidiagonalize L in U */
02870 /*                 (CWorkspace: need 3*M, prefer 2*M+2*M*NB) */
02871 /*                 (RWorkspace: need M) */
02872 
02873                     i__2 = *lwork - iwork + 1;
02874                     cgebrd_(m, m, &u[u_offset], ldu, &s[1], &rwork[ie], &work[
02875                             itauq], &work[itaup], &work[iwork], &i__2, &ierr);
02876 
02877 /*                 Multiply right vectors bidiagonalizing L by Q in A */
02878 /*                 (CWorkspace: need 2*M+N, prefer 2*M+N*NB) */
02879 /*                 (RWorkspace: 0) */
02880 
02881                     i__2 = *lwork - iwork + 1;
02882                     cunmbr_("P", "L", "C", m, n, m, &u[u_offset], ldu, &work[
02883                             itaup], &a[a_offset], lda, &work[iwork], &i__2, &
02884                             ierr);
02885 
02886 /*                 Generate left vectors bidiagonalizing L in U */
02887 /*                 (CWorkspace: need 3*M, prefer 2*M+M*NB) */
02888 /*                 (RWorkspace: 0) */
02889 
02890                     i__2 = *lwork - iwork + 1;
02891                     cungbr_("Q", m, m, m, &u[u_offset], ldu, &work[itauq], &
02892                             work[iwork], &i__2, &ierr);
02893                     irwork = ie + *m;
02894 
02895 /*                 Perform bidiagonal QR iteration, computing left */
02896 /*                 singular vectors of A in U and computing right */
02897 /*                 singular vectors of A in A */
02898 /*                 (CWorkspace: 0) */
02899 /*                 (RWorkspace: need BDSPAC) */
02900 
02901                     cbdsqr_("U", m, n, m, &c__0, &s[1], &rwork[ie], &a[
02902                             a_offset], lda, &u[u_offset], ldu, cdum, &c__1, &
02903                             rwork[irwork], info);
02904 
02905                 }
02906 
02907             } else if (wntvs) {
02908 
02909                 if (wntun) {
02910 
02911 /*                 Path 4t(N much larger than M, JOBU='N', JOBVT='S') */
02912 /*                 M right singular vectors to be computed in VT and */
02913 /*                 no left singular vectors to be computed */
02914 
02915                     if (*lwork >= *m * *m + *m * 3) {
02916 
02917 /*                    Sufficient workspace for a fast algorithm */
02918 
02919                         ir = 1;
02920                         if (*lwork >= wrkbl + *lda * *m) {
02921 
02922 /*                       WORK(IR) is LDA by M */
02923 
02924                             ldwrkr = *lda;
02925                         } else {
02926 
02927 /*                       WORK(IR) is M by M */
02928 
02929                             ldwrkr = *m;
02930                         }
02931                         itau = ir + ldwrkr * *m;
02932                         iwork = itau + *m;
02933 
02934 /*                    Compute A=L*Q */
02935 /*                    (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
02936 /*                    (RWorkspace: 0) */
02937 
02938                         i__2 = *lwork - iwork + 1;
02939                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
02940                                 iwork], &i__2, &ierr);
02941 
02942 /*                    Copy L to WORK(IR), zeroing out above it */
02943 
02944                         clacpy_("L", m, m, &a[a_offset], lda, &work[ir], &
02945                                 ldwrkr);
02946                         i__2 = *m - 1;
02947                         i__3 = *m - 1;
02948                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &work[ir + 
02949                                 ldwrkr], &ldwrkr);
02950 
02951 /*                    Generate Q in A */
02952 /*                    (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
02953 /*                    (RWorkspace: 0) */
02954 
02955                         i__2 = *lwork - iwork + 1;
02956                         cunglq_(m, n, m, &a[a_offset], lda, &work[itau], &
02957                                 work[iwork], &i__2, &ierr);
02958                         ie = 1;
02959                         itauq = itau;
02960                         itaup = itauq + *m;
02961                         iwork = itaup + *m;
02962 
02963 /*                    Bidiagonalize L in WORK(IR) */
02964 /*                    (CWorkspace: need M*M+3*M, prefer M*M+2*M+2*M*NB) */
02965 /*                    (RWorkspace: need M) */
02966 
02967                         i__2 = *lwork - iwork + 1;
02968                         cgebrd_(m, m, &work[ir], &ldwrkr, &s[1], &rwork[ie], &
02969                                 work[itauq], &work[itaup], &work[iwork], &
02970                                 i__2, &ierr);
02971 
02972 /*                    Generate right vectors bidiagonalizing L in */
02973 /*                    WORK(IR) */
02974 /*                    (CWorkspace: need M*M+3*M, prefer M*M+2*M+(M-1)*NB) */
02975 /*                    (RWorkspace: 0) */
02976 
02977                         i__2 = *lwork - iwork + 1;
02978                         cungbr_("P", m, m, m, &work[ir], &ldwrkr, &work[itaup]
02979 , &work[iwork], &i__2, &ierr);
02980                         irwork = ie + *m;
02981 
02982 /*                    Perform bidiagonal QR iteration, computing right */
02983 /*                    singular vectors of L in WORK(IR) */
02984 /*                    (CWorkspace: need M*M) */
02985 /*                    (RWorkspace: need BDSPAC) */
02986 
02987                         cbdsqr_("U", m, m, &c__0, &c__0, &s[1], &rwork[ie], &
02988                                 work[ir], &ldwrkr, cdum, &c__1, cdum, &c__1, &
02989                                 rwork[irwork], info);
02990 
02991 /*                    Multiply right singular vectors of L in WORK(IR) by */
02992 /*                    Q in A, storing result in VT */
02993 /*                    (CWorkspace: need M*M) */
02994 /*                    (RWorkspace: 0) */
02995 
02996                         cgemm_("N", "N", m, n, m, &c_b2, &work[ir], &ldwrkr, &
02997                                 a[a_offset], lda, &c_b1, &vt[vt_offset], ldvt);
02998 
02999                     } else {
03000 
03001 /*                    Insufficient workspace for a fast algorithm */
03002 
03003                         itau = 1;
03004                         iwork = itau + *m;
03005 
03006 /*                    Compute A=L*Q */
03007 /*                    (CWorkspace: need 2*M, prefer M+M*NB) */
03008 /*                    (RWorkspace: 0) */
03009 
03010                         i__2 = *lwork - iwork + 1;
03011                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03012                                 iwork], &i__2, &ierr);
03013 
03014 /*                    Copy result to VT */
03015 
03016                         clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], 
03017                                 ldvt);
03018 
03019 /*                    Generate Q in VT */
03020 /*                    (CWorkspace: need 2*M, prefer M+M*NB) */
03021 /*                    (RWorkspace: 0) */
03022 
03023                         i__2 = *lwork - iwork + 1;
03024                         cunglq_(m, n, m, &vt[vt_offset], ldvt, &work[itau], &
03025                                 work[iwork], &i__2, &ierr);
03026                         ie = 1;
03027                         itauq = itau;
03028                         itaup = itauq + *m;
03029                         iwork = itaup + *m;
03030 
03031 /*                    Zero out above L in A */
03032 
03033                         i__2 = *m - 1;
03034                         i__3 = *m - 1;
03035                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &a[(a_dim1 <<
03036                                  1) + 1], lda);
03037 
03038 /*                    Bidiagonalize L in A */
03039 /*                    (CWorkspace: need 3*M, prefer 2*M+2*M*NB) */
03040 /*                    (RWorkspace: need M) */
03041 
03042                         i__2 = *lwork - iwork + 1;
03043                         cgebrd_(m, m, &a[a_offset], lda, &s[1], &rwork[ie], &
03044                                 work[itauq], &work[itaup], &work[iwork], &
03045                                 i__2, &ierr);
03046 
03047 /*                    Multiply right vectors bidiagonalizing L by Q in VT */
03048 /*                    (CWorkspace: need 2*M+N, prefer 2*M+N*NB) */
03049 /*                    (RWorkspace: 0) */
03050 
03051                         i__2 = *lwork - iwork + 1;
03052                         cunmbr_("P", "L", "C", m, n, m, &a[a_offset], lda, &
03053                                 work[itaup], &vt[vt_offset], ldvt, &work[
03054                                 iwork], &i__2, &ierr);
03055                         irwork = ie + *m;
03056 
03057 /*                    Perform bidiagonal QR iteration, computing right */
03058 /*                    singular vectors of A in VT */
03059 /*                    (CWorkspace: 0) */
03060 /*                    (RWorkspace: need BDSPAC) */
03061 
03062                         cbdsqr_("U", m, n, &c__0, &c__0, &s[1], &rwork[ie], &
03063                                 vt[vt_offset], ldvt, cdum, &c__1, cdum, &c__1, 
03064                                  &rwork[irwork], info);
03065 
03066                     }
03067 
03068                 } else if (wntuo) {
03069 
03070 /*                 Path 5t(N much larger than M, JOBU='O', JOBVT='S') */
03071 /*                 M right singular vectors to be computed in VT and */
03072 /*                 M left singular vectors to be overwritten on A */
03073 
03074                     if (*lwork >= (*m << 1) * *m + *m * 3) {
03075 
03076 /*                    Sufficient workspace for a fast algorithm */
03077 
03078                         iu = 1;
03079                         if (*lwork >= wrkbl + (*lda << 1) * *m) {
03080 
03081 /*                       WORK(IU) is LDA by M and WORK(IR) is LDA by M */
03082 
03083                             ldwrku = *lda;
03084                             ir = iu + ldwrku * *m;
03085                             ldwrkr = *lda;
03086                         } else if (*lwork >= wrkbl + (*lda + *m) * *m) {
03087 
03088 /*                       WORK(IU) is LDA by M and WORK(IR) is M by M */
03089 
03090                             ldwrku = *lda;
03091                             ir = iu + ldwrku * *m;
03092                             ldwrkr = *m;
03093                         } else {
03094 
03095 /*                       WORK(IU) is M by M and WORK(IR) is M by M */
03096 
03097                             ldwrku = *m;
03098                             ir = iu + ldwrku * *m;
03099                             ldwrkr = *m;
03100                         }
03101                         itau = ir + ldwrkr * *m;
03102                         iwork = itau + *m;
03103 
03104 /*                    Compute A=L*Q */
03105 /*                    (CWorkspace: need 2*M*M+2*M, prefer 2*M*M+M+M*NB) */
03106 /*                    (RWorkspace: 0) */
03107 
03108                         i__2 = *lwork - iwork + 1;
03109                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03110                                 iwork], &i__2, &ierr);
03111 
03112 /*                    Copy L to WORK(IU), zeroing out below it */
03113 
03114                         clacpy_("L", m, m, &a[a_offset], lda, &work[iu], &
03115                                 ldwrku);
03116                         i__2 = *m - 1;
03117                         i__3 = *m - 1;
03118                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &work[iu + 
03119                                 ldwrku], &ldwrku);
03120 
03121 /*                    Generate Q in A */
03122 /*                    (CWorkspace: need 2*M*M+2*M, prefer 2*M*M+M+M*NB) */
03123 /*                    (RWorkspace: 0) */
03124 
03125                         i__2 = *lwork - iwork + 1;
03126                         cunglq_(m, n, m, &a[a_offset], lda, &work[itau], &
03127                                 work[iwork], &i__2, &ierr);
03128                         ie = 1;
03129                         itauq = itau;
03130                         itaup = itauq + *m;
03131                         iwork = itaup + *m;
03132 
03133 /*                    Bidiagonalize L in WORK(IU), copying result to */
03134 /*                    WORK(IR) */
03135 /*                    (CWorkspace: need   2*M*M+3*M, */
03136 /*                                 prefer 2*M*M+2*M+2*M*NB) */
03137 /*                    (RWorkspace: need   M) */
03138 
03139                         i__2 = *lwork - iwork + 1;
03140                         cgebrd_(m, m, &work[iu], &ldwrku, &s[1], &rwork[ie], &
03141                                 work[itauq], &work[itaup], &work[iwork], &
03142                                 i__2, &ierr);
03143                         clacpy_("L", m, m, &work[iu], &ldwrku, &work[ir], &
03144                                 ldwrkr);
03145 
03146 /*                    Generate right bidiagonalizing vectors in WORK(IU) */
03147 /*                    (CWorkspace: need   2*M*M+3*M-1, */
03148 /*                                 prefer 2*M*M+2*M+(M-1)*NB) */
03149 /*                    (RWorkspace: 0) */
03150 
03151                         i__2 = *lwork - iwork + 1;
03152                         cungbr_("P", m, m, m, &work[iu], &ldwrku, &work[itaup]
03153 , &work[iwork], &i__2, &ierr);
03154 
03155 /*                    Generate left bidiagonalizing vectors in WORK(IR) */
03156 /*                    (CWorkspace: need 2*M*M+3*M, prefer 2*M*M+2*M+M*NB) */
03157 /*                    (RWorkspace: 0) */
03158 
03159                         i__2 = *lwork - iwork + 1;
03160                         cungbr_("Q", m, m, m, &work[ir], &ldwrkr, &work[itauq]
03161 , &work[iwork], &i__2, &ierr);
03162                         irwork = ie + *m;
03163 
03164 /*                    Perform bidiagonal QR iteration, computing left */
03165 /*                    singular vectors of L in WORK(IR) and computing */
03166 /*                    right singular vectors of L in WORK(IU) */
03167 /*                    (CWorkspace: need 2*M*M) */
03168 /*                    (RWorkspace: need BDSPAC) */
03169 
03170                         cbdsqr_("U", m, m, m, &c__0, &s[1], &rwork[ie], &work[
03171                                 iu], &ldwrku, &work[ir], &ldwrkr, cdum, &c__1, 
03172                                  &rwork[irwork], info);
03173 
03174 /*                    Multiply right singular vectors of L in WORK(IU) by */
03175 /*                    Q in A, storing result in VT */
03176 /*                    (CWorkspace: need M*M) */
03177 /*                    (RWorkspace: 0) */
03178 
03179                         cgemm_("N", "N", m, n, m, &c_b2, &work[iu], &ldwrku, &
03180                                 a[a_offset], lda, &c_b1, &vt[vt_offset], ldvt);
03181 
03182 /*                    Copy left singular vectors of L to A */
03183 /*                    (CWorkspace: need M*M) */
03184 /*                    (RWorkspace: 0) */
03185 
03186                         clacpy_("F", m, m, &work[ir], &ldwrkr, &a[a_offset], 
03187                                 lda);
03188 
03189                     } else {
03190 
03191 /*                    Insufficient workspace for a fast algorithm */
03192 
03193                         itau = 1;
03194                         iwork = itau + *m;
03195 
03196 /*                    Compute A=L*Q, copying result to VT */
03197 /*                    (CWorkspace: need 2*M, prefer M+M*NB) */
03198 /*                    (RWorkspace: 0) */
03199 
03200                         i__2 = *lwork - iwork + 1;
03201                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03202                                 iwork], &i__2, &ierr);
03203                         clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], 
03204                                 ldvt);
03205 
03206 /*                    Generate Q in VT */
03207 /*                    (CWorkspace: need 2*M, prefer M+M*NB) */
03208 /*                    (RWorkspace: 0) */
03209 
03210                         i__2 = *lwork - iwork + 1;
03211                         cunglq_(m, n, m, &vt[vt_offset], ldvt, &work[itau], &
03212                                 work[iwork], &i__2, &ierr);
03213                         ie = 1;
03214                         itauq = itau;
03215                         itaup = itauq + *m;
03216                         iwork = itaup + *m;
03217 
03218 /*                    Zero out above L in A */
03219 
03220                         i__2 = *m - 1;
03221                         i__3 = *m - 1;
03222                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &a[(a_dim1 <<
03223                                  1) + 1], lda);
03224 
03225 /*                    Bidiagonalize L in A */
03226 /*                    (CWorkspace: need 3*M, prefer 2*M+2*M*NB) */
03227 /*                    (RWorkspace: need M) */
03228 
03229                         i__2 = *lwork - iwork + 1;
03230                         cgebrd_(m, m, &a[a_offset], lda, &s[1], &rwork[ie], &
03231                                 work[itauq], &work[itaup], &work[iwork], &
03232                                 i__2, &ierr);
03233 
03234 /*                    Multiply right vectors bidiagonalizing L by Q in VT */
03235 /*                    (CWorkspace: need 2*M+N, prefer 2*M+N*NB) */
03236 /*                    (RWorkspace: 0) */
03237 
03238                         i__2 = *lwork - iwork + 1;
03239                         cunmbr_("P", "L", "C", m, n, m, &a[a_offset], lda, &
03240                                 work[itaup], &vt[vt_offset], ldvt, &work[
03241                                 iwork], &i__2, &ierr);
03242 
03243 /*                    Generate left bidiagonalizing vectors of L in A */
03244 /*                    (CWorkspace: need 3*M, prefer 2*M+M*NB) */
03245 /*                    (RWorkspace: 0) */
03246 
03247                         i__2 = *lwork - iwork + 1;
03248                         cungbr_("Q", m, m, m, &a[a_offset], lda, &work[itauq], 
03249                                  &work[iwork], &i__2, &ierr);
03250                         irwork = ie + *m;
03251 
03252 /*                    Perform bidiagonal QR iteration, computing left */
03253 /*                    singular vectors of A in A and computing right */
03254 /*                    singular vectors of A in VT */
03255 /*                    (CWorkspace: 0) */
03256 /*                    (RWorkspace: need BDSPAC) */
03257 
03258                         cbdsqr_("U", m, n, m, &c__0, &s[1], &rwork[ie], &vt[
03259                                 vt_offset], ldvt, &a[a_offset], lda, cdum, &
03260                                 c__1, &rwork[irwork], info);
03261 
03262                     }
03263 
03264                 } else if (wntuas) {
03265 
03266 /*                 Path 6t(N much larger than M, JOBU='S' or 'A', */
03267 /*                         JOBVT='S') */
03268 /*                 M right singular vectors to be computed in VT and */
03269 /*                 M left singular vectors to be computed in U */
03270 
03271                     if (*lwork >= *m * *m + *m * 3) {
03272 
03273 /*                    Sufficient workspace for a fast algorithm */
03274 
03275                         iu = 1;
03276                         if (*lwork >= wrkbl + *lda * *m) {
03277 
03278 /*                       WORK(IU) is LDA by N */
03279 
03280                             ldwrku = *lda;
03281                         } else {
03282 
03283 /*                       WORK(IU) is LDA by M */
03284 
03285                             ldwrku = *m;
03286                         }
03287                         itau = iu + ldwrku * *m;
03288                         iwork = itau + *m;
03289 
03290 /*                    Compute A=L*Q */
03291 /*                    (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
03292 /*                    (RWorkspace: 0) */
03293 
03294                         i__2 = *lwork - iwork + 1;
03295                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03296                                 iwork], &i__2, &ierr);
03297 
03298 /*                    Copy L to WORK(IU), zeroing out above it */
03299 
03300                         clacpy_("L", m, m, &a[a_offset], lda, &work[iu], &
03301                                 ldwrku);
03302                         i__2 = *m - 1;
03303                         i__3 = *m - 1;
03304                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &work[iu + 
03305                                 ldwrku], &ldwrku);
03306 
03307 /*                    Generate Q in A */
03308 /*                    (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
03309 /*                    (RWorkspace: 0) */
03310 
03311                         i__2 = *lwork - iwork + 1;
03312                         cunglq_(m, n, m, &a[a_offset], lda, &work[itau], &
03313                                 work[iwork], &i__2, &ierr);
03314                         ie = 1;
03315                         itauq = itau;
03316                         itaup = itauq + *m;
03317                         iwork = itaup + *m;
03318 
03319 /*                    Bidiagonalize L in WORK(IU), copying result to U */
03320 /*                    (CWorkspace: need M*M+3*M, prefer M*M+2*M+2*M*NB) */
03321 /*                    (RWorkspace: need M) */
03322 
03323                         i__2 = *lwork - iwork + 1;
03324                         cgebrd_(m, m, &work[iu], &ldwrku, &s[1], &rwork[ie], &
03325                                 work[itauq], &work[itaup], &work[iwork], &
03326                                 i__2, &ierr);
03327                         clacpy_("L", m, m, &work[iu], &ldwrku, &u[u_offset], 
03328                                 ldu);
03329 
03330 /*                    Generate right bidiagonalizing vectors in WORK(IU) */
03331 /*                    (CWorkspace: need   M*M+3*M-1, */
03332 /*                                 prefer M*M+2*M+(M-1)*NB) */
03333 /*                    (RWorkspace: 0) */
03334 
03335                         i__2 = *lwork - iwork + 1;
03336                         cungbr_("P", m, m, m, &work[iu], &ldwrku, &work[itaup]
03337 , &work[iwork], &i__2, &ierr);
03338 
03339 /*                    Generate left bidiagonalizing vectors in U */
03340 /*                    (CWorkspace: need M*M+3*M, prefer M*M+2*M+M*NB) */
03341 /*                    (RWorkspace: 0) */
03342 
03343                         i__2 = *lwork - iwork + 1;
03344                         cungbr_("Q", m, m, m, &u[u_offset], ldu, &work[itauq], 
03345                                  &work[iwork], &i__2, &ierr);
03346                         irwork = ie + *m;
03347 
03348 /*                    Perform bidiagonal QR iteration, computing left */
03349 /*                    singular vectors of L in U and computing right */
03350 /*                    singular vectors of L in WORK(IU) */
03351 /*                    (CWorkspace: need M*M) */
03352 /*                    (RWorkspace: need BDSPAC) */
03353 
03354                         cbdsqr_("U", m, m, m, &c__0, &s[1], &rwork[ie], &work[
03355                                 iu], &ldwrku, &u[u_offset], ldu, cdum, &c__1, 
03356                                 &rwork[irwork], info);
03357 
03358 /*                    Multiply right singular vectors of L in WORK(IU) by */
03359 /*                    Q in A, storing result in VT */
03360 /*                    (CWorkspace: need M*M) */
03361 /*                    (RWorkspace: 0) */
03362 
03363                         cgemm_("N", "N", m, n, m, &c_b2, &work[iu], &ldwrku, &
03364                                 a[a_offset], lda, &c_b1, &vt[vt_offset], ldvt);
03365 
03366                     } else {
03367 
03368 /*                    Insufficient workspace for a fast algorithm */
03369 
03370                         itau = 1;
03371                         iwork = itau + *m;
03372 
03373 /*                    Compute A=L*Q, copying result to VT */
03374 /*                    (CWorkspace: need 2*M, prefer M+M*NB) */
03375 /*                    (RWorkspace: 0) */
03376 
03377                         i__2 = *lwork - iwork + 1;
03378                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03379                                 iwork], &i__2, &ierr);
03380                         clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], 
03381                                 ldvt);
03382 
03383 /*                    Generate Q in VT */
03384 /*                    (CWorkspace: need 2*M, prefer M+M*NB) */
03385 /*                    (RWorkspace: 0) */
03386 
03387                         i__2 = *lwork - iwork + 1;
03388                         cunglq_(m, n, m, &vt[vt_offset], ldvt, &work[itau], &
03389                                 work[iwork], &i__2, &ierr);
03390 
03391 /*                    Copy L to U, zeroing out above it */
03392 
03393                         clacpy_("L", m, m, &a[a_offset], lda, &u[u_offset], 
03394                                 ldu);
03395                         i__2 = *m - 1;
03396                         i__3 = *m - 1;
03397                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &u[(u_dim1 <<
03398                                  1) + 1], ldu);
03399                         ie = 1;
03400                         itauq = itau;
03401                         itaup = itauq + *m;
03402                         iwork = itaup + *m;
03403 
03404 /*                    Bidiagonalize L in U */
03405 /*                    (CWorkspace: need 3*M, prefer 2*M+2*M*NB) */
03406 /*                    (RWorkspace: need M) */
03407 
03408                         i__2 = *lwork - iwork + 1;
03409                         cgebrd_(m, m, &u[u_offset], ldu, &s[1], &rwork[ie], &
03410                                 work[itauq], &work[itaup], &work[iwork], &
03411                                 i__2, &ierr);
03412 
03413 /*                    Multiply right bidiagonalizing vectors in U by Q */
03414 /*                    in VT */
03415 /*                    (CWorkspace: need 2*M+N, prefer 2*M+N*NB) */
03416 /*                    (RWorkspace: 0) */
03417 
03418                         i__2 = *lwork - iwork + 1;
03419                         cunmbr_("P", "L", "C", m, n, m, &u[u_offset], ldu, &
03420                                 work[itaup], &vt[vt_offset], ldvt, &work[
03421                                 iwork], &i__2, &ierr);
03422 
03423 /*                    Generate left bidiagonalizing vectors in U */
03424 /*                    (CWorkspace: need 3*M, prefer 2*M+M*NB) */
03425 /*                    (RWorkspace: 0) */
03426 
03427                         i__2 = *lwork - iwork + 1;
03428                         cungbr_("Q", m, m, m, &u[u_offset], ldu, &work[itauq], 
03429                                  &work[iwork], &i__2, &ierr);
03430                         irwork = ie + *m;
03431 
03432 /*                    Perform bidiagonal QR iteration, computing left */
03433 /*                    singular vectors of A in U and computing right */
03434 /*                    singular vectors of A in VT */
03435 /*                    (CWorkspace: 0) */
03436 /*                    (RWorkspace: need BDSPAC) */
03437 
03438                         cbdsqr_("U", m, n, m, &c__0, &s[1], &rwork[ie], &vt[
03439                                 vt_offset], ldvt, &u[u_offset], ldu, cdum, &
03440                                 c__1, &rwork[irwork], info);
03441 
03442                     }
03443 
03444                 }
03445 
03446             } else if (wntva) {
03447 
03448                 if (wntun) {
03449 
03450 /*                 Path 7t(N much larger than M, JOBU='N', JOBVT='A') */
03451 /*                 N right singular vectors to be computed in VT and */
03452 /*                 no left singular vectors to be computed */
03453 
03454 /* Computing MAX */
03455                     i__2 = *n + *m, i__3 = *m * 3;
03456                     if (*lwork >= *m * *m + max(i__2,i__3)) {
03457 
03458 /*                    Sufficient workspace for a fast algorithm */
03459 
03460                         ir = 1;
03461                         if (*lwork >= wrkbl + *lda * *m) {
03462 
03463 /*                       WORK(IR) is LDA by M */
03464 
03465                             ldwrkr = *lda;
03466                         } else {
03467 
03468 /*                       WORK(IR) is M by M */
03469 
03470                             ldwrkr = *m;
03471                         }
03472                         itau = ir + ldwrkr * *m;
03473                         iwork = itau + *m;
03474 
03475 /*                    Compute A=L*Q, copying result to VT */
03476 /*                    (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
03477 /*                    (RWorkspace: 0) */
03478 
03479                         i__2 = *lwork - iwork + 1;
03480                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03481                                 iwork], &i__2, &ierr);
03482                         clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], 
03483                                 ldvt);
03484 
03485 /*                    Copy L to WORK(IR), zeroing out above it */
03486 
03487                         clacpy_("L", m, m, &a[a_offset], lda, &work[ir], &
03488                                 ldwrkr);
03489                         i__2 = *m - 1;
03490                         i__3 = *m - 1;
03491                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &work[ir + 
03492                                 ldwrkr], &ldwrkr);
03493 
03494 /*                    Generate Q in VT */
03495 /*                    (CWorkspace: need M*M+M+N, prefer M*M+M+N*NB) */
03496 /*                    (RWorkspace: 0) */
03497 
03498                         i__2 = *lwork - iwork + 1;
03499                         cunglq_(n, n, m, &vt[vt_offset], ldvt, &work[itau], &
03500                                 work[iwork], &i__2, &ierr);
03501                         ie = 1;
03502                         itauq = itau;
03503                         itaup = itauq + *m;
03504                         iwork = itaup + *m;
03505 
03506 /*                    Bidiagonalize L in WORK(IR) */
03507 /*                    (CWorkspace: need M*M+3*M, prefer M*M+2*M+2*M*NB) */
03508 /*                    (RWorkspace: need M) */
03509 
03510                         i__2 = *lwork - iwork + 1;
03511                         cgebrd_(m, m, &work[ir], &ldwrkr, &s[1], &rwork[ie], &
03512                                 work[itauq], &work[itaup], &work[iwork], &
03513                                 i__2, &ierr);
03514 
03515 /*                    Generate right bidiagonalizing vectors in WORK(IR) */
03516 /*                    (CWorkspace: need   M*M+3*M-1, */
03517 /*                                 prefer M*M+2*M+(M-1)*NB) */
03518 /*                    (RWorkspace: 0) */
03519 
03520                         i__2 = *lwork - iwork + 1;
03521                         cungbr_("P", m, m, m, &work[ir], &ldwrkr, &work[itaup]
03522 , &work[iwork], &i__2, &ierr);
03523                         irwork = ie + *m;
03524 
03525 /*                    Perform bidiagonal QR iteration, computing right */
03526 /*                    singular vectors of L in WORK(IR) */
03527 /*                    (CWorkspace: need M*M) */
03528 /*                    (RWorkspace: need BDSPAC) */
03529 
03530                         cbdsqr_("U", m, m, &c__0, &c__0, &s[1], &rwork[ie], &
03531                                 work[ir], &ldwrkr, cdum, &c__1, cdum, &c__1, &
03532                                 rwork[irwork], info);
03533 
03534 /*                    Multiply right singular vectors of L in WORK(IR) by */
03535 /*                    Q in VT, storing result in A */
03536 /*                    (CWorkspace: need M*M) */
03537 /*                    (RWorkspace: 0) */
03538 
03539                         cgemm_("N", "N", m, n, m, &c_b2, &work[ir], &ldwrkr, &
03540                                 vt[vt_offset], ldvt, &c_b1, &a[a_offset], lda);
03541 
03542 /*                    Copy right singular vectors of A from A to VT */
03543 
03544                         clacpy_("F", m, n, &a[a_offset], lda, &vt[vt_offset], 
03545                                 ldvt);
03546 
03547                     } else {
03548 
03549 /*                    Insufficient workspace for a fast algorithm */
03550 
03551                         itau = 1;
03552                         iwork = itau + *m;
03553 
03554 /*                    Compute A=L*Q, copying result to VT */
03555 /*                    (CWorkspace: need 2*M, prefer M+M*NB) */
03556 /*                    (RWorkspace: 0) */
03557 
03558                         i__2 = *lwork - iwork + 1;
03559                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03560                                 iwork], &i__2, &ierr);
03561                         clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], 
03562                                 ldvt);
03563 
03564 /*                    Generate Q in VT */
03565 /*                    (CWorkspace: need M+N, prefer M+N*NB) */
03566 /*                    (RWorkspace: 0) */
03567 
03568                         i__2 = *lwork - iwork + 1;
03569                         cunglq_(n, n, m, &vt[vt_offset], ldvt, &work[itau], &
03570                                 work[iwork], &i__2, &ierr);
03571                         ie = 1;
03572                         itauq = itau;
03573                         itaup = itauq + *m;
03574                         iwork = itaup + *m;
03575 
03576 /*                    Zero out above L in A */
03577 
03578                         i__2 = *m - 1;
03579                         i__3 = *m - 1;
03580                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &a[(a_dim1 <<
03581                                  1) + 1], lda);
03582 
03583 /*                    Bidiagonalize L in A */
03584 /*                    (CWorkspace: need 3*M, prefer 2*M+2*M*NB) */
03585 /*                    (RWorkspace: need M) */
03586 
03587                         i__2 = *lwork - iwork + 1;
03588                         cgebrd_(m, m, &a[a_offset], lda, &s[1], &rwork[ie], &
03589                                 work[itauq], &work[itaup], &work[iwork], &
03590                                 i__2, &ierr);
03591 
03592 /*                    Multiply right bidiagonalizing vectors in A by Q */
03593 /*                    in VT */
03594 /*                    (CWorkspace: need 2*M+N, prefer 2*M+N*NB) */
03595 /*                    (RWorkspace: 0) */
03596 
03597                         i__2 = *lwork - iwork + 1;
03598                         cunmbr_("P", "L", "C", m, n, m, &a[a_offset], lda, &
03599                                 work[itaup], &vt[vt_offset], ldvt, &work[
03600                                 iwork], &i__2, &ierr);
03601                         irwork = ie + *m;
03602 
03603 /*                    Perform bidiagonal QR iteration, computing right */
03604 /*                    singular vectors of A in VT */
03605 /*                    (CWorkspace: 0) */
03606 /*                    (RWorkspace: need BDSPAC) */
03607 
03608                         cbdsqr_("U", m, n, &c__0, &c__0, &s[1], &rwork[ie], &
03609                                 vt[vt_offset], ldvt, cdum, &c__1, cdum, &c__1, 
03610                                  &rwork[irwork], info);
03611 
03612                     }
03613 
03614                 } else if (wntuo) {
03615 
03616 /*                 Path 8t(N much larger than M, JOBU='O', JOBVT='A') */
03617 /*                 N right singular vectors to be computed in VT and */
03618 /*                 M left singular vectors to be overwritten on A */
03619 
03620 /* Computing MAX */
03621                     i__2 = *n + *m, i__3 = *m * 3;
03622                     if (*lwork >= (*m << 1) * *m + max(i__2,i__3)) {
03623 
03624 /*                    Sufficient workspace for a fast algorithm */
03625 
03626                         iu = 1;
03627                         if (*lwork >= wrkbl + (*lda << 1) * *m) {
03628 
03629 /*                       WORK(IU) is LDA by M and WORK(IR) is LDA by M */
03630 
03631                             ldwrku = *lda;
03632                             ir = iu + ldwrku * *m;
03633                             ldwrkr = *lda;
03634                         } else if (*lwork >= wrkbl + (*lda + *m) * *m) {
03635 
03636 /*                       WORK(IU) is LDA by M and WORK(IR) is M by M */
03637 
03638                             ldwrku = *lda;
03639                             ir = iu + ldwrku * *m;
03640                             ldwrkr = *m;
03641                         } else {
03642 
03643 /*                       WORK(IU) is M by M and WORK(IR) is M by M */
03644 
03645                             ldwrku = *m;
03646                             ir = iu + ldwrku * *m;
03647                             ldwrkr = *m;
03648                         }
03649                         itau = ir + ldwrkr * *m;
03650                         iwork = itau + *m;
03651 
03652 /*                    Compute A=L*Q, copying result to VT */
03653 /*                    (CWorkspace: need 2*M*M+2*M, prefer 2*M*M+M+M*NB) */
03654 /*                    (RWorkspace: 0) */
03655 
03656                         i__2 = *lwork - iwork + 1;
03657                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03658                                 iwork], &i__2, &ierr);
03659                         clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], 
03660                                 ldvt);
03661 
03662 /*                    Generate Q in VT */
03663 /*                    (CWorkspace: need 2*M*M+M+N, prefer 2*M*M+M+N*NB) */
03664 /*                    (RWorkspace: 0) */
03665 
03666                         i__2 = *lwork - iwork + 1;
03667                         cunglq_(n, n, m, &vt[vt_offset], ldvt, &work[itau], &
03668                                 work[iwork], &i__2, &ierr);
03669 
03670 /*                    Copy L to WORK(IU), zeroing out above it */
03671 
03672                         clacpy_("L", m, m, &a[a_offset], lda, &work[iu], &
03673                                 ldwrku);
03674                         i__2 = *m - 1;
03675                         i__3 = *m - 1;
03676                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &work[iu + 
03677                                 ldwrku], &ldwrku);
03678                         ie = 1;
03679                         itauq = itau;
03680                         itaup = itauq + *m;
03681                         iwork = itaup + *m;
03682 
03683 /*                    Bidiagonalize L in WORK(IU), copying result to */
03684 /*                    WORK(IR) */
03685 /*                    (CWorkspace: need   2*M*M+3*M, */
03686 /*                                 prefer 2*M*M+2*M+2*M*NB) */
03687 /*                    (RWorkspace: need   M) */
03688 
03689                         i__2 = *lwork - iwork + 1;
03690                         cgebrd_(m, m, &work[iu], &ldwrku, &s[1], &rwork[ie], &
03691                                 work[itauq], &work[itaup], &work[iwork], &
03692                                 i__2, &ierr);
03693                         clacpy_("L", m, m, &work[iu], &ldwrku, &work[ir], &
03694                                 ldwrkr);
03695 
03696 /*                    Generate right bidiagonalizing vectors in WORK(IU) */
03697 /*                    (CWorkspace: need   2*M*M+3*M-1, */
03698 /*                                 prefer 2*M*M+2*M+(M-1)*NB) */
03699 /*                    (RWorkspace: 0) */
03700 
03701                         i__2 = *lwork - iwork + 1;
03702                         cungbr_("P", m, m, m, &work[iu], &ldwrku, &work[itaup]
03703 , &work[iwork], &i__2, &ierr);
03704 
03705 /*                    Generate left bidiagonalizing vectors in WORK(IR) */
03706 /*                    (CWorkspace: need 2*M*M+3*M, prefer 2*M*M+2*M+M*NB) */
03707 /*                    (RWorkspace: 0) */
03708 
03709                         i__2 = *lwork - iwork + 1;
03710                         cungbr_("Q", m, m, m, &work[ir], &ldwrkr, &work[itauq]
03711 , &work[iwork], &i__2, &ierr);
03712                         irwork = ie + *m;
03713 
03714 /*                    Perform bidiagonal QR iteration, computing left */
03715 /*                    singular vectors of L in WORK(IR) and computing */
03716 /*                    right singular vectors of L in WORK(IU) */
03717 /*                    (CWorkspace: need 2*M*M) */
03718 /*                    (RWorkspace: need BDSPAC) */
03719 
03720                         cbdsqr_("U", m, m, m, &c__0, &s[1], &rwork[ie], &work[
03721                                 iu], &ldwrku, &work[ir], &ldwrkr, cdum, &c__1, 
03722                                  &rwork[irwork], info);
03723 
03724 /*                    Multiply right singular vectors of L in WORK(IU) by */
03725 /*                    Q in VT, storing result in A */
03726 /*                    (CWorkspace: need M*M) */
03727 /*                    (RWorkspace: 0) */
03728 
03729                         cgemm_("N", "N", m, n, m, &c_b2, &work[iu], &ldwrku, &
03730                                 vt[vt_offset], ldvt, &c_b1, &a[a_offset], lda);
03731 
03732 /*                    Copy right singular vectors of A from A to VT */
03733 
03734                         clacpy_("F", m, n, &a[a_offset], lda, &vt[vt_offset], 
03735                                 ldvt);
03736 
03737 /*                    Copy left singular vectors of A from WORK(IR) to A */
03738 
03739                         clacpy_("F", m, m, &work[ir], &ldwrkr, &a[a_offset], 
03740                                 lda);
03741 
03742                     } else {
03743 
03744 /*                    Insufficient workspace for a fast algorithm */
03745 
03746                         itau = 1;
03747                         iwork = itau + *m;
03748 
03749 /*                    Compute A=L*Q, copying result to VT */
03750 /*                    (CWorkspace: need 2*M, prefer M+M*NB) */
03751 /*                    (RWorkspace: 0) */
03752 
03753                         i__2 = *lwork - iwork + 1;
03754                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03755                                 iwork], &i__2, &ierr);
03756                         clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], 
03757                                 ldvt);
03758 
03759 /*                    Generate Q in VT */
03760 /*                    (CWorkspace: need M+N, prefer M+N*NB) */
03761 /*                    (RWorkspace: 0) */
03762 
03763                         i__2 = *lwork - iwork + 1;
03764                         cunglq_(n, n, m, &vt[vt_offset], ldvt, &work[itau], &
03765                                 work[iwork], &i__2, &ierr);
03766                         ie = 1;
03767                         itauq = itau;
03768                         itaup = itauq + *m;
03769                         iwork = itaup + *m;
03770 
03771 /*                    Zero out above L in A */
03772 
03773                         i__2 = *m - 1;
03774                         i__3 = *m - 1;
03775                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &a[(a_dim1 <<
03776                                  1) + 1], lda);
03777 
03778 /*                    Bidiagonalize L in A */
03779 /*                    (CWorkspace: need 3*M, prefer 2*M+2*M*NB) */
03780 /*                    (RWorkspace: need M) */
03781 
03782                         i__2 = *lwork - iwork + 1;
03783                         cgebrd_(m, m, &a[a_offset], lda, &s[1], &rwork[ie], &
03784                                 work[itauq], &work[itaup], &work[iwork], &
03785                                 i__2, &ierr);
03786 
03787 /*                    Multiply right bidiagonalizing vectors in A by Q */
03788 /*                    in VT */
03789 /*                    (CWorkspace: need 2*M+N, prefer 2*M+N*NB) */
03790 /*                    (RWorkspace: 0) */
03791 
03792                         i__2 = *lwork - iwork + 1;
03793                         cunmbr_("P", "L", "C", m, n, m, &a[a_offset], lda, &
03794                                 work[itaup], &vt[vt_offset], ldvt, &work[
03795                                 iwork], &i__2, &ierr);
03796 
03797 /*                    Generate left bidiagonalizing vectors in A */
03798 /*                    (CWorkspace: need 3*M, prefer 2*M+M*NB) */
03799 /*                    (RWorkspace: 0) */
03800 
03801                         i__2 = *lwork - iwork + 1;
03802                         cungbr_("Q", m, m, m, &a[a_offset], lda, &work[itauq], 
03803                                  &work[iwork], &i__2, &ierr);
03804                         irwork = ie + *m;
03805 
03806 /*                    Perform bidiagonal QR iteration, computing left */
03807 /*                    singular vectors of A in A and computing right */
03808 /*                    singular vectors of A in VT */
03809 /*                    (CWorkspace: 0) */
03810 /*                    (RWorkspace: need BDSPAC) */
03811 
03812                         cbdsqr_("U", m, n, m, &c__0, &s[1], &rwork[ie], &vt[
03813                                 vt_offset], ldvt, &a[a_offset], lda, cdum, &
03814                                 c__1, &rwork[irwork], info);
03815 
03816                     }
03817 
03818                 } else if (wntuas) {
03819 
03820 /*                 Path 9t(N much larger than M, JOBU='S' or 'A', */
03821 /*                         JOBVT='A') */
03822 /*                 N right singular vectors to be computed in VT and */
03823 /*                 M left singular vectors to be computed in U */
03824 
03825 /* Computing MAX */
03826                     i__2 = *n + *m, i__3 = *m * 3;
03827                     if (*lwork >= *m * *m + max(i__2,i__3)) {
03828 
03829 /*                    Sufficient workspace for a fast algorithm */
03830 
03831                         iu = 1;
03832                         if (*lwork >= wrkbl + *lda * *m) {
03833 
03834 /*                       WORK(IU) is LDA by M */
03835 
03836                             ldwrku = *lda;
03837                         } else {
03838 
03839 /*                       WORK(IU) is M by M */
03840 
03841                             ldwrku = *m;
03842                         }
03843                         itau = iu + ldwrku * *m;
03844                         iwork = itau + *m;
03845 
03846 /*                    Compute A=L*Q, copying result to VT */
03847 /*                    (CWorkspace: need M*M+2*M, prefer M*M+M+M*NB) */
03848 /*                    (RWorkspace: 0) */
03849 
03850                         i__2 = *lwork - iwork + 1;
03851                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03852                                 iwork], &i__2, &ierr);
03853                         clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], 
03854                                 ldvt);
03855 
03856 /*                    Generate Q in VT */
03857 /*                    (CWorkspace: need M*M+M+N, prefer M*M+M+N*NB) */
03858 /*                    (RWorkspace: 0) */
03859 
03860                         i__2 = *lwork - iwork + 1;
03861                         cunglq_(n, n, m, &vt[vt_offset], ldvt, &work[itau], &
03862                                 work[iwork], &i__2, &ierr);
03863 
03864 /*                    Copy L to WORK(IU), zeroing out above it */
03865 
03866                         clacpy_("L", m, m, &a[a_offset], lda, &work[iu], &
03867                                 ldwrku);
03868                         i__2 = *m - 1;
03869                         i__3 = *m - 1;
03870                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &work[iu + 
03871                                 ldwrku], &ldwrku);
03872                         ie = 1;
03873                         itauq = itau;
03874                         itaup = itauq + *m;
03875                         iwork = itaup + *m;
03876 
03877 /*                    Bidiagonalize L in WORK(IU), copying result to U */
03878 /*                    (CWorkspace: need M*M+3*M, prefer M*M+2*M+2*M*NB) */
03879 /*                    (RWorkspace: need M) */
03880 
03881                         i__2 = *lwork - iwork + 1;
03882                         cgebrd_(m, m, &work[iu], &ldwrku, &s[1], &rwork[ie], &
03883                                 work[itauq], &work[itaup], &work[iwork], &
03884                                 i__2, &ierr);
03885                         clacpy_("L", m, m, &work[iu], &ldwrku, &u[u_offset], 
03886                                 ldu);
03887 
03888 /*                    Generate right bidiagonalizing vectors in WORK(IU) */
03889 /*                    (CWorkspace: need M*M+3*M, prefer M*M+2*M+(M-1)*NB) */
03890 /*                    (RWorkspace: 0) */
03891 
03892                         i__2 = *lwork - iwork + 1;
03893                         cungbr_("P", m, m, m, &work[iu], &ldwrku, &work[itaup]
03894 , &work[iwork], &i__2, &ierr);
03895 
03896 /*                    Generate left bidiagonalizing vectors in U */
03897 /*                    (CWorkspace: need M*M+3*M, prefer M*M+2*M+M*NB) */
03898 /*                    (RWorkspace: 0) */
03899 
03900                         i__2 = *lwork - iwork + 1;
03901                         cungbr_("Q", m, m, m, &u[u_offset], ldu, &work[itauq], 
03902                                  &work[iwork], &i__2, &ierr);
03903                         irwork = ie + *m;
03904 
03905 /*                    Perform bidiagonal QR iteration, computing left */
03906 /*                    singular vectors of L in U and computing right */
03907 /*                    singular vectors of L in WORK(IU) */
03908 /*                    (CWorkspace: need M*M) */
03909 /*                    (RWorkspace: need BDSPAC) */
03910 
03911                         cbdsqr_("U", m, m, m, &c__0, &s[1], &rwork[ie], &work[
03912                                 iu], &ldwrku, &u[u_offset], ldu, cdum, &c__1, 
03913                                 &rwork[irwork], info);
03914 
03915 /*                    Multiply right singular vectors of L in WORK(IU) by */
03916 /*                    Q in VT, storing result in A */
03917 /*                    (CWorkspace: need M*M) */
03918 /*                    (RWorkspace: 0) */
03919 
03920                         cgemm_("N", "N", m, n, m, &c_b2, &work[iu], &ldwrku, &
03921                                 vt[vt_offset], ldvt, &c_b1, &a[a_offset], lda);
03922 
03923 /*                    Copy right singular vectors of A from A to VT */
03924 
03925                         clacpy_("F", m, n, &a[a_offset], lda, &vt[vt_offset], 
03926                                 ldvt);
03927 
03928                     } else {
03929 
03930 /*                    Insufficient workspace for a fast algorithm */
03931 
03932                         itau = 1;
03933                         iwork = itau + *m;
03934 
03935 /*                    Compute A=L*Q, copying result to VT */
03936 /*                    (CWorkspace: need 2*M, prefer M+M*NB) */
03937 /*                    (RWorkspace: 0) */
03938 
03939                         i__2 = *lwork - iwork + 1;
03940                         cgelqf_(m, n, &a[a_offset], lda, &work[itau], &work[
03941                                 iwork], &i__2, &ierr);
03942                         clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], 
03943                                 ldvt);
03944 
03945 /*                    Generate Q in VT */
03946 /*                    (CWorkspace: need M+N, prefer M+N*NB) */
03947 /*                    (RWorkspace: 0) */
03948 
03949                         i__2 = *lwork - iwork + 1;
03950                         cunglq_(n, n, m, &vt[vt_offset], ldvt, &work[itau], &
03951                                 work[iwork], &i__2, &ierr);
03952 
03953 /*                    Copy L to U, zeroing out above it */
03954 
03955                         clacpy_("L", m, m, &a[a_offset], lda, &u[u_offset], 
03956                                 ldu);
03957                         i__2 = *m - 1;
03958                         i__3 = *m - 1;
03959                         claset_("U", &i__2, &i__3, &c_b1, &c_b1, &u[(u_dim1 <<
03960                                  1) + 1], ldu);
03961                         ie = 1;
03962                         itauq = itau;
03963                         itaup = itauq + *m;
03964                         iwork = itaup + *m;
03965 
03966 /*                    Bidiagonalize L in U */
03967 /*                    (CWorkspace: need 3*M, prefer 2*M+2*M*NB) */
03968 /*                    (RWorkspace: need M) */
03969 
03970                         i__2 = *lwork - iwork + 1;
03971                         cgebrd_(m, m, &u[u_offset], ldu, &s[1], &rwork[ie], &
03972                                 work[itauq], &work[itaup], &work[iwork], &
03973                                 i__2, &ierr);
03974 
03975 /*                    Multiply right bidiagonalizing vectors in U by Q */
03976 /*                    in VT */
03977 /*                    (CWorkspace: need 2*M+N, prefer 2*M+N*NB) */
03978 /*                    (RWorkspace: 0) */
03979 
03980                         i__2 = *lwork - iwork + 1;
03981                         cunmbr_("P", "L", "C", m, n, m, &u[u_offset], ldu, &
03982                                 work[itaup], &vt[vt_offset], ldvt, &work[
03983                                 iwork], &i__2, &ierr);
03984 
03985 /*                    Generate left bidiagonalizing vectors in U */
03986 /*                    (CWorkspace: need 3*M, prefer 2*M+M*NB) */
03987 /*                    (RWorkspace: 0) */
03988 
03989                         i__2 = *lwork - iwork + 1;
03990                         cungbr_("Q", m, m, m, &u[u_offset], ldu, &work[itauq], 
03991                                  &work[iwork], &i__2, &ierr);
03992                         irwork = ie + *m;
03993 
03994 /*                    Perform bidiagonal QR iteration, computing left */
03995 /*                    singular vectors of A in U and computing right */
03996 /*                    singular vectors of A in VT */
03997 /*                    (CWorkspace: 0) */
03998 /*                    (RWorkspace: need BDSPAC) */
03999 
04000                         cbdsqr_("U", m, n, m, &c__0, &s[1], &rwork[ie], &vt[
04001                                 vt_offset], ldvt, &u[u_offset], ldu, cdum, &
04002                                 c__1, &rwork[irwork], info);
04003 
04004                     }
04005 
04006                 }
04007 
04008             }
04009 
04010         } else {
04011 
04012 /*           N .LT. MNTHR */
04013 
04014 /*           Path 10t(N greater than M, but not much larger) */
04015 /*           Reduce to bidiagonal form without LQ decomposition */
04016 
04017             ie = 1;
04018             itauq = 1;
04019             itaup = itauq + *m;
04020             iwork = itaup + *m;
04021 
04022 /*           Bidiagonalize A */
04023 /*           (CWorkspace: need 2*M+N, prefer 2*M+(M+N)*NB) */
04024 /*           (RWorkspace: M) */
04025 
04026             i__2 = *lwork - iwork + 1;
04027             cgebrd_(m, n, &a[a_offset], lda, &s[1], &rwork[ie], &work[itauq], 
04028                     &work[itaup], &work[iwork], &i__2, &ierr);
04029             if (wntuas) {
04030 
04031 /*              If left singular vectors desired in U, copy result to U */
04032 /*              and generate left bidiagonalizing vectors in U */
04033 /*              (CWorkspace: need 3*M-1, prefer 2*M+(M-1)*NB) */
04034 /*              (RWorkspace: 0) */
04035 
04036                 clacpy_("L", m, m, &a[a_offset], lda, &u[u_offset], ldu);
04037                 i__2 = *lwork - iwork + 1;
04038                 cungbr_("Q", m, m, n, &u[u_offset], ldu, &work[itauq], &work[
04039                         iwork], &i__2, &ierr);
04040             }
04041             if (wntvas) {
04042 
04043 /*              If right singular vectors desired in VT, copy result to */
04044 /*              VT and generate right bidiagonalizing vectors in VT */
04045 /*              (CWorkspace: need 2*M+NRVT, prefer 2*M+NRVT*NB) */
04046 /*              (RWorkspace: 0) */
04047 
04048                 clacpy_("U", m, n, &a[a_offset], lda, &vt[vt_offset], ldvt);
04049                 if (wntva) {
04050                     nrvt = *n;
04051                 }
04052                 if (wntvs) {
04053                     nrvt = *m;
04054                 }
04055                 i__2 = *lwork - iwork + 1;
04056                 cungbr_("P", &nrvt, n, m, &vt[vt_offset], ldvt, &work[itaup], 
04057                         &work[iwork], &i__2, &ierr);
04058             }
04059             if (wntuo) {
04060 
04061 /*              If left singular vectors desired in A, generate left */
04062 /*              bidiagonalizing vectors in A */
04063 /*              (CWorkspace: need 3*M-1, prefer 2*M+(M-1)*NB) */
04064 /*              (RWorkspace: 0) */
04065 
04066                 i__2 = *lwork - iwork + 1;
04067                 cungbr_("Q", m, m, n, &a[a_offset], lda, &work[itauq], &work[
04068                         iwork], &i__2, &ierr);
04069             }
04070             if (wntvo) {
04071 
04072 /*              If right singular vectors desired in A, generate right */
04073 /*              bidiagonalizing vectors in A */
04074 /*              (CWorkspace: need 3*M, prefer 2*M+M*NB) */
04075 /*              (RWorkspace: 0) */
04076 
04077                 i__2 = *lwork - iwork + 1;
04078                 cungbr_("P", m, n, m, &a[a_offset], lda, &work[itaup], &work[
04079                         iwork], &i__2, &ierr);
04080             }
04081             irwork = ie + *m;
04082             if (wntuas || wntuo) {
04083                 nru = *m;
04084             }
04085             if (wntun) {
04086                 nru = 0;
04087             }
04088             if (wntvas || wntvo) {
04089                 ncvt = *n;
04090             }
04091             if (wntvn) {
04092                 ncvt = 0;
04093             }
04094             if (! wntuo && ! wntvo) {
04095 
04096 /*              Perform bidiagonal QR iteration, if desired, computing */
04097 /*              left singular vectors in U and computing right singular */
04098 /*              vectors in VT */
04099 /*              (CWorkspace: 0) */
04100 /*              (RWorkspace: need BDSPAC) */
04101 
04102                 cbdsqr_("L", m, &ncvt, &nru, &c__0, &s[1], &rwork[ie], &vt[
04103                         vt_offset], ldvt, &u[u_offset], ldu, cdum, &c__1, &
04104                         rwork[irwork], info);
04105             } else if (! wntuo && wntvo) {
04106 
04107 /*              Perform bidiagonal QR iteration, if desired, computing */
04108 /*              left singular vectors in U and computing right singular */
04109 /*              vectors in A */
04110 /*              (CWorkspace: 0) */
04111 /*              (RWorkspace: need BDSPAC) */
04112 
04113                 cbdsqr_("L", m, &ncvt, &nru, &c__0, &s[1], &rwork[ie], &a[
04114                         a_offset], lda, &u[u_offset], ldu, cdum, &c__1, &
04115                         rwork[irwork], info);
04116             } else {
04117 
04118 /*              Perform bidiagonal QR iteration, if desired, computing */
04119 /*              left singular vectors in A and computing right singular */
04120 /*              vectors in VT */
04121 /*              (CWorkspace: 0) */
04122 /*              (RWorkspace: need BDSPAC) */
04123 
04124                 cbdsqr_("L", m, &ncvt, &nru, &c__0, &s[1], &rwork[ie], &vt[
04125                         vt_offset], ldvt, &a[a_offset], lda, cdum, &c__1, &
04126                         rwork[irwork], info);
04127             }
04128 
04129         }
04130 
04131     }
04132 
04133 /*     Undo scaling if necessary */
04134 
04135     if (iscl == 1) {
04136         if (anrm > bignum) {
04137             slascl_("G", &c__0, &c__0, &bignum, &anrm, &minmn, &c__1, &s[1], &
04138                     minmn, &ierr);
04139         }
04140         if (*info != 0 && anrm > bignum) {
04141             i__2 = minmn - 1;
04142             slascl_("G", &c__0, &c__0, &bignum, &anrm, &i__2, &c__1, &rwork[
04143                     ie], &minmn, &ierr);
04144         }
04145         if (anrm < smlnum) {
04146             slascl_("G", &c__0, &c__0, &smlnum, &anrm, &minmn, &c__1, &s[1], &
04147                     minmn, &ierr);
04148         }
04149         if (*info != 0 && anrm < smlnum) {
04150             i__2 = minmn - 1;
04151             slascl_("G", &c__0, &c__0, &smlnum, &anrm, &i__2, &c__1, &rwork[
04152                     ie], &minmn, &ierr);
04153         }
04154     }
04155 
04156 /*     Return optimal workspace in WORK(1) */
04157 
04158     work[1].r = (real) maxwrk, work[1].i = 0.f;
04159 
04160     return 0;
04161 
04162 /*     End of CGESVD */
04163 
04164 } /* cgesvd_ */