188 SUBROUTINE sstedc( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
198 INTEGER INFO, LDZ, LIWORK, LWORK, N
202 REAL D( * ), E( * ), WORK( * ), Z( ldz, * )
209 parameter( zero = 0.0e0, one = 1.0e0, two = 2.0e0 )
213 INTEGER FINISH, I, ICOMPZ, II, J, K, LGN, LIWMIN,
214 $ lwmin, m, smlsiz, start, storez, strtrw
215 REAL EPS, ORGNRM, P, TINY
221 EXTERNAL ilaenv, lsame, slamch, slanst
228 INTRINSIC abs, int, log, max, mod,
REAL, SQRT
235 lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 )
237 IF( lsame( compz,
'N' ) )
THEN
239 ELSE IF( lsame( compz,
'V' ) )
THEN
241 ELSE IF( lsame( compz,
'I' ) )
THEN
246 IF( icompz.LT.0 )
THEN
248 ELSE IF( n.LT.0 )
THEN
250 ELSE IF( ( ldz.LT.1 ) .OR.
251 $ ( icompz.GT.0 .AND. ldz.LT.max( 1, n ) ) )
THEN
259 smlsiz = ilaenv( 9,
'SSTEDC',
' ', 0, 0, 0, 0 )
260 IF( n.LE.1 .OR. icompz.EQ.0 )
THEN
263 ELSE IF( n.LE.smlsiz )
THEN
267 lgn = int( log(
REAL( N ) )/log( TWO ) )
272 IF( icompz.EQ.1 )
THEN
273 lwmin = 1 + 3*n + 2*n*lgn + 4*n**2
274 liwmin = 6 + 6*n + 5*n*lgn
275 ELSE IF( icompz.EQ.2 )
THEN
276 lwmin = 1 + 4*n + n**2
283 IF( lwork.LT.lwmin .AND. .NOT. lquery )
THEN
285 ELSE IF( liwork.LT.liwmin .AND. .NOT. lquery )
THEN
291 CALL xerbla(
'SSTEDC', -info )
293 ELSE IF (lquery)
THEN
318 IF( icompz.EQ.0 )
THEN
319 CALL ssterf( n, d, e, info )
326 IF( n.LE.smlsiz )
THEN
328 CALL ssteqr( compz, n, d, e, z, ldz, work, info )
335 IF( icompz.EQ.1 )
THEN
341 IF( icompz.EQ.2 )
THEN
342 CALL slaset(
'Full', n, n, zero, one, z, ldz )
347 orgnrm = slanst(
'M', n, d, e )
351 eps = slamch(
'Epsilon' )
358 IF( start.LE.n )
THEN
368 IF( finish.LT.n )
THEN
369 tiny = eps*sqrt( abs( d( finish ) ) )*
370 $ sqrt( abs( d( finish+1 ) ) )
371 IF( abs( e( finish ) ).GT.tiny )
THEN
379 m = finish - start + 1
384 IF( m.GT.smlsiz )
THEN
388 orgnrm = slanst(
'M', m, d( start ), e( start ) )
389 CALL slascl(
'G', 0, 0, orgnrm, one, m, 1, d( start ), m,
391 CALL slascl(
'G', 0, 0, orgnrm, one, m-1, 1, e( start ),
394 IF( icompz.EQ.1 )
THEN
399 CALL slaed0( icompz, n, m, d( start ), e( start ),
400 $ z( strtrw, start ), ldz, work( 1 ), n,
401 $ work( storez ), iwork, info )
403 info = ( info / ( m+1 )+start-1 )*( n+1 ) +
404 $ mod( info, ( m+1 ) ) + start - 1
410 CALL slascl(
'G', 0, 0, one, orgnrm, m, 1, d( start ), m,
414 IF( icompz.EQ.1 )
THEN
420 CALL ssteqr(
'I', m, d( start ), e( start ), work, m,
421 $ work( m*m+1 ), info )
422 CALL slacpy(
'A', n, m, z( 1, start ), ldz,
423 $ work( storez ), n )
424 CALL sgemm(
'N',
'N', n, m, m, one,
425 $ work( storez ), n, work, m, zero,
426 $ z( 1, start ), ldz )
427 ELSE IF( icompz.EQ.2 )
THEN
428 CALL ssteqr(
'I', m, d( start ), e( start ),
429 $ z( start, start ), ldz, work, info )
431 CALL ssterf( m, d( start ), e( start ), info )
434 info = start*( n+1 ) + finish
450 IF( icompz.EQ.0 )
THEN
454 CALL slasrt(
'I', n, d, info )
465 IF( d( j ).LT.p )
THEN
473 CALL sswap( n, z( 1, i ), 1, z( 1, k ), 1 )
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
subroutine sswap(N, SX, INCX, SY, INCY)
SSWAP
subroutine slaed0(ICOMPQ, QSIZ, N, D, E, Q, LDQ, QSTORE, LDQS, WORK, IWORK, INFO)
SLAED0 used by sstedc. Computes all eigenvalues and corresponding eigenvectors of an unreduced symmet...
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.
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
subroutine sstedc(COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK, LIWORK, INFO)
SSTEBZ
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
subroutine slasrt(ID, N, D, INFO)
SLASRT sorts numbers in increasing or decreasing order.
subroutine ssterf(N, D, E, INFO)
SSTERF
subroutine ssteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
SSTEQR