function deigvals(A) result(lam)
real(8), intent(in) :: A(:, :) ! matrix for eigenvalue compuation
complex(8) :: lam(size(A,1)) ! eigenvalues: A c = lam c
real(8), allocatable :: At(:,:),vl(:,: ),vr(:,:),wi(:),work(:),wr(:)
real(8) :: lwork_guess(1)
integer :: info, lda, ldvl, ldvr, lwork, n
lda = size(A(:,1))
n = size(A(1,:))
call assert_shape(A,[n,n],"solve","A")
ldvl = n
ldvr = n
allocate(At(n,n),wr(n),wi(n),vl(ldvl,n),vr(ldvr,n))
!Copy the input Matrix
At = A
!1st Call: Query the right size for the working array.
call dgeev('N', 'N', n, At, lda, wr, wi, vl, ldvl, vr, ldvr, lwork_guess, -1, info)
if(info /= 0) then
print*, "dgeev returned info = ",info
if (info < 0) then
print*, "the",-info,"-th argument had an illegal value"
else
print*, "the QR algorithm failed to compute all the"
print*, "eigenvalues, and no eigenvectors have been computed;"
print*, "elements ", info+1, ":", n, "of WR and WI contain eigenvalues which"
print*, "have converged."
end if
stop 'deig error: 1st call dgeev'
endif
lwork = int(lwork_guess(1)) ;allocate(work(lwork))
!2nd Call: Actual solution of the eigenproblem
call dgeev('N', 'N', n, At, lda, wr, wi, vl, ldvl, vr, ldvr, work, lwork, info)
if(info /= 0) then
print *, "dgeev returned info = ",info
if (info < 0) then
print*, "the",-info,"-th argument had an illegal value"
else
print*, "the QR algorithm failed to compute all the"
print*, "eigenvalues, and no eigenvectors have been computed;"
print*, "elements ", info+1, ":", n, "of WR and WI contain eigenvalues which"
print*, "have converged."
end if
stop 'deig error: 2nd call dgeev'
end if
lam = wr + xi*wi
deallocate(At,wr,wi,vl,vr)
end function deigvals