subroutine store_krylov_vec(ivec)
use FciMCData, only: TotWalkers, CurrentDets
use global_det_data, only: det_diagH
use hash, only: hash_table_lookup, add_hash_table_entry
use Loggingdata, only: tPrintDataTables
use semi_stoch_procs, only: check_determ_flag
use SystemData, only: nel
use util_mod, only: int_fmt
integer, intent(in) :: ivec
integer :: idet, iamp, sign_ind, flag_ind, hash_val, det_ind
integer :: nI(nel)
integer(n_int) :: temp, int_sign(lenof_sign_kp)
logical :: tDetFound, tCoreDet
real(dp) :: amp_fraction, real_sign(lenof_sign_kp)
character(len=*), parameter :: t_r = "store_krylov_vec"
if (tPrintDataTables) then
write(stdout, '(a71)', advance='no') "# Adding the current walker configuration to the Krylov vector array..."
call neci_flush(stdout)
end if
! The index of the first element referring to the sign, for this ivec.
sign_ind = nifd + lenof_sign_kp * (ivec - 1) + 1
flag_ind = nifd + lenof_all_signs + 1
! Loop over all occupied determinants for this new Krylov vector.
do idet = 1, int(TotWalkers)
int_sign = CurrentDets(IlutBits%ind_pop:IlutBits%ind_pop + lenof_sign_kp - 1, idet)
call extract_sign(CurrentDets(:, idet), real_sign)
tCoreDet = check_determ_flag(CurrentDets(:, idet))
! Don't add unoccpied determinants, unless they are core determinants.
if (IsUnoccDet(real_sign) .and. (.not. tCoreDet)) cycle
tDetFound = .false.
call decode_bit_det(nI, CurrentDets(:, idet))
! Search the hash table for this determinant.
call hash_table_lookup(nI, CurrentDets(:, idet), nifd, krylov_vecs_ht, krylov_vecs, det_ind, hash_val, tDetFound)
if (tDetFound) then
! In this case the determinant is already in the Krylov vector
! array.
krylov_vecs(sign_ind:sign_ind + lenof_sign_kp - 1, det_ind) = int_sign
else
! A new determiant needs to be added.
TotWalkersKP = TotWalkersKP + 1
det_ind = TotWalkersKP
if (TotWalkersKP > krylov_vecs_length) then
call stop_all(t_r, "There are no slots left in the krylov_vecs array for the next determinant. &
&You can increase the size of this array using the memory-factor option in &
&the kp-fciqmc block of the input file.")
end if
! Add the determinant's index to the hash table.
call add_hash_table_entry(krylov_vecs_ht, det_ind, hash_val)
! Copy determinant data across.
krylov_vecs(0:nifd, det_ind) = CurrentDets(0:nifd, idet)
krylov_vecs(sign_ind:sign_ind + lenof_sign_kp - 1, det_ind) = int_sign
krylov_helems(det_ind) = det_diagH(idet)
krylov_vecs(flag_ind, det_ind) = CurrentDets(IlutBits%ind_flag, idet)
end if
! Update information about how much of the hash table is filled.
do iamp = 1, lenof_sign_kp
if (.not. near_zero(real_sign(iamp))) then
nkrylov_amp_elems_used = nkrylov_amp_elems_used + 1
end if
end do
end do
if (tPrintDataTables) then
write(stdout, '(1x,a5)', advance='yes') "Done."
write(stdout, '(a56,'//int_fmt(TotWalkersKP, 1)//',1x,a17,'//int_fmt(krylov_vecs_length, 1)//')') &
"# Number unique determinants in the Krylov vector array:", TotWalkersKP, "out of a possible", krylov_vecs_length
amp_fraction = real(nkrylov_amp_elems_used, dp) / real(nkrylov_amp_elems_tot, dp)
write(stdout, '(a69,1x,es10.4)') "# Fraction of the amplitude elements used in the Krylov vector array:", amp_fraction
call neci_flush(stdout)
end if
end subroutine store_krylov_vec