From 55b58268f872405e13f522d3dc093e042477151d Mon Sep 17 00:00:00 2001
From: Anthony Scemama <scemama@irsamc.ups-tlse.fr>
Date: Wed, 29 May 2019 09:04:43 +0200
Subject: [PATCH] F12 without singles

---
 CI/F12CI.ml | 18 +++++++++---------
 1 file changed, 9 insertions(+), 9 deletions(-)

diff --git a/CI/F12CI.ml b/CI/F12CI.ml
index f727fd6..db70656 100644
--- a/CI/F12CI.ml
+++ b/CI/F12CI.ml
@@ -61,7 +61,7 @@ let hf_ij mo_basis ki kj =
   [ h_ij mo_basis ki kj ; f_ij mo_basis ki kj ]
 
 
-let is_internal det_space =
+let is_a_double det_space =
     let mo_class      = DeterminantSpace.mo_class det_space in
     let mo_num = Array.length @@ MOClass.mo_class_array mo_class in
     let m l =
@@ -83,8 +83,8 @@ let is_internal det_space =
       and b = Bitstring.logand aux_mask beta
       in
       match Bitstring.popcount a + Bitstring.popcount b with
-      | 1 | 2 -> false
-      | _ -> true
+      | 2 -> true
+      | _ -> false
 
 
 let dressing_vector ~frozen_core aux_basis f12_amplitudes ci =
@@ -107,13 +107,12 @@ let dressing_vector ~frozen_core aux_basis f12_amplitudes ci =
         |> DeterminantSpace.determinant_stream
 
       in
-      (* Select only singly and doubly excited determinants
-          wrt FCI space *)
+      (* Select only doubly excited determinants wrt FCI space *)
       Stream.from (fun _ ->
         try
           let rec result () = 
             let ki = Stream.next s in
-            if is_internal ci.CI.det_space ki then
+            if not (is_a_double ci.CI.det_space ki) then
                 result ()
             else
                 Some ki
@@ -227,15 +226,12 @@ Printf.printf "Add aux basis\n%!";
     MOBasis.of_mo_basis s mo_basis
   in
   let () = 
-Printf.printf "F12 ints\n%!";
     ignore @@ MOBasis.f12_ints aux_basis
   in
   let () =
-Printf.printf "2e ints\n%!";
     ignore @@ MOBasis.two_e_ints aux_basis
   in
 
-Printf.printf "det space\n%!";
   let det_space = 
     DeterminantSpace.fci_f12_of_mo_basis  aux_basis  ~frozen_core  mo_num
   in
@@ -269,12 +265,16 @@ Printf.printf "det space\n%!";
          |> Matrix.to_mat
       in
 
+      Printf.printf "Cmax : %e\n" psi.{column_idx,state};
+      Printf.printf "Norm : %e\n" (sqrt (gemm ~transa:`T delta delta).{state,state});
+
       let f = 1.0 /. psi.{column_idx,state} in
       let delta_00 =
         (* Delta_00 = {% $\sum_{j \ne x} delta_j c_j / c_x$ %} *)
         f *. ( (gemm ~transa:`T delta psi).{state,state} -.
           delta.{column_idx,state} *. psi.{column_idx,state} )
       in
+      Printf.printf "Delta_00 : %e %e\n" delta.{column_idx,state} delta_00;
       delta.{column_idx,state} <- delta.{column_idx,state} -. delta_00;
 
       let diagonal =