From 9a8f440f7a945a5595774254eeb94bfe52574d5c Mon Sep 17 00:00:00 2001
From: Anthony Scemama <scemama@irsamc.ups-tlse.fr>
Date: Fri, 2 Feb 2018 18:51:17 +0100
Subject: [PATCH] Transposed 0m

---
 Basis/TwoElectronRRVectorized.ml | 27 ++++++++++++++++++++-------
 1 file changed, 20 insertions(+), 7 deletions(-)

diff --git a/Basis/TwoElectronRRVectorized.ml b/Basis/TwoElectronRRVectorized.ml
index d971aec..91b7a76 100644
--- a/Basis/TwoElectronRRVectorized.ml
+++ b/Basis/TwoElectronRRVectorized.ml
@@ -43,10 +43,10 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
       in
       let f = expo_b *. (Coordinate.coord center_ab xyz) in
       Array.mapi (fun k c -> c *.  expo_inv_p *.
-        ( (Coordinate.coord center_pq.(k) xyz) *. zero_m_array.(k).(m+1)
-          -. f *. zero_m_array.(k).(m) ) ) coef_prod
+        ( (Coordinate.coord center_pq.(k) xyz) *. zero_m_array.(m+1).(k)
+          -. f *. zero_m_array.(m).(k) ) ) coef_prod
 
-    | 0 -> Array.mapi (fun k c -> c *. zero_m_array.(k).(m)) coef_prod
+    | 0 -> Array.mapi (fun k c -> c *. zero_m_array.(m).(k)) coef_prod
     | totAngMom_a ->
       let key = Zkey.of_int_tuple (Zkey.Three angMom_a)
       in
@@ -93,8 +93,7 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
     | (i,0) -> if (i>0) then
                  vrr0_v m angMom_a totAngMom_a
                else
-                (* TODO : TRANSPOSE zero_m_array *)
-                 Array.mapi (fun k c -> c *. zero_m_array.(k).(m)) coef_prod
+                 Array.mapi (fun k c -> c *. zero_m_array.(m).(k)) coef_prod
     | (_,_) ->
 
       let key =  Zkey.of_int_tuple (Zkey.Six (angMom_a, angMom_c))
@@ -196,7 +195,7 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
     | 0 ->
       begin
         match (totAngMom_a, totAngMom_c) with
-        | (0,0) -> Array.mapi (fun k c -> c *. zero_m_array.(k).(0)) coef_prod
+        | (0,0) -> Array.mapi (fun k c -> c *. zero_m_array.(0).(k)) coef_prod
         | (_,0) -> vrr0_v 0 angMom_a totAngMom_a
         | (_,_) -> vrr_v 0 angMom_a angMom_c totAngMom_a totAngMom_c
       end
@@ -336,7 +335,7 @@ let contracted_class_shell_pairs ~zero_m ?schwartz_p ?schwartz_q shell_p shell_q
                 in
 
                 let zero_m_array =
-                  zero_m ~maxm ~expo_pq_inv ~norm_pq_sq
+                  zero_m ~maxm:0 ~expo_pq_inv ~norm_pq_sq
                 in
 
                 accu +. coef_prod *. zero_m_array.(0)
@@ -391,6 +390,20 @@ let contracted_class_shell_pairs ~zero_m ?schwartz_p ?schwartz_q shell_p shell_q
           and coef_prod = Array.map (fun (zero_m_array, expo_inv, d, center_cd,
                                           center_pq,coef_prod) -> coef_prod) common
           in
+          (* Transpose zero_m_array 
+          *)
+          let zero_m_array =
+            let result = Array.init (maxm+1) (fun _ ->
+              Array.make (Array.length coef_prod) 0.)
+            in
+            for m=0 to maxm do
+              for k=0 to (Array.length coef_prod-1) do
+               result.(m).(k) <- zero_m_array.(k).(m)
+              done;
+            done;
+            result
+          in
+
             (* Compute the integral class from the primitive shell quartet *)
           let map_1d = Array.init maxm (fun _ -> Zmap.create (4*maxm)) in
           let map_2d = Array.init maxm (fun _ -> Zmap.create (Array.length class_indices)) in