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 doublereal c_b60 = 0.;
00037 static integer c__2 = 2;
00038
00039 int ddrvpp_(logical *dotype, integer *nn, integer *nval,
00040 integer *nrhs, doublereal *thresh, logical *tsterr, integer *nmax,
00041 doublereal *a, doublereal *afac, doublereal *asav, doublereal *b,
00042 doublereal *bsav, doublereal *x, doublereal *xact, doublereal *s,
00043 doublereal *work, doublereal *rwork, integer *iwork, integer *nout)
00044 {
00045
00046
00047 static integer iseedy[4] = { 1988,1989,1990,1991 };
00048 static char uplos[1*2] = "U" "L";
00049 static char facts[1*3] = "F" "N" "E";
00050 static char packs[1*2] = "C" "R";
00051 static char equeds[1*2] = "N" "Y";
00052
00053
00054 static char fmt_9999[] = "(1x,a,\002, UPLO='\002,a1,\002', N =\002,i5"
00055 ",\002, type \002,i1,\002, test(\002,i1,\002)=\002,g12.5)";
00056 static char fmt_9997[] = "(1x,a,\002, FACT='\002,a1,\002', UPLO='\002,"
00057 "a1,\002', N=\002,i5,\002, EQUED='\002,a1,\002', type \002,i1,"
00058 "\002, test(\002,i1,\002)=\002,g12.5)";
00059 static char fmt_9998[] = "(1x,a,\002, FACT='\002,a1,\002', UPLO='\002,"
00060 "a1,\002', N=\002,i5,\002, type \002,i1,\002, test(\002,i1,\002)"
00061 "=\002,g12.5)";
00062
00063
00064 address a__1[2];
00065 integer i__1, i__2, i__3, i__4, i__5[2];
00066 char ch__1[2];
00067
00068
00069 int s_copy(char *, char *, ftnlen, ftnlen);
00070 integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00071 int s_cat(char *, char **, integer *, integer *, ftnlen);
00072
00073
00074 integer i__, k, n, k1, in, kl, ku, nt, lda, npp;
00075 char fact[1];
00076 integer ioff, mode;
00077 doublereal amax;
00078 char path[3];
00079 integer imat, info;
00080 char dist[1], uplo[1], type__[1];
00081 integer nrun, ifact;
00082 extern int dget04_(integer *, integer *, doublereal *,
00083 integer *, doublereal *, integer *, doublereal *, doublereal *);
00084 integer nfail, iseed[4], nfact;
00085 extern doublereal dget06_(doublereal *, doublereal *);
00086 extern logical lsame_(char *, char *);
00087 char equed[1];
00088 doublereal roldc, rcond, scond;
00089 integer nimat;
00090 extern int dppt01_(char *, integer *, doublereal *,
00091 doublereal *, doublereal *, doublereal *), dppt02_(char *,
00092 integer *, integer *, doublereal *, doublereal *, integer *,
00093 doublereal *, integer *, doublereal *, doublereal *);
00094 doublereal anorm;
00095 extern int dppt05_(char *, integer *, integer *,
00096 doublereal *, doublereal *, integer *, doublereal *, integer *,
00097 doublereal *, integer *, doublereal *, doublereal *, doublereal *), dcopy_(integer *, doublereal *, integer *, doublereal *,
00098 integer *);
00099 logical equil;
00100 integer iuplo, izero, nerrs;
00101 extern int dppsv_(char *, integer *, integer *,
00102 doublereal *, doublereal *, integer *, integer *);
00103 logical zerot;
00104 char xtype[1];
00105 extern int dlatb4_(char *, integer *, integer *, integer
00106 *, char *, integer *, integer *, doublereal *, integer *,
00107 doublereal *, char *), aladhd_(integer *,
00108 char *), alaerh_(char *, char *, integer *, integer *,
00109 char *, integer *, integer *, integer *, integer *, integer *,
00110 integer *, integer *, integer *, integer *);
00111 logical prefac;
00112 doublereal rcondc;
00113 logical nofact;
00114 char packit[1];
00115 integer iequed;
00116 extern int dlacpy_(char *, integer *, integer *,
00117 doublereal *, integer *, doublereal *, integer *),
00118 dlarhs_(char *, char *, char *, char *, integer *, integer *,
00119 integer *, integer *, integer *, doublereal *, integer *,
00120 doublereal *, integer *, doublereal *, integer *, integer *,
00121 integer *), dlaset_(char *,
00122 integer *, integer *, doublereal *, doublereal *, doublereal *,
00123 integer *);
00124 extern doublereal dlansp_(char *, char *, integer *, doublereal *,
00125 doublereal *);
00126 extern int alasvm_(char *, integer *, integer *, integer
00127 *, integer *);
00128 doublereal cndnum;
00129 extern int dlaqsp_(char *, integer *, doublereal *,
00130 doublereal *, doublereal *, doublereal *, char *),
00131 dlatms_(integer *, integer *, char *, integer *, char *,
00132 doublereal *, integer *, doublereal *, doublereal *, integer *,
00133 integer *, char *, doublereal *, integer *, doublereal *, integer
00134 *);
00135 doublereal ainvnm;
00136 extern int dppequ_(char *, integer *, doublereal *,
00137 doublereal *, doublereal *, doublereal *, integer *),
00138 dpptrf_(char *, integer *, doublereal *, integer *),
00139 dpptri_(char *, integer *, doublereal *, integer *),
00140 derrvx_(char *, integer *);
00141 doublereal result[6];
00142 extern int dppsvx_(char *, char *, integer *, integer *,
00143 doublereal *, doublereal *, char *, doublereal *, doublereal *,
00144 integer *, doublereal *, integer *, doublereal *, doublereal *,
00145 doublereal *, doublereal *, integer *, integer *);
00146
00147
00148 static cilist io___49 = { 0, 0, 0, fmt_9999, 0 };
00149 static cilist io___52 = { 0, 0, 0, fmt_9997, 0 };
00150 static cilist io___53 = { 0, 0, 0, fmt_9998, 0 };
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 --iwork;
00248 --rwork;
00249 --work;
00250 --s;
00251 --xact;
00252 --x;
00253 --bsav;
00254 --b;
00255 --asav;
00256 --afac;
00257 --a;
00258 --nval;
00259 --dotype;
00260
00261
00262
00263
00264
00265
00266
00267 s_copy(path, "Double precision", (ftnlen)1, (ftnlen)16);
00268 s_copy(path + 1, "PP", (ftnlen)2, (ftnlen)2);
00269 nrun = 0;
00270 nfail = 0;
00271 nerrs = 0;
00272 for (i__ = 1; i__ <= 4; ++i__) {
00273 iseed[i__ - 1] = iseedy[i__ - 1];
00274
00275 }
00276
00277
00278
00279 if (*tsterr) {
00280 derrvx_(path, nout);
00281 }
00282 infoc_1.infot = 0;
00283
00284
00285
00286 i__1 = *nn;
00287 for (in = 1; in <= i__1; ++in) {
00288 n = nval[in];
00289 lda = max(n,1);
00290 npp = n * (n + 1) / 2;
00291 *(unsigned char *)xtype = 'N';
00292 nimat = 9;
00293 if (n <= 0) {
00294 nimat = 1;
00295 }
00296
00297 i__2 = nimat;
00298 for (imat = 1; imat <= i__2; ++imat) {
00299
00300
00301
00302 if (! dotype[imat]) {
00303 goto L130;
00304 }
00305
00306
00307
00308 zerot = imat >= 3 && imat <= 5;
00309 if (zerot && n < imat - 2) {
00310 goto L130;
00311 }
00312
00313
00314
00315 for (iuplo = 1; iuplo <= 2; ++iuplo) {
00316 *(unsigned char *)uplo = *(unsigned char *)&uplos[iuplo - 1];
00317 *(unsigned char *)packit = *(unsigned char *)&packs[iuplo - 1]
00318 ;
00319
00320
00321
00322
00323 dlatb4_(path, &imat, &n, &n, type__, &kl, &ku, &anorm, &mode,
00324 &cndnum, dist);
00325 rcondc = 1. / cndnum;
00326
00327 s_copy(srnamc_1.srnamt, "DLATMS", (ftnlen)32, (ftnlen)6);
00328 dlatms_(&n, &n, dist, iseed, type__, &rwork[1], &mode, &
00329 cndnum, &anorm, &kl, &ku, packit, &a[1], &lda, &work[
00330 1], &info);
00331
00332
00333
00334 if (info != 0) {
00335 alaerh_(path, "DLATMS", &info, &c__0, uplo, &n, &n, &c_n1,
00336 &c_n1, &c_n1, &imat, &nfail, &nerrs, nout);
00337 goto L120;
00338 }
00339
00340
00341
00342
00343 if (zerot) {
00344 if (imat == 3) {
00345 izero = 1;
00346 } else if (imat == 4) {
00347 izero = n;
00348 } else {
00349 izero = n / 2 + 1;
00350 }
00351
00352
00353
00354 if (iuplo == 1) {
00355 ioff = (izero - 1) * izero / 2;
00356 i__3 = izero - 1;
00357 for (i__ = 1; i__ <= i__3; ++i__) {
00358 a[ioff + i__] = 0.;
00359
00360 }
00361 ioff += izero;
00362 i__3 = n;
00363 for (i__ = izero; i__ <= i__3; ++i__) {
00364 a[ioff] = 0.;
00365 ioff += i__;
00366
00367 }
00368 } else {
00369 ioff = izero;
00370 i__3 = izero - 1;
00371 for (i__ = 1; i__ <= i__3; ++i__) {
00372 a[ioff] = 0.;
00373 ioff = ioff + n - i__;
00374
00375 }
00376 ioff -= izero;
00377 i__3 = n;
00378 for (i__ = izero; i__ <= i__3; ++i__) {
00379 a[ioff + i__] = 0.;
00380
00381 }
00382 }
00383 } else {
00384 izero = 0;
00385 }
00386
00387
00388
00389 dcopy_(&npp, &a[1], &c__1, &asav[1], &c__1);
00390
00391 for (iequed = 1; iequed <= 2; ++iequed) {
00392 *(unsigned char *)equed = *(unsigned char *)&equeds[
00393 iequed - 1];
00394 if (iequed == 1) {
00395 nfact = 3;
00396 } else {
00397 nfact = 1;
00398 }
00399
00400 i__3 = nfact;
00401 for (ifact = 1; ifact <= i__3; ++ifact) {
00402 *(unsigned char *)fact = *(unsigned char *)&facts[
00403 ifact - 1];
00404 prefac = lsame_(fact, "F");
00405 nofact = lsame_(fact, "N");
00406 equil = lsame_(fact, "E");
00407
00408 if (zerot) {
00409 if (prefac) {
00410 goto L100;
00411 }
00412 rcondc = 0.;
00413
00414 } else if (! lsame_(fact, "N"))
00415 {
00416
00417
00418
00419
00420
00421
00422 dcopy_(&npp, &asav[1], &c__1, &afac[1], &c__1);
00423 if (equil || iequed > 1) {
00424
00425
00426
00427
00428 dppequ_(uplo, &n, &afac[1], &s[1], &scond, &
00429 amax, &info);
00430 if (info == 0 && n > 0) {
00431 if (iequed > 1) {
00432 scond = 0.;
00433 }
00434
00435
00436
00437 dlaqsp_(uplo, &n, &afac[1], &s[1], &scond,
00438 &amax, equed);
00439 }
00440 }
00441
00442
00443
00444
00445 if (equil) {
00446 roldc = rcondc;
00447 }
00448
00449
00450
00451 anorm = dlansp_("1", uplo, &n, &afac[1], &rwork[1]
00452 );
00453
00454
00455
00456 dpptrf_(uplo, &n, &afac[1], &info);
00457
00458
00459
00460 dcopy_(&npp, &afac[1], &c__1, &a[1], &c__1);
00461 dpptri_(uplo, &n, &a[1], &info);
00462
00463
00464
00465 ainvnm = dlansp_("1", uplo, &n, &a[1], &rwork[1]);
00466 if (anorm <= 0. || ainvnm <= 0.) {
00467 rcondc = 1.;
00468 } else {
00469 rcondc = 1. / anorm / ainvnm;
00470 }
00471 }
00472
00473
00474
00475 dcopy_(&npp, &asav[1], &c__1, &a[1], &c__1);
00476
00477
00478
00479 s_copy(srnamc_1.srnamt, "DLARHS", (ftnlen)32, (ftnlen)
00480 6);
00481 dlarhs_(path, xtype, uplo, " ", &n, &n, &kl, &ku,
00482 nrhs, &a[1], &lda, &xact[1], &lda, &b[1], &
00483 lda, iseed, &info);
00484 *(unsigned char *)xtype = 'C';
00485 dlacpy_("Full", &n, nrhs, &b[1], &lda, &bsav[1], &lda);
00486
00487 if (nofact) {
00488
00489
00490
00491
00492
00493
00494 dcopy_(&npp, &a[1], &c__1, &afac[1], &c__1);
00495 dlacpy_("Full", &n, nrhs, &b[1], &lda, &x[1], &
00496 lda);
00497
00498 s_copy(srnamc_1.srnamt, "DPPSV ", (ftnlen)32, (
00499 ftnlen)6);
00500 dppsv_(uplo, &n, nrhs, &afac[1], &x[1], &lda, &
00501 info);
00502
00503
00504
00505 if (info != izero) {
00506 alaerh_(path, "DPPSV ", &info, &izero, uplo, &
00507 n, &n, &c_n1, &c_n1, nrhs, &imat, &
00508 nfail, &nerrs, nout);
00509 goto L70;
00510 } else if (info != 0) {
00511 goto L70;
00512 }
00513
00514
00515
00516
00517 dppt01_(uplo, &n, &a[1], &afac[1], &rwork[1],
00518 result);
00519
00520
00521
00522 dlacpy_("Full", &n, nrhs, &b[1], &lda, &work[1], &
00523 lda);
00524 dppt02_(uplo, &n, nrhs, &a[1], &x[1], &lda, &work[
00525 1], &lda, &rwork[1], &result[1]);
00526
00527
00528
00529 dget04_(&n, nrhs, &x[1], &lda, &xact[1], &lda, &
00530 rcondc, &result[2]);
00531 nt = 3;
00532
00533
00534
00535
00536 i__4 = nt;
00537 for (k = 1; k <= i__4; ++k) {
00538 if (result[k - 1] >= *thresh) {
00539 if (nfail == 0 && nerrs == 0) {
00540 aladhd_(nout, path);
00541 }
00542 io___49.ciunit = *nout;
00543 s_wsfe(&io___49);
00544 do_fio(&c__1, "DPPSV ", (ftnlen)6);
00545 do_fio(&c__1, uplo, (ftnlen)1);
00546 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(
00547 integer));
00548 do_fio(&c__1, (char *)&imat, (ftnlen)
00549 sizeof(integer));
00550 do_fio(&c__1, (char *)&k, (ftnlen)sizeof(
00551 integer));
00552 do_fio(&c__1, (char *)&result[k - 1], (
00553 ftnlen)sizeof(doublereal));
00554 e_wsfe();
00555 ++nfail;
00556 }
00557
00558 }
00559 nrun += nt;
00560 L70:
00561 ;
00562 }
00563
00564
00565
00566 if (! prefac && npp > 0) {
00567 dlaset_("Full", &npp, &c__1, &c_b60, &c_b60, &
00568 afac[1], &npp);
00569 }
00570 dlaset_("Full", &n, nrhs, &c_b60, &c_b60, &x[1], &lda);
00571 if (iequed > 1 && n > 0) {
00572
00573
00574
00575
00576 dlaqsp_(uplo, &n, &a[1], &s[1], &scond, &amax,
00577 equed);
00578 }
00579
00580
00581
00582
00583 s_copy(srnamc_1.srnamt, "DPPSVX", (ftnlen)32, (ftnlen)
00584 6);
00585 dppsvx_(fact, uplo, &n, nrhs, &a[1], &afac[1], equed,
00586 &s[1], &b[1], &lda, &x[1], &lda, &rcond, &
00587 rwork[1], &rwork[*nrhs + 1], &work[1], &iwork[
00588 1], &info);
00589
00590
00591
00592 if (info != izero) {
00593
00594 i__5[0] = 1, a__1[0] = fact;
00595 i__5[1] = 1, a__1[1] = uplo;
00596 s_cat(ch__1, a__1, i__5, &c__2, (ftnlen)2);
00597 alaerh_(path, "DPPSVX", &info, &izero, ch__1, &n,
00598 &n, &c_n1, &c_n1, nrhs, &imat, &nfail, &
00599 nerrs, nout);
00600 goto L90;
00601 }
00602
00603 if (info == 0) {
00604 if (! prefac) {
00605
00606
00607
00608
00609 dppt01_(uplo, &n, &a[1], &afac[1], &rwork[(*
00610 nrhs << 1) + 1], result);
00611 k1 = 1;
00612 } else {
00613 k1 = 2;
00614 }
00615
00616
00617
00618 dlacpy_("Full", &n, nrhs, &bsav[1], &lda, &work[1]
00619 , &lda);
00620 dppt02_(uplo, &n, nrhs, &asav[1], &x[1], &lda, &
00621 work[1], &lda, &rwork[(*nrhs << 1) + 1], &
00622 result[1]);
00623
00624
00625
00626 if (nofact || prefac && lsame_(equed, "N")) {
00627 dget04_(&n, nrhs, &x[1], &lda, &xact[1], &lda,
00628 &rcondc, &result[2]);
00629 } else {
00630 dget04_(&n, nrhs, &x[1], &lda, &xact[1], &lda,
00631 &roldc, &result[2]);
00632 }
00633
00634
00635
00636
00637 dppt05_(uplo, &n, nrhs, &asav[1], &b[1], &lda, &x[
00638 1], &lda, &xact[1], &lda, &rwork[1], &
00639 rwork[*nrhs + 1], &result[3]);
00640 } else {
00641 k1 = 6;
00642 }
00643
00644
00645
00646
00647 result[5] = dget06_(&rcond, &rcondc);
00648
00649
00650
00651
00652 for (k = k1; k <= 6; ++k) {
00653 if (result[k - 1] >= *thresh) {
00654 if (nfail == 0 && nerrs == 0) {
00655 aladhd_(nout, path);
00656 }
00657 if (prefac) {
00658 io___52.ciunit = *nout;
00659 s_wsfe(&io___52);
00660 do_fio(&c__1, "DPPSVX", (ftnlen)6);
00661 do_fio(&c__1, fact, (ftnlen)1);
00662 do_fio(&c__1, uplo, (ftnlen)1);
00663 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(
00664 integer));
00665 do_fio(&c__1, equed, (ftnlen)1);
00666 do_fio(&c__1, (char *)&imat, (ftnlen)
00667 sizeof(integer));
00668 do_fio(&c__1, (char *)&k, (ftnlen)sizeof(
00669 integer));
00670 do_fio(&c__1, (char *)&result[k - 1], (
00671 ftnlen)sizeof(doublereal));
00672 e_wsfe();
00673 } else {
00674 io___53.ciunit = *nout;
00675 s_wsfe(&io___53);
00676 do_fio(&c__1, "DPPSVX", (ftnlen)6);
00677 do_fio(&c__1, fact, (ftnlen)1);
00678 do_fio(&c__1, uplo, (ftnlen)1);
00679 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(
00680 integer));
00681 do_fio(&c__1, (char *)&imat, (ftnlen)
00682 sizeof(integer));
00683 do_fio(&c__1, (char *)&k, (ftnlen)sizeof(
00684 integer));
00685 do_fio(&c__1, (char *)&result[k - 1], (
00686 ftnlen)sizeof(doublereal));
00687 e_wsfe();
00688 }
00689 ++nfail;
00690 }
00691
00692 }
00693 nrun = nrun + 7 - k1;
00694 L90:
00695 L100:
00696 ;
00697 }
00698
00699 }
00700 L120:
00701 ;
00702 }
00703 L130:
00704 ;
00705 }
00706
00707 }
00708
00709
00710
00711 alasvm_(path, nout, &nfail, &nrun, &nerrs);
00712
00713 return 0;
00714
00715
00716
00717 }