subroutine second_real_time_spawn(err)
! routine for the second spawning step in the 2nd order RK method
! to create the k2 spawning list. An important change:
! this "spawning" list also has to contain the diagonal death/cloning
! step influence to combine it with the original y(n)
! integer, intent(inout) :: n_determ_states
integer, intent(out) :: err
character(*), parameter :: this_routine = "second_real_time_spawn"
! mimic the most of this routine to the already written first
! spawning step, but with a different death step and without the
! annihilation step at the end!
integer :: idet, parent_flags, nI_parent(nel), unused_flags, ex_level_to_ref, &
ireplica, nspawn, ispawn, nI_child(nel), ic, ex(2, maxExcit), &
ex_level_to_hf
integer(n_int) :: ilut_child(0:niftot)
real(dp) :: parent_sign(lenof_sign), parent_hdiag, prob, child_sign(lenof_sign), &
unused_sign(lenof_sign), unused_rdm_real, diag_sign(lenof_sign)
logical :: tParentIsDeterm, tParentUnoccupied, tParity, break
HElement_t(dp) :: HelGen, parent_hoffdiag
integer :: determ_index
real(dp) :: prefactor, unused_fac, unused_avEx
! declare this is the second runge kutta step
runge_kutta_step = 2
! index for counting deterministic states
determ_index = 1
! prefactor for unbiasing if the number of spawns is cut off
prefactor = 1.0_dp
! use part of nicks code, and remove the parts, that dont matter
do idet = 1, int(TotWalkers)
parent_flags = 0_n_int
! Indicate that the scratch storage used for excitation generation from the
! same walker has not been filled (it is filled when we excite from the first
! particle on a determinant).
fcimc_excit_gen_store%tFilled = .false.
call extract_bit_rep(CurrentDets(:, idet), nI_parent, parent_sign, unused_flags, &
idet, fcimc_excit_gen_store)
ex_level_to_ref = FindBitExcitLevel(iLutRef, CurrentDets(:, idet), max_calc_ex_level)
if (tRef_Not_HF) then
ex_level_to_hf = FindBitExcitLevel(iLutHF_true, CurrentDets(:, idet), &
max_calc_ex_level)
else
ex_level_to_hf = ex_level_to_ref
end if
tParentIsDeterm = check_determ_flag(CurrentDets(:, idet), core_run)
tParentUnoccupied = IsUnoccDet(parent_sign)
if (tParentIsDeterm) then
associate(rep => cs_replicas(core_run))
rep%indices_of_determ_states(determ_index) = idet
rep%partial_determ_vecs(:, determ_index) = parent_sign
end associate
determ_index = determ_index + 1
if (IsUnoccDet(parent_sign)) cycle
end if
if (tTruncInitiator) call CalcParentFlag(idet, parent_flags)
! If this slot is unoccupied (and also not a core determinant) then add it to
! the list of free slots and cycle.
! actually dont need to update this here since nothing gets merged
! at the end.. only k2 gets created
if (tParentUnoccupied) then
! this is unnecessary in the end, the merge is done
! with the original ensemble, with the original FreeSlots
! iEndFreeSlot = iEndFreeSlot + 1
! FreeSlot(iEndFreeSlot) = idet
cycle
end if
! The current diagonal matrix element is stored persistently.
parent_hdiag = det_diagH(idet)
parent_hoffdiag = det_offdiagH(idet)
! UPDATE: call this routine anyway to update info on noathf
! and noatdoubs, for the intermediate step
call SumEContrib(nI_parent, ex_level_to_ref, parent_sign, CurrentDets(:, idet), &
parent_hdiag, parent_hoffdiag, 1.0_dp, tPairedReplicas, idet)
! If we're on the Hartree-Fock, and all singles and
! doubles are in the core space, then there will be
! no stochastic spawning from this determinant, so
! we can the rest of this loop.
if (ss_space_in%tDoubles .and. ex_level_to_hf == 0 .and. tDetermHFSpawning) cycle
! If this condition is not met (if all electrons have spin up or all have spin down)
! then there will be no determinants to spawn to, so don't attempt spawning.
! thats a really specific condition.. shouldnt be checked each
! cycle.. since this is only input dependent..
if (.not. tGZero) then ! skip this if we only want the corespace-evolution
do ireplica = 1, lenof_sign
call decide_num_to_spawn(parent_sign(ireplica), AvMCExcits, nspawn)
!call merge_spawn(nspawn,prefactor)
do ispawn = 1, nspawn
! Zero the bit representation, to ensure no extraneous data gets through.
ilut_child = 0_n_int
call generate_excitation(nI_parent, CurrentDets(:, idet), nI_child, &
ilut_child, exFlag, ic, ex, tParity, prob, &
HElGen, fcimc_excit_gen_store)
! If a valid excitation.
if (.not. IsNullDet(nI_child)) then
call encode_child(CurrentDets(:, idet), ilut_child, ic, ex)
ilut_child(IlutBits%ind_flag) = 0_n_int
if (tSemiStochastic) then
break = check_semistoch_flags(ilut_child, nI_child, part_type_to_run(ireplica), tParentIsDeterm)
if (break) cycle
end if
! unbias if the number of spawns was truncated
child_sign = attempt_create(nI_parent, CurrentDets(:, idet), parent_sign, &
nI_child, ilut_child, prob, HElGen, ic, ex, tParity, &
ex_level_to_ref, ireplica, unused_sign, AvMCExcits, unused_rdm_real, 1.0_dp)
child_sign = prefactor * child_sign
else
child_sign = 0.0_dp
end if
! If any (valid) children have been spawned.
if ((any(abs(child_sign) > EPS)) .and. (ic /= 0) .and. (ic <= 2)) then
! not quite sure about how to collect the child
! stats in the new rt-fciqmc ..
call new_child_stats_normal(second_spawn_iter_data, CurrentDets(:, idet), &
nI_child, ilut_child, ic, ex_level_to_ref, &
child_sign, parent_flags, ireplica)
call create_particle_with_hash_table(nI_child, ilut_child, child_sign, &
ireplica, CurrentDets(:, idet), &
second_spawn_iter_data, err)
end if ! If a child was spawned.
end do ! Over mulitple particles on same determinant.
end do ! Over the replicas on the same determinant.
end if
! If this is a core-space determinant then the death step is done in
! determ_projection.
if (.not. tParentIsDeterm) then
! essentially have to treat the diagonal clone/death step like
! a spawning step and store the result into the spawned array..
! call walker_death_spawn ()!iter_data_fciqmc, nI_parent, &
! ilut_parent, parent_hdiag, parent_sign, idet, ex_level_to_ref)
! also not quite sure how to do the child_stats here ...
! or in general for now in the rt-fciqmc
! have to write all new book-keeping routines i guess..
! i also have to change the sign convention here, since in
! the annihilation(and in the spawing routine above) the
! particles get merged with a + instead of the - in the
! original death routine for a "normal" death
diag_sign = -attempt_die_realtime(parent_hdiag, &
parent_sign, ex_level_to_ref)
if (any(abs(diag_sign) > EPS)) then
call create_diagonal_as_spawn(CurrentDets(:, idet), &
diag_sign, second_spawn_iter_data)
end if
end if
end do ! Over all determinants.
! the deterministic time evoulution is performed here according to the
! RK scheme (annihilation is done separately as deterministic_annihilation)
if (tSemiStochastic) call real_time_determ_projection()
end subroutine second_real_time_spawn