| Loic Grenie on Sun, 19 Nov 2006 17:27:29 +0100 |
[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]
| rnfkummer modification |
Hi again,
I've asked in the pari-users mailing list help for determining a
"basis" of linearly indenpendant Kummer extensions of a base field
(or bnr). I've devised a modification of _rnfkummer and rndkummersimple
that does the following:
- rnfkummer(bnr, subgroup, -ell) will output a list of linearly
independant extensions for the given subgroup.
--> - rnfkummer(bnr, subgroup, -1) will output the same thing as
--> before, rnfkummer(bnr, subgroup)) i.e. the equation of a Kummer
--> extension for the given subgroup but using an algorithm which
--> should be, most of the time,
--> QUADRATIC INSTEAD OF EXPONENTIAL
--> (sorry for shouting), when going through rnfkummersimple() (it
--> might apply also to _rnfkummer but I don't know and I've
--> preferred not to do something that is beyond my limited knowledge).
While the second feature is not very useful to me, the first one is
nearly vital (well, I can always apply the patch myself !). The second
is a incentive for the developpers to look benevolently at my patch !
(I confess that I'm very proud of it, though).
Full set of disclaimers:
- I don't know the theorem that are behind the algorithm. I've just
used the F_ell-vector space structure of the
Gal(compositum(Kummer extensions of deg ell of K))/K (when K
contains all ell-th roots of unity) + a trick to improve the
probability that ok_congruence() will return 1.
- I've tried to check it only on a small set of examples.
- I'm a very bad programmer.
Note: I've not indented two do {...} while (firstpass--); loops
to minimize the set of differences and make clearer what I've changed.
LoïcIndex: src/modules/kummer.c
===================================================================
RCS file: /home/cvs/pari/src/modules/kummer.c,v
retrieving revision 1.130
diff -c -r1.130 kummer.c
*** src/modules/kummer.c 30 Sep 2006 23:43:03 -0000 1.130
--- src/modules/kummer.c 19 Nov 2006 15:47:15 -0000
***************
*** 340,348 ****
get_gell(GEN bnr, GEN subgp, long all)
{
GEN gell;
! if (all) gell = stoi(all);
! else if (subgp) gell = det(subgp);
! else gell = det(diagonal_i(gmael(bnr,5,2)));
if (typ(gell) != t_INT) pari_err(arither1);
return gell;
}
--- 340,348 ----
get_gell(GEN bnr, GEN subgp, long all)
{
GEN gell;
! if (all && all != -1) gell = stoi(all);
! else if (subgp) gell = det(subgp);
! else gell = det(diagonal_i(gmael(bnr,5,2)));
if (typ(gell) != t_INT) pari_err(arither1);
return gell;
}
***************
*** 509,514 ****
--- 509,537 ----
return x;
}
+ /* A column vector representing a subgroup of prime index */
+ /* !! Not Stack Clean: return pointers to subgroup's elements !! */
+ static GEN
+ grptocol(GEN subgroup)
+ {
+ long l, i, j;
+ GEN col;
+
+ l = lg(subgroup);
+ col = cgetg(l, t_COL);
+ for (i = 1; i < l; i++)
+ if (gcmp1(gcoeff(subgroup, i, i)))
+ gel(col, i) = gen_0;
+ else
+ {
+ gel(col, i) = gen_m1;
+ break;
+ }
+ for (j=i; ++j < l; )
+ gel(col, j) = gcoeff(subgroup, i, j);
+ return col;
+ }
+
/* if all!=0, give all equations of degree 'all'. Assume bnr modulus is the
* conductor */
***************
*** 518,528 ****
long ell, i, j, degK, dK;
long lSml2, lSl2, lSp, rc, lW;
long prec;
GEN bnf,nf,bid,ideal,arch,cycgen;
GEN clgp,cyc;
GEN Sp,listprSp,matP;
! GEN res,u,M,K,y,vecMsup,vecW,vecWB,vecBp,msign;
primlist L;
bnf = gel(bnr,1);
--- 541,554 ----
long ell, i, j, degK, dK;
long lSml2, lSl2, lSp, rc, lW;
long prec;
+ long rk=0, ncyc=0;
+ GEN mat=NULL, matgrp=NULL, xell;
+ long firstpass = all<0;
GEN bnf,nf,bid,ideal,arch,cycgen;
GEN clgp,cyc;
GEN Sp,listprSp,matP;
! GEN res=NULL,u,M,K,y,vecMsup,vecW,vecWB,vecBp,msign;
primlist L;
bnf = gel(bnr,1);
***************
*** 551,560 ****
matP = cgetg(lSp+1, t_MAT);
for (j=1; j<=lSp; j++)
{
! GEN e, a, L;
L = isprincipalell(bnf,gel(Sp,j), cycgen,u,gell,rc);
! e = gel(L,1); gel(matP,j) = e;
! a = gel(L,2); gel(vecBp,j) = a;
}
vecWB = shallowconcat(vecW, vecBp);
--- 577,586 ----
matP = cgetg(lSp+1, t_MAT);
for (j=1; j<=lSp; j++)
{
! GEN L;
L = isprincipalell(bnf,gel(Sp,j), cycgen,u,gell,rc);
! gel( matP,j) = gel(L,1);
! gel(vecBp,j) = gel(L,2);
}
vecWB = shallowconcat(vecW, vecBp);
***************
*** 580,615 ****
lW = lg(vecW);
M = vconcat(M, shallowconcat(zeromat(rc,lW-1), matP));
! K = FpM_ker(M, gell);
! dK = lg(K)-1;
y = cgetg(dK+1,t_VECSMALL);
! res = cgetg(1,t_VEC); /* in case all = 1 */
while (dK)
{
for (i=1; i<dK; i++) y[i] = 0;
y[i] = 1; /* y = [0,...,0,1,0,...,0], 1 at dK'th position */
do
{
pari_sp av = avma;
! GEN be, P, X = FpC_red(ZM_zc_mul(K, y), gell);
if (ok_congruence(X, gell, lW, vecMsup) && ok_sign(X, msign, arch))
{/* be satisfies all congruences, x^ell - be is irreducible, signature
! * and relative discriminant are correct */
! be = compute_beta(X, vecWB, gell, bnf);
! be = lift_if_rational(coltoalg(nf, be));
! P = gsub(monomial(gen_1, ell, 0), be);
! if (all) res = shallowconcat(res, gerepileupto(av, P));
! else
! {
! if (gequal(rnfnormgroup(bnr,P),subgroup)) return P; /*DONE*/
! avma = av;
! }
}
else avma = av;
} while (increment(y, dK, ell));
y[dK--] = 0;
}
! if (all) return res;
return gen_0;
}
--- 606,854 ----
lW = lg(vecW);
M = vconcat(M, shallowconcat(zeromat(rc,lW-1), matP));
! if (all < 0)
! {
! pari_sp av2 = avma;
! GEN Kidx;
! long idx, ffree;
! /* Reorganize the kernel so that the tests of ok_congruence can be ok
! * for y[ncyc]=1 and y[1..ncyc]=1 */
! K = ker(gmodulo(M, gell));
! dK = lg(K)-1;
! Kidx = cgetg(dK+1, t_VECSMALL);
! /* First step: Gauss elimination on vectors lW...lg(M) */
! for (idx = lg(K), i=lg(M); --i >= lW; )
! {
! GEN tmp, Ki;
! for (j=dK; j > 0; j--) if (!gcmp0(gcoeff(K, i, j))) break;
! if (!j)
! {
! /* Do our best to ensure that K[dK,i]!=0 */
! if (gcmp0(gcoeff(K, i, dK)))
! {
! for (j = idx; j < dK; j++)
! if (!gcmp0(gcoeff(K, i, j)) &&
! itos(gcoeff(K, Kidx[j], dK)) < ell - 1)
! gel(K, dK) = gadd(gel(K, dK), gel(K, j));
! }
! continue;
! }
! if (j != --idx)
! {
! tmp = gel(K, j);
! gel(K, j) = gel(K, idx);
! gel(K, idx) = tmp;
! }
! Kidx[idx] = i;
! if (!gcmp1(gcoeff(K, idx, idx)))
! gel(K, idx) = gdiv(gel(K, idx), gcoeff(K, i, idx));
! Ki = gel(K, idx);
! if (!gcmp1(gcoeff(K, i, dK)))
! gel(K, dK) = gsub(gel(K, dK), gmul(gsubgs(gcoeff(K, i, dK), 1), Ki));
! for (j = dK; --j > 0; )
! {
! if (j == idx) continue;
! if (!gcmp0(gcoeff(K, i, j)))
! gel(K, j) = gsub(gel(K, j), gmul(gcoeff(K, i, j), Ki));
! }
! }
! /* ffree = first vector that is not "free" for the scalar products */
! ffree = idx;
! /* Second step: for each hyperplace equation in vecMsup, do the same
! * thing as before. */
! for (i=1; i < lg(vecMsup); i++)
! {
! GEN tmp, Ki;
! long scalarprod;
! for (j=ffree; --j > 0; )
! if (lg(gel(vecMsup, i)) == 2 &&
! (scalarprod = itos(lift(gmul(gel(vecMsup, i), gel(K, j))))))
! break;
! if (!j)
! {
! /* Do our best to ensure that vecMsup.K[dK] != 0 */
! if (!signe(lift(gmul(gel(vecMsup, i), gel(K, dK)))))
! {
! for (j = ffree; j < dK; j++)
! if (itos(lift(gmul(gel(vecMsup, i), gel(K, j)))) &&
! itos(gcoeff(K, Kidx[j], dK)) < ell - 1)
! gel(K, dK) = gadd(gel(K, dK), gel(K, j));
! }
! continue;
! }
! if (j != --ffree)
! {
! tmp = gel(K, j);
! gel(K, j) = gel(K, ffree);
! gel(K, ffree) = tmp;
! }
! if (!signe(lift(gmul(gel(vecMsup, i), gel(K, dK)))))
! gel(K, ffree) = gdiv(gel(K, ffree), gcoeff(K, i, ffree));
! Ki = gel(K, ffree);
! if (!gcmp1(gcoeff(K, i, dK)))
! gel(K, dK) = gsub(gel(K, dK), gmul(gsubgs(gcoeff(K, i, dK), 1), Ki));
! for (j = dK; --j > 0; )
! {
! if (j == ffree) continue;
! if (!gcmp0(gcoeff(K, i, j)))
! gel(K, j) = gsub(gel(K, j), gmul(gcoeff(K, i, j), Ki));
! }
! }
! if (ell == 2)
! {
! for (i = ffree, j = ffree - 1; i < dK; i++, j--)
! {
! GEN tmp = gel(K, i);
! gel(K, i) = gel(K, j);
! gel(K, j) = tmp;
! }
! }
! K = gerepileupto(av2, lift(K));
! }
! else
! {
! K = FpM_ker(M, gell);
! dK = lg(K)-1;
! }
y = cgetg(dK+1,t_VECSMALL);
! if (all) res = cgetg(1,t_VEC); /* in case all = 1 */
! if (all<0)
! {
! GEN col;
! ncyc = dK;
! col = cgetg(lg(M)+1, t_COL);
! for (i=1; i <= lg(M); i++)
! gel(col, i) = gen_0;
! mat = cgetg(ncyc+1, t_MAT);
! for (i = 1; i <= ncyc; i++)
! gel(mat, i) = col;
! if (all == -1)
! {
! long l = lg(gmael(bnr, 5, 2));
! col = cgetg(l, t_COL);
! for (i = 1; i < l; i++)
! gel(col, i) = gen_0;
! matgrp = cgetg(ncyc+2, t_MAT);
! gel(matgrp, 1) = grptocol(subgroup);
! for (i = 2; i <= ncyc+1; i++)
! gel(matgrp, i) = col;
! }
! rk = 0;
! }
! xell = monomial(gen_1, ell, 0);
! do {
! dK = lg(K)-1;
while (dK)
{
for (i=1; i<dK; i++) y[i] = 0;
y[i] = 1; /* y = [0,...,0,1,0,...,0], 1 at dK'th position */
+ if (firstpass) y[ncyc] = 1;
do
{
pari_sp av = avma;
! GEN be, P=NULL, X = FpC_red(ZM_zc_mul(K, y), gell);
if (ok_congruence(X, gell, lW, vecMsup) && ok_sign(X, msign, arch))
{/* be satisfies all congruences, x^ell - be is irreducible, signature
! * and relative discriminant are correct */
! if (all<0)
! {
! gel(mat, rk+1) = X;
! if (rank(mat) > rk)
! {
! rk++;
! gel(mat, rk) = gclone(X);
! }
! else
! continue;
! }
! be = compute_beta(X, vecWB, gell, bnf);
! be = lift_if_rational(coltoalg(nf, be));
! if (all != -1) P = gsub(xell, be);
! if (all)
! {
! res = shallowconcat(res, gerepileupto(av, (all+1) ? P:gcopy(be)));
! }
! else
! {
! if (gequal(rnfnormgroup(bnr,P),subgroup)) return P; /*DONE*/
! avma = av;
! continue;
! }
! if (all == -1)
! {
! pari_sp av2 = avma;
! GEN Kgrp, colgrp;
! long dKgrp;
! colgrp = grptocol(rnfnormgroup(bnr, gsub(xell, gel(res, rk))));
! if (ell != 2 && rk != 1)
! {
! /* Compute the pesky scalar */
! GEN mat2, K2, multiplicator;
! mat2 = cgetg(4, t_MAT);
! gel(mat2, 1) = gel(matgrp, 2);
! gel(mat2, 2) = colgrp;
! gel(mat2, 3) = grptocol(rnfnormgroup(bnr, gsub(xell, gmul(gel(res, 1), gel(res, rk)))));
! K2 = FpM_ker(mat2, gell);
! if (lg(K2) != 2)
! pari_err(talker, "L'algèbre linéaire m'a tuer");
! K2 = gel(K2, 1);
! if (!gequal(gel(K2, 1), gel(K2, 2)))
! {
! multiplicator = gdiv(gel(K2, 2), gmodulo(gel(K2, 1), gell));
! colgrp = lift(gmul(multiplicator, colgrp));
! }
! }
! gel(matgrp, rk+1) = gclone(colgrp);
! setlg(matgrp, rk+2);
! Kgrp = FpM_ker(matgrp, gell);
! setlg(matgrp, ncyc+2);
! dKgrp = lg(Kgrp)-1;
! if (dKgrp /* == 1 */)
! {
! long l;
! Kgrp = gel(Kgrp, 1);
! l = lg(Kgrp);
! for (i = 2; i < l; i++) {
! GEN pow = gel(Kgrp, i);
!
! if (signe(pow))
! {
! if (pow[2] != 1)
! gel(res, i-1) = gpowgs(gel(res, i-1), pow[2]);
! }
! else
! gel(res, i-1) = gen_1;
! }
! be = divide_conquer_prod(res, gmul);
! be = reducebeta(bnf, be, gell);
! be = coltoalg(nf, be);
! res = gsub(xell, be);
! goto DONE2;
! }
! avma = av2;
! }
! if (all < 0 && rk == ncyc)
! goto DONE2;
! if (all < 0) break;
}
else avma = av;
} while (increment(y, dK, ell));
y[dK--] = 0;
}
! } while (firstpass--);
! if (all)
! {
! if (all<0)
! {
! DONE2:
! for (i=1;i<=rk;i++)
! gunclone(gel(mat, i));
! if (all == -1)
! for (i=1;i<=rk;i++)
! gunclone(gel(matgrp, i+1));
! }
! return res;
! }
return gen_0;
}
***************
*** 882,894 ****
GEN polnf,bnf,nf,bnfz,nfz,bid,ideal,cycgen,gell,p1,p2,wk,U,vselmer;
GEN clgp,cyc,gen;
GEN Q,idealz,gothf;
! GEN res,u,M,K,y,vecMsup,vecW,vecWA,vecWB,vecB,vecC,vecAp,vecBp;
GEN matP,Sp,listprSp,Tc,Tv,P;
primlist L;
toK_s T;
tau_s _tau, *tau;
compo_s COMPO;
pari_timer t;
if (DEBUGLEVEL) TIMERstart(&t);
checkbnrgen(bnr);
--- 1121,1136 ----
GEN polnf,bnf,nf,bnfz,nfz,bid,ideal,cycgen,gell,p1,p2,wk,U,vselmer;
GEN clgp,cyc,gen;
GEN Q,idealz,gothf;
! GEN res=NULL,u,M,K,y,vecMsup,vecW,vecWA,vecWB,vecB,vecC,vecAp,vecBp;
GEN matP,Sp,listprSp,Tc,Tv,P;
primlist L;
toK_s T;
tau_s _tau, *tau;
compo_s COMPO;
pari_timer t;
+ long rk=0, ncyc=0;
+ GEN mat=NULL;
+ long firstpass = all<0;
if (DEBUGLEVEL) TIMERstart(&t);
checkbnrgen(bnr);
***************
*** 902,913 ****
if (DEBUGLEVEL) msgTIMER(&t, "[rnfkummer] conductor");
bnr = gel(p1,2);
subgroup = gel(p1,3);
! gell = get_gell(bnr,subgroup,all);
ell = itos(gell);
if (ell == 1) return pol_x(vnf);
if (!uisprime(ell)) pari_err(impl,"kummer for composite relative degree");
if (!umodiu(wk,ell)) return rnfkummersimple(bnr, subgroup, gell, all);
bid = gel(bnr,2);
ideal = gmael(bid,1,1);
/* step 1 of alg 5.3.5. */
--- 1144,1156 ----
if (DEBUGLEVEL) msgTIMER(&t, "[rnfkummer] conductor");
bnr = gel(p1,2);
subgroup = gel(p1,3);
! gell = get_gell(bnr,subgroup,all<-1?-all:all);
ell = itos(gell);
if (ell == 1) return pol_x(vnf);
if (!uisprime(ell)) pari_err(impl,"kummer for composite relative degree");
if (!umodiu(wk,ell)) return rnfkummersimple(bnr, subgroup, gell, all);
+ if (all == -1) all = 0;
bid = gel(bnr,2);
ideal = gmael(bid,1,1);
/* step 1 of alg 5.3.5. */
***************
*** 1073,1080 ****
if (DEBUGLEVEL>2) fprintferr("Step 18\n");
dK = lg(K)-1;
y = cgetg(dK+1,t_VECSMALL);
! res = cgetg(1, t_VEC); /* in case all = 1 */
if (DEBUGLEVEL) msgTIMER(&t, "[rnfkummer] candidate list");
while (dK)
{
for (i=1; i<dK; i++) y[i] = 0;
--- 1316,1337 ----
if (DEBUGLEVEL>2) fprintferr("Step 18\n");
dK = lg(K)-1;
y = cgetg(dK+1,t_VECSMALL);
! if (all) res = cgetg(1, t_VEC);
if (DEBUGLEVEL) msgTIMER(&t, "[rnfkummer] candidate list");
+ if (all<0)
+ {
+ GEN col;
+ ncyc = dK;
+ col = cgetg(lg(M)+1, t_COL);
+ for (i = 1; i <= lg(M); i++)
+ gel(col, i) = gen_0;
+ mat = cgetg(ncyc+1, t_MAT);
+ for (i = 1; i <= ncyc; i++)
+ gel(mat, i) = col;
+ rk = 0;
+ }
+ do {
+ dK = lg(K)-1;
while (dK)
{
for (i=1; i<dK; i++) y[i] = 0;
***************
*** 1084,1105 ****
GEN H, be, P, X = FpC_red(ZM_zc_mul(K, y), gell);
if (ok_congruence(X, gell, lW, vecMsup))
{
! be = compute_beta(X, vecWB, gell, bnfz);
! P = compute_polrel(nfz, &T, be, g, ell);
! P = lift_if_rational(P);
! if (DEBUGLEVEL>1) fprintferr("polrel(beta) = %Z\n", P);
! H = rnfnormgroup(bnr, P);
! if (!all) {
! if (gequal(subgroup, H)) return P; /* DONE */
! } else {
! if (!gequal(subgroup,H) && conductor(bnr, H, -1) == gen_0) continue;
! }
! res = shallowconcat(res, P);
}
} while (increment(y, dK, ell));
y[dK--] = 0;
}
! if (all) return res;
return gen_0; /* FAIL */
}
--- 1341,1387 ----
GEN H, be, P, X = FpC_red(ZM_zc_mul(K, y), gell);
if (ok_congruence(X, gell, lW, vecMsup))
{
! if (all<0)
! {
! gel(mat, rk+1) = gclone(X);
! if (rank(mat) > rk)
! {
! rk++;
! gel(mat, rk) = gclone(X);
! }
! else
! continue;
! }
! be = compute_beta(X, vecWB, gell, bnfz);
! P = compute_polrel(nfz, &T, be, g, ell);
! P = lift_if_rational(P);
! if (DEBUGLEVEL>1) fprintferr("polrel(beta) = %Z\n", P);
! H = rnfnormgroup(bnr, P);
! if (!all) {
! if (gequal(subgroup, H)) return P; /* DONE */
! continue;
! } else {
! if (!gequal(subgroup,H) && conductor(bnr, H, -1) == gen_0) {
! continue;
! }
! }
! res = shallowconcat(res, P);
! if (all < 0 && rk == ncyc)
! goto DONE2;
! if (all < 0) break;
}
} while (increment(y, dK, ell));
y[dK--] = 0;
}
! } while (firstpass--);
! if (all)
! {
! if (all<0)
! DONE2:
! for (i=1;i<=rk;i++)
! gunclone(gel(mat, i));
! return res;
! }
return gen_0; /* FAIL */
}