Karim BELABAS on Thu, 21 Jan 1999 15:36:41 +0100 (MET) |
[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]
Re: poldisc |
[Ilya:] > b) I tried to do > > p=a*x^3+b*x^2*y+c*x*y^2+d*y^3 + e*x^2+f*x*y+g*y^2 + h*x+i*y + h > di=poldisc(p) > > and it will not work with 16M stack, and does not finish in half an hour > anyway. > > Should it be *that* resource-intensive? I think I may be able to calculate > given a sheet of paper... PARI's internal representation for polynomials... (highly non-symetrical, optimized for 1 or 2 variables). It is instantaneous with default stack, *IF* you make sure the important variables have high priority: gp> y; p=a*x^3+b*x^2*y+c*x*y^2+d*y^3 + e*x^2+f*x*y+g*y^2 + h*x+i*y + h ^^^ (now x and y are priviledged) > a) It looks like you cannot advice poldisc() which variable you are > interested in. True; easy to add, though. The patch below does what you want (see online help). Karim. *** src/basemath/polarit2.c.orig Tue Jan 19 13:51:09 1999 --- src/basemath/polarit2.c Thu Jan 21 15:20:43 1999 *************** *** 2049,2055 **** } GEN ! discsr(GEN x) { long av=avma,tetpil,tx=typ(x),i; GEN z,p1,p2; --- 2062,2068 ---- } GEN ! poldisc0(GEN x, long v) { long av=avma,tetpil,tx=typ(x),i; GEN z,p1,p2; *************** *** 2058,2064 **** { case t_POL: if (gcmp0(x)) return gzero; ! p1 = subres(x, derivpol(x)); p2 = leading_term(x); if (!gcmp1(p2)) p1 = gdiv(p1,p2); if ((lgef(x)-3) & 2) p1 = gneg(p1); return gerepileupto(av,p1); --- 2071,2080 ---- { case t_POL: if (gcmp0(x)) return gzero; ! if (v < 0) ! p1 = subres(x, derivpol(x)); ! else ! p1 = polresultant0(x, deriv(x,v), v, 0); p2 = leading_term(x); if (!gcmp1(p2)) p1 = gdiv(p1,p2); if ((lgef(x)-3) & 2) p1 = gneg(p1); return gerepileupto(av,p1); *************** *** 2067,2073 **** return stoi(-4); case t_QUAD: case t_POLMOD: ! return discsr((GEN)x[1]); case t_QFR: case t_QFI: p1=sqri((GEN)x[2]); --- 2083,2089 ---- return stoi(-4); case t_QUAD: case t_POLMOD: ! return poldisc0((GEN)x[1], v); case t_QFR: case t_QFI: p1=sqri((GEN)x[2]); *************** *** 2076,2082 **** case t_VEC: case t_COL: case t_MAT: i=lg(x); z=cgetg(i,tx); ! for (i--; i; i--) z[i]=(long)discsr((GEN)x[i]); return z; } err(typeer,"discsr"); --- 2092,2098 ---- case t_VEC: case t_COL: case t_MAT: i=lg(x); z=cgetg(i,tx); ! for (i--; i; i--) z[i]=(long)poldisc0((GEN)x[i], v); return z; } err(typeer,"discsr"); *************** *** 2084,2089 **** --- 2100,2111 ---- } GEN + discsr(GEN x) + { + return poldisc0(x, -1); + } + + GEN reduceddiscsmith(GEN pol) { long av=avma,tetpil,i,j,n; *** src/headers/paridecl.h.orig Tue Jan 19 18:25:58 1999 --- src/headers/paridecl.h Thu Jan 21 15:17:03 1999 *************** *** 997,1002 **** --- 997,1003 ---- GEN nfisincl0(GEN a, GEN b,long flag); GEN nfisisom0(GEN a, GEN b,long flag); GEN nfiso(GEN a, GEN b); + GEN poldisc0(GEN x, long v); GEN polfnf(GEN a, GEN t); GEN polresultant0(GEN x, GEN y,long v,long flag); GEN polsym(GEN x, long n); *** src/language/init.c.orig Fri Jan 15 13:59:44 1999 --- src/language/init.c Thu Jan 21 15:16:44 1999 *************** *** 1495,1501 **** {"polcompositum",25,(void*)polcompositum0,6,"GGD0,L,"}, {"polcyclo",99,(void*)cyclo,7,"LDn"}, {"poldegree",99,(void*)poldegree,7,"lGDn"}, ! {"poldisc",18,(void*)discsr,7,"G"}, {"poldiscreduced",18,(void*)reduceddiscsmith,7,"G"}, {"polgalois",99,(void*)galois,6,"Gp"}, {"polinterpolate",31,(void*)polint,7,"GGDGD&"}, --- 1495,1501 ---- {"polcompositum",25,(void*)polcompositum0,6,"GGD0,L,"}, {"polcyclo",99,(void*)cyclo,7,"LDn"}, {"poldegree",99,(void*)poldegree,7,"lGDn"}, ! {"poldisc",99,(void*)poldisc0,7,"GDn"}, {"poldiscreduced",18,(void*)reduceddiscsmith,7,"G"}, {"polgalois",99,(void*)galois,6,"Gp"}, {"polinterpolate",31,(void*)polint,7,"GGDGD&"}, *** src/language/helpmsg.c.orig Tue Dec 15 16:30:18 1998 --- src/language/helpmsg.c Thu Jan 21 15:18:10 1999 *************** *** 294,300 **** "polcompositum(pol1,pol2,{flag=0}): vector of all possible compositums of the number fields defined by the polynomials pol1 and pol2. If (optional) flag is set (i.e non-null), output for each compositum, not only the compositum polynomial pol, but a vector [pol,al1,al2,k] where al1 (resp. al2) is a root of pol1 (resp. pol2) expressed as a polynomial modulo pol, and a small integer k such that al2+k*al1 is the chosen root of pol", "polcyclo(n,{v=x}): n-th cyclotomic polynomial (in variable v)", "poldegree(x,{v}): degree of the polynomial or rational function x with respect to main variable if v is omitted, with respect to v otherwise. Return -1 if x = 0, and 0 if it's a non-zero scalar", ! "poldisc(x): discriminant of the polynomial x", "poldiscreduced(f): vector of elementary divisors of Z[a]/f'(a)Z[a], where a is a root of the polynomial f", "polgalois(x): Galois group of the polynomial x (see manual for group coding)", "polinterpolate(xa,ya,{x}): polynomial interpolation at x according to data vectors xa, ya", --- 294,300 ---- "polcompositum(pol1,pol2,{flag=0}): vector of all possible compositums of the number fields defined by the polynomials pol1 and pol2. If (optional) flag is set (i.e non-null), output for each compositum, not only the compositum polynomial pol, but a vector [pol,al1,al2,k] where al1 (resp. al2) is a root of pol1 (resp. pol2) expressed as a polynomial modulo pol, and a small integer k such that al2+k*al1 is the chosen root of pol", "polcyclo(n,{v=x}): n-th cyclotomic polynomial (in variable v)", "poldegree(x,{v}): degree of the polynomial or rational function x with respect to main variable if v is omitted, with respect to v otherwise. Return -1 if x = 0, and 0 if it's a non-zero scalar", ! "poldisc(x,{v}): discriminant of the polynomial x, with respect to main variable if v is omitted, with respect to v otherwise", "poldiscreduced(f): vector of elementary divisors of Z[a]/f'(a)Z[a], where a is a root of the polynomial f", "polgalois(x): Galois group of the polynomial x (see manual for group coding)", "polinterpolate(xa,ya,{x}): polynomial interpolation at x according to data vectors xa, ya", -- Karim Belabas email: Karim.Belabas@math.u-psud.fr Dep. de Mathematiques, Bat. 425 Universite Paris-Sud Tel: (00 33) 1 69 15 57 48 F-91405 Orsay (France) Fax: (00 33) 1 69 15 60 19 -- PARI/GP Home Page: http://hasse.mathematik.-tu-muenchen.de/ntsw/pari/