Determinants check the nb of alpha and beta electrons

ocaml/Determinant.ml

@@ -6,11 +6,6 @@ type t = int64 array with sexp
 let to_int64_array (x:t) = (x:int64 array)
 ;;
 
-let of_int64_array n_int x =
-   assert ((Array.length x) = (N_int_number.to_int n_int)*2) ;
-   x
-;;
-
 let to_alpha_beta x = 
   let x = to_int64_array x in
   let n_int = (Array.length x)/2 in
@@ -29,11 +24,20 @@ let to_bitlist_couple x =
   in (xa,xb)
 ;;
 
-let of_bitlist_couple (xa,xb) =
+
+let of_int64_array ~n_int ~alpha ~beta x =
+   assert ((Array.length x) = (N_int_number.to_int n_int)*2) ;
+   let (a,b) = to_bitlist_couple x in
+   assert (Bitlist.popcnt a = Elec_alpha_number.to_int alpha);
+   assert (Bitlist.popcnt b = Elec_beta_number.to_int  beta );
+   x
+;;
+
+let of_bitlist_couple ~alpha ~beta (xa,xb) =
   let ba = Bitlist.to_int64_list xa in
   let bb = Bitlist.to_int64_list xb in
-  let n_int = Bitlist.n_int_of_mo_tot_num (List.length ba) in
-  of_int64_array n_int (Array.of_list (ba@bb))
+  let n_int = Bitlist.n_int_of_mo_tot_num (List.length xa) in
+  of_int64_array ~n_int:n_int ~alpha:alpha ~beta:beta (Array.of_list (ba@bb))
 ;;
 
 let bitlist_to_string ~mo_tot_num x =

ocaml/Determinant.mli

@@ -10,8 +10,12 @@ type t = int64 array with sexp
 (** Transform to an int64 array *)
 val to_int64_array : t -> int64 array
 
-(** Create from an int64 array *)
-val of_int64_array : Qptypes.N_int_number.t -> int64 array -> t
+(** Create from an int64 array, checking the number of alpha 
+  * and beta electrons *)
+val of_int64_array : n_int:Qptypes.N_int_number.t -> 
+  alpha:Qptypes.Elec_alpha_number.t ->
+  beta:Qptypes.Elec_beta_number.t -> 
+  int64 array -> t
 
 (** Split into an alpha-only and a beta-only determinant *)
 val to_alpha_beta : t -> (int64 array)*(int64 array)
@@ -20,7 +24,9 @@ val to_alpha_beta : t -> (int64 array)*(int64 array)
 val to_bitlist_couple : t -> Bitlist.t * Bitlist.t
 
 (** Create from a bit list *)
-val of_bitlist_couple : Bitlist.t * Bitlist.t -> t
+val of_bitlist_couple : alpha:Qptypes.Elec_alpha_number.t ->
+  beta:Qptypes.Elec_beta_number.t -> 
+  Bitlist.t * Bitlist.t -> t
 
 (** String representation *)
 val to_string : mo_tot_num:Qptypes.MO_number.t -> t -> string

ocaml/Input_determinants.ml

@@ -246,11 +246,17 @@ end = struct
             |> List.map ~f:Det_coef.to_float
     in
     Ezfio.ezfio_array_of_list ~rank:1 ~dim:[| n_det |] ~data:c
+    |> Ezfio.set_determinants_psi_coef 
   ;;
 
 
   let read_psi_det () =
-    let n_int = read_n_int () in
+    let n_int = read_n_int () 
+    and n_alpha = Ezfio.get_electrons_elec_alpha_num ()
+        |> Elec_alpha_number.of_int 
+    and n_beta = Ezfio.get_electrons_elec_beta_num ()
+        |> Elec_beta_number.of_int 
+    in
     if not (Ezfio.has_determinants_psi_det ()) then
       begin
         let mo_tot_num = MO_number.get_max () in
@@ -258,9 +264,9 @@ end = struct
           | 0 -> accu
           | n -> build_data ((MO_number.of_int ~max:mo_tot_num n)::accu) (n-1)
         in
-        let det_a = build_data [] (Ezfio.get_electrons_elec_alpha_num ())
+        let det_a = build_data [] (Elec_alpha_number.to_int n_alpha)
           |> Bitlist.of_mo_number_list n_int
-        and det_b = build_data [] (Ezfio.get_electrons_elec_beta_num  ())
+        and det_b = build_data [] (Elec_beta_number.to_int n_beta)
           |> Bitlist.of_mo_number_list n_int
         in
         let data = ( (Bitlist.to_int64_list det_a) @ 
@@ -280,7 +286,8 @@ end = struct
     List.init dim.(2) ~f:(fun i ->
       Array.sub ~pos:(2*n_int*i) ~len:(2*n_int) data)
     |> List.map ~f:(Determinant.of_int64_array
-      (N_int_number.of_int n_int))
+      ~n_int:(N_int_number.of_int n_int)
+      ~alpha:n_alpha ~beta:n_beta )
     |> Array.of_list
   ;;
 
@@ -500,6 +507,11 @@ psi_det                = %s
 
     (* Handle determinants *)
     let psi_det = 
+      let n_alpha = Ezfio.get_electrons_elec_alpha_num ()
+        |> Elec_alpha_number.of_int 
+      and n_beta = Ezfio.get_electrons_elec_beta_num ()
+        |> Elec_beta_number.of_int 
+      in
       let rec read_dets accu = function
       | [] -> List.rev accu
       | ""::c::alpha::beta::tail -> 
@@ -507,7 +519,8 @@ psi_det                = %s
             let alpha = String.rev alpha |> Bitlist.of_string ~zero:'-' ~one:'+'
             and beta  = String.rev beta  |> Bitlist.of_string ~zero:'-' ~one:'+'
             in
-            let newdet = Determinant.of_bitlist_couple (alpha,beta) 
+            let newdet = Determinant.of_bitlist_couple 
+              ~alpha:n_alpha ~beta:n_beta (alpha,beta) 
             |> Determinant.sexp_of_t |> Sexplib.Sexp.to_string
             in
             read_dets (newdet::accu) tail
@@ -525,8 +538,6 @@ psi_det                = %s
     and n_int = Printf.sprintf "(n_int %d)" (N_int_number.get_max ()) in
     let s = String.concat [ header ; bitkind ; n_int ; psi_coef ; psi_det]
     in
-    ("("^s^")")
-    |> print_endline ;
     Sexp.of_string ("("^s^")")
     |> t_of_sexp
   ;;

ocaml/test_determinants.ml

@@ -4,7 +4,10 @@ let test_module () =
   let mo_tot_num = MO_number.of_int 10 in
   let det =
    [| 15L ; 7L |]
-   |> Determinant.of_int64_array (N_int_number.of_int 1)
+   |> Determinant.of_int64_array 
+      ~n_int:(N_int_number.of_int 1)
+      ~alpha:(Elec_alpha_number.of_int 4)
+      ~beta:(Elec_beta_number.of_int 3)
   in
   Printf.printf "%s\n" (Determinant.to_string (~mo_tot_num:mo_tot_num) det)
 ;;