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

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

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

Purpose:
 DSBEVX computes selected eigenvalues and, optionally, eigenvectors
 of a real symmetric 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 DOUBLE PRECISION array, dimension (LDAB, N)
          On entry, the upper or lower triangle of the symmetric 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.  If UPLO = 'U', the first
          superdiagonal and the diagonal of the tridiagonal matrix T
          are returned in rows KD and KD+1 of AB, and if UPLO = 'L',
          the diagonal and first subdiagonal of T are returned in the
          first two rows of AB.
[in]LDAB
          LDAB is INTEGER
          The leading dimension of the array AB.  LDAB >= KD + 1.
[out]Q
          Q is DOUBLE PRECISION array, dimension (LDQ, N)
          If JOBZ = 'V', the N-by-N orthogonal 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 DOUBLE PRECISION
          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 DOUBLE PRECISION
          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 DOUBLE PRECISION
          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*DLAMCH('S'), not zero.
          If this routine returns with INFO>0, indicating that some
          eigenvectors did not converge, try setting ABSTOL to
          2*DLAMCH('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 DOUBLE PRECISION array, dimension (N)
          The first M elements contain the selected eigenvalues in
          ascending order.
[out]Z
          Z is DOUBLE PRECISION 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 DOUBLE PRECISION 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 267 of file dsbevx.f.

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

Here is the call graph for this function:

Here is the caller graph for this function: