version 1.1.1.1, 2001/10/02 11:17:01 |
version 1.2, 2002/09/11 07:26:48 |
Line 20 Foundation, Inc., 59 Temple Place - Suite 330, Boston, |
|
Line 20 Foundation, Inc., 59 Temple Place - Suite 330, Boston, |
|
/**************************************************************/ |
/**************************************************************/ |
#include "pari.h" |
#include "pari.h" |
#include "parinf.h" |
#include "parinf.h" |
GEN idealhermite_aux(GEN nf, GEN x); |
|
|
|
|
int new_galois_format = 0; |
|
|
|
extern GEN R_from_QR(GEN x, long prec); |
|
extern GEN to_polmod(GEN x, GEN mod); |
|
extern GEN roots_to_pol_r1r2(GEN a, long r1, long v); |
|
extern GEN idealhermite_aux(GEN nf, GEN x); |
|
extern GEN cauchy_bound(GEN p); |
|
extern GEN galoisbig(GEN x, long prec); |
|
extern GEN lllfp_marked(int M, GEN x, long D, long flag, long prec, int gram); |
|
extern GEN lllint_marked(long M, GEN x, long D, int g, GEN *h, GEN *fl, GEN *B); |
|
extern GEN mulmat_pol(GEN A, GEN x); |
|
extern ulong smulss(ulong x, ulong y, ulong *rem); |
|
|
void |
void |
checkrnf(GEN rnf) |
checkrnf(GEN rnf) |
{ |
{ |
if (typ(rnf)!=t_VEC) err(idealer1); |
if (typ(rnf)!=t_VEC || lg(rnf)!=12) err(idealer1); |
if (lg(rnf)!=12) err(idealer1); |
|
} |
} |
|
|
GEN |
GEN |
checkbnf(GEN bnf) |
_checkbnf(GEN bnf) |
{ |
{ |
if (typ(bnf)!=t_VEC) err(idealer1); |
if (typ(bnf) == t_VEC) |
switch (lg(bnf)) |
switch (lg(bnf)) |
|
{ |
|
case 11: return bnf; |
|
case 7: return checkbnf((GEN)bnf[1]); |
|
|
|
case 3: |
|
if (typ(bnf[2])==t_POLMOD) |
|
return checkbnf((GEN)bnf[1]); |
|
} |
|
return NULL; |
|
} |
|
|
|
GEN |
|
_checknf(GEN nf) |
|
{ |
|
if (typ(nf)==t_VEC) |
|
switch(lg(nf)) |
|
{ |
|
case 10: return nf; |
|
case 11: return checknf((GEN)nf[7]); |
|
case 7: return checknf((GEN)nf[1]); |
|
case 3: if (typ(nf[2]) == t_POLMOD) return checknf((GEN)nf[1]); |
|
} |
|
return NULL; |
|
} |
|
|
|
GEN |
|
checkbnf(GEN x) |
|
{ |
|
GEN bnf = _checkbnf(x); |
|
if (!bnf) |
{ |
{ |
case 11: return bnf; |
if (_checknf(x)) err(talker,"please apply bnfinit first"); |
case 10: |
err(idealer1); |
if (typ(bnf[1])==t_POL) |
|
err(talker,"please apply bnfinit first"); |
|
break; |
|
case 7: |
|
return checkbnf((GEN)bnf[1]); |
|
|
|
case 3: |
|
if (typ(bnf[2])==t_POLMOD) |
|
return checkbnf((GEN)bnf[1]); |
|
} |
} |
err(idealer1); |
return bnf; |
return NULL; /* not reached */ |
|
} |
} |
|
|
GEN |
GEN |
Line 61 checkbnf_discard(GEN bnf) |
|
Line 92 checkbnf_discard(GEN bnf) |
|
} |
} |
|
|
GEN |
GEN |
checknf(GEN nf) |
checknf(GEN x) |
{ |
{ |
if (typ(nf)==t_POL) err(talker,"please apply nfinit first"); |
GEN nf = _checknf(x); |
if (typ(nf)!=t_VEC) err(idealer1); |
if (!nf) |
switch(lg(nf)) |
|
{ |
{ |
case 10: return nf; |
if (typ(x)==t_POL) err(talker,"please apply nfinit first"); |
case 11: return checknf((GEN)nf[7]); |
err(idealer1); |
case 7: return checknf((GEN)nf[1]); |
|
case 3: if (typ(nf[2]) == t_POLMOD) return checknf((GEN)nf[1]); |
|
} |
} |
err(idealer1); |
return nf; |
return NULL; /* not reached */ |
|
} |
} |
|
|
void |
void |
Line 126 checkprimeid(GEN id) |
|
Line 153 checkprimeid(GEN id) |
|
} |
} |
|
|
void |
void |
checkprhall(GEN prhall) |
|
{ |
|
if (typ(prhall) != t_VEC || lg(prhall) != 3 || typ(prhall[1]) != t_MAT) |
|
err(talker,"incorrect prhall format"); |
|
} |
|
|
|
void |
|
check_pol_int(GEN x, char *s) |
check_pol_int(GEN x, char *s) |
{ |
{ |
long k = lgef(x)-1; |
long k = lgef(x)-1; |
Line 248 transroot(GEN x, int i, int j) |
|
Line 268 transroot(GEN x, int i, int j) |
|
GEN |
GEN |
tschirnhaus(GEN x) |
tschirnhaus(GEN x) |
{ |
{ |
long a, v = varn(x), av = avma; |
gpmem_t av = avma, av2; |
|
long a, v = varn(x); |
GEN u, p1 = cgetg(5,t_POL); |
GEN u, p1 = cgetg(5,t_POL); |
|
|
if (typ(x)!=t_POL) err(notpoler,"tschirnhaus"); |
if (typ(x)!=t_POL) err(notpoler,"tschirnhaus"); |
Line 260 tschirnhaus(GEN x) |
|
Line 281 tschirnhaus(GEN x) |
|
a = mymyrand() >> randshift; if (a==0) a =1; p1[4]=lstoi(a); |
a = mymyrand() >> randshift; if (a==0) a =1; p1[4]=lstoi(a); |
a = mymyrand() >> (randshift-1); if (a>=4) a-=8; p1[3]=lstoi(a); |
a = mymyrand() >> (randshift-1); if (a>=4) a-=8; p1[3]=lstoi(a); |
a = mymyrand() >> (randshift-1); if (a>=4) a-=8; p1[2]=lstoi(a); |
a = mymyrand() >> (randshift-1); if (a>=4) a-=8; p1[2]=lstoi(a); |
u = caract2(x,p1,v); a=avma; |
u = caract2(x,p1,v); av2=avma; |
} |
} |
while (lgef(srgcd(u,derivpol(u))) > 3); /* while u not separable */ |
while (lgef(srgcd(u,derivpol(u))) > 3); /* while u not separable */ |
if (DEBUGLEVEL>1) |
if (DEBUGLEVEL>1) |
fprintferr("Tschirnhaus transform. New pol: %Z",u); |
fprintferr("Tschirnhaus transform. New pol: %Z",u); |
avma=a; return gerepileupto(av,u); |
avma=av2; return gerepileupto(av,u); |
} |
} |
#undef randshift |
#undef randshift |
|
|
Line 284 gpolcomp(GEN p1, GEN p2) |
|
Line 305 gpolcomp(GEN p1, GEN p2) |
|
return 0; |
return 0; |
} |
} |
|
|
/* assume pol in Z[X] */ |
/* assume pol in Z[X]. Find C, L in Z such that POL = C pol(x / L) monic in Z[X]. |
|
* Return POL and set *ptlead = L */ |
GEN |
GEN |
primitive_pol_to_monic(GEN pol, GEN *ptlead) |
primitive_pol_to_monic(GEN pol, GEN *ptlead) |
{ |
{ |
|
|
return p1; |
return p1; |
} |
} |
|
|
GEN galoisbig(GEN x, long prec); |
static GEN |
GEN roots_to_pol(GEN a, long v); |
roots_to_ZX(GEN z, long *e) |
|
{ |
|
GEN a = roots_to_pol(z,0); |
|
GEN b = grndtoi(greal(a),e); |
|
long e1 = gexpo(gimag(a)); |
|
if (e1 > *e) *e = e1; |
|
return b; |
|
} |
|
|
|
static GEN |
|
_res(long n, long s, long k) |
|
{ |
|
GEN y = cgetg(4, t_VEC); |
|
y[1] = lstoi(n); |
|
y[2] = lstoi(s); |
|
if (!new_galois_format) k = (n == 6 && (k == 6 || k == 2))? 2: 1; |
|
y[3] = lstoi(k); return y; |
|
} |
|
|
GEN |
GEN |
galois(GEN x, long prec) |
galois(GEN x, long prec) |
{ |
{ |
long av=avma,av1,i,j,k,n,f,l,l2,e,e1,pr,ind; |
gpmem_t av = avma, av1; |
GEN x1,p1,p2,p3,p4,p5,p6,w,y,z,ee; |
long i,j,k,n,f,l,l2,e,e1,pr,ind; |
|
GEN x1,p1,p2,p3,p4,p5,w,z,ee; |
static int ind5[20]={2,5,3,4, 1,3,4,5, 1,5,2,4, 1,2,3,5, 1,4,2,3}; |
static int ind5[20]={2,5,3,4, 1,3,4,5, 1,5,2,4, 1,2,3,5, 1,4,2,3}; |
static int ind6[60]={3,5,4,6, 2,6,4,5, 2,3,5,6, 2,4,3,6, 2,5,3,4, |
static int ind6[60]={3,5,4,6, 2,6,4,5, 2,3,5,6, 2,4,3,6, 2,5,3,4, |
1,4,5,6, 1,5,3,6, 1,6,3,4, 1,3,4,5, 1,6,2,5, |
1,4,5,6, 1,5,3,6, 1,6,3,4, 1,3,4,5, 1,6,2,5, |
Line 357 galois(GEN x, long prec) |
|
Line 397 galois(GEN x, long prec) |
|
if (typ(x)!=t_POL) err(notpoler,"galois"); |
if (typ(x)!=t_POL) err(notpoler,"galois"); |
n=degpol(x); if (n<=0) err(constpoler,"galois"); |
n=degpol(x); if (n<=0) err(constpoler,"galois"); |
if (n>11) err(impl,"galois of degree higher than 11"); |
if (n>11) err(impl,"galois of degree higher than 11"); |
x = gdiv(x,content(x)); |
x = primpart(x); |
check_pol_int(x, "galois"); |
check_pol_int(x, "galois"); |
if (gisirreducible(x) != gun) |
if (gisirreducible(x) != gun) |
err(impl,"galois of reducible polynomial"); |
err(impl,"galois of reducible polynomial"); |
|
|
if (n<4) |
if (n<4) |
{ |
{ |
if (n<3) |
if (n == 1) { avma = av; return _res(1, 1,1); } |
{ |
if (n == 2) { avma = av; return _res(2,-1,1); } |
avma=av; y=cgetg(4,t_VEC); |
/* n = 3 */ |
y[1] = (n==1)? un: deux; |
f = carreparfait(ZX_disc(x)); |
y[2]=lnegi(gun); |
avma = av; |
} |
return f? _res(3,1,1): _res(6,-1,2); |
else /* n=3 */ |
|
{ |
|
f=carreparfait(discsr(x)); |
|
avma=av; y=cgetg(4,t_VEC); |
|
if (f) { y[1]=lstoi(3); y[2]=un; } |
|
else { y[1]=lstoi(6); y[2]=lnegi(gun); } |
|
} |
|
y[3]=un; return y; |
|
} |
} |
x1 = x = primitive_pol_to_monic(x,NULL); av1=avma; |
x1 = x = primitive_pol_to_monic(x,NULL); av1=avma; |
if (n>7) return galoisbig(x,prec); |
if (n > 7) return galoisbig(x,prec); |
for(;;) |
for(;;) |
{ |
{ |
|
GEN cb = cauchy_bound(x); |
switch(n) |
switch(n) |
{ |
{ |
case 4: z = cgetg(7,t_VEC); |
case 4: z = cgetg(7,t_VEC); |
|
prec = DEFAULTPREC + (long)((gexpo(cb)*18. / BITS_IN_LONG)); |
for(;;) |
for(;;) |
{ |
{ |
p1=roots(x,prec); |
p1=roots(x,prec); |
Line 395 galois(GEN x, long prec) |
|
Line 429 galois(GEN x, long prec) |
|
z[4] = (long)get_F4(transroot(p1,1,4)); |
z[4] = (long)get_F4(transroot(p1,1,4)); |
z[5] = (long)get_F4(transroot(p1,2,3)); |
z[5] = (long)get_F4(transroot(p1,2,3)); |
z[6] = (long)get_F4(transroot(p1,3,4)); |
z[6] = (long)get_F4(transroot(p1,3,4)); |
p4 = roots_to_pol(z,0); |
p5 = roots_to_ZX(z, &e); |
p5=grndtoi(greal(p4),&e); |
|
e1=gexpo(gimag(p4)); if (e1>e) e=e1; |
|
if (e <= -10) break; |
if (e <= -10) break; |
prec = (prec<<1)-2; |
prec = (prec<<1)-2; |
} |
} |
p6 = modulargcd(p5, derivpol(p5)); |
if (!ZX_is_squarefree(p5)) goto tchi; |
if (typ(p6)==t_POL && lgef(p6)>3) goto tchi; |
|
p2 = (GEN)factor(p5)[1]; |
p2 = (GEN)factor(p5)[1]; |
switch(lg(p2)-1) |
switch(lg(p2)-1) |
{ |
{ |
case 1: f=carreparfait(discsr(x)); avma=av; y=cgetg(4,t_VEC); |
case 1: f = carreparfait(ZX_disc(x)); avma = av; |
y[3]=un; |
return f? _res(12,1,4): _res(24,-1,5); |
if (f) { y[2]=un; y[1]=lstoi(12); return y; } |
|
y[2]=lnegi(gun); y[1]=lstoi(24); return y; |
|
|
|
case 2: avma=av; y=cgetg(4,t_VEC); |
case 2: avma = av; return _res(8,-1,3); |
y[3]=un; y[2]=lnegi(gun); y[1]=lstoi(8); return y; |
|
|
|
case 3: avma=av; y=cgetg(4,t_VEC); |
case 3: avma = av; |
y[1]=lstoi(4); y[3]=un; |
return (degpol(p2[1])==2)? _res(4,1,2): _res(4,-1,1); |
y[2] = (lgef(p2[1])==5)? un: lnegi(gun); |
|
return y; |
|
|
|
default: err(bugparier,"galois (bug1)"); |
default: err(bugparier,"galois (bug1)"); |
} |
} |
|
|
case 5: z = cgetg(7,t_VEC); |
case 5: z = cgetg(7,t_VEC); |
ee = new_chunk(7); w = new_chunk(7); |
ee= cgetg(7,t_VECSMALL); |
|
w = cgetg(7,t_VECSMALL); |
|
prec = DEFAULTPREC + (long)((gexpo(cb)*21. / BITS_IN_LONG)); |
for(;;) |
for(;;) |
{ |
{ |
for(;;) |
for(;;) |
Line 443 galois(GEN x, long prec) |
|
Line 471 galois(GEN x, long prec) |
|
e1 = gexpo(gimag(p3)); if (e1>e) e=e1; |
e1 = gexpo(gimag(p3)); if (e1>e) e=e1; |
ee[l]=e; z[l] = (long)p3; |
ee[l]=e; z[l] = (long)p3; |
} |
} |
p4 = roots_to_pol(z,0); |
p5 = roots_to_ZX(z, &e); |
p5=grndtoi(greal(p4),&e); |
|
e1 = gexpo(gimag(p4)); if (e1>e) e=e1; |
|
if (e <= -10) break; |
if (e <= -10) break; |
prec = (prec<<1)-2; |
prec = (prec<<1)-2; |
} |
} |
p6 = modulargcd(p5,derivpol(p5)); |
if (!ZX_is_squarefree(p5)) goto tchi; |
if (typ(p6)==t_POL && lgef(p6)>3) goto tchi; |
|
p3=(GEN)factor(p5)[1]; |
p3=(GEN)factor(p5)[1]; |
f=carreparfait(discsr(x)); |
f=carreparfait(ZX_disc(x)); |
if (lg(p3)-1==1) |
if (lg(p3)-1==1) |
{ |
{ |
avma=av; y=cgetg(4,t_VEC); y[3]=un; |
avma = av; |
if (f) { y[2]=un; y[1]=lstoi(60); return y; } |
return f? _res(60,1,4): _res(120,-1,5); |
else { y[2]=lneg(gun); y[1]=lstoi(120); return y; } |
|
} |
} |
if (!f) |
if (!f) { avma = av; return _res(20,-1,3); } |
{ |
|
avma=av; y=cgetg(4,t_VEC); |
|
y[3]=un; y[2]=lneg(gun); y[1]=lstoi(20); return y; |
|
} |
|
pr = - (bit_accuracy(prec) >> 1); |
pr = - (bit_accuracy(prec) >> 1); |
for (l=1; l<=6; l++) |
for (l=1; l<=6; l++) |
if (ee[l] <= pr && gcmp0(poleval(p5,(GEN)w[l]))) break; |
if (ee[l] <= pr && gcmp0(poleval(p5,(GEN)w[l]))) break; |
Line 472 galois(GEN x, long prec) |
|
Line 493 galois(GEN x, long prec) |
|
p3=gzero; |
p3=gzero; |
for (i=1; i<=5; i++) |
for (i=1; i<=5; i++) |
{ |
{ |
j=(i%5)+1; |
j = i+1; if (j>5) j -= 5; |
p3=gadd(p3,gmul(gmul((GEN)p2[i],(GEN)p2[j]), |
p3 = gadd(p3,gmul(gmul((GEN)p2[i],(GEN)p2[j]), |
gsub((GEN)p2[j],(GEN)p2[i]))); |
gsub((GEN)p2[j],(GEN)p2[i]))); |
} |
} |
p5=gsqr(p3); p4=grndtoi(greal(p5),&e); |
p5=gsqr(p3); p4=grndtoi(greal(p5),&e); |
e1 = gexpo(gimag(p5)); if (e1>e) e=e1; |
e1 = gexpo(gimag(p5)); if (e1>e) e=e1; |
if (e <= -10) |
if (e <= -10) |
{ |
{ |
if (gcmp0(p4)) goto tchi; |
if (gcmp0(p4)) goto tchi; |
f=carreparfait(p4); avma=av; y=cgetg(4,t_VEC); |
f = carreparfait(p4); avma = av; |
y[3]=y[2]=un; y[1]=lstoi(f?5:10); |
return f? _res(5,1,1): _res(10,1,2); |
return y; |
|
} |
} |
prec=(prec<<1)-2; |
prec=(prec<<1)-2; |
} |
} |
|
|
case 6: z = cgetg(7, t_VEC); |
case 6: z = cgetg(7, t_VEC); |
|
prec = DEFAULTPREC + (long)((gexpo(cb)*42. / BITS_IN_LONG)); |
for(;;) |
for(;;) |
{ |
{ |
for(;;) |
for(;;) |
Line 502 galois(GEN x, long prec) |
|
Line 523 galois(GEN x, long prec) |
|
{ |
{ |
p5=gadd(gmul((GEN)p2[ind6[k]],(GEN)p2[ind6[k+1]]), |
p5=gadd(gmul((GEN)p2[ind6[k]],(GEN)p2[ind6[k+1]]), |
gmul((GEN)p2[ind6[k+2]],(GEN)p2[ind6[k+3]])); |
gmul((GEN)p2[ind6[k+2]],(GEN)p2[ind6[k+3]])); |
p3=gadd(p3,gmul(gsqr(gmul((GEN)p2[i],(GEN)p2[j])),p5)); k+=4; |
p3=gadd(p3,gmul(gsqr(gmul((GEN)p2[i],(GEN)p2[j])),p5)); |
|
k += 4; |
} |
} |
z[l] = (long)p3; |
z[l] = (long)p3; |
} |
} |
p4=roots_to_pol(z,0); |
p5 = roots_to_ZX(z, &e); |
p5=grndtoi(greal(p4),&e); |
|
e1 = gexpo(gimag(p4)); if (e1>e) e=e1; |
|
if (e <= -10) break; |
if (e <= -10) break; |
prec=(prec<<1)-2; |
prec=(prec<<1)-2; |
} |
} |
p6 = modulargcd(p5,derivpol(p5)); |
if (!ZX_is_squarefree(p5)) goto tchi; |
if (typ(p6)==t_POL && lgef(p6)>3) goto tchi; |
|
p2=(GEN)factor(p5)[1]; |
p2=(GEN)factor(p5)[1]; |
switch(lg(p2)-1) |
switch(lg(p2)-1) |
{ |
{ |
Line 530 galois(GEN x, long prec) |
|
Line 549 galois(GEN x, long prec) |
|
gmul(gmul((GEN)p2[4],(GEN)p2[5]),(GEN)p2[6])); |
gmul(gmul((GEN)p2[4],(GEN)p2[5]),(GEN)p2[6])); |
z[++ind] = (long)p3; |
z[++ind] = (long)p3; |
} |
} |
p4 = roots_to_pol(z,0); |
p5 = roots_to_ZX(z, &e); |
p5=grndtoi(greal(p4),&e); |
|
e1 = gexpo(gimag(p4)); if (e1>e) e=e1; |
|
if (e <= -10) |
if (e <= -10) |
{ |
{ |
p6 = modulargcd(p5,derivpol(p5)); |
if (!ZX_is_squarefree(p5)) goto tchi; |
if (typ(p6)==t_POL && lgef(p6)>3) goto tchi; |
p2 = (GEN)factor(p5)[1]; |
p2=(GEN)factor(p5)[1]; |
f = carreparfait(ZX_disc(x)); |
f=carreparfait(discsr(x)); |
avma = av; |
avma=av; y=cgetg(4,t_VEC); y[3]=un; |
|
if (lg(p2)-1==1) |
if (lg(p2)-1==1) |
{ |
return f? _res(360,1,15): _res(720,-1,16); |
if (f) { y[2]=un; y[1]=lstoi(360); } |
|
else { y[2]=lnegi(gun); y[1]=lstoi(720); } |
|
} |
|
else |
else |
{ |
return f? _res(36,1,10): _res(72,-1,13); |
if (f) { y[2]=un; y[1]=lstoi(36); } |
|
else { y[2]=lnegi(gun); y[1]=lstoi(72); } |
|
} |
|
return y; |
|
} |
} |
prec=(prec<<1)-2; break; |
prec=(prec<<1)-2; break; |
|
|
case 2: l2=degpol(p2[1]); if (l2>3) l2=6-l2; |
case 2: l2=degpol(p2[1]); if (l2>3) l2=6-l2; |
switch(l2) |
switch(l2) |
{ |
{ |
case 1: f=carreparfait(discsr(x)); |
case 1: f = carreparfait(ZX_disc(x)); |
avma=av; y=cgetg(4,t_VEC); y[3]=un; |
avma = av; |
if (f) { y[2]=un; y[1]=lstoi(60); } |
return f? _res(60,1,12): _res(120,-1,14); |
else { y[2]=lneg(gun); y[1]=lstoi(120); } |
case 2: f = carreparfait(ZX_disc(x)); |
return y; |
if (f) { avma = av; return _res(24,1,7); } |
case 2: f=carreparfait(discsr(x)); |
p3 = (degpol(p2[1])==2)? (GEN)p2[2]: (GEN)p2[1]; |
if (f) |
f = carreparfait(ZX_disc(p3)); |
{ |
avma = av; |
avma=av; y=cgetg(4,t_VEC); |
return f? _res(24,-1,6): _res(48,-1,11); |
y[3]=y[2]=un; y[1]=lstoi(24); |
case 3: f = carreparfait(ZX_disc((GEN)p2[1])) |
} |
|| carreparfait(ZX_disc((GEN)p2[2])); |
else |
avma = av; |
{ |
return f? _res(18,-1,5): _res(36,-1,9); |
p3=(lgef(p2[1])==5) ? (GEN)p2[2]:(GEN)p2[1]; |
|
f=carreparfait(discsr(p3)); |
|
avma=av; y=cgetg(4,t_VEC); y[2]=lneg(gun); |
|
if (f) { y[1]=lstoi(24); y[3]=deux; } |
|
else { y[1]=lstoi(48); y[3]=un; } |
|
} |
|
return y; |
|
case 3: f=carreparfait(discsr((GEN)p2[1])) |
|
|| carreparfait(discsr((GEN)p2[2])); |
|
avma=av; y=cgetg(4,t_VEC); |
|
y[3]=un; y[2]=lneg(gun); y[1]=lstoi(f? 18: 36); |
|
return y; |
|
} |
} |
case 3: |
case 3: |
for (l2=1; l2<=3; l2++) |
for (l2=1; l2<=3; l2++) |
if (lgef(p2[l2])>=6) p3=(GEN)p2[l2]; |
if (degpol(p2[l2]) >= 3) p3 = (GEN)p2[l2]; |
if (lgef(p3)==6) |
if (degpol(p3) == 3) |
{ |
{ |
f=carreparfait(discsr(p3)); avma=av; y=cgetg(4,t_VEC); |
f = carreparfait(ZX_disc(p3)); avma = av; |
y[2]=lneg(gun); y[1]=lstoi(f? 6: 12); |
return f? _res(6,-1,1): _res(12,-1,3); |
} |
} |
else |
else |
{ |
{ |
f=carreparfait(discsr(x)); avma=av; y=cgetg(4,t_VEC); |
f = carreparfait(ZX_disc(x)); avma = av; |
if (f) { y[2]=un; y[1]=lstoi(12); } |
return f? _res(12,1,4): _res(24,-1,8); |
else { y[2]=lneg(gun); y[1]=lstoi(24); } |
|
} |
} |
y[3]=un; return y; |
case 4: avma = av; return _res(6,-1,2); |
case 4: avma=av; y=cgetg(4,t_VEC); |
|
y[1]=lstoi(6); y[2]=lneg(gun); y[3]=deux; return y; |
|
default: err(bugparier,"galois (bug3)"); |
default: err(bugparier,"galois (bug3)"); |
} |
} |
} |
} |
|
|
case 7: z = cgetg(36,t_VEC); |
case 7: z = cgetg(36,t_VEC); |
|
prec = DEFAULTPREC + (long)((gexpo(cb)*7. / BITS_IN_LONG)); |
for(;;) |
for(;;) |
{ |
{ |
ind = 0; p1=roots(x,prec); p4=gun; |
ind = 0; p1=roots(x,prec); |
for (i=1; i<=5; i++) |
for (i=1; i<=5; i++) |
for (j=i+1; j<=6; j++) |
for (j=i+1; j<=6; j++) |
{ |
{ |
p6 = gadd((GEN)p1[i],(GEN)p1[j]); |
GEN t = gadd((GEN)p1[i],(GEN)p1[j]); |
for (k=j+1; k<=7; k++) |
for (k=j+1; k<=7; k++) z[++ind] = ladd(t, (GEN)p1[k]); |
z[++ind] = (long) gadd(p6,(GEN)p1[k]); |
|
} |
} |
p4 = roots_to_pol(z,0); |
p5 = roots_to_ZX(z, &e); |
p5 = grndtoi(greal(p4),&e); |
|
e1 = gexpo(gimag(p4)); if (e1>e) e=e1; |
|
if (e <= -10) break; |
if (e <= -10) break; |
prec = (prec<<1)-2; |
prec = (prec<<1)-2; |
} |
} |
p6 = modulargcd(p5,derivpol(p5)); |
if (!ZX_is_squarefree(p5)) goto tchi; |
if (typ(p6)==t_POL && lgef(p6)>3) goto tchi; |
|
p2=(GEN)factor(p5)[1]; |
p2=(GEN)factor(p5)[1]; |
switch(lg(p2)-1) |
switch(lg(p2)-1) |
{ |
{ |
case 1: f=carreparfait(discsr(x)); avma=av; y=cgetg(4,t_VEC); y[3]=un; |
case 1: f = carreparfait(ZX_disc(x)); avma = av; |
if (f) { y[2]=un; y[1]=lstoi(2520); } |
return f? _res(2520,1,6): _res(5040,-1,7); |
else { y[2]=lneg(gun); y[1]=lstoi(5040); } |
case 2: f = degpol(p2[1]); avma = av; |
return y; |
return (f==7 || f==28)? _res(168,1,5): _res(42,-1,4); |
case 2: f=degpol(p2[1]); avma=av; y=cgetg(4,t_VEC); y[3]=un; |
case 3: avma = av; return _res(21,1,3); |
if (f==7 || f==28) { y[2]=un; y[1]=lstoi(168); } |
case 4: avma = av; return _res(14,-1,2); |
else { y[2]=lneg(gun); y[1]=lstoi(42); } |
case 5: avma = av; return _res(7,1,1); |
return y; |
|
case 3: avma=av; y=cgetg(4,t_VEC); |
|
y[3]=y[2]=un; y[1]=lstoi(21); return y; |
|
case 4: avma=av; y=cgetg(4,t_VEC); |
|
y[3]=un; y[2]=lneg(gun); y[1]=lstoi(14); return y; |
|
case 5: avma=av; y=cgetg(4,t_VEC); |
|
y[3]=y[2]=un; y[1]=lstoi(7); return y; |
|
default: err(talker,"galois (bug2)"); |
default: err(talker,"galois (bug2)"); |
} |
} |
} |
} |
tchi: avma=av1; x=tschirnhaus(x1); |
tchi: avma = av1; x = tschirnhaus(x1); |
} |
} |
} |
} |
|
|
GEN |
GEN |
galoisapply(GEN nf, GEN aut, GEN x) |
galoisapply(GEN nf, GEN aut, GEN x) |
{ |
{ |
long av=avma,tetpil,lx,j,N; |
gpmem_t av=avma,tetpil; |
|
long lx,j,N; |
GEN p,p1,y,pol; |
GEN p,p1,y,pol; |
|
|
nf=checknf(nf); pol=(GEN)nf[1]; |
nf=checknf(nf); pol=(GEN)nf[1]; |
Line 713 GEN pol_to_monic(GEN pol, GEN *lead); |
|
Line 698 GEN pol_to_monic(GEN pol, GEN *lead); |
|
int cmp_pol(GEN x, GEN y); |
int cmp_pol(GEN x, GEN y); |
|
|
GEN |
GEN |
|
get_bnfpol(GEN x, GEN *bnf, GEN *nf) |
|
{ |
|
*bnf = _checkbnf(x); |
|
*nf = _checknf(x); |
|
if (*nf) return (GEN)(*nf)[1]; |
|
if (typ(x) != t_POL) err(idealer1); |
|
return x; |
|
} |
|
|
|
GEN |
get_nfpol(GEN x, GEN *nf) |
get_nfpol(GEN x, GEN *nf) |
{ |
{ |
if (typ(x) == t_POL) { *nf = NULL; return x; } |
if (typ(x) == t_POL) { *nf = NULL; return x; } |
Line 723 get_nfpol(GEN x, GEN *nf) |
|
Line 718 get_nfpol(GEN x, GEN *nf) |
|
static GEN |
static GEN |
nfiso0(GEN a, GEN b, long fliso) |
nfiso0(GEN a, GEN b, long fliso) |
{ |
{ |
ulong av = avma; |
gpmem_t av = avma; |
long n,m,i,vb,lx; |
long n,m,i,vb,lx; |
GEN nfa,nfb,p1,y,la,lb; |
GEN nfa,nfb,p1,y,la,lb; |
|
|
Line 751 nfiso0(GEN a, GEN b, long fliso) |
|
Line 746 nfiso0(GEN a, GEN b, long fliso) |
|
} |
} |
else |
else |
{ |
{ |
GEN da = nfa? (GEN)nfa[3]: discsr(a); |
GEN da = nfa? (GEN)nfa[3]: ZX_disc(a); |
GEN db = nfb? (GEN)nfb[3]: discsr(b); |
GEN db = nfb? (GEN)nfb[3]: ZX_disc(b); |
if (fliso) |
if (fliso) |
{ |
{ |
p1=gdiv(da,db); |
p1=gdiv(da,db); |
Line 787 nfiso0(GEN a, GEN b, long fliso) |
|
Line 782 nfiso0(GEN a, GEN b, long fliso) |
|
} |
} |
y = gen_sort(y, 0, cmp_pol); |
y = gen_sort(y, 0, cmp_pol); |
settyp(y, t_VEC); |
settyp(y, t_VEC); |
if (vb == 0) delete_var(); |
if (vb == 0) (void)delete_var(); |
} |
} |
lx = lg(y); if (lx==1) { avma=av; return gzero; } |
lx = lg(y); if (lx==1) { avma=av; return gzero; } |
for (i=1; i<lx; i++) |
for (i=1; i<lx; i++) |
Line 801 nfiso0(GEN a, GEN b, long fliso) |
|
Line 796 nfiso0(GEN a, GEN b, long fliso) |
|
} |
} |
|
|
GEN |
GEN |
nfisisom(GEN a, GEN b) |
nfisisom(GEN a, GEN b) { return nfiso0(a,b,1); } |
{ |
|
return nfiso0(a,b,1); |
|
} |
|
|
|
GEN |
GEN |
nfisincl(GEN a, GEN b) |
nfisincl(GEN a, GEN b) { return nfiso0(a,b,0); } |
{ |
|
return nfiso0(a,b,0); |
|
} |
|
|
|
/*************************************************************************/ |
/*************************************************************************/ |
/** **/ |
/** **/ |
/** INITALG **/ |
/** INITALG **/ |
/** **/ |
/** **/ |
/*************************************************************************/ |
/*************************************************************************/ |
extern GEN LLL_nfbasis(GEN *x, GEN polr, GEN base, long prec); |
|
extern GEN eleval(GEN f,GEN h,GEN a); |
|
extern int canon_pol(GEN z); |
|
extern GEN mat_to_vecpol(GEN x, long v); |
|
extern GEN vecpol_to_mat(GEN v, long n); |
|
|
|
|
GEN |
|
get_roots(GEN x,long r1,long prec) |
|
{ |
|
GEN roo = (typ(x)==t_VEC)? dummycopy(x): roots(x,prec); |
|
long i, ru = (lg(roo)-1 + r1) >> 1; |
|
|
|
for (i=1; i<=r1; i++) roo[i]=lreal((GEN)roo[i]); |
|
for ( ; i<=ru; i++) roo[i]=roo[(i<<1)-r1]; |
|
roo[0]=evaltyp(t_VEC)|evallg(ru+1); return roo; |
|
} |
|
|
/* For internal use. compute trace(x mod pol), sym=polsym(pol,deg(pol)-1) */ |
/* For internal use. compute trace(x mod pol), sym=polsym(pol,deg(pol)-1) */ |
GEN |
GEN |
quicktrace(GEN x, GEN sym) |
quicktrace(GEN x, GEN sym) |
Line 843 quicktrace(GEN x, GEN sym) |
|
Line 838 quicktrace(GEN x, GEN sym) |
|
static GEN |
static GEN |
trace_col(GEN x, GEN T) |
trace_col(GEN x, GEN T) |
{ |
{ |
ulong av = avma; |
gpmem_t av = avma; |
GEN t = gzero; |
GEN t = gzero; |
long i, l = lg(x); |
long i, l = lg(x); |
|
|
Line 853 trace_col(GEN x, GEN T) |
|
Line 848 trace_col(GEN x, GEN T) |
|
return gerepileuptoint(av, t); |
return gerepileuptoint(av, t); |
} |
} |
|
|
/* Seek a new, simpler, polynomial pol defining the same number field as |
|
* x (assumed to be monic at this point). |
|
* *ptx receives pol |
|
* *ptdx receives disc(pol) |
|
* *ptrev expresses the new root in terms of the old one. |
|
* base updated in place. |
|
* prec = 0 <==> field is totally real |
|
*/ |
|
static void |
|
nfinit_reduce(long flag, GEN *px, GEN *pdx, GEN *prev, GEN *pbase, long prec) |
|
{ |
|
GEN chbas,a,phimax,dxn,s,sn,p1,p2,p3,polmax,rev,polr; |
|
GEN x = *px, dx = *pdx, base = *pbase; |
|
long i,j,nmax,numb,flc,v=varn(x), n=lg(base)-1; |
|
|
|
if (n == 1) |
|
{ |
|
*px = gsub(polx[v],gun); *pdx = gun; |
|
*prev = polymodrecip(gmodulcp(gun, x)); return; |
|
} |
|
polr = prec? roots(x, prec): NULL; |
|
if (polr) |
|
for (s=gzero,i=1; i<=n; i++) s = gadd(s,gnorm((GEN)polr[i])); |
|
else |
|
s = subii(sqri((GEN)x[n+1]), shifti((GEN)x[n],1)); |
|
|
|
chbas = LLL_nfbasis(&x,polr,base,prec); |
|
if (DEBUGLEVEL) msgtimer("LLL basis"); |
|
|
|
phimax=polx[v]; polmax=dummycopy(x); |
|
nmax=(flag & nf_PARTIAL)?min(n,3):n; |
|
|
|
for (numb=0,i=2; i<=nmax || (!numb && i<=n); i++) |
|
{ /* cf pols_for_polred */ |
|
if (DEBUGLEVEL>2) { fprintferr("i = %ld\n",i); flusherr(); } |
|
p1 = a = gmul(base,(GEN)chbas[i]); p3=content(p1); |
|
if (gcmp1(p3)) p3 = NULL; else p1 = gdiv(p1,p3); |
|
p1 = caract2(x,p1,v); |
|
if (p3) |
|
for (p2=gun, j=lgef(p1)-2; j>=2; j--) |
|
{ |
|
p2 = gmul(p2,p3); p1[j] = lmul((GEN)p1[j], p2); |
|
} |
|
p2 = modulargcd(derivpol(p1),p1); |
|
if (lgef(p2) > 3) continue; |
|
|
|
if (DEBUGLEVEL>3) outerr(p1); |
|
dxn=discsr(p1); flc=absi_cmp(dxn,dx); numb++; |
|
if (flc>0) continue; |
|
|
|
if (polr) |
|
for (sn=gzero,j=1; j<=n; j++) |
|
sn = gadd(sn,gnorm(poleval(a,(GEN)polr[j]))); |
|
else |
|
sn = subii(sqri((GEN)p1[n+1]), shifti((GEN)p1[n],1)); |
|
if (flc>=0) |
|
{ |
|
flc=gcmp(sn,s); |
|
if (flc>0 || (!flc && gpolcomp(p1,polmax) >= 0)) continue; |
|
} |
|
dx=dxn; s=sn; polmax=p1; phimax=a; |
|
} |
|
if (!numb) |
|
{ |
|
if (gisirreducible(x)!=gun) err(redpoler,"nfinit_reduce"); |
|
err(talker,"you have found a counter-example to a conjecture, " |
|
"please send us\nthe polynomial as soon as possible"); |
|
} |
|
if (phimax == polx[v]) /* no improvement */ |
|
rev = gmodulcp(phimax,x); |
|
else |
|
{ |
|
if (canon_pol(polmax) < 0) phimax = gneg_i(phimax); |
|
if (DEBUGLEVEL>1) fprintferr("polmax = %Z\n",polmax); |
|
p2 = gmodulcp(phimax,x); rev = polymodrecip(p2); |
|
a = base; base = cgetg(n+1,t_VEC); |
|
p1 = (GEN)rev[2]; |
|
for (i=1; i<=n; i++) |
|
base[i] = (long)eleval(polmax, (GEN)a[i], p1); |
|
p1 = vecpol_to_mat(base,n); p2=denom(p1); p1=gmul(p2,p1); |
|
p1 = gdiv(hnfmodid(p1,p2), p2); |
|
base = mat_to_vecpol(p1,v); |
|
} |
|
*px=polmax; *pdx=dx; *prev=rev; *pbase = base; |
|
} |
|
|
|
/* pol belonging to Z[x], return a monic polynomial generating the same field |
/* pol belonging to Z[x], return a monic polynomial generating the same field |
* as pol (x-> ax+b)) set lead = NULL if pol was monic (after dividing |
* as pol (x-> ax+b)) set lead = NULL if pol was monic (after dividing |
* by the content), and to to leading coeff otherwise. |
* by the content), and to to leading coeff otherwise. |
|
|
pol_to_monic(GEN pol, GEN *lead) |
pol_to_monic(GEN pol, GEN *lead) |
{ |
{ |
long n = lgef(pol)-1; |
long n = lgef(pol)-1; |
GEN p1; |
|
|
|
if (n==1 || gcmp1((GEN)pol[n])) { *lead = NULL; return pol; } |
if (n==1 || gcmp1((GEN)pol[n])) { *lead = NULL; return pol; } |
|
return primitive_pol_to_monic(primpart(pol), lead); |
p1=content(pol); if (!gcmp1(p1)) pol = gdiv(pol,p1); |
|
return primitive_pol_to_monic(pol,lead); |
|
} |
} |
|
|
/* NB: TI integral, det(TI) = d_K, ideal/dideal = codifferent */ |
/* Let T = (Tr(wi wj)), then T^(-1) Z^n = codifferent = coD/dcoD |
GEN |
* NB: det(T) = dK, TI = dK T^(-1) integral */ |
make_MDI(GEN nf, GEN TI, GEN *ideal, GEN *dideal) |
static GEN |
|
make_MDI(GEN nf, GEN TI, GEN *coD, GEN *dcoD) |
{ |
{ |
GEN c = content(TI); |
GEN c, MDI, dK = (GEN)nf[3]; |
GEN p1; |
|
|
|
*dideal = divii((GEN)nf[3], c); |
TI = Q_primitive_part(TI, &c); |
*ideal = p1 = hnfmodid(gdiv(TI,c), *dideal); |
*dcoD = c? diviiexact(dK, c): dK; |
return gmul(c, ideal_two_elt(nf, p1)); |
*coD = hnfmodid(TI, *dcoD); |
|
MDI = ideal_two_elt(nf, *coD); |
|
return c? gmul(c, MDI): MDI; |
} |
} |
|
|
/* [bas[i]/den[i]]= integer basis. roo = real part of the roots */ |
|
GEN |
|
make_M(GEN basden,GEN roo) |
|
{ |
|
GEN p1,d,M, bas=(GEN)basden[1], den=(GEN)basden[2]; |
|
long i,j, ru = lg(roo), n = lg(bas); |
|
M = cgetg(n,t_MAT); |
|
for (j=1; j<n; j++) |
|
{ |
|
p1=cgetg(ru,t_COL); M[j]=(long)p1; |
|
for (i=1; i<ru; i++) |
|
p1[i]=(long)poleval((GEN)bas[j],(GEN)roo[i]); |
|
} |
|
if (den) |
|
{ |
|
long prec = precision((GEN)roo[1]); |
|
GEN invd,rd = cgetr(prec+1); |
|
for (j=1; j<n; j++) |
|
{ |
|
d = (GEN)den[j]; if (!d) continue; |
|
p1 = (GEN)M[j]; affir(d,rd); invd = ginv(rd); |
|
for (i=1; i<ru; i++) p1[i] = lmul((GEN)p1[i], invd); |
|
} |
|
} |
|
if (DEBUGLEVEL>4) msgtimer("matrix M"); |
|
return M; |
|
} |
|
|
|
GEN |
|
make_MC(long r1,GEN M) |
|
{ |
|
long i,j,av,tetpil, n = lg(M), ru = lg(M[1]); |
|
GEN p1,p2,MC=cgetg(ru,t_MAT); |
|
|
|
for (j=1; j<ru; j++) |
|
{ |
|
p1=cgetg(n,t_COL); MC[j]=(long)p1; |
|
for (i=1; i<n; i++) |
|
{ |
|
av=avma; p2=gconj(gcoeff(M,j,i)); tetpil=avma; |
|
p1[i] = (j<=r1)? (long)p2: lpile(av,tetpil,gmul2n(p2,1)); |
|
} |
|
} |
|
if (DEBUGLEVEL>4) msgtimer("matrix MC"); |
|
return MC; |
|
} |
|
|
|
GEN |
|
get_roots(GEN x,long r1,long ru,long prec) |
|
{ |
|
GEN roo = (typ(x)==t_VEC)? dummycopy(x): roots(x,prec); |
|
long i; |
|
|
|
for (i=1; i<=r1; i++) roo[i]=lreal((GEN)roo[i]); |
|
for ( ; i<=ru; i++) roo[i]=roo[(i<<1)-r1]; |
|
roo[0]=evaltyp(t_VEC)|evallg(ru+1); return roo; |
|
} |
|
|
|
/* assume lg(nf) > 3 && typ(nf) = container [hopefully a genuine nf] */ |
/* assume lg(nf) > 3 && typ(nf) = container [hopefully a genuine nf] */ |
long |
long |
nf_get_r1(GEN nf) |
nf_get_r1(GEN nf) |
Line 1043 nf_get_r2(GEN nf) |
|
Line 893 nf_get_r2(GEN nf) |
|
err(talker,"false nf in nf_get_r2"); |
err(talker,"false nf in nf_get_r2"); |
return itos((GEN)x[2]); |
return itos((GEN)x[2]); |
} |
} |
|
void |
|
nf_get_sign(GEN nf, long *r1, long *r2) |
|
{ |
|
GEN x = (GEN)nf[2]; |
|
if (typ(x) != t_VEC || lg(x) != 3 |
|
|| typ(x[1]) != t_INT || typ(x[2]) != t_INT) |
|
err(talker,"false nf in nf_get_sign"); |
|
*r1 = itos((GEN)x[1]); |
|
*r2 = (degpol(nf[1]) - *r1) >> 1; |
|
} |
|
|
extern GEN mulmat_pol(GEN A, GEN x); |
|
|
|
static GEN |
static GEN |
get_T(GEN mul, GEN x, GEN bas, GEN den) |
get_Tr(GEN mul, GEN x, GEN basden) |
{ |
{ |
|
GEN tr,T,T1,sym, bas = (GEN)basden[1], den = (GEN)basden[2]; |
long i,j,n = lg(bas)-1; |
long i,j,n = lg(bas)-1; |
GEN tr,T,T1,sym; |
|
T = cgetg(n+1,t_MAT); |
T = cgetg(n+1,t_MAT); |
T1 = cgetg(n+1,t_COL); |
T1 = cgetg(n+1,t_COL); |
sym = polsym(x, n-1); |
sym = polsym(x, n-1); |
Line 1076 get_T(GEN mul, GEN x, GEN bas, GEN den) |
|
Line 934 get_T(GEN mul, GEN x, GEN bas, GEN den) |
|
GEN |
GEN |
get_bas_den(GEN bas) |
get_bas_den(GEN bas) |
{ |
{ |
GEN z,d,den, dbas = dummycopy(bas); |
GEN z,b,d,den, dbas = dummycopy(bas); |
long i, c = 0, l = lg(bas); |
long i, l = lg(bas); |
|
int power = 1; |
den = cgetg(l,t_VEC); |
den = cgetg(l,t_VEC); |
for (i=1; i<l; i++) |
for (i=1; i<l; i++) |
{ |
{ |
d = denom(content((GEN)dbas[i])); |
b = Q_remove_denom((GEN)bas[i], &d); |
if (is_pm1(d)) d = NULL; else { dbas[i] = lmul((GEN)dbas[i],d); c++; } |
dbas[i]= (long)b; |
den[i] = (long)d; |
den[i] = (long)d; if (d) power = 0; |
} |
} |
if (!c) den = NULL; /* power basis */ |
if (power) den = NULL; /* power basis */ |
z = cgetg(3,t_VEC); |
z = cgetg(3,t_VEC); |
z[1] = (long)dbas; |
z[1] = (long)dbas; |
z[2] = (long)den; return z; |
z[2] = (long)den; return z; |
Line 1101 _mulii(GEN x, GEN y) |
|
Line 960 _mulii(GEN x, GEN y) |
|
} |
} |
|
|
GEN |
GEN |
get_mul_table(GEN x,GEN basden,GEN invbas,GEN *T) |
get_mul_table(GEN x,GEN basden,GEN invbas) |
{ |
{ |
long i,j, n = degpol(x); |
long i,j, n = degpol(x); |
GEN z,d, mul = cgetg(n*n+1,t_MAT), bas=(GEN)basden[1], den=(GEN)basden[2]; |
GEN z,d, mul = cgetg(n*n+1,t_MAT), bas=(GEN)basden[1], den=(GEN)basden[2]; |
|
|
for (j=1; j<=n*n; j++) mul[j]=lgetg(n+1,t_COL); |
for (j=1; j<=n*n; j++) mul[j]=lgetg(n+1,t_COL); |
for (i=1; i<=n; i++) |
for (i=1; i<=n; i++) |
for (j=i; j<=n; j++) |
for (j=i; j<=n; j++) |
{ |
{ |
ulong av = avma; |
gpmem_t av = avma; |
z = gres(gmul((GEN)bas[j],(GEN)bas[i]), x); |
z = gres(gmul((GEN)bas[j],(GEN)bas[i]), x); |
z = mulmat_pol(invbas, z); /* integral column */ |
z = mulmat_pol(invbas, z); /* integral column */ |
if (den) |
if (den) |
Line 1120 get_mul_table(GEN x,GEN basden,GEN invbas,GEN *T) |
|
Line 979 get_mul_table(GEN x,GEN basden,GEN invbas,GEN *T) |
|
} |
} |
mul[j+(i-1)*n] = mul[i+(j-1)*n] = lpileupto(av,z); |
mul[j+(i-1)*n] = mul[i+(j-1)*n] = lpileupto(av,z); |
} |
} |
if (T) *T = get_T(mul,x,bas,den); |
|
return mul; |
return mul; |
} |
} |
|
|
/* fill mat = nf[5], as well as nf[8] and nf[9] |
/* as get_Tr, mul_table not precomputed */ |
* If (small) only compute a subset (use dummy 0s for the rest) */ |
static GEN |
|
make_Tr(GEN x, GEN w) |
|
{ |
|
long i,j, n = degpol(x); |
|
GEN p1,p2,t,d; |
|
GEN sym = cgetg(n+2,t_VEC); |
|
GEN den = cgetg(n+1,t_VEC); |
|
GEN T = cgetg(n+1,t_MAT); |
|
gpmem_t av; |
|
|
|
sym = polsym(x, n-1); |
|
p1 = get_bas_den(w); |
|
w = (GEN)p1[1]; |
|
den = (GEN)p1[2]; |
|
for (i=1; i<=n; i++) |
|
{ |
|
p1 = cgetg(n+1,t_COL); T[i] = (long)p1; |
|
for (j=1; j<i ; j++) p1[j] = coeff(T,i,j); |
|
for ( ; j<=n; j++) |
|
{ |
|
av = avma; |
|
p2 = gres(gmul((GEN)w[i],(GEN)w[j]), x); |
|
t = quicktrace(p2, sym); |
|
if (den) |
|
{ |
|
d = _mulii((GEN)den[i],(GEN)den[j]); |
|
if (d) t = diviiexact(t, d); |
|
} |
|
p1[j] = lpileuptoint(av, t); |
|
} |
|
} |
|
return T; |
|
} |
|
|
|
/* compute roots so that _absolute_ accuracy of M >= prec [also holds for G] */ |
|
static void |
|
get_roots_for_M(nffp_t *F) |
|
{ |
|
long n, eBD, er, PREC; |
|
|
|
if (F->extraprec < 0) |
|
{ /* not initialized yet */ |
|
n = degpol(F->x); |
|
eBD = 1 + gexpo((GEN)F->basden[1]); |
|
er = 1 + gexpo(F->ro? F->ro: cauchy_bound(F->x)); |
|
if (er < 0) er = 0; |
|
F->extraprec = ((n*er + eBD + (long)log2(n)) >> TWOPOTBITS_IN_LONG); |
|
} |
|
|
|
PREC = F->prec + F->extraprec; |
|
if (F->ro && gprecision((GEN)F->ro[1]) >= PREC) return; |
|
F->ro = get_roots(F->x, F->r1, PREC); |
|
} |
|
|
|
/* [bas[i]/den[i]]= integer basis. roo = real part of the roots */ |
|
static void |
|
make_M(nffp_t *F, int trunc) |
|
{ |
|
GEN bas = (GEN)F->basden[1], den = (GEN)F->basden[2], ro = F->ro; |
|
GEN m, d, M; |
|
long i, j, l = lg(ro), n = lg(bas); |
|
M = cgetg(n,t_MAT); |
|
m = cgetg(l, t_COL); M[1] = (long)m; |
|
for (i=1; i<l; i++) m[i] = un; /* bas[1] = 1 */ |
|
for (j=2; j<n; j++) |
|
{ |
|
m = cgetg(l,t_COL); M[j] = (long)m; |
|
for (i=1; i<l; i++) m[i] = (long)poleval((GEN)bas[j], (GEN)ro[i]); |
|
} |
|
if (den) |
|
{ |
|
GEN invd, rd = cgetr(F->prec + F->extraprec); |
|
for (j=2; j<n; j++) |
|
{ |
|
d = (GEN)den[j]; if (!d) continue; |
|
m = (GEN)M[j]; affir(d,rd); invd = ginv(rd); |
|
for (i=1; i<l; i++) m[i] = lmul((GEN)m[i], invd); |
|
} |
|
} |
|
|
|
if (trunc && gprecision(M) > F->prec) |
|
{ |
|
M = gprec_w(M, F->prec); |
|
F->ro = gprec_w(ro,F->prec); |
|
} |
|
if (DEBUGLEVEL>4) msgtimer("matrix M"); |
|
F->M = M; |
|
} |
|
|
|
/* return G real such that G~ * G = T_2 */ |
|
static void |
|
make_G(nffp_t *F) |
|
{ |
|
GEN G, m, g, r, M = F->M, sqrt2 = gsqrt(gdeux, F->prec + F->extraprec); |
|
long i, j, k, r1 = F->r1, l = lg(M); |
|
|
|
G = cgetg(l, t_MAT); |
|
for (j=1; j<l; j++) |
|
{ |
|
g = cgetg(l, t_COL); |
|
G[j] = (long)g; m = (GEN)M[j]; |
|
for (k=i=1; i<=r1; i++) g[k++] = m[i]; |
|
for ( ; k < l; i++) |
|
{ |
|
r = (GEN)m[i]; |
|
if (typ(r) == t_COMPLEX) |
|
{ |
|
g[k++] = lmpmul(sqrt2, (GEN)r[1]); |
|
g[k++] = lmpmul(sqrt2, (GEN)r[2]); |
|
} |
|
else |
|
{ |
|
g[k++] = lmpmul(sqrt2, r); |
|
g[k++] = zero; |
|
} |
|
} |
|
} |
|
F->G = G; |
|
} |
|
|
|
static void |
|
make_M_G(nffp_t *F, int trunc) |
|
{ |
|
get_roots_for_M(F); |
|
make_M(F, trunc); |
|
make_G(F); |
|
} |
|
|
void |
void |
get_nf_matrices(GEN nf, long small) |
remake_GM(GEN nf, nffp_t *F, long prec) |
{ |
{ |
GEN x=(GEN)nf[1],dK=(GEN)nf[3],index=(GEN)nf[4],ro=(GEN)nf[6],bas=(GEN)nf[7]; |
F->x = (GEN)nf[1]; |
GEN basden,mul,invbas,M,MC,T,MDI,D,TI,A,dA,mat; |
F->ro = NULL; |
long r1 = nf_get_r1(nf), n = lg(bas)-1; |
F->r1 = nf_get_r1(nf); |
|
F->basden = get_bas_den((GEN)nf[7]); |
|
F->extraprec = -1; |
|
F->prec = prec; make_M_G(F, 1); |
|
} |
|
|
mat = cgetg(small? 4: 8,t_VEC); nf[5] = (long)mat; |
static void |
basden = get_bas_den(bas); |
nffp_init(nffp_t *F, nfbasic_t *T, GEN ro, long prec) |
M = make_M(basden,ro); |
{ |
MC = make_MC(r1,M); |
F->x = T->x; |
mat[1]=(long)M; |
F->ro = ro; |
mat[3]=(long)mulmat_real(MC,M); |
F->r1 = T->r1; |
if (small) { nf[8]=nf[9]=mat[2]=zero; return; } |
if (!T->basden) T->basden = get_bas_den(T->bas); |
|
F->basden = T->basden; |
|
F->extraprec = -1; |
|
F->prec = prec; |
|
} |
|
|
invbas = QM_inv(vecpol_to_mat(bas,n), gun); |
static void |
mul = get_mul_table(x,basden,invbas,&T); |
get_nf_fp_compo(nfbasic_t *T, nffp_t *F, GEN ro, long prec) |
|
{ |
|
nffp_init(F,T,ro,prec); |
|
make_M_G(F, 1); |
|
} |
|
|
|
static GEN |
|
get_sign(long r1, long n) |
|
{ |
|
GEN s = cgetg(3, t_VEC); |
|
s[1] = lstoi(r1); |
|
s[2] = lstoi((n - r1) >> 1); return s; |
|
} |
|
|
|
GEN |
|
nfbasic_to_nf(nfbasic_t *T, GEN ro, long prec) |
|
{ |
|
GEN nf = cgetg(10,t_VEC); |
|
GEN x = T->x; |
|
GEN invbas,Tr,D,TI,A,dA, mat = cgetg(8,t_VEC); |
|
nffp_t F; |
|
get_nf_fp_compo(T, &F, ro, prec); |
|
|
|
nf[1] = (long)T->x; |
|
nf[2] = (long)get_sign(T->r1, degpol(T->x)); |
|
nf[3] = (long)T->dK; |
|
nf[4] = (long)T->index; |
|
nf[6] = (long)F.ro; |
|
nf[5] = (long)mat; |
|
nf[7] = (long)T->bas; |
|
|
|
mat[1] = (long)F.M; |
|
mat[2] = (long)F.G; |
|
|
|
invbas = QM_inv(vecpol_to_mat(T->bas, lg(T->bas)-1), gun); |
|
nf[8] = (long)invbas; |
|
nf[9] = (long)get_mul_table(x, F.basden, invbas); |
if (DEBUGLEVEL) msgtimer("mult. table"); |
if (DEBUGLEVEL) msgtimer("mult. table"); |
|
|
nf[8]=(long)invbas; |
|
nf[9]=(long)mul; |
|
|
|
TI = ZM_inv(T, dK); /* dK T^-1 */ |
Tr = get_Tr((GEN)nf[9], x, F.basden); |
MDI = make_MDI(nf,TI, &A, &dA); |
TI = ZM_inv(Tr, T->dK); /* dK T^-1 */ |
mat[6]=(long)TI; |
mat[6] = (long)TI; |
mat[7]=(long)MDI; /* needed in idealinv below */ |
mat[7] = (long)make_MDI(nf,TI, &A, &dA); /* needed in idealinv below */ |
if (is_pm1(index)) |
if (is_pm1(T->index)) |
D = idealhermite_aux(nf, derivpol(x)); |
D = idealhermite_aux(nf, derivpol(x)); |
else |
else |
D = gmul(dA, idealinv(nf, A)); |
D = gmul(dA, idealinv(nf, A)); |
mat[2]=(long)MC; |
mat[3] = zero; /* FIXME: was gram matrix of current mat[2]. Useless */ |
mat[4]=(long)T; |
mat[4] = (long)Tr; |
mat[5]=(long)D; |
mat[5] = (long)D; |
if (DEBUGLEVEL) msgtimer("matrices"); |
if (DEBUGLEVEL) msgtimer("matrices"); |
|
return nf; |
} |
} |
|
|
/* Initialize the number field defined by the polynomial x (in variable v) |
static GEN |
* flag & nf_REGULAR |
hnffromLLL(GEN nf) |
* regular behaviour. |
{ |
* flag & nf_SMALL |
GEN d, x; |
* compute only nf[1] (pol), nf[2] (signature), nf[5][3] (T2) and |
x = vecpol_to_mat((GEN)nf[7], degpol(nf[1])); |
* nf[7] (integer basis), the other components are filled with gzero. |
x = Q_remove_denom(x, &d); |
* flag & nf_REDUCE |
if (!d) return x; /* power basis */ |
* try a polred first. |
return gauss(hnfmodid(x, d), x); |
* flag & nf_PARTIAL |
} |
* do a partial polred, not a polredabs |
|
* flag & nf_ORIG |
|
* do a polred and return [nfinit(x),Mod(a,red)], where |
|
* Mod(a,red)=Mod(v,x) (i.e return the base change). |
|
*/ |
|
|
|
/* here x can be a polynomial, an nf or a bnf */ |
static GEN |
|
nfbasechange(GEN u, GEN x) |
|
{ |
|
long i,lx; |
|
GEN y; |
|
switch(typ(x)) |
|
{ |
|
case t_COL: /* nfelt */ |
|
return gmul(u, x); |
|
|
|
case t_MAT: /* ideal */ |
|
y = dummycopy(x); |
|
lx = lg(x); |
|
for (i=1; i<lx; i++) y[i] = lmul(u, (GEN)y[i]); |
|
break; |
|
|
|
case t_VEC: /* pr */ |
|
checkprimeid(x); |
|
y = dummycopy(x); |
|
y[2] = lmul(u, (GEN)y[2]); |
|
y[5] = lmul(u, (GEN)y[5]); |
|
break; |
|
default: y = x; |
|
} |
|
return y; |
|
} |
|
|
GEN |
GEN |
initalgall0(GEN x, long flag, long prec) |
nffromhnfbasis(GEN nf, GEN x) |
{ |
{ |
GEN lead = NULL,nf,ro,bas,mat,rev,dK,dx,index,res; |
long tx = typ(x); |
long av=avma,n,i,r1,r2,ru,PRECREG; |
gpmem_t av = avma; |
|
GEN u; |
|
if (!is_vec_t(tx)) return gcopy(x); |
|
nf = checknf(nf); |
|
u = hnffromLLL(nf); |
|
return gerepilecopy(av, nfbasechange(u, x)); |
|
} |
|
|
if (DEBUGLEVEL) timer2(); |
GEN |
if (typ(x)==t_POL) |
nftohnfbasis(GEN nf, GEN x) |
|
{ |
|
long tx = typ(x); |
|
gpmem_t av = avma; |
|
GEN u; |
|
if (!is_vec_t(tx)) return gcopy(x); |
|
nf = checknf(nf); |
|
u = ZM_inv(hnffromLLL(nf), gun); |
|
return gerepilecopy(av, nfbasechange(u, x)); |
|
} |
|
|
|
static GEN |
|
get_red_G(nfbasic_t *T, GEN *pro) |
|
{ |
|
GEN G, u, u0 = NULL; |
|
gpmem_t av; |
|
long i, prec, extraprec, n = degpol(T->x); |
|
nffp_t F; |
|
|
|
extraprec = (long) (0.5 * n * sizeof(long) / 8.); |
|
prec = DEFAULTPREC + extraprec; |
|
nffp_init(&F, T, *pro, prec); |
|
av = avma; |
|
for (i=1; ; i++) |
{ |
{ |
n=degpol(x); if (n<=0) err(constpoler,"initalgall0"); |
F.prec = prec; make_M_G(&F, 0); G = F.G; |
check_pol_int(x,"initalg"); |
if (u0) G = gmul(G, u0); |
if (gisirreducible(x) == gzero) err(redpoler,"nfinit"); |
if (DEBUGLEVEL) |
if (!gcmp1((GEN)x[n+2])) |
fprintferr("get_red_G: starting LLL, prec = %ld (%ld + %ld)\n", |
|
prec + F.extraprec, prec, F.extraprec); |
|
if ((u = lllfp_marked(1, G, 100, 2, prec, 0))) |
{ |
{ |
x = pol_to_monic(x,&lead); |
if (typ(u) == t_MAT) break; |
if (!(flag & nf_SMALL)) |
u = (GEN)u[1]; |
{ |
if (u0) u0 = gerepileupto(av, gmul(u0,u)); |
if (!(flag & nf_REDUCE)) |
else u0 = gerepilecopy(av, u); |
err(warner,"non-monic polynomial. Result of the form [nf,c]"); |
|
flag = flag | nf_REDUCE | nf_ORIG; |
|
} |
|
} |
} |
bas = allbase4(x,0,&dK,NULL); |
if (i == MAXITERPOL) err(accurer,"red_T2"); |
if (DEBUGLEVEL) msgtimer("round4"); |
prec = (prec<<1)-2 + (gexpo(u0) >> TWOPOTBITS_IN_LONG); |
dx = discsr(x); r1 = sturm(x); |
F.ro = NULL; |
|
if (DEBUGLEVEL) err(warnprec,"get_red_G", prec); |
} |
} |
|
*pro = F.ro; return u0? gmul(u0,u): u; |
|
} |
|
|
|
/* Return the base change matrix giving an LLL-reduced basis for the |
|
* integer basis of nf(T). |
|
* *ro = roots of x, computed to precision prec [or NULL] */ |
|
static GEN |
|
get_LLL_basis(nfbasic_t *T, GEN *pro) |
|
{ |
|
GEN u; |
|
if (T->r1 != degpol(T->x)) u = get_red_G(T, pro); |
else |
else |
{ |
{ |
i = lg(x); |
u = lllint_marked(1, make_Tr(T->x, T->bas), 100, 1, &u,NULL,NULL); |
if (typ(x) == t_VEC && i<=4 && i>=3 && typ(x[1])==t_POL) |
if (!u) u = idmat(1); |
{ /* polynomial + integer basis */ |
|
bas=(GEN)x[2]; x=(GEN)x[1]; n=degpol(x); |
|
if (typ(bas) == t_MAT) { mat = bas; bas = mat_to_vecpol(mat,varn(x)); } |
|
else mat = vecpol_to_mat(bas,n); |
|
dx = discsr(x); r1 = sturm(x); |
|
dK = gmul(dx, gsqr(det2(mat))); |
|
} |
|
else |
|
{ |
|
nf=checknf(x); |
|
bas=(GEN)nf[7]; x=(GEN)nf[1]; n=degpol(x); |
|
dK=(GEN)nf[3]; dx=mulii(dK, sqri((GEN)nf[4])); |
|
r1 = nf_get_r1(nf); |
|
} |
|
bas[1]=lpolun[varn(x)]; /* it may be gun => SEGV later */ |
|
} |
} |
r2=(n-r1)>>1; ru=r1+r2; |
return u; |
PRECREG = prec + (expi(dK)>>(TWOPOTBITS_IN_LONG+1)) |
} |
+ (long)((sqrt((double)n)+3) / sizeof(long) * 4); |
|
rev = NULL; |
static void |
if (flag & nf_REDUCE) |
set_LLL_basis(nfbasic_t *T, GEN *pro) |
|
{ |
|
T->bas = gmul(T->bas, get_LLL_basis(T, pro)); |
|
T->basden = NULL; /* recompute */ |
|
} |
|
|
|
typedef struct { |
|
GEN xbest, dxbest; |
|
long ind, indmax, indbest; |
|
} ok_pol_t; |
|
|
|
/* is xn better than x ? */ |
|
static int |
|
better_pol(GEN xn, GEN dxn, GEN x, GEN dx) |
|
{ |
|
int fl; |
|
if (!x) return 1; |
|
fl = absi_cmp(dxn, dx); |
|
return (fl < 0 || (!fl && gpolcomp(xn, x) < 0)); |
|
} |
|
|
|
static GEN |
|
ok_pol(void *TT, GEN xn) |
|
{ |
|
ok_pol_t *T = (ok_pol_t*)TT; |
|
GEN dxn; |
|
|
|
if (++T->ind > T->indmax && T->xbest) return T->xbest; |
|
|
|
if (!ZX_is_squarefree(xn)) return (T->ind == T->indmax)? T->xbest: NULL; |
|
if (DEBUGLEVEL>3) outerr(xn); |
|
dxn = ZX_disc(xn); |
|
if (better_pol(xn, dxn, T->xbest, T->dxbest)) |
{ |
{ |
nfinit_reduce(flag, &x, &dx, &rev, &bas, r1==n? 0: prec); |
T->dxbest = dxn; T->xbest = xn; T->indbest = T->ind; |
if (DEBUGLEVEL) msgtimer("polred"); |
|
} |
} |
if (!carrecomplet(divii(dx,dK),&index)) |
if (T->ind >= T->indmax) return T->xbest; |
err(bugparier,"nfinit (incorrect discriminant)"); |
return NULL; |
|
} |
|
|
ro=get_roots(x,r1,ru,PRECREG); |
/* z in Z[X] with positive leading coeff. Set z := z(-X) or z(X) so that |
if (DEBUGLEVEL) msgtimer("roots"); |
* second coeff is > 0. In place. */ |
|
static int |
|
canon_pol(GEN z) |
|
{ |
|
long i,s; |
|
|
nf=cgetg(10,t_VEC); |
for (i=lgef(z)-2; i>=2; i-=2) |
nf[1]=(long)x; |
|
nf[2]=lgetg(3,t_VEC); |
|
mael(nf,2,1)=lstoi(r1); |
|
mael(nf,2,2)=lstoi(r2); |
|
nf[3]=(long)dK; |
|
nf[4]=(long)index; |
|
nf[6]=(long)ro; |
|
nf[7]=(long)bas; |
|
get_nf_matrices(nf, flag & nf_SMALL); |
|
|
|
if (flag & nf_ORIG) |
|
{ |
{ |
if (!rev) err(talker,"bad flag in initalgall0"); |
s = signe(z[i]); |
res = cgetg(3,t_VEC); |
if (s > 0) |
res[1]=(long)nf; nf = res; |
{ |
res[2]=lead? ldiv(rev,lead): (long)rev; |
for (; i>=2; i-=2) z[i]=lnegi((GEN)z[i]); |
|
return -1; |
|
} |
|
if (s) return 1; |
} |
} |
return gerepilecopy(av, nf); |
return 0; |
} |
} |
|
|
|
static GEN _polred(GEN x, GEN a, GEN *pta, FP_chk_fun *CHECK); |
|
|
|
/* Seek a simpler, polynomial pol defining the same number field as |
|
* x (assumed to be monic at this point) */ |
|
static GEN |
|
nfpolred(int part, nfbasic_t *T) |
|
{ |
|
GEN x = T->x, dx = T->dx, a = T->bas; |
|
GEN phi, xbest, dxbest, mat, d, rev; |
|
long i, n = lg(a)-1, v = varn(x); |
|
ok_pol_t O; |
|
FP_chk_fun chk; |
|
|
|
if (degpol(x) == 1) { T->x = gsub(polx[v],gun); return gun; } |
|
|
|
if (!dx) dx = mulii(T->dK, sqri(T->index)); |
|
|
|
O.ind = 0; |
|
O.indmax = part? min(n,3): n; |
|
O.xbest = NULL; |
|
chk.f = &ok_pol; |
|
chk.data = (void*)&O; |
|
if (!_polred(x, a, NULL, &chk)) |
|
err(talker,"you found a counter-example to a conjecture, please report!"); |
|
xbest = O.xbest; dxbest = O.dxbest; |
|
if (!better_pol(xbest, dxbest, x, dx)) return NULL; /* no improvement */ |
|
|
|
/* update T */ |
|
phi = (GEN)a[O.indbest]; |
|
if (canon_pol(xbest) < 0) phi = gneg_i(phi); |
|
if (DEBUGLEVEL>1) fprintferr("xbest = %Z\n",xbest); |
|
rev = modreverse_i(phi, x); |
|
for (i=1; i<=n; i++) a[i] = (long)RX_RXQ_compo((GEN)a[i], rev, xbest); |
|
mat = vecpol_to_mat(Q_remove_denom(a, &d), n); |
|
if (!is_pm1(d)) mat = gdiv(hnfmodid(mat,d), d); else mat = idmat(n); |
|
|
|
(void)carrecomplet(diviiexact(dxbest,T->dK), &(T->index)); |
|
T->bas= mat_to_vecpol(mat, v); |
|
T->dx = dxbest; |
|
T->x = xbest; return rev; |
|
} |
|
|
GEN |
GEN |
initalgred(GEN x, long prec) |
get_nfindex(GEN bas) |
{ |
{ |
return initalgall0(x,nf_REDUCE,prec); |
gpmem_t av = avma; |
|
long n = lg(bas)-1; |
|
GEN d, mat = vecpol_to_mat(Q_remove_denom(bas, &d), n); |
|
if (!d) { avma = av; return gun; } |
|
return gerepileuptoint(av, diviiexact(gpowgs(d, n), det(mat))); |
} |
} |
|
|
|
void |
|
nfbasic_init(GEN x, long flag, GEN fa, nfbasic_t *T) |
|
{ |
|
GEN bas, dK, dx, index; |
|
long r1; |
|
|
|
T->basden = NULL; |
|
T->lead = NULL; |
|
if (DEBUGLEVEL) (void)timer2(); |
|
if (typ(x) == t_POL) |
|
{ |
|
check_pol_int(x, "nfinit"); |
|
if (gisirreducible(x) == gzero) err(redpoler, "nfinit"); |
|
x = pol_to_monic(x, &(T->lead)); |
|
bas = allbase(x, flag, &dx, &dK, &index, &fa); |
|
if (DEBUGLEVEL) msgtimer("round4"); |
|
r1 = sturm(x); |
|
} |
|
else if (typ(x) == t_VEC && lg(x)==3 && typ(x[1])==t_POL) |
|
{ /* monic integral polynomial + integer basis */ |
|
GEN mat; |
|
bas = (GEN)x[2]; x = (GEN)x[1]; |
|
if (typ(bas) == t_MAT) |
|
{ mat = bas; bas = mat_to_vecpol(mat,varn(x)); } |
|
else |
|
mat = vecpol_to_mat(bas, lg(bas)-1); |
|
index = get_nfindex(bas); |
|
dx = ZX_disc(x); |
|
dK = diviiexact(dx, sqri(index)); |
|
r1 = sturm(x); |
|
} |
|
else |
|
{ /* nf, bnf, bnr */ |
|
GEN nf = checknf(x); |
|
x = (GEN)nf[1]; |
|
dK = (GEN)nf[3]; |
|
index = (GEN)nf[4]; |
|
bas = (GEN)nf[7]; |
|
r1 = nf_get_r1(nf); dx = NULL; |
|
} |
|
|
|
T->x = x; |
|
T->r1 = r1; |
|
T->dx = dx; |
|
T->dK = dK; |
|
T->bas = bas; |
|
T->index = index; |
|
} |
|
|
|
/* Initialize the number field defined by the polynomial x (in variable v) |
|
* flag & nf_RED: try a polred first. |
|
* flag & nf_PARTRED: as nf_RED but check the first two elements only. |
|
* flag & nf_ORIG |
|
* do a polred and return [nfinit(x), Mod(a,red)], where |
|
* Mod(a,red) = Mod(v,x) (i.e return the base change). */ |
GEN |
GEN |
initalgred2(GEN x, long prec) |
_initalg(GEN x, long flag, long prec) |
{ |
{ |
return initalgall0(x,nf_REDUCE|nf_ORIG,prec); |
const gpmem_t av = avma; |
|
GEN nf, rev = NULL, ro = NULL; |
|
nfbasic_t T; |
|
|
|
nfbasic_init(x, flag, NULL, &T); |
|
|
|
if (T.lead && !(flag & (nf_RED|nf_PARTRED))) |
|
{ |
|
err(warner,"non-monic polynomial. Result of the form [nf,c]"); |
|
flag |= nf_PARTRED | nf_ORIG; |
|
} |
|
if (flag & (nf_RED|nf_PARTRED)) |
|
{ |
|
set_LLL_basis(&T, &ro); |
|
rev = nfpolred(flag & nf_PARTRED, &T); |
|
if (rev) ro = NULL; /* changed T.x */ |
|
if (flag & nf_ORIG) |
|
{ |
|
if (!rev) rev = polx[varn(T.x)]; /* no improvement */ |
|
if (T.lead) rev = gdiv(rev, T.lead); |
|
rev = to_polmod(rev, T.x); |
|
} |
|
if (DEBUGLEVEL) msgtimer("polred"); |
|
} |
|
|
|
set_LLL_basis(&T, &ro); |
|
if (DEBUGLEVEL) msgtimer("LLL basis"); |
|
|
|
nf = nfbasic_to_nf(&T, ro, prec); |
|
if (flag & nf_ORIG) |
|
{ |
|
GEN res = cgetg(3,t_VEC); |
|
res[1] = (long)nf; |
|
res[2] = (long)rev; nf = res; |
|
} |
|
return gerepilecopy(av, nf); |
} |
} |
|
|
GEN |
GEN |
|
initalgred(GEN x, long prec) { return _initalg(x, nf_RED, prec); } |
|
GEN |
|
initalgred2(GEN x, long prec) { return _initalg(x, nf_RED|nf_ORIG, prec); } |
|
GEN |
|
initalg(GEN x, long prec) { return _initalg(x, 0, prec); } |
|
|
|
GEN |
nfinit0(GEN x, long flag,long prec) |
nfinit0(GEN x, long flag,long prec) |
{ |
{ |
switch(flag) |
switch(flag) |
{ |
{ |
case 0: |
case 0: |
case 1: return initalgall0(x,nf_REGULAR,prec); |
case 1: return _initalg(x,0,prec); |
case 2: return initalgall0(x,nf_REDUCE,prec); |
case 2: return _initalg(x,nf_RED,prec); |
case 3: return initalgall0(x,nf_REDUCE|nf_ORIG,prec); |
case 3: return _initalg(x,nf_RED|nf_ORIG,prec); |
case 4: return initalgall0(x,nf_REDUCE|nf_PARTIAL,prec); |
case 4: return _initalg(x,nf_PARTRED,prec); |
case 5: return initalgall0(x,nf_REDUCE|nf_ORIG|nf_PARTIAL,prec); |
case 5: return _initalg(x,nf_PARTRED|nf_ORIG,prec); |
case 6: return initalgall0(x,nf_SMALL,prec); |
|
default: err(flagerr,"nfinit"); |
default: err(flagerr,"nfinit"); |
} |
} |
return NULL; /* not reached */ |
return NULL; /* not reached */ |
} |
} |
|
|
GEN |
|
initalg(GEN x, long prec) |
|
{ |
|
return initalgall0(x,nf_REGULAR,prec); |
|
} |
|
|
|
/* assume x a bnr/bnf/nf */ |
/* assume x a bnr/bnf/nf */ |
long |
long |
nfgetprec(GEN x) |
nfgetprec(GEN x) |
{ |
{ |
GEN nf = checknf(x), ro = (GEN)nf[6]; |
GEN nf = checknf(x), ro = (GEN)nf[6]; |
return (typ(ro)==t_VEC)?precision((GEN)ro[1]): DEFAULTPREC; |
return (typ(ro)==t_VEC)? precision((GEN)ro[1]): DEFAULTPREC; |
} |
} |
|
|
GEN |
GEN |
nfnewprec(GEN nf, long prec) |
nfnewprec(GEN nf, long prec) |
{ |
{ |
long av=avma,r1,r2,ru,n; |
gpmem_t av = avma; |
GEN y,pol,ro,basden,MC,mat,M; |
GEN NF; |
|
nffp_t F; |
|
|
if (typ(nf) != t_VEC) err(talker,"incorrect nf in nfnewprec"); |
if (typ(nf) != t_VEC) err(talker,"incorrect nf in nfnewprec"); |
if (lg(nf) == 11) return bnfnewprec(nf,prec); |
if (lg(nf) == 11) return bnfnewprec(nf,prec); |
if (lg(nf) == 7) return bnrnewprec(nf,prec); |
if (lg(nf) == 7) return bnrnewprec(nf,prec); |
(void)checknf(nf); |
(void)checknf(nf); |
y = dummycopy(nf); |
NF = dummycopy(nf); |
pol=(GEN)nf[1]; n=degpol(pol); |
NF[5] = (long)dummycopy((GEN)NF[5]); |
r1 = nf_get_r1(nf); |
remake_GM(NF, &F, prec); |
r2 = (n - r1) >> 1; ru = r1+r2; |
NF[6] = (long)F.ro; |
mat=dummycopy((GEN)nf[5]); |
mael(NF,5,1) = (long)F.M; |
ro=get_roots(pol,r1,ru,prec); |
mael(NF,5,2) = (long)F.G; |
y[5]=(long)mat; |
return gerepilecopy(av, NF); |
y[6]=(long)ro; |
} |
basden = get_bas_den((GEN)nf[7]); |
|
M = make_M(basden,ro); |
/********************************************************************/ |
MC = make_MC(r1,M); |
/** **/ |
mat[1]=(long)M; |
/** POLRED **/ |
if (typ(nf[8]) != t_INT) mat[2]=(long)MC; /* not a small nf */ |
/** **/ |
mat[3]=(long)mulmat_real(MC,M); |
/********************************************************************/ |
|
/* remove duplicate polynomials in y, updating a (same indexes), in place */ |
|
static void |
|
remove_duplicates(GEN y, GEN a) |
|
{ |
|
long k, i, l = lg(y); |
|
gpmem_t av = avma; |
|
GEN z; |
|
|
|
if (l < 2) return; |
|
z = new_chunk(3); |
|
z[1] = (long)y; |
|
z[2] = (long)a; (void)sort_factor(z, gcmp); |
|
for (k=1, i=2; i<l; i++) |
|
if (!gegal((GEN)y[i], (GEN)y[k])) |
|
{ |
|
k++; |
|
a[k] = a[i]; |
|
y[k] = y[i]; |
|
} |
|
l = k+1; setlg(a,l); setlg(y,l); |
|
avma = av; |
|
} |
|
|
|
/* if CHECK != NULL, return the first polynomial pol found such that |
|
* CHECK->f(CHECK->data, pol) != 0; return NULL if there are no such pol */ |
|
static GEN |
|
_polred(GEN x, GEN a, GEN *pta, FP_chk_fun *CHECK) |
|
{ |
|
long i, v = varn(x), l = lg(a); |
|
GEN ch, d, y = cgetg(l,t_VEC); |
|
|
|
for (i=1; i<l; i++) |
|
{ |
|
if (DEBUGLEVEL>2) { fprintferr("i = %ld\n",i); flusherr(); } |
|
ch = QX_caract(x, (GEN)a[i], v); |
|
if (CHECK) |
|
{ |
|
ch = CHECK->f(CHECK->data, ch); |
|
if (!ch) continue; |
|
return ch; |
|
} |
|
d = modulargcd(derivpol(ch), ch); |
|
if (degpol(d)) ch = gdivexact(ch,d); |
|
|
|
if (canon_pol(ch) < 0 && pta) a[i] = (long)gneg_i((GEN)a[i]); |
|
if (DEBUGLEVEL > 3) outerr(ch); |
|
y[i] = (long)ch; |
|
} |
|
if (CHECK) return NULL; /* no suitable polynomial found */ |
|
remove_duplicates(y,a); |
|
if (pta) *pta = a; |
|
return y; |
|
} |
|
|
|
static GEN |
|
allpolred(GEN x, long flag, GEN fa, GEN *pta, FP_chk_fun *CHECK) |
|
{ |
|
GEN ro = NULL; |
|
nfbasic_t T; nfbasic_init(x, flag, fa, &T); |
|
if (T.lead) err(impl,"polred for non-monic polynomial"); |
|
set_LLL_basis(&T, &ro); |
|
return _polred(T.x, T.bas, pta, CHECK); |
|
} |
|
|
|
/* FIXME: backward compatibility */ |
|
#define red_PARTIAL 1 |
|
#define red_ORIG 2 |
|
|
|
GEN |
|
polred0(GEN x, long flag, GEN fa) |
|
{ |
|
gpmem_t av = avma; |
|
GEN y, a; |
|
int fl = 0; |
|
|
|
if (fa && gcmp0(fa)) fa = NULL; /* compatibility */ |
|
if (flag & red_PARTIAL) fl |= nf_PARTIALFACT; |
|
if (flag & red_ORIG) fl |= nf_ORIG; |
|
y = allpolred(x, fl, fa, &a, NULL); |
|
if (fl & nf_ORIG) { |
|
GEN z = cgetg(3,t_MAT); |
|
z[1] = (long)a; |
|
z[2] = (long)y; y = z; |
|
} |
return gerepilecopy(av, y); |
return gerepilecopy(av, y); |
} |
} |
|
|
|
GEN |
|
ordred(GEN x) |
|
{ |
|
gpmem_t av = avma; |
|
GEN y; |
|
|
|
if (typ(x) != t_POL) err(typeer,"ordred"); |
|
if (!gcmp1(leading_term(x))) err(impl,"ordred"); |
|
if (!signe(x)) return gcopy(x); |
|
y = cgetg(3,t_VEC); |
|
y[1] = (long)x; |
|
y[2] = (long)idmat(degpol(x)); |
|
return gerepileupto(av, polred(y)); |
|
} |
|
|
|
GEN |
|
polredfirstpol(GEN x, long flag, FP_chk_fun *CHECK) |
|
{ return allpolred(x,flag,NULL,NULL,CHECK); } |
|
GEN |
|
polred(GEN x) |
|
{ return polred0(x, 0, NULL); } |
|
GEN |
|
smallpolred(GEN x) |
|
{ return polred0(x, red_PARTIAL, NULL); } |
|
GEN |
|
factoredpolred(GEN x, GEN fa) |
|
{ return polred0(x, 0, fa); } |
|
GEN |
|
polred2(GEN x) |
|
{ return polred0(x, red_ORIG, NULL); } |
|
GEN |
|
smallpolred2(GEN x) |
|
{ return polred0(x, red_PARTIAL|red_ORIG, NULL); } |
|
GEN |
|
factoredpolred2(GEN x, GEN fa) |
|
{ return polred0(x, red_PARTIAL, fa); } |
|
|
|
/********************************************************************/ |
|
/** **/ |
|
/** POLREDABS **/ |
|
/** **/ |
|
/********************************************************************/ |
|
|
|
GEN |
|
T2_from_embed_norm(GEN x, long r1) |
|
{ |
|
GEN p = sum(x, 1, r1); |
|
GEN q = sum(x, r1+1, lg(x)-1); |
|
if (q != gzero) p = gadd(p, gmul2n(q,1)); |
|
return p; |
|
} |
|
|
|
GEN |
|
T2_from_embed(GEN x, long r1) |
|
{ |
|
return T2_from_embed_norm(gnorm(x), r1); |
|
} |
|
|
|
typedef struct { |
|
long r1, v, prec; |
|
GEN ZKembed; /* embeddings of fincke-pohst-reduced Zk basis */ |
|
GEN u; /* matrix giving fincke-pohst-reduced Zk basis */ |
|
GEN M; /* embeddings of initial (LLL-reduced) Zk basis */ |
|
GEN bound; /* T2 norm of the polynomial defining nf */ |
|
} CG_data; |
|
|
|
/* characteristic pol of x */ |
|
static GEN |
|
get_polchar(CG_data *d, GEN x) |
|
{ |
|
GEN g = gmul(d->ZKembed,x); |
|
long e; |
|
g = grndtoi(roots_to_pol_r1r2(g, d->r1, d->v), &e); |
|
if (e > -5) err(precer, "get_polchar"); |
|
return g; |
|
} |
|
|
|
/* return a defining polynomial for Q(x) */ |
|
static GEN |
|
get_polmin(CG_data *d, GEN x) |
|
{ |
|
GEN g = get_polchar(d,x); |
|
GEN h = modulargcd(g,derivpol(g)); |
|
if (degpol(h)) g = gdivexact(g,h); |
|
return g; |
|
} |
|
|
|
/* does x generate the correct field ? */ |
|
static GEN |
|
chk_gen(void *data, GEN x) |
|
{ |
|
gpmem_t av = avma; |
|
GEN g = get_polchar((CG_data*)data,x); |
|
GEN h = modulargcd(g,derivpol(g)); |
|
if (degpol(h)) { avma = av; return NULL; } |
|
if (DEBUGLEVEL>3) fprintferr(" generator: %Z\n",g); |
|
return g; |
|
} |
|
|
|
/* mat = base change matrix, gram = LLL-reduced T2 matrix */ |
|
static GEN |
|
chk_gen_init(FP_chk_fun *chk, GEN gram, GEN mat) |
|
{ |
|
CG_data *d = (CG_data*)chk->data; |
|
GEN P,bound,prev,x,B; |
|
long l = lg(gram), N = l-1,i,prec,prec2; |
|
int skipfirst = 0; |
|
|
|
d->u = mat; |
|
d->ZKembed = gmul(d->M, mat); |
|
|
|
bound = d->bound; |
|
prev = NULL; |
|
x = zerocol(N); |
|
for (i=1; i<l; i++) |
|
{ |
|
B = gcoeff(gram,i,i); |
|
if (gcmp(B,bound) >= 0) continue; /* don't assume increasing norms */ |
|
|
|
x[i] = un; P = get_polmin(d,x); |
|
x[i] = zero; |
|
if (degpol(P) == N) |
|
{ |
|
if (gcmp(B,bound) < 0) bound = B ; |
|
continue; |
|
} |
|
if (DEBUGLEVEL>2) fprintferr("chk_gen_init: subfield %Z\n",P); |
|
if (skipfirst != i-1) continue; |
|
|
|
/* Q(w[1],...,w[i-1]) = Q[X]/(prev) is a strict subfield of nf */ |
|
if (prev && degpol(prev) < N && !gegal(prev,P)) |
|
{ |
|
if (degpol(prev) * degpol(P) > 64) continue; /* too expensive */ |
|
P = compositum(prev,P); |
|
P = (GEN)P[lg(P)-1]; |
|
if (degpol(P) == N) continue; |
|
if (DEBUGLEVEL>2 && degpol(P) != degpol(prev)) |
|
fprintferr("chk_gen_init: subfield %Z\n",P); |
|
} |
|
skipfirst++; prev = P; |
|
} |
|
/* x_1,...,x_skipfirst generate a strict subfield [unless N=skipfirst=1] */ |
|
chk->skipfirst = skipfirst; |
|
if (DEBUGLEVEL>2) fprintferr("chk_gen_init: skipfirst = %ld\n",skipfirst); |
|
|
|
/* should be gexpo( max_k C^n_k (bound/k)^(k/2) ) */ |
|
prec2 = (1 + (((gexpo(bound)*N)/2) >> TWOPOTBITS_IN_LONG)); |
|
if (prec2 < 0) prec2 = 0; |
|
prec = 3 + prec2; |
|
if (DEBUGLEVEL) |
|
fprintferr("chk_gen_init: new prec = %ld (initially %ld)\n", prec, d->prec); |
|
if (prec > d->prec) err(bugparier, "precision problem in polredabs"); |
|
if (prec < d->prec) d->ZKembed = gprec_w(d->ZKembed, prec); |
|
return bound; |
|
} |
|
|
|
/* store phi(beta mod z). Suitable for gerepileupto */ |
|
static GEN |
|
storeeval(GEN beta, GEN z, GEN lead) |
|
{ |
|
GEN y = cgetg(3,t_VEC); |
|
z = gcopy(z); |
|
beta = lead? gdiv(beta, lead): gcopy(beta); |
|
y[1] = (long)z; |
|
y[2] = (long)to_polmod(beta,z); |
|
return y; |
|
} |
|
|
|
static GEN |
|
storeraw(GEN beta, GEN z) |
|
{ |
|
GEN y = cgetg(3,t_VEC); |
|
y[1] = lcopy(z); |
|
y[2] = lcopy(beta); return y; |
|
} |
|
|
|
static void |
|
findmindisc(GEN *py, GEN *pa) |
|
{ |
|
GEN v, dmin, z, b, discs, y = *py, a = *pa; |
|
long i,k, l = lg(y); |
|
|
|
if (l == 2) { *py = (GEN)y[1]; *pa = (GEN)a[1]; return; } |
|
|
|
discs = cgetg(l,t_VEC); |
|
for (i=1; i<l; i++) discs[i] = labsi(ZX_disc((GEN)y[i])); |
|
v = sindexsort(discs); |
|
k = v[1]; |
|
dmin = (GEN)discs[k]; |
|
z = (GEN)y[k]; |
|
b = (GEN)a[k]; |
|
for (i=2; i<l; i++) |
|
{ |
|
k = v[i]; |
|
if (!egalii((GEN)discs[k], dmin)) break; |
|
if (gpolcomp((GEN)y[k],z) < 0) { z = (GEN)y[k]; b = (GEN)a[k]; } |
|
} |
|
*py = z; *pa = b; |
|
} |
|
|
|
static GEN |
|
storepol(GEN x, GEN z, GEN a, GEN lead, long flag) |
|
{ |
|
GEN y, b; |
|
if (flag & nf_RAW) |
|
y = storeraw(a, z); |
|
else if (flag & nf_ORIG) |
|
{ |
|
b = modreverse_i(a, x); |
|
y = storeeval(b,z, lead); |
|
} |
|
else |
|
y = gcopy(z); |
|
return y; |
|
} |
|
|
|
static GEN |
|
storeallpol(GEN x, GEN z, GEN a, GEN lead, long flag) |
|
{ |
|
GEN y, b; |
|
|
|
if (flag & nf_RAW) |
|
{ |
|
long i, c = lg(z); |
|
y = cgetg(c,t_VEC); |
|
for (i=1; i<c; i++) y[i] = (long)storeraw((GEN)a[i], (GEN)z[i]); |
|
} |
|
else if (flag & nf_ORIG) |
|
{ |
|
long i, c = lg(z); |
|
b = cgetg(c, t_VEC); |
|
for (i=1; i<c; i++) |
|
b[i] = (long)modreverse_i((GEN)a[i], x); |
|
|
|
y = cgetg(c,t_VEC); |
|
for (i=1; i<c; i++) y[i] = (long)storeeval((GEN)b[i], (GEN)z[i], lead); |
|
} |
|
else |
|
y = gcopy(z); |
|
return y; |
|
} |
|
|
|
GEN |
|
_polredabs(nfbasic_t *T, GEN *u) |
|
{ |
|
long i, prec, e, n = degpol(T->x); |
|
GEN v, ro = NULL; |
|
FP_chk_fun chk = { &chk_gen, &chk_gen_init, NULL, 0 }; |
|
nffp_t F; |
|
CG_data d; chk.data = (void*)&d; |
|
|
|
set_LLL_basis(T, &ro); |
|
if (DEBUGLEVEL) msgtimer("LLL basis"); |
|
|
|
/* || polchar ||_oo < 2^e */ |
|
e = n * (gexpo(gmulsg(n, cauchy_bound(T->x))) + 1); |
|
prec = DEFAULTPREC + (e >> TWOPOTBITS_IN_LONG); |
|
|
|
get_nf_fp_compo(T, &F, ro, prec); |
|
|
|
d.v = varn(T->x); |
|
d.r1= T->r1; |
|
d.bound = gmul(T2_from_embed(F.ro, d.r1), dbltor(1.00000001)); |
|
for (i=1; ; i++) |
|
{ |
|
GEN R = R_from_QR(F.G, prec); |
|
d.prec = prec; |
|
d.M = F.M; |
|
if (R) |
|
{ |
|
v = fincke_pohst(_vec(R),NULL,5000,1, 0, &chk); |
|
if (v) break; |
|
} |
|
if (i==MAXITERPOL) err(accurer,"polredabs0"); |
|
prec = (prec<<1)-2; |
|
get_nf_fp_compo(T, &F, NULL, prec); |
|
if (DEBUGLEVEL) err(warnprec,"polredabs0",prec); |
|
} |
|
*u = d.u; return v; |
|
} |
|
|
|
GEN |
|
polredabs0(GEN x, long flag) |
|
{ |
|
long i, l, vx; |
|
gpmem_t av = avma; |
|
GEN y, a, u; |
|
nfbasic_t T; |
|
|
|
nfbasic_init(x, flag & nf_PARTIALFACT, NULL, &T); |
|
x = T.x; vx = varn(x); |
|
|
|
if (degpol(x) == 1) |
|
{ |
|
u = NULL; |
|
y = _vec(polx[vx]); |
|
a = _vec(gsub((GEN)y[1], (GEN)x[2])); |
|
} |
|
else |
|
{ |
|
GEN v = _polredabs(&T, &u); |
|
y = (GEN)v[1]; |
|
a = (GEN)v[2]; |
|
} |
|
l = lg(a); |
|
for (i=1; i<l; i++) |
|
if (canon_pol((GEN)y[i]) < 0) a[i] = (long)gneg_i((GEN)a[i]); |
|
remove_duplicates(y,a); |
|
l = lg(a); |
|
if (l == 1) |
|
{ |
|
y = _vec(x); |
|
a = _vec(polx[vx]); |
|
} |
|
if (DEBUGLEVEL) fprintferr("%ld minimal vectors found.\n",l-1); |
|
if (flag & nf_ALL) |
|
{ |
|
if (u) for (i=1; i<l; i++) a[i] = lmul(T.bas, gmul(u, (GEN)a[i])); |
|
y = storeallpol(x,y,a,T.lead,flag); |
|
} |
|
else |
|
{ |
|
findmindisc(&y, &a); |
|
if (u) a = gmul(T.bas, gmul(u, a)); |
|
y = storepol(x,y,a,T.lead,flag); |
|
} |
|
return gerepileupto(av, y); |
|
} |
|
|
|
GEN |
|
polredabsall(GEN x, long flun) { return polredabs0(x, flun | nf_ALL); } |
|
GEN |
|
polredabs(GEN x) { return polredabs0(x,0); } |
|
GEN |
|
polredabs2(GEN x) { return polredabs0(x,nf_ORIG); } |
|
|
static long |
static long |
nf_pm1(GEN y) |
nf_pm1(GEN y) |
{ |
{ |
|
|
for (i=2; i<l; i++) |
for (i=2; i<l; i++) |
if (signe(y[i])) return 0; |
if (signe(y[i])) return 0; |
return signe(y[1]); |
return signe(y[1]); |
|
|
} |
} |
|
|
static GEN |
static GEN |
Line 1362 is_primitive_root(GEN nf, GEN fa, GEN x, long w) |
|
Line 2019 is_primitive_root(GEN nf, GEN fa, GEN x, long w) |
|
return x; |
return x; |
} |
} |
|
|
|
/* only roots of 1 are +/- 1 */ |
|
static GEN |
|
trivroots(GEN nf) |
|
{ |
|
GEN y = cgetg(3, t_VEC); |
|
y[1] = deux; |
|
y[2] = (long)gscalcol_i(negi(gun), degpol(nf[1])); |
|
return y; |
|
} |
|
|
GEN |
GEN |
rootsof1(GEN nf) |
rootsof1(GEN nf) |
{ |
{ |
ulong av; |
gpmem_t av = avma; |
long N,k,i,ws,prec; |
long N,k,i,ws,prec; |
GEN algun,p1,y,R1,d,list,w; |
GEN p1,y,d,list,w; |
|
|
y=cgetg(3,t_VEC); av=avma; nf=checknf(nf); |
nf = checknf(nf); |
R1=gmael(nf,2,1); algun=gmael(nf,8,1); |
if ( nf_get_r1(nf) ) return trivroots(nf); |
if (signe(R1)) |
|
{ |
N = degpol(nf[1]); prec = nfgetprec(nf); |
y[1]=deux; |
|
y[2]=lneg(algun); return y; |
|
} |
|
N=degpol(nf[1]); prec=gprecision((GEN)nf[6]); |
|
#ifdef LONG_IS_32BIT |
|
if (prec < 10) prec = 10; |
|
#else |
|
if (prec < 6) prec = 6; |
|
#endif |
|
for (i=1; ; i++) |
for (i=1; ; i++) |
{ |
{ |
p1 = fincke_pohst(nf,stoi(N),1000,1,prec,NULL); |
GEN R = R_from_QR(gmael(nf,5,2), prec); |
if (p1) break; |
if (R) |
|
{ |
|
p1 = fincke_pohst(_vec(R),stoi(N),1000,1, 0, NULL); |
|
if (p1) break; |
|
} |
if (i == MAXITERPOL) err(accurer,"rootsof1"); |
if (i == MAXITERPOL) err(accurer,"rootsof1"); |
prec=(prec<<1)-2; |
prec = (prec<<1)-2; |
if (DEBUGLEVEL) err(warnprec,"rootsof1",prec); |
if (DEBUGLEVEL) err(warnprec,"rootsof1",prec); |
nf=nfnewprec(nf,prec); |
nf = nfnewprec(nf,prec); |
} |
} |
if (itos(ground((GEN)p1[2])) != N) err(bugparier,"rootsof1 (bug1)"); |
if (itos(ground((GEN)p1[2])) != N) err(bugparier,"rootsof1 (bug1)"); |
w=(GEN)p1[1]; ws = itos(w); |
w = (GEN)p1[1]; ws = itos(w); |
if (ws == 2) |
if (ws == 2) { avma = av; return trivroots(nf); } |
{ |
|
y[1]=deux; avma=av; |
|
y[2]=lneg(algun); return y; |
|
} |
|
|
|
d = decomp(w); list = (GEN)p1[3]; k = lg(list); |
d = decomp(w); list = (GEN)p1[3]; k = lg(list); |
for (i=1; i<k; i++) |
for (i=1; i<k; i++) |
{ |
{ |
p1 = (GEN)list[i]; |
p1 = (GEN)list[i]; |
p1 = is_primitive_root(nf,d,p1,ws); |
p1 = is_primitive_root(nf,d,p1,ws); |
if (p1) |
if (p1) break; |
{ |
|
y[2]=lpilecopy(av,p1); |
|
y[1]=lstoi(ws); return y; |
|
} |
|
} |
} |
err(bugparier,"rootsof1"); |
if (!p1) err(bugparier,"rootsof1"); |
return NULL; /* not reached */ |
y = cgetg(3, t_VEC); |
|
y[1] = lstoi(ws); |
|
y[2] = (long)p1; return gerepilecopy(av, y); |
} |
} |
|
|
/*******************************************************************/ |
/*******************************************************************/ |
Line 1419 rootsof1(GEN nf) |
|
Line 2075 rootsof1(GEN nf) |
|
/* DEDEKIND ZETA FUNCTION */ |
/* DEDEKIND ZETA FUNCTION */ |
/* */ |
/* */ |
/*******************************************************************/ |
/*******************************************************************/ |
|
|
ulong smulss(ulong x, ulong y, ulong *rem); |
|
|
|
static GEN |
static GEN |
dirzetak0(GEN nf, long N0) |
dirzetak0(GEN nf, long N0) |
{ |
{ |
GEN vect,p1,pol,disc,c,c2; |
GEN vect,p1,pol,disc,c,c2; |
long av=avma,i,j,k,limk,lx; |
gpmem_t av=avma; |
|
long i,j,k,limk,lx; |
ulong q,p,rem; |
ulong q,p,rem; |
byteptr d=diffptr; |
byteptr d=diffptr; |
long court[] = {evaltyp(t_INT)|m_evallg(3), evalsigne(1)|evallgefint(3),0}; |
long court[] = {evaltyp(t_INT)|_evallg(3), evalsigne(1)|evallgefint(3),0}; |
|
|
pol=(GEN)nf[1]; disc=(GEN)nf[4]; |
pol=(GEN)nf[1]; disc=(GEN)nf[4]; |
c = (GEN) gpmalloc((N0+1)*sizeof(long)); |
c = (GEN) gpmalloc((N0+1)*sizeof(long)); |
Line 1440 dirzetak0(GEN nf, long N0) |
|
Line 2094 dirzetak0(GEN nf, long N0) |
|
|
|
while (court[2]<=N0) |
while (court[2]<=N0) |
{ |
{ |
court[2] += *d++; if (! *d) err(primer1); |
NEXT_PRIME_VIADIFF_CHECK(court[2], d); |
if (smodis(disc,court[2])) /* court does not divide index */ |
if (smodis(disc,court[2])) /* court does not divide index */ |
{ vect = (GEN) simplefactmod(pol,court)[1]; lx=lg(vect); } |
{ vect = (GEN) simplefactmod(pol,court)[1]; lx=lg(vect); } |
else |
else |
Line 1494 initzeta(GEN pol, long prec) |
|
Line 2148 initzeta(GEN pol, long prec) |
|
GEN nfz,nf,alpha,beta,mu,gr1,gr2,gru,p1,p2,cst,A0,c0,c1,c2,eps,coef; |
GEN nfz,nf,alpha,beta,mu,gr1,gr2,gru,p1,p2,cst,A0,c0,c1,c2,eps,coef; |
GEN limx,bnf,resi,zet,C,coeflog,racpi,aij,tabj,colzero, *tabcstn, *tabcstni; |
GEN limx,bnf,resi,zet,C,coeflog,racpi,aij,tabj,colzero, *tabcstn, *tabcstni; |
GEN c_even,ck_even,c_odd,ck_odd,serie_even,serie_odd,serie_exp,Pi; |
GEN c_even,ck_even,c_odd,ck_odd,serie_even,serie_odd,serie_exp,Pi; |
long N0,imin,imax,r1,r2,ru,R,N,i,j,k,n, av,av2,tetpil; |
long N0,imin,imax,r1,r2,ru,R,N,i,j,k,n; |
long court[] = {evaltyp(t_INT)|m_evallg(3), evalsigne(1)|evallgefint(3),0}; |
gpmem_t av,av2,tetpil; |
|
long court[] = {evaltyp(t_INT)|_evallg(3), evalsigne(1)|evallgefint(3),0}; |
stackzone *zone, *zone0, *zone1; |
stackzone *zone, *zone0, *zone1; |
|
|
/*************** Calcul du residu et des constantes ***************/ |
/*************** Calcul du residu et des constantes ***************/ |
eps=gmul2n(gun,-bit_accuracy(prec)-6); p1=dbltor(0.5); |
eps=gmul2n(gun,-bit_accuracy(prec)-6); p1=dbltor(0.5); |
nfz=cgetg(10,t_VEC); |
nfz=cgetg(10,t_VEC); |
bnf=buchinit(pol,p1,p1,prec+1); prec=(prec<<1)-1; |
bnf = bnfinit0(pol,2,NULL,prec+1); prec=(prec<<1)-1; |
bnf = checkbnf_discard(bnf); |
bnf = checkbnf_discard(bnf); |
Pi = mppi(prec); racpi=gsqrt(Pi,prec); |
Pi = mppi(prec); racpi=gsqrt(Pi,prec); |
|
|
Line 1540 initzeta(GEN pol, long prec) |
|
Line 2195 initzeta(GEN pol, long prec) |
|
if (is_bigint(p1) || N0 > 10000000) |
if (is_bigint(p1) || N0 > 10000000) |
err(talker,"discriminant too large for initzeta, sorry"); |
err(talker,"discriminant too large for initzeta, sorry"); |
if (DEBUGLEVEL>=2) |
if (DEBUGLEVEL>=2) |
{ fprintferr("\ninitzeta:\nN0 = %ld\n",N0); flusherr(); timer2(); } |
{ fprintferr("\ninitzeta:\nN0 = %ld\n",N0); flusherr(); (void)timer2(); } |
|
|
/* Calcul de imax */ |
/* Calcul de imax */ |
|
|
Line 1563 initzeta(GEN pol, long prec) |
|
Line 2218 initzeta(GEN pol, long prec) |
|
zone = switch_stack(NULL,i + 2*(N0+1) + 6*prec); |
zone = switch_stack(NULL,i + 2*(N0+1) + 6*prec); |
zone1 = switch_stack(NULL,2*i); |
zone1 = switch_stack(NULL,2*i); |
zone0 = switch_stack(NULL,2*i); |
zone0 = switch_stack(NULL,2*i); |
switch_stack(zone,1); |
(void)switch_stack(zone,1); |
tabcstn = (GEN*) cgetg(N0+1,t_VEC); |
tabcstn = (GEN*) cgetg(N0+1,t_VEC); |
tabcstni = (GEN*) cgetg(N0+1,t_VEC); |
tabcstni = (GEN*) cgetg(N0+1,t_VEC); |
p1 = ginv(cst); |
p1 = ginv(cst); |
for (i=1; i<=N0; i++) { tabcstni[i] = gun; tabcstn[i] = mulsr(i,p1); } |
for (i=1; i<=N0; i++) { tabcstni[i] = gun; tabcstn[i] = mulsr(i,p1); } |
switch_stack(zone,0); |
(void)switch_stack(zone,0); |
|
|
/********** Calcul des coefficients a(i,j) independants de s **********/ |
/********** Calcul des coefficients a(i,j) independants de s **********/ |
|
|
Line 1595 initzeta(GEN pol, long prec) |
|
Line 2250 initzeta(GEN pol, long prec) |
|
if (k&1) ck_odd[k]=lneg((GEN)ck_odd[k]); |
if (k&1) ck_odd[k]=lneg((GEN)ck_odd[k]); |
ck_odd[k]=ladd(gru,(GEN)ck_odd[k]); |
ck_odd[k]=ladd(gru,(GEN)ck_odd[k]); |
} |
} |
ck_odd[1]=lsub((GEN)ck_odd[1],gmul(gr1,glog(gdeux,prec))); |
ck_odd[1]=lsub((GEN)ck_odd[1],gmul(gr1,mplog2(prec))); |
serie_even =cgetg(ru+3,t_SER); serie_odd=cgetg(r2+3,t_SER); |
serie_even =cgetg(ru+3,t_SER); serie_odd=cgetg(r2+3,t_SER); |
serie_even[1] = serie_odd[1] = evalsigne(1)+evalvalp(1); |
serie_even[1] = serie_odd[1] = evalsigne(1)+evalvalp(1); |
i=0; |
i=0; |
Line 1656 initzeta(GEN pol, long prec) |
|
Line 2311 initzeta(GEN pol, long prec) |
|
t[1] = lstoi(coef[i]); t[2] = lmul((GEN)t[1],p1); |
t[1] = lstoi(coef[i]); t[2] = lmul((GEN)t[1],p1); |
for (j=2; j<=ru; j++) |
for (j=2; j<=ru; j++) |
{ |
{ |
long av2 = avma; p2 = gmul((GEN)t[j],p1); |
gpmem_t av2 = avma; |
|
p2 = gmul((GEN)t[j], p1); |
t[j+1] = (long)gerepileuptoleaf(av2, divrs(p2,j)); |
t[j+1] = (long)gerepileuptoleaf(av2, divrs(p2,j)); |
} |
} |
/* tabj[n,j]=coef(n)*ln(c/n)^(j-1)/(j-1)! */ |
/* tabj[n,j]=coef(n)*ln(c/n)^(j-1)/(j-1)! */ |
Line 1705 initzeta(GEN pol, long prec) |
|
Line 2361 initzeta(GEN pol, long prec) |
|
} |
} |
/* use a parallel stack */ |
/* use a parallel stack */ |
z = i&1? zone1: zone0; |
z = i&1? zone1: zone0; |
switch_stack(z, 1); |
(void)switch_stack(z, 1); |
for (n=1; n<=N0; n++) |
for (n=1; n<=N0; n++) |
if (coef[n]) tabcstni[n] = mpmul(tabcstni[n],tabcstn[n]); |
if (coef[n]) tabcstni[n] = mpmul(tabcstni[n],tabcstn[n]); |
/* come back */ |
/* come back */ |
switch_stack(z, 0); |
(void)switch_stack(z, 0); |
} |
} |
nfz[4] = (long) C; |
nfz[4] = (long) C; |
if (DEBUGLEVEL>=2) msgtimer("Cik"); |
if (DEBUGLEVEL>=2) msgtimer("Cik"); |
Line 1724 gzetakall(GEN nfz, GEN s, long flag, long prec2) |
|
Line 2380 gzetakall(GEN nfz, GEN s, long flag, long prec2) |
|
GEN resi,C,cst,cstlog,coeflog,cs,coef; |
GEN resi,C,cst,cstlog,coeflog,cs,coef; |
GEN lambd,gammas,gammaunmoins,gammas2,gammaunmoins2,var1,var2; |
GEN lambd,gammas,gammaunmoins,gammas2,gammaunmoins2,var1,var2; |
GEN p1,unmoins,gexpro,gar,val,valm,valk,valkm; |
GEN p1,unmoins,gexpro,gar,val,valm,valk,valkm; |
long ts = typ(s), r1,r2,ru,imax,i,j,k,N0,sl,prec,bigprec, av = avma; |
long ts = typ(s), r1,r2,ru,imax,i,j,k,N0,sl,prec,bigprec; |
|
gpmem_t av = avma; |
|
|
if (typ(nfz)!=t_VEC || lg(nfz)!=10 || typ(nfz[1]) != t_VEC) |
if (typ(nfz)!=t_VEC || lg(nfz)!=10 || typ(nfz[1]) != t_VEC) |
err(talker,"not a zeta number field in zetakall"); |
err(talker,"not a zeta number field in zetakall"); |
Line 1824 gzetakall(GEN nfz, GEN s, long flag, long prec2) |
|
Line 2481 gzetakall(GEN nfz, GEN s, long flag, long prec2) |
|
for (k=1; k<=ru; k++) |
for (k=1; k<=ru; k++) |
{ |
{ |
p1 = gcoeff(C,i,k); |
p1 = gcoeff(C,i,k); |
lambd = gsub(lambd,gdiv(p1,valk )); valk = mulii(val,valk); |
lambd = gsub(lambd,gdiv(p1,valk )); valk = mulii(val,valk); |
lambd = gsub(lambd,gdiv(p1,valkm)); valkm = mulii(valm,valkm); |
lambd = gsub(lambd,gdiv(p1,valkm)); valkm = mulii(valm,valkm); |
} |
} |
val = addis(val,1); |
val = addis(val,1); |