PC_SinglesWeighted_gen_exc Subroutine

private subroutine PC_SinglesWeighted_gen_exc(this, nI, ilutI, nJ, ilutJ, exFlag, ic, ex, tParity, pGen, hel, store, part_type)

Type Bound

PC_SinglesWeighted_t

Arguments

Type IntentOptional Attributes Name
class(PC_SinglesWeighted_t), intent(inout) :: this
integer, intent(in) :: nI(nel)
integer(kind=n_int), intent(in) :: ilutI(0:NIfTot)
integer, intent(out) :: nJ(nel)
integer(kind=n_int), intent(out) :: ilutJ(0:NifTot)
integer, intent(in) :: exFlag
integer, intent(out) :: ic
integer, intent(out) :: ex(2,maxExcit)
logical, intent(out) :: tParity
real(kind=dp), intent(out) :: pGen
real(kind=dp), intent(out) :: hel
type(excit_gen_store_type), intent(inout), target :: store
integer, intent(in), optional :: part_type

Contents


Source Code

    subroutine PC_SinglesWeighted_gen_exc(this, nI, ilutI, nJ, ilutJ, exFlag, ic, &
                                     ex, tParity, pGen, hel, store, part_type)
        class(PC_SinglesWeighted_t), intent(inout) :: this
        integer, intent(in) :: nI(nel), exFlag
        integer(n_int), intent(in) :: ilutI(0:NIfTot)
        integer, intent(out) :: nJ(nel), ic, ex(2, maxExcit)
        integer(n_int), intent(out) :: ilutJ(0:NifTot)
        real(dp), intent(out) :: pGen
        logical, intent(out) :: tParity
        HElement_t(dp), intent(out) :: hel
        type(excit_gen_store_type), intent(inout), target :: store
        integer, intent(in), optional :: part_type

        integer :: elec, src, tgt, i_sg
        real(dp) :: p_src, p_tgt

#ifdef WARNING_WORKAROUND_
        associate(exFlag => exFlag); end associate
        associate(part_type => part_type); end associate
#endif
#ifdef WARNING_WORKAROUND_
        hel = 0.0_dp
#endif
        ic = 1

        if (this%use_lookup) then
            i_sg = this%indexer%lookup_supergroup_idx(store%idx_curr_dets, nI)
        else
            i_sg = this%indexer%idx_nI(nI)
        end if

        select_particle: block
            elec = int(genrand_real2_dSFMT() * nel) + 1
            src = nI(elec)
            p_src = 1._dp / real(nEl, dp)
        end block select_particle

        select_hole: block
            real(dp) :: renorm_tgt
            integer :: dummy, unoccupied(nBasis - nEl)
            integer(n_int) :: ilut_unoccupied(0 : nIfD)
            call this%get_unoccupied(ilutI(0 : nIfD), ilut_unoccupied, unoccupied)
            renorm_tgt = 1._dp - sum(this%A_sampler%get_prob(src, i_sg, nI))
            if (do_direct_calculation(renorm_tgt)) then
                renorm_tgt = sum(this%A_sampler%get_prob(src, i_sg, unoccupied))
            end if

            call this%A_sampler%constrained_sample(&
                 src, i_sg, unoccupied, ilut_unoccupied, renorm_tgt, dummy, tgt, p_tgt)
            if (tgt == 0) then
                call make_invalid()
                return
            end if
        end block select_hole

        pGen = p_src * p_tgt
        call make_single(nI, nJ, elec, tgt, ex, tParity)
        ilutJ = ilutI
        clr_orb(ilutJ, src)
        set_orb(ilutJ, tgt)
        contains

            subroutine make_invalid()
                nJ(1) = 0
                ilutJ = 0_n_int
            end subroutine
    end subroutine PC_SinglesWeighted_gen_exc