[libmath-prime-util-perl] 36/55: More aggressive binomial reduction

[Partha P. Mukherjee]

This is an automated email from the git hooks/post-receive script.

ppm-guest pushed a commit to annotated tag v0.41
in repository libmath-prime-util-perl.

commit cc14faa3dae65d9ca75c9688c472f0f95697fe63
Author: Dana Jacobsen <dana at acm.org>
Date:   Mon May 12 17:48:09 2014 -0700

    More aggressive binomial reduction
---
 t/19-moebius.t |  7 ++++---
 util.c         | 14 ++++++++------
 2 files changed, 12 insertions(+), 9 deletions(-)

diff --git a/t/19-moebius.t b/t/19-moebius.t
index cba20f2..5be67ca 100644
--- a/t/19-moebius.t
+++ b/t/19-moebius.t
@@ -405,9 +405,10 @@ my @binomials = (
  [ 1,2, 0 ],
  [ 13,13, 1 ],
  [ 13,14, 0 ],
- [ 40,19, "131282408400" ],
- [ 67,31, "11923179284862717872" ],
- [ 228,12, "30689926618143230620" ],
+ [ 35,16, 4059928950 ],             # We can do this natively even in 32-bit
+ [ 40,19, "131282408400" ],         # We can do this in 64-bit
+ [ 67,31, "11923179284862717872" ], # ...and this
+ [ 228,12, "30689926618143230620" ],# But the result of this is too big.
  [ 177,78, "3314450882216440395106465322941753788648564665022000" ],
 );
 
diff --git a/util.c b/util.c
index 1702d7a..0563105 100644
--- a/util.c
+++ b/util.c
@@ -1326,16 +1326,18 @@ int kronecker_ss(IV a, IV b) {
 
 /* Thanks to MJD and RosettaCode */
 UV binomial(UV n, UV k) {
-  UV d, r = 1;
+  UV d, g, r = 1;
   if (k >= n) return (k == n);
   if (k > n/2) k = n-k;
   for (d = 1; d <= k; d++) {
     if (r >= UV_MAX/n) {  /* Possible overflow */
-      UV g = gcd_ui(r, d);
-      r /= g;
-      if (r >= UV_MAX/n) return 0;  /* Unavoidable overflow */
-      r *= n--;
-      r /= (d/g);
+      UV nr, dr;  /* reduced numerator / denominator */
+      g = gcd_ui(n, d);  nr = n/g;  dr = d/g;
+      g = gcd_ui(r, dr);  r = r/g;  dr = dr/g;
+      if (r >= UV_MAX/nr) return 0;  /* Unavoidable overflow */
+      r *= nr;
+      r /= dr;
+      n--;
     } else {
       r *= n--;
       r /= d;

-- 
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-perl/packages/libmath-prime-util-perl.git