subroutine Gen_Hist_ExcDjs(iLutnI)
! This uses GenExcitations3 in symexcit3.F90 to generate all the
! possible either single or double excitations from D_i, finds the
! processor they would be on if occupied, and puts them in the
! SingExcDjs array according to that processor.
use DetBitOps, only: EncodeBitDet
use hist_data, only: AllHistogram
use load_balance_calcnodes, only: DetermineDetNode
use LoggingData, only: RDMExcitLevel
use Parallel_neci, only: nProcessors
use rdm_data, only: Sing_ExcList, Doub_ExcList, ExcNorm, OneEl_Gap, TwoEl_Gap
use rdm_data, only: Sing_ExcDjs, Doub_ExcDjs
use rdm_data, only: one_rdms, two_rdm_spawn
use rdm_filling, only: fill_diag_1rdm, fill_spawn_rdm_diag
use SymExcit3, only: GenExcitations3
use SymExcit4, only: GenExcitations4, ExcitGenSessionType
use SystemData, only: nel
integer(n_int), intent(in) :: iLutnI(0:NIfTot)
integer(n_int) :: iLutnJ(0:NIfTot)
integer, dimension(lenof_sign) :: HistPos
real(dp), dimension(lenof_sign) :: RealHistPos
integer :: ExcitMat3(2, 2), nI(NEl), nJ(NEl), Proc, FlagsDi
integer :: a, b, exflag
logical :: tAllExcitFound, tParity
real(dp) :: SignDi(lenof_sign), full_sign(1)
type(ExcitGenSessionType) :: session
! Unfortunately uses the decoded determinant - might want to look at this.
call extract_bit_rep(iLutnI, nI, RealHistPos, FlagsDi)
HistPos = int(RealHistPos)
SignDi(1) = AllHistogram(1, HistPos(1)) / ExcNorm
SignDi(lenof_sign) = AllHistogram(1, HistPos(1)) / ExcNorm
if (RDMExcitLevel == 1) then
call fill_diag_1rdm(one_rdms, nI, SignDi)
else
full_sign = SignDi(1) * SignDi(lenof_sign)
call fill_spawn_rdm_diag(two_rdm_spawn, nI, full_sign)
end if
! Zeros in ExcitMat3 starts off at the first single excitation.
ExcitMat3(:, :) = 0
! This becomes true when all the excitations have been found.
tAllExcitFound = .false.
do while (.not. tAllExcitFound)
exflag = 1
! Passed out of here is the singly excited determinant, nJ.
! Information such as the orbitals involved in the excitation and
! the parity is also found in this step, we are not currently
! storing this, and it is re-calculated later on (after the
! determinants are passed to the relevant processor) - but the speed
! of sending this information vs recalculating it will be tested.
! RDMExcitLevel is passed through, if this is 1, only singles are
! generated, if it is 2 only doubles are found.
if (tReltvy) then
call GenExcitations4(session, nI, nJ, exFlag, ExcitMat3(:, :), tParity, tAllExcitFound, .true.)
else
call GenExcitations3(nI, iLutnI, nJ, exflag, ExcitMat3(:, :), tParity, tAllExcitFound, .true.)
end if
if (tAllExcitFound) exit
iLutnJ(:) = 0
call EncodeBitDet(nJ, iLutnJ)
Proc = DetermineDetNode(nel, nJ, 0)
! This will return a value between 0 -> nProcessors-1.
Sing_ExcDjs(:, Sing_ExcList(Proc)) = iLutnJ(:)
Sing_ExcList(Proc) = Sing_ExcList(Proc) + 1
! Want a quick test to see if arrays are getting full.
if (Sing_ExcList(Proc) > nint(OneEl_Gap * (Proc + 1))) then
write(stdout, *) 'Proc', Proc
write(stdout, *) 'Sing_ExcList', Sing_ExcList
write(stdout, *) 'No. spaces for each proc', nint(OneEl_Gap)
call Stop_All('GenExcDjs', 'Too many excitations for space available.')
end if
end do
if (RDMExcitLevel /= 1) then
ExcitMat3(:, :) = 0
! Zeros in ExcitMat3 starts off at the first single excitation.
tAllExcitFound = .false.
! This becomes true when all the excitations have been found.
do while (.not. tAllExcitFound)
exflag = 2
! Passed out of here is the doubly excited determinant, nJ.
! Information such as the orbitals involved in the excitation
! and the parity is also found in this step, we are not currently
! storing this, and it is re-calculated later on (after the
! determinants are passed to the relevant processor) - but the
! speed of sending this information vs recalculating it will be
! tested. RDMExcitLevel is passed through, if this is 1, only
! singles are generated, if it is 2 only doubles are found.
if (tReltvy) then
call GenExcitations4(session, nI, nJ, exFlag, ExcitMat3(:, :), tParity, tAllExcitFound, .true.)
else
call GenExcitations3(nI, iLutnI, nJ, exflag, ExcitMat3(:, :), tParity, tAllExcitFound, .true.)
end if
if (tAllExcitFound) exit
iLutnJ(:) = 0
call EncodeBitDet(nJ, iLutnJ)
Proc = DetermineDetNode(nel, nJ, 0)
! This will return a value between 0 -> nProcessors-1.
Doub_ExcDjs(:, Doub_ExcList(Proc)) = iLutnJ(:)
! All the double excitations from this particular nI are being
! stored in Doub_ExcDjs.
Doub_ExcList(Proc) = Doub_ExcList(Proc) + 1
! Want a quick test to see if arrays are getting full.
if (Doub_ExcList(Proc) > nint(TwoEl_Gap * (Proc + 1))) then
write(stdout, *) 'Proc', Proc
write(stdout, *) 'Doub_ExcList', Doub_ExcList
write(stdout, *) 'No. spaces for each proc', nint(TwoEl_Gap)
call Stop_All('GenExcDjs', 'Too many excitations for space available.')
end if
end do
end if
end subroutine Gen_Hist_ExcDjs