subroutine generate_using_mp1_criterion(target_ndets, ilut_list, space_size)
! In: target_ndets - The number of determinants to keep, unless there are less
! singles and doubles than this value, in which case all singles and doubles
! will be kept.
! In/Out: ilut_list - List of determinants generated.
! In/Out: space_size - Number of determinants in the generated space.
! If ilut_list is not empty on input and you want to keep
! the states already in it, then on input space_size should
! be equal to the number of states to be kept in ilut_list,
! and new states will be added in from space_size+1.
! Otherwise, space_size must equal 0 on input.
! On output space_size will equal the total number of
! generated plus what space_size was on input.
integer, intent(in) :: target_ndets
integer(n_int), intent(inout) :: ilut_list(0:, :)
integer, intent(inout) :: space_size
integer(n_int), allocatable :: temp_list(:, :)
real(dp), allocatable :: amp_list(:)
integer(n_int) :: ilut(0:NIfTot)
integer :: nI(nel)
integer :: ex(2, maxExcit), ex_flag, ndets
integer :: pos, i
real(dp) :: amp, energy_contrib
logical :: tAllExcitFound, tParity
integer(n_int), allocatable :: excitations(:, :)
integer(n_int) :: ilutG(0:nifguga)
allocate(amp_list(target_ndets))
allocate(temp_list(0:NIfD, target_ndets))
amp_list = 0.0_dp
temp_list = 0_n_int
! Should we generate just singles (1), just doubles (2), or both (3)?
if (tUEG .or. tNoSingExcits) then
ex_flag = 2
else if (tHub) then
if (tReal) then
ex_flag = 1
else
ex_flag = 2
end if
else
! Generate both the single and double excitations.
ex_flag = 3
end if
tAllExcitFound = .false.
! Count the HF determinant.
ndets = 1
ex = 0
ilut = 0_n_int
! Start by adding the HF state.
temp_list(0:NIfD, 1) = ilutHF(0:NIfD)
amp_list = huge(amp)
! Set this amplitude to be the lowest possible number so that it doesn't get removed from
! the list in the following selection - we always want to keep the HF determinant.
amp_list(1) = -huge(amp)
! Loop through all connections to the HF determinant and keep the required number which
! have the largest MP1 weights.
if (tGUGA) then
! in guga, create all excitations at once and then check for the
! MP1 amplitude in an additional loop
call convert_ilut_toGUGA(ilutHF, ilutG)
call actHamiltonian(ilutG, CSF_Info_t(ilutG), excitations, ndets)
do i = 1, ndets
call convert_ilut_toNECI(excitations(:, i), ilut)
call decode_bit_det(nI, ilut)
call return_mp1_amp_and_mp2_energy(nI, ilut, ex, tParity, amp, &
energy_contrib)
pos = binary_search_real(amp_list, -abs(amp), 1.0e-8_dp)
if (pos < 0) pos = -pos
if (pos > 0 .and. pos <= target_ndets) then
! Shuffle all less significant determinants down one slot, and throw away the
! previous least significant determinant.
temp_list(0:NIfD, pos + 1:target_ndets) = temp_list(0:NIfD, pos:target_ndets - 1)
amp_list(pos + 1:target_ndets) = amp_list(pos:target_ndets - 1)
! Add in the new ilut and amplitude in the correct position.
temp_list(0:NIfD, pos) = ilut(0:NIfD)
! Store the negative absolute value, because the binary search sorts from
! lowest (most negative) to highest.
amp_list(pos) = -abs(amp)
end if
end do
else
do while (.true.)
call GenExcitations3(HFDet, ilutHF, nI, ex_flag, ex, tParity, tAllExcitFound, .false.)
! When no more basis functions are found, this value is returned and the loop is exited.
if (tAllExcitFound) exit
call EncodeBitDet(nI, ilut)
if (tHPHF) then
if (.not. IsAllowedHPHF(ilut(0:NIfD))) cycle
end if
ndets = ndets + 1
! If a determinant is returned (if we did not find the final one last time.)
if (.not. tAllExcitFound) then
call return_mp1_amp_and_mp2_energy(nI, ilut, ex, tParity, amp, energy_contrib)
pos = binary_search_real(amp_list, -abs(amp), 1.0e-8_dp)
! If pos is less then there isn't another determinant with the same amplitude
! (which will be common), but -pos specifies where in the list it should be
! inserted to keep amp_list in order.
if (pos < 0) pos = -pos
if (pos > 0 .and. pos <= target_ndets) then
! Shuffle all less significant determinants down one slot, and throw away the
! previous least significant determinant.
temp_list(0:NIfD, pos + 1:target_ndets) = temp_list(0:NIfD, pos:target_ndets - 1)
amp_list(pos + 1:target_ndets) = amp_list(pos:target_ndets - 1)
! Add in the new ilut and amplitude in the correct position.
temp_list(0:NIfD, pos) = ilut(0:NIfD)
! Store the negative absolute value, because the binary search sorts from
! lowest (most negative) to highest.
amp_list(pos) = -abs(amp)
end if
end if
end do
end if
call warning_neci("generate_using_mp1_criterion", &
"Note that there are less connections to the Hartree-Fock than the requested &
&size of the space. The space will therefore be smaller than requested, &
&containing all connections.")
! Now that the correct determinants have been selected, add them to the desired space.
do i = 1, min(ndets, target_ndets)
ilut(0:NIfD) = temp_list(0:NIfD, i)
call add_state_to_space(ilut, ilut_list, space_size)
end do
end subroutine generate_using_mp1_criterion