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 integer c__1 = 1;
00019 static integer c_n1 = -1;
00020 static integer c__0 = 0;
00021 static integer c__4 = 4;
00022 static doublereal c_b36 = 0.;
00023 static integer c__2 = 2;
00024 static doublereal c_b42 = 1.;
00025 static integer c__3 = 3;
00026 static logical c_true = TRUE_;
00027 static logical c_false = FALSE_;
00028
00029 int ddrvgg_(integer *nsizes, integer *nn, integer *ntypes,
00030 logical *dotype, integer *iseed, doublereal *thresh, doublereal *
00031 thrshn, integer *nounit, doublereal *a, integer *lda, doublereal *b,
00032 doublereal *s, doublereal *t, doublereal *s2, doublereal *t2,
00033 doublereal *q, integer *ldq, doublereal *z__, doublereal *alphr1,
00034 doublereal *alphi1, doublereal *beta1, doublereal *alphr2, doublereal
00035 *alphi2, doublereal *beta2, doublereal *vl, doublereal *vr,
00036 doublereal *work, integer *lwork, doublereal *result, integer *info)
00037 {
00038
00039
00040 static integer kclass[26] = { 1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,
00041 2,2,2,3 };
00042 static integer kbmagn[26] = { 1,1,1,1,1,1,1,1,3,2,3,2,2,3,1,1,1,1,1,1,1,3,
00043 2,3,2,1 };
00044 static integer ktrian[26] = { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,
00045 1,1,1,1 };
00046 static integer iasign[26] = { 0,0,0,0,0,0,2,0,2,2,0,0,2,2,2,0,2,0,0,0,2,2,
00047 2,2,2,0 };
00048 static integer ibsign[26] = { 0,0,0,0,0,0,0,2,0,0,2,2,0,0,2,0,2,0,0,0,0,0,
00049 0,0,0,0 };
00050 static integer kz1[6] = { 0,1,2,1,3,3 };
00051 static integer kz2[6] = { 0,0,1,2,1,1 };
00052 static integer kadd[6] = { 0,0,0,0,3,2 };
00053 static integer katype[26] = { 0,1,0,1,2,3,4,1,4,4,1,1,4,4,4,2,4,5,8,7,9,4,
00054 4,4,4,0 };
00055 static integer kbtype[26] = { 0,0,1,1,2,-3,1,4,1,1,4,4,1,1,-4,2,-4,8,8,8,
00056 8,8,8,8,8,0 };
00057 static integer kazero[26] = { 1,1,1,1,1,1,2,1,2,2,1,1,2,2,3,1,3,5,5,5,5,3,
00058 3,3,3,1 };
00059 static integer kbzero[26] = { 1,1,1,1,1,1,1,2,1,1,2,2,1,1,4,1,4,6,6,6,6,4,
00060 4,4,4,1 };
00061 static integer kamagn[26] = { 1,1,1,1,1,1,1,1,2,3,2,3,2,3,1,1,1,1,1,1,1,2,
00062 3,3,2,1 };
00063
00064
00065 static char fmt_9999[] = "(\002 DDRVGG: \002,a,\002 returned INFO=\002,i"
00066 "6,\002.\002,/9x,\002N=\002,i6,\002, JTYPE=\002,i6,\002, ISEED="
00067 "(\002,3(i5,\002,\002),i5,\002)\002)";
00068 static char fmt_9997[] = "(\002 DDRVGG: DGET53 returned INFO=\002,i1,"
00069 "\002 for eigenvalue \002,i6,\002.\002,/9x,\002N=\002,i6,\002, JT"
00070 "YPE=\002,i6,\002, ISEED=(\002,3(i5,\002,\002),i5,\002)\002)";
00071 static char fmt_9996[] = "(\002 DDRVGG: S not in Schur form at eigenvalu"
00072 "e \002,i6,\002.\002,/9x,\002N=\002,i6,\002, JTYPE=\002,i6,\002, "
00073 "ISEED=(\002,3(i5,\002,\002),i5,\002)\002)";
00074 static char fmt_9998[] = "(\002 DDRVGG: \002,a,\002 Eigenvectors from"
00075 " \002,a,\002 incorrectly \002,\002normalized.\002,/\002 Bits of "
00076 "error=\002,0p,g10.3,\002,\002,9x,\002N=\002,i6,\002, JTYPE=\002,"
00077 "i6,\002, ISEED=(\002,3(i5,\002,\002),i5,\002)\002)";
00078 static char fmt_9995[] = "(/1x,a3,\002 -- Real Generalized eigenvalue pr"
00079 "oblem driver\002)";
00080 static char fmt_9994[] = "(\002 Matrix types (see DDRVGG for details):"
00081 " \002)";
00082 static char fmt_9993[] = "(\002 Special Matrices:\002,23x,\002(J'=transp"
00083 "osed Jordan block)\002,/\002 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I"
00084 ") 5=(J',J') \002,\0026=(diag(J',I), diag(I,J'))\002,/\002 Diag"
00085 "onal Matrices: ( \002,\002D=diag(0,1,2,...) )\002,/\002 7=(D,"
00086 "I) 9=(large*D, small*I\002,\002) 11=(large*I, small*D) 13=(l"
00087 "arge*D, large*I)\002,/\002 8=(I,D) 10=(small*D, large*I) 12="
00088 "(small*I, large*D) \002,\002 14=(small*D, small*I)\002,/\002 15"
00089 "=(D, reversed D)\002)";
00090 static char fmt_9992[] = "(\002 Matrices Rotated by Random \002,a,\002 M"
00091 "atrices U, V:\002,/\002 16=Transposed Jordan Blocks "
00092 " 19=geometric \002,\002alpha, beta=0,1\002,/\002 17=arithm. alp"
00093 "ha&beta \002,\002 20=arithmetic alpha, beta=0,"
00094 "1\002,/\002 18=clustered \002,\002alpha, beta=0,1 21"
00095 "=random alpha, beta=0,1\002,/\002 Large & Small Matrices:\002,"
00096 "/\002 22=(large, small) \002,\00223=(small,large) 24=(smal"
00097 "l,small) 25=(large,large)\002,/\002 26=random O(1) matrices"
00098 ".\002)";
00099 static char fmt_9991[] = "(/\002 Tests performed: (S is Schur, T is tri"
00100 "angular, \002,\002Q and Z are \002,a,\002,\002,/20x,\002l and r "
00101 "are the appropriate left and right\002,/19x,\002eigenvectors, re"
00102 "sp., a is alpha, b is beta, and\002,/19x,a,\002 means \002,a,"
00103 "\002.)\002,/\002 1 = | A - Q S Z\002,a,\002 | / ( |A| n ulp ) "
00104 " 2 = | B - Q T Z\002,a,\002 | / ( |B| n ulp )\002,/\002 3 = | "
00105 "I - QQ\002,a,\002 | / ( n ulp ) 4 = | I - ZZ\002,a"
00106 ",\002 | / ( n ulp )\002,/\002 5 = difference between (alpha,beta"
00107 ") and diagonals of\002,\002 (S,T)\002,/\002 6 = max | ( b A - a "
00108 "B )\002,a,\002 l | / const. 7 = max | ( b A - a B ) r | / cons"
00109 "t.\002,/1x)";
00110 static char fmt_9990[] = "(\002 Matrix order=\002,i5,\002, type=\002,i2"
00111 ",\002, seed=\002,4(i4,\002,\002),\002 result \002,i3,\002 is\002"
00112 ",0p,f8.2)";
00113 static char fmt_9989[] = "(\002 Matrix order=\002,i5,\002, type=\002,i2"
00114 ",\002, seed=\002,4(i4,\002,\002),\002 result \002,i3,\002 is\002"
00115 ",1p,d10.3)";
00116
00117
00118 integer a_dim1, a_offset, b_dim1, b_offset, q_dim1, q_offset, s_dim1,
00119 s_offset, s2_dim1, s2_offset, t_dim1, t_offset, t2_dim1,
00120 t2_offset, vl_dim1, vl_offset, vr_dim1, vr_offset, z_dim1,
00121 z_offset, i__1, i__2, i__3, i__4;
00122 doublereal d__1, d__2, d__3, d__4, d__5, d__6, d__7, d__8, d__9, d__10;
00123
00124
00125 double d_sign(doublereal *, doublereal *);
00126 integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00127
00128
00129 integer j, n, i1, n1, jc, nb, in, jr, ns, nbz;
00130 doublereal ulp;
00131 integer iadd, nmax;
00132 doublereal temp1, temp2;
00133 logical badnn;
00134 extern int dgegs_(char *, char *, integer *, doublereal *
00135 , integer *, doublereal *, integer *, doublereal *, doublereal *,
00136 doublereal *, doublereal *, integer *, doublereal *, integer *,
00137 doublereal *, integer *, integer *), dget51_(
00138 integer *, integer *, doublereal *, integer *, doublereal *,
00139 integer *, doublereal *, integer *, doublereal *, integer *,
00140 doublereal *, doublereal *), dgegv_(char *, char *, integer *,
00141 doublereal *, integer *, doublereal *, integer *, doublereal *,
00142 doublereal *, doublereal *, doublereal *, integer *, doublereal *,
00143 integer *, doublereal *, integer *, integer *),
00144 dget52_(logical *, integer *, doublereal *, integer *, doublereal
00145 *, integer *, doublereal *, integer *, doublereal *, doublereal *,
00146 doublereal *, doublereal *, doublereal *), dget53_(doublereal *,
00147 integer *, doublereal *, integer *, doublereal *, doublereal *,
00148 doublereal *, doublereal *, integer *);
00149 doublereal dumma[4];
00150 integer iinfo;
00151 doublereal rmagn[4];
00152 integer nmats, jsize, nerrs, jtype, ntest;
00153 extern int dlatm4_(integer *, integer *, integer *,
00154 integer *, integer *, doublereal *, doublereal *, doublereal *,
00155 integer *, integer *, doublereal *, integer *), dorm2r_(char *,
00156 char *, integer *, integer *, integer *, doublereal *, integer *,
00157 doublereal *, doublereal *, integer *, doublereal *, integer *), dlabad_(doublereal *, doublereal *);
00158 logical ilabad;
00159 extern doublereal dlamch_(char *);
00160 extern int dlarfg_(integer *, doublereal *, doublereal *,
00161 integer *, doublereal *);
00162 extern doublereal dlarnd_(integer *, integer *);
00163 extern int dlacpy_(char *, integer *, integer *,
00164 doublereal *, integer *, doublereal *, integer *);
00165 doublereal safmin, safmax;
00166 integer ioldsd[4];
00167 extern integer ilaenv_(integer *, char *, char *, integer *, integer *,
00168 integer *, integer *);
00169 extern int alasvm_(char *, integer *, integer *, integer
00170 *, integer *), dlaset_(char *, integer *, integer *,
00171 doublereal *, doublereal *, doublereal *, integer *),
00172 xerbla_(char *, integer *);
00173 doublereal ulpinv;
00174 integer lwkopt, mtypes, ntestt;
00175
00176
00177 static cilist io___42 = { 0, 0, 0, fmt_9999, 0 };
00178 static cilist io___43 = { 0, 0, 0, fmt_9999, 0 };
00179 static cilist io___47 = { 0, 0, 0, fmt_9997, 0 };
00180 static cilist io___48 = { 0, 0, 0, fmt_9996, 0 };
00181 static cilist io___49 = { 0, 0, 0, fmt_9999, 0 };
00182 static cilist io___51 = { 0, 0, 0, fmt_9998, 0 };
00183 static cilist io___52 = { 0, 0, 0, fmt_9998, 0 };
00184 static cilist io___53 = { 0, 0, 0, fmt_9996, 0 };
00185 static cilist io___54 = { 0, 0, 0, fmt_9995, 0 };
00186 static cilist io___55 = { 0, 0, 0, fmt_9994, 0 };
00187 static cilist io___56 = { 0, 0, 0, fmt_9993, 0 };
00188 static cilist io___57 = { 0, 0, 0, fmt_9992, 0 };
00189 static cilist io___58 = { 0, 0, 0, fmt_9991, 0 };
00190 static cilist io___59 = { 0, 0, 0, fmt_9990, 0 };
00191 static cilist io___60 = { 0, 0, 0, fmt_9989, 0 };
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 --nn;
00529 --dotype;
00530 --iseed;
00531 t2_dim1 = *lda;
00532 t2_offset = 1 + t2_dim1;
00533 t2 -= t2_offset;
00534 s2_dim1 = *lda;
00535 s2_offset = 1 + s2_dim1;
00536 s2 -= s2_offset;
00537 t_dim1 = *lda;
00538 t_offset = 1 + t_dim1;
00539 t -= t_offset;
00540 s_dim1 = *lda;
00541 s_offset = 1 + s_dim1;
00542 s -= s_offset;
00543 b_dim1 = *lda;
00544 b_offset = 1 + b_dim1;
00545 b -= b_offset;
00546 a_dim1 = *lda;
00547 a_offset = 1 + a_dim1;
00548 a -= a_offset;
00549 vr_dim1 = *ldq;
00550 vr_offset = 1 + vr_dim1;
00551 vr -= vr_offset;
00552 vl_dim1 = *ldq;
00553 vl_offset = 1 + vl_dim1;
00554 vl -= vl_offset;
00555 z_dim1 = *ldq;
00556 z_offset = 1 + z_dim1;
00557 z__ -= z_offset;
00558 q_dim1 = *ldq;
00559 q_offset = 1 + q_dim1;
00560 q -= q_offset;
00561 --alphr1;
00562 --alphi1;
00563 --beta1;
00564 --alphr2;
00565 --alphi2;
00566 --beta2;
00567 --work;
00568 --result;
00569
00570
00571
00572
00573
00574
00575
00576 *info = 0;
00577
00578 badnn = FALSE_;
00579 nmax = 1;
00580 i__1 = *nsizes;
00581 for (j = 1; j <= i__1; ++j) {
00582
00583 i__2 = nmax, i__3 = nn[j];
00584 nmax = max(i__2,i__3);
00585 if (nn[j] < 0) {
00586 badnn = TRUE_;
00587 }
00588
00589 }
00590
00591
00592
00593
00594
00595 i__1 = 1, i__2 = ilaenv_(&c__1, "DGEQRF", " ", &nmax, &nmax, &c_n1, &c_n1), i__1 = max(i__1,i__2), i__2 = ilaenv_(&
00596 c__1, "DORMQR", "LT", &nmax, &nmax, &nmax, &c_n1), i__1 = max(i__1,i__2), i__2 = ilaenv_(&c__1, "DORGQR",
00597 " ", &nmax, &nmax, &nmax, &c_n1);
00598 nb = max(i__1,i__2);
00599 nbz = ilaenv_(&c__1, "DHGEQZ", "SII", &nmax, &c__1, &nmax, &c__0);
00600 ns = ilaenv_(&c__4, "DHGEQZ", "SII", &nmax, &c__1, &nmax, &c__0);
00601 i1 = nbz + ns;
00602
00603 i__1 = nmax * 6, i__2 = nmax * (nb + 1), i__1 = max(i__1,i__2), i__2 = ((
00604 i1 << 1) + nmax + 1) * (i1 + 1);
00605 lwkopt = (nmax << 1) + max(i__1,i__2);
00606
00607
00608
00609 if (*nsizes < 0) {
00610 *info = -1;
00611 } else if (badnn) {
00612 *info = -2;
00613 } else if (*ntypes < 0) {
00614 *info = -3;
00615 } else if (*thresh < 0.) {
00616 *info = -6;
00617 } else if (*lda <= 1 || *lda < nmax) {
00618 *info = -10;
00619 } else if (*ldq <= 1 || *ldq < nmax) {
00620 *info = -19;
00621 } else if (lwkopt > *lwork) {
00622 *info = -30;
00623 }
00624
00625 if (*info != 0) {
00626 i__1 = -(*info);
00627 xerbla_("DDRVGG", &i__1);
00628 return 0;
00629 }
00630
00631
00632
00633 if (*nsizes == 0 || *ntypes == 0) {
00634 return 0;
00635 }
00636
00637 safmin = dlamch_("Safe minimum");
00638 ulp = dlamch_("Epsilon") * dlamch_("Base");
00639 safmin /= ulp;
00640 safmax = 1. / safmin;
00641 dlabad_(&safmin, &safmax);
00642 ulpinv = 1. / ulp;
00643
00644
00645
00646 rmagn[0] = 0.;
00647 rmagn[1] = 1.;
00648
00649
00650
00651 ntestt = 0;
00652 nerrs = 0;
00653 nmats = 0;
00654
00655 i__1 = *nsizes;
00656 for (jsize = 1; jsize <= i__1; ++jsize) {
00657 n = nn[jsize];
00658 n1 = max(1,n);
00659 rmagn[2] = safmax * ulp / (doublereal) n1;
00660 rmagn[3] = safmin * ulpinv * n1;
00661
00662 if (*nsizes != 1) {
00663 mtypes = min(26,*ntypes);
00664 } else {
00665 mtypes = min(27,*ntypes);
00666 }
00667
00668 i__2 = mtypes;
00669 for (jtype = 1; jtype <= i__2; ++jtype) {
00670 if (! dotype[jtype]) {
00671 goto L160;
00672 }
00673 ++nmats;
00674 ntest = 0;
00675
00676
00677
00678 for (j = 1; j <= 4; ++j) {
00679 ioldsd[j - 1] = iseed[j];
00680
00681 }
00682
00683
00684
00685 for (j = 1; j <= 15; ++j) {
00686 result[j] = 0.;
00687
00688 }
00689
00690
00691
00692
00693
00694
00695
00696
00697
00698
00699
00700
00701
00702
00703
00704
00705
00706
00707
00708
00709
00710
00711
00712
00713 if (mtypes > 26) {
00714 goto L110;
00715 }
00716 iinfo = 0;
00717 if (kclass[jtype - 1] < 3) {
00718
00719
00720
00721 if ((i__3 = katype[jtype - 1], abs(i__3)) == 3) {
00722 in = ((n - 1) / 2 << 1) + 1;
00723 if (in != n) {
00724 dlaset_("Full", &n, &n, &c_b36, &c_b36, &a[a_offset],
00725 lda);
00726 }
00727 } else {
00728 in = n;
00729 }
00730 dlatm4_(&katype[jtype - 1], &in, &kz1[kazero[jtype - 1] - 1],
00731 &kz2[kazero[jtype - 1] - 1], &iasign[jtype - 1], &
00732 rmagn[kamagn[jtype - 1]], &ulp, &rmagn[ktrian[jtype -
00733 1] * kamagn[jtype - 1]], &c__2, &iseed[1], &a[
00734 a_offset], lda);
00735 iadd = kadd[kazero[jtype - 1] - 1];
00736 if (iadd > 0 && iadd <= n) {
00737 a[iadd + iadd * a_dim1] = 1.;
00738 }
00739
00740
00741
00742 if ((i__3 = kbtype[jtype - 1], abs(i__3)) == 3) {
00743 in = ((n - 1) / 2 << 1) + 1;
00744 if (in != n) {
00745 dlaset_("Full", &n, &n, &c_b36, &c_b36, &b[b_offset],
00746 lda);
00747 }
00748 } else {
00749 in = n;
00750 }
00751 dlatm4_(&kbtype[jtype - 1], &in, &kz1[kbzero[jtype - 1] - 1],
00752 &kz2[kbzero[jtype - 1] - 1], &ibsign[jtype - 1], &
00753 rmagn[kbmagn[jtype - 1]], &c_b42, &rmagn[ktrian[jtype
00754 - 1] * kbmagn[jtype - 1]], &c__2, &iseed[1], &b[
00755 b_offset], lda);
00756 iadd = kadd[kbzero[jtype - 1] - 1];
00757 if (iadd != 0 && iadd <= n) {
00758 b[iadd + iadd * b_dim1] = 1.;
00759 }
00760
00761 if (kclass[jtype - 1] == 2 && n > 0) {
00762
00763
00764
00765
00766
00767
00768 i__3 = n - 1;
00769 for (jc = 1; jc <= i__3; ++jc) {
00770 i__4 = n;
00771 for (jr = jc; jr <= i__4; ++jr) {
00772 q[jr + jc * q_dim1] = dlarnd_(&c__3, &iseed[1]);
00773 z__[jr + jc * z_dim1] = dlarnd_(&c__3, &iseed[1]);
00774
00775 }
00776 i__4 = n + 1 - jc;
00777 dlarfg_(&i__4, &q[jc + jc * q_dim1], &q[jc + 1 + jc *
00778 q_dim1], &c__1, &work[jc]);
00779 work[(n << 1) + jc] = d_sign(&c_b42, &q[jc + jc *
00780 q_dim1]);
00781 q[jc + jc * q_dim1] = 1.;
00782 i__4 = n + 1 - jc;
00783 dlarfg_(&i__4, &z__[jc + jc * z_dim1], &z__[jc + 1 +
00784 jc * z_dim1], &c__1, &work[n + jc]);
00785 work[n * 3 + jc] = d_sign(&c_b42, &z__[jc + jc *
00786 z_dim1]);
00787 z__[jc + jc * z_dim1] = 1.;
00788
00789 }
00790 q[n + n * q_dim1] = 1.;
00791 work[n] = 0.;
00792 d__1 = dlarnd_(&c__2, &iseed[1]);
00793 work[n * 3] = d_sign(&c_b42, &d__1);
00794 z__[n + n * z_dim1] = 1.;
00795 work[n * 2] = 0.;
00796 d__1 = dlarnd_(&c__2, &iseed[1]);
00797 work[n * 4] = d_sign(&c_b42, &d__1);
00798
00799
00800
00801 i__3 = n;
00802 for (jc = 1; jc <= i__3; ++jc) {
00803 i__4 = n;
00804 for (jr = 1; jr <= i__4; ++jr) {
00805 a[jr + jc * a_dim1] = work[(n << 1) + jr] * work[
00806 n * 3 + jc] * a[jr + jc * a_dim1];
00807 b[jr + jc * b_dim1] = work[(n << 1) + jr] * work[
00808 n * 3 + jc] * b[jr + jc * b_dim1];
00809
00810 }
00811
00812 }
00813 i__3 = n - 1;
00814 dorm2r_("L", "N", &n, &n, &i__3, &q[q_offset], ldq, &work[
00815 1], &a[a_offset], lda, &work[(n << 1) + 1], &
00816 iinfo);
00817 if (iinfo != 0) {
00818 goto L100;
00819 }
00820 i__3 = n - 1;
00821 dorm2r_("R", "T", &n, &n, &i__3, &z__[z_offset], ldq, &
00822 work[n + 1], &a[a_offset], lda, &work[(n << 1) +
00823 1], &iinfo);
00824 if (iinfo != 0) {
00825 goto L100;
00826 }
00827 i__3 = n - 1;
00828 dorm2r_("L", "N", &n, &n, &i__3, &q[q_offset], ldq, &work[
00829 1], &b[b_offset], lda, &work[(n << 1) + 1], &
00830 iinfo);
00831 if (iinfo != 0) {
00832 goto L100;
00833 }
00834 i__3 = n - 1;
00835 dorm2r_("R", "T", &n, &n, &i__3, &z__[z_offset], ldq, &
00836 work[n + 1], &b[b_offset], lda, &work[(n << 1) +
00837 1], &iinfo);
00838 if (iinfo != 0) {
00839 goto L100;
00840 }
00841 }
00842 } else {
00843
00844
00845
00846 i__3 = n;
00847 for (jc = 1; jc <= i__3; ++jc) {
00848 i__4 = n;
00849 for (jr = 1; jr <= i__4; ++jr) {
00850 a[jr + jc * a_dim1] = rmagn[kamagn[jtype - 1]] *
00851 dlarnd_(&c__2, &iseed[1]);
00852 b[jr + jc * b_dim1] = rmagn[kbmagn[jtype - 1]] *
00853 dlarnd_(&c__2, &iseed[1]);
00854
00855 }
00856
00857 }
00858 }
00859
00860 L100:
00861
00862 if (iinfo != 0) {
00863 io___42.ciunit = *nounit;
00864 s_wsfe(&io___42);
00865 do_fio(&c__1, "Generator", (ftnlen)9);
00866 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00867 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00868 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00869 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00870 e_wsfe();
00871 *info = abs(iinfo);
00872 return 0;
00873 }
00874
00875 L110:
00876
00877
00878
00879 dlacpy_(" ", &n, &n, &a[a_offset], lda, &s[s_offset], lda);
00880 dlacpy_(" ", &n, &n, &b[b_offset], lda, &t[t_offset], lda);
00881 ntest = 1;
00882 result[1] = ulpinv;
00883
00884 dgegs_("V", "V", &n, &s[s_offset], lda, &t[t_offset], lda, &
00885 alphr1[1], &alphi1[1], &beta1[1], &q[q_offset], ldq, &z__[
00886 z_offset], ldq, &work[1], lwork, &iinfo);
00887 if (iinfo != 0) {
00888 io___43.ciunit = *nounit;
00889 s_wsfe(&io___43);
00890 do_fio(&c__1, "DGEGS", (ftnlen)5);
00891 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00892 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00893 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00894 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00895 e_wsfe();
00896 *info = abs(iinfo);
00897 goto L140;
00898 }
00899
00900 ntest = 4;
00901
00902
00903
00904 dget51_(&c__1, &n, &a[a_offset], lda, &s[s_offset], lda, &q[
00905 q_offset], ldq, &z__[z_offset], ldq, &work[1], &result[1])
00906 ;
00907 dget51_(&c__1, &n, &b[b_offset], lda, &t[t_offset], lda, &q[
00908 q_offset], ldq, &z__[z_offset], ldq, &work[1], &result[2])
00909 ;
00910 dget51_(&c__3, &n, &b[b_offset], lda, &t[t_offset], lda, &q[
00911 q_offset], ldq, &q[q_offset], ldq, &work[1], &result[3]);
00912 dget51_(&c__3, &n, &b[b_offset], lda, &t[t_offset], lda, &z__[
00913 z_offset], ldq, &z__[z_offset], ldq, &work[1], &result[4])
00914 ;
00915
00916
00917
00918
00919 temp1 = 0.;
00920
00921 i__3 = n;
00922 for (j = 1; j <= i__3; ++j) {
00923 ilabad = FALSE_;
00924 if (alphi1[j] == 0.) {
00925
00926 d__7 = safmin, d__8 = (d__2 = alphr1[j], abs(d__2)), d__7
00927 = max(d__7,d__8), d__8 = (d__3 = s[j + j * s_dim1]
00928 , abs(d__3));
00929
00930 d__9 = safmin, d__10 = (d__5 = beta1[j], abs(d__5)), d__9
00931 = max(d__9,d__10), d__10 = (d__6 = t[j + j *
00932 t_dim1], abs(d__6));
00933 temp2 = ((d__1 = alphr1[j] - s[j + j * s_dim1], abs(d__1))
00934 / max(d__7,d__8) + (d__4 = beta1[j] - t[j + j *
00935 t_dim1], abs(d__4)) / max(d__9,d__10)) / ulp;
00936 if (j < n) {
00937 if (s[j + 1 + j * s_dim1] != 0.) {
00938 ilabad = TRUE_;
00939 }
00940 }
00941 if (j > 1) {
00942 if (s[j + (j - 1) * s_dim1] != 0.) {
00943 ilabad = TRUE_;
00944 }
00945 }
00946 } else {
00947 if (alphi1[j] > 0.) {
00948 i1 = j;
00949 } else {
00950 i1 = j - 1;
00951 }
00952 if (i1 <= 0 || i1 >= n) {
00953 ilabad = TRUE_;
00954 } else if (i1 < n - 1) {
00955 if (s[i1 + 2 + (i1 + 1) * s_dim1] != 0.) {
00956 ilabad = TRUE_;
00957 }
00958 } else if (i1 > 1) {
00959 if (s[i1 + (i1 - 1) * s_dim1] != 0.) {
00960 ilabad = TRUE_;
00961 }
00962 }
00963 if (! ilabad) {
00964 dget53_(&s[i1 + i1 * s_dim1], lda, &t[i1 + i1 *
00965 t_dim1], lda, &beta1[j], &alphr1[j], &alphi1[
00966 j], &temp2, &iinfo);
00967 if (iinfo >= 3) {
00968 io___47.ciunit = *nounit;
00969 s_wsfe(&io___47);
00970 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(
00971 integer));
00972 do_fio(&c__1, (char *)&j, (ftnlen)sizeof(integer))
00973 ;
00974 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer))
00975 ;
00976 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(
00977 integer));
00978 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(
00979 integer));
00980 e_wsfe();
00981 *info = abs(iinfo);
00982 }
00983 } else {
00984 temp2 = ulpinv;
00985 }
00986 }
00987 temp1 = max(temp1,temp2);
00988 if (ilabad) {
00989 io___48.ciunit = *nounit;
00990 s_wsfe(&io___48);
00991 do_fio(&c__1, (char *)&j, (ftnlen)sizeof(integer));
00992 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00993 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00994 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
00995 ;
00996 e_wsfe();
00997 }
00998
00999 }
01000 result[5] = temp1;
01001
01002
01003
01004
01005
01006 dlacpy_(" ", &n, &n, &a[a_offset], lda, &s2[s2_offset], lda);
01007 dlacpy_(" ", &n, &n, &b[b_offset], lda, &t2[t2_offset], lda);
01008 ntest = 6;
01009 result[6] = ulpinv;
01010
01011 dgegv_("V", "V", &n, &s2[s2_offset], lda, &t2[t2_offset], lda, &
01012 alphr2[1], &alphi2[1], &beta2[1], &vl[vl_offset], ldq, &
01013 vr[vr_offset], ldq, &work[1], lwork, &iinfo);
01014 if (iinfo != 0) {
01015 io___49.ciunit = *nounit;
01016 s_wsfe(&io___49);
01017 do_fio(&c__1, "DGEGV", (ftnlen)5);
01018 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
01019 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01020 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01021 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01022 e_wsfe();
01023 *info = abs(iinfo);
01024 goto L140;
01025 }
01026
01027 ntest = 7;
01028
01029
01030
01031 dget52_(&c_true, &n, &a[a_offset], lda, &b[b_offset], lda, &vl[
01032 vl_offset], ldq, &alphr2[1], &alphi2[1], &beta2[1], &work[
01033 1], dumma);
01034 result[6] = dumma[0];
01035 if (dumma[1] > *thrshn) {
01036 io___51.ciunit = *nounit;
01037 s_wsfe(&io___51);
01038 do_fio(&c__1, "Left", (ftnlen)4);
01039 do_fio(&c__1, "DGEGV", (ftnlen)5);
01040 do_fio(&c__1, (char *)&dumma[1], (ftnlen)sizeof(doublereal));
01041 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01042 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01043 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01044 e_wsfe();
01045 }
01046
01047 dget52_(&c_false, &n, &a[a_offset], lda, &b[b_offset], lda, &vr[
01048 vr_offset], ldq, &alphr2[1], &alphi2[1], &beta2[1], &work[
01049 1], dumma);
01050 result[7] = dumma[0];
01051 if (dumma[1] > *thresh) {
01052 io___52.ciunit = *nounit;
01053 s_wsfe(&io___52);
01054 do_fio(&c__1, "Right", (ftnlen)5);
01055 do_fio(&c__1, "DGEGV", (ftnlen)5);
01056 do_fio(&c__1, (char *)&dumma[1], (ftnlen)sizeof(doublereal));
01057 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01058 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01059 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01060 e_wsfe();
01061 }
01062
01063
01064
01065 i__3 = n;
01066 for (j = 1; j <= i__3; ++j) {
01067 ilabad = FALSE_;
01068 if (alphi2[j] > 0.) {
01069 if (j == n) {
01070 ilabad = TRUE_;
01071 } else if (alphi2[j + 1] >= 0.) {
01072 ilabad = TRUE_;
01073 }
01074 } else if (alphi2[j] < 0.) {
01075 if (j == 1) {
01076 ilabad = TRUE_;
01077 } else if (alphi2[j - 1] <= 0.) {
01078 ilabad = TRUE_;
01079 }
01080 }
01081 if (ilabad) {
01082 io___53.ciunit = *nounit;
01083 s_wsfe(&io___53);
01084 do_fio(&c__1, (char *)&j, (ftnlen)sizeof(integer));
01085 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01086 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01087 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
01088 ;
01089 e_wsfe();
01090 }
01091
01092 }
01093
01094
01095
01096 L140:
01097
01098 ntestt += ntest;
01099
01100
01101
01102 i__3 = ntest;
01103 for (jr = 1; jr <= i__3; ++jr) {
01104 if (result[jr] >= *thresh) {
01105
01106
01107
01108
01109 if (nerrs == 0) {
01110 io___54.ciunit = *nounit;
01111 s_wsfe(&io___54);
01112 do_fio(&c__1, "DGG", (ftnlen)3);
01113 e_wsfe();
01114
01115
01116
01117 io___55.ciunit = *nounit;
01118 s_wsfe(&io___55);
01119 e_wsfe();
01120 io___56.ciunit = *nounit;
01121 s_wsfe(&io___56);
01122 e_wsfe();
01123 io___57.ciunit = *nounit;
01124 s_wsfe(&io___57);
01125 do_fio(&c__1, "Orthogonal", (ftnlen)10);
01126 e_wsfe();
01127
01128
01129
01130 io___58.ciunit = *nounit;
01131 s_wsfe(&io___58);
01132 do_fio(&c__1, "orthogonal", (ftnlen)10);
01133 do_fio(&c__1, "'", (ftnlen)1);
01134 do_fio(&c__1, "transpose", (ftnlen)9);
01135 for (j = 1; j <= 5; ++j) {
01136 do_fio(&c__1, "'", (ftnlen)1);
01137 }
01138 e_wsfe();
01139
01140 }
01141 ++nerrs;
01142 if (result[jr] < 1e4) {
01143 io___59.ciunit = *nounit;
01144 s_wsfe(&io___59);
01145 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01146 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer))
01147 ;
01148 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(
01149 integer));
01150 do_fio(&c__1, (char *)&jr, (ftnlen)sizeof(integer));
01151 do_fio(&c__1, (char *)&result[jr], (ftnlen)sizeof(
01152 doublereal));
01153 e_wsfe();
01154 } else {
01155 io___60.ciunit = *nounit;
01156 s_wsfe(&io___60);
01157 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01158 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer))
01159 ;
01160 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(
01161 integer));
01162 do_fio(&c__1, (char *)&jr, (ftnlen)sizeof(integer));
01163 do_fio(&c__1, (char *)&result[jr], (ftnlen)sizeof(
01164 doublereal));
01165 e_wsfe();
01166 }
01167 }
01168
01169 }
01170
01171 L160:
01172 ;
01173 }
01174
01175 }
01176
01177
01178
01179 alasvm_("DGG", nounit, &nerrs, &ntestt, &c__0);
01180 return 0;
01181
01182
01183
01184
01185
01186
01187
01188
01189
01190
01191
01192 }