132 SUBROUTINE dsteqr( COMPZ, N, D, E, Z, LDZ, WORK, INFO )
144 DOUBLE PRECISION D( * ), E( * ), WORK( * ), Z( ldz, * )
150 DOUBLE PRECISION ZERO, ONE, TWO, THREE
151 parameter( zero = 0.0d0, one = 1.0d0, two = 2.0d0,
154 parameter( maxit = 30 )
157 INTEGER I, ICOMPZ, II, ISCALE, J, JTOT, K, L, L1, LEND,
158 $ lendm1, lendp1, lendsv, lm1, lsv, m, mm, mm1,
160 DOUBLE PRECISION ANORM, B, C, EPS, EPS2, F, G, P, R, RT1, RT2,
161 $ s, safmax, safmin, ssfmax, ssfmin, tst
165 DOUBLE PRECISION DLAMCH, DLANST, DLAPY2
166 EXTERNAL lsame, dlamch, dlanst, dlapy2
173 INTRINSIC abs, max, sign, sqrt
181 IF( lsame( compz,
'N' ) )
THEN
183 ELSE IF( lsame( compz,
'V' ) )
THEN
185 ELSE IF( lsame( compz,
'I' ) )
THEN
190 IF( icompz.LT.0 )
THEN
192 ELSE IF( n.LT.0 )
THEN
194 ELSE IF( ( ldz.LT.1 ) .OR. ( icompz.GT.0 .AND. ldz.LT.max( 1,
199 CALL
xerbla(
'DSTEQR', -info )
218 safmin = dlamch(
'S' )
219 safmax = one / safmin
220 ssfmax = sqrt( safmax ) / three
221 ssfmin = sqrt( safmin ) / eps2
227 $ CALL
dlaset(
'Full', n, n, zero, one, z, ldz )
249 IF( tst.LE.( sqrt( abs( d( m ) ) )*sqrt( abs( d( m+
250 $ 1 ) ) ) )*eps )
THEN
269 anorm = dlanst(
'M', lend-l+1, d( l ), e( l ) )
273 IF( anorm.GT.ssfmax )
THEN
275 CALL
dlascl(
'G', 0, 0, anorm, ssfmax, lend-l+1, 1, d( l ), n,
277 CALL
dlascl(
'G', 0, 0, anorm, ssfmax, lend-l, 1, e( l ), n,
279 ELSE IF( anorm.LT.ssfmin )
THEN
281 CALL
dlascl(
'G', 0, 0, anorm, ssfmin, lend-l+1, 1, d( l ), n,
283 CALL
dlascl(
'G', 0, 0, anorm, ssfmin, lend-l, 1, e( l ), n,
289 IF( abs( d( lend ) ).LT.abs( d( l ) ) )
THEN
304 tst = abs( e( m ) )**2
305 IF( tst.LE.( eps2*abs( d( m ) ) )*abs( d( m+1 ) )+
323 IF( icompz.GT.0 )
THEN
324 CALL
dlaev2( d( l ), e( l ), d( l+1 ), rt1, rt2, c, s )
327 CALL
dlasr(
'R',
'V',
'B', n, 2, work( l ),
328 $ work( n-1+l ), z( 1, l ), ldz )
330 CALL
dlae2( d( l ), e( l ), d( l+1 ), rt1, rt2 )
347 g = ( d( l+1 )-p ) / ( two*e( l ) )
349 g = d( m ) - p + ( e( l ) / ( g+sign( r, g ) ) )
361 CALL
dlartg( g, f, c, s, r )
365 r = ( d( i )-g )*s + two*c*b
372 IF( icompz.GT.0 )
THEN
381 IF( icompz.GT.0 )
THEN
383 CALL
dlasr(
'R',
'V',
'B', n, mm, work( l ), work( n-1+l ),
410 DO 100 m = l, lendp1, -1
411 tst = abs( e( m-1 ) )**2
412 IF( tst.LE.( eps2*abs( d( m ) ) )*abs( d( m-1 ) )+
430 IF( icompz.GT.0 )
THEN
431 CALL
dlaev2( d( l-1 ), e( l-1 ), d( l ), rt1, rt2, c, s )
434 CALL
dlasr(
'R',
'V',
'F', n, 2, work( m ),
435 $ work( n-1+m ), z( 1, l-1 ), ldz )
437 CALL
dlae2( d( l-1 ), e( l-1 ), d( l ), rt1, rt2 )
454 g = ( d( l-1 )-p ) / ( two*e( l-1 ) )
456 g = d( m ) - p + ( e( l-1 ) / ( g+sign( r, g ) ) )
468 CALL
dlartg( g, f, c, s, r )
472 r = ( d( i+1 )-g )*s + two*c*b
479 IF( icompz.GT.0 )
THEN
488 IF( icompz.GT.0 )
THEN
490 CALL
dlasr(
'R',
'V',
'F', n, mm, work( m ), work( n-1+m ),
513 IF( iscale.EQ.1 )
THEN
514 CALL
dlascl(
'G', 0, 0, ssfmax, anorm, lendsv-lsv+1, 1,
515 $ d( lsv ), n, info )
516 CALL
dlascl(
'G', 0, 0, ssfmax, anorm, lendsv-lsv, 1, e( lsv ),
518 ELSE IF( iscale.EQ.2 )
THEN
519 CALL
dlascl(
'G', 0, 0, ssfmin, anorm, lendsv-lsv+1, 1,
520 $ d( lsv ), n, info )
521 CALL
dlascl(
'G', 0, 0, ssfmin, anorm, lendsv-lsv, 1, e( lsv ),
539 IF( icompz.EQ.0 )
THEN
543 CALL
dlasrt(
'I', n, d, info )
554 IF( d( j ).LT.p )
THEN
562 CALL
dswap( n, z( 1, i ), 1, z( 1, k ), 1 )
subroutine dsteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
DSTEQR
subroutine dlasrt(ID, N, D, INFO)
DLASRT sorts numbers in increasing or decreasing order.
subroutine dlartg(F, G, CS, SN, R)
DLARTG generates a plane rotation with real cosine and real sine.
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.
subroutine dlasr(SIDE, PIVOT, DIRECT, M, N, C, S, A, LDA)
DLASR applies a sequence of plane rotations to a general rectangular matrix.
subroutine dlaev2(A, B, C, RT1, RT2, CS1, SN1)
DLAEV2 computes the eigenvalues and eigenvectors of a 2-by-2 symmetric/Hermitian matrix.
subroutine dswap(N, DX, INCX, DY, INCY)
DSWAP
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...