brintos

brintos / llvm-project-archived public Read only

0
0
Text · 22.1 KiB · 620c4ca Raw
824 lines · c
1/*2  Name:     gmp_compat.c3  Purpose:  Provide GMP compatiable routines for imath library4  Author:   David Peixotto5 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 copy9  of this software and associated documentation files (the "Software"), to deal10  in the Software without restriction, including without limitation the rights11  to use, copy, modify, merge, publish, distribute, sublicense, and/or sell12  copies of the Software, and to permit persons to whom the Software is13  furnished to do so, subject to the following conditions:14 15  The above copyright notice and this permission notice shall be included in16  all copies or substantial portions of the Software.17 18  THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR19  IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,20  FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.  IN NO EVENT SHALL THE21  AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER22  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 THE24  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#else37#include <sys/types.h>38#endif39 40#ifdef NDEBUG41#define CHECK(res) (res)42#else43#define CHECK(res) assert(((res) == MP_OK) && "expected MP_OK")44#endif45 46/* *(signed char *)&endian_test will thus either be:47 *     0b00000001 =  1 on big-endian48 *     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 translations55 *56 *************************************************************************/57/* gmp: mpq_clear */58void GMPQAPI(clear)(mp_rat x) { mp_rat_clear(x); }59 60/* gmp: mpq_cmp */61int GMPQAPI(cmp)(mp_rat op1, mp_rat op2) { return mp_rat_compare(op1, 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(x); }84 85/* gmp: mpz_cmp_si */86int GMPZAPI(cmp_si)(mp_int op1, long op2) {87  return mp_int_compare_value(op1, op2);88}89 90/* gmp: mpz_cmpabs */91int GMPZAPI(cmpabs)(mp_int op1, mp_int op2) {92  return mp_int_compare_unsigned(op1, op2);93}94 95/* gmp: mpz_cmp */96int GMPZAPI(cmp)(mp_int op1, mp_int op2) { return mp_int_compare(op1, 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(rop1, rop2); }122 123/* gmp: mpq_sgn */124int GMPQAPI(sgn)(mp_rat op) { return mp_rat_compare_zero(op); }125 126/* gmp: mpz_sgn */127int GMPZAPI(sgn)(mp_int op) { return mp_int_compare_zero(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(op); }141 142/* gmp: mpq_num_ref */143mp_int GMPQAPI(numref)(mp_rat op) { return mp_rat_numer_ref(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 functions150 */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(temp);158 159  CHECK(mp_int_mul(op1, op2, temp));160  CHECK(mp_int_add(rop, temp, rop));161  mp_int_clear(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(n) == 0;181  int d_is_zero = mp_int_compare_zero(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(r) == 0;188  mp_int_clear(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(temp);199 200  CHECK(mp_int_mul(op1, op2, temp));201  CHECK(mp_int_sub(rop, temp, rop));202 203  mp_int_clear(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(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(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(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(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(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(temp);267}268 269/*************************************************************************270 *271 * Functions with different behavior in corner cases272 *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(op1) == 0;278  int op2_is_zero = mp_int_compare_zero(op2) == 0;279 280  if (op1_is_zero && op2_is_zero) {281    mp_int_zero(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(op, r);298  if (str == NULL) {299    str = malloc(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(str[i]);309    } else {310      str[i] = tolower(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(mp_rat_denom_ref(op), 1) == 0)322    return GMPZAPI(get_str)(str, radix, mp_rat_numer_ref(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(op, r);330  if (str == NULL) {331    str = malloc(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(str[i]);341    } else {342      str[i] = tolower(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(rop, 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(strlen(s) + 1);365  strcpy(str, s);366 367  /* Properly format the string as an int by terminating at the / */368  slash = strchr(str, '/');369  if (slash) *slash = '\0';370 371  /* Parse numerator */372  resN = mp_int_read_string(mp_rat_numer_ref(rop), base, str);373 374  /* Parse denominator if given or set to 1 if not */375  if (slash) {376    resD = mp_int_read_string(mp_rat_denom_ref(rop), base, slash + 1);377  } else {378    resD = mp_int_set_uvalue(mp_rat_denom_ref(rop), 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(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 grab392   * the least significant digits that will fit into the long.  Read the digits393   * into the long starting at the most significant digit that fits into a394   * 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 the397   * rest of the imath library) to accommodate architectures that don't deal398   * 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(op)) {403    digits_to_copy = MP_USED(op);404  }405 406  mp_digit *digits = MP_DIGITS(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(op, &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(op, &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(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(op1) == 0;468  int op2_is_zero = mp_int_compare_zero(op2) == 0;469 470  if (op1_is_zero || op2_is_zero) {471    mp_int_zero(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(op1) == 0)483    mp_int_zero(rop);484  else485    CHECK(mp_int_mul_pow2(op1, op2, rop));486}487 488/*489 * Functions needing expanded functionality490 */491/* [Note]Overview of division implementation492 493    All division operations (N / D) compute q and r such that494 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 q498    and r values should be used, GMP implements three different rounding modes499    for integer division:500 501      ceiling  - round q twords +infinity, r has opposite sign as d502      floor    - round q twords -infinity, r has same sign as d503      truncate - round q twords zero,      r has same sign as n504 505    The imath library only supports truncate as a rounding mode. We need to506    implement the other rounding modes in terms of truncating division. We first507    perform the division in trucate mode and then adjust q accordingly. Once we508    know q, we can easily compute the correct r according the the formula above509    by computing:510 511      r = N - q * D512 513    The main task is to compute q. We can compute the correct q from a trucated514    version as follows.515 516    For ceiling rounding mode, if q is less than 0 then the truncated rounding517    mode is the same as the ceiling rounding mode.  If q is greater than zero518    then we need to round q up by one because the truncated version was rounded519    down to zero. If q equals zero then check to see if the result of the520    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 rounding523    mode is the same as the floor rounding mode. If q is less than zero then we524    need to round q down by one because the trucated mode rounded q up by one525    twords zero. If q is zero then we need to check to see if the result of the526    division is negative. A negative result needs to decrement q to negative527    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(n);539  dsign = mp_int_compare_zero(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(q);546  rsign = mp_int_compare_zero(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(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(n);570  dsign = mp_int_compare_zero(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(q);577  rsign = mp_int_compare_zero(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(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(q, orig_d, temp);614  mp_int_sub(orig_n, temp, r);615 616  mp_int_clear(q);617  mp_int_clear(temp);618  mp_int_clear(orig_d);619  mp_int_clear(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, temp);644  GMPZAPI(fdiv_r)(r, orig_n, temp);645  CHECK(mp_int_to_uint(r, &rl));646 647  mp_int_clear(temp);648  mp_int_clear(r);649  mp_int_clear(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 our666   * 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(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(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 is686     greater than the size of words used internally by imath. The number of687     missing bytes is the number of bytes needed to fill out the last word. If688     this number is greater than the size of a single mp_digit, then we need to689     pad the word with extra zeros. Otherwise, the missing bytes can be filled690     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(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(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 our748   * 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(tmp, 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(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(tmp) + uz_ - 1;792    while (uz_ > 1 && (*dz_-- == 0)) --uz_;793    tmp->used = uz_;794  }795 796  /* Copy to destination */797  mp_int_copy(tmp, rop);798  mp_int_clear(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(op) == 0) return 1;808 809  /* Compute string length in base */810  res = mp_int_string_len(op, 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(op) < 0) size -= 1;821 822  return size;823}824