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

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

commit 06eadfb62748fc57f10b8ffbf4798edc5071e31b
Author: Dana Jacobsen <d...@acm.org>
Date:   Sun Jan 5 03:42:40 2014 -0800

    Tweaks to SQUFOF, pbrent, and strategy
---
 factor.c | 117 ++++++++++++++++++++++++++++++++-------------------------------
 1 file changed, 60 insertions(+), 57 deletions(-)

diff --git a/factor.c b/factor.c
index ed61fea..b6ed3e0 100644
--- a/factor.c
+++ b/factor.c
@@ -95,8 +95,8 @@ int factor(UV n, UV *factors)
     while ( (n >= f*f) && (!_XS_is_prime(n)) ) {
       int split_success = 0;
       /* Adjust the number of rounds based on the number size */
-      UV const br_rounds = ((n>>29) < 100000) ?  1500 :  1500;
-      UV const sq_rounds = 80000; /* 20k 91%, 40k 98%, 80k 99.9%, 120k 99.99% 
*/
+      UV const br_rounds = ((n>>29) < 100000) ?  1500 :  2000;
+      UV const sq_rounds =100000; /* 20k 91%, 40k 98%, 80k 99.9%, 120k 99.99% 
*/
 
       /* 99.7% of 32-bit, 94% of 64-bit random inputs factored here */
       if (!split_success) {
@@ -111,7 +111,7 @@ int factor(UV n, UV *factors)
       /* At this point we should only have 16+ digit semiprimes. */
       /* This p-1 gets about 2/3 of what makes it through the above */
       if (!split_success) {
-        split_success = pminus1_factor(n, tofac_stack+ntofac, 5000, 80000)-1;
+        split_success = pminus1_factor(n, tofac_stack+ntofac, 5000, 100000)-1;
         if (verbose) printf("pminus1 %d\n", split_success);
       }
       /* Some rounds of HOLF, good for close to perfect squares which are
@@ -446,7 +446,7 @@ int holf_factor(UV n, UV *factors, UV rounds)
 /* Pollard / Brent.  Brent's modifications to Pollard's Rho.  Maybe faster. */
 int pbrent_factor(UV n, UV *factors, UV rounds, UV a)
 {
-  UV f, i, r;
+  UV f, m, r;
   UV Xi = 2;
   UV Xm = 2;
   const UV inner = 64;
@@ -460,15 +460,16 @@ int pbrent_factor(UV n, UV *factors, UV rounds, UV a)
     /* Do rleft rounds, inner at a time */
     while (rleft > 0) {
       UV dorounds = (rleft > inner) ? inner : rleft;
-      UV m = 1;
       saveXi = Xi;
-      for (i = 0; i < dorounds; i++) {
+      rleft -= dorounds;
+      rounds -= dorounds;
+      Xi = sqraddmod(Xi, a, n);        /* First iteration, no mulmod needed */
+      m = (Xi>Xm) ? Xi-Xm : Xm-Xi;
+      while (--dorounds > 0) {         /* Now do inner-1=63 more iterations */
         Xi = sqraddmod(Xi, a, n);
         f = (Xi>Xm) ? Xi-Xm : Xm-Xi;
         m = mulmod(m, f, n);
       }
-      rleft -= dorounds;
-      rounds -= dorounds;
       f = gcd_ui(m, n);
       if (f != 1)
         break;
@@ -739,8 +740,7 @@ int pplus1_factor(UV n, UV *factors, UV B1)
 
 typedef struct
 {
-  UV mult;
-  UV valid;
+  int valid;
   UV P;
   UV bn;
   UV Qn;
@@ -750,14 +750,12 @@ typedef struct
   UV imax;
 } mult_t;
 
-/* N < 2^63 (or 2^31).  *f == 1 if no factor found */
-static void squfof_unit(UV n, mult_t* mult_save, UV* f)
+/* N < 2^63 (or 2^31).  Returns 0 or a factor */
+static UV squfof_unit(UV n, mult_t* mult_save)
 {
   UV imax,i,Q0,b0,Qn,bn,P,bbn,Ro,S,So,t1,t2;
   int j;
 
-  *f = 0;
-
   P = mult_save->P;
   bn = mult_save->bn;
   Qn = mult_save->Qn;
@@ -790,8 +788,7 @@ static void squfof_unit(UV n, mult_t* mult_save, UV* f)
         mult_save->Qn = Qn;
         mult_save->Q0 = Q0;
         mult_save->it = i;
-        *f = 0;
-        return;
+        return 0;
       }
 
       SQUARE_SEARCH_ITERATION;
@@ -830,25 +827,34 @@ static void squfof_unit(UV n, mult_t* mult_save, UV* f)
       SYMMETRY_POINT_ITERATION;
       if (j++ > 2000000) {
          mult_save->valid = 0;
-         *f = 0;
-         return;
+         return 0;
       }
     }
 
-    *f = gcd_ui(Ro, n);
-    if (*f > 1)
-      return;
+    t1 = gcd_ui(Ro, n);
+    if (t1 > 1)
+      return t1;
   }
 }
 
+/* Gower and Wagstaff 2008:
+ *    http://www.ams.org/journals/mcom/2008-77-261/S0025-5718-07-02010-8/
+ * Section 5.3.  I've added some with 13,17,19.  Sorted by F(). */
 static const UV squfof_multipliers[] =
-  { 3*5*7*11, 3*5*7, 3*5*11, 3*5, 3*7*11, 3*7, 5*7*11, 5*7,
-    3*11,     3,     5*11,   5,   7*11,   7,   11,     1     };
+  /* { 3*5*7*11, 3*5*7, 3*5*11, 3*5, 3*7*11, 3*7, 5*7*11, 5*7,
+       3*11,     3,     5*11,   5,   7*11,   7,   11,     1   }; */
+  { 3*5*7*11, 3*5*7,  3*5*7*11*13, 3*5*7*13, 3*5*7*11*17, 3*5*11,
+    3*5*7*17, 3*5,    3*5*7*11*19, 3*5*11*13,3*5*7*19,    3*5*7*13*17,
+    3*5*13,   3*7*11, 3*7,         5*7*11,   3*7*13,      5*7,
+    3*5*17,   5*7*13, 3*5*19,      3*11,     3*7*17,      3,
+    3*11*13,  5*11,   3*7*19,      3*13,     5,           5*11*13,
+    5*7*19,   5*13,   7*11,        7,        3*17,        7*13,
+    11,       1 };
 #define NSQUFOF_MULT (sizeof(squfof_multipliers)/sizeof(squfof_multipliers[0]))
 
 int squfof_factor(UV n, UV *factors, UV rounds)
 {
-  const UV big2 = UV_MAX >> 2;
+  const UV big2 = UV_MAX;
   mult_t mult_save[NSQUFOF_MULT];
   int still_racing;
   UV i, nn64, mult, f64;
@@ -862,45 +868,42 @@ int squfof_factor(UV n, UV *factors, UV rounds)
     factors[0] = n;  return 1;
   }
 
-  for (i = 0; i < NSQUFOF_MULT; i++) {
-    mult = squfof_multipliers[i];
-    nn64 = n * mult;
-    mult_save[i].mult = mult;
-    if ((big2 / mult) < n) {
-      mult_save[i].valid = 0; /* This multiplier would overflow 64-bit */
-      continue;
-    }
-    mult_save[i].valid = 1;
-
-    mult_save[i].b0 = isqrt(nn64);
-    mult_save[i].imax = (UV) (sqrt(mult_save[i].b0) / 3);
-    if (mult_save[i].imax < 20)     mult_save[i].imax = 20;
-    if (mult_save[i].imax > rounds) mult_save[i].imax = rounds;
-
-    mult_save[i].Q0 = 1;
-    mult_save[i].P  = mult_save[i].b0;
-    mult_save[i].Qn = nn64 - (mult_save[i].b0 * mult_save[i].b0);
-    if (mult_save[i].Qn == 0) {
-      factors[0] = mult_save[i].b0;
-      factors[1] = n / mult_save[i].b0;
-      MPUassert( factors[0] * factors[1] == n , "incorrect factoring");
-      return 2;
-    }
-    mult_save[i].bn = (mult_save[i].b0 + mult_save[i].P) / mult_save[i].Qn;
-    mult_save[i].it = 0;
-  }
+  for (i = 0; i < NSQUFOF_MULT; i++)
+    mult_save[i].valid = -1;
 
   /* Process the multipliers a little at a time: 0.33*(n*mult)^1/4: 20-20k */
   do {
     still_racing = 0;
     for (i = 0; i < NSQUFOF_MULT; i++) {
-      if (!mult_save[i].valid)
-        continue;
-      nn64 = n * squfof_multipliers[i];
-      squfof_unit(nn64, &mult_save[i], &f64);
+      if (mult_save[i].valid == 0)  continue;
+      mult = squfof_multipliers[i];
+      nn64 = n * mult;
+      if (mult_save[i].valid == -1) {
+        if ((big2 / mult) < n) {
+          mult_save[i].valid = 0; /* This multiplier would overflow 64-bit */
+          continue;
+        }
+        mult_save[i].valid = 1;
+        mult_save[i].b0 = isqrt(nn64);
+        mult_save[i].imax = (UV) (sqrt(mult_save[i].b0) / 16);
+        if (mult_save[i].imax < 20)     mult_save[i].imax = 20;
+        if (mult_save[i].imax > rounds) mult_save[i].imax = rounds;
+        mult_save[i].Q0 = 1;
+        mult_save[i].P  = mult_save[i].b0;
+        mult_save[i].Qn = nn64 - (mult_save[i].b0 * mult_save[i].b0);
+        if (mult_save[i].Qn == 0) {
+          factors[0] = mult_save[i].b0;
+          factors[1] = n / mult_save[i].b0;
+          MPUassert( factors[0] * factors[1] == n , "incorrect factoring");
+          return 2;
+        }
+        mult_save[i].bn = (mult_save[i].b0 + mult_save[i].P) / mult_save[i].Qn;
+        mult_save[i].it = 0;
+      }
+      f64 = squfof_unit(nn64, &mult_save[i]);
       if (f64 > 1) {
-        if (f64 != squfof_multipliers[i]) {
-          f64 /= gcd_ui(f64, squfof_multipliers[i]);
+        if (f64 != mult) {
+          f64 /= gcd_ui(f64, mult);
           if (f64 != 1) {
             factors[0] = f64;
             factors[1] = n / f64;

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

_______________________________________________
Pkg-perl-cvs-commits mailing list
Pkg-perl-cvs-commits@lists.alioth.debian.org
http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/pkg-perl-cvs-commits

Reply via email to