QCaml/CI/Determinant.ml

type phase_type =
| Plus
| Minus

type spin_determinant_type =
{
  bitstring : Z.t ;
  phase  : phase_type;
}

type determinant_type =
{
  alpha : Z.t ;
  beta  : Z.t ;
  phase  : phase_type;
}


type t = determinant_type


let spindet_of_list l = 
  let rec aux accu nperm = function
  | [] -> accu, if nperm mod 2 = 0 then Plus else Minus
  | i :: rest -> 
    let i = pred i in
    let x = Z.(shift_left one i) in
    let accu = Z.logor accu x in
    let nperm =
      let mask = Z.(x-one) in
      let r = Z.logand accu mask in
      if r = Z.zero then
        nperm
      else
        nperm + (Z.popcount r)
    in
    aux accu nperm rest
  in
  List.rev l
  |> aux Z.zero 0



(* TODO Phase *)
let rec to_list spindet =
  let rec aux accu z =
    if z <> Z.zero then
      let element = (Z.(trailing_zeros z)+1) in
      aux (element::accu) Z.(z land (pred z))
    else List.rev accu
  in aux [] spindet
  

let of_lists a b =
  let alpha, phase_a = 
    spindet_of_list a
  in
  let beta, phase_b = 
    spindet_of_list b
  in
  let phase =
    match phase_a, phase_b with
    | Plus , Plus -> Plus
    | Minus, Minus-> Plus
    | _ -> Minus
  in
  { alpha ; beta ; phase }

let alpha det = det.alpha
let beta  det = det.beta

let sgn det =
  match det.phase with
  | Plus  ->  1.
  | Minus -> -1.

let pp_spindet ppf spindet =
  String.init (Z.numbits spindet) (fun i -> if (Z.testbit spindet i) then '+' else '-')
  |> Format.fprintf ppf "@[<v> @[<h> %s @]@]"


let pp_det ppf det =
  Format.fprintf ppf "@[<v> a: %a @; b: %a @]@." pp_spindet det.alpha pp_spindet det.beta

let test_case () =

  let test_creation () =
      let l_a = [ 1 ; 2 ; 3 ; 5 ; 64 ]
      and l_b = [ 2 ; 3 ; 5 ; 65 ] in
      let det = of_lists l_a l_b in
      let z_a = alpha det
      and z_b = beta  det in
      Alcotest.(check (list int )) "alpha" (to_list z_a) l_a;
      Alcotest.(check (list int )) "beta"  (to_list z_b) l_b;
      Alcotest.(check bool) "phase" (sgn det = 1.)  true;
  in
  let test_phase () =
      let l_a = [ 1 ; 2 ; 3 ; 64 ; 5 ]
      and l_b = [ 2 ; 3 ; 5 ; 65 ] in
      let det = of_lists l_a l_b in
      Alcotest.(check bool) "phase" (sgn det = -1.)  true;
      let l_a = [ 1 ; 2 ; 3 ; 64 ; 5 ]
      and l_b = [ 3 ; 2 ; 5 ; 65 ] in
      let det = of_lists l_a l_b in
      Alcotest.(check bool) "phase" (sgn det = 1.)  true;
      let l_a = [ 1 ; 3 ; 2 ; 64 ; 5 ]
      and l_b = [ 3 ; 2 ; 5 ; 65 ] in
      let det = of_lists l_a l_b in
      Alcotest.(check bool) "phase" (sgn det = -1.)  true;
      let l_a = [ 1 ; 3 ; 2 ; 64 ; 5 ]
      and l_b = [ 3 ; 2 ; 65 ; 5 ] in
      let det = of_lists l_a l_b in
      Alcotest.(check bool) "phase" (sgn det = 1.)  true;
  in
  [
    "Creation", `Quick, test_creation;
    "Phase"   , `Quick, test_phase;
  ]