From c9c039c2fe800cbd1f867242e450c84f2c939600 Mon Sep 17 00:00:00 2001
From: Anthony Scemama <scemama@irsamc.ups-tlse.fr>
Date: Tue, 5 Mar 2019 00:39:23 +0100
Subject: [PATCH] Improved Huckel guess

---
 CI/DeterminantSpace.ml |  4 ++--
 SCF/Guess.ml           | 34 ++++++++++++++++++----------------
 2 files changed, 20 insertions(+), 18 deletions(-)

diff --git a/CI/DeterminantSpace.ml b/CI/DeterminantSpace.ml
index 7345e46..bbc2586 100644
--- a/CI/DeterminantSpace.ml
+++ b/CI/DeterminantSpace.ml
@@ -153,13 +153,12 @@ let fci_of_mo_basis ?(frozen_core=true) mo_basis =
     let n_det_beta = Array.length det_beta in
     let n_det_alfa = Array.length det_alfa in
 
-    (*
     let ndet = n_det_alfa * n_det_beta in
     Format.printf "Number of determinants : %d %d %d\n%!" 
       n_det_alfa n_det_beta ndet;
     Spin (det_alfa, det_beta)
-    *)
 
+    (*
     let det = Array.make n_det_alfa
         (Array.init n_det_beta (fun i -> i))
     in
@@ -171,6 +170,7 @@ let fci_of_mo_basis ?(frozen_core=true) mo_basis =
     Arbitrary {
       det_alfa ; det_beta ; det ; index_start
     }
+    *)
   in
   { n_alfa ; n_beta ; mo_class ; mo_basis ; determinants }
 
diff --git a/SCF/Guess.ml b/SCF/Guess.ml
index aa8d1e4..7813554 100644
--- a/SCF/Guess.ml
+++ b/SCF/Guess.ml
@@ -21,29 +21,31 @@ let hcore_guess ao_basis =
 
 let huckel_guess ao_basis =
   let c = 0.5 *. 1.75 in
-  let ao_num   = Ao.basis    ao_basis |> Basis.size in
   let eN_ints  = Ao.eN_ints  ao_basis |> NucInt.matrix
   and kin_ints = Ao.kin_ints ao_basis |> KinInt.matrix
-  and overlap  = Ao.overlap  ao_basis |> Ov.matrix
-  and m_X      = Ao.ortho    ao_basis
   in
-  let diag = Array.init (ao_num+1) (fun i -> if i=0 then 0. else
-    eN_ints.{i,i} +. kin_ints.{i,i})
+  let m_F = 
+    Mat.add eN_ints kin_ints
+  in
+  let ao_num   = Ao.basis    ao_basis |> Basis.size 
+  and overlap  = Ao.overlap  ao_basis |> Ov.matrix
+  in
+  let diag = Vec.init ao_num (fun i -> 
+    m_F.{i,i} )
   in
 
   function 
   | 0 -> invalid_arg "Huckel guess needs a non-zero number of occupied MOs."
-  | nocc ->
-    let density = gemm ~alpha:2. ~transb:`T ~k:nocc m_X m_X in
-    let f = Fock.make_rhf ~density ao_basis in
-    let m_F = Fock.fock f in
-    for j=1 to ao_num do
-      for i=1 to ao_num do
-        if (i <> j) then
-          m_F.{i,j} <- c *. overlap.{i,j} *. (diag.(i) +. diag.(j)) (*TODO Pseudo *)
-      done;
-    done;
-    m_F
+  | nocc -> 
+    Mat.init_cols ao_num ao_num (fun i j ->
+        if (i<>j) then
+          if (diag.{i} +. diag.{j}) < 5. then
+            c *. overlap.{i,j} *. (diag.{i} +. diag.{j}) +. m_F.{i,j} (*TODO Pseudo *)
+          else
+            m_F.{i,j} (*TODO Pseudo *)
+        else
+          diag.{i}
+    )
 
 
 let make ?(nocc=0) ~guess ao_basis =