LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine chbevx ( character  JOBZ,
character  RANGE,
character  UPLO,
integer  N,
integer  KD,
complex, dimension( ldab, * )  AB,
integer  LDAB,
complex, dimension( ldq, * )  Q,
integer  LDQ,
real  VL,
real  VU,
integer  IL,
integer  IU,
real  ABSTOL,
integer  M,
real, dimension( * )  W,
complex, dimension( ldz, * )  Z,
integer  LDZ,
complex, dimension( * )  WORK,
real, dimension( * )  RWORK,
integer, dimension( * )  IWORK,
integer, dimension( * )  IFAIL,
integer  INFO 
)

CHBEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices

Download CHBEVX + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 CHBEVX computes selected eigenvalues and, optionally, eigenvectors
 of a complex Hermitian band matrix A.  Eigenvalues and eigenvectors
 can be selected by specifying either a range of values or a range of
 indices for the desired eigenvalues.
Parameters
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]RANGE
          RANGE is CHARACTER*1
          = 'A': all eigenvalues will be found;
          = 'V': all eigenvalues in the half-open interval (VL,VU]
                 will be found;
          = 'I': the IL-th through IU-th eigenvalues will be found.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in]KD
          KD is INTEGER
          The number of superdiagonals of the matrix A if UPLO = 'U',
          or the number of subdiagonals if UPLO = 'L'.  KD >= 0.
[in,out]AB
          AB is COMPLEX array, dimension (LDAB, N)
          On entry, the upper or lower triangle of the Hermitian band
          matrix A, stored in the first KD+1 rows of the array.  The
          j-th column of A is stored in the j-th column of the array AB
          as follows:
          if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(n,j+kd).

          On exit, AB is overwritten by values generated during the
          reduction to tridiagonal form.
[in]LDAB
          LDAB is INTEGER
          The leading dimension of the array AB.  LDAB >= KD + 1.
[out]Q
          Q is COMPLEX array, dimension (LDQ, N)
          If JOBZ = 'V', the N-by-N unitary matrix used in the
                          reduction to tridiagonal form.
          If JOBZ = 'N', the array Q is not referenced.
[in]LDQ
          LDQ is INTEGER
          The leading dimension of the array Q.  If JOBZ = 'V', then
          LDQ >= max(1,N).
[in]VL
          VL is REAL
          If RANGE='V', the lower bound of the interval to
          be searched for eigenvalues. VL < VU.
          Not referenced if RANGE = 'A' or 'I'.
[in]VU
          VU is REAL
          If RANGE='V', the upper bound of the interval to
          be searched for eigenvalues. VL < VU.
          Not referenced if RANGE = 'A' or 'I'.
[in]IL
          IL is INTEGER
          If RANGE='I', the index of the
          smallest eigenvalue to be returned.
          1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
          Not referenced if RANGE = 'A' or 'V'.
[in]IU
          IU is INTEGER
          If RANGE='I', the index of the
          largest eigenvalue to be returned.
          1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
          Not referenced if RANGE = 'A' or 'V'.
[in]ABSTOL
          ABSTOL is REAL
          The absolute error tolerance for the eigenvalues.
          An approximate eigenvalue is accepted as converged
          when it is determined to lie in an interval [a,b]
          of width less than or equal to

                  ABSTOL + EPS *   max( |a|,|b| ) ,

          where EPS is the machine precision.  If ABSTOL is less than
          or equal to zero, then  EPS*|T|  will be used in its place,
          where |T| is the 1-norm of the tridiagonal matrix obtained
          by reducing AB to tridiagonal form.

          Eigenvalues will be computed most accurately when ABSTOL is
          set to twice the underflow threshold 2*SLAMCH('S'), not zero.
          If this routine returns with INFO>0, indicating that some
          eigenvectors did not converge, try setting ABSTOL to
          2*SLAMCH('S').

          See "Computing Small Singular Values of Bidiagonal Matrices
          with Guaranteed High Relative Accuracy," by Demmel and
          Kahan, LAPACK Working Note #3.
[out]M
          M is INTEGER
          The total number of eigenvalues found.  0 <= M <= N.
          If RANGE = 'A', M = N, and if RANGE = 'I', M = IU-IL+1.
[out]W
          W is REAL array, dimension (N)
          The first M elements contain the selected eigenvalues in
          ascending order.
[out]Z
          Z is COMPLEX array, dimension (LDZ, max(1,M))
          If JOBZ = 'V', then if INFO = 0, the first M columns of Z
          contain the orthonormal eigenvectors of the matrix A
          corresponding to the selected eigenvalues, with the i-th
          column of Z holding the eigenvector associated with W(i).
          If an eigenvector fails to converge, then that column of Z
          contains the latest approximation to the eigenvector, and the
          index of the eigenvector is returned in IFAIL.
          If JOBZ = 'N', then Z is not referenced.
          Note: the user must ensure that at least max(1,M) columns are
          supplied in the array Z; if RANGE = 'V', the exact value of M
          is not known in advance and an upper bound must be used.
[in]LDZ
          LDZ is INTEGER
          The leading dimension of the array Z.  LDZ >= 1, and if
          JOBZ = 'V', LDZ >= max(1,N).
[out]WORK
          WORK is COMPLEX array, dimension (N)
[out]RWORK
          RWORK is REAL array, dimension (7*N)
[out]IWORK
          IWORK is INTEGER array, dimension (5*N)
[out]IFAIL
          IFAIL is INTEGER array, dimension (N)
          If JOBZ = 'V', then if INFO = 0, the first M elements of
          IFAIL are zero.  If INFO > 0, then IFAIL contains the
          indices of the eigenvectors that failed to converge.
          If JOBZ = 'N', then IFAIL is not referenced.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, then i eigenvectors failed to converge.
                Their indices are stored in array IFAIL.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
June 2016

Definition at line 269 of file chbevx.f.

269 *
270 * -- LAPACK driver routine (version 3.6.1) --
271 * -- LAPACK is a software package provided by Univ. of Tennessee, --
272 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
273 * June 2016
274 *
275 * .. Scalar Arguments ..
276  CHARACTER jobz, range, uplo
277  INTEGER il, info, iu, kd, ldab, ldq, ldz, m, n
278  REAL abstol, vl, vu
279 * ..
280 * .. Array Arguments ..
281  INTEGER ifail( * ), iwork( * )
282  REAL rwork( * ), w( * )
283  COMPLEX ab( ldab, * ), q( ldq, * ), work( * ),
284  $ z( ldz, * )
285 * ..
286 *
287 * =====================================================================
288 *
289 * .. Parameters ..
290  REAL zero, one
291  parameter ( zero = 0.0e0, one = 1.0e0 )
292  COMPLEX czero, cone
293  parameter ( czero = ( 0.0e0, 0.0e0 ),
294  $ cone = ( 1.0e0, 0.0e0 ) )
295 * ..
296 * .. Local Scalars ..
297  LOGICAL alleig, indeig, lower, test, valeig, wantz
298  CHARACTER order
299  INTEGER i, iinfo, imax, indd, inde, indee, indibl,
300  $ indisp, indiwk, indrwk, indwrk, iscale, itmp1,
301  $ j, jj, nsplit
302  REAL abstll, anrm, bignum, eps, rmax, rmin, safmin,
303  $ sigma, smlnum, tmp1, vll, vuu
304  COMPLEX ctmp1
305 * ..
306 * .. External Functions ..
307  LOGICAL lsame
308  REAL clanhb, slamch
309  EXTERNAL lsame, clanhb, slamch
310 * ..
311 * .. External Subroutines ..
312  EXTERNAL ccopy, cgemv, chbtrd, clacpy, clascl, cstein,
314  $ xerbla
315 * ..
316 * .. Intrinsic Functions ..
317  INTRINSIC max, min, REAL, sqrt
318 * ..
319 * .. Executable Statements ..
320 *
321 * Test the input parameters.
322 *
323  wantz = lsame( jobz, 'V' )
324  alleig = lsame( range, 'A' )
325  valeig = lsame( range, 'V' )
326  indeig = lsame( range, 'I' )
327  lower = lsame( uplo, 'L' )
328 *
329  info = 0
330  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
331  info = -1
332  ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) ) THEN
333  info = -2
334  ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
335  info = -3
336  ELSE IF( n.LT.0 ) THEN
337  info = -4
338  ELSE IF( kd.LT.0 ) THEN
339  info = -5
340  ELSE IF( ldab.LT.kd+1 ) THEN
341  info = -7
342  ELSE IF( wantz .AND. ldq.LT.max( 1, n ) ) THEN
343  info = -9
344  ELSE
345  IF( valeig ) THEN
346  IF( n.GT.0 .AND. vu.LE.vl )
347  $ info = -11
348  ELSE IF( indeig ) THEN
349  IF( il.LT.1 .OR. il.GT.max( 1, n ) ) THEN
350  info = -12
351  ELSE IF( iu.LT.min( n, il ) .OR. iu.GT.n ) THEN
352  info = -13
353  END IF
354  END IF
355  END IF
356  IF( info.EQ.0 ) THEN
357  IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) )
358  $ info = -18
359  END IF
360 *
361  IF( info.NE.0 ) THEN
362  CALL xerbla( 'CHBEVX', -info )
363  RETURN
364  END IF
365 *
366 * Quick return if possible
367 *
368  m = 0
369  IF( n.EQ.0 )
370  $ RETURN
371 *
372  IF( n.EQ.1 ) THEN
373  m = 1
374  IF( lower ) THEN
375  ctmp1 = ab( 1, 1 )
376  ELSE
377  ctmp1 = ab( kd+1, 1 )
378  END IF
379  tmp1 = REAL( ctmp1 )
380  IF( valeig ) THEN
381  IF( .NOT.( vl.LT.tmp1 .AND. vu.GE.tmp1 ) )
382  $ m = 0
383  END IF
384  IF( m.EQ.1 ) THEN
385  w( 1 ) = ctmp1
386  IF( wantz )
387  $ z( 1, 1 ) = cone
388  END IF
389  RETURN
390  END IF
391 *
392 * Get machine constants.
393 *
394  safmin = slamch( 'Safe minimum' )
395  eps = slamch( 'Precision' )
396  smlnum = safmin / eps
397  bignum = one / smlnum
398  rmin = sqrt( smlnum )
399  rmax = min( sqrt( bignum ), one / sqrt( sqrt( safmin ) ) )
400 *
401 * Scale matrix to allowable range, if necessary.
402 *
403  iscale = 0
404  abstll = abstol
405  IF ( valeig ) THEN
406  vll = vl
407  vuu = vu
408  ELSE
409  vll = zero
410  vuu = zero
411  ENDIF
412  anrm = clanhb( 'M', uplo, n, kd, ab, ldab, rwork )
413  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
414  iscale = 1
415  sigma = rmin / anrm
416  ELSE IF( anrm.GT.rmax ) THEN
417  iscale = 1
418  sigma = rmax / anrm
419  END IF
420  IF( iscale.EQ.1 ) THEN
421  IF( lower ) THEN
422  CALL clascl( 'B', kd, kd, one, sigma, n, n, ab, ldab, info )
423  ELSE
424  CALL clascl( 'Q', kd, kd, one, sigma, n, n, ab, ldab, info )
425  END IF
426  IF( abstol.GT.0 )
427  $ abstll = abstol*sigma
428  IF( valeig ) THEN
429  vll = vl*sigma
430  vuu = vu*sigma
431  END IF
432  END IF
433 *
434 * Call CHBTRD to reduce Hermitian band matrix to tridiagonal form.
435 *
436  indd = 1
437  inde = indd + n
438  indrwk = inde + n
439  indwrk = 1
440  CALL chbtrd( jobz, uplo, n, kd, ab, ldab, rwork( indd ),
441  $ rwork( inde ), q, ldq, work( indwrk ), iinfo )
442 *
443 * If all eigenvalues are desired and ABSTOL is less than or equal
444 * to zero, then call SSTERF or CSTEQR. If this fails for some
445 * eigenvalue, then try SSTEBZ.
446 *
447  test = .false.
448  IF (indeig) THEN
449  IF (il.EQ.1 .AND. iu.EQ.n) THEN
450  test = .true.
451  END IF
452  END IF
453  IF ((alleig .OR. test) .AND. (abstol.LE.zero)) THEN
454  CALL scopy( n, rwork( indd ), 1, w, 1 )
455  indee = indrwk + 2*n
456  IF( .NOT.wantz ) THEN
457  CALL scopy( n-1, rwork( inde ), 1, rwork( indee ), 1 )
458  CALL ssterf( n, w, rwork( indee ), info )
459  ELSE
460  CALL clacpy( 'A', n, n, q, ldq, z, ldz )
461  CALL scopy( n-1, rwork( inde ), 1, rwork( indee ), 1 )
462  CALL csteqr( jobz, n, w, rwork( indee ), z, ldz,
463  $ rwork( indrwk ), info )
464  IF( info.EQ.0 ) THEN
465  DO 10 i = 1, n
466  ifail( i ) = 0
467  10 CONTINUE
468  END IF
469  END IF
470  IF( info.EQ.0 ) THEN
471  m = n
472  GO TO 30
473  END IF
474  info = 0
475  END IF
476 *
477 * Otherwise, call SSTEBZ and, if eigenvectors are desired, CSTEIN.
478 *
479  IF( wantz ) THEN
480  order = 'B'
481  ELSE
482  order = 'E'
483  END IF
484  indibl = 1
485  indisp = indibl + n
486  indiwk = indisp + n
487  CALL sstebz( range, order, n, vll, vuu, il, iu, abstll,
488  $ rwork( indd ), rwork( inde ), m, nsplit, w,
489  $ iwork( indibl ), iwork( indisp ), rwork( indrwk ),
490  $ iwork( indiwk ), info )
491 *
492  IF( wantz ) THEN
493  CALL cstein( n, rwork( indd ), rwork( inde ), m, w,
494  $ iwork( indibl ), iwork( indisp ), z, ldz,
495  $ rwork( indrwk ), iwork( indiwk ), ifail, info )
496 *
497 * Apply unitary matrix used in reduction to tridiagonal
498 * form to eigenvectors returned by CSTEIN.
499 *
500  DO 20 j = 1, m
501  CALL ccopy( n, z( 1, j ), 1, work( 1 ), 1 )
502  CALL cgemv( 'N', n, n, cone, q, ldq, work, 1, czero,
503  $ z( 1, j ), 1 )
504  20 CONTINUE
505  END IF
506 *
507 * If matrix was scaled, then rescale eigenvalues appropriately.
508 *
509  30 CONTINUE
510  IF( iscale.EQ.1 ) THEN
511  IF( info.EQ.0 ) THEN
512  imax = m
513  ELSE
514  imax = info - 1
515  END IF
516  CALL sscal( imax, one / sigma, w, 1 )
517  END IF
518 *
519 * If eigenvalues are not in order, then sort them, along with
520 * eigenvectors.
521 *
522  IF( wantz ) THEN
523  DO 50 j = 1, m - 1
524  i = 0
525  tmp1 = w( j )
526  DO 40 jj = j + 1, m
527  IF( w( jj ).LT.tmp1 ) THEN
528  i = jj
529  tmp1 = w( jj )
530  END IF
531  40 CONTINUE
532 *
533  IF( i.NE.0 ) THEN
534  itmp1 = iwork( indibl+i-1 )
535  w( i ) = w( j )
536  iwork( indibl+i-1 ) = iwork( indibl+j-1 )
537  w( j ) = tmp1
538  iwork( indibl+j-1 ) = itmp1
539  CALL cswap( n, z( 1, i ), 1, z( 1, j ), 1 )
540  IF( info.NE.0 ) THEN
541  itmp1 = ifail( i )
542  ifail( i ) = ifail( j )
543  ifail( j ) = itmp1
544  END IF
545  END IF
546  50 CONTINUE
547  END IF
548 *
549  RETURN
550 *
551 * End of CHBEVX
552 *
subroutine sstebz(RANGE, ORDER, N, VL, VU, IL, IU, ABSTOL, D, E, M, NSPLIT, W, IBLOCK, ISPLIT, WORK, IWORK, INFO)
SSTEBZ
Definition: sstebz.f:275
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:145
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine chbtrd(VECT, UPLO, N, KD, AB, LDAB, D, E, Q, LDQ, WORK, INFO)
CHBTRD
Definition: chbtrd.f:165
subroutine csteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
CSTEQR
Definition: csteqr.f:134
subroutine cgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
CGEMV
Definition: cgemv.f:160
subroutine cstein(N, D, E, M, W, IBLOCK, ISPLIT, Z, LDZ, WORK, IWORK, IFAIL, INFO)
CSTEIN
Definition: cstein.f:184
real function clanhb(NORM, UPLO, N, K, AB, LDAB, WORK)
CLANHB returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a Hermitian band matrix.
Definition: clanhb.f:134
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:55
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
Definition: ccopy.f:52
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
Definition: cswap.f:52
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine ssterf(N, D, E, INFO)
SSTERF
Definition: ssterf.f:88
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:53

Here is the call graph for this function:

Here is the caller graph for this function: