OpenXM_contrib/pari-2.2/src/basemath/gen3.c - diff

Return to gen3.c CVS log

Up to [local] / OpenXM_contrib / pari-2.2 / src / basemath

Diff for /OpenXM_contrib/pari-2.2/src/basemath/Attic/gen3.c between version 1.1 and 1.2

version 1.1, 2001/10/02 11:17:04

version 1.2, 2002/09/11 07:26:51

Line 197 padicprec(GEN x, GEN p)

return 0; /* not reached */

}

/* Degre de x par rapport a la variable v si v>=0, par rapport a la variable

/* Degree of x (scalar, t_POL, t_RFRAC) wrt variable v if v >= 0,

* principale si v<0. On suppose x fraction rationnelle ou polynome.

* wrt to main variable if v < 0. Convention: deg(0) = -1.

* Convention deg(0)=-1.

long

poldegree(GEN x, long v)

{

long tx=typ(x), av, w, d;

long tx = typ(x), lx,w,i,d;

if (is_scalar_t(tx)) return gcmp0(x)? -1: 0;

switch(tx)

Line 213 poldegree(GEN x, long v)

Line 212 poldegree(GEN x, long v)

w = varn(x);

if (v < 0 || v == w) return degpol(x);

if (v < w) return signe(x)? 0: -1;

av = avma; x = gsubst(gsubst(x,w,polx[MAXVARN]),v,polx[0]);

lx = lgef(x); d = -1;

if (gvar(x)) { d = gcmp0(x)? -1: 0; } else d = degpol(x);

for (i=2; i<lx; i++)

avma = av; return d;

{

long e = poldegree((GEN)x[i], v);

if (e > d) d = e;

}

return d;

case t_RFRAC: case t_RFRACN:

if (gcmp0((GEN)x[1])) return -1;

Line 237 degree(GEN x)

Line 240 degree(GEN x)

GEN

pollead(GEN x, long v)

{

long l,tx = typ(x),av,tetpil,w;

long l, tx = typ(x), w;

gpmem_t av, tetpil;

GEN xinit;

if (is_scalar_t(tx)) return gcopy(x);

Line 375 ismonome(GEN x)

Line 379 ismonome(GEN x)

/* assume z[1] was created last */

#define fix_frac_if_int(z) if (is_pm1(z[2]))\

z = gerepileuptoint((long)(z+3), (GEN)z[1]);

z = gerepileuptoint((gpmem_t)(z+3), (GEN)z[1]);

GEN

gmulsg(long s, GEN y)

{

long ty=typ(y),ly=lg(y),i,av,tetpil;

long ty=typ(y), ly=lg(y), i;

gpmem_t av, tetpil;

GEN z,p1,p2;

switch(ty)

Line 393 gmulsg(long s, GEN y)

Line 398 gmulsg(long s, GEN y)

case t_INTMOD: z=cgetg(3,t_INTMOD); p2=(GEN)y[1];

(void)new_chunk(lgefint(p2)<<2); /* HACK */

p1=mulsi(s,(GEN)y[2]); avma=(long)z;

p1=mulsi(s,(GEN)y[2]); avma=(gpmem_t)z;

z[2]=lmodii(p1,p2); icopyifstack(p2,z[1]); return z;

case t_FRAC:

Line 453 gmulsg(long s, GEN y)

Line 458 gmulsg(long s, GEN y)

if (!s) return zeropol(gvar(y));

z = cgetg(3, t_RFRAC);

i = ggcd(stoi(s),(GEN)y[2])[2];

avma = (long)z;

avma = (gpmem_t)z;

if (i == 1)

{

z[1]=lmulgs((GEN)y[1], s);

Line 490 gmulsg(long s, GEN y)

Line 495 gmulsg(long s, GEN y)

GEN

gdivgs(GEN x, long s)

{

static long court[] = { evaltyp(t_INT) | m_evallg(3),0,0 };

static long court[] = { evaltyp(t_INT) | _evallg(3),0,0 };

long tx=typ(x),lx=lg(x),i,av;

long tx=typ(x), lx=lg(x), i;

gpmem_t av;

GEN z,p1;

if (!s) err(gdiver2);

Line 560 gdivgs(GEN x, long s)

Line 566 gdivgs(GEN x, long s)

z[1]=x[1]; return z;

case t_RFRAC:

z = cgetg(3, t_RFRAC);

av = avma;

i = ggcd(stoi(s),(GEN)x[1])[2];

p1 = ggcd(stoi(s),(GEN)x[1]);

avma = (long)z;

if (typ(p1) == t_INT)

if (i == 1)

{

z[2]=lmulsg(s,(GEN)x[2]);

avma = av;

z[1]=lcopy((GEN)x[1]); return z;

z = cgetg(3, t_RFRAC);

i = p1[2];

if (i == 1)

{

z[1] = lcopy((GEN)x[1]);

z[2] = lmulsg(s,(GEN)x[2]);

}

else

{

z[1] = ldivgs((GEN)x[1], i);

z[2] = lmulgs((GEN)x[2], s/i);

}

z[1] = ldivgs((GEN)x[1], i);

else /* t_FRAC */

z[2] = lmulgs((GEN)x[2], s/i); return z;

{

z = cgetg(3, t_RFRAC);

z[1] = ldiv((GEN)x[1], p1);

z[2] = lmul((GEN)x[2], gdivsg(s,p1));

z = gerepilecopy(av, z);

}

return z;

case t_RFRACN: z=cgetg(3,t_RFRACN);

z[2]=lmulsg(s,(GEN)x[2]);

z[1]=lcopy((GEN)x[1]);

z[1]=lcopy((GEN)x[1]); return z;

z[2]=lmulsg(s,(GEN)x[2]); return z;

case t_VEC: case t_COL: case t_MAT: z=cgetg(lx,tx);

for (i=1; i<lx; i++) z[i]=ldivgs((GEN)x[i],s);

Line 584 gdivgs(GEN x, long s)

Line 606 gdivgs(GEN x, long s)

/*******************************************************************/

/* */

/* MODULO GENERAL */

/* GENERIC REMAINDER */

/* */

/*******************************************************************/

/* euclidean quotient for scalars of admissible types */

static GEN

_quot(GEN x, GEN y)

{

GEN q = gdiv(x,y), f = gfloor(q);

if (gsigne(y) < 0 && !gegal(f,q)) f = gadd(f,gun);

return f;

}

static GEN

quot(GEN x, GEN y)

{

gpmem_t av = avma;

return gerepileupto(av, _quot(x, y));

}

GEN

gmod(GEN x, GEN y)

{

long av,tetpil,i, tx=typ(x), ty = typ(y);

gpmem_t av, tetpil;

long i,lx,ty, tx = typ(x);

GEN z,p1;

if (is_matvec_t(tx))

{

i=lg(x); z=cgetg(i,tx);

lx=lg(x); z=cgetg(lx,tx);

for (i--; i; i--) z[i]=lmod((GEN)x[i],y);

for (i=1; i<lx; i++) z[i] = lmod((GEN)x[i],y);

return z;

}

ty = typ(y);

switch(ty)

{

case t_INT:

Line 612 gmod(GEN x, GEN y)

Line 650 gmod(GEN x, GEN y)

z[1]=lmppgcd((GEN)x[1],y);

z[2]=lmodii((GEN)x[2],(GEN)z[1]); return z;

case t_FRAC: case t_FRACN:

case t_FRACN:

av=avma; if (tx==t_FRACN) x=gred(x);

av=avma; x = gred(x);

return gerepileupto(av, gmod(x,y));

case t_FRAC:

av=avma;

p1=mulii((GEN)x[1],mpinvmod((GEN)x[2],y));

tetpil=avma; return gerepile(av,tetpil,modii(p1,y));

Line 624 gmod(GEN x, GEN y)

Line 665 gmod(GEN x, GEN y)

case t_PADIC:

{

long p1[] = {evaltyp(t_INTMOD)|m_evallg(3),0,0};

long p1[] = {evaltyp(t_INTMOD)|_evallg(3),0,0};

p1[1] = (long)y;

p1[2] = lgeti(lgefint(y));

gaffect(x,p1); return (GEN)p1[2];

}

case t_POLMOD: case t_POL:

return gzero;

/* case t_REAL could be defined as below, but conlicting semantic

* with lift(x * Mod(1,y)). */

default: err(operf,"%",tx,ty);

default: err(operf,"%",x,y);

}

case t_REAL: case t_FRAC: case t_FRACN:

switch(tx)

{

case t_INT: case t_REAL: case t_FRAC: case t_FRACN:

av=avma; p1 = gfloor(gdiv(x,y)); p1 = gneg_i(gmul(p1,y));

av = avma;

tetpil=avma; return gerepile(av,tetpil,gadd(x,p1));

return gerepileupto(av, gsub(x, gmul(_quot(x,y),y)));

case t_POLMOD: case t_POL:

return gzero;

default: err(operf,"%",tx,ty);

default: err(operf,"%",x,y);

}

case t_POL:

if (is_scalar_t(tx))

{

if (tx!=t_POLMOD || varn(x[1]) > varn(y))

return (lgef(y) < 4)? gzero: gcopy(x);

return degpol(y)? gcopy(x): gzero;

if (varn(x[1])!=varn(y)) return gzero;

z=cgetg(3,t_POLMOD);

z[1]=lgcd((GEN)x[1],y);

Line 663 gmod(GEN x, GEN y)

Line 706 gmod(GEN x, GEN y)

case t_POL:

return gres(x,y);

case t_RFRAC: case t_RFRACN:

case t_RFRACN:

av=avma; if (tx==t_RFRACN) x=gred_rfrac(x);

av=avma; return gerepileupto(av, gmod(gred_rfrac(x), y));

case t_RFRAC:

av=avma;

p1=gmul((GEN)x[1],ginvmod((GEN)x[2],y)); tetpil=avma;

return gerepile(av,tetpil,gres(p1,y));

Line 672 gmod(GEN x, GEN y)

Line 717 gmod(GEN x, GEN y)

if (ismonome(y) && varn(x) == varn(y))

{

long d = degpol(y);

if (lg(x)-2 + valp(x) < d) err(operi,"%",tx,ty);

if (lg(x)-2 + valp(x) < d) err(operi,"%",x,y);

av = avma;

return gerepileupto(av, gmod(gtrunc(x), y));

}

default: err(operf,"%",tx,ty);

default: err(operf,"%",x,y);

}

err(operf,"%",tx,ty);

err(operf,"%",x,y);

return NULL; /* not reached */

}

Line 696 gmodulsg(long x, GEN y)

Line 741 gmodulsg(long x, GEN y)

case t_POL: z=cgetg(3,t_POLMOD);

z[1]=lcopy(y); z[2]=lstoi(x); return z;

}

err(operf,"%",t_INT,typ(y)); return NULL; /* not reached */

err(operf,"%",stoi(x),y); return NULL; /* not reached */

}

GEN

Line 743 gmodulo(GEN x,GEN y)

Line 788 gmodulo(GEN x,GEN y)

if (tx!=t_POL && !is_rfrac_t(tx) && tx!=t_SER) break;

z[2]=(long)specialmod(x,y); return z;

}

err(operf,"%",tx,typ(y)); return NULL; /* not reached */

err(operf,"%",x,y); return NULL; /* not reached */

}

GEN

Line 776 gmodulcp(GEN x,GEN y)

Line 821 gmodulcp(GEN x,GEN y)

if (tx!=t_POL && !is_rfrac_t(tx) && tx!=t_SER) break;

z[2]=(long)specialmod(x,y); return z;

}

err(operf,"%",tx,typ(y)); return NULL; /* not reached */

err(operf,"%",x,y); return NULL; /* not reached */

}

GEN

Line 793 Mod0(GEN x,GEN y,long flag)

Line 838 Mod0(GEN x,GEN y,long flag)

/*******************************************************************/

/* */

/* DIVISION ENTIERE GENERALE */

/* GENERIC EUCLIDEAN DIVISION */

/* DIVISION ENTIERE AVEC RESTE GENERALE */

/* */

/*******************************************************************/

GEN

gdivent(GEN x, GEN y)

{

long tx=typ(x),ty=typ(y);

long tx = typ(x);

if (tx == t_INT)

if (is_matvec_t(tx))

{

if (ty == t_INT) return truedvmdii(x,y,NULL);

long i, lx = lg(x);

if (ty!=t_POL) err(typeer,"gdivent");

GEN z = cgetg(lx,tx);

return gzero;

for (i=1; i<lx; i++) z[i] = (long)gdivent((GEN)x[i],y);

return z;

}

if (ty != t_POL) err(typeer,"gdivent");

switch(typ(y))

if (tx == t_POL) return gdeuc(x,y);

{

if (! is_scalar_t(tx)) err(typeer,"gdivent");

case t_INT:

return gzero;

switch(tx)

{ /* equal to, but more efficient than, quot(x,y) */

case t_INT: return truedvmdii(x,y,NULL);

case t_REAL:

case t_FRAC:

case t_FRACN: return quot(x,y);

case t_POL: return gdiv(x,y);

}

break;

case t_REAL:

case t_FRAC:

case t_FRACN: return quot(x,y);

case t_POL:

if (is_scalar_t(tx))

{

if (tx == t_POLMOD) break;

return degpol(y)? gzero: gdiv(x,y);

}

if (tx == t_POL) return gdeuc(x,y);

}

err(operf,"\\",x,y);

return NULL; /* not reached */

}

/* with remainder */

static GEN

quotrem(GEN x, GEN y, GEN *r)

{

gpmem_t av;

GEN q = quot(x,y);

av = avma;

*r = gerepileupto(av, gsub(x, gmul(q,y)));

return q;

}

GEN

gdiventres(GEN x, GEN y)

{

long tx=typ(x),ty=typ(y);

long tx = typ(x);

GEN z = cgetg(3,t_COL);

GEN z,q,r;

if (tx==t_INT)

if (is_matvec_t(tx))

{

if (ty==t_INT)

long i, lx = lg(x);

z = cgetg(lx,tx);

for (i=1; i<lx; i++) z[i] = (long)gdiventres((GEN)x[i],y);

return z;

}

z = cgetg(3,t_COL);

switch(typ(y))

{

case t_INT:

switch(tx)

{ /* equal to, but more efficient than next case */

case t_INT:

z[1] = (long)truedvmdii(x,y,(GEN*)(z+2));

return z;

case t_REAL:

case t_FRAC:

case t_FRACN:

q = quotrem(x,y,&r);

z[1] = (long)q;

z[2] = (long)r; return z;

case t_POL: return gdiv(x,y);

}

break;

case t_REAL:

case t_FRAC:

case t_FRACN:

q = quotrem(x,y,&r);

z[1] = (long)q;

z[2] = (long)r; return z;

case t_POL:

if (is_scalar_t(tx))

{

if (tx == t_POLMOD) break;

if (degpol(y))

{

q = gzero;

r = gcopy(x);

}

else

{

q = gdiv(x,y);

r = gzero;

}

return q;

}

if (tx == t_POL)

{

z[1]=ldivres(x,y,(GEN *)(z+2));

return z;

}

err(operf,"\\",x,y);

return NULL; /* not reached */

}

extern GEN swap_vars(GEN b0, long v);

GEN

divrem(GEN x, GEN y, long v)

{

gpmem_t av = avma;

long vx,vy;

GEN z,q,r;

if (v < 0 || typ(y) != t_POL || typ(x) != t_POL) return gdiventres(x,y);

vx = varn(x); if (vx != v) x = swap_vars(x,v);

vy = varn(y); if (vy != v) y = swap_vars(y,v);

q = poldivres(x,y, &r);

if (v && (vx != v || vy != v))

{

q = poleval(q, polx[v]);

r = poleval(r, polx[v]);

}

z = cgetg(3,t_COL);

z[1] = (long)q;

z[2] = (long)r; return gerepilecopy(av, z);

}

static int

is_scal(long t) { return t==t_INT || t==t_FRAC || t==t_FRACN; }

GEN

diviiround(GEN x, GEN y)

{

gpmem_t av1, av = avma;

GEN q,r;

int fl;

q = dvmdii(x,y,&r); /* q = x/y rounded towards 0, sgn(r)=sgn(x) */

if (r==gzero) return q;

av1 = avma;

fl = absi_cmp(shifti(r,1),y);

avma = av1; cgiv(r);

if (fl >= 0) /* If 2*|r| >= |y| */

{

long sz = signe(x)*signe(y);

if (fl || sz > 0) q = gerepileuptoint(av, addis(q,sz));

}

return q;

}

/* If x and y are not both scalars, same as gdivent.

* Otherwise, compute the quotient x/y, rounded to the nearest integer

* (towards +oo in case of tie). */

GEN

gdivround(GEN x, GEN y)

{

gpmem_t av1, av;

long tx=typ(x),ty=typ(y);

GEN q,r;

int fl;

if (tx==t_INT && ty==t_INT) return diviiround(x,y);

av = avma;

if (is_scal(tx) && is_scal(ty))

{ /* same as diviiround but less efficient */

q = quotrem(x,y,&r);

av1 = avma;

fl = gcmp(gmul2n(gabs(r,0),1), gabs(y,0));

avma = av1; cgiv(r);

if (fl >= 0) /* If 2*|r| >= |y| */

{

z[1]=(long)truedvmdii(x,y,(GEN*)(z+2));

long sz = gsigne(y);

return z;

if (fl || sz > 0) q = gerepileupto(av, gaddgs(q, sz));

}

if (ty!=t_POL) err(typeer,"gdiventres");

return q;

z[1]=zero; z[2]=licopy(x); return z;

}

if (ty != t_POL) err(typeer,"gdiventres");

if (is_matvec_t(tx))

if (tx == t_POL)

{

z[1]=ldivres(x,y,(GEN *)(z+2));

long i,lx = lg(x);

GEN z = cgetg(lx,tx);

for (i=1; i<lx; i++) z[i] = (long)gdivround((GEN)x[i],y);

return z;

}

if (! is_scalar_t(tx)) err(typeer,"gdiventres");

return gdivent(x,y);

z[1]=zero; z[2]=lcopy(x); return z;

}

GEN

Line 857 gdivmod(GEN x, GEN y, GEN *pr)

Line 1052 gdivmod(GEN x, GEN y, GEN *pr)

return poldivres(x,y,pr);

}

/* When x and y are integers, compute the quotient x/y, rounded to the

* nearest integer. If there is a tie, the quotient is rounded towards zero.

* If x and y are not both integers, same as gdivent.

GEN

gdivround(GEN x, GEN y)

{

long tx=typ(x),ty=typ(y);

if (tx==t_INT)

{

if (ty==t_INT)

{

long av = avma, av1,fl;

GEN r, q=dvmdii(x,y,&r); /* q = x/y rounded towards 0, sqn(r)=sgn(x) */

if (r==gzero) return q;

av1 = avma;

fl = absi_cmp(shifti(r,1),y);

avma = av1; cgiv(r);

if (fl >= 0) /* If 2*|r| >= |y| */

{

long sz = signe(x)*signe(y);

if (fl || sz > 0)

{ av1=avma; q = gerepile(av,av1,addis(q,sz)); }

}

return q;

}

if (ty!=t_POL) err(typeer,"gdivround");

return gzero;

}

if (ty != t_POL) err(typeer,"gdivround");

if (tx == t_POL) return gdeuc(x,y);

if (! is_scalar_t(tx)) err(typeer,"gdivround");

return gzero;

}

/*******************************************************************/

/* */

/* SHIFT D'UN GEN */

Line 904 gdivround(GEN x, GEN y)

Line 1062 gdivround(GEN x, GEN y)

GEN

gshift(GEN x, long n)

{

return gshift3(x, n, 0);

}

GEN

gshift3(GEN x, long n, long flag)

{

long i,l,lx, tx = typ(x);

GEN y;

switch(tx)

{

case t_INT:

return shifti(x,n);

return shifti3(x,n,flag);

case t_REAL:

return shiftr(x,n);

Line 929 extern GEN mulscalrfrac(GEN x, GEN y);

Line 1093 extern GEN mulscalrfrac(GEN x, GEN y);

GEN

gmul2n(GEN x, long n)

{

long tx,lx,i,k,l,av,tetpil;

long tx, lx, i, k, l;

gpmem_t av, tetpil;

GEN p2,p1,y;

tx=typ(x);

Line 951 gmul2n(GEN x, long n)

Line 1116 gmul2n(GEN x, long n)

if (n > 0)

{

y=cgetg(3,t_INTMOD); p2=(GEN)x[1];

av=avma; new_chunk(2 * lgefint(p2) + n); /* HACK */

av=avma; (void)new_chunk(2 * lgefint(p2) + n); /* HACK */

p1 = shifti((GEN)x[2],n); avma=av;

y[2]=lmodii(p1,p2); icopyifstack(p2,y[1]); return y;

}

l=avma; y=gmul2n(gun,n); tetpil=avma;

av=avma; y=gmul2n(gun,n); tetpil=avma;

return gerepile(l,tetpil,gmul(y,x));

return gerepile(av,tetpil,gmul(y,x));

case t_FRAC: case t_FRACN:

l = vali((GEN)x[1]);

Line 1007 gmul2n(GEN x, long n)

Line 1172 gmul2n(GEN x, long n)

return y;

case t_PADIC:

l=avma; y=gmul2n(gun,n); tetpil=avma;

av=avma; y=gmul2n(gun,n); tetpil=avma;

return gerepile(l,tetpil,gmul(y,x));

return gerepile(av,tetpil,gmul(y,x));

}

err(typeer,"gmul2n");

return NULL; /* not reached */

}

GEN

real2n(long n, long prec) { GEN z = realun(prec); setexpo(z, n); return z; }

/*******************************************************************/

/* */

/* INVERSE D'UN GEN */

Line 1024 extern GEN fix_rfrac_if_pol(GEN x, GEN y);

Line 1192 extern GEN fix_rfrac_if_pol(GEN x, GEN y);

GEN

ginv(GEN x)

{

long tx=typ(x),av,tetpil,s;

long tx=typ(x), s;

gpmem_t av, tetpil;

GEN z,y,p1,p2;

switch(tx)

Line 1098 ginv(GEN x)

Line 1267 ginv(GEN x)

return y;

case t_MAT:

return (lg(x)==1)? cgetg(1,t_MAT): invmat(x);

case t_VECSMALL:

{

long i,lx = lg(x);

y = cgetg(lx,t_VECSMALL);

for (i=1; i<lx; i++)

{

long xi=x[i];

if (xi<1 || xi>=lx) err(talker,"incorrect permtuation to inverse");

y[xi] = i;

}

return y;

}

err(typeer,"inverse");

return NULL; /* not reached */

Line 1113 ginv(GEN x)

Line 1294 ginv(GEN x)

static GEN

gconvsp(GEN x, int flpile)

{

long v = varn(x), av,tetpil,i;

long v = varn(x), i;

gpmem_t av, tetpil;

GEN p1,y;

if (gcmp0(x)) return zeropol(v);

av=avma; y=dummycopy(x); settyp(y,t_POL);

i=lg(x)-1; while (i>1 && gcmp0((GEN)y[i])) i--;

setlgef(y,i+1);

p1=gpuigs(polx[v],valp(x));

p1=gpowgs(polx[v],valp(x));

tetpil=avma; p1=gmul(p1,y);

return flpile? gerepile(av,tetpil,p1): p1;

}

Line 1128 gconvsp(GEN x, int flpile)

Line 1310 gconvsp(GEN x, int flpile)

GEN

gsubst0(GEN x, GEN T, GEN y)

{

ulong av;

gpmem_t av;

long d, v;

if (typ(T) != t_POL || !ismonome(T) || !gcmp1(leading_term(T)))

err(talker,"variable number expected in subst");

Line 1142 GEN

Line 1324 GEN

gsubst(GEN x, long v, GEN y)

{

long tx = typ(x), ty = typ(y), lx = lg(x), ly = lg(y);

long l,vx,vy,e,ex,ey,tetpil,av,i,j,k,jb,limite;

long l, vx, vy, e, ex, ey, i, j, k, jb;

gpmem_t tetpil, av, limite;

GEN t,p1,p2,z;

if (ty==t_MAT)

Line 1178 gsubst(GEN x, long v, GEN y)

Line 1361 gsubst(GEN x, long v, GEN y)

switch(tx)

{

case t_POL:

l=lgef(x);

l = lgef(x);

if (l==2)

return (ty==t_MAT)? gscalmat(gzero,ly-1): gzero;

vx=varn(x);

vx = varn(x);

if (vx<v)

if (vx < v)

{

if (gvar(y) > vx)

{ /* easy special case */

z = cgetg(l, t_POL); z[1] = x[1];

for (i=2; i<l; i++) z[i] = lsubst((GEN)x[i],v,y);

return normalizepol_i(z,l);

}

/* general case */

av=avma; p1=polx[vx]; z= gsubst((GEN)x[l-1],v,y);

for (i=l-1; i>2; i--) z=gadd(gmul(z,p1),gsubst((GEN)x[i-1],v,y));

return gerepileupto(av,z);

}

/* v <= vx */

if (ty!=t_MAT)

return (vx>v)? gcopy(x): poleval(x,y);

return (vx > v)? gcopy(x): poleval(x,y);

if (vx>v) return gscalmat(x,ly-1);

if (vx > v) return gscalmat(x,ly-1);

if (l==3) return gscalmat((GEN)x[2],ly-1);

av=avma; z=(GEN)x[l-1];

for (i=l-1; i>2; i--) z=gaddmat((GEN)x[i-1],gmul(z,y));

Line 1211 gsubst(GEN x, long v, GEN y)

Line 1402 gsubst(GEN x, long v, GEN y)

p2=gsubst(p1,v,y); tetpil=avma; z=tayl(p2,vx,lx-2);

if (ex)

{

p1=gpuigs(polx[vx],ex); tetpil=avma; z=gmul(z,p1);

p1=gpowgs(polx[vx],ex); tetpil=avma; z=gmul(z,p1);

}

return gerepile(av,tetpil,z);

}

Line 1232 gsubst(GEN x, long v, GEN y)

Line 1423 gsubst(GEN x, long v, GEN y)

av=avma; z = zeroser(vy,0);

for (i=lx-1; i>=2; i--)

z = gadd((GEN)x[i], gmul(y,z));

if (ex) z = gmul(z, gpuigs(y,ex));

if (ex) z = gmul(z, gpowgs(y,ex));

return gerepileupto(av,z);

}

Line 1265 gsubst(GEN x, long v, GEN y)

Line 1456 gsubst(GEN x, long v, GEN y)

}

if (!ex) return gerepilecopy(av,z);

if (l<ly) { setlg(y,l); p1=gpuigs(y,ex); setlg(y,ly); }

if (l<ly) { setlg(y,l); p1=gpowgs(y,ex); setlg(y,ly); }

else p1=gpuigs(y,ex);

else p1=gpowgs(y,ex);

tetpil=avma; return gerepile(av,tetpil,gmul(z,p1));

case t_POL: case t_RFRAC: case t_RFRACN:

Line 1289 gsubst(GEN x, long v, GEN y)

Line 1480 gsubst(GEN x, long v, GEN y)

case t_VEC: case t_COL: case t_MAT: z=cgetg(lx,tx);

for (i=1; i<lx; i++) z[i]=lsubst((GEN)x[i],v,y);

return z;

}

return z;

return gcopy(x);

}

/*******************************************************************/

Line 1302 gsubst(GEN x, long v, GEN y)

Line 1494 gsubst(GEN x, long v, GEN y)

GEN

recip(GEN x)

{

long tetpil, av=avma, v=varn(x);

long v=varn(x), lx = lg(x);

gpmem_t tetpil, av=avma;

GEN p1,p2,a,y,u;

if (typ(x)!=t_SER) err(talker,"not a series in serreverse");

if (valp(x)!=1) err(talker,"valuation not equal to 1 in serreverse");

if (valp(x)!=1 || lx < 3)

err(talker,"valuation not equal to 1 in serreverse");

a=(GEN)x[2];

if (gcmp1(a))

{

long i,j,k, lx=lg(x), lim=stack_lim(av,2);

long i, j, k, mi;

gpmem_t lim=stack_lim(av, 2);

u=cgetg(lx,t_SER); y=cgetg(lx,t_SER);

mi = lx-1; while (mi>=3 && gcmp0((GEN)x[mi])) mi--;

u[1]=y[1]=evalsigne(1) | evalvalp(1) | evalvarn(v);

u = cgetg(lx,t_SER);

u[2]=un; u[3]=lmulsg(-2,(GEN)x[3]);

y = cgetg(lx,t_SER);

y[2]=un; y[3]=lneg((GEN)x[3]);

u[1] = y[1] = evalsigne(1) | evalvalp(1) | evalvarn(v);

u[2] = y[2] = un;

if (lx > 3)

{

u[3] = lmulsg(-2,(GEN)x[3]);

y[3] = lneg((GEN)x[3]);

}

for (i=3; i<lx-1; )

{

for (j=3; j<i+1; j++)

{

p1=(GEN)u[j];

p1 = (GEN)x[j];

for (k=j-1; k>2; k--)

for (k=max(3,j+2-mi); k<j; k++)

p1=gsub(p1,gmul((GEN)u[k],(GEN)x[j-k+2]));

p1 = gadd(p1, gmul((GEN)u[k],(GEN)x[j-k+2]));

u[j]=lsub(p1,(GEN)x[j]);

u[j] = lsub((GEN)u[j], p1);

}

p1=gmulsg(i,(GEN)x[i+1]);

p1 = gmulsg(i,(GEN)x[i+1]);

for (k=2; k<i; k++)

for (k=2; k<min(i,mi); k++)

{

p2=gmul((GEN)x[k+1],(GEN)u[i-k+2]);

p2 = gmul((GEN)x[k+1],(GEN)u[i-k+2]);

p1=gadd(p1,gmulsg(k,p2));

p1 = gadd(p1, gmulsg(k,p2));

}

i++; u[i]=lneg(p1); y[i]=ldivgs((GEN)u[i],i-1);

i++;

u[i] = lneg(p1);

y[i] = ldivgs((GEN)u[i],i-1);

if (low_stack(lim, stack_lim(av,2)))

{

GEN *gptr[2];

Line 1343 recip(GEN x)

Line 1546 recip(GEN x)

}

return gerepilecopy(av,y);

}

y=gdiv(x,a); y[2]=un; y=recip(y);

y = gdiv(x,a); y[2] = un; y = recip(y);

a=gdiv(polx[v],a); tetpil=avma;

a = gdiv(polx[v],a); tetpil = avma;

return gerepile(av,tetpil,gsubst(y,v,a));

return gerepile(av,tetpil, gsubst(y,v,a));

}

/*******************************************************************/

Line 1359 derivpol(GEN x)

Line 1562 derivpol(GEN x)

long i,lx = lgef(x)-1;

GEN y;

if (lx<3) return gzero;

if (lx<3) return zeropol(varn(x));

y = cgetg(lx,t_POL);

for (i=2; i<lx ; i++) y[i] = lmulsg(i-1,(GEN)x[i+1]);

y[1] = x[1]; return normalizepol_i(y,i);

Line 1370 derivser(GEN x)

Line 1573 derivser(GEN x)

{

long i,j,vx = varn(x), e = valp(x), lx = lg(x);

GEN y;

if (gcmp0(x)) return zeroser(vx,e-1);

if (gcmp0(x)) return zeroser(vx,e? e-1: 0);

if (e)

{

y=cgetg(lx,t_SER); y[1] = evalsigne(1) | evalvalp(e-1) | evalvarn(vx);

Line 1389 derivser(GEN x)

Line 1592 derivser(GEN x)

GEN

deriv(GEN x, long v)

{

long lx,vx,tx,e,i,j,l,av,tetpil;

long lx, vx, tx, e, i, j;

gpmem_t av, tetpil;

GEN y,p1,p2;

tx=typ(x); if (is_const_t(tx)) return gzero;

Line 1419 deriv(GEN x, long v)

Line 1623 deriv(GEN x, long v)

{

if (!signe(x)) return gcopy(x);

lx=lg(x); e=valp(x);

l=avma;

av=avma;

for (i=2; i<lx; i++)

{

if (!gcmp0(deriv((GEN)x[i],v))) break;

avma=l;

avma=av;

}

if (i==lx) return ggrando(polx[vx],e+lx-2);

y=cgetg(lx-i+2,t_SER);

Line 1434 deriv(GEN x, long v)

Line 1638 deriv(GEN x, long v)

return derivser(x);

case t_RFRAC: case t_RFRACN: av=avma; y=cgetg(3,tx);

y[2]=lsqr((GEN)x[2]); l=avma;

y[2]=lsqr((GEN)x[2]); av=avma;

p1=gmul((GEN)x[2],deriv((GEN)x[1],v));

p2=gmul(gneg_i((GEN)x[1]),deriv((GEN)x[2],v));

tetpil=avma; p1=gadd(p1,p2);

if (tx==t_RFRACN) { y[1]=lpile(l,tetpil,p1); return y; }

if (tx==t_RFRACN) { y[1]=lpile(av,tetpil,p1); return y; }

y[1]=(long)p1; return gerepileupto(av,gred_rfrac(y));

case t_VEC: case t_COL: case t_MAT: lx=lg(x); y=cgetg(lx,tx);

Line 1469 triv_integ(GEN x, long v, long tx, long lx)

Line 1673 triv_integ(GEN x, long v, long tx, long lx)

GEN

integ(GEN x, long v)

{

long lx,tx,e,i,j,vx,n,av=avma,tetpil;

long lx, tx, e, i, j, vx, n;

gpmem_t av=avma, tetpil;

GEN y,p1;

tx = typ(x);

Line 1614 gfloor(GEN x)

Line 1819 gfloor(GEN x)

GEN

gfrac(GEN x)

{

long av = avma,tetpil;

gpmem_t av = avma, tetpil;

GEN p1 = gneg_i(gfloor(x));

tetpil=avma; return gerepile(av,tetpil,gadd(x,p1));

Line 1624 GEN

Line 1829 GEN

gceil(GEN x)

{

GEN y,p1;

long i,lx,tx=typ(x),av,tetpil;

long i, lx, tx=typ(x);

gpmem_t av, tetpil;

switch(tx)

{

Line 1671 GEN

Line 1877 GEN

ground(GEN x)

{

GEN y,p1;

long i,lx,tx=typ(x),av,tetpil;

long i, lx, tx=typ(x);

gpmem_t av;

switch(tx)

{

Line 1683 ground(GEN x)

Line 1890 ground(GEN x)

long ex, s = signe(x);

if (s==0 || (ex=expo(x)) < -1) return gzero;

if (ex < 0) return s>0? gun: negi(gun);

av=avma; p1 = realun(3 + (ex>>TWOPOTBITS_IN_LONG));

av = avma;

setexpo(p1,-1); /* p1 = 0.5 */

p1 = real2n(-1, 3 + (ex>>TWOPOTBITS_IN_LONG)); /* = 0.5 */

p1 = addrr(x,p1); tetpil = avma;

p1 = addrr(x,p1);

return gerepile(av,tetpil,mpent(p1));

return gerepileuptoint(av, mpent(p1));

}

case t_FRAC: case t_FRACN:

av=avma; p1 = addii(shifti((GEN)x[2], -1), (GEN)x[1]);

return diviiround((GEN)x[1], (GEN)x[2]);

return gerepileuptoint(av, truedvmdii(p1, (GEN)x[2], NULL));

case t_POLMOD: y=cgetg(3,t_POLMOD);

copyifstack(x[1],y[1]);

Line 1715 GEN

Line 1921 GEN

grndtoi(GEN x, long *e)

{

GEN y,p1;

long i,tx=typ(x), lx,av,ex,e1;

long i, tx=typ(x), lx, ex, e1;

gpmem_t av;

*e = -HIGHEXPOBIT;

switch(tx)

Line 1734 grndtoi(GEN x, long *e)

Line 1941 grndtoi(GEN x, long *e)

lx=lg(x); e1 = ex - bit_accuracy(lx) + 1;

settyp(p1,t_INT); setlgefint(p1,lx);

y=shifti(p1,e1); if (signe(x)<0) y=addsi(-1,y);

y = gerepileupto(av,y);

y = gerepileuptoint(av,y);

if (e1<=0) { av=avma; e1=expo(subri(x,y)); avma=av; }

*e=e1; return y;

Line 1763 grndtoi(GEN x, long *e)

Line 1970 grndtoi(GEN x, long *e)

GEN

gcvtoi(GEN x, long *e)

{

long tx=typ(x), lx,i,ex,av,e1;

long tx=typ(x), lx, i, ex, e1;

gpmem_t av;

GEN y;

*e = -HIGHEXPOBIT;

Line 1795 gcvtoi(GEN x, long *e)

Line 2003 gcvtoi(GEN x, long *e)

GEN

ceil_safe(GEN x)

{

ulong av = avma;

gpmem_t av = avma;

long e;

GEN y = gcvtoi(x,&e);

if (e < 0) e = 0;

Line 1806 ceil_safe(GEN x)

Line 2014 ceil_safe(GEN x)

GEN

gtrunc(GEN x)

{

long tx=typ(x),av,tetpil,i,v;

long tx=typ(x), i, v;

gpmem_t av;

GEN y;

switch(tx)

Line 1826 gtrunc(GEN x)

Line 2035 gtrunc(GEN x)

if (!v) return gcopy((GEN)x[4]);

if (v>0)

{ /* here p^v is an integer */

av=avma; y=gpuigs((GEN)x[2],v); tetpil=avma;

av = avma; y = gpowgs((GEN)x[2],v);

return gerepile(av,tetpil, mulii(y,(GEN)x[4]));

return gerepileuptoint(av, mulii(y,(GEN)x[4]));

}

y=cgetg(3,t_FRAC);

y[1]=licopy((GEN)x[4]);

y[1] = licopy((GEN)x[4]);

y[2]=lpuigs((GEN)x[2],-v);

y[2] = lpowgs((GEN)x[2],-v);

return y;

case t_RFRAC: case t_RFRACN:

Line 2032 scalarser(GEN x, long v, long prec)

Line 2241 scalarser(GEN x, long v, long prec)

GEN

gtoser(GEN x, long v)

{

long tx=typ(x),lx,i,j,l,av,tetpil;

long tx=typ(x), lx, i, j, l;

gpmem_t av, tetpil;

GEN y,p1,p2;

if (v<0) v = 0;

Line 2112 gtovec(GEN x)

Line 2322 gtovec(GEN x)

for (i=1; i<=lx-2; i++) y[i]=lcopy((GEN)x[i]);

return y;

}

if (tx == t_VECSMALL)

if (tx == t_VECSMALL) return small_to_vec(x);

{

lx=lg(x); y=cgetg(lx,t_VEC);

for (i=1; i<lx; i++) y[i]=lstoi(x[i]);

return y;

}

if (!signe(x)) { y=cgetg(2,t_VEC); y[1]=zero; return y; }

lx=lg(x); y=cgetg(lx-1,t_VEC); x++;

for (i=1; i<=lx-2; i++) y[i]=lcopy((GEN)x[i]);

Line 2160 compo(GEN x, long n)

Line 2365 compo(GEN x, long n)

return gcopy((GEN)x[l]);

}

/* assume v > varn(x), extract coeff of polx[v]^n */

static GEN

multi_coeff(GEN x, long n, long v, long dx)

{

long i;

GEN z;

if (dx == 0) return polcoeff_i((GEN)x[2], n, v);

z = cgetg(dx+3, t_POL);

z[1] = x[1];

for (i=2; i<dx+3; i++) z[i] = (long)polcoeff_i((GEN)x[i], n, v);

return normalizepol(z);

}

/* assume x a t_POL */

static GEN

_polcoeff(GEN x, long n, long v)

{

long w, dx;

dx = degpol(x);

if (dx < 0) return gzero;

if (v < 0 || v == (w=varn(x)))

return (n < 0 || n > dx)? gzero: (GEN)x[n+2];

if (w > v) return n? gzero: x;

/* w < v */

return multi_coeff(x, n, v, dx);

}

/* assume x a t_SER */

static GEN

_sercoeff(GEN x, long n, long v)

{

long w, dx, ex;

if (!signe(x)) return gzero;

dx = lg(x)-3; ex = valp(x); n -= ex;

if (v < 0 || v == (w=varn(x)))

{

if (n > dx) err(talker,"non existent component in truecoeff");

return (n < 0)? gzero: (GEN)x[n+2];

}

if (w > v) return n? gzero: x;

/* w < v */

return multi_coeff(x, n, v, dx);

}

/* assume x a t_RFRAC(n) */

static GEN

_rfraccoeff(GEN x, long n, long v)

{

GEN P,Q, p = (GEN)x[1], q = (GEN)x[2];

long vp = gvar(p), vq = gvar(q);

if (v < 0) v = min(vp, vq);

P = (vp == v)? p: swap_vars(p, v);

Q = (vq == v)? q: swap_vars(q, v);

if (!ismonome(Q)) err(typeer, "polcoeff");

n += degpol(Q);

return gdiv(_polcoeff(P, n, v), leading_term(Q));

}

GEN

polcoeff_i(GEN x, long n, long v)

{

switch(typ(x))

{

case t_POL: return _polcoeff(x,n,v);

case t_SER: return _sercoeff(x,n,v);

case t_RFRACN: err(typeer, "polcoeff"); /* too expensive to reduce */

case t_RFRAC: return _rfraccoeff(x,n,v);

default: return n? gzero: x;

}

/* with respect to the main variable if v<0, with respect to the variable v

otherwise. v ignored if x is not a polynomial/series. */

GEN

polcoeff0(GEN x, long n, long v)

{

long tx=typ(x),lx,ex,w,av,tetpil;

long tx=typ(x);

GEN xinit;

gpmem_t av;

if (is_scalar_t(tx)) return n? gzero: gcopy(x);

av = avma;

switch(tx)

{

case t_POL: x = _polcoeff(x,n,v); break;

case t_SER: x = _sercoeff(x,n,v); break;

case t_RFRAC:

case t_RFRACN: x = _rfraccoeff(x,n,v); break;

case t_QFR: case t_QFI: case t_VEC: case t_COL: case t_MAT:

if (n<1 || n>=lg(x)) break;

if (n>=1 && n<lg(x)) return gcopy((GEN)x[n]);

return gcopy((GEN)x[n]);

/* fall through */

case t_POL:

default: err(talker,"nonexistent component in truecoeff");

if (n<0) return gzero;

w=varn(x);

if (v < 0 || v == w)

return (n>=lgef(x)-2)? gzero: gcopy((GEN)x[n+2]);

if (v < w) return n? gzero: gcopy(x);

av=avma; xinit=x;

x=gsubst(gsubst(x,w,polx[MAXVARN]),v,polx[0]);

if (gvar(x)) { avma=av; return n? gzero: gcopy(xinit); }

if (typ(x) == t_POL)

{

if (n>=lgef(x)-2) { avma=av; return gzero; }

x = (GEN)x[n+2];

}

else

x = polcoeff0(x, n, 0);

tetpil=avma; return gerepile(av,tetpil,gsubst(x,MAXVARN,polx[w]));

case t_SER:

w=varn(x);

if (v < 0 || v == w)

{

if (!signe(x)) return gzero;

lx=lg(x); ex=valp(x); if (n<ex) return gzero;

if (n>=ex+lx-2) break;

return gcopy((GEN)x[n-ex+2]);

}

if (v < w) return n? gzero: gcopy(x);

av=avma; xinit=x;

x=gsubst(gsubst(x,w,polx[MAXVARN]),v,polx[0]);

if (gvar(x)) { avma=av; return n? gzero: gcopy(xinit); }

if (gcmp0(x)) { avma=av; return gzero; }

if (typ(x) == t_SER)

{

lx=lg(x); ex=valp(x); if (n<ex) return gzero;

if (n>=ex+lx-2) break;

x = (GEN)x[n-ex+2];

}

else

x = polcoeff0(x, n, 0);

tetpil=avma; return gerepile(av,tetpil,gsubst(x,MAXVARN,polx[w]));

case t_RFRAC: case t_RFRACN:

w = precdl; av = avma;

if (v<0) v = gvar(x);

ex = ggval((GEN)x[2], polx[v]);

precdl = n + ex + 1; x = gtoser(x, v); precdl = w;

return gerepileupto(av, polcoeff0(x, n, v));

}

err(talker,"nonexistent component in truecoeff");

if (x == gzero) return x;

return NULL; /* not reached */

if (avma == av) return gcopy(x);

return gerepilecopy(av, x);

}

GEN

Line 2239 truecoeff(GEN x, long n)

Line 2474 truecoeff(GEN x, long n)

GEN

denom(GEN x)

{

long lx,i,av,tetpil;

long lx, i;

gpmem_t av, tetpil;

GEN s,t;

switch(typ(x))

Line 2285 denom(GEN x)

Line 2521 denom(GEN x)

GEN

numer(GEN x)

{

long av,tetpil;

gpmem_t av, tetpil;

GEN s;

switch(typ(x))

Line 2407 lift_intern0(GEN x, long v)

Line 2643 lift_intern0(GEN x, long v)

GEN

centerlift0(GEN x, long v)

{

long lx,tx=typ(x),i,av;

long lx, tx=typ(x), i;

gpmem_t av;

GEN y;

switch(tx)

Line 2463 centerlift(GEN x)

Line 2700 centerlift(GEN x)

static GEN

op_ReIm(GEN f(GEN), GEN x)

{

long lx,i,j,av, tx = typ(x);

long lx, i, j, tx = typ(x);

gpmem_t av;

GEN z;

switch(tx)

Line 2549 gimag(GEN x)

Line 2787 gimag(GEN x)

static long

_egal(GEN x, GEN y)

{

long av = avma;

gpmem_t av = avma;

long r = gegal(simplify_i(x), simplify_i(y));

avma = av; return r;

}

Line 2590 gnot(GEN x) { return gcmp0(x)? gun: gzero; }

Line 2828 gnot(GEN x) { return gcmp0(x)? gun: gzero; }

GEN

geval(GEN x)

{

long av,tetpil,lx,i, tx = typ(x);

long lx, i, tx = typ(x);

gpmem_t av, tetpil;

GEN y,z;

if (is_const_t(tx)) return gcopy(x);

Line 2710 simplify_i(GEN x)

Line 2949 simplify_i(GEN x)

GEN

simplify(GEN x)

{

ulong av = avma;

gpmem_t av = avma;

return gerepilecopy(av, simplify_i(x));

}

Line 2722 simplify(GEN x)

Line 2961 simplify(GEN x)

static GEN

qfeval0_i(GEN q, GEN x, long n)

{

long i,j,av=avma;

long i, j;

gpmem_t av=avma;

GEN res=gzero;

for (i=2;i<n;i++)

Line 2738 qfeval0_i(GEN q, GEN x, long n)

Line 2978 qfeval0_i(GEN q, GEN x, long n)

static GEN

qfeval0(GEN q, GEN x, long n)

{

long i,j,av=avma;

long i, j;

gpmem_t av=avma;

GEN res=gzero;

for (i=2;i<n;i++)

Line 2753 qfeval0(GEN q, GEN x, long n)

Line 2994 qfeval0(GEN q, GEN x, long n)

static GEN

qfeval0(GEN q, GEN x, long n)

{

long i,j,av=avma;

long i, j;

gpmem_t av=avma;

GEN res = gmul(gcoeff(q,1,1), gsqr((GEN)x[1]));

for (i=2; i<n; i++)

Line 2796 qfeval(GEN q, GEN x)

Line 3038 qfeval(GEN q, GEN x)

static GEN

qfbeval0_i(GEN q, GEN x, GEN y, long n)

{

long i,j,av=avma;

long i, j;

gpmem_t av=avma;

GEN res = gmul(gcoeff(q,1,1), mulii((GEN)x[1],(GEN)y[1]));

for (i=2;i<n;i++)

Line 2815 qfbeval0_i(GEN q, GEN x, GEN y, long n)

Line 3058 qfbeval0_i(GEN q, GEN x, GEN y, long n)

static GEN

qfbeval0(GEN q, GEN x, GEN y, long n)

{

long i,j,av=avma;

long i, j;

gpmem_t av=avma;

GEN res = gmul(gcoeff(q,1,1), gmul((GEN)x[1],(GEN)y[1]));

for (i=2;i<n;i++)

Line 2833 qfbeval0(GEN q, GEN x, GEN y, long n)

Line 3077 qfbeval0(GEN q, GEN x, GEN y, long n)

static GEN

qfbeval0(GEN q, GEN x, GEN y, long n)

{

long i,j,av=avma;

long i, j;

gpmem_t av=avma;

GEN res = gmul(gcoeff(q,1,1), gmul((GEN)x[1], (GEN)y[1]));

for (i=2; i<n; i++)

Line 2908 qf_base_change(GEN q, GEN M, int flag)

Line 3153 qf_base_change(GEN q, GEN M, int flag)

return res;

}

/* compute M'.M */

GEN

gram_matrix(GEN M)

{

long n=lg(M),i,j,k, av;

GEN res = cgetg(n,t_MAT),p1;

if (n==1)

{

if (typ(M) != t_MAT)

err(talker,"invalid data in gram_matrix");

return res;

}

if (typ(M) != t_MAT || lg(M[1]) != n)

err(talker,"not a square matrix in gram_matrix");

for (i=1;i<n;i++) res[i]=lgetg(n,t_COL);

av=avma;

for (i=1;i<n;i++)

{

p1 = gzero;

for (k=1;k<n;k++)

p1 = gadd(p1, gsqr(gcoeff(M,k,i)));

coeff(res,i,i) = (long) gerepileupto(av,p1);

av=avma;

}

for (i=2;i<n;i++)

for (j=1;j<i;j++)

{

p1=gzero;

for (k=1;k<n;k++)

p1 = gadd(p1, gmul(gcoeff(M,k,i),gcoeff(M,k,j)));

coeff(res,i,j)=coeff(res,j,i) = lpileupto(av,p1);

av=avma;

}

return res;

}

/* return Re(x * y), x and y scalars */

GEN

mul_real(GEN x, GEN y)

Line 2954 mul_real(GEN x, GEN y)

Line 3161 mul_real(GEN x, GEN y)

{

if (typ(y) == t_COMPLEX)

{

long av=avma, tetpil;

gpmem_t av=avma, tetpil;

GEN p1 = gmul((GEN)x[1], (GEN) y[1]);

GEN p2 = gneg(gmul((GEN)x[2], (GEN) y[2]));

tetpil=avma; return gerepile(av,tetpil,gadd(p1,p2));

Line 2970 mul_real(GEN x, GEN y)

Line 3177 mul_real(GEN x, GEN y)

GEN

mulmat_real(GEN x, GEN y)

{

long av,i,j,k, lx = lg(x), ly = lg(y), l = lg(x[1]);

long i, j, k, lx = lg(x), ly = lg(y), l = lg(x[1]);

gpmem_t av;

GEN p1, z = cgetg(ly,t_MAT);

for (j=1; j<ly; j++)

Line 2990 mulmat_real(GEN x, GEN y)

Line 3198 mulmat_real(GEN x, GEN y)

static GEN

hqfeval0(GEN q, GEN x, long n)

{

long i,j, av=avma;

long i, j;

gpmem_t av=avma;

GEN res=gzero;

for (i=2;i<n;i++)

Line 3028 hqfeval(GEN q, GEN x)

Line 3237 hqfeval(GEN q, GEN x)

GEN

poleval(GEN x, GEN y)

{

long av,tetpil,i,j,imin,tx=typ(x);

long i, j, imin, tx = typ(x);

GEN p1,p2,p3,r,s;

gpmem_t av0 = avma, av, lim;

GEN p1, p2, r, s;

if (is_scalar_t(tx)) return gcopy(x);

switch(tx)

{

case t_POL:

i=lgef(x)-1; imin=2; break;

i = lgef(x)-1; imin = 2; break;

case t_RFRAC: case t_RFRACN: av=avma;

case t_RFRAC: case t_RFRACN:

p1=poleval((GEN)x[1],y);

p1 = poleval((GEN)x[1],y);

p2=poleval((GEN)x[2],y); tetpil=avma;

p2 = poleval((GEN)x[2],y);

return gerepile(av,tetpil,gdiv(p1,p2));

return gerepileupto(av0, gdiv(p1,p2));

case t_VEC: case t_COL:

i=lg(x)-1; imin=1; break;

i = lg(x)-1; imin = 1; break;

default: err(typeer,"poleval");

return NULL; /* not reached */

}

if (i<=imin)

return (i==imin)? gcopy((GEN)x[imin]): gzero;

av=avma; p1=(GEN)x[i]; i--;

lim = stack_lim(av0,2);

p1 = (GEN)x[i]; i--;

if (typ(y)!=t_COMPLEX)

{

#if 0 /* standard Horner's rule */

Line 3063 poleval(GEN x, GEN y)

Line 3274 poleval(GEN x, GEN y)

for (j=i; gcmp0((GEN)x[j]); j--)

if (j==imin)

{

if (i!=j) y = gpuigs(y,i-j+1);

if (i!=j) y = gpowgs(y, i-j+1);

tetpil=avma; return gerepile(av,tetpil,gmul(p1,y));

return gerepileupto(av0, gmul(p1,y));

}

r = (i==j)? y: gpuigs(y,i-j+1);

r = (i==j)? y: gpowgs(y, i-j+1);

p1 = gadd(gmul(p1,r), (GEN)x[j]);

if (low_stack(lim, stack_lim(av0,2)))

{

if (DEBUGMEM>1) err(warnmem,"poleval: i = %ld",i);

p1 = gerepileupto(av0, p1);

}

return gerepileupto(av,p1);

return gerepileupto(av0,p1);

}

p2=(GEN)x[i]; i--; r=gtrace(y); s=gneg_i(gnorm(y));

p2 = (GEN)x[i]; i--; r = gtrace(y); s = gneg_i(gnorm(y));

av = avma;

for ( ; i>=imin; i--)

{

p3=gadd(p2,gmul(r,p1));

GEN p3 = gadd(p2, gmul(r, p1));

p2=gadd((GEN)x[i],gmul(s,p1)); p1=p3;

p2 = gadd((GEN)x[i], gmul(s, p1)); p1 = p3;

if (low_stack(lim, stack_lim(av0,2)))

{

if (DEBUGMEM>1) err(warnmem,"poleval: i = %ld",i);

gerepileall(av, 2, &p1, &p2);

}

p1=gmul(y,p1); tetpil=avma;

return gerepileupto(av0, gadd(p2, gmul(y,p1)));

return gerepile(av,tetpil,gadd(p1,p2));

}

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>