![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DDISNA 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download DDISNA + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ddisna.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ddisna.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ddisna.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE DDISNA( JOB, M, N, D, SEP, INFO ) 00022 * 00023 * .. Scalar Arguments .. 00024 * CHARACTER JOB 00025 * INTEGER INFO, M, N 00026 * .. 00027 * .. Array Arguments .. 00028 * DOUBLE PRECISION D( * ), SEP( * ) 00029 * .. 00030 * 00031 * 00032 *> \par Purpose: 00033 * ============= 00034 *> 00035 *> \verbatim 00036 *> 00037 *> DDISNA computes the reciprocal condition numbers for the eigenvectors 00038 *> of a real symmetric or complex Hermitian matrix or for the left or 00039 *> right singular vectors of a general m-by-n matrix. The reciprocal 00040 *> condition number is the 'gap' between the corresponding eigenvalue or 00041 *> singular value and the nearest other one. 00042 *> 00043 *> The bound on the error, measured by angle in radians, in the I-th 00044 *> computed vector is given by 00045 *> 00046 *> DLAMCH( 'E' ) * ( ANORM / SEP( I ) ) 00047 *> 00048 *> where ANORM = 2-norm(A) = max( abs( D(j) ) ). SEP(I) is not allowed 00049 *> to be smaller than DLAMCH( 'E' )*ANORM in order to limit the size of 00050 *> the error bound. 00051 *> 00052 *> DDISNA may also be used to compute error bounds for eigenvectors of 00053 *> the generalized symmetric definite eigenproblem. 00054 *> \endverbatim 00055 * 00056 * Arguments: 00057 * ========== 00058 * 00059 *> \param[in] JOB 00060 *> \verbatim 00061 *> JOB is CHARACTER*1 00062 *> Specifies for which problem the reciprocal condition numbers 00063 *> should be computed: 00064 *> = 'E': the eigenvectors of a symmetric/Hermitian matrix; 00065 *> = 'L': the left singular vectors of a general matrix; 00066 *> = 'R': the right singular vectors of a general matrix. 00067 *> \endverbatim 00068 *> 00069 *> \param[in] M 00070 *> \verbatim 00071 *> M is INTEGER 00072 *> The number of rows of the matrix. M >= 0. 00073 *> \endverbatim 00074 *> 00075 *> \param[in] N 00076 *> \verbatim 00077 *> N is INTEGER 00078 *> If JOB = 'L' or 'R', the number of columns of the matrix, 00079 *> in which case N >= 0. Ignored if JOB = 'E'. 00080 *> \endverbatim 00081 *> 00082 *> \param[in] D 00083 *> \verbatim 00084 *> D is DOUBLE PRECISION array, dimension (M) if JOB = 'E' 00085 *> dimension (min(M,N)) if JOB = 'L' or 'R' 00086 *> The eigenvalues (if JOB = 'E') or singular values (if JOB = 00087 *> 'L' or 'R') of the matrix, in either increasing or decreasing 00088 *> order. If singular values, they must be non-negative. 00089 *> \endverbatim 00090 *> 00091 *> \param[out] SEP 00092 *> \verbatim 00093 *> SEP is DOUBLE PRECISION array, dimension (M) if JOB = 'E' 00094 *> dimension (min(M,N)) if JOB = 'L' or 'R' 00095 *> The reciprocal condition numbers of the vectors. 00096 *> \endverbatim 00097 *> 00098 *> \param[out] INFO 00099 *> \verbatim 00100 *> INFO is INTEGER 00101 *> = 0: successful exit. 00102 *> < 0: if INFO = -i, the i-th argument had an illegal value. 00103 *> \endverbatim 00104 * 00105 * Authors: 00106 * ======== 00107 * 00108 *> \author Univ. of Tennessee 00109 *> \author Univ. of California Berkeley 00110 *> \author Univ. of Colorado Denver 00111 *> \author NAG Ltd. 00112 * 00113 *> \date November 2011 00114 * 00115 *> \ingroup auxOTHERcomputational 00116 * 00117 * ===================================================================== 00118 SUBROUTINE DDISNA( JOB, M, N, D, SEP, INFO ) 00119 * 00120 * -- LAPACK computational routine (version 3.4.0) -- 00121 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00122 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00123 * November 2011 00124 * 00125 * .. Scalar Arguments .. 00126 CHARACTER JOB 00127 INTEGER INFO, M, N 00128 * .. 00129 * .. Array Arguments .. 00130 DOUBLE PRECISION D( * ), SEP( * ) 00131 * .. 00132 * 00133 * ===================================================================== 00134 * 00135 * .. Parameters .. 00136 DOUBLE PRECISION ZERO 00137 PARAMETER ( ZERO = 0.0D+0 ) 00138 * .. 00139 * .. Local Scalars .. 00140 LOGICAL DECR, EIGEN, INCR, LEFT, RIGHT, SING 00141 INTEGER I, K 00142 DOUBLE PRECISION ANORM, EPS, NEWGAP, OLDGAP, SAFMIN, THRESH 00143 * .. 00144 * .. External Functions .. 00145 LOGICAL LSAME 00146 DOUBLE PRECISION DLAMCH 00147 EXTERNAL LSAME, DLAMCH 00148 * .. 00149 * .. Intrinsic Functions .. 00150 INTRINSIC ABS, MAX, MIN 00151 * .. 00152 * .. External Subroutines .. 00153 EXTERNAL XERBLA 00154 * .. 00155 * .. Executable Statements .. 00156 * 00157 * Test the input arguments 00158 * 00159 INFO = 0 00160 EIGEN = LSAME( JOB, 'E' ) 00161 LEFT = LSAME( JOB, 'L' ) 00162 RIGHT = LSAME( JOB, 'R' ) 00163 SING = LEFT .OR. RIGHT 00164 IF( EIGEN ) THEN 00165 K = M 00166 ELSE IF( SING ) THEN 00167 K = MIN( M, N ) 00168 END IF 00169 IF( .NOT.EIGEN .AND. .NOT.SING ) THEN 00170 INFO = -1 00171 ELSE IF( M.LT.0 ) THEN 00172 INFO = -2 00173 ELSE IF( K.LT.0 ) THEN 00174 INFO = -3 00175 ELSE 00176 INCR = .TRUE. 00177 DECR = .TRUE. 00178 DO 10 I = 1, K - 1 00179 IF( INCR ) 00180 $ INCR = INCR .AND. D( I ).LE.D( I+1 ) 00181 IF( DECR ) 00182 $ DECR = DECR .AND. D( I ).GE.D( I+1 ) 00183 10 CONTINUE 00184 IF( SING .AND. K.GT.0 ) THEN 00185 IF( INCR ) 00186 $ INCR = INCR .AND. ZERO.LE.D( 1 ) 00187 IF( DECR ) 00188 $ DECR = DECR .AND. D( K ).GE.ZERO 00189 END IF 00190 IF( .NOT.( INCR .OR. DECR ) ) 00191 $ INFO = -4 00192 END IF 00193 IF( INFO.NE.0 ) THEN 00194 CALL XERBLA( 'DDISNA', -INFO ) 00195 RETURN 00196 END IF 00197 * 00198 * Quick return if possible 00199 * 00200 IF( K.EQ.0 ) 00201 $ RETURN 00202 * 00203 * Compute reciprocal condition numbers 00204 * 00205 IF( K.EQ.1 ) THEN 00206 SEP( 1 ) = DLAMCH( 'O' ) 00207 ELSE 00208 OLDGAP = ABS( D( 2 )-D( 1 ) ) 00209 SEP( 1 ) = OLDGAP 00210 DO 20 I = 2, K - 1 00211 NEWGAP = ABS( D( I+1 )-D( I ) ) 00212 SEP( I ) = MIN( OLDGAP, NEWGAP ) 00213 OLDGAP = NEWGAP 00214 20 CONTINUE 00215 SEP( K ) = OLDGAP 00216 END IF 00217 IF( SING ) THEN 00218 IF( ( LEFT .AND. M.GT.N ) .OR. ( RIGHT .AND. M.LT.N ) ) THEN 00219 IF( INCR ) 00220 $ SEP( 1 ) = MIN( SEP( 1 ), D( 1 ) ) 00221 IF( DECR ) 00222 $ SEP( K ) = MIN( SEP( K ), D( K ) ) 00223 END IF 00224 END IF 00225 * 00226 * Ensure that reciprocal condition numbers are not less than 00227 * threshold, in order to limit the size of the error bound 00228 * 00229 EPS = DLAMCH( 'E' ) 00230 SAFMIN = DLAMCH( 'S' ) 00231 ANORM = MAX( ABS( D( 1 ) ), ABS( D( K ) ) ) 00232 IF( ANORM.EQ.ZERO ) THEN 00233 THRESH = EPS 00234 ELSE 00235 THRESH = MAX( EPS*ANORM, SAFMIN ) 00236 END IF 00237 DO 30 I = 1, K 00238 SEP( I ) = MAX( SEP( I ), THRESH ) 00239 30 CONTINUE 00240 * 00241 RETURN 00242 * 00243 * End of DDISNA 00244 * 00245 END