cchkbb.c
Go to the documentation of this file.
00001 /* cchkbb.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__0 = 0;
00021 static integer c__6 = 6;
00022 static real c_b33 = 1.f;
00023 static integer c__1 = 1;
00024 static real c_b41 = 0.f;
00025 static integer c__4 = 4;
00026 static integer c_n1 = -1;
00027 
00028 /* Subroutine */ int cchkbb_(integer *nsizes, integer *mval, integer *nval, 
00029         integer *nwdths, integer *kk, integer *ntypes, logical *dotype, 
00030         integer *nrhs, integer *iseed, real *thresh, integer *nounit, complex 
00031         *a, integer *lda, complex *ab, integer *ldab, real *bd, real *be, 
00032         complex *q, integer *ldq, complex *p, integer *ldp, complex *c__, 
00033         integer *ldc, complex *cc, complex *work, integer *lwork, real *rwork, 
00034          real *result, integer *info)
00035 {
00036     /* Initialized data */
00037 
00038     static integer ktype[15] = { 1,2,4,4,4,4,4,6,6,6,6,6,9,9,9 };
00039     static integer kmagn[15] = { 1,1,1,1,1,2,3,1,1,1,2,3,1,2,3 };
00040     static integer kmode[15] = { 0,0,4,3,1,4,4,4,3,1,4,4,0,0,0 };
00041 
00042     /* Format strings */
00043     static char fmt_9999[] = "(\002 CCHKBB: \002,a,\002 returned INFO=\002,i"
00044             "5,\002.\002,/9x,\002M=\002,i5,\002 N=\002,i5,\002 K=\002,i5,\002"
00045             ", JTYPE=\002,i5,\002, ISEED=(\002,3(i5,\002,\002),i5,\002)\002)";
00046     static char fmt_9998[] = "(\002 M =\002,i4,\002 N=\002,i4,\002, K=\002,i"
00047             "3,\002, seed=\002,4(i4,\002,\002),\002 type \002,i2,\002, test"
00048             "(\002,i2,\002)=\002,g10.3)";
00049 
00050     /* System generated locals */
00051     integer a_dim1, a_offset, ab_dim1, ab_offset, c_dim1, c_offset, cc_dim1, 
00052             cc_offset, p_dim1, p_offset, q_dim1, q_offset, i__1, i__2, i__3, 
00053             i__4, i__5, i__6, i__7, i__8, i__9;
00054 
00055     /* Builtin functions */
00056     double sqrt(doublereal);
00057     integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00058 
00059     /* Local variables */
00060     integer i__, j, k, m, n, kl, jr, ku;
00061     real ulp, cond;
00062     integer jcol, kmax, mmax, nmax;
00063     real unfl, ovfl;
00064     extern /* Subroutine */ int cbdt01_(integer *, integer *, integer *, 
00065             complex *, integer *, complex *, integer *, real *, real *, 
00066             complex *, integer *, complex *, real *, real *), cbdt02_(integer 
00067             *, integer *, complex *, integer *, complex *, integer *, complex 
00068             *, integer *, complex *, real *, real *);
00069     logical badmm, badnn;
00070     integer imode, iinfo;
00071     extern /* Subroutine */ int cunt01_(char *, integer *, integer *, complex 
00072             *, integer *, complex *, integer *, real *, real *);
00073     real anorm;
00074     integer mnmin, mnmax, nmats, jsize, nerrs, itype, jtype, ntest;
00075     extern /* Subroutine */ int slahd2_(integer *, char *), cgbbrd_(
00076             char *, integer *, integer *, integer *, integer *, integer *, 
00077             complex *, integer *, real *, real *, complex *, integer *, 
00078             complex *, integer *, complex *, integer *, complex *, real *, 
00079             integer *);
00080     logical badnnb;
00081     extern doublereal slamch_(char *);
00082     extern /* Subroutine */ int clacpy_(char *, integer *, integer *, complex 
00083             *, integer *, complex *, integer *);
00084     integer idumma[1];
00085     extern /* Subroutine */ int claset_(char *, integer *, integer *, complex 
00086             *, complex *, complex *, integer *);
00087     integer ioldsd[4];
00088     extern /* Subroutine */ int xerbla_(char *, integer *), clatmr_(
00089             integer *, integer *, char *, integer *, char *, complex *, 
00090             integer *, real *, complex *, char *, char *, complex *, integer *
00091 , real *, complex *, integer *, real *, char *, integer *, 
00092             integer *, integer *, real *, real *, char *, complex *, integer *
00093 , integer *, integer *), clatms_(integer *, integer *, char *, integer *, char *, 
00094             real *, integer *, real *, real *, integer *, integer *, char *, 
00095             complex *, integer *, complex *, integer *);
00096     real amninv;
00097     integer jwidth;
00098     extern /* Subroutine */ int slasum_(char *, integer *, integer *, integer 
00099             *);
00100     real rtunfl, rtovfl, ulpinv;
00101     integer mtypes, ntestt;
00102 
00103     /* Fortran I/O blocks */
00104     static cilist io___41 = { 0, 0, 0, fmt_9999, 0 };
00105     static cilist io___43 = { 0, 0, 0, fmt_9999, 0 };
00106     static cilist io___45 = { 0, 0, 0, fmt_9998, 0 };
00107 
00108 
00109 
00110 /*  -- LAPACK test routine (new routine for release 2.0) -- */
00111 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00112 /*     November 2006 */
00113 
00114 /*     .. Scalar Arguments .. */
00115 /*     .. */
00116 /*     .. Array Arguments .. */
00117 /*     .. */
00118 
00119 /*  Purpose */
00120 /*  ======= */
00121 
00122 /*  CCHKBB tests the reduction of a general complex rectangular band */
00123 /*  matrix to real bidiagonal form. */
00124 
00125 /*  CGBBRD factors a general band matrix A as  Q B P* , where * means */
00126 /*  conjugate transpose, B is upper bidiagonal, and Q and P are unitary; */
00127 /*  CGBBRD can also overwrite a given matrix C with Q* C . */
00128 
00129 /*  For each pair of matrix dimensions (M,N) and each selected matrix */
00130 /*  type, an M by N matrix A and an M by NRHS matrix C are generated. */
00131 /*  The problem dimensions are as follows */
00132 /*     A:          M x N */
00133 /*     Q:          M x M */
00134 /*     P:          N x N */
00135 /*     B:          min(M,N) x min(M,N) */
00136 /*     C:          M x NRHS */
00137 
00138 /*  For each generated matrix, 4 tests are performed: */
00139 
00140 /*  (1)   | A - Q B PT | / ( |A| max(M,N) ulp ), PT = P' */
00141 
00142 /*  (2)   | I - Q' Q | / ( M ulp ) */
00143 
00144 /*  (3)   | I - PT PT' | / ( N ulp ) */
00145 
00146 /*  (4)   | Y - Q' C | / ( |Y| max(M,NRHS) ulp ), where Y = Q' C. */
00147 
00148 /*  The "types" are specified by a logical array DOTYPE( 1:NTYPES ); */
00149 /*  if DOTYPE(j) is .TRUE., then matrix type "j" will be generated. */
00150 /*  Currently, the list of possible types is: */
00151 
00152 /*  The possible matrix types are */
00153 
00154 /*  (1)  The zero matrix. */
00155 /*  (2)  The identity matrix. */
00156 
00157 /*  (3)  A diagonal matrix with evenly spaced entries */
00158 /*       1, ..., ULP  and random signs. */
00159 /*       (ULP = (first number larger than 1) - 1 ) */
00160 /*  (4)  A diagonal matrix with geometrically spaced entries */
00161 /*       1, ..., ULP  and random signs. */
00162 /*  (5)  A diagonal matrix with "clustered" entries 1, ULP, ..., ULP */
00163 /*       and random signs. */
00164 
00165 /*  (6)  Same as (3), but multiplied by SQRT( overflow threshold ) */
00166 /*  (7)  Same as (3), but multiplied by SQRT( underflow threshold ) */
00167 
00168 /*  (8)  A matrix of the form  U D V, where U and V are orthogonal and */
00169 /*       D has evenly spaced entries 1, ..., ULP with random signs */
00170 /*       on the diagonal. */
00171 
00172 /*  (9)  A matrix of the form  U D V, where U and V are orthogonal and */
00173 /*       D has geometrically spaced entries 1, ..., ULP with random */
00174 /*       signs on the diagonal. */
00175 
00176 /*  (10) A matrix of the form  U D V, where U and V are orthogonal and */
00177 /*       D has "clustered" entries 1, ULP,..., ULP with random */
00178 /*       signs on the diagonal. */
00179 
00180 /*  (11) Same as (8), but multiplied by SQRT( overflow threshold ) */
00181 /*  (12) Same as (8), but multiplied by SQRT( underflow threshold ) */
00182 
00183 /*  (13) Rectangular matrix with random entries chosen from (-1,1). */
00184 /*  (14) Same as (13), but multiplied by SQRT( overflow threshold ) */
00185 /*  (15) Same as (13), but multiplied by SQRT( underflow threshold ) */
00186 
00187 /*  Arguments */
00188 /*  ========= */
00189 
00190 /*  NSIZES  (input) INTEGER */
00191 /*          The number of values of M and N contained in the vectors */
00192 /*          MVAL and NVAL.  The matrix sizes are used in pairs (M,N). */
00193 /*          If NSIZES is zero, CCHKBB does nothing.  NSIZES must be at */
00194 /*          least zero. */
00195 
00196 /*  MVAL    (input) INTEGER array, dimension (NSIZES) */
00197 /*          The values of the matrix row dimension M. */
00198 
00199 /*  NVAL    (input) INTEGER array, dimension (NSIZES) */
00200 /*          The values of the matrix column dimension N. */
00201 
00202 /*  NWDTHS  (input) INTEGER */
00203 /*          The number of bandwidths to use.  If it is zero, */
00204 /*          CCHKBB does nothing.  It must be at least zero. */
00205 
00206 /*  KK      (input) INTEGER array, dimension (NWDTHS) */
00207 /*          An array containing the bandwidths to be used for the band */
00208 /*          matrices.  The values must be at least zero. */
00209 
00210 /*  NTYPES  (input) INTEGER */
00211 /*          The number of elements in DOTYPE.   If it is zero, CCHKBB */
00212 /*          does nothing.  It must be at least zero.  If it is MAXTYP+1 */
00213 /*          and NSIZES is 1, then an additional type, MAXTYP+1 is */
00214 /*          defined, which is to use whatever matrix is in A.  This */
00215 /*          is only useful if DOTYPE(1:MAXTYP) is .FALSE. and */
00216 /*          DOTYPE(MAXTYP+1) is .TRUE. . */
00217 
00218 /*  DOTYPE  (input) LOGICAL array, dimension (NTYPES) */
00219 /*          If DOTYPE(j) is .TRUE., then for each size in NN a */
00220 /*          matrix of that size and of type j will be generated. */
00221 /*          If NTYPES is smaller than the maximum number of types */
00222 /*          defined (PARAMETER MAXTYP), then types NTYPES+1 through */
00223 /*          MAXTYP will not be generated.  If NTYPES is larger */
00224 /*          than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES) */
00225 /*          will be ignored. */
00226 
00227 /*  NRHS    (input) INTEGER */
00228 /*          The number of columns in the "right-hand side" matrix C. */
00229 /*          If NRHS = 0, then the operations on the right-hand side will */
00230 /*          not be tested. NRHS must be at least 0. */
00231 
00232 /*  ISEED   (input/output) INTEGER array, dimension (4) */
00233 /*          On entry ISEED specifies the seed of the random number */
00234 /*          generator. The array elements should be between 0 and 4095; */
00235 /*          if not they will be reduced mod 4096.  Also, ISEED(4) must */
00236 /*          be odd.  The random number generator uses a linear */
00237 /*          congruential sequence limited to small integers, and so */
00238 /*          should produce machine independent random numbers. The */
00239 /*          values of ISEED are changed on exit, and can be used in the */
00240 /*          next call to CCHKBB to continue the same random number */
00241 /*          sequence. */
00242 
00243 /*  THRESH  (input) REAL */
00244 /*          A test will count as "failed" if the "error", computed as */
00245 /*          described above, exceeds THRESH.  Note that the error */
00246 /*          is scaled to be O(1), so THRESH should be a reasonably */
00247 /*          small multiple of 1, e.g., 10 or 100.  In particular, */
00248 /*          it should not depend on the precision (single vs. double) */
00249 /*          or the size of the matrix.  It must be at least zero. */
00250 
00251 /*  NOUNIT  (input) INTEGER */
00252 /*          The FORTRAN unit number for printing out error messages */
00253 /*          (e.g., if a routine returns IINFO not equal to 0.) */
00254 
00255 /*  A       (input/workspace) REAL array, dimension */
00256 /*                            (LDA, max(NN)) */
00257 /*          Used to hold the matrix A. */
00258 
00259 /*  LDA     (input) INTEGER */
00260 /*          The leading dimension of A.  It must be at least 1 */
00261 /*          and at least max( NN ). */
00262 
00263 /*  AB      (workspace) REAL array, dimension (LDAB, max(NN)) */
00264 /*          Used to hold A in band storage format. */
00265 
00266 /*  LDAB    (input) INTEGER */
00267 /*          The leading dimension of AB.  It must be at least 2 (not 1!) */
00268 /*          and at least max( KK )+1. */
00269 
00270 /*  BD      (workspace) REAL array, dimension (max(NN)) */
00271 /*          Used to hold the diagonal of the bidiagonal matrix computed */
00272 /*          by CGBBRD. */
00273 
00274 /*  BE      (workspace) REAL array, dimension (max(NN)) */
00275 /*          Used to hold the off-diagonal of the bidiagonal matrix */
00276 /*          computed by CGBBRD. */
00277 
00278 /*  Q       (workspace) COMPLEX array, dimension (LDQ, max(NN)) */
00279 /*          Used to hold the unitary matrix Q computed by CGBBRD. */
00280 
00281 /*  LDQ     (input) INTEGER */
00282 /*          The leading dimension of Q.  It must be at least 1 */
00283 /*          and at least max( NN ). */
00284 
00285 /*  P       (workspace) COMPLEX array, dimension (LDP, max(NN)) */
00286 /*          Used to hold the unitary matrix P computed by CGBBRD. */
00287 
00288 /*  LDP     (input) INTEGER */
00289 /*          The leading dimension of P.  It must be at least 1 */
00290 /*          and at least max( NN ). */
00291 
00292 /*  C       (workspace) COMPLEX array, dimension (LDC, max(NN)) */
00293 /*          Used to hold the matrix C updated by CGBBRD. */
00294 
00295 /*  LDC     (input) INTEGER */
00296 /*          The leading dimension of U.  It must be at least 1 */
00297 /*          and at least max( NN ). */
00298 
00299 /*  CC      (workspace) COMPLEX array, dimension (LDC, max(NN)) */
00300 /*          Used to hold a copy of the matrix C. */
00301 
00302 /*  WORK    (workspace) COMPLEX array, dimension (LWORK) */
00303 
00304 /*  LWORK   (input) INTEGER */
00305 /*          The number of entries in WORK.  This must be at least */
00306 /*          max( LDA+1, max(NN)+1 )*max(NN). */
00307 
00308 /*  RWORK   (workspace) REAL array, dimension (max(NN)) */
00309 
00310 /*  RESULT  (output) REAL array, dimension (4) */
00311 /*          The values computed by the tests described above. */
00312 /*          The values are currently limited to 1/ulp, to avoid */
00313 /*          overflow. */
00314 
00315 /*  INFO    (output) INTEGER */
00316 /*          If 0, then everything ran OK. */
00317 
00318 /* ----------------------------------------------------------------------- */
00319 
00320 /*       Some Local Variables and Parameters: */
00321 /*       ---- ----- --------- --- ---------- */
00322 /*       ZERO, ONE       Real 0 and 1. */
00323 /*       MAXTYP          The number of types defined. */
00324 /*       NTEST           The number of tests performed, or which can */
00325 /*                       be performed so far, for the current matrix. */
00326 /*       NTESTT          The total number of tests performed so far. */
00327 /*       NMAX            Largest value in NN. */
00328 /*       NMATS           The number of matrices generated so far. */
00329 /*       NERRS           The number of tests which have exceeded THRESH */
00330 /*                       so far. */
00331 /*       COND, IMODE     Values to be passed to the matrix generators. */
00332 /*       ANORM           Norm of A; passed to matrix generators. */
00333 
00334 /*       OVFL, UNFL      Overflow and underflow thresholds. */
00335 /*       ULP, ULPINV     Finest relative precision and its inverse. */
00336 /*       RTOVFL, RTUNFL  Square roots of the previous 2 values. */
00337 /*               The following four arrays decode JTYPE: */
00338 /*       KTYPE(j)        The general type (1-10) for type "j". */
00339 /*       KMODE(j)        The MODE value to be passed to the matrix */
00340 /*                       generator for type "j". */
00341 /*       KMAGN(j)        The order of magnitude ( O(1), */
00342 /*                       O(overflow^(1/2) ), O(underflow^(1/2) ) */
00343 
00344 /*  ===================================================================== */
00345 
00346 /*     .. Parameters .. */
00347 /*     .. */
00348 /*     .. Local Scalars .. */
00349 /*     .. */
00350 /*     .. Local Arrays .. */
00351 /*     .. */
00352 /*     .. External Functions .. */
00353 /*     .. */
00354 /*     .. External Subroutines .. */
00355 /*     .. */
00356 /*     .. Intrinsic Functions .. */
00357 /*     .. */
00358 /*     .. Data statements .. */
00359     /* Parameter adjustments */
00360     --mval;
00361     --nval;
00362     --kk;
00363     --dotype;
00364     --iseed;
00365     a_dim1 = *lda;
00366     a_offset = 1 + a_dim1;
00367     a -= a_offset;
00368     ab_dim1 = *ldab;
00369     ab_offset = 1 + ab_dim1;
00370     ab -= ab_offset;
00371     --bd;
00372     --be;
00373     q_dim1 = *ldq;
00374     q_offset = 1 + q_dim1;
00375     q -= q_offset;
00376     p_dim1 = *ldp;
00377     p_offset = 1 + p_dim1;
00378     p -= p_offset;
00379     cc_dim1 = *ldc;
00380     cc_offset = 1 + cc_dim1;
00381     cc -= cc_offset;
00382     c_dim1 = *ldc;
00383     c_offset = 1 + c_dim1;
00384     c__ -= c_offset;
00385     --work;
00386     --rwork;
00387     --result;
00388 
00389     /* Function Body */
00390 /*     .. */
00391 /*     .. Executable Statements .. */
00392 
00393 /*     Check for errors */
00394 
00395     ntestt = 0;
00396     *info = 0;
00397 
00398 /*     Important constants */
00399 
00400     badmm = FALSE_;
00401     badnn = FALSE_;
00402     mmax = 1;
00403     nmax = 1;
00404     mnmax = 1;
00405     i__1 = *nsizes;
00406     for (j = 1; j <= i__1; ++j) {
00407 /* Computing MAX */
00408         i__2 = mmax, i__3 = mval[j];
00409         mmax = max(i__2,i__3);
00410         if (mval[j] < 0) {
00411             badmm = TRUE_;
00412         }
00413 /* Computing MAX */
00414         i__2 = nmax, i__3 = nval[j];
00415         nmax = max(i__2,i__3);
00416         if (nval[j] < 0) {
00417             badnn = TRUE_;
00418         }
00419 /* Computing MAX */
00420 /* Computing MIN */
00421         i__4 = mval[j], i__5 = nval[j];
00422         i__2 = mnmax, i__3 = min(i__4,i__5);
00423         mnmax = max(i__2,i__3);
00424 /* L10: */
00425     }
00426 
00427     badnnb = FALSE_;
00428     kmax = 0;
00429     i__1 = *nwdths;
00430     for (j = 1; j <= i__1; ++j) {
00431 /* Computing MAX */
00432         i__2 = kmax, i__3 = kk[j];
00433         kmax = max(i__2,i__3);
00434         if (kk[j] < 0) {
00435             badnnb = TRUE_;
00436         }
00437 /* L20: */
00438     }
00439 
00440 /*     Check for errors */
00441 
00442     if (*nsizes < 0) {
00443         *info = -1;
00444     } else if (badmm) {
00445         *info = -2;
00446     } else if (badnn) {
00447         *info = -3;
00448     } else if (*nwdths < 0) {
00449         *info = -4;
00450     } else if (badnnb) {
00451         *info = -5;
00452     } else if (*ntypes < 0) {
00453         *info = -6;
00454     } else if (*nrhs < 0) {
00455         *info = -8;
00456     } else if (*lda < nmax) {
00457         *info = -13;
00458     } else if (*ldab < (kmax << 1) + 1) {
00459         *info = -15;
00460     } else if (*ldq < nmax) {
00461         *info = -19;
00462     } else if (*ldp < nmax) {
00463         *info = -21;
00464     } else if (*ldc < nmax) {
00465         *info = -23;
00466     } else if ((max(*lda,nmax) + 1) * nmax > *lwork) {
00467         *info = -26;
00468     }
00469 
00470     if (*info != 0) {
00471         i__1 = -(*info);
00472         xerbla_("CCHKBB", &i__1);
00473         return 0;
00474     }
00475 
00476 /*     Quick return if possible */
00477 
00478     if (*nsizes == 0 || *ntypes == 0 || *nwdths == 0) {
00479         return 0;
00480     }
00481 
00482 /*     More Important constants */
00483 
00484     unfl = slamch_("Safe minimum");
00485     ovfl = 1.f / unfl;
00486     ulp = slamch_("Epsilon") * slamch_("Base");
00487     ulpinv = 1.f / ulp;
00488     rtunfl = sqrt(unfl);
00489     rtovfl = sqrt(ovfl);
00490 
00491 /*     Loop over sizes, widths, types */
00492 
00493     nerrs = 0;
00494     nmats = 0;
00495 
00496     i__1 = *nsizes;
00497     for (jsize = 1; jsize <= i__1; ++jsize) {
00498         m = mval[jsize];
00499         n = nval[jsize];
00500         mnmin = min(m,n);
00501 /* Computing MAX */
00502         i__2 = max(1,m);
00503         amninv = 1.f / (real) max(i__2,n);
00504 
00505         i__2 = *nwdths;
00506         for (jwidth = 1; jwidth <= i__2; ++jwidth) {
00507             k = kk[jwidth];
00508             if (k >= m && k >= n) {
00509                 goto L150;
00510             }
00511 /* Computing MAX */
00512 /* Computing MIN */
00513             i__5 = m - 1;
00514             i__3 = 0, i__4 = min(i__5,k);
00515             kl = max(i__3,i__4);
00516 /* Computing MAX */
00517 /* Computing MIN */
00518             i__5 = n - 1;
00519             i__3 = 0, i__4 = min(i__5,k);
00520             ku = max(i__3,i__4);
00521 
00522             if (*nsizes != 1) {
00523                 mtypes = min(15,*ntypes);
00524             } else {
00525                 mtypes = min(16,*ntypes);
00526             }
00527 
00528             i__3 = mtypes;
00529             for (jtype = 1; jtype <= i__3; ++jtype) {
00530                 if (! dotype[jtype]) {
00531                     goto L140;
00532                 }
00533                 ++nmats;
00534                 ntest = 0;
00535 
00536                 for (j = 1; j <= 4; ++j) {
00537                     ioldsd[j - 1] = iseed[j];
00538 /* L30: */
00539                 }
00540 
00541 /*              Compute "A". */
00542 
00543 /*              Control parameters: */
00544 
00545 /*                  KMAGN  KMODE        KTYPE */
00546 /*              =1  O(1)   clustered 1  zero */
00547 /*              =2  large  clustered 2  identity */
00548 /*              =3  small  exponential  (none) */
00549 /*              =4         arithmetic   diagonal, (w/ singular values) */
00550 /*              =5         random log   (none) */
00551 /*              =6         random       nonhermitian, w/ singular values */
00552 /*              =7                      (none) */
00553 /*              =8                      (none) */
00554 /*              =9                      random nonhermitian */
00555 
00556                 if (mtypes > 15) {
00557                     goto L90;
00558                 }
00559 
00560                 itype = ktype[jtype - 1];
00561                 imode = kmode[jtype - 1];
00562 
00563 /*              Compute norm */
00564 
00565                 switch (kmagn[jtype - 1]) {
00566                     case 1:  goto L40;
00567                     case 2:  goto L50;
00568                     case 3:  goto L60;
00569                 }
00570 
00571 L40:
00572                 anorm = 1.f;
00573                 goto L70;
00574 
00575 L50:
00576                 anorm = rtovfl * ulp * amninv;
00577                 goto L70;
00578 
00579 L60:
00580                 anorm = rtunfl * max(m,n) * ulpinv;
00581                 goto L70;
00582 
00583 L70:
00584 
00585                 claset_("Full", lda, &n, &c_b1, &c_b1, &a[a_offset], lda);
00586                 claset_("Full", ldab, &n, &c_b1, &c_b1, &ab[ab_offset], ldab);
00587                 iinfo = 0;
00588                 cond = ulpinv;
00589 
00590 /*              Special Matrices -- Identity & Jordan block */
00591 
00592 /*                 Zero */
00593 
00594                 if (itype == 1) {
00595                     iinfo = 0;
00596 
00597                 } else if (itype == 2) {
00598 
00599 /*                 Identity */
00600 
00601                     i__4 = n;
00602                     for (jcol = 1; jcol <= i__4; ++jcol) {
00603                         i__5 = jcol + jcol * a_dim1;
00604                         a[i__5].r = anorm, a[i__5].i = 0.f;
00605 /* L80: */
00606                     }
00607 
00608                 } else if (itype == 4) {
00609 
00610 /*                 Diagonal Matrix, singular values specified */
00611 
00612                     clatms_(&m, &n, "S", &iseed[1], "N", &rwork[1], &imode, &
00613                             cond, &anorm, &c__0, &c__0, "N", &a[a_offset], 
00614                             lda, &work[1], &iinfo);
00615 
00616                 } else if (itype == 6) {
00617 
00618 /*                 Nonhermitian, singular values specified */
00619 
00620                     clatms_(&m, &n, "S", &iseed[1], "N", &rwork[1], &imode, &
00621                             cond, &anorm, &kl, &ku, "N", &a[a_offset], lda, &
00622                             work[1], &iinfo);
00623 
00624                 } else if (itype == 9) {
00625 
00626 /*                 Nonhermitian, random entries */
00627 
00628                     clatmr_(&m, &n, "S", &iseed[1], "N", &work[1], &c__6, &
00629                             c_b33, &c_b2, "T", "N", &work[n + 1], &c__1, &
00630                             c_b33, &work[(n << 1) + 1], &c__1, &c_b33, "N", 
00631                             idumma, &kl, &ku, &c_b41, &anorm, "N", &a[
00632                             a_offset], lda, idumma, &iinfo);
00633 
00634                 } else {
00635 
00636                     iinfo = 1;
00637                 }
00638 
00639 /*              Generate Right-Hand Side */
00640 
00641                 clatmr_(&m, nrhs, "S", &iseed[1], "N", &work[1], &c__6, &
00642                         c_b33, &c_b2, "T", "N", &work[m + 1], &c__1, &c_b33, &
00643                         work[(m << 1) + 1], &c__1, &c_b33, "N", idumma, &m, 
00644                         nrhs, &c_b41, &c_b33, "NO", &c__[c_offset], ldc, 
00645                         idumma, &iinfo);
00646 
00647                 if (iinfo != 0) {
00648                     io___41.ciunit = *nounit;
00649                     s_wsfe(&io___41);
00650                     do_fio(&c__1, "Generator", (ftnlen)9);
00651                     do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00652                     do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00653                     do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00654                     do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
00655                             ;
00656                     e_wsfe();
00657                     *info = abs(iinfo);
00658                     return 0;
00659                 }
00660 
00661 L90:
00662 
00663 /*              Copy A to band storage. */
00664 
00665                 i__4 = n;
00666                 for (j = 1; j <= i__4; ++j) {
00667 /* Computing MAX */
00668                     i__5 = 1, i__6 = j - ku;
00669 /* Computing MIN */
00670                     i__8 = m, i__9 = j + kl;
00671                     i__7 = min(i__8,i__9);
00672                     for (i__ = max(i__5,i__6); i__ <= i__7; ++i__) {
00673                         i__5 = ku + 1 + i__ - j + j * ab_dim1;
00674                         i__6 = i__ + j * a_dim1;
00675                         ab[i__5].r = a[i__6].r, ab[i__5].i = a[i__6].i;
00676 /* L100: */
00677                     }
00678 /* L110: */
00679                 }
00680 
00681 /*              Copy C */
00682 
00683                 clacpy_("Full", &m, nrhs, &c__[c_offset], ldc, &cc[cc_offset], 
00684                          ldc);
00685 
00686 /*              Call CGBBRD to compute B, Q and P, and to update C. */
00687 
00688                 cgbbrd_("B", &m, &n, nrhs, &kl, &ku, &ab[ab_offset], ldab, &
00689                         bd[1], &be[1], &q[q_offset], ldq, &p[p_offset], ldp, &
00690                         cc[cc_offset], ldc, &work[1], &rwork[1], &iinfo);
00691 
00692                 if (iinfo != 0) {
00693                     io___43.ciunit = *nounit;
00694                     s_wsfe(&io___43);
00695                     do_fio(&c__1, "CGBBRD", (ftnlen)6);
00696                     do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00697                     do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00698                     do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00699                     do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
00700                             ;
00701                     e_wsfe();
00702                     *info = abs(iinfo);
00703                     if (iinfo < 0) {
00704                         return 0;
00705                     } else {
00706                         result[1] = ulpinv;
00707                         goto L120;
00708                     }
00709                 }
00710 
00711 /*              Test 1:  Check the decomposition A := Q * B * P' */
00712 /*                   2:  Check the orthogonality of Q */
00713 /*                   3:  Check the orthogonality of P */
00714 /*                   4:  Check the computation of Q' * C */
00715 
00716                 cbdt01_(&m, &n, &c_n1, &a[a_offset], lda, &q[q_offset], ldq, &
00717                         bd[1], &be[1], &p[p_offset], ldp, &work[1], &rwork[1], 
00718                          &result[1]);
00719                 cunt01_("Columns", &m, &m, &q[q_offset], ldq, &work[1], lwork, 
00720                          &rwork[1], &result[2]);
00721                 cunt01_("Rows", &n, &n, &p[p_offset], ldp, &work[1], lwork, &
00722                         rwork[1], &result[3]);
00723                 cbdt02_(&m, nrhs, &c__[c_offset], ldc, &cc[cc_offset], ldc, &
00724                         q[q_offset], ldq, &work[1], &rwork[1], &result[4]);
00725 
00726 /*              End of Loop -- Check for RESULT(j) > THRESH */
00727 
00728                 ntest = 4;
00729 L120:
00730                 ntestt += ntest;
00731 
00732 /*              Print out tests which fail. */
00733 
00734                 i__4 = ntest;
00735                 for (jr = 1; jr <= i__4; ++jr) {
00736                     if (result[jr] >= *thresh) {
00737                         if (nerrs == 0) {
00738                             slahd2_(nounit, "CBB");
00739                         }
00740                         ++nerrs;
00741                         io___45.ciunit = *nounit;
00742                         s_wsfe(&io___45);
00743                         do_fio(&c__1, (char *)&m, (ftnlen)sizeof(integer));
00744                         do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00745                         do_fio(&c__1, (char *)&k, (ftnlen)sizeof(integer));
00746                         do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(
00747                                 integer));
00748                         do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer))
00749                                 ;
00750                         do_fio(&c__1, (char *)&jr, (ftnlen)sizeof(integer));
00751                         do_fio(&c__1, (char *)&result[jr], (ftnlen)sizeof(
00752                                 real));
00753                         e_wsfe();
00754                     }
00755 /* L130: */
00756                 }
00757 
00758 L140:
00759                 ;
00760             }
00761 L150:
00762             ;
00763         }
00764 /* L160: */
00765     }
00766 
00767 /*     Summary */
00768 
00769     slasum_("CBB", nounit, &nerrs, &ntestt);
00770     return 0;
00771 
00772 
00773 /*     End of CCHKBB */
00774 
00775 } /* cchkbb_ */


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