BinSearchListHPHF Subroutine

public subroutine BinSearchListHPHF(iLut, List, Length, MinInd, MaxInd, PartInd, tSuccess)

Arguments

Type IntentOptional Attributes Name
integer(kind=n_int) :: iLut(0:NIfTot)
integer(kind=n_int) :: List(0:NIfTot,Length)
integer :: Length
integer :: MinInd
integer :: MaxInd
integer :: PartInd
logical :: tSuccess

Contents

Source Code


Source Code

    SUBROUTINE BinSearchListHPHF(iLut, List, Length, MinInd, MaxInd, PartInd, tSuccess)
        INTEGER :: Length, MinInd, MaxInd, PartInd
        INTEGER(KIND=n_int) :: iLut(0:NIfTot), List(0:NIfTot, Length)
        INTEGER :: i, j, N, Comp
        LOGICAL :: tSuccess

        i = MinInd
        j = MaxInd
        IF(i - j == 0) THEN
            Comp = DetBitLT(List(:, MaxInd), iLut(:), nifd)
            IF(Comp == 0) THEN
                tSuccess = .true.
                PartInd = MaxInd
                RETURN
            ELSE
                tSuccess = .false.
                PartInd = MinInd
            end if
        end if
        do while(j - i > 0)  !End when the upper and lower bound are the same.
            N = (i + j) / 2       !Find the midpoint of the two indices

!Comp is 1 if CyrrebtDets(N) is "less" than iLut, and -1 if it is more or 0 if they are the same
            Comp = DetBitLT(List(:, N), iLut(:), nifd)

            IF(Comp == 0) THEN
!Praise the lord, we've found it!
                tSuccess = .true.
                PartInd = N
                RETURN
            else if((Comp == 1) .and. (i /= N)) THEN
!The value of the determinant at N is LESS than the determinant we're looking for.
!Therefore, move the lower bound of the search up to N.
!However, if the lower bound is already equal to N then the two bounds are consecutive and we have failed...
                i = N
            else if(i == N) THEN

                IF(i == MaxInd - 1) THEN
!This deals with the case where we are interested in the final/first entry in the list. Check the final entry of the list and leave
!We need to check the last index.
                    Comp = DetBitLT(List(:, i + 1), iLut(:), nifd)
                    IF(Comp == 0) THEN
                        tSuccess = .true.
                        PartInd = i + 1
                        RETURN
                    else if(Comp == 1) THEN
!final entry is less than the one we want.
                        tSuccess = .false.
                        PartInd = i + 1
                        RETURN
                    ELSE
                        tSuccess = .false.
                        PartInd = i
                        RETURN
                    end if

                else if(i == MinInd) THEN
                    tSuccess = .false.
                    PartInd = i
                    RETURN
                ELSE
                    i = j
                end if

            else if(Comp == -1) THEN
!The value of the determinant at N is MORE than the determinant we're looking for. Move the upper bound of the search down to N.
                j = N
            ELSE
!We have failed - exit loop
                i = j
            end if

        end do

!If we have failed, then we want to find the index that is one less than where the particle would have been.
        tSuccess = .false.
        PartInd = MAX(MinInd, i - 1)

    END SUBROUTINE BinSearchListHPHF