Annotation of OpenXM_contrib/pari/src/basemath/galconj.c, Revision 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>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to galconj.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / pari / src / basemath