added new timing test

src/tc_bi_ortho/slater_tc_opt_double.irp.f

@@ -32,19 +32,23 @@ subroutine double_htilde_mu_mat_fock_bi_ortho(Nint, key_j, key_i, hmono, htwoe,
   if(degree.ne.2)then
    return
   endif
-
+  integer :: degree_i,degree_j
-  call bitstring_to_list_ab(key_i, occ, Ne, Nint)
+  call get_excitation_degree(ref_bitmask,key_i,degree_i,N_int)
+  call get_excitation_degree(ref_bitmask,key_j,degree_j,N_int)
   call get_double_excitation(key_i, key_j, exc, phase, Nint)
   call decode_exc(exc, 2, h1, p1, h2, p2, s1, s2)
 
   if(s1.ne.s2)then
    ! opposite spin two-body 
-!   key_j, key_i
     htwoe  = mo_bi_ortho_tc_two_e(p2,p1,h2,h1) 
     if(three_body_h_tc)then
      if(.not.double_normal_ord)then
-      call three_comp_two_e_elem(key_i,h1,h2,p1,p2,s1,s2,hthree)
+      if(degree_i>degree_j)then
-     elseif(double_normal_ord.and.+Ne(1).gt.2)then
+       call three_comp_two_e_elem(key_j,h1,h2,p1,p2,s1,s2,hthree)
+      else
+       call three_comp_two_e_elem(key_i,h1,h2,p1,p2,s1,s2,hthree)
+      endif
+     elseif(double_normal_ord.and.elec_num+elec_num.gt.2)then
       htwoe += normal_two_body_bi_orth(p2,h2,p1,h1)!!! WTF ???
      endif
     endif
@@ -56,8 +60,12 @@ subroutine double_htilde_mu_mat_fock_bi_ortho(Nint, key_j, key_i, hmono, htwoe,
    htwoe -= mo_bi_ortho_tc_two_e(p1,p2,h2,h1) 
    if(three_body_h_tc)then
     if(.not.double_normal_ord)then
-     call three_comp_two_e_elem(key_i,h1,h2,p1,p2,s1,s2,hthree)
+     if(degree_i>degree_j)then
-    elseif(double_normal_ord.and.+Ne(1).gt.2)then
+      call three_comp_two_e_elem(key_j,h1,h2,p1,p2,s1,s2,hthree)
+     else
+      call three_comp_two_e_elem(key_i,h1,h2,p1,p2,s1,s2,hthree)
+     endif
+    elseif(double_normal_ord.and.elec_num+elec_num.gt.2)then
       htwoe -= normal_two_body_bi_orth(h2,p1,h1,p2)!!! WTF ???
       htwoe += normal_two_body_bi_orth(h1,p1,h2,p2)!!! WTF ???
     endif

src/tc_bi_ortho/test_tc_bi_ortho.irp.f

@@ -11,13 +11,16 @@ program tc_bi_ortho
   touch read_wf
   touch  my_grid_becke my_n_pt_r_grid my_n_pt_a_grid
 
-!  call test_slater_tc_opt
+! call test_slater_tc_opt
- call timing_hij
+ call timing_tot
+! call timing_diag
+! call timing_single
+! call timing_double
 end
 
 subroutine test_slater_tc_opt
  implicit none
- integer :: i,j
+ integer :: i,j,degree
  double precision :: hmono, htwoe, htot, hthree 
  double precision :: hnewmono, hnewtwoe, hnewthree, hnewtot
  double precision :: accu_d ,i_count, accu
@@ -25,84 +28,193 @@ subroutine test_slater_tc_opt
  accu_d = 0.d0
  i_count = 0.d0
  do i = 1, N_det
-! do i = 1,1
-  call diag_htilde_mu_mat_bi_ortho(N_int, psi_det(1,1,i), hmono, htwoe, htot)
-  call htilde_mu_mat_bi_ortho(psi_det(1,1,i), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
-  call diag_htilde_mu_mat_fock_bi_ortho(N_int, psi_det(1,1,i), hnewmono, hnewtwoe, hnewthree, hnewtot)
-!  print*,hthree,hnewthree
-!  print*,htot,hnewtot,dabs(hnewtot-htot)
-  accu_d += dabs(htot-hnewtot) 
-  if(dabs(htot-hnewtot).gt.1.d-8)then
-   print*,i
-   print*,htot,hnewtot,dabs(htot-hnewtot)
-  endif
-!  do j = 319,319
   do j = 1,N_det
-   if(i==j)cycle
-   integer :: degree 
-   call get_excitation_degree(psi_det(1,1,j), psi_det(1,1,i),degree,N_int)
-!   if(degree .ne. 1)cycle
    call htilde_mu_mat_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
    call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hnewmono, hnewtwoe, hnewthree, hnewtot)
-!   if(dabs(hthree).gt.1.d-15)then
    if(dabs(htot).gt.1.d-15)then
-!    if(dabs(htot-hnewtot).gt.1.d-8.or.dabs(htot-hnewtot).gt.dabs(htot))then
      i_count += 1.D0
      accu += dabs(htot-hnewtot) 
-!    if(dabs(hthree-hnewthree).gt.1.d-8.or.dabs(hthree-hnewthree).gt.dabs(hthree))then
+     if(dabs(htot-hnewtot).gt.1.d-8.or.dabs(htot-hnewtot).gt.dabs(htot))then
-    if(dabs(htot-hnewtot).gt.1.d-8.or.dabs(htot-hnewtot).gt.dabs(htot))then
+      call get_excitation_degree(psi_det(1,1,j), psi_det(1,1,i),degree,N_int)
-     print*,j,i,degree
+      print*,j,i,degree
-     call debug_det(psi_det(1,1,i),N_int)
+      call debug_det(psi_det(1,1,i),N_int)
-     call debug_det(psi_det(1,1,j),N_int)
+      call debug_det(psi_det(1,1,j),N_int)
-     print*,htot,hnewtot,dabs(htot-hnewtot) 
+      print*,htot,hnewtot,dabs(htot-hnewtot)
-!     print*,hthree,hnewthree,dabs(hthree-hnewthree) 
+      print*,hthree,hnewthree,dabs(hthree-hnewthree)
-     stop
+      stop
-    endif
+     endif
-!    print*,htot,hnewtot,dabs(htot-hnewtot) 
    endif
   enddo
  enddo
- print*,'accu_d = ',accu_d/N_det
  print*,'accu   = ',accu/i_count
 
 end
 
-
+subroutine timing_tot
-subroutine timing_hij
  implicit none
  integer :: i,j
  double precision :: wall0, wall1
  double precision, allocatable :: mat_old(:,:),mat_new(:,:)
- double precision :: hmono, htwoe, hthree, htot
+ double precision :: hmono, htwoe, hthree, htot, i_count
-! allocate(mat_old(N_det,N_det))
+ integer :: degree 
+ call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,1), psi_det(1,1,2), N_int, hmono, htwoe, hthree, htot)
+ call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,1), psi_det(1,1,2), N_int, hmono, htwoe, hthree, htot)
+ call wall_time(wall0)
+ i_count = 0.d0
+ do i = 1, N_det
+  do j = 1, N_det
+!   call get_excitation_degree(psi_det(1,1,j), psi_det(1,1,i),degree,N_int)
+   i_count += 1.d0
+   call htilde_mu_mat_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
+  enddo
+ enddo
+ call wall_time(wall1)
+ print*,'i_count = ',i_count
+ print*,'time for old hij for total   = ',wall1 - wall0
+
+ call wall_time(wall0)
+ i_count = 0.d0
+ do i = 1, N_det
+  do j = 1, N_det
+!   call get_excitation_degree(psi_det(1,1,j), psi_det(1,1,i),degree,N_int)
+   i_count += 1.d0
+   call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
+  enddo
+ enddo
+ call wall_time(wall1)
+ print*,'i_count = ',i_count
+ print*,'time for new hij for total   = ',wall1 - wall0
+ call i_H_j(psi_det(1,1,1), psi_det(1,1,2),N_int,htot)
+ call wall_time(wall0)
+ i_count = 0.d0
+ do i = 1, N_det
+  do j = 1, N_det
+   call i_H_j(psi_det(1,1,j), psi_det(1,1,i),N_int,htot)
+   i_count += 1.d0
+  enddo
+ enddo
+ call wall_time(wall1)
+ print*,'i_count = ',i_count
+ print*,'time for new hij STANDARD    = ',wall1 - wall0
+
+end
+
+subroutine timing_diag
+ implicit none
+ integer :: i,j
+ double precision :: wall0, wall1
+ double precision, allocatable :: mat_old(:,:),mat_new(:,:)
+ double precision :: hmono, htwoe, hthree, htot, i_count
+ integer :: degree 
  call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,1), psi_det(1,1,1), N_int, hmono, htwoe, hthree, htot)
  call wall_time(wall0)
+ i_count = 0.d0
  do i = 1, N_det
-  do j = 1, N_det
+  do j = i,i 
+   i_count += 1.d0
    call htilde_mu_mat_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
-!   mat_old(j,i) = htot
   enddo
  enddo
  call wall_time(wall1)
- print*,'time for old hij = ',wall1 - wall0
+ print*,'i_count = ',i_count
+ print*,'time for old hij for diagonal= ',wall1 - wall0
 
-! allocate(mat_new(N_det,N_det))
  call wall_time(wall0)
+ i_count = 0.d0
  do i = 1, N_det
-  do j = 1, N_det
+  do j = i,i
+   i_count += 1.d0
    call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
-!   mat_new(j,i) = htot
   enddo
  enddo
  call wall_time(wall1)
- print*,'time for new hij = ',wall1 - wall0
+ print*,'i_count = ',i_count
- double precision :: accu
+ print*,'time for new hij for diagonal= ',wall1 - wall0
+
+end
+
+subroutine timing_single
+ implicit none
+ integer :: i,j
+ double precision :: wall0, wall1,accu
+ double precision, allocatable :: mat_old(:,:),mat_new(:,:)
+ double precision :: hmono, htwoe, hthree, htot, i_count
+ integer :: degree 
+ call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,1), psi_det(1,1,1), N_int, hmono, htwoe, hthree, htot)
+ i_count = 0.d0
  accu = 0.d0
  do i = 1, N_det
   do j = 1, N_det
-!   accu += dabs(mat_new(j,i) - mat_old(j,i))
+   call get_excitation_degree(psi_det(1,1,j), psi_det(1,1,i),degree,N_int)
+   if(degree.ne.1)cycle
+   i_count += 1.d0
+   call wall_time(wall0)
+   call htilde_mu_mat_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
+   call wall_time(wall1)
+   accu += wall1 - wall0
   enddo
  enddo
- print*,'accu = ',accu
+ print*,'i_count = ',i_count
+ print*,'time for old hij for singles = ',accu
+
+ i_count = 0.d0
+ accu = 0.d0
+ do i = 1, N_det
+  do j = 1, N_det
+   call get_excitation_degree(psi_det(1,1,j), psi_det(1,1,i),degree,N_int)
+   if(degree.ne.1)cycle
+   i_count += 1.d0
+   call wall_time(wall0)
+   call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
+   call wall_time(wall1)
+   accu += wall1 - wall0
+  enddo
+ enddo
+ print*,'i_count = ',i_count
+ print*,'time for new hij for singles = ',accu
 
 end
+
+subroutine timing_double
+ implicit none
+ integer :: i,j
+ double precision :: wall0, wall1,accu
+ double precision, allocatable :: mat_old(:,:),mat_new(:,:)
+ double precision :: hmono, htwoe, hthree, htot, i_count
+ integer :: degree 
+ call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,1), psi_det(1,1,1), N_int, hmono, htwoe, hthree, htot)
+ i_count = 0.d0
+ accu = 0.d0
+ do i = 1, N_det
+  do j = 1, N_det
+   call get_excitation_degree(psi_det(1,1,j), psi_det(1,1,i),degree,N_int)
+   if(degree.ne.2)cycle
+   i_count += 1.d0
+   call wall_time(wall0)
+   call htilde_mu_mat_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
+   call wall_time(wall1)
+   accu += wall1 - wall0
+  enddo
+ enddo
+ print*,'i_count = ',i_count
+ print*,'time for old hij for doubles = ',accu
+
+ i_count = 0.d0
+ accu = 0.d0
+ do i = 1, N_det
+  do j = 1, N_det
+   call get_excitation_degree(psi_det(1,1,j), psi_det(1,1,i),degree,N_int)
+   if(degree.ne.2)cycle
+   i_count += 1.d0
+   call wall_time(wall0)
+   call htilde_mu_mat_opt_bi_ortho(psi_det(1,1,j), psi_det(1,1,i), N_int, hmono, htwoe, hthree, htot)
+   call wall_time(wall1)
+   accu += wall1 - wall0
+  enddo
+ enddo
+ call wall_time(wall1)
+ print*,'i_count = ',i_count
+ print*,'time for new hij for doubles = ',accu
+
+end
+