stdlib_lapack_eigv_gen.fypp Source File


Source Code

#:include "common.fypp" 
submodule(stdlib_lapack_eig_svd_lsq) stdlib_lapack_eigv_gen
  implicit none


  contains
#:for ik,it,ii in LINALG_INT_KINDS_TYPES

     module subroutine stdlib${ii}$_sgeev( jobvl, jobvr, n, a, lda, wr, wi, vl, ldvl, vr,ldvr, work, lwork, &
     !! SGEEV computes for an N-by-N real nonsymmetric matrix A, the
     !! eigenvalues and, optionally, the left and/or right eigenvectors.
     !! The right eigenvector v(j) of A satisfies
     !! A * v(j) = lambda(j) * v(j)
     !! where lambda(j) is its eigenvalue.
     !! The left eigenvector u(j) of A satisfies
     !! u(j)**H * A = lambda(j) * u(j)**H
     !! where u(j)**H denotes the conjugate-transpose of u(j).
     !! The computed eigenvectors are normalized to have Euclidean norm
     !! equal to 1 and largest component real.
               info )
        ! -- lapack driver routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: jobvl, jobvr
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: lda, ldvl, ldvr, lwork, n
           ! Array Arguments 
           real(sp), intent(inout) :: a(lda,*)
           real(sp), intent(out) :: vl(ldvl,*), vr(ldvr,*), wi(*), work(*), wr(*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lquery, scalea, wantvl, wantvr
           character :: side
           integer(${ik}$) :: hswork, i, ibal, ierr, ihi, ilo, itau, iwrk, k, lwork_trevc, maxwrk, &
                     minwrk, nout
           real(sp) :: anrm, bignum, cs, cscale, eps, r, scl, smlnum, sn
           ! Local Arrays 
           logical(lk) :: select(1_${ik}$)
           real(sp) :: dum(1_${ik}$)
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input arguments
           info = 0_${ik}$
           lquery = ( lwork==-1_${ik}$ )
           wantvl = stdlib_lsame( jobvl, 'V' )
           wantvr = stdlib_lsame( jobvr, 'V' )
           if( ( .not.wantvl ) .and. ( .not.stdlib_lsame( jobvl, 'N' ) ) ) then
              info = -1_${ik}$
           else if( ( .not.wantvr ) .and. ( .not.stdlib_lsame( jobvr, 'N' ) ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           else if( lda<max( 1_${ik}$, n ) ) then
              info = -5_${ik}$
           else if( ldvl<1_${ik}$ .or. ( wantvl .and. ldvl<n ) ) then
              info = -9_${ik}$
           else if( ldvr<1_${ik}$ .or. ( wantvr .and. ldvr<n ) ) then
              info = -11_${ik}$
           end if
           ! compute workspace
            ! (note: comments in the code beginning "workspace:" describe the
             ! minimal amount of workspace needed at that point in the code,
             ! as well as the preferred amount for good performance.
             ! nb refers to the optimal block size for the immediately
             ! following subroutine, as returned by stdlib${ii}$_ilaenv.
             ! hswork refers to the workspace preferred by stdlib${ii}$_shseqr, as
             ! calculated below. hswork is computed assuming ilo=1 and ihi=n,
             ! the worst case.)
           if( info==0_${ik}$ ) then
              if( n==0_${ik}$ ) then
                 minwrk = 1_${ik}$
                 maxwrk = 1_${ik}$
              else
                 maxwrk = 2_${ik}$*n + n*stdlib${ii}$_ilaenv( 1_${ik}$, 'SGEHRD', ' ', n, 1_${ik}$, n, 0_${ik}$ )
                 if( wantvl ) then
                    minwrk = 4_${ik}$*n
                    maxwrk = max( maxwrk, 2_${ik}$*n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$,'SORGHR', ' ', n, 1_${ik}$, n,&
                               -1_${ik}$ ) )
                    call stdlib${ii}$_shseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vl, ldvl,work, -1_${ik}$, &
                              info )
                    hswork = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + 1_${ik}$, n + hswork )
                    call stdlib${ii}$_strevc3( 'L', 'B', select, n, a, lda,vl, ldvl, vr, ldvr, n, nout,&
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    maxwrk = max( maxwrk, 4_${ik}$*n )
                 else if( wantvr ) then
                    minwrk = 4_${ik}$*n
                    maxwrk = max( maxwrk, 2_${ik}$*n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$,'SORGHR', ' ', n, 1_${ik}$, n,&
                               -1_${ik}$ ) )
                    call stdlib${ii}$_shseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vr, ldvr,work, -1_${ik}$, &
                              info )
                    hswork = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + 1_${ik}$, n + hswork )
                    call stdlib${ii}$_strevc3( 'R', 'B', select, n, a, lda,vl, ldvl, vr, ldvr, n, nout,&
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    maxwrk = max( maxwrk, 4_${ik}$*n )
                 else
                    minwrk = 3_${ik}$*n
                    call stdlib${ii}$_shseqr( 'E', 'N', n, 1_${ik}$, n, a, lda, wr, wi, vr, ldvr,work, -1_${ik}$, &
                              info )
                    hswork = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + 1_${ik}$, n + hswork )
                 end if
                 maxwrk = max( maxwrk, minwrk )
              end if
              work( 1_${ik}$ ) = maxwrk
              if( lwork<minwrk .and. .not.lquery ) then
                 info = -13_${ik}$
              end if
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'SGEEV ', -info )
              return
           else if( lquery ) then
              return
           end if
           ! quick return if possible
           if( n==0 )return
           ! get machine constants
           eps = stdlib${ii}$_slamch( 'P' )
           smlnum = stdlib${ii}$_slamch( 'S' )
           bignum = one / smlnum
           call stdlib${ii}$_slabad( smlnum, bignum )
           smlnum = sqrt( smlnum ) / eps
           bignum = one / smlnum
           ! scale a if max element outside range [smlnum,bignum]
           anrm = stdlib${ii}$_slange( 'M', n, n, a, lda, dum )
           scalea = .false.
           if( anrm>zero .and. anrm<smlnum ) then
              scalea = .true.
              cscale = smlnum
           else if( anrm>bignum ) then
              scalea = .true.
              cscale = bignum
           end if
           if( scalea )call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, anrm, cscale, n, n, a, lda, ierr )
           ! balance the matrix
           ! (workspace: need n)
           ibal = 1_${ik}$
           call stdlib${ii}$_sgebal( 'B', n, a, lda, ilo, ihi, work( ibal ), ierr )
           ! reduce to upper hessenberg form
           ! (workspace: need 3*n, prefer 2*n+n*nb)
           itau = ibal + n
           iwrk = itau + n
           call stdlib${ii}$_sgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),lwork-iwrk+1, ierr &
                     )
           if( wantvl ) then
              ! want left eigenvectors
              ! copy householder vectors to vl
              side = 'L'
              call stdlib${ii}$_slacpy( 'L', n, n, a, lda, vl, ldvl )
              ! generate orthogonal matrix in vl
              ! (workspace: need 3*n-1, prefer 2*n+(n-1)*nb)
              call stdlib${ii}$_sorghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vl
              ! (workspace: need n+1, prefer n+hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_shseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vl, ldvl,work( iwrk ), &
                        lwork-iwrk+1, info )
              if( wantvr ) then
                 ! want left and right eigenvectors
                 ! copy schur vectors to vr
                 side = 'B'
                 call stdlib${ii}$_slacpy( 'F', n, n, vl, ldvl, vr, ldvr )
              end if
           else if( wantvr ) then
              ! want right eigenvectors
              ! copy householder vectors to vr
              side = 'R'
              call stdlib${ii}$_slacpy( 'L', n, n, a, lda, vr, ldvr )
              ! generate orthogonal matrix in vr
              ! (workspace: need 3*n-1, prefer 2*n+(n-1)*nb)
              call stdlib${ii}$_sorghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vr
              ! (workspace: need n+1, prefer n+hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_shseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           else
              ! compute eigenvalues only
              ! (workspace: need n+1, prefer n+hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_shseqr( 'E', 'N', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           end if
           ! if info /= 0 from stdlib${ii}$_shseqr, then quit
           if( info/=0 )go to 50
           if( wantvl .or. wantvr ) then
              ! compute left and/or right eigenvectors
              ! (workspace: need 4*n, prefer n + n + 2*n*nb)
              call stdlib${ii}$_strevc3( side, 'B', select, n, a, lda, vl, ldvl, vr, ldvr,n, nout, work(&
                         iwrk ), lwork-iwrk+1, ierr )
           end if
           if( wantvl ) then
              ! undo balancing of left eigenvectors
              ! (workspace: need n)
              call stdlib${ii}$_sgebak( 'B', 'L', n, ilo, ihi, work( ibal ), n, vl, ldvl,ierr )
              ! normalize left eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_snrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_sscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_slapy2( stdlib${ii}$_snrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_snrm2( n, &
                              vl( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_sscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_sscal( n, scl, vl( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( iwrk+k-1 ) = vl( k, i )**2_${ik}$ + vl( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_isamax( n, work( iwrk ), 1_${ik}$ )
                    call stdlib${ii}$_slartg( vl( k, i ), vl( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_srot( n, vl( 1_${ik}$, i ), 1_${ik}$, vl( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vl( k, i+1 ) = zero
                 end if
              end do
           end if
           if( wantvr ) then
              ! undo balancing of right eigenvectors
              ! (workspace: need n)
              call stdlib${ii}$_sgebak( 'B', 'R', n, ilo, ihi, work( ibal ), n, vr, ldvr,ierr )
              ! normalize right eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_snrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_sscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_slapy2( stdlib${ii}$_snrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_snrm2( n, &
                              vr( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_sscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_sscal( n, scl, vr( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( iwrk+k-1 ) = vr( k, i )**2_${ik}$ + vr( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_isamax( n, work( iwrk ), 1_${ik}$ )
                    call stdlib${ii}$_slartg( vr( k, i ), vr( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_srot( n, vr( 1_${ik}$, i ), 1_${ik}$, vr( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vr( k, i+1 ) = zero
                 end if
              end do
           end if
           ! undo scaling if necessary
           50 continue
           if( scalea ) then
              call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wr( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wi( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              if( info>0_${ik}$ ) then
                 call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wr, n,ierr )
                 call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wi, n,ierr )
              end if
           end if
           work( 1_${ik}$ ) = maxwrk
           return
     end subroutine stdlib${ii}$_sgeev

     module subroutine stdlib${ii}$_dgeev( jobvl, jobvr, n, a, lda, wr, wi, vl, ldvl, vr,ldvr, work, lwork, &
     !! DGEEV computes for an N-by-N real nonsymmetric matrix A, the
     !! eigenvalues and, optionally, the left and/or right eigenvectors.
     !! The right eigenvector v(j) of A satisfies
     !! A * v(j) = lambda(j) * v(j)
     !! where lambda(j) is its eigenvalue.
     !! The left eigenvector u(j) of A satisfies
     !! u(j)**H * A = lambda(j) * u(j)**H
     !! where u(j)**H denotes the conjugate-transpose of u(j).
     !! The computed eigenvectors are normalized to have Euclidean norm
     !! equal to 1 and largest component real.
               info )
        ! -- lapack driver routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: jobvl, jobvr
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: lda, ldvl, ldvr, lwork, n
           ! Array Arguments 
           real(dp), intent(inout) :: a(lda,*)
           real(dp), intent(out) :: vl(ldvl,*), vr(ldvr,*), wi(*), work(*), wr(*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lquery, scalea, wantvl, wantvr
           character :: side
           integer(${ik}$) :: hswork, i, ibal, ierr, ihi, ilo, itau, iwrk, k, lwork_trevc, maxwrk, &
                     minwrk, nout
           real(dp) :: anrm, bignum, cs, cscale, eps, r, scl, smlnum, sn
           ! Local Arrays 
           logical(lk) :: select(1_${ik}$)
           real(dp) :: dum(1_${ik}$)
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input arguments
           info = 0_${ik}$
           lquery = ( lwork==-1_${ik}$ )
           wantvl = stdlib_lsame( jobvl, 'V' )
           wantvr = stdlib_lsame( jobvr, 'V' )
           if( ( .not.wantvl ) .and. ( .not.stdlib_lsame( jobvl, 'N' ) ) ) then
              info = -1_${ik}$
           else if( ( .not.wantvr ) .and. ( .not.stdlib_lsame( jobvr, 'N' ) ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           else if( lda<max( 1_${ik}$, n ) ) then
              info = -5_${ik}$
           else if( ldvl<1_${ik}$ .or. ( wantvl .and. ldvl<n ) ) then
              info = -9_${ik}$
           else if( ldvr<1_${ik}$ .or. ( wantvr .and. ldvr<n ) ) then
              info = -11_${ik}$
           end if
           ! compute workspace
            ! (note: comments in the code beginning "workspace:" describe the
             ! minimal amount of workspace needed at that point in the code,
             ! as well as the preferred amount for good performance.
             ! nb refers to the optimal block size for the immediately
             ! following subroutine, as returned by stdlib${ii}$_ilaenv.
             ! hswork refers to the workspace preferred by stdlib${ii}$_dhseqr, as
             ! calculated below. hswork is computed assuming ilo=1 and ihi=n,
             ! the worst case.)
           if( info==0_${ik}$ ) then
              if( n==0_${ik}$ ) then
                 minwrk = 1_${ik}$
                 maxwrk = 1_${ik}$
              else
                 maxwrk = 2_${ik}$*n + n*stdlib${ii}$_ilaenv( 1_${ik}$, 'DGEHRD', ' ', n, 1_${ik}$, n, 0_${ik}$ )
                 if( wantvl ) then
                    minwrk = 4_${ik}$*n
                    maxwrk = max( maxwrk, 2_${ik}$*n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$,'DORGHR', ' ', n, 1_${ik}$, n,&
                               -1_${ik}$ ) )
                    call stdlib${ii}$_dhseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vl, ldvl,work, -1_${ik}$, &
                              info )
                    hswork = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + 1_${ik}$, n + hswork )
                    call stdlib${ii}$_dtrevc3( 'L', 'B', select, n, a, lda,vl, ldvl, vr, ldvr, n, nout,&
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    maxwrk = max( maxwrk, 4_${ik}$*n )
                 else if( wantvr ) then
                    minwrk = 4_${ik}$*n
                    maxwrk = max( maxwrk, 2_${ik}$*n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$,'DORGHR', ' ', n, 1_${ik}$, n,&
                               -1_${ik}$ ) )
                    call stdlib${ii}$_dhseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vr, ldvr,work, -1_${ik}$, &
                              info )
                    hswork = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + 1_${ik}$, n + hswork )
                    call stdlib${ii}$_dtrevc3( 'R', 'B', select, n, a, lda,vl, ldvl, vr, ldvr, n, nout,&
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    maxwrk = max( maxwrk, 4_${ik}$*n )
                 else
                    minwrk = 3_${ik}$*n
                    call stdlib${ii}$_dhseqr( 'E', 'N', n, 1_${ik}$, n, a, lda, wr, wi, vr, ldvr,work, -1_${ik}$, &
                              info )
                    hswork = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + 1_${ik}$, n + hswork )
                 end if
                 maxwrk = max( maxwrk, minwrk )
              end if
              work( 1_${ik}$ ) = maxwrk
              if( lwork<minwrk .and. .not.lquery ) then
                 info = -13_${ik}$
              end if
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'DGEEV ', -info )
              return
           else if( lquery ) then
              return
           end if
           ! quick return if possible
           if( n==0 )return
           ! get machine constants
           eps = stdlib${ii}$_dlamch( 'P' )
           smlnum = stdlib${ii}$_dlamch( 'S' )
           bignum = one / smlnum
           call stdlib${ii}$_dlabad( smlnum, bignum )
           smlnum = sqrt( smlnum ) / eps
           bignum = one / smlnum
           ! scale a if max element outside range [smlnum,bignum]
           anrm = stdlib${ii}$_dlange( 'M', n, n, a, lda, dum )
           scalea = .false.
           if( anrm>zero .and. anrm<smlnum ) then
              scalea = .true.
              cscale = smlnum
           else if( anrm>bignum ) then
              scalea = .true.
              cscale = bignum
           end if
           if( scalea )call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, anrm, cscale, n, n, a, lda, ierr )
           ! balance the matrix
           ! (workspace: need n)
           ibal = 1_${ik}$
           call stdlib${ii}$_dgebal( 'B', n, a, lda, ilo, ihi, work( ibal ), ierr )
           ! reduce to upper hessenberg form
           ! (workspace: need 3*n, prefer 2*n+n*nb)
           itau = ibal + n
           iwrk = itau + n
           call stdlib${ii}$_dgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),lwork-iwrk+1, ierr &
                     )
           if( wantvl ) then
              ! want left eigenvectors
              ! copy householder vectors to vl
              side = 'L'
              call stdlib${ii}$_dlacpy( 'L', n, n, a, lda, vl, ldvl )
              ! generate orthogonal matrix in vl
              ! (workspace: need 3*n-1, prefer 2*n+(n-1)*nb)
              call stdlib${ii}$_dorghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vl
              ! (workspace: need n+1, prefer n+hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_dhseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vl, ldvl,work( iwrk ), &
                        lwork-iwrk+1, info )
              if( wantvr ) then
                 ! want left and right eigenvectors
                 ! copy schur vectors to vr
                 side = 'B'
                 call stdlib${ii}$_dlacpy( 'F', n, n, vl, ldvl, vr, ldvr )
              end if
           else if( wantvr ) then
              ! want right eigenvectors
              ! copy householder vectors to vr
              side = 'R'
              call stdlib${ii}$_dlacpy( 'L', n, n, a, lda, vr, ldvr )
              ! generate orthogonal matrix in vr
              ! (workspace: need 3*n-1, prefer 2*n+(n-1)*nb)
              call stdlib${ii}$_dorghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vr
              ! (workspace: need n+1, prefer n+hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_dhseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           else
              ! compute eigenvalues only
              ! (workspace: need n+1, prefer n+hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_dhseqr( 'E', 'N', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           end if
           ! if info /= 0 from stdlib${ii}$_dhseqr, then quit
           if( info/=0 )go to 50
           if( wantvl .or. wantvr ) then
              ! compute left and/or right eigenvectors
              ! (workspace: need 4*n, prefer n + n + 2*n*nb)
              call stdlib${ii}$_dtrevc3( side, 'B', select, n, a, lda, vl, ldvl, vr, ldvr,n, nout, work(&
                         iwrk ), lwork-iwrk+1, ierr )
           end if
           if( wantvl ) then
              ! undo balancing of left eigenvectors
              ! (workspace: need n)
              call stdlib${ii}$_dgebak( 'B', 'L', n, ilo, ihi, work( ibal ), n, vl, ldvl,ierr )
              ! normalize left eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_dnrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_dscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_dlapy2( stdlib${ii}$_dnrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_dnrm2( n, &
                              vl( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_dscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_dscal( n, scl, vl( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( iwrk+k-1 ) = vl( k, i )**2_${ik}$ + vl( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_idamax( n, work( iwrk ), 1_${ik}$ )
                    call stdlib${ii}$_dlartg( vl( k, i ), vl( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_drot( n, vl( 1_${ik}$, i ), 1_${ik}$, vl( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vl( k, i+1 ) = zero
                 end if
              end do
           end if
           if( wantvr ) then
              ! undo balancing of right eigenvectors
              ! (workspace: need n)
              call stdlib${ii}$_dgebak( 'B', 'R', n, ilo, ihi, work( ibal ), n, vr, ldvr,ierr )
              ! normalize right eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_dnrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_dscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_dlapy2( stdlib${ii}$_dnrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_dnrm2( n, &
                              vr( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_dscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_dscal( n, scl, vr( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( iwrk+k-1 ) = vr( k, i )**2_${ik}$ + vr( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_idamax( n, work( iwrk ), 1_${ik}$ )
                    call stdlib${ii}$_dlartg( vr( k, i ), vr( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_drot( n, vr( 1_${ik}$, i ), 1_${ik}$, vr( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vr( k, i+1 ) = zero
                 end if
              end do
           end if
           ! undo scaling if necessary
           50 continue
           if( scalea ) then
              call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wr( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wi( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              if( info>0_${ik}$ ) then
                 call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wr, n,ierr )
                 call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wi, n,ierr )
              end if
           end if
           work( 1_${ik}$ ) = maxwrk
           return
     end subroutine stdlib${ii}$_dgeev

#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
     module subroutine stdlib${ii}$_${ri}$geev( jobvl, jobvr, n, a, lda, wr, wi, vl, ldvl, vr,ldvr, work, lwork, &
     !! DGEEV: computes for an N-by-N real nonsymmetric matrix A, the
     !! eigenvalues and, optionally, the left and/or right eigenvectors.
     !! The right eigenvector v(j) of A satisfies
     !! A * v(j) = lambda(j) * v(j)
     !! where lambda(j) is its eigenvalue.
     !! The left eigenvector u(j) of A satisfies
     !! u(j)**H * A = lambda(j) * u(j)**H
     !! where u(j)**H denotes the conjugate-transpose of u(j).
     !! The computed eigenvectors are normalized to have Euclidean norm
     !! equal to 1 and largest component real.
               info )
        ! -- lapack driver routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_${rk}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: jobvl, jobvr
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: lda, ldvl, ldvr, lwork, n
           ! Array Arguments 
           real(${rk}$), intent(inout) :: a(lda,*)
           real(${rk}$), intent(out) :: vl(ldvl,*), vr(ldvr,*), wi(*), work(*), wr(*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lquery, scalea, wantvl, wantvr
           character :: side
           integer(${ik}$) :: hswork, i, ibal, ierr, ihi, ilo, itau, iwrk, k, lwork_trevc, maxwrk, &
                     minwrk, nout
           real(${rk}$) :: anrm, bignum, cs, cscale, eps, r, scl, smlnum, sn
           ! Local Arrays 
           logical(lk) :: select(1_${ik}$)
           real(${rk}$) :: dum(1_${ik}$)
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input arguments
           info = 0_${ik}$
           lquery = ( lwork==-1_${ik}$ )
           wantvl = stdlib_lsame( jobvl, 'V' )
           wantvr = stdlib_lsame( jobvr, 'V' )
           if( ( .not.wantvl ) .and. ( .not.stdlib_lsame( jobvl, 'N' ) ) ) then
              info = -1_${ik}$
           else if( ( .not.wantvr ) .and. ( .not.stdlib_lsame( jobvr, 'N' ) ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           else if( lda<max( 1_${ik}$, n ) ) then
              info = -5_${ik}$
           else if( ldvl<1_${ik}$ .or. ( wantvl .and. ldvl<n ) ) then
              info = -9_${ik}$
           else if( ldvr<1_${ik}$ .or. ( wantvr .and. ldvr<n ) ) then
              info = -11_${ik}$
           end if
           ! compute workspace
            ! (note: comments in the code beginning "workspace:" describe the
             ! minimal amount of workspace needed at that point in the code,
             ! as well as the preferred amount for good performance.
             ! nb refers to the optimal block size for the immediately
             ! following subroutine, as returned by stdlib${ii}$_ilaenv.
             ! hswork refers to the workspace preferred by stdlib${ii}$_${ri}$hseqr, as
             ! calculated below. hswork is computed assuming ilo=1 and ihi=n,
             ! the worst case.)
           if( info==0_${ik}$ ) then
              if( n==0_${ik}$ ) then
                 minwrk = 1_${ik}$
                 maxwrk = 1_${ik}$
              else
                 maxwrk = 2_${ik}$*n + n*stdlib${ii}$_ilaenv( 1_${ik}$, 'DGEHRD', ' ', n, 1_${ik}$, n, 0_${ik}$ )
                 if( wantvl ) then
                    minwrk = 4_${ik}$*n
                    maxwrk = max( maxwrk, 2_${ik}$*n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$,'DORGHR', ' ', n, 1_${ik}$, n,&
                               -1_${ik}$ ) )
                    call stdlib${ii}$_${ri}$hseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vl, ldvl,work, -1_${ik}$, &
                              info )
                    hswork = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + 1_${ik}$, n + hswork )
                    call stdlib${ii}$_${ri}$trevc3( 'L', 'B', select, n, a, lda,vl, ldvl, vr, ldvr, n, nout,&
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    maxwrk = max( maxwrk, 4_${ik}$*n )
                 else if( wantvr ) then
                    minwrk = 4_${ik}$*n
                    maxwrk = max( maxwrk, 2_${ik}$*n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$,'DORGHR', ' ', n, 1_${ik}$, n,&
                               -1_${ik}$ ) )
                    call stdlib${ii}$_${ri}$hseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vr, ldvr,work, -1_${ik}$, &
                              info )
                    hswork = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + 1_${ik}$, n + hswork )
                    call stdlib${ii}$_${ri}$trevc3( 'R', 'B', select, n, a, lda,vl, ldvl, vr, ldvr, n, nout,&
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    maxwrk = max( maxwrk, 4_${ik}$*n )
                 else
                    minwrk = 3_${ik}$*n
                    call stdlib${ii}$_${ri}$hseqr( 'E', 'N', n, 1_${ik}$, n, a, lda, wr, wi, vr, ldvr,work, -1_${ik}$, &
                              info )
                    hswork = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + 1_${ik}$, n + hswork )
                 end if
                 maxwrk = max( maxwrk, minwrk )
              end if
              work( 1_${ik}$ ) = maxwrk
              if( lwork<minwrk .and. .not.lquery ) then
                 info = -13_${ik}$
              end if
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'DGEEV ', -info )
              return
           else if( lquery ) then
              return
           end if
           ! quick return if possible
           if( n==0 )return
           ! get machine constants
           eps = stdlib${ii}$_${ri}$lamch( 'P' )
           smlnum = stdlib${ii}$_${ri}$lamch( 'S' )
           bignum = one / smlnum
           call stdlib${ii}$_${ri}$labad( smlnum, bignum )
           smlnum = sqrt( smlnum ) / eps
           bignum = one / smlnum
           ! scale a if max element outside range [smlnum,bignum]
           anrm = stdlib${ii}$_${ri}$lange( 'M', n, n, a, lda, dum )
           scalea = .false.
           if( anrm>zero .and. anrm<smlnum ) then
              scalea = .true.
              cscale = smlnum
           else if( anrm>bignum ) then
              scalea = .true.
              cscale = bignum
           end if
           if( scalea )call stdlib${ii}$_${ri}$lascl( 'G', 0_${ik}$, 0_${ik}$, anrm, cscale, n, n, a, lda, ierr )
           ! balance the matrix
           ! (workspace: need n)
           ibal = 1_${ik}$
           call stdlib${ii}$_${ri}$gebal( 'B', n, a, lda, ilo, ihi, work( ibal ), ierr )
           ! reduce to upper hessenberg form
           ! (workspace: need 3*n, prefer 2*n+n*nb)
           itau = ibal + n
           iwrk = itau + n
           call stdlib${ii}$_${ri}$gehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),lwork-iwrk+1, ierr &
                     )
           if( wantvl ) then
              ! want left eigenvectors
              ! copy householder vectors to vl
              side = 'L'
              call stdlib${ii}$_${ri}$lacpy( 'L', n, n, a, lda, vl, ldvl )
              ! generate orthogonal matrix in vl
              ! (workspace: need 3*n-1, prefer 2*n+(n-1)*nb)
              call stdlib${ii}$_${ri}$orghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vl
              ! (workspace: need n+1, prefer n+hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_${ri}$hseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vl, ldvl,work( iwrk ), &
                        lwork-iwrk+1, info )
              if( wantvr ) then
                 ! want left and right eigenvectors
                 ! copy schur vectors to vr
                 side = 'B'
                 call stdlib${ii}$_${ri}$lacpy( 'F', n, n, vl, ldvl, vr, ldvr )
              end if
           else if( wantvr ) then
              ! want right eigenvectors
              ! copy householder vectors to vr
              side = 'R'
              call stdlib${ii}$_${ri}$lacpy( 'L', n, n, a, lda, vr, ldvr )
              ! generate orthogonal matrix in vr
              ! (workspace: need 3*n-1, prefer 2*n+(n-1)*nb)
              call stdlib${ii}$_${ri}$orghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vr
              ! (workspace: need n+1, prefer n+hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_${ri}$hseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           else
              ! compute eigenvalues only
              ! (workspace: need n+1, prefer n+hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_${ri}$hseqr( 'E', 'N', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           end if
           ! if info /= 0 from stdlib${ii}$_${ri}$hseqr, then quit
           if( info/=0 )go to 50
           if( wantvl .or. wantvr ) then
              ! compute left and/or right eigenvectors
              ! (workspace: need 4*n, prefer n + n + 2*n*nb)
              call stdlib${ii}$_${ri}$trevc3( side, 'B', select, n, a, lda, vl, ldvl, vr, ldvr,n, nout, work(&
                         iwrk ), lwork-iwrk+1, ierr )
           end if
           if( wantvl ) then
              ! undo balancing of left eigenvectors
              ! (workspace: need n)
              call stdlib${ii}$_${ri}$gebak( 'B', 'L', n, ilo, ihi, work( ibal ), n, vl, ldvl,ierr )
              ! normalize left eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_${ri}$nrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_${ri}$scal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_${ri}$lapy2( stdlib${ii}$_${ri}$nrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_${ri}$nrm2( n, &
                              vl( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_${ri}$scal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_${ri}$scal( n, scl, vl( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( iwrk+k-1 ) = vl( k, i )**2_${ik}$ + vl( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_i${ri}$amax( n, work( iwrk ), 1_${ik}$ )
                    call stdlib${ii}$_${ri}$lartg( vl( k, i ), vl( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_${ri}$rot( n, vl( 1_${ik}$, i ), 1_${ik}$, vl( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vl( k, i+1 ) = zero
                 end if
              end do
           end if
           if( wantvr ) then
              ! undo balancing of right eigenvectors
              ! (workspace: need n)
              call stdlib${ii}$_${ri}$gebak( 'B', 'R', n, ilo, ihi, work( ibal ), n, vr, ldvr,ierr )
              ! normalize right eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_${ri}$nrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_${ri}$scal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_${ri}$lapy2( stdlib${ii}$_${ri}$nrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_${ri}$nrm2( n, &
                              vr( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_${ri}$scal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_${ri}$scal( n, scl, vr( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( iwrk+k-1 ) = vr( k, i )**2_${ik}$ + vr( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_i${ri}$amax( n, work( iwrk ), 1_${ik}$ )
                    call stdlib${ii}$_${ri}$lartg( vr( k, i ), vr( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_${ri}$rot( n, vr( 1_${ik}$, i ), 1_${ik}$, vr( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vr( k, i+1 ) = zero
                 end if
              end do
           end if
           ! undo scaling if necessary
           50 continue
           if( scalea ) then
              call stdlib${ii}$_${ri}$lascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wr( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              call stdlib${ii}$_${ri}$lascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wi( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              if( info>0_${ik}$ ) then
                 call stdlib${ii}$_${ri}$lascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wr, n,ierr )
                 call stdlib${ii}$_${ri}$lascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wi, n,ierr )
              end if
           end if
           work( 1_${ik}$ ) = maxwrk
           return
     end subroutine stdlib${ii}$_${ri}$geev

#:endif
#:endfor

     module subroutine stdlib${ii}$_cgeev( jobvl, jobvr, n, a, lda, w, vl, ldvl, vr, ldvr,work, lwork, rwork, &
     !! CGEEV computes for an N-by-N complex nonsymmetric matrix A, the
     !! eigenvalues and, optionally, the left and/or right eigenvectors.
     !! The right eigenvector v(j) of A satisfies
     !! A * v(j) = lambda(j) * v(j)
     !! where lambda(j) is its eigenvalue.
     !! The left eigenvector u(j) of A satisfies
     !! u(j)**H * A = lambda(j) * u(j)**H
     !! where u(j)**H denotes the conjugate transpose of u(j).
     !! The computed eigenvectors are normalized to have Euclidean norm
     !! equal to 1 and largest component real.
               info )
        ! -- lapack driver routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: jobvl, jobvr
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: lda, ldvl, ldvr, lwork, n
           ! Array Arguments 
           real(sp), intent(out) :: rwork(*)
           complex(sp), intent(inout) :: a(lda,*)
           complex(sp), intent(out) :: vl(ldvl,*), vr(ldvr,*), w(*), work(*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lquery, scalea, wantvl, wantvr
           character :: side
           integer(${ik}$) :: hswork, i, ibal, ierr, ihi, ilo, irwork, itau, iwrk, k, lwork_trevc, &
                     maxwrk, minwrk, nout
           real(sp) :: anrm, bignum, cscale, eps, scl, smlnum
           complex(sp) :: tmp
           ! Local Arrays 
           logical(lk) :: select(1_${ik}$)
           real(sp) :: dum(1_${ik}$)
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input arguments
           info = 0_${ik}$
           lquery = ( lwork==-1_${ik}$ )
           wantvl = stdlib_lsame( jobvl, 'V' )
           wantvr = stdlib_lsame( jobvr, 'V' )
           if( ( .not.wantvl ) .and. ( .not.stdlib_lsame( jobvl, 'N' ) ) ) then
              info = -1_${ik}$
           else if( ( .not.wantvr ) .and. ( .not.stdlib_lsame( jobvr, 'N' ) ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           else if( lda<max( 1_${ik}$, n ) ) then
              info = -5_${ik}$
           else if( ldvl<1_${ik}$ .or. ( wantvl .and. ldvl<n ) ) then
              info = -8_${ik}$
           else if( ldvr<1_${ik}$ .or. ( wantvr .and. ldvr<n ) ) then
              info = -10_${ik}$
           end if
           ! compute workspace
            ! (note: comments in the code beginning "workspace:" describe the
             ! minimal amount of workspace needed at that point in the code,
             ! as well as the preferred amount for good performance.
             ! cworkspace refers to complex workspace, and rworkspace to real
             ! workspace. nb refers to the optimal block size for the
             ! immediately following subroutine, as returned by stdlib${ii}$_ilaenv.
             ! hswork refers to the workspace preferred by stdlib${ii}$_chseqr, as
             ! calculated below. hswork is computed assuming ilo=1 and ihi=n,
             ! the worst case.)
           if( info==0_${ik}$ ) then
              if( n==0_${ik}$ ) then
                 minwrk = 1_${ik}$
                 maxwrk = 1_${ik}$
              else
                 maxwrk = n + n*stdlib${ii}$_ilaenv( 1_${ik}$, 'CGEHRD', ' ', n, 1_${ik}$, n, 0_${ik}$ )
                 minwrk = 2_${ik}$*n
                 if( wantvl ) then
                    maxwrk = max( maxwrk, n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$, 'CUNGHR',' ', n, 1_${ik}$, n, -&
                              1_${ik}$ ) )
                    call stdlib${ii}$_ctrevc3( 'L', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, rwork, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_chseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, w, vl, ldvl,work, -1_${ik}$, info )
                              
                 else if( wantvr ) then
                    maxwrk = max( maxwrk, n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$, 'CUNGHR',' ', n, 1_${ik}$, n, -&
                              1_${ik}$ ) )
                    call stdlib${ii}$_ctrevc3( 'R', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, rwork, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_chseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, w, vr, ldvr,work, -1_${ik}$, info )
                              
                 else
                    call stdlib${ii}$_chseqr( 'E', 'N', n, 1_${ik}$, n, a, lda, w, vr, ldvr,work, -1_${ik}$, info )
                              
                 end if
                 hswork = int( work(1_${ik}$),KIND=${ik}$)
                 maxwrk = max( maxwrk, hswork, minwrk )
              end if
              work( 1_${ik}$ ) = maxwrk
              if( lwork<minwrk .and. .not.lquery ) then
                 info = -12_${ik}$
              end if
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'CGEEV ', -info )
              return
           else if( lquery ) then
              return
           end if
           ! quick return if possible
           if( n==0 )return
           ! get machine constants
           eps = stdlib${ii}$_slamch( 'P' )
           smlnum = stdlib${ii}$_slamch( 'S' )
           bignum = one / smlnum
           call stdlib${ii}$_slabad( smlnum, bignum )
           smlnum = sqrt( smlnum ) / eps
           bignum = one / smlnum
           ! scale a if max element outside range [smlnum,bignum]
           anrm = stdlib${ii}$_clange( 'M', n, n, a, lda, dum )
           scalea = .false.
           if( anrm>zero .and. anrm<smlnum ) then
              scalea = .true.
              cscale = smlnum
           else if( anrm>bignum ) then
              scalea = .true.
              cscale = bignum
           end if
           if( scalea )call stdlib${ii}$_clascl( 'G', 0_${ik}$, 0_${ik}$, anrm, cscale, n, n, a, lda, ierr )
           ! balance the matrix
           ! (cworkspace: none)
           ! (rworkspace: need n)
           ibal = 1_${ik}$
           call stdlib${ii}$_cgebal( 'B', n, a, lda, ilo, ihi, rwork( ibal ), ierr )
           ! reduce to upper hessenberg form
           ! (cworkspace: need 2*n, prefer n+n*nb)
           ! (rworkspace: none)
           itau = 1_${ik}$
           iwrk = itau + n
           call stdlib${ii}$_cgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),lwork-iwrk+1, ierr &
                     )
           if( wantvl ) then
              ! want left eigenvectors
              ! copy householder vectors to vl
              side = 'L'
              call stdlib${ii}$_clacpy( 'L', n, n, a, lda, vl, ldvl )
              ! generate unitary matrix in vl
              ! (cworkspace: need 2*n-1, prefer n+(n-1)*nb)
              ! (rworkspace: none)
              call stdlib${ii}$_cunghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vl
              ! (cworkspace: need 1, prefer hswork (see comments) )
              ! (rworkspace: none)
              iwrk = itau
              call stdlib${ii}$_chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vl, ldvl,work( iwrk ), lwork-&
                        iwrk+1, info )
              if( wantvr ) then
                 ! want left and right eigenvectors
                 ! copy schur vectors to vr
                 side = 'B'
                 call stdlib${ii}$_clacpy( 'F', n, n, vl, ldvl, vr, ldvr )
              end if
           else if( wantvr ) then
              ! want right eigenvectors
              ! copy householder vectors to vr
              side = 'R'
              call stdlib${ii}$_clacpy( 'L', n, n, a, lda, vr, ldvr )
              ! generate unitary matrix in vr
              ! (cworkspace: need 2*n-1, prefer n+(n-1)*nb)
              ! (rworkspace: none)
              call stdlib${ii}$_cunghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vr
              ! (cworkspace: need 1, prefer hswork (see comments) )
              ! (rworkspace: none)
              iwrk = itau
              call stdlib${ii}$_chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vr, ldvr,work( iwrk ), lwork-&
                        iwrk+1, info )
           else
              ! compute eigenvalues only
              ! (cworkspace: need 1, prefer hswork (see comments) )
              ! (rworkspace: none)
              iwrk = itau
              call stdlib${ii}$_chseqr( 'E', 'N', n, ilo, ihi, a, lda, w, vr, ldvr,work( iwrk ), lwork-&
                        iwrk+1, info )
           end if
           ! if info /= 0 from stdlib${ii}$_chseqr, then quit
           if( info/=0 )go to 50
           if( wantvl .or. wantvr ) then
              ! compute left and/or right eigenvectors
              ! (cworkspace: need 2*n, prefer n + 2*n*nb)
              ! (rworkspace: need 2*n)
              irwork = ibal + n
              call stdlib${ii}$_ctrevc3( side, 'B', select, n, a, lda, vl, ldvl, vr, ldvr,n, nout, work(&
                         iwrk ), lwork-iwrk+1,rwork( irwork ), n, ierr )
           end if
           if( wantvl ) then
              ! undo balancing of left eigenvectors
              ! (cworkspace: none)
              ! (rworkspace: need n)
              call stdlib${ii}$_cgebak( 'B', 'L', n, ilo, ihi, rwork( ibal ), n, vl, ldvl,ierr )
              ! normalize left eigenvectors and make largest component real
              do i = 1, n
                 scl = one / stdlib${ii}$_scnrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ )
                 call stdlib${ii}$_csscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                 do k = 1, n
                    rwork( irwork+k-1 ) = real( vl( k, i ),KIND=sp)**2_${ik}$ +aimag( vl( k, i ) )&
                              **2_${ik}$
                 end do
                 k = stdlib${ii}$_isamax( n, rwork( irwork ), 1_${ik}$ )
                 tmp = conjg( vl( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
                 call stdlib${ii}$_cscal( n, tmp, vl( 1_${ik}$, i ), 1_${ik}$ )
                 vl( k, i ) = cmplx( real( vl( k, i ),KIND=sp), zero,KIND=sp)
              end do
           end if
           if( wantvr ) then
              ! undo balancing of right eigenvectors
              ! (cworkspace: none)
              ! (rworkspace: need n)
              call stdlib${ii}$_cgebak( 'B', 'R', n, ilo, ihi, rwork( ibal ), n, vr, ldvr,ierr )
              ! normalize right eigenvectors and make largest component real
              do i = 1, n
                 scl = one / stdlib${ii}$_scnrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ )
                 call stdlib${ii}$_csscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                 do k = 1, n
                    rwork( irwork+k-1 ) = real( vr( k, i ),KIND=sp)**2_${ik}$ +aimag( vr( k, i ) )&
                              **2_${ik}$
                 end do
                 k = stdlib${ii}$_isamax( n, rwork( irwork ), 1_${ik}$ )
                 tmp = conjg( vr( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
                 call stdlib${ii}$_cscal( n, tmp, vr( 1_${ik}$, i ), 1_${ik}$ )
                 vr( k, i ) = cmplx( real( vr( k, i ),KIND=sp), zero,KIND=sp)
              end do
           end if
           ! undo scaling if necessary
           50 continue
           if( scalea ) then
              call stdlib${ii}$_clascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, w( info+1 ),max( n-info, 1_${ik}$ )&
                        , ierr )
              if( info>0_${ik}$ ) then
                 call stdlib${ii}$_clascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, w, n, ierr )
              end if
           end if
           work( 1_${ik}$ ) = maxwrk
           return
     end subroutine stdlib${ii}$_cgeev

     module subroutine stdlib${ii}$_zgeev( jobvl, jobvr, n, a, lda, w, vl, ldvl, vr, ldvr,work, lwork, rwork, &
     !! ZGEEV computes for an N-by-N complex nonsymmetric matrix A, the
     !! eigenvalues and, optionally, the left and/or right eigenvectors.
     !! The right eigenvector v(j) of A satisfies
     !! A * v(j) = lambda(j) * v(j)
     !! where lambda(j) is its eigenvalue.
     !! The left eigenvector u(j) of A satisfies
     !! u(j)**H * A = lambda(j) * u(j)**H
     !! where u(j)**H denotes the conjugate transpose of u(j).
     !! The computed eigenvectors are normalized to have Euclidean norm
     !! equal to 1 and largest component real.
               info )
        ! -- lapack driver routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: jobvl, jobvr
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: lda, ldvl, ldvr, lwork, n
           ! Array Arguments 
           real(dp), intent(out) :: rwork(*)
           complex(dp), intent(inout) :: a(lda,*)
           complex(dp), intent(out) :: vl(ldvl,*), vr(ldvr,*), w(*), work(*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lquery, scalea, wantvl, wantvr
           character :: side
           integer(${ik}$) :: hswork, i, ibal, ierr, ihi, ilo, irwork, itau, iwrk, k, lwork_trevc, &
                     maxwrk, minwrk, nout
           real(dp) :: anrm, bignum, cscale, eps, scl, smlnum
           complex(dp) :: tmp
           ! Local Arrays 
           logical(lk) :: select(1_${ik}$)
           real(dp) :: dum(1_${ik}$)
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input arguments
           info = 0_${ik}$
           lquery = ( lwork==-1_${ik}$ )
           wantvl = stdlib_lsame( jobvl, 'V' )
           wantvr = stdlib_lsame( jobvr, 'V' )
           if( ( .not.wantvl ) .and. ( .not.stdlib_lsame( jobvl, 'N' ) ) ) then
              info = -1_${ik}$
           else if( ( .not.wantvr ) .and. ( .not.stdlib_lsame( jobvr, 'N' ) ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           else if( lda<max( 1_${ik}$, n ) ) then
              info = -5_${ik}$
           else if( ldvl<1_${ik}$ .or. ( wantvl .and. ldvl<n ) ) then
              info = -8_${ik}$
           else if( ldvr<1_${ik}$ .or. ( wantvr .and. ldvr<n ) ) then
              info = -10_${ik}$
           end if
           ! compute workspace
            ! (note: comments in the code beginning "workspace:" describe the
             ! minimal amount of workspace needed at that point in the code,
             ! as well as the preferred amount for good performance.
             ! cworkspace refers to complex workspace, and rworkspace to real
             ! workspace. nb refers to the optimal block size for the
             ! immediately following subroutine, as returned by stdlib${ii}$_ilaenv.
             ! hswork refers to the workspace preferred by stdlib${ii}$_zhseqr, as
             ! calculated below. hswork is computed assuming ilo=1 and ihi=n,
             ! the worst case.)
           if( info==0_${ik}$ ) then
              if( n==0_${ik}$ ) then
                 minwrk = 1_${ik}$
                 maxwrk = 1_${ik}$
              else
                 maxwrk = n + n*stdlib${ii}$_ilaenv( 1_${ik}$, 'ZGEHRD', ' ', n, 1_${ik}$, n, 0_${ik}$ )
                 minwrk = 2_${ik}$*n
                 if( wantvl ) then
                    maxwrk = max( maxwrk, n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$, 'ZUNGHR',' ', n, 1_${ik}$, n, -&
                              1_${ik}$ ) )
                    call stdlib${ii}$_ztrevc3( 'L', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, rwork, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_zhseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, w, vl, ldvl,work, -1_${ik}$, info )
                              
                 else if( wantvr ) then
                    maxwrk = max( maxwrk, n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$, 'ZUNGHR',' ', n, 1_${ik}$, n, -&
                              1_${ik}$ ) )
                    call stdlib${ii}$_ztrevc3( 'R', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, rwork, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_zhseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, w, vr, ldvr,work, -1_${ik}$, info )
                              
                 else
                    call stdlib${ii}$_zhseqr( 'E', 'N', n, 1_${ik}$, n, a, lda, w, vr, ldvr,work, -1_${ik}$, info )
                              
                 end if
                 hswork = int( work(1_${ik}$),KIND=${ik}$)
                 maxwrk = max( maxwrk, hswork, minwrk )
              end if
              work( 1_${ik}$ ) = maxwrk
              if( lwork<minwrk .and. .not.lquery ) then
                 info = -12_${ik}$
              end if
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'ZGEEV ', -info )
              return
           else if( lquery ) then
              return
           end if
           ! quick return if possible
           if( n==0 )return
           ! get machine constants
           eps = stdlib${ii}$_dlamch( 'P' )
           smlnum = stdlib${ii}$_dlamch( 'S' )
           bignum = one / smlnum
           call stdlib${ii}$_dlabad( smlnum, bignum )
           smlnum = sqrt( smlnum ) / eps
           bignum = one / smlnum
           ! scale a if max element outside range [smlnum,bignum]
           anrm = stdlib${ii}$_zlange( 'M', n, n, a, lda, dum )
           scalea = .false.
           if( anrm>zero .and. anrm<smlnum ) then
              scalea = .true.
              cscale = smlnum
           else if( anrm>bignum ) then
              scalea = .true.
              cscale = bignum
           end if
           if( scalea )call stdlib${ii}$_zlascl( 'G', 0_${ik}$, 0_${ik}$, anrm, cscale, n, n, a, lda, ierr )
           ! balance the matrix
           ! (cworkspace: none)
           ! (rworkspace: need n)
           ibal = 1_${ik}$
           call stdlib${ii}$_zgebal( 'B', n, a, lda, ilo, ihi, rwork( ibal ), ierr )
           ! reduce to upper hessenberg form
           ! (cworkspace: need 2*n, prefer n+n*nb)
           ! (rworkspace: none)
           itau = 1_${ik}$
           iwrk = itau + n
           call stdlib${ii}$_zgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),lwork-iwrk+1, ierr &
                     )
           if( wantvl ) then
              ! want left eigenvectors
              ! copy householder vectors to vl
              side = 'L'
              call stdlib${ii}$_zlacpy( 'L', n, n, a, lda, vl, ldvl )
              ! generate unitary matrix in vl
              ! (cworkspace: need 2*n-1, prefer n+(n-1)*nb)
              ! (rworkspace: none)
              call stdlib${ii}$_zunghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vl
              ! (cworkspace: need 1, prefer hswork (see comments) )
              ! (rworkspace: none)
              iwrk = itau
              call stdlib${ii}$_zhseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vl, ldvl,work( iwrk ), lwork-&
                        iwrk+1, info )
              if( wantvr ) then
                 ! want left and right eigenvectors
                 ! copy schur vectors to vr
                 side = 'B'
                 call stdlib${ii}$_zlacpy( 'F', n, n, vl, ldvl, vr, ldvr )
              end if
           else if( wantvr ) then
              ! want right eigenvectors
              ! copy householder vectors to vr
              side = 'R'
              call stdlib${ii}$_zlacpy( 'L', n, n, a, lda, vr, ldvr )
              ! generate unitary matrix in vr
              ! (cworkspace: need 2*n-1, prefer n+(n-1)*nb)
              ! (rworkspace: none)
              call stdlib${ii}$_zunghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vr
              ! (cworkspace: need 1, prefer hswork (see comments) )
              ! (rworkspace: none)
              iwrk = itau
              call stdlib${ii}$_zhseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vr, ldvr,work( iwrk ), lwork-&
                        iwrk+1, info )
           else
              ! compute eigenvalues only
              ! (cworkspace: need 1, prefer hswork (see comments) )
              ! (rworkspace: none)
              iwrk = itau
              call stdlib${ii}$_zhseqr( 'E', 'N', n, ilo, ihi, a, lda, w, vr, ldvr,work( iwrk ), lwork-&
                        iwrk+1, info )
           end if
           ! if info /= 0 from stdlib${ii}$_zhseqr, then quit
           if( info/=0 )go to 50
           if( wantvl .or. wantvr ) then
              ! compute left and/or right eigenvectors
              ! (cworkspace: need 2*n, prefer n + 2*n*nb)
              ! (rworkspace: need 2*n)
              irwork = ibal + n
              call stdlib${ii}$_ztrevc3( side, 'B', select, n, a, lda, vl, ldvl, vr, ldvr,n, nout, work(&
                         iwrk ), lwork-iwrk+1,rwork( irwork ), n, ierr )
           end if
           if( wantvl ) then
              ! undo balancing of left eigenvectors
              ! (cworkspace: none)
              ! (rworkspace: need n)
              call stdlib${ii}$_zgebak( 'B', 'L', n, ilo, ihi, rwork( ibal ), n, vl, ldvl,ierr )
              ! normalize left eigenvectors and make largest component real
              do i = 1, n
                 scl = one / stdlib${ii}$_dznrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ )
                 call stdlib${ii}$_zdscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                 do k = 1, n
                    rwork( irwork+k-1 ) = real( vl( k, i ),KIND=dp)**2_${ik}$ +aimag( vl( k, i ) )&
                              **2_${ik}$
                 end do
                 k = stdlib${ii}$_idamax( n, rwork( irwork ), 1_${ik}$ )
                 tmp = conjg( vl( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
                 call stdlib${ii}$_zscal( n, tmp, vl( 1_${ik}$, i ), 1_${ik}$ )
                 vl( k, i ) = cmplx( real( vl( k, i ),KIND=dp), zero,KIND=dp)
              end do
           end if
           if( wantvr ) then
              ! undo balancing of right eigenvectors
              ! (cworkspace: none)
              ! (rworkspace: need n)
              call stdlib${ii}$_zgebak( 'B', 'R', n, ilo, ihi, rwork( ibal ), n, vr, ldvr,ierr )
              ! normalize right eigenvectors and make largest component real
              do i = 1, n
                 scl = one / stdlib${ii}$_dznrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ )
                 call stdlib${ii}$_zdscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                 do k = 1, n
                    rwork( irwork+k-1 ) = real( vr( k, i ),KIND=dp)**2_${ik}$ +aimag( vr( k, i ) )&
                              **2_${ik}$
                 end do
                 k = stdlib${ii}$_idamax( n, rwork( irwork ), 1_${ik}$ )
                 tmp = conjg( vr( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
                 call stdlib${ii}$_zscal( n, tmp, vr( 1_${ik}$, i ), 1_${ik}$ )
                 vr( k, i ) = cmplx( real( vr( k, i ),KIND=dp), zero,KIND=dp)
              end do
           end if
           ! undo scaling if necessary
           50 continue
           if( scalea ) then
              call stdlib${ii}$_zlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, w( info+1 ),max( n-info, 1_${ik}$ )&
                        , ierr )
              if( info>0_${ik}$ ) then
                 call stdlib${ii}$_zlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, w, n, ierr )
              end if
           end if
           work( 1_${ik}$ ) = maxwrk
           return
     end subroutine stdlib${ii}$_zgeev

#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
     module subroutine stdlib${ii}$_${ci}$geev( jobvl, jobvr, n, a, lda, w, vl, ldvl, vr, ldvr,work, lwork, rwork, &
     !! ZGEEV: computes for an N-by-N complex nonsymmetric matrix A, the
     !! eigenvalues and, optionally, the left and/or right eigenvectors.
     !! The right eigenvector v(j) of A satisfies
     !! A * v(j) = lambda(j) * v(j)
     !! where lambda(j) is its eigenvalue.
     !! The left eigenvector u(j) of A satisfies
     !! u(j)**H * A = lambda(j) * u(j)**H
     !! where u(j)**H denotes the conjugate transpose of u(j).
     !! The computed eigenvectors are normalized to have Euclidean norm
     !! equal to 1 and largest component real.
               info )
        ! -- lapack driver routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_${ck}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: jobvl, jobvr
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: lda, ldvl, ldvr, lwork, n
           ! Array Arguments 
           real(${ck}$), intent(out) :: rwork(*)
           complex(${ck}$), intent(inout) :: a(lda,*)
           complex(${ck}$), intent(out) :: vl(ldvl,*), vr(ldvr,*), w(*), work(*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lquery, scalea, wantvl, wantvr
           character :: side
           integer(${ik}$) :: hswork, i, ibal, ierr, ihi, ilo, irwork, itau, iwrk, k, lwork_trevc, &
                     maxwrk, minwrk, nout
           real(${ck}$) :: anrm, bignum, cscale, eps, scl, smlnum
           complex(${ck}$) :: tmp
           ! Local Arrays 
           logical(lk) :: select(1_${ik}$)
           real(${ck}$) :: dum(1_${ik}$)
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input arguments
           info = 0_${ik}$
           lquery = ( lwork==-1_${ik}$ )
           wantvl = stdlib_lsame( jobvl, 'V' )
           wantvr = stdlib_lsame( jobvr, 'V' )
           if( ( .not.wantvl ) .and. ( .not.stdlib_lsame( jobvl, 'N' ) ) ) then
              info = -1_${ik}$
           else if( ( .not.wantvr ) .and. ( .not.stdlib_lsame( jobvr, 'N' ) ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           else if( lda<max( 1_${ik}$, n ) ) then
              info = -5_${ik}$
           else if( ldvl<1_${ik}$ .or. ( wantvl .and. ldvl<n ) ) then
              info = -8_${ik}$
           else if( ldvr<1_${ik}$ .or. ( wantvr .and. ldvr<n ) ) then
              info = -10_${ik}$
           end if
           ! compute workspace
            ! (note: comments in the code beginning "workspace:" describe the
             ! minimal amount of workspace needed at that point in the code,
             ! as well as the preferred amount for good performance.
             ! cworkspace refers to complex workspace, and rworkspace to real
             ! workspace. nb refers to the optimal block size for the
             ! immediately following subroutine, as returned by stdlib${ii}$_ilaenv.
             ! hswork refers to the workspace preferred by stdlib${ii}$_${ci}$hseqr, as
             ! calculated below. hswork is computed assuming ilo=1 and ihi=n,
             ! the worst case.)
           if( info==0_${ik}$ ) then
              if( n==0_${ik}$ ) then
                 minwrk = 1_${ik}$
                 maxwrk = 1_${ik}$
              else
                 maxwrk = n + n*stdlib${ii}$_ilaenv( 1_${ik}$, 'ZGEHRD', ' ', n, 1_${ik}$, n, 0_${ik}$ )
                 minwrk = 2_${ik}$*n
                 if( wantvl ) then
                    maxwrk = max( maxwrk, n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$, 'ZUNGHR',' ', n, 1_${ik}$, n, -&
                              1_${ik}$ ) )
                    call stdlib${ii}$_${ci}$trevc3( 'L', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, rwork, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_${ci}$hseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, w, vl, ldvl,work, -1_${ik}$, info )
                              
                 else if( wantvr ) then
                    maxwrk = max( maxwrk, n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$, 'ZUNGHR',' ', n, 1_${ik}$, n, -&
                              1_${ik}$ ) )
                    call stdlib${ii}$_${ci}$trevc3( 'R', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, rwork, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_${ci}$hseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, w, vr, ldvr,work, -1_${ik}$, info )
                              
                 else
                    call stdlib${ii}$_${ci}$hseqr( 'E', 'N', n, 1_${ik}$, n, a, lda, w, vr, ldvr,work, -1_${ik}$, info )
                              
                 end if
                 hswork = int( work(1_${ik}$),KIND=${ik}$)
                 maxwrk = max( maxwrk, hswork, minwrk )
              end if
              work( 1_${ik}$ ) = maxwrk
              if( lwork<minwrk .and. .not.lquery ) then
                 info = -12_${ik}$
              end if
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'ZGEEV ', -info )
              return
           else if( lquery ) then
              return
           end if
           ! quick return if possible
           if( n==0 )return
           ! get machine constants
           eps = stdlib${ii}$_${c2ri(ci)}$lamch( 'P' )
           smlnum = stdlib${ii}$_${c2ri(ci)}$lamch( 'S' )
           bignum = one / smlnum
           call stdlib${ii}$_${c2ri(ci)}$labad( smlnum, bignum )
           smlnum = sqrt( smlnum ) / eps
           bignum = one / smlnum
           ! scale a if max element outside range [smlnum,bignum]
           anrm = stdlib${ii}$_${ci}$lange( 'M', n, n, a, lda, dum )
           scalea = .false.
           if( anrm>zero .and. anrm<smlnum ) then
              scalea = .true.
              cscale = smlnum
           else if( anrm>bignum ) then
              scalea = .true.
              cscale = bignum
           end if
           if( scalea )call stdlib${ii}$_${ci}$lascl( 'G', 0_${ik}$, 0_${ik}$, anrm, cscale, n, n, a, lda, ierr )
           ! balance the matrix
           ! (cworkspace: none)
           ! (rworkspace: need n)
           ibal = 1_${ik}$
           call stdlib${ii}$_${ci}$gebal( 'B', n, a, lda, ilo, ihi, rwork( ibal ), ierr )
           ! reduce to upper hessenberg form
           ! (cworkspace: need 2*n, prefer n+n*nb)
           ! (rworkspace: none)
           itau = 1_${ik}$
           iwrk = itau + n
           call stdlib${ii}$_${ci}$gehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),lwork-iwrk+1, ierr &
                     )
           if( wantvl ) then
              ! want left eigenvectors
              ! copy householder vectors to vl
              side = 'L'
              call stdlib${ii}$_${ci}$lacpy( 'L', n, n, a, lda, vl, ldvl )
              ! generate unitary matrix in vl
              ! (cworkspace: need 2*n-1, prefer n+(n-1)*nb)
              ! (rworkspace: none)
              call stdlib${ii}$_${ci}$unghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vl
              ! (cworkspace: need 1, prefer hswork (see comments) )
              ! (rworkspace: none)
              iwrk = itau
              call stdlib${ii}$_${ci}$hseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vl, ldvl,work( iwrk ), lwork-&
                        iwrk+1, info )
              if( wantvr ) then
                 ! want left and right eigenvectors
                 ! copy schur vectors to vr
                 side = 'B'
                 call stdlib${ii}$_${ci}$lacpy( 'F', n, n, vl, ldvl, vr, ldvr )
              end if
           else if( wantvr ) then
              ! want right eigenvectors
              ! copy householder vectors to vr
              side = 'R'
              call stdlib${ii}$_${ci}$lacpy( 'L', n, n, a, lda, vr, ldvr )
              ! generate unitary matrix in vr
              ! (cworkspace: need 2*n-1, prefer n+(n-1)*nb)
              ! (rworkspace: none)
              call stdlib${ii}$_${ci}$unghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vr
              ! (cworkspace: need 1, prefer hswork (see comments) )
              ! (rworkspace: none)
              iwrk = itau
              call stdlib${ii}$_${ci}$hseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vr, ldvr,work( iwrk ), lwork-&
                        iwrk+1, info )
           else
              ! compute eigenvalues only
              ! (cworkspace: need 1, prefer hswork (see comments) )
              ! (rworkspace: none)
              iwrk = itau
              call stdlib${ii}$_${ci}$hseqr( 'E', 'N', n, ilo, ihi, a, lda, w, vr, ldvr,work( iwrk ), lwork-&
                        iwrk+1, info )
           end if
           ! if info /= 0 from stdlib${ii}$_${ci}$hseqr, then quit
           if( info/=0 )go to 50
           if( wantvl .or. wantvr ) then
              ! compute left and/or right eigenvectors
              ! (cworkspace: need 2*n, prefer n + 2*n*nb)
              ! (rworkspace: need 2*n)
              irwork = ibal + n
              call stdlib${ii}$_${ci}$trevc3( side, 'B', select, n, a, lda, vl, ldvl, vr, ldvr,n, nout, work(&
                         iwrk ), lwork-iwrk+1,rwork( irwork ), n, ierr )
           end if
           if( wantvl ) then
              ! undo balancing of left eigenvectors
              ! (cworkspace: none)
              ! (rworkspace: need n)
              call stdlib${ii}$_${ci}$gebak( 'B', 'L', n, ilo, ihi, rwork( ibal ), n, vl, ldvl,ierr )
              ! normalize left eigenvectors and make largest component real
              do i = 1, n
                 scl = one / stdlib${ii}$_${c2ri(ci)}$znrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ )
                 call stdlib${ii}$_${ci}$dscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                 do k = 1, n
                    rwork( irwork+k-1 ) = real( vl( k, i ),KIND=${ck}$)**2_${ik}$ +aimag( vl( k, i ) )&
                              **2_${ik}$
                 end do
                 k = stdlib${ii}$_i${c2ri(ci)}$amax( n, rwork( irwork ), 1_${ik}$ )
                 tmp = conjg( vl( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
                 call stdlib${ii}$_${ci}$scal( n, tmp, vl( 1_${ik}$, i ), 1_${ik}$ )
                 vl( k, i ) = cmplx( real( vl( k, i ),KIND=${ck}$), zero,KIND=${ck}$)
              end do
           end if
           if( wantvr ) then
              ! undo balancing of right eigenvectors
              ! (cworkspace: none)
              ! (rworkspace: need n)
              call stdlib${ii}$_${ci}$gebak( 'B', 'R', n, ilo, ihi, rwork( ibal ), n, vr, ldvr,ierr )
              ! normalize right eigenvectors and make largest component real
              do i = 1, n
                 scl = one / stdlib${ii}$_${c2ri(ci)}$znrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ )
                 call stdlib${ii}$_${ci}$dscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                 do k = 1, n
                    rwork( irwork+k-1 ) = real( vr( k, i ),KIND=${ck}$)**2_${ik}$ +aimag( vr( k, i ) )&
                              **2_${ik}$
                 end do
                 k = stdlib${ii}$_i${c2ri(ci)}$amax( n, rwork( irwork ), 1_${ik}$ )
                 tmp = conjg( vr( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
                 call stdlib${ii}$_${ci}$scal( n, tmp, vr( 1_${ik}$, i ), 1_${ik}$ )
                 vr( k, i ) = cmplx( real( vr( k, i ),KIND=${ck}$), zero,KIND=${ck}$)
              end do
           end if
           ! undo scaling if necessary
           50 continue
           if( scalea ) then
              call stdlib${ii}$_${ci}$lascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, w( info+1 ),max( n-info, 1_${ik}$ )&
                        , ierr )
              if( info>0_${ik}$ ) then
                 call stdlib${ii}$_${ci}$lascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, w, n, ierr )
              end if
           end if
           work( 1_${ik}$ ) = maxwrk
           return
     end subroutine stdlib${ii}$_${ci}$geev

#:endif
#:endfor



     module subroutine stdlib${ii}$_sgeevx( balanc, jobvl, jobvr, sense, n, a, lda, wr, wi,vl, ldvl, vr, ldvr, &
     !! SGEEVX computes for an N-by-N real nonsymmetric matrix A, the
     !! eigenvalues and, optionally, the left and/or right eigenvectors.
     !! Optionally also, it computes a balancing transformation to improve
     !! the conditioning of the eigenvalues and eigenvectors (ILO, IHI,
     !! SCALE, and ABNRM), reciprocal condition numbers for the eigenvalues
     !! (RCONDE), and reciprocal condition numbers for the right
     !! eigenvectors (RCONDV).
     !! The right eigenvector v(j) of A satisfies
     !! A * v(j) = lambda(j) * v(j)
     !! where lambda(j) is its eigenvalue.
     !! The left eigenvector u(j) of A satisfies
     !! u(j)**H * A = lambda(j) * u(j)**H
     !! where u(j)**H denotes the conjugate-transpose of u(j).
     !! The computed eigenvectors are normalized to have Euclidean norm
     !! equal to 1 and largest component real.
     !! Balancing a matrix means permuting the rows and columns to make it
     !! more nearly upper triangular, and applying a diagonal similarity
     !! transformation D * A * D**(-1), where D is a diagonal matrix, to
     !! make its rows and columns closer in norm and the condition numbers
     !! of its eigenvalues and eigenvectors smaller.  The computed
     !! reciprocal condition numbers correspond to the balanced matrix.
     !! Permuting rows and columns will not change the condition numbers
     !! (in exact arithmetic) but diagonal scaling will.  For further
     !! explanation of balancing, see section 4.10.2_sp of the LAPACK
     !! Users' Guide.
               ilo, ihi, scale, abnrm,rconde, rcondv, work, lwork, iwork, info )
        ! -- lapack driver routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: balanc, jobvl, jobvr, sense
           integer(${ik}$), intent(out) :: ihi, ilo, info
           integer(${ik}$), intent(in) :: lda, ldvl, ldvr, lwork, n
           real(sp), intent(out) :: abnrm
           ! Array Arguments 
           integer(${ik}$), intent(out) :: iwork(*)
           real(sp), intent(inout) :: a(lda,*)
           real(sp), intent(out) :: rconde(*), rcondv(*), scale(*), vl(ldvl,*), vr(ldvr,*), wi(*),&
                      work(*), wr(*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lquery, scalea, wantvl, wantvr, wntsnb, wntsne, wntsnn, wntsnv
           character :: job, side
           integer(${ik}$) :: hswork, i, icond, ierr, itau, iwrk, k, lwork_trevc, maxwrk, minwrk, &
                     nout
           real(sp) :: anrm, bignum, cs, cscale, eps, r, scl, smlnum, sn
           ! Local Arrays 
           logical(lk) :: select(1_${ik}$)
           real(sp) :: dum(1_${ik}$)
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input arguments
           info = 0_${ik}$
           lquery = ( lwork==-1_${ik}$ )
           wantvl = stdlib_lsame( jobvl, 'V' )
           wantvr = stdlib_lsame( jobvr, 'V' )
           wntsnn = stdlib_lsame( sense, 'N' )
           wntsne = stdlib_lsame( sense, 'E' )
           wntsnv = stdlib_lsame( sense, 'V' )
           wntsnb = stdlib_lsame( sense, 'B' )
           if( .not.( stdlib_lsame( balanc, 'N' ) .or. stdlib_lsame( balanc, 'S' ).or. &
                     stdlib_lsame( balanc, 'P' ) .or. stdlib_lsame( balanc, 'B' ) ) )then
              info = -1_${ik}$
           else if( ( .not.wantvl ) .and. ( .not.stdlib_lsame( jobvl, 'N' ) ) ) then
              info = -2_${ik}$
           else if( ( .not.wantvr ) .and. ( .not.stdlib_lsame( jobvr, 'N' ) ) ) then
              info = -3_${ik}$
           else if( .not.( wntsnn .or. wntsne .or. wntsnb .or. wntsnv ) .or.( ( wntsne .or. &
                     wntsnb ) .and. .not.( wantvl .and.wantvr ) ) ) then
              info = -4_${ik}$
           else if( n<0_${ik}$ ) then
              info = -5_${ik}$
           else if( lda<max( 1_${ik}$, n ) ) then
              info = -7_${ik}$
           else if( ldvl<1_${ik}$ .or. ( wantvl .and. ldvl<n ) ) then
              info = -11_${ik}$
           else if( ldvr<1_${ik}$ .or. ( wantvr .and. ldvr<n ) ) then
              info = -13_${ik}$
           end if
           ! compute workspace
            ! (note: comments in the code beginning "workspace:" describe the
             ! minimal amount of workspace needed at that point in the code,
             ! as well as the preferred amount for good performance.
             ! nb refers to the optimal block size for the immediately
             ! following subroutine, as returned by stdlib${ii}$_ilaenv.
             ! hswork refers to the workspace preferred by stdlib${ii}$_shseqr, as
             ! calculated below. hswork is computed assuming ilo=1 and ihi=n,
             ! the worst case.)
           if( info==0_${ik}$ ) then
              if( n==0_${ik}$ ) then
                 minwrk = 1_${ik}$
                 maxwrk = 1_${ik}$
              else
                 maxwrk = n + n*stdlib${ii}$_ilaenv( 1_${ik}$, 'SGEHRD', ' ', n, 1_${ik}$, n, 0_${ik}$ )
                 if( wantvl ) then
                    call stdlib${ii}$_strevc3( 'L', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_shseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vl, ldvl,work, -1_${ik}$, &
                              info )
                 else if( wantvr ) then
                    call stdlib${ii}$_strevc3( 'R', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_shseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vr, ldvr,work, -1_${ik}$, &
                              info )
                 else
                    if( wntsnn ) then
                       call stdlib${ii}$_shseqr( 'E', 'N', n, 1_${ik}$, n, a, lda, wr, wi, vr,ldvr, work, -1_${ik}$, &
                                 info )
                    else
                       call stdlib${ii}$_shseqr( 'S', 'N', n, 1_${ik}$, n, a, lda, wr, wi, vr,ldvr, work, -1_${ik}$, &
                                 info )
                    end if
                 end if
                 hswork = int( work(1_${ik}$),KIND=${ik}$)
                 if( ( .not.wantvl ) .and. ( .not.wantvr ) ) then
                    minwrk = 2_${ik}$*n
                    if( .not.wntsnn )minwrk = max( minwrk, n*n+6*n )
                    maxwrk = max( maxwrk, hswork )
                    if( .not.wntsnn )maxwrk = max( maxwrk, n*n + 6_${ik}$*n )
                 else
                    minwrk = 3_${ik}$*n
                    if( ( .not.wntsnn ) .and. ( .not.wntsne ) )minwrk = max( minwrk, n*n + 6_${ik}$*n )
                              
                    maxwrk = max( maxwrk, hswork )
                    maxwrk = max( maxwrk, n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$, 'SORGHR',' ', n, 1_${ik}$, n, -&
                              1_${ik}$ ) )
                    if( ( .not.wntsnn ) .and. ( .not.wntsne ) )maxwrk = max( maxwrk, n*n + 6_${ik}$*n )
                              
                    maxwrk = max( maxwrk, 3_${ik}$*n )
                 end if
                 maxwrk = max( maxwrk, minwrk )
              end if
              work( 1_${ik}$ ) = maxwrk
              if( lwork<minwrk .and. .not.lquery ) then
                 info = -21_${ik}$
              end if
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'SGEEVX', -info )
              return
           else if( lquery ) then
              return
           end if
           ! quick return if possible
           if( n==0 )return
           ! get machine constants
           eps = stdlib${ii}$_slamch( 'P' )
           smlnum = stdlib${ii}$_slamch( 'S' )
           bignum = one / smlnum
           call stdlib${ii}$_slabad( smlnum, bignum )
           smlnum = sqrt( smlnum ) / eps
           bignum = one / smlnum
           ! scale a if max element outside range [smlnum,bignum]
           icond = 0_${ik}$
           anrm = stdlib${ii}$_slange( 'M', n, n, a, lda, dum )
           scalea = .false.
           if( anrm>zero .and. anrm<smlnum ) then
              scalea = .true.
              cscale = smlnum
           else if( anrm>bignum ) then
              scalea = .true.
              cscale = bignum
           end if
           if( scalea )call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, anrm, cscale, n, n, a, lda, ierr )
           ! balance the matrix and compute abnrm
           call stdlib${ii}$_sgebal( balanc, n, a, lda, ilo, ihi, scale, ierr )
           abnrm = stdlib${ii}$_slange( '1', n, n, a, lda, dum )
           if( scalea ) then
              dum( 1_${ik}$ ) = abnrm
              call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, 1_${ik}$, 1_${ik}$, dum, 1_${ik}$, ierr )
              abnrm = dum( 1_${ik}$ )
           end if
           ! reduce to upper hessenberg form
           ! (workspace: need 2*n, prefer n+n*nb)
           itau = 1_${ik}$
           iwrk = itau + n
           call stdlib${ii}$_sgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),lwork-iwrk+1, ierr &
                     )
           if( wantvl ) then
              ! want left eigenvectors
              ! copy householder vectors to vl
              side = 'L'
              call stdlib${ii}$_slacpy( 'L', n, n, a, lda, vl, ldvl )
              ! generate orthogonal matrix in vl
              ! (workspace: need 2*n-1, prefer n+(n-1)*nb)
              call stdlib${ii}$_sorghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vl
              ! (workspace: need 1, prefer hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_shseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vl, ldvl,work( iwrk ), &
                        lwork-iwrk+1, info )
              if( wantvr ) then
                 ! want left and right eigenvectors
                 ! copy schur vectors to vr
                 side = 'B'
                 call stdlib${ii}$_slacpy( 'F', n, n, vl, ldvl, vr, ldvr )
              end if
           else if( wantvr ) then
              ! want right eigenvectors
              ! copy householder vectors to vr
              side = 'R'
              call stdlib${ii}$_slacpy( 'L', n, n, a, lda, vr, ldvr )
              ! generate orthogonal matrix in vr
              ! (workspace: need 2*n-1, prefer n+(n-1)*nb)
              call stdlib${ii}$_sorghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vr
              ! (workspace: need 1, prefer hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_shseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           else
              ! compute eigenvalues only
              ! if condition numbers desired, compute schur form
              if( wntsnn ) then
                 job = 'E'
              else
                 job = 'S'
              end if
              ! (workspace: need 1, prefer hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_shseqr( job, 'N', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           end if
           ! if info /= 0 from stdlib${ii}$_shseqr, then quit
           if( info/=0 )go to 50
           if( wantvl .or. wantvr ) then
              ! compute left and/or right eigenvectors
              ! (workspace: need 3*n, prefer n + 2*n*nb)
              call stdlib${ii}$_strevc3( side, 'B', select, n, a, lda, vl, ldvl, vr, ldvr,n, nout, work(&
                         iwrk ), lwork-iwrk+1, ierr )
           end if
           ! compute condition numbers if desired
           ! (workspace: need n*n+6*n unless sense = 'e')
           if( .not.wntsnn ) then
              call stdlib${ii}$_strsna( sense, 'A', select, n, a, lda, vl, ldvl, vr, ldvr,rconde, &
                        rcondv, n, nout, work( iwrk ), n, iwork,icond )
           end if
           if( wantvl ) then
              ! undo balancing of left eigenvectors
              call stdlib${ii}$_sgebak( balanc, 'L', n, ilo, ihi, scale, n, vl, ldvl,ierr )
              ! normalize left eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_snrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_sscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_slapy2( stdlib${ii}$_snrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_snrm2( n, &
                              vl( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_sscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_sscal( n, scl, vl( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( k ) = vl( k, i )**2_${ik}$ + vl( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_isamax( n, work, 1_${ik}$ )
                    call stdlib${ii}$_slartg( vl( k, i ), vl( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_srot( n, vl( 1_${ik}$, i ), 1_${ik}$, vl( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vl( k, i+1 ) = zero
                 end if
              end do
           end if
           if( wantvr ) then
              ! undo balancing of right eigenvectors
              call stdlib${ii}$_sgebak( balanc, 'R', n, ilo, ihi, scale, n, vr, ldvr,ierr )
              ! normalize right eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_snrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_sscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_slapy2( stdlib${ii}$_snrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_snrm2( n, &
                              vr( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_sscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_sscal( n, scl, vr( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( k ) = vr( k, i )**2_${ik}$ + vr( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_isamax( n, work, 1_${ik}$ )
                    call stdlib${ii}$_slartg( vr( k, i ), vr( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_srot( n, vr( 1_${ik}$, i ), 1_${ik}$, vr( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vr( k, i+1 ) = zero
                 end if
              end do
           end if
           ! undo scaling if necessary
           50 continue
           if( scalea ) then
              call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wr( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wi( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              if( info==0_${ik}$ ) then
                 if( ( wntsnv .or. wntsnb ) .and. icond==0_${ik}$ )call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale,&
                            anrm, n, 1_${ik}$, rcondv, n,ierr )
              else
                 call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wr, n,ierr )
                 call stdlib${ii}$_slascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wi, n,ierr )
              end if
           end if
           work( 1_${ik}$ ) = maxwrk
           return
     end subroutine stdlib${ii}$_sgeevx

     module subroutine stdlib${ii}$_dgeevx( balanc, jobvl, jobvr, sense, n, a, lda, wr, wi,vl, ldvl, vr, ldvr, &
     !! DGEEVX computes for an N-by-N real nonsymmetric matrix A, the
     !! eigenvalues and, optionally, the left and/or right eigenvectors.
     !! Optionally also, it computes a balancing transformation to improve
     !! the conditioning of the eigenvalues and eigenvectors (ILO, IHI,
     !! SCALE, and ABNRM), reciprocal condition numbers for the eigenvalues
     !! (RCONDE), and reciprocal condition numbers for the right
     !! eigenvectors (RCONDV).
     !! The right eigenvector v(j) of A satisfies
     !! A * v(j) = lambda(j) * v(j)
     !! where lambda(j) is its eigenvalue.
     !! The left eigenvector u(j) of A satisfies
     !! u(j)**H * A = lambda(j) * u(j)**H
     !! where u(j)**H denotes the conjugate-transpose of u(j).
     !! The computed eigenvectors are normalized to have Euclidean norm
     !! equal to 1 and largest component real.
     !! Balancing a matrix means permuting the rows and columns to make it
     !! more nearly upper triangular, and applying a diagonal similarity
     !! transformation D * A * D**(-1), where D is a diagonal matrix, to
     !! make its rows and columns closer in norm and the condition numbers
     !! of its eigenvalues and eigenvectors smaller.  The computed
     !! reciprocal condition numbers correspond to the balanced matrix.
     !! Permuting rows and columns will not change the condition numbers
     !! (in exact arithmetic) but diagonal scaling will.  For further
     !! explanation of balancing, see section 4.10.2_dp of the LAPACK
     !! Users' Guide.
               ilo, ihi, scale, abnrm,rconde, rcondv, work, lwork, iwork, info )
        ! -- lapack driver routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: balanc, jobvl, jobvr, sense
           integer(${ik}$), intent(out) :: ihi, ilo, info
           integer(${ik}$), intent(in) :: lda, ldvl, ldvr, lwork, n
           real(dp), intent(out) :: abnrm
           ! Array Arguments 
           integer(${ik}$), intent(out) :: iwork(*)
           real(dp), intent(inout) :: a(lda,*)
           real(dp), intent(out) :: rconde(*), rcondv(*), scale(*), vl(ldvl,*), vr(ldvr,*), wi(*),&
                      work(*), wr(*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lquery, scalea, wantvl, wantvr, wntsnb, wntsne, wntsnn, wntsnv
           character :: job, side
           integer(${ik}$) :: hswork, i, icond, ierr, itau, iwrk, k, lwork_trevc, maxwrk, minwrk, &
                     nout
           real(dp) :: anrm, bignum, cs, cscale, eps, r, scl, smlnum, sn
           ! Local Arrays 
           logical(lk) :: select(1_${ik}$)
           real(dp) :: dum(1_${ik}$)
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input arguments
           info = 0_${ik}$
           lquery = ( lwork==-1_${ik}$ )
           wantvl = stdlib_lsame( jobvl, 'V' )
           wantvr = stdlib_lsame( jobvr, 'V' )
           wntsnn = stdlib_lsame( sense, 'N' )
           wntsne = stdlib_lsame( sense, 'E' )
           wntsnv = stdlib_lsame( sense, 'V' )
           wntsnb = stdlib_lsame( sense, 'B' )
           if( .not.( stdlib_lsame( balanc, 'N' ) .or. stdlib_lsame( balanc, 'S' ).or. &
                     stdlib_lsame( balanc, 'P' ) .or. stdlib_lsame( balanc, 'B' ) ) )then
              info = -1_${ik}$
           else if( ( .not.wantvl ) .and. ( .not.stdlib_lsame( jobvl, 'N' ) ) ) then
              info = -2_${ik}$
           else if( ( .not.wantvr ) .and. ( .not.stdlib_lsame( jobvr, 'N' ) ) ) then
              info = -3_${ik}$
           else if( .not.( wntsnn .or. wntsne .or. wntsnb .or. wntsnv ) .or.( ( wntsne .or. &
                     wntsnb ) .and. .not.( wantvl .and.wantvr ) ) ) then
              info = -4_${ik}$
           else if( n<0_${ik}$ ) then
              info = -5_${ik}$
           else if( lda<max( 1_${ik}$, n ) ) then
              info = -7_${ik}$
           else if( ldvl<1_${ik}$ .or. ( wantvl .and. ldvl<n ) ) then
              info = -11_${ik}$
           else if( ldvr<1_${ik}$ .or. ( wantvr .and. ldvr<n ) ) then
              info = -13_${ik}$
           end if
           ! compute workspace
            ! (note: comments in the code beginning "workspace:" describe the
             ! minimal amount of workspace needed at that point in the code,
             ! as well as the preferred amount for good performance.
             ! nb refers to the optimal block size for the immediately
             ! following subroutine, as returned by stdlib${ii}$_ilaenv.
             ! hswork refers to the workspace preferred by stdlib${ii}$_dhseqr, as
             ! calculated below. hswork is computed assuming ilo=1 and ihi=n,
             ! the worst case.)
           if( info==0_${ik}$ ) then
              if( n==0_${ik}$ ) then
                 minwrk = 1_${ik}$
                 maxwrk = 1_${ik}$
              else
                 maxwrk = n + n*stdlib${ii}$_ilaenv( 1_${ik}$, 'DGEHRD', ' ', n, 1_${ik}$, n, 0_${ik}$ )
                 if( wantvl ) then
                    call stdlib${ii}$_dtrevc3( 'L', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_dhseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vl, ldvl,work, -1_${ik}$, &
                              info )
                 else if( wantvr ) then
                    call stdlib${ii}$_dtrevc3( 'R', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_dhseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vr, ldvr,work, -1_${ik}$, &
                              info )
                 else
                    if( wntsnn ) then
                       call stdlib${ii}$_dhseqr( 'E', 'N', n, 1_${ik}$, n, a, lda, wr, wi, vr,ldvr, work, -1_${ik}$, &
                                 info )
                    else
                       call stdlib${ii}$_dhseqr( 'S', 'N', n, 1_${ik}$, n, a, lda, wr, wi, vr,ldvr, work, -1_${ik}$, &
                                 info )
                    end if
                 end if
                 hswork = int( work(1_${ik}$),KIND=${ik}$)
                 if( ( .not.wantvl ) .and. ( .not.wantvr ) ) then
                    minwrk = 2_${ik}$*n
                    if( .not.wntsnn )minwrk = max( minwrk, n*n+6*n )
                    maxwrk = max( maxwrk, hswork )
                    if( .not.wntsnn )maxwrk = max( maxwrk, n*n + 6_${ik}$*n )
                 else
                    minwrk = 3_${ik}$*n
                    if( ( .not.wntsnn ) .and. ( .not.wntsne ) )minwrk = max( minwrk, n*n + 6_${ik}$*n )
                              
                    maxwrk = max( maxwrk, hswork )
                    maxwrk = max( maxwrk, n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$, 'DORGHR',' ', n, 1_${ik}$, n, -&
                              1_${ik}$ ) )
                    if( ( .not.wntsnn ) .and. ( .not.wntsne ) )maxwrk = max( maxwrk, n*n + 6_${ik}$*n )
                              
                    maxwrk = max( maxwrk, 3_${ik}$*n )
                 end if
                 maxwrk = max( maxwrk, minwrk )
              end if
              work( 1_${ik}$ ) = maxwrk
              if( lwork<minwrk .and. .not.lquery ) then
                 info = -21_${ik}$
              end if
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'DGEEVX', -info )
              return
           else if( lquery ) then
              return
           end if
           ! quick return if possible
           if( n==0 )return
           ! get machine constants
           eps = stdlib${ii}$_dlamch( 'P' )
           smlnum = stdlib${ii}$_dlamch( 'S' )
           bignum = one / smlnum
           call stdlib${ii}$_dlabad( smlnum, bignum )
           smlnum = sqrt( smlnum ) / eps
           bignum = one / smlnum
           ! scale a if max element outside range [smlnum,bignum]
           icond = 0_${ik}$
           anrm = stdlib${ii}$_dlange( 'M', n, n, a, lda, dum )
           scalea = .false.
           if( anrm>zero .and. anrm<smlnum ) then
              scalea = .true.
              cscale = smlnum
           else if( anrm>bignum ) then
              scalea = .true.
              cscale = bignum
           end if
           if( scalea )call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, anrm, cscale, n, n, a, lda, ierr )
           ! balance the matrix and compute abnrm
           call stdlib${ii}$_dgebal( balanc, n, a, lda, ilo, ihi, scale, ierr )
           abnrm = stdlib${ii}$_dlange( '1', n, n, a, lda, dum )
           if( scalea ) then
              dum( 1_${ik}$ ) = abnrm
              call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, 1_${ik}$, 1_${ik}$, dum, 1_${ik}$, ierr )
              abnrm = dum( 1_${ik}$ )
           end if
           ! reduce to upper hessenberg form
           ! (workspace: need 2*n, prefer n+n*nb)
           itau = 1_${ik}$
           iwrk = itau + n
           call stdlib${ii}$_dgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),lwork-iwrk+1, ierr &
                     )
           if( wantvl ) then
              ! want left eigenvectors
              ! copy householder vectors to vl
              side = 'L'
              call stdlib${ii}$_dlacpy( 'L', n, n, a, lda, vl, ldvl )
              ! generate orthogonal matrix in vl
              ! (workspace: need 2*n-1, prefer n+(n-1)*nb)
              call stdlib${ii}$_dorghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vl
              ! (workspace: need 1, prefer hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_dhseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vl, ldvl,work( iwrk ), &
                        lwork-iwrk+1, info )
              if( wantvr ) then
                 ! want left and right eigenvectors
                 ! copy schur vectors to vr
                 side = 'B'
                 call stdlib${ii}$_dlacpy( 'F', n, n, vl, ldvl, vr, ldvr )
              end if
           else if( wantvr ) then
              ! want right eigenvectors
              ! copy householder vectors to vr
              side = 'R'
              call stdlib${ii}$_dlacpy( 'L', n, n, a, lda, vr, ldvr )
              ! generate orthogonal matrix in vr
              ! (workspace: need 2*n-1, prefer n+(n-1)*nb)
              call stdlib${ii}$_dorghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vr
              ! (workspace: need 1, prefer hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_dhseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           else
              ! compute eigenvalues only
              ! if condition numbers desired, compute schur form
              if( wntsnn ) then
                 job = 'E'
              else
                 job = 'S'
              end if
              ! (workspace: need 1, prefer hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_dhseqr( job, 'N', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           end if
           ! if info /= 0 from stdlib${ii}$_dhseqr, then quit
           if( info/=0 )go to 50
           if( wantvl .or. wantvr ) then
              ! compute left and/or right eigenvectors
              ! (workspace: need 3*n, prefer n + 2*n*nb)
              call stdlib${ii}$_dtrevc3( side, 'B', select, n, a, lda, vl, ldvl, vr, ldvr,n, nout, work(&
                         iwrk ), lwork-iwrk+1, ierr )
           end if
           ! compute condition numbers if desired
           ! (workspace: need n*n+6*n unless sense = 'e')
           if( .not.wntsnn ) then
              call stdlib${ii}$_dtrsna( sense, 'A', select, n, a, lda, vl, ldvl, vr, ldvr,rconde, &
                        rcondv, n, nout, work( iwrk ), n, iwork,icond )
           end if
           if( wantvl ) then
              ! undo balancing of left eigenvectors
              call stdlib${ii}$_dgebak( balanc, 'L', n, ilo, ihi, scale, n, vl, ldvl,ierr )
              ! normalize left eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_dnrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_dscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_dlapy2( stdlib${ii}$_dnrm2( n, vl( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_dnrm2( n, &
                              vl( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_dscal( n, scl, vl( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_dscal( n, scl, vl( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( k ) = vl( k, i )**2_${ik}$ + vl( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_idamax( n, work, 1_${ik}$ )
                    call stdlib${ii}$_dlartg( vl( k, i ), vl( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_drot( n, vl( 1_${ik}$, i ), 1_${ik}$, vl( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vl( k, i+1 ) = zero
                 end if
              end do
           end if
           if( wantvr ) then
              ! undo balancing of right eigenvectors
              call stdlib${ii}$_dgebak( balanc, 'R', n, ilo, ihi, scale, n, vr, ldvr,ierr )
              ! normalize right eigenvectors and make largest component real
              do i = 1, n
                 if( wi( i )==zero ) then
                    scl = one / stdlib${ii}$_dnrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_dscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                 else if( wi( i )>zero ) then
                    scl = one / stdlib${ii}$_dlapy2( stdlib${ii}$_dnrm2( n, vr( 1_${ik}$, i ), 1_${ik}$ ),stdlib${ii}$_dnrm2( n, &
                              vr( 1_${ik}$, i+1 ), 1_${ik}$ ) )
                    call stdlib${ii}$_dscal( n, scl, vr( 1_${ik}$, i ), 1_${ik}$ )
                    call stdlib${ii}$_dscal( n, scl, vr( 1_${ik}$, i+1 ), 1_${ik}$ )
                    do k = 1, n
                       work( k ) = vr( k, i )**2_${ik}$ + vr( k, i+1 )**2_${ik}$
                    end do
                    k = stdlib${ii}$_idamax( n, work, 1_${ik}$ )
                    call stdlib${ii}$_dlartg( vr( k, i ), vr( k, i+1 ), cs, sn, r )
                    call stdlib${ii}$_drot( n, vr( 1_${ik}$, i ), 1_${ik}$, vr( 1_${ik}$, i+1 ), 1_${ik}$, cs, sn )
                    vr( k, i+1 ) = zero
                 end if
              end do
           end if
           ! undo scaling if necessary
           50 continue
           if( scalea ) then
              call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wr( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, n-info, 1_${ik}$, wi( info+1 ),max( n-info, 1_${ik}$ &
                        ), ierr )
              if( info==0_${ik}$ ) then
                 if( ( wntsnv .or. wntsnb ) .and. icond==0_${ik}$ )call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale,&
                            anrm, n, 1_${ik}$, rcondv, n,ierr )
              else
                 call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wr, n,ierr )
                 call stdlib${ii}$_dlascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, ilo-1, 1_${ik}$, wi, n,ierr )
              end if
           end if
           work( 1_${ik}$ ) = maxwrk
           return
     end subroutine stdlib${ii}$_dgeevx

#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
     module subroutine stdlib${ii}$_${ri}$geevx( balanc, jobvl, jobvr, sense, n, a, lda, wr, wi,vl, ldvl, vr, ldvr, &
     !! DGEEVX: computes for an N-by-N real nonsymmetric matrix A, the
     !! eigenvalues and, optionally, the left and/or right eigenvectors.
     !! Optionally also, it computes a balancing transformation to improve
     !! the conditioning of the eigenvalues and eigenvectors (ILO, IHI,
     !! SCALE, and ABNRM), reciprocal condition numbers for the eigenvalues
     !! (RCONDE), and reciprocal condition numbers for the right
     !! eigenvectors (RCONDV).
     !! The right eigenvector v(j) of A satisfies
     !! A * v(j) = lambda(j) * v(j)
     !! where lambda(j) is its eigenvalue.
     !! The left eigenvector u(j) of A satisfies
     !! u(j)**H * A = lambda(j) * u(j)**H
     !! where u(j)**H denotes the conjugate-transpose of u(j).
     !! The computed eigenvectors are normalized to have Euclidean norm
     !! equal to 1 and largest component real.
     !! Balancing a matrix means permuting the rows and columns to make it
     !! more nearly upper triangular, and applying a diagonal similarity
     !! transformation D * A * D**(-1), where D is a diagonal matrix, to
     !! make its rows and columns closer in norm and the condition numbers
     !! of its eigenvalues and eigenvectors smaller.  The computed
     !! reciprocal condition numbers correspond to the balanced matrix.
     !! Permuting rows and columns will not change the condition numbers
     !! (in exact arithmetic) but diagonal scaling will.  For further
     !! explanation of balancing, see section 4.10.2_${rk}$ of the LAPACK
     !! Users' Guide.
               ilo, ihi, scale, abnrm,rconde, rcondv, work, lwork, iwork, info )
        ! -- lapack driver routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_${rk}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: balanc, jobvl, jobvr, sense
           integer(${ik}$), intent(out) :: ihi, ilo, info
           integer(${ik}$), intent(in) :: lda, ldvl, ldvr, lwork, n
           real(${rk}$), intent(out) :: abnrm
           ! Array Arguments 
           integer(${ik}$), intent(out) :: iwork(*)
           real(${rk}$), intent(inout) :: a(lda,*)
           real(${rk}$), intent(out) :: rconde(*), rcondv(*), scale(*), vl(ldvl,*), vr(ldvr,*), wi(*),&
                      work(*), wr(*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lquery, scalea, wantvl, wantvr, wntsnb, wntsne, wntsnn, wntsnv
           character :: job, side
           integer(${ik}$) :: hswork, i, icond, ierr, itau, iwrk, k, lwork_trevc, maxwrk, minwrk, &
                     nout
           real(${rk}$) :: anrm, bignum, cs, cscale, eps, r, scl, smlnum, sn
           ! Local Arrays 
           logical(lk) :: select(1_${ik}$)
           real(${rk}$) :: dum(1_${ik}$)
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input arguments
           info = 0_${ik}$
           lquery = ( lwork==-1_${ik}$ )
           wantvl = stdlib_lsame( jobvl, 'V' )
           wantvr = stdlib_lsame( jobvr, 'V' )
           wntsnn = stdlib_lsame( sense, 'N' )
           wntsne = stdlib_lsame( sense, 'E' )
           wntsnv = stdlib_lsame( sense, 'V' )
           wntsnb = stdlib_lsame( sense, 'B' )
           if( .not.( stdlib_lsame( balanc, 'N' ) .or. stdlib_lsame( balanc, 'S' ).or. &
                     stdlib_lsame( balanc, 'P' ) .or. stdlib_lsame( balanc, 'B' ) ) )then
              info = -1_${ik}$
           else if( ( .not.wantvl ) .and. ( .not.stdlib_lsame( jobvl, 'N' ) ) ) then
              info = -2_${ik}$
           else if( ( .not.wantvr ) .and. ( .not.stdlib_lsame( jobvr, 'N' ) ) ) then
              info = -3_${ik}$
           else if( .not.( wntsnn .or. wntsne .or. wntsnb .or. wntsnv ) .or.( ( wntsne .or. &
                     wntsnb ) .and. .not.( wantvl .and.wantvr ) ) ) then
              info = -4_${ik}$
           else if( n<0_${ik}$ ) then
              info = -5_${ik}$
           else if( lda<max( 1_${ik}$, n ) ) then
              info = -7_${ik}$
           else if( ldvl<1_${ik}$ .or. ( wantvl .and. ldvl<n ) ) then
              info = -11_${ik}$
           else if( ldvr<1_${ik}$ .or. ( wantvr .and. ldvr<n ) ) then
              info = -13_${ik}$
           end if
           ! compute workspace
            ! (note: comments in the code beginning "workspace:" describe the
             ! minimal amount of workspace needed at that point in the code,
             ! as well as the preferred amount for good performance.
             ! nb refers to the optimal block size for the immediately
             ! following subroutine, as returned by stdlib${ii}$_ilaenv.
             ! hswork refers to the workspace preferred by stdlib${ii}$_${ri}$hseqr, as
             ! calculated below. hswork is computed assuming ilo=1 and ihi=n,
             ! the worst case.)
           if( info==0_${ik}$ ) then
              if( n==0_${ik}$ ) then
                 minwrk = 1_${ik}$
                 maxwrk = 1_${ik}$
              else
                 maxwrk = n + n*stdlib${ii}$_ilaenv( 1_${ik}$, 'DGEHRD', ' ', n, 1_${ik}$, n, 0_${ik}$ )
                 if( wantvl ) then
                    call stdlib${ii}$_${ri}$trevc3( 'L', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_${ri}$hseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vl, ldvl,work, -1_${ik}$, &
                              info )
                 else if( wantvr ) then
                    call stdlib${ii}$_${ri}$trevc3( 'R', 'B', select, n, a, lda,vl, ldvl, vr, ldvr,n, nout, &
                              work, -1_${ik}$, ierr )
                    lwork_trevc = int( work(1_${ik}$),KIND=${ik}$)
                    maxwrk = max( maxwrk, n + lwork_trevc )
                    call stdlib${ii}$_${ri}$hseqr( 'S', 'V', n, 1_${ik}$, n, a, lda, wr, wi, vr, ldvr,work, -1_${ik}$, &
                              info )
                 else
                    if( wntsnn ) then
                       call stdlib${ii}$_${ri}$hseqr( 'E', 'N', n, 1_${ik}$, n, a, lda, wr, wi, vr,ldvr, work, -1_${ik}$, &
                                 info )
                    else
                       call stdlib${ii}$_${ri}$hseqr( 'S', 'N', n, 1_${ik}$, n, a, lda, wr, wi, vr,ldvr, work, -1_${ik}$, &
                                 info )
                    end if
                 end if
                 hswork = int( work(1_${ik}$),KIND=${ik}$)
                 if( ( .not.wantvl ) .and. ( .not.wantvr ) ) then
                    minwrk = 2_${ik}$*n
                    if( .not.wntsnn )minwrk = max( minwrk, n*n+6*n )
                    maxwrk = max( maxwrk, hswork )
                    if( .not.wntsnn )maxwrk = max( maxwrk, n*n + 6_${ik}$*n )
                 else
                    minwrk = 3_${ik}$*n
                    if( ( .not.wntsnn ) .and. ( .not.wntsne ) )minwrk = max( minwrk, n*n + 6_${ik}$*n )
                              
                    maxwrk = max( maxwrk, hswork )
                    maxwrk = max( maxwrk, n + ( n - 1_${ik}$ )*stdlib${ii}$_ilaenv( 1_${ik}$, 'DORGHR',' ', n, 1_${ik}$, n, -&
                              1_${ik}$ ) )
                    if( ( .not.wntsnn ) .and. ( .not.wntsne ) )maxwrk = max( maxwrk, n*n + 6_${ik}$*n )
                              
                    maxwrk = max( maxwrk, 3_${ik}$*n )
                 end if
                 maxwrk = max( maxwrk, minwrk )
              end if
              work( 1_${ik}$ ) = maxwrk
              if( lwork<minwrk .and. .not.lquery ) then
                 info = -21_${ik}$
              end if
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'DGEEVX', -info )
              return
           else if( lquery ) then
              return
           end if
           ! quick return if possible
           if( n==0 )return
           ! get machine constants
           eps = stdlib${ii}$_${ri}$lamch( 'P' )
           smlnum = stdlib${ii}$_${ri}$lamch( 'S' )
           bignum = one / smlnum
           call stdlib${ii}$_${ri}$labad( smlnum, bignum )
           smlnum = sqrt( smlnum ) / eps
           bignum = one / smlnum
           ! scale a if max element outside range [smlnum,bignum]
           icond = 0_${ik}$
           anrm = stdlib${ii}$_${ri}$lange( 'M', n, n, a, lda, dum )
           scalea = .false.
           if( anrm>zero .and. anrm<smlnum ) then
              scalea = .true.
              cscale = smlnum
           else if( anrm>bignum ) then
              scalea = .true.
              cscale = bignum
           end if
           if( scalea )call stdlib${ii}$_${ri}$lascl( 'G', 0_${ik}$, 0_${ik}$, anrm, cscale, n, n, a, lda, ierr )
           ! balance the matrix and compute abnrm
           call stdlib${ii}$_${ri}$gebal( balanc, n, a, lda, ilo, ihi, scale, ierr )
           abnrm = stdlib${ii}$_${ri}$lange( '1', n, n, a, lda, dum )
           if( scalea ) then
              dum( 1_${ik}$ ) = abnrm
              call stdlib${ii}$_${ri}$lascl( 'G', 0_${ik}$, 0_${ik}$, cscale, anrm, 1_${ik}$, 1_${ik}$, dum, 1_${ik}$, ierr )
              abnrm = dum( 1_${ik}$ )
           end if
           ! reduce to upper hessenberg form
           ! (workspace: need 2*n, prefer n+n*nb)
           itau = 1_${ik}$
           iwrk = itau + n
           call stdlib${ii}$_${ri}$gehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),lwork-iwrk+1, ierr &
                     )
           if( wantvl ) then
              ! want left eigenvectors
              ! copy householder vectors to vl
              side = 'L'
              call stdlib${ii}$_${ri}$lacpy( 'L', n, n, a, lda, vl, ldvl )
              ! generate orthogonal matrix in vl
              ! (workspace: need 2*n-1, prefer n+(n-1)*nb)
              call stdlib${ii}$_${ri}$orghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vl
              ! (workspace: need 1, prefer hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_${ri}$hseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vl, ldvl,work( iwrk ), &
                        lwork-iwrk+1, info )
              if( wantvr ) then
                 ! want left and right eigenvectors
                 ! copy schur vectors to vr
                 side = 'B'
                 call stdlib${ii}$_${ri}$lacpy( 'F', n, n, vl, ldvl, vr, ldvr )
              end if
           else if( wantvr ) then
              ! want right eigenvectors
              ! copy householder vectors to vr
              side = 'R'
              call stdlib${ii}$_${ri}$lacpy( 'L', n, n, a, lda, vr, ldvr )
              ! generate orthogonal matrix in vr
              ! (workspace: need 2*n-1, prefer n+(n-1)*nb)
              call stdlib${ii}$_${ri}$orghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),lwork-iwrk+1, &
                        ierr )
              ! perform qr iteration, accumulating schur vectors in vr
              ! (workspace: need 1, prefer hswork (see comments) )
              iwrk = itau
              call stdlib${ii}$_${ri}$hseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,work( iwrk ), &
                        lwork-iwrk+1, info )
           else
              ! compute eigenvalues only
              ! if condition numbers desired, compute schur form
              if(