1 /* Copyright (C) 2016 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 #include "pari.h"
16 #include "paripriv.h"
17 /*******************************************************************/
18 /* LOGARITHMIC CLASS GROUP */
19 /*******************************************************************/
20 /* min(v, v(Log_p Norm_{F_\p/Q_p}(x))) */
21 static long
vlognorm(GEN nf,GEN T,GEN x,GEN p,long v)22 vlognorm(GEN nf, GEN T, GEN x, GEN p, long v)
23 {
24 GEN a = nf_to_scalar_or_alg(nf, x);
25 GEN N = RgXQ_norm(a, T);
26 if (typ(N) != t_PADIC) N = cvtop(N, p, v);
27 return minss(v, valp( Qp_log(N) ));
28 }
29 /* K number field, pr a maximal ideal, let K_pr be the attached local
30 * field, K_pr = Q_p[X] / (T), T irreducible. Return \tilde{e}(K_pr/Q_p) */
31 static long
etilde(GEN nf,GEN pr,GEN T)32 etilde(GEN nf, GEN pr, GEN T)
33 {
34 GEN gp = pr_get_p(pr);
35 ulong e = pr_get_e(pr);
36 long v, voo, vmin, p, k;
37
38 if (!u_pval(e, gp))
39 {
40 v = u_pval(pr_get_f(pr), gp);
41 return itou( mului(e, powiu(gp, v)) );
42 }
43 nf = checknf(nf);
44 p = itou(gp);
45 k = e / (p-1) + 1;
46 /* log Norm_{F_P/Q_p} (1 + P^k) = Tr(P^k) = p^[(k + v(Diff))/ e] Z_p */
47 voo = (k + idealval(nf, nf_get_diff(nf), pr)) / e;
48 vmin = vlognorm(nf, T, pr_get_gen(pr), gp, voo);
49 if (k > 1)
50 {
51 GEN U = idealprincipalunits(nf, pr, k);
52 GEN gen = abgrp_get_gen(U), cyc = abgrp_get_cyc(U);
53 long i, l = lg(cyc);
54 for (i = 1; i < l; i++)
55 {
56 if (voo - Z_lval(gel(cyc,i), p) >= vmin) break;
57 vmin = vlognorm(nf, T, gel(gen,i), gp, vmin);
58 }
59 }
60 v = u_lval(degpol(T), p) + (p == 2UL? 2 : 1) - vmin;
61 (void)u_lvalrem(e, p, &e);
62 return e * upowuu(p,v);
63 }
64 static long
ftilde_from_e(GEN pr,long e)65 ftilde_from_e(GEN pr, long e) { return pr_get_e(pr) * pr_get_f(pr) / e; }
66 static long
ftilde(GEN K,GEN pr,GEN T)67 ftilde(GEN K, GEN pr, GEN T) { return ftilde_from_e(pr, etilde(K,pr, T)); }
68
69 static long
get_ZpX_index(GEN K,GEN pr,GEN T)70 get_ZpX_index(GEN K, GEN pr, GEN T)
71 {
72 GEN p, pi;
73 long j, l = lg(T);
74 if (l == 2) return 1;
75 p = pr_get_p(pr); pi = nf_to_scalar_or_alg(K, pr_get_gen(pr));
76 for (j = 1; j < l; j++)
77 {
78 GEN t = gel(T,j);
79 if (t && gvaluation(RgXQ_norm(pi, t), p)) return j;
80 }
81 return 0;
82 }
83
84 /* Given a number field K and a prime p, return
85 * S = places of K above p [primedec]
86 * R = corresponding p-adic factors of K.pol (mod p^k), in the same order */
87 static GEN
padicfact(GEN K,GEN S,long k)88 padicfact(GEN K, GEN S, long k)
89 {
90 GEN R, p = pr_get_p(gel(S,1));
91 GEN T = gel(factorpadic(nf_get_pol(K), p, k), 1);
92 long l, i;
93 S = idealprimedec(K, p);
94 R = cgetg_copy(S, &l);
95 for (i = 1; i < l; i++)
96 {
97 long j = get_ZpX_index(K, gel(S,i), T);
98 gel(R,i) = gel(T,j);
99 gel(T,j) = NULL;
100 }
101 return R;
102 }
103
104 /* K a bnf, compute Cl'(K) = ell-Sylow of Cl(K) / (places above ell).
105 * Return [D, u, R0, U0, ordS]
106 * - D: cyclic factors for Cl'(K)
107 * - u: generators of cyclic factors (all coprime to ell)
108 * - R0: subgroup isprincipal(<S>) (divides K.cyc)
109 * - U0: generators of R0 are of the form S . U0
110 * - ordS[i] = order of S[i] in CL(K) */
111 static GEN
CL_prime(GEN K,GEN ell,GEN Sell)112 CL_prime(GEN K, GEN ell, GEN Sell)
113 {
114 GEN g, ordS, R0, U0, U, D, u, cyc = bnf_get_cyc(K);
115 long i, l, lD, lS = lg(Sell);
116
117 g = leafcopy(bnf_get_gen(K));
118 l = lg(g);
119 for (i = 1; i < l; i++)
120 {
121 GEN A = gel(g,i), a = gcoeff(A,1,1);
122 long v = Z_pvalrem(a, ell, &a);
123 if (v) gel(g,i) = hnfmodid(A, a); /* make coprime to ell */
124 }
125 R0 = cgetg(lS, t_MAT);
126 ordS = cgetg(lS, t_VEC);
127 for (i = 1; i < lS; i++)
128 {
129 gel(R0,i) = isprincipal(K, gel(Sell,i));
130 gel(ordS,i) = charorder(cyc, gel(R0,i)); /* order of Sell[i] */
131 }
132 R0 = shallowconcat(R0, diagonal_shallow(cyc));
133 /* R0 = subgroup generated by S in Cl(K) [ divides diagonal(K.cyc) ]*/
134 R0 = ZM_hnfall(R0, &U0, 2); /* [S | cyc] * U0 = R0 in HNF */
135 D = ZM_snfall(R0, &U,NULL);
136 D = RgM_diagonal_shallow(D);
137 lD = lg(D);
138 u = ZM_inv(U, NULL); settyp(u, t_VEC);
139 for (i = 1; i < lD; i++) gel(u,i) = idealfactorback(K,g,gel(u,i),1);
140 setlg(U0, l);
141 U0 = rowslice(U0,1,lS-1); /* restrict to 'S' part */
142 return mkvec5(D, u, R0, U0, ordS);
143 }
144
145 static GEN
ell1(GEN ell)146 ell1(GEN ell) { return equaliu(ell,2)? utoipos(5): addiu(ell,1); }
147
148 /* log N_{F_P/Q_p}(x) */
149 static GEN
vtilde_i(GEN K,GEN x,GEN T,GEN ell,long prec)150 vtilde_i(GEN K, GEN x, GEN T, GEN ell, long prec)
151 {
152 GEN N, cx;
153 if (typ(x) != t_POL) x = nf_to_scalar_or_alg(K, x);
154 if (typ(x) != t_POL) { cx = x; N = NULL; }
155 else
156 {
157 x = Q_primitive_part(x,&cx);
158 N = resultant(RgX_rem(x,T), T);
159 N = cvtop(N,ell,prec);
160 }
161 if (cx)
162 {
163 (void)Q_pvalrem(cx, ell, &cx);
164 if (!isint1(cx))
165 {
166 cx = gpowgs(cvtop(cx,ell,prec), degpol(T));
167 N = N? gmul(N, cx): cx;
168 }
169 }
170 return N? Qp_log(N): gen_0;
171 }
172 static GEN
vecvtilde_i(GEN K,GEN x,GEN T,GEN ell,long prec)173 vecvtilde_i(GEN K, GEN x, GEN T, GEN ell, long prec)
174 { pari_APPLY_same(vtilde_i(K, gel(x,i), T, ell, prec)); }
175 static GEN
vtilde(GEN K,GEN x,GEN T,GEN deg,GEN ell,long prec)176 vtilde(GEN K, GEN x, GEN T, GEN deg, GEN ell, long prec)
177 {
178 pari_sp av;
179 GEN v, G, E;
180 if (typ(x) != t_MAT) return gdiv(vtilde_i(K,x,T,ell,prec), deg);
181 G = gel(x,1);
182 E = gel(x,2); av = avma; v = vecvtilde_i(K,G,T,ell,prec);
183 return gerepileupto(av, gdiv(RgV_dotproduct(E, v), deg));
184 }
185
186 /* v[i] = deg S[i] mod p^prec */
187 static GEN
get_vdegS(GEN Ftilde,GEN ell,long prec)188 get_vdegS(GEN Ftilde, GEN ell, long prec)
189 {
190 long i, l = lg(Ftilde);
191 GEN v = cgetg(l, t_VEC), degell = Qp_log( cvtop(ell1(ell), ell, prec) );
192 for (i = 1; i < l; i++) gel(v,i) = gmulsg(Ftilde[i], degell);
193 return v;
194 }
195 /* K a bnf. Compute kernel \tilde{Cl}_K(ell); return cyclic factors.
196 * Set *pM to (vtilde_S[i](US[j]))_{i,j} */
197 static GEN
CL_tilde(GEN K,GEN US,GEN ell,GEN T,long imin,GEN vdegS,GEN * pM,long prec)198 CL_tilde(GEN K, GEN US, GEN ell, GEN T, long imin, GEN vdegS,
199 GEN *pM, long prec)
200 {
201 long i, j, k, lD, l = lg(T), lU = lg(US);
202 GEN D, M, ellk;
203
204 /* p = P^e: \tilde{Cl}(l) = (1) */
205 if (l == 2) { *pM = cgetg(1, t_MAT); return cgetg(1, t_VEC); }
206 M = cgetg(lU, t_MAT);
207 for (j = 1; j < lU; j++)
208 {
209 GEN c = cgetg(l, t_COL), a = gel(US,j);
210 for (i = 1; i < l; i++)
211 gel(c,i) = vtilde(K, a, gel(T,i), gel(vdegS,i), ell, prec);
212 gel(M,j) = c;
213 }
214 k = padicprec(M, ell); ellk = powiu(ell, k);
215 *pM = M = gmod(M, ellk);
216 M = ZM_hnfmodid(rowsplice(M, imin), ellk);
217 D = matsnf0(M, 4); lD = lg(D);
218 if (lD > 1 && Z_pval(gel(D,1), ell) >= k) return NULL;
219 return D;
220 }
221
222 /* [L:K] = ell^k; return 1 if L/K is locally cyclotomic at ell, 0 otherwise */
223 long
rnfislocalcyclo(GEN rnf)224 rnfislocalcyclo(GEN rnf)
225 {
226 pari_sp av = avma;
227 GEN K, L, S, SK, TK, SLs, SL2, TL, ell;
228 ulong ll;
229 long i, j, k, lk, lSK;
230 checkrnf(rnf);
231 lk = rnf_get_degree(rnf);
232 if (lk == 1) return 1;
233 k = uisprimepower(lk, &ll);
234 if (!k) pari_err_IMPL("rnfislocalcyclo for non-l-extensions");
235 ell = utoi(ll);
236 K = rnf_get_nf(rnf);
237 L = rnf_build_nfabs(rnf, nf_get_prec(K));
238 S = rnfidealprimedec(rnf, ell);
239 SK = gel(S,1);
240 SLs = gel(S,2);
241 SL2 = shallowconcat1(SLs);
242 TK = padicfact(K, SK, 100); lSK = lg(SK);
243 TL = padicfact(L, SL2, 100);
244 for (i = 1; i < lSK; i++)
245 {
246 long eK = etilde(K, gel(SK,i), gel(TK,i));
247 GEN SL = gel(SLs,i);
248 long lSL = lg(SL);
249 for (j = 1; j < lSL; j++)
250 {
251 long iS = gen_search(SL2, gel(SL,j), 0, (void*)&cmp_prime_over_p,
252 &cmp_nodata);
253 long eL = etilde(L, gel(SL,j), gel(TL,iS));
254 if (dvdui(eL/eK, ell)) return gc_long(av,0);
255 }
256 };
257 return gc_long(av,1);
258 }
259
260 #if 0
261 /* Return 1 if L/Q is locally cyclotomic at ell */
262 static int
263 islocalcycloQ(GEN L, GEN ell)
264 {
265 GEN SL = idealprimedec(L,ell), TL;
266 long i, lSL = lg(SL);
267 TL = padicfact(L, SL, 100);
268 for (i = 1; i < lSL; i++)
269 {
270 long eL = etilde(L, gel(SL,i), gel(TL,i));
271 if (dvdui(eL,ell)) return 0;
272 }
273 return 1;
274 }
275 #endif
276
277 /* true nf, pr a prid */
278 static long
nfislocalpower_i(GEN nf,GEN pr,GEN a,GEN n)279 nfislocalpower_i(GEN nf, GEN pr, GEN a, GEN n)
280 {
281 long v, e, t;
282 GEN p, G, L;
283 a = nf_to_scalar_or_basis(nf,a);
284 if (!signe(n)) return isint1(a);
285 v = nfvalrem(nf, a, pr, &a); if (!dvdsi(v, n)) return 0;
286 p = pr_get_p(pr);
287 v = Z_pvalrem(n, p, &n);
288 if (!equali1(n))
289 {
290 GEN T, modpr = zk_to_Fq_init(nf, &pr, &T, &p);
291 GEN ap = nf_to_Fq(nf, a, modpr);
292 if (!Fq_ispower(ap, n, T, p)) return 0;
293 }
294 if (!v) return 1;
295 e = pr_get_e(pr);
296 if (v == 1) /* optimal formula */
297 t = itos( divii(mului(e,p), subiu(p,1)) ) + 1;
298 else /* straight Hensel */
299 t = 2 * e * v + 1;
300 G = Idealstarprk(nf, pr, t, nf_INIT);
301 L = ideallogmod(nf, a, G, powiu(p, v));
302 return ZV_equal0(L) || ZV_pval(L, p) >= v;
303 }
304 long
nfislocalpower(GEN nf,GEN pr,GEN a,GEN n)305 nfislocalpower(GEN nf, GEN pr, GEN a, GEN n)
306 {
307 pari_sp av = avma;
308 if (typ(n) != t_INT) pari_err_TYPE("nfislocalpower",n);
309 nf = checknf(nf); checkprid(pr);
310 return gc_long(av, nfislocalpower_i(nf, pr, a, n));
311 }
312
313 /* v_ell( exponent(D) ) */
314 static long
ellexpo(GEN D,GEN ell)315 ellexpo(GEN D, GEN ell) { return lg(D) == 1? 0: Z_pval(gel(D,1), ell); }
316
317 static GEN
ellsylow(GEN cyc,GEN ell)318 ellsylow(GEN cyc, GEN ell)
319 {
320 long i, l;
321 GEN d = cgetg_copy(cyc, &l);
322 for (i = 1; i < l; i++)
323 {
324 GEN c = gel(cyc,i), a;
325 if (!Z_pvalrem(c, ell, &a)) break;
326 gel(d,i) = diviiexact(c, a);
327 }
328 setlg(d, i); return d;
329 }
330
331 static long
vnorm_x(GEN nf,GEN x,GEN ell)332 vnorm_x(GEN nf, GEN x, GEN ell)
333 {
334 x = nf_to_scalar_or_alg(nf,x);
335 if (typ(x) != t_POL) return 0;
336 x = Q_primpart(x);
337 return Q_pval(nfnorm(nf,x), ell);
338 }
339 static long
vtilde_prec_x(GEN nf,GEN x,GEN ell)340 vtilde_prec_x(GEN nf, GEN x, GEN ell)
341 {
342 long i, l, v;
343 GEN G;
344 if (typ(x) != t_MAT) return vnorm_x(nf,x,ell);
345 G = gel(x,1); l = lg(G); v = 0;
346 for (i = 1; i < l; i++) v = maxss(v, vnorm_x(nf,gel(G,i),ell));
347 return v;
348 }
349 /* upper bound for \delta(vec): estimate loss of accuracy when evaluating
350 * \tilde{v} on the vec[i] */
351 static long
vtilde_prec(GEN nf,GEN vec,GEN ell)352 vtilde_prec(GEN nf, GEN vec, GEN ell)
353 {
354 long v0 = 0, i, l = lg(vec);
355 for (i = 1; i < l; i++)
356 v0 = maxss(v0, vtilde_prec_x(nf, gel(vec,i), ell));
357 return 3 + v0 + z_pval(nf_get_degree(nf), ell);
358 }
359 static GEN
get_Ftilde(GEN nf,GEN S,GEN T,GEN ell,long * pimin)360 get_Ftilde(GEN nf, GEN S, GEN T, GEN ell, long *pimin)
361 {
362 long j, lS = lg(S), vmin = lS;
363 GEN Ftilde = cgetg(lS, t_VECSMALL);
364 *pimin = 1;
365 for (j = 1; j < lS; j++)
366 {
367 long f = ftilde(nf, gel(S,j), gel(T,j)), v = z_pval(f, ell);
368 Ftilde[j] = f; if (v < vmin) { vmin = v; *pimin = j; }
369 }
370 return Ftilde;
371 }
372 static GEN
bnflog_i(GEN bnf,GEN ell)373 bnflog_i(GEN bnf, GEN ell)
374 {
375 long prec0, prec;
376 GEN nf, US, vdegS, S, T, M, CLp, CLt, Ftilde, vtG, ellk;
377 GEN D, Ap, cycAp, fu;
378 long imin, i, j, lvAp;
379
380 bnf = checkbnf(bnf); nf = bnf_get_nf(bnf);
381 S = idealprimedec(nf, ell);
382 US = sunits_mod_units(bnf, S);
383 prec0 = maxss(30, vtilde_prec(nf, US, ell));
384 if (!(fu = bnf_build_cheapfu(bnf)) && !(fu = bnf_compactfu(bnf)))
385 bnf_build_units(bnf);
386 US = shallowconcat(fu, US);
387 settyp(US, t_COL);
388 T = padicfact(nf, S, prec0);
389 Ftilde = get_Ftilde(nf, S, T, ell, &imin);
390 CLp = CL_prime(bnf, ell, S);
391 cycAp = gel(CLp,1);
392 Ap = gel(CLp,2);
393 for(;;)
394 {
395 vdegS = get_vdegS(Ftilde, ell, prec0);
396 CLt = CL_tilde(nf, US, ell, T, imin, vdegS, &vtG, prec0);
397 if (CLt) break;
398 prec0 <<= 1;
399 T = padicfact(nf, S, prec0);
400 }
401 prec = ellexpo(cycAp, ell) + ellexpo(CLt,ell) + 1;
402 if (prec == 1) return mkvec3(cgetg(1,t_VEC), cgetg(1,t_VEC), cgetg(1,t_VEC));
403
404 ellk = powiu(ell, prec);
405 lvAp = lg(Ap);
406 if (lvAp > 1)
407 {
408 long lS = lg(S);
409 GEN Kcyc = bnf_get_cyc(bnf);
410 GEN C = zeromatcopy(lvAp-1, lS-1);
411 GEN Rell = gel(CLp,3), Uell = gel(CLp,4), ordS = gel(CLp,5);
412 for (i = 1; i < lvAp; i++)
413 {
414 GEN a, b, bi, A = gel(Ap,i), d = gel(cycAp,i);
415 bi = isprincipal(bnf, A);
416 a = vecmodii(ZC_Z_mul(bi,d), Kcyc);
417 /* a in subgroup generated by S = Rell; hence b integral */
418 b = hnf_invimage(Rell, a);
419 b = vecmodii(ZM_ZC_mul(Uell, ZC_neg(b)), ordS);
420 A = mkvec2(A, trivial_fact());
421 A = idealpowred(nf, A, d);
422 /* find a principal representative of A_i^cycA_i up to elements of S */
423 a = isprincipalfact(bnf,gel(A,1),S,b,nf_GENMAT|nf_FORCE);
424 if (!gequal0(gel(a,1))) pari_err_BUG("bnflog");
425 a = famat_mul_shallow(gel(A,2), gel(a,2)); /* principal part */
426 if (lg(a) == 1) continue;
427 for (j = 1; j < lS; j++)
428 gcoeff(C,i,j) = vtilde(nf, a, gel(T,j), gel(vdegS,j), ell, prec0);
429 }
430 C = gmod(gneg(C),ellk);
431 C = shallowtrans(C);
432 M = mkmat2(mkcol2(diagonal_shallow(cycAp), C), mkcol2(gen_0, vtG));
433 M = shallowmatconcat(M); /* relation matrix */
434 }
435 else
436 M = vtG;
437 M = ZM_hnfmodid(M, ellk);
438 D = matsnf0(M, 4);
439 if (lg(D) == 1 || !dvdii(gel(D,1), ellk))
440 pari_err_BUG("bnflog [missing Z_l component]");
441 D = vecslice(D,2,lg(D)-1);
442 return mkvec3(D, CLt, ellsylow(cycAp, ell));
443 }
444 GEN
bnflog(GEN bnf,GEN ell)445 bnflog(GEN bnf, GEN ell)
446 {
447 pari_sp av = avma;
448 return gerepilecopy(av, bnflog_i(bnf, ell));
449 }
450
451 GEN
bnflogef(GEN nf,GEN pr)452 bnflogef(GEN nf, GEN pr)
453 {
454 pari_sp av = avma;
455 long e, f, ef;
456 GEN p;
457 checkprid(pr); p = pr_get_p(pr);
458 nf = checknf(nf);
459 e = pr_get_e(pr);
460 f = pr_get_f(pr); ef = e*f;
461 if (u_pval(ef, p))
462 {
463 GEN T = gel(factorpadic(nf_get_pol(nf), p, 100), 1);
464 long j = get_ZpX_index(nf, pr, T);
465 e = etilde(nf, pr, gel(T,j));
466 f = ef / e;
467 }
468 set_avma(av); return mkvec2s(e,f);
469 }
470
471 GEN
bnflogdegree(GEN nf,GEN A,GEN ell)472 bnflogdegree(GEN nf, GEN A, GEN ell)
473 {
474 pari_sp av = avma;
475 GEN AZ, A0Z, NA0;
476 long vAZ;
477
478 if (typ(ell) != t_INT) pari_err_TYPE("bnflogdegree", ell);
479 nf = checknf(nf);
480 A = idealhnf(nf, A);
481 AZ = gcoeff(A,1,1);
482 vAZ = Z_pvalrem(AZ, ell, &A0Z);
483 if (is_pm1(A0Z))
484 NA0 = gen_1;
485 else
486 (void)Z_pvalrem(idealnorm(nf,A), ell, &NA0);
487 if (vAZ)
488 {
489 GEN Aell = ZM_hnfmodid(A, powiu(ell,vAZ));
490 GEN S = idealprimedec(nf, ell), T;
491 long l, i, s = 0;
492 T = padicfact(nf, S, 100);
493 l = lg(S);
494 for (i = 1; i < l; i++)
495 {
496 GEN P = gel(S,i);
497 long v = idealval(nf, Aell, P);
498 if (v) s += v * ftilde(nf, P, gel(T,i));
499 }
500 if (s) NA0 = gmul(NA0, gpowgs(ell1(ell), s));
501 }
502 return gerepileupto(av, NA0);
503 }
504