diff --git a/CI/CI.ml b/CI/CI.ml index bef7328..d7476be 100644 --- a/CI/CI.ml +++ b/CI/CI.ml @@ -637,12 +637,14 @@ let second_order_sum2 { det_space ; m_H ; m_S2 ; eigensystem ; n_states } let is_internal det_space = + let mo_class = DeterminantSpace.mo_class det_space in + let numbits = Array.length @@ MOClass.mo_class_array mo_class in let m l = List.fold_left (fun accu i -> - let j = i-1 in Bitstring.(logor accu (shift_left_one j)) - ) Bitstring.zero l + let j = i-1 in + Bitstring.logor accu (Bitstring.shift_left_one numbits j) + ) (Bitstring.zero numbits) l in - let mo_class = DeterminantSpace.mo_class det_space in let active_mask = m (MOClass.active_mos mo_class) in let occ_mask = m (MOClass.core_mos mo_class) in let inactive_mask = m (MOClass.inactive_mos mo_class) in diff --git a/CI/Determinant.ml b/CI/Determinant.ml index 04bce94..d932db1 100644 --- a/CI/Determinant.ml +++ b/CI/Determinant.ml @@ -14,10 +14,10 @@ let alfa t = t.alfa let beta t = t.beta -let vac = +let vac n = { - alfa = Spindeterminant.vac ; - beta = Spindeterminant.vac ; + alfa = Spindeterminant.vac n; + beta = Spindeterminant.vac n; } let phase t = @@ -55,9 +55,9 @@ let degree t t' = (degree_alfa t t') + (degree_beta t t') -let of_lists a b = - let alfa = Spindeterminant.of_list a - and beta = Spindeterminant.of_list b +let of_lists n a b = + let alfa = Spindeterminant.of_list n a + and beta = Spindeterminant.of_list n b in of_spindeterminants alfa beta @@ -103,7 +103,7 @@ 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 det = of_lists 66 l_a l_b in let z_a = alfa det and z_b = beta det in Alcotest.(check (list int )) "alfa" (Spindeterminant.to_list z_a) l_a; @@ -114,19 +114,19 @@ let test_case () = 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 + let det = of_lists 66 l_a l_b in Alcotest.(check bool) "phase" (phase det = Phase.Neg) 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 + let det = of_lists 66 l_a l_b in Alcotest.(check bool) "phase" (phase det = Phase.Pos) 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 + let det = of_lists 66 l_a l_b in Alcotest.(check bool) "phase" (phase det = Phase.Neg) 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 + let det = of_lists 66 l_a l_b in Alcotest.(check bool) "phase" (phase det = Phase.Pos) true; in @@ -134,22 +134,22 @@ let test_case () = let det = let open Spin in creation Alfa 1 @@ creation Alfa 3 @@ creation Alfa 2 @@ creation Alfa 5 @@ - creation Beta 1 @@ creation Beta 3 @@ creation Beta 4 @@ creation Beta 5 @@ vac + creation Beta 1 @@ creation Beta 3 @@ creation Beta 4 @@ creation Beta 5 @@ vac 10 in Alcotest.(check bool) "creation 1" true - (det = of_lists [ 1 ; 3 ; 2 ; 5 ] [1 ; 3 ; 4 ; 5 ] ); + (det = of_lists 10 [ 1 ; 3 ; 2 ; 5 ] [1 ; 3 ; 4 ; 5 ] ); let det' = single_excitation Spin.Alfa 3 6 det in Alcotest.(check bool) "single_exc 1" true - (det' = of_lists [ 1 ; 6 ; 2 ; 5 ] [1 ; 3 ; 4 ; 5 ] ); + (det' = of_lists 10 [ 1 ; 6 ; 2 ; 5 ] [1 ; 3 ; 4 ; 5 ] ); let det' = single_excitation Spin.Beta 3 6 det in Alcotest.(check bool) "single_exc 2" true - (det' = of_lists [ 1 ; 3 ; 2 ; 5 ] [1 ; 6 ; 4 ; 5 ] ); + (det' = of_lists 10 [ 1 ; 3 ; 2 ; 5 ] [1 ; 6 ; 4 ; 5 ] ); let det' = single_excitation Spin.Alfa 4 6 det @@ -164,26 +164,26 @@ let test_case () = let det' = double_excitation Spin.Alfa 3 6 Spin.Alfa 2 7 det in - let det'' = of_lists [ 1 ; 6 ; 7 ; 5 ] [1 ; 3 ; 4 ; 5 ] in + let det'' = of_lists 10 [ 1 ; 6 ; 7 ; 5 ] [1 ; 3 ; 4 ; 5 ] in Alcotest.(check bool) "double_exc 1" true (det' = det''); let det' = double_excitation Spin.Beta 3 6 Spin.Beta 5 7 det in Alcotest.(check bool) "double_exc 2" true - (det' = of_lists [ 1 ; 3 ; 2 ; 5 ] [1 ; 6 ; 4 ; 7 ] ); + (det' = of_lists 10 [ 1 ; 3 ; 2 ; 5 ] [1 ; 6 ; 4 ; 7 ] ); let det' = double_excitation Spin.Alfa 3 6 Spin.Beta 5 7 det in Alcotest.(check bool) "double_exc 3" true - (det' = of_lists [ 1 ; 6 ; 2 ; 5 ] [1 ; 3 ; 4 ; 7 ] ); + (det' = of_lists 10 [ 1 ; 6 ; 2 ; 5 ] [1 ; 3 ; 4 ; 7 ] ); let det' = double_excitation Spin.Beta 5 7 Spin.Alfa 3 6 det in Alcotest.(check bool) "double_exc 4" true - (det' = of_lists [ 1 ; 6 ; 2 ; 5 ] [1 ; 3 ; 4 ; 7 ] ); + (det' = of_lists 10 [ 1 ; 6 ; 2 ; 5 ] [1 ; 3 ; 4 ; 7 ] ); let det' = double_excitation Spin.Alfa 4 6 Spin.Alfa 2 7 det diff --git a/CI/Determinant.mli b/CI/Determinant.mli index 2c75204..aa0836f 100644 --- a/CI/Determinant.mli +++ b/CI/Determinant.mli @@ -29,8 +29,9 @@ val is_none : t -> bool (** {1 Second quantization operators} *) -val vac : t -(** Vacuum state, [vac = Some ]{% $|\rangle$ %} *) +val vac : int -> t +(** Vacuum state, [vac = Some ]{% $|\rangle$ %}. The integer parameter is the size of the + MO basis set. *) val creation : Spin.t -> particle -> t -> t (** [creation spin p] is the creation operator {% $a^\dagger_p$ %}. *) @@ -65,8 +66,9 @@ val of_spindeterminants : Spindeterminant.t -> Spindeterminant.t -> t [Spindeterminant.t]. *) -val of_lists : int list -> int list -> t -(** Creates a Slater determinant from a two lists of orbital indices. *) +val of_lists : int -> int list -> int list -> t +(** Creates a Slater determinant from a two lists of orbital indices. + The integer parameter is the size of the MO basis set. *) val negate_phase : t -> t (** Returns the same determinant with the phase negated. *) diff --git a/CI/Excitation.ml b/CI/Excitation.ml index e9854b8..6b154da 100644 --- a/CI/Excitation.ml +++ b/CI/Excitation.ml @@ -18,6 +18,7 @@ let single_of_spindet t t' = and d' = Spindeterminant.bitstring t' in let tmp = Bitstring.logxor d d' in + let shift_left_one = Bitstring.(shift_left_one (numbits tmp)) in let hole_z = Bitstring.logand (Spindeterminant.bitstring t ) tmp and particle_z = Bitstring.logand (Spindeterminant.bitstring t') tmp in @@ -32,8 +33,8 @@ let single_of_spindet t t' = let mask = let h = high-1 in let l = low in - let mask_up = Bitstring.shift_left_one h |> Bitstring.minus_one - and mask_dn = Bitstring.plus_one @@ Bitstring.lognot (Bitstring.shift_left_one l) + let mask_up = shift_left_one h |> Bitstring.minus_one + and mask_dn = Bitstring.plus_one @@ Bitstring.lognot (shift_left_one l) in Bitstring.logand mask_up mask_dn in let phase = @@ -150,7 +151,7 @@ let test_case () = let test_single () = let l_a = [ 1 ; 2 ; 3 ; 5 ; 64 ] and l_b = [ 2 ; 3 ; 5 ; 65 ] in - let det1 = Determinant.of_lists l_a l_b in + let det1 = Determinant.of_lists 66 l_a l_b in let det2 = Determinant.single_excitation Spin.Alfa 3 7 det1 in let t = single_of_det det1 det2 in Alcotest.(check bool) "single 1" true (t = Single (Phase.Pos, { hole=3 ; particle=7 ; spin=Spin.Alfa}) ); @@ -171,7 +172,7 @@ let test_case () = let test_double () = let l_a = [ 1 ; 2 ; 3 ; 5 ; 64 ] and l_b = [ 2 ; 3 ; 5 ; 65 ] in - let det1 = Determinant.of_lists l_a l_b in + let det1 = Determinant.of_lists 66 l_a l_b in let det2 = Determinant.double_excitation Spin.Alfa 3 7 Spin.Alfa 2 6 det1 in let t = double_of_det det1 det2 in Alcotest.(check bool) "double 1" true diff --git a/CI/F12CI.ml b/CI/F12CI.ml index 9143dc7..781ce40 100644 --- a/CI/F12CI.ml +++ b/CI/F12CI.ml @@ -64,18 +64,29 @@ let f_ij mo_basis ki kj = let is_internal det_space = + let mo_class = DeterminantSpace.mo_class det_space in + let mo_num = Array.length @@ MOClass.mo_class_array mo_class in let m l = List.fold_left (fun accu i -> - let j = i-1 in Bitstring.(logor accu (shift_left_one j)) - ) Bitstring.zero l + let j = i-1 in Bitstring.logor accu (Bitstring.shift_left_one mo_num j) + ) (Bitstring.zero mo_num) l in - let mo_class = DeterminantSpace.mo_class det_space in - let aux_mask = m (MOClass.auxiliary_mos mo_class) in + let aux_mask = m (MOClass.auxiliary_mos mo_class) in fun a -> let alfa = Determinant.alfa a |> Spindeterminant.bitstring in + (* + let beta = + Determinant.beta a + |> Spindeterminant.bitstring + in + let a = Bitstring.logand aux_mask alfa + and b = Bitstring.logand aux_mask beta + in Bitstring.popcount a + Bitstring.popcount b < 2 + *) + if not (Bitstring.logand aux_mask alfa |> Bitstring.is_zero ) then false else diff --git a/CI/Spindeterminant.ml b/CI/Spindeterminant.ml index 241ae0b..1f6aff2 100644 --- a/CI/Spindeterminant.ml +++ b/CI/Spindeterminant.ml @@ -23,8 +23,8 @@ let bitstring = function | None -> invalid_arg "Spindeterminant is None" -let vac = - Some { bitstring = Bitstring.zero; +let vac n = + Some { bitstring = Bitstring.zero n; phase = Phase.Pos; } @@ -36,7 +36,8 @@ let creation p = function None else begin - let x = Bitstring.shift_left_one i in + let numbits = Bitstring.numbits spindet.bitstring in + let x = Bitstring.shift_left_one numbits i in let bitstring = Bitstring.logor spindet.bitstring x in let mask = Bitstring.minus_one x in let r = Bitstring.logand bitstring mask in @@ -53,7 +54,8 @@ let annihilation h = function None else begin - let x = Bitstring.shift_left_one i in + let numbits = Bitstring.numbits spindet.bitstring in + let x = Bitstring.shift_left_one numbits i in let mask = Bitstring.minus_one x in let r = Bitstring.logand spindet.bitstring mask in let phase = Phase.add_nperm spindet.phase (Bitstring.popcount r) in @@ -81,16 +83,16 @@ let degree t = let holes_of t t' = Bitstring.logand (bitstring t) (Bitstring.logxor (bitstring t) (bitstring t')) - |> Bitstring.to_list [] + |> Bitstring.to_list let particles_of t t' = Bitstring.logand (bitstring t') (Bitstring.logxor (bitstring t) (bitstring t')) - |> Bitstring.to_list [] + |> Bitstring.to_list let holes_particles_of t t' = let x = Bitstring.logxor (bitstring t) (bitstring t') in - let holes = Bitstring.logand (bitstring t) x |> Bitstring.to_list [] - and particles = Bitstring.logand (bitstring t') x |> Bitstring.to_list [] + let holes = Bitstring.logand (bitstring t) x |> Bitstring.to_list + and particles = Bitstring.logand (bitstring t') x |> Bitstring.to_list in List.map2 (fun h p -> (h,p)) holes particles @@ -103,9 +105,9 @@ let negate_phase = function let of_bitstring ?(phase=Phase.Pos) bitstring = Some { bitstring ; phase } -let of_list l = +let of_list n l = List.rev l - |> List.fold_left (fun accu p -> creation p accu) vac + |> List.fold_left (fun accu p -> creation p accu) (vac n) let rec to_list = function @@ -127,7 +129,7 @@ let n_electrons = function let pp_spindet n ppf = function | None -> Format.fprintf ppf "@[None@]" | Some s -> - Format.fprintf ppf "@[%a %a@]" Phase.pp_phase s.phase (Bitstring.pp_bitstring n) + Format.fprintf ppf "@[%a %a@]" Phase.pp_phase s.phase Bitstring.pp s.bitstring @@ -140,36 +142,36 @@ let test_case () = let test_creation () = let l_a = [ 1 ; 2 ; 3 ; 5 ] in - let det = of_list l_a in + let det = of_list 10 l_a in Alcotest.(check (list int )) "bitstring 1" l_a (to_list det); Alcotest.(check bool) "phase 2" true (phase det = Phase.Pos); let l_b = [ 1 ; 3 ; 2 ; 5 ] in - let det = of_list l_b in + let det = of_list 10 l_b in Alcotest.(check (list int )) "bitstring 2" l_a (to_list det); Alcotest.(check bool) "phase 2" true (phase det = Phase.Neg); in let test_a_operators () = let det = - creation 5 @@ creation 2 @@ creation 2 @@ creation 1 @@ vac + creation 5 @@ creation 2 @@ creation 2 @@ creation 1 @@ (vac 10) in Alcotest.(check bool) "none 1" true (is_none det); let det = - creation 5 @@ creation 3 @@ creation 2 @@ creation 1 @@ vac + creation 5 @@ creation 3 @@ creation 2 @@ creation 1 @@ (vac 10) in let l_a = [ 1 ; 2 ; 3 ; 5 ] in Alcotest.(check (list int )) "bitstring 1" l_a (to_list det); Alcotest.(check bool) "phase 1" true (phase det = Phase.Pos); let det = - creation 1 @@ creation 3 @@ creation 2 @@ creation 5 @@ vac + creation 1 @@ creation 3 @@ creation 2 @@ creation 5 @@ (vac 10) in Alcotest.(check (list int )) "bitstring 2" l_a (to_list det); Alcotest.(check bool) "phase 2" true (phase det = Phase.Neg); let l_b = [ 1 ; 3 ; 2 ; 5 ] in - let det = of_list l_b in + let det = of_list 10 l_b in Alcotest.(check (list int )) "bitstring 3" l_a (to_list det); Alcotest.(check bool) "phase 3" true (phase det = Phase.Neg); @@ -197,16 +199,19 @@ let test_case () = let test_exc_operators () = let l_a = [ 1 ; 2 ; 3 ; 5 ] in - let det = of_list l_a in + let det = of_list 10 l_a in let l_b = [ 1 ; 7 ; 3 ; 5 ] in - let det2 = of_list l_b in + let det2 = of_list 10 l_b in + Format.printf "%a@." (pp_spindet 7) det; + Format.printf "%a@." (pp_spindet 7) det2; + Format.printf "%a@." (pp_spindet 7) (single_excitation_reference 2 7 det); Alcotest.(check bool) "single 1" true (single_excitation_reference 2 7 det = det2); Alcotest.(check bool) "single 2" true (single_excitation 2 7 det = single_excitation_reference 2 7 det); Alcotest.(check bool) "single 3" true (single_excitation_reference 4 7 det |> is_none); Alcotest.(check bool) "single 4" true (single_excitation 4 7 det |> is_none); let l_c = [ 1 ; 7 ; 6 ; 5 ] in - let det3 = of_list l_c in + let det3 = of_list 10 l_c in Alcotest.(check bool) "double 1" true (double_excitation_reference 2 7 3 6 det = det3); Alcotest.(check bool) "double 2" true (double_excitation 2 7 3 6 det = double_excitation_reference 2 7 3 6 det); Alcotest.(check bool) "double 3" true (double_excitation_reference 4 7 3 6 det |> is_none); @@ -215,14 +220,14 @@ let test_case () = let test_exc_spindet () = let l_a = [ 1 ; 2 ; 3 ; 5 ] in - let det = of_list l_a in + let det = of_list 10 l_a in let l_b = [ 1 ; 7 ; 3 ; 5 ] in - let det2 = of_list l_b in + let det2 = of_list 10 l_b in Alcotest.(check int) "single" 1 (degree det det2); Alcotest.(check (list int)) "holes" [2] (holes_of det det2); Alcotest.(check (list int)) "particles" [7] (particles_of det det2); let l_b = [ 1 ; 7 ; 3 ; 6 ] in - let det2 = of_list l_b in + let det2 = of_list 10 l_b in Alcotest.(check int) "double" 2 (degree det det2); Alcotest.(check (list int)) "holes" [2 ; 5] (holes_of det det2); Alcotest.(check (list int)) "particles" [6 ; 7] (particles_of det det2); diff --git a/CI/Spindeterminant.mli b/CI/Spindeterminant.mli index e9df5a1..d87df4f 100644 --- a/CI/Spindeterminant.mli +++ b/CI/Spindeterminant.mli @@ -29,8 +29,9 @@ val negate_phase : t -> t (** {1 Second quantization operators} *) -val vac : t -(** Vacuum state, [vac = Some ]{% $|\rangle$ %} *) +val vac : int -> t +(** Vacuum state, [vac = Some ]{% $|\rangle$ %}. The integer parameter contains the + number of orbitals in the basis set. *) val creation : particle -> t -> t (** [creation p] is the creation operator {% $a^\dagger_p$ %}. *) @@ -66,10 +67,11 @@ val n_electrons : t -> int val of_bitstring : ?phase:Phase.t -> Bitstring.t -> t (** Creates from a bitstring and an optional phase.*) -val of_list : int list -> t +val of_list : int -> int list -> t (** Builds a spin-determinant from a list of orbital indices. If the creation of the spin-determinant is not possible because of Pauli's exclusion principle, a [None] - spin-determinant is returned. *) + spin-determinant is returned. + The first integer is the size of the MO basis set. *) val to_list : t -> int list (** Transforms a spin-determinant into a list of orbital indices. *) diff --git a/CI/SpindeterminantSpace.ml b/CI/SpindeterminantSpace.ml index d5da05c..589687b 100644 --- a/CI/SpindeterminantSpace.ml +++ b/CI/SpindeterminantSpace.ml @@ -17,8 +17,9 @@ let fci_of_mo_basis ~frozen_core mo_basis elec_num = let mo_num = MOBasis.size mo_basis in let mo_class = MOClass.fci ~frozen_core mo_basis in let m l = - List.fold_left (fun accu i -> let j = i-1 in Bitstring.(logor accu (shift_left_one j)) - ) Bitstring.zero l + List.fold_left (fun accu i -> let j = i-1 in + Bitstring.logor accu (Bitstring.shift_left_one mo_num j) + ) (Bitstring.zero mo_num) l in let occ_mask = m (MOClass.core_mos mo_class) and active_mask = m (MOClass.active_mos mo_class) @@ -39,8 +40,9 @@ let cas_of_mo_basis mo_basis ~frozen_core elec_num n m = let mo_num = MOBasis.size mo_basis in let mo_class = MOClass.cas_sd ~frozen_core mo_basis n m in let m l = - List.fold_left (fun accu i -> let j = i-1 in Bitstring.(logor accu (Bitstring.shift_left_one j)) - ) Bitstring.zero l + List.fold_left (fun accu i -> let j = i-1 in + Bitstring.logor accu (Bitstring.shift_left_one mo_num j) + ) (Bitstring.zero mo_num) l in let active_mask = m (MOClass.active_mos mo_class) in let occ_mask = m (MOClass.core_mos mo_class) in diff --git a/Utils/Bitstring.ml b/Utils/Bitstring.ml index 46fe8e6..3e6382c 100644 --- a/Utils/Bitstring.ml +++ b/Utils/Bitstring.ml @@ -1,131 +1,171 @@ +module One = struct + type t = int + + + let of_int x = + assert (x > 0); x + + let numbits _ = 63 + let zero = 0 + let is_zero x = x = 0 + let shift_left x i = x lsl i + let shift_right x i = x lsr i + let shift_left_one i = 1 lsl i + let testbit x i = ( (x lsr i) land 1 ) = 1 + let logor a b = a lor b + let logxor a b = a lxor b + let logand a b = a land b + let lognot a = lnot a + let minus_one a = a - 1 + let plus_one a = a + 1 + + let popcount = function + | 0 -> 0 + | r -> Util.popcnt (Int64.of_int r) + + + let trailing_zeros r = + Util.trailz (Int64.of_int r) + + + let hamdist a b = + a lxor b + |> Int64.of_int + |> Util.popcnt + + + let pp ppf s = + Format.fprintf ppf "@[@[%a@]i@]" (Util.pp_bitstring 64) + (Z.of_int s) + +end + + + + +module Many = struct + + type t = Z.t + + let of_int = Z.of_int + let of_z x = x + let zero = Z.zero + let is_zero x = x = Z.zero + let shift_left = Z.shift_left + let shift_right = Z.shift_right + let shift_left_one i = Z.shift_left Z.one i + let testbit = Z.testbit + let logor = Z.logor + let logxor = Z.logxor + let logand = Z.logand + let lognot = Z.lognot + let minus_one = Z.pred + let plus_one = Z.succ + let trailing_zeros = Z.trailing_zeros + let hamdist = Z.hamdist + let numbits i = max (Z.numbits i) 64 + + let popcount z = + if z = Z.zero then 0 else Z.popcount z + + let pp ppf s = + Format.fprintf ppf "@[@[%a@]m@]" (Util.pp_bitstring (Z.numbits s)) s + +end + + type t = -| One of int +| One of int | Many of Z.t +let of_int x = + One (One.of_int x) -(* -let of_int x = One x -*) -let of_int x = Many (Z.of_int x) +let of_z x = + if Z.numbits x < 64 then One (Z.to_int x) else Many (Many.of_z x) -let of_z x = if (Z.lt x (Z.of_int max_int)) then One (Z.to_int x) else Many x - -(* -let zero = One 0 -*) -let zero = Many Z.zero +let zero = function +| n when n < 64 -> One (One.zero) +| _ -> Many (Many.zero) +let numbits = function +| One x -> One.numbits x +| Many x -> Many.numbits x let is_zero = function -| One x -> assert false (* x = 0 *) -| Many x -> x = Z.zero +| One x -> One.is_zero x +| Many x -> Many.is_zero x +let shift_left x i = match x with +| One x -> One (One.shift_left x i) +| Many x -> Many (Many.shift_left x i) -let shift_left x i = - match x with - | One x -> assert false (* - let y = x lsl i in - if y lsr i = x then - One y - else - Many (Z.shift_left (Z.of_int x) i) - *) - | Many x -> Many (Z.shift_left x i) +let shift_right x i = match x with +| One x -> One (One.shift_right x i) +| Many x -> Many (Many.shift_right x i) +let shift_left_one = function +| n when n < 64 -> fun i -> One (One.shift_left_one i) +| _ -> fun i -> Many (Many.shift_left_one i) -let shift_right x i = - match x with - | One x -> assert false (* One (x lsr i) *) - | Many x -> Z.shift_right x i |> of_z - - -let shift_left_one = - let memo = - Array.init 512 (fun i -> - (* - if i < 63 then - One (1 lsl i) - else - *) - Many (Z.(shift_left one i))) - in - fun i -> - if i < 512 then - memo.(i) - else - Many (Z.(shift_left one i)) - - -let testbit bs i = - match bs with - | One one -> assert false (* (one lsr i) land 1 = 1 *) - | Many z -> Z.testbit z i - +let testbit = function +| One x -> One.testbit x +| Many x -> Many.testbit x let logor a b = match a,b with - | One a, One b -> assert false (* One (a lor b) *) - | One a, Many b -> Many (Z.logor (Z.of_int a) b) - | Many a, One b -> Many (Z.logor a (Z.of_int b)) - | Many a, Many b -> Many (Z.logor a b) - + | One a, One b -> One (One.logor a b) + | Many a, Many b -> Many (Many.logor a b) + | _ -> invalid_arg "Bitstring.logor" let logxor a b = match a,b with - | One a, One b -> assert false (* One (a lxor b) *) - | One a, Many b -> Many (Z.logxor (Z.of_int a) b) - | Many a, One b -> Many (Z.logxor a (Z.of_int b)) - | Many a, Many b -> Many (Z.logxor a b) - + | One a, One b -> One (One.logxor a b) + | Many a, Many b -> Many (Many.logxor a b) + | _ -> invalid_arg "Bitstring.logxor" let logand a b = match a,b with - | One a, One b -> assert false (* One (a land b) *) - | One a, Many b -> Many (Z.logand (Z.of_int a) b) - | Many a, One b -> Many (Z.logand a (Z.of_int b)) - | Many a, Many b -> Many (Z.logand a b) - + | One a, One b -> One (One.logand a b) + | Many a, Many b -> Many (Many.logand a b) + | _ -> invalid_arg "Bitstring.logand" let lognot = function - | One a -> Many ( Z.(lognot @@ of_int a) ) - | Many a -> Many ( Z.lognot a ) - +| One x -> One (One.lognot x) +| Many x -> Many (Many.lognot x) let minus_one = function - | One a -> assert false (* One ( a-1 ) *) - | Many a -> Many ( Z.(a-one) ) - +| One x -> One (One.minus_one x) +| Many x -> Many (Many.minus_one x) let plus_one = function - | One a -> assert false (* One ( a+1 ) *) - | Many a -> Many ( Z.(a+one) ) - - -let popcount = function - | One r -> assert false (* Util.popcnt (Int64.of_int r) *) - | Many r when r = Z.zero -> 0 - | Many r -> Z.popcount r - +| One x -> One (One.plus_one x) +| Many x -> Many (Many.plus_one x) let trailing_zeros = function - | One r -> assert false (* Util.trailz (Int64.of_int r) *) - | Many r -> Z.trailing_zeros r +| One x -> One.trailing_zeros x +| Many x -> Many.trailing_zeros x + +let hamdist a b = match a, b with +| One a, One b -> One.hamdist a b +| Many a, Many b -> Many.hamdist a b +| _ -> invalid_arg "Bitstring.hamdist" + +let popcount = function +| One x -> One.popcount x +| Many x -> Many.popcount x + +let pp ppf = function +| One x -> One.pp ppf x +| Many x -> Many.pp ppf x -let hamdist a b = - match a,b with - | One a, One b -> assert false (* a lxor b |> Int64.of_int |> Util.popcnt *) - | One a, Many b -> Z.hamdist (Z.of_int a) b - | Many a, One b -> Z.hamdist a (Z.of_int b) - | Many a, Many b -> Z.hamdist a b - - -let rec to_list accu = function - | t when (t = Many Z.zero || t = One 0) -> List.rev accu +let rec to_list ?(accu=[]) = function + | t when (is_zero t) -> List.rev accu | t -> let newlist = (trailing_zeros t + 1)::accu in - to_list newlist (logand t (minus_one t)) + to_list ~accu:newlist (logand t (minus_one t)) (** [permtutations m n] generates the list of all possible [n]-bit @@ -149,118 +189,46 @@ let permtutations m n = let t'' = shift_right (minus_one (logand (lognot t) t')) (trailing_zeros u + 1) in aux (k-1) (logor t' t'') (u :: rest) in - aux (Util.binom n m) (minus_one (shift_left_one m)) [] - - - -let pp_bitstring n ppf s = - Format.fprintf ppf "@[%a@]" (Util.pp_bitstring n) ( - match s with - | Many b -> b - | One b -> Z.of_int b ) - - + aux (Util.binom n m) (minus_one (shift_left_one m m)) [] (*-----------------------------------------------------------------------------------*) -(* TODO let test_case () = - let test_creation () = - let l_a = [ 1 ; 2 ; 3 ; 5 ] in - let det = of_list l_a in - Alcotest.(check (list int )) "bitstring 1" l_a (to_list det); - Alcotest.(check bool) "phase 2" true (phase det = Phase.Pos); - let l_b = [ 1 ; 3 ; 2 ; 5 ] in - let det = of_list l_b in - Alcotest.(check (list int )) "bitstring 2" l_a (to_list det); - Alcotest.(check bool) "phase 2" true (phase det = Phase.Neg); + let test_one_many () = + let x = 8745687 in + let z = Z.of_int 8745687 in + let one_x = One x in + let many_x = Many z in + Alcotest.(check bool) "of_x" true (one_x = (of_int x)); + Alcotest.(check bool) "of_z" true (one_x = (of_z z)); + Alcotest.(check bool) "shift_left1" true (One (x lsl 3) = shift_left one_x 3); + Alcotest.(check bool) "shift_left2" true (Many (Z.shift_left z 3) = shift_left many_x 3); + Alcotest.(check bool) "shift_left3" true (Many (Z.shift_left z 100) = shift_left many_x 100); + Alcotest.(check bool) "shift_right1" true (One (x lsr 3) = shift_right one_x 3); + Alcotest.(check bool) "shift_right2" true (Many (Z.shift_right z 3) = shift_right many_x 3); + Alcotest.(check bool) "shift_left_one1" true (One (1 lsl 3) = shift_left_one 4 3); + Alcotest.(check bool) "shift_left_one2" true (Many (Z.shift_left Z.one 200) = shift_left_one 300 200); + Alcotest.(check bool) "testbit1" true (testbit (One 8) 3); + Alcotest.(check bool) "testbit2" false (testbit (One 8) 2); + Alcotest.(check bool) "testbit3" false (testbit (One 8) 4); + Alcotest.(check bool) "testbit4" true (testbit (Many (Z.of_int 8)) 3); + Alcotest.(check bool) "testbit5" false (testbit (Many (Z.of_int 8)) 2); + Alcotest.(check bool) "testbit6" false (testbit (Many (Z.of_int 8)) 4); + Alcotest.(check bool) "logor1" true (One (1 lor 2) = logor (One 1) (One 2)); + Alcotest.(check bool) "logor2" true (Many (Z.of_int (1 lor 2)) = logor (Many Z.one) (Many (Z.of_int 2))); + Alcotest.(check bool) "logxor1" true (One (1 lxor 2) = logxor (One 1) (One 2)); + Alcotest.(check bool) "logxor2" true (Many (Z.of_int (1 lxor 2)) = logxor (Many Z.one) (Many (Z.of_int 2))); + Alcotest.(check bool) "logand1" true (One (1 land 2) = logand (One 1) (One 2)); + Alcotest.(check bool) "logand2" true (Many (Z.of_int (1 land 2)) = logand (Many Z.one) (Many (Z.of_int 2))); + Alcotest.(check bool) "to_list" true ([ 1 ; 3 ; 4 ; 6 ] = (to_list (One 45))); in - let test_a_operators () = - let det = - creation 5 @@ creation 2 @@ creation 2 @@ creation 1 @@ vac - in - Alcotest.(check bool) "none 1" true (is_none det); - - let det = - creation 5 @@ creation 3 @@ creation 2 @@ creation 1 @@ vac - in - let l_a = [ 1 ; 2 ; 3 ; 5 ] in - Alcotest.(check (list int )) "bitstring 1" l_a (to_list det); - Alcotest.(check bool) "phase 1" true (phase det = Phase.Pos); - - let det = - creation 1 @@ creation 3 @@ creation 2 @@ creation 5 @@ vac - in - Alcotest.(check (list int )) "bitstring 2" l_a (to_list det); - Alcotest.(check bool) "phase 2" true (phase det = Phase.Neg); - - let l_b = [ 1 ; 3 ; 2 ; 5 ] in - let det = of_list l_b in - Alcotest.(check (list int )) "bitstring 3" l_a (to_list det); - Alcotest.(check bool) "phase 3" true (phase det = Phase.Neg); - - Alcotest.(check bool) "none 1" true (annihilation 4 det |> is_none); - - let det = - annihilation 1 det - in - Alcotest.(check (list int )) "bitstring 4" (List.tl l_a) (to_list det); - Alcotest.(check bool) "phase 4" true (phase det = Phase.Neg); - - let det = - annihilation 3 det - in - Alcotest.(check (list int )) "bitstring 5" [ 2 ; 5 ] (to_list det); - Alcotest.(check bool) "phase 5" true (phase det = Phase.Pos); - - let det = - annihilation 5 @@ annihilation 2 det - in - Alcotest.(check (list int )) "bitstring 6" [] (to_list det); - Alcotest.(check bool) "phase 6" true (phase det = Phase.Pos); - - in - - let test_exc_operators () = - let l_a = [ 1 ; 2 ; 3 ; 5 ] in - let det = of_list l_a in - let l_b = [ 1 ; 7 ; 3 ; 5 ] in - let det2 = of_list l_b in - Alcotest.(check bool) "single 1" true (single_excitation_reference 2 7 det = det2); - Alcotest.(check bool) "single 2" true (single_excitation 2 7 det = single_excitation_reference 2 7 det); - Alcotest.(check bool) "single 3" true (single_excitation_reference 4 7 det |> is_none); - Alcotest.(check bool) "single 4" true (single_excitation 4 7 det |> is_none); - - let l_c = [ 1 ; 7 ; 6 ; 5 ] in - let det3 = of_list l_c in - Alcotest.(check bool) "double 1" true (double_excitation_reference 2 7 3 6 det = det3); - Alcotest.(check bool) "double 2" true (double_excitation 2 7 3 6 det = double_excitation_reference 2 7 3 6 det); - Alcotest.(check bool) "double 3" true (double_excitation_reference 4 7 3 6 det |> is_none); - Alcotest.(check bool) "double 4" true (double_excitation 4 7 3 6 det |> is_none); - in - - let test_exc_spindet () = - let l_a = [ 1 ; 2 ; 3 ; 5 ] in - let det = of_list l_a in - let l_b = [ 1 ; 7 ; 3 ; 5 ] in - let det2 = of_list l_b in - Alcotest.(check int) "single" 1 (degree det det2); - Alcotest.(check (list int)) "holes" [2] (holes_of det det2); - Alcotest.(check (list int)) "particles" [7] (particles_of det det2); - let l_b = [ 1 ; 7 ; 3 ; 6 ] in - let det2 = of_list l_b in - Alcotest.(check int) "double" 2 (degree det det2); - Alcotest.(check (list int)) "holes" [2 ; 5] (holes_of det det2); - Alcotest.(check (list int)) "particles" [6 ; 7] (particles_of det det2); - in [ - "Creation", `Quick, test_creation; - "Creation/Annihilation Operators", `Quick, test_a_operators; - "Excitation Operators", `Quick, test_exc_operators; - "Excitation of spindet", `Quick, test_exc_spindet; + "One-many", `Quick, test_one_many; ] -*) + + + diff --git a/Utils/Util.mli b/Utils/Util.mli index c716ec4..53212d7 100644 --- a/Utils/Util.mli +++ b/Utils/Util.mli @@ -32,6 +32,9 @@ val popcnt : int64 -> int val trailz : int64 -> int (** ctz instruction *) +val leadz : int64 -> int +(** ctz instruction *) + (** {2 General functions} *) diff --git a/run_tests.ml b/run_tests.ml index 9e66476..5bdeebb 100644 --- a/run_tests.ml +++ b/run_tests.ml @@ -14,6 +14,7 @@ let test_water_dz () = in Alcotest.run "Unit tests" [ "Util", Util.test_case (); + "Bitstring", Bitstring.test_case (); "Spindeterminant", Spindeterminant.test_case (); "Determinant", Determinant.test_case (); "Excitation", Excitation.test_case ();