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 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 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.
real function slamch(CMACH)
SLAMCH
subroutine slasrt(ID, N, D, INFO)
SLASRT sorts numbers in increasing or decreasing order.
subroutine ssterf(N, D, E, INFO)
SSTERF
real function slapy2(X, Y)
SLAPY2 returns sqrt(x2+y2).