Purpose
 =======
 
    LA_SPGVX and LA_HPGVX compute selected eigenvalues and, optionally, 
 the corresponding eigenvectors of generalized eigenvalue problems of 
 the form
       A*z = lambda*B*z, A*B*z = lambda*z, and B*A*z = lambda*z,
 where A and B are real symmetric in the case of LA_SPGVX and complex
 Hermitian in the case of LA_HPGVX. In both cases B is positive
 definite. Eigenvalues and eigenvectors can be selected by specifying
 either a range of values or a range of indices for the desired 
 eigenvalues. Matrices A and B are stored in a packed format.
 
 =========
 
       SUBROUTINE LA_SPGVX / LA_HPGVX( AP, BP, W, ITYPE= itype, &
             UPLO= uplo, Z= z, VL= vl, VU= vu, IL= il, IU= iu, M= m, &
             IFAIL= ifail, ABSTOL= abstol, INFO= info )
         <type>(<wp>), INTENT(INOUT) :: AP(:), BP(:)
         REAL(<wp>), INTENT(OUT) :: W(:)
         INTEGER, INTENT(IN), OPTIONAL :: ITYPE
         CHARACTER(LEN=1), INTENT(IN), OPTIONAL :: UPLO
         <type>(<wp>), INTENT(OUT), OPTIONAL :: Z(:,:)
         REAL(<wp>), INTENT(IN), OPTIONAL :: VL, VU
         INTEGER, INTENT(IN), OPTIONAL :: IL, IU
         INTEGER, INTENT(OUT), OPTIONAL :: M
         INTEGER, INTENT(OUT), OPTIONAL :: IFAIL(:)
         REAL(<wp>), INTENT(IN), OPTIONAL :: ABSTOL
         INTEGER, INTENT(OUT), OPTIONAL :: INFO
       where
         <type> ::= REAL | COMPLEX
         <wp>   ::= KIND(1.0) | KIND(1.0D0)
 
 Arguments
 =========
 
 AP       (input/output) REAL or COMPLEX array, shape (:) with 
          size(AP) = n*(n + 1)/2, where n is the order of A and B.
          On entry, the upper or lower triangle of matrix A in packed
 	  storage. The elements are stored columnwise as follows:
 	  if UPLO = 'U', AP(i +(j-1)*j/2) = A(i,j) for 1<=i<=j<=n;
 	  if UPLO = 'L', AP(i +(j-1)*(2*n-j)/2) = A(i,j) for 1<=j<=i<=n.
          On exit, the contents of AP are destroyed.
 BP       (input/output) REAL or COMPLEX array, shape (:) with 
          size(BP) = size(AP).
          On entry, the upper or lower triangle of matrix B in packed 
 	  storage. The elements are stored columnwise as follows:
 	  if UPLO = 'U', BP(i +(j-1)*j/2) = B(i,j) for 1<=i<=j<=n;
 	  if UPLO = 'L', BP(i +(j-1)*(2*n-j)/2) = B(i,j) for 1<=j<=i<=n.
          On exit, the triangular factor U or L of the Cholesky 
 	  factorization B = U^T*U or B = L*L^T, in the same storage 
 	  format as B.
 W        (output) REAL array, shape (:) with size(W) = n.
          The eigenvalues in ascending order.
 ITYPE    Optional (input) INTEGER.
          Specifies the problem type to be solved:
             = 1: A*z = lambda*B*z
             = 2: A*B*z = lambda*z
             = 3: B*A*z = lambda*z
          Default value: 1.
 UPLO     Optional (input) CHARACTER(LEN=1).
             = 'U': Upper triangles of A and B are stored;
             = 'L': Lower triangles of A and B are stored.
          Default value: 'U'.
 Z        Optional (output) REAL or COMPLEX rectangular array, shape 
          (:,:) with size(Z,1) = n and size(Z,2) = M.
          The first M columns of Z contain the orthonormal eigenvectors
 	  corresponding to the selected eigenvalues, with the i-th 
 	  column of Z holding the eigenvector associated with the 
 	  eigenvalue in W(i). The eigenvectors are normalized as 
 	  follows:
            if ITYPE = 1 or 2: Z^H * B * Z = I ,
            if ITYPE = 3: Z^H * B^-1 * 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.
 VL,VU    Optional (input) REAL.
          The lower and upper bounds of the interval to be searched for
 	  eigenvalues. VL < VU.
          Default values: VL = -HUGE(<wp>) and VU = HUGE(<wp>), where 
 	  <wp> ::= KIND(1.0) | KIND(1.0D0).
          Note: Neither VL nor VU may be present if IL and/or IU is 
 	  present.
 IL,IU    Optional (input) INTEGER.
          The indices of the smallest and largest eigenvalues to be 
 	  returned. The IL-th through IU-th eigenvalues will be found. 
 	  1 <= IL <= IU <= size(A,1).
          Default values: IL = 1 and IU = size(A,1).
          Note: Neither IL nor IU may be present if VL and/or VU is 
 	  present.
          Note: All eigenvalues are calculated if none of the arguments
 	  VL, VU, IL and IU are present.
 M        Optional (output) INTEGER.
          The total number of eigenvalues found. 0 <= M <= size(A,1).
          Note: If IL and IU are present then M = IU - IL + 1.
 IFAIL    Optional (output) INTEGER array, shape (:) with size(IFAIL) =
          size(A,1).
          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.
          Note: If Z is present then IFAIL should also be present.
 ABSTOL   Optional (input) 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 + EPSILON(1.0_<wp>) * max(|a|,|b|),
          where <wp> is the working precision. If ABSTOL <= 0, then 
 	  EPSILON(1.0_<wp>) * ||T||1 will be used in its place, where 
 	  ||T||1 is the l1 norm of the tridiagonal matrix obtained by
 	  reducing the generalized eigenvalue problem to tridiagonal
 	  form. Eigenvalues will be computed most accurately when
 	  ABSTOL is set to twice the underflow threshold 
 	  2 * LA_LAMCH(1.0_<wp>, 'S'), not zero.
          Default value: 0.0_<wp>.
          Note: If this routine returns with 0 < INFO <= n, then some 
 	  eigenvectors did not converge.
          Try setting ABSTOL to 2 * LA_LAMCH(1.0_<wp>, 'S').
 INFO     Optional (output) INTEGER.
          = 0: successful exit.
          < 0: if INFO = -i, the i-th argument had an illegal value
          > 0: the algorithm failed to converge or matrix B is not 
 	  positive definite:
             <= n: the algorithm failed to converge; if INFO = i, then 
 	       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 
 	       factorization of B could not be completed and no
               eigenvalues or eigenvectors were computed.
          If INFO is not present and an error occurs, then the program
          is terminated with an error message.