version 1.1, 2001/10/02 11:17:11 |
version 1.2, 2002/09/11 07:27:06 |
Line 21 Foundation, Inc., 59 Temple Place - Suite 330, Boston, |
|
Line 21 Foundation, Inc., 59 Temple Place - Suite 330, Boston, |
|
#include "pari.h" |
#include "pari.h" |
#include "parinf.h" |
#include "parinf.h" |
|
|
|
extern void init_dalloc(); |
|
extern double *dalloc(size_t n); |
|
extern GEN gmul_mati_smallvec(GEN x, GEN y); |
|
extern GEN GS_norms(GEN B, long prec); |
|
extern GEN RXQX_divrem(GEN x, GEN y, GEN T, GEN *pr); |
|
extern GEN RXQX_red(GEN P, GEN T); |
|
extern GEN apply_kummer(GEN nf,GEN pol,GEN e,GEN p); |
|
extern GEN centermod_i(GEN x, GEN p, GEN ps2); |
|
extern GEN dim1proj(GEN prh); |
extern GEN hensel_lift_fact(GEN pol, GEN fact, GEN T, GEN p, GEN pev, long e); |
extern GEN hensel_lift_fact(GEN pol, GEN fact, GEN T, GEN p, GEN pev, long e); |
extern GEN nf_get_T2(GEN base, GEN polr); |
extern GEN initgaloisborne(GEN T, GEN dn, long prec, GEN *ptL, GEN *ptprep, GEN *ptdis); |
extern GEN nfreducemodpr_i(GEN x, GEN prh); |
extern GEN max_modulus(GEN p, double tau); |
|
extern GEN mulmat_pol(GEN A, GEN x); |
|
extern GEN nfgcd(GEN P, GEN Q, GEN nf, GEN den); |
|
extern GEN polsym_gen(GEN P, GEN y0, long n, GEN T, GEN N); |
extern GEN sort_factor(GEN y, int (*cmp)(GEN,GEN)); |
extern GEN sort_factor(GEN y, int (*cmp)(GEN,GEN)); |
extern GEN pidealprimeinv(GEN nf, GEN x); |
extern GEN special_pivot(GEN x); |
|
extern GEN trivfact(void); |
|
extern GEN vandermondeinverse(GEN L, GEN T, GEN den, GEN prep); |
|
extern GEN vconcat(GEN A, GEN B); |
|
extern int cmbf_precs(GEN q, GEN A, GEN B, long *a, long *b, GEN *qa, GEN *qb); |
|
extern int isrational(GEN x); |
|
extern GEN LLL_check_progress(GEN Bnorm, long n0, GEN m, int final, pari_timer *T, long *ti_LLL); |
|
extern void remake_GM(GEN nf, nffp_t *F, long prec); |
|
#define RXQX_div(x,y,T) RXQX_divrem((x),(y),(T),NULL) |
|
#define RXQX_rem(x,y,T) RXQX_divrem((x),(y),(T),ONLY_REM) |
|
|
static GEN nffactormod2(GEN nf,GEN pol,GEN pr); |
|
static GEN nfmod_split2(GEN nf, GEN prhall, GEN polb, GEN v, GEN q); |
|
static GEN nf_pol_mul(GEN nf,GEN pol1,GEN pol2); |
|
static GEN nf_pol_divres(GEN nf,GEN pol1,GEN pol2, GEN *pr); |
|
static GEN nf_pol_subres(GEN nf,GEN pol1,GEN pol2); |
|
static GEN nfmod_pol_reduce(GEN nf,GEN prhall,GEN pol); |
|
static GEN nfmod_pol_divres(GEN nf,GEN prhall,GEN pol1,GEN pol2, GEN *pr); |
|
static GEN nfmod_pol_gcd(GEN nf,GEN prhall,GEN pol1,GEN pol2); |
|
static GEN nf_bestlift(GEN id,GEN idinv,GEN den,GEN elt); |
|
static GEN nf_pol_lift(GEN id,GEN idinv,GEN den,GEN pol); |
|
static GEN nfsqff(GEN nf,GEN pol,long fl); |
static GEN nfsqff(GEN nf,GEN pol,long fl); |
static int nf_combine_factors(GEN nf,long fxn,GEN psf,long dlim,long hint); |
|
|
|
typedef struct Nfcmbf /* for use in nfsqff */ |
/* for nf_bestlift: reconstruction of algebraic integers known mod \wp^k */ |
|
/* FIXME: assume \wp has degree 1 for now */ |
|
typedef struct { |
|
long k; /* input known mod \wp^k */ |
|
GEN pk; /* p^k = Norm(\wp^k)*/ |
|
GEN den; /* denom(prk^-1) = p^k */ |
|
GEN prk; /* |.|^2 LLL-reduced basis (b_i) of \wp^k (NOT T2-reduced) */ |
|
GEN iprk; /* den * prk^-1 */ |
|
GEN GSmin; /* min |b_i^*|^2 */ |
|
GEN ZpProj;/* projector to Zp / \wp^k = Z/p^k (\wp unramified, degree 1) */ |
|
GEN tozk; |
|
GEN topow; |
|
} nflift_t; |
|
|
|
typedef struct /* for use in nfsqff */ |
{ |
{ |
GEN pol, h, hinv, fact, res, lt, den; |
nflift_t *L; |
long nfact, nfactmod; |
GEN nf, pol, fact, res, bound, pr, pk, Br, ZC, dn, polbase, BS_2; |
} Nfcmbf; |
long hint; |
static Nfcmbf nfcmbf; |
} nfcmbf_t; |
|
|
|
GEN |
|
RXQX_mul(GEN x, GEN y, GEN T) |
|
{ |
|
return RXQX_red(gmul(x,y), T); |
|
} |
|
|
static GEN |
static GEN |
unifpol0(GEN nf,GEN pol,long flag) |
unifpol0(GEN nf,GEN pol,long flag) |
{ |
{ |
static long n = 0; |
static long n = 0; |
static GEN vun = NULL; |
static GEN vun = NULL; |
GEN f = (GEN) nf[1]; |
GEN f = (GEN) nf[1]; |
long i = degpol(f), av; |
long i = degpol(f); |
|
gpmem_t av; |
|
|
if (i != n) |
if (i != n) |
{ |
{ |
Line 67 unifpol0(GEN nf,GEN pol,long flag) |
|
Line 98 unifpol0(GEN nf,GEN pol,long flag) |
|
switch(typ(pol)) |
switch(typ(pol)) |
{ |
{ |
case t_INT: case t_FRAC: case t_RFRAC: |
case t_INT: case t_FRAC: case t_RFRAC: |
|
if (flag) return gcopy(pol); |
pol = gmul(pol,vun); break; |
pol = gmul(pol,vun); break; |
|
|
case t_POL: |
case t_POL: |
|
if (flag && !degpol(pol)) return gcopy(constant_term(pol)); |
pol = gmodulcp(pol,f); /* fall through */ |
pol = gmodulcp(pol,f); /* fall through */ |
case t_POLMOD: |
case t_POLMOD: |
pol = algtobasis(nf,pol); |
pol = algtobasis(nf,pol); |
Line 100 unifpol(GEN nf,GEN pol,long flag) |
|
Line 133 unifpol(GEN nf,GEN pol,long flag) |
|
return unifpol0(nf,(GEN) pol, flag); |
return unifpol0(nf,(GEN) pol, flag); |
} |
} |
|
|
#if 0 |
/* factorization of x modulo pr. Assume x already reduced */ |
/* return a monic polynomial of degree d with random coefficients in Z_nf */ |
GEN |
static GEN |
FqX_factor(GEN x, GEN T, GEN p) |
random_pol(GEN nf,long d) |
|
{ |
{ |
long i,j, n = degpol(nf[1]); |
GEN rep; |
GEN pl,p; |
if (!T) |
|
|
pl=cgetg(d+3,t_POL); |
|
for (i=2; i<d+2; i++) |
|
{ |
{ |
p=cgetg(n+1,t_COL); pl[i]=(long)p; |
rep = factmod0(x, p); |
for (j=1; j<=n; j++) |
rep[2] = (long)small_to_vec((GEN)rep[2]); |
p[j] = lstoi(mymyrand()%101 - 50); |
|
} |
} |
p=cgetg(n+1,t_COL); pl[i]=(long)p; |
else |
p[1]=un; for (i=2; i<=n; i++) p[i]=zero; |
{ |
|
rep = factmod9(x, p, T); |
pl[1] = evalsigne(1) | evallgef(d+3) | evalvarn(0); |
rep = lift_intern(lift_intern(rep)); |
return pl; |
} |
|
return rep; |
} |
} |
#endif |
|
|
|
/* multiplication of x by y */ |
GEN |
static GEN |
nffactormod(GEN nf, GEN x, GEN pr) |
nf_pol_mul(GEN nf,GEN x,GEN y) |
|
{ |
{ |
long tetpil,av=avma; |
long j, l, vx = varn(x), vn; |
GEN res = gmul(unifpol(nf,x,1), unifpol(nf,y,1)); |
gpmem_t av = avma; |
|
GEN z, rep, xrd, modpr, T, p; |
|
|
tetpil = avma; |
nf = checknf(nf); |
return gerepile(av,tetpil,unifpol(nf,res,0)); |
vn = varn((GEN)nf[1]); |
|
if (typ(x)!=t_POL) err(typeer,"nffactormod"); |
|
if (vx >= vn) |
|
err(talker,"polynomial variable must have highest priority in nffactormod"); |
|
|
|
modpr = nf_to_ff_init(nf, &pr, &T, &p); |
|
xrd = modprX(x, nf, modpr); |
|
rep = FqX_factor(xrd,T,p); |
|
z = (GEN)rep[1]; l = lg(z); |
|
for (j = 1; j < l; j++) z[j] = (long)modprX_lift((GEN)z[j], modpr); |
|
return gerepilecopy(av, rep); |
} |
} |
|
|
/* compute x^2 in nf */ |
/* If p doesn't divide bad and has a divisor of degree 1, return the latter. */ |
static GEN |
static GEN |
nf_pol_sqr(GEN nf,GEN x) |
choose_prime(GEN nf, GEN bad, GEN *p, byteptr *PT) |
{ |
{ |
long tetpil,av=avma; |
GEN q = icopy(gun), r, x = (GEN)nf[1]; |
GEN res = gsqr(unifpol(nf,x,1)); |
ulong pp = *p? itou(*p): 0; |
|
byteptr pt = *PT; |
tetpil = avma; |
gpmem_t av = avma; |
return gerepile(av,tetpil,unifpol(nf,res,0)); |
for (;;) |
|
{ |
|
NEXT_PRIME_VIADIFF_CHECK(pp, pt); |
|
if (! smodis(bad,pp)) continue; |
|
affui(pp, q); |
|
r = rootmod(x, q); if (lg(r) > 1) break; |
|
avma = av; |
|
} |
|
r = gsub(polx[varn(x)], lift_intern((GEN)r[1])); |
|
*p = q; |
|
*PT = pt; return apply_kummer(nf,r,gun,q); |
} |
} |
|
|
/* reduce the coefficients of pol modulo prhall */ |
|
static GEN |
static GEN |
nfmod_pol_reduce(GEN nf,GEN prhall,GEN pol) |
QXQ_normalize(GEN P, GEN T) |
{ |
{ |
long av=avma,tetpil,i; |
GEN t = leading_term(P); |
GEN p1; |
if (!gcmp1(t)) |
|
{ |
if (typ(pol)!=t_POL) return nfreducemodpr(nf,pol,prhall); |
if (is_rational_t(typ(t))) |
pol=unifpol(nf,pol,0); |
P = gdiv(P, t); |
|
else |
tetpil=avma; i=lgef(pol); |
P = RXQX_mul(QX_invmod(t,T), P, T); |
p1=cgetg(i,t_POL); p1[1]=pol[1]; |
} |
for (i--; i>=2; i--) |
return P; |
p1[i] = (long) nfreducemodpr(nf,(GEN)pol[i],prhall); |
|
return gerepile(av,tetpil, normalizepol(p1)); |
|
} |
} |
|
|
/* x^2 modulo prhall */ |
/* return the roots of pol in nf */ |
static GEN |
GEN |
nfmod_pol_sqr(GEN nf,GEN prhall,GEN x) |
nfroots(GEN nf,GEN pol) |
{ |
{ |
long av=avma,tetpil; |
gpmem_t av = avma; |
GEN px; |
int d = degpol(pol); |
|
GEN A,g, T; |
|
|
px = nfmod_pol_reduce(nf,prhall,x); |
nf = checknf(nf); T = (GEN)nf[1]; |
px = unifpol(nf,lift(px),1); |
if (typ(pol) != t_POL) err(notpoler,"nfroots"); |
px = unifpol(nf,nf_pol_sqr(nf,px),0); |
if (varn(pol) >= varn(T)) |
tetpil=avma; |
err(talker,"polynomial variable must have highest priority in nfroots"); |
return gerepile(av,tetpil,nfmod_pol_reduce(nf,prhall,px)); |
if (d == 0) return cgetg(1,t_VEC); |
|
if (d == 1) |
|
{ |
|
A = gneg_i(gdiv((GEN)pol[2],(GEN)pol[3])); |
|
return gerepilecopy(av, _vec( basistoalg(nf,A) )); |
|
} |
|
A = fix_relative_pol(nf,pol,0); |
|
A = primpart( lift_intern(A) ); |
|
if (DEBUGLEVEL>3) fprintferr("test if polynomial is square-free\n"); |
|
g = nfgcd(A, derivpol(A), T, NULL); |
|
|
|
if (degpol(g)) |
|
{ /* not squarefree */ |
|
g = QXQ_normalize(g, T); |
|
A = RXQX_div(A,g,T); |
|
} |
|
A = QXQ_normalize(A, T); |
|
A = primpart(A); |
|
A = nfsqff(nf,A,1); |
|
return gerepileupto(av, gen_sort(A, 0, cmp_pol)); |
} |
} |
|
|
/* multiplication of x by y modulo prhall */ |
int |
static GEN |
nfissplit(GEN nf, GEN x) |
nfmod_pol_mul(GEN nf,GEN prhall,GEN x,GEN y) |
|
{ |
{ |
long av=avma,tetpil; |
gpmem_t av = avma; |
GEN px,py; |
long l = lg(nfsqff(checknf(nf), x, 2)); |
|
avma = av; return l != 1; |
px = nfmod_pol_reduce(nf,prhall,x); px = unifpol(nf,lift(px),1); |
|
py = nfmod_pol_reduce(nf,prhall,y); py = unifpol(nf,lift(py),1); |
|
px = unifpol(nf,nf_pol_mul(nf,px,py),0); |
|
tetpil=avma; |
|
return gerepile(av,tetpil,nfmod_pol_reduce(nf,prhall,px)); |
|
} |
} |
|
|
/* Euclidan division of x by y */ |
/* nf = K[y] / (P). Galois over K ? */ |
static GEN |
int |
nf_pol_divres(GEN nf,GEN x,GEN y,GEN *pr) |
nfisgalois(GEN nf, GEN P) |
{ |
{ |
long av = avma,tetpil; |
return degpol(P) <= 2 || nfissplit(nf, P); |
GEN nq = poldivres(unifpol(nf,x,1),unifpol(nf,y,1),pr); |
|
GEN *gptr[2]; |
|
|
|
tetpil=avma; nq=unifpol(nf,nq,0); |
|
if (pr) *pr = unifpol(nf,*pr,0); |
|
gptr[0]=&nq; gptr[1]=pr; |
|
gerepilemanysp(av,tetpil,gptr,pr ? 2:1); |
|
return nq; |
|
} |
} |
|
|
/* Euclidan division of x by y modulo prhall */ |
/* return a minimal lift of elt modulo id */ |
static GEN |
static GEN |
nfmod_pol_divres(GEN nf,GEN prhall,GEN x,GEN y, GEN *pr) |
nf_bestlift(GEN elt, GEN bound, nflift_t *T) |
{ |
{ |
long av=avma,dx,dy,dz,i,j,k,l,n,tetpil; |
GEN u; |
GEN z,p1,p3,px,py; |
long i,l = lg(T->prk), t = typ(elt); |
|
if (t != t_INT) |
py = nfmod_pol_reduce(nf,prhall,y); |
|
if (gcmp0(py)) |
|
err(talker, "division by zero in nfmod_pol_divres"); |
|
|
|
tetpil=avma; |
|
px=nfmod_pol_reduce(nf,prhall,x); |
|
dx=degpol(px); dy=degpol(py); dz=dx-dy; |
|
if (dx<dy) |
|
{ |
{ |
GEN vzero; |
if (t == t_POL) elt = mulmat_pol(T->tozk, elt); |
|
u = gmul(T->iprk,elt); |
if (pr) *pr = gerepile(av,tetpil,px); |
for (i=1; i<l; i++) u[i] = (long)diviiround((GEN)u[i], T->pk); |
else avma = av; |
|
|
|
n=degpol(nf[1]); |
|
vzero = cgetg(n+1,t_COL); |
|
n=degpol(nf[1]); |
|
for (i=1; i<=n; i++) vzero[i]=zero; |
|
|
|
z=cgetg(3,t_POL); z[2]=(long)vzero; |
|
z[1]=evallgef(2) | evalvarn(varn(px)); |
|
return z; |
|
} |
} |
p1 = NULL; /* gcc -Wall */ |
else |
|
|
z=cgetg(dz+3,t_POL); z[1]=evalsigne(1) | evallgef(3+dz); |
|
setvarn(z,varn(px)); |
|
z[dz+2] = (long) element_divmodpr(nf,(GEN)px[dx+2],(GEN)py[dy+2],prhall); |
|
for (i=dx-1; i>=dy; --i) |
|
{ |
{ |
l=avma; p1=(GEN)px[i+2]; |
u = gmul(elt, (GEN)T->iprk[1]); |
for (j=i-dy+1; j<=i && j<=dz; j++) |
for (i=1; i<l; i++) u[i] = (long)diviiround((GEN)u[i], T->pk); |
p1 = gsub(p1, element_mul(nf,(GEN)z[j+2],(GEN)py[i-j+2])); |
elt = gscalcol(elt, l-1); |
p1 = nfreducemodpr(nf,p1,prhall); |
|
tetpil=avma; p3=element_divmodpr(nf,p1,(GEN)py[dy+2],prhall); |
|
z[i-dy+2]=lpile(l,tetpil,p3); |
|
z[i-dy+2]=(long)nfreducemodpr(nf,(GEN)z[i-dy+2],prhall); |
|
} |
} |
l=avma; |
u = gsub(elt, gmul(T->prk, u)); |
for (i=dy-1; i>=0; --i) |
if (bound && gcmp(QuickNormL2(u,DEFAULTPREC), bound) > 0) u = NULL; |
{ |
return u; |
l=avma; p1=((GEN)px[i+2]); |
|
for (j=0; j<=i && j<=dz; j++) |
|
p1 = gsub(p1, element_mul(nf,(GEN)z[j+2],(GEN)py[i-j+2])); |
|
p1 = gerepileupto(l, nfreducemodpr(nf,p1,prhall)); |
|
if (!gcmp0(p1)) break; |
|
} |
|
|
|
if (!pr) { avma = l; return z; } |
|
|
|
if (i<0) |
|
{ |
|
avma=l; |
|
p3 = cgetg(3,t_POL); p3[2]=zero; |
|
p3[1] = evallgef(2) | evalvarn(varn(px)); |
|
*pr=p3; return z; |
|
} |
|
|
|
p3=cgetg(i+3,t_POL); |
|
p3[1]=evalsigne(1) | evallgef(i+3) | evalvarn(varn(px)); |
|
p3[i+2]=(long)nfreducemodpr(nf,p1,prhall); |
|
for (k=i-1; k>=0; --k) |
|
{ |
|
l=avma; p1=((GEN)px[k+2]); |
|
for (j=0; j<=k && j<=dz; j++) |
|
p1 = gsub(p1, element_mul(nf,(GEN)z[j+2],(GEN)py[k-j+2])); |
|
p3[k+2]=lpileupto(l,nfreducemodpr(nf,p1,prhall)); |
|
} |
|
*pr=p3; return z; |
|
} |
} |
|
|
/* GCD of x and y */ |
|
static GEN |
static GEN |
nf_pol_subres(GEN nf,GEN x,GEN y) |
nf_bestlift_to_pol(GEN elt, GEN bound, nflift_t *T) |
{ |
{ |
long av=avma,tetpil; |
GEN u = nf_bestlift(elt,bound,T); |
GEN s = srgcd(unifpol(nf,x,1), unifpol(nf,y,1)); |
if (u) u = gmul(T->topow, u); |
|
return u; |
tetpil=avma; return gerepile(av,tetpil,unifpol(nf,s,1)); |
|
} |
} |
|
|
/* GCD of x and y modulo prhall */ |
/* return the lift of pol with coefficients of T2-norm <= C (if possible) */ |
static GEN |
static GEN |
nfmod_pol_gcd(GEN nf,GEN prhall,GEN x,GEN y) |
nf_pol_lift(GEN pol, GEN bound, nfcmbf_t *T) |
{ |
{ |
long av=avma; |
long i, d = lgef(pol); |
GEN p1,p2; |
GEN x = cgetg(d,t_POL); |
|
nflift_t *L = T->L; |
|
|
if (lgef(x)<lgef(y)) { p1=y; y=x; x=p1; } |
x[1] = pol[1]; |
p1=nfmod_pol_reduce(nf,prhall,x); |
for (i=2; i<d; i++) |
p2=nfmod_pol_reduce(nf,prhall,y); |
|
while (!isexactzero(p2)) |
|
{ |
{ |
GEN p3; |
x[i] = (long)nf_bestlift_to_pol((GEN)pol[i], bound, L); |
|
if (!x[i]) return NULL; |
nfmod_pol_divres(nf,prhall,p1,p2,&p3); |
|
p1=p2; p2=p3; |
|
} |
} |
return gerepileupto(av,p1); |
return x; |
} |
} |
|
|
/* return pol^e modulo prhall and the polynomial pmod */ |
/* return the factorization of x in nf */ |
static GEN |
GEN |
nfmod_pol_pow(GEN nf,GEN prhall,GEN pmod,GEN pol,GEN e) |
nffactor(GEN nf,GEN pol) |
{ |
{ |
long i, av = avma, n = degpol(nf[1]); |
GEN A,g,y,p1,rep,T; |
GEN p1,p2,vun; |
long l, j, d = degpol(pol); |
|
gpmem_t av; |
|
if (DEBUGLEVEL>3) (void)timer2(); |
|
|
vun=cgetg(n+1,t_COL); vun[1]=un; for (i=2; i<=n; i++) vun[i]=zero; |
nf = checknf(nf); T = (GEN)nf[1]; |
p1=gcopy(polun[varn(pol)]); p1[2]=(long)vun; |
if (typ(pol) != t_POL) err(notpoler,"nffactor"); |
if (gcmp0(e)) return p1; |
if (varn(pol) >= varn(T)) |
|
err(talker,"polynomial variable must have highest priority in nffactor"); |
|
|
p2=nfmod_pol_reduce(nf,prhall,pol); |
if (d == 0) return trivfact(); |
for(;;) |
rep = cgetg(3, t_MAT); av = avma; |
|
if (d == 1) |
{ |
{ |
if (!vali(e)) |
rep[1] = (long)_col( gcopy(pol) ); |
{ |
rep[2] = (long)_col( gun ); |
p1=nfmod_pol_mul(nf,prhall,p1,p2); |
return rep; |
nfmod_pol_divres(nf,prhall,p1,pmod,&p1); |
|
} |
|
if (gcmp1(e)) break; |
|
|
|
e=shifti(e,-1); |
|
p2=nfmod_pol_sqr(nf,prhall,p2); |
|
nfmod_pol_divres(nf,prhall,p2,pmod,&p2); |
|
} |
} |
return gerepileupto(av,p1); |
|
} |
|
|
|
static long |
A = fix_relative_pol(nf,pol,0); |
isdivbyprime(GEN nf, GEN x, GEN pr) |
if (degpol(nf[1]) == 1) |
{ |
return gerepileupto(av, factpol(simplify(pol), 0)); |
GEN elt, p = (GEN)pr[1], tau = (GEN)pr[5]; |
|
|
|
elt = element_mul(nf, x, tau); |
A = primpart( lift_intern(A) ); |
if (divise(content(elt), p)) return 1; |
if (DEBUGLEVEL>3) fprintferr("test if polynomial is square-free\n"); |
|
g = nfgcd(A, derivpol(A), T, NULL); |
|
|
return 0; |
A = QXQ_normalize(A, T); |
} |
A = primpart(A); |
|
|
/* return the factor of nf.pol modulo p corresponding to pr */ |
if (degpol(g)) |
static GEN |
{ /* not squarefree */ |
localpol(GEN nf, GEN pr) |
gpmem_t av1; |
{ |
GEN ex; |
long i, l; |
g = QXQ_normalize(g, T); |
GEN fct, pol = (GEN)nf[1], p = (GEN)pr[1]; |
A = RXQX_div(A,g, T); |
|
|
fct = lift((GEN)factmod(pol, p)[1]); |
y = nfsqff(nf,A,0); av1 = avma; |
l = lg(fct) - 1; |
l = lg(y); |
for (i = 1; i <= l; i++) |
ex=(GEN)gpmalloc(l * sizeof(long)); |
if (isdivbyprime(nf, (GEN)fct[i], pr)) return (GEN)fct[i]; |
for (j=l-1; j>=1; j--) |
|
{ |
err(talker, "cannot find a suitable factor in localpol"); |
GEN fact = lift((GEN)y[j]), quo = g, q; |
return NULL; /* not reached */ |
long e = 0; |
|
for(e = 1;; e++) |
|
{ |
|
q = RXQX_divrem(quo,fact,T, ONLY_DIVIDES); |
|
if (!q) break; |
|
quo = q; |
|
} |
|
ex[j] = e; |
|
} |
|
avma = av1; y = gerepileupto(av, y); |
|
p1 = cgetg(l, t_COL); for (j=l-1; j>=1; j--) p1[j] = lstoi(ex[j]); |
|
free(ex); |
|
} |
|
else |
|
{ |
|
y = gerepileupto(av, nfsqff(nf,A,0)); |
|
l = lg(y); |
|
p1 = cgetg(l, t_COL); for (j=l-1; j>=1; j--) p1[j] = un; |
|
} |
|
if (DEBUGLEVEL>3) |
|
fprintferr("number of factor(s) found: %ld\n", lg(y)-1); |
|
rep[1] = (long)y; |
|
rep[2] = (long)p1; return sort_factor(rep, cmp_pol); |
} |
} |
|
|
/* factorization of x modulo pr */ |
/* return a bound for T_2(P), P | polbase in C[X] |
|
* NB: Mignotte bound: A | S ==> |
|
* |a_i| <= binom(d-1, i-1) || S ||_2 + binom(d-1, i) lc(S) |
|
* |
|
* Apply to sigma(S) for all embeddings sigma, then take the L_2 norm over |
|
* sigma, then take the sup over i. |
|
**/ |
static GEN |
static GEN |
nffactormod0(GEN nf, GEN x, GEN pr) |
nf_Mignotte_bound(GEN nf, GEN polbase) |
{ |
{ |
long av = avma, j, l, vx = varn(x), vn; |
GEN G = gmael(nf,5,2), lS = leading_term(polbase); /* t_INT */ |
GEN rep, pol, xrd, prh, p1; |
GEN p1, C, N2, matGS, binlS, bin; |
|
long prec, i, j, d = degpol(polbase), n = degpol(nf[1]), r1 = nf_get_r1(nf); |
|
|
nf=checknf(nf); |
if (typ(lS) == t_COL) lS = (GEN)lS[1]; |
vn = varn((GEN)nf[1]); |
binlS = bin = vecbinome(d-1); |
if (typ(x)!=t_POL) err(typeer,"nffactormod"); |
if (!gcmp1(lS)) binlS = gmul(lS, bin); |
if (vx >= vn) |
|
err(talker,"polynomial variable must have highest priority in nffactormod"); |
|
|
|
prh = nfmodprinit(nf, pr); |
N2 = cgetg(n+1, t_VEC); |
|
prec = gprecision(G); |
|
for (;;) |
|
{ |
|
nffp_t F; |
|
|
if (divise((GEN)nf[4], (GEN)pr[1])) |
matGS = cgetg(d+2, t_MAT); |
return gerepileupto(av, nffactormod2(nf,x,pr)); |
for (j=0; j<=d; j++) matGS[j+1] = lmul(G, (GEN)polbase[j+2]); |
|
matGS = gtrans_i(matGS); |
|
for (j=1; j <= r1; j++) /* N2[j] = || sigma_j(S) ||_2 */ |
|
{ |
|
N2[j] = lsqrt( QuickNormL2((GEN)matGS[j], DEFAULTPREC), DEFAULTPREC ); |
|
if (lg(N2[j]) < DEFAULTPREC) goto PRECPB; |
|
} |
|
for ( ; j <= n; j+=2) |
|
{ |
|
GEN q1 = QuickNormL2((GEN)matGS[j ], DEFAULTPREC); |
|
GEN q2 = QuickNormL2((GEN)matGS[j+1], DEFAULTPREC); |
|
p1 = gmul2n(mpadd(q1, q2), -1); |
|
N2[j] = N2[j+1] = lsqrt( p1, DEFAULTPREC ); |
|
if (lg(N2[j]) < DEFAULTPREC) goto PRECPB; |
|
} |
|
if (j > n) break; /* done */ |
|
PRECPB: |
|
prec = (prec<<1)-2; |
|
remake_GM(nf, &F, prec); G = F.G; |
|
if (DEBUGLEVEL>1) err(warnprec, "nf_factor_bound", prec); |
|
} |
|
|
xrd = nfmod_pol_reduce(nf, prh, x); |
/* Take sup over 0 <= i <= d of |
if (gcmp1((GEN)pr[4])) |
* sum_sigma | binom(d-1, i-1) ||sigma(S)||_2 + binom(d-1,i) lc(S) |^2 */ |
|
|
|
/* i = 0: n lc(S)^2 */ |
|
C = mulsi(n, sqri(lS)); |
|
/* i = d: sum_sigma ||sigma(S)||_2^2 */ |
|
p1 = gnorml2(N2); if (gcmp(C, p1) < 0) C = p1; |
|
for (i = 1; i < d; i++) |
{ |
{ |
pol = gun; /* dummy */ |
GEN s = gzero; |
for (j = 2; j < lg(xrd); j++) |
for (j = 1; j <= n; j++) |
xrd[j] = mael(xrd, j, 1); |
|
rep = factmod(xrd, (GEN)pr[1]); |
|
rep = lift(rep); |
|
} |
|
else |
|
{ |
|
pol = localpol(nf, pr); |
|
xrd = lift(unifpol(nf, xrd, 1)); |
|
p1 = gun; |
|
for (j = 2; j < lg(xrd); j++) |
|
{ |
{ |
xrd[j] = lmod((GEN)xrd[j], pol); |
p1 = mpadd( mpmul((GEN)bin[i], (GEN)N2[j]), (GEN)binlS[i+1] ); |
p1 = mpppcm(p1, denom(content((GEN)xrd[j]))); |
s = mpadd(s, gsqr(p1)); |
} |
} |
rep = factmod9(gmul(xrd, p1), (GEN)pr[1], pol); |
if (gcmp(C, s) < 0) C = s; |
rep = lift(lift(rep)); |
|
} |
} |
|
return C; |
|
} |
|
|
l = lg((GEN)rep[1]); |
/* return a bound for T_2(P), P | polbase |
for (j = 1; j < l; j++) |
* max |b_i|^2 <= 3^{3/2 + d} / (4 \pi d) [P]_2, |
mael(rep, 1, j) = (long)unifpol(nf, gmael(rep, 1, j), 1); |
* where [P]_2 is Bombieri's 2-norm |
|
* Sum over conjugates |
|
*/ |
|
static GEN |
|
nf_Beauzamy_bound(GEN nf, GEN polbase) |
|
{ |
|
GEN lt,C,run,s, G = gmael(nf,5,2), POL, bin; |
|
long i,prec,precnf, d = degpol(polbase), n = degpol(nf[1]); |
|
|
return gerepilecopy(av, rep); |
precnf = gprecision(G); |
|
prec = MEDDEFAULTPREC; |
|
bin = vecbinome(d); |
|
POL = polbase + 2; |
|
/* compute [POL]_2 */ |
|
for (;;) |
|
{ |
|
run= realun(prec); |
|
s = realzero(prec); |
|
for (i=0; i<=d; i++) |
|
{ |
|
GEN p1 = gnorml2(gmul(G, gmul(run, (GEN)POL[i]))); /* T2(POL[i]) */ |
|
if (!signe(p1)) continue; |
|
if (lg(p1) == 3) break; |
|
/* s += T2(POL[i]) / binomial(d,i) */ |
|
s = addrr(s, gdiv(p1, (GEN)bin[i+1])); |
|
} |
|
if (i > d) break; |
|
|
|
prec = (prec<<1)-2; |
|
if (prec > precnf) |
|
{ |
|
nffp_t F; remake_GM(nf, &F, prec); G = F.G; |
|
if (DEBUGLEVEL>1) err(warnprec, "nf_factor_bound", prec); |
|
} |
|
} |
|
lt = (GEN)leading_term(polbase)[1]; |
|
s = gmul(s, mulis(sqri(lt), n)); |
|
C = gpow(stoi(3), dbltor(1.5 + d), DEFAULTPREC); /* 3^{3/2 + d} */ |
|
return gdiv(gmul(C, s), gmulsg(d, mppi(DEFAULTPREC))); |
} |
} |
|
|
GEN |
static GEN |
nffactormod(GEN nf, GEN x, GEN pr) |
nf_factor_bound(GEN nf, GEN polbase) |
{ |
{ |
return nffactormod0(nf, x, pr); |
gpmem_t av = avma; |
|
GEN a = nf_Mignotte_bound(nf, polbase); |
|
GEN b = nf_Beauzamy_bound(nf, polbase); |
|
if (DEBUGLEVEL>2) |
|
{ |
|
fprintferr("Mignotte bound: %Z\n",a); |
|
fprintferr("Beauzamy bound: %Z\n",b); |
|
} |
|
return gerepileupto(av, gmin(a, b)); |
} |
} |
|
|
extern GEN trivfact(void); |
static long |
|
ZXY_get_prec(GEN P) |
/* factorization of x modulo pr */ |
|
GEN |
|
nffactormod2(GEN nf,GEN pol,GEN pr) |
|
{ |
{ |
long av = avma, tetpil,lb,nbfact,psim,N,n,i,j,k,d,e,vf,r,kk; |
long i, j, z, prec = 0; |
GEN y,ex,*t,f1,f2,f3,df1,g1,polb,pold,polu,vker; |
for (i=2; i<lgef(P); i++) |
GEN Q,f,x,u,v,v2,v3,vz,q,vun,vzero,prhall; |
|
|
|
nf=checknf(nf); |
|
if (typ(pol)!=t_POL) err(typeer,"nffactormod"); |
|
if (varn(pol) >= varn(nf[1])) |
|
err(talker,"polynomial variable must have highest priority in nffactormod"); |
|
|
|
prhall=nfmodprinit(nf,pr); n=degpol(nf[1]); |
|
vun = gscalcol_i(gun, n); |
|
vzero = gscalcol_i(gzero, n); |
|
q = vker = NULL; /* gcc -Wall */ |
|
|
|
f=unifpol(nf,pol,0); f=nfmod_pol_reduce(nf,prhall,f); |
|
if (!signe(f)) err(zeropoler,"nffactormod"); |
|
d=degpol(f); vf=varn(f); |
|
if (d == 0) return trivfact(); |
|
t = (GEN*)new_chunk(d+1); ex = cgetg(d+1, t_VECSMALL); |
|
x=gcopy(polx[vf]); x[3]=(long)vun; x[2]=(long)vzero; |
|
if (d<=1) { nbfact=2; t[1] = f; ex[1]=1; } |
|
else |
|
{ |
{ |
q = (GEN)pr[1]; psim = VERYBIGINT; |
GEN p = (GEN)P[i]; |
if (!is_bigint(q)) psim = itos(q); |
if (typ(p) == t_INT) |
/* psim has an effect only when p is small. If too big, set it to a huge |
|
* number (i.e ignore it) to avoid an error in itos on next line. |
|
*/ |
|
q=gpuigs(q, itos((GEN)pr[4])); |
|
f1=f; e=1; nbfact=1; |
|
while (lgef(f1)>3) |
|
{ |
{ |
df1=deriv(f1,vf); f2=nfmod_pol_gcd(nf,prhall,f1,df1); |
z = lgefint(p); |
g1=nfmod_pol_divres(nf,prhall,f1,f2,NULL); k=0; |
if (z > prec) prec = z; |
while (lgef(g1)>3) |
} |
|
else |
|
{ |
|
for (j=2; j<lgef(p); j++) |
{ |
{ |
k++; |
z = lgefint(p[j]); |
if (k%psim == 0) |
if (z > prec) prec = z; |
{ |
|
k++; f2=nfmod_pol_divres(nf,prhall,f2,g1,NULL); |
|
} |
|
f3=nfmod_pol_gcd(nf,prhall,f2,g1); |
|
u = nfmod_pol_divres(nf,prhall,g1,f3,NULL); |
|
f2= nfmod_pol_divres(nf,prhall,f2,f3,NULL); |
|
g1=f3; |
|
if (lgef(u)>3) |
|
{ |
|
N=degpol(u); Q=cgetg(N+1,t_MAT); |
|
v3=cgetg(N+1,t_COL); Q[1]=(long)v3; |
|
v3[1]=(long)vun; for (i=2; i<=N; i++) v3[i]=(long)vzero; |
|
|
|
v2 = v = nfmod_pol_pow(nf,prhall,u,x,q); |
|
for (j=2; j<=N; j++) |
|
{ |
|
v3=cgetg(N+1,t_COL); Q[j]=(long)v3; |
|
for (i=1; i<=lgef(v2)-2; i++) v3[i]=v2[i+1]; |
|
for (; i<=N; i++) v3[i]=(long)vzero; |
|
if (j<N) |
|
{ |
|
v2=nfmod_pol_mul(nf,prhall,v2,v); |
|
nfmod_pol_divres(nf,prhall,v2,u,&v2); |
|
} |
|
} |
|
for (i=1; i<=N; i++) |
|
coeff(Q,i,i)=lsub((GEN)coeff(Q,i,i),vun); |
|
v2=nfkermodpr(nf,Q,prhall); r=lg(v2)-1; t[nbfact]=gcopy(u); kk=1; |
|
if (r>1) |
|
{ |
|
vker=cgetg(r+1,t_COL); |
|
for (i=1; i<=r; i++) |
|
{ |
|
v3=cgetg(N+2,t_POL); |
|
v3[1]=evalsigne(1)+evallgef(2+N); setvarn(v3,vf); |
|
vker[i]=(long)v3; for (j=1; j<=N; j++) v3[j+1]=coeff(v2,j,i); |
|
normalizepol(v3); |
|
} |
|
} |
|
while (kk<r) |
|
{ |
|
v=gcopy(polun[vf]); v[2]=(long)vzero; |
|
for (i=1; i<=r; i++) |
|
{ |
|
vz=cgetg(n+1,t_COL); |
|
for (j=1; j<=n; j++) |
|
vz[j] = lmodsi(mymyrand()>>8, q); |
|
vz=nfreducemodpr(nf,vz,prhall); |
|
v=gadd(v,nfmod_pol_mul(nf,prhall,vz,(GEN)vker[i])); |
|
} |
|
for (i=1; i<=kk && kk<r; i++) |
|
{ |
|
polb=t[nbfact+i-1]; lb=lgef(polb); |
|
if (lb>4) |
|
{ |
|
if(psim==2) |
|
polu=nfmod_split2(nf,prhall,polb,v,q); |
|
else |
|
{ |
|
polu=nfmod_pol_pow(nf,prhall,polb,v,shifti(q,-1)); |
|
polu=gsub(polu,vun); |
|
} |
|
pold=nfmod_pol_gcd(nf,prhall,polb,polu); |
|
if (lgef(pold)>3 && lgef(pold)<lb) |
|
{ |
|
t[nbfact+i-1]=pold; kk++; |
|
t[nbfact+kk-1]=nfmod_pol_divres(nf,prhall,polb,pold,NULL); |
|
} |
|
} |
|
} |
|
} |
|
for (i=nbfact; i<nbfact+r; i++) ex[i]=e*k; |
|
nbfact+=r; |
|
} |
|
} |
} |
e*=psim; j=(degpol(f2))/psim+3; f1=cgetg(j,t_POL); |
|
f1[1] = evalsigne(1) | evallgef(j) | evalvarn(vf); |
|
for (i=2; i<j; i++) |
|
f1[i]=(long)element_powmodpr(nf,(GEN)f2[psim*(i-2)+2], |
|
gdiv(q,(GEN)pr[1]),prhall); |
|
} |
} |
} |
} |
if (nbfact < 2) |
return prec + 1; |
err(talker, "%Z not a prime (nffactormod)", q); |
} |
for (j=1; j<nbfact; j++) |
|
|
long |
|
ZM_get_prec(GEN x) |
|
{ |
|
long i, j, l, k = 2, lx = lg(x); |
|
|
|
for (j=1; j<lx; j++) |
{ |
{ |
v = element_divmodpr(nf,vun,gmael(t,j,lgef(t[j])-1),prhall); |
GEN c = (GEN)x[j]; |
t[j] = unifpol(nf,nfmod_pol_mul(nf,prhall,v,(GEN)t[j]),1); |
for (i=1; i<lx; i++) { l = lgefint(c[i]); if (l > k) k = l; } |
} |
} |
|
return k; |
tetpil=avma; y=cgetg(3,t_MAT); |
|
u=cgetg(nbfact,t_COL); y[1]=(long)u; |
|
v=cgetg(nbfact,t_COL); y[2]=(long)v; |
|
for (j=1,k=0; j<nbfact; j++) |
|
if (ex[j]) |
|
{ |
|
k++; |
|
u[k]=lcopy((GEN)t[j]); |
|
v[k]=lstoi(ex[j]); |
|
} |
|
return gerepile(av,tetpil,y); |
|
} |
} |
|
long |
/* return pol + pol^2 + ... + pol^(q/2) modulo prhall and |
ZX_get_prec(GEN x) |
the polynomial pmod */ |
|
static GEN |
|
nfmod_split2(GEN nf,GEN prhall,GEN pmod,GEN pol,GEN exp) |
|
{ |
{ |
long av = avma; |
long j, l, k = 2, lx = lgef(x); |
GEN p1,p2,q; |
|
|
|
if (cmpis(exp,2)<=0) return pol; |
for (j=2; j<lx; j++) |
p2=p1=pol; q=shifti(exp,-1); |
|
while (!gcmp1(q)) |
|
{ |
{ |
p2=nfmod_pol_sqr(nf,prhall,p2); |
l = lgefint(x[j]); if (l > k) k = l; |
nfmod_pol_divres(nf,prhall,p2,pmod,&p2); |
|
q=shifti(q,-1); p1=gadd(p1,p2); |
|
} |
} |
return gerepileupto(av,p1); |
return k; |
} |
} |
|
|
/* If p doesn't divide either a or b and has a divisor of degree 1, return it. |
|
* Return NULL otherwise. |
|
*/ |
|
static GEN |
static GEN |
p_ok(GEN nf, GEN p, GEN a) |
complex_bound(GEN P) |
{ |
{ |
long av,m,i; |
return gmul(max_modulus(P, 0.01), dbltor(1.0101)); /* exp(0.01) = 1.0100 */ |
GEN dec; |
|
|
|
if (divise(a,p)) return NULL; |
|
av = avma; dec = primedec(nf,p); m=lg(dec); |
|
for (i=1; i<m; i++) |
|
{ |
|
GEN pr = (GEN)dec[i]; |
|
if (is_pm1(pr[4])) |
|
return pr; |
|
} |
|
avma = av; return NULL; |
|
} |
} |
|
|
/* for each new prime ct--, if ct = 0, return NULL */ |
/* return Bs: if r a root of sigma_i(P), |r| < Bs[i] */ |
static GEN |
static GEN |
choose_prime(GEN nf, GEN dk, GEN lim, long ct) |
nf_root_bounds(GEN P, GEN T) |
{ |
{ |
GEN p, pr; |
long lR, i, j, l, prec; |
|
GEN Ps, R, V, nf; |
|
|
p = nextprime(lim); |
if (isrational(P)) return complex_bound(P); |
for (;;) |
|
|
T = get_nfpol(T, &nf); |
|
|
|
prec = ZXY_get_prec(P); |
|
l = lgef(P); |
|
if (nf && nfgetprec(nf) >= prec) |
|
R = (GEN)nf[6]; |
|
else |
|
R = roots(T, prec); |
|
lR = lg(R); |
|
V = cgetg(lR, t_VEC); |
|
Ps = cgetg(l, t_POL); /* sigma (P) */ |
|
Ps[1] = P[1]; |
|
for (j=1; j<lg(R); j++) |
{ |
{ |
if ((pr = p_ok(nf,p,dk))) break; |
GEN r = (GEN)R[j]; |
ct--; |
for (i=2; i<l; i++) Ps[i] = (long)poleval((GEN)P[i], r); |
if (!ct) return NULL; |
V[j] = (long)complex_bound(Ps); |
p = nextprime(addis(p,2)); |
|
} |
} |
|
return V; |
return pr; |
|
} |
} |
|
|
/* test if the discriminant of polbase modulo some few primes |
/* return B such that if x in O_K, K = Z[X]/(T), then the L2-norm of the |
is non-zero. Return 1 if it is so (=> polbase is square-free) |
* coordinates of the numerator of x [on the power, resp. integral, basis if T |
and 0 otherwise (=> polbase may or may not be square-free) */ |
* is a polynomial, resp. an nf] is <= B T_2(x) |
static int |
* *ptden = multiplicative bound for denom(x) */ |
is_sqf(GEN nf, GEN polbase) |
static GEN |
|
L2_bound(GEN T, GEN tozk, GEN *ptden) |
{ |
{ |
GEN lt, pr, prh, p2, p; |
GEN nf, M, L, prep, den, u; |
long i, d = lgef(polbase), ct = 5; |
long prec; |
|
|
lt = (GEN)leading_term(polbase)[1]; |
T = get_nfpol(T, &nf); |
p = stoi(101); |
u = NULL; /* gcc -Wall */ |
|
|
while (ct > 0) |
prec = ZX_get_prec(T) + ZM_get_prec(tozk); |
{ |
den = nf? gun: NULL; |
/* small primes tend to divide discriminants more often |
den = initgaloisborne(T, den, prec, &L, &prep, NULL); |
than large ones so we look at primes >= 101 */ |
M = vandermondeinverse(L, gmul(T, realun(prec)), den, prep); |
pr = choose_prime(nf,lt,p,30); |
if (nf) M = gmul(tozk, M); |
if (!pr) break; |
if (gcmp1(den)) den = NULL; |
|
*ptden = den; |
|
return QuickNormL2(M, DEFAULTPREC); |
|
} |
|
|
p=(GEN)pr[1]; |
/* || L ||_p^p in dimension n (L may be a scalar) */ |
prh=prime_to_ideal(nf,pr); |
static GEN |
|
normlp(GEN L, long p, long n) |
|
{ |
|
long i,l, t = typ(L); |
|
GEN z; |
|
|
p2=gcopy(polbase); |
if (!is_vec_t(t)) return gmulsg(n, gpowgs(L, p)); |
lt=mpinvmod(lt,p); |
|
|
|
for (i=2; i<d; i++) |
l = lg(L); z = gzero; |
p2[i] = nfreducemodpr_i(gmul(lt,(GEN)p2[i]), prh)[1]; |
/* assert(n == l-1); */ |
p2 = normalizepol(p2); |
for (i=1; i<l; i++) |
|
z = gadd(z, gpowgs((GEN)L[i], p)); |
/* discriminant is non-zero => polynomial is square-free */ |
return z; |
if (!gcmp0(p2) && !divise(discsr(p2),p)) { return 1; } |
|
|
|
ct--; |
|
p=addis(p,1); |
|
} |
|
|
|
return 0; |
|
} |
} |
|
|
/* rescale p in K[X] (coeffs in algtobasis form) --> primitive in O_K[X] */ |
/* S = S0 + qS1, P = P0 + qP1 (Euclidean div. by q integer). For a true |
GEN |
* factor (vS, vP), we have: |
nf_pol_to_int(GEN p, GEN *den) |
* | S vS + P vP |^2 < Btra |
|
* This implies | S1 vS + P1 vP |^2 < Blow, assuming q > sqrt(Btra). |
|
* d = dimension of low part (= [nf:Q]) |
|
* n0 = bound for |vS|^2 |
|
* */ |
|
static double |
|
get_Blow(long n0, long d, GEN q) |
{ |
{ |
long i, l = lgef(p); |
#if 0 |
GEN d = gun; |
double sqrtn0 = sqrt((double)n0); |
for (i=2; i<l; i++) |
double t = sqrtn0 + sqrt((double)d) * (sqrtn0 + 1); |
d = mpppcm(d,denom((GEN)p[i])); |
#else |
if (! gcmp1(d)) p = gmul(p, d); |
double sqrtd = sqrt((double)d); |
*den = d; return p; |
double t, aux; |
|
|
|
if (gexpo(q)>30) |
|
aux = 1.0001; |
|
else |
|
{ |
|
aux = gtodouble(q); |
|
aux = sqrt(1 + 4/(aux*aux)); |
|
} |
|
t = n0*sqrtd + sqrtd/2 * aux * (sqrtd * (n0+1)); /* assume pr degree 1 */ |
|
#endif |
|
t = 1. + 0.5 * t; |
|
return t * t; |
} |
} |
|
|
/* return the roots of pol in nf */ |
typedef struct { |
GEN |
GEN d; |
nfroots(GEN nf,GEN pol) |
GEN dPinvS; /* d P^(-1) S [ integral ] */ |
|
double **PinvSdbl; /* P^(-1) S as double */ |
|
GEN S1, P1; /* S = S0 + S1 q, idem P */ |
|
} trace_data; |
|
|
|
/* S1 * u - P1 * round(P^-1 S u). K non-zero coords in u given by ind */ |
|
static GEN |
|
get_trace(GEN ind, trace_data *T) |
{ |
{ |
long av=avma,tetpil,d=lgef(pol),fl; |
long i, j, l, K = lg(ind)-1; |
GEN p1,p2,polbase,polmod,den; |
GEN z, s, v; |
|
|
p2=NULL; /* gcc -Wall */ |
s = (GEN)T->S1[ ind[1] ]; |
nf=checknf(nf); |
if (K == 1) return s; |
if (typ(pol)!=t_POL) err(talker,"not a polynomial in nfroots"); |
|
if (varn(pol) >= varn(nf[1])) |
|
err(talker,"polynomial variable must have highest priority in nfroots"); |
|
|
|
polbase=unifpol(nf,pol,0); |
/* compute s = S1 u */ |
|
for (j=2; j<=K; j++) s = gadd(s, (GEN)T->S1[ ind[j] ]); |
|
|
if (d==3) |
/* compute v := - round(P^1 S u) */ |
|
l = lg(s); |
|
v = cgetg(l, t_VECSMALL); |
|
for (i=1; i<l; i++) |
{ |
{ |
tetpil=avma; p1=cgetg(1,t_VEC); |
double r, t = 0.; |
return gerepile(av,tetpil,p1); |
/* quick approximate computation */ |
|
for (j=1; j<=K; j++) t += T->PinvSdbl[ ind[j] ][i]; |
|
r = floor(t + 0.5); |
|
if (fabs(t + 0.5 - r) < 0.01) |
|
{ /* dubious, compute exactly */ |
|
z = gzero; |
|
for (j=1; j<=K; j++) z = addii(z, ((GEN**)T->dPinvS)[ ind[j] ][i]); |
|
v[i] = - itos( diviiround(z, T->d) ); |
|
} |
|
else |
|
v[i] = - (long)r; |
} |
} |
|
return gadd(s, gmul_mati_smallvec(T->P1, v)); |
|
} |
|
|
if (d==4) |
static trace_data * |
{ |
init_trace(trace_data *T, GEN S, nflift_t *L, GEN q) |
tetpil=avma; p1=cgetg(2,t_VEC); |
{ |
p1[1] = (long)basistoalg(nf,gneg_i( |
long e = gexpo((GEN)S), i,j, l,h; |
element_div(nf,(GEN)polbase[2],(GEN)polbase[3]))); |
GEN qgood, S1, invd; |
return gerepile(av,tetpil,p1); |
|
} |
|
|
|
p1=element_inv(nf,leading_term(polbase)); |
if (e < 0) return NULL; /* S = 0 */ |
polbase=nf_pol_mul(nf,p1,polbase); |
|
|
|
polbase = nf_pol_to_int(polbase, &den); |
qgood = shifti(gun, e - 32); /* single precision check */ |
polmod=unifpol(nf,polbase,1); |
if (cmpii(qgood, q) > 0) q = qgood; |
|
|
if (DEBUGLEVEL>=4) |
S1 = gdivround(S, q); |
fprintferr("test if the polynomial is square-free\n"); |
if (gcmp0(S1)) return NULL; |
|
|
fl = is_sqf(nf, polbase); |
invd = ginv(itor(L->den, DEFAULTPREC)); |
|
|
/* the polynomial may not be square-free ... */ |
T->dPinvS = gmul(L->iprk, S); |
if (!fl) |
l = lg(S); |
|
h = lg(T->dPinvS[1]); |
|
T->PinvSdbl = (double**)cgetg(l, t_MAT); |
|
init_dalloc(); |
|
for (j = 1; j < l; j++) |
{ |
{ |
p1=derivpol(polmod); |
double *t = dalloc(h * sizeof(double)); |
p2=nf_pol_subres(nf,polmod,p1); |
GEN c = (GEN)T->dPinvS[j]; |
if (degpol(p2) == 0) fl = 1; |
T->PinvSdbl[j] = t; |
|
for (i=1; i < h; i++) t[i] = rtodbl(gmul(invd, (GEN)c[i])); |
} |
} |
|
|
if (!fl) |
T->d = L->den; |
{ |
T->P1 = gdivround(L->prk, q); |
p1=element_inv(nf,leading_term(p2)); |
T->S1 = S1; return T; |
p2=nf_pol_mul(nf,p1,p2); |
|
polmod=nf_pol_divres(nf,polmod,p2,NULL); |
|
|
|
p1=element_inv(nf,leading_term(polmod)); |
|
polmod=nf_pol_mul(nf,p1,polmod); |
|
|
|
polmod = nf_pol_to_int(polmod, &den); |
|
polmod=unifpol(nf,polmod,1); |
|
} |
|
|
|
p1 = nfsqff(nf,polmod,1); |
|
tetpil=avma; return gerepile(av, tetpil, gen_sort(p1, 0, cmp_pol)); |
|
} |
} |
|
|
/* return a minimal lift of elt modulo id */ |
static void |
static GEN |
update_trace(trace_data *T, long k, long i) |
nf_bestlift(GEN id,GEN idinv,GEN den,GEN elt) |
|
{ |
{ |
GEN u = gmul(idinv,elt); |
if (!T) return; |
long i, l = lg(u); |
T->S1[k] = T->S1[i]; |
for (i=1; i<l; i++) u[i] = (long)gdivround((GEN)u[i], den); |
T->dPinvS[k] = T->dPinvS[i]; |
return gsub(elt, gmul(id, u)); |
T->PinvSdbl[k] = T->PinvSdbl[i]; |
} |
} |
|
|
/* return the lift of pol with coefficients of T2-norm <= C (if possible) */ |
/* Naive recombination of modular factors: combine up to maxK modular |
|
* factors, degree <= klim and divisible by hint |
|
* |
|
* target = polynomial we want to factor |
|
* famod = array of modular factors. Product should be congruent to |
|
* target/lc(target) modulo p^a |
|
* For true factors: S1,S2 <= p^b, with b <= a and p^(b-a) < 2^31 */ |
static GEN |
static GEN |
nf_pol_lift(GEN id,GEN idinv,GEN den,GEN pol) |
nfcmbf(nfcmbf_t *T, GEN p, long a, long maxK, long klim) |
{ |
{ |
long i, d = lgef(pol); |
long Sbound; |
GEN x = cgetg(d,t_POL); |
GEN pol = T->pol, nf = T->nf, famod = T->fact; |
x[1] = pol[1]; |
GEN bound = T->bound; |
for (i=2; i<d; i++) |
GEN nfpol = (GEN)nf[1]; |
x[i] = (long) nf_bestlift(id,idinv,den,(GEN)pol[i]); |
long K = 1, cnt = 1, i,j,k, curdeg, lfamod = lg(famod)-1, dnf = degpol(nfpol); |
return x; |
GEN lc, lcpol; |
} |
GEN pk = gpowgs(p,a), pas2 = shifti(pk,-1); |
|
|
#if 0 |
GEN trace1 = cgetg(lfamod+1, t_MAT); |
/* return pol(elt) */ |
GEN trace2 = cgetg(lfamod+1, t_MAT); |
static GEN |
GEN ind = cgetg(lfamod+1, t_VECSMALL); |
nf_pol_eval(GEN nf,GEN pol,GEN elt) |
GEN degpol = cgetg(lfamod+1, t_VECSMALL); |
{ |
GEN degsofar = cgetg(lfamod+1, t_VECSMALL); |
long av=avma,tetpil,i; |
GEN listmod = cgetg(lfamod+1, t_COL); |
GEN p1; |
GEN fa = cgetg(lfamod+1, t_COL); |
|
GEN res = cgetg(3, t_VEC); |
|
GEN q = ceil_safe(mpsqrt(T->BS_2)); |
|
const double Blow = get_Blow(lfamod, dnf, q); |
|
trace_data _T1, _T2, *T1, *T2; |
|
|
i=lgef(pol)-1; if (i==2) return gcopy((GEN)pol[2]); |
if (maxK < 0) maxK = lfamod-1; |
|
|
p1=element_mul(nf,(GEN)pol[i],elt); |
lc = absi(leading_term(pol)); |
for (i--; i>=3; i--) |
if (gcmp1(lc)) lc = NULL; |
p1=element_mul(nf,elt,gadd((GEN)pol[i],p1)); |
lcpol = lc? gmul(lc,pol): pol; |
tetpil=avma; return gerepile(av,tetpil,gadd(p1,(GEN)pol[2])); |
|
} |
|
#endif |
|
|
|
/* return the factorization of x in nf */ |
|
GEN |
|
nffactor(GEN nf,GEN pol) |
|
{ |
|
GEN y,p1,p2,den,p3,p4,quot,rep=cgetg(3,t_MAT); |
|
long av = avma,tetpil,i,j,d,fl; |
|
|
|
if (DEBUGLEVEL >= 4) timer2(); |
|
|
|
p3=NULL; /* gcc -Wall */ |
|
nf=checknf(nf); |
|
if (typ(pol)!=t_POL) err(typeer,"nffactor"); |
|
if (varn(pol) >= varn(nf[1])) |
|
err(talker,"polynomial variable must have highest priority in nffactor"); |
|
|
|
d=lgef(pol); |
|
if (d==3) |
|
{ |
{ |
rep[1]=lgetg(1,t_COL); |
GEN t1,t2, lc2 = lc? sqri(lc): NULL; |
rep[2]=lgetg(1,t_COL); |
|
return rep; |
|
} |
|
if (d==4) |
|
{ |
|
p1=cgetg(2,t_COL); rep[1]=(long)p1; p1[1]=lcopy(pol); |
|
p1=cgetg(2,t_COL); rep[2]=(long)p1; p1[1]=un; |
|
return rep; |
|
} |
|
|
|
p1=element_inv(nf,leading_term(pol)); |
for (i=1; i <= lfamod; i++) |
pol=nf_pol_mul(nf,p1,pol); |
{ |
|
GEN P = (GEN)famod[i]; |
|
long d = degpol(P); |
|
|
p1=unifpol(nf,pol,0); |
degpol[i] = d; P += 2; |
p1 = nf_pol_to_int(p1, &den); |
t1 = (GEN)P[d-1];/* = - S_1 */ |
|
t2 = sqri(t1); |
|
if (d > 1) t2 = subii(t2, shifti((GEN)P[d-2],1)); |
|
t2 = modii(t2, pk); /* = S_2 Newton sum */ |
|
if (lc) |
|
{ |
|
t1 = modii(mulii(lc, t1), pk); |
|
t2 = modii(mulii(lc2,t2), pk); |
|
} |
|
trace1[i] = (long)nf_bestlift(t1, NULL, T->L); |
|
trace2[i] = (long)nf_bestlift(t2, NULL, T->L); |
|
} |
|
|
if (DEBUGLEVEL>=4) |
T1 = init_trace(&_T1, trace1, T->L, q); |
fprintferr("test if the polynomial is square-free\n"); |
T2 = init_trace(&_T2, trace2, T->L, q); |
|
|
fl = is_sqf(nf, p1); |
|
|
|
/* polynomial may not be square-free ... */ |
|
if (!fl) |
|
{ |
|
p2=derivpol(p1); |
|
p3=nf_pol_subres(nf,p1,p2); |
|
if (degpol(p3) == 0) fl = 1; |
|
} |
} |
|
degsofar[0] = 0; /* sentinel */ |
|
|
if (!fl) |
/* ind runs through strictly increasing sequences of length K, |
{ |
* 1 <= ind[i] <= lfamod */ |
p4=element_inv(nf,leading_term(p3)); |
nextK: |
p3=nf_pol_mul(nf,p4,p3); |
if (K > maxK || 2*K > lfamod) goto END; |
|
if (DEBUGLEVEL > 3) |
|
fprintferr("\n### K = %d, %Z combinations\n", K,binome(stoi(lfamod), K)); |
|
setlg(ind, K+1); ind[1] = 1; |
|
Sbound = ((K+1)>>1); |
|
i = 1; curdeg = degpol[ind[1]]; |
|
for(;;) |
|
{ /* try all combinations of K factors */ |
|
for (j = i; j < K; j++) |
|
{ |
|
degsofar[j] = curdeg; |
|
ind[j+1] = ind[j]+1; curdeg += degpol[ind[j+1]]; |
|
} |
|
if (curdeg <= klim && curdeg % T->hint == 0) /* trial divide */ |
|
{ |
|
GEN t, y, q, list; |
|
gpmem_t av; |
|
|
p2=nf_pol_divres(nf,p1,p3,NULL); |
av = avma; |
p4=element_inv(nf,leading_term(p2)); |
/* d - 1 test */ |
p2=nf_pol_mul(nf,p4,p2); |
if (T1) |
|
{ |
|
t = get_trace(ind, T1); |
|
if (rtodbl(QuickNormL2(t,DEFAULTPREC)) > Blow) |
|
{ |
|
if (DEBUGLEVEL>6) fprintferr("."); |
|
avma = av; goto NEXT; |
|
} |
|
} |
|
/* d - 2 test */ |
|
if (T2) |
|
{ |
|
t = get_trace(ind, T2); |
|
if (rtodbl(QuickNormL2(t,DEFAULTPREC)) > Blow) |
|
{ |
|
if (DEBUGLEVEL>3) fprintferr("|"); |
|
avma = av; goto NEXT; |
|
} |
|
} |
|
avma = av; |
|
y = lc; /* full computation */ |
|
for (i=1; i<=K; i++) |
|
{ |
|
GEN q = (GEN)famod[ind[i]]; |
|
if (y) q = gmul(y, q); |
|
y = centermod_i(q, pk, pas2); |
|
} |
|
y = nf_pol_lift(y, bound, T); |
|
if (!y) |
|
{ |
|
if (DEBUGLEVEL>3) fprintferr("@"); |
|
avma = av; goto NEXT; |
|
} |
|
/* try out the new combination: y is the candidate factor */ |
|
q = RXQX_divrem(lcpol,y, nfpol, ONLY_DIVIDES); |
|
if (!q) |
|
{ |
|
if (DEBUGLEVEL>3) fprintferr("*"); |
|
avma = av; goto NEXT; |
|
} |
|
|
p2 = nf_pol_to_int(p2, &den); |
/* found a factor */ |
|
list = cgetg(K+1, t_VEC); |
|
listmod[cnt] = (long)list; |
|
for (i=1; i<=K; i++) list[i] = famod[ind[i]]; |
|
|
p2=unifpol(nf,p2,1); |
y = primpart(y); |
tetpil = avma; y = nfsqff(nf,p2,0); |
fa[cnt++] = (long)QXQ_normalize(y, nfpol); |
i = nfcmbf.nfact; |
/* fix up pol */ |
|
pol = q; |
|
if (lc) pol = primpart(pol); |
|
for (i=j=k=1; i <= lfamod; i++) |
|
{ /* remove used factors */ |
|
if (j <= K && i == ind[j]) j++; |
|
else |
|
{ |
|
famod[k] = famod[i]; |
|
update_trace(T1, k, i); |
|
update_trace(T2, k, i); |
|
degpol[k] = degpol[i]; k++; |
|
} |
|
} |
|
lfamod -= K; |
|
if (lfamod < 2*K) goto END; |
|
i = 1; curdeg = degpol[ind[1]]; |
|
if (lc) lc = absi(leading_term(pol)); |
|
lcpol = lc? gmul(lc,pol): pol; |
|
if (DEBUGLEVEL > 2) |
|
{ |
|
fprintferr("\n"); msgtimer("to find factor %Z",y); |
|
fprintferr("remaining modular factor(s): %ld\n", lfamod); |
|
} |
|
continue; |
|
} |
|
|
quot=nf_pol_divres(nf,p1,p2,NULL); |
NEXT: |
p3=(GEN)gpmalloc((i+1) * sizeof(long)); |
for (i = K+1;;) |
for (j=i; j>=1; j--) |
|
{ |
{ |
GEN fact=(GEN)y[j], quo = quot, rem; |
if (--i == 0) { K++; goto nextK; } |
long e = 0; |
if (++ind[i] <= lfamod - K + i) |
|
|
do |
|
{ |
{ |
quo = nf_pol_divres(nf,quo,fact,&rem); |
curdeg = degsofar[i-1] + degpol[ind[i]]; |
e++; |
if (curdeg <= klim) break; |
} |
} |
while (gcmp0(rem)); |
|
p3[j]=lstoi(e); |
|
} |
} |
avma = (long)y; y = gerepile(av, tetpil, y); |
|
p2=cgetg(i+1, t_COL); for (; i>=1; i--) p2[i]=lcopy((GEN)p3[i]); |
|
free(p3); |
|
} |
} |
else |
END: |
{ |
if (degpol(pol) > 0) |
tetpil=avma; y = gerepile(av,tetpil,nfsqff(nf,p1,0)); |
{ /* leftover factor */ |
i = nfcmbf.nfact; |
if (signe(leading_term(pol)) < 0) pol = gneg_i(pol); |
p2=cgetg(i+1, t_COL); for (; i>=1; i--) p2[i]=un; |
|
} |
|
if (DEBUGLEVEL>=4) |
|
fprintferr("number of factor(s) found: %ld\n", nfcmbf.nfact); |
|
rep[1]=(long)y; |
|
rep[2]=(long)p2; return sort_factor(rep, cmp_pol); |
|
} |
|
|
|
/* test if the matrix M is suitable */ |
if (lc && lfamod < 2*K) pol = QXQ_normalize(primpart(pol), nfpol); |
static long |
setlg(famod, lfamod+1); |
test_mat(GEN M, GEN p, GEN C2, long k) |
listmod[cnt] = (long)dummycopy(famod); |
{ |
fa[cnt++] = (long)pol; |
long av = avma, i, N = lg(M); |
|
GEN min, prod, L2, R; |
|
|
|
min = prod = gcoeff(M,1,1); |
|
for (i = 2; i < N; i++) |
|
{ |
|
L2 = gcoeff(M,i,i); prod = mpmul(prod,L2); |
|
if (mpcmp(L2,min) < 0) min = L2; |
|
} |
} |
R = mpmul(min, gpowgs(p, k<<1)); |
if (DEBUGLEVEL>6) fprintferr("\n"); |
i = mpcmp(mpmul(C2,prod), R); avma = av; |
setlg(listmod, cnt); setlg(fa, cnt); |
return (i < 0); |
res[1] = (long)fa; |
|
res[2] = (long)listmod; return res; |
} |
} |
|
|
/* return the matrix corresponding to pr^e with R(pr^e) > C */ |
|
static GEN |
static GEN |
T2_matrix_pow(GEN nf, GEN pre, GEN p, GEN C, long *kmax, long prec) |
nf_check_factors(nfcmbf_t *T, GEN P, GEN M_L, GEN famod, GEN pk) |
{ |
{ |
long av=avma,av1,lim, k = *kmax, N = degpol((GEN)nf[1]); |
GEN nf = T->nf, bound = T->bound; |
int tot_real = !signe(gmael(nf,2,2)); |
GEN nfT = (GEN)nf[1]; |
GEN p1,p2,p3,u,C2,T2, x = (GEN)nf[1]; |
long i, j, r, n0; |
|
GEN pol = P, lcpol, lc, list, piv, y, pas2; |
|
|
C2 = gdiv(gmul2n(C,2), absi((GEN)nf[3])); |
piv = special_pivot(M_L); |
p1 = gmul(pre,lllintpartial(pre)); av1 = avma; |
if (!piv) return NULL; |
T2 = tot_real? gmael(nf,5,4) |
if (DEBUGLEVEL>3) fprintferr("special_pivot output:\n%Z\n",piv); |
: nf_get_T2((GEN) nf[7], roots(x,prec)); |
|
p3 = qf_base_change(T2,p1,1); |
|
|
|
if (N <= 6 && test_mat(p3,p,C2,k)) |
pas2 = shifti(pk,-1); |
|
r = lg(piv)-1; |
|
n0 = lg(piv[1])-1; |
|
list = cgetg(r+1, t_COL); |
|
lc = absi(leading_term(pol)); |
|
if (is_pm1(lc)) lc = NULL; |
|
lcpol = lc? gmul(lc, pol): pol; |
|
for (i=1; i<r; i++) |
{ |
{ |
avma = av1; return gerepileupto(av,p1); |
GEN c = (GEN)piv[i]; |
} |
if (DEBUGLEVEL) fprintferr("nf_LLL_cmbf: checking factor %ld\n",i); |
|
|
av1=avma; lim = stack_lim(av1,2); |
y = lc; |
for (;;) |
for (j=1; j<=n0; j++) |
{ |
if (signe(c[j])) |
if (DEBUGLEVEL>2) fprintferr("exponent: %ld\n",k); |
{ |
|
GEN q = (GEN)famod[j]; |
|
if (y) q = gmul(y, q); |
|
y = centermod_i(q, pk, pas2); |
|
} |
|
y = nf_pol_lift(y, bound, T); |
|
if (!y) return NULL; |
|
|
for(;;) |
/* y is the candidate factor */ |
{ |
pol = RXQX_divrem(lcpol,y,nfT, ONLY_DIVIDES); |
u = tot_real? lllgramall(p3,2,lll_IM) : lllgramintern(p3,2,1,prec); |
if (!pol) return NULL; |
if (u) break; |
|
|
|
prec=(prec<<1)-2; |
y = primpart(y); |
if (DEBUGLEVEL > 1) err(warnprec,"nffactor[1]",prec); |
if (lc) |
T2 = nf_get_T2((GEN) nf[7], roots(x,prec)); |
|
p3 = qf_base_change(T2,p1,1); |
|
} |
|
if (DEBUGLEVEL>2) msgtimer("lllgram + base change"); |
|
p3 = qf_base_change(p3,u,1); |
|
if (test_mat(p3,p,C2,k)) |
|
{ |
{ |
*kmax = k; |
pol = primpart(pol); |
return gerepileupto(av,gmul(p1,u)); |
lc = absi(leading_term(pol)); |
} |
} |
|
lcpol = lc? gmul(lc, pol): pol; |
/* we also need to increase the precision */ |
list[i] = (long)QXQ_normalize(y, nfT); |
p2=shifti(gceil(mulsr(k, glog(p,DEFAULTPREC))),-1); |
|
prec += (long)(itos(p2)*pariK1); |
|
if (DEBUGLEVEL > 1) err(warnprec,"nffactor[2]",prec); |
|
k = k<<1; p1 = idealmullll(nf,p1,p1); |
|
if (low_stack(lim, stack_lim(av1,2))) |
|
{ |
|
if (DEBUGMEM>1) err(warnmem,"T2_matrix_pow"); |
|
p1 = gerepileupto(av1,p1); |
|
} |
|
if (!tot_real) T2 = nf_get_T2((GEN) nf[7], roots(x,prec)); |
|
p3 = qf_base_change(T2,p1,1); |
|
} |
} |
|
y = primpart(pol); |
|
list[i] = (long)QXQ_normalize(y, nfT); return list; |
} |
} |
|
|
/* return the factorization of the square-free polynomial x. |
|
The coeff of x are in Z_nf and its leading term is a rational |
|
integer. If fl = 1,return only the roots of x in nf */ |
|
static GEN |
static GEN |
nfsqff(GEN nf,GEN pol, long fl) |
nf_to_Zp(GEN x, GEN pk, GEN pas2, GEN ffproj) |
{ |
{ |
long d=lgef(pol),i,k,m,n,av=avma,tetpil,newprec,prec,nbf=BIGINT,anbf,ct=5; |
return centermodii(gmul(ffproj, x), pk, pas2); |
GEN p1,pr,p2,rep,k2,C,h,dk,dki,p,prh,p3,T2,polbase,fact,pk,ap,apr; |
} |
GEN polmod,polred,hinv,lt,minp,den,maxp=shifti(gun,32),run,aprh; |
|
|
|
if (DEBUGLEVEL>=4) msgtimer("square-free"); |
/* assume lc(pol) != 0 mod prh. Reduce P to Zp[X] mod p^a */ |
|
static GEN |
|
ZpX(GEN P, GEN pk, GEN ffproj) |
|
{ |
|
long i, l; |
|
GEN z, pas2 = shifti(pk,-1); |
|
|
dk=absi((GEN)nf[3]); |
if (typ(P)!=t_POL) return nf_to_Zp(P,pk,pas2,ffproj); |
dki=mulii(dk,(GEN)nf[4]); |
l = lgef(P); |
n=degpol(nf[1]); |
z = cgetg(l,t_POL); z[1] = P[1]; |
|
for (i=2; i<l; i++) z[i] = (long)nf_to_Zp((GEN)P[i],pk,pas2,ffproj); |
|
return normalizepol(z); |
|
} |
|
|
polbase = unifpol(nf,pol,0); |
/* We want to be able to reconstruct x, |x|^2 < C, from x mod pr^k |
polmod = unifpol(nf,pol,1); |
* We want B := min B_i >= C. Let alpha the LLL parameter, |
dki=mulii(dki,gnorm((GEN)polmod[d-1])); |
* y = 1/(alpha-1/4) > 4/3, the theoretical guaranteed bound follows from |
|
* (Npr)^(2k/n) = det(L)^(2/n) <= 1/n sum B_i |
|
* <= 1/n B sum y^i |
|
* <= 1/n B y^(n-1) / (y-1) |
|
* |
|
* Hence log(B/n) >= 2k/n log(Npr) - (n-1)log(y) + log(y-1) |
|
* >= log(C/n), provided |
|
* k >= [ log(C/n) + (n-1)log(y) - log(y-1) ] * n/(2log(Npr)) |
|
*/ |
|
static double |
|
bestlift_bound(GEN C, long n, double alpha, GEN Npr) |
|
{ |
|
const double y = 1 / (alpha - 0.25); /* = 2 if alpha = 3/4 */ |
|
double t; |
|
if (typ(C) != t_REAL) C = gmul(C, realun(DEFAULTPREC)); |
|
setlg(C, DEFAULTPREC); |
|
t = rtodbl(mplog(divrs(C,n))) + (n-1) * log(y) - log(y-1); |
|
return ceil((t * n) / (2.* log(gtodouble(Npr)))); |
|
} |
|
|
prec = DEFAULTPREC; |
static void |
for (;;) |
bestlift_init(long a, GEN nf, GEN pr, GEN C, nflift_t *T) |
{ |
{ |
if (prec <= gprecision(nf)) |
const int D = 4; |
T2 = gprec_w(gmael(nf,5,3), prec); |
const double alpha = ((double)D-1) / D; /* LLL parameter */ |
else |
gpmem_t av = avma; |
T2 = nf_get_T2((GEN)nf[7], roots((GEN)nf[1], prec)); |
GEN prk, PRK, B, GSmin, pk; |
|
|
run=realun(prec); |
if (!a) a = (long)bestlift_bound(C, degpol(nf[1]), alpha, idealnorm(nf,pr)); |
p1=realzero(prec); |
|
for (i=2; i<d; i++) |
|
{ |
|
p2 = gmul(run, (GEN)polbase[i]); |
|
p2 = qfeval(T2, p2); |
|
p1 = addrr(p1, gdiv(p2, binome(stoi(d), i-2))); |
|
if (signe(p1) < 0) break; |
|
} |
|
|
|
if (signe(p1) > 0) break; |
for (;; avma = av, a<<=1) |
|
{ |
|
if (DEBUGLEVEL>2) fprintferr("exponent: %ld\n",a); |
|
prk = idealpows(nf, pr, a); |
|
pk = gcoeff(prk,1,1); |
|
PRK = hnfmodid(prk, pk); |
|
|
prec = (prec<<1)-2; |
PRK = lllint_i(PRK, D, 0, NULL, NULL, &B); |
if (DEBUGLEVEL > 1) err(warnprec, "nfsqff", prec); |
GSmin = vecmin(GS_norms(B, DEFAULTPREC)); |
|
if (gcmp(GSmin, C) >= 0) break; |
} |
} |
|
if (DEBUGLEVEL>2) fprintferr("for this exponent, GSmin = %Z\n",GSmin); |
|
T->k = a; |
|
T->pk = T->den = pk; |
|
T->prk = PRK; |
|
T->iprk = ZM_inv(PRK, pk); |
|
T->GSmin= GSmin; |
|
T->ZpProj = dim1proj(prk); /* nf --> Zp */ |
|
} |
|
|
lt = (GEN)leading_term(polbase)[1]; |
/* assume pr has degree 1 */ |
p1 = gmul(p1, mulis(sqri(lt), n)); |
static GEN |
C = gpow(stoi(3), gadd(gmulsg(3, ghalf), stoi(d)), prec); |
nf_LLL_cmbf(nfcmbf_t *T, GEN p, long k, long rec) |
C = gdiv(gmul(C, p1), gmulsg(d, mppi(prec))); |
{ |
|
nflift_t *L = T->L; |
if (DEBUGLEVEL>=4) |
GEN pk = L->pk, PRK = L->prk, PRKinv = L->iprk, GSmin = L->GSmin; |
fprintferr("the bound on the T2-norm of the coeff. is: %Z\n", C); |
|
|
|
/* the theoretical bound for the exponent should be: |
GEN famod = T->fact, nf = T->nf, ZC = T->ZC, Br = T->Br; |
k2=gadd(glog(gdivgs(C,n),DEFAULTPREC), mulsr(n*(n-1), dbltor(0.347))); */ |
GEN Pbase = T->polbase, P = T->pol, dn = T->dn; |
k2=gadd(glog(gdivgs(C,n),DEFAULTPREC), mulsr(n*(n-1), dbltor(0.15))); |
GEN nfT = (GEN)nf[1]; |
k2=gmul2n(gmulgs(k2,n),-1); |
GEN Btra; |
|
long dnf = degpol(nfT), dP = degpol(P); |
|
|
minp=gmin(gexp(gmul2n(k2,-6),BIGDEFAULTPREC), maxp); |
double BitPerFactor = 0.5; /* nb bits / modular factor */ |
minp=gceil(minp); |
long i, C, tmax, n0; |
|
GEN lP, Bnorm, Tra, T2, TT, CM_L, m, list, ZERO; |
if (DEBUGLEVEL>=4) |
double Blow; |
{ |
gpmem_t av, av2, lim; |
fprintferr("lower bound for the prime numbers: %Z\n", minp); |
long ti_LLL = 0, ti_CF = 0; |
msgtimer("bounds computation"); |
pari_timer ti, ti2, TI; |
} |
|
|
|
p = rep = polred = NULL; /* gcc -Wall */ |
if (DEBUGLEVEL>2) (void)TIMER(&ti); |
pr=NULL; |
|
for (;;) |
lP = absi(leading_term(P)); |
|
if (is_pm1(lP)) lP = NULL; |
|
|
|
n0 = lg(famod) - 1; |
|
/* Lattice: (S PRK), small vector (vS vP). To find k bound for the image, |
|
* write S = S1 q + S0, P = P1 q + P0 |
|
* |S1 vS + P1 vP|^2 <= Blow for all (vS,vP) assoc. to true factors */ |
|
Btra = mulrr(ZC, mulsr(dP*dP, normlp(Br, 2, dnf))); |
|
Blow = get_Blow(n0, dnf, gceil(Btra)); |
|
C = (long)ceil(sqrt(Blow/n0)) + 1; /* C^2 n0 ~ Blow */ |
|
Bnorm = dbltor( n0 * C * C + Blow ); |
|
ZERO = zeromat(n0, dnf); |
|
|
|
av = avma; lim = stack_lim(av, 1); |
|
TT = cgetg(n0+1, t_VEC); |
|
Tra = cgetg(n0+1, t_MAT); |
|
for (i=1; i<=n0; i++) TT[i] = 0; |
|
CM_L = gscalsmat(C, n0); |
|
/* tmax = current number of traces used (and computed so far) */ |
|
for(tmax = 0;; tmax++) |
{ |
{ |
apr = choose_prime(nf,dki,minp, pr?30:0); |
long a, b, bmin, bgood, delta, tnew = tmax + 1, r = lg(CM_L)-1; |
if (!apr) break; |
GEN oldCM_L, M_L, q, S1, P1, VV; |
|
int first = 1; |
|
|
ap=(GEN)apr[1]; |
/* bound for f . S_k(genuine factor) = ZC * bound for T_2(S_tnew) */ |
aprh=prime_to_ideal(nf,apr); |
Btra = mulrr(ZC, mulsr(dP*dP, normlp(Br, 2*tnew, dnf))); |
|
bmin = logint(ceil_safe(mpsqrt(Btra)), gdeux, NULL); |
|
if (DEBUGLEVEL>2) |
|
fprintferr("\nLLL_cmbf: %ld potential factors (tmax = %ld, bmin = %ld)\n", |
|
r, tmax, bmin); |
|
|
polred=gcopy(polbase); |
/* compute Newton sums (possibly relifting first) */ |
lt=(GEN)leading_term(polbase)[1]; |
if (gcmp(GSmin, Btra) < 0) |
lt=mpinvmod(lt,ap); |
{ |
|
nflift_t L1; |
|
GEN polred; |
|
|
for (i=2; i<d; i++) |
bestlift_init(k<<1, nf, T->pr, Btra, &L1); |
polred[i] = nfreducemodpr_i(gmul(lt,(GEN)polbase[i]), aprh)[1]; |
k = L1.k; |
|
pk = L1.pk; |
|
PRK = L1.prk; |
|
PRKinv = L1.iprk; |
|
GSmin = L1.GSmin; |
|
|
if (!divise(discsr(polred),ap)) |
polred = ZpX(Pbase, pk, L1.ZpProj); |
|
famod = hensel_lift_fact(polred,famod,NULL,p,pk,k); |
|
for (i=1; i<=n0; i++) TT[i] = 0; |
|
} |
|
for (i=1; i<=n0; i++) |
{ |
{ |
rep=(GEN)simplefactmod(polred,ap)[1]; |
GEN p1, lPpow; |
anbf=lg(rep)-1; |
GEN p2 = polsym_gen((GEN)famod[i], (GEN)TT[i], tnew, NULL, pk); |
ct--; |
TT[i] = (long)p2; |
if (anbf < nbf) |
p1 = (GEN)p2[tnew+1]; |
|
/* make Newton sums integral */ |
|
if (lP) |
{ |
{ |
nbf=anbf; |
lPpow = gpowgs(lP, tnew); |
pr=gcopy(apr); |
p1 = mulii(p1, lPpow); /* assume pr has degree 1 */ |
p=gcopy(ap); |
|
if (DEBUGLEVEL>=4) |
|
{ |
|
fprintferr("prime ideal considered: %Z\n", pr); |
|
fprintferr("number of irreducible factors: %ld\n", nbf); |
|
} |
|
if (nbf == 1) break; |
|
} |
} |
|
if (dn) p1 = mulii(p1,dn); |
|
if (dn || lP) p1 = modii(p1, pk); |
|
Tra[i] = (long)nf_bestlift(p1, NULL, L); /* S_tnew(famod) */ |
} |
} |
if (pr && !ct) break; |
|
|
|
minp = addis(ap,1); |
/* compute truncation parameter */ |
} |
if (DEBUGLEVEL>2) { TIMERstart(&ti2); TIMERstart(&TI); } |
|
oldCM_L = CM_L; |
|
av2 = avma; |
|
b = delta = 0; /* -Wall */ |
|
AGAIN: |
|
M_L = gdivexact(CM_L, stoi(C)); |
|
T2 = gmul(Tra, M_L); |
|
VV = gdivround(gmul(PRKinv, T2), pk); |
|
T2 = gsub(T2, gmul(PRK, VV)); |
|
|
k = itos(gceil(gdiv(k2,glog(p,BIGDEFAULTPREC)))); |
if (first) |
|
{ /* initialize lattice, using few p-adic digits for traces */ |
|
a = gexpo(T2); |
|
bgood = (long)(a - max(32, BitPerFactor * r)); |
|
b = max(bmin, bgood); |
|
delta = a - b; |
|
} |
|
else |
|
{ /* add more p-adic digits and continue reduction */ |
|
long b0 = gexpo(T2); |
|
if (b0 < b) b = b0; |
|
b = max(b-delta, bmin); |
|
if (b - delta/2 < bmin) b = bmin; /* near there. Go all the way */ |
|
} |
|
|
if (DEBUGLEVEL>=4) |
/* restart with truncated entries */ |
{ |
q = shifti(gun, b); |
fprintferr("prime ideal chosen: %Z\n", pr); |
P1 = gdivround(PRK, q); |
msgtimer("choice of the prime ideal"); |
S1 = gdivround(Tra, q); |
} |
T2 = gsub(gmul(S1, M_L), gmul(P1, VV)); |
|
m = vconcat( CM_L, T2 ); |
|
if (first) |
|
{ |
|
first = 0; |
|
m = concatsp( m, vconcat(ZERO, P1) ); |
|
/* [ C M_L 0 ] |
|
* m = [ ] square matrix |
|
* [ T2' PRK ] T2' = Tra * M_L truncated |
|
*/ |
|
} |
|
CM_L = LLL_check_progress(Bnorm, n0, m, b == bmin, /*dbg:*/ &ti, &ti_LLL); |
|
if (DEBUGLEVEL>2) |
|
fprintferr("LLL_cmbf: b =%4ld; r =%3ld -->%3ld, time = %ld\n", |
|
b, lg(m)-1, CM_L? lg(CM_L)-1: 1, TIMER(&TI)); |
|
if (!CM_L) { list = _col(QXQ_normalize(P,nfT)); break; } |
|
if (b > bmin) |
|
{ |
|
CM_L = gerepilecopy(av2, CM_L); |
|
goto AGAIN; |
|
} |
|
if (DEBUGLEVEL>2) msgTIMER(&ti2, "for this trace"); |
|
|
if (lg(rep)==2) |
i = lg(CM_L) - 1; |
{ |
if (i == r && gegal(CM_L, oldCM_L)) |
if (fl) { avma=av; return cgetg(1,t_VEC); } |
{ |
rep=cgetg(2,t_VEC); rep[1]=lcopy(polmod); |
CM_L = oldCM_L; |
nfcmbf.nfact=1; return gerepileupto(av, rep); |
avma = av2; continue; |
|
} |
|
|
|
if (i <= r && i*rec < n0) |
|
{ |
|
if (DEBUGLEVEL>2) (void)TIMER(&ti); |
|
list = nf_check_factors(T, P, M_L, famod, pk); |
|
if (DEBUGLEVEL>2) ti_CF += TIMER(&ti); |
|
if (list) break; |
|
CM_L = gerepilecopy(av2, CM_L); |
|
} |
|
if (low_stack(lim, stack_lim(av,1))) |
|
{ |
|
if(DEBUGMEM>1) err(warnmem,"nf_LLL_cmbf"); |
|
gerepileall(av, 8, &CM_L,&TT,&Tra,&famod,&pk,&GSmin,&PRK,&PRKinv); |
|
} |
} |
} |
|
if (DEBUGLEVEL>2) |
|
fprintferr("* Time LLL: %ld\n* Time Check Factor: %ld\n",ti_LLL,ti_CF); |
|
return list; |
|
} |
|
|
p2=cgetr(DEFAULTPREC); |
static GEN |
affir(mulii(absi(dk),gpowgs(p,k)),p2); |
nf_combine_factors(nfcmbf_t *T, GEN polred, GEN p, long a, long klim) |
p2=shifti(gceil(mplog(p2)),-1); |
{ |
|
GEN z, res, L, listmod, famod = T->fact, nf = T->nf; |
|
long i, m, l, maxK = 3, nft = lg(famod)-1; |
|
|
newprec = MEDDEFAULTPREC + (long)(itos(p2)*pariK1); |
T->fact = hensel_lift_fact(polred,famod,NULL,p,T->pk,a); |
if (DEBUGLEVEL>=4) |
if (nft < 11) maxK = -1; /* few modular factors: try all posibilities */ |
fprintferr("new precision: %ld\n",newprec); |
if (DEBUGLEVEL>3) msgtimer("Hensel lift"); |
|
|
prh = idealpows(nf,pr,k); m = k; |
/* FIXME: neither nfcmbf nor LLL_cmbf can handle the non-nf case */ |
h = T2_matrix_pow(nf,prh,p,C, &k, newprec); |
|
if (m != k) prh=idealpows(nf,pr,k); |
|
|
|
if (DEBUGLEVEL>=4) |
T->res = cgetg(nft+1,t_VEC); |
|
L = nfcmbf(T, p, a, maxK, klim); |
|
|
|
res = (GEN)L[1]; |
|
listmod = (GEN)L[2]; l = lg(listmod)-1; |
|
famod = (GEN)listmod[l]; |
|
if (maxK >= 0 && lg(famod)-1 > 2*maxK) |
{ |
{ |
fprintferr("a suitable exponent is: %ld\n",(long)k); |
if (l > 1) T->polbase = unifpol(nf, (GEN)res[l], 0); |
msgtimer("computation of H"); |
L = nf_LLL_cmbf(T, p, a, maxK); |
|
/* remove last elt, possibly unfactored. Add all new ones. */ |
|
setlg(res, l); res = concatsp(res, L); |
} |
} |
|
if (DEBUGLEVEL>3) msgtimer("computation of the factors"); |
|
|
pk = gcoeff(prh,1,1); |
m = lg(res); z = cgetg(m, t_VEC); |
lt=(GEN)leading_term(polbase)[1]; |
for (i=1;i<m; i++) z[i] = (long)unifpol(nf,(GEN)res[i],1); |
lt=mpinvmod(lt,pk); |
return z; |
|
} |
|
|
polred[1] = polbase[1]; |
/* return the factorization of the square-free polynomial x. |
for (i=2; i<d; i++) |
The coeff of x are in Z_nf and its leading term is a rational integer. |
polred[i] = nfreducemodpr_i(gmul(lt,(GEN)polbase[i]), prh)[1]; |
deg(x) > 1, deg(nfpol) > 1 |
|
If fl = 1, return only the roots of x in nf |
|
If fl = 2, as fl=1 if pol splits, [] otherwise */ |
|
static GEN |
|
nfsqff(GEN nf, GEN pol, long fl) |
|
{ |
|
long i, m, n, nbf, ct, dpol = degpol(pol); |
|
gpmem_t av = avma; |
|
GEN pr, C0, C, dk, bad, polbase, pk; |
|
GEN N2, rep, p, ap, polmod, polred, lt, nfpol; |
|
byteptr pt=diffptr; |
|
nfcmbf_t T; |
|
nflift_t L; |
|
|
fact = lift_intern((GEN)factmod(polred,p)[1]); |
if (DEBUGLEVEL>3) msgtimer("square-free"); |
rep = hensel_lift_fact(polred,fact,NULL,p,pk,k); |
nfpol = (GEN)nf[1]; n = degpol(nfpol); |
|
polbase = unifpol(nf,pol,0); |
|
polmod = unifpol(nf,pol,1); |
|
/* heuristic */ |
|
#if 1 |
|
if (dpol*4 < n) |
|
{ |
|
if (DEBUGLEVEL>2) fprintferr("Using Trager's method\n"); |
|
return gerepilecopy(av, (GEN)polfnf(polmod, nfpol)[1]); |
|
} |
|
#endif |
|
|
if (DEBUGLEVEL >= 4) msgtimer("computation of the p-adic factorization"); |
pol = simplify_i(lift(polmod)); |
|
lt = (GEN)leading_term(polbase)[1]; /* t_INT */ |
|
|
den = det(h); /* |den| = NP^k */ |
dk = absi((GEN)nf[3]); |
hinv= adj(h); |
bad = mulii(mulii(dk,(GEN)nf[4]), lt); |
lt=(GEN)leading_term(polbase)[1]; |
|
|
|
if (fl) |
p = polred = pr = NULL; /* gcc -Wall */ |
|
nbf = 0; ap = NULL; |
|
for (ct = 5;;) |
{ |
{ |
long x_a[] = { evaltyp(t_POL)|m_evallg(4), 0,0,0 }; |
GEN apr = choose_prime(nf, bad, &ap, &pt); |
GEN mlt = gneg_i(lt), rem; |
GEN modpr = zkmodprinit(nf, apr); |
x_a[1] = polbase[1]; setlgef(x_a, 4); |
long anbf; |
x_a[3] = un; |
|
p1=cgetg(lg(rep)+1,t_VEC); |
|
polbase = unifpol(nf,polbase,1); |
|
for (m=1,i=1; i<lg(rep); i++) |
|
{ |
|
p2=(GEN)rep[i]; if (lgef(p2)!=4) break; |
|
|
|
p3 = algtobasis(nf, gmul(mlt,(GEN)p2[2])); |
polred = modprX(polbase, nf, modpr); |
p3 = nf_bestlift(h,hinv,den,p3); |
if (!FpX_is_squarefree(polred,ap)) continue; |
p3 = gdiv(p3,lt); |
|
x_a[2] = lneg(p3); /* check P(p3) == 0 */ |
anbf = fl? FpX_nbroots(polred,ap): FpX_nbfact(polred,ap); |
p2 = poldivres(polbase, unifpol(nf,x_a,1), &rem); |
if (!nbf || anbf < nbf) |
if (!signe(rem)) { polbase = p2; p1[m++] = (long)p3; } |
{ |
|
nbf = anbf; pr = apr; p = ap; |
|
if (DEBUGLEVEL>3) |
|
fprintferr("%3ld %s at prime ideal above %Z\n", |
|
nbf, fl?"roots": "factors", p); |
|
if (fl == 2 && nbf < dpol) return cgetg(1,t_VEC); |
|
if (nbf <= 1) |
|
{ |
|
if (!fl) /* irreducible */ |
|
return gerepilecopy(av, _vec(QXQ_normalize(polmod, nfpol))); |
|
if (!nbf) return cgetg(1,t_VEC); /* no root */ |
|
} |
} |
} |
tetpil=avma; rep=cgetg(m,t_VEC); |
if (--ct <= 0) break; |
for (i=1; i<m; i++) rep[i]=(long)basistoalg(nf,(GEN)p1[i]); |
|
return gerepile(av,tetpil,rep); |
|
} |
} |
|
if (DEBUGLEVEL>3) { |
|
msgtimer("choice of a prime ideal"); |
|
fprintferr("Prime ideal chosen: %Z\n", pr); |
|
} |
|
|
for (i=1; i<lg(rep); i++) |
L.tozk = (GEN)nf[8]; |
rep[i] = (long)unifpol(nf,(GEN)rep[i],0); |
L.topow= (GEN)nf[7]; |
|
T.ZC = L2_bound(nf, L.tozk, &(T.dn)); |
|
T.Br = gmul(lt, nf_root_bounds(pol, nf)); |
|
|
nfcmbf.pol = polmod; |
if (fl) C0 = normlp(T.Br, 2, n); |
nfcmbf.lt = lt; |
else C0 = nf_factor_bound(nf, polbase); /* bound for T_2(Q_i), Q | P */ |
nfcmbf.h = h; |
C = mulrr(T.ZC, C0); /* bound for |Q_i|^2 in Z^n on chosen Z-basis */ |
nfcmbf.hinv = hinv; |
T.bound = C; |
nfcmbf.den = den; |
|
nfcmbf.fact = rep; |
|
nfcmbf.res = cgetg(lg(rep)+1,t_VEC); |
|
nfcmbf.nfact = 0; |
|
nfcmbf.nfactmod = lg(rep)-1; |
|
nf_combine_factors(nf,1,NULL,d-3,1); |
|
|
|
if (DEBUGLEVEL >= 4) msgtimer("computation of the factors"); |
N2 = mulsr(dpol*dpol, normlp(T.Br, 4, n)); /* bound for T_2(lt * S_2) */ |
|
T.BS_2 = mulrr(T.ZC, N2); /* bound for |S_2|^2 on chosen Z-basis */ |
|
|
i = nfcmbf.nfact; |
if (DEBUGLEVEL>2) { |
if (lgef(nfcmbf.pol)>3) |
msgtimer("bound computation"); |
{ |
fprintferr(" 1) T_2 bound for %s: %Z\n", fl?"root":"factor", C0); |
nfcmbf.res[++i] = (long) nf_pol_divres(nf,nfcmbf.pol,nfcmbf.lt,NULL); |
fprintferr(" 2) Conversion from T_2 --> | |^2 bound : %Z\n", T.ZC); |
nfcmbf.nfact = i; |
fprintferr(" 3) Final bound: %Z\n", C); |
} |
} |
|
|
tetpil=avma; rep=cgetg(i+1,t_VEC); |
if (fl) |
for ( ; i>=1; i--) |
(void)logint(C, p, &pk); |
rep[i]=(long)unifpol(nf,(GEN)nfcmbf.res[i],1); |
#if 0 /* overkill */ |
return gerepile(av,tetpil,rep); |
else |
} |
{ /* overlift needed for d-1/d-2 tests? */ |
|
GEN pb; long b; /* junk */ |
static int |
if (cmbf_precs(p, C, T.BS_2, &a, &b, &pk, &pb)) |
nf_combine_factors(GEN nf,long fxn,GEN psf,long dlim,long hint) |
{ /* Rare */ |
{ |
if (DEBUGLEVEL) err(warner,"nffactor: overlift for d-1/d-2 test"); |
int val = 0; /* assume failure */ |
C = itor(pk, DEFAULTPREC); |
GEN newf, newpsf = NULL; |
} |
long av,newd,ltop,i; |
|
|
|
/* Assertion: fxn <= nfcmbf.nfactmod && dlim > 0 */ |
|
|
|
/* first, try deeper factors without considering the current one */ |
|
if (fxn != nfcmbf.nfactmod) |
|
{ |
|
val = nf_combine_factors(nf,fxn+1,psf,dlim,hint); |
|
if (val && psf) return 1; |
|
} |
} |
|
#endif |
|
|
/* second, try including the current modular factor in the product */ |
bestlift_init(0, nf, pr, C, &L); |
newf = (GEN)nfcmbf.fact[fxn]; |
T.pr = pr; |
if (!newf) return val; /* modular factor already used */ |
T.L = &L; |
newd = degpol(newf); |
|
if (newd > dlim) return val; /* degree of new factor is too large */ |
|
|
|
av = avma; |
/* polred is monic */ |
if (newd % hint == 0) |
pk = L.pk; |
{ |
polred = ZpX(gmul(mpinvmod(lt,pk), polbase), pk, L.ZpProj); |
GEN p, quot,rem; |
|
|
|
newpsf = nf_pol_mul(nf, psf? psf: nfcmbf.lt, newf); |
if (fl) |
newpsf = nf_pol_lift(nfcmbf.h,nfcmbf.hinv,nfcmbf.den,newpsf); |
{ /* roots only */ |
/* try out the new combination */ |
long x_r[] = { evaltyp(t_POL)|_evallg(4), 0,0,0 }; |
ltop=avma; |
rep = rootpadicfast(polred, p, L.k); |
quot=nf_pol_divres(nf,nfcmbf.pol,newpsf,&rem); |
x_r[1] = evalsigne(1) | evalvarn(varn(pol)) | evallgef(4); |
if (gcmp0(rem)) /* found a factor */ |
x_r[3] = un; |
|
for (m=1,i=1; i<lg(rep); i++) |
{ |
{ |
p = nf_pol_mul(nf,element_inv(nf,leading_term(newpsf)),newpsf); |
GEN q, r = (GEN)rep[i]; |
nfcmbf.res[++nfcmbf.nfact] = (long) p; /* store factor */ |
|
nfcmbf.fact[fxn]=0; /* remove used modular factor */ |
|
|
|
/* fix up target */ |
r = nf_bestlift_to_pol(gmul(lt,r), NULL, &L); |
p=gun; quot=unifpol(nf,quot,0); |
r = gdiv(r,lt); |
for (i=2; i<lgef(quot); i++) |
x_r[2] = lneg(r); /* check P(r) == 0 */ |
p = mpppcm(p, denom((GEN)quot[i])); |
q = RXQX_divrem(pol, x_r, nfpol, ONLY_DIVIDES); |
|
if (q) { pol = q; rep[m++] = (long)r; } |
nfcmbf.pol = nf_pol_mul(nf,p,quot); |
else if (fl == 2) return cgetg(1, t_VEC); |
nfcmbf.lt = leading_term(nfcmbf.pol); |
|
return 1; |
|
} |
} |
avma=ltop; |
rep[0] = evaltyp(t_VEC) | evallg(m); |
|
return gerepilecopy(av, rep); |
} |
} |
|
|
/* If room in degree limit + more modular factors to try, add more factors to |
T.polbase = polbase; |
* newpsf */ |
T.pol = pol; |
if (newd < dlim && fxn < nfcmbf.nfactmod |
T.nf = nf; |
&& nf_combine_factors(nf,fxn+1,newpsf,dlim-newd,hint)) |
T.fact = (GEN)factmod0(polred,p)[1]; |
{ |
T.hint = 1; /* useless */ |
nfcmbf.fact[fxn]=0; /* remove used modular factor */ |
|
return 1; |
rep = nf_combine_factors(&T, polred, p, L.k, dpol-1); |
} |
return gerepileupto(av, rep); |
avma = av; return val; |
|
} |
} |
|
|
/* return the characteristic polynomial of alpha over nf, where alpha |
/* return the characteristic polynomial of alpha over nf, where alpha |
is an element of the algebra nf[X]/(T) given as a polynomial in X */ |
is an element of the algebra nf[X]/(T) given as a polynomial in X */ |
GEN |
GEN |
rnfcharpoly(GEN nf,GEN T,GEN alpha,int v) |
rnfcharpoly(GEN nf,GEN T,GEN alpha,int v) |
{ |
{ |
long av = avma, vnf, vT, lT; |
long vnf, vT, lT; |
|
gpmem_t av = avma; |
GEN p1; |
GEN p1; |
|
|
nf=checknf(nf); vnf = varn(nf[1]); |
nf=checknf(nf); vnf = varn(nf[1]); |
Line 1243 rnfcharpoly(GEN nf,GEN T,GEN alpha,int v) |
|
Line 1414 rnfcharpoly(GEN nf,GEN T,GEN alpha,int v) |
|
return gerepileupto(av, gsub(polx[v], alpha)); |
return gerepileupto(av, gsub(polx[v], alpha)); |
p1 = caract2(unifpol(nf,T,1), unifpol(nf,alpha,1), v); |
p1 = caract2(unifpol(nf,T,1), unifpol(nf,alpha,1), v); |
return gerepileupto(av, unifpol(nf,p1,1)); |
return gerepileupto(av, unifpol(nf,p1,1)); |
} |
|
|
|
#if 0 |
|
/* return the minimal polynomial of alpha over nf, where alpha is |
|
an element of the algebra nf[X]/(T) given as a polynomial in X */ |
|
GEN |
|
rnfminpoly(GEN nf,GEN T,GEN alpha,int n) |
|
{ |
|
long av=avma,tetpil; |
|
GEN p1,p2; |
|
|
|
nf=checknf(nf); p1=rnfcharpoly(nf,T,alpha,n); |
|
tetpil=avma; p2=nf_pol_subres(nf,p1,deriv(p1,varn(T))); |
|
if (lgef(p2)==3) { avma=tetpil; return p1; } |
|
|
|
p1 = nf_pol_divres(nf,p1,p2,NULL); |
|
p2 = element_inv(nf,leading_term(p1)); |
|
tetpil=avma; return gerepile(av,tetpil,unifpol(nf,nf_pol_mul(nf,p2,p1),1)); |
|
} |
|
#endif |
|
|
|
/* relative Dedekind criterion over nf, applied to the order defined by a |
|
* root of irreducible polynomial T, modulo the prime ideal pr. Returns |
|
* [flag,basis,val], where basis is a pseudo-basis of the enlarged order, |
|
* flag is 1 iff this order is pr-maximal, and val is the valuation in pr of |
|
* the order discriminant |
|
*/ |
|
GEN |
|
rnfdedekind(GEN nf,GEN T,GEN pr) |
|
{ |
|
long av=avma,vt,r,d,da,n,m,i,j; |
|
GEN p1,p2,p,tau,g,vecun,veczero,matid; |
|
GEN prhall,res,h,k,base,Ca; |
|
|
|
nf=checknf(nf); Ca=unifpol(nf,T,0); |
|
res=cgetg(4,t_VEC); |
|
if (typ(pr)==t_VEC && lg(pr)==3) |
|
{ prhall = (GEN)pr[2]; pr = (GEN)pr[1]; } |
|
else |
|
prhall=nfmodprinit(nf,pr); |
|
p=(GEN)pr[1]; tau=gdiv((GEN)pr[5],p); |
|
n=degpol(nf[1]); m=degpol(T); |
|
|
|
vecun = gscalcol_i(gun,n); |
|
veczero = zerocol(n); |
|
|
|
p1=(GEN)nffactormod(nf,Ca,pr)[1]; |
|
r=lg(p1); if (r < 2) err(constpoler,"rnfdedekind"); |
|
g=lift((GEN)p1[1]); |
|
for (i=2; i<r; i++) |
|
g = nf_pol_mul(nf,g, lift((GEN)p1[i])); |
|
h=nfmod_pol_divres(nf,prhall,Ca,g,NULL); |
|
k=nf_pol_mul(nf,tau,gsub(Ca, nf_pol_mul(nf,lift(g),lift(h)))); |
|
p2=nfmod_pol_gcd(nf,prhall,g,h); |
|
k= nfmod_pol_gcd(nf,prhall,p2,k); |
|
|
|
d=degpol(k); /* <= m */ |
|
vt = idealval(nf,discsr(T),pr) - 2*d; |
|
res[3]=lstoi(vt); |
|
if (!d || vt<=1) res[1]=un; else res[1]=zero; |
|
|
|
base=cgetg(3,t_VEC); |
|
p1=cgetg(m+d+1,t_MAT); base[1]=(long)p1; |
|
p2=cgetg(m+d+1,t_VEC); base[2]=(long)p2; |
|
/* if d > 0, base[2] temporarily multiplied by p, for the final nfhermitemod |
|
* [ which requires integral ideals ] */ |
|
matid = gscalmat(d? p: gun, n); |
|
for (j=1; j<=m; j++) |
|
{ |
|
p2[j]=(long)matid; |
|
p1[j]=lgetg(m+1,t_COL); |
|
for (i=1; i<=m; i++) |
|
coeff(p1,i,j) = (i==j)?(long)vecun:(long)veczero; |
|
} |
|
if (d) |
|
{ |
|
GEN pal = lift(nfmod_pol_divres(nf,prhall,Ca,k,NULL)); |
|
GEN prinvp = pidealprimeinv(nf,pr); /* again multiplied by p */ |
|
GEN nfx = unifpol(nf,polx[varn(T)],0); |
|
|
|
for ( ; j<=m+d; j++) |
|
{ |
|
p1[j]=lgetg(m+1,t_COL); |
|
da=degpol(pal); |
|
for (i=1; i<=da+1; i++) coeff(p1,i,j)=pal[i+1]; |
|
for ( ; i<=m; i++) coeff(p1,i,j)=(long)veczero; |
|
p2[j]=(long)prinvp; |
|
nf_pol_divres(nf,nf_pol_mul(nf,pal,nfx),T,&pal); |
|
} |
|
/* the modulus is integral */ |
|
base = nfhermitemod(nf,base, gmul(gpowgs(p, m-d), |
|
idealpows(nf, prinvp, d))); |
|
base[2] = ldiv((GEN)base[2], p); /* cancel the factor p */ |
|
} |
|
res[2]=(long)base; return gerepilecopy(av, res); |
|
} |
} |