subroutine singleEnd(ilut, csf_i, excitInfo, tempExcits, nExcits, excitations)
! end function to calculate all single excitations of a given CSF
! ilut
integer(n_int), intent(in) :: ilut(0:nifguga)
type(CSF_Info_t), intent(in) :: csf_i
type(ExcitationInformation_t), intent(in) :: excitInfo
integer(n_int), intent(inout), allocatable :: tempExcits(:, :)
integer, intent(inout) :: nExcits
integer(n_int), intent(out), allocatable :: excitations(:, :)
character(*), parameter :: this_routine = "singleEnd"
integer :: iEx, cnt, ende, ierr, gen, deltaB
integer(n_int) :: t(0:nifguga)
real(dp) :: bVal, tempWeight
! with the correct and consistent use of the probabilistic weight
! functions during the excitation creation, and the assertment that
! the indices and provided CSF are compatible to provide possible
! excitation, no wrong value deltaB branches should arrive here..
! nevertheless do some asserts in debug mode
#ifdef DEBUG_
ASSERT(excitInfo%currentGen /= 0)
ASSERT(isProperCSF_ilut(ilut))
if (excitInfo%currentGen == gen_type%R) then
ASSERT(.not. isZero(ilut, excitInfo%fullEnd))
else if (excitInfo%currentGen == gen_type%L) then
ASSERT(.not. isThree(ilut, excitInfo%fullEnd))
end if
#endif
! for easier readability
ende = excitInfo%fullEnd
! hm... not sure how to treat the end b value... essentially need
! it from one orbital more than -> maybe jsut add or distract one
! if 2/1 step
bVal = csf_i%B_real(ende)
gen = excitInfo%currentGen
! although I do not think I have to check if deltaB values fit,
! still do it for now and check if its really always correct...
select case (csf_i%stepvector(ende))
case (0)
! if it is zero it implies i have a lowering generator, otherwise
! i wouldnt even be here.
do iEx = 1, nExcits
! here depending on deltab b set to 1 or 2
t = tempExcits(:, iEx)
deltaB = getDeltaB(t)
if (deltaB == -1) then
! set alpha bit
set_orb(t, 2 * ende)
tempWeight = getSingleMatrixElement(2, 0, deltaB, gen, bVal)
else
! set beta bit
set_orb(t, 2 * ende - 1)
tempWeight = getSingleMatrixElement(1, 0, deltaB, gen, bVal)
end if
! and fill it in with matrix element and maybe even store
! the calculated matrix element in the signed weight
! entry here... or update weights() -> design decision todo
call update_matrix_element(t, tempWeight, 1)
tempExcits(:, iEx) = t
end do
case (3)
! if d = 3, it implies a raising generator or otherwise we
! would not be here ..
do iEx = 1, nExcits
t = tempExcits(:, iEx)
deltaB = getDeltaB(t)
if (deltaB == -1) then
! clear beta bit
clr_orb(t, 2 * ende - 1)
tempWeight = getSingleMatrixElement(2, 3, deltaB, gen, bVal)
else
! clear alpha bit
clr_orb(t, 2 * ende)
tempWeight = getSingleMatrixElement(1, 3, deltaB, gen, bVal)
end if
call update_matrix_element(t, tempWeight, 1)
tempExcits(:, iEx) = t
end do
case (1)
! d=1 needs a deltaB -1 branch, so delete +1 excitations if they
! happen to get here. alhough that shouldnt happen with the correct
! use of probabilistic weight functions.
cnt = 1
! have to check generator type in this case, as both are possible
! and should do that outside of the loop as always the same...
if (gen == gen_type%R) then
do iEx = 1, nExcits
deltaB = getDeltaB(tempExcits(:, cnt))
if (deltaB == 1) then
! insert last excitation in current place and delete one
! excitation then cycle
tempExcits(:, cnt) = tempExcits(:, nExcits)
nExcits = nExcits - 1
! also make noise for now to indicate smth didnt work
! as expected
cycle
end if
! otherwise deal normally with them
t = tempExcits(:, cnt)
! for raising gen 1->0 : clear beta bit
clr_orb(t, 2 * ende - 1)
! matrix element is just 1 so dont do anything
tempExcits(:, cnt) = t
cnt = cnt + 1
end do
else ! lowering generator
do iEx = 1, nExcits
deltaB = getDeltaB(tempExcits(:, cnt))
if (deltaB == 1) then
tempExcits(:, cnt) = tempExcits(:, nExcits)
nExcits = nExcits - 1
cycle
end if
t = tempExcits(:, cnt)
! for lowering gen: 1 -> 3: set alpha bit
set_orb(t, 2 * ende)
tempWeight = getSingleMatrixElement(3, 1, deltaB, gen, bVal)
call update_matrix_element(t, tempWeight, 1)
tempExcits(:, cnt) = t
cnt = cnt + 1
end do
end if
case (2)
cnt = 1
if (gen == gen_type%R) then
do iEx = 1, nExcits
deltaB = getDeltaB(tempExcits(:, cnt))
if (deltaB == -1) then
tempExcits(:, cnt) = tempExcits(:, nExcits)
nExcits = nExcits - 1
cycle
end if
t = tempExcits(:, cnt)
! for raising: 2 -> 0: clear alpha bit
clr_orb(t, 2 * ende)
! matrix element is just 1 -> so dont do anything
tempExcits(:, cnt) = t
cnt = cnt + 1
end do
else ! lowering gen
do iEx = 1, nExcits
deltaB = getDeltaB(tempExcits(:, cnt))
if (deltaB == -1) then
tempExcits(:, cnt) = tempExcits(:, nExcits)
nExcits = nExcits - 1
cycle
end if
t = tempExcits(:, cnt)
! for lowering gen: 2 -> 3 : set beta bit
set_orb(t, 2 * ende - 1)
tempWeight = getSingleMatrixElement(3, 2, deltaB, gen, bVal)
call update_matrix_element(t, tempWeight, 1)
tempExcits(:, cnt) = t
cnt = cnt + 1
end do
end if
end select
! now have to cut tempExcits to only necessary
! and used entries and output it!
cnt = 1
! fill up the excitations and store matrix element
do iEx = 1, nExcits
! maybe check here again to not have excitations with zero matrix
! elements
if (near_zero(extract_matrix_element(tempExcits(:, iEx), 1))) cycle
tempExcits(:, cnt) = tempExcits(:, iEx)
cnt = cnt + 1
end do
! also update nExcits one final time if something gets cut
nExcits = cnt - 1
allocate(excitations(0:nifguga, nExcits), stat=ierr)
excitations = tempExcits(:, 1:nExcits)
! and get rid of container variables
deallocate(tempExcits)
end subroutine singleEnd