src/bignum.c

   1 /* Big numbers for Emacs.
   2
   3 Copyright 2018-2019 Free Software Foundation, Inc.
   4
   5 This file is part of GNU Emacs.
   6
   7 GNU Emacs is free software: you can redistribute it and/or modify
   8 it under the terms of the GNU General Public License as published by
   9 the Free Software Foundation, either version 3 of the License, or (at
  10 your option) any later version.
  11
  12 GNU Emacs is distributed in the hope that it will be useful,
  13 but WITHOUT ANY WARRANTY; without even the implied warranty of
  14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15 GNU General Public License for more details.
  16
  17 You should have received a copy of the GNU General Public License
  18 along with GNU Emacs.  If not, see <https://www.gnu.org/licenses/>.  */
  19
  20 #include <config.h>
  21
  22 #include "bignum.h"
  23
  24 #include "lisp.h"
  25
  26 #include <math.h>
  27 #include <stdlib.h>
  28
  29 /* mpz global temporaries.  Making them global saves the trouble of
  30    properly using mpz_init and mpz_clear on temporaries even when
  31    storage is exhausted.  Admittedly this is not ideal.  An mpz value
  32    in a temporary is made permanent by mpz_swapping it with a bignum's
  33    value.  Although typically at most two temporaries are needed,
  34    time_arith, rounddiv_q and rounding_driver each need four.  */
  35
  36 mpz_t mpz[4];
  37
  38 static void *
  39 xrealloc_for_gmp (void *ptr, size_t ignore, size_t size)
  40 {
  41   return xrealloc (ptr, size);
  42 }
  43
  44 static void
  45 xfree_for_gmp (void *ptr, size_t ignore)
  46 {
  47   xfree (ptr);
  48 }
  49
  50 void
  51 init_bignum (void)
  52 {
  53   eassert (mp_bits_per_limb == GMP_NUMB_BITS);
  54   integer_width = 1 << 16;
  55   mp_set_memory_functions (xmalloc, xrealloc_for_gmp, xfree_for_gmp);
  56
  57   for (int i = 0; i < ARRAYELTS (mpz); i++)
  58     mpz_init (mpz[i]);
  59 }
  60
  61 /* Return the value of the Lisp bignum N, as a double.  */
  62 double
  63 bignum_to_double (Lisp_Object n)
  64 {
  65   return mpz_get_d_rounded (XBIGNUM (n)->value);
  66 }
  67
  68 /* Return D, converted to a Lisp integer.  Discard any fraction.
  69    Signal an error if D cannot be converted.  */
  70 Lisp_Object
  71 double_to_integer (double d)
  72 {
  73   if (!isfinite (d))
  74     overflow_error ();
  75   mpz_set_d (mpz[0], d);
  76   return make_integer_mpz ();
  77 }
  78
  79 /* Return a Lisp integer equal to mpz[0], which has BITS bits and which
  80    must not be in fixnum range.  Set mpz[0] to a junk value.  */
  81 static Lisp_Object
  82 make_bignum_bits (size_t bits)
  83 {
  84   /* The documentation says integer-width should be nonnegative, so
  85      comparing it to BITS works even though BITS is unsigned.  Treat
  86      integer-width as if it were at least twice the machine integer width,
  87      so that timefns.c can safely use bignums for double-precision
  88      timestamps.  */
  89   if (integer_width < bits && 2 * max (INTMAX_WIDTH, UINTMAX_WIDTH) < bits)
  90     overflow_error ();
  91
  92   struct Lisp_Bignum *b = ALLOCATE_PLAIN_PSEUDOVECTOR (struct Lisp_Bignum,
  93                                                        PVEC_BIGNUM);
  94   mpz_init (b->value);
  95   mpz_swap (b->value, mpz[0]);
  96   return make_lisp_ptr (b, Lisp_Vectorlike);
  97 }
  98
  99 /* Return a Lisp integer equal to mpz[0], which must not be in fixnum range.
 100    Set mpz[0] to a junk value.  */
 101 static Lisp_Object
 102 make_bignum (void)
 103 {
 104   return make_bignum_bits (mpz_sizeinbase (mpz[0], 2));
 105 }
 106
 107 /* Return a Lisp integer equal to N, which must not be in fixnum range.  */
 108 Lisp_Object
 109 make_bigint (intmax_t n)
 110 {
 111   eassert (FIXNUM_OVERFLOW_P (n));
 112   mpz_set_intmax (mpz[0], n);
 113   return make_bignum ();
 114 }
 115 Lisp_Object
 116 make_biguint (uintmax_t n)
 117 {
 118   eassert (FIXNUM_OVERFLOW_P (n));
 119   mpz_set_uintmax (mpz[0], n);
 120   return make_bignum ();
 121 }
 122
 123 /* Return a Lisp integer equal to -N, which must not be in fixnum range.  */
 124 Lisp_Object
 125 make_neg_biguint (uintmax_t n)
 126 {
 127   eassert (-MOST_NEGATIVE_FIXNUM < n);
 128   mpz_set_uintmax (mpz[0], n);
 129   mpz_neg (mpz[0], mpz[0]);
 130   return make_bignum ();
 131 }
 132
 133 /* Return a Lisp integer with value taken from mpz[0].
 134    Set mpz[0] to a junk value.  */
 135 Lisp_Object
 136 make_integer_mpz (void)
 137 {
 138   size_t bits = mpz_sizeinbase (mpz[0], 2);
 139
 140   if (bits <= FIXNUM_BITS)
 141     {
 142       EMACS_INT v = 0;
 143       int i = 0, shift = 0;
 144
 145       do
 146         {
 147           EMACS_INT limb = mpz_getlimbn (mpz[0], i++);
 148           v += limb << shift;
 149           shift += GMP_NUMB_BITS;
 150         }
 151       while (shift < bits);
 152
 153       if (mpz_sgn (mpz[0]) < 0)
 154         v = -v;
 155
 156       if (!FIXNUM_OVERFLOW_P (v))
 157         return make_fixnum (v);
 158     }
 159
 160   return make_bignum_bits (bits);
 161 }
 162
 163 /* Set RESULT to V.  This code is for when intmax_t is wider than long.  */
 164 void
 165 mpz_set_intmax_slow (mpz_t result, intmax_t v)
 166 {
 167   int maxlimbs = (INTMAX_WIDTH + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS;
 168   mp_limb_t *limb = mpz_limbs_write (result, maxlimbs);
 169   int n = 0;
 170   uintmax_t u = v;
 171   bool negative = v < 0;
 172   if (negative)
 173     {
 174       uintmax_t two = 2;
 175       u = -u & ((two << (UINTMAX_WIDTH - 1)) - 1);
 176     }
 177
 178   do
 179     {
 180       limb[n++] = u;
 181       u = GMP_NUMB_BITS < UINTMAX_WIDTH ? u >> GMP_NUMB_BITS : 0;
 182     }
 183   while (u != 0);
 184
 185   mpz_limbs_finish (result, negative ? -n : n);
 186 }
 187 void
 188 mpz_set_uintmax_slow (mpz_t result, uintmax_t v)
 189 {
 190   int maxlimbs = (UINTMAX_WIDTH + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS;
 191   mp_limb_t *limb = mpz_limbs_write (result, maxlimbs);
 192   int n = 0;
 193
 194   do
 195     {
 196       limb[n++] = v;
 197       v = GMP_NUMB_BITS < INTMAX_WIDTH ? v >> GMP_NUMB_BITS : 0;
 198     }
 199   while (v != 0);
 200
 201   mpz_limbs_finish (result, n);
 202 }
 203
 204 /* If Z fits into *PI, store its value there and return true.
 205    Return false otherwise.  */
 206 bool
 207 mpz_to_intmax (mpz_t const z, intmax_t *pi)
 208 {
 209   ptrdiff_t bits = mpz_sizeinbase (z, 2);
 210   bool negative = mpz_sgn (z) < 0;
 211
 212   if (bits < INTMAX_WIDTH)
 213     {
 214       intmax_t v = 0;
 215       int i = 0, shift = 0;
 216
 217       do
 218         {
 219           intmax_t limb = mpz_getlimbn (z, i++);
 220           v += limb << shift;
 221           shift += GMP_NUMB_BITS;
 222         }
 223       while (shift < bits);
 224
 225       *pi = negative ? -v : v;
 226       return true;
 227     }
 228   if (bits == INTMAX_WIDTH && INTMAX_MIN < -INTMAX_MAX && negative
 229       && mpz_scan1 (z, 0) == INTMAX_WIDTH - 1)
 230     {
 231       *pi = INTMAX_MIN;
 232       return true;
 233     }
 234   return false;
 235 }
 236 bool
 237 mpz_to_uintmax (mpz_t const z, uintmax_t *pi)
 238 {
 239   if (mpz_sgn (z) < 0)
 240     return false;
 241   ptrdiff_t bits = mpz_sizeinbase (z, 2);
 242   if (UINTMAX_WIDTH < bits)
 243     return false;
 244
 245   uintmax_t v = 0;
 246   int i = 0, shift = 0;
 247
 248   do
 249     {
 250       uintmax_t limb = mpz_getlimbn (z, i++);
 251       v += limb << shift;
 252       shift += GMP_NUMB_BITS;
 253     }
 254   while (shift < bits);
 255
 256   *pi = v;
 257   return true;
 258 }
 259
 260 /* Return the value of the bignum X if it fits, 0 otherwise.
 261    A bignum cannot be zero, so 0 indicates failure reliably.  */
 262 intmax_t
 263 bignum_to_intmax (Lisp_Object x)
 264 {
 265   intmax_t i;
 266   return mpz_to_intmax (XBIGNUM (x)->value, &i) ? i : 0;
 267 }
 268 uintmax_t
 269 bignum_to_uintmax (Lisp_Object x)
 270 {
 271   uintmax_t i;
 272   return mpz_to_uintmax (XBIGNUM (x)->value, &i) ? i : 0;
 273 }
 274
 275 /* Yield an upper bound on the buffer size needed to contain a C
 276    string representing the NUM in base BASE.  This includes any
 277    preceding '-' and the terminating NUL.  */
 278 static ptrdiff_t
 279 mpz_bufsize (mpz_t const num, int base)
 280 {
 281   return mpz_sizeinbase (num, base) + 2;
 282 }
 283 ptrdiff_t
 284 bignum_bufsize (Lisp_Object num, int base)
 285 {
 286   return mpz_bufsize (XBIGNUM (num)->value, base);
 287 }
 288
 289 /* Convert NUM to a nearest double, as opposed to mpz_get_d which
 290    truncates toward zero.  */
 291 double
 292 mpz_get_d_rounded (mpz_t const num)
 293 {
 294   ptrdiff_t size = mpz_bufsize (num, 10);
 295
 296   /* Use mpz_get_d as a shortcut for a bignum so small that rounding
 297      errors cannot occur, which is possible if EMACS_INT (not counting
 298      sign) has fewer bits than a double significand.  */
 299   if (! ((FLT_RADIX == 2 && DBL_MANT_DIG <= FIXNUM_BITS - 1)
 300          || (FLT_RADIX == 16 && DBL_MANT_DIG * 4 <= FIXNUM_BITS - 1))
 301       && size <= DBL_DIG + 2)
 302     return mpz_get_d (num);
 303
 304   USE_SAFE_ALLOCA;
 305   char *buf = SAFE_ALLOCA (size);
 306   mpz_get_str (buf, 10, num);
 307   double result = strtod (buf, NULL);
 308   SAFE_FREE ();
 309   return result;
 310 }
 311
 312 /* Store into BUF (of size SIZE) the value of NUM as a base-BASE string.
 313    If BASE is negative, use upper-case digits in base -BASE.
 314    Return the string's length.
 315    SIZE must equal bignum_bufsize (NUM, abs (BASE)).  */
 316 ptrdiff_t
 317 bignum_to_c_string (char *buf, ptrdiff_t size, Lisp_Object num, int base)
 318 {
 319   eassert (bignum_bufsize (num, abs (base)) == size);
 320   mpz_get_str (buf, base, XBIGNUM (num)->value);
 321   ptrdiff_t n = size - 2;
 322   return !buf[n - 1] ? n - 1 : n + !!buf[n];
 323 }
 324
 325 /* Convert NUM to a base-BASE Lisp string.
 326    If BASE is negative, use upper-case digits in base -BASE.  */
 327
 328 Lisp_Object
 329 bignum_to_string (Lisp_Object num, int base)
 330 {
 331   ptrdiff_t size = bignum_bufsize (num, abs (base));
 332   USE_SAFE_ALLOCA;
 333   char *str = SAFE_ALLOCA (size);
 334   ptrdiff_t len = bignum_to_c_string (str, size, num, base);
 335   Lisp_Object result = make_unibyte_string (str, len);
 336   SAFE_FREE ();
 337   return result;
 338 }
 339
 340 /* Create a bignum by scanning NUM, with digits in BASE.
 341    NUM must consist of an optional '-', a nonempty sequence
 342    of base-BASE digits, and a terminating NUL byte, and
 343    the represented number must not be in fixnum range.  */
 344
 345 Lisp_Object
 346 make_bignum_str (char const *num, int base)
 347 {
 348   struct Lisp_Bignum *b = ALLOCATE_PLAIN_PSEUDOVECTOR (struct Lisp_Bignum,
 349                                                        PVEC_BIGNUM);
 350   mpz_init (b->value);
 351   int check = mpz_set_str (b->value, num, base);
 352   eassert (check == 0);
 353   return make_lisp_ptr (b, Lisp_Vectorlike);
 354 }