few mom opt

plugins/local/bi_ort_ints/total_twoe_pot.irp.f

@@ -40,16 +40,72 @@ BEGIN_PROVIDER [double precision, mo_bi_ortho_tc_two_e_chemist, (mo_num, mo_num,
   END_DOC
 
   implicit none
-  integer                       :: i, j, k, l, m, n, p, q
+  integer                       :: i, j, k, l, m, n, p, q, s, r
+  double precision              :: t1, t2
   double precision, allocatable :: a1(:,:,:,:), a2(:,:,:,:)
+  double precision, allocatable :: a_jkp(:,:,:), a_kpq(:,:,:), a_pqr(:,:,:)
+
+  print *, ' PROVIDING mo_bi_ortho_tc_two_e_chemist ...'
+  call wall_time(t1)
+  call print_memory_usage()
 
   PROVIDE mo_r_coef mo_l_coef
+  PROVIDe ao_two_e_tc_tot
+
+  if(ao_to_mo_tc_n3) then
+
+    print*, ' memory scale of TC ao -> mo: O(N3) '
+
+    allocate(a_jkp(ao_num,ao_num,mo_num))
+    allocate(a_kpq(ao_num,mo_num,mo_num))
+    allocate(a_pqr(mo_num,mo_num,mo_num))
+
+    do s = 1, mo_num
+      mo_bi_ortho_tc_two_e_chemist(:,:,:,s) = 0.d0
+
+      do l = 1, ao_num
+
+        call dgemm( 'T', 'N', ao_num*ao_num, mo_num, ao_num, 1.d0 &
+                  , ao_two_e_tc_tot(1,1,1,l), ao_num, mo_l_coef(1,1), ao_num  &
+                  , 0.d0, a_jkp(1,1,1), ao_num*ao_num)
+  
+        call dgemm( 'T', 'N', ao_num*mo_num, mo_num, ao_num, 1.d0 &
+                  , a_jkp(1,1,1), ao_num, mo_r_coef(1,1), ao_num  &
+                  , 0.d0, a_kpq(1,1,1), ao_num*mo_num)
+  
+        call dgemm( 'T', 'N', mo_num*mo_num, mo_num, ao_num, 1.d0 &
+                  , a_kpq(1,1,1), ao_num, mo_l_coef(1,1), ao_num  &
+                  , 0.d0, a_pqr(1,1,1), mo_num*mo_num)
+
+        !$OMP PARALLEL         &
+        !$OMP DEFAULT(NONE)    &
+        !$OMP PRIVATE(p, q, r) &
+        !$OMP SHARED(s, l, mo_num, mo_bi_ortho_tc_two_e_chemist, mo_r_coef, a_pqr)
+        !$OMP DO COLLAPSE(2)
+        do p = 1, mo_num
+          do q = 1, mo_num
+            do r = 1, mo_num
+              mo_bi_ortho_tc_two_e_chemist(p,q,r,s) = mo_bi_ortho_tc_two_e_chemist(p,q,r,s) + mo_r_coef(l,s) * a_pqr(p,q,r)
+            enddo
+          enddo
+        enddo
+        !$OMP END DO
+        !$OMP END PARALLEL
+
+      enddo ! l
+    enddo ! s
+
+    deallocate(a_jkp, a_kpq, a_pqr)
+
+  else
+
+    print*, ' memory scale of TC ao -> mo: O(N4) '
 
     allocate(a2(ao_num,ao_num,ao_num,mo_num))
   
     call dgemm( 'T', 'N', ao_num*ao_num*ao_num, mo_num, ao_num, 1.d0     &
               , ao_two_e_tc_tot(1,1,1,1), ao_num, mo_l_coef(1,1), ao_num &
-            , 0.d0 , a2(1,1,1,1), ao_num*ao_num*ao_num)
+              , 0.d0, a2(1,1,1,1), ao_num*ao_num*ao_num)
   
     allocate(a1(ao_num,ao_num,mo_num,mo_num))
   
@@ -72,6 +128,7 @@ BEGIN_PROVIDER [double precision, mo_bi_ortho_tc_two_e_chemist, (mo_num, mo_num,
   
     deallocate(a2)
   
+  endif
 
   !allocate(a1(mo_num,ao_num,ao_num,ao_num))
   !a1 = 0.d0
@@ -135,6 +192,10 @@ BEGIN_PROVIDER [double precision, mo_bi_ortho_tc_two_e_chemist, (mo_num, mo_num,
   !enddo
   !deallocate(a1)
 
+  call wall_time(t2)
+  print *, ' WALL TIME for PROVIDING mo_bi_ortho_tc_two_e_chemist (min)', (t2-t1)/60.d0
+  call print_memory_usage()
+
 END_PROVIDER 
 
 ! ---

plugins/local/non_h_ints_mu/total_tc_int.irp.f

@@ -201,6 +201,8 @@ BEGIN_PROVIDER [double precision, ao_two_e_tc_tot, (ao_num, ao_num, ao_num, ao_n
     !$OMP END DO
     !$OMP END PARALLEL
 
+    call clear_ao_map()
+
     if(tc_integ_type .eq. "numeric") then
       FREE int2_grad1_u12_ao_num int2_grad1_u12_square_ao_num
     endif

src/tools/print_detweights.irp.f

@@ -0,0 +1,66 @@
+program print_detweights
+
+  implicit none
+
+  read_wf = .True.
+  touch read_wf
+
+  call main()
+
+end
+
+! ---
+
+subroutine main()
+
+  implicit none
+  integer                       :: i
+  integer                       :: degree
+  integer                       :: ios
+  integer,          allocatable :: deg(:), ii(:), deg_sorted(:)
+  double precision, allocatable :: c(:)
+
+  PROVIDE N_int
+  PROVIDE N_det
+  PROVIDE psi_det
+  PROVIDe psi_coef
+
+  allocate(deg(N_det), ii(N_det), deg_sorted(N_det), c(N_det))
+
+  do i = 1, N_det
+
+    call debug_det(psi_det(1,1,i), N_int)
+    call get_excitation_degree(psi_det(1,1,i), psi_det(1,1,1), degree, N_int)
+
+    ii (i) = i
+    deg(i) = degree
+    c  (i) = dabs(psi_coef(i,1))
+  enddo
+
+  call dsort(c, ii, N_det)
+
+  do i = 1, N_det
+    deg_sorted(i) = deg(ii(i))
+  enddo
+
+  print *, ' saving psi'
+
+  ! Writing output in binary format
+  open(unit=10, file="coef.bin", status="replace", action="write", iostat=ios, form="unformatted")
+
+    if(ios /= 0) then
+      print *, ' Error opening file!'
+      stop
+    endif
+    
+    write(10) N_det
+    write(10) deg_sorted
+    write(10) c
+
+  close(10)
+
+  deallocate(deg, ii, deg_sorted, c)
+
+end
+
+