subroutine FindResidualParticle(cum_det, new_det, part_type, iter_data, &
Spawned_No, Parent_Array_Ind)
! This routine is called whilst compressing the spawned list during
! annihilation. It considers the sign and flags from two particles
! on the same determinant, and calculates the correct sign/flags
! for the compressed particle.
!
! --> The information is stored within the first particle in a block
! --> Should be called for real/imaginary particles seperately
integer(n_int), intent(inout) :: cum_det(0:nIfTot)
integer(n_int), intent(in) :: new_det(0:IlutBits%len_bcast)
integer, intent(in) :: part_type, Spawned_No
integer, intent(inout) :: Parent_Array_Ind
type(fcimc_iter_data), intent(inout) :: iter_data
real(dp) :: new_sgn, cum_sgn, updated_sign, sgn_prod
integer :: run
new_sgn = extract_part_sign(new_det, part_type)
cum_sgn = extract_part_sign(cum_det, part_type)
! If the cumulative and new signs for this replica are both non-zero
! then there have been at least two spawning events to this site, so
! set the initiator flag.
! Also set the initiator flag if the new walker has its initiator flag
! set.
if (tTruncInitiator) then
if (tInitCoherentRule) then
if ((abs(cum_sgn) > 1.e-12_dp .and. abs(new_sgn) > 1.e-12_dp) .or. &
test_flag(new_det, get_initiator_flag(part_type))) &
call set_flag(cum_det, get_initiator_flag(part_type))
else
if (test_flag(new_det, get_initiator_flag(part_type))) &
call set_flag(cum_det, get_initiator_flag(part_type))
end if
end if
sgn_prod = cum_sgn * new_sgn
! Update annihilation statistics.
if (.not. tPrecond) then
if (sgn_prod < 0.0_dp) then
run = part_type_to_run(part_type)
if (.not. t_real_time_fciqmc .or. runge_kutta_step == 2) then
Annihilated(run) = Annihilated(run) + 2 * min(abs(cum_sgn), abs(new_sgn))
iter_data%nannihil(part_type) = iter_data%nannihil(part_type) &
+ 2 * min(abs(cum_sgn), abs(new_sgn))
end if
end if
end if
! Update the cumulative sign count.
updated_sign = cum_sgn + new_sgn
call encode_part_sign(cum_det, updated_sign, part_type)
! Obviously only add the parent determinant into the parent array if it is
! actually being stored - and is therefore not zero.
if (tFillingStochRDMonFly .and. (.not. tNoNewRDMContrib)) then
if (.not. DetBitZero( &
new_det(IlutBits%ind_parent:IlutBits%ind_parent + IlutBits%len_orb))) then
if (abs(new_sgn) > 1.e-12_dp) then
! Add parent (Di) stored in SpawnedParts to the parent array.
Spawned_Parents(0:IlutBitsParent%ind_flag, Parent_Array_Ind) = &
new_det(IlutBits%ind_parent:IlutBits%ind_parent_flag)
Spawned_Parents(IlutBitsParent%ind_source, Parent_Array_Ind) = &
part_type
! in the guga implementation we also need to transfer rdm information
if (tGUGA) then
call transfer_stochastic_rdm_info(new_det, &
Spawned_Parents(:, Parent_Array_Ind), &
BitIndex_from=IlutBits, BitIndex_to=IlutBitsParent)
end if
Parent_Array_Ind = Parent_Array_Ind + 1
Spawned_Parents_Index(2, Spawned_No) = &
Spawned_Parents_Index(2, Spawned_No) + 1
end if
end if
end if
end subroutine FindResidualParticle