update_gf_overlap Subroutine

    subroutine update_gf_overlap()
        ! subroutine to calculate the overlap of the current y(t) = a_j(a^+_j)(t)y(0)>
        ! time evolved wavefunction to the saved <y(0)|a^+_i(a_i)
        use timing_neci, only: timer, get_total_time
        implicit none
        integer :: idet, nI(nel), det_ind, hash_val, runA, runB, iGf
        real(dp) :: real_sign_1(lenof_sign), real_sign_2(lenof_sign)
        complex(dp) :: overlap(normsize)
        logical :: tDetFound
        real(dp) :: gf_time

        call set_timer(calc_gf_time)

        do iGf = 1, gf_count
            overlap = cmplx(0.0_dp, 0.0_dp, dp)
            do idet = 1, overlap_states(iGf)%nDets

                call extract_sign(overlap_states(iGf)%dets(:, idet), real_sign_1)

                if (IsUnoccDet(real_sign_1)) cycle

                call decode_bit_det(nI, overlap_states(iGf)%dets(:, idet))

                ! search for the hash table associated with the time evolved
                ! wavefunction -> is this already initialized correctly?
                call hash_table_lookup(nI, overlap_states(iGf)%dets(:, idet), nifd, &
                                       HashIndex, CurrentDets, det_ind, hash_val, tDetFound)

                if (tDetFound) then
                    ! both real and imaginary part of the time-evolved wf are required
                    call extract_sign(CurrentDets(:, det_ind), real_sign_2)

                    do runA = 1, inum_runs
                        do runB = 1, inum_runs
                            ! overlap is now treated as complex type
                            ! this only works for the complex code
                            overlap(overlap_index(runA, runB)) = overlap(overlap_index(runA, runB)) &
                                                                 + conjg(cmplx(real_sign_1(min_part_type(runA)), &
                                                                               real_sign_1(max_part_type(runA)), dp)) &
                                                             * cmplx(real_sign_2(min_part_type(runB)), real_sign_2(max_part_type(runB)), dp)
                        end do
                    end do
                end if
            end do

            ! rmneci_setup: the overlap has to be reduced as each proc
            ! only computes its own part
            call MPIReduce(overlap, MPI_SUM, gf_overlap(:, iGf))
        end do

        call halt_timer(calc_gf_time)
        gf_time = get_total_time(calc_gf_time)

    end subroutine update_gf_overlap