104 REAL zero, one, two, three
105 parameter ( zero = 0.0e0, one = 1.0e0, two = 2.0e0,
108 parameter ( maxit = 30 )
111 INTEGER i, iscale, jtot, l, l1, lend, lendsv, lsv, m,
113 REAL alpha, anorm, bb, c, eps, eps2, gamma, oldc,
114 $ oldgam, p, r, rt1, rt2, rte, s, safmax, safmin,
115 $ sigma, ssfmax, ssfmin
125 INTRINSIC abs, sign, sqrt
137 CALL xerbla(
'SSTERF', -info )
148 safmax = one / safmin
149 ssfmax = sqrt( safmax ) / three
150 ssfmin = sqrt( safmin ) / eps2
170 IF( abs( e( m ) ).LE.( sqrt( abs( d( m ) ) )*
171 $ sqrt( abs( d( m+1 ) ) ) )*eps )
THEN
189 anorm =
slanst(
'M', lend-l+1, d( l ), e( l ) )
193 IF( anorm.GT.ssfmax )
THEN
195 CALL slascl(
'G', 0, 0, anorm, ssfmax, lend-l+1, 1, d( l ), n,
197 CALL slascl(
'G', 0, 0, anorm, ssfmax, lend-l, 1, e( l ), n,
199 ELSE IF( anorm.LT.ssfmin )
THEN
201 CALL slascl(
'G', 0, 0, anorm, ssfmin, lend-l+1, 1, d( l ), n,
203 CALL slascl(
'G', 0, 0, anorm, ssfmin, lend-l, 1, e( l ), n,
207 DO 40 i = l, lend - 1
213 IF( abs( d( lend ) ).LT.abs( d( l ) ) )
THEN
226 DO 60 m = l, lend - 1
227 IF( abs( e( m ) ).LE.eps2*abs( d( m )*d( m+1 ) ) )
245 CALL slae2( d( l ), rte, d( l+1 ), rt1, rt2 )
262 sigma = ( d( l+1 )-p ) / ( two*rte )
264 sigma = p - ( rte / ( sigma+sign( r, sigma ) ) )
268 gamma = d( m ) - sigma
273 DO 80 i = m - 1, l, -1
283 gamma = c*( alpha-sigma ) - s*oldgam
284 d( i+1 ) = oldgam + ( alpha-gamma )
286 p = ( gamma*gamma ) / c
293 d( l ) = sigma + gamma
313 DO 110 m = l, lend + 1, -1
314 IF( abs( e( m-1 ) ).LE.eps2*abs( d( m )*d( m-1 ) ) )
330 rte = sqrt( e( l-1 ) )
331 CALL slae2( d( l ), rte, d( l-1 ), rt1, rt2 )
347 rte = sqrt( e( l-1 ) )
348 sigma = ( d( l-1 )-p ) / ( two*rte )
350 sigma = p - ( rte / ( sigma+sign( r, sigma ) ) )
354 gamma = d( m ) - sigma
369 gamma = c*( alpha-sigma ) - s*oldgam
370 d( i ) = oldgam + ( alpha-gamma )
372 p = ( gamma*gamma ) / c
379 d( l ) = sigma + gamma
398 $
CALL slascl(
'G', 0, 0, ssfmax, anorm, lendsv-lsv+1, 1,
399 $ d( lsv ), n, info )
401 $
CALL slascl(
'G', 0, 0, ssfmin, anorm, lendsv-lsv+1, 1,
402 $ d( lsv ), n, info )
418 CALL slasrt(
'I', n, d, info )
subroutine slae2(A, B, C, RT1, RT2)
SLAE2 computes the eigenvalues of a 2-by-2 symmetric matrix.
subroutine slasrt(ID, N, D, INFO)
SLASRT sorts numbers in increasing or decreasing order.
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine slascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
real function slamch(CMACH)
SLAMCH
real function slapy2(X, Y)
SLAPY2 returns sqrt(x2+y2).
real function slanst(NORM, N, D, E)
SLANST 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.