Annotation of OpenXM_contrib/pari/src/modules/nffactor.c, Revision 1.1
1.1 ! maekawa 1: /*******************************************************************/
! 2: /* */
! 3: /* Factorisation dans un corps de nombres */
! 4: /* et modulo un ideal premier */
! 5: /* */
! 6: /*******************************************************************/
! 7: /* $Id: nffactor.c,v 1.5 1999/09/23 18:48:42 xavier Exp $ */
! 8: #include "pari.h"
! 9:
! 10: GEN hensel_lift(GEN pol,GEN fk,GEN fkk,GEN p,long e);
! 11: GEN hensel_lift_fact(GEN pol, GEN fact, GEN p, GEN pev, long e);
! 12: GEN nf_get_T2(GEN base, GEN polr);
! 13: GEN nfreducemodpr_i(GEN x, GEN prh);
! 14: GEN sort_factor(GEN y, int (*cmp)(GEN,GEN));
! 15:
! 16: static GEN nf_pol_mul(GEN nf,GEN pol1,GEN pol2);
! 17: static GEN nf_pol_sqr(GEN nf,GEN pol1);
! 18: static GEN nf_pol_divres(GEN nf,GEN pol1,GEN pol2, GEN *pr);
! 19: static GEN nf_pol_subres(GEN nf,GEN pol1,GEN pol2);
! 20: static GEN nfmod_pol_reduce(GEN nf,GEN prhall,GEN pol);
! 21: static GEN nfmod_pol_mul(GEN nf,GEN prhall,GEN pol1,GEN pol2);
! 22: static GEN nfmod_pol_sqr(GEN nf,GEN prhall,GEN pol1);
! 23: static GEN nfmod_pol_divres(GEN nf,GEN prhall,GEN pol1,GEN pol2, GEN *pr);
! 24: static GEN nfmod_split2(GEN nf,GEN prhall,GEN pmod,GEN pol,GEN exp);
! 25: static GEN nfmod_pol_gcd(GEN nf,GEN prhall,GEN pol1,GEN pol2);
! 26: static GEN nfmod_pol_pow(GEN nf,GEN prhall,GEN pmod,GEN pol,GEN exp);
! 27: static GEN nf_bestlift(GEN id,GEN idinv,GEN den,GEN elt);
! 28: static GEN nf_pol_lift(GEN id,GEN idinv,GEN den,GEN pol);
! 29: static GEN nfsqff(GEN nf,GEN pol,long fl);
! 30: static int nf_combine_factors(GEN nf,long fxn,GEN psf,long dlim,long hint);
! 31:
! 32: typedef struct Nfcmbf /* for use in nfsqff */
! 33: {
! 34: GEN pol, h, hinv, fact, res, lt, den;
! 35: long nfact, nfactmod;
! 36: } Nfcmbf;
! 37: static Nfcmbf nfcmbf;
! 38:
! 39: static GEN
! 40: unifpol0(GEN nf,GEN pol,long flag)
! 41: {
! 42: static long n = 0;
! 43: static GEN vun = NULL;
! 44: GEN f = (GEN) nf[1];
! 45: long i = lgef(f)-3, av;
! 46:
! 47: if (i != n)
! 48: {
! 49: n=i; if (vun) free(vun);
! 50: vun = (GEN) gpmalloc((n+1)*sizeof(long));
! 51: vun[0] = evaltyp(t_COL) | evallg(n+1); vun[1]=un;
! 52: for ( ; i>=2; i--) vun[i]=zero;
! 53: }
! 54:
! 55: av = avma;
! 56: switch(typ(pol))
! 57: {
! 58: case t_INT: case t_FRAC: case t_RFRAC:
! 59: pol = gmul(pol,vun); break;
! 60:
! 61: case t_POL:
! 62: pol = gmodulcp(pol,f); /* fall through */
! 63: case t_POLMOD:
! 64: pol = algtobasis(nf,pol);
! 65: }
! 66: if (flag) pol = basistoalg(nf, lift(pol));
! 67: return gerepileupto(av,pol);
! 68: }
! 69:
! 70: /* Pour pol un polynome a coefficients dans Z, nf ou l'algebre correspondant
! 71: * renvoie le meme polynome avec pour coefficients:
! 72: * si flag=0: des vecteurs sur la base d'entiers.
! 73: * si flag=1: des polymods.
! 74: */
! 75: GEN
! 76: unifpol(GEN nf,GEN pol,long flag)
! 77: {
! 78: if (typ(pol)==t_POL && varn(pol) < varn(nf[1]))
! 79: {
! 80: long i, d = lgef(pol);
! 81: GEN p1 = pol;
! 82:
! 83: pol=cgetg(d,t_POL); pol[1]=p1[1];
! 84: for (i=2; i<d; i++)
! 85: pol[i] = (long) unifpol0(nf,(GEN) p1[i],flag);
! 86:
! 87: return pol;
! 88: }
! 89: return unifpol0(nf,(GEN) pol, flag);
! 90: }
! 91:
! 92: /* cree un polynome unitaire de degre d a entrees au hazard dans Z_nf */
! 93: GEN
! 94: random_pol(GEN nf,long d)
! 95: {
! 96: long i,j, n = lgef(nf[1])-3;
! 97: GEN pl,p;
! 98:
! 99: pl=cgetg(d+3,t_POL);
! 100: for (i=2; i<d+2; i++)
! 101: {
! 102: p=cgetg(n+1,t_COL); pl[i]=(long)p;
! 103: for (j=1; j<=n; j++)
! 104: p[j] = lstoi(mymyrand()%101 - 50);
! 105: }
! 106: p=cgetg(n+1,t_COL); pl[i]=(long)p;
! 107: p[1]=un; for (i=2; i<=n; i++) p[i]=zero;
! 108:
! 109: pl[1] = evalsigne(1) | evallgef(d+3) | evalvarn(0);
! 110: return pl;
! 111: }
! 112:
! 113: /* multiplication de x par y dans nf (x element accepte) */
! 114: static GEN
! 115: nf_pol_mul(GEN nf,GEN x,GEN y)
! 116: {
! 117: long tetpil,av=avma;
! 118: GEN res = gmul(unifpol(nf,x,1), unifpol(nf,y,1));
! 119:
! 120: tetpil = avma;
! 121: return gerepile(av,tetpil,unifpol(nf,res,0));
! 122: }
! 123:
! 124: /* calcule x^2 dans nf (x element accepte) */
! 125: static GEN
! 126: nf_pol_sqr(GEN nf,GEN x)
! 127: {
! 128: long tetpil,av=avma;
! 129: GEN res = gsqr(unifpol(nf,x,1));
! 130:
! 131: tetpil = avma;
! 132: return gerepile(av,tetpil,unifpol(nf,res,0));
! 133: }
! 134:
! 135: /* Reduction modulo phall des coefficients de pol (pol element accepte) */
! 136: static GEN
! 137: nfmod_pol_reduce(GEN nf,GEN prhall,GEN pol)
! 138: {
! 139: long av=avma,tetpil,i;
! 140: GEN p1;
! 141:
! 142: if (typ(pol)!=t_POL) return nfreducemodpr(nf,pol,prhall);
! 143: pol=unifpol(nf,pol,0);
! 144:
! 145: tetpil=avma; i=lgef(pol);
! 146: p1=cgetg(i,t_POL); p1[1]=pol[1];
! 147: for (i--; i>=2; i--)
! 148: p1[i] = (long) nfreducemodpr(nf,(GEN)pol[i],prhall);
! 149: return gerepile(av,tetpil, normalizepol(p1));
! 150: }
! 151:
! 152: /* x^2 modulo prhall ds nf (x elt accepte) */
! 153: static GEN
! 154: nfmod_pol_sqr(GEN nf,GEN prhall,GEN x)
! 155: {
! 156: long av=avma,tetpil;
! 157: GEN px;
! 158:
! 159: px = nfmod_pol_reduce(nf,prhall,x);
! 160: px = unifpol(nf,lift(px),1);
! 161: px = unifpol(nf,nf_pol_sqr(nf,px),0);
! 162: tetpil=avma;
! 163: return gerepile(av,tetpil,nfmod_pol_reduce(nf,prhall,px));
! 164: }
! 165:
! 166: /* multiplication des polynomes x et y modulo prhall ds nf (x elt accepte) */
! 167: static GEN
! 168: nfmod_pol_mul(GEN nf,GEN prhall,GEN x,GEN y)
! 169: {
! 170: long av=avma,tetpil;
! 171: GEN px,py;
! 172:
! 173: px = nfmod_pol_reduce(nf,prhall,x); px = unifpol(nf,lift(px),1);
! 174: py = nfmod_pol_reduce(nf,prhall,y); py = unifpol(nf,lift(py),1);
! 175: px = unifpol(nf,nf_pol_mul(nf,px,py),0);
! 176: tetpil=avma;
! 177: return gerepile(av,tetpil,nfmod_pol_reduce(nf,prhall,px));
! 178: }
! 179:
! 180: /* division euclidienne du polynome x par le polynome y */
! 181: static GEN
! 182: nf_pol_divres(GEN nf,GEN x,GEN y,GEN *pr)
! 183: {
! 184: long av = avma,tetpil;
! 185: GEN nq = poldivres(unifpol(nf,x,1),unifpol(nf,y,1),pr);
! 186: GEN *gptr[2];
! 187:
! 188: tetpil=avma; nq=unifpol(nf,nq,0);
! 189: if (pr) *pr = unifpol(nf,*pr,0);
! 190: gptr[0]=&nq; gptr[1]=pr;
! 191: gerepilemanysp(av,tetpil,gptr,pr ? 2:1);
! 192: return nq;
! 193: }
! 194:
! 195: /* division euclidienne du polynome x par le polynome y modulo prhall */
! 196: static GEN
! 197: nfmod_pol_divres(GEN nf,GEN prhall,GEN x,GEN y, GEN *pr)
! 198: {
! 199: long av=avma,dx,dy,dz,i,j,k,l,n,tetpil;
! 200: GEN z,p1,p2,p3,px,py;
! 201:
! 202: py = nfmod_pol_reduce(nf,prhall,y);
! 203: if (gcmp0(py))
! 204: err(talker, "division by zero in nfmod_pol_divres");
! 205:
! 206: tetpil=avma;
! 207: px=nfmod_pol_reduce(nf,prhall,x);
! 208: dx=lgef(px)-3; dy=lgef(py)-3; dz=dx-dy;
! 209: if (dx<dy)
! 210: {
! 211: GEN vzero;
! 212:
! 213: if (pr) *pr = gerepile(av,tetpil,px);
! 214: else avma = av;
! 215:
! 216: n=lgef(nf[1])-3;
! 217: vzero = cgetg(n+1,t_COL);
! 218: n=lgef(nf[1])-3;
! 219: for (i=1; i<=n; i++) vzero[i]=zero;
! 220:
! 221: z=cgetg(3,t_POL); z[2]=(long)vzero;
! 222: z[1]=evallgef(2) | evalvarn(varn(px));
! 223: return z;
! 224: }
! 225:
! 226: z=cgetg(dz+3,t_POL); z[1]=evalsigne(1) | evallgef(3+dz);
! 227: setvarn(z,varn(px));
! 228: z[dz+2] = (long) element_divmodpr(nf,(GEN)px[dx+2],(GEN)py[dy+2],prhall);
! 229: for (i=dx-1; i>=dy; --i)
! 230: {
! 231: l=avma; p1=nfreducemodpr(nf,(GEN)px[i+2],prhall);
! 232: for (j=i-dy+1; j<=i && j<=dz; j++)
! 233: {
! 234: p2=element_mul(nf,(GEN)z[j+2],(GEN)py[i-j+2]);
! 235: p2=nfreducemodpr(nf,p2,prhall); p1=gsub(p1,p2);
! 236: }
! 237: tetpil=avma; p3=element_divmodpr(nf,p1,(GEN)py[dy+2],prhall);
! 238: z[i-dy+2]=lpile(l,tetpil,p3);
! 239: z[i-dy+2]=(long)nfreducemodpr(nf,(GEN)z[i-dy+2],prhall);
! 240: }
! 241: l=avma;
! 242: for (i=dy-1; i>=0; --i)
! 243: {
! 244: l=avma; p1=((GEN)px[i+2]);
! 245: for (j=0; j<=i && j<=dz; j++)
! 246: {
! 247: p2=element_mul(nf,(GEN)z[j+2],(GEN)py[i-j+2]);
! 248: p2=nfreducemodpr(nf,p2,prhall); tetpil=avma; p1=gsub(p1,p2);
! 249: }
! 250: p1=gerepile(l,tetpil,p1);
! 251: if (!gcmp0(nfreducemodpr(nf,p1,prhall))) break;
! 252: }
! 253:
! 254: if (!pr) { avma = l; return z; }
! 255:
! 256: if (i<0)
! 257: {
! 258: avma=l;
! 259: p3 = cgetg(3,t_POL); p3[2]=zero;
! 260: p3[1] = evallgef(2) | evalvarn(varn(px));
! 261: *pr=p3; return z;
! 262: }
! 263:
! 264: p3=cgetg(i+3,t_POL);
! 265: p3[1]=evalsigne(1) | evallgef(i+3) | evalvarn(varn(px));
! 266: p3[i+2]=(long)nfreducemodpr(nf,p1,prhall);
! 267: for (k=i-1; k>=0; --k)
! 268: {
! 269: l=avma; p1=((GEN)px[k+2]);
! 270: for (j=0; j<=k && j<=dz; j++)
! 271: {
! 272: p2=element_mul(nf,(GEN)z[j+2],(GEN)py[k-j+2]);
! 273: p2=nfreducemodpr(nf,p2,prhall); tetpil=avma; p1=gsub(p1,p2);
! 274: }
! 275: p3[k+2]=lpile(l,tetpil,nfreducemodpr(nf,p1,prhall));
! 276: }
! 277: *pr=p3; return z;
! 278: }
! 279:
! 280: /* PGCD des polynomes x et y, par l'algorithme du sub-resultant */
! 281: static GEN
! 282: nf_pol_subres(GEN nf,GEN x,GEN y)
! 283: {
! 284: long av=avma,tetpil;
! 285: GEN s = srgcd(unifpol(nf,x,1), unifpol(nf,y,1));
! 286:
! 287: tetpil=avma; return gerepile(av,tetpil,unifpol(nf,s,1));
! 288: }
! 289:
! 290: /* PGCD des polynomes x et y modulo prhall */
! 291: static GEN
! 292: nfmod_pol_gcd(GEN nf,GEN prhall,GEN x,GEN y)
! 293: {
! 294: long av=avma;
! 295: GEN p1,p2;
! 296:
! 297: if (lgef(x)<lgef(y)) { p1=y; y=x; x=p1; }
! 298: p1=nfmod_pol_reduce(nf,prhall,x);
! 299: p2=nfmod_pol_reduce(nf,prhall,y);
! 300: while (!isexactzero(p2))
! 301: {
! 302: GEN p3;
! 303:
! 304: nfmod_pol_divres(nf,prhall,p1,p2,&p3);
! 305: p1=p2; p2=p3;
! 306: }
! 307: return gerepileupto(av,p1);
! 308: }
! 309:
! 310: /* Calcul pol^e modulo prhall et le polynome pmod */
! 311: static GEN
! 312: nfmod_pol_pow(GEN nf,GEN prhall,GEN pmod,GEN pol,GEN e)
! 313: {
! 314: long i, av = avma, n = lgef(nf[1])-3;
! 315: GEN p1,p2,vun;
! 316:
! 317: vun=cgetg(n+1,t_COL); vun[1]=un; for (i=2; i<=n; i++) vun[i]=zero;
! 318: p1=gcopy(polun[varn(pol)]); p1[2]=(long)vun;
! 319: if (gcmp0(e)) return p1;
! 320:
! 321: p2=nfmod_pol_reduce(nf,prhall,pol);
! 322: for(;;)
! 323: {
! 324: if (!vali(e))
! 325: {
! 326: p1=nfmod_pol_mul(nf,prhall,p1,p2);
! 327: nfmod_pol_divres(nf,prhall,p1,pmod,&p1);
! 328: }
! 329: if (gcmp1(e)) break;
! 330:
! 331: e=shifti(e,-1);
! 332: p2=nfmod_pol_sqr(nf,prhall,p2);
! 333: nfmod_pol_divres(nf,prhall,p2,pmod,&p2);
! 334: }
! 335: return gerepileupto(av,p1);
! 336: }
! 337:
! 338: /* factorisation du polynome x modulo pr */
! 339: GEN
! 340: nffactormod(GEN nf,GEN pol,GEN pr)
! 341: {
! 342: long av = avma, tetpil,lb,nbfact,psim,N,n,i,j,k,d,e,vf,r,kk;
! 343: GEN y,ex,*t,f1,f2,f3,df1,g1,polb,pold,polu,vker;
! 344: GEN Q,f,x,u,v,v2,v3,vz,q,vun,vzero,prhall;
! 345:
! 346: nf=checknf(nf);
! 347: if (typ(pol)!=t_POL) err(typeer,"nffactormod");
! 348: if (varn(pol) >= varn(nf[1]))
! 349: err(talker,"polynomial variable must have highest priority in nffactormod");
! 350:
! 351: prhall=nfmodprinit(nf,pr); n=lgef(nf[1])-3;
! 352: vun = gscalcol_i(gun, n);
! 353: vzero = gscalcol_i(gzero, n);
! 354:
! 355: f=unifpol(nf,pol,0); f=nfmod_pol_reduce(nf,prhall,f);
! 356: d=lgef(f)-3; vf=varn(f);
! 357: t = (GEN*)new_chunk(d+1);
! 358: ex= new_chunk(d+1);
! 359: x=gcopy(polx[vf]); x[3]=(long)vun; x[2]=(long)vzero;
! 360: if (d<=1) { nbfact=2; t[1] = f; ex[1]=1; }
! 361: else
! 362: {
! 363: q = (GEN)pr[1]; psim = VERYBIGINT;
! 364: if (!is_bigint(q)) psim = itos(q);
! 365: /* psim has an effect only when p is small. If too big, set it to a huge
! 366: * number (i.e ignore it) to avoid an error in itos on next line.
! 367: */
! 368: q=gpuigs(q, itos((GEN)pr[4]));
! 369: f1=f; e=1; nbfact=1;
! 370: while (lgef(f1)>3)
! 371: {
! 372: df1=deriv(f1,vf); f2=nfmod_pol_gcd(nf,prhall,f1,df1);
! 373: g1=nfmod_pol_divres(nf,prhall,f1,f2,NULL); k=0;
! 374: while (lgef(g1)>3)
! 375: {
! 376: k++;
! 377: if (k%psim == 0)
! 378: {
! 379: k++; f2=nfmod_pol_divres(nf,prhall,f2,g1,NULL);
! 380: }
! 381: f3=nfmod_pol_gcd(nf,prhall,f2,g1);
! 382: u = nfmod_pol_divres(nf,prhall,g1,f3,NULL);
! 383: f2= nfmod_pol_divres(nf,prhall,f2,f3,NULL);
! 384: g1=f3;
! 385: if (lgef(u)>3)
! 386: {
! 387: N=lgef(u)-3; Q=cgetg(N+1,t_MAT);
! 388: v3=cgetg(N+1,t_COL); Q[1]=(long)v3;
! 389: v3[1]=(long)vun; for (i=2; i<=N; i++) v3[i]=(long)vzero;
! 390:
! 391: v2 = v = nfmod_pol_pow(nf,prhall,u,x,q);
! 392: for (j=2; j<=N; j++)
! 393: {
! 394: v3=cgetg(N+1,t_COL); Q[j]=(long)v3;
! 395: for (i=1; i<=lgef(v2)-2; i++) v3[i]=v2[i+1];
! 396: for (; i<=N; i++) v3[i]=(long)vzero;
! 397: if (j<N)
! 398: {
! 399: v2=nfmod_pol_mul(nf,prhall,v2,v);
! 400: nfmod_pol_divres(nf,prhall,v2,u,&v2);
! 401: }
! 402: }
! 403: for (i=1; i<=N; i++)
! 404: coeff(Q,i,i)=lsub((GEN)coeff(Q,i,i),vun);
! 405: v2=nfkermodpr(nf,Q,prhall); r=lg(v2)-1; t[nbfact]=gcopy(u); kk=1;
! 406: if (r>1)
! 407: {
! 408: vker=cgetg(r+1,t_COL);
! 409: for (i=1; i<=r; i++)
! 410: {
! 411: v3=cgetg(N+2,t_POL);
! 412: v3[1]=evalsigne(1)+evallgef(2+N); setvarn(v3,vf);
! 413: vker[i]=(long)v3; for (j=1; j<=N; j++) v3[j+1]=coeff(v2,j,i);
! 414: normalizepol(v3);
! 415: }
! 416: }
! 417: while (kk<r)
! 418: {
! 419: v=gcopy(polun[vf]); v[2]=(long)vzero;
! 420: for (i=1; i<=r; i++)
! 421: {
! 422: vz=cgetg(n+1,t_COL);
! 423: for (j=1; j<=n; j++)
! 424: vz[j] = lmodsi(mymyrand()>>8, q);
! 425: vz=nfreducemodpr(nf,vz,prhall);
! 426: v=gadd(v,nfmod_pol_mul(nf,prhall,vz,(GEN)vker[i]));
! 427: }
! 428: for (i=1; i<=kk && kk<r; i++)
! 429: {
! 430: polb=t[nbfact+i-1]; lb=lgef(polb);
! 431: if (lb>4)
! 432: {
! 433: if(psim==2)
! 434: polu=nfmod_split2(nf,prhall,polb,v,q);
! 435: else
! 436: {
! 437: polu=nfmod_pol_pow(nf,prhall,polb,v,shifti(q,-1));
! 438: polu=gsub(polu,vun);
! 439: }
! 440: pold=nfmod_pol_gcd(nf,prhall,polb,polu);
! 441: if (lgef(pold)>3 && lgef(pold)<lb)
! 442: {
! 443: t[nbfact+i-1]=pold; kk++;
! 444: t[nbfact+kk-1]=nfmod_pol_divres(nf,prhall,polb,pold,NULL);
! 445: }
! 446: }
! 447: }
! 448: }
! 449: for (i=nbfact; i<nbfact+r; i++) ex[i]=e*k;
! 450: nbfact+=r;
! 451: }
! 452: }
! 453: e*=psim; j=(lgef(f2)-3)/psim+3; f1=cgetg(j,t_POL);
! 454: f1[1] = evalsigne(1) | evallgef(j) | evalvarn(vf);
! 455: for (i=2; i<j; i++)
! 456: f1[i]=(long)element_powmodpr(nf,(GEN)f2[psim*(i-2)+2],
! 457: gdiv(q,(GEN)pr[1]),prhall);
! 458: }
! 459: }
! 460: if (nbfact < 2)
! 461: err(talker, "%Z not a prime (nffactormod)", q);
! 462: v=element_divmodpr(nf,vun,gmael(t,1,lgef(t[1])-1),prhall);
! 463: t[1]=unifpol(nf,nfmod_pol_mul(nf,prhall,v,(GEN)t[1]),1);
! 464: for (j=2; j<nbfact; j++)
! 465: if (ex[j])
! 466: {
! 467: v=element_divmodpr(nf,vun,gmael(t,j,lgef(t[j])-1),prhall);
! 468: t[j]=unifpol(nf,nfmod_pol_mul(nf,prhall,v,(GEN)t[j]),1);
! 469: }
! 470:
! 471: tetpil=avma; y=cgetg(3,t_MAT);
! 472: u=cgetg(nbfact,t_COL); y[1]=(long)u;
! 473: v=cgetg(nbfact,t_COL); y[2]=(long)v;
! 474: for (j=1,k=0; j<nbfact; j++)
! 475: if (ex[j])
! 476: {
! 477: k++;
! 478: u[k]=lcopy((GEN)t[j]);
! 479: v[k]=lstoi(ex[j]);
! 480: }
! 481: return gerepile(av,tetpil,y);
! 482: }
! 483:
! 484: /* Calcule pol + pol^2 + ... + pol^(q/2) modulo prhall et le polynome pmod */
! 485: static GEN
! 486: nfmod_split2(GEN nf,GEN prhall,GEN pmod,GEN pol,GEN exp)
! 487: {
! 488: long av = avma;
! 489: GEN p1,p2,q;
! 490:
! 491: if (cmpis(exp,2)<=0) return pol;
! 492: p2=p1=pol; q=shifti(exp,-1);
! 493: while (!gcmp1(q))
! 494: {
! 495: p2=nfmod_pol_sqr(nf,prhall,p2);
! 496: nfmod_pol_divres(nf,prhall,p2,pmod,&p2);
! 497: q=shifti(q,-1); p1=gadd(p1,p2);
! 498: }
! 499: return gerepileupto(av,p1);
! 500: }
! 501:
! 502: /* If p doesn't divide either a or b and has a divisor of degree 1, return it.
! 503: * Return NULL otherwise.
! 504: */
! 505: static GEN
! 506: p_ok(GEN nf, GEN p, GEN a)
! 507: {
! 508: long av,m,i;
! 509: GEN dec;
! 510:
! 511: if (divise(a,p)) return NULL;
! 512: av = avma; dec = primedec(nf,p); m=lg(dec);
! 513: for (i=1; i<m; i++)
! 514: {
! 515: GEN pr = (GEN)dec[i];
! 516: if (is_pm1(pr[4]))
! 517: {
! 518: if (DEBUGLEVEL>=4)
! 519: fprintferr("Premier choisi pour decomposition: %Z\n", pr);
! 520: return pr;
! 521: }
! 522: }
! 523: avma = av; return NULL;
! 524: }
! 525:
! 526: static GEN
! 527: choose_prime(GEN nf, GEN dk, GEN lim)
! 528: {
! 529: GEN p, pr;
! 530:
! 531: p = nextprime(lim);
! 532: for (;;)
! 533: {
! 534: if ((pr = p_ok(nf,p,dk))) break;
! 535: p = nextprime(addis(p,2));
! 536: }
! 537:
! 538: return pr;
! 539: }
! 540:
! 541: /* Renvoie les racines de pol contenues dans nf */
! 542: GEN
! 543: nfroots(GEN nf,GEN pol)
! 544: {
! 545: long av=avma,tetpil,i,d=lgef(pol);
! 546: GEN p1,p2,polbase,polmod,den;
! 547:
! 548: nf=checknf(nf);
! 549: if (typ(pol)!=t_POL) err(talker,"not a polynomial in nfroots");
! 550: if (varn(pol) >= varn(nf[1]))
! 551: err(talker,"polynomial variable must have highest priority in nfroots");
! 552:
! 553: polbase=unifpol(nf,pol,0);
! 554:
! 555: if (d==3)
! 556: {
! 557: tetpil=avma; p1=cgetg(1,t_VEC);
! 558: return gerepile(av,tetpil,p1);
! 559: }
! 560:
! 561: if (d==4)
! 562: {
! 563: tetpil=avma; p1=cgetg(2,t_VEC);
! 564: p1[1] = (long)basistoalg(nf,gneg_i(
! 565: element_div(nf,(GEN)polbase[2],(GEN)polbase[3])));
! 566: return gerepile(av,tetpil,p1);
! 567: }
! 568:
! 569: p1=element_inv(nf,leading_term(polbase));
! 570: polbase=nf_pol_mul(nf,p1,polbase);
! 571:
! 572: den=gun;
! 573: for (i=2; i<d; i++)
! 574: if (! gcmp0((GEN)polbase[i]))
! 575: den = glcm(den,denom((GEN)polbase[i]));
! 576: if (! gcmp1(absi(den)))
! 577: for (i=2; i<d; i++)
! 578: polbase[i] = lmul(den,(GEN)polbase[i]);
! 579:
! 580: polmod=unifpol(nf,polbase,1);
! 581:
! 582: if (DEBUGLEVEL>=4)
! 583: fprintferr("On teste si le polynome est square-free\n");
! 584:
! 585: p1=derivpol(polmod);
! 586: p2=nf_pol_subres(nf,polmod,p1);
! 587:
! 588: if (degree(p2) > 0)
! 589: {
! 590: p1=element_inv(nf,leading_term(p2));
! 591: p2=nf_pol_mul(nf,p1,p2);
! 592: polmod=nf_pol_divres(nf,polmod,p2,NULL);
! 593:
! 594: p1=element_inv(nf,leading_term(polmod));
! 595: polmod=nf_pol_mul(nf,p1,polmod);
! 596:
! 597: d=lgef(polmod); den=gun;
! 598: for (i=2; i<d; i++)
! 599: if (!gcmp0((GEN)polmod[i]))
! 600: den = glcm(den,denom((GEN)polmod[i]));
! 601: if (!gcmp1(absi(den)))
! 602: for (i=2; i<d; i++)
! 603: polmod[i] = lmul(den,(GEN)polmod[i]);
! 604:
! 605: polmod=unifpol(nf,polmod,1);
! 606: }
! 607:
! 608: p1 = nfsqff(nf,polmod,1);
! 609: tetpil=avma; return gerepile(av, tetpil, gen_sort(p1, 0, cmp_pol));
! 610: }
! 611:
! 612: /* Relevement de elt modulo id minimal */
! 613: static GEN
! 614: nf_bestlift(GEN id,GEN idinv,GEN den,GEN elt)
! 615: {
! 616: return gsub(elt,gmul(id,ground(gmul(den,gmul(idinv,elt)))));
! 617: }
! 618:
! 619: /* Releve le polynome pol avec des coeff de norme t2 <= C si possible */
! 620: static GEN
! 621: nf_pol_lift(GEN id,GEN idinv,GEN den,GEN pol)
! 622: {
! 623: long i, d = lgef(pol);
! 624: GEN p1 = pol;
! 625:
! 626: pol=cgetg(d,t_POL); pol[1]=p1[1];
! 627: for (i=2; i<d; i++)
! 628: pol[i] = (long) nf_bestlift(id,idinv,den,(GEN)p1[i]);
! 629: return pol;
! 630: }
! 631:
! 632: #if 0
! 633: /* Evalue le polynome pol en elt */
! 634: static GEN
! 635: nf_pol_eval(GEN nf,GEN pol,GEN elt)
! 636: {
! 637: long av=avma,tetpil,i;
! 638: GEN p1;
! 639:
! 640: i=lgef(pol)-1; if (i==2) return gcopy((GEN)pol[2]);
! 641:
! 642: p1=element_mul(nf,(GEN)pol[i],elt);
! 643: for (i--; i>=3; i--)
! 644: p1=element_mul(nf,elt,gadd((GEN)pol[i],p1));
! 645: tetpil=avma; return gerepile(av,tetpil,gadd(p1,(GEN)pol[2]));
! 646: }
! 647: #endif
! 648:
! 649: /* Calcule la factorisation du polynome x dans nf */
! 650: GEN
! 651: nffactor(GEN nf,GEN pol)
! 652: { GEN y,p1,p2,den,p3,p4,quot, rep = cgetg(3,t_MAT);
! 653: long av = avma,tetpil,i,d;
! 654:
! 655: if (DEBUGLEVEL >= 4) timer2();
! 656:
! 657: nf=checknf(nf);
! 658: if (typ(pol)!=t_POL) err(typeer,"nffactor");
! 659: if (varn(pol) >= varn(nf[1]))
! 660: err(talker,"polynomial variable must have highest priority in nffactor");
! 661:
! 662: d=lgef(pol);
! 663: if (d==3)
! 664: {
! 665: rep[1]=lgetg(1,t_COL);
! 666: rep[2]=lgetg(1,t_COL);
! 667: return rep;
! 668: }
! 669: if (d==4)
! 670: {
! 671: p1=cgetg(2,t_COL); rep[1]=(long)p1; p1[1]=lcopy(pol);
! 672: p1=cgetg(2,t_COL); rep[2]=(long)p1; p1[1]=un;
! 673: return rep;
! 674: }
! 675:
! 676: p1=element_inv(nf,leading_term(pol));
! 677: pol=nf_pol_mul(nf,p1,pol);
! 678:
! 679: p1=unifpol(nf,pol,0); den=gun;
! 680: for (i=2; i<d; i++)
! 681: if (! gcmp0((GEN)p1[i]))
! 682: den = glcm(den,denom((GEN)p1[i]));
! 683: if (! gcmp1(absi(den)))
! 684: for (i=2; i<d; i++)
! 685: p1[i] = lmul(den,(GEN)p1[i]);
! 686:
! 687: if (DEBUGLEVEL>=4)
! 688: fprintferr("On teste si le polynome est square-free\n");
! 689:
! 690: p2=derivpol(p1);
! 691: p3=nf_pol_subres(nf,p1,p2);
! 692:
! 693: if (degree(p3) > 0)
! 694: {
! 695: p4=element_inv(nf,leading_term(p3));
! 696: p3=nf_pol_mul(nf,p4,p3);
! 697:
! 698: p2=nf_pol_divres(nf,p1,p3,NULL);
! 699: p4=element_inv(nf,leading_term(p2));
! 700: p2=nf_pol_mul(nf,p4,p2);
! 701:
! 702: d=lgef(p2); den=gun;
! 703: for (i=2; i<d; i++)
! 704: if (!gcmp0((GEN)p2[i]))
! 705: den = glcm(den,denom((GEN)p2[i]));
! 706: if (!gcmp1(absi(den)))
! 707: for (i=2; i<d; i++)
! 708: p2[i] = lmul(den,(GEN)p2[i]);
! 709:
! 710: p2=unifpol(nf,p2,1);
! 711: tetpil = avma; y = nfsqff(nf,p2,0);
! 712: i = nfcmbf.nfact;
! 713:
! 714: quot=nf_pol_divres(nf,p1,p2,NULL);
! 715: p3=(GEN)gpmalloc((i+1) * sizeof(long));
! 716: for ( ; i>=1; i--)
! 717: {
! 718: GEN fact=(GEN)y[i], quo = quot, rem;
! 719: long e = 0;
! 720:
! 721: do
! 722: {
! 723: quo = nf_pol_divres(nf,quo,fact,&rem);
! 724: e++;
! 725: }
! 726: while (gcmp0(rem));
! 727: p3[i]=lstoi(e);
! 728: }
! 729: avma = (long)y; y = gerepile(av, tetpil, y);
! 730: p2=cgetg(i+1, t_COL); for (; i>=1; i--) p2[i]=lstoi(p3[i]);
! 731: free(p3);
! 732: }
! 733: else
! 734: {
! 735: tetpil=avma; y = gerepile(av,tetpil,nfsqff(nf,p1,0));
! 736: i = nfcmbf.nfact;
! 737: p2=cgetg(i+1, t_COL); for (; i>=1; i--) p2[i]=un;
! 738: }
! 739: if (DEBUGLEVEL>=4)
! 740: fprintferr("Nombre de facteur(s) trouve(s) : %ld\n", nfcmbf.nfact);
! 741: rep[1]=(long)y;
! 742: rep[2]=(long)p2; return sort_factor(rep, cmp_pol);
! 743: }
! 744:
! 745: /* teste si la matrice M est suffisament orthogonale */
! 746: static long
! 747: test_mat(GEN M, GEN p, GEN C2, long k)
! 748: {
! 749: long av = avma, i, N = lg(M);
! 750: GEN min, prod, L2, R;
! 751:
! 752: min = prod = gcoeff(M,1,1);
! 753: for (i = 2; i < N; i++)
! 754: {
! 755: L2 = gcoeff(M,i,i); prod = mpmul(prod,L2);
! 756: if (mpcmp(L2,min) < 0) min = L2;
! 757: }
! 758: R = mpmul(min, gpowgs(p, k<<1));
! 759: i = mpcmp(mpmul(C2,prod), R); avma = av;
! 760: return (i < 0);
! 761: }
! 762:
! 763: /* calcule la matrice correspondant a pr^e jusqu'a ce que R(pr^e) > C */
! 764: static GEN
! 765: T2_matrix_pow(GEN nf, GEN pre, GEN p, GEN C, long *kmax, long prec)
! 766: {
! 767: long av=avma,av1,lim, k = *kmax, N = lgef((GEN)nf[1])-3;
! 768: int tot_real = !signe(gmael(nf,2,2));
! 769: GEN p1,p2,p3,u,C2,T2, x = (GEN)nf[1];
! 770:
! 771: C2 = gdiv(gmul2n(C,2), absi((GEN)nf[3]));
! 772: p1 = gmul(pre,lllintpartial(pre)); av1 = avma;
! 773: T2 = tot_real? gmael(nf,5,4)
! 774: : nf_get_T2((GEN) nf[7], roots(x,prec));
! 775: p3 = qf_base_change(T2,p1,1);
! 776:
! 777: if (N <= 6 && test_mat(p3,p,C2,k))
! 778: {
! 779: avma = av1; return gerepileupto(av,p1);
! 780: }
! 781:
! 782: av1=avma; lim = stack_lim(av1,2);
! 783: for (;;)
! 784: {
! 785: if (DEBUGLEVEL>2) fprintferr("exponent: %ld\n",k);
! 786:
! 787: for(;;)
! 788: {
! 789: u = tot_real? lllgramint(p3): lllgramintern(p3,100,1,prec);
! 790: if (u) break;
! 791:
! 792: prec=(prec<<1)-2;
! 793: if (DEBUGLEVEL > 1) err(warnprec,"nffactor[1]",prec);
! 794: T2 = nf_get_T2((GEN) nf[7], roots(x,prec));
! 795: p3 = qf_base_change(T2,p1,1);
! 796: }
! 797: if (DEBUGLEVEL>2) msgtimer("lllgram + base change");
! 798: p3 = qf_base_change(p3,u,1);
! 799: if (test_mat(p3,p,C2,k))
! 800: {
! 801: *kmax = k;
! 802: return gerepileupto(av,gmul(p1,u));
! 803: }
! 804:
! 805: /* il faut augmenter la precision en meme temps */
! 806: p2=shifti(gceil(mulsr(k, glog(p,DEFAULTPREC))),-1);
! 807: prec += (long)(itos(p2)*pariK1);
! 808: if (DEBUGLEVEL > 1) err(warnprec,"nffactor[2]",prec);
! 809: k = k<<1; p1 = idealmullll(nf,p1,p1);
! 810: if (low_stack(lim, stack_lim(av1,2)))
! 811: {
! 812: if (DEBUGMEM>1) err(warnmem,"T2_matrix_pow");
! 813: p1 = gerepileupto(av1,p1);
! 814: }
! 815: if (!tot_real) T2 = nf_get_T2((GEN) nf[7], roots(x,prec));
! 816: p3 = qf_base_change(T2,p1,1);
! 817: }
! 818: }
! 819:
! 820: /* Calcule la factorisation du polynome x qui est sans facteurs carres dans nf,
! 821: Le polynome est a coefficients dans Z_k et son terme dominant est un entier
! 822: rationnel, si fl=1 renvoie seulement les racines du polynome contenues
! 823: dans le corps */
! 824: static GEN
! 825: nfsqff(GEN nf,GEN pol, long fl)
! 826: {
! 827: long d=lgef(pol),i,k,m,n,av=avma,tetpil,newprec,prec;
! 828: GEN p1,pr,p2,rep,k2,C,h,dk,dki,p,prh,p3,T2,polbase,fact,pk;
! 829: GEN polmod,polred,hinv,lt,minp,den,maxp=shifti(gun,32),run;
! 830:
! 831: dk=absi((GEN)nf[3]);
! 832: dki=mulii(dk,(GEN)nf[4]);
! 833: n=lgef(nf[1])-3;
! 834:
! 835: polbase = unifpol(nf,pol,0);
! 836: polmod = unifpol(nf,pol,1);
! 837: dki=mulii(dki,gnorm((GEN)polmod[d-1]));
! 838:
! 839: prec = DEFAULTPREC;
! 840: for (;;)
! 841: {
! 842: if (prec <= gprecision(nf))
! 843: T2 = gprec_w(gmael(nf,5,3), prec);
! 844: else
! 845: T2 = nf_get_T2((GEN)nf[7], roots((GEN)nf[1], prec));
! 846:
! 847: run=realun(prec);
! 848: p1=realzero(prec);
! 849: for (i=2; i<d; i++)
! 850: {
! 851: p2 = gmul(run, (GEN)polbase[i]);
! 852: p1 = addrr(p1, qfeval(T2, p2));
! 853: if (signe(p1) < 0) break;
! 854: }
! 855:
! 856: if (signe(p1) > 0) break;
! 857:
! 858: prec = (prec<<1)-2;
! 859: if (DEBUGLEVEL > 1) err(warnprec, "nfsqff", prec);
! 860: }
! 861:
! 862: if (DEBUGLEVEL>=4)
! 863: fprintferr("La norme de ce polynome est : %Z\n", p1);
! 864:
! 865: lt=(GEN)leading_term(polbase)[1];
! 866: p2=mulis(sqri(lt),n);
! 867: C=mpadd(p1,p2);
! 868: C=mpadd(C,gsqrt(gmul(p1,p2),prec));
! 869:
! 870: if (!fl)
! 871: C=mpmul(C,sqri(binome(shifti(stoi(n-1),-1),(n-1)>>2)));
! 872:
! 873: C=gmul(C,sqri(lt));
! 874:
! 875: if (DEBUGLEVEL>=4)
! 876: fprintferr("La borne de la norme des coeff du diviseur est : %Z\n", C);
! 877:
! 878: k2=gmul2n(gmulgs(glog(gdivgs(gmul2n(C,2),n),DEFAULTPREC),n),-1);
! 879: if (!fl) k2=mulrr(k2,dbltor(1.1));
! 880:
! 881: minp=gmin(gexp(gmul2n(k2,-6),BIGDEFAULTPREC), maxp);
! 882: minp=gceil(minp);
! 883:
! 884: if (DEBUGLEVEL>=4)
! 885: {
! 886: fprintferr("borne inf. sur les nombres premiers : %Z\n", minp);
! 887: msgtimer("Calcul des bornes");
! 888: }
! 889: for (;;)
! 890: {
! 891: pr=choose_prime(nf,dki,minp); p=(GEN)pr[1];
! 892: prh=prime_to_ideal(nf,pr);
! 893:
! 894: polred=gcopy(polbase);
! 895: lt=(GEN)leading_term(polbase)[1];
! 896: lt=mpinvmod(lt,p);
! 897:
! 898: for (i=2; i<d; i++)
! 899: polred[i] = nfreducemodpr_i(gmul(lt,(GEN)polbase[i]), prh)[1];
! 900: if (!divise(discsr(polred),p)) break;
! 901: minp = addis(p,1);
! 902: }
! 903:
! 904: k = itos(gceil(gdiv(k2,glog(p,BIGDEFAULTPREC))));
! 905: rep=(GEN)simplefactmod(polred,p)[1];
! 906:
! 907: if (DEBUGLEVEL>=4)
! 908: {
! 909: msgtimer("Choix de l'ideal premier");
! 910: fprintferr("Nombre de facteurs irreductibles modulo pr = %ld\n", lg(rep)-1);
! 911: }
! 912: if (lg(rep)==2)
! 913: {
! 914: if (fl) { avma=av; return cgetg(1,t_VEC); }
! 915: rep=cgetg(2,t_VEC); rep[1]=lcopy(polmod);
! 916: nfcmbf.nfact=1; return gerepileupto(av, rep);
! 917: }
! 918:
! 919: p2=cgetr(DEFAULTPREC);
! 920: affir(mulii(absi(dk),gpowgs(p,k)),p2);
! 921: p2=shifti(gceil(mplog(p2)),-1);
! 922:
! 923: newprec = MEDDEFAULTPREC + (long)(itos(p2)*pariK1);
! 924: if (DEBUGLEVEL>=4)
! 925: fprintferr("nouvelle precision : %ld\n",newprec);
! 926:
! 927: prh = idealpows(nf,pr,k); m = k;
! 928: h = T2_matrix_pow(nf,prh,p,C, &k, newprec);
! 929: if (m != k) prh=idealpows(nf,pr,k);
! 930:
! 931: if (DEBUGLEVEL>=4)
! 932: {
! 933: fprintferr("un exposant convenable est : %ld\n",(long)k);
! 934: msgtimer("Calcul de H");
! 935: }
! 936:
! 937: pk = gcoeff(prh,1,1);
! 938: lt=(GEN)leading_term(polbase)[1];
! 939: lt=mpinvmod(lt,pk);
! 940:
! 941: polred[1] = polbase[1];
! 942: for (i=2; i<d; i++)
! 943: polred[i] = nfreducemodpr_i(gmul(lt,(GEN)polbase[i]), prh)[1];
! 944:
! 945: fact = lift_intern((GEN)factmod(polred,p)[1]);
! 946: rep = hensel_lift_fact(polred,fact,p,pk,k);
! 947:
! 948: if (DEBUGLEVEL >= 4) msgtimer("Calcul de la factorisation pr-adique");
! 949:
! 950: den=ginv(det(h));
! 951: hinv=adj(h);
! 952: lt=(GEN)leading_term(polbase)[1];
! 953:
! 954: if (fl)
! 955: {
! 956: long x_a[] = { evaltyp(t_POL)|m_evallg(4), 0,0,0 };
! 957: GEN mlt = negi(lt), rem;
! 958: x_a[1] = polbase[1]; setlgef(x_a, 4);
! 959: x_a[3] = un;
! 960: p1=cgetg(lg(rep)+1,t_VEC);
! 961: polbase = unifpol(nf,polbase,1);
! 962: for (m=1,i=1; i<lg(rep); i++)
! 963: {
! 964: p2=(GEN)rep[i]; if (lgef(p2)!=4) break;
! 965:
! 966: p3 = algtobasis(nf, gmul(mlt,(GEN)p2[2]));
! 967: p3 = nf_bestlift(h,hinv,den,p3);
! 968: p3 = gdiv(p3,lt);
! 969: x_a[2] = lneg(p3); /* check P(p3) == 0 */
! 970: p2 = poldivres(polbase, unifpol(nf,x_a,1), &rem);
! 971: if (!signe(rem)) { polbase = p2; p1[m++] = (long)p3; }
! 972: }
! 973: tetpil=avma; rep=cgetg(m,t_VEC);
! 974: for (i=1; i<m; i++) rep[i]=(long)basistoalg(nf,(GEN)p1[i]);
! 975: return gerepile(av,tetpil,rep);
! 976: }
! 977:
! 978: for (i=1; i<lg(rep); i++)
! 979: rep[i] = (long)unifpol(nf,(GEN)rep[i],0);
! 980:
! 981: nfcmbf.pol = polmod;
! 982: nfcmbf.lt = lt;
! 983: nfcmbf.h = h;
! 984: nfcmbf.hinv = hinv;
! 985: nfcmbf.den = den;
! 986: nfcmbf.fact = rep;
! 987: nfcmbf.res = cgetg(lg(rep)+1,t_VEC);
! 988: nfcmbf.nfact = 0;
! 989: nfcmbf.nfactmod = lg(rep)-1;
! 990: nf_combine_factors(nf,1,NULL,d-3,1);
! 991:
! 992: if (DEBUGLEVEL >= 4) msgtimer("Reconnaissance des facteurs");
! 993:
! 994: i = nfcmbf.nfact;
! 995: if (lgef(nfcmbf.pol)>3)
! 996: {
! 997: nfcmbf.res[++i] = (long) nf_pol_divres(nf,nfcmbf.pol,nfcmbf.lt,NULL);
! 998: nfcmbf.nfact = i;
! 999: }
! 1000:
! 1001: tetpil=avma; rep=cgetg(i+1,t_VEC);
! 1002: for ( ; i>=1; i--)
! 1003: rep[i]=(long)unifpol(nf,(GEN)nfcmbf.res[i],1);
! 1004: return gerepile(av,tetpil,rep);
! 1005: }
! 1006:
! 1007: static int
! 1008: nf_combine_factors(GEN nf,long fxn,GEN psf,long dlim,long hint)
! 1009: {
! 1010: int val=0; /* assume failure */
! 1011: GEN newf, newpsf = NULL;
! 1012: long newd,ltop,i;
! 1013:
! 1014: if (dlim<=0) return 0;
! 1015: if (fxn > nfcmbf.nfactmod) return 0;
! 1016: /* first, try deeper factors without considering the current one */
! 1017: if (fxn != nfcmbf.nfactmod)
! 1018: {
! 1019: val=nf_combine_factors(nf,fxn+1,psf,dlim,hint);
! 1020: if (val && psf) return 1;
! 1021: }
! 1022:
! 1023: /* second, try including the current modular factor in the product */
! 1024: newf=(GEN)nfcmbf.fact[fxn];
! 1025: if (!newf) return val; /* modular factor already used */
! 1026: newd=lgef(newf)-3;
! 1027: if (newd>dlim) return val; /* degree of new factor is too large */
! 1028:
! 1029: if (newd%hint == 0)
! 1030: {
! 1031: GEN p, quot,rem;
! 1032:
! 1033: newpsf = nf_pol_mul(nf, (psf)? psf: nfcmbf.lt, newf);
! 1034: newpsf = nf_pol_lift(nfcmbf.h,nfcmbf.hinv,nfcmbf.den,newpsf);
! 1035: /* try out the new combination */
! 1036: ltop=avma;
! 1037: quot=nf_pol_divres(nf,nfcmbf.pol,newpsf,&rem);
! 1038: if (gcmp0(rem)) /* found a factor */
! 1039: {
! 1040: p = nf_pol_mul(nf,element_inv(nf,leading_term(newpsf)),newpsf);
! 1041: nfcmbf.res[++nfcmbf.nfact] = (long) p; /* store factor */
! 1042: nfcmbf.fact[fxn]=0; /* remove used modular factor */
! 1043:
! 1044: /* fix up target */
! 1045: p=gun; quot=unifpol(nf,quot,0);
! 1046: for (i=2; i<lgef(quot); i++)
! 1047: if (!gcmp0((GEN)quot[i]))
! 1048: p = glcm(p, denom((GEN)quot[i]));
! 1049:
! 1050: nfcmbf.pol = nf_pol_mul(nf,p,quot);
! 1051: nfcmbf.lt = leading_term(nfcmbf.pol);
! 1052: return 1;
! 1053: }
! 1054: avma=ltop;
! 1055: }
! 1056:
! 1057: /* newpsf needs more; try for it */
! 1058: if (newd==dlim) return val; /* no more room in degree limit */
! 1059: if (fxn==nfcmbf.nfactmod) return val; /* no more modular factors to try */
! 1060:
! 1061: if (nf_combine_factors(nf,fxn+1,newpsf,dlim-newd,hint))
! 1062: {
! 1063: nfcmbf.fact[fxn]=0; /* remove used modular factor */
! 1064: return 1;
! 1065: }
! 1066: return val;
! 1067: }
! 1068:
! 1069: /* Calcule le polynome caracteristique de alpha sur nf ou alpha est un
! 1070: element de l'algebre nf[X]/(T) exprime comme polynome en x */
! 1071: GEN
! 1072: rnfcharpoly(GEN nf,GEN T,GEN alpha,int v)
! 1073: {
! 1074: long av=avma,tetpil;
! 1075: GEN p1;
! 1076:
! 1077: nf=checknf(nf); if (v<0) v = 0;
! 1078: if (typ(alpha) == t_POLMOD) alpha = lift_to_pol(alpha);
! 1079: p1 = caract2(unifpol(nf,T,1), unifpol(nf,alpha,1), v);
! 1080: tetpil=avma; return gerepile(av,tetpil,unifpol(nf,p1,1));
! 1081: }
! 1082:
! 1083: #if 0
! 1084: /* Calcule le polynome minimal de alpha sur nf ou alpha est un
! 1085: element de l'algebre nf[X]/(T) exprime comme polynome en x */
! 1086: GEN
! 1087: rnfminpoly(GEN nf,GEN T,GEN alpha,int n)
! 1088: {
! 1089: long av=avma,tetpil;
! 1090: GEN p1,p2;
! 1091:
! 1092: nf=checknf(nf); p1=rnfcharpoly(nf,T,alpha,n);
! 1093: tetpil=avma; p2=nf_pol_subres(nf,p1,deriv(p1,varn(T)));
! 1094: if (lgef(p2)==3) { avma=tetpil; return p1; }
! 1095:
! 1096: p1 = nf_pol_divres(nf,p1,p2,NULL);
! 1097: p2 = element_inv(nf,leading_term(p1));
! 1098: tetpil=avma; return gerepile(av,tetpil,unifpol(nf,nf_pol_mul(nf,p2,p1),1));
! 1099: }
! 1100: #endif
! 1101:
! 1102: /* relative Dedekind criterion over nf, applied to the order defined by a
! 1103: * root of irreducible polynomial T, modulo the prime ideal pr. Returns
! 1104: * [flag,basis,val], where basis is a pseudo-basis of the enlarged order,
! 1105: * flag is 1 iff this order is pr-maximal, and val is the valuation in pr of
! 1106: * the order discriminant
! 1107: */
! 1108: GEN
! 1109: rnfdedekind(GEN nf,GEN T,GEN pr)
! 1110: {
! 1111: long av=avma,vt,tetpil,r,d,da,n,m,i,j;
! 1112: GEN p1,p2,p,tau,g,vecun,veczero,matid;
! 1113: GEN prhall,res,h,k,base,Ca;
! 1114:
! 1115: nf=checknf(nf); Ca=unifpol(nf,T,0);
! 1116: res=cgetg(4,t_VEC);
! 1117: if (typ(pr)==t_VEC && lg(pr)==3)
! 1118: { prhall = (GEN)pr[2]; pr = (GEN)pr[1]; }
! 1119: else
! 1120: prhall=nfmodprinit(nf,pr);
! 1121: p=(GEN)pr[1]; tau=gdiv((GEN)pr[5],p);
! 1122: n=lgef(nf[1])-3; m=lgef(T)-3;
! 1123:
! 1124: vecun=cgetg(n+1,t_COL); vecun[1]=un;
! 1125: veczero=cgetg(n+1,t_COL); veczero[1]=zero;
! 1126: for (i=2; i<=n; i++) vecun[i] = veczero[i] = zero;
! 1127: matid=idmat(n);
! 1128:
! 1129: p1=(GEN)nffactormod(nf,Ca,pr)[1];
! 1130: g=lift((GEN)p1[1]); r=lg(p1);
! 1131: for (i=2; i<r; i++)
! 1132: g = nf_pol_mul(nf,g, lift((GEN)p1[i]));
! 1133: h=nfmod_pol_divres(nf,prhall,Ca,g,NULL);
! 1134: k=nf_pol_mul(nf,tau,gsub(Ca, nf_pol_mul(nf,lift(g),lift(h))));
! 1135: p2=nfmod_pol_gcd(nf,prhall,g,h);
! 1136: k= nfmod_pol_gcd(nf,prhall,p2,k);
! 1137:
! 1138: d=lgef(k)-3;
! 1139: vt = idealval(nf,discsr(T),pr) - 2*d;
! 1140: res[3]=lstoi(vt);
! 1141: if (!d || vt<=1) res[1]=un; else res[1]=zero;
! 1142:
! 1143: base=cgetg(3,t_VEC);
! 1144: p1=cgetg(m+d+1,t_MAT); base[1]=(long)p1;
! 1145: p2=cgetg(m+d+1,t_VEC); base[2]=(long)p2;
! 1146: for (j=1; j<=m; j++)
! 1147: {
! 1148: p2[j]=(long)matid;
! 1149: p1[j]=lgetg(m+1,t_COL);
! 1150: for (i=1; i<=m; i++)
! 1151: coeff(p1,i,j) = (i==j)?(long)vecun:(long)veczero;
! 1152: }
! 1153:
! 1154: if (d)
! 1155: {
! 1156: GEN pal = lift(nfmod_pol_divres(nf,prhall,Ca,k,NULL));
! 1157: GEN prinv=idealinv(nf,pr);
! 1158: GEN nfx=unifpol(nf,polx[varn(T)],0);
! 1159:
! 1160: for (j=m+1; j<=m+d; j++)
! 1161: {
! 1162: p1[j]=lgetg(m+1,t_COL);
! 1163: da=lgef(pal)-3;
! 1164: for (i=1; i<=da+1; i++) coeff(p1,i,j)=pal[i+1];
! 1165: for ( ; i<=m; i++) coeff(p1,i,j)=(long)veczero;
! 1166: p2[j]=(long)prinv;
! 1167: nf_pol_divres(nf,nf_pol_mul(nf,pal,nfx),T,&pal);
! 1168: }
! 1169: base=nfhermite(nf,base);
! 1170: }
! 1171: res[2]=(long)base; tetpil=avma; return gerepile(av,tetpil,gcopy(res));
! 1172: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to nffactor.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