259 SUBROUTINE zhbevx( JOBZ, RANGE, UPLO, N, KD, AB, LDAB, Q, LDQ, VL,
260 $ vu, il, iu, abstol, m, w, z, ldz, work, rwork,
261 $ iwork, ifail, info )
269 CHARACTER JOBZ, RANGE, UPLO
270 INTEGER IL, INFO, IU, KD, LDAB, LDQ, LDZ, M, N
271 DOUBLE PRECISION ABSTOL, VL, VU
274 INTEGER IFAIL( * ), IWORK( * )
275 DOUBLE PRECISION RWORK( * ), W( * )
276 COMPLEX*16 AB( ldab, * ), Q( ldq, * ), WORK( * ),
283 DOUBLE PRECISION ZERO, ONE
284 parameter( zero = 0.0d0, one = 1.0d0 )
285 COMPLEX*16 CZERO, CONE
286 parameter( czero = ( 0.0d0, 0.0d0 ),
287 $ cone = ( 1.0d0, 0.0d0 ) )
290 LOGICAL ALLEIG, INDEIG, LOWER, TEST, VALEIG, WANTZ
292 INTEGER I, IINFO, IMAX, INDD, INDE, INDEE, INDIBL,
293 $ indisp, indiwk, indrwk, indwrk, iscale, itmp1,
295 DOUBLE PRECISION ABSTLL, ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN,
296 $ sigma, smlnum, tmp1, vll, vuu
301 DOUBLE PRECISION DLAMCH, ZLANHB
302 EXTERNAL lsame, dlamch, zlanhb
310 INTRINSIC dble, max, min, sqrt
316 wantz = lsame( jobz,
'V' )
317 alleig = lsame( range,
'A' )
318 valeig = lsame( range,
'V' )
319 indeig = lsame( range,
'I' )
320 lower = lsame( uplo,
'L' )
323 IF( .NOT.( wantz .OR. lsame( jobz,
'N' ) ) )
THEN
325 ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) )
THEN
327 ELSE IF( .NOT.( lower .OR. lsame( uplo,
'U' ) ) )
THEN
329 ELSE IF( n.LT.0 )
THEN
331 ELSE IF( kd.LT.0 )
THEN
333 ELSE IF( ldab.LT.kd+1 )
THEN
335 ELSE IF( wantz .AND. ldq.LT.max( 1, n ) )
THEN
339 IF( n.GT.0 .AND. vu.LE.vl )
341 ELSE IF( indeig )
THEN
342 IF( il.LT.1 .OR. il.GT.max( 1, n ) )
THEN
344 ELSE IF( iu.LT.min( n, il ) .OR. iu.GT.n )
THEN
350 IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) )
355 CALL
xerbla(
'ZHBEVX', -info )
370 ctmp1 = ab( kd+1, 1 )
374 IF( .NOT.( vl.LT.tmp1 .AND. vu.GE.tmp1 ) )
387 safmin = dlamch(
'Safe minimum' )
388 eps = dlamch(
'Precision' )
389 smlnum = safmin / eps
390 bignum = one / smlnum
391 rmin = sqrt( smlnum )
392 rmax = min( sqrt( bignum ), one / sqrt( sqrt( safmin ) ) )
405 anrm = zlanhb(
'M', uplo, n, kd, ab, ldab, rwork )
406 IF( anrm.GT.zero .AND. anrm.LT.rmin )
THEN
409 ELSE IF( anrm.GT.rmax )
THEN
413 IF( iscale.EQ.1 )
THEN
415 CALL
zlascl(
'B', kd, kd, one, sigma, n, n, ab, ldab, info )
417 CALL
zlascl(
'Q', kd, kd, one, sigma, n, n, ab, ldab, info )
420 $ abstll = abstol*sigma
433 CALL
zhbtrd( jobz, uplo, n, kd, ab, ldab, rwork( indd ),
434 $ rwork( inde ), q, ldq, work( indwrk ), iinfo )
442 IF (il.EQ.1 .AND. iu.EQ.n)
THEN
446 IF ((alleig .OR. test) .AND. (abstol.LE.zero))
THEN
447 CALL
dcopy( n, rwork( indd ), 1, w, 1 )
449 IF( .NOT.wantz )
THEN
450 CALL
dcopy( n-1, rwork( inde ), 1, rwork( indee ), 1 )
451 CALL
dsterf( n, w, rwork( indee ), info )
453 CALL
zlacpy(
'A', n, n, q, ldq, z, ldz )
454 CALL
dcopy( n-1, rwork( inde ), 1, rwork( indee ), 1 )
455 CALL
zsteqr( jobz, n, w, rwork( indee ), z, ldz,
456 $ rwork( indrwk ), info )
480 CALL
dstebz( range, order, n, vll, vuu, il, iu, abstll,
481 $ rwork( indd ), rwork( inde ), m, nsplit, w,
482 $ iwork( indibl ), iwork( indisp ), rwork( indrwk ),
483 $ iwork( indiwk ), info )
486 CALL
zstein( n, rwork( indd ), rwork( inde ), m, w,
487 $ iwork( indibl ), iwork( indisp ), z, ldz,
488 $ rwork( indrwk ), iwork( indiwk ), ifail, info )
494 CALL
zcopy( n, z( 1, j ), 1, work( 1 ), 1 )
495 CALL
zgemv(
'N', n, n, cone, q, ldq, work, 1, czero,
503 IF( iscale.EQ.1 )
THEN
509 CALL
dscal( imax, one / sigma, w, 1 )
520 IF( w( jj ).LT.tmp1 )
THEN
527 itmp1 = iwork( indibl+i-1 )
529 iwork( indibl+i-1 ) = iwork( indibl+j-1 )
531 iwork( indibl+j-1 ) = itmp1
532 CALL zswap( n, z( 1, i ), 1, z( 1, j ), 1 )
535 ifail( i ) = ifail( j )
subroutine zgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
ZGEMV
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
subroutine zstein(N, D, E, M, W, IBLOCK, ISPLIT, Z, LDZ, WORK, IWORK, IFAIL, INFO)
ZSTEIN
subroutine dsterf(N, D, E, INFO)
DSTERF
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine zhbtrd(VECT, UPLO, N, KD, AB, LDAB, D, E, Q, LDQ, WORK, INFO)
ZHBTRD
subroutine zlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
ZLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
subroutine dscal(N, DA, DX, INCX)
DSCAL
subroutine zsteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
ZSTEQR
subroutine zhbevx(JOBZ, RANGE, UPLO, N, KD, AB, LDAB, Q, LDQ, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, WORK, RWORK, IWORK, IFAIL, INFO)
ZHBEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matric...
subroutine dstebz(RANGE, ORDER, N, VL, VU, IL, IU, ABSTOL, D, E, M, NSPLIT, W, IBLOCK, ISPLIT, WORK, IWORK, INFO)
DSTEBZ