CHEGVX  -  compute  selected  eigenvalues, and optionally,
       eigenvectors of a complex  generalized  Hermitian-definite
       eigenproblem,     of     the     form    A*x=(lambda)*B*x,
       A*Bx=(lambda)*x, or B*A*x=(lambda)*x


SYNOPSIS

       SUBROUTINE CHEGVX( ITYPE, JOBZ, RANGE, UPLO, N, A, LDA, B,
                          LDB,  VL,  VU, IL, IU, ABSTOL, M, W, Z,
                          LDZ, WORK, LWORK, RWORK, IWORK,  IFAIL,
                          INFO )

           CHARACTER      JOBZ, RANGE, UPLO

           INTEGER        IL,  INFO,  ITYPE,  IU,  LDA, LDB, LDZ,
                          LWORK, M, N

           REAL           ABSTOL, VL, VU

           INTEGER        IFAIL( * ), IWORK( * )

           REAL           RWORK( * ), W( * )

           COMPLEX        A( LDA, * ), B( LDB, * ), WORK( * ), Z(
                          LDZ, * )


PURPOSE

       CHEGVX  computes  selected  eigenvalues,  and  optionally,
       eigenvectors of a complex  generalized  Hermitian-definite
       eigenproblem,     of     the     form    A*x=(lambda)*B*x,
       A*Bx=(lambda)*x, or B*A*x=(lambda)*x. Here  A  and  B  are
       assumed  to  be Hermitian and B is also positive definite.
       Eigenvalues and eigenvectors can be selected by specifying
       either  a  range  of  values or a range of indices for the
       desired eigenvalues.


ARGUMENTS

       ITYPE   (input) INTEGER
               Specifies the problem type to be solved:
               = 1:  A*x = (lambda)*B*x
               = 2:  A*B*x = (lambda)*x
               = 3:  B*A*x = (lambda)*x

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

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

               = 'U':  Upper triangles of A and B are stored;
               = 'L':  Lower triangles of A and B are stored.

       N       (input) INTEGER
               The order of the matrices A and B.  N >= 0.

       A       (input/output) COMPLEX array, dimension (LDA, N)
               On entry, the Hermitian 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) INTEGER
               The  leading  dimension  of  the  array A.  LDA >=
               max(1,N).

       B       (input/output) COMPLEX array, dimension (LDB, N)
               On entry, the Hermitian matrix B.  If UPLO =  'U',
               the leading N-by-N upper triangular part of B con­
               tains the upper triangular part of the  matrix  B.
               If UPLO = 'L', the leading N-by-N lower triangular
               part of B contains the lower  triangular  part  of
               the matrix B.

               On  exit,  if  INFO <= N, the part of B containing
               the matrix is overwritten by the triangular factor
               U  or L from the Cholesky factorization B = U**H*U
               or B = L*L**H.

       LDB     (input) INTEGER
               The leading dimension of  the  array  B.   LDB  >=
               max(1,N).

       VL      (input) REAL
               VU       (input)  REAL 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'.

       IL      (input) INTEGER
               IU      (input) INTEGER 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 refer­
               enced if RANGE = 'A' or 'V'.

               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
               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').

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

       W       (output) REAL array, dimension (N)
               The first M elements contain the  selected  eigen­
               values in ascending order.

       Z       (output) COMPLEX array, dimension (LDZ, max(1,M))
               If  JOBZ = 'N', then Z is not referenced.  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 asso­
               ciated with W(i).  The eigenvectors are normalized
               as  follows:  if  ITYPE = 1 or 2, Z**T*B*Z = I; if
               ITYPE = 3, Z**T*inv(B)*Z = 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.  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) INTEGER
               The  leading  dimension of the array Z.  LDZ >= 1,
               and if JOBZ = 'V', LDZ >= max(1,N).

       WORK    (workspace/output) COMPLEX array, dimension
               (LWORK)

               LWORK.

       LWORK   (input) INTEGER
               The  length  of  the   array   WORK.    LWORK   >=
               max(1,2*N-1).   For  optimal  efficiency, LWORK >=
               (NB+1)*N, where NB is  the  blocksize  for  CHETRD
               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.

       RWORK   (workspace) REAL array, dimension (7*N)

       IWORK   (workspace) INTEGER array, dimension (5*N)

       IFAIL   (output) INTEGER array, dimension (N)
               If  JOBZ = 'V', then if INFO = 0, the first M ele­
               ments 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.

       INFO    (output) INTEGER
               = 0:  successful exit
               < 0:  if INFO = -i, the i-th argument had an ille­
               gal value
               > 0:  CPOTRF or CHEEVX returned an error code:
               <= N:  if INFO = i, CHEEVX failed to  converge;  i
               eigenvectors  failed  to  converge.  Their indices
               are stored in array IFAIL.  > N:   if INFO =  N  +
               i,  for  1  <=  i  <= N, then the leading minor of
               order i of B is not positive definite.   The  fac­
               torization  of  B  could  not  be completed and no
               eigenvalues or eigenvectors were computed.


FURTHER DETAILS

       Based on contributions by
          Mark Fahey, Department of Mathematics,  Univ.  of  Ken­
       tucky, USA


Man(1) output converted with man2html