calcSingleOverlapRaising Subroutine

public subroutine calcSingleOverlapRaising(ilut, csf_i, excitInfo, excitations, nExcits, posSwitches, negSwitches)

Arguments

Type IntentOptional Attributes Name
integer(kind=n_int), intent(in) :: ilut(0:nifguga)
type(CSF_Info_t), intent(in) :: csf_i
type(ExcitationInformation_t), intent(inout) :: excitInfo
integer(kind=n_int), intent(out), allocatable :: excitations(:,:)
integer, intent(out) :: nExcits
real(kind=dp), intent(in) :: posSwitches(nSpatOrbs)
real(kind=dp), intent(in) :: negSwitches(nSpatOrbs)

Contents


Source Code

    subroutine calcSingleOverlapRaising(ilut, csf_i, excitInfo, excitations, &
                                        nExcits, posSwitches, negSwitches)
        ! special function to calculate all excitations for a single overlap
        ! double excitations with only raising generators
        integer(n_int), intent(in) :: ilut(0:nifguga)
        type(CSF_Info_t), intent(in) :: csf_i
        type(ExcitationInformation_t), intent(inout) :: excitInfo
        integer(n_int), intent(out), allocatable :: excitations(:, :)
        integer, intent(out) :: nExcits
        real(dp), intent(in) :: posSwitches(nSpatOrbs), negSwitches(nSpatOrbs)
        character(*), parameter :: this_routine = "calcSingleOverlapRaising"

        integer(n_int), allocatable :: tempExcits(:, :)
        integer(n_int) :: t(0:nifguga)
        integer :: i, iEx, deltaB, se, en
        type(WeightObj_t) :: weightObj
        real(dp) :: plusWeight, minusWeight

        ASSERT(isProperCSF_ilut(ilut))
        ASSERT(excitInfo%typ == excit_type%single_overlap_raising)
        ASSERT(excitInfo%firstGen == gen_type%R)
        ASSERT(excitInfo%lastGen == gen_type%R)

        se = excitInfo%secondStart
        en = excitInfo%fullEnd
        ! create weight object here
        weightObj = init_singleOverlapRaising(csf_i, se, en, negSwitches(se), &
                                              posSwitches(se), csf_i%B_real(se))

        excitInfo%currentGen = excitInfo%firstGen
        ! have to make specific single start to correctly adress the weights
        call createSingleStart(ilut, csf_i, excitInfo, posSwitches, negSwitches, &
                               weightObj, tempExcits, nExcits)

        ! loop until overlap site
        do i = excitInfo%fullStart + 1, se - 1
            call singleUpdate(ilut, csf_i, i, excitInfo, posSwitches, &
                              negSwitches, weightObj, tempExcits, nExcits)
        end do

        i = excitInfo%secondStart

        ! update weights
        weightObj = init_singleWeight(csf_i, en)
        plusWeight = weightObj%proc%plus(posSwitches(se), &
                                         csf_i%B_real(se), weightObj%dat)
        minusWeight = weightObj%proc%minus(negSwitches(se), &
                                           csf_i%B_real(se), weightObj%dat)

        ASSERT(plusWeight + minusWeight > 0.0_dp)

        ! do special stuff at single overlap raising site
        if (csf_i%stepvector(se) == 1) then
            ! in this case there should come a deltaB=+1 branch, nevertheless
            ! check for now..
            ! have to include probabilistic weights too... which are the normal
            ! single excitations weights at this point
            do iEx = 1, nExcits
                deltaB = getDeltaB(tempExcits(:, iEx))
                if (deltaB == 1) then
                    call stop_all(this_routine, &
                                  "there should not be a DeltaB= +1 at a RR(1) single overlap site")
                end if
                ! only need to check switching weight, since on track has to
                ! be greater than 0, to be even here TODO: think about that!
                ! not quite sure yeah.. since this is a switch possib..
                ! jsut to be save to a check like above
                ! update: to be save change that to check probs

                if (near_zero(minusWeight)) then
                    ! do only switch
                    t = tempExcits(:, iEx)
                    ! change 1 -> 2
                    set_orb(t, 2 * se)
                    clr_orb(t, 2 * se - 1)

                    call setDeltaB(1, t)

                    ! no change in matrix elements
                    tempExcits(:, iEx) = t

                else
                    ! staying on -1 doesnt change anything, so just check if
                    ! a switch is also possible in this case: and add at end

                    if (plusWeight > 0.0_dp) then
                        ! and all deltaB=-1 are branch possibs. so add an excitations
                        nExcits = nExcits + 1
                        t = tempExcits(:, iEx)
                        ! change 1 -> 2
                        set_orb(t, 2 * se)
                        clr_orb(t, 2 * se - 1)

                        call setDeltaB(1, t)

                        ! no change in matrix elements
                        tempExcits(:, nExcits) = t
                    end if
                end if
            end do

        else if (csf_i%stepvector(se) == 2) then
            ! in this case always a switch, and i b allows also a stay
            ! how are the probs here...

            if (near_zero(plusWeight)) then
                ! just switch
                ! only switches in place
                do iEx = 1, nExcits
                    if (getDeltaB(tempExcits(:, iEx)) == -1) then
                        call stop_all(this_routine, &
                                      "there should be no deltaB=-1 at RR(2) single overlap site")
                    end if

                    t = tempExcits(:, iEx)
                    ! change 2 -> 1
                    clr_orb(t, 2 * se)
                    set_orb(t, 2 * se - 1)

                    call setDeltaB(-1, t)

                    tempExcits(:, iEx) = t
                end do

            else
                ! staying doesnt change anything just check if switch is
                ! additionally possible

                if (minusWeight > 0.0_dp) then
                    ! stay and add switch
                    do iEx = 1, nExcits
                        ! still check delta B just to be sure
                        deltaB = getDeltaB(tempExcits(:, iEx))
                        if (deltaB == -1) then
                            call stop_all(this_routine, &
                                          "there should be no deltaB=-1 at RR(2) single overlap site")
                        end if
                        nExcits = nExcits + 1

                        t = tempExcits(:, iEx)
                        ! change 2 -> 1
                        clr_orb(t, 2 * se)
                        set_orb(t, 2 * se - 1)

                        call setDeltaB(-1, t)

                        tempExcits(:, nExcits) = t
                    end do
                end if
            end if
        end if

        ! continue with secon region normally
        ! have to reset weight object here to use "normal" single excitation
        ! also no change in generator type since alike generators
        excitInfo%currentGen = excitInfo%lastGen

        do i = excitInfo%secondStart + 1, excitInfo%fullEnd - 1
            call singleUpdate(ilut, csf_i, i, excitInfo, posSwitches, negSwitches, &
                              weightObj, tempExcits, nExcits)
        end do

        ! normal end step
        call singleEnd(ilut, csf_i, excitInfo, tempExcits, nExcits, excitations)

    end subroutine calcSingleOverlapRaising