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
00037 int ddrvac_(logical *dotype, integer *nm, integer *mval,
00038 integer *nns, integer *nsval, doublereal *thresh, integer *nmax,
00039 doublereal *a, doublereal *afac, doublereal *b, doublereal *x,
00040 doublereal *work, doublereal *rwork, real *swork, integer *nout)
00041 {
00042
00043
00044 static integer iseedy[4] = { 1988,1989,1990,1991 };
00045 static char uplos[1*2] = "U" "L";
00046
00047
00048 static char fmt_9988[] = "(\002 *** \002,a6,\002 returned with INFO ="
00049 "\002,i5,\002 instead of \002,i5,/\002 ==> N =\002,i5,\002, type"
00050 " \002,i2)";
00051 static char fmt_9975[] = "(\002 *** Error code from \002,a6,\002=\002,"
00052 "i5,\002 for M=\002,i5,\002, type \002,i2)";
00053 static char fmt_8999[] = "(/1x,a3,\002: positive definite dense matri"
00054 "ces\002)";
00055 static char fmt_8979[] = "(4x,\0021. Diagonal\002,24x,\0027. Last n/2 co"
00056 "lumns zero\002,/4x,\0022. Upper triangular\002,16x,\0028. Random"
00057 ", CNDNUM = sqrt(0.1/EPS)\002,/4x,\0023. Lower triangular\002,16x,"
00058 "\0029. Random, CNDNUM = 0.1/EPS\002,/4x,\0024. Random, CNDNUM = 2"
00059 "\002,13x,\00210. Scaled near underflow\002,/4x,\0025. First colu"
00060 "mn zero\002,14x,\00211. Scaled near overflow\002,/4x,\0026. Last"
00061 " column zero\002)";
00062 static char fmt_8960[] = "(3x,i2,\002: norm_1( B - A * X ) / \002,\002("
00063 " norm_1(A) * norm_1(X) * EPS * SQRT(N) ) > 1 if ITERREF\002,/4x"
00064 ",\002or norm_1( B - A * X ) / \002,\002( norm_1(A) * norm_1(X) "
00065 "* EPS ) > THRES if DPOTRF\002)";
00066 static char fmt_9998[] = "(\002 UPLO='\002,a1,\002', N =\002,i5,\002, NR"
00067 "HS=\002,i3,\002, type \002,i2,\002, test(\002,i2,\002) =\002,g12"
00068 ".5)";
00069 static char fmt_9996[] = "(1x,a6,\002: \002,i6,\002 out of \002,i6,\002 "
00070 "tests failed to pass the threshold\002)";
00071 static char fmt_9995[] = "(/1x,\002All tests for \002,a6,\002 routines p"
00072 "assed the threshold (\002,i6,\002 tests run)\002)";
00073 static char fmt_9994[] = "(6x,i6,\002 error messages recorded\002)";
00074
00075
00076 integer i__1, i__2, i__3;
00077 cilist ci__1;
00078
00079
00080 int s_copy(char *, char *, ftnlen, ftnlen);
00081 integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00082 double sqrt(doublereal);
00083
00084
00085 integer i__, n, im, kl, ku, lda, ioff, mode, kase, imat, info;
00086 char path[3], dist[1];
00087 integer irhs, iter, nrhs;
00088 char uplo[1], type__[1];
00089 integer nrun;
00090 extern int alahd_(integer *, char *);
00091 integer nfail, iseed[4], nimat;
00092 doublereal anorm;
00093 extern int dpot06_(char *, integer *, integer *,
00094 doublereal *, integer *, doublereal *, integer *, doublereal *,
00095 integer *, doublereal *, doublereal *);
00096 integer iuplo, izero, nerrs;
00097 logical zerot;
00098 char xtype[1];
00099 extern int dlatb4_(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 *), dlacpy_(char *, integer *,
00105 integer *, doublereal *, integer *, doublereal *, integer *), dlarhs_(char *, char *, char *, char *, integer *,
00106 integer *, integer *, integer *, integer *, doublereal *, integer
00107 *, doublereal *, integer *, doublereal *, integer *, integer *,
00108 integer *);
00109 doublereal cndnum;
00110 extern int dlatms_(integer *, integer *, char *, integer
00111 *, char *, doublereal *, integer *, doublereal *, doublereal *,
00112 integer *, integer *, char *, doublereal *, integer *, doublereal
00113 *, integer *);
00114 doublereal result[1];
00115 extern int dsposv_(char *, integer *, integer *,
00116 doublereal *, integer *, doublereal *, integer *, doublereal *,
00117 integer *, doublereal *, real *, integer *, integer *);
00118
00119
00120 static cilist io___32 = { 0, 0, 0, fmt_9988, 0 };
00121 static cilist io___33 = { 0, 0, 0, fmt_9975, 0 };
00122 static cilist io___35 = { 0, 0, 0, fmt_8999, 0 };
00123 static cilist io___36 = { 0, 0, 0, fmt_8979, 0 };
00124 static cilist io___37 = { 0, 0, 0, fmt_8960, 0 };
00125 static cilist io___38 = { 0, 0, 0, fmt_9998, 0 };
00126 static cilist io___39 = { 0, 0, 0, fmt_9996, 0 };
00127 static cilist io___40 = { 0, 0, 0, fmt_9995, 0 };
00128 static cilist io___41 = { 0, 0, 0, fmt_9994, 0 };
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138
00139
00140
00141
00142
00143
00144
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 --swork;
00218 --rwork;
00219 --work;
00220 --x;
00221 --b;
00222 --afac;
00223 --a;
00224 --nsval;
00225 --mval;
00226 --dotype;
00227
00228
00229
00230
00231
00232
00233
00234 kase = 0;
00235 s_copy(path, "Double precision", (ftnlen)1, (ftnlen)16);
00236 s_copy(path + 1, "PO", (ftnlen)2, (ftnlen)2);
00237 nrun = 0;
00238 nfail = 0;
00239 nerrs = 0;
00240 for (i__ = 1; i__ <= 4; ++i__) {
00241 iseed[i__ - 1] = iseedy[i__ - 1];
00242
00243 }
00244
00245 infoc_1.infot = 0;
00246
00247
00248
00249 i__1 = *nm;
00250 for (im = 1; im <= i__1; ++im) {
00251 n = mval[im];
00252 lda = max(n,1);
00253 nimat = 9;
00254 if (n <= 0) {
00255 nimat = 1;
00256 }
00257
00258 i__2 = nimat;
00259 for (imat = 1; imat <= i__2; ++imat) {
00260
00261
00262
00263 if (! dotype[imat]) {
00264 goto L110;
00265 }
00266
00267
00268
00269 zerot = imat >= 3 && imat <= 5;
00270 if (zerot && n < imat - 2) {
00271 goto L110;
00272 }
00273
00274
00275
00276 for (iuplo = 1; iuplo <= 2; ++iuplo) {
00277 *(unsigned char *)uplo = *(unsigned char *)&uplos[iuplo - 1];
00278
00279
00280
00281
00282 dlatb4_(path, &imat, &n, &n, type__, &kl, &ku, &anorm, &mode,
00283 &cndnum, dist);
00284
00285 s_copy(srnamc_1.srnamt, "DLATMS", (ftnlen)32, (ftnlen)6);
00286 dlatms_(&n, &n, dist, iseed, type__, &rwork[1], &mode, &
00287 cndnum, &anorm, &kl, &ku, uplo, &a[1], &lda, &work[1],
00288 &info);
00289
00290
00291
00292 if (info != 0) {
00293 alaerh_(path, "DLATMS", &info, &c__0, uplo, &n, &n, &c_n1,
00294 &c_n1, &c_n1, &imat, &nfail, &nerrs, nout);
00295 goto L100;
00296 }
00297
00298
00299
00300
00301 if (zerot) {
00302 if (imat == 3) {
00303 izero = 1;
00304 } else if (imat == 4) {
00305 izero = n;
00306 } else {
00307 izero = n / 2 + 1;
00308 }
00309 ioff = (izero - 1) * lda;
00310
00311
00312
00313 if (iuplo == 1) {
00314 i__3 = izero - 1;
00315 for (i__ = 1; i__ <= i__3; ++i__) {
00316 a[ioff + i__] = 0.;
00317
00318 }
00319 ioff += izero;
00320 i__3 = n;
00321 for (i__ = izero; i__ <= i__3; ++i__) {
00322 a[ioff] = 0.;
00323 ioff += lda;
00324
00325 }
00326 } else {
00327 ioff = izero;
00328 i__3 = izero - 1;
00329 for (i__ = 1; i__ <= i__3; ++i__) {
00330 a[ioff] = 0.;
00331 ioff += lda;
00332
00333 }
00334 ioff -= izero;
00335 i__3 = n;
00336 for (i__ = izero; i__ <= i__3; ++i__) {
00337 a[ioff + i__] = 0.;
00338
00339 }
00340 }
00341 } else {
00342 izero = 0;
00343 }
00344
00345 i__3 = *nns;
00346 for (irhs = 1; irhs <= i__3; ++irhs) {
00347 nrhs = nsval[irhs];
00348 *(unsigned char *)xtype = 'N';
00349
00350
00351
00352 s_copy(srnamc_1.srnamt, "DLARHS", (ftnlen)32, (ftnlen)6);
00353 dlarhs_(path, xtype, uplo, " ", &n, &n, &kl, &ku, &nrhs, &
00354 a[1], &lda, &x[1], &lda, &b[1], &lda, iseed, &
00355 info);
00356
00357
00358
00359
00360 s_copy(srnamc_1.srnamt, "DSPOSV ", (ftnlen)32, (ftnlen)7);
00361 ++kase;
00362
00363 dlacpy_("All", &n, &n, &a[1], &lda, &afac[1], &lda);
00364
00365 dsposv_(uplo, &n, &nrhs, &afac[1], &lda, &b[1], &lda, &x[
00366 1], &lda, &work[1], &swork[1], &iter, &info);
00367 if (iter < 0) {
00368 dlacpy_("All", &n, &n, &a[1], &lda, &afac[1], &lda);
00369 }
00370
00371
00372
00373 if (info != izero) {
00374
00375 if (nfail == 0 && nerrs == 0) {
00376 alahd_(nout, path);
00377 }
00378 ++nerrs;
00379
00380 if (info != izero && izero != 0) {
00381 io___32.ciunit = *nout;
00382 s_wsfe(&io___32);
00383 do_fio(&c__1, "DSPOSV", (ftnlen)6);
00384 do_fio(&c__1, (char *)&info, (ftnlen)sizeof(
00385 integer));
00386 do_fio(&c__1, (char *)&izero, (ftnlen)sizeof(
00387 integer));
00388 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer))
00389 ;
00390 do_fio(&c__1, (char *)&imat, (ftnlen)sizeof(
00391 integer));
00392 e_wsfe();
00393 } else {
00394 io___33.ciunit = *nout;
00395 s_wsfe(&io___33);
00396 do_fio(&c__1, "DSPOSV", (ftnlen)6);
00397 do_fio(&c__1, (char *)&info, (ftnlen)sizeof(
00398 integer));
00399 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer))
00400 ;
00401 do_fio(&c__1, (char *)&imat, (ftnlen)sizeof(
00402 integer));
00403 e_wsfe();
00404 }
00405 }
00406
00407
00408
00409 if (info != 0) {
00410 goto L110;
00411 }
00412
00413
00414
00415 dlacpy_("All", &n, &nrhs, &b[1], &lda, &work[1], &lda);
00416
00417 dpot06_(uplo, &n, &nrhs, &a[1], &lda, &x[1], &lda, &work[
00418 1], &lda, &rwork[1], result);
00419
00420
00421
00422
00423
00424
00425
00426
00427
00428
00429
00430
00431
00432 if (*thresh <= 0.f || iter >= 0 && n > 0 && result[0] >=
00433 sqrt((doublereal) n) || iter < 0 && result[0] >= *
00434 thresh) {
00435
00436 if (nfail == 0 && nerrs == 0) {
00437 io___35.ciunit = *nout;
00438 s_wsfe(&io___35);
00439 do_fio(&c__1, "DPO", (ftnlen)3);
00440 e_wsfe();
00441 ci__1.cierr = 0;
00442 ci__1.ciunit = *nout;
00443 ci__1.cifmt = "( ' Matrix types:' )";
00444 s_wsfe(&ci__1);
00445 e_wsfe();
00446 io___36.ciunit = *nout;
00447 s_wsfe(&io___36);
00448 e_wsfe();
00449 ci__1.cierr = 0;
00450 ci__1.ciunit = *nout;
00451 ci__1.cifmt = "( ' Test ratios:' )";
00452 s_wsfe(&ci__1);
00453 e_wsfe();
00454 io___37.ciunit = *nout;
00455 s_wsfe(&io___37);
00456 do_fio(&c__1, (char *)&c__1, (ftnlen)sizeof(
00457 integer));
00458 e_wsfe();
00459 ci__1.cierr = 0;
00460 ci__1.ciunit = *nout;
00461 ci__1.cifmt = "( ' Messages:' )";
00462 s_wsfe(&ci__1);
00463 e_wsfe();
00464 }
00465
00466 io___38.ciunit = *nout;
00467 s_wsfe(&io___38);
00468 do_fio(&c__1, uplo, (ftnlen)1);
00469 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00470 do_fio(&c__1, (char *)&nrhs, (ftnlen)sizeof(integer));
00471 do_fio(&c__1, (char *)&imat, (ftnlen)sizeof(integer));
00472 do_fio(&c__1, (char *)&c__1, (ftnlen)sizeof(integer));
00473 do_fio(&c__1, (char *)&result[0], (ftnlen)sizeof(
00474 doublereal));
00475 e_wsfe();
00476
00477 ++nfail;
00478
00479 }
00480
00481 ++nrun;
00482
00483
00484 }
00485 L100:
00486 ;
00487 }
00488 L110:
00489 ;
00490 }
00491
00492 }
00493
00494
00495
00496
00497
00498 if (nfail > 0) {
00499 io___39.ciunit = *nout;
00500 s_wsfe(&io___39);
00501 do_fio(&c__1, "DSPOSV", (ftnlen)6);
00502 do_fio(&c__1, (char *)&nfail, (ftnlen)sizeof(integer));
00503 do_fio(&c__1, (char *)&nrun, (ftnlen)sizeof(integer));
00504 e_wsfe();
00505 } else {
00506 io___40.ciunit = *nout;
00507 s_wsfe(&io___40);
00508 do_fio(&c__1, "DSPOSV", (ftnlen)6);
00509 do_fio(&c__1, (char *)&nrun, (ftnlen)sizeof(integer));
00510 e_wsfe();
00511 }
00512 if (nerrs > 0) {
00513 io___41.ciunit = *nout;
00514 s_wsfe(&io___41);
00515 do_fio(&c__1, (char *)&nerrs, (ftnlen)sizeof(integer));
00516 e_wsfe();
00517 }
00518
00519
00520
00521
00522
00523
00524
00525 return 0;
00526
00527
00528
00529 }