Annotation of OpenXM_contrib/pari/src/basemath/galconj.c, Revision 1.1.1.1
1.1 maekawa 1: /*************************************************************************/
2: /** **/
3: /** GALOIS CONJUGATES **/
4: /** **/
5: /*************************************************************************/
6: /* $Id: galconj.c,v 1.8 1999/09/23 17:50:56 karim Exp $ */
7: #include "pari.h"
8: GEN
9: galoisconj(GEN nf)
10: {
11: GEN x, y, z;
12: long i, lz, lx, v, av = avma;
13: nf = checknf(nf);
14: x = (GEN) nf[1];
15: v = varn(x);
16: lx = lgef(x) - 2;
17: if (v == 0)
18: nf = gsubst(nf, 0, polx[MAXVARN]);
19: else
20: {
21: x = dummycopy(x);
22: setvarn(x, 0);
23: }
24: z = nfroots(nf, x);
25: lz = lg(z);
26: y = cgetg(lz, t_COL);
27: for (i = 1; i < lz; i++)
28: {
29: GEN p1 = lift((GEN) z[i]);
30: setvarn(p1, v);
31: y[i] = (long) p1;
32: }
33: return gerepileupto(av, y);
34: }
35: /* nbmax: maximum number of possible conjugates */
36: GEN
37: galoisconj2pol(GEN x, long nbmax, long prec)
38: {
39: long av = avma, i, n, v, nbauto;
40: GEN y, w, polr, p1, p2;
41: n = lgef(x) - 3;
42: if (n <= 0)
43: return cgetg(1, t_VEC);
44: if (gisirreducible(x) == gzero)
45: err(redpoler, "galoisconj2pol");
46: polr = roots(x, prec);
47: p1 = (GEN) polr[1];
48: nbauto = 1;
49: prec = (long) (bit_accuracy(prec) * L2SL10 * 0.75);
50: w = cgetg(n + 2, t_VEC);
51: w[1] = un;
52: for (i = 2; i <= n; i++)
53: w[i] = lmul(p1, (GEN) w[i - 1]);
54: v = varn(x);
55: y = cgetg(nbmax + 1, t_COL);
56: y[1] = (long) polx[v];
57: for (i = 2; i <= n && nbauto < nbmax; i++)
58: {
59: w[n + 1] = polr[i];
60: p1 = lindep2(w, prec);
61: if (signe(p1[n + 1]))
62: {
63: setlg(p1, n + 1);
64: p2 = gdiv(gtopolyrev(p1, v), negi((GEN) p1[n + 1]));
65: if (gdivise(poleval(x, p2), x))
66: {
67: y[++nbauto] = (long) p2;
68: if (DEBUGLEVEL > 1)
69: fprintferr("conjugate %ld: %Z\n", i, y[nbauto]);
70: }
71: }
72: }
73: setlg(y, 1 + nbauto);
74: return gerepileupto(av, gen_sort(y, 0, cmp_pol));
75: }
76: GEN
77: galoisconj2(GEN nf, long nbmax, long prec)
78: {
79: long av = avma, i, j, n, r1, ru, nbauto;
80: GEN x, y, w, polr, p1, p2;
81: if (typ(nf) == t_POL)
82: return galoisconj2pol(nf, nbmax, prec);
83: nf = checknf(nf);
84: x = (GEN) nf[1];
85: n = lgef(x) - 3;
86: if (n <= 0)
87: return cgetg(1, t_VEC);
88: r1 = itos(gmael(nf, 2, 1));
89: p1 = (GEN) nf[6];
90: prec = precision((GEN) p1[1]);
91: /* accuracy in decimal digits */
92: prec = (long) (bit_accuracy(prec) * L2SL10 * 0.75);
93: ru = (n + r1) >> 1;
94: nbauto = 1;
95: polr = cgetg(n + 1, t_VEC);
96: for (i = 1; i <= r1; i++)
97: polr[i] = p1[i];
98: for (j = i; i <= ru; i++)
99: {
100: polr[j++] = p1[i];
101: polr[j++] = lconj((GEN) p1[i]);
102: }
103: p2 = gmael(nf, 5, 1);
104: w = cgetg(n + 2, t_VEC);
105: for (i = 1; i <= n; i++)
106: w[i] = coeff(p2, 1, i);
107: y = cgetg(nbmax + 1, t_COL);
108: y[1] = (long) polx[varn(x)];
109: for (i = 2; i <= n && nbauto < nbmax; i++)
110: {
111: w[n + 1] = polr[i];
112: p1 = lindep2(w, prec);
113: if (signe(p1[n + 1]))
114: {
115: setlg(p1, n + 1);
116: settyp(p1, t_COL);
117: p2 = gdiv(gmul((GEN) nf[7], p1), negi((GEN) p1[n + 1]));
118: if (gdivise(poleval(x, p2), x))
119: {
120: y[++nbauto] = (long) p2;
121: if (DEBUGLEVEL > 1)
122: fprintferr("conjugate %ld: %Z\n", i, y[nbauto]);
123: }
124: }
125: }
126: setlg(y, 1 + nbauto);
127: return gerepileupto(av, gen_sort(y, 0, cmp_pol));
128: }
129: /*************************************************************************/
130: /** **/
131: /** GALOIS4 **/
132: /** **/
133: /** **/
134: /** **/
135: /*************************************************************************/
136: struct galois_lift
137: {
138: GEN T;
139: GEN den;
140: GEN p;
141: GEN borne;
142: GEN L;
143: long e;
144: GEN Q;
145: };
146: /* Initialise la structure galois_lift */
147: void
148: initlift(GEN T, GEN den, GEN p, GEN borne, GEN L, struct galois_lift * gl)
149: {
150: ulong ltop, lbot;
151: gl->T = T;
152: gl->den = den;
153: gl->p = p;
154: gl->borne = borne;
155: gl->L = L;
156: ltop = avma;
157: gl->e = itos(gmax(gdeux, gceil(gdiv(glog(gmul2n(borne, 1), DEFAULTPREC),
158: glog(p, DEFAULTPREC)))));
159: lbot = avma;
160: gl->Q = gpowgs(p, gl->e);
161: gl->Q = gerepile(ltop, lbot, gl->Q);
162: }
163: /*
164: * T est le polynome \sum t_i*X^i S est un Polmod T Retourne un vecteur Spow
165: * verifiant la condition: i>=1 et t_i!=0 ==> Spow[i]=S^(i-1)*t_i Essaie
166: * d'etre efficace sur les polynomes lacunaires
167: */
168: GEN
169: compoTS(GEN T, GEN S)
170: {
171: GEN Spow;
172: int i;
173: Spow = cgetg(lgef(T) - 2, t_VEC);
174: for (i = 3; i < lg(Spow); i++)
175: Spow[i] = 0;
176: Spow[1] = un;
177: Spow[2] = (long) S;
178: for (i = 2; i < lg(Spow) - 1; i++)
179: {
180: if (!gcmp0((GEN) T[i + 3]))
181: for (;;)
182: {
183: int k, l;
184: for (k = 1; k <= (i >> 1); k++)
185: if (Spow[k + 1] && Spow[i - k + 1])
186: break;
187: if ((k << 1) < i)
188: {
189: Spow[i + 1] = lmul((GEN) Spow[k + 1], (GEN) Spow[i - k + 1]);
190: break;
191: } else if ((k << 1) == i) /* En esperant que sqr est plus
192: * rapide... */
193: {
194: Spow[i + 1] = lsqr((GEN) Spow[k + 1]);
195: break;
196: }
197: for (k = i - 1; k > 0; k--)
198: if (Spow[k + 1])
199: break;
200: if ((k << 1) < i)
201: {
202: Spow[(k << 1) + 1] = lsqr((GEN) Spow[k + 1]);
203: continue;
204: }
205: for (l = i - k; l > 0; l--)
206: if (Spow[l + 1])
207: break;
208: if (Spow[i - l - k + 1])
209: Spow[i - k + 1] = lmul((GEN) Spow[i - l - k + 1], (GEN) Spow[l + 1]);
210: else
211: Spow[l + k + 1] = lmul((GEN) Spow[k + 1], (GEN) Spow[l + 1]);
212: }
213: }
214: for (i = 1; i < lg(Spow); i++)
215: if (!gcmp0((GEN) T[i + 2]))
216: Spow[i] = lmul((GEN) Spow[i], (GEN) T[i + 2]);
217: return Spow;
218: }
219: /*
220: * Calcule T(S) a l'aide du vecteur Spow
221: */
222: static GEN
223: calcTS(GEN Spow, GEN T, GEN S)
224: {
225: GEN res = gzero;
226: int i;
227: for (i = 1; i < lg(Spow); i++)
228: if (!gcmp0((GEN) T[i + 2]))
229: res = gadd(res, (GEN) Spow[i]);
230: res = gmul(res, S);
231: return gadd(res, (GEN) T[2]);
232: }
233: /*
234: * Calcule T'(S) a l'aide du vecteur Spow
235: */
236: static GEN
237: calcderivTS(GEN Spow, GEN T, GEN mod)
238: {
239: GEN res = gzero;
240: int i;
241: for (i = 1; i < lg(Spow); i++)
242: if (!gcmp0((GEN) T[i + 2]))
243: res = gadd(res, gmul((GEN) Spow[i], stoi(i)));
244: res = gmul(lift(lift(res)), mod);
245: return gmodulcp(res, T);
246: }
247: /*
248: * Verifie que S est une solution presque surement et calcule sa permutation
249: */
250: static int
251: poltopermtest(GEN f, GEN L, GEN pf)
252: {
253: ulong ltop;
254: GEN fx, fp;
255: int i, j;
256: fp = cgetg(lg(L), t_VECSMALL);
257: ltop = avma;
258: for (i = 1; i < lg(fp); i++)
259: fp[i] = 1;
260: for (i = 1; i < lg(L); i++)
261: {
262: fx = gsubst(f, varn(f), (GEN) L[i]);
263: for (j = 1; j < lg(L); j++)
264: if (fp[j] && gegal(fx, (GEN) L[j]))
265: {
266: pf[i] = j;
267: fp[j] = 0;
268: break;
269: }
270: if (j == lg(L))
271: return 0;
272: avma = ltop;
273: }
274: return 1;
275: }
276: /*
277: * Soit T un polynome de \ZZ[X] , p un nombre premier , S\in\FF_p[X]/(T) tel
278: * que T(S)=0 [p,T] Relever S en S_0 tel que T(S_0)=0 [T,p^e]
279: */
280: GEN
281: automorphismlift(GEN S, struct galois_lift * gl)
282: {
283: ulong ltop, lbot;
284: long x;
285: GEN q, mod, modold;
286: GEN W, Tr, Sr, Wr = gzero, Trold, Spow;
287: int flag, init, ex;
288: GEN *gptr[2];
289: if (DEBUGLEVEL >= 1)
290: timer2();
291: x = varn(gl->T);
292: Tr = gmul(gl->T, gmodulcp(gun, gl->p));
293: W = ginv(gsubst(deriv(Tr, x), x, S));
294: q = gl->p;
295: modold = gl->p;
296: Trold = Tr;
297: ex = 1;
298: flag = 1;
299: init = 0;
300: gptr[0] = &S;
301: gptr[1] = &Wr;
302: while (flag)
303: {
304: if (DEBUGLEVEL >= 1)
305: timer2();
306: q = gsqr(q);
307: ex <<= 1;
308: if (ex >= gl->e)
309: {
310: flag = 0;
311: q = gl->Q;
312: ex = gl->e;
313: }
314: mod = gmodulcp(gun, q);
315: Tr = gmul(gl->T, mod);
316: ltop = avma;
317: Sr = gmodulcp(gmul(lift(lift(S)), mod), Tr);
318: Spow = compoTS(Tr, Sr);
319: if (init)
320: W = gmul(Wr, gsub(gdeux, gmul(Wr, calcderivTS(Spow, Trold, modold))));
321: else
322: init = 1;
323: Wr = gmodulcp(gmul(lift(lift(W)), mod), Tr);
324: S = gmul(Wr, calcTS(Spow, Tr, Sr));
325: lbot = avma;
326: Wr = gcopy(Wr);
327: S = gsub(Sr, S);
328: gerepilemanysp(ltop, lbot, gptr, 2);
329: modold = mod;
330: Trold = Tr;
331: }
332: if (DEBUGLEVEL >= 1)
333: msgtimer("automorphismlift()");
334: return S;
335: }
336: struct galois_testlift
337: {
338: long n;
339: long f;
340: long g;
341: GEN bezoutcoeff;
342: GEN pauto;
343: };
344: /*
345: * Si polb est un polynomes de Z[X] et pola un facteur modulo p, retourne b*v
346: * telqu' il existe u et a tel que: a=pola [mod p] a*b=polb [mod p^e]
347: * b*v+a*u=1 [mod p^e]
348: */
349: GEN
350: bezout_lift_fact(GEN pola, GEN polb, GEN p, GEN pev, long e)
351: {
352: long ev;
353: GEN ae, be, u, v, ae2, be2, s, t, pe, pe2, z, g;
354: long ltop = avma, lbot;
355: if (DEBUGLEVEL >= 1)
356: timer2();
357: ae = lift(pola);
358: be = Fp_poldivres(polb, ae, p, NULL);
359: g = (GEN) Fp_pol_extgcd(ae, be, p, &u, &v)[2]; /* deg g = 0 */
360: if (!gcmp1(g))
361: {
362: g = mpinvmod(g, p);
363: u = gmul(u, g);
364: v = gmul(v, g);
365: }
366: for (pe = p, ev = 1; ev < e;)
367: {
368: ev <<= 1;
369: pe2 = (ev >= e) ? pev : sqri(pe);
370: g = gadd(polb, gneg_i(gmul(ae, be)));
371: g = Fp_pol_red(g, pe2);
372: g = gdivexact(g, pe);
373: z = Fp_pol_red(gmul(v, g), pe);
374: t = Fp_poldivres(z, ae, pe, &s);
375: t = gadd(gmul(u, g), gmul(t, be));
376: t = Fp_pol_red(t, pe);
377: be2 = gadd(be, gmul(t, pe));
378: ae2 = gadd(ae, gmul(s, pe));/* already reduced mod pe2 */
379: g = gadd(gun, gneg_i(gadd(gmul(u, ae2), gmul(v, be2))));
380: g = Fp_pol_red(g, pe2);
381: g = gdivexact(g, pe);
382: z = Fp_pol_red(gmul(v, g), pe);
383: t = Fp_poldivres(z, ae, pe, &s);
384: t = gadd(gmul(u, g), gmul(t, be));
385: t = Fp_pol_red(t, pe);
386: u = gadd(u, gmul(t, pe));
387: v = gadd(v, gmul(s, pe));
388: pe = pe2;
389: ae = ae2;
390: be = be2;
391: }
392: lbot = avma;
393: g = gmul(v, be);
394: g = gerepile(ltop, lbot, g);
395: if (DEBUGLEVEL >= 1)
396: msgtimer("bezout_lift_fact()");
397: return g;
398: }
399: static long
400: inittestlift(GEN Tmod, long elift, struct galois_lift * gl, struct galois_testlift * gt, GEN frob)
401: {
402: ulong ltop = avma;
403: int i, j;
404: long v;
405: GEN pe, autpow, plift;
406: GEN Tmodp, xmodp, modQ, TmodQ, xmodQ;
407: GEN *gptr[2];
408: v = varn(gl->T);
409: gt->n = lg(gl->L) - 1;
410: gt->g = lg(Tmod) - 1;
411: gt->f = gt->n / gt->g;
412: Tmodp = gmul(gl->T, gmodulcp(gun, gl->p));
413: xmodp = gmodulcp(gmul(polx[v], gmodulcp(gun, gl->p)), Tmodp);
414: pe = gpowgs(gl->p, elift);
415: plift = automorphismlift(powgi(xmodp, pe), gl);
416: if (frob)
417: {
418: GEN tlift;
419: tlift = gdiv(centerlift(gmul((GEN) plift[2], gl->den)), gl->den);
420: if (poltopermtest(tlift, gl->L, frob))
421: {
422: avma = ltop;
423: return 1;
424: }
425: }
426: modQ = gmodulcp(gun, gl->Q);
427: TmodQ = gmul(gl->T, modQ);
428: xmodQ = gmodulcp(gmul(polx[v], modQ), TmodQ);
429: gt->bezoutcoeff = cgetg(gt->g + 1, t_VEC);
430: for (i = 1; i <= gt->g; i++)
431: {
432: GEN blift;
433: blift = bezout_lift_fact((GEN) Tmod[i], gl->T, gl->p, gl->Q, gl->e);
434: gt->bezoutcoeff[i] = (long) gmodulcp(gmul(blift, modQ), TmodQ);
435: }
436: gt->pauto = cgetg(gt->f + 1, t_VEC);
437: gt->pauto[1] = (long) xmodQ;
438: gt->pauto[2] = (long) plift;
439: if (gt->f > 2)
440: {
441: autpow = cgetg(gt->n, t_VEC);
442: autpow[1] = (long) plift;
443: for (i = 2; i < gt->n; i++)
444: autpow[i] = lmul((GEN) autpow[i - 1], plift);
445: for (i = 3; i <= gt->f; i++)
446: {
447: GEN s, P;
448: P = ((GEN **) gt->pauto)[i - 1][2];
449: s = (GEN) P[2];
450: for (j = 1; j < lgef(P) - 2; j++)
451: s = gadd(s, gmul((GEN) autpow[j], (GEN) P[j + 2]));
452: gt->pauto[i] = (long) s;
453: }
454: }
455: gptr[0] = >->bezoutcoeff;
456: gptr[1] = >->pauto;
457: gerepilemany(ltop, gptr, 2);
458: return 0;
459: }
460: /*
461: *
462: */
463: long
464: frobeniusliftall(GEN sg, GEN * psi, struct galois_lift * gl, struct galois_testlift * gt, GEN frob)
465: {
466: ulong ltop = avma, av, ltop2;
467: long d, N, z, m, c;
468: int i, j, k;
469: GEN pf, u, v;
470: m = gt->g;
471: c = lg(sg) - 1;
472: d = m / c;
473: pf = cgetg(m, t_VECSMALL);
474: *psi = pf;
475: ltop2 = avma;
476: N = itos(gdiv(mpfact(m), gmul(stoi(c), gpowgs(mpfact(d), c))));
477: avma = ltop2;
478: v = gmul((GEN) gt->pauto[2], (GEN) gt->bezoutcoeff[m]);
479: for (i = 1; i < m; i++)
480: pf[i] = 1 + i / d;
481: av = avma;
482: for (i = 0;; i++)
483: {
484: u = v;
485: if (DEBUGLEVEL >= 4)
486: {
487: fprintferr("GaloisConj:Testing %Z:%d:%Z:", sg, i, pf);
488: timer2();
489: }
490: for (j = 1; j < m; j++)
491: u = gadd(u, gmul((GEN) gt->pauto[sg[pf[j]] + 1], (GEN) gt->bezoutcoeff[j]));
492: u = gdiv(centerlift(gmul((GEN) u[2], gl->den)), gl->den);
493: if (poltopermtest(u, gl->L, frob))
494: {
495: if (DEBUGLEVEL >= 4)
496: msgtimer("");
497: avma = ltop2;
498: return 1;
499: }
500: if (DEBUGLEVEL >= 4)
501: msgtimer("");
502: avma = av;
503: if (i == N - 1)
504: break;
505: for (j = 2; j < m && pf[j - 1] >= pf[j]; j++);
506: for (k = 1; k < j - k && pf[k] != pf[j - k]; k++)
507: {
508: z = pf[k];
509: pf[k] = pf[j - k];
510: pf[j - k] = z;
511: }
512: for (k = j - 1; pf[k] >= pf[j]; k--);
513: z = pf[j];
514: pf[j] = pf[k];
515: pf[k] = z;
516: }
517: avma = ltop;
518: *psi = NULL;
519: return 0;
520: }
521: /*
522: * applique une permutation t a un vecteur s sans duplication
523: */
524: static GEN
525: applyperm(GEN s, GEN t)
526: {
527: GEN u;
528: int i;
529: if (lg(s) < lg(t))
530: err(talker, "First permutation shorter than second in applyperm");
531: u = cgetg(lg(s), typ(s));
532: for (i = 1; i < lg(s); i++)
533: u[i] = s[t[i]];
534: return u;
535: }
536: /*
537: * alloue une ardoise pour n entiers de longueurs pour les test de
538: * permutation
539: */
540: GEN
541: alloue_ardoise(long n, long s)
542: {
543: GEN ar;
544: long i;
545: ar = cgetg(n + 1, t_VEC);
546: for (i = 1; i <= n; i++)
547: ar[i] = lgetg(s, t_INT);
548: return ar;
549: }
550: /*
551: * structure contenant toutes les données pour le tests des permutations:
552: *
553: * ordre :ordre des tests pour verifie_test ordre[lg(ordre)]: numero du test
554: * principal borne : borne sur les coefficients a trouver ladic: modulo
555: * l-adique des racines lborne:ladic-borne TM:vecteur des ligne de M
556: * PV:vecteur des clones des matrices de test (Vmatrix) (ou NULL si non
557: * calcule) L,M comme d'habitude (voir plus bas)
558: */
559: struct galois_test
560: {
561: GEN ordre;
562: GEN borne, lborne, ladic;
563: GEN PV, TM;
564: GEN L, M;
565: };
566: /* Calcule la matrice de tests correspondant a la n-ieme ligne de V */
567: static GEN
568: Vmatrix(long n, struct galois_test * td)
569: {
570: ulong ltop = avma, lbot;
571: GEN V;
572: long i;
573: V = cgetg(lg(td->L), t_VEC);
574: for (i = 1; i < lg(V); i++)
575: V[i] = ((GEN **) td->M)[i][n][2];
576: V = centerlift(gmul(td->L, V));
577: lbot = avma;
578: V = gmod(V, td->ladic);
579: return gerepile(ltop, lbot, V);
580: }
581: /*
582: * Initialise la structure galois_test
583: */
584: static void
585: inittest(GEN L, GEN M, GEN borne, GEN ladic, struct galois_test * td)
586: {
587: ulong ltop;
588: long i;
589: int n = lg(L) - 1;
590: if (DEBUGLEVEL >= 8)
591: fprintferr("GaloisConj:Entree Init Test\n");
592: td->ordre = cgetg(n + 1, t_VECSMALL);
593: for (i = 1; i <= n - 2; i++)
594: td->ordre[i] = i + 2;
595: for (; i <= n; i++)
596: td->ordre[i] = i - n + 2;
597: td->borne = borne;
598: td->lborne = gsub(ladic, borne);
599: td->ladic = ladic;
600: td->L = L;
601: td->M = M;
602: td->PV = cgetg(n + 1, t_VECSMALL); /* peut-etre t_VEC est correct ici */
603: for (i = 1; i <= n; i++)
604: td->PV[i] = 0;
605: ltop = avma;
606: td->PV[td->ordre[n]] = (long) gclone(Vmatrix(td->ordre[n], td));
607: avma = ltop;
608: td->TM = gtrans(M);
609: settyp(td->TM, t_VEC);
610: for (i = 1; i < lg(td->TM); i++)
611: settyp(td->TM[i], t_VEC);
612: if (DEBUGLEVEL >= 8)
613: fprintferr("GaloisConj:Sortie Init Test\n");
614: }
615: /*
616: * liberer les clones de la structure galois_test
617: *
618: * Reservé a l'accreditation ultra-violet:Liberez les clones! Paranoia(tm)
619: */
620: static void
621: freetest(struct galois_test * td)
622: {
623: int i;
624: for (i = 1; i < lg(td->PV); i++)
625: if (td->PV[i])
626: {
627: gunclone((GEN) td->PV[i]);
628: td->PV[i] = 0;
629: }
630: }
631: /*
632: * Test si le nombre padique P est proche d'un entier inferieur a td->borne
633: * en valeur absolue.
634: */
635: long
636: padicisint(GEN P, struct galois_test * td)
637: {
638: long r;
639: ulong ltop = avma;
640: GEN U;
641: if (typ(P) != t_INT)
642: err(typeer, "padicisint");
643: U = modii(P, td->ladic);
644: r = gcmp(U, (GEN) td->borne) <= 0 || gcmp(U, (GEN) td->lborne) >= 0;
645: avma = ltop;
646: return r;
647: }
648: /*
649: * Verifie si pf est une vrai solution et non pas un "hop"
650: */
651: static long
652: verifietest(GEN pf, struct galois_test * td)
653: {
654: ulong av = avma;
655: GEN P, V;
656: int i, j;
657: int n = lg(td->L) - 1;
658: if (DEBUGLEVEL >= 8)
659: fprintferr("GaloisConj:Entree Verifie Test\n");
660: P = applyperm(td->L, pf);
661: for (i = 1; i < n; i++)
662: {
663: long ord;
664: GEN PW;
665: ord = td->ordre[i];
666: PW = (GEN) td->PV[ord];
667: if (PW)
668: {
669: V = ((GEN **) PW)[1][pf[1]];
670: for (j = 2; j <= n; j++)
671: V = gadd(V, ((GEN **) PW)[j][pf[j]]);
672: } else
673: V = centerlift(gmul((GEN) td->TM[ord], P));
674: if (!padicisint(V, td))
675: break;
676: }
677: if (i == n)
678: {
679: if (DEBUGLEVEL >= 8)
680: fprintferr("GaloisConj:Sortie Verifie Test:1\n");
681: avma = av;
682: return 1;
683: }
684: if (!td->PV[td->ordre[i]])
685: {
686: td->PV[td->ordre[i]] = (long) gclone(Vmatrix(td->ordre[i], td));
687: if (DEBUGLEVEL >= 2)
688: fprintferr("M");
689: }
690: if (DEBUGLEVEL >= 2)
691: fprintferr("%d.", i);
692: if (i > 1)
693: {
694: long z;
695: z = td->ordre[i];
696: for (j = i; j > 1; j--)
697: td->ordre[j] = td->ordre[j - 1];
698: td->ordre[1] = z;
699: if (DEBUGLEVEL >= 6)
700: fprintferr("%Z", td->ordre);
701: }
702: if (DEBUGLEVEL >= 8)
703: fprintferr("GaloisConj:Sortie Verifie Test:0\n");
704: avma = av;
705: return 0;
706: }
707: /*
708: * F est la decomposition en cycle de sigma, B est la decomposition en cycle
709: * de cl(tau) Teste toutes les permutations pf pouvant correspondre a tau tel
710: * que : tau*sigma*tau^-1=sigma^s et tau^d=sigma^t ou d est l'ordre de
711: * cl(tau) --------- x: vecteur des choix y: vecteur des mises a jour G:
712: * vecteur d'acces a F linéaire
713: */
714: GEN
715: testpermutation(GEN F, GEN B, long s, long t, long cut, struct galois_test * td)
716: {
717: ulong av, avm = avma, av2;
718: long a, b, c, d, n;
719: GEN pf, x, ar, y, *G;
720: int N, cx, p1, p2, p3, p4, p5, p6;
721: int i, j, hop = 0;
722: GEN V, W;
723: if (DEBUGLEVEL >= 1)
724: timer2();
725: a = lg(F) - 1;
726: b = lg(F[1]) - 1;
727: c = lg(B) - 1;
728: d = lg(B[1]) - 1;
729: n = a * b;
730: s = (b + s) % b;
731: pf = cgetg(n + 1, t_VECSMALL);
732: av = avma;
733: ar = alloue_ardoise(a, 1 + lg(td->ladic));
734: x = cgetg(a + 1, t_VECSMALL);
735: y = cgetg(a + 1, t_VECSMALL);
736: G = (GEN *) cgetg(a + 1, t_VECSMALL); /* Don't worry */
737: av2 = avma;
738: W = (GEN) td->PV[td->ordre[n]];
739: for (cx = 1, i = 1, j = 1; cx <= a; cx++, i++)
740: {
741: x[cx] = (i != d) ? 0 : t;
742: y[cx] = 1;
743: G[cx] = (GEN) F[((long **) B)[j][i]]; /* Be happy */
744: if (i == d)
745: {
746: i = 0;
747: j++;
748: }
749: } /* Be happy now! */
750: N = itos(gpowgs(stoi(b), c * (d - 1))) / cut;
751: avma = av2;
752: if (DEBUGLEVEL >= 4)
753: fprintferr("GaloisConj:I will try %d permutations\n", N);
754: for (cx = 0; cx < N; cx++)
755: {
756: if (DEBUGLEVEL >= 2 && cx % 1000 == 999)
757: fprintferr("%d%% ", (100 * cx) / N);
758: if (cx)
759: {
760: for (i = 1, j = d; i < a;)
761: {
762: y[i] = 1;
763: if ((++(x[i])) != b)
764: break;
765: x[i++] = 0;
766: if (i == j)
767: {
768: y[i++] = 1;
769: j += d;
770: }
771: }
772: y[i + 1] = 1;
773: }
774: for (p1 = 1, p5 = d; p1 <= a; p1++, p5++)
775: if (y[p1])
776: {
777: if (p5 == d)
778: {
779: p5 = 0;
780: p4 = 0;
781: } else
782: p4 = x[p1 - 1];
783: if (p5 == d - 1)
784: p6 = p1 + 1 - d;
785: else
786: p6 = p1 + 1;
787: y[p1] = 0;
788: V = gzero;
789: for (p2 = 1 + p4, p3 = 1 + x[p1]; p2 <= b; p2++)
790: {
791: V = gadd(V, ((GEN **) W)[G[p6][p3]][G[p1][p2]]);
792: p3 += s;
793: if (p3 > b)
794: p3 -= b;
795: }
796: p3 = 1 + x[p1] - s;
797: if (p3 <= 0)
798: p3 += b;
799: for (p2 = p4; p2 >= 1; p2--)
800: {
801: V = gadd(V, ((GEN **) W)[G[p6][p3]][G[p1][p2]]);
802: p3 -= s;
803: if (p3 <= 0)
804: p3 += b;
805: }
806: gaffect((GEN) V, (GEN) ar[p1]);
807: }
808: V = gzero;
809: for (p1 = 1; p1 <= a; p1++)
810: V = gadd(V, (GEN) ar[p1]);
811: if (padicisint(V, td))
812: {
813: for (p1 = 1, p5 = d; p1 <= a; p1++, p5++)
814: {
815: if (p5 == d)
816: {
817: p5 = 0;
818: p4 = 0;
819: } else
820: p4 = x[p1 - 1];
821: if (p5 == d - 1)
822: p6 = p1 + 1 - d;
823: else
824: p6 = p1 + 1;
825: for (p2 = 1 + p4, p3 = 1 + x[p1]; p2 <= b; p2++)
826: {
827: pf[G[p1][p2]] = G[p6][p3];
828: p3 += s;
829: if (p3 > b)
830: p3 -= b;
831: }
832: p3 = 1 + x[p1] - s;
833: if (p3 <= 0)
834: p3 += b;
835: for (p2 = p4; p2 >= 1; p2--)
836: {
837: pf[G[p1][p2]] = G[p6][p3];
838: p3 -= s;
839: if (p3 <= 0)
840: p3 += b;
841: }
842: }
843: if (verifietest(pf, td))
844: {
845: avma = av;
846: if (DEBUGLEVEL >= 1)
847: msgtimer("testpermutation(%d)", cx + 1);
848: if (DEBUGLEVEL >= 2 && hop)
849: fprintferr("GaloisConj:%d hop sur %d iterations\n", hop, cx + 1);
850: return pf;
851: } else
852: hop++;
853: }
854: avma = av2;
855: }
856: avma = avm;
857: if (DEBUGLEVEL >= 1)
858: msgtimer("testpermutation(%d)", N);
859: if (DEBUGLEVEL >= 2 && hop)
860: fprintferr("GaloisConj:%d hop sur %d iterations\n", hop, N);
861: return gzero;
862: }
863: /*
864: * Calcule la liste des sous groupes de \ZZ/m\ZZ d'ordre divisant p et
865: * retourne la liste de leurs elements
866: */
867: GEN
868: listsousgroupes(long m, long p)
869: {
870: ulong ltop = avma, lbot;
871: GEN zns, lss, res, sg, gen, ord, res1;
872: int k, card, i, j, l, n, phi, h;
873: if (m == 2)
874: {
875: res = cgetg(2, t_VEC);
876: sg = cgetg(2, t_VECSMALL);
877: res[1] = (long) sg;
878: sg[1] = 1;
879: return res;
880: }
881: zns = znstar(stoi(m));
882: phi = itos((GEN) zns[1]);
883: lss = subgrouplist((GEN) zns[2], 0);
884: gen = cgetg(lg(zns[3]), t_VECSMALL);
885: ord = cgetg(lg(zns[3]), t_VECSMALL);
886: res1 = cgetg(lg(lss), t_VECSMALL);
887: lbot = avma;
888: for (k = 1, i = lg(lss) - 1; i >= 1; i--)
889: {
890: long av;
891: av = avma;
892: card = phi / itos(det((GEN) lss[i]));
893: avma = av;
894: if (p % card == 0)
895: {
896: sg = cgetg(1 + card, t_VECSMALL);
897: sg[1] = 1;
898: av = avma;
899: for (j = 1; j < lg(gen); j++)
900: {
901: gen[j] = 1;
902: for (h = 1; h < lg(lss[i]); h++)
903: gen[j] = (gen[j] * itos(lift(powgi(((GEN **) zns)[3][h],
904: ((GEN ***) lss)[i][j][h])))) % m;
905: ord[j] = itos(((GEN **) zns)[2][j]) / itos(((GEN ***) lss)[i][j][j]);
906: }
907: avma = av;
908: for (j = 1, l = 1; j < lg(gen); j++)
909: {
910: int c = l * (ord[j] - 1);
911: for (n = 1; n <= c; n++)/* I like it */
912: sg[++l] = (sg[n] * gen[j]) % m;
913: }
914: res1[k++] = (long) sg;
915: }
916: }
917: res = cgetg(k, t_VEC);
918: for (i = 1; i < k; i++)
919: res[i] = res1[i];
920: return gerepile(ltop, lbot, res);
921: }
922: /* retourne la permutation identite */
923: GEN
924: permidentity(long l)
925: {
926: GEN perm;
927: int i;
928: perm = cgetg(l + 1, t_VECSMALL);
929: for (i = 1; i <= l; i++)
930: perm[i] = i;
931: return perm;
932: }
933: /* retourne la decomposition en cycle */
934: GEN
935: permtocycle(GEN p)
936: {
937: long ltop = avma, lbot;
938: int i, j, k, l, m;
939: GEN bit, cycle, cy;
940: cycle = cgetg(lg(p), t_VEC);
941: bit = cgetg(lg(p), t_VECSMALL);
942: for (i = 1; i < lg(p); i++)
943: bit[i] = 0;
944: for (k = 1, l = 1; k < lg(p);)
945: {
946: for (j = 1; bit[j]; j++);
947: cy = cgetg(lg(p), t_VECSMALL);
948: m = 1;
949: do
950: {
951: bit[j] = 1;
952: k++;
953: cy[m++] = j;
954: j = p[j];
955: } while (!bit[j]);
956: setlg(cy, m);
957: cycle[l++] = (long) cy;
958: }
959: setlg(cycle, l);
960: lbot = avma;
961: cycle = gcopy(cycle);
962: return gerepile(ltop, lbot, cycle);
963: }
964: /* Calcule les orbites d'un ensemble de permutations */
965: GEN
966: vecpermorbite(GEN v)
967: {
968: ulong ltop = avma, lbot;
969: int i, j, k, l, m, n, o, p, flag;
970: GEN bit, cycle, cy;
971: n = lg(v[1]);
972: cycle = cgetg(n, t_VEC);
973: bit = cgetg(n, t_VECSMALL);
974: for (i = 1; i < n; i++)
975: bit[i] = 0;
976: for (k = 1, l = 1; k < n;)
977: {
978: for (j = 1; bit[j]; j++);
979: cy = cgetg(n, t_VECSMALL);
980: m = 1;
981: cy[m++] = j;
982: bit[j] = 1;
983: k++;
984: do
985: {
986: flag = 0;
987: for (o = 1; o < lg(v); o++)
988: {
989: for (p = 1; p < m; p++) /* m varie! */
990: {
991: j = ((long **) v)[o][cy[p]];
992: if (!bit[j])
993: {
994: flag = 1;
995: bit[j] = 1;
996: k++;
997: cy[m++] = j;
998: }
999: }
1000: }
1001: } while (flag);
1002: setlg(cy, m);
1003: cycle[l++] = (long) cy;
1004: }
1005: setlg(cycle, l);
1006: lbot = avma;
1007: cycle = gcopy(cycle);
1008: return gerepile(ltop, lbot, cycle);
1009: }
1010: /*
1011: * Calcule la permutation associe a un polynome f des racines L
1012: */
1013: GEN
1014: poltoperm(GEN f, GEN L)
1015: {
1016: ulong ltop, ltop2;
1017: GEN pf, fx, fp;
1018: int i, j;
1019: pf = cgetg(lg(L), t_VECSMALL);
1020: ltop = avma;
1021: fp = cgetg(lg(L), t_VECSMALL);
1022: ltop2 = avma;
1023: for (i = 1; i < lg(fp); i++)
1024: fp[i] = 1;
1025: for (i = 1; i < lg(L); i++)
1026: {
1027: fx = gsubst(f, varn(f), (GEN) L[i]);
1028: for (j = 1; j < lg(L); j++)
1029: if (fp[j] && gegal(fx, (GEN) L[j]))
1030: {
1031: pf[i] = j;
1032: fp[j] = 0;
1033: break;
1034: }
1035: avma = ltop2;
1036: }
1037: avma = ltop;
1038: return pf;
1039: }
1040: /*
1041: * Calcule un polynome R definissant le corps fixe de T pour les orbites O
1042: * tel que R est sans-facteur carre modulo mod et l retourne dans U les
1043: * racines de R
1044: */
1045: GEN
1046: corpsfixeorbitemod(GEN L, GEN O, long x, GEN mod, GEN l, GEN * U)
1047: {
1048: ulong ltop = avma, lbot, av, av2;
1049: GEN g, p, PL, *gptr[2], gmod, gl, modl;
1050: GEN id;
1051: int i, j, d, dmax;
1052: PL = cgetg(lg(O), t_COL);
1053: modl = gmodulcp(gun, l);
1054: av2 = avma;
1055: dmax = lg(L) + 1;
1056: d = 0;
1057: do
1058: {
1059: avma = av2;
1060: id = stoi(d);
1061: g = polun[x];
1062: for (i = 1; i < lg(O); i++)
1063: {
1064: p = gadd(id, (GEN) L[((GEN *) O)[i][1]]);
1065: for (j = 2; j < lg(O[i]); j++)
1066: p = gmul(p, gadd(id, (GEN) L[((GEN *) O)[i][j]]));
1067: PL[i] = (long) p;
1068: g = gmul(g, gsub(polx[x], p));
1069: }
1070: lbot = avma;
1071: g = centerlift(g);
1072: av = avma;
1073: gmod = gmul(g, mod);
1074: gl = gmul(g, modl);
1075: if (DEBUGLEVEL >= 6)
1076: fprintferr("GaloisConj:corps fixe:%d:%Z\n", d, g);
1077: d = (d > 0 ? -d : 1 - d);
1078: if (d > dmax)
1079: err(talker, "prime too small in corpsfixeorbitemod");
1080: } while (degree(ggcd(deriv(gl, x), gl)) || degree(ggcd(deriv(gmod, x), gmod)));
1081: avma = av;
1082: *U = gcopy(PL);
1083: gptr[0] = &g;
1084: gptr[1] = U;
1085: gerepilemanysp(ltop, lbot, gptr, 2);
1086: return g;
1087: }
1088: /*
1089: * Calcule l'inclusion de R dans T i.e. S telque T|R\compo S
1090: */
1091: GEN
1092: corpsfixeinclusion(GEN O, GEN PL)
1093: {
1094: GEN S;
1095: int i, j;
1096: S = cgetg((lg(O) - 1) * (lg(O[1]) - 1) + 1, t_COL);
1097: for (i = 1; i < lg(O); i++)
1098: for (j = 1; j < lg(O[i]); j++)
1099: S[((GEN *) O)[i][j]] = lcopy((GEN) PL[i]);
1100: return S;
1101: }
1102: /* Calcule l'inverse de la matrice de van der Monde de T multiplie par den */
1103: GEN
1104: vandermondeinverse(GEN L, GEN T, GEN den)
1105: {
1106: ulong ltop = avma, lbot;
1107: int i, j, n = lg(L);
1108: long x = varn(T);
1109: GEN M, P, Tp;
1110: if (DEBUGLEVEL >= 1)
1111: timer2();
1112: M = cgetg(n, t_MAT);
1113: Tp = deriv(T, x);
1114: for (i = 1; i < n; i++)
1115: {
1116: M[i] = lgetg(n, t_COL);
1117: P = gdiv(gdeuc(T, gsub(polx[x], (GEN) L[i])), gsubst(Tp, x, (GEN) L[i]));
1118: for (j = 1; j < n; j++)
1119: ((GEN *) M)[i][j] = P[1 + j];
1120: }
1121: if (DEBUGLEVEL >= 1)
1122: msgtimer("vandermondeinverse");
1123: lbot = avma;
1124: M = gmul(den, M);
1125: return gerepile(ltop, lbot, M);
1126: }
1127: /* Calcule le polynome associe a un vecteur conjugue v */
1128: static GEN
1129: vectopol(GEN v, GEN L, GEN M, GEN den, long x)
1130: {
1131: ulong ltop = avma, lbot;
1132: GEN res;
1133: res = gmul(M, v);
1134: res = gtopolyrev(centerlift(res), x);
1135: lbot = avma;
1136: res = gdiv(res, den);
1137: return gerepile(ltop, lbot, res);
1138: }
1139: /* Calcule le polynome associe a une permuation p */
1140: static GEN
1141: permtopol(GEN p, GEN L, GEN M, GEN den, long x)
1142: {
1143: ulong ltop = avma, lbot;
1144: GEN res;
1145: res = gmul(M, applyperm(L, p));
1146: res = gtopolyrev(centerlift(res), x);
1147: lbot = avma;
1148: res = gdiv(res, den);
1149: return gerepile(ltop, lbot, res);
1150: }
1151: /*
1152: * factorise partiellement n sous la forme n=d*u*f^2 ou d est sans facteur
1153: * carre et u n'est pas un carre parfait et retourne u*f
1154: */
1155: GEN
1156: mycoredisc(GEN n)
1157: {
1158: ulong av = avma, tetpil, r;
1159: GEN y, p1, p2, ret;
1160: {
1161: long av = avma, tetpil, i;
1162: GEN fa, p1, p2, p3, res = gun, res2 = gun, y;
1163: fa = auxdecomp(n, 0);
1164: p1 = (GEN) fa[1];
1165: p2 = (GEN) fa[2];
1166: for (i = 1; i < lg(p1) - 1; i++)
1167: {
1168: p3 = (GEN) p2[i];
1169: if (mod2(p3))
1170: res = mulii(res, (GEN) p1[i]);
1171: if (!gcmp1(p3))
1172: res2 = mulii(res2, gpui((GEN) p1[i], shifti(p3, -1), 0));
1173: }
1174: p3 = (GEN) p2[i];
1175: if (mod2(p3)) /* impaire: verifions */
1176: {
1177: if (!gcmp1(p3))
1178: res2 = mulii(res2, gpui((GEN) p1[i], shifti(p3, -1), 0));
1179: if (isprime((GEN) p1[i]))
1180: res = mulii(res, (GEN) p1[i]);
1181: else
1182: res2 = mulii(res2, (GEN) p1[i]);
1183: } else /* paire:OK */
1184: res2 = mulii(res2, gpui((GEN) p1[i], shifti(p3, -1), 0));
1185: tetpil = avma;
1186: y = cgetg(3, t_VEC);
1187: y[1] = (long) icopy(res);
1188: y[2] = (long) icopy(res2);
1189: ret = gerepile(av, tetpil, y);
1190: }
1191: p2 = ret;
1192: p1 = (GEN) p2[1];
1193: r = mod4(p1);
1194: if (signe(p1) < 0)
1195: r = 4 - r;
1196: if (r == 1 || r == 4)
1197: return (GEN) p2[2];
1198: tetpil = avma;
1199: y = gmul2n((GEN) p2[2], -1);
1200: return gerepile(av, tetpil, y);
1201: }
1202: /* Calcule la puissance exp d'une permutation decompose en cycle */
1203: GEN
1204: permcyclepow(GEN O, long exp)
1205: {
1206: int j, k, n;
1207: GEN p;
1208: for (n = 0, j = 1; j < lg(O); j++)
1209: n += lg(O[j]) - 1;
1210: p = cgetg(n + 1, t_VECSMALL);
1211: for (j = 1; j < lg(O); j++)
1212: {
1213: n = lg(O[j]) - 1;
1214: for (k = 1; k <= n; k++)
1215: p[((long **) O)[j][k]] = ((long **) O)[j][1 + (k + exp - 1) % n];
1216: }
1217: return p;
1218: }
1219: /*
1220: * Casse l'orbite O en sous-orbite d'ordre premier correspondant a des
1221: * puissance de l'element
1222: */
1223: GEN
1224: splitorbite(GEN O)
1225: {
1226: ulong ltop = avma, lbot;
1227: int i, n;
1228: GEN F, fc, res;
1229: n = lg(O[1]) - 1;
1230: F = factor(stoi(n));
1231: fc = cgetg(lg(F[1]), t_VECSMALL);
1232: for (i = 1; i < lg(fc); i++)
1233: fc[i] = itos(gpow(((GEN **) F)[1][i], ((GEN **) F)[2][i], DEFAULTPREC));
1234: lbot = avma;
1235: res = cgetg(lg(fc), t_VEC);
1236: for (i = 1; i < lg(fc); i++)
1237: {
1238: GEN v;
1239: v = cgetg(3, t_VEC);
1240: res[i] = (long) v;
1241: v[1] = (long) permcyclepow(O, n / fc[i]);
1242: v[2] = (long) stoi(fc[i]);
1243: }
1244: return gerepile(ltop, lbot, res);
1245: }
1246: /*
1247: * structure qui contient l'analyse du groupe de Galois par etude des degres
1248: * de Frobenius:
1249: *
1250: * p: nombre premier a relever deg: degre du lift à effectuer pgp: plus grand
1251: * diviseur premier du degre de T. exception: 1 si le pgp-Sylow est
1252: * non-cyclique. l: un nombre premier tel que T est totalement decompose
1253: * modulo l ppp: plus petit diviseur premier du degre de T. primepointer:
1254: * permet de calculer les premiers suivant p.
1255: */
1256: struct galois_analysis
1257: {
1258: long p;
1259: long deg;
1260: long exception;
1261: long ord;
1262: long l;
1263: long ppp;
1264: long p4;
1265: byteptr primepointer;
1266: };
1267: /* peut-etre on peut accelerer(distinct degree factorisation */
1268: void
1269: galoisanalysis(GEN T, struct galois_analysis * ga, long calcul_l)
1270: {
1271: ulong ltop = avma;
1272: long n, p, ex, plift, nbmax, nbtest, exist6 = 0, p4 = 0;
1273: GEN F, fc;
1274: byteptr primepointer, pp = 0;
1275: long pha, ord, deg, ppp, pgp, ordmax = 0, l = 0;
1276: int i;
1277: /* Etude du cardinal: */
1278: /* Calcul de l'ordre garanti d'un sous-groupe cyclique */
1279: /* Tout les groupes d'ordre n sont cyclique ssi gcd(n,phi(n))==1 */
1280: if (DEBUGLEVEL >= 1)
1281: timer2();
1282: n = degree(T);
1283: F = factor(stoi(n));
1284: fc = cgetg(lg(F[1]), t_VECSMALL);
1285: for (i = 1; i < lg(fc); i++)
1286: fc[i] = itos(gpow(((GEN **) F)[1][i], ((GEN **) F)[2][i], DEFAULTPREC));
1287: ppp = itos(((GEN **) F)[1][1]); /* Plus Petit diviseur Premier */
1288: pgp = itos(((GEN **) F)[1][lg(F[1]) - 1]); /* Plus Grand diviseur
1289: * Premier */
1290: pha = 1;
1291: ord = 1;
1292: ex = 0;
1293: for (i = lg(F[1]) - 1; i > 0; i--)
1294: {
1295: p = itos(((GEN **) F)[1][i]);
1296: if (pha % p != 0)
1297: {
1298: ord *= p;
1299: pha *= p - 1;
1300: } else
1301: {
1302: ex = 1;
1303: break;
1304: }
1305: if (!gcmp1(((GEN **) F)[2][i]))
1306: break;
1307: }
1308: plift = 0;
1309: /* Etude des ordres des Frobenius */
1310: nbmax = max(4, n / 2); /* Nombre maxi d'éléments à tester */
1311: nbtest = 0;
1312: deg = 0;
1313: for (p = 0, primepointer = diffptr; plift == 0 || (nbtest < nbmax && ord != n) || (n == 24 && exist6 == 0 && p4 == 0);)
1314: {
1315: long u, s, c;
1316: long isram;
1317: GEN S;
1318: c = *primepointer++;
1319: if (!c)
1320: err(primer1);
1321: p += c;
1322: if (p <= (n << 1))
1323: continue;
1324: S = (GEN) simplefactmod(T, stoi(p));
1325: isram = 0;
1326: for (i = 1; i < lg(S[2]) && !isram; i++)
1327: if (!gcmp1(((GEN **) S)[2][i]))
1328: isram = 1;
1329: if (isram == 0)
1330: {
1331: s = itos(((GEN **) S)[1][1]);
1332: for (i = 2; i < lg(S[1]) && !isram; i++)
1333: if (itos(((GEN **) S)[1][i]) != s)
1334: {
1335: avma = ltop;
1336: if (DEBUGLEVEL >= 2)
1337: fprintferr("GaloisAnalysis:non Galois for p=%d\n", p);
1338: ga->p = p;
1339: ga->deg = 0;
1340: return; /* Not a Galois polynomial */
1341: }
1342: if (l == 0 && s == 1)
1343: l = p;
1344: nbtest++;
1345: if (nbtest > (3 * nbmax) && (n == 60 || n == 75))
1346: break;
1347: if (s % 6 == 0)
1348: exist6 = 1;
1349: if (p4 == 0 && s == 4)
1350: p4 = p;
1351: if (DEBUGLEVEL >= 6)
1352: fprintferr("GaloisAnalysis:Nbtest=%d,plift=%d,p=%d,s=%d,ord=%d\n", nbtest, plift, p, s, ord);
1353: if (s > ordmax)
1354: ordmax = s;
1355: if (s >= ord) /* Calcul de l'exposant distinguant minimal
1356: * garanti */
1357: {
1358: if (s * ppp == n) /* Tout sous groupe d'indice ppp est distingué */
1359: u = s;
1360: else /* Théorème de structure des groupes
1361: * hyper-résoluble */
1362: {
1363: u = 1;
1364: for (i = lg(fc) - 1; i > 0; i--)
1365: {
1366: if (s % fc[i] == 0)
1367: u *= fc[i];
1368: else
1369: break;
1370: }
1371: }
1372: if (u != 1)
1373: {
1374: if (!ex || s > ord || (s == ord && (plift == 0 || u > deg)))
1375: {
1376: deg = u;
1377: ord = s;
1378: plift = p;
1379: pp = primepointer;
1380: ex = 1;
1381: }
1382: } else if (!ex && (plift == 0 || s > ord))
1383: /*
1384: * J'ai besoin de plus de connaissance sur les p-groupes, surtout
1385: * pour p=2;
1386: */
1387: {
1388: deg = pgp;
1389: ord = s;
1390: plift = p;
1391: pp = primepointer;
1392: ex = 0;
1393: }
1394: }
1395: }
1396: }
1397: if (plift == 0 || ((n == 36 || n == 48) && !exist6) || (n == 56 && ordmax == 7) || (n == 60 && ordmax == 5) || (n == 72 && !exist6) || (n == 80 && ordmax == 5))
1398: {
1399: deg = 0;
1400: err(warner, "Galois group almost certainly not weakly super solvable");
1401: }
1402: avma = ltop;
1403: if (calcul_l)
1404: {
1405: ulong av;
1406: GEN x;
1407: long c;
1408: x = cgetg(3, t_POLMOD);
1409: x[2] = lpolx[varn(T)];
1410: av = avma;
1411: while (l == 0)
1412: {
1413: c = *primepointer++;
1414: if (!c)
1415: err(primer1);
1416: p += c;
1417: x[1] = lmul(T, gmodulss(1, p));
1418: if (gegal(gpuigs(x, p), x))
1419: l = p;
1420: avma = av;
1421: }
1422: avma = ltop;
1423: }
1424: ga->p = plift;
1425: ga->exception = ex;
1426: ga->deg = deg;
1427: ga->ord = ord;
1428: ga->l = l;
1429: ga->primepointer = pp;
1430: ga->ppp = ppp;
1431: ga->p4 = p4;
1432: if (DEBUGLEVEL >= 4)
1433: fprintferr("GaloisAnalysis:p=%d l=%d exc=%d deg=%d ord=%d ppp=%d\n", p, l, ex, deg, ord, ppp);
1434: if (DEBUGLEVEL >= 1)
1435: msgtimer("galoisanalysis()");
1436: }
1437: /* Calcule les bornes sur les coefficients a chercher */
1438: struct galois_borne
1439: {
1440: GEN l;
1441: long valsol;
1442: long valabs;
1443: GEN bornesol;
1444: GEN ladicsol;
1445: GEN ladicabs;
1446: };
1447: /*
1448: * recalcule L avec une precision superieur
1449: */
1450: GEN
1451: recalculeL(GEN T, GEN Tp, GEN L, struct galois_borne * gb, struct galois_borne * ngb)
1452: {
1453: ulong ltop = avma, lbot = avma;
1454: GEN L2, y, z, lad;
1455: long i, val;
1456: if (DEBUGLEVEL >= 1)
1457: timer2();
1458: L2 = cgetg(lg(L), typ(L));
1459: for (i = 1; i < lg(L); i++)
1460: {
1461: ltop = avma;
1462: val = gb->valabs;
1463: lad = gb->ladicabs;
1464: z = (GEN) L[i];
1465: while (val < ngb->valabs)
1466: {
1467: val <<= 1;
1468: if (val >= ngb->valabs)
1469: {
1470: val = ngb->valabs;
1471: lad = ngb->ladicabs;
1472: } else
1473: lad = gsqr(lad);
1474: y = gmodulcp((GEN) z[2], lad);
1475: z = gdiv(poleval(T, y), poleval(Tp, y));
1476: lbot = avma;
1477: z = gsub(y, z);
1478: }
1479: L2[i] = (long) gerepile(ltop, lbot, z);
1480: }
1481: return L2;
1482: }
1483: void
1484: initborne(GEN T, GEN disc, struct galois_borne * gb, long ppp)
1485: {
1486: ulong ltop = avma, lbot, av2;
1487: GEN borne, borneroots, borneabs;
1488: int i, j;
1489: int n;
1490: GEN L, M, z;
1491: L = roots(T, DEFAULTPREC);
1492: n = lg(L) - 1;
1493: for (i = 1; i <= n; i++)
1494: {
1495: z = (GEN) L[i];
1496: if (signe(z[2]))
1497: break;
1498: L[i] = z[1];
1499: }
1500: M = vandermondeinverse(L, gmul(T, realun(DEFAULTPREC)), disc);
1501: borne = gzero;
1502: for (i = 1; i <= n; i++)
1503: {
1504: z = gzero;
1505: for (j = 1; j <= n; j++)
1506: z = gadd(z, gabs(((GEN **) M)[j][i], DEFAULTPREC));
1507: if (gcmp(z, borne) > 0)
1508: borne = z;
1509: }
1510: borneroots = gzero;
1511: for (i = 1; i <= n; i++)
1512: {
1513: z = gabs((GEN) L[i], DEFAULTPREC);
1514: if (gcmp(z, borneroots) > 0)
1515: borneroots = z;
1516: }
1517: borneabs = addsr(1, gpowgs(addsr(n, borneroots), n / ppp));
1518: lbot = avma;
1519: borneroots = addsr(1, gmul(borne, borneroots));
1520: av2 = avma;
1521: borneabs = gmul2n(gmul(borne, borneabs), 4);
1522: gb->valsol = itos(gceil(gdiv(glog(gmul2n(borneroots, 4 + (n >> 1)), DEFAULTPREC), glog(gb->l, DEFAULTPREC))));
1523: if (DEBUGLEVEL >= 4)
1524: fprintferr("GaloisConj:val1=%d\n", gb->valsol);
1525: gb->valabs = max(gb->valsol, itos(gceil(gdiv(glog(borneabs, DEFAULTPREC), glog(gb->l, DEFAULTPREC)))));
1526: if (DEBUGLEVEL >= 4)
1527: fprintferr("GaloisConj:val2=%d\n", gb->valabs);
1528: avma = av2;
1529: gb->bornesol = gerepile(ltop, lbot, borneroots);
1530: gb->ladicsol = gpowgs(gb->l, gb->valsol);
1531: gb->ladicabs = gpowgs(gb->l, gb->valabs);
1532: }
1533: /* Groupe A4 */
1534: GEN
1535: a4galoisgen(GEN T, struct galois_test * td)
1536: {
1537: int ltop = avma, av, av2;
1538: int i, j, k;
1539: int n;
1540: int N, hop = 0;
1541: GEN *ar, **mt; /* tired of casting */
1542: GEN t, u, O;
1543: GEN res, orb, ry;
1544: GEN pft, pfu, pfv;
1545: n = degree(T);
1546: res = cgetg(4, t_VEC);
1547: ry = cgetg(3, t_VEC);
1548: res[1] = (long) ry;
1549: pft = cgetg(n + 1, t_VECSMALL);
1550: ry[1] = (long) pft;
1551: ry[2] = (long) stoi(2);
1552: ry = cgetg(3, t_VEC);
1553: pfu = cgetg(n + 1, t_VECSMALL);
1554: ry[1] = (long) pfu;
1555: ry[2] = (long) stoi(2);
1556: res[2] = (long) ry;
1557: ry = cgetg(3, t_VEC);
1558: pfv = cgetg(n + 1, t_VECSMALL);
1559: ry[1] = (long) pfv;
1560: ry[2] = (long) stoi(3);
1561: res[3] = (long) ry;
1562: av = avma;
1563: ar = (GEN *) alloue_ardoise(n, 1 + lg(td->ladic));
1564: mt = (GEN **) td->PV[td->ordre[n]];
1565: t = cgetg(n + 1, t_VECSMALL) + 1; /* Sorry for this hack */
1566: u = cgetg(n + 1, t_VECSMALL) + 1; /* too lazy to correct */
1567: av2 = avma;
1568: N = itos(gdiv(mpfact(n), mpfact(n >> 1))) >> (n >> 1);
1569: if (DEBUGLEVEL >= 4)
1570: fprintferr("A4GaloisConj:I will test %d permutations\n", N);
1571: avma = av2;
1572: for (i = 0; i < n; i++)
1573: t[i] = i + 1;
1574: for (i = 0; i < N; i++)
1575: {
1576: GEN g;
1577: int a, x, y;
1578: if (i == 0)
1579: {
1580: gaffect(gzero, ar[(n - 2) >> 1]);
1581: for (k = n - 2; k > 2; k -= 2)
1582: gaddz(ar[k >> 1], gadd(mt[k + 1][k + 2], mt[k + 2][k + 1]), ar[(k >> 1) - 1]);
1583: } else
1584: {
1585: x = i;
1586: y = 1;
1587: do
1588: {
1589: hiremainder = 0;
1590: y += 2;
1591: x = divll(x, y);
1592: a = hiremainder;
1593: }
1594: while (!a);
1595: switch (y)
1596: {
1597: case 3:
1598: x = t[2];
1599: if (a == 1)
1600: {
1601: t[2] = t[1];
1602: t[1] = x;
1603: } else
1604: {
1605: t[2] = t[0];
1606: t[0] = x;
1607: }
1608: break;
1609: case 5:
1610: x = t[0];
1611: t[0] = t[2];
1612: t[2] = t[1];
1613: t[1] = x;
1614: x = t[4];
1615: t[4] = t[4 - a];
1616: t[4 - a] = x;
1617: gaddz(ar[2], gadd(mt[t[4]][t[5]], mt[t[5]][t[4]]), ar[1]);
1618: break;
1619: case 7:
1620: x = t[0];
1621: t[0] = t[4];
1622: t[4] = t[3];
1623: t[3] = t[1];
1624: t[1] = t[2];
1625: t[2] = x;
1626: x = t[6];
1627: t[6] = t[6 - a];
1628: t[6 - a] = x;
1629: gaddz(ar[3], gadd(mt[t[6]][t[7]], mt[t[7]][t[6]]), ar[2]);
1630: gaddz(ar[2], gadd(mt[t[4]][t[5]], mt[t[5]][t[4]]), ar[1]);
1631: break;
1632: case 9:
1633: x = t[0];
1634: t[0] = t[6];
1635: t[6] = t[5];
1636: t[5] = t[3];
1637: t[3] = x;
1638: x = t[4];
1639: t[4] = t[1];
1640: t[1] = x;
1641: x = t[8];
1642: t[8] = t[8 - a];
1643: t[8 - a] = x;
1644: gaddz(ar[4], gadd(mt[t[8]][t[9]], mt[t[9]][t[8]]), ar[3]);
1645: gaddz(ar[3], gadd(mt[t[6]][t[7]], mt[t[7]][t[6]]), ar[2]);
1646: gaddz(ar[2], gadd(mt[t[4]][t[5]], mt[t[5]][t[4]]), ar[1]);
1647: break;
1648: default:
1649: y--;
1650: x = t[0];
1651: t[0] = t[2];
1652: t[2] = t[1];
1653: t[1] = x;
1654: for (k = 4; k < y; k += 2)
1655: {
1656: int j;
1657: x = t[k];
1658: for (j = k; j > 0; j--)
1659: t[j] = t[j - 1];
1660: t[0] = x;
1661: }
1662: x = t[y];
1663: t[y] = t[y - a];
1664: t[y - a] = x;
1665: for (k = y; k > 2; k -= 2)
1666: gaddz(ar[k >> 1], gadd(mt[t[k]][t[k + 1]], mt[t[k + 1]][t[k]]), ar[(k >> 1) - 1]);
1667: }
1668: }
1669: g = gadd(ar[1], gadd(gadd(mt[t[0]][t[1]], mt[t[1]][t[0]]),
1670: gadd(mt[t[2]][t[3]], mt[t[3]][t[2]])));
1671: if (padicisint(g, td))
1672: {
1673: for (k = 0; k < n; k += 2)
1674: {
1675: pft[t[k]] = t[k + 1];
1676: pft[t[k + 1]] = t[k];
1677: }
1678: if (verifietest(pft, td))
1679: break;
1680: else
1681: hop++;
1682: }
1683: avma = av2;
1684: }
1685: if (i == N)
1686: {
1687: avma = ltop;
1688: if (DEBUGLEVEL >= 1 && hop)
1689: fprintferr("A4GaloisConj: %d hop sur %d iterations\n", hop, N);
1690: return gzero;
1691: }
1692: if (DEBUGLEVEL >= 1 && hop)
1693: fprintferr("A4GaloisConj: %d hop sur %d iterations\n", hop, N);
1694: N = itos(gdiv(mpfact(n >> 1), mpfact(n >> 2))) >> 1;
1695: avma = av2;
1696: if (DEBUGLEVEL >= 4)
1697: fprintferr("A4GaloisConj:sigma=%Z \n", pft);
1698: for (i = 0; i < N; i++)
1699: {
1700: GEN g;
1701: int a, x, y;
1702: if (i == 0)
1703: {
1704: for (k = 0; k < n; k += 4)
1705: {
1706: u[k + 3] = t[k + 3];
1707: u[k + 2] = t[k + 1];
1708: u[k + 1] = t[k + 2];
1709: u[k] = t[k];
1710: }
1711: } else
1712: {
1713: x = i;
1714: y = -2;
1715: do
1716: {
1717: hiremainder = 0;
1718: y += 4;
1719: x = divll(x, y);
1720: a = hiremainder;
1721: } while (!a);
1722: x = u[2];
1723: u[2] = u[0];
1724: u[0] = x;
1725: switch (y)
1726: {
1727: case 2:
1728: break;
1729: case 6:
1730: x = u[4];
1731: u[4] = u[6];
1732: u[6] = x;
1733: if (a % 2 == 0)
1734: {
1735: a = 4 - a / 2;
1736: x = u[6];
1737: u[6] = u[a];
1738: u[a] = x;
1739: x = u[4];
1740: u[4] = u[a - 2];
1741: u[a - 2] = x;
1742: }
1743: break;
1744: case 10:
1745: x = u[6];
1746: u[6] = u[3];
1747: u[3] = u[2];
1748: u[2] = u[4];
1749: u[4] = u[1];
1750: u[1] = u[0];
1751: u[0] = x;
1752: if (a >= 3)
1753: a += 2;
1754: a = 8 - a;
1755: x = u[10];
1756: u[10] = u[a];
1757: u[a] = x;
1758: x = u[8];
1759: u[8] = u[a - 2];
1760: u[a - 2] = x;
1761: break;
1762: }
1763: }
1764: g = gzero;
1765: for (k = 0; k < n; k += 2)
1766: g = gadd(g, gadd(mt[u[k]][u[k + 1]], mt[u[k + 1]][u[k]]));
1767: if (padicisint(g, td))
1768: {
1769: for (k = 0; k < n; k += 2)
1770: {
1771: pfu[u[k]] = u[k + 1];
1772: pfu[u[k + 1]] = u[k];
1773: }
1774: if (verifietest(pfu, td))
1775: break;
1776: else
1777: hop++;
1778: }
1779: avma = av2;
1780: }
1781: if (i == N)
1782: {
1783: avma = ltop;
1784: return gzero;
1785: }
1786: if (DEBUGLEVEL >= 1 && hop)
1787: fprintferr("A4GaloisConj: %d hop sur %d iterations\n", hop, N);
1788: if (DEBUGLEVEL >= 4)
1789: fprintferr("A4GaloisConj:tau=%Z \n", u);
1790: avma = av2;
1791: orb = cgetg(3, t_VEC);
1792: orb[1] = (long) pft;
1793: orb[2] = (long) pfu;
1794: if (DEBUGLEVEL >= 4)
1795: fprintferr("A4GaloisConj:orb=%Z \n", orb);
1796: O = vecpermorbite(orb);
1797: if (DEBUGLEVEL >= 4)
1798: fprintferr("A4GaloisConj:O=%Z \n", O);
1799: av2 = avma;
1800: for (j = 0; j < 2; j++)
1801: {
1802: pfv[((long **) O)[1][1]] = ((long **) O)[2][1];
1803: pfv[((long **) O)[1][2]] = ((long **) O)[2][3 + j];
1804: pfv[((long **) O)[1][3]] = ((long **) O)[2][4 - (j << 1)];
1805: pfv[((long **) O)[1][4]] = ((long **) O)[2][2 + j];
1806: for (i = 0; i < 4; i++)
1807: {
1808: long x;
1809: GEN g;
1810: switch (i)
1811: {
1812: case 0:
1813: break;
1814: case 1:
1815: x = ((long **) O)[3][1];
1816: ((long **) O)[3][1] = ((long **) O)[3][2];
1817: ((long **) O)[3][2] = x;
1818: x = ((long **) O)[3][3];
1819: ((long **) O)[3][3] = ((long **) O)[3][4];
1820: ((long **) O)[3][4] = x;
1821: break;
1822: case 2:
1823: x = ((long **) O)[3][1];
1824: ((long **) O)[3][1] = ((long **) O)[3][4];
1825: ((long **) O)[3][4] = x;
1826: x = ((long **) O)[3][2];
1827: ((long **) O)[3][2] = ((long **) O)[3][3];
1828: ((long **) O)[3][3] = x;
1829: break;
1830: case 3:
1831: x = ((long **) O)[3][1];
1832: ((long **) O)[3][1] = ((long **) O)[3][2];
1833: ((long **) O)[3][2] = x;
1834: x = ((long **) O)[3][3];
1835: ((long **) O)[3][3] = ((long **) O)[3][4];
1836: ((long **) O)[3][4] = x;
1837: }
1838: pfv[((long **) O)[2][1]] = ((long **) O)[3][1];
1839: pfv[((long **) O)[2][3 + j]] = ((long **) O)[3][4 - j];
1840: pfv[((long **) O)[2][4 - (j << 1)]] = ((long **) O)[3][2 + (j << 1)];
1841: pfv[((long **) O)[2][2 + j]] = ((long **) O)[3][3 - j];
1842: pfv[((long **) O)[3][1]] = ((long **) O)[1][1];
1843: pfv[((long **) O)[3][4 - j]] = ((long **) O)[1][2];
1844: pfv[((long **) O)[3][2 + (j << 1)]] = ((long **) O)[1][3];
1845: pfv[((long **) O)[3][3 - j]] = ((long **) O)[1][4];
1846: g = gzero;
1847: for (k = 1; k <= n; k++)
1848: g = gadd(g, mt[k][pfv[k]]);
1849: if (padicisint(g, td) && verifietest(pfv, td))
1850: {
1851: avma = av;
1852: if (DEBUGLEVEL >= 1)
1853: fprintferr("A4GaloisConj:%d hop sur %d iterations max\n", hop, 10395 + 68);
1854: return res;
1855: } else
1856: hop++;
1857: avma = av2;
1858: }
1859: }
1860: /* Echec? */
1861: avma = ltop;
1862: return gzero;
1863: }
1864: /* Groupe S4 */
1865: GEN
1866: Fpisom(GEN p, GEN T1, GEN T2)
1867: {
1868: ulong ltop = avma, lbot;
1869: GEN T, res;
1870: long v;
1871: if (T1 == T2)
1872: return gmodulcp(polx[varn(T1)], T1);
1873: v = varn(T2);
1874: setvarn(T2, MAXVARN);
1875: T = (GEN) factmod9(T1, p, T2)[1];
1876: setvarn(T2, v);
1877: lbot = avma;
1878: res = gneg(((GEN **) T)[1][2]);
1879: setvarn(res[1], v);
1880: setvarn(res[2], v);
1881: return gerepile(ltop, lbot, res);
1882: }
1883: GEN
1884: Fpinvisom(GEN S)
1885: {
1886: ulong ltop = avma, lbot;
1887: GEN M, U, V;
1888: int n, i, j;
1889: long x;
1890: x = varn(S[1]);
1891: U = gmodulcp(polun[x], (GEN) S[1]);
1892: n = degree((GEN) S[1]);
1893: M = cgetg(n + 1, t_MAT);
1894: for (i = 1; i <= n; i++)
1895: {
1896: int d;
1897: if (i > 1)
1898: U = gmul(U, S);
1899: M[i] = lgetg(n + 1, t_COL);
1900: d = degree((GEN) U[2]);
1901: for (j = 1; j <= d + 1; j++)
1902: ((GEN *) M)[i][j] = ((GEN *) U)[2][1 + j];
1903: for (j = d + 2; j <= n; j++)
1904: ((GEN *) M)[i][j] = zero;
1905: }
1906: V = cgetg(n + 1, t_COL);
1907: for (i = 1; i <= n; i++)
1908: V[i] = zero;
1909: V[2] = un;
1910: V = gauss(M, V);
1911: lbot = avma;
1912: V = gtopolyrev(V, x);
1913: return gerepile(ltop, lbot, V);
1914: }
1915: static void
1916: makelift(GEN u, struct galois_lift * gl, GEN liftpow)
1917: {
1918: int i;
1919: liftpow[1] = (long) automorphismlift(u, gl);
1920: for (i = 2; i < lg(liftpow); i++)
1921: liftpow[i] = lmul((GEN) liftpow[i - 1], (GEN) liftpow[1]);
1922: }
1923: static long
1924: s4test(GEN u, GEN liftpow, struct galois_lift * gl, GEN phi)
1925: {
1926: GEN res;
1927: int i;
1928: if (DEBUGLEVEL >= 6)
1929: timer2();
1930: u = (GEN) u[2];
1931: res = (GEN) u[2];
1932: for (i = 1; i < lgef(u) - 2; i++)
1933: res = gadd(res, gmul((GEN) liftpow[i], (GEN) u[i + 2]));
1934: res = gdiv(centerlift(gmul((GEN) res[2], gl->den)), gl->den);
1935: if (DEBUGLEVEL >= 6)
1936: msgtimer("s4test()");
1937: return poltopermtest(res, gl->L, phi);
1938: }
1939: GEN
1940: s4galoisgen(struct galois_lift * gl)
1941: {
1942: struct galois_testlift gt;
1943: ulong ltop = avma, av, ltop2;
1944: GEN Tmod, isom, isominv, misom;
1945: int i, j;
1946: GEN sg;
1947: GEN sigma, tau, phi;
1948: GEN res, ry;
1949: GEN pj;
1950: GEN p;
1951: GEN bezoutcoeff, pauto, liftpow;
1952: long v;
1953: p = gl->p;
1954: v = varn(gl->T);
1955: res = cgetg(5, t_VEC);
1956: ry = cgetg(3, t_VEC);
1957: res[1] = (long) ry;
1958: ry[1] = lgetg(lg(gl->L), t_VECSMALL);
1959: ry[2] = (long) stoi(2);
1960: ry = cgetg(3, t_VEC);
1961: res[2] = (long) ry;
1962: ry[1] = lgetg(lg(gl->L), t_VECSMALL);
1963: ry[2] = (long) stoi(2);
1964: ry = cgetg(3, t_VEC);
1965: res[3] = (long) ry;
1966: ry[1] = lgetg(lg(gl->L), t_VECSMALL);
1967: ry[2] = (long) stoi(3);
1968: ry = cgetg(3, t_VEC);
1969: res[4] = (long) ry;
1970: ry[1] = lgetg(lg(gl->L), t_VECSMALL);
1971: ry[2] = (long) stoi(2);
1972: ltop2 = avma;
1973: sg = cgetg(7, t_VECSMALL);
1974: pj = cgetg(7, t_VECSMALL);
1975: sigma = cgetg(lg(gl->L), t_VECSMALL);
1976: tau = cgetg(lg(gl->L), t_VECSMALL);
1977: phi = cgetg(lg(gl->L), t_VECSMALL);
1978: for (i = 1; i < lg(sg); i++)
1979: sg[i] = i;
1980: Tmod = (GEN) factmod(gl->T, p)[1];
1981: isom = cgetg(lg(Tmod), t_VEC);
1982: isominv = cgetg(lg(Tmod), t_VEC);
1983: misom = cgetg(lg(Tmod), t_MAT);
1984: inittestlift(Tmod, 1, gl, >, NULL);
1985: bezoutcoeff = gt.bezoutcoeff;
1986: pauto = gt.pauto;
1987: for (i = 1; i < lg(pj); i++)
1988: pj[i] = 0;
1989: for (i = 1; i < lg(isom); i++)
1990: {
1991: misom[i] = lgetg(lg(Tmod), t_COL);
1992: isom[i] = (long) Fpisom(p, (GEN) Tmod[1], (GEN) Tmod[i]);
1993: if (DEBUGLEVEL >= 4)
1994: fprintferr("S4GaloisConj:Computing isomorphisms %d:%Z\n", i, (GEN) isom[i]);
1995: isominv[i] = (long) lift(Fpinvisom((GEN) isom[i]));
1996: }
1997: for (i = 1; i < lg(isom); i++)
1998: for (j = 1; j < lg(isom); j++)
1999: ((GEN **) misom)[i][j] = gsubst((GEN) isominv[i], v, (GEN) isom[j]);
2000: liftpow = cgetg(24, t_VEC);
2001: av = avma;
2002: for (i = 0; i < 3; i++)
2003: {
2004: ulong av2;
2005: GEN u;
2006: int j1, j2, j3;
2007: ulong avm1, avm2;
2008: GEN u1, u2, u3;
2009: if (i)
2010: {
2011: int x;
2012: x = sg[3];
2013: if (i == 1)
2014: {
2015: sg[3] = sg[2];
2016: sg[2] = x;
2017: } else
2018: {
2019: sg[3] = sg[1];
2020: sg[1] = x;
2021: }
2022: }
2023: u = chinois(((GEN **) misom)[sg[1]][sg[2]],
2024: ((GEN **) misom)[sg[2]][sg[1]]);
2025: u = chinois(u, ((GEN **) misom)[sg[3]][sg[4]]);
2026: u = chinois(u, ((GEN **) misom)[sg[4]][sg[3]]);
2027: u = chinois(u, ((GEN **) misom)[sg[5]][sg[6]]);
2028: u = chinois(u, ((GEN **) misom)[sg[6]][sg[5]]);
2029: makelift(u, gl, liftpow);
2030: av2 = avma;
2031: for (j1 = 0; j1 < 4; j1++)
2032: {
2033: u1 = gadd(gmul((GEN) bezoutcoeff[sg[5]], (GEN) pauto[1 + j1]),
2034: gmul((GEN) bezoutcoeff[sg[6]], (GEN) pauto[((-j1) & 3) + 1]));
2035: avm1 = avma;
2036: for (j2 = 0; j2 < 4; j2++)
2037: {
2038: u2 = gadd(u1, gmul((GEN) bezoutcoeff[sg[3]], (GEN) pauto[1 + j2]));
2039: u2 = gadd(u2, gmul((GEN) bezoutcoeff[sg[4]], (GEN) pauto[((-j2) & 3) + 1]));
2040: avm2 = avma;
2041: for (j3 = 0; j3 < 4; j3++)
2042: {
2043: u3 = gadd(u2, gmul((GEN) bezoutcoeff[sg[1]], (GEN) pauto[1 + j3]));
2044: u3 = gadd(u3, gmul((GEN) bezoutcoeff[sg[2]], (GEN) pauto[((-j3) & 3) + 1]));
2045: if (DEBUGLEVEL >= 4)
2046: fprintferr("S4GaloisConj:Testing %d/3:%d/4:%d/4:%d/4:%Z\n", i, j1, j2, j3, sg);
2047: if (s4test(u3, liftpow, gl, sigma))
2048: {
2049: pj[1] = j3;
2050: pj[2] = j2, pj[3] = j1;
2051: goto suites4;
2052: }
2053: avma = avm2;
2054: }
2055: avma = avm1;
2056: }
2057: avma = av2;
2058: }
2059: avma = av;
2060: }
2061: avma = ltop;
2062: return gzero;
2063: suites4:
2064: if (DEBUGLEVEL >= 4)
2065: fprintferr("S4GaloisConj:sigma=%Z\n", sigma);
2066: if (DEBUGLEVEL >= 4)
2067: fprintferr("S4GaloisConj:pj=%Z\n", pj);
2068: avma = av;
2069: for (j = 1; j <= 3; j++)
2070: {
2071: ulong av2;
2072: GEN u;
2073: int w, l;
2074: int z;
2075: z = sg[1];
2076: sg[1] = sg[3];
2077: sg[3] = sg[5];
2078: sg[5] = z;
2079: z = sg[2];
2080: sg[2] = sg[4];
2081: sg[4] = sg[6];
2082: sg[6] = z;
2083: z = pj[1];
2084: pj[1] = pj[2];
2085: pj[2] = pj[3];
2086: pj[3] = z;
2087: for (l = 0; l < 2; l++)
2088: {
2089: u = chinois(((GEN **) misom)[sg[1]][sg[3]],
2090: ((GEN **) misom)[sg[3]][sg[1]]);
2091: u = chinois(u, ((GEN **) misom)[sg[2]][sg[4]]);
2092: u = chinois(u, ((GEN **) misom)[sg[4]][sg[2]]);
2093: u = chinois(u, ((GEN **) misom)[sg[5]][sg[6]]);
2094: u = chinois(u, ((GEN **) misom)[sg[6]][sg[5]]);
2095: makelift(u, gl, liftpow);
2096: av2 = avma;
2097: for (w = 0; w < 4; w += 2)
2098: {
2099: GEN uu;
2100: long av3;
2101: pj[6] = (w + pj[3]) & 3;
2102: uu = gadd(gmul((GEN) bezoutcoeff[sg[5]], (GEN) pauto[(pj[6] & 3) + 1]),
2103: gmul((GEN) bezoutcoeff[sg[6]], (GEN) pauto[((-pj[6]) & 3) + 1]));
2104: av3 = avma;
2105: for (i = 0; i < 4; i++)
2106: {
2107: GEN u;
2108: pj[4] = i;
2109: pj[5] = (i + pj[2] - pj[1]) & 3;
2110: if (DEBUGLEVEL >= 4)
2111: fprintferr("S4GaloisConj:Testing %d/3:%d/2:%d/2:%d/4:%Z:%Z\n", j - 1, w >> 1, l, i, sg, pj);
2112: u = gadd(uu, gmul((GEN) pauto[(pj[4] & 3) + 1], (GEN) bezoutcoeff[sg[1]]));
2113: u = gadd(u, gmul((GEN) pauto[((-pj[4]) & 3) + 1], (GEN) bezoutcoeff[sg[3]]));
2114: u = gadd(u, gmul((GEN) pauto[(pj[5] & 3) + 1], (GEN) bezoutcoeff[sg[2]]));
2115: u = gadd(u, gmul((GEN) pauto[((-pj[5]) & 3) + 1], (GEN) bezoutcoeff[sg[4]]));
2116: if (s4test(u, liftpow, gl, tau))
2117: goto suites4_2;
2118: avma = av3;
2119: }
2120: avma = av2;
2121: }
2122: z = sg[4];
2123: sg[4] = sg[3];
2124: sg[3] = z;
2125: pj[2] = (-pj[2]) & 3;
2126: avma = av;
2127: }
2128: }
2129: avma = ltop;
2130: return gzero;
2131: suites4_2:
2132: avma = av;
2133: {
2134: int abc, abcdef;
2135: GEN u;
2136: ulong av2;
2137: abc = (pj[1] + pj[2] + pj[3]) & 3;
2138: abcdef = (((abc + pj[4] + pj[5] - pj[6]) & 3) >> 1);
2139: u = chinois(((GEN **) misom)[sg[1]][sg[4]],
2140: ((GEN **) misom)[sg[4]][sg[1]]);
2141: u = chinois(u, ((GEN **) misom)[sg[2]][sg[5]]);
2142: u = chinois(u, ((GEN **) misom)[sg[5]][sg[2]]);
2143: u = chinois(u, ((GEN **) misom)[sg[3]][sg[6]]);
2144: u = chinois(u, ((GEN **) misom)[sg[6]][sg[3]]);
2145: makelift(u, gl, liftpow);
2146: av2 = avma;
2147: for (j = 0; j < 8; j++)
2148: {
2149: int h, g, i;
2150: h = j & 3;
2151: g = abcdef + ((j & 4) >> 1);
2152: i = h + abc - g;
2153: u = gmul((GEN) pauto[(g & 3) + 1], (GEN) bezoutcoeff[sg[1]]);
2154: u = gadd(u, gmul((GEN) pauto[((-g) & 3) + 1], (GEN) bezoutcoeff[sg[4]]));
2155: u = gadd(u, gmul((GEN) pauto[(h & 3) + 1], (GEN) bezoutcoeff[sg[2]]));
2156: u = gadd(u, gmul((GEN) pauto[((-h) & 3) + 1], (GEN) bezoutcoeff[sg[5]]));
2157: u = gadd(u, gmul((GEN) pauto[(i & 3) + 1], (GEN) bezoutcoeff[sg[3]]));
2158: u = gadd(u, gmul((GEN) pauto[((-i) & 3) + 1], (GEN) bezoutcoeff[sg[6]]));
2159: if (DEBUGLEVEL >= 4)
2160: fprintferr("S4GaloisConj:Testing %d/8 %d:%d:%d\n", j, g & 3, h & 3, i & 3);
2161: if (s4test(u, liftpow, gl, phi))
2162: break;
2163: avma = av2;
2164: }
2165: }
2166: if (j == 8)
2167: {
2168: avma = ltop;
2169: return gzero;
2170: }
2171: for (i = 1; i < lg(gl->L); i++)
2172: {
2173: ((GEN **) res)[1][1][i] = sigma[tau[i]];
2174: ((GEN **) res)[2][1][i] = phi[sigma[tau[phi[i]]]];
2175: ((GEN **) res)[3][1][i] = phi[sigma[i]];
2176: ((GEN **) res)[4][1][i] = sigma[i];
2177: }
2178: avma = ltop2;
2179: return res;
2180: }
2181: GEN
2182: galoisgen(GEN T, GEN L, GEN M, GEN den, struct galois_borne * gb, const struct galois_analysis * ga)
2183: {
2184: struct galois_analysis Pga;
2185: struct galois_borne Pgb;
2186: struct galois_test td;
2187: struct galois_testlift gt;
2188: struct galois_lift gl;
2189: ulong ltop = avma, lbot;
2190: long n, p, deg, ex;
2191: byteptr primepointer;
2192: long sg, Pm = 0, fp;
2193: long x = varn(T);
2194: int i, j;
2195: GEN Tmod, res, pf = gzero, split, psi, ip, mod, ppsi;
2196: GEN frob;
2197: GEN O;
2198: GEN P, PG, PL, Pden, PM, Pmod, Pp;
2199: GEN *lo; /* tired of casting */
2200: n = degree(T);
2201: if (!ga->deg)
2202: return gzero;
2203: p = ga->p;
2204: ex = ga->exception;
2205: deg = ga->deg;
2206: primepointer = ga->primepointer;
2207: if (DEBUGLEVEL >= 9)
2208: fprintferr("GaloisConj:denominator:%Z\n", den);
2209: if (n == 12 && ga->ord == 3) /* A4 is very probable,so test it first */
2210: {
2211: long av = avma;
2212: if (DEBUGLEVEL >= 4)
2213: fprintferr("GaloisConj:Testing A4 first\n");
2214: inittest(L, M, gb->bornesol, gb->ladicsol, &td);
2215: lbot = avma;
2216: PG = a4galoisgen(T, &td);
2217: freetest(&td);
2218: if (PG != gzero)
2219: {
2220: return gerepile(ltop, lbot, PG);
2221: }
2222: avma = av;
2223: }
2224: if (n == 24 && ga->ord == 3) /* S4 is very probable,so test it first */
2225: {
2226: long av = avma;
2227: if (DEBUGLEVEL >= 4)
2228: fprintferr("GaloisConj:Testing S4 first\n");
2229: lbot = avma;
2230: initlift(T, den, stoi(ga->p4), gb->bornesol, L, &gl);
2231: PG = s4galoisgen(&gl);
2232: if (PG != gzero)
2233: {
2234: return gerepile(ltop, lbot, PG);
2235: }
2236: avma = av;
2237: }
2238: frob = cgetg(lg(L), t_VECSMALL);
2239: for (;;)
2240: {
2241: ulong av = avma;
2242: long isram;
2243: long c;
2244: ip = stoi(p);
2245: Tmod = (GEN) factmod(T, ip);
2246: isram = 0;
2247: for (i = 1; i < lg(Tmod[2]) && !isram; i++)
2248: if (!gcmp1(((GEN **) Tmod)[2][i]))
2249: isram = 1;
2250: if (isram == 0)
2251: {
2252: fp = degree(((GEN **) Tmod)[1][1]);
2253: for (i = 2; i < lg(Tmod[1]); i++)
2254: if (degree(((GEN **) Tmod)[1][i]) != fp)
2255: {
2256: avma = ltop;
2257: return gzero; /* Not Galois polynomial */
2258: }
2259: if (fp % deg == 0)
2260: {
2261: if (DEBUGLEVEL >= 4)
2262: fprintferr("Galoisconj:p=%d deg=%d fp=%d\n", p, deg, fp);
2263: lo = (GEN *) listsousgroupes(deg, n / fp);
2264: initlift(T, den, ip, gb->bornesol, L, &gl);
2265: if (inittestlift((GEN) Tmod[1], fp / deg, &gl, >, frob))
2266: {
2267: int k;
2268: sg = lg(lo) - 1;
2269: psi = cgetg(gt.g + 1, t_VECSMALL);
2270: for (k = 1; k <= gt.g; k++)
2271: psi[k] = 1;
2272: goto suite;
2273: }
2274: if (DEBUGLEVEL >= 4)
2275: fprintferr("Galoisconj:Subgroups list:%Z\n", (GEN) lo);
2276: if (deg == fp && cgcd(deg, n / deg) == 1) /* Pretty sure it's the
2277: * biggest ? */
2278: {
2279: for (sg = 1; sg < lg(lo) - 1; sg++)
2280: if (frobeniusliftall(lo[sg], &psi, &gl, >, frob))
2281: goto suite;
2282: } else
2283: {
2284: for (sg = lg(lo) - 2; sg >= 1; sg--) /* else start with the
2285: * fastest! */
2286: if (frobeniusliftall(lo[sg], &psi, &gl, >, frob))
2287: goto suite;
2288: }
2289: if (ex == 1 && (n == 12 || n % 12 != 0))
2290: {
2291: avma = ltop;
2292: return gzero;
2293: } else
2294: {
2295: ex--;
2296: if (n % 12 == 0)
2297: {
2298: if (n % 36 == 0)
2299: deg = 5 - deg;
2300: else
2301: deg = 2; /* For group like Z/2ZxA4 */
2302: }
2303: if (n % 12 == 0 && ex == -5)
2304: err(warner, "galoisconj _may_ hang up for this polynomial");
2305: }
2306: }
2307: }
2308: c = *primepointer++;
2309: if (!c)
2310: err(primer1);
2311: p += c;
2312: if (DEBUGLEVEL >= 4)
2313: fprintferr("GaloisConj:next p=%d\n", p);
2314: avma = av;
2315: }
2316: suite: /* Djikstra probably hates me. (Linus
2317: * Torvalds linux/kernel/sched.c) */
2318: if (deg == n) /* Cyclique */
2319: {
2320: lbot = avma;
2321: res = cgetg(2, t_VEC);
2322: res[1] = lgetg(3, t_VEC);
2323: ((GEN **) res)[1][1] = gcopy(frob);
2324: ((GEN **) res)[1][2] = stoi(deg);
2325: return gerepile(ltop, lbot, res);
2326: }
2327: if (DEBUGLEVEL >= 9)
2328: fprintferr("GaloisConj:Frobenius:%Z\n", frob);
2329: O = permtocycle(frob);
2330: if (DEBUGLEVEL >= 9)
2331: fprintferr("GaloisConj:Orbite:%Z\n", O);
2332: {
2333: GEN S, Tp, Fi, Sp;
2334: long gp = n / fp;
2335: ppsi = cgetg(gp + 1, t_VECSMALL);
2336: mod = gmodulcp(gun, ip);
2337: P = corpsfixeorbitemod(L, O, x, mod, gb->l, &PL);
2338: S = corpsfixeinclusion(O, PL);
2339: S = vectopol(S, L, M, den, x);
2340: if (DEBUGLEVEL >= 6)
2341: fprintferr("GaloisConj:Inclusion:%Z\n", S);
2342: Pmod = (GEN) factmod(P, ip)[1];
2343: Tp = gmul(T, mod);
2344: Sp = gmul(S, mod);
2345: Pp = gmul(P, mod);
2346: if (DEBUGLEVEL >= 4)
2347: fprintferr("GaloisConj:psi=%Z\n", psi);
2348: if (DEBUGLEVEL >= 4)
2349: fprintferr("GaloisConj:Pmod=%Z\n", Pmod);
2350: if (DEBUGLEVEL >= 4)
2351: fprintferr("GaloisConj:Tmod=%Z\n", Tmod);
2352: for (i = 1; i <= gp; i++)
2353: {
2354: Fi = ggcd(Tp, gsubst((GEN) Pmod[i], x, Sp));
2355: Fi = gdiv(Fi, (GEN) Fi[lgef(Fi) - 1]);
2356: if (DEBUGLEVEL >= 6)
2357: fprintferr("GaloisConj:Fi=%Z %d", Fi, i);
2358: for (j = 1; j <= gp; j++)
2359: if (gegal(Fi, ((GEN **) Tmod)[1][j]))
2360: break;
2361: if (DEBUGLEVEL >= 6)
2362: fprintferr("-->%d\n", j);
2363: if (j == gp + 1)
2364: {
2365: avma = ltop;
2366: return gzero;
2367: }
2368: if (j == gp)
2369: {
2370: Pm = i;
2371: ppsi[i] = 1;
2372: } else
2373: ppsi[i] = psi[j];
2374: }
2375: if (DEBUGLEVEL >= 4)
2376: fprintferr("GaloisConj:Pm=%d ppsi=%Z\n", Pm, ppsi);
2377: galoisanalysis(P, &Pga, 0);
2378: if (Pga.deg == 0)
2379: {
2380: avma = ltop;
2381: return gzero; /* Avoid computing the discriminant */
2382: }
2383: Pden = gabs(mycoredisc(discsr(P)), DEFAULTPREC);
2384: Pgb.l = gb->l;
2385: initborne(P, Pden, &Pgb, Pga.ppp);
2386: if (Pgb.valabs > gb->valabs)
2387: PL = recalculeL(P, deriv(P, x), PL, gb, &Pgb);
2388: PM = vandermondeinverse(PL, P, Pden);
2389: PG = galoisgen(P, PL, PM, Pden, &Pgb, &Pga);
2390: }
2391: if (DEBUGLEVEL >= 4)
2392: fprintferr("GaloisConj:Retour sur Terre:%Z\n", PG);
2393: if (PG == gzero)
2394: {
2395: avma = ltop;
2396: return gzero;
2397: }
2398: inittest(L, M, gb->bornesol, gb->ladicsol, &td);
2399: split = splitorbite(O);
2400: lbot = avma;
2401: res = cgetg(lg(PG) + lg(split) - 1, t_VEC);
2402: for (i = 1; i < lg(split); i++)
2403: res[i] = lcopy((GEN) split[i]);
2404: for (j = 1; j < lg(PG); j++)
2405: {
2406: ulong lbot = 0, av = avma;
2407: GEN B, tau;
2408: long t, g;
2409: long w, sr, dss;
2410: if (DEBUGLEVEL >= 6)
2411: fprintferr("GaloisConj:G[%d][1]=%Z\n", j, ((GEN **) PG)[j][1]);
2412: B = permtocycle(((GEN **) PG)[j][1]);
2413: if (DEBUGLEVEL >= 6)
2414: fprintferr("GaloisConj:B=%Z\n", B);
2415: tau = gmul(mod, permtopol(((GEN **) PG)[j][1], PL, PM, Pden, x));
2416: tau = gsubst((GEN) Pmod[Pm], x, tau);
2417: tau = ggcd(Pp, tau);
2418: tau = gdiv(tau, (GEN) tau[lgef(tau) - 1]);
2419: if (DEBUGLEVEL >= 6)
2420: fprintferr("GaloisConj:tau=%Z\n", tau);
2421: for (g = 1; g < lg(Pmod); g++)
2422: if (gegal(tau, (GEN) Pmod[g]))
2423: break;
2424: if (g == lg(Pmod))
2425: {
2426: freetest(&td);
2427: avma = ltop;
2428: return gzero;
2429: }
2430: w = ((long **) lo)[sg][ppsi[g]];
2431: dss = deg / cgcd(deg, w - 1);
2432: sr = 1;
2433: for (i = 1; i < lg(B[1]) - 1; i++)
2434: sr = (1 + sr * w) % deg;
2435: sr = cgcd(sr, deg);
2436: if (DEBUGLEVEL >= 6)
2437: fprintferr("GaloisConj:w=%d [%d] sr=%d dss=%d\n", w, deg, sr, dss);
2438: for (t = 0; t < sr; t += dss)
2439: {
2440: lbot = avma;
2441: pf = testpermutation(O, B, w, t, sr, &td);
2442: if (pf != gzero)
2443: break;
2444: }
2445: if (pf == gzero)
2446: {
2447: freetest(&td);
2448: avma = ltop;
2449: return gzero;
2450: } else
2451: {
2452: GEN f;
2453: f = cgetg(3, t_VEC);
2454: f[1] = (long) pf;
2455: f[2] = (long) gcopy(((GEN **) PG)[j][2]);
2456: res[lg(split) + j - 1] = (long) gerepile(av, lbot, f);
2457: }
2458: }
2459: if (DEBUGLEVEL >= 4)
2460: fprintferr("GaloisConj:Fini!\n");
2461: freetest(&td);
2462: return gerepile(ltop, lbot, res);
2463: }
2464: /*
2465: * T: polynome ou nf den optionnel: si présent doit etre un multiple du
2466: * dénominateur commun des solutions Si T est un nf et den n'est pas présent,
2467: * le denominateur de la base d'entiers est pris pour den
2468: */
2469: GEN
2470: galoisconj4(GEN T, GEN den, long flag)
2471: {
2472: ulong ltop = avma, lbot;
2473: GEN G, L, M, res, aut, L2, M2, grp;
2474: struct galois_analysis ga;
2475: struct galois_borne gb;
2476: int n, i, j, k;
2477: if (typ(T) != t_POL)
2478: {
2479: T = checknf(T);
2480: if (!den)
2481: den = gcoeff((GEN) T[8], degree((GEN) T[1]), degree((GEN) T[1]));
2482: T = (GEN) T[1];
2483: }
2484: n = lgef(T) - 3;
2485: if (n <= 0)
2486: err(constpoler, "galoisconj4");
2487: for (k = 2; k <= n + 2; k++)
2488: if (typ(T[k]) != t_INT)
2489: err(talker, "polynomial not in Z[X] in galoisconj4");
2490: if (!gcmp1((GEN) T[n + 2]))
2491: err(talker, "non-monic polynomial in galoisconj4");
2492: n = degree(T);
2493: if (n == 1) /* Too easy! */
2494: {
2495: res = cgetg(2, t_COL);
2496: res[1] = (long) polx[varn(T)];
2497: return res;
2498: }
2499: galoisanalysis(T, &ga, 1);
2500: if (ga.deg == 0)
2501: {
2502: avma = ltop;
2503: return stoi(ga.p); /* Avoid computing the discriminant */
2504: }
2505: if (!den)
2506: den = absi(mycoredisc(discsr(T)));
2507: else
2508: {
2509: if (typ(den) != t_INT)
2510: err(talker, "Second arg. must be integer in galoisconj4");
2511: den = absi(den);
2512: }
2513: gb.l = stoi(ga.l);
2514: initborne(T, den, &gb, ga.ppp);
2515: if (DEBUGLEVEL >= 1)
2516: timer2();
2517: {
2518: GEN f = rootpadic(T, gb.l, gb.valabs);
2519: GEN _p = gmael(f, 1, 2);
2520: L = gmul(gtrunc(f), gmodulcp(gun, gb.ladicabs));
2521: gunclone(_p);
2522: }
2523: if (DEBUGLEVEL >= 1)
2524: msgtimer("rootpadic()");
2525: M = vandermondeinverse(L, T, den);
2526: G = galoisgen(T, L, M, den, &gb, &ga);
2527: if (DEBUGLEVEL >= 6)
2528: fprintferr("GaloisConj:%Z\n", G);
2529: if (G == gzero)
2530: {
2531: avma = ltop;
2532: return gzero;
2533: }
2534: if (DEBUGLEVEL >= 1)
2535: timer2();
2536: L2 = gmul(L, gmodulcp(gun, gb.ladicsol));
2537: M2 = gmul(M, gmodulcp(gun, gb.ladicsol));
2538: if (flag)
2539: {
2540: lbot = avma;
2541: grp = cgetg(7, t_VEC);
2542: grp[1] = (long) gcopy(T);
2543: grp[2] = (long) stoi(ga.l);
2544: grp[3] = (long) gcopy(L);
2545: grp[4] = (long) gcopy(M);
2546: grp[5] = (long) gcopy(den);
2547: }
2548: res = cgetg(n + 1, t_VEC);
2549: res[1] = (long) permidentity(n);
2550: k = 1;
2551: for (i = 1; i < lg(G); i++)
2552: {
2553: int c = k * (itos(((GEN **) G)[i][2]) - 1);
2554: for (j = 1; j <= c; j++) /* I like it */
2555: res[++k] = (long) applyperm((GEN) res[j], ((GEN **) G)[i][1]);
2556: }
2557: if (flag)
2558: {
2559: grp[6] = (long) res;
2560: return gerepile(ltop, lbot, grp);
2561: }
2562: aut = cgetg(n + 1, t_COL);
2563: for (i = 1; i <= n; i++)
2564: aut[i] = (long) permtopol((GEN) res[i], L2, M2, den, varn(T));
2565: if (DEBUGLEVEL >= 1)
2566: msgtimer("Calcul polynomes");
2567: return gerepileupto(ltop, aut);
2568: }
2569: /* Calcule le nombre d'automorphisme de T avec forte probabilité */
2570: /* pdepart premier nombre premier a tester */
2571: long
2572: numberofconjugates(GEN T, long pdepart)
2573: {
2574: ulong ltop = avma, ltop2;
2575: long n, p, nbmax, nbtest;
2576: long card;
2577: byteptr primepointer;
2578: int i;
2579: GEN L;
2580: n = degree(T);
2581: card = sturm(T);
2582: card = cgcd(card, n - card);
2583: nbmax = n + 1;
2584: nbtest = 0;
2585: L = cgetg(n + 1, t_VECSMALL);
2586: ltop2 = avma;
2587: for (p = 0, primepointer = diffptr; nbtest < nbmax && card > 1;)
2588: {
2589: long s, c;
2590: long isram;
2591: GEN S;
2592: c = *primepointer++;
2593: if (!c)
2594: err(primer1);
2595: p += c;
2596: if (p < pdepart)
2597: continue;
2598: S = (GEN) simplefactmod(T, stoi(p));
2599: isram = 0;
2600: for (i = 1; i < lg(S[2]) && !isram; i++)
2601: if (!gcmp1(((GEN **) S)[2][i]))
2602: isram = 1;
2603: if (isram == 0)
2604: {
2605: for (i = 1; i <= n; i++)
2606: L[i] = 0;
2607: for (i = 1; i < lg(S[1]) && !isram; i++)
2608: L[itos(((GEN **) S)[1][i])]++;
2609: s = L[1];
2610: for (i = 2; i <= n; i++)
2611: s = cgcd(s, L[i] * i);
2612: card = cgcd(s, card);
2613: }
2614: if (DEBUGLEVEL >= 6)
2615: fprintferr("NumberOfConjugates:Nbtest=%d,card=%d,p=%d\n", nbtest, card, p);
2616: nbtest++;
2617: avma = ltop2;
2618: }
2619: if (DEBUGLEVEL >= 2)
2620: fprintferr("NumberOfConjugates:card=%d,p=%d\n", card, p);
2621: avma = ltop;
2622: return card;
2623: }
2624: GEN
2625: galoisconj0(GEN nf, long flag, GEN d, long prec)
2626: {
2627: GEN G, T;
2628: long card;
2629: if (typ(nf) != t_POL)
2630: {
2631: nf = checknf(nf);
2632: T = (GEN) nf[1];
2633: } else
2634: T = nf;
2635: switch (flag)
2636: {
2637: case 0:
2638: G = galoisconj4(nf, d, 0);
2639: if (typ(G) != t_INT) /* Success */
2640: return G;
2641: else
2642: {
2643: card = G == gzero ? degree(T) : numberofconjugates(T, itos(G));
2644: if (card != 1)
2645: {
2646: if (typ(nf) == t_POL)
2647: return galoisconj2pol(nf, card, prec);
2648: else
2649: return galoisconj(nf);
2650: }
2651: }
2652: break; /* Failure */
2653: case 1:
2654: return galoisconj(nf);
2655: case 2:
2656: return galoisconj2(nf, degree(T), prec);
2657: case 4:
2658: G = galoisconj4(nf, d, 0);
2659: if (typ(G) != t_INT)
2660: return G;
2661: break; /* Failure */
2662: default:
2663: err(flagerr, "nfgaloisconj");
2664: }
2665: G = cgetg(2, t_COL);
2666: G[1] = (long) polx[varn(T)];
2667: return G; /* Failure */
2668: }
2669: /******************************************************************************/
2670: /* Galois theory related algorithms */
2671: /******************************************************************************/
2672: GEN
2673: checkgal(GEN gal)
2674: {
2675: if (typ(gal) == t_POL)
2676: err(talker, "please apply galoisinit first");
2677: if (typ(gal) != t_VEC || lg(gal) != 7)
2678: err(talker, "Not a Galois field in a Galois related function");
2679: return gal;
2680: }
2681: GEN
2682: galoisinit(GEN nf, GEN den)
2683: {
2684: GEN G;
2685: G = galoisconj4(nf, den, 1);
2686: if (typ(G) == t_INT)
2687: err(talker, "galoisinit: field not Galois or Galois group not weakly super solvable");
2688: return G;
2689: }
2690: GEN
2691: galoispermtopol(GEN gal, GEN perm)
2692: {
2693: gal = checkgal(gal);
2694: if (typ(perm) != t_VEC)
2695: err(typeer, "galoispermtopol:");
2696: return permtopol(perm, (GEN) gal[3], (GEN) gal[4], (GEN) gal[5], varn((GEN) gal[1]));
2697: }
2698: GEN
2699: galoisfixedfield(GEN gal, GEN perm, GEN p)
2700: {
2701: ulong ltop = avma, lbot;
2702: GEN P, S, PL, O, res, mod;
2703: long x;
2704: x = varn((GEN) gal[1]);
2705: gal = checkgal(gal);
2706: if (typ(perm) != t_VEC)
2707: err(typeer, "galoisfixedfield:");
2708: if (p)
2709: {
2710: if (typ(p) != t_INT)
2711: err(talker, "galoisfixedfield: optionnal argument must be an prime integer");
2712: mod = gmodulcp(gun, p);
2713: } else
2714: mod = gun;
2715: O = vecpermorbite(perm);
2716: P = corpsfixeorbitemod((GEN) gal[3], O, x, mod, (GEN) gal[2], &PL);
2717: S = corpsfixeinclusion(O, PL);
2718: S = vectopol(S, (GEN) gal[3], (GEN) gal[4], (GEN) gal[5], x);
2719: lbot = avma;
2720: res = cgetg(3, t_VEC);
2721: res[1] = (long) gcopy(P);
2722: res[2] = (long) gmodulcp(S, (GEN) gal[1]);
2723: return gerepile(ltop, lbot, res);
2724: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>