[Ecm-commits] r2614 - trunk

cvs commits ecm-commits at lists.gforge.inria.fr
Sat Feb 14 13:03:40 CET 2015


Author: zimmerma
Date: 2015-02-14 13:03:40 +0100 (Sat, 14 Feb 2015)
New Revision: 2614

Modified:
   trunk/ecm-impl.h
   trunk/listz.c
   trunk/polyeval.c
   trunk/stage2.c
Log:
removed old polyeval() code [not tested any more, make check did fail with it]


Modified: trunk/ecm-impl.h
===================================================================
--- trunk/ecm-impl.h	2015-02-14 11:38:55 UTC (rev 2613)
+++ trunk/ecm-impl.h	2015-02-14 12:03:40 UTC (rev 2614)
@@ -107,10 +107,6 @@
 #define PM1FS2_COST 1.0 / 4.0
 #define PP1FS2_COST 1.0 / 4.0
 
-/* if POLYEVALTELLEGEN is defined, use polyeval_tellegen(),
-   otherwise use polyeval() */
-#define POLYEVALTELLEGEN
-
 /* define top-level multiplication */
 #define KARA 2
 #define TOOM3 3
@@ -123,10 +119,6 @@
 
 #include "sp.h"
 
-#ifdef POLYEVALTELLEGEN
-#define USE_SHORT_PRODUCT
-#endif
-
 #include <assert.h>
 #define ASSERT_ALWAYS(expr)   assert (expr)
 #ifdef WANT_ASSERT

Modified: trunk/listz.c
===================================================================
--- trunk/listz.c	2015-02-14 11:38:55 UTC (rev 2613)
+++ trunk/listz.c	2015-02-14 12:03:40 UTC (rev 2614)
@@ -218,21 +218,6 @@
     mpz_sub (p[i], q[i], r[i]);
 }
 
-#ifndef POLYEVALTELLEGEN
-/* p[i] <- q[i] * r mod m */
-void
-list_mul_z (listz_t p, listz_t q, mpz_t r, unsigned int n, mpz_t m)
-{
-  unsigned int i;
-
-  for (i = 0; i < n; i++)
-    {
-      mpz_mul (p[i], q[i], r);
-      mpz_mod (p[i], p[i], m);
-    }
-}
-#endif
-
 /* Multiply up the integers in l, modulo n. Each entry becomes the
    product (mod n) of itself and all previous entries */
    
@@ -530,85 +515,7 @@
     }
 }
 
-#ifndef POLYEVALTELLEGEN
 /*
-  divides a[0]+a[1]*x+...+a[2K-1]*x^(2K-1)
-  By b[0]+b[1]*x+...+b[K-1]*x^(K-1)+x^K
-  i.e. a polynomial of 2K coefficients divided by a monic polynomial
-  with K+1 coefficients (b[K]=1 is implicit).
-  Puts the quotient in q[0]+q[1]*x+...+q[K-1]*x^(K-1)
-  and the remainder in a[0]+a[1]*x+...+a[K-1]*x^(K-1)
-  Needs space for list_mul_mem(K) coefficients in t.
-  If top is non-zero, a[0]..a[K-1] are reduced mod n.
-*/
-void
-RecursiveDivision (listz_t q, listz_t a, listz_t b, unsigned int K,
-                   listz_t t, mpz_t n, int top)
-{
-  if (K == 1) /* a0+a1*x = a1*(b0+x) + a0-a1*b0 */
-    {
-      mpz_mod (a[1], a[1], n);
-      mpz_mul (q[0], a[1], b[0]);
-      mpz_mod (q[0], q[0], n);
-      mpz_sub (a[0], a[0], q[0]);
-      if (top)
-        mpz_mod (a[0], a[0], n);
-      mpz_set (q[0], a[1]);
-    }
-  else
-    {
-      unsigned int k, l, i, po2;
-
-      k = K / 2;
-      l = K - k;
-      for (po2 = K; (po2 && 1) == 0; po2 >>= 1);
-      po2 = (po2 == 1);
-
-      /* first perform a (2l) / l division */
-      RecursiveDivision (q + k, a + 2 * k, b + k, l, t, n, 0);
-      /* subtract q[k..k+l-1] * b[0..k-1] */
-      ASSERTD(list_check(q+l,k,n) && list_check(b,k,n));
-      if (po2 && Fermat)
-        F_mul (t, q + l, b, k, DEFAULT, Fermat, t + K); /* sets t[0..2*k-2]*/
-      else
-        list_mult_n (t, q + l, b, k); /* sets t[0..2*k-2] */
-      list_sub (a + l, a + l, t, 2 * k - 1);
-      if (k < l) /* don't forget to subtract q[k] * b[0..k-1] */
-        {
-	  for (i=0; i<k; i++)
-	    {
-	      mpz_mul (t[0], q[k], b[i]); /* TODO: need to reduce t[0]? */
-	      mpz_sub (a[k+i], a[k+i], t[0]);
-	    }
-        }
-      /* remainder is in a[0..K+k-1] */
-
-      /* then perform a (2k) / k division */
-      RecursiveDivision (q, a + l, b + l, k, t, n, 0);
-      /* subtract q[0..k-1] * b[0..l-1] */
-      ASSERTD(list_check(q,k,n) && list_check(b,k,n));
-      if (po2 && Fermat)
-        F_mul (t, q, b, k, DEFAULT, Fermat, t + K);
-      else
-        list_mult_n (t, q, b, k);
-      list_sub (a, a, t, 2 * k - 1);
-      if (k < l) /* don't forget to subtract q[0..k-1] * b[k] */
-        {
-          for (i=0; i<k; i++)
-            {
-              mpz_mul (t[0], q[i], b[k]); /* TODO: need to reduce t[0]? */
-              mpz_sub (a[k+i], a[k+i], t[0]);
-            }
-        }
-
-      /* normalizes the remainder wrt n */
-      if (top)
-        list_mod (a, a, K, n);
-    }
-}
-#endif
-
-/*
   Returns in a[0]+a[1]*x+...+a[K-1]*x^(K-1)
   the remainder of the division of
   A = a[0]+a[1]*x+...+a[2K-2]*x^(2K-2)

Modified: trunk/polyeval.c
===================================================================
--- trunk/polyeval.c	2015-02-14 11:38:55 UTC (rev 2613)
+++ trunk/polyeval.c	2015-02-14 12:03:40 UTC (rev 2614)
@@ -37,83 +37,6 @@
 
 extern unsigned int Fermat;
 
-#ifndef POLYEVALTELLEGEN
-/* algorithm polyeval from section 3.7 of Peter Montgomery's dissertation.
-Input: 
-   G - an array of k elements of R, G[i], 0 <= i < k
-       representing the coefficients of a polynomial G(x) of degree < k
-   Tree - the product tree produced by PolyFromRoots
-   Tree[0][0..k-1] (degree k/2)
-   Tree[1][0..k-1] (degree k/4), ...,
-   Tree[lgk-1][0..k-1] (degree 1)
-Output: the sequence of values of G(a[i]) are stored in G[i] for 0 <= i < k
-Remark: we need an auxiliary (k+1)-th cell G[k] in G.
-The memory used is M(k) = max(3*floor(k/2)+list_mul_mem(floor(k/2)),
-                              k+list_mul_mem(ceil(k/2)),
-                              floor(k/2) + M(ceil(k/2))).
-Since list_mul_mem(k) >= 2*k, the maximum is the 1st.
-*/
-void
-polyeval (listz_t G, unsigned int k, listz_t *Tree, listz_t T, mpz_t n,
-          unsigned int sh)
-{
-  unsigned int l, m;
-  listz_t T0;
-
-  if (k == 1)
-    return;
-
-  T0 = Tree[0] + sh;
-  
-  m = k / 2;
-  l = k - m;
-
-  /* divide G[0]+G[1]*x+...+G[k-1]*x^(k-1) by
-            T0[l]+...+T0[k-1]*x^(m-1)+x^m,
-            quotient in {T+m,l-1}, remainder in {T,m} */
-
-  if (k == 2 * m)
-    {
-      /* FIXME: avoid the copy here by giving different 2nd and 3rd arguments
-         to RecursiveDivision */
-      list_set (T, G, k);
-      /* the following needs k+m+list_mul_mem(m) in T */
-      RecursiveDivision (T + k, T, T0 + l, m, T + k + m, n, 1);
-    }
-  else /* k = 2m+1: subtract G[k-1]*x^(l-1) * T0 from G */
-    {
-      /* G - G[k-1] * (x^m + {T0+l,m}) * x^m */
-      list_set (T, G, m);
-      list_mul_z (T + m, T0 + l, G[k - 1], m, n);
-      list_sub (T + m, G + m, T + m, m);
-      /* the following needs 3m+list_mul_mem(m) in T */
-      RecursiveDivision (T + 2 * m, T, T0 + l, m, T + 3 * m, n, 1);
-    }
-  /* in both cases we need 3*(k/2)+list_mul_mem(k/2) */
-
-  /* right remainder is in {T,m} */
-
-  /* k = 2l or k = 2l-1 */
-  
-  /* divide G[0]+G[1]*x+...+G[k-1]*x^(k-1) by
-            T0[0]+...+T0[l-1]*x^(l-1)+x^l:
-            quotient in {T+m,m-1}, remainder in {G,l} */
-
-  if (k < 2 * l)
-    mpz_set_ui (G[k], 0);
-  /* the following needs k+list_mul_mem(l) in T */
-  RecursiveDivision (T + m, G, T0, l, T + k, n, 1);
-
-  /* left remainder is in {G,l} */
-  
-  polyeval (G, l, Tree + 1, T + m, n, sh);
-
-  /* copy right remainder in {G+l,m} */
-  list_set (G + l, T, m);
-  polyeval (G + l, m, Tree + 1, T, n, sh + l);
-}
-#endif
-
 #if defined(DEBUG) || defined(DEBUG_TREEDATA)
 void
 print_vect (listz_t t, unsigned int l)
@@ -252,6 +175,23 @@
     return r1;
 }
 
+/* This is the documentation of the (now removed) polyeval() function.
+   Algorithm polyeval from section 3.7 of Peter Montgomery's dissertation.
+Input: 
+   G - an array of k elements of R, G[i], 0 <= i < k
+       representing the coefficients of a polynomial G(x) of degree < k
+   Tree - the product tree produced by PolyFromRoots
+   Tree[0][0..k-1] (degree k/2)
+   Tree[1][0..k-1] (degree k/4), ...,
+   Tree[lgk-1][0..k-1] (degree 1)
+Output: the sequence of values of G(a[i]) are stored in G[i] for 0 <= i < k
+Remark: we need an auxiliary (k+1)-th cell G[k] in G.
+The memory used is M(k) = max(3*floor(k/2)+list_mul_mem(floor(k/2)),
+                              k+list_mul_mem(ceil(k/2)),
+                              floor(k/2) + M(ceil(k/2))).
+Since list_mul_mem(k) >= 2*k, the maximum is the 1st.
+*/
+
 /* Same as polyeval. Needs invF as extra argument.
    Return non-zero iff an error occurred.
 */
@@ -269,11 +209,7 @@
     ASSERT(Tree != NULL || TreeFilename != NULL);
     
     tupspace = TUpTree_space (k) + k;
-#ifndef USE_SHORT_PRODUCT
-    tkspace = TMulGen_space (k - 1, k - 1, k - 1) + k;
-#else
     tkspace = 2 * k - 1 + list_mul_mem (k);
-#endif
 
     tupspace = MAX (tupspace, tkspace);
     
@@ -306,15 +242,9 @@
       }
     else
       {
-#ifdef USE_SHORT_PRODUCT
         /* need space 2k-1+list_mul_mem(k) in T */
         list_mul_high (T, invF, b, k);
         list_mod (T, T + k - 1, k, n);
-#else
-        /* revert invF for call to TMulGen below */
-        list_revert (invF, k);
-        TMulGen (T, k - 1, invF, k - 1, b, k - 1, T + k, n);
-#endif
       }
     list_revert (T, k);
     if (TreeFilename != NULL)

Modified: trunk/stage2.c
===================================================================
--- trunk/stage2.c	2015-02-14 11:38:55 UTC (rev 2613)
+++ trunk/stage2.c	2015-02-14 12:03:40 UTC (rev 2614)
@@ -761,7 +761,6 @@
   clear_list (G, dF);
   G = NULL;
   st = cputime ();
-#ifdef POLYEVALTELLEGEN
   if (use_ntt)
     youpi = ntt_polyevalT (T, dF, Tree, T + dF + 1, sp_invF,
 	mpzspm, TreeFilename);
@@ -774,11 +773,6 @@
       outputf (OUTPUT_ERROR, "Error, not enough memory\n");
       goto clear_fd;
     }
-#else
-  clear_list (invF, dF + 1);
-  invF = NULL;
-  polyeval (T, dF, Tree, T + dF + 1, n, 0);
-#endif
   treefiles_used = 0; /* Polyeval deletes treefiles by itself */
 
   if (test_verbose (OUTPUT_TRACE))



More information about the Ecm-commits mailing list