subroutine gen_ab_cum_list_1_1(csf_i, occ_orbs, cum_arr, excit_arr, orb_arr)
! specific routine when the occupaton of the already picked orbitals
! is 2 -> still think if i really want to outpout a cum_arr of lenght
! nBasis -> nSpatOrbs would be better and doable... !! todo
! only difference to 3_3 case is, that (a) can be (i,j), but which also
! forces (b) to be the other of (i,j), but there is the additional
! contraint then, that there has to be a possible switch between them!
! so, probably a good idea to check if there is a possible switch
! between i and j first, and only then allow a = (i,j)
type(CSF_Info_t), intent(in) :: csf_i
integer, intent(in) :: occ_orbs(2)
real(dp), intent(out) :: cum_arr(nSpatOrbs)
type(ExcitationInformation_t), intent(out), optional :: excit_arr(nSpatOrbs)
integer, intent(out), optional :: orb_arr(nSpatOrbs)
character(*), parameter :: this_routine = "gen_ab_cum_list_1_1"
integer :: ki(3), kj(3), n_id(2), orb_a, ka(3), kb(3), orb_b, st, en
integer :: z_ind
real(dp) :: cum_sum, contrib
logical :: tSwitch
type(ExcitationInformation_t) :: excitInfo
! for legacy compatibility
z_ind = ubound(kPointToBasisFn, 3)
ki = G1(occ_orbs(1))%k
kj = G1(occ_orbs(2))%k
cum_sum = 0.0_dp
orb_arr = 0
cum_arr = 0.0_dp
n_id = gtID(occ_orbs)
tSwitch = .true.
! determine if there is a possible switch between i and j, to check if
! a = i/j should be allowed
! if the symmetry assumptions below(to be tested!) hold, i do not
! need to care about possible switches even here!
! meh except that a fullstart into fullstop mixed is possible..
if (csf_i%stepvector(n_id(1)) == csf_i%stepvector(n_id(2))) then
if (csf_i%stepvector(n_id(1)) == 1) then
if (count_alpha_orbs_ij(csf_i, n_id(1), n_id(2)) == 0) tSwitch = .false.
else if (csf_i%stepvector(n_id(1)) == 2) then
if (count_beta_orbs_ij(csf_i, n_id(1), n_id(2)) == 0) tSwitch = .false.
end if
end if
! due to k-point symmetry alot of excitation types:
! eg. if i and j are seperate a and b also have to be!
! not true! ki + kj = 2*ka can defo be!
do orb_a = 1, n_id(1) - 1
contrib = 0.0_dp
excitInfo%valid = .false.
! avoid doubly occupied orbitals
if (csf_i%stepvector(orb_a) /= 3) then
ka = G1(2 * orb_a)%k
kb = ki + kj - ka
call MomPbcSym(kb, nBasisMax)
orb_b = gtID(KPointToBasisFn(kb(1), kb(2), z_ind, 1))
! just check if k-point restriction works as i thought
ASSERT(orb_b /= n_id(1) .and. orb_b /= n_id(2))
if (orb_a == orb_b) then
! (a) must be empty and there must be a switch possible!
if (csf_i%stepvector(orb_a) == 0 .and. tSwitch) then
! have to make contrib twice as high here since
! there is no second chance to pick it the other way
! around..
contrib = 2.0_dp
! _RR_(ab) > ^RR(i) > ^R(j)
excitInfo = assign_excitInfo_values_double( &
excit_type%fullstart_raising, &
gen_type%R, gen_type%R, gen_type%R, gen_type%R, gen_type%R, &
orb_a, n_id(1), orb_a, n_id(2), orb_a, orb_a, n_id(1), n_id(2), &
0, 2, 1.0_dp, 1.0_dp)
end if
else
! check guga restrictions todo
if (csf_i%stepvector(orb_b) /= 3) then
! no additional restrictions or?...
contrib = 1.0_dp
if (orb_b > n_id(2)) then
! _R(a) > _LR(i) > ^RL(j) > ^L(b)
excitInfo = assign_excitInfo_values_double( &
excit_type%double_R_to_L, &
gen_type%L, gen_type%R, gen_type%R, gen_type%R, gen_type%L, &
orb_a, n_id(2), orb_b, n_id(1), orb_a, n_id(1), n_id(2), orb_b, &
0, 4, 1.0_dp, 1.0_dp)
else if (orb_b < n_id(1)) then
! check if a == b
! check maxima
st = min(orb_a, orb_b)
en = max(orb_a, orb_b)
! _R(min) > _RR(max) > ^RR(i) > ^R(j)
excitInfo = assign_excitInfo_values_double( &
excit_type%double_raising, &
gen_type%R, gen_type%R, gen_type%R, gen_type%R, gen_type%R, &
en, n_id(2), st, n_id(1), st, en, n_id(1), n_id(2), &
0, 4, 1.0_dp, 1.0_dp)
else
! b is between i and j
! _R(a) > _LR(i) > ^LR(b) > ^R(j)
excitInfo = assign_excitInfo_values_double( &
excit_type%double_R_to_L_to_R, &
gen_type%L, gen_type%R, gen_type%R, gen_type%R, gen_type%R, &
orb_a, n_id(2), orb_b, n_id(1), orb_a, n_id(1), orb_b, n_id(2), &
0, 4, 1.0_dp, 1.0_dp)
end if
end if
end if
end if
cum_sum = cum_sum + contrib
cum_arr(orb_a) = cum_sum
excit_arr(orb_a) = excitInfo
orb_arr(orb_a) = orb_b
end do
! check if orbital (i) is valid to choose from
if (tSwitch) then
! do not have to deal specifically with orbital 1 ...or?
! do i need to consider kb restriction... yes i do think so..
kb = kj
! and also need no mapping back here?
! hm..
call MomPbcSym(kb, nBasisMax)
orb_b = gtID(KPointToBasisFn(kb(1), kb(2), z_ind, 1))
! b has to be (j)!
ASSERT(orb_b == n_id(2))
contrib = 1.0_dp
excitInfo = assign_excitInfo_values_double( &
excit_type%fullstart_stop_mixed, &
gen_type%L, gen_type%R, gen_type%R, gen_type%R, gen_type%R, &
n_id(1), n_id(2), n_id(2), n_id(1), n_id(1), n_id(1), n_id(2), n_id(2), &
0, 2, 1.0_dp, 1.0_dp)
cum_sum = cum_sum + contrib
cum_arr(n_id(1)) = cum_sum
excit_arr(n_id(1)) = excitInfo
orb_arr(n_id(1)) = n_id(2)
else
! otherwise orb (i) off-limits
if (n_id(1) == 1) then
cum_arr(1) = 0.0_dp
else
cum_arr(n_id(1)) = cum_arr(n_id(1) - 1)
end if
end if
do orb_a = n_id(1) + 1, n_id(2) - 1
contrib = 0.0_dp
excitInfo%valid = .false.
! avoid doubly occupied orbitals
if (csf_i%stepvector(orb_a) /= 3) then
ka = G1(2 * orb_a)%k
kb = ki + kj - ka
call MomPbcSym(kb, nBasisMax)
orb_b = gtID(KPointToBasisFn(kb(1), kb(2), z_ind, 1))
ASSERT(orb_b /= n_id(1) .and. orb_b /= n_id(2))
if (orb_a == orb_b) then
if (csf_i%stepvector(orb_a) == 0 .and. tSwitch) then
contrib = 2.0_dp
! _L(i) > ^LR_(ab) > ^R(j)
excitInfo = assign_excitInfo_values_double( &
excit_type%single_overlap_L_to_R, &
gen_type%L, gen_type%R, gen_type%L, gen_type%L, gen_type%R, &
orb_a, n_id(1), orb_a, n_id(2), n_id(1), orb_a, orb_a, n_id(2), &
0, 2, 1.0_dp, 1.0_dp, 1)
end if
else
! check guga restrictions todo
if (csf_i%stepvector(orb_b) /= 3) then
! no additional restrictions or?...
contrib = 1.0_dp
if (orb_b < n_id(1)) then
! _R(b) > _LR(i) > ^LR(a) > ^R(j)
excitInfo = assign_excitInfo_values_double( &
excit_type%double_R_to_L_to_R, &
gen_type%L, gen_type%R, gen_type%R, gen_type%R, gen_type%R, &
orb_b, n_id(2), orb_a, n_id(1), orb_b, n_id(1), orb_a, n_id(2), &
0, 4, 1.0_dp, 1.0_dp)
else if (orb_b > n_id(2)) then
! _L(i) > _RL(a) > ^RL(j) > ^L(b)
excitInfo = assign_excitInfo_values_double( &
excit_type%double_L_to_R_to_L, &
gen_type%L, gen_type%R, gen_type%L, gen_type%L, gen_type%L, &
orb_b, n_id(1), orb_a, n_id(2), n_id(1), orb_a, n_id(2), orb_b, &
0, 4, 1.0_dp, 1.0_dp)
else
st = min(orb_a, orb_b)
en = max(orb_a, orb_b)
! _L(i) > _RL(min) > ^LR(max) > ^R(j)
excitInfo = assign_excitInfo_values_double( &
excit_type%double_L_to_R, &
gen_type%L, gen_type%R, gen_type%L, gen_type%L, gen_type%R, &
en, n_id(1), st, n_id(2), n_id(1), st, en, n_id(2), &
0, 4, 1.0_dp, 1.0_dp)
end if
end if
end if
end if
cum_sum = cum_sum + contrib
cum_arr(orb_a) = cum_sum
excit_arr(orb_a) = excitInfo
orb_arr(orb_a) = orb_b
end do
! fill in orbital j, depending if a switch is possible
! check if orbital (i) is valid to choose from
if (tSwitch) then
!todo: think about, that this is actually the same as the
! RL > RL considered above, just picked in different order..
! do not have to deal specifically with orbital 1 ...or?
! do i need to consider kb restriction... yes i do think so..
kb = ki
call MomPbcSym(kb, nBasisMax)
orb_b = gtID(KPointToBasisFn(kb(1), kb(2), z_ind, 1))
ASSERT(orb_b == n_id(1))
! b has to be (j)!
contrib = 1.0_dp
excitInfo = assign_excitInfo_values_double( &
excit_type%fullstart_stop_mixed, &
gen_type%L, gen_type%R, gen_type%R, gen_type%R, gen_type%R, &
n_id(1), n_id(2), n_id(2), n_id(1), n_id(1), n_id(1), n_id(2), n_id(2), &
0, 2, 1.0_dp, 1.0_dp)
cum_sum = cum_sum + contrib
cum_arr(n_id(2)) = cum_sum
excit_arr(n_id(2)) = excitInfo
orb_arr(n_id(2)) = n_id(1)
else
cum_arr(n_id(2)) = cum_arr(n_id(2) - 1)
excitInfo%valid = .false.
excit_arr(n_id(2)) = excitInfo
end if
do orb_a = n_id(2) + 1, nSpatOrbs
contrib = 0.0_dp
excitInfo%valid = .false.
if (csf_i%stepvector(orb_a) /= 3) then
ka = G1(2 * orb_a)%k
kb = ki + kj - ka
call MomPbcSym(kb, nBasisMax)
! there is no "spin" restriction in the guga case
! so both possible b orbs have to be checked
! its actually NOT possible that ka = ki !! so do not have
! to check that case!
orb_b = gtID(KPointToBasisFn(kb(1), kb(2), z_ind, 1))
! this should work as we pick beta orbital first
ASSERT(orb_b /= n_id(1) .and. orb_b /= n_id(2))
if (orb_a == orb_b) then
if (csf_i%stepvector(orb_a) == 0 .and. tSwitch) then
contrib = 2.0_dp
! _L(i) > _LL(j) > ^LL^(ab)
excitInfo = assign_excitInfo_values_double( &
excit_type%fullstop_lowering, &
gen_type%L, gen_type%L, gen_type%L, gen_type%L, gen_type%L, &
orb_a, n_id(1), orb_b, n_id(2), n_id(1), n_id(2), orb_a, orb_a, &
0, 2, 1.0_dp, 1.0_dp)
end if
else
if (csf_i%stepvector(orb_b) /= 3) then
! only then its a possible excitation
contrib = 1.0_dp
! check type of excitation
if (orb_b < n_id(1)) then
! _R(b) > _LR(i) > ^RL(j) > ^L(a)
excitInfo = assign_excitInfo_values_double( &
excit_type%double_R_to_L, &
gen_type%L, gen_type%R, gen_type%R, gen_type%R, gen_type%L, &
orb_b, n_id(2), orb_a, n_id(1), orb_b, n_id(1), n_id(2), orb_a, &
0, 4, 1.0_dp, 1.0_dp)
else if (orb_b > n_id(2)) then
st = min(orb_a, orb_b)
en = max(orb_a, orb_b)
! _L(i) > _LL(j) > ^LL(min) > ^L(max)
excitInfo = assign_excitInfo_values_double( &
excit_type%double_lowering, &
gen_type%L, gen_type%L, gen_type%L, gen_type%L, gen_type%L, &
st, n_id(1), en, n_id(2), n_id(1), n_id(2), st, en, &
0, 4, 1.0_dp, 1.0_dp)
else
! _L(i) > _RL(b) > ^RL(j) > ^L(a)
excitInfo = assign_excitInfo_values_double( &
excit_type%double_L_to_R_to_L, &
gen_type%L, gen_type%R, gen_type%L, gen_type%L, gen_type%L, &
orb_a, n_id(1), orb_b, n_id(2), n_id(1), orb_b, n_id(2), orb_a, &
0, 4, 1.0_dp, 1.0_dp)
end if
end if
end if
end if
cum_sum = cum_sum + contrib
cum_arr(orb_a) = cum_sum
excit_arr(orb_a) = excitInfo
orb_arr(orb_a) = orb_b
end do
end subroutine gen_ab_cum_list_1_1