diff --git a/Basis/KinInt.ml b/Basis/KinInt.ml
index e10a1ae..0200ddb 100644
--- a/Basis/KinInt.ml
+++ b/Basis/KinInt.ml
@@ -143,7 +143,8 @@ let of_basis basis =
                 try Zmap.find cls key
                 with Not_found -> failwith "Bug in kinetic integrals"
               in
-              result.{i_c,j_c} <- value
+              result.{i_c,j_c} <- value;
+              result.{j_c,i_c} <- value;
             ) (Contracted_shell.powers shell.(i));
         ) (Contracted_shell.powers shell.(j))
     done;
diff --git a/Basis/NucInt.ml b/Basis/NucInt.ml
index cd85b00..f5e0e69 100644
--- a/Basis/NucInt.ml
+++ b/Basis/NucInt.ml
@@ -83,6 +83,7 @@ let of_basis_nuclei basis nuclei =
                 Zmap.find cls key 
               in
               eni_array.{j_c,i_c} <- value;
+              eni_array.{i_c,j_c} <- value;
             ) (Contracted_shell.powers shell.(j))
         ) (Contracted_shell.powers shell.(i));
     done;
diff --git a/Basis/Overlap.ml b/Basis/Overlap.ml
index 22ef83b..cb363b8 100644
--- a/Basis/Overlap.ml
+++ b/Basis/Overlap.ml
@@ -116,7 +116,8 @@ let of_basis basis =
                 try Zmap.find cls key 
                 with Not_found -> failwith "Bug in overlap integrals"
               in
-              result.{i_c,j_c} <- value
+              result.{i_c,j_c} <- value;
+              result.{j_c,i_c} <- value;
             ) (Contracted_shell.powers shell.(i));
         ) (Contracted_shell.powers shell.(j))
     done;
@@ -143,3 +144,6 @@ let to_file ~filename overlap =
   done;
   close_out oc
 
+
+
+
diff --git a/Simulation.ml b/Simulation.ml
index f87723c..9809c80 100644
--- a/Simulation.ml
+++ b/Simulation.ml
@@ -2,6 +2,7 @@ type t = {
   basis : Basis.t;
   nuclei : Nuclei.t;
   overlap : Overlap.t lazy_t;
+  overlap_ortho : Orthonormalization.t lazy_t;
   eN_ints : NucInt.t lazy_t;
   kin_ints : KinInt.t lazy_t;
   ee_ints : ERI.t lazy_t;
@@ -9,9 +10,12 @@ type t = {
 }
 
 let make ~nuclei ~basis =
+  let overlap = 
+    lazy (Overlap.of_basis basis)
+  in
   {
-    basis ; nuclei ;
-    overlap  = lazy (Overlap.of_basis basis);
+    basis ; nuclei ; overlap ;
+    overlap_ortho = lazy (Orthonormalization.make (Lazy.force overlap));
     eN_ints  = lazy (NucInt.of_basis_nuclei basis nuclei);
     kin_ints = lazy (KinInt.of_basis basis);
     ee_ints  = lazy (ERI.of_basis basis);
diff --git a/Utils/Util.ml b/Utils/Util.ml
index b04c0f2..5c67741 100644
--- a/Utils/Util.ml
+++ b/Utils/Util.ml
@@ -3,6 +3,7 @@ external erf_float : float -> float = "erf_float_bytecode" "erf_float" [@@unboxe
 external erfc_float : float -> float = "erfc_float_bytecode" "erfc_float" [@@unboxed]  [@@noalloc]
 external gamma_float : float -> float = "gamma_float_bytecode" "gamma_float" [@@unboxed]  [@@noalloc]
 
+open Lacaml.D
 open Constants
 
 let factmax = 150
@@ -152,3 +153,11 @@ let array_sum a =
 let array_product a = 
    Array.fold_left ( *. ) 0. a
 
+
+let diagonalize_symm h = 
+  let v = lacpy h in
+  let w = Vec.create (Mat.dim1 h) in
+  let result = 
+    syevd ~vectors:true ~w v
+  in
+  v, result