Bill Allombert on Wed, 31 Aug 2005 14:52:44 +0200


[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

patch for a new divide_conquer_assoc function


Hello PARI-dev,

Here a patch that add a new function divide_conquer_assoc:

divide_conquer_assoc(GEN x, GEN (*mul)(void*,GEN,GEN), void *data);

This function is essentially the same than divide_conquer_prod but for
the extra data parameter. This avoid the need of static variables,
which are not thread-safe.

A better name for this function would be much welcome.

I elected to keep divide_conquer_prod() interface unchanged to preserve
backward compatibility and because there are several direct use of it
in libpari with arguments like gmul or mulii, though maybe they could
be replaced by a dedicated function ZV_prod or RgV_prod.

Cheers,
Bill.

Index: pari/src/basemath/Flx.c
===================================================================
--- pari.orig/src/basemath/Flx.c	2005-08-30 16:18:31.000000000 +0200
+++ pari/src/basemath/Flx.c	2005-08-30 23:23:21.000000000 +0200
@@ -1211,11 +1211,10 @@
   return p1;
 }
 
-static ulong global_pp;
 static GEN
-_Flx_mul(GEN a, GEN b)
+_Flx_mul(void *p, GEN a, GEN b)
 {
-  return Flx_mul(a,b, global_pp);
+  return Flx_mul(a,b, (ulong)p);
 }
 
 /* compute prod (x - a[i]) */
@@ -1225,7 +1224,7 @@
   long i,k,lx = lg(a);
   GEN p1,p2;
   if (lx == 1) return Fl_to_Flx(1,vs);
-  p1 = cgetg(lx, t_VEC); global_pp = p;
+  p1 = cgetg(lx, t_VEC); 
   for (k=1,i=1; i<lx-1; i+=2)
   {
     p2 = cgetg(5,t_VECSMALL); gel(p1,k++) = p2;
@@ -1243,7 +1242,7 @@
     p2[2] = a[i]?p - a[i]:0;
     p2[3] = 1;
   }
-  setlg(p1, k); return divide_conquer_prod(p1, _Flx_mul);
+  setlg(p1, k); return divide_conquer_assoc(p1, _Flx_mul,(void *)p);
 }
 
 GEN
@@ -2125,15 +2124,18 @@
   return gerepileupto(av0, y);
 }
 
-static GEN Tmodulo;
-static ulong modulo; 
-static GEN _FlxqX_mul(GEN a,GEN b){return FlxqX_mul(a,b,Tmodulo,modulo);}
+struct _FlxqX {ulong p; GEN T;};
+static GEN _FlxqX_mul(void *data,GEN a,GEN b)
+{
+  struct _FlxqX *d=(struct _FlxqX*)data;
+  return FlxqX_mul(a,b,d->T,d->p);
+}
 
 GEN 
 FlxqXV_prod(GEN V, GEN T, ulong p)
 {
-  modulo = p; Tmodulo = T;
-  return divide_conquer_prod(V, &_FlxqX_mul);
+  struct _FlxqX d; d.p=p; d.T=T;
+  return divide_conquer_assoc(V, &_FlxqX_mul, (void*)&d);
 }
 
 GEN
Index: pari/src/basemath/polarit2.c
===================================================================
--- pari.orig/src/basemath/polarit2.c	2005-08-30 16:18:31.000000000 +0200
+++ pari/src/basemath/polarit2.c	2005-08-30 23:23:21.000000000 +0200
@@ -2321,7 +2321,7 @@
 }
 
 GEN
-divide_conquer_prod(GEN x, GEN (*mul)(GEN,GEN))
+divide_conquer_assoc(GEN x, GEN (*mul)(void *,GEN,GEN),void *data)
 {
   pari_sp ltop, lim;
   long i,k,lx = lg(x);
@@ -2336,7 +2336,7 @@
       fprintferr("prod: remaining objects %ld\n",k-1);
     lx = k; k = 1;
     for (i=1; i<lx-1; i+=2)
-      gel(x,k++) = mul(gel(x,i),gel(x,i+1));
+      gel(x,k++) = mul(data,gel(x,i),gel(x,i+1));
     if (i < lx) x[k++] = x[i];
     if (low_stack(lim,stack_lim(av,1)))
       gerepilecoeffs(ltop,x+1,k-1);
@@ -2344,30 +2344,41 @@
   return gel(x,1);
 }
 
-static GEN static_OBJ; /* nf or ell */
+static GEN
+_domul(void *data, GEN x, GEN y)
+{
+  GEN (*mul)(GEN,GEN)=(GEN (*)(GEN,GEN)) data;
+  return mul(x,y);
+}
+
+GEN
+divide_conquer_prod(GEN x, GEN (*mul)(GEN,GEN))
+{
+  return divide_conquer_assoc(x, _domul, (void *)mul);
+}
 
 static GEN
-idmulred(GEN x, GEN y) { return idealmulred(static_OBJ, x, y, 0); }
+idmulred(void *nf, GEN x, GEN y) { return idealmulred((GEN) nf, x, y, 0); }
 static GEN
-idpowred(GEN x, GEN n) { return idealpowred(static_OBJ, x, n, 0); }
+idpowred(void *nf, GEN x, GEN n) { return idealpowred((GEN) nf, x, n, 0); }
 static GEN
-idmul(GEN x, GEN y) { return idealmul(static_OBJ, x, y); }
+idmul(void *nf, GEN x, GEN y) { return idealmul((GEN) nf, x, y); }
 static GEN
-idpow(GEN x, GEN n) { return idealpow(static_OBJ, x, n); }
+idpow(void *nf, GEN x, GEN n) { return idealpow((GEN) nf, x, n); }
 static GEN
-eltmul(GEN x, GEN y) { return element_mul(static_OBJ, x, y); }
+eltmul(void *nf, GEN x, GEN y) { return element_mul((GEN) nf, x, y); }
 static GEN
-eltpow(GEN x, GEN n) { return element_pow(static_OBJ, x, n); }
+eltpow(void *nf, GEN x, GEN n) { return element_pow((GEN) nf, x, n); }
 
 #if 0
 static GEN
-ellmul(GEN x, GEN y) { return powell(static_OBJ, x, y); }
+ellmul(void *ell, GEN x, GEN y) { return addell((GEN) ell, x, y); }
 static GEN
-ellpow(GEN x, GEN n) { return idealpow(static_OBJ, x, n); }
+ellpow(void *GEN x, GEN n) { return powell((GEN) ell, x, n); }
 #endif
 
 GEN
-_factorback(GEN fa, GEN e, GEN (*_mul)(GEN,GEN), GEN (*_pow)(GEN,GEN))
+_factorback(GEN fa, GEN e, GEN (*_mul)(void*,GEN,GEN), GEN (*_pow)(void*,GEN,GEN), void *data)
 {
   pari_sp av = avma;
   long k,l,lx,t = typ(fa);
@@ -2383,7 +2394,7 @@
       if (l != 3) err(talker,"not a factorisation in factorback");
     } else {
       if (!is_vec_t(t)) err(talker,"not a factorisation in factorback");
-      return gerepileupto(av, divide_conquer_prod(fa, _mul));
+      return gerepileupto(av, divide_conquer_assoc(fa, _mul,data));
     }
     p = gel(fa,1);
     e = gel(fa,2);
@@ -2402,11 +2413,14 @@
   x = cgetg(lx,t_VEC);
   for (l=1,k=1; k<lx; k++)
     if (signe(e[k]))
-      gel(x,l++) = _pow(gel(p,k),gel(e,k));
+      gel(x,l++) = _pow(data,gel(p,k),gel(e,k));
   setlg(x,l);
-  return gerepileupto(av, divide_conquer_prod(x, _mul));
+  return gerepileupto(av, divide_conquer_assoc(x, _mul,data));
 }
 
+static GEN _agmul(void *a, GEN x, GEN y) { return gmul(x,y);}
+static GEN _apowgi(void *a, GEN x, GEN y) { return powgi(x,y);}
+
 GEN
 factorback_i(GEN fa, GEN e, GEN OBJ, int red)
 {
@@ -2415,11 +2429,10 @@
     if (e) {
       OBJ = _checknf(e); if (OBJ) e = NULL;
     }
-    if (!OBJ) return _factorback(fa, e, &gmul, &powgi);
+    if (!OBJ) return _factorback(fa, e, &_agmul, &_apowgi, NULL);
   }
-  static_OBJ = OBJ;
-  if (red) return _factorback(fa, e, &idmulred, &idpowred);
-  else     return _factorback(fa, e, &idmul,    &idpow);
+  if (red) return _factorback(fa, e, &idmulred, &idpowred, OBJ);
+  else     return _factorback(fa, e, &idmul,    &idpow, OBJ);
 }
 
 GEN
@@ -2428,8 +2441,8 @@
   if (!nf && e && lg(e) > 1 && typ(e[1]) != t_INT) { nf = e; e = NULL; }
   if (!nf) err(talker, "missing nf in factorbackelt");
 
-  static_OBJ = checknf(nf);
-  return _factorback(fa, e, &eltmul, &eltpow);
+  nf = checknf(nf);
+  return _factorback(fa, e, &eltmul, &eltpow, nf);
 }
 
 GEN
Index: pari/src/basemath/polarit3.c
===================================================================
--- pari.orig/src/basemath/polarit3.c	2005-08-30 16:18:31.000000000 +0200
+++ pari/src/basemath/polarit3.c	2005-08-30 23:23:21.000000000 +0200
@@ -669,13 +669,11 @@
   return gerepileupto(av, y);
 }
 
-static GEN modulo;
-static GEN _FpX_mul(GEN a,GEN b){return FpX_mul(a,b,modulo);}
+static GEN _FpX_mul(void *p,GEN a,GEN b){return FpX_mul(a,b,(GEN)p);}
 GEN 
 FpXV_prod(GEN V, GEN p)
 {
-  modulo = p; 
-  return divide_conquer_prod(V, &_FpX_mul);
+  return divide_conquer_assoc(V, &_FpX_mul,(void *)p);
 }
 
 GEN
@@ -1121,8 +1119,12 @@
   return T? FpXQX_divrem(x,y,T,p,z): FpX_divrem(x,y,p,z);
 }
 
-static GEN Tmodulo;
-static GEN _FpXQX_mul(GEN a,GEN b){return FpXQX_mul(a,b,Tmodulo,modulo);}
+struct _FpXQX { GEN T,p; };
+static GEN _FpXQX_mul(void *data, GEN a,GEN b)
+{
+  struct _FpXQX *d=(struct _FpXQX*)data;
+  return FpXQX_mul(a,b,d->T,d->p);
+}
 GEN 
 FpXQXV_prod(GEN V, GEN T, GEN p)
 {
@@ -1135,8 +1137,13 @@
     Tl = FlxqXV_prod(Vl, Tl, pp);
     return gerepileupto(av, FlxX_to_ZXX(Tl));
   }
-  modulo = p; Tmodulo = T;
-  return divide_conquer_prod(V, &_FpXQX_mul);
+  else
+  {
+    struct _FpXQX d;
+    d.p=p; 
+    d.T=T;
+    return divide_conquer_assoc(V, &_FpXQX_mul,(void*)&d);
+  }
 }
 
 GEN
Index: pari/src/headers/paridecl.h
===================================================================
--- pari.orig/src/headers/paridecl.h	2005-08-30 19:52:03.000000000 +0200
+++ pari/src/headers/paridecl.h	2005-08-30 23:23:21.000000000 +0200
@@ -1491,6 +1491,7 @@
 GEN     deg1_from_roots(GEN L, long v);
 GEN     discsr(GEN x);
 GEN     divide_conquer_prod(GEN x, GEN (*mul)(GEN,GEN));
+GEN     divide_conquer_assoc(GEN x, GEN (*mul)(void*,GEN,GEN), void *data);
 GEN     factor(GEN x);
 GEN     factor0(GEN x,long flag);
 GEN     factorback(GEN fa,GEN nf);