/* mini-gmp, a minimalistic implementation of a GNU GMP subset.
Contributed to the GNU project by Niels Möller
Copyright 1991-1997, 1999-2019 Free Software Foundation, Inc.
This file is part of the GNU MP Library.
The GNU MP Library is free software; you can redistribute it and/or modify
it under the terms of the GNU Lesser General Public License as published by
the Free Software Foundation; either version 3 of the License, or (at your
option) any later version.
The GNU MP Library is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
License for more details.
You should have received a copy of the GNU Lesser General Public License
along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */
/* NOTE: All functions in this file which are not declared in
mini-gmp.h are internal, and are not intended to be compatible
neither with GMP nor with future versions of mini-gmp. */
/* Much of the material copied from GMP files, including: gmp-impl.h,
longlong.h, mpn/generic/add_n.c, mpn/generic/addmul_1.c,
mpn/generic/lshift.c, mpn/generic/mul_1.c,
mpn/generic/mul_basecase.c, mpn/generic/rshift.c,
mpn/generic/sbpi1_div_qr.c, mpn/generic/sub_n.c,
mpn/generic/submul_1.c. */
#include <assert.h>
#include <ctype.h>
#include <limits.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "mini-gmp.h"
#if !defined(MINI_GMP_DONT_USE_FLOAT_H)
#include <float.h>
#endif
/* Macros */
#define GMP_LIMB_BITS (sizeof(mp_limb_t) * CHAR_BIT)
#define GMP_LIMB_MAX ((mp_limb_t) ~ (mp_limb_t) 0)
#define GMP_LIMB_HIGHBIT ((mp_limb_t) 1 << (GMP_LIMB_BITS - 1))
#define GMP_HLIMB_BIT ((mp_limb_t) 1 << (GMP_LIMB_BITS / 2))
#define GMP_LLIMB_MASK (GMP_HLIMB_BIT - 1)
#define GMP_ULONG_BITS (sizeof(unsigned long) * CHAR_BIT)
#define GMP_ULONG_HIGHBIT ((unsigned long) 1 << (GMP_ULONG_BITS - 1))
#define GMP_ABS(x) ((x) >= 0 ? (x) : -(x))
#define GMP_NEG_CAST(T,x) (-((T)((x) + 1) - 1))
#define GMP_MIN(a, b) ((a) < (b) ? (a) : (b))
#define GMP_MAX(a, b) ((a) > (b) ? (a) : (b))
#define GMP_CMP(a,b) (((a) > (b)) - ((a) < (b)))
#if defined(DBL_MANT_DIG) && FLT_RADIX == 2
#define GMP_DBL_MANT_BITS DBL_MANT_DIG
#else
#define GMP_DBL_MANT_BITS (53)
#endif
/* Return non-zero if xp,xsize and yp,ysize overlap.
If xp+xsize<=yp there's no overlap, or if yp+ysize<=xp there's no
overlap. If both these are false, there's an overlap. */
#define GMP_MPN_OVERLAP_P(xp, xsize, yp, ysize) \
((xp) + (xsize) > (yp) && (yp) + (ysize) > (xp))
#define gmp_assert_nocarry(x) do { \
mp_limb_t __cy = (x); \
assert (__cy == 0); \
} while (0)
#define gmp_clz(count, x) do { \
mp_limb_t __clz_x = (x); \
unsigned __clz_c = 0; \
int LOCAL_SHIFT_BITS = 8; \
if (GMP_LIMB_BITS > LOCAL_SHIFT_BITS) \
for (; \
(__clz_x & ((mp_limb_t) 0xff << (GMP_LIMB_BITS - 8))) == 0; \
__clz_c += 8) \
{ __clz_x <<= LOCAL_SHIFT_BITS; } \
for (; (__clz_x & GMP_LIMB_HIGHBIT) == 0; __clz_c++) \
__clz_x <<= 1; \
(count) = __clz_c; \
} while (0)
#define gmp_ctz(count, x) do { \
mp_limb_t __ctz_x = (x); \
unsigned __ctz_c = 0; \
gmp_clz (__ctz_c, __ctz_x & - __ctz_x); \
(count) = GMP_LIMB_BITS - 1 - __ctz_c; \
} while (0)
#define gmp_add_ssaaaa(sh, sl, ah, al, bh, bl) \
do { \
mp_limb_t __x; \
__x = (al) + (bl); \
(sh) = (ah) + (bh) + (__x < (al)); \
(sl) = __x; \
} while (0)
#define gmp_sub_ddmmss(sh, sl, ah, al, bh, bl) \
do { \
mp_limb_t __x; \
__x = (al) - (bl); \
(sh) = (ah) - (bh) - ((al) < (bl)); \
(sl) = __x; \
} while (0)
#define gmp_umul_ppmm(w1, w0, u, v) \
do { \
int LOCAL_GMP_LIMB_BITS = GMP_LIMB_BITS; \
if (sizeof(unsigned int) * CHAR_BIT >= 2 * GMP_LIMB_BITS) \
{ \
unsigned int __ww = (unsigned int) (u) * (v); \
w0 = (mp_limb_t) __ww; \
w1 = (mp_limb_t) (__ww >> LOCAL_GMP_LIMB_BITS); \
} \
else if (GMP_ULONG_BITS >= 2 * GMP_LIMB_BITS) \
{ \
unsigned long int __ww = (unsigned long int) (u) * (v); \
w0 = (mp_limb_t) __ww; \
w1 = (mp_limb_t) (__ww >> LOCAL_GMP_LIMB_BITS); \
} \
else { \
mp_limb_t __x0, __x1, __x2, __x3; \
unsigned __ul, __vl, __uh, __vh; \
mp_limb_t __u = (u), __v = (v); \
\
__ul = __u & GMP_LLIMB_MASK; \
__uh = __u >> (GMP_LIMB_BITS / 2); \
__vl = __v & GMP_LLIMB_MASK; \
__vh = __v >> (GMP_LIMB_BITS / 2); \
\
__x0 = (mp_limb_t) __ul * __vl; \
__x1 = (mp_limb_t) __ul * __vh; \
__x2 = (mp_limb_t) __uh * __vl; \
__x3 = (mp_limb_t) __uh * __vh; \
\
__x1 += __x0 >> (GMP_LIMB_BITS / 2);/* this can't give carry */ \
__x1 += __x2; /* but this indeed can */ \
if (__x1 < __x2) /* did we get it? */ \
__x3 += GMP_HLIMB_BIT; /* yes, add it in the proper pos. */ \
\
(w1) = __x3 + (__x1 >> (GMP_LIMB_BITS / 2)); \
(w0) = (__x1 << (GMP_LIMB_BITS / 2)) + (__x0 & GMP_LLIMB_MASK); \
} \
} while (0)
#define gmp_udiv_qrnnd_preinv(q, r, nh, nl, d, di) \
do { \
mp_limb_t _qh, _ql, _r, _mask; \
gmp_umul_ppmm (_qh, _ql, (nh), (di)); \
gmp_add_ssaaaa (_qh, _ql, _qh, _ql, (nh) + 1, (nl)); \
_r = (nl) - _qh * (d); \
_mask = -(mp_limb_t) (_r > _ql); /* both > and >= are OK */ \
_qh += _mask; \
_r += _mask & (d); \
if (_r >= (d)) \
{ \
_r -= (d); \
_qh++; \
} \
\
(r) = _r; \
(q) = _qh; \
} while (0)
#define gmp_udiv_qr_3by2(q, r1, r0, n2, n1, n0, d1, d0, dinv) \
do { \
mp_limb_t _q0, _t1, _t0, _mask; \
gmp_umul_ppmm ((q), _q0, (n2), (dinv)); \
gmp_add_ssaaaa ((q), _q0, (q), _q0, (n2), (n1)); \
\
/* Compute the two most significant limbs of n - q'd */ \
(r1) = (n1) - (d1) * (q); \
gmp_sub_ddmmss ((r1), (r0), (r1), (n0), (d1), (d0)); \
gmp_umul_ppmm (_t1, _t0, (d0), (q)); \
gmp_sub_ddmmss ((r1), (r0), (r1), (r0), _t1, _t0); \
(q)++; \
\
/* Conditionally adjust q and the remainders */ \
_mask = - (mp_limb_t) ((r1) >= _q0); \
(q) += _mask; \
gmp_add_ssaaaa ((r1), (r0), (r1), (r0), _mask & (d1), _mask & (d0)); \
if ((r1) >= (d1)) \
{ \
if ((r1) > (d1) || (r0) >= (d0)) \
{ \
(q)++; \
gmp_sub_ddmmss ((r1), (r0), (r1), (r0), (d1), (d0)); \
} \
| 