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