Merge branch 'master' of github.com:scemama/quantum_package

TODO

@@ -13,6 +13,7 @@
   * dft_utils_one_body
   * dm_for_dft
   * integrals_erf
+* prendre des bouts du src/README.rst  et en mettre partout 
 
 
 

data/module_gitignore

@@ -45,3 +45,4 @@ slave
 tools
 utils
 zmq
+

src/ao_two_e_erf_integrals/EZFIO.cfg

@@ -0,0 +1,13 @@
+[disk_access_ao_integrals_erf]
+type: Disk_access
+doc: Read/Write |AO| integrals with the long range interaction from/to disk [ Write | Read | None ]
+interface: ezfio,provider,ocaml
+default: None
+
+[mu_erf]
+type: double precision
+doc: cutting of the interaction in the range separated model
+interface: ezfio,provider,ocaml
+default: 0.5
+ezfio_name: mu_erf
+

src/ao_two_e_erf_integrals/NEED

@@ -0,0 +1 @@
+ao_two_e_integrals

src/ao_two_e_erf_integrals/providers_ao_erf.irp.f

@@ -90,7 +90,7 @@ BEGIN_PROVIDER [ logical, ao_bielec_integrals_erf_in_map ]
   if (write_ao_integrals_erf) then
     call ezfio_set_work_empty(.False.)
     call map_save_to_disk(trim(ezfio_filename)//'/work/ao_ints_erf',ao_integrals_erf_map)
-    call ezfio_set_ao_two_e_integrals_disk_access_ao_integrals_erf("Read")
+    call ezfio_set_ao_two_e_erf_integrals_disk_access_ao_integrals_erf("Read")
   endif
   
 END_PROVIDER

src/ao_two_e_erf_integrals/read_write_erf.irp.f

@@ -0,0 +1,27 @@
+ BEGIN_PROVIDER [ logical, read_ao_integrals_erf ]
+&BEGIN_PROVIDER [ logical, write_ao_integrals_erf ]
+   implicit none
+   
+   BEGIN_DOC
+   ! Flag to read or write the |AO| erf integrals
+   END_DOC
+   
+   if (disk_access_ao_integrals_erf.EQ.'Read') then
+     read_ao_integrals_erf =  .True.
+     write_ao_integrals_erf = .False.
+     
+   else if  (disk_access_ao_integrals_erf.EQ.'Write') then
+     read_ao_integrals_erf = .False.
+     write_ao_integrals_erf =  .True.
+     
+   else if (disk_access_ao_integrals_erf.EQ.'None') then
+     read_ao_integrals_erf = .False.
+     write_ao_integrals_erf = .False.
+     
+   else
+     print *, 'disk_access_ao_integrals_erf has a wrong type'
+     stop 1
+     
+   endif
+   
+END_PROVIDER

src/ao_two_e_erf_integrals/routines_save_integrals_erf.irp.f

@@ -4,7 +4,7 @@ subroutine save_erf_bi_elec_integrals_ao
  PROVIDE ao_bielec_integrals_erf_in_map
  call ezfio_set_work_empty(.False.)
  call map_save_to_disk(trim(ezfio_filename)//'/work/ao_ints_erf',ao_integrals_erf_map)
- call ezfio_set_ao_two_e_integrals_disk_access_ao_integrals('Read')
+ call ezfio_set_ao_two_e_erf_integrals_disk_access_ao_integrals_erf('Read')
 end
 
 subroutine save_erf_bielec_ints_ao_into_ints_ao

src/ao_two_e_integrals/EZFIO.cfg

@@ -4,12 +4,6 @@ doc: Read/Write |AO| integrals from/to disk [ Write | Read | None ]
 interface: ezfio,provider,ocaml
 default: None
 
-[disk_access_ao_integrals_erf]
-type: Disk_access
-doc: Read/Write |AO| integrals with the long range interaction from/to disk [ Write | Read | None ]
-interface: ezfio,provider,ocaml
-default: None
-
 [ao_integrals_threshold]
 type: Threshold
 doc: If | (pq|rs) | < `ao_integrals_threshold` then (pq|rs) is zero
@@ -24,10 +18,3 @@ interface: ezfio,provider,ocaml
 default: False
 ezfio_name: direct
 
-[mu_erf]
-type: double precision
-doc: cutting of the interaction in the range separated model
-interface: ezfio,provider,ocaml
-default: 0.5
-ezfio_name: mu_erf
-

src/ao_two_e_integrals/read_write.irp.f

@@ -24,32 +24,4 @@
      
    endif
    
-END_PROVIDER
- 
- BEGIN_PROVIDER [ logical, read_ao_integrals_erf ]
-&BEGIN_PROVIDER [ logical, write_ao_integrals_erf ]
-   implicit none
-   
-   BEGIN_DOC
-   ! Flag to read or write the |AO| erf integrals
-   END_DOC
-   
-   if (disk_access_ao_integrals_erf.EQ.'Read') then
-     read_ao_integrals_erf =  .True.
-     write_ao_integrals_erf = .False.
-     
-   else if  (disk_access_ao_integrals_erf.EQ.'Write') then
-     read_ao_integrals_erf = .False.
-     write_ao_integrals_erf =  .True.
-     
-   else if (disk_access_ao_integrals_erf.EQ.'None') then
-     read_ao_integrals_erf = .False.
-     write_ao_integrals_erf = .False.
-     
-   else
-     print *, 'disk_access_ao_integrals_erf has a wrong type'
-     stop 1
-     
-   endif
-   
 END_PROVIDER

src/aux_quantities/save_one_body_dm.irp.f

@@ -1,7 +1,7 @@
 program save_one_body_dm
   implicit none
  BEGIN_DOC 
-! programs that computes the one body density on the mo basis for alpha and beta electrons from the wave function stored in the EZFIO folder, and then save it into the EZFIO folder data_energy_and_density. 
+! programs that computes the one body density on the mo basis for alpha and beta electrons from the wave function stored in the EZFIO folder, and then save it into the EZFIO folder aux_quantities. 
 ! 
 ! Then, the global variable data_one_body_alpha_dm_mo and data_one_body_beta_dm_mo will automatically read the density in a further calculation. 
 !
@@ -15,6 +15,6 @@ end
 
 subroutine routine
  
- call ezfio_set_data_energy_and_density_data_one_body_alpha_dm_mo(one_body_dm_mo_alpha)
+ call ezfio_set_aux_quantities_data_one_body_alpha_dm_mo(one_body_dm_mo_alpha)
- call ezfio_set_data_energy_and_density_data_one_body_beta_dm_mo(one_body_dm_mo_beta)
+ call ezfio_set_aux_quantities_data_one_body_beta_dm_mo(one_body_dm_mo_beta)
 end

src/davidson/u0_wee_u0.irp.f

@@ -304,7 +304,7 @@ subroutine H_S2_u_0_bielec_nstates_openmp_work_$N_int(v_t,s_t,u_t,N_st,sze,istar
       ASSERT (lrow <= N_det_alpha_unique)
 
       tmp_det2(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, lrow)
-      call i_H_j_mono_spin_bielec( tmp_det, tmp_det2, $N_int, 1, hij) 
+      call i_Wee_j_mono( tmp_det, tmp_det2, $N_int, 1, hij) 
 
       do l=1,N_st
         v_t(l,k_a) = v_t(l,k_a) + hij * u_t(l,l_a)
@@ -384,7 +384,7 @@ subroutine H_S2_u_0_bielec_nstates_openmp_work_$N_int(v_t,s_t,u_t,N_st,sze,istar
       ASSERT (lcol <= N_det_beta_unique)
 
       tmp_det2(1:$N_int,2) = psi_det_beta_unique (1:$N_int, lcol)
-      call i_H_j_mono_spin_bielec( tmp_det, tmp_det2, $N_int, 2, hij)
+      call i_Wee_j_mono( tmp_det, tmp_det2, $N_int, 2, hij)
       l_a = psi_bilinear_matrix_transp_order(l_b)
       ASSERT (l_a <= N_det)
       do l=1,N_st
@@ -430,9 +430,9 @@ subroutine H_S2_u_0_bielec_nstates_openmp_work_$N_int(v_t,s_t,u_t,N_st,sze,istar
     tmp_det(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, krow)
     tmp_det(1:$N_int,2) = psi_det_beta_unique (1:$N_int, kcol)
     
-    double precision, external :: diag_H_mat_elem_bielec, diag_S_mat_elem
+    double precision, external :: diag_wee_mat_elem, diag_S_mat_elem
   
-    hij = diag_H_mat_elem_bielec(tmp_det,$N_int) 
+    hij = diag_wee_mat_elem(tmp_det,$N_int) 
     sij = diag_S_mat_elem(tmp_det,$N_int)
     do l=1,N_st
       v_t(l,k_a) = v_t(l,k_a) + hij * u_t(l,k_a)

src/density_for_dft/NEED

@@ -1,3 +1,3 @@
 determinants 
 dft_keywords 
-data_energy_and_density
+aux_quantities

src/determinants/mono_excitations_bielec.irp.f

@@ -1,5 +1,5 @@
   use bitmasks
-subroutine get_mono_excitation_from_fock_bielec(det_1,det_2,h,p,spin,phase,hij)
+subroutine mono_excitation_wee(det_1,det_2,h,p,spin,phase,hij)
  use bitmasks
  implicit none
  integer,intent(in) :: h,p,spin
@@ -23,7 +23,7 @@ subroutine get_mono_excitation_from_fock_bielec(det_1,det_2,h,p,spin,phase,hij)
  enddo
  call bitstring_to_list_ab(hole, occ_hole, n_occ_ab_hole, N_int)
  call bitstring_to_list_ab(partcl, occ_partcl, n_occ_ab_partcl, N_int)
- hij = fock_operator_bielec_closed_shell_ref_bitmask(h,p)
+ hij = fock_wee_closed_shell(h,p)
  ! holes :: direct terms
  do i0 = 1, n_occ_ab_hole(1)
   i = occ_hole(i0,1)
@@ -60,7 +60,7 @@ subroutine get_mono_excitation_from_fock_bielec(det_1,det_2,h,p,spin,phase,hij)
 end
 
 
-BEGIN_PROVIDER [double precision, fock_operator_bielec_closed_shell_ref_bitmask, (mo_tot_num, mo_tot_num) ]
+BEGIN_PROVIDER [double precision, fock_wee_closed_shell, (mo_tot_num, mo_tot_num) ]
  implicit none
  integer :: i0,j0,i,j,k0,k
  integer :: n_occ_ab(2)
@@ -92,8 +92,8 @@ BEGIN_PROVIDER [double precision, fock_operator_bielec_closed_shell_ref_bitmask,
     k = occ(k0,1)
     accu += 2.d0 * array_coulomb(k) - array_exchange(k)
    enddo
-   fock_operator_bielec_closed_shell_ref_bitmask(i,j) = accu 
+   fock_wee_closed_shell(i,j) = accu 
-   fock_operator_bielec_closed_shell_ref_bitmask(j,i) = accu 
+   fock_wee_closed_shell(j,i) = accu 
   enddo
  enddo
 
@@ -109,8 +109,8 @@ BEGIN_PROVIDER [double precision, fock_operator_bielec_closed_shell_ref_bitmask,
     k = occ(k0,1)
     accu += 2.d0 * array_coulomb(k) - array_exchange(k)
    enddo
-   fock_operator_bielec_closed_shell_ref_bitmask(i,j) = accu 
+   fock_wee_closed_shell(i,j) = accu 
-   fock_operator_bielec_closed_shell_ref_bitmask(j,i) = accu 
+   fock_wee_closed_shell(j,i) = accu 
   enddo
  enddo
 
@@ -126,8 +126,8 @@ BEGIN_PROVIDER [double precision, fock_operator_bielec_closed_shell_ref_bitmask,
     k = occ(k0,1)
     accu += 2.d0 * array_coulomb(k) - array_exchange(k)
    enddo
-   fock_operator_bielec_closed_shell_ref_bitmask(i,j) = accu 
+   fock_wee_closed_shell(i,j) = accu 
-   fock_operator_bielec_closed_shell_ref_bitmask(j,i) = accu 
+   fock_wee_closed_shell(j,i) = accu 
   enddo
  enddo
 

src/determinants/slater_rules_wee_mono.irp.f

@@ -1,5 +1,5 @@
 
-subroutine i_H_j_mono_spin_bielec(key_i,key_j,Nint,spin,hij)
+subroutine i_Wee_j_mono(key_i,key_j,Nint,spin,hij)
   use bitmasks
   implicit none
   BEGIN_DOC
@@ -15,11 +15,11 @@ subroutine i_H_j_mono_spin_bielec(key_i,key_j,Nint,spin,hij)
   PROVIDE big_array_exchange_integrals mo_bielec_integrals_in_map
 
   call get_mono_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
-  call get_mono_excitation_from_fock_bielec(key_i,key_j,exc(1,1),exc(1,2),spin,phase,hij)
+  call mono_excitation_wee(key_i,key_j,exc(1,1),exc(1,2),spin,phase,hij)
 end
 
 
-double precision function diag_H_mat_elem_bielec(det_in,Nint)
+double precision function diag_wee_mat_elem(det_in,Nint)
   implicit none
   BEGIN_DOC
   ! Computes <i|H|i>
@@ -52,7 +52,7 @@ double precision function diag_H_mat_elem_bielec(det_in,Nint)
     nexc(2)       = nexc(2) + popcnt(hole(i,2))
   enddo
   
-  diag_H_mat_elem_bielec = bi_elec_ref_bitmask_energy
+  diag_wee_mat_elem = bi_elec_ref_bitmask_energy
   if (nexc(1)+nexc(2) == 0) then
     return
   endif
@@ -74,9 +74,9 @@ double precision function diag_H_mat_elem_bielec(det_in,Nint)
     nb = elec_num_tab(iand(ispin,1)+1)
     do i=1,nexc(ispin)
       !DIR$ FORCEINLINE
-      call ac_operator_bielec( occ_particle(i,ispin), ispin, det_tmp, diag_H_mat_elem_bielec, Nint,na,nb)
+      call ac_operator_bielec( occ_particle(i,ispin), ispin, det_tmp, diag_wee_mat_elem, Nint,na,nb)
       !DIR$ FORCEINLINE
-      call a_operator_bielec ( occ_hole    (i,ispin), ispin, det_tmp, diag_H_mat_elem_bielec, Nint,na,nb)
+      call a_operator_bielec ( occ_hole    (i,ispin), ispin, det_tmp, diag_wee_mat_elem, Nint,na,nb)
     enddo
   enddo
 end
@@ -352,10 +352,10 @@ subroutine i_H_j_bielec(key_i,key_j,Nint,hij)
         p = exc(1,2,2)
         spin = 2
       endif
-      call get_mono_excitation_from_fock_bielec(key_i,key_j,p,m,spin,phase,hij)
+      call mono_excitation_wee(key_i,key_j,p,m,spin,phase,hij)
     case (0)
-      double precision :: diag_H_mat_elem_bielec
+      double precision :: diag_wee_mat_elem
-      hij = diag_H_mat_elem_bielec(key_i,Nint)
+      hij = diag_wee_mat_elem(key_i,Nint)
   end select
 end
 

src/dft_utils_one_e/NEED

@@ -1,6 +1,8 @@
 density_for_dft 
-dft_utils_on_grid 
+dft_utils_in_r
 mo_one_e_integrals
 mo_two_e_integrals
 ao_one_e_integrals
 ao_two_e_integrals
+mo_two_e_erf_integrals
+ao_two_e_erf_integrals

src/dft_utils_one_e/e_xc_general.irp.f

@@ -0,0 +1,38 @@
+
+  BEGIN_PROVIDER [double precision, energy_x, (N_states)]
+ &BEGIN_PROVIDER [double precision, energy_c, (N_states)]
+ implicit none
+ BEGIN_DOC
+ ! correlation and exchange energies general providers.
+ END_DOC
+  if(trim(exchange_functional)=="short_range_LDA")then
+   energy_x = energy_sr_x_LDA
+   energy_x = energy_sr_x_LDA
+  else if(exchange_functional.EQ."short_range_PBE")then
+   energy_x = energy_sr_x_PBE
+   energy_x = energy_sr_x_PBE
+  else if(exchange_functional.EQ."None")then
+   energy_x = 0.d0 
+   energy_x = 0.d0 
+  else 
+   print*, 'Exchange functional required does not exist ...'
+   print*,'exchange_functional',exchange_functional
+   stop
+  endif
+
+  if(trim(correlation_functional)=="short_range_LDA")then
+   energy_c = energy_sr_c_LDA
+   energy_c = energy_sr_c_LDA
+  else if(correlation_functional.EQ."short_range_PBE")then
+   energy_c = energy_sr_c_PBE
+   energy_c = energy_sr_c_PBE
+  else if(correlation_functional.EQ."None")then
+   energy_c = 0.d0 
+   energy_c = 0.d0 
+  else 
+   print*, 'Correlation functional required does not ecist ...'
+   print*,'correlation_functional',correlation_functional
+   stop
+  endif
+ 
+END_PROVIDER 

src/dft_utils_one_e/exc_sr_pbe.irp.f

@@ -1,19 +1,42 @@
 subroutine ec_pbe_sr(mu,rhoc,rhoo,sigmacc,sigmaco,sigmaoo,ec,vrhoc,vrhoo,vsigmacc,vsigmaco,vsigmaoo)
-!************************************************************************
+ BEGIN_DOC 
-!     Short-range PBE correlation energy functional for erf interaction
+! Short-range PBE correlation energy functional for erf interaction
 !
+! input : ==========
+! 
+! mu = range separated parameter 
 !
-!************************************************************************
+! rhoc, rhoo = total density and spin density
+!
+! sigmacc    = square of the gradient of the total density 
+!
+! sigmaco    = square of the gradient of the spin density
+!
+! sigmaoo    = scalar product between the gradient of the total density and the one of the spin density
+!
+! output: ==========
+! 
+! ec         = correlation energy 
+!
+! all variables v** are energy derivatives with respect to components of the density 
+! 
+! vrhoc      = derivative with respect to the total density 
+!
+! vrhoo      = derivative with respect to spin density
+!
+! vsigmacc   = derivative with respect to the square of the gradient of the total density
+!
+! vsigmaco   = derivative with respect to scalar product between the gradients of total and spin densities 
+!
+! vsigmaoo   = derivative with respect to the square of the gradient of the psin density
+ END_DOC
 include 'constants.include.F'
       implicit none
-! input
       double precision, intent(in) ::  rhoc,rhoo,mu
       double precision, intent(in) ::  sigmacc,sigmaco,sigmaoo
-! output
       double precision, intent(out) ::  ec
       double precision, intent(out) ::  vrhoc,vrhoo
       double precision, intent(out) ::  vsigmacc,vsigmaco,vsigmaoo
-! local
       double precision tol
       parameter(tol=1d-12)
 
@@ -82,8 +105,6 @@ include 'constants.include.F'
      double precision :: vc_a_lda,vc_b_lda
      call ec_lda(rhoa,rhob,ecclda,vc_a_lda,vc_b_lda)
      eclda = ecclda
-    !decldadrho = 0.5d0 * (vc_a_lda+vc_b_lda)
-    !decldadrho = 0.5d0 * (vc_a_lda-vc_b_lda)
 
      if ((ecerflda/eclda).le.0d0) then
         beta=0d0
@@ -108,7 +129,6 @@ include 'constants.include.F'
 
      ec = ecerfpbe
 
-!     if(ldebug) write(*,*)"ecerfpbe=",ecerfpbe
 
 ! Derive
 
@@ -168,12 +188,13 @@ end
 
 subroutine ex_pbe_sr(mu,rho_a,rho_b,grd_rho_a_2,grd_rho_b_2,grd_rho_a_b,ex,vx_rho_a,vx_rho_b,vx_grd_rho_a_2,vx_grd_rho_b_2,vx_grd_rho_a_b)
 BEGIN_DOC
+!mu    = range separation parameter
 !rho_a = density alpha
 !rho_b = density beta
 !grd_rho_a_2 = (gradient rho_a)^2
 !grd_rho_b_2 = (gradient rho_b)^2
 !grd_rho_a_b = (gradient rho_a).(gradient rho_b)
-!ex = exchange energy density at point r
+!ex = exchange energy density at the density and corresponding gradients of the density
 !vx_rho_a = d ex / d rho_a
 !vx_rho_b = d ex / d rho_b
 !vx_grd_rho_a_2 = d ex / d grd_rho_a_2
@@ -374,11 +395,26 @@ END_DOC
 
 
 subroutine ec_pbe_only(mu,rhoc,rhoo,sigmacc,sigmaco,sigmaoo,ec)
-!************************************************************************
+ BEGIN_DOC 
-!     Short-range PBE correlation energy functional for erf interaction
+! Short-range PBE correlation energy functional for erf interaction
 !
+! input : ==========
+! 
+! mu = range separated parameter 
 !
-!************************************************************************
+! rhoc, rhoo = total density and spin density
+!
+! sigmacc    = square of the gradient of the total density 
+!
+! sigmaco    = square of the gradient of the spin density
+!
+! sigmaoo    = scalar product between the gradient of the total density and the one of the spin density
+!
+! output: ==========
+! 
+! ec         = correlation energy 
+!
+ END_DOC
 include 'constants.include.F'
       implicit none
 ! input

src/dft_utils_one_e/pot_general.irp.f

@@ -88,46 +88,6 @@ END_PROVIDER
 
 END_PROVIDER 
 
-
-  BEGIN_PROVIDER [double precision, energy_x, (N_states)]
- &BEGIN_PROVIDER [double precision, energy_c, (N_states)]
- implicit none
- BEGIN_DOC
- ! correlation and exchange energies general providers.
- END_DOC
-  if(trim(exchange_functional)=="short_range_LDA")then
-   energy_x = energy_sr_x_LDA
-   energy_x = energy_sr_x_LDA
-  else if(exchange_functional.EQ."short_range_PBE")then
-   energy_x = energy_sr_x_PBE
-   energy_x = energy_sr_x_PBE
-  else if(exchange_functional.EQ."None")then
-   energy_x = 0.d0 
-   energy_x = 0.d0 
-  else 
-   print*, 'Exchange functional required does not exist ...'
-   print*,'exchange_functional',exchange_functional
-   stop
-  endif
-
-  if(trim(correlation_functional)=="short_range_LDA")then
-   energy_c = energy_sr_c_LDA
-   energy_c = energy_sr_c_LDA
-  else if(correlation_functional.EQ."short_range_PBE")then
-   energy_c = energy_sr_c_PBE
-   energy_c = energy_sr_c_PBE
-  else if(correlation_functional.EQ."None")then
-   energy_c = 0.d0 
-   energy_c = 0.d0 
-  else 
-   print*, 'Correlation functional required does not ecist ...'
-   print*,'correlation_functional',correlation_functional
-   stop
-  endif
- 
-END_PROVIDER 
-
-
  BEGIN_PROVIDER [double precision, Trace_v_xc, (N_states)]
 &BEGIN_PROVIDER [double precision, Trace_v_H, (N_states)]
 &BEGIN_PROVIDER [double precision, Trace_v_Hxc, (N_states)]

src/dft_utils_two_body/NEED

@@ -1,3 +1,3 @@
-dft_utils_one_body 
+dft_utils_one_e 
 determinants
 davidson_undressed

src/dft_utils_two_body/psi_energy_erf.irp.f

@@ -61,19 +61,19 @@ END_PROVIDER
   
   do j= 1, elec_alpha_num
     do i = j+1, elec_alpha_num
-      bi_elec_ref_bitmask_energy_erf += mo_bielec_integral_erf_jj_anti(occ(i,1),occ(j,1))
+      bi_elec_ref_bitmask_energy_erf += mo_two_e_int_erf_jj_anti(occ(i,1),occ(j,1))
-      ref_bitmask_energy_erf += mo_bielec_integral_erf_jj_anti(occ(i,1),occ(j,1))
+      ref_bitmask_energy_erf += mo_two_e_int_erf_jj_anti(occ(i,1),occ(j,1))
     enddo
   enddo
   
   do j= 1, elec_beta_num
     do i = j+1, elec_beta_num
-      bi_elec_ref_bitmask_energy_erf += mo_bielec_integral_erf_jj_anti(occ(i,2),occ(j,2))
+      bi_elec_ref_bitmask_energy_erf += mo_two_e_int_erf_jj_anti(occ(i,2),occ(j,2))
-      ref_bitmask_energy_erf += mo_bielec_integral_erf_jj_anti(occ(i,2),occ(j,2))
+      ref_bitmask_energy_erf += mo_two_e_int_erf_jj_anti(occ(i,2),occ(j,2))
     enddo
     do i= 1, elec_alpha_num
-      bi_elec_ref_bitmask_energy_erf += mo_bielec_integral_erf_jj(occ(i,1),occ(j,2))
+      bi_elec_ref_bitmask_energy_erf += mo_two_e_int_erf_jj(occ(i,1),occ(j,2))
-      ref_bitmask_energy_erf += mo_bielec_integral_erf_jj(occ(i,1),occ(j,2))
+      ref_bitmask_energy_erf += mo_two_e_int_erf_jj(occ(i,1),occ(j,2))
     enddo
   enddo
   

src/dft_utils_two_body/slater_rules_erf.irp.f

@@ -20,7 +20,7 @@ subroutine i_H_j_erf(key_i,key_j,Nint,hij)
   integer                        :: occ(Nint*bit_kind_size,2)
   double precision               :: diag_H_mat_elem_erf, phase,phase_2
   integer                        :: n_occ_ab(2)
-  PROVIDE mo_bielec_integrals_erf_in_map mo_integrals_erf_map big_array_exchange_integrals_erf
+  PROVIDE mo_bielec_integrals_erf_in_map mo_integrals_erf_map int_erf_3_index_exc
   
   ASSERT (Nint > 0)
   ASSERT (Nint == N_int)
@@ -39,9 +39,9 @@ subroutine i_H_j_erf(key_i,key_j,Nint,hij)
       if (exc(0,1,1) == 1) then
         ! Mono alpha, mono beta
         if(exc(1,1,1) == exc(1,2,2) )then
-         hij = phase * big_array_exchange_integrals_erf(exc(1,1,1),exc(1,1,2),exc(1,2,1))
+         hij = phase * int_erf_3_index_exc(exc(1,1,1),exc(1,1,2),exc(1,2,1))
         else if (exc(1,2,1) ==exc(1,1,2))then
-         hij = phase * big_array_exchange_integrals_erf(exc(1,2,1),exc(1,1,1),exc(1,2,2))
+         hij = phase * int_erf_3_index_exc(exc(1,2,1),exc(1,1,1),exc(1,2,2))
         else
          hij = phase*get_mo_bielec_integral_erf(                          &
              exc(1,1,1),                                              &
@@ -84,10 +84,10 @@ subroutine i_H_j_erf(key_i,key_j,Nint,hij)
         p = exc(1,2,1)
         spin = 1
         do i = 1, n_occ_ab(1)
-         hij += -big_array_exchange_integrals_erf(occ(i,1),m,p) + big_array_coulomb_integrals_erf(occ(i,1),m,p)
+         hij += -int_erf_3_index_exc(occ(i,1),m,p) + int_erf_3_index(occ(i,1),m,p)
         enddo
         do i = 1, n_occ_ab(2)
-         hij += big_array_coulomb_integrals_erf(occ(i,2),m,p)
+         hij += int_erf_3_index(occ(i,2),m,p)
         enddo
       else
         ! Mono beta
@@ -95,10 +95,10 @@ subroutine i_H_j_erf(key_i,key_j,Nint,hij)
         p = exc(1,2,2)
         spin = 2
         do i = 1, n_occ_ab(2)
-         hij += -big_array_exchange_integrals_erf(occ(i,2),m,p) + big_array_coulomb_integrals_erf(occ(i,2),m,p)
+         hij += -int_erf_3_index_exc(occ(i,2),m,p) + int_erf_3_index(occ(i,2),m,p)
         enddo
         do i = 1, n_occ_ab(1)
-         hij += big_array_coulomb_integrals_erf(occ(i,1),m,p)
+         hij += int_erf_3_index(occ(i,1),m,p)
         enddo
       endif
       hij = hij * phase
@@ -124,21 +124,21 @@ double precision function diag_H_mat_elem_erf(key_i,Nint)
  ! alpha - alpha
  do i = 1, n_occ_ab(1)
   do j = i+1, n_occ_ab(1)
-   diag_H_mat_elem_erf += mo_bielec_integral_erf_jj_anti(occ(i,1),occ(j,1))
+   diag_H_mat_elem_erf += mo_two_e_int_erf_jj_anti(occ(i,1),occ(j,1))
   enddo
  enddo
 
  ! beta - beta 
  do i = 1, n_occ_ab(2)
   do j = i+1, n_occ_ab(2)
-   diag_H_mat_elem_erf += mo_bielec_integral_erf_jj_anti(occ(i,2),occ(j,2))
+   diag_H_mat_elem_erf += mo_two_e_int_erf_jj_anti(occ(i,2),occ(j,2))
   enddo
  enddo
 
  ! alpha - beta 
  do i = 1, n_occ_ab(1)
   do j = 1, n_occ_ab(2)
-   diag_H_mat_elem_erf += mo_bielec_integral_erf_jj(occ(i,1),occ(j,2))
+   diag_H_mat_elem_erf += mo_two_e_int_erf_jj(occ(i,1),occ(j,2))
   enddo
  enddo
 end
@@ -155,7 +155,7 @@ subroutine i_H_j_mono_spin_erf(key_i,key_j,Nint,spin,hij)
   integer                        :: exc(0:2,2)
   double precision               :: phase
 
-  PROVIDE big_array_exchange_integrals_erf mo_bielec_integrals_erf_in_map
+  PROVIDE int_erf_3_index_exc mo_bielec_integrals_erf_in_map
 
   call i_H_j_erf(key_i,key_j,Nint,hij)
 end
@@ -176,7 +176,7 @@ subroutine i_H_j_double_spin_erf(key_i,key_j,Nint,hij)
   double precision               :: phase
   double precision, external     :: get_mo_bielec_integral_erf
 
-  PROVIDE big_array_exchange_integrals_erf mo_bielec_integrals_erf_in_map
+  PROVIDE int_erf_3_index_exc mo_bielec_integrals_erf_in_map
 
   call get_double_excitation_spin(key_i,key_j,exc,phase,Nint)
   hij = phase*(get_mo_bielec_integral_erf(                               &
@@ -205,15 +205,15 @@ subroutine i_H_j_double_alpha_beta_erf(key_i,key_j,Nint,hij)
   double precision               :: phase, phase2
   double precision, external     :: get_mo_bielec_integral_erf
 
-  PROVIDE big_array_exchange_integrals_erf mo_bielec_integrals_erf_in_map
+  PROVIDE int_erf_3_index_exc mo_bielec_integrals_erf_in_map
 
   call get_mono_excitation_spin(key_i(1,1),key_j(1,1),exc(0,1,1),phase,Nint)
   call get_mono_excitation_spin(key_i(1,2),key_j(1,2),exc(0,1,2),phase2,Nint)
   phase = phase*phase2
   if (exc(1,1,1) == exc(1,2,2)) then
-    hij = phase * big_array_exchange_integrals_erf(exc(1,1,1),exc(1,1,2),exc(1,2,1))
+    hij = phase * int_erf_3_index_exc(exc(1,1,1),exc(1,1,2),exc(1,2,1))
   else if (exc(1,2,1) == exc(1,1,2)) then
-    hij = phase * big_array_exchange_integrals_erf(exc(1,2,1),exc(1,1,1),exc(1,2,2))
+    hij = phase * int_erf_3_index_exc(exc(1,2,1),exc(1,1,1),exc(1,2,2))
   else
     hij = phase*get_mo_bielec_integral_erf(                              &
         exc(1,1,1),                                                  &

src/dft_utils_two_body/write_effective_rsdft_hamiltonian.irp.f

@@ -26,8 +26,8 @@ end
 subroutine routines_compute_energy
  implicit none
  call print_variational_energy_dft
- call ezfio_set_data_energy_and_density_data_one_body_alpha_dm_mo(one_body_dm_mo_alpha)
+ call ezfio_set_aux_quantities_data_one_body_alpha_dm_mo(one_body_dm_mo_alpha)
- call ezfio_set_data_energy_and_density_data_one_body_beta_dm_mo(one_body_dm_mo_beta)
+ call ezfio_set_aux_quantities_data_one_body_beta_dm_mo(one_body_dm_mo_beta)
 
 end
 

src/dummy/NEED

@@ -1,6 +1,8 @@
 ao_basis
 ao_one_e_integrals
+ao_two_e_erf_integrals
 ao_two_e_integrals
+aux_quantities
 becke_numerical_grid
 bitmask
 cis
@@ -8,8 +10,12 @@ cisd
 davidson
 davidson_dressed
 davidson_undressed
+density_for_dft
 determinants
-dft_utils_one_body
+dft_keywords
+dft_utils_in_r
+dft_utils_one_e
+dft_utils_two_body
 dressing
 electrons
 ezfio_files
@@ -17,9 +23,13 @@ fci
 generators_cas
 generators_full
 hartree_fock
+iterations
+kohn_sham
+kohn_sham_rs
 mo_basis
 mo_guess
 mo_one_e_integrals
+mo_two_e_erf_integrals
 mo_two_e_integrals
 mpi
 mrpt_utils
@@ -28,9 +38,12 @@ perturbation
 pseudo
 psiref_cas
 psiref_utils
+scf_utils
 selectors_cassd
 selectors_full
 selectors_utils
 single_ref_method
+slave
+tools
 utils
 zmq

src/hartree_fock/README.rst

@@ -0,0 +1,33 @@
+============
+Hartree-Fock
+============
+
+
+The Hartree-Fock module performs *Restricted* Hartree-Fock calculations (the
+spatial part of the |MOs| is common for alpha and beta spinorbitals).
+
+The Hartree-Fock program does the following:
+
+#. Compute/Read all the one- and two-electron integrals, and store them in memory
+#. Check in the |EZFIO| database if there is a set of |MOs|. If there is, it
+   will read them as initial guess. Otherwise, it will create a guess.
+#. Perform the |SCF| iterations
+
+At each iteration, the |MOs| are saved in the |EZFIO| database. Hence, if the calculation
+crashes for any unexpected reason, the calculation can be restarted by running again
+the |SCF| with the same |EZFIO| database.
+
+The `DIIS`_ algorithm is implemented, as well as the `level-shifting`_ method.
+If the |SCF| does not converge, try again with a higher value of :option:`level_shift`.
+
+To start a calculation from scratch, the simplest way is to remove the
+``mo_basis`` directory from the |EZFIO| database, and run the |SCF| again.
+
+
+
+
+.. _DIIS: https://en.wikipedia.org/w/index.php?title=DIIS
+.. _level-shifting: https://doi.org/10.1002/qua.560070407
+
+
+

src/kohn_sham/NEED

@@ -1,2 +1,2 @@
-dft_utils_one_body 
+dft_utils_one_e 
 scf_utils

src/kohn_sham/README.rst

@@ -0,0 +1,33 @@
+============
+Hartree-Fock
+============
+
+
+The Hartree-Fock module performs *Restricted* Hartree-Fock calculations (the
+spatial part of the |MOs| is common for alpha and beta spinorbitals).
+
+The Hartree-Fock program does the following:
+
+#. Compute/Read all the one- and two-electron integrals, and store them in memory
+#. Check in the |EZFIO| database if there is a set of |MOs|. If there is, it
+   will read them as initial guess. Otherwise, it will create a guess.
+#. Perform the |SCF| iterations
+
+At each iteration, the |MOs| are saved in the |EZFIO| database. Hence, if the calculation
+crashes for any unexpected reason, the calculation can be restarted by running again
+the |SCF| with the same |EZFIO| database.
+
+The `DIIS`_ algorithm is implemented, as well as the `level-shifting`_ method.
+If the |SCF| does not converge, try again with a higher value of :option:`level_shift`.
+
+To start a calculation from scratch, the simplest way is to remove the
+``mo_basis`` directory from the |EZFIO| database, and run the |SCF| again.
+
+
+
+
+.. _DIIS: https://en.wikipedia.org/w/index.php?title=DIIS
+.. _level-shifting: https://doi.org/10.1002/qua.560070407
+
+
+

src/kohn_sham_range_separated/NEED

@@ -1,2 +0,0 @@
-dft_utils_one_body 
-scf_utils 

src/kohn_sham_rs/NEED

@@ -0,0 +1,2 @@
+dft_utils_one_e 
+scf_utils 

src/kohn_sham_rs/README.rst

@@ -0,0 +1,33 @@
+============
+Hartree-Fock
+============
+
+
+The Hartree-Fock module performs *Restricted* Hartree-Fock calculations (the
+spatial part of the |MOs| is common for alpha and beta spinorbitals).
+
+The Hartree-Fock program does the following:
+
+#. Compute/Read all the one- and two-electron integrals, and store them in memory
+#. Check in the |EZFIO| database if there is a set of |MOs|. If there is, it
+   will read them as initial guess. Otherwise, it will create a guess.
+#. Perform the |SCF| iterations
+
+At each iteration, the |MOs| are saved in the |EZFIO| database. Hence, if the calculation
+crashes for any unexpected reason, the calculation can be restarted by running again
+the |SCF| with the same |EZFIO| database.
+
+The `DIIS`_ algorithm is implemented, as well as the `level-shifting`_ method.
+If the |SCF| does not converge, try again with a higher value of :option:`level_shift`.
+
+To start a calculation from scratch, the simplest way is to remove the
+``mo_basis`` directory from the |EZFIO| database, and run the |SCF| again.
+
+
+
+
+.. _DIIS: https://en.wikipedia.org/w/index.php?title=DIIS
+.. _level-shifting: https://doi.org/10.1002/qua.560070407
+
+
+

src/mo_two_e_erf_integrals/EZFIO.cfg

@@ -0,0 +1,6 @@
+[disk_access_mo_integrals_erf]
+type: Disk_access
+doc: Read/Write MO integrals with the long range interaction from/to disk [    Write | Read | None ]
+interface: ezfio,provider,ocaml
+default: None
+

src/mo_two_e_erf_integrals/NEED

@@ -0,0 +1,3 @@
+ao_two_e_erf_integrals
+mo_two_e_integrals
+mo_basis

src/mo_two_e_erf_integrals/README.rst

@@ -0,0 +1,20 @@
+======================
+mo_two_e_erf_integrals
+======================
+
+Here, all two-electron integrals (:math:`erf({\mu}_{erf} * r_{12})/r_{12}`) are computed.
+As they have 4 indices and many are zero, they are stored in a map, as defined
+in :file:`Utils/map_module.f90`. 
+
+The range separation parameter :math:`{\mu}_{erf}` is the variable :option:`ao_two_e_erf_integrals mu_erf`. 
+
+To fetch an |MO| integral, use
+`get_mo_bielec_integral_erf(i,j,k,l,mo_integrals_map_erf)`
+
+The conventions are: 
+
+* For |MO| integrals : <ij|kl> = <12|12>
+
+Be aware that it might not be the same conventions for |MO| and |AO| integrals. 
+
+

src/mo_two_e_erf_integrals/core_quantities_erf.irp.f

@@ -0,0 +1,38 @@
+BEGIN_PROVIDER [double precision, core_energy_erf]
+ implicit none
+ BEGIN_DOC 
+! energy from the core : contains all core-core contributionswith the erf interaction
+ END_DOC
+ integer :: i,j,k,l
+ core_energy_erf = 0.d0
+ do i = 1, n_core_orb
+  j = list_core(i)
+  core_energy_erf += 2.d0 * mo_mono_elec_integral(j,j) + mo_two_e_int_erf_jj(j,j)
+  do k = i+1, n_core_orb
+   l = list_core(k)
+   core_energy_erf += 2.d0 * (2.d0 * mo_two_e_int_erf_jj(j,l) - mo_two_e_int_erf_jj_exchange(j,l))
+  enddo
+ enddo
+ core_energy_erf += nuclear_repulsion
+
+END_PROVIDER
+
+BEGIN_PROVIDER [double precision, core_fock_operator_erf, (mo_tot_num,mo_tot_num)]
+ implicit none
+ integer :: i,j,k,l,m,n
+ double precision :: get_mo_bielec_integral_erf
+ BEGIN_DOC 
+! this is the contribution to the Fock operator from the core electrons with the erf interaction
+ END_DOC
+ core_fock_operator_erf = 0.d0
+ do i = 1, n_act_orb
+  j = list_act(i)
+  do k = 1, n_act_orb
+   l = list_act(k)
+   do m = 1, n_core_orb
+    n = list_core(m)
+    core_fock_operator_erf(j,l) += 2.d0 * get_mo_bielec_integral_erf(j,n,l,n,mo_integrals_erf_map) - get_mo_bielec_integral_erf(j,n,n,l,mo_integrals_erf_map)
+   enddo
+  enddo
+ enddo
+END_PROVIDER

src/mo_two_e_erf_integrals/ints_erf_3_index.irp.f

@@ -0,0 +1,28 @@
+ BEGIN_PROVIDER [double precision, int_erf_3_index, (mo_tot_num,mo_tot_num, mo_tot_num)]
+&BEGIN_PROVIDER [double precision, int_erf_3_index_exc,(mo_tot_num,mo_tot_num, mo_tot_num)]
+ implicit none
+ BEGIN_DOC
+ ! int_erf_3_index(i,j)     = <ij|ij> = (ii|jj) with the erf interaction
+ !
+ ! int_erf_3_index_exc(i,j) = <ij|ji> = (ij|ij) with the erf interaction
+ END_DOC
+ integer :: i,j,k,l
+ double precision :: get_mo_bielec_integral_erf
+ double precision :: integral
+
+ do k = 1, mo_tot_num
+  do i = 1, mo_tot_num
+   do j = 1, mo_tot_num
+     l = j
+     integral = get_mo_bielec_integral_erf(i,j,k,l,mo_integrals_erf_map)
+     int_erf_3_index(j,i,k) = integral
+     l = j
+     integral = get_mo_bielec_integral_erf(i,j,l,k,mo_integrals_erf_map)
+     int_erf_3_index_exc(j,i,k) = integral
+   enddo
+  enddo
+ enddo
+
+
+END_PROVIDER 
+

src/mo_two_e_erf_integrals/mo_bi_integrals_erf.irp.f

@@ -55,15 +55,15 @@ BEGIN_PROVIDER [ logical, mo_bielec_integrals_erf_in_map ]
   if (write_mo_integrals_erf) then
     call ezfio_set_work_empty(.False.)
     call map_save_to_disk(trim(ezfio_filename)//'/work/mo_ints_erf',mo_integrals_erf_map)
-    call ezfio_set_mo_two_e_integrals_disk_access_mo_integrals_erf("Read")
+    call ezfio_set_mo_two_e_erf_integrals_disk_access_mo_integrals_erf("Read")
   endif
   
 END_PROVIDER
 
 
- BEGIN_PROVIDER [ double precision, mo_bielec_integral_erf_jj_from_ao, (mo_tot_num,mo_tot_num) ]
+ BEGIN_PROVIDER [ double precision, mo_two_e_int_erf_jj_from_ao, (mo_tot_num,mo_tot_num) ]
-&BEGIN_PROVIDER [ double precision, mo_bielec_integral_erf_jj_exchange_from_ao, (mo_tot_num,mo_tot_num) ]
+&BEGIN_PROVIDER [ double precision, mo_two_e_int_erf_jj_exchange_from_ao, (mo_tot_num,mo_tot_num) ]
-&BEGIN_PROVIDER [ double precision, mo_bielec_integral_erf_jj_anti_from_ao, (mo_tot_num,mo_tot_num) ]
+&BEGIN_PROVIDER [ double precision, mo_two_e_int_erf_jj_anti_from_ao, (mo_tot_num,mo_tot_num) ]
   BEGIN_DOC
   ! mo_bielec_integral_jj_from_ao(i,j) = J_ij
   ! mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij
@@ -83,8 +83,8 @@ END_PROVIDER
     PROVIDE ao_bielec_integrals_erf_in_map mo_coef
   endif
   
-  mo_bielec_integral_erf_jj_from_ao = 0.d0
+  mo_two_e_int_erf_jj_from_ao = 0.d0
-  mo_bielec_integral_erf_jj_exchange_from_ao = 0.d0
+  mo_two_e_int_erf_jj_exchange_from_ao = 0.d0
   
   !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: iqrs, iqsr
   
@@ -94,7 +94,7 @@ END_PROVIDER
       !$OMP  iqrs, iqsr,iqri,iqis)                                   &
       !$OMP SHARED(mo_tot_num,mo_coef_transp,ao_num,&
       !$OMP  ao_integrals_threshold,do_direct_integrals)             &
-      !$OMP REDUCTION(+:mo_bielec_integral_erf_jj_from_ao,mo_bielec_integral_erf_jj_exchange_from_ao)
+      !$OMP REDUCTION(+:mo_two_e_int_erf_jj_from_ao,mo_two_e_int_erf_jj_exchange_from_ao)
   
   allocate( int_value(ao_num), int_idx(ao_num),                      &
       iqrs(mo_tot_num,ao_num), iqis(mo_tot_num), iqri(mo_tot_num),&
@@ -177,8 +177,8 @@ END_PROVIDER
         !DIR$ VECTOR ALIGNED
         do j=1,mo_tot_num
           c = mo_coef_transp(j,q)*mo_coef_transp(j,s)
-          mo_bielec_integral_erf_jj_from_ao(j,i) += c * iqis(i)
+          mo_two_e_int_erf_jj_from_ao(j,i) += c * iqis(i)
-          mo_bielec_integral_erf_jj_exchange_from_ao(j,i) += c * iqri(i)
+          mo_two_e_int_erf_jj_exchange_from_ao(j,i) += c * iqri(i)
         enddo
       enddo
       
@@ -188,16 +188,16 @@ END_PROVIDER
   deallocate(iqrs,iqsr,int_value,int_idx)
   !$OMP END PARALLEL
   
-  mo_bielec_integral_erf_jj_anti_from_ao = mo_bielec_integral_erf_jj_from_ao - mo_bielec_integral_erf_jj_exchange_from_ao
+  mo_two_e_int_erf_jj_anti_from_ao = mo_two_e_int_erf_jj_from_ao - mo_two_e_int_erf_jj_exchange_from_ao
   
   
 ! end
 END_PROVIDER 
 
 
- BEGIN_PROVIDER [ double precision, mo_bielec_integral_erf_jj, (mo_tot_num,mo_tot_num) ]
+ BEGIN_PROVIDER [ double precision, mo_two_e_int_erf_jj, (mo_tot_num,mo_tot_num) ]
-&BEGIN_PROVIDER [ double precision, mo_bielec_integral_erf_jj_exchange, (mo_tot_num,mo_tot_num) ]
+&BEGIN_PROVIDER [ double precision, mo_two_e_int_erf_jj_exchange, (mo_tot_num,mo_tot_num) ]
-&BEGIN_PROVIDER [ double precision, mo_bielec_integral_erf_jj_anti, (mo_tot_num,mo_tot_num) ]
+&BEGIN_PROVIDER [ double precision, mo_two_e_int_erf_jj_anti, (mo_tot_num,mo_tot_num) ]
   implicit none
   BEGIN_DOC
   ! mo_bielec_integral_jj(i,j) = J_ij
@@ -209,14 +209,14 @@ END_PROVIDER
   double precision               :: get_mo_bielec_integral_erf
   
   PROVIDE mo_bielec_integrals_erf_in_map
-  mo_bielec_integral_erf_jj = 0.d0
+  mo_two_e_int_erf_jj = 0.d0
-  mo_bielec_integral_erf_jj_exchange = 0.d0
+  mo_two_e_int_erf_jj_exchange = 0.d0
   
   do j=1,mo_tot_num
     do i=1,mo_tot_num
-      mo_bielec_integral_erf_jj(i,j) = get_mo_bielec_integral_erf(i,j,i,j,mo_integrals_erf_map)
+      mo_two_e_int_erf_jj(i,j) = get_mo_bielec_integral_erf(i,j,i,j,mo_integrals_erf_map)
-      mo_bielec_integral_erf_jj_exchange(i,j) = get_mo_bielec_integral_erf(i,j,j,i,mo_integrals_erf_map)
+      mo_two_e_int_erf_jj_exchange(i,j) = get_mo_bielec_integral_erf(i,j,j,i,mo_integrals_erf_map)
-      mo_bielec_integral_erf_jj_anti(i,j) = mo_bielec_integral_erf_jj(i,j) - mo_bielec_integral_erf_jj_exchange(i,j)
+      mo_two_e_int_erf_jj_anti(i,j) = mo_two_e_int_erf_jj(i,j) - mo_two_e_int_erf_jj_exchange(i,j)
     enddo
   enddo
   
@@ -229,16 +229,16 @@ subroutine clear_mo_erf_map
   ! Frees the memory of the MO map
   END_DOC
   call map_deinit(mo_integrals_erf_map)
-  FREE mo_integrals_erf_map mo_bielec_integral_erf_jj mo_bielec_integral_erf_jj_anti
+  FREE mo_integrals_erf_map mo_two_e_int_erf_jj mo_two_e_int_erf_jj_anti
-  FREE mo_bielec_integral_Erf_jj_exchange mo_bielec_integrals_erf_in_map
+  FREE mo_two_e_int_erf_jj_exchange mo_bielec_integrals_erf_in_map
   
   
 end
 
 subroutine provide_all_mo_integrals_erf
   implicit none
-  provide mo_integrals_erf_map mo_bielec_integral_erf_jj mo_bielec_integral_erf_jj_anti
+  provide mo_integrals_erf_map mo_two_e_int_erf_jj mo_two_e_int_erf_jj_anti
-  provide mo_bielec_integral_erf_jj_exchange mo_bielec_integrals_erf_in_map
+  provide mo_two_e_int_erf_jj_exchange mo_bielec_integrals_erf_in_map
   
 end
 

src/mo_two_e_erf_integrals/read_write_erf.irp.f

@@ -0,0 +1,27 @@
+ BEGIN_PROVIDER [ logical, read_mo_integrals_erf ]
+&BEGIN_PROVIDER [ logical, write_mo_integrals_erf ]
+   implicit none
+   
+   BEGIN_DOC
+   ! Flag to read or write the |MO| erf integrals
+   END_DOC
+   
+   if (disk_access_mo_integrals_erf.EQ.'Read') then
+     read_mo_integrals_erf =  .True.
+     write_mo_integrals_erf = .False.
+     
+   else if  (disk_access_mo_integrals_erf.EQ.'Write') then
+     read_mo_integrals_erf = .False.
+     write_mo_integrals_erf =  .True.
+     
+   else if (disk_access_mo_integrals_erf.EQ.'None') then
+     read_mo_integrals_erf = .False.
+     write_mo_integrals_erf = .False.
+     
+   else
+     print *, 'disk_access_mo_integrals_erf has a wrong type'
+     stop 1
+     
+   endif
+   
+END_PROVIDER

src/mo_two_e_erf_integrals/routines_save_integrals_erf.irp.f

@@ -4,7 +4,7 @@ subroutine save_erf_bi_elec_integrals_mo
  PROVIDE mo_bielec_integrals_erf_in_map
  call ezfio_set_work_empty(.False.)
  call map_save_to_disk(trim(ezfio_filename)//'/work/mo_ints_erf',mo_integrals_erf_map)
- call ezfio_set_mo_two_e_integrals_disk_access_mo_integrals('Read')
+ call ezfio_set_mo_two_e_erf_integrals_disk_access_mo_integrals_erf('Read')
 end
 
 subroutine save_erf_bielec_ints_mo_into_ints_mo

src/mo_two_e_integrals/EZFIO.cfg

@@ -4,12 +4,6 @@ doc: Read/Write |MO| integrals from/to disk [ Write | Read | None ]
 interface: ezfio,provider,ocaml
 default: None
 
-[disk_access_mo_integrals_erf]
-type: Disk_access
-doc: Read/Write MO integrals with the long range interaction from/to disk [    Write | Read | None ]
-interface: ezfio,provider,ocaml
-default: None
-
 [mo_integrals_threshold]
 type: Threshold
 doc: If | <ij|kl> | < `mo_integrals_threshold` then <ij|kl> is zero

src/mo_two_e_integrals/core_quantities.irp.f

@@ -1,5 +1,8 @@
 BEGIN_PROVIDER [double precision, core_energy]
  implicit none
+ BEGIN_DOC 
+! energy from the core : contains all core-core contributions
+ END_DOC
  integer :: i,j,k,l
  core_energy = 0.d0
  do i = 1, n_core_orb

src/mo_two_e_integrals/integrals_3_index.irp.f

@@ -1,6 +1,11 @@
  BEGIN_PROVIDER [double precision, big_array_coulomb_integrals, (mo_tot_num,mo_tot_num, mo_tot_num)]
 &BEGIN_PROVIDER [double precision, big_array_exchange_integrals,(mo_tot_num,mo_tot_num, mo_tot_num)]
  implicit none
+ BEGIN_DOC
+ ! big_array_coulomb_integrals(i,j)  = <ij|ij> = (ii|jj) 
+ !
+ ! big_array_exchange_integrals(i,j) = <ij|ji> = (ij|ij) 
+ END_DOC
  integer :: i,j,k,l
  double precision :: get_mo_bielec_integral
  double precision :: integral
@@ -20,27 +25,3 @@
 
 END_PROVIDER 
 
-
- BEGIN_PROVIDER [double precision, big_array_coulomb_integrals_erf, (mo_tot_num,mo_tot_num, mo_tot_num)]
-&BEGIN_PROVIDER [double precision, big_array_exchange_integrals_erf,(mo_tot_num,mo_tot_num, mo_tot_num)]
- implicit none
- integer :: i,j,k,l
- double precision :: get_mo_bielec_integral_erf
- double precision :: integral
-
- do k = 1, mo_tot_num
-  do i = 1, mo_tot_num
-   do j = 1, mo_tot_num
-     l = j
-     integral = get_mo_bielec_integral_erf(i,j,k,l,mo_integrals_erf_map)
-     big_array_coulomb_integrals_erf(j,i,k) = integral
-     l = j
-     integral = get_mo_bielec_integral_erf(i,j,l,k,mo_integrals_erf_map)
-     big_array_exchange_integrals_erf(j,i,k) = integral
-   enddo
-  enddo
- enddo
-
-
-END_PROVIDER 
-

src/mo_two_e_integrals/mo_bi_integrals.irp.f

@@ -152,7 +152,7 @@ BEGIN_PROVIDER [ logical, mo_bielec_integrals_in_map ]
   if (write_mo_integrals.and.mpi_master) then
     call ezfio_set_work_empty(.False.)
     call map_save_to_disk(trim(ezfio_filename)//'/work/mo_ints',mo_integrals_map)
-    call ezfio_set_mo_two_e_integrals_disk_access_mo_integrals_erf('Read')
+    call ezfio_set_mo_two_e_integrals_disk_access_mo_integrals('Read')
   endif
   
 END_PROVIDER

src/mo_two_e_integrals/read_write.irp.f

@@ -24,32 +24,4 @@
      
    endif
    
-END_PROVIDER
- 
- BEGIN_PROVIDER [ logical, read_mo_integrals_erf ]
-&BEGIN_PROVIDER [ logical, write_mo_integrals_erf ]
-   implicit none
-   
-   BEGIN_DOC
-   ! Flag to read or write the |MO| erf integrals
-   END_DOC
-   
-   if (disk_access_mo_integrals_erf.EQ.'Read') then
-     read_mo_integrals_erf =  .True.
-     write_mo_integrals_erf = .False.
-     
-   else if  (disk_access_mo_integrals_erf.EQ.'Write') then
-     read_mo_integrals_erf = .False.
-     write_mo_integrals_erf =  .True.
-     
-   else if (disk_access_mo_integrals_erf.EQ.'None') then
-     read_mo_integrals_erf = .False.
-     write_mo_integrals_erf = .False.
-     
-   else
-     print *, 'disk_access_mo_integrals_erf has a wrong type'
-     stop 1
-     
-   endif
-   
 END_PROVIDER

src/tools/NEED

@@ -1 +1,2 @@
 fci
+mo_two_e_erf_integrals