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