#: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(