subroutine calc_projected_hamil(nvecs, krylov_array, krylov_ht, ndets, h_matrix, partial_vecs, full_vecs, h_diag)
use bit_rep_data, only: IlutBits
use bit_reps, only: decode_bit_det
use CalcData, only: ss_space_in, tDetermHFSpawning
use constants
use DetBitOps, only: return_ms, FindBitExcitLevel
use enumerate_excitations, only: generate_connection_normal, generate_connection_kpnt
use enumerate_excitations, only: init_generate_connected_space
use FciMCData, only: fcimc_excit_gen_store, exFlag, SpawnVecKP, SpawnVecKP2
use FciMCData, only: SpawnVec, SpawnVec2, SpawnedParts
use core_space_util, only: cs_replicas
use FciMCData, only: SpawnedParts2, InitialSpawnedSlots, ValidSpawnedList, ll_node
use FciMCData, only: spawn_ht, subspace_hamil_time, iLutHF_True, max_calc_ex_level
use FCIMCData, only: HFConn, Hii
use global_det_data, only: det_diagH
use hash, only: clear_hash_table
use kp_fciqmc_data_mod, only: tSemiStochasticKPHamil, tExcitedStateKP, av_mc_excits_kp
use kp_fciqmc_data_mod, only: tExactHamilSpawning, kp_hamil_exact_frac
use procedure_pointers, only: generate_excitation, encode_child, get_spawn_helement
use semi_stoch_procs, only: is_core_state, check_determ_flag
use semi_stoch_procs, only: determ_projection_kp_hamil
use SystemData, only: nbasis, tAllSymSectors, nOccAlpha, nOccBeta, nel
use SystemData, only: tKPntSym, tUseBrillouin
use timing_neci, only: set_timer, halt_timer
use util_mod, only: stochastic_round, near_zero
integer, intent(in) :: nvecs
integer(n_int), intent(in) :: krylov_array(0:, :)
type(ll_node), pointer, intent(in) :: krylov_ht(:)
integer, intent(in) :: ndets
real(dp), intent(out) :: h_matrix(:, :)
real(dp), allocatable, intent(inout), optional :: partial_vecs(:, :), full_vecs(:, :)
real(dp), intent(in), optional :: h_diag(:)
integer :: idet, ispawn, nspawn, i, j, ex_level_to_hf
integer :: determ_ind, flag_ind, ic, ex(2, maxExcit), ms_parent
integer :: ex_flag, nStore(6)
integer, allocatable :: excit_gen(:)
integer :: nI_parent(nel), nI_child(nel)
integer(n_int) :: ilut_child(0:NIfTot), ilut_parent(0:NIfTot)
real(dp) :: prob, tot_pop, h_diag_elem
real(dp), dimension(lenof_all_signs) :: child_sign, parent_sign
integer(n_int) :: int_sign(lenof_all_signs)
logical :: tChildIsDeterm, tParentIsDeterm, tParentUnoccupied, tParity
logical :: tNearlyFull, tAllFinished, all_excits_found, tTempUseBrill
logical :: tFinished
HElement_t(dp) :: HElGen, HEl
call set_timer(subspace_hamil_time)
h_matrix(:, :) = 0.0_dp
ilut_parent = 0_n_int
flag_ind = nifd + lenof_all_signs + 1
ValidSpawnedList = InitialSpawnedSlots
call clear_hash_table(spawn_ht)
tNearlyFull = .false.
tFinished = .false.
determ_ind = 1
if (.not. tExcitedStateKP) then
! We want SpawnedParts and SpawnedParts2 to point to the 'wider' spawning arrays which
! hold all the signs of *all* the Krylov vectors. This is because SendProcNewParts uses
! SpawnedParts and SpawnedParts2.
SpawnedParts => SpawnVecKP
SpawnedParts2 => SpawnVecKP2
end if
do idet = 1, ndets
! The 'parent' determinant from which spawning is to be attempted.
ilut_parent(0:nifd) = krylov_array(0:nifd, idet)
ilut_parent(IlutBits%ind_flag) = krylov_array(flag_ind, idet)
! Indicate that the scratch storage used for excitation generation from the
! same walker has not been filled (it is filled when we excite from the first
! particle on a determinant).
fcimc_excit_gen_store%tFilled = .false.
call decode_bit_det(nI_parent, ilut_parent)
int_sign = krylov_array(IlutBits%ind_pop:IlutBits%ind_pop + lenof_all_signs - 1, idet)
parent_sign = transfer(int_sign, parent_sign)
tot_pop = sum(abs(parent_sign))
tParentIsDeterm = check_determ_flag(ilut_parent)
tParentUnoccupied = all(abs(parent_sign) < 1.e-16_dp)
! If this determinant is in the deterministic space then store the relevant
! data in arrays for later use.
if (tParentIsDeterm) then
! Add the amplitude to the deterministic vector.
partial_vecs(:, determ_ind) = parent_sign
determ_ind = determ_ind + 1
end if
if (tParentUnoccupied) cycle
if (ss_space_in%tDoubles .and. tDetermHFSpawning) then
ex_level_to_hf = FindBitExcitLevel(iLutHF_true, ilut_parent, max_calc_ex_level)
if (ex_level_to_hf == 0) cycle
end if
if (tAllSymSectors) then
ms_parent = return_ms(ilut_parent)
! If this condition is met (if all electrons have spin up or all have spin
! down) then there will be no determinants to spawn to, so don't attempt spawning.
if (abs(ms_parent) == nbasis / 2) cycle
nOccAlpha = (nel + ms_parent) / 2
nOccBeta = (nel - ms_parent) / 2
end if
if (tExactHamilSpawning .and. av_mc_excits_kp * tot_pop >= real(HFConn, dp) * kp_hamil_exact_frac) then
! We no longer want to treat these deterministic states by the
! usual deterministic projection. So set the corresponding
! amplitude in the deterministic vector to 0 and add in the
! diagonal element separately.
if (tParentIsDeterm) then
partial_vecs(:, determ_ind - 1) = 0.0_dp
if (present(h_diag)) then
h_diag_elem = h_diag(idet) + Hii
else
h_diag_elem = det_diagH(idet) + Hii
end if
child_sign = calc_amp_kp_hamil(parent_sign, 1.0_dp, 1.0_dp, h_diag_elem)
call create_particle_kp_estimates(nI_parent, ilut_parent, child_sign, tNearlyFull)
end if
call init_generate_connected_space(nI_parent, ex_flag, all_excits_found, ex, excit_gen, &
nstore, tTempUseBrill)
do while (.true.)
if (tKPntSym) then
call generate_connection_kpnt(nI_parent, ilut_parent, nI_child, ilut_child, ex_flag, &
all_excits_found, nstore, excit_gen, hel=HEl)
else
call generate_connection_normal(nI_parent, ilut_parent, nI_child, ilut_child, ex_flag, &
ex, all_excits_found, hel=HEl)
end if
if (all_excits_found) exit
child_sign = calc_amp_kp_hamil(parent_sign, 1.0_dp, 1.0_dp, real(HEl, dp))
call create_particle_kp_estimates(nI_child, ilut_child, child_sign, tNearlyFull)
if (tNearlyFull) then
call add_in_contribs(nvecs, krylov_array, krylov_ht, tFinished, tAllFinished, 1.0_dp, h_matrix)
call clear_hash_table(spawn_ht)
tNearlyFull = .false.
end if
end do
if (tKPntSym) then
tUseBrillouin = tTempUseBrill
deallocate(excit_gen)
end if
else
nspawn = stochastic_round(av_mc_excits_kp * tot_pop)
do ispawn = 1, nspawn
! Zero the bit representation, to ensure no extraneous data gets through.
ilut_child = 0_n_int
call generate_excitation(nI_parent, ilut_parent, nI_child, &
ilut_child, exFlag, ic, ex, tParity, prob, &
HElGen, fcimc_excit_gen_store)
! If a valid excitation.
if (.not. IsNullDet(nI_child)) then
call encode_child(ilut_parent, ilut_child, ic, ex)
! If spawning is both too and from the core space, abort it.
if (tSemiStochasticKPHamil .and. tParentisDeterm) then
if (is_core_state(ilut_child, nI_child)) cycle
end if
HEl = get_spawn_helement(nI_parent, nI_child, ilut_parent, ilut_child, ic, ex, &
tParity, HElGen)
child_sign = calc_amp_kp_hamil(parent_sign, prob, av_mc_excits_kp * tot_pop, real(HEl, dp))
else
child_sign = 0.0_dp
end if
! If any (valid) children have been spawned.
if (.not. (all(near_zero(child_sign)) .or. ic == 0 .or. ic > 2)) then
call create_particle_kp_estimates(nI_child, ilut_child, child_sign, tNearlyFull)
if (tNearlyFull) then
call add_in_contribs(nvecs, krylov_array, krylov_ht, tFinished, tAllFinished, 1.0_dp, h_matrix)
call clear_hash_table(spawn_ht)
tNearlyFull = .false.
end if
end if ! If a child was spawned.
end do ! Over mulitple spawns from the same determinant.
end if
end do ! Over all determinants.
tFinished = .true.
do
call add_in_contribs(nvecs, krylov_array, krylov_ht, tFinished, tAllFinished, 1.0_dp, h_matrix)
if (tAllFinished) exit
end do
call calc_hamil_contribs_diag(nvecs, krylov_array, ndets, h_matrix, h_diag)
if (tSemiStochasticKPHamil) then
associate(rep => cs_replicas(core_run))
call determ_projection_kp_hamil(partial_vecs, full_vecs, rep)
end associate
call calc_hamil_contribs_semistoch(nvecs, krylov_array, h_matrix, partial_vecs)
end if
! Symmetrise the projected Hamiltonian.
do i = 1, nvecs
do j = 1, i - 1
h_matrix(i, j) = h_matrix(j, i)
end do
end do
if (.not. tExcitedStateKP) then
! Now let SpawnedParts and SpawnedParts2 point back to their
! original arrays.
SpawnedParts => SpawnVec
SpawnedParts2 => SpawnVec2
end if
call halt_timer(subspace_hamil_time)
end subroutine calc_projected_hamil