*     F08UPF Example Program Text
*     Mark 21.  NAG Copyright 2004.
*     .. Parameters ..
      INTEGER          NIN, NOUT
      PARAMETER        (NIN=5,NOUT=6)
      INTEGER          KAMAX, KBMAX, NMAX, MMAX
      PARAMETER        (KAMAX=5,KBMAX=5,NMAX=20,MMAX=10)
      INTEGER          LDAB, LDBB, LDQ, LDZ
      PARAMETER        (LDAB=KAMAX+1,LDBB=KBMAX+1,LDQ=NMAX,LDZ=NMAX)
      CHARACTER        UPLO
      PARAMETER        (UPLO='U')
      DOUBLE PRECISION ZERO
      PARAMETER        (ZERO=0.0D+0)
*     .. Local Scalars ..
      DOUBLE PRECISION ABSTOL, VL, VU
      INTEGER          I, IFAIL, IL, INFO, IU, J, KA, KB, M, N
*     .. Local Arrays ..
      COMPLEX *16      AB(LDAB,NMAX), BB(LDBB,NMAX), Q(LDQ,NMAX),
     +                 WORK(NMAX), Z(LDZ,MMAX)
      DOUBLE PRECISION RWORK(7*NMAX), W(NMAX)
      INTEGER          INDEX(NMAX), IWORK(5*NMAX)
*     .. External Subroutines ..
      EXTERNAL         X04DAF, ZHBGVX
*     .. Intrinsic Functions ..
      INTRINSIC        MAX, MIN
*     .. Executable Statements ..
      WRITE (NOUT,*) 'F08UPF Example Program Results'
      WRITE (NOUT,*)
*     Skip heading in data file
      READ (NIN,*)
      READ (NIN,*) N, KA, KB
      IF (N.LE.NMAX .AND. KA.LE.KAMAX .AND. KB.LE.KBMAX) THEN
*
*        Read the lower and upper bounds of the interval to be searched,
*        and read the upper or lower triangular parts of the matrices A
*        and B from data file
*
         READ (NIN,*) VL, VU
         IF (UPLO.EQ.'U') THEN
            READ (NIN,*) ((AB(KA+1+I-J,J),J=I,MIN(N,I+KA)),I=1,N)
            READ (NIN,*) ((BB(KB+1+I-J,J),J=I,MIN(N,I+KB)),I=1,N)
         ELSE IF (UPLO.EQ.'L') THEN
            READ (NIN,*) ((AB(1+I-J,J),J=MAX(1,I-KA),I),I=1,N)
            READ (NIN,*) ((BB(1+I-J,J),J=MAX(1,I-KB),I),I=1,N)
         END IF
*
*        Set the absolute error tolerance for eigenvalues.  With ABSTOL
*        set to zero, the default value is used instead
*
         ABSTOL = ZERO
*
*        Solve the generalized symmetric eigenvalue problem
*        A*x = lambda*B*x
*
         CALL ZHBGVX('Vectors','Values in range',UPLO,N,KA,KB,AB,LDAB,
     +               BB,LDBB,Q,LDQ,VL,VU,IL,IU,ABSTOL,M,W,Z,LDZ,WORK,
     +               RWORK,IWORK,INDEX,INFO)
*
         IF (INFO.GE.0 .AND. INFO.LE.N .AND. M.LE.MMAX) THEN
*
*           Print solution
*
            WRITE (NOUT,99999) 'Number of eigenvalues found =', M
            WRITE (NOUT,*)
            WRITE (NOUT,*) 'Eigenvalues'
            WRITE (NOUT,99998) (W(J),J=1,M)
*
            IFAIL = 0
            CALL X04DAF('General',' ',N,M,Z,LDZ,'Selected eigenvectors',
     +                  IFAIL)
            IF (INFO.GT.0) THEN
               WRITE (NOUT,99999)
     +           'INFO eigenvectors failed to converge, INFO =', INFO
               WRITE (NOUT,*)
     +           'Indices of eigenvectors that did not converge'
               WRITE (NOUT,99997) (INDEX(J),J=1,M)
            END IF
         ELSE IF (INFO.GT.N .AND. INFO.LE.2*N) THEN
            I = INFO - N
            WRITE (NOUT,99996) 'The leading minor of order ', I,
     +        ' of B is not positive definite'
         ELSE IF (M.GT.MMAX) THEN
            WRITE (NOUT,99995) 'M greater than MMAX, M =', M,
     +        ', MMAX =', MMAX
         ELSE
            WRITE (NOUT,99999) 'Failure in ZHBGVX. INFO =', INFO
         END IF
      ELSE
         WRITE (NOUT,*)
         WRITE (NOUT,*) 'NMAX too small'
      END IF
      STOP
*
99999 FORMAT (1X,A,I5)
99998 FORMAT (3X,(8F8.4))
99997 FORMAT (3X,(8I8))
99996 FORMAT (1X,A,I4,A)
99995 FORMAT (1X,A,I5,A,I5)
      END