gen_cum_list_guga_single_1 Subroutine

private subroutine gen_cum_list_guga_single_1(nI, csf_i, orb_i, cc_i, cum_arr)

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: nI(nel)
type(CSF_Info_t), intent(in) :: csf_i
integer, intent(in) :: orb_i
integer, intent(in) :: cc_i
real(kind=dp), intent(out) :: cum_arr(OrbClassCount(cc_i))

Contents


Source Code

    subroutine gen_cum_list_guga_single_1(nI, csf_i, orb_i, cc_i, cum_arr)
        ! specific single orbital picker if stepvector of electron (i) is 1
        integer, intent(in) :: nI(nel)
        type(CSF_Info_t), intent(in) :: csf_i
        integer, intent(in) :: orb_i, cc_i
        real(dp), intent(out) :: cum_arr(OrbClassCount(cc_i))
        character(*), parameter :: this_routine = "gen_cum_list_guga_single_1"

        integer :: nOrb, i, label_index, j, n_id(nEl), id_i, &
                   lower, upper, s_orb, spin
        real(dp) :: cum_sum, hel

        ! if d(i) = 1 -> i can only pick d(a) = 1 if there is a switch possib
        ! d(j) = 2 inbetween! -> include that in cummulative probabilities!
        ! still not quite sure how to effectively include that...
        nOrb = OrbClassCount(cc_i)
        label_index = SymLabelCounts2(1, cc_i)

        n_id = gtID(nI)
        id_i = gtID(orb_i)

        cum_sum = 0.0_dp

        ! do some precalcing here to determine, which orbitals to exclude due
        ! to GUGA restrictions?..
        call find_switches(csf_i, id_i, lower, upper)

        if (is_beta(orb_i)) then
            spin = 1
        else
            spin = 0
        end if

        do i = 1, nOrb
            s_orb = sym_label_list_spat(label_index + i - 1)

            if (s_orb == id_i) then
                cum_arr(i) = cum_sum
                cycle
            end if

            hel = 0.0_dp

            select case (csf_i%stepvector(s_orb))

                ! include here the guga restrictions...
            case (0)
                ! no restrictions if 0 since both branches allowed

                ! single particle matrix element:
                hel = hel + abs(GetTMatEl(orb_i, 2 * s_orb - spin))

                ! do the loop over all the other electrons
                ! (is this always symmetrie allowed?..)

                ! i know both occupation numbers! of start and end!
                ! -> calc. the specific U(i,i,j,i)n(i) and U(i,j,j,j)n(j)
                ! specifically
                ! in case d=0 -> no influence for WR_(s_orb) and WL^(s_orb)
                ! and it can only be one of those
                ! also n(id_i) = 1 which is also 0 -> so exclude this from
                ! loop below

                ! if depending on the type of generator
                ! either the topCont sign (if L)
                ! or botCont sign (if R)
                ! is unknown! -> so for every k < st if R
                ! or k > en if L the sign of the two-particle matrix
                ! elements is unkown! and has to be added in terms of
                ! absolute value!

                ! problem is, if i do not know, all the signs correctly i
                ! have to add up all matrix element contributions with
                ! there absolute value.. otherwise the order of summing
                ! influences the result, and thus cannot be true!
                ! "good" thing is, i do not need to consider, so many
                ! different possibilities

                if (.not. t_mixed_hubbard) then
                    do j = 1, nEl

                        ! todo: finish all contributions later for now only do
                        ! those which are the same for all
                        ! exclude initial orbital, since this case gets
                        ! contributed already outside of loop over electrons!
                        ! but only spin-orbital or spatial??
                        if (n_id(j) == id_i) cycle
                        hel = hel + abs(get_umat_el(id_i, n_id(j), s_orb, n_id(j)))

                        ! now depending on generator and relation of j to
                        ! st and en -> i know sign or don't

                        hel = hel + abs(get_umat_el(id_i, n_id(j), n_id(j), s_orb))

                    end do
                end if

            case (1)
                ! here i have to somehow find out if there is a
                ! (2) between s_orb and id_i
                ! how to do that?
                ! could also use the loop over nEl to check if there
                ! is a switch between (i) and (a) ->  and set matrix
                ! element to 0 otherwise... -> would make effor O(n) again
                if (s_orb < lower .or. s_orb > upper) then
                    ! only allowed if possible switch

                    hel = hel + abs(GetTMatEl(orb_i, 2 * s_orb - spin))

                    ! also need the weight contribution at start

                    if (.not. t_mixed_hubbard) then
                        hel = hel + abs(get_umat_el(id_i, s_orb, s_orb, s_orb))

                        ! do the loop over all the other electrons
                        ! (is this always symmetrie allowed?..)

                        do j = 1, nEl

                            ! todo: finish all contributions later for now only do
                            ! those which are the same for all
                            if (n_id(j) == id_i .or. n_id(j) == s_orb) cycle
                            hel = hel + abs(get_umat_el(id_i, n_id(j), s_orb, n_id(j)))

                            hel = hel + abs(get_umat_el(id_i, n_id(j), n_id(j), s_orb))

                        end do
                    end if
                end if

            case (2)
                ! no restrictions for 2 -> 1 excitations
                hel = hel + abs(GetTMatEl(orb_i, 2 * s_orb - spin))
                ! do the loop over all the other electrons
                ! (is this always symmetrie allowed?..)

                if (.not. t_mixed_hubbard) then
                    hel = hel + abs(get_umat_el(id_i, s_orb, s_orb, s_orb))

                    do j = 1, nEl
                        ! todo: finish all contributions later for now only do
                        ! those which are the same for all
                        if (n_id(j) == id_i .or. n_id(j) == s_orb) cycle
                        hel = hel + abs(get_umat_el(id_i, n_id(j), s_orb, n_id(j)))

                        hel = hel + abs(get_umat_el(id_i, n_id(j), n_id(j), s_orb))

                    end do
                end if

            case (3)
                ! do nothing in this case!

            case default
                ! should not be here!
                call stop_all(this_routine, "stepvalue /= {0,1,2,3}! something is wrong!")

            end select
            cum_sum = cum_sum + abs_l1(hel)
            cum_arr(i) = cum_sum

        end do

    end subroutine gen_cum_list_guga_single_1