![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b SLARRK 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download SLARRK + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slarrk.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slarrk.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slarrk.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE SLARRK( N, IW, GL, GU, 00022 * D, E2, PIVMIN, RELTOL, W, WERR, INFO) 00023 * 00024 * .. Scalar Arguments .. 00025 * INTEGER INFO, IW, N 00026 * REAL PIVMIN, RELTOL, GL, GU, W, WERR 00027 * .. 00028 * .. Array Arguments .. 00029 * REAL D( * ), E2( * ) 00030 * .. 00031 * 00032 * 00033 *> \par Purpose: 00034 * ============= 00035 *> 00036 *> \verbatim 00037 *> 00038 *> SLARRK computes one eigenvalue of a symmetric tridiagonal 00039 *> matrix T to suitable accuracy. This is an auxiliary code to be 00040 *> called from SSTEMR. 00041 *> 00042 *> To avoid overflow, the matrix must be scaled so that its 00043 *> largest element is no greater than overflow**(1/2) * underflow**(1/4) in absolute value, and for greatest 00044 *> accuracy, it should not be much smaller than that. 00045 *> 00046 *> See W. Kahan "Accurate Eigenvalues of a Symmetric Tridiagonal 00047 *> Matrix", Report CS41, Computer Science Dept., Stanford 00048 *> University, July 21, 1966. 00049 *> \endverbatim 00050 * 00051 * Arguments: 00052 * ========== 00053 * 00054 *> \param[in] N 00055 *> \verbatim 00056 *> N is INTEGER 00057 *> The order of the tridiagonal matrix T. N >= 0. 00058 *> \endverbatim 00059 *> 00060 *> \param[in] IW 00061 *> \verbatim 00062 *> IW is INTEGER 00063 *> The index of the eigenvalues to be returned. 00064 *> \endverbatim 00065 *> 00066 *> \param[in] GL 00067 *> \verbatim 00068 *> GL is REAL 00069 *> \endverbatim 00070 *> 00071 *> \param[in] GU 00072 *> \verbatim 00073 *> GU is REAL 00074 *> An upper and a lower bound on the eigenvalue. 00075 *> \endverbatim 00076 *> 00077 *> \param[in] D 00078 *> \verbatim 00079 *> D is REAL array, dimension (N) 00080 *> The n diagonal elements of the tridiagonal matrix T. 00081 *> \endverbatim 00082 *> 00083 *> \param[in] E2 00084 *> \verbatim 00085 *> E2 is REAL array, dimension (N-1) 00086 *> The (n-1) squared off-diagonal elements of the tridiagonal matrix T. 00087 *> \endverbatim 00088 *> 00089 *> \param[in] PIVMIN 00090 *> \verbatim 00091 *> PIVMIN is REAL 00092 *> The minimum pivot allowed in the Sturm sequence for T. 00093 *> \endverbatim 00094 *> 00095 *> \param[in] RELTOL 00096 *> \verbatim 00097 *> RELTOL is REAL 00098 *> The minimum relative width of an interval. When an interval 00099 *> is narrower than RELTOL times the larger (in 00100 *> magnitude) endpoint, then it is considered to be 00101 *> sufficiently small, i.e., converged. Note: this should 00102 *> always be at least radix*machine epsilon. 00103 *> \endverbatim 00104 *> 00105 *> \param[out] W 00106 *> \verbatim 00107 *> W is REAL 00108 *> \endverbatim 00109 *> 00110 *> \param[out] WERR 00111 *> \verbatim 00112 *> WERR is REAL 00113 *> The error bound on the corresponding eigenvalue approximation 00114 *> in W. 00115 *> \endverbatim 00116 *> 00117 *> \param[out] INFO 00118 *> \verbatim 00119 *> INFO is INTEGER 00120 *> = 0: Eigenvalue converged 00121 *> = -1: Eigenvalue did NOT converge 00122 *> \endverbatim 00123 * 00124 *> \par Internal Parameters: 00125 * ========================= 00126 *> 00127 *> \verbatim 00128 *> FUDGE REAL , default = 2 00129 *> A "fudge factor" to widen the Gershgorin intervals. 00130 *> \endverbatim 00131 * 00132 * Authors: 00133 * ======== 00134 * 00135 *> \author Univ. of Tennessee 00136 *> \author Univ. of California Berkeley 00137 *> \author Univ. of Colorado Denver 00138 *> \author NAG Ltd. 00139 * 00140 *> \date November 2011 00141 * 00142 *> \ingroup auxOTHERauxiliary 00143 * 00144 * ===================================================================== 00145 SUBROUTINE SLARRK( N, IW, GL, GU, 00146 $ D, E2, PIVMIN, RELTOL, W, WERR, INFO) 00147 * 00148 * -- LAPACK auxiliary routine (version 3.4.0) -- 00149 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00150 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00151 * November 2011 00152 * 00153 * .. Scalar Arguments .. 00154 INTEGER INFO, IW, N 00155 REAL PIVMIN, RELTOL, GL, GU, W, WERR 00156 * .. 00157 * .. Array Arguments .. 00158 REAL D( * ), E2( * ) 00159 * .. 00160 * 00161 * ===================================================================== 00162 * 00163 * .. Parameters .. 00164 REAL FUDGE, HALF, TWO, ZERO 00165 PARAMETER ( HALF = 0.5E0, TWO = 2.0E0, 00166 $ FUDGE = TWO, ZERO = 0.0E0 ) 00167 * .. 00168 * .. Local Scalars .. 00169 INTEGER I, IT, ITMAX, NEGCNT 00170 REAL ATOLI, EPS, LEFT, MID, RIGHT, RTOLI, TMP1, 00171 $ TMP2, TNORM 00172 * .. 00173 * .. External Functions .. 00174 REAL SLAMCH 00175 EXTERNAL SLAMCH 00176 * .. 00177 * .. Intrinsic Functions .. 00178 INTRINSIC ABS, INT, LOG, MAX 00179 * .. 00180 * .. Executable Statements .. 00181 * 00182 * Get machine constants 00183 EPS = SLAMCH( 'P' ) 00184 00185 TNORM = MAX( ABS( GL ), ABS( GU ) ) 00186 RTOLI = RELTOL 00187 ATOLI = FUDGE*TWO*PIVMIN 00188 00189 ITMAX = INT( ( LOG( TNORM+PIVMIN )-LOG( PIVMIN ) ) / 00190 $ LOG( TWO ) ) + 2 00191 00192 INFO = -1 00193 00194 LEFT = GL - FUDGE*TNORM*EPS*N - FUDGE*TWO*PIVMIN 00195 RIGHT = GU + FUDGE*TNORM*EPS*N + FUDGE*TWO*PIVMIN 00196 IT = 0 00197 00198 10 CONTINUE 00199 * 00200 * Check if interval converged or maximum number of iterations reached 00201 * 00202 TMP1 = ABS( RIGHT - LEFT ) 00203 TMP2 = MAX( ABS(RIGHT), ABS(LEFT) ) 00204 IF( TMP1.LT.MAX( ATOLI, PIVMIN, RTOLI*TMP2 ) ) THEN 00205 INFO = 0 00206 GOTO 30 00207 ENDIF 00208 IF(IT.GT.ITMAX) 00209 $ GOTO 30 00210 00211 * 00212 * Count number of negative pivots for mid-point 00213 * 00214 IT = IT + 1 00215 MID = HALF * (LEFT + RIGHT) 00216 NEGCNT = 0 00217 TMP1 = D( 1 ) - MID 00218 IF( ABS( TMP1 ).LT.PIVMIN ) 00219 $ TMP1 = -PIVMIN 00220 IF( TMP1.LE.ZERO ) 00221 $ NEGCNT = NEGCNT + 1 00222 * 00223 DO 20 I = 2, N 00224 TMP1 = D( I ) - E2( I-1 ) / TMP1 - MID 00225 IF( ABS( TMP1 ).LT.PIVMIN ) 00226 $ TMP1 = -PIVMIN 00227 IF( TMP1.LE.ZERO ) 00228 $ NEGCNT = NEGCNT + 1 00229 20 CONTINUE 00230 00231 IF(NEGCNT.GE.IW) THEN 00232 RIGHT = MID 00233 ELSE 00234 LEFT = MID 00235 ENDIF 00236 GOTO 10 00237 00238 30 CONTINUE 00239 * 00240 * Converged or maximum number of iterations reached 00241 * 00242 W = HALF * (LEFT + RIGHT) 00243 WERR = HALF * ABS( RIGHT - LEFT ) 00244 00245 RETURN 00246 * 00247 * End of SLARRK 00248 * 00249 END