subroutine singleUpdate(ilut, csf_i, sOrb, excitInfo, posSwitches, negSwitches, &
weightObj, tempExcits, nExcits)
! update function for calculation of all single excitations of a
! given CSF ilut. this implementation takes a general
! probabilistic weight obj as input, to calculate the specific
! probabilisitc weights, determining the different needed ones for
! specific excitations
integer(n_int), intent(in) :: ilut(0:nifguga)
type(CSF_Info_t), intent(in) :: csf_i
integer, intent(in) :: sOrb
type(ExcitationInformation_t), intent(in) :: excitInfo
real(dp), intent(in) :: posSwitches(nSpatOrbs), negSwitches(nSpatOrbs)
type(WeightObj_t) :: weightObj
integer(n_int), intent(inout) :: tempExcits(:, :)
integer, intent(inout) :: nExcits
character(*), parameter :: this_routine = "singleUpdate"
integer(n_int) :: t(0:nifguga)
integer :: iEx, deltaB, gen
real(dp) :: plusWeight, minusWeight, tempWeight, bVal
ASSERT(isProperCSF_ilut(ilut))
ASSERT(sOrb > 0 .and. sOrb < nSpatOrbs)
select case (csf_i%stepvector(sOrb))
! if (csf_i%stepvector(sOrb) == 0) then
! if (isZero(ilut, sOrb)) then
! do nothin actually.. not even change matrix elements
case (0)
return
! else if (csf_i%stepvector(sOrb) == 3) then
! else if (isThree(ilut, sOrb)) then
! only change matrix element to negative one
case (3)
do iEx = 1, nExcits
call update_matrix_element(tempExcits(:, iEx), -1.0_dp, 1)
end do
return
end select
! end if
! to generally use this function i need to define a current generator...
gen = excitInfo%currentGen
bVal = csf_i%B_real(sOrb)
! if 1 or 2 as stepvalue need probWeight
! here the change to the regular excitations is made:
! smth like: have yet to determine, when weight object gets initialized
! but have to check b-value here to avoid division by 0
! or include that in the probWeight calculation?
! b can only be 0 after a 2. which is only a switch possib. for +1
! in the single excitation atleast. and if i set the plus weight to
! 0 and the -1 weight depending on the switch possibilities and the
! end value i could somehow include that?!
! have the bValue restriction now included in the weight calc.
! atleast for pure single excitations... have to additionally do that
! for the more complicated excitation types too
plusWeight = weightObj%proc%plus(posSwitches(sOrb), bVal, weightObj%dat)
minusWeight = weightObj%proc%minus(negSwitches(sOrb), bVal, weightObj%dat)
! have to do some sort of abort_excitations funciton here too
ASSERT(.not. near_zero(plusWeight + minusWeight))
if (csf_i%stepvector(sOrb) == 1) then
! if its a deltaB = -1 this is a switch possib, indepentent
! of the b-vector..
! have to loop over already created excitations
! maybe check if weight is positive before loop, as it is always
! the same and then do an according update..
! is positive weight is 0, negative weight has to be >0
! or else we wouldn be here -> no switches just update -1 branches
if (near_zero(plusWeight)) then
! no switches lead to a nonzero excitation, just update
! matrix element and stay on track
do iEx = 1, nExcits
! need deltaB value
#ifdef DEBUG_
deltaB = getDeltaB(tempExcits(:, iEx))
! check if a -1 is encountert
ASSERT(deltaB == -1)
#endif
call update_matrix_element(tempExcits(:, iEx), &
getSingleMatrixElement(1, 1, -1, gen, bVal), 1)
end do
! when negative weight is 0, positiv weight has to be > 0
! so update positive branches and switch negative ones
else if (near_zero(minusWeight)) then
do iEx = 1, nExcits
t = tempExcits(:, iEx)
deltaB = getDeltaB(t)
if (deltaB == 1) then
! only encode new matrix element in this case
call update_matrix_element(t, &
getSingleMatrixElement(1, 1, deltaB, gen, bVal), 1)
else
! switch also 1 -> 2
clr_orb(t, 2 * sOrb - 1)
set_orb(t, 2 * sOrb)
call setDeltaB(1, t)
call update_matrix_element(t, &
getSingleMatrixElement(2, 1, deltaB, gen, bVal), 1)
end if
tempExcits(:, iEx) = t
end do
! both weights are positiv, staying on track and switching possib
else
! check if deltaB=-1
do iEx = 1, nExcits
! staying on track possible for all excitations. calculate
deltaB = getDeltaB(tempExcits(:, iEx))
tempWeight = extract_matrix_element(tempExcits(:, iEx), 1)
call encode_matrix_element(tempExcits(:, iEx), tempWeight * &
getSingleMatrixElement(1, 1, deltaB, gen, bVal), 1)
if (deltaB == -1) then
! for deltaB=-1 branch a 1 is a switch possibility, if
! weight of positive branch is non-zero.
! new excitations get appended at the end of tempExcits
nExcits = nExcits + 1
t = tempExcits(:, iEx)
! have to change corresponding bits
clr_orb(t, 2 * sOrb - 1)
set_orb(t, 2 * sOrb)
! update matrix element and deltaB flag on new branch
call setDeltaB(1, t)
! forgot the decision with which deltaB i should access
call encode_matrix_element(t, tempWeight * &
getSingleMatrixElement(2, 1, deltaB, gen, bVal), 1)
! and fill into list with correct deltaB
tempExcits(:, nExcits) = t
end if
end do
end if
else
! if it is a 2, it gets a bit more complicated, with the forced
! switches and non-possible stayings...
! if it is a forcesd switch but the negative weight is zero
! i should set the matrix element to zero and move on
! if the weight is not zero but its still a forced switch
! i should update it on the spot, and not add any additional
! excitations. but if switch is possible and weights is
! bigger then zero do it as above
! if b < 1 there is always a forced switch for +1 branch,
! but this situation is not covered by my probWeights...
! hence here i have to additionally check and it could lead to
! no excitations at all.
! did some bullshit thinking here....
! can change down here something now, since i changed the bValue
! inclusion in the weight functions..
if (near_zero(plusWeight)) then
! update on spot and switch
! in this case all +1 branches HAVE to switch, but leave them
! in same position
do iEx = 1, nExcits
t = tempExcits(:, iEx)
if (getDeltaB(tempExcits(:, iEx)) == 1) then
! change stepvectors and branch
set_orb(t, 2 * sOrb - 1)
clr_orb(t, 2 * sOrb)
call setDeltaB(-1, t)
tempWeight = getSingleMatrixElement(1, 2, 1, gen, bVal)
else
! else just update the matrix elements :
tempWeight = getSingleMatrixElement(2, 2, -1, gen, bVal)
end if
call update_matrix_element(t, tempWeight, 1)
tempExcits(:, iEx) = t
end do
else if (near_zero(minusWeight)) then
! update on spot stay
! in this case staying on branch is possible for +1 and a
! -1 branch would have a zero weight, so just update matrix
! but no additional excitations...
do iEx = 1, nExcits
#ifdef DEBUG_
deltaB = getDeltaB(tempExcits(:, iEx))
! ASSERT(deltaB == +1)
#endif
call update_matrix_element(tempExcits(:, iEx), &
getSingleMatrixElement(2, 2, 1, gen, bVal), 1)
end do
else if (minusWeight > 0.0_dp .and. plusWeight > 0.0_dp) then
! update stay on spot and add additional excitations
do iEx = 1, nExcits
! first update matrix elements from staying excitations
deltaB = getDeltaB(tempExcits(:, iEx))
tempWeight = extract_matrix_element(tempExcits(:, iEx), 1)
call update_matrix_element(tempExcits(:, iEx), &
getSingleMatrixElement(2, 2, deltaB, gen, bVal), 1)
if (deltaB == 1) then
! for deltaB +1 branches switch possibility
nExcits = nExcits + 1
! change bits and branch and store at end
t = tempExcits(:, iEx)
set_orb(t, 2 * sOrb - 1)
clr_orb(t, 2 * sOrb)
call setDeltaB(-1, t)
call encode_matrix_element(t, tempWeight * &
getSingleMatrixElement(1, 2, deltaB, gen, bVal), 1)
tempExcits(:, nExcits) = t
end if
end do
else
! something went wrong in this case
call stop_all(this_routine, &
"something went wrong in the excitation generation, shouldnt be here!")
end if
end if
end subroutine singleUpdate