*> \brief \b DGET53 * * =========== DOCUMENTATION =========== * * Online html documentation available at * http://www.netlib.org/lapack/explore-html/ * * Definition: * =========== * * SUBROUTINE DGET53( A, LDA, B, LDB, SCALE, WR, WI, RESULT, INFO ) * * .. Scalar Arguments .. * INTEGER INFO, LDA, LDB * DOUBLE PRECISION RESULT, SCALE, WI, WR * .. * .. Array Arguments .. * DOUBLE PRECISION A( LDA, * ), B( LDB, * ) * .. * * *> \par Purpose: * ============= *> *> \verbatim *> *> DGET53 checks the generalized eigenvalues computed by DLAG2. *> *> The basic test for an eigenvalue is: *> *> | det( s A - w B ) | *> RESULT = --------------------------------------------------- *> ulp max( s norm(A), |w| norm(B) )*norm( s A - w B ) *> *> Two "safety checks" are performed: *> *> (1) ulp*max( s*norm(A), |w|*norm(B) ) must be at least *> safe_minimum. This insures that the test performed is *> not essentially det(0*A + 0*B)=0. *> *> (2) s*norm(A) + |w|*norm(B) must be less than 1/safe_minimum. *> This insures that s*A - w*B will not overflow. *> *> If these tests are not passed, then s and w are scaled and *> tested anyway, if this is possible. *> \endverbatim * * Arguments: * ========== * *> \param[in] A *> \verbatim *> A is DOUBLE PRECISION array, dimension (LDA, 2) *> The 2x2 matrix A. *> \endverbatim *> *> \param[in] LDA *> \verbatim *> LDA is INTEGER *> The leading dimension of A. It must be at least 2. *> \endverbatim *> *> \param[in] B *> \verbatim *> B is DOUBLE PRECISION array, dimension (LDB, N) *> The 2x2 upper-triangular matrix B. *> \endverbatim *> *> \param[in] LDB *> \verbatim *> LDB is INTEGER *> The leading dimension of B. It must be at least 2. *> \endverbatim *> *> \param[in] SCALE *> \verbatim *> SCALE is DOUBLE PRECISION *> The "scale factor" s in the formula s A - w B . It is *> assumed to be non-negative. *> \endverbatim *> *> \param[in] WR *> \verbatim *> WR is DOUBLE PRECISION *> The real part of the eigenvalue w in the formula *> s A - w B . *> \endverbatim *> *> \param[in] WI *> \verbatim *> WI is DOUBLE PRECISION *> The imaginary part of the eigenvalue w in the formula *> s A - w B . *> \endverbatim *> *> \param[out] RESULT *> \verbatim *> RESULT is DOUBLE PRECISION *> If INFO is 2 or less, the value computed by the test *> described above. *> If INFO=3, this will just be 1/ulp. *> \endverbatim *> *> \param[out] INFO *> \verbatim *> INFO is INTEGER *> =0: The input data pass the "safety checks". *> =1: s*norm(A) + |w|*norm(B) > 1/safe_minimum. *> =2: ulp*max( s*norm(A), |w|*norm(B) ) < safe_minimum *> =3: same as INFO=2, but s and w could not be scaled so *> as to compute the test. *> \endverbatim * * Authors: * ======== * *> \author Univ. of Tennessee *> \author Univ. of California Berkeley *> \author Univ. of Colorado Denver *> \author NAG Ltd. * *> \date November 2011 * *> \ingroup double_eig * * ===================================================================== SUBROUTINE DGET53( A, LDA, B, LDB, SCALE, WR, WI, RESULT, INFO ) * * -- LAPACK test routine (version 3.4.0) -- * -- LAPACK is a software package provided by Univ. of Tennessee, -- * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- * November 2011 * * .. Scalar Arguments .. INTEGER INFO, LDA, LDB DOUBLE PRECISION RESULT, SCALE, WI, WR * .. * .. Array Arguments .. DOUBLE PRECISION A( LDA, * ), B( LDB, * ) * .. * * ===================================================================== * * .. Parameters .. DOUBLE PRECISION ZERO, ONE PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 ) * .. * .. Local Scalars .. DOUBLE PRECISION ABSW, ANORM, BNORM, CI11, CI12, CI22, CNORM, $ CR11, CR12, CR21, CR22, CSCALE, DETI, DETR, S1, $ SAFMIN, SCALES, SIGMIN, TEMP, ULP, WIS, WRS * .. * .. External Functions .. DOUBLE PRECISION DLAMCH EXTERNAL DLAMCH * .. * .. Intrinsic Functions .. INTRINSIC ABS, MAX, SQRT * .. * .. Executable Statements .. * * Initialize * INFO = 0 RESULT = ZERO SCALES = SCALE WRS = WR WIS = WI * * Machine constants and norms * SAFMIN = DLAMCH( 'Safe minimum' ) ULP = DLAMCH( 'Epsilon' )*DLAMCH( 'Base' ) ABSW = ABS( WRS ) + ABS( WIS ) ANORM = MAX( ABS( A( 1, 1 ) )+ABS( A( 2, 1 ) ), $ ABS( A( 1, 2 ) )+ABS( A( 2, 2 ) ), SAFMIN ) BNORM = MAX( ABS( B( 1, 1 ) ), ABS( B( 1, 2 ) )+ABS( B( 2, 2 ) ), $ SAFMIN ) * * Check for possible overflow. * TEMP = ( SAFMIN*BNORM )*ABSW + ( SAFMIN*ANORM )*SCALES IF( TEMP.GE.ONE ) THEN * * Scale down to avoid overflow * INFO = 1 TEMP = ONE / TEMP SCALES = SCALES*TEMP WRS = WRS*TEMP WIS = WIS*TEMP ABSW = ABS( WRS ) + ABS( WIS ) END IF S1 = MAX( ULP*MAX( SCALES*ANORM, ABSW*BNORM ), $ SAFMIN*MAX( SCALES, ABSW ) ) * * Check for W and SCALE essentially zero. * IF( S1.LT.SAFMIN ) THEN INFO = 2 IF( SCALES.LT.SAFMIN .AND. ABSW.LT.SAFMIN ) THEN INFO = 3 RESULT = ONE / ULP RETURN END IF * * Scale up to avoid underflow * TEMP = ONE / MAX( SCALES*ANORM+ABSW*BNORM, SAFMIN ) SCALES = SCALES*TEMP WRS = WRS*TEMP WIS = WIS*TEMP ABSW = ABS( WRS ) + ABS( WIS ) S1 = MAX( ULP*MAX( SCALES*ANORM, ABSW*BNORM ), $ SAFMIN*MAX( SCALES, ABSW ) ) IF( S1.LT.SAFMIN ) THEN INFO = 3 RESULT = ONE / ULP RETURN END IF END IF * * Compute C = s A - w B * CR11 = SCALES*A( 1, 1 ) - WRS*B( 1, 1 ) CI11 = -WIS*B( 1, 1 ) CR21 = SCALES*A( 2, 1 ) CR12 = SCALES*A( 1, 2 ) - WRS*B( 1, 2 ) CI12 = -WIS*B( 1, 2 ) CR22 = SCALES*A( 2, 2 ) - WRS*B( 2, 2 ) CI22 = -WIS*B( 2, 2 ) * * Compute the smallest singular value of s A - w B: * * |det( s A - w B )| * sigma_min = ------------------ * norm( s A - w B ) * CNORM = MAX( ABS( CR11 )+ABS( CI11 )+ABS( CR21 ), $ ABS( CR12 )+ABS( CI12 )+ABS( CR22 )+ABS( CI22 ), SAFMIN ) CSCALE = ONE / SQRT( CNORM ) DETR = ( CSCALE*CR11 )*( CSCALE*CR22 ) - $ ( CSCALE*CI11 )*( CSCALE*CI22 ) - $ ( CSCALE*CR12 )*( CSCALE*CR21 ) DETI = ( CSCALE*CR11 )*( CSCALE*CI22 ) + $ ( CSCALE*CI11 )*( CSCALE*CR22 ) - $ ( CSCALE*CI12 )*( CSCALE*CR21 ) SIGMIN = ABS( DETR ) + ABS( DETI ) RESULT = SIGMIN / S1 RETURN * * End of DGET53 * END