Annotation of OpenXM_contrib/pari/src/modules/nffactor.c, Revision 1.1.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