00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015
00016
00017
00018 static real c_b18 = 0.f;
00019 static integer c__0 = 0;
00020 static real c_b32 = 1.f;
00021 static integer c__4 = 4;
00022 static integer c__6 = 6;
00023 static integer c__1 = 1;
00024
00025 int schkhs_(integer *nsizes, integer *nn, integer *ntypes,
00026 logical *dotype, integer *iseed, real *thresh, integer *nounit, real *
00027 a, integer *lda, real *h__, real *t1, real *t2, real *u, integer *ldu,
00028 real *z__, real *uz, real *wr1, real *wi1, real *wr3, real *wi3,
00029 real *evectl, real *evectr, real *evecty, real *evectx, real *uu,
00030 real *tau, real *work, integer *nwork, integer *iwork, logical *
00031 select, real *result, integer *info)
00032 {
00033
00034
00035 static integer ktype[21] = { 1,2,3,4,4,4,4,4,6,6,6,6,6,6,6,6,6,6,9,9,9 };
00036 static integer kmagn[21] = { 1,1,1,1,1,1,2,3,1,1,1,1,1,1,1,1,2,3,1,2,3 };
00037 static integer kmode[21] = { 0,0,0,4,3,1,4,4,4,3,1,5,4,3,1,5,5,5,4,3,1 };
00038 static integer kconds[21] = { 0,0,0,0,0,0,0,0,1,1,1,1,2,2,2,2,2,2,0,0,0 };
00039
00040
00041 static char fmt_9999[] = "(\002 SCHKHS: \002,a,\002 returned INFO=\002,i"
00042 "6,\002.\002,/9x,\002N=\002,i6,\002, JTYPE=\002,i6,\002, ISEED="
00043 "(\002,3(i5,\002,\002),i5,\002)\002)";
00044 static char fmt_9998[] = "(\002 SCHKHS: \002,a,\002 Eigenvectors from"
00045 " \002,a,\002 incorrectly \002,\002normalized.\002,/\002 Bits of "
00046 "error=\002,0p,g10.3,\002,\002,9x,\002N=\002,i6,\002, JTYPE=\002,"
00047 "i6,\002, ISEED=(\002,3(i5,\002,\002),i5,\002)\002)";
00048 static char fmt_9997[] = "(\002 SCHKHS: Selected \002,a,\002 Eigenvector"
00049 "s from \002,a,\002 do not match other eigenvectors \002,9x,\002N="
00050 "\002,i6,\002, JTYPE=\002,i6,\002, ISEED=(\002,3(i5,\002,\002),i5,"
00051 "\002)\002)";
00052
00053
00054 integer a_dim1, a_offset, evectl_dim1, evectl_offset, evectr_dim1,
00055 evectr_offset, evectx_dim1, evectx_offset, evecty_dim1,
00056 evecty_offset, h_dim1, h_offset, t1_dim1, t1_offset, t2_dim1,
00057 t2_offset, u_dim1, u_offset, uu_dim1, uu_offset, uz_dim1,
00058 uz_offset, z_dim1, z_offset, i__1, i__2, i__3, i__4;
00059 real r__1, r__2, r__3, r__4, r__5, r__6;
00060
00061
00062 double sqrt(doublereal);
00063 integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00064
00065
00066 integer i__, j, k, n, n1, jj, in, ihi, ilo;
00067 real ulp, cond;
00068 integer jcol, nmax;
00069 real unfl, ovfl, temp1, temp2;
00070 logical badnn, match;
00071 integer imode;
00072 real dumma[6];
00073 integer iinfo, nselc;
00074 real conds;
00075 extern int sget10_(integer *, integer *, real *, integer
00076 *, real *, integer *, real *, real *), sgemm_(char *, char *,
00077 integer *, integer *, integer *, real *, real *, integer *, real *
00078 , integer *, real *, real *, integer *), sget22_(
00079 char *, char *, char *, integer *, real *, integer *, real *,
00080 integer *, real *, real *, real *, real *)
00081 ;
00082 real aninv, anorm;
00083 integer nmats, nselr, jsize;
00084 extern int shst01_(integer *, integer *, integer *, real
00085 *, integer *, real *, integer *, real *, integer *, real *,
00086 integer *, real *);
00087 integer nerrs, itype, jtype, ntest;
00088 extern int scopy_(integer *, real *, integer *, real *,
00089 integer *);
00090 real rtulp;
00091 extern int slabad_(real *, real *);
00092 char adumma[1*1];
00093 extern doublereal slamch_(char *);
00094 integer idumma[1];
00095 extern int sgehrd_(integer *, integer *, integer *, real
00096 *, integer *, real *, real *, integer *, integer *);
00097 integer ioldsd[4];
00098 extern int xerbla_(char *, integer *), slatme_(
00099 integer *, char *, integer *, real *, integer *, real *, real *,
00100 char *, char *, char *, char *, real *, integer *, real *,
00101 integer *, integer *, real *, real *, integer *, real *, integer *
00102 ), shsein_(char *, char *,
00103 char *, logical *, integer *, real *, integer *, real *, real *,
00104 real *, integer *, real *, integer *, integer *, integer *, real *
00105 , integer *, integer *, integer *),
00106 slacpy_(char *, integer *, integer *, real *, integer *, real *,
00107 integer *), slafts_(char *, integer *, integer *, integer
00108 *, integer *, real *, integer *, real *, integer *, integer *), slaset_(char *, integer *, integer *, real *, real *,
00109 real *, integer *), slatmr_(integer *, integer *, char *,
00110 integer *, char *, real *, integer *, real *, real *, char *,
00111 char *, real *, integer *, real *, real *, integer *, real *,
00112 char *, integer *, integer *, integer *, real *, real *, char *,
00113 real *, integer *, integer *, integer *), slatms_(integer *, integer *, char *,
00114 integer *, char *, real *, integer *, real *, real *, integer *,
00115 integer *, char *, real *, integer *, real *, integer *), shseqr_(char *, char *, integer *, integer *,
00116 integer *, real *, integer *, real *, real *, real *, integer *,
00117 real *, integer *, integer *), slasum_(char *,
00118 integer *, integer *, integer *), sorghr_(integer *,
00119 integer *, integer *, real *, integer *, real *, real *, integer *
00120 , integer *), strevc_(char *, char *, logical *, integer *, real *
00121 , integer *, real *, integer *, real *, integer *, integer *,
00122 integer *, real *, integer *);
00123 real rtunfl, rtovfl, rtulpi, ulpinv;
00124 integer mtypes, ntestt;
00125 extern int sormhr_(char *, char *, integer *, integer *,
00126 integer *, integer *, real *, integer *, real *, real *, integer *
00127 , real *, integer *, integer *);
00128
00129
00130 static cilist io___36 = { 0, 0, 0, fmt_9999, 0 };
00131 static cilist io___39 = { 0, 0, 0, fmt_9999, 0 };
00132 static cilist io___41 = { 0, 0, 0, fmt_9999, 0 };
00133 static cilist io___42 = { 0, 0, 0, fmt_9999, 0 };
00134 static cilist io___43 = { 0, 0, 0, fmt_9999, 0 };
00135 static cilist io___50 = { 0, 0, 0, fmt_9999, 0 };
00136 static cilist io___51 = { 0, 0, 0, fmt_9998, 0 };
00137 static cilist io___52 = { 0, 0, 0, fmt_9999, 0 };
00138 static cilist io___56 = { 0, 0, 0, fmt_9997, 0 };
00139 static cilist io___57 = { 0, 0, 0, fmt_9999, 0 };
00140 static cilist io___58 = { 0, 0, 0, fmt_9998, 0 };
00141 static cilist io___59 = { 0, 0, 0, fmt_9999, 0 };
00142 static cilist io___60 = { 0, 0, 0, fmt_9997, 0 };
00143 static cilist io___61 = { 0, 0, 0, fmt_9999, 0 };
00144 static cilist io___62 = { 0, 0, 0, fmt_9998, 0 };
00145 static cilist io___63 = { 0, 0, 0, fmt_9999, 0 };
00146 static cilist io___64 = { 0, 0, 0, fmt_9998, 0 };
00147 static cilist io___65 = { 0, 0, 0, fmt_9999, 0 };
00148 static cilist io___66 = { 0, 0, 0, fmt_9999, 0 };
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194
00195
00196
00197
00198
00199
00200
00201
00202
00203
00204
00205
00206
00207
00208
00209
00210
00211
00212
00213
00214
00215
00216
00217
00218
00219
00220
00221
00222
00223
00224
00225
00226
00227
00228
00229
00230
00231
00232
00233
00234
00235
00236
00237
00238
00239
00240
00241
00242
00243
00244
00245
00246
00247
00248
00249
00250
00251
00252
00253
00254
00255
00256
00257
00258
00259
00260
00261
00262
00263
00264
00265
00266
00267
00268
00269
00270
00271
00272
00273
00274
00275
00276
00277
00278
00279
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294
00295
00296
00297
00298
00299
00300
00301
00302
00303
00304
00305
00306
00307
00308
00309
00310
00311
00312
00313
00314
00315
00316
00317
00318
00319
00320
00321
00322
00323
00324
00325
00326
00327
00328
00329
00330
00331
00332
00333
00334
00335
00336
00337
00338
00339
00340
00341
00342
00343
00344
00345
00346
00347
00348
00349
00350
00351
00352
00353
00354
00355
00356
00357
00358
00359
00360
00361
00362
00363
00364
00365
00366
00367
00368
00369
00370
00371
00372
00373
00374
00375
00376
00377
00378
00379
00380
00381
00382
00383
00384
00385
00386
00387
00388
00389
00390
00391
00392
00393
00394
00395
00396
00397
00398
00399
00400
00401
00402
00403
00404
00405
00406
00407
00408
00409
00410
00411
00412
00413
00414
00415
00416
00417
00418
00419
00420
00421
00422
00423
00424
00425
00426
00427
00428
00429
00430
00431
00432
00433
00434
00435
00436
00437
00438
00439
00440
00441
00442
00443
00444
00445
00446
00447
00448
00449
00450
00451
00452
00453
00454
00455
00456
00457
00458
00459
00460
00461
00462
00463
00464
00465
00466
00467
00468
00469
00470
00471
00472
00473
00474
00475
00476
00477
00478
00479
00480
00481
00482
00483
00484
00485
00486
00487
00488
00489
00490
00491
00492
00493
00494
00495
00496
00497
00498
00499
00500
00501
00502
00503
00504
00505
00506
00507
00508
00509
00510
00511
00512
00513
00514
00515
00516
00517
00518
00519
00520
00521
00522
00523
00524
00525
00526
00527 --nn;
00528 --dotype;
00529 --iseed;
00530 t2_dim1 = *lda;
00531 t2_offset = 1 + t2_dim1;
00532 t2 -= t2_offset;
00533 t1_dim1 = *lda;
00534 t1_offset = 1 + t1_dim1;
00535 t1 -= t1_offset;
00536 h_dim1 = *lda;
00537 h_offset = 1 + h_dim1;
00538 h__ -= h_offset;
00539 a_dim1 = *lda;
00540 a_offset = 1 + a_dim1;
00541 a -= a_offset;
00542 uu_dim1 = *ldu;
00543 uu_offset = 1 + uu_dim1;
00544 uu -= uu_offset;
00545 evectx_dim1 = *ldu;
00546 evectx_offset = 1 + evectx_dim1;
00547 evectx -= evectx_offset;
00548 evecty_dim1 = *ldu;
00549 evecty_offset = 1 + evecty_dim1;
00550 evecty -= evecty_offset;
00551 evectr_dim1 = *ldu;
00552 evectr_offset = 1 + evectr_dim1;
00553 evectr -= evectr_offset;
00554 evectl_dim1 = *ldu;
00555 evectl_offset = 1 + evectl_dim1;
00556 evectl -= evectl_offset;
00557 uz_dim1 = *ldu;
00558 uz_offset = 1 + uz_dim1;
00559 uz -= uz_offset;
00560 z_dim1 = *ldu;
00561 z_offset = 1 + z_dim1;
00562 z__ -= z_offset;
00563 u_dim1 = *ldu;
00564 u_offset = 1 + u_dim1;
00565 u -= u_offset;
00566 --wr1;
00567 --wi1;
00568 --wr3;
00569 --wi3;
00570 --tau;
00571 --work;
00572 --iwork;
00573 --select;
00574 --result;
00575
00576
00577
00578
00579
00580
00581
00582 ntestt = 0;
00583 *info = 0;
00584
00585 badnn = FALSE_;
00586 nmax = 0;
00587 i__1 = *nsizes;
00588 for (j = 1; j <= i__1; ++j) {
00589
00590 i__2 = nmax, i__3 = nn[j];
00591 nmax = max(i__2,i__3);
00592 if (nn[j] < 0) {
00593 badnn = TRUE_;
00594 }
00595
00596 }
00597
00598
00599
00600 if (*nsizes < 0) {
00601 *info = -1;
00602 } else if (badnn) {
00603 *info = -2;
00604 } else if (*ntypes < 0) {
00605 *info = -3;
00606 } else if (*thresh < 0.f) {
00607 *info = -6;
00608 } else if (*lda <= 1 || *lda < nmax) {
00609 *info = -9;
00610 } else if (*ldu <= 1 || *ldu < nmax) {
00611 *info = -14;
00612 } else if ((nmax << 2) * nmax + 2 > *nwork) {
00613 *info = -28;
00614 }
00615
00616 if (*info != 0) {
00617 i__1 = -(*info);
00618 xerbla_("SCHKHS", &i__1);
00619 return 0;
00620 }
00621
00622
00623
00624 if (*nsizes == 0 || *ntypes == 0) {
00625 return 0;
00626 }
00627
00628
00629
00630 unfl = slamch_("Safe minimum");
00631 ovfl = slamch_("Overflow");
00632 slabad_(&unfl, &ovfl);
00633 ulp = slamch_("Epsilon") * slamch_("Base");
00634 ulpinv = 1.f / ulp;
00635 rtunfl = sqrt(unfl);
00636 rtovfl = sqrt(ovfl);
00637 rtulp = sqrt(ulp);
00638 rtulpi = 1.f / rtulp;
00639
00640
00641
00642 nerrs = 0;
00643 nmats = 0;
00644
00645 i__1 = *nsizes;
00646 for (jsize = 1; jsize <= i__1; ++jsize) {
00647 n = nn[jsize];
00648 if (n == 0) {
00649 goto L270;
00650 }
00651 n1 = max(1,n);
00652 aninv = 1.f / (real) n1;
00653
00654 if (*nsizes != 1) {
00655 mtypes = min(21,*ntypes);
00656 } else {
00657 mtypes = min(22,*ntypes);
00658 }
00659
00660 i__2 = mtypes;
00661 for (jtype = 1; jtype <= i__2; ++jtype) {
00662 if (! dotype[jtype]) {
00663 goto L260;
00664 }
00665 ++nmats;
00666 ntest = 0;
00667
00668
00669
00670 for (j = 1; j <= 4; ++j) {
00671 ioldsd[j - 1] = iseed[j];
00672
00673 }
00674
00675
00676
00677 for (j = 1; j <= 14; ++j) {
00678 result[j] = 0.f;
00679
00680 }
00681
00682
00683
00684
00685
00686
00687
00688
00689
00690
00691
00692
00693
00694
00695
00696
00697
00698 if (mtypes > 21) {
00699 goto L100;
00700 }
00701
00702 itype = ktype[jtype - 1];
00703 imode = kmode[jtype - 1];
00704
00705
00706
00707 switch (kmagn[jtype - 1]) {
00708 case 1: goto L40;
00709 case 2: goto L50;
00710 case 3: goto L60;
00711 }
00712
00713 L40:
00714 anorm = 1.f;
00715 goto L70;
00716
00717 L50:
00718 anorm = rtovfl * ulp * aninv;
00719 goto L70;
00720
00721 L60:
00722 anorm = rtunfl * n * ulpinv;
00723 goto L70;
00724
00725 L70:
00726
00727 slaset_("Full", lda, &n, &c_b18, &c_b18, &a[a_offset], lda);
00728 iinfo = 0;
00729 cond = ulpinv;
00730
00731
00732
00733 if (itype == 1) {
00734
00735
00736
00737 iinfo = 0;
00738
00739 } else if (itype == 2) {
00740
00741
00742
00743 i__3 = n;
00744 for (jcol = 1; jcol <= i__3; ++jcol) {
00745 a[jcol + jcol * a_dim1] = anorm;
00746
00747 }
00748
00749 } else if (itype == 3) {
00750
00751
00752
00753 i__3 = n;
00754 for (jcol = 1; jcol <= i__3; ++jcol) {
00755 a[jcol + jcol * a_dim1] = anorm;
00756 if (jcol > 1) {
00757 a[jcol + (jcol - 1) * a_dim1] = 1.f;
00758 }
00759
00760 }
00761
00762 } else if (itype == 4) {
00763
00764
00765
00766 slatms_(&n, &n, "S", &iseed[1], "S", &work[1], &imode, &cond,
00767 &anorm, &c__0, &c__0, "N", &a[a_offset], lda, &work[n
00768 + 1], &iinfo);
00769
00770 } else if (itype == 5) {
00771
00772
00773
00774 slatms_(&n, &n, "S", &iseed[1], "S", &work[1], &imode, &cond,
00775 &anorm, &n, &n, "N", &a[a_offset], lda, &work[n + 1],
00776 &iinfo);
00777
00778 } else if (itype == 6) {
00779
00780
00781
00782 if (kconds[jtype - 1] == 1) {
00783 conds = 1.f;
00784 } else if (kconds[jtype - 1] == 2) {
00785 conds = rtulpi;
00786 } else {
00787 conds = 0.f;
00788 }
00789
00790 *(unsigned char *)&adumma[0] = ' ';
00791 slatme_(&n, "S", &iseed[1], &work[1], &imode, &cond, &c_b32,
00792 adumma, "T", "T", "T", &work[n + 1], &c__4, &conds, &
00793 n, &n, &anorm, &a[a_offset], lda, &work[(n << 1) + 1],
00794 &iinfo);
00795
00796 } else if (itype == 7) {
00797
00798
00799
00800 slatmr_(&n, &n, "S", &iseed[1], "S", &work[1], &c__6, &c_b32,
00801 &c_b32, "T", "N", &work[n + 1], &c__1, &c_b32, &work[(
00802 n << 1) + 1], &c__1, &c_b32, "N", idumma, &c__0, &
00803 c__0, &c_b18, &anorm, "NO", &a[a_offset], lda, &iwork[
00804 1], &iinfo);
00805
00806 } else if (itype == 8) {
00807
00808
00809
00810 slatmr_(&n, &n, "S", &iseed[1], "S", &work[1], &c__6, &c_b32,
00811 &c_b32, "T", "N", &work[n + 1], &c__1, &c_b32, &work[(
00812 n << 1) + 1], &c__1, &c_b32, "N", idumma, &n, &n, &
00813 c_b18, &anorm, "NO", &a[a_offset], lda, &iwork[1], &
00814 iinfo);
00815
00816 } else if (itype == 9) {
00817
00818
00819
00820 slatmr_(&n, &n, "S", &iseed[1], "N", &work[1], &c__6, &c_b32,
00821 &c_b32, "T", "N", &work[n + 1], &c__1, &c_b32, &work[(
00822 n << 1) + 1], &c__1, &c_b32, "N", idumma, &n, &n, &
00823 c_b18, &anorm, "NO", &a[a_offset], lda, &iwork[1], &
00824 iinfo);
00825
00826 } else if (itype == 10) {
00827
00828
00829
00830 slatmr_(&n, &n, "S", &iseed[1], "N", &work[1], &c__6, &c_b32,
00831 &c_b32, "T", "N", &work[n + 1], &c__1, &c_b32, &work[(
00832 n << 1) + 1], &c__1, &c_b32, "N", idumma, &n, &c__0, &
00833 c_b18, &anorm, "NO", &a[a_offset], lda, &iwork[1], &
00834 iinfo);
00835
00836 } else {
00837
00838 iinfo = 1;
00839 }
00840
00841 if (iinfo != 0) {
00842 io___36.ciunit = *nounit;
00843 s_wsfe(&io___36);
00844 do_fio(&c__1, "Generator", (ftnlen)9);
00845 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00846 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00847 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00848 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00849 e_wsfe();
00850 *info = abs(iinfo);
00851 return 0;
00852 }
00853
00854 L100:
00855
00856
00857
00858 slacpy_(" ", &n, &n, &a[a_offset], lda, &h__[h_offset], lda);
00859
00860 ntest = 1;
00861
00862 ilo = 1;
00863 ihi = n;
00864
00865 i__3 = *nwork - n;
00866 sgehrd_(&n, &ilo, &ihi, &h__[h_offset], lda, &work[1], &work[n +
00867 1], &i__3, &iinfo);
00868
00869 if (iinfo != 0) {
00870 result[1] = ulpinv;
00871 io___39.ciunit = *nounit;
00872 s_wsfe(&io___39);
00873 do_fio(&c__1, "SGEHRD", (ftnlen)6);
00874 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00875 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00876 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00877 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00878 e_wsfe();
00879 *info = abs(iinfo);
00880 goto L250;
00881 }
00882
00883 i__3 = n - 1;
00884 for (j = 1; j <= i__3; ++j) {
00885 uu[j + 1 + j * uu_dim1] = 0.f;
00886 i__4 = n;
00887 for (i__ = j + 2; i__ <= i__4; ++i__) {
00888 u[i__ + j * u_dim1] = h__[i__ + j * h_dim1];
00889 uu[i__ + j * uu_dim1] = h__[i__ + j * h_dim1];
00890 h__[i__ + j * h_dim1] = 0.f;
00891
00892 }
00893
00894 }
00895 i__3 = n - 1;
00896 scopy_(&i__3, &work[1], &c__1, &tau[1], &c__1);
00897 i__3 = *nwork - n;
00898 sorghr_(&n, &ilo, &ihi, &u[u_offset], ldu, &work[1], &work[n + 1],
00899 &i__3, &iinfo);
00900 ntest = 2;
00901
00902 shst01_(&n, &ilo, &ihi, &a[a_offset], lda, &h__[h_offset], lda, &
00903 u[u_offset], ldu, &work[1], nwork, &result[1]);
00904
00905
00906
00907
00908
00909 slacpy_(" ", &n, &n, &h__[h_offset], lda, &t2[t2_offset], lda);
00910 ntest = 3;
00911 result[3] = ulpinv;
00912
00913 shseqr_("E", "N", &n, &ilo, &ihi, &t2[t2_offset], lda, &wr3[1], &
00914 wi3[1], &uz[uz_offset], ldu, &work[1], nwork, &iinfo);
00915 if (iinfo != 0) {
00916 io___41.ciunit = *nounit;
00917 s_wsfe(&io___41);
00918 do_fio(&c__1, "SHSEQR(E)", (ftnlen)9);
00919 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00920 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00921 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00922 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00923 e_wsfe();
00924 if (iinfo <= n + 2) {
00925 *info = abs(iinfo);
00926 goto L250;
00927 }
00928 }
00929
00930
00931
00932 slacpy_(" ", &n, &n, &h__[h_offset], lda, &t2[t2_offset], lda);
00933
00934 shseqr_("S", "N", &n, &ilo, &ihi, &t2[t2_offset], lda, &wr1[1], &
00935 wi1[1], &uz[uz_offset], ldu, &work[1], nwork, &iinfo);
00936 if (iinfo != 0 && iinfo <= n + 2) {
00937 io___42.ciunit = *nounit;
00938 s_wsfe(&io___42);
00939 do_fio(&c__1, "SHSEQR(S)", (ftnlen)9);
00940 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00941 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00942 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00943 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00944 e_wsfe();
00945 *info = abs(iinfo);
00946 goto L250;
00947 }
00948
00949
00950
00951
00952 slacpy_(" ", &n, &n, &h__[h_offset], lda, &t1[t1_offset], lda);
00953 slacpy_(" ", &n, &n, &u[u_offset], ldu, &uz[uz_offset], lda);
00954
00955 shseqr_("S", "V", &n, &ilo, &ihi, &t1[t1_offset], lda, &wr1[1], &
00956 wi1[1], &uz[uz_offset], ldu, &work[1], nwork, &iinfo);
00957 if (iinfo != 0 && iinfo <= n + 2) {
00958 io___43.ciunit = *nounit;
00959 s_wsfe(&io___43);
00960 do_fio(&c__1, "SHSEQR(V)", (ftnlen)9);
00961 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00962 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00963 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00964 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00965 e_wsfe();
00966 *info = abs(iinfo);
00967 goto L250;
00968 }
00969
00970
00971
00972 sgemm_("T", "N", &n, &n, &n, &c_b32, &u[u_offset], ldu, &uz[
00973 uz_offset], ldu, &c_b18, &z__[z_offset], ldu);
00974 ntest = 8;
00975
00976
00977
00978
00979 shst01_(&n, &ilo, &ihi, &h__[h_offset], lda, &t1[t1_offset], lda,
00980 &z__[z_offset], ldu, &work[1], nwork, &result[3]);
00981
00982
00983
00984
00985 shst01_(&n, &ilo, &ihi, &a[a_offset], lda, &t1[t1_offset], lda, &
00986 uz[uz_offset], ldu, &work[1], nwork, &result[5]);
00987
00988
00989
00990 sget10_(&n, &n, &t2[t2_offset], lda, &t1[t1_offset], lda, &work[1]
00991 , &result[7]);
00992
00993
00994
00995 temp1 = 0.f;
00996 temp2 = 0.f;
00997 i__3 = n;
00998 for (j = 1; j <= i__3; ++j) {
00999
01000 r__5 = temp1, r__6 = (r__1 = wr1[j], dabs(r__1)) + (r__2 =
01001 wi1[j], dabs(r__2)), r__5 = max(r__5,r__6), r__6 = (
01002 r__3 = wr3[j], dabs(r__3)) + (r__4 = wi3[j], dabs(
01003 r__4));
01004 temp1 = dmax(r__5,r__6);
01005
01006 r__3 = temp2, r__4 = (r__1 = wr1[j] - wr3[j], dabs(r__1)) + (
01007 r__2 = wr1[j] - wr3[j], dabs(r__2));
01008 temp2 = dmax(r__3,r__4);
01009
01010 }
01011
01012
01013 r__1 = unfl, r__2 = ulp * dmax(temp1,temp2);
01014 result[8] = temp2 / dmax(r__1,r__2);
01015
01016
01017
01018
01019
01020 ntest = 9;
01021 result[9] = ulpinv;
01022
01023
01024
01025 nselc = 0;
01026 nselr = 0;
01027 j = n;
01028 L140:
01029 if (wi1[j] == 0.f) {
01030
01031 i__3 = n / 4;
01032 if (nselr < max(i__3,1)) {
01033 ++nselr;
01034 select[j] = TRUE_;
01035 } else {
01036 select[j] = FALSE_;
01037 }
01038 --j;
01039 } else {
01040
01041 i__3 = n / 4;
01042 if (nselc < max(i__3,1)) {
01043 ++nselc;
01044 select[j] = TRUE_;
01045 select[j - 1] = FALSE_;
01046 } else {
01047 select[j] = FALSE_;
01048 select[j - 1] = FALSE_;
01049 }
01050 j += -2;
01051 }
01052 if (j > 0) {
01053 goto L140;
01054 }
01055
01056 strevc_("Right", "All", &select[1], &n, &t1[t1_offset], lda,
01057 dumma, ldu, &evectr[evectr_offset], ldu, &n, &in, &work[1]
01058 , &iinfo);
01059 if (iinfo != 0) {
01060 io___50.ciunit = *nounit;
01061 s_wsfe(&io___50);
01062 do_fio(&c__1, "STREVC(R,A)", (ftnlen)11);
01063 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01064 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01065 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01066 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01067 e_wsfe();
01068 *info = abs(iinfo);
01069 goto L250;
01070 }
01071
01072
01073
01074 sget22_("N", "N", "N", &n, &t1[t1_offset], lda, &evectr[
01075 evectr_offset], ldu, &wr1[1], &wi1[1], &work[1], dumma);
01076 result[9] = dumma[0];
01077 if (dumma[1] > *thresh) {
01078 io___51.ciunit = *nounit;
01079 s_wsfe(&io___51);
01080 do_fio(&c__1, "Right", (ftnlen)5);
01081 do_fio(&c__1, "STREVC", (ftnlen)6);
01082 do_fio(&c__1, (char *)&dumma[1], (ftnlen)sizeof(real));
01083 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01084 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01085 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01086 e_wsfe();
01087 }
01088
01089
01090
01091
01092 strevc_("Right", "Some", &select[1], &n, &t1[t1_offset], lda,
01093 dumma, ldu, &evectl[evectl_offset], ldu, &n, &in, &work[1]
01094 , &iinfo);
01095 if (iinfo != 0) {
01096 io___52.ciunit = *nounit;
01097 s_wsfe(&io___52);
01098 do_fio(&c__1, "STREVC(R,S)", (ftnlen)11);
01099 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01100 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01101 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01102 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01103 e_wsfe();
01104 *info = abs(iinfo);
01105 goto L250;
01106 }
01107
01108 k = 1;
01109 match = TRUE_;
01110 i__3 = n;
01111 for (j = 1; j <= i__3; ++j) {
01112 if (select[j] && wi1[j] == 0.f) {
01113 i__4 = n;
01114 for (jj = 1; jj <= i__4; ++jj) {
01115 if (evectr[jj + j * evectr_dim1] != evectl[jj + k *
01116 evectl_dim1]) {
01117 match = FALSE_;
01118 goto L180;
01119 }
01120
01121 }
01122 ++k;
01123 } else if (select[j] && wi1[j] != 0.f) {
01124 i__4 = n;
01125 for (jj = 1; jj <= i__4; ++jj) {
01126 if (evectr[jj + j * evectr_dim1] != evectl[jj + k *
01127 evectl_dim1] || evectr[jj + (j + 1) *
01128 evectr_dim1] != evectl[jj + (k + 1) *
01129 evectl_dim1]) {
01130 match = FALSE_;
01131 goto L180;
01132 }
01133
01134 }
01135 k += 2;
01136 }
01137
01138 }
01139 L180:
01140 if (! match) {
01141 io___56.ciunit = *nounit;
01142 s_wsfe(&io___56);
01143 do_fio(&c__1, "Right", (ftnlen)5);
01144 do_fio(&c__1, "STREVC", (ftnlen)6);
01145 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01146 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01147 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01148 e_wsfe();
01149 }
01150
01151
01152
01153 ntest = 10;
01154 result[10] = ulpinv;
01155 strevc_("Left", "All", &select[1], &n, &t1[t1_offset], lda, &
01156 evectl[evectl_offset], ldu, dumma, ldu, &n, &in, &work[1],
01157 &iinfo);
01158 if (iinfo != 0) {
01159 io___57.ciunit = *nounit;
01160 s_wsfe(&io___57);
01161 do_fio(&c__1, "STREVC(L,A)", (ftnlen)11);
01162 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01163 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01164 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01165 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01166 e_wsfe();
01167 *info = abs(iinfo);
01168 goto L250;
01169 }
01170
01171
01172
01173 sget22_("Trans", "N", "Conj", &n, &t1[t1_offset], lda, &evectl[
01174 evectl_offset], ldu, &wr1[1], &wi1[1], &work[1], &dumma[2]
01175 );
01176 result[10] = dumma[2];
01177 if (dumma[3] > *thresh) {
01178 io___58.ciunit = *nounit;
01179 s_wsfe(&io___58);
01180 do_fio(&c__1, "Left", (ftnlen)4);
01181 do_fio(&c__1, "STREVC", (ftnlen)6);
01182 do_fio(&c__1, (char *)&dumma[3], (ftnlen)sizeof(real));
01183 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01184 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01185 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01186 e_wsfe();
01187 }
01188
01189
01190
01191
01192 strevc_("Left", "Some", &select[1], &n, &t1[t1_offset], lda, &
01193 evectr[evectr_offset], ldu, dumma, ldu, &n, &in, &work[1],
01194 &iinfo);
01195 if (iinfo != 0) {
01196 io___59.ciunit = *nounit;
01197 s_wsfe(&io___59);
01198 do_fio(&c__1, "STREVC(L,S)", (ftnlen)11);
01199 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01200 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01201 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01202 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01203 e_wsfe();
01204 *info = abs(iinfo);
01205 goto L250;
01206 }
01207
01208 k = 1;
01209 match = TRUE_;
01210 i__3 = n;
01211 for (j = 1; j <= i__3; ++j) {
01212 if (select[j] && wi1[j] == 0.f) {
01213 i__4 = n;
01214 for (jj = 1; jj <= i__4; ++jj) {
01215 if (evectl[jj + j * evectl_dim1] != evectr[jj + k *
01216 evectr_dim1]) {
01217 match = FALSE_;
01218 goto L220;
01219 }
01220
01221 }
01222 ++k;
01223 } else if (select[j] && wi1[j] != 0.f) {
01224 i__4 = n;
01225 for (jj = 1; jj <= i__4; ++jj) {
01226 if (evectl[jj + j * evectl_dim1] != evectr[jj + k *
01227 evectr_dim1] || evectl[jj + (j + 1) *
01228 evectl_dim1] != evectr[jj + (k + 1) *
01229 evectr_dim1]) {
01230 match = FALSE_;
01231 goto L220;
01232 }
01233
01234 }
01235 k += 2;
01236 }
01237
01238 }
01239 L220:
01240 if (! match) {
01241 io___60.ciunit = *nounit;
01242 s_wsfe(&io___60);
01243 do_fio(&c__1, "Left", (ftnlen)4);
01244 do_fio(&c__1, "STREVC", (ftnlen)6);
01245 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01246 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01247 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01248 e_wsfe();
01249 }
01250
01251
01252
01253 ntest = 11;
01254 result[11] = ulpinv;
01255 i__3 = n;
01256 for (j = 1; j <= i__3; ++j) {
01257 select[j] = TRUE_;
01258
01259 }
01260
01261 shsein_("Right", "Qr", "Ninitv", &select[1], &n, &h__[h_offset],
01262 lda, &wr3[1], &wi3[1], dumma, ldu, &evectx[evectx_offset],
01263 ldu, &n1, &in, &work[1], &iwork[1], &iwork[1], &iinfo);
01264 if (iinfo != 0) {
01265 io___61.ciunit = *nounit;
01266 s_wsfe(&io___61);
01267 do_fio(&c__1, "SHSEIN(R)", (ftnlen)9);
01268 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01269 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01270 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01271 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01272 e_wsfe();
01273 *info = abs(iinfo);
01274 if (iinfo < 0) {
01275 goto L250;
01276 }
01277 } else {
01278
01279
01280
01281
01282
01283 sget22_("N", "N", "N", &n, &h__[h_offset], lda, &evectx[
01284 evectx_offset], ldu, &wr3[1], &wi3[1], &work[1],
01285 dumma);
01286 if (dumma[0] < ulpinv) {
01287 result[11] = dumma[0] * aninv;
01288 }
01289 if (dumma[1] > *thresh) {
01290 io___62.ciunit = *nounit;
01291 s_wsfe(&io___62);
01292 do_fio(&c__1, "Right", (ftnlen)5);
01293 do_fio(&c__1, "SHSEIN", (ftnlen)6);
01294 do_fio(&c__1, (char *)&dumma[1], (ftnlen)sizeof(real));
01295 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01296 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01297 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
01298 ;
01299 e_wsfe();
01300 }
01301 }
01302
01303
01304
01305 ntest = 12;
01306 result[12] = ulpinv;
01307 i__3 = n;
01308 for (j = 1; j <= i__3; ++j) {
01309 select[j] = TRUE_;
01310
01311 }
01312
01313 shsein_("Left", "Qr", "Ninitv", &select[1], &n, &h__[h_offset],
01314 lda, &wr3[1], &wi3[1], &evecty[evecty_offset], ldu, dumma,
01315 ldu, &n1, &in, &work[1], &iwork[1], &iwork[1], &iinfo);
01316 if (iinfo != 0) {
01317 io___63.ciunit = *nounit;
01318 s_wsfe(&io___63);
01319 do_fio(&c__1, "SHSEIN(L)", (ftnlen)9);
01320 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01321 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01322 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01323 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01324 e_wsfe();
01325 *info = abs(iinfo);
01326 if (iinfo < 0) {
01327 goto L250;
01328 }
01329 } else {
01330
01331
01332
01333
01334
01335 sget22_("C", "N", "C", &n, &h__[h_offset], lda, &evecty[
01336 evecty_offset], ldu, &wr3[1], &wi3[1], &work[1], &
01337 dumma[2]);
01338 if (dumma[2] < ulpinv) {
01339 result[12] = dumma[2] * aninv;
01340 }
01341 if (dumma[3] > *thresh) {
01342 io___64.ciunit = *nounit;
01343 s_wsfe(&io___64);
01344 do_fio(&c__1, "Left", (ftnlen)4);
01345 do_fio(&c__1, "SHSEIN", (ftnlen)6);
01346 do_fio(&c__1, (char *)&dumma[3], (ftnlen)sizeof(real));
01347 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01348 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01349 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
01350 ;
01351 e_wsfe();
01352 }
01353 }
01354
01355
01356
01357 ntest = 13;
01358 result[13] = ulpinv;
01359
01360 sormhr_("Left", "No transpose", &n, &n, &ilo, &ihi, &uu[uu_offset]
01361 , ldu, &tau[1], &evectx[evectx_offset], ldu, &work[1],
01362 nwork, &iinfo);
01363 if (iinfo != 0) {
01364 io___65.ciunit = *nounit;
01365 s_wsfe(&io___65);
01366 do_fio(&c__1, "SORMHR(R)", (ftnlen)9);
01367 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01368 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01369 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01370 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01371 e_wsfe();
01372 *info = abs(iinfo);
01373 if (iinfo < 0) {
01374 goto L250;
01375 }
01376 } else {
01377
01378
01379
01380
01381
01382 sget22_("N", "N", "N", &n, &a[a_offset], lda, &evectx[
01383 evectx_offset], ldu, &wr3[1], &wi3[1], &work[1],
01384 dumma);
01385 if (dumma[0] < ulpinv) {
01386 result[13] = dumma[0] * aninv;
01387 }
01388 }
01389
01390
01391
01392 ntest = 14;
01393 result[14] = ulpinv;
01394
01395 sormhr_("Left", "No transpose", &n, &n, &ilo, &ihi, &uu[uu_offset]
01396 , ldu, &tau[1], &evecty[evecty_offset], ldu, &work[1],
01397 nwork, &iinfo);
01398 if (iinfo != 0) {
01399 io___66.ciunit = *nounit;
01400 s_wsfe(&io___66);
01401 do_fio(&c__1, "SORMHR(L)", (ftnlen)9);
01402 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01403 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01404 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01405 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01406 e_wsfe();
01407 *info = abs(iinfo);
01408 if (iinfo < 0) {
01409 goto L250;
01410 }
01411 } else {
01412
01413
01414
01415
01416
01417 sget22_("C", "N", "C", &n, &a[a_offset], lda, &evecty[
01418 evecty_offset], ldu, &wr3[1], &wi3[1], &work[1], &
01419 dumma[2]);
01420 if (dumma[2] < ulpinv) {
01421 result[14] = dumma[2] * aninv;
01422 }
01423 }
01424
01425
01426
01427 L250:
01428
01429 ntestt += ntest;
01430 slafts_("SHS", &n, &n, &jtype, &ntest, &result[1], ioldsd, thresh,
01431 nounit, &nerrs);
01432
01433 L260:
01434 ;
01435 }
01436 L270:
01437 ;
01438 }
01439
01440
01441
01442 slasum_("SHS", nounit, &nerrs, &ntestt);
01443
01444 return 0;
01445
01446
01447
01448
01449 }