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 doublereal c_b18 = 0.;
00019 static integer c__0 = 0;
00020 static doublereal c_b32 = 1.;
00021 static integer c__4 = 4;
00022 static integer c__6 = 6;
00023 static integer c__1 = 1;
00024
00025 int dchkhs_(integer *nsizes, integer *nn, integer *ntypes,
00026 logical *dotype, integer *iseed, doublereal *thresh, integer *nounit,
00027 doublereal *a, integer *lda, doublereal *h__, doublereal *t1,
00028 doublereal *t2, doublereal *u, integer *ldu, doublereal *z__,
00029 doublereal *uz, doublereal *wr1, doublereal *wi1, doublereal *wr3,
00030 doublereal *wi3, doublereal *evectl, doublereal *evectr, doublereal *
00031 evecty, doublereal *evectx, doublereal *uu, doublereal *tau,
00032 doublereal *work, integer *nwork, integer *iwork, logical *select,
00033 doublereal *result, integer *info)
00034 {
00035
00036
00037 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 };
00038 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 };
00039 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 };
00040 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 };
00041
00042
00043 static char fmt_9999[] = "(\002 DCHKHS: \002,a,\002 returned INFO=\002,i"
00044 "6,\002.\002,/9x,\002N=\002,i6,\002, JTYPE=\002,i6,\002, ISEED="
00045 "(\002,3(i5,\002,\002),i5,\002)\002)";
00046 static char fmt_9998[] = "(\002 DCHKHS: \002,a,\002 Eigenvectors from"
00047 " \002,a,\002 incorrectly \002,\002normalized.\002,/\002 Bits of "
00048 "error=\002,0p,g10.3,\002,\002,9x,\002N=\002,i6,\002, JTYPE=\002,"
00049 "i6,\002, ISEED=(\002,3(i5,\002,\002),i5,\002)\002)";
00050 static char fmt_9997[] = "(\002 DCHKHS: Selected \002,a,\002 Eigenvector"
00051 "s from \002,a,\002 do not match other eigenvectors \002,9x,\002N="
00052 "\002,i6,\002, JTYPE=\002,i6,\002, ISEED=(\002,3(i5,\002,\002),i5,"
00053 "\002)\002)";
00054
00055
00056 integer a_dim1, a_offset, evectl_dim1, evectl_offset, evectr_dim1,
00057 evectr_offset, evectx_dim1, evectx_offset, evecty_dim1,
00058 evecty_offset, h_dim1, h_offset, t1_dim1, t1_offset, t2_dim1,
00059 t2_offset, u_dim1, u_offset, uu_dim1, uu_offset, uz_dim1,
00060 uz_offset, z_dim1, z_offset, i__1, i__2, i__3, i__4;
00061 doublereal d__1, d__2, d__3, d__4, d__5, d__6;
00062
00063
00064 double sqrt(doublereal);
00065 integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00066
00067
00068 integer i__, j, k, n, n1, jj, in, ihi, ilo;
00069 doublereal ulp, cond;
00070 integer jcol, nmax;
00071 doublereal unfl, ovfl, temp1, temp2;
00072 logical badnn;
00073 extern int dget10_(integer *, integer *, doublereal *,
00074 integer *, doublereal *, integer *, doublereal *, doublereal *),
00075 dget22_(char *, char *, char *, integer *, doublereal *, integer *
00076 , doublereal *, integer *, doublereal *, doublereal *, doublereal
00077 *, doublereal *), dgemm_(char *, char *,
00078 integer *, integer *, integer *, doublereal *, doublereal *,
00079 integer *, doublereal *, integer *, doublereal *, doublereal *,
00080 integer *);
00081 logical match;
00082 integer imode;
00083 doublereal dumma[6];
00084 integer iinfo, nselc;
00085 doublereal conds;
00086 extern int dhst01_(integer *, integer *, integer *,
00087 doublereal *, integer *, doublereal *, integer *, doublereal *,
00088 integer *, doublereal *, integer *, doublereal *);
00089 doublereal aninv, anorm;
00090 extern int dcopy_(integer *, doublereal *, integer *,
00091 doublereal *, integer *);
00092 integer nmats, nselr, jsize, nerrs, itype, jtype, ntest;
00093 doublereal rtulp;
00094 extern int dlabad_(doublereal *, doublereal *);
00095 extern doublereal dlamch_(char *);
00096 extern int dgehrd_(integer *, integer *, integer *,
00097 doublereal *, integer *, doublereal *, doublereal *, integer *,
00098 integer *);
00099 char adumma[1*1];
00100 extern int dlatme_(integer *, char *, integer *,
00101 doublereal *, integer *, doublereal *, doublereal *, char *, char
00102 *, char *, char *, doublereal *, integer *, doublereal *, integer
00103 *, integer *, doublereal *, doublereal *, integer *, doublereal *,
00104 integer *), dhsein_(char
00105 *, char *, char *, logical *, integer *, doublereal *, integer *,
00106 doublereal *, doublereal *, doublereal *, integer *, doublereal *,
00107 integer *, integer *, integer *, doublereal *, integer *,
00108 integer *, integer *);
00109 integer idumma[1];
00110 extern int dlacpy_(char *, integer *, integer *,
00111 doublereal *, integer *, doublereal *, integer *);
00112 integer ioldsd[4];
00113 extern int dlafts_(char *, integer *, integer *, integer
00114 *, integer *, doublereal *, integer *, doublereal *, integer *,
00115 integer *), dlaset_(char *, integer *, integer *,
00116 doublereal *, doublereal *, doublereal *, integer *),
00117 dlatmr_(integer *, integer *, char *, integer *, char *,
00118 doublereal *, integer *, doublereal *, doublereal *, char *, char
00119 *, doublereal *, integer *, doublereal *, doublereal *, integer *,
00120 doublereal *, char *, integer *, integer *, integer *,
00121 doublereal *, doublereal *, char *, doublereal *, integer *,
00122 integer *, integer *), dlasum_(char *, integer *, integer *, integer *),
00123 dhseqr_(char *, char *, integer *, integer *, integer *,
00124 doublereal *, integer *, doublereal *, doublereal *, doublereal *,
00125 integer *, doublereal *, integer *, integer *),
00126 dlatms_(integer *, integer *, char *, integer *, char *,
00127 doublereal *, integer *, doublereal *, doublereal *, integer *,
00128 integer *, char *, doublereal *, integer *, doublereal *, integer
00129 *), dorghr_(integer *, integer *, integer
00130 *, doublereal *, integer *, doublereal *, doublereal *, integer *,
00131 integer *), dtrevc_(char *, char *, logical *, integer *,
00132 doublereal *, integer *, doublereal *, integer *, doublereal *,
00133 integer *, integer *, integer *, doublereal *, integer *), dormhr_(char *, char *, integer *, integer *, integer *,
00134 integer *, doublereal *, integer *, doublereal *, doublereal *,
00135 integer *, doublereal *, integer *, integer *),
00136 xerbla_(char *, integer *);
00137 doublereal rtunfl, rtovfl, rtulpi, ulpinv;
00138 integer mtypes, ntestt;
00139
00140
00141 static cilist io___36 = { 0, 0, 0, fmt_9999, 0 };
00142 static cilist io___39 = { 0, 0, 0, fmt_9999, 0 };
00143 static cilist io___41 = { 0, 0, 0, fmt_9999, 0 };
00144 static cilist io___42 = { 0, 0, 0, fmt_9999, 0 };
00145 static cilist io___43 = { 0, 0, 0, fmt_9999, 0 };
00146 static cilist io___50 = { 0, 0, 0, fmt_9999, 0 };
00147 static cilist io___51 = { 0, 0, 0, fmt_9998, 0 };
00148 static cilist io___52 = { 0, 0, 0, fmt_9999, 0 };
00149 static cilist io___56 = { 0, 0, 0, fmt_9997, 0 };
00150 static cilist io___57 = { 0, 0, 0, fmt_9999, 0 };
00151 static cilist io___58 = { 0, 0, 0, fmt_9998, 0 };
00152 static cilist io___59 = { 0, 0, 0, fmt_9999, 0 };
00153 static cilist io___60 = { 0, 0, 0, fmt_9997, 0 };
00154 static cilist io___61 = { 0, 0, 0, fmt_9999, 0 };
00155 static cilist io___62 = { 0, 0, 0, fmt_9998, 0 };
00156 static cilist io___63 = { 0, 0, 0, fmt_9999, 0 };
00157 static cilist io___64 = { 0, 0, 0, fmt_9998, 0 };
00158 static cilist io___65 = { 0, 0, 0, fmt_9999, 0 };
00159 static cilist io___66 = { 0, 0, 0, fmt_9999, 0 };
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
00528
00529
00530
00531
00532
00533
00534
00535
00536
00537
00538 --nn;
00539 --dotype;
00540 --iseed;
00541 t2_dim1 = *lda;
00542 t2_offset = 1 + t2_dim1;
00543 t2 -= t2_offset;
00544 t1_dim1 = *lda;
00545 t1_offset = 1 + t1_dim1;
00546 t1 -= t1_offset;
00547 h_dim1 = *lda;
00548 h_offset = 1 + h_dim1;
00549 h__ -= h_offset;
00550 a_dim1 = *lda;
00551 a_offset = 1 + a_dim1;
00552 a -= a_offset;
00553 uu_dim1 = *ldu;
00554 uu_offset = 1 + uu_dim1;
00555 uu -= uu_offset;
00556 evectx_dim1 = *ldu;
00557 evectx_offset = 1 + evectx_dim1;
00558 evectx -= evectx_offset;
00559 evecty_dim1 = *ldu;
00560 evecty_offset = 1 + evecty_dim1;
00561 evecty -= evecty_offset;
00562 evectr_dim1 = *ldu;
00563 evectr_offset = 1 + evectr_dim1;
00564 evectr -= evectr_offset;
00565 evectl_dim1 = *ldu;
00566 evectl_offset = 1 + evectl_dim1;
00567 evectl -= evectl_offset;
00568 uz_dim1 = *ldu;
00569 uz_offset = 1 + uz_dim1;
00570 uz -= uz_offset;
00571 z_dim1 = *ldu;
00572 z_offset = 1 + z_dim1;
00573 z__ -= z_offset;
00574 u_dim1 = *ldu;
00575 u_offset = 1 + u_dim1;
00576 u -= u_offset;
00577 --wr1;
00578 --wi1;
00579 --wr3;
00580 --wi3;
00581 --tau;
00582 --work;
00583 --iwork;
00584 --select;
00585 --result;
00586
00587
00588
00589
00590
00591
00592
00593 ntestt = 0;
00594 *info = 0;
00595
00596 badnn = FALSE_;
00597 nmax = 0;
00598 i__1 = *nsizes;
00599 for (j = 1; j <= i__1; ++j) {
00600
00601 i__2 = nmax, i__3 = nn[j];
00602 nmax = max(i__2,i__3);
00603 if (nn[j] < 0) {
00604 badnn = TRUE_;
00605 }
00606
00607 }
00608
00609
00610
00611 if (*nsizes < 0) {
00612 *info = -1;
00613 } else if (badnn) {
00614 *info = -2;
00615 } else if (*ntypes < 0) {
00616 *info = -3;
00617 } else if (*thresh < 0.) {
00618 *info = -6;
00619 } else if (*lda <= 1 || *lda < nmax) {
00620 *info = -9;
00621 } else if (*ldu <= 1 || *ldu < nmax) {
00622 *info = -14;
00623 } else if ((nmax << 2) * nmax + 2 > *nwork) {
00624 *info = -28;
00625 }
00626
00627 if (*info != 0) {
00628 i__1 = -(*info);
00629 xerbla_("DCHKHS", &i__1);
00630 return 0;
00631 }
00632
00633
00634
00635 if (*nsizes == 0 || *ntypes == 0) {
00636 return 0;
00637 }
00638
00639
00640
00641 unfl = dlamch_("Safe minimum");
00642 ovfl = dlamch_("Overflow");
00643 dlabad_(&unfl, &ovfl);
00644 ulp = dlamch_("Epsilon") * dlamch_("Base");
00645 ulpinv = 1. / ulp;
00646 rtunfl = sqrt(unfl);
00647 rtovfl = sqrt(ovfl);
00648 rtulp = sqrt(ulp);
00649 rtulpi = 1. / rtulp;
00650
00651
00652
00653 nerrs = 0;
00654 nmats = 0;
00655
00656 i__1 = *nsizes;
00657 for (jsize = 1; jsize <= i__1; ++jsize) {
00658 n = nn[jsize];
00659 if (n == 0) {
00660 goto L270;
00661 }
00662 n1 = max(1,n);
00663 aninv = 1. / (doublereal) n1;
00664
00665 if (*nsizes != 1) {
00666 mtypes = min(21,*ntypes);
00667 } else {
00668 mtypes = min(22,*ntypes);
00669 }
00670
00671 i__2 = mtypes;
00672 for (jtype = 1; jtype <= i__2; ++jtype) {
00673 if (! dotype[jtype]) {
00674 goto L260;
00675 }
00676 ++nmats;
00677 ntest = 0;
00678
00679
00680
00681 for (j = 1; j <= 4; ++j) {
00682 ioldsd[j - 1] = iseed[j];
00683
00684 }
00685
00686
00687
00688 for (j = 1; j <= 14; ++j) {
00689 result[j] = 0.;
00690
00691 }
00692
00693
00694
00695
00696
00697
00698
00699
00700
00701
00702
00703
00704
00705
00706
00707
00708
00709 if (mtypes > 21) {
00710 goto L100;
00711 }
00712
00713 itype = ktype[jtype - 1];
00714 imode = kmode[jtype - 1];
00715
00716
00717
00718 switch (kmagn[jtype - 1]) {
00719 case 1: goto L40;
00720 case 2: goto L50;
00721 case 3: goto L60;
00722 }
00723
00724 L40:
00725 anorm = 1.;
00726 goto L70;
00727
00728 L50:
00729 anorm = rtovfl * ulp * aninv;
00730 goto L70;
00731
00732 L60:
00733 anorm = rtunfl * n * ulpinv;
00734 goto L70;
00735
00736 L70:
00737
00738 dlaset_("Full", lda, &n, &c_b18, &c_b18, &a[a_offset], lda);
00739 iinfo = 0;
00740 cond = ulpinv;
00741
00742
00743
00744 if (itype == 1) {
00745
00746
00747
00748 iinfo = 0;
00749
00750 } else if (itype == 2) {
00751
00752
00753
00754 i__3 = n;
00755 for (jcol = 1; jcol <= i__3; ++jcol) {
00756 a[jcol + jcol * a_dim1] = anorm;
00757
00758 }
00759
00760 } else if (itype == 3) {
00761
00762
00763
00764 i__3 = n;
00765 for (jcol = 1; jcol <= i__3; ++jcol) {
00766 a[jcol + jcol * a_dim1] = anorm;
00767 if (jcol > 1) {
00768 a[jcol + (jcol - 1) * a_dim1] = 1.;
00769 }
00770
00771 }
00772
00773 } else if (itype == 4) {
00774
00775
00776
00777 dlatms_(&n, &n, "S", &iseed[1], "S", &work[1], &imode, &cond,
00778 &anorm, &c__0, &c__0, "N", &a[a_offset], lda, &work[n
00779 + 1], &iinfo);
00780
00781 } else if (itype == 5) {
00782
00783
00784
00785 dlatms_(&n, &n, "S", &iseed[1], "S", &work[1], &imode, &cond,
00786 &anorm, &n, &n, "N", &a[a_offset], lda, &work[n + 1],
00787 &iinfo);
00788
00789 } else if (itype == 6) {
00790
00791
00792
00793 if (kconds[jtype - 1] == 1) {
00794 conds = 1.;
00795 } else if (kconds[jtype - 1] == 2) {
00796 conds = rtulpi;
00797 } else {
00798 conds = 0.;
00799 }
00800
00801 *(unsigned char *)&adumma[0] = ' ';
00802 dlatme_(&n, "S", &iseed[1], &work[1], &imode, &cond, &c_b32,
00803 adumma, "T", "T", "T", &work[n + 1], &c__4, &conds, &
00804 n, &n, &anorm, &a[a_offset], lda, &work[(n << 1) + 1],
00805 &iinfo);
00806
00807 } else if (itype == 7) {
00808
00809
00810
00811 dlatmr_(&n, &n, "S", &iseed[1], "S", &work[1], &c__6, &c_b32,
00812 &c_b32, "T", "N", &work[n + 1], &c__1, &c_b32, &work[(
00813 n << 1) + 1], &c__1, &c_b32, "N", idumma, &c__0, &
00814 c__0, &c_b18, &anorm, "NO", &a[a_offset], lda, &iwork[
00815 1], &iinfo);
00816
00817 } else if (itype == 8) {
00818
00819
00820
00821 dlatmr_(&n, &n, "S", &iseed[1], "S", &work[1], &c__6, &c_b32,
00822 &c_b32, "T", "N", &work[n + 1], &c__1, &c_b32, &work[(
00823 n << 1) + 1], &c__1, &c_b32, "N", idumma, &n, &n, &
00824 c_b18, &anorm, "NO", &a[a_offset], lda, &iwork[1], &
00825 iinfo);
00826
00827 } else if (itype == 9) {
00828
00829
00830
00831 dlatmr_(&n, &n, "S", &iseed[1], "N", &work[1], &c__6, &c_b32,
00832 &c_b32, "T", "N", &work[n + 1], &c__1, &c_b32, &work[(
00833 n << 1) + 1], &c__1, &c_b32, "N", idumma, &n, &n, &
00834 c_b18, &anorm, "NO", &a[a_offset], lda, &iwork[1], &
00835 iinfo);
00836
00837 } else if (itype == 10) {
00838
00839
00840
00841 dlatmr_(&n, &n, "S", &iseed[1], "N", &work[1], &c__6, &c_b32,
00842 &c_b32, "T", "N", &work[n + 1], &c__1, &c_b32, &work[(
00843 n << 1) + 1], &c__1, &c_b32, "N", idumma, &n, &c__0, &
00844 c_b18, &anorm, "NO", &a[a_offset], lda, &iwork[1], &
00845 iinfo);
00846
00847 } else {
00848
00849 iinfo = 1;
00850 }
00851
00852 if (iinfo != 0) {
00853 io___36.ciunit = *nounit;
00854 s_wsfe(&io___36);
00855 do_fio(&c__1, "Generator", (ftnlen)9);
00856 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00857 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00858 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00859 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00860 e_wsfe();
00861 *info = abs(iinfo);
00862 return 0;
00863 }
00864
00865 L100:
00866
00867
00868
00869 dlacpy_(" ", &n, &n, &a[a_offset], lda, &h__[h_offset], lda);
00870
00871 ntest = 1;
00872
00873 ilo = 1;
00874 ihi = n;
00875
00876 i__3 = *nwork - n;
00877 dgehrd_(&n, &ilo, &ihi, &h__[h_offset], lda, &work[1], &work[n +
00878 1], &i__3, &iinfo);
00879
00880 if (iinfo != 0) {
00881 result[1] = ulpinv;
00882 io___39.ciunit = *nounit;
00883 s_wsfe(&io___39);
00884 do_fio(&c__1, "DGEHRD", (ftnlen)6);
00885 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00886 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00887 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00888 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00889 e_wsfe();
00890 *info = abs(iinfo);
00891 goto L250;
00892 }
00893
00894 i__3 = n - 1;
00895 for (j = 1; j <= i__3; ++j) {
00896 uu[j + 1 + j * uu_dim1] = 0.;
00897 i__4 = n;
00898 for (i__ = j + 2; i__ <= i__4; ++i__) {
00899 u[i__ + j * u_dim1] = h__[i__ + j * h_dim1];
00900 uu[i__ + j * uu_dim1] = h__[i__ + j * h_dim1];
00901 h__[i__ + j * h_dim1] = 0.;
00902
00903 }
00904
00905 }
00906 i__3 = n - 1;
00907 dcopy_(&i__3, &work[1], &c__1, &tau[1], &c__1);
00908 i__3 = *nwork - n;
00909 dorghr_(&n, &ilo, &ihi, &u[u_offset], ldu, &work[1], &work[n + 1],
00910 &i__3, &iinfo);
00911 ntest = 2;
00912
00913 dhst01_(&n, &ilo, &ihi, &a[a_offset], lda, &h__[h_offset], lda, &
00914 u[u_offset], ldu, &work[1], nwork, &result[1]);
00915
00916
00917
00918
00919
00920 dlacpy_(" ", &n, &n, &h__[h_offset], lda, &t2[t2_offset], lda);
00921 ntest = 3;
00922 result[3] = ulpinv;
00923
00924 dhseqr_("E", "N", &n, &ilo, &ihi, &t2[t2_offset], lda, &wr3[1], &
00925 wi3[1], &uz[uz_offset], ldu, &work[1], nwork, &iinfo);
00926 if (iinfo != 0) {
00927 io___41.ciunit = *nounit;
00928 s_wsfe(&io___41);
00929 do_fio(&c__1, "DHSEQR(E)", (ftnlen)9);
00930 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00931 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00932 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00933 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00934 e_wsfe();
00935 if (iinfo <= n + 2) {
00936 *info = abs(iinfo);
00937 goto L250;
00938 }
00939 }
00940
00941
00942
00943 dlacpy_(" ", &n, &n, &h__[h_offset], lda, &t2[t2_offset], lda);
00944
00945 dhseqr_("S", "N", &n, &ilo, &ihi, &t2[t2_offset], lda, &wr1[1], &
00946 wi1[1], &uz[uz_offset], ldu, &work[1], nwork, &iinfo);
00947 if (iinfo != 0 && iinfo <= n + 2) {
00948 io___42.ciunit = *nounit;
00949 s_wsfe(&io___42);
00950 do_fio(&c__1, "DHSEQR(S)", (ftnlen)9);
00951 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00952 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00953 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00954 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00955 e_wsfe();
00956 *info = abs(iinfo);
00957 goto L250;
00958 }
00959
00960
00961
00962
00963 dlacpy_(" ", &n, &n, &h__[h_offset], lda, &t1[t1_offset], lda);
00964 dlacpy_(" ", &n, &n, &u[u_offset], ldu, &uz[uz_offset], lda);
00965
00966 dhseqr_("S", "V", &n, &ilo, &ihi, &t1[t1_offset], lda, &wr1[1], &
00967 wi1[1], &uz[uz_offset], ldu, &work[1], nwork, &iinfo);
00968 if (iinfo != 0 && iinfo <= n + 2) {
00969 io___43.ciunit = *nounit;
00970 s_wsfe(&io___43);
00971 do_fio(&c__1, "DHSEQR(V)", (ftnlen)9);
00972 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00973 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00974 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00975 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00976 e_wsfe();
00977 *info = abs(iinfo);
00978 goto L250;
00979 }
00980
00981
00982
00983 dgemm_("T", "N", &n, &n, &n, &c_b32, &u[u_offset], ldu, &uz[
00984 uz_offset], ldu, &c_b18, &z__[z_offset], ldu);
00985 ntest = 8;
00986
00987
00988
00989
00990 dhst01_(&n, &ilo, &ihi, &h__[h_offset], lda, &t1[t1_offset], lda,
00991 &z__[z_offset], ldu, &work[1], nwork, &result[3]);
00992
00993
00994
00995
00996 dhst01_(&n, &ilo, &ihi, &a[a_offset], lda, &t1[t1_offset], lda, &
00997 uz[uz_offset], ldu, &work[1], nwork, &result[5]);
00998
00999
01000
01001 dget10_(&n, &n, &t2[t2_offset], lda, &t1[t1_offset], lda, &work[1]
01002 , &result[7]);
01003
01004
01005
01006 temp1 = 0.;
01007 temp2 = 0.;
01008 i__3 = n;
01009 for (j = 1; j <= i__3; ++j) {
01010
01011 d__5 = temp1, d__6 = (d__1 = wr1[j], abs(d__1)) + (d__2 = wi1[
01012 j], abs(d__2)), d__5 = max(d__5,d__6), d__6 = (d__3 =
01013 wr3[j], abs(d__3)) + (d__4 = wi3[j], abs(d__4));
01014 temp1 = max(d__5,d__6);
01015
01016 d__3 = temp2, d__4 = (d__1 = wr1[j] - wr3[j], abs(d__1)) + (
01017 d__2 = wr1[j] - wr3[j], abs(d__2));
01018 temp2 = max(d__3,d__4);
01019
01020 }
01021
01022
01023 d__1 = unfl, d__2 = ulp * max(temp1,temp2);
01024 result[8] = temp2 / max(d__1,d__2);
01025
01026
01027
01028
01029
01030 ntest = 9;
01031 result[9] = ulpinv;
01032
01033
01034
01035 nselc = 0;
01036 nselr = 0;
01037 j = n;
01038 L140:
01039 if (wi1[j] == 0.) {
01040
01041 i__3 = n / 4;
01042 if (nselr < max(i__3,1)) {
01043 ++nselr;
01044 select[j] = TRUE_;
01045 } else {
01046 select[j] = FALSE_;
01047 }
01048 --j;
01049 } else {
01050
01051 i__3 = n / 4;
01052 if (nselc < max(i__3,1)) {
01053 ++nselc;
01054 select[j] = TRUE_;
01055 select[j - 1] = FALSE_;
01056 } else {
01057 select[j] = FALSE_;
01058 select[j - 1] = FALSE_;
01059 }
01060 j += -2;
01061 }
01062 if (j > 0) {
01063 goto L140;
01064 }
01065
01066 dtrevc_("Right", "All", &select[1], &n, &t1[t1_offset], lda,
01067 dumma, ldu, &evectr[evectr_offset], ldu, &n, &in, &work[1]
01068 , &iinfo);
01069 if (iinfo != 0) {
01070 io___50.ciunit = *nounit;
01071 s_wsfe(&io___50);
01072 do_fio(&c__1, "DTREVC(R,A)", (ftnlen)11);
01073 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01074 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01075 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01076 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01077 e_wsfe();
01078 *info = abs(iinfo);
01079 goto L250;
01080 }
01081
01082
01083
01084 dget22_("N", "N", "N", &n, &t1[t1_offset], lda, &evectr[
01085 evectr_offset], ldu, &wr1[1], &wi1[1], &work[1], dumma);
01086 result[9] = dumma[0];
01087 if (dumma[1] > *thresh) {
01088 io___51.ciunit = *nounit;
01089 s_wsfe(&io___51);
01090 do_fio(&c__1, "Right", (ftnlen)5);
01091 do_fio(&c__1, "DTREVC", (ftnlen)6);
01092 do_fio(&c__1, (char *)&dumma[1], (ftnlen)sizeof(doublereal));
01093 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01094 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01095 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01096 e_wsfe();
01097 }
01098
01099
01100
01101
01102 dtrevc_("Right", "Some", &select[1], &n, &t1[t1_offset], lda,
01103 dumma, ldu, &evectl[evectl_offset], ldu, &n, &in, &work[1]
01104 , &iinfo);
01105 if (iinfo != 0) {
01106 io___52.ciunit = *nounit;
01107 s_wsfe(&io___52);
01108 do_fio(&c__1, "DTREVC(R,S)", (ftnlen)11);
01109 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01110 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01111 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01112 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01113 e_wsfe();
01114 *info = abs(iinfo);
01115 goto L250;
01116 }
01117
01118 k = 1;
01119 match = TRUE_;
01120 i__3 = n;
01121 for (j = 1; j <= i__3; ++j) {
01122 if (select[j] && wi1[j] == 0.) {
01123 i__4 = n;
01124 for (jj = 1; jj <= i__4; ++jj) {
01125 if (evectr[jj + j * evectr_dim1] != evectl[jj + k *
01126 evectl_dim1]) {
01127 match = FALSE_;
01128 goto L180;
01129 }
01130
01131 }
01132 ++k;
01133 } else if (select[j] && wi1[j] != 0.) {
01134 i__4 = n;
01135 for (jj = 1; jj <= i__4; ++jj) {
01136 if (evectr[jj + j * evectr_dim1] != evectl[jj + k *
01137 evectl_dim1] || evectr[jj + (j + 1) *
01138 evectr_dim1] != evectl[jj + (k + 1) *
01139 evectl_dim1]) {
01140 match = FALSE_;
01141 goto L180;
01142 }
01143
01144 }
01145 k += 2;
01146 }
01147
01148 }
01149 L180:
01150 if (! match) {
01151 io___56.ciunit = *nounit;
01152 s_wsfe(&io___56);
01153 do_fio(&c__1, "Right", (ftnlen)5);
01154 do_fio(&c__1, "DTREVC", (ftnlen)6);
01155 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01156 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01157 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01158 e_wsfe();
01159 }
01160
01161
01162
01163 ntest = 10;
01164 result[10] = ulpinv;
01165 dtrevc_("Left", "All", &select[1], &n, &t1[t1_offset], lda, &
01166 evectl[evectl_offset], ldu, dumma, ldu, &n, &in, &work[1],
01167 &iinfo);
01168 if (iinfo != 0) {
01169 io___57.ciunit = *nounit;
01170 s_wsfe(&io___57);
01171 do_fio(&c__1, "DTREVC(L,A)", (ftnlen)11);
01172 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01173 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01174 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01175 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01176 e_wsfe();
01177 *info = abs(iinfo);
01178 goto L250;
01179 }
01180
01181
01182
01183 dget22_("Trans", "N", "Conj", &n, &t1[t1_offset], lda, &evectl[
01184 evectl_offset], ldu, &wr1[1], &wi1[1], &work[1], &dumma[2]
01185 );
01186 result[10] = dumma[2];
01187 if (dumma[3] > *thresh) {
01188 io___58.ciunit = *nounit;
01189 s_wsfe(&io___58);
01190 do_fio(&c__1, "Left", (ftnlen)4);
01191 do_fio(&c__1, "DTREVC", (ftnlen)6);
01192 do_fio(&c__1, (char *)&dumma[3], (ftnlen)sizeof(doublereal));
01193 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01194 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01195 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01196 e_wsfe();
01197 }
01198
01199
01200
01201
01202 dtrevc_("Left", "Some", &select[1], &n, &t1[t1_offset], lda, &
01203 evectr[evectr_offset], ldu, dumma, ldu, &n, &in, &work[1],
01204 &iinfo);
01205 if (iinfo != 0) {
01206 io___59.ciunit = *nounit;
01207 s_wsfe(&io___59);
01208 do_fio(&c__1, "DTREVC(L,S)", (ftnlen)11);
01209 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01210 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01211 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01212 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01213 e_wsfe();
01214 *info = abs(iinfo);
01215 goto L250;
01216 }
01217
01218 k = 1;
01219 match = TRUE_;
01220 i__3 = n;
01221 for (j = 1; j <= i__3; ++j) {
01222 if (select[j] && wi1[j] == 0.) {
01223 i__4 = n;
01224 for (jj = 1; jj <= i__4; ++jj) {
01225 if (evectl[jj + j * evectl_dim1] != evectr[jj + k *
01226 evectr_dim1]) {
01227 match = FALSE_;
01228 goto L220;
01229 }
01230
01231 }
01232 ++k;
01233 } else if (select[j] && wi1[j] != 0.) {
01234 i__4 = n;
01235 for (jj = 1; jj <= i__4; ++jj) {
01236 if (evectl[jj + j * evectl_dim1] != evectr[jj + k *
01237 evectr_dim1] || evectl[jj + (j + 1) *
01238 evectl_dim1] != evectr[jj + (k + 1) *
01239 evectr_dim1]) {
01240 match = FALSE_;
01241 goto L220;
01242 }
01243
01244 }
01245 k += 2;
01246 }
01247
01248 }
01249 L220:
01250 if (! match) {
01251 io___60.ciunit = *nounit;
01252 s_wsfe(&io___60);
01253 do_fio(&c__1, "Left", (ftnlen)4);
01254 do_fio(&c__1, "DTREVC", (ftnlen)6);
01255 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01256 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01257 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01258 e_wsfe();
01259 }
01260
01261
01262
01263 ntest = 11;
01264 result[11] = ulpinv;
01265 i__3 = n;
01266 for (j = 1; j <= i__3; ++j) {
01267 select[j] = TRUE_;
01268
01269 }
01270
01271 dhsein_("Right", "Qr", "Ninitv", &select[1], &n, &h__[h_offset],
01272 lda, &wr3[1], &wi3[1], dumma, ldu, &evectx[evectx_offset],
01273 ldu, &n1, &in, &work[1], &iwork[1], &iwork[1], &iinfo);
01274 if (iinfo != 0) {
01275 io___61.ciunit = *nounit;
01276 s_wsfe(&io___61);
01277 do_fio(&c__1, "DHSEIN(R)", (ftnlen)9);
01278 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01279 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01280 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01281 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01282 e_wsfe();
01283 *info = abs(iinfo);
01284 if (iinfo < 0) {
01285 goto L250;
01286 }
01287 } else {
01288
01289
01290
01291
01292
01293 dget22_("N", "N", "N", &n, &h__[h_offset], lda, &evectx[
01294 evectx_offset], ldu, &wr3[1], &wi3[1], &work[1],
01295 dumma);
01296 if (dumma[0] < ulpinv) {
01297 result[11] = dumma[0] * aninv;
01298 }
01299 if (dumma[1] > *thresh) {
01300 io___62.ciunit = *nounit;
01301 s_wsfe(&io___62);
01302 do_fio(&c__1, "Right", (ftnlen)5);
01303 do_fio(&c__1, "DHSEIN", (ftnlen)6);
01304 do_fio(&c__1, (char *)&dumma[1], (ftnlen)sizeof(
01305 doublereal));
01306 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01307 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01308 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
01309 ;
01310 e_wsfe();
01311 }
01312 }
01313
01314
01315
01316 ntest = 12;
01317 result[12] = ulpinv;
01318 i__3 = n;
01319 for (j = 1; j <= i__3; ++j) {
01320 select[j] = TRUE_;
01321
01322 }
01323
01324 dhsein_("Left", "Qr", "Ninitv", &select[1], &n, &h__[h_offset],
01325 lda, &wr3[1], &wi3[1], &evecty[evecty_offset], ldu, dumma,
01326 ldu, &n1, &in, &work[1], &iwork[1], &iwork[1], &iinfo);
01327 if (iinfo != 0) {
01328 io___63.ciunit = *nounit;
01329 s_wsfe(&io___63);
01330 do_fio(&c__1, "DHSEIN(L)", (ftnlen)9);
01331 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01332 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01333 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01334 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01335 e_wsfe();
01336 *info = abs(iinfo);
01337 if (iinfo < 0) {
01338 goto L250;
01339 }
01340 } else {
01341
01342
01343
01344
01345
01346 dget22_("C", "N", "C", &n, &h__[h_offset], lda, &evecty[
01347 evecty_offset], ldu, &wr3[1], &wi3[1], &work[1], &
01348 dumma[2]);
01349 if (dumma[2] < ulpinv) {
01350 result[12] = dumma[2] * aninv;
01351 }
01352 if (dumma[3] > *thresh) {
01353 io___64.ciunit = *nounit;
01354 s_wsfe(&io___64);
01355 do_fio(&c__1, "Left", (ftnlen)4);
01356 do_fio(&c__1, "DHSEIN", (ftnlen)6);
01357 do_fio(&c__1, (char *)&dumma[3], (ftnlen)sizeof(
01358 doublereal));
01359 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01360 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01361 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
01362 ;
01363 e_wsfe();
01364 }
01365 }
01366
01367
01368
01369 ntest = 13;
01370 result[13] = ulpinv;
01371
01372 dormhr_("Left", "No transpose", &n, &n, &ilo, &ihi, &uu[uu_offset]
01373 , ldu, &tau[1], &evectx[evectx_offset], ldu, &work[1],
01374 nwork, &iinfo);
01375 if (iinfo != 0) {
01376 io___65.ciunit = *nounit;
01377 s_wsfe(&io___65);
01378 do_fio(&c__1, "DORMHR(R)", (ftnlen)9);
01379 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01380 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01381 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01382 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01383 e_wsfe();
01384 *info = abs(iinfo);
01385 if (iinfo < 0) {
01386 goto L250;
01387 }
01388 } else {
01389
01390
01391
01392
01393
01394 dget22_("N", "N", "N", &n, &a[a_offset], lda, &evectx[
01395 evectx_offset], ldu, &wr3[1], &wi3[1], &work[1],
01396 dumma);
01397 if (dumma[0] < ulpinv) {
01398 result[13] = dumma[0] * aninv;
01399 }
01400 }
01401
01402
01403
01404 ntest = 14;
01405 result[14] = ulpinv;
01406
01407 dormhr_("Left", "No transpose", &n, &n, &ilo, &ihi, &uu[uu_offset]
01408 , ldu, &tau[1], &evecty[evecty_offset], ldu, &work[1],
01409 nwork, &iinfo);
01410 if (iinfo != 0) {
01411 io___66.ciunit = *nounit;
01412 s_wsfe(&io___66);
01413 do_fio(&c__1, "DORMHR(L)", (ftnlen)9);
01414 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01415 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01416 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01417 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01418 e_wsfe();
01419 *info = abs(iinfo);
01420 if (iinfo < 0) {
01421 goto L250;
01422 }
01423 } else {
01424
01425
01426
01427
01428
01429 dget22_("C", "N", "C", &n, &a[a_offset], lda, &evecty[
01430 evecty_offset], ldu, &wr3[1], &wi3[1], &work[1], &
01431 dumma[2]);
01432 if (dumma[2] < ulpinv) {
01433 result[14] = dumma[2] * aninv;
01434 }
01435 }
01436
01437
01438
01439 L250:
01440
01441 ntestt += ntest;
01442 dlafts_("DHS", &n, &n, &jtype, &ntest, &result[1], ioldsd, thresh,
01443 nounit, &nerrs);
01444
01445 L260:
01446 ;
01447 }
01448 L270:
01449 ;
01450 }
01451
01452
01453
01454 dlasum_("DHS", nounit, &nerrs, &ntestt);
01455
01456 return 0;
01457
01458
01459
01460
01461 }