Line data Source code
1 : /* Copyright (C) 2020 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 : #include "pari.h"
15 : #include "paripriv.h"
16 :
17 : #define DEBUGLEVEL DEBUGLEVEL_bnf
18 :
19 : /* if x a famat, assume it is a true unit (very costly to check even that
20 : * it's an algebraic integer) */
21 : GEN
22 1203 : bnfisunit(GEN bnf, GEN x)
23 : {
24 1203 : long tx = typ(x), i, r1, RU, e, n, prec;
25 1203 : pari_sp av = avma;
26 : GEN t, logunit, ex, nf, pi2_sur_w, rx, emb, arg;
27 :
28 1203 : bnf = checkbnf(bnf); nf = bnf_get_nf(bnf);
29 1203 : RU = lg(bnf_get_logfu(bnf));
30 1203 : n = bnf_get_tuN(bnf); /* # { roots of 1 } */
31 1203 : if (tx == t_MAT)
32 : { /* famat, assumed OK */
33 954 : if (lg(x) != 3) pari_err_TYPE("bnfisunit [not a factorization]", x);
34 : } else {
35 249 : x = nf_to_scalar_or_basis(nf,x);
36 249 : if (typ(x) != t_COL)
37 : { /* rational unit ? */
38 : GEN v;
39 : long s;
40 102 : if (typ(x) != t_INT || !is_pm1(x)) return cgetg(1,t_COL);
41 96 : s = signe(x); set_avma(av); v = zerocol(RU);
42 96 : gel(v,RU) = utoi((s > 0)? 0: n>>1);
43 96 : return v;
44 : }
45 147 : if (!isint1(Q_denom(x))) retgc_const(av, cgetg(1, t_COL));
46 : }
47 :
48 1101 : r1 = nf_get_r1(nf);
49 1101 : prec = nf_get_prec(nf);
50 1101 : for (i = 1;; i++)
51 0 : {
52 : GEN rlog;
53 1101 : nf = bnf_get_nf(bnf);
54 1101 : logunit = bnf_get_logfu(bnf);
55 1101 : rlog = real_i(logunit);
56 1101 : rx = nflogembed(nf,x,&emb, prec);
57 1101 : if (rx)
58 : {
59 1101 : GEN logN = RgV_sum(rx); /* log(Nx), should be ~ 0 */
60 1101 : if (gexpo(logN) > -20)
61 : { /* precision problem ? */
62 12 : if (typ(logN) != t_REAL) retgc_const(av, cgetg(1, t_COL)); /*no*/
63 24 : if (i == 1 && typ(x) != t_MAT &&
64 24 : !is_pm1(nfnorm(nf, x))) retgc_const(av, cgetg(1, t_COL));
65 : }
66 : else
67 : {
68 360 : ex = RU == 1? cgetg(1,t_COL)
69 1089 : : RgM_solve(rlog, rx); /* ~ fundamental units exponents */
70 1089 : if (ex) { ex = grndtoi(ex, &e); if (e < -4) break; }
71 : }
72 : }
73 0 : if (i == 1)
74 0 : prec = nbits2prec(gexpo(x) + 128);
75 : else
76 : {
77 0 : if (i > 4) pari_err_PREC("bnfisunit");
78 0 : prec = precdbl(prec);
79 : }
80 0 : if (DEBUGLEVEL) pari_warn(warnprec,"bnfisunit",prec);
81 0 : bnf = bnfnewprec_shallow(bnf, prec);
82 : }
83 : /* choose a large embedding => small relative error */
84 1489 : for (i = 1; i < RU; i++)
85 1027 : if (signe(gel(rx,i)) > -1) break;
86 1089 : if (RU == 1) t = gen_0;
87 : else
88 : {
89 729 : t = imag_i( row_i(logunit,i, 1,RU-1) );
90 729 : t = RgV_dotproduct(t, ex);
91 729 : if (i > r1) t = gmul2n(t, -1);
92 : }
93 1089 : if (typ(emb) != t_MAT) arg = garg(gel(emb,i), prec);
94 : else
95 : {
96 954 : GEN p = gel(emb,1), e = gel(emb,2);
97 954 : long j, l = lg(p);
98 954 : arg = NULL;
99 52690 : for (j = 1; j < l; j++)
100 : {
101 51736 : GEN a = gmul(gel(e,j), garg(gel(gel(p,j),i), prec));
102 51736 : arg = arg? gadd(arg, a): a;
103 : }
104 : }
105 1089 : t = gsub(arg, t); /* = arg(the missing root of 1) */
106 1089 : pi2_sur_w = divru(mppi(prec), n>>1); /* 2pi / n */
107 1089 : e = umodiu(roundr(divrr(t, pi2_sur_w)), n);
108 1089 : if (n > 2)
109 : {
110 6 : GEN z = algtobasis(nf, bnf_get_tuU(bnf)); /* primitive root of 1 */
111 6 : GEN ro = RgV_dotproduct(row(nf_get_M(nf), i), z);
112 6 : GEN p2 = roundr(divrr(garg(ro, prec), pi2_sur_w));
113 6 : e *= Fl_inv(umodiu(p2,n), n);
114 6 : e %= n;
115 : }
116 1089 : gel(ex,RU) = utoi(e); setlg(ex, RU+1); return gc_GEN(av, ex);
117 : }
118 :
119 : /* split M a square ZM in HNF as [H, B; 0, Id], H in HNF without 1-eigenvalue */
120 : static GEN
121 850 : hnfsplit(GEN M, GEN *pB)
122 : {
123 850 : long i, l = lg(M);
124 2757 : for (i = l-1; i; i--)
125 2513 : if (!equali1(gcoeff(M,i,i))) break;
126 850 : if (!i) { *pB = zeromat(0, l-1); return cgetg(1, t_MAT); }
127 606 : *pB = matslice(M, 1, i, i+1, l-1); return matslice(M, 1, i, 1, i);
128 : }
129 :
130 : /* S a list of prime ideal in idealprimedec format. If pH != NULL, set it to
131 : * the HNF of the S-class group and return bnfsunit, else return bnfunits */
132 : static GEN
133 1228 : bnfsunit_i(GEN bnf, GEN S, GEN *pH, GEN *pA, GEN *pden)
134 : {
135 1228 : long FLAG, i, lS = lg(S);
136 : GEN M, U1, U2, U, V, H, Sunit, B, g;
137 :
138 1228 : if (!is_vec_t(typ(S))) pari_err_TYPE("bnfsunit",S);
139 1228 : bnf = checkbnf(bnf);
140 1228 : if (lS == 1)
141 : {
142 378 : *pA = cgetg(1,t_MAT);
143 378 : *pden = gen_1; return cgetg(1,t_VEC);
144 : }
145 850 : M = cgetg(lS,t_MAT); /* relation matrix for the S class group */
146 850 : g = cgetg(lS,t_MAT); /* principal part */
147 850 : FLAG = pH ? 0: nf_GENMAT;
148 4963 : for (i = 1; i < lS; i++)
149 : {
150 4113 : GEN pr = gel(S,i);
151 4113 : checkprid(pr);
152 4113 : if (pH)
153 1335 : gel(M,i) = isprincipal(bnf, pr);
154 : else
155 : {
156 2778 : GEN v = bnfisprincipal0(bnf, pr, FLAG);
157 2778 : gel(M,i) = gel(v,1);
158 2778 : gel(g,i) = gel(v,2);
159 : }
160 : }
161 : /* S class group and S units, use ZM_hnflll to get small 'U' */
162 850 : M = shallowconcat(M, diagonal_shallow(bnf_get_cyc(bnf)));
163 850 : H = ZM_hnflll(M, &U, 1); setlg(U, lS); if (pH) *pH = H;
164 850 : U1 = rowslice(U,1, lS-1);
165 850 : U2 = rowslice(U,lS, lg(M)-1); /* (M | cyc) [U1; U2] = 0 */
166 850 : H = ZM_hnflll(U1, pH? NULL: &V, 0);
167 : /* U1 = upper left corner of U, invertible. S * U1 = principal ideals
168 : * whose generators generate the S-units */
169 850 : H = hnfsplit(H, &B);
170 : /* [ H B ] [ H^-1 - H^-1 B ]
171 : * U1 * V = HNF(U1) = [ 0 Id ], inverse = [ 0 Id ]
172 : * S * HNF(U1) = integral generators for S-units = Sunit */
173 850 : Sunit = cgetg(lS, t_VEC);
174 850 : if (pH)
175 : {
176 124 : long nH = lg(H) - 1;
177 124 : FLAG = nf_FORCE | nf_GEN;
178 1459 : for (i = 1; i < lS; i++)
179 : {
180 1335 : GEN C = NULL, E;
181 1335 : if (i <= nH) E = gel(H,i); else { C = gel(S,i), E = gel(B,i-nH); }
182 1335 : gel(Sunit,i) = gel(isprincipalfact(bnf, C, S, E, FLAG), 2);
183 : }
184 : }
185 : else
186 : {
187 726 : GEN cycgen = bnf_build_cycgen(bnf);
188 726 : U1 = ZM_mul(U1, V);
189 726 : U2 = ZM_mul(U2, V);
190 726 : FLAG = nf_FORCE | nf_GENMAT;
191 3504 : for (i = 1; i < lS; i++)
192 : {
193 2778 : GEN a = famatV_factorback(g, gel(U1,i));
194 2778 : GEN b = famatV_factorback(cycgen, ZC_neg(gel(U2,i)));
195 2778 : gel(Sunit,i) = famat_reduce(famat_mul(a, b));
196 : }
197 : }
198 850 : H = ZM_inv(H, pden);
199 850 : *pA = shallowconcat(H, ZM_neg(ZM_mul(H,B))); /* top inverse * den */
200 850 : return Sunit;
201 : }
202 : GEN
203 130 : bnfsunit(GEN bnf,GEN S,long prec)
204 : {
205 130 : pari_sp av = avma;
206 130 : long i, l = lg(S);
207 130 : GEN v, R, nf, A, den, U, cl, H = NULL;
208 130 : bnf = checkbnf(bnf); nf = bnf_get_nf(bnf);
209 130 : v = cgetg(7, t_VEC);
210 130 : gel(v,1) = U = bnfsunit_i(bnf, S, &H, &A, &den);
211 130 : gel(v,2) = mkvec2(A, den);
212 130 : gel(v,3) = cgetg(1,t_VEC); /* dummy */
213 130 : R = bnf_get_reg(bnf);
214 130 : cl = bnf_get_clgp(bnf);
215 130 : if (l != 1)
216 : {
217 124 : GEN u,A, G = bnf_get_gen(bnf), D = ZM_snf_group(H,NULL,&u), h = ZV_prod(D);
218 124 : long lD = lg(D);
219 124 : A = cgetg(lD, t_VEC);
220 140 : for(i = 1; i < lD; i++) gel(A,i) = idealfactorback(nf, G, gel(u,i), 1);
221 124 : cl = mkvec3(h, D, A);
222 124 : R = mpmul(R, h);
223 1459 : for (i = 1; i < l; i++)
224 : {
225 1335 : GEN pr = gel(S,i), p = pr_get_p(pr);
226 1335 : long f = pr_get_f(pr);
227 1335 : R = mpmul(R, logr_abs(itor(p,prec)));
228 1335 : if (f != 1) R = mulru(R, f);
229 1335 : gel(U,i) = nf_to_scalar_or_alg(nf, gel(U,i));
230 : }
231 : }
232 130 : gel(v,4) = R;
233 130 : gel(v,5) = cl;
234 130 : gel(v,6) = S; return gc_GEN(av, v);
235 : }
236 : GEN
237 942 : bnfunits(GEN bnf, GEN S)
238 : {
239 942 : pari_sp av = avma;
240 : GEN A, den, U, fu, tu;
241 942 : bnf = checkbnf(bnf);
242 942 : U = bnfsunit_i(bnf, S? S: cgetg(1,t_VEC), NULL, &A, &den);
243 942 : if (!S) S = cgetg(1,t_VEC);
244 942 : fu = bnf_compactfu(bnf);
245 942 : if (!fu)
246 : {
247 : long i, l;
248 204 : fu = bnf_has_fu(bnf); if (!fu) bnf_build_units(bnf);
249 204 : fu = shallowcopy(fu); l = lg(fu);
250 594 : for (i = 1; i < l; i++) gel(fu,i) = to_famat_shallow(gel(fu,i),gen_1);
251 : }
252 942 : tu = nf_to_scalar_or_basis(bnf_get_nf(bnf), bnf_get_tuU(bnf));
253 942 : U = shallowconcat(U, vec_append(fu, to_famat_shallow(tu,gen_1)));
254 942 : return gc_GEN(av, mkvec4(U, S, A, den));
255 : }
256 : GEN
257 156 : sunits_mod_units(GEN bnf, GEN S)
258 : {
259 156 : pari_sp av = avma;
260 : GEN A, den;
261 156 : bnf = checkbnf(bnf);
262 156 : return gc_GEN(av, bnfsunit_i(bnf, S, NULL, &A, &den));
263 : }
264 :
265 : /* v_S(x), x in famat form */
266 : static GEN
267 90 : sunit_famat_val(GEN nf, GEN S, GEN x)
268 : {
269 90 : long i, l = lg(S);
270 90 : GEN v = cgetg(l, t_VEC);
271 930 : for (i = 1; i < l; i++) gel(v,i) = famat_nfvalrem(nf, x, gel(S,i), NULL);
272 90 : return v;
273 : }
274 : /* v_S(x) */
275 : static GEN
276 1077 : sunit_val(GEN nf, GEN S, GEN x)
277 : {
278 1077 : long i, l = lg(S);
279 1077 : GEN v, Nx, dx, Ndx = gen_1;
280 1077 : if (isintzero(x)) return NULL;
281 1077 : v = zerocol(l-1);
282 1077 : x = Q_remove_denom(x, &dx);
283 1077 : Nx = idealnorm(nf, x);
284 1077 : if (dx) Ndx = powiu(dx, nf_get_degree(nf)); else if (is_pm1(Nx)) return v;
285 26552 : for (i = 1; i < l; i++)
286 : {
287 25670 : GEN P = gel(S,i), p = pr_get_p(P);
288 25670 : long f = pr_get_f(P), vn = 0, vd = 0;
289 25670 : if (dvdii(Nx, p))
290 : {
291 1564 : vn = nfval(nf, x, P);
292 1564 : if (vn) Nx = diviiexact(Nx, powiu(p, vn * f));
293 : }
294 25670 : if (dx && dvdii(dx, p))
295 : {
296 28 : vd = pr_get_e(P) * Z_pval(dx, p);
297 28 : if (vd) Ndx = diviiexact(Ndx, powiu(p, vd * f));
298 : }
299 25670 : vn -= vd; if (vn) gel(v,i) = stoi(vn);
300 : }
301 882 : return is_pm1(Nx) && is_pm1(Ndx)? v: NULL;
302 : }
303 :
304 : /* if *px a famat, assume it's an S-unit */
305 : static GEN
306 1191 : make_unit(GEN nf, GEN U, GEN *px)
307 : {
308 1191 : GEN den, v, w, A, gen = gel(U,1), S = gel(U,2), x = *px;
309 1191 : long cH, i, l = lg(S);
310 :
311 1191 : if (l == 1) return cgetg(1, t_COL);
312 1167 : A = gel(U,3); den = gel(U,4);
313 1167 : cH = nbrows(A);
314 1167 : if (typ(x) == t_MAT && lg(x) == 3)
315 90 : w = sunit_famat_val(nf, S, x); /* x = S v */
316 : else
317 : {
318 1077 : x = nf_to_scalar_or_basis(nf,x);
319 1077 : w = sunit_val(nf, S, x);
320 1077 : if (!w) return NULL;
321 : }
322 1149 : v = ZM_ZC_mul(A, w);
323 1149 : w += cH; w[0] = evaltyp(t_COL) | _evallg(lg(A) - cH);
324 2163 : if (!is_pm1(den)) for (i = 1; i <= cH; i++)
325 : {
326 : GEN r;
327 1014 : gel(v,i) = dvmdii(gel(v,i), den, &r);
328 1014 : if (r != gen_0) return NULL;
329 : }
330 1149 : v = shallowconcat(v, w); /* append bottom of w (= [0 Id] part) */
331 29606 : for (i = 1; i < l; i++)
332 : {
333 28457 : GEN e = gel(v,i);
334 28457 : if (signe(e)) x = famat_mulpow_shallow(x, gel(gen,i), negi(e));
335 : }
336 1149 : *px = x; return v;
337 : }
338 :
339 : static GEN
340 1191 : bnfissunit_i(GEN bnf, GEN x, GEN U)
341 : {
342 1191 : GEN w, nf, v = NULL;
343 1191 : bnf = checkbnf(bnf);
344 1191 : nf = bnf_get_nf(bnf);
345 1191 : if ( (w = make_unit(nf, U, &x)) ) v = bnfisunit(bnf, x);
346 1191 : if (!v || lg(v) == 1) return NULL;
347 1161 : return mkvec2(v, w);
348 : }
349 : static int
350 108 : checkU(GEN U)
351 : {
352 108 : if (typ(U) != t_VEC || lg(U) != 5) return 0;
353 108 : return typ(gel(U,1)) == t_VEC && is_vec_t(typ(gel(U,2)))
354 216 : && typ(gel(U,4))==t_INT;
355 : }
356 : GEN
357 138 : bnfisunit0(GEN bnf, GEN x, GEN U)
358 : {
359 138 : pari_sp av = avma;
360 : GEN z;
361 138 : if (!U) return bnfisunit(bnf, x);
362 108 : if (!checkU(U)) pari_err_TYPE("bnfisunit",U);
363 108 : z = bnfissunit_i(bnf, x, U);
364 108 : if (!z) retgc_const(av, cgetg(1, t_COL));
365 96 : return gc_GEN(av, shallowconcat(gel(z,2), gel(z,1)));
366 : }
367 :
368 : /* OBSOLETE */
369 : static int
370 1083 : checkbnfS_i(GEN v)
371 : {
372 : GEN S, g, w;
373 1083 : if (typ(v) != t_VEC || lg(v) != 7) return 0;
374 1083 : g = gel(v,1); w = gel(v,2); S = gel(v,6);
375 1083 : if (typ(g) != t_VEC || !is_vec_t(typ(S)) || lg(S) != lg(g)) return 0;
376 1083 : return typ(w) == t_VEC && lg(w) == 3;
377 : }
378 : /* OBSOLETE */
379 : GEN
380 1083 : bnfissunit(GEN bnf, GEN bnfS, GEN x)
381 : {
382 1083 : pari_sp av = avma;
383 : GEN z, U;
384 1083 : if (!checkbnfS_i(bnfS)) pari_err_TYPE("bnfissunit",bnfS);
385 1083 : U = mkvec4(gel(bnfS,1), gel(bnfS,6), gmael(bnfS,2,1), gmael(bnfS,2,2));
386 1083 : z = bnfissunit_i(bnf, x, U);
387 1083 : if (!z) retgc_const(av, cgetg(1, t_COL));
388 1065 : return gc_GEN(av, shallowconcat(gel(z,1), gel(z,2)));
389 : }
|
