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 39652 : reduce_mod_Qell(GEN nf, GEN b, ulong ell)
49 : {
50 : GEN c;
51 39652 : b = nf_to_scalar_or_basis(nf, b);
52 39652 : b = Q_primitive_part(b, &c);
53 39652 : if (c)
54 : {
55 11748 : GEN d, fa = Q_factor_limit(c, 1000000);
56 11748 : d = factorback2(gel(fa,1), ZV_to_Flv(gel(fa,2), ell));
57 11748 : b = typ(b) == t_INT? mulii(b,d): ZC_Z_mul(b, d);
58 : }
59 39652 : return b;
60 : }
61 :
62 : static GEN
63 39652 : reducebeta(GEN bnfz, GEN b, long ell)
64 : {
65 39652 : GEN t, cb, fu, nf = bnf_get_nf(bnfz);
66 :
67 39652 : if (DEBUGLEVEL>1) err_printf("reducing beta = %Ps\n",b);
68 39652 : b = reduce_mod_Qell(nf, b, ell);
69 39652 : t = idealredmodpower(nf, b, ell, 1000000);
70 39652 : if (!isint1(t)) b = nfmul(nf, b, nfpow_u(nf, t, ell));
71 39652 : if (DEBUGLEVEL>1) err_printf("beta reduced via ell-th root = %Ps\n",b);
72 39652 : b = Q_primitive_part(b, &cb);
73 39652 : if (cb && nfispower(nf, b, ell, NULL)) return cb;
74 38831 : if ((fu = bnf_build_cheapfu(bnfz)))
75 : { /* log. embeddings of fu^ell */
76 38807 : GEN elllogfu = gmulgs(real_i(bnf_get_logfu(bnfz)), ell);
77 38807 : long prec = nf_get_prec(nf);
78 : for (;;)
79 0 : {
80 38807 : GEN ex, y, z = nflogembed(nf, b, NULL, prec);
81 38807 : if (z && (ex = RgM_Babai(elllogfu, z)))
82 : {
83 38807 : if (ZV_equal0(ex)) break;
84 5434 : y = nffactorback(nf, fu, ZC_z_mul(ex,ell));
85 5434 : 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 38831 : 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 2431 : get_tau(struct rnfkummer *kum)
108 : { /* compute action of tau: q^g + kp */
109 2431 : compo_s *C = &kum->COMPO;
110 2431 : GEN U = RgX_add(RgXQ_powu(C->q, kum->g, C->R), RgX_muls(C->p, C->k));
111 2431 : kum->tau.x = RgX_RgXQ_eval(C->rev, U, C->R);
112 2431 : kum->tau.R = C->R;
113 2431 : kum->tau.zk = nfgaloismatrix(bnf_get_nf(kum->bnfz), kum->tau.x);
114 2431 : }
115 :
116 : static GEN RgV_tau(GEN x, tau_s *tau);
117 : static GEN
118 237488 : Rg_tau(GEN x, tau_s *tau)
119 : {
120 237488 : switch(typ(x))
121 : {
122 14537 : case t_INT: case t_FRAC: return x;
123 208665 : case t_COL: return RgM_RgC_mul(tau->zk, x);
124 14286 : 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 15485 : RgV_tau(GEN x, tau_s *tau)
130 210515 : { pari_APPLY_same(Rg_tau(gel(x,i), tau)); }
131 : /* [x, tau(x), ..., tau^(m-1)(x)] */
132 : static GEN
133 5078 : powtau(GEN x, long m, tau_s *tau)
134 : {
135 5078 : GEN y = cgetg(m+1, t_VEC);
136 : long i;
137 5078 : gel(y,1) = x;
138 12301 : for (i=2; i<=m; i++) gel(y,i) = Rg_tau(gel(y,i-1), tau);
139 5078 : return y;
140 : }
141 : /* x^lambda */
142 : static GEN
143 7862 : Rg_lambda(GEN x, toK_s *T)
144 : {
145 7862 : tau_s *tau = T->tau;
146 7862 : long i, m = T->m;
147 7862 : GEN y = trivial_fact(), powg = T->powg; /* powg[i] = g^i */
148 20772 : for (i=1; i<m; i++)
149 : {
150 12910 : y = famat_mulpows_shallow(y, x, uel(powg,m-i+1));
151 12910 : x = Rg_tau(x, tau);
152 : }
153 7862 : return famat_mul_shallow(y, x);
154 : }
155 : static GEN
156 2539 : RgV_lambda(GEN x, toK_s *T)
157 8244 : { pari_APPLY_same(Rg_lambda(gel(x,i), T)); }
158 :
159 : static int
160 5688 : prconj(GEN P, GEN Q, tau_s *tau)
161 : {
162 5688 : GEN p = pr_get_p(P), x = pr_get_gen(P);
163 : for(;;)
164 : {
165 17026 : if (ZC_prdvd(x,Q)) return 1;
166 12998 : x = FpC_red(Rg_tau(x, tau), p);
167 12998 : if (ZC_prdvd(x,P)) return 0;
168 : }
169 : }
170 : static int
171 99404 : prconj_in_list(GEN S, GEN P, tau_s *tau)
172 : {
173 : long i, l, e, f;
174 : GEN p, x;
175 99404 : if (!tau) return 0;
176 9213 : p = pr_get_p(P); x = pr_get_gen(P);
177 9213 : e = pr_get_e(P); f = pr_get_f(P); l = lg(S);
178 11205 : for (i = 1; i < l; i++)
179 : {
180 6020 : GEN Q = gel(S, i);
181 6020 : if (equalii(p, pr_get_p(Q)) && e == pr_get_e(Q) && f == pr_get_f(Q))
182 5688 : if (ZV_equal(x, pr_get_gen(Q)) || prconj(gel(S,i), P, tau)) return 1;
183 : }
184 5185 : return 0;
185 : }
186 :
187 : /* >= ell */
188 : static long
189 53075 : get_z(GEN pr, long ell) { return ell * (pr_get_e(pr) / (ell-1)); }
190 : /* zeta_ell in nfz */
191 : static void
192 39652 : list_Hecke(GEN *pSp, GEN *pvsprk, GEN nfz, GEN fa, GEN gell, tau_s *tau)
193 : {
194 39652 : GEN P = gel(fa,1), E = gel(fa,2), faell, Sl, S, Sl1, Sl2, Vl, Vl2;
195 39652 : long i, l = lg(P), ell = gell[2];
196 :
197 39652 : S = vectrunc_init(l);
198 39652 : Sl1= vectrunc_init(l);
199 39652 : Sl2= vectrunc_init(l);
200 39652 : Vl2= vectrunc_init(l);
201 115859 : for (i = 1; i < l; i++)
202 : {
203 76207 : GEN pr = gel(P,i);
204 76207 : if (!equaliu(pr_get_p(pr), ell))
205 54925 : { if (!prconj_in_list(S,pr,tau)) vectrunc_append(S,pr); }
206 : else
207 : { /* pr | ell */
208 21282 : long a = get_z(pr, ell) + 1 - itou(gel(E,i));
209 21282 : if (!a)
210 10430 : { if (!prconj_in_list(Sl1,pr,tau)) vectrunc_append(Sl1, pr); }
211 10852 : else if (a != 1 && !prconj_in_list(Sl2,pr,tau))
212 : {
213 2245 : vectrunc_append(Sl2, pr);
214 2245 : vectrunc_append(Vl2, log_prk_init(nfz, pr, a, gell));
215 : }
216 : }
217 : }
218 39652 : faell = idealprimedec(nfz, gell); l = lg(faell);
219 39652 : Vl = vectrunc_init(l);
220 39652 : Sl = vectrunc_init(l);
221 92732 : for (i = 1; i < l; i++)
222 : {
223 53080 : GEN pr = gel(faell,i);
224 53080 : if (!tablesearch(P, pr, cmp_prime_ideal) && !prconj_in_list(Sl, pr, tau))
225 : {
226 31793 : vectrunc_append(Sl, pr);
227 31793 : vectrunc_append(Vl, log_prk_init(nfz, pr, get_z(pr,ell), gell));
228 : }
229 : }
230 39652 : *pvsprk = shallowconcat(Vl2, Vl); /* divide ell */
231 39652 : *pSp = shallowconcat(S, Sl1);
232 39652 : }
233 :
234 : /* Return a Flm, sprk mod pr^k, pr | ell, k >= 2 */
235 : static GEN
236 34038 : logall(GEN nf, GEN v, long lW, long mgi, GEN gell, GEN sprk)
237 : {
238 34038 : long i, ell = gell[2], l = lg(v);
239 34038 : GEN M = cgetg(l,t_MAT);
240 140323 : for (i = 1; i < l; i++)
241 : {
242 106285 : GEN c = log_prk(nf, gel(v,i), sprk, gell); /* ell-rank = #c */
243 106285 : c = ZV_to_Flv(c, ell);
244 106285 : if (i < lW) c = Flv_Fl_mul(c, mgi, ell);
245 106285 : gel(M,i) = c;
246 : }
247 34038 : return M;
248 : }
249 : static GEN
250 39652 : matlogall(GEN nf, GEN v, long lW, long mgi, GEN gell, GEN vsprk)
251 : {
252 39652 : GEN M = NULL;
253 39652 : long i, l = lg(vsprk);
254 73690 : for (i = 1; i < l; i++)
255 34038 : M = vconcat(M, logall(nf, v, lW, mgi, gell, gel(vsprk,i)));
256 39652 : 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 62223 : isprincipalell(GEN bnfz, GEN id, GEN cycgenmod, ulong ell, long rc,
264 : GEN *pv, GEN *pb)
265 : {
266 62223 : long i, l = lg(cycgenmod);
267 62223 : GEN y = bnfisprincipal0(bnfz, id, nf_FORCE|nf_GENMAT);
268 62223 : GEN v = ZV_to_Flv(gel(y,1), ell), b = gel(y,2);
269 63010 : for (i = rc+1; i < l; i++)
270 787 : b = famat_mulpows_shallow(b, gel(cycgenmod,i), v[i]);
271 62223 : setlg(v,rc+1); *pv = v; *pb = b;
272 62223 : }
273 :
274 : static GEN
275 39652 : compute_beta(GEN X, GEN vecWB, GEN ell, GEN bnfz)
276 : {
277 39652 : GEN be = famat_reduce(famatV_zv_factorback(vecWB, X));
278 39652 : if (typ(be) == t_MAT)
279 : {
280 39473 : gel(be,2) = centermod(gel(be,2), ell);
281 39473 : be = nffactorback(bnfz, be, NULL);
282 : }
283 39652 : be = reducebeta(bnfz, be, itou(ell));
284 39652 : if (DEBUGLEVEL>1) err_printf("beta reduced = %Ps\n",be);
285 39652 : return be;
286 : }
287 :
288 : GEN
289 63855 : lift_if_rational(GEN x)
290 : {
291 : long lx, i;
292 : GEN y;
293 :
294 63855 : switch(typ(x))
295 : {
296 8758 : default: break;
297 :
298 35302 : case t_POLMOD:
299 35302 : y = gel(x,2);
300 35302 : if (typ(y) == t_POL)
301 : {
302 10436 : long d = degpol(y);
303 10436 : if (d > 0) return x;
304 2321 : return (d < 0)? gen_0: gel(y,2);
305 : }
306 24866 : return y;
307 :
308 8155 : case t_POL: lx = lg(x);
309 30699 : for (i=2; i<lx; i++) gel(x,i) = lift_if_rational(gel(x,i));
310 8155 : break;
311 11640 : case t_VEC: case t_COL: case t_MAT: lx = lg(x);
312 46042 : for (i=1; i<lx; i++) gel(x,i) = lift_if_rational(gel(x,i));
313 : }
314 28553 : return x;
315 : }
316 :
317 : /* lift elt t in nf to nfz, algebraic form */
318 : static GEN
319 732 : lifttoKz(GEN nf, GEN t, compo_s *C)
320 : {
321 732 : GEN x = nf_to_scalar_or_alg(nf, t);
322 732 : if (typ(x) != t_POL) return x;
323 732 : return RgX_RgXQ_eval(x, C->p, C->R);
324 : }
325 : /* lift ideal id in nf to nfz */
326 : static GEN
327 2539 : ideallifttoKz(GEN nfz, GEN nf, GEN id, compo_s *C)
328 : {
329 2539 : GEN I = idealtwoelt(nf,id);
330 2539 : GEN x = nf_to_scalar_or_alg(nf, gel(I,2));
331 2539 : if (typ(x) != t_POL) return gel(I,1);
332 1794 : gel(I,2) = algtobasis(nfz, RgX_RgXQ_eval(x, C->p, C->R));
333 1794 : return idealhnf_two(nfz,I);
334 : }
335 :
336 : static GEN
337 791 : prlifttoKz_i(GEN nfz, GEN nf, GEN pr, compo_s *C)
338 : {
339 791 : GEN p = pr_get_p(pr), T = nf_get_pol(nfz);
340 791 : if (nf_get_degree(nf) != 1)
341 : { /* restrict to primes above pr */
342 732 : GEN t = pr_get_gen(pr);
343 732 : t = Q_primpart( lifttoKz(nf,t,C) );
344 732 : T = FpX_gcd(FpX_red(T,p), FpX_red(t,p), p);
345 732 : T = FpX_normalize(T, p);
346 : }
347 791 : 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 369 : prlifttoKz(GEN nfz, GEN nf, GEN pr, compo_s *C)
354 : {
355 369 : GEN P = prlifttoKz_i(nfz, nf, pr, C);
356 369 : return idealprimedec_kummer(nfz, gel(P,1), pr_get_e(pr), pr_get_p(pr));
357 : }
358 : static GEN
359 422 : prlifttoKzall(GEN nfz, GEN nf, GEN pr, compo_s *C)
360 : {
361 422 : GEN P = prlifttoKz_i(nfz, nf, pr, C), p = pr_get_p(pr), vP;
362 422 : long l = lg(P), e = pr_get_e(pr), i;
363 422 : vP = cgetg(l, t_VEC);
364 1598 : for (i = 1; i < l; i++)
365 1176 : gel(vP,i) = idealprimedec_kummer(nfz,gel(P,i), e, p);
366 422 : return vP;
367 : }
368 :
369 : static GEN
370 42191 : get_badbnf(GEN bnf)
371 : {
372 : long i, l;
373 42191 : GEN bad = gen_1, gen = bnf_get_gen(bnf);
374 42191 : l = lg(gen);
375 46498 : for (i = 1; i < l; i++)
376 : {
377 4307 : GEN g = gel(gen,i);
378 4307 : bad = lcmii(bad, gcoeff(g,1,1));
379 : }
380 42191 : return bad;
381 : }
382 : /* test whether H has index p */
383 : static int
384 66355 : H_is_good(GEN H, GEN p)
385 : {
386 66355 : long l = lg(H), status = 0, i;
387 161961 : for (i = 1; i < l; i++)
388 122139 : if (equalii(gcoeff(H,i,i), p) && ++status > 1) return 0;
389 39822 : return status == 1;
390 : }
391 : static GEN
392 35590 : 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 39652 : get_prlist(GEN bnr, GEN H, GEN gell, GEN *pfa, struct rnfkummer *kum)
407 : {
408 39652 : pari_sp av0 = avma;
409 39652 : GEN Hz = NULL, bnrz = NULL, cycz = NULL, nfz = NULL;
410 39652 : GEN cyc = bnr_get_cyc(bnr), nf = bnr_get_nf(bnr), F = gel(bnr_get_mod(bnr),1);
411 39652 : GEN bad, Hsofar, L = cgetg(1, t_VEC);
412 : forprime_t T;
413 39652 : ulong p, ell = gell[2];
414 39652 : int Ldone = 0;
415 :
416 39652 : bad = get_badbnf(bnr_get_bnf(bnr));
417 39652 : if (kum)
418 : {
419 2539 : GEN bnfz = kum->bnfz, ideal = gel(bnr_get_mod(bnr), 1);
420 2539 : GEN badz = lcmii(get_badbnf(bnfz), nf_get_index(bnf_get_nf(bnfz)));
421 2539 : bad = lcmii(bad,badz);
422 2539 : nfz = bnf_get_nf(bnfz);
423 2539 : ideal = ideallifttoKz(nfz, nf, ideal, &kum->COMPO);
424 2539 : *pfa = idealfactor_partial(nfz, ideal, bid_primes(bnr_get_bid(bnr)));
425 2539 : 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 290 : bnrz = Buchraymod(bnfz, *pfa, nf_INIT, gell);
429 290 : cycz = bnr_get_cyc(bnrz);
430 290 : Hz = diagonal_shallow(ZV_snf_gcd(cycz, gell));
431 290 : if (H_is_good(Hz, gell))
432 : {
433 96 : *pfa = gel(bnrconductor_factored(bnrz, Hz), 2);
434 96 : return gc_all(av0, 2, &L, pfa);
435 : }
436 : }
437 : }
438 39556 : bad = lcmii(gcoeff(F,1,1), bad);
439 39556 : cyc = ZV_snf_gcd(cyc, gell);
440 39556 : Hsofar = diagonal_shallow(cyc);
441 39556 : if (H_is_good(Hsofar, gell))
442 : {
443 23502 : if (!Hz) return gc_all(av0, pfa? 2: 1, &L, pfa);
444 98 : Ldone = 1;
445 : }
446 : /* restrict to primes not dividing bad and 1 mod ell */
447 16152 : u_forprime_arith_init(&T, 2, ULONG_MAX, 1, ell);
448 113327 : while ((p = u_forprime_next(&T)))
449 : {
450 : GEN LP;
451 : long i, l;
452 113327 : if (!umodiu(bad, p)) continue;
453 100626 : LP = idealprimedec_limit_f(nf, utoipos(p), 1);
454 100626 : l = lg(LP);
455 144042 : for (i = 1; i < l; i++)
456 : {
457 59568 : pari_sp av = avma;
458 59568 : GEN M, P = gel(LP,i), v = bnrisprincipalmod(bnr, P, gell, 0);
459 59568 : if (!hnf_invimage(H, v)) { set_avma(av); continue; }
460 : /* P in H */
461 30478 : M = ZM_hnfmodid(shallowconcat(Hsofar, v), cyc);
462 30478 : if (Hz)
463 : { /* N_{Kz/K} P in H hence P in Hz */
464 422 : GEN vP = prlifttoKzall(nfz, nf, P, &kum->COMPO);
465 422 : long j, lv = lg(vP);
466 1118 : for (j = 1; j < lv; j++)
467 : {
468 890 : v = bnrisprincipalmod(bnrz, gel(vP,j), gell, 0);
469 890 : if (!ZV_equal0(v))
470 : {
471 890 : Hz = ZM_hnfmodid(shallowconcat(Hz,v), cycz);
472 890 : if (H_is_good(Hz, gell))
473 : {
474 194 : *pfa = gel(bnrconductor_factored(bnrz, Hz), 2);
475 194 : if (!Ldone) L = vec_append(L, gel(vP,1));
476 194 : return gc_all(av0, 2, &L, pfa);
477 : }
478 : }
479 : }
480 228 : P = gel(vP,1);
481 : }
482 30056 : else if (kum) P = prlifttoKz(nfz, nf, P, &kum->COMPO);
483 30284 : if (Ldone || ZM_equal(M, Hsofar)) continue;
484 25619 : L = vec_append(L, P);
485 25619 : Hsofar = M;
486 25619 : if (H_is_good(Hsofar, gell))
487 : {
488 16030 : if (!Hz) return gc_all(av0, pfa? 2: 1, &L, pfa);
489 72 : 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 39652 : subgroup_info(GEN bnfz, GEN Lprz, GEN gell, GEN vecWA)
501 : {
502 39652 : GEN M, nfz = bnf_get_nf(bnfz), Lell = mkvec(gell);
503 39652 : long i, j, ell = gell[2], l = lg(vecWA), lz = lg(Lprz);
504 39652 : M = cgetg(l, t_MAT);
505 166827 : for (j=1; j<l; j++) gel(M,j) = cgetg(lz, t_VECSMALL);
506 65295 : for (i=1; i < lz; i++)
507 : {
508 25643 : GEN pr = gel(Lprz,i), EX = subiu(pr_norm(pr), 1);
509 25643 : GEN N, g,T,p, prM = idealhnf_shallow(nfz, pr);
510 25643 : GEN modpr = zk_to_Fq_init(nfz, &pr,&T,&p);
511 25643 : long v = Z_lvalrem(divis(EX,ell), ell, &N) + 1; /* Norm(pr)-1 = N * ell^v */
512 25643 : GEN ellv = powuu(ell, v);
513 25643 : g = gener_Fq_local(T,p, Lell);
514 25643 : g = Fq_pow(g,N, T,p); /* order ell^v */
515 129549 : for (j=1; j < l; j++)
516 : {
517 103906 : GEN logc, c = gel(vecWA,j);
518 103906 : if (typ(c) == t_MAT) /* famat */
519 65990 : c = famat_makecoprime(nfz, gel(c,1), gel(c,2), pr, prM, EX);
520 103906 : c = nf_to_Fq(nfz, c, modpr);
521 103906 : c = Fq_pow(c, N, T,p);
522 103906 : logc = Fq_log(c, g, ellv, T,p);
523 103906 : ucoeff(M, i,j) = umodiu(logc, ell);
524 : }
525 : }
526 39652 : return M;
527 : }
528 :
529 : static GEN
530 39544 : futu(GEN bnf)
531 : {
532 39544 : GEN fu, tu, SUnits = bnf_get_sunits(bnf);
533 39544 : if (SUnits)
534 : {
535 31091 : tu = nf_to_scalar_or_basis(bnf_get_nf(bnf), bnf_get_tuU(bnf));
536 31091 : fu = bnf_compactfu(bnf);
537 : }
538 : else
539 : {
540 8453 : GEN U = bnf_build_units(bnf);
541 8453 : tu = gel(U,1); fu = vecslice(U, 2, lg(U)-1);
542 : }
543 39544 : return vec_append(fu, tu);
544 : }
545 : static GEN
546 39544 : bnf_cycgenmod(GEN bnf, long ell, GEN *pselmer, long *prc)
547 : {
548 39544 : GEN gen = bnf_get_gen(bnf), cyc = bnf_get_cyc(bnf), nf = bnf_get_nf(bnf);
549 39544 : GEN B, r = ZV_to_Flv(cyc, ell);
550 39544 : long i, rc, l = lg(gen);
551 39544 : B = cgetg(l, t_VEC);
552 41836 : for (i = 1; i < l && !r[i]; i++);
553 39544 : *prc = rc = i-1; /* ell-rank */
554 42583 : for (i = 1; i < l; i++)
555 : {
556 3039 : GEN G, q, c = gel(cyc,i), g = gel(gen,i);
557 3039 : if (i > rc && r[i] != 1) c = muliu(c, Fl_inv(r[i], ell));
558 3039 : q = divis(c, ell); /* remainder = 0 (i <= rc) or 1 */
559 : /* compute (b) = g^c mod ell-th powers */
560 3039 : G = equali1(q)? g: idealpowred(nf, g, q); /* lose principal part */
561 3039 : G = idealpows(nf, G, ell);
562 3039 : if (i > rc) G = idealmul(nf, G, g);
563 3039 : gel(B,i) = gel(bnfisprincipal0(bnf, G, nf_GENMAT|nf_FORCE), 2);
564 : }
565 39544 : *pselmer = shallowconcat(futu(bnf), vecslice(B,1,rc));
566 39544 : return B;
567 : }
568 :
569 : static GEN
570 37113 : rnfkummersimple(GEN bnr, GEN H, long ell)
571 : {
572 : long j, lSp, rc;
573 : GEN bnf, nf,bid, cycgenmod, Sp, vsprk, matP;
574 37113 : GEN be, M, K, vecW, vecWB, vecBp, gell = utoipos(ell);
575 : /* primes landing in H must be totally split */
576 37113 : GEN Lpr = get_prlist(bnr, H, gell, NULL, NULL);
577 :
578 37113 : bnf = bnr_get_bnf(bnr); if (!bnf_get_sunits(bnf)) bnf_build_units(bnf);
579 37113 : nf = bnf_get_nf(bnf);
580 37113 : bid = bnr_get_bid(bnr);
581 37113 : list_Hecke(&Sp, &vsprk, nf, bid_get_fact2(bid), gell, NULL);
582 :
583 37113 : cycgenmod = bnf_cycgenmod(bnf, ell, &vecW, &rc);
584 37113 : lSp = lg(Sp);
585 37113 : vecBp = cgetg(lSp, t_VEC);
586 37113 : matP = cgetg(lSp, t_MAT);
587 96294 : for (j = 1; j < lSp; j++)
588 59181 : isprincipalell(bnf,gel(Sp,j), cycgenmod,ell,rc, &gel(matP,j),&gel(vecBp,j));
589 37113 : vecWB = shallowconcat(vecW, vecBp);
590 :
591 37113 : M = matlogall(nf, vecWB, 0, 0, gell, vsprk);
592 37113 : M = vconcat(M, shallowconcat(zero_Flm(rc,lg(vecW)-1), matP));
593 37113 : M = vconcat(M, subgroup_info(bnf, Lpr, gell, vecWB));
594 37113 : K = Flm_ker(M, ell);
595 37113 : if (ell == 2)
596 : {
597 35007 : GEN msign = nfsign(nf, vecWB), y;
598 35007 : GEN arch = ZV_to_zv(bid_get_arch(bid)); /* the conductor */
599 35007 : msign = Flm_mul(msign, K, 2);
600 35007 : y = Flm_ker(msign, 2);
601 35007 : y = zv_equal0(arch)? gel(y,1): Flm_Flc_invimage(msign, arch, 2);
602 35007 : K = Flm_Flc_mul(K, y, 2);
603 : }
604 : else
605 2106 : K = gel(K,1);
606 37113 : be = compute_beta(K, vecWB, gell, bnf);
607 37113 : be = nf_to_scalar_or_alg(nf, be);
608 37113 : if (typ(be) == t_POL) be = mkpolmod(be, nf_get_pol(nf));
609 37113 : return gsub(pol_xn(ell, 0), be);
610 : }
611 :
612 : static ulong
613 97542 : nf_to_logFl(GEN nf, GEN x, GEN modpr, ulong g, ulong q, ulong ell, ulong p)
614 : {
615 97542 : x = nf_to_Fp_coprime(nf, x, modpr);
616 97542 : return Fl_log(Fl_powu(umodiu(x, p), q, p), g, ell, p);
617 : }
618 : static GEN
619 30020 : nfV_to_logFlv(GEN nf, GEN x, GEN modpr, ulong g, ulong q, ulong ell, ulong p)
620 127562 : { 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 2431 : kervirtualunit(struct rnfkummer *kum, GEN vselmer)
630 : {
631 2431 : GEN bnf = kum->bnfz, cyc = bnf_get_cyc(bnf), nf = bnf_get_nf(bnf);
632 2431 : GEN W, B, vy, vz, M, U1, U2, vtau, vell, SUnits = bnf_get_sunits(bnf);
633 2431 : long i, j, r, l = lg(vselmer), rc = kum->rc, ru = l-1 - rc, ell = kum->ell;
634 2431 : long LIMC = SUnits? itou(gel(SUnits,4)): 1;
635 : ulong p;
636 : forprime_t T;
637 :
638 2431 : vtau = cgetg(l, t_VEC);
639 2431 : vell = cgetg(l, t_VEC);
640 11758 : for (j = 1; j < l; j++)
641 : {
642 9327 : GEN t = gel(vselmer,j);
643 9327 : if (typ(t) == t_MAT)
644 : {
645 : GEN ct;
646 6896 : t = nffactorback(bnf, gel(t,1), ZV_to_Flv(gel(t,2), ell));
647 6896 : t = Q_primitive_part(t, &ct);
648 6896 : if (ct)
649 : {
650 3185 : GEN F = Q_factor(ct);
651 3185 : ct = factorback2(gel(F,1), ZV_to_Flv(gel(F,2), ell));
652 3185 : t = (typ(t) == t_INT)? ct: ZC_Z_mul(t, ct);
653 : }
654 : }
655 9327 : gel(vell,j) = t; /* integral, not too far from primitive */
656 9327 : gel(vtau,j) = Rg_tau(t, &kum->tau);
657 : }
658 2431 : U1 = vecslice(vell, 1, ru); /* units */
659 2431 : U2 = vecslice(vell, ru+1, ru+rc); /* cycgen (mod ell-th powers) */
660 2431 : B = nf_get_index(nf); /* bad primes; from 1 to ru are LIMC-units */
661 3316 : for (i = 1; i <= rc; i++) B = mulii(B, nfnorm(nf, gel(U2,i)));
662 2431 : if (LIMC > 1)
663 : {
664 2431 : GEN U, fa = absZ_factor_limit_strict(B, LIMC, &U), P = gel(fa,1);
665 2431 : long lP = lg(P);
666 2431 : B = U? gel(U,1): gen_1;
667 2431 : if (lP > 1 && cmpiu(gel(P,lP-1), LIMC) >= 0) B = mulii(B, gel(P,lP-1));
668 : }
669 2431 : if (is_pm1(B)) B = NULL;
670 2431 : vy = cgetg(l, t_MAT);
671 10873 : for (j = 1; j <= ru; j++) gel(vy,j) = zero_Flv(rc); /* units */
672 3316 : for ( ; j < l; j++)
673 : {
674 885 : GEN y, w, u = gel(vtau, j); /* virtual unit */
675 885 : if (!idealispower(nf, u, ell, &w)) pari_err_BUG("kervirtualunit");
676 885 : y = isprincipal(bnf, w); setlg(y, rc+1);
677 885 : if (!ZV_equal0(y))
678 2294 : for (i = 1; i <= rc; i++)
679 1409 : gel(y,i) = diviiexact(mului(ell,gel(y,i)), gel(cyc,i));
680 885 : gel(vy,j) = ZV_to_Flv(y, ell);
681 : }
682 2431 : u_forprime_arith_init(&T, LIMC+1, ULONG_MAX, 1, ell);
683 2431 : M = cgetg(ru+1, t_MAT); r = 1; setlg(M,2);
684 2431 : vz = cgetg(ru+1, t_MAT);
685 10420 : while ((p = u_forprime_next(&T))) if (!B || umodiu(B,p))
686 : {
687 10372 : GEN P = idealprimedec_limit_f(nf, utoipos(p), 1);
688 10372 : long nP = lg(P)-1;
689 10372 : ulong g = rootsof1_Fl(ell, p), q = p / ell; /* (p-1) / ell */
690 21077 : for (i = 1; i <= nP; i++)
691 : {
692 13136 : GEN modpr = zkmodprinit(nf, gel(P,i));
693 : GEN z, v2;
694 13136 : gel(M, r) = nfV_to_logFlv(nf, U1, modpr, g, q, ell, p); /* log futu */
695 13136 : if (Flm_rank(M, ell) < r) continue; /* discard */
696 :
697 8442 : v2 = nfV_to_logFlv(nf, U2, modpr, g, q, ell, p); /* log alpha[1..rc] */
698 8442 : gel(vz, r) = z = nfV_to_logFlv(nf, vtau, modpr, g, q, ell, p);
699 11593 : for (j = ru+1; j < l; j++)
700 3151 : uel(z,j) = Fl_sub(uel(z,j), Flv_dotproduct(v2, gel(vy,j), ell), ell);
701 8442 : if (r == ru) break;
702 6011 : r++; setlg(M, r+1);
703 : }
704 10372 : if (i < nP) break;
705 : }
706 2431 : 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 2431 : M = Flm_inv(Flm_transpose(M), ell);
710 2431 : vz = Flm_transpose(vz); /* now ru x #vtau */
711 11758 : for (j = 1; j < l; j++)
712 9327 : gel(vy,j) = shallowconcat(Flm_Flc_mul(M, gel(vz,j), ell), gel(vy,j));
713 2431 : W = Flm_ker(Flm_Fl_sub(vy, kum->g, ell), ell); l = lg(W);
714 7842 : for (j = 1; j < l; j++)
715 5411 : gel(W,j) = famat_reduce(famatV_zv_factorback(vselmer, gel(W,j)));
716 2431 : 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 4881 : get_mmu(long b, GEN r, long ell)
723 : {
724 4881 : long i, m = lg(r)-1;
725 4881 : GEN M = cgetg(m, t_VECSMALL);
726 23828 : for (i = 0; i < m-1; i++) M[i+1] = (r[b + 1] * r[m - i]) / ell;
727 4881 : return M;
728 : }
729 : /* max_b zv_sum(mu_b) < m ell */
730 : static long
731 2155 : max_smu(GEN r, long ell)
732 : {
733 2155 : long i, s = 0, z = vecsmall_max(r), l = lg(r);
734 6374 : for (i = 2; i < l; i++) s += (z * r[i]) / ell;
735 2155 : return s;
736 : }
737 :
738 : /* coeffs(x, a..b) in variable 0 >= varn(x) */
739 : static GEN
740 15216 : split_pol(GEN x, long a, long b)
741 : {
742 15216 : long i, l = degpol(x);
743 15216 : GEN y = x + a, z;
744 :
745 15216 : if (l < b) b = l;
746 15216 : if (a > b || varn(x) != 0) return pol_0(0);
747 15216 : l = b-a + 3;
748 15216 : z = cgetg(l, t_POL); z[1] = x[1];
749 87614 : for (i = 2; i < l; i++) gel(z,i) = gel(y,i);
750 15216 : 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 7608 : mod_Xell_a(GEN z, long ell, GEN an, GEN ad, GEN T)
757 : {
758 7608 : GEN z1 = split_pol(z, ell, degpol(z));
759 7608 : GEN z0 = split_pol(z, 0, ell-1); /* z = v^ell z1 + z0*/
760 7608 : if (ad) z0 = ZXX_Z_mul(z0, ad);
761 7608 : 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 7420 : to_alg(GEN nfz, GEN c, GEN D)
766 : {
767 7420 : if (typ(c) != t_COL) return D? mulii(D,c): c;
768 7420 : return RgV_dotproduct(nf_get_zkprimpart(nfz), c);
769 : }
770 : /* assume x in alg form */
771 : static GEN
772 7608 : downtoK(toK_s *T, GEN x)
773 : {
774 7608 : if (typ(x) != t_POL) return x;
775 7608 : x = RgM_RgC_mul(T->invexpoteta1, RgX_to_RgC(x, lg(T->invexpoteta1) - 1));
776 7608 : 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 2539 : compute_polrel(struct rnfkummer *kum, GEN be)
782 : {
783 2539 : toK_s *T = &kum->T;
784 2539 : long i, k, MU = 0, ell = kum->ell, m = T->m;
785 2539 : 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 2539 : GEN C, prim_Rk, C_Rk, prim_root, C_root, mell = utoineg(ell);
788 2539 : GEN nfz = bnf_get_nf(kum->bnfz), Tz = nf_get_pol(nfz), Dz = nf_get_zkden(nfz);
789 : pari_timer ti;
790 :
791 2539 : if (DEBUGLEVEL>1) { err_printf("Computing Newton sums: "); timer_start(&ti); }
792 2539 : if (equali1(Dz)) Dz = NULL;
793 2539 : D = Dz;
794 2539 : Ninvbe = Q_remove_denom(nfinv(nfz, be), &Dinvbe);
795 2539 : powtau_Ninvbe = powtau(Ninvbe, m-1, T->tau);
796 2539 : if (Dinvbe)
797 : {
798 2155 : MU = max_smu(r, ell);
799 2155 : D = mul_denom(Dz, powiu(Dinvbe, MU));
800 : }
801 :
802 2539 : root = cgetg(ell + 2, t_POL); /* compute D*root, will correct at the end */
803 2539 : root[1] = evalsigne(1) | evalvarn(0);
804 2539 : gel(root,2) = gen_0;
805 2539 : gel(root,3) = D? D: gen_1;
806 7608 : for (i = 2; i < ell; i++) gel(root,2+i) = gen_0;
807 7420 : for (i = 1; i < m; i++)
808 : { /* compute (1/be) ^ (-mu) instead of be^mu [mu < 0].
809 : * 1/be = Ninvbe / Dinvbe */
810 4881 : GEN mmu = get_mmu(i, r, ell), t;
811 4881 : t = to_alg(nfz, nffactorback(nfz, powtau_Ninvbe, mmu), Dz);/* Ninvbe^-mu */
812 4881 : if (Dinvbe)
813 : {
814 4219 : long a = MU - zv_sum(mmu);
815 4219 : if (a) t = gmul(t, powiu(Dinvbe, a));
816 : }
817 4881 : gel(root, 2 + r[i+1]) = t; /* root += D * (z_ell*T)^{r_i} be^mu_i */
818 : }
819 2539 : 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 2539 : prim_root = Q_primitive_part(root, &C_root);
823 2539 : C_root = div_content(C_root, D);
824 :
825 : /* theta^ell = be^( sum tau^a r_{d-1-a} ) = a = numa / Dz */
826 2539 : numa = to_alg(nfz, nffactorback(nfz, powtau(be, m, T->tau),
827 : vecsmall_reverse(r)), Dz);
828 2539 : if (DEBUGLEVEL>1) err_printf("root(%ld) ", timer_delay(&ti));
829 :
830 : /* Compute mod (X^ell - t, nfz.pol) */
831 2539 : C_Rk = C_root; prim_Rk = prim_root;
832 2539 : C = div_content(C_root, Dz);
833 2539 : S = cgetg(ell+3, t_POL); /* Newton sums */
834 2539 : S[1] = evalsigne(1) | evalvarn(0);
835 2539 : gel(S,2) = gen_0;
836 10147 : for (k = 2; k <= ell; k++)
837 : { /* compute the k-th Newton sum; here C_Rk ~ C_root */
838 7608 : pari_sp av = avma;
839 7608 : GEN z, C_z, d, Rk = ZXQX_mul(prim_Rk, prim_root, Tz);
840 7608 : Rk = mod_Xell_a(Rk, ell, numa, Dz, Tz); /* (mod X^ell - a, nfz.pol) */
841 7608 : prim_Rk = Q_primitive_part(Rk, &d); /* d C_root ~ 1 */
842 7608 : C_Rk = mul_content(C_Rk, mul_content(d, C));
843 : /* root^k = prim_Rk * C_Rk */
844 7608 : z = Q_primitive_part(gel(prim_Rk,2), &C_z); /* C_z ~ 1/C_root ~ 1/C_Rk */
845 7608 : z = downtoK(T, z);
846 7608 : C_z = mul_content(mul_content(C_z, C_Rk), mell);
847 7608 : z = gmul(z, C_z); /* C_z ~ 1 */
848 7608 : (void)gc_all(av, C_Rk? 3: 2, &z, &prim_Rk, &C_Rk);
849 7608 : if (DEBUGLEVEL>1) err_printf("%ld(%ld) ", k, timer_delay(&ti));
850 7608 : gel(S,k+1) = z; /* - Newton sum */
851 : }
852 2539 : gel(S,ell+2) = gen_m1; if (DEBUGLEVEL>1) err_printf("\n");
853 2539 : return RgX_recip(RgXn_expint(S,ell+1));
854 : }
855 :
856 : static void
857 2431 : compositum_red(compo_s *C, GEN P, GEN Q)
858 : {
859 2431 : GEN p, q, a, z = gel(compositum2(P, Q),1);
860 2431 : a = gel(z,1);
861 2431 : p = gel(gel(z,2), 2);
862 2431 : q = gel(gel(z,3), 2);
863 2431 : C->k = itos( gel(z,4) );
864 2431 : z = polredbest(a, nf_ORIG);
865 2431 : C->R = gel(z,1);
866 2431 : a = gel(gel(z,2), 2);
867 2431 : C->p = RgX_RgXQ_eval(p, a, C->R);
868 2431 : C->q = RgX_RgXQ_eval(q, a, C->R);
869 2431 : C->rev = QXQ_reverse(a, C->R);
870 2431 : }
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 2539 : nfX_Z_normalize(GEN nf, GEN P)
876 : {
877 : long i, l;
878 2539 : GEN C, Cj, PZ = cgetg_copy(P, &l);
879 2539 : PZ[1] = P[1];
880 15225 : for (i = 2; i < l; i++) /* minor variation on RgX_to_nfX (create PZ) */
881 : {
882 12686 : GEN z = nf_to_scalar_or_basis(nf, gel(P,i));
883 12686 : if (typ(z) == t_INT)
884 9224 : gel(PZ,i) = gel(P,i) = z;
885 : else
886 3462 : gel(PZ,i) = ZV_content(z);
887 : }
888 2539 : (void)ZX_Z_normalize(PZ, &C);
889 :
890 2539 : if (C == gen_1) return;
891 394 : Cj = C;
892 1840 : for (i = l-2; i > 1; i--)
893 : {
894 1446 : if (i != l-2) Cj = mulii(Cj, C);
895 1446 : gel(P,i) = gdiv(gel(P,i), Cj);
896 : }
897 : }
898 :
899 : /* set kum->vecC, kum->tQ */
900 : static void
901 701 : _rnfkummer_step4(struct rnfkummer *kum, long d, long m)
902 : {
903 701 : long i, j, rc = kum->rc;
904 701 : GEN Q, vT, vB, vC, vz, B = cgetg(rc+1,t_VEC), T = cgetg(rc+1,t_MAT);
905 701 : GEN gen = bnf_get_gen(kum->bnfz), cycgenmod = kum->cycgenmod;
906 701 : ulong ell = kum->ell;
907 :
908 1586 : for (j = 1; j <= rc; j++)
909 : {
910 885 : GEN t = gel(gen,j);
911 885 : t = ZM_hnfmodid(RgM_mul(kum->tau.zk, t), gcoeff(t, 1,1)); /* tau(t) */
912 885 : isprincipalell(kum->bnfz, t, cycgenmod,ell,rc, &gel(T,j), &gel(B,j));
913 : }
914 701 : Q = Flm_ker(Flm_Fl_sub(Flm_transpose(T), kum->g, ell), ell);
915 701 : kum->tQ = lg(Q) == 1? NULL: Flm_transpose(Q);
916 701 : kum->vecC = vC = cgetg(rc+1, t_VEC);
917 : /* T = rc x rc matrix */
918 701 : vT = Flm_powers(T, m-2, ell);
919 701 : vB = cgetg(m, t_VEC);
920 701 : vz = cgetg(rc+1, t_VEC);
921 1586 : for (i = 1; i <= rc; i++) gel(vz, i) = cgetg(m, t_VEC);
922 1900 : for (j = 1; j < m; j++)
923 : {
924 1199 : GEN Tj = Flm_Fl_mul(gel(vT,m-j), Fl_mul(j,d,ell), ell);
925 1199 : gel(vB, j) = RgV_tau(j == 1? B: gel(vB, j-1), &kum->tau);
926 2688 : for (i = 1; i <= rc; i++) gmael(vz, i, j) = gel(Tj, i);
927 : }
928 701 : vB = shallowconcat1(vB);
929 1586 : for (i = 1; i <= rc; i++)
930 : {
931 885 : GEN z = shallowconcat1(gel(vz,i));
932 885 : gel(vC,i) = famat_reduce(famatV_zv_factorback(vB, z));
933 : }
934 701 : }
935 :
936 : /* alg 5.3.5 */
937 : static void
938 2431 : rnfkummer_init(struct rnfkummer *kum, GEN bnf, GEN P, ulong ell, long prec)
939 : {
940 2431 : compo_s *COMPO = &kum->COMPO;
941 2431 : toK_s *T = &kum->T;
942 2431 : 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 2431 : if (DEBUGLEVEL>2) err_printf("Step 1\n");
947 2431 : if (DEBUGLEVEL) timer_start(&ti);
948 2431 : compositum_red(COMPO, polnf, polcyclo(ell, varn(polnf)));
949 2431 : 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 2431 : if (DEBUGLEVEL>2) err_printf("Step 2\n");
955 2431 : degK = degpol(polnf);
956 2431 : degKz = degpol(COMPO->R);
957 2431 : m = degKz / degK; /* > 1 */
958 2431 : d = (ell-1) / m;
959 2431 : g = Fl_powu(pgener_Fl(ell), d, ell);
960 2431 : if (Fl_powu(g, m, ell*ell) == 1) g += ell;
961 : /* ord(g) = m in all (Z/ell^k)^* */
962 2431 : if (DEBUGLEVEL>2) err_printf("Step 3\n");
963 2431 : nfz = nfinit(mkvec2(COMPO->R, P), prec);
964 2431 : if (lg(nfcertify(nfz)) > 1) nfz = nfinit(COMPO->R, prec); /* paranoia */
965 2431 : kum->bnfz = bnfz = Buchall(nfz, nf_FORCE, prec);
966 2431 : if (DEBUGLEVEL) timer_printf(&ti, "[rnfkummer] bnfinit(Kz)");
967 2431 : kum->cycgenmod = bnf_cycgenmod(bnfz, ell, &vselmer, &kum->rc);
968 2431 : kum->ell = ell;
969 2431 : kum->g = g;
970 2431 : kum->mgi = Fl_div(m, g, ell);
971 2431 : get_tau(kum);
972 2431 : if (DEBUGLEVEL>2) err_printf("Step 4\n");
973 2431 : if (kum->rc)
974 701 : _rnfkummer_step4(kum, d, m);
975 : else
976 1730 : { kum->vecC = cgetg(1, t_VEC); kum->tQ = NULL; }
977 2431 : if (DEBUGLEVEL>2) err_printf("Step 5\n");
978 2431 : kum->vecW = kervirtualunit(kum, vselmer);
979 2431 : if (DEBUGLEVEL>2) err_printf("Step 8\n");
980 : /* left inverse */
981 2431 : T->invexpoteta1 = QM_inv(RgXQ_matrix_pow(COMPO->p, degKz, degK, COMPO->R));
982 2431 : T->polnf = polnf;
983 2431 : T->tau = &kum->tau;
984 2431 : T->m = m;
985 2431 : T->powg = Fl_powers(g, m, ell);
986 2431 : }
987 :
988 : static GEN
989 2539 : rnfkummer_ell(struct rnfkummer *kum, GEN bnr, GEN H)
990 : {
991 2539 : ulong ell = kum->ell;
992 2539 : GEN bnfz = kum->bnfz, nfz = bnf_get_nf(bnfz), gell = utoipos(ell);
993 2539 : GEN vecC = kum->vecC, vecW = kum->vecW, cycgenmod = kum->cycgenmod;
994 2539 : long lW = lg(vecW), rc = kum->rc, j, lSp;
995 2539 : 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 2539 : GEN Lpr = get_prlist(bnr, H, gell, &faFz, kum);
999 :
1000 2539 : if (DEBUGLEVEL>2) err_printf("Step 9, 10 and 11\n");
1001 2539 : list_Hecke(&Sp, &vsprk, nfz, faFz, gell, T->tau);
1002 :
1003 2539 : if (DEBUGLEVEL>2) err_printf("Step 12\n");
1004 2539 : lSp = lg(Sp);
1005 2539 : vecAp = cgetg(lSp, t_VEC);
1006 2539 : vecBp = cgetg(lSp, t_VEC);
1007 2539 : matP = cgetg(lSp, t_MAT);
1008 4696 : for (j = 1; j < lSp; j++)
1009 : {
1010 : GEN e, a;
1011 2157 : isprincipalell(bnfz, gel(Sp,j), cycgenmod,ell,rc, &e, &a);
1012 2157 : gel(matP,j) = e;
1013 2157 : gel(vecBp,j) = famat_mul_shallow(famatV_zv_factorback(vecC, zv_neg(e)), a);
1014 2157 : gel(vecAp,j) = Rg_lambda(gel(vecBp,j), T);
1015 : }
1016 2539 : if (DEBUGLEVEL>2) err_printf("Step 13\n");
1017 2539 : vecWA = shallowconcat(vecW, vecAp);
1018 2539 : vecWB = shallowconcat(vecW, vecBp);
1019 :
1020 2539 : if (DEBUGLEVEL>2) err_printf("Step 14, 15 and 17\n");
1021 2539 : M = matlogall(nfz, vecWA, lW, kum->mgi, gell, vsprk);
1022 2539 : if (kum->tQ)
1023 : {
1024 274 : GEN QtP = Flm_mul(kum->tQ, matP, ell);
1025 274 : M = vconcat(M, shallowconcat(zero_Flm(lgcols(kum->tQ)-1,lW-1), QtP));
1026 : }
1027 2539 : lambdaWB = shallowconcat(RgV_lambda(vecW, T), vecAp);/*vecWB^lambda*/
1028 2539 : M = vconcat(M, subgroup_info(bnfz, Lpr, gell, lambdaWB));
1029 2539 : if (DEBUGLEVEL>2) err_printf("Step 16\n");
1030 2539 : K = Flm_ker(M, ell);
1031 2539 : if (DEBUGLEVEL>2) err_printf("Step 18\n");
1032 2539 : be = compute_beta(gel(K,1), vecWB, gell, kum->bnfz);
1033 2539 : P = compute_polrel(kum, be);
1034 2539 : nfX_Z_normalize(bnr_get_nf(bnr), P);
1035 2539 : if (DEBUGLEVEL>1) err_printf("polrel(beta) = %Ps\n", P);
1036 2539 : return P;
1037 : }
1038 :
1039 : static void
1040 34131 : bnr_check_var(GEN bnr)
1041 : {
1042 34131 : GEN T = nf_get_pol(bnr_get_nf(bnr));
1043 34131 : if (!varn(T)) pari_err_PRIORITY("bnrclassfield", T, "=", 0);
1044 34114 : }
1045 :
1046 : static void
1047 7608 : bnrclassfield_sanitize(GEN *pbnr, GEN *pH)
1048 : {
1049 7608 : bnr_subgroup_sanitize(pbnr, pH);
1050 7578 : bnr_check_var(*pbnr);
1051 7567 : }
1052 :
1053 : static GEN
1054 1002 : _rnfkummer(GEN bnr, GEN H, long prec)
1055 : {
1056 : ulong ell;
1057 : GEN gell, bnf, nf, P;
1058 : struct rnfkummer kum;
1059 :
1060 1002 : bnrclassfield_sanitize(&bnr, &H);
1061 997 : gell = H? ZM_det(H): ZV_prod(bnr_get_cyc(bnr));
1062 997 : ell = itou(gell);
1063 997 : if (ell == 1) return pol_x(0);
1064 997 : if (!uisprime(ell)) pari_err_IMPL("rnfkummer for composite relative degree");
1065 997 : if (bnf_get_tuN(bnr_get_bnf(bnr)) % ell == 0)
1066 639 : return rnfkummersimple(bnr, H, ell);
1067 358 : bnf = bnr_get_bnf(bnr); nf = bnf_get_nf(bnf);
1068 358 : P = ZV_union_shallow(nf_get_ramified_primes(nf), mkvec(gell));
1069 358 : rnfkummer_init(&kum, bnf, P, ell, maxss(prec,BIGDEFAULTPREC));
1070 358 : return rnfkummer_ell(&kum, bnr, H);
1071 : }
1072 :
1073 : GEN
1074 582 : rnfkummer(GEN bnr, GEN H, long prec)
1075 582 : { pari_sp av = avma; return gc_GEN(av, _rnfkummer(bnr, H, prec)); }
1076 :
1077 : /*******************************************************************/
1078 : /* bnrclassfield */
1079 : /*******************************************************************/
1080 :
1081 : /* TODO: could be exported */
1082 : static void
1083 248640 : gsetvarn(GEN x, long v)
1084 : {
1085 : long i;
1086 248640 : switch(typ(x))
1087 : {
1088 2538 : case t_POL: case t_SER:
1089 2538 : setvarn(x, v); return;
1090 0 : case t_LIST:
1091 0 : x = list_data(x); if (!x) return;
1092 : /* fall through t_VEC */
1093 : case t_VEC: case t_COL: case t_MAT:
1094 284874 : for (i = lg(x)-1; i > 0; i--) gsetvarn(gel(x,i), v);
1095 : }
1096 : }
1097 :
1098 : /* emb root of pol as polmod modulo pol2, return relative polynomial */
1099 : static GEN
1100 402 : relative_pol(GEN pol, GEN emb, GEN pol2)
1101 : {
1102 : GEN eqn, polrel;
1103 402 : if (degree(pol)==1) return pol2;
1104 360 : eqn = gsub(liftpol_shallow(emb), pol_x(varn(pol)));
1105 360 : eqn = Q_remove_denom(eqn, NULL);
1106 360 : polrel = nfgcd(pol2, eqn, pol, NULL);
1107 360 : return RgX_Rg_div(polrel, leading_coeff(polrel));
1108 : }
1109 :
1110 : /* pol defines K/nf */
1111 : static GEN
1112 420 : bnrclassfield_tower(GEN bnr, GEN subgroup, GEN TB, GEN p, long finaldeg, long absolute, long prec)
1113 : {
1114 420 : pari_sp av = avma;
1115 420 : GEN nf = bnr_get_nf(bnr), nf2, rnf, bnf2, bnr2;
1116 : GEN H, cyc, pk, sgpk, pol2, emb, emb2, famod, fa, Lbad;
1117 : long i, r1, ell, sp, spk, last;
1118 : forprime_t iter;
1119 :
1120 420 : rnf = rnfinit0(nf, TB, 1);
1121 420 : nf2 = rnf_build_nfabs(rnf, prec);
1122 420 : if (lg(nfcertify(nf2)) > 1)
1123 : {
1124 0 : rnf = rnfinit0(nf, gel(TB,1), 1);
1125 0 : nf2 = rnf_build_nfabs(rnf, prec);
1126 : }
1127 420 : gsetvarn(nf2, varn(nf_get_pol(nf)));
1128 :
1129 420 : r1 = nf_get_r1(nf2);
1130 420 : bnf2 = Buchall(nf2, nf_FORCE, prec);
1131 :
1132 420 : sp = itos(p);
1133 420 : spk = sp * rnf_get_degree(rnf);
1134 420 : pk = stoi(spk);
1135 420 : sgpk = hnfmodid(subgroup,pk);
1136 420 : last = spk==finaldeg;
1137 :
1138 : /* compute conductor */
1139 420 : famod = gel(bid_get_fact2(bnr_get_bid(bnr)), 1);
1140 420 : if (lg(famod)==1)
1141 : {
1142 120 : fa = trivial_fact();
1143 120 : Lbad = cgetg(1, t_VECSMALL);
1144 : }
1145 : else
1146 : {
1147 300 : long j, l = lg(famod);
1148 : GEN P, E;
1149 300 : P = cgetg(l, t_COL);
1150 300 : Lbad = cgetg(l, t_VEC);
1151 1116 : for(i = j = 1; i < l; i++)
1152 : {
1153 816 : GEN q, pr = gel(famod,i);
1154 816 : gel(P,i) = rnfidealprimedec(rnf, pr);
1155 816 : q = pr_get_p(pr);
1156 816 : if (lgefint(q) == 3) gel(Lbad,j++) = q;
1157 : }
1158 300 : setlg(Lbad,j);
1159 300 : Lbad = ZV_to_zv(ZV_sort_uniq_shallow(Lbad));
1160 300 : P = shallowconcat1(P);
1161 300 : P = gen_sort(P, (void*)&cmp_prime_ideal, &cmp_nodata);
1162 300 : settyp(P, t_COL); l = lg(P);
1163 300 : E = cgetg(l, t_COL);
1164 1200 : for (i = 1; i < l; i++)
1165 : {
1166 900 : GEN pr = gel(P,i);
1167 900 : long e = equalii(p, pr_get_p(pr))? 1 + (pr_get_e(pr)*sp) / (sp-1): 1;
1168 900 : gel(E,i) = utoipos(e);
1169 : }
1170 300 : fa = mkmat2(P, E);
1171 : }
1172 420 : bnr2 = Buchraymod(bnf2, mkvec2(fa, const_vec(r1,gen_1)), nf_INIT, pk);
1173 :
1174 : /* compute subgroup */
1175 420 : cyc = bnr_get_cyc(bnr2);
1176 420 : H = Flm_image(zv_diagonal(ZV_to_Flv(cyc,sp)), sp);
1177 420 : u_forprime_init(&iter, 2, ULONG_MAX);
1178 26766 : while ((ell = u_forprime_next(&iter))) if (!zv_search(Lbad, ell))
1179 : {
1180 26268 : GEN dec = idealprimedec_limit_f(nf, utoi(ell), 1);
1181 26268 : long l = lg(dec);
1182 50316 : for (i = 1; i < l; i++)
1183 : {
1184 24048 : GEN pr = gel(dec,i), Pr = gel(rnfidealprimedec(rnf, pr), 1);
1185 24048 : long f = pr_get_f(Pr) / pr_get_f(pr);
1186 24048 : GEN vpr = FpC_Fp_mul(bnrisprincipalmod(bnr,pr,pk,0), utoi(f), pk);
1187 24048 : if (gequal0(ZC_hnfrem(vpr,sgpk)))
1188 3372 : H = vec_append(H, ZV_to_Flv(bnrisprincipalmod(bnr2,Pr,p,0), sp));
1189 : }
1190 26268 : if (lg(H) > lg(cyc)+3)
1191 : {
1192 420 : H = Flm_image(H, sp);
1193 420 : if (lg(cyc)-lg(H) == 1) break;
1194 : }
1195 : }
1196 420 : H = hnfmodid(shallowconcat(zm_to_ZM(H), diagonal_shallow(cyc)), p);
1197 :
1198 : /* polynomial over nf2 */
1199 420 : pol2 = _rnfkummer(bnr2, H, prec);
1200 : /* absolute polynomial */
1201 420 : pol2 = rnfequation2(nf2, pol2);
1202 420 : emb2 = gel(pol2,2); /* generator of nf2 as polmod modulo pol2 */
1203 420 : pol2 = gel(pol2,1);
1204 : /* polynomial over nf */
1205 420 : if (!absolute || !last)
1206 : {
1207 402 : emb = rnf_get_alpha(rnf); /* generator of nf as polynomial in nf2 */
1208 402 : emb = poleval(emb, emb2); /* generator of nf as polmod modulo pol2 */
1209 402 : pol2 = relative_pol(nf_get_pol(nf), emb, pol2);
1210 : }
1211 420 : if (!last) pol2 = rnfpolredbest(nf, pol2, 0);
1212 :
1213 420 : obj_free(rnf);
1214 420 : pol2 = gc_GEN(av, pol2);
1215 420 : if (last) return pol2;
1216 48 : TB = mkvec2(pol2, gel(TB,2));
1217 48 : return bnrclassfield_tower(bnr, subgroup, TB, p, finaldeg, absolute, prec);
1218 : }
1219 :
1220 : /* subgroups H_i of bnr s.t. bnr/H_i is cyclic and inter_i H_i = subgroup */
1221 : static GEN
1222 34401 : cyclic_compos(GEN subgroup)
1223 : {
1224 34401 : pari_sp av = avma;
1225 34401 : GEN Ui, L, pe, D = ZM_snf_group(subgroup, NULL, &Ui);
1226 34401 : long i, l = lg(D);
1227 :
1228 34401 : L = cgetg(l, t_VEC);
1229 34401 : if (l == 1) return L;
1230 34401 : pe = gel(D,1);
1231 73056 : for (i = 1; i < l; i++)
1232 38655 : gel(L,i) = hnfmodid(shallowconcat(subgroup, vecsplice(Ui,i)),pe);
1233 34401 : return gc_GEN(av, L);
1234 : }
1235 :
1236 : /* p prime; set pkum=NULL if p-th root of unity in base field
1237 : * absolute=1 allowed if extension is cyclic with exponent>1 */
1238 : static GEN
1239 34401 : bnrclassfield_primepower(struct rnfkummer *pkum, GEN bnr, GEN subgroup, GEN p,
1240 : GEN P, long absolute, long prec)
1241 : {
1242 34401 : GEN res, subs = cyclic_compos(subgroup);
1243 34401 : long i, l = lg(subs);
1244 :
1245 34401 : res = cgetg(l,t_VEC);
1246 73056 : for (i = 1; i < l; i++)
1247 : {
1248 38655 : GEN H = gel(subs,i), Hp = hnfmodid(H,p);
1249 38655 : GEN pol, pe, bnr2 = bnrtoprimitive(bnr, Hp, p);
1250 38655 : if (!bnr2) bnr2 = bnr;
1251 : else
1252 : {
1253 7494 : GEN map = bnrsurjection(bnr,bnr2);
1254 7494 : Hp = abmap_subgroup_image(map, Hp);
1255 : }
1256 38655 : if (pkum) pol = rnfkummer_ell(pkum, bnr2, Hp);
1257 36474 : else pol = rnfkummersimple(bnr2, Hp, itos(p));
1258 38655 : pe = ZM_det_triangular(H);
1259 38655 : if (!equalii(p,pe))
1260 372 : pol = bnrclassfield_tower(bnr, H, mkvec2(pol,P), p, itos(pe), absolute, prec);
1261 38655 : gel(res,i) = pol;
1262 : }
1263 34401 : return res;
1264 : }
1265 :
1266 : /* partition of v into two subsets whose products are as balanced as possible */
1267 : /* assume v sorted */
1268 : static GEN
1269 114 : vecsmall_balance(GEN v)
1270 : {
1271 : forvec_t T;
1272 114 : GEN xbounds, x, vuniq, mult, ind, prod, prodbest = gen_0, bound,
1273 114 : xbest = NULL, res1, res2;
1274 114 : long i=1, j, k1, k2;
1275 114 : if (lg(v) == 3) return mkvec2(mkvecsmall(1), mkvecsmall(2));
1276 36 : vuniq = cgetg(lg(v), t_VECSMALL);
1277 36 : mult = cgetg(lg(v), t_VECSMALL);
1278 36 : ind = cgetg(lg(v), t_VECSMALL);
1279 36 : vuniq[1] = v[1];
1280 36 : mult[1] = 1;
1281 36 : ind[1] = 1;
1282 138 : for (j=2; j<lg(v); j++)
1283 : {
1284 102 : if (v[j] == vuniq[i]) mult[i]++;
1285 : else
1286 : {
1287 12 : i++;
1288 12 : vuniq[i] = v[j];
1289 12 : mult[i] = 1;
1290 12 : ind[i] = j;
1291 : }
1292 : }
1293 36 : setlg(vuniq, ++i);
1294 36 : setlg(mult, i);
1295 36 : setlg(ind, i);
1296 :
1297 36 : vuniq = zv_to_ZV(vuniq);
1298 36 : prod = factorback2(vuniq, mult);
1299 36 : bound = sqrti(prod);
1300 36 : xbounds = cgetg(lg(mult), t_VEC);
1301 84 : for (i=1; i<lg(mult); i++) gel(xbounds,i) = mkvec2s(0,mult[i]);
1302 :
1303 36 : forvec_init(&T, xbounds, 0);
1304 234 : while ((x = forvec_next(&T)))
1305 : {
1306 198 : prod = factorback2(vuniq, x);
1307 198 : if (cmpii(prod,bound)<=0 && cmpii(prod,prodbest)>0)
1308 : {
1309 90 : prodbest = prod;
1310 90 : xbest = gcopy(x);
1311 : }
1312 : }
1313 36 : res1 = cgetg(lg(v), t_VECSMALL);
1314 36 : res2 = cgetg(lg(v), t_VECSMALL);
1315 84 : for (i=1,k1=1,k2=1; i<lg(xbest); i++)
1316 : {
1317 102 : for (j=0; j<itos(gel(xbest,i)); j++) res1[k1++] = ind[i]+j;
1318 132 : for (; j<mult[i]; j++) res2[k2++] = ind[i]+j;
1319 : }
1320 36 : setlg(res1, k1);
1321 36 : setlg(res2, k2); return mkvec2(res1, res2);
1322 : }
1323 :
1324 : /* TODO nfcompositum should accept vectors of pols */
1325 : /* assume all fields are linearly disjoint */
1326 : /* assume the polynomials are sorted by degree */
1327 : static GEN
1328 384 : nfcompositumall(GEN nf, GEN L)
1329 : {
1330 : GEN pol, vdeg, part;
1331 : long i;
1332 384 : if (lg(L)==2) return gel(L,1);
1333 114 : vdeg = cgetg(lg(L), t_VECSMALL);
1334 408 : for (i=1; i<lg(L); i++) vdeg[i] = degree(gel(L,i));
1335 114 : part = vecsmall_balance(vdeg);
1336 114 : pol = cgetg(3, t_VEC);
1337 342 : for (i = 1; i < 3; i++)
1338 : {
1339 228 : GEN L2 = vecpermute(L, gel(part,i)), T = nfcompositumall(nf, L2);
1340 228 : gel(pol,i) = rnfpolredbest(nf, T, 0);
1341 : }
1342 114 : return nfcompositum(nf, gel(pol,1), gel(pol,2), 2);
1343 : }
1344 :
1345 : static GEN
1346 33051 : disc_primes(GEN bnr)
1347 : {
1348 33051 : GEN v = bid_primes(bnr_get_bid(bnr));
1349 33051 : GEN r = nf_get_ramified_primes(bnr_get_nf(bnr));
1350 33051 : return ZV_union_shallow(r, v);
1351 : }
1352 : static struct rnfkummer **
1353 33033 : rnfkummer_initall(GEN bnr, GEN vP, GEN P, long prec)
1354 : {
1355 33033 : long i, w, l = lg(vP), S = vP[l-1] + 1;
1356 33033 : GEN bnf = bnr_get_bnf(bnr);
1357 : struct rnfkummer **vkum;
1358 :
1359 33033 : w = bnf_get_tuN(bnf);
1360 33033 : vkum = (struct rnfkummer **)stack_malloc(S * sizeof(struct rnfkummer*));
1361 33033 : if (prec < BIGDEFAULTPREC) prec = BIGDEFAULTPREC;
1362 66138 : for (i = 1; i < l; i++)
1363 : {
1364 33105 : long p = vP[i];
1365 33105 : if (w % p == 0) { vkum[p] = NULL; continue; }
1366 2073 : vkum[p] = (struct rnfkummer *)stack_malloc(sizeof(struct rnfkummer));
1367 2073 : rnfkummer_init(vkum[p], bnf, P, p, prec);
1368 : }
1369 33033 : return vkum;
1370 : }
1371 : /* fully handle a single congruence subgroup H */
1372 : static GEN
1373 34365 : bnrclassfield_H(struct rnfkummer **vkum, GEN bnr, GEN bad, GEN H0, GEN fa, long flag,
1374 : long prec)
1375 : {
1376 34365 : GEN PN = gel(fa,1), EN = gel(fa,2), res;
1377 34365 : long i, lPN = lg(PN), absolute;
1378 :
1379 34365 : if (lPN == 1) switch(flag)
1380 : {
1381 12 : case 0:
1382 12 : return mkvec(pol_x(0));
1383 12 : case 1:
1384 12 : return pol_x(0);
1385 12 : default: /* 2 */
1386 12 : res = shallowcopy(nf_get_pol(bnr_get_nf(bnr)));
1387 12 : setvarn(res,0); return res;
1388 : }
1389 34329 : absolute = flag==2 && lPN==2 && !equali1(gel(EN,1)); /* one prime, exponent > 1 */
1390 34329 : res = cgetg(lPN, t_VEC);
1391 68730 : for (i = 1; i < lPN; i++)
1392 : {
1393 34401 : GEN p = gel(PN,i), H = hnfmodid(H0, powii(p, gel(EN,i)));
1394 34401 : long sp = itos(p);
1395 34401 : if (absolute) absolute = FpM_rank(H,p)==lg(H)-2; /* cyclic */
1396 34401 : gel(res,i) = bnrclassfield_primepower(vkum[sp], bnr, H, p, bad, absolute, prec);
1397 : }
1398 34329 : res = liftpol_shallow(shallowconcat1(res));
1399 34329 : res = gen_sort_shallow(res, (void*)cmp_RgX, gen_cmp_RgX);
1400 34329 : if (flag)
1401 : {
1402 156 : GEN nf = bnr_get_nf(bnr);
1403 156 : res = nfcompositumall(nf, res);
1404 156 : if (flag==2 && !absolute) res = rnfequation(nf, res);
1405 : }
1406 34329 : return res;
1407 : }
1408 : /* for a vector of subgroups */
1409 : static GEN
1410 26541 : bnrclassfieldvec(GEN bnr, GEN v, long flag, long prec)
1411 : {
1412 26541 : long j, lv = lg(v);
1413 26541 : GEN vH, vfa, vP, P, w = cgetg(lv, t_VEC);
1414 26541 : struct rnfkummer **vkum = NULL;
1415 :
1416 26541 : if (lv == 1) return w;
1417 26535 : vH = cgetg(lv, t_VEC);
1418 26535 : vP = cgetg(lv, t_VEC);
1419 26535 : vfa = cgetg(lv, t_VEC);
1420 54378 : for (j = 1; j < lv; j++)
1421 : {
1422 27849 : GEN N, fa, H = bnr_subgroup_check(bnr, gel(v,j), &N);
1423 27849 : if (is_bigint(N)) pari_err_OVERFLOW("bnrclassfield [too large degree]");
1424 27843 : if (!H) H = diagonal_shallow(bnr_get_cyc(bnr));
1425 27843 : gel(vH,j) = H;
1426 27843 : gel(vfa,j) = fa = Z_factor(N);
1427 27843 : gel(vP,j) = ZV_to_zv(gel(fa, 1));
1428 : }
1429 26529 : vP = shallowconcat1(vP); vecsmall_sort(vP);
1430 26529 : vP = vecsmall_uniq_sorted(vP);
1431 26529 : P = disc_primes(bnr);
1432 26529 : if (lg(vP) > 1) vkum = rnfkummer_initall(bnr, vP, P, prec);
1433 54372 : for (j = 1; j < lv; j++)
1434 27843 : gel(w,j) = bnrclassfield_H(vkum, bnr, P, gel(vH,j), gel(vfa,j), flag, prec);
1435 26529 : return w;
1436 : }
1437 : /* flag:
1438 : * 0 t_VEC of polynomials whose compositum is the extension
1439 : * 1 single polynomial
1440 : * 2 single absolute polynomial */
1441 : GEN
1442 33165 : bnrclassfield(GEN bnr, GEN subgroup, long flag, long prec)
1443 : {
1444 33165 : pari_sp av = avma;
1445 : GEN N, fa, P;
1446 : struct rnfkummer **vkum;
1447 :
1448 33165 : if (flag<0 || flag>2) pari_err_FLAG("bnrclassfield [must be 0,1 or 2]");
1449 33153 : if (subgroup && typ(subgroup) == t_VEC)
1450 : {
1451 26553 : bnr_sanitize(&bnr);
1452 26553 : bnr_check_var(bnr);
1453 26547 : if (!char_check(bnr_get_cyc(bnr), subgroup))
1454 26541 : return gc_GEN(av, bnrclassfieldvec(bnr, subgroup, flag, prec));
1455 : }
1456 6606 : bnrclassfield_sanitize(&bnr, &subgroup);
1457 6570 : N = ZM_det_triangular(subgroup);
1458 6570 : if (equali1(N)) switch(flag)
1459 : {
1460 24 : case 0: set_avma(av); retmkvec(pol_x(0));
1461 6 : case 1: set_avma(av); return pol_x(0);
1462 6 : default: /* 2 */
1463 6 : P = shallowcopy(nf_get_pol(bnr_get_nf(bnr)));
1464 6 : setvarn(P,0); return gc_GEN(av,P);
1465 : }
1466 6534 : if (is_bigint(N)) pari_err_OVERFLOW("bnrclassfield [too large degree]");
1467 6522 : fa = Z_factor(N); P = disc_primes(bnr);
1468 6522 : vkum = rnfkummer_initall(bnr, ZV_to_zv(gel(fa,1)), P, prec);
1469 6522 : return gc_GEN(av, bnrclassfield_H(vkum, bnr, P, subgroup, fa, flag, prec));
1470 : }
|
