subroutine calcAllExcitations_single(ilut, csf_i, i, j, excitations, nExcits, t_full)
! function to calculate all possible single excitation for a CSF
! given in (ilut) format and indices (i,j). used to calculated
! H|D> to calculate the projected energy.
integer(n_int), intent(in) :: ilut(0:nifguga)
type(CSF_Info_t), intent(in) :: csf_i
integer, intent(in) :: i, j
integer(n_int), intent(out), allocatable :: excitations(:, :)
integer, intent(out) :: nExcits
logical, intent(in), optional :: t_full
character(*), parameter :: this_routine = "calcAllExcitations_single_ilut"
type(ExcitationInformation_t) :: excitInfo
integer(n_int), allocatable :: tempExcits(:, :)
integer :: ierr, iOrb, iEx, st
real(dp) :: posSwitches(nSpatOrbs), negSwitches(nSpatOrbs), &
plusWeight, minusWeight
HElement_t(dp) :: tmat
type(WeightObj_t) :: weights
logical :: t_full_
def_default(t_full_, t_full, .true.)
ASSERT(i > 0 .and. i <= nSpatOrbs)
ASSERT(j > 0 .and. j <= nSpatOrbs)
! set nExcits to zero by default and only change when excitations
! happen
nExcits = 0
! first check the single particle matrix element, it it is zero leave
! have index it with spin orbitals: assume non-UHF basis
if (t_full_) then
tmat = GetTMatEl(2 * i, 2 * j)
else
tmat = h_cast(1.0_dp)
end if
if (near_zero(tmat)) then
allocate(excitations(0, 0), stat=ierr)
return
end if
! first determine excitation info:
excitInfo = excitationIdentifier(i, j)
! determine the main restrictions of stepvalues depending on generator
! type -> no raising start at d = 3, no lowering start at d = 0 etc
if (excitInfo%gen1 /= 0) then
if (csf_i%stepvector(i) == 3 .or. csf_i%stepvector(j) == 0) then
allocate(excitations(0, 0), stat=ierr)
return
end if
if (t_tJ_model) then
! restrict hops to singly occupied orbitals in t-J model
if (csf_i%Occ_int(i) == 1) then
allocate(excitations(0, 0), stat=ierr)
return
end if
end if
end if
! also check necessary ending restrictions and probability weights
! to check if no excitation is possible
call calcRemainingSwitches_excitInfo_single(csf_i, excitInfo, posSwitches, negSwitches)
! change it here to also use the functions involving
! the WeightObj_t objects.. to efficiently determine
! if excitations have to aborted
! and to make the whole code more general
weights = init_singleWeight(csf_i, excitInfo%fullEnd)
plusWeight = weights%proc%plus(posSwitches(excitInfo%fullStart), &
csf_i%B_real(excitInfo%fullStart), weights%dat)
minusWeight = weights%proc%minus(negSwitches(excitInfo%fullStart), &
csf_i%B_real(excitInfo%fullStart), weights%dat)
! calc total weight functions of all possible branches
! if that is zero no excitation possible..
! do not need exact weight, but only knowledge if it is zero...
! so I do not yet need the bVector
! now we have to really calculate the excitations
! but also use the stochastic weight functions to not
! calculate unnecessary excitations...
! maybe here allocate arrays with the maximum number of
! possible excitations, and fill it up step after step
! 2^|i-j| is more then the maximum number of excitations..
! this could be a problem with this kind of implementation
! actually... -> since for every possible number of
! i,j index pair( and for double even i,j,k,l) this gets
! out of hand for big systems. maybe i have to switch back
! to determine the kind of excitations between to arbitrarily
! given CSFs and then calc. the overlap between them...
! wait and talk to Ali about that...
! should calculate bVector beforehand, and maybe store whole
! b vector for all excitation calculation (i,j) for a given
! CSF, since always the same and needed...
! when i use the remaining switches and end restrictions
! correctly i should not need to delete already created
! excitations. -> so i can fill up the excitation list step
! by step. -> and maybe use top most value as indication of
! delta b value... so i dont need two lists or some more
! advanced data-structure... but when i do it ilut format
! maybe i can use the already provided flag structure which
! is usually used for the sign of the walkers on a given CSF.
! have to figure out how to access and effectively adress
! them
! do it again in this kind of fashion:
st = excitInfo%fullStart
! check compatibility of chosen indices
if ((csf_i%stepvector(st) == 1 .and. near_zero(plusWeight)) .or. &
(csf_i%stepvector(st) == 2 .and. near_zero(minusWeight)) .or. &
near_zero(minusWeight + plusWeight)) then
allocate(excitations(0, 0), stat=ierr)
return
end if
! have to give probabilistic weight object as input, to deal
call createSingleStart(ilut, csf_i, excitInfo, posSwitches, &
negSwitches, weights, tempExcits, nExcits)
! to not call getTmatEl again in createSingleStart loop over
! the atmost two excitations here and multiply with tmat
do iOrb = excitInfo%fullStart + 1, excitInfo%fullEnd - 1
call singleUpdate(ilut, csf_i, iOrb, excitInfo, posSwitches, &
negSwitches, weights, tempExcits, nExcits)
end do
call singleEnd(ilut, csf_i, excitInfo, tempExcits, &
nExcits, excitations)
! encode IC = 1 in the deltB information of the GUGA
! excitation to handle it in the remaining NECI code
! correctly
do iEx = 1, nExcits
call encode_matrix_element(excitations(:, iEx), 0.0_dp, 2)
call update_matrix_element(excitations(:, iEx), tmat, 1)
call setDeltaB(1, excitations(:, iEx))
end do
end subroutine calcAllExcitations_single