stdlib/gmp-impl.h

   1 /* Include file for internal GNU MP types and definitions.
   2
   3 Copyright (C) 1991, 1993, 1994, 1995, 1996 Free Software Foundation, Inc.
   4
   5 This file is part of the GNU MP Library.
   6
   7 The GNU MP Library is free software; you can redistribute it and/or modify
   8 it under the terms of the GNU Library General Public License as published by
   9 the Free Software Foundation; either version 2 of the License, or (at your
  10 option) any later version.
  11
  12 The GNU MP Library is distributed in the hope that it will be useful, but
  13 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
  14 or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Library General Public
  15 License for more details.
  16
  17 You should have received a copy of the GNU Library General Public License
  18 along with the GNU MP Library; see the file COPYING.LIB.  If not, write to
  19 the Free Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. */
  20
  21 /* When using gcc, make sure to use its builtin alloca.  */
  22 #if ! defined (alloca) && defined (__GNUC__)
  23 #define alloca __builtin_alloca
  24 #define HAVE_ALLOCA
  25 #endif
  26
  27 /* When using cc, do whatever necessary to allow use of alloca.  For many
  28    machines, this means including alloca.h.  IBM's compilers need a #pragma
  29    in "each module that needs to use alloca".  */
  30 #if ! defined (alloca)
  31 /* We need lots of variants for MIPS, to cover all versions and perversions
  32    of OSes for MIPS.  */
  33 #if defined (__mips) || defined (MIPSEL) || defined (MIPSEB) \
  34  || defined (_MIPSEL) || defined (_MIPSEB) || defined (__sgi) \
  35  || defined (__alpha) || defined (__sparc) || defined (sparc) \
  36  || defined (__ksr__)
  37 #include <alloca.h>
  38 #define HAVE_ALLOCA
  39 #endif
  40 #if defined (_IBMR2)
  41 #pragma alloca
  42 #define HAVE_ALLOCA
  43 #endif
  44 #if defined (__DECC)
  45 #define alloca(x) __ALLOCA(x)
  46 #define HAVE_ALLOCA
  47 #endif
  48 #endif
  49
  50 #if ! defined (HAVE_ALLOCA) || USE_STACK_ALLOC
  51 #include "stack-alloc.h"
  52 #else
  53 #define TMP_DECL(m)
  54 #define TMP_ALLOC(x) alloca(x)
  55 #define TMP_MARK(m)
  56 #define TMP_FREE(m)
  57 #endif
  58
  59 #ifndef NULL
  60 #define NULL ((void *) 0)
  61 #endif
  62
  63 #if ! defined (__GNUC__)
  64 #define inline                  /* Empty */
  65 #endif
  66
  67 #define ABS(x) (x >= 0 ? x : -x)
  68 #define MIN(l,o) ((l) < (o) ? (l) : (o))
  69 #define MAX(h,i) ((h) > (i) ? (h) : (i))
  70
  71 #include "gmp-mparam.h"
  72 /* #include "longlong.h" */
  73
  74 #if defined (__STDC__)  || defined (__cplusplus)
  75 void *malloc (size_t);
  76 void *realloc (void *, size_t);
  77 void free (void *);
  78
  79 extern void *   (*_mp_allocate_func) (size_t);
  80 extern void *   (*_mp_reallocate_func) (void *, size_t, size_t);
  81 extern void     (*_mp_free_func) (void *, size_t);
  82
  83 void *_mp_default_allocate (size_t);
  84 void *_mp_default_reallocate (void *, size_t, size_t);
  85 void _mp_default_free (void *, size_t);
  86
  87 #else
  88
  89 #define const                   /* Empty */
  90 #define signed                  /* Empty */
  91
  92 void *malloc ();
  93 void *realloc ();
  94 void free ();
  95
  96 extern void *   (*_mp_allocate_func) ();
  97 extern void *   (*_mp_reallocate_func) ();
  98 extern void     (*_mp_free_func) ();
  99
 100 void *_mp_default_allocate ();
 101 void *_mp_default_reallocate ();
 102 void _mp_default_free ();
 103 #endif
 104
 105 /* Copy NLIMBS *limbs* from SRC to DST.  */
 106 #define MPN_COPY_INCR(DST, SRC, NLIMBS) \
 107   do {                                                                  \
 108     mp_size_t __i;                                                      \
 109     for (__i = 0; __i < (NLIMBS); __i++)                                \
 110       (DST)[__i] = (SRC)[__i];                                          \
 111   } while (0)
 112 #define MPN_COPY_DECR(DST, SRC, NLIMBS) \
 113   do {                                                                  \
 114     mp_size_t __i;                                                      \
 115     for (__i = (NLIMBS) - 1; __i >= 0; __i--)                           \
 116       (DST)[__i] = (SRC)[__i];                                          \
 117   } while (0)
 118 #define MPN_COPY MPN_COPY_INCR
 119
 120 /* Zero NLIMBS *limbs* AT DST.  */
 121 #define MPN_ZERO(DST, NLIMBS) \
 122   do {                                                                  \
 123     mp_size_t __i;                                                      \
 124     for (__i = 0; __i < (NLIMBS); __i++)                                \
 125       (DST)[__i] = 0;                                                   \
 126   } while (0)
 127
 128 #define MPN_NORMALIZE(DST, NLIMBS) \
 129   do {                                                                  \
 130     while (NLIMBS > 0)                                                  \
 131       {                                                                 \
 132         if ((DST)[(NLIMBS) - 1] != 0)                                   \
 133           break;                                                        \
 134         NLIMBS--;                                                       \
 135       }                                                                 \
 136   } while (0)
 137 #define MPN_NORMALIZE_NOT_ZERO(DST, NLIMBS) \
 138   do {                                                                  \
 139     while (1)                                                           \
 140       {                                                                 \
 141         if ((DST)[(NLIMBS) - 1] != 0)                                   \
 142           break;                                                        \
 143         NLIMBS--;                                                       \
 144       }                                                                 \
 145   } while (0)
 146
 147 /* Initialize the MP_INT X with space for NLIMBS limbs.
 148    X should be a temporary variable, and it will be automatically
 149    cleared out when the running function returns.
 150    We use __x here to make it possible to accept both mpz_ptr and mpz_t
 151    arguments.  */
 152 #define MPZ_TMP_INIT(X, NLIMBS) \
 153   do {                                                                  \
 154     mpz_ptr __x = (X);                                                  \
 155     __x->_mp_alloc = (NLIMBS);                                          \
 156     __x->_mp_d = (mp_ptr) TMP_ALLOC ((NLIMBS) * BYTES_PER_MP_LIMB);     \
 157   } while (0)
 158
 159 #define MPN_MUL_N_RECURSE(prodp, up, vp, size, tspace) \
 160   do {                                                                  \
 161     if ((size) < KARATSUBA_THRESHOLD)                                   \
 162       impn_mul_n_basecase (prodp, up, vp, size);                        \
 163     else                                                                \
 164       impn_mul_n (prodp, up, vp, size, tspace);                 \
 165   } while (0);
 166 #define MPN_SQR_N_RECURSE(prodp, up, size, tspace) \
 167   do {                                                                  \
 168     if ((size) < KARATSUBA_THRESHOLD)                                   \
 169       impn_sqr_n_basecase (prodp, up, size);                            \
 170     else                                                                \
 171       impn_sqr_n (prodp, up, size, tspace);                             \
 172   } while (0);
 173
 174 /* Structure for conversion between internal binary format and
 175    strings in base 2..36.  */
 176 struct bases
 177 {
 178   /* Number of digits in the conversion base that always fits in an mp_limb.
 179      For example, for base 10 on a machine where a mp_limb has 32 bits this
 180      is 9, since 10**9 is the largest number that fits into a mp_limb.  */
 181   int chars_per_limb;
 182
 183   /* log(2)/log(conversion_base) */
 184   float chars_per_bit_exactly;
 185
 186   /* base**chars_per_limb, i.e. the biggest number that fits a word, built by
 187      factors of base.  Exception: For 2, 4, 8, etc, big_base is log2(base),
 188      i.e. the number of bits used to represent each digit in the base.  */
 189   mp_limb big_base;
 190
 191   /* A BITS_PER_MP_LIMB bit approximation to 1/big_base, represented as a
 192      fixed-point number.  Instead of dividing by big_base an application can
 193      choose to multiply by big_base_inverted.  */
 194   mp_limb big_base_inverted;
 195 };
 196
 197 /* Access macros for structure fields for user-visible structures with
 198    hidden fields.  */
 199 #define size(X) (X)._mp_size
 200 #define alloc(X) (X)._mp_alloc
 201 #define prec(X) (X)._mp_prec
 202 #define limbs(X) (X)._mp_d
 203
 204 extern const struct bases __mp_bases[];
 205 extern mp_size_t __gmp_default_fp_limb_precision;
 206
 207 /* Divide the two-limb number in (NH,,NL) by D, with DI being the largest
 208    limb not larger than (2**(2*BITS_PER_MP_LIMB))/D - (2**BITS_PER_MP_LIMB).
 209    If this would yield overflow, DI should be the largest possible number
 210    (i.e., only ones).  For correct operation, the most significant bit of D
 211    has to be set.  Put the quotient in Q and the remainder in R.  */
 212 #define udiv_qrnnd_preinv(q, r, nh, nl, d, di) \
 213   do {                                                                  \
 214     mp_limb _q, _ql, _r;                                                \
 215     mp_limb _xh, _xl;                                                   \
 216     umul_ppmm (_q, _ql, (nh), (di));                                    \
 217     _q += (nh);                 /* DI is 2**BITS_PER_MP_LIMB too small */\
 218     umul_ppmm (_xh, _xl, _q, (d));                                      \
 219     sub_ddmmss (_xh, _r, (nh), (nl), _xh, _xl);                         \
 220     if (_xh != 0)                                                       \
 221       {                                                                 \
 222         sub_ddmmss (_xh, _r, _xh, _r, 0, (d));                          \
 223         _q += 1;                                                        \
 224         if (_xh != 0)                                                   \
 225           {                                                             \
 226             sub_ddmmss (_xh, _r, _xh, _r, 0, (d));                      \
 227             _q += 1;                                                    \
 228           }                                                             \
 229       }                                                                 \
 230     if (_r >= (d))                                                      \
 231       {                                                                 \
 232         _r -= (d);                                                      \
 233         _q += 1;                                                        \
 234       }                                                                 \
 235     (r) = _r;                                                           \
 236     (q) = _q;                                                           \
 237   } while (0)
 238 /* Like udiv_qrnnd_preinv, but for for any value D.  DNORM is D shifted left
 239    so that its most significant bit is set.  LGUP is ceil(log2(D)).  */
 240 #define udiv_qrnnd_preinv2gen(q, r, nh, nl, d, di, dnorm, lgup) \
 241   do {                                                                  \
 242     mp_limb n2, n10, n1, nadj, q1;                                      \
 243     mp_limb _xh, _xl;                                                   \
 244     n2 = ((nh) << (BITS_PER_MP_LIMB - (lgup))) + ((nl) >> 1 >> (l - 1));\
 245     n10 = (nl) << (BITS_PER_MP_LIMB - (lgup));                          \
 246     n1 = ((mp_limb_signed) n10 >> (BITS_PER_MP_LIMB - 1));              \
 247     nadj = n10 + (n1 & (dnorm));                                        \
 248     umul_ppmm (_xh, _xl, di, n2 - n1);                                  \
 249     add_ssaaaa (_xh, _xl, _xh, _xl, 0, nadj);                           \
 250     q1 = ~(n2 + _xh);                                                   \
 251     umul_ppmm (_xh, _xl, q1, d);                                        \
 252     add_ssaaaa (_xh, _xl, _xh, _xl, nh, nl);                            \
 253     _xh -= (d);                                                         \
 254     (r) = _xl + ((d) & _xh);                                            \
 255     (q) = _xh - q1;                                                     \
 256   } while (0)
 257 /* Exactly like udiv_qrnnd_preinv, but branch-free.  It is not clear which
 258    version to use.  */
 259 #define udiv_qrnnd_preinv2norm(q, r, nh, nl, d, di) \
 260   do {                                                                  \
 261     mp_limb n2, n10, n1, nadj, q1;                                      \
 262     mp_limb _xh, _xl;                                                   \
 263     n2 = (nh);                                                          \
 264     n10 = (nl);                                                         \
 265     n1 = ((mp_limb_signed) n10 >> (BITS_PER_MP_LIMB - 1));              \
 266     nadj = n10 + (n1 & (d));                                            \
 267     umul_ppmm (_xh, _xl, di, n2 - n1);                                  \
 268     add_ssaaaa (_xh, _xl, _xh, _xl, 0, nadj);                           \
 269     q1 = ~(n2 + _xh);                                                   \
 270     umul_ppmm (_xh, _xl, q1, d);                                        \
 271     add_ssaaaa (_xh, _xl, _xh, _xl, nh, nl);                            \
 272     _xh -= (d);                                                         \
 273     (r) = _xl + ((d) & _xh);                                            \
 274     (q) = _xh - q1;                                                     \
 275   } while (0)
 276
 277 #if defined (__GNUC__)
 278 /* Define stuff for longlong.h.  */
 279 typedef unsigned int UQItype    __attribute__ ((mode (QI)));
 280 typedef          int SItype     __attribute__ ((mode (SI)));
 281 typedef unsigned int USItype    __attribute__ ((mode (SI)));
 282 typedef          int DItype     __attribute__ ((mode (DI)));
 283 typedef unsigned int UDItype    __attribute__ ((mode (DI)));
 284 #else
 285 typedef unsigned char UQItype;
 286 typedef          long SItype;
 287 typedef unsigned long USItype;
 288 #endif
 289
 290 typedef mp_limb UWtype;
 291 typedef unsigned int UHWtype;
 292 #define W_TYPE_SIZE BITS_PER_MP_LIMB
 293
 294 /* Internal mpn calls */
 295 #define impn_mul_n_basecase     __MPN(impn_mul_n_basecase)
 296 #define impn_mul_n              __MPN(impn_mul_n)
 297 #define impn_sqr_n_basecase     __MPN(impn_sqr_n_basecase)
 298 #define impn_sqr_n              __MPN(impn_sqr_n)
 299
 300 #ifndef IEEE_DOUBLE_BIG_ENDIAN
 301 #define IEEE_DOUBLE_BIG_ENDIAN 1
 302 #endif
 303
 304 #if IEEE_DOUBLE_BIG_ENDIAN
 305 union ieee_double_extract
 306 {
 307   struct
 308     {
 309       unsigned int sig:1;
 310       unsigned int exp:11;
 311       unsigned int manh:20;
 312       unsigned int manl:32;
 313     } s;
 314   double d;
 315 };
 316 #else
 317 union ieee_double_extract
 318 {
 319   struct
 320     {
 321       unsigned int manl:32;
 322       unsigned int manh:20;
 323       unsigned int exp:11;
 324       unsigned int sig:1;
 325     } s;
 326   double d;
 327 };
 328 #endif