From 85862811ad398ad5e51337e75a12332df694870f Mon Sep 17 00:00:00 2001
From: Anthony Scemama <scemama@irsamc.ups-tlse.fr>
Date: Fri, 14 Jul 2017 01:24:51 +0000
Subject: [PATCH] Truncation in QMC and 4-idx

---
 plugins/QMC/EZFIO.cfg                         |  6 ++++++
 plugins/QMC/truncate_wf_spin.irp.f            | 12 ++++++++----
 src/Integrals_Bielec/four_idx_transform.irp.f | 12 ++++++++++++
 3 files changed, 26 insertions(+), 4 deletions(-)
 create mode 100644 plugins/QMC/EZFIO.cfg
 create mode 100644 src/Integrals_Bielec/four_idx_transform.irp.f

diff --git a/plugins/QMC/EZFIO.cfg b/plugins/QMC/EZFIO.cfg
new file mode 100644
index 00000000..8960a0fa
--- /dev/null
+++ b/plugins/QMC/EZFIO.cfg
@@ -0,0 +1,6 @@
+[ci_threshold]
+type: Threshold
+doc: Threshold on the CI coefficients as in QMCChem
+interface: ezfio,provider,ocaml
+default: 1.e-8
+
diff --git a/plugins/QMC/truncate_wf_spin.irp.f b/plugins/QMC/truncate_wf_spin.irp.f
index e8cde79e..12de19ca 100644
--- a/plugins/QMC/truncate_wf_spin.irp.f
+++ b/plugins/QMC/truncate_wf_spin.irp.f
@@ -15,15 +15,13 @@ program e_curve
   nab = n_det_alpha_unique+n_det_beta_unique
   allocate ( norm_sort(0:nab), iorder(0:nab) )
 
-  double precision :: thresh
   integer(bit_kind), allocatable :: det_i(:,:), det_j(:,:)
   double precision, allocatable  :: u_t(:,:), v_t(:,:), s_t(:,:)
   double precision, allocatable  :: u_0(:,:), v_0(:,:)
   allocate(u_t(N_states,N_det),v_t(N_states,N_det),s_t(N_states,N_det))
   allocate(u_0(N_det,N_states),v_0(N_det,N_states))
 
-  print *,  'Threshold?'
-  read(*,*) thresh
+  read(*,*) ci_threshold
 
   norm_sort(0) = 0.d0
   iorder(0) = 0
@@ -45,19 +43,23 @@ program e_curve
   do j=0,nab
     i = iorder(j)
     if (i<0) then
+      !$OMP PARALLEL DO PRIVATE(k)
       do k=1,n_det
         if (psi_bilinear_matrix_columns(k) == -i) then
           psi_bilinear_matrix_values(k,1) = 0.d0
         endif
       enddo
+      !$OMP END PARALLEL DO
     else
+      !$OMP PARALLEL DO PRIVATE(k)
       do k=1,n_det
         if (psi_bilinear_matrix_rows(k) ==  i) then
           psi_bilinear_matrix_values(k,1) = 0.d0
         endif
       enddo
+      !$OMP END PARALLEL DO
     endif
-    if (thresh > norm_sort(j)) then
+    if (ci_threshold > norm_sort(j)) then
       cycle
     endif
 
@@ -70,7 +72,9 @@ program e_curve
         u_t,                                                           &
         size(u_t, 1),                                                  &
         N_det, N_states)
+    print *, 'Computing H|Psi> ...'
     call H_S2_u_0_nstates_openmp_work(v_t,s_t,u_t,N_states,N_det,1,N_det,0,1)
+    print *, 'Done'
     call dtranspose(                                                   &
         v_t,                                                           &
         size(v_t, 1),                                                  &
diff --git a/src/Integrals_Bielec/four_idx_transform.irp.f b/src/Integrals_Bielec/four_idx_transform.irp.f
new file mode 100644
index 00000000..a98ce2aa
--- /dev/null
+++ b/src/Integrals_Bielec/four_idx_transform.irp.f
@@ -0,0 +1,12 @@
+program four_idx
+  implicit none
+  BEGIN_DOC
+! 4-index transformation from AO to MO integrals
+  END_DOC
+  
+  disk_access_mo_integrals = 'Write'
+  SOFT_TOUCH disk_access_mo_integrals
+  if (.true.) then
+    PROVIDE mo_bielec_integrals_in_map
+  endif
+end