generate_core_connections Subroutine

    subroutine generate_core_connections(rep)

        use DetBitOps, only: FindBitExcitLevel, GetBitExcitation
        use Parallel_neci, only: MPIAllGatherV
        type(core_space_t), intent(inout) :: rep
        integer :: i, j, ic, counter, ierr
        integer :: Ex(2, nel)
        logical :: tSign
        integer(n_int), allocatable, dimension(:, :) :: temp_store
        integer(TagIntType) :: TempStoreTag
        character(len=*), parameter :: t_r = "calculate_det_hamiltonian_sparse"

        allocate(rep%core_connections(rep%determ_sizes(iProcIndex)))

        allocate(temp_store(0:NIfTot, rep%determ_space_size), stat=ierr)
        call LogMemAlloc('temp_store', maxval(rep%determ_sizes) * (NIfTot + 1), 8, t_r, &
                         TempStoreTag, ierr)

        ! Stick together the deterministic states from all processors, on all processors.
        call MPIAllGatherV(SpawnedParts(0:NIfTot, 1:rep%determ_sizes(iProcIndex)), temp_store, &
                           rep%determ_sizes, rep%determ_displs)

        ! Over all core states on this processor.
        do i = 1, rep%determ_sizes(iProcIndex)

            ! The number of non-zero elements in this array will be almost the same as in
            ! the core Hamiltonian array, except the diagonal element is not considered,
            ! so there will actually be one less.
            allocate(rep%core_connections(i)%elements(rep%sparse_core_ham(i)%num_elements - 1))
            allocate(rep%core_connections(i)%positions(rep%sparse_core_ham(i)%num_elements - 1))

            ! The total number of non-zero elements in row i.
            ! thats a problem for the RDM sampling! I need not only
            ! connected states by the Hamiltonian.. because RDMs are more general...
            rep%core_connections(i)%num_elements = rep%sparse_core_ham(i)%num_elements - 1

            counter = 0
            do j = 1, rep%sparse_core_ham(i)%num_elements
                ! If not the diagonal element.
                if (rep%sparse_core_ham(i)%positions(j) /= i + rep%determ_displs(iProcIndex)) then
                    Ex = 0
                    Ex(1, 1) = nel
                    counter = counter + 1
                    ! The positions of the non-zero and non-diagonal elements in this row i.
                    rep%core_connections(i)%positions(counter) = rep%sparse_core_ham(i)%positions(j)

                    ! for the GUGA implementation this has to be changed in the
                    ! future. but since this routine is only called if we calc.
                    ! RDMs on the fly, i can postpone that until then.. todo
                    ic = FindBitExcitLevel(SpawnedParts(:, i), &
                        temp_store(:, rep%sparse_core_ham(i)%positions(j)))
                    call GetBitExcitation(SpawnedParts(0:NIfD, i), temp_store(0:NIfD, &
                                        rep%sparse_core_ham(i)%positions(j)), Ex, tSign)
                    if (tSign) then
                        ! Odd number of permutations. Minus the excitation level.
                        rep%core_connections(i)%elements(counter) = -ic
                    else
                        ! Even number of permutations. The excitation level.
                        rep%core_connections(i)%elements(counter) = ic
                    end if
                end if
            end do

        end do

        deallocate(temp_store, stat=ierr)
        call LogMemDealloc(t_r, TempStoreTag, ierr)

    end subroutine generate_core_connections