subroutine Add_RDM_From_IJ_Pair_GUGA(spawn, one_rdms, nI, nJ, sign_i, &
sign_j, t_bra_to_ket, t_fast, rdm_ind_in, x0, x1)
! corresponding GUGA routine from function found in rdm_filling.F90
type(rdm_spawn_t), intent(inout) :: spawn
type(one_rdm_t), intent(inout) :: one_rdms(:)
integer, intent(in) :: nI(nel), nJ(nel)
real(dp), intent(in) :: sign_i(:), sign_j(:)
logical, intent(in) :: t_bra_to_ket
logical, intent(in), optional :: t_fast
integer(int_rdm), intent(in), optional :: rdm_ind_in
real(dp), intent(in), optional :: x0, x1
#ifdef DEBUG_
character(*), parameter :: this_routine = "Add_RDM_From_IJ_Pair_GUGA"
#endif
type(ExcitationInformation_t) :: excitInfo
HElement_t(dp) :: mat_ele
integer(int_rdm), allocatable :: rdm_ind(:)
real(dp), allocatable :: rdm_mat(:)
integer :: i, j, k, l, n, ex_lvl, ex_typ
real(dp) :: full_sign(spawn%rdm_send%sign_length)
integer(n_int) :: ilutGi(0:GugaBits%len_tot), ilutGj(0:GugaBits%len_tot)
integer(int_rdm) :: pure_ind
logical :: t_fast_
def_default(t_fast_, t_fast, .true.)
if (t_fast_) then
ASSERT(present(rdm_ind_in))
ASSERT(present(x0))
ASSERT(present(x1))
ASSERT(rdm_ind_in /= 0_int_rdm)
ex_lvl = extract_excit_lvl_rdm(rdm_ind_in)
ex_typ = extract_excit_type_rdm(rdm_ind_in)
ASSERT(ex_lvl == 1 .or. ex_lvl == 2)
pure_ind = pure_rdm_ind(rdm_ind_in)
if (ex_lvl == 1) then
if (RDMExcitLevel == 1) then
call fill_sings_1rdm_guga(one_rdms, sign_I, sign_J, &
x0, pure_ind)
else
if (t_bra_to_ket) then
call EncodeBitDet_guga(nI, ilutGi)
call EncodeBitDet_guga(nJ, ilutGj)
else
call EncodeBitDet_guga(nI, ilutGJ)
call EncodeBitDet_guga(nJ, ilutGi)
end if
call fill_sings_2rdm_guga(&
spawn, IlutGi, ilutGj, sign_i, sign_j, x0, pure_ind)
end if
else if (ex_lvl == 2 .and. RDMExcitLevel /= 1) then
if (t_bra_to_ket) then
call EncodeBitDet_guga(nI, ilutGi)
call EncodeBitDet_guga(nJ, ilutGj)
else
call EncodeBitDet_guga(nI, ilutGJ)
call EncodeBitDet_guga(nJ, ilutGi)
end if
call create_all_rdm_contribs(rdm_ind_in, x0, x1, rdm_ind, rdm_mat)
do n = 1, size(rdm_ind)
if (.not. near_zero(rdm_mat(n))) then
call extract_2_rdm_ind(rdm_ind(n), i, j, k, l)
full_sign = sign_I * sign_J * rdm_mat(n)
if (.not. &
(ex_typ == excit_type%fullstart_stop_alike)) then
full_sign = 2.0_dp * full_sign
end if
call add_to_rdm_spawn_t(spawn, i, j, k, l, &
full_sign, .true.)
end if
end do
end if
else
if (t_bra_to_ket) then
call EncodeBitDet_guga(nI, ilutGi)
call EncodeBitDet_guga(nJ, ilutGj)
else
call EncodeBitDet_guga(nI, ilutGJ)
call EncodeBitDet_guga(nJ, ilutGi)
end if
block
type(CSF_Info_t) :: csf_i, csf_j
csf_i = CSF_Info_t(ilutGi)
csf_j = CSF_Info_t(ilutGj)
call calc_guga_matrix_element(IlutGi, csf_i, ilutGj, csf_j, excitInfo, mat_ele, &
t_hamil=.false., rdm_ind=rdm_ind, &
rdm_mat=rdm_mat)
end block
! i assume sign_i and sign_j are not 0 if we end up here..
do n = 1, size(rdm_ind)
if (.not. near_zero(rdm_mat(n))) then
if (excitInfo%excitLvl == 1) then
if (RDMExcitLevel == 1) then
call fill_sings_1rdm_guga(one_rdms, sign_I, sign_J, &
rdm_mat(n), rdm_ind(n))
else
call fill_sings_2rdm_guga(spawn, IlutGi, &
ilutGj, sign_i, sign_j, rdm_mat(n), rdm_ind(n))
end if
else if (excitInfo%excitLvl == 2 .and. RDMExcitLevel /= 1) then
call extract_2_rdm_ind(rdm_ind(n), i, j, k, l)
full_sign = sign_I * sign_J * rdm_mat(n)
! here in the 'exact' filling (coming from HF or
! within the semistochastic space I think it makes
! sense to fill symmetrically.. since here no
! stochastic spawning is happening and this does not
! give us information about the hermiticity error!
call add_to_rdm_spawn_t(spawn, i, j, k, l, &
full_sign, .true.)
if (.not. &
(excitInfo%typ == excit_type%fullstart_stop_alike)) then
call add_to_rdm_spawn_t(spawn, k, l, i, j, &
full_sign, .true.)
end if
end if
end if
end do
end if
end subroutine Add_RDM_From_IJ_Pair_GUGA