Add_ExplicitRDM_Contrib Subroutine

    subroutine Add_ExplicitRDM_Contrib(iLutnI, blank_det)

        ! This is the general routine for taking a particular determinant in
        ! the spawned list, D_i and adding it's contribution to the reduced
        ! density matrix.

        use LoggingData, only: RDMExcitLevel
        use Parallel_neci, only: nProcessors
        use rdm_data, only: Sing_ExcList, Doub_ExcList, Sing_InitExcSlots, Doub_InitExcSlots
        use rdm_data, only: Sing_ExcDjs, Doub_ExcDjs
        use bit_rep_data, only: nifguga

        integer(n_int), intent(in) :: iLutnI(0:NIfTot)
        logical, intent(in) :: blank_det
        integer :: i, ni(nel), FlagsDi
        integer(n_int) :: ilutG(0:nifguga)
        real(dp) :: SignDi(lenof_sign)

        ! Set up excitation arrays.
        ! These are blocked according to the processor the excitation would be
        ! on if occupied. In each block, the first entry is the sign of
        ! determinant D_i and the second the bit string of the determinant
        ! (these need to be sent along with the excitations). Each processor
        ! will have a different Di.

        Sing_ExcDjs(:, :) = 0
        Sing_ExcList(:) = 0
        Sing_ExcList(:) = Sing_InitExcSlots(:)

        if (tGUGA) then
            call convert_ilut_toGUGA(ilutNi, ilutG)

            call extract_bit_rep(iLutnI, nI, SignDi, FlagsDi)
            call encode_matrix_element(ilutG, SignDi(1), 2)

            do i = 0, nProcessors - 1
                Sing_ExcDjs(:, Sing_ExcList(i)) = ilutG
                Sing_ExcList(i) = Sing_ExcList(i) + 1
            end do

        else

            do i = 0, nProcessors - 1
                Sing_ExcDjs(:, Sing_ExcList(i)) = iLutnI(:)
                Sing_ExcList(i) = Sing_ExcList(i) + 1
            end do
        end if

        if (RDMExcitLevel /= 1) then
            Doub_ExcDjs(:, :) = 0
            Doub_ExcList(:) = 0
            Doub_ExcList(:) = Doub_InitExcSlots(:)

            if (tGUGA) then
                do i = 0, nProcessors - 1
                    Doub_ExcDjs(:, Doub_ExcList(i)) = ilutG
                    Doub_ExcList(i) = Doub_ExcList(i) + 1
                end do
            else
                do i = 0, nProcessors - 1
                    Doub_ExcDjs(:, Doub_ExcList(i)) = iLutnI(:)
                    Doub_ExcList(i) = Doub_ExcList(i) + 1
                end do
            end if
        end if

        ! Out of here we will get a filled ExcDjs array with all the single or
        ! double excitations from Dj, this will be done for each proc.
        if (.not. blank_det) then
            if (tGUGA) then
                call gen_exc_djs_guga(ilutnI, CSF_Info_t(ilutnI))
            else
                call GenExcDjs(iLutnI)
            end if
        end if

        ! We then need to send the excitations to the relevant processors.
        ! This routine then calls SearchOccDets which takes each excitation
        ! and and binary searches the occupied determinants for this. If found,
        ! we re-find the orbitals and parity involved in the excitation, and
        ! add the c_i*c_j contributions to the corresponding matrix element.
        if (tGUGA) then
            call send_proc_ex_djs()
        else
            call SendProcExcDjs()
        end if

    end subroutine Add_ExplicitRDM_Contrib