contrib/gmp/gmp-impl.h

   1 /* Include file for internal GNU MP types and definitions.
   2
   3    THE CONTENTS OF THIS FILE ARE FOR INTERNAL USE AND ARE ALMOST CERTAIN TO
   4    BE SUBJECT TO INCOMPATIBLE CHANGES IN FUTURE GNU MP RELEASES.
   5
   6 Copyright 1991, 1993, 1994, 1995, 1996, 1997, 1999, 2000, 2001, 2002, 2003,
   7 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012 Free Software Foundation,
   8 Inc.
   9
  10 This file is part of the GNU MP Library.
  11
  12 The GNU MP Library is free software; you can redistribute it and/or modify
  13 it under the terms of the GNU Lesser General Public License as published by
  14 the Free Software Foundation; either version 3 of the License, or (at your
  15 option) any later version.
  16
  17 The GNU MP Library is distributed in the hope that it will be useful, but
  18 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
  19 or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
  20 License for more details.
  21
  22 You should have received a copy of the GNU Lesser General Public License
  23 along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.  */
  24
  25
  26 /* __GMP_DECLSPEC must be given on any global data that will be accessed
  27    from outside libgmp, meaning from the test or development programs, or
  28    from libgmpxx.  Failing to do this will result in an incorrect address
  29    being used for the accesses.  On functions __GMP_DECLSPEC makes calls
  30    from outside libgmp more efficient, but they'll still work fine without
  31    it.  */
  32
  33
  34 #ifndef __GMP_IMPL_H__
  35 #define __GMP_IMPL_H__
  36
  37 #if defined _CRAY
  38 #include <intrinsics.h>  /* for _popcnt */
  39 #endif
  40
  41 /* limits.h is not used in general, since it's an ANSI-ism, and since on
  42    solaris gcc 2.95 under -mcpu=ultrasparc in ABI=32 ends up getting wrong
  43    values (the ABI=64 values).
  44
  45    On Cray vector systems, however, we need the system limits.h since sizes
  46    of signed and unsigned types can differ there, depending on compiler
  47    options (eg. -hnofastmd), making our SHRT_MAX etc expressions fail.  For
  48    reference, int can be 46 or 64 bits, whereas uint is always 64 bits; and
  49    short can be 24, 32, 46 or 64 bits, and different for ushort.  */
  50
  51 #if defined _CRAY
  52 #include <limits.h>
  53 #endif
  54
  55 /* For fat.h and other fat binary stuff.
  56    No need for __GMP_ATTRIBUTE_PURE or __GMP_NOTHROW, since functions
  57    declared this way are only used to set function pointers in __gmp_cpuvec,
  58    they're not called directly.  */
  59 #define DECL_add_n(name) \
  60   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t))
  61 #define DECL_addmul_1(name) \
  62   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t))
  63 #define DECL_copyd(name) \
  64   __GMP_DECLSPEC void name __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t))
  65 #define DECL_copyi(name) \
  66   DECL_copyd (name)
  67 #define DECL_divexact_1(name) \
  68   __GMP_DECLSPEC void name __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t))
  69 #define DECL_divexact_by3c(name) \
  70   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t))
  71 #define DECL_divrem_1(name) \
  72   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t))
  73 #define DECL_gcd_1(name) \
  74   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t))
  75 #define DECL_lshift(name) \
  76   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, unsigned))
  77 #define DECL_mod_1(name) \
  78   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t))
  79 #define DECL_mod_34lsub1(name) \
  80   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_srcptr, mp_size_t))
  81 #define DECL_modexact_1c_odd(name) \
  82   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t))
  83 #define DECL_mul_1(name) \
  84   DECL_addmul_1 (name)
  85 #define DECL_mul_basecase(name) \
  86   __GMP_DECLSPEC void name __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t))
  87 #define DECL_preinv_divrem_1(name) \
  88   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t, int))
  89 #define DECL_preinv_mod_1(name) \
  90   __GMP_DECLSPEC mp_limb_t name __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t))
  91 #define DECL_rshift(name) \
  92   DECL_lshift (name)
  93 #define DECL_sqr_basecase(name) \
  94   __GMP_DECLSPEC void name __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t))
  95 #define DECL_sub_n(name) \
  96   DECL_add_n (name)
  97 #define DECL_submul_1(name) \
  98   DECL_addmul_1 (name)
  99
 100 #if ! __GMP_WITHIN_CONFIGURE
 101 #include "config.h"
 102 #include "gmp-mparam.h"
 103 #include "fib_table.h"
 104 #include "mp_bases.h"
 105 #if WANT_FAT_BINARY
 106 #include "fat.h"
 107 #endif
 108 #endif
 109
 110 #if HAVE_INTTYPES_H      /* for uint_least32_t */
 111 # include <inttypes.h>
 112 #else
 113 # if HAVE_STDINT_H
 114 #  include <stdint.h>
 115 # endif
 116 #endif
 117
 118 #ifdef __cplusplus
 119 #include <cstring>  /* for strlen */
 120 #include <string>   /* for std::string */
 121 #endif
 122
 123
 124 #ifndef WANT_TMP_DEBUG  /* for TMP_ALLOC_LIMBS_2 and others */
 125 #define WANT_TMP_DEBUG 0
 126 #endif
 127
 128 /* The following tries to get a good version of alloca.  The tests are
 129    adapted from autoconf AC_FUNC_ALLOCA, with a couple of additions.
 130    Whether this succeeds is tested by GMP_FUNC_ALLOCA and HAVE_ALLOCA will
 131    be setup appropriately.
 132
 133    ifndef alloca - a cpp define might already exist.
 134        glibc <stdlib.h> includes <alloca.h> which uses GCC __builtin_alloca.
 135        HP cc +Olibcalls adds a #define of alloca to __builtin_alloca.
 136
 137    GCC __builtin_alloca - preferred whenever available.
 138
 139    _AIX pragma - IBM compilers need a #pragma in "each module that needs to
 140        use alloca".  Pragma indented to protect pre-ANSI cpp's.  _IBMR2 was
 141        used in past versions of GMP, retained still in case it matters.
 142
 143        The autoconf manual says this pragma needs to be at the start of a C
 144        file, apart from comments and preprocessor directives.  Is that true?
 145        xlc on aix 4.xxx doesn't seem to mind it being after prototypes etc
 146        from gmp.h.
 147 */
 148
 149 #ifndef alloca
 150 # ifdef __GNUC__
 151 #  define alloca __builtin_alloca
 152 # else
 153 #  ifdef __DECC
 154 #   define alloca(x) __ALLOCA(x)
 155 #  else
 156 #   ifdef _MSC_VER
 157 #    include <malloc.h>
 158 #    define alloca _alloca
 159 #   else
 160 #    if HAVE_ALLOCA_H
 161 #     include <alloca.h>
 162 #    else
 163 #     if defined (_AIX) || defined (_IBMR2)
 164  #pragma alloca
 165 #     else
 166        char *alloca ();
 167 #     endif
 168 #    endif
 169 #   endif
 170 #  endif
 171 # endif
 172 #endif
 173
 174
 175 /* if not provided by gmp-mparam.h */
 176 #ifndef BYTES_PER_MP_LIMB
 177 #define BYTES_PER_MP_LIMB  SIZEOF_MP_LIMB_T
 178 #endif
 179 #define GMP_LIMB_BYTES  BYTES_PER_MP_LIMB
 180 #ifndef GMP_LIMB_BITS
 181 #define GMP_LIMB_BITS  (8 * SIZEOF_MP_LIMB_T)
 182 #endif
 183
 184 #define BITS_PER_ULONG  (8 * SIZEOF_UNSIGNED_LONG)
 185
 186
 187 /* gmp_uint_least32_t is an unsigned integer type with at least 32 bits. */
 188 #if HAVE_UINT_LEAST32_T
 189 typedef uint_least32_t      gmp_uint_least32_t;
 190 #else
 191 #if SIZEOF_UNSIGNED_SHORT >= 4
 192 typedef unsigned short      gmp_uint_least32_t;
 193 #else
 194 #if SIZEOF_UNSIGNED >= 4
 195 typedef unsigned            gmp_uint_least32_t;
 196 #else
 197 typedef unsigned long       gmp_uint_least32_t;
 198 #endif
 199 #endif
 200 #endif
 201
 202
 203 /* gmp_intptr_t, for pointer to integer casts */
 204 #if HAVE_INTPTR_T
 205 typedef intptr_t            gmp_intptr_t;
 206 #else /* fallback */
 207 typedef size_t              gmp_intptr_t;
 208 #endif
 209
 210
 211 /* pre-inverse types for truncating division and modulo */
 212 typedef struct {mp_limb_t inv32;} gmp_pi1_t;
 213 typedef struct {mp_limb_t inv21, inv32, inv53;} gmp_pi2_t;
 214
 215
 216 /* const and signed must match __gmp_const and __gmp_signed, so follow the
 217    decision made for those in gmp.h.    */
 218 #if ! __GMP_HAVE_CONST
 219 #define const   /* empty */
 220 #define signed  /* empty */
 221 #endif
 222
 223 /* "const" basically means a function does nothing but examine its arguments
 224    and give a return value, it doesn't read or write any memory (neither
 225    global nor pointed to by arguments), and has no other side-effects.  This
 226    is more restrictive than "pure".  See info node "(gcc)Function
 227    Attributes".  __GMP_NO_ATTRIBUTE_CONST_PURE lets tune/common.c etc turn
 228    this off when trying to write timing loops.  */
 229 #if HAVE_ATTRIBUTE_CONST && ! defined (__GMP_NO_ATTRIBUTE_CONST_PURE)
 230 #define ATTRIBUTE_CONST  __attribute__ ((const))
 231 #else
 232 #define ATTRIBUTE_CONST
 233 #endif
 234
 235 #if HAVE_ATTRIBUTE_NORETURN
 236 #define ATTRIBUTE_NORETURN  __attribute__ ((noreturn))
 237 #else
 238 #define ATTRIBUTE_NORETURN
 239 #endif
 240
 241 /* "malloc" means a function behaves like malloc in that the pointer it
 242    returns doesn't alias anything.  */
 243 #if HAVE_ATTRIBUTE_MALLOC
 244 #define ATTRIBUTE_MALLOC  __attribute__ ((malloc))
 245 #else
 246 #define ATTRIBUTE_MALLOC
 247 #endif
 248
 249
 250 #if ! HAVE_STRCHR
 251 #define strchr(s,c)  index(s,c)
 252 #endif
 253
 254 #if ! HAVE_MEMSET
 255 #define memset(p, c, n)                 \
 256   do {                                  \
 257     ASSERT ((n) >= 0);                  \
 258     char *__memset__p = (p);            \
 259     int  __i;                           \
 260     for (__i = 0; __i < (n); __i++)     \
 261       __memset__p[__i] = (c);           \
 262   } while (0)
 263 #endif
 264
 265 /* va_copy is standard in C99, and gcc provides __va_copy when in strict C89
 266    mode.  Falling back to a memcpy will give maximum portability, since it
 267    works no matter whether va_list is a pointer, struct or array.  */
 268 #if ! defined (va_copy) && defined (__va_copy)
 269 #define va_copy(dst,src)  __va_copy(dst,src)
 270 #endif
 271 #if ! defined (va_copy)
 272 #define va_copy(dst,src) \
 273   do { memcpy (&(dst), &(src), sizeof (va_list)); } while (0)
 274 #endif
 275
 276
 277 /* HAVE_HOST_CPU_alpha_CIX is 1 on an alpha with the CIX instructions
 278    (ie. ctlz, ctpop, cttz).  */
 279 #if HAVE_HOST_CPU_alphaev67 || HAVE_HOST_CPU_alphaev68  \
 280   || HAVE_HOST_CPU_alphaev7
 281 #define HAVE_HOST_CPU_alpha_CIX 1
 282 #endif
 283
 284
 285 #if defined (__cplusplus)
 286 extern "C" {
 287 #endif
 288
 289
 290 /* Usage: TMP_DECL;
 291           TMP_MARK;
 292           ptr = TMP_ALLOC (bytes);
 293           TMP_FREE;
 294
 295    Small allocations should use TMP_SALLOC, big allocations should use
 296    TMP_BALLOC.  Allocations that might be small or big should use TMP_ALLOC.
 297
 298    Functions that use just TMP_SALLOC should use TMP_SDECL, TMP_SMARK, and
 299    TMP_SFREE.
 300
 301    TMP_DECL just declares a variable, but might be empty and so must be last
 302    in a list of variables.  TMP_MARK must be done before any TMP_ALLOC.
 303    TMP_ALLOC(0) is not allowed.  TMP_FREE doesn't need to be done if a
 304    TMP_MARK was made, but then no TMP_ALLOCs.  */
 305
 306 /* The alignment in bytes, used for TMP_ALLOCed blocks, when alloca or
 307    __gmp_allocate_func doesn't already determine it.  Currently TMP_ALLOC
 308    isn't used for "double"s, so that's not in the union.  */
 309 union tmp_align_t {
 310   mp_limb_t  l;
 311   char       *p;
 312 };
 313 #define __TMP_ALIGN  sizeof (union tmp_align_t)
 314
 315 /* Return "a" rounded upwards to a multiple of "m", if it isn't already.
 316    "a" must be an unsigned type.
 317    This is designed for use with a compile-time constant "m".
 318    The POW2 case is expected to be usual, and gcc 3.0 and up recognises
 319    "(-(8*n))%8" or the like is always zero, which means the rounding up in
 320    the WANT_TMP_NOTREENTRANT version of TMP_ALLOC below will be a noop.  */
 321 #define ROUND_UP_MULTIPLE(a,m)          \
 322   (POW2_P(m) ? (a) + (-(a))%(m)         \
 323    : (a)+(m)-1 - (((a)+(m)-1) % (m)))
 324
 325 #if defined (WANT_TMP_ALLOCA) || defined (WANT_TMP_REENTRANT)
 326 struct tmp_reentrant_t {
 327   struct tmp_reentrant_t  *next;
 328   size_t                  size;   /* bytes, including header */
 329 };
 330 __GMP_DECLSPEC void *__gmp_tmp_reentrant_alloc __GMP_PROTO ((struct tmp_reentrant_t **, size_t)) ATTRIBUTE_MALLOC;
 331 __GMP_DECLSPEC void  __gmp_tmp_reentrant_free __GMP_PROTO ((struct tmp_reentrant_t *));
 332 #endif
 333
 334 #if WANT_TMP_ALLOCA
 335 #define TMP_SDECL
 336 #define TMP_DECL                struct tmp_reentrant_t *__tmp_marker
 337 #define TMP_SMARK
 338 #define TMP_MARK                __tmp_marker = 0
 339 #define TMP_SALLOC(n)           alloca(n)
 340 #define TMP_BALLOC(n)           __gmp_tmp_reentrant_alloc (&__tmp_marker, n)
 341 #define TMP_ALLOC(n)                                                    \
 342   (LIKELY ((n) < 65536) ? TMP_SALLOC(n) : TMP_BALLOC(n))
 343 #define TMP_SFREE
 344 #define TMP_FREE                                                           \
 345   do {                                                                     \
 346     if (UNLIKELY (__tmp_marker != 0)) __gmp_tmp_reentrant_free (__tmp_marker); \
 347   } while (0)
 348 #endif
 349
 350 #if WANT_TMP_REENTRANT
 351 #define TMP_SDECL               TMP_DECL
 352 #define TMP_DECL                struct tmp_reentrant_t *__tmp_marker
 353 #define TMP_SMARK               TMP_MARK
 354 #define TMP_MARK                __tmp_marker = 0
 355 #define TMP_SALLOC(n)           TMP_ALLOC(n)
 356 #define TMP_BALLOC(n)           TMP_ALLOC(n)
 357 #define TMP_ALLOC(n)            __gmp_tmp_reentrant_alloc (&__tmp_marker, n)
 358 #define TMP_SFREE               TMP_FREE
 359 #define TMP_FREE                __gmp_tmp_reentrant_free (__tmp_marker)
 360 #endif
 361
 362 #if WANT_TMP_NOTREENTRANT
 363 struct tmp_marker
 364 {
 365   struct tmp_stack *which_chunk;
 366   void *alloc_point;
 367 };
 368 __GMP_DECLSPEC void *__gmp_tmp_alloc __GMP_PROTO ((unsigned long)) ATTRIBUTE_MALLOC;
 369 __GMP_DECLSPEC void __gmp_tmp_mark __GMP_PROTO ((struct tmp_marker *));
 370 __GMP_DECLSPEC void __gmp_tmp_free __GMP_PROTO ((struct tmp_marker *));
 371 #define TMP_SDECL               TMP_DECL
 372 #define TMP_DECL                struct tmp_marker __tmp_marker
 373 #define TMP_SMARK               TMP_MARK
 374 #define TMP_MARK                __gmp_tmp_mark (&__tmp_marker)
 375 #define TMP_SALLOC(n)           TMP_ALLOC(n)
 376 #define TMP_BALLOC(n)           TMP_ALLOC(n)
 377 #define TMP_ALLOC(n)                                                    \
 378   __gmp_tmp_alloc (ROUND_UP_MULTIPLE ((unsigned long) (n), __TMP_ALIGN))
 379 #define TMP_SFREE               TMP_FREE
 380 #define TMP_FREE                __gmp_tmp_free (&__tmp_marker)
 381 #endif
 382
 383 #if WANT_TMP_DEBUG
 384 /* See tal-debug.c for some comments. */
 385 struct tmp_debug_t {
 386   struct tmp_debug_entry_t  *list;
 387   const char                *file;
 388   int                       line;
 389 };
 390 struct tmp_debug_entry_t {
 391   struct tmp_debug_entry_t  *next;
 392   char                      *block;
 393   size_t                    size;
 394 };
 395 __GMP_DECLSPEC void  __gmp_tmp_debug_mark  __GMP_PROTO ((const char *, int, struct tmp_debug_t **,
 396                                                          struct tmp_debug_t *,
 397                                                          const char *, const char *));
 398 __GMP_DECLSPEC void *__gmp_tmp_debug_alloc __GMP_PROTO ((const char *, int, int,
 399                                                          struct tmp_debug_t **, const char *,
 400                                                          size_t)) ATTRIBUTE_MALLOC;
 401 __GMP_DECLSPEC void  __gmp_tmp_debug_free  __GMP_PROTO ((const char *, int, int,
 402                                                          struct tmp_debug_t **,
 403                                                          const char *, const char *));
 404 #define TMP_SDECL TMP_DECL_NAME(__tmp_xmarker, "__tmp_marker")
 405 #define TMP_DECL TMP_DECL_NAME(__tmp_xmarker, "__tmp_marker")
 406 #define TMP_SMARK TMP_MARK_NAME(__tmp_xmarker, "__tmp_marker")
 407 #define TMP_MARK TMP_MARK_NAME(__tmp_xmarker, "__tmp_marker")
 408 #define TMP_SFREE TMP_FREE_NAME(__tmp_xmarker, "__tmp_marker")
 409 #define TMP_FREE TMP_FREE_NAME(__tmp_xmarker, "__tmp_marker")
 410 /* The marker variable is designed to provoke an uninitialized variable
 411    warning from the compiler if TMP_FREE is used without a TMP_MARK.
 412    __tmp_marker_inscope does the same for TMP_ALLOC.  Runtime tests pick
 413    these things up too.  */
 414 #define TMP_DECL_NAME(marker, marker_name)                      \
 415   int marker;                                                   \
 416   int __tmp_marker_inscope;                                     \
 417   const char *__tmp_marker_name = marker_name;                  \
 418   struct tmp_debug_t  __tmp_marker_struct;                      \
 419   /* don't demand NULL, just cast a zero */                     \
 420   struct tmp_debug_t  *__tmp_marker = (struct tmp_debug_t *) 0
 421 #define TMP_MARK_NAME(marker, marker_name)                      \
 422   do {                                                          \
 423     marker = 1;                                                 \
 424     __tmp_marker_inscope = 1;                                   \
 425     __gmp_tmp_debug_mark  (ASSERT_FILE, ASSERT_LINE,            \
 426                            &__tmp_marker, &__tmp_marker_struct, \
 427                            __tmp_marker_name, marker_name);     \
 428   } while (0)
 429 #define TMP_SALLOC(n)           TMP_ALLOC(n)
 430 #define TMP_BALLOC(n)           TMP_ALLOC(n)
 431 #define TMP_ALLOC(size)                                                 \
 432   __gmp_tmp_debug_alloc (ASSERT_FILE, ASSERT_LINE,                      \
 433                          __tmp_marker_inscope,                          \
 434                          &__tmp_marker, __tmp_marker_name, size)
 435 #define TMP_FREE_NAME(marker, marker_name)                      \
 436   do {                                                          \
 437     __gmp_tmp_debug_free  (ASSERT_FILE, ASSERT_LINE,            \
 438                            marker, &__tmp_marker,               \
 439                            __tmp_marker_name, marker_name);     \
 440   } while (0)
 441 #endif /* WANT_TMP_DEBUG */
 442
 443
 444 /* Allocating various types. */
 445 #define TMP_ALLOC_TYPE(n,type)  ((type *) TMP_ALLOC ((n) * sizeof (type)))
 446 #define TMP_SALLOC_TYPE(n,type) ((type *) TMP_SALLOC ((n) * sizeof (type)))
 447 #define TMP_BALLOC_TYPE(n,type) ((type *) TMP_BALLOC ((n) * sizeof (type)))
 448 #define TMP_ALLOC_LIMBS(n)      TMP_ALLOC_TYPE(n,mp_limb_t)
 449 #define TMP_SALLOC_LIMBS(n)     TMP_SALLOC_TYPE(n,mp_limb_t)
 450 #define TMP_BALLOC_LIMBS(n)     TMP_BALLOC_TYPE(n,mp_limb_t)
 451 #define TMP_ALLOC_MP_PTRS(n)    TMP_ALLOC_TYPE(n,mp_ptr)
 452 #define TMP_SALLOC_MP_PTRS(n)   TMP_SALLOC_TYPE(n,mp_ptr)
 453 #define TMP_BALLOC_MP_PTRS(n)   TMP_BALLOC_TYPE(n,mp_ptr)
 454
 455 /* It's more efficient to allocate one block than two.  This is certainly
 456    true of the malloc methods, but it can even be true of alloca if that
 457    involves copying a chunk of stack (various RISCs), or a call to a stack
 458    bounds check (mingw).  In any case, when debugging keep separate blocks
 459    so a redzoning malloc debugger can protect each individually.  */
 460 #define TMP_ALLOC_LIMBS_2(xp,xsize, yp,ysize)           \
 461   do {                                                  \
 462     if (WANT_TMP_DEBUG)                                 \
 463       {                                                 \
 464         (xp) = TMP_ALLOC_LIMBS (xsize);                 \
 465         (yp) = TMP_ALLOC_LIMBS (ysize);                 \
 466       }                                                 \
 467     else                                                \
 468       {                                                 \
 469         (xp) = TMP_ALLOC_LIMBS ((xsize) + (ysize));     \
 470         (yp) = (xp) + (xsize);                          \
 471       }                                                 \
 472   } while (0)
 473
 474
 475 /* From gmp.h, nicer names for internal use. */
 476 #define CRAY_Pragma(str)               __GMP_CRAY_Pragma(str)
 477 #define MPN_CMP(result, xp, yp, size)  __GMPN_CMP(result, xp, yp, size)
 478 #define LIKELY(cond)                   __GMP_LIKELY(cond)
 479 #define UNLIKELY(cond)                 __GMP_UNLIKELY(cond)
 480
 481 #define ABS(x) ((x) >= 0 ? (x) : -(x))
 482 #define ABS_CAST(T,x) ((x) >= 0 ? (T)(x) : -((T)((x) + 1) - 1))
 483 #undef MIN
 484 #define MIN(l,o) ((l) < (o) ? (l) : (o))
 485 #undef MAX
 486 #define MAX(h,i) ((h) > (i) ? (h) : (i))
 487 #define numberof(x)  (sizeof (x) / sizeof ((x)[0]))
 488
 489 /* Field access macros.  */
 490 #define SIZ(x) ((x)->_mp_size)
 491 #define ABSIZ(x) ABS (SIZ (x))
 492 #define PTR(x) ((x)->_mp_d)
 493 #define LIMBS(x) ((x)->_mp_d)
 494 #define EXP(x) ((x)->_mp_exp)
 495 #define PREC(x) ((x)->_mp_prec)
 496 #define ALLOC(x) ((x)->_mp_alloc)
 497
 498 /* n-1 inverts any low zeros and the lowest one bit.  If n&(n-1) leaves zero
 499    then that lowest one bit must have been the only bit set.  n==0 will
 500    return true though, so avoid that.  */
 501 #define POW2_P(n)  (((n) & ((n) - 1)) == 0)
 502
 503
 504 /* The "short" defines are a bit different because shorts are promoted to
 505    ints by ~ or >> etc.
 506
 507    #ifndef's are used since on some systems (HP?) header files other than
 508    limits.h setup these defines.  We could forcibly #undef in that case, but
 509    there seems no need to worry about that.  */
 510
 511 #ifndef ULONG_MAX
 512 #define ULONG_MAX   __GMP_ULONG_MAX
 513 #endif
 514 #ifndef UINT_MAX
 515 #define UINT_MAX    __GMP_UINT_MAX
 516 #endif
 517 #ifndef USHRT_MAX
 518 #define USHRT_MAX   __GMP_USHRT_MAX
 519 #endif
 520 #define MP_LIMB_T_MAX      (~ (mp_limb_t) 0)
 521
 522 /* Must cast ULONG_MAX etc to unsigned long etc, since they might not be
 523    unsigned on a K&R compiler.  In particular the HP-UX 10 bundled K&R cc
 524    treats the plain decimal values in <limits.h> as signed.  */
 525 #define ULONG_HIGHBIT      (ULONG_MAX ^ ((unsigned long) ULONG_MAX >> 1))
 526 #define UINT_HIGHBIT       (UINT_MAX ^ ((unsigned) UINT_MAX >> 1))
 527 #define USHRT_HIGHBIT      ((unsigned short) (USHRT_MAX ^ ((unsigned short) USHRT_MAX >> 1)))
 528 #define GMP_LIMB_HIGHBIT  (MP_LIMB_T_MAX ^ (MP_LIMB_T_MAX >> 1))
 529
 530 #ifndef LONG_MIN
 531 #define LONG_MIN           ((long) ULONG_HIGHBIT)
 532 #endif
 533 #ifndef LONG_MAX
 534 #define LONG_MAX           (-(LONG_MIN+1))
 535 #endif
 536
 537 #ifndef INT_MIN
 538 #define INT_MIN            ((int) UINT_HIGHBIT)
 539 #endif
 540 #ifndef INT_MAX
 541 #define INT_MAX            (-(INT_MIN+1))
 542 #endif
 543
 544 #ifndef SHRT_MIN
 545 #define SHRT_MIN           ((short) USHRT_HIGHBIT)
 546 #endif
 547 #ifndef SHRT_MAX
 548 #define SHRT_MAX           ((short) (-(SHRT_MIN+1)))
 549 #endif
 550
 551 #if __GMP_MP_SIZE_T_INT
 552 #define MP_SIZE_T_MAX      INT_MAX
 553 #define MP_SIZE_T_MIN      INT_MIN
 554 #else
 555 #define MP_SIZE_T_MAX      LONG_MAX
 556 #define MP_SIZE_T_MIN      LONG_MIN
 557 #endif
 558
 559 /* mp_exp_t is the same as mp_size_t */
 560 #define MP_EXP_T_MAX   MP_SIZE_T_MAX
 561 #define MP_EXP_T_MIN   MP_SIZE_T_MIN
 562
 563 #define LONG_HIGHBIT       LONG_MIN
 564 #define INT_HIGHBIT        INT_MIN
 565 #define SHRT_HIGHBIT       SHRT_MIN
 566
 567
 568 #define GMP_NUMB_HIGHBIT  (CNST_LIMB(1) << (GMP_NUMB_BITS-1))
 569
 570 #if GMP_NAIL_BITS == 0
 571 #define GMP_NAIL_LOWBIT   CNST_LIMB(0)
 572 #else
 573 #define GMP_NAIL_LOWBIT   (CNST_LIMB(1) << GMP_NUMB_BITS)
 574 #endif
 575
 576 #if GMP_NAIL_BITS != 0
 577 /* Set various *_THRESHOLD values to be used for nails.  Thus we avoid using
 578    code that has not yet been qualified.  */
 579
 580 #undef  DC_DIV_QR_THRESHOLD
 581 #define DC_DIV_QR_THRESHOLD              50
 582
 583 #undef DIVREM_1_NORM_THRESHOLD
 584 #undef DIVREM_1_UNNORM_THRESHOLD
 585 #undef MOD_1_NORM_THRESHOLD
 586 #undef MOD_1_UNNORM_THRESHOLD
 587 #undef USE_PREINV_DIVREM_1
 588 #undef DIVREM_2_THRESHOLD
 589 #undef DIVEXACT_1_THRESHOLD
 590 #define DIVREM_1_NORM_THRESHOLD           MP_SIZE_T_MAX  /* no preinv */
 591 #define DIVREM_1_UNNORM_THRESHOLD         MP_SIZE_T_MAX  /* no preinv */
 592 #define MOD_1_NORM_THRESHOLD              MP_SIZE_T_MAX  /* no preinv */
 593 #define MOD_1_UNNORM_THRESHOLD            MP_SIZE_T_MAX  /* no preinv */
 594 #define USE_PREINV_DIVREM_1               0  /* no preinv */
 595 #define DIVREM_2_THRESHOLD                MP_SIZE_T_MAX  /* no preinv */
 596
 597 /* mpn/generic/mul_fft.c is not nails-capable. */
 598 #undef  MUL_FFT_THRESHOLD
 599 #undef  SQR_FFT_THRESHOLD
 600 #define MUL_FFT_THRESHOLD                MP_SIZE_T_MAX
 601 #define SQR_FFT_THRESHOLD                MP_SIZE_T_MAX
 602 #endif
 603
 604 /* Swap macros. */
 605
 606 #define MP_LIMB_T_SWAP(x, y)                    \
 607   do {                                          \
 608     mp_limb_t __mp_limb_t_swap__tmp = (x);      \
 609     (x) = (y);                                  \
 610     (y) = __mp_limb_t_swap__tmp;                \
 611   } while (0)
 612 #define MP_SIZE_T_SWAP(x, y)                    \
 613   do {                                          \
 614     mp_size_t __mp_size_t_swap__tmp = (x);      \
 615     (x) = (y);                                  \
 616     (y) = __mp_size_t_swap__tmp;                \
 617   } while (0)
 618
 619 #define MP_PTR_SWAP(x, y)               \
 620   do {                                  \
 621     mp_ptr __mp_ptr_swap__tmp = (x);    \
 622     (x) = (y);                          \
 623     (y) = __mp_ptr_swap__tmp;           \
 624   } while (0)
 625 #define MP_SRCPTR_SWAP(x, y)                    \
 626   do {                                          \
 627     mp_srcptr __mp_srcptr_swap__tmp = (x);      \
 628     (x) = (y);                                  \
 629     (y) = __mp_srcptr_swap__tmp;                \
 630   } while (0)
 631
 632 #define MPN_PTR_SWAP(xp,xs, yp,ys)      \
 633   do {                                  \
 634     MP_PTR_SWAP (xp, yp);               \
 635     MP_SIZE_T_SWAP (xs, ys);            \
 636   } while(0)
 637 #define MPN_SRCPTR_SWAP(xp,xs, yp,ys)   \
 638   do {                                  \
 639     MP_SRCPTR_SWAP (xp, yp);            \
 640     MP_SIZE_T_SWAP (xs, ys);            \
 641   } while(0)
 642
 643 #define MPZ_PTR_SWAP(x, y)              \
 644   do {                                  \
 645     mpz_ptr __mpz_ptr_swap__tmp = (x);  \
 646     (x) = (y);                          \
 647     (y) = __mpz_ptr_swap__tmp;          \
 648   } while (0)
 649 #define MPZ_SRCPTR_SWAP(x, y)                   \
 650   do {                                          \
 651     mpz_srcptr __mpz_srcptr_swap__tmp = (x);    \
 652     (x) = (y);                                  \
 653     (y) = __mpz_srcptr_swap__tmp;               \
 654   } while (0)
 655
 656
 657 /* Enhancement: __gmp_allocate_func could have "__attribute__ ((malloc))",
 658    but current gcc (3.0) doesn't seem to support that.  */
 659 __GMP_DECLSPEC extern void * (*__gmp_allocate_func) __GMP_PROTO ((size_t));
 660 __GMP_DECLSPEC extern void * (*__gmp_reallocate_func) __GMP_PROTO ((void *, size_t, size_t));
 661 __GMP_DECLSPEC extern void   (*__gmp_free_func) __GMP_PROTO ((void *, size_t));
 662
 663 __GMP_DECLSPEC void *__gmp_default_allocate __GMP_PROTO ((size_t));
 664 __GMP_DECLSPEC void *__gmp_default_reallocate __GMP_PROTO ((void *, size_t, size_t));
 665 __GMP_DECLSPEC void __gmp_default_free __GMP_PROTO ((void *, size_t));
 666
 667 #define __GMP_ALLOCATE_FUNC_TYPE(n,type) \
 668   ((type *) (*__gmp_allocate_func) ((n) * sizeof (type)))
 669 #define __GMP_ALLOCATE_FUNC_LIMBS(n)   __GMP_ALLOCATE_FUNC_TYPE (n, mp_limb_t)
 670
 671 #define __GMP_REALLOCATE_FUNC_TYPE(p, old_size, new_size, type) \
 672   ((type *) (*__gmp_reallocate_func)                            \
 673    (p, (old_size) * sizeof (type), (new_size) * sizeof (type)))
 674 #define __GMP_REALLOCATE_FUNC_LIMBS(p, old_size, new_size) \
 675   __GMP_REALLOCATE_FUNC_TYPE(p, old_size, new_size, mp_limb_t)
 676
 677 #define __GMP_FREE_FUNC_TYPE(p,n,type) (*__gmp_free_func) (p, (n) * sizeof (type))
 678 #define __GMP_FREE_FUNC_LIMBS(p,n)     __GMP_FREE_FUNC_TYPE (p, n, mp_limb_t)
 679
 680 #define __GMP_REALLOCATE_FUNC_MAYBE(ptr, oldsize, newsize)      \
 681   do {                                                          \
 682     if ((oldsize) != (newsize))                                 \
 683       (ptr) = (*__gmp_reallocate_func) (ptr, oldsize, newsize); \
 684   } while (0)
 685
 686 #define __GMP_REALLOCATE_FUNC_MAYBE_TYPE(ptr, oldsize, newsize, type)   \
 687   do {                                                                  \
 688     if ((oldsize) != (newsize))                                         \
 689       (ptr) = (type *) (*__gmp_reallocate_func)                         \
 690         (ptr, (oldsize) * sizeof (type), (newsize) * sizeof (type));    \
 691   } while (0)
 692
 693
 694 /* Dummy for non-gcc, code involving it will go dead. */
 695 #if ! defined (__GNUC__) || __GNUC__ < 2
 696 #define __builtin_constant_p(x)   0
 697 #endif
 698
 699
 700 /* In gcc 2.96 and up on i386, tail calls are optimized to jumps if the
 701    stack usage is compatible.  __attribute__ ((regparm (N))) helps by
 702    putting leading parameters in registers, avoiding extra stack.
 703
 704    regparm cannot be used with calls going through the PLT, because the
 705    binding code there may clobber the registers (%eax, %edx, %ecx) used for
 706    the regparm parameters.  Calls to local (ie. static) functions could
 707    still use this, if we cared to differentiate locals and globals.
 708
 709    On athlon-unknown-freebsd4.9 with gcc 3.3.3, regparm cannot be used with
 710    -p or -pg profiling, since that version of gcc doesn't realize the
 711    .mcount calls will clobber the parameter registers.  Other systems are
 712    ok, like debian with glibc 2.3.2 (mcount doesn't clobber), but we don't
 713    bother to try to detect this.  regparm is only an optimization so we just
 714    disable it when profiling (profiling being a slowdown anyway).  */
 715
 716 #if HAVE_HOST_CPU_FAMILY_x86 && __GMP_GNUC_PREREQ (2,96) && ! defined (PIC) \
 717   && ! WANT_PROFILING_PROF && ! WANT_PROFILING_GPROF
 718 #define USE_LEADING_REGPARM 1
 719 #else
 720 #define USE_LEADING_REGPARM 0
 721 #endif
 722
 723 /* Macros for altering parameter order according to regparm usage. */
 724 #if USE_LEADING_REGPARM
 725 #define REGPARM_2_1(a,b,x)    x,a,b
 726 #define REGPARM_3_1(a,b,c,x)  x,a,b,c
 727 #define REGPARM_ATTR(n) __attribute__ ((regparm (n)))
 728 #else
 729 #define REGPARM_2_1(a,b,x)    a,b,x
 730 #define REGPARM_3_1(a,b,c,x)  a,b,c,x
 731 #define REGPARM_ATTR(n)
 732 #endif
 733
 734
 735 /* ASM_L gives a local label for a gcc asm block, for use when temporary
 736    local labels like "1:" might not be available, which is the case for
 737    instance on the x86s (the SCO assembler doesn't support them).
 738
 739    The label generated is made unique by including "%=" which is a unique
 740    number for each insn.  This ensures the same name can be used in multiple
 741    asm blocks, perhaps via a macro.  Since jumps between asm blocks are not
 742    allowed there's no need for a label to be usable outside a single
 743    block.  */
 744
 745 #define ASM_L(name)  LSYM_PREFIX "asm_%=_" #name
 746
 747
 748 #if defined (__GNUC__) && HAVE_HOST_CPU_FAMILY_x86
 749 #if 0
 750 /* FIXME: Check that these actually improve things.
 751    FIXME: Need a cld after each std.
 752    FIXME: Can't have inputs in clobbered registers, must describe them as
 753    dummy outputs, and add volatile. */
 754 #define MPN_COPY_INCR(DST, SRC, N)                                      \
 755   __asm__ ("cld\n\trep\n\tmovsl" : :                                    \
 756            "D" (DST), "S" (SRC), "c" (N) :                              \
 757            "cx", "di", "si", "memory")
 758 #define MPN_COPY_DECR(DST, SRC, N)                                      \
 759   __asm__ ("std\n\trep\n\tmovsl" : :                                    \
 760            "D" ((DST) + (N) - 1), "S" ((SRC) + (N) - 1), "c" (N) :      \
 761            "cx", "di", "si", "memory")
 762 #endif
 763 #endif
 764
 765
 766 __GMP_DECLSPEC void __gmpz_aorsmul_1 __GMP_PROTO ((REGPARM_3_1 (mpz_ptr, mpz_srcptr, mp_limb_t, mp_size_t))) REGPARM_ATTR(1);
 767 #define mpz_aorsmul_1(w,u,v,sub)  __gmpz_aorsmul_1 (REGPARM_3_1 (w, u, v, sub))
 768
 769 #define mpz_n_pow_ui __gmpz_n_pow_ui
 770 __GMP_DECLSPEC void    mpz_n_pow_ui __GMP_PROTO ((mpz_ptr, mp_srcptr, mp_size_t, unsigned long));
 771
 772
 773 #define mpn_addmul_1c __MPN(addmul_1c)
 774 __GMP_DECLSPEC mp_limb_t mpn_addmul_1c __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t));
 775
 776 #define mpn_addmul_2 __MPN(addmul_2)
 777 __GMP_DECLSPEC mp_limb_t mpn_addmul_2 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 778
 779 #define mpn_addmul_3 __MPN(addmul_3)
 780 __GMP_DECLSPEC mp_limb_t mpn_addmul_3 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 781
 782 #define mpn_addmul_4 __MPN(addmul_4)
 783 __GMP_DECLSPEC mp_limb_t mpn_addmul_4 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 784
 785 #define mpn_addmul_5 __MPN(addmul_5)
 786 __GMP_DECLSPEC mp_limb_t mpn_addmul_5 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 787
 788 #define mpn_addmul_6 __MPN(addmul_6)
 789 __GMP_DECLSPEC mp_limb_t mpn_addmul_6 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 790
 791 #define mpn_addmul_7 __MPN(addmul_7)
 792 __GMP_DECLSPEC mp_limb_t mpn_addmul_7 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 793
 794 #define mpn_addmul_8 __MPN(addmul_8)
 795 __GMP_DECLSPEC mp_limb_t mpn_addmul_8 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 796
 797 /* mpn_addlsh1_n(c,a,b,n), when it exists, sets {c,n} to {a,n}+2*{b,n}, and
 798    returns the carry out (0, 1 or 2).  */
 799 #define mpn_addlsh1_n __MPN(addlsh1_n)
 800 __GMP_DECLSPEC mp_limb_t mpn_addlsh1_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 801
 802 /* mpn_addlsh2_n(c,a,b,n), when it exists, sets {c,n} to {a,n}+4*{b,n}, and
 803    returns the carry out (0, ..., 4).  */
 804 #define mpn_addlsh2_n __MPN(addlsh2_n)
 805 __GMP_DECLSPEC mp_limb_t mpn_addlsh2_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 806
 807 /* mpn_addlsh_n(c,a,b,n,k), when it exists, sets {c,n} to {a,n}+2^k*{b,n}, and
 808    returns the carry out (0, ..., 2^k).  */
 809 #define mpn_addlsh_n __MPN(addlsh_n)
 810   __GMP_DECLSPEC mp_limb_t mpn_addlsh_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, unsigned int));
 811
 812 /* mpn_sublsh1_n(c,a,b,n), when it exists, sets {c,n} to {a,n}-2*{b,n}, and
 813    returns the borrow out (0, 1 or 2).  */
 814 #define mpn_sublsh1_n __MPN(sublsh1_n)
 815 __GMP_DECLSPEC mp_limb_t mpn_sublsh1_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 816
 817 /* mpn_rsblsh1_n(c,a,b,n), when it exists, sets {c,n} to 2*{b,n}-{a,n}, and
 818    returns the carry out (-1, 0, 1).  */
 819 #define mpn_rsblsh1_n __MPN(rsblsh1_n)
 820 __GMP_DECLSPEC mp_limb_signed_t mpn_rsblsh1_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 821
 822 /* mpn_sublsh2_n(c,a,b,n), when it exists, sets {c,n} to {a,n}-4*{b,n}, and
 823    returns the borrow out (FIXME 0, 1, 2 or 3).  */
 824 #define mpn_sublsh2_n __MPN(sublsh2_n)
 825 __GMP_DECLSPEC mp_limb_t mpn_sublsh2_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 826
 827 /* mpn_rsblsh2_n(c,a,b,n), when it exists, sets {c,n} to 4*{b,n}-{a,n}, and
 828    returns the carry out (-1, ..., 3).  */
 829 #define mpn_rsblsh2_n __MPN(rsblsh2_n)
 830 __GMP_DECLSPEC mp_limb_signed_t mpn_rsblsh2_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 831
 832 /* mpn_rsblsh_n(c,a,b,n,k), when it exists, sets {c,n} to 2^k*{b,n}-{a,n}, and
 833    returns the carry out (-1, 0, ..., 2^k-1).  */
 834 #define mpn_rsblsh_n __MPN(rsblsh_n)
 835 __GMP_DECLSPEC mp_limb_signed_t mpn_rsblsh_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, unsigned int));
 836
 837 /* mpn_rsh1add_n(c,a,b,n), when it exists, sets {c,n} to ({a,n} + {b,n}) >> 1,
 838    and returns the bit rshifted out (0 or 1).  */
 839 #define mpn_rsh1add_n __MPN(rsh1add_n)
 840 __GMP_DECLSPEC mp_limb_t mpn_rsh1add_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 841 #define mpn_rsh1add_nc __MPN(rsh1add_nc)
 842 __GMP_DECLSPEC mp_limb_t mpn_rsh1add_nc __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, mp_limb_t));
 843
 844 /* mpn_rsh1sub_n(c,a,b,n), when it exists, sets {c,n} to ({a,n} - {b,n}) >> 1,
 845    and returns the bit rshifted out (0 or 1).  If there's a borrow from the
 846    subtract, it's stored as a 1 in the high bit of c[n-1], like a twos
 847    complement negative.  */
 848 #define mpn_rsh1sub_n __MPN(rsh1sub_n)
 849 __GMP_DECLSPEC mp_limb_t mpn_rsh1sub_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 850 #define mpn_rsh1sub_nc __MPN(rsh1sub_nc)
 851 __GMP_DECLSPEC mp_limb_t mpn_rsh1sub_nc __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, mp_limb_t));
 852
 853 #define mpn_lshiftc __MPN(lshiftc)
 854 __GMP_DECLSPEC mp_limb_t mpn_lshiftc __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, unsigned int));
 855
 856 #define mpn_add_n_sub_n __MPN(add_n_sub_n)
 857 __GMP_DECLSPEC mp_limb_t mpn_add_n_sub_n __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 858
 859 #define mpn_add_n_sub_nc __MPN(add_n_sub_nc)
 860 __GMP_DECLSPEC mp_limb_t mpn_add_n_sub_nc __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, mp_limb_t));
 861
 862 #define mpn_addaddmul_1msb0 __MPN(addaddmul_1msb0)
 863 __GMP_DECLSPEC mp_limb_t mpn_addaddmul_1msb0 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t));
 864
 865 #define mpn_divrem_1c __MPN(divrem_1c)
 866 __GMP_DECLSPEC mp_limb_t mpn_divrem_1c __GMP_PROTO ((mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t));
 867
 868 #define mpn_dump __MPN(dump)
 869 __GMP_DECLSPEC void mpn_dump __GMP_PROTO ((mp_srcptr, mp_size_t));
 870
 871 #define mpn_fib2_ui __MPN(fib2_ui)
 872 __GMP_DECLSPEC mp_size_t mpn_fib2_ui __GMP_PROTO ((mp_ptr, mp_ptr, unsigned long));
 873
 874 /* Remap names of internal mpn functions.  */
 875 #define __clz_tab               __MPN(clz_tab)
 876 #define mpn_udiv_w_sdiv         __MPN(udiv_w_sdiv)
 877
 878 #define mpn_jacobi_base __MPN(jacobi_base)
 879 __GMP_DECLSPEC int mpn_jacobi_base __GMP_PROTO ((mp_limb_t, mp_limb_t, int)) ATTRIBUTE_CONST;
 880
 881 #define mpn_mod_1c __MPN(mod_1c)
 882 __GMP_DECLSPEC mp_limb_t mpn_mod_1c __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t)) __GMP_ATTRIBUTE_PURE;
 883
 884 #define mpn_mul_1c __MPN(mul_1c)
 885 __GMP_DECLSPEC mp_limb_t mpn_mul_1c __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t));
 886
 887 #define mpn_mul_2 __MPN(mul_2)
 888 __GMP_DECLSPEC mp_limb_t mpn_mul_2 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 889
 890 #define mpn_mul_3 __MPN(mul_3)
 891 __GMP_DECLSPEC mp_limb_t mpn_mul_3 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 892
 893 #define mpn_mul_4 __MPN(mul_4)
 894 __GMP_DECLSPEC mp_limb_t mpn_mul_4 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 895
 896 #ifndef mpn_mul_basecase  /* if not done with cpuvec in a fat binary */
 897 #define mpn_mul_basecase __MPN(mul_basecase)
 898 __GMP_DECLSPEC void mpn_mul_basecase __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t));
 899 #endif
 900
 901 #define mpn_mullo_n __MPN(mullo_n)
 902 __GMP_DECLSPEC void mpn_mullo_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 903
 904 #define mpn_mullo_basecase __MPN(mullo_basecase)
 905 __GMP_DECLSPEC void mpn_mullo_basecase __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t));
 906
 907 #define mpn_sqr __MPN(sqr)
 908 __GMP_DECLSPEC void mpn_sqr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t));
 909
 910 #ifndef mpn_sqr_basecase  /* if not done with cpuvec in a fat binary */
 911 #define mpn_sqr_basecase __MPN(sqr_basecase)
 912 __GMP_DECLSPEC void mpn_sqr_basecase __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t));
 913 #endif
 914
 915 #define mpn_submul_1c __MPN(submul_1c)
 916 __GMP_DECLSPEC mp_limb_t mpn_submul_1c __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t));
 917
 918 #define mpn_redc_1 __MPN(redc_1)
 919 __GMP_DECLSPEC void mpn_redc_1 __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_limb_t));
 920
 921 #define mpn_redc_2 __MPN(redc_2)
 922 __GMP_DECLSPEC void mpn_redc_2 __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 923 #define mpn_redc_n __MPN(redc_n)
 924 __GMP_DECLSPEC void mpn_redc_n __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_srcptr));
 925
 926
 927 #define mpn_mod_1_1p_cps __MPN(mod_1_1p_cps)
 928 __GMP_DECLSPEC void mpn_mod_1_1p_cps __GMP_PROTO ((mp_limb_t [4], mp_limb_t));
 929 #define mpn_mod_1_1p __MPN(mod_1_1p)
 930 __GMP_DECLSPEC mp_limb_t mpn_mod_1_1p __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t [4])) __GMP_ATTRIBUTE_PURE;
 931
 932 #define mpn_mod_1s_2p_cps __MPN(mod_1s_2p_cps)
 933 __GMP_DECLSPEC void mpn_mod_1s_2p_cps __GMP_PROTO ((mp_limb_t [5], mp_limb_t));
 934 #define mpn_mod_1s_2p __MPN(mod_1s_2p)
 935 __GMP_DECLSPEC mp_limb_t mpn_mod_1s_2p __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t [5])) __GMP_ATTRIBUTE_PURE;
 936
 937 #define mpn_mod_1s_3p_cps __MPN(mod_1s_3p_cps)
 938 __GMP_DECLSPEC void mpn_mod_1s_3p_cps __GMP_PROTO ((mp_limb_t [6], mp_limb_t));
 939 #define mpn_mod_1s_3p __MPN(mod_1s_3p)
 940 __GMP_DECLSPEC mp_limb_t mpn_mod_1s_3p __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t [6])) __GMP_ATTRIBUTE_PURE;
 941
 942 #define mpn_mod_1s_4p_cps __MPN(mod_1s_4p_cps)
 943 __GMP_DECLSPEC void mpn_mod_1s_4p_cps __GMP_PROTO ((mp_limb_t [7], mp_limb_t));
 944 #define mpn_mod_1s_4p __MPN(mod_1s_4p)
 945 __GMP_DECLSPEC mp_limb_t mpn_mod_1s_4p __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t [7])) __GMP_ATTRIBUTE_PURE;
 946
 947 #define mpn_bc_mulmod_bnm1 __MPN(bc_mulmod_bnm1)
 948 __GMP_DECLSPEC void mpn_bc_mulmod_bnm1 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, mp_ptr));
 949 #define mpn_mulmod_bnm1 __MPN(mulmod_bnm1)
 950 __GMP_DECLSPEC void mpn_mulmod_bnm1 __GMP_PROTO ((mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
 951 #define mpn_mulmod_bnm1_next_size __MPN(mulmod_bnm1_next_size)
 952 __GMP_DECLSPEC mp_size_t mpn_mulmod_bnm1_next_size __GMP_PROTO ((mp_size_t)) ATTRIBUTE_CONST;
 953 static inline mp_size_t
 954 mpn_mulmod_bnm1_itch (mp_size_t rn, mp_size_t an, mp_size_t bn) {
 955   mp_size_t n, itch;
 956   n = rn >> 1;
 957   itch = rn + 4 +
 958     (an > n ? (bn > n ? rn : n) : 0);
 959   return itch;
 960 }
 961
 962 #define mpn_sqrmod_bnm1 __MPN(sqrmod_bnm1)
 963 __GMP_DECLSPEC void mpn_sqrmod_bnm1 __GMP_PROTO ((mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
 964 #define mpn_sqrmod_bnm1_next_size __MPN(sqrmod_bnm1_next_size)
 965 __GMP_DECLSPEC mp_size_t mpn_sqrmod_bnm1_next_size __GMP_PROTO ((mp_size_t)) ATTRIBUTE_CONST;
 966 static inline mp_size_t
 967 mpn_sqrmod_bnm1_itch (mp_size_t rn, mp_size_t an) {
 968   mp_size_t n, itch;
 969   n = rn >> 1;
 970   itch = rn + 3 +
 971     (an > n ? an : 0);
 972   return itch;
 973 }
 974
 975 typedef __gmp_randstate_struct *gmp_randstate_ptr;
 976 typedef const __gmp_randstate_struct *gmp_randstate_srcptr;
 977
 978 /* Pseudo-random number generator function pointers structure.  */
 979 typedef struct {
 980   void (*randseed_fn) __GMP_PROTO ((gmp_randstate_t, mpz_srcptr));
 981   void (*randget_fn) __GMP_PROTO ((gmp_randstate_t, mp_ptr, unsigned long int));
 982   void (*randclear_fn) __GMP_PROTO ((gmp_randstate_t));
 983   void (*randiset_fn) __GMP_PROTO ((gmp_randstate_ptr, gmp_randstate_srcptr));
 984 } gmp_randfnptr_t;
 985
 986 /* Macro to obtain a void pointer to the function pointers structure.  */
 987 #define RNG_FNPTR(rstate) ((rstate)->_mp_algdata._mp_lc)
 988
 989 /* Macro to obtain a pointer to the generator's state.
 990    When used as a lvalue the rvalue needs to be cast to mp_ptr.  */
 991 #define RNG_STATE(rstate) ((rstate)->_mp_seed->_mp_d)
 992
 993 /* Write a given number of random bits to rp.  */
 994 #define _gmp_rand(rp, state, bits)                              \
 995   do {                                                          \
 996     gmp_randstate_ptr  __rstate = (state);                      \
 997     (*((gmp_randfnptr_t *) RNG_FNPTR (__rstate))->randget_fn)   \
 998        (__rstate, rp, bits);                                    \
 999   } while (0)
1000
1001 __GMP_DECLSPEC void __gmp_randinit_mt_noseed __GMP_PROTO ((gmp_randstate_t));
1002
1003
1004 /* __gmp_rands is the global state for the old-style random functions, and
1005    is also used in the test programs (hence the __GMP_DECLSPEC).
1006
1007    There's no seeding here, so mpz_random etc will generate the same
1008    sequence every time.  This is not unlike the C library random functions
1009    if you don't seed them, so perhaps it's acceptable.  Digging up a seed
1010    from /dev/random or the like would work on many systems, but might
1011    encourage a false confidence, since it'd be pretty much impossible to do
1012    something that would work reliably everywhere.  In any case the new style
1013    functions are recommended to applications which care about randomness, so
1014    the old functions aren't too important.  */
1015
1016 __GMP_DECLSPEC extern char             __gmp_rands_initialized;
1017 __GMP_DECLSPEC extern gmp_randstate_t  __gmp_rands;
1018
1019 #define RANDS                                       \
1020   ((__gmp_rands_initialized ? 0                     \
1021     : (__gmp_rands_initialized = 1,                 \
1022        __gmp_randinit_mt_noseed (__gmp_rands), 0)), \
1023    __gmp_rands)
1024
1025 /* this is used by the test programs, to free memory */
1026 #define RANDS_CLEAR()                   \
1027   do {                                  \
1028     if (__gmp_rands_initialized)        \
1029       {                                 \
1030         __gmp_rands_initialized = 0;    \
1031         gmp_randclear (__gmp_rands);    \
1032       }                                 \
1033   } while (0)
1034
1035
1036 /* For a threshold between algorithms A and B, size>=thresh is where B
1037    should be used.  Special value MP_SIZE_T_MAX means only ever use A, or
1038    value 0 means only ever use B.  The tests for these special values will
1039    be compile-time constants, so the compiler should be able to eliminate
1040    the code for the unwanted algorithm.  */
1041
1042 #define ABOVE_THRESHOLD(size,thresh)    \
1043   ((thresh) == 0                        \
1044    || ((thresh) != MP_SIZE_T_MAX        \
1045        && (size) >= (thresh)))
1046 #define BELOW_THRESHOLD(size,thresh)  (! ABOVE_THRESHOLD (size, thresh))
1047
1048 #define MPN_TOOM22_MUL_MINSIZE    4
1049 #define MPN_TOOM2_SQR_MINSIZE     4
1050
1051 #define MPN_TOOM33_MUL_MINSIZE   17
1052 #define MPN_TOOM3_SQR_MINSIZE    17
1053
1054 #define MPN_TOOM44_MUL_MINSIZE   30
1055 #define MPN_TOOM4_SQR_MINSIZE    30
1056
1057 #define MPN_TOOM6H_MUL_MINSIZE   46
1058 #define MPN_TOOM6_SQR_MINSIZE    46
1059
1060 #define MPN_TOOM8H_MUL_MINSIZE   86
1061 #define MPN_TOOM8_SQR_MINSIZE    86
1062
1063 #define MPN_TOOM32_MUL_MINSIZE   10
1064 #define MPN_TOOM42_MUL_MINSIZE   10
1065 #define MPN_TOOM43_MUL_MINSIZE   49 /* ??? */
1066 #define MPN_TOOM53_MUL_MINSIZE   49 /* ??? */
1067 #define MPN_TOOM63_MUL_MINSIZE   49
1068
1069 #define   mpn_sqr_diagonal __MPN(sqr_diagonal)
1070 __GMP_DECLSPEC void      mpn_sqr_diagonal __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t));
1071
1072 #define   mpn_toom_interpolate_5pts __MPN(toom_interpolate_5pts)
1073 __GMP_DECLSPEC void      mpn_toom_interpolate_5pts __GMP_PROTO ((mp_ptr, mp_ptr, mp_ptr, mp_size_t, mp_size_t, int, mp_limb_t));
1074
1075 enum toom6_flags {toom6_all_pos = 0, toom6_vm1_neg = 1, toom6_vm2_neg = 2};
1076 #define   mpn_toom_interpolate_6pts __MPN(toom_interpolate_6pts)
1077 __GMP_DECLSPEC void      mpn_toom_interpolate_6pts __GMP_PROTO ((mp_ptr, mp_size_t, enum toom6_flags, mp_ptr, mp_ptr, mp_ptr, mp_size_t));
1078
1079 enum toom7_flags { toom7_w1_neg = 1, toom7_w3_neg = 2 };
1080 #define   mpn_toom_interpolate_7pts __MPN(toom_interpolate_7pts)
1081 __GMP_DECLSPEC void      mpn_toom_interpolate_7pts __GMP_PROTO ((mp_ptr, mp_size_t, enum toom7_flags, mp_ptr, mp_ptr, mp_ptr, mp_ptr, mp_size_t, mp_ptr));
1082
1083 #define mpn_toom_interpolate_8pts __MPN(toom_interpolate_8pts)
1084 __GMP_DECLSPEC void      mpn_toom_interpolate_8pts __GMP_PROTO ((mp_ptr, mp_size_t, mp_ptr, mp_ptr, mp_size_t, mp_ptr));
1085
1086 #define mpn_toom_interpolate_12pts __MPN(toom_interpolate_12pts)
1087 __GMP_DECLSPEC void      mpn_toom_interpolate_12pts __GMP_PROTO ((mp_ptr, mp_ptr, mp_ptr, mp_ptr, mp_size_t, mp_size_t, int, mp_ptr));
1088
1089 #define mpn_toom_interpolate_16pts __MPN(toom_interpolate_16pts)
1090 __GMP_DECLSPEC void      mpn_toom_interpolate_16pts __GMP_PROTO ((mp_ptr, mp_ptr, mp_ptr, mp_ptr, mp_ptr, mp_size_t, mp_size_t, int, mp_ptr));
1091
1092 #define   mpn_toom_couple_handling __MPN(toom_couple_handling)
1093 __GMP_DECLSPEC void mpn_toom_couple_handling __GMP_PROTO ((mp_ptr, mp_size_t, mp_ptr, int, mp_size_t, int, int));
1094
1095 #define   mpn_toom_eval_dgr3_pm1 __MPN(toom_eval_dgr3_pm1)
1096 __GMP_DECLSPEC int mpn_toom_eval_dgr3_pm1 __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_size_t, mp_ptr));
1097
1098 #define   mpn_toom_eval_dgr3_pm2 __MPN(toom_eval_dgr3_pm2)
1099 __GMP_DECLSPEC int mpn_toom_eval_dgr3_pm2 __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_size_t, mp_ptr));
1100
1101 #define   mpn_toom_eval_pm1 __MPN(toom_eval_pm1)
1102 __GMP_DECLSPEC int mpn_toom_eval_pm1 __GMP_PROTO ((mp_ptr, mp_ptr, unsigned, mp_srcptr, mp_size_t, mp_size_t, mp_ptr));
1103
1104 #define   mpn_toom_eval_pm2 __MPN(toom_eval_pm2)
1105 __GMP_DECLSPEC int mpn_toom_eval_pm2 __GMP_PROTO ((mp_ptr, mp_ptr, unsigned, mp_srcptr, mp_size_t, mp_size_t, mp_ptr));
1106
1107 #define   mpn_toom_eval_pm2exp __MPN(toom_eval_pm2exp)
1108 __GMP_DECLSPEC int mpn_toom_eval_pm2exp __GMP_PROTO ((mp_ptr, mp_ptr, unsigned, mp_srcptr, mp_size_t, mp_size_t, unsigned, mp_ptr));
1109
1110 #define   mpn_toom_eval_pm2rexp __MPN(toom_eval_pm2rexp)
1111 __GMP_DECLSPEC int mpn_toom_eval_pm2rexp __GMP_PROTO ((mp_ptr, mp_ptr, unsigned, mp_srcptr, mp_size_t, mp_size_t, unsigned, mp_ptr));
1112
1113 #define   mpn_toom22_mul __MPN(toom22_mul)
1114 __GMP_DECLSPEC void      mpn_toom22_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1115
1116 #define   mpn_toom32_mul __MPN(toom32_mul)
1117 __GMP_DECLSPEC void      mpn_toom32_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1118
1119 #define   mpn_toom42_mul __MPN(toom42_mul)
1120 __GMP_DECLSPEC void      mpn_toom42_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1121
1122 #define   mpn_toom52_mul __MPN(toom52_mul)
1123 __GMP_DECLSPEC void      mpn_toom52_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1124
1125 #define   mpn_toom62_mul __MPN(toom62_mul)
1126 __GMP_DECLSPEC void      mpn_toom62_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1127
1128 #define   mpn_toom2_sqr __MPN(toom2_sqr)
1129 __GMP_DECLSPEC void      mpn_toom2_sqr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
1130
1131 #define   mpn_toom33_mul __MPN(toom33_mul)
1132 __GMP_DECLSPEC void      mpn_toom33_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1133
1134 #define   mpn_toom43_mul __MPN(toom43_mul)
1135 __GMP_DECLSPEC void      mpn_toom43_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1136
1137 #define   mpn_toom53_mul __MPN(toom53_mul)
1138 __GMP_DECLSPEC void      mpn_toom53_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1139
1140 #define   mpn_toom63_mul __MPN(toom63_mul)
1141 __GMP_DECLSPEC void      mpn_toom63_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1142
1143 #define   mpn_toom3_sqr __MPN(toom3_sqr)
1144 __GMP_DECLSPEC void      mpn_toom3_sqr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
1145
1146 #define   mpn_toom44_mul __MPN(toom44_mul)
1147 __GMP_DECLSPEC void      mpn_toom44_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1148
1149 #define   mpn_toom4_sqr __MPN(toom4_sqr)
1150 __GMP_DECLSPEC void      mpn_toom4_sqr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
1151
1152 #define   mpn_toom6h_mul __MPN(toom6h_mul)
1153 __GMP_DECLSPEC void      mpn_toom6h_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1154
1155 #define   mpn_toom6_sqr __MPN(toom6_sqr)
1156 __GMP_DECLSPEC void      mpn_toom6_sqr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
1157
1158 #define   mpn_toom8h_mul __MPN(toom8h_mul)
1159 __GMP_DECLSPEC void      mpn_toom8h_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1160
1161 #define   mpn_toom8_sqr __MPN(toom8_sqr)
1162 __GMP_DECLSPEC void      mpn_toom8_sqr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
1163
1164 #define   mpn_fft_best_k __MPN(fft_best_k)
1165 __GMP_DECLSPEC int       mpn_fft_best_k __GMP_PROTO ((mp_size_t, int)) ATTRIBUTE_CONST;
1166
1167 #define   mpn_mul_fft __MPN(mul_fft)
1168 __GMP_DECLSPEC mp_limb_t mpn_mul_fft __GMP_PROTO ((mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, int));
1169
1170 #define   mpn_mul_fft_full __MPN(mul_fft_full)
1171 __GMP_DECLSPEC void      mpn_mul_fft_full __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t));
1172
1173 #define   mpn_nussbaumer_mul __MPN(nussbaumer_mul)
1174 __GMP_DECLSPEC void      mpn_nussbaumer_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t));
1175
1176 #define   mpn_fft_next_size __MPN(fft_next_size)
1177 __GMP_DECLSPEC mp_size_t mpn_fft_next_size __GMP_PROTO ((mp_size_t, int)) ATTRIBUTE_CONST;
1178
1179 #define   mpn_sbpi1_div_qr __MPN(sbpi1_div_qr)
1180 __GMP_DECLSPEC mp_limb_t mpn_sbpi1_div_qr __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t));
1181
1182 #define   mpn_sbpi1_div_q __MPN(sbpi1_div_q)
1183 __GMP_DECLSPEC mp_limb_t mpn_sbpi1_div_q __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t));
1184
1185 #define   mpn_sbpi1_divappr_q __MPN(sbpi1_divappr_q)
1186 __GMP_DECLSPEC mp_limb_t mpn_sbpi1_divappr_q __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t));
1187
1188 #define   mpn_dcpi1_div_qr __MPN(dcpi1_div_qr)
1189 __GMP_DECLSPEC mp_limb_t mpn_dcpi1_div_qr __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, gmp_pi1_t *));
1190 #define   mpn_dcpi1_div_qr_n __MPN(dcpi1_div_qr_n)
1191 __GMP_DECLSPEC mp_limb_t mpn_dcpi1_div_qr_n __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, gmp_pi1_t *, mp_ptr));
1192
1193 #define   mpn_dcpi1_div_q __MPN(dcpi1_div_q)
1194 __GMP_DECLSPEC mp_limb_t mpn_dcpi1_div_q __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, gmp_pi1_t *));
1195
1196 #define   mpn_dcpi1_divappr_q __MPN(dcpi1_divappr_q)
1197 __GMP_DECLSPEC mp_limb_t mpn_dcpi1_divappr_q __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, gmp_pi1_t *));
1198 #define   mpn_dcpi1_divappr_q_n __MPN(dcpi1_divappr_q_n)
1199 __GMP_DECLSPEC mp_limb_t mpn_dcpi1_divappr_q_n __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, gmp_pi1_t *, mp_ptr));
1200
1201 #define   mpn_mu_div_qr __MPN(mu_div_qr)
1202 __GMP_DECLSPEC mp_limb_t mpn_mu_div_qr __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1203 #define   mpn_mu_div_qr_itch __MPN(mu_div_qr_itch)
1204 __GMP_DECLSPEC mp_size_t mpn_mu_div_qr_itch __GMP_PROTO ((mp_size_t, mp_size_t, int));
1205 #define   mpn_mu_div_qr_choose_in __MPN(mu_div_qr_choose_in)
1206 __GMP_DECLSPEC mp_size_t mpn_mu_div_qr_choose_in __GMP_PROTO ((mp_size_t, mp_size_t, int));
1207
1208 #define   mpn_preinv_mu_div_qr __MPN(preinv_mu_div_qr)
1209 __GMP_DECLSPEC mp_limb_t mpn_preinv_mu_div_qr __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1210 #define   mpn_preinv_mu_div_qr_itch __MPN(preinv_mu_div_qr_itch)
1211 __GMP_DECLSPEC mp_size_t mpn_preinv_mu_div_qr_itch __GMP_PROTO ((mp_size_t, mp_size_t, mp_size_t));
1212
1213 #define   mpn_mu_divappr_q __MPN(mu_divappr_q)
1214 __GMP_DECLSPEC mp_limb_t mpn_mu_divappr_q __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1215 #define   mpn_mu_divappr_q_itch __MPN(mu_divappr_q_itch)
1216 __GMP_DECLSPEC mp_size_t mpn_mu_divappr_q_itch __GMP_PROTO ((mp_size_t, mp_size_t, int));
1217 #define   mpn_mu_divappr_q_choose_in __MPN(mu_divappr_q_choose_in)
1218 __GMP_DECLSPEC mp_size_t mpn_mu_divappr_q_choose_in __GMP_PROTO ((mp_size_t, mp_size_t, int));
1219
1220 #define   mpn_preinv_mu_divappr_q __MPN(preinv_mu_divappr_q)
1221 __GMP_DECLSPEC mp_limb_t mpn_preinv_mu_divappr_q __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1222
1223 #define   mpn_mu_div_q __MPN(mu_div_q)
1224 __GMP_DECLSPEC mp_limb_t mpn_mu_div_q __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1225 #define   mpn_mu_div_q_itch __MPN(mu_div_q_itch)
1226 __GMP_DECLSPEC mp_size_t mpn_mu_div_q_itch __GMP_PROTO ((mp_size_t, mp_size_t, int));
1227
1228 #define  mpn_div_q __MPN(div_q)
1229 __GMP_DECLSPEC void mpn_div_q __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1230
1231 #define   mpn_invert __MPN(invert)
1232 __GMP_DECLSPEC void      mpn_invert __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
1233 #define mpn_invert_itch(n)  mpn_invertappr_itch(n)
1234
1235 #define   mpn_ni_invertappr __MPN(ni_invertappr)
1236 __GMP_DECLSPEC mp_limb_t mpn_ni_invertappr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
1237 #define   mpn_invertappr __MPN(invertappr)
1238 __GMP_DECLSPEC mp_limb_t mpn_invertappr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
1239 #define mpn_invertappr_itch(n)  (3 * (n) + 2)
1240
1241 #define   mpn_binvert __MPN(binvert)
1242 __GMP_DECLSPEC void      mpn_binvert __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
1243 #define   mpn_binvert_itch __MPN(binvert_itch)
1244 __GMP_DECLSPEC mp_size_t mpn_binvert_itch __GMP_PROTO ((mp_size_t));
1245
1246 #define mpn_bdiv_q_1 __MPN(bdiv_q_1)
1247 __GMP_DECLSPEC mp_limb_t mpn_bdiv_q_1 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t));
1248
1249 #define mpn_pi1_bdiv_q_1 __MPN(pi1_bdiv_q_1)
1250 __GMP_DECLSPEC mp_limb_t mpn_pi1_bdiv_q_1 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t, int));
1251
1252 #define   mpn_sbpi1_bdiv_qr __MPN(sbpi1_bdiv_qr)
1253 __GMP_DECLSPEC mp_limb_t mpn_sbpi1_bdiv_qr __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t));
1254
1255 #define   mpn_sbpi1_bdiv_q __MPN(sbpi1_bdiv_q)
1256 __GMP_DECLSPEC void      mpn_sbpi1_bdiv_q __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t));
1257
1258 #define   mpn_dcpi1_bdiv_qr __MPN(dcpi1_bdiv_qr)
1259 __GMP_DECLSPEC mp_limb_t mpn_dcpi1_bdiv_qr __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t));
1260 #define   mpn_dcpi1_bdiv_qr_n_itch __MPN(dcpi1_bdiv_qr_n_itch)
1261 __GMP_DECLSPEC mp_size_t mpn_dcpi1_bdiv_qr_n_itch __GMP_PROTO ((mp_size_t));
1262
1263 #define   mpn_dcpi1_bdiv_qr_n __MPN(dcpi1_bdiv_qr_n)
1264 __GMP_DECLSPEC mp_limb_t mpn_dcpi1_bdiv_qr_n __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_limb_t, mp_ptr));
1265 #define   mpn_dcpi1_bdiv_q __MPN(dcpi1_bdiv_q)
1266 __GMP_DECLSPEC void      mpn_dcpi1_bdiv_q __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t));
1267
1268 #define   mpn_dcpi1_bdiv_q_n_itch __MPN(dcpi1_bdiv_q_n_itch)
1269 __GMP_DECLSPEC mp_size_t mpn_dcpi1_bdiv_q_n_itch __GMP_PROTO ((mp_size_t));
1270 #define   mpn_dcpi1_bdiv_q_n __MPN(dcpi1_bdiv_q_n)
1271 __GMP_DECLSPEC void      mpn_dcpi1_bdiv_q_n __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_limb_t, mp_ptr));
1272
1273 #define   mpn_mu_bdiv_qr __MPN(mu_bdiv_qr)
1274 __GMP_DECLSPEC mp_limb_t mpn_mu_bdiv_qr __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1275 #define   mpn_mu_bdiv_qr_itch __MPN(mu_bdiv_qr_itch)
1276 __GMP_DECLSPEC mp_size_t mpn_mu_bdiv_qr_itch __GMP_PROTO ((mp_size_t, mp_size_t));
1277
1278 #define   mpn_mu_bdiv_q __MPN(mu_bdiv_q)
1279 __GMP_DECLSPEC void      mpn_mu_bdiv_q __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1280 #define   mpn_mu_bdiv_q_itch __MPN(mu_bdiv_q_itch)
1281 __GMP_DECLSPEC mp_size_t mpn_mu_bdiv_q_itch __GMP_PROTO ((mp_size_t, mp_size_t));
1282
1283 #define   mpn_bdiv_qr __MPN(bdiv_qr)
1284 __GMP_DECLSPEC mp_limb_t mpn_bdiv_qr __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1285 #define   mpn_bdiv_qr_itch __MPN(bdiv_qr_itch)
1286 __GMP_DECLSPEC mp_size_t mpn_bdiv_qr_itch __GMP_PROTO ((mp_size_t, mp_size_t));
1287
1288 #define   mpn_bdiv_q __MPN(bdiv_q)
1289 __GMP_DECLSPEC void      mpn_bdiv_q __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1290 #define   mpn_bdiv_q_itch __MPN(bdiv_q_itch)
1291 __GMP_DECLSPEC mp_size_t mpn_bdiv_q_itch __GMP_PROTO ((mp_size_t, mp_size_t));
1292
1293 #define   mpn_divexact __MPN(divexact)
1294 __GMP_DECLSPEC void      mpn_divexact __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t));
1295 #define   mpn_divexact_itch __MPN(divexact_itch)
1296 __GMP_DECLSPEC mp_size_t mpn_divexact_itch __GMP_PROTO ((mp_size_t, mp_size_t));
1297
1298 #define   mpn_bdiv_dbm1c __MPN(bdiv_dbm1c)
1299 __GMP_DECLSPEC mp_limb_t mpn_bdiv_dbm1c __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t));
1300 #define   mpn_bdiv_dbm1(dst, src, size, divisor) \
1301   mpn_bdiv_dbm1c (dst, src, size, divisor, __GMP_CAST (mp_limb_t, 0))
1302
1303 #define   mpn_powm __MPN(powm)
1304 __GMP_DECLSPEC void      mpn_powm __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1305 #define   mpn_powlo __MPN(powlo)
1306 __GMP_DECLSPEC void      mpn_powlo __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, mp_size_t, mp_ptr));
1307 #define   mpn_powm_sec __MPN(powm_sec)
1308 __GMP_DECLSPEC void      mpn_powm_sec __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
1309 #define   mpn_powm_sec_itch __MPN(powm_sec_itch)
1310 __GMP_DECLSPEC mp_size_t mpn_powm_sec_itch __GMP_PROTO ((mp_size_t, mp_size_t, mp_size_t));
1311 #define   mpn_subcnd_n __MPN(subcnd_n)
1312 __GMP_DECLSPEC mp_limb_t mpn_subcnd_n __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, mp_limb_t));
1313 #define   mpn_tabselect __MPN(tabselect)
1314 __GMP_DECLSPEC void      mpn_tabselect __GMP_PROTO ((volatile mp_limb_t *, volatile mp_limb_t *, mp_size_t, mp_size_t, mp_size_t));
1315 #define mpn_redc_1_sec __MPN(redc_1_sec)
1316 __GMP_DECLSPEC void mpn_redc_1_sec __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_limb_t));
1317
1318 #ifndef DIVEXACT_BY3_METHOD
1319 #if GMP_NUMB_BITS % 2 == 0 && ! defined (HAVE_NATIVE_mpn_divexact_by3c)
1320 #define DIVEXACT_BY3_METHOD 0   /* default to using mpn_bdiv_dbm1c */
1321 #else
1322 #define DIVEXACT_BY3_METHOD 1
1323 #endif
1324 #endif
1325
1326 #if DIVEXACT_BY3_METHOD == 0
1327 #undef mpn_divexact_by3
1328 #define mpn_divexact_by3(dst,src,size) \
1329   (3 & mpn_bdiv_dbm1 (dst, src, size, __GMP_CAST (mp_limb_t, GMP_NUMB_MASK / 3)))
1330 /* override mpn_divexact_by3c defined in gmp.h */
1331 /*
1332 #undef mpn_divexact_by3c
1333 #define mpn_divexact_by3c(dst,src,size,cy) \
1334   (3 & mpn_bdiv_dbm1c (dst, src, size, __GMP_CAST (mp_limb_t, GMP_NUMB_MASK / 3, GMP_NUMB_MASK / 3 * cy)))
1335 */
1336 #endif
1337
1338 #if GMP_NUMB_BITS % 4 == 0
1339 #define mpn_divexact_by5(dst,src,size) \
1340   (7 & 3 * mpn_bdiv_dbm1 (dst, src, size, __GMP_CAST (mp_limb_t, GMP_NUMB_MASK / 5)))
1341 #endif
1342
1343 #if GMP_NUMB_BITS % 6 == 0
1344 #define mpn_divexact_by7(dst,src,size) \
1345   (7 & 1 * mpn_bdiv_dbm1 (dst, src, size, __GMP_CAST (mp_limb_t, GMP_NUMB_MASK / 7)))
1346 #endif
1347
1348 #if GMP_NUMB_BITS % 6 == 0
1349 #define mpn_divexact_by9(dst,src,size) \
1350   (15 & 7 * mpn_bdiv_dbm1 (dst, src, size, __GMP_CAST (mp_limb_t, GMP_NUMB_MASK / 9)))
1351 #endif
1352
1353 #if GMP_NUMB_BITS % 10 == 0
1354 #define mpn_divexact_by11(dst,src,size) \
1355   (15 & 5 * mpn_bdiv_dbm1 (dst, src, size, __GMP_CAST (mp_limb_t, GMP_NUMB_MASK / 11)))
1356 #endif
1357
1358 #if GMP_NUMB_BITS % 12 == 0
1359 #define mpn_divexact_by13(dst,src,size) \
1360   (15 & 3 * mpn_bdiv_dbm1 (dst, src, size, __GMP_CAST (mp_limb_t, GMP_NUMB_MASK / 13)))
1361 #endif
1362
1363 #if GMP_NUMB_BITS % 4 == 0
1364 #define mpn_divexact_by15(dst,src,size) \
1365   (15 & 1 * mpn_bdiv_dbm1 (dst, src, size, __GMP_CAST (mp_limb_t, GMP_NUMB_MASK / 15)))
1366 #endif
1367
1368 #define mpz_divexact_gcd  __gmpz_divexact_gcd
1369 __GMP_DECLSPEC void    mpz_divexact_gcd __GMP_PROTO ((mpz_ptr, mpz_srcptr, mpz_srcptr));
1370
1371 #define mpz_inp_str_nowhite __gmpz_inp_str_nowhite
1372 #ifdef _GMP_H_HAVE_FILE
1373 __GMP_DECLSPEC size_t  mpz_inp_str_nowhite __GMP_PROTO ((mpz_ptr, FILE *, int, int, size_t));
1374 #endif
1375
1376 #define mpn_divisible_p __MPN(divisible_p)
1377 __GMP_DECLSPEC int     mpn_divisible_p __GMP_PROTO ((mp_srcptr, mp_size_t, mp_srcptr, mp_size_t)) __GMP_ATTRIBUTE_PURE;
1378
1379 #define   mpn_rootrem __MPN(rootrem)
1380 __GMP_DECLSPEC mp_size_t mpn_rootrem __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_limb_t));
1381
1382
1383 #if defined (_CRAY)
1384 #define MPN_COPY_INCR(dst, src, n)                                      \
1385   do {                                                                  \
1386     int __i;            /* Faster on some Crays with plain int */       \
1387     _Pragma ("_CRI ivdep");                                             \
1388     for (__i = 0; __i < (n); __i++)                                     \
1389       (dst)[__i] = (src)[__i];                                          \
1390   } while (0)
1391 #endif
1392
1393 /* used by test programs, hence __GMP_DECLSPEC */
1394 #ifndef mpn_copyi  /* if not done with cpuvec in a fat binary */
1395 #define mpn_copyi __MPN(copyi)
1396 __GMP_DECLSPEC void mpn_copyi __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t));
1397 #endif
1398
1399 #if ! defined (MPN_COPY_INCR) && HAVE_NATIVE_mpn_copyi
1400 #define MPN_COPY_INCR(dst, src, size)                   \
1401   do {                                                  \
1402     ASSERT ((size) >= 0);                               \
1403     ASSERT (MPN_SAME_OR_INCR_P (dst, src, size));       \
1404     mpn_copyi (dst, src, size);                         \
1405   } while (0)
1406 #endif
1407
1408 /* Copy N limbs from SRC to DST incrementing, N==0 allowed.  */
1409 #if ! defined (MPN_COPY_INCR)
1410 #define MPN_COPY_INCR(dst, src, n)                      \
1411   do {                                                  \
1412     ASSERT ((n) >= 0);                                  \
1413     ASSERT (MPN_SAME_OR_INCR_P (dst, src, n));          \
1414     if ((n) != 0)                                       \
1415       {                                                 \
1416         mp_size_t __n = (n) - 1;                        \
1417         mp_ptr __dst = (dst);                           \
1418         mp_srcptr __src = (src);                        \
1419         mp_limb_t __x;                                  \
1420         __x = *__src++;                                 \
1421         if (__n != 0)                                   \
1422           {                                             \
1423             do                                          \
1424               {                                         \
1425                 *__dst++ = __x;                         \
1426                 __x = *__src++;                         \
1427               }                                         \
1428             while (--__n);                              \
1429           }                                             \
1430         *__dst++ = __x;                                 \
1431       }                                                 \
1432   } while (0)
1433 #endif
1434
1435
1436 #if defined (_CRAY)
1437 #define MPN_COPY_DECR(dst, src, n)                                      \
1438   do {                                                                  \
1439     int __i;            /* Faster on some Crays with plain int */       \
1440     _Pragma ("_CRI ivdep");                                             \
1441     for (__i = (n) - 1; __i >= 0; __i--)                                \
1442       (dst)[__i] = (src)[__i];                                          \
1443   } while (0)
1444 #endif
1445
1446 /* used by test programs, hence __GMP_DECLSPEC */
1447 #ifndef mpn_copyd  /* if not done with cpuvec in a fat binary */
1448 #define mpn_copyd __MPN(copyd)
1449 __GMP_DECLSPEC void mpn_copyd __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t));
1450 #endif
1451
1452 #if ! defined (MPN_COPY_DECR) && HAVE_NATIVE_mpn_copyd
1453 #define MPN_COPY_DECR(dst, src, size)                   \
1454   do {                                                  \
1455     ASSERT ((size) >= 0);                               \
1456     ASSERT (MPN_SAME_OR_DECR_P (dst, src, size));       \
1457     mpn_copyd (dst, src, size);                         \
1458   } while (0)
1459 #endif
1460
1461 /* Copy N limbs from SRC to DST decrementing, N==0 allowed.  */
1462 #if ! defined (MPN_COPY_DECR)
1463 #define MPN_COPY_DECR(dst, src, n)                      \
1464   do {                                                  \
1465     ASSERT ((n) >= 0);                                  \
1466     ASSERT (MPN_SAME_OR_DECR_P (dst, src, n));          \
1467     if ((n) != 0)                                       \
1468       {                                                 \
1469         mp_size_t __n = (n) - 1;                        \
1470         mp_ptr __dst = (dst) + __n;                     \
1471         mp_srcptr __src = (src) + __n;                  \
1472         mp_limb_t __x;                                  \
1473         __x = *__src--;                                 \
1474         if (__n != 0)                                   \
1475           {                                             \
1476             do                                          \
1477               {                                         \
1478                 *__dst-- = __x;                         \
1479                 __x = *__src--;                         \
1480               }                                         \
1481             while (--__n);                              \
1482           }                                             \
1483         *__dst-- = __x;                                 \
1484       }                                                 \
1485   } while (0)
1486 #endif
1487
1488
1489 #ifndef MPN_COPY
1490 #define MPN_COPY(d,s,n)                         \
1491   do {                                          \
1492     ASSERT (MPN_SAME_OR_SEPARATE_P (d, s, n));  \
1493     MPN_COPY_INCR (d, s, n);                    \
1494   } while (0)
1495 #endif
1496
1497
1498 /* Set {dst,size} to the limbs of {src,size} in reverse order. */
1499 #define MPN_REVERSE(dst, src, size)                     \
1500   do {                                                  \
1501     mp_ptr     __dst = (dst);                           \
1502     mp_size_t  __size = (size);                         \
1503     mp_srcptr  __src = (src) + __size - 1;              \
1504     mp_size_t  __i;                                     \
1505     ASSERT ((size) >= 0);                               \
1506     ASSERT (! MPN_OVERLAP_P (dst, size, src, size));    \
1507     CRAY_Pragma ("_CRI ivdep");                         \
1508     for (__i = 0; __i < __size; __i++)                  \
1509       {                                                 \
1510         *__dst = *__src;                                \
1511         __dst++;                                        \
1512         __src--;                                        \
1513       }                                                 \
1514   } while (0)
1515
1516
1517 /* Zero n limbs at dst.
1518
1519    For power and powerpc we want an inline stu/bdnz loop for zeroing.  On
1520    ppc630 for instance this is optimal since it can sustain only 1 store per
1521    cycle.
1522
1523    gcc 2.95.x (for powerpc64 -maix64, or powerpc32) doesn't recognise the
1524    "for" loop in the generic code below can become stu/bdnz.  The do/while
1525    here helps it get to that.  The same caveat about plain -mpowerpc64 mode
1526    applies here as to __GMPN_COPY_INCR in gmp.h.
1527
1528    xlc 3.1 already generates stu/bdnz from the generic C, and does so from
1529    this loop too.
1530
1531    Enhancement: GLIBC does some trickery with dcbz to zero whole cache lines
1532    at a time.  MPN_ZERO isn't all that important in GMP, so it might be more
1533    trouble than it's worth to do the same, though perhaps a call to memset
1534    would be good when on a GNU system.  */
1535
1536 #if HAVE_HOST_CPU_FAMILY_power || HAVE_HOST_CPU_FAMILY_powerpc
1537 #define MPN_ZERO(dst, n)                        \
1538   do {                                          \
1539     ASSERT ((n) >= 0);                          \
1540     if ((n) != 0)                               \
1541       {                                         \
1542         mp_ptr __dst = (dst) - 1;               \
1543         mp_size_t __n = (n);                    \
1544         do                                      \
1545           *++__dst = 0;                         \
1546         while (--__n);                          \
1547       }                                         \
1548   } while (0)
1549 #endif
1550
1551 #ifndef MPN_ZERO
1552 #define MPN_ZERO(dst, n)                        \
1553   do {                                          \
1554     ASSERT ((n) >= 0);                          \
1555     if ((n) != 0)                               \
1556       {                                         \
1557         mp_ptr __dst = (dst);                   \
1558         mp_size_t __n = (n);                    \
1559         do                                      \
1560           *__dst++ = 0;                         \
1561         while (--__n);                          \
1562       }                                         \
1563   } while (0)
1564 #endif
1565
1566
1567 /* On the x86s repe/scasl doesn't seem useful, since it takes many cycles to
1568    start up and would need to strip a lot of zeros before it'd be faster
1569    than a simple cmpl loop.  Here are some times in cycles for
1570    std/repe/scasl/cld and cld/repe/scasl (the latter would be for stripping
1571    low zeros).
1572
1573                 std   cld
1574            P5    18    16
1575            P6    46    38
1576            K6    36    13
1577            K7    21    20
1578 */
1579 #ifndef MPN_NORMALIZE
1580 #define MPN_NORMALIZE(DST, NLIMBS) \
1581   do {                                                                  \
1582     while ((NLIMBS) > 0)                                                \
1583       {                                                                 \
1584         if ((DST)[(NLIMBS) - 1] != 0)                                   \
1585           break;                                                        \
1586         (NLIMBS)--;                                                     \
1587       }                                                                 \
1588   } while (0)
1589 #endif
1590 #ifndef MPN_NORMALIZE_NOT_ZERO
1591 #define MPN_NORMALIZE_NOT_ZERO(DST, NLIMBS)     \
1592   do {                                          \
1593     ASSERT ((NLIMBS) >= 1);                     \
1594     while (1)                                   \
1595       {                                         \
1596         if ((DST)[(NLIMBS) - 1] != 0)           \
1597           break;                                \
1598         (NLIMBS)--;                             \
1599       }                                         \
1600   } while (0)
1601 #endif
1602
1603 /* Strip least significant zero limbs from {ptr,size} by incrementing ptr
1604    and decrementing size.  low should be ptr[0], and will be the new ptr[0]
1605    on returning.  The number in {ptr,size} must be non-zero, ie. size!=0 and
1606    somewhere a non-zero limb.  */
1607 #define MPN_STRIP_LOW_ZEROS_NOT_ZERO(ptr, size, low)    \
1608   do {                                                  \
1609     ASSERT ((size) >= 1);                               \
1610     ASSERT ((low) == (ptr)[0]);                         \
1611                                                         \
1612     while ((low) == 0)                                  \
1613       {                                                 \
1614         (size)--;                                       \
1615         ASSERT ((size) >= 1);                           \
1616         (ptr)++;                                        \
1617         (low) = *(ptr);                                 \
1618       }                                                 \
1619   } while (0)
1620
1621 /* Initialize X of type mpz_t with space for NLIMBS limbs.  X should be a
1622    temporary variable; it will be automatically cleared out at function
1623    return.  We use __x here to make it possible to accept both mpz_ptr and
1624    mpz_t arguments.  */
1625 #define MPZ_TMP_INIT(X, NLIMBS)                                         \
1626   do {                                                                  \
1627     mpz_ptr __x = (X);                                                  \
1628     ASSERT ((NLIMBS) >= 1);                                             \
1629     __x->_mp_alloc = (NLIMBS);                                          \
1630     __x->_mp_d = TMP_ALLOC_LIMBS (NLIMBS);                              \
1631   } while (0)
1632
1633 /* Realloc for an mpz_t WHAT if it has less than NEEDED limbs.  */
1634 #define MPZ_REALLOC(z,n) (UNLIKELY ((n) > ALLOC(z))     \
1635                           ? (mp_ptr) _mpz_realloc(z,n)  \
1636                           : PTR(z))
1637
1638 #define MPZ_EQUAL_1_P(z)  (SIZ(z)==1 && PTR(z)[0] == 1)
1639
1640
1641 /* MPN_FIB2_SIZE(n) is the size in limbs required by mpn_fib2_ui for fp and
1642    f1p.
1643
1644    From Knuth vol 1 section 1.2.8, F[n] = phi^n/sqrt(5) rounded to the
1645    nearest integer, where phi=(1+sqrt(5))/2 is the golden ratio.  So the
1646    number of bits required is n*log_2((1+sqrt(5))/2) = n*0.6942419.
1647
1648    The multiplier used is 23/32=0.71875 for efficient calculation on CPUs
1649    without good floating point.  There's +2 for rounding up, and a further
1650    +2 since at the last step x limbs are doubled into a 2x+1 limb region
1651    whereas the actual F[2k] value might be only 2x-1 limbs.
1652
1653    Note that a division is done first, since on a 32-bit system it's at
1654    least conceivable to go right up to n==ULONG_MAX.  (F[2^32-1] would be
1655    about 380Mbytes, plus temporary workspace of about 1.2Gbytes here and
1656    whatever a multiply of two 190Mbyte numbers takes.)
1657
1658    Enhancement: When GMP_NUMB_BITS is not a power of 2 the division could be
1659    worked into the multiplier.  */
1660
1661 #define MPN_FIB2_SIZE(n) \
1662   ((mp_size_t) ((n) / 32 * 23 / GMP_NUMB_BITS) + 4)
1663
1664
1665 /* FIB_TABLE(n) returns the Fibonacci number F[n].  Must have n in the range
1666    -1 <= n <= FIB_TABLE_LIMIT (that constant in fib_table.h).
1667
1668    FIB_TABLE_LUCNUM_LIMIT (in fib_table.h) is the largest n for which L[n] =
1669    F[n] + 2*F[n-1] fits in a limb.  */
1670
1671 __GMP_DECLSPEC extern const mp_limb_t __gmp_fib_table[];
1672 #define FIB_TABLE(n)  (__gmp_fib_table[(n)+1])
1673
1674 #define SIEVESIZE 512           /* FIXME: Allow gmp_init_primesieve to choose */
1675 typedef struct
1676 {
1677   unsigned long d;                 /* current index in s[] */
1678   unsigned long s0;                /* number corresponding to s[0] */
1679   unsigned long sqrt_s0;           /* misnomer for sqrt(s[SIEVESIZE-1]) */
1680   unsigned char s[SIEVESIZE + 1];  /* sieve table */
1681 } gmp_primesieve_t;
1682
1683 #define gmp_init_primesieve __gmp_init_primesieve
1684 __GMP_DECLSPEC void gmp_init_primesieve (gmp_primesieve_t *);
1685
1686 #define gmp_nextprime __gmp_nextprime
1687 __GMP_DECLSPEC unsigned long int gmp_nextprime (gmp_primesieve_t *);
1688
1689
1690 #ifndef MUL_TOOM22_THRESHOLD
1691 #define MUL_TOOM22_THRESHOLD             30
1692 #endif
1693
1694 #ifndef MUL_TOOM33_THRESHOLD
1695 #define MUL_TOOM33_THRESHOLD            100
1696 #endif
1697
1698 #ifndef MUL_TOOM44_THRESHOLD
1699 #define MUL_TOOM44_THRESHOLD            300
1700 #endif
1701
1702 #ifndef MUL_TOOM6H_THRESHOLD
1703 #define MUL_TOOM6H_THRESHOLD            350
1704 #endif
1705
1706 #ifndef SQR_TOOM6_THRESHOLD
1707 #define SQR_TOOM6_THRESHOLD MUL_TOOM6H_THRESHOLD
1708 #endif
1709
1710 #ifndef MUL_TOOM8H_THRESHOLD
1711 #define MUL_TOOM8H_THRESHOLD            450
1712 #endif
1713
1714 #ifndef SQR_TOOM8_THRESHOLD
1715 #define SQR_TOOM8_THRESHOLD MUL_TOOM8H_THRESHOLD
1716 #endif
1717
1718 #ifndef MUL_TOOM32_TO_TOOM43_THRESHOLD
1719 #define MUL_TOOM32_TO_TOOM43_THRESHOLD  100
1720 #endif
1721
1722 #ifndef MUL_TOOM32_TO_TOOM53_THRESHOLD
1723 #define MUL_TOOM32_TO_TOOM53_THRESHOLD  110
1724 #endif
1725
1726 #ifndef MUL_TOOM42_TO_TOOM53_THRESHOLD
1727 #define MUL_TOOM42_TO_TOOM53_THRESHOLD  100
1728 #endif
1729
1730 #ifndef MUL_TOOM42_TO_TOOM63_THRESHOLD
1731 #define MUL_TOOM42_TO_TOOM63_THRESHOLD  110
1732 #endif
1733
1734 /* MUL_TOOM22_THRESHOLD_LIMIT is the maximum for MUL_TOOM22_THRESHOLD.  In a
1735    normal build MUL_TOOM22_THRESHOLD is a constant and we use that.  In a fat
1736    binary or tune program build MUL_TOOM22_THRESHOLD is a variable and a
1737    separate hard limit will have been defined.  Similarly for TOOM3.  */
1738 #ifndef MUL_TOOM22_THRESHOLD_LIMIT
1739 #define MUL_TOOM22_THRESHOLD_LIMIT  MUL_TOOM22_THRESHOLD
1740 #endif
1741 #ifndef MUL_TOOM33_THRESHOLD_LIMIT
1742 #define MUL_TOOM33_THRESHOLD_LIMIT  MUL_TOOM33_THRESHOLD
1743 #endif
1744 #ifndef MULLO_BASECASE_THRESHOLD_LIMIT
1745 #define MULLO_BASECASE_THRESHOLD_LIMIT  MULLO_BASECASE_THRESHOLD
1746 #endif
1747
1748 /* SQR_BASECASE_THRESHOLD is where mpn_sqr_basecase should take over from
1749    mpn_mul_basecase.  Default is to use mpn_sqr_basecase from 0.  (Note that we
1750    certainly always want it if there's a native assembler mpn_sqr_basecase.)
1751
1752    If it turns out that mpn_toom2_sqr becomes faster than mpn_mul_basecase
1753    before mpn_sqr_basecase does, then SQR_BASECASE_THRESHOLD is the toom2
1754    threshold and SQR_TOOM2_THRESHOLD is 0.  This oddity arises more or less
1755    because SQR_TOOM2_THRESHOLD represents the size up to which mpn_sqr_basecase
1756    should be used, and that may be never.  */
1757
1758 #ifndef SQR_BASECASE_THRESHOLD
1759 #define SQR_BASECASE_THRESHOLD            0
1760 #endif
1761
1762 #ifndef SQR_TOOM2_THRESHOLD
1763 #define SQR_TOOM2_THRESHOLD              50
1764 #endif
1765
1766 #ifndef SQR_TOOM3_THRESHOLD
1767 #define SQR_TOOM3_THRESHOLD             120
1768 #endif
1769
1770 #ifndef SQR_TOOM4_THRESHOLD
1771 #define SQR_TOOM4_THRESHOLD             400
1772 #endif
1773
1774 /* See comments above about MUL_TOOM33_THRESHOLD_LIMIT.  */
1775 #ifndef SQR_TOOM3_THRESHOLD_LIMIT
1776 #define SQR_TOOM3_THRESHOLD_LIMIT  SQR_TOOM3_THRESHOLD
1777 #endif
1778
1779 #ifndef DC_DIV_QR_THRESHOLD
1780 #define DC_DIV_QR_THRESHOLD              50
1781 #endif
1782
1783 #ifndef DC_DIVAPPR_Q_THRESHOLD
1784 #define DC_DIVAPPR_Q_THRESHOLD          200
1785 #endif
1786
1787 #ifndef DC_BDIV_QR_THRESHOLD
1788 #define DC_BDIV_QR_THRESHOLD             50
1789 #endif
1790
1791 #ifndef DC_BDIV_Q_THRESHOLD
1792 #define DC_BDIV_Q_THRESHOLD             180
1793 #endif
1794
1795 #ifndef DIVEXACT_JEB_THRESHOLD
1796 #define DIVEXACT_JEB_THRESHOLD           25
1797 #endif
1798
1799 #ifndef INV_MULMOD_BNM1_THRESHOLD
1800 #define INV_MULMOD_BNM1_THRESHOLD  (5*MULMOD_BNM1_THRESHOLD)
1801 #endif
1802
1803 #ifndef INV_APPR_THRESHOLD
1804 #define INV_APPR_THRESHOLD         INV_NEWTON_THRESHOLD
1805 #endif
1806
1807 #ifndef INV_NEWTON_THRESHOLD
1808 #define INV_NEWTON_THRESHOLD            200
1809 #endif
1810
1811 #ifndef BINV_NEWTON_THRESHOLD
1812 #define BINV_NEWTON_THRESHOLD           300
1813 #endif
1814
1815 #ifndef MU_DIVAPPR_Q_THRESHOLD
1816 #define MU_DIVAPPR_Q_THRESHOLD         2000
1817 #endif
1818
1819 #ifndef MU_DIV_QR_THRESHOLD
1820 #define MU_DIV_QR_THRESHOLD            2000
1821 #endif
1822
1823 #ifndef MUPI_DIV_QR_THRESHOLD
1824 #define MUPI_DIV_QR_THRESHOLD           200
1825 #endif
1826
1827 #ifndef MU_BDIV_Q_THRESHOLD
1828 #define MU_BDIV_Q_THRESHOLD            2000
1829 #endif
1830
1831 #ifndef MU_BDIV_QR_THRESHOLD
1832 #define MU_BDIV_QR_THRESHOLD           2000
1833 #endif
1834
1835 #ifndef MULMOD_BNM1_THRESHOLD
1836 #define MULMOD_BNM1_THRESHOLD            16
1837 #endif
1838
1839 #ifndef SQRMOD_BNM1_THRESHOLD
1840 #define SQRMOD_BNM1_THRESHOLD            16
1841 #endif
1842
1843 #ifndef MUL_TO_MULMOD_BNM1_FOR_2NXN_THRESHOLD
1844 #define MUL_TO_MULMOD_BNM1_FOR_2NXN_THRESHOLD  (INV_MULMOD_BNM1_THRESHOLD/2)
1845 #endif
1846
1847 #if HAVE_NATIVE_mpn_addmul_2 || HAVE_NATIVE_mpn_redc_2
1848
1849 #ifndef REDC_1_TO_REDC_2_THRESHOLD
1850 #define REDC_1_TO_REDC_2_THRESHOLD       15
1851 #endif
1852 #ifndef REDC_2_TO_REDC_N_THRESHOLD
1853 #define REDC_2_TO_REDC_N_THRESHOLD      100
1854 #endif
1855
1856 #else
1857
1858 #ifndef REDC_1_TO_REDC_N_THRESHOLD
1859 #define REDC_1_TO_REDC_N_THRESHOLD      100
1860 #endif
1861
1862 #endif /* HAVE_NATIVE_mpn_addmul_2 || HAVE_NATIVE_mpn_redc_2 */
1863
1864
1865 /* First k to use for an FFT modF multiply.  A modF FFT is an order
1866    log(2^k)/log(2^(k-1)) algorithm, so k=3 is merely 1.5 like karatsuba,
1867    whereas k=4 is 1.33 which is faster than toom3 at 1.485.    */
1868 #define FFT_FIRST_K  4
1869
1870 /* Threshold at which FFT should be used to do a modF NxN -> N multiply. */
1871 #ifndef MUL_FFT_MODF_THRESHOLD
1872 #define MUL_FFT_MODF_THRESHOLD   (MUL_TOOM33_THRESHOLD * 3)
1873 #endif
1874 #ifndef SQR_FFT_MODF_THRESHOLD
1875 #define SQR_FFT_MODF_THRESHOLD   (SQR_TOOM3_THRESHOLD * 3)
1876 #endif
1877
1878 /* Threshold at which FFT should be used to do an NxN -> 2N multiply.  This
1879    will be a size where FFT is using k=7 or k=8, since an FFT-k used for an
1880    NxN->2N multiply and not recursing into itself is an order
1881    log(2^k)/log(2^(k-2)) algorithm, so it'll be at least k=7 at 1.39 which
1882    is the first better than toom3.  */
1883 #ifndef MUL_FFT_THRESHOLD
1884 #define MUL_FFT_THRESHOLD   (MUL_FFT_MODF_THRESHOLD * 10)
1885 #endif
1886 #ifndef SQR_FFT_THRESHOLD
1887 #define SQR_FFT_THRESHOLD   (SQR_FFT_MODF_THRESHOLD * 10)
1888 #endif
1889
1890 /* Table of thresholds for successive modF FFT "k"s.  The first entry is
1891    where FFT_FIRST_K+1 should be used, the second FFT_FIRST_K+2,
1892    etc.  See mpn_fft_best_k(). */
1893 #ifndef MUL_FFT_TABLE
1894 #define MUL_FFT_TABLE                           \
1895   { MUL_TOOM33_THRESHOLD * 4,   /* k=5 */        \
1896     MUL_TOOM33_THRESHOLD * 8,   /* k=6 */        \
1897     MUL_TOOM33_THRESHOLD * 16,  /* k=7 */        \
1898     MUL_TOOM33_THRESHOLD * 32,  /* k=8 */        \
1899     MUL_TOOM33_THRESHOLD * 96,  /* k=9 */        \
1900     MUL_TOOM33_THRESHOLD * 288, /* k=10 */       \
1901     0 }
1902 #endif
1903 #ifndef SQR_FFT_TABLE
1904 #define SQR_FFT_TABLE                           \
1905   { SQR_TOOM3_THRESHOLD * 4,   /* k=5 */        \
1906     SQR_TOOM3_THRESHOLD * 8,   /* k=6 */        \
1907     SQR_TOOM3_THRESHOLD * 16,  /* k=7 */        \
1908     SQR_TOOM3_THRESHOLD * 32,  /* k=8 */        \
1909     SQR_TOOM3_THRESHOLD * 96,  /* k=9 */        \
1910     SQR_TOOM3_THRESHOLD * 288, /* k=10 */       \
1911     0 }
1912 #endif
1913
1914 struct fft_table_nk
1915 {
1916   unsigned int n:27;
1917   unsigned int k:5;
1918 };
1919
1920 #ifndef FFT_TABLE_ATTRS
1921 #define FFT_TABLE_ATTRS   static const
1922 #endif
1923
1924 #define MPN_FFT_TABLE_SIZE  16
1925
1926
1927 #ifndef DC_DIV_QR_THRESHOLD
1928 #define DC_DIV_QR_THRESHOLD    (3 * MUL_TOOM22_THRESHOLD)
1929 #endif
1930
1931 #ifndef GET_STR_DC_THRESHOLD
1932 #define GET_STR_DC_THRESHOLD             18
1933 #endif
1934
1935 #ifndef GET_STR_PRECOMPUTE_THRESHOLD
1936 #define GET_STR_PRECOMPUTE_THRESHOLD     35
1937 #endif
1938
1939 #ifndef SET_STR_DC_THRESHOLD
1940 #define SET_STR_DC_THRESHOLD            750
1941 #endif
1942
1943 #ifndef SET_STR_PRECOMPUTE_THRESHOLD
1944 #define SET_STR_PRECOMPUTE_THRESHOLD   2000
1945 #endif
1946
1947 /* Return non-zero if xp,xsize and yp,ysize overlap.
1948    If xp+xsize<=yp there's no overlap, or if yp+ysize<=xp there's no
1949    overlap.  If both these are false, there's an overlap. */
1950 #define MPN_OVERLAP_P(xp, xsize, yp, ysize) \
1951   ((xp) + (xsize) > (yp) && (yp) + (ysize) > (xp))
1952 #define MEM_OVERLAP_P(xp, xsize, yp, ysize)     \
1953   (   (char *) (xp) + (xsize) > (char *) (yp)   \
1954    && (char *) (yp) + (ysize) > (char *) (xp))
1955
1956 /* Return non-zero if xp,xsize and yp,ysize are either identical or not
1957    overlapping.  Return zero if they're partially overlapping. */
1958 #define MPN_SAME_OR_SEPARATE_P(xp, yp, size)    \
1959   MPN_SAME_OR_SEPARATE2_P(xp, size, yp, size)
1960 #define MPN_SAME_OR_SEPARATE2_P(xp, xsize, yp, ysize)           \
1961   ((xp) == (yp) || ! MPN_OVERLAP_P (xp, xsize, yp, ysize))
1962
1963 /* Return non-zero if dst,dsize and src,ssize are either identical or
1964    overlapping in a way suitable for an incrementing/decrementing algorithm.
1965    Return zero if they're partially overlapping in an unsuitable fashion. */
1966 #define MPN_SAME_OR_INCR2_P(dst, dsize, src, ssize)             \
1967   ((dst) <= (src) || ! MPN_OVERLAP_P (dst, dsize, src, ssize))
1968 #define MPN_SAME_OR_INCR_P(dst, src, size)      \
1969   MPN_SAME_OR_INCR2_P(dst, size, src, size)
1970 #define MPN_SAME_OR_DECR2_P(dst, dsize, src, ssize)             \
1971   ((dst) >= (src) || ! MPN_OVERLAP_P (dst, dsize, src, ssize))
1972 #define MPN_SAME_OR_DECR_P(dst, src, size)      \
1973   MPN_SAME_OR_DECR2_P(dst, size, src, size)
1974
1975
1976 /* ASSERT() is a private assertion checking scheme, similar to <assert.h>.
1977    ASSERT() does the check only if WANT_ASSERT is selected, ASSERT_ALWAYS()
1978    does it always.  Generally assertions are meant for development, but
1979    might help when looking for a problem later too.
1980
1981    Note that strings shouldn't be used within the ASSERT expression,
1982    eg. ASSERT(strcmp(s,"notgood")!=0), since the quotes upset the "expr"
1983    used in the !HAVE_STRINGIZE case (ie. K&R).  */
1984
1985 #ifdef __LINE__
1986 #define ASSERT_LINE  __LINE__
1987 #else
1988 #define ASSERT_LINE  -1
1989 #endif
1990
1991 #ifdef __FILE__
1992 #define ASSERT_FILE  __FILE__
1993 #else
1994 #define ASSERT_FILE  ""
1995 #endif
1996
1997 __GMP_DECLSPEC void __gmp_assert_header __GMP_PROTO ((const char *, int));
1998 __GMP_DECLSPEC void __gmp_assert_fail __GMP_PROTO ((const char *, int, const char *)) ATTRIBUTE_NORETURN;
1999
2000 #if HAVE_STRINGIZE
2001 #define ASSERT_FAIL(expr)  __gmp_assert_fail (ASSERT_FILE, ASSERT_LINE, #expr)
2002 #else
2003 #define ASSERT_FAIL(expr)  __gmp_assert_fail (ASSERT_FILE, ASSERT_LINE, "expr")
2004 #endif
2005
2006 #define ASSERT_ALWAYS(expr)     \
2007   do {                          \
2008     if (!(expr))                \
2009       ASSERT_FAIL (expr);       \
2010   } while (0)
2011
2012 #if WANT_ASSERT
2013 #define ASSERT(expr)   ASSERT_ALWAYS (expr)
2014 #else
2015 #define ASSERT(expr)   do {} while (0)
2016 #endif
2017
2018
2019 /* ASSERT_CARRY checks the expression is non-zero, and ASSERT_NOCARRY checks
2020    that it's zero.  In both cases if assertion checking is disabled the
2021    expression is still evaluated.  These macros are meant for use with
2022    routines like mpn_add_n() where the return value represents a carry or
2023    whatever that should or shouldn't occur in some context.  For example,
2024    ASSERT_NOCARRY (mpn_add_n (rp, s1p, s2p, size)); */
2025 #if WANT_ASSERT
2026 #define ASSERT_CARRY(expr)     ASSERT_ALWAYS ((expr) != 0)
2027 #define ASSERT_NOCARRY(expr)   ASSERT_ALWAYS ((expr) == 0)
2028 #else
2029 #define ASSERT_CARRY(expr)     (expr)
2030 #define ASSERT_NOCARRY(expr)   (expr)
2031 #endif
2032
2033
2034 /* ASSERT_CODE includes code when assertion checking is wanted.  This is the
2035    same as writing "#if WANT_ASSERT", but more compact.  */
2036 #if WANT_ASSERT
2037 #define ASSERT_CODE(expr)  expr
2038 #else
2039 #define ASSERT_CODE(expr)
2040 #endif
2041
2042
2043 /* Test that an mpq_t is in fully canonical form.  This can be used as
2044    protection on routines like mpq_equal which give wrong results on
2045    non-canonical inputs.  */
2046 #if WANT_ASSERT
2047 #define ASSERT_MPQ_CANONICAL(q)                         \
2048   do {                                                  \
2049     ASSERT (q->_mp_den._mp_size > 0);                   \
2050     if (q->_mp_num._mp_size == 0)                       \
2051       {                                                 \
2052         /* zero should be 0/1 */                        \
2053         ASSERT (mpz_cmp_ui (mpq_denref(q), 1L) == 0);   \
2054       }                                                 \
2055     else                                                \
2056       {                                                 \
2057         /* no common factors */                         \
2058         mpz_t  __g;                                     \
2059         mpz_init (__g);                                 \
2060         mpz_gcd (__g, mpq_numref(q), mpq_denref(q));    \
2061         ASSERT (mpz_cmp_ui (__g, 1) == 0);              \
2062         mpz_clear (__g);                                \
2063       }                                                 \
2064   } while (0)
2065 #else
2066 #define ASSERT_MPQ_CANONICAL(q)  do {} while (0)
2067 #endif
2068
2069 /* Check that the nail parts are zero. */
2070 #define ASSERT_ALWAYS_LIMB(limb)                \
2071   do {                                          \
2072     mp_limb_t  __nail = (limb) & GMP_NAIL_MASK; \
2073     ASSERT_ALWAYS (__nail == 0);                \
2074   } while (0)
2075 #define ASSERT_ALWAYS_MPN(ptr, size)            \
2076   do {                                          \
2077     /* let whole loop go dead when no nails */  \
2078     if (GMP_NAIL_BITS != 0)                     \
2079       {                                         \
2080         mp_size_t  __i;                         \
2081         for (__i = 0; __i < (size); __i++)      \
2082           ASSERT_ALWAYS_LIMB ((ptr)[__i]);      \
2083       }                                         \
2084   } while (0)
2085 #if WANT_ASSERT
2086 #define ASSERT_LIMB(limb)       ASSERT_ALWAYS_LIMB (limb)
2087 #define ASSERT_MPN(ptr, size)   ASSERT_ALWAYS_MPN (ptr, size)
2088 #else
2089 #define ASSERT_LIMB(limb)       do {} while (0)
2090 #define ASSERT_MPN(ptr, size)   do {} while (0)
2091 #endif
2092
2093
2094 /* Assert that an mpn region {ptr,size} is zero, or non-zero.
2095    size==0 is allowed, and in that case {ptr,size} considered to be zero.  */
2096 #if WANT_ASSERT
2097 #define ASSERT_MPN_ZERO_P(ptr,size)     \
2098   do {                                  \
2099     mp_size_t  __i;                     \
2100     ASSERT ((size) >= 0);               \
2101     for (__i = 0; __i < (size); __i++)  \
2102       ASSERT ((ptr)[__i] == 0);         \
2103   } while (0)
2104 #define ASSERT_MPN_NONZERO_P(ptr,size)  \
2105   do {                                  \
2106     mp_size_t  __i;                     \
2107     int        __nonzero = 0;           \
2108     ASSERT ((size) >= 0);               \
2109     for (__i = 0; __i < (size); __i++)  \
2110       if ((ptr)[__i] != 0)              \
2111         {                               \
2112           __nonzero = 1;                \
2113           break;                        \
2114         }                               \
2115     ASSERT (__nonzero);                 \
2116   } while (0)
2117 #else
2118 #define ASSERT_MPN_ZERO_P(ptr,size)     do {} while (0)
2119 #define ASSERT_MPN_NONZERO_P(ptr,size)  do {} while (0)
2120 #endif
2121
2122
2123 #if ! HAVE_NATIVE_mpn_com
2124 #undef mpn_com
2125 #define mpn_com(d,s,n)                                  \
2126   do {                                                  \
2127     mp_ptr     __d = (d);                               \
2128     mp_srcptr  __s = (s);                               \
2129     mp_size_t  __n = (n);                               \
2130     ASSERT (__n >= 1);                                  \
2131     ASSERT (MPN_SAME_OR_SEPARATE_P (__d, __s, __n));    \
2132     do                                                  \
2133       *__d++ = (~ *__s++) & GMP_NUMB_MASK;              \
2134     while (--__n);                                      \
2135   } while (0)
2136 #endif
2137
2138 #define MPN_LOGOPS_N_INLINE(rp, up, vp, n, operation)                   \
2139   do {                                                                  \
2140     mp_srcptr   __up = (up);                                            \
2141     mp_srcptr   __vp = (vp);                                            \
2142     mp_ptr      __rp = (rp);                                            \
2143     mp_size_t   __n = (n);                                              \
2144     mp_limb_t __a, __b;                                                 \
2145     ASSERT (__n > 0);                                                   \
2146     ASSERT (MPN_SAME_OR_SEPARATE_P (__rp, __up, __n));                  \
2147     ASSERT (MPN_SAME_OR_SEPARATE_P (__rp, __vp, __n));                  \
2148     __up += __n;                                                        \
2149     __vp += __n;                                                        \
2150     __rp += __n;                                                        \
2151     __n = -__n;                                                         \
2152     do {                                                                \
2153       __a = __up[__n];                                                  \
2154       __b = __vp[__n];                                                  \
2155       __rp[__n] = operation;                                            \
2156     } while (++__n);                                                    \
2157   } while (0)
2158
2159
2160 #if ! HAVE_NATIVE_mpn_and_n
2161 #undef mpn_and_n
2162 #define mpn_and_n(rp, up, vp, n) \
2163   MPN_LOGOPS_N_INLINE (rp, up, vp, n, __a & __b)
2164 #endif
2165
2166 #if ! HAVE_NATIVE_mpn_andn_n
2167 #undef mpn_andn_n
2168 #define mpn_andn_n(rp, up, vp, n) \
2169   MPN_LOGOPS_N_INLINE (rp, up, vp, n, __a & ~__b)
2170 #endif
2171
2172 #if ! HAVE_NATIVE_mpn_nand_n
2173 #undef mpn_nand_n
2174 #define mpn_nand_n(rp, up, vp, n) \
2175   MPN_LOGOPS_N_INLINE (rp, up, vp, n, ~(__a & __b) & GMP_NUMB_MASK)
2176 #endif
2177
2178 #if ! HAVE_NATIVE_mpn_ior_n
2179 #undef mpn_ior_n
2180 #define mpn_ior_n(rp, up, vp, n) \
2181   MPN_LOGOPS_N_INLINE (rp, up, vp, n, __a | __b)
2182 #endif
2183
2184 #if ! HAVE_NATIVE_mpn_iorn_n
2185 #undef mpn_iorn_n
2186 #define mpn_iorn_n(rp, up, vp, n) \
2187   MPN_LOGOPS_N_INLINE (rp, up, vp, n, (__a | ~__b) & GMP_NUMB_MASK)
2188 #endif
2189
2190 #if ! HAVE_NATIVE_mpn_nior_n
2191 #undef mpn_nior_n
2192 #define mpn_nior_n(rp, up, vp, n) \
2193   MPN_LOGOPS_N_INLINE (rp, up, vp, n, ~(__a | __b) & GMP_NUMB_MASK)
2194 #endif
2195
2196 #if ! HAVE_NATIVE_mpn_xor_n
2197 #undef mpn_xor_n
2198 #define mpn_xor_n(rp, up, vp, n) \
2199   MPN_LOGOPS_N_INLINE (rp, up, vp, n, __a ^ __b)
2200 #endif
2201
2202 #if ! HAVE_NATIVE_mpn_xnor_n
2203 #undef mpn_xnor_n
2204 #define mpn_xnor_n(rp, up, vp, n) \
2205   MPN_LOGOPS_N_INLINE (rp, up, vp, n, ~(__a ^ __b) & GMP_NUMB_MASK)
2206 #endif
2207
2208 #define mpn_trialdiv __MPN(trialdiv)
2209 __GMP_DECLSPEC mp_limb_t mpn_trialdiv __GMP_PROTO ((mp_srcptr, mp_size_t, mp_size_t, int *));
2210
2211 #define mpn_remove __MPN(remove)
2212 __GMP_DECLSPEC mp_bitcnt_t mpn_remove __GMP_PROTO ((mp_ptr, mp_size_t *, mp_ptr, mp_size_t, mp_ptr, mp_size_t, mp_bitcnt_t));
2213
2214
2215 /* ADDC_LIMB sets w=x+y and cout to 0 or 1 for a carry from that addition. */
2216 #if GMP_NAIL_BITS == 0
2217 #define ADDC_LIMB(cout, w, x, y)        \
2218   do {                                  \
2219     mp_limb_t  __x = (x);               \
2220     mp_limb_t  __y = (y);               \
2221     mp_limb_t  __w = __x + __y;         \
2222     (w) = __w;                          \
2223     (cout) = __w < __x;                 \
2224   } while (0)
2225 #else
2226 #define ADDC_LIMB(cout, w, x, y)        \
2227   do {                                  \
2228     mp_limb_t  __w;                     \
2229     ASSERT_LIMB (x);                    \
2230     ASSERT_LIMB (y);                    \
2231     __w = (x) + (y);                    \
2232     (w) = __w & GMP_NUMB_MASK;          \
2233     (cout) = __w >> GMP_NUMB_BITS;      \
2234   } while (0)
2235 #endif
2236
2237 /* SUBC_LIMB sets w=x-y and cout to 0 or 1 for a borrow from that
2238    subtract.  */
2239 #if GMP_NAIL_BITS == 0
2240 #define SUBC_LIMB(cout, w, x, y)        \
2241   do {                                  \
2242     mp_limb_t  __x = (x);               \
2243     mp_limb_t  __y = (y);               \
2244     mp_limb_t  __w = __x - __y;         \
2245     (w) = __w;                          \
2246     (cout) = __w > __x;                 \
2247   } while (0)
2248 #else
2249 #define SUBC_LIMB(cout, w, x, y)        \
2250   do {                                  \
2251     mp_limb_t  __w = (x) - (y);         \
2252     (w) = __w & GMP_NUMB_MASK;          \
2253     (cout) = __w >> (GMP_LIMB_BITS-1);  \
2254   } while (0)
2255 #endif
2256
2257
2258 /* MPN_INCR_U does {ptr,size} += n, MPN_DECR_U does {ptr,size} -= n, both
2259    expecting no carry (or borrow) from that.
2260
2261    The size parameter is only for the benefit of assertion checking.  In a
2262    normal build it's unused and the carry/borrow is just propagated as far
2263    as it needs to go.
2264
2265    On random data, usually only one or two limbs of {ptr,size} get updated,
2266    so there's no need for any sophisticated looping, just something compact
2267    and sensible.
2268
2269    FIXME: Switch all code from mpn_{incr,decr}_u to MPN_{INCR,DECR}_U,
2270    declaring their operand sizes, then remove the former.  This is purely
2271    for the benefit of assertion checking.  */
2272
2273 #if defined (__GNUC__) && HAVE_HOST_CPU_FAMILY_x86 && GMP_NAIL_BITS == 0      \
2274   && GMP_LIMB_BITS == 32 && ! defined (NO_ASM) && ! WANT_ASSERT
2275 /* Better flags handling than the generic C gives on i386, saving a few
2276    bytes of code and maybe a cycle or two.  */
2277
2278 #define MPN_IORD_U(ptr, incr, aors)                                     \
2279   do {                                                                  \
2280     mp_ptr  __ptr_dummy;                                                \
2281     if (__builtin_constant_p (incr) && (incr) == 1)                     \
2282       {                                                                 \
2283         __asm__ __volatile__                                            \
2284           ("\n" ASM_L(top) ":\n"                                        \
2285            "\t" aors " $1, (%0)\n"                                      \
2286            "\tleal 4(%0),%0\n"                                          \
2287            "\tjc " ASM_L(top)                                           \
2288            : "=r" (__ptr_dummy)                                         \
2289            : "0"  (ptr)                                                 \
2290            : "memory");                                                 \
2291       }                                                                 \
2292     else                                                                \
2293       {                                                                 \
2294         __asm__ __volatile__                                            \
2295           (   aors  " %2,(%0)\n"                                        \
2296            "\tjnc " ASM_L(done) "\n"                                    \
2297            ASM_L(top) ":\n"                                             \
2298            "\t" aors " $1,4(%0)\n"                                      \
2299            "\tleal 4(%0),%0\n"                                          \
2300            "\tjc " ASM_L(top) "\n"                                      \
2301            ASM_L(done) ":\n"                                            \
2302            : "=r" (__ptr_dummy)                                         \
2303            : "0"  (ptr),                                                \
2304              "ri" (incr)                                                \
2305            : "memory");                                                 \
2306       }                                                                 \
2307   } while (0)
2308
2309 #define MPN_INCR_U(ptr, size, incr)  MPN_IORD_U (ptr, incr, "addl")
2310 #define MPN_DECR_U(ptr, size, incr)  MPN_IORD_U (ptr, incr, "subl")
2311 #define mpn_incr_u(ptr, incr)  MPN_INCR_U (ptr, 0, incr)
2312 #define mpn_decr_u(ptr, incr)  MPN_DECR_U (ptr, 0, incr)
2313 #endif
2314
2315 #if GMP_NAIL_BITS == 0
2316 #ifndef mpn_incr_u
2317 #define mpn_incr_u(p,incr)                              \
2318   do {                                                  \
2319     mp_limb_t __x;                                      \
2320     mp_ptr __p = (p);                                   \
2321     if (__builtin_constant_p (incr) && (incr) == 1)     \
2322       {                                                 \
2323         while (++(*(__p++)) == 0)                       \
2324           ;                                             \
2325       }                                                 \
2326     else                                                \
2327       {                                                 \
2328         __x = *__p + (incr);                            \
2329         *__p = __x;                                     \
2330         if (__x < (incr))                               \
2331           while (++(*(++__p)) == 0)                     \
2332             ;                                           \
2333       }                                                 \
2334   } while (0)
2335 #endif
2336 #ifndef mpn_decr_u
2337 #define mpn_decr_u(p,incr)                              \
2338   do {                                                  \
2339     mp_limb_t __x;                                      \
2340     mp_ptr __p = (p);                                   \
2341     if (__builtin_constant_p (incr) && (incr) == 1)     \
2342       {                                                 \
2343         while ((*(__p++))-- == 0)                       \
2344           ;                                             \
2345       }                                                 \
2346     else                                                \
2347       {                                                 \
2348         __x = *__p;                                     \
2349         *__p = __x - (incr);                            \
2350         if (__x < (incr))                               \
2351           while ((*(++__p))-- == 0)                     \
2352             ;                                           \
2353       }                                                 \
2354   } while (0)
2355 #endif
2356 #endif
2357
2358 #if GMP_NAIL_BITS >= 1
2359 #ifndef mpn_incr_u
2360 #define mpn_incr_u(p,incr)                              \
2361   do {                                                  \
2362     mp_limb_t __x;                                      \
2363     mp_ptr __p = (p);                                   \
2364     if (__builtin_constant_p (incr) && (incr) == 1)     \
2365       {                                                 \
2366         do                                              \
2367           {                                             \
2368             __x = (*__p + 1) & GMP_NUMB_MASK;           \
2369             *__p++ = __x;                               \
2370           }                                             \
2371         while (__x == 0);                               \
2372       }                                                 \
2373     else                                                \
2374       {                                                 \
2375         __x = (*__p + (incr));                          \
2376         *__p++ = __x & GMP_NUMB_MASK;                   \
2377         if (__x >> GMP_NUMB_BITS != 0)                  \
2378           {                                             \
2379             do                                          \
2380               {                                         \
2381                 __x = (*__p + 1) & GMP_NUMB_MASK;       \
2382                 *__p++ = __x;                           \
2383               }                                         \
2384             while (__x == 0);                           \
2385           }                                             \
2386       }                                                 \
2387   } while (0)
2388 #endif
2389 #ifndef mpn_decr_u
2390 #define mpn_decr_u(p,incr)                              \
2391   do {                                                  \
2392     mp_limb_t __x;                                      \
2393     mp_ptr __p = (p);                                   \
2394     if (__builtin_constant_p (incr) && (incr) == 1)     \
2395       {                                                 \
2396         do                                              \
2397           {                                             \
2398             __x = *__p;                                 \
2399             *__p++ = (__x - 1) & GMP_NUMB_MASK;         \
2400           }                                             \
2401         while (__x == 0);                               \
2402       }                                                 \
2403     else                                                \
2404       {                                                 \
2405         __x = *__p - (incr);                            \
2406         *__p++ = __x & GMP_NUMB_MASK;                   \
2407         if (__x >> GMP_NUMB_BITS != 0)                  \
2408           {                                             \
2409             do                                          \
2410               {                                         \
2411                 __x = *__p;                             \
2412                 *__p++ = (__x - 1) & GMP_NUMB_MASK;     \
2413               }                                         \
2414             while (__x == 0);                           \
2415           }                                             \
2416       }                                                 \
2417   } while (0)
2418 #endif
2419 #endif
2420
2421 #ifndef MPN_INCR_U
2422 #if WANT_ASSERT
2423 #define MPN_INCR_U(ptr, size, n)                        \
2424   do {                                                  \
2425     ASSERT ((size) >= 1);                               \
2426     ASSERT_NOCARRY (mpn_add_1 (ptr, ptr, size, n));     \
2427   } while (0)
2428 #else
2429 #define MPN_INCR_U(ptr, size, n)   mpn_incr_u (ptr, n)
2430 #endif
2431 #endif
2432
2433 #ifndef MPN_DECR_U
2434 #if WANT_ASSERT
2435 #define MPN_DECR_U(ptr, size, n)                        \
2436   do {                                                  \
2437     ASSERT ((size) >= 1);                               \
2438     ASSERT_NOCARRY (mpn_sub_1 (ptr, ptr, size, n));     \
2439   } while (0)
2440 #else
2441 #define MPN_DECR_U(ptr, size, n)   mpn_decr_u (ptr, n)
2442 #endif
2443 #endif
2444
2445
2446 /* Structure for conversion between internal binary format and
2447    strings in base 2..36.  */
2448 struct bases
2449 {
2450   /* Number of digits in the conversion base that always fits in an mp_limb_t.
2451      For example, for base 10 on a machine where a mp_limb_t has 32 bits this
2452      is 9, since 10**9 is the largest number that fits into a mp_limb_t.  */
2453   int chars_per_limb;
2454
2455   /* log(2)/log(conversion_base) */
2456   double chars_per_bit_exactly;
2457
2458   /* base**chars_per_limb, i.e. the biggest number that fits a word, built by
2459      factors of base.  Exception: For 2, 4, 8, etc, big_base is log2(base),
2460      i.e. the number of bits used to represent each digit in the base.  */
2461   mp_limb_t big_base;
2462
2463   /* A GMP_LIMB_BITS bit approximation to 1/big_base, represented as a
2464      fixed-point number.  Instead of dividing by big_base an application can
2465      choose to multiply by big_base_inverted.  */
2466   mp_limb_t big_base_inverted;
2467 };
2468
2469 #define   mp_bases __MPN(bases)
2470 __GMP_DECLSPEC extern const struct bases mp_bases[257];
2471
2472
2473 /* For power of 2 bases this is exact.  For other bases the result is either
2474    exact or one too big.
2475
2476    To be exact always it'd be necessary to examine all the limbs of the
2477    operand, since numbers like 100..000 and 99...999 generally differ only
2478    in the lowest limb.  It'd be possible to examine just a couple of high
2479    limbs to increase the probability of being exact, but that doesn't seem
2480    worth bothering with.  */
2481
2482 #define MPN_SIZEINBASE(result, ptr, size, base)                         \
2483   do {                                                                  \
2484     int       __lb_base, __cnt;                                         \
2485     size_t __totbits;                                                   \
2486                                                                         \
2487     ASSERT ((size) >= 0);                                               \
2488     ASSERT ((base) >= 2);                                               \
2489     ASSERT ((base) < numberof (mp_bases));                              \
2490                                                                         \
2491     /* Special case for X == 0.  */                                     \
2492     if ((size) == 0)                                                    \
2493       (result) = 1;                                                     \
2494     else                                                                \
2495       {                                                                 \
2496         /* Calculate the total number of significant bits of X.  */     \
2497         count_leading_zeros (__cnt, (ptr)[(size)-1]);                   \
2498         __totbits = (size_t) (size) * GMP_NUMB_BITS - (__cnt - GMP_NAIL_BITS);\
2499                                                                         \
2500         if (POW2_P (base))                                              \
2501           {                                                             \
2502             __lb_base = mp_bases[base].big_base;                        \
2503             (result) = (__totbits + __lb_base - 1) / __lb_base;         \
2504           }                                                             \
2505         else                                                            \
2506           (result) = (size_t)                                           \
2507             (__totbits * mp_bases[base].chars_per_bit_exactly) + 1;     \
2508       }                                                                 \
2509   } while (0)
2510
2511 /* eliminate mp_bases lookups for base==16 */
2512 #define MPN_SIZEINBASE_16(result, ptr, size)                            \
2513   do {                                                                  \
2514     int       __cnt;                                                    \
2515     mp_size_t __totbits;                                                \
2516                                                                         \
2517     ASSERT ((size) >= 0);                                               \
2518                                                                         \
2519     /* Special case for X == 0.  */                                     \
2520     if ((size) == 0)                                                    \
2521       (result) = 1;                                                     \
2522     else                                                                \
2523       {                                                                 \
2524         /* Calculate the total number of significant bits of X.  */     \
2525         count_leading_zeros (__cnt, (ptr)[(size)-1]);                   \
2526         __totbits = (size_t) (size) * GMP_NUMB_BITS - (__cnt - GMP_NAIL_BITS);\
2527         (result) = (__totbits + 4 - 1) / 4;                             \
2528       }                                                                 \
2529   } while (0)
2530
2531 /* bit count to limb count, rounding up */
2532 #define BITS_TO_LIMBS(n)  (((n) + (GMP_NUMB_BITS - 1)) / GMP_NUMB_BITS)
2533
2534 /* MPN_SET_UI sets an mpn (ptr, cnt) to given ui.  MPZ_FAKE_UI creates fake
2535    mpz_t from ui.  The zp argument must have room for LIMBS_PER_ULONG limbs
2536    in both cases (LIMBS_PER_ULONG is also defined here.) */
2537 #if BITS_PER_ULONG <= GMP_NUMB_BITS /* need one limb per ulong */
2538
2539 #define LIMBS_PER_ULONG 1
2540 #define MPN_SET_UI(zp, zn, u)   \
2541   (zp)[0] = (u);                \
2542   (zn) = ((zp)[0] != 0);
2543 #define MPZ_FAKE_UI(z, zp, u)   \
2544   (zp)[0] = (u);                \
2545   PTR (z) = (zp);               \
2546   SIZ (z) = ((zp)[0] != 0);     \
2547   ASSERT_CODE (ALLOC (z) = 1);
2548
2549 #else /* need two limbs per ulong */
2550
2551 #define LIMBS_PER_ULONG 2
2552 #define MPN_SET_UI(zp, zn, u)                          \
2553   (zp)[0] = (u) & GMP_NUMB_MASK;                       \
2554   (zp)[1] = (u) >> GMP_NUMB_BITS;                      \
2555   (zn) = ((zp)[1] != 0 ? 2 : (zp)[0] != 0 ? 1 : 0);
2556 #define MPZ_FAKE_UI(z, zp, u)                          \
2557   (zp)[0] = (u) & GMP_NUMB_MASK;                       \
2558   (zp)[1] = (u) >> GMP_NUMB_BITS;                      \
2559   SIZ (z) = ((zp)[1] != 0 ? 2 : (zp)[0] != 0 ? 1 : 0); \
2560   PTR (z) = (zp);                                      \
2561   ASSERT_CODE (ALLOC (z) = 2);
2562
2563 #endif
2564
2565
2566 #if HAVE_HOST_CPU_FAMILY_x86
2567 #define TARGET_REGISTER_STARVED 1
2568 #else
2569 #define TARGET_REGISTER_STARVED 0
2570 #endif
2571
2572
2573 /* LIMB_HIGHBIT_TO_MASK(n) examines the high bit of a limb value and turns 1
2574    or 0 there into a limb 0xFF..FF or 0 respectively.
2575
2576    On most CPUs this is just an arithmetic right shift by GMP_LIMB_BITS-1,
2577    but C99 doesn't guarantee signed right shifts are arithmetic, so we have
2578    a little compile-time test and a fallback to a "? :" form.  The latter is
2579    necessary for instance on Cray vector systems.
2580
2581    Recent versions of gcc (eg. 3.3) will in fact optimize a "? :" like this
2582    to an arithmetic right shift anyway, but it's good to get the desired
2583    shift on past versions too (in particular since an important use of
2584    LIMB_HIGHBIT_TO_MASK is in udiv_qrnnd_preinv).  */
2585
2586 #define LIMB_HIGHBIT_TO_MASK(n)                                 \
2587   (((mp_limb_signed_t) -1 >> 1) < 0                             \
2588    ? (mp_limb_signed_t) (n) >> (GMP_LIMB_BITS - 1)              \
2589    : (n) & GMP_LIMB_HIGHBIT ? MP_LIMB_T_MAX : CNST_LIMB(0))
2590
2591
2592 /* Use a library function for invert_limb, if available. */
2593 #define   mpn_invert_limb __MPN(invert_limb)
2594 __GMP_DECLSPEC mp_limb_t mpn_invert_limb __GMP_PROTO ((mp_limb_t)) ATTRIBUTE_CONST;
2595 #if ! defined (invert_limb) && HAVE_NATIVE_mpn_invert_limb
2596 #define invert_limb(invxl,xl)           \
2597   do {                                  \
2598     (invxl) = mpn_invert_limb (xl);     \
2599   } while (0)
2600 #endif
2601
2602 #ifndef invert_limb
2603 #define invert_limb(invxl,xl)                   \
2604   do {                                          \
2605     mp_limb_t dummy;                            \
2606     ASSERT ((xl) != 0);                         \
2607     udiv_qrnnd (invxl, dummy, ~(xl), ~CNST_LIMB(0), xl);  \
2608   } while (0)
2609 #endif
2610
2611 #define invert_pi1(dinv, d1, d0)                                \
2612   do {                                                          \
2613     mp_limb_t v, p, t1, t0, mask;                               \
2614     invert_limb (v, d1);                                        \
2615     p = d1 * v;                                                 \
2616     p += d0;                                                    \
2617     if (p < d0)                                                 \
2618       {                                                         \
2619         v--;                                                    \
2620         mask = -(p >= d1);                                      \
2621         p -= d1;                                                \
2622         v += mask;                                              \
2623         p -= mask & d1;                                         \
2624       }                                                         \
2625     umul_ppmm (t1, t0, d0, v);                                  \
2626     p += t1;                                                    \
2627     if (p < t1)                                                 \
2628       {                                                         \
2629         v--;                                                    \
2630         if (UNLIKELY (p >= d1))                                 \
2631           {                                                     \
2632             if (p > d1 || t0 >= d0)                             \
2633               v--;                                              \
2634           }                                                     \
2635       }                                                         \
2636     (dinv).inv32 = v;                                           \
2637   } while (0)
2638
2639
2640 #ifndef udiv_qrnnd_preinv
2641 #define udiv_qrnnd_preinv udiv_qrnnd_preinv3
2642 #endif
2643
2644 /* Divide the two-limb number in (NH,,NL) by D, with DI being the largest
2645    limb not larger than (2**(2*GMP_LIMB_BITS))/D - (2**GMP_LIMB_BITS).
2646    If this would yield overflow, DI should be the largest possible number
2647    (i.e., only ones).  For correct operation, the most significant bit of D
2648    has to be set.  Put the quotient in Q and the remainder in R.  */
2649 #define udiv_qrnnd_preinv1(q, r, nh, nl, d, di)                         \
2650   do {                                                                  \
2651     mp_limb_t _q, _ql, _r;                                              \
2652     mp_limb_t _xh, _xl;                                                 \
2653     ASSERT ((d) != 0);                                                  \
2654     umul_ppmm (_q, _ql, (nh), (di));                                    \
2655     _q += (nh); /* Compensate, di is 2**GMP_LIMB_BITS too small */      \
2656     umul_ppmm (_xh, _xl, _q, (d));                                      \
2657     sub_ddmmss (_xh, _r, (nh), (nl), _xh, _xl);                         \
2658     if (_xh != 0)                                                       \
2659       {                                                                 \
2660         sub_ddmmss (_xh, _r, _xh, _r, 0, (d));                          \
2661         _q += 1;                                                        \
2662         if (_xh != 0)                                                   \
2663           {                                                             \
2664             _r -= (d);                                                  \
2665             _q += 1;                                                    \
2666           }                                                             \
2667       }                                                                 \
2668     if (_r >= (d))                                                      \
2669       {                                                                 \
2670         _r -= (d);                                                      \
2671         _q += 1;                                                        \
2672       }                                                                 \
2673     (r) = _r;                                                           \
2674     (q) = _q;                                                           \
2675   } while (0)
2676
2677 /* Like udiv_qrnnd_preinv, but branch-free. */
2678 #define udiv_qrnnd_preinv2(q, r, nh, nl, d, di)                         \
2679   do {                                                                  \
2680     mp_limb_t _n2, _n10, _nmask, _nadj, _q1;                            \
2681     mp_limb_t _xh, _xl;                                                 \
2682     _n2 = (nh);                                                         \
2683     _n10 = (nl);                                                        \
2684     _nmask = LIMB_HIGHBIT_TO_MASK (_n10);                               \
2685     _nadj = _n10 + (_nmask & (d));                                      \
2686     umul_ppmm (_xh, _xl, di, _n2 - _nmask);                             \
2687     add_ssaaaa (_xh, _xl, _xh, _xl, _n2, _nadj);                        \
2688     _q1 = ~_xh;                                                         \
2689     umul_ppmm (_xh, _xl, _q1, d);                                       \
2690     add_ssaaaa (_xh, _xl, _xh, _xl, nh, nl);                            \
2691     _xh -= (d);                                 /* xh = 0 or -1 */      \
2692     (r) = _xl + ((d) & _xh);                                            \
2693     (q) = _xh - _q1;                                                    \
2694   } while (0)
2695
2696 /* Like udiv_qrnnd_preinv2, but for for any value D.  DNORM is D shifted left
2697    so that its most significant bit is set.  LGUP is ceil(log2(D)).  */
2698 #define udiv_qrnnd_preinv2gen(q, r, nh, nl, d, di, dnorm, lgup) \
2699   do {                                                                  \
2700     mp_limb_t _n2, _n10, _nmask, _nadj, _q1;                            \
2701     mp_limb_t _xh, _xl;                                                 \
2702     _n2 = ((nh) << (GMP_LIMB_BITS - (lgup))) + ((nl) >> 1 >> (l - 1));  \
2703     _n10 = (nl) << (GMP_LIMB_BITS - (lgup));                            \
2704     _nmask = LIMB_HIGHBIT_TO_MASK (_n10);                               \
2705     _nadj = _n10 + (_nmask & (dnorm));                                  \
2706     umul_ppmm (_xh, _xl, di, _n2 - _nmask);                             \
2707     add_ssaaaa (_xh, _xl, _xh, _xl, _n2, _nadj);                        \
2708     _q1 = ~_xh;                                                         \
2709     umul_ppmm (_xh, _xl, _q1, d);                                       \
2710     add_ssaaaa (_xh, _xl, _xh, _xl, nh, nl);                            \
2711     _xh -= (d);                                                         \
2712     (r) = _xl + ((d) & _xh);                                            \
2713     (q) = _xh - _q1;                                                    \
2714   } while (0)
2715
2716 /* udiv_qrnnd_preinv3 -- Based on work by Niels Möller and Torbjörn Granlund.
2717
2718    We write things strangely below, to help gcc.  A more straightforward
2719    version:
2720
2721    _r = (nl) - _qh * (d);
2722    _t = _r + (d);
2723    if (_r >= _ql)
2724      {
2725        _qh--;
2726        _r = _t;
2727      }
2728
2729    For one operation shorter critical path, one may want to use this form:
2730
2731    _p = _qh * (d)
2732    _s = (nl) + (d);
2733    _r = (nl) - _p;
2734    _t = _s - _p;
2735    if (_r >= _ql)
2736      {
2737        _qh--;
2738        _r = _t;
2739      }
2740 */
2741 #define udiv_qrnnd_preinv3(q, r, nh, nl, d, di)                         \
2742   do {                                                                  \
2743     mp_limb_t _qh, _ql, _r;                                             \
2744     umul_ppmm (_qh, _ql, (nh), (di));                                   \
2745     if (__builtin_constant_p (nl) && (nl) == 0)                         \
2746       _qh += (nh) + 1;                                                  \
2747     else                                                                \
2748       add_ssaaaa (_qh, _ql, _qh, _ql, (nh) + 1, (nl));                  \
2749     _r = (nl) - _qh * (d);                                              \
2750     if (_r > _ql)       /* both > and >= should be OK */                \
2751       {                                                                 \
2752         _r += (d);                                                      \
2753         _qh--;                                                          \
2754       }                                                                 \
2755     if (UNLIKELY (_r >= (d)))                                           \
2756       {                                                                 \
2757         _r -= (d);                                                      \
2758         _qh++;                                                          \
2759       }                                                                 \
2760     (r) = _r;                                                           \
2761     (q) = _qh;                                                          \
2762   } while (0)
2763
2764 /* Compute r = nh*B mod d, where di is the inverse of d.  */
2765 #define udiv_rnd_preinv(r, nh, d, di)                                   \
2766   do {                                                                  \
2767     mp_limb_t _qh, _ql, _r;                                             \
2768     umul_ppmm (_qh, _ql, (nh), (di));                                   \
2769     _qh += (nh) + 1;                                                    \
2770     _r = - _qh * (d);                                                   \
2771     if (_r > _ql)                                                       \
2772       _r += (d);                                                        \
2773     (r) = _r;                                                           \
2774   } while (0)
2775
2776 /* Compute quotient the quotient and remainder for n / d. Requires d
2777    >= B^2 / 2 and n < d B. di is the inverse
2778
2779      floor ((B^3 - 1) / (d0 + d1 B)) - B.
2780
2781    NOTE: Output variables are updated multiple times. Only some inputs
2782    and outputs may overlap.
2783 */
2784 #define udiv_qr_3by2(q, r1, r0, n2, n1, n0, d1, d0, dinv)               \
2785   do {                                                                  \
2786     mp_limb_t _q0, _t1, _t0, _mask;                                     \
2787     umul_ppmm ((q), _q0, (n2), (dinv));                                 \
2788     add_ssaaaa ((q), _q0, (q), _q0, (n2), (n1));                        \
2789                                                                         \
2790     /* Compute the two most significant limbs of n - q'd */             \
2791     (r1) = (n1) - (d1) * (q);                                           \
2792     (r0) = (n0);                                                        \
2793     sub_ddmmss ((r1), (r0), (r1), (r0), (d1), (d0));                    \
2794     umul_ppmm (_t1, _t0, (d0), (q));                                    \
2795     sub_ddmmss ((r1), (r0), (r1), (r0), _t1, _t0);                      \
2796     (q)++;                                                              \
2797                                                                         \
2798     /* Conditionally adjust q and the remainders */                     \
2799     _mask = - (mp_limb_t) ((r1) >= _q0);                                \
2800     (q) += _mask;                                                       \
2801     add_ssaaaa ((r1), (r0), (r1), (r0), _mask & (d1), _mask & (d0));    \
2802     if (UNLIKELY ((r1) >= (d1)))                                        \
2803       {                                                                 \
2804         if ((r1) > (d1) || (r0) >= (d0))                                \
2805           {                                                             \
2806             (q)++;                                                      \
2807             sub_ddmmss ((r1), (r0), (r1), (r0), (d1), (d0));            \
2808           }                                                             \
2809       }                                                                 \
2810   } while (0)
2811
2812 #ifndef mpn_preinv_divrem_1  /* if not done with cpuvec in a fat binary */
2813 #define   mpn_preinv_divrem_1 __MPN(preinv_divrem_1)
2814 __GMP_DECLSPEC mp_limb_t mpn_preinv_divrem_1 __GMP_PROTO ((mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t, int));
2815 #endif
2816
2817
2818 /* USE_PREINV_DIVREM_1 is whether to use mpn_preinv_divrem_1, as opposed to the
2819    plain mpn_divrem_1.  The default is yes, since the few CISC chips where
2820    preinv is not good have defines saying so.  */
2821 #ifndef USE_PREINV_DIVREM_1
2822 #define USE_PREINV_DIVREM_1   1
2823 #endif
2824
2825 #if USE_PREINV_DIVREM_1
2826 #define MPN_DIVREM_OR_PREINV_DIVREM_1(qp,xsize,ap,size,d,dinv,shift)    \
2827   mpn_preinv_divrem_1 (qp, xsize, ap, size, d, dinv, shift)
2828 #else
2829 #define MPN_DIVREM_OR_PREINV_DIVREM_1(qp,xsize,ap,size,d,dinv,shift)    \
2830   mpn_divrem_1 (qp, xsize, ap, size, d)
2831 #endif
2832
2833 #ifndef PREINV_MOD_1_TO_MOD_1_THRESHOLD
2834 #define PREINV_MOD_1_TO_MOD_1_THRESHOLD 10
2835 #endif
2836
2837 /* This selection may seem backwards.  The reason mpn_mod_1 typically takes
2838    over for larger sizes is that it uses the mod_1_1 function.  */
2839 #define MPN_MOD_OR_PREINV_MOD_1(src,size,divisor,inverse)               \
2840   (BELOW_THRESHOLD (size, PREINV_MOD_1_TO_MOD_1_THRESHOLD)              \
2841    ? mpn_preinv_mod_1 (src, size, divisor, inverse)                     \
2842    : mpn_mod_1 (src, size, divisor))
2843
2844
2845 #ifndef mpn_mod_34lsub1  /* if not done with cpuvec in a fat binary */
2846 #define   mpn_mod_34lsub1 __MPN(mod_34lsub1)
2847 __GMP_DECLSPEC mp_limb_t mpn_mod_34lsub1 __GMP_PROTO ((mp_srcptr, mp_size_t)) __GMP_ATTRIBUTE_PURE;
2848 #endif
2849
2850
2851 /* DIVEXACT_1_THRESHOLD is at what size to use mpn_divexact_1, as opposed to
2852    plain mpn_divrem_1.  Likewise BMOD_1_TO_MOD_1_THRESHOLD for
2853    mpn_modexact_1_odd against plain mpn_mod_1.  On most CPUs divexact and
2854    modexact are faster at all sizes, so the defaults are 0.  Those CPUs
2855    where this is not right have a tuned threshold.  */
2856 #ifndef DIVEXACT_1_THRESHOLD
2857 #define DIVEXACT_1_THRESHOLD  0
2858 #endif
2859 #ifndef BMOD_1_TO_MOD_1_THRESHOLD
2860 #define BMOD_1_TO_MOD_1_THRESHOLD  10
2861 #endif
2862
2863 #ifndef mpn_divexact_1  /* if not done with cpuvec in a fat binary */
2864 #define mpn_divexact_1 __MPN(divexact_1)
2865 __GMP_DECLSPEC void    mpn_divexact_1 __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_limb_t));
2866 #endif
2867
2868 #define MPN_DIVREM_OR_DIVEXACT_1(dst, src, size, divisor)                     \
2869   do {                                                                        \
2870     if (BELOW_THRESHOLD (size, DIVEXACT_1_THRESHOLD))                         \
2871       ASSERT_NOCARRY (mpn_divrem_1 (dst, (mp_size_t) 0, src, size, divisor)); \
2872     else                                                                      \
2873       {                                                                       \
2874         ASSERT (mpn_mod_1 (src, size, divisor) == 0);                         \
2875         mpn_divexact_1 (dst, src, size, divisor);                             \
2876       }                                                                       \
2877   } while (0)
2878
2879 #ifndef mpn_modexact_1c_odd  /* if not done with cpuvec in a fat binary */
2880 #define   mpn_modexact_1c_odd __MPN(modexact_1c_odd)
2881 __GMP_DECLSPEC mp_limb_t mpn_modexact_1c_odd __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t, mp_limb_t)) __GMP_ATTRIBUTE_PURE;
2882 #endif
2883
2884 #if HAVE_NATIVE_mpn_modexact_1_odd
2885 #define   mpn_modexact_1_odd  __MPN(modexact_1_odd)
2886 __GMP_DECLSPEC mp_limb_t mpn_modexact_1_odd __GMP_PROTO ((mp_srcptr, mp_size_t, mp_limb_t)) __GMP_ATTRIBUTE_PURE;
2887 #else
2888 #define mpn_modexact_1_odd(src,size,divisor) \
2889   mpn_modexact_1c_odd (src, size, divisor, CNST_LIMB(0))
2890 #endif
2891
2892 #define MPN_MOD_OR_MODEXACT_1_ODD(src,size,divisor)                     \
2893   (BELOW_THRESHOLD (size, BMOD_1_TO_MOD_1_THRESHOLD)                    \
2894    ? mpn_modexact_1_odd (src, size, divisor)                            \
2895    : mpn_mod_1 (src, size, divisor))
2896
2897 /* binvert_limb() sets inv to the multiplicative inverse of n modulo
2898    2^GMP_NUMB_BITS, ie. satisfying inv*n == 1 mod 2^GMP_NUMB_BITS.
2899    n must be odd (otherwise such an inverse doesn't exist).
2900
2901    This is not to be confused with invert_limb(), which is completely
2902    different.
2903
2904    The table lookup gives an inverse with the low 8 bits valid, and each
2905    multiply step doubles the number of bits.  See Jebelean "An algorithm for
2906    exact division" end of section 4 (reference in gmp.texi).
2907
2908    Possible enhancement: Could use UHWtype until the last step, if half-size
2909    multiplies are faster (might help under _LONG_LONG_LIMB).
2910
2911    Alternative: As noted in Granlund and Montgomery "Division by Invariant
2912    Integers using Multiplication" (reference in gmp.texi), n itself gives a
2913    3-bit inverse immediately, and could be used instead of a table lookup.
2914    A 4-bit inverse can be obtained effectively from xoring bits 1 and 2 into
2915    bit 3, for instance with (((n + 2) & 4) << 1) ^ n.  */
2916
2917 #define binvert_limb_table  __gmp_binvert_limb_table
2918 __GMP_DECLSPEC extern const unsigned char  binvert_limb_table[128];
2919
2920 #define binvert_limb(inv,n)                                             \
2921   do {                                                                  \
2922     mp_limb_t  __n = (n);                                               \
2923     mp_limb_t  __inv;                                                   \
2924     ASSERT ((__n & 1) == 1);                                            \
2925                                                                         \
2926     __inv = binvert_limb_table[(__n/2) & 0x7F]; /*  8 */                \
2927     if (GMP_NUMB_BITS > 8)   __inv = 2 * __inv - __inv * __inv * __n;   \
2928     if (GMP_NUMB_BITS > 16)  __inv = 2 * __inv - __inv * __inv * __n;   \
2929     if (GMP_NUMB_BITS > 32)  __inv = 2 * __inv - __inv * __inv * __n;   \
2930                                                                         \
2931     if (GMP_NUMB_BITS > 64)                                             \
2932       {                                                                 \
2933         int  __invbits = 64;                                            \
2934         do {                                                            \
2935           __inv = 2 * __inv - __inv * __inv * __n;                      \
2936           __invbits *= 2;                                               \
2937         } while (__invbits < GMP_NUMB_BITS);                            \
2938       }                                                                 \
2939                                                                         \
2940     ASSERT ((__inv * __n & GMP_NUMB_MASK) == 1);                        \
2941     (inv) = __inv & GMP_NUMB_MASK;                                      \
2942   } while (0)
2943 #define modlimb_invert binvert_limb  /* backward compatibility */
2944
2945 /* Multiplicative inverse of 3, modulo 2^GMP_NUMB_BITS.
2946    Eg. 0xAAAAAAAB for 32 bits, 0xAAAAAAAAAAAAAAAB for 64 bits.
2947    GMP_NUMB_MAX/3*2+1 is right when GMP_NUMB_BITS is even, but when it's odd
2948    we need to start from GMP_NUMB_MAX>>1. */
2949 #define MODLIMB_INVERSE_3 (((GMP_NUMB_MAX >> (GMP_NUMB_BITS % 2)) / 3) * 2 + 1)
2950
2951 /* ceil(GMP_NUMB_MAX/3) and ceil(2*GMP_NUMB_MAX/3).
2952    These expressions work because GMP_NUMB_MAX%3 != 0 for all GMP_NUMB_BITS. */
2953 #define GMP_NUMB_CEIL_MAX_DIV3   (GMP_NUMB_MAX / 3 + 1)
2954 #define GMP_NUMB_CEIL_2MAX_DIV3  ((GMP_NUMB_MAX>>1) / 3 + 1 + GMP_NUMB_HIGHBIT)
2955
2956
2957 /* Set r to -a mod d.  a>=d is allowed.  Can give r>d.  All should be limbs.
2958
2959    It's not clear whether this is the best way to do this calculation.
2960    Anything congruent to -a would be fine for the one limb congruence
2961    tests.  */
2962
2963 #define NEG_MOD(r, a, d)                                                \
2964   do {                                                                  \
2965     ASSERT ((d) != 0);                                                  \
2966     ASSERT_LIMB (a);                                                    \
2967     ASSERT_LIMB (d);                                                    \
2968                                                                         \
2969     if ((a) <= (d))                                                     \
2970       {                                                                 \
2971         /* small a is reasonably likely */                              \
2972         (r) = (d) - (a);                                                \
2973       }                                                                 \
2974     else                                                                \
2975       {                                                                 \
2976         unsigned   __twos;                                              \
2977         mp_limb_t  __dnorm;                                             \
2978         count_leading_zeros (__twos, d);                                \
2979         __twos -= GMP_NAIL_BITS;                                        \
2980         __dnorm = (d) << __twos;                                        \
2981         (r) = ((a) <= __dnorm ? __dnorm : 2*__dnorm) - (a);             \
2982       }                                                                 \
2983                                                                         \
2984     ASSERT_LIMB (r);                                                    \
2985   } while (0)
2986
2987 /* A bit mask of all the least significant zero bits of n, or -1 if n==0. */
2988 #define LOW_ZEROS_MASK(n)  (((n) & -(n)) - 1)
2989
2990
2991 /* ULONG_PARITY sets "p" to 1 if there's an odd number of 1 bits in "n", or
2992    to 0 if there's an even number.  "n" should be an unsigned long and "p"
2993    an int.  */
2994
2995 #if defined (__GNUC__) && ! defined (NO_ASM) && HAVE_HOST_CPU_alpha_CIX
2996 #define ULONG_PARITY(p, n)                                              \
2997   do {                                                                  \
2998     int __p;                                                            \
2999     __asm__ ("ctpop %1, %0" : "=r" (__p) : "r" (n));                    \
3000     (p) = __p & 1;                                                      \
3001   } while (0)
3002 #endif
3003
3004 /* Cray intrinsic _popcnt. */
3005 #ifdef _CRAY
3006 #define ULONG_PARITY(p, n)      \
3007   do {                          \
3008     (p) = _popcnt (n) & 1;      \
3009   } while (0)
3010 #endif
3011
3012 #if defined (__GNUC__) && ! defined (__INTEL_COMPILER)                  \
3013     && ! defined (NO_ASM) && defined (__ia64)
3014 /* unsigned long is either 32 or 64 bits depending on the ABI, zero extend
3015    to a 64 bit unsigned long long for popcnt */
3016 #define ULONG_PARITY(p, n)                                              \
3017   do {                                                                  \
3018     unsigned long long  __n = (unsigned long) (n);                      \
3019     int  __p;                                                           \
3020     __asm__ ("popcnt %0 = %1" : "=r" (__p) : "r" (__n));                \
3021     (p) = __p & 1;                                                      \
3022   } while (0)
3023 #endif
3024
3025 #if defined (__GNUC__) && ! defined (__INTEL_COMPILER)                  \
3026     && ! defined (NO_ASM) && HAVE_HOST_CPU_FAMILY_x86
3027 #if __GMP_GNUC_PREREQ (3,1)
3028 #define __GMP_qm "=Qm"
3029 #define __GMP_q "=Q"
3030 #else
3031 #define __GMP_qm "=qm"
3032 #define __GMP_q "=q"
3033 #endif
3034 #define ULONG_PARITY(p, n)                                              \
3035   do {                                                                  \
3036     char           __p;                                                 \
3037     unsigned long  __n = (n);                                           \
3038     __n ^= (__n >> 16);                                                 \
3039     __asm__ ("xorb %h1, %b1\n\t"                                        \
3040              "setpo %0"                                                 \
3041          : __GMP_qm (__p), __GMP_q (__n)                                \
3042          : "1" (__n));                                                  \
3043     (p) = __p;                                                          \
3044   } while (0)
3045 #endif
3046
3047 #if ! defined (ULONG_PARITY)
3048 #define ULONG_PARITY(p, n)                                              \
3049   do {                                                                  \
3050     unsigned long  __n = (n);                                           \
3051     int  __p = 0;                                                       \
3052     do                                                                  \
3053       {                                                                 \
3054         __p ^= 0x96696996L >> (__n & 0x1F);                             \
3055         __n >>= 5;                                                      \
3056       }                                                                 \
3057     while (__n != 0);                                                   \
3058                                                                         \
3059     (p) = __p & 1;                                                      \
3060   } while (0)
3061 #endif
3062
3063
3064 /* 3 cycles on 604 or 750 since shifts and rlwimi's can pair.  gcc (as of
3065    version 3.1 at least) doesn't seem to know how to generate rlwimi for
3066    anything other than bit-fields, so use "asm".  */
3067 #if defined (__GNUC__) && ! defined (NO_ASM)                    \
3068   && HAVE_HOST_CPU_FAMILY_powerpc && GMP_LIMB_BITS == 32
3069 #define BSWAP_LIMB(dst, src)                                            \
3070   do {                                                                  \
3071     mp_limb_t  __bswapl_src = (src);                                    \
3072     mp_limb_t  __tmp1 = __bswapl_src >> 24;             /* low byte */  \
3073     mp_limb_t  __tmp2 = __bswapl_src << 24;             /* high byte */ \
3074     __asm__ ("rlwimi %0, %2, 24, 16, 23"                /* 2nd low */   \
3075          : "=r" (__tmp1) : "0" (__tmp1), "r" (__bswapl_src));           \
3076     __asm__ ("rlwimi %0, %2,  8,  8, 15"                /* 3nd high */  \
3077          : "=r" (__tmp2) : "0" (__tmp2), "r" (__bswapl_src));           \
3078     (dst) = __tmp1 | __tmp2;                            /* whole */     \
3079   } while (0)
3080 #endif
3081
3082 /* bswap is available on i486 and up and is fast.  A combination rorw $8 /
3083    roll $16 / rorw $8 is used in glibc for plain i386 (and in the linux
3084    kernel with xchgb instead of rorw), but this is not done here, because
3085    i386 means generic x86 and mixing word and dword operations will cause
3086    partial register stalls on P6 chips.  */
3087 #if defined (__GNUC__) && ! defined (NO_ASM)            \
3088   && HAVE_HOST_CPU_FAMILY_x86 && ! HAVE_HOST_CPU_i386   \
3089   && GMP_LIMB_BITS == 32
3090 #define BSWAP_LIMB(dst, src)                                            \
3091   do {                                                                  \
3092     __asm__ ("bswap %0" : "=r" (dst) : "0" (src));                      \
3093   } while (0)
3094 #endif
3095
3096 #if defined (__GNUC__) && ! defined (NO_ASM)            \
3097   && defined (__amd64__) && GMP_LIMB_BITS == 64
3098 #define BSWAP_LIMB(dst, src)                                            \
3099   do {                                                                  \
3100     __asm__ ("bswap %q0" : "=r" (dst) : "0" (src));                     \
3101   } while (0)
3102 #endif
3103
3104 #if defined (__GNUC__) && ! defined (__INTEL_COMPILER)                  \
3105     && ! defined (NO_ASM) && defined (__ia64) && GMP_LIMB_BITS == 64
3106 #define BSWAP_LIMB(dst, src)                                            \
3107   do {                                                                  \
3108     __asm__ ("mux1 %0 = %1, @rev" : "=r" (dst) :  "r" (src));           \
3109   } while (0)
3110 #endif
3111
3112 /* As per glibc. */
3113 #if defined (__GNUC__) && ! defined (NO_ASM)                    \
3114   && HAVE_HOST_CPU_FAMILY_m68k && GMP_LIMB_BITS == 32
3115 #define BSWAP_LIMB(dst, src)                                            \
3116   do {                                                                  \
3117     mp_limb_t  __bswapl_src = (src);                                    \
3118     __asm__ ("ror%.w %#8, %0\n\t"                                       \
3119              "swap   %0\n\t"                                            \
3120              "ror%.w %#8, %0"                                           \
3121              : "=d" (dst)                                               \
3122              : "0" (__bswapl_src));                                     \
3123   } while (0)
3124 #endif
3125
3126 #if ! defined (BSWAP_LIMB)
3127 #if GMP_LIMB_BITS == 8
3128 #define BSWAP_LIMB(dst, src)            \
3129   do { (dst) = (src); } while (0)
3130 #endif
3131 #if GMP_LIMB_BITS == 16
3132 #define BSWAP_LIMB(dst, src)                    \
3133   do {                                          \
3134     (dst) = ((src) << 8) + ((src) >> 8);        \
3135   } while (0)
3136 #endif
3137 #if GMP_LIMB_BITS == 32
3138 #define BSWAP_LIMB(dst, src)    \
3139   do {                          \
3140     (dst) =                     \
3141       ((src) << 24)             \
3142       + (((src) & 0xFF00) << 8) \
3143       + (((src) >> 8) & 0xFF00) \
3144       + ((src) >> 24);          \
3145   } while (0)
3146 #endif
3147 #if GMP_LIMB_BITS == 64
3148 #define BSWAP_LIMB(dst, src)            \
3149   do {                                  \
3150     (dst) =                             \
3151       ((src) << 56)                     \
3152       + (((src) & 0xFF00) << 40)        \
3153       + (((src) & 0xFF0000) << 24)      \
3154       + (((src) & 0xFF000000) << 8)     \
3155       + (((src) >> 8) & 0xFF000000)     \
3156       + (((src) >> 24) & 0xFF0000)      \
3157       + (((src) >> 40) & 0xFF00)        \
3158       + ((src) >> 56);                  \
3159   } while (0)
3160 #endif
3161 #endif
3162
3163 #if ! defined (BSWAP_LIMB)
3164 #define BSWAP_LIMB(dst, src)                            \
3165   do {                                                  \
3166     mp_limb_t  __bswapl_src = (src);                    \
3167     mp_limb_t  __dst = 0;                               \
3168     int        __i;                                     \
3169     for (__i = 0; __i < BYTES_PER_MP_LIMB; __i++)       \
3170       {                                                 \
3171         __dst = (__dst << 8) | (__bswapl_src & 0xFF);   \
3172         __bswapl_src >>= 8;                             \
3173       }                                                 \
3174     (dst) = __dst;                                      \
3175   } while (0)
3176 #endif
3177
3178
3179 /* Apparently lwbrx might be slow on some PowerPC chips, so restrict it to
3180    those we know are fast.  */
3181 #if defined (__GNUC__) && ! defined (NO_ASM)                            \
3182   && GMP_LIMB_BITS == 32 && HAVE_LIMB_BIG_ENDIAN                        \
3183   && (HAVE_HOST_CPU_powerpc604                                          \
3184       || HAVE_HOST_CPU_powerpc604e                                      \
3185       || HAVE_HOST_CPU_powerpc750                                       \
3186       || HAVE_HOST_CPU_powerpc7400)
3187 #define BSWAP_LIMB_FETCH(limb, src)                                     \
3188   do {                                                                  \
3189     mp_srcptr  __blf_src = (src);                                       \
3190     mp_limb_t  __limb;                                                  \
3191     __asm__ ("lwbrx %0, 0, %1"                                          \
3192              : "=r" (__limb)                                            \
3193              : "r" (__blf_src),                                         \
3194                "m" (*__blf_src));                                       \
3195     (limb) = __limb;                                                    \
3196   } while (0)
3197 #endif
3198
3199 #if ! defined (BSWAP_LIMB_FETCH)
3200 #define BSWAP_LIMB_FETCH(limb, src)  BSWAP_LIMB (limb, *(src))
3201 #endif
3202
3203
3204 /* On the same basis that lwbrx might be slow, restrict stwbrx to those we
3205    know are fast.  FIXME: Is this necessary?  */
3206 #if defined (__GNUC__) && ! defined (NO_ASM)                            \
3207   && GMP_LIMB_BITS == 32 && HAVE_LIMB_BIG_ENDIAN                        \
3208   && (HAVE_HOST_CPU_powerpc604                                          \
3209       || HAVE_HOST_CPU_powerpc604e                                      \
3210       || HAVE_HOST_CPU_powerpc750                                       \
3211       || HAVE_HOST_CPU_powerpc7400)
3212 #define BSWAP_LIMB_STORE(dst, limb)                                     \
3213   do {                                                                  \
3214     mp_ptr     __dst = (dst);                                           \
3215     mp_limb_t  __limb = (limb);                                         \
3216     __asm__ ("stwbrx %1, 0, %2"                                         \
3217              : "=m" (*__dst)                                            \
3218              : "r" (__limb),                                            \
3219                "r" (__dst));                                            \
3220   } while (0)
3221 #endif
3222
3223 #if ! defined (BSWAP_LIMB_STORE)
3224 #define BSWAP_LIMB_STORE(dst, limb)  BSWAP_LIMB (*(dst), limb)
3225 #endif
3226
3227
3228 /* Byte swap limbs from {src,size} and store at {dst,size}. */
3229 #define MPN_BSWAP(dst, src, size)                       \
3230   do {                                                  \
3231     mp_ptr     __dst = (dst);                           \
3232     mp_srcptr  __src = (src);                           \
3233     mp_size_t  __size = (size);                         \
3234     mp_size_t  __i;                                     \
3235     ASSERT ((size) >= 0);                               \
3236     ASSERT (MPN_SAME_OR_SEPARATE_P (dst, src, size));   \
3237     CRAY_Pragma ("_CRI ivdep");                         \
3238     for (__i = 0; __i < __size; __i++)                  \
3239       {                                                 \
3240         BSWAP_LIMB_FETCH (*__dst, __src);               \
3241         __dst++;                                        \
3242         __src++;                                        \
3243       }                                                 \
3244   } while (0)
3245
3246 /* Byte swap limbs from {dst,size} and store in reverse order at {src,size}. */
3247 #define MPN_BSWAP_REVERSE(dst, src, size)               \
3248   do {                                                  \
3249     mp_ptr     __dst = (dst);                           \
3250     mp_size_t  __size = (size);                         \
3251     mp_srcptr  __src = (src) + __size - 1;              \
3252     mp_size_t  __i;                                     \
3253     ASSERT ((size) >= 0);                               \
3254     ASSERT (! MPN_OVERLAP_P (dst, size, src, size));    \
3255     CRAY_Pragma ("_CRI ivdep");                         \
3256     for (__i = 0; __i < __size; __i++)                  \
3257       {                                                 \
3258         BSWAP_LIMB_FETCH (*__dst, __src);               \
3259         __dst++;                                        \
3260         __src--;                                        \
3261       }                                                 \
3262   } while (0)
3263
3264
3265 /* No processor claiming to be SPARC v9 compliant seems to
3266    implement the POPC instruction.  Disable pattern for now.  */
3267 #if 0
3268 #if defined __GNUC__ && defined __sparc_v9__ && GMP_LIMB_BITS == 64
3269 #define popc_limb(result, input)                                        \
3270   do {                                                                  \
3271     DItype __res;                                                       \
3272     __asm__ ("popc %1,%0" : "=r" (result) : "rI" (input));              \
3273   } while (0)
3274 #endif
3275 #endif
3276
3277 #if defined (__GNUC__) && ! defined (NO_ASM) && HAVE_HOST_CPU_alpha_CIX
3278 #define popc_limb(result, input)                                        \
3279   do {                                                                  \
3280     __asm__ ("ctpop %1, %0" : "=r" (result) : "r" (input));             \
3281   } while (0)
3282 #endif
3283
3284 /* Cray intrinsic. */
3285 #ifdef _CRAY
3286 #define popc_limb(result, input)        \
3287   do {                                  \
3288     (result) = _popcnt (input);         \
3289   } while (0)
3290 #endif
3291
3292 #if defined (__GNUC__) && ! defined (__INTEL_COMPILER)                  \
3293     && ! defined (NO_ASM) && defined (__ia64) && GMP_LIMB_BITS == 64
3294 #define popc_limb(result, input)                                        \
3295   do {                                                                  \
3296     __asm__ ("popcnt %0 = %1" : "=r" (result) : "r" (input));           \
3297   } while (0)
3298 #endif
3299
3300 /* Cool population count of an mp_limb_t.
3301    You have to figure out how this works, We won't tell you!
3302
3303    The constants could also be expressed as:
3304      0x55... = [2^N / 3]     = [(2^N-1)/3]
3305      0x33... = [2^N / 5]     = [(2^N-1)/5]
3306      0x0f... = [2^N / 17]    = [(2^N-1)/17]
3307      (N is GMP_LIMB_BITS, [] denotes truncation.) */
3308
3309 #if ! defined (popc_limb) && GMP_LIMB_BITS == 8
3310 #define popc_limb(result, input)                                        \
3311   do {                                                                  \
3312     mp_limb_t  __x = (input);                                           \
3313     __x -= (__x >> 1) & MP_LIMB_T_MAX/3;                                \
3314     __x = ((__x >> 2) & MP_LIMB_T_MAX/5) + (__x & MP_LIMB_T_MAX/5);     \
3315     __x = ((__x >> 4) + __x) & MP_LIMB_T_MAX/17;                        \
3316     (result) = __x & 0xff;                                              \
3317   } while (0)
3318 #endif
3319
3320 #if ! defined (popc_limb) && GMP_LIMB_BITS == 16
3321 #define popc_limb(result, input)                                        \
3322   do {                                                                  \
3323     mp_limb_t  __x = (input);                                           \
3324     __x -= (__x >> 1) & MP_LIMB_T_MAX/3;                                \
3325     __x = ((__x >> 2) & MP_LIMB_T_MAX/5) + (__x & MP_LIMB_T_MAX/5);     \
3326     __x = ((__x >> 4) + __x) & MP_LIMB_T_MAX/17;                        \
3327     __x = ((__x >> 8) + __x);                                           \
3328     (result) = __x & 0xff;                                              \
3329   } while (0)
3330 #endif
3331
3332 #if ! defined (popc_limb) && GMP_LIMB_BITS == 32
3333 #define popc_limb(result, input)                                        \
3334   do {                                                                  \
3335     mp_limb_t  __x = (input);                                           \
3336     __x -= (__x >> 1) & MP_LIMB_T_MAX/3;                                \
3337     __x = ((__x >> 2) & MP_LIMB_T_MAX/5) + (__x & MP_LIMB_T_MAX/5);     \
3338     __x = ((__x >> 4) + __x) & MP_LIMB_T_MAX/17;                        \
3339     __x = ((__x >> 8) + __x);                                           \
3340     __x = ((__x >> 16) + __x);                                          \
3341     (result) = __x & 0xff;                                              \
3342   } while (0)
3343 #endif
3344
3345 #if ! defined (popc_limb) && GMP_LIMB_BITS == 64
3346 #define popc_limb(result, input)                                        \
3347   do {                                                                  \
3348     mp_limb_t  __x = (input);                                           \
3349     __x -= (__x >> 1) & MP_LIMB_T_MAX/3;                                \
3350     __x = ((__x >> 2) & MP_LIMB_T_MAX/5) + (__x & MP_LIMB_T_MAX/5);     \
3351     __x = ((__x >> 4) + __x) & MP_LIMB_T_MAX/17;                        \
3352     __x = ((__x >> 8) + __x);                                           \
3353     __x = ((__x >> 16) + __x);                                          \
3354     __x = ((__x >> 32) + __x);                                          \
3355     (result) = __x & 0xff;                                              \
3356   } while (0)
3357 #endif
3358
3359
3360 /* Define stuff for longlong.h.  */
3361 #if HAVE_ATTRIBUTE_MODE
3362 typedef unsigned int UQItype    __attribute__ ((mode (QI)));
3363 typedef          int SItype     __attribute__ ((mode (SI)));
3364 typedef unsigned int USItype    __attribute__ ((mode (SI)));
3365 typedef          int DItype     __attribute__ ((mode (DI)));
3366 typedef unsigned int UDItype    __attribute__ ((mode (DI)));
3367 #else
3368 typedef unsigned char UQItype;
3369 typedef          long SItype;
3370 typedef unsigned long USItype;
3371 #if HAVE_LONG_LONG
3372 typedef long long int DItype;
3373 typedef unsigned long long int UDItype;
3374 #else /* Assume `long' gives us a wide enough type.  Needed for hppa2.0w.  */
3375 typedef long int DItype;
3376 typedef unsigned long int UDItype;
3377 #endif
3378 #endif
3379
3380 typedef mp_limb_t UWtype;
3381 typedef unsigned int UHWtype;
3382 #define W_TYPE_SIZE GMP_LIMB_BITS
3383
3384 /* Define ieee_double_extract and _GMP_IEEE_FLOATS.
3385
3386    Bit field packing is "implementation defined" according to C99, which
3387    leaves us at the compiler's mercy here.  For some systems packing is
3388    defined in the ABI (eg. x86).  In any case so far it seems universal that
3389    little endian systems pack from low to high, and big endian from high to
3390    low within the given type.
3391
3392    Within the fields we rely on the integer endianness being the same as the
3393    float endianness, this is true everywhere we know of and it'd be a fairly
3394    strange system that did anything else.  */
3395
3396 #if HAVE_DOUBLE_IEEE_LITTLE_SWAPPED
3397 #define _GMP_IEEE_FLOATS 1
3398 union ieee_double_extract
3399 {
3400   struct
3401     {
3402       gmp_uint_least32_t manh:20;
3403       gmp_uint_least32_t exp:11;
3404       gmp_uint_least32_t sig:1;
3405       gmp_uint_least32_t manl:32;
3406     } s;
3407   double d;
3408 };
3409 #endif
3410
3411 #if HAVE_DOUBLE_IEEE_LITTLE_ENDIAN
3412 #define _GMP_IEEE_FLOATS 1
3413 union ieee_double_extract
3414 {
3415   struct
3416     {
3417       gmp_uint_least32_t manl:32;
3418       gmp_uint_least32_t manh:20;
3419       gmp_uint_least32_t exp:11;
3420       gmp_uint_least32_t sig:1;
3421     } s;
3422   double d;
3423 };
3424 #endif
3425
3426 #if HAVE_DOUBLE_IEEE_BIG_ENDIAN
3427 #define _GMP_IEEE_FLOATS 1
3428 union ieee_double_extract
3429 {
3430   struct
3431     {
3432       gmp_uint_least32_t sig:1;
3433       gmp_uint_least32_t exp:11;
3434       gmp_uint_least32_t manh:20;
3435       gmp_uint_least32_t manl:32;
3436     } s;
3437   double d;
3438 };
3439 #endif
3440
3441
3442 /* Use (4.0 * ...) instead of (2.0 * ...) to work around buggy compilers
3443    that don't convert ulong->double correctly (eg. SunOS 4 native cc).  */
3444 #define MP_BASE_AS_DOUBLE (4.0 * ((mp_limb_t) 1 << (GMP_NUMB_BITS - 2)))
3445 /* Maximum number of limbs it will take to store any `double'.
3446    We assume doubles have 53 mantissa bits.  */
3447 #define LIMBS_PER_DOUBLE ((53 + GMP_NUMB_BITS - 2) / GMP_NUMB_BITS + 1)
3448
3449 __GMP_DECLSPEC int __gmp_extract_double __GMP_PROTO ((mp_ptr, double));
3450
3451 #define mpn_get_d __gmpn_get_d
3452 __GMP_DECLSPEC double mpn_get_d __GMP_PROTO ((mp_srcptr, mp_size_t, mp_size_t, long)) __GMP_ATTRIBUTE_PURE;
3453
3454
3455 /* DOUBLE_NAN_INF_ACTION executes code a_nan if x is a NaN, or executes
3456    a_inf if x is an infinity.  Both are considered unlikely values, for
3457    branch prediction.  */
3458
3459 #if _GMP_IEEE_FLOATS
3460 #define DOUBLE_NAN_INF_ACTION(x, a_nan, a_inf)  \
3461   do {                                          \
3462     union ieee_double_extract  u;               \
3463     u.d = (x);                                  \
3464     if (UNLIKELY (u.s.exp == 0x7FF))            \
3465       {                                         \
3466         if (u.s.manl == 0 && u.s.manh == 0)     \
3467           { a_inf; }                            \
3468         else                                    \
3469           { a_nan; }                            \
3470       }                                         \
3471   } while (0)
3472 #endif
3473
3474 #if HAVE_DOUBLE_VAX_D || HAVE_DOUBLE_VAX_G || HAVE_DOUBLE_CRAY_CFP
3475 /* no nans or infs in these formats */
3476 #define DOUBLE_NAN_INF_ACTION(x, a_nan, a_inf)  \
3477   do { } while (0)
3478 #endif
3479
3480 #ifndef DOUBLE_NAN_INF_ACTION
3481 /* Unknown format, try something generic.
3482    NaN should be "unordered", so x!=x.
3483    Inf should be bigger than DBL_MAX.  */
3484 #define DOUBLE_NAN_INF_ACTION(x, a_nan, a_inf)                  \
3485   do {                                                          \
3486     {                                                           \
3487       if (UNLIKELY ((x) != (x)))                                \
3488         { a_nan; }                                              \
3489       else if (UNLIKELY ((x) > DBL_MAX || (x) < -DBL_MAX))      \
3490         { a_inf; }                                              \
3491     }                                                           \
3492   } while (0)
3493 #endif
3494
3495 /* On m68k, x86 and amd64, gcc (and maybe other compilers) can hold doubles
3496    in the coprocessor, which means a bigger exponent range than normal, and
3497    depending on the rounding mode, a bigger mantissa than normal.  (See
3498    "Disappointments" in the gcc manual.)  FORCE_DOUBLE stores and fetches
3499    "d" through memory to force any rounding and overflows to occur.
3500
3501    On amd64, and on x86s with SSE2, gcc (depending on options) uses the xmm
3502    registers, where there's no such extra precision and no need for the
3503    FORCE_DOUBLE.  We don't bother to detect this since the present uses for
3504    FORCE_DOUBLE are only in test programs and default generic C code.
3505
3506    Not quite sure that an "automatic volatile" will use memory, but it does
3507    in gcc.  An asm("":"=m"(d):"0"(d)) can't be used to trick gcc, since
3508    apparently matching operands like "0" are only allowed on a register
3509    output.  gcc 3.4 warns about this, though in fact it and past versions
3510    seem to put the operand through memory as hoped.  */
3511
3512 #if (HAVE_HOST_CPU_FAMILY_m68k || HAVE_HOST_CPU_FAMILY_x86      \
3513      || defined (__amd64__))
3514 #define FORCE_DOUBLE(d) \
3515   do { volatile double __gmp_force = (d); (d) = __gmp_force; } while (0)
3516 #else
3517 #define FORCE_DOUBLE(d)  do { } while (0)
3518 #endif
3519
3520
3521 __GMP_DECLSPEC extern int __gmp_junk;
3522 __GMP_DECLSPEC extern const int __gmp_0;
3523 __GMP_DECLSPEC void __gmp_exception __GMP_PROTO ((int)) ATTRIBUTE_NORETURN;
3524 __GMP_DECLSPEC void __gmp_divide_by_zero __GMP_PROTO ((void)) ATTRIBUTE_NORETURN;
3525 __GMP_DECLSPEC void __gmp_sqrt_of_negative __GMP_PROTO ((void)) ATTRIBUTE_NORETURN;
3526 __GMP_DECLSPEC void __gmp_invalid_operation __GMP_PROTO ((void)) ATTRIBUTE_NORETURN;
3527 #define GMP_ERROR(code)   __gmp_exception (code)
3528 #define DIVIDE_BY_ZERO    __gmp_divide_by_zero ()
3529 #define SQRT_OF_NEGATIVE  __gmp_sqrt_of_negative ()
3530
3531 #if defined _LONG_LONG_LIMB
3532 #if __GMP_HAVE_TOKEN_PASTE
3533 #define CNST_LIMB(C) ((mp_limb_t) C##LL)
3534 #else
3535 #define CNST_LIMB(C) ((mp_limb_t) C/**/LL)
3536 #endif
3537 #else /* not _LONG_LONG_LIMB */
3538 #if __GMP_HAVE_TOKEN_PASTE
3539 #define CNST_LIMB(C) ((mp_limb_t) C##L)
3540 #else
3541 #define CNST_LIMB(C) ((mp_limb_t) C/**/L)
3542 #endif
3543 #endif /* _LONG_LONG_LIMB */
3544
3545 /* Stuff used by mpn/generic/perfsqr.c and mpz/prime_p.c */
3546 #if GMP_NUMB_BITS == 2
3547 #define PP 0x3                                  /* 3 */
3548 #define PP_FIRST_OMITTED 5
3549 #endif
3550 #if GMP_NUMB_BITS == 4
3551 #define PP 0xF                                  /* 3 x 5 */
3552 #define PP_FIRST_OMITTED 7
3553 #endif
3554 #if GMP_NUMB_BITS == 8
3555 #define PP 0x69                                 /* 3 x 5 x 7 */
3556 #define PP_FIRST_OMITTED 11
3557 #endif
3558 #if GMP_NUMB_BITS == 16
3559 #define PP 0x3AA7                               /* 3 x 5 x 7 x 11 x 13 */
3560 #define PP_FIRST_OMITTED 17
3561 #endif
3562 #if GMP_NUMB_BITS == 32
3563 #define PP 0xC0CFD797L                          /* 3 x 5 x 7 x 11 x ... x 29 */
3564 #define PP_INVERTED 0x53E5645CL
3565 #define PP_FIRST_OMITTED 31
3566 #endif
3567 #if GMP_NUMB_BITS == 64
3568 #define PP CNST_LIMB(0xE221F97C30E94E1D)        /* 3 x 5 x 7 x 11 x ... x 53 */
3569 #define PP_INVERTED CNST_LIMB(0x21CFE6CFC938B36B)
3570 #define PP_FIRST_OMITTED 59
3571 #endif
3572 #ifndef PP_FIRST_OMITTED
3573 #define PP_FIRST_OMITTED 3
3574 #endif
3575
3576
3577
3578 /* BIT1 means a result value in bit 1 (second least significant bit), with a
3579    zero bit representing +1 and a one bit representing -1.  Bits other than
3580    bit 1 are garbage.  These are meant to be kept in "int"s, and casts are
3581    used to ensure the expressions are "int"s even if a and/or b might be
3582    other types.
3583
3584    JACOBI_TWOS_U_BIT1 and JACOBI_RECIP_UU_BIT1 are used in mpn_jacobi_base
3585    and their speed is important.  Expressions are used rather than
3586    conditionals to accumulate sign changes, which effectively means XORs
3587    instead of conditional JUMPs. */
3588
3589 /* (a/0), with a signed; is 1 if a=+/-1, 0 otherwise */
3590 #define JACOBI_S0(a)   (((a) == 1) | ((a) == -1))
3591
3592 /* (a/0), with a unsigned; is 1 if a=+/-1, 0 otherwise */
3593 #define JACOBI_U0(a)   ((a) == 1)
3594
3595 /* (a/0), with a given by low and size;
3596    is 1 if a=+/-1, 0 otherwise */
3597 #define JACOBI_LS0(alow,asize) \
3598   (((asize) == 1 || (asize) == -1) && (alow) == 1)
3599
3600 /* (a/0), with a an mpz_t;
3601    fetch of low limb always valid, even if size is zero */
3602 #define JACOBI_Z0(a)   JACOBI_LS0 (PTR(a)[0], SIZ(a))
3603
3604 /* (0/b), with b unsigned; is 1 if b=1, 0 otherwise */
3605 #define JACOBI_0U(b)   ((b) == 1)
3606
3607 /* (0/b), with b unsigned; is 1 if b=+/-1, 0 otherwise */
3608 #define JACOBI_0S(b)   ((b) == 1 || (b) == -1)
3609
3610 /* (0/b), with b given by low and size; is 1 if b=+/-1, 0 otherwise */
3611 #define JACOBI_0LS(blow,bsize) \
3612   (((bsize) == 1 || (bsize) == -1) && (blow) == 1)
3613
3614 /* Convert a bit1 to +1 or -1. */
3615 #define JACOBI_BIT1_TO_PN(result_bit1) \
3616   (1 - ((int) (result_bit1) & 2))
3617
3618 /* (2/b), with b unsigned and odd;
3619    is (-1)^((b^2-1)/8) which is 1 if b==1,7mod8 or -1 if b==3,5mod8 and
3620    hence obtained from (b>>1)^b */
3621 #define JACOBI_TWO_U_BIT1(b) \
3622   ((int) (((b) >> 1) ^ (b)))
3623
3624 /* (2/b)^twos, with b unsigned and odd */
3625 #define JACOBI_TWOS_U_BIT1(twos, b) \
3626   ((int) ((twos) << 1) & JACOBI_TWO_U_BIT1 (b))
3627
3628 /* (2/b)^twos, with b unsigned and odd */
3629 #define JACOBI_TWOS_U(twos, b) \
3630   (JACOBI_BIT1_TO_PN (JACOBI_TWOS_U_BIT1 (twos, b)))
3631
3632 /* (-1/b), with b odd (signed or unsigned);
3633    is (-1)^((b-1)/2) */
3634 #define JACOBI_N1B_BIT1(b) \
3635   ((int) (b))
3636
3637 /* (a/b) effect due to sign of a: signed/unsigned, b odd;
3638    is (-1/b) if a<0, or +1 if a>=0 */
3639 #define JACOBI_ASGN_SU_BIT1(a, b) \
3640   ((((a) < 0) << 1) & JACOBI_N1B_BIT1(b))
3641
3642 /* (a/b) effect due to sign of b: signed/signed;
3643    is -1 if a and b both negative, +1 otherwise */
3644 #define JACOBI_BSGN_SS_BIT1(a, b) \
3645   ((((a)<0) & ((b)<0)) << 1)
3646
3647 /* (a/b) effect due to sign of b: signed/mpz;
3648    is -1 if a and b both negative, +1 otherwise */
3649 #define JACOBI_BSGN_SZ_BIT1(a, b) \
3650   JACOBI_BSGN_SS_BIT1 (a, SIZ(b))
3651
3652 /* (a/b) effect due to sign of b: mpz/signed;
3653    is -1 if a and b both negative, +1 otherwise */
3654 #define JACOBI_BSGN_ZS_BIT1(a, b) \
3655   JACOBI_BSGN_SZ_BIT1 (b, a)
3656
3657 /* (a/b) reciprocity to switch to (b/a), a,b both unsigned and odd;
3658    is (-1)^((a-1)*(b-1)/4), which means +1 if either a,b==1mod4, or -1 if
3659    both a,b==3mod4, achieved in bit 1 by a&b.  No ASSERT()s about a,b odd
3660    because this is used in a couple of places with only bit 1 of a or b
3661    valid. */
3662 #define JACOBI_RECIP_UU_BIT1(a, b) \
3663   ((int) ((a) & (b)))
3664
3665 /* Strip low zero limbs from {b_ptr,b_size} by incrementing b_ptr and
3666    decrementing b_size.  b_low should be b_ptr[0] on entry, and will be
3667    updated for the new b_ptr.  result_bit1 is updated according to the
3668    factors of 2 stripped, as per (a/2).  */
3669 #define JACOBI_STRIP_LOW_ZEROS(result_bit1, a, b_ptr, b_size, b_low)    \
3670   do {                                                                  \
3671     ASSERT ((b_size) >= 1);                                             \
3672     ASSERT ((b_low) == (b_ptr)[0]);                                     \
3673                                                                         \
3674     while (UNLIKELY ((b_low) == 0))                                     \
3675       {                                                                 \
3676         (b_size)--;                                                     \
3677         ASSERT ((b_size) >= 1);                                         \
3678         (b_ptr)++;                                                      \
3679         (b_low) = *(b_ptr);                                             \
3680                                                                         \
3681         ASSERT (((a) & 1) != 0);                                        \
3682         if ((GMP_NUMB_BITS % 2) == 1)                                   \
3683           (result_bit1) ^= JACOBI_TWO_U_BIT1(a);                        \
3684       }                                                                 \
3685   } while (0)
3686
3687 /* Set a_rem to {a_ptr,a_size} reduced modulo b, either using mod_1 or
3688    modexact_1_odd, but in either case leaving a_rem<b.  b must be odd and
3689    unsigned.  modexact_1_odd effectively calculates -a mod b, and
3690    result_bit1 is adjusted for the factor of -1.
3691
3692    The way mpn_modexact_1_odd sometimes bases its remainder on a_size and
3693    sometimes on a_size-1 means if GMP_NUMB_BITS is odd we can't know what
3694    factor to introduce into result_bit1, so for that case use mpn_mod_1
3695    unconditionally.
3696
3697    FIXME: mpn_modexact_1_odd is more efficient, so some way to get it used
3698    for odd GMP_NUMB_BITS would be good.  Perhaps it could mung its result,
3699    or not skip a divide step, or something. */
3700
3701 #define JACOBI_MOD_OR_MODEXACT_1_ODD(result_bit1, a_rem, a_ptr, a_size, b) \
3702   do {                                                                     \
3703     mp_srcptr  __a_ptr  = (a_ptr);                                         \
3704     mp_size_t  __a_size = (a_size);                                        \
3705     mp_limb_t  __b      = (b);                                             \
3706                                                                            \
3707     ASSERT (__a_size >= 1);                                                \
3708     ASSERT (__b & 1);                                                      \
3709                                                                            \
3710     if ((GMP_NUMB_BITS % 2) != 0                                           \
3711         || ABOVE_THRESHOLD (__a_size, BMOD_1_TO_MOD_1_THRESHOLD))          \
3712       {                                                                    \
3713         (a_rem) = mpn_mod_1 (__a_ptr, __a_size, __b);                      \
3714       }                                                                    \
3715     else                                                                   \
3716       {                                                                    \
3717         (result_bit1) ^= JACOBI_N1B_BIT1 (__b);                            \
3718         (a_rem) = mpn_modexact_1_odd (__a_ptr, __a_size, __b);             \
3719       }                                                                    \
3720   } while (0)
3721
3722 /* Matrix multiplication */
3723 #define   mpn_matrix22_mul __MPN(matrix22_mul)
3724 __GMP_DECLSPEC void      mpn_matrix22_mul __GMP_PROTO ((mp_ptr, mp_ptr, mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_srcptr, mp_srcptr, mp_srcptr, mp_size_t, mp_ptr));
3725 #define   mpn_matrix22_mul_strassen __MPN(matrix22_mul_strassen)
3726 __GMP_DECLSPEC void      mpn_matrix22_mul_strassen __GMP_PROTO ((mp_ptr, mp_ptr, mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_srcptr, mp_srcptr, mp_srcptr, mp_size_t, mp_ptr));
3727 #define   mpn_matrix22_mul_itch __MPN(matrix22_mul_itch)
3728 __GMP_DECLSPEC mp_size_t mpn_matrix22_mul_itch __GMP_PROTO ((mp_size_t, mp_size_t));
3729
3730 #ifndef MATRIX22_STRASSEN_THRESHOLD
3731 #define MATRIX22_STRASSEN_THRESHOLD 30
3732 #endif
3733
3734 /* HGCD definitions */
3735
3736 /* Extract one numb, shifting count bits left
3737     ________  ________
3738    |___xh___||___xl___|
3739           |____r____|
3740    >count <
3741
3742    The count includes any nail bits, so it should work fine if count
3743    is computed using count_leading_zeros. If GMP_NAIL_BITS > 0, all of
3744    xh, xl and r include nail bits. Must have 0 < count < GMP_LIMB_BITS.
3745
3746    FIXME: Omit masking with GMP_NUMB_MASK, and let callers do that for
3747    those calls where the count high bits of xh may be non-zero.
3748 */
3749
3750 #define MPN_EXTRACT_NUMB(count, xh, xl)                         \
3751   ((((xh) << ((count) - GMP_NAIL_BITS)) & GMP_NUMB_MASK) |      \
3752    ((xl) >> (GMP_LIMB_BITS - (count))))
3753
3754
3755 /* The matrix non-negative M = (u, u'; v,v') keeps track of the
3756    reduction (a;b) = M (alpha; beta) where alpha, beta are smaller
3757    than a, b. The determinant must always be one, so that M has an
3758    inverse (v', -u'; -v, u). Elements always fit in GMP_NUMB_BITS - 1
3759    bits. */
3760 struct hgcd_matrix1
3761 {
3762   mp_limb_t u[2][2];
3763 };
3764
3765 #define mpn_hgcd2 __MPN (hgcd2)
3766 __GMP_DECLSPEC int mpn_hgcd2 __GMP_PROTO ((mp_limb_t, mp_limb_t, mp_limb_t, mp_limb_t,  struct hgcd_matrix1 *));
3767
3768 #define mpn_hgcd_mul_matrix1_vector __MPN (hgcd_mul_matrix1_vector)
3769 __GMP_DECLSPEC mp_size_t mpn_hgcd_mul_matrix1_vector __GMP_PROTO ((const struct hgcd_matrix1 *, mp_ptr, mp_srcptr, mp_ptr, mp_size_t));
3770
3771 #define mpn_hgcd_mul_matrix1_inverse_vector __MPN (hgcd_mul_matrix1_inverse_vector)
3772 __GMP_DECLSPEC mp_size_t mpn_hgcd_mul_matrix1_inverse_vector __GMP_PROTO ((const struct hgcd_matrix1 *, mp_ptr, mp_srcptr, mp_ptr, mp_size_t));
3773
3774 struct hgcd_matrix
3775 {
3776   mp_size_t alloc;              /* for sanity checking only */
3777   mp_size_t n;
3778   mp_ptr p[2][2];
3779 };
3780
3781 #define MPN_HGCD_MATRIX_INIT_ITCH(n) (4 * ((n+1)/2 + 1))
3782
3783 #define mpn_hgcd_matrix_init __MPN (hgcd_matrix_init)
3784 __GMP_DECLSPEC void mpn_hgcd_matrix_init __GMP_PROTO ((struct hgcd_matrix *, mp_size_t, mp_ptr));
3785
3786 #define mpn_hgcd_matrix_mul __MPN (hgcd_matrix_mul)
3787 __GMP_DECLSPEC void mpn_hgcd_matrix_mul __GMP_PROTO ((struct hgcd_matrix *, const struct hgcd_matrix *, mp_ptr));
3788
3789 #define mpn_hgcd_matrix_adjust __MPN (hgcd_matrix_adjust)
3790 __GMP_DECLSPEC mp_size_t mpn_hgcd_matrix_adjust __GMP_PROTO ((struct hgcd_matrix *, mp_size_t, mp_ptr, mp_ptr, mp_size_t, mp_ptr));
3791
3792 #define mpn_hgcd_itch __MPN (hgcd_itch)
3793 __GMP_DECLSPEC mp_size_t mpn_hgcd_itch __GMP_PROTO ((mp_size_t));
3794
3795 #define mpn_hgcd __MPN (hgcd)
3796 __GMP_DECLSPEC mp_size_t mpn_hgcd __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, struct hgcd_matrix *, mp_ptr));
3797
3798 #define MPN_HGCD_LEHMER_ITCH(n) (n)
3799
3800 #define mpn_hgcd_lehmer __MPN (hgcd_lehmer)
3801 __GMP_DECLSPEC mp_size_t mpn_hgcd_lehmer __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, struct hgcd_matrix *, mp_ptr));
3802
3803 /* Needs storage for the quotient */
3804 #define MPN_GCD_SUBDIV_STEP_ITCH(n) (n)
3805
3806 #define mpn_gcd_subdiv_step __MPN(gcd_subdiv_step)
3807 __GMP_DECLSPEC mp_size_t mpn_gcd_subdiv_step __GMP_PROTO ((mp_ptr, mp_size_t *, mp_ptr, mp_ptr, mp_size_t, mp_ptr));
3808
3809 #define MPN_GCD_LEHMER_N_ITCH(n) (n)
3810
3811 #define mpn_gcd_lehmer_n __MPN(gcd_lehmer_n)
3812 __GMP_DECLSPEC mp_size_t mpn_gcd_lehmer_n __GMP_PROTO ((mp_ptr, mp_ptr, mp_ptr, mp_size_t, mp_ptr));
3813
3814 #define mpn_gcdext_subdiv_step __MPN(gcdext_subdiv_step)
3815 __GMP_DECLSPEC mp_size_t mpn_gcdext_subdiv_step __GMP_PROTO ((mp_ptr, mp_size_t *, mp_ptr, mp_size_t *, mp_ptr, mp_ptr, mp_size_t, mp_ptr, mp_ptr, mp_size_t *, mp_ptr, mp_ptr));
3816
3817 #define MPN_GCDEXT_LEHMER_N_ITCH(n) (4*(n) + 3)
3818
3819 #define mpn_gcdext_lehmer_n __MPN(gcdext_lehmer_n)
3820 __GMP_DECLSPEC mp_size_t mpn_gcdext_lehmer_n __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t *, mp_ptr, mp_ptr, mp_size_t, mp_ptr));
3821
3822 /* 4*(an + 1) + 4*(bn + 1) + an */
3823 #define MPN_GCDEXT_LEHMER_ITCH(an, bn) (5*(an) + 4*(bn) + 8)
3824
3825 #ifndef HGCD_THRESHOLD
3826 #define HGCD_THRESHOLD 400
3827 #endif
3828
3829 #ifndef GCD_DC_THRESHOLD
3830 #define GCD_DC_THRESHOLD 1000
3831 #endif
3832
3833 #ifndef GCDEXT_DC_THRESHOLD
3834 #define GCDEXT_DC_THRESHOLD 600
3835 #endif
3836
3837 /* Definitions for mpn_set_str and mpn_get_str */
3838 struct powers
3839 {
3840   mp_ptr p;                     /* actual power value */
3841   mp_size_t n;                  /* # of limbs at p */
3842   mp_size_t shift;              /* weight of lowest limb, in limb base B */
3843   size_t digits_in_base;        /* number of corresponding digits */
3844   int base;
3845 };
3846 typedef struct powers powers_t;
3847 #define mpn_dc_set_str_powtab_alloc(n) ((n) + GMP_LIMB_BITS)
3848 #define mpn_dc_set_str_itch(n) ((n) + GMP_LIMB_BITS)
3849 #define mpn_dc_get_str_powtab_alloc(n) ((n) + 2 * GMP_LIMB_BITS)
3850 #define mpn_dc_get_str_itch(n) ((n) + GMP_LIMB_BITS)
3851
3852 #define   mpn_dc_set_str __MPN(dc_set_str)
3853 __GMP_DECLSPEC mp_size_t mpn_dc_set_str __GMP_PROTO ((mp_ptr, const unsigned char *, size_t, const powers_t *, mp_ptr));
3854 #define   mpn_bc_set_str __MPN(bc_set_str)
3855 __GMP_DECLSPEC mp_size_t mpn_bc_set_str __GMP_PROTO ((mp_ptr, const unsigned char *, size_t, int));
3856 #define   mpn_set_str_compute_powtab __MPN(set_str_compute_powtab)
3857 __GMP_DECLSPEC void      mpn_set_str_compute_powtab __GMP_PROTO ((powers_t *, mp_ptr, mp_size_t, int));
3858
3859
3860 /* __GMPF_BITS_TO_PREC applies a minimum 53 bits, rounds upwards to a whole
3861    limb and adds an extra limb.  __GMPF_PREC_TO_BITS drops that extra limb,
3862    hence giving back the user's size in bits rounded up.  Notice that
3863    converting prec->bits->prec gives an unchanged value.  */
3864 #define __GMPF_BITS_TO_PREC(n)                                          \
3865   ((mp_size_t) ((__GMP_MAX (53, n) + 2 * GMP_NUMB_BITS - 1) / GMP_NUMB_BITS))
3866 #define __GMPF_PREC_TO_BITS(n) \
3867   ((mp_bitcnt_t) (n) * GMP_NUMB_BITS - GMP_NUMB_BITS)
3868
3869 __GMP_DECLSPEC extern mp_size_t __gmp_default_fp_limb_precision;
3870
3871
3872 /* Set n to the number of significant digits an mpf of the given _mp_prec
3873    field, in the given base.  This is a rounded up value, designed to ensure
3874    there's enough digits to reproduce all the guaranteed part of the value.
3875
3876    There are prec many limbs, but the high might be only "1" so forget it
3877    and just count prec-1 limbs into chars.  +1 rounds that upwards, and a
3878    further +1 is because the limbs usually won't fall on digit boundaries.
3879
3880    FIXME: If base is a power of 2 and the bits per digit divides
3881    GMP_LIMB_BITS then the +2 is unnecessary.  This happens always for
3882    base==2, and in base==16 with the current 32 or 64 bit limb sizes. */
3883
3884 #define MPF_SIGNIFICANT_DIGITS(n, base, prec)                           \
3885   do {                                                                  \
3886     ASSERT (base >= 2 && base < numberof (mp_bases));                   \
3887     (n) = 2 + (size_t) ((((size_t) (prec) - 1) * GMP_NUMB_BITS)         \
3888                         * mp_bases[(base)].chars_per_bit_exactly);      \
3889   } while (0)
3890
3891
3892 /* Decimal point string, from the current C locale.  Needs <langinfo.h> for
3893    nl_langinfo and constants, preferably with _GNU_SOURCE defined to get
3894    DECIMAL_POINT from glibc, and needs <locale.h> for localeconv, each under
3895    their respective #if HAVE_FOO_H.
3896
3897    GLIBC recommends nl_langinfo because getting only one facet can be
3898    faster, apparently. */
3899
3900 /* DECIMAL_POINT seems to need _GNU_SOURCE defined to get it from glibc. */
3901 #if HAVE_NL_LANGINFO && defined (DECIMAL_POINT)
3902 #define GMP_DECIMAL_POINT  (nl_langinfo (DECIMAL_POINT))
3903 #endif
3904 /* RADIXCHAR is deprecated, still in unix98 or some such. */
3905 #if HAVE_NL_LANGINFO && defined (RADIXCHAR) && ! defined (GMP_DECIMAL_POINT)
3906 #define GMP_DECIMAL_POINT  (nl_langinfo (RADIXCHAR))
3907 #endif
3908 /* localeconv is slower since it returns all locale stuff */
3909 #if HAVE_LOCALECONV && ! defined (GMP_DECIMAL_POINT)
3910 #define GMP_DECIMAL_POINT  (localeconv()->decimal_point)
3911 #endif
3912 #if ! defined (GMP_DECIMAL_POINT)
3913 #define GMP_DECIMAL_POINT  (".")
3914 #endif
3915
3916
3917 #define DOPRNT_CONV_FIXED        1
3918 #define DOPRNT_CONV_SCIENTIFIC   2
3919 #define DOPRNT_CONV_GENERAL      3
3920
3921 #define DOPRNT_JUSTIFY_NONE      0
3922 #define DOPRNT_JUSTIFY_LEFT      1
3923 #define DOPRNT_JUSTIFY_RIGHT     2
3924 #define DOPRNT_JUSTIFY_INTERNAL  3
3925
3926 #define DOPRNT_SHOWBASE_YES      1
3927 #define DOPRNT_SHOWBASE_NO       2
3928 #define DOPRNT_SHOWBASE_NONZERO  3
3929
3930 struct doprnt_params_t {
3931   int         base;          /* negative for upper case */
3932   int         conv;          /* choices above */
3933   const char  *expfmt;       /* exponent format */
3934   int         exptimes4;     /* exponent multiply by 4 */
3935   char        fill;          /* character */
3936   int         justify;       /* choices above */
3937   int         prec;          /* prec field, or -1 for all digits */
3938   int         showbase;      /* choices above */
3939   int         showpoint;     /* if radix point always shown */
3940   int         showtrailing;  /* if trailing zeros wanted */
3941   char        sign;          /* '+', ' ', or '\0' */
3942   int         width;         /* width field */
3943 };
3944
3945 #if _GMP_H_HAVE_VA_LIST
3946
3947 __GMP_DECLSPEC typedef int (*doprnt_format_t) __GMP_PROTO ((void *, const char *, va_list));
3948 __GMP_DECLSPEC typedef int (*doprnt_memory_t) __GMP_PROTO ((void *, const char *, size_t));
3949 __GMP_DECLSPEC typedef int (*doprnt_reps_t)   __GMP_PROTO ((void *, int, int));
3950 __GMP_DECLSPEC typedef int (*doprnt_final_t)  __GMP_PROTO ((void *));
3951
3952 struct doprnt_funs_t {
3953   doprnt_format_t  format;
3954   doprnt_memory_t  memory;
3955   doprnt_reps_t    reps;
3956   doprnt_final_t   final;   /* NULL if not required */
3957 };
3958
3959 extern const struct doprnt_funs_t  __gmp_fprintf_funs;
3960 extern const struct doprnt_funs_t  __gmp_sprintf_funs;
3961 extern const struct doprnt_funs_t  __gmp_snprintf_funs;
3962 extern const struct doprnt_funs_t  __gmp_obstack_printf_funs;
3963 extern const struct doprnt_funs_t  __gmp_ostream_funs;
3964
3965 /* "buf" is a __gmp_allocate_func block of "alloc" many bytes.  The first
3966    "size" of these have been written.  "alloc > size" is maintained, so
3967    there's room to store a '\0' at the end.  "result" is where the
3968    application wants the final block pointer.  */
3969 struct gmp_asprintf_t {
3970   char    **result;
3971   char    *buf;
3972   size_t  size;
3973   size_t  alloc;
3974 };
3975
3976 #define GMP_ASPRINTF_T_INIT(d, output)                          \
3977   do {                                                          \
3978     (d).result = (output);                                      \
3979     (d).alloc = 256;                                            \
3980     (d).buf = (char *) (*__gmp_allocate_func) ((d).alloc);      \
3981     (d).size = 0;                                               \
3982   } while (0)
3983
3984 /* If a realloc is necessary, use twice the size actually required, so as to
3985    avoid repeated small reallocs.  */
3986 #define GMP_ASPRINTF_T_NEED(d, n)                                       \
3987   do {                                                                  \
3988     size_t  alloc, newsize, newalloc;                                   \
3989     ASSERT ((d)->alloc >= (d)->size + 1);                               \
3990                                                                         \
3991     alloc = (d)->alloc;                                                 \
3992     newsize = (d)->size + (n);                                          \
3993     if (alloc <= newsize)                                               \
3994       {                                                                 \
3995         newalloc = 2*newsize;                                           \
3996         (d)->alloc = newalloc;                                          \
3997         (d)->buf = __GMP_REALLOCATE_FUNC_TYPE ((d)->buf,                \
3998                                                alloc, newalloc, char);  \
3999       }                                                                 \
4000   } while (0)
4001
4002 __GMP_DECLSPEC int __gmp_asprintf_memory __GMP_PROTO ((struct gmp_asprintf_t *, const char *, size_t));
4003 __GMP_DECLSPEC int __gmp_asprintf_reps __GMP_PROTO ((struct gmp_asprintf_t *, int, int));
4004 __GMP_DECLSPEC int __gmp_asprintf_final __GMP_PROTO ((struct gmp_asprintf_t *));
4005
4006 /* buf is where to write the next output, and size is how much space is left
4007    there.  If the application passed size==0 then that's what we'll have
4008    here, and nothing at all should be written.  */
4009 struct gmp_snprintf_t {
4010   char    *buf;
4011   size_t  size;
4012 };
4013
4014 /* Add the bytes printed by the call to the total retval, or bail out on an
4015    error.  */
4016 #define DOPRNT_ACCUMULATE(call) \
4017   do {                          \
4018     int  __ret;                 \
4019     __ret = call;               \
4020     if (__ret == -1)            \
4021       goto error;               \
4022     retval += __ret;            \
4023   } while (0)
4024 #define DOPRNT_ACCUMULATE_FUN(fun, params)      \
4025   do {                                          \
4026     ASSERT ((fun) != NULL);                     \
4027     DOPRNT_ACCUMULATE ((*(fun)) params);        \
4028   } while (0)
4029
4030 #define DOPRNT_FORMAT(fmt, ap)                          \
4031   DOPRNT_ACCUMULATE_FUN (funs->format, (data, fmt, ap))
4032 #define DOPRNT_MEMORY(ptr, len)                                 \
4033   DOPRNT_ACCUMULATE_FUN (funs->memory, (data, ptr, len))
4034 #define DOPRNT_REPS(c, n)                               \
4035   DOPRNT_ACCUMULATE_FUN (funs->reps, (data, c, n))
4036
4037 #define DOPRNT_STRING(str)      DOPRNT_MEMORY (str, strlen (str))
4038
4039 #define DOPRNT_REPS_MAYBE(c, n) \
4040   do {                          \
4041     if ((n) != 0)               \
4042       DOPRNT_REPS (c, n);       \
4043   } while (0)
4044 #define DOPRNT_MEMORY_MAYBE(ptr, len)   \
4045   do {                                  \
4046     if ((len) != 0)                     \
4047       DOPRNT_MEMORY (ptr, len);         \
4048   } while (0)
4049
4050 __GMP_DECLSPEC int __gmp_doprnt __GMP_PROTO ((const struct doprnt_funs_t *, void *, const char *, va_list));
4051 __GMP_DECLSPEC int __gmp_doprnt_integer __GMP_PROTO ((const struct doprnt_funs_t *, void *, const struct doprnt_params_t *, const char *));
4052
4053 #define __gmp_doprnt_mpf __gmp_doprnt_mpf2
4054 __GMP_DECLSPEC int __gmp_doprnt_mpf __GMP_PROTO ((const struct doprnt_funs_t *, void *, const struct doprnt_params_t *, const char *, mpf_srcptr));
4055
4056 __GMP_DECLSPEC int __gmp_replacement_vsnprintf __GMP_PROTO ((char *, size_t, const char *, va_list));
4057 #endif /* _GMP_H_HAVE_VA_LIST */
4058
4059
4060 typedef int (*gmp_doscan_scan_t)  __GMP_PROTO ((void *, const char *, ...));
4061 typedef void *(*gmp_doscan_step_t) __GMP_PROTO ((void *, int));
4062 typedef int (*gmp_doscan_get_t)   __GMP_PROTO ((void *));
4063 typedef int (*gmp_doscan_unget_t) __GMP_PROTO ((int, void *));
4064
4065 struct gmp_doscan_funs_t {
4066   gmp_doscan_scan_t   scan;
4067   gmp_doscan_step_t   step;
4068   gmp_doscan_get_t    get;
4069   gmp_doscan_unget_t  unget;
4070 };
4071 extern const struct gmp_doscan_funs_t  __gmp_fscanf_funs;
4072 extern const struct gmp_doscan_funs_t  __gmp_sscanf_funs;
4073
4074 #if _GMP_H_HAVE_VA_LIST
4075 __GMP_DECLSPEC int __gmp_doscan __GMP_PROTO ((const struct gmp_doscan_funs_t *, void *, const char *, va_list));
4076 #endif
4077
4078
4079 /* For testing and debugging.  */
4080 #define MPZ_CHECK_FORMAT(z)                                     \
4081   do {                                                          \
4082     ASSERT_ALWAYS (SIZ(z) == 0 || PTR(z)[ABSIZ(z) - 1] != 0);   \
4083     ASSERT_ALWAYS (ALLOC(z) >= ABSIZ(z));                       \
4084     ASSERT_ALWAYS_MPN (PTR(z), ABSIZ(z));                       \
4085   } while (0)
4086
4087 #define MPQ_CHECK_FORMAT(q)                             \
4088   do {                                                  \
4089     MPZ_CHECK_FORMAT (mpq_numref (q));                  \
4090     MPZ_CHECK_FORMAT (mpq_denref (q));                  \
4091     ASSERT_ALWAYS (SIZ(mpq_denref(q)) >= 1);            \
4092                                                         \
4093     if (SIZ(mpq_numref(q)) == 0)                        \
4094       {                                                 \
4095         /* should have zero as 0/1 */                   \
4096         ASSERT_ALWAYS (SIZ(mpq_denref(q)) == 1          \
4097                        && PTR(mpq_denref(q))[0] == 1);  \
4098       }                                                 \
4099     else                                                \
4100       {                                                 \
4101         /* should have no common factors */             \
4102         mpz_t  g;                                       \
4103         mpz_init (g);                                   \
4104         mpz_gcd (g, mpq_numref(q), mpq_denref(q));      \
4105         ASSERT_ALWAYS (mpz_cmp_ui (g, 1) == 0);         \
4106         mpz_clear (g);                                  \
4107       }                                                 \
4108   } while (0)
4109
4110 #define MPF_CHECK_FORMAT(f)                             \
4111   do {                                                  \
4112     ASSERT_ALWAYS (PREC(f) >= __GMPF_BITS_TO_PREC(53)); \
4113     ASSERT_ALWAYS (ABSIZ(f) <= PREC(f)+1);              \
4114     if (SIZ(f) == 0)                                    \
4115       ASSERT_ALWAYS (EXP(f) == 0);                      \
4116     if (SIZ(f) != 0)                                    \
4117       ASSERT_ALWAYS (PTR(f)[ABSIZ(f) - 1] != 0);        \
4118   } while (0)
4119
4120
4121 #define MPZ_PROVOKE_REALLOC(z)                                  \
4122   do { ALLOC(z) = ABSIZ(z); } while (0)
4123
4124
4125 /* Enhancement: The "mod" and "gcd_1" functions below could have
4126    __GMP_ATTRIBUTE_PURE, but currently (gcc 3.3) that's not supported on
4127    function pointers, only actual functions.  It probably doesn't make much
4128    difference to the gmp code, since hopefully we arrange calls so there's
4129    no great need for the compiler to move things around.  */
4130
4131 #if WANT_FAT_BINARY && (HAVE_HOST_CPU_FAMILY_x86 || HAVE_HOST_CPU_FAMILY_x86_64)
4132 /* NOTE: The function pointers in this struct are also in CPUVEC_FUNCS_LIST
4133    in mpn/x86/x86-defs.m4.  Be sure to update that when changing here.  */
4134 struct cpuvec_t {
4135   DECL_add_n           ((*add_n));
4136   DECL_addmul_1        ((*addmul_1));
4137   DECL_copyd           ((*copyd));
4138   DECL_copyi           ((*copyi));
4139   DECL_divexact_1      ((*divexact_1));
4140   DECL_divexact_by3c   ((*divexact_by3c));
4141   DECL_divrem_1        ((*divrem_1));
4142   DECL_gcd_1           ((*gcd_1));
4143   DECL_lshift          ((*lshift));
4144   DECL_mod_1           ((*mod_1));
4145   DECL_mod_34lsub1     ((*mod_34lsub1));
4146   DECL_modexact_1c_odd ((*modexact_1c_odd));
4147   DECL_mul_1           ((*mul_1));
4148   DECL_mul_basecase    ((*mul_basecase));
4149   DECL_preinv_divrem_1 ((*preinv_divrem_1));
4150   DECL_preinv_mod_1    ((*preinv_mod_1));
4151   DECL_rshift          ((*rshift));
4152   DECL_sqr_basecase    ((*sqr_basecase));
4153   DECL_sub_n           ((*sub_n));
4154   DECL_submul_1        ((*submul_1));
4155   int                  initialized;
4156   mp_size_t            mul_toom22_threshold;
4157   mp_size_t            mul_toom33_threshold;
4158   mp_size_t            sqr_toom2_threshold;
4159   mp_size_t            sqr_toom3_threshold;
4160 };
4161 __GMP_DECLSPEC extern struct cpuvec_t __gmpn_cpuvec;
4162 #endif /* x86 fat binary */
4163
4164 __GMP_DECLSPEC void __gmpn_cpuvec_init __GMP_PROTO ((void));
4165
4166 /* Get a threshold "field" from __gmpn_cpuvec, running __gmpn_cpuvec_init()
4167    if that hasn't yet been done (to establish the right values).  */
4168 #define CPUVEC_THRESHOLD(field)                                               \
4169   ((LIKELY (__gmpn_cpuvec.initialized) ? 0 : (__gmpn_cpuvec_init (), 0)),     \
4170    __gmpn_cpuvec.field)
4171
4172
4173 #if HAVE_NATIVE_mpn_add_nc
4174 #define mpn_add_nc __MPN(add_nc)
4175 __GMP_DECLSPEC mp_limb_t mpn_add_nc __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, mp_limb_t));
4176 #else
4177 static inline
4178 mp_limb_t
4179 mpn_add_nc (mp_ptr rp, mp_srcptr up, mp_srcptr vp, mp_size_t n, mp_limb_t ci)
4180 {
4181   mp_limb_t co;
4182   co = mpn_add_n (rp, up, vp, n);
4183   co += mpn_add_1 (rp, rp, n, ci);
4184   return co;
4185 }
4186 #endif
4187
4188 #if HAVE_NATIVE_mpn_sub_nc
4189 #define mpn_sub_nc __MPN(sub_nc)
4190 __GMP_DECLSPEC mp_limb_t mpn_sub_nc __GMP_PROTO ((mp_ptr, mp_srcptr, mp_srcptr, mp_size_t, mp_limb_t));
4191 #else
4192 static inline mp_limb_t
4193 mpn_sub_nc (mp_ptr rp, mp_srcptr up, mp_srcptr vp, mp_size_t n, mp_limb_t ci)
4194 {
4195   mp_limb_t co;
4196   co = mpn_sub_n (rp, up, vp, n);
4197   co += mpn_sub_1 (rp, rp, n, ci);
4198   return co;
4199 }
4200 #endif
4201
4202 static inline int
4203 mpn_zero_p (mp_srcptr ap, mp_size_t n)
4204 {
4205   mp_size_t i;
4206   for (i = n - 1; i >= 0; i--)
4207     {
4208       if (ap[i] != 0)
4209         return 0;
4210     }
4211   return 1;
4212 }
4213
4214 #if TUNE_PROGRAM_BUILD
4215 /* Some extras wanted when recompiling some .c files for use by the tune
4216    program.  Not part of a normal build.
4217
4218    It's necessary to keep these thresholds as #defines (just to an
4219    identically named variable), since various defaults are established based
4220    on #ifdef in the .c files.  For some this is not so (the defaults are
4221    instead established above), but all are done this way for consistency. */
4222
4223 #undef  MUL_TOOM22_THRESHOLD
4224 #define MUL_TOOM22_THRESHOLD            mul_toom22_threshold
4225 extern mp_size_t                        mul_toom22_threshold;
4226
4227 #undef  MUL_TOOM33_THRESHOLD
4228 #define MUL_TOOM33_THRESHOLD            mul_toom33_threshold
4229 extern mp_size_t                        mul_toom33_threshold;
4230
4231 #undef  MUL_TOOM44_THRESHOLD
4232 #define MUL_TOOM44_THRESHOLD            mul_toom44_threshold
4233 extern mp_size_t                        mul_toom44_threshold;
4234
4235 #undef  MUL_TOOM6H_THRESHOLD
4236 #define MUL_TOOM6H_THRESHOLD            mul_toom6h_threshold
4237 extern mp_size_t                        mul_toom6h_threshold;
4238
4239 #undef  MUL_TOOM8H_THRESHOLD
4240 #define MUL_TOOM8H_THRESHOLD            mul_toom8h_threshold
4241 extern mp_size_t                        mul_toom8h_threshold;
4242
4243 #undef  MUL_TOOM32_TO_TOOM43_THRESHOLD
4244 #define MUL_TOOM32_TO_TOOM43_THRESHOLD  mul_toom32_to_toom43_threshold
4245 extern mp_size_t                        mul_toom32_to_toom43_threshold;
4246
4247 #undef  MUL_TOOM32_TO_TOOM53_THRESHOLD
4248 #define MUL_TOOM32_TO_TOOM53_THRESHOLD  mul_toom32_to_toom53_threshold
4249 extern mp_size_t                        mul_toom32_to_toom53_threshold;
4250
4251 #undef  MUL_TOOM42_TO_TOOM53_THRESHOLD
4252 #define MUL_TOOM42_TO_TOOM53_THRESHOLD  mul_toom42_to_toom53_threshold
4253 extern mp_size_t                        mul_toom42_to_toom53_threshold;
4254
4255 #undef  MUL_TOOM42_TO_TOOM63_THRESHOLD
4256 #define MUL_TOOM42_TO_TOOM63_THRESHOLD  mul_toom42_to_toom63_threshold
4257 extern mp_size_t                        mul_toom42_to_toom63_threshold;
4258
4259 #undef  MUL_FFT_THRESHOLD
4260 #define MUL_FFT_THRESHOLD               mul_fft_threshold
4261 extern mp_size_t                        mul_fft_threshold;
4262
4263 #undef  MUL_FFT_MODF_THRESHOLD
4264 #define MUL_FFT_MODF_THRESHOLD          mul_fft_modf_threshold
4265 extern mp_size_t                        mul_fft_modf_threshold;
4266
4267 #undef  MUL_FFT_TABLE
4268 #define MUL_FFT_TABLE                   { 0 }
4269
4270 #undef  MUL_FFT_TABLE3
4271 #define MUL_FFT_TABLE3                  { {0,0} }
4272
4273 /* A native mpn_sqr_basecase is not tuned and SQR_BASECASE_THRESHOLD should
4274    remain as zero (always use it). */
4275 #if ! HAVE_NATIVE_mpn_sqr_basecase
4276 #undef  SQR_BASECASE_THRESHOLD
4277 #define SQR_BASECASE_THRESHOLD          sqr_basecase_threshold
4278 extern mp_size_t                        sqr_basecase_threshold;
4279 #endif
4280
4281 #if TUNE_PROGRAM_BUILD_SQR
4282 #undef  SQR_TOOM2_THRESHOLD
4283 #define SQR_TOOM2_THRESHOLD             SQR_TOOM2_MAX_GENERIC
4284 #else
4285 #undef  SQR_TOOM2_THRESHOLD
4286 #define SQR_TOOM2_THRESHOLD             sqr_toom2_threshold
4287 extern mp_size_t                        sqr_toom2_threshold;
4288 #endif
4289
4290 #undef  SQR_TOOM3_THRESHOLD
4291 #define SQR_TOOM3_THRESHOLD             sqr_toom3_threshold
4292 extern mp_size_t                        sqr_toom3_threshold;
4293
4294 #undef  SQR_TOOM4_THRESHOLD
4295 #define SQR_TOOM4_THRESHOLD             sqr_toom4_threshold
4296 extern mp_size_t                        sqr_toom4_threshold;
4297
4298 #undef  SQR_TOOM6_THRESHOLD
4299 #define SQR_TOOM6_THRESHOLD             sqr_toom6_threshold
4300 extern mp_size_t                        sqr_toom6_threshold;
4301
4302 #undef  SQR_TOOM8_THRESHOLD
4303 #define SQR_TOOM8_THRESHOLD             sqr_toom8_threshold
4304 extern mp_size_t                        sqr_toom8_threshold;
4305
4306 #undef  SQR_FFT_THRESHOLD
4307 #define SQR_FFT_THRESHOLD               sqr_fft_threshold
4308 extern mp_size_t                        sqr_fft_threshold;
4309
4310 #undef  SQR_FFT_MODF_THRESHOLD
4311 #define SQR_FFT_MODF_THRESHOLD          sqr_fft_modf_threshold
4312 extern mp_size_t                        sqr_fft_modf_threshold;
4313
4314 #undef  SQR_FFT_TABLE
4315 #define SQR_FFT_TABLE                   { 0 }
4316
4317 #undef  SQR_FFT_TABLE3
4318 #define SQR_FFT_TABLE3                  { {0,0} }
4319
4320 #undef  MULLO_BASECASE_THRESHOLD
4321 #define MULLO_BASECASE_THRESHOLD        mullo_basecase_threshold
4322 extern mp_size_t                        mullo_basecase_threshold;
4323
4324 #undef  MULLO_DC_THRESHOLD
4325 #define MULLO_DC_THRESHOLD              mullo_dc_threshold
4326 extern mp_size_t                        mullo_dc_threshold;
4327
4328 #undef  MULLO_MUL_N_THRESHOLD
4329 #define MULLO_MUL_N_THRESHOLD           mullo_mul_n_threshold
4330 extern mp_size_t                        mullo_mul_n_threshold;
4331
4332 #undef  DC_DIV_QR_THRESHOLD
4333 #define DC_DIV_QR_THRESHOLD             dc_div_qr_threshold
4334 extern mp_size_t                        dc_div_qr_threshold;
4335
4336 #undef  DC_DIVAPPR_Q_THRESHOLD
4337 #define DC_DIVAPPR_Q_THRESHOLD          dc_divappr_q_threshold
4338 extern mp_size_t                        dc_divappr_q_threshold;
4339
4340 #undef  DC_BDIV_Q_THRESHOLD
4341 #define DC_BDIV_Q_THRESHOLD             dc_bdiv_q_threshold
4342 extern mp_size_t                        dc_bdiv_q_threshold;
4343
4344 #undef  DC_BDIV_QR_THRESHOLD
4345 #define DC_BDIV_QR_THRESHOLD            dc_bdiv_qr_threshold
4346 extern mp_size_t                        dc_bdiv_qr_threshold;
4347
4348 #undef  MU_DIV_QR_THRESHOLD
4349 #define MU_DIV_QR_THRESHOLD             mu_div_qr_threshold
4350 extern mp_size_t                        mu_div_qr_threshold;
4351
4352 #undef  MU_DIVAPPR_Q_THRESHOLD
4353 #define MU_DIVAPPR_Q_THRESHOLD          mu_divappr_q_threshold
4354 extern mp_size_t                        mu_divappr_q_threshold;
4355
4356 #undef  MUPI_DIV_QR_THRESHOLD
4357 #define MUPI_DIV_QR_THRESHOLD           mupi_div_qr_threshold
4358 extern mp_size_t                        mupi_div_qr_threshold;
4359
4360 #undef  MU_BDIV_QR_THRESHOLD
4361 #define MU_BDIV_QR_THRESHOLD            mu_bdiv_qr_threshold
4362 extern mp_size_t                        mu_bdiv_qr_threshold;
4363
4364 #undef  MU_BDIV_Q_THRESHOLD
4365 #define MU_BDIV_Q_THRESHOLD             mu_bdiv_q_threshold
4366 extern mp_size_t                        mu_bdiv_q_threshold;
4367
4368 #undef  INV_MULMOD_BNM1_THRESHOLD
4369 #define INV_MULMOD_BNM1_THRESHOLD       inv_mulmod_bnm1_threshold
4370 extern mp_size_t                        inv_mulmod_bnm1_threshold;
4371
4372 #undef  INV_NEWTON_THRESHOLD
4373 #define INV_NEWTON_THRESHOLD            inv_newton_threshold
4374 extern mp_size_t                        inv_newton_threshold;
4375
4376 #undef  INV_APPR_THRESHOLD
4377 #define INV_APPR_THRESHOLD              inv_appr_threshold
4378 extern mp_size_t                        inv_appr_threshold;
4379
4380 #undef  BINV_NEWTON_THRESHOLD
4381 #define BINV_NEWTON_THRESHOLD           binv_newton_threshold
4382 extern mp_size_t                        binv_newton_threshold;
4383
4384 #undef  REDC_1_TO_REDC_2_THRESHOLD
4385 #define REDC_1_TO_REDC_2_THRESHOLD      redc_1_to_redc_2_threshold
4386 extern mp_size_t                        redc_1_to_redc_2_threshold;
4387
4388 #undef  REDC_2_TO_REDC_N_THRESHOLD
4389 #define REDC_2_TO_REDC_N_THRESHOLD      redc_2_to_redc_n_threshold
4390 extern mp_size_t                        redc_2_to_redc_n_threshold;
4391
4392 #undef  REDC_1_TO_REDC_N_THRESHOLD
4393 #define REDC_1_TO_REDC_N_THRESHOLD      redc_1_to_redc_n_threshold
4394 extern mp_size_t                        redc_1_to_redc_n_threshold;
4395
4396 #undef  MATRIX22_STRASSEN_THRESHOLD
4397 #define MATRIX22_STRASSEN_THRESHOLD     matrix22_strassen_threshold
4398 extern mp_size_t                        matrix22_strassen_threshold;
4399
4400 #undef  HGCD_THRESHOLD
4401 #define HGCD_THRESHOLD                  hgcd_threshold
4402 extern mp_size_t                        hgcd_threshold;
4403
4404 #undef  GCD_DC_THRESHOLD
4405 #define GCD_DC_THRESHOLD                gcd_dc_threshold
4406 extern mp_size_t                        gcd_dc_threshold;
4407
4408 #undef  GCDEXT_DC_THRESHOLD
4409 #define GCDEXT_DC_THRESHOLD             gcdext_dc_threshold
4410 extern mp_size_t                        gcdext_dc_threshold;
4411
4412 #undef  DIVREM_1_NORM_THRESHOLD
4413 #define DIVREM_1_NORM_THRESHOLD         divrem_1_norm_threshold
4414 extern mp_size_t                        divrem_1_norm_threshold;
4415
4416 #undef  DIVREM_1_UNNORM_THRESHOLD
4417 #define DIVREM_1_UNNORM_THRESHOLD       divrem_1_unnorm_threshold
4418 extern mp_size_t                        divrem_1_unnorm_threshold;
4419
4420 #undef  MOD_1_NORM_THRESHOLD
4421 #define MOD_1_NORM_THRESHOLD            mod_1_norm_threshold
4422 extern mp_size_t                        mod_1_norm_threshold;
4423
4424 #undef  MOD_1_UNNORM_THRESHOLD
4425 #define MOD_1_UNNORM_THRESHOLD          mod_1_unnorm_threshold
4426 extern mp_size_t                        mod_1_unnorm_threshold;
4427
4428 #undef  MOD_1N_TO_MOD_1_1_THRESHOLD
4429 #define MOD_1N_TO_MOD_1_1_THRESHOLD     mod_1n_to_mod_1_1_threshold
4430 extern mp_size_t                        mod_1n_to_mod_1_1_threshold;
4431
4432 #undef  MOD_1U_TO_MOD_1_1_THRESHOLD
4433 #define MOD_1U_TO_MOD_1_1_THRESHOLD     mod_1u_to_mod_1_1_threshold
4434 extern mp_size_t                        mod_1u_to_mod_1_1_threshold;
4435
4436 #undef  MOD_1_1_TO_MOD_1_2_THRESHOLD
4437 #define MOD_1_1_TO_MOD_1_2_THRESHOLD    mod_1_1_to_mod_1_2_threshold
4438 extern mp_size_t                        mod_1_1_to_mod_1_2_threshold;
4439
4440 #undef  MOD_1_2_TO_MOD_1_4_THRESHOLD
4441 #define MOD_1_2_TO_MOD_1_4_THRESHOLD    mod_1_2_to_mod_1_4_threshold
4442 extern mp_size_t                        mod_1_2_to_mod_1_4_threshold;
4443
4444 #undef  PREINV_MOD_1_TO_MOD_1_THRESHOLD
4445 #define PREINV_MOD_1_TO_MOD_1_THRESHOLD preinv_mod_1_to_mod_1_threshold
4446 extern mp_size_t                        preinv_mod_1_to_mod_1_threshold;
4447
4448 #if ! UDIV_PREINV_ALWAYS
4449 #undef  DIVREM_2_THRESHOLD
4450 #define DIVREM_2_THRESHOLD              divrem_2_threshold
4451 extern mp_size_t                        divrem_2_threshold;
4452 #endif
4453
4454 #undef  MULMOD_BNM1_THRESHOLD
4455 #define MULMOD_BNM1_THRESHOLD           mulmod_bnm1_threshold
4456 extern mp_size_t                        mulmod_bnm1_threshold;
4457
4458 #undef  SQRMOD_BNM1_THRESHOLD
4459 #define SQRMOD_BNM1_THRESHOLD           sqrmod_bnm1_threshold
4460 extern mp_size_t                        sqrmod_bnm1_threshold;
4461
4462 #undef  GET_STR_DC_THRESHOLD
4463 #define GET_STR_DC_THRESHOLD            get_str_dc_threshold
4464 extern mp_size_t                        get_str_dc_threshold;
4465
4466 #undef  GET_STR_PRECOMPUTE_THRESHOLD
4467 #define GET_STR_PRECOMPUTE_THRESHOLD    get_str_precompute_threshold
4468 extern mp_size_t                        get_str_precompute_threshold;
4469
4470 #undef  SET_STR_DC_THRESHOLD
4471 #define SET_STR_DC_THRESHOLD            set_str_dc_threshold
4472 extern mp_size_t                        set_str_dc_threshold;
4473
4474 #undef  SET_STR_PRECOMPUTE_THRESHOLD
4475 #define SET_STR_PRECOMPUTE_THRESHOLD    set_str_precompute_threshold
4476 extern mp_size_t                        set_str_precompute_threshold;
4477
4478 #undef  FFT_TABLE_ATTRS
4479 #define FFT_TABLE_ATTRS
4480 extern mp_size_t  mpn_fft_table[2][MPN_FFT_TABLE_SIZE];
4481 #define FFT_TABLE3_SIZE 2000    /* generous space for tuning */
4482 extern struct fft_table_nk mpn_fft_table3[2][FFT_TABLE3_SIZE];
4483
4484 /* Sizes the tune program tests up to, used in a couple of recompilations. */
4485 #undef MUL_TOOM22_THRESHOLD_LIMIT
4486 #undef MUL_TOOM33_THRESHOLD_LIMIT
4487 #undef MULLO_BASECASE_THRESHOLD_LIMIT
4488 #undef SQR_TOOM3_THRESHOLD_LIMIT
4489 #define SQR_TOOM2_MAX_GENERIC           200
4490 #define MUL_TOOM22_THRESHOLD_LIMIT      700
4491 #define MUL_TOOM33_THRESHOLD_LIMIT      700
4492 #define SQR_TOOM3_THRESHOLD_LIMIT       400
4493 #define MUL_TOOM44_THRESHOLD_LIMIT     1000
4494 #define SQR_TOOM4_THRESHOLD_LIMIT      1000
4495 #define MUL_TOOM6H_THRESHOLD_LIMIT     1100
4496 #define SQR_TOOM6_THRESHOLD_LIMIT      1100
4497 #define MUL_TOOM8H_THRESHOLD_LIMIT     1200
4498 #define SQR_TOOM8_THRESHOLD_LIMIT      1200
4499 #define MULLO_BASECASE_THRESHOLD_LIMIT  200
4500 #define GET_STR_THRESHOLD_LIMIT         150
4501
4502 #endif /* TUNE_PROGRAM_BUILD */
4503
4504 #if defined (__cplusplus)
4505 }
4506 #endif
4507
4508 /* FIXME: Make these itch functions less conservative.  Also consider making
4509    them dependent on just 'an', and compute the allocation directly from 'an'
4510    instead of via n.  */
4511
4512 /* toom22/toom2: Scratch need is 2*(an + k), k is the recursion depth.
4513    k is ths smallest k such that
4514      ceil(an/2^k) < MUL_TOOM22_THRESHOLD.
4515    which implies that
4516      k = bitsize of floor ((an-1)/(MUL_TOOM22_THRESHOLD-1))
4517        = 1 + floor (log_2 (floor ((an-1)/(MUL_TOOM22_THRESHOLD-1))))
4518 */
4519 #define mpn_toom22_mul_itch(an, bn) \
4520   (2 * ((an) + GMP_NUMB_BITS))
4521 #define mpn_toom2_sqr_itch(an) \
4522   (2 * ((an) + GMP_NUMB_BITS))
4523
4524 /* toom33/toom3: Scratch need is 5an/2 + 10k, k is the recursion depth.
4525    We use 3an + C, so that we can use a smaller constant.
4526  */
4527 #define mpn_toom33_mul_itch(an, bn) \
4528   (3 * (an) + GMP_NUMB_BITS)
4529 #define mpn_toom3_sqr_itch(an) \
4530   (3 * (an) + GMP_NUMB_BITS)
4531
4532 /* toom33/toom3: Scratch need is 8an/3 + 13k, k is the recursion depth.
4533    We use 3an + C, so that we can use a smaller constant.
4534  */
4535 #define mpn_toom44_mul_itch(an, bn) \
4536   (3 * (an) + GMP_NUMB_BITS)
4537 #define mpn_toom4_sqr_itch(an) \
4538   (3 * (an) + GMP_NUMB_BITS)
4539
4540 #define mpn_toom6_sqr_itch(n)                                           \
4541 ( ((n) - SQR_TOOM6_THRESHOLD)*2 +                                       \
4542    MAX(SQR_TOOM6_THRESHOLD*2 + GMP_NUMB_BITS*6,                 \
4543        mpn_toom4_sqr_itch(SQR_TOOM6_THRESHOLD)) )
4544
4545 #define mpn_toom6_mul_n_itch(n)                                         \
4546 ( ((n) - MUL_TOOM6H_THRESHOLD)*2 +                                      \
4547    MAX(MUL_TOOM6H_THRESHOLD*2 + GMP_NUMB_BITS*6,                        \
4548        mpn_toom44_mul_itch(MUL_TOOM6H_THRESHOLD,MUL_TOOM6H_THRESHOLD)) )
4549
4550 static inline mp_size_t
4551 mpn_toom6h_mul_itch (mp_size_t an, mp_size_t bn) {
4552   mp_size_t estimatedN;
4553   estimatedN = (an + bn) / (size_t) 10 + 1;
4554   return mpn_toom6_mul_n_itch (estimatedN * 6);
4555 }
4556
4557 #define mpn_toom8_sqr_itch(n)                                           \
4558 ( (((n)*15)>>3) - ((SQR_TOOM8_THRESHOLD*15)>>3) +                       \
4559    MAX(((SQR_TOOM8_THRESHOLD*15)>>3) + GMP_NUMB_BITS*6,         \
4560        mpn_toom6_sqr_itch(SQR_TOOM8_THRESHOLD)) )
4561
4562 #define mpn_toom8_mul_n_itch(n)                                         \
4563 ( (((n)*15)>>3) - ((MUL_TOOM8H_THRESHOLD*15)>>3) +                      \
4564    MAX(((MUL_TOOM8H_THRESHOLD*15)>>3) + GMP_NUMB_BITS*6,                \
4565        mpn_toom6_mul_n_itch(MUL_TOOM8H_THRESHOLD)) )
4566
4567 static inline mp_size_t
4568 mpn_toom8h_mul_itch (mp_size_t an, mp_size_t bn) {
4569   mp_size_t estimatedN;
4570   estimatedN = (an + bn) / (size_t) 14 + 1;
4571   return mpn_toom8_mul_n_itch (estimatedN * 8);
4572 }
4573
4574 static inline mp_size_t
4575 mpn_toom32_mul_itch (mp_size_t an, mp_size_t bn)
4576 {
4577   mp_size_t n = 1 + (2 * an >= 3 * bn ? (an - 1) / (size_t) 3 : (bn - 1) >> 1);
4578   mp_size_t itch = 2 * n + 1;
4579
4580   return itch;
4581 }
4582
4583 static inline mp_size_t
4584 mpn_toom42_mul_itch (mp_size_t an, mp_size_t bn)
4585 {
4586   mp_size_t n = an >= 2 * bn ? (an + 3) >> 2 : (bn + 1) >> 1;
4587   return 6 * n + 3;
4588 }
4589
4590 static inline mp_size_t
4591 mpn_toom43_mul_itch (mp_size_t an, mp_size_t bn)
4592 {
4593   mp_size_t n = 1 + (3 * an >= 4 * bn ? (an - 1) >> 2 : (bn - 1) / (size_t) 3);
4594
4595   return 6*n + 4;
4596 }
4597
4598 static inline mp_size_t
4599 mpn_toom52_mul_itch (mp_size_t an, mp_size_t bn)
4600 {
4601   mp_size_t n = 1 + (2 * an >= 5 * bn ? (an - 1) / (size_t) 5 : (bn - 1) >> 1);
4602   return 6*n + 4;
4603 }
4604
4605 static inline mp_size_t
4606 mpn_toom53_mul_itch (mp_size_t an, mp_size_t bn)
4607 {
4608   mp_size_t n = 1 + (3 * an >= 5 * bn ? (an - 1) / (size_t) 5 : (bn - 1) / (size_t) 3);
4609   return 10 * n + 10;
4610 }
4611
4612 static inline mp_size_t
4613 mpn_toom62_mul_itch (mp_size_t an, mp_size_t bn)
4614 {
4615   mp_size_t n = 1 + (an >= 3 * bn ? (an - 1) / (size_t) 6 : (bn - 1) >> 1);
4616   return 10 * n + 10;
4617 }
4618
4619 static inline mp_size_t
4620 mpn_toom63_mul_itch (mp_size_t an, mp_size_t bn)
4621 {
4622   mp_size_t n = 1 + (an >= 2 * bn ? (an - 1) / (size_t) 6 : (bn - 1) / (size_t) 3);
4623   return 9 * n + 3;
4624 }
4625
4626 #if 0
4627 #define mpn_fft_mul mpn_mul_fft_full
4628 #else
4629 #define mpn_fft_mul mpn_nussbaumer_mul
4630 #endif
4631
4632 #ifdef __cplusplus
4633
4634 /* A little helper for a null-terminated __gmp_allocate_func string.
4635    The destructor ensures it's freed even if an exception is thrown.
4636    The len field is needed by the destructor, and can be used by anyone else
4637    to avoid a second strlen pass over the data.
4638
4639    Since our input is a C string, using strlen is correct.  Perhaps it'd be
4640    more C++-ish style to use std::char_traits<char>::length, but char_traits
4641    isn't available in gcc 2.95.4.  */
4642
4643 class gmp_allocated_string {
4644  public:
4645   char *str;
4646   size_t len;
4647   gmp_allocated_string(char *arg)
4648   {
4649     str = arg;
4650     len = std::strlen (str);
4651   }
4652   ~gmp_allocated_string()
4653   {
4654     (*__gmp_free_func) (str, len+1);
4655   }
4656 };
4657
4658 std::istream &__gmpz_operator_in_nowhite (std::istream &, mpz_ptr, char);
4659 int __gmp_istream_set_base (std::istream &, char &, bool &, bool &);
4660 void __gmp_istream_set_digits (std::string &, std::istream &, char &, bool &, int);
4661 void __gmp_doprnt_params_from_ios (struct doprnt_params_t *p, std::ios &o);
4662 std::ostream& __gmp_doprnt_integer_ostream (std::ostream &o, struct doprnt_params_t *p, char *s);
4663 extern const struct doprnt_funs_t  __gmp_asprintf_funs_noformat;
4664
4665 #endif /* __cplusplus */
4666
4667 #endif /* __GMP_IMPL_H__ */