diagonalise_tridiagonal_non_hermitian Subroutine

public subroutine diagonalise_tridiagonal_non_hermitian(this, N, t_calc_eigenvectors)

Arguments

Type IntentOptional Attributes Name
type(LanczosCalcType), intent(inout) :: this
integer, intent(in) :: N
logical, intent(in) :: t_calc_eigenvectors

Contents


Source Code

    subroutine diagonalise_tridiagonal_non_hermitian(this, N, t_calc_eigenvectors)
        type(LanczosCalcType), intent(inout) :: this
        integer, intent(in) :: N
        logical, intent(in) :: t_calc_eigenvectors
        integer :: lwork, info
        real(dp) :: rwork(4 * N)
        character :: jobz
        real(dp) :: dummy_ev(N), dummy_vec(1, N), e_vectors(N, N)
        logical :: t_sort = .true.
        integer :: sort_ind(N)

        lwork = 4 * N
        if (t_calc_eigenvectors) then
            jobz = 'V'
        else
            jobz = 'N'
        end if

        if (.not. t_calc_eigenvectors) then
            this%ritz_values_old(1:this%n_states) = this%ritz_values(1:this%n_states)
            this%super%projected_hamil_work(1:N, 1:N) = this%super%projected_hamil(1:N, 1:N)
        end if

        call dgeev( &
            'N', &
            jobz, &
            N, &
            this%super%projected_hamil_work(1:N, 1:N), &
            N, &
            this%ritz_values(1:N), &
            dummy_ev, &
            dummy_vec, &
            1, &
            e_vectors, &
            N, &
            rwork, &
            lwork, &
            info)

        if (t_calc_eigenvectors) then
            ! move the eigenvectors out of the working array
#ifdef CMPLX_
            this%T_eigenvectors(1:N, 1:this%n_states) = &
                cmplx(e_vectors(1:N, 1:this%n_states), kind=dp)
#else
            this%T_eigenvectors(1:N, 1:this%n_states) = &
                e_vectors(1:N, 1:this%n_states)
#endif
        end if

    end subroutine diagonalise_tridiagonal_non_hermitian