PTEROSOR/QuantumPackage
10
0
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-12-22 20:34:58 +01:00
Code Releases Activity

tc two rdm works for O CISD in biorthonormal basis

Browse Source
This commit is contained in:
eginer 2023-08-16 18:07:50 +02:00
parent 251a71ce23
commit 0ffaf820a2
2 changed files with 19 additions and 14 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

19
src/tc_bi_ortho/test_tc_two_rdm.irp.f
View File

@ -35,6 +35,10 @@ subroutine test
do p2 = 1, mo_num do p2 = 1, mo_num
integral = mo_bi_ortho_tc_two_e(p2,p1,h2,h1) integral = mo_bi_ortho_tc_two_e(p2,p1,h2,h1)
rdm = tc_two_rdm(p1,h1,p2,h2) rdm = tc_two_rdm(p1,h1,p2,h2)
! if(dabs(rdm).gt.1.d-10)then
! print*,h1,p1,h2,p2
! print*,rdm,integral,rdm*integral
! endif
accu += integral * rdm accu += integral * rdm
enddo enddo
enddo enddo
@ -42,19 +46,28 @@ subroutine test
enddo enddo
accu *= 0.5d0 accu *= 0.5d0
print*,'accu = ',accu print*,'accu = ',accu
! print*,tc_two_rdm(1,1,1,1),mo_bi_ortho_tc_two_e(1,1,1,1) ! print*,mo_bi_ortho_tc_two_e(2,15,2,1)
! print*,mo_bi_ortho_tc_two_e(15,2,2,1)
! print*,mo_bi_ortho_tc_two_e(2,1,2,15)
! print*,mo_bi_ortho_tc_two_e(2,1,15,2)
ref = 0.d0 ref = 0.d0
do i = 1, N_det do i = 1, N_det
do j = 1, N_det do j = 1, N_det
! if(i.ne.j)cycle ! if(i.eq.j)cycle
call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,i), psi_det(1,1,j), N_int, hmono, htwoe, hthree, htot) call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,i), psi_det(1,1,j), N_int, hmono, htwoe, hthree, htot)
do istate = 1,N_states do istate = 1,N_states
! print*,'i,j',i,j,psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) * state_average_weight(istate) * htwoe ! print*,'i,j',i,j
! print*,psi_l_coef_bi_ortho(i,istate) , psi_r_coef_bi_ortho(j,istate) , htwoe
! print*,psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) * htwoe
! if(i.ne.j)then
! print*,psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) , htwoe ! print*,psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) , htwoe
! print*,psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) * htwoe
! endif
ref += psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) * state_average_weight(istate) * htwoe ref += psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) * state_average_weight(istate) * htwoe
enddo enddo
enddo enddo
enddo enddo
print*,' ref = ',ref print*,' ref = ',ref
print*,'delta= ',ref-accu
end end

14
src/tc_bi_ortho/two_rdm_naive.irp.f
View File

@ -30,7 +30,6 @@ BEGIN_PROVIDER [ double precision, tc_two_rdm, (mo_num, mo_num, mo_num, mo_num)]
if(degree == 2)then if(degree == 2)then
call update_tc_rdm(h1,p1,h2,p2,s1,s2,tc_two_rdm,mo_num,contrib) call update_tc_rdm(h1,p1,h2,p2,s1,s2,tc_two_rdm,mo_num,contrib)
else if(degree==1)then else if(degree==1)then
! cycle
! occupation of the determinant psi_det(j) ! occupation of the determinant psi_det(j)
call bitstring_to_list_ab(psi_det(1,1,j), occ, n_occ_ab, N_int) call bitstring_to_list_ab(psi_det(1,1,j), occ, n_occ_ab, N_int)
@ -48,13 +47,12 @@ BEGIN_PROVIDER [ double precision, tc_two_rdm, (mo_num, mo_num, mo_num, mo_num)]
m = occ(mm,s2) m = occ(mm,s2)
h2 = m h2 = m
p2 = m p2 = m
if(h2.le.h1)cycle
call update_tc_rdm(h1,p1,h2,p2,s1,s2,tc_two_rdm,mo_num,contrib) call update_tc_rdm(h1,p1,h2,p2,s1,s2,tc_two_rdm,mo_num,contrib)
enddo enddo
endif endif
else if(degree == 0)then else if(degree == 0)then
! cycle
contrib = psi_l_coef_bi_ortho(i,1) * psi_r_coef_bi_ortho(j,1) * state_average_weight(1) contrib = psi_l_coef_bi_ortho(i,1) * psi_r_coef_bi_ortho(j,1) * state_average_weight(1)
! print*,'contrib',contrib
do istate = 2, N_states do istate = 2, N_states
contrib += psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) * state_average_weight(istate) contrib += psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) * state_average_weight(istate)
enddo enddo
@ -115,18 +113,12 @@ subroutine update_tc_rdm(h1,p1,h2,p2,s1,s2,array,sze,contrib)
else ! same spin double excitation else ! same spin double excitation
array(p1,h1,p2,h2) += contrib array(p1,h1,p2,h2) += contrib
! exchange ! exchange
! exchanging the holes
array(p2,h1,p1,h2) -= contrib
! exchanging the particles ! exchanging the particles
array(p2,h1,p1,h2) -= contrib
! exchanging the
array(p1,h2,p2,h1) -= contrib array(p1,h2,p2,h1) -= contrib
! permutation for particle symmetry ! permutation for particle symmetry
array(p2,h2,p1,h1) += contrib array(p2,h2,p1,h1) += contrib
! exchange
! exchanging the holes
array(p1,h2,p2,h1) -= contrib
! exchanging the particles
array(p2,h1,p1,h2) -= contrib
endif endif
end end