1*> \brief \b DLAGTF computes an LU factorization of a matrix T-λI, where T is a general tridiagonal matrix, and λ a scalar, using partial pivoting with row interchanges. 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8*> \htmlonly 9*> Download DLAGTF + dependencies 10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlagtf.f"> 11*> [TGZ]</a> 12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlagtf.f"> 13*> [ZIP]</a> 14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlagtf.f"> 15*> [TXT]</a> 16*> \endhtmlonly 17* 18* Definition: 19* =========== 20* 21* SUBROUTINE DLAGTF( N, A, LAMBDA, B, C, TOL, D, IN, INFO ) 22* 23* .. Scalar Arguments .. 24* INTEGER INFO, N 25* DOUBLE PRECISION LAMBDA, TOL 26* .. 27* .. Array Arguments .. 28* INTEGER IN( * ) 29* DOUBLE PRECISION A( * ), B( * ), C( * ), D( * ) 30* .. 31* 32* 33*> \par Purpose: 34* ============= 35*> 36*> \verbatim 37*> 38*> DLAGTF factorizes the matrix (T - lambda*I), where T is an n by n 39*> tridiagonal matrix and lambda is a scalar, as 40*> 41*> T - lambda*I = PLU, 42*> 43*> where P is a permutation matrix, L is a unit lower tridiagonal matrix 44*> with at most one non-zero sub-diagonal elements per column and U is 45*> an upper triangular matrix with at most two non-zero super-diagonal 46*> elements per column. 47*> 48*> The factorization is obtained by Gaussian elimination with partial 49*> pivoting and implicit row scaling. 50*> 51*> The parameter LAMBDA is included in the routine so that DLAGTF may 52*> be used, in conjunction with DLAGTS, to obtain eigenvectors of T by 53*> inverse iteration. 54*> \endverbatim 55* 56* Arguments: 57* ========== 58* 59*> \param[in] N 60*> \verbatim 61*> N is INTEGER 62*> The order of the matrix T. 63*> \endverbatim 64*> 65*> \param[in,out] A 66*> \verbatim 67*> A is DOUBLE PRECISION array, dimension (N) 68*> On entry, A must contain the diagonal elements of T. 69*> 70*> On exit, A is overwritten by the n diagonal elements of the 71*> upper triangular matrix U of the factorization of T. 72*> \endverbatim 73*> 74*> \param[in] LAMBDA 75*> \verbatim 76*> LAMBDA is DOUBLE PRECISION 77*> On entry, the scalar lambda. 78*> \endverbatim 79*> 80*> \param[in,out] B 81*> \verbatim 82*> B is DOUBLE PRECISION array, dimension (N-1) 83*> On entry, B must contain the (n-1) super-diagonal elements of 84*> T. 85*> 86*> On exit, B is overwritten by the (n-1) super-diagonal 87*> elements of the matrix U of the factorization of T. 88*> \endverbatim 89*> 90*> \param[in,out] C 91*> \verbatim 92*> C is DOUBLE PRECISION array, dimension (N-1) 93*> On entry, C must contain the (n-1) sub-diagonal elements of 94*> T. 95*> 96*> On exit, C is overwritten by the (n-1) sub-diagonal elements 97*> of the matrix L of the factorization of T. 98*> \endverbatim 99*> 100*> \param[in] TOL 101*> \verbatim 102*> TOL is DOUBLE PRECISION 103*> On entry, a relative tolerance used to indicate whether or 104*> not the matrix (T - lambda*I) is nearly singular. TOL should 105*> normally be chose as approximately the largest relative error 106*> in the elements of T. For example, if the elements of T are 107*> correct to about 4 significant figures, then TOL should be 108*> set to about 5*10**(-4). If TOL is supplied as less than eps, 109*> where eps is the relative machine precision, then the value 110*> eps is used in place of TOL. 111*> \endverbatim 112*> 113*> \param[out] D 114*> \verbatim 115*> D is DOUBLE PRECISION array, dimension (N-2) 116*> On exit, D is overwritten by the (n-2) second super-diagonal 117*> elements of the matrix U of the factorization of T. 118*> \endverbatim 119*> 120*> \param[out] IN 121*> \verbatim 122*> IN is INTEGER array, dimension (N) 123*> On exit, IN contains details of the permutation matrix P. If 124*> an interchange occurred at the kth step of the elimination, 125*> then IN(k) = 1, otherwise IN(k) = 0. The element IN(n) 126*> returns the smallest positive integer j such that 127*> 128*> abs( u(j,j) ) <= norm( (T - lambda*I)(j) )*TOL, 129*> 130*> where norm( A(j) ) denotes the sum of the absolute values of 131*> the jth row of the matrix A. If no such j exists then IN(n) 132*> is returned as zero. If IN(n) is returned as positive, then a 133*> diagonal element of U is small, indicating that 134*> (T - lambda*I) is singular or nearly singular, 135*> \endverbatim 136*> 137*> \param[out] INFO 138*> \verbatim 139*> INFO is INTEGER 140*> = 0: successful exit 141*> < 0: if INFO = -k, the kth argument had an illegal value 142*> \endverbatim 143* 144* Authors: 145* ======== 146* 147*> \author Univ. of Tennessee 148*> \author Univ. of California Berkeley 149*> \author Univ. of Colorado Denver 150*> \author NAG Ltd. 151* 152*> \date December 2016 153* 154*> \ingroup auxOTHERcomputational 155* 156* ===================================================================== 157 SUBROUTINE DLAGTF( N, A, LAMBDA, B, C, TOL, D, IN, INFO ) 158* 159* -- LAPACK computational routine (version 3.7.0) -- 160* -- LAPACK is a software package provided by Univ. of Tennessee, -- 161* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 162* December 2016 163* 164* .. Scalar Arguments .. 165 INTEGER INFO, N 166 DOUBLE PRECISION LAMBDA, TOL 167* .. 168* .. Array Arguments .. 169 INTEGER IN( * ) 170 DOUBLE PRECISION A( * ), B( * ), C( * ), D( * ) 171* .. 172* 173* ===================================================================== 174* 175* .. Parameters .. 176 DOUBLE PRECISION ZERO 177 PARAMETER ( ZERO = 0.0D+0 ) 178* .. 179* .. Local Scalars .. 180 INTEGER K 181 DOUBLE PRECISION EPS, MULT, PIV1, PIV2, SCALE1, SCALE2, TEMP, TL 182* .. 183* .. Intrinsic Functions .. 184 INTRINSIC ABS, MAX 185* .. 186* .. External Functions .. 187 DOUBLE PRECISION DLAMCH 188 EXTERNAL DLAMCH 189* .. 190* .. External Subroutines .. 191 EXTERNAL XERBLA 192* .. 193* .. Executable Statements .. 194* 195 INFO = 0 196 IF( N.LT.0 ) THEN 197 INFO = -1 198 CALL XERBLA( 'DLAGTF', -INFO ) 199 RETURN 200 END IF 201* 202 IF( N.EQ.0 ) 203 $ RETURN 204* 205 A( 1 ) = A( 1 ) - LAMBDA 206 IN( N ) = 0 207 IF( N.EQ.1 ) THEN 208 IF( A( 1 ).EQ.ZERO ) 209 $ IN( 1 ) = 1 210 RETURN 211 END IF 212* 213 EPS = DLAMCH( 'Epsilon' ) 214* 215 TL = MAX( TOL, EPS ) 216 SCALE1 = ABS( A( 1 ) ) + ABS( B( 1 ) ) 217 DO 10 K = 1, N - 1 218 A( K+1 ) = A( K+1 ) - LAMBDA 219 SCALE2 = ABS( C( K ) ) + ABS( A( K+1 ) ) 220 IF( K.LT.( N-1 ) ) 221 $ SCALE2 = SCALE2 + ABS( B( K+1 ) ) 222 IF( A( K ).EQ.ZERO ) THEN 223 PIV1 = ZERO 224 ELSE 225 PIV1 = ABS( A( K ) ) / SCALE1 226 END IF 227 IF( C( K ).EQ.ZERO ) THEN 228 IN( K ) = 0 229 PIV2 = ZERO 230 SCALE1 = SCALE2 231 IF( K.LT.( N-1 ) ) 232 $ D( K ) = ZERO 233 ELSE 234 PIV2 = ABS( C( K ) ) / SCALE2 235 IF( PIV2.LE.PIV1 ) THEN 236 IN( K ) = 0 237 SCALE1 = SCALE2 238 C( K ) = C( K ) / A( K ) 239 A( K+1 ) = A( K+1 ) - C( K )*B( K ) 240 IF( K.LT.( N-1 ) ) 241 $ D( K ) = ZERO 242 ELSE 243 IN( K ) = 1 244 MULT = A( K ) / C( K ) 245 A( K ) = C( K ) 246 TEMP = A( K+1 ) 247 A( K+1 ) = B( K ) - MULT*TEMP 248 IF( K.LT.( N-1 ) ) THEN 249 D( K ) = B( K+1 ) 250 B( K+1 ) = -MULT*D( K ) 251 END IF 252 B( K ) = TEMP 253 C( K ) = MULT 254 END IF 255 END IF 256 IF( ( MAX( PIV1, PIV2 ).LE.TL ) .AND. ( IN( N ).EQ.0 ) ) 257 $ IN( N ) = K 258 10 CONTINUE 259 IF( ( ABS( A( N ) ).LE.SCALE1*TL ) .AND. ( IN( N ).EQ.0 ) ) 260 $ IN( N ) = N 261* 262 RETURN 263* 264* End of DLAGTF 265* 266 END 267