File: [local] / OpenXM_contrib / pari / src / modules / Attic / subfield.c (download)
Revision 1.1.1.1 (vendor branch), Sun Jan 9 17:35:33 2000 UTC (24 years, 5 months ago) by maekawa
Branch: PARI_GP
CVS Tags: maekawa-ipv6, VERSION_2_0_17_BETA, RELEASE_20000124, RELEASE_1_2_3, RELEASE_1_2_2_KNOPPIX_b, RELEASE_1_2_2_KNOPPIX, RELEASE_1_2_2, RELEASE_1_2_1, RELEASE_1_1_3, RELEASE_1_1_2
Changes since 1.1: +0 -0
lines
Import PARI/GP 2.0.17 beta.
|
/*******************************************************************/
/* */
/* SUBFIELDS OF A NUMBER FIELD */
/* */
/* J. Klueners and M. Pohst, J. Symb. Comp. (1996), vol. 11 */
/* */
/*******************************************************************/
/* $Id: subfield.c,v 1.1.1.1 1999/09/16 13:48:20 karim Exp $ */
#include "pari.h"
#ifdef __WIN32
# include <io.h> /* for open, read, close */
#endif
GEN roots_to_pol(GEN a, long v);
static long TR; /* nombre de changements de polynomes (degre fixe) */
static GEN FACTORDL; /* factorisation of |disc(L)| */
static GEN print_block_system(long N,GEN Y,long d, GEN vbs);
static GEN compute_data(GEN nf,GEN ff,GEN p,long m,long nn);
/* Computation of potential block systems of given size d associated to a
* rational prime p: give a row vector of row vectors containing the
* potential block systems of imprimitivity; a potential block system is a
* vector of row vectors (enumeration of the roots).
*/
#define BIL 32 /* for 64bit machines also */
static GEN
calc_block(long N,GEN Z,long d,GEN Y,GEN vbs)
{
long r,lK,i,j,k,t,tp,T,lpn,u,nn,lnon,lY;
GEN K,n,non,pn,pnon,e,Yp,Zp,Zpp;
if (DEBUGLEVEL>3)
{
long l = vbs? lg(vbs): 0;
fprintferr("avma = %ld, lg(Z) = %ld, lg(Y) = %ld, lg(vbs) = %ld\n",
avma,lg(Z),lg(Y),l);
if (DEBUGLEVEL > 5)
{
fprintferr("Z = %Z\n",Z);
fprintferr("Y = %Z\n",Y);
if (DEBUGLEVEL > 7) fprintferr("vbs = %Z\n",vbs);
}
}
r=lg(Z); lnon = min(BIL, r);
e = new_chunk(BIL);
n = new_chunk(r);
non = new_chunk(lnon);
pnon = new_chunk(lnon);
pn = new_chunk(lnon);
Zp = cgetg(lnon,t_VEC);
Zpp = cgetg(lnon,t_VEC);
for (i=1; i<r; i++) n[i] = lg(Z[i])-1;
K=divisors(stoi(n[1])); lK=lg(K);
for (i=1; i<lK; i++)
{
lpn=0; k = itos((GEN)K[i]);
for (j=2; j<r; j++)
if (n[j]%k == 0)
{
if (++lpn >= BIL) err(talker,"overflow in calc_block");
pn[lpn]=n[j]; pnon[lpn]=j;
}
if (!lpn)
{
if (d*k != n[1]) continue;
T=1; lpn=1;
}
else
T = 1<<lpn;
for (t=0; t<T; t++)
{
for (nn=n[1],tp=t, u=1; u<=lpn; u++,tp>>=1)
{
if (tp&1) { nn += pn[u]; e[u]=1; } else e[u]=0;
}
if (d*k == nn)
{
long av=avma;
int Z_equal_Zp = 1;
for (j=1; j<lnon; j++) non[j]=0;
Zp[1]=Z[1];
for (u=2,j=1; j<=lpn; j++)
if (e[j])
{
Zp[u]=Z[pnon[j]]; non[pnon[j]]=1;
if (Zp[u] != Z[u]) Z_equal_Zp = 0;
u++;
}
setlg(Zp, u);
lY=lg(Y); Yp = cgetg(lY+1,t_VEC);
for (j=1; j<lY; j++) Yp[j]=Y[j];
Yp[lY]=(long)Zp;
if (r == u && Z_equal_Zp)
vbs = print_block_system(N,Yp,d,vbs);
else
{
for (u=1,j=2; j<r; j++)
if (!non[j]) Zpp[u++] = Z[j];
setlg(Zpp, u);
vbs = calc_block(N,Zpp,d,Yp,vbs);
}
avma=av;
}
}
}
return vbs;
}
static GEN
potential_block_systems(long N, long d,GEN ff,long *n)
{
long av=avma,r,i,j,k;
GEN p1,vbs,Z;
r=lg(ff); Z=cgetg(r,t_VEC);
for (k=0,i=1; i<r; i++)
{
p1=cgetg(n[i]+1,t_VECSMALL); Z[i]=(long)p1;
for (j=1; j<=n[i]; j++) p1[j]= ++k;
}
vbs=calc_block(N,Z,d, cgetg(1,t_VEC), NULL);
avma=av; return vbs;
}
/* product of permutations. Put the result in perm1. */
static void
perm_mul(GEN perm1,GEN perm2)
{
long av = avma,i, N = lg(perm1);
GEN perm=new_chunk(N);
for (i=1; i<N; i++) perm[i]=perm1[perm2[i]];
for (i=1; i<N; i++) perm1[i]=perm[i];
avma=av;
}
/* cy is a cycle; compute cy^l as a permutation */
static GEN
cycle_power_to_perm(GEN perm,GEN cy,long l)
{
long lp,i,j,b, N = lg(perm), lcy = lg(cy)-1;
lp = l % lcy;
for (i=1; i<N; i++) perm[i] = i;
if (lp)
{
long av = avma;
GEN p1 = new_chunk(N);
b = cy[1];
for (i=1; i<lcy; i++) b = (perm[b] = cy[i+1]);
perm[b] = cy[1];
for (i=1; i<N; i++) p1[i] = perm[i];
for (j=2; j<=lp; j++) perm_mul(perm,p1);
avma = av;
}
return perm;
}
/* image du block system D par la permutation perm */
static GEN
im_block_by_perm(GEN D,GEN perm)
{
long i,j,lb,lcy;
GEN Dn,cy,p1;
lb=lg(D); Dn=cgetg(lb,t_VEC);
for (i=1; i<lb; i++)
{
cy=(GEN)D[i]; lcy=lg(cy);
Dn[i]=lgetg(lcy,t_VECSMALL); p1=(GEN)Dn[i];
for (j=1; j<lcy; j++) p1[j] = perm[cy[j]];
}
return Dn;
}
/* cy is a cycle; recturn cy(a) */
static long
im_by_cy(long a,GEN cy)
{
long k, l = lg(cy);
k=1; while (k<l && cy[k] != a) k++;
if (k == l) return a;
k++; if (k == l) k = 1;
return cy[k];
}
/* renvoie 0 si l'un des coefficients de g[i] est de module > M[i]; 1 sinon */
static long
ok_coeffs(GEN g,GEN M)
{
long i, lg = lgef(g)-1; /* g is monic, and cst term is ok */
for (i=3; i<lg; i++)
if (absi_cmp((GEN)g[i], (GEN)M[i]) > 0) return 0;
return 1;
}
/* vbs[0] = current cardinality+1, vbs[1] = max number of elts */
static GEN
append_vbs(GEN vbs, GEN D)
{
long l,maxl,i,j,n, lD = lg(D);
GEN Dn, last;
n = 0; for (i=1; i<lD; i++) n += lg(D[i]);
Dn = (GEN)gpmalloc((lD + n) * sizeof(long));
last = Dn + lD; Dn[0] = D[0];
for (i=1; i<lD; i++)
{
GEN cy = (GEN)D[i], cn = last;
for (j=0; j<lg(cy); j++) cn[j] = cy[j];
Dn[i] = (long)cn; last = cn + j;
}
if (!vbs)
{
maxl = 1024;
vbs = (GEN)gpmalloc((2 + maxl)*sizeof(GEN));
*vbs = maxl; vbs++; setlg(vbs, 1);
}
l = lg(vbs); maxl = vbs[-1];
if (l == maxl)
{
vbs = (GEN)gprealloc((void*)(vbs-1), (2 + (maxl<<1))*sizeof(GEN),
(2 + maxl)*sizeof(GEN));
*vbs = maxl<<1; vbs++; setlg(vbs, 1);
}
if (DEBUGLEVEL>5) fprintferr("appending D = %Z\n",D);
vbs[l] = (long)Dn; setlg(vbs, l+1); return vbs;
}
GEN
myconcat(GEN D, long a)
{
long i,l = lg(D);
GEN x = cgetg(l+1,t_VECSMALL);
for (i=1; i<l; i++) x[i]=D[i];
x[l] = a; return x;
}
void
myconcat2(GEN D, GEN a)
{
long i,l = lg(D), m = lg(a);
GEN x = D + (l-1);
for (i=1; i<m; i++) x[i]=a[i];
setlg(D, l+m-1);
}
static GEN
print_block_system(long N,GEN Y,long d, GEN vbs)
{
long a,i,j,l,ll,*k,*n,lp,**e,u,v,*t,ns, r = lg(Y);
GEN D,De,Z,cyperm,perm,p1,empty;
if (DEBUGLEVEL>5) fprintferr("Y = %Z\n",Y);
empty = cgetg(1,t_VEC);
n = new_chunk(N+1);
D = cgetg(N+1, t_VEC); setlg(D,1);
t=new_chunk(r+1); k=new_chunk(r+1); Z=cgetg(r+1,t_VEC);
for (ns=0,i=1; i<r; i++)
{
GEN Yi = (GEN)Y[i], cy;
long ki = 0, si = lg(Yi)-1;
for (j=1; j<=si; j++) { n[j]=lg(Yi[j])-1; ki += n[j]; }
ki /= d;
De=cgetg(ki+1,t_VEC);
for (j=1; j<=ki; j++) De[j]=(long)empty;
for (j=1; j<=si; j++)
{
a = mael(Yi,j,1); cy = (GEN)Yi[j];
for (l=1,lp=0; l<=n[j]; l++)
{
lp++; if (lp>ki) lp = 1;
a = im_by_cy(a, cy);
De[lp] = (long)myconcat((GEN)De[lp], a);
}
}
if (si>1 && ki>1)
{
ns++; t[ns]=si-1; k[ns]=ki;
Z[ns]=lgetg(si,t_VEC); p1=(GEN)Z[ns];
for (j=2; j<=si; j++) p1[j-1]=Yi[j];
}
myconcat2(D,De);
}
if (DEBUGLEVEL>2) { fprintferr("\nns = %ld\n",ns); flusherr(); }
if (!ns) return append_vbs(vbs,D);
setlg(Z, ns+1);
e=(long**)new_chunk(ns+1);
for (i=1; i<=ns; i++)
{
e[i]=new_chunk(t[i]+1);
for (j=1; j<=t[i]; j++) e[i][j]=0;
}
cyperm = cgetg(N+1,t_VEC);
perm = cgetg(N+1,t_VEC); i=ns;
do
{
long av = avma;
if (DEBUGLEVEL>5)
{
for (l=1; l<=ns; l++)
{
for (ll=1; ll<=t[l]; ll++)
fprintferr("e[%ld][%ld] = %ld, ",l,ll,e[l][ll]);
fprintferr("\n");
}
fprintferr("\n"); flusherr();
}
for (u=1; u<=N; u++) perm[u]=u;
for (u=1; u<=ns; u++)
for (v=1; v<=t[u]; v++)
perm_mul(perm, cycle_power_to_perm(cyperm,gmael(Z,u,v),e[u][v]));
vbs = append_vbs(vbs, im_block_by_perm(D,perm));
avma = av;
e[ns][t[ns]]++;
if (e[ns][t[ns]] >= k[ns])
{
j=t[ns]-1;
while (j>=1 && e[ns][j] == k[ns]-1) j--;
if (j>=1) { e[ns][j]++; for (l=j+1; l<=t[ns]; l++) e[ns][l]=0; }
else
{
i=ns-1;
while (i>=1)
{
j=t[i];
while (j>=1 && e[i][j] == k[i]-1) j--;
if (j<1) i--;
else
{
e[i][j]++;
for (l=j+1; l<=t[i]; l++) e[i][l]=0;
for (ll=i+1; ll<=ns; ll++)
for (l=1; l<=t[ll]; l++) e[ll][l]=0;
break;
}
}
}
}
}
while (i>0);
return vbs;
}
/* rend le numero du cycle (a1,...,an) dans le support duquel se trouve a */
/* met dans *pt l'indice i tq ai = a */
static long
in_what_cycle(long a,GEN cys,long *pt)
{
long i,k,nk, lcys=lg(cys);
for (k=1; k<lcys; k++)
{
GEN c = (GEN)cys[k]; nk = lg(c);
for (i=1; i<nk; i++)
if (a == c[i]) { *pt = i; return k; }
}
err(talker,"impossible to find %d in in_what_cycle",a);
return 0; /* not reached */
}
/* Return common factors to A and B + s the prime to A part of B */
static GEN
commonfactor(GEN A, GEN B)
{
GEN f,p1,p2,s, y = cgetg(3,t_MAT);
long lf,i;
s = absi(B); f=(GEN)A[1]; lf=lg(f);
p1=cgetg(lf+1,t_COL); y[1]=(long)p1;
p2=cgetg(lf+1,t_COL); y[2]=(long)p2;
for (i=1; i<lf; i++)
{
p1[i] = f[i];
p2[i] = lstoi(pvaluation(s,(GEN)f[i], &s));
}
p1[i] = (long)s;
p2[i] = un; return y;
}
static void
polsimplify(GEN x)
{
long i,lx = lgef(x);
for (i=2; i<lx; i++)
if (typ(x[i]) == t_POL) x[i] = signe(x[i])? mael(x,i,2): zero;
}
/* Renvoie un polynome g definissant un sous-corps potentiel, ou
* 0: si le polynome trouve n'est pas separable,
* 1: si les coefficients du polynome trouve sont plus grands que la borne M,
* 2: si p divise le discriminant de g,
* 3: si le discriminant de g est nul,
* 4: si la partie s de d(g) premiere avec d(L) n'est pas un carre,
* 5: si s est un carre et si un des facteurs premiers communs a d(g) et d(L)
* a un exposant impair dans d(g) et un exposant plus petit que d dans d(L),
* 6: si le discriminant du corps defini par g a la puissance d ne divise pas
* le discriminant du corps nf (soit L).
*/
static GEN
cand_for_subfields(GEN A,GEN DATA,GEN *ptdelta,GEN *ptrootsA)
{
long av=avma,N,m,i,j,d,lf;
GEN P,pe,p,pol,cys,tabroots,delta,g,dg,unmodpe,tabrA;
GEN factcommon,ff1,ff2,p1;
GEN *gptr[3];
pol=(GEN)DATA[1]; N=lgef(pol)-3; m=lg(A)-1; d=N/m;
if (N%m) err(talker,"incompatible block system in cand_for_subfields");
p = (GEN)DATA[2];
cys=(GEN)DATA[5];
tabroots=(GEN)DATA[10];
pe = gclone((GEN)DATA[9]);
unmodpe = cgetg(3,t_INTMOD); unmodpe[1]=(long)pe; unmodpe[2]=un;
delta = cgetg(m+1,t_VEC);
tabrA = cgetg(m+1,t_VEC);
for (i=1; i<=m; i++)
{
GEN Ai=(GEN)A[i], col = cgetg(d+1,t_VEC);
long l,k;
tabrA[i]=(long)col; p1 = unmodpe;
for (j=1; j<=d; j++)
{
l=in_what_cycle(Ai[j],cys,&k);
col[j] = mael(tabroots, l, k);
p1 = gmul(p1, (GEN)col[j]);
}
p1 = lift_intern((GEN)p1[2]);
for (j=1; j<i; j++)
if (gegal(p1,(GEN)delta[j])) { avma=av; return gzero; }
if (DEBUGLEVEL>2) fprintferr("delta[%ld] = %Z\n",i,p1);
delta[i] = (long)p1;
}
P = gmael3(tabroots,1,1,1);
for (i=1; i<=m; i++)
{
p1 = cgetg(3,t_POLMOD); p1[1]=(long)P; p1[2]=delta[i];
delta[i] = (long)p1;
}
g = roots_to_pol(gmul(unmodpe,delta),0);
g=centerlift(lift_intern(g)); polsimplify(g);
if (DEBUGLEVEL>2) fprintferr("pol. found = %Z\n",g);
if (!ok_coeffs(g,(GEN)DATA[8])) return gun;
dg=discsr(g);
if (!signe(dg)) return stoi(3);
if (!signe(resii(dg,p))) return gdeux;
factcommon=commonfactor(FACTORDL,dg);
ff1=(GEN)factcommon[1]; lf=lg(ff1)-1;
if (!carreparfait((GEN)ff1[lf])) return stoi(4);
ff2=(GEN)factcommon[2];
for (i=1; i<lf; i++)
if (mod2((GEN)ff2[i]) && itos(gmael(FACTORDL,2,i)) < d) return stoi(5);
gunclone(pe);
*ptdelta=delta; *ptrootsA=tabrA;
gptr[0]=&g; gptr[1]=ptdelta; gptr[2]=ptrootsA;
gerepilemany(av,gptr,3); return g;
}
/* a partir d'un polynome h(x) dont les coefficients sont definis mod p^k,
* on construit un polynome a coefficients dans Q dont les coefficients ont
* pour approximation p-adique les coefficients de h */
static GEN
retrieve_p_adique_polynomial_in_Q(GEN ind,GEN h)
{
return gdiv(centerlift(gmul(h,ind)), ind);
}
/* interpolation polynomial P(x) s.t P(T[j][i]) = delta[i] mod p */
static GEN
interpolation_polynomial(GEN T, GEN delta)
{
long i,j,i1,j1, m = lg(T), d = lg(T[1]);
GEN P = NULL, x0 = gneg(polx[0]);
for (j=1; j<m; j++)
{
GEN p3 = NULL;
for (i=1; i<d; i++)
{
GEN p1=gun, p2=gun, a = gneg(gmael(T,j,i));
for (j1=1; j1<m; j1++)
for (i1=1; i1<d; i1++)
if (i1 != i || j1 != j)
{
p1 = gmul(p1,gadd(gmael(T,j1,i1), x0));
p2 = gmul(p2,gadd(gmael(T,j1,i1), a));
}
p1 = gdiv(p1,p2);
p3 = p3? gadd(p3, p1): p1;
}
p3 = gmul((GEN)delta[j],p3);
P = P? gadd(P,p3): p3;
}
return P;
}
/* nf est le corps de nombres, g un polynome de Z[x] candidat
* pour definir un sous-corps, p le nombre premier ayant servi a definir le
* potential block system rootsA donne par les racines avec une approximation
* convenable, e est la precision p-adique des elements de rootsA et delta la
* liste des racines de g dans une extension convenable en precision p^e.
* Renvoie un polynome h de Q[x] tel que f divise g o h et donc tel que le
* couple (g,h) definisse un sous-corps, ou bien gzero si rootsA n'est pas un
* block system
*/
static GEN
embedding_of_potential_subfields(GEN nf,GEN g,GEN DATA,GEN rootsA,GEN delta)
{
GEN w0_inQ,w0,w1,h0,gp,p2,f,unmodp,p,ind, maxp;
long av = avma, av1;
f=(GEN)nf[1]; ind=(GEN)nf[4]; p=(GEN)DATA[2];
maxp=mulii((GEN)DATA[11],ind);
gp=deriv(g,varn(g)); unmodp=gmodulsg(1,p);
av1 = avma;
w0 = interpolation_polynomial(gmul(rootsA,unmodp), delta);
w0 = lift_intern(w0); /* in Fp[x] */
polsimplify(w0);
w0_inQ = retrieve_p_adique_polynomial_in_Q(ind,w0);
(void)gbezout(poleval(gp,w0), gmul(unmodp,f), &h0, &p2);
w0 = lift_intern(w0); /* in Z[x] */
h0 = lift_intern(lift_intern(h0));
for(;;)
{
GEN p1;
/* Given g in Z[x], gp its derivative, p a prime, [w0,h0] in Z[x] s.t.
* h0(x).gp(w0(x)) = 1 and g(w0(x)) = 0 (mod f,mod p), return
* [w1,h1] satisfying the same condition mod p^2. Moreover,
* [w1,h1] = [w0,h0] (mod p)
* (cf. Dixon: J. Austral. Math. Soc., Series A, vol.49, 1990, p.445) */
if (DEBUGLEVEL>2)
{
fprintferr("w = "); outerr(w0);
fprintferr("h = "); outerr(h0);
}
p = sqri(p); unmodp[1] = (long)p;
p1 = gneg(gmul(h0, poleval(g,w0)));
w1 = gres(gmul(unmodp,gadd(w0,p1)), f);
p2 = retrieve_p_adique_polynomial_in_Q(ind,w1);
if ((gegal(p2, w0_inQ) || cmpii(p,maxp)) && gdivise(poleval(g,p2), f))
return gerepileupto(av, poleval(p2, gadd(polx[0],stoi(TR))));
if (DEBUGLEVEL>2)
{
fprintferr("Old Q-polynomial: "); outerr(w0_inQ);
fprintferr("New Q-polynomial: "); outerr(p2);
}
if (cmpii(p, maxp) > 0)
{
if (DEBUGLEVEL) fprintferr("coeff too big for embedding\n");
avma=av; return gzero;
}
w1 = lift_intern(w1);
p1 = gneg(gmul(h0, poleval(gp,w1)));
p1 = gmul(h0, gadd(gdeux,p1));
h0 = lift_intern(gres(gmul(unmodp,p1), f));
w0 = w1; w0_inQ = p2;
{
GEN *gptr[4]; gptr[0]=&w0; gptr[1]=&h0; gptr[2]=&w0_inQ; gptr[3]=&p;
gerepilemany(av1,gptr,4);
}
}
}
static long
choose_prime(GEN pol,GEN dpol,long d,GEN *ptff,GEN *ptlistpotbl, long *ptnn)
{
long j,k,oldllist,llist,r,nn,oldnn,*n,N,pp;
GEN p,listpotbl,oldlistpotbl,ff,oldff,p3;
byteptr di=diffptr;
if (DEBUGLEVEL) timer2();
di++; p = stoi(2); N = lgef(pol)-3;
while (p[2]<=N) p[2] += *di++;
oldllist = oldnn = BIGINT;
n = new_chunk(N+1);
for(k=1; k<11 || oldnn == BIGINT; k++,p[2]+= *di++)
{
long av=avma;
while (!smodis(dpol,p[2])) p[2] += *di++;
ff=(GEN)factmod(pol,p)[1]; r=lg(ff)-1;
if (r>1 && r<N)
{
for (j=1; j<=r; j++) n[j]=lgef(ff[j])-3;
p3 = stoi(n[1]);
for (j=2; j<=r; j++) p3 = glcm(p3,stoi(n[j]));
nn=itos(p3);
if (nn > oldnn)
{
if (DEBUGLEVEL)
{
fprintferr("p = %ld,\tr = %ld,\tnn = %ld,\t#pbs = skipped\n",
p[2],r,nn);
}
continue;
}
listpotbl=potential_block_systems(N,d,ff,n);
if (!listpotbl) { oldlistpotbl = NULL; pp = p[2]; break; }
llist=lg(listpotbl)-1;
if (DEBUGLEVEL)
{
fprintferr("Time: %ldms,\tp = %ld,\tr = %ld,\tnn = %ld,\t#pbs = %ld\n",
timer2(),p[2],r,nn,llist);
flusherr();
}
if (nn<oldnn || llist<oldllist)
{
oldllist=llist; oldlistpotbl=listpotbl;
pp=p[2]; oldff=ff; oldnn=nn; continue;
}
for (j=1; j<llist; j++) free((void*)listpotbl[j]);
free((void*)(listpotbl-1));
}
avma = av;
}
if (DEBUGLEVEL)
{
fprintferr("Chosen prime: p = %ld\n",pp);
if (DEBUGLEVEL>2)
fprintferr("List of potential block systems of size %ld: %Z\n",
d,oldlistpotbl);
flusherr();
}
*ptlistpotbl=oldlistpotbl; *ptff=oldff; *ptnn=oldnn; return pp;
}
static GEN
change_pol(GEN nf, GEN ff)
{
long i,l;
GEN pol = (GEN)nf[1], p1 = gsub(polx[0],gun);
TR++; pol=poleval(pol, p1);
nf = dummycopy(nf);
nf[6] = (long)dummycopy((GEN)nf[6]);
l=lg(ff);
for (i=1; i<l; i++) ff[i]=(long)poleval((GEN)ff[i], p1);
l=lg(nf[6]); p1=(GEN)nf[6];
for (i=1; i<l; i++) p1[i]=ladd(gun,(GEN)p1[i]);
nf[1]=(long)pol; return nf;
}
static GEN
bound_for_coeff(long m,GEN rr,long r1, GEN *maxroot)
{
long i, lrr=lg(rr);
GEN p1,b1,b2,B,M, C = matpascal(m-1);
rr = gabs(rr,DEFAULTPREC); *maxroot = vecmax(rr);
for (i=1; i<lrr; i++)
if (gcmp((GEN)rr[i], gun) < 0) rr[i] = un;
for (b1=gun,i=1; i<=r1; i++) b1 = gmul(b1, (GEN)rr[i]);
for (b2=gun ; i<lrr; i++) b2 = gmul(b2, (GEN)rr[i]);
B = gmul(b1, gsqr(b2));
M = cgetg(m+2, t_VEC); M[1]=M[2]=zero; /* unused */
for (i=1; i<m; i++)
{
p1 = gadd(gmul(gcoeff(C, m, i), B),
gcoeff(C, m, i+1));
M[i+2] = lceil(p1);
}
return M;
}
/* liste des sous corps de degre d du corps de nombres nf */
static GEN
subfields_of_given_degree(GEN nf,GEN dpol,long d)
{
long av,av1,av2,tetpil,pp,llist,i,nn,N;
GEN listpotbl,ff,A,delta,rootsA,CSF,ESF,p1,p2,LSB;
GEN DATA, pol = (GEN)nf[1];
av=avma;
N = lgef(pol)-3;
pp=choose_prime(pol,dpol,N/d,&ff,&listpotbl,&nn);
if (!listpotbl) { avma=av; return cgetg(1,t_VEC); }
llist=lg(listpotbl);
LAB0:
av1=avma; LSB=cgetg(1,t_VEC);
DATA=compute_data(nf,ff,stoi(pp),d,nn);
for (i=1; i<llist; i++)
{
av2=avma; A=(GEN)listpotbl[i];
if (DEBUGLEVEL > 1)
fprintferr("\n* Potential block # %ld: %Z\n",i,A);
CSF=cand_for_subfields(A,DATA,&delta,&rootsA);
if (typ(CSF)==t_INT)
{
if (DEBUGLEVEL > 1) switch(itos(CSF))
{
case 0: fprintferr("changing f(x): non separable g(x)\n"); break;
case 1: fprintferr("coeff too big for pol g(x)\n"); break;
case 2: fprintferr("changing f(x): p divides disc(g(x))\n"); break;
case 3: fprintferr("non irreducible polynomial g(x)\n"); break;
case 4: fprintferr("prime to d(L) part of d(g) not a square\n"); break;
case 5: fprintferr("too small exponent of a prime factor in d(L)\n"); break;
case 6: fprintferr("the d-th power of d(K) does not divide d(L)\n");
}
switch(itos(CSF))
{
case 0: case 2:
avma=av1; nf = change_pol(nf,ff); pol = (GEN)nf[1];
if (DEBUGLEVEL) fprintferr("new f = %Z\n",pol);
goto LAB0;
}
avma=av2;
}
else
{
if (DEBUGLEVEL) fprintferr("candidate = %Z\n",CSF);
ESF=embedding_of_potential_subfields(nf,CSF,DATA,rootsA,delta);
if (ESF == gzero) avma=av2;
else
{
if (DEBUGLEVEL) fprintferr("embedding = %Z\n",ESF);
p1=cgetg(3,t_VEC); p2=cgetg(2,t_VEC); p2[1]=(long)p1;
p1[1]=(long)CSF;
p1[2]=(long)ESF; tetpil=avma;
LSB=gerepile(av2,tetpil, concat(LSB,p2));
}
}
}
for (i=1; i<llist; i++) free((void*)listpotbl[i]);
free((void*)(listpotbl-1)); tetpil=avma;
return gerepile(av,tetpil,gcopy(LSB));
}
GEN
subfields(GEN nf,GEN d)
{
long av=avma,di,N,v0,lp1,i;
GEN dpol,p1,LSB,p2,pol;
nf=checknf(nf); pol = (GEN)nf[1];
v0=varn(pol); N=lgef(pol)-3; di=itos(d);
if (di==N)
{
LSB=cgetg(2,t_VEC); p1=cgetg(3,t_VEC); LSB[1]=(long)p1;
p1[1]=lcopy(pol); p1[2]=lpolx[v0]; return LSB;
}
if (di==1)
{
LSB=cgetg(2,t_VEC); p1=cgetg(3,t_VEC); LSB[1]=(long)p1;
p1[1]=lpolx[v0]; p1[2]=lcopy(pol); return LSB;
}
if (di<=0 || di>N || N%di) return cgetg(1,t_VEC);
TR=0; dpol=mulii((GEN)nf[3],sqri((GEN)nf[4]));
if (v0) nf=gsubst(nf,v0,polx[0]);
FACTORDL=factor(absi((GEN)nf[3]));
p1=subfields_of_given_degree(nf,dpol,di); lp1=lg(p1)-1;
if (v0)
for (i=1; i<=lp1; i++)
{ p2=(GEN)p1[i]; setvarn(p2[1],v0); setvarn(p2[2],v0); }
return gerepileupto(av,p1);
}
static GEN
subfieldsall(GEN nf)
{
long av=avma,av1,N,ld,d,i,j,lNLSB,v0,lp1;
GEN pol,dpol,dg,LSB,NLSB,p1,p2;
nf=checknf(nf); pol = (GEN)nf[1];
v0=varn(pol); N=lgef(pol)-3;
if (isprime(stoi(N)))
{
avma=av; LSB=cgetg(3,t_VEC);
LSB[1]=lgetg(3,t_VEC); LSB[2]=lgetg(3,t_VEC);
p1=(GEN)LSB[1]; p1[1]=lcopy(pol); p1[2]=lpolx[v0];
p2=(GEN)LSB[2]; p2[1]=p1[2]; p2[2]=p1[1];
return LSB;
}
FACTORDL=factor(absi((GEN)nf[3])); dg=divisors(stoi(N));
dpol=mulii(sqri((GEN)nf[4]),(GEN)nf[3]);
if (DEBUGLEVEL>0)
{
fprintferr("\n***** Entering subfields\n\n");
fprintferr("pol = "); outerr(pol);
fprintferr("dpol = "); outerr(dpol);
fprintferr("divisors = "); outerr(dg);
}
ld=lg(dg)-1; LSB=cgetg(2,t_VEC); LSB[1]=lgetg(3,t_VEC);
p1=(GEN)LSB[1]; p1[1]=(long)pol; p1[2]=(long)polx[0];
if (v0) nf=gsubst(nf,v0,polx[0]);
for (i=2; i<ld; i++)
{
TR=0; av1=avma; d=itos((GEN)dg[i]);
if (DEBUGLEVEL>0)
{
fprintferr("\n*** Looking for subfields of degree %ld\n\n",N/d);
flusherr();
}
NLSB=subfields_of_given_degree(nf,dpol,N/d);
if (DEBUGLEVEL)
{
fprintferr("\nSubfields of degree %ld:\n",N/d);
lNLSB=lg(NLSB)-1; for (j=1; j<=lNLSB; j++) outerr((GEN)NLSB[j]);
}
if (lg(NLSB)>1) LSB = concatsp(LSB,NLSB); else avma=av1;
}
p1=cgetg(2,t_VEC); p1[1]=lgetg(3,t_VEC);p2=(GEN)p1[1];
p2[1]=(long)polx[0]; p2[2]=(long)pol;
LSB=concatsp(LSB,p1); lp1=lg(LSB)-1;
LSB = gerepileupto(av, gcopy(LSB));
if (v0)
for (i=1; i<=lp1; i++)
{ p2=(GEN)LSB[i]; setvarn(p2[1],v0); setvarn(p2[2],v0); }
if (DEBUGLEVEL>0) fprintferr("\n***** Leaving subfields\n\n");
return LSB;
}
GEN
subfields0(GEN nf,GEN d)
{
return d? subfields(nf,d): subfieldsall(nf);
}
/* irreducible (unitary) polynomial of degree n over Fp[v] */
GEN
ffinit(GEN p,long n,long v)
{
long av,av1,tetpil,i,*a,j,l,pp;
GEN pol,fpol;
if (n<=0) err(talker,"non positive degree in ffinit");
if (is_bigint(p)) err(talker,"prime field too big in ffinit");
if (v<0) v = 0;
av=avma; pp=itos(p); pol = cgetg(n+3,t_POL);
pol[1] = evalsigne(1)|evalvarn(v)|evallgef(n+3);
a=new_chunk(n+2);
a[1]=1; for (i=2; i<=n+1; i++) a[i]=0;
pol[n+2]=un; av1=avma;
for(;;)
{
a[n+1]++;
if (a[n+1]>=pp)
{
j=n; while (j>=2 && a[j]==pp-1) j--;
if (j>=2) { a[j]++; for (l=j+1; l<=n+1; l++) a[l]=0; }
}
for (i=2; i<=n+1; i++) pol[i]=lstoi(a[n+3-i]);
fpol=simplefactmod(pol,p);
if (lg(fpol[1])==2 && gcmp1(gmael(fpol,2,1))) break;
avma=av1;
}
tetpil=avma; return gerepile(av,tetpil,Fp_pol(pol,p));
}
static GEN
lift_coeff(GEN x, GEN fq)
{
GEN r;
if (typ(x) == t_POLMOD) { r = x; x = (GEN)x[2]; }
else r = cgetg(3,t_POLMOD);
r[1]=(long)fq; r[2]=(long)lift_intern(x); return r;
}
/* a is a polynomial whose coeffs are in Fq (= (Z/p)[y] / (fqbar), where
* fqbar is the reduction of fq mod p).
* Lift _in place_ the coeffs so that they belong to Z[y] / (fq)
*/
static GEN
special_lift(GEN a,GEN fq)
{
long la,i;
GEN c;
if (typ(a)==t_POL)
{
la=lgef(a); c=cgetg(la,t_POL); c[1]=a[1];
for (i=2; i<la; i++) c[i]=(long)lift_coeff((GEN)a[i],fq);
return c;
}
return lift_coeff(a,fq);
}
/* Hensel lift: fk = vector of factors of pol (unramified) in finite field
* Fp / fkk. Lift it to the precision p^e. This is equivalent to working
* in precision pi^e in the unramified extension of Qp given by fkk.
*/
GEN
hensel_lift(GEN pol,GEN fk,GEN fkk,GEN p,long e)
{
long av = avma, i, r = lg(fk)-1;
GEN p1,A,B,C,R,U,V,fklift,fklift2,fk2;
GEN unmodp = gmodulsg(1,p), fq = lift(fkk);
fk2=cgetg(r+1,t_VEC);
fklift=cgetg(r+1,t_VEC);
fklift2=cgetg(r+1,t_VEC);
fk2[r] = fklift2[r] = un;
for (i=r; i>1; i--)
{
fk2[i-1] = lmul((GEN)fk2[i],(GEN)fk[i]);
fklift[i] = (long)special_lift(gcopy((GEN)fk[i]),fq);
fklift2[i-1] = lmul((GEN)fklift2[i],(GEN)fklift[i]);
}
fklift[1] = (long)special_lift(gcopy((GEN)fk[1]),fq);
R=cgetg(r+1,t_VEC); C=pol;
for (i=1; i<r; i++)
{ /* treat factors two by two: fk[i] and fk2[i] = product fk[i+1..] */
long av1 = avma,tetpil1, ex = 1;
GEN pp;
(void)gbezout((GEN)fk[i],(GEN)fk2[i],&U,&V);
A = (GEN)fklift[i]; U = special_lift(U,fq);
B = (GEN)fklift2[i]; V = special_lift(V,fq);
for (pp=p;; pp=sqri(pp))
{ /* Algorithm 3.5.[5,6] H. Cohen page 137 (1995) */
GEN f,t,A0,B0,U0,V0;
unmodp[1] = (long)pp;
p1 = gneg_i(gmul(A,B));
p1=gdiv(gadd(C,p1),pp);
f=gmul(p1,unmodp);
t=poldivres(gmul(V,f),A, &A0);
A0=special_lift(A0,fq);
B0=special_lift(gadd(gmul(U,f),gmul(B,t)),fq);
A0 = gmul(A0,pp);
B0 = gmul(B0,pp); tetpil1 = avma;
A = gadd(A, A0);
B = gadd(B, B0); ex <<= 1;
if (ex>=e)
{
GEN *gptr[2]; gptr[0]=&A; gptr[1]=&B;
gerepilemanysp(av1,tetpil1,gptr,2);
C = B; R[i] = (long)A; break;
}
p1 = gneg_i(gadd(gmul(U,A),gmul(V,B)));
p1=gdiv(gadd(gun,p1),pp);
f=gmul(p1,unmodp);
t=poldivres(gmul(V,f),A, &V0);
U0=special_lift(gadd(gmul(U,f),gmul(B,t)),fq);
V0=special_lift(V0,fq);
U = gadd(U, gmul(U0,pp));
V = gadd(V, gmul(V0,pp));
}
}
if (r==1) C = gcopy(C);
R[r] = (long)C; return gerepileupto(av,R);
}
/* etant donne nf et p et la factorisation de nf[1] mod p, et le degre m des
* sous corps cherches, cree un vecteur ligne a 13 composantes:
* 1 : le polynome nf[1],
* 2 : le premier p,
* 3 : la factorisation ff,
* 4 : la longeur des cycles associes (n_1,...,n_r),
* 5 : les cycles associes,
* 6 : le corps F_(p^q),
* 7 : les racines de f dans F_(p^q) par facteur de ff,
* 8 : la borne M pour les sous-corps,
* 9 : l'exposant e telle que la precision des lifts soit p^e>2.M,
* 10: le lift de Hensel a la precision p^e de la factorisation en facteurs
* lineaires de nf[1] dans F_(p^q),
* 11: la borne de Hadamard pour les coefficients de h(x) tel que g o h = 0
* mod nf[1].
* ces donnees sont valides pour nf, p et m (d) donnes...
*/
static GEN
compute_data(GEN nf, GEN ff, GEN p, long m, long nn)
{
long i,j,l,r,*n,e,N,pp,d,r1;
GEN DATA,p1,p2,cys,fhk,tabroots,MM,fk,dpol,maxroot,maxMM,pol;
if (DEBUGLEVEL>1) { fprintferr("Entering compute_data()\n\n"); flusherr(); }
pol = (GEN)nf[1]; N = lgef(pol)-3;
DATA=cgetg(14,t_VEC);
DATA[1]=(long)pol;
DATA[2]=(long)p; r=lg(ff)-1;
DATA[3]=(long)ff;
n = cgetg(r+1, t_VECSMALL);
DATA[4]= (long)n;
for (j=1; j<=r; j++) n[j]=lgef(ff[j])-3;
cys=cgetg(r+1,t_VEC); l=0;
for (i=1; i<=r; i++)
{
p1 = cgetg(n[i]+1, t_VECSMALL);
cys[i] = (long)p1; for (j=1; j<=n[i]; j++) p1[j]=++l;
}
DATA[5]=(long)cys;
DATA[6]=(long)ffinit(p,nn,MAXVARN);
tabroots=cgetg(r+1,t_VEC);
for (j=1; j<=r; j++)
{
p1=(GEN)factmod9((GEN)ff[j],p,(GEN)DATA[6])[1];
p2=cgetg(n[j]+1,t_VEC); tabroots[j]=(long)p2;
p2[1]=lneg(gmael(p1,1,2));
for (i=2; i<=n[j]; i++) p2[i]=(long)powgi((GEN)p2[i-1],p);
}
DATA[7]=(long)tabroots;
r1=itos(gmael(nf,2,1));
MM = bound_for_coeff(m, (GEN)nf[6], r1, &maxroot);
MM = gmul2n(MM,1);
DATA[8]=(long)MM;
pp=itos(p); maxMM = vecmax(MM);
for (e=1,p1=p; cmpii(p1, maxMM) < 0; ) { p1 = mulis(p1,pp); e++; }
DATA[9]=lpuigs(p,e); fk=cgetg(N+1,t_VEC);
for (l=1,j=1; j<=r; j++)
for (i=1; i<=n[j]; i++)
fk[l++] = lsub(polx[0],gmael(tabroots,j,i));
fhk = hensel_lift(pol,fk,(GEN)DATA[6],p,e);
tabroots=cgetg(r+1,t_VEC);
for (l=1,j=1; j<=r; j++)
{
p1 = cgetg(n[j]+1,t_VEC); tabroots[j]=(long)p1;
for (i=1; i<=n[j]; i++,l++) p1[i] = lneg(gmael(fhk,l,2));
}
DATA[10]=(long)tabroots;
d=N/m; p1=gmul(stoi(N), gsqrt(gpuigs(stoi(N-1),N-1),DEFAULTPREC));
p2 = gpuigs(maxroot, d + N*(N-1)/2);
dpol=mulii(sqri((GEN)nf[4]),(GEN)nf[3]);
p1 = gdiv(gmul(p1,p2), gsqrt(absi(dpol),DEFAULTPREC));
p1 = grndtoi(p1, &e);
if (e>=0) p1 = addii(p1, shifti(gun, e));
p1 = shifti(p1, 1);
DATA[11]=(long)p1;
if (DEBUGLEVEL>1)
{
fprintferr("DATA =\n");
fprintferr("f = "); outerr((GEN)DATA[1]);
fprintferr("p = "); outerr((GEN)DATA[2]);
fprintferr("ff = "); outerr((GEN)DATA[3]);
fprintferr("lcy = "); outerr((GEN)DATA[4]);
fprintferr("cys = "); outerr((GEN)DATA[5]);
fprintferr("bigfq = "); outerr((GEN)DATA[6]);
fprintferr("roots = "); outerr((GEN)DATA[7]);
fprintferr("2 * M = "); outerr((GEN)DATA[8]);
fprintferr("p^e = "); outerr((GEN)DATA[9]);
fprintferr("lifted roots = "); outerr((GEN)DATA[10]);
fprintferr("2 * Hadamard bound = "); outerr((GEN)DATA[11]);
}
return DATA;
}
/*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
/*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
/* */
/* AUTOMORPHISMS OF AN ABELIAN NUMBER FIELD */
/* */
/* V. Acciaro and J. Klueners (1996) */
/* */
/*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
/*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
/* calcul du frobenius en p pour le corps abelien defini par le polynome pol,
* par relevement de hensel du frobenius frobp de l'extension des corps
* residuels (frobp est un polynome mod pol a coefficients dans F_p)
*/
static GEN
frobenius(GEN pol,GEN frobp,GEN p,GEN B,GEN d)
{
long av=avma,v0,deg,i,depas;
GEN b0,b1,pold,polp,poldp,w0,w1,g0,g1,unmodp,polpp,v,pp,unmodpp,poldpp,bl0,bl1;
v0=varn(pol); unmodp=gmodulsg(1,p); pold=deriv(pol,v0);
b0=frobp; polp=gmul(unmodp,pol);
poldp=gsubst(deriv(polp,v0),v0,frobp);
w0=ginv(poldp);
bl0=lift(b0); deg=lgef(bl0)-3;
v=cgetg(deg+2,t_VEC);
for (i=1; i<=deg+1; i++)
v[i]=ldiv(centerlift(gmul(d,compo(bl0,deg+2-i))),d);
g0=gtopoly(v,v0);
if (DEBUGLEVEL>2)
{
fprintferr("val. initiales:\n");
fprintferr("b0 = "); outerr(b0);
fprintferr("w0 = "); outerr(w0);
fprintferr("g0 = "); outerr(g0);
}
depas=1; pp=gsqr(p);
for(;;)
{
if (gcmp(pp,B)>0) depas=0;
unmodpp=gmodulsg(1,pp);
polpp=gmul(unmodpp,pol); poldpp=gmul(unmodpp,pold);
b0=gmodulcp(gmul(unmodpp,lift_intern(lift_intern(b0))),polpp);
w0=gmodulcp(gmul(unmodpp,lift_intern(lift_intern(w0))),polpp);
b1=gsub(b0,gmul(w0,gsubst(polpp,v0,b0)));
w1=gmul(w0,gsub(gdeux,gmul(w0,gsubst(poldpp,v0,b1))));
bl1=lift(b1); deg=lgef(bl1)-3;
v=cgetg(deg+2,t_VEC);
for (i=1; i<=deg+1; i++)
v[i]=ldiv(centerlift(gmul(d,compo(bl1,deg+2-i))),d);
g1=gtopoly(v,v0);
if (DEBUGLEVEL>2)
{
fprintferr("pp = "); outerr(pp);
fprintferr("b1 = "); outerr(b1);
fprintferr("w1 = "); outerr(w1);
fprintferr("g1 = "); outerr(g1);
}
if (gegal(g0,g1)) return gerepileupto(av,g1);
pp=gsqr(pp); b0=b1; w0=w1; g0=g1;
if (!depas) err(talker,"the number field is not an Abelian number field");
}
}
static GEN
compute_denom(GEN dpol)
{
long av=avma,lf,i,a;
GEN d,f1,f2, f = decomp(dpol);
f1=(GEN)f[1];
f2=(GEN)f[2]; lf=lg(f1);
for (d=gun,i=1; i<lf; i++)
{
a = itos((GEN)f2[i]) >> 1;
d = mulii(d, gpuigs((GEN)f1[i],a));
}
return gerepileupto(av,d);
}
static GEN
compute_bound_for_lift(GEN pol,GEN dpol,GEN d)
{
long av=avma,n,i;
GEN p1,p2,p3;
n=lgef(pol)-3;
p1=gdiv(gmul(stoi(n),gpui(stoi(n-1),gdivgs(stoi(n-1),2),DEFAULTPREC)),
gsqrt(dpol,DEFAULTPREC));
p2=gzero;
for (i=2; i<=n+2; i++) p2=gadd(p2,gsqr((GEN)pol[i]));
p2=gpuigs(gsqrt(p2,DEFAULTPREC),n-1);
p1=gmul(p1,p2); p2=gzero;
for (i=2; i<=n+2; i++)
{
p3 = gabs((GEN)pol[i],DEFAULTPREC);
if (gcmp(p3,p2)>0) p2 = p3;
}
p2=gmul(d,gadd(gun,p2));
return gerepileupto(av, gmul2n(gsqr(gmul(p1,p2)),1));
/* Borne heuristique de P. S. Wang, Math. Comp. 30, 1976, p. 332
p2=gzero; for (i=2; i<=n+2; i++) p2=gadd(p2,gsqr((GEN)pol[i]));
p1=gzero;
for (i=2; i<=n+2; i++){ if (gcmp(gabs((GEN)pol[i],4),p1)>0) p1=gabs((GEN)pol[i],4); }
if (gcmp(p2,p1)>0) p1=p2;
p2=gmul(gdiv(mpfactr(n,4),gsqr(mpfactr(n/2,4))),d);
B=gmul(p1,p2);
tetpil=avma; return gerepile(av,tetpil,gcopy(B));
*/
}
static long
isinlist(GEN T,long longT,GEN x)
{
long i;
for (i=1; i<=longT; i++)
if (gegal(x,(GEN)T[i])) return i;
return 0;
}
/* renvoie 0 si frobp n'est pas dans la liste T; sinon le no de frobp dans T */
static long
isinlistmodp(GEN T,long longT,GEN frobp,GEN p)
{
long av=avma,i;
GEN p1,p2,unmodp;
p1=lift_intern(lift_intern(frobp)); unmodp=gmodulsg(1,p);
for (i=1; i<=longT; i++)
{
p2=lift_intern(gmul(unmodp,(GEN)T[i]));
if (gegal(p2,p1)) { avma=av; return i; }
}
avma=av; return 0;
}
/* renvoie le plus petit f tel que frobp^f est dans la liste T */
static long
minimalexponent(GEN T,long longT,GEN frobp,GEN p,long N)
{
long av=avma,i;
GEN p1 = frobp;
for (i=1; i<=N; i++)
{
if (isinlistmodp(T,longT,p1,p)) {avma=av; return i;}
p1 = gpui(p1,p,DEFAULTPREC);
}
err(talker,"missing frobenius (field not abelian ?)");
return 0; /* not reached */
}
/* Computation of all the automorphisms of the abelian number field
defined by the monic irreducible polynomial pol with integral coefficients */
GEN
conjugates(GEN pol)
{
long av,tetpil,N,i,j,pp,bound_primes,nbprimes,longT,v0,flL,f,longTnew,*tab,nop,flnf;
GEN T,S,p1,p2,p,dpol,modunp,polp,xbar,frobp,frob,d,B,nf;
byteptr di;
if (DEBUGLEVEL>2){ fprintferr("** Entree dans conjugates\n"); flusherr(); }
flnf=0; if (typ(pol)!=t_POL){ nf=checknf(pol); flnf=1; pol=(GEN)nf[1]; }
av=avma; N=lgef(pol)-3; v0=varn(pol);
if (N==1) { S=cgetg(2,t_VEC); S[1]=(long)polx[v0]; return S; }
if (N==2)
{
S=cgetg(3,t_VEC); S[1]=(long)polx[v0];
S[2]=lsub(gneg(polx[v0]),(GEN)pol[3]);
tetpil=avma; return gerepile(av,tetpil,gcopy(S));
}
dpol=absi(discsr(pol));
if (DEBUGLEVEL>2)
{ fprintferr("discriminant du polynome: "); outerr(dpol); }
d = flnf? (GEN)nf[4]: compute_denom(dpol);
if (DEBUGLEVEL>2)
{ fprintferr("facteur carre du discriminant: "); outerr(d); }
B=compute_bound_for_lift(pol,dpol,d);
if (DEBUGLEVEL>2) { fprintferr("borne pour les lifts: "); outerr(B); }
/* sous GRH il faut en fait 3.47*log(dpol) */
p1=gfloor(glog(dpol,DEFAULTPREC));
bound_primes=itos(p1);
if (DEBUGLEVEL>2)
{ fprintferr("borne pour les premiers: %ld\n",bound_primes); flusherr(); }
nbprimes=itos(gfloor(gmul(dbltor(1.25506),
gdiv(p1,glog(p1,DEFAULTPREC)))));
if (DEBUGLEVEL>2)
{ fprintferr("borne pour le nombre de premiers: %ld\n",nbprimes); flusherr(); }
S=cgetg(nbprimes+1,t_VEC);
di=diffptr; pp=*di; i=0;
while (pp<=bound_primes)
{
if (smodis(dpol,pp)) { i++; S[i]=lstoi(pp); }
pp = pp + (*(++di));
}
for (j=i+1; j<=nbprimes; j++) S[j]=zero;
nbprimes=i; tab=new_chunk(nbprimes+1);
for (i=1; i<=nbprimes; i++) tab[i]=0;
if (DEBUGLEVEL>2)
{
fprintferr("nombre de premiers: %ld\n",nbprimes);
fprintferr("table des premiers: "); outerr(S);
}
T=cgetg(N+1,t_VEC); T[1]=(long)polx[v0];
for (i=2; i<=N; i++) T[i]=zero; longT=1;
if (DEBUGLEVEL>2) { fprintferr("table initiale: "); outerr(T); }
for(;;)
{
do
{
do
{
nop = 1+itos(shifti(mulss(mymyrand(),nbprimes),-(BITS_IN_RANDOM-1)));
}
while (tab[nop]);
tab[nop]=1; p=(GEN)S[nop];
if (DEBUGLEVEL>2) { fprintferr("\nnombre premier: "); outerr(p); }
modunp=gmodulsg(1,p);
polp=gmul(modunp,pol);
xbar=gmodulcp(gmul(polx[v0],modunp),polp);
frobp=gpui(xbar,p,4);
if (DEBUGLEVEL>2) { fprintferr("frobenius mod p: "); outerr(frobp); }
flL=isinlistmodp(T,longT,frobp,p);
if (DEBUGLEVEL>2){ fprintferr("flL: %ld\n",flL); flusherr(); }
}
while (flL);
f=minimalexponent(T,longT,frobp,p,N);
if (DEBUGLEVEL>2){ fprintferr("exposant minimum: %ld\n",f); flusherr(); }
frob=frobenius(pol,frobp,p,B,d);
if (DEBUGLEVEL>2) { fprintferr("frobenius: "); outerr(frob); }
/* Ce passage n'est vrai que si le corps est abelien !! */
longTnew=longT;
p2=gmodulcp(frob,pol);
for (i=1; i<=longTnew; i++)
for (j=1; j<f; j++)
{
p1=lift(gsubst((GEN)T[i],v0,gpuigs(p2,j)));
if (DEBUGLEVEL>2)
{
fprintferr("test de la puissance (%ld,%ld): ",i,j); outerr(p1);
}
if (!isinlist(T,longTnew,p1))
{
longT++; T[longT]=(long)p1;
if (longT==N)
{
if (DEBUGLEVEL>2)
{ fprintferr("** Sortie de conjugates\n"); flusherr(); }
tetpil=avma; return gerepile(av,tetpil,gcopy(T));
}
}
}
if (DEBUGLEVEL>2) { fprintferr("nouvelle table: "); outerr(T); }
}
}
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to subfield.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / pari / src / modules