subroutine CompressSpawnedList(ValidSpawned, iter_data)
! This sorts and compresses the spawned list to make it easier for the
! rest of the annihilation process. This is not essential, but should
! prove worthwhile.
type(fcimc_iter_data), intent(inout) :: iter_data
integer :: VecInd, ValidSpawned, DetsMerged, i, BeginningBlockDet, FirstInitIndex, CurrentBlockDet
real(dp) :: SpawnedSign(lenof_sign), Temp_Sign(lenof_sign)
integer :: EndBlockDet, part_type, Parent_Array_Ind
integer :: No_Spawned_Parents
integer(kind=n_int), pointer :: PointTemp(:, :)
integer(n_int) :: cum_det(0:IlutBits%len_bcast), temp_det(0:IlutBits%len_bcast)
character(len=*), parameter :: t_r = 'CompressSpawnedList'
type(timer), save :: Sort_time
! We want to sort the list of newly spawned particles, in order for
! quicker binary searching later on. They should remain sorted after
! annihilation between spawned.
if (.not. bNodeRoot) return
Sort_time%timer_name = 'Compress Sort interface'
call set_timer(Sort_time, 20)
call sort(SpawnedParts(0:IlutBits%len_bcast, 1:ValidSpawned), ilut_lt, ilut_gt)
call halt_timer(Sort_time)
if (tHistSpawn) HistMinInd2(1:NEl) = FCIDetIndex(1:NEl)
! First, we compress the list of spawned particles, so that they are
! only specified at most once in each processors list. During this, we
! transfer the particles from SpawnedParts to SpawnedParts2. If we are
! working with complex walkers, we essentially do the same thing twice,
! annihilating real and imaginary particles seperately.
! This is the index in the SpawnedParts2 array to copy the compressed
! walkers into.
VecInd = 1
! BeginningBlockDet will indicate the index of the first entry for a
! given determinant in SpawnedParts.
BeginningBlockDet = 1
DetsMerged = 0
Parent_Array_Ind = 1
Spawned_Parts_Zero = 0
do while (BeginningBlockDet <= ValidSpawned)
! Loop in blocks of the same determinant to the end of the list of
! walkers.
FirstInitIndex = 0
CurrentBlockDet = BeginningBlockDet + 1
do while (CurrentBlockDet <= ValidSpawned)
if (.not. (DetBitEQ(SpawnedParts(:, BeginningBlockDet), &
SpawnedParts(:, CurrentBlockDet)))) exit
! Loop over walkers on the same determinant in SpawnedParts.
CurrentBlockDet = CurrentBlockDet + 1
end do
! EndBlockDet indicates that we have reached the end of the block
! of similar dets.
EndBlockDet = CurrentBlockDet - 1
if (EndBlockDet == BeginningBlockDet) then
! Optimisation: This block only consists of one entry. Simply
! copy it across rather than explicitly searching
! the list.
! If this one entry has no amplitude then don't add it to the
! compressed list, but just cycle.
call extract_sign(SpawnedParts(:, BeginningBlockDet), temp_sign)
if ((sum(abs(temp_sign)) < 1.e-12_dp) .and. &
(.not. (tFillingStochRDMonFly .and. (.not. tNoNewRDMContrib)))) then
DetsMerged = DetsMerged + 1
BeginningBlockDet = CurrentBlockDet
cycle
end if
! Transfer all info to the other array.
SpawnedParts2(:, VecInd) = SpawnedParts(:, BeginningBlockDet)
if (tFillingStochRDMonFly .and. (.not. tNoNewRDMContrib)) then
! SpawnedParts contains the determinants spawned on (Dj),
! and it's parent (Di) plus it's sign (Cj).
! As in | Dj | Di | Ci |
! We then compress multiple occurances of Dj, but these may
! have come from different parents, and we want to keep
! track of all Di's. As we compress SpawnedParts, we
! therefore move all the parents (Di's) into Spawned_Parents.
! If the compressed Dj is at position VecInd in SpawnedParts,
! then Spawned_Parents_Index(1,VecInd) is the starting point
! of it's parents (Di) in Spawned_Parents, and there are
! Spawned_Parents_Index(2,VecInd) entries corresponding to
! this Dj.
if (.not. bit_parent_zero(SpawnedParts(:, BeginningBlockDet))) then
! If the parent determinant is null, the contribution to
! the RDM is zero. No point in doing anything more with it.
Spawned_Parents(0:IlutBitsParent%ind_flag, Parent_Array_Ind) = &
SpawnedParts(IlutBits%ind_parent:IlutBits%ind_parent_flag, BeginningBlockDet)
call extract_sign(SpawnedParts(:, BeginningBlockDet), temp_sign)
! Search to see which sign is non-zero, and therefore
! find which simulation the spawning occured from and to.
! NOTE: it is safe to compare against zero exactly here,
! because all other components will have been set to zero
! exactly and can't have changed at all.
Spawned_Parents(IlutBitsParent%ind_source, Parent_Array_Ind) = 0
do part_type = 1, lenof_sign
if (abs(temp_sign(part_type)) > 1.0e-12_dp) then
Spawned_Parents(IlutBitsParent%ind_source, Parent_Array_Ind) &
= part_type
exit
end if
end do
! in the guga case we also need to transfer the rdm information
if (tGUGA) then
call transfer_stochastic_rdm_info( &
SpawnedParts(:, BeginningBlockDet), &
Spawned_Parents(:, Parent_Array_Ind), &
BitIndex_from=IlutBits, &
BitIndex_to=IlutBitsParent)
end if
! The first nifd of the Spawned_Parents entry is the
! parent determinant, the nifd + 1 entry is the Ci.
! Parent_Array_Ind keeps track of the position in
! Spawned_Parents.
Spawned_Parents_Index(1, VecInd) = Parent_Array_Ind
Spawned_Parents_Index(2, VecInd) = 1
! In this case there is only one instance of Dj - so
! therefore only 1 parent Di.
Parent_Array_Ind = Parent_Array_Ind + 1
else
Spawned_Parents_Index(1, VecInd) = Parent_Array_Ind
Spawned_Parents_Index(2, VecInd) = 0
end if
call extract_sign(SpawnedParts(:, BeginningBlockDet), temp_sign)
if (IsUnoccDet(temp_sign)) then
Spawned_Parts_Zero = Spawned_Parts_Zero + 1
end if
end if
VecInd = VecInd + 1
! Move onto the next block of determinants.
BeginningBlockDet = CurrentBlockDet
cycle ! Skip the rest of this block.
end if
! Reset the cumulative determinant
cum_det = 0_n_int
cum_det(0:nifd) = SpawnedParts(0:nifd, BeginningBlockDet)
if (tPreCond .or. tReplicaEstimates) then
cum_det(IlutBits%ind_hdiag) = SpawnedParts(IlutBits%ind_hdiag, BeginningBlockDet)
end if
if (tFillingStochRDMonFly .and. (.not. tNoNewRDMContrib)) then
! This is the first Dj determinant - set the index for the
! beginning of where the parents for this Dj can be found in
! Spawned_Parents.
Spawned_Parents_Index(1, VecInd) = Parent_Array_Ind
! In this case, multiple Dj's must be compressed, and therefore
! the Di's dealt with as described above. We first just
! initialise the position in the Spawned_Parents array to enter
! the Di's.
Spawned_Parents_Index(2, VecInd) = 0
end if
do i = BeginningBlockDet, EndBlockDet
! if logged, accumulate the number of spawn events
if (tLogNumSpawns) then
cum_det(IlutBits%ind_spawn) = cum_det(IlutBits%ind_spawn) + &
SpawnedParts(IlutBits%ind_spawn, i)
end if
! Annihilate in this block seperately for walkers of different types.
do part_type = 1, lenof_sign
if (tHistSpawn) then
call extract_sign(SpawnedParts(:, i), SpawnedSign)
call extract_sign(SpawnedParts(:, BeginningBlockDet), &
temp_sign)
call HistAnnihilEvent(SpawnedParts, SpawnedSign, &
temp_sign, part_type)
end if
call FindResidualParticle(cum_det, SpawnedParts(:, i), &
part_type, iter_data, VecInd, Parent_Array_Ind)
end do
end do ! Loop over particle type.
! Copy details into the final array.
call extract_sign(cum_det, temp_sign)
if ((sum(abs(temp_sign)) > 1.e-12_dp) .or. (tFillingStochRDMonFly .and. (.not. tNoNewRDMContrib))) then
! Transfer all info into the other array.
! Usually this is only done if the final sign on the compressed
! Dj is not equal to zero. But in the case of the stochastic RDM,
! we are concerned with the sign of Dj in the CurrentDets array,
! not the newly spawned sign. We still want to check if Dj has
! a non-zero Cj in CurrentDets, so we need to carry this Dj
! through to the stage of checking CurrentDets regardless of
! the sign here. Also getting rid of them here would make the
! biased sign of Ci slightly wrong.
! for the stochastic GUGA RDMs we lose the information about
! the RDMs here in the SpawnedParts array
! So we have to use the Spawned_Parents array..
SpawnedParts2(0:NIfTot, VecInd) = cum_det(0:NIfTot)
if (tPreCond .or. tReplicaEstimates) then
SpawnedParts2(IlutBits%ind_hdiag, VecInd) = cum_det(IlutBits%ind_hdiag)
end if
VecInd = VecInd + 1
DetsMerged = DetsMerged + EndBlockDet - BeginningBlockDet
! Spawned_Parts_Zero is the number of spawned parts that are
! zero after compression of the spawned_parts list - and should
! have been removed from SpawnedParts if we weren't calculating
! the RDM - need this for a check later.
if (IsUnoccDet(temp_sign)) Spawned_Parts_Zero = Spawned_Parts_Zero + 1
else
! All particles from block have been annihilated.
DetsMerged = DetsMerged + EndBlockDet - BeginningBlockDet + 1
! This spawned entry will be removed - don't want to store any
! parents. Reset Parent_Array_Ind so that the parents will be
! written over.
! For semi-stochastic simulations, store this state so that we
! can set the flags of the next state to be the same, if it is
! the same state.
temp_det(0:NIfTot) = cum_det(0:NIfTot)
end if
! Move onto the next block of determinants.
BeginningBlockDet = CurrentBlockDet
end do
if (tFillingStochRDMonFly .and. (.not. tNoNewRDMContrib)) No_Spawned_Parents = Parent_Array_Ind - 1
! This is the new number of unique spawned determinants on the processor.
ValidSpawned = ValidSpawned - DetsMerged
if (ValidSpawned /= (VecInd - 1)) then
call stop_all(t_r, "Error in compression of spawned list")
end if
! Want the compressed list in spawnedparts at the end of it - swap
! pointers around.
PointTemp => SpawnedParts2
SpawnedParts2 => SpawnedParts
SpawnedParts => PointTemp
if (tAutoAdaptiveShift) then
PointTemp => SpawnInfo2
SpawnInfo2 => SpawnInfo
SpawnInfo => PointTemp
end if
end subroutine CompressSpawnedList