LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
sdrvev.f
Go to the documentation of this file.
00001 *> \brief \b SDRVEV
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *  Definition:
00009 *  ===========
00010 *
00011 *       SUBROUTINE SDRVEV( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
00012 *                          NOUNIT, A, LDA, H, WR, WI, WR1, WI1, VL, LDVL,
00013 *                          VR, LDVR, LRE, LDLRE, RESULT, WORK, NWORK,
00014 *                          IWORK, INFO )
00015 * 
00016 *       .. Scalar Arguments ..
00017 *       INTEGER            INFO, LDA, LDLRE, LDVL, LDVR, NOUNIT, NSIZES,
00018 *      $                   NTYPES, NWORK
00019 *       REAL               THRESH
00020 *       ..
00021 *       .. Array Arguments ..
00022 *       LOGICAL            DOTYPE( * )
00023 *       INTEGER            ISEED( 4 ), IWORK( * ), NN( * )
00024 *       REAL               A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
00025 *      $                   RESULT( 7 ), VL( LDVL, * ), VR( LDVR, * ),
00026 *      $                   WI( * ), WI1( * ), WORK( * ), WR( * ), WR1( * )
00027 *       ..
00028 *  
00029 *
00030 *> \par Purpose:
00031 *  =============
00032 *>
00033 *> \verbatim
00034 *>
00035 *>    SDRVEV  checks the nonsymmetric eigenvalue problem driver SGEEV.
00036 *>
00037 *>    When SDRVEV is called, a number of matrix "sizes" ("n's") and a
00038 *>    number of matrix "types" are specified.  For each size ("n")
00039 *>    and each type of matrix, one matrix will be generated and used
00040 *>    to test the nonsymmetric eigenroutines.  For each matrix, 7
00041 *>    tests will be performed:
00042 *>
00043 *>    (1)     | A * VR - VR * W | / ( n |A| ulp )
00044 *>
00045 *>      Here VR is the matrix of unit right eigenvectors.
00046 *>      W is a block diagonal matrix, with a 1x1 block for each
00047 *>      real eigenvalue and a 2x2 block for each complex conjugate
00048 *>      pair.  If eigenvalues j and j+1 are a complex conjugate pair,
00049 *>      so WR(j) = WR(j+1) = wr and WI(j) = - WI(j+1) = wi, then the
00050 *>      2 x 2 block corresponding to the pair will be:
00051 *>
00052 *>              (  wr  wi  )
00053 *>              ( -wi  wr  )
00054 *>
00055 *>      Such a block multiplying an n x 2 matrix  ( ur ui ) on the
00056 *>      right will be the same as multiplying  ur + i*ui  by  wr + i*wi.
00057 *>
00058 *>    (2)     | A**H * VL - VL * W**H | / ( n |A| ulp )
00059 *>
00060 *>      Here VL is the matrix of unit left eigenvectors, A**H is the
00061 *>      conjugate transpose of A, and W is as above.
00062 *>
00063 *>    (3)     | |VR(i)| - 1 | / ulp and whether largest component real
00064 *>
00065 *>      VR(i) denotes the i-th column of VR.
00066 *>
00067 *>    (4)     | |VL(i)| - 1 | / ulp and whether largest component real
00068 *>
00069 *>      VL(i) denotes the i-th column of VL.
00070 *>
00071 *>    (5)     W(full) = W(partial)
00072 *>
00073 *>      W(full) denotes the eigenvalues computed when both VR and VL
00074 *>      are also computed, and W(partial) denotes the eigenvalues
00075 *>      computed when only W, only W and VR, or only W and VL are
00076 *>      computed.
00077 *>
00078 *>    (6)     VR(full) = VR(partial)
00079 *>
00080 *>      VR(full) denotes the right eigenvectors computed when both VR
00081 *>      and VL are computed, and VR(partial) denotes the result
00082 *>      when only VR is computed.
00083 *>
00084 *>     (7)     VL(full) = VL(partial)
00085 *>
00086 *>      VL(full) denotes the left eigenvectors computed when both VR
00087 *>      and VL are also computed, and VL(partial) denotes the result
00088 *>      when only VL is computed.
00089 *>
00090 *>    The "sizes" are specified by an array NN(1:NSIZES); the value of
00091 *>    each element NN(j) specifies one size.
00092 *>    The "types" are specified by a logical array DOTYPE( 1:NTYPES );
00093 *>    if DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
00094 *>    Currently, the list of possible types is:
00095 *>
00096 *>    (1)  The zero matrix.
00097 *>    (2)  The identity matrix.
00098 *>    (3)  A (transposed) Jordan block, with 1's on the diagonal.
00099 *>
00100 *>    (4)  A diagonal matrix with evenly spaced entries
00101 *>         1, ..., ULP  and random signs.
00102 *>         (ULP = (first number larger than 1) - 1 )
00103 *>    (5)  A diagonal matrix with geometrically spaced entries
00104 *>         1, ..., ULP  and random signs.
00105 *>    (6)  A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
00106 *>         and random signs.
00107 *>
00108 *>    (7)  Same as (4), but multiplied by a constant near
00109 *>         the overflow threshold
00110 *>    (8)  Same as (4), but multiplied by a constant near
00111 *>         the underflow threshold
00112 *>
00113 *>    (9)  A matrix of the form  U' T U, where U is orthogonal and
00114 *>         T has evenly spaced entries 1, ..., ULP with random signs
00115 *>         on the diagonal and random O(1) entries in the upper
00116 *>         triangle.
00117 *>
00118 *>    (10) A matrix of the form  U' T U, where U is orthogonal and
00119 *>         T has geometrically spaced entries 1, ..., ULP with random
00120 *>         signs on the diagonal and random O(1) entries in the upper
00121 *>         triangle.
00122 *>
00123 *>    (11) A matrix of the form  U' T U, where U is orthogonal and
00124 *>         T has "clustered" entries 1, ULP,..., ULP with random
00125 *>         signs on the diagonal and random O(1) entries in the upper
00126 *>         triangle.
00127 *>
00128 *>    (12) A matrix of the form  U' T U, where U is orthogonal and
00129 *>         T has real or complex conjugate paired eigenvalues randomly
00130 *>         chosen from ( ULP, 1 ) and random O(1) entries in the upper
00131 *>         triangle.
00132 *>
00133 *>    (13) A matrix of the form  X' T X, where X has condition
00134 *>         SQRT( ULP ) and T has evenly spaced entries 1, ..., ULP
00135 *>         with random signs on the diagonal and random O(1) entries
00136 *>         in the upper triangle.
00137 *>
00138 *>    (14) A matrix of the form  X' T X, where X has condition
00139 *>         SQRT( ULP ) and T has geometrically spaced entries
00140 *>         1, ..., ULP with random signs on the diagonal and random
00141 *>         O(1) entries in the upper triangle.
00142 *>
00143 *>    (15) A matrix of the form  X' T X, where X has condition
00144 *>         SQRT( ULP ) and T has "clustered" entries 1, ULP,..., ULP
00145 *>         with random signs on the diagonal and random O(1) entries
00146 *>         in the upper triangle.
00147 *>
00148 *>    (16) A matrix of the form  X' T X, where X has condition
00149 *>         SQRT( ULP ) and T has real or complex conjugate paired
00150 *>         eigenvalues randomly chosen from ( ULP, 1 ) and random
00151 *>         O(1) entries in the upper triangle.
00152 *>
00153 *>    (17) Same as (16), but multiplied by a constant
00154 *>         near the overflow threshold
00155 *>    (18) Same as (16), but multiplied by a constant
00156 *>         near the underflow threshold
00157 *>
00158 *>    (19) Nonsymmetric matrix with random entries chosen from (-1,1).
00159 *>         If N is at least 4, all entries in first two rows and last
00160 *>         row, and first column and last two columns are zero.
00161 *>    (20) Same as (19), but multiplied by a constant
00162 *>         near the overflow threshold
00163 *>    (21) Same as (19), but multiplied by a constant
00164 *>         near the underflow threshold
00165 *> \endverbatim
00166 *
00167 *  Arguments:
00168 *  ==========
00169 *
00170 *> \param[in] NSIZES
00171 *> \verbatim
00172 *>          NSIZES is INTEGER
00173 *>          The number of sizes of matrices to use.  If it is zero,
00174 *>          SDRVEV does nothing.  It must be at least zero.
00175 *> \endverbatim
00176 *>
00177 *> \param[in] NN
00178 *> \verbatim
00179 *>          NN is INTEGER array, dimension (NSIZES)
00180 *>          An array containing the sizes to be used for the matrices.
00181 *>          Zero values will be skipped.  The values must be at least
00182 *>          zero.
00183 *> \endverbatim
00184 *>
00185 *> \param[in] NTYPES
00186 *> \verbatim
00187 *>          NTYPES is INTEGER
00188 *>          The number of elements in DOTYPE.   If it is zero, SDRVEV
00189 *>          does nothing.  It must be at least zero.  If it is MAXTYP+1
00190 *>          and NSIZES is 1, then an additional type, MAXTYP+1 is
00191 *>          defined, which is to use whatever matrix is in A.  This
00192 *>          is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
00193 *>          DOTYPE(MAXTYP+1) is .TRUE. .
00194 *> \endverbatim
00195 *>
00196 *> \param[in] DOTYPE
00197 *> \verbatim
00198 *>          DOTYPE is LOGICAL array, dimension (NTYPES)
00199 *>          If DOTYPE(j) is .TRUE., then for each size in NN a
00200 *>          matrix of that size and of type j will be generated.
00201 *>          If NTYPES is smaller than the maximum number of types
00202 *>          defined (PARAMETER MAXTYP), then types NTYPES+1 through
00203 *>          MAXTYP will not be generated.  If NTYPES is larger
00204 *>          than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
00205 *>          will be ignored.
00206 *> \endverbatim
00207 *>
00208 *> \param[in,out] ISEED
00209 *> \verbatim
00210 *>          ISEED is INTEGER array, dimension (4)
00211 *>          On entry ISEED specifies the seed of the random number
00212 *>          generator. The array elements should be between 0 and 4095;
00213 *>          if not they will be reduced mod 4096.  Also, ISEED(4) must
00214 *>          be odd.  The random number generator uses a linear
00215 *>          congruential sequence limited to small integers, and so
00216 *>          should produce machine independent random numbers. The
00217 *>          values of ISEED are changed on exit, and can be used in the
00218 *>          next call to SDRVEV to continue the same random number
00219 *>          sequence.
00220 *> \endverbatim
00221 *>
00222 *> \param[in] THRESH
00223 *> \verbatim
00224 *>          THRESH is REAL
00225 *>          A test will count as "failed" if the "error", computed as
00226 *>          described above, exceeds THRESH.  Note that the error
00227 *>          is scaled to be O(1), so THRESH should be a reasonably
00228 *>          small multiple of 1, e.g., 10 or 100.  In particular,
00229 *>          it should not depend on the precision (single vs. double)
00230 *>          or the size of the matrix.  It must be at least zero.
00231 *> \endverbatim
00232 *>
00233 *> \param[in] NOUNIT
00234 *> \verbatim
00235 *>          NOUNIT is INTEGER
00236 *>          The FORTRAN unit number for printing out error messages
00237 *>          (e.g., if a routine returns INFO not equal to 0.)
00238 *> \endverbatim
00239 *>
00240 *> \param[out] A
00241 *> \verbatim
00242 *>          A is REAL array, dimension (LDA, max(NN))
00243 *>          Used to hold the matrix whose eigenvalues are to be
00244 *>          computed.  On exit, A contains the last matrix actually used.
00245 *> \endverbatim
00246 *>
00247 *> \param[in] LDA
00248 *> \verbatim
00249 *>          LDA is INTEGER
00250 *>          The leading dimension of A, and H. LDA must be at
00251 *>          least 1 and at least max(NN).
00252 *> \endverbatim
00253 *>
00254 *> \param[out] H
00255 *> \verbatim
00256 *>          H is REAL array, dimension (LDA, max(NN))
00257 *>          Another copy of the test matrix A, modified by SGEEV.
00258 *> \endverbatim
00259 *>
00260 *> \param[out] WR
00261 *> \verbatim
00262 *>          WR is REAL array, dimension (max(NN))
00263 *> \endverbatim
00264 *>
00265 *> \param[out] WI
00266 *> \verbatim
00267 *>          WI is REAL array, dimension (max(NN))
00268 *>
00269 *>          The real and imaginary parts of the eigenvalues of A.
00270 *>          On exit, WR + WI*i are the eigenvalues of the matrix in A.
00271 *> \endverbatim
00272 *>
00273 *> \param[out] WR1
00274 *> \verbatim
00275 *>          WR1 is REAL array, dimension (max(NN))
00276 *> \endverbatim
00277 *>
00278 *> \param[out] WI1
00279 *> \verbatim
00280 *>          WI1 is REAL array, dimension (max(NN))
00281 *>
00282 *>          Like WR, WI, these arrays contain the eigenvalues of A,
00283 *>          but those computed when SGEEV only computes a partial
00284 *>          eigendecomposition, i.e. not the eigenvalues and left
00285 *>          and right eigenvectors.
00286 *> \endverbatim
00287 *>
00288 *> \param[out] VL
00289 *> \verbatim
00290 *>          VL is REAL array, dimension (LDVL, max(NN))
00291 *>          VL holds the computed left eigenvectors.
00292 *> \endverbatim
00293 *>
00294 *> \param[in] LDVL
00295 *> \verbatim
00296 *>          LDVL is INTEGER
00297 *>          Leading dimension of VL. Must be at least max(1,max(NN)).
00298 *> \endverbatim
00299 *>
00300 *> \param[out] VR
00301 *> \verbatim
00302 *>          VR is REAL array, dimension (LDVR, max(NN))
00303 *>          VR holds the computed right eigenvectors.
00304 *> \endverbatim
00305 *>
00306 *> \param[in] LDVR
00307 *> \verbatim
00308 *>          LDVR is INTEGER
00309 *>          Leading dimension of VR. Must be at least max(1,max(NN)).
00310 *> \endverbatim
00311 *>
00312 *> \param[out] LRE
00313 *> \verbatim
00314 *>          LRE is REAL array, dimension (LDLRE,max(NN))
00315 *>          LRE holds the computed right or left eigenvectors.
00316 *> \endverbatim
00317 *>
00318 *> \param[in] LDLRE
00319 *> \verbatim
00320 *>          LDLRE is INTEGER
00321 *>          Leading dimension of LRE. Must be at least max(1,max(NN)).
00322 *> \endverbatim
00323 *>
00324 *> \param[out] RESULT
00325 *> \verbatim
00326 *>          RESULT is REAL array, dimension (7)
00327 *>          The values computed by the seven tests described above.
00328 *>          The values are currently limited to 1/ulp, to avoid overflow.
00329 *> \endverbatim
00330 *>
00331 *> \param[out] WORK
00332 *> \verbatim
00333 *>          WORK is REAL array, dimension (NWORK)
00334 *> \endverbatim
00335 *>
00336 *> \param[in] NWORK
00337 *> \verbatim
00338 *>          NWORK is INTEGER
00339 *>          The number of entries in WORK.  This must be at least
00340 *>          5*NN(j)+2*NN(j)**2 for all j.
00341 *> \endverbatim
00342 *>
00343 *> \param[out] IWORK
00344 *> \verbatim
00345 *>          IWORK is INTEGER array, dimension (max(NN))
00346 *> \endverbatim
00347 *>
00348 *> \param[out] INFO
00349 *> \verbatim
00350 *>          INFO is INTEGER
00351 *>          If 0, then everything ran OK.
00352 *>           -1: NSIZES < 0
00353 *>           -2: Some NN(j) < 0
00354 *>           -3: NTYPES < 0
00355 *>           -6: THRESH < 0
00356 *>           -9: LDA < 1 or LDA < NMAX, where NMAX is max( NN(j) ).
00357 *>          -16: LDVL < 1 or LDVL < NMAX, where NMAX is max( NN(j) ).
00358 *>          -18: LDVR < 1 or LDVR < NMAX, where NMAX is max( NN(j) ).
00359 *>          -20: LDLRE < 1 or LDLRE < NMAX, where NMAX is max( NN(j) ).
00360 *>          -23: NWORK too small.
00361 *>          If  SLATMR, SLATMS, SLATME or SGEEV returns an error code,
00362 *>              the absolute value of it is returned.
00363 *>
00364 *>-----------------------------------------------------------------------
00365 *>
00366 *>     Some Local Variables and Parameters:
00367 *>     ---- ----- --------- --- ----------
00368 *>
00369 *>     ZERO, ONE       Real 0 and 1.
00370 *>     MAXTYP          The number of types defined.
00371 *>     NMAX            Largest value in NN.
00372 *>     NERRS           The number of tests which have exceeded THRESH
00373 *>     COND, CONDS,
00374 *>     IMODE           Values to be passed to the matrix generators.
00375 *>     ANORM           Norm of A; passed to matrix generators.
00376 *>
00377 *>     OVFL, UNFL      Overflow and underflow thresholds.
00378 *>     ULP, ULPINV     Finest relative precision and its inverse.
00379 *>     RTULP, RTULPI   Square roots of the previous 4 values.
00380 *>
00381 *>             The following four arrays decode JTYPE:
00382 *>     KTYPE(j)        The general type (1-10) for type "j".
00383 *>     KMODE(j)        The MODE value to be passed to the matrix
00384 *>                     generator for type "j".
00385 *>     KMAGN(j)        The order of magnitude ( O(1),
00386 *>                     O(overflow^(1/2) ), O(underflow^(1/2) )
00387 *>     KCONDS(j)       Selectw whether CONDS is to be 1 or
00388 *>                     1/sqrt(ulp).  (0 means irrelevant.)
00389 *> \endverbatim
00390 *
00391 *  Authors:
00392 *  ========
00393 *
00394 *> \author Univ. of Tennessee 
00395 *> \author Univ. of California Berkeley 
00396 *> \author Univ. of Colorado Denver 
00397 *> \author NAG Ltd. 
00398 *
00399 *> \date November 2011
00400 *
00401 *> \ingroup single_eig
00402 *
00403 *  =====================================================================
00404       SUBROUTINE SDRVEV( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
00405      $                   NOUNIT, A, LDA, H, WR, WI, WR1, WI1, VL, LDVL,
00406      $                   VR, LDVR, LRE, LDLRE, RESULT, WORK, NWORK,
00407      $                   IWORK, INFO )
00408 *
00409 *  -- LAPACK test routine (version 3.4.0) --
00410 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00411 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00412 *     November 2011
00413 *
00414 *     .. Scalar Arguments ..
00415       INTEGER            INFO, LDA, LDLRE, LDVL, LDVR, NOUNIT, NSIZES,
00416      $                   NTYPES, NWORK
00417       REAL               THRESH
00418 *     ..
00419 *     .. Array Arguments ..
00420       LOGICAL            DOTYPE( * )
00421       INTEGER            ISEED( 4 ), IWORK( * ), NN( * )
00422       REAL               A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
00423      $                   RESULT( 7 ), VL( LDVL, * ), VR( LDVR, * ),
00424      $                   WI( * ), WI1( * ), WORK( * ), WR( * ), WR1( * )
00425 *     ..
00426 *
00427 *  =====================================================================
00428 *
00429 *     .. Parameters ..
00430       REAL               ZERO, ONE
00431       PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0 )
00432       REAL               TWO
00433       PARAMETER          ( TWO = 2.0E0 )
00434       INTEGER            MAXTYP
00435       PARAMETER          ( MAXTYP = 21 )
00436 *     ..
00437 *     .. Local Scalars ..
00438       LOGICAL            BADNN
00439       CHARACTER*3        PATH
00440       INTEGER            IINFO, IMODE, ITYPE, IWK, J, JCOL, JJ, JSIZE,
00441      $                   JTYPE, MTYPES, N, NERRS, NFAIL, NMAX,
00442      $                   NNWORK, NTEST, NTESTF, NTESTT
00443       REAL               ANORM, COND, CONDS, OVFL, RTULP, RTULPI, TNRM,
00444      $                   ULP, ULPINV, UNFL, VMX, VRMX, VTST
00445 *     ..
00446 *     .. Local Arrays ..
00447       CHARACTER          ADUMMA( 1 )
00448       INTEGER            IDUMMA( 1 ), IOLDSD( 4 ), KCONDS( MAXTYP ),
00449      $                   KMAGN( MAXTYP ), KMODE( MAXTYP ),
00450      $                   KTYPE( MAXTYP )
00451       REAL               DUM( 1 ), RES( 2 )
00452 *     ..
00453 *     .. External Functions ..
00454       REAL               SLAMCH, SLAPY2, SNRM2
00455       EXTERNAL           SLAMCH, SLAPY2, SNRM2
00456 *     ..
00457 *     .. External Subroutines ..
00458       EXTERNAL           SGEEV, SGET22, SLABAD, SLACPY, SLASUM, SLATME,
00459      $                   SLATMR, SLATMS, SLASET, XERBLA
00460 *     ..
00461 *     .. Intrinsic Functions ..
00462       INTRINSIC          ABS, MAX, MIN, SQRT
00463 *     ..
00464 *     .. Data statements ..
00465       DATA               KTYPE / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
00466       DATA               KMAGN / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
00467      $                   3, 1, 2, 3 /
00468       DATA               KMODE / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
00469      $                   1, 5, 5, 5, 4, 3, 1 /
00470       DATA               KCONDS / 3*0, 5*0, 4*1, 6*2, 3*0 /
00471 *     ..
00472 *     .. Executable Statements ..
00473 *
00474       PATH( 1: 1 ) = 'Single precision'
00475       PATH( 2: 3 ) = 'EV'
00476 *
00477 *     Check for errors
00478 *
00479       NTESTT = 0
00480       NTESTF = 0
00481       INFO = 0
00482 *
00483 *     Important constants
00484 *
00485       BADNN = .FALSE.
00486       NMAX = 0
00487       DO 10 J = 1, NSIZES
00488          NMAX = MAX( NMAX, NN( J ) )
00489          IF( NN( J ).LT.0 )
00490      $      BADNN = .TRUE.
00491    10 CONTINUE
00492 *
00493 *     Check for errors
00494 *
00495       IF( NSIZES.LT.0 ) THEN
00496          INFO = -1
00497       ELSE IF( BADNN ) THEN
00498          INFO = -2
00499       ELSE IF( NTYPES.LT.0 ) THEN
00500          INFO = -3
00501       ELSE IF( THRESH.LT.ZERO ) THEN
00502          INFO = -6
00503       ELSE IF( NOUNIT.LE.0 ) THEN
00504          INFO = -7
00505       ELSE IF( LDA.LT.1 .OR. LDA.LT.NMAX ) THEN
00506          INFO = -9
00507       ELSE IF( LDVL.LT.1 .OR. LDVL.LT.NMAX ) THEN
00508          INFO = -16
00509       ELSE IF( LDVR.LT.1 .OR. LDVR.LT.NMAX ) THEN
00510          INFO = -18
00511       ELSE IF( LDLRE.LT.1 .OR. LDLRE.LT.NMAX ) THEN
00512          INFO = -20
00513       ELSE IF( 5*NMAX+2*NMAX**2.GT.NWORK ) THEN
00514          INFO = -23
00515       END IF
00516 *
00517       IF( INFO.NE.0 ) THEN
00518          CALL XERBLA( 'SDRVEV', -INFO )
00519          RETURN
00520       END IF
00521 *
00522 *     Quick return if nothing to do
00523 *
00524       IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 )
00525      $   RETURN
00526 *
00527 *     More Important constants
00528 *
00529       UNFL = SLAMCH( 'Safe minimum' )
00530       OVFL = ONE / UNFL
00531       CALL SLABAD( UNFL, OVFL )
00532       ULP = SLAMCH( 'Precision' )
00533       ULPINV = ONE / ULP
00534       RTULP = SQRT( ULP )
00535       RTULPI = ONE / RTULP
00536 *
00537 *     Loop over sizes, types
00538 *
00539       NERRS = 0
00540 *
00541       DO 270 JSIZE = 1, NSIZES
00542          N = NN( JSIZE )
00543          IF( NSIZES.NE.1 ) THEN
00544             MTYPES = MIN( MAXTYP, NTYPES )
00545          ELSE
00546             MTYPES = MIN( MAXTYP+1, NTYPES )
00547          END IF
00548 *
00549          DO 260 JTYPE = 1, MTYPES
00550             IF( .NOT.DOTYPE( JTYPE ) )
00551      $         GO TO 260
00552 *
00553 *           Save ISEED in case of an error.
00554 *
00555             DO 20 J = 1, 4
00556                IOLDSD( J ) = ISEED( J )
00557    20       CONTINUE
00558 *
00559 *           Compute "A"
00560 *
00561 *           Control parameters:
00562 *
00563 *           KMAGN  KCONDS  KMODE        KTYPE
00564 *       =1  O(1)   1       clustered 1  zero
00565 *       =2  large  large   clustered 2  identity
00566 *       =3  small          exponential  Jordan
00567 *       =4                 arithmetic   diagonal, (w/ eigenvalues)
00568 *       =5                 random log   symmetric, w/ eigenvalues
00569 *       =6                 random       general, w/ eigenvalues
00570 *       =7                              random diagonal
00571 *       =8                              random symmetric
00572 *       =9                              random general
00573 *       =10                             random triangular
00574 *
00575             IF( MTYPES.GT.MAXTYP )
00576      $         GO TO 90
00577 *
00578             ITYPE = KTYPE( JTYPE )
00579             IMODE = KMODE( JTYPE )
00580 *
00581 *           Compute norm
00582 *
00583             GO TO ( 30, 40, 50 )KMAGN( JTYPE )
00584 *
00585    30       CONTINUE
00586             ANORM = ONE
00587             GO TO 60
00588 *
00589    40       CONTINUE
00590             ANORM = OVFL*ULP
00591             GO TO 60
00592 *
00593    50       CONTINUE
00594             ANORM = UNFL*ULPINV
00595             GO TO 60
00596 *
00597    60       CONTINUE
00598 *
00599             CALL SLASET( 'Full', LDA, N, ZERO, ZERO, A, LDA )
00600             IINFO = 0
00601             COND = ULPINV
00602 *
00603 *           Special Matrices -- Identity & Jordan block
00604 *
00605 *              Zero
00606 *
00607             IF( ITYPE.EQ.1 ) THEN
00608                IINFO = 0
00609 *
00610             ELSE IF( ITYPE.EQ.2 ) THEN
00611 *
00612 *              Identity
00613 *
00614                DO 70 JCOL = 1, N
00615                   A( JCOL, JCOL ) = ANORM
00616    70          CONTINUE
00617 *
00618             ELSE IF( ITYPE.EQ.3 ) THEN
00619 *
00620 *              Jordan Block
00621 *
00622                DO 80 JCOL = 1, N
00623                   A( JCOL, JCOL ) = ANORM
00624                   IF( JCOL.GT.1 )
00625      $               A( JCOL, JCOL-1 ) = ONE
00626    80          CONTINUE
00627 *
00628             ELSE IF( ITYPE.EQ.4 ) THEN
00629 *
00630 *              Diagonal Matrix, [Eigen]values Specified
00631 *
00632                CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND,
00633      $                      ANORM, 0, 0, 'N', A, LDA, WORK( N+1 ),
00634      $                      IINFO )
00635 *
00636             ELSE IF( ITYPE.EQ.5 ) THEN
00637 *
00638 *              Symmetric, eigenvalues specified
00639 *
00640                CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND,
00641      $                      ANORM, N, N, 'N', A, LDA, WORK( N+1 ),
00642      $                      IINFO )
00643 *
00644             ELSE IF( ITYPE.EQ.6 ) THEN
00645 *
00646 *              General, eigenvalues specified
00647 *
00648                IF( KCONDS( JTYPE ).EQ.1 ) THEN
00649                   CONDS = ONE
00650                ELSE IF( KCONDS( JTYPE ).EQ.2 ) THEN
00651                   CONDS = RTULPI
00652                ELSE
00653                   CONDS = ZERO
00654                END IF
00655 *
00656                ADUMMA( 1 ) = ' '
00657                CALL SLATME( N, 'S', ISEED, WORK, IMODE, COND, ONE,
00658      $                      ADUMMA, 'T', 'T', 'T', WORK( N+1 ), 4,
00659      $                      CONDS, N, N, ANORM, A, LDA, WORK( 2*N+1 ),
00660      $                      IINFO )
00661 *
00662             ELSE IF( ITYPE.EQ.7 ) THEN
00663 *
00664 *              Diagonal, random eigenvalues
00665 *
00666                CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE,
00667      $                      'T', 'N', WORK( N+1 ), 1, ONE,
00668      $                      WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0,
00669      $                      ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
00670 *
00671             ELSE IF( ITYPE.EQ.8 ) THEN
00672 *
00673 *              Symmetric, random eigenvalues
00674 *
00675                CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE,
00676      $                      'T', 'N', WORK( N+1 ), 1, ONE,
00677      $                      WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N,
00678      $                      ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
00679 *
00680             ELSE IF( ITYPE.EQ.9 ) THEN
00681 *
00682 *              General, random eigenvalues
00683 *
00684                CALL SLATMR( N, N, 'S', ISEED, 'N', WORK, 6, ONE, ONE,
00685      $                      'T', 'N', WORK( N+1 ), 1, ONE,
00686      $                      WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N,
00687      $                      ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
00688                IF( N.GE.4 ) THEN
00689                   CALL SLASET( 'Full', 2, N, ZERO, ZERO, A, LDA )
00690                   CALL SLASET( 'Full', N-3, 1, ZERO, ZERO, A( 3, 1 ),
00691      $                         LDA )
00692                   CALL SLASET( 'Full', N-3, 2, ZERO, ZERO, A( 3, N-1 ),
00693      $                         LDA )
00694                   CALL SLASET( 'Full', 1, N, ZERO, ZERO, A( N, 1 ),
00695      $                         LDA )
00696                END IF
00697 *
00698             ELSE IF( ITYPE.EQ.10 ) THEN
00699 *
00700 *              Triangular, random eigenvalues
00701 *
00702                CALL SLATMR( N, N, 'S', ISEED, 'N', WORK, 6, ONE, ONE,
00703      $                      'T', 'N', WORK( N+1 ), 1, ONE,
00704      $                      WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, 0,
00705      $                      ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
00706 *
00707             ELSE
00708 *
00709                IINFO = 1
00710             END IF
00711 *
00712             IF( IINFO.NE.0 ) THEN
00713                WRITE( NOUNIT, FMT = 9993 )'Generator', IINFO, N, JTYPE,
00714      $            IOLDSD
00715                INFO = ABS( IINFO )
00716                RETURN
00717             END IF
00718 *
00719    90       CONTINUE
00720 *
00721 *           Test for minimal and generous workspace
00722 *
00723             DO 250 IWK = 1, 2
00724                IF( IWK.EQ.1 ) THEN
00725                   NNWORK = 4*N
00726                ELSE
00727                   NNWORK = 5*N + 2*N**2
00728                END IF
00729                NNWORK = MAX( NNWORK, 1 )
00730 *
00731 *              Initialize RESULT
00732 *
00733                DO 100 J = 1, 7
00734                   RESULT( J ) = -ONE
00735   100          CONTINUE
00736 *
00737 *              Compute eigenvalues and eigenvectors, and test them
00738 *
00739                CALL SLACPY( 'F', N, N, A, LDA, H, LDA )
00740                CALL SGEEV( 'V', 'V', N, H, LDA, WR, WI, VL, LDVL, VR,
00741      $                     LDVR, WORK, NNWORK, IINFO )
00742                IF( IINFO.NE.0 ) THEN
00743                   RESULT( 1 ) = ULPINV
00744                   WRITE( NOUNIT, FMT = 9993 )'SGEEV1', IINFO, N, JTYPE,
00745      $               IOLDSD
00746                   INFO = ABS( IINFO )
00747                   GO TO 220
00748                END IF
00749 *
00750 *              Do Test (1)
00751 *
00752                CALL SGET22( 'N', 'N', 'N', N, A, LDA, VR, LDVR, WR, WI,
00753      $                      WORK, RES )
00754                RESULT( 1 ) = RES( 1 )
00755 *
00756 *              Do Test (2)
00757 *
00758                CALL SGET22( 'T', 'N', 'T', N, A, LDA, VL, LDVL, WR, WI,
00759      $                      WORK, RES )
00760                RESULT( 2 ) = RES( 1 )
00761 *
00762 *              Do Test (3)
00763 *
00764                DO 120 J = 1, N
00765                   TNRM = ONE
00766                   IF( WI( J ).EQ.ZERO ) THEN
00767                      TNRM = SNRM2( N, VR( 1, J ), 1 )
00768                   ELSE IF( WI( J ).GT.ZERO ) THEN
00769                      TNRM = SLAPY2( SNRM2( N, VR( 1, J ), 1 ),
00770      $                      SNRM2( N, VR( 1, J+1 ), 1 ) )
00771                   END IF
00772                   RESULT( 3 ) = MAX( RESULT( 3 ),
00773      $                          MIN( ULPINV, ABS( TNRM-ONE ) / ULP ) )
00774                   IF( WI( J ).GT.ZERO ) THEN
00775                      VMX = ZERO
00776                      VRMX = ZERO
00777                      DO 110 JJ = 1, N
00778                         VTST = SLAPY2( VR( JJ, J ), VR( JJ, J+1 ) )
00779                         IF( VTST.GT.VMX )
00780      $                     VMX = VTST
00781                         IF( VR( JJ, J+1 ).EQ.ZERO .AND.
00782      $                      ABS( VR( JJ, J ) ).GT.VRMX )
00783      $                      VRMX = ABS( VR( JJ, J ) )
00784   110                CONTINUE
00785                      IF( VRMX / VMX.LT.ONE-TWO*ULP )
00786      $                  RESULT( 3 ) = ULPINV
00787                   END IF
00788   120          CONTINUE
00789 *
00790 *              Do Test (4)
00791 *
00792                DO 140 J = 1, N
00793                   TNRM = ONE
00794                   IF( WI( J ).EQ.ZERO ) THEN
00795                      TNRM = SNRM2( N, VL( 1, J ), 1 )
00796                   ELSE IF( WI( J ).GT.ZERO ) THEN
00797                      TNRM = SLAPY2( SNRM2( N, VL( 1, J ), 1 ),
00798      $                      SNRM2( N, VL( 1, J+1 ), 1 ) )
00799                   END IF
00800                   RESULT( 4 ) = MAX( RESULT( 4 ),
00801      $                          MIN( ULPINV, ABS( TNRM-ONE ) / ULP ) )
00802                   IF( WI( J ).GT.ZERO ) THEN
00803                      VMX = ZERO
00804                      VRMX = ZERO
00805                      DO 130 JJ = 1, N
00806                         VTST = SLAPY2( VL( JJ, J ), VL( JJ, J+1 ) )
00807                         IF( VTST.GT.VMX )
00808      $                     VMX = VTST
00809                         IF( VL( JJ, J+1 ).EQ.ZERO .AND.
00810      $                      ABS( VL( JJ, J ) ).GT.VRMX )
00811      $                      VRMX = ABS( VL( JJ, J ) )
00812   130                CONTINUE
00813                      IF( VRMX / VMX.LT.ONE-TWO*ULP )
00814      $                  RESULT( 4 ) = ULPINV
00815                   END IF
00816   140          CONTINUE
00817 *
00818 *              Compute eigenvalues only, and test them
00819 *
00820                CALL SLACPY( 'F', N, N, A, LDA, H, LDA )
00821                CALL SGEEV( 'N', 'N', N, H, LDA, WR1, WI1, DUM, 1, DUM,
00822      $                     1, WORK, NNWORK, IINFO )
00823                IF( IINFO.NE.0 ) THEN
00824                   RESULT( 1 ) = ULPINV
00825                   WRITE( NOUNIT, FMT = 9993 )'SGEEV2', IINFO, N, JTYPE,
00826      $               IOLDSD
00827                   INFO = ABS( IINFO )
00828                   GO TO 220
00829                END IF
00830 *
00831 *              Do Test (5)
00832 *
00833                DO 150 J = 1, N
00834                   IF( WR( J ).NE.WR1( J ) .OR. WI( J ).NE.WI1( J ) )
00835      $               RESULT( 5 ) = ULPINV
00836   150          CONTINUE
00837 *
00838 *              Compute eigenvalues and right eigenvectors, and test them
00839 *
00840                CALL SLACPY( 'F', N, N, A, LDA, H, LDA )
00841                CALL SGEEV( 'N', 'V', N, H, LDA, WR1, WI1, DUM, 1, LRE,
00842      $                     LDLRE, WORK, NNWORK, IINFO )
00843                IF( IINFO.NE.0 ) THEN
00844                   RESULT( 1 ) = ULPINV
00845                   WRITE( NOUNIT, FMT = 9993 )'SGEEV3', IINFO, N, JTYPE,
00846      $               IOLDSD
00847                   INFO = ABS( IINFO )
00848                   GO TO 220
00849                END IF
00850 *
00851 *              Do Test (5) again
00852 *
00853                DO 160 J = 1, N
00854                   IF( WR( J ).NE.WR1( J ) .OR. WI( J ).NE.WI1( J ) )
00855      $               RESULT( 5 ) = ULPINV
00856   160          CONTINUE
00857 *
00858 *              Do Test (6)
00859 *
00860                DO 180 J = 1, N
00861                   DO 170 JJ = 1, N
00862                      IF( VR( J, JJ ).NE.LRE( J, JJ ) )
00863      $                  RESULT( 6 ) = ULPINV
00864   170             CONTINUE
00865   180          CONTINUE
00866 *
00867 *              Compute eigenvalues and left eigenvectors, and test them
00868 *
00869                CALL SLACPY( 'F', N, N, A, LDA, H, LDA )
00870                CALL SGEEV( 'V', 'N', N, H, LDA, WR1, WI1, LRE, LDLRE,
00871      $                     DUM, 1, WORK, NNWORK, IINFO )
00872                IF( IINFO.NE.0 ) THEN
00873                   RESULT( 1 ) = ULPINV
00874                   WRITE( NOUNIT, FMT = 9993 )'SGEEV4', IINFO, N, JTYPE,
00875      $               IOLDSD
00876                   INFO = ABS( IINFO )
00877                   GO TO 220
00878                END IF
00879 *
00880 *              Do Test (5) again
00881 *
00882                DO 190 J = 1, N
00883                   IF( WR( J ).NE.WR1( J ) .OR. WI( J ).NE.WI1( J ) )
00884      $               RESULT( 5 ) = ULPINV
00885   190          CONTINUE
00886 *
00887 *              Do Test (7)
00888 *
00889                DO 210 J = 1, N
00890                   DO 200 JJ = 1, N
00891                      IF( VL( J, JJ ).NE.LRE( J, JJ ) )
00892      $                  RESULT( 7 ) = ULPINV
00893   200             CONTINUE
00894   210          CONTINUE
00895 *
00896 *              End of Loop -- Check for RESULT(j) > THRESH
00897 *
00898   220          CONTINUE
00899 *
00900                NTEST = 0
00901                NFAIL = 0
00902                DO 230 J = 1, 7
00903                   IF( RESULT( J ).GE.ZERO )
00904      $               NTEST = NTEST + 1
00905                   IF( RESULT( J ).GE.THRESH )
00906      $               NFAIL = NFAIL + 1
00907   230          CONTINUE
00908 *
00909                IF( NFAIL.GT.0 )
00910      $            NTESTF = NTESTF + 1
00911                IF( NTESTF.EQ.1 ) THEN
00912                   WRITE( NOUNIT, FMT = 9999 )PATH
00913                   WRITE( NOUNIT, FMT = 9998 )
00914                   WRITE( NOUNIT, FMT = 9997 )
00915                   WRITE( NOUNIT, FMT = 9996 )
00916                   WRITE( NOUNIT, FMT = 9995 )THRESH
00917                   NTESTF = 2
00918                END IF
00919 *
00920                DO 240 J = 1, 7
00921                   IF( RESULT( J ).GE.THRESH ) THEN
00922                      WRITE( NOUNIT, FMT = 9994 )N, IWK, IOLDSD, JTYPE,
00923      $                  J, RESULT( J )
00924                   END IF
00925   240          CONTINUE
00926 *
00927                NERRS = NERRS + NFAIL
00928                NTESTT = NTESTT + NTEST
00929 *
00930   250       CONTINUE
00931   260    CONTINUE
00932   270 CONTINUE
00933 *
00934 *     Summary
00935 *
00936       CALL SLASUM( PATH, NOUNIT, NERRS, NTESTT )
00937 *
00938  9999 FORMAT( / 1X, A3, ' -- Real Eigenvalue-Eigenvector Decomposition',
00939      $      ' Driver', / ' Matrix types (see SDRVEV for details): ' )
00940 *
00941  9998 FORMAT( / ' Special Matrices:', / '  1=Zero matrix.             ',
00942      $      '           ', '  5=Diagonal: geometr. spaced entries.',
00943      $      / '  2=Identity matrix.                    ', '  6=Diagona',
00944      $      'l: clustered entries.', / '  3=Transposed Jordan block.  ',
00945      $      '          ', '  7=Diagonal: large, evenly spaced.', / '  ',
00946      $      '4=Diagonal: evenly spaced entries.    ', '  8=Diagonal: s',
00947      $      'mall, evenly spaced.' )
00948  9997 FORMAT( ' Dense, Non-Symmetric Matrices:', / '  9=Well-cond., ev',
00949      $      'enly spaced eigenvals.', ' 14=Ill-cond., geomet. spaced e',
00950      $      'igenals.', / ' 10=Well-cond., geom. spaced eigenvals. ',
00951      $      ' 15=Ill-conditioned, clustered e.vals.', / ' 11=Well-cond',
00952      $      'itioned, clustered e.vals. ', ' 16=Ill-cond., random comp',
00953      $      'lex ', / ' 12=Well-cond., random complex ', 6X, '   ',
00954      $      ' 17=Ill-cond., large rand. complx ', / ' 13=Ill-condi',
00955      $      'tioned, evenly spaced.     ', ' 18=Ill-cond., small rand.',
00956      $      ' complx ' )
00957  9996 FORMAT( ' 19=Matrix with random O(1) entries.    ', ' 21=Matrix ',
00958      $      'with small random entries.', / ' 20=Matrix with large ran',
00959      $      'dom entries.   ', / )
00960  9995 FORMAT( ' Tests performed with test threshold =', F8.2,
00961      $      / / ' 1 = | A VR - VR W | / ( n |A| ulp ) ',
00962      $      / ' 2 = | transpose(A) VL - VL W | / ( n |A| ulp ) ',
00963      $      / ' 3 = | |VR(i)| - 1 | / ulp ',
00964      $      / ' 4 = | |VL(i)| - 1 | / ulp ',
00965      $      / ' 5 = 0 if W same no matter if VR or VL computed,',
00966      $      ' 1/ulp otherwise', /
00967      $      ' 6 = 0 if VR same no matter if VL computed,',
00968      $      '  1/ulp otherwise', /
00969      $      ' 7 = 0 if VL same no matter if VR computed,',
00970      $      '  1/ulp otherwise', / )
00971  9994 FORMAT( ' N=', I5, ', IWK=', I2, ', seed=', 4( I4, ',' ),
00972      $      ' type ', I2, ', test(', I2, ')=', G10.3 )
00973  9993 FORMAT( ' SDRVEV: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
00974      $      I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
00975 *
00976       RETURN
00977 *
00978 *     End of SDRVEV
00979 *
00980       END
 All Files Functions