diff --git a/src/tc_bi_ortho/test_tc_fock.irp.f b/src/tc_bi_ortho/test_tc_fock.irp.f index 26446daf..ebd43a7a 100644 --- a/src/tc_bi_ortho/test_tc_fock.irp.f +++ b/src/tc_bi_ortho/test_tc_fock.irp.f @@ -15,7 +15,8 @@ program test_tc_fock !call routine_2 ! call routine_3() - call test_3e +! call test_3e + call routine_tot end ! --- @@ -84,8 +85,8 @@ subroutine routine_3() print*, i, a stop endif - !print*, ' excited det' - !call debug_det(det_i, N_int) + print*, ' excited det' + call debug_det(det_i, N_int) call htilde_mu_mat_bi_ortho(det_i, ref_bitmask, N_int, hmono, htwoe, hthree, htilde_ij) if(dabs(hthree).lt.1.d-10)cycle @@ -116,3 +117,78 @@ subroutine routine_3() end subroutine routine_3 ! --- +subroutine routine_tot() + + use bitmasks ! you need to include the bitmasks_module.f90 features + + implicit none + integer :: i, a, i_ok, s1,other_spin(2) + double precision :: hmono, htwoe, hthree, htilde_ij + double precision :: err_ai, err_tot, ref, new + integer(bit_kind), allocatable :: det_i(:,:) + + allocate(det_i(N_int,2)) + other_spin(1) = 2 + other_spin(2) = 1 + + err_tot = 0.d0 + +! do s1 = 1, 2 + s1 = 2 + det_i = ref_bitmask + call debug_det(det_i, N_int) + print*, ' HF det' + call debug_det(det_i, N_int) + +! do i = 1, elec_num_tab(s1) +! do a = elec_num_tab(s1)+1, mo_num ! virtual + do i = 1, elec_beta_num + do a = elec_beta_num+1, elec_alpha_num! virtual +! do i = elec_beta_num+1, elec_alpha_num +! do a = elec_alpha_num+1, mo_num! virtual + print*,i,a + + det_i = ref_bitmask + call do_single_excitation(det_i, i, a, s1, i_ok) + if(i_ok == -1) then + print*, 'PB !!' + print*, i, a + stop + endif + + call htilde_mu_mat_bi_ortho(det_i, ref_bitmask, N_int, hmono, htwoe, hthree, htilde_ij) + print*,htilde_ij + if(dabs(htilde_ij).lt.1.d-10)cycle + print*, ' excited det' + call debug_det(det_i, N_int) + + if(s1 == 1)then + new = Fock_matrix_tc_mo_alpha(a,i) + else + new = Fock_matrix_tc_mo_beta(a,i) + endif + ref = htilde_ij +! if(s1 == 1)then +! new = fock_a_tot_3e_bi_orth(a,i) +! else if(s1 == 2)then +! new = fock_b_tot_3e_bi_orth(a,i) +! endif + err_ai = dabs(dabs(ref) - dabs(new)) + if(err_ai .gt. 1d-7) then + print*,'s1 = ',s1 + print*, ' warning on', i, a + print*, ref,new,err_ai + endif + print*, ref,new,err_ai + err_tot += err_ai + + write(22, *) htilde_ij + enddo + enddo +! enddo + + print *, ' err_tot = ', err_tot + + deallocate(det_i) + +end subroutine routine_3 diff --git a/src/tc_scf/fock_tc.irp.f b/src/tc_scf/fock_tc.irp.f index c3642a7e..2a08e469 100644 --- a/src/tc_scf/fock_tc.irp.f +++ b/src/tc_scf/fock_tc.irp.f @@ -89,7 +89,7 @@ BEGIN_PROVIDER [ double precision, Fock_matrix_tc_mo_alpha, (mo_num, mo_num) ] call ao_to_mo_bi_ortho( Fock_matrix_tc_ao_alpha, size(Fock_matrix_tc_ao_alpha, 1) & , Fock_matrix_tc_mo_alpha, size(Fock_matrix_tc_mo_alpha, 1) ) - if(three_body_h_tc) then + if(three_body_h_tc.and.elec_alpha_num == elec_beta_num) then Fock_matrix_tc_mo_alpha += fock_a_tot_3e_bi_orth endif @@ -116,7 +116,7 @@ BEGIN_PROVIDER [ double precision, Fock_matrix_tc_mo_beta, (mo_num,mo_num) ] call ao_to_mo_bi_ortho( Fock_matrix_tc_ao_beta, size(Fock_matrix_tc_ao_beta, 1) & , Fock_matrix_tc_mo_beta, size(Fock_matrix_tc_mo_beta, 1) ) - if(three_body_h_tc) then + if(three_body_h_tc.and.elec_alpha_num == elec_beta_num) then Fock_matrix_tc_mo_beta += fock_b_tot_3e_bi_orth endif