subroutine init_semi_stochastic(core_in, tStartedFromCoreGround)
! Initialise the semi-stochastic information. This includes enumerating a list of all
! determinants or CSFs in the deterministic space and calculating and storing the resulting
! Hamiltonian matrix elements. The lists which will store the walker amplitude vectors in
! the deterministic space are also allocated.
type(subspace_in) :: core_in
logical, intent(out) :: tStartedFromCoreGround
integer :: i, ierr, run, num_core_runs
integer :: nI(nel)
integer(MPIArg) :: mpi_temp
character(len=*), parameter :: t_r = "init_semi_stochastic"
#ifndef CMPLX_
real(dp), allocatable :: e_values(:)
HElement_t(dp), allocatable :: e_vectors(:, :), gs_vector(:)
real(dp) :: gs_energy
#endif
! If we are load balancing, this gets disabled once semi stochastic
! has been initialised. Therefore we should do a last-gasp load
! adjustment at this point.
if (tLoadBalanceBlocks .and. .not. tFillingStochRDMOnFly) then
call adjust_load_balance(iter_data_fciqmc)
end if
call MPIBarrier(ierr, tTimeIn=.false.)
call set_timer(SemiStoch_Init_Time)
if (t_global_core_space) then
num_core_runs = 1
else
num_core_runs = inum_runs
end if
! Allocate the corespace replicas
allocate(cs_replicas(num_core_runs))
write(stdout, '(/,12("="),1x,a30,1x,12("="))') "Semi-stochastic initialisation"
call neci_flush(stdout)
do run = 1, size(cs_replicas)
associate(rep => cs_replicas(run))
allocate(rep%determ_sizes(0:nProcessors - 1))
allocate(rep%determ_displs(0:nProcessors - 1))
allocate(rep%determ_last(0:nProcessors - 1))
call rep%associate_run(t_global_core_space, run)
rep%determ_sizes = 0_MPIArg
rep%determ_last = 0_MPIArg
if (.not. (tStartCAS .or. core_in%tPops .or. core_in%tDoubles .or. core_in%tCAS .or. core_in%tRAS .or. &
core_in%tOptimised .or. core_in%tLowE .or. core_in%tRead .or. core_in%tMP1 .or. &
core_in%tFCI .or. core_in%tHeisenbergFCI .or. core_in%tHF .or. &
core_in%tPopsAuto .or. core_in%tPopsProportion)) then
call stop_all("init_semi_stochastic", "You have not selected a semi-stochastic core space to use.")
end if
if (.not. tUseRealCoeffs) call stop_all(t_r, "To use semi-stochastic you must also use real coefficients.")
! Call the enumerating subroutines to create all excitations and add these states to
! SpawnedParts on the correct processor. As they do this, they count the size of the
! deterministic space (on their own processor only).
write(stdout, '("Generating the deterministic space...")'); call neci_flush(stdout)
if (core_in%tPopsAuto) then
write(stdout, '("Choosing 10% of initiator space as core space, if this number is larger than 50000, then use 50000")')
call neci_flush(stdout)
! from my understanding npops refers to the total core space size
write(stdout, '("Estimated size of core space:",1X,i5)') int(AllNoInitDets(run) * 0.1_dp)
call neci_flush(stdout)
if (int(AllNoInitDets(run) * 0.1_dp) > 50000) then
core_in%npops = 50000
else
core_in%npops = int(AllNoInitDets(run) * 0.1_dp)
end if
! might also need to check if the total space is too large so that
! tApproxSpace should be used instead...
end if
if (core_in%tPopsProportion) then
write(stdout, '("Choosing ",F7.2,"% of initiator space as core space")') 100 * core_in%npops_proportion
write(stdout, *) "Estimated size of core space: ", int(AllNoInitDets(run) * core_in%npops_proportion)
core_in%npops = max(1, int(AllNoInitDets(run) * core_in%npops_proportion))
end if
if (core_in%tApproxSpace) write(stdout, '(" ... approximately using the factor of",1X,i5)') core_in%nApproxSpace
call generate_space(core_in, run)
! So that all procs store the size of the deterministic spaces on all procs.
mpi_temp = rep%determ_sizes(iProcIndex)
call MPIAllGather(mpi_temp, rep%determ_sizes, ierr)
rep%determ_space_size = sum(rep%determ_sizes)
rep%determ_space_size_int = int(rep%determ_space_size)
! This is now the total size on the replica with the largest space
! Typically, all replicas will have either similar or the same space size
write(stdout, '("Total size of deterministic space:",1X,i8)') rep%determ_space_size
if (rep%determ_space_size > (AllNoInitDets(run) .div. 2_int64) .and. tTruncInitiator) then
write(stdout, *)"WARNING: Total size of deterministic space is greater than&
& 50% of the initiator space."
write(stdout, *)" Reducing the size of the deterministic space is&
& encouraged."
if (iProcIndex == 0) then
write(stderr, *)"WARNING: Total size of deterministic space is greater than&
& 50% of the initiator space."
write(stderr, *)" Reducing the size of the deterministic space is&
& encouraged."
end if
end if
write(stdout, '("Size of deterministic space on this processor:",1X,i8)') rep%determ_sizes(iProcIndex)
call neci_flush(stdout)
call rep%alloc_wf()
! This array will hold the positions of the deterministic states in CurrentDets.
allocate(rep%indices_of_determ_states(rep%determ_sizes(iProcIndex)), stat=ierr)
call LogMemAlloc('indices_of_determ_states', int(rep%determ_sizes(iProcIndex)), &
sizeof_int, t_r, rep%IDetermTag, ierr)
! Calculate the indices in the full vector at which the various processors end.
rep%determ_last(0) = rep%determ_sizes(0)
do i = 1, nProcessors - 1
rep%determ_last(i) = rep%determ_last(i - 1) + rep%determ_sizes(i)
end do
call sort(spawnedparts(0:NIfTot, 1:rep%determ_sizes(iprocindex)), ilut_lt, ilut_gt)
! Do a check that no states are in the deterministic space twice. The list is sorted
! already so simply check states next to each other in the list.
do i = 2, rep%determ_sizes(iProcIndex)
if (all(SpawnedParts(0:nifd, i - 1) == SpawnedParts(0:nifd, i))) then
call decode_bit_det(nI, SpawnedParts(:, i))
write(stdout, '("State found twice:")')
write(stdout, *) SpawnedParts(:, i)
call write_det(stdout, nI, .true.)
call stop_all("init_semi_stochastic", &
"The same state has been found twice in the deterministic space.")
end if
end do
! Store every core determinant from all processors on all processors, in core_space.
call store_whole_core_space(rep)
! Create the hash table to address the core determinants.
call initialise_shared_rht(rep%core_space, rep%determ_space_size_int, rep%core_ht)
if (tWriteCore) call write_core_space(rep)
write(stdout, '("Generating the Hamiltonian in the deterministic space...")'); call neci_flush(stdout)
if (tAllSymSectors .or. tReltvy .or. nOccAlpha <= 1 .or. nOccBeta <= 1 &
.or. tGUGA .or. t_mol_3_body .or. t_non_hermitian_2_body) then
! In the above cases the faster generation is not implemented, so
! use the original algorithm.
call set_timer(SemiStoch_Hamil_Time)
if (tHPHF) then
call calc_determ_hamil_sparse_hphf(rep)
else
call calc_determ_hamil_sparse(rep)
end if
call halt_timer(SemiStoch_Hamil_Time)
write(stdout, '("Total time (seconds) taken for Hamiltonian generation:", f9.3)') &
get_total_time(SemiStoch_Hamil_Time)
else
if (tHPHF) then
call calc_determ_hamil_opt_hphf(rep)
else
call calc_determ_hamil_opt(rep)
end if
end if
#ifndef CMPLX_
if (t_print_core_info) then
! i think i also want information, like the energy and the
! eigenvectors of the core-space
if (t_print_core_hamil) then
call print_basis(rep)
call print_hamiltonian(rep)
end if
root_print "I am before the diagonalization step with", t_non_hermitian_2_body
if (t_non_hermitian_2_body) then
call diagonalize_core_non_hermitian(e_values, e_vectors, rep)
if (t_choose_trial_state) then
gs_energy = e_values(trial_excit_choice(1))
else
gs_energy = e_values(1)
end if
else
call diagonalize_core(gs_energy, gs_vector, rep)
end if
root_print "semi-stochastic space GS energy: ", gs_energy
if (t_print_core_vec) then
root_print "semi-stochastic gs vector:"
call print_vec(gs_vector, "semistoch-vec")
end if
end if
#endif
if (tRDMonFly) call generate_core_connections(rep)
if (tDetermProjApproxHamil) call init_var_space(rep)
! Move the states to CurrentDets.
call add_core_states_currentdet_hash(run)
! If using a trial wavefunction, and that initialisation has already
! been performed, then the current_trial_amps array needs correcting
! after the core states were added and sorted into CurrentDets.
call reinit_current_trial_amps()
! If starting from a popsfile then global_determinant_data will not
! have been initialised, or if in the middle of a calculation then new
! determinants may have been added.
if (tReadPops .or. (semistoch_shift_iter /= 0)) call fill_in_diag_helements()
SpawnedParts = 0_n_int
TotWalkersOld = TotWalkers
tStartedFromCoreGround = .false.
if (tStartCoreGroundState .and. (.not. tReadPops) .and. tStaticCore .and. (.not. tTrialInit)) then
if (t_non_hermitian_2_body) then
call set_timer(SemiStoch_nonhermit_Time)
call start_walkers_from_core_ground_nonhermit(tPrintInfo=.true., run=run)
call halt_timer(SemiStoch_nonhermit_Time)
tStartedFromCoreGround = .true.
write(stdout, '("Total time (seconds) taken for non-hermitian diagonalization:", f9.3)') &
get_total_time(SemiStoch_nonhermit_Time)
else
call set_timer(SemiStoch_Davidson_Time)
call start_walkers_from_core_ground(tPrintInfo=.true., run=run)
call halt_timer(SemiStoch_Davidson_Time)
tStartedFromCoreGround = .true.
write(stdout, '("Total time (seconds) taken for Davidson calculation:", f9.3)') &
get_total_time(SemiStoch_Davidson_Time)
end if
end if
end associate
end do
! Call MPIBarrier here so that Semistoch_Init_Time will give the
! initialisation time for all processors to finish.
call MPIBarrier(ierr, tTimeIn=.false.)
call halt_timer(SemiStoch_Init_Time)
write(stdout, '("Semi-stochastic initialisation complete.")')
write(stdout, '("Time (seconds) taken for semi-stochastic initialisation:", f9.3, /)') &
get_total_time(SemiStoch_Init_Time)
call neci_flush(stdout)
end subroutine init_semi_stochastic