98 DOUBLE PRECISION d( * ), e( * )
104 DOUBLE PRECISION zero, one, two, three
105 parameter( zero = 0.0d0, one = 1.0d0, two = 2.0d0,
108 parameter( maxit = 30 )
111 INTEGER i, iscale, jtot, l, l1, lend, lendsv, lsv, m,
113 DOUBLE PRECISION alpha, anorm, bb, c, eps, eps2, gamma, oldc,
114 $ oldgam, p, r, rt1, rt2, rte, s, safmax, safmin,
115 $ sigma, ssfmax, ssfmin, rmax
125 INTRINSIC abs, sign, sqrt
137 CALL
xerbla(
'DSTERF', -info )
148 safmax = one / safmin
149 ssfmax = sqrt( safmax ) / three
150 ssfmin = sqrt( safmin ) / eps2
171 IF( abs( e( m ) ).LE.( sqrt( abs( d( m ) ) )*sqrt( abs( d( m+
172 $ 1 ) ) ) )*eps )
THEN
190 anorm =
dlanst(
'M', lend-l+1, d( l ), e( l ) )
194 IF( (anorm.GT.ssfmax) )
THEN
196 CALL
dlascl(
'G', 0, 0, anorm, ssfmax, lend-l+1, 1, d( l ), n,
198 CALL
dlascl(
'G', 0, 0, anorm, ssfmax, lend-l, 1, e( l ), n,
200 ELSE IF( anorm.LT.ssfmin )
THEN
202 CALL
dlascl(
'G', 0, 0, anorm, ssfmin, lend-l+1, 1, d( l ), n,
204 CALL
dlascl(
'G', 0, 0, anorm, ssfmin, lend-l, 1, e( l ), n,
208 DO 40 i = l, lend - 1
214 IF( abs( d( lend ) ).LT.abs( d( l ) ) )
THEN
227 DO 60 m = l, lend - 1
228 IF( abs( e( m ) ).LE.eps2*abs( d( m )*d( m+1 ) ) )
246 CALL
dlae2( d( l ), rte, d( l+1 ), rt1, rt2 )
263 sigma = ( d( l+1 )-p ) / ( two*rte )
265 sigma = p - ( rte / ( sigma+sign( r, sigma ) ) )
269 gamma = d( m ) - sigma
274 DO 80 i = m - 1, l, -1
284 gamma = c*( alpha-sigma ) - s*oldgam
285 d( i+1 ) = oldgam + ( alpha-gamma )
287 p = ( gamma*gamma ) / c
294 d( l ) = sigma + gamma
314 DO 110 m = l, lend + 1, -1
315 IF( abs( e( m-1 ) ).LE.eps2*abs( d( m )*d( m-1 ) ) )
331 rte = sqrt( e( l-1 ) )
332 CALL
dlae2( d( l ), rte, d( l-1 ), rt1, rt2 )
348 rte = sqrt( e( l-1 ) )
349 sigma = ( d( l-1 )-p ) / ( two*rte )
351 sigma = p - ( rte / ( sigma+sign( r, sigma ) ) )
355 gamma = d( m ) - sigma
370 gamma = c*( alpha-sigma ) - s*oldgam
371 d( i ) = oldgam + ( alpha-gamma )
373 p = ( gamma*gamma ) / c
380 d( l ) = sigma + gamma
399 $ CALL
dlascl(
'G', 0, 0, ssfmax, anorm, lendsv-lsv+1, 1,
400 $ d( lsv ), n, info )
402 $ CALL
dlascl(
'G', 0, 0, ssfmin, anorm, lendsv-lsv+1, 1,
403 $ d( lsv ), n, info )
419 CALL
dlasrt(
'I', n, d, info )
subroutine dlasrt(ID, N, D, INFO)
DLASRT sorts numbers in increasing or decreasing order.
subroutine dsterf(N, D, E, INFO)
DSTERF
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine dlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
subroutine dlae2(A, B, C, RT1, RT2)
DLAE2 computes the eigenvalues of a 2-by-2 symmetric matrix.
double precision function dlapy2(X, Y)
DLAPY2 returns sqrt(x2+y2).
double precision function dlanst(NORM, N, D, E)
DLANST returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a real symmetric tridiagonal matrix.
double precision function dlamch(CMACH)
DLAMCH