Purpose
 =======
 
      LA_SYEVX / LA_HEEVX compute selected eigenvalues and, optionally,
 the corresponding eigenvectors of a real symmetric/complex Hermitian 
 matrix A. Eigenvalues and eigenvectors can be selected by specifying
 either a range of values or a range of indices for the desired 
 eigenvalues.
 
 =========
 
         SUBROUTINE LA_SYEVX / LA_HEEVX ( A, W, JOBZ=jobz, UPLO=uplo, &
                        VL=vl, VU=vu, IL=il, IU=iu, M=m, IFAIL=ifail, &
                        ABSTOL=abstol, INFO=info )
              <type>(<wp>), INTENT(INOUT) :: A(:,:)
              REAL(<wp>), INTENT(OUT) :: W(:)
              CHARACTER(LEN=1), INTENT(IN), OPTIONAL :: JOBZ, UPLO
              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
 =========
 
 A      (input/output) REAL or COMPLEX square array, shape (:,:).
        On entry, the matrix A.
        If UPLO = 'U', the upper triangular part of A contains the upper
        triangular part of the matrix A. If UPLO = 'L', the lower 
        triangular part of A contains the lower triangular part of the
        matrix A.
        On exit:
        If JOBZ = 'V', then the first M columns of A contain the 
        orthonormal eigenvectors of the matrix A corresponding to the 
        selected eigenvalues, with the i-th column of A containing the
        eigenvector associated with the eigenvalue in W(i) . If an 
        eigenvector fails to converge, then that column of A contains the
        latest approximation to the eigenvector and the index of the 
        eigenvector is returned in IFAIL.
        If JOBZ = 'N', then the upper triangle (if UPLO = 'U') or the 
        lower triangle (if UPLO = 'L') of A, including the diagonal, is 
        destroyed.
 W      (output) REAL array, shape (:) with size(W) = size(A,1).
        The first M elements contain the selected eigenvalues in 
        ascending order.
 JOBZ   Optional (input) CHARACTER(LEN=1).
        = 'N': Computes eigenvalues only;
        = 'V': Computes eigenvalues and eigenvectors.
        Default value: 'N'.
 UPLO   Optional (input) CHARACTER(LEN=1).
        = 'U': Upper triangle of A is stored;
        = 'L': Lower triangle of A is stored.
        Default value: 'U'.
 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: IFAIL must be absent if JOBZ = 'N'.
 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 A
        to tridiagonal form. Eigenvalues will be computed most accurately 
        when ABSTOL is set to twice the underflow threshold
        2 * LA_LAMCH(1.0_<wp>, 'Save minimum'), not zero.
        Default value: 0.0_<wp>.
        Note: If this routine returns with INFO > 0, then some 
        eigenvectors did not converge. Try setting ABSTOL to 
        2 * LA_LAMCH(1.0_<wp>, 'Save minimum').
 INFO   Optional (output) 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.
        If INFO is not present and an error occurs, then the program is 
        terminated with an error message.