Contents


NAME

     ssyevr  -  compute  selected  eigenvalues  and,  optionally,
     eigenvectors of a real symmetric tridiagonal matrix T

SYNOPSIS

     SUBROUTINE SSYEVR(JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU,
           ABSTOL, M, W, Z, LDZ, ISUPPZ, WORK, LWORK, IWORK, LIWORK, INFO)

     CHARACTER * 1 JOBZ, RANGE, UPLO
     INTEGER N, LDA, IL, IU, M, LDZ, LWORK, LIWORK, INFO
     INTEGER ISUPPZ(*), IWORK(*)
     REAL VL, VU, ABSTOL
     REAL A(LDA,*), W(*), Z(LDZ,*), WORK(*)

     SUBROUTINE SSYEVR_64(JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU,
           ABSTOL, M, W, Z, LDZ, ISUPPZ, WORK, LWORK, IWORK, LIWORK, INFO)

     CHARACTER * 1 JOBZ, RANGE, UPLO
     INTEGER*8 N, LDA, IL, IU, M, LDZ, LWORK, LIWORK, INFO
     INTEGER*8 ISUPPZ(*), IWORK(*)
     REAL VL, VU, ABSTOL
     REAL A(LDA,*), W(*), Z(LDZ,*), WORK(*)

  F95 INTERFACE
     SUBROUTINE SYEVR(JOBZ, RANGE, UPLO, [N], A, [LDA], VL, VU, IL, IU,
            ABSTOL, M, W, Z, [LDZ], ISUPPZ, [WORK], [LWORK], [IWORK], [LIWORK],
            [INFO])

     CHARACTER(LEN=1) :: JOBZ, RANGE, UPLO
     INTEGER :: N, LDA, IL, IU, M, LDZ, LWORK, LIWORK, INFO
     INTEGER, DIMENSION(:) :: ISUPPZ, IWORK
     REAL :: VL, VU, ABSTOL
     REAL, DIMENSION(:) :: W, WORK
     REAL, DIMENSION(:,:) :: A, Z

     SUBROUTINE SYEVR_64(JOBZ, RANGE, UPLO, [N], A, [LDA], VL, VU, IL, IU,
            ABSTOL, M, W, Z, [LDZ], ISUPPZ, [WORK], [LWORK], [IWORK], [LIWORK],
            [INFO])

     CHARACTER(LEN=1) :: JOBZ, RANGE, UPLO
     INTEGER(8) :: N, LDA, IL, IU, M, LDZ, LWORK, LIWORK, INFO
     INTEGER(8), DIMENSION(:) :: ISUPPZ, IWORK
     REAL :: VL, VU, ABSTOL
     REAL, DIMENSION(:) :: W, WORK
     REAL, DIMENSION(:,:) :: A, Z
  C INTERFACE
     #include <sunperf.h>

     void ssyevr(char jobz, char range, char uplo, int  n,  float
               *a,  int  lda, float vl, float vu, int il, int iu,
               float abstol, int *m, float *w, float *z, int ldz,
               int *isuppz, int *info);

     void ssyevr_64(char jobz, char range,  char  uplo,  long  n,
               float  *a,  long lda, float vl, float vu, long il,
               long iu, float abstol, long *m,  float  *w,  float
               *z, long ldz, long *isuppz, long *info);

PURPOSE

     ssyevr computes selected eigenvalues and, optionally, eigen-
     vectors  of  a  real symmetric tridiagonal matrix T.  Eigen-
     values and eigenvectors can be selected by specifying either
     a  range  of  values  or  a range of indices for the desired
     eigenvalues.

     Whenever possible, SSYEVR calls SSTEGR to compute the
     eigenspectrum  using  Relatively   Robust   Representations.
     SSTEGR  computes  eigenvalues  by  the dqds algorithm, while
     orthogonal eigenvectors are computed from various "good" L D
     L^T   representations   (also  known  as  Relatively  Robust
     Representations). Gram-Schmidt orthogonalization is  avoided
     as  far as possible. More specifically, the various steps of
     the algorithm are as follows. For the i-th  unreduced  block
     of T,
        (a) Compute T - sigma_i = L_i D_i L_i^T,  such  that  L_i
     D_i L_i^T
             is a relatively robust representation,
        (b) Compute the eigenvalues, lambda_j, of L_i  D_i  L_i^T
     to high
            relative accuracy by the dqds algorithm,
        (c) If there is a cluster of close eigenvalues,  "choose"
     sigma_i
            close to the cluster, and go to step (a),
        (d) Given the approximate eigenvalue lambda_j of L_i  D_i
     L_i^T,
            compute the corresponding eigenvector by forming a
            rank-revealing twisted factorization.
     The desired accuracy of the output can be specified  by  the
     input parameter ABSTOL.

     For more details, see "A new O(n^2) algorithm for  the  sym-
     metric   tridiagonal   eigenvalue/eigenvector  problem",  by
     Inderjit Dhillon, Computer Science Division Technical Report
     No. UCB//CSD-97-971, UC Berkeley, May 1997.
     Note 1 : SSYEVR calls  SSTEGR  when  the  full  spectrum  is
     requested on machines which conform to the ieee-754 floating
     point standard.  SSYEVR calls SSTEBZ and SSTEIN on  non-ieee
     machines and
     when partial spectrum requests are made.

     Normal execution of SSTEGR may create  NaNs  and  infinities
     and  hence  may  abort  due to a floating point exception in
     environments which do not handle NaNs and infinities in  the
     ieee standard default manner.

ARGUMENTS

     JOBZ (input)
               = 'N':  Compute eigenvalues only;
               = 'V':  Compute eigenvalues and eigenvectors.

     RANGE (input)
               = '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.

     UPLO (input)
               = 'U':  Upper triangle of A is stored;
               = 'L':  Lower triangle of A is stored.

     N (input) The order of the matrix A.  N >= 0.

     A (input/output)
               On entry, the symmetric matrix A.  If UPLO =  'U',
               the leading N-by-N upper triangular part of A con-
               tains the upper triangular part of the  matrix  A.
               If UPLO = 'L', the leading N-by-N lower triangular
               part of A contains the lower  triangular  part  of
               the  matrix  A.   On  exit, the lower triangle (if
               UPLO='L') or the upper triangle (if  UPLO='U')  of
               A, including the diagonal, is destroyed.

     LDA (input)
               The leading dimension of  the  array  A.   LDA  >=
               max(1,N).

     VL (input)
               If RANGE='V', the lower and upper  bounds  of  the
               interval  to be searched for eigenvalues. VL < VU.
               Not referenced if RANGE = 'A' or 'I'.

     VU (input)
               See the description of VL.

     IL (input)
               If RANGE='I', the indices (in ascending order)  of
               the   smallest   and  largest  eigenvalues  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'.

     IU (input)
               See the description of IL.

     ABSTOL (input)
               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 A to tri-
               diagonal form.

               See "Computing Small Singular Values of Bidiagonal
               Matrices  with Guaranteed High Relative Accuracy,"
               by Demmel and Kahan, LAPACK Working Note #3.

               If high relative accuracy is important, set ABSTOL
               to  SLAMCH(  'Safe  minimum'  ).   Doing  so  will
               guarantee that eigenvalues are  computed  to  high
               relative   accuracy   when   possible   in  future
               releases.  The current  code  does  not  make  any
               guarantees about high relative accuracy, but furu-
               tre releases will. See J. Barlow  and  J.  Demmel,
               "Computing  Accurate Eigensystems of Scaled Diago-
               nally Dominant Matrices", LAPACK Working Note  #7,
               for  a  discussion  of which matrices define their
               eigenvalues to high relative accuracy.

     M (output)
               The total number of eigenvalues found.  0 <= M  <=
               N.  If RANGE = 'A', M = N, and if RANGE = 'I', M =
               IU-IL+1.

     W (output)
               The first M elements contain the  selected  eigen-
               values in ascending order.

     Z (input) 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 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.

     LDZ (input)
               The leading dimension of the array Z.  LDZ  >=  1,
               and if JOBZ = 'V', LDZ >= max(1,N).

     ISUPPZ (output)
               The support of the eigenvectors in  Z,  i.e.,  the
               indices  indicating the nonzero elements in Z. The
               i-th  eigenvector  is  nonzero  only  in  elements
               ISUPPZ( 2*i-1 ) through ISUPPZ( 2*i ).

     WORK (workspace)
               On exit, if INFO = 0, WORK(1) returns the  optimal
               LWORK.

     LWORK (input)
               The  dimension  of  the  array  WORK.   LWORK   >=
               max(1,26*N).   For  optimal  efficiency,  LWORK >=
               (NB+6)*N, where NB is the max of the blocksize for
               SSYTRD and SORMTR returned by ILAENV.

               If LWORK = -1, then a workspace query is  assumed;
               the  routine  only  calculates the optimal size of
               the WORK array, returns this value  as  the  first
               entry  of  the  WORK  array,  and no error message
               related to LWORK is issued by XERBLA.
     IWORK (workspace/output)
               On exit, if INFO = 0, IWORK(1) returns the optimal
               LWORK.

     LIWORK (input)
               The dimension  of  the  array  IWORK.   LIWORK  >=
               max(1,10*N).

               If LIWORK = -1, then a workspace query is assumed;
               the  routine  only  calculates the optimal size of
               the IWORK array, returns this value as  the  first
               entry  of  the  IWORK  array, and no error message
               related to LIWORK is issued by XERBLA.

     INFO (output)
               = 0:  successful exit
               < 0:  if INFO = -i, the i-th argument had an ille-
               gal value
               > 0:  Internal error

FURTHER DETAILS

     Based on contributions by
        Inderjit Dhillon, IBM Almaden, USA
        Osni Marques, LBNL/NERSC, USA
        Ken Stanley, Computer Science Division, University of
          California at Berkeley, USA