From 0e31cfee7fa4eb5aaf8df9d751dfcd2a79ab978a Mon Sep 17 00:00:00 2001 From: Kevin Gasperich Date: Sun, 23 Feb 2020 16:40:26 -0600 Subject: [PATCH] complex slater_rules_wee_mono --- src/determinants/slater_rules_wee_mono.irp.f | 192 +++++++++++++++++++ src/utils_complex/qp2-pbc-diff.txt | 6 +- 2 files changed, 195 insertions(+), 3 deletions(-) diff --git a/src/determinants/slater_rules_wee_mono.irp.f b/src/determinants/slater_rules_wee_mono.irp.f index 3a8c9075..92754104 100644 --- a/src/determinants/slater_rules_wee_mono.irp.f +++ b/src/determinants/slater_rules_wee_mono.irp.f @@ -361,3 +361,195 @@ subroutine i_H_j_two_e(key_i,key_j,Nint,hij) end select end +!==============================================================================! +! ! +! Complex ! +! ! +!==============================================================================! + +subroutine i_Wee_j_single_complex(key_i,key_j,Nint,spin,hij) + use bitmasks + implicit none + BEGIN_DOC + ! Returns $\langle i|H|j \rangle$ where $i$ and $j$ are determinants differing by a + ! single excitation. + END_DOC + integer, intent(in) :: Nint, spin + integer(bit_kind), intent(in) :: key_i(Nint,2), key_j(Nint,2) + complex*16, intent(out) :: hij + + integer :: exc(0:2,2) + double precision :: phase + + PROVIDE big_array_exchange_integrals_complex mo_two_e_integrals_in_map + + call get_single_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint) + call single_excitation_wee_complex(key_i,key_j,exc(1,1),exc(1,2),spin,phase,hij) +end + + +subroutine i_H_j_mono_spin_one_e_complex(key_i,key_j,Nint,spin,hij) + !todo: check hole/particle m/p ordering? + if (is_complex) then + print*,irp_here,' not implemented for complex' + stop -1 + endif + use bitmasks + implicit none + BEGIN_DOC + ! Returns $\langle i|H|j \rangle$ where $i$ and $j$ are determinants differing by + ! a single excitation. + END_DOC + integer, intent(in) :: Nint, spin + integer(bit_kind), intent(in) :: key_i(Nint,2), key_j(Nint,2) + complex*16, intent(out) :: hij + + integer :: exc(0:2,2) + double precision :: phase + + call get_single_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint) + integer :: m,p + m = exc(1,1) + p = exc(1,2) + hij = phase * mo_one_e_integrals_complex(m,p) +end + +subroutine i_H_j_one_e_complex(key_i,key_j,Nint,hij) + !todo: check hole/particle m/p ordering? + if (is_complex) then + print*,irp_here,' not implemented for complex' + stop -1 + endif + use bitmasks + implicit none + BEGIN_DOC + ! Returns $\langle i|H|j \rangle$ where $i$ and $j$ are determinants. + END_DOC + integer, intent(in) :: Nint + integer(bit_kind), intent(in) :: key_i(Nint,2), key_j(Nint,2) + complex*16, intent(out) :: hij + + integer :: degree,m,p + double precision :: diag_h_mat_elem_one_e,phase + integer :: exc(0:2,2,2) + call get_excitation_degree(key_i,key_j,degree,Nint) + hij = (0.d0,0.d0) + if(degree>1)then + return + endif + if(degree==0)then + hij = dcmplx(diag_H_mat_elem_one_e(key_i,N_int),0.d0) + else + call get_single_excitation(key_i,key_j,exc,phase,Nint) + if (exc(0,1,1) == 1) then + ! Mono alpha + m = exc(1,1,1) + p = exc(1,2,1) + else + ! Mono beta + m = exc(1,1,2) + p = exc(1,2,2) + endif + hij = phase * mo_one_e_integrals_complex(m,p) + endif + +end + +subroutine i_H_j_two_e_complex(key_i,key_j,Nint,hij) + !todo: check hole/particle m/p ordering? + if (is_complex) then + print*,irp_here,' not implemented for complex' + stop -1 + endif + use bitmasks + implicit none + BEGIN_DOC + ! Returns $\langle i|H|j \rangle$ where $i$ and $j$ are determinants. + END_DOC + integer, intent(in) :: Nint + integer(bit_kind), intent(in) :: key_i(Nint,2), key_j(Nint,2) + complex*16, intent(out) :: hij + + integer :: exc(0:2,2,2) + integer :: degree + complex*16 :: get_two_e_integral_complex + integer :: m,n,p,q + integer :: i,j,k + integer :: occ(Nint*bit_kind_size,2) + double precision :: diag_H_mat_elem, phase,phase_2 + integer :: n_occ_ab(2) + PROVIDE mo_two_e_integrals_in_map mo_integrals_map big_array_exchange_integrals_complex ref_bitmask_two_e_energy + + ASSERT (Nint > 0) + ASSERT (Nint == N_int) + ASSERT (sum(popcnt(key_i(:,1))) == elec_alpha_num) + ASSERT (sum(popcnt(key_i(:,2))) == elec_beta_num) + ASSERT (sum(popcnt(key_j(:,1))) == elec_alpha_num) + ASSERT (sum(popcnt(key_j(:,2))) == elec_beta_num) + + hij = (0.d0,0.d0) + !DIR$ FORCEINLINE + call get_excitation_degree(key_i,key_j,degree,Nint) + integer :: spin + select case (degree) + case (2) + call get_double_excitation(key_i,key_j,exc,phase,Nint) + if (exc(0,1,1) == 1) then + ! Mono alpha, mono beta + if(exc(1,1,1) == exc(1,2,2) )then + hij = phase * big_array_exchange_integrals_complex(exc(1,1,1),exc(1,1,2),exc(1,2,1)) + else if (exc(1,2,1) ==exc(1,1,2))then + hij = phase * big_array_exchange_integrals_complex(exc(1,2,1),exc(1,1,1),exc(1,2,2)) + else + hij = phase*get_two_e_integral_complex( & + exc(1,1,1), & + exc(1,1,2), & + exc(1,2,1), & + exc(1,2,2) ,mo_integrals_map,mo_integrals_map_2) + endif + else if (exc(0,1,1) == 2) then + ! Double alpha + hij = phase*(get_two_e_integral_complex( & + exc(1,1,1), & + exc(2,1,1), & + exc(1,2,1), & + exc(2,2,1) ,mo_integrals_map,mo_integrals_map_2) - & + get_two_e_integral_complex( & + exc(1,1,1), & + exc(2,1,1), & + exc(2,2,1), & + exc(1,2,1) ,mo_integrals_map,mo_integrals_map_2) ) + else if (exc(0,1,2) == 2) then + ! Double beta + hij = phase*(get_two_e_integral_complex( & + exc(1,1,2), & + exc(2,1,2), & + exc(1,2,2), & + exc(2,2,2) ,mo_integrals_map,mo_integrals_map_2) - & + get_two_e_integral_complex( & + exc(1,1,2), & + exc(2,1,2), & + exc(2,2,2), & + exc(1,2,2) ,mo_integrals_map,mo_integrals_map_2) ) + endif + case (1) + call get_single_excitation(key_i,key_j,exc,phase,Nint) + !DIR$ FORCEINLINE + call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint) + if (exc(0,1,1) == 1) then + ! Mono alpha + m = exc(1,1,1) + p = exc(1,2,1) + spin = 1 + else + ! Mono beta + m = exc(1,1,2) + p = exc(1,2,2) + spin = 2 + endif + call single_excitation_wee_complex(key_i,key_j,p,m,spin,phase,hij) + case (0) + double precision :: diag_wee_mat_elem + hij = dcmplx(diag_wee_mat_elem(key_i,Nint),0.d0) + end select +end diff --git a/src/utils_complex/qp2-pbc-diff.txt b/src/utils_complex/qp2-pbc-diff.txt index 46667e68..f2b3ed38 100644 --- a/src/utils_complex/qp2-pbc-diff.txt +++ b/src/utils_complex/qp2-pbc-diff.txt @@ -8,7 +8,7 @@ general: determinants: (done) connected_to_ref.irp.f (done) create_excitations.irp.f - (done?)density_matrix{,_complex}.irp.f + (done?) density_matrix{,_complex}.irp.f no one_e_dm_dagger_mo_spin_index_complex need to test for complex (done) determinants_bitmasks.irp.f @@ -45,11 +45,11 @@ determinants: made copies of needed functions for complex still need to do implementation (done) single_excitations.irp.f - (done?)single_excitation_two_e.irp.f + (done?) single_excitation_two_e.irp.f (****) slater_rules.irp.f made copies of needed functions for complex still need to do implementation - (****) slater_rules_wee_mono.irp.f + (done?) slater_rules_wee_mono.irp.f (done) sort_dets_ab.irp.f spindeterminants.ezfio_config need svd complex?