function read_pops_general(iunit, PopNel, det_tmp, binary_pops, det_list, &
max_dets, ReadBatch, EndPopsList, &
PopNIfSgn, gdata_read_handler, trimmed_parts) &
result(CurrWalkers)
! General all-purpose pops reading routine.
!
! In: iunit - The popsfile being read from
! ReadBatch - The size of the buffer array to use. Normally
! MaxSpawned, unless specified manually.
integer, intent(in) :: iunit, PopNel, ReadBatch, max_dets
integer, intent(out) :: det_tmp(PopNel)
integer, intent(in) :: PopNIfSgn
logical, intent(in) :: binary_pops
integer(n_int), intent(out) :: det_list(0:NIfTot, max_dets)
integer(int64), intent(in) :: EndPopsList
type(gdata_io_t), intent(in) :: gdata_read_handler
integer(int64) :: CurrWalkers
logical, intent(inout) :: trimmed_parts
character(*), parameter :: this_routine = 'read_pops_general'
logical :: tEOF
logical :: tReadAllPops
integer(MPIArg) :: sendcounts(nNodes), disps(nNodes), recvcount
integer :: PopsInitialSlots(0:nNodes - 1), PopsSendList(0:nNodes - 1)
integer :: batch_size, MaxSendIndex, i, j, nBatches, err, proc
integer(n_int) :: ilut_tmp(0:NIfTot)
real(dp), allocatable :: gdataRead(:, :), gdata(:, :)
real(dp), allocatable :: gdata_tmp(:)
integer :: gdata_size
integer(int64) :: det_attempt, nread
integer :: nelem
real(dp) :: gdata_factor
integer(n_int), allocatable :: BatchRead(:, :)
! If we are on the root processor, we need a buffer array which stores
! the particles as they are read in.
!
! If we are on another processor, allocate a one-element array to avoid
! segfaults in MPIScatter.
allocate (BatchRead(0:NIfTot, merge(ReadBatch, 1, iProcIndex == root)))
! global determinant data (if available)
gdata_size = gdata_read_handler%entry_size()
allocate (gdataRead(gdata_size, ReadBatch))
allocate (gdata_tmp(gdata_size))
if (iProcIndex == root) then
! This is the size of the batches in the above array (only has
! meaning on the root processor)
batch_size = int(real(ReadBatch, dp) / real(nNodes, dp))
! What is the first element in the buffer array that the particles
! for ! each processor are placed in?
forall (i=0:nNodes - 1) PopsInitialSlots(i) = batch_size * i + 1
write (stdout, '(a,i12,a)') "Reading in a maximum of ", ReadBatch, &
" determinants at a time from POPSFILE"
call neci_flush(stdout)
end if
! Keep reading until all of the particles have been read in!
det_attempt = 1
tReadAllPops = .false.
CurrWalkers = 0
pops_norm = 0.0_dp
nBatches = 0
sendcounts = 0
disps = 0
MaxSendIndex = 1
r_loop: do while (.not. tReadAllPops)
! We read in the particles on the root node.
! --> Only do particle receiving on teh other nodes
if (iProcIndex == root) then
! Get ready for reading in the next batch of walkers
nBatches = nBatches + 1
BatchRead(:, :) = 0
PopsSendList(:) = PopsInitialSlots(:)
do while (.true.)
! Read the next entry, and store the walker in ilut_tmp.
! The decoded form is placed in det_tmp
tEOF = read_popsfile_det(iunit, PopNel, binary_pops, &
ilut_tmp, det_tmp, PopNIfSgn, &
.true., nread, gdata_tmp, &
trimmed_parts=trimmed_parts)
det_attempt = det_attempt + nread
! When we have got to the end of the file, we are done.
if (tEOF) then
if (det_attempt <= EndPopsList) then
call stop_all(this_routine, &
"Too few determinants found.")
end if
exit
end if
! Add the contribution from this determinant to the
! norm of the popsfile wave function.
call add_pops_norm_contrib(ilut_tmp)
! Where should this particle be going?
proc = DetermineDetNode(PopNel, det_tmp, 0)
! Store the found determinant in the temporary list,
! and if we have filled up the slot in the list then
! distribute it when it is full.
BatchRead(:, PopsSendList(proc)) = ilut_tmp(:)
gdataRead(:, PopsSendList(proc)) = gdata_tmp(:)
PopsSendList(proc) = PopsSendList(proc) + 1
! If we have filled up the lists, exit the loop so that the
! particles get distributed
if (proc /= nNodes - 1) then
if (PopsInitialSlots(proc + 1) - PopsSendList(proc) < 2) &
exit
else
if (ReadBatch - PopsSendList(proc) < 2) exit
endif
enddo
! Have we read in all of the particles?
if (det_attempt > EndPopsList) tReadAllPops = .true.
! Prep the counts for transmitting the particles to all of
! the nodes.
do j = 0, nNodes - 1
sendcounts(j + 1) = &
int((PopsSendList(j) - PopsInitialSlots(j)) * &
(NIfTot + 1), MPIArg)
disps(j + 1) = &
int((PopsInitialSlots(j) - 1) * (NIfTot + 1), MPIArg)
enddo
MaxSendIndex = (disps(nNodes) + sendcounts(nNodes)) &
/ (nIfTot + 1)
endif
! Now scatter the particles read in to their correct processors.
if (bNodeRoot) then
! How much data goes to each processor?
call MPIScatter(sendcounts, recvcount, err, roots)
! allocate the buffer for the acc/tot spawns
allocate (gdata(gdata_size, recvcount / (NIfTot + 1)))
if (err /= 0) &
call stop_all(this_routine, "MPI scatter error")
! Transfer the walkers.
call MPIScatterV(BatchRead(:, 1:MaxSendIndex), sendcounts, &
disps, &
det_list(:, CurrWalkers + 1:CurrWalkers + 1 + (recvcount / (NIfTot + 1))), &
recvcount, err, Roots)
if (err /= 0) &
call stop_all(this_routine, "MPI scatterV error")
! number of communicated elements
nelem = int(recvcount) .div.int(NIfTot + 1)
! in auto-adaptive shift mode, also communicate the accumulated
! acc/tot spawns so far
! same for the Hij/pgen ratios, if available
! aas data has a different first dimension -> rescale
gdata_factor = real(gdata_size, dp) / real(NIfTot + 1, dp)
recvcount = int(recvcount * gdata_factor, MPIArg)
sendcounts = int(sendcounts * gdata_factor, MPIArg)
call MPIScatterV(gdataRead(1:gdata_size, 1:MaxSendIndex), sendcounts, disps, &
gdata(1:gdata_size, 1:nelem), recvcount, err, Roots)
call gdata_read_handler%read_gdata(gdata, nelem, int(CurrWalkers) + 1)
deallocate (gdata)
CurrWalkers = CurrWalkers + nelem
end if
! Are we done?
call MPIBCast(tReadAllPops)
end do r_loop
write (stdout, "(a,i8)") "Number of batches required to distribute all &
&determinants in POPSFILE: ", nBatches
write (stdout, *) "Number of configurations read in to this process: ", &
CurrWalkers
deallocate (gdataRead)
deallocate (BatchRead)
end function read_pops_general