subroutine create_crude_guga_double(ilut, nI, csf_i, exc, pgen, excitInfo_in)
integer(n_int), intent(in) :: ilut(0:nifguga)
integer, intent(in) :: nI(nel)
type(CSF_Info_t), intent(in) :: csf_i
integer(n_int), intent(out) :: exc(0:nifguga)
real(dp), intent(out) :: pgen
type(ExcitationInformation_t), intent(in), optional :: excitInfo_in
type(ExcitationInformation_t) :: excitInfo
logical :: compFlag
real(dp) :: posSwitches(nSpatOrbs), negSwitches(nSpatOrbs)
real(dp) :: branch_pgen, orb_pgen
HElement_t(dp) :: mat_ele
integer(n_int) :: ilutI(0:niftot), ilutJ(0:niftot)
type(WeightObj_t) :: weights
integer :: elecs(2), orbs(2)
if (present(excitInfo_in)) then
excitInfo = excitInfo_in
else
call pickOrbitals_double(ilut, nI, csf_i, excitInfo, orb_pgen)
if (.not. excitInfo%valid) then
! if no valid indices were picked, return 0 excitation and return
exc = 0_n_int
pgen = 0.0_dp
return
end if
call checkCompatibility(&
csf_i, excitInfo, compFlag, posSwitches, negSwitches, weights)
if (.not. compFlag) then
exc = 0_n_int
pgen = 0.0_dp
return
end if
end if
! here I have to do the actual crude double excitation..
! my idea for now is to create pseudo random spin-orbital from the
! picked spatial orbitals
call create_random_spin_orbs(ilut, csf_i, excitInfo, elecs, orbs, branch_pgen)
if (any(elecs == 0) .or. any(orbs == 0)) then
exc = 0_n_int
pgen = 0.0_dp
return
end if
if (near_zero(branch_pgen)) then
exc = 0_n_int
pgen = 0.0_dp
return
end if
exc = ilut
clr_orb(exc, elecs(1))
clr_orb(exc, elecs(2))
set_orb(exc, orbs(1))
set_orb(exc, orbs(2))
! this check is the same at the end of singles:
! maybe I also need to do this only if no excitInfo_in is provided..
! we also want to check if we produced a valid CSF here
if (.not. isProperCSF_ilut(exc, .true.)) then
exc = 0_n_int
pgen = 0.0_dp
return
end if
! we also need to calculate the matrix element here!
call convert_ilut_toNECI(ilut, ilutI)
call convert_ilut_toNECI(exc, ilutJ)
call calc_guga_matrix_element(ilutI, csf_i, ilutJ, CSF_Info_t(ilutJ), excitInfo, mat_ele, .true.)
if (near_zero(mat_ele)) then
exc = 0_n_int
pgen = 0.0_dp
return
end if
call encode_matrix_element(exc, 0.0_dp, 2)
call encode_matrix_element(exc, mat_ele, 1)
global_excitInfo = excitInfo
if (present(excitInfo_in)) then
! then the orbitals were already picked before and we only want
! to give the branch_pgen here
pgen = branch_pgen
else
! otherwise we want the full pgen
pgen = orb_pgen * branch_pgen
end if
end subroutine create_crude_guga_double