00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015
00016
00017
00018 struct {
00019 integer infot, nunit;
00020 logical ok, lerr;
00021 } infoc_;
00022
00023 #define infoc_1 infoc_
00024
00025 struct {
00026 char srnamt[32];
00027 } srnamc_;
00028
00029 #define srnamc_1 srnamc_
00030
00031
00032
00033 static integer c__0 = 0;
00034 static integer c_n1 = -1;
00035 static integer c__1 = 1;
00036 static doublecomplex c_b50 = {0.,0.};
00037 static doublecomplex c_b51 = {1.,0.};
00038 static integer c__7 = 7;
00039
00040 int zchkpb_(logical *dotype, integer *nn, integer *nval,
00041 integer *nnb, integer *nbval, integer *nns, integer *nsval,
00042 doublereal *thresh, logical *tsterr, integer *nmax, doublecomplex *a,
00043 doublecomplex *afac, doublecomplex *ainv, doublecomplex *b,
00044 doublecomplex *x, doublecomplex *xact, doublecomplex *work,
00045 doublereal *rwork, integer *nout)
00046 {
00047
00048
00049 static integer iseedy[4] = { 1988,1989,1990,1991 };
00050
00051
00052 static char fmt_9999[] = "(\002 UPLO='\002,a1,\002', N=\002,i5,\002, KD"
00053 "=\002,i5,\002, NB=\002,i4,\002, type \002,i2,\002, test \002,i2"
00054 ",\002, ratio= \002,g12.5)";
00055 static char fmt_9998[] = "(\002 UPLO='\002,a1,\002', N=\002,i5,\002, KD"
00056 "=\002,i5,\002, NRHS=\002,i3,\002, type \002,i2,\002, test(\002,i"
00057 "2,\002) = \002,g12.5)";
00058 static char fmt_9997[] = "(\002 UPLO='\002,a1,\002', N=\002,i5,\002, KD"
00059 "=\002,i5,\002,\002,10x,\002 type \002,i2,\002, test(\002,i2,\002"
00060 ") = \002,g12.5)";
00061
00062
00063 integer i__1, i__2, i__3, i__4, i__5, i__6;
00064
00065
00066 int s_copy(char *, char *, ftnlen, ftnlen);
00067 integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00068
00069
00070 integer i__, k, n, i1, i2, kd, nb, in, kl, iw, ku, lda, ikd, inb, nkd,
00071 ldab, ioff, mode, koff, imat, info;
00072 char path[3], dist[1];
00073 integer irhs, nrhs;
00074 char uplo[1], type__[1];
00075 integer nrun;
00076 extern int alahd_(integer *, char *);
00077 integer nfail, iseed[4];
00078 extern doublereal dget06_(doublereal *, doublereal *);
00079 integer kdval[4];
00080 doublereal rcond;
00081 integer nimat;
00082 doublereal anorm;
00083 extern int zget04_(integer *, integer *, doublecomplex *,
00084 integer *, doublecomplex *, integer *, doublereal *, doublereal *
00085 ), zpbt01_(char *, integer *, integer *, doublecomplex *, integer
00086 *, doublecomplex *, integer *, doublereal *, doublereal *)
00087 , zpbt02_(char *, integer *, integer *, integer *, doublecomplex *
00088 , integer *, doublecomplex *, integer *, doublecomplex *, integer
00089 *, doublereal *, doublereal *), zpbt05_(char *, integer *,
00090 integer *, integer *, doublecomplex *, integer *, doublecomplex *
00091 , integer *, doublecomplex *, integer *, doublecomplex *, integer
00092 *, doublereal *, doublereal *, doublereal *);
00093 integer iuplo, izero, nerrs;
00094 logical zerot;
00095 extern int zcopy_(integer *, doublecomplex *, integer *,
00096 doublecomplex *, integer *), zswap_(integer *, doublecomplex *,
00097 integer *, doublecomplex *, integer *);
00098 char xtype[1];
00099 extern int zlatb4_(char *, integer *, integer *, integer
00100 *, char *, integer *, integer *, doublereal *, integer *,
00101 doublereal *, char *), alaerh_(char *,
00102 char *, integer *, integer *, char *, integer *, integer *,
00103 integer *, integer *, integer *, integer *, integer *, integer *,
00104 integer *);
00105 doublereal rcondc;
00106 char packit[1];
00107 extern doublereal zlanhb_(char *, char *, integer *, integer *,
00108 doublecomplex *, integer *, doublereal *),
00109 zlange_(char *, integer *, integer *, doublecomplex *, integer *,
00110 doublereal *);
00111 extern int alasum_(char *, integer *, integer *, integer
00112 *, integer *);
00113 doublereal cndnum;
00114 extern int zlaipd_(integer *, doublecomplex *, integer *,
00115 integer *);
00116 doublereal ainvnm;
00117 extern int zpbcon_(char *, integer *, integer *,
00118 doublecomplex *, integer *, doublereal *, doublereal *,
00119 doublecomplex *, doublereal *, integer *), xlaenv_(
00120 integer *, integer *), zlacpy_(char *, integer *, integer *,
00121 doublecomplex *, integer *, doublecomplex *, integer *),
00122 zlarhs_(char *, char *, char *, char *, integer *, integer *,
00123 integer *, integer *, integer *, doublecomplex *, integer *,
00124 doublecomplex *, integer *, doublecomplex *, integer *, integer *,
00125 integer *), zlaset_(char *,
00126 integer *, integer *, doublecomplex *, doublecomplex *,
00127 doublecomplex *, integer *), zpbrfs_(char *, integer *,
00128 integer *, integer *, doublecomplex *, integer *, doublecomplex *,
00129 integer *, doublecomplex *, integer *, doublecomplex *, integer *
00130 , doublereal *, doublereal *, doublecomplex *, doublereal *,
00131 integer *), zpbtrf_(char *, integer *, integer *,
00132 doublecomplex *, integer *, integer *), zlatms_(integer *,
00133 integer *, char *, integer *, char *, doublereal *, integer *,
00134 doublereal *, doublereal *, integer *, integer *, char *,
00135 doublecomplex *, integer *, doublecomplex *, integer *);
00136 doublereal result[7];
00137 extern int zerrpo_(char *, integer *), zpbtrs_(
00138 char *, integer *, integer *, integer *, doublecomplex *, integer
00139 *, doublecomplex *, integer *, integer *);
00140
00141
00142 static cilist io___40 = { 0, 0, 0, fmt_9999, 0 };
00143 static cilist io___46 = { 0, 0, 0, fmt_9998, 0 };
00144 static cilist io___48 = { 0, 0, 0, fmt_9997, 0 };
00145
00146
00147
00148
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 --rwork;
00243 --work;
00244 --xact;
00245 --x;
00246 --b;
00247 --ainv;
00248 --afac;
00249 --a;
00250 --nsval;
00251 --nbval;
00252 --nval;
00253 --dotype;
00254
00255
00256
00257
00258
00259
00260
00261 s_copy(path, "Zomplex precision", (ftnlen)1, (ftnlen)17);
00262 s_copy(path + 1, "PB", (ftnlen)2, (ftnlen)2);
00263 nrun = 0;
00264 nfail = 0;
00265 nerrs = 0;
00266 for (i__ = 1; i__ <= 4; ++i__) {
00267 iseed[i__ - 1] = iseedy[i__ - 1];
00268
00269 }
00270
00271
00272
00273 if (*tsterr) {
00274 zerrpo_(path, nout);
00275 }
00276 infoc_1.infot = 0;
00277 kdval[0] = 0;
00278
00279
00280
00281 i__1 = *nn;
00282 for (in = 1; in <= i__1; ++in) {
00283 n = nval[in];
00284 lda = max(n,1);
00285 *(unsigned char *)xtype = 'N';
00286
00287
00288
00289
00290 i__2 = 1, i__3 = min(n,4);
00291 nkd = max(i__2,i__3);
00292 nimat = 8;
00293 if (n == 0) {
00294 nimat = 1;
00295 }
00296
00297 kdval[1] = n + (n + 1) / 4;
00298 kdval[2] = (n * 3 - 1) / 4;
00299 kdval[3] = (n + 1) / 4;
00300
00301 i__2 = nkd;
00302 for (ikd = 1; ikd <= i__2; ++ikd) {
00303
00304
00305
00306
00307
00308 kd = kdval[ikd - 1];
00309 ldab = kd + 1;
00310
00311
00312
00313 for (iuplo = 1; iuplo <= 2; ++iuplo) {
00314 koff = 1;
00315 if (iuplo == 1) {
00316 *(unsigned char *)uplo = 'U';
00317
00318 i__3 = 1, i__4 = kd + 2 - n;
00319 koff = max(i__3,i__4);
00320 *(unsigned char *)packit = 'Q';
00321 } else {
00322 *(unsigned char *)uplo = 'L';
00323 *(unsigned char *)packit = 'B';
00324 }
00325
00326 i__3 = nimat;
00327 for (imat = 1; imat <= i__3; ++imat) {
00328
00329
00330
00331 if (! dotype[imat]) {
00332 goto L60;
00333 }
00334
00335
00336
00337 zerot = imat >= 2 && imat <= 4;
00338 if (zerot && n < imat - 1) {
00339 goto L60;
00340 }
00341
00342 if (! zerot || ! dotype[1]) {
00343
00344
00345
00346
00347 zlatb4_(path, &imat, &n, &n, type__, &kl, &ku, &anorm,
00348 &mode, &cndnum, dist);
00349
00350 s_copy(srnamc_1.srnamt, "ZLATMS", (ftnlen)32, (ftnlen)
00351 6);
00352 zlatms_(&n, &n, dist, iseed, type__, &rwork[1], &mode,
00353 &cndnum, &anorm, &kd, &kd, packit, &a[koff],
00354 &ldab, &work[1], &info);
00355
00356
00357
00358 if (info != 0) {
00359 alaerh_(path, "ZLATMS", &info, &c__0, uplo, &n, &
00360 n, &kd, &kd, &c_n1, &imat, &nfail, &nerrs,
00361 nout);
00362 goto L60;
00363 }
00364 } else if (izero > 0) {
00365
00366
00367
00368
00369 iw = (lda << 1) + 1;
00370 if (iuplo == 1) {
00371 ioff = (izero - 1) * ldab + kd + 1;
00372 i__4 = izero - i1;
00373 zcopy_(&i__4, &work[iw], &c__1, &a[ioff - izero +
00374 i1], &c__1);
00375 iw = iw + izero - i1;
00376 i__4 = i2 - izero + 1;
00377
00378 i__6 = ldab - 1;
00379 i__5 = max(i__6,1);
00380 zcopy_(&i__4, &work[iw], &c__1, &a[ioff], &i__5);
00381 } else {
00382 ioff = (i1 - 1) * ldab + 1;
00383 i__4 = izero - i1;
00384
00385 i__6 = ldab - 1;
00386 i__5 = max(i__6,1);
00387 zcopy_(&i__4, &work[iw], &c__1, &a[ioff + izero -
00388 i1], &i__5);
00389 ioff = (izero - 1) * ldab + 1;
00390 iw = iw + izero - i1;
00391 i__4 = i2 - izero + 1;
00392 zcopy_(&i__4, &work[iw], &c__1, &a[ioff], &c__1);
00393 }
00394 }
00395
00396
00397
00398
00399 izero = 0;
00400 if (zerot) {
00401 if (imat == 2) {
00402 izero = 1;
00403 } else if (imat == 3) {
00404 izero = n;
00405 } else {
00406 izero = n / 2 + 1;
00407 }
00408
00409
00410
00411 iw = lda << 1;
00412
00413 i__5 = (kd << 1) + 1;
00414 i__4 = min(i__5,n);
00415 for (i__ = 1; i__ <= i__4; ++i__) {
00416 i__5 = iw + i__;
00417 work[i__5].r = 0., work[i__5].i = 0.;
00418
00419 }
00420 ++iw;
00421
00422 i__4 = izero - kd;
00423 i1 = max(i__4,1);
00424
00425 i__4 = izero + kd;
00426 i2 = min(i__4,n);
00427
00428 if (iuplo == 1) {
00429 ioff = (izero - 1) * ldab + kd + 1;
00430 i__4 = izero - i1;
00431 zswap_(&i__4, &a[ioff - izero + i1], &c__1, &work[
00432 iw], &c__1);
00433 iw = iw + izero - i1;
00434 i__4 = i2 - izero + 1;
00435
00436 i__6 = ldab - 1;
00437 i__5 = max(i__6,1);
00438 zswap_(&i__4, &a[ioff], &i__5, &work[iw], &c__1);
00439 } else {
00440 ioff = (i1 - 1) * ldab + 1;
00441 i__4 = izero - i1;
00442
00443 i__6 = ldab - 1;
00444 i__5 = max(i__6,1);
00445 zswap_(&i__4, &a[ioff + izero - i1], &i__5, &work[
00446 iw], &c__1);
00447 ioff = (izero - 1) * ldab + 1;
00448 iw = iw + izero - i1;
00449 i__4 = i2 - izero + 1;
00450 zswap_(&i__4, &a[ioff], &c__1, &work[iw], &c__1);
00451 }
00452 }
00453
00454
00455
00456 if (iuplo == 1) {
00457 zlaipd_(&n, &a[kd + 1], &ldab, &c__0);
00458 } else {
00459 zlaipd_(&n, &a[1], &ldab, &c__0);
00460 }
00461
00462
00463
00464 i__4 = *nnb;
00465 for (inb = 1; inb <= i__4; ++inb) {
00466 nb = nbval[inb];
00467 xlaenv_(&c__1, &nb);
00468
00469
00470
00471
00472 i__5 = kd + 1;
00473 zlacpy_("Full", &i__5, &n, &a[1], &ldab, &afac[1], &
00474 ldab);
00475 s_copy(srnamc_1.srnamt, "ZPBTRF", (ftnlen)32, (ftnlen)
00476 6);
00477 zpbtrf_(uplo, &n, &kd, &afac[1], &ldab, &info);
00478
00479
00480
00481 if (info != izero) {
00482 alaerh_(path, "ZPBTRF", &info, &izero, uplo, &n, &
00483 n, &kd, &kd, &nb, &imat, &nfail, &nerrs,
00484 nout);
00485 goto L50;
00486 }
00487
00488
00489
00490 if (info != 0) {
00491 goto L50;
00492 }
00493
00494
00495
00496
00497
00498 i__5 = kd + 1;
00499 zlacpy_("Full", &i__5, &n, &afac[1], &ldab, &ainv[1],
00500 &ldab);
00501 zpbt01_(uplo, &n, &kd, &a[1], &ldab, &ainv[1], &ldab,
00502 &rwork[1], result);
00503
00504
00505
00506 if (result[0] >= *thresh) {
00507 if (nfail == 0 && nerrs == 0) {
00508 alahd_(nout, path);
00509 }
00510 io___40.ciunit = *nout;
00511 s_wsfe(&io___40);
00512 do_fio(&c__1, uplo, (ftnlen)1);
00513 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer))
00514 ;
00515 do_fio(&c__1, (char *)&kd, (ftnlen)sizeof(integer)
00516 );
00517 do_fio(&c__1, (char *)&nb, (ftnlen)sizeof(integer)
00518 );
00519 do_fio(&c__1, (char *)&imat, (ftnlen)sizeof(
00520 integer));
00521 do_fio(&c__1, (char *)&c__1, (ftnlen)sizeof(
00522 integer));
00523 do_fio(&c__1, (char *)&result[0], (ftnlen)sizeof(
00524 doublereal));
00525 e_wsfe();
00526 ++nfail;
00527 }
00528 ++nrun;
00529
00530
00531
00532 if (inb > 1) {
00533 goto L50;
00534 }
00535
00536
00537
00538
00539 zlaset_("Full", &n, &n, &c_b50, &c_b51, &ainv[1], &
00540 lda);
00541 s_copy(srnamc_1.srnamt, "ZPBTRS", (ftnlen)32, (ftnlen)
00542 6);
00543 zpbtrs_(uplo, &n, &kd, &n, &afac[1], &ldab, &ainv[1],
00544 &lda, &info);
00545
00546
00547
00548 anorm = zlanhb_("1", uplo, &n, &kd, &a[1], &ldab, &
00549 rwork[1]);
00550 ainvnm = zlange_("1", &n, &n, &ainv[1], &lda, &rwork[
00551 1]);
00552 if (anorm <= 0. || ainvnm <= 0.) {
00553 rcondc = 1.;
00554 } else {
00555 rcondc = 1. / anorm / ainvnm;
00556 }
00557
00558 i__5 = *nns;
00559 for (irhs = 1; irhs <= i__5; ++irhs) {
00560 nrhs = nsval[irhs];
00561
00562
00563
00564
00565 s_copy(srnamc_1.srnamt, "ZLARHS", (ftnlen)32, (
00566 ftnlen)6);
00567 zlarhs_(path, xtype, uplo, " ", &n, &n, &kd, &kd,
00568 &nrhs, &a[1], &ldab, &xact[1], &lda, &b[1]
00569 , &lda, iseed, &info);
00570 zlacpy_("Full", &n, &nrhs, &b[1], &lda, &x[1], &
00571 lda);
00572
00573 s_copy(srnamc_1.srnamt, "ZPBTRS", (ftnlen)32, (
00574 ftnlen)6);
00575 zpbtrs_(uplo, &n, &kd, &nrhs, &afac[1], &ldab, &x[
00576 1], &lda, &info);
00577
00578
00579
00580 if (info != 0) {
00581 alaerh_(path, "ZPBTRS", &info, &c__0, uplo, &
00582 n, &n, &kd, &kd, &nrhs, &imat, &nfail,
00583 &nerrs, nout);
00584 }
00585
00586 zlacpy_("Full", &n, &nrhs, &b[1], &lda, &work[1],
00587 &lda);
00588 zpbt02_(uplo, &n, &kd, &nrhs, &a[1], &ldab, &x[1],
00589 &lda, &work[1], &lda, &rwork[1], &result[
00590 1]);
00591
00592
00593
00594
00595 zget04_(&n, &nrhs, &x[1], &lda, &xact[1], &lda, &
00596 rcondc, &result[2]);
00597
00598
00599
00600
00601 s_copy(srnamc_1.srnamt, "ZPBRFS", (ftnlen)32, (
00602 ftnlen)6);
00603 zpbrfs_(uplo, &n, &kd, &nrhs, &a[1], &ldab, &afac[
00604 1], &ldab, &b[1], &lda, &x[1], &lda, &
00605 rwork[1], &rwork[nrhs + 1], &work[1], &
00606 rwork[(nrhs << 1) + 1], &info);
00607
00608
00609
00610 if (info != 0) {
00611 alaerh_(path, "ZPBRFS", &info, &c__0, uplo, &
00612 n, &n, &kd, &kd, &nrhs, &imat, &nfail,
00613 &nerrs, nout);
00614 }
00615
00616 zget04_(&n, &nrhs, &x[1], &lda, &xact[1], &lda, &
00617 rcondc, &result[3]);
00618 zpbt05_(uplo, &n, &kd, &nrhs, &a[1], &ldab, &b[1],
00619 &lda, &x[1], &lda, &xact[1], &lda, &
00620 rwork[1], &rwork[nrhs + 1], &result[4]);
00621
00622
00623
00624
00625 for (k = 2; k <= 6; ++k) {
00626 if (result[k - 1] >= *thresh) {
00627 if (nfail == 0 && nerrs == 0) {
00628 alahd_(nout, path);
00629 }
00630 io___46.ciunit = *nout;
00631 s_wsfe(&io___46);
00632 do_fio(&c__1, uplo, (ftnlen)1);
00633 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(
00634 integer));
00635 do_fio(&c__1, (char *)&kd, (ftnlen)sizeof(
00636 integer));
00637 do_fio(&c__1, (char *)&nrhs, (ftnlen)
00638 sizeof(integer));
00639 do_fio(&c__1, (char *)&imat, (ftnlen)
00640 sizeof(integer));
00641 do_fio(&c__1, (char *)&k, (ftnlen)sizeof(
00642 integer));
00643 do_fio(&c__1, (char *)&result[k - 1], (
00644 ftnlen)sizeof(doublereal));
00645 e_wsfe();
00646 ++nfail;
00647 }
00648
00649 }
00650 nrun += 5;
00651
00652 }
00653
00654
00655
00656
00657 s_copy(srnamc_1.srnamt, "ZPBCON", (ftnlen)32, (ftnlen)
00658 6);
00659 zpbcon_(uplo, &n, &kd, &afac[1], &ldab, &anorm, &
00660 rcond, &work[1], &rwork[1], &info);
00661
00662
00663
00664 if (info != 0) {
00665 alaerh_(path, "ZPBCON", &info, &c__0, uplo, &n, &
00666 n, &kd, &kd, &c_n1, &imat, &nfail, &nerrs,
00667 nout);
00668 }
00669
00670 result[6] = dget06_(&rcond, &rcondc);
00671
00672
00673
00674 if (result[6] >= *thresh) {
00675 if (nfail == 0 && nerrs == 0) {
00676 alahd_(nout, path);
00677 }
00678 io___48.ciunit = *nout;
00679 s_wsfe(&io___48);
00680 do_fio(&c__1, uplo, (ftnlen)1);
00681 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer))
00682 ;
00683 do_fio(&c__1, (char *)&kd, (ftnlen)sizeof(integer)
00684 );
00685 do_fio(&c__1, (char *)&imat, (ftnlen)sizeof(
00686 integer));
00687 do_fio(&c__1, (char *)&c__7, (ftnlen)sizeof(
00688 integer));
00689 do_fio(&c__1, (char *)&result[6], (ftnlen)sizeof(
00690 doublereal));
00691 e_wsfe();
00692 ++nfail;
00693 }
00694 ++nrun;
00695 L50:
00696 ;
00697 }
00698 L60:
00699 ;
00700 }
00701
00702 }
00703
00704 }
00705
00706 }
00707
00708
00709
00710 alasum_(path, nout, &nfail, &nrun, &nerrs);
00711
00712 return 0;
00713
00714
00715
00716 }