subroutine ReadPopsHeadv4(iunithead, tPop64Bit, tPopHPHF, tPopLz, iPopLenof_Sign, iPopNel, &
iPopAllTotWalkers, PopDiagSft, PopSumNoatHF_out, PopAllSumENum, iPopIter, &
PopNIfD, PopNIfSgn, Popinum_runs, PopNIfFlag, PopNIfTot, read_tau, &
PopBlockingIter, PopRandomHash, read_psingles, read_pparallel, read_ptriples, &
read_nnodes, read_walkers_on_nodes, PopBalanceBlocks)
integer, intent(in) :: iunithead
logical, intent(out) :: tPop64Bit, tPopHPHF, tPopLz
integer, intent(out) :: iPopLenof_sign, iPopNel, iPopIter, PopNIfD
integer, intent(out) :: PopNIfSgn, PopNIfFlag, PopNIfTot
integer, intent(out) :: PopBlockingIter, read_nnodes, Popinum_runs
integer, intent(out) :: PopRandomHash(2056), PopBalanceBlocks
integer(int64), intent(out) :: read_walkers_on_nodes(0:nProcessors - 1)
integer(int64), intent(out) :: iPopAllTotWalkers
real(dp), intent(out) :: PopDiagSft(inum_runs), read_tau, read_psingles, read_ptriples
real(dp), intent(out) :: PopSumNoatHF_out(lenof_sign)
real(dp) :: PopSumNoatHF(2056), PopMultiSft(2056)
real(dp) :: PopMultiSumENum(2056), PopMultiSumNoatHF(2056)
real(dp), intent(out) :: read_pparallel
HElement_t(dp), intent(out) :: PopAllSumENum(inum_runs)
integer :: PopsVersion
!Variables for the namelist
logical :: Pop64Bit, PopLz, PopHPHF
integer :: PopLensign, PopNEl, PopCyc, PopiBlockingIter
integer, parameter :: max_nodes = 30000
integer(int64) :: PopTotwalk
integer(int64), allocatable :: PopWalkersOnNodes(:)
integer :: PopNNodes
real(dp) :: PopSft, PopTau, PopPSingles, PopPParallel, PopGammaSing, PopPTriples
real(dp) :: PopPDoubles, PopPSing_spindiff1, PopPDoub_spindiff1, PopPDoub_spindiff2
real(dp) :: PopGammaDoub, PopGammaTrip, PopGammaOpp, PopGammaPar, PopMaxDeathCpt
real(dp) :: PopTotImagTime, PopSft2, PopParBias
real(dp) :: PopGammaSing_spindiff1, PopGammaDoub_spindiff1, PopGammaDoub_spindiff2
integer :: PopAccumPopsCounter
character(*), parameter :: this_routine = 'ReadPopsHeadv4'
! need dummy read-in variable, since we start from a converged real
! calculation usually! atleast thats the default for now!
HElement_t(dp) :: PopSumENum
integer :: PopNifY
namelist /POPSHEAD/ Pop64Bit, PopHPHF, PopLz, PopLensign, PopNEl, &
PopTotwalk, PopSft, PopSft2, PopSumNoatHF, PopSumENum, &
PopCyc, PopNIfD, PopNifY, PopNIfSgn, PopNIfFlag, PopNIfTot, &
PopTau, PopiBlockingIter, PopRandomHash, &
PopPSingles, PopPSing_spindiff1, PopPDoubles, PopPDoub_spindiff1, PopPDoub_spindiff2, PopPTriples, &
PopPParallel, PopNNodes, PopWalkersOnNodes, PopGammaSing, &
PopGammaDoub, PopGammaTrip, PopGammaOpp, PopGammaPar, PopMaxDeathCpt, &
PopGammaSing_spindiff1, PopGammaDoub_spindiff1, PopGammaDoub_spindiff2, &
PopTotImagTime, Popinum_runs, PopParBias, PopMultiSft, &
PopMultiSumNoatHF, PopMultiSumENum, PopBalanceBlocks, &
tPopAutoAdaptiveShift, tPopScaleBlooms, &
tPopAccumPops, PopAccumPopsCounter
PopsVersion = FindPopsfileVersion(iunithead)
if (PopsVersion /= 4) call stop_all("ReadPopsfileHeadv4", "Wrong popsfile version for this routine.")
PopParBias = 0.0_dp
PopPParallel = 0.0_dp
PopPSingles = 0.0_dp
PopPTriples = 0.0_dp
PopMultiSft = 0.0_dp
PopMaxDeathCpt = 0.0_dp
PopTotImagTime = 0.0_dp
allocate(PopWalkersOnNodes(max_nodes), source=0_int64)
PopBalanceBlocks = -1
PopNNodes = 0
tPopAutoAdaptiveShift = .false.
tPopScaleBlooms = .false.
tPopAccumPops = .false.
if (iProcIndex == root) then
read (iunithead, POPSHEAD)
endif
!Broadcast the read in values from the header to all nodes.
call MPIBCast(Pop64Bit)
call MPIBCast(PopHPHF)
call MPIBCast(PopLz)
call MPIBCast(PopLensign)
call MPIBCast(PopNEl)
call MPIBCast(PopTotwalk)
call MPIBCast(PopSft)
call MPIBCast(PopSft2)
call MPIBCast(PopSumENum)
call MPIBCast(PopCyc)
call MPIBCast(PopNIfD)
call MPIBCast(PopNIfSgn)
call MPIBCast(PopNIfFlag)
call MPIBCast(PopNIfTot)
call MPIBCast(Popinum_runs)
call MPIBCast(PopTau)
call MPIBCast(PopiBlockingIter)
call MPIBCast(PopPSingles)
call MPIBCast(PopPTriples)
if (tReltvy) then
call MPIBCast(PopPDoubles)
call MPIBCast(PopPSing_spindiff1)
call MPIBCast(PopPDoub_spindiff1)
call MPIBCast(PopPDoub_spindiff2)
endif
call MPIBCast(PopPParallel)
call MPIBCast(PopParBias)
call MPIBCast(PopNNodes)
call MPIBcast(PopTotImagTime)
call MPIBCast(PopMultiSft)
call MPIBCast(PopMultiSumNoatHF)
call MPIBCast(PopMultiSumENum)
if (PopNNodes == nProcessors) then
! What is the maximum number of nodes currently supported. We might
! need to update this...
if (PopNNodes > max_nodes) &
call stop_all(this_routine, "Too many processors in POPSFILE. Update &
&max_nodes")
call MPIBCast(PopWalkersOnNodes(1:PopNNodes))
read_walkers_on_nodes(0:PopNNodes - 1) = &
PopWalkersOnNodes(1:PopNNodes)
end if
call MPIBcast(PopGammaSing)
call MPIBcast(PopGammaDoub)
call MPIBcast(PopGammaTrip)
if (tReltvy) then
call MPIBcast(PopGammaSing_spindiff1)
call MPIBcast(PopGammaDoub_spindiff1)
call MPIBcast(PopGammaDoub_spindiff2)
endif
call MPIBcast(PopGammaOpp)
call MPIBCast(PopGammaPar)
call MPIBcast(PopMaxDeathCpt)
call MPIBcast(PopRandomHash)
call MPIBcast(PopBalanceBlocks)
call MPIBCast(tPopAutoAdaptiveShift)
call MPIBCast(tPopScaleBlooms)
call MPIBCast(tPopAccumPops)
tPop64Bit = Pop64Bit
tPopHPHF = PopHPHF
tPopLz = PopLz
iPopLenof_sign = PopLensign
iPopNel = PopNel
iPopAllTotWalkers = PopTotwalk
iPopIter = PopCyc
read_tau = PopTau
PopBlockingIter = PopiBlockingIter
read_psingles = PopPSingles
read_pTriples = PopPTriples
if (abs(PopParBias) > 1.0e-12_dp .and. abs(PopPParallel) < 1.0e-12_dp) then
popPParallel = (PopParBias * par_elec_pairs) &
/ (PopParBias * par_elec_pairs + AB_elec_pairs)
end if
read_pParallel = PopPParallel
read_nnodes = PopNNodes
TotImagTime = PopTotImagTime
! in output generation, these fields are used when tMultiReplicas is set, so this should be
! used here, too (not tReplicaReferencesDiffer), given that the number of runs did not
! change (this would break the read)
if (tMultiReplicas) then
PopDiagSft = PopMultiSft(1:inum_runs)
PopAllSumENum = PopMultiSumENum(1:inum_runs)
PopSumNoatHF_out = PopMultiSumNoatHF(1:lenof_sign)
else
! If we want the walker number to be stable, take the shift
! from the POPSFILE, otherwise, keep the input value.
if (inum_runs == 2) then
if (Popinum_runs == 2) then
PopDiagSft(1) = PopSft
PopDiagSft(inum_runs) = PopSft2
else
! Previously we only had a single run, now we are
! restarting with double run
PopDiagSft(1) = PopSft
PopDiagSft(inum_runs) = PopSft
! I do not think this works, because PopSumNoatHF will not be of size lenof_sign
endif
PopSumNoatHF_out = PopSumNoatHF(1)
else
PopDiagSft(1:inum_runs) = PopSft
! Here, lenof_sign = lenof_sign/inum_runs
PopSumNoatHF_out = PopSumNoatHF(1:lenof_sign)
endif
PopAllSumENum(1:inum_runs) = PopSumENum
end if
call MPIBCast(PopSumNoatHF_out)
! Fill the tau-searching accumulators, to avoid blips in tau etc.
tau_search_stats%gamma_sing = PopGammaSing
tau_search_stats%gamma_doub = PopGammaDoub
tau_search_stats%gamma_trip = PopGammaTrip
if (tReltvy) then
tau_search_stats%gamma_sing_spindiff1 = PopGammaSing_spindiff1
tau_search_stats%gamma_doub_spindiff1 = PopGammaDoub_spindiff1
tau_search_stats%gamma_doub_spindiff2 = PopGammaDoub_spindiff2
endif
tau_search_stats%gamma_opp = PopGammaOpp
tau_search_stats%gamma_par = PopGammaPar
max_death_cpt = PopMaxDeathCpt
end subroutine ReadPopsHeadv4