subroutine CheckPopsParams(tPop64Bit, tPopHPHF, tPopLz, iPopLenof_Sign, iPopNel, &
iPopAllTotWalkers, PopDiagSft, PopSumNoatHF, PopAllSumENum, iPopIter, &
PopNIfD, PopNIfSgn, PopNIfTot, &
WalkerListSize, read_tau, PopBlockingIter, read_psingles, &
read_pparallel, read_ptriples, perturb_ncreate, perturb_nann)
use LoggingData, only: tZeroProjE
logical, intent(in) :: tPop64Bit, tPopHPHF, tPopLz
integer, intent(in) :: iPopLenof_sign, iPopNel, iPopIter, PopNIfD, PopNIfSgn, PopNIfTot
integer, intent(in) :: PopBlockingIter
integer(int64), intent(in) :: iPopAllTotWalkers
real(dp), intent(in) :: PopDiagSft(inum_runs), read_tau
real(dp), intent(in) :: read_psingles, read_pparallel, read_ptriples
real(dp), dimension(lenof_sign), intent(in) :: PopSumNoatHF
HElement_t(dp), intent(in) :: PopAllSumENum(inum_runs)
integer, intent(in) :: perturb_ncreate, perturb_nann
integer, intent(out) :: WalkerListSize
character(len=*), parameter :: this_routine = 'CheckPopsParams'
!Ensure all NIF and symmetry options the same as when popsfile was written out.
#ifdef INT64_
if (.not. tPop64Bit) call stop_all(this_routine, "Popsfile created with 32 bit walkers, but now using 64 bit.")
#else
if (tPop64Bit) call stop_all(this_routine, "Popsfile created with 64 bit walkers, but now using 32 bit.")
#endif
if (tPopHPHF .neqv. tHPHF) call stop_all(this_routine, "Popsfile HPHF and input HPHF not same")
if (tPopLz .neqv. tFixLz) call stop_all(this_routine, "Popsfile Lz and input Lz not same")
if (iPopNEl + perturb_ncreate - perturb_nann /= NEl) call stop_all(this_routine, "The number of electrons &
&in the POPSFILE is not consistent with the number you have asked to run with.")
if (PopNIfD /= NIfD) call stop_all(this_routine, "Popsfile NIfD and calculated NIfD not same")
if (inum_runs == 1) then
!We allow these two values to be different if we're reading in a popsfile fine inum_runs=1 and we want to
!continue the calculation with inum_runs=2
! also allow these values to differ for a real-time fciqmc calc.
! started from a converged groundstate, which is assumed to be
! real valued only
if (iPopLenof_sign /= lenof_sign) then
call stop_all(this_routine, "Popsfile lenof_sign and input lenof_sign not same")
end if
if (PopNIfSgn /= IlutBits%len_pop) then
call stop_all(this_routine, "Popsfile NIfSgn and calculated NIfSgn not same")
end if
endif
if (inum_runs == 1 .and. .not. t_real_time_fciqmc) then
if (PopNIfTot /= NIfTot) then
call stop_all(this_routine, &
"Popsfile NIfTot and calculated NIfTot not same")
end if
endif
if (.not. tWalkContGrow) then
DiagSft = PopDiagSft
end if
if (abs(PopDiagSft(1)) < 1.0e-12_dp) then
!If the popsfile has a shift of zero, continue letting the population grow
tWalkContGrow = .true.
if (inum_runs == 2) then
DiagSft(1) = PopDiagSft(1)
DiagSft(inum_runs) = PopDiagSft(1)
else
DiagSft = PopDiagSft(1)
endif
endif
if (tWalkContGrow) then
tSinglePartPhase = .true.
!If continuing to grow, ensure we can allocate enough memory for what we hope to get the walker population to,
!rather than the average number of determinants in the popsfile.
WalkerListSize = int(max(int(initwalkers, int64), NINT(real(iPopAllTotWalkers, dp) / real(nNodes, dp), int64)))
else
tSinglePartPhase = .false.
WalkerListSize = NINT(real(iPopAllTotWalkers, dp) / real(nNodes, dp))
endif
AllSumENum(1:inum_runs) = PopAllSumENum
AllSumNoatHF(1:lenof_sign) = PopSumNoatHF
PreviousCycles = iPopIter
write (stdout, *) "Number of iterations in previous simulation: ", PreviousCycles
IF (NEquilSteps > 0) THEN
write (stdout, *) "Removing equilibration steps since reading in from POPSFILE."
NEquilSteps = 0
ENDIF
IF (TZeroProjE) THEN
!Reset energy estimator
write (stdout, *) "Resetting projected energy counters to zero..."
AllSumENum = 0.0_dp
AllSumNoatHF = 0
ENDIF
if (near_zero(read_tau)) then
!Using popsfile v.3, where tau is not written out.
!Exit if trying to dynamically search for timestep
if (tau_start_val == possible_tau_start%from_popsfile) then
call stop_all(this_routine, "Cannot read tau from popsfile version <= 3.x")
endif
write (stdout, *) "Old popsfile detected."
write (stdout, *) "Therefore automatic blocking will only start from current run"
iBlockingIter = PreviousCycles
else
if (t_real_time_fciqmc) then
! if reading from a real-time popsfile, also read in tau
! this works because the real-time popsfile is read last
if (tau_start_val == possible_tau_start%from_popsfile) then
call assign_value_to_tau(clamp(read_tau, min_tau, max_tau), 'Initialization from popsfile.')
endif
! also use the adjusted pSingle etc. if provided
if (.not. near_zero(read_psingles)) then
pSingles = read_psingles
pDoubles = 1.0_dp - pSingles
write (stdout, *) " use pSingles and pDoubles provided by POPSFILE!"
write (stdout, *) " pSingles set to: ", pSingles
write (stdout, *) " pDoubles set to: ", pDoubles
end if
tReadPTriples = .false.
if (.not. near_zero(read_ptriples)) then
pTriples = read_ptriples
tReadPTriples = .true.
write (stdout, *) "Using pTriples provided by POPSFILE"
end if
! also do that for pParallel
if (.not. near_zero(read_pparallel)) then
write (stdout, *) " use pParallel provided by POPSFILE: ", read_pparallel
pParallel = read_pparallel
end if
else
!Using popsfile v.4, where tau is written out and read in
which_tau_to_use: if (tau_start_val == possible_tau_start%from_popsfile) then
call assign_value_to_tau(clamp(read_tau, min_tau, max_tau), 'Initialization from popsfile.')
! If we have been searching for tau, we may have been searching
! for psingles (it is done at the same time).
if (allocated(pSinglesIn) .or. allocated(pDoublesIn)) then
write (stdout, *) "using pSingles/pDoubles specified in input file!"
else
if (.not. near_zero(read_psingles)) then
pSingles = read_psingles
if (.not. tReltvy) then
pDoubles = 1.0_dp - pSingles
end if
write (stdout, "(A)") "Using pSingles and pDoubles from POPSFILE: "
write (stdout, "(A,F12.8)") " pSingles: ", pSingles
write (stdout, "(A,F12.8)") " pDoubles: ", pDoubles
end if
end if
if (allocated(pParallelIn)) then
write (stdout, "(A)") "Using pParallel specified in input file!"
else
if (.not. near_zero(read_pparallel)) then
pParallel = read_pparallel
write (stdout, "(A)") "Using pParallel from POPSFILE: "
write (stdout, "(A,F12.8)") " pParallel: ", pParallel
end if
end if
end if which_tau_to_use
if (allocated(pSinglesIn) .or. allocated(pDoublesIn)) then
write (stdout, *) "using pSingles/pDoubles specified in input file!"
else
if (.not. near_zero(read_psingles)) then
pSingles = read_psingles
if (.not. tReltvy) then
pDoubles = 1.0_dp - pSingles
end if
write (stdout, *) "Using read-in pSingles=", pSingles
write (stdout, *) "Using read-in pDoubles=", pDoubles
end if
end if
tReadPTriples = .false.
if (allocated(pTriplesIn)) then
write (stdout, "(A)") "Using pTriples specified in input file!"
else
if (.not. near_zero(read_pTriples)) then
pTriples = read_pTriples
tReadPTriples = .true.
end if
end if
endif
iBlockingIter = PopBlockingIter
endif
end subroutine CheckPopsParams