Line data Source code
1 : /* Copyright (C) 2000 The PARI group.
2 :
3 : This file is part of the PARI/GP package.
4 :
5 : PARI/GP is free software; you can redistribute it and/or modify it under the
6 : terms of the GNU General Public License as published by the Free Software
7 : Foundation; either version 2 of the License, or (at your option) any later
8 : version. It is distributed in the hope that it will be useful, but WITHOUT
9 : ANY WARRANTY WHATSOEVER.
10 :
11 : Check the License for details. You should have received a copy of it, along
12 : with the package; see the file 'COPYING'. If not, write to the Free Software
13 : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
14 :
15 : /*******************************************************************/
16 : /* */
17 : /* KUMMER EXTENSIONS */
18 : /* */
19 : /*******************************************************************/
20 : #include "pari.h"
21 : #include "paripriv.h"
22 :
23 : #define DEBUGLEVEL DEBUGLEVEL_bnrclassfield
24 :
25 : typedef struct {
26 : GEN R; /* nf.pol */
27 : GEN x; /* tau ( Mod(x, R) ) */
28 : GEN zk;/* action of tau on nf.zk (as t_MAT) */
29 : } tau_s;
30 :
31 : typedef struct {
32 : GEN polnf, invexpoteta1, powg;
33 : tau_s *tau;
34 : long m;
35 : } toK_s;
36 :
37 : typedef struct {
38 : GEN R; /* ZX, compositum(P,Q) */
39 : GEN p; /* QX, Mod(p,R) root of P */
40 : GEN q; /* QX, Mod(q,R) root of Q */
41 : long k; /* Q[X]/R generated by q + k p */
42 : GEN rev;
43 : } compo_s;
44 :
45 : /* REDUCTION MOD ell-TH POWERS */
46 : /* make b integral by multiplying by t in (Q^*)^ell */
47 : static GEN
48 36743 : reduce_mod_Qell(GEN nf, GEN b, ulong ell)
49 : {
50 : GEN c;
51 36743 : b = nf_to_scalar_or_basis(nf, b);
52 36743 : b = Q_primitive_part(b, &c);
53 36742 : if (c)
54 : {
55 9578 : GEN d, fa = Q_factor_limit(c, 1000000);
56 9578 : d = factorback2(gel(fa,1), ZV_to_Flv(gel(fa,2), ell));
57 9578 : b = typ(b) == t_INT? mulii(b,d): ZC_Z_mul(b, d);
58 : }
59 36742 : return b;
60 : }
61 :
62 : static GEN
63 36743 : reducebeta(GEN bnfz, GEN b, long ell)
64 : {
65 36743 : GEN t, cb, fu, nf = bnf_get_nf(bnfz);
66 :
67 36743 : if (DEBUGLEVEL>1) err_printf("reducing beta = %Ps\n",b);
68 36743 : b = reduce_mod_Qell(nf, b, ell);
69 36742 : t = idealredmodpower(nf, b, ell, 1000000);
70 36743 : if (!isint1(t)) b = nfmul(nf, b, nfpow_u(nf, t, ell));
71 36743 : if (DEBUGLEVEL>1) err_printf("beta reduced via ell-th root = %Ps\n",b);
72 36743 : b = Q_primitive_part(b, &cb);
73 36743 : if (cb && nfispower(nf, b, ell, NULL)) return cb;
74 36407 : if ((fu = bnf_build_cheapfu(bnfz)))
75 : { /* log. embeddings of fu^ell */
76 36379 : GEN elllogfu = gmulgs(real_i(bnf_get_logfu(bnfz)), ell);
77 36379 : long prec = nf_get_prec(nf);
78 : for (;;)
79 0 : {
80 36379 : GEN ex, y, z = nflogembed(nf, b, NULL, prec);
81 36379 : if (z && (ex = RgM_Babai(elllogfu, z)))
82 : {
83 36379 : if (ZV_equal0(ex)) break;
84 4028 : y = nffactorback(nf, fu, ZC_z_mul(ex,ell));
85 4028 : b = nfdiv(nf, b, y); break;
86 : }
87 0 : prec = precdbl(prec);
88 0 : if (DEBUGLEVEL) pari_warn(warnprec,"reducebeta",prec);
89 0 : nf = nfnewprec_shallow(nf,prec);
90 : }
91 : }
92 36407 : return cb? gmul(b, cb): b;
93 : }
94 :
95 : struct rnfkummer
96 : {
97 : GEN bnfz, cycgenmod, u, vecC, tQ, vecW;
98 : ulong mgi, g, ell;
99 : long rc;
100 : compo_s COMPO;
101 : tau_s tau;
102 : toK_s T;
103 : };
104 :
105 : /* set kum->tau; compute Gal(K(\zeta_l)/K) */
106 : static void
107 2877 : get_tau(struct rnfkummer *kum)
108 : { /* compute action of tau: q^g + kp */
109 2877 : compo_s *C = &kum->COMPO;
110 2877 : GEN U = RgX_add(RgXQ_powu(C->q, kum->g, C->R), RgX_muls(C->p, C->k));
111 2877 : kum->tau.x = RgX_RgXQ_eval(C->rev, U, C->R);
112 2877 : kum->tau.R = C->R;
113 2877 : kum->tau.zk = nfgaloismatrix(bnf_get_nf(kum->bnfz), kum->tau.x);
114 2877 : }
115 :
116 : static GEN tauofvec(GEN x, tau_s *tau);
117 : static GEN
118 282185 : tauofelt(GEN x, tau_s *tau)
119 : {
120 282185 : switch(typ(x))
121 : {
122 16910 : case t_INT: case t_FRAC: return x;
123 248398 : case t_COL: return RgM_RgC_mul(tau->zk, x);
124 16877 : case t_MAT: return mkmat2(tauofvec(gel(x,1), tau), gel(x,2));
125 : default: pari_err_TYPE("tauofelt",x); return NULL;/*LCOV_EXCL_LINE*/
126 : }
127 : }
128 : static GEN
129 18333 : tauofvec(GEN x, tau_s *tau)
130 : {
131 : long i, l;
132 18333 : GEN y = cgetg_copy(x, &l);
133 250305 : for (i=1; i<l; i++) gel(y,i) = tauofelt(gel(x,i), tau);
134 18333 : return y;
135 : }
136 : /* [x, tau(x), ..., tau^(m-1)(x)] */
137 : static GEN
138 5992 : powtau(GEN x, long m, tau_s *tau)
139 : {
140 5992 : GEN y = cgetg(m+1, t_VEC);
141 : long i;
142 5992 : gel(y,1) = x;
143 14532 : for (i=2; i<=m; i++) gel(y,i) = tauofelt(gel(y,i-1), tau);
144 5992 : return y;
145 : }
146 : /* x^lambda */
147 : static GEN
148 9233 : lambdaofelt(GEN x, toK_s *T)
149 : {
150 9233 : tau_s *tau = T->tau;
151 9233 : long i, m = T->m;
152 9233 : GEN y = trivial_fact(), powg = T->powg; /* powg[i] = g^i */
153 24444 : for (i=1; i<m; i++)
154 : {
155 15211 : y = famat_mulpows_shallow(y, x, uel(powg,m-i+1));
156 15211 : x = tauofelt(x, tau);
157 : }
158 9233 : return famat_mul_shallow(y, x);
159 : }
160 : static GEN
161 2996 : lambdaofvec(GEN x, toK_s *T)
162 : {
163 : long i, l;
164 2996 : GEN y = cgetg_copy(x, &l);
165 9723 : for (i=1; i<l; i++) gel(y,i) = lambdaofelt(gel(x,i), T);
166 2996 : return y;
167 : }
168 :
169 : static GEN
170 1071 : tauofideal(GEN id, tau_s *tau)
171 : {
172 1071 : return ZM_hnfmodid(RgM_mul(tau->zk, id), gcoeff(id, 1,1));
173 : }
174 :
175 : static int
176 6699 : prconj(GEN P, GEN Q, tau_s *tau)
177 : {
178 6699 : GEN p = pr_get_p(P), x = pr_get_gen(P);
179 : for(;;)
180 : {
181 20128 : if (ZC_prdvd(x,Q)) return 1;
182 15395 : x = FpC_red(tauofelt(x, tau), p);
183 15396 : if (ZC_prdvd(x,P)) return 0;
184 : }
185 : }
186 : static int
187 91761 : prconj_in_list(GEN S, GEN P, tau_s *tau)
188 : {
189 : long i, l, e, f;
190 : GEN p, x;
191 91761 : if (!tau) return 0;
192 10808 : p = pr_get_p(P); x = pr_get_gen(P);
193 10808 : e = pr_get_e(P); f = pr_get_f(P); l = lg(S);
194 13153 : for (i = 1; i < l; i++)
195 : {
196 7077 : GEN Q = gel(S, i);
197 7077 : if (equalii(p, pr_get_p(Q)) && e == pr_get_e(Q) && f == pr_get_f(Q))
198 6699 : if (ZV_equal(x, pr_get_gen(Q)) || prconj(gel(S,i), P, tau)) return 1;
199 : }
200 6076 : return 0;
201 : }
202 :
203 : /* >= ell */
204 : static long
205 42866 : get_z(GEN pr, long ell) { return ell * (pr_get_e(pr) / (ell-1)); }
206 : /* zeta_ell in nfz */
207 : static void
208 36743 : list_Hecke(GEN *pSp, GEN *pvsprk, GEN nfz, GEN fa, GEN gell, tau_s *tau)
209 : {
210 36743 : GEN P = gel(fa,1), E = gel(fa,2), faell, Sl, S, Sl1, Sl2, Vl, Vl2;
211 36743 : long i, l = lg(P), ell = gell[2];
212 :
213 36743 : S = vectrunc_init(l);
214 36743 : Sl1= vectrunc_init(l);
215 36743 : Sl2= vectrunc_init(l);
216 36743 : Vl2= vectrunc_init(l);
217 101647 : for (i = 1; i < l; i++)
218 : {
219 64904 : GEN pr = gel(P,i);
220 64904 : if (!equaliu(pr_get_p(pr), ell))
221 53991 : { if (!prconj_in_list(S,pr,tau)) vectrunc_append(S,pr); }
222 : else
223 : { /* pr | ell */
224 10913 : long a = get_z(pr, ell) + 1 - itou(gel(E,i));
225 10913 : if (!a)
226 3178 : { if (!prconj_in_list(Sl1,pr,tau)) vectrunc_append(Sl1, pr); }
227 7735 : else if (a != 1 && !prconj_in_list(Sl2,pr,tau))
228 : {
229 2625 : vectrunc_append(Sl2, pr);
230 2625 : vectrunc_append(Vl2, log_prk_init(nfz, pr, a, gell));
231 : }
232 : }
233 : }
234 36743 : faell = idealprimedec(nfz, gell); l = lg(faell);
235 36743 : Vl = vectrunc_init(l);
236 36743 : Sl = vectrunc_init(l);
237 79618 : for (i = 1; i < l; i++)
238 : {
239 42875 : GEN pr = gel(faell,i);
240 42875 : if (!tablesearch(P, pr, cmp_prime_ideal) && !prconj_in_list(Sl, pr, tau))
241 : {
242 31955 : vectrunc_append(Sl, pr);
243 31954 : vectrunc_append(Vl, log_prk_init(nfz, pr, get_z(pr,ell), gell));
244 : }
245 : }
246 36743 : *pvsprk = shallowconcat(Vl2, Vl); /* divide ell */
247 36742 : *pSp = shallowconcat(S, Sl1);
248 36742 : }
249 :
250 : /* Return a Flm, sprk mod pr^k, pr | ell, k >= 2 */
251 : static GEN
252 34580 : logall(GEN nf, GEN v, long lW, long mgi, GEN gell, GEN sprk)
253 : {
254 34580 : long i, ell = gell[2], l = lg(v);
255 34580 : GEN M = cgetg(l,t_MAT);
256 139827 : for (i = 1; i < l; i++)
257 : {
258 105249 : GEN c = log_prk(nf, gel(v,i), sprk, gell); /* ell-rank = #c */
259 105240 : c = ZV_to_Flv(c, ell);
260 105247 : if (i < lW) c = Flv_Fl_mul(c, mgi, ell);
261 105247 : gel(M,i) = c;
262 : }
263 34578 : return M;
264 : }
265 : static GEN
266 36743 : matlogall(GEN nf, GEN v, long lW, long mgi, GEN gell, GEN vsprk)
267 : {
268 36743 : GEN M = NULL;
269 36743 : long i, l = lg(vsprk);
270 71323 : for (i = 1; i < l; i++)
271 34580 : M = vconcat(M, logall(nf, v, lW, mgi, gell, gel(vsprk,i)));
272 36743 : return M;
273 : }
274 :
275 : /* id = (b) prod_{i <= rc} bnfz.gen[i]^v[i] (mod K^*)^ell,
276 : * - i <= rc: gen[i]^cyc[i] = (cycgenmod[i]); ell | cyc[i]
277 : * - i > rc: gen[i]^(u[i]*cyc[i]) = (cycgenmod[i]); u[i] cyc[i] = 1 mod ell */
278 : static void
279 53522 : isprincipalell(GEN bnfz, GEN id, GEN cycgenmod, ulong ell, long rc,
280 : GEN *pv, GEN *pb)
281 : {
282 53522 : long i, l = lg(cycgenmod);
283 53522 : GEN y = bnfisprincipal0(bnfz, id, nf_FORCE|nf_GENMAT);
284 53522 : GEN v = ZV_to_Flv(gel(y,1), ell), b = gel(y,2);
285 54467 : for (i = rc+1; i < l; i++)
286 945 : b = famat_mulpows_shallow(b, gel(cycgenmod,i), v[i]);
287 53522 : setlg(v,rc+1); *pv = v; *pb = b;
288 53522 : }
289 :
290 : static GEN
291 36743 : compute_beta(GEN X, GEN vecWB, GEN ell, GEN bnfz)
292 : {
293 36743 : GEN be = famat_reduce(famatV_zv_factorback(vecWB, X));
294 36743 : if (typ(be) == t_MAT)
295 : {
296 36722 : gel(be,2) = centermod(gel(be,2), ell);
297 36722 : be = nffactorback(bnfz, be, NULL);
298 : }
299 36743 : be = reducebeta(bnfz, be, itou(ell));
300 36742 : if (DEBUGLEVEL>1) err_printf("beta reduced = %Ps\n",be);
301 36742 : return be;
302 : }
303 :
304 : GEN
305 67949 : lift_if_rational(GEN x)
306 : {
307 : long lx, i;
308 : GEN y;
309 :
310 67949 : switch(typ(x))
311 : {
312 9611 : default: break;
313 :
314 38409 : case t_POLMOD:
315 38409 : y = gel(x,2);
316 38409 : if (typ(y) == t_POL)
317 : {
318 12570 : long d = degpol(y);
319 12570 : if (d > 0) return x;
320 2758 : return (d < 0)? gen_0: gel(y,2);
321 : }
322 25839 : return y;
323 :
324 8477 : case t_POL: lx = lg(x);
325 33460 : for (i=2; i<lx; i++) gel(x,i) = lift_if_rational(gel(x,i));
326 8477 : break;
327 11452 : case t_VEC: case t_COL: case t_MAT: lx = lg(x);
328 47187 : for (i=1; i<lx; i++) gel(x,i) = lift_if_rational(gel(x,i));
329 : }
330 29540 : return x;
331 : }
332 :
333 : /* lift elt t in nf to nfz, algebraic form */
334 : static GEN
335 889 : lifttoKz(GEN nf, GEN t, compo_s *C)
336 : {
337 889 : GEN x = nf_to_scalar_or_alg(nf, t);
338 889 : if (typ(x) != t_POL) return x;
339 889 : return RgX_RgXQ_eval(x, C->p, C->R);
340 : }
341 : /* lift ideal id in nf to nfz */
342 : static GEN
343 2996 : ideallifttoKz(GEN nfz, GEN nf, GEN id, compo_s *C)
344 : {
345 2996 : GEN I = idealtwoelt(nf,id);
346 2996 : GEN x = nf_to_scalar_or_alg(nf, gel(I,2));
347 2996 : if (typ(x) != t_POL) return gel(I,1);
348 2128 : gel(I,2) = algtobasis(nfz, RgX_RgXQ_eval(x, C->p, C->R));
349 2128 : return idealhnf_two(nfz,I);
350 : }
351 :
352 : static GEN
353 959 : prlifttoKz_i(GEN nfz, GEN nf, GEN pr, compo_s *C)
354 : {
355 959 : GEN p = pr_get_p(pr), T = nf_get_pol(nfz);
356 959 : if (nf_get_degree(nf) != 1)
357 : { /* restrict to primes above pr */
358 889 : GEN t = pr_get_gen(pr);
359 889 : t = Q_primpart( lifttoKz(nf,t,C) );
360 889 : T = FpX_gcd(FpX_red(T,p), FpX_red(t,p), p);
361 889 : T = FpX_normalize(T, p);
362 : }
363 959 : return gel(FpX_factor(T, p), 1);
364 : }
365 : /* lift ideal pr in nf to ONE prime in nfz (the others are conjugate under tau
366 : * and bring no further information on e_1 W). Assume pr coprime to
367 : * index of both nf and nfz, and unramified in Kz/K (minor simplification) */
368 : static GEN
369 420 : prlifttoKz(GEN nfz, GEN nf, GEN pr, compo_s *C)
370 : {
371 420 : GEN P = prlifttoKz_i(nfz, nf, pr, C);
372 420 : return idealprimedec_kummer(nfz, gel(P,1), pr_get_e(pr), pr_get_p(pr));
373 : }
374 : static GEN
375 539 : prlifttoKzall(GEN nfz, GEN nf, GEN pr, compo_s *C)
376 : {
377 539 : GEN P = prlifttoKz_i(nfz, nf, pr, C), p = pr_get_p(pr), vP;
378 539 : long l = lg(P), e = pr_get_e(pr), i;
379 539 : vP = cgetg(l, t_VEC);
380 2037 : for (i = 1; i < l; i++)
381 1498 : gel(vP,i) = idealprimedec_kummer(nfz,gel(P,i), e, p);
382 539 : return vP;
383 : }
384 :
385 : static GEN
386 39739 : get_badbnf(GEN bnf)
387 : {
388 : long i, l;
389 39739 : GEN bad = gen_1, gen = bnf_get_gen(bnf);
390 39739 : l = lg(gen);
391 44842 : for (i = 1; i < l; i++)
392 : {
393 5103 : GEN g = gel(gen,i);
394 5103 : bad = lcmii(bad, gcoeff(g,1,1));
395 : }
396 39739 : return bad;
397 : }
398 : /* test whether H has index p */
399 : static int
400 52702 : H_is_good(GEN H, GEN p)
401 : {
402 52702 : long l = lg(H), status = 0, i;
403 121927 : for (i = 1; i < l; i++)
404 84968 : if (equalii(gcoeff(H,i,i), p) && ++status > 1) return 0;
405 36959 : return status == 1;
406 : }
407 : static GEN
408 36806 : bid_primes(GEN bid) { return prV_primes(gel(bid_get_fact(bid),1)); }
409 : /* Let K base field, L/K described by bnr (conductor F) + H. Return a list of
410 : * primes coprime to f*ell of degree 1 in K whose images in Cl_f(K) together
411 : * with ell*Cl_f(K), generate H:
412 : * thus they all split in Lz/Kz; t in Kz is such that
413 : * t^(1/p) generates Lz => t is an ell-th power in k(pr) for all such primes.
414 : * Restrict to primes not dividing
415 : * - the index of the polynomial defining Kz,
416 : * - the modulus,
417 : * - ell,
418 : * - a generator in bnf.gen or bnfz.gen.
419 : * If ell | F and Kz != K, set compute the congruence group Hz over Kz
420 : * and set *pfa to the conductor factorization. */
421 : static GEN
422 36743 : get_prlist(GEN bnr, GEN H, GEN gell, GEN *pfa, struct rnfkummer *kum)
423 : {
424 36743 : pari_sp av0 = avma;
425 36743 : GEN Hz = NULL, bnrz = NULL, cycz = NULL, nfz = NULL;
426 36743 : GEN cyc = bnr_get_cyc(bnr), nf = bnr_get_nf(bnr), F = gel(bnr_get_mod(bnr),1);
427 36743 : GEN bad, Hsofar, L = cgetg(1, t_VEC);
428 : forprime_t T;
429 36743 : ulong p, ell = gell[2];
430 36743 : int Ldone = 0;
431 :
432 36743 : bad = get_badbnf(bnr_get_bnf(bnr));
433 36743 : if (kum)
434 : {
435 2996 : GEN bnfz = kum->bnfz, ideal = gel(bnr_get_mod(bnr), 1);
436 2996 : GEN badz = lcmii(get_badbnf(bnfz), nf_get_index(bnf_get_nf(bnfz)));
437 2996 : bad = lcmii(bad,badz);
438 2996 : nfz = bnf_get_nf(bnfz);
439 2996 : ideal = ideallifttoKz(nfz, nf, ideal, &kum->COMPO);
440 2996 : *pfa = idealfactor_partial(nfz, ideal, bid_primes(bnr_get_bid(bnr)));
441 2996 : if (dvdiu(idealdown(nf, ideal), ell))
442 : { /* ell | N(ideal), need Hz = Ker N: Cl_Kz(bothz) -> Cl_K(ideal)/H
443 : * to update conductor */
444 357 : bnrz = Buchraymod(bnfz, *pfa, nf_INIT, gell);
445 357 : cycz = bnr_get_cyc(bnrz);
446 357 : Hz = diagonal_shallow(ZV_snf_gcd(cycz, gell));
447 357 : if (H_is_good(Hz, gell))
448 : {
449 112 : *pfa = gel(bnrconductor_factored(bnrz, Hz), 2);
450 112 : return gc_all(av0, 2, &L, pfa);
451 : }
452 : }
453 : }
454 36631 : bad = lcmii(gcoeff(F,1,1), bad);
455 36631 : cyc = ZV_snf_gcd(cyc, gell);
456 36629 : Hsofar = diagonal_shallow(cyc);
457 36631 : if (H_is_good(Hsofar, gell))
458 : {
459 25067 : if (!Hz) return gc_all(av0, pfa? 2: 1, &L, pfa);
460 119 : Ldone = 1;
461 : }
462 : /* restrict to primes not dividing bad and 1 mod ell */
463 11683 : u_forprime_arith_init(&T, 2, ULONG_MAX, 1, ell);
464 60108 : while ((p = u_forprime_next(&T)))
465 : {
466 : GEN LP;
467 : long i, l;
468 60107 : if (!umodiu(bad, p)) continue;
469 50896 : LP = idealprimedec_limit_f(nf, utoipos(p), 1);
470 50897 : l = lg(LP);
471 74499 : for (i = 1; i < l; i++)
472 : {
473 35285 : pari_sp av = avma;
474 35285 : GEN M, P = gel(LP,i), v = bnrisprincipalmod(bnr, P, gell, 0);
475 35284 : if (!hnf_invimage(H, v)) { set_avma(av); continue; }
476 : /* P in H */
477 17465 : M = ZM_hnfmodid(shallowconcat(Hsofar, v), cyc);
478 17465 : if (Hz)
479 : { /* N_{Kz/K} P in H hence P in Hz */
480 539 : GEN vP = prlifttoKzall(nfz, nf, P, &kum->COMPO);
481 539 : long j, lv = lg(vP);
482 1435 : for (j = 1; j < lv; j++)
483 : {
484 1141 : v = bnrisprincipalmod(bnrz, gel(vP,j), gell, 0);
485 1141 : if (!ZV_equal0(v))
486 : {
487 1141 : Hz = ZM_hnfmodid(shallowconcat(Hz,v), cycz);
488 1141 : if (H_is_good(Hz, gell))
489 : {
490 245 : *pfa = gel(bnrconductor_factored(bnrz, Hz), 2);
491 245 : if (!Ldone) L = vec_append(L, gel(vP,1));
492 245 : return gc_all(av0, 2, &L, pfa);
493 : }
494 : }
495 : }
496 294 : P = gel(vP,1);
497 : }
498 16926 : else if (kum) P = prlifttoKz(nfz, nf, P, &kum->COMPO);
499 17220 : if (Ldone || ZM_equal(M, Hsofar)) continue;
500 14574 : L = vec_append(L, P);
501 14574 : Hsofar = M;
502 14574 : if (H_is_good(Hsofar, gell))
503 : {
504 11536 : if (!Hz) return gc_all(av0, pfa? 2: 1, &L, pfa);
505 98 : Ldone = 1;
506 : }
507 : }
508 : }
509 : pari_err_BUG("rnfkummer [get_prlist]"); return NULL;/*LCOV_EXCL_LINE*/
510 : }
511 : /*Lprz list of prime ideals in Kz that must split completely in Lz/Kz, vecWA
512 : * generators for the S-units used to build the Kummer generators. Return
513 : * matsmall M such that \prod WA[j]^x[j] ell-th power mod pr[i] iff
514 : * \sum M[i,j] x[j] = 0 (mod ell) */
515 : static GEN
516 36743 : subgroup_info(GEN bnfz, GEN Lprz, GEN gell, GEN vecWA)
517 : {
518 36743 : GEN M, nfz = bnf_get_nf(bnfz), Lell = mkvec(gell);
519 36743 : long i, j, ell = gell[2], l = lg(vecWA), lz = lg(Lprz);
520 36743 : M = cgetg(l, t_MAT);
521 146544 : for (j=1; j<l; j++) gel(M,j) = cgetg(lz, t_VECSMALL);
522 51344 : for (i=1; i < lz; i++)
523 : {
524 14602 : GEN pr = gel(Lprz,i), EX = subiu(pr_norm(pr), 1);
525 14602 : GEN N, g,T,p, prM = idealhnf_shallow(nfz, pr);
526 14602 : GEN modpr = zk_to_Fq_init(nfz, &pr,&T,&p);
527 14602 : long v = Z_lvalrem(divis(EX,ell), ell, &N) + 1; /* Norm(pr)-1 = N * ell^v */
528 14602 : GEN ellv = powuu(ell, v);
529 14602 : g = gener_Fq_local(T,p, Lell);
530 14602 : g = Fq_pow(g,N, T,p); /* order ell^v */
531 62950 : for (j=1; j < l; j++)
532 : {
533 48348 : GEN logc, c = gel(vecWA,j);
534 48348 : if (typ(c) == t_MAT) /* famat */
535 34019 : c = famat_makecoprime(nfz, gel(c,1), gel(c,2), pr, prM, EX);
536 48349 : c = nf_to_Fq(nfz, c, modpr);
537 48348 : c = Fq_pow(c, N, T,p);
538 48348 : logc = Fq_log(c, g, ellv, T,p);
539 48348 : ucoeff(M, i,j) = umodiu(logc, ell);
540 : }
541 : }
542 36742 : return M;
543 : }
544 :
545 : static GEN
546 36623 : futu(GEN bnf)
547 : {
548 36623 : GEN fu, tu, SUnits = bnf_get_sunits(bnf);
549 36623 : if (SUnits)
550 : {
551 36098 : tu = nf_to_scalar_or_basis(bnf_get_nf(bnf), bnf_get_tuU(bnf));
552 36099 : fu = bnf_compactfu(bnf);
553 : }
554 : else
555 : {
556 525 : GEN U = bnf_build_units(bnf);
557 525 : tu = gel(U,1); fu = vecslice(U, 2, lg(U)-1);
558 : }
559 36624 : return vec_append(fu, tu);
560 : }
561 : static GEN
562 36623 : bnf_cycgenmod(GEN bnf, long ell, GEN *pselmer, long *prc)
563 : {
564 36623 : GEN gen = bnf_get_gen(bnf), cyc = bnf_get_cyc(bnf), nf = bnf_get_nf(bnf);
565 36623 : GEN B, r = ZV_to_Flv(cyc, ell);
566 36623 : long i, rc, l = lg(gen);
567 36623 : B = cgetg(l, t_VEC);
568 39339 : for (i = 1; i < l && !r[i]; i++);
569 36623 : *prc = rc = i-1; /* ell-rank */
570 40235 : for (i = 1; i < l; i++)
571 : {
572 3612 : GEN G, q, c = gel(cyc,i), g = gel(gen,i);
573 3612 : if (i > rc && r[i] != 1) c = muliu(c, Fl_inv(r[i], ell));
574 3612 : q = divis(c, ell); /* remainder = 0 (i <= rc) or 1 */
575 : /* compute (b) = g^c mod ell-th powers */
576 3612 : G = equali1(q)? g: idealpowred(nf, g, q); /* lose principal part */
577 3612 : G = idealpows(nf, G, ell);
578 3612 : if (i > rc) G = idealmul(nf, G, g);
579 3612 : gel(B,i) = gel(bnfisprincipal0(bnf, G, nf_GENMAT|nf_FORCE), 2);
580 : }
581 36623 : *pselmer = shallowconcat(futu(bnf), vecslice(B,1,rc));
582 36624 : return B;
583 : }
584 :
585 : static GEN
586 33747 : rnfkummersimple(GEN bnr, GEN H, long ell)
587 : {
588 : long j, lSp, rc;
589 : GEN bnf, nf,bid, cycgenmod, Sp, vsprk, matP;
590 33747 : GEN be, M, K, vecW, vecWB, vecBp, gell = utoipos(ell);
591 : /* primes landing in H must be totally split */
592 33747 : GEN Lpr = get_prlist(bnr, H, gell, NULL, NULL);
593 :
594 33747 : bnf = bnr_get_bnf(bnr); if (!bnf_get_sunits(bnf)) bnf_build_units(bnf);
595 33747 : nf = bnf_get_nf(bnf);
596 33747 : bid = bnr_get_bid(bnr);
597 33747 : list_Hecke(&Sp, &vsprk, nf, bid_get_fact2(bid), gell, NULL);
598 :
599 33746 : cycgenmod = bnf_cycgenmod(bnf, ell, &vecW, &rc);
600 33747 : lSp = lg(Sp);
601 33747 : vecBp = cgetg(lSp, t_VEC);
602 33747 : matP = cgetg(lSp, t_MAT);
603 83692 : for (j = 1; j < lSp; j++)
604 49945 : isprincipalell(bnf,gel(Sp,j), cycgenmod,ell,rc, &gel(matP,j),&gel(vecBp,j));
605 33747 : vecWB = shallowconcat(vecW, vecBp);
606 :
607 33747 : M = matlogall(nf, vecWB, 0, 0, gell, vsprk);
608 33747 : M = vconcat(M, shallowconcat(zero_Flm(rc,lg(vecW)-1), matP));
609 33747 : M = vconcat(M, subgroup_info(bnf, Lpr, gell, vecWB));
610 33747 : K = Flm_ker(M, ell);
611 33747 : if (ell == 2)
612 : {
613 31290 : GEN msign = nfsign(nf, vecWB), y;
614 31290 : GEN arch = ZV_to_zv(bid_get_arch(bid)); /* the conductor */
615 31290 : msign = Flm_mul(msign, K, 2);
616 31290 : y = Flm_ker(msign, 2);
617 31290 : y = zv_equal0(arch)? gel(y,1): Flm_Flc_invimage(msign, arch, 2);
618 31290 : K = Flm_Flc_mul(K, y, 2);
619 : }
620 : else
621 2457 : K = gel(K,1);
622 33747 : be = compute_beta(K, vecWB, gell, bnf);
623 33746 : be = nf_to_scalar_or_alg(nf, be);
624 33747 : if (typ(be) == t_POL) be = mkpolmod(be, nf_get_pol(nf));
625 33746 : return gsub(pol_xn(ell, 0), be);
626 : }
627 :
628 : static ulong
629 115556 : nf_to_logFl(GEN nf, GEN x, GEN modpr, ulong g, ulong q, ulong ell, ulong p)
630 : {
631 115556 : x = nf_to_Fp_coprime(nf, x, modpr);
632 115556 : return Fl_log(Fl_powu(umodiu(x, p), q, p), g, ell, p);
633 : }
634 : static GEN
635 35560 : nfV_to_logFlv(GEN nf, GEN x, GEN modpr, ulong g, ulong q, ulong ell, ulong p)
636 151116 : { pari_APPLY_long(nf_to_logFl(nf, gel(x,i), modpr, g, q, ell, p)); }
637 :
638 : /* Compute e_1 Cl(K)/Cl(K)^ell. If u = w^ell a virtual unit, compute
639 : * discrete log mod ell on units.gen + bnf.gen (efficient variant of algo
640 : * 5.3.11) by finding ru degree 1 primes Pj coprime to everything, and gj
641 : * in k(Pj)^* of order ell such that
642 : * log_gj(u_i^((Pj.p - 1) / ell) mod Pj), j = 1..ru
643 : * has maximal F_ell rank ru then solve linear system */
644 : static GEN
645 2877 : kervirtualunit(struct rnfkummer *kum, GEN vselmer)
646 : {
647 2877 : GEN bnf = kum->bnfz, cyc = bnf_get_cyc(bnf), nf = bnf_get_nf(bnf);
648 2877 : GEN W, B, vy, vz, M, U1, U2, vtau, vell, SUnits = bnf_get_sunits(bnf);
649 2877 : long i, j, r, l = lg(vselmer), rc = kum->rc, ru = l-1 - rc, ell = kum->ell;
650 2877 : long LIMC = SUnits? itou(gel(SUnits,4)): 1;
651 : ulong p;
652 : forprime_t T;
653 :
654 2877 : vtau = cgetg(l, t_VEC);
655 2877 : vell = cgetg(l, t_VEC);
656 13944 : for (j = 1; j < l; j++)
657 : {
658 11067 : GEN t = gel(vselmer,j);
659 11067 : if (typ(t) == t_MAT)
660 : {
661 : GEN ct;
662 8190 : t = nffactorback(bnf, gel(t,1), ZV_to_Flv(gel(t,2), ell));
663 8190 : t = Q_primitive_part(t, &ct);
664 8190 : if (ct)
665 : {
666 3814 : GEN F = Q_factor(ct);
667 3814 : ct = factorback2(gel(F,1), ZV_to_Flv(gel(F,2), ell));
668 3814 : t = (typ(t) == t_INT)? ct: ZC_Z_mul(t, ct);
669 : }
670 : }
671 11067 : gel(vell,j) = t; /* integral, not too far from primitive */
672 11067 : gel(vtau,j) = tauofelt(t, &kum->tau);
673 : }
674 2877 : U1 = vecslice(vell, 1, ru); /* units */
675 2877 : U2 = vecslice(vell, ru+1, ru+rc); /* cycgen (mod ell-th powers) */
676 2877 : B = nf_get_index(nf); /* bad primes; from 1 to ru are LIMC-units */
677 3948 : for (i = 1; i <= rc; i++) B = mulii(B, nfnorm(nf, gel(U2,i)));
678 2877 : if (LIMC > 1)
679 : {
680 2877 : GEN U, fa = absZ_factor_limit_strict(B, LIMC, &U), P = gel(fa,1);
681 2877 : long lP = lg(P);
682 2877 : B = U? gel(U,1): gen_1;
683 2877 : if (lP > 1 && cmpiu(gel(P,lP-1), LIMC) >= 0) B = mulii(B, gel(P,lP-1));
684 : }
685 2877 : if (is_pm1(B)) B = NULL;
686 2877 : vy = cgetg(l, t_MAT);
687 12873 : for (j = 1; j <= ru; j++) gel(vy,j) = zero_Flv(rc); /* units */
688 3948 : for ( ; j < l; j++)
689 : {
690 1071 : GEN y, w, u = gel(vtau, j); /* virtual unit */
691 1071 : if (!idealispower(nf, u, ell, &w)) pari_err_BUG("kervirtualunit");
692 1071 : y = isprincipal(bnf, w); setlg(y, rc+1);
693 1071 : if (!ZV_equal0(y))
694 2716 : for (i = 1; i <= rc; i++)
695 1645 : gel(y,i) = diviiexact(mului(ell,gel(y,i)), gel(cyc,i));
696 1071 : gel(vy,j) = ZV_to_Flv(y, ell);
697 : }
698 2877 : u_forprime_arith_init(&T, LIMC+1, ULONG_MAX, 1, ell);
699 2877 : M = cgetg(ru+1, t_MAT); r = 1; setlg(M,2);
700 2877 : vz = cgetg(ru+1, t_MAT);
701 12355 : while ((p = u_forprime_next(&T))) if (!B || umodiu(B,p))
702 : {
703 12292 : GEN P = idealprimedec_limit_f(nf, utoipos(p), 1);
704 12292 : long nP = lg(P)-1;
705 12292 : ulong g = rootsof1_Fl(ell, p), q = p / ell; /* (p-1) / ell */
706 24983 : for (i = 1; i <= nP; i++)
707 : {
708 15568 : GEN modpr = zkmodprinit(nf, gel(P,i));
709 : GEN z, v2;
710 15568 : gel(M, r) = nfV_to_logFlv(nf, U1, modpr, g, q, ell, p); /* log futu */
711 15568 : if (Flm_rank(M, ell) < r) continue; /* discard */
712 :
713 9996 : v2 = nfV_to_logFlv(nf, U2, modpr, g, q, ell, p); /* log alpha[1..rc] */
714 9996 : gel(vz, r) = z = nfV_to_logFlv(nf, vtau, modpr, g, q, ell, p);
715 13776 : for (j = ru+1; j < l; j++)
716 3780 : uel(z,j) = Fl_sub(uel(z,j), Flv_dotproduct(v2, gel(vy,j), ell), ell);
717 9996 : if (r == ru) break;
718 7119 : r++; setlg(M, r+1);
719 : }
720 12292 : if (i < nP) break;
721 : }
722 2877 : if (r != ru) pari_err_BUG("kervirtualunit");
723 : /* Solve prod_k U[k]^x[j,k] = vtau[j] / prod_i alpha[i]^vy[j,i] mod (K^*)^ell
724 : * for 1 <= j <= #vtau. I.e. for a fixed j: M x[j] = vz[j] (mod ell) */
725 2877 : M = Flm_inv(Flm_transpose(M), ell);
726 2877 : vz = Flm_transpose(vz); /* now ru x #vtau */
727 13944 : for (j = 1; j < l; j++)
728 11067 : gel(vy,j) = shallowconcat(Flm_Flc_mul(M, gel(vz,j), ell), gel(vy,j));
729 2877 : W = Flm_ker(Flm_Fl_sub(vy, kum->g, ell), ell); l = lg(W);
730 9282 : for (j = 1; j < l; j++)
731 6405 : gel(W,j) = famat_reduce(famatV_zv_factorback(vselmer, gel(W,j)));
732 2877 : settyp(W, t_VEC); return W;
733 : }
734 :
735 : /* - mu_b = sum_{0 <= i < m} floor(r_b r_{m-1-i} / ell) tau^i.
736 : * Note that i is in fact restricted to i < m-1 */
737 : static GEN
738 5768 : get_mmu(long b, GEN r, long ell)
739 : {
740 5768 : long i, m = lg(r)-1;
741 5768 : GEN M = cgetg(m, t_VECSMALL);
742 28224 : for (i = 0; i < m-1; i++) M[i+1] = (r[b + 1] * r[m - i]) / ell;
743 5768 : return M;
744 : }
745 : /* max_b zv_sum(mu_b) < m ell */
746 : static long
747 2548 : max_smu(GEN r, long ell)
748 : {
749 2548 : long i, s = 0, z = vecsmall_max(r), l = lg(r);
750 7534 : for (i = 2; i < l; i++) s += (z * r[i]) / ell;
751 2548 : return s;
752 : }
753 :
754 : /* coeffs(x, a..b) in variable 0 >= varn(x) */
755 : static GEN
756 17976 : split_pol(GEN x, long a, long b)
757 : {
758 17976 : long i, l = degpol(x);
759 17976 : GEN y = x + a, z;
760 :
761 17976 : if (l < b) b = l;
762 17976 : if (a > b || varn(x) != 0) return pol_0(0);
763 17976 : l = b-a + 3;
764 17976 : z = cgetg(l, t_POL); z[1] = x[1];
765 103670 : for (i = 2; i < l; i++) gel(z,i) = gel(y,i);
766 17976 : return normalizepol_lg(z, l);
767 : }
768 :
769 : /* return (ad * z) mod (T^ell - an/ad), assuming deg_T(z) < 2*ell
770 : * allow ad to be NULL (= 1) */
771 : static GEN
772 8988 : mod_Xell_a(GEN z, long ell, GEN an, GEN ad, GEN T)
773 : {
774 8988 : GEN z1 = split_pol(z, ell, degpol(z));
775 8988 : GEN z0 = split_pol(z, 0, ell-1); /* z = v^ell z1 + z0*/
776 8988 : if (ad) z0 = ZXX_Z_mul(z0, ad);
777 8988 : return gadd(z0, ZXQX_ZXQ_mul(z1, an, T));
778 : }
779 : /* D*basistoalg(nfz, c), in variable v. Result is integral */
780 : static GEN
781 8764 : to_alg(GEN nfz, GEN c, GEN D)
782 : {
783 8764 : if (typ(c) != t_COL) return D? mulii(D,c): c;
784 8764 : return RgV_dotproduct(nf_get_zkprimpart(nfz), c);
785 : }
786 : /* assume x in alg form */
787 : static GEN
788 8988 : downtoK(toK_s *T, GEN x)
789 : {
790 8988 : if (typ(x) != t_POL) return x;
791 8988 : x = RgM_RgC_mul(T->invexpoteta1, RgX_to_RgC(x, lg(T->invexpoteta1) - 1));
792 8988 : return mkpolmod(RgV_to_RgX(x, varn(T->polnf)), T->polnf);
793 : }
794 :
795 : /* th. 5.3.5. and prop. 5.3.9. */
796 : static GEN
797 2996 : compute_polrel(struct rnfkummer *kum, GEN be)
798 : {
799 2996 : toK_s *T = &kum->T;
800 2996 : long i, k, MU = 0, ell = kum->ell, m = T->m;
801 2996 : GEN r = Fl_powers(kum->g, m-1, ell); /* r[i+1] = g^i mod ell */
802 : GEN D, S, root, numa, powtau_Ninvbe, Ninvbe, Dinvbe;
803 2996 : GEN C, prim_Rk, C_Rk, prim_root, C_root, mell = utoineg(ell);
804 2996 : GEN nfz = bnf_get_nf(kum->bnfz), Tz = nf_get_pol(nfz), Dz = nf_get_zkden(nfz);
805 : pari_timer ti;
806 :
807 2996 : if (DEBUGLEVEL>1) { err_printf("Computing Newton sums: "); timer_start(&ti); }
808 2996 : if (equali1(Dz)) Dz = NULL;
809 2996 : D = Dz;
810 2996 : Ninvbe = Q_remove_denom(nfinv(nfz, be), &Dinvbe);
811 2996 : powtau_Ninvbe = powtau(Ninvbe, m-1, T->tau);
812 2996 : if (Dinvbe)
813 : {
814 2548 : MU = max_smu(r, ell);
815 2548 : D = mul_denom(Dz, powiu(Dinvbe, MU));
816 : }
817 :
818 2996 : root = cgetg(ell + 2, t_POL); /* compute D*root, will correct at the end */
819 2996 : root[1] = evalsigne(1) | evalvarn(0);
820 2996 : gel(root,2) = gen_0;
821 2996 : gel(root,3) = D? D: gen_1;
822 8988 : for (i = 2; i < ell; i++) gel(root,2+i) = gen_0;
823 8764 : for (i = 1; i < m; i++)
824 : { /* compute (1/be) ^ (-mu) instead of be^mu [mu < 0].
825 : * 1/be = Ninvbe / Dinvbe */
826 5768 : GEN mmu = get_mmu(i, r, ell), t;
827 5768 : t = to_alg(nfz, nffactorback(nfz, powtau_Ninvbe, mmu), Dz);/* Ninvbe^-mu */
828 5768 : if (Dinvbe)
829 : {
830 4986 : long a = MU - zv_sum(mmu);
831 4986 : if (a) t = gmul(t, powiu(Dinvbe, a));
832 : }
833 5768 : gel(root, 2 + r[i+1]) = t; /* root += D * (z_ell*T)^{r_i} be^mu_i */
834 : }
835 2996 : root = ZXX_renormalize(root, ell+2);
836 : /* Other roots are as above with z_ell -> z_ell^j.
837 : * Treat all contents (C_*) and principal parts (prim_*) separately */
838 2996 : prim_root = Q_primitive_part(root, &C_root);
839 2996 : C_root = div_content(C_root, D);
840 :
841 : /* theta^ell = be^( sum tau^a r_{d-1-a} ) = a = numa / Dz */
842 2996 : numa = to_alg(nfz, nffactorback(nfz, powtau(be, m, T->tau),
843 : vecsmall_reverse(r)), Dz);
844 2996 : if (DEBUGLEVEL>1) err_printf("root(%ld) ", timer_delay(&ti));
845 :
846 : /* Compute mod (X^ell - t, nfz.pol) */
847 2996 : C_Rk = C_root; prim_Rk = prim_root;
848 2996 : C = div_content(C_root, Dz);
849 2996 : S = cgetg(ell+3, t_POL); /* Newton sums */
850 2996 : S[1] = evalsigne(1) | evalvarn(0);
851 2996 : gel(S,2) = gen_0;
852 11984 : for (k = 2; k <= ell; k++)
853 : { /* compute the k-th Newton sum; here C_Rk ~ C_root */
854 8988 : pari_sp av = avma;
855 8988 : GEN z, C_z, d, Rk = ZXQX_mul(prim_Rk, prim_root, Tz);
856 8988 : Rk = mod_Xell_a(Rk, ell, numa, Dz, Tz); /* (mod X^ell - a, nfz.pol) */
857 8988 : prim_Rk = Q_primitive_part(Rk, &d); /* d C_root ~ 1 */
858 8988 : C_Rk = mul_content(C_Rk, mul_content(d, C));
859 : /* root^k = prim_Rk * C_Rk */
860 8988 : z = Q_primitive_part(gel(prim_Rk,2), &C_z); /* C_z ~ 1/C_root ~ 1/C_Rk */
861 8988 : z = downtoK(T, z);
862 8988 : C_z = mul_content(mul_content(C_z, C_Rk), mell);
863 8988 : z = gmul(z, C_z); /* C_z ~ 1 */
864 8988 : gerepileall(av, C_Rk? 3: 2, &z, &prim_Rk, &C_Rk);
865 8988 : if (DEBUGLEVEL>1) err_printf("%ld(%ld) ", k, timer_delay(&ti));
866 8988 : gel(S,k+1) = z; /* - Newton sum */
867 : }
868 2996 : gel(S,ell+2) = gen_m1; if (DEBUGLEVEL>1) err_printf("\n");
869 2996 : return RgX_recip(RgXn_expint(S,ell+1));
870 : }
871 :
872 : static void
873 2877 : compositum_red(compo_s *C, GEN P, GEN Q)
874 : {
875 2877 : GEN p, q, a, z = gel(compositum2(P, Q),1);
876 2877 : a = gel(z,1);
877 2877 : p = gel(gel(z,2), 2);
878 2877 : q = gel(gel(z,3), 2);
879 2877 : C->k = itos( gel(z,4) );
880 2877 : z = polredbest(a, nf_ORIG);
881 2877 : C->R = gel(z,1);
882 2877 : a = gel(gel(z,2), 2);
883 2877 : C->p = RgX_RgXQ_eval(p, a, C->R);
884 2877 : C->q = RgX_RgXQ_eval(q, a, C->R);
885 2877 : C->rev = QXQ_reverse(a, C->R);
886 2877 : }
887 :
888 : /* replace P->C^(-deg P) P(xC) for the largest integer C such that coefficients
889 : * remain algebraic integers. Lift *rational* coefficients */
890 : static void
891 2996 : nfX_Z_normalize(GEN nf, GEN P)
892 : {
893 : long i, l;
894 2996 : GEN C, Cj, PZ = cgetg_copy(P, &l);
895 2996 : PZ[1] = P[1];
896 17976 : for (i = 2; i < l; i++) /* minor variation on RgX_to_nfX (create PZ) */
897 : {
898 14980 : GEN z = nf_to_scalar_or_basis(nf, gel(P,i));
899 14980 : if (typ(z) == t_INT)
900 10939 : gel(PZ,i) = gel(P,i) = z;
901 : else
902 4041 : gel(PZ,i) = ZV_content(z);
903 : }
904 2996 : (void)ZX_Z_normalize(PZ, &C);
905 :
906 2996 : if (C == gen_1) return;
907 493 : Cj = C;
908 2342 : for (i = l-2; i > 1; i--)
909 : {
910 1849 : if (i != l-2) Cj = mulii(Cj, C);
911 1849 : gel(P,i) = gdiv(gel(P,i), Cj);
912 : }
913 : }
914 :
915 : /* set kum->vecC, kum->tQ */
916 : static void
917 868 : _rnfkummer_step4(struct rnfkummer *kum, long d, long m)
918 : {
919 868 : long i, j, rc = kum->rc;
920 868 : GEN Q, vT, vB, vC, vz, B = cgetg(rc+1,t_VEC), T = cgetg(rc+1,t_MAT);
921 868 : GEN gen = bnf_get_gen(kum->bnfz), cycgenmod = kum->cycgenmod;
922 868 : ulong ell = kum->ell;
923 :
924 1939 : for (j = 1; j <= rc; j++)
925 : {
926 1071 : GEN t = tauofideal(gel(gen,j), &kum->tau);
927 1071 : isprincipalell(kum->bnfz, t, cycgenmod,ell,rc, &gel(T,j), &gel(B,j));
928 : }
929 868 : Q = Flm_ker(Flm_Fl_sub(Flm_transpose(T), kum->g, ell), ell);
930 868 : kum->tQ = lg(Q) == 1? NULL: Flm_transpose(Q);
931 868 : kum->vecC = vC = cgetg(rc+1, t_VEC);
932 : /* T = rc x rc matrix */
933 868 : vT = Flm_powers(T, m-2, ell);
934 868 : vB = cgetg(m, t_VEC);
935 868 : vz = cgetg(rc+1, t_VEC);
936 1939 : for (i = 1; i <= rc; i++) gel(vz, i) = cgetg(m, t_VEC);
937 2324 : for (j = 1; j < m; j++)
938 : {
939 1456 : GEN Tj = Flm_Fl_mul(gel(vT,m-j), Fl_mul(j,d,ell), ell);
940 1456 : gel(vB, j) = tauofvec(j == 1? B: gel(vB, j-1), &kum->tau);
941 3241 : for (i = 1; i <= rc; i++) gmael(vz, i, j) = gel(Tj, i);
942 : }
943 868 : vB = shallowconcat1(vB);
944 1939 : for (i = 1; i <= rc; i++)
945 : {
946 1071 : GEN z = shallowconcat1(gel(vz,i));
947 1071 : gel(vC,i) = famat_reduce(famatV_zv_factorback(vB, z));
948 : }
949 868 : }
950 :
951 : /* alg 5.3.5 */
952 : static void
953 2877 : rnfkummer_init(struct rnfkummer *kum, GEN bnf, GEN P, ulong ell, long prec)
954 : {
955 2877 : compo_s *COMPO = &kum->COMPO;
956 2877 : toK_s *T = &kum->T;
957 2877 : GEN nf = bnf_get_nf(bnf), polnf = nf_get_pol(nf), vselmer, bnfz, nfz;
958 : long degK, degKz, m, d;
959 : ulong g;
960 : pari_timer ti;
961 2877 : if (DEBUGLEVEL>2) err_printf("Step 1\n");
962 2877 : if (DEBUGLEVEL) timer_start(&ti);
963 2877 : compositum_red(COMPO, polnf, polcyclo(ell, varn(polnf)));
964 2877 : if (DEBUGLEVEL)
965 : {
966 0 : timer_printf(&ti, "[rnfkummer] compositum");
967 0 : if (DEBUGLEVEL>1) err_printf("polred(compositum) = %Ps\n",COMPO->R);
968 : }
969 2877 : if (DEBUGLEVEL>2) err_printf("Step 2\n");
970 2877 : degK = degpol(polnf);
971 2877 : degKz = degpol(COMPO->R);
972 2877 : m = degKz / degK; /* > 1 */
973 2877 : d = (ell-1) / m;
974 2877 : g = Fl_powu(pgener_Fl(ell), d, ell);
975 2877 : if (Fl_powu(g, m, ell*ell) == 1) g += ell;
976 : /* ord(g) = m in all (Z/ell^k)^* */
977 2877 : if (DEBUGLEVEL>2) err_printf("Step 3\n");
978 2877 : nfz = nfinit(mkvec2(COMPO->R, P), prec);
979 2877 : if (lg(nfcertify(nfz)) > 1) nfz = nfinit(COMPO->R, prec); /* paranoia */
980 2877 : kum->bnfz = bnfz = Buchall(nfz, nf_FORCE, prec);
981 2877 : if (DEBUGLEVEL) timer_printf(&ti, "[rnfkummer] bnfinit(Kz)");
982 2877 : kum->cycgenmod = bnf_cycgenmod(bnfz, ell, &vselmer, &kum->rc);
983 2877 : kum->ell = ell;
984 2877 : kum->g = g;
985 2877 : kum->mgi = Fl_div(m, g, ell);
986 2877 : get_tau(kum);
987 2877 : if (DEBUGLEVEL>2) err_printf("Step 4\n");
988 2877 : if (kum->rc)
989 868 : _rnfkummer_step4(kum, d, m);
990 : else
991 2009 : { kum->vecC = cgetg(1, t_VEC); kum->tQ = NULL; }
992 2877 : if (DEBUGLEVEL>2) err_printf("Step 5\n");
993 2877 : kum->vecW = kervirtualunit(kum, vselmer);
994 2877 : if (DEBUGLEVEL>2) err_printf("Step 8\n");
995 : /* left inverse */
996 2877 : T->invexpoteta1 = QM_inv(RgXQ_matrix_pow(COMPO->p, degKz, degK, COMPO->R));
997 2877 : T->polnf = polnf;
998 2877 : T->tau = &kum->tau;
999 2877 : T->m = m;
1000 2877 : T->powg = Fl_powers(g, m, ell);
1001 2877 : }
1002 :
1003 : static GEN
1004 2996 : rnfkummer_ell(struct rnfkummer *kum, GEN bnr, GEN H)
1005 : {
1006 2996 : ulong ell = kum->ell;
1007 2996 : GEN bnfz = kum->bnfz, nfz = bnf_get_nf(bnfz), gell = utoipos(ell);
1008 2996 : GEN vecC = kum->vecC, vecW = kum->vecW, cycgenmod = kum->cycgenmod;
1009 2996 : long lW = lg(vecW), rc = kum->rc, j, lSp;
1010 2996 : toK_s *T = &kum->T;
1011 : GEN K, be, P, faFz, vsprk, Sp, vecAp, vecBp, matP, vecWA, vecWB, M, lambdaWB;
1012 : /* primes landing in H must be totally split */
1013 2996 : GEN Lpr = get_prlist(bnr, H, gell, &faFz, kum);
1014 :
1015 2996 : if (DEBUGLEVEL>2) err_printf("Step 9, 10 and 11\n");
1016 2996 : list_Hecke(&Sp, &vsprk, nfz, faFz, gell, T->tau);
1017 :
1018 2996 : if (DEBUGLEVEL>2) err_printf("Step 12\n");
1019 2996 : lSp = lg(Sp);
1020 2996 : vecAp = cgetg(lSp, t_VEC);
1021 2996 : vecBp = cgetg(lSp, t_VEC);
1022 2996 : matP = cgetg(lSp, t_MAT);
1023 5502 : for (j = 1; j < lSp; j++)
1024 : {
1025 : GEN e, a;
1026 2506 : isprincipalell(bnfz, gel(Sp,j), cycgenmod,ell,rc, &e, &a);
1027 2506 : gel(matP,j) = e;
1028 2506 : gel(vecBp,j) = famat_mul_shallow(famatV_zv_factorback(vecC, zv_neg(e)), a);
1029 2506 : gel(vecAp,j) = lambdaofelt(gel(vecBp,j), T);
1030 : }
1031 2996 : if (DEBUGLEVEL>2) err_printf("Step 13\n");
1032 2996 : vecWA = shallowconcat(vecW, vecAp);
1033 2996 : vecWB = shallowconcat(vecW, vecBp);
1034 :
1035 2996 : if (DEBUGLEVEL>2) err_printf("Step 14, 15 and 17\n");
1036 2996 : M = matlogall(nfz, vecWA, lW, kum->mgi, gell, vsprk);
1037 2996 : if (kum->tQ)
1038 : {
1039 322 : GEN QtP = Flm_mul(kum->tQ, matP, ell);
1040 322 : M = vconcat(M, shallowconcat(zero_Flm(lgcols(kum->tQ)-1,lW-1), QtP));
1041 : }
1042 2996 : lambdaWB = shallowconcat(lambdaofvec(vecW, T), vecAp);/*vecWB^lambda*/
1043 2996 : M = vconcat(M, subgroup_info(bnfz, Lpr, gell, lambdaWB));
1044 2996 : if (DEBUGLEVEL>2) err_printf("Step 16\n");
1045 2996 : K = Flm_ker(M, ell);
1046 2996 : if (DEBUGLEVEL>2) err_printf("Step 18\n");
1047 2996 : be = compute_beta(gel(K,1), vecWB, gell, kum->bnfz);
1048 2996 : P = compute_polrel(kum, be);
1049 2996 : nfX_Z_normalize(bnr_get_nf(bnr), P);
1050 2996 : if (DEBUGLEVEL>1) err_printf("polrel(beta) = %Ps\n", P);
1051 2996 : return P;
1052 : }
1053 :
1054 : static void
1055 3920 : bnrclassfield_sanitize(GEN *pbnr, GEN *pH)
1056 : {
1057 : GEN T;
1058 3920 : bnr_subgroup_sanitize(pbnr, pH);
1059 3885 : T = nf_get_pol(bnr_get_nf(*pbnr));
1060 3885 : if (!varn(T)) pari_err_PRIORITY("bnrclassfield", T, "=", 0);
1061 3871 : }
1062 :
1063 : static GEN
1064 966 : _rnfkummer(GEN bnr, GEN H, long prec)
1065 : {
1066 : ulong ell;
1067 : GEN gell, bnf, nf, P;
1068 : struct rnfkummer kum;
1069 :
1070 966 : bnrclassfield_sanitize(&bnr, &H);
1071 959 : gell = H? ZM_det(H): ZV_prod(bnr_get_cyc(bnr));
1072 959 : ell = itou(gell);
1073 959 : if (ell == 1) return pol_x(0);
1074 959 : if (!uisprime(ell)) pari_err_IMPL("rnfkummer for composite relative degree");
1075 959 : if (bnf_get_tuN(bnr_get_bnf(bnr)) % ell == 0)
1076 532 : return rnfkummersimple(bnr, H, ell);
1077 427 : bnf = bnr_get_bnf(bnr); nf = bnf_get_nf(bnf);
1078 427 : P = ZV_union_shallow(nf_get_ramified_primes(nf), mkvec(gell));
1079 427 : rnfkummer_init(&kum, bnf, P, ell, maxss(prec,BIGDEFAULTPREC));
1080 427 : return rnfkummer_ell(&kum, bnr, H);
1081 : }
1082 :
1083 : GEN
1084 700 : rnfkummer(GEN bnr, GEN H, long prec)
1085 700 : { pari_sp av = avma; return gerepilecopy(av, _rnfkummer(bnr, H, prec)); }
1086 :
1087 : /*******************************************************************/
1088 : /* bnrclassfield */
1089 : /*******************************************************************/
1090 :
1091 : /* TODO: could be exported */
1092 : static void
1093 229152 : gsetvarn(GEN x, long v)
1094 : {
1095 : long i;
1096 229152 : switch(typ(x))
1097 : {
1098 1841 : case t_POL: case t_SER:
1099 1841 : setvarn(x, v); return;
1100 0 : case t_LIST:
1101 0 : x = list_data(x); if (!x) return;
1102 : /* fall through t_VEC */
1103 : case t_VEC: case t_COL: case t_MAT:
1104 257313 : for (i = lg(x)-1; i > 0; i--) gsetvarn(gel(x,i), v);
1105 : }
1106 : }
1107 :
1108 : /* emb root of pol as polmod modulo pol2, return relative polynomial */
1109 : static GEN
1110 245 : relative_pol(GEN pol, GEN emb, GEN pol2)
1111 : {
1112 : GEN eqn, polrel;
1113 245 : if (degree(pol)==1) return pol2;
1114 196 : eqn = gsub(liftpol_shallow(emb), pol_x(varn(pol)));
1115 196 : eqn = Q_remove_denom(eqn, NULL);
1116 196 : polrel = nfgcd(pol2, eqn, pol, NULL);
1117 196 : return RgX_Rg_div(polrel, leading_coeff(polrel));
1118 : }
1119 :
1120 : /* pol defines K/nf */
1121 : static GEN
1122 266 : bnrclassfield_tower(GEN bnr, GEN subgroup, GEN TB, GEN p, long finaldeg, long absolute, long prec)
1123 : {
1124 266 : pari_sp av = avma;
1125 : GEN nf, nf2, rnf, bnf, bnf2, bnr2, q, H, dec, cyc, pk, sgpk, pol2, emb, emb2, famod, fa, Lbad;
1126 : long i, r1, ell, sp, spk, last;
1127 : forprime_t iter;
1128 :
1129 266 : bnf = bnr_get_bnf(bnr);
1130 266 : nf = bnf_get_nf(bnf);
1131 266 : rnf = rnfinit0(nf, TB, 1);
1132 266 : nf2 = rnf_build_nfabs(rnf, prec);
1133 266 : gsetvarn(nf2, varn(nf_get_pol(nf)));
1134 :
1135 266 : if (lg(nfcertify(nf2)) > 1)
1136 : {
1137 0 : rnf = rnfinit0(nf, gel(TB,1), 1);
1138 0 : nf2 = rnf_build_nfabs(rnf, prec);
1139 0 : gsetvarn(nf2, varn(nf_get_pol(nf)));
1140 : }
1141 :
1142 266 : r1 = nf_get_r1(nf2);
1143 266 : bnf2 = Buchall(nf2, nf_FORCE, prec);
1144 :
1145 266 : sp = itos(p);
1146 266 : spk = sp * rnf_get_degree(rnf);
1147 266 : pk = stoi(spk);
1148 266 : sgpk = hnfmodid(subgroup,pk);
1149 266 : last = spk==finaldeg;
1150 :
1151 : /* compute conductor */
1152 266 : famod = gel(bid_get_fact2(bnr_get_bid(bnr)),1);
1153 266 : if (lg(famod)==1)
1154 : {
1155 140 : fa = trivial_fact();
1156 140 : Lbad = cgetg(1, t_VECSMALL);
1157 : }
1158 : else
1159 : {
1160 126 : long j=1;
1161 126 : fa = cgetg(3, t_MAT);
1162 126 : gel(fa,1) = cgetg(lg(famod), t_VEC);
1163 126 : Lbad = cgetg(lg(famod), t_VEC);
1164 294 : for(i=1; i<lg(famod); i++)
1165 : {
1166 168 : GEN pr = gel(famod,i);
1167 168 : gmael(fa,1,i) = rnfidealprimedec(rnf, pr);
1168 168 : q = pr_get_p(pr);
1169 168 : if (lgefint(q) == 3) gel(Lbad,j++) = q;
1170 : }
1171 126 : setlg(Lbad,j);
1172 126 : Lbad = ZV_to_zv(ZV_sort_uniq_shallow(Lbad));
1173 126 : gel(fa,1) = shallowconcat1(gel(fa,1));
1174 126 : settyp(gel(fa,1), t_COL);
1175 126 : gel(fa,2) = cgetg(lg(gel(fa,1)), t_COL);
1176 392 : for (i=1; i<lg(gel(fa,1)); i++)
1177 : {
1178 266 : GEN pr = gcoeff(fa,i,1);
1179 266 : long e = equalii(p, pr_get_p(pr))? 1 + (pr_get_e(pr)*sp) / (sp-1): 1;
1180 266 : gcoeff(fa,i,2) = utoipos(e);
1181 : }
1182 : }
1183 266 : bnr2 = Buchraymod(bnf2, mkvec2(fa, const_vec(r1,gen_1)), nf_INIT, pk);
1184 :
1185 : /* compute subgroup */
1186 266 : cyc = bnr_get_cyc(bnr2);
1187 266 : H = Flm_image(zv_diagonal(ZV_to_Flv(cyc,sp)), sp);
1188 266 : u_forprime_init(&iter, 2, ULONG_MAX);
1189 16513 : while ((ell = u_forprime_next(&iter))) if (!zv_search(Lbad, ell))
1190 : {
1191 16366 : dec = idealprimedec_limit_f(nf, utoi(ell), 1);
1192 32298 : for (i=1; i<lg(dec); i++)
1193 : {
1194 15932 : GEN pr = gel(dec,i), Pr = gel(rnfidealprimedec(rnf, pr), 1);
1195 15932 : long f = pr_get_f(Pr) / pr_get_f(pr);
1196 15932 : GEN vpr = FpC_Fp_mul(bnrisprincipalmod(bnr,pr,pk,0), utoi(f), pk);
1197 15932 : if (gequal0(ZC_hnfrem(vpr,sgpk)))
1198 1512 : H = vec_append(H, ZV_to_Flv(bnrisprincipalmod(bnr2,Pr,p,0), sp));
1199 : }
1200 16366 : if (lg(H) > lg(cyc)+3)
1201 : {
1202 266 : H = Flm_image(H, sp);
1203 266 : if (lg(cyc)-lg(H) == 1) break;
1204 : }
1205 : }
1206 266 : H = hnfmodid(shallowconcat(zm_to_ZM(H), diagonal_shallow(cyc)), p);
1207 :
1208 : /* polynomial over nf2 */
1209 266 : pol2 = _rnfkummer(bnr2, H, prec);
1210 : /* absolute polynomial */
1211 266 : pol2 = rnfequation2(nf2, pol2);
1212 266 : emb2 = gel(pol2,2); /* generator of nf2 as polmod modulo pol2 */
1213 266 : pol2 = gel(pol2,1);
1214 : /* polynomial over nf */
1215 266 : if (!absolute || !last)
1216 : {
1217 245 : emb = rnf_get_alpha(rnf); /* generator of nf as polynomial in nf2 */
1218 245 : emb = poleval(emb, emb2); /* generator of nf as polmod modulo pol2 */
1219 245 : pol2 = relative_pol(nf_get_pol(nf), emb, pol2);
1220 : }
1221 266 : if (!last) pol2 = rnfpolredbest(nf, pol2, 0);
1222 :
1223 266 : obj_free(rnf);
1224 266 : pol2 = gerepilecopy(av, pol2);
1225 266 : if (last) return pol2;
1226 56 : TB = mkvec2(pol2, gel(TB,2));
1227 56 : return bnrclassfield_tower(bnr, subgroup, TB, p, finaldeg, absolute, prec);
1228 : }
1229 :
1230 : /* subgroups H_i of bnr s.t. bnr/H_i is cyclic and inter_i H_i = subgroup */
1231 : static GEN
1232 35385 : cyclic_compos(GEN subgroup)
1233 : {
1234 35385 : pari_sp av = avma;
1235 35385 : GEN Ui, L, pe, D = ZM_snf_group(subgroup, NULL, &Ui);
1236 35384 : long i, l = lg(D);
1237 :
1238 35384 : L = cgetg(l, t_VEC);
1239 35384 : if (l == 1) return L;
1240 35384 : pe = gel(D,1);
1241 71167 : for (i = 1; i < l; i++)
1242 35783 : gel(L,i) = hnfmodid(shallowconcat(subgroup, vecsplice(Ui,i)),pe);
1243 35384 : return gerepilecopy(av, L);
1244 : }
1245 :
1246 : /* p prime; set pkum=NULL if p-th root of unity in base field
1247 : * absolute=1 allowed if extension is cyclic with exponent>1 */
1248 : static GEN
1249 35385 : bnrclassfield_primepower(struct rnfkummer *pkum, GEN bnr, GEN subgroup, GEN p,
1250 : GEN P, long absolute, long prec)
1251 : {
1252 35385 : GEN res, subs = cyclic_compos(subgroup);
1253 35384 : long i, l = lg(subs);
1254 :
1255 35384 : res = cgetg(l,t_VEC);
1256 71168 : for (i = 1; i < l; i++)
1257 : {
1258 35783 : GEN H = gel(subs,i), cnd = bnrconductormod(bnr, hnfmodid(H,p), p);
1259 35784 : GEN pol, pe, bnr2 = gel(cnd,2), Hp = gel(cnd,3);
1260 35784 : if (pkum) pol = rnfkummer_ell(pkum, bnr2, Hp);
1261 33215 : else pol = rnfkummersimple(bnr2, Hp, itos(p));
1262 35784 : pe = ZM_det_triangular(H);
1263 35784 : if (!equalii(p,pe))
1264 210 : pol = bnrclassfield_tower(bnr, H, mkvec2(pol,P), p, itos(pe), absolute, prec);
1265 35784 : gel(res,i) = pol;
1266 : }
1267 35385 : return res;
1268 : }
1269 :
1270 : /* partition of v into two subsets whose products are as balanced as possible */
1271 : /* assume v sorted */
1272 : static GEN
1273 133 : vecsmall_balance(GEN v)
1274 : {
1275 : forvec_t T;
1276 133 : GEN xbounds, x, vuniq, mult, ind, prod, prodbest = gen_0, bound,
1277 133 : xbest = NULL, res1, res2;
1278 133 : long i=1, j, k1, k2;
1279 133 : if (lg(v) == 3) return mkvec2(mkvecsmall(1), mkvecsmall(2));
1280 42 : vuniq = cgetg(lg(v), t_VECSMALL);
1281 42 : mult = cgetg(lg(v), t_VECSMALL);
1282 42 : ind = cgetg(lg(v), t_VECSMALL);
1283 42 : vuniq[1] = v[1];
1284 42 : mult[1] = 1;
1285 42 : ind[1] = 1;
1286 161 : for (j=2; j<lg(v); j++)
1287 : {
1288 119 : if (v[j] == vuniq[i]) mult[i]++;
1289 : else
1290 : {
1291 14 : i++;
1292 14 : vuniq[i] = v[j];
1293 14 : mult[i] = 1;
1294 14 : ind[i] = j;
1295 : }
1296 : }
1297 42 : setlg(vuniq, ++i);
1298 42 : setlg(mult, i);
1299 42 : setlg(ind, i);
1300 :
1301 42 : vuniq = zv_to_ZV(vuniq);
1302 42 : prod = factorback2(vuniq, mult);
1303 42 : bound = sqrti(prod);
1304 42 : xbounds = cgetg(lg(mult), t_VEC);
1305 98 : for (i=1; i<lg(mult); i++) gel(xbounds,i) = mkvec2s(0,mult[i]);
1306 :
1307 42 : forvec_init(&T, xbounds, 0);
1308 273 : while ((x = forvec_next(&T)))
1309 : {
1310 231 : prod = factorback2(vuniq, x);
1311 231 : if (cmpii(prod,bound)<=0 && cmpii(prod,prodbest)>0)
1312 : {
1313 105 : prodbest = prod;
1314 105 : xbest = gcopy(x);
1315 : }
1316 : }
1317 42 : res1 = cgetg(lg(v), t_VECSMALL);
1318 42 : res2 = cgetg(lg(v), t_VECSMALL);
1319 98 : for (i=1,k1=1,k2=1; i<lg(xbest); i++)
1320 : {
1321 119 : for (j=0; j<itos(gel(xbest,i)); j++) res1[k1++] = ind[i]+j;
1322 154 : for (; j<mult[i]; j++) res2[k2++] = ind[i]+j;
1323 : }
1324 42 : setlg(res1, k1);
1325 42 : setlg(res2, k2); return mkvec2(res1, res2);
1326 : }
1327 :
1328 : /* TODO nfcompositum should accept vectors of pols */
1329 : /* assume all fields are linearly disjoint */
1330 : /* assume the polynomials are sorted by degree */
1331 : static GEN
1332 448 : nfcompositumall(GEN nf, GEN L)
1333 : {
1334 : GEN pol, vdeg, part;
1335 : long i;
1336 448 : if (lg(L)==2) return gel(L,1);
1337 133 : vdeg = cgetg(lg(L), t_VECSMALL);
1338 476 : for (i=1; i<lg(L); i++) vdeg[i] = degree(gel(L,i));
1339 133 : part = vecsmall_balance(vdeg);
1340 133 : pol = cgetg(3, t_VEC);
1341 399 : for (i = 1; i < 3; i++)
1342 : {
1343 266 : GEN L2 = vecpermute(L, gel(part,i)), T = nfcompositumall(nf, L2);
1344 266 : gel(pol,i) = rnfpolredbest(nf, T, 0);
1345 : }
1346 133 : return nfcompositum(nf, gel(pol,1), gel(pol,2), 2);
1347 : }
1348 :
1349 : static GEN
1350 33810 : disc_primes(GEN bnr)
1351 : {
1352 33810 : GEN v = bid_primes(bnr_get_bid(bnr));
1353 33809 : GEN r = nf_get_ramified_primes(bnr_get_nf(bnr));
1354 33809 : return ZV_union_shallow(r, v);
1355 : }
1356 : static struct rnfkummer **
1357 33788 : rnfkummer_initall(GEN bnr, GEN vP, GEN P, long prec)
1358 : {
1359 33788 : long i, w, l = lg(vP), S = vP[l-1] + 1;
1360 33788 : GEN bnf = bnr_get_bnf(bnr);
1361 : struct rnfkummer **vkum;
1362 :
1363 33788 : w = bnf_get_tuN(bnf);
1364 33788 : vkum = (struct rnfkummer **)stack_malloc(S * sizeof(struct rnfkummer*));
1365 33788 : if (prec < BIGDEFAULTPREC) prec = BIGDEFAULTPREC;
1366 67660 : for (i = 1; i < l; i++)
1367 : {
1368 33872 : long p = vP[i];
1369 33872 : if (w % p == 0) { vkum[p] = NULL; continue; }
1370 2450 : vkum[p] = (struct rnfkummer *)stack_malloc(sizeof(struct rnfkummer));
1371 2450 : rnfkummer_init(vkum[p], bnf, P, p, prec);
1372 : }
1373 33788 : return vkum;
1374 : }
1375 : /* fully handle a single congruence subgroup H */
1376 : static GEN
1377 35342 : bnrclassfield_H(struct rnfkummer **vkum, GEN bnr, GEN bad, GEN H0, GEN fa, long flag,
1378 : long prec)
1379 : {
1380 35342 : GEN PN = gel(fa,1), EN = gel(fa,2), res;
1381 35342 : long i, lPN = lg(PN), absolute;
1382 :
1383 35342 : if (lPN == 1) switch(flag)
1384 : {
1385 14 : case 0:
1386 14 : return mkvec(pol_x(0));
1387 14 : case 1:
1388 14 : return pol_x(0);
1389 14 : default: /* 2 */
1390 14 : res = shallowcopy(nf_get_pol(bnr_get_nf(bnr)));
1391 14 : setvarn(res,0); return res;
1392 : }
1393 35300 : absolute = flag==2 && lPN==2 && !equali1(gel(EN,1)); /* one prime, exponent > 1 */
1394 35300 : res = cgetg(lPN, t_VEC);
1395 70685 : for (i = 1; i < lPN; i++)
1396 : {
1397 35384 : GEN p = gel(PN,i), H = hnfmodid(H0, powii(p, gel(EN,i)));
1398 35385 : long sp = itos(p);
1399 35385 : if (absolute) absolute = FpM_rank(H,p)==lg(H)-2; /* cyclic */
1400 35385 : gel(res,i) = bnrclassfield_primepower(vkum[sp], bnr, H, p, bad, absolute, prec);
1401 : }
1402 35301 : res = liftpol_shallow(shallowconcat1(res));
1403 35301 : res = gen_sort_shallow(res, (void*)cmp_RgX, gen_cmp_RgX);
1404 35301 : if (flag)
1405 : {
1406 182 : GEN nf = bnr_get_nf(bnr);
1407 182 : res = nfcompositumall(nf, res);
1408 182 : if (flag==2 && !absolute) res = rnfequation(nf, res);
1409 : }
1410 35301 : return res;
1411 : }
1412 : /* for a vector of subgroups */
1413 : static GEN
1414 30968 : bnrclassfieldvec(GEN bnr, GEN v, long flag, long prec)
1415 : {
1416 30968 : long j, lv = lg(v);
1417 30968 : GEN vH, vfa, vP, P, w = cgetg(lv, t_VEC);
1418 30967 : struct rnfkummer **vkum = NULL;
1419 :
1420 30967 : if (lv == 1) return w;
1421 30960 : vH = cgetg(lv, t_VEC);
1422 30960 : vP = cgetg(lv, t_VEC);
1423 30960 : vfa = cgetg(lv, t_VEC);
1424 63447 : for (j = 1; j < lv; j++)
1425 : {
1426 32493 : GEN N, fa, H = bnr_subgroup_check(bnr, gel(v,j), &N);
1427 32494 : if (is_bigint(N)) pari_err_OVERFLOW("bnrclassfield [too large degree]");
1428 32487 : if (!H) H = diagonal_shallow(bnr_get_cyc(bnr));
1429 32487 : gel(vH,j) = H;
1430 32487 : gel(vfa,j) = fa = Z_factor(N);
1431 32487 : gel(vP,j) = ZV_to_zv(gel(fa, 1));
1432 : }
1433 30954 : vP = shallowconcat1(vP); vecsmall_sort(vP);
1434 30954 : vP = vecsmall_uniq_sorted(vP);
1435 30954 : P = disc_primes(bnr);
1436 30953 : if (lg(vP) > 1) vkum = rnfkummer_initall(bnr, vP, P, prec);
1437 63440 : for (j = 1; j < lv; j++)
1438 32486 : gel(w,j) = bnrclassfield_H(vkum, bnr, P, gel(vH,j), gel(vfa,j), flag, prec);
1439 30954 : return w;
1440 : }
1441 : /* flag:
1442 : * 0 t_VEC of polynomials whose compositum is the extension
1443 : * 1 single polynomial
1444 : * 2 single absolute polynomial */
1445 : GEN
1446 33936 : bnrclassfield(GEN bnr, GEN subgroup, long flag, long prec)
1447 : {
1448 33936 : pari_sp av = avma;
1449 : GEN N, fa, P;
1450 : struct rnfkummer **vkum;
1451 :
1452 33936 : if (flag<0 || flag>2) pari_err_FLAG("bnrclassfield [must be 0,1 or 2]");
1453 33922 : if (subgroup && typ(subgroup) == t_VEC)
1454 : {
1455 30975 : if (nftyp(bnr)==typ_BNF) bnr = Buchray(bnr, gen_1, nf_INIT);
1456 30940 : else checkbnr(bnr);
1457 30975 : if (!char_check(bnr_get_cyc(bnr), subgroup))
1458 30968 : return gerepilecopy(av, bnrclassfieldvec(bnr, subgroup, flag, prec));
1459 : }
1460 2954 : bnrclassfield_sanitize(&bnr, &subgroup);
1461 2912 : N = ZM_det_triangular(subgroup);
1462 2912 : if (equali1(N)) switch(flag)
1463 : {
1464 28 : case 0: set_avma(av); retmkvec(pol_x(0));
1465 7 : case 1: set_avma(av); return pol_x(0);
1466 7 : default: /* 2 */
1467 7 : P = shallowcopy(nf_get_pol(bnr_get_nf(bnr)));
1468 7 : setvarn(P,0); return gerepilecopy(av,P);
1469 : }
1470 2870 : if (is_bigint(N)) pari_err_OVERFLOW("bnrclassfield [too large degree]");
1471 2856 : fa = Z_factor(N); P = disc_primes(bnr);
1472 2856 : vkum = rnfkummer_initall(bnr, ZV_to_zv(gel(fa,1)), P, prec);
1473 2856 : return gerepilecopy(av, bnrclassfield_H(vkum, bnr, P, subgroup, fa, flag, prec));
1474 : }
|
