![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b SGET23 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 SGET23( COMP, BALANC, JTYPE, THRESH, ISEED, NOUNIT, N, 00012 * A, LDA, H, WR, WI, WR1, WI1, VL, LDVL, VR, 00013 * LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN, 00014 * RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT, 00015 * WORK, LWORK, IWORK, INFO ) 00016 * 00017 * .. Scalar Arguments .. 00018 * LOGICAL COMP 00019 * CHARACTER BALANC 00020 * INTEGER INFO, JTYPE, LDA, LDLRE, LDVL, LDVR, LWORK, N, 00021 * $ NOUNIT 00022 * REAL THRESH 00023 * .. 00024 * .. Array Arguments .. 00025 * INTEGER ISEED( 4 ), IWORK( * ) 00026 * REAL A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ), 00027 * $ RCDEIN( * ), RCDVIN( * ), RCNDE1( * ), 00028 * $ RCNDV1( * ), RCONDE( * ), RCONDV( * ), 00029 * $ RESULT( 11 ), SCALE( * ), SCALE1( * ), 00030 * $ VL( LDVL, * ), VR( LDVR, * ), WI( * ), 00031 * $ WI1( * ), WORK( * ), WR( * ), WR1( * ) 00032 * .. 00033 * 00034 * 00035 *> \par Purpose: 00036 * ============= 00037 *> 00038 *> \verbatim 00039 *> 00040 *> SGET23 checks the nonsymmetric eigenvalue problem driver SGEEVX. 00041 *> If COMP = .FALSE., the first 8 of the following tests will be 00042 *> performed on the input matrix A, and also test 9 if LWORK is 00043 *> sufficiently large. 00044 *> if COMP is .TRUE. all 11 tests will be performed. 00045 *> 00046 *> (1) | A * VR - VR * W | / ( n |A| ulp ) 00047 *> 00048 *> Here VR is the matrix of unit right eigenvectors. 00049 *> W is a block diagonal matrix, with a 1x1 block for each 00050 *> real eigenvalue and a 2x2 block for each complex conjugate 00051 *> pair. If eigenvalues j and j+1 are a complex conjugate pair, 00052 *> so WR(j) = WR(j+1) = wr and WI(j) = - WI(j+1) = wi, then the 00053 *> 2 x 2 block corresponding to the pair will be: 00054 *> 00055 *> ( wr wi ) 00056 *> ( -wi wr ) 00057 *> 00058 *> Such a block multiplying an n x 2 matrix ( ur ui ) on the 00059 *> right will be the same as multiplying ur + i*ui by wr + i*wi. 00060 *> 00061 *> (2) | A**H * VL - VL * W**H | / ( n |A| ulp ) 00062 *> 00063 *> Here VL is the matrix of unit left eigenvectors, A**H is the 00064 *> conjugate transpose of A, and W is as above. 00065 *> 00066 *> (3) | |VR(i)| - 1 | / ulp and largest component real 00067 *> 00068 *> VR(i) denotes the i-th column of VR. 00069 *> 00070 *> (4) | |VL(i)| - 1 | / ulp and largest component real 00071 *> 00072 *> VL(i) denotes the i-th column of VL. 00073 *> 00074 *> (5) 0 if W(full) = W(partial), 1/ulp otherwise 00075 *> 00076 *> W(full) denotes the eigenvalues computed when VR, VL, RCONDV 00077 *> and RCONDE are also computed, and W(partial) denotes the 00078 *> eigenvalues computed when only some of VR, VL, RCONDV, and 00079 *> RCONDE are computed. 00080 *> 00081 *> (6) 0 if VR(full) = VR(partial), 1/ulp otherwise 00082 *> 00083 *> VR(full) denotes the right eigenvectors computed when VL, RCONDV 00084 *> and RCONDE are computed, and VR(partial) denotes the result 00085 *> when only some of VL and RCONDV are computed. 00086 *> 00087 *> (7) 0 if VL(full) = VL(partial), 1/ulp otherwise 00088 *> 00089 *> VL(full) denotes the left eigenvectors computed when VR, RCONDV 00090 *> and RCONDE are computed, and VL(partial) denotes the result 00091 *> when only some of VR and RCONDV are computed. 00092 *> 00093 *> (8) 0 if SCALE, ILO, IHI, ABNRM (full) = 00094 *> SCALE, ILO, IHI, ABNRM (partial) 00095 *> 1/ulp otherwise 00096 *> 00097 *> SCALE, ILO, IHI and ABNRM describe how the matrix is balanced. 00098 *> (full) is when VR, VL, RCONDE and RCONDV are also computed, and 00099 *> (partial) is when some are not computed. 00100 *> 00101 *> (9) 0 if RCONDV(full) = RCONDV(partial), 1/ulp otherwise 00102 *> 00103 *> RCONDV(full) denotes the reciprocal condition numbers of the 00104 *> right eigenvectors computed when VR, VL and RCONDE are also 00105 *> computed. RCONDV(partial) denotes the reciprocal condition 00106 *> numbers when only some of VR, VL and RCONDE are computed. 00107 *> 00108 *> (10) |RCONDV - RCDVIN| / cond(RCONDV) 00109 *> 00110 *> RCONDV is the reciprocal right eigenvector condition number 00111 *> computed by SGEEVX and RCDVIN (the precomputed true value) 00112 *> is supplied as input. cond(RCONDV) is the condition number of 00113 *> RCONDV, and takes errors in computing RCONDV into account, so 00114 *> that the resulting quantity should be O(ULP). cond(RCONDV) is 00115 *> essentially given by norm(A)/RCONDE. 00116 *> 00117 *> (11) |RCONDE - RCDEIN| / cond(RCONDE) 00118 *> 00119 *> RCONDE is the reciprocal eigenvalue condition number 00120 *> computed by SGEEVX and RCDEIN (the precomputed true value) 00121 *> is supplied as input. cond(RCONDE) is the condition number 00122 *> of RCONDE, and takes errors in computing RCONDE into account, 00123 *> so that the resulting quantity should be O(ULP). cond(RCONDE) 00124 *> is essentially given by norm(A)/RCONDV. 00125 *> \endverbatim 00126 * 00127 * Arguments: 00128 * ========== 00129 * 00130 *> \param[in] COMP 00131 *> \verbatim 00132 *> COMP is LOGICAL 00133 *> COMP describes which input tests to perform: 00134 *> = .FALSE. if the computed condition numbers are not to 00135 *> be tested against RCDVIN and RCDEIN 00136 *> = .TRUE. if they are to be compared 00137 *> \endverbatim 00138 *> 00139 *> \param[in] BALANC 00140 *> \verbatim 00141 *> BALANC is CHARACTER 00142 *> Describes the balancing option to be tested. 00143 *> = 'N' for no permuting or diagonal scaling 00144 *> = 'P' for permuting but no diagonal scaling 00145 *> = 'S' for no permuting but diagonal scaling 00146 *> = 'B' for permuting and diagonal scaling 00147 *> \endverbatim 00148 *> 00149 *> \param[in] JTYPE 00150 *> \verbatim 00151 *> JTYPE is INTEGER 00152 *> Type of input matrix. Used to label output if error occurs. 00153 *> \endverbatim 00154 *> 00155 *> \param[in] THRESH 00156 *> \verbatim 00157 *> THRESH is REAL 00158 *> A test will count as "failed" if the "error", computed as 00159 *> described above, exceeds THRESH. Note that the error 00160 *> is scaled to be O(1), so THRESH should be a reasonably 00161 *> small multiple of 1, e.g., 10 or 100. In particular, 00162 *> it should not depend on the precision (single vs. double) 00163 *> or the size of the matrix. It must be at least zero. 00164 *> \endverbatim 00165 *> 00166 *> \param[in] ISEED 00167 *> \verbatim 00168 *> ISEED is INTEGER array, dimension (4) 00169 *> If COMP = .FALSE., the random number generator seed 00170 *> used to produce matrix. 00171 *> If COMP = .TRUE., ISEED(1) = the number of the example. 00172 *> Used to label output if error occurs. 00173 *> \endverbatim 00174 *> 00175 *> \param[in] NOUNIT 00176 *> \verbatim 00177 *> NOUNIT is INTEGER 00178 *> The FORTRAN unit number for printing out error messages 00179 *> (e.g., if a routine returns INFO not equal to 0.) 00180 *> \endverbatim 00181 *> 00182 *> \param[in] N 00183 *> \verbatim 00184 *> N is INTEGER 00185 *> The dimension of A. N must be at least 0. 00186 *> \endverbatim 00187 *> 00188 *> \param[in,out] A 00189 *> \verbatim 00190 *> A is REAL array, dimension (LDA,N) 00191 *> Used to hold the matrix whose eigenvalues are to be 00192 *> computed. 00193 *> \endverbatim 00194 *> 00195 *> \param[in] LDA 00196 *> \verbatim 00197 *> LDA is INTEGER 00198 *> The leading dimension of A, and H. LDA must be at 00199 *> least 1 and at least N. 00200 *> \endverbatim 00201 *> 00202 *> \param[out] H 00203 *> \verbatim 00204 *> H is REAL array, dimension (LDA,N) 00205 *> Another copy of the test matrix A, modified by SGEEVX. 00206 *> \endverbatim 00207 *> 00208 *> \param[out] WR 00209 *> \verbatim 00210 *> WR is REAL array, dimension (N) 00211 *> \endverbatim 00212 *> 00213 *> \param[out] WI 00214 *> \verbatim 00215 *> WI is REAL array, dimension (N) 00216 *> 00217 *> The real and imaginary parts of the eigenvalues of A. 00218 *> On exit, WR + WI*i are the eigenvalues of the matrix in A. 00219 *> \endverbatim 00220 *> 00221 *> \param[out] WR1 00222 *> \verbatim 00223 *> WR1 is REAL array, dimension (N) 00224 *> \endverbatim 00225 *> 00226 *> \param[out] WI1 00227 *> \verbatim 00228 *> WI1 is REAL array, dimension (N) 00229 *> 00230 *> Like WR, WI, these arrays contain the eigenvalues of A, 00231 *> but those computed when SGEEVX only computes a partial 00232 *> eigendecomposition, i.e. not the eigenvalues and left 00233 *> and right eigenvectors. 00234 *> \endverbatim 00235 *> 00236 *> \param[out] VL 00237 *> \verbatim 00238 *> VL is REAL array, dimension (LDVL,N) 00239 *> VL holds the computed left eigenvectors. 00240 *> \endverbatim 00241 *> 00242 *> \param[in] LDVL 00243 *> \verbatim 00244 *> LDVL is INTEGER 00245 *> Leading dimension of VL. Must be at least max(1,N). 00246 *> \endverbatim 00247 *> 00248 *> \param[out] VR 00249 *> \verbatim 00250 *> VR is REAL array, dimension (LDVR,N) 00251 *> VR holds the computed right eigenvectors. 00252 *> \endverbatim 00253 *> 00254 *> \param[in] LDVR 00255 *> \verbatim 00256 *> LDVR is INTEGER 00257 *> Leading dimension of VR. Must be at least max(1,N). 00258 *> \endverbatim 00259 *> 00260 *> \param[out] LRE 00261 *> \verbatim 00262 *> LRE is REAL array, dimension (LDLRE,N) 00263 *> LRE holds the computed right or left eigenvectors. 00264 *> \endverbatim 00265 *> 00266 *> \param[in] LDLRE 00267 *> \verbatim 00268 *> LDLRE is INTEGER 00269 *> Leading dimension of LRE. Must be at least max(1,N). 00270 *> \endverbatim 00271 *> 00272 *> \param[out] RCONDV 00273 *> \verbatim 00274 *> RCONDV is REAL array, dimension (N) 00275 *> RCONDV holds the computed reciprocal condition numbers 00276 *> for eigenvectors. 00277 *> \endverbatim 00278 *> 00279 *> \param[out] RCNDV1 00280 *> \verbatim 00281 *> RCNDV1 is REAL array, dimension (N) 00282 *> RCNDV1 holds more computed reciprocal condition numbers 00283 *> for eigenvectors. 00284 *> \endverbatim 00285 *> 00286 *> \param[in] RCDVIN 00287 *> \verbatim 00288 *> RCDVIN is REAL array, dimension (N) 00289 *> When COMP = .TRUE. RCDVIN holds the precomputed reciprocal 00290 *> condition numbers for eigenvectors to be compared with 00291 *> RCONDV. 00292 *> \endverbatim 00293 *> 00294 *> \param[out] RCONDE 00295 *> \verbatim 00296 *> RCONDE is REAL array, dimension (N) 00297 *> RCONDE holds the computed reciprocal condition numbers 00298 *> for eigenvalues. 00299 *> \endverbatim 00300 *> 00301 *> \param[out] RCNDE1 00302 *> \verbatim 00303 *> RCNDE1 is REAL array, dimension (N) 00304 *> RCNDE1 holds more computed reciprocal condition numbers 00305 *> for eigenvalues. 00306 *> \endverbatim 00307 *> 00308 *> \param[in] RCDEIN 00309 *> \verbatim 00310 *> RCDEIN is REAL array, dimension (N) 00311 *> When COMP = .TRUE. RCDEIN holds the precomputed reciprocal 00312 *> condition numbers for eigenvalues to be compared with 00313 *> RCONDE. 00314 *> \endverbatim 00315 *> 00316 *> \param[out] SCALE 00317 *> \verbatim 00318 *> SCALE is REAL array, dimension (N) 00319 *> Holds information describing balancing of matrix. 00320 *> \endverbatim 00321 *> 00322 *> \param[out] SCALE1 00323 *> \verbatim 00324 *> SCALE1 is REAL array, dimension (N) 00325 *> Holds information describing balancing of matrix. 00326 *> \endverbatim 00327 *> 00328 *> \param[out] RESULT 00329 *> \verbatim 00330 *> RESULT is REAL array, dimension (11) 00331 *> The values computed by the 11 tests described above. 00332 *> The values are currently limited to 1/ulp, to avoid 00333 *> overflow. 00334 *> \endverbatim 00335 *> 00336 *> \param[out] WORK 00337 *> \verbatim 00338 *> WORK is REAL array, dimension (LWORK) 00339 *> \endverbatim 00340 *> 00341 *> \param[in] LWORK 00342 *> \verbatim 00343 *> LWORK is INTEGER 00344 *> The number of entries in WORK. This must be at least 00345 *> 3*N, and 6*N+N**2 if tests 9, 10 or 11 are to be performed. 00346 *> \endverbatim 00347 *> 00348 *> \param[out] IWORK 00349 *> \verbatim 00350 *> IWORK is INTEGER array, dimension (2*N) 00351 *> \endverbatim 00352 *> 00353 *> \param[out] INFO 00354 *> \verbatim 00355 *> INFO is INTEGER 00356 *> If 0, successful exit. 00357 *> If <0, input parameter -INFO had an incorrect value. 00358 *> If >0, SGEEVX returned an error code, the absolute 00359 *> value of which is returned. 00360 *> \endverbatim 00361 * 00362 * Authors: 00363 * ======== 00364 * 00365 *> \author Univ. of Tennessee 00366 *> \author Univ. of California Berkeley 00367 *> \author Univ. of Colorado Denver 00368 *> \author NAG Ltd. 00369 * 00370 *> \date November 2011 00371 * 00372 *> \ingroup single_eig 00373 * 00374 * ===================================================================== 00375 SUBROUTINE SGET23( COMP, BALANC, JTYPE, THRESH, ISEED, NOUNIT, N, 00376 $ A, LDA, H, WR, WI, WR1, WI1, VL, LDVL, VR, 00377 $ LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN, 00378 $ RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT, 00379 $ WORK, LWORK, IWORK, INFO ) 00380 * 00381 * -- LAPACK test routine (version 3.4.0) -- 00382 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00383 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00384 * November 2011 00385 * 00386 * .. Scalar Arguments .. 00387 LOGICAL COMP 00388 CHARACTER BALANC 00389 INTEGER INFO, JTYPE, LDA, LDLRE, LDVL, LDVR, LWORK, N, 00390 $ NOUNIT 00391 REAL THRESH 00392 * .. 00393 * .. Array Arguments .. 00394 INTEGER ISEED( 4 ), IWORK( * ) 00395 REAL A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ), 00396 $ RCDEIN( * ), RCDVIN( * ), RCNDE1( * ), 00397 $ RCNDV1( * ), RCONDE( * ), RCONDV( * ), 00398 $ RESULT( 11 ), SCALE( * ), SCALE1( * ), 00399 $ VL( LDVL, * ), VR( LDVR, * ), WI( * ), 00400 $ WI1( * ), WORK( * ), WR( * ), WR1( * ) 00401 * .. 00402 * 00403 * ===================================================================== 00404 * 00405 * 00406 * .. Parameters .. 00407 REAL ZERO, ONE, TWO 00408 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0 ) 00409 REAL EPSIN 00410 PARAMETER ( EPSIN = 5.9605E-8 ) 00411 * .. 00412 * .. Local Scalars .. 00413 LOGICAL BALOK, NOBAL 00414 CHARACTER SENSE 00415 INTEGER I, IHI, IHI1, IINFO, ILO, ILO1, ISENS, ISENSM, 00416 $ J, JJ, KMIN 00417 REAL ABNRM, ABNRM1, EPS, SMLNUM, TNRM, TOL, TOLIN, 00418 $ ULP, ULPINV, V, VIMIN, VMAX, VMX, VRMIN, VRMX, 00419 $ VTST 00420 * .. 00421 * .. Local Arrays .. 00422 CHARACTER SENS( 2 ) 00423 REAL DUM( 1 ), RES( 2 ) 00424 * .. 00425 * .. External Functions .. 00426 LOGICAL LSAME 00427 REAL SLAMCH, SLAPY2, SNRM2 00428 EXTERNAL LSAME, SLAMCH, SLAPY2, SNRM2 00429 * .. 00430 * .. External Subroutines .. 00431 EXTERNAL SGEEVX, SGET22, SLACPY, XERBLA 00432 * .. 00433 * .. Intrinsic Functions .. 00434 INTRINSIC ABS, MAX, MIN, REAL 00435 * .. 00436 * .. Data statements .. 00437 DATA SENS / 'N', 'V' / 00438 * .. 00439 * .. Executable Statements .. 00440 * 00441 * Check for errors 00442 * 00443 NOBAL = LSAME( BALANC, 'N' ) 00444 BALOK = NOBAL .OR. LSAME( BALANC, 'P' ) .OR. 00445 $ LSAME( BALANC, 'S' ) .OR. LSAME( BALANC, 'B' ) 00446 INFO = 0 00447 IF( .NOT.BALOK ) THEN 00448 INFO = -2 00449 ELSE IF( THRESH.LT.ZERO ) THEN 00450 INFO = -4 00451 ELSE IF( NOUNIT.LE.0 ) THEN 00452 INFO = -6 00453 ELSE IF( N.LT.0 ) THEN 00454 INFO = -7 00455 ELSE IF( LDA.LT.1 .OR. LDA.LT.N ) THEN 00456 INFO = -9 00457 ELSE IF( LDVL.LT.1 .OR. LDVL.LT.N ) THEN 00458 INFO = -16 00459 ELSE IF( LDVR.LT.1 .OR. LDVR.LT.N ) THEN 00460 INFO = -18 00461 ELSE IF( LDLRE.LT.1 .OR. LDLRE.LT.N ) THEN 00462 INFO = -20 00463 ELSE IF( LWORK.LT.3*N .OR. ( COMP .AND. LWORK.LT.6*N+N*N ) ) THEN 00464 INFO = -31 00465 END IF 00466 * 00467 IF( INFO.NE.0 ) THEN 00468 CALL XERBLA( 'SGET23', -INFO ) 00469 RETURN 00470 END IF 00471 * 00472 * Quick return if nothing to do 00473 * 00474 DO 10 I = 1, 11 00475 RESULT( I ) = -ONE 00476 10 CONTINUE 00477 * 00478 IF( N.EQ.0 ) 00479 $ RETURN 00480 * 00481 * More Important constants 00482 * 00483 ULP = SLAMCH( 'Precision' ) 00484 SMLNUM = SLAMCH( 'S' ) 00485 ULPINV = ONE / ULP 00486 * 00487 * Compute eigenvalues and eigenvectors, and test them 00488 * 00489 IF( LWORK.GE.6*N+N*N ) THEN 00490 SENSE = 'B' 00491 ISENSM = 2 00492 ELSE 00493 SENSE = 'E' 00494 ISENSM = 1 00495 END IF 00496 CALL SLACPY( 'F', N, N, A, LDA, H, LDA ) 00497 CALL SGEEVX( BALANC, 'V', 'V', SENSE, N, H, LDA, WR, WI, VL, LDVL, 00498 $ VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV, 00499 $ WORK, LWORK, IWORK, IINFO ) 00500 IF( IINFO.NE.0 ) THEN 00501 RESULT( 1 ) = ULPINV 00502 IF( JTYPE.NE.22 ) THEN 00503 WRITE( NOUNIT, FMT = 9998 )'SGEEVX1', IINFO, N, JTYPE, 00504 $ BALANC, ISEED 00505 ELSE 00506 WRITE( NOUNIT, FMT = 9999 )'SGEEVX1', IINFO, N, ISEED( 1 ) 00507 END IF 00508 INFO = ABS( IINFO ) 00509 RETURN 00510 END IF 00511 * 00512 * Do Test (1) 00513 * 00514 CALL SGET22( 'N', 'N', 'N', N, A, LDA, VR, LDVR, WR, WI, WORK, 00515 $ RES ) 00516 RESULT( 1 ) = RES( 1 ) 00517 * 00518 * Do Test (2) 00519 * 00520 CALL SGET22( 'T', 'N', 'T', N, A, LDA, VL, LDVL, WR, WI, WORK, 00521 $ RES ) 00522 RESULT( 2 ) = RES( 1 ) 00523 * 00524 * Do Test (3) 00525 * 00526 DO 30 J = 1, N 00527 TNRM = ONE 00528 IF( WI( J ).EQ.ZERO ) THEN 00529 TNRM = SNRM2( N, VR( 1, J ), 1 ) 00530 ELSE IF( WI( J ).GT.ZERO ) THEN 00531 TNRM = SLAPY2( SNRM2( N, VR( 1, J ), 1 ), 00532 $ SNRM2( N, VR( 1, J+1 ), 1 ) ) 00533 END IF 00534 RESULT( 3 ) = MAX( RESULT( 3 ), 00535 $ MIN( ULPINV, ABS( TNRM-ONE ) / ULP ) ) 00536 IF( WI( J ).GT.ZERO ) THEN 00537 VMX = ZERO 00538 VRMX = ZERO 00539 DO 20 JJ = 1, N 00540 VTST = SLAPY2( VR( JJ, J ), VR( JJ, J+1 ) ) 00541 IF( VTST.GT.VMX ) 00542 $ VMX = VTST 00543 IF( VR( JJ, J+1 ).EQ.ZERO .AND. ABS( VR( JJ, J ) ).GT. 00544 $ VRMX )VRMX = ABS( VR( JJ, J ) ) 00545 20 CONTINUE 00546 IF( VRMX / VMX.LT.ONE-TWO*ULP ) 00547 $ RESULT( 3 ) = ULPINV 00548 END IF 00549 30 CONTINUE 00550 * 00551 * Do Test (4) 00552 * 00553 DO 50 J = 1, N 00554 TNRM = ONE 00555 IF( WI( J ).EQ.ZERO ) THEN 00556 TNRM = SNRM2( N, VL( 1, J ), 1 ) 00557 ELSE IF( WI( J ).GT.ZERO ) THEN 00558 TNRM = SLAPY2( SNRM2( N, VL( 1, J ), 1 ), 00559 $ SNRM2( N, VL( 1, J+1 ), 1 ) ) 00560 END IF 00561 RESULT( 4 ) = MAX( RESULT( 4 ), 00562 $ MIN( ULPINV, ABS( TNRM-ONE ) / ULP ) ) 00563 IF( WI( J ).GT.ZERO ) THEN 00564 VMX = ZERO 00565 VRMX = ZERO 00566 DO 40 JJ = 1, N 00567 VTST = SLAPY2( VL( JJ, J ), VL( JJ, J+1 ) ) 00568 IF( VTST.GT.VMX ) 00569 $ VMX = VTST 00570 IF( VL( JJ, J+1 ).EQ.ZERO .AND. ABS( VL( JJ, J ) ).GT. 00571 $ VRMX )VRMX = ABS( VL( JJ, J ) ) 00572 40 CONTINUE 00573 IF( VRMX / VMX.LT.ONE-TWO*ULP ) 00574 $ RESULT( 4 ) = ULPINV 00575 END IF 00576 50 CONTINUE 00577 * 00578 * Test for all options of computing condition numbers 00579 * 00580 DO 200 ISENS = 1, ISENSM 00581 * 00582 SENSE = SENS( ISENS ) 00583 * 00584 * Compute eigenvalues only, and test them 00585 * 00586 CALL SLACPY( 'F', N, N, A, LDA, H, LDA ) 00587 CALL SGEEVX( BALANC, 'N', 'N', SENSE, N, H, LDA, WR1, WI1, DUM, 00588 $ 1, DUM, 1, ILO1, IHI1, SCALE1, ABNRM1, RCNDE1, 00589 $ RCNDV1, WORK, LWORK, IWORK, IINFO ) 00590 IF( IINFO.NE.0 ) THEN 00591 RESULT( 1 ) = ULPINV 00592 IF( JTYPE.NE.22 ) THEN 00593 WRITE( NOUNIT, FMT = 9998 )'SGEEVX2', IINFO, N, JTYPE, 00594 $ BALANC, ISEED 00595 ELSE 00596 WRITE( NOUNIT, FMT = 9999 )'SGEEVX2', IINFO, N, 00597 $ ISEED( 1 ) 00598 END IF 00599 INFO = ABS( IINFO ) 00600 GO TO 190 00601 END IF 00602 * 00603 * Do Test (5) 00604 * 00605 DO 60 J = 1, N 00606 IF( WR( J ).NE.WR1( J ) .OR. WI( J ).NE.WI1( J ) ) 00607 $ RESULT( 5 ) = ULPINV 00608 60 CONTINUE 00609 * 00610 * Do Test (8) 00611 * 00612 IF( .NOT.NOBAL ) THEN 00613 DO 70 J = 1, N 00614 IF( SCALE( J ).NE.SCALE1( J ) ) 00615 $ RESULT( 8 ) = ULPINV 00616 70 CONTINUE 00617 IF( ILO.NE.ILO1 ) 00618 $ RESULT( 8 ) = ULPINV 00619 IF( IHI.NE.IHI1 ) 00620 $ RESULT( 8 ) = ULPINV 00621 IF( ABNRM.NE.ABNRM1 ) 00622 $ RESULT( 8 ) = ULPINV 00623 END IF 00624 * 00625 * Do Test (9) 00626 * 00627 IF( ISENS.EQ.2 .AND. N.GT.1 ) THEN 00628 DO 80 J = 1, N 00629 IF( RCONDV( J ).NE.RCNDV1( J ) ) 00630 $ RESULT( 9 ) = ULPINV 00631 80 CONTINUE 00632 END IF 00633 * 00634 * Compute eigenvalues and right eigenvectors, and test them 00635 * 00636 CALL SLACPY( 'F', N, N, A, LDA, H, LDA ) 00637 CALL SGEEVX( BALANC, 'N', 'V', SENSE, N, H, LDA, WR1, WI1, DUM, 00638 $ 1, LRE, LDLRE, ILO1, IHI1, SCALE1, ABNRM1, RCNDE1, 00639 $ RCNDV1, WORK, LWORK, IWORK, IINFO ) 00640 IF( IINFO.NE.0 ) THEN 00641 RESULT( 1 ) = ULPINV 00642 IF( JTYPE.NE.22 ) THEN 00643 WRITE( NOUNIT, FMT = 9998 )'SGEEVX3', IINFO, N, JTYPE, 00644 $ BALANC, ISEED 00645 ELSE 00646 WRITE( NOUNIT, FMT = 9999 )'SGEEVX3', IINFO, N, 00647 $ ISEED( 1 ) 00648 END IF 00649 INFO = ABS( IINFO ) 00650 GO TO 190 00651 END IF 00652 * 00653 * Do Test (5) again 00654 * 00655 DO 90 J = 1, N 00656 IF( WR( J ).NE.WR1( J ) .OR. WI( J ).NE.WI1( J ) ) 00657 $ RESULT( 5 ) = ULPINV 00658 90 CONTINUE 00659 * 00660 * Do Test (6) 00661 * 00662 DO 110 J = 1, N 00663 DO 100 JJ = 1, N 00664 IF( VR( J, JJ ).NE.LRE( J, JJ ) ) 00665 $ RESULT( 6 ) = ULPINV 00666 100 CONTINUE 00667 110 CONTINUE 00668 * 00669 * Do Test (8) again 00670 * 00671 IF( .NOT.NOBAL ) THEN 00672 DO 120 J = 1, N 00673 IF( SCALE( J ).NE.SCALE1( J ) ) 00674 $ RESULT( 8 ) = ULPINV 00675 120 CONTINUE 00676 IF( ILO.NE.ILO1 ) 00677 $ RESULT( 8 ) = ULPINV 00678 IF( IHI.NE.IHI1 ) 00679 $ RESULT( 8 ) = ULPINV 00680 IF( ABNRM.NE.ABNRM1 ) 00681 $ RESULT( 8 ) = ULPINV 00682 END IF 00683 * 00684 * Do Test (9) again 00685 * 00686 IF( ISENS.EQ.2 .AND. N.GT.1 ) THEN 00687 DO 130 J = 1, N 00688 IF( RCONDV( J ).NE.RCNDV1( J ) ) 00689 $ RESULT( 9 ) = ULPINV 00690 130 CONTINUE 00691 END IF 00692 * 00693 * Compute eigenvalues and left eigenvectors, and test them 00694 * 00695 CALL SLACPY( 'F', N, N, A, LDA, H, LDA ) 00696 CALL SGEEVX( BALANC, 'V', 'N', SENSE, N, H, LDA, WR1, WI1, LRE, 00697 $ LDLRE, DUM, 1, ILO1, IHI1, SCALE1, ABNRM1, RCNDE1, 00698 $ RCNDV1, WORK, LWORK, IWORK, IINFO ) 00699 IF( IINFO.NE.0 ) THEN 00700 RESULT( 1 ) = ULPINV 00701 IF( JTYPE.NE.22 ) THEN 00702 WRITE( NOUNIT, FMT = 9998 )'SGEEVX4', IINFO, N, JTYPE, 00703 $ BALANC, ISEED 00704 ELSE 00705 WRITE( NOUNIT, FMT = 9999 )'SGEEVX4', IINFO, N, 00706 $ ISEED( 1 ) 00707 END IF 00708 INFO = ABS( IINFO ) 00709 GO TO 190 00710 END IF 00711 * 00712 * Do Test (5) again 00713 * 00714 DO 140 J = 1, N 00715 IF( WR( J ).NE.WR1( J ) .OR. WI( J ).NE.WI1( J ) ) 00716 $ RESULT( 5 ) = ULPINV 00717 140 CONTINUE 00718 * 00719 * Do Test (7) 00720 * 00721 DO 160 J = 1, N 00722 DO 150 JJ = 1, N 00723 IF( VL( J, JJ ).NE.LRE( J, JJ ) ) 00724 $ RESULT( 7 ) = ULPINV 00725 150 CONTINUE 00726 160 CONTINUE 00727 * 00728 * Do Test (8) again 00729 * 00730 IF( .NOT.NOBAL ) THEN 00731 DO 170 J = 1, N 00732 IF( SCALE( J ).NE.SCALE1( J ) ) 00733 $ RESULT( 8 ) = ULPINV 00734 170 CONTINUE 00735 IF( ILO.NE.ILO1 ) 00736 $ RESULT( 8 ) = ULPINV 00737 IF( IHI.NE.IHI1 ) 00738 $ RESULT( 8 ) = ULPINV 00739 IF( ABNRM.NE.ABNRM1 ) 00740 $ RESULT( 8 ) = ULPINV 00741 END IF 00742 * 00743 * Do Test (9) again 00744 * 00745 IF( ISENS.EQ.2 .AND. N.GT.1 ) THEN 00746 DO 180 J = 1, N 00747 IF( RCONDV( J ).NE.RCNDV1( J ) ) 00748 $ RESULT( 9 ) = ULPINV 00749 180 CONTINUE 00750 END IF 00751 * 00752 190 CONTINUE 00753 * 00754 200 CONTINUE 00755 * 00756 * If COMP, compare condition numbers to precomputed ones 00757 * 00758 IF( COMP ) THEN 00759 CALL SLACPY( 'F', N, N, A, LDA, H, LDA ) 00760 CALL SGEEVX( 'N', 'V', 'V', 'B', N, H, LDA, WR, WI, VL, LDVL, 00761 $ VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV, 00762 $ WORK, LWORK, IWORK, IINFO ) 00763 IF( IINFO.NE.0 ) THEN 00764 RESULT( 1 ) = ULPINV 00765 WRITE( NOUNIT, FMT = 9999 )'SGEEVX5', IINFO, N, ISEED( 1 ) 00766 INFO = ABS( IINFO ) 00767 GO TO 250 00768 END IF 00769 * 00770 * Sort eigenvalues and condition numbers lexicographically 00771 * to compare with inputs 00772 * 00773 DO 220 I = 1, N - 1 00774 KMIN = I 00775 VRMIN = WR( I ) 00776 VIMIN = WI( I ) 00777 DO 210 J = I + 1, N 00778 IF( WR( J ).LT.VRMIN ) THEN 00779 KMIN = J 00780 VRMIN = WR( J ) 00781 VIMIN = WI( J ) 00782 END IF 00783 210 CONTINUE 00784 WR( KMIN ) = WR( I ) 00785 WI( KMIN ) = WI( I ) 00786 WR( I ) = VRMIN 00787 WI( I ) = VIMIN 00788 VRMIN = RCONDE( KMIN ) 00789 RCONDE( KMIN ) = RCONDE( I ) 00790 RCONDE( I ) = VRMIN 00791 VRMIN = RCONDV( KMIN ) 00792 RCONDV( KMIN ) = RCONDV( I ) 00793 RCONDV( I ) = VRMIN 00794 220 CONTINUE 00795 * 00796 * Compare condition numbers for eigenvectors 00797 * taking their condition numbers into account 00798 * 00799 RESULT( 10 ) = ZERO 00800 EPS = MAX( EPSIN, ULP ) 00801 V = MAX( REAL( N )*EPS*ABNRM, SMLNUM ) 00802 IF( ABNRM.EQ.ZERO ) 00803 $ V = ONE 00804 DO 230 I = 1, N 00805 IF( V.GT.RCONDV( I )*RCONDE( I ) ) THEN 00806 TOL = RCONDV( I ) 00807 ELSE 00808 TOL = V / RCONDE( I ) 00809 END IF 00810 IF( V.GT.RCDVIN( I )*RCDEIN( I ) ) THEN 00811 TOLIN = RCDVIN( I ) 00812 ELSE 00813 TOLIN = V / RCDEIN( I ) 00814 END IF 00815 TOL = MAX( TOL, SMLNUM / EPS ) 00816 TOLIN = MAX( TOLIN, SMLNUM / EPS ) 00817 IF( EPS*( RCDVIN( I )-TOLIN ).GT.RCONDV( I )+TOL ) THEN 00818 VMAX = ONE / EPS 00819 ELSE IF( RCDVIN( I )-TOLIN.GT.RCONDV( I )+TOL ) THEN 00820 VMAX = ( RCDVIN( I )-TOLIN ) / ( RCONDV( I )+TOL ) 00821 ELSE IF( RCDVIN( I )+TOLIN.LT.EPS*( RCONDV( I )-TOL ) ) THEN 00822 VMAX = ONE / EPS 00823 ELSE IF( RCDVIN( I )+TOLIN.LT.RCONDV( I )-TOL ) THEN 00824 VMAX = ( RCONDV( I )-TOL ) / ( RCDVIN( I )+TOLIN ) 00825 ELSE 00826 VMAX = ONE 00827 END IF 00828 RESULT( 10 ) = MAX( RESULT( 10 ), VMAX ) 00829 230 CONTINUE 00830 * 00831 * Compare condition numbers for eigenvalues 00832 * taking their condition numbers into account 00833 * 00834 RESULT( 11 ) = ZERO 00835 DO 240 I = 1, N 00836 IF( V.GT.RCONDV( I ) ) THEN 00837 TOL = ONE 00838 ELSE 00839 TOL = V / RCONDV( I ) 00840 END IF 00841 IF( V.GT.RCDVIN( I ) ) THEN 00842 TOLIN = ONE 00843 ELSE 00844 TOLIN = V / RCDVIN( I ) 00845 END IF 00846 TOL = MAX( TOL, SMLNUM / EPS ) 00847 TOLIN = MAX( TOLIN, SMLNUM / EPS ) 00848 IF( EPS*( RCDEIN( I )-TOLIN ).GT.RCONDE( I )+TOL ) THEN 00849 VMAX = ONE / EPS 00850 ELSE IF( RCDEIN( I )-TOLIN.GT.RCONDE( I )+TOL ) THEN 00851 VMAX = ( RCDEIN( I )-TOLIN ) / ( RCONDE( I )+TOL ) 00852 ELSE IF( RCDEIN( I )+TOLIN.LT.EPS*( RCONDE( I )-TOL ) ) THEN 00853 VMAX = ONE / EPS 00854 ELSE IF( RCDEIN( I )+TOLIN.LT.RCONDE( I )-TOL ) THEN 00855 VMAX = ( RCONDE( I )-TOL ) / ( RCDEIN( I )+TOLIN ) 00856 ELSE 00857 VMAX = ONE 00858 END IF 00859 RESULT( 11 ) = MAX( RESULT( 11 ), VMAX ) 00860 240 CONTINUE 00861 250 CONTINUE 00862 * 00863 END IF 00864 * 00865 9999 FORMAT( ' SGET23: ', A, ' returned INFO=', I6, '.', / 9X, 'N=', 00866 $ I6, ', INPUT EXAMPLE NUMBER = ', I4 ) 00867 9998 FORMAT( ' SGET23: ', A, ' returned INFO=', I6, '.', / 9X, 'N=', 00868 $ I6, ', JTYPE=', I6, ', BALANC = ', A, ', ISEED=(', 00869 $ 3( I5, ',' ), I5, ')' ) 00870 * 00871 RETURN 00872 * 00873 * End of SGET23 00874 * 00875 END