added providers for the totally symmetrized integrals

src/tc_bi_ortho/slater_tc_opt_diag.irp.f

@@ -113,7 +113,7 @@ subroutine ac_tc_operator(iorb,ispin,key,hmono,htwoe,hthree,Nint,na,nb)
   integer                        :: occ(Nint*bit_kind_size,2)
   integer                        :: other_spin
   integer                        :: k,l,i,jj,mm,j,m
-  double precision :: three_e_diag_parrallel_spin, direct_int, exchange_int
+  double precision ::  direct_int, exchange_int
   
 
   if (iorb < 1) then
@@ -163,7 +163,7 @@ subroutine ac_tc_operator(iorb,ispin,key,hmono,htwoe,hthree,Nint,na,nb)
     jj = occ(j,ispin)
     do m = j+1, na
      mm = occ(m,ispin)
-     hthree += three_e_diag_parrallel_spin(mm,jj,iorb)
+     hthree += three_e_diag_parrallel_spin_prov(mm,jj,iorb)
     enddo
    enddo
    !! same-spin/oposite-spin
@@ -210,7 +210,7 @@ subroutine a_tc_operator(iorb,ispin,key,hmono,htwoe,hthree,Nint,na,nb)
   integer(bit_kind), intent(inout) :: key(Nint,2)
   double precision, intent(inout) :: hmono,htwoe,hthree
   
-  double precision  :: direct_int, exchange_int, three_e_diag_parrallel_spin
+  double precision  :: direct_int, exchange_int
   integer                        :: occ(Nint*bit_kind_size,2)
   integer                        :: other_spin
   integer                        :: k,l,i,jj,mm,j,m
@@ -250,7 +250,7 @@ subroutine a_tc_operator(iorb,ispin,key,hmono,htwoe,hthree,Nint,na,nb)
     jj = occ(j,ispin)
     do m = j+1, na
      mm = occ(m,ispin)
-     hthree -= three_e_diag_parrallel_spin(mm,jj,iorb)
+     hthree -= three_e_diag_parrallel_spin_prov(mm,jj,iorb)
     enddo
    enddo
    !! same-spin/oposite-spin

src/tc_bi_ortho/slater_tc_opt_double.irp.f

@@ -84,7 +84,7 @@ subroutine three_comp_two_e_elem(key_i,h1,h2,p1,p2,s1,s2,hthree)
  integer :: n_occ_ab_hole(2),n_occ_ab_particle(2)
  integer(bit_kind)              :: det_tmp(N_int,2)
  integer :: ipart, ihole
- double precision :: direct_int, exchange_int, three_e_double_parrallel_spin
+ double precision :: direct_int, exchange_int
 
   nexc(1) = 0
   nexc(2) = 0
@@ -118,9 +118,9 @@ subroutine three_comp_two_e_elem(key_i,h1,h2,p1,p2,s1,s2,hthree)
     ispin = 1 ! i==alpha ==> pure same spin terms
     do i = 1, nexc(ispin) ! number of couple of holes/particles 
      ipart=occ_particle(i,ispin)
-     hthree += three_e_double_parrallel_spin(ipart,p2,h2,p1,h1)
+     hthree += three_e_double_parrallel_spin_prov(ipart,p2,h2,p1,h1)
      ihole=occ_hole(i,ispin)
-     hthree -= three_e_double_parrallel_spin(ihole,p2,h2,p1,h1)
+     hthree -= three_e_double_parrallel_spin_prov(ihole,p2,h2,p1,h1)
     enddo
     ispin = 2 ! i==beta ==> alpha/alpha/beta terms
     do i = 1, nexc(ispin) ! number of couple of holes/particles 
@@ -145,9 +145,9 @@ subroutine three_comp_two_e_elem(key_i,h1,h2,p1,p2,s1,s2,hthree)
    ispin = 2 ! i==beta  ==> pure same spin terms
    do i = 1, nexc(ispin) ! number of couple of holes/particles 
     ipart=occ_particle(i,ispin)
-    hthree += three_e_double_parrallel_spin(ipart,p2,h2,p1,h1)
+    hthree += three_e_double_parrallel_spin_prov(ipart,p2,h2,p1,h1)
     ihole=occ_hole(i,ispin)
-    hthree -= three_e_double_parrallel_spin(ihole,p2,h2,p1,h1)
+    hthree -= three_e_double_parrallel_spin_prov(ihole,p2,h2,p1,h1)
    enddo
    ispin = 1 ! i==alpha==> beta/beta/alpha terms
    do i = 1, nexc(ispin) ! number of couple of holes/particles 
@@ -305,13 +305,12 @@ subroutine give_contrib_for_aaaa(h1,h2,p1,p2,occ,Ne,contrib)
  double precision, intent(out) :: contrib
  integer :: mm,m 
  double precision :: direct_int, exchange_int
- double precision :: three_e_double_parrallel_spin
  !! h1,p1 == alpha 
  !! h2,p2 == alpha
  contrib = 0.d0
  do mm = 1, Ne(1) !! alpha ==> pure parallele spin contribution
    m = occ(mm,1)
-   contrib += three_e_double_parrallel_spin(m,p2,h2,p1,h1)
+   contrib += three_e_double_parrallel_spin_prov(m,p2,h2,p1,h1)
  enddo
 
  do mm = 1, Ne(2) !! beta
@@ -371,13 +370,12 @@ subroutine give_contrib_for_bbbb(h1,h2,p1,p2,occ,Ne,contrib)
  double precision, intent(out) :: contrib
  integer :: mm,m 
  double precision :: direct_int, exchange_int
- double precision :: three_e_double_parrallel_spin
  !! h1,p1 == beta
  !! h2,p2 == beta
  contrib = 0.d0
  do mm = 1, Ne(2) !! beta ==> pure parallele spin contribution
    m = occ(mm,1)
-   contrib += three_e_double_parrallel_spin(m,p2,h2,p1,h1)
+   contrib += three_e_double_parrallel_spin_prov(m,p2,h2,p1,h1)
  enddo
 
  do mm = 1, Ne(1) !! alpha

src/tc_bi_ortho/slater_tc_opt_single.irp.f

@@ -196,7 +196,7 @@ subroutine fock_ac_tc_operator(iorb,ispin,key, h_fock,p_fock, ispin_fock,hthree,
   integer                        :: occ(Nint*bit_kind_size,2)
   integer                        :: other_spin
   integer                        :: k,l,i,jj,j
-  double precision :: three_e_single_parrallel_spin, direct_int, exchange_int
+  double precision :: direct_int, exchange_int
   
 
   if (iorb < 1) then
@@ -236,7 +236,7 @@ subroutine fock_ac_tc_operator(iorb,ispin,key, h_fock,p_fock, ispin_fock,hthree,
    !! jj = ispin = ispin_fock >> pure parallel spin
    do j = 1, na
     jj = occ(j,ispin)
-    hthree += three_e_single_parrallel_spin(jj,iorb,p_fock,h_fock)
+    hthree += three_e_single_parrallel_spin_prov(jj,iorb,p_fock,h_fock)
    enddo
    !! spin of jj == other spin than ispin AND ispin_fock
    !! exchange between the iorb and (h_fock, p_fock)
@@ -287,7 +287,7 @@ subroutine fock_a_tc_operator(iorb,ispin,key, h_fock,p_fock, ispin_fock,hthree,N
   integer(bit_kind), intent(inout) :: key(Nint,2)
   double precision, intent(inout) :: hthree
   
-  double precision  :: direct_int, exchange_int, three_e_single_parrallel_spin
+  double precision  :: direct_int, exchange_int
   integer                        :: occ(Nint*bit_kind_size,2)
   integer                        :: other_spin
   integer                        :: k,l,i,jj,mm,j,m
@@ -315,7 +315,7 @@ subroutine fock_a_tc_operator(iorb,ispin,key, h_fock,p_fock, ispin_fock,hthree,N
    !! jj = ispin = ispin_fock >> pure parallel spin
    do j = 1, na
     jj = occ(j,ispin)
-    hthree -= three_e_single_parrallel_spin(jj,iorb,p_fock,h_fock)
+    hthree -= three_e_single_parrallel_spin_prov(jj,iorb,p_fock,h_fock)
    enddo
    !! spin of jj == other spin than ispin AND ispin_fock
    !! exchange between the iorb and (h_fock, p_fock)

src/tc_bi_ortho/symmetrized_3_e_int_prov.irp.f

@@ -0,0 +1,140 @@
+
+BEGIN_PROVIDER [ double precision, three_e_diag_parrallel_spin_prov, (mo_num, mo_num, mo_num)]
+
+  BEGIN_DOC
+  !
+  ! matrix element of the -L  three-body operator ON A BI ORTHONORMAL BASIS 
+  !
+  ! three_e_diag_parrallel_spin_prov(m,j,i) = All combinations of the form <mji|-L|mji> for same spin matrix elements  
+  ! 
+  ! notice the -1 sign: in this way three_e_diag_parrallel_spin_prov can be directly used to compute Slater rules with a + sign
+  !
+  END_DOC
+
+  implicit none
+  integer          :: i, j, m
+  double precision :: integral, wall1, wall0, three_e_diag_parrallel_spin
+
+  three_e_diag_parrallel_spin_prov = 0.d0
+  print *, ' Providing the three_e_diag_parrallel_spin_prov ...'
+
+ integral = three_e_diag_parrallel_spin(1,1,1) ! to provide all stuffs
+  call wall_time(wall0)
+ !$OMP PARALLEL                 &
+ !$OMP DEFAULT (NONE)           &
+ !$OMP PRIVATE (i,j,m,integral) & 
+ !$OMP SHARED (mo_num,three_e_diag_parrallel_spin_prov)
+ !$OMP DO SCHEDULE (dynamic)
+  do i = 1, mo_num
+    do j = 1, mo_num
+      do m = j, mo_num
+        three_e_diag_parrallel_spin_prov(m,j,i) =  three_e_diag_parrallel_spin(m,j,i)
+      enddo
+    enddo
+  enddo
+ !$OMP END DO
+ !$OMP END PARALLEL
+
+  do i = 1, mo_num
+    do j = 1, mo_num
+      do m = 1, j
+        three_e_diag_parrallel_spin_prov(m,j,i) = three_e_diag_parrallel_spin_prov(j,m,i)
+      enddo
+    enddo
+  enddo
+
+  call wall_time(wall1)
+  print *, ' wall time for three_e_diag_parrallel_spin_prov', wall1 - wall0
+
+END_PROVIDER 
+
+BEGIN_PROVIDER [ double precision, three_e_single_parrallel_spin_prov, (mo_num, mo_num, mo_num, mo_num)]
+
+  BEGIN_DOC
+  !
+  ! matrix element of the -L  three-body operator FOR THE DIRECT TERMS OF SINGLE EXCITATIONS AND BI ORTHO MOs
+  !
+  ! three_e_single_parrallel_spin_prov(m,j,k,i) = All combination of <mjk|-L|mji> for same spin matrix elements 
+  !
+  ! notice the -1 sign: in this way three_e_3_idx_direct_bi_ort can be directly used to compute Slater rules with a + sign
+  !
+  END_DOC
+
+ implicit none
+ integer          :: i, j, k, m
+ double precision :: integral, wall1, wall0, three_e_single_parrallel_spin
+
+  three_e_single_parrallel_spin_prov = 0.d0
+  print *, ' Providing the three_e_single_parrallel_spin_prov ...'
+
+  integral = three_e_single_parrallel_spin(1,1,1,1)
+  call wall_time(wall0)
+ !$OMP PARALLEL                   &
+ !$OMP DEFAULT (NONE)             &
+ !$OMP PRIVATE (i,j,k,m,integral) & 
+ !$OMP SHARED (mo_num,three_e_single_parrallel_spin_prov)
+ !$OMP DO SCHEDULE (dynamic)
+  do i = 1, mo_num
+    do k = 1, mo_num
+      do j = 1, mo_num
+        do m = 1, mo_num
+          three_e_single_parrallel_spin_prov(m,j,k,i) = three_e_single_parrallel_spin(m,j,k,i)
+        enddo
+      enddo
+    enddo
+  enddo
+ !$OMP END DO
+ !$OMP END PARALLEL
+
+  call wall_time(wall1)
+  print *, ' wall time for three_e_single_parrallel_spin_prov', wall1 - wall0
+
+END_PROVIDER 
+
+
+! ---
+
+BEGIN_PROVIDER [ double precision, three_e_double_parrallel_spin_prov, (mo_num, mo_num, mo_num, mo_num, mo_num)]
+
+  BEGIN_DOC
+  !
+  ! matrix element of the -L  three-body operator FOR THE DIRECT TERMS OF DOUBLE EXCITATIONS AND BI ORTHO MOs
+  !
+  ! three_e_double_parrallel_spin_prov(m,l,j,k,i) = <mlk|-L|mji> ::: notice that i is the RIGHT MO and k is the LEFT MO
+  !
+  ! notice the -1 sign: in this way three_e_3_idx_direct_bi_ort can be directly used to compute Slater rules with a + sign
+  END_DOC
+
+  implicit none
+  integer          :: i, j, k, m, l
+  double precision :: integral, wall1, wall0, three_e_double_parrallel_spin
+
+  three_e_double_parrallel_spin_prov = 0.d0
+  print *, ' Providing the three_e_double_parrallel_spin_prov ...'
+  call wall_time(wall0)
+
+ integral = three_e_double_parrallel_spin(1,1,1,1,1)
+ !$OMP PARALLEL                     &
+ !$OMP DEFAULT (NONE)               &
+ !$OMP PRIVATE (i,j,k,m,l,integral) & 
+ !$OMP SHARED (mo_num,three_e_double_parrallel_spin_prov)
+ !$OMP DO SCHEDULE (dynamic)
+  do i = 1, mo_num
+    do k = 1, mo_num
+      do j = 1, mo_num
+        do l = 1, mo_num
+          do m = 1, mo_num
+            three_e_double_parrallel_spin_prov(m,l,j,k,i) = three_e_double_parrallel_spin(m,l,j,k,i)
+          enddo
+        enddo
+      enddo
+    enddo
+  enddo
+ !$OMP END DO
+ !$OMP END PARALLEL
+
+  call wall_time(wall1)
+  print *, ' wall time for three_e_double_parrallel_spin_prov', wall1 - wall0
+
+END_PROVIDER 
+

src/tc_bi_ortho/test_tc_bi_ortho.irp.f

@@ -57,7 +57,7 @@ subroutine test_slater_tc_opt
 !     print*,hthree,hnewthree,dabs(hthree-hnewthree) 
      stop
     endif
-    print*,htot,hnewtot,dabs(htot-hnewtot) 
+!    print*,htot,hnewtot,dabs(htot-hnewtot) 
    endif
   enddo
  enddo