subroutine create_particle_with_hash_table(nI_child, ilut_child, child_sign, &
part_type, ilut_parent, iter_data, err)
use hash, only: hash_table_lookup, add_hash_table_entry
integer, intent(in) :: nI_child(nel), part_type
integer(n_int), intent(in) :: ilut_child(0:NIfTot), ilut_parent(0:NIfTot)
real(dp), intent(in) :: child_sign(lenof_sign)
type(fcimc_iter_data), intent(inout) :: iter_data
integer, intent(out) :: err
integer :: proc, ind, hash_val_cd, hash_val, i, run
integer(n_int) :: int_sign(lenof_sign)
real(dp) :: real_sign_old(lenof_sign), real_sign_new(lenof_sign)
real(dp) :: sgn_prod(lenof_sign)
logical :: list_full, tSuccess, allowed_child
integer :: global_position
integer, parameter :: flags = 0
character(*), parameter :: this_routine = 'create_particle_with_hash_table'
!Only one element of child should be non-zero except for real-time evolution
if (.not. (t_rotated_time .or. tVerletSweep)) then
ASSERT((sum(abs(child_sign)) - maxval(abs(child_sign))) < 1.0e-12_dp)
end if
err = 0
! Only one element of child should be non-zero
if (tSimpleInit) then
call stop_all(this_routine, "Cannot use a hash table to the spawned list when using the &
&simple-initiator option.")
end if
call hash_table_lookup(nI_child, ilut_child, nifd, spawn_ht, &
SpawnedParts, ind, hash_val, tSuccess)
if (tSuccess) then
! If the spawned child is already in the spawning array.
! Extract the old sign.
call extract_sign(SpawnedParts(:, ind), real_sign_old)
! If the new child has an opposite sign to that of walkers already
! on the site, then annihilation occurs. The stats for this need
! accumulating.
! in the second real-time spawn loop, i can spawn also to
! determinants, which are actually diagonal particles
! hence i have to update the diag_spawn flag if i annihilate all
! particles eg. and maybe also update the ndied and nborn
! quantities, as this then is not done in the Annihilation if
! no info about the diagonal particles remain ...
! but to know if its a death or a cloning event i have to know
! the original sign in the stored y(n) array... but i dont
! wanna do a lookup in this original list..
! hm: an idea, maybe in the end create a new SpawnedPartsDiag
! array to store the "spawning" from the diagonal step which
! where annhilations in the DirectAnnihilation routine gets
! treated as deaths/births.. -> yes! thats a good idea!
! also there i would be sure to not find the determinants if
! i loop over them, since it essentially is only a copy of the
! worked on y(n) + k1/2 list
! and it would be nicer to seperate those 2 steps as they are
! essentially smth different...
! UPDATE! decided to store the diagonal particles in the 2nd
! RK loop in a seperate DiagParts array -> so no need to
! distinguish here, as only "proper" spawns are treated here!
if (.not. tPrecond) then
sgn_prod = real_sign_old * child_sign
do i = 1, lenof_sign
if (sgn_prod(i) < 0.0_dp) then
iter_data%nannihil(i) = iter_data%nannihil(i) + 2 * min(abs(real_sign_old(i)), abs(child_sign(i)))
end if
end do
end if
! Find the total new sign.
real_sign_new = real_sign_old + child_sign
! Encode the new sign.
call encode_sign(SpawnedParts(:, ind), real_sign_new)
! this is not correctly considered for the real-time or complex
! code .. probably nobody thought about using this in the __cmplx
! implementation..
! Set the initiator flags appropriately.
! If this determinant (on this replica) has already been spawned to
! then set the initiator flag. Also if this child was spawned from
! an initiator, set the initiator flag.
! (There is now an option (tInitCoherentRule = .false.) to turn this
! coherent spawning rule off, mainly for testing purposes).
if (tTruncInitiator) then
if (tInitCoherentRule) then
if (abs(real_sign_old(part_type)) > 1.e-12_dp .or. &
test_flag(ilut_parent, get_initiator_flag(part_type))) then
call set_flag(SpawnedParts(:, ind), get_initiator_flag(part_type))
end if
else
if (test_flag(ilut_parent, get_initiator_flag(part_type))) then
call set_flag(SpawnedParts(:, ind), get_initiator_flag(part_type))
end if
end if
end if
! log the spawn
global_position = ind
else
! Determine which processor the particle should end up on in the
! DirectAnnihilation algorithm.
proc = DetermineDetNode(nel, nI_child, 0)
! Check that the position described by ValidSpawnedList is acceptable.
! If we have filled up the memory that would be acceptable, then
! kill the calculation hard (i.e. stop_all) with a descriptive
! error message.
if (checkValidSpawnedList(proc)) then
err = 1
return
end if
call encode_bit_rep(SpawnedParts(:, ValidSpawnedList(proc)), &
ilut_child(0:nifd), child_sign, flags)
! If the parent was an initiator then set the initiator flag for the
! child, to allow it to survive.
if (t_real_time_fciqmc) then
! for real-time purpose: if the spawn is already populated, also
! set the initiator flag to prevent abort due to the RK reset
if (tTruncInitiator .and. runge_kutta_step == 2) then
! check whether the target is already in CurrentDets
call hash_table_lookup(nI_child, ilut_child, nifd, HashIndex, &
CurrentDets, ind, hash_val_cd, tSuccess)
if (tSuccess) then
call extract_sign(CurrentDets(:, ind), sgn_prod)
! check whether the target is populated in this run
if (.not. is_run_unnocc(sgn_prod, part_type_to_run(part_type))) then
call set_flag(SpawnedParts(:, ValidSpawnedList(proc)), &
get_initiator_flag(part_type))
end if
end if
end if
end if
if (tTruncInitiator) then
if (test_flag(ilut_parent, get_initiator_flag(part_type))) then
call set_flag(SpawnedParts(:, ValidSpawnedList(proc)), &
get_initiator_flag(part_type))
end if
end if
! where to store the global data
global_position = ValidSpawnedList(proc)
call add_hash_table_entry(spawn_ht, ValidSpawnedList(proc), hash_val)
ValidSpawnedList(proc) = ValidSpawnedList(proc) + 1
end if
! store global data
if (tLogNumSpawns) call increase_spawn_counter(SpawnedParts(:, global_position))
! Sum the number of created children to use in acceptance ratio.
! in the rt-fciqmc i have to track the stats of the 2 RK steps
! seperately
! RT_M_Merge: Merge
! rmneci_setup: introduced multirun support, fixed issue in non
! real-time scheme
run = part_type_to_run(part_type)
if (.not. t_real_time_fciqmc .or. runge_kutta_step == 2) then
acceptances(run) = &
acceptances(run) + maxval(abs(child_sign))
end if
end subroutine create_particle_with_hash_table