setup_pchb_sampler_conditional Subroutine

private subroutine setup_pchb_sampler_conditional(this)

Type Bound

GugaAliasSampler_t

Arguments

Type IntentOptional Attributes Name
class(GugaAliasSampler_t), intent(inout) :: this

Contents


Source Code

    subroutine setup_pchb_sampler_conditional(this)
        class(GugaAliasSampler_t), intent(inout) :: this
        integer :: i, j, ij, a, b, ab
        integer(int64) :: ijMax, abMax
        integer(int64), allocatable :: excit_info(:)
        real(dp), allocatable :: w(:)

        ijMax = fuseIndex(nSpatOrbs, nSpatOrbs)
        abMax = fuseIndex(nSpatOrbs, nSpatOrbs)
        ! just to be save also code up a setup function with the
        ! necessary if statements, so i can check if my optimized
        ! sampler is correct (and if it is even faster..)
        ! weights per pair
        allocate(w(abMax), source = 0.0_dp)
        ! excitation information per pair:
        allocate(excit_info(abMax))

        ! allocate: all samplers have the same size
        call this%alias_sampler%shared_alloc(int(ijMax), int(abMax), "GUGA-pchb")
        ! todo: do the same for the excit_info array!
        call this%new_info_table(ijMax, abMax)

        ! the encode_excit_info function encodes as: E_{ai}E_{bj)
        ! but for some index combinations we have to change the input so
        ! actually E_{aj}E_{bi} is encoded! convention to use only
        ! certain type of excit-types
        do i = 1, nSpatOrbs
            do j = i, nSpatOrbs
                w = 0.0_dp
                ij = fuseIndex(i,j)
                excit_info = 0_int64
                do a = 1, nSpatOrbs
                    do b = a, nSpatOrbs
                        ! shoud i point group symmetry restrctions here?
                        ! this would avoid unnecessary if statements..

                        if (RandExcitSymLabelProd(&
                                        SpinOrbSymLabel(2*i), SpinOrbSymLabel(2*j)) &
                                /= RandExcitSymLabelProd(&
                                        SpinOrbSymLabel(2*a), SpinOrbSymLabel(2*b))) then
                            cycle
                        end if

                        ab = fuseIndex(a,b)
                        call get_weight_and_info(i, j, a, b, w(ab), excit_info(ab))
                    end do
                end do
                call this%alias_sampler%setup_entry(ij, root, w)
                call this%set_info_entry(ij, excit_info)

            end do
        end do
    end subroutine setup_pchb_sampler_conditional