Branch data Line data Source code
1 : : /*
2 : : * Copyright (c) 2017 Andrew Kelley
3 : : *
4 : : * This file is part of zig, which is MIT licensed.
5 : : * See http://opensource.org/licenses/MIT
6 : : */
7 : :
8 : : #include "bigfloat.hpp"
9 : : #include "bigint.hpp"
10 : : #include "buffer.hpp"
11 : : #include "list.hpp"
12 : : #include "os.hpp"
13 : : #include "softfloat.hpp"
14 : :
15 : : #include <limits>
16 : : #include <algorithm>
17 : :
18 : : static uint64_t bigint_as_unsigned(const BigInt *bigint);
19 : :
20 : 867849 : static void bigint_normalize(BigInt *dest) {
21 : 867849 : const uint64_t *digits = bigint_ptr(dest);
22 : :
23 : 867849 : size_t last_nonzero_digit = SIZE_MAX;
24 [ + + ]: 1991448 : for (size_t i = 0; i < dest->digit_count; i += 1) {
25 : 1123599 : uint64_t digit = digits[i];
26 [ + + ]: 1123599 : if (digit != 0) {
27 : 953657 : last_nonzero_digit = i;
28 : : }
29 : : }
30 [ + + ]: 867849 : if (last_nonzero_digit == SIZE_MAX) {
31 : 2737 : dest->is_negative = false;
32 : 2737 : dest->digit_count = 0;
33 : : } else {
34 : 865112 : dest->digit_count = last_nonzero_digit + 1;
35 [ + + ]: 865112 : if (last_nonzero_digit == 0) {
36 : 742589 : dest->data.digit = digits[0];
37 : : }
38 : : }
39 : 867849 : }
40 : :
41 : 0 : static uint8_t digit_to_char(uint8_t digit, bool uppercase) {
42 [ # # ]: 0 : if (digit <= 9) {
43 : 0 : return digit + '0';
44 [ # # ]: 0 : } else if (digit <= 35) {
45 [ # # ]: 0 : return (digit - 10) + (uppercase ? 'A' : 'a');
46 : : } else {
47 : 0 : zig_unreachable();
48 : : }
49 : : }
50 : :
51 : 1163 : size_t bigint_bits_needed(const BigInt *op) {
52 : 1163 : size_t full_bits = op->digit_count * 64;
53 : 1163 : size_t leading_zero_count = bigint_clz(op, full_bits);
54 : 1163 : size_t bits_needed = full_bits - leading_zero_count;
55 : 1163 : return bits_needed + op->is_negative;
56 : : }
57 : :
58 : 18156 : static void to_twos_complement(BigInt *dest, const BigInt *op, size_t bit_count) {
59 [ + + ][ + + ]: 18156 : if (bit_count == 0 || op->digit_count == 0) {
60 : 172 : bigint_init_unsigned(dest, 0);
61 : 172 : return;
62 : : }
63 [ + + ]: 17984 : if (op->is_negative) {
64 : 376 : BigInt negated = {0};
65 : 376 : bigint_negate(&negated, op);
66 : :
67 : 376 : BigInt inverted = {0};
68 : 376 : bigint_not(&inverted, &negated, bit_count, false);
69 : :
70 : 376 : BigInt one = {0};
71 : 376 : bigint_init_unsigned(&one, 1);
72 : :
73 : 376 : bigint_add(dest, &inverted, &one);
74 : 376 : return;
75 : : }
76 : :
77 : 17608 : dest->is_negative = false;
78 : 17608 : const uint64_t *op_digits = bigint_ptr(op);
79 [ + + ]: 17608 : if (op->digit_count == 1) {
80 : 17160 : dest->data.digit = op_digits[0];
81 [ + + ]: 17160 : if (bit_count < 64) {
82 : 4424 : dest->data.digit &= (1ULL << bit_count) - 1;
83 : : }
84 : 17160 : dest->digit_count = 1;
85 : 17160 : bigint_normalize(dest);
86 : 17160 : return;
87 : : }
88 : 448 : size_t digits_to_copy = bit_count / 64;
89 : 448 : size_t leftover_bits = bit_count % 64;
90 : 448 : dest->digit_count = digits_to_copy + ((leftover_bits == 0) ? 0 : 1);
91 [ + + ][ + - ]: 448 : if (dest->digit_count == 1 && leftover_bits == 0) {
92 : 12 : dest->data.digit = op_digits[0];
93 [ - + ]: 12 : if (dest->data.digit == 0) dest->digit_count = 0;
94 : 12 : return;
95 : : }
96 : 436 : dest->data.digits = allocate_nonzero<uint64_t>(dest->digit_count);
97 [ + + ]: 1308 : for (size_t i = 0; i < digits_to_copy; i += 1) {
98 [ + - ]: 872 : uint64_t digit = (i < op->digit_count) ? op_digits[i] : 0;
99 : 872 : dest->data.digits[i] = digit;
100 : : }
101 [ - + ]: 436 : if (leftover_bits != 0) {
102 [ # # ]: 0 : uint64_t digit = (digits_to_copy < op->digit_count) ? op_digits[digits_to_copy] : 0;
103 : 0 : dest->data.digits[digits_to_copy] = digit & ((1ULL << leftover_bits) - 1);
104 : : }
105 : 436 : bigint_normalize(dest);
106 : : }
107 : :
108 : 19358126 : static bool bit_at_index(const BigInt *bi, size_t index) {
109 : 19358126 : size_t digit_index = index / 64;
110 [ + + ]: 19358126 : if (digit_index >= bi->digit_count)
111 : 1296 : return false;
112 : 19356830 : size_t digit_bit_index = index % 64;
113 : 19356830 : const uint64_t *digits = bigint_ptr(bi);
114 : 19356830 : uint64_t digit = digits[digit_index];
115 : 19356830 : return ((digit >> digit_bit_index) & 0x1) == 0x1;
116 : : }
117 : :
118 : 8562 : static void from_twos_complement(BigInt *dest, const BigInt *src, size_t bit_count, bool is_signed) {
119 : 8562 : assert(!src->is_negative);
120 : :
121 [ + + ][ + + ]: 8562 : if (bit_count == 0 || src->digit_count == 0) {
122 : 64 : bigint_init_unsigned(dest, 0);
123 : 64 : return;
124 : : }
125 : :
126 [ + + ][ + + ]: 8498 : if (is_signed && bit_at_index(src, bit_count - 1)) {
[ + + ]
127 : 548 : BigInt negative_one = {0};
128 : 548 : bigint_init_signed(&negative_one, -1);
129 : :
130 : 548 : BigInt minus_one = {0};
131 : 548 : bigint_add(&minus_one, src, &negative_one);
132 : :
133 : 548 : BigInt inverted = {0};
134 : 548 : bigint_not(&inverted, &minus_one, bit_count, false);
135 : :
136 : 548 : bigint_negate(dest, &inverted);
137 : 548 : return;
138 : :
139 : : }
140 : :
141 : 7950 : bigint_init_bigint(dest, src);
142 : : }
143 : :
144 : 2410399 : void bigint_init_unsigned(BigInt *dest, uint64_t x) {
145 [ + + ]: 2410399 : if (x == 0) {
146 : 747765 : dest->digit_count = 0;
147 : 747765 : dest->is_negative = false;
148 : 747765 : return;
149 : : }
150 : 1662634 : dest->digit_count = 1;
151 : 1662634 : dest->data.digit = x;
152 : 1662634 : dest->is_negative = false;
153 : : }
154 : :
155 : 580 : void bigint_init_signed(BigInt *dest, int64_t x) {
156 [ - + ]: 580 : if (x >= 0) {
157 : 0 : return bigint_init_unsigned(dest, x);
158 : : }
159 : 580 : dest->is_negative = true;
160 : 580 : dest->digit_count = 1;
161 : 580 : dest->data.digit = ((uint64_t)(-(x + 1))) + 1;
162 : : }
163 : :
164 : 0 : void bigint_init_data(BigInt *dest, const uint64_t *digits, size_t digit_count, bool is_negative) {
165 [ # # ]: 0 : if (digit_count == 0) {
166 : 0 : return bigint_init_unsigned(dest, 0);
167 [ # # ]: 0 : } else if (digit_count == 1) {
168 : 0 : dest->digit_count = 1;
169 : 0 : dest->data.digit = digits[0];
170 : 0 : dest->is_negative = is_negative;
171 : 0 : bigint_normalize(dest);
172 : 0 : return;
173 : : }
174 : :
175 : 0 : dest->digit_count = digit_count;
176 : 0 : dest->is_negative = is_negative;
177 : 0 : dest->data.digits = allocate_nonzero<uint64_t>(digit_count);
178 : 0 : memcpy(dest->data.digits, digits, sizeof(uint64_t) * digit_count);
179 : :
180 : 0 : bigint_normalize(dest);
181 : : }
182 : :
183 : 1195597 : void bigint_init_bigint(BigInt *dest, const BigInt *src) {
184 [ + + ]: 1195597 : if (src->digit_count == 0) {
185 : 184427 : return bigint_init_unsigned(dest, 0);
186 [ + + ]: 1011170 : } else if (src->digit_count == 1) {
187 : 979340 : dest->digit_count = 1;
188 : 979340 : dest->data.digit = src->data.digit;
189 : 979340 : dest->is_negative = src->is_negative;
190 : 979340 : return;
191 : : }
192 : 31830 : dest->is_negative = src->is_negative;
193 : 31830 : dest->digit_count = src->digit_count;
194 : 31830 : dest->data.digits = allocate_nonzero<uint64_t>(dest->digit_count);
195 : 31830 : memcpy(dest->data.digits, src->data.digits, sizeof(uint64_t) * dest->digit_count);
196 : : }
197 : :
198 : 48 : void bigint_init_bigfloat(BigInt *dest, const BigFloat *op) {
199 : : float128_t zero;
200 : 48 : ui32_to_f128M(0, &zero);
201 : :
202 : 48 : dest->is_negative = f128M_lt(&op->value, &zero);
203 : : float128_t abs_val;
204 [ - + ]: 48 : if (dest->is_negative) {
205 : 0 : f128M_sub(&zero, &op->value, &abs_val);
206 : : } else {
207 : 48 : memcpy(&abs_val, &op->value, sizeof(float128_t));
208 : : }
209 : :
210 : : float128_t max_u64;
211 : 48 : ui64_to_f128M(UINT64_MAX, &max_u64);
212 [ + - ]: 48 : if (f128M_le(&abs_val, &max_u64)) {
213 : 48 : dest->digit_count = 1;
214 : 48 : dest->data.digit = f128M_to_ui64(&op->value, softfloat_round_minMag, false);
215 : 48 : bigint_normalize(dest);
216 : 48 : return;
217 : : }
218 : :
219 : : float128_t amt;
220 : 0 : f128M_div(&abs_val, &max_u64, &amt);
221 : : float128_t remainder;
222 : 0 : f128M_rem(&abs_val, &max_u64, &remainder);
223 : :
224 : 0 : dest->digit_count = 2;
225 : 0 : dest->data.digits = allocate_nonzero<uint64_t>(dest->digit_count);
226 : 0 : dest->data.digits[0] = f128M_to_ui64(&remainder, softfloat_round_minMag, false);
227 : 0 : dest->data.digits[1] = f128M_to_ui64(&amt, softfloat_round_minMag, false);
228 : 0 : bigint_normalize(dest);
229 : : }
230 : :
231 : 452444 : bool bigint_fits_in_bits(const BigInt *bn, size_t bit_count, bool is_signed) {
232 [ + + ][ + - ]: 452444 : assert(bn->digit_count != 1 || bn->data.digit != 0);
233 [ + + ]: 452444 : if (bit_count == 0) {
234 : 472 : return bigint_cmp_zero(bn) == CmpEQ;
235 : : }
236 [ + + ]: 451972 : if (bn->digit_count == 0) {
237 : 38810 : return true;
238 : : }
239 : :
240 [ + + ]: 413162 : if (!is_signed) {
241 : 361456 : size_t full_bits = bn->digit_count * 64;
242 : 361456 : size_t leading_zero_count = bigint_clz(bn, full_bits);
243 : 361456 : return bit_count >= full_bits - leading_zero_count;
244 : : }
245 : :
246 : 51706 : BigInt one = {0};
247 : 51706 : bigint_init_unsigned(&one, 1);
248 : :
249 : 51706 : BigInt shl_amt = {0};
250 : 51706 : bigint_init_unsigned(&shl_amt, bit_count - 1);
251 : :
252 : 51706 : BigInt max_value_plus_one = {0};
253 : 51706 : bigint_shl(&max_value_plus_one, &one, &shl_amt);
254 : :
255 : 51706 : BigInt max_value = {0};
256 : 51706 : bigint_sub(&max_value, &max_value_plus_one, &one);
257 : :
258 : 51706 : BigInt min_value = {0};
259 : 51706 : bigint_negate(&min_value, &max_value_plus_one);
260 : :
261 : 51706 : Cmp min_cmp = bigint_cmp(bn, &min_value);
262 : 51706 : Cmp max_cmp = bigint_cmp(bn, &max_value);
263 : :
264 [ + - ][ + + ]: 452444 : return (min_cmp == CmpGT || min_cmp == CmpEQ) && (max_cmp == CmpLT || max_cmp == CmpEQ);
[ + - ][ + + ]
265 : : }
266 : :
267 : 11038 : void bigint_write_twos_complement(const BigInt *big_int, uint8_t *buf, size_t bit_count, bool is_big_endian) {
268 [ - + ]: 11038 : if (bit_count == 0)
269 : 0 : return;
270 : :
271 : 11038 : BigInt twos_comp = {0};
272 : 11038 : to_twos_complement(&twos_comp, big_int, bit_count);
273 : :
274 : 11038 : const uint64_t *twos_comp_digits = bigint_ptr(&twos_comp);
275 : :
276 : 11038 : size_t bits_in_last_digit = bit_count % 64;
277 [ + + ]: 11038 : if (bits_in_last_digit == 0) bits_in_last_digit = 64;
278 : 11038 : size_t bytes_in_last_digit = (bits_in_last_digit + 7) / 8;
279 : 11038 : size_t unwritten_byte_count = 8 - bytes_in_last_digit;
280 : :
281 [ + + ]: 11038 : if (is_big_endian) {
282 : 256 : size_t last_digit_index = (bit_count - 1) / 64;
283 : 256 : size_t digit_index = last_digit_index;
284 : 256 : size_t buf_index = 0;
285 : : for (;;) {
286 [ + - ]: 288 : uint64_t x = (digit_index < twos_comp.digit_count) ? twos_comp_digits[digit_index] : 0;
287 : :
288 : 288 : for (size_t byte_index = 7;;) {
289 : 1632 : uint8_t byte = x & 0xff;
290 [ + + ]: 1632 : if (digit_index == last_digit_index) {
291 : 1376 : buf[buf_index + byte_index - unwritten_byte_count] = byte;
292 [ + + ]: 1376 : if (byte_index == unwritten_byte_count) break;
293 : : } else {
294 : 256 : buf[buf_index + byte_index] = byte;
295 : : }
296 : :
297 [ + + ]: 1376 : if (byte_index == 0) break;
298 : 1344 : byte_index -= 1;
299 : 1344 : x >>= 8;
300 : 1344 : }
301 : :
302 [ + + ]: 288 : if (digit_index == 0) break;
303 : 32 : digit_index -= 1;
304 [ - + ]: 32 : if (digit_index == last_digit_index) {
305 : 0 : buf_index += bytes_in_last_digit;
306 : : } else {
307 : 32 : buf_index += 8;
308 : : }
309 : 288 : }
310 : : } else {
311 : 10782 : size_t digit_count = (bit_count + 63) / 64;
312 : 10782 : size_t buf_index = 0;
313 [ + + ]: 22084 : for (size_t digit_index = 0; digit_index < digit_count; digit_index += 1) {
314 [ + + ]: 11046 : uint64_t x = (digit_index < twos_comp.digit_count) ? twos_comp_digits[digit_index] : 0;
315 : :
316 : 89742 : for (size_t byte_index = 0;
317 [ + + ][ + + ]: 89742 : byte_index < 8 && (digit_index + 1 < digit_count || byte_index < bytes_in_last_digit);
[ + + ]
318 : 78696 : byte_index += 1)
319 : : {
320 : 78696 : uint8_t byte = x & 0xff;
321 : 78696 : buf[buf_index] = byte;
322 : 78696 : buf_index += 1;
323 : 78696 : x >>= 8;
324 : : }
325 : : }
326 : : }
327 : : }
328 : :
329 : :
330 : 13234 : void bigint_read_twos_complement(BigInt *dest, const uint8_t *buf, size_t bit_count, bool is_big_endian,
331 : : bool is_signed)
332 : : {
333 [ - + ]: 13234 : if (bit_count == 0) {
334 : 0 : bigint_init_unsigned(dest, 0);
335 : 0 : return;
336 : : }
337 : :
338 : 13234 : dest->digit_count = (bit_count + 63) / 64;
339 : : uint64_t *digits;
340 [ + + ]: 13234 : if (dest->digit_count == 1) {
341 : 12914 : digits = &dest->data.digit;
342 : : } else {
343 : 320 : digits = allocate_nonzero<uint64_t>(dest->digit_count);
344 : 320 : dest->data.digits = digits;
345 : : }
346 : :
347 : 13234 : size_t bits_in_last_digit = bit_count % 64;
348 [ + + ]: 13234 : if (bits_in_last_digit == 0) {
349 : 11198 : bits_in_last_digit = 64;
350 : : }
351 : 13234 : size_t bytes_in_last_digit = (bits_in_last_digit + 7) / 8;
352 : 13234 : size_t unread_byte_count = 8 - bytes_in_last_digit;
353 : :
354 [ - + ]: 13234 : if (is_big_endian) {
355 : 0 : size_t buf_index = 0;
356 : 0 : uint64_t digit = 0;
357 [ # # ]: 0 : for (size_t byte_index = unread_byte_count; byte_index < 8; byte_index += 1) {
358 : 0 : uint8_t byte = buf[buf_index];
359 : 0 : buf_index += 1;
360 : 0 : digit <<= 8;
361 : 0 : digit |= byte;
362 : : }
363 : 0 : digits[dest->digit_count - 1] = digit;
364 [ # # ]: 0 : for (size_t digit_index = 1; digit_index < dest->digit_count; digit_index += 1) {
365 : 0 : digit = 0;
366 [ # # ]: 0 : for (size_t byte_index = 0; byte_index < 8; byte_index += 1) {
367 : 0 : uint8_t byte = buf[buf_index];
368 : 0 : buf_index += 1;
369 : 0 : digit <<= 8;
370 : 0 : digit |= byte;
371 : : }
372 : 0 : digits[dest->digit_count - 1 - digit_index] = digit;
373 : : }
374 : : } else {
375 : 13234 : size_t buf_index = 0;
376 [ + + ]: 26788 : for (size_t digit_index = 0; digit_index < dest->digit_count; digit_index += 1) {
377 : 13554 : uint64_t digit = 0;
378 [ + + ]: 13554 : size_t end_byte_index = (digit_index == dest->digit_count - 1) ? bytes_in_last_digit : 8;
379 [ + + ]: 112618 : for (size_t byte_index = 0; byte_index < end_byte_index; byte_index += 1) {
380 : 99064 : uint64_t byte = buf[buf_index];
381 : 99064 : buf_index += 1;
382 : :
383 : 99064 : digit |= byte << (8 * byte_index);
384 : : }
385 : 13554 : digits[digit_index] = digit;
386 : : }
387 : : }
388 : :
389 [ + + ]: 13234 : if (is_signed) {
390 : 1548 : bigint_normalize(dest);
391 : 1548 : BigInt tmp = {0};
392 : 1548 : bigint_init_bigint(&tmp, dest);
393 : 1548 : from_twos_complement(dest, &tmp, bit_count, true);
394 : : } else {
395 : 11686 : dest->is_negative = false;
396 : 11686 : bigint_normalize(dest);
397 : : }
398 : : }
399 : :
400 : : #if defined(_MSC_VER)
401 : : static bool add_u64_overflow(uint64_t op1, uint64_t op2, uint64_t *result) {
402 : : *result = op1 + op2;
403 : : return *result < op1 || *result < op2;
404 : : }
405 : :
406 : : static bool sub_u64_overflow(uint64_t op1, uint64_t op2, uint64_t *result) {
407 : : *result = op1 - op2;
408 : : return *result > op1;
409 : : }
410 : :
411 : : bool mul_u64_overflow(uint64_t op1, uint64_t op2, uint64_t *result) {
412 : : *result = op1 * op2;
413 : :
414 : : if (op1 == 0 || op2 == 0)
415 : : return false;
416 : :
417 : : if (op1 > UINT64_MAX / op2)
418 : : return true;
419 : :
420 : : if (op2 > UINT64_MAX / op1)
421 : : return true;
422 : :
423 : : return false;
424 : : }
425 : : #else
426 : 390084 : static bool add_u64_overflow(uint64_t op1, uint64_t op2, uint64_t *result) {
427 : 390084 : return __builtin_uaddll_overflow((unsigned long long)op1, (unsigned long long)op2,
428 : 390084 : (unsigned long long *)result);
429 : : }
430 : :
431 : 62831 : static bool sub_u64_overflow(uint64_t op1, uint64_t op2, uint64_t *result) {
432 : 62831 : return __builtin_usubll_overflow((unsigned long long)op1, (unsigned long long)op2,
433 : 62831 : (unsigned long long *)result);
434 : : }
435 : :
436 : 162 : bool mul_u64_overflow(uint64_t op1, uint64_t op2, uint64_t *result) {
437 : 162 : return __builtin_umulll_overflow((unsigned long long)op1, (unsigned long long)op2,
438 : 162 : (unsigned long long *)result);
439 : : }
440 : : #endif
441 : :
442 : 524701 : void bigint_add(BigInt *dest, const BigInt *op1, const BigInt *op2) {
443 [ + + ]: 524701 : if (op1->digit_count == 0) {
444 : 524197 : return bigint_init_bigint(dest, op2);
445 : : }
446 [ + + ]: 428869 : if (op2->digit_count == 0) {
447 : 70709 : return bigint_init_bigint(dest, op1);
448 : : }
449 [ + + ]: 358160 : if (op1->is_negative == op2->is_negative) {
450 : 294653 : dest->is_negative = op1->is_negative;
451 : :
452 : 294653 : const uint64_t *op1_digits = bigint_ptr(op1);
453 : 294653 : const uint64_t *op2_digits = bigint_ptr(op2);
454 : 294653 : bool overflow = add_u64_overflow(op1_digits[0], op2_digits[0], &dest->data.digit);
455 [ + + ][ + - ]: 294653 : if (overflow == 0 && op1->digit_count == 1 && op2->digit_count == 1) {
[ + + ]
456 : 237463 : dest->digit_count = 1;
457 : 237463 : bigint_normalize(dest);
458 : 237463 : return;
459 : : }
460 : 57190 : size_t i = 1;
461 : 57190 : uint64_t first_digit = dest->data.digit;
462 : 57190 : dest->data.digits = allocate_nonzero<uint64_t>(max(op1->digit_count, op2->digit_count) + 1);
463 : 57190 : dest->data.digits[0] = first_digit;
464 : :
465 : : for (;;) {
466 : 133268 : bool found_digit = false;
467 : 133268 : uint64_t x = overflow;
468 : 133268 : overflow = 0;
469 : :
470 [ + + ]: 133268 : if (i < op1->digit_count) {
471 : 76046 : found_digit = true;
472 : 76046 : uint64_t digit = op1_digits[i];
473 : 76046 : overflow += add_u64_overflow(x, digit, &x);
474 : : }
475 : :
476 [ + + ]: 133268 : if (i < op2->digit_count) {
477 : 19385 : found_digit = true;
478 : 19385 : uint64_t digit = op2_digits[i];
479 : 19385 : overflow += add_u64_overflow(x, digit, &x);
480 : : }
481 : :
482 : 133268 : dest->data.digits[i] = x;
483 : 133268 : i += 1;
484 : :
485 [ + + ]: 133268 : if (!found_digit) {
486 : 57190 : dest->digit_count = i;
487 : 57190 : bigint_normalize(dest);
488 : 57190 : return;
489 : : }
490 : 76078 : }
491 : : }
492 : : const BigInt *op_pos;
493 : : const BigInt *op_neg;
494 [ + + ]: 63507 : if (op1->is_negative) {
495 : 2056 : op_neg = op1;
496 : 2056 : op_pos = op2;
497 : : } else {
498 : 61451 : op_pos = op1;
499 : 61451 : op_neg = op2;
500 : : }
501 : :
502 : 63507 : BigInt op_neg_abs = {0};
503 : 63507 : bigint_negate(&op_neg_abs, op_neg);
504 : : const BigInt *bigger_op;
505 : : const BigInt *smaller_op;
506 [ + + + - ]: 63507 : switch (bigint_cmp(op_pos, &op_neg_abs)) {
507 : 1192 : case CmpEQ:
508 : 1192 : bigint_init_unsigned(dest, 0);
509 : 1192 : return;
510 : 1784 : case CmpLT:
511 : 1784 : bigger_op = &op_neg_abs;
512 : 1784 : smaller_op = op_pos;
513 : 1784 : dest->is_negative = true;
514 : 1784 : break;
515 : 60531 : case CmpGT:
516 : 60531 : bigger_op = op_pos;
517 : 60531 : smaller_op = &op_neg_abs;
518 : 60531 : dest->is_negative = false;
519 : 60531 : break;
520 : : }
521 : 62315 : const uint64_t *bigger_op_digits = bigint_ptr(bigger_op);
522 : 62315 : const uint64_t *smaller_op_digits = bigint_ptr(smaller_op);
523 : 62315 : uint64_t overflow = sub_u64_overflow(bigger_op_digits[0], smaller_op_digits[0], &dest->data.digit);
524 [ + + ][ + - ]: 62315 : if (overflow == 0 && bigger_op->digit_count == 1 && smaller_op->digit_count == 1) {
[ + + ]
525 : 61811 : dest->digit_count = 1;
526 : 61811 : bigint_normalize(dest);
527 : 61811 : return;
528 : : }
529 : 504 : uint64_t first_digit = dest->data.digit;
530 : 504 : dest->data.digits = allocate_nonzero<uint64_t>(bigger_op->digit_count);
531 : 504 : dest->data.digits[0] = first_digit;
532 : 504 : size_t i = 1;
533 : :
534 : : for (;;) {
535 : 516 : bool found_digit = false;
536 : 516 : uint64_t x = bigger_op_digits[i];
537 : 516 : uint64_t prev_overflow = overflow;
538 : 516 : overflow = 0;
539 : :
540 [ - + ]: 516 : if (i < smaller_op->digit_count) {
541 : 0 : found_digit = true;
542 : 0 : uint64_t digit = smaller_op_digits[i];
543 : 0 : overflow += sub_u64_overflow(x, digit, &x);
544 : : }
545 [ + + ]: 516 : if (sub_u64_overflow(x, prev_overflow, &x)) {
546 : 12 : found_digit = true;
547 : 12 : overflow += 1;
548 : : }
549 : 516 : dest->data.digits[i] = x;
550 : 516 : i += 1;
551 : :
552 [ + + ][ + - ]: 516 : if (!found_digit || i >= bigger_op->digit_count)
553 : : break;
554 : 12 : }
555 : 504 : assert(overflow == 0);
556 : 504 : dest->digit_count = i;
557 : 504 : bigint_normalize(dest);
558 : : }
559 : :
560 : 2184 : void bigint_add_wrap(BigInt *dest, const BigInt *op1, const BigInt *op2, size_t bit_count, bool is_signed) {
561 : 2184 : BigInt unwrapped = {0};
562 : 2184 : bigint_add(&unwrapped, op1, op2);
563 : 2184 : bigint_truncate(dest, &unwrapped, bit_count, is_signed);
564 : 2184 : }
565 : :
566 : 62417 : void bigint_sub(BigInt *dest, const BigInt *op1, const BigInt *op2) {
567 : 62417 : BigInt op2_negated = {0};
568 : 62417 : bigint_negate(&op2_negated, op2);
569 : 62417 : return bigint_add(dest, op1, &op2_negated);
570 : : }
571 : :
572 : 88 : void bigint_sub_wrap(BigInt *dest, const BigInt *op1, const BigInt *op2, size_t bit_count, bool is_signed) {
573 : 88 : BigInt op2_negated = {0};
574 : 88 : bigint_negate(&op2_negated, op2);
575 : 88 : return bigint_add_wrap(dest, op1, &op2_negated, bit_count, is_signed);
576 : : }
577 : :
578 : 271532 : static void mul_overflow(uint64_t op1, uint64_t op2, uint64_t *lo, uint64_t *hi) {
579 : 271532 : uint64_t u1 = (op1 & 0xffffffff);
580 : 271532 : uint64_t v1 = (op2 & 0xffffffff);
581 : 271532 : uint64_t t = (u1 * v1);
582 : 271532 : uint64_t w3 = (t & 0xffffffff);
583 : 271532 : uint64_t k = (t >> 32);
584 : :
585 : 271532 : op1 >>= 32;
586 : 271532 : t = (op1 * v1) + k;
587 : 271532 : k = (t & 0xffffffff);
588 : 271532 : uint64_t w1 = (t >> 32);
589 : :
590 : 271532 : op2 >>= 32;
591 : 271532 : t = (u1 * op2) + k;
592 : 271532 : k = (t >> 32);
593 : :
594 : 271532 : *hi = (op1 * op2) + w1 + k;
595 : 271532 : *lo = (t << 32) + w3;
596 : 271532 : }
597 : :
598 : 19159 : static void mul_scalar(BigInt *dest, const BigInt *op, uint64_t scalar) {
599 : 19159 : bigint_init_unsigned(dest, 0);
600 : :
601 : : BigInt bi_64;
602 : 19159 : bigint_init_unsigned(&bi_64, 64);
603 : :
604 : 19159 : const uint64_t *op_digits = bigint_ptr(op);
605 : 19159 : size_t i = op->digit_count - 1;
606 : :
607 : : for (;;) {
608 : : BigInt shifted;
609 : 42142 : bigint_shl(&shifted, dest, &bi_64);
610 : :
611 : : uint64_t result_scalar;
612 : : uint64_t carry_scalar;
613 : 42142 : mul_overflow(scalar, op_digits[i], &result_scalar, &carry_scalar);
614 : :
615 : : BigInt result;
616 : 42142 : bigint_init_unsigned(&result, result_scalar);
617 : :
618 : : BigInt carry;
619 : 42142 : bigint_init_unsigned(&carry, carry_scalar);
620 : :
621 : : BigInt carry_shifted;
622 : 42142 : bigint_shl(&carry_shifted, &carry, &bi_64);
623 : :
624 : : BigInt tmp;
625 : 42142 : bigint_add(&tmp, &shifted, &carry_shifted);
626 : :
627 : 42142 : bigint_add(dest, &tmp, &result);
628 : :
629 [ + + ]: 42142 : if (i == 0) {
630 : 19159 : break;
631 : : }
632 : 22983 : i -= 1;
633 : 22983 : }
634 : 19159 : }
635 : :
636 : 281111 : void bigint_mul(BigInt *dest, const BigInt *op1, const BigInt *op2) {
637 [ + + ][ - + ]: 281111 : if (op1->digit_count == 0 || op2->digit_count == 0) {
638 : 261984 : return bigint_init_unsigned(dest, 0);
639 : : }
640 : 229390 : const uint64_t *op1_digits = bigint_ptr(op1);
641 : 229390 : const uint64_t *op2_digits = bigint_ptr(op2);
642 : :
643 : : uint64_t carry;
644 : 229390 : mul_overflow(op1_digits[0], op2_digits[0], &dest->data.digit, &carry);
645 [ + + ][ + - ]: 229390 : if (carry == 0 && op1->digit_count == 1 && op2->digit_count == 1) {
[ + + ]
646 : 210263 : dest->is_negative = (op1->is_negative != op2->is_negative);
647 : 210263 : dest->digit_count = 1;
648 : 210263 : bigint_normalize(dest);
649 : 210263 : return;
650 : : }
651 : :
652 : 19127 : bigint_init_unsigned(dest, 0);
653 : :
654 : : BigInt bi_64;
655 : 19127 : bigint_init_unsigned(&bi_64, 64);
656 : :
657 : 19127 : size_t i = op2->digit_count - 1;
658 : : for (;;) {
659 : : BigInt shifted;
660 : 19159 : bigint_shl(&shifted, dest, &bi_64);
661 : :
662 : : BigInt scalar_result;
663 : 19159 : mul_scalar(&scalar_result, op1, op2_digits[i]);
664 : :
665 : 19159 : bigint_add(dest, &scalar_result, &shifted);
666 : :
667 [ + + ]: 19159 : if (i == 0) {
668 : 19127 : break;
669 : : }
670 : 32 : i -= 1;
671 : 32 : }
672 : :
673 : 19127 : dest->is_negative = (op1->is_negative != op2->is_negative);
674 : 19127 : bigint_normalize(dest);
675 : : }
676 : :
677 : 32 : void bigint_mul_wrap(BigInt *dest, const BigInt *op1, const BigInt *op2, size_t bit_count, bool is_signed) {
678 : 32 : BigInt unwrapped = {0};
679 : 32 : bigint_mul(&unwrapped, op1, op2);
680 : 32 : bigint_truncate(dest, &unwrapped, bit_count, is_signed);
681 : 32 : }
682 : :
683 : : enum ZeroBehavior {
684 : : /// \brief The returned value is undefined.
685 : : ZB_Undefined,
686 : : /// \brief The returned value is numeric_limits<T>::max()
687 : : ZB_Max,
688 : : /// \brief The returned value is numeric_limits<T>::digits
689 : : ZB_Width
690 : : };
691 : :
692 : : template <typename T, std::size_t SizeOfT> struct LeadingZerosCounter {
693 : : static std::size_t count(T Val, ZeroBehavior) {
694 : : if (!Val)
695 : : return std::numeric_limits<T>::digits;
696 : :
697 : : // Bisection method.
698 : : std::size_t ZeroBits = 0;
699 : : for (T Shift = std::numeric_limits<T>::digits >> 1; Shift; Shift >>= 1) {
700 : : T Tmp = Val >> Shift;
701 : : if (Tmp)
702 : : Val = Tmp;
703 : : else
704 : : ZeroBits |= Shift;
705 : : }
706 : : return ZeroBits;
707 : : }
708 : : };
709 : :
710 : : #if __GNUC__ >= 4 || defined(_MSC_VER)
711 : : template <typename T> struct LeadingZerosCounter<T, 4> {
712 : 96 : static std::size_t count(T Val, ZeroBehavior ZB) {
713 [ + - ][ - + ]: 96 : if (ZB != ZB_Undefined && Val == 0)
714 : 0 : return 32;
715 : :
716 : : #if defined(_MSC_VER)
717 : : unsigned long Index;
718 : : _BitScanReverse(&Index, Val);
719 : : return Index ^ 31;
720 : : #else
721 : 96 : return __builtin_clz(Val);
722 : : #endif
723 : : }
724 : : };
725 : :
726 : : #if !defined(_MSC_VER) || defined(_M_X64)
727 : : template <typename T> struct LeadingZerosCounter<T, 8> {
728 : : static std::size_t count(T Val, ZeroBehavior ZB) {
729 : : if (ZB != ZB_Undefined && Val == 0)
730 : : return 64;
731 : :
732 : : #if defined(_MSC_VER)
733 : : unsigned long Index;
734 : : _BitScanReverse64(&Index, Val);
735 : : return Index ^ 63;
736 : : #else
737 : : return __builtin_clzll(Val);
738 : : #endif
739 : : }
740 : : };
741 : : #endif
742 : : #endif
743 : :
744 : : /// \brief Count number of 0's from the most significant bit to the least
745 : : /// stopping at the first 1.
746 : : ///
747 : : /// Only unsigned integral types are allowed.
748 : : ///
749 : : /// \param ZB the behavior on an input of 0. Only ZB_Width and ZB_Undefined are
750 : : /// valid arguments.
751 : : template <typename T>
752 : 96 : std::size_t countLeadingZeros(T Val, ZeroBehavior ZB = ZB_Width) {
753 : : static_assert(std::numeric_limits<T>::is_integer &&
754 : : !std::numeric_limits<T>::is_signed,
755 : : "Only unsigned integral types are allowed.");
756 : 96 : return LeadingZerosCounter<T, sizeof(T)>::count(Val, ZB);
757 : : }
758 : :
759 : : /// Make a 64-bit integer from a high / low pair of 32-bit integers.
760 : 496 : constexpr inline uint64_t Make_64(uint32_t High, uint32_t Low) {
761 : 496 : return ((uint64_t)High << 32) | (uint64_t)Low;
762 : : }
763 : :
764 : : /// Return the high 32 bits of a 64 bit value.
765 : 1600 : constexpr inline uint32_t Hi_32(uint64_t Value) {
766 : 1600 : return static_cast<uint32_t>(Value >> 32);
767 : : }
768 : :
769 : : /// Return the low 32 bits of a 64 bit value.
770 : 1984 : constexpr inline uint32_t Lo_32(uint64_t Value) {
771 : 1984 : return static_cast<uint32_t>(Value);
772 : : }
773 : :
774 : : /// Implementation of Knuth's Algorithm D (Division of nonnegative integers)
775 : : /// from "Art of Computer Programming, Volume 2", section 4.3.1, p. 272. The
776 : : /// variables here have the same names as in the algorithm. Comments explain
777 : : /// the algorithm and any deviation from it.
778 : 96 : static void KnuthDiv(uint32_t *u, uint32_t *v, uint32_t *q, uint32_t* r,
779 : : unsigned m, unsigned n)
780 : : {
781 : 96 : assert(u && "Must provide dividend");
782 : 96 : assert(v && "Must provide divisor");
783 : 96 : assert(q && "Must provide quotient");
784 [ + - ][ + - ]: 96 : assert(u != v && u != q && v != q && "Must use different memory");
[ + - ]
785 : 96 : assert(n>1 && "n must be > 1");
786 : :
787 : : // b denotes the base of the number system. In our case b is 2^32.
788 : 96 : const uint64_t b = uint64_t(1) << 32;
789 : :
790 : : // D1. [Normalize.] Set d = b / (v[n-1] + 1) and multiply all the digits of
791 : : // u and v by d. Note that we have taken Knuth's advice here to use a power
792 : : // of 2 value for d such that d * v[n-1] >= b/2 (b is the base). A power of
793 : : // 2 allows us to shift instead of multiply and it is easy to determine the
794 : : // shift amount from the leading zeros. We are basically normalizing the u
795 : : // and v so that its high bits are shifted to the top of v's range without
796 : : // overflow. Note that this can require an extra word in u so that u must
797 : : // be of length m+n+1.
798 : 96 : unsigned shift = countLeadingZeros(v[n-1]);
799 : 96 : uint32_t v_carry = 0;
800 : 96 : uint32_t u_carry = 0;
801 [ + - ]: 96 : if (shift) {
802 [ + + ]: 576 : for (unsigned i = 0; i < m+n; ++i) {
803 : 480 : uint32_t u_tmp = u[i] >> (32 - shift);
804 : 480 : u[i] = (u[i] << shift) | u_carry;
805 : 480 : u_carry = u_tmp;
806 : : }
807 [ + + ]: 576 : for (unsigned i = 0; i < n; ++i) {
808 : 480 : uint32_t v_tmp = v[i] >> (32 - shift);
809 : 480 : v[i] = (v[i] << shift) | v_carry;
810 : 480 : v_carry = v_tmp;
811 : : }
812 : : }
813 : 96 : u[m+n] = u_carry;
814 : :
815 : : // D2. [Initialize j.] Set j to m. This is the loop counter over the places.
816 : 96 : int j = m;
817 : 0 : do {
818 : : // D3. [Calculate q'.].
819 : : // Set qp = (u[j+n]*b + u[j+n-1]) / v[n-1]. (qp=qprime=q')
820 : : // Set rp = (u[j+n]*b + u[j+n-1]) % v[n-1]. (rp=rprime=r')
821 : : // Now test if qp == b or qp*v[n-2] > b*rp + u[j+n-2]; if so, decrease
822 : : // qp by 1, increase rp by v[n-1], and repeat this test if rp < b. The test
823 : : // on v[n-2] determines at high speed most of the cases in which the trial
824 : : // value qp is one too large, and it eliminates all cases where qp is two
825 : : // too large.
826 : 96 : uint64_t dividend = Make_64(u[j+n], u[j+n-1]);
827 : 96 : uint64_t qp = dividend / v[n-1];
828 : 96 : uint64_t rp = dividend % v[n-1];
829 [ - + ][ + - ]: 96 : if (qp == b || qp*v[n-2] > b*rp + u[j+n-2]) {
830 : 0 : qp--;
831 : 0 : rp += v[n-1];
832 [ # # ][ # # ]: 0 : if (rp < b && (qp == b || qp*v[n-2] > b*rp + u[j+n-2]))
[ # # ]
833 : 0 : qp--;
834 : : }
835 : :
836 : : // D4. [Multiply and subtract.] Replace (u[j+n]u[j+n-1]...u[j]) with
837 : : // (u[j+n]u[j+n-1]..u[j]) - qp * (v[n-1]...v[1]v[0]). This computation
838 : : // consists of a simple multiplication by a one-place number, combined with
839 : : // a subtraction.
840 : : // The digits (u[j+n]...u[j]) should be kept positive; if the result of
841 : : // this step is actually negative, (u[j+n]...u[j]) should be left as the
842 : : // true value plus b**(n+1), namely as the b's complement of
843 : : // the true value, and a "borrow" to the left should be remembered.
844 : 96 : int64_t borrow = 0;
845 [ + + ]: 576 : for (unsigned i = 0; i < n; ++i) {
846 : 480 : uint64_t p = uint64_t(qp) * uint64_t(v[i]);
847 : 480 : int64_t subres = int64_t(u[j+i]) - borrow - Lo_32(p);
848 : 480 : u[j+i] = Lo_32(subres);
849 : 480 : borrow = Hi_32(p) - Hi_32(subres);
850 : : }
851 : 96 : bool isNeg = u[j+n] < borrow;
852 : 96 : u[j+n] -= Lo_32(borrow);
853 : :
854 : : // D5. [Test remainder.] Set q[j] = qp. If the result of step D4 was
855 : : // negative, go to step D6; otherwise go on to step D7.
856 : 96 : q[j] = Lo_32(qp);
857 [ - + ]: 96 : if (isNeg) {
858 : : // D6. [Add back]. The probability that this step is necessary is very
859 : : // small, on the order of only 2/b. Make sure that test data accounts for
860 : : // this possibility. Decrease q[j] by 1
861 : 0 : q[j]--;
862 : : // and add (0v[n-1]...v[1]v[0]) to (u[j+n]u[j+n-1]...u[j+1]u[j]).
863 : : // A carry will occur to the left of u[j+n], and it should be ignored
864 : : // since it cancels with the borrow that occurred in D4.
865 : 0 : bool carry = false;
866 [ # # ]: 0 : for (unsigned i = 0; i < n; i++) {
867 : 0 : uint32_t limit = std::min(u[j+i],v[i]);
868 : 0 : u[j+i] += v[i] + carry;
869 [ # # ][ # # ]: 0 : carry = u[j+i] < limit || (carry && u[j+i] == limit);
[ # # ]
870 : : }
871 : 0 : u[j+n] += carry;
872 : : }
873 : :
874 : : // D7. [Loop on j.] Decrease j by one. Now if j >= 0, go back to D3.
875 [ - + ]: 96 : } while (--j >= 0);
876 : :
877 : : // D8. [Unnormalize]. Now q[...] is the desired quotient, and the desired
878 : : // remainder may be obtained by dividing u[...] by d. If r is non-null we
879 : : // compute the remainder (urem uses this).
880 [ + + ]: 96 : if (r) {
881 : : // The value d is expressed by the "shift" value above since we avoided
882 : : // multiplication by d by using a shift left. So, all we have to do is
883 : : // shift right here.
884 [ + - ]: 32 : if (shift) {
885 : 32 : uint32_t carry = 0;
886 [ + + ]: 192 : for (int i = n-1; i >= 0; i--) {
887 : 160 : r[i] = (u[i] >> shift) | carry;
888 : 160 : carry = u[i] << (32 - shift);
889 : : }
890 : : } else {
891 [ # # ]: 32 : for (int i = n-1; i >= 0; i--) {
892 : 0 : r[i] = u[i];
893 : : }
894 : : }
895 : : }
896 : 96 : }
897 : :
898 : : // Implementation ported from LLVM/lib/Support/APInt.cpp
899 : 112 : static void bigint_unsigned_division(const BigInt *op1, const BigInt *op2, BigInt *Quotient, BigInt *Remainder) {
900 : 112 : Cmp cmp = bigint_cmp(op1, op2);
901 [ - + ]: 112 : if (cmp == CmpLT) {
902 [ # # ]: 0 : if (Quotient != nullptr) {
903 : 0 : bigint_init_unsigned(Quotient, 0);
904 : : }
905 [ # # ]: 0 : if (Remainder != nullptr) {
906 : 0 : bigint_init_bigint(Remainder, op1);
907 : : }
908 : 0 : return;
909 : : }
910 [ - + ]: 112 : if (cmp == CmpEQ) {
911 [ # # ]: 0 : if (Quotient != nullptr) {
912 : 0 : bigint_init_unsigned(Quotient, 1);
913 : : }
914 [ # # ]: 0 : if (Remainder != nullptr) {
915 : 0 : bigint_init_unsigned(Remainder, 0);
916 : : }
917 : 0 : return;
918 : : }
919 : :
920 : 112 : const uint64_t *LHS = bigint_ptr(op1);
921 : 112 : const uint64_t *RHS = bigint_ptr(op2);
922 : 112 : unsigned lhsWords = op1->digit_count;
923 : 112 : unsigned rhsWords = op2->digit_count;
924 : :
925 : : // First, compose the values into an array of 32-bit words instead of
926 : : // 64-bit words. This is a necessity of both the "short division" algorithm
927 : : // and the Knuth "classical algorithm" which requires there to be native
928 : : // operations for +, -, and * on an m bit value with an m*2 bit result. We
929 : : // can't use 64-bit operands here because we don't have native results of
930 : : // 128-bits. Furthermore, casting the 64-bit values to 32-bit values won't
931 : : // work on large-endian machines.
932 : 112 : unsigned n = rhsWords * 2;
933 : 112 : unsigned m = (lhsWords * 2) - n;
934 : :
935 : : // Allocate space for the temporary values we need either on the stack, if
936 : : // it will fit, or on the heap if it won't.
937 : : uint32_t SPACE[128];
938 : 112 : uint32_t *U = nullptr;
939 : 112 : uint32_t *V = nullptr;
940 : 112 : uint32_t *Q = nullptr;
941 : 112 : uint32_t *R = nullptr;
942 [ + - ][ + + ]: 112 : if ((Remainder?4:3)*n+2*m+1 <= 128) {
943 : 112 : U = &SPACE[0];
944 : 112 : V = &SPACE[m+n+1];
945 : 112 : Q = &SPACE[(m+n+1) + n];
946 [ + + ]: 112 : if (Remainder)
947 : 48 : R = &SPACE[(m+n+1) + n + (m+n)];
948 : : } else {
949 [ # # ]: 0 : U = new uint32_t[m + n + 1];
950 [ # # ]: 0 : V = new uint32_t[n];
951 [ # # ]: 0 : Q = new uint32_t[m+n];
952 [ # # ]: 0 : if (Remainder)
953 [ # # ]: 0 : R = new uint32_t[n];
954 : : }
955 : :
956 : : // Initialize the dividend
957 : 112 : memset(U, 0, (m+n+1)*sizeof(uint32_t));
958 [ + + ]: 448 : for (unsigned i = 0; i < lhsWords; ++i) {
959 : 336 : uint64_t tmp = LHS[i];
960 : 336 : U[i * 2] = Lo_32(tmp);
961 : 336 : U[i * 2 + 1] = Hi_32(tmp);
962 : : }
963 : 112 : U[m+n] = 0; // this extra word is for "spill" in the Knuth algorithm.
964 : :
965 : : // Initialize the divisor
966 : 112 : memset(V, 0, (n)*sizeof(uint32_t));
967 [ + + ]: 416 : for (unsigned i = 0; i < rhsWords; ++i) {
968 : 304 : uint64_t tmp = RHS[i];
969 : 304 : V[i * 2] = Lo_32(tmp);
970 : 304 : V[i * 2 + 1] = Hi_32(tmp);
971 : : }
972 : :
973 : : // initialize the quotient and remainder
974 : 112 : memset(Q, 0, (m+n) * sizeof(uint32_t));
975 [ + + ]: 112 : if (Remainder)
976 : 48 : memset(R, 0, n * sizeof(uint32_t));
977 : :
978 : : // Now, adjust m and n for the Knuth division. n is the number of words in
979 : : // the divisor. m is the number of words by which the dividend exceeds the
980 : : // divisor (i.e. m+n is the length of the dividend). These sizes must not
981 : : // contain any zero words or the Knuth algorithm fails.
982 [ + - ][ + + ]: 224 : for (unsigned i = n; i > 0 && V[i-1] == 0; i--) {
983 : 112 : n--;
984 : 112 : m++;
985 : : }
986 [ + - ][ + + ]: 208 : for (unsigned i = m+n; i > 0 && U[i-1] == 0; i--)
987 : 96 : m--;
988 : :
989 : : // If we're left with only a single word for the divisor, Knuth doesn't work
990 : : // so we implement the short division algorithm here. This is much simpler
991 : : // and faster because we are certain that we can divide a 64-bit quantity
992 : : // by a 32-bit quantity at hardware speed and short division is simply a
993 : : // series of such operations. This is just like doing short division but we
994 : : // are using base 2^32 instead of base 10.
995 : 112 : assert(n != 0 && "Divide by zero?");
996 [ + + ]: 112 : if (n == 1) {
997 : 16 : uint32_t divisor = V[0];
998 : 16 : uint32_t remainder = 0;
999 [ + + ]: 112 : for (int i = m; i >= 0; i--) {
1000 : 96 : uint64_t partial_dividend = Make_64(remainder, U[i]);
1001 [ - + ]: 96 : if (partial_dividend == 0) {
1002 : 0 : Q[i] = 0;
1003 : 0 : remainder = 0;
1004 [ - + ]: 96 : } else if (partial_dividend < divisor) {
1005 : 0 : Q[i] = 0;
1006 : 0 : remainder = Lo_32(partial_dividend);
1007 [ - + ]: 96 : } else if (partial_dividend == divisor) {
1008 : 0 : Q[i] = 1;
1009 : 0 : remainder = 0;
1010 : : } else {
1011 : 96 : Q[i] = Lo_32(partial_dividend / divisor);
1012 : 96 : remainder = Lo_32(partial_dividend - (Q[i] * divisor));
1013 : : }
1014 : : }
1015 [ + - ]: 16 : if (R)
1016 : 16 : R[0] = remainder;
1017 : : } else {
1018 : : // Now we're ready to invoke the Knuth classical divide algorithm. In this
1019 : : // case n > 1.
1020 : 96 : KnuthDiv(U, V, Q, R, m, n);
1021 : : }
1022 : :
1023 : : // If the caller wants the quotient
1024 [ + + ]: 112 : if (Quotient) {
1025 : 64 : Quotient->is_negative = false;
1026 : 64 : Quotient->digit_count = lhsWords;
1027 [ - + ]: 64 : if (lhsWords == 1) {
1028 : 0 : Quotient->data.digit = Make_64(Q[1], Q[0]);
1029 : : } else {
1030 : 64 : Quotient->data.digits = allocate<uint64_t>(lhsWords);
1031 [ + + ]: 256 : for (size_t i = 0; i < lhsWords; i += 1) {
1032 : 192 : Quotient->data.digits[i] = Make_64(Q[i*2+1], Q[i*2]);
1033 : : }
1034 : : }
1035 : : }
1036 : :
1037 : : // If the caller wants the remainder
1038 [ + + ]: 112 : if (Remainder) {
1039 : 48 : Remainder->is_negative = false;
1040 : 48 : Remainder->digit_count = rhsWords;
1041 [ + + ]: 48 : if (rhsWords == 1) {
1042 : 16 : Remainder->data.digit = Make_64(R[1], R[0]);
1043 : : } else {
1044 : 32 : Remainder->data.digits = allocate<uint64_t>(rhsWords);
1045 [ + + ]: 208 : for (size_t i = 0; i < rhsWords; i += 1) {
1046 : 96 : Remainder->data.digits[i] = Make_64(R[i*2+1], R[i*2]);
1047 : : }
1048 : : }
1049 : : }
1050 : : }
1051 : :
1052 : 602 : void bigint_div_trunc(BigInt *dest, const BigInt *op1, const BigInt *op2) {
1053 : 602 : assert(op2->digit_count != 0); // division by zero
1054 [ + + ]: 602 : if (op1->digit_count == 0) {
1055 : 8 : bigint_init_unsigned(dest, 0);
1056 : 538 : return;
1057 : : }
1058 : 594 : const uint64_t *op1_digits = bigint_ptr(op1);
1059 : 594 : const uint64_t *op2_digits = bigint_ptr(op2);
1060 [ + - ][ + + ]: 594 : if (op1->digit_count == 1 && op2->digit_count == 1) {
1061 : 530 : dest->data.digit = op1_digits[0] / op2_digits[0];
1062 : 530 : dest->digit_count = 1;
1063 : 530 : dest->is_negative = op1->is_negative != op2->is_negative;
1064 : 530 : bigint_normalize(dest);
1065 : 530 : return;
1066 : : }
1067 [ - + ][ # # ]: 64 : if (op2->digit_count == 1 && op2_digits[0] == 1) {
1068 : : // X / 1 == X
1069 : 0 : bigint_init_bigint(dest, op1);
1070 : 0 : dest->is_negative = op1->is_negative != op2->is_negative;
1071 : 0 : bigint_normalize(dest);
1072 : 0 : return;
1073 : : }
1074 : :
1075 : : const BigInt *op1_positive;
1076 : : BigInt op1_positive_data;
1077 [ + + ]: 64 : if (op1->is_negative) {
1078 : 32 : bigint_negate(&op1_positive_data, op1);
1079 : 32 : op1_positive = &op1_positive_data;
1080 : : } else {
1081 : 32 : op1_positive = op1;
1082 : : }
1083 : :
1084 : : const BigInt *op2_positive;
1085 : : BigInt op2_positive_data;
1086 [ + + ]: 64 : if (op2->is_negative) {
1087 : 32 : bigint_negate(&op2_positive_data, op2);
1088 : 32 : op2_positive = &op2_positive_data;
1089 : : } else {
1090 : 32 : op2_positive = op2;
1091 : : }
1092 : :
1093 : 64 : bigint_unsigned_division(op1_positive, op2_positive, dest, nullptr);
1094 : 64 : dest->is_negative = op1->is_negative != op2->is_negative;
1095 : 64 : bigint_normalize(dest);
1096 : : }
1097 : :
1098 : 80 : void bigint_div_floor(BigInt *dest, const BigInt *op1, const BigInt *op2) {
1099 [ + + ]: 80 : if (op1->is_negative != op2->is_negative) {
1100 : 56 : bigint_div_trunc(dest, op1, op2);
1101 : 56 : BigInt mult_again = {0};
1102 : 56 : bigint_mul(&mult_again, dest, op2);
1103 : 56 : mult_again.is_negative = op1->is_negative;
1104 [ + + ]: 56 : if (bigint_cmp(&mult_again, op1) != CmpEQ) {
1105 : 32 : BigInt tmp = {0};
1106 : 32 : bigint_init_bigint(&tmp, dest);
1107 : 32 : BigInt neg_one = {0};
1108 : 32 : bigint_init_signed(&neg_one, -1);
1109 : 32 : bigint_add(dest, &tmp, &neg_one);
1110 : : }
1111 : 56 : bigint_normalize(dest);
1112 : : } else {
1113 : 24 : bigint_div_trunc(dest, op1, op2);
1114 : : }
1115 : 80 : }
1116 : :
1117 : 408 : void bigint_rem(BigInt *dest, const BigInt *op1, const BigInt *op2) {
1118 : 408 : assert(op2->digit_count != 0); // division by zero
1119 [ - + ]: 408 : if (op1->digit_count == 0) {
1120 : 0 : bigint_init_unsigned(dest, 0);
1121 : 360 : return;
1122 : : }
1123 : 408 : const uint64_t *op1_digits = bigint_ptr(op1);
1124 : 408 : const uint64_t *op2_digits = bigint_ptr(op2);
1125 : :
1126 [ + - ][ + + ]: 408 : if (op1->digit_count == 1 && op2->digit_count == 1) {
1127 : 360 : dest->data.digit = op1_digits[0] % op2_digits[0];
1128 : 360 : dest->digit_count = 1;
1129 : 360 : dest->is_negative = op1->is_negative;
1130 : 360 : bigint_normalize(dest);
1131 : 360 : return;
1132 : : }
1133 [ - + ][ # # ]: 48 : if (op2->digit_count == 2 && op2_digits[0] == 0 && op2_digits[1] == 1) {
[ # # ]
1134 : : // special case this divisor
1135 : 0 : bigint_init_unsigned(dest, op1_digits[0]);
1136 : 0 : dest->is_negative = op1->is_negative;
1137 : 0 : bigint_normalize(dest);
1138 : 0 : return;
1139 : : }
1140 : :
1141 [ + + ][ - + ]: 48 : if (op2->digit_count == 1 && op2_digits[0] == 1) {
1142 : : // X % 1 == 0
1143 : 0 : bigint_init_unsigned(dest, 0);
1144 : 0 : return;
1145 : : }
1146 : :
1147 : : const BigInt *op1_positive;
1148 : : BigInt op1_positive_data;
1149 [ + + ]: 48 : if (op1->is_negative) {
1150 : 16 : bigint_negate(&op1_positive_data, op1);
1151 : 16 : op1_positive = &op1_positive_data;
1152 : : } else {
1153 : 32 : op1_positive = op1;
1154 : : }
1155 : :
1156 : : const BigInt *op2_positive;
1157 : : BigInt op2_positive_data;
1158 [ - + ]: 48 : if (op2->is_negative) {
1159 : 0 : bigint_negate(&op2_positive_data, op2);
1160 : 0 : op2_positive = &op2_positive_data;
1161 : : } else {
1162 : 48 : op2_positive = op2;
1163 : : }
1164 : :
1165 : 48 : bigint_unsigned_division(op1_positive, op2_positive, nullptr, dest);
1166 : 48 : dest->is_negative = op1->is_negative;
1167 : 48 : bigint_normalize(dest);
1168 : : }
1169 : :
1170 : 40 : void bigint_mod(BigInt *dest, const BigInt *op1, const BigInt *op2) {
1171 [ + + ]: 40 : if (op1->is_negative) {
1172 : : BigInt first_rem;
1173 : 24 : bigint_rem(&first_rem, op1, op2);
1174 : 24 : first_rem.is_negative = !op2->is_negative;
1175 : : BigInt op2_minus_rem;
1176 : 24 : bigint_add(&op2_minus_rem, op2, &first_rem);
1177 : 24 : bigint_rem(dest, &op2_minus_rem, op2);
1178 : 24 : dest->is_negative = false;
1179 : : } else {
1180 : 16 : bigint_rem(dest, op1, op2);
1181 : 16 : dest->is_negative = false;
1182 : : }
1183 : 40 : }
1184 : :
1185 : 3107 : void bigint_or(BigInt *dest, const BigInt *op1, const BigInt *op2) {
1186 [ + + ]: 3107 : if (op1->digit_count == 0) {
1187 : 475 : return bigint_init_bigint(dest, op2);
1188 : : }
1189 [ + + ]: 2632 : if (op2->digit_count == 0) {
1190 : 8 : return bigint_init_bigint(dest, op1);
1191 : : }
1192 [ + + ][ + + ]: 2624 : if (op1->is_negative || op2->is_negative) {
1193 : 16 : size_t big_bit_count = max(bigint_bits_needed(op1), bigint_bits_needed(op2));
1194 : :
1195 : 16 : BigInt twos_comp_op1 = {0};
1196 : 16 : to_twos_complement(&twos_comp_op1, op1, big_bit_count);
1197 : :
1198 : 16 : BigInt twos_comp_op2 = {0};
1199 : 16 : to_twos_complement(&twos_comp_op2, op2, big_bit_count);
1200 : :
1201 : 16 : BigInt twos_comp_dest = {0};
1202 : 16 : bigint_or(&twos_comp_dest, &twos_comp_op1, &twos_comp_op2);
1203 : :
1204 : 16 : from_twos_complement(dest, &twos_comp_dest, big_bit_count, true);
1205 : : } else {
1206 : 2608 : dest->is_negative = false;
1207 : 2608 : const uint64_t *op1_digits = bigint_ptr(op1);
1208 : 2608 : const uint64_t *op2_digits = bigint_ptr(op2);
1209 [ + + ][ + + ]: 2608 : if (op1->digit_count == 1 && op2->digit_count == 1) {
1210 : 2584 : dest->digit_count = 1;
1211 : 2584 : dest->data.digit = op1_digits[0] | op2_digits[0];
1212 : 2584 : bigint_normalize(dest);
1213 : 2584 : return;
1214 : : }
1215 : 24 : dest->digit_count = max(op1->digit_count, op2->digit_count);
1216 : 24 : dest->data.digits = allocate_nonzero<uint64_t>(dest->digit_count);
1217 [ + + ]: 72 : for (size_t i = 0; i < dest->digit_count; i += 1) {
1218 : 48 : uint64_t digit = 0;
1219 [ + + ]: 48 : if (i < op1->digit_count) {
1220 : 40 : digit |= op1_digits[i];
1221 : : }
1222 [ + - ]: 48 : if (i < op2->digit_count) {
1223 : 48 : digit |= op2_digits[i];
1224 : : }
1225 : 48 : dest->data.digits[i] = digit;
1226 : : }
1227 : 24 : bigint_normalize(dest);
1228 : : }
1229 : : }
1230 : :
1231 : 6280 : void bigint_and(BigInt *dest, const BigInt *op1, const BigInt *op2) {
1232 [ + - ][ - + ]: 6280 : if (op1->digit_count == 0 || op2->digit_count == 0) {
1233 : 0 : return bigint_init_unsigned(dest, 0);
1234 : : }
1235 [ + + ][ + + ]: 6280 : if (op1->is_negative || op2->is_negative) {
1236 : 24 : size_t big_bit_count = max(bigint_bits_needed(op1), bigint_bits_needed(op2));
1237 : :
1238 : 24 : BigInt twos_comp_op1 = {0};
1239 : 24 : to_twos_complement(&twos_comp_op1, op1, big_bit_count);
1240 : :
1241 : 24 : BigInt twos_comp_op2 = {0};
1242 : 24 : to_twos_complement(&twos_comp_op2, op2, big_bit_count);
1243 : :
1244 : 24 : BigInt twos_comp_dest = {0};
1245 : 24 : bigint_and(&twos_comp_dest, &twos_comp_op1, &twos_comp_op2);
1246 : :
1247 : 24 : from_twos_complement(dest, &twos_comp_dest, big_bit_count, true);
1248 : : } else {
1249 : 6256 : dest->is_negative = false;
1250 : 6256 : const uint64_t *op1_digits = bigint_ptr(op1);
1251 : 6256 : const uint64_t *op2_digits = bigint_ptr(op2);
1252 [ + + ][ + + ]: 6256 : if (op1->digit_count == 1 && op2->digit_count == 1) {
1253 : 6216 : dest->digit_count = 1;
1254 : 6216 : dest->data.digit = op1_digits[0] & op2_digits[0];
1255 : 6216 : bigint_normalize(dest);
1256 : 6216 : return;
1257 : : }
1258 : :
1259 : 40 : dest->digit_count = max(op1->digit_count, op2->digit_count);
1260 : 40 : dest->data.digits = allocate_nonzero<uint64_t>(dest->digit_count);
1261 : :
1262 : 40 : size_t i = 0;
1263 [ + + ][ + + ]: 96 : for (; i < op1->digit_count && i < op2->digit_count; i += 1) {
1264 : 56 : dest->data.digits[i] = op1_digits[i] & op2_digits[i];
1265 : : }
1266 [ + + ]: 64 : for (; i < dest->digit_count; i += 1) {
1267 : 24 : dest->data.digits[i] = 0;
1268 : : }
1269 : 40 : bigint_normalize(dest);
1270 : : }
1271 : : }
1272 : :
1273 : 276 : void bigint_xor(BigInt *dest, const BigInt *op1, const BigInt *op2) {
1274 [ - + ]: 276 : if (op1->digit_count == 0) {
1275 : 0 : return bigint_init_bigint(dest, op2);
1276 : : }
1277 [ + + ]: 276 : if (op2->digit_count == 0) {
1278 : 24 : return bigint_init_bigint(dest, op1);
1279 : : }
1280 [ + + ][ + + ]: 252 : if (op1->is_negative || op2->is_negative) {
1281 : 16 : size_t big_bit_count = max(bigint_bits_needed(op1), bigint_bits_needed(op2));
1282 : :
1283 : 16 : BigInt twos_comp_op1 = {0};
1284 : 16 : to_twos_complement(&twos_comp_op1, op1, big_bit_count);
1285 : :
1286 : 16 : BigInt twos_comp_op2 = {0};
1287 : 16 : to_twos_complement(&twos_comp_op2, op2, big_bit_count);
1288 : :
1289 : 16 : BigInt twos_comp_dest = {0};
1290 : 16 : bigint_xor(&twos_comp_dest, &twos_comp_op1, &twos_comp_op2);
1291 : :
1292 : 16 : from_twos_complement(dest, &twos_comp_dest, big_bit_count, true);
1293 : : } else {
1294 : 236 : dest->is_negative = false;
1295 : 236 : const uint64_t *op1_digits = bigint_ptr(op1);
1296 : 236 : const uint64_t *op2_digits = bigint_ptr(op2);
1297 : :
1298 [ + - ][ + - ]: 236 : assert(op1->digit_count > 0 && op2->digit_count > 0);
1299 [ + + ][ + + ]: 236 : if (op1->digit_count == 1 && op2->digit_count == 1) {
1300 : 152 : dest->digit_count = 1;
1301 : 152 : dest->data.digit = op1_digits[0] ^ op2_digits[0];
1302 : 152 : bigint_normalize(dest);
1303 : 152 : return;
1304 : : }
1305 : 84 : dest->digit_count = max(op1->digit_count, op2->digit_count);
1306 : 84 : dest->data.digits = allocate_nonzero<uint64_t>(dest->digit_count);
1307 : 84 : size_t i = 0;
1308 [ + + ][ + + ]: 228 : for (; i < op1->digit_count && i < op2->digit_count; i += 1) {
1309 : 144 : dest->data.digits[i] = op1_digits[i] ^ op2_digits[i];
1310 : : }
1311 [ + + ]: 108 : for (; i < dest->digit_count; i += 1) {
1312 [ + + ]: 24 : if (i < op1->digit_count) {
1313 : 8 : dest->data.digits[i] = op1_digits[i];
1314 [ + - ]: 16 : } else if (i < op2->digit_count) {
1315 : 16 : dest->data.digits[i] = op2_digits[i];
1316 : : } else {
1317 : 0 : zig_unreachable();
1318 : : }
1319 : : }
1320 : 84 : bigint_normalize(dest);
1321 : : }
1322 : : }
1323 : :
1324 : 161624 : void bigint_shl(BigInt *dest, const BigInt *op1, const BigInt *op2) {
1325 : 161624 : assert(!op2->is_negative);
1326 : :
1327 [ + + ]: 161624 : if (op2->digit_count == 0) {
1328 : 525 : bigint_init_bigint(dest, op1);
1329 : 525 : return;
1330 : : }
1331 : :
1332 [ + + ]: 161099 : if (op1->digit_count == 0) {
1333 : 59889 : bigint_init_unsigned(dest, 0);
1334 : 59889 : return;
1335 : : }
1336 : :
1337 [ - + ]: 101210 : if (op2->digit_count != 1) {
1338 : 0 : zig_panic("TODO shift left by amount greater than 64 bit integer");
1339 : : }
1340 : :
1341 : 101210 : const uint64_t *op1_digits = bigint_ptr(op1);
1342 : 101210 : uint64_t shift_amt = bigint_as_unsigned(op2);
1343 : :
1344 [ + + ][ + + ]: 101210 : if (op1->digit_count == 1 && shift_amt < 64) {
1345 : 57156 : dest->data.digit = op1_digits[0] << shift_amt;
1346 [ + + ]: 57156 : if (dest->data.digit > op1_digits[0]) {
1347 : 57144 : dest->digit_count = 1;
1348 : 57144 : dest->is_negative = op1->is_negative;
1349 : 57144 : return;
1350 : : }
1351 : : }
1352 : :
1353 : 44066 : uint64_t digit_shift_count = shift_amt / 64;
1354 : 44066 : uint64_t leftover_shift_count = shift_amt % 64;
1355 : :
1356 : 44066 : dest->data.digits = allocate<uint64_t>(op1->digit_count + digit_shift_count + 1);
1357 : 44066 : dest->digit_count = digit_shift_count;
1358 : 44066 : uint64_t carry = 0;
1359 [ + + ]: 95714 : for (size_t i = 0; i < op1->digit_count; i += 1) {
1360 : 51648 : uint64_t digit = op1_digits[i];
1361 : 51648 : dest->data.digits[dest->digit_count] = carry | (digit << leftover_shift_count);
1362 : 51648 : dest->digit_count += 1;
1363 [ + + ]: 51648 : if (leftover_shift_count > 0) {
1364 : 492 : carry = digit >> (64 - leftover_shift_count);
1365 : : } else {
1366 : 51156 : carry = 0;
1367 : : }
1368 : : }
1369 : 44066 : dest->data.digits[dest->digit_count] = carry;
1370 : 44066 : dest->digit_count += 1;
1371 : 44066 : dest->is_negative = op1->is_negative;
1372 : 44066 : bigint_normalize(dest);
1373 : : }
1374 : :
1375 : 4564 : void bigint_shl_trunc(BigInt *dest, const BigInt *op1, const BigInt *op2, size_t bit_count, bool is_signed) {
1376 : 4564 : BigInt unwrapped = {0};
1377 : 4564 : bigint_shl(&unwrapped, op1, op2);
1378 : 4564 : bigint_truncate(dest, &unwrapped, bit_count, is_signed);
1379 : 4564 : }
1380 : :
1381 : 10927 : void bigint_shr(BigInt *dest, const BigInt *op1, const BigInt *op2) {
1382 : 10927 : assert(!op2->is_negative);
1383 : :
1384 [ + + ]: 10927 : if (op1->digit_count == 0) {
1385 : 8 : return bigint_init_unsigned(dest, 0);
1386 : : }
1387 : :
1388 [ - + ]: 10919 : if (op2->digit_count == 0) {
1389 : 0 : return bigint_init_bigint(dest, op1);
1390 : : }
1391 : :
1392 [ - + ]: 10919 : if (op2->digit_count != 1) {
1393 : 0 : zig_panic("TODO shift right by amount greater than 64 bit integer");
1394 : : }
1395 : :
1396 : 10919 : const uint64_t *op1_digits = bigint_ptr(op1);
1397 : 10919 : uint64_t shift_amt = bigint_as_unsigned(op2);
1398 : :
1399 [ + + ]: 10919 : if (op1->digit_count == 1) {
1400 [ + + ]: 10871 : dest->data.digit = (shift_amt < 64) ? op1_digits[0] >> shift_amt : 0;
1401 : 10871 : dest->digit_count = 1;
1402 : 10871 : dest->is_negative = op1->is_negative;
1403 : 10871 : bigint_normalize(dest);
1404 : 10871 : return;
1405 : : }
1406 : :
1407 : 48 : size_t digit_shift_count = shift_amt / 64;
1408 : 48 : size_t leftover_shift_count = shift_amt % 64;
1409 : :
1410 [ - + ]: 48 : if (digit_shift_count >= op1->digit_count) {
1411 : 0 : return bigint_init_unsigned(dest, 0);
1412 : : }
1413 : :
1414 : 48 : dest->digit_count = op1->digit_count - digit_shift_count;
1415 : : uint64_t *digits;
1416 [ + + ]: 48 : if (dest->digit_count == 1) {
1417 : 16 : digits = &dest->data.digit;
1418 : : } else {
1419 : 32 : digits = allocate<uint64_t>(dest->digit_count);
1420 : 32 : dest->data.digits = digits;
1421 : : }
1422 : :
1423 : 48 : uint64_t carry = 0;
1424 : 48 : for (size_t op_digit_index = op1->digit_count - 1;;) {
1425 : 80 : uint64_t digit = op1_digits[op_digit_index];
1426 : 80 : size_t dest_digit_index = op_digit_index - digit_shift_count;
1427 : 80 : digits[dest_digit_index] = carry | (digit >> leftover_shift_count);
1428 : 80 : carry = digit << (64 - leftover_shift_count);
1429 : :
1430 [ + + ]: 80 : if (dest_digit_index == 0) { break; }
1431 : 32 : op_digit_index -= 1;
1432 : 32 : }
1433 : 48 : dest->is_negative = op1->is_negative;
1434 : 48 : bigint_normalize(dest);
1435 : : }
1436 : :
1437 : 184422 : void bigint_negate(BigInt *dest, const BigInt *op) {
1438 : 184422 : bigint_init_bigint(dest, op);
1439 : 184422 : dest->is_negative = !dest->is_negative;
1440 : 184422 : bigint_normalize(dest);
1441 : 184422 : }
1442 : :
1443 : 40 : void bigint_negate_wrap(BigInt *dest, const BigInt *op, size_t bit_count) {
1444 : : BigInt zero;
1445 : 40 : bigint_init_unsigned(&zero, 0);
1446 : 40 : bigint_sub_wrap(dest, &zero, op, bit_count, true);
1447 : 40 : }
1448 : :
1449 : 1074 : void bigint_not(BigInt *dest, const BigInt *op, size_t bit_count, bool is_signed) {
1450 [ - + ]: 1074 : if (bit_count == 0) {
1451 : 0 : bigint_init_unsigned(dest, 0);
1452 : 0 : return;
1453 : : }
1454 : :
1455 [ + + ]: 1074 : if (is_signed) {
1456 : 26 : BigInt twos_comp = {0};
1457 : 26 : to_twos_complement(&twos_comp, op, bit_count);
1458 : :
1459 : 26 : BigInt inverted = {0};
1460 : 26 : bigint_not(&inverted, &twos_comp, bit_count, false);
1461 : :
1462 : 26 : from_twos_complement(dest, &inverted, bit_count, true);
1463 : 26 : return;
1464 : : }
1465 : :
1466 : 1048 : assert(!op->is_negative);
1467 : :
1468 : 1048 : dest->is_negative = false;
1469 : 1048 : const uint64_t *op_digits = bigint_ptr(op);
1470 [ + + ]: 1048 : if (bit_count <= 64) {
1471 : 772 : dest->digit_count = 1;
1472 [ + + ]: 772 : if (op->digit_count == 0) {
1473 [ + + ]: 48 : if (bit_count == 64) {
1474 : 32 : dest->data.digit = UINT64_MAX;
1475 : : } else {
1476 : 48 : dest->data.digit = (1ULL << bit_count) - 1;
1477 : : }
1478 [ + - ]: 724 : } else if (op->digit_count == 1) {
1479 : 724 : dest->data.digit = ~op_digits[0];
1480 [ + + ]: 724 : if (bit_count != 64) {
1481 : 570 : uint64_t mask = (1ULL << bit_count) - 1;
1482 : 724 : dest->data.digit &= mask;
1483 : : }
1484 : : }
1485 : 772 : bigint_normalize(dest);
1486 : 772 : return;
1487 : : }
1488 : 276 : dest->digit_count = (bit_count + 63) / 64;
1489 : 276 : assert(dest->digit_count >= op->digit_count);
1490 : 276 : dest->data.digits = allocate_nonzero<uint64_t>(dest->digit_count);
1491 : 276 : size_t i = 0;
1492 [ + + ]: 692 : for (; i < op->digit_count; i += 1) {
1493 : 416 : dest->data.digits[i] = ~op_digits[i];
1494 : : }
1495 [ + + ]: 412 : for (; i < dest->digit_count; i += 1) {
1496 : 136 : dest->data.digits[i] = 0xffffffffffffffffULL;
1497 : : }
1498 : 276 : size_t digit_index = dest->digit_count - 1;
1499 : 276 : size_t digit_bit_index = bit_count % 64;
1500 [ + + ]: 276 : if (digit_bit_index != 0) {
1501 : 128 : uint64_t mask = (1ULL << digit_bit_index) - 1;
1502 : 128 : dest->data.digits[digit_index] &= mask;
1503 : : }
1504 : 276 : bigint_normalize(dest);
1505 : : }
1506 : :
1507 : 6932 : void bigint_truncate(BigInt *dest, const BigInt *op, size_t bit_count, bool is_signed) {
1508 : : BigInt twos_comp;
1509 : 6932 : to_twos_complement(&twos_comp, op, bit_count);
1510 : 6932 : from_twos_complement(dest, &twos_comp, bit_count, is_signed);
1511 : 6932 : }
1512 : :
1513 : 772066 : Cmp bigint_cmp(const BigInt *op1, const BigInt *op2) {
1514 [ + + ][ + + ]: 772066 : if (op1->is_negative && !op2->is_negative) {
1515 : 14050 : return CmpLT;
1516 [ + + ][ + + ]: 758016 : } else if (!op1->is_negative && op2->is_negative) {
1517 : 39260 : return CmpGT;
1518 [ + + ]: 718756 : } else if (op1->digit_count > op2->digit_count) {
1519 [ - + ]: 19115 : return op1->is_negative ? CmpLT : CmpGT;
1520 [ + + ]: 699641 : } else if (op2->digit_count > op1->digit_count) {
1521 [ + + ]: 3282 : return op1->is_negative ? CmpGT : CmpLT;
1522 [ + + ]: 696359 : } else if (op1->digit_count == 0) {
1523 : 37875 : return CmpEQ;
1524 : : }
1525 : 658484 : const uint64_t *op1_digits = bigint_ptr(op1);
1526 : 658484 : const uint64_t *op2_digits = bigint_ptr(op2);
1527 : 658484 : for (size_t i = op1->digit_count - 1; ;) {
1528 : 658972 : uint64_t op1_digit = op1_digits[i];
1529 : 658972 : uint64_t op2_digit = op2_digits[i];
1530 : :
1531 [ + + ]: 658972 : if (op1_digit > op2_digit) {
1532 [ + + ]: 272291 : return op1->is_negative ? CmpLT : CmpGT;
1533 : : }
1534 [ + + ]: 386681 : if (op1_digit < op2_digit) {
1535 [ + + ]: 166996 : return op1->is_negative ? CmpGT : CmpLT;
1536 : : }
1537 : :
1538 [ + + ]: 219685 : if (i == 0) {
1539 : 219197 : return CmpEQ;
1540 : : }
1541 : 488 : i -= 1;
1542 : 488 : }
1543 : : }
1544 : :
1545 : 88 : void bigint_append_buf(Buf *buf, const BigInt *op, uint64_t base) {
1546 [ + + ]: 88 : if (op->digit_count == 0) {
1547 : 24 : buf_append_char(buf, '0');
1548 : 88 : return;
1549 : : }
1550 [ - + ]: 64 : if (op->is_negative) {
1551 : 0 : buf_append_char(buf, '-');
1552 : : }
1553 [ + - ][ + - ]: 64 : if (op->digit_count == 1 && base == 10) {
1554 : 64 : buf_appendf(buf, "%" ZIG_PRI_u64, op->data.digit);
1555 : 64 : return;
1556 : : }
1557 [ # # ][ # # ]: 0 : if (op->digit_count == 1 && base == 16) {
1558 : 0 : buf_appendf(buf, "%" ZIG_PRI_x64, op->data.digit);
1559 : 0 : return;
1560 : : }
1561 : 0 : size_t first_digit_index = buf_len(buf);
1562 : :
1563 : 0 : BigInt digit_bi = {0};
1564 : 0 : BigInt a1 = {0};
1565 : 0 : BigInt a2 = {0};
1566 : :
1567 : 0 : BigInt *a = &a1;
1568 : 0 : BigInt *other_a = &a2;
1569 : 0 : bigint_init_bigint(a, op);
1570 : :
1571 : 0 : BigInt base_bi = {0};
1572 : 0 : bigint_init_unsigned(&base_bi, base);
1573 : :
1574 : : for (;;) {
1575 : 0 : bigint_rem(&digit_bi, a, &base_bi);
1576 : 0 : uint8_t digit = bigint_as_unsigned(&digit_bi);
1577 : 0 : buf_append_char(buf, digit_to_char(digit, false));
1578 : 0 : bigint_div_trunc(other_a, a, &base_bi);
1579 : : {
1580 : 0 : BigInt *tmp = a;
1581 : 0 : a = other_a;
1582 : 0 : other_a = tmp;
1583 : : }
1584 [ # # ]: 0 : if (bigint_cmp_zero(a) == CmpEQ) {
1585 : 0 : break;
1586 : : }
1587 : 0 : }
1588 : :
1589 : : // reverse
1590 [ # # ]: 0 : for (size_t i = first_digit_index; i < buf_len(buf) / 2; i += 1) {
1591 : 0 : size_t other_i = buf_len(buf) + first_digit_index - i - 1;
1592 : 0 : uint8_t tmp = buf_ptr(buf)[i];
1593 : 0 : buf_ptr(buf)[i] = buf_ptr(buf)[other_i];
1594 : 0 : buf_ptr(buf)[other_i] = tmp;
1595 : : }
1596 : : }
1597 : :
1598 : 538 : size_t bigint_popcount_unsigned(const BigInt *bi) {
1599 : 538 : assert(!bi->is_negative);
1600 [ + + ]: 538 : if (bi->digit_count == 0)
1601 : 61 : return 0;
1602 : :
1603 : 477 : size_t count = 0;
1604 : 477 : size_t bit_count = bi->digit_count * 64;
1605 [ + + ]: 31005 : for (size_t i = 0; i < bit_count; i += 1) {
1606 [ + + ]: 30528 : if (bit_at_index(bi, i))
1607 : 1779 : count += 1;
1608 : : }
1609 : 477 : return count;
1610 : : }
1611 : :
1612 : 16 : size_t bigint_popcount_signed(const BigInt *bi, size_t bit_count) {
1613 [ - + ]: 16 : if (bit_count == 0)
1614 : 0 : return 0;
1615 [ - + ]: 16 : if (bi->digit_count == 0)
1616 : 0 : return 0;
1617 : :
1618 : 16 : BigInt twos_comp = {0};
1619 : 16 : to_twos_complement(&twos_comp, bi, bit_count);
1620 : :
1621 : 16 : size_t count = 0;
1622 [ + + ]: 208 : for (size_t i = 0; i < bit_count; i += 1) {
1623 [ + + ]: 192 : if (bit_at_index(&twos_comp, i))
1624 : 144 : count += 1;
1625 : : }
1626 : 16 : return count;
1627 : : }
1628 : :
1629 : 48 : size_t bigint_ctz(const BigInt *bi, size_t bit_count) {
1630 [ - + ]: 48 : if (bit_count == 0)
1631 : 0 : return 0;
1632 [ + + ]: 48 : if (bi->digit_count == 0)
1633 : 16 : return bit_count;
1634 : :
1635 : 32 : BigInt twos_comp = {0};
1636 : 32 : to_twos_complement(&twos_comp, bi, bit_count);
1637 : :
1638 : 32 : size_t count = 0;
1639 [ + - ]: 2104 : for (size_t i = 0; i < bit_count; i += 1) {
1640 [ + + ]: 2104 : if (bit_at_index(&twos_comp, i))
1641 : 32 : return count;
1642 : 2072 : count += 1;
1643 : : }
1644 : 48 : return count;
1645 : : }
1646 : :
1647 : 362763 : size_t bigint_clz(const BigInt *bi, size_t bit_count) {
1648 [ + + ][ + + ]: 362763 : if (bi->is_negative || bit_count == 0)
1649 : 153 : return 0;
1650 [ + + ]: 362610 : if (bi->digit_count == 0)
1651 : 8 : return bit_count;
1652 : :
1653 : 362602 : size_t count = 0;
1654 : 362602 : for (size_t i = bit_count - 1;;) {
1655 [ + + ]: 19321480 : if (bit_at_index(bi, i))
1656 : 362602 : return count;
1657 : 18958878 : count += 1;
1658 : :
1659 [ - + ]: 18958878 : if (i == 0) break;
1660 : 18958878 : i -= 1;
1661 : : }
1662 : 0 : return count;
1663 : : }
1664 : :
1665 : 221326 : static uint64_t bigint_as_unsigned(const BigInt *bigint) {
1666 : 221326 : assert(!bigint->is_negative);
1667 [ + + ]: 221326 : if (bigint->digit_count == 0) {
1668 : 13381 : return 0;
1669 [ + - ]: 207945 : } else if (bigint->digit_count == 1) {
1670 : 207945 : return bigint->data.digit;
1671 : : } else {
1672 : 0 : zig_unreachable();
1673 : : }
1674 : : }
1675 : :
1676 : 101536 : uint64_t bigint_as_u64(const BigInt *bigint)
1677 : : {
1678 : 101536 : return bigint_as_unsigned(bigint);
1679 : : }
1680 : :
1681 : 4637 : uint32_t bigint_as_u32(const BigInt *bigint) {
1682 : 4637 : uint64_t value64 = bigint_as_unsigned(bigint);
1683 : 4637 : uint32_t value32 = (uint32_t)value64;
1684 : 4637 : assert (value64 == value32);
1685 : 4637 : return value32;
1686 : : }
1687 : :
1688 : 3024 : size_t bigint_as_usize(const BigInt *bigint) {
1689 : 3024 : uint64_t value64 = bigint_as_unsigned(bigint);
1690 : 3024 : size_t valueUsize = (size_t)value64;
1691 : 3024 : assert (value64 == valueUsize);
1692 : 3024 : return valueUsize;
1693 : : }
1694 : :
1695 : 3429 : int64_t bigint_as_signed(const BigInt *bigint) {
1696 [ + + ]: 3429 : if (bigint->digit_count == 0) {
1697 : 547 : return 0;
1698 [ + - ]: 2882 : } else if (bigint->digit_count == 1) {
1699 [ + + ]: 2882 : if (bigint->is_negative) {
1700 [ + - ]: 8 : if (bigint->data.digit <= 9223372036854775808ULL) {
1701 : 8 : return (-((int64_t)(bigint->data.digit - 1))) - 1;
1702 : : } else {
1703 : 0 : zig_unreachable();
1704 : : }
1705 : : } else {
1706 : 2874 : return bigint->data.digit;
1707 : : }
1708 : : } else {
1709 : 0 : zig_unreachable();
1710 : : }
1711 : : }
1712 : :
1713 : 342808 : Cmp bigint_cmp_zero(const BigInt *op) {
1714 [ + + ]: 342808 : if (op->digit_count == 0) {
1715 : 18796 : return CmpEQ;
1716 : : }
1717 [ + + ]: 324012 : return op->is_negative ? CmpLT : CmpGT;
1718 : : }
1719 : :
1720 : 406035 : uint32_t bigint_hash(BigInt x) {
1721 [ + + ]: 406035 : if (x.digit_count == 0) {
1722 : 53695 : return 0;
1723 : : } else {
1724 : 352340 : return bigint_ptr(&x)[0];
1725 : : }
1726 : : }
1727 : :
1728 : 90929 : bool bigint_eql(BigInt a, BigInt b) {
1729 : 90929 : return bigint_cmp(&a, &b) == CmpEQ;
1730 : : }
1731 : :
1732 : 113751 : void bigint_incr(BigInt *x) {
1733 [ + + ]: 113751 : if (x->digit_count == 0) {
1734 : 12709 : bigint_init_unsigned(x, 1);
1735 : 113751 : return;
1736 : : }
1737 : :
1738 [ + - ]: 101042 : if (x->digit_count == 1) {
1739 [ + + ][ + - ]: 101042 : if (x->is_negative && x->data.digit != 0) {
1740 : 8 : x->data.digit -= 1;
1741 : 8 : return;
1742 [ + - ][ + - ]: 101034 : } else if (!x->is_negative && x->data.digit != UINT64_MAX) {
1743 : 101034 : x->data.digit += 1;
1744 : 101034 : return;
1745 : : }
1746 : : }
1747 : :
1748 : : BigInt copy;
1749 : 0 : bigint_init_bigint(©, x);
1750 : :
1751 : : BigInt one;
1752 : 0 : bigint_init_unsigned(&one, 1);
1753 : :
1754 : 0 : bigint_add(x, ©, &one);
1755 : : }
1756 : :
|
