fill_in_hash_table Subroutine

public subroutine fill_in_hash_table(hash_table, table_length, walker_list, list_length, ignore_unocc)

Uses

Arguments

Type IntentOptional Attributes Name
type(ll_node), intent(inout), pointer :: hash_table(:)
integer, intent(in) :: table_length
integer(kind=n_int), intent(in) :: walker_list(0:,:)
integer, intent(in) :: list_length
logical, intent(in) :: ignore_unocc

Contents

Source Code


Source Code

    subroutine fill_in_hash_table(hash_table, table_length, walker_list, list_length, ignore_unocc)

        ! This assumes that the input hash table is clear (use clear_hash_table)
        ! and that there are no repeats in walker_list.

        ! If ignore_unocc dets is input as .true. then unoccupied determinants
        ! will not be included in the hash table, unless they are core
        ! determinants.

        use DetBitOps, only: tAccumEmptyDet
        integer, intent(in) :: table_length, list_length
        type(ll_node), pointer, intent(inout) :: hash_table(:)
        integer(n_int), intent(in) :: walker_list(0:, :)
        logical, intent(in) :: ignore_unocc

        type(ll_node), pointer :: temp_node
        integer :: i, hash_val, nI(nel), run
        real(dp) :: real_sign(lenof_sign)
        logical :: tCoreDet
#ifdef DEBUG_
        character(*), parameter :: this_routine = "fill_in_hash_table"
#endif

        tCoreDet = .false.

        do i = 1, list_length
            ! If the ignore_unocc option is true then we don't want to include
            ! unoccupied determinants in the hash table, unless they're
            ! deterministic states.
            if (ignore_unocc) then
                if (tSemiStochastic) then
                    tCoreDet = test_flag_multi(walker_list(:, i), flag_deterministic)
                end if
                call extract_sign(walker_list(:, i), real_sign)
                if (IsUnoccDet(real_sign) .and. (.not. tCoreDet) .and. (.not. tAccumEmptyDet(walker_list(:, i)))) cycle
            end if

            call decode_bit_det(nI, walker_list(:, i))
            ! Find the hash value corresponding to this determinant.
            ASSERT(all(nI <= nBasis))
            ASSERT(all(nI > 0))

            hash_val = FindWalkerHash(nI, table_length)

            call add_hash_table_entry(hash_table, i, hash_val)
        end do

    end subroutine fill_in_hash_table