zlatsy.c
Go to the documentation of this file.
00001 /* zlatsy.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 
00021 /* Subroutine */ int zlatsy_(char *uplo, integer *n, doublecomplex *x, 
00022         integer *ldx, integer *iseed)
00023 {
00024     /* System generated locals */
00025     integer x_dim1, x_offset, i__1, i__2, i__3;
00026     doublecomplex z__1, z__2, z__3;
00027 
00028     /* Builtin functions */
00029     double sqrt(doublereal), z_abs(doublecomplex *);
00030 
00031     /* Local variables */
00032     doublecomplex a, b, c__;
00033     integer i__, j;
00034     doublecomplex r__;
00035     integer n5;
00036     doublereal beta, alpha, alpha3;
00037     extern /* Double Complex */ VOID zlarnd_(doublecomplex *, integer *, 
00038             integer *);
00039 
00040 
00041 /*  -- LAPACK auxiliary test routine (version 3.1) -- */
00042 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00043 /*     November 2006 */
00044 
00045 /*     .. Scalar Arguments .. */
00046 /*     .. */
00047 /*     .. Array Arguments .. */
00048 /*     .. */
00049 
00050 /*  Purpose */
00051 /*  ======= */
00052 
00053 /*  ZLATSY generates a special test matrix for the complex symmetric */
00054 /*  (indefinite) factorization.  The pivot blocks of the generated matrix */
00055 /*  will be in the following order: */
00056 /*     2x2 pivot block, non diagonalizable */
00057 /*     1x1 pivot block */
00058 /*     2x2 pivot block, diagonalizable */
00059 /*     (cycle repeats) */
00060 /*  A row interchange is required for each non-diagonalizable 2x2 block. */
00061 
00062 /*  Arguments */
00063 /*  ========= */
00064 
00065 /*  UPLO    (input) CHARACTER */
00066 /*          Specifies whether the generated matrix is to be upper or */
00067 /*          lower triangular. */
00068 /*          = 'U':  Upper triangular */
00069 /*          = 'L':  Lower triangular */
00070 
00071 /*  N       (input) INTEGER */
00072 /*          The dimension of the matrix to be generated. */
00073 
00074 /*  X       (output) COMPLEX*16 array, dimension (LDX,N) */
00075 /*          The generated matrix, consisting of 3x3 and 2x2 diagonal */
00076 /*          blocks which result in the pivot sequence given above. */
00077 /*          The matrix outside of these diagonal blocks is zero. */
00078 
00079 /*  LDX     (input) INTEGER */
00080 /*          The leading dimension of the array X. */
00081 
00082 /*  ISEED   (input/output) INTEGER array, dimension (4) */
00083 /*          On entry, the seed for the random number generator.  The last */
00084 /*          of the four integers must be odd.  (modified on exit) */
00085 
00086 /*  ===================================================================== */
00087 
00088 /*     .. Parameters .. */
00089 /*     .. */
00090 /*     .. Local Scalars .. */
00091 /*     .. */
00092 /*     .. External Functions .. */
00093 /*     .. */
00094 /*     .. Intrinsic Functions .. */
00095 /*     .. */
00096 /*     .. Executable Statements .. */
00097 
00098 /*     Initialize constants */
00099 
00100     /* Parameter adjustments */
00101     x_dim1 = *ldx;
00102     x_offset = 1 + x_dim1;
00103     x -= x_offset;
00104     --iseed;
00105 
00106     /* Function Body */
00107     alpha = (sqrt(17.) + 1.) / 8.;
00108     beta = alpha - .001;
00109     alpha3 = alpha * alpha * alpha;
00110 
00111 /*     UPLO = 'U':  Upper triangular storage */
00112 
00113     if (*(unsigned char *)uplo == 'U') {
00114 
00115 /*        Fill the upper triangle of the matrix with zeros. */
00116 
00117         i__1 = *n;
00118         for (j = 1; j <= i__1; ++j) {
00119             i__2 = j;
00120             for (i__ = 1; i__ <= i__2; ++i__) {
00121                 i__3 = i__ + j * x_dim1;
00122                 x[i__3].r = 0., x[i__3].i = 0.;
00123 /* L10: */
00124             }
00125 /* L20: */
00126         }
00127         n5 = *n / 5;
00128         n5 = *n - n5 * 5 + 1;
00129 
00130         i__1 = n5;
00131         for (i__ = *n; i__ >= i__1; i__ += -5) {
00132             zlarnd_(&z__2, &c__5, &iseed[1]);
00133             z__1.r = alpha3 * z__2.r, z__1.i = alpha3 * z__2.i;
00134             a.r = z__1.r, a.i = z__1.i;
00135             zlarnd_(&z__2, &c__5, &iseed[1]);
00136             z__1.r = z__2.r / alpha, z__1.i = z__2.i / alpha;
00137             b.r = z__1.r, b.i = z__1.i;
00138             z__3.r = b.r * 2., z__3.i = b.i * 2.;
00139             z__2.r = z__3.r * 0. - z__3.i * 1., z__2.i = z__3.r * 1. + z__3.i 
00140                     * 0.;
00141             z__1.r = a.r - z__2.r, z__1.i = a.i - z__2.i;
00142             c__.r = z__1.r, c__.i = z__1.i;
00143             z__1.r = c__.r / beta, z__1.i = c__.i / beta;
00144             r__.r = z__1.r, r__.i = z__1.i;
00145             i__2 = i__ + i__ * x_dim1;
00146             x[i__2].r = a.r, x[i__2].i = a.i;
00147             i__2 = i__ - 2 + i__ * x_dim1;
00148             x[i__2].r = b.r, x[i__2].i = b.i;
00149             i__2 = i__ - 2 + (i__ - 1) * x_dim1;
00150             x[i__2].r = r__.r, x[i__2].i = r__.i;
00151             i__2 = i__ - 2 + (i__ - 2) * x_dim1;
00152             x[i__2].r = c__.r, x[i__2].i = c__.i;
00153             i__2 = i__ - 1 + (i__ - 1) * x_dim1;
00154             zlarnd_(&z__1, &c__2, &iseed[1]);
00155             x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00156             i__2 = i__ - 3 + (i__ - 3) * x_dim1;
00157             zlarnd_(&z__1, &c__2, &iseed[1]);
00158             x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00159             i__2 = i__ - 4 + (i__ - 4) * x_dim1;
00160             zlarnd_(&z__1, &c__2, &iseed[1]);
00161             x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00162             if (z_abs(&x[i__ - 3 + (i__ - 3) * x_dim1]) > z_abs(&x[i__ - 4 + (
00163                     i__ - 4) * x_dim1])) {
00164                 i__2 = i__ - 4 + (i__ - 3) * x_dim1;
00165                 i__3 = i__ - 3 + (i__ - 3) * x_dim1;
00166                 z__1.r = x[i__3].r * 2., z__1.i = x[i__3].i * 2.;
00167                 x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00168             } else {
00169                 i__2 = i__ - 4 + (i__ - 3) * x_dim1;
00170                 i__3 = i__ - 4 + (i__ - 4) * x_dim1;
00171                 z__1.r = x[i__3].r * 2., z__1.i = x[i__3].i * 2.;
00172                 x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00173             }
00174 /* L30: */
00175         }
00176 
00177 /*        Clean-up for N not a multiple of 5. */
00178 
00179         i__ = n5 - 1;
00180         if (i__ > 2) {
00181             zlarnd_(&z__2, &c__5, &iseed[1]);
00182             z__1.r = alpha3 * z__2.r, z__1.i = alpha3 * z__2.i;
00183             a.r = z__1.r, a.i = z__1.i;
00184             zlarnd_(&z__2, &c__5, &iseed[1]);
00185             z__1.r = z__2.r / alpha, z__1.i = z__2.i / alpha;
00186             b.r = z__1.r, b.i = z__1.i;
00187             z__3.r = b.r * 2., z__3.i = b.i * 2.;
00188             z__2.r = z__3.r * 0. - z__3.i * 1., z__2.i = z__3.r * 1. + z__3.i 
00189                     * 0.;
00190             z__1.r = a.r - z__2.r, z__1.i = a.i - z__2.i;
00191             c__.r = z__1.r, c__.i = z__1.i;
00192             z__1.r = c__.r / beta, z__1.i = c__.i / beta;
00193             r__.r = z__1.r, r__.i = z__1.i;
00194             i__1 = i__ + i__ * x_dim1;
00195             x[i__1].r = a.r, x[i__1].i = a.i;
00196             i__1 = i__ - 2 + i__ * x_dim1;
00197             x[i__1].r = b.r, x[i__1].i = b.i;
00198             i__1 = i__ - 2 + (i__ - 1) * x_dim1;
00199             x[i__1].r = r__.r, x[i__1].i = r__.i;
00200             i__1 = i__ - 2 + (i__ - 2) * x_dim1;
00201             x[i__1].r = c__.r, x[i__1].i = c__.i;
00202             i__1 = i__ - 1 + (i__ - 1) * x_dim1;
00203             zlarnd_(&z__1, &c__2, &iseed[1]);
00204             x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00205             i__ += -3;
00206         }
00207         if (i__ > 1) {
00208             i__1 = i__ + i__ * x_dim1;
00209             zlarnd_(&z__1, &c__2, &iseed[1]);
00210             x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00211             i__1 = i__ - 1 + (i__ - 1) * x_dim1;
00212             zlarnd_(&z__1, &c__2, &iseed[1]);
00213             x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00214             if (z_abs(&x[i__ + i__ * x_dim1]) > z_abs(&x[i__ - 1 + (i__ - 1) *
00215                      x_dim1])) {
00216                 i__1 = i__ - 1 + i__ * x_dim1;
00217                 i__2 = i__ + i__ * x_dim1;
00218                 z__1.r = x[i__2].r * 2., z__1.i = x[i__2].i * 2.;
00219                 x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00220             } else {
00221                 i__1 = i__ - 1 + i__ * x_dim1;
00222                 i__2 = i__ - 1 + (i__ - 1) * x_dim1;
00223                 z__1.r = x[i__2].r * 2., z__1.i = x[i__2].i * 2.;
00224                 x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00225             }
00226             i__ += -2;
00227         } else if (i__ == 1) {
00228             i__1 = i__ + i__ * x_dim1;
00229             zlarnd_(&z__1, &c__2, &iseed[1]);
00230             x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00231             --i__;
00232         }
00233 
00234 /*     UPLO = 'L':  Lower triangular storage */
00235 
00236     } else {
00237 
00238 /*        Fill the lower triangle of the matrix with zeros. */
00239 
00240         i__1 = *n;
00241         for (j = 1; j <= i__1; ++j) {
00242             i__2 = *n;
00243             for (i__ = j; i__ <= i__2; ++i__) {
00244                 i__3 = i__ + j * x_dim1;
00245                 x[i__3].r = 0., x[i__3].i = 0.;
00246 /* L40: */
00247             }
00248 /* L50: */
00249         }
00250         n5 = *n / 5;
00251         n5 *= 5;
00252 
00253         i__1 = n5;
00254         for (i__ = 1; i__ <= i__1; i__ += 5) {
00255             zlarnd_(&z__2, &c__5, &iseed[1]);
00256             z__1.r = alpha3 * z__2.r, z__1.i = alpha3 * z__2.i;
00257             a.r = z__1.r, a.i = z__1.i;
00258             zlarnd_(&z__2, &c__5, &iseed[1]);
00259             z__1.r = z__2.r / alpha, z__1.i = z__2.i / alpha;
00260             b.r = z__1.r, b.i = z__1.i;
00261             z__3.r = b.r * 2., z__3.i = b.i * 2.;
00262             z__2.r = z__3.r * 0. - z__3.i * 1., z__2.i = z__3.r * 1. + z__3.i 
00263                     * 0.;
00264             z__1.r = a.r - z__2.r, z__1.i = a.i - z__2.i;
00265             c__.r = z__1.r, c__.i = z__1.i;
00266             z__1.r = c__.r / beta, z__1.i = c__.i / beta;
00267             r__.r = z__1.r, r__.i = z__1.i;
00268             i__2 = i__ + i__ * x_dim1;
00269             x[i__2].r = a.r, x[i__2].i = a.i;
00270             i__2 = i__ + 2 + i__ * x_dim1;
00271             x[i__2].r = b.r, x[i__2].i = b.i;
00272             i__2 = i__ + 2 + (i__ + 1) * x_dim1;
00273             x[i__2].r = r__.r, x[i__2].i = r__.i;
00274             i__2 = i__ + 2 + (i__ + 2) * x_dim1;
00275             x[i__2].r = c__.r, x[i__2].i = c__.i;
00276             i__2 = i__ + 1 + (i__ + 1) * x_dim1;
00277             zlarnd_(&z__1, &c__2, &iseed[1]);
00278             x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00279             i__2 = i__ + 3 + (i__ + 3) * x_dim1;
00280             zlarnd_(&z__1, &c__2, &iseed[1]);
00281             x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00282             i__2 = i__ + 4 + (i__ + 4) * x_dim1;
00283             zlarnd_(&z__1, &c__2, &iseed[1]);
00284             x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00285             if (z_abs(&x[i__ + 3 + (i__ + 3) * x_dim1]) > z_abs(&x[i__ + 4 + (
00286                     i__ + 4) * x_dim1])) {
00287                 i__2 = i__ + 4 + (i__ + 3) * x_dim1;
00288                 i__3 = i__ + 3 + (i__ + 3) * x_dim1;
00289                 z__1.r = x[i__3].r * 2., z__1.i = x[i__3].i * 2.;
00290                 x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00291             } else {
00292                 i__2 = i__ + 4 + (i__ + 3) * x_dim1;
00293                 i__3 = i__ + 4 + (i__ + 4) * x_dim1;
00294                 z__1.r = x[i__3].r * 2., z__1.i = x[i__3].i * 2.;
00295                 x[i__2].r = z__1.r, x[i__2].i = z__1.i;
00296             }
00297 /* L60: */
00298         }
00299 
00300 /*        Clean-up for N not a multiple of 5. */
00301 
00302         i__ = n5 + 1;
00303         if (i__ < *n - 1) {
00304             zlarnd_(&z__2, &c__5, &iseed[1]);
00305             z__1.r = alpha3 * z__2.r, z__1.i = alpha3 * z__2.i;
00306             a.r = z__1.r, a.i = z__1.i;
00307             zlarnd_(&z__2, &c__5, &iseed[1]);
00308             z__1.r = z__2.r / alpha, z__1.i = z__2.i / alpha;
00309             b.r = z__1.r, b.i = z__1.i;
00310             z__3.r = b.r * 2., z__3.i = b.i * 2.;
00311             z__2.r = z__3.r * 0. - z__3.i * 1., z__2.i = z__3.r * 1. + z__3.i 
00312                     * 0.;
00313             z__1.r = a.r - z__2.r, z__1.i = a.i - z__2.i;
00314             c__.r = z__1.r, c__.i = z__1.i;
00315             z__1.r = c__.r / beta, z__1.i = c__.i / beta;
00316             r__.r = z__1.r, r__.i = z__1.i;
00317             i__1 = i__ + i__ * x_dim1;
00318             x[i__1].r = a.r, x[i__1].i = a.i;
00319             i__1 = i__ + 2 + i__ * x_dim1;
00320             x[i__1].r = b.r, x[i__1].i = b.i;
00321             i__1 = i__ + 2 + (i__ + 1) * x_dim1;
00322             x[i__1].r = r__.r, x[i__1].i = r__.i;
00323             i__1 = i__ + 2 + (i__ + 2) * x_dim1;
00324             x[i__1].r = c__.r, x[i__1].i = c__.i;
00325             i__1 = i__ + 1 + (i__ + 1) * x_dim1;
00326             zlarnd_(&z__1, &c__2, &iseed[1]);
00327             x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00328             i__ += 3;
00329         }
00330         if (i__ < *n) {
00331             i__1 = i__ + i__ * x_dim1;
00332             zlarnd_(&z__1, &c__2, &iseed[1]);
00333             x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00334             i__1 = i__ + 1 + (i__ + 1) * x_dim1;
00335             zlarnd_(&z__1, &c__2, &iseed[1]);
00336             x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00337             if (z_abs(&x[i__ + i__ * x_dim1]) > z_abs(&x[i__ + 1 + (i__ + 1) *
00338                      x_dim1])) {
00339                 i__1 = i__ + 1 + i__ * x_dim1;
00340                 i__2 = i__ + i__ * x_dim1;
00341                 z__1.r = x[i__2].r * 2., z__1.i = x[i__2].i * 2.;
00342                 x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00343             } else {
00344                 i__1 = i__ + 1 + i__ * x_dim1;
00345                 i__2 = i__ + 1 + (i__ + 1) * x_dim1;
00346                 z__1.r = x[i__2].r * 2., z__1.i = x[i__2].i * 2.;
00347                 x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00348             }
00349             i__ += 2;
00350         } else if (i__ == *n) {
00351             i__1 = i__ + i__ * x_dim1;
00352             zlarnd_(&z__1, &c__2, &iseed[1]);
00353             x[i__1].r = z__1.r, x[i__1].i = z__1.i;
00354             ++i__;
00355         }
00356     }
00357 
00358     return 0;
00359 
00360 /*     End of ZLATSY */
00361 
00362 } /* zlatsy_ */

swiftnav
Author(s):
autogenerated on Sat Jun 8 2019 18:56:42