assign_excits_to_proc_guga Subroutine

private subroutine assign_excits_to_proc_guga(n_tot, excits, excit_lvl)

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: n_tot
integer(kind=n_int), intent(in), allocatable :: excits(:,:)
integer, intent(in) :: excit_lvl

Contents


Source Code

    subroutine assign_excits_to_proc_guga(n_tot, excits, excit_lvl)
        integer, intent(in) :: n_tot, excit_lvl
        integer(n_int), intent(in), allocatable :: excits(:, :)
        character(*), parameter :: this_routine = "assign_excits_to_proc_guga"

        integer :: i, proc, nJ(nel)
        integer(n_int) :: ilut(0:niftot)

        ASSERT(excit_lvl == 1 .or. excit_lvl == 2)

        if (excit_lvl == 1) then
            do i = 1, n_tot

                call convert_ilut_toNECI(excits(:, i), ilut)

                call decode_bit_det(nJ, ilut)

                proc = DetermineDetNode(nel, nJ, 0)

                Sing_ExcDjs(:, Sing_ExcList(proc)) = excits(:, i)
                Sing_ExcList(proc) = Sing_ExcList(proc) + 1

#ifdef DEBUG_
                if (Sing_ExcList(Proc) > nint(OneEl_Gap * (Proc + 1))) then
                    write(stdout, *) 'Proc', Proc
                    write(stdout, *) 'Sing_ExcList', Sing_ExcList
                    write(stdout, *) 'No. spaces for each proc', nint(OneEl_Gap)
                    call Stop_All('GenExcDjs', 'Too many excitations for space available.')
                end if
#endif
            end do

        else if (excit_lvl == 2) then
            do i = 1, n_tot

                call convert_ilut_toNECI(excits(:, i), ilut)

                call decode_bit_det(nJ, ilut)
                proc = DetermineDetNode(nel, nJ, 0)

                Doub_ExcDjs(:, Doub_ExcList(proc)) = excits(:, i)
                Doub_ExcList(proc) = Doub_ExcList(proc) + 1

#ifdef DEBUG_
                if (Doub_ExcList(Proc) > nint(TwoEl_Gap * (Proc + 1))) then
                    write(stdout, *) 'Proc', Proc
                    write(stdout, *) 'Doub_ExcList', Doub_ExcList
                    write(stdout, *) 'No. spaces for each proc', nint(TwoEl_Gap)
                    call Stop_All('GenExcDjs', 'Too many excitations for space available.')
                end if
#endif
            end do
        end if

    end subroutine assign_excits_to_proc_guga