mirror of
https://github.com/QuantumPackage/qp2.git
synced 2024-09-01 13:43:40 +02:00
added optimization for Slater_tc in two-e elements
This commit is contained in:
eginer
parent
ec05b8c329
commit
2e45413f44
@ -199,3 +199,29 @@ END_PROVIDER
|
||||
|
||||
! ---
|
||||
|
||||
|
||||
BEGIN_PROVIDER [ double precision, mo_bi_ortho_tc_two_e_jj, (mo_num,mo_num) ]
|
||||
&BEGIN_PROVIDER [ double precision, mo_bi_ortho_tc_two_e_jj_exchange, (mo_num,mo_num) ]
|
||||
&BEGIN_PROVIDER [ double precision, mo_bi_ortho_tc_two_e_jj_anti, (mo_num,mo_num) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! mo_bi_ortho_tc_two_e_jj(i,j) = J_ij = <ji|W-K|ji>
|
||||
! mo_bi_ortho_tc_two_e_jj_exchange(i,j) = K_ij = <ij|W-K|ji>
|
||||
! mo_bi_ortho_tc_two_e_jj_anti(i,j) = J_ij - K_ij
|
||||
END_DOC
|
||||
|
||||
integer :: i,j
|
||||
double precision :: get_two_e_integral
|
||||
|
||||
mo_bi_ortho_tc_two_e_jj = 0.d0
|
||||
mo_bi_ortho_tc_two_e_jj_exchange = 0.d0
|
||||
|
||||
do i=1,mo_num
|
||||
do j=1,mo_num
|
||||
mo_bi_ortho_tc_two_e_jj(i,j) = mo_bi_ortho_tc_two_e(j,i,j,i)
|
||||
mo_bi_ortho_tc_two_e_jj_exchange(i,j) = mo_bi_ortho_tc_two_e(i,j,j,i)
|
||||
mo_bi_ortho_tc_two_e_jj_anti(i,j) = mo_bi_ortho_tc_two_e_jj(i,j) - mo_bi_ortho_tc_two_e_jj_exchange(i,j)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
@ -1790,12 +1790,12 @@ double precision function diag_H_mat_elem(det_in,Nint)
|
||||
integer :: tmp(2)
|
||||
!DIR$ FORCEINLINE
|
||||
call bitstring_to_list_ab(particle, occ_particle, tmp, Nint)
|
||||
ASSERT (tmp(1) == nexc(1))
|
||||
ASSERT (tmp(2) == nexc(2))
|
||||
ASSERT (tmp(1) == nexc(1)) ! Number of particles alpha
|
||||
ASSERT (tmp(2) == nexc(2)) ! Number of particle beta
|
||||
!DIR$ FORCEINLINE
|
||||
call bitstring_to_list_ab(hole, occ_hole, tmp, Nint)
|
||||
ASSERT (tmp(1) == nexc(1))
|
||||
ASSERT (tmp(2) == nexc(2))
|
||||
ASSERT (tmp(1) == nexc(1)) ! Number of holes alpha
|
||||
ASSERT (tmp(2) == nexc(2)) ! Number of holes beta
|
||||
|
||||
det_tmp = ref_bitmask
|
||||
do ispin=1,2
|
||||
|
||||
208
src/tc_bi_ortho/slater_tc_opt.irp.f
Normal file
208
src/tc_bi_ortho/slater_tc_opt.irp.f
Normal file
@ -0,0 +1,208 @@
|
||||
BEGIN_PROVIDER [ double precision, ref_tc_energy_tot]
|
||||
&BEGIN_PROVIDER [ double precision, ref_tc_energy_1e]
|
||||
&BEGIN_PROVIDER [ double precision, ref_tc_energy_2e]
|
||||
&BEGIN_PROVIDER [ double precision, ref_tc_energy_3e]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Various component of the TC energy for the reference "HF" Slater determinant
|
||||
END_DOC
|
||||
double precision :: hmono, htwoe, htot, hthree
|
||||
call diag_htilde_mu_mat_bi_ortho(N_int,HF_bitmask , hmono, htwoe, htot)
|
||||
ref_tc_energy_1e = hmono
|
||||
ref_tc_energy_2e = htwoe
|
||||
if(three_body_h_tc)then
|
||||
call diag_htilde_three_body_ints_bi_ort(N_int, HF_bitmask, hthree)
|
||||
ref_tc_energy_3e = hthree
|
||||
else
|
||||
ref_tc_energy_3e = 0.d0
|
||||
endif
|
||||
ref_tc_energy_tot = ref_tc_energy_1e + ref_tc_energy_2e + ref_tc_energy_3e
|
||||
END_PROVIDER
|
||||
|
||||
subroutine diag_htilde_mu_mat_fock_bi_ortho(Nint, det_in, hmono, htwoe, hthree, htot)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Computes $\langle i|H|i \rangle$.
|
||||
END_DOC
|
||||
integer,intent(in) :: Nint
|
||||
integer(bit_kind),intent(in) :: det_in(Nint,2)
|
||||
double precision, intent(out) :: hmono,htwoe,htot,hthree
|
||||
|
||||
integer(bit_kind) :: hole(Nint,2)
|
||||
integer(bit_kind) :: particle(Nint,2)
|
||||
integer :: i, nexc(2), ispin
|
||||
integer :: occ_particle(Nint*bit_kind_size,2)
|
||||
integer :: occ_hole(Nint*bit_kind_size,2)
|
||||
integer(bit_kind) :: det_tmp(Nint,2)
|
||||
integer :: na, nb
|
||||
|
||||
ASSERT (Nint > 0)
|
||||
ASSERT (sum(popcnt(det_in(:,1))) == elec_alpha_num)
|
||||
ASSERT (sum(popcnt(det_in(:,2))) == elec_beta_num)
|
||||
|
||||
|
||||
nexc(1) = 0
|
||||
nexc(2) = 0
|
||||
do i=1,Nint
|
||||
hole(i,1) = xor(det_in(i,1),ref_bitmask(i,1))
|
||||
hole(i,2) = xor(det_in(i,2),ref_bitmask(i,2))
|
||||
particle(i,1) = iand(hole(i,1),det_in(i,1))
|
||||
particle(i,2) = iand(hole(i,2),det_in(i,2))
|
||||
hole(i,1) = iand(hole(i,1),ref_bitmask(i,1))
|
||||
hole(i,2) = iand(hole(i,2),ref_bitmask(i,2))
|
||||
nexc(1) = nexc(1) + popcnt(hole(i,1))
|
||||
nexc(2) = nexc(2) + popcnt(hole(i,2))
|
||||
enddo
|
||||
|
||||
if (nexc(1)+nexc(2) == 0) then
|
||||
htot = ref_tc_energy_tot
|
||||
return
|
||||
endif
|
||||
|
||||
!call debug_det(det_in,Nint)
|
||||
integer :: tmp(2)
|
||||
!DIR$ FORCEINLINE
|
||||
call bitstring_to_list_ab(particle, occ_particle, tmp, Nint)
|
||||
ASSERT (tmp(1) == nexc(1)) ! Number of particles alpha
|
||||
ASSERT (tmp(2) == nexc(2)) ! Number of particle beta
|
||||
!DIR$ FORCEINLINE
|
||||
call bitstring_to_list_ab(hole, occ_hole, tmp, Nint)
|
||||
ASSERT (tmp(1) == nexc(1)) ! Number of holes alpha
|
||||
ASSERT (tmp(2) == nexc(2)) ! Number of holes beta
|
||||
|
||||
|
||||
det_tmp = ref_bitmask
|
||||
hmono = ref_tc_energy_1e
|
||||
htwoe = ref_tc_energy_2e
|
||||
hthree= ref_tc_energy_3e
|
||||
|
||||
do ispin=1,2
|
||||
na = elec_num_tab(ispin)
|
||||
nb = elec_num_tab(iand(ispin,1)+1)
|
||||
do i=1,nexc(ispin)
|
||||
!DIR$ FORCEINLINE
|
||||
call ac_tc_operator( occ_particle(i,ispin), ispin, det_tmp, hmono,htwoe,hthree, Nint,na,nb)
|
||||
!DIR$ FORCEINLINE
|
||||
call a_tc_operator ( occ_hole (i,ispin), ispin, det_tmp, hmono,htwoe,hthree, Nint,na,nb)
|
||||
enddo
|
||||
enddo
|
||||
htot = hmono+htwoe+hthree
|
||||
end
|
||||
|
||||
subroutine ac_tc_operator(iorb,ispin,key,hmono,htwoe,hthree,Nint,na,nb)
|
||||
use bitmasks
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Routine that computes one- and two-body energy corresponding
|
||||
!
|
||||
! to the ADDITION of an electron in an orbital 'iorb' of spin 'ispin'
|
||||
!
|
||||
! onto a determinant 'key'.
|
||||
!
|
||||
! in output, the determinant key is changed by the ADDITION of that electron
|
||||
!
|
||||
! and the quantities hmono,htwoe,hthree are INCREMENTED
|
||||
END_DOC
|
||||
integer, intent(in) :: iorb, ispin, Nint
|
||||
integer, intent(inout) :: na, nb
|
||||
integer(bit_kind), intent(inout) :: key(Nint,2)
|
||||
double precision, intent(inout) :: hmono,htwoe,hthree
|
||||
|
||||
integer :: occ(Nint*bit_kind_size,2)
|
||||
integer :: other_spin
|
||||
integer :: k,l,i
|
||||
|
||||
if (iorb < 1) then
|
||||
print *, irp_here, ': iorb < 1'
|
||||
print *, iorb, mo_num
|
||||
stop -1
|
||||
endif
|
||||
if (iorb > mo_num) then
|
||||
print *, irp_here, ': iorb > mo_num'
|
||||
print *, iorb, mo_num
|
||||
stop -1
|
||||
endif
|
||||
|
||||
ASSERT (ispin > 0)
|
||||
ASSERT (ispin < 3)
|
||||
ASSERT (Nint > 0)
|
||||
|
||||
integer :: tmp(2)
|
||||
!DIR$ FORCEINLINE
|
||||
call bitstring_to_list_ab(key, occ, tmp, Nint)
|
||||
ASSERT (tmp(1) == elec_alpha_num)
|
||||
ASSERT (tmp(2) == elec_beta_num)
|
||||
|
||||
k = shiftr(iorb-1,bit_kind_shift)+1
|
||||
ASSERT (k >0)
|
||||
l = iorb - shiftl(k-1,bit_kind_shift)-1
|
||||
ASSERT (l >= 0)
|
||||
key(k,ispin) = ibset(key(k,ispin),l)
|
||||
other_spin = iand(ispin,1)+1
|
||||
|
||||
hmono = hmono + mo_bi_ortho_tc_one_e(iorb,iorb)
|
||||
|
||||
! Same spin
|
||||
do i=1,na
|
||||
htwoe = htwoe + mo_bi_ortho_tc_two_e_jj_anti(occ(i,ispin),iorb)
|
||||
enddo
|
||||
|
||||
! Opposite spin
|
||||
do i=1,nb
|
||||
htwoe = htwoe + mo_bi_ortho_tc_two_e_jj(occ(i,other_spin),iorb)
|
||||
enddo
|
||||
na = na+1
|
||||
end
|
||||
|
||||
subroutine a_tc_operator(iorb,ispin,key,hmono,htwoe,hthree,Nint,na,nb)
|
||||
use bitmasks
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Routine that computes one- and two-body energy corresponding
|
||||
!
|
||||
! to the REMOVAL of an electron in an orbital 'iorb' of spin 'ispin'
|
||||
!
|
||||
! onto a determinant 'key'.
|
||||
!
|
||||
! in output, the determinant key is changed by the REMOVAL of that electron
|
||||
!
|
||||
! and the quantities hmono,htwoe,hthree are INCREMENTED
|
||||
END_DOC
|
||||
integer, intent(in) :: iorb, ispin, Nint
|
||||
integer, intent(inout) :: na, nb
|
||||
integer(bit_kind), intent(inout) :: key(Nint,2)
|
||||
double precision, intent(inout) :: hmono,htwoe,hthree
|
||||
|
||||
integer :: occ(Nint*bit_kind_size,2)
|
||||
integer :: other_spin
|
||||
integer :: k,l,i
|
||||
integer :: tmp(2)
|
||||
|
||||
ASSERT (iorb > 0)
|
||||
ASSERT (ispin > 0)
|
||||
ASSERT (ispin < 3)
|
||||
ASSERT (Nint > 0)
|
||||
|
||||
k = shiftr(iorb-1,bit_kind_shift)+1
|
||||
ASSERT (k>0)
|
||||
l = iorb - shiftl(k-1,bit_kind_shift)-1
|
||||
key(k,ispin) = ibclr(key(k,ispin),l)
|
||||
other_spin = iand(ispin,1)+1
|
||||
|
||||
!DIR$ FORCEINLINE
|
||||
call bitstring_to_list_ab(key, occ, tmp, Nint)
|
||||
na = na-1
|
||||
|
||||
hmono = hmono - mo_bi_ortho_tc_one_e(iorb,iorb)
|
||||
|
||||
! Same spin
|
||||
do i=1,na
|
||||
htwoe= htwoe- mo_bi_ortho_tc_two_e_jj_anti(occ(i,ispin),iorb)
|
||||
enddo
|
||||
|
||||
! Opposite spin
|
||||
do i=1,nb
|
||||
htwoe= htwoe- mo_bi_ortho_tc_two_e_jj(occ(i,other_spin),iorb)
|
||||
enddo
|
||||
|
||||
end
|
||||
@ -11,121 +11,18 @@ program tc_bi_ortho
|
||||
touch read_wf
|
||||
touch my_grid_becke my_n_pt_r_grid my_n_pt_a_grid
|
||||
|
||||
! call routine_2
|
||||
call test_rout
|
||||
call test_slater_tc_opt
|
||||
end
|
||||
|
||||
subroutine test_rout
|
||||
subroutine test_slater_tc_opt
|
||||
implicit none
|
||||
integer :: i,j,ii,jj
|
||||
use bitmasks ! you need to include the bitmasks_module.f90 features
|
||||
integer(bit_kind), allocatable :: det_i(:,:)
|
||||
allocate(det_i(N_int,2))
|
||||
det_i(:,:)= psi_det(:,:,1)
|
||||
call debug_det(det_i,N_int)
|
||||
integer, allocatable :: occ(:,:)
|
||||
integer :: n_occ_ab(2)
|
||||
allocate(occ(N_int*bit_kind_size,2))
|
||||
call bitstring_to_list_ab(det_i, occ, n_occ_ab, N_int)
|
||||
double precision :: hmono, htwoe, htot
|
||||
call diag_htilde_mu_mat_bi_ortho(N_int, det_i, hmono, htwoe, htot)
|
||||
print*,'hmono, htwoe, htot'
|
||||
print*, hmono, htwoe, htot
|
||||
print*,'alpha electrons orbital occupancy'
|
||||
do i = 1, n_occ_ab(1) ! browsing the alpha electrons
|
||||
j = occ(i,1)
|
||||
print*,j,mo_bi_ortho_tc_one_e(j,j)
|
||||
enddo
|
||||
print*,'beta electrons orbital occupancy'
|
||||
do i = 1, n_occ_ab(2) ! browsing the beta electrons
|
||||
j = occ(i,2)
|
||||
print*,j,mo_bi_ortho_tc_one_e(j,j)
|
||||
enddo
|
||||
print*,'alpha beta'
|
||||
do i = 1, n_occ_ab(1)
|
||||
ii = occ(i,1)
|
||||
do j = 1, n_occ_ab(2)
|
||||
jj = occ(j,2)
|
||||
print*,ii,jj,mo_bi_ortho_tc_two_e(jj,ii,jj,ii)
|
||||
enddo
|
||||
enddo
|
||||
print*,'alpha alpha'
|
||||
do i = 1, n_occ_ab(1)
|
||||
ii = occ(i,1)
|
||||
do j = 1, n_occ_ab(1)
|
||||
jj = occ(j,1)
|
||||
print*,ii,jj,mo_bi_ortho_tc_two_e(jj,ii,jj,ii), mo_bi_ortho_tc_two_e(ii,jj,jj,ii)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
print*,'beta beta'
|
||||
do i = 1, n_occ_ab(2)
|
||||
ii = occ(i,2)
|
||||
do j = 1, n_occ_ab(2)
|
||||
jj = occ(j,2)
|
||||
print*,ii,jj,mo_bi_ortho_tc_two_e(jj,ii,jj,ii), mo_bi_ortho_tc_two_e(ii,jj,jj,ii)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
|
||||
end
|
||||
|
||||
subroutine routine_2
|
||||
implicit none
|
||||
integer :: i
|
||||
double precision :: bi_ortho_mo_ints
|
||||
print*,'H matrix'
|
||||
integer :: i,j
|
||||
double precision :: hmono, htwoe, htot, hthree
|
||||
double precision :: hnewmono, hnewtwoe, hnewthnewree, hnewtot
|
||||
do i = 1, N_det
|
||||
write(*,'(1000(F16.5,X))')htilde_matrix_elmt_bi_ortho(:,i)
|
||||
enddo
|
||||
i = 1
|
||||
double precision :: phase
|
||||
integer :: degree,h1, p1, h2, p2, s1, s2, exc(0:2,2,2)
|
||||
call get_excitation_degree(ref_bitmask, psi_det(1,1,i), degree, N_int)
|
||||
if(degree==2)then
|
||||
call get_double_excitation(ref_bitmask, psi_det(1,1,i), exc, phase, N_int)
|
||||
call decode_exc(exc, 2, h1, p1, h2, p2, s1, s2)
|
||||
print*,'h1,h2,p1,p2'
|
||||
print*, h1,h2,p1,p2
|
||||
print*,mo_bi_ortho_tc_two_e(p1,p2,h1,h2),mo_bi_ortho_tc_two_e(h1,h2,p1,p2)
|
||||
endif
|
||||
|
||||
|
||||
print*,'coef'
|
||||
do i = 1, ao_num
|
||||
print*,i,mo_l_coef(i,8),mo_r_coef(i,8)
|
||||
enddo
|
||||
! print*,'mdlqfmlqgmqglj'
|
||||
! print*,'mo_bi_ortho_tc_two_e()',mo_bi_ortho_tc_two_e(2,2,3,3)
|
||||
! print*,'bi_ortho_mo_ints ',bi_ortho_mo_ints(2,2,3,3)
|
||||
print*,'Overlap'
|
||||
do i = 1, mo_num
|
||||
write(*,'(100(F16.10,X))')overlap_bi_ortho(:,i)
|
||||
call diag_htilde_mu_mat_bi_ortho(N_int, psi_det(1,1,i), hmono, htwoe, htot)
|
||||
call diag_htilde_mu_mat_fock_bi_ortho(N_int, psi_det(1,1,i), hnewmono, hnewtwoe, hnewthnewree, hnewtot)
|
||||
print*,htot,hnewtot,dabs(htot-hnewtot)
|
||||
enddo
|
||||
|
||||
end
|
||||
|
||||
subroutine routine
|
||||
implicit none
|
||||
double precision :: hmono,htwoe,hthree,htot
|
||||
integer(bit_kind), allocatable :: key1(:,:)
|
||||
integer(bit_kind), allocatable :: key2(:,:)
|
||||
allocate(key1(N_int,2),key2(N_int,2))
|
||||
use bitmasks
|
||||
key1 = ref_bitmask
|
||||
call htilde_mu_mat_bi_ortho(key1,key1, N_int, hmono,htwoe,hthree,htot)
|
||||
key2 = key1
|
||||
integer :: h,p,i_ok
|
||||
h = 1
|
||||
p = 8
|
||||
call do_single_excitation(key2,h,p,1,i_ok)
|
||||
call debug_det(key2,N_int)
|
||||
call htilde_mu_mat_bi_ortho(key2,key1, N_int, hmono,htwoe,hthree,htot)
|
||||
! print*,'fock_matrix_tc_mo_alpha(p,h) = ',fock_matrix_tc_mo_alpha(p,h)
|
||||
print*,'htot = ',htot
|
||||
print*,'hmono = ',hmono
|
||||
print*,'htwoe = ',htwoe
|
||||
double precision :: bi_ortho_mo_ints
|
||||
print*,'bi_ortho_mo_ints(1,p,1,h)',bi_ortho_mo_ints(1,p,1,h)
|
||||
|
||||
end
|
||||
|
||||
Loading…
Reference in New Issue
Block a user