From 0ffaf820a2bc156e021e0f68bddf119e04ce82b5 Mon Sep 17 00:00:00 2001 From: eginer Date: Wed, 16 Aug 2023 18:07:50 +0200 Subject: [PATCH] tc two rdm works for O CISD in biorthonormal basis --- src/tc_bi_ortho/test_tc_two_rdm.irp.f | 19 ++++++++++++++++--- src/tc_bi_ortho/two_rdm_naive.irp.f | 14 +++----------- 2 files changed, 19 insertions(+), 14 deletions(-) diff --git a/src/tc_bi_ortho/test_tc_two_rdm.irp.f b/src/tc_bi_ortho/test_tc_two_rdm.irp.f index ecdeef43..3e556312 100644 --- a/src/tc_bi_ortho/test_tc_two_rdm.irp.f +++ b/src/tc_bi_ortho/test_tc_two_rdm.irp.f @@ -35,6 +35,10 @@ subroutine test do p2 = 1, mo_num integral = mo_bi_ortho_tc_two_e(p2,p1,h2,h1) 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 enddo enddo @@ -42,19 +46,28 @@ subroutine test enddo accu *= 0.5d0 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 do i = 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) 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 +! endif ref += psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) * state_average_weight(istate) * htwoe enddo enddo enddo print*,' ref = ',ref + print*,'delta= ',ref-accu end diff --git a/src/tc_bi_ortho/two_rdm_naive.irp.f b/src/tc_bi_ortho/two_rdm_naive.irp.f index 9694c653..8fd34975 100644 --- a/src/tc_bi_ortho/two_rdm_naive.irp.f +++ b/src/tc_bi_ortho/two_rdm_naive.irp.f @@ -30,7 +30,6 @@ BEGIN_PROVIDER [ double precision, tc_two_rdm, (mo_num, mo_num, mo_num, mo_num)] if(degree == 2)then call update_tc_rdm(h1,p1,h2,p2,s1,s2,tc_two_rdm,mo_num,contrib) else if(degree==1)then -! cycle ! occupation of the determinant psi_det(j) 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) h2 = m p2 = m - if(h2.le.h1)cycle call update_tc_rdm(h1,p1,h2,p2,s1,s2,tc_two_rdm,mo_num,contrib) enddo endif 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) -! print*,'contrib',contrib do istate = 2, N_states contrib += psi_l_coef_bi_ortho(i,istate) * psi_r_coef_bi_ortho(j,istate) * state_average_weight(istate) enddo @@ -115,18 +113,12 @@ subroutine update_tc_rdm(h1,p1,h2,p2,s1,s2,array,sze,contrib) else ! same spin double excitation array(p1,h1,p2,h2) += contrib ! exchange - ! exchanging the holes - array(p2,h1,p1,h2) -= contrib ! exchanging the particles + array(p2,h1,p1,h2) -= contrib + ! exchanging the array(p1,h2,p2,h1) -= contrib - ! permutation for particle symmetry 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 end