Re: two bugs in 2.2.0

[Karim BELABAS <Karim.Belabas@math.u-psud.fr>]

On Sun, 23 Dec 2001, Ilya Zakharevich wrote:
[...]
> is the following working as expected?  I do not think that
> such simple polynomials would require something as drastic as this:
>
>   ones(N) = (t^N - 1)/(t - 1)
>   allocatemem(132*1024^2)
>   p = 3456; ones(3*p-1) % ones(p)
>
> runs out of memory (and is quite slow if not).

OK, I made the fix I mentionned in my previous email plus a trivial one to
take advantage of polynomials with few non-zero entries [ benchmarks said the
penalty was negligible for small degrees ]

On my old Ultra, I now get [with default stack]

(18:49) gp > ones(N) = (t^N - 1)/(t - 1)
(18:49) gp > p = 3456; ones(3*p-1) % ones(p);
time = 720 ms.

which improves a bit on the situation for largish degrees.

Cheers,

    Karim.

P.S: CVS server for the development version updated.

P.S2: I've attached a patch for pari-2.2.0 for which the original bug report
was made [ but if you're using development versions anyway, you should rather
upgrade: 2.2.0 is about 8 months old ]
-- 
Karim Belabas                    Tel: (+33) (0)1 69 15 57 48
Dép. de Mathematiques, Bat. 425  Fax: (+33) (0)1 69 15 60 19
Université Paris-Sud             Email: Karim.Belabas@math.u-psud.fr
F-91405 Orsay (France)           http://www.math.u-psud.fr/~belabas
--
PARI/GP Home Page: http://www.parigp-home.de/

RCS file: /home/megrez/cvsroot/pari/src/basemath/polarit1.c,v
retrieving revision 1.47
retrieving revision 1.77
diff -c -r1.47 -r1.77
*** src/basemath/polarit1.c	2001/04/13 09:38:19	1.47
--- src/basemath/polarit1.c	2002/01/08 17:47:14	1.77
***************
*** 207,213 ****
  {
    ulong avy,av,av1;
    long ty=typ(y),tx,vx,vy,dx,dy,dz,i,j,sx,lrem;
-   int remainder;
    GEN z,p1,rem,y_lead,mod;
    GEN (*f)(GEN,GEN);
  
--- 118,123 ----
***************
*** 294,308 ****
    z[1]=evalsigne(1) | evallgef(dz+3) | evalvarn(vx);
    x += 2; y += 2; z += 2;
  
!   p1 = (GEN)x[dx]; remainder = (pr == ONLY_REM);
    z[dz]=y_lead? (long)f(p1,y_lead): lcopy(p1);
    for (i=dx-1; i>=dy; i--)
    {
      av1=avma; p1=(GEN)x[i];
      for (j=i-dy+1; j<=i && j<=dz; j++)
!       if (y[i-j]) p1 = gadd(p1, gmul((GEN)z[j],(GEN)y[i-j]));
      if (y_lead) p1 = f(p1,y_lead);
!     if (!remainder) p1 = avma==av1? gcopy(p1): gerepileupto(av1,p1);
      z[i-dy] = (long)p1;
    }
    if (!pr) return gerepileupto(av,z-2);
--- 204,221 ----
    z[1]=evalsigne(1) | evallgef(dz+3) | evalvarn(vx);
    x += 2; y += 2; z += 2;
  
!   p1 = (GEN)x[dx];
    z[dz]=y_lead? (long)f(p1,y_lead): lcopy(p1);
    for (i=dx-1; i>=dy; i--)
    {
      av1=avma; p1=(GEN)x[i];
      for (j=i-dy+1; j<=i && j<=dz; j++)
!       if (y[i-j] && z[j] != zero) p1 = gadd(p1, gmul((GEN)z[j],(GEN)y[i-j]));
      if (y_lead) p1 = f(p1,y_lead);
!    
!     if (isexactzero(p1)) { avma=av1; p1 = gzero; }
!     else
!       p1 = avma==av1? gcopy(p1): gerepileupto(av1,p1);
      z[i-dy] = (long)p1;
    }
    if (!pr) return gerepileupto(av,z-2);
***************
*** 313,319 ****
      p1 = (GEN)x[i];
      /* we always enter this loop at least once */
      for (j=0; j<=i && j<=dz; j++)
!       if (y[i-j]) p1 = gadd(p1, gmul((GEN)z[j],(GEN)y[i-j]));
      if (mod && avma==av1) p1 = gmul(p1,mod);
      if (!gcmp0(p1)) { sx = 1; break; } /* remainder is non-zero */
      if (!isinexactreal(p1) && !isexactzero(p1)) break;
--- 226,232 ----
      p1 = (GEN)x[i];
      /* we always enter this loop at least once */
      for (j=0; j<=i && j<=dz; j++)
!       if (y[i-j] && z[j] != zero) p1 = gadd(p1, gmul((GEN)z[j],(GEN)y[i-j]));
      if (mod && avma==av1) p1 = gmul(p1,mod);
      if (!gcmp0(p1)) { sx = 1; break; } /* remainder is non-zero */
      if (!isinexactreal(p1) && !isexactzero(p1)) break;
***************
*** 337,349 ****
    {
      av1=avma; p1 = (GEN)x[i];
      for (j=0; j<=i && j<=dz; j++)
!       if (y[i-j]) p1 = gadd(p1, gmul((GEN)z[j],(GEN)y[i-j]));
      if (mod && avma==av1) p1 = gmul(p1,mod);
      rem[i]=avma==av1? lcopy(p1):lpileupto(av1,p1);
    }
    rem -= 2;
    if (!sx) normalizepol_i(rem, lrem);
!   if (remainder) return gerepileupto(av,rem);
    z -= 2;
    {
      GEN *gptr[2]; gptr[0]=&z; gptr[1]=&rem;
--- 250,262 ----
    {
      av1=avma; p1 = (GEN)x[i];
      for (j=0; j<=i && j<=dz; j++)
!       if (y[i-j] && z[j] != zero) p1 = gadd(p1, gmul((GEN)z[j],(GEN)y[i-j]));
      if (mod && avma==av1) p1 = gmul(p1,mod);
      rem[i]=avma==av1? lcopy(p1):lpileupto(av1,p1);
    }
    rem -= 2;
    if (!sx) normalizepol_i(rem, lrem);
!   if (pr == ONLY_REM) return gerepileupto(av,rem);
    z -= 2;
    {
      GEN *gptr[2]; gptr[0]=&z; gptr[1]=&rem;