Working libint

plugins/Hartree_Fock/debug_libinit.irp.f

@@ -80,10 +80,11 @@
                            !Compute the norm
                            double precision:: coef1, coef2, coef3, coef4, norm
 
-                           coef1 = ao_coef_normalization_factor(bf1)
+                           coef1 = ao_coef_normalization_libint_factor(bf1)
-                           coef2 = ao_coef_normalization_factor(bf2)
+                           coef2 = ao_coef_normalization_libint_factor(bf2) 
-                           coef3 = ao_coef_normalization_factor(bf3)
+                           coef3 = ao_coef_normalization_libint_factor(bf3) 
-                           coef4 = ao_coef_normalization_factor(bf4)
+                           coef4 = ao_coef_normalization_libint_factor(bf4)
+
                            norm = coef1*coef2*coef3*coef4
    
                            double precision:: libint, ref
@@ -95,9 +96,14 @@
 !                           ref = ao_bielec_integral(bf1,bf2,bf3,bf4)
 !   
 !                           if ( (ABS(ABS(ref) - ABS(libint)).GE.1e-6) )  THEN
+!
 !                              print*, bf1,bf2,bf3,bf4
-!                              print*, ref
+!                              print*,"r", ref
-!                              print*, libint
+!                              print*,"l", libint
+!                              print*,"r/l", ref/libint
+!                              print*,"l/r", libint/ref
+!                              print*,"n", norm
+!
 !                              call exit(1)
 !                           end if
 

src/AO_Basis/aos.irp.f

@@ -25,7 +25,7 @@ END_PROVIDER
   BEGIN_DOC
   ! Coefficients including the AO normalization
   END_DOC
-  double precision               :: norm, norm2,overlap_x,overlap_y,overlap_z,C_A(3), c
+  double precision               :: norm,overlap_x,overlap_y,overlap_z,C_A(3), c
   integer                        :: l, powA(3), nz
   integer                        :: i,j,k
   nz=100
@@ -34,9 +34,11 @@ END_PROVIDER
   C_A(3) = 0.d0
   ao_coef_normalized = 0.d0
   do i=1,ao_num
-    powA(1) = ao_l(i)
+
-    powA(2) = 0
+    powA(1) = ao_power(i,1)
-    powA(3) = 0
+    powA(2) = ao_power(i,2)
+    powA(3) = ao_power(i,3)
+
     do j=1,ao_prim_num(i)
       call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j),powA,powA,overlap_x,overlap_y,overlap_z,norm,nz)
       ao_coef_normalized(i,j) = ao_coef(i,j)/sqrt(norm)
@@ -57,6 +59,39 @@ END_PROVIDER
 
 END_PROVIDER
 
+BEGIN_PROVIDER [ double precision, ao_coef_normalization_libint_factor, (ao_num) ]
+  implicit none
+  BEGIN_DOC
+  ! Coefficients including the AO normalization
+  END_DOC
+  double precision               :: norm,overlap_x,overlap_y,overlap_z,C_A(3), c
+  integer                        :: l, powA(3), nz
+  integer                        :: i,j,k
+  nz=100
+  C_A(1) = 0.d0
+  C_A(2) = 0.d0
+  C_A(3) = 0.d0
+
+  do i=1,ao_num
+    powA(1) = ao_l(i)
+    powA(2) = 0
+    powA(3) = 0
+ 
+    ! Normalization of the contracted basis functions
+    norm = 0.d0
+    do j=1,ao_prim_num(i)
+     do k=1,ao_prim_num(i)
+      call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,k),powA,powA,overlap_x,overlap_y,overlap_z,c,nz)
+      norm = norm+c*ao_coef_normalized(i,j)*ao_coef_normalized(i,k)
+     enddo
+    enddo
+    ao_coef_normalization_libint_factor(i) = ao_coef_normalization_factor(i) * sqrt(norm)
+
+  enddo
+
+END_PROVIDER
+
+
  BEGIN_PROVIDER [ double precision, ao_coef_normalized_ordered, (ao_num_align,ao_prim_num_max) ]
 &BEGIN_PROVIDER [ double precision, ao_expo_ordered, (ao_num_align,ao_prim_num_max) ]
    implicit none

src/Integrals_Bielec/ao_bi_integrals.irp.f

@@ -5,11 +5,6 @@ double precision function ao_bielec_integral(i,j,k,l)
   !     i(r1) j(r1) 1/r12 k(r2) l(r2)
   END_DOC
 
- ! use libint_module
- ! call init_libint(trim(ezfio_filename)//char(0))
- ! ao_bielec_integral = ao_bielec_integral_libint(i,j,k,l)
- ! return
-
   integer,intent(in)             :: i,j,k,l
   integer                        :: p,q,r,s
   double precision               :: I_center(3),J_center(3),K_center(3),L_center(3)