clattp.c
Go to the documentation of this file.
00001 /* clattp.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 integer c__5 = 5;
00019 static integer c__2 = 2;
00020 static integer c__1 = 1;
00021 static integer c__4 = 4;
00022 static real c_b93 = 2.f;
00023 
00024 /* Subroutine */ int clattp_(integer *imat, char *uplo, char *trans, char *
00025         diag, integer *iseed, integer *n, complex *ap, complex *b, complex *
00026         work, real *rwork, integer *info)
00027 {
00028     /* System generated locals */
00029     integer i__1, i__2, i__3, i__4, i__5;
00030     real r__1, r__2;
00031     doublereal d__1, d__2;
00032     complex q__1, q__2, q__3, q__4, q__5;
00033 
00034     /* Builtin functions */
00035     /* Subroutine */ int s_copy(char *, char *, ftnlen, ftnlen);
00036     void c_div(complex *, complex *, complex *);
00037     double pow_dd(doublereal *, doublereal *), sqrt(doublereal);
00038     void r_cnjg(complex *, complex *);
00039     double c_abs(complex *);
00040 
00041     /* Local variables */
00042     real c__;
00043     integer i__, j;
00044     complex s;
00045     real t, x, y, z__;
00046     integer jc;
00047     complex ra;
00048     integer jj;
00049     complex rb;
00050     integer jl, kl, jr, ku, iy, jx;
00051     real ulp, sfac;
00052     integer mode;
00053     char path[3], dist[1];
00054     real unfl;
00055     extern /* Subroutine */ int crot_(integer *, complex *, integer *, 
00056             complex *, integer *, real *, complex *);
00057     real rexp;
00058     char type__[1];
00059     real texp;
00060     complex star1, plus1, plus2;
00061     real bscal;
00062     extern logical lsame_(char *, char *);
00063     real tscal;
00064     complex ctemp;
00065     real anorm, bnorm, tleft;
00066     extern /* Subroutine */ int crotg_(complex *, complex *, real *, complex *
00067 );
00068     logical upper;
00069     extern /* Subroutine */ int clatb4_(char *, integer *, integer *, integer 
00070             *, char *, integer *, integer *, real *, integer *, real *, char *
00071 ), slabad_(real *, real *);
00072     extern integer icamax_(integer *, complex *, integer *);
00073     extern /* Complex */ VOID clarnd_(complex *, integer *, integer *);
00074     extern doublereal slamch_(char *);
00075     extern /* Subroutine */ int csscal_(integer *, real *, complex *, integer 
00076             *);
00077     char packit[1];
00078     real bignum;
00079     extern /* Subroutine */ int clatms_(integer *, integer *, char *, integer 
00080             *, char *, real *, integer *, real *, real *, integer *, integer *
00081 , char *, complex *, integer *, complex *, integer *);
00082     real cndnum;
00083     extern /* Subroutine */ int clarnv_(integer *, integer *, integer *, 
00084             complex *);
00085     integer jcnext, jcount;
00086     extern /* Subroutine */ int slarnv_(integer *, integer *, integer *, real 
00087             *);
00088     real smlnum;
00089 
00090 
00091 /*  -- LAPACK test routine (version 3.1) -- */
00092 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00093 /*     November 2006 */
00094 
00095 /*     .. Scalar Arguments .. */
00096 /*     .. */
00097 /*     .. Array Arguments .. */
00098 /*     .. */
00099 
00100 /*  Purpose */
00101 /*  ======= */
00102 
00103 /*  CLATTP generates a triangular test matrix in packed storage. */
00104 /*  IMAT and UPLO uniquely specify the properties of the test matrix, */
00105 /*  which is returned in the array AP. */
00106 
00107 /*  Arguments */
00108 /*  ========= */
00109 
00110 /*  IMAT    (input) INTEGER */
00111 /*          An integer key describing which matrix to generate for this */
00112 /*          path. */
00113 
00114 /*  UPLO    (input) CHARACTER*1 */
00115 /*          Specifies whether the matrix A will be upper or lower */
00116 /*          triangular. */
00117 /*          = 'U':  Upper triangular */
00118 /*          = 'L':  Lower triangular */
00119 
00120 /*  TRANS   (input) CHARACTER*1 */
00121 /*          Specifies whether the matrix or its transpose will be used. */
00122 /*          = 'N':  No transpose */
00123 /*          = 'T':  Transpose */
00124 /*          = 'C':  Conjugate transpose */
00125 
00126 /*  DIAG    (output) CHARACTER*1 */
00127 /*          Specifies whether or not the matrix A is unit triangular. */
00128 /*          = 'N':  Non-unit triangular */
00129 /*          = 'U':  Unit triangular */
00130 
00131 /*  ISEED   (input/output) INTEGER array, dimension (4) */
00132 /*          The seed vector for the random number generator (used in */
00133 /*          CLATMS).  Modified on exit. */
00134 
00135 /*  N       (input) INTEGER */
00136 /*          The order of the matrix to be generated. */
00137 
00138 /*  AP      (output) COMPLEX array, dimension (N*(N+1)/2) */
00139 /*          The upper or lower triangular matrix A, packed columnwise in */
00140 /*          a linear array.  The j-th column of A is stored in the array */
00141 /*          AP as follows: */
00142 /*          if UPLO = 'U', AP((j-1)*j/2 + i) = A(i,j) for 1<=i<=j; */
00143 /*          if UPLO = 'L', */
00144 /*             AP((j-1)*(n-j) + j*(j+1)/2 + i-j) = A(i,j) for j<=i<=n. */
00145 
00146 /*  B       (output) COMPLEX array, dimension (N) */
00147 /*          The right hand side vector, if IMAT > 10. */
00148 
00149 /*  WORK    (workspace) COMPLEX array, dimension (2*N) */
00150 
00151 /*  RWORK   (workspace) REAL array, dimension (N) */
00152 
00153 /*  INFO    (output) INTEGER */
00154 /*          = 0:  successful exit */
00155 /*          < 0:  if INFO = -i, the i-th argument had an illegal value */
00156 
00157 /*  ===================================================================== */
00158 
00159 /*     .. Parameters .. */
00160 /*     .. */
00161 /*     .. Local Scalars .. */
00162 /*     .. */
00163 /*     .. External Functions .. */
00164 /*     .. */
00165 /*     .. External Subroutines .. */
00166 /*     .. */
00167 /*     .. Intrinsic Functions .. */
00168 /*     .. */
00169 /*     .. Executable Statements .. */
00170 
00171     /* Parameter adjustments */
00172     --rwork;
00173     --work;
00174     --b;
00175     --ap;
00176     --iseed;
00177 
00178     /* Function Body */
00179     s_copy(path, "Complex precision", (ftnlen)1, (ftnlen)17);
00180     s_copy(path + 1, "TP", (ftnlen)2, (ftnlen)2);
00181     unfl = slamch_("Safe minimum");
00182     ulp = slamch_("Epsilon") * slamch_("Base");
00183     smlnum = unfl;
00184     bignum = (1.f - ulp) / smlnum;
00185     slabad_(&smlnum, &bignum);
00186     if (*imat >= 7 && *imat <= 10 || *imat == 18) {
00187         *(unsigned char *)diag = 'U';
00188     } else {
00189         *(unsigned char *)diag = 'N';
00190     }
00191     *info = 0;
00192 
00193 /*     Quick return if N.LE.0. */
00194 
00195     if (*n <= 0) {
00196         return 0;
00197     }
00198 
00199 /*     Call CLATB4 to set parameters for CLATMS. */
00200 
00201     upper = lsame_(uplo, "U");
00202     if (upper) {
00203         clatb4_(path, imat, n, n, type__, &kl, &ku, &anorm, &mode, &cndnum, 
00204                 dist);
00205         *(unsigned char *)packit = 'C';
00206     } else {
00207         i__1 = -(*imat);
00208         clatb4_(path, &i__1, n, n, type__, &kl, &ku, &anorm, &mode, &cndnum, 
00209                 dist);
00210         *(unsigned char *)packit = 'R';
00211     }
00212 
00213 /*     IMAT <= 6:  Non-unit triangular matrix */
00214 
00215     if (*imat <= 6) {
00216         clatms_(n, n, dist, &iseed[1], type__, &rwork[1], &mode, &cndnum, &
00217                 anorm, &kl, &ku, packit, &ap[1], n, &work[1], info);
00218 
00219 /*     IMAT > 6:  Unit triangular matrix */
00220 /*     The diagonal is deliberately set to something other than 1. */
00221 
00222 /*     IMAT = 7:  Matrix is the identity */
00223 
00224     } else if (*imat == 7) {
00225         if (upper) {
00226             jc = 1;
00227             i__1 = *n;
00228             for (j = 1; j <= i__1; ++j) {
00229                 i__2 = j - 1;
00230                 for (i__ = 1; i__ <= i__2; ++i__) {
00231                     i__3 = jc + i__ - 1;
00232                     ap[i__3].r = 0.f, ap[i__3].i = 0.f;
00233 /* L10: */
00234                 }
00235                 i__2 = jc + j - 1;
00236                 ap[i__2].r = (real) j, ap[i__2].i = 0.f;
00237                 jc += j;
00238 /* L20: */
00239             }
00240         } else {
00241             jc = 1;
00242             i__1 = *n;
00243             for (j = 1; j <= i__1; ++j) {
00244                 i__2 = jc;
00245                 ap[i__2].r = (real) j, ap[i__2].i = 0.f;
00246                 i__2 = *n;
00247                 for (i__ = j + 1; i__ <= i__2; ++i__) {
00248                     i__3 = jc + i__ - j;
00249                     ap[i__3].r = 0.f, ap[i__3].i = 0.f;
00250 /* L30: */
00251                 }
00252                 jc = jc + *n - j + 1;
00253 /* L40: */
00254             }
00255         }
00256 
00257 /*     IMAT > 7:  Non-trivial unit triangular matrix */
00258 
00259 /*     Generate a unit triangular matrix T with condition CNDNUM by */
00260 /*     forming a triangular matrix with known singular values and */
00261 /*     filling in the zero entries with Givens rotations. */
00262 
00263     } else if (*imat <= 10) {
00264         if (upper) {
00265             jc = 0;
00266             i__1 = *n;
00267             for (j = 1; j <= i__1; ++j) {
00268                 i__2 = j - 1;
00269                 for (i__ = 1; i__ <= i__2; ++i__) {
00270                     i__3 = jc + i__;
00271                     ap[i__3].r = 0.f, ap[i__3].i = 0.f;
00272 /* L50: */
00273                 }
00274                 i__2 = jc + j;
00275                 ap[i__2].r = (real) j, ap[i__2].i = 0.f;
00276                 jc += j;
00277 /* L60: */
00278             }
00279         } else {
00280             jc = 1;
00281             i__1 = *n;
00282             for (j = 1; j <= i__1; ++j) {
00283                 i__2 = jc;
00284                 ap[i__2].r = (real) j, ap[i__2].i = 0.f;
00285                 i__2 = *n;
00286                 for (i__ = j + 1; i__ <= i__2; ++i__) {
00287                     i__3 = jc + i__ - j;
00288                     ap[i__3].r = 0.f, ap[i__3].i = 0.f;
00289 /* L70: */
00290                 }
00291                 jc = jc + *n - j + 1;
00292 /* L80: */
00293             }
00294         }
00295 
00296 /*        Since the trace of a unit triangular matrix is 1, the product */
00297 /*        of its singular values must be 1.  Let s = sqrt(CNDNUM), */
00298 /*        x = sqrt(s) - 1/sqrt(s), y = sqrt(2/(n-2))*x, and z = x**2. */
00299 /*        The following triangular matrix has singular values s, 1, 1, */
00300 /*        ..., 1, 1/s: */
00301 
00302 /*        1  y  y  y  ...  y  y  z */
00303 /*           1  0  0  ...  0  0  y */
00304 /*              1  0  ...  0  0  y */
00305 /*                 .  ...  .  .  . */
00306 /*                     .   .  .  . */
00307 /*                         1  0  y */
00308 /*                            1  y */
00309 /*                               1 */
00310 
00311 /*        To fill in the zeros, we first multiply by a matrix with small */
00312 /*        condition number of the form */
00313 
00314 /*        1  0  0  0  0  ... */
00315 /*           1  +  *  0  0  ... */
00316 /*              1  +  0  0  0 */
00317 /*                 1  +  *  0  0 */
00318 /*                    1  +  0  0 */
00319 /*                       ... */
00320 /*                          1  +  0 */
00321 /*                             1  0 */
00322 /*                                1 */
00323 
00324 /*        Each element marked with a '*' is formed by taking the product */
00325 /*        of the adjacent elements marked with '+'.  The '*'s can be */
00326 /*        chosen freely, and the '+'s are chosen so that the inverse of */
00327 /*        T will have elements of the same magnitude as T.  If the *'s in */
00328 /*        both T and inv(T) have small magnitude, T is well conditioned. */
00329 /*        The two offdiagonals of T are stored in WORK. */
00330 
00331 /*        The product of these two matrices has the form */
00332 
00333 /*        1  y  y  y  y  y  .  y  y  z */
00334 /*           1  +  *  0  0  .  0  0  y */
00335 /*              1  +  0  0  .  0  0  y */
00336 /*                 1  +  *  .  .  .  . */
00337 /*                    1  +  .  .  .  . */
00338 /*                       .  .  .  .  . */
00339 /*                          .  .  .  . */
00340 /*                             1  +  y */
00341 /*                                1  y */
00342 /*                                   1 */
00343 
00344 /*        Now we multiply by Givens rotations, using the fact that */
00345 
00346 /*              [  c   s ] [  1   w ] [ -c  -s ] =  [  1  -w ] */
00347 /*              [ -s   c ] [  0   1 ] [  s  -c ]    [  0   1 ] */
00348 /*        and */
00349 /*              [ -c  -s ] [  1   0 ] [  c   s ] =  [  1   0 ] */
00350 /*              [  s  -c ] [  w   1 ] [ -s   c ]    [ -w   1 ] */
00351 
00352 /*        where c = w / sqrt(w**2+4) and s = 2 / sqrt(w**2+4). */
00353 
00354         clarnd_(&q__2, &c__5, &iseed[1]);
00355         q__1.r = q__2.r * .25f, q__1.i = q__2.i * .25f;
00356         star1.r = q__1.r, star1.i = q__1.i;
00357         sfac = .5f;
00358         clarnd_(&q__2, &c__5, &iseed[1]);
00359         q__1.r = sfac * q__2.r, q__1.i = sfac * q__2.i;
00360         plus1.r = q__1.r, plus1.i = q__1.i;
00361         i__1 = *n;
00362         for (j = 1; j <= i__1; j += 2) {
00363             c_div(&q__1, &star1, &plus1);
00364             plus2.r = q__1.r, plus2.i = q__1.i;
00365             i__2 = j;
00366             work[i__2].r = plus1.r, work[i__2].i = plus1.i;
00367             i__2 = *n + j;
00368             work[i__2].r = star1.r, work[i__2].i = star1.i;
00369             if (j + 1 <= *n) {
00370                 i__2 = j + 1;
00371                 work[i__2].r = plus2.r, work[i__2].i = plus2.i;
00372                 i__2 = *n + j + 1;
00373                 work[i__2].r = 0.f, work[i__2].i = 0.f;
00374                 c_div(&q__1, &star1, &plus2);
00375                 plus1.r = q__1.r, plus1.i = q__1.i;
00376                 clarnd_(&q__1, &c__2, &iseed[1]);
00377                 rexp = q__1.r;
00378                 if (rexp < 0.f) {
00379                     d__1 = (doublereal) sfac;
00380                     d__2 = (doublereal) (1.f - rexp);
00381                     r__1 = -pow_dd(&d__1, &d__2);
00382                     clarnd_(&q__2, &c__5, &iseed[1]);
00383                     q__1.r = r__1 * q__2.r, q__1.i = r__1 * q__2.i;
00384                     star1.r = q__1.r, star1.i = q__1.i;
00385                 } else {
00386                     d__1 = (doublereal) sfac;
00387                     d__2 = (doublereal) (rexp + 1.f);
00388                     r__1 = pow_dd(&d__1, &d__2);
00389                     clarnd_(&q__2, &c__5, &iseed[1]);
00390                     q__1.r = r__1 * q__2.r, q__1.i = r__1 * q__2.i;
00391                     star1.r = q__1.r, star1.i = q__1.i;
00392                 }
00393             }
00394 /* L90: */
00395         }
00396 
00397         x = sqrt(cndnum) - 1.f / sqrt(cndnum);
00398         if (*n > 2) {
00399             y = sqrt(2.f / (real) (*n - 2)) * x;
00400         } else {
00401             y = 0.f;
00402         }
00403         z__ = x * x;
00404 
00405         if (upper) {
00406 
00407 /*           Set the upper triangle of A with a unit triangular matrix */
00408 /*           of known condition number. */
00409 
00410             jc = 1;
00411             i__1 = *n;
00412             for (j = 2; j <= i__1; ++j) {
00413                 i__2 = jc + 1;
00414                 ap[i__2].r = y, ap[i__2].i = 0.f;
00415                 if (j > 2) {
00416                     i__2 = jc + j - 1;
00417                     i__3 = j - 2;
00418                     ap[i__2].r = work[i__3].r, ap[i__2].i = work[i__3].i;
00419                 }
00420                 if (j > 3) {
00421                     i__2 = jc + j - 2;
00422                     i__3 = *n + j - 3;
00423                     ap[i__2].r = work[i__3].r, ap[i__2].i = work[i__3].i;
00424                 }
00425                 jc += j;
00426 /* L100: */
00427             }
00428             jc -= *n;
00429             i__1 = jc + 1;
00430             ap[i__1].r = z__, ap[i__1].i = 0.f;
00431             i__1 = *n - 1;
00432             for (j = 2; j <= i__1; ++j) {
00433                 i__2 = jc + j;
00434                 ap[i__2].r = y, ap[i__2].i = 0.f;
00435 /* L110: */
00436             }
00437         } else {
00438 
00439 /*           Set the lower triangle of A with a unit triangular matrix */
00440 /*           of known condition number. */
00441 
00442             i__1 = *n - 1;
00443             for (i__ = 2; i__ <= i__1; ++i__) {
00444                 i__2 = i__;
00445                 ap[i__2].r = y, ap[i__2].i = 0.f;
00446 /* L120: */
00447             }
00448             i__1 = *n;
00449             ap[i__1].r = z__, ap[i__1].i = 0.f;
00450             jc = *n + 1;
00451             i__1 = *n - 1;
00452             for (j = 2; j <= i__1; ++j) {
00453                 i__2 = jc + 1;
00454                 i__3 = j - 1;
00455                 ap[i__2].r = work[i__3].r, ap[i__2].i = work[i__3].i;
00456                 if (j < *n - 1) {
00457                     i__2 = jc + 2;
00458                     i__3 = *n + j - 1;
00459                     ap[i__2].r = work[i__3].r, ap[i__2].i = work[i__3].i;
00460                 }
00461                 i__2 = jc + *n - j;
00462                 ap[i__2].r = y, ap[i__2].i = 0.f;
00463                 jc = jc + *n - j + 1;
00464 /* L130: */
00465             }
00466         }
00467 
00468 /*        Fill in the zeros using Givens rotations */
00469 
00470         if (upper) {
00471             jc = 1;
00472             i__1 = *n - 1;
00473             for (j = 1; j <= i__1; ++j) {
00474                 jcnext = jc + j;
00475                 i__2 = jcnext + j - 1;
00476                 ra.r = ap[i__2].r, ra.i = ap[i__2].i;
00477                 rb.r = 2.f, rb.i = 0.f;
00478                 crotg_(&ra, &rb, &c__, &s);
00479 
00480 /*              Multiply by [ c  s; -conjg(s)  c] on the left. */
00481 
00482                 if (*n > j + 1) {
00483                     jx = jcnext + j;
00484                     i__2 = *n;
00485                     for (i__ = j + 2; i__ <= i__2; ++i__) {
00486                         i__3 = jx + j;
00487                         q__2.r = c__ * ap[i__3].r, q__2.i = c__ * ap[i__3].i;
00488                         i__4 = jx + j + 1;
00489                         q__3.r = s.r * ap[i__4].r - s.i * ap[i__4].i, q__3.i =
00490                                  s.r * ap[i__4].i + s.i * ap[i__4].r;
00491                         q__1.r = q__2.r + q__3.r, q__1.i = q__2.i + q__3.i;
00492                         ctemp.r = q__1.r, ctemp.i = q__1.i;
00493                         i__3 = jx + j + 1;
00494                         r_cnjg(&q__4, &s);
00495                         q__3.r = -q__4.r, q__3.i = -q__4.i;
00496                         i__4 = jx + j;
00497                         q__2.r = q__3.r * ap[i__4].r - q__3.i * ap[i__4].i, 
00498                                 q__2.i = q__3.r * ap[i__4].i + q__3.i * ap[
00499                                 i__4].r;
00500                         i__5 = jx + j + 1;
00501                         q__5.r = c__ * ap[i__5].r, q__5.i = c__ * ap[i__5].i;
00502                         q__1.r = q__2.r + q__5.r, q__1.i = q__2.i + q__5.i;
00503                         ap[i__3].r = q__1.r, ap[i__3].i = q__1.i;
00504                         i__3 = jx + j;
00505                         ap[i__3].r = ctemp.r, ap[i__3].i = ctemp.i;
00506                         jx += i__;
00507 /* L140: */
00508                     }
00509                 }
00510 
00511 /*              Multiply by [-c -s;  conjg(s) -c] on the right. */
00512 
00513                 if (j > 1) {
00514                     i__2 = j - 1;
00515                     r__1 = -c__;
00516                     q__1.r = -s.r, q__1.i = -s.i;
00517                     crot_(&i__2, &ap[jcnext], &c__1, &ap[jc], &c__1, &r__1, &
00518                             q__1);
00519                 }
00520 
00521 /*              Negate A(J,J+1). */
00522 
00523                 i__2 = jcnext + j - 1;
00524                 i__3 = jcnext + j - 1;
00525                 q__1.r = -ap[i__3].r, q__1.i = -ap[i__3].i;
00526                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00527                 jc = jcnext;
00528 /* L150: */
00529             }
00530         } else {
00531             jc = 1;
00532             i__1 = *n - 1;
00533             for (j = 1; j <= i__1; ++j) {
00534                 jcnext = jc + *n - j + 1;
00535                 i__2 = jc + 1;
00536                 ra.r = ap[i__2].r, ra.i = ap[i__2].i;
00537                 rb.r = 2.f, rb.i = 0.f;
00538                 crotg_(&ra, &rb, &c__, &s);
00539                 r_cnjg(&q__1, &s);
00540                 s.r = q__1.r, s.i = q__1.i;
00541 
00542 /*              Multiply by [ c -s;  conjg(s) c] on the right. */
00543 
00544                 if (*n > j + 1) {
00545                     i__2 = *n - j - 1;
00546                     q__1.r = -s.r, q__1.i = -s.i;
00547                     crot_(&i__2, &ap[jcnext + 1], &c__1, &ap[jc + 2], &c__1, &
00548                             c__, &q__1);
00549                 }
00550 
00551 /*              Multiply by [-c  s; -conjg(s) -c] on the left. */
00552 
00553                 if (j > 1) {
00554                     jx = 1;
00555                     i__2 = j - 1;
00556                     for (i__ = 1; i__ <= i__2; ++i__) {
00557                         r__1 = -c__;
00558                         i__3 = jx + j - i__;
00559                         q__2.r = r__1 * ap[i__3].r, q__2.i = r__1 * ap[i__3]
00560                                 .i;
00561                         i__4 = jx + j - i__ + 1;
00562                         q__3.r = s.r * ap[i__4].r - s.i * ap[i__4].i, q__3.i =
00563                                  s.r * ap[i__4].i + s.i * ap[i__4].r;
00564                         q__1.r = q__2.r + q__3.r, q__1.i = q__2.i + q__3.i;
00565                         ctemp.r = q__1.r, ctemp.i = q__1.i;
00566                         i__3 = jx + j - i__ + 1;
00567                         r_cnjg(&q__4, &s);
00568                         q__3.r = -q__4.r, q__3.i = -q__4.i;
00569                         i__4 = jx + j - i__;
00570                         q__2.r = q__3.r * ap[i__4].r - q__3.i * ap[i__4].i, 
00571                                 q__2.i = q__3.r * ap[i__4].i + q__3.i * ap[
00572                                 i__4].r;
00573                         i__5 = jx + j - i__ + 1;
00574                         q__5.r = c__ * ap[i__5].r, q__5.i = c__ * ap[i__5].i;
00575                         q__1.r = q__2.r - q__5.r, q__1.i = q__2.i - q__5.i;
00576                         ap[i__3].r = q__1.r, ap[i__3].i = q__1.i;
00577                         i__3 = jx + j - i__;
00578                         ap[i__3].r = ctemp.r, ap[i__3].i = ctemp.i;
00579                         jx = jx + *n - i__ + 1;
00580 /* L160: */
00581                     }
00582                 }
00583 
00584 /*              Negate A(J+1,J). */
00585 
00586                 i__2 = jc + 1;
00587                 i__3 = jc + 1;
00588                 q__1.r = -ap[i__3].r, q__1.i = -ap[i__3].i;
00589                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00590                 jc = jcnext;
00591 /* L170: */
00592             }
00593         }
00594 
00595 /*     IMAT > 10:  Pathological test cases.  These triangular matrices */
00596 /*     are badly scaled or badly conditioned, so when used in solving a */
00597 /*     triangular system they may cause overflow in the solution vector. */
00598 
00599     } else if (*imat == 11) {
00600 
00601 /*        Type 11:  Generate a triangular matrix with elements between */
00602 /*        -1 and 1. Give the diagonal norm 2 to make it well-conditioned. */
00603 /*        Make the right hand side large so that it requires scaling. */
00604 
00605         if (upper) {
00606             jc = 1;
00607             i__1 = *n;
00608             for (j = 1; j <= i__1; ++j) {
00609                 i__2 = j - 1;
00610                 clarnv_(&c__4, &iseed[1], &i__2, &ap[jc]);
00611                 i__2 = jc + j - 1;
00612                 clarnd_(&q__2, &c__5, &iseed[1]);
00613                 q__1.r = q__2.r * 2.f, q__1.i = q__2.i * 2.f;
00614                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00615                 jc += j;
00616 /* L180: */
00617             }
00618         } else {
00619             jc = 1;
00620             i__1 = *n;
00621             for (j = 1; j <= i__1; ++j) {
00622                 if (j < *n) {
00623                     i__2 = *n - j;
00624                     clarnv_(&c__4, &iseed[1], &i__2, &ap[jc + 1]);
00625                 }
00626                 i__2 = jc;
00627                 clarnd_(&q__2, &c__5, &iseed[1]);
00628                 q__1.r = q__2.r * 2.f, q__1.i = q__2.i * 2.f;
00629                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00630                 jc = jc + *n - j + 1;
00631 /* L190: */
00632             }
00633         }
00634 
00635 /*        Set the right hand side so that the largest value is BIGNUM. */
00636 
00637         clarnv_(&c__2, &iseed[1], n, &b[1]);
00638         iy = icamax_(n, &b[1], &c__1);
00639         bnorm = c_abs(&b[iy]);
00640         bscal = bignum / dmax(1.f,bnorm);
00641         csscal_(n, &bscal, &b[1], &c__1);
00642 
00643     } else if (*imat == 12) {
00644 
00645 /*        Type 12:  Make the first diagonal element in the solve small to */
00646 /*        cause immediate overflow when dividing by T(j,j). */
00647 /*        In type 12, the offdiagonal elements are small (CNORM(j) < 1). */
00648 
00649         clarnv_(&c__2, &iseed[1], n, &b[1]);
00650 /* Computing MAX */
00651         r__1 = 1.f, r__2 = (real) (*n - 1);
00652         tscal = 1.f / dmax(r__1,r__2);
00653         if (upper) {
00654             jc = 1;
00655             i__1 = *n;
00656             for (j = 1; j <= i__1; ++j) {
00657                 i__2 = j - 1;
00658                 clarnv_(&c__4, &iseed[1], &i__2, &ap[jc]);
00659                 i__2 = j - 1;
00660                 csscal_(&i__2, &tscal, &ap[jc], &c__1);
00661                 i__2 = jc + j - 1;
00662                 clarnd_(&q__1, &c__5, &iseed[1]);
00663                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00664                 jc += j;
00665 /* L200: */
00666             }
00667             i__1 = *n * (*n + 1) / 2;
00668             i__2 = *n * (*n + 1) / 2;
00669             q__1.r = smlnum * ap[i__2].r, q__1.i = smlnum * ap[i__2].i;
00670             ap[i__1].r = q__1.r, ap[i__1].i = q__1.i;
00671         } else {
00672             jc = 1;
00673             i__1 = *n;
00674             for (j = 1; j <= i__1; ++j) {
00675                 i__2 = *n - j;
00676                 clarnv_(&c__2, &iseed[1], &i__2, &ap[jc + 1]);
00677                 i__2 = *n - j;
00678                 csscal_(&i__2, &tscal, &ap[jc + 1], &c__1);
00679                 i__2 = jc;
00680                 clarnd_(&q__1, &c__5, &iseed[1]);
00681                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00682                 jc = jc + *n - j + 1;
00683 /* L210: */
00684             }
00685             q__1.r = smlnum * ap[1].r, q__1.i = smlnum * ap[1].i;
00686             ap[1].r = q__1.r, ap[1].i = q__1.i;
00687         }
00688 
00689     } else if (*imat == 13) {
00690 
00691 /*        Type 13:  Make the first diagonal element in the solve small to */
00692 /*        cause immediate overflow when dividing by T(j,j). */
00693 /*        In type 13, the offdiagonal elements are O(1) (CNORM(j) > 1). */
00694 
00695         clarnv_(&c__2, &iseed[1], n, &b[1]);
00696         if (upper) {
00697             jc = 1;
00698             i__1 = *n;
00699             for (j = 1; j <= i__1; ++j) {
00700                 i__2 = j - 1;
00701                 clarnv_(&c__4, &iseed[1], &i__2, &ap[jc]);
00702                 i__2 = jc + j - 1;
00703                 clarnd_(&q__1, &c__5, &iseed[1]);
00704                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00705                 jc += j;
00706 /* L220: */
00707             }
00708             i__1 = *n * (*n + 1) / 2;
00709             i__2 = *n * (*n + 1) / 2;
00710             q__1.r = smlnum * ap[i__2].r, q__1.i = smlnum * ap[i__2].i;
00711             ap[i__1].r = q__1.r, ap[i__1].i = q__1.i;
00712         } else {
00713             jc = 1;
00714             i__1 = *n;
00715             for (j = 1; j <= i__1; ++j) {
00716                 i__2 = *n - j;
00717                 clarnv_(&c__4, &iseed[1], &i__2, &ap[jc + 1]);
00718                 i__2 = jc;
00719                 clarnd_(&q__1, &c__5, &iseed[1]);
00720                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00721                 jc = jc + *n - j + 1;
00722 /* L230: */
00723             }
00724             q__1.r = smlnum * ap[1].r, q__1.i = smlnum * ap[1].i;
00725             ap[1].r = q__1.r, ap[1].i = q__1.i;
00726         }
00727 
00728     } else if (*imat == 14) {
00729 
00730 /*        Type 14:  T is diagonal with small numbers on the diagonal to */
00731 /*        make the growth factor underflow, but a small right hand side */
00732 /*        chosen so that the solution does not overflow. */
00733 
00734         if (upper) {
00735             jcount = 1;
00736             jc = (*n - 1) * *n / 2 + 1;
00737             for (j = *n; j >= 1; --j) {
00738                 i__1 = j - 1;
00739                 for (i__ = 1; i__ <= i__1; ++i__) {
00740                     i__2 = jc + i__ - 1;
00741                     ap[i__2].r = 0.f, ap[i__2].i = 0.f;
00742 /* L240: */
00743                 }
00744                 if (jcount <= 2) {
00745                     i__1 = jc + j - 1;
00746                     clarnd_(&q__2, &c__5, &iseed[1]);
00747                     q__1.r = smlnum * q__2.r, q__1.i = smlnum * q__2.i;
00748                     ap[i__1].r = q__1.r, ap[i__1].i = q__1.i;
00749                 } else {
00750                     i__1 = jc + j - 1;
00751                     clarnd_(&q__1, &c__5, &iseed[1]);
00752                     ap[i__1].r = q__1.r, ap[i__1].i = q__1.i;
00753                 }
00754                 ++jcount;
00755                 if (jcount > 4) {
00756                     jcount = 1;
00757                 }
00758                 jc = jc - j + 1;
00759 /* L250: */
00760             }
00761         } else {
00762             jcount = 1;
00763             jc = 1;
00764             i__1 = *n;
00765             for (j = 1; j <= i__1; ++j) {
00766                 i__2 = *n;
00767                 for (i__ = j + 1; i__ <= i__2; ++i__) {
00768                     i__3 = jc + i__ - j;
00769                     ap[i__3].r = 0.f, ap[i__3].i = 0.f;
00770 /* L260: */
00771                 }
00772                 if (jcount <= 2) {
00773                     i__2 = jc;
00774                     clarnd_(&q__2, &c__5, &iseed[1]);
00775                     q__1.r = smlnum * q__2.r, q__1.i = smlnum * q__2.i;
00776                     ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00777                 } else {
00778                     i__2 = jc;
00779                     clarnd_(&q__1, &c__5, &iseed[1]);
00780                     ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00781                 }
00782                 ++jcount;
00783                 if (jcount > 4) {
00784                     jcount = 1;
00785                 }
00786                 jc = jc + *n - j + 1;
00787 /* L270: */
00788             }
00789         }
00790 
00791 /*        Set the right hand side alternately zero and small. */
00792 
00793         if (upper) {
00794             b[1].r = 0.f, b[1].i = 0.f;
00795             for (i__ = *n; i__ >= 2; i__ += -2) {
00796                 i__1 = i__;
00797                 b[i__1].r = 0.f, b[i__1].i = 0.f;
00798                 i__1 = i__ - 1;
00799                 clarnd_(&q__2, &c__5, &iseed[1]);
00800                 q__1.r = smlnum * q__2.r, q__1.i = smlnum * q__2.i;
00801                 b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00802 /* L280: */
00803             }
00804         } else {
00805             i__1 = *n;
00806             b[i__1].r = 0.f, b[i__1].i = 0.f;
00807             i__1 = *n - 1;
00808             for (i__ = 1; i__ <= i__1; i__ += 2) {
00809                 i__2 = i__;
00810                 b[i__2].r = 0.f, b[i__2].i = 0.f;
00811                 i__2 = i__ + 1;
00812                 clarnd_(&q__2, &c__5, &iseed[1]);
00813                 q__1.r = smlnum * q__2.r, q__1.i = smlnum * q__2.i;
00814                 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00815 /* L290: */
00816             }
00817         }
00818 
00819     } else if (*imat == 15) {
00820 
00821 /*        Type 15:  Make the diagonal elements small to cause gradual */
00822 /*        overflow when dividing by T(j,j).  To control the amount of */
00823 /*        scaling needed, the matrix is bidiagonal. */
00824 
00825 /* Computing MAX */
00826         r__1 = 1.f, r__2 = (real) (*n - 1);
00827         texp = 1.f / dmax(r__1,r__2);
00828         d__1 = (doublereal) smlnum;
00829         d__2 = (doublereal) texp;
00830         tscal = pow_dd(&d__1, &d__2);
00831         clarnv_(&c__4, &iseed[1], n, &b[1]);
00832         if (upper) {
00833             jc = 1;
00834             i__1 = *n;
00835             for (j = 1; j <= i__1; ++j) {
00836                 i__2 = j - 2;
00837                 for (i__ = 1; i__ <= i__2; ++i__) {
00838                     i__3 = jc + i__ - 1;
00839                     ap[i__3].r = 0.f, ap[i__3].i = 0.f;
00840 /* L300: */
00841                 }
00842                 if (j > 1) {
00843                     i__2 = jc + j - 2;
00844                     ap[i__2].r = -1.f, ap[i__2].i = -1.f;
00845                 }
00846                 i__2 = jc + j - 1;
00847                 clarnd_(&q__2, &c__5, &iseed[1]);
00848                 q__1.r = tscal * q__2.r, q__1.i = tscal * q__2.i;
00849                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00850                 jc += j;
00851 /* L310: */
00852             }
00853             i__1 = *n;
00854             b[i__1].r = 1.f, b[i__1].i = 1.f;
00855         } else {
00856             jc = 1;
00857             i__1 = *n;
00858             for (j = 1; j <= i__1; ++j) {
00859                 i__2 = *n;
00860                 for (i__ = j + 2; i__ <= i__2; ++i__) {
00861                     i__3 = jc + i__ - j;
00862                     ap[i__3].r = 0.f, ap[i__3].i = 0.f;
00863 /* L320: */
00864                 }
00865                 if (j < *n) {
00866                     i__2 = jc + 1;
00867                     ap[i__2].r = -1.f, ap[i__2].i = -1.f;
00868                 }
00869                 i__2 = jc;
00870                 clarnd_(&q__2, &c__5, &iseed[1]);
00871                 q__1.r = tscal * q__2.r, q__1.i = tscal * q__2.i;
00872                 ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00873                 jc = jc + *n - j + 1;
00874 /* L330: */
00875             }
00876             b[1].r = 1.f, b[1].i = 1.f;
00877         }
00878 
00879     } else if (*imat == 16) {
00880 
00881 /*        Type 16:  One zero diagonal element. */
00882 
00883         iy = *n / 2 + 1;
00884         if (upper) {
00885             jc = 1;
00886             i__1 = *n;
00887             for (j = 1; j <= i__1; ++j) {
00888                 clarnv_(&c__4, &iseed[1], &j, &ap[jc]);
00889                 if (j != iy) {
00890                     i__2 = jc + j - 1;
00891                     clarnd_(&q__2, &c__5, &iseed[1]);
00892                     q__1.r = q__2.r * 2.f, q__1.i = q__2.i * 2.f;
00893                     ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00894                 } else {
00895                     i__2 = jc + j - 1;
00896                     ap[i__2].r = 0.f, ap[i__2].i = 0.f;
00897                 }
00898                 jc += j;
00899 /* L340: */
00900             }
00901         } else {
00902             jc = 1;
00903             i__1 = *n;
00904             for (j = 1; j <= i__1; ++j) {
00905                 i__2 = *n - j + 1;
00906                 clarnv_(&c__4, &iseed[1], &i__2, &ap[jc]);
00907                 if (j != iy) {
00908                     i__2 = jc;
00909                     clarnd_(&q__2, &c__5, &iseed[1]);
00910                     q__1.r = q__2.r * 2.f, q__1.i = q__2.i * 2.f;
00911                     ap[i__2].r = q__1.r, ap[i__2].i = q__1.i;
00912                 } else {
00913                     i__2 = jc;
00914                     ap[i__2].r = 0.f, ap[i__2].i = 0.f;
00915                 }
00916                 jc = jc + *n - j + 1;
00917 /* L350: */
00918             }
00919         }
00920         clarnv_(&c__2, &iseed[1], n, &b[1]);
00921         csscal_(n, &c_b93, &b[1], &c__1);
00922 
00923     } else if (*imat == 17) {
00924 
00925 /*        Type 17:  Make the offdiagonal elements large to cause overflow */
00926 /*        when adding a column of T.  In the non-transposed case, the */
00927 /*        matrix is constructed to cause overflow when adding a column in */
00928 /*        every other step. */
00929 
00930         tscal = unfl / ulp;
00931         tscal = (1.f - ulp) / tscal;
00932         i__1 = *n * (*n + 1) / 2;
00933         for (j = 1; j <= i__1; ++j) {
00934             i__2 = j;
00935             ap[i__2].r = 0.f, ap[i__2].i = 0.f;
00936 /* L360: */
00937         }
00938         texp = 1.f;
00939         if (upper) {
00940             jc = (*n - 1) * *n / 2 + 1;
00941             for (j = *n; j >= 2; j += -2) {
00942                 i__1 = jc;
00943                 r__1 = -tscal / (real) (*n + 1);
00944                 ap[i__1].r = r__1, ap[i__1].i = 0.f;
00945                 i__1 = jc + j - 1;
00946                 ap[i__1].r = 1.f, ap[i__1].i = 0.f;
00947                 i__1 = j;
00948                 r__1 = texp * (1.f - ulp);
00949                 b[i__1].r = r__1, b[i__1].i = 0.f;
00950                 jc = jc - j + 1;
00951                 i__1 = jc;
00952                 r__1 = -(tscal / (real) (*n + 1)) / (real) (*n + 2);
00953                 ap[i__1].r = r__1, ap[i__1].i = 0.f;
00954                 i__1 = jc + j - 2;
00955                 ap[i__1].r = 1.f, ap[i__1].i = 0.f;
00956                 i__1 = j - 1;
00957                 r__1 = texp * (real) (*n * *n + *n - 1);
00958                 b[i__1].r = r__1, b[i__1].i = 0.f;
00959                 texp *= 2.f;
00960                 jc = jc - j + 2;
00961 /* L370: */
00962             }
00963             r__1 = (real) (*n + 1) / (real) (*n + 2) * tscal;
00964             b[1].r = r__1, b[1].i = 0.f;
00965         } else {
00966             jc = 1;
00967             i__1 = *n - 1;
00968             for (j = 1; j <= i__1; j += 2) {
00969                 i__2 = jc + *n - j;
00970                 r__1 = -tscal / (real) (*n + 1);
00971                 ap[i__2].r = r__1, ap[i__2].i = 0.f;
00972                 i__2 = jc;
00973                 ap[i__2].r = 1.f, ap[i__2].i = 0.f;
00974                 i__2 = j;
00975                 r__1 = texp * (1.f - ulp);
00976                 b[i__2].r = r__1, b[i__2].i = 0.f;
00977                 jc = jc + *n - j + 1;
00978                 i__2 = jc + *n - j - 1;
00979                 r__1 = -(tscal / (real) (*n + 1)) / (real) (*n + 2);
00980                 ap[i__2].r = r__1, ap[i__2].i = 0.f;
00981                 i__2 = jc;
00982                 ap[i__2].r = 1.f, ap[i__2].i = 0.f;
00983                 i__2 = j + 1;
00984                 r__1 = texp * (real) (*n * *n + *n - 1);
00985                 b[i__2].r = r__1, b[i__2].i = 0.f;
00986                 texp *= 2.f;
00987                 jc = jc + *n - j;
00988 /* L380: */
00989             }
00990             i__1 = *n;
00991             r__1 = (real) (*n + 1) / (real) (*n + 2) * tscal;
00992             b[i__1].r = r__1, b[i__1].i = 0.f;
00993         }
00994 
00995     } else if (*imat == 18) {
00996 
00997 /*        Type 18:  Generate a unit triangular matrix with elements */
00998 /*        between -1 and 1, and make the right hand side large so that it */
00999 /*        requires scaling. */
01000 
01001         if (upper) {
01002             jc = 1;
01003             i__1 = *n;
01004             for (j = 1; j <= i__1; ++j) {
01005                 i__2 = j - 1;
01006                 clarnv_(&c__4, &iseed[1], &i__2, &ap[jc]);
01007                 i__2 = jc + j - 1;
01008                 ap[i__2].r = 0.f, ap[i__2].i = 0.f;
01009                 jc += j;
01010 /* L390: */
01011             }
01012         } else {
01013             jc = 1;
01014             i__1 = *n;
01015             for (j = 1; j <= i__1; ++j) {
01016                 if (j < *n) {
01017                     i__2 = *n - j;
01018                     clarnv_(&c__4, &iseed[1], &i__2, &ap[jc + 1]);
01019                 }
01020                 i__2 = jc;
01021                 ap[i__2].r = 0.f, ap[i__2].i = 0.f;
01022                 jc = jc + *n - j + 1;
01023 /* L400: */
01024             }
01025         }
01026 
01027 /*        Set the right hand side so that the largest value is BIGNUM. */
01028 
01029         clarnv_(&c__2, &iseed[1], n, &b[1]);
01030         iy = icamax_(n, &b[1], &c__1);
01031         bnorm = c_abs(&b[iy]);
01032         bscal = bignum / dmax(1.f,bnorm);
01033         csscal_(n, &bscal, &b[1], &c__1);
01034 
01035     } else if (*imat == 19) {
01036 
01037 /*        Type 19:  Generate a triangular matrix with elements between */
01038 /*        BIGNUM/(n-1) and BIGNUM so that at least one of the column */
01039 /*        norms will exceed BIGNUM. */
01040 /*        1/3/91:  CLATPS no longer can handle this case */
01041 
01042 /* Computing MAX */
01043         r__1 = 1.f, r__2 = (real) (*n - 1);
01044         tleft = bignum / dmax(r__1,r__2);
01045 /* Computing MAX */
01046         r__1 = 1.f, r__2 = (real) (*n);
01047         tscal = bignum * ((real) (*n - 1) / dmax(r__1,r__2));
01048         if (upper) {
01049             jc = 1;
01050             i__1 = *n;
01051             for (j = 1; j <= i__1; ++j) {
01052                 clarnv_(&c__5, &iseed[1], &j, &ap[jc]);
01053                 slarnv_(&c__1, &iseed[1], &j, &rwork[1]);
01054                 i__2 = j;
01055                 for (i__ = 1; i__ <= i__2; ++i__) {
01056                     i__3 = jc + i__ - 1;
01057                     i__4 = jc + i__ - 1;
01058                     r__1 = tleft + rwork[i__] * tscal;
01059                     q__1.r = r__1 * ap[i__4].r, q__1.i = r__1 * ap[i__4].i;
01060                     ap[i__3].r = q__1.r, ap[i__3].i = q__1.i;
01061 /* L410: */
01062                 }
01063                 jc += j;
01064 /* L420: */
01065             }
01066         } else {
01067             jc = 1;
01068             i__1 = *n;
01069             for (j = 1; j <= i__1; ++j) {
01070                 i__2 = *n - j + 1;
01071                 clarnv_(&c__5, &iseed[1], &i__2, &ap[jc]);
01072                 i__2 = *n - j + 1;
01073                 slarnv_(&c__1, &iseed[1], &i__2, &rwork[1]);
01074                 i__2 = *n;
01075                 for (i__ = j; i__ <= i__2; ++i__) {
01076                     i__3 = jc + i__ - j;
01077                     i__4 = jc + i__ - j;
01078                     r__1 = tleft + rwork[i__ - j + 1] * tscal;
01079                     q__1.r = r__1 * ap[i__4].r, q__1.i = r__1 * ap[i__4].i;
01080                     ap[i__3].r = q__1.r, ap[i__3].i = q__1.i;
01081 /* L430: */
01082                 }
01083                 jc = jc + *n - j + 1;
01084 /* L440: */
01085             }
01086         }
01087         clarnv_(&c__2, &iseed[1], n, &b[1]);
01088         csscal_(n, &c_b93, &b[1], &c__1);
01089     }
01090 
01091 /*     Flip the matrix across its counter-diagonal if the transpose will */
01092 /*     be used. */
01093 
01094     if (! lsame_(trans, "N")) {
01095         if (upper) {
01096             jj = 1;
01097             jr = *n * (*n + 1) / 2;
01098             i__1 = *n / 2;
01099             for (j = 1; j <= i__1; ++j) {
01100                 jl = jj;
01101                 i__2 = *n - j;
01102                 for (i__ = j; i__ <= i__2; ++i__) {
01103                     i__3 = jr - i__ + j;
01104                     t = ap[i__3].r;
01105                     i__3 = jr - i__ + j;
01106                     i__4 = jl;
01107                     ap[i__3].r = ap[i__4].r, ap[i__3].i = ap[i__4].i;
01108                     i__3 = jl;
01109                     ap[i__3].r = t, ap[i__3].i = 0.f;
01110                     jl += i__;
01111 /* L450: */
01112                 }
01113                 jj = jj + j + 1;
01114                 jr -= *n - j + 1;
01115 /* L460: */
01116             }
01117         } else {
01118             jl = 1;
01119             jj = *n * (*n + 1) / 2;
01120             i__1 = *n / 2;
01121             for (j = 1; j <= i__1; ++j) {
01122                 jr = jj;
01123                 i__2 = *n - j;
01124                 for (i__ = j; i__ <= i__2; ++i__) {
01125                     i__3 = jl + i__ - j;
01126                     t = ap[i__3].r;
01127                     i__3 = jl + i__ - j;
01128                     i__4 = jr;
01129                     ap[i__3].r = ap[i__4].r, ap[i__3].i = ap[i__4].i;
01130                     i__3 = jr;
01131                     ap[i__3].r = t, ap[i__3].i = 0.f;
01132                     jr -= i__;
01133 /* L470: */
01134                 }
01135                 jl = jl + *n - j + 1;
01136                 jj = jj - j - 1;
01137 /* L480: */
01138             }
01139         }
01140     }
01141 
01142     return 0;
01143 
01144 /*     End of CLATTP */
01145 
01146 } /* clattp_ */


swiftnav
Author(s):
autogenerated on Sat Jun 8 2019 18:55:32