=================================================================== RCS file: /home/cvs/OpenXM_contrib/pari-2.2/src/basemath/Attic/gen1.c,v retrieving revision 1.1 retrieving revision 1.2 diff -u -p -r1.1 -r1.2 --- OpenXM_contrib/pari-2.2/src/basemath/Attic/gen1.c 2001/10/02 11:17:04 1.1 +++ OpenXM_contrib/pari-2.2/src/basemath/Attic/gen1.c 2002/09/11 07:26:50 1.2 @@ -1,4 +1,4 @@ -/* $Id: gen1.c,v 1.1 2001/10/02 11:17:04 noro Exp $ +/* $Id: gen1.c,v 1.2 2002/09/11 07:26:50 noro Exp $ Copyright (C) 2000 The PARI group. @@ -30,13 +30,13 @@ Foundation, Inc., 59 Temple Place - Suite 330, Boston, /* assume z[1] was created last */ #define fix_frac_if_int(z) if (is_pm1(z[2]))\ - z = gerepileupto((long)(z+3), (GEN)z[1]); + z = gerepileupto((gpmem_t)(z+3), (GEN)z[1]); /* assume z[1] was created last */ #define fix_frac_if_int_GC(z,tetpil) { if (is_pm1(z[2]))\ - z = gerepileupto((long)(z+3), (GEN)z[1]);\ + z = gerepileupto((gpmem_t)(z+3), (GEN)z[1]);\ else\ - gerepilemanyvec((long)z, tetpil, z+1, 2); } + gerepilemanyvec((gpmem_t)z, tetpil, z+1, 2); } GEN quickmul(GEN a, GEN b, long na, long nb); @@ -46,7 +46,7 @@ static GEN op_polmod(GEN f(GEN,GEN), GEN x, GEN y, long tx) { GEN mod,k,l, z=cgetg(3,t_POLMOD); - long av,tetpil; + gpmem_t av, tetpil; l=(GEN)y[1]; if (tx==t_POLMOD) @@ -105,7 +105,7 @@ gred_rfrac_simple(GEN x1, GEN x2) c = denom(x1); y = cgetg(3,t_RFRAC); - y[1] = (long)numer(x1); + y[1] = lmul(x1,c); y[2] = lmul(x2,c); return y; } @@ -116,8 +116,8 @@ gred_rfrac2_i(GEN x1, GEN x2) long tx,ty; if (gcmp0(x1)) return gcopy(x1); - - tx=typ(x1); ty=typ(x2); + x1 = simplify_i(x1); tx = typ(x1); + x2 = simplify_i(x2); ty = typ(x2); if (ty!=t_POL) { if (tx!=t_POL) return gred_rfrac_copy(x1,x2); @@ -165,7 +165,7 @@ gred_rfrac_i(GEN x) GEN gred_rfrac2(GEN x1, GEN x2) { - ulong av = avma; + gpmem_t av = avma; return gerepileupto(av, gred_rfrac2_i(x1, x2)); } @@ -180,7 +180,7 @@ GEN gred_frac2(GEN x1, GEN x2) { GEN p1, y = dvmdii(x1,x2,&p1); - ulong av; + gpmem_t av; if (p1 == gzero) return y; /* gzero intended */ av = avma; @@ -221,7 +221,7 @@ gred(GEN x) GEN gsub(GEN x, GEN y) { - long tetpil, av = avma; + gpmem_t tetpil, av = avma; y=gneg_i(y); tetpil=avma; return gerepile(av,tetpil,gadd(x,y)); } @@ -235,90 +235,108 @@ gsub(GEN x, GEN y) static GEN addpadic(GEN x, GEN y) { - long c,e,r,d,r1,r2,av,tetpil; - GEN z,p1,p2, p = (GEN)x[2]; + gpmem_t av = avma; + long c,d,e,r,rx,ry; + GEN u,z,p,mod; - z=cgetg(5,t_PADIC); icopyifstack(p, z[2]); av=avma; - e=valp(x); r=valp(y); d = r-e; - if (d<0) { p1=x; x=y; y=p1; e=r; d = -d; } - r1=precp(x); r2=precp(y); - if (d) + (void)new_chunk(5+lgefint(x[3])+lgefint(y[3])); + e = valp(x); + r = valp(y); d = r-e; + if (d < 0) { GEN p1=x; x=y; y=p1; e = r; d = -d; } + rx = precp(x); p = (GEN)x[2]; + ry = precp(y); + if (d) /* v(x) < v(y) */ { - r = d+r2; - p1 = (d==1)? p: gclone(gpuigs(p,d)); - avma=av; - if (r=r) + else { - avma=av; z[4]=zero; z[3]=un; - z[1]=evalvalp(e+r); return z; - } - if (c) - { - p2=gclone(p2); avma=av; - if (c==1) - z[3] = ldivii((GEN)x[3], p); - else + if (ry < rx) { r=ry; mod = (GEN)x[3]; } else { r=rx; mod = (GEN)y[3]; } + u = addii((GEN)x[4], (GEN)y[4]); + if (!signe(u) || (c = pvaluation(u,p,&u)) >= r) { - p1 = gpuigs(p,c); tetpil=avma; - z[3] = lpile(av,tetpil, divii((GEN)x[3], p1)); + avma = av; return padiczero(p, e+r); } - z[4]=lmodii(p2,(GEN)z[3]); gunclone(p2); - z[1]=evalprecp(r-c) | evalvalp(e+c); return z; + if (c) + { + mod = divii(mod, gpowgs(p,c)); + r -= c; + e += c; + } } - tetpil=avma; - z[4]=lpile(av,tetpil,modii(p1,(GEN)x[3])); - z[3]=licopy((GEN)x[3]); - z[1]=evalprecp(r) | evalvalp(e); return z; + avma = av; z = cgetg(5,t_PADIC); + z[1] = evalprecp(r) | evalvalp(e); + z[3] = licopy(mod); + z[4] = lmodii(u,(GEN)z[3]); + icopyifstack(p, z[2]); return z; } /* return x + y, where x is t_INT or t_FRAC(N), y t_PADIC */ static GEN gaddpex(GEN x, GEN y) { - long tx,e1,e2,e3,av,tetpil; - GEN z,p,p1,p2; + gpmem_t av; + long tx,vy,py,d,r,e; + GEN z,q,p,p1,p2,mod,u; if (gcmp0(x)) return gcopy(y); - av=avma; p=(GEN)y[2]; tx=typ(x); - z=cgetg(5,t_PADIC); z[2]=(long)p; - e3 = (tx == t_INT)? pvaluation(x,p,&p1) + av = avma; p = (GEN)y[2]; tx = typ(x); + e = (tx == t_INT)? pvaluation(x,p,&p1) : pvaluation((GEN)x[1],p,&p1) - pvaluation((GEN)x[2],p,&p2); - e1 = valp(y)-e3; e2 = signe(y[4])? e1+precp(y): e1; - if (e2<=0) + vy = valp(y); d = vy - e; py = precp(y); r = d + py; + if (r <= 0) { avma = av; return gcopy(y); } + mod = (GEN)y[3]; + u = (GEN)y[4]; + (void)new_chunk(5 + lgefint(mod) + lgefint(p)*labs(d)); + + if (d > 0) { - z[1] = evalprecp(0) | evalvalp(e3); - z[3] = un; - z[4] = zero; + q = gpowgs(p,d); + mod = mulii(mod, q); + u = mulii(u, q); + if (tx != t_INT && !is_pm1(p2)) p1 = mulii(p1, mpinvmod(p2,mod)); + u = addii(u, p1); } + else if (d < 0) + { + q = gpowgs(p,-d); + if (tx != t_INT && !is_pm1(p2)) p1 = mulii(p1, mpinvmod(p2,mod)); + p1 = mulii(p1, q); + u = addii(u, p1); + r = py; e = vy; + } else { - if (tx != t_INT && !is_pm1(p2)) p1 = gdiv(p1,p2); - z[1] = evalprecp(e2) | evalvalp(e3); - z[3] = e1? lmul((GEN)y[3], gpuigs(p,e1)): y[3]; - z[4] = lmod(p1,(GEN)z[3]); + long c; + if (tx != t_INT && !is_pm1(p2)) p1 = mulii(p1, mpinvmod(p2,mod)); + u = addii(u, p1); + if (!signe(u) || (c = pvaluation(u,p,&u)) >= r) + { + avma = av; return padiczero(p,e+r); + } + if (c) + { + mod = divii(mod, gpowgs(p,c)); + r -= c; + e += c; + } } - tetpil=avma; return gerepile(av,tetpil,addpadic(z,y)); + avma = av; z = cgetg(5,t_PADIC); + z[1] = evalprecp(r) | evalvalp(e); + z[3] = licopy(mod); + z[4] = lmodii(u,(GEN)z[3]); + icopyifstack(p, z[2]); return z; } static long kro_quad(GEN x, GEN y) { - long k, av=avma; + long k; + gpmem_t av=avma; x = subii(sqri((GEN)x[3]), shifti((GEN)x[2],2)); k = kronecker(x,y); avma=av; return k; @@ -337,20 +355,20 @@ addfrac(GEN x, GEN y) if (is_pm1(delta)) { p1 = mulii(x1,y2); - p2 = mulii(y1,x2); avma = (long)z; + p2 = mulii(y1,x2); avma = (gpmem_t)z; z[1] = laddii(p1,p2); z[2] = lmulii(x2,y2); return z; } x2 = divii(x2,delta); y2 = divii(y2,delta); n = addii(mulii(x1,y2), mulii(y1,x2)); - if (!signe(n)) { avma = (long)(z+3); return gzero; } + if (!signe(n)) { avma = (gpmem_t)(z+3); return gzero; } d = mulii(x2, y2); p1 = dvmdii(n, delta, &p2); if (p2 == gzero) { - if (is_pm1(d)) { avma = (long)(z+3); return icopy(p1); } - avma = (long)z; + if (is_pm1(d)) { avma = (gpmem_t)(z+3); return icopy(p1); } + avma = (gpmem_t)z; z[1] = licopy(p1); z[2] = licopy(d); return z; } @@ -361,7 +379,7 @@ addfrac(GEN x, GEN y) n = divii(n, p1); } d = mulii(d,delta); - avma = (long)z; + avma = (gpmem_t)z; z[1] = licopy(n); z[2] = licopy(d); return z; } @@ -372,7 +390,7 @@ addrfrac(GEN x, GEN y) GEN z = cgetg(3,t_RFRAC); GEN x1 = (GEN)x[1], x2 = (GEN)x[2]; GEN y1 = (GEN)y[1], y2 = (GEN)y[2], p1,p2,n,d,delta; - long tetpil; + gpmem_t tetpil; delta = ggcd(x2,y2); if (gcmp1(delta)) @@ -380,26 +398,26 @@ addrfrac(GEN x, GEN y) p1 = gmul(x1,y2); p2 = gmul(y1,x2); tetpil = avma; /* numerator is non-zero */ - z[1] = lpile((long)z,tetpil, gadd(p1,p2)); + z[1] = lpile((gpmem_t)z,tetpil, gadd(p1,p2)); z[2] = lmul(x2, y2); return z; } x2 = gdeuc(x2,delta); y2 = gdeuc(y2,delta); n = gadd(gmul(x1,y2), gmul(y1,x2)); - if (!signe(n)) return gerepileupto((long)(z+3), n); + if (gcmp0(n)) return gerepileupto((gpmem_t)(z+3), n); tetpil = avma; d = gmul(x2, y2); p1 = poldivres(n, delta, &p2); /* we want gcd(n,delta) */ - if (!signe(p2)) + if (gcmp0(p2)) { if (lgef(d) == 3) /* "constant" denominator */ { d = (GEN)d[2]; if (gcmp_1(d)) p1 = gneg(p1); else if (!gcmp1(d)) p1 = gdiv(p1, d); - return gerepileupto((long)(z+3), p1); + return gerepileupto((gpmem_t)(z+3), p1); } z[1]=(long)p1; z[2]=(long)d; - gerepilemanyvec((long)z,tetpil,z+1,2); return z; + gerepilemanyvec((gpmem_t)z,tetpil,z+1,2); return z; } p1 = ggcd(delta, p2); if (gcmp1(p1)) @@ -414,14 +432,14 @@ addrfrac(GEN x, GEN y) z[1] = ldeuc(n,p1); } z[2] = lmul(d,delta); - gerepilemanyvec((long)z,tetpil,z+1,2); return z; + gerepilemanyvec((gpmem_t)z,tetpil,z+1,2); return z; } static GEN addscalrfrac(GEN x, GEN y) { GEN p1,num, z = cgetg(3,t_RFRAC); - long tetpil, av; + gpmem_t tetpil, av; p1 = gmul(x,(GEN)y[2]); tetpil = avma; num = gadd(p1,(GEN)y[1]); @@ -435,11 +453,11 @@ addscalrfrac(GEN x, GEN y) tetpil = avma; z[1] = ldiv(num, p1); z[2] = ldiv((GEN)y[2], p1); - gerepilemanyvec((long)z,tetpil,z+1,2); return z; + gerepilemanyvec((gpmem_t)z,tetpil,z+1,2); return z; } } avma = av; - z[1]=lpile((long)z,tetpil, num); + z[1]=lpile((gpmem_t)z,tetpil, num); z[2]=lcopy((GEN)y[2]); return z; } @@ -460,7 +478,8 @@ to_polmod(GEN x, GEN mod) GEN gadd(GEN x, GEN y) { - long tx = typ(x), ty = typ(y), vx,vy,lx,ly,i,j,k,l,av,tetpil; + long tx = typ(x), ty = typ(y), vx, vy, lx, ly, i, j, k, l; + gpmem_t av, tetpil; GEN z,p1,p2; if (is_const_t(tx) && is_const_t(ty)) @@ -476,13 +495,13 @@ gadd(GEN x, GEN y) case t_INTMOD: z=cgetg(3,t_INTMOD); p2=(GEN)y[1]; (void)new_chunk(lgefint(p2)+1); /* HACK */ - p1 = addii(modii(x,p2),(GEN)y[2]); avma = (long)z; + p1 = addii(modii(x,p2),(GEN)y[2]); avma = (gpmem_t)z; z[2] = (cmpii(p1,p2) >=0)? lsubii(p1,p2): licopy(p1); icopyifstack(p2,z[1]); return z; case t_FRAC: case t_FRACN: z=cgetg(3,ty); (void)new_chunk(lgefint(x)+lgefint(y[1])+lgefint(y[2])+1); /*HACK*/ - p1 = mulii((GEN)y[2],x); avma = (long)z; + p1 = mulii((GEN)y[2],x); avma = (gpmem_t)z; z[1] = laddii((GEN)y[1], p1); z[2] = licopy((GEN)y[2]); return z; @@ -529,7 +548,7 @@ gadd(GEN x, GEN y) p1=co8(y,lg(x)+i); tetpil=avma; return gerepile(av,tetpil,gadd(p1,x)); - case t_INTMOD: case t_PADIC: err(operf,"+",tx,ty); + case t_INTMOD: case t_PADIC: err(operf,"+",x,y); } case t_INTMOD: @@ -554,7 +573,7 @@ gadd(GEN x, GEN y) case t_FRAC: case t_FRACN: z=cgetg(3,t_INTMOD); p2=(GEN)x[1]; (void)new_chunk(lgefint(p2)<<2); /* HACK */ p1 = mulii((GEN)y[1], mpinvmod((GEN)y[2],p2)); - p1 = addii(modii(p1,p2), (GEN)x[2]); avma=(long)z; + p1 = addii(modii(p1,p2), (GEN)x[2]); avma=(gpmem_t)z; z[2]=lmodii(p1,p2); icopyifstack(p2,z[1]); return z; case t_COMPLEX: z=cgetg(3,t_COMPLEX); @@ -562,10 +581,10 @@ gadd(GEN x, GEN y) z[2]=lcopy((GEN)y[2]); return z; case t_PADIC: - l=avma; p1=cgetg(3,t_INTMOD); + av=avma; p1=cgetg(3,t_INTMOD); p1[1]=x[1]; p1[2]=lgeti(lgefint(x[1])); gaffect(y,p1); tetpil=avma; - return gerepile(l,tetpil,gadd(p1,x)); + return gerepile(av,tetpil,gadd(p1,x)); case t_QUAD: z=cgetg(4,t_QUAD); copyifstack(y[1], z[1]); @@ -577,10 +596,10 @@ gadd(GEN x, GEN y) switch (ty) { case t_FRAC: return addfrac(x,y); - case t_FRACN: z=cgetg(3,t_FRACN); l=avma; + case t_FRACN: z=cgetg(3,t_FRACN); av=avma; p1=mulii((GEN)x[1],(GEN)y[2]); p2=mulii((GEN)x[2],(GEN)y[1]); - tetpil=avma; z[1]=lpile(l,tetpil,addii(p1,p2)); + tetpil=avma; z[1]=lpile(av,tetpil,addii(p1,p2)); z[2]=lmulii((GEN)x[2],(GEN)y[2]); return z; @@ -616,7 +635,7 @@ gadd(GEN x, GEN y) return gerepile(av,tetpil,gadd(p1,y)); case t_QUAD: - lx=precision(x); if (!lx) err(operi,"+",tx,ty); + lx=precision(x); if (!lx) err(operi,"+",x,y); if (gcmp0(y)) return gcopy(x); av=avma; i=gexpo(y)-gexpo(x); @@ -629,7 +648,7 @@ gadd(GEN x, GEN y) switch(ty) { case t_PADIC: - if (!egalii((GEN)x[2],(GEN)y[2])) err(operi,"+",tx,ty); + if (!egalii((GEN)x[2],(GEN)y[2])) err(operi,"+",x,y); return addpadic(x,y); case t_QUAD: @@ -646,7 +665,7 @@ gadd(GEN x, GEN y) } case t_QUAD: z=cgetg(4,t_QUAD); k=x[1]; l=y[1]; - if (!gegal((GEN)k,(GEN)l)) err(operi,"+",tx,ty); + if (!gegal((GEN)k,(GEN)l)) err(operi,"+",x,y); copyifstack(l, z[1]); z[2]=ladd((GEN)x[2],(GEN)y[2]); z[3]=ladd((GEN)x[3],(GEN)y[3]); return z; @@ -733,9 +752,9 @@ gadd(GEN x, GEN y) if (isexactzero(x)) return gcopy(y); if (ty == t_MAT) return gaddmat(x,y); /* fall through */ - case t_QFR: case t_QFI: err(operf,"+",tx,ty); + case t_QFR: case t_QFI: err(operf,"+",x,y); } - err(operf,"+",tx,ty); + err(operf,"+",x,y); } /* here !isscalar(x) && isscalar(y) && (vx=vy || ismatvec(x and y)) */ @@ -751,7 +770,7 @@ gadd(GEN x, GEN y) for (i=2; i varn(x)) return gdiv(x,y); } else if (varn(y) > varn(x)) return gdiv(x,y); - return NULL; + avma = av; return NULL; } static long @@ -895,7 +916,8 @@ GEN mulscalrfrac(GEN x, GEN y) { GEN p1,z,y1,y2,cx,cy1,cy2; - long tetpil,tx; + long tx; + gpmem_t tetpil; if (gcmp0(x)) return gcopy(x); @@ -926,8 +948,8 @@ mulscalrfrac(GEN x, GEN y) z[2] = lmul(y2, cy2); z[1] = lmul(y1, cy1); p1 = fix_rfrac_if_pol((GEN)z[1],(GEN)z[2]); - if (p1) return gerepileupto((long)(z+3), p1); - gerepilemanyvec((long)z,tetpil,z+1,2); return z; + if (p1) return gerepileupto((gpmem_t)(z+3), p1); + gerepilemanyvec((gpmem_t)z,tetpil,z+1,2); return z; } static GEN @@ -936,7 +958,7 @@ mulrfrac(GEN x, GEN y) GEN z = cgetg(3,t_RFRAC), p1; GEN x1 = (GEN)x[1], x2 = (GEN)x[2]; GEN y1 = (GEN)y[1], y2 = (GEN)y[2]; - long tetpil; + gpmem_t tetpil; p1 = ggcd(x1, y2); if (!gcmp1(p1)) { x1 = gdiv(x1,p1); y2 = gdiv(y2,p1); } p1 = ggcd(x2, y1); if (!gcmp1(p1)) { x2 = gdiv(x2,p1); y1 = gdiv(y1,p1); } @@ -944,8 +966,8 @@ mulrfrac(GEN x, GEN y) z[2] = lmul(x2,y2); z[1] = lmul(x1,y1); p1 = fix_rfrac_if_pol((GEN)z[1],(GEN)z[2]); - if (p1) return gerepileupto((long)(z+3), p1); - gerepilemanyvec((long)z,tetpil,z+1,2); return z; + if (p1) return gerepileupto((gpmem_t)(z+3), p1); + gerepilemanyvec((gpmem_t)z,tetpil,z+1,2); return z; } GEN @@ -1040,18 +1062,24 @@ gmul_err(GEN x, GEN y, long tx, long ty) case t_QFR: return compreal(x,y); case t_VECSMALL: l = lg(x); z = cgetg(l, t_VECSMALL); - if (l != lg(y)) err(operf,"*",tx,ty); - for (i=1; i=l) err(operf,"*",x,y); + z[i]=x[y[i]]; + } return z; } - err(operf,"*",tx,ty); + err(operf,"*",x,y); return NULL; /* not reached */ } GEN gmul(GEN x, GEN y) { - long tx,ty,lx,ly,vx,vy,i,j,k,l,av,tetpil; + long tx, ty, lx, ly, vx, vy, i, j, k, l; + gpmem_t av, tetpil; GEN z,p1,p2,p3,p4; if (x == y) return gsqr(x); @@ -1071,7 +1099,7 @@ gmul(GEN x, GEN y) case t_INTMOD: z=cgetg(3,t_INTMOD); p2=(GEN)y[1]; (void)new_chunk(lgefint(p2)<<2); /* HACK */ - p1=mulii(modii(x,p2),(GEN)y[2]); avma=(long)z; + p1=mulii(modii(x,p2),(GEN)y[2]); avma=(gpmem_t)z; z[2]=lmodii(p1,p2); icopyifstack(p2,z[1]); return z; case t_FRAC: @@ -1080,7 +1108,7 @@ gmul(GEN x, GEN y) p1 = mppgcd(x,(GEN)y[2]); if (is_pm1(p1)) { - avma = (long)z; + avma = (gpmem_t)z; z[2] = licopy((GEN)y[2]); z[1] = lmulii((GEN)y[1], x); } @@ -1103,8 +1131,8 @@ gmul(GEN x, GEN y) case t_PADIC: if (!signe(x)) return gzero; - l=avma; p1=cgetp(y); gaffect(x,p1); tetpil=avma; - return gerepile(l,tetpil,gmul(p1,y)); + av=avma; p1=cgetp(y); gaffect(x,p1); tetpil=avma; + return gerepile(av,tetpil,gmul(p1,y)); case t_QUAD: z=cgetg(4,t_QUAD); copyifstack(y[1], z[1]); @@ -1118,18 +1146,18 @@ gmul(GEN x, GEN y) case t_REAL: return mulrr(x,y); case t_FRAC: case t_FRACN: - l=avma; p1=cgetr(lg(x)); tetpil=avma; gaffect(y,p1); - p2=mulrr(p1,x); return gerepile(l,tetpil,p2); + av=avma; p1=mulri(x,(GEN)y[1]); tetpil=avma; + return gerepile(av, tetpil, divri(p1, (GEN)y[2])); case t_COMPLEX: z=cgetg(3,t_COMPLEX); z[1]=lmul(x,(GEN)y[1]); z[2]=lmul(x,(GEN)y[2]); return z; case t_QUAD: - l=avma; p1=co8(y,lg(x)); tetpil=avma; - return gerepile(l,tetpil,gmul(p1,x)); + av=avma; p1=co8(y,lg(x)); tetpil=avma; + return gerepile(av,tetpil,gmul(p1,x)); - default: err(operf,"*",tx,ty); + default: err(operf,"*",x,y); } case t_INTMOD: @@ -1155,7 +1183,7 @@ gmul(GEN x, GEN y) case t_FRAC: case t_FRACN: z=cgetg(3,t_INTMOD); p2=(GEN)x[1]; (void)new_chunk(lgefint(p2)<<2); /* HACK */ p1 = mulii((GEN)y[1], mpinvmod((GEN)y[2],p2)); - p1 = mulii(modii(p1,p2),(GEN)x[2]); avma=(long)z; + p1 = mulii(modii(p1,p2),(GEN)x[2]); avma=(gpmem_t)z; z[2] = lmodii(p1,p2); icopyifstack(p2,z[1]); return z; case t_COMPLEX: z=cgetg(3,t_COMPLEX); @@ -1163,10 +1191,10 @@ gmul(GEN x, GEN y) z[2]=lmul(x,(GEN)y[2]); return z; case t_PADIC: - l=avma; p1=cgetg(3,t_INTMOD); + av=avma; p1=cgetg(3,t_INTMOD); p1[1]=x[1]; p1[2]=lgeti(lg(x[1])); gaffect(y,p1); tetpil=avma; - return gerepile(l,tetpil,gmul(x,p1)); + return gerepile(av,tetpil,gmul(x,p1)); case t_QUAD: z=cgetg(4,t_QUAD); copyifstack(y[1], z[1]); @@ -1201,8 +1229,8 @@ gmul(GEN x, GEN y) case t_PADIC: if (!signe(x[1])) return gzero; - l=avma; p1=cgetp(y); gaffect(x,p1); tetpil=avma; - return gerepile(l,tetpil,gmul(p1,y)); + av=avma; p1=cgetp(y); gaffect(x,p1); tetpil=avma; + return gerepile(av,tetpil,gmul(p1,y)); case t_QUAD: z=cgetg(4,t_QUAD); copyifstack(y[1], z[1]); @@ -1213,14 +1241,14 @@ gmul(GEN x, GEN y) case t_COMPLEX: switch(ty) { - case t_COMPLEX: z=cgetg(3,t_COMPLEX); l=avma; + case t_COMPLEX: z=cgetg(3,t_COMPLEX); av=avma; p1=gmul((GEN)x[1],(GEN)y[1]); p2=gmul((GEN)x[2],(GEN)y[2]); x=gadd((GEN)x[1],(GEN)x[2]); y=gadd((GEN)y[1],(GEN)y[2]); y=gmul(x,y); x=gadd(p1,p2); tetpil=avma; z[1]=lsub(p1,p2); z[2]=lsub(y,x); - gerepilemanyvec(l,tetpil,z+1,2); return z; + gerepilemanyvec(av,tetpil,z+1,2); return z; case t_PADIC: if (krosg(-1,(GEN)y[2])) @@ -1239,16 +1267,16 @@ gmul(GEN x, GEN y) return gerepile(av,tetpil,gmul(p1,y)); case t_QUAD: - lx=precision(x); if (!lx) err(operi,"*",tx,ty); - l=avma; p1=co8(y,lx); tetpil=avma; - return gerepile(l,tetpil,gmul(p1,x)); + lx=precision(x); if (!lx) err(operi,"*",x,y); + av=avma; p1=co8(y,lx); tetpil=avma; + return gerepile(av,tetpil,gmul(p1,x)); } case t_PADIC: switch(ty) { case t_PADIC: - if (!egalii((GEN)x[2],(GEN)y[2])) err(operi,"*",tx,ty); + if (!egalii((GEN)x[2],(GEN)y[2])) err(operi,"*",x,y); l = valp(x)+valp(y); if (!signe(x[4])) { z=gcopy(x); setvalp(z,l); return z; } if (!signe(y[4])) { z=gcopy(y); setvalp(z,l); return z; } @@ -1273,9 +1301,9 @@ gmul(GEN x, GEN y) case t_QUAD: z=cgetg(4,t_QUAD); p1=(GEN)x[1]; p2=(GEN)y[1]; - if (!gegal(p1,p2)) err(operi,"*",tx,ty); + if (!gegal(p1,p2)) err(operi,"*",x,y); - copyifstack(p2, z[1]); l=avma; + copyifstack(p2, z[1]); av=avma; p2=gmul((GEN)x[2],(GEN)y[2]); p3=gmul((GEN)x[3],(GEN)y[3]); p4=gmul(gneg_i((GEN)p1[2]),p3); @@ -1283,17 +1311,17 @@ gmul(GEN x, GEN y) if (gcmp0((GEN)p1[3])) { tetpil=avma; - z[2]=lpile(l,tetpil,gadd(p4,p2)); l=avma; + z[2]=lpile(av,tetpil,gadd(p4,p2)); av=avma; p2=gmul((GEN)x[2],(GEN)y[3]); p3=gmul((GEN)x[3],(GEN)y[2]); tetpil=avma; - z[3]=lpile(l,tetpil,gadd(p2,p3)); return z; + z[3]=lpile(av,tetpil,gadd(p2,p3)); return z; } p1 = gadd(gmul((GEN)x[2],(GEN)y[3]), gmul((GEN)x[3],(GEN)y[2])); tetpil=avma; z[2]=ladd(p2,p4); z[3]=ladd(p1,p3); - gerepilemanyvec(l,tetpil,z+2,2); return z; + gerepilemanyvec(av,tetpil,z+2,2); return z; } err(bugparier,"multiplication"); } @@ -1327,18 +1355,18 @@ gmul(GEN x, GEN y) switch(ty) { case t_COL: - if (lx!=ly) err(operi,"*",tx,ty); - z=gzero; l=avma; + if (lx!=ly) err(operi,"*",x,y); + z=gzero; av=avma; for (i=1; ily) { k=ly; ly=lx; lx=k; p1=y; y=x; x=p1; } @@ -1536,31 +1566,35 @@ gmul(GEN x, GEN y) z[1] = evalvalp(valp(x)+valp(y)) | evalvarn(vx) | evalsigne(1); x += 2; y += 2; z += 2; lx -= 3; p2 = (GEN)gpmalloc((lx+1)*sizeof(long)); + p3 = (GEN)gpmalloc((ly+1)*sizeof(long)); + mix = miy = 0; for (i=0; i<=lx; i++) { - p2[i] = !isexactzero((GEN)y[i]); + p2[i] = !isexactzero((GEN)y[i]); if (p2[i]) miy = i; + p3[i] = !isexactzero((GEN)x[i]); if (p3[i]) mix = i; p1 = gzero; av = avma; - for (j=0; j<=i; j++) - if (p2[j]) - p1 = gadd(p1, gmul((GEN)y[j],(GEN)x[i-j])); + for (j=i-mix; j<=min(i,miy); j++) + if (p2[j]) p1 = gadd(p1, gmul((GEN)y[j],(GEN)x[i-j])); z[i] = lpileupto(av,p1); } z -= 2; /* back to normalcy */ - free(p2); return normalize(z); + free(p2); + free(p3); return normalize(z); + } case t_RFRAC: case t_RFRACN: if (gcmp0(y)) return zeropol(vx); if (gcmp0(x)) return zeroser(vx, valp(x)+gval(y,vx)); - l=avma; p1=gmul((GEN)y[1],x); tetpil=avma; - return gerepile(l,tetpil,gdiv(p1,(GEN)y[2])); + av=avma; p1=gmul((GEN)y[1],x); tetpil=avma; + return gerepile(av,tetpil,gdiv(p1,(GEN)y[2])); - default: err(operf,"*",tx,ty); + default: err(operf,"*",x,y); } /* (tx,ty) == t_RFRAC <==> ty == t_RFRAC */ case t_RFRAC: return mulrfrac(x,y); case t_RFRACN: - if (!is_rfrac_t(ty)) err(operf,"*",tx,ty); + if (!is_rfrac_t(ty)) err(operf,"*",x,y); av=avma; z=cgetg(3,ty); z[1]=lmul((GEN)x[1],(GEN)y[1]); z[2]=lmul((GEN)x[2],(GEN)y[2]); return z; @@ -1571,7 +1605,8 @@ gmul(GEN x, GEN y) GEN gsqr(GEN x) { - long tx=typ(x),lx,i,j,k,l,av,tetpil; + long tx=typ(x), lx, i, j, k, l; + gpmem_t av, tetpil; GEN z,p1,p2,p3,p4; if (is_scalar_t(tx)) @@ -1585,7 +1620,7 @@ gsqr(GEN x) case t_INTMOD: z=cgetg(3,t_INTMOD); p2=(GEN)x[1]; (void)new_chunk(lgefint(p2)<<2); /* HACK */ - p1=sqri((GEN)x[2]); avma=(long)z; + p1=sqri((GEN)x[2]); avma=(gpmem_t)z; z[2]=lmodii(p1,p2); icopyifstack(p2,z[1]); return z; case t_FRAC: case t_FRACN: @@ -1595,13 +1630,13 @@ gsqr(GEN x) return z; case t_COMPLEX: - z=cgetg(lg(x),tx); l=avma; + z=cgetg(3,t_COMPLEX); av=avma; p1=gadd((GEN)x[1],(GEN)x[2]); p2=gadd((GEN)x[1],gneg_i((GEN)x[2])); p3=gmul((GEN)x[1],(GEN)x[2]); tetpil=avma; z[1]=lmul(p1,p2); z[2]=lshift(p3,1); - gerepilemanyvec(l,tetpil,z+1,2); + gerepilemanyvec(av,tetpil,z+1,2); return z; case t_PADIC: @@ -1611,33 +1646,33 @@ gsqr(GEN x) z[1] = evalprecp(precp(x)+i) | evalvalp(2*valp(x)); icopyifstack(x[2], z[2]); z[3] = lshifti((GEN)x[3], i); av = avma; - z[4] = (long)gerepileuptoint(av, modii(sqri((GEN)x[4]), (GEN)z[3])); + z[4] = lpileuptoint(av, modii(sqri((GEN)x[4]), (GEN)z[3])); return z; case t_QUAD: - p1=(GEN)x[1]; z=cgetg(lg(x),tx); l=avma; + p1=(GEN)x[1]; z=cgetg(4,t_QUAD); av=avma; p2=gsqr((GEN)x[2]); p3=gsqr((GEN)x[3]); p4=gmul(gneg_i((GEN)p1[2]),p3); if (gcmp0((GEN)p1[3])) { tetpil=avma; - z[2]=lpile(l,tetpil,gadd(p4,p2)); - l=avma; p2=gmul((GEN)x[2],(GEN)x[3]); tetpil=avma; - z[3]=lpile(l,tetpil,gmul2n(p2,1)); + z[2]=lpile(av,tetpil,gadd(p4,p2)); + av=avma; p2=gmul((GEN)x[2],(GEN)x[3]); tetpil=avma; + z[3]=lpile(av,tetpil,gmul2n(p2,1)); copyifstack(p1,z[1]); return z; } p1=gmul((GEN)x[2],(GEN)x[3]); p1=gmul2n(p1,1); tetpil=avma; z[2]=ladd(p2,p4); z[3]=ladd(p1,p3); - gerepilemanyvec(l,tetpil,z+2,2); + gerepilemanyvec(av,tetpil,z+2,2); copyifstack(x[1],z[1]); return z; case t_POLMOD: - z=cgetg(lg(x),tx); copyifstack(x[1],z[1]); - l=avma; p1=gsqr((GEN)x[2]); tetpil=avma; - z[2]=lpile(l,tetpil, gres(p1,(GEN)z[1])); + z=cgetg(3,t_POLMOD); copyifstack(x[1],z[1]); + av=avma; p1=gsqr((GEN)x[2]); tetpil=avma; + z[2]=lpile(av,tetpil, gres(p1,(GEN)z[1])); return z; } @@ -1664,25 +1699,27 @@ gsqr(GEN x) } case t_SER: + { + long mi; if (gcmp0(x)) return zeroser(varn(x), 2*valp(x)); lx = lg(x); z = cgetg(lx,tx); z[1] = evalsigne(1) | evalvalp(2*valp(x)) | evalvarn(varn(x)); x += 2; z += 2; lx -= 3; p2 = (GEN)gpmalloc((lx+1)*sizeof(long)); + mi = 0; for (i=0; i<=lx; i++) { - p2[i] = !isexactzero((GEN)x[i]); - p1=gzero; av=avma; l=(i+1)>>1; - for (j=0; j>1) - 1; + for (j=i-mi; j<=min(l,mi); j++) + if (p2[j] && p2[i-j]) p1 = gadd(p1, gmul((GEN)x[j],(GEN)x[i-j])); p1 = gshift(p1,1); if ((i&1) == 0 && p2[i>>1]) p1 = gadd(p1, gsqr((GEN)x[i>>1])); z[i] = lpileupto(av,p1); } z -= 2; free(p2); return normalize(z); - + } case t_RFRAC: case t_RFRACN: z=cgetg(3,tx); z[1]=lsqr((GEN)x[1]); @@ -1691,17 +1728,17 @@ gsqr(GEN x) case t_MAT: lx=lg(x); if (lx==1) return cgetg(1,tx); - if (lx != lg(x[1])) err(operi,"*",tx,tx); + if (lx != lg(x[1])) err(operi,"*",x,x); z=cgetg(lx,tx); for (j=1; j0)?valp(x):1); - tetpil=avma; return gerepile(l,tetpil,gdiv(x,p1)); + tetpil=avma; return gerepile(av,tetpil,gdiv(x,p1)); case t_INTMOD: - l=avma; p1=cgetg(3,t_INTMOD); + av=avma; p1=cgetg(3,t_INTMOD); p1[1]=y[1]; p1[2]=lgeti(lg(y[1])); gaffect(x,p1); tetpil=avma; - return gerepile(l,tetpil,gdiv(p1,y)); + return gerepile(av,tetpil,gdiv(p1,y)); case t_PADIC: - if (!egalii((GEN)x[2],(GEN)y[2])) err(operi,"/",tx,ty); + if (!egalii((GEN)x[2],(GEN)y[2])) err(operi,"/",x,y); if (!signe(x[4])) { z=gcopy(x); setvalp(z,valp(x)-valp(y)); @@ -2020,18 +2058,18 @@ gdiv(GEN x, GEN y) } p1=(precp(x)>precp(y)) ? y : x; - z=cgetp(p1); l=avma; + z=cgetp(p1); av=avma; setvalp(z,valp(x)-valp(y)); p2=mpinvmod((GEN)y[4],(GEN)p1[3]); modiiz(mulii((GEN)x[4],p2),(GEN)p1[3],(GEN)z[4]); - avma=l; return z; + avma=av; return z; case t_COMPLEX: case t_QUAD: - l=avma; p1=gmul(x,gconj(y)); p2=gnorm(y); tetpil=avma; - return gerepile(l,tetpil,gdiv(p1,p2)); + av=avma; p1=gmul(x,gconj(y)); p2=gnorm(y); tetpil=avma; + return gerepile(av,tetpil,gdiv(p1,p2)); case t_REAL: - err(operf,"/",tx,ty); + err(operf,"/",x,y); } case t_QUAD: @@ -2044,23 +2082,23 @@ gdiv(GEN x, GEN y) return z; case t_REAL: - l=avma; p1=co8(x,lg(y)); tetpil=avma; - return gerepile(l,tetpil,gdiv(p1,y)); + av=avma; p1=co8(x,lg(y)); tetpil=avma; + return gerepile(av,tetpil,gdiv(p1,y)); case t_PADIC: - l=avma; p1=cvtop(x,(GEN)y[2],precp(y)); - tetpil=avma; return gerepile(l,tetpil,gdiv(p1,y)); + av=avma; p1=cvtop(x,(GEN)y[2],precp(y)); + tetpil=avma; return gerepile(av,tetpil,gdiv(p1,y)); case t_COMPLEX: - ly=precision(y); if (!ly) err(operi,"/",tx,ty); - l=avma; p1=co8(x,ly); tetpil=avma; - return gerepile(l,tetpil,gdiv(p1,y)); + ly=precision(y); if (!ly) err(operi,"/",x,y); + av=avma; p1=co8(x,ly); tetpil=avma; + return gerepile(av,tetpil,gdiv(p1,y)); case t_QUAD: k=x[1]; l=y[1]; - if (!gegal((GEN)k,(GEN)l)) err(operi,"/",tx,ty); - l=avma; p1=gnorm(y); p2=gmul(x,gconj(y)); tetpil=avma; - return gerepile(l,tetpil,gdiv(p2,p1)); + if (!gegal((GEN)k,(GEN)l)) err(operi,"/",x,y); + av=avma; p1=gnorm(y); p2=gmul(x,gconj(y)); tetpil=avma; + return gerepile(av,tetpil,gdiv(p2,p1)); } } err(bugparier,"division"); @@ -2085,7 +2123,7 @@ gdiv(GEN x, GEN y) } return gerepileupto(av, gmul(x, ginv(y))); } - if (is_noncalc_t(tx) || is_noncalc_t(ty)) err(operf,"/",tx, ty); + if (is_noncalc_t(tx) || is_noncalc_t(ty)) err(operf,"/",x,y); /* now x and y are not both is_scalar_t */ lx = lg(x); @@ -2107,7 +2145,7 @@ gdiv(GEN x, GEN y) for (i=lontyp[tx]; i=10 et ty>=10*/ @@ -2167,7 +2205,7 @@ gdiv(GEN x, GEN y) z[1]=lmul(x,(GEN)y[2]); z[2]=lcopy((GEN)y[1]); return z; - default: err(operf,"/",tx,ty); + default: err(operf,"/",x,y); } case t_SER: @@ -2227,10 +2265,10 @@ gdiv(GEN x, GEN y) } case t_RFRAC: case t_RFRACN: - l=avma; p2=gmul(x,(GEN)y[2]); tetpil=avma; - return gerepile(l,tetpil,gdiv(p2,(GEN)y[1])); + av=avma; p2=gmul(x,(GEN)y[2]); tetpil=avma; + return gerepile(av,tetpil,gdiv(p2,(GEN)y[1])); - default: err(operf,"/",tx,ty); + default: err(operf,"/",x,y); } case t_RFRAC: case t_RFRACN: @@ -2243,8 +2281,8 @@ gdiv(GEN x, GEN y) z[1]=lcopy((GEN)x[1]); return z; case t_SER: - l=avma; p2=gmul((GEN)x[2],y); tetpil=avma; - return gerepile(l,tetpil, gdiv((GEN)x[1],p2)); + av=avma; p2=gmul((GEN)x[2],y); tetpil=avma; + return gerepile(av,tetpil, gdiv((GEN)x[1],p2)); case t_RFRAC: case t_RFRACN: if (tx == t_RFRACN) ty=t_RFRACN; @@ -2253,7 +2291,7 @@ gdiv(GEN x, GEN y) z[1]=lmul((GEN)x[1],(GEN)y[2]); z[2]=lmul((GEN)x[2],(GEN)y[1]); return z; - default: err(operf,"/",tx,ty); + default: err(operf,"/",x,y); } case t_VEC: case t_COL: case t_MAT: @@ -2263,12 +2301,12 @@ gdiv(GEN x, GEN y) for (i=1; i