calc_pgen_k_space_hubbard_par Function

    function calc_pgen_k_space_hubbard_par(nI, ilutI, ex, ic) result(pgen)
        integer, intent(in) :: nI(nel), ex(:, :), ic
        integer(n_int), intent(in) :: ilutI(0:niftot)
        real(dp) :: pgen
#ifdef DEBUG_
        real(dp) :: test
#endif
        real(dp) :: p_elec, p_orb, cum_arr(nbasis / 2), cum_sum
        integer :: orb_list(nbasis / 2, 2)

        ! check ic:
        if (ic /= 2) then
            pgen = 0.0_dp
            return
        end if

        ! check spin:
        if (.not. (same_spin(ex(1, 1), ex(1, 2)) .and. same_spin(ex(2, 1), ex(2, 2)) .and. &
                   same_spin(ex(1, 1), ex(2, 1)))) then
            pgen = 0.0_dp
            return
        end if

        if (get_ispn(ex(1, 1:2)) == 1) then
            p_elec = 1.0_dp / real(nbasis / 2 - nOccBeta, dp)
        else
            p_elec = 1.0_dp / real(nbasis / 2 - nOccAlpha, dp)
        end if

        call create_ab_list_par_hubbard(nI, ilutI, ex(1, 1:2), orb_list, cum_arr, &
                                        cum_sum, ex(2, 1), p_orb)

        pgen = p_elec * p_orb * 2.0_dp

#ifdef DEBUG_
        call create_ab_list_par_hubbard(nI, ilutI, ex(1, 1:2), orb_list, cum_arr, &
                                        cum_sum, ex(2, 1), test)

        if (abs(test - p_orb) > 1.e-8) then
            print *, "pgen assumption wrong:!"
            print *, "p_orb: ", p_orb
            print *, "test: ", test
            print *, "ex(2,:): ", ex(2, :)
        end if

#endif

    end function calc_pgen_k_space_hubbard_par