1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
# 1 "kernel/uint63_63.ml"
(************************************************************************)
(*         *   The Coq Proof Assistant / The Coq Development Team       *)
(*  v      *   INRIA, CNRS and contributors - Copyright 1999-2019       *)
(* <O___,, *       (see CREDITS file for the list of authors)           *)
(*   \VV/  **************************************************************)
(*    //   *    This file is distributed under the terms of the         *)
(*         *     GNU Lesser General Public License Version 2.1          *)
(*         *     (see LICENSE file for the text of the license)         *)
(************************************************************************)

type t = int

let _ = assert (Sys.word_size = 64)

let uint_size = 63

let maxuint63 = Int64.of_string "0x7FFFFFFFFFFFFFFF"
let maxuint31 = 0x7FFFFFFF

    (* conversion from an int *)
let to_uint64 i = Int64.logand (Int64.of_int i) maxuint63

let of_int i = i
[@@ocaml.inline always]

let to_int2 i = (0,i)

let of_int64 _i = assert false

let hash i = i
[@@ocaml.inline always]

    (* conversion of an uint63 to a string *)
let to_string i = Int64.to_string (to_uint64 i)

let of_string s =
  let i64 = Int64.of_string s in
  if Int64.compare Int64.zero i64 <= 0
      && Int64.compare i64 maxuint63 <= 0
  then Int64.to_int i64
  else raise (Failure "Int64.of_string")

(* Compiles an unsigned int to OCaml code *)
let compile i = Printf.sprintf "Uint63.of_int (%i)" i

let zero = 0
let one = 1

    (* logical shift *)
let l_sl x y =
  if 0 <= y && y < 63 then x lsl y else 0

let l_sr x y =
  if 0 <= y && y < 63 then x lsr y else 0

let l_and x y = x land y
[@@ocaml.inline always]

let l_or x y = x lor y
[@@ocaml.inline always]

let l_xor x y = x lxor y
[@@ocaml.inline always]

    (* addition of int63 *)
let add x y = x + y
[@@ocaml.inline always]

    (* subtraction *)
let sub x y = x - y
[@@ocaml.inline always]

    (* multiplication *)
let mul x y = x * y
[@@ocaml.inline always]

    (* division *)
let div (x : int) (y : int) =
  if y = 0 then 0 else Int64.to_int (Int64.div (to_uint64 x) (to_uint64 y))

    (* modulo *)
let rem (x : int) (y : int) =
  if y = 0 then 0 else Int64.to_int (Int64.rem (to_uint64 x) (to_uint64 y))

let diveucl x y = (div x y, rem x y)

let addmuldiv p x y =
  l_or (l_sl x p) (l_sr y (uint_size - p))

    (* comparison *)
let lt (x : int) (y : int) =
  (x lxor 0x4000000000000000) < (y lxor 0x4000000000000000)
[@@ocaml.inline always]

let le (x : int) (y : int) =
  (x lxor 0x4000000000000000) <= (y lxor 0x4000000000000000)
[@@ocaml.inline always]

    (* division of two numbers by one *)
(* precondition: xh < y *)
(* outputs: q, r s.t. x = q * y + r, r < y *)
let div21 xh xl y =
  (* nh might temporarily grow as large as 2*y - 1 in the loop body,
     so we store it as a 64-bit unsigned integer *)
  let nh = ref xh in
  let nl = ref xl in
  let q = ref 0 in
  for _i = 0 to 62 do
    (* invariants: 0 <= nh < y, nl = (xl*2^i) % 2^63,
       (q*y + nh) * 2^(63-i) + (xl % 2^(63-i)) = (xh%y) * 2^63 + xl *)
    nh := Int64.logor (Int64.shift_left !nh 1) (Int64.of_int (!nl lsr 62));
    nl := !nl lsl 1;
    q := !q lsl 1;
    (* TODO: use "Int64.unsigned_compare !nh y >= 0",
       once OCaml 4.08 becomes the minimal required version *)
    if Int64.compare !nh 0L < 0 || Int64.compare !nh y >= 0 then
      begin q := !q lor 1; nh := Int64.sub !nh y; end
  done;
  !q, Int64.to_int !nh

let div21 xh xl y =
  let xh = to_uint64 xh in
  let y = to_uint64 y in
  if Int64.compare y xh <= 0 then 0, 0 else div21 xh xl y

     (* exact multiplication *)
(* TODO: check that none of these additions could be a logical or *)


(* size = 32 + 31
   (hx << 31 + lx) * (hy << 31 + ly) =
   hxhy << 62 + (hxly + lxhy) << 31 + lxly
*)

(* precondition : (x lsr 62 = 0 || y lsr 62 = 0) *)
let mulc_aux x y =
  let lx = x land maxuint31 in
  let ly = y land maxuint31 in
  let hx = x lsr 31  in
  let hy = y lsr 31 in
    (* hx and hy are 32 bits value but at most one is 32 *)
  let hxy  = hx * hy in (* 63 bits *)
  let hxly = hx * ly in (* 63 bits *)
  let lxhy = lx * hy in (* 63 bits *)
  let lxy  = lx * ly in (* 62 bits *)
  let l  = lxy lor (hxy lsl 62) (* 63 bits *) in
  let h  = hxy lsr 1 in (* 62 bits *)
  let hl = hxly + lxhy in (* We can have a carry *)
  let h  = if lt hl hxly then h + (1 lsl 31) else h in
  let hl'= hl lsl 31 in
  let l  = l + hl' in
  let h  = if lt l hl' then h + 1 else h in
  let h  = h + (hl lsr 32) in
  (h,l)

let mulc x y =
  if (x lsr 62 = 0 || y lsr 62 = 0) then mulc_aux x y
  else
    let yl = y lxor (1 lsl 62) in
    let (h,l) = mulc_aux x yl in
    (* h << 63 + l = x * yl
       x * y = x * (1 << 62 + yl)  =
       x << 62 + x*yl = x << 62 + h << 63 + l *)
    let l' = l + (x lsl 62) in
    let h = if lt l' l then h + 1 else h in
    (h + (x lsr 1), l')

let equal (x : int) (y : int) = x = y
[@@ocaml.inline always]

let compare (x:int) (y:int) =
  let x = x lxor 0x4000000000000000 in
  let y = y lxor 0x4000000000000000 in
  if x > y then 1
  else if y > x then -1
  else 0

    (* head tail *)

let head0 x =
  let r = ref 0 in
  let x = ref x in
  if !x land 0x7FFFFFFF00000000 = 0 then r := !r + 31
  else x := !x lsr 31;
  if !x land 0xFFFF0000 = 0 then (x := !x lsl 16; r := !r + 16);
  if !x land 0xFF000000 = 0 then (x := !x lsl 8; r := !r + 8);
  if !x land 0xF0000000 = 0 then (x := !x lsl 4; r := !r + 4);
  if !x land 0xC0000000 = 0 then (x := !x lsl 2; r := !r + 2);
  if !x land 0x80000000 = 0 then (x := !x lsl 1; r := !r + 1);
  if !x land 0x80000000 = 0 then (               r := !r + 1);
  !r;;

let tail0 x =
  let r = ref 0 in
  let x = ref x in
  if !x land 0xFFFFFFFF = 0 then (x := !x lsr 32; r := !r + 32);
  if !x land 0xFFFF = 0 then (x := !x lsr 16; r := !r + 16);
  if !x land 0xFF = 0   then (x := !x lsr 8;  r := !r + 8);
  if !x land 0xF = 0    then (x := !x lsr 4;  r := !r + 4);
  if !x land 0x3 = 0    then (x := !x lsr 2;  r := !r + 2);
  if !x land 0x1 = 0    then (                r := !r + 1);
  !r

let is_uint63 t =
  Obj.is_int t
[@@ocaml.inline always]

(* Arithmetic with explicit carries *)

(* Analog of Numbers.Abstract.Cyclic.carry *)
type 'a carry = C0 of 'a | C1 of 'a

let addc x y =
  let r = x + y in
  if lt r x then C1 r else C0 r

let addcarryc x y =
  let r = x + y + 1 in
  if le r x then C1 r else C0 r

let subc x y =
  let r = x - y in
  if le y x then C0 r else C1 r

let subcarryc x y =
  let r = x - y - 1 in
  if lt y x then C0 r else C1 r