1 SUBROUTINE SLAGTF( N, A, LAMBDA, B, C, TOL, D, IN, INFO ) 2* 3* -- LAPACK routine (version 3.0) -- 4* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., 5* Courant Institute, Argonne National Lab, and Rice University 6* June 30, 1999 7* 8* .. Scalar Arguments .. 9 INTEGER INFO, N 10 REAL LAMBDA, TOL 11* .. 12* .. Array Arguments .. 13 INTEGER IN( * ) 14 REAL A( * ), B( * ), C( * ), D( * ) 15* .. 16* 17* Purpose 18* ======= 19* 20* SLAGTF factorizes the matrix (T - lambda*I), where T is an n by n 21* tridiagonal matrix and lambda is a scalar, as 22* 23* T - lambda*I = PLU, 24* 25* where P is a permutation matrix, L is a unit lower tridiagonal matrix 26* with at most one non-zero sub-diagonal elements per column and U is 27* an upper triangular matrix with at most two non-zero super-diagonal 28* elements per column. 29* 30* The factorization is obtained by Gaussian elimination with partial 31* pivoting and implicit row scaling. 32* 33* The parameter LAMBDA is included in the routine so that SLAGTF may 34* be used, in conjunction with SLAGTS, to obtain eigenvectors of T by 35* inverse iteration. 36* 37* Arguments 38* ========= 39* 40* N (input) INTEGER 41* The order of the matrix T. 42* 43* A (input/output) REAL array, dimension (N) 44* On entry, A must contain the diagonal elements of T. 45* 46* On exit, A is overwritten by the n diagonal elements of the 47* upper triangular matrix U of the factorization of T. 48* 49* LAMBDA (input) REAL 50* On entry, the scalar lambda. 51* 52* B (input/output) REAL array, dimension (N-1) 53* On entry, B must contain the (n-1) super-diagonal elements of 54* T. 55* 56* On exit, B is overwritten by the (n-1) super-diagonal 57* elements of the matrix U of the factorization of T. 58* 59* C (input/output) REAL array, dimension (N-1) 60* On entry, C must contain the (n-1) sub-diagonal elements of 61* T. 62* 63* On exit, C is overwritten by the (n-1) sub-diagonal elements 64* of the matrix L of the factorization of T. 65* 66* TOL (input) REAL 67* On entry, a relative tolerance used to indicate whether or 68* not the matrix (T - lambda*I) is nearly singular. TOL should 69* normally be chose as approximately the largest relative error 70* in the elements of T. For example, if the elements of T are 71* correct to about 4 significant figures, then TOL should be 72* set to about 5*10**(-4). If TOL is supplied as less than eps, 73* where eps is the relative machine precision, then the value 74* eps is used in place of TOL. 75* 76* D (output) REAL array, dimension (N-2) 77* On exit, D is overwritten by the (n-2) second super-diagonal 78* elements of the matrix U of the factorization of T. 79* 80* IN (output) INTEGER array, dimension (N) 81* On exit, IN contains details of the permutation matrix P. If 82* an interchange occurred at the kth step of the elimination, 83* then IN(k) = 1, otherwise IN(k) = 0. The element IN(n) 84* returns the smallest positive integer j such that 85* 86* abs( u(j,j) ).le. norm( (T - lambda*I)(j) )*TOL, 87* 88* where norm( A(j) ) denotes the sum of the absolute values of 89* the jth row of the matrix A. If no such j exists then IN(n) 90* is returned as zero. If IN(n) is returned as positive, then a 91* diagonal element of U is small, indicating that 92* (T - lambda*I) is singular or nearly singular, 93* 94* INFO (output) INTEGER 95* = 0 : successful exit 96* .lt. 0: if INFO = -k, the kth argument had an illegal value 97* 98* ===================================================================== 99* 100* .. Parameters .. 101 REAL ZERO 102 PARAMETER ( ZERO = 0.0E+0 ) 103* .. 104* .. Local Scalars .. 105 INTEGER K 106 REAL EPS, MULT, PIV1, PIV2, SCALE1, SCALE2, TEMP, TL 107* .. 108* .. Intrinsic Functions .. 109 INTRINSIC ABS, MAX 110* .. 111* .. External Functions .. 112 REAL SLAMCH 113 EXTERNAL SLAMCH 114* .. 115* .. External Subroutines .. 116 EXTERNAL XERBLA 117* .. 118* .. Executable Statements .. 119* 120 INFO = 0 121 IF( N.LT.0 ) THEN 122 INFO = -1 123 CALL XERBLA( 'SLAGTF', -INFO ) 124 RETURN 125 END IF 126* 127 IF( N.EQ.0 ) 128 $ RETURN 129* 130 A( 1 ) = A( 1 ) - LAMBDA 131 IN( N ) = 0 132 IF( N.EQ.1 ) THEN 133 IF( A( 1 ).EQ.ZERO ) 134 $ IN( 1 ) = 1 135 RETURN 136 END IF 137* 138 EPS = SLAMCH( 'Epsilon' ) 139* 140 TL = MAX( TOL, EPS ) 141 SCALE1 = ABS( A( 1 ) ) + ABS( B( 1 ) ) 142 DO 10 K = 1, N - 1 143 A( K+1 ) = A( K+1 ) - LAMBDA 144 SCALE2 = ABS( C( K ) ) + ABS( A( K+1 ) ) 145 IF( K.LT.( N-1 ) ) 146 $ SCALE2 = SCALE2 + ABS( B( K+1 ) ) 147 IF( A( K ).EQ.ZERO ) THEN 148 PIV1 = ZERO 149 ELSE 150 PIV1 = ABS( A( K ) ) / SCALE1 151 END IF 152 IF( C( K ).EQ.ZERO ) THEN 153 IN( K ) = 0 154 PIV2 = ZERO 155 SCALE1 = SCALE2 156 IF( K.LT.( N-1 ) ) 157 $ D( K ) = ZERO 158 ELSE 159 PIV2 = ABS( C( K ) ) / SCALE2 160 IF( PIV2.LE.PIV1 ) THEN 161 IN( K ) = 0 162 SCALE1 = SCALE2 163 C( K ) = C( K ) / A( K ) 164 A( K+1 ) = A( K+1 ) - C( K )*B( K ) 165 IF( K.LT.( N-1 ) ) 166 $ D( K ) = ZERO 167 ELSE 168 IN( K ) = 1 169 MULT = A( K ) / C( K ) 170 A( K ) = C( K ) 171 TEMP = A( K+1 ) 172 A( K+1 ) = B( K ) - MULT*TEMP 173 IF( K.LT.( N-1 ) ) THEN 174 D( K ) = B( K+1 ) 175 B( K+1 ) = -MULT*D( K ) 176 END IF 177 B( K ) = TEMP 178 C( K ) = MULT 179 END IF 180 END IF 181 IF( ( MAX( PIV1, PIV2 ).LE.TL ) .AND. ( IN( N ).EQ.0 ) ) 182 $ IN( N ) = K 183 10 CONTINUE 184 IF( ( ABS( A( N ) ).LE.SCALE1*TL ) .AND. ( IN( N ).EQ.0 ) ) 185 $ IN( N ) = N 186* 187 RETURN 188* 189* End of SLAGTF 190* 191 END 192