cdrvpb.c
Go to the documentation of this file.
00001 /* cdrvpb.f -- translated by f2c (version 20061008).
00002    You must link the resulting object file with libf2c:
00003         on Microsoft Windows system, link with libf2c.lib;
00004         on Linux or Unix systems, link with .../path/to/libf2c.a -lm
00005         or, if you install libf2c.a in a standard place, with -lf2c -lm
00006         -- in that order, at the end of the command line, as in
00007                 cc *.o -lf2c -lm
00008         Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
00009 
00010                 http://www.netlib.org/f2c/libf2c.zip
00011 */
00012 
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015 
00016 /* Common Block Declarations */
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 /* Table of constant values */
00032 
00033 static integer c__1 = 1;
00034 static integer c__2 = 2;
00035 static integer c__0 = 0;
00036 static integer c_n1 = -1;
00037 static complex c_b47 = {0.f,0.f};
00038 static complex c_b48 = {1.f,0.f};
00039 
00040 /* Subroutine */ int cdrvpb_(logical *dotype, integer *nn, integer *nval, 
00041         integer *nrhs, real *thresh, logical *tsterr, integer *nmax, complex *
00042         a, complex *afac, complex *asav, complex *b, complex *bsav, complex *
00043         x, complex *xact, real *s, complex *work, real *rwork, integer *nout)
00044 {
00045     /* Initialized data */
00046 
00047     static integer iseedy[4] = { 1988,1989,1990,1991 };
00048     static char facts[1*3] = "F" "N" "E";
00049     static char equeds[1*2] = "N" "Y";
00050 
00051     /* Format strings */
00052     static char fmt_9999[] = "(1x,a,\002, UPLO='\002,a1,\002', N =\002,i5"
00053             ",\002, KD =\002,i5,\002, type \002,i1,\002, test(\002,i1,\002)"
00054             "=\002,g12.5)";
00055     static char fmt_9997[] = "(1x,a,\002( '\002,a1,\002', '\002,a1,\002',"
00056             " \002,i5,\002, \002,i5,\002, ... ), EQUED='\002,a1,\002', type"
00057             " \002,i1,\002, test(\002,i1,\002)=\002,g12.5)";
00058     static char fmt_9998[] = "(1x,a,\002( '\002,a1,\002', '\002,a1,\002',"
00059             " \002,i5,\002, \002,i5,\002, ... ), type \002,i1,\002, test(\002"
00060             ",i1,\002)=\002,g12.5)";
00061 
00062     /* System generated locals */
00063     address a__1[2];
00064     integer i__1, i__2, i__3, i__4, i__5, i__6, i__7[2];
00065     char ch__1[2];
00066 
00067     /* Builtin functions */
00068     /* Subroutine */ int s_copy(char *, char *, ftnlen, ftnlen);
00069     integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00070     /* Subroutine */ int s_cat(char *, char **, integer *, integer *, ftnlen);
00071 
00072     /* Local variables */
00073     integer i__, k, n, i1, i2, k1, kd, nb, in, kl, iw, ku, nt, lda, ikd, nkd, 
00074             ldab;
00075     char fact[1];
00076     integer ioff, mode, koff;
00077     real amax;
00078     char path[3];
00079     integer imat, info;
00080     char dist[1], uplo[1], type__[1];
00081     integer nrun, ifact;
00082     extern /* Subroutine */ int cget04_(integer *, integer *, complex *, 
00083             integer *, complex *, integer *, real *, real *);
00084     integer nfail, iseed[4], nfact;
00085     extern /* Subroutine */ int cpbt01_(char *, integer *, integer *, complex 
00086             *, integer *, complex *, integer *, real *, real *), 
00087             cpbt02_(char *, integer *, integer *, integer *, complex *, 
00088             integer *, complex *, integer *, complex *, integer *, real *, 
00089             real *), cpbt05_(char *, integer *, integer *, integer *, 
00090             complex *, integer *, complex *, integer *, complex *, integer *, 
00091             complex *, integer *, real *, real *, real *);
00092     integer kdval[4];
00093     extern logical lsame_(char *, char *);
00094     char equed[1];
00095     integer nbmin;
00096     real rcond, roldc, scond;
00097     integer nimat;
00098     extern doublereal sget06_(real *, real *);
00099     real anorm;
00100     extern /* Subroutine */ int ccopy_(integer *, complex *, integer *, 
00101             complex *, integer *), cpbsv_(char *, integer *, integer *, 
00102             integer *, complex *, integer *, complex *, integer *, integer *);
00103     logical equil;
00104     extern /* Subroutine */ int cswap_(integer *, complex *, integer *, 
00105             complex *, integer *);
00106     integer iuplo, izero, nerrs;
00107     logical zerot;
00108     char xtype[1];
00109     extern /* Subroutine */ int clatb4_(char *, integer *, integer *, integer 
00110             *, char *, integer *, integer *, real *, integer *, real *, char *
00111 ), aladhd_(integer *, char *);
00112     extern doublereal clanhb_(char *, char *, integer *, integer *, complex *, 
00113              integer *, real *), clange_(char *, integer *, 
00114             integer *, complex *, integer *, real *);
00115     extern /* Subroutine */ int claqhb_(char *, integer *, integer *, complex 
00116             *, integer *, real *, real *, real *, char *), 
00117             alaerh_(char *, char *, integer *, integer *, char *, integer *, 
00118             integer *, integer *, integer *, integer *, integer *, integer *, 
00119             integer *, integer *), claipd_(integer *, 
00120             complex *, integer *, integer *);
00121     logical prefac;
00122     real rcondc;
00123     logical nofact;
00124     char packit[1];
00125     integer iequed;
00126     extern /* Subroutine */ int clacpy_(char *, integer *, integer *, complex 
00127             *, integer *, complex *, integer *), clarhs_(char *, char 
00128             *, char *, char *, integer *, integer *, integer *, integer *, 
00129             integer *, complex *, integer *, complex *, integer *, complex *, 
00130             integer *, integer *, integer *), 
00131             claset_(char *, integer *, integer *, complex *, complex *, 
00132             complex *, integer *), cpbequ_(char *, integer *, integer 
00133             *, complex *, integer *, real *, real *, real *, integer *), alasvm_(char *, integer *, integer *, integer *, integer 
00134             *);
00135     real cndnum;
00136     extern /* Subroutine */ int clatms_(integer *, integer *, char *, integer 
00137             *, char *, real *, integer *, real *, real *, integer *, integer *
00138 , char *, complex *, integer *, complex *, integer *), cpbtrf_(char *, integer *, integer *, complex *, 
00139             integer *, integer *);
00140     real ainvnm;
00141     extern /* Subroutine */ int cpbtrs_(char *, integer *, integer *, integer 
00142             *, complex *, integer *, complex *, integer *, integer *),
00143              xlaenv_(integer *, integer *), cpbsvx_(char *, char *, integer *, 
00144              integer *, integer *, complex *, integer *, complex *, integer *, 
00145              char *, real *, complex *, integer *, complex *, integer *, real 
00146             *, real *, real *, complex *, real *, integer *), cerrvx_(char *, integer *);
00147     real result[6];
00148 
00149     /* Fortran I/O blocks */
00150     static cilist io___57 = { 0, 0, 0, fmt_9999, 0 };
00151     static cilist io___60 = { 0, 0, 0, fmt_9997, 0 };
00152     static cilist io___61 = { 0, 0, 0, fmt_9998, 0 };
00153 
00154 
00155 
00156 /*  -- LAPACK test routine (version 3.1) -- */
00157 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00158 /*     November 2006 */
00159 
00160 /*     .. Scalar Arguments .. */
00161 /*     .. */
00162 /*     .. Array Arguments .. */
00163 /*     .. */
00164 
00165 /*  Purpose */
00166 /*  ======= */
00167 
00168 /*  CDRVPB tests the driver routines CPBSV and -SVX. */
00169 
00170 /*  Arguments */
00171 /*  ========= */
00172 
00173 /*  DOTYPE  (input) LOGICAL array, dimension (NTYPES) */
00174 /*          The matrix types to be used for testing.  Matrices of type j */
00175 /*          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) = */
00176 /*          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used. */
00177 
00178 /*  NN      (input) INTEGER */
00179 /*          The number of values of N contained in the vector NVAL. */
00180 
00181 /*  NVAL    (input) INTEGER array, dimension (NN) */
00182 /*          The values of the matrix dimension N. */
00183 
00184 /*  NRHS    (input) INTEGER */
00185 /*          The number of right hand side vectors to be generated for */
00186 /*          each linear system. */
00187 
00188 /*  THRESH  (input) REAL */
00189 /*          The threshold value for the test ratios.  A result is */
00190 /*          included in the output file if RESULT >= THRESH.  To have */
00191 /*          every test ratio printed, use THRESH = 0. */
00192 
00193 /*  TSTERR  (input) LOGICAL */
00194 /*          Flag that indicates whether error exits are to be tested. */
00195 
00196 /*  NMAX    (input) INTEGER */
00197 /*          The maximum value permitted for N, used in dimensioning the */
00198 /*          work arrays. */
00199 
00200 /*  A       (workspace) COMPLEX array, dimension (NMAX*NMAX) */
00201 
00202 /*  AFAC    (workspace) COMPLEX array, dimension (NMAX*NMAX) */
00203 
00204 /*  ASAV    (workspace) COMPLEX array, dimension (NMAX*NMAX) */
00205 
00206 /*  B       (workspace) COMPLEX array, dimension (NMAX*NRHS) */
00207 
00208 /*  BSAV    (workspace) COMPLEX array, dimension (NMAX*NRHS) */
00209 
00210 /*  X       (workspace) COMPLEX array, dimension (NMAX*NRHS) */
00211 
00212 /*  XACT    (workspace) COMPLEX array, dimension (NMAX*NRHS) */
00213 
00214 /*  S       (workspace) REAL array, dimension (NMAX) */
00215 
00216 /*  WORK    (workspace) COMPLEX array, dimension */
00217 /*                      (NMAX*max(3,NRHS)) */
00218 
00219 /*  RWORK   (workspace) REAL array, dimension (NMAX+2*NRHS) */
00220 
00221 /*  NOUT    (input) INTEGER */
00222 /*          The unit number for output. */
00223 
00224 /*  ===================================================================== */
00225 
00226 /*     .. Parameters .. */
00227 /*     .. */
00228 /*     .. Local Scalars .. */
00229 /*     .. */
00230 /*     .. Local Arrays .. */
00231 /*     .. */
00232 /*     .. External Functions .. */
00233 /*     .. */
00234 /*     .. External Subroutines .. */
00235 /*     .. */
00236 /*     .. Intrinsic Functions .. */
00237 /*     .. */
00238 /*     .. Scalars in Common .. */
00239 /*     .. */
00240 /*     .. Common blocks .. */
00241 /*     .. */
00242 /*     .. Data statements .. */
00243     /* Parameter adjustments */
00244     --rwork;
00245     --work;
00246     --s;
00247     --xact;
00248     --x;
00249     --bsav;
00250     --b;
00251     --asav;
00252     --afac;
00253     --a;
00254     --nval;
00255     --dotype;
00256 
00257     /* Function Body */
00258 /*     .. */
00259 /*     .. Executable Statements .. */
00260 
00261 /*     Initialize constants and the random number seed. */
00262 
00263     s_copy(path, "Complex precision", (ftnlen)1, (ftnlen)17);
00264     s_copy(path + 1, "PB", (ftnlen)2, (ftnlen)2);
00265     nrun = 0;
00266     nfail = 0;
00267     nerrs = 0;
00268     for (i__ = 1; i__ <= 4; ++i__) {
00269         iseed[i__ - 1] = iseedy[i__ - 1];
00270 /* L10: */
00271     }
00272 
00273 /*     Test the error exits */
00274 
00275     if (*tsterr) {
00276         cerrvx_(path, nout);
00277     }
00278     infoc_1.infot = 0;
00279     kdval[0] = 0;
00280 
00281 /*     Set the block size and minimum block size for testing. */
00282 
00283     nb = 1;
00284     nbmin = 2;
00285     xlaenv_(&c__1, &nb);
00286     xlaenv_(&c__2, &nbmin);
00287 
00288 /*     Do for each value of N in NVAL */
00289 
00290     i__1 = *nn;
00291     for (in = 1; in <= i__1; ++in) {
00292         n = nval[in];
00293         lda = max(n,1);
00294         *(unsigned char *)xtype = 'N';
00295 
00296 /*        Set limits on the number of loop iterations. */
00297 
00298 /* Computing MAX */
00299         i__2 = 1, i__3 = min(n,4);
00300         nkd = max(i__2,i__3);
00301         nimat = 8;
00302         if (n == 0) {
00303             nimat = 1;
00304         }
00305 
00306         kdval[1] = n + (n + 1) / 4;
00307         kdval[2] = (n * 3 - 1) / 4;
00308         kdval[3] = (n + 1) / 4;
00309 
00310         i__2 = nkd;
00311         for (ikd = 1; ikd <= i__2; ++ikd) {
00312 
00313 /*           Do for KD = 0, (5*N+1)/4, (3N-1)/4, and (N+1)/4. This order */
00314 /*           makes it easier to skip redundant values for small values */
00315 /*           of N. */
00316 
00317             kd = kdval[ikd - 1];
00318             ldab = kd + 1;
00319 
00320 /*           Do first for UPLO = 'U', then for UPLO = 'L' */
00321 
00322             for (iuplo = 1; iuplo <= 2; ++iuplo) {
00323                 koff = 1;
00324                 if (iuplo == 1) {
00325                     *(unsigned char *)uplo = 'U';
00326                     *(unsigned char *)packit = 'Q';
00327 /* Computing MAX */
00328                     i__3 = 1, i__4 = kd + 2 - n;
00329                     koff = max(i__3,i__4);
00330                 } else {
00331                     *(unsigned char *)uplo = 'L';
00332                     *(unsigned char *)packit = 'B';
00333                 }
00334 
00335                 i__3 = nimat;
00336                 for (imat = 1; imat <= i__3; ++imat) {
00337 
00338 /*                 Do the tests only if DOTYPE( IMAT ) is true. */
00339 
00340                     if (! dotype[imat]) {
00341                         goto L80;
00342                     }
00343 
00344 /*                 Skip types 2, 3, or 4 if the matrix size is too small. */
00345 
00346                     zerot = imat >= 2 && imat <= 4;
00347                     if (zerot && n < imat - 1) {
00348                         goto L80;
00349                     }
00350 
00351                     if (! zerot || ! dotype[1]) {
00352 
00353 /*                    Set up parameters with CLATB4 and generate a test */
00354 /*                    matrix with CLATMS. */
00355 
00356                         clatb4_(path, &imat, &n, &n, type__, &kl, &ku, &anorm, 
00357                                  &mode, &cndnum, dist);
00358 
00359                         s_copy(srnamc_1.srnamt, "CLATMS", (ftnlen)32, (ftnlen)
00360                                 6);
00361                         clatms_(&n, &n, dist, iseed, type__, &rwork[1], &mode, 
00362                                  &cndnum, &anorm, &kd, &kd, packit, &a[koff], 
00363                                 &ldab, &work[1], &info);
00364 
00365 /*                    Check error code from CLATMS. */
00366 
00367                         if (info != 0) {
00368                             alaerh_(path, "CLATMS", &info, &c__0, uplo, &n, &
00369                                     n, &c_n1, &c_n1, &c_n1, &imat, &nfail, &
00370                                     nerrs, nout);
00371                             goto L80;
00372                         }
00373                     } else if (izero > 0) {
00374 
00375 /*                    Use the same matrix for types 3 and 4 as for type */
00376 /*                    2 by copying back the zeroed out column, */
00377 
00378                         iw = (lda << 1) + 1;
00379                         if (iuplo == 1) {
00380                             ioff = (izero - 1) * ldab + kd + 1;
00381                             i__4 = izero - i1;
00382                             ccopy_(&i__4, &work[iw], &c__1, &a[ioff - izero + 
00383                                     i1], &c__1);
00384                             iw = iw + izero - i1;
00385                             i__4 = i2 - izero + 1;
00386 /* Computing MAX */
00387                             i__6 = ldab - 1;
00388                             i__5 = max(i__6,1);
00389                             ccopy_(&i__4, &work[iw], &c__1, &a[ioff], &i__5);
00390                         } else {
00391                             ioff = (i1 - 1) * ldab + 1;
00392                             i__4 = izero - i1;
00393 /* Computing MAX */
00394                             i__6 = ldab - 1;
00395                             i__5 = max(i__6,1);
00396                             ccopy_(&i__4, &work[iw], &c__1, &a[ioff + izero - 
00397                                     i1], &i__5);
00398                             ioff = (izero - 1) * ldab + 1;
00399                             iw = iw + izero - i1;
00400                             i__4 = i2 - izero + 1;
00401                             ccopy_(&i__4, &work[iw], &c__1, &a[ioff], &c__1);
00402                         }
00403                     }
00404 
00405 /*                 For types 2-4, zero one row and column of the matrix */
00406 /*                 to test that INFO is returned correctly. */
00407 
00408                     izero = 0;
00409                     if (zerot) {
00410                         if (imat == 2) {
00411                             izero = 1;
00412                         } else if (imat == 3) {
00413                             izero = n;
00414                         } else {
00415                             izero = n / 2 + 1;
00416                         }
00417 
00418 /*                    Save the zeroed out row and column in WORK(*,3) */
00419 
00420                         iw = lda << 1;
00421 /* Computing MIN */
00422                         i__5 = (kd << 1) + 1;
00423                         i__4 = min(i__5,n);
00424                         for (i__ = 1; i__ <= i__4; ++i__) {
00425                             i__5 = iw + i__;
00426                             work[i__5].r = 0.f, work[i__5].i = 0.f;
00427 /* L20: */
00428                         }
00429                         ++iw;
00430 /* Computing MAX */
00431                         i__4 = izero - kd;
00432                         i1 = max(i__4,1);
00433 /* Computing MIN */
00434                         i__4 = izero + kd;
00435                         i2 = min(i__4,n);
00436 
00437                         if (iuplo == 1) {
00438                             ioff = (izero - 1) * ldab + kd + 1;
00439                             i__4 = izero - i1;
00440                             cswap_(&i__4, &a[ioff - izero + i1], &c__1, &work[
00441                                     iw], &c__1);
00442                             iw = iw + izero - i1;
00443                             i__4 = i2 - izero + 1;
00444 /* Computing MAX */
00445                             i__6 = ldab - 1;
00446                             i__5 = max(i__6,1);
00447                             cswap_(&i__4, &a[ioff], &i__5, &work[iw], &c__1);
00448                         } else {
00449                             ioff = (i1 - 1) * ldab + 1;
00450                             i__4 = izero - i1;
00451 /* Computing MAX */
00452                             i__6 = ldab - 1;
00453                             i__5 = max(i__6,1);
00454                             cswap_(&i__4, &a[ioff + izero - i1], &i__5, &work[
00455                                     iw], &c__1);
00456                             ioff = (izero - 1) * ldab + 1;
00457                             iw = iw + izero - i1;
00458                             i__4 = i2 - izero + 1;
00459                             cswap_(&i__4, &a[ioff], &c__1, &work[iw], &c__1);
00460                         }
00461                     }
00462 
00463 /*                 Set the imaginary part of the diagonals. */
00464 
00465                     if (iuplo == 1) {
00466                         claipd_(&n, &a[kd + 1], &ldab, &c__0);
00467                     } else {
00468                         claipd_(&n, &a[1], &ldab, &c__0);
00469                     }
00470 
00471 /*                 Save a copy of the matrix A in ASAV. */
00472 
00473                     i__4 = kd + 1;
00474                     clacpy_("Full", &i__4, &n, &a[1], &ldab, &asav[1], &ldab);
00475 
00476                     for (iequed = 1; iequed <= 2; ++iequed) {
00477                         *(unsigned char *)equed = *(unsigned char *)&equeds[
00478                                 iequed - 1];
00479                         if (iequed == 1) {
00480                             nfact = 3;
00481                         } else {
00482                             nfact = 1;
00483                         }
00484 
00485                         i__4 = nfact;
00486                         for (ifact = 1; ifact <= i__4; ++ifact) {
00487                             *(unsigned char *)fact = *(unsigned char *)&facts[
00488                                     ifact - 1];
00489                             prefac = lsame_(fact, "F");
00490                             nofact = lsame_(fact, "N");
00491                             equil = lsame_(fact, "E");
00492 
00493                             if (zerot) {
00494                                 if (prefac) {
00495                                     goto L60;
00496                                 }
00497                                 rcondc = 0.f;
00498 
00499                             } else if (! lsame_(fact, "N")) {
00500 
00501 /*                          Compute the condition number for comparison */
00502 /*                          with the value returned by CPBSVX (FACT = */
00503 /*                          'N' reuses the condition number from the */
00504 /*                          previous iteration with FACT = 'F'). */
00505 
00506                                 i__5 = kd + 1;
00507                                 clacpy_("Full", &i__5, &n, &asav[1], &ldab, &
00508                                         afac[1], &ldab);
00509                                 if (equil || iequed > 1) {
00510 
00511 /*                             Compute row and column scale factors to */
00512 /*                             equilibrate the matrix A. */
00513 
00514                                     cpbequ_(uplo, &n, &kd, &afac[1], &ldab, &
00515                                             s[1], &scond, &amax, &info);
00516                                     if (info == 0 && n > 0) {
00517                                         if (iequed > 1) {
00518                                             scond = 0.f;
00519                                         }
00520 
00521 /*                                Equilibrate the matrix. */
00522 
00523                                         claqhb_(uplo, &n, &kd, &afac[1], &
00524                                                 ldab, &s[1], &scond, &amax, 
00525                                                 equed);
00526                                     }
00527                                 }
00528 
00529 /*                          Save the condition number of the */
00530 /*                          non-equilibrated system for use in CGET04. */
00531 
00532                                 if (equil) {
00533                                     roldc = rcondc;
00534                                 }
00535 
00536 /*                          Compute the 1-norm of A. */
00537 
00538                                 anorm = clanhb_("1", uplo, &n, &kd, &afac[1], 
00539                                         &ldab, &rwork[1]);
00540 
00541 /*                          Factor the matrix A. */
00542 
00543                                 cpbtrf_(uplo, &n, &kd, &afac[1], &ldab, &info);
00544 
00545 /*                          Form the inverse of A. */
00546 
00547                                 claset_("Full", &n, &n, &c_b47, &c_b48, &a[1], 
00548                                          &lda);
00549                                 s_copy(srnamc_1.srnamt, "CPBTRS", (ftnlen)32, 
00550                                         (ftnlen)6);
00551                                 cpbtrs_(uplo, &n, &kd, &n, &afac[1], &ldab, &
00552                                         a[1], &lda, &info);
00553 
00554 /*                          Compute the 1-norm condition number of A. */
00555 
00556                                 ainvnm = clange_("1", &n, &n, &a[1], &lda, &
00557                                         rwork[1]);
00558                                 if (anorm <= 0.f || ainvnm <= 0.f) {
00559                                     rcondc = 1.f;
00560                                 } else {
00561                                     rcondc = 1.f / anorm / ainvnm;
00562                                 }
00563                             }
00564 
00565 /*                       Restore the matrix A. */
00566 
00567                             i__5 = kd + 1;
00568                             clacpy_("Full", &i__5, &n, &asav[1], &ldab, &a[1], 
00569                                      &ldab);
00570 
00571 /*                       Form an exact solution and set the right hand */
00572 /*                       side. */
00573 
00574                             s_copy(srnamc_1.srnamt, "CLARHS", (ftnlen)32, (
00575                                     ftnlen)6);
00576                             clarhs_(path, xtype, uplo, " ", &n, &n, &kd, &kd, 
00577                                     nrhs, &a[1], &ldab, &xact[1], &lda, &b[1], 
00578                                      &lda, iseed, &info);
00579                             *(unsigned char *)xtype = 'C';
00580                             clacpy_("Full", &n, nrhs, &b[1], &lda, &bsav[1], &
00581                                     lda);
00582 
00583                             if (nofact) {
00584 
00585 /*                          --- Test CPBSV  --- */
00586 
00587 /*                          Compute the L*L' or U'*U factorization of the */
00588 /*                          matrix and solve the system. */
00589 
00590                                 i__5 = kd + 1;
00591                                 clacpy_("Full", &i__5, &n, &a[1], &ldab, &
00592                                         afac[1], &ldab);
00593                                 clacpy_("Full", &n, nrhs, &b[1], &lda, &x[1], 
00594                                         &lda);
00595 
00596                                 s_copy(srnamc_1.srnamt, "CPBSV ", (ftnlen)32, 
00597                                         (ftnlen)6);
00598                                 cpbsv_(uplo, &n, &kd, nrhs, &afac[1], &ldab, &
00599                                         x[1], &lda, &info);
00600 
00601 /*                          Check error code from CPBSV . */
00602 
00603                                 if (info != izero) {
00604                                     alaerh_(path, "CPBSV ", &info, &izero, 
00605                                             uplo, &n, &n, &kd, &kd, nrhs, &
00606                                             imat, &nfail, &nerrs, nout);
00607                                     goto L40;
00608                                 } else if (info != 0) {
00609                                     goto L40;
00610                                 }
00611 
00612 /*                          Reconstruct matrix from factors and compute */
00613 /*                          residual. */
00614 
00615                                 cpbt01_(uplo, &n, &kd, &a[1], &ldab, &afac[1], 
00616                                          &ldab, &rwork[1], result);
00617 
00618 /*                          Compute residual of the computed solution. */
00619 
00620                                 clacpy_("Full", &n, nrhs, &b[1], &lda, &work[
00621                                         1], &lda);
00622                                 cpbt02_(uplo, &n, &kd, nrhs, &a[1], &ldab, &x[
00623                                         1], &lda, &work[1], &lda, &rwork[1], &
00624                                         result[1]);
00625 
00626 /*                          Check solution from generated exact solution. */
00627 
00628                                 cget04_(&n, nrhs, &x[1], &lda, &xact[1], &lda, 
00629                                          &rcondc, &result[2]);
00630                                 nt = 3;
00631 
00632 /*                          Print information about the tests that did */
00633 /*                          not pass the threshold. */
00634 
00635                                 i__5 = nt;
00636                                 for (k = 1; k <= i__5; ++k) {
00637                                     if (result[k - 1] >= *thresh) {
00638                                         if (nfail == 0 && nerrs == 0) {
00639                                             aladhd_(nout, path);
00640                                         }
00641                                         io___57.ciunit = *nout;
00642                                         s_wsfe(&io___57);
00643                                         do_fio(&c__1, "CPBSV ", (ftnlen)6);
00644                                         do_fio(&c__1, uplo, (ftnlen)1);
00645                                         do_fio(&c__1, (char *)&n, (ftnlen)
00646                                                 sizeof(integer));
00647                                         do_fio(&c__1, (char *)&kd, (ftnlen)
00648                                                 sizeof(integer));
00649                                         do_fio(&c__1, (char *)&imat, (ftnlen)
00650                                                 sizeof(integer));
00651                                         do_fio(&c__1, (char *)&k, (ftnlen)
00652                                                 sizeof(integer));
00653                                         do_fio(&c__1, (char *)&result[k - 1], 
00654                                                 (ftnlen)sizeof(real));
00655                                         e_wsfe();
00656                                         ++nfail;
00657                                     }
00658 /* L30: */
00659                                 }
00660                                 nrun += nt;
00661 L40:
00662                                 ;
00663                             }
00664 
00665 /*                       --- Test CPBSVX --- */
00666 
00667                             if (! prefac) {
00668                                 i__5 = kd + 1;
00669                                 claset_("Full", &i__5, &n, &c_b47, &c_b47, &
00670                                         afac[1], &ldab);
00671                             }
00672                             claset_("Full", &n, nrhs, &c_b47, &c_b47, &x[1], &
00673                                     lda);
00674                             if (iequed > 1 && n > 0) {
00675 
00676 /*                          Equilibrate the matrix if FACT='F' and */
00677 /*                          EQUED='Y' */
00678 
00679                                 claqhb_(uplo, &n, &kd, &a[1], &ldab, &s[1], &
00680                                         scond, &amax, equed);
00681                             }
00682 
00683 /*                       Solve the system and compute the condition */
00684 /*                       number and error bounds using CPBSVX. */
00685 
00686                             s_copy(srnamc_1.srnamt, "CPBSVX", (ftnlen)32, (
00687                                     ftnlen)6);
00688                             cpbsvx_(fact, uplo, &n, &kd, nrhs, &a[1], &ldab, &
00689                                     afac[1], &ldab, equed, &s[1], &b[1], &lda, 
00690                                      &x[1], &lda, &rcond, &rwork[1], &rwork[*
00691                                     nrhs + 1], &work[1], &rwork[(*nrhs << 1) 
00692                                     + 1], &info);
00693 
00694 /*                       Check the error code from CPBSVX. */
00695 
00696                             if (info != izero) {
00697 /* Writing concatenation */
00698                                 i__7[0] = 1, a__1[0] = fact;
00699                                 i__7[1] = 1, a__1[1] = uplo;
00700                                 s_cat(ch__1, a__1, i__7, &c__2, (ftnlen)2);
00701                                 alaerh_(path, "CPBSVX", &info, &izero, ch__1, 
00702                                         &n, &n, &kd, &kd, nrhs, &imat, &nfail, 
00703                                          &nerrs, nout);
00704                                 goto L60;
00705                             }
00706 
00707                             if (info == 0) {
00708                                 if (! prefac) {
00709 
00710 /*                             Reconstruct matrix from factors and */
00711 /*                             compute residual. */
00712 
00713                                     cpbt01_(uplo, &n, &kd, &a[1], &ldab, &
00714                                             afac[1], &ldab, &rwork[(*nrhs << 
00715                                             1) + 1], result);
00716                                     k1 = 1;
00717                                 } else {
00718                                     k1 = 2;
00719                                 }
00720 
00721 /*                          Compute residual of the computed solution. */
00722 
00723                                 clacpy_("Full", &n, nrhs, &bsav[1], &lda, &
00724                                         work[1], &lda);
00725                                 cpbt02_(uplo, &n, &kd, nrhs, &asav[1], &ldab, 
00726                                         &x[1], &lda, &work[1], &lda, &rwork[(*
00727                                         nrhs << 1) + 1], &result[1]);
00728 
00729 /*                          Check solution from generated exact solution. */
00730 
00731                                 if (nofact || prefac && lsame_(equed, "N")) {
00732                                     cget04_(&n, nrhs, &x[1], &lda, &xact[1], &
00733                                             lda, &rcondc, &result[2]);
00734                                 } else {
00735                                     cget04_(&n, nrhs, &x[1], &lda, &xact[1], &
00736                                             lda, &roldc, &result[2]);
00737                                 }
00738 
00739 /*                          Check the error bounds from iterative */
00740 /*                          refinement. */
00741 
00742                                 cpbt05_(uplo, &n, &kd, nrhs, &asav[1], &ldab, 
00743                                         &b[1], &lda, &x[1], &lda, &xact[1], &
00744                                         lda, &rwork[1], &rwork[*nrhs + 1], &
00745                                         result[3]);
00746                             } else {
00747                                 k1 = 6;
00748                             }
00749 
00750 /*                       Compare RCOND from CPBSVX with the computed */
00751 /*                       value in RCONDC. */
00752 
00753                             result[5] = sget06_(&rcond, &rcondc);
00754 
00755 /*                       Print information about the tests that did not */
00756 /*                       pass the threshold. */
00757 
00758                             for (k = k1; k <= 6; ++k) {
00759                                 if (result[k - 1] >= *thresh) {
00760                                     if (nfail == 0 && nerrs == 0) {
00761                                         aladhd_(nout, path);
00762                                     }
00763                                     if (prefac) {
00764                                         io___60.ciunit = *nout;
00765                                         s_wsfe(&io___60);
00766                                         do_fio(&c__1, "CPBSVX", (ftnlen)6);
00767                                         do_fio(&c__1, fact, (ftnlen)1);
00768                                         do_fio(&c__1, uplo, (ftnlen)1);
00769                                         do_fio(&c__1, (char *)&n, (ftnlen)
00770                                                 sizeof(integer));
00771                                         do_fio(&c__1, (char *)&kd, (ftnlen)
00772                                                 sizeof(integer));
00773                                         do_fio(&c__1, equed, (ftnlen)1);
00774                                         do_fio(&c__1, (char *)&imat, (ftnlen)
00775                                                 sizeof(integer));
00776                                         do_fio(&c__1, (char *)&k, (ftnlen)
00777                                                 sizeof(integer));
00778                                         do_fio(&c__1, (char *)&result[k - 1], 
00779                                                 (ftnlen)sizeof(real));
00780                                         e_wsfe();
00781                                     } else {
00782                                         io___61.ciunit = *nout;
00783                                         s_wsfe(&io___61);
00784                                         do_fio(&c__1, "CPBSVX", (ftnlen)6);
00785                                         do_fio(&c__1, fact, (ftnlen)1);
00786                                         do_fio(&c__1, uplo, (ftnlen)1);
00787                                         do_fio(&c__1, (char *)&n, (ftnlen)
00788                                                 sizeof(integer));
00789                                         do_fio(&c__1, (char *)&kd, (ftnlen)
00790                                                 sizeof(integer));
00791                                         do_fio(&c__1, (char *)&imat, (ftnlen)
00792                                                 sizeof(integer));
00793                                         do_fio(&c__1, (char *)&k, (ftnlen)
00794                                                 sizeof(integer));
00795                                         do_fio(&c__1, (char *)&result[k - 1], 
00796                                                 (ftnlen)sizeof(real));
00797                                         e_wsfe();
00798                                     }
00799                                     ++nfail;
00800                                 }
00801 /* L50: */
00802                             }
00803                             nrun = nrun + 7 - k1;
00804 L60:
00805                             ;
00806                         }
00807 /* L70: */
00808                     }
00809 L80:
00810                     ;
00811                 }
00812 /* L90: */
00813             }
00814 /* L100: */
00815         }
00816 /* L110: */
00817     }
00818 
00819 /*     Print a summary of the results. */
00820 
00821     alasvm_(path, nout, &nfail, &nrun, &nerrs);
00822 
00823     return 0;
00824 
00825 /*     End of CDRVPB */
00826 
00827 } /* cdrvpb_ */


swiftnav
Author(s):
autogenerated on Sat Jun 8 2019 18:55:21