1/*
2 Name: gmp_compat.c
3 Purpose: Provide GMP compatiable routines for imath library
4 Author: David Peixotto
5
6 Copyright (c) 2012 Qualcomm Innovation Center, Inc. All rights reserved.
7
8 Permission is hereby granted, free of charge, to any person obtaining a copy
9 of this software and associated documentation files (the "Software"), to deal
10 in the Software without restriction, including without limitation the rights
11 to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
12 copies of the Software, and to permit persons to whom the Software is
13 furnished to do so, subject to the following conditions:
14
15 The above copyright notice and this permission notice shall be included in
16 all copies or substantial portions of the Software.
17
18 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
19 IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
20 FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
21 AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
22 LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
23 OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
24 SOFTWARE.
25 */
26#include "gmp_compat.h"
27#include <assert.h>
28#include <ctype.h>
29#include <stdio.h>
30#include <stdlib.h>
31#include <string.h>
32
33#if defined(_MSC_VER)
34#include <BaseTsd.h>
35typedef SSIZE_T ssize_t;
36#else
37#include <sys/types.h>
38#endif
39
40#ifdef NDEBUG
41#define CHECK(res) (res)
42#else
43#define CHECK(res) assert(((res) == MP_OK) && "expected MP_OK")
44#endif
45
46/* *(signed char *)&endian_test will thus either be:
47 * 0b00000001 = 1 on big-endian
48 * 0b11111111 = -1 on little-endian */
49static const uint16_t endian_test = 0x1FF;
50#define HOST_ENDIAN (*(signed char *)&endian_test)
51
52/*************************************************************************
53 *
54 * Functions with direct translations
55 *
56 *************************************************************************/
57/* gmp: mpq_clear */
58void GMPQAPI(clear)(mp_rat x) { mp_rat_clear(r: x); }
59
60/* gmp: mpq_cmp */
61int GMPQAPI(cmp)(mp_rat op1, mp_rat op2) { return mp_rat_compare(a: op1, b: op2); }
62
63/* gmp: mpq_init */
64void GMPQAPI(init)(mp_rat x) { CHECK(mp_rat_init(x)); }
65
66/* gmp: mpq_mul */
67void GMPQAPI(mul)(mp_rat product, mp_rat multiplier, mp_rat multiplicand) {
68 CHECK(mp_rat_mul(multiplier, multiplicand, product));
69}
70
71/* gmp: mpq_set */
72void GMPQAPI(set)(mp_rat rop, mp_rat op) { CHECK(mp_rat_copy(op, rop)); }
73
74/* gmp: mpz_abs */
75void GMPZAPI(abs)(mp_int rop, mp_int op) { CHECK(mp_int_abs(op, rop)); }
76
77/* gmp: mpz_add */
78void GMPZAPI(add)(mp_int rop, mp_int op1, mp_int op2) {
79 CHECK(mp_int_add(op1, op2, rop));
80}
81
82/* gmp: mpz_clear */
83void GMPZAPI(clear)(mp_int x) { mp_int_clear(z: x); }
84
85/* gmp: mpz_cmp_si */
86int GMPZAPI(cmp_si)(mp_int op1, long op2) {
87 return mp_int_compare_value(z: op1, v: op2);
88}
89
90/* gmp: mpz_cmpabs */
91int GMPZAPI(cmpabs)(mp_int op1, mp_int op2) {
92 return mp_int_compare_unsigned(a: op1, b: op2);
93}
94
95/* gmp: mpz_cmp */
96int GMPZAPI(cmp)(mp_int op1, mp_int op2) { return mp_int_compare(a: op1, b: op2); }
97
98/* gmp: mpz_init */
99void GMPZAPI(init)(mp_int x) { CHECK(mp_int_init(x)); }
100
101/* gmp: mpz_mul */
102void GMPZAPI(mul)(mp_int rop, mp_int op1, mp_int op2) {
103 CHECK(mp_int_mul(op1, op2, rop));
104}
105
106/* gmp: mpz_neg */
107void GMPZAPI(neg)(mp_int rop, mp_int op) { CHECK(mp_int_neg(op, rop)); }
108
109/* gmp: mpz_set_si */
110void GMPZAPI(set_si)(mp_int rop, long op) { CHECK(mp_int_set_value(rop, op)); }
111
112/* gmp: mpz_set */
113void GMPZAPI(set)(mp_int rop, mp_int op) { CHECK(mp_int_copy(op, rop)); }
114
115/* gmp: mpz_sub */
116void GMPZAPI(sub)(mp_int rop, mp_int op1, mp_int op2) {
117 CHECK(mp_int_sub(op1, op2, rop));
118}
119
120/* gmp: mpz_swap */
121void GMPZAPI(swap)(mp_int rop1, mp_int rop2) { mp_int_swap(a: rop1, c: rop2); }
122
123/* gmp: mpq_sgn */
124int GMPQAPI(sgn)(mp_rat op) { return mp_rat_compare_zero(r: op); }
125
126/* gmp: mpz_sgn */
127int GMPZAPI(sgn)(mp_int op) { return mp_int_compare_zero(z: op); }
128
129/* gmp: mpq_set_ui */
130void GMPQAPI(set_ui)(mp_rat rop, unsigned long op1, unsigned long op2) {
131 CHECK(mp_rat_set_uvalue(rop, op1, op2));
132}
133
134/* gmp: mpz_set_ui */
135void GMPZAPI(set_ui)(mp_int rop, unsigned long op) {
136 CHECK(mp_int_set_uvalue(rop, op));
137}
138
139/* gmp: mpq_den_ref */
140mp_int GMPQAPI(denref)(mp_rat op) { return mp_rat_denom_ref(r: op); }
141
142/* gmp: mpq_num_ref */
143mp_int GMPQAPI(numref)(mp_rat op) { return mp_rat_numer_ref(r: op); }
144
145/* gmp: mpq_canonicalize */
146void GMPQAPI(canonicalize)(mp_rat op) { CHECK(mp_rat_reduce(op)); }
147
148/*
149 * Functions that can be implemented as a combination of imath functions
150 */
151
152/* gmp: mpz_addmul */
153/* gmp: rop = rop + (op1 * op2) */
154void GMPZAPI(addmul)(mp_int rop, mp_int op1, mp_int op2) {
155 mpz_t tempz;
156 mp_int temp = &tempz;
157 mp_int_init(z: temp);
158
159 CHECK(mp_int_mul(op1, op2, temp));
160 CHECK(mp_int_add(rop, temp, rop));
161 mp_int_clear(z: temp);
162}
163
164/* gmp: mpz_divexact */
165/* gmp: only produces correct results when d divides n */
166void GMPZAPI(divexact)(mp_int q, mp_int n, mp_int d) {
167 CHECK(mp_int_div(n, d, q, NULL));
168}
169
170/* gmp: mpz_divisible_p */
171/* gmp: return 1 if d divides n, 0 otherwise */
172/* gmp: 0 is considered to divide only 0 */
173int GMPZAPI(divisible_p)(mp_int n, mp_int d) {
174 /* variables to hold remainder */
175 mpz_t rz;
176 mp_int r = &rz;
177 int r_is_zero;
178
179 /* check for d = 0 */
180 int n_is_zero = mp_int_compare_zero(z: n) == 0;
181 int d_is_zero = mp_int_compare_zero(z: d) == 0;
182 if (d_is_zero) return n_is_zero;
183
184 /* return true if remainder is 0 */
185 CHECK(mp_int_init(r));
186 CHECK(mp_int_div(n, d, NULL, r));
187 r_is_zero = mp_int_compare_zero(z: r) == 0;
188 mp_int_clear(z: r);
189
190 return r_is_zero;
191}
192
193/* gmp: mpz_submul */
194/* gmp: rop = rop - (op1 * op2) */
195void GMPZAPI(submul)(mp_int rop, mp_int op1, mp_int op2) {
196 mpz_t tempz;
197 mp_int temp = &tempz;
198 mp_int_init(z: temp);
199
200 CHECK(mp_int_mul(op1, op2, temp));
201 CHECK(mp_int_sub(rop, temp, rop));
202
203 mp_int_clear(z: temp);
204}
205
206/* gmp: mpz_add_ui */
207void GMPZAPI(add_ui)(mp_int rop, mp_int op1, unsigned long op2) {
208 mpz_t tempz;
209 mp_int temp = &tempz;
210 CHECK(mp_int_init_uvalue(temp, op2));
211
212 CHECK(mp_int_add(op1, temp, rop));
213
214 mp_int_clear(z: temp);
215}
216
217/* gmp: mpz_divexact_ui */
218/* gmp: only produces correct results when d divides n */
219void GMPZAPI(divexact_ui)(mp_int q, mp_int n, unsigned long d) {
220 mpz_t tempz;
221 mp_int temp = &tempz;
222 CHECK(mp_int_init_uvalue(temp, d));
223
224 CHECK(mp_int_div(n, temp, q, NULL));
225
226 mp_int_clear(z: temp);
227}
228
229/* gmp: mpz_mul_ui */
230void GMPZAPI(mul_ui)(mp_int rop, mp_int op1, unsigned long op2) {
231 mpz_t tempz;
232 mp_int temp = &tempz;
233 CHECK(mp_int_init_uvalue(temp, op2));
234
235 CHECK(mp_int_mul(op1, temp, rop));
236
237 mp_int_clear(z: temp);
238}
239
240/* gmp: mpz_pow_ui */
241/* gmp: 0^0 = 1 */
242void GMPZAPI(pow_ui)(mp_int rop, mp_int base, unsigned long exp) {
243 mpz_t tempz;
244 mp_int temp = &tempz;
245
246 /* check for 0^0 */
247 if (exp == 0 && mp_int_compare_zero(z: base) == 0) {
248 CHECK(mp_int_set_value(rop, 1));
249 return;
250 }
251
252 /* rop = base^exp */
253 CHECK(mp_int_init_uvalue(temp, exp));
254 CHECK(mp_int_expt_full(base, temp, rop));
255 mp_int_clear(z: temp);
256}
257
258/* gmp: mpz_sub_ui */
259void GMPZAPI(sub_ui)(mp_int rop, mp_int op1, unsigned long op2) {
260 mpz_t tempz;
261 mp_int temp = &tempz;
262 CHECK(mp_int_init_uvalue(temp, op2));
263
264 CHECK(mp_int_sub(op1, temp, rop));
265
266 mp_int_clear(z: temp);
267}
268
269/*************************************************************************
270 *
271 * Functions with different behavior in corner cases
272 *
273 *************************************************************************/
274
275/* gmp: mpz_gcd */
276void GMPZAPI(gcd)(mp_int rop, mp_int op1, mp_int op2) {
277 int op1_is_zero = mp_int_compare_zero(z: op1) == 0;
278 int op2_is_zero = mp_int_compare_zero(z: op2) == 0;
279
280 if (op1_is_zero && op2_is_zero) {
281 mp_int_zero(z: rop);
282 return;
283 }
284
285 CHECK(mp_int_gcd(op1, op2, rop));
286}
287
288/* gmp: mpz_get_str */
289char *GMPZAPI(get_str)(char *str, int radix, mp_int op) {
290 int i, r, len;
291
292 /* Support negative radix like gmp */
293 r = radix;
294 if (r < 0) r = -r;
295
296 /* Compute the length of the string needed to hold the int */
297 len = mp_int_string_len(z: op, radix: r);
298 if (str == NULL) {
299 str = malloc(size: len);
300 }
301
302 /* Convert to string using imath function */
303 CHECK(mp_int_to_string(op, r, str, len));
304
305 /* Change case to match gmp */
306 for (i = 0; i < len - 1; i++) {
307 if (radix < 0) {
308 str[i] = toupper(c: str[i]);
309 } else {
310 str[i] = tolower(c: str[i]);
311 }
312 }
313 return str;
314}
315
316/* gmp: mpq_get_str */
317char *GMPQAPI(get_str)(char *str, int radix, mp_rat op) {
318 int i, r, len;
319
320 /* Only print numerator if it is a whole number */
321 if (mp_int_compare_value(z: mp_rat_denom_ref(r: op), v: 1) == 0)
322 return GMPZAPI(get_str)(str, radix, op: mp_rat_numer_ref(r: op));
323
324 /* Support negative radix like gmp */
325 r = radix;
326 if (r < 0) r = -r;
327
328 /* Compute the length of the string needed to hold the int */
329 len = mp_rat_string_len(r: op, radix: r);
330 if (str == NULL) {
331 str = malloc(size: len);
332 }
333
334 /* Convert to string using imath function */
335 CHECK(mp_rat_to_string(op, r, str, len));
336
337 /* Change case to match gmp */
338 for (i = 0; i < len; i++) {
339 if (radix < 0) {
340 str[i] = toupper(c: str[i]);
341 } else {
342 str[i] = tolower(c: str[i]);
343 }
344 }
345
346 return str;
347}
348
349/* gmp: mpz_set_str */
350int GMPZAPI(set_str)(mp_int rop, char *str, int base) {
351 mp_result res = mp_int_read_string(z: rop, radix: base, str);
352 return ((res == MP_OK) ? 0 : -1);
353}
354
355/* gmp: mpq_set_str */
356int GMPQAPI(set_str)(mp_rat rop, char *s, int base) {
357 char *slash;
358 char *str;
359 mp_result resN;
360 mp_result resD;
361 int res = 0;
362
363 /* Copy string to temporary storage so we can modify it below */
364 str = malloc(size: strlen(s: s) + 1);
365 strcpy(dest: str, src: s);
366
367 /* Properly format the string as an int by terminating at the / */
368 slash = strchr(s: str, c: '/');
369 if (slash) *slash = '\0';
370
371 /* Parse numerator */
372 resN = mp_int_read_string(z: mp_rat_numer_ref(r: rop), radix: base, str);
373
374 /* Parse denominator if given or set to 1 if not */
375 if (slash) {
376 resD = mp_int_read_string(z: mp_rat_denom_ref(r: rop), radix: base, str: slash + 1);
377 } else {
378 resD = mp_int_set_uvalue(z: mp_rat_denom_ref(r: rop), uvalue: 1);
379 }
380
381 /* Return failure if either parse failed */
382 if (resN != MP_OK || resD != MP_OK) {
383 res = -1;
384 }
385
386 free(ptr: str);
387 return res;
388}
389
390static unsigned long get_long_bits(mp_int op) {
391 /* Deal with integer that does not fit into unsigned long. We want to grab
392 * the least significant digits that will fit into the long. Read the digits
393 * into the long starting at the most significant digit that fits into a
394 * long. The long is shifted over by MP_DIGIT_BIT before each digit is added.
395 *
396 * The shift is decomposed into two steps (following the pattern used in the
397 * rest of the imath library) to accommodate architectures that don't deal
398 * well with 32-bit shifts.
399 */
400 mp_size digits_to_copy =
401 (sizeof(unsigned long) + sizeof(mp_digit) - 1) / sizeof(mp_digit);
402 if (digits_to_copy > MP_USED(Z: op)) {
403 digits_to_copy = MP_USED(Z: op);
404 }
405
406 mp_digit *digits = MP_DIGITS(Z: op);
407 unsigned long out = 0;
408
409 for (int i = digits_to_copy - 1; i >= 0; i--) {
410 out <<= (MP_DIGIT_BIT / 2);
411 out <<= (MP_DIGIT_BIT / 2);
412 out |= digits[i];
413 }
414
415 return out;
416}
417
418/* gmp: mpz_get_ui */
419unsigned long GMPZAPI(get_ui)(mp_int op) {
420 unsigned long out;
421
422 /* Try a standard conversion that fits into an unsigned long */
423 mp_result res = mp_int_to_uint(z: op, out: &out);
424 if (res == MP_OK) return out;
425
426 /* Abort the try if we don't have a range error in the conversion.
427 * The range error indicates that the value cannot fit into a long. */
428 CHECK(res == MP_RANGE ? MP_OK : MP_RANGE);
429 if (res != MP_RANGE) return 0;
430
431 return get_long_bits(op);
432}
433
434/* gmp: mpz_get_si */
435long GMPZAPI(get_si)(mp_int op) {
436 long out;
437 unsigned long uout;
438 int long_msb;
439
440 /* Try a standard conversion that fits into a long */
441 mp_result res = mp_int_to_int(z: op, out: &out);
442 if (res == MP_OK) return out;
443
444 /* Abort the try if we don't have a range error in the conversion.
445 * The range error indicates that the value cannot fit into a long. */
446 CHECK(res == MP_RANGE ? MP_OK : MP_RANGE);
447 if (res != MP_RANGE) return 0;
448
449 /* get least significant bits into an unsigned long */
450 uout = get_long_bits(op);
451
452 /* clear the top bit */
453 long_msb = (sizeof(unsigned long) * 8) - 1;
454 uout &= (~(1UL << long_msb));
455
456 /* convert to negative if needed based on sign of op */
457 if (MP_SIGN(Z: op) == MP_NEG) {
458 uout = 0 - uout;
459 }
460
461 out = (long)uout;
462 return out;
463}
464
465/* gmp: mpz_lcm */
466void GMPZAPI(lcm)(mp_int rop, mp_int op1, mp_int op2) {
467 int op1_is_zero = mp_int_compare_zero(z: op1) == 0;
468 int op2_is_zero = mp_int_compare_zero(z: op2) == 0;
469
470 if (op1_is_zero || op2_is_zero) {
471 mp_int_zero(z: rop);
472 return;
473 }
474
475 CHECK(mp_int_lcm(op1, op2, rop));
476 CHECK(mp_int_abs(rop, rop));
477}
478
479/* gmp: mpz_mul_2exp */
480/* gmp: allow big values for op2 when op1 == 0 */
481void GMPZAPI(mul_2exp)(mp_int rop, mp_int op1, unsigned long op2) {
482 if (mp_int_compare_zero(z: op1) == 0)
483 mp_int_zero(z: rop);
484 else
485 CHECK(mp_int_mul_pow2(op1, op2, rop));
486}
487
488/*
489 * Functions needing expanded functionality
490 */
491/* [Note]Overview of division implementation
492
493 All division operations (N / D) compute q and r such that
494
495 N = q * D + r, with 0 <= abs(r) < abs(d)
496
497 The q and r values are not uniquely specified by N and D. To specify which q
498 and r values should be used, GMP implements three different rounding modes
499 for integer division:
500
501 ceiling - round q twords +infinity, r has opposite sign as d
502 floor - round q twords -infinity, r has same sign as d
503 truncate - round q twords zero, r has same sign as n
504
505 The imath library only supports truncate as a rounding mode. We need to
506 implement the other rounding modes in terms of truncating division. We first
507 perform the division in trucate mode and then adjust q accordingly. Once we
508 know q, we can easily compute the correct r according the the formula above
509 by computing:
510
511 r = N - q * D
512
513 The main task is to compute q. We can compute the correct q from a trucated
514 version as follows.
515
516 For ceiling rounding mode, if q is less than 0 then the truncated rounding
517 mode is the same as the ceiling rounding mode. If q is greater than zero
518 then we need to round q up by one because the truncated version was rounded
519 down to zero. If q equals zero then check to see if the result of the
520 divison is positive. A positive result needs to increment q to one.
521
522 For floor rounding mode, if q is greater than 0 then the trucated rounding
523 mode is the same as the floor rounding mode. If q is less than zero then we
524 need to round q down by one because the trucated mode rounded q up by one
525 twords zero. If q is zero then we need to check to see if the result of the
526 division is negative. A negative result needs to decrement q to negative
527 one.
528 */
529
530/* gmp: mpz_cdiv_q */
531void GMPZAPI(cdiv_q)(mp_int q, mp_int n, mp_int d) {
532 mpz_t rz;
533 mp_int r = &rz;
534 int qsign, rsign, nsign, dsign;
535 CHECK(mp_int_init(r));
536
537 /* save signs before division because q can alias with n or d */
538 nsign = mp_int_compare_zero(z: n);
539 dsign = mp_int_compare_zero(z: d);
540
541 /* truncating division */
542 CHECK(mp_int_div(n, d, q, r));
543
544 /* see: [Note]Overview of division implementation */
545 qsign = mp_int_compare_zero(z: q);
546 rsign = mp_int_compare_zero(z: r);
547 if (qsign > 0) { /* q > 0 */
548 if (rsign != 0) { /* r != 0 */
549 CHECK(mp_int_add_value(q, 1, q));
550 }
551 } else if (qsign == 0) { /* q == 0 */
552 if (rsign != 0) { /* r != 0 */
553 if ((nsign > 0 && dsign > 0) || (nsign < 0 && dsign < 0)) {
554 CHECK(mp_int_set_value(q, 1));
555 }
556 }
557 }
558 mp_int_clear(z: r);
559}
560
561/* gmp: mpz_fdiv_q */
562void GMPZAPI(fdiv_q)(mp_int q, mp_int n, mp_int d) {
563 mpz_t rz;
564 mp_int r = &rz;
565 int qsign, rsign, nsign, dsign;
566 CHECK(mp_int_init(r));
567
568 /* save signs before division because q can alias with n or d */
569 nsign = mp_int_compare_zero(z: n);
570 dsign = mp_int_compare_zero(z: d);
571
572 /* truncating division */
573 CHECK(mp_int_div(n, d, q, r));
574
575 /* see: [Note]Overview of division implementation */
576 qsign = mp_int_compare_zero(z: q);
577 rsign = mp_int_compare_zero(z: r);
578 if (qsign < 0) { /* q < 0 */
579 if (rsign != 0) { /* r != 0 */
580 CHECK(mp_int_sub_value(q, 1, q));
581 }
582 } else if (qsign == 0) { /* q == 0 */
583 if (rsign != 0) { /* r != 0 */
584 if ((nsign < 0 && dsign > 0) || (nsign > 0 && dsign < 0)) {
585 CHECK(mp_int_set_value(q, -1));
586 }
587 }
588 }
589 mp_int_clear(z: r);
590}
591
592/* gmp: mpz_fdiv_r */
593void GMPZAPI(fdiv_r)(mp_int r, mp_int n, mp_int d) {
594 mpz_t qz;
595 mpz_t tempz;
596 mpz_t orig_dz;
597 mpz_t orig_nz;
598 mp_int q = &qz;
599 mp_int temp = &tempz;
600 mp_int orig_d = &orig_dz;
601 mp_int orig_n = &orig_nz;
602 CHECK(mp_int_init(q));
603 CHECK(mp_int_init(temp));
604 /* Make a copy of n in case n and d in case they overlap with q */
605 CHECK(mp_int_init_copy(orig_d, d));
606 CHECK(mp_int_init_copy(orig_n, n));
607
608 /* floor division */
609 GMPZAPI(fdiv_q)(q, n, d);
610
611 /* see: [Note]Overview of division implementation */
612 /* n = q * d + r ==> r = n - q * d */
613 mp_int_mul(a: q, b: orig_d, c: temp);
614 mp_int_sub(a: orig_n, b: temp, c: r);
615
616 mp_int_clear(z: q);
617 mp_int_clear(z: temp);
618 mp_int_clear(z: orig_d);
619 mp_int_clear(z: orig_n);
620}
621
622/* gmp: mpz_tdiv_q */
623void GMPZAPI(tdiv_q)(mp_int q, mp_int n, mp_int d) {
624 /* truncating division*/
625 CHECK(mp_int_div(n, d, q, NULL));
626}
627
628/* gmp: mpz_fdiv_q_ui */
629unsigned long GMPZAPI(fdiv_q_ui)(mp_int q, mp_int n, unsigned long d) {
630 mpz_t tempz;
631 mp_int temp = &tempz;
632 mpz_t rz;
633 mp_int r = &rz;
634 mpz_t orig_nz;
635 mp_int orig_n = &orig_nz;
636 unsigned long rl;
637 CHECK(mp_int_init_uvalue(temp, d));
638 CHECK(mp_int_init(r));
639 /* Make a copy of n in case n and q overlap */
640 CHECK(mp_int_init_copy(orig_n, n));
641
642 /* use floor division mode to compute q and r */
643 GMPZAPI(fdiv_q)(q, n, d: temp);
644 GMPZAPI(fdiv_r)(r, n: orig_n, d: temp);
645 CHECK(mp_int_to_uint(r, &rl));
646
647 mp_int_clear(z: temp);
648 mp_int_clear(z: r);
649 mp_int_clear(z: orig_n);
650
651 return rl;
652}
653
654/* gmp: mpz_export */
655void *GMPZAPI(export)(void *rop, size_t *countp, int order, size_t size,
656 int endian, size_t nails, mp_int op) {
657 size_t i, j;
658 size_t num_used_bytes;
659 size_t num_words, num_missing_bytes;
660 ssize_t word_offset;
661 unsigned char *dst;
662 mp_digit *src;
663 int src_bits;
664
665 /* We do not have a complete implementation. Assert to ensure our
666 * restrictions are in place.
667 */
668 assert(nails == 0 && "Do not support non-full words");
669 assert(endian == 1 || endian == 0 || endian == -1);
670 assert(order == 1 || order == -1);
671
672 /* Test for zero */
673 if (mp_int_compare_zero(z: op) == 0) {
674 if (countp) *countp = 0;
675 return rop;
676 }
677
678 /* Calculate how many words we need */
679 num_used_bytes = mp_int_unsigned_len(z: op);
680 num_words = (num_used_bytes + (size - 1)) / size; /* ceil division */
681 assert(num_used_bytes > 0);
682
683 /* Check to see if we will have missing bytes in the last word.
684
685 Missing bytes can only occur when the size of words we output is
686 greater than the size of words used internally by imath. The number of
687 missing bytes is the number of bytes needed to fill out the last word. If
688 this number is greater than the size of a single mp_digit, then we need to
689 pad the word with extra zeros. Otherwise, the missing bytes can be filled
690 directly from the zeros in the last digit in the number.
691 */
692 num_missing_bytes = (size * num_words) - num_used_bytes;
693 assert(num_missing_bytes < size);
694
695 /* Allocate space for the result if needed */
696 if (rop == NULL) {
697 rop = malloc(size: num_words * size);
698 }
699
700 if (endian == 0) {
701 endian = HOST_ENDIAN;
702 }
703
704 /* Initialize dst and src pointers */
705 dst = (unsigned char *)rop + (order >= 0 ? (num_words - 1) * size : 0) +
706 (endian >= 0 ? size - 1 : 0);
707 src = MP_DIGITS(Z: op);
708 src_bits = MP_DIGIT_BIT;
709
710 word_offset = (endian >= 0 ? size : -size) + (order < 0 ? size : -size);
711
712 for (i = 0; i < num_words; i++) {
713 for (j = 0; j < size && i * size + j < num_used_bytes; j++) {
714 if (src_bits == 0) {
715 ++src;
716 src_bits = MP_DIGIT_BIT;
717 }
718 *dst = (*src >> (MP_DIGIT_BIT - src_bits)) & 0xFF;
719 src_bits -= 8;
720 dst -= endian;
721 }
722 for (; j < size; j++) {
723 *dst = 0;
724 dst -= endian;
725 }
726 dst += word_offset;
727 }
728
729 if (countp) *countp = num_words;
730 return rop;
731}
732
733/* gmp: mpz_import */
734void GMPZAPI(import)(mp_int rop, size_t count, int order, size_t size,
735 int endian, size_t nails, const void *op) {
736 mpz_t tmpz;
737 mp_int tmp = &tmpz;
738 size_t total_size;
739 size_t num_digits;
740 ssize_t word_offset;
741 const unsigned char *src;
742 mp_digit *dst;
743 int dst_bits;
744 size_t i, j;
745 if (count == 0 || op == NULL) return;
746
747 /* We do not have a complete implementation. Assert to ensure our
748 * restrictions are in place. */
749 assert(nails == 0 && "Do not support non-full words");
750 assert(endian == 1 || endian == 0 || endian == -1);
751 assert(order == 1 || order == -1);
752
753 if (endian == 0) {
754 endian = HOST_ENDIAN;
755 }
756
757 /* Compute number of needed digits by ceil division */
758 total_size = count * size;
759 num_digits = (total_size + sizeof(mp_digit) - 1) / sizeof(mp_digit);
760
761 /* Init temporary */
762 mp_int_init_size(z: tmp, prec: num_digits);
763 for (i = 0; i < num_digits; i++) tmp->digits[i] = 0;
764
765 /* Copy bytes */
766 src = (const unsigned char *)op + (order >= 0 ? (count - 1) * size : 0) +
767 (endian >= 0 ? size - 1 : 0);
768 dst = MP_DIGITS(Z: tmp);
769 dst_bits = 0;
770
771 word_offset = (endian >= 0 ? size : -size) + (order < 0 ? size : -size);
772
773 for (i = 0; i < count; i++) {
774 for (j = 0; j < size; j++) {
775 if (dst_bits == MP_DIGIT_BIT) {
776 ++dst;
777 dst_bits = 0;
778 }
779 *dst |= ((mp_digit)*src) << dst_bits;
780 dst_bits += 8;
781 src -= endian;
782 }
783 src += word_offset;
784 }
785
786 tmp->used = num_digits;
787
788 /* Remove leading zeros from number */
789 {
790 mp_size uz_ = tmp->used;
791 mp_digit *dz_ = MP_DIGITS(Z: tmp) + uz_ - 1;
792 while (uz_ > 1 && (*dz_-- == 0)) --uz_;
793 tmp->used = uz_;
794 }
795
796 /* Copy to destination */
797 mp_int_copy(a: tmp, c: rop);
798 mp_int_clear(z: tmp);
799}
800
801/* gmp: mpz_sizeinbase */
802size_t GMPZAPI(sizeinbase)(mp_int op, int base) {
803 mp_result res;
804 size_t size;
805
806 /* If op == 0, return 1 */
807 if (mp_int_compare_zero(z: op) == 0) return 1;
808
809 /* Compute string length in base */
810 res = mp_int_string_len(z: op, radix: base);
811 CHECK((res > 0) == MP_OK);
812
813 /* Now adjust the final size by getting rid of string artifacts */
814 size = res;
815
816 /* subtract one for the null terminator */
817 size -= 1;
818
819 /* subtract one for the negative sign */
820 if (mp_int_compare_zero(z: op) < 0) size -= 1;
821
822 return size;
823}
824

source code of polly/lib/External/isl/imath/gmp_compat.c