[Gmp-commit] /var/hg/gmp: 4 new changesets
mercurial at gmplib.org
mercurial at gmplib.org
Tue Nov 19 15:16:15 UTC 2019
details: /var/hg/gmp/rev/500a0b8a9ff1
changeset: 17972:500a0b8a9ff1
user: Marco Bodrato <bodrato at mail.dm.unipi.it>
date: Tue Nov 19 15:22:19 2019 +0100
description:
mini-gmp/mini-gmp.c (mpn_invert_3by2): Add a new shortcut.
details: /var/hg/gmp/rev/f68dad45e793
changeset: 17973:f68dad45e793
user: Marco Bodrato <bodrato at mail.dm.unipi.it>
date: Tue Nov 19 15:47:07 2019 +0100
description:
mini-gmp/mini-gmp.c (mpn_invert_3by2): Move an assert earlier.
details: /var/hg/gmp/rev/0d1224f34ab7
changeset: 17974:0d1224f34ab7
user: Marco Bodrato <bodrato at mail.dm.unipi.it>
date: Tue Nov 19 15:47:28 2019 +0100
description:
mini-gmp/mini-gmp.c: Indent
details: /var/hg/gmp/rev/be2862e85ed4
changeset: 17975:be2862e85ed4
user: Marco Bodrato <bodrato at mail.dm.unipi.it>
date: Tue Nov 19 16:16:05 2019 +0100
description:
mini-gmp/tests: Remove a couple of warnings
diffstat:
mini-gmp/mini-gmp.c | 239 +++++++++++++++++++-----------------
mini-gmp/tests/t-mpq_muldiv_2exp.c | 2 +-
mini-gmp/tests/t-mpq_str.c | 2 +-
3 files changed, 127 insertions(+), 116 deletions(-)
diffs (truncated from 303 to 300 lines):
diff -r 29b2a74384c8 -r be2862e85ed4 mini-gmp/mini-gmp.c
--- a/mini-gmp/mini-gmp.c Mon Nov 18 00:23:03 2019 +0100
+++ b/mini-gmp/mini-gmp.c Tue Nov 19 16:16:05 2019 +0100
@@ -145,27 +145,27 @@
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); \
+ 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); \
+ __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; \
+ __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. */ \
+ __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); \
+ (w1) = __x3 + (__x1 >> (GMP_LIMB_BITS / 2)); \
+ (w0) = (__x1 << (GMP_LIMB_BITS / 2)) + (__x0 & GMP_LLIMB_MASK); \
} \
} while (0)
@@ -771,6 +771,8 @@
mpn_invert_3by2 (mp_limb_t u1, mp_limb_t u0)
{
int GMP_LIMB_BITS_MUL_3 = GMP_LIMB_BITS * 3;
+ assert (u1 >= GMP_LIMB_HIGHBIT);
+
if (sizeof (unsigned) * CHAR_BIT > GMP_LIMB_BITS * 3)
{
return (((unsigned) 1 << GMP_LIMB_BITS_MUL_3) - 1) /
@@ -782,107 +784,116 @@
(((unsigned long) u1 << GMP_LIMB_BITS_MUL_3 / 3) + u0);
}
else {
- mp_limb_t r, p, m, ql;
- unsigned ul, uh, qh;
-
- assert (u1 >= GMP_LIMB_HIGHBIT);
-
- /* For notation, let b denote the half-limb base, so that B = b^2.
- Split u1 = b uh + ul. */
- ul = u1 & GMP_LLIMB_MASK;
- uh = u1 >> (GMP_LIMB_BITS / 2);
-
- /* Approximation of the high half of quotient. Differs from the 2/1
- inverse of the half limb uh, since we have already subtracted
- u0. */
- qh = ~u1 / uh;
-
- /* Adjust to get a half-limb 3/2 inverse, i.e., we want
-
- qh' = floor( (b^3 - 1) / u) - b = floor ((b^3 - b u - 1) / u
- = floor( (b (~u) + b-1) / u),
-
- and the remainder
-
- r = b (~u) + b-1 - qh (b uh + ul)
- = b (~u - qh uh) + b-1 - qh ul
-
- Subtraction of qh ul may underflow, which implies adjustments.
- But by normalization, 2 u >= B > qh ul, so we need to adjust by
- at most 2.
- */
-
- r = ((~u1 - (mp_limb_t) qh * uh) << (GMP_LIMB_BITS / 2)) | GMP_LLIMB_MASK;
-
- p = (mp_limb_t) qh * ul;
- /* Adjustment steps taken from udiv_qrnnd_c */
- if (r < p)
- {
- qh--;
- r += u1;
- if (r >= u1) /* i.e. we didn't get carry when adding to r */
+ mp_limb_t r, m;
+
+ if (GMP_ULONG_BITS >= GMP_LIMB_BITS * 2)
+ {
+ /* Set m to the 2/1 inverse of u1. */
+ m = ~((unsigned long) u1 << GMP_LIMB_BITS_MUL_3 / 3) / u1;
+ r = ~(m * u1);
+ }
+ else
+ {
+ mp_limb_t p, ql;
+ unsigned ul, uh, qh;
+
+ /* For notation, let b denote the half-limb base, so that B = b^2.
+ Split u1 = b uh + ul. */
+ ul = u1 & GMP_LLIMB_MASK;
+ uh = u1 >> (GMP_LIMB_BITS / 2);
+
+ /* Approximation of the high half of quotient. Differs from the 2/1
+ inverse of the half limb uh, since we have already subtracted
+ u0. */
+ qh = ~u1 / uh;
+
+ /* Adjust to get a half-limb 3/2 inverse, i.e., we want
+
+ qh' = floor( (b^3 - 1) / u) - b = floor ((b^3 - b u - 1) / u
+ = floor( (b (~u) + b-1) / u),
+
+ and the remainder
+
+ r = b (~u) + b-1 - qh (b uh + ul)
+ = b (~u - qh uh) + b-1 - qh ul
+
+ Subtraction of qh ul may underflow, which implies adjustments.
+ But by normalization, 2 u >= B > qh ul, so we need to adjust by
+ at most 2.
+ */
+
+ r = ((~u1 - (mp_limb_t) qh * uh) << (GMP_LIMB_BITS / 2)) | GMP_LLIMB_MASK;
+
+ p = (mp_limb_t) qh * ul;
+ /* Adjustment steps taken from udiv_qrnnd_c */
if (r < p)
{
qh--;
r += u1;
+ if (r >= u1) /* i.e. we didn't get carry when adding to r */
+ if (r < p)
+ {
+ qh--;
+ r += u1;
+ }
}
- }
- r -= p;
-
- /* Low half of the quotient is
-
- ql = floor ( (b r + b-1) / u1).
-
- This is a 3/2 division (on half-limbs), for which qh is a
- suitable inverse. */
-
- p = (r >> (GMP_LIMB_BITS / 2)) * qh + r;
- /* Unlike full-limb 3/2, we can add 1 without overflow. For this to
- work, it is essential that ql is a full mp_limb_t. */
- ql = (p >> (GMP_LIMB_BITS / 2)) + 1;
-
- /* By the 3/2 trick, we don't need the high half limb. */
- r = (r << (GMP_LIMB_BITS / 2)) + GMP_LLIMB_MASK - ql * u1;
-
- if (r >= (p << (GMP_LIMB_BITS / 2)))
- {
- ql--;
- r += u1;
- }
- m = ((mp_limb_t) qh << (GMP_LIMB_BITS / 2)) + ql;
- if (r >= u1)
- {
- m++;
- r -= u1;
- }
-
- /* Now m is the 2/1 inverse of u1. If u0 > 0, adjust it to become a
- 3/2 inverse. */
- if (u0 > 0)
- {
- mp_limb_t th, tl;
- r = ~r;
- r += u0;
- if (r < u0)
- {
- m--;
- if (r >= u1)
- {
- m--;
- r -= u1;
- }
- r -= u1;
- }
- gmp_umul_ppmm (th, tl, u0, m);
- r += th;
- if (r < th)
- {
- m--;
- m -= ((r > u1) | ((r == u1) & (tl > u0)));
- }
- }
-
- return m;
+ r -= p;
+
+ /* Low half of the quotient is
+
+ ql = floor ( (b r + b-1) / u1).
+
+ This is a 3/2 division (on half-limbs), for which qh is a
+ suitable inverse. */
+
+ p = (r >> (GMP_LIMB_BITS / 2)) * qh + r;
+ /* Unlike full-limb 3/2, we can add 1 without overflow. For this to
+ work, it is essential that ql is a full mp_limb_t. */
+ ql = (p >> (GMP_LIMB_BITS / 2)) + 1;
+
+ /* By the 3/2 trick, we don't need the high half limb. */
+ r = (r << (GMP_LIMB_BITS / 2)) + GMP_LLIMB_MASK - ql * u1;
+
+ if (r >= (p << (GMP_LIMB_BITS / 2)))
+ {
+ ql--;
+ r += u1;
+ }
+ m = ((mp_limb_t) qh << (GMP_LIMB_BITS / 2)) + ql;
+ if (r >= u1)
+ {
+ m++;
+ r -= u1;
+ }
+ }
+
+ /* Now m is the 2/1 inverse of u1. If u0 > 0, adjust it to become a
+ 3/2 inverse. */
+ if (u0 > 0)
+ {
+ mp_limb_t th, tl;
+ r = ~r;
+ r += u0;
+ if (r < u0)
+ {
+ m--;
+ if (r >= u1)
+ {
+ m--;
+ r -= u1;
+ }
+ r -= u1;
+ }
+ gmp_umul_ppmm (th, tl, u0, m);
+ r += th;
+ if (r < th)
+ {
+ m--;
+ m -= ((r > u1) | ((r == u1) & (tl > u0)));
+ }
+ }
+
+ return m;
}
}
@@ -3334,7 +3345,7 @@
mpz_fac_ui (t, k);
for (; k > 0; --k)
- mpz_mul_ui (r, r, n--);
+ mpz_mul_ui (r, r, n--);
mpz_divexact (r, r, t);
mpz_clear (t);
diff -r 29b2a74384c8 -r be2862e85ed4 mini-gmp/tests/t-mpq_muldiv_2exp.c
--- a/mini-gmp/tests/t-mpq_muldiv_2exp.c Mon Nov 18 00:23:03 2019 +0100
+++ b/mini-gmp/tests/t-mpq_muldiv_2exp.c Tue Nov 19 16:16:05 2019 +0100
@@ -100,7 +100,7 @@
|| mpz_sizeinbase (t, 2) - 1 != e || mpz_cmp_ui (mpq_denref (aq), 1) != 0)
{
fprintf (stderr, "mpq_div_2exp failed: %lu\n", e);
- fprintf (stderr, "%li %li %lu %li\n", e2, t2, mpz_scan1 (t, 0), mpz_sizeinbase (t, 2));
+ fprintf (stderr, "%li %li %lu %zu\n", e2, t2, mpz_scan1 (t, 0), mpz_sizeinbase (t, 2));
dump ("na", a);
dump ("da", b);
dump ("nr", mpq_numref (rq));
diff -r 29b2a74384c8 -r be2862e85ed4 mini-gmp/tests/t-mpq_str.c
--- a/mini-gmp/tests/t-mpq_str.c Mon Nov 18 00:23:03 2019 +0100
+++ b/mini-gmp/tests/t-mpq_str.c Tue Nov 19 16:16:05 2019 +0100
@@ -144,7 +144,7 @@
char *ap;
char *bp;
char *rp;
- size_t bn, rn, arn;
+ size_t rn, arn;
More information about the gmp-commit
mailing list