LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
dchkee.f
Go to the documentation of this file.
00001 *> \brief \b DCHKEE
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *  Definition:
00009 *  ===========
00010 *
00011 *       PROGRAM DCHKEE
00012 * 
00013 *
00014 *> \par Purpose:
00015 *  =============
00016 *>
00017 *> \verbatim
00018 *>
00019 *> DCHKEE tests the DOUBLE PRECISION LAPACK subroutines for the matrix
00020 *> eigenvalue problem.  The test paths in this version are
00021 *>
00022 *> NEP (Nonsymmetric Eigenvalue Problem):
00023 *>     Test DGEHRD, DORGHR, DHSEQR, DTREVC, DHSEIN, and DORMHR
00024 *>
00025 *> SEP (Symmetric Eigenvalue Problem):
00026 *>     Test DSYTRD, DORGTR, DSTEQR, DSTERF, DSTEIN, DSTEDC,
00027 *>     and drivers DSYEV(X), DSBEV(X), DSPEV(X), DSTEV(X),
00028 *>                 DSYEVD,   DSBEVD,   DSPEVD,   DSTEVD
00029 *>
00030 *> SVD (Singular Value Decomposition):
00031 *>     Test DGEBRD, DORGBR, DBDSQR, DBDSDC
00032 *>     and the drivers DGESVD, DGESDD
00033 *>
00034 *> DEV (Nonsymmetric Eigenvalue/eigenvector Driver):
00035 *>     Test DGEEV
00036 *>
00037 *> DES (Nonsymmetric Schur form Driver):
00038 *>     Test DGEES
00039 *>
00040 *> DVX (Nonsymmetric Eigenvalue/eigenvector Expert Driver):
00041 *>     Test DGEEVX
00042 *>
00043 *> DSX (Nonsymmetric Schur form Expert Driver):
00044 *>     Test DGEESX
00045 *>
00046 *> DGG (Generalized Nonsymmetric Eigenvalue Problem):
00047 *>     Test DGGHRD, DGGBAL, DGGBAK, DHGEQZ, and DTGEVC
00048 *>     and the driver routines DGEGS and DGEGV
00049 *>
00050 *> DGS (Generalized Nonsymmetric Schur form Driver):
00051 *>     Test DGGES
00052 *>
00053 *> DGV (Generalized Nonsymmetric Eigenvalue/eigenvector Driver):
00054 *>     Test DGGEV
00055 *>
00056 *> DGX (Generalized Nonsymmetric Schur form Expert Driver):
00057 *>     Test DGGESX
00058 *>
00059 *> DXV (Generalized Nonsymmetric Eigenvalue/eigenvector Expert Driver):
00060 *>     Test DGGEVX
00061 *>
00062 *> DSG (Symmetric Generalized Eigenvalue Problem):
00063 *>     Test DSYGST, DSYGV, DSYGVD, DSYGVX, DSPGST, DSPGV, DSPGVD,
00064 *>     DSPGVX, DSBGST, DSBGV, DSBGVD, and DSBGVX
00065 *>
00066 *> DSB (Symmetric Band Eigenvalue Problem):
00067 *>     Test DSBTRD
00068 *>
00069 *> DBB (Band Singular Value Decomposition):
00070 *>     Test DGBBRD
00071 *>
00072 *> DEC (Eigencondition estimation):
00073 *>     Test DLALN2, DLASY2, DLAEQU, DLAEXC, DTRSYL, DTREXC, DTRSNA,
00074 *>     DTRSEN, and DLAQTR
00075 *>
00076 *> DBL (Balancing a general matrix)
00077 *>     Test DGEBAL
00078 *>
00079 *> DBK (Back transformation on a balanced matrix)
00080 *>     Test DGEBAK
00081 *>
00082 *> DGL (Balancing a matrix pair)
00083 *>     Test DGGBAL
00084 *>
00085 *> DGK (Back transformation on a matrix pair)
00086 *>     Test DGGBAK
00087 *>
00088 *> GLM (Generalized Linear Regression Model):
00089 *>     Tests DGGGLM
00090 *>
00091 *> GQR (Generalized QR and RQ factorizations):
00092 *>     Tests DGGQRF and DGGRQF
00093 *>
00094 *> GSV (Generalized Singular Value Decomposition):
00095 *>     Tests DGGSVD, DGGSVP, DTGSJA, DLAGS2, DLAPLL, and DLAPMT
00096 *>
00097 *> CSD (CS decomposition):
00098 *>     Tests DORCSD
00099 *>
00100 *> LSE (Constrained Linear Least Squares):
00101 *>     Tests DGGLSE
00102 *>
00103 *> Each test path has a different set of inputs, but the data sets for
00104 *> the driver routines xEV, xES, xVX, and xSX can be concatenated in a
00105 *> single input file.  The first line of input should contain one of the
00106 *> 3-character path names in columns 1-3.  The number of remaining lines
00107 *> depends on what is found on the first line.
00108 *>
00109 *> The number of matrix types used in testing is often controllable from
00110 *> the input file.  The number of matrix types for each path, and the
00111 *> test routine that describes them, is as follows:
00112 *>
00113 *> Path name(s)  Types    Test routine
00114 *>
00115 *> DHS or NEP      21     DCHKHS
00116 *> DST or SEP      21     DCHKST (routines)
00117 *>                 18     DDRVST (drivers)
00118 *> DBD or SVD      16     DCHKBD (routines)
00119 *>                  5     DDRVBD (drivers)
00120 *> DEV             21     DDRVEV
00121 *> DES             21     DDRVES
00122 *> DVX             21     DDRVVX
00123 *> DSX             21     DDRVSX
00124 *> DGG             26     DCHKGG (routines)
00125 *>                 26     DDRVGG (drivers)
00126 *> DGS             26     DDRGES
00127 *> DGX              5     DDRGSX
00128 *> DGV             26     DDRGEV
00129 *> DXV              2     DDRGVX
00130 *> DSG             21     DDRVSG
00131 *> DSB             15     DCHKSB
00132 *> DBB             15     DCHKBB
00133 *> DEC              -     DCHKEC
00134 *> DBL              -     DCHKBL
00135 *> DBK              -     DCHKBK
00136 *> DGL              -     DCHKGL
00137 *> DGK              -     DCHKGK
00138 *> GLM              8     DCKGLM
00139 *> GQR              8     DCKGQR
00140 *> GSV              8     DCKGSV
00141 *> CSD              3     DCKCSD
00142 *> LSE              8     DCKLSE
00143 *>
00144 *>-----------------------------------------------------------------------
00145 *>
00146 *> NEP input file:
00147 *>
00148 *> line 2:  NN, INTEGER
00149 *>          Number of values of N.
00150 *>
00151 *> line 3:  NVAL, INTEGER array, dimension (NN)
00152 *>          The values for the matrix dimension N.
00153 *>
00154 *> line 4:  NPARMS, INTEGER
00155 *>          Number of values of the parameters NB, NBMIN, NX, NS, and
00156 *>          MAXB.
00157 *>
00158 *> line 5:  NBVAL, INTEGER array, dimension (NPARMS)
00159 *>          The values for the blocksize NB.
00160 *>
00161 *> line 6:  NBMIN, INTEGER array, dimension (NPARMS)
00162 *>          The values for the minimum blocksize NBMIN.
00163 *>
00164 *> line 7:  NXVAL, INTEGER array, dimension (NPARMS)
00165 *>          The values for the crossover point NX.
00166 *>
00167 *> line 8:  INMIN, INTEGER array, dimension (NPARMS)
00168 *>          LAHQR vs TTQRE crossover point, >= 11
00169 *>
00170 *> line 9:  INWIN, INTEGER array, dimension (NPARMS)
00171 *>          recommended deflation window size
00172 *>
00173 *> line 10: INIBL, INTEGER array, dimension (NPARMS)
00174 *>          nibble crossover point
00175 *>
00176 *> line 11: ISHFTS, INTEGER array, dimension (NPARMS)
00177 *>          number of simultaneous shifts)
00178 *>
00179 *> line 12: IACC22, INTEGER array, dimension (NPARMS)
00180 *>          select structured matrix multiply: 0, 1 or 2)
00181 *>
00182 *> line 13: THRESH
00183 *>          Threshold value for the test ratios.  Information will be
00184 *>          printed about each test for which the test ratio is greater
00185 *>          than or equal to the threshold.  To have all of the test
00186 *>          ratios printed, use THRESH = 0.0 .
00187 *>
00188 *> line 14: NEWSD, INTEGER
00189 *>          A code indicating how to set the random number seed.
00190 *>          = 0:  Set the seed to a default value before each run
00191 *>          = 1:  Initialize the seed to a default value only before the
00192 *>                first run
00193 *>          = 2:  Like 1, but use the seed values on the next line
00194 *>
00195 *> If line 14 was 2:
00196 *>
00197 *> line 15: INTEGER array, dimension (4)
00198 *>          Four integer values for the random number seed.
00199 *>
00200 *> lines 15-EOF:  The remaining lines occur in sets of 1 or 2 and allow
00201 *>          the user to specify the matrix types.  Each line contains
00202 *>          a 3-character path name in columns 1-3, and the number
00203 *>          of matrix types must be the first nonblank item in columns
00204 *>          4-80.  If the number of matrix types is at least 1 but is
00205 *>          less than the maximum number of possible types, a second
00206 *>          line will be read to get the numbers of the matrix types to
00207 *>          be used.  For example,
00208 *> NEP 21
00209 *>          requests all of the matrix types for the nonsymmetric
00210 *>          eigenvalue problem, while
00211 *> NEP  4
00212 *> 9 10 11 12
00213 *>          requests only matrices of type 9, 10, 11, and 12.
00214 *>
00215 *>          The valid 3-character path names are 'NEP' or 'SHS' for the
00216 *>          nonsymmetric eigenvalue routines.
00217 *>
00218 *>-----------------------------------------------------------------------
00219 *>
00220 *> SEP or DSG input file:
00221 *>
00222 *> line 2:  NN, INTEGER
00223 *>          Number of values of N.
00224 *>
00225 *> line 3:  NVAL, INTEGER array, dimension (NN)
00226 *>          The values for the matrix dimension N.
00227 *>
00228 *> line 4:  NPARMS, INTEGER
00229 *>          Number of values of the parameters NB, NBMIN, and NX.
00230 *>
00231 *> line 5:  NBVAL, INTEGER array, dimension (NPARMS)
00232 *>          The values for the blocksize NB.
00233 *>
00234 *> line 6:  NBMIN, INTEGER array, dimension (NPARMS)
00235 *>          The values for the minimum blocksize NBMIN.
00236 *>
00237 *> line 7:  NXVAL, INTEGER array, dimension (NPARMS)
00238 *>          The values for the crossover point NX.
00239 *>
00240 *> line 8:  THRESH
00241 *>          Threshold value for the test ratios.  Information will be
00242 *>          printed about each test for which the test ratio is greater
00243 *>          than or equal to the threshold.
00244 *>
00245 *> line 9:  TSTCHK, LOGICAL
00246 *>          Flag indicating whether or not to test the LAPACK routines.
00247 *>
00248 *> line 10: TSTDRV, LOGICAL
00249 *>          Flag indicating whether or not to test the driver routines.
00250 *>
00251 *> line 11: TSTERR, LOGICAL
00252 *>          Flag indicating whether or not to test the error exits for
00253 *>          the LAPACK routines and driver routines.
00254 *>
00255 *> line 12: NEWSD, INTEGER
00256 *>          A code indicating how to set the random number seed.
00257 *>          = 0:  Set the seed to a default value before each run
00258 *>          = 1:  Initialize the seed to a default value only before the
00259 *>                first run
00260 *>          = 2:  Like 1, but use the seed values on the next line
00261 *>
00262 *> If line 12 was 2:
00263 *>
00264 *> line 13: INTEGER array, dimension (4)
00265 *>          Four integer values for the random number seed.
00266 *>
00267 *> lines 13-EOF:  Lines specifying matrix types, as for NEP.
00268 *>          The 3-character path names are 'SEP' or 'SST' for the
00269 *>          symmetric eigenvalue routines and driver routines, and
00270 *>          'DSG' for the routines for the symmetric generalized
00271 *>          eigenvalue problem.
00272 *>
00273 *>-----------------------------------------------------------------------
00274 *>
00275 *> SVD input file:
00276 *>
00277 *> line 2:  NN, INTEGER
00278 *>          Number of values of M and N.
00279 *>
00280 *> line 3:  MVAL, INTEGER array, dimension (NN)
00281 *>          The values for the matrix row dimension M.
00282 *>
00283 *> line 4:  NVAL, INTEGER array, dimension (NN)
00284 *>          The values for the matrix column dimension N.
00285 *>
00286 *> line 5:  NPARMS, INTEGER
00287 *>          Number of values of the parameter NB, NBMIN, NX, and NRHS.
00288 *>
00289 *> line 6:  NBVAL, INTEGER array, dimension (NPARMS)
00290 *>          The values for the blocksize NB.
00291 *>
00292 *> line 7:  NBMIN, INTEGER array, dimension (NPARMS)
00293 *>          The values for the minimum blocksize NBMIN.
00294 *>
00295 *> line 8:  NXVAL, INTEGER array, dimension (NPARMS)
00296 *>          The values for the crossover point NX.
00297 *>
00298 *> line 9:  NSVAL, INTEGER array, dimension (NPARMS)
00299 *>          The values for the number of right hand sides NRHS.
00300 *>
00301 *> line 10: THRESH
00302 *>          Threshold value for the test ratios.  Information will be
00303 *>          printed about each test for which the test ratio is greater
00304 *>          than or equal to the threshold.
00305 *>
00306 *> line 11: TSTCHK, LOGICAL
00307 *>          Flag indicating whether or not to test the LAPACK routines.
00308 *>
00309 *> line 12: TSTDRV, LOGICAL
00310 *>          Flag indicating whether or not to test the driver routines.
00311 *>
00312 *> line 13: TSTERR, LOGICAL
00313 *>          Flag indicating whether or not to test the error exits for
00314 *>          the LAPACK routines and driver routines.
00315 *>
00316 *> line 14: NEWSD, INTEGER
00317 *>          A code indicating how to set the random number seed.
00318 *>          = 0:  Set the seed to a default value before each run
00319 *>          = 1:  Initialize the seed to a default value only before the
00320 *>                first run
00321 *>          = 2:  Like 1, but use the seed values on the next line
00322 *>
00323 *> If line 14 was 2:
00324 *>
00325 *> line 15: INTEGER array, dimension (4)
00326 *>          Four integer values for the random number seed.
00327 *>
00328 *> lines 15-EOF:  Lines specifying matrix types, as for NEP.
00329 *>          The 3-character path names are 'SVD' or 'SBD' for both the
00330 *>          SVD routines and the SVD driver routines.
00331 *>
00332 *>-----------------------------------------------------------------------
00333 *>
00334 *> DEV and DES data files:
00335 *>
00336 *> line 1:  'DEV' or 'DES' in columns 1 to 3.
00337 *>
00338 *> line 2:  NSIZES, INTEGER
00339 *>          Number of sizes of matrices to use. Should be at least 0
00340 *>          and at most 20. If NSIZES = 0, no testing is done
00341 *>          (although the remaining  3 lines are still read).
00342 *>
00343 *> line 3:  NN, INTEGER array, dimension(NSIZES)
00344 *>          Dimensions of matrices to be tested.
00345 *>
00346 *> line 4:  NB, NBMIN, NX, NS, NBCOL, INTEGERs
00347 *>          These integer parameters determine how blocking is done
00348 *>          (see ILAENV for details)
00349 *>          NB     : block size
00350 *>          NBMIN  : minimum block size
00351 *>          NX     : minimum dimension for blocking
00352 *>          NS     : number of shifts in xHSEQR
00353 *>          NBCOL  : minimum column dimension for blocking
00354 *>
00355 *> line 5:  THRESH, REAL
00356 *>          The test threshold against which computed residuals are
00357 *>          compared. Should generally be in the range from 10. to 20.
00358 *>          If it is 0., all test case data will be printed.
00359 *>
00360 *> line 6:  TSTERR, LOGICAL
00361 *>          Flag indicating whether or not to test the error exits.
00362 *>
00363 *> line 7:  NEWSD, INTEGER
00364 *>          A code indicating how to set the random number seed.
00365 *>          = 0:  Set the seed to a default value before each run
00366 *>          = 1:  Initialize the seed to a default value only before the
00367 *>                first run
00368 *>          = 2:  Like 1, but use the seed values on the next line
00369 *>
00370 *> If line 7 was 2:
00371 *>
00372 *> line 8:  INTEGER array, dimension (4)
00373 *>          Four integer values for the random number seed.
00374 *>
00375 *> lines 9 and following:  Lines specifying matrix types, as for NEP.
00376 *>          The 3-character path name is 'DEV' to test SGEEV, or
00377 *>          'DES' to test SGEES.
00378 *>
00379 *>-----------------------------------------------------------------------
00380 *>
00381 *> The DVX data has two parts. The first part is identical to DEV,
00382 *> and the second part consists of test matrices with precomputed
00383 *> solutions.
00384 *>
00385 *> line 1:  'DVX' in columns 1-3.
00386 *>
00387 *> line 2:  NSIZES, INTEGER
00388 *>          If NSIZES = 0, no testing of randomly generated examples
00389 *>          is done, but any precomputed examples are tested.
00390 *>
00391 *> line 3:  NN, INTEGER array, dimension(NSIZES)
00392 *>
00393 *> line 4:  NB, NBMIN, NX, NS, NBCOL, INTEGERs
00394 *>
00395 *> line 5:  THRESH, REAL
00396 *>
00397 *> line 6:  TSTERR, LOGICAL
00398 *>
00399 *> line 7:  NEWSD, INTEGER
00400 *>
00401 *> If line 7 was 2:
00402 *>
00403 *> line 8:  INTEGER array, dimension (4)
00404 *>
00405 *> lines 9 and following: The first line contains 'DVX' in columns 1-3
00406 *>          followed by the number of matrix types, possibly with
00407 *>          a second line to specify certain matrix types.
00408 *>          If the number of matrix types = 0, no testing of randomly
00409 *>          generated examples is done, but any precomputed examples
00410 *>          are tested.
00411 *>
00412 *> remaining lines : Each matrix is stored on 1+2*N lines, where N is
00413 *>          its dimension. The first line contains the dimension (a
00414 *>          single integer). The next N lines contain the matrix, one
00415 *>          row per line. The last N lines correspond to each
00416 *>          eigenvalue. Each of these last N lines contains 4 real
00417 *>          values: the real part of the eigenvalue, the imaginary
00418 *>          part of the eigenvalue, the reciprocal condition number of
00419 *>          the eigenvalues, and the reciprocal condition number of the
00420 *>          eigenvector.  The end of data is indicated by dimension N=0.
00421 *>          Even if no data is to be tested, there must be at least one
00422 *>          line containing N=0.
00423 *>
00424 *>-----------------------------------------------------------------------
00425 *>
00426 *> The DSX data is like DVX. The first part is identical to DEV, and the
00427 *> second part consists of test matrices with precomputed solutions.
00428 *>
00429 *> line 1:  'DSX' in columns 1-3.
00430 *>
00431 *> line 2:  NSIZES, INTEGER
00432 *>          If NSIZES = 0, no testing of randomly generated examples
00433 *>          is done, but any precomputed examples are tested.
00434 *>
00435 *> line 3:  NN, INTEGER array, dimension(NSIZES)
00436 *>
00437 *> line 4:  NB, NBMIN, NX, NS, NBCOL, INTEGERs
00438 *>
00439 *> line 5:  THRESH, REAL
00440 *>
00441 *> line 6:  TSTERR, LOGICAL
00442 *>
00443 *> line 7:  NEWSD, INTEGER
00444 *>
00445 *> If line 7 was 2:
00446 *>
00447 *> line 8:  INTEGER array, dimension (4)
00448 *>
00449 *> lines 9 and following: The first line contains 'DSX' in columns 1-3
00450 *>          followed by the number of matrix types, possibly with
00451 *>          a second line to specify certain matrix types.
00452 *>          If the number of matrix types = 0, no testing of randomly
00453 *>          generated examples is done, but any precomputed examples
00454 *>          are tested.
00455 *>
00456 *> remaining lines : Each matrix is stored on 3+N lines, where N is its
00457 *>          dimension. The first line contains the dimension N and the
00458 *>          dimension M of an invariant subspace. The second line
00459 *>          contains M integers, identifying the eigenvalues in the
00460 *>          invariant subspace (by their position in a list of
00461 *>          eigenvalues ordered by increasing real part). The next N
00462 *>          lines contain the matrix. The last line contains the
00463 *>          reciprocal condition number for the average of the selected
00464 *>          eigenvalues, and the reciprocal condition number for the
00465 *>          corresponding right invariant subspace. The end of data is
00466 *>          indicated by a line containing N=0 and M=0. Even if no data
00467 *>          is to be tested, there must be at least one line containing
00468 *>          N=0 and M=0.
00469 *>
00470 *>-----------------------------------------------------------------------
00471 *>
00472 *> DGG input file:
00473 *>
00474 *> line 2:  NN, INTEGER
00475 *>          Number of values of N.
00476 *>
00477 *> line 3:  NVAL, INTEGER array, dimension (NN)
00478 *>          The values for the matrix dimension N.
00479 *>
00480 *> line 4:  NPARMS, INTEGER
00481 *>          Number of values of the parameters NB, NBMIN, NS, MAXB, and
00482 *>          NBCOL.
00483 *>
00484 *> line 5:  NBVAL, INTEGER array, dimension (NPARMS)
00485 *>          The values for the blocksize NB.
00486 *>
00487 *> line 6:  NBMIN, INTEGER array, dimension (NPARMS)
00488 *>          The values for NBMIN, the minimum row dimension for blocks.
00489 *>
00490 *> line 7:  NSVAL, INTEGER array, dimension (NPARMS)
00491 *>          The values for the number of shifts.
00492 *>
00493 *> line 8:  MXBVAL, INTEGER array, dimension (NPARMS)
00494 *>          The values for MAXB, used in determining minimum blocksize.
00495 *>
00496 *> line 9:  NBCOL, INTEGER array, dimension (NPARMS)
00497 *>          The values for NBCOL, the minimum column dimension for
00498 *>          blocks.
00499 *>
00500 *> line 10: THRESH
00501 *>          Threshold value for the test ratios.  Information will be
00502 *>          printed about each test for which the test ratio is greater
00503 *>          than or equal to the threshold.
00504 *>
00505 *> line 11: TSTCHK, LOGICAL
00506 *>          Flag indicating whether or not to test the LAPACK routines.
00507 *>
00508 *> line 12: TSTDRV, LOGICAL
00509 *>          Flag indicating whether or not to test the driver routines.
00510 *>
00511 *> line 13: TSTERR, LOGICAL
00512 *>          Flag indicating whether or not to test the error exits for
00513 *>          the LAPACK routines and driver routines.
00514 *>
00515 *> line 14: NEWSD, INTEGER
00516 *>          A code indicating how to set the random number seed.
00517 *>          = 0:  Set the seed to a default value before each run
00518 *>          = 1:  Initialize the seed to a default value only before the
00519 *>                first run
00520 *>          = 2:  Like 1, but use the seed values on the next line
00521 *>
00522 *> If line 14 was 2:
00523 *>
00524 *> line 15: INTEGER array, dimension (4)
00525 *>          Four integer values for the random number seed.
00526 *>
00527 *> lines 15-EOF:  Lines specifying matrix types, as for NEP.
00528 *>          The 3-character path name is 'DGG' for the generalized
00529 *>          eigenvalue problem routines and driver routines.
00530 *>
00531 *>-----------------------------------------------------------------------
00532 *>
00533 *> DGS and DGV input files:
00534 *>
00535 *> line 1:  'DGS' or 'DGV' in columns 1 to 3.
00536 *>
00537 *> line 2:  NN, INTEGER
00538 *>          Number of values of N.
00539 *>
00540 *> line 3:  NVAL, INTEGER array, dimension(NN)
00541 *>          Dimensions of matrices to be tested.
00542 *>
00543 *> line 4:  NB, NBMIN, NX, NS, NBCOL, INTEGERs
00544 *>          These integer parameters determine how blocking is done
00545 *>          (see ILAENV for details)
00546 *>          NB     : block size
00547 *>          NBMIN  : minimum block size
00548 *>          NX     : minimum dimension for blocking
00549 *>          NS     : number of shifts in xHGEQR
00550 *>          NBCOL  : minimum column dimension for blocking
00551 *>
00552 *> line 5:  THRESH, REAL
00553 *>          The test threshold against which computed residuals are
00554 *>          compared. Should generally be in the range from 10. to 20.
00555 *>          If it is 0., all test case data will be printed.
00556 *>
00557 *> line 6:  TSTERR, LOGICAL
00558 *>          Flag indicating whether or not to test the error exits.
00559 *>
00560 *> line 7:  NEWSD, INTEGER
00561 *>          A code indicating how to set the random number seed.
00562 *>          = 0:  Set the seed to a default value before each run
00563 *>          = 1:  Initialize the seed to a default value only before the
00564 *>                first run
00565 *>          = 2:  Like 1, but use the seed values on the next line
00566 *>
00567 *> If line 17 was 2:
00568 *>
00569 *> line 7:  INTEGER array, dimension (4)
00570 *>          Four integer values for the random number seed.
00571 *>
00572 *> lines 7-EOF:  Lines specifying matrix types, as for NEP.
00573 *>          The 3-character path name is 'DGS' for the generalized
00574 *>          eigenvalue problem routines and driver routines.
00575 *>
00576 *>-----------------------------------------------------------------------
00577 *>
00578 *> DXV input files:
00579 *>
00580 *> line 1:  'DXV' in columns 1 to 3.
00581 *>
00582 *> line 2:  N, INTEGER
00583 *>          Value of N.
00584 *>
00585 *> line 3:  NB, NBMIN, NX, NS, NBCOL, INTEGERs
00586 *>          These integer parameters determine how blocking is done
00587 *>          (see ILAENV for details)
00588 *>          NB     : block size
00589 *>          NBMIN  : minimum block size
00590 *>          NX     : minimum dimension for blocking
00591 *>          NS     : number of shifts in xHGEQR
00592 *>          NBCOL  : minimum column dimension for blocking
00593 *>
00594 *> line 4:  THRESH, REAL
00595 *>          The test threshold against which computed residuals are
00596 *>          compared. Should generally be in the range from 10. to 20.
00597 *>          Information will be printed about each test for which the
00598 *>          test ratio is greater than or equal to the threshold.
00599 *>
00600 *> line 5:  TSTERR, LOGICAL
00601 *>          Flag indicating whether or not to test the error exits for
00602 *>          the LAPACK routines and driver routines.
00603 *>
00604 *> line 6:  NEWSD, INTEGER
00605 *>          A code indicating how to set the random number seed.
00606 *>          = 0:  Set the seed to a default value before each run
00607 *>          = 1:  Initialize the seed to a default value only before the
00608 *>                first run
00609 *>          = 2:  Like 1, but use the seed values on the next line
00610 *>
00611 *> If line 6 was 2:
00612 *>
00613 *> line 7: INTEGER array, dimension (4)
00614 *>          Four integer values for the random number seed.
00615 *>
00616 *> If line 2 was 0:
00617 *>
00618 *> line 7-EOF: Precomputed examples are tested.
00619 *>
00620 *> remaining lines : Each example is stored on 3+2*N lines, where N is
00621 *>          its dimension. The first line contains the dimension (a
00622 *>          single integer). The next N lines contain the matrix A, one
00623 *>          row per line. The next N lines contain the matrix B.  The
00624 *>          next line contains the reciprocals of the eigenvalue
00625 *>          condition numbers.  The last line contains the reciprocals of
00626 *>          the eigenvector condition numbers.  The end of data is
00627 *>          indicated by dimension N=0.  Even if no data is to be tested,
00628 *>          there must be at least one line containing N=0.
00629 *>
00630 *>-----------------------------------------------------------------------
00631 *>
00632 *> DGX input files:
00633 *>
00634 *> line 1:  'DGX' in columns 1 to 3.
00635 *>
00636 *> line 2:  N, INTEGER
00637 *>          Value of N.
00638 *>
00639 *> line 3:  NB, NBMIN, NX, NS, NBCOL, INTEGERs
00640 *>          These integer parameters determine how blocking is done
00641 *>          (see ILAENV for details)
00642 *>          NB     : block size
00643 *>          NBMIN  : minimum block size
00644 *>          NX     : minimum dimension for blocking
00645 *>          NS     : number of shifts in xHGEQR
00646 *>          NBCOL  : minimum column dimension for blocking
00647 *>
00648 *> line 4:  THRESH, REAL
00649 *>          The test threshold against which computed residuals are
00650 *>          compared. Should generally be in the range from 10. to 20.
00651 *>          Information will be printed about each test for which the
00652 *>          test ratio is greater than or equal to the threshold.
00653 *>
00654 *> line 5:  TSTERR, LOGICAL
00655 *>          Flag indicating whether or not to test the error exits for
00656 *>          the LAPACK routines and driver routines.
00657 *>
00658 *> line 6:  NEWSD, INTEGER
00659 *>          A code indicating how to set the random number seed.
00660 *>          = 0:  Set the seed to a default value before each run
00661 *>          = 1:  Initialize the seed to a default value only before the
00662 *>                first run
00663 *>          = 2:  Like 1, but use the seed values on the next line
00664 *>
00665 *> If line 6 was 2:
00666 *>
00667 *> line 7: INTEGER array, dimension (4)
00668 *>          Four integer values for the random number seed.
00669 *>
00670 *> If line 2 was 0:
00671 *>
00672 *> line 7-EOF: Precomputed examples are tested.
00673 *>
00674 *> remaining lines : Each example is stored on 3+2*N lines, where N is
00675 *>          its dimension. The first line contains the dimension (a
00676 *>          single integer).  The next line contains an integer k such
00677 *>          that only the last k eigenvalues will be selected and appear
00678 *>          in the leading diagonal blocks of $A$ and $B$. The next N
00679 *>          lines contain the matrix A, one row per line.  The next N
00680 *>          lines contain the matrix B.  The last line contains the
00681 *>          reciprocal of the eigenvalue cluster condition number and the
00682 *>          reciprocal of the deflating subspace (associated with the
00683 *>          selected eigencluster) condition number.  The end of data is
00684 *>          indicated by dimension N=0.  Even if no data is to be tested,
00685 *>          there must be at least one line containing N=0.
00686 *>
00687 *>-----------------------------------------------------------------------
00688 *>
00689 *> DSB input file:
00690 *>
00691 *> line 2:  NN, INTEGER
00692 *>          Number of values of N.
00693 *>
00694 *> line 3:  NVAL, INTEGER array, dimension (NN)
00695 *>          The values for the matrix dimension N.
00696 *>
00697 *> line 4:  NK, INTEGER
00698 *>          Number of values of K.
00699 *>
00700 *> line 5:  KVAL, INTEGER array, dimension (NK)
00701 *>          The values for the matrix dimension K.
00702 *>
00703 *> line 6:  THRESH
00704 *>          Threshold value for the test ratios.  Information will be
00705 *>          printed about each test for which the test ratio is greater
00706 *>          than or equal to the threshold.
00707 *>
00708 *> line 7:  NEWSD, INTEGER
00709 *>          A code indicating how to set the random number seed.
00710 *>          = 0:  Set the seed to a default value before each run
00711 *>          = 1:  Initialize the seed to a default value only before the
00712 *>                first run
00713 *>          = 2:  Like 1, but use the seed values on the next line
00714 *>
00715 *> If line 7 was 2:
00716 *>
00717 *> line 8:  INTEGER array, dimension (4)
00718 *>          Four integer values for the random number seed.
00719 *>
00720 *> lines 8-EOF:  Lines specifying matrix types, as for NEP.
00721 *>          The 3-character path name is 'DSB'.
00722 *>
00723 *>-----------------------------------------------------------------------
00724 *>
00725 *> DBB input file:
00726 *>
00727 *> line 2:  NN, INTEGER
00728 *>          Number of values of M and N.
00729 *>
00730 *> line 3:  MVAL, INTEGER array, dimension (NN)
00731 *>          The values for the matrix row dimension M.
00732 *>
00733 *> line 4:  NVAL, INTEGER array, dimension (NN)
00734 *>          The values for the matrix column dimension N.
00735 *>
00736 *> line 4:  NK, INTEGER
00737 *>          Number of values of K.
00738 *>
00739 *> line 5:  KVAL, INTEGER array, dimension (NK)
00740 *>          The values for the matrix bandwidth K.
00741 *>
00742 *> line 6:  NPARMS, INTEGER
00743 *>          Number of values of the parameter NRHS
00744 *>
00745 *> line 7:  NSVAL, INTEGER array, dimension (NPARMS)
00746 *>          The values for the number of right hand sides NRHS.
00747 *>
00748 *> line 8:  THRESH
00749 *>          Threshold value for the test ratios.  Information will be
00750 *>          printed about each test for which the test ratio is greater
00751 *>          than or equal to the threshold.
00752 *>
00753 *> line 9:  NEWSD, INTEGER
00754 *>          A code indicating how to set the random number seed.
00755 *>          = 0:  Set the seed to a default value before each run
00756 *>          = 1:  Initialize the seed to a default value only before the
00757 *>                first run
00758 *>          = 2:  Like 1, but use the seed values on the next line
00759 *>
00760 *> If line 9 was 2:
00761 *>
00762 *> line 10: INTEGER array, dimension (4)
00763 *>          Four integer values for the random number seed.
00764 *>
00765 *> lines 10-EOF:  Lines specifying matrix types, as for SVD.
00766 *>          The 3-character path name is 'DBB'.
00767 *>
00768 *>-----------------------------------------------------------------------
00769 *>
00770 *> DEC input file:
00771 *>
00772 *> line  2: THRESH, REAL
00773 *>          Threshold value for the test ratios.  Information will be
00774 *>          printed about each test for which the test ratio is greater
00775 *>          than or equal to the threshold.
00776 *>
00777 *> lines  3-EOF:
00778 *>
00779 *> Input for testing the eigencondition routines consists of a set of
00780 *> specially constructed test cases and their solutions.  The data
00781 *> format is not intended to be modified by the user.
00782 *>
00783 *>-----------------------------------------------------------------------
00784 *>
00785 *> DBL and DBK input files:
00786 *>
00787 *> line 1:  'DBL' in columns 1-3 to test SGEBAL, or 'DBK' in
00788 *>          columns 1-3 to test SGEBAK.
00789 *>
00790 *> The remaining lines consist of specially constructed test cases.
00791 *>
00792 *>-----------------------------------------------------------------------
00793 *>
00794 *> DGL and DGK input files:
00795 *>
00796 *> line 1:  'DGL' in columns 1-3 to test DGGBAL, or 'DGK' in
00797 *>          columns 1-3 to test DGGBAK.
00798 *>
00799 *> The remaining lines consist of specially constructed test cases.
00800 *>
00801 *>-----------------------------------------------------------------------
00802 *>
00803 *> GLM data file:
00804 *>
00805 *> line 1:  'GLM' in columns 1 to 3.
00806 *>
00807 *> line 2:  NN, INTEGER
00808 *>          Number of values of M, P, and N.
00809 *>
00810 *> line 3:  MVAL, INTEGER array, dimension(NN)
00811 *>          Values of M (row dimension).
00812 *>
00813 *> line 4:  PVAL, INTEGER array, dimension(NN)
00814 *>          Values of P (row dimension).
00815 *>
00816 *> line 5:  NVAL, INTEGER array, dimension(NN)
00817 *>          Values of N (column dimension), note M <= N <= M+P.
00818 *>
00819 *> line 6:  THRESH, REAL
00820 *>          Threshold value for the test ratios.  Information will be
00821 *>          printed about each test for which the test ratio is greater
00822 *>          than or equal to the threshold.
00823 *>
00824 *> line 7:  TSTERR, LOGICAL
00825 *>          Flag indicating whether or not to test the error exits for
00826 *>          the LAPACK routines and driver routines.
00827 *>
00828 *> line 8:  NEWSD, INTEGER
00829 *>          A code indicating how to set the random number seed.
00830 *>          = 0:  Set the seed to a default value before each run
00831 *>          = 1:  Initialize the seed to a default value only before the
00832 *>                first run
00833 *>          = 2:  Like 1, but use the seed values on the next line
00834 *>
00835 *> If line 8 was 2:
00836 *>
00837 *> line 9:  INTEGER array, dimension (4)
00838 *>          Four integer values for the random number seed.
00839 *>
00840 *> lines 9-EOF:  Lines specifying matrix types, as for NEP.
00841 *>          The 3-character path name is 'GLM' for the generalized
00842 *>          linear regression model routines.
00843 *>
00844 *>-----------------------------------------------------------------------
00845 *>
00846 *> GQR data file:
00847 *>
00848 *> line 1:  'GQR' in columns 1 to 3.
00849 *>
00850 *> line 2:  NN, INTEGER
00851 *>          Number of values of M, P, and N.
00852 *>
00853 *> line 3:  MVAL, INTEGER array, dimension(NN)
00854 *>          Values of M.
00855 *>
00856 *> line 4:  PVAL, INTEGER array, dimension(NN)
00857 *>          Values of P.
00858 *>
00859 *> line 5:  NVAL, INTEGER array, dimension(NN)
00860 *>          Values of N.
00861 *>
00862 *> line 6:  THRESH, REAL
00863 *>          Threshold value for the test ratios.  Information will be
00864 *>          printed about each test for which the test ratio is greater
00865 *>          than or equal to the threshold.
00866 *>
00867 *> line 7:  TSTERR, LOGICAL
00868 *>          Flag indicating whether or not to test the error exits for
00869 *>          the LAPACK routines and driver routines.
00870 *>
00871 *> line 8:  NEWSD, INTEGER
00872 *>          A code indicating how to set the random number seed.
00873 *>          = 0:  Set the seed to a default value before each run
00874 *>          = 1:  Initialize the seed to a default value only before the
00875 *>                first run
00876 *>          = 2:  Like 1, but use the seed values on the next line
00877 *>
00878 *> If line 8 was 2:
00879 *>
00880 *> line 9:  INTEGER array, dimension (4)
00881 *>          Four integer values for the random number seed.
00882 *>
00883 *> lines 9-EOF:  Lines specifying matrix types, as for NEP.
00884 *>          The 3-character path name is 'GQR' for the generalized
00885 *>          QR and RQ routines.
00886 *>
00887 *>-----------------------------------------------------------------------
00888 *>
00889 *> GSV data file:
00890 *>
00891 *> line 1:  'GSV' in columns 1 to 3.
00892 *>
00893 *> line 2:  NN, INTEGER
00894 *>          Number of values of M, P, and N.
00895 *>
00896 *> line 3:  MVAL, INTEGER array, dimension(NN)
00897 *>          Values of M (row dimension).
00898 *>
00899 *> line 4:  PVAL, INTEGER array, dimension(NN)
00900 *>          Values of P (row dimension).
00901 *>
00902 *> line 5:  NVAL, INTEGER array, dimension(NN)
00903 *>          Values of N (column dimension).
00904 *>
00905 *> line 6:  THRESH, REAL
00906 *>          Threshold value for the test ratios.  Information will be
00907 *>          printed about each test for which the test ratio is greater
00908 *>          than or equal to the threshold.
00909 *>
00910 *> line 7:  TSTERR, LOGICAL
00911 *>          Flag indicating whether or not to test the error exits for
00912 *>          the LAPACK routines and driver routines.
00913 *>
00914 *> line 8:  NEWSD, INTEGER
00915 *>          A code indicating how to set the random number seed.
00916 *>          = 0:  Set the seed to a default value before each run
00917 *>          = 1:  Initialize the seed to a default value only before the
00918 *>                first run
00919 *>          = 2:  Like 1, but use the seed values on the next line
00920 *>
00921 *> If line 8 was 2:
00922 *>
00923 *> line 9:  INTEGER array, dimension (4)
00924 *>          Four integer values for the random number seed.
00925 *>
00926 *> lines 9-EOF:  Lines specifying matrix types, as for NEP.
00927 *>          The 3-character path name is 'GSV' for the generalized
00928 *>          SVD routines.
00929 *>
00930 *>-----------------------------------------------------------------------
00931 *>
00932 *> CSD data file:
00933 *>
00934 *> line 1:  'CSD' in columns 1 to 3.
00935 *>
00936 *> line 2:  NM, INTEGER
00937 *>          Number of values of M, P, and N.
00938 *>
00939 *> line 3:  MVAL, INTEGER array, dimension(NM)
00940 *>          Values of M (row and column dimension of orthogonal matrix).
00941 *>
00942 *> line 4:  PVAL, INTEGER array, dimension(NM)
00943 *>          Values of P (row dimension of top-left block).
00944 *>
00945 *> line 5:  NVAL, INTEGER array, dimension(NM)
00946 *>          Values of N (column dimension of top-left block).
00947 *>
00948 *> line 6:  THRESH, REAL
00949 *>          Threshold value for the test ratios.  Information will be
00950 *>          printed about each test for which the test ratio is greater
00951 *>          than or equal to the threshold.
00952 *>
00953 *> line 7:  TSTERR, LOGICAL
00954 *>          Flag indicating whether or not to test the error exits for
00955 *>          the LAPACK routines and driver routines.
00956 *>
00957 *> line 8:  NEWSD, INTEGER
00958 *>          A code indicating how to set the random number seed.
00959 *>          = 0:  Set the seed to a default value before each run
00960 *>          = 1:  Initialize the seed to a default value only before the
00961 *>                first run
00962 *>          = 2:  Like 1, but use the seed values on the next line
00963 *>
00964 *> If line 8 was 2:
00965 *>
00966 *> line 9:  INTEGER array, dimension (4)
00967 *>          Four integer values for the random number seed.
00968 *>
00969 *> lines 9-EOF:  Lines specifying matrix types, as for NEP.
00970 *>          The 3-character path name is 'CSD' for the CSD routine.
00971 *>
00972 *>-----------------------------------------------------------------------
00973 *>
00974 *> LSE data file:
00975 *>
00976 *> line 1:  'LSE' in columns 1 to 3.
00977 *>
00978 *> line 2:  NN, INTEGER
00979 *>          Number of values of M, P, and N.
00980 *>
00981 *> line 3:  MVAL, INTEGER array, dimension(NN)
00982 *>          Values of M.
00983 *>
00984 *> line 4:  PVAL, INTEGER array, dimension(NN)
00985 *>          Values of P.
00986 *>
00987 *> line 5:  NVAL, INTEGER array, dimension(NN)
00988 *>          Values of N, note P <= N <= P+M.
00989 *>
00990 *> line 6:  THRESH, REAL
00991 *>          Threshold value for the test ratios.  Information will be
00992 *>          printed about each test for which the test ratio is greater
00993 *>          than or equal to the threshold.
00994 *>
00995 *> line 7:  TSTERR, LOGICAL
00996 *>          Flag indicating whether or not to test the error exits for
00997 *>          the LAPACK routines and driver routines.
00998 *>
00999 *> line 8:  NEWSD, INTEGER
01000 *>          A code indicating how to set the random number seed.
01001 *>          = 0:  Set the seed to a default value before each run
01002 *>          = 1:  Initialize the seed to a default value only before the
01003 *>                first run
01004 *>          = 2:  Like 1, but use the seed values on the next line
01005 *>
01006 *> If line 8 was 2:
01007 *>
01008 *> line 9:  INTEGER array, dimension (4)
01009 *>          Four integer values for the random number seed.
01010 *>
01011 *> lines 9-EOF:  Lines specifying matrix types, as for NEP.
01012 *>          The 3-character path name is 'GSV' for the generalized
01013 *>          SVD routines.
01014 *>
01015 *>-----------------------------------------------------------------------
01016 *>
01017 *> NMAX is currently set to 132 and must be at least 12 for some of the
01018 *> precomputed examples, and LWORK = NMAX*(5*NMAX+5)+1 in the parameter
01019 *> statements below.  For SVD, we assume NRHS may be as big as N.  The
01020 *> parameter NEED is set to 14 to allow for 14 N-by-N matrices for DGG.
01021 *> \endverbatim
01022 *
01023 *  Arguments:
01024 *  ==========
01025 *
01026 *
01027 *  Authors:
01028 *  ========
01029 *
01030 *> \author Univ. of Tennessee 
01031 *> \author Univ. of California Berkeley 
01032 *> \author Univ. of Colorado Denver 
01033 *> \author NAG Ltd. 
01034 *
01035 *> \date April 2012
01036 *
01037 *> \ingroup double_eig
01038 *
01039 *  =====================================================================
01040       PROGRAM DCHKEE
01041 *
01042 *  -- LAPACK test routine (version 3.4.1) --
01043 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
01044 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
01045 *     April 2012
01046 *
01047 *  =====================================================================
01048 *
01049 *     .. Parameters ..
01050       INTEGER            NMAX
01051       PARAMETER          ( NMAX = 132 )
01052       INTEGER            NCMAX
01053       PARAMETER          ( NCMAX = 20 )
01054       INTEGER            NEED
01055       PARAMETER          ( NEED = 14 )
01056       INTEGER            LWORK
01057       PARAMETER          ( LWORK = NMAX*( 5*NMAX+5 )+1 )
01058       INTEGER            LIWORK
01059       PARAMETER          ( LIWORK = NMAX*( 5*NMAX+20 ) )
01060       INTEGER            MAXIN
01061       PARAMETER          ( MAXIN = 20 )
01062       INTEGER            MAXT
01063       PARAMETER          ( MAXT = 30 )
01064       INTEGER            NIN, NOUT
01065       PARAMETER          ( NIN = 5, NOUT = 6 )
01066 *     ..
01067 *     .. Local Scalars ..
01068       LOGICAL            CSD, DBB, DGG, DSB, FATAL, GLM, GQR, GSV, LSE,
01069      $                   NEP, DBK, DBL, SEP, DES, DEV, DGK, DGL, DGS,
01070      $                   DGV, DGX, DSX, SVD, DVX, DXV, TSTCHK, TSTDIF,
01071      $                   TSTDRV, TSTERR
01072       CHARACTER          C1
01073       CHARACTER*3        C3, PATH
01074       CHARACTER*32       VNAME
01075       CHARACTER*10       INTSTR
01076       CHARACTER*80       LINE
01077       INTEGER            I, I1, IC, INFO, ITMP, K, LENP, MAXTYP, NEWSD,
01078      $                   NK, NN, NPARMS, NRHS, NTYPES,
01079      $                   VERS_MAJOR, VERS_MINOR, VERS_PATCH 
01080       DOUBLE PRECISION   EPS, S1, S2, THRESH, THRSHN
01081 *     ..
01082 *     .. Local Arrays ..
01083       LOGICAL            DOTYPE( MAXT ), LOGWRK( NMAX )
01084       INTEGER            IOLDSD( 4 ), ISEED( 4 ), IWORK( LIWORK ),
01085      $                   KVAL( MAXIN ), MVAL( MAXIN ), MXBVAL( MAXIN ),
01086      $                   NBCOL( MAXIN ), NBMIN( MAXIN ), NBVAL( MAXIN ),
01087      $                   NSVAL( MAXIN ), NVAL( MAXIN ), NXVAL( MAXIN ),
01088      $                   PVAL( MAXIN )
01089       INTEGER            INMIN( MAXIN ), INWIN( MAXIN ), INIBL( MAXIN ),
01090      $                   ISHFTS( MAXIN ), IACC22( MAXIN )
01091       DOUBLE PRECISION   A( NMAX*NMAX, NEED ), B( NMAX*NMAX, 5 ),
01092      $                   C( NCMAX*NCMAX, NCMAX*NCMAX ), D( NMAX, 12 ),
01093      $                   RESULT( 500 ), TAUA( NMAX ), TAUB( NMAX ),
01094      $                   WORK( LWORK ), X( 5*NMAX )
01095 *     ..
01096 *     .. External Functions ..
01097       LOGICAL            LSAMEN
01098       DOUBLE PRECISION   DLAMCH, DSECND
01099       EXTERNAL           LSAMEN, DLAMCH, DSECND
01100 *     ..
01101 *     .. External Subroutines ..
01102       EXTERNAL           ALAREQ, DCHKBB, DCHKBD, DCHKBK, DCHKBL, DCHKEC,
01103      $                   DCHKGG, DCHKGK, DCHKGL, DCHKHS, DCHKSB, DCHKST,
01104      $                   DCKCSD, DCKGLM, DCKGQR, DCKGSV, DCKLSE, DDRGES,
01105      $                   DDRGEV, DDRGSX, DDRGVX, DDRVBD, DDRVES, DDRVEV,
01106      $                   DDRVGG, DDRVSG, DDRVST, DDRVSX, DDRVVX, DERRBD,
01107      $                   DERRED, DERRGG, DERRHS, DERRST, ILAVER, XLAENV
01108 *     ..
01109 *     .. Intrinsic Functions ..
01110       INTRINSIC          LEN, MIN
01111 *     ..
01112 *     .. Scalars in Common ..
01113       LOGICAL            LERR, OK
01114       CHARACTER*32       SRNAMT
01115       INTEGER            INFOT, MAXB, NPROC, NSHIFT, NUNIT, SELDIM,
01116      $                   SELOPT
01117 *     ..
01118 *     .. Arrays in Common ..
01119       LOGICAL            SELVAL( 20 )
01120       INTEGER            IPARMS( 100 )
01121       DOUBLE PRECISION   SELWI( 20 ), SELWR( 20 )
01122 *     ..
01123 *     .. Common blocks ..
01124       COMMON             / CENVIR / NPROC, NSHIFT, MAXB
01125       COMMON             / INFOC / INFOT, NUNIT, OK, LERR
01126       COMMON             / SRNAMC / SRNAMT
01127       COMMON             / SSLCT / SELOPT, SELDIM, SELVAL, SELWR, SELWI
01128       COMMON             / CLAENV / IPARMS
01129 *     ..
01130 *     .. Data statements ..
01131       DATA               INTSTR / '0123456789' /
01132       DATA               IOLDSD / 0, 0, 0, 1 /
01133 *     ..
01134 *     .. Executable Statements ..
01135 *
01136       S1 = DSECND( )
01137       FATAL = .FALSE.
01138       NUNIT = NOUT
01139 *
01140 *     Return to here to read multiple sets of data
01141 *
01142    10 CONTINUE
01143 *
01144 *     Read the first line and set the 3-character test path
01145 *
01146       READ( NIN, FMT = '(A80)', END = 380 )LINE
01147       PATH = LINE( 1: 3 )
01148       NEP = LSAMEN( 3, PATH, 'NEP' ) .OR. LSAMEN( 3, PATH, 'DHS' )
01149       SEP = LSAMEN( 3, PATH, 'SEP' ) .OR. LSAMEN( 3, PATH, 'DST' ) .OR.
01150      $      LSAMEN( 3, PATH, 'DSG' )
01151       SVD = LSAMEN( 3, PATH, 'SVD' ) .OR. LSAMEN( 3, PATH, 'DBD' )
01152       DEV = LSAMEN( 3, PATH, 'DEV' )
01153       DES = LSAMEN( 3, PATH, 'DES' )
01154       DVX = LSAMEN( 3, PATH, 'DVX' )
01155       DSX = LSAMEN( 3, PATH, 'DSX' )
01156       DGG = LSAMEN( 3, PATH, 'DGG' )
01157       DGS = LSAMEN( 3, PATH, 'DGS' )
01158       DGX = LSAMEN( 3, PATH, 'DGX' )
01159       DGV = LSAMEN( 3, PATH, 'DGV' )
01160       DXV = LSAMEN( 3, PATH, 'DXV' )
01161       DSB = LSAMEN( 3, PATH, 'DSB' )
01162       DBB = LSAMEN( 3, PATH, 'DBB' )
01163       GLM = LSAMEN( 3, PATH, 'GLM' )
01164       GQR = LSAMEN( 3, PATH, 'GQR' ) .OR. LSAMEN( 3, PATH, 'GRQ' )
01165       GSV = LSAMEN( 3, PATH, 'GSV' )
01166       CSD = LSAMEN( 3, PATH, 'CSD' )
01167       LSE = LSAMEN( 3, PATH, 'LSE' )
01168       DBL = LSAMEN( 3, PATH, 'DBL' )
01169       DBK = LSAMEN( 3, PATH, 'DBK' )
01170       DGL = LSAMEN( 3, PATH, 'DGL' )
01171       DGK = LSAMEN( 3, PATH, 'DGK' )
01172 *
01173 *     Report values of parameters.
01174 *
01175       IF( PATH.EQ.'   ' ) THEN
01176          GO TO 10
01177       ELSE IF( NEP ) THEN
01178          WRITE( NOUT, FMT = 9987 )
01179       ELSE IF( SEP ) THEN
01180          WRITE( NOUT, FMT = 9986 )
01181       ELSE IF( SVD ) THEN
01182          WRITE( NOUT, FMT = 9985 )
01183       ELSE IF( DEV ) THEN
01184          WRITE( NOUT, FMT = 9979 )
01185       ELSE IF( DES ) THEN
01186          WRITE( NOUT, FMT = 9978 )
01187       ELSE IF( DVX ) THEN
01188          WRITE( NOUT, FMT = 9977 )
01189       ELSE IF( DSX ) THEN
01190          WRITE( NOUT, FMT = 9976 )
01191       ELSE IF( DGG ) THEN
01192          WRITE( NOUT, FMT = 9975 )
01193       ELSE IF( DGS ) THEN
01194          WRITE( NOUT, FMT = 9964 )
01195       ELSE IF( DGX ) THEN
01196          WRITE( NOUT, FMT = 9965 )
01197       ELSE IF( DGV ) THEN
01198          WRITE( NOUT, FMT = 9963 )
01199       ELSE IF( DXV ) THEN
01200          WRITE( NOUT, FMT = 9962 )
01201       ELSE IF( DSB ) THEN
01202          WRITE( NOUT, FMT = 9974 )
01203       ELSE IF( DBB ) THEN
01204          WRITE( NOUT, FMT = 9967 )
01205       ELSE IF( GLM ) THEN
01206          WRITE( NOUT, FMT = 9971 )
01207       ELSE IF( GQR ) THEN
01208          WRITE( NOUT, FMT = 9970 )
01209       ELSE IF( GSV ) THEN
01210          WRITE( NOUT, FMT = 9969 )
01211       ELSE IF( CSD ) THEN
01212          WRITE( NOUT, FMT = 9960 )
01213       ELSE IF( LSE ) THEN
01214          WRITE( NOUT, FMT = 9968 )
01215       ELSE IF( DBL ) THEN
01216 *
01217 *        DGEBAL:  Balancing
01218 *
01219          CALL DCHKBL( NIN, NOUT )
01220          GO TO 10
01221       ELSE IF( DBK ) THEN
01222 *
01223 *        DGEBAK:  Back transformation
01224 *
01225          CALL DCHKBK( NIN, NOUT )
01226          GO TO 10
01227       ELSE IF( DGL ) THEN
01228 *
01229 *        DGGBAL:  Balancing
01230 *
01231          CALL DCHKGL( NIN, NOUT )
01232          GO TO 10
01233       ELSE IF( DGK ) THEN
01234 *
01235 *        DGGBAK:  Back transformation
01236 *
01237          CALL DCHKGK( NIN, NOUT )
01238          GO TO 10
01239       ELSE IF( LSAMEN( 3, PATH, 'DEC' ) ) THEN
01240 *
01241 *        DEC:  Eigencondition estimation
01242 *
01243          READ( NIN, FMT = * )THRESH
01244          CALL XLAENV( 1, 1 )
01245          CALL XLAENV( 12, 11 )
01246          CALL XLAENV( 13, 2 )
01247          CALL XLAENV( 14, 0 )
01248          CALL XLAENV( 15, 2 )
01249          CALL XLAENV( 16, 2 )
01250          TSTERR = .TRUE.
01251          CALL DCHKEC( THRESH, TSTERR, NIN, NOUT )
01252          GO TO 10
01253       ELSE
01254          WRITE( NOUT, FMT = 9992 )PATH
01255          GO TO 10
01256       END IF
01257       CALL ILAVER( VERS_MAJOR, VERS_MINOR, VERS_PATCH )
01258       WRITE( NOUT, FMT = 9972 ) VERS_MAJOR, VERS_MINOR, VERS_PATCH
01259       WRITE( NOUT, FMT = 9984 )
01260 *
01261 *     Read the number of values of M, P, and N.
01262 *
01263       READ( NIN, FMT = * )NN
01264       IF( NN.LT.0 ) THEN
01265          WRITE( NOUT, FMT = 9989 )'   NN ', NN, 1
01266          NN = 0
01267          FATAL = .TRUE.
01268       ELSE IF( NN.GT.MAXIN ) THEN
01269          WRITE( NOUT, FMT = 9988 )'   NN ', NN, MAXIN
01270          NN = 0
01271          FATAL = .TRUE.
01272       END IF
01273 *
01274 *     Read the values of M
01275 *
01276       IF( .NOT.( DGX .OR. DXV ) ) THEN
01277          READ( NIN, FMT = * )( MVAL( I ), I = 1, NN )
01278          IF( SVD ) THEN
01279             VNAME = '    M '
01280          ELSE
01281             VNAME = '    N '
01282          END IF
01283          DO 20 I = 1, NN
01284             IF( MVAL( I ).LT.0 ) THEN
01285                WRITE( NOUT, FMT = 9989 )VNAME, MVAL( I ), 0
01286                FATAL = .TRUE.
01287             ELSE IF( MVAL( I ).GT.NMAX ) THEN
01288                WRITE( NOUT, FMT = 9988 )VNAME, MVAL( I ), NMAX
01289                FATAL = .TRUE.
01290             END IF
01291    20    CONTINUE
01292          WRITE( NOUT, FMT = 9983 )'M:    ', ( MVAL( I ), I = 1, NN )
01293       END IF
01294 *
01295 *     Read the values of P
01296 *
01297       IF( GLM .OR. GQR .OR. GSV .OR. CSD .OR. LSE ) THEN
01298          READ( NIN, FMT = * )( PVAL( I ), I = 1, NN )
01299          DO 30 I = 1, NN
01300             IF( PVAL( I ).LT.0 ) THEN
01301                WRITE( NOUT, FMT = 9989 )' P  ', PVAL( I ), 0
01302                FATAL = .TRUE.
01303             ELSE IF( PVAL( I ).GT.NMAX ) THEN
01304                WRITE( NOUT, FMT = 9988 )' P  ', PVAL( I ), NMAX
01305                FATAL = .TRUE.
01306             END IF
01307    30    CONTINUE
01308          WRITE( NOUT, FMT = 9983 )'P:    ', ( PVAL( I ), I = 1, NN )
01309       END IF
01310 *
01311 *     Read the values of N
01312 *
01313       IF( SVD .OR. DBB .OR. GLM .OR. GQR .OR. GSV .OR. CSD .OR.
01314      $    LSE ) THEN
01315          READ( NIN, FMT = * )( NVAL( I ), I = 1, NN )
01316          DO 40 I = 1, NN
01317             IF( NVAL( I ).LT.0 ) THEN
01318                WRITE( NOUT, FMT = 9989 )'    N ', NVAL( I ), 0
01319                FATAL = .TRUE.
01320             ELSE IF( NVAL( I ).GT.NMAX ) THEN
01321                WRITE( NOUT, FMT = 9988 )'    N ', NVAL( I ), NMAX
01322                FATAL = .TRUE.
01323             END IF
01324    40    CONTINUE
01325       ELSE
01326          DO 50 I = 1, NN
01327             NVAL( I ) = MVAL( I )
01328    50    CONTINUE
01329       END IF
01330       IF( .NOT.( DGX .OR. DXV ) ) THEN
01331          WRITE( NOUT, FMT = 9983 )'N:    ', ( NVAL( I ), I = 1, NN )
01332       ELSE
01333          WRITE( NOUT, FMT = 9983 )'N:    ', NN
01334       END IF
01335 *
01336 *     Read the number of values of K, followed by the values of K
01337 *
01338       IF( DSB .OR. DBB ) THEN
01339          READ( NIN, FMT = * )NK
01340          READ( NIN, FMT = * )( KVAL( I ), I = 1, NK )
01341          DO 60 I = 1, NK
01342             IF( KVAL( I ).LT.0 ) THEN
01343                WRITE( NOUT, FMT = 9989 )'    K ', KVAL( I ), 0
01344                FATAL = .TRUE.
01345             ELSE IF( KVAL( I ).GT.NMAX ) THEN
01346                WRITE( NOUT, FMT = 9988 )'    K ', KVAL( I ), NMAX
01347                FATAL = .TRUE.
01348             END IF
01349    60    CONTINUE
01350          WRITE( NOUT, FMT = 9983 )'K:    ', ( KVAL( I ), I = 1, NK )
01351       END IF
01352 *
01353       IF( DEV .OR. DES .OR. DVX .OR. DSX ) THEN
01354 *
01355 *        For the nonsymmetric QR driver routines, only one set of
01356 *        parameters is allowed.
01357 *
01358          READ( NIN, FMT = * )NBVAL( 1 ), NBMIN( 1 ), NXVAL( 1 ),
01359      $      INMIN( 1 ), INWIN( 1 ), INIBL(1), ISHFTS(1), IACC22(1)
01360          IF( NBVAL( 1 ).LT.1 ) THEN
01361             WRITE( NOUT, FMT = 9989 )'   NB ', NBVAL( 1 ), 1
01362             FATAL = .TRUE.
01363          ELSE IF( NBMIN( 1 ).LT.1 ) THEN
01364             WRITE( NOUT, FMT = 9989 )'NBMIN ', NBMIN( 1 ), 1
01365             FATAL = .TRUE.
01366          ELSE IF( NXVAL( 1 ).LT.1 ) THEN
01367             WRITE( NOUT, FMT = 9989 )'   NX ', NXVAL( 1 ), 1
01368             FATAL = .TRUE.
01369          ELSE IF( INMIN( 1 ).LT.1 ) THEN
01370             WRITE( NOUT, FMT = 9989 )'   INMIN ', INMIN( 1 ), 1
01371             FATAL = .TRUE.
01372          ELSE IF( INWIN( 1 ).LT.1 ) THEN
01373             WRITE( NOUT, FMT = 9989 )'   INWIN ', INWIN( 1 ), 1
01374             FATAL = .TRUE.
01375          ELSE IF( INIBL( 1 ).LT.1 ) THEN
01376             WRITE( NOUT, FMT = 9989 )'   INIBL ', INIBL( 1 ), 1
01377             FATAL = .TRUE.
01378          ELSE IF( ISHFTS( 1 ).LT.1 ) THEN
01379             WRITE( NOUT, FMT = 9989 )'   ISHFTS ', ISHFTS( 1 ), 1
01380             FATAL = .TRUE.
01381          ELSE IF( IACC22( 1 ).LT.0 ) THEN
01382             WRITE( NOUT, FMT = 9989 )'   IACC22 ', IACC22( 1 ), 0
01383             FATAL = .TRUE.
01384          END IF
01385          CALL XLAENV( 1, NBVAL( 1 ) )
01386          CALL XLAENV( 2, NBMIN( 1 ) )
01387          CALL XLAENV( 3, NXVAL( 1 ) )
01388          CALL XLAENV(12, MAX( 11, INMIN( 1 ) ) )
01389          CALL XLAENV(13, INWIN( 1 ) )
01390          CALL XLAENV(14, INIBL( 1 ) )
01391          CALL XLAENV(15, ISHFTS( 1 ) )
01392          CALL XLAENV(16, IACC22( 1 ) )
01393          WRITE( NOUT, FMT = 9983 )'NB:   ', NBVAL( 1 )
01394          WRITE( NOUT, FMT = 9983 )'NBMIN:', NBMIN( 1 )
01395          WRITE( NOUT, FMT = 9983 )'NX:   ', NXVAL( 1 )
01396          WRITE( NOUT, FMT = 9983 )'INMIN:   ', INMIN( 1 )
01397          WRITE( NOUT, FMT = 9983 )'INWIN: ', INWIN( 1 )
01398          WRITE( NOUT, FMT = 9983 )'INIBL: ', INIBL( 1 )
01399          WRITE( NOUT, FMT = 9983 )'ISHFTS: ', ISHFTS( 1 )
01400          WRITE( NOUT, FMT = 9983 )'IACC22: ', IACC22( 1 )
01401 *
01402       ELSEIF( DGS .OR. DGX .OR. DGV .OR.  DXV ) THEN
01403 *
01404 *        For the nonsymmetric generalized driver routines, only one set
01405 *        of parameters is allowed.
01406 *
01407          READ( NIN, FMT = * )NBVAL( 1 ), NBMIN( 1 ), NXVAL( 1 ),
01408      $      NSVAL( 1 ), MXBVAL( 1 )
01409          IF( NBVAL( 1 ).LT.1 ) THEN
01410             WRITE( NOUT, FMT = 9989 )'   NB ', NBVAL( 1 ), 1
01411             FATAL = .TRUE.
01412          ELSE IF( NBMIN( 1 ).LT.1 ) THEN
01413             WRITE( NOUT, FMT = 9989 )'NBMIN ', NBMIN( 1 ), 1
01414             FATAL = .TRUE.
01415          ELSE IF( NXVAL( 1 ).LT.1 ) THEN
01416             WRITE( NOUT, FMT = 9989 )'   NX ', NXVAL( 1 ), 1
01417             FATAL = .TRUE.
01418          ELSE IF( NSVAL( 1 ).LT.2 ) THEN
01419             WRITE( NOUT, FMT = 9989 )'   NS ', NSVAL( 1 ), 2
01420             FATAL = .TRUE.
01421          ELSE IF( MXBVAL( 1 ).LT.1 ) THEN
01422             WRITE( NOUT, FMT = 9989 )' MAXB ', MXBVAL( 1 ), 1
01423             FATAL = .TRUE.
01424          END IF
01425          CALL XLAENV( 1, NBVAL( 1 ) )
01426          CALL XLAENV( 2, NBMIN( 1 ) )
01427          CALL XLAENV( 3, NXVAL( 1 ) )
01428          CALL XLAENV( 4, NSVAL( 1 ) )
01429          CALL XLAENV( 8, MXBVAL( 1 ) )
01430          WRITE( NOUT, FMT = 9983 )'NB:   ', NBVAL( 1 )
01431          WRITE( NOUT, FMT = 9983 )'NBMIN:', NBMIN( 1 )
01432          WRITE( NOUT, FMT = 9983 )'NX:   ', NXVAL( 1 )
01433          WRITE( NOUT, FMT = 9983 )'NS:   ', NSVAL( 1 )
01434          WRITE( NOUT, FMT = 9983 )'MAXB: ', MXBVAL( 1 )
01435 *
01436       ELSE IF( .NOT.DSB .AND. .NOT.GLM .AND. .NOT.GQR .AND. .NOT.
01437      $         GSV .AND. .NOT.CSD .AND. .NOT.LSE ) THEN
01438 *
01439 *        For the other paths, the number of parameters can be varied
01440 *        from the input file.  Read the number of parameter values.
01441 *
01442          READ( NIN, FMT = * )NPARMS
01443          IF( NPARMS.LT.1 ) THEN
01444             WRITE( NOUT, FMT = 9989 )'NPARMS', NPARMS, 1
01445             NPARMS = 0
01446             FATAL = .TRUE.
01447          ELSE IF( NPARMS.GT.MAXIN ) THEN
01448             WRITE( NOUT, FMT = 9988 )'NPARMS', NPARMS, MAXIN
01449             NPARMS = 0
01450             FATAL = .TRUE.
01451          END IF
01452 *
01453 *        Read the values of NB
01454 *
01455          IF( .NOT.DBB ) THEN
01456             READ( NIN, FMT = * )( NBVAL( I ), I = 1, NPARMS )
01457             DO 70 I = 1, NPARMS
01458                IF( NBVAL( I ).LT.0 ) THEN
01459                   WRITE( NOUT, FMT = 9989 )'   NB ', NBVAL( I ), 0
01460                   FATAL = .TRUE.
01461                ELSE IF( NBVAL( I ).GT.NMAX ) THEN
01462                   WRITE( NOUT, FMT = 9988 )'   NB ', NBVAL( I ), NMAX
01463                   FATAL = .TRUE.
01464                END IF
01465    70       CONTINUE
01466             WRITE( NOUT, FMT = 9983 )'NB:   ',
01467      $         ( NBVAL( I ), I = 1, NPARMS )
01468          END IF
01469 *
01470 *        Read the values of NBMIN
01471 *
01472          IF( NEP .OR. SEP .OR. SVD .OR. DGG ) THEN
01473             READ( NIN, FMT = * )( NBMIN( I ), I = 1, NPARMS )
01474             DO 80 I = 1, NPARMS
01475                IF( NBMIN( I ).LT.0 ) THEN
01476                   WRITE( NOUT, FMT = 9989 )'NBMIN ', NBMIN( I ), 0
01477                   FATAL = .TRUE.
01478                ELSE IF( NBMIN( I ).GT.NMAX ) THEN
01479                   WRITE( NOUT, FMT = 9988 )'NBMIN ', NBMIN( I ), NMAX
01480                   FATAL = .TRUE.
01481                END IF
01482    80       CONTINUE
01483             WRITE( NOUT, FMT = 9983 )'NBMIN:',
01484      $         ( NBMIN( I ), I = 1, NPARMS )
01485          ELSE
01486             DO 90 I = 1, NPARMS
01487                NBMIN( I ) = 1
01488    90       CONTINUE
01489          END IF
01490 *
01491 *        Read the values of NX
01492 *
01493          IF( NEP .OR. SEP .OR. SVD ) THEN
01494             READ( NIN, FMT = * )( NXVAL( I ), I = 1, NPARMS )
01495             DO 100 I = 1, NPARMS
01496                IF( NXVAL( I ).LT.0 ) THEN
01497                   WRITE( NOUT, FMT = 9989 )'   NX ', NXVAL( I ), 0
01498                   FATAL = .TRUE.
01499                ELSE IF( NXVAL( I ).GT.NMAX ) THEN
01500                   WRITE( NOUT, FMT = 9988 )'   NX ', NXVAL( I ), NMAX
01501                   FATAL = .TRUE.
01502                END IF
01503   100       CONTINUE
01504             WRITE( NOUT, FMT = 9983 )'NX:   ',
01505      $         ( NXVAL( I ), I = 1, NPARMS )
01506          ELSE
01507             DO 110 I = 1, NPARMS
01508                NXVAL( I ) = 1
01509   110       CONTINUE
01510          END IF
01511 *
01512 *        Read the values of NSHIFT (if DGG) or NRHS (if SVD
01513 *        or DBB).
01514 *
01515          IF( SVD .OR. DBB .OR. DGG ) THEN
01516             READ( NIN, FMT = * )( NSVAL( I ), I = 1, NPARMS )
01517             DO 120 I = 1, NPARMS
01518                IF( NSVAL( I ).LT.0 ) THEN
01519                   WRITE( NOUT, FMT = 9989 )'   NS ', NSVAL( I ), 0
01520                   FATAL = .TRUE.
01521                ELSE IF( NSVAL( I ).GT.NMAX ) THEN
01522                   WRITE( NOUT, FMT = 9988 )'   NS ', NSVAL( I ), NMAX
01523                   FATAL = .TRUE.
01524                END IF
01525   120       CONTINUE
01526             WRITE( NOUT, FMT = 9983 )'NS:   ',
01527      $         ( NSVAL( I ), I = 1, NPARMS )
01528          ELSE
01529             DO 130 I = 1, NPARMS
01530                NSVAL( I ) = 1
01531   130       CONTINUE
01532          END IF
01533 *
01534 *        Read the values for MAXB.
01535 *
01536          IF( DGG ) THEN
01537             READ( NIN, FMT = * )( MXBVAL( I ), I = 1, NPARMS )
01538             DO 140 I = 1, NPARMS
01539                IF( MXBVAL( I ).LT.0 ) THEN
01540                   WRITE( NOUT, FMT = 9989 )' MAXB ', MXBVAL( I ), 0
01541                   FATAL = .TRUE.
01542                ELSE IF( MXBVAL( I ).GT.NMAX ) THEN
01543                   WRITE( NOUT, FMT = 9988 )' MAXB ', MXBVAL( I ), NMAX
01544                   FATAL = .TRUE.
01545                END IF
01546   140       CONTINUE
01547             WRITE( NOUT, FMT = 9983 )'MAXB: ',
01548      $         ( MXBVAL( I ), I = 1, NPARMS )
01549          ELSE
01550             DO 150 I = 1, NPARMS
01551                MXBVAL( I ) = 1
01552   150       CONTINUE
01553          END IF
01554 *
01555 *        Read the values for INMIN.
01556 *
01557          IF( NEP ) THEN
01558             READ( NIN, FMT = * )( INMIN( I ), I = 1, NPARMS )
01559             DO 540 I = 1, NPARMS
01560                IF( INMIN( I ).LT.0 ) THEN
01561                   WRITE( NOUT, FMT = 9989 )' INMIN ', INMIN( I ), 0
01562                   FATAL = .TRUE.
01563                END IF
01564   540       CONTINUE
01565             WRITE( NOUT, FMT = 9983 )'INMIN: ',
01566      $         ( INMIN( I ), I = 1, NPARMS )
01567          ELSE
01568             DO 550 I = 1, NPARMS
01569                INMIN( I ) = 1
01570   550       CONTINUE
01571          END IF
01572 *
01573 *        Read the values for INWIN.
01574 *
01575          IF( NEP ) THEN
01576             READ( NIN, FMT = * )( INWIN( I ), I = 1, NPARMS )
01577             DO 560 I = 1, NPARMS
01578                IF( INWIN( I ).LT.0 ) THEN
01579                   WRITE( NOUT, FMT = 9989 )' INWIN ', INWIN( I ), 0
01580                   FATAL = .TRUE.
01581                END IF
01582   560       CONTINUE
01583             WRITE( NOUT, FMT = 9983 )'INWIN: ',
01584      $         ( INWIN( I ), I = 1, NPARMS )
01585          ELSE
01586             DO 570 I = 1, NPARMS
01587                INWIN( I ) = 1
01588   570       CONTINUE
01589          END IF
01590 *
01591 *        Read the values for INIBL.
01592 *
01593          IF( NEP ) THEN
01594             READ( NIN, FMT = * )( INIBL( I ), I = 1, NPARMS )
01595             DO 580 I = 1, NPARMS
01596                IF( INIBL( I ).LT.0 ) THEN
01597                   WRITE( NOUT, FMT = 9989 )' INIBL ', INIBL( I ), 0
01598                   FATAL = .TRUE.
01599                END IF
01600   580       CONTINUE
01601             WRITE( NOUT, FMT = 9983 )'INIBL: ',
01602      $         ( INIBL( I ), I = 1, NPARMS )
01603          ELSE
01604             DO 590 I = 1, NPARMS
01605                INIBL( I ) = 1
01606   590       CONTINUE
01607          END IF
01608 *
01609 *        Read the values for ISHFTS.
01610 *
01611          IF( NEP ) THEN
01612             READ( NIN, FMT = * )( ISHFTS( I ), I = 1, NPARMS )
01613             DO 600 I = 1, NPARMS
01614                IF( ISHFTS( I ).LT.0 ) THEN
01615                   WRITE( NOUT, FMT = 9989 )' ISHFTS ', ISHFTS( I ), 0
01616                   FATAL = .TRUE.
01617                END IF
01618   600       CONTINUE
01619             WRITE( NOUT, FMT = 9983 )'ISHFTS: ',
01620      $         ( ISHFTS( I ), I = 1, NPARMS )
01621          ELSE
01622             DO 610 I = 1, NPARMS
01623                ISHFTS( I ) = 1
01624   610       CONTINUE
01625          END IF
01626 *
01627 *        Read the values for IACC22.
01628 *
01629          IF( NEP ) THEN
01630             READ( NIN, FMT = * )( IACC22( I ), I = 1, NPARMS )
01631             DO 620 I = 1, NPARMS
01632                IF( IACC22( I ).LT.0 ) THEN
01633                   WRITE( NOUT, FMT = 9989 )' IACC22 ', IACC22( I ), 0
01634                   FATAL = .TRUE.
01635                END IF
01636   620       CONTINUE
01637             WRITE( NOUT, FMT = 9983 )'IACC22: ',
01638      $         ( IACC22( I ), I = 1, NPARMS )
01639          ELSE
01640             DO 630 I = 1, NPARMS
01641                IACC22( I ) = 1
01642   630       CONTINUE
01643          END IF
01644 *
01645 *        Read the values for NBCOL.
01646 *
01647          IF( DGG ) THEN
01648             READ( NIN, FMT = * )( NBCOL( I ), I = 1, NPARMS )
01649             DO 160 I = 1, NPARMS
01650                IF( NBCOL( I ).LT.0 ) THEN
01651                   WRITE( NOUT, FMT = 9989 )'NBCOL ', NBCOL( I ), 0
01652                   FATAL = .TRUE.
01653                ELSE IF( NBCOL( I ).GT.NMAX ) THEN
01654                   WRITE( NOUT, FMT = 9988 )'NBCOL ', NBCOL( I ), NMAX
01655                   FATAL = .TRUE.
01656                END IF
01657   160       CONTINUE
01658             WRITE( NOUT, FMT = 9983 )'NBCOL:',
01659      $         ( NBCOL( I ), I = 1, NPARMS )
01660          ELSE
01661             DO 170 I = 1, NPARMS
01662                NBCOL( I ) = 1
01663   170       CONTINUE
01664          END IF
01665       END IF
01666 *
01667 *     Calculate and print the machine dependent constants.
01668 *
01669       WRITE( NOUT, FMT = * )
01670       EPS = DLAMCH( 'Underflow threshold' )
01671       WRITE( NOUT, FMT = 9981 )'underflow', EPS
01672       EPS = DLAMCH( 'Overflow threshold' )
01673       WRITE( NOUT, FMT = 9981 )'overflow ', EPS
01674       EPS = DLAMCH( 'Epsilon' )
01675       WRITE( NOUT, FMT = 9981 )'precision', EPS
01676 *
01677 *     Read the threshold value for the test ratios.
01678 *
01679       READ( NIN, FMT = * )THRESH
01680       WRITE( NOUT, FMT = 9982 )THRESH
01681       IF( SEP .OR. SVD .OR. DGG ) THEN
01682 *
01683 *        Read the flag that indicates whether to test LAPACK routines.
01684 *
01685          READ( NIN, FMT = * )TSTCHK
01686 *
01687 *        Read the flag that indicates whether to test driver routines.
01688 *
01689          READ( NIN, FMT = * )TSTDRV
01690       END IF
01691 *
01692 *     Read the flag that indicates whether to test the error exits.
01693 *
01694       READ( NIN, FMT = * )TSTERR
01695 *
01696 *     Read the code describing how to set the random number seed.
01697 *
01698       READ( NIN, FMT = * )NEWSD
01699 *
01700 *     If NEWSD = 2, read another line with 4 integers for the seed.
01701 *
01702       IF( NEWSD.EQ.2 )
01703      $   READ( NIN, FMT = * )( IOLDSD( I ), I = 1, 4 )
01704 *
01705       DO 180 I = 1, 4
01706          ISEED( I ) = IOLDSD( I )
01707   180 CONTINUE
01708 *
01709       IF( FATAL ) THEN
01710          WRITE( NOUT, FMT = 9999 )
01711          STOP
01712       END IF
01713 *
01714 *     Read the input lines indicating the test path and its parameters.
01715 *     The first three characters indicate the test path, and the number
01716 *     of test matrix types must be the first nonblank item in columns
01717 *     4-80.
01718 *
01719   190 CONTINUE
01720 *
01721       IF( .NOT.( DGX .OR. DXV ) ) THEN
01722 *
01723   200    CONTINUE
01724          READ( NIN, FMT = '(A80)', END = 380 )LINE
01725          C3 = LINE( 1: 3 )
01726          LENP = LEN( LINE )
01727          I = 3
01728          ITMP = 0
01729          I1 = 0
01730   210    CONTINUE
01731          I = I + 1
01732          IF( I.GT.LENP ) THEN
01733             IF( I1.GT.0 ) THEN
01734                GO TO 240
01735             ELSE
01736                NTYPES = MAXT
01737                GO TO 240
01738             END IF
01739          END IF
01740          IF( LINE( I: I ).NE.' ' .AND. LINE( I: I ).NE.',' ) THEN
01741             I1 = I
01742             C1 = LINE( I1: I1 )
01743 *
01744 *        Check that a valid integer was read
01745 *
01746             DO 220 K = 1, 10
01747                IF( C1.EQ.INTSTR( K: K ) ) THEN
01748                   IC = K - 1
01749                   GO TO 230
01750                END IF
01751   220       CONTINUE
01752             WRITE( NOUT, FMT = 9991 )I, LINE
01753             GO TO 200
01754   230       CONTINUE
01755             ITMP = 10*ITMP + IC
01756             GO TO 210
01757          ELSE IF( I1.GT.0 ) THEN
01758             GO TO 240
01759          ELSE
01760             GO TO 210
01761          END IF
01762   240    CONTINUE
01763          NTYPES = ITMP
01764 *
01765 *     Skip the tests if NTYPES is <= 0.
01766 *
01767          IF( .NOT.( DEV .OR. DES .OR. DVX .OR. DSX .OR. DGV .OR.
01768      $       DGS ) .AND. NTYPES.LE.0 ) THEN
01769             WRITE( NOUT, FMT = 9990 )C3
01770             GO TO 200
01771          END IF
01772 *
01773       ELSE
01774          IF( DXV )
01775      $      C3 = 'DXV'
01776          IF( DGX )
01777      $      C3 = 'DGX'
01778       END IF
01779 *
01780 *     Reset the random number seed.
01781 *
01782       IF( NEWSD.EQ.0 ) THEN
01783          DO 250 K = 1, 4
01784             ISEED( K ) = IOLDSD( K )
01785   250    CONTINUE
01786       END IF
01787 *
01788       IF( LSAMEN( 3, C3, 'DHS' ) .OR. LSAMEN( 3, C3, 'NEP' ) ) THEN
01789 *
01790 *        -------------------------------------
01791 *        NEP:  Nonsymmetric Eigenvalue Problem
01792 *        -------------------------------------
01793 *        Vary the parameters
01794 *           NB    = block size
01795 *           NBMIN = minimum block size
01796 *           NX    = crossover point
01797 *           NS    = number of shifts
01798 *           MAXB  = minimum submatrix size
01799 *
01800          MAXTYP = 21
01801          NTYPES = MIN( MAXTYP, NTYPES )
01802          CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
01803          CALL XLAENV( 1, 1 )
01804          IF( TSTERR )
01805      $      CALL DERRHS( 'DHSEQR', NOUT )
01806          DO 270 I = 1, NPARMS
01807             CALL XLAENV( 1, NBVAL( I ) )
01808             CALL XLAENV( 2, NBMIN( I ) )
01809             CALL XLAENV( 3, NXVAL( I ) )
01810             CALL XLAENV(12, MAX( 11, INMIN( I ) ) )
01811             CALL XLAENV(13, INWIN( I ) )
01812             CALL XLAENV(14, INIBL( I ) )
01813             CALL XLAENV(15, ISHFTS( I ) )
01814             CALL XLAENV(16, IACC22( I ) )
01815 *
01816             IF( NEWSD.EQ.0 ) THEN
01817                DO 260 K = 1, 4
01818                   ISEED( K ) = IOLDSD( K )
01819   260          CONTINUE
01820             END IF
01821             WRITE( NOUT, FMT = 9961 )C3, NBVAL( I ), NBMIN( I ),
01822      $         NXVAL( I ), MAX( 11, INMIN(I)),
01823      $         INWIN( I ), INIBL( I ), ISHFTS( I ), IACC22( I )
01824             CALL DCHKHS( NN, NVAL, MAXTYP, DOTYPE, ISEED, THRESH, NOUT,
01825      $                   A( 1, 1 ), NMAX, A( 1, 2 ), A( 1, 3 ),
01826      $                   A( 1, 4 ), A( 1, 5 ), NMAX, A( 1, 6 ),
01827      $                   A( 1, 7 ), D( 1, 1 ), D( 1, 2 ), D( 1, 3 ),
01828      $                   D( 1, 4 ), A( 1, 8 ), A( 1, 9 ), A( 1, 10 ),
01829      $                   A( 1, 11 ), A( 1, 12 ), D( 1, 5 ), WORK, LWORK,
01830      $                   IWORK, LOGWRK, RESULT, INFO )
01831             IF( INFO.NE.0 )
01832      $         WRITE( NOUT, FMT = 9980 )'DCHKHS', INFO
01833   270    CONTINUE
01834 *
01835       ELSE IF( LSAMEN( 3, C3, 'DST' ) .OR. LSAMEN( 3, C3, 'SEP' ) ) THEN
01836 *
01837 *        ----------------------------------
01838 *        SEP:  Symmetric Eigenvalue Problem
01839 *        ----------------------------------
01840 *        Vary the parameters
01841 *           NB    = block size
01842 *           NBMIN = minimum block size
01843 *           NX    = crossover point
01844 *
01845          MAXTYP = 21
01846          NTYPES = MIN( MAXTYP, NTYPES )
01847          CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
01848          CALL XLAENV( 1, 1 )
01849          CALL XLAENV( 9, 25 )
01850          IF( TSTERR )
01851      $      CALL DERRST( 'DST', NOUT )
01852          DO 290 I = 1, NPARMS
01853             CALL XLAENV( 1, NBVAL( I ) )
01854             CALL XLAENV( 2, NBMIN( I ) )
01855             CALL XLAENV( 3, NXVAL( I ) )
01856 *
01857             IF( NEWSD.EQ.0 ) THEN
01858                DO 280 K = 1, 4
01859                   ISEED( K ) = IOLDSD( K )
01860   280          CONTINUE
01861             END IF
01862             WRITE( NOUT, FMT = 9997 )C3, NBVAL( I ), NBMIN( I ),
01863      $         NXVAL( I )
01864             IF( TSTCHK ) THEN
01865                CALL DCHKST( NN, NVAL, MAXTYP, DOTYPE, ISEED, THRESH,
01866      $                      NOUT, A( 1, 1 ), NMAX, A( 1, 2 ), D( 1, 1 ),
01867      $                      D( 1, 2 ), D( 1, 3 ), D( 1, 4 ), D( 1, 5 ),
01868      $                      D( 1, 6 ), D( 1, 7 ), D( 1, 8 ), D( 1, 9 ),
01869      $                      D( 1, 10 ), D( 1, 11 ), A( 1, 3 ), NMAX,
01870      $                      A( 1, 4 ), A( 1, 5 ), D( 1, 12 ), A( 1, 6 ),
01871      $                      WORK, LWORK, IWORK, LIWORK, RESULT, INFO )
01872                IF( INFO.NE.0 )
01873      $            WRITE( NOUT, FMT = 9980 )'DCHKST', INFO
01874             END IF
01875             IF( TSTDRV ) THEN
01876                CALL DDRVST( NN, NVAL, 18, DOTYPE, ISEED, THRESH, NOUT,
01877      $                      A( 1, 1 ), NMAX, D( 1, 3 ), D( 1, 4 ),
01878      $                      D( 1, 5 ), D( 1, 6 ), D( 1, 8 ), D( 1, 9 ),
01879      $                      D( 1, 10 ), D( 1, 11 ), A( 1, 2 ), NMAX,
01880      $                      A( 1, 3 ), D( 1, 12 ), A( 1, 4 ), WORK,
01881      $                      LWORK, IWORK, LIWORK, RESULT, INFO )
01882                IF( INFO.NE.0 )
01883      $            WRITE( NOUT, FMT = 9980 )'DDRVST', INFO
01884             END IF
01885   290    CONTINUE
01886 *
01887       ELSE IF( LSAMEN( 3, C3, 'DSG' ) ) THEN
01888 *
01889 *        ----------------------------------------------
01890 *        DSG:  Symmetric Generalized Eigenvalue Problem
01891 *        ----------------------------------------------
01892 *        Vary the parameters
01893 *           NB    = block size
01894 *           NBMIN = minimum block size
01895 *           NX    = crossover point
01896 *
01897          MAXTYP = 21
01898          NTYPES = MIN( MAXTYP, NTYPES )
01899          CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
01900          CALL XLAENV( 9, 25 )
01901          DO 310 I = 1, NPARMS
01902             CALL XLAENV( 1, NBVAL( I ) )
01903             CALL XLAENV( 2, NBMIN( I ) )
01904             CALL XLAENV( 3, NXVAL( I ) )
01905 *
01906             IF( NEWSD.EQ.0 ) THEN
01907                DO 300 K = 1, 4
01908                   ISEED( K ) = IOLDSD( K )
01909   300          CONTINUE
01910             END IF
01911             WRITE( NOUT, FMT = 9997 )C3, NBVAL( I ), NBMIN( I ),
01912      $         NXVAL( I )
01913             IF( TSTCHK ) THEN
01914                CALL DDRVSG( NN, NVAL, MAXTYP, DOTYPE, ISEED, THRESH,
01915      $                      NOUT, A( 1, 1 ), NMAX, A( 1, 2 ), NMAX,
01916      $                      D( 1, 3 ), A( 1, 3 ), NMAX, A( 1, 4 ),
01917      $                      A( 1, 5 ), A( 1, 6 ), A( 1, 7 ), WORK,
01918      $                      LWORK, IWORK, LIWORK, RESULT, INFO )
01919                IF( INFO.NE.0 )
01920      $            WRITE( NOUT, FMT = 9980 )'DDRVSG', INFO
01921             END IF
01922   310    CONTINUE
01923 *
01924       ELSE IF( LSAMEN( 3, C3, 'DBD' ) .OR. LSAMEN( 3, C3, 'SVD' ) ) THEN
01925 *
01926 *        ----------------------------------
01927 *        SVD:  Singular Value Decomposition
01928 *        ----------------------------------
01929 *        Vary the parameters
01930 *           NB    = block size
01931 *           NBMIN = minimum block size
01932 *           NX    = crossover point
01933 *           NRHS  = number of right hand sides
01934 *
01935          MAXTYP = 16
01936          NTYPES = MIN( MAXTYP, NTYPES )
01937          CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
01938          CALL XLAENV( 1, 1 )
01939          CALL XLAENV( 9, 25 )
01940 *
01941 *        Test the error exits
01942 *
01943          IF( TSTERR .AND. TSTCHK )
01944      $      CALL DERRBD( 'DBD', NOUT )
01945          IF( TSTERR .AND. TSTDRV )
01946      $      CALL DERRED( 'DBD', NOUT )
01947 *
01948          DO 330 I = 1, NPARMS
01949             NRHS = NSVAL( I )
01950             CALL XLAENV( 1, NBVAL( I ) )
01951             CALL XLAENV( 2, NBMIN( I ) )
01952             CALL XLAENV( 3, NXVAL( I ) )
01953             IF( NEWSD.EQ.0 ) THEN
01954                DO 320 K = 1, 4
01955                   ISEED( K ) = IOLDSD( K )
01956   320          CONTINUE
01957             END IF
01958             WRITE( NOUT, FMT = 9995 )C3, NBVAL( I ), NBMIN( I ),
01959      $         NXVAL( I ), NRHS
01960             IF( TSTCHK ) THEN
01961                CALL DCHKBD( NN, MVAL, NVAL, MAXTYP, DOTYPE, NRHS, ISEED,
01962      $                      THRESH, A( 1, 1 ), NMAX, D( 1, 1 ),
01963      $                      D( 1, 2 ), D( 1, 3 ), D( 1, 4 ), A( 1, 2 ),
01964      $                      NMAX, A( 1, 3 ), A( 1, 4 ), A( 1, 5 ), NMAX,
01965      $                      A( 1, 6 ), NMAX, A( 1, 7 ), A( 1, 8 ), WORK,
01966      $                      LWORK, IWORK, NOUT, INFO )
01967                IF( INFO.NE.0 )
01968      $            WRITE( NOUT, FMT = 9980 )'DCHKBD', INFO
01969             END IF
01970             IF( TSTDRV )
01971      $         CALL DDRVBD( NN, MVAL, NVAL, MAXTYP, DOTYPE, ISEED,
01972      $                      THRESH, A( 1, 1 ), NMAX, A( 1, 2 ), NMAX,
01973      $                      A( 1, 3 ), NMAX, A( 1, 4 ), A( 1, 5 ),
01974      $                      A( 1, 6 ), D( 1, 1 ), D( 1, 2 ), D( 1, 3 ),
01975      $                      WORK, LWORK, IWORK, NOUT, INFO )
01976   330    CONTINUE
01977 *
01978       ELSE IF( LSAMEN( 3, C3, 'DEV' ) ) THEN
01979 *
01980 *        --------------------------------------------
01981 *        DEV:  Nonsymmetric Eigenvalue Problem Driver
01982 *              DGEEV (eigenvalues and eigenvectors)
01983 *        --------------------------------------------
01984 *
01985          MAXTYP = 21
01986          NTYPES = MIN( MAXTYP, NTYPES )
01987          IF( NTYPES.LE.0 ) THEN
01988             WRITE( NOUT, FMT = 9990 )C3
01989          ELSE
01990             IF( TSTERR )
01991      $         CALL DERRED( C3, NOUT )
01992             CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
01993             CALL DDRVEV( NN, NVAL, NTYPES, DOTYPE, ISEED, THRESH, NOUT,
01994      $                   A( 1, 1 ), NMAX, A( 1, 2 ), D( 1, 1 ),
01995      $                   D( 1, 2 ), D( 1, 3 ), D( 1, 4 ), A( 1, 3 ),
01996      $                   NMAX, A( 1, 4 ), NMAX, A( 1, 5 ), NMAX, RESULT,
01997      $                   WORK, LWORK, IWORK, INFO )
01998             IF( INFO.NE.0 )
01999      $         WRITE( NOUT, FMT = 9980 )'DGEEV', INFO
02000          END IF
02001          WRITE( NOUT, FMT = 9973 )
02002          GO TO 10
02003 *
02004       ELSE IF( LSAMEN( 3, C3, 'DES' ) ) THEN
02005 *
02006 *        --------------------------------------------
02007 *        DES:  Nonsymmetric Eigenvalue Problem Driver
02008 *              DGEES (Schur form)
02009 *        --------------------------------------------
02010 *
02011          MAXTYP = 21
02012          NTYPES = MIN( MAXTYP, NTYPES )
02013          IF( NTYPES.LE.0 ) THEN
02014             WRITE( NOUT, FMT = 9990 )C3
02015          ELSE
02016             IF( TSTERR )
02017      $         CALL DERRED( C3, NOUT )
02018             CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02019             CALL DDRVES( NN, NVAL, NTYPES, DOTYPE, ISEED, THRESH, NOUT,
02020      $                   A( 1, 1 ), NMAX, A( 1, 2 ), A( 1, 3 ),
02021      $                   D( 1, 1 ), D( 1, 2 ), D( 1, 3 ), D( 1, 4 ),
02022      $                   A( 1, 4 ), NMAX, RESULT, WORK, LWORK, IWORK,
02023      $                   LOGWRK, INFO )
02024             IF( INFO.NE.0 )
02025      $         WRITE( NOUT, FMT = 9980 )'DGEES', INFO
02026          END IF
02027          WRITE( NOUT, FMT = 9973 )
02028          GO TO 10
02029 *
02030       ELSE IF( LSAMEN( 3, C3, 'DVX' ) ) THEN
02031 *
02032 *        --------------------------------------------------------------
02033 *        DVX:  Nonsymmetric Eigenvalue Problem Expert Driver
02034 *              DGEEVX (eigenvalues, eigenvectors and condition numbers)
02035 *        --------------------------------------------------------------
02036 *
02037          MAXTYP = 21
02038          NTYPES = MIN( MAXTYP, NTYPES )
02039          IF( NTYPES.LT.0 ) THEN
02040             WRITE( NOUT, FMT = 9990 )C3
02041          ELSE
02042             IF( TSTERR )
02043      $         CALL DERRED( C3, NOUT )
02044             CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02045             CALL DDRVVX( NN, NVAL, NTYPES, DOTYPE, ISEED, THRESH, NIN,
02046      $                   NOUT, A( 1, 1 ), NMAX, A( 1, 2 ), D( 1, 1 ),
02047      $                   D( 1, 2 ), D( 1, 3 ), D( 1, 4 ), A( 1, 3 ),
02048      $                   NMAX, A( 1, 4 ), NMAX, A( 1, 5 ), NMAX,
02049      $                   D( 1, 5 ), D( 1, 6 ), D( 1, 7 ), D( 1, 8 ),
02050      $                   D( 1, 9 ), D( 1, 10 ), D( 1, 11 ), D( 1, 12 ),
02051      $                   RESULT, WORK, LWORK, IWORK, INFO )
02052             IF( INFO.NE.0 )
02053      $         WRITE( NOUT, FMT = 9980 )'DGEEVX', INFO
02054          END IF
02055          WRITE( NOUT, FMT = 9973 )
02056          GO TO 10
02057 *
02058       ELSE IF( LSAMEN( 3, C3, 'DSX' ) ) THEN
02059 *
02060 *        ---------------------------------------------------
02061 *        DSX:  Nonsymmetric Eigenvalue Problem Expert Driver
02062 *              DGEESX (Schur form and condition numbers)
02063 *        ---------------------------------------------------
02064 *
02065          MAXTYP = 21
02066          NTYPES = MIN( MAXTYP, NTYPES )
02067          IF( NTYPES.LT.0 ) THEN
02068             WRITE( NOUT, FMT = 9990 )C3
02069          ELSE
02070             IF( TSTERR )
02071      $         CALL DERRED( C3, NOUT )
02072             CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02073             CALL DDRVSX( NN, NVAL, NTYPES, DOTYPE, ISEED, THRESH, NIN,
02074      $                   NOUT, A( 1, 1 ), NMAX, A( 1, 2 ), A( 1, 3 ),
02075      $                   D( 1, 1 ), D( 1, 2 ), D( 1, 3 ), D( 1, 4 ),
02076      $                   D( 1, 5 ), D( 1, 6 ), A( 1, 4 ), NMAX,
02077      $                   A( 1, 5 ), RESULT, WORK, LWORK, IWORK, LOGWRK,
02078      $                   INFO )
02079             IF( INFO.NE.0 )
02080      $         WRITE( NOUT, FMT = 9980 )'DGEESX', INFO
02081          END IF
02082          WRITE( NOUT, FMT = 9973 )
02083          GO TO 10
02084 *
02085       ELSE IF( LSAMEN( 3, C3, 'DGG' ) ) THEN
02086 *
02087 *        -------------------------------------------------
02088 *        DGG:  Generalized Nonsymmetric Eigenvalue Problem
02089 *        -------------------------------------------------
02090 *        Vary the parameters
02091 *           NB    = block size
02092 *           NBMIN = minimum block size
02093 *           NS    = number of shifts
02094 *           MAXB  = minimum submatrix size
02095 *           NBCOL = minimum column dimension for blocks
02096 *
02097          MAXTYP = 26
02098          NTYPES = MIN( MAXTYP, NTYPES )
02099          CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02100          IF( TSTCHK .AND. TSTERR )
02101      $      CALL DERRGG( C3, NOUT )
02102          DO 350 I = 1, NPARMS
02103             CALL XLAENV( 1, NBVAL( I ) )
02104             CALL XLAENV( 2, NBMIN( I ) )
02105             CALL XLAENV( 4, NSVAL( I ) )
02106             CALL XLAENV( 8, MXBVAL( I ) )
02107             CALL XLAENV( 5, NBCOL( I ) )
02108 *
02109             IF( NEWSD.EQ.0 ) THEN
02110                DO 340 K = 1, 4
02111                   ISEED( K ) = IOLDSD( K )
02112   340          CONTINUE
02113             END IF
02114             WRITE( NOUT, FMT = 9996 )C3, NBVAL( I ), NBMIN( I ),
02115      $         NSVAL( I ), MXBVAL( I ), NBCOL( I )
02116             TSTDIF = .FALSE.
02117             THRSHN = 10.D0
02118             IF( TSTCHK ) THEN
02119                CALL DCHKGG( NN, NVAL, MAXTYP, DOTYPE, ISEED, THRESH,
02120      $                      TSTDIF, THRSHN, NOUT, A( 1, 1 ), NMAX,
02121      $                      A( 1, 2 ), A( 1, 3 ), A( 1, 4 ), A( 1, 5 ),
02122      $                      A( 1, 6 ), A( 1, 7 ), A( 1, 8 ), A( 1, 9 ),
02123      $                      NMAX, A( 1, 10 ), A( 1, 11 ), A( 1, 12 ),
02124      $                      D( 1, 1 ), D( 1, 2 ), D( 1, 3 ), D( 1, 4 ),
02125      $                      D( 1, 5 ), D( 1, 6 ), A( 1, 13 ),
02126      $                      A( 1, 14 ), WORK, LWORK, LOGWRK, RESULT,
02127      $                      INFO )
02128                IF( INFO.NE.0 )
02129      $            WRITE( NOUT, FMT = 9980 )'DCHKGG', INFO
02130             END IF
02131             CALL XLAENV( 1, 1 )
02132             IF( TSTDRV ) THEN
02133                CALL DDRVGG( NN, NVAL, MAXTYP, DOTYPE, ISEED, THRESH,
02134      $                      THRSHN, NOUT, A( 1, 1 ), NMAX, A( 1, 2 ),
02135      $                      A( 1, 3 ), A( 1, 4 ), A( 1, 5 ), A( 1, 6 ),
02136      $                      A( 1, 7 ), NMAX, A( 1, 8 ), D( 1, 1 ),
02137      $                      D( 1, 2 ), D( 1, 3 ), D( 1, 4 ), D( 1, 5 ),
02138      $                      D( 1, 6 ), A( 1, 13 ), A( 1, 14 ), WORK,
02139      $                      LWORK, RESULT, INFO )
02140                IF( INFO.NE.0 )
02141      $            WRITE( NOUT, FMT = 9980 )'DDRVGG', INFO
02142             END IF
02143   350    CONTINUE
02144 *
02145       ELSE IF( LSAMEN( 3, C3, 'DGS' ) ) THEN
02146 *
02147 *        -------------------------------------------------
02148 *        DGS:  Generalized Nonsymmetric Eigenvalue Problem
02149 *              DGGES (Schur form)
02150 *        -------------------------------------------------
02151 *
02152          MAXTYP = 26
02153          NTYPES = MIN( MAXTYP, NTYPES )
02154          IF( NTYPES.LE.0 ) THEN
02155             WRITE( NOUT, FMT = 9990 )C3
02156          ELSE
02157             IF( TSTERR )
02158      $         CALL DERRGG( C3, NOUT )
02159             CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02160             CALL DDRGES( NN, NVAL, MAXTYP, DOTYPE, ISEED, THRESH, NOUT,
02161      $                   A( 1, 1 ), NMAX, A( 1, 2 ), A( 1, 3 ),
02162      $                   A( 1, 4 ), A( 1, 7 ), NMAX, A( 1, 8 ),
02163      $                   D( 1, 1 ), D( 1, 2 ), D( 1, 3 ), WORK, LWORK,
02164      $                   RESULT, LOGWRK, INFO )
02165 *
02166             IF( INFO.NE.0 )
02167      $         WRITE( NOUT, FMT = 9980 )'DDRGES', INFO
02168          END IF
02169          WRITE( NOUT, FMT = 9973 )
02170          GO TO 10
02171 *
02172       ELSE IF( DGX ) THEN
02173 *
02174 *        -------------------------------------------------
02175 *        DGX:  Generalized Nonsymmetric Eigenvalue Problem
02176 *              DGGESX (Schur form and condition numbers)
02177 *        -------------------------------------------------
02178 *
02179          MAXTYP = 5
02180          NTYPES = MAXTYP
02181          IF( NN.LT.0 ) THEN
02182             WRITE( NOUT, FMT = 9990 )C3
02183          ELSE
02184             IF( TSTERR )
02185      $         CALL DERRGG( C3, NOUT )
02186             CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02187             CALL XLAENV( 5, 2 )
02188             CALL DDRGSX( NN, NCMAX, THRESH, NIN, NOUT, A( 1, 1 ), NMAX,
02189      $                   A( 1, 2 ), A( 1, 3 ), A( 1, 4 ), A( 1, 5 ),
02190      $                   A( 1, 6 ), D( 1, 1 ), D( 1, 2 ), D( 1, 3 ),
02191      $                   C( 1, 1 ), NCMAX*NCMAX, A( 1, 12 ), WORK,
02192      $                   LWORK, IWORK, LIWORK, LOGWRK, INFO )
02193             IF( INFO.NE.0 )
02194      $         WRITE( NOUT, FMT = 9980 )'DDRGSX', INFO
02195          END IF
02196          WRITE( NOUT, FMT = 9973 )
02197          GO TO 10
02198 *
02199       ELSE IF( LSAMEN( 3, C3, 'DGV' ) ) THEN
02200 *
02201 *        -------------------------------------------------
02202 *        DGV:  Generalized Nonsymmetric Eigenvalue Problem
02203 *              DGGEV (Eigenvalue/vector form)
02204 *        -------------------------------------------------
02205 *
02206          MAXTYP = 26
02207          NTYPES = MIN( MAXTYP, NTYPES )
02208          IF( NTYPES.LE.0 ) THEN
02209             WRITE( NOUT, FMT = 9990 )C3
02210          ELSE
02211             IF( TSTERR )
02212      $         CALL DERRGG( C3, NOUT )
02213             CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02214             CALL DDRGEV( NN, NVAL, MAXTYP, DOTYPE, ISEED, THRESH, NOUT,
02215      $                   A( 1, 1 ), NMAX, A( 1, 2 ), A( 1, 3 ),
02216      $                   A( 1, 4 ), A( 1, 7 ), NMAX, A( 1, 8 ),
02217      $                   A( 1, 9 ), NMAX, D( 1, 1 ), D( 1, 2 ),
02218      $                   D( 1, 3 ), D( 1, 4 ), D( 1, 5 ), D( 1, 6 ),
02219      $                   WORK, LWORK, RESULT, INFO )
02220             IF( INFO.NE.0 )
02221      $         WRITE( NOUT, FMT = 9980 )'DDRGEV', INFO
02222          END IF
02223          WRITE( NOUT, FMT = 9973 )
02224          GO TO 10
02225 *
02226       ELSE IF( DXV ) THEN
02227 *
02228 *        -------------------------------------------------
02229 *        DXV:  Generalized Nonsymmetric Eigenvalue Problem
02230 *              DGGEVX (eigenvalue/vector with condition numbers)
02231 *        -------------------------------------------------
02232 *
02233          MAXTYP = 2
02234          NTYPES = MAXTYP
02235          IF( NN.LT.0 ) THEN
02236             WRITE( NOUT, FMT = 9990 )C3
02237          ELSE
02238             IF( TSTERR )
02239      $         CALL DERRGG( C3, NOUT )
02240             CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02241             CALL DDRGVX( NN, THRESH, NIN, NOUT, A( 1, 1 ), NMAX,
02242      $                   A( 1, 2 ), A( 1, 3 ), A( 1, 4 ), D( 1, 1 ),
02243      $                   D( 1, 2 ), D( 1, 3 ), A( 1, 5 ), A( 1, 6 ),
02244      $                   IWORK( 1 ), IWORK( 2 ), D( 1, 4 ), D( 1, 5 ),
02245      $                   D( 1, 6 ), D( 1, 7 ), D( 1, 8 ), D( 1, 9 ),
02246      $                   WORK, LWORK, IWORK( 3 ), LIWORK-2, RESULT,
02247      $                   LOGWRK, INFO )
02248 *
02249             IF( INFO.NE.0 )
02250      $         WRITE( NOUT, FMT = 9980 )'DDRGVX', INFO
02251          END IF
02252          WRITE( NOUT, FMT = 9973 )
02253          GO TO 10
02254 *
02255       ELSE IF( LSAMEN( 3, C3, 'DSB' ) ) THEN
02256 *
02257 *        ------------------------------
02258 *        DSB:  Symmetric Band Reduction
02259 *        ------------------------------
02260 *
02261          MAXTYP = 15
02262          NTYPES = MIN( MAXTYP, NTYPES )
02263          CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02264          IF( TSTERR )
02265      $      CALL DERRST( 'DSB', NOUT )
02266          CALL DCHKSB( NN, NVAL, NK, KVAL, MAXTYP, DOTYPE, ISEED, THRESH,
02267      $                NOUT, A( 1, 1 ), NMAX, D( 1, 1 ), D( 1, 2 ),
02268      $                A( 1, 2 ), NMAX, WORK, LWORK, RESULT, INFO )
02269          IF( INFO.NE.0 )
02270      $      WRITE( NOUT, FMT = 9980 )'DCHKSB', INFO
02271 *
02272       ELSE IF( LSAMEN( 3, C3, 'DBB' ) ) THEN
02273 *
02274 *        ------------------------------
02275 *        DBB:  General Band Reduction
02276 *        ------------------------------
02277 *
02278          MAXTYP = 15
02279          NTYPES = MIN( MAXTYP, NTYPES )
02280          CALL ALAREQ( C3, NTYPES, DOTYPE, MAXTYP, NIN, NOUT )
02281          DO 370 I = 1, NPARMS
02282             NRHS = NSVAL( I )
02283 *
02284             IF( NEWSD.EQ.0 ) THEN
02285                DO 360 K = 1, 4
02286                   ISEED( K ) = IOLDSD( K )
02287   360          CONTINUE
02288             END IF
02289             WRITE( NOUT, FMT = 9966 )C3, NRHS
02290             CALL DCHKBB( NN, MVAL, NVAL, NK, KVAL, MAXTYP, DOTYPE, NRHS,
02291      $                   ISEED, THRESH, NOUT, A( 1, 1 ), NMAX,
02292      $                   A( 1, 2 ), 2*NMAX, D( 1, 1 ), D( 1, 2 ),
02293      $                   A( 1, 4 ), NMAX, A( 1, 5 ), NMAX, A( 1, 6 ),
02294      $                   NMAX, A( 1, 7 ), WORK, LWORK, RESULT, INFO )
02295             IF( INFO.NE.0 )
02296      $         WRITE( NOUT, FMT = 9980 )'DCHKBB', INFO
02297   370    CONTINUE
02298 *
02299       ELSE IF( LSAMEN( 3, C3, 'GLM' ) ) THEN
02300 *
02301 *        -----------------------------------------
02302 *        GLM:  Generalized Linear Regression Model
02303 *        -----------------------------------------
02304 *
02305          CALL XLAENV( 1, 1 )
02306          IF( TSTERR )
02307      $      CALL DERRGG( 'GLM', NOUT )
02308          CALL DCKGLM( NN, MVAL, PVAL, NVAL, NTYPES, ISEED, THRESH, NMAX,
02309      $                A( 1, 1 ), A( 1, 2 ), B( 1, 1 ), B( 1, 2 ), X,
02310      $                WORK, D( 1, 1 ), NIN, NOUT, INFO )
02311          IF( INFO.NE.0 )
02312      $      WRITE( NOUT, FMT = 9980 )'DCKGLM', INFO
02313 *
02314       ELSE IF( LSAMEN( 3, C3, 'GQR' ) ) THEN
02315 *
02316 *        ------------------------------------------
02317 *        GQR:  Generalized QR and RQ factorizations
02318 *        ------------------------------------------
02319 *
02320          CALL XLAENV( 1, 1 )
02321          IF( TSTERR )
02322      $      CALL DERRGG( 'GQR', NOUT )
02323          CALL DCKGQR( NN, MVAL, NN, PVAL, NN, NVAL, NTYPES, ISEED,
02324      $                THRESH, NMAX, A( 1, 1 ), A( 1, 2 ), A( 1, 3 ),
02325      $                A( 1, 4 ), TAUA, B( 1, 1 ), B( 1, 2 ), B( 1, 3 ),
02326      $                B( 1, 4 ), B( 1, 5 ), TAUB, WORK, D( 1, 1 ), NIN,
02327      $                NOUT, INFO )
02328          IF( INFO.NE.0 )
02329      $      WRITE( NOUT, FMT = 9980 )'DCKGQR', INFO
02330 *
02331       ELSE IF( LSAMEN( 3, C3, 'GSV' ) ) THEN
02332 *
02333 *        ----------------------------------------------
02334 *        GSV:  Generalized Singular Value Decomposition
02335 *        ----------------------------------------------
02336 *
02337          IF( TSTERR )
02338      $      CALL DERRGG( 'GSV', NOUT )
02339          CALL DCKGSV( NN, MVAL, PVAL, NVAL, NTYPES, ISEED, THRESH, NMAX,
02340      $                A( 1, 1 ), A( 1, 2 ), B( 1, 1 ), B( 1, 2 ),
02341      $                A( 1, 3 ), B( 1, 3 ), A( 1, 4 ), TAUA, TAUB,
02342      $                B( 1, 4 ), IWORK, WORK, D( 1, 1 ), NIN, NOUT,
02343      $                INFO )
02344          IF( INFO.NE.0 )
02345      $      WRITE( NOUT, FMT = 9980 )'DCKGSV', INFO
02346 *
02347       ELSE IF( LSAMEN( 3, C3, 'CSD' ) ) THEN
02348 *
02349 *        ----------------------------------------------
02350 *        CSD:  CS Decomposition
02351 *        ----------------------------------------------
02352 *
02353          CALL XLAENV(1,1)
02354          IF( TSTERR )
02355      $      CALL DERRGG( 'CSD', NOUT )
02356          CALL DCKCSD( NN, MVAL, PVAL, NVAL, NTYPES, ISEED, THRESH, NMAX,
02357      $                A( 1, 1 ), A( 1, 2 ), A( 1, 3 ), A( 1, 4 ),
02358      $                A( 1, 5 ), A( 1, 6 ), A( 1, 7 ), IWORK, WORK,
02359      $                D( 1, 1 ), NIN, NOUT, INFO )
02360          IF( INFO.NE.0 )
02361      $      WRITE( NOUT, FMT = 9980 )'DCKCSD', INFO
02362 *
02363       ELSE IF( LSAMEN( 3, C3, 'LSE' ) ) THEN
02364 *
02365 *        --------------------------------------
02366 *        LSE:  Constrained Linear Least Squares
02367 *        --------------------------------------
02368 *
02369          CALL XLAENV( 1, 1 )
02370          IF( TSTERR )
02371      $      CALL DERRGG( 'LSE', NOUT )
02372          CALL DCKLSE( NN, MVAL, PVAL, NVAL, NTYPES, ISEED, THRESH, NMAX,
02373      $                A( 1, 1 ), A( 1, 2 ), B( 1, 1 ), B( 1, 2 ), X,
02374      $                WORK, D( 1, 1 ), NIN, NOUT, INFO )
02375          IF( INFO.NE.0 )
02376      $      WRITE( NOUT, FMT = 9980 )'DCKLSE', INFO
02377 *
02378       ELSE
02379          WRITE( NOUT, FMT = * )
02380          WRITE( NOUT, FMT = * )
02381          WRITE( NOUT, FMT = 9992 )C3
02382       END IF
02383       IF( .NOT.( DGX .OR. DXV ) )
02384      $   GO TO 190
02385   380 CONTINUE
02386       WRITE( NOUT, FMT = 9994 )
02387       S2 = DSECND( )
02388       WRITE( NOUT, FMT = 9993 )S2 - S1
02389 *
02390  9999 FORMAT( / ' Execution not attempted due to input errors' )
02391  9997 FORMAT( / / 1X, A3, ':  NB =', I4, ', NBMIN =', I4, ', NX =', I4 )
02392  9996 FORMAT( / / 1X, A3, ':  NB =', I4, ', NBMIN =', I4, ', NS =', I4,
02393      $      ', MAXB =', I4, ', NBCOL =', I4 )
02394  9995 FORMAT( / / 1X, A3, ':  NB =', I4, ', NBMIN =', I4, ', NX =', I4,
02395      $      ', NRHS =', I4 )
02396  9994 FORMAT( / / ' End of tests' )
02397  9993 FORMAT( ' Total time used = ', F12.2, ' seconds', / )
02398  9992 FORMAT( 1X, A3, ':  Unrecognized path name' )
02399  9991 FORMAT( / / ' *** Invalid integer value in column ', I2,
02400      $      ' of input', ' line:', / A79 )
02401  9990 FORMAT( / / 1X, A3, ' routines were not tested' )
02402  9989 FORMAT( ' Invalid input value: ', A, '=', I6, '; must be >=',
02403      $      I6 )
02404  9988 FORMAT( ' Invalid input value: ', A, '=', I6, '; must be <=',
02405      $      I6 )
02406  9987 FORMAT( ' Tests of the Nonsymmetric Eigenvalue Problem routines' )
02407  9986 FORMAT( ' Tests of the Symmetric Eigenvalue Problem routines' )
02408  9985 FORMAT( ' Tests of the Singular Value Decomposition routines' )
02409  9984 FORMAT( / ' The following parameter values will be used:' )
02410  9983 FORMAT( 4X, A, 10I6, / 10X, 10I6 )
02411  9982 FORMAT( / ' Routines pass computational tests if test ratio is ',
02412      $      'less than', F8.2, / )
02413  9981 FORMAT( ' Relative machine ', A, ' is taken to be', D16.6 )
02414  9980 FORMAT( ' *** Error code from ', A, ' = ', I4 )
02415  9979 FORMAT( / ' Tests of the Nonsymmetric Eigenvalue Problem Driver',
02416      $      / '    DGEEV (eigenvalues and eigevectors)' )
02417  9978 FORMAT( / ' Tests of the Nonsymmetric Eigenvalue Problem Driver',
02418      $      / '    DGEES (Schur form)' )
02419  9977 FORMAT( / ' Tests of the Nonsymmetric Eigenvalue Problem Expert',
02420      $      ' Driver', / '    DGEEVX (eigenvalues, eigenvectors and',
02421      $      ' condition numbers)' )
02422  9976 FORMAT( / ' Tests of the Nonsymmetric Eigenvalue Problem Expert',
02423      $      ' Driver', / '    DGEESX (Schur form and condition',
02424      $      ' numbers)' )
02425  9975 FORMAT( / ' Tests of the Generalized Nonsymmetric Eigenvalue ',
02426      $      'Problem routines' )
02427  9974 FORMAT( ' Tests of DSBTRD', / ' (reduction of a symmetric band ',
02428      $      'matrix to tridiagonal form)' )
02429  9973 FORMAT( / 1X, 71( '-' ) )
02430  9972 FORMAT( / ' LAPACK VERSION ', I1, '.', I1, '.', I1 )
02431  9971 FORMAT( / ' Tests of the Generalized Linear Regression Model ',
02432      $      'routines' )
02433  9970 FORMAT( / ' Tests of the Generalized QR and RQ routines' )
02434  9969 FORMAT( / ' Tests of the Generalized Singular Value',
02435      $      ' Decomposition routines' )
02436  9968 FORMAT( / ' Tests of the Linear Least Squares routines' )
02437  9967 FORMAT( ' Tests of DGBBRD', / ' (reduction of a general band ',
02438      $      'matrix to real bidiagonal form)' )
02439  9966 FORMAT( / / 1X, A3, ':  NRHS =', I4 )
02440  9965 FORMAT( / ' Tests of the Generalized Nonsymmetric Eigenvalue ',
02441      $      'Problem Expert Driver DGGESX' )
02442  9964 FORMAT( / ' Tests of the Generalized Nonsymmetric Eigenvalue ',
02443      $      'Problem Driver DGGES' )
02444  9963 FORMAT( / ' Tests of the Generalized Nonsymmetric Eigenvalue ',
02445      $      'Problem Driver DGGEV' )
02446  9962 FORMAT( / ' Tests of the Generalized Nonsymmetric Eigenvalue ',
02447      $      'Problem Expert Driver DGGEVX' )
02448  9961 FORMAT( / / 1X, A3, ':  NB =', I4, ', NBMIN =', I4, ', NX =', I4,
02449      $      ', INMIN=', I4, 
02450      $      ', INWIN =', I4, ', INIBL =', I4, ', ISHFTS =', I4,
02451      $      ', IACC22 =', I4)
02452  9960 FORMAT( / ' Tests of the CS Decomposition routines' )
02453 *
02454 *     End of DCHKEE
02455 *
02456       END
 All Files Functions