HDIAG_neci Subroutine

public subroutine HDIAG_neci(NDET, HAMIL, LAB, NROW, CK, W, WORK2, WORK, NBLOCKSTARTS, NBLOCKS)


.. Diagonalize

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: NDET
real(kind=dp), intent(in) :: HAMIL(*)
integer, intent(in) :: LAB(*)
integer, intent(in) :: NROW(NDET)
real(kind=dp), intent(out) :: CK(NDET,NDET)
real(kind=dp), intent(inout) :: W(NDET)
real(kind=dp), intent(out) :: WORK2(3*NDET)
real(kind=dp), intent(out) :: WORK(*)
integer, intent(in) :: NBLOCKSTARTS(NBLOCKS+1)
integer, intent(in) :: NBLOCKS

Contents

Source Code


Source Code

    SUBROUTINE HDIAG_neci(NDET, HAMIL, LAB, NROW, CK, W, WORK2, WORK, NBLOCKSTARTS, NBLOCKS)
        INTEGER, intent(in) :: NDET, NROW(NDET), LAB(*)
        HElement_t(dp), intent(in) :: HAMIL(*)
        HElement_t(dp), intent(out) :: CK(NDET, NDET), WORK(*), WORK2(3 * NDET)
        real(dp), intent(inout) :: W(NDET)
        INTEGER, intent(in) :: NBLOCKS, NBLOCKSTARTS(NBLOCKS + 1)
        integer :: i, j, ind, indz, nbs
        INTEGER(int32) :: INFO
        type(timer), save :: proc_timer
        real(dp) :: GSEN
        character(len=*), parameter :: t_r = "HDIAG_neci"
#ifndef CMPLX_
        associate(tmp => WORK(1:1)); end associate
#endif

        proc_timer%timer_name = 'HDIAG     '
        call set_timer(proc_timer)
        GSEN = 1.D100
        CK(:, :) = 0._dp
!.. Now we fill the RIJ array
        IND = 1
        INDZ = 1
        DO I = 1, NDET
            INDZ = INDZ + NROW(I)
            DO WHILE (IND < INDZ)
                CK(I, LAB(IND)) = HAMIL(IND)
                IND = IND + 1
            END DO
        END DO
        IF (HElement_t_size /= 1) THEN
!ensure hermitian
            do i = 1, ndet
                do j = 1, i - 1
#ifdef CMPLX_
                    !To avoid compiler warnings
                    CK(i, j) = CONJG(CMPLX(CK(j, i), kind=dp))
#endif
                end do
            end do
        END IF

!****************************
!.. Diagonalize
        DO I = 1, NBLOCKS
            NBS = NBLOCKSTARTS(I)
            IF (HElement_t_size == 1) THEN
#ifndef CMPLX_
                CALL DSYEV('V', 'U', NBLOCKSTARTS(I + 1) - NBS, CK(NBS, NBS), NDET, W(NBS), WORK2, 3 * NDET, INFO)
#else
                call stop_all(t_r, "This block for real calcs only")
#endif
            ELSE
#ifdef CMPLX_
                CALL ZHEEV('V', 'U', NBLOCKSTARTS(I + 1) - NBS, CK(NBS, NBS), NDET, W(NBS), WORK, 4 * NDET, WORK2, INFO)
#else
                call stop_all(t_r, "This block for complex calcs only")
#endif
            END IF
            IF (INFO /= 0) THEN
                WRITE(stdout, *) 'DYSEV error: ', INFO
                call stop_all(t_r, "DSYEV error")
            END IF
            IF (W(NBS) < GSEN) GSEN = W(NBS)
        END DO

!.. CK now contains the eigenvectors, and W the eigenvalues
        WRITE(stdout, "(A,F19.11,I4)") "GROUND E=", GSEN
        call halt_timer(proc_timer)
    end subroutine hdiag_neci