Annotation of rpl/lapack/lapack/dlagtf.f, revision 1.1
1.1 ! bertrand 1: SUBROUTINE DLAGTF( N, A, LAMBDA, B, C, TOL, D, IN, INFO )
! 2: *
! 3: * -- LAPACK routine (version 3.2) --
! 4: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 6: * November 2006
! 7: *
! 8: * .. Scalar Arguments ..
! 9: INTEGER INFO, N
! 10: DOUBLE PRECISION LAMBDA, TOL
! 11: * ..
! 12: * .. Array Arguments ..
! 13: INTEGER IN( * )
! 14: DOUBLE PRECISION A( * ), B( * ), C( * ), D( * )
! 15: * ..
! 16: *
! 17: * Purpose
! 18: * =======
! 19: *
! 20: * DLAGTF 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 DLAGTF may
! 34: * be used, in conjunction with DLAGTS, 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) DOUBLE PRECISION 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) DOUBLE PRECISION
! 50: * On entry, the scalar lambda.
! 51: *
! 52: * B (input/output) DOUBLE PRECISION 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) DOUBLE PRECISION 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) DOUBLE PRECISION
! 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) DOUBLE PRECISION 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: DOUBLE PRECISION ZERO
! 102: PARAMETER ( ZERO = 0.0D+0 )
! 103: * ..
! 104: * .. Local Scalars ..
! 105: INTEGER K
! 106: DOUBLE PRECISION EPS, MULT, PIV1, PIV2, SCALE1, SCALE2, TEMP, TL
! 107: * ..
! 108: * .. Intrinsic Functions ..
! 109: INTRINSIC ABS, MAX
! 110: * ..
! 111: * .. External Functions ..
! 112: DOUBLE PRECISION DLAMCH
! 113: EXTERNAL DLAMCH
! 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( 'DLAGTF', -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 = DLAMCH( '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 DLAGTF
! 190: *
! 191: END
CVSweb interface <joel.bertrand@systella.fr>