LCPQ/qp2
9
1
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2025-01-03 09:05:39 +01:00
Code Releases Activity

two body part up to single excitations work with fock

Browse Source
This commit is contained in:
eginer 2023-01-19 19:29:26 +01:00
parent 1651242fba
commit 4ee0802150
5 changed files with 165 additions and 97 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

23
src/bi_ort_ints/total_twoe_pot.irp.f
View File

@ -224,4 +224,27 @@ END_PROVIDER
enddo enddo
enddo enddo
END_PROVIDER
BEGIN_PROVIDER [double precision, tc_2e_3idx_coulomb_integrals, (mo_num,mo_num, mo_num)]
&BEGIN_PROVIDER [double precision, tc_2e_3idx_exchange_integrals,(mo_num,mo_num, mo_num)]
implicit none
BEGIN_DOC
! tc_2e_3idx_coulomb_integrals(j,k,i) = <jk|ji>
!
! tc_2e_3idx_exchange_integrals(j,k,i) = <kj|ji>
END_DOC
integer :: i,j,k,l
double precision :: get_two_e_integral
double precision :: integral
do i = 1, mo_num
do k = 1, mo_num
do j = 1, mo_num
tc_2e_3idx_coulomb_integrals(j, k,i) = mo_bi_ortho_tc_two_e(j ,k ,j ,i )
tc_2e_3idx_exchange_integrals(j,k,i) = mo_bi_ortho_tc_two_e(k ,j ,j ,i )
enddo
enddo
enddo
END_PROVIDER END_PROVIDER

2
src/determinants/single_excitations.irp.f
View File

@ -28,7 +28,7 @@ BEGIN_PROVIDER [double precision, fock_operator_closed_shell_ref_bitmask, (mo_nu
integer :: occ_virt(N_int*bit_kind_size,2) integer :: occ_virt(N_int*bit_kind_size,2)
integer(bit_kind) :: key_test(N_int) integer(bit_kind) :: key_test(N_int)
integer(bit_kind) :: key_virt(N_int,2) integer(bit_kind) :: key_virt(N_int,2)
fock_operator_closed_shell_ref_bitmask = 0.d0
call bitstring_to_list_ab(ref_closed_shell_bitmask, occ, n_occ_ab, N_int) call bitstring_to_list_ab(ref_closed_shell_bitmask, occ, n_occ_ab, N_int)
do i = 1, N_int do i = 1, N_int
key_virt(i,1) = full_ijkl_bitmask(i) key_virt(i,1) = full_ijkl_bitmask(i)

3
src/tc_bi_ortho/slater_tc.irp.f
View File

@ -324,6 +324,9 @@ subroutine single_htilde_mu_mat_bi_ortho(Nint, key_j, key_i, hmono, htwoe, htot)
call get_single_excitation(key_i, key_j, exc, phase, Nint) call get_single_excitation(key_i, key_j, exc, phase, Nint)
call decode_exc(exc,1,h1,p1,h2,p2,s1,s2) call decode_exc(exc,1,h1,p1,h2,p2,s1,s2)
! if(h1==14.and.p1==2)then
! print*,'h1,p1 old = ',h1,p1
! endif
hmono = mo_bi_ortho_tc_one_e(p1,h1) * phase hmono = mo_bi_ortho_tc_one_e(p1,h1) * phase

214
src/tc_bi_ortho/slater_tc_opt_single.irp.f
View File

@ -1,6 +1,6 @@
subroutine single_htilde_mu_mat_fock_bi_ortho (Nint, key_j, key_i, hmono, htwoe, htot) subroutine single_htilde_mu_mat_fock_bi_ortho (Nint, key_j, key_i, hmono, htwoe, hthree, htot)
BEGIN_DOC BEGIN_DOC
! <key_j | H_tilde | key_i> for single excitation ONLY FOR ONE- AND TWO-BODY TERMS ! <key_j | H_tilde | key_i> for single excitation ONLY FOR ONE- AND TWO-BODY TERMS
!! !!
@ -14,7 +14,7 @@ subroutine single_htilde_mu_mat_fock_bi_ortho (Nint, key_j, key_i, hmono, htwoe
implicit none implicit none
integer, intent(in) :: Nint integer, intent(in) :: Nint
integer(bit_kind), intent(in) :: key_j(Nint,2), key_i(Nint,2) integer(bit_kind), intent(in) :: key_j(Nint,2), key_i(Nint,2)
double precision, intent(out) :: hmono, htwoe, htot double precision, intent(out) :: hmono, htwoe, hthree, htot
integer :: occ(Nint*bit_kind_size,2) integer :: occ(Nint*bit_kind_size,2)
integer :: Ne(2), i, j, ii, jj, ispin, jspin, k, kk integer :: Ne(2), i, j, ii, jj, ispin, jspin, k, kk
integer :: degree,exc(0:2,2,2) integer :: degree,exc(0:2,2,2)
@ -27,109 +27,21 @@ subroutine single_htilde_mu_mat_fock_bi_ortho (Nint, key_j, key_i, hmono, htwoe
other_spin(1) = 2 other_spin(1) = 2
other_spin(2) = 1 other_spin(2) = 1
hmono = 0.d0 hmono = 0.d0
htwoe= 0.d0 htwoe = 0.d0
htot = 0.d0 hthree = 0.d0
htot = 0.d0
call get_excitation_degree(key_i, key_j, degree, Nint) call get_excitation_degree(key_i, key_j, degree, Nint)
if(degree.ne.1)then if(degree.ne.1)then
return return
endif endif
call bitstring_to_list_ab(key_i, occ, Ne, Nint) call bitstring_to_list_ab(key_i, occ, Ne, Nint)
call get_single_excitation(key_i, key_j, exc, phase, Nint) call get_single_excitation(key_i, key_j, exc, phase, Nint)
call decode_exc(exc,1,h1,p1,h2,p2,s1,s2) call decode_exc(exc,1,h1,p1,h2,p2,s1,s2)
call get_single_excitation_from_fock_tc(key_i,key_j,h1,p1,s1,phase,hmono,htwoe,hthree,htot) call get_single_excitation_from_fock_tc(key_i,key_j,h1,p1,s1,phase,hmono,htwoe,hthree,htot)
end end
BEGIN_PROVIDER [double precision, tc_2e_3idx_coulomb_integrals, (mo_num,mo_num, mo_num)]
&BEGIN_PROVIDER [double precision, tc_2e_3idx_exchange_integrals,(mo_num,mo_num, mo_num)]
implicit none
BEGIN_DOC
! tc_2e_3idx_coulomb_integrals(j,k,i) = <jk|ji>
!
! tc_2e_3idx_exchange_integrals(j,k,i) = <kj|ji>
END_DOC
integer :: i,j,k,l
double precision :: get_two_e_integral
double precision :: integral
do k = 1, mo_num
do i = 1, mo_num
do j = 1, mo_num
tc_2e_3idx_coulomb_integrals(j, k,i) = mo_bi_ortho_tc_two_e(j ,k ,j ,i )
tc_2e_3idx_exchange_integrals(j,k,i) = mo_bi_ortho_tc_two_e(k ,j ,j ,i )
enddo
enddo
enddo
END_PROVIDER
BEGIN_PROVIDER [double precision, fock_op_2_e_tc_closed_shell, (mo_num, mo_num) ]
implicit none
BEGIN_DOC
! Closed-shell part of the Fock operator for the TC operator
END_DOC
integer :: h0,p0,h,p,k0,k,i
integer :: n_occ_ab(2)
integer :: occ(N_int*bit_kind_size,2)
integer :: n_occ_ab_virt(2)
integer :: occ_virt(N_int*bit_kind_size,2)
integer(bit_kind) :: key_test(N_int)
integer(bit_kind) :: key_virt(N_int,2)
call bitstring_to_list_ab(ref_closed_shell_bitmask, occ, n_occ_ab, N_int)
do i = 1, N_int
key_virt(i,1) = full_ijkl_bitmask(i)
key_virt(i,2) = full_ijkl_bitmask(i)
key_virt(i,1) = xor(key_virt(i,1),ref_closed_shell_bitmask(i,1))
key_virt(i,2) = xor(key_virt(i,2),ref_closed_shell_bitmask(i,2))
enddo
call bitstring_to_list_ab(key_virt, occ_virt, n_occ_ab_virt, N_int)
! docc ---> virt single excitations
do h0 = 1, n_occ_ab(1)
h=occ(h0,1)
do p0 = 1, n_occ_ab_virt(1)
p = occ_virt(p0,1)
accu = 0.d0
do k0 = 1, n_occ_ab(1)
k = occ(k0,1)
accu += 2.d0 * tc_2e_3idx_coulomb_integrals(k,p,h) - array_exchange(k,p,h)
enddo
fock_op_2_e_tc_closed_shell(p,h) = accu
enddo
enddo
! virt ---> virt single excitations
do h0 = 1, n_occ_ab_virt(1)
h=occ_virt(h0,1)
do p0 = 1, n_occ_ab_virt(1)
p = occ_virt(p0,1)
accu = 0.d0
do k0 = 1, n_occ_ab(1)
k = occ(k0,1)
accu += 2.d0 * tc_2e_3idx_coulomb_integrals(k,p,h) - array_exchange(k,p,h)
enddo
fock_op_2_e_tc_closed_shell(p,h) = accu
enddo
enddo
! docc ---> docc single excitations
do h0 = 1, n_occ_ab(1)
h=occ(h0,1)
do p0 = 1, n_occ_ab(1)
p = occ(p0,1)
accu = 0.d0
do k0 = 1, n_occ_ab(1)
k = occ(k0,1)
accu += 2.d0 * tc_2e_3idx_coulomb_integrals(k,p,h) - array_exchange(k,p,h)
enddo
fock_op_2_e_tc_closed_shell(p,h) = accu
enddo
enddo
END_PROVIDER
subroutine get_single_excitation_from_fock_tc(key_i,key_j,h,p,spin,phase,hmono,htwoe,hthree,htot) subroutine get_single_excitation_from_fock_tc(key_i,key_j,h,p,spin,phase,hmono,htwoe,hthree,htot)
use bitmasks use bitmasks
@ -200,3 +112,115 @@ subroutine get_single_excitation_from_fock_tc(key_i,key_j,h,p,spin,phase,hmono,h
end end
BEGIN_PROVIDER [double precision, fock_op_2_e_tc_closed_shell, (mo_num, mo_num) ]
implicit none
BEGIN_DOC
! Closed-shell part of the Fock operator for the TC operator
END_DOC
integer :: h0,p0,h,p,k0,k,i
integer :: n_occ_ab(2)
integer :: occ(N_int*bit_kind_size,2)
integer :: n_occ_ab_virt(2)
integer :: occ_virt(N_int*bit_kind_size,2)
integer(bit_kind) :: key_test(N_int)
integer(bit_kind) :: key_virt(N_int,2)
double precision :: accu
fock_op_2_e_tc_closed_shell = -1000.d0
call bitstring_to_list_ab(ref_closed_shell_bitmask, occ, n_occ_ab, N_int)
do i = 1, N_int
key_virt(i,1) = full_ijkl_bitmask(i)
key_virt(i,2) = full_ijkl_bitmask(i)
key_virt(i,1) = xor(key_virt(i,1),ref_closed_shell_bitmask(i,1))
key_virt(i,2) = xor(key_virt(i,2),ref_closed_shell_bitmask(i,2))
enddo
call bitstring_to_list_ab(key_virt, occ_virt, n_occ_ab_virt, N_int)
! docc ---> virt single excitations
do h0 = 1, n_occ_ab(1)
h=occ(h0,1)
do p0 = 1, n_occ_ab_virt(1)
p = occ_virt(p0,1)
accu = 0.d0
do k0 = 1, n_occ_ab(1)
k = occ(k0,1)
accu += 2.d0 * tc_2e_3idx_coulomb_integrals(k,p,h) - tc_2e_3idx_exchange_integrals(k,p,h)
enddo
fock_op_2_e_tc_closed_shell(p,h) = accu
enddo
enddo
do h0 = 1, n_occ_ab_virt(1)
h = occ_virt(h0,1)
do p0 = 1, n_occ_ab(1)
p=occ(p0,1)
accu = 0.d0
do k0 = 1, n_occ_ab(1)
k = occ(k0,1)
accu += 2.d0 * tc_2e_3idx_coulomb_integrals(k,p,h) - tc_2e_3idx_exchange_integrals(k,p,h)
enddo
fock_op_2_e_tc_closed_shell(p,h) = accu
enddo
enddo
! virt ---> virt single excitations
do h0 = 1, n_occ_ab_virt(1)
h=occ_virt(h0,1)
do p0 = 1, n_occ_ab_virt(1)
p = occ_virt(p0,1)
accu = 0.d0
do k0 = 1, n_occ_ab(1)
k = occ(k0,1)
accu += 2.d0 * tc_2e_3idx_coulomb_integrals(k,p,h) - tc_2e_3idx_exchange_integrals(k,p,h)
enddo
fock_op_2_e_tc_closed_shell(p,h) = accu
enddo
enddo
do h0 = 1, n_occ_ab_virt(1)
h = occ_virt(h0,1)
do p0 = 1, n_occ_ab_virt(1)
p=occ_virt(p0,1)
accu = 0.d0
do k0 = 1, n_occ_ab(1)
k = occ(k0,1)
accu += 2.d0 * tc_2e_3idx_coulomb_integrals(k,p,h) - tc_2e_3idx_exchange_integrals(k,p,h)
enddo
fock_op_2_e_tc_closed_shell(p,h) = accu
enddo
enddo
! docc ---> docc single excitations
do h0 = 1, n_occ_ab(1)
h=occ(h0,1)
do p0 = 1, n_occ_ab(1)
p = occ(p0,1)
accu = 0.d0
do k0 = 1, n_occ_ab(1)
k = occ(k0,1)
accu += 2.d0 * tc_2e_3idx_coulomb_integrals(k,p,h) - tc_2e_3idx_exchange_integrals(k,p,h)
enddo
fock_op_2_e_tc_closed_shell(p,h) = accu
enddo
enddo
do h0 = 1, n_occ_ab(1)
h = occ(h0,1)
do p0 = 1, n_occ_ab(1)
p=occ(p0,1)
accu = 0.d0
do k0 = 1, n_occ_ab(1)
k = occ(k0,1)
accu += 2.d0 * tc_2e_3idx_coulomb_integrals(k,p,h) - tc_2e_3idx_exchange_integrals(k,p,h)
enddo
fock_op_2_e_tc_closed_shell(p,h) = accu
enddo
enddo
! do i = 1, mo_num
! write(*,'(100(F10.5,X))')fock_op_2_e_tc_closed_shell(:,i)
! enddo
END_PROVIDER

20
src/tc_bi_ortho/test_tc_bi_ortho.irp.f
View File

@ -19,10 +19,28 @@ subroutine test_slater_tc_opt
integer :: i,j integer :: i,j
double precision :: hmono, htwoe, htot, hthree double precision :: hmono, htwoe, htot, hthree
double precision :: hnewmono, hnewtwoe, hnewthnewree, hnewtot double precision :: hnewmono, hnewtwoe, hnewthnewree, hnewtot
double precision :: accu ,i_count
accu = 0.d0
i_count = 0.d0
do i = 1, N_det do i = 1, N_det
! do i = 14,14
call diag_htilde_mu_mat_bi_ortho(N_int, psi_det(1,1,i), hmono, htwoe, htot) 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) 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) do j = 1, N_det
! do j = 1, 1
if(i==j)cycle
call single_htilde_mu_mat_bi_ortho(N_int, psi_det(1,1,j), psi_det(1,1,i), hmono, htwoe, htot)
call single_htilde_mu_mat_fock_bi_ortho (N_int, psi_det(1,1,j), psi_det(1,1,i), hnewmono, hnewtwoe, hnewthnewree, hnewtot)
if(dabs(htot).gt.1.d-10)then
! if(dabs(htot-hnewtot).gt.1.d-8.or.dabs(htot-hnewtot).gt.dabs(htot))then
print*,j,i
i_count += 1.D0
print*,htot,hnewtot,dabs(htot-hnewtot)
accu += dabs(htot-hnewtot)
! endif
endif
enddo
enddo enddo
print*,'accu = ',accu/i_count
end end