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 : /* BASIC NF OPERATIONS */
18 : /* */
19 : /*******************************************************************/
20 : #include "pari.h"
21 : #include "paripriv.h"
22 :
23 : #define DEBUGLEVEL DEBUGLEVEL_nf
24 :
25 : /*******************************************************************/
26 : /* */
27 : /* OPERATIONS OVER NUMBER FIELD ELEMENTS. */
28 : /* represented as column vectors over the integral basis */
29 : /* */
30 : /*******************************************************************/
31 : static GEN
32 39661790 : get_tab(GEN nf, long *N)
33 : {
34 39661790 : GEN tab = (typ(nf) == t_MAT)? nf: gel(nf,9);
35 39661790 : *N = nbrows(tab); return tab;
36 : }
37 :
38 : /* x != 0, y t_INT. Return x * y (not memory clean if x = 1) */
39 : static GEN
40 1075840270 : _mulii(GEN x, GEN y) {
41 1737741701 : return is_pm1(x)? (signe(x) < 0)? negi(y): y
42 1737577769 : : mulii(x, y);
43 : }
44 :
45 : GEN
46 21818 : tablemul_ei_ej(GEN M, long i, long j)
47 : {
48 : long N;
49 21818 : GEN tab = get_tab(M, &N);
50 21818 : tab += (i-1)*N; return gel(tab,j);
51 : }
52 :
53 : /* Outputs x.ei, where ei is the i-th elt of the algebra basis.
54 : * x an RgV of correct length and arbitrary content (polynomials, scalars...).
55 : * M is the multiplication table ei ej = sum_k M_k^(i,j) ek */
56 : GEN
57 11137 : tablemul_ei(GEN M, GEN x, long i)
58 : {
59 : long j, k, N;
60 : GEN v, tab;
61 :
62 11137 : if (i==1) return gcopy(x);
63 11137 : tab = get_tab(M, &N);
64 11137 : if (typ(x) != t_COL) { v = zerocol(N); gel(v,i) = gcopy(x); return v; }
65 11137 : tab += (i-1)*N; v = cgetg(N+1,t_COL);
66 : /* wi . x = [ sum_j tab[k,j] x[j] ]_k */
67 76139 : for (k=1; k<=N; k++)
68 : {
69 65002 : pari_sp av = avma;
70 65002 : GEN s = gen_0;
71 462126 : for (j=1; j<=N; j++)
72 : {
73 397124 : GEN c = gcoeff(tab,k,j);
74 397124 : if (!gequal0(c)) s = gadd(s, gmul(c, gel(x,j)));
75 : }
76 65002 : gel(v,k) = gerepileupto(av,s);
77 : }
78 11137 : return v;
79 : }
80 : /* as tablemul_ei, assume x a ZV of correct length */
81 : GEN
82 23891452 : zk_ei_mul(GEN nf, GEN x, long i)
83 : {
84 : long j, k, N;
85 : GEN v, tab;
86 :
87 23891452 : if (i==1) return ZC_copy(x);
88 23891452 : tab = get_tab(nf, &N); tab += (i-1)*N;
89 23891789 : v = cgetg(N+1,t_COL);
90 169385176 : for (k=1; k<=N; k++)
91 : {
92 145497330 : pari_sp av = avma;
93 145497330 : GEN s = gen_0;
94 2139579567 : for (j=1; j<=N; j++)
95 : {
96 1994241852 : GEN c = gcoeff(tab,k,j);
97 1994241852 : if (signe(c)) s = addii(s, _mulii(c, gel(x,j)));
98 : }
99 145337715 : gel(v,k) = gerepileuptoint(av, s);
100 : }
101 23887846 : return v;
102 : }
103 :
104 : /* table of multiplication by wi in R[w1,..., wN] */
105 : GEN
106 39133 : ei_multable(GEN TAB, long i)
107 : {
108 : long k,N;
109 39133 : GEN m, tab = get_tab(TAB, &N);
110 39133 : tab += (i-1)*N;
111 39133 : m = cgetg(N+1,t_MAT);
112 153553 : for (k=1; k<=N; k++) gel(m,k) = gel(tab,k);
113 39133 : return m;
114 : }
115 :
116 : GEN
117 10816860 : zk_multable(GEN nf, GEN x)
118 : {
119 10816860 : long i, l = lg(x);
120 10816860 : GEN mul = cgetg(l,t_MAT);
121 10816772 : gel(mul,1) = x; /* assume w_1 = 1 */
122 34346803 : for (i=2; i<l; i++) gel(mul,i) = zk_ei_mul(nf,x,i);
123 10813017 : return mul;
124 : }
125 : GEN
126 2611 : multable(GEN M, GEN x)
127 : {
128 : long i, N;
129 : GEN mul;
130 2611 : if (typ(x) == t_MAT) return x;
131 0 : M = get_tab(M, &N);
132 0 : if (typ(x) != t_COL) return scalarmat(x, N);
133 0 : mul = cgetg(N+1,t_MAT);
134 0 : gel(mul,1) = x; /* assume w_1 = 1 */
135 0 : for (i=2; i<=N; i++) gel(mul,i) = tablemul_ei(M,x,i);
136 0 : return mul;
137 : }
138 :
139 : /* x integral in nf; table of multiplication by x in ZK = Z[w1,..., wN].
140 : * Return a t_INT if x is scalar, and a ZM otherwise */
141 : GEN
142 4992220 : zk_scalar_or_multable(GEN nf, GEN x)
143 : {
144 4992220 : long tx = typ(x);
145 4992220 : if (tx == t_MAT || tx == t_INT) return x;
146 4831003 : x = nf_to_scalar_or_basis(nf, x);
147 4830889 : return (typ(x) == t_COL)? zk_multable(nf, x): x;
148 : }
149 :
150 : GEN
151 21300 : nftrace(GEN nf, GEN x)
152 : {
153 21300 : pari_sp av = avma;
154 21300 : nf = checknf(nf);
155 21300 : x = nf_to_scalar_or_basis(nf, x);
156 21285 : x = (typ(x) == t_COL)? RgV_dotproduct(x, gel(nf_get_Tr(nf),1))
157 21306 : : gmulgu(x, nf_get_degree(nf));
158 21298 : return gerepileupto(av, x);
159 : }
160 : GEN
161 1015 : rnfelttrace(GEN rnf, GEN x)
162 : {
163 1015 : pari_sp av = avma;
164 1015 : checkrnf(rnf);
165 : /* avoid rnfabstorel special t_POL case misinterpretation */
166 1008 : if (typ(x) == t_POL && varn(x) == rnf_get_varn(rnf))
167 63 : x = gmodulo(x, rnf_get_pol(rnf));
168 1008 : x = rnfeltabstorel(rnf, x);
169 693 : x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gtrace(x))
170 798 : : gmulgu(x, rnf_get_degree(rnf));
171 798 : return gerepileupto(av, x);
172 : }
173 :
174 : static GEN
175 35 : famatQ_to_famatZ(GEN fa)
176 : {
177 35 : GEN E, F, Q, P = gel(fa,1);
178 35 : long i, j, l = lg(P);
179 35 : if (l == 1 || RgV_is_ZV(P)) return fa;
180 7 : Q = cgetg(2*l, t_COL);
181 7 : F = cgetg(2*l, t_COL); E = gel(fa, 2);
182 35 : for (i = j = 1; i < l; i++)
183 : {
184 28 : GEN p = gel(P,i);
185 28 : if (typ(p) == t_INT)
186 14 : { gel(Q, j) = p; gel(F, j) = gel(E, i); j++; }
187 : else
188 : {
189 14 : gel(Q, j) = gel(p,1); gel(F, j) = gel(E, i); j++;
190 14 : gel(Q, j) = gel(p,2); gel(F, j) = negi(gel(E, i)); j++;
191 : }
192 : }
193 7 : setlg(Q, j); setlg(F, j); return mkmat2(Q, F);
194 : }
195 : static GEN
196 35 : famat_cba(GEN fa)
197 : {
198 35 : GEN Q, F, P = gel(fa, 1), E = gel(fa, 2);
199 35 : long i, j, lQ, l = lg(P);
200 35 : if (l == 1) return fa;
201 28 : Q = ZV_cba(P); lQ = lg(Q); settyp(Q, t_COL);
202 28 : F = cgetg(lQ, t_COL);
203 77 : for (j = 1; j < lQ; j++)
204 : {
205 49 : GEN v = gen_0, q = gel(Q,j);
206 49 : if (!equali1(q))
207 203 : for (i = 1; i < l; i++)
208 : {
209 161 : long e = Z_pval(gel(P,i), q);
210 161 : if (e) v = addii(v, muliu(gel(E,i), e));
211 : }
212 49 : gel(F, j) = v;
213 : }
214 28 : return mkmat2(Q, F);
215 : }
216 : static long
217 35 : famat_sign(GEN fa)
218 : {
219 35 : GEN P = gel(fa,1), E = gel(fa,2);
220 35 : long i, l = lg(P), s = 1;
221 126 : for (i = 1; i < l; i++)
222 91 : if (signe(gel(P,i)) < 0 && mpodd(gel(E,i))) s = -s;
223 35 : return s;
224 : }
225 : static GEN
226 35 : famat_abs(GEN fa)
227 : {
228 35 : GEN Q, P = gel(fa,1);
229 : long i, l;
230 35 : Q = cgetg_copy(P, &l);
231 126 : for (i = 1; i < l; i++) gel(Q,i) = absi_shallow(gel(P,i));
232 35 : return mkmat2(Q, gel(fa,2));
233 : }
234 :
235 : /* assume nf is a genuine nf, fa a famat */
236 : static GEN
237 35 : famat_norm(GEN nf, GEN fa)
238 : {
239 35 : pari_sp av = avma;
240 35 : GEN G, g = gel(fa,1);
241 : long i, l, s;
242 :
243 35 : G = cgetg_copy(g, &l);
244 112 : for (i = 1; i < l; i++) gel(G,i) = nfnorm(nf, gel(g,i));
245 35 : fa = mkmat2(G, gel(fa,2));
246 35 : fa = famatQ_to_famatZ(fa);
247 35 : s = famat_sign(fa);
248 35 : fa = famat_reduce(famat_abs(fa));
249 35 : fa = famat_cba(fa);
250 35 : g = factorback(fa);
251 35 : return gerepileupto(av, s < 0? gneg(g): g);
252 : }
253 : GEN
254 222793 : nfnorm(GEN nf, GEN x)
255 : {
256 222793 : pari_sp av = avma;
257 : GEN c, den;
258 : long n;
259 222793 : nf = checknf(nf);
260 222793 : n = nf_get_degree(nf);
261 222793 : if (typ(x) == t_MAT) return famat_norm(nf, x);
262 222758 : x = nf_to_scalar_or_basis(nf, x);
263 222757 : if (typ(x)!=t_COL)
264 126728 : return gerepileupto(av, gpowgs(x, n));
265 96029 : x = nf_to_scalar_or_alg(nf, Q_primitive_part(x, &c));
266 96029 : x = Q_remove_denom(x, &den);
267 96030 : x = ZX_resultant_all(nf_get_pol(nf), x, den, 0);
268 96029 : return gerepileupto(av, c ? gmul(x, gpowgs(c, n)): x);
269 : }
270 :
271 : static GEN
272 119 : to_RgX(GEN P, long vx)
273 : {
274 119 : return varn(P) == vx ? P: scalarpol_shallow(P, vx);
275 : }
276 :
277 : GEN
278 462 : rnfeltnorm(GEN rnf, GEN x)
279 : {
280 462 : pari_sp av = avma;
281 : GEN nf, pol;
282 : long v;
283 462 : checkrnf(rnf);
284 455 : v = rnf_get_varn(rnf);
285 : /* avoid rnfabstorel special t_POL case misinterpretation */
286 455 : if (typ(x) == t_POL && varn(x) == v) x = gmodulo(x, rnf_get_pol(rnf));
287 455 : x = liftpol_shallow(rnfeltabstorel(rnf, x));
288 245 : nf = rnf_get_nf(rnf); pol = rnf_get_pol(rnf);
289 490 : x = (typ(x) == t_POL)
290 119 : ? rnfeltdown(rnf, nfX_resultant(nf,pol,to_RgX(x,v)))
291 245 : : gpowgs(x, rnf_get_degree(rnf));
292 245 : return gerepileupto(av, x);
293 : }
294 :
295 : /* x + y in nf */
296 : GEN
297 23401929 : nfadd(GEN nf, GEN x, GEN y)
298 : {
299 23401929 : pari_sp av = avma;
300 : GEN z;
301 :
302 23401929 : nf = checknf(nf);
303 23401929 : x = nf_to_scalar_or_basis(nf, x);
304 23401929 : y = nf_to_scalar_or_basis(nf, y);
305 23401929 : if (typ(x) != t_COL)
306 17642336 : { z = (typ(y) == t_COL)? RgC_Rg_add(y, x): gadd(x,y); }
307 : else
308 5759593 : { z = (typ(y) == t_COL)? RgC_add(x, y): RgC_Rg_add(x, y); }
309 23401929 : return gerepileupto(av, z);
310 : }
311 : /* x - y in nf */
312 : GEN
313 1809883 : nfsub(GEN nf, GEN x, GEN y)
314 : {
315 1809883 : pari_sp av = avma;
316 : GEN z;
317 :
318 1809883 : nf = checknf(nf);
319 1809883 : x = nf_to_scalar_or_basis(nf, x);
320 1809883 : y = nf_to_scalar_or_basis(nf, y);
321 1809883 : if (typ(x) != t_COL)
322 1278354 : { z = (typ(y) == t_COL)? Rg_RgC_sub(x,y): gsub(x,y); }
323 : else
324 531529 : { z = (typ(y) == t_COL)? RgC_sub(x,y): RgC_Rg_sub(x,y); }
325 1809883 : return gerepileupto(av, z);
326 : }
327 :
328 : /* product of ZC x,y in (true) nf; ( sum_i x_i sum_j y_j m^{i,j}_k )_k */
329 : static GEN
330 8604550 : nfmuli_ZC(GEN nf, GEN x, GEN y)
331 : {
332 : long i, j, k, N;
333 8604550 : GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
334 :
335 41996987 : for (k = 1; k <= N; k++)
336 : {
337 33392526 : pari_sp av = avma;
338 33392526 : GEN s, TABi = TAB;
339 33392526 : if (k == 1)
340 8604532 : s = mulii(gel(x,1),gel(y,1));
341 : else
342 24787838 : s = addii(mulii(gel(x,1),gel(y,k)),
343 24787994 : mulii(gel(x,k),gel(y,1)));
344 217870705 : for (i=2; i<=N; i++)
345 : {
346 184482172 : GEN t, xi = gel(x,i);
347 184482172 : TABi += N;
348 184482172 : if (!signe(xi)) continue;
349 :
350 93467458 : t = NULL;
351 1063171783 : for (j=2; j<=N; j++)
352 : {
353 969705258 : GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
354 969705258 : if (!signe(c)) continue;
355 281215372 : p1 = _mulii(c, gel(y,j));
356 281220719 : t = t? addii(t, p1): p1;
357 : }
358 93466525 : if (t) s = addii(s, mulii(xi, t));
359 : }
360 33388533 : gel(v,k) = gerepileuptoint(av,s);
361 : }
362 8604461 : return v;
363 : }
364 : static int
365 74531399 : is_famat(GEN x) { return typ(x) == t_MAT && lg(x) == 3; }
366 : /* product of x and y in nf */
367 : GEN
368 36262880 : nfmul(GEN nf, GEN x, GEN y)
369 : {
370 : GEN z;
371 36262880 : pari_sp av = avma;
372 :
373 36262880 : if (x == y) return nfsqr(nf,x);
374 :
375 32184502 : nf = checknf(nf);
376 32184501 : if (is_famat(x) || is_famat(y)) return famat_mul(x, y);
377 32184191 : x = nf_to_scalar_or_basis(nf, x);
378 32184193 : y = nf_to_scalar_or_basis(nf, y);
379 32184196 : if (typ(x) != t_COL)
380 : {
381 21772238 : if (isintzero(x)) return gen_0;
382 15729753 : z = (typ(y) == t_COL)? RgC_Rg_mul(y, x): gmul(x,y); }
383 : else
384 : {
385 10411958 : if (typ(y) != t_COL)
386 : {
387 4533560 : if (isintzero(y)) return gen_0;
388 1607667 : z = RgC_Rg_mul(x, y);
389 : }
390 : else
391 : {
392 : GEN dx, dy;
393 5878398 : x = Q_remove_denom(x, &dx);
394 5878398 : y = Q_remove_denom(y, &dy);
395 5878398 : z = nfmuli_ZC(nf,x,y);
396 5878399 : dx = mul_denom(dx,dy);
397 5878399 : if (dx) z = ZC_Z_div(z, dx);
398 : }
399 : }
400 23215812 : return gerepileupto(av, z);
401 : }
402 : /* square of ZC x in nf */
403 : static GEN
404 7096127 : nfsqri_ZC(GEN nf, GEN x)
405 : {
406 : long i, j, k, N;
407 7096127 : GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
408 :
409 38890353 : for (k = 1; k <= N; k++)
410 : {
411 31794246 : pari_sp av = avma;
412 31794246 : GEN s, TABi = TAB;
413 31794246 : if (k == 1)
414 7096311 : s = sqri(gel(x,1));
415 : else
416 24697935 : s = shifti(mulii(gel(x,1),gel(x,k)), 1);
417 253509925 : for (i=2; i<=N; i++)
418 : {
419 221732338 : GEN p1, c, t, xi = gel(x,i);
420 221732338 : TABi += N;
421 221732338 : if (!signe(xi)) continue;
422 :
423 79856580 : c = gcoeff(TABi, k, i);
424 79856580 : t = signe(c)? _mulii(c,xi): NULL;
425 675787763 : for (j=i+1; j<=N; j++)
426 : {
427 595930195 : c = gcoeff(TABi, k, j);
428 595930195 : if (!signe(c)) continue;
429 231907241 : p1 = _mulii(c, shifti(gel(x,j),1));
430 231912365 : t = t? addii(t, p1): p1;
431 : }
432 79857568 : if (t) s = addii(s, mulii(xi, t));
433 : }
434 31777587 : gel(v,k) = gerepileuptoint(av,s);
435 : }
436 7096107 : return v;
437 : }
438 : /* square of x in nf */
439 : GEN
440 8896169 : nfsqr(GEN nf, GEN x)
441 : {
442 8896169 : pari_sp av = avma;
443 : GEN z;
444 :
445 8896169 : nf = checknf(nf);
446 8896169 : if (is_famat(x)) return famat_sqr(x);
447 8896171 : x = nf_to_scalar_or_basis(nf, x);
448 8896172 : if (typ(x) != t_COL) z = gsqr(x);
449 : else
450 : {
451 : GEN dx;
452 2631242 : x = Q_remove_denom(x, &dx);
453 2631240 : z = nfsqri_ZC(nf,x);
454 2631246 : if (dx) z = RgC_Rg_div(z, sqri(dx));
455 : }
456 8896176 : return gerepileupto(av, z);
457 : }
458 :
459 : /* x a ZC, v a t_COL of ZC/Z */
460 : GEN
461 205488 : zkC_multable_mul(GEN v, GEN x)
462 : {
463 205488 : long i, l = lg(v);
464 205488 : GEN y = cgetg(l, t_COL);
465 799570 : for (i = 1; i < l; i++)
466 : {
467 594082 : GEN c = gel(v,i);
468 594082 : if (typ(c)!=t_COL) {
469 0 : if (!isintzero(c)) c = ZC_Z_mul(gel(x,1), c);
470 : } else {
471 594082 : c = ZM_ZC_mul(x,c);
472 594082 : if (ZV_isscalar(c)) c = gel(c,1);
473 : }
474 594082 : gel(y,i) = c;
475 : }
476 205488 : return y;
477 : }
478 :
479 : GEN
480 56807 : nfC_multable_mul(GEN v, GEN x)
481 : {
482 56807 : long i, l = lg(v);
483 56807 : GEN y = cgetg(l, t_COL);
484 383543 : for (i = 1; i < l; i++)
485 : {
486 326736 : GEN c = gel(v,i);
487 326736 : if (typ(c)!=t_COL) {
488 272231 : if (!isintzero(c)) c = RgC_Rg_mul(gel(x,1), c);
489 : } else {
490 54505 : c = RgM_RgC_mul(x,c);
491 54505 : if (QV_isscalar(c)) c = gel(c,1);
492 : }
493 326736 : gel(y,i) = c;
494 : }
495 56807 : return y;
496 : }
497 :
498 : GEN
499 197908 : nfC_nf_mul(GEN nf, GEN v, GEN x)
500 : {
501 : long tx;
502 : GEN y;
503 :
504 197908 : x = nf_to_scalar_or_basis(nf, x);
505 197908 : tx = typ(x);
506 197908 : if (tx != t_COL)
507 : {
508 : long l, i;
509 149584 : if (tx == t_INT)
510 : {
511 140519 : long s = signe(x);
512 140519 : if (!s) return zerocol(lg(v)-1);
513 133157 : if (is_pm1(x)) return s > 0? leafcopy(v): RgC_neg(v);
514 : }
515 48468 : l = lg(v); y = cgetg(l, t_COL);
516 347242 : for (i=1; i < l; i++)
517 : {
518 298774 : GEN c = gel(v,i);
519 298774 : if (typ(c) != t_COL) c = gmul(c, x); else c = RgC_Rg_mul(c, x);
520 298774 : gel(y,i) = c;
521 : }
522 48468 : return y;
523 : }
524 : else
525 : {
526 : GEN dx;
527 48324 : x = zk_multable(nf, Q_remove_denom(x,&dx));
528 48324 : y = nfC_multable_mul(v, x);
529 48324 : return dx? RgC_Rg_div(y, dx): y;
530 : }
531 : }
532 : static GEN
533 10961 : mulbytab(GEN M, GEN c)
534 10961 : { return typ(c) == t_COL? RgM_RgC_mul(M,c): RgC_Rg_mul(gel(M,1), c); }
535 : GEN
536 2611 : tablemulvec(GEN M, GEN x, GEN v)
537 : {
538 : long l, i;
539 : GEN y;
540 :
541 2611 : if (typ(x) == t_COL && RgV_isscalar(x))
542 : {
543 0 : x = gel(x,1);
544 0 : return typ(v) == t_POL? RgX_Rg_mul(v,x): RgV_Rg_mul(v,x);
545 : }
546 2611 : x = multable(M, x); /* multiplication table by x */
547 2611 : y = cgetg_copy(v, &l);
548 2611 : if (typ(v) == t_POL)
549 : {
550 2611 : y[1] = v[1];
551 13572 : for (i=2; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
552 2611 : y = normalizepol(y);
553 : }
554 : else
555 : {
556 0 : for (i=1; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
557 : }
558 2611 : return y;
559 : }
560 :
561 : GEN
562 1260340 : zkmultable_capZ(GEN mx) { return Q_denom(zkmultable_inv(mx)); }
563 : GEN
564 1579066 : zkmultable_inv(GEN mx) { return ZM_gauss(mx, col_ei(lg(mx)-1,1)); }
565 : /* nf a true nf, x a ZC */
566 : GEN
567 318726 : zk_inv(GEN nf, GEN x) { return zkmultable_inv(zk_multable(nf,x)); }
568 :
569 : /* inverse of x in nf */
570 : GEN
571 238924 : nfinv(GEN nf, GEN x)
572 : {
573 238924 : pari_sp av = avma;
574 : GEN z;
575 :
576 238924 : nf = checknf(nf);
577 238924 : if (is_famat(x)) return famat_inv(x);
578 238924 : x = nf_to_scalar_or_basis(nf, x);
579 238924 : if (typ(x) == t_COL)
580 : {
581 : GEN d;
582 190543 : x = Q_remove_denom(x, &d);
583 190543 : z = zk_inv(nf, x);
584 190543 : if (d) z = RgC_Rg_mul(z, d);
585 : }
586 : else
587 48381 : z = ginv(x);
588 238924 : return gerepileupto(av, z);
589 : }
590 :
591 : /* quotient of x and y in nf */
592 : GEN
593 36279 : nfdiv(GEN nf, GEN x, GEN y)
594 : {
595 36279 : pari_sp av = avma;
596 : GEN z;
597 :
598 36279 : nf = checknf(nf);
599 36279 : if (is_famat(x) || is_famat(y)) return famat_div(x,y);
600 36188 : y = nf_to_scalar_or_basis(nf, y);
601 36188 : if (typ(y) != t_COL)
602 : {
603 22085 : x = nf_to_scalar_or_basis(nf, x);
604 22085 : z = (typ(x) == t_COL)? RgC_Rg_div(x, y): gdiv(x,y);
605 : }
606 : else
607 : {
608 : GEN d;
609 14103 : y = Q_remove_denom(y, &d);
610 14103 : z = nfmul(nf, x, zk_inv(nf,y));
611 14103 : if (d) z = typ(z) == t_COL? RgC_Rg_mul(z, d): gmul(z, d);
612 : }
613 36188 : return gerepileupto(av, z);
614 : }
615 :
616 : /* product of INTEGERS (t_INT or ZC) x and y in (true) nf */
617 : GEN
618 4097329 : nfmuli(GEN nf, GEN x, GEN y)
619 : {
620 4097329 : if (typ(x) == t_INT) return (typ(y) == t_COL)? ZC_Z_mul(y, x): mulii(x,y);
621 2958736 : if (typ(y) == t_INT) return ZC_Z_mul(x, y);
622 2726120 : return nfmuli_ZC(nf, x, y);
623 : }
624 : GEN
625 4464935 : nfsqri(GEN nf, GEN x)
626 4464935 : { return (typ(x) == t_INT)? sqri(x): nfsqri_ZC(nf, x); }
627 :
628 : /* both x and y are RgV */
629 : GEN
630 0 : tablemul(GEN TAB, GEN x, GEN y)
631 : {
632 : long i, j, k, N;
633 : GEN s, v;
634 0 : if (typ(x) != t_COL) return gmul(x, y);
635 0 : if (typ(y) != t_COL) return gmul(y, x);
636 0 : N = lg(x)-1;
637 0 : v = cgetg(N+1,t_COL);
638 0 : for (k=1; k<=N; k++)
639 : {
640 0 : pari_sp av = avma;
641 0 : GEN TABi = TAB;
642 0 : if (k == 1)
643 0 : s = gmul(gel(x,1),gel(y,1));
644 : else
645 0 : s = gadd(gmul(gel(x,1),gel(y,k)),
646 0 : gmul(gel(x,k),gel(y,1)));
647 0 : for (i=2; i<=N; i++)
648 : {
649 0 : GEN t, xi = gel(x,i);
650 0 : TABi += N;
651 0 : if (gequal0(xi)) continue;
652 :
653 0 : t = NULL;
654 0 : for (j=2; j<=N; j++)
655 : {
656 0 : GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
657 0 : if (gequal0(c)) continue;
658 0 : p1 = gmul(c, gel(y,j));
659 0 : t = t? gadd(t, p1): p1;
660 : }
661 0 : if (t) s = gadd(s, gmul(xi, t));
662 : }
663 0 : gel(v,k) = gerepileupto(av,s);
664 : }
665 0 : return v;
666 : }
667 : GEN
668 48614 : tablesqr(GEN TAB, GEN x)
669 : {
670 : long i, j, k, N;
671 : GEN s, v;
672 :
673 48614 : if (typ(x) != t_COL) return gsqr(x);
674 48614 : N = lg(x)-1;
675 48614 : v = cgetg(N+1,t_COL);
676 :
677 346630 : for (k=1; k<=N; k++)
678 : {
679 298016 : pari_sp av = avma;
680 298016 : GEN TABi = TAB;
681 298016 : if (k == 1)
682 48614 : s = gsqr(gel(x,1));
683 : else
684 249402 : s = gmul2n(gmul(gel(x,1),gel(x,k)), 1);
685 1898410 : for (i=2; i<=N; i++)
686 : {
687 1600394 : GEN p1, c, t, xi = gel(x,i);
688 1600394 : TABi += N;
689 1600394 : if (gequal0(xi)) continue;
690 :
691 416801 : c = gcoeff(TABi, k, i);
692 416801 : t = !gequal0(c)? gmul(c,xi): NULL;
693 1668394 : for (j=i+1; j<=N; j++)
694 : {
695 1251593 : c = gcoeff(TABi, k, j);
696 1251593 : if (gequal0(c)) continue;
697 643118 : p1 = gmul(gmul2n(c,1), gel(x,j));
698 643118 : t = t? gadd(t, p1): p1;
699 : }
700 416801 : if (t) s = gadd(s, gmul(xi, t));
701 : }
702 298016 : gel(v,k) = gerepileupto(av,s);
703 : }
704 48614 : return v;
705 : }
706 :
707 : static GEN
708 354916 : _mul(void *data, GEN x, GEN y) { return nfmuli((GEN)data,x,y); }
709 : static GEN
710 967083 : _sqr(void *data, GEN x) { return nfsqri((GEN)data,x); }
711 :
712 : /* Compute z^n in nf, left-shift binary powering */
713 : GEN
714 940094 : nfpow(GEN nf, GEN z, GEN n)
715 : {
716 940094 : pari_sp av = avma;
717 : long s;
718 : GEN x, cx;
719 :
720 940094 : if (typ(n)!=t_INT) pari_err_TYPE("nfpow",n);
721 940094 : nf = checknf(nf);
722 940093 : s = signe(n); if (!s) return gen_1;
723 940093 : if (is_famat(z)) return famat_pow(z, n);
724 879448 : x = nf_to_scalar_or_basis(nf, z);
725 879449 : if (typ(x) != t_COL) return powgi(x,n);
726 760150 : if (s < 0)
727 : { /* simplified nfinv */
728 : GEN d;
729 45717 : x = Q_remove_denom(x, &d);
730 45717 : x = zk_inv(nf, x);
731 45717 : x = primitive_part(x, &cx);
732 45717 : cx = mul_content(cx, d);
733 45717 : n = negi(n);
734 : }
735 : else
736 714433 : x = primitive_part(x, &cx);
737 760139 : x = gen_pow_i(x, n, (void*)nf, _sqr, _mul);
738 760147 : if (cx)
739 46733 : x = gerepileupto(av, gmul(x, powgi(cx, n)));
740 : else
741 713414 : x = gerepilecopy(av, x);
742 760161 : return x;
743 : }
744 : /* Compute z^n in nf, left-shift binary powering */
745 : GEN
746 345100 : nfpow_u(GEN nf, GEN z, ulong n)
747 : {
748 345100 : pari_sp av = avma;
749 : GEN x, cx;
750 :
751 345100 : if (!n) return gen_1;
752 345100 : x = nf_to_scalar_or_basis(nf, z);
753 345100 : if (typ(x) != t_COL) return gpowgs(x,n);
754 309207 : x = primitive_part(x, &cx);
755 309205 : x = gen_powu_i(x, n, (void*)nf, _sqr, _mul);
756 309205 : if (cx)
757 : {
758 114514 : x = gmul(x, powgi(cx, utoipos(n)));
759 114514 : return gerepileupto(av,x);
760 : }
761 194691 : return gerepilecopy(av, x);
762 : }
763 :
764 : long
765 588 : nfissquare(GEN nf, GEN z, GEN *px)
766 : {
767 588 : pari_sp av = avma;
768 588 : long v = fetch_var_higher();
769 : GEN R;
770 588 : nf = checknf(nf);
771 588 : if (nf_get_degree(nf) == 1)
772 : {
773 21 : z = algtobasis(nf, z);
774 21 : if (!issquareall(gel(z,1), px)) return gc_long(av, 0);
775 14 : if (px) *px = gerepileupto(av, *px); else set_avma(av);
776 14 : return 1;
777 : }
778 567 : z = nf_to_scalar_or_alg(nf, z);
779 567 : R = nfroots(nf, deg2pol_shallow(gen_m1, gen_0, z, v));
780 567 : delete_var(); if (lg(R) == 1) return gc_long(av, 0);
781 546 : if (px) *px = gerepilecopy(av, nf_to_scalar_or_basis(nf, gel(R,1)));
782 14 : else set_avma(av);
783 546 : return 1;
784 : }
785 :
786 : long
787 7709 : nfispower(GEN nf, GEN z, long n, GEN *px)
788 : {
789 7709 : pari_sp av = avma;
790 7709 : long v = fetch_var_higher();
791 : GEN R;
792 7709 : nf = checknf(nf);
793 7709 : if (nf_get_degree(nf) == 1)
794 : {
795 329 : z = algtobasis(nf, z);
796 329 : if (!ispower(gel(z,1), stoi(n), px)) return gc_long(av, 0);
797 147 : if (px) *px = gerepileupto(av, *px); else set_avma(av);
798 147 : return 1;
799 : }
800 7380 : if (n <= 0)
801 0 : pari_err_DOMAIN("nfeltispower","exponent","<=",gen_0,stoi(n));
802 7380 : z = nf_to_scalar_or_alg(nf, z);
803 7380 : if (n==1)
804 : {
805 0 : if (px) *px = gerepilecopy(av, z);
806 0 : return 1;
807 : }
808 7380 : R = nfroots(nf, gsub(pol_xn(n, v), z));
809 7380 : delete_var(); if (lg(R) == 1) return gc_long(av, 0);
810 3157 : if (px) *px = gerepilecopy(av, nf_to_scalar_or_basis(nf, gel(R,1)));
811 3143 : else set_avma(av);
812 3157 : return 1;
813 : }
814 :
815 : static GEN
816 56 : idmulred(void *nf, GEN x, GEN y) { return idealmulred((GEN) nf, x, y); }
817 : static GEN
818 413 : idpowred(void *nf, GEN x, GEN n) { return idealpowred((GEN) nf, x, n); }
819 : static GEN
820 70361 : idmul(void *nf, GEN x, GEN y) { return idealmul((GEN) nf, x, y); }
821 : static GEN
822 86123 : idpow(void *nf, GEN x, GEN n) { return idealpow((GEN) nf, x, n); }
823 : GEN
824 85177 : idealfactorback(GEN nf, GEN L, GEN e, int red)
825 : {
826 85177 : nf = checknf(nf);
827 85177 : if (red) return gen_factorback(L, e, (void*)nf, &idmulred, &idpowred, NULL);
828 84820 : if (!e && typ(L) == t_MAT && lg(L) == 3) { e = gel(L,2); L = gel(L,1); }
829 84820 : if (is_vec_t(typ(L)) && RgV_is_prV(L))
830 : { /* don't use gen_factorback since *= pr^v can be done more efficiently */
831 64396 : pari_sp av = avma;
832 64396 : long i, l = lg(L);
833 : GEN a;
834 64396 : if (!e) e = const_vec(l-1, gen_1);
835 61540 : else switch(typ(e))
836 : {
837 7264 : case t_VECSMALL: e = zv_to_ZV(e); break;
838 54276 : case t_VEC: case t_COL:
839 54276 : if (!RgV_is_ZV(e))
840 0 : pari_err_TYPE("factorback [not an exponent vector]", e);
841 54276 : break;
842 0 : default: pari_err_TYPE("idealfactorback", e);
843 : }
844 64396 : if (l != lg(e))
845 0 : pari_err_TYPE("factorback [not an exponent vector]", e);
846 64396 : if (l == 1 || ZV_equal0(e)) return gc_const(av, gen_1);
847 23289 : a = idealpow(nf, gel(L,1), gel(e,1));
848 242883 : for (i = 2; i < l; i++)
849 219594 : if (signe(gel(e,i))) a = idealmulpowprime(nf, a, gel(L,i), gel(e,i));
850 23289 : return gerepileupto(av, a);
851 : }
852 20424 : return gen_factorback(L, e, (void*)nf, &idmul, &idpow, NULL);
853 : }
854 : static GEN
855 327448 : eltmul(void *nf, GEN x, GEN y) { return nfmul((GEN) nf, x, y); }
856 : static GEN
857 464370 : eltpow(void *nf, GEN x, GEN n) { return nfpow((GEN) nf, x, n); }
858 : GEN
859 264986 : nffactorback(GEN nf, GEN L, GEN e)
860 264986 : { return gen_factorback(L, e, (void*)checknf(nf), &eltmul, &eltpow, NULL); }
861 :
862 : static GEN
863 3071451 : _nf_red(void *E, GEN x) { (void)E; return gcopy(x); }
864 :
865 : static GEN
866 12603562 : _nf_add(void *E, GEN x, GEN y) { return nfadd((GEN)E,x,y); }
867 :
868 : static GEN
869 745719 : _nf_neg(void *E, GEN x) { (void)E; return gneg(x); }
870 :
871 : static GEN
872 15127605 : _nf_mul(void *E, GEN x, GEN y) { return nfmul((GEN)E,x,y); }
873 :
874 : static GEN
875 52951 : _nf_inv(void *E, GEN x) { return nfinv((GEN)E,x); }
876 :
877 : static GEN
878 10799 : _nf_s(void *E, long x) { (void)E; return stoi(x); }
879 :
880 : static const struct bb_field nf_field={_nf_red,_nf_add,_nf_mul,_nf_neg,
881 : _nf_inv,&gequal0,_nf_s };
882 :
883 225006 : const struct bb_field *get_nf_field(void **E, GEN nf)
884 225006 : { *E = (void*)nf; return &nf_field; }
885 :
886 : GEN
887 14 : nfM_det(GEN nf, GEN M)
888 : {
889 : void *E;
890 14 : const struct bb_field *S = get_nf_field(&E, nf);
891 14 : return gen_det(M, E, S);
892 : }
893 : GEN
894 10785 : nfM_inv(GEN nf, GEN M)
895 : {
896 : void *E;
897 10785 : const struct bb_field *S = get_nf_field(&E, nf);
898 10785 : return gen_Gauss(M, matid(lg(M)-1), E, S);
899 : }
900 :
901 : GEN
902 0 : nfM_ker(GEN nf, GEN M)
903 : {
904 : void *E;
905 0 : const struct bb_field *S = get_nf_field(&E, nf);
906 0 : return gen_ker(M, 0, E, S);
907 : }
908 :
909 : GEN
910 10316 : nfM_mul(GEN nf, GEN A, GEN B)
911 : {
912 : void *E;
913 10316 : const struct bb_field *S = get_nf_field(&E, nf);
914 10316 : return gen_matmul(A, B, E, S);
915 : }
916 : GEN
917 203891 : nfM_nfC_mul(GEN nf, GEN A, GEN B)
918 : {
919 : void *E;
920 203891 : const struct bb_field *S = get_nf_field(&E, nf);
921 203891 : return gen_matcolmul(A, B, E, S);
922 : }
923 :
924 : /* valuation of integral x (ZV), with resp. to prime ideal pr */
925 : long
926 24017020 : ZC_nfvalrem(GEN x, GEN pr, GEN *newx)
927 : {
928 24017020 : pari_sp av = avma;
929 : long i, v, l;
930 24017020 : GEN r, y, p = pr_get_p(pr), mul = pr_get_tau(pr);
931 :
932 : /* p inert */
933 24017034 : if (typ(mul) == t_INT) return newx? ZV_pvalrem(x, p, newx):ZV_pval(x, p);
934 23013284 : y = cgetg_copy(x, &l); /* will hold the new x */
935 23013636 : x = leafcopy(x);
936 37176372 : for(v=0;; v++)
937 : {
938 143008539 : for (i=1; i<l; i++)
939 : { /* is (x.b)[i] divisible by p ? */
940 128840291 : gel(y,i) = dvmdii(ZMrow_ZC_mul(mul,x,i),p,&r);
941 128843807 : if (r != gen_0) { if (newx) *newx = x; return v; }
942 : }
943 14168248 : swap(x, y);
944 14168248 : if (!newx && (v & 0xf) == 0xf) v += pr_get_e(pr) * ZV_pvalrem(x, p, &x);
945 14168248 : if (gc_needed(av,1))
946 : {
947 0 : if(DEBUGMEM>1) pari_warn(warnmem,"ZC_nfvalrem, v >= %ld", v);
948 0 : gerepileall(av, 2, &x, &y);
949 : }
950 : }
951 : }
952 : long
953 19746188 : ZC_nfval(GEN x, GEN P)
954 19746188 : { return ZC_nfvalrem(x, P, NULL); }
955 :
956 : /* v_P(x) != 0, x a ZV. Simpler version of ZC_nfvalrem */
957 : int
958 1246559 : ZC_prdvd(GEN x, GEN P)
959 : {
960 1246559 : pari_sp av = avma;
961 : long i, l;
962 1246559 : GEN p = pr_get_p(P), mul = pr_get_tau(P);
963 1246593 : if (typ(mul) == t_INT) return ZV_Z_dvd(x, p);
964 1246054 : l = lg(x);
965 5050343 : for (i=1; i<l; i++)
966 4534617 : if (!dvdii(ZMrow_ZC_mul(mul,x,i), p)) return gc_bool(av,0);
967 515726 : return gc_bool(av,1);
968 : }
969 :
970 : int
971 357 : pr_equal(GEN P, GEN Q)
972 : {
973 357 : GEN gQ, p = pr_get_p(P);
974 357 : long e = pr_get_e(P), f = pr_get_f(P), n;
975 357 : if (!equalii(p, pr_get_p(Q)) || e != pr_get_e(Q) || f != pr_get_f(Q))
976 336 : return 0;
977 21 : gQ = pr_get_gen(Q); n = lg(gQ)-1;
978 21 : if (2*e*f > n) return 1; /* room for only one such pr */
979 14 : return ZV_equal(pr_get_gen(P), gQ) || ZC_prdvd(gQ, P);
980 : }
981 :
982 : GEN
983 420721 : famat_nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
984 : {
985 420721 : pari_sp av = avma;
986 420721 : GEN P = gel(x,1), E = gel(x,2), V = gen_0, y = NULL;
987 420721 : long l = lg(P), simplify = 0, i;
988 420721 : if (py) { *py = gen_1; y = cgetg(l, t_COL); }
989 :
990 2258510 : for (i = 1; i < l; i++)
991 : {
992 1837789 : GEN e = gel(E,i);
993 : long v;
994 1837789 : if (!signe(e))
995 : {
996 7 : if (py) gel(y,i) = gen_1;
997 7 : simplify = 1; continue;
998 : }
999 1837782 : v = nfvalrem(nf, gel(P,i), pr, py? &gel(y,i): NULL);
1000 1837782 : if (v == LONG_MAX) { set_avma(av); if (py) *py = gen_0; return mkoo(); }
1001 1837782 : V = addmulii(V, stoi(v), e);
1002 : }
1003 420721 : if (!py) V = gerepileuptoint(av, V);
1004 : else
1005 : {
1006 42 : y = mkmat2(y, gel(x,2));
1007 42 : if (simplify) y = famat_remove_trivial(y);
1008 42 : gerepileall(av, 2, &V, &y); *py = y;
1009 : }
1010 420721 : return V;
1011 : }
1012 : long
1013 5621848 : nfval(GEN nf, GEN x, GEN pr)
1014 : {
1015 5621848 : pari_sp av = avma;
1016 : long w, e;
1017 : GEN cx, p;
1018 :
1019 5621848 : if (gequal0(x)) return LONG_MAX;
1020 5608530 : nf = checknf(nf);
1021 5608526 : checkprid(pr);
1022 5608518 : p = pr_get_p(pr);
1023 5608517 : e = pr_get_e(pr);
1024 5608512 : x = nf_to_scalar_or_basis(nf, x);
1025 5608428 : if (typ(x) != t_COL) return e*Q_pval(x,p);
1026 2376562 : x = Q_primitive_part(x, &cx);
1027 2376617 : w = ZC_nfval(x,pr);
1028 2376564 : if (cx) w += e*Q_pval(cx,p);
1029 2376562 : return gc_long(av,w);
1030 : }
1031 :
1032 : /* want to write p^v = uniformizer^(e*v) * z^v, z coprime to pr */
1033 : /* z := tau^e / p^(e-1), algebraic integer coprime to pr; return z^v */
1034 : static GEN
1035 951048 : powp(GEN nf, GEN pr, long v)
1036 : {
1037 : GEN b, z;
1038 : long e;
1039 951048 : if (!v) return gen_1;
1040 424473 : b = pr_get_tau(pr);
1041 424473 : if (typ(b) == t_INT) return gen_1;
1042 121898 : e = pr_get_e(pr);
1043 121898 : z = gel(b,1);
1044 121898 : if (e != 1) z = gdiv(nfpow_u(nf, z, e), powiu(pr_get_p(pr),e-1));
1045 121898 : if (v < 0) { v = -v; z = nfinv(nf, z); }
1046 121898 : if (v != 1) z = nfpow_u(nf, z, v);
1047 121898 : return z;
1048 : }
1049 : long
1050 3639128 : nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
1051 : {
1052 3639128 : pari_sp av = avma;
1053 : long w, e;
1054 : GEN cx, p, t;
1055 :
1056 3639128 : if (!py) return nfval(nf,x,pr);
1057 1787955 : if (gequal0(x)) { *py = gen_0; return LONG_MAX; }
1058 1787899 : nf = checknf(nf);
1059 1787899 : checkprid(pr);
1060 1787899 : p = pr_get_p(pr);
1061 1787898 : e = pr_get_e(pr);
1062 1787898 : x = nf_to_scalar_or_basis(nf, x);
1063 1787900 : if (typ(x) != t_COL) {
1064 538531 : w = Q_pvalrem(x,p, py);
1065 538531 : if (!w) { *py = gerepilecopy(av, x); return 0; }
1066 330239 : *py = gerepileupto(av, gmul(powp(nf, pr, w), *py));
1067 330239 : return e*w;
1068 : }
1069 1249369 : x = Q_primitive_part(x, &cx);
1070 1249366 : w = ZC_nfvalrem(x,pr, py);
1071 1249359 : if (cx)
1072 : {
1073 620809 : long v = Q_pvalrem(cx,p, &t);
1074 620809 : *py = nfmul(nf, *py, gmul(powp(nf,pr,v), t));
1075 620809 : *py = gerepileupto(av, *py);
1076 620809 : w += e*v;
1077 : }
1078 : else
1079 628550 : *py = gerepilecopy(av, *py);
1080 1249373 : return w;
1081 : }
1082 : GEN
1083 15015 : gpnfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
1084 : {
1085 : long v;
1086 15015 : if (is_famat(x)) return famat_nfvalrem(nf, x, pr, py);
1087 15008 : v = nfvalrem(nf,x,pr,py);
1088 15008 : return v == LONG_MAX? mkoo(): stoi(v);
1089 : }
1090 :
1091 : /* true nf */
1092 : GEN
1093 305620 : coltoalg(GEN nf, GEN x)
1094 : {
1095 305620 : return mkpolmod( nf_to_scalar_or_alg(nf, x), nf_get_pol(nf) );
1096 : }
1097 :
1098 : GEN
1099 358712 : basistoalg(GEN nf, GEN x)
1100 : {
1101 : GEN T;
1102 :
1103 358712 : nf = checknf(nf);
1104 358712 : switch(typ(x))
1105 : {
1106 299355 : case t_COL: {
1107 299355 : pari_sp av = avma;
1108 299355 : return gerepilecopy(av, coltoalg(nf, x));
1109 : }
1110 33390 : case t_POLMOD:
1111 33390 : T = nf_get_pol(nf);
1112 33390 : if (!RgX_equal_var(T,gel(x,1)))
1113 0 : pari_err_MODULUS("basistoalg", T,gel(x,1));
1114 33390 : return gcopy(x);
1115 2359 : case t_POL:
1116 2359 : T = nf_get_pol(nf);
1117 2359 : if (varn(T) != varn(x)) pari_err_VAR("basistoalg",x,T);
1118 2359 : retmkpolmod(RgX_rem(x, T), ZX_copy(T));
1119 23608 : case t_INT:
1120 : case t_FRAC:
1121 23608 : T = nf_get_pol(nf);
1122 23608 : retmkpolmod(gcopy(x), ZX_copy(T));
1123 0 : default:
1124 0 : pari_err_TYPE("basistoalg",x);
1125 : return NULL; /* LCOV_EXCL_LINE */
1126 : }
1127 : }
1128 :
1129 : /* true nf, x a t_POL */
1130 : static GEN
1131 4543382 : pol_to_scalar_or_basis(GEN nf, GEN x)
1132 : {
1133 4543382 : GEN T = nf_get_pol(nf);
1134 4543382 : long l = lg(x);
1135 4543382 : if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_basis", x,T);
1136 4543277 : if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
1137 4543277 : if (l == 2) return gen_0;
1138 3539190 : if (l == 3)
1139 : {
1140 816312 : x = gel(x,2);
1141 816312 : if (!is_rational_t(typ(x))) pari_err_TYPE("nf_to_scalar_or_basis",x);
1142 816305 : return x;
1143 : }
1144 2722878 : return poltobasis(nf,x);
1145 : }
1146 : /* Assume nf is a genuine nf. */
1147 : GEN
1148 161334673 : nf_to_scalar_or_basis(GEN nf, GEN x)
1149 : {
1150 161334673 : switch(typ(x))
1151 : {
1152 96904560 : case t_INT: case t_FRAC:
1153 96904560 : return x;
1154 556274 : case t_POLMOD:
1155 556274 : x = checknfelt_mod(nf,x,"nf_to_scalar_or_basis");
1156 556138 : switch(typ(x))
1157 : {
1158 85428 : case t_INT: case t_FRAC: return x;
1159 470710 : case t_POL: return pol_to_scalar_or_basis(nf,x);
1160 : }
1161 0 : break;
1162 4072674 : case t_POL: return pol_to_scalar_or_basis(nf,x);
1163 59805347 : case t_COL:
1164 59805347 : if (lg(x)-1 != nf_get_degree(nf)) break;
1165 59805051 : return QV_isscalar(x)? gel(x,1): x;
1166 : }
1167 96 : pari_err_TYPE("nf_to_scalar_or_basis",x);
1168 : return NULL; /* LCOV_EXCL_LINE */
1169 : }
1170 : /* Let x be a polynomial with coefficients in Q or nf. Return the same
1171 : * polynomial with coefficients expressed as vectors (on the integral basis).
1172 : * No consistency checks, not memory-clean. */
1173 : GEN
1174 28737 : RgX_to_nfX(GEN nf, GEN x)
1175 : {
1176 : long i, l;
1177 28737 : GEN y = cgetg_copy(x, &l); y[1] = x[1];
1178 235369 : for (i=2; i<l; i++) gel(y,i) = nf_to_scalar_or_basis(nf, gel(x,i));
1179 28737 : return y;
1180 : }
1181 :
1182 : /* Assume nf is a genuine nf. */
1183 : GEN
1184 3881534 : nf_to_scalar_or_alg(GEN nf, GEN x)
1185 : {
1186 3881534 : switch(typ(x))
1187 : {
1188 84795 : case t_INT: case t_FRAC:
1189 84795 : return x;
1190 420 : case t_POLMOD:
1191 420 : x = checknfelt_mod(nf,x,"nf_to_scalar_or_alg");
1192 420 : if (typ(x) != t_POL) return x;
1193 : /* fall through */
1194 : case t_POL:
1195 : {
1196 5124 : GEN T = nf_get_pol(nf);
1197 5124 : long l = lg(x);
1198 5124 : if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_alg", x,T);
1199 5124 : if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
1200 5124 : if (l == 2) return gen_0;
1201 5124 : if (l == 3) return gel(x,2);
1202 3612 : return x;
1203 : }
1204 3791568 : case t_COL:
1205 : {
1206 : GEN dx;
1207 3791568 : if (lg(x)-1 != nf_get_degree(nf)) break;
1208 7506334 : if (QV_isscalar(x)) return gel(x,1);
1209 3714574 : x = Q_remove_denom(x, &dx);
1210 3714617 : x = RgV_RgC_mul(nf_get_zkprimpart(nf), x);
1211 3714763 : dx = mul_denom(dx, nf_get_zkden(nf));
1212 3714747 : return gdiv(x,dx);
1213 : }
1214 : }
1215 48 : pari_err_TYPE("nf_to_scalar_or_alg",x);
1216 : return NULL; /* LCOV_EXCL_LINE */
1217 : }
1218 :
1219 : /* gmul(A, RgX_to_RgC(x)), A t_MAT of compatible dimensions */
1220 : GEN
1221 1365 : RgM_RgX_mul(GEN A, GEN x)
1222 : {
1223 1365 : long i, l = lg(x)-1;
1224 : GEN z;
1225 1365 : if (l == 1) return zerocol(nbrows(A));
1226 1351 : z = gmul(gel(x,2), gel(A,1));
1227 2555 : for (i = 2; i < l; i++)
1228 1204 : if (!gequal0(gel(x,i+1))) z = gadd(z, gmul(gel(x,i+1), gel(A,i)));
1229 1351 : return z;
1230 : }
1231 : GEN
1232 10315968 : ZM_ZX_mul(GEN A, GEN x)
1233 : {
1234 10315968 : long i, l = lg(x)-1;
1235 : GEN z;
1236 10315968 : if (l == 1) return zerocol(nbrows(A));
1237 10314834 : z = ZC_Z_mul(gel(A,1), gel(x,2));
1238 32201782 : for (i = 2; i < l ; i++)
1239 21889392 : if (signe(gel(x,i+1))) z = ZC_add(z, ZC_Z_mul(gel(A,i), gel(x,i+1)));
1240 10312390 : return z;
1241 : }
1242 : /* x a t_POL, nf a genuine nf. No garbage collecting. No check. */
1243 : GEN
1244 9716526 : poltobasis(GEN nf, GEN x)
1245 : {
1246 9716526 : GEN d, T = nf_get_pol(nf);
1247 9716487 : if (varn(x) != varn(T)) pari_err_VAR( "poltobasis", x,T);
1248 9716354 : if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
1249 9716352 : x = Q_remove_denom(x, &d);
1250 9716543 : if (!RgX_is_ZX(x)) pari_err_TYPE("poltobasis",x);
1251 9716457 : x = ZM_ZX_mul(nf_get_invzk(nf), x);
1252 9714525 : if (d) x = RgC_Rg_div(x, d);
1253 9714566 : return x;
1254 : }
1255 :
1256 : GEN
1257 921437 : algtobasis(GEN nf, GEN x)
1258 : {
1259 : pari_sp av;
1260 :
1261 921437 : nf = checknf(nf);
1262 921437 : switch(typ(x))
1263 : {
1264 113254 : case t_POLMOD:
1265 113254 : if (!RgX_equal_var(nf_get_pol(nf),gel(x,1)))
1266 7 : pari_err_MODULUS("algtobasis", nf_get_pol(nf),gel(x,1));
1267 113247 : x = gel(x,2);
1268 113247 : switch(typ(x))
1269 : {
1270 8197 : case t_INT:
1271 8197 : case t_FRAC: return scalarcol(x, nf_get_degree(nf));
1272 105050 : case t_POL:
1273 105050 : av = avma;
1274 105050 : return gerepileupto(av,poltobasis(nf,x));
1275 : }
1276 0 : break;
1277 :
1278 249760 : case t_POL:
1279 249760 : av = avma;
1280 249760 : return gerepileupto(av,poltobasis(nf,x));
1281 :
1282 83178 : case t_COL:
1283 83178 : if (!RgV_is_QV(x)) pari_err_TYPE("nfalgtobasis",x);
1284 83171 : if (lg(x)-1 != nf_get_degree(nf)) pari_err_DIM("nfalgtobasis");
1285 83171 : return gcopy(x);
1286 :
1287 475249 : case t_INT:
1288 475249 : case t_FRAC: return scalarcol(x, nf_get_degree(nf));
1289 : }
1290 0 : pari_err_TYPE("algtobasis",x);
1291 : return NULL; /* LCOV_EXCL_LINE */
1292 : }
1293 :
1294 : GEN
1295 47488 : rnfbasistoalg(GEN rnf,GEN x)
1296 : {
1297 47488 : const char *f = "rnfbasistoalg";
1298 : long lx, i;
1299 47488 : pari_sp av = avma;
1300 : GEN z, nf, R, T;
1301 :
1302 47488 : checkrnf(rnf);
1303 47488 : nf = rnf_get_nf(rnf);
1304 47488 : T = nf_get_pol(nf);
1305 47488 : R = QXQX_to_mod_shallow(rnf_get_pol(rnf), T);
1306 47488 : switch(typ(x))
1307 : {
1308 875 : case t_COL:
1309 875 : z = cgetg_copy(x, &lx);
1310 2597 : for (i=1; i<lx; i++)
1311 : {
1312 1778 : GEN c = nf_to_scalar_or_alg(nf, gel(x,i));
1313 1722 : if (typ(c) == t_POL) c = mkpolmod(c,T);
1314 1722 : gel(z,i) = c;
1315 : }
1316 819 : z = RgV_RgC_mul(gel(rnf_get_zk(rnf),1), z);
1317 735 : return gerepileupto(av, gmodulo(z,R));
1318 :
1319 31227 : case t_POLMOD:
1320 31227 : x = polmod_nffix(f, rnf, x, 0);
1321 30954 : if (typ(x) != t_POL) break;
1322 14261 : retmkpolmod(RgX_copy(x), RgX_copy(R));
1323 1274 : case t_POL:
1324 1274 : if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = gmodulo(x,T); break; }
1325 1029 : if (varn(x) == varn(R))
1326 : {
1327 973 : x = RgX_nffix(f,nf_get_pol(nf),x,0);
1328 973 : return gmodulo(x, R);
1329 : }
1330 56 : pari_err_VAR(f, x,R);
1331 : }
1332 30994 : retmkpolmod(scalarpol(x, varn(R)), RgX_copy(R));
1333 : }
1334 :
1335 : GEN
1336 2275 : matbasistoalg(GEN nf,GEN x)
1337 : {
1338 : long i, j, li, lx;
1339 2275 : GEN z = cgetg_copy(x, &lx);
1340 :
1341 2275 : if (lx == 1) return z;
1342 2268 : switch(typ(x))
1343 : {
1344 77 : case t_VEC: case t_COL:
1345 273 : for (i=1; i<lx; i++) gel(z,i) = basistoalg(nf, gel(x,i));
1346 77 : return z;
1347 2191 : case t_MAT: break;
1348 0 : default: pari_err_TYPE("matbasistoalg",x);
1349 : }
1350 2191 : li = lgcols(x);
1351 8183 : for (j=1; j<lx; j++)
1352 : {
1353 5992 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1354 5992 : gel(z,j) = c;
1355 28077 : for (i=1; i<li; i++) gel(c,i) = basistoalg(nf, gel(xj,i));
1356 : }
1357 2191 : return z;
1358 : }
1359 :
1360 : GEN
1361 30736 : matalgtobasis(GEN nf,GEN x)
1362 : {
1363 : long i, j, li, lx;
1364 30736 : GEN z = cgetg_copy(x, &lx);
1365 :
1366 30736 : if (lx == 1) return z;
1367 30274 : switch(typ(x))
1368 : {
1369 30267 : case t_VEC: case t_COL:
1370 79584 : for (i=1; i<lx; i++) gel(z,i) = algtobasis(nf, gel(x,i));
1371 30267 : return z;
1372 7 : case t_MAT: break;
1373 0 : default: pari_err_TYPE("matalgtobasis",x);
1374 : }
1375 7 : li = lgcols(x);
1376 14 : for (j=1; j<lx; j++)
1377 : {
1378 7 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1379 7 : gel(z,j) = c;
1380 21 : for (i=1; i<li; i++) gel(c,i) = algtobasis(nf, gel(xj,i));
1381 : }
1382 7 : return z;
1383 : }
1384 : GEN
1385 10932 : RgM_to_nfM(GEN nf,GEN x)
1386 : {
1387 : long i, j, li, lx;
1388 10932 : GEN z = cgetg_copy(x, &lx);
1389 :
1390 10932 : if (lx == 1) return z;
1391 10932 : li = lgcols(x);
1392 81564 : for (j=1; j<lx; j++)
1393 : {
1394 70632 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1395 70632 : gel(z,j) = c;
1396 462195 : for (i=1; i<li; i++) gel(c,i) = nf_to_scalar_or_basis(nf, gel(xj,i));
1397 : }
1398 10932 : return z;
1399 : }
1400 : GEN
1401 148496 : RgC_to_nfC(GEN nf, GEN x)
1402 908609 : { pari_APPLY_type(t_COL, nf_to_scalar_or_basis(nf, gel(x,i))) }
1403 :
1404 : /* x a t_POLMOD, supposedly in rnf = K[z]/(T), K = Q[y]/(Tnf) */
1405 : GEN
1406 141436 : polmod_nffix(const char *f, GEN rnf, GEN x, int lift)
1407 141436 : { return polmod_nffix2(f, rnf_get_nfpol(rnf), rnf_get_pol(rnf), x,lift); }
1408 : GEN
1409 141527 : polmod_nffix2(const char *f, GEN T, GEN R, GEN x, int lift)
1410 : {
1411 141527 : if (RgX_equal_var(gel(x,1), R))
1412 : {
1413 129087 : x = gel(x,2);
1414 129087 : if (typ(x) == t_POL && varn(x) == varn(R))
1415 : {
1416 98349 : x = RgX_nffix(f, T, x, lift);
1417 98349 : switch(lg(x))
1418 : {
1419 5817 : case 2: return gen_0;
1420 12197 : case 3: return gel(x,2);
1421 : }
1422 80335 : return x;
1423 : }
1424 : }
1425 43178 : return Rg_nffix(f, T, x, lift);
1426 : }
1427 : GEN
1428 1428 : rnfalgtobasis(GEN rnf,GEN x)
1429 : {
1430 1428 : const char *f = "rnfalgtobasis";
1431 1428 : pari_sp av = avma;
1432 : GEN T, R;
1433 :
1434 1428 : checkrnf(rnf);
1435 1428 : R = rnf_get_pol(rnf);
1436 1428 : T = rnf_get_nfpol(rnf);
1437 1428 : switch(typ(x))
1438 : {
1439 98 : case t_COL:
1440 98 : if (lg(x)-1 != rnf_get_degree(rnf)) pari_err_DIM(f);
1441 49 : x = RgV_nffix(f, T, x, 0);
1442 42 : return gerepilecopy(av, x);
1443 :
1444 1162 : case t_POLMOD:
1445 1162 : x = polmod_nffix(f, rnf, x, 0);
1446 1057 : if (typ(x) != t_POL) break;
1447 714 : return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
1448 112 : case t_POL:
1449 112 : if (varn(x) == varn(T))
1450 : {
1451 42 : RgX_check_QX(x,f);
1452 28 : if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
1453 28 : x = mkpolmod(x,T); break;
1454 : }
1455 70 : x = RgX_nffix(f, T, x, 0);
1456 56 : if (degpol(x) >= degpol(R)) x = RgX_rem(x, R);
1457 56 : return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
1458 : }
1459 427 : return gerepileupto(av, scalarcol(x, rnf_get_degree(rnf)));
1460 : }
1461 :
1462 : /* Given a and b in nf, gives an algebraic integer y in nf such that a-b.y
1463 : * is "small" */
1464 : GEN
1465 259 : nfdiveuc(GEN nf, GEN a, GEN b)
1466 : {
1467 259 : pari_sp av = avma;
1468 259 : a = nfdiv(nf,a,b);
1469 259 : return gerepileupto(av, ground(a));
1470 : }
1471 :
1472 : /* Given a and b in nf, gives a "small" algebraic integer r in nf
1473 : * of the form a-b.y */
1474 : GEN
1475 259 : nfmod(GEN nf, GEN a, GEN b)
1476 : {
1477 259 : pari_sp av = avma;
1478 259 : GEN p1 = gneg_i(nfmul(nf,b,ground(nfdiv(nf,a,b))));
1479 259 : return gerepileupto(av, nfadd(nf,a,p1));
1480 : }
1481 :
1482 : /* Given a and b in nf, gives a two-component vector [y,r] in nf such
1483 : * that r=a-b.y is "small". */
1484 : GEN
1485 259 : nfdivrem(GEN nf, GEN a, GEN b)
1486 : {
1487 259 : pari_sp av = avma;
1488 259 : GEN p1,z, y = ground(nfdiv(nf,a,b));
1489 :
1490 259 : p1 = gneg_i(nfmul(nf,b,y));
1491 259 : z = cgetg(3,t_VEC);
1492 259 : gel(z,1) = gcopy(y);
1493 259 : gel(z,2) = nfadd(nf,a,p1); return gerepileupto(av, z);
1494 : }
1495 :
1496 : /*************************************************************************/
1497 : /** **/
1498 : /** LOGARITHMIC EMBEDDINGS **/
1499 : /** **/
1500 : /*************************************************************************/
1501 :
1502 : static int
1503 4611246 : low_prec(GEN x)
1504 : {
1505 4611246 : switch(typ(x))
1506 : {
1507 0 : case t_INT: return !signe(x);
1508 4611246 : case t_REAL: return !signe(x) || realprec(x) <= DEFAULTPREC;
1509 0 : default: return 0;
1510 : }
1511 : }
1512 :
1513 : static GEN
1514 23102 : cxlog_1(GEN nf) { return zerocol(lg(nf_get_roots(nf))-1); }
1515 : static GEN
1516 545 : cxlog_m1(GEN nf, long prec)
1517 : {
1518 545 : long i, l = lg(nf_get_roots(nf)), r1 = nf_get_r1(nf);
1519 545 : GEN v = cgetg(l, t_COL), p, P;
1520 545 : p = mppi(prec); P = mkcomplex(gen_0, p);
1521 1224 : for (i = 1; i <= r1; i++) gel(v,i) = P; /* IPi*/
1522 545 : if (i < l) P = gmul2n(P,1);
1523 1160 : for ( ; i < l; i++) gel(v,i) = P; /* 2IPi */
1524 545 : return v;
1525 : }
1526 : static GEN
1527 1714781 : ZC_cxlog(GEN nf, GEN x, long prec)
1528 : {
1529 : long i, l, r1;
1530 : GEN v;
1531 1714781 : x = RgM_RgC_mul(nf_get_M(nf), Q_primpart(x));
1532 1714782 : l = lg(x); r1 = nf_get_r1(nf);
1533 4330070 : for (i = 1; i <= r1; i++)
1534 2615288 : if (low_prec(gel(x,i))) return NULL;
1535 3513920 : for ( ; i < l; i++)
1536 1799140 : if (low_prec(gnorm(gel(x,i)))) return NULL;
1537 1714780 : v = cgetg(l,t_COL);
1538 4330068 : for (i = 1; i <= r1; i++) gel(v,i) = glog(gel(x,i),prec);
1539 3513919 : for ( ; i < l; i++) gel(v,i) = gmul2n(glog(gel(x,i),prec),1);
1540 1714781 : return v;
1541 : }
1542 : static GEN
1543 223273 : famat_cxlog(GEN nf, GEN fa, long prec)
1544 : {
1545 223273 : GEN G, E, y = NULL;
1546 : long i, l;
1547 :
1548 223273 : if (typ(fa) != t_MAT) pari_err_TYPE("famat_cxlog",fa);
1549 223273 : if (lg(fa) == 1) return cxlog_1(nf);
1550 223273 : G = gel(fa,1);
1551 223273 : E = gel(fa,2); l = lg(E);
1552 1119327 : for (i = 1; i < l; i++)
1553 : {
1554 896054 : GEN t, e = gel(E,i), x = nf_to_scalar_or_basis(nf, gel(G,i));
1555 : /* multiplicative arch would be better (save logs), but exponents overflow
1556 : * [ could keep track of expo separately, but not worth it ] */
1557 896054 : switch(typ(x))
1558 : { /* ignore positive rationals */
1559 16415 : case t_FRAC: x = gel(x,1); /* fall through */
1560 266282 : case t_INT: if (signe(x) > 0) continue;
1561 97 : if (!mpodd(e)) continue;
1562 41 : t = cxlog_m1(nf, prec); /* we probably should not reach this line */
1563 41 : break;
1564 629772 : default: /* t_COL */
1565 629772 : t = ZC_cxlog(nf,x,prec); if (!t) return NULL;
1566 629772 : t = RgC_Rg_mul(t, e);
1567 : }
1568 629813 : y = y? RgV_add(y,t): t;
1569 : }
1570 223273 : return y ? y: cxlog_1(nf);
1571 : }
1572 : /* Archimedean components: [e_i Log( sigma_i(X) )], where X = primpart(x),
1573 : * and e_i = 1 (resp 2.) for i <= R1 (resp. > R1) */
1574 : GEN
1575 1309430 : nf_cxlog(GEN nf, GEN x, long prec)
1576 : {
1577 1309430 : if (typ(x) == t_MAT) return famat_cxlog(nf,x,prec);
1578 1086157 : x = nf_to_scalar_or_basis(nf,x);
1579 1086157 : switch(typ(x))
1580 : {
1581 0 : case t_FRAC: x = gel(x,1); /* fall through */
1582 1148 : case t_INT:
1583 1148 : return signe(x) > 0? cxlog_1(nf): cxlog_m1(nf, prec);
1584 1085009 : default:
1585 1085009 : return ZC_cxlog(nf, x, prec);
1586 : }
1587 : }
1588 : GEN
1589 97 : nfV_cxlog(GEN nf, GEN x, long prec)
1590 : {
1591 : long i, l;
1592 97 : GEN v = cgetg_copy(x, &l);
1593 167 : for (i = 1; i < l; i++)
1594 70 : if (!(gel(v,i) = nf_cxlog(nf, gel(x,i), prec))) return NULL;
1595 97 : return v;
1596 : }
1597 :
1598 : static GEN
1599 15232 : scalar_logembed(GEN nf, GEN u, GEN *emb)
1600 : {
1601 : GEN v, logu;
1602 15232 : long i, s = signe(u), RU = lg(nf_get_roots(nf))-1, R1 = nf_get_r1(nf);
1603 :
1604 15232 : if (!s) pari_err_DOMAIN("nflogembed","argument","=",gen_0,u);
1605 15232 : v = cgetg(RU+1, t_COL); logu = logr_abs(u);
1606 17213 : for (i = 1; i <= R1; i++) gel(v,i) = logu;
1607 15232 : if (i <= RU)
1608 : {
1609 14350 : GEN logu2 = shiftr(logu,1);
1610 55839 : for ( ; i <= RU; i++) gel(v,i) = logu2;
1611 : }
1612 15232 : if (emb) *emb = const_col(RU, u);
1613 15232 : return v;
1614 : }
1615 :
1616 : static GEN
1617 1309 : famat_logembed(GEN nf,GEN x,GEN *emb,long prec)
1618 : {
1619 1309 : GEN A, M, T, a, t, g = gel(x,1), e = gel(x,2);
1620 1309 : long i, l = lg(e);
1621 :
1622 1309 : if (l == 1) return scalar_logembed(nf, real_1(prec), emb);
1623 1309 : A = NULL; T = emb? cgetg(l, t_COL): NULL;
1624 1309 : if (emb) *emb = M = mkmat2(T, e);
1625 62132 : for (i = 1; i < l; i++)
1626 : {
1627 60823 : a = nflogembed(nf, gel(g,i), &t, prec);
1628 60823 : if (!a) return NULL;
1629 60823 : a = RgC_Rg_mul(a, gel(e,i));
1630 60823 : A = A? RgC_add(A, a): a;
1631 60823 : if (emb) gel(T,i) = t;
1632 : }
1633 1309 : return A;
1634 : }
1635 :
1636 : /* Get archimedean components: [e_i log( | sigma_i(x) | )], with e_i = 1
1637 : * (resp 2.) for i <= R1 (resp. > R1) and set emb to the embeddings of x.
1638 : * Return NULL if precision problem */
1639 : GEN
1640 98658 : nflogembed(GEN nf, GEN x, GEN *emb, long prec)
1641 : {
1642 : long i, l, r1;
1643 : GEN v, t;
1644 :
1645 98658 : if (typ(x) == t_MAT) return famat_logembed(nf,x,emb,prec);
1646 97349 : x = nf_to_scalar_or_basis(nf,x);
1647 97349 : if (typ(x) != t_COL) return scalar_logembed(nf, gtofp(x,prec), emb);
1648 82117 : x = RgM_RgC_mul(nf_get_M(nf), x);
1649 82117 : l = lg(x); r1 = nf_get_r1(nf); v = cgetg(l,t_COL);
1650 109046 : for (i = 1; i <= r1; i++)
1651 : {
1652 26929 : t = gabs(gel(x,i),prec); if (low_prec(t)) return NULL;
1653 26929 : gel(v,i) = glog(t,prec);
1654 : }
1655 252006 : for ( ; i < l; i++)
1656 : {
1657 169890 : t = gnorm(gel(x,i)); if (low_prec(t)) return NULL;
1658 169890 : gel(v,i) = glog(t,prec);
1659 : }
1660 82116 : if (emb) *emb = x;
1661 82116 : return v;
1662 : }
1663 :
1664 : /*************************************************************************/
1665 : /** **/
1666 : /** REAL EMBEDDINGS **/
1667 : /** **/
1668 : /*************************************************************************/
1669 : static GEN
1670 486240 : sarch_get_cyc(GEN sarch) { return gel(sarch,1); }
1671 : static GEN
1672 665696 : sarch_get_archp(GEN sarch) { return gel(sarch,2); }
1673 : static GEN
1674 163632 : sarch_get_MI(GEN sarch) { return gel(sarch,3); }
1675 : static GEN
1676 163632 : sarch_get_lambda(GEN sarch) { return gel(sarch,4); }
1677 : static GEN
1678 163632 : sarch_get_F(GEN sarch) { return gel(sarch,5); }
1679 :
1680 : /* x not a scalar, true nf, return number of positive roots of char_x */
1681 : static long
1682 1290 : num_positive(GEN nf, GEN x)
1683 : {
1684 1290 : GEN T = nf_get_pol(nf), B, charx;
1685 : long dnf, vnf, N;
1686 1290 : x = nf_to_scalar_or_alg(nf, x); /* not a scalar */
1687 1290 : charx = ZXQ_charpoly(x, T, 0);
1688 1290 : charx = ZX_radical(charx);
1689 1290 : N = degpol(T) / degpol(charx);
1690 : /* real places are unramified ? */
1691 1290 : if (N == 1 || ZX_sturm(charx) * N == nf_get_r1(nf))
1692 1283 : return ZX_sturmpart(charx, mkvec2(gen_0,mkoo())) * N;
1693 : /* painful case, multiply by random square until primitive */
1694 7 : dnf = nf_get_degree(nf);
1695 7 : vnf = varn(T);
1696 7 : B = int2n(10);
1697 : for(;;)
1698 0 : {
1699 7 : GEN y = RgXQ_sqr(random_FpX(dnf, vnf, B), T);
1700 7 : y = RgXQ_mul(x, y, T);
1701 7 : charx = ZXQ_charpoly(y, T, 0);
1702 7 : if (ZX_is_squarefree(charx))
1703 7 : return ZX_sturmpart(charx, mkvec2(gen_0,mkoo())) * N;
1704 : }
1705 : }
1706 :
1707 : /* x a QC: return sigma_k(x) where 1 <= k <= r1+r2; correct but inefficient
1708 : * if x in Q. M = nf_get_M(nf) */
1709 : static GEN
1710 2134 : nfembed_i(GEN M, GEN x, long k)
1711 : {
1712 2134 : long i, l = lg(M);
1713 2134 : GEN z = gel(x,1);
1714 24356 : for (i = 2; i < l; i++) z = gadd(z, gmul(gcoeff(M,k,i), gel(x,i)));
1715 2134 : return z;
1716 : }
1717 : GEN
1718 0 : nfembed(GEN nf, GEN x, long k)
1719 : {
1720 0 : pari_sp av = avma;
1721 0 : nf = checknf(nf);
1722 0 : x = nf_to_scalar_or_basis(nf,x);
1723 0 : if (typ(x) != t_COL) return gerepilecopy(av, x);
1724 0 : return gerepileupto(av, nfembed_i(nf_get_M(nf),x,k));
1725 : }
1726 :
1727 : /* x a ZC */
1728 : static GEN
1729 905988 : zk_embed(GEN M, GEN x, long k)
1730 : {
1731 905988 : long i, l = lg(x);
1732 905988 : GEN z = gel(x,1); /* times M[k,1], which is 1 */
1733 3009444 : for (i = 2; i < l; i++) z = mpadd(z, mpmul(gcoeff(M,k,i), gel(x,i)));
1734 905982 : return z;
1735 : }
1736 :
1737 : /* Given floating point approximation z of sigma_k(x), decide its sign
1738 : * [0/+, 1/- and -1 for FAIL] */
1739 : static long
1740 887641 : eval_sign_embed(GEN z)
1741 : {
1742 887641 : if (typ(z) == t_REAL)
1743 : {
1744 887641 : long l = realprec(z);
1745 887641 : if (l <= LOWDEFAULTPREC
1746 887641 : || (l == LOWDEFAULTPREC + 1 && !z[l-1])) return -1; /* dubious, fail */
1747 886884 : if (expo(z) < 16 - l) return -1; /* same */
1748 : }
1749 886841 : return (signe(z) < 1)? 1: 0;
1750 : }
1751 : /* return v such that (-1)^v = sign(sigma_k(x)), x primitive ZC */
1752 : static long
1753 791125 : eval_sign(GEN M, GEN x, long k)
1754 791125 : { return eval_sign_embed( zk_embed(M, x, k) ); }
1755 :
1756 : /* check that signs[i..#signs] == s; signs = NULL encodes "totally positive" */
1757 : static int
1758 0 : oksigns(long l, GEN signs, long i, long s)
1759 : {
1760 0 : if (!signs) return s == 0;
1761 0 : for (; i < l; i++)
1762 0 : if (signs[i] != s) return 0;
1763 0 : return 1;
1764 : }
1765 : /* check that signs[i] = s and signs[i+1..#signs] = 1-s */
1766 : static int
1767 0 : oksigns2(long l, GEN signs, long i, long s)
1768 : {
1769 0 : if (!signs) return s == 0 && i == l-1;
1770 0 : return signs[i] == s && oksigns(l, signs, i+1, 1-s);
1771 : }
1772 :
1773 : /* true nf, x a ZC (primitive for efficiency) which is not a scalar; embx its
1774 : * embeddings or NULL */
1775 : static int
1776 80246 : nfchecksigns_i(GEN nf, GEN x, GEN embx, GEN signs, GEN archp)
1777 : {
1778 80246 : long i, l = lg(archp), np = -1;
1779 80246 : long bigx = embx? 0: gexpo(x) >= nf_get_prec(nf);
1780 80246 : GEN M = nf_get_M(nf), sarch = NULL;
1781 126373 : for (i = 1; i < l; i++)
1782 : {
1783 97892 : long s = -1;
1784 97892 : if (embx)
1785 96528 : s = eval_sign_embed(gel(embx,i));
1786 1364 : else if (!bigx)
1787 1364 : s = eval_sign(M, x, archp[i]);
1788 : /* 0 / + or 1 / -; -1 for FAIL */
1789 97892 : if (s < 0) /* failure */
1790 : {
1791 0 : long ni, r1 = nf_get_r1(nf);
1792 : GEN xi;
1793 0 : if (np < 0)
1794 : {
1795 0 : np = num_positive(nf, x);
1796 0 : if (np == 0) return oksigns(l, signs, i, 1);
1797 0 : if (np == r1) return oksigns(l, signs, i, 0);
1798 0 : sarch = nfarchstar(nf, NULL, identity_perm(r1));
1799 : }
1800 0 : xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
1801 0 : xi = Q_primpart(xi);
1802 0 : ni = num_positive(nf, nfmuli(nf,x,xi));
1803 0 : if (ni == 0) return oksigns2(l, signs, i, 0);
1804 0 : if (ni == r1) return oksigns2(l, signs, i, 1);
1805 0 : s = ni < np? 0: 1;
1806 : }
1807 97892 : if (s != (signs? signs[i]: 0)) return 0;
1808 : }
1809 28481 : return 1;
1810 : }
1811 : static void
1812 775 : pl_convert(GEN pl, GEN *psigns, GEN *parchp)
1813 : {
1814 775 : long i, j, l = lg(pl);
1815 775 : GEN signs = cgetg(l, t_VECSMALL);
1816 775 : GEN archp = cgetg(l, t_VECSMALL);
1817 2576 : for (i = j = 1; i < l; i++)
1818 : {
1819 1801 : if (!pl[i]) continue;
1820 1403 : archp[j] = i;
1821 1403 : signs[j] = (pl[i] < 0)? 1: 0;
1822 1403 : j++;
1823 : }
1824 775 : setlg(archp, j); *parchp = archp;
1825 775 : setlg(signs, j); *psigns = signs;
1826 775 : }
1827 : /* pl : requested signs for real embeddings, 0 = no sign constraint */
1828 : int
1829 14719 : nfchecksigns(GEN nf, GEN x, GEN pl)
1830 : {
1831 14719 : pari_sp av = avma;
1832 : GEN signs, archp;
1833 14719 : nf = checknf(nf);
1834 14719 : x = nf_to_scalar_or_basis(nf,x);
1835 14719 : if (typ(x) != t_COL)
1836 : {
1837 13944 : long i, l = lg(pl), s = gsigne(x);
1838 27853 : for (i = 1; i < l; i++)
1839 13909 : if (pl[i] && pl[i] != s) return gc_bool(av,0);
1840 13944 : return gc_bool(av,1);
1841 : }
1842 775 : pl_convert(pl, &signs, &archp);
1843 775 : return gc_bool(av, nfchecksigns_i(nf, x, NULL, signs, archp));
1844 : }
1845 :
1846 : /* signs = NULL: totally positive, else sign[i] = 0 (+) or 1 (-) */
1847 : static GEN
1848 163632 : get_C(GEN lambda, long l, GEN signs)
1849 : {
1850 : long i;
1851 : GEN C, mlambda;
1852 163632 : if (!signs) return const_vec(l-1, lambda);
1853 133882 : C = cgetg(l, t_COL); mlambda = gneg(lambda);
1854 342945 : for (i = 1; i < l; i++) gel(C,i) = signs[i]? mlambda: lambda;
1855 133885 : return C;
1856 : }
1857 : /* signs = NULL: totally positive at archp.
1858 : * Assume that a t_COL x is not a scalar */
1859 : static GEN
1860 277347 : nfsetsigns(GEN nf, GEN signs, GEN x, GEN sarch)
1861 : {
1862 277347 : long i, l = lg(sarch_get_archp(sarch));
1863 : GEN ex;
1864 : /* Is signature already correct ? */
1865 277346 : if (typ(x) != t_COL)
1866 : {
1867 197882 : long s = gsigne(x);
1868 197882 : if (!s) i = 1;
1869 197861 : else if (!signs)
1870 7427 : i = (s < 0)? 1: l;
1871 : else
1872 : {
1873 190434 : s = s < 0? 1: 0;
1874 324137 : for (i = 1; i < l; i++)
1875 245441 : if (signs[i] != s) break;
1876 : }
1877 197882 : ex = (i < l)? const_col(l-1, x): NULL;
1878 : }
1879 : else
1880 : { /* inefficient if x scalar, wrong if x = 0 */
1881 79464 : pari_sp av = avma;
1882 79464 : GEN cex, M = nf_get_M(nf), archp = sarch_get_archp(sarch);
1883 79468 : GEN xp = Q_primitive_part(x,&cex);
1884 79468 : ex = cgetg(l,t_COL);
1885 194332 : for (i = 1; i < l; i++) gel(ex,i) = zk_embed(M,xp,archp[i]);
1886 79471 : if (nfchecksigns_i(nf, xp, ex, signs, archp)) { ex = NULL; set_avma(av); }
1887 51728 : else if (cex) ex = RgC_Rg_mul(ex, cex); /* put back content */
1888 : }
1889 277353 : if (ex)
1890 : { /* If no, fix it */
1891 163632 : GEN MI = sarch_get_MI(sarch), F = sarch_get_F(sarch);
1892 163632 : GEN lambda = sarch_get_lambda(sarch);
1893 163632 : GEN t = RgC_sub(get_C(lambda, l, signs), ex);
1894 163625 : t = grndtoi(RgM_RgC_mul(MI,t), NULL);
1895 163617 : if (lg(F) != 1) t = ZM_ZC_mul(F, t);
1896 163623 : x = typ(x) == t_COL? RgC_add(t, x): RgC_Rg_add(t, x);
1897 : }
1898 277341 : return x;
1899 : }
1900 : /* - true nf
1901 : * - sarch = nfarchstar(nf, F);
1902 : * - x encodes a vector of signs at arch.archp: either a t_VECSMALL
1903 : * (vector of signs as {0,1}-vector), NULL (totally positive at archp),
1904 : * or a nonzero number field element (replaced by its signature at archp);
1905 : * - y is a nonzero number field element
1906 : * Return z = y (mod F) with signs(y, archp) = signs(x) (a {0,1}-vector).
1907 : * Not stack-clean */
1908 : GEN
1909 308890 : set_sign_mod_divisor(GEN nf, GEN x, GEN y, GEN sarch)
1910 : {
1911 308890 : GEN archp = sarch_get_archp(sarch);
1912 308888 : if (lg(archp) == 1) return y;
1913 275582 : if (x && typ(x) != t_VECSMALL) x = nfsign_arch(nf, x, archp);
1914 275582 : return nfsetsigns(nf, x, nf_to_scalar_or_basis(nf,y), sarch);
1915 : }
1916 :
1917 : static GEN
1918 83452 : setsigns_init(GEN nf, GEN archp, GEN F, GEN DATA)
1919 : {
1920 83452 : GEN lambda, Mr = rowpermute(nf_get_M(nf), archp), MI = F? RgM_mul(Mr,F): Mr;
1921 83457 : lambda = gmul2n(matrixnorm(MI,DEFAULTPREC), -1);
1922 83458 : if (typ(lambda) != t_REAL) lambda = gmul(lambda, uutoQ(1001,1000));
1923 83458 : if (lg(archp) < lg(MI))
1924 : {
1925 58924 : GEN perm = gel(indexrank(MI), 2);
1926 58926 : if (!F) F = matid(nf_get_degree(nf));
1927 58926 : MI = vecpermute(MI, perm);
1928 58925 : F = vecpermute(F, perm);
1929 : }
1930 83459 : if (!F) F = cgetg(1,t_MAT);
1931 83459 : MI = RgM_inv(MI);
1932 83459 : return mkvec5(DATA, archp, MI, lambda, F);
1933 : }
1934 : /* F nonzero integral ideal in HNF (or NULL: Z_K), compute elements in 1+F
1935 : * whose sign matrix at archp is identity; archp in 'indices' format */
1936 : GEN
1937 259447 : nfarchstar(GEN nf, GEN F, GEN archp)
1938 : {
1939 259447 : long nba = lg(archp) - 1;
1940 259447 : if (!nba) return mkvec2(cgetg(1,t_VEC), archp);
1941 81697 : if (F && equali1(gcoeff(F,1,1))) F = NULL;
1942 81697 : if (F) F = idealpseudored(F, nf_get_roundG(nf));
1943 81684 : return setsigns_init(nf, archp, F, const_vec(nba, gen_2));
1944 : }
1945 :
1946 : /*************************************************************************/
1947 : /** **/
1948 : /** IDEALCHINESE **/
1949 : /** **/
1950 : /*************************************************************************/
1951 : static int
1952 4206 : isprfact(GEN x)
1953 : {
1954 : long i, l;
1955 : GEN L, E;
1956 4206 : if (typ(x) != t_MAT || lg(x) != 3) return 0;
1957 4206 : L = gel(x,1); l = lg(L);
1958 4206 : E = gel(x,2);
1959 13993 : for(i=1; i<l; i++)
1960 : {
1961 9787 : checkprid(gel(L,i));
1962 9787 : if (typ(gel(E,i)) != t_INT) return 0;
1963 : }
1964 4206 : return 1;
1965 : }
1966 :
1967 : /* initialize projectors mod pr[i]^e[i] for idealchinese */
1968 : static GEN
1969 4206 : pr_init(GEN nf, GEN fa, GEN w, GEN dw)
1970 : {
1971 4206 : GEN U, E, F, FZ, L = gel(fa,1), E0 = gel(fa,2);
1972 4206 : long i, r = lg(L);
1973 :
1974 4206 : if (w && lg(w) != r) pari_err_TYPE("idealchinese", w);
1975 4206 : if (r == 1 && !dw) return cgetg(1,t_VEC);
1976 4192 : E = leafcopy(E0); /* do not destroy fa[2] */
1977 13979 : for (i = 1; i < r; i++)
1978 9787 : if (signe(gel(E,i)) < 0) gel(E,i) = gen_0;
1979 4192 : F = factorbackprime(nf, L, E);
1980 4192 : if (dw)
1981 : {
1982 693 : F = ZM_Z_mul(F, dw);
1983 1596 : for (i = 1; i < r; i++)
1984 : {
1985 903 : GEN pr = gel(L,i);
1986 903 : long e = itos(gel(E0,i)), v = idealval(nf, dw, pr);
1987 903 : if (e >= 0)
1988 896 : gel(E,i) = addiu(gel(E,i), v);
1989 7 : else if (v + e <= 0)
1990 0 : F = idealmulpowprime(nf, F, pr, stoi(-v)); /* coprime to pr */
1991 : else
1992 : {
1993 7 : F = idealmulpowprime(nf, F, pr, stoi(e));
1994 7 : gel(E,i) = stoi(v + e);
1995 : }
1996 : }
1997 : }
1998 4192 : U = cgetg(r, t_VEC);
1999 13979 : for (i = 1; i < r; i++)
2000 : {
2001 : GEN u;
2002 9787 : if (w && gequal0(gel(w,i))) u = gen_0; /* unused */
2003 : else
2004 : {
2005 9710 : GEN pr = gel(L,i), e = gel(E,i), t;
2006 9710 : t = idealdivpowprime(nf,F, pr, e);
2007 9710 : u = hnfmerge_get_1(t, idealpow(nf, pr, e));
2008 9710 : if (!u) pari_err_COPRIME("idealchinese", t,pr);
2009 : }
2010 9787 : gel(U,i) = u;
2011 : }
2012 4192 : FZ = gcoeff(F, 1, 1);
2013 4192 : F = idealpseudored(F, nf_get_roundG(nf));
2014 4192 : return mkvec2(mkvec2(F, FZ), U);
2015 : }
2016 :
2017 : static GEN
2018 2261 : pl_normalize(GEN nf, GEN pl)
2019 : {
2020 2261 : const char *fun = "idealchinese";
2021 2261 : if (lg(pl)-1 != nf_get_r1(nf)) pari_err_TYPE(fun,pl);
2022 2261 : switch(typ(pl))
2023 : {
2024 707 : case t_VEC: RgV_check_ZV(pl,fun); pl = ZV_to_zv(pl);
2025 : /* fall through */
2026 2261 : case t_VECSMALL: break;
2027 0 : default: pari_err_TYPE(fun,pl);
2028 : }
2029 2261 : return pl;
2030 : }
2031 :
2032 : static int
2033 9443 : is_chineseinit(GEN x)
2034 : {
2035 : GEN fa, pl;
2036 : long l;
2037 9443 : if (typ(x) != t_VEC || lg(x)!=3) return 0;
2038 7602 : fa = gel(x,1);
2039 7602 : pl = gel(x,2);
2040 7602 : if (typ(fa) != t_VEC || typ(pl) != t_VEC) return 0;
2041 4207 : l = lg(fa);
2042 4207 : if (l != 1)
2043 : {
2044 : GEN z;
2045 4165 : if (l != 3) return 0;
2046 4165 : z = gel(fa, 1);
2047 4165 : if (typ(z) != t_VEC || lg(z) != 3 || typ(gel(z,1)) != t_MAT
2048 4158 : || typ(gel(z,2)) != t_INT
2049 4158 : || typ(gel(fa,2)) != t_VEC)
2050 7 : return 0;
2051 : }
2052 4200 : l = lg(pl);
2053 4200 : if (l != 1)
2054 : {
2055 665 : if (l != 6 || typ(gel(pl,3)) != t_MAT || typ(gel(pl,1)) != t_VECSMALL
2056 665 : || typ(gel(pl,2)) != t_VECSMALL)
2057 0 : return 0;
2058 : }
2059 4200 : return 1;
2060 : }
2061 :
2062 : /* nf a true 'nf' */
2063 : static GEN
2064 4661 : chineseinit_i(GEN nf, GEN fa, GEN w, GEN dw)
2065 : {
2066 4661 : const char *fun = "idealchineseinit";
2067 4661 : GEN archp = NULL, pl = NULL;
2068 4661 : switch(typ(fa))
2069 : {
2070 2261 : case t_VEC:
2071 2261 : if (is_chineseinit(fa))
2072 : {
2073 0 : if (dw) pari_err_DOMAIN(fun, "denom(y)", "!=", gen_1, w);
2074 0 : return fa;
2075 : }
2076 2261 : if (lg(fa) != 3) pari_err_TYPE(fun, fa);
2077 : /* of the form [x,s] */
2078 2261 : pl = pl_normalize(nf, gel(fa,2));
2079 2261 : fa = gel(fa,1);
2080 2261 : archp = vecsmall01_to_indices(pl);
2081 : /* keep pr_init, reset pl */
2082 2261 : if (is_chineseinit(fa)) { fa = gel(fa,1); break; }
2083 : /* fall through */
2084 : case t_MAT: /* factorization? */
2085 4206 : if (isprfact(fa)) { fa = pr_init(nf, fa, w, dw); break; }
2086 0 : default: pari_err_TYPE(fun,fa);
2087 : }
2088 :
2089 4661 : if (!pl) pl = cgetg(1,t_VEC);
2090 : else
2091 : {
2092 2261 : long r = lg(archp);
2093 2261 : if (r == 1) pl = cgetg(1, t_VEC);
2094 : else
2095 : {
2096 1757 : GEN F = (lg(fa) == 1)? NULL: gmael(fa,1,1), signs = cgetg(r, t_VECSMALL);
2097 : long i;
2098 5082 : for (i = 1; i < r; i++) signs[i] = (pl[archp[i]] < 0)? 1: 0;
2099 1757 : pl = setsigns_init(nf, archp, F, signs);
2100 : }
2101 : }
2102 4661 : return mkvec2(fa, pl);
2103 : }
2104 :
2105 : /* Given a prime ideal factorization x, possibly with 0 or negative exponents,
2106 : * and a vector w of elements of nf, gives b such that
2107 : * v_p(b-w_p)>=v_p(x) for all prime ideals p in the ideal factorization
2108 : * and v_p(b)>=0 for all other p, using the standard proof given in GTM 138. */
2109 : GEN
2110 8406 : idealchinese(GEN nf, GEN x0, GEN w)
2111 : {
2112 8406 : const char *fun = "idealchinese";
2113 8406 : pari_sp av = avma;
2114 8406 : GEN x = x0, x1, x2, s, dw, F;
2115 :
2116 8406 : nf = checknf(nf);
2117 8406 : if (!w) return gerepilecopy(av, chineseinit_i(nf,x,NULL,NULL));
2118 :
2119 4921 : if (typ(w) != t_VEC) pari_err_TYPE(fun,w);
2120 4921 : w = Q_remove_denom(matalgtobasis(nf,w), &dw);
2121 4921 : if (!is_chineseinit(x)) x = chineseinit_i(nf,x,w,dw);
2122 : /* x is a 'chineseinit' */
2123 4921 : x1 = gel(x,1); s = NULL;
2124 4921 : x2 = gel(x,2);
2125 4921 : if (lg(x1) == 1) { F = NULL; dw = NULL; }
2126 : else
2127 : {
2128 4879 : GEN U = gel(x1,2), FZ;
2129 4879 : long i, r = lg(w);
2130 4879 : F = gmael(x1,1,1); FZ = gmael(x1,1,2);
2131 17624 : for (i=1; i<r; i++)
2132 12745 : if (!ZV_equal0(gel(w,i)))
2133 : {
2134 9640 : GEN t = nfmuli(nf, gel(U,i), gel(w,i));
2135 9640 : s = s? ZC_add(s,t): t;
2136 : }
2137 4879 : if (s)
2138 : {
2139 4858 : s = ZC_reducemodmatrix(s, F);
2140 4858 : if (dw && x == x0) /* input was a chineseinit */
2141 : {
2142 7 : dw = modii(dw, FZ);
2143 7 : s = FpC_Fp_mul(s, Fp_inv(dw, FZ), FZ);
2144 7 : dw = NULL;
2145 : }
2146 4858 : if (ZV_isscalar(s)) s = icopy(gel(s,1));
2147 : }
2148 : }
2149 4921 : if (lg(x2) != 1)
2150 : {
2151 1764 : s = nfsetsigns(nf, gel(x2,1), s? s: gen_0, x2);
2152 1764 : if (typ(s) == t_COL && QV_isscalar(s))
2153 : {
2154 294 : s = gel(s,1); if (!dw) s = gcopy(s);
2155 : }
2156 : }
2157 3157 : else if (!s) return gc_const(av, gen_0);
2158 4872 : return gerepileupto(av, dw? gdiv(s, dw): s);
2159 : }
2160 :
2161 : /*************************************************************************/
2162 : /** **/
2163 : /** (Z_K/I)^* **/
2164 : /** **/
2165 : /*************************************************************************/
2166 : GEN
2167 2261 : vecsmall01_to_indices(GEN v)
2168 : {
2169 2261 : long i, k, l = lg(v);
2170 2261 : GEN p = new_chunk(l) + l;
2171 6636 : for (k=1, i=l-1; i; i--)
2172 4375 : if (v[i]) { *--p = i; k++; }
2173 2261 : *--p = _evallg(k) | evaltyp(t_VECSMALL);
2174 2261 : set_avma((pari_sp)p); return p;
2175 : }
2176 : GEN
2177 1092461 : vec01_to_indices(GEN v)
2178 : {
2179 : long i, k, l;
2180 : GEN p;
2181 :
2182 1092461 : switch (typ(v))
2183 : {
2184 1045729 : case t_VECSMALL: return v;
2185 46732 : case t_VEC: break;
2186 0 : default: pari_err_TYPE("vec01_to_indices",v);
2187 : }
2188 46732 : l = lg(v);
2189 46732 : p = new_chunk(l) + l;
2190 140553 : for (k=1, i=l-1; i; i--)
2191 93821 : if (signe(gel(v,i))) { *--p = i; k++; }
2192 46732 : *--p = _evallg(k) | evaltyp(t_VECSMALL);
2193 46732 : set_avma((pari_sp)p); return p;
2194 : }
2195 : GEN
2196 136882 : indices_to_vec01(GEN p, long r)
2197 : {
2198 136882 : long i, l = lg(p);
2199 136882 : GEN v = zerovec(r);
2200 206620 : for (i = 1; i < l; i++) gel(v, p[i]) = gen_1;
2201 136880 : return v;
2202 : }
2203 :
2204 : /* return (column) vector of R1 signatures of x (0 or 1) */
2205 : GEN
2206 1045730 : nfsign_arch(GEN nf, GEN x, GEN arch)
2207 : {
2208 1045730 : GEN sarch, M, V, archp = vec01_to_indices(arch);
2209 1045729 : long i, s, np, bigx, n = lg(archp)-1;
2210 : pari_sp av;
2211 :
2212 1045729 : if (!n) return cgetg(1,t_VECSMALL);
2213 844129 : if (typ(x) == t_MAT)
2214 : { /* factorisation */
2215 276147 : GEN g = gel(x,1), e = gel(x,2);
2216 276147 : long l = lg(g);
2217 276147 : V = zero_zv(n);
2218 831600 : for (i = 1; i < l; i++)
2219 555452 : if (mpodd(gel(e,i)))
2220 435803 : Flv_add_inplace(V, nfsign_arch(nf,gel(g,i),archp), 2);
2221 276148 : set_avma((pari_sp)V); return V;
2222 : }
2223 567982 : av = avma; V = cgetg(n+1,t_VECSMALL);
2224 567981 : x = nf_to_scalar_or_basis(nf, x);
2225 567981 : switch(typ(x))
2226 : {
2227 182987 : case t_INT:
2228 182987 : s = signe(x);
2229 182987 : if (!s) pari_err_DOMAIN("nfsign_arch","element","=",gen_0,x);
2230 182987 : set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
2231 644 : case t_FRAC:
2232 644 : s = signe(gel(x,1));
2233 644 : set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
2234 : }
2235 384350 : x = Q_primpart(x); M = nf_get_M(nf); sarch = NULL;
2236 384353 : np = -1; bigx = gexpo(x) >= nf_get_prec(nf);
2237 1173347 : for (i = 1; i <= n; i++)
2238 : {
2239 789831 : long s = bigx ? -1: eval_sign(M, x, archp[i]);
2240 789827 : if (s < 0) /* failure */
2241 : {
2242 870 : long ni, r1 = nf_get_r1(nf);
2243 : GEN xi;
2244 870 : if (np < 0)
2245 : {
2246 870 : np = num_positive(nf, x);
2247 870 : if (np == 0) { set_avma(av); return const_vecsmall(n, 1); }
2248 820 : if (np == r1){ set_avma(av); return const_vecsmall(n, 0); }
2249 420 : sarch = nfarchstar(nf, NULL, identity_perm(r1));
2250 : }
2251 420 : xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
2252 420 : xi = Q_primpart(xi);
2253 420 : ni = num_positive(nf, nfmuli(nf,x,xi));
2254 420 : if (ni == 0) { set_avma(av); V = const_vecsmall(n, 1); V[i] = 0; return V; }
2255 413 : if (ni == r1){ set_avma(av); V = const_vecsmall(n, 0); V[i] = 1; return V; }
2256 35 : s = ni < np? 0: 1;
2257 : }
2258 788992 : V[i] = s;
2259 : }
2260 383516 : set_avma((pari_sp)V); return V;
2261 : }
2262 : static void
2263 35483 : chk_ind(const char *s, long i, long r1)
2264 : {
2265 35483 : if (i <= 0) pari_err_DOMAIN(s, "index", "<=", gen_0, stoi(i));
2266 35469 : if (i > r1) pari_err_DOMAIN(s, "index", ">", utoi(r1), utoi(i));
2267 35434 : }
2268 : static GEN
2269 126385 : parse_embed(GEN ind, long r, const char *f)
2270 : {
2271 : long l, i;
2272 126385 : if (!ind) return identity_perm(r);
2273 33418 : switch(typ(ind))
2274 : {
2275 70 : case t_INT: ind = mkvecsmall(itos(ind)); break;
2276 84 : case t_VEC: case t_COL: ind = vec_to_vecsmall(ind); break;
2277 33264 : case t_VECSMALL: break;
2278 0 : default: pari_err_TYPE(f, ind);
2279 : }
2280 33418 : l = lg(ind);
2281 68852 : for (i = 1; i < l; i++) chk_ind(f, ind[i], r);
2282 33369 : return ind;
2283 : }
2284 : GEN
2285 124068 : nfeltsign(GEN nf, GEN x, GEN ind0)
2286 : {
2287 124068 : pari_sp av = avma;
2288 : long i, l;
2289 : GEN v, ind;
2290 124068 : nf = checknf(nf);
2291 124068 : ind = parse_embed(ind0, nf_get_r1(nf), "nfeltsign");
2292 124047 : l = lg(ind);
2293 124047 : if (is_rational_t(typ(x)))
2294 : { /* nfsign_arch would test this, but avoid converting t_VECSMALL -> t_VEC */
2295 : GEN s;
2296 30975 : switch(gsigne(x))
2297 : {
2298 16366 : case -1:s = gen_m1; break;
2299 14602 : case 1: s = gen_1; break;
2300 7 : default: s = gen_0; break;
2301 : }
2302 30975 : set_avma(av);
2303 30975 : return (ind0 && typ(ind0) == t_INT)? s: const_vec(l-1, s);
2304 : }
2305 93072 : v = nfsign_arch(nf, x, ind);
2306 93072 : if (ind0 && typ(ind0) == t_INT) { set_avma(av); return v[1]? gen_m1: gen_1; }
2307 93058 : settyp(v, t_VEC);
2308 262325 : for (i = 1; i < l; i++) gel(v,i) = v[i]? gen_m1: gen_1;
2309 93058 : return gerepileupto(av, v);
2310 : }
2311 :
2312 : /* true nf */
2313 : GEN
2314 728 : nfeltembed_i(GEN *pnf, GEN x, GEN ind0, long prec0)
2315 : {
2316 : long i, e, l, r1, r2, prec, prec1;
2317 728 : GEN v, ind, cx, nf = *pnf;
2318 728 : nf_get_sign(nf,&r1,&r2);
2319 728 : x = nf_to_scalar_or_basis(nf, x);
2320 721 : ind = parse_embed(ind0, r1+r2, "nfeltembed");
2321 714 : l = lg(ind);
2322 714 : if (typ(x) != t_COL)
2323 : {
2324 224 : if (!(ind0 && typ(ind0) == t_INT)) x = const_vec(l-1, x);
2325 224 : return x;
2326 : }
2327 490 : x = Q_primitive_part(x, &cx);
2328 490 : prec1 = prec0; e = gexpo(x);
2329 490 : if (e > 8) prec1 += nbits2extraprec(e);
2330 490 : prec = prec1;
2331 490 : if (nf_get_prec(nf) < prec) nf = nfnewprec_shallow(nf, prec);
2332 490 : v = cgetg(l, t_VEC);
2333 : for(;;)
2334 132 : {
2335 622 : GEN M = nf_get_M(nf);
2336 2624 : for (i = 1; i < l; i++)
2337 : {
2338 2134 : GEN t = nfembed_i(M, x, ind[i]);
2339 2134 : long e = gexpo(t);
2340 2134 : if (gequal0(t) || precision(t) < prec0
2341 2134 : || (e < 0 && prec < prec1 + nbits2extraprec(-e)) ) break;
2342 2002 : if (cx) t = gmul(t, cx);
2343 2002 : gel(v,i) = t;
2344 : }
2345 622 : if (i == l) break;
2346 132 : prec = precdbl(prec);
2347 132 : if (DEBUGLEVEL>1) pari_warn(warnprec,"eltnfembed", prec);
2348 132 : *pnf = nf = nfnewprec_shallow(nf, prec);
2349 : }
2350 490 : if (ind0 && typ(ind0) == t_INT) v = gel(v,1);
2351 490 : return v;
2352 : }
2353 : GEN
2354 728 : nfeltembed(GEN nf, GEN x, GEN ind0, long prec0)
2355 : {
2356 728 : pari_sp av = avma; nf = checknf(nf);
2357 728 : return gerepilecopy(av, nfeltembed_i(&nf, x, ind0, prec0));
2358 : }
2359 :
2360 : /* number of distinct roots of sigma(f) */
2361 : GEN
2362 1596 : nfpolsturm(GEN nf, GEN f, GEN ind0)
2363 : {
2364 1596 : pari_sp av = avma;
2365 : long d, l, r1, single;
2366 : GEN ind, u, v, vr1, T, s, t;
2367 :
2368 1596 : nf = checknf(nf); T = nf_get_pol(nf); r1 = nf_get_r1(nf);
2369 1596 : ind = parse_embed(ind0, r1, "nfpolsturm");
2370 1575 : single = ind0 && typ(ind0) == t_INT;
2371 1575 : l = lg(ind);
2372 :
2373 1575 : if (gequal0(f)) pari_err_ROOTS0("nfpolsturm");
2374 1568 : if (typ(f) == t_POL && varn(f) != varn(T))
2375 : {
2376 1547 : f = RgX_nffix("nfpolsturm", T, f,1);
2377 1547 : if (lg(f) == 3) f = NULL;
2378 : }
2379 : else
2380 : {
2381 21 : (void)Rg_nffix("nfpolsturm", T, f, 0);
2382 21 : f = NULL;
2383 : }
2384 1568 : if (!f) { set_avma(av); return single? gen_0: zerovec(l-1); }
2385 1547 : d = degpol(f);
2386 1547 : if (d == 1) { set_avma(av); return single? gen_1: const_vec(l-1,gen_1); }
2387 :
2388 1505 : vr1 = const_vecsmall(l-1, 1);
2389 1505 : u = Q_primpart(f); s = ZV_to_zv(nfeltsign(nf, gel(u,d+2), ind));
2390 1505 : v = RgX_deriv(u); t = odd(d)? leafcopy(s): zv_neg(s);
2391 : for(;;)
2392 182 : {
2393 1687 : GEN r = RgX_neg( Q_primpart(RgX_pseudorem(u, v)) ), sr;
2394 1687 : long i, dr = degpol(r);
2395 1687 : if (dr < 0) break;
2396 1687 : sr = ZV_to_zv(nfeltsign(nf, gel(r,dr+2), ind));
2397 4144 : for (i = 1; i < l; i++)
2398 2457 : if (sr[i] != s[i]) { s[i] = sr[i], vr1[i]--; }
2399 1687 : if (odd(dr)) sr = zv_neg(sr);
2400 4144 : for (i = 1; i < l; i++)
2401 2457 : if (sr[i] != t[i]) { t[i] = sr[i], vr1[i]++; }
2402 1687 : if (!dr) break;
2403 182 : u = v; v = r;
2404 : }
2405 1505 : if (single) return gc_stoi(av,vr1[1]);
2406 1498 : return gerepileupto(av, zv_to_ZV(vr1));
2407 : }
2408 :
2409 : /* True nf; return the vector of signs of x; the matrix of such if x is a vector
2410 : * of nf elements */
2411 : GEN
2412 43959 : nfsign(GEN nf, GEN x)
2413 : {
2414 : long i, l;
2415 : GEN archp, S;
2416 :
2417 43959 : archp = identity_perm( nf_get_r1(nf) );
2418 43959 : if (typ(x) != t_VEC) return nfsign_arch(nf, x, archp);
2419 35937 : l = lg(x); S = cgetg(l, t_MAT);
2420 148055 : for (i=1; i<l; i++) gel(S,i) = nfsign_arch(nf, gel(x,i), archp);
2421 35936 : return S;
2422 : }
2423 :
2424 : /* x integral elt, A integral ideal in HNF; reduce x mod A */
2425 : static GEN
2426 7806804 : zk_modHNF(GEN x, GEN A)
2427 7806804 : { return (typ(x) == t_COL)? ZC_hnfrem(x, A): modii(x, gcoeff(A,1,1)); }
2428 :
2429 : /* given an element x in Z_K and an integral ideal y in HNF, coprime with x,
2430 : outputs an element inverse of x modulo y */
2431 : GEN
2432 189 : nfinvmodideal(GEN nf, GEN x, GEN y)
2433 : {
2434 189 : pari_sp av = avma;
2435 189 : GEN a, yZ = gcoeff(y,1,1);
2436 :
2437 189 : if (equali1(yZ)) return gen_0;
2438 189 : x = nf_to_scalar_or_basis(nf, x);
2439 189 : if (typ(x) == t_INT) return gerepileupto(av, Fp_inv(x, yZ));
2440 :
2441 79 : a = hnfmerge_get_1(idealhnf_principal(nf,x), y);
2442 79 : if (!a) pari_err_INV("nfinvmodideal", x);
2443 79 : return gerepileupto(av, zk_modHNF(nfdiv(nf,a,x), y));
2444 : }
2445 :
2446 : static GEN
2447 2682291 : nfsqrmodideal(GEN nf, GEN x, GEN id)
2448 2682291 : { return zk_modHNF(nfsqri(nf,x), id); }
2449 : static GEN
2450 7284813 : nfmulmodideal(GEN nf, GEN x, GEN y, GEN id)
2451 7284813 : { return x? zk_modHNF(nfmuli(nf,x,y), id): y; }
2452 : /* assume x integral, k integer, A in HNF */
2453 : GEN
2454 5841252 : nfpowmodideal(GEN nf,GEN x,GEN k,GEN A)
2455 : {
2456 5841252 : long s = signe(k);
2457 : pari_sp av;
2458 : GEN y;
2459 :
2460 5841252 : if (!s) return gen_1;
2461 5841252 : av = avma;
2462 5841252 : x = nf_to_scalar_or_basis(nf, x);
2463 5841548 : if (typ(x) != t_COL) return Fp_pow(x, k, gcoeff(A,1,1));
2464 2625729 : if (s < 0) { k = negi(k); x = nfinvmodideal(nf, x,A); }
2465 2625729 : if (equali1(k)) return gerepileupto(av, s > 0? zk_modHNF(x, A): x);
2466 1148734 : for(y = NULL;;)
2467 : {
2468 3831080 : if (mpodd(k)) y = nfmulmodideal(nf,y,x,A);
2469 3831080 : k = shifti(k,-1); if (!signe(k)) break;
2470 2681947 : x = nfsqrmodideal(nf,x,A);
2471 : }
2472 1148725 : return gerepileupto(av, y);
2473 : }
2474 :
2475 : /* a * g^n mod id */
2476 : static GEN
2477 4691934 : nfmulpowmodideal(GEN nf, GEN a, GEN g, GEN n, GEN id)
2478 : {
2479 4691934 : return nfmulmodideal(nf, a, nfpowmodideal(nf,g,n,id), id);
2480 : }
2481 :
2482 : /* assume (num(g[i]), id) = 1 for all i. Return prod g[i]^e[i] mod id.
2483 : * EX = multiple of exponent of (O_K/id)^* */
2484 : GEN
2485 2620654 : famat_to_nf_modideal_coprime(GEN nf, GEN g, GEN e, GEN id, GEN EX)
2486 : {
2487 2620654 : GEN EXo2, plus = NULL, minus = NULL, idZ = gcoeff(id,1,1);
2488 2620654 : long i, lx = lg(g);
2489 :
2490 2620654 : if (equali1(idZ)) return gen_1; /* id = Z_K */
2491 2620160 : EXo2 = (expi(EX) > 10)? shifti(EX,-1): NULL;
2492 8328695 : for (i = 1; i < lx; i++)
2493 : {
2494 5708642 : GEN h, n = centermodii(gel(e,i), EX, EXo2);
2495 5708136 : long sn = signe(n);
2496 5708136 : if (!sn) continue;
2497 :
2498 4039634 : h = nf_to_scalar_or_basis(nf, gel(g,i));
2499 4040065 : switch(typ(h))
2500 : {
2501 2383256 : case t_INT: break;
2502 0 : case t_FRAC:
2503 0 : h = Fp_div(gel(h,1), gel(h,2), idZ); break;
2504 1656809 : default:
2505 : {
2506 : GEN dh;
2507 1656809 : h = Q_remove_denom(h, &dh);
2508 1656968 : if (dh) h = FpC_Fp_mul(h, Fp_inv(dh,idZ), idZ);
2509 : }
2510 : }
2511 4040111 : if (sn > 0)
2512 4038276 : plus = nfmulpowmodideal(nf, plus, h, n, id);
2513 : else /* sn < 0 */
2514 1835 : minus = nfmulpowmodideal(nf, minus, h, negi(n), id);
2515 : }
2516 2620053 : if (minus) plus = nfmulmodideal(nf, plus, nfinvmodideal(nf,minus,id), id);
2517 2620167 : return plus? plus: gen_1;
2518 : }
2519 :
2520 : /* given 2 integral ideals x, y in HNF s.t x | y | x^2, compute (1+x)/(1+y) in
2521 : * the form [[cyc],[gen], U], where U := ux^-1 as a pair [ZM, denom(U)] */
2522 : static GEN
2523 236929 : zidealij(GEN x, GEN y)
2524 : {
2525 236929 : GEN U, G, cyc, xp = gcoeff(x,1,1), xi = hnf_invscale(x, xp);
2526 : long j, N;
2527 :
2528 : /* x^(-1) y = relations between the 1 + x_i (HNF) */
2529 236925 : cyc = ZM_snf_group(ZM_Z_divexact(ZM_mul(xi, y), xp), &U, &G);
2530 236930 : N = lg(cyc); G = ZM_mul(x,G); settyp(G, t_VEC); /* new generators */
2531 573857 : for (j=1; j<N; j++)
2532 : {
2533 336956 : GEN c = gel(G,j);
2534 336956 : gel(c,1) = addiu(gel(c,1), 1); /* 1 + g_j */
2535 336940 : if (ZV_isscalar(c)) gel(G,j) = gel(c,1);
2536 : }
2537 236901 : return mkvec4(cyc, G, ZM_mul(U,xi), xp);
2538 : }
2539 :
2540 : /* lg(x) > 1, x + 1; shallow */
2541 : static GEN
2542 169703 : ZC_add1(GEN x)
2543 : {
2544 169703 : long i, l = lg(x);
2545 169703 : GEN y = cgetg(l, t_COL);
2546 396222 : for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
2547 169702 : gel(y,1) = addiu(gel(x,1), 1); return y;
2548 : }
2549 : /* lg(x) > 1, x - 1; shallow */
2550 : static GEN
2551 70497 : ZC_sub1(GEN x)
2552 : {
2553 70497 : long i, l = lg(x);
2554 70497 : GEN y = cgetg(l, t_COL);
2555 176939 : for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
2556 70497 : gel(y,1) = subiu(gel(x,1), 1); return y;
2557 : }
2558 :
2559 : /* x,y are t_INT or ZC */
2560 : static GEN
2561 0 : zkadd(GEN x, GEN y)
2562 : {
2563 0 : long tx = typ(x);
2564 0 : if (tx == typ(y))
2565 0 : return tx == t_INT? addii(x,y): ZC_add(x,y);
2566 : else
2567 0 : return tx == t_INT? ZC_Z_add(y,x): ZC_Z_add(x,y);
2568 : }
2569 : /* x a t_INT or ZC, x+1; shallow */
2570 : static GEN
2571 255364 : zkadd1(GEN x)
2572 : {
2573 255364 : long tx = typ(x);
2574 255364 : return tx == t_INT? addiu(x,1): ZC_add1(x);
2575 : }
2576 : /* x a t_INT or ZC, x-1; shallow */
2577 : static GEN
2578 255410 : zksub1(GEN x)
2579 : {
2580 255410 : long tx = typ(x);
2581 255410 : return tx == t_INT? subiu(x,1): ZC_sub1(x);
2582 : }
2583 : /* x,y are t_INT or ZC; x - y */
2584 : static GEN
2585 0 : zksub(GEN x, GEN y)
2586 : {
2587 0 : long tx = typ(x), ty = typ(y);
2588 0 : if (tx == ty)
2589 0 : return tx == t_INT? subii(x,y): ZC_sub(x,y);
2590 : else
2591 0 : return tx == t_INT? Z_ZC_sub(x,y): ZC_Z_sub(x,y);
2592 : }
2593 : /* x is t_INT or ZM (mult. map), y is t_INT or ZC; x * y */
2594 : static GEN
2595 255366 : zkmul(GEN x, GEN y)
2596 : {
2597 255366 : long tx = typ(x), ty = typ(y);
2598 255366 : if (ty == t_INT)
2599 184903 : return tx == t_INT? mulii(x,y): ZC_Z_mul(gel(x,1),y);
2600 : else
2601 70463 : return tx == t_INT? ZC_Z_mul(y,x): ZM_ZC_mul(x,y);
2602 : }
2603 :
2604 : /* (U,V) = 1 coprime ideals. Want z = x mod U, = y mod V; namely
2605 : * z =vx + uy = v(x-y) + y, where u + v = 1, u in U, v in V.
2606 : * zkc = [v, UV], v a t_INT or ZM (mult. by v map), UV a ZM (ideal in HNF);
2607 : * shallow */
2608 : GEN
2609 0 : zkchinese(GEN zkc, GEN x, GEN y)
2610 : {
2611 0 : GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd(zkmul(v, zksub(x,y)), y);
2612 0 : return zk_modHNF(z, UV);
2613 : }
2614 : /* special case z = x mod U, = 1 mod V; shallow */
2615 : GEN
2616 255411 : zkchinese1(GEN zkc, GEN x)
2617 : {
2618 255411 : GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd1(zkmul(v, zksub1(x)));
2619 255392 : return (typ(z) == t_INT)? z: ZC_hnfrem(z, UV);
2620 : }
2621 : static GEN
2622 237392 : zkVchinese1(GEN zkc, GEN v)
2623 : {
2624 : long i, ly;
2625 237392 : GEN y = cgetg_copy(v, &ly);
2626 492764 : for (i=1; i<ly; i++) gel(y,i) = zkchinese1(zkc, gel(v,i));
2627 237354 : return y;
2628 : }
2629 :
2630 : /* prepare to solve z = x (mod A), z = y mod (B) [zkchinese or zkchinese1] */
2631 : GEN
2632 237152 : zkchineseinit(GEN nf, GEN A, GEN B, GEN AB)
2633 : {
2634 237152 : GEN v = idealaddtoone_raw(nf, A, B);
2635 : long e;
2636 237145 : if ((e = gexpo(v)) > 5)
2637 : {
2638 83279 : GEN b = (typ(v) == t_COL)? v: scalarcol_shallow(v, nf_get_degree(nf));
2639 83279 : b= ZC_reducemodlll(b, AB);
2640 83285 : if (gexpo(b) < e) v = b;
2641 : }
2642 237142 : return mkvec2(zk_scalar_or_multable(nf,v), AB);
2643 : }
2644 : /* prepare to solve z = x (mod A), z = 1 mod (B)
2645 : * and then z = 1 (mod A), z = y mod (B) [zkchinese1 twice] */
2646 : static GEN
2647 259 : zkchinese1init2(GEN nf, GEN A, GEN B, GEN AB)
2648 : {
2649 259 : GEN zkc = zkchineseinit(nf, A, B, AB);
2650 259 : GEN mv = gel(zkc,1), mu;
2651 259 : if (typ(mv) == t_INT) return mkvec2(zkc, mkvec2(subui(1,mv),AB));
2652 35 : mu = RgM_Rg_add_shallow(ZM_neg(mv), gen_1);
2653 35 : return mkvec2(mkvec2(mv,AB), mkvec2(mu,AB));
2654 : }
2655 :
2656 : static GEN
2657 2151558 : apply_U(GEN L, GEN a)
2658 : {
2659 2151558 : GEN e, U = gel(L,3), dU = gel(L,4);
2660 2151558 : if (typ(a) == t_INT)
2661 671620 : e = ZC_Z_mul(gel(U,1), subiu(a, 1));
2662 : else
2663 : { /* t_COL */
2664 1479938 : GEN t = shallowcopy(a);
2665 1479987 : gel(t,1) = subiu(gel(t,1), 1); /* t = a - 1 */
2666 1479879 : e = ZM_ZC_mul(U, t);
2667 : }
2668 2151452 : return gdiv(e, dU);
2669 : }
2670 :
2671 : /* true nf; vectors of [[cyc],[g],U.X^-1]. Assume k > 1. */
2672 : static GEN
2673 169130 : principal_units(GEN nf, GEN pr, long k, GEN prk)
2674 : {
2675 : GEN list, prb;
2676 169130 : ulong mask = quadratic_prec_mask(k);
2677 169130 : long a = 1;
2678 :
2679 169130 : prb = pr_hnf(nf,pr);
2680 169126 : list = vectrunc_init(k);
2681 406052 : while (mask > 1)
2682 : {
2683 236930 : GEN pra = prb;
2684 236930 : long b = a << 1;
2685 :
2686 236930 : if (mask & 1) b--;
2687 236930 : mask >>= 1;
2688 : /* compute 1 + pr^a / 1 + pr^b, 2a <= b */
2689 236930 : prb = (b >= k)? prk: idealpows(nf,pr,b);
2690 236930 : vectrunc_append(list, zidealij(pra, prb));
2691 236929 : a = b;
2692 : }
2693 169122 : return list;
2694 : }
2695 : /* a = 1 mod (pr) return log(a) on local-gens of 1+pr/1+pr^k */
2696 : static GEN
2697 1329349 : log_prk1(GEN nf, GEN a, long nh, GEN L2, GEN prk)
2698 : {
2699 1329349 : GEN y = cgetg(nh+1, t_COL);
2700 1329361 : long j, iy, c = lg(L2)-1;
2701 3480845 : for (j = iy = 1; j <= c; j++)
2702 : {
2703 2151554 : GEN L = gel(L2,j), cyc = gel(L,1), gen = gel(L,2), E = apply_U(L,a);
2704 2151357 : long i, nc = lg(cyc)-1;
2705 2151357 : int last = (j == c);
2706 5814250 : for (i = 1; i <= nc; i++, iy++)
2707 : {
2708 3662766 : GEN t, e = gel(E,i);
2709 3662766 : if (typ(e) != t_INT) pari_err_COPRIME("zlog_prk1", a, prk);
2710 3662759 : t = Fp_neg(e, gel(cyc,i));
2711 3662750 : gel(y,iy) = negi(t);
2712 3662903 : if (!last && signe(t)) a = nfmulpowmodideal(nf, a, gel(gen,i), t, prk);
2713 : }
2714 : }
2715 1329291 : return y;
2716 : }
2717 : /* true nf */
2718 : static GEN
2719 56657 : principal_units_relations(GEN nf, GEN L2, GEN prk, long nh)
2720 : {
2721 56657 : GEN h = cgetg(nh+1,t_MAT);
2722 56656 : long ih, j, c = lg(L2)-1;
2723 181121 : for (j = ih = 1; j <= c; j++)
2724 : {
2725 124464 : GEN L = gel(L2,j), F = gel(L,1), G = gel(L,2);
2726 124464 : long k, lG = lg(G);
2727 304290 : for (k = 1; k < lG; k++,ih++)
2728 : { /* log(g^f) mod pr^e */
2729 179825 : GEN a = nfpowmodideal(nf,gel(G,k),gel(F,k),prk);
2730 179824 : gel(h,ih) = ZC_neg(log_prk1(nf, a, nh, L2, prk));
2731 179826 : gcoeff(h,ih,ih) = gel(F,k);
2732 : }
2733 : }
2734 56657 : return h;
2735 : }
2736 : /* true nf; k > 1; multiplicative group (1 + pr) / (1 + pr^k) */
2737 : static GEN
2738 169124 : idealprincipalunits_i(GEN nf, GEN pr, long k, GEN *pU)
2739 : {
2740 169124 : GEN cyc, gen, L2, prk = idealpows(nf, pr, k);
2741 :
2742 169130 : L2 = principal_units(nf, pr, k, prk);
2743 169127 : if (k == 2)
2744 : {
2745 112470 : GEN L = gel(L2,1);
2746 112470 : cyc = gel(L,1);
2747 112470 : gen = gel(L,2);
2748 112470 : if (pU) *pU = matid(lg(gen)-1);
2749 : }
2750 : else
2751 : {
2752 56657 : long c = lg(L2), j;
2753 56657 : GEN EX, h, Ui, vg = cgetg(c, t_VEC);
2754 181122 : for (j = 1; j < c; j++) gel(vg, j) = gmael(L2,j,2);
2755 56657 : vg = shallowconcat1(vg);
2756 56657 : h = principal_units_relations(nf, L2, prk, lg(vg)-1);
2757 56657 : h = ZM_hnfall_i(h, NULL, 0);
2758 56657 : cyc = ZM_snf_group(h, pU, &Ui);
2759 56657 : c = lg(Ui); gen = cgetg(c, t_VEC); EX = cyc_get_expo(cyc);
2760 188297 : for (j = 1; j < c; j++)
2761 131640 : gel(gen,j) = famat_to_nf_modideal_coprime(nf, vg, gel(Ui,j), prk, EX);
2762 : }
2763 169126 : return mkvec4(cyc, gen, prk, L2);
2764 : }
2765 : GEN
2766 154 : idealprincipalunits(GEN nf, GEN pr, long k)
2767 : {
2768 : pari_sp av;
2769 : GEN v;
2770 154 : nf = checknf(nf);
2771 154 : if (k == 1) { checkprid(pr); retmkvec3(gen_1,cgetg(1,t_VEC),cgetg(1,t_VEC)); }
2772 147 : av = avma; v = idealprincipalunits_i(nf, pr, k, NULL);
2773 147 : return gerepilecopy(av, mkvec3(powiu(pr_norm(pr), k-1), gel(v,1), gel(v,2)));
2774 : }
2775 :
2776 : /* true nf; given an ideal pr^k dividing an integral ideal x (in HNF form)
2777 : * compute an 'sprk', the structure of G = (Z_K/pr^k)^* [ x = NULL for x=pr^k ]
2778 : * Return a vector with at least 4 components [cyc],[gen],[HNF pr^k,pr,k],ff,
2779 : * where
2780 : * cyc : type of G as abelian group (SNF)
2781 : * gen : generators of G, coprime to x
2782 : * pr^k: in HNF
2783 : * ff : data for log_g in (Z_K/pr)^*
2784 : * Two extra components are present iff k > 1: L2, U
2785 : * L2 : list of data structures to compute local DL in (Z_K/pr)^*,
2786 : * and 1 + pr^a/ 1 + pr^b for various a < b <= min(2a, k)
2787 : * U : base change matrices to convert a vector of local DL to DL wrt gen
2788 : * If MOD is not NULL, initialize G / G^MOD instead */
2789 : static GEN
2790 425855 : sprkinit(GEN nf, GEN pr, long k, GEN x, GEN MOD)
2791 : {
2792 425855 : GEN T, p, Ld, modpr, cyc, gen, g, g0, A, prk, U, L2, ord0 = NULL;
2793 425855 : long f = pr_get_f(pr);
2794 :
2795 425855 : if(DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",k,pr);
2796 425855 : modpr = nf_to_Fq_init(nf, &pr,&T,&p);
2797 425879 : if (MOD)
2798 : {
2799 378335 : GEN o = subiu(powiu(p,f), 1), d = gcdii(o, MOD), fa = Z_factor(d);
2800 378294 : ord0 = mkvec2(o, fa); /* true order, factorization of order in G/G^MOD */
2801 378293 : Ld = gel(fa,1);
2802 378293 : if (lg(Ld) > 1 && equaliu(gel(Ld,1),2)) Ld = vecslice(Ld,2,lg(Ld)-1);
2803 : }
2804 : /* (Z_K / pr)^* */
2805 425840 : if (f == 1)
2806 : {
2807 336695 : g0 = g = MOD? pgener_Fp_local(p, Ld): pgener_Fp(p);
2808 336714 : if (!ord0) ord0 = get_arith_ZZM(subiu(p,1));
2809 : }
2810 : else
2811 : {
2812 89145 : g0 = g = MOD? gener_FpXQ_local(T, p, Ld): gener_FpXQ(T,p, &ord0);
2813 89144 : g = Fq_to_nf(g, modpr);
2814 89144 : if (typ(g) == t_POL) g = poltobasis(nf, g);
2815 : }
2816 425878 : A = gel(ord0, 1); /* Norm(pr)-1 */
2817 : /* If MOD != NULL, d = gcd(A, MOD): g^(A/d) has order d */
2818 425878 : if (k == 1)
2819 : {
2820 256901 : cyc = mkvec(A);
2821 256901 : gen = mkvec(g);
2822 256896 : prk = pr_hnf(nf,pr);
2823 256910 : L2 = U = NULL;
2824 : }
2825 : else
2826 : { /* local-gens of (1 + pr)/(1 + pr^k) = SNF-gens * U */
2827 : GEN AB, B, u, v, w;
2828 : long j, l;
2829 168977 : w = idealprincipalunits_i(nf, pr, k, &U);
2830 : /* incorporate (Z_K/pr)^*, order A coprime to B = expo(1+pr/1+pr^k)*/
2831 168978 : cyc = leafcopy(gel(w,1)); B = cyc_get_expo(cyc); AB = mulii(A,B);
2832 168963 : gen = leafcopy(gel(w,2));
2833 168962 : prk = gel(w,3);
2834 168962 : g = nfpowmodideal(nf, g, B, prk);
2835 168976 : g0 = Fq_pow(g0, modii(B,A), T, p); /* update primitive root */
2836 168976 : L2 = mkvec3(A, g, gel(w,4));
2837 168976 : gel(cyc,1) = AB;
2838 168976 : gel(gen,1) = nfmulmodideal(nf, gel(gen,1), g, prk);
2839 168973 : u = mulii(Fp_inv(A,B), A);
2840 168959 : v = subui(1, u); l = lg(U);
2841 505459 : for (j = 1; j < l; j++) gcoeff(U,1,j) = Fp_mul(u, gcoeff(U,1,j), AB);
2842 168967 : U = mkvec2(Rg_col_ei(v, lg(gen)-1, 1), U);
2843 : }
2844 : /* local-gens of (Z_K/pr^k)^* = SNF-gens * U */
2845 425886 : if (x)
2846 : {
2847 236898 : GEN uv = zkchineseinit(nf, idealmulpowprime(nf,x,pr,utoineg(k)), prk, x);
2848 236874 : gen = zkVchinese1(uv, gen);
2849 : }
2850 425823 : return mkvecn(U? 6: 4, cyc, gen, prk, mkvec3(modpr,g0,ord0), L2, U);
2851 : }
2852 : GEN
2853 3980747 : sprk_get_cyc(GEN s) { return gel(s,1); }
2854 : GEN
2855 1968678 : sprk_get_expo(GEN s) { return cyc_get_expo(sprk_get_cyc(s)); }
2856 : GEN
2857 335792 : sprk_get_gen(GEN s) { return gel(s,2); }
2858 : GEN
2859 4913733 : sprk_get_prk(GEN s) { return gel(s,3); }
2860 : GEN
2861 2542292 : sprk_get_ff(GEN s) { return gel(s,4); }
2862 : GEN
2863 2602999 : sprk_get_pr(GEN s) { GEN ff = gel(s,4); return modpr_get_pr(gel(ff,1)); }
2864 : /* L2 to 1 + pr / 1 + pr^k */
2865 : static GEN
2866 1211783 : sprk_get_L2(GEN s) { return gmael(s,5,3); }
2867 : /* lift to nf of primitive root of k(pr) */
2868 : static GEN
2869 318199 : sprk_get_gnf(GEN s) { return gmael(s,5,2); }
2870 : /* A = Npr-1, <g> = (Z_K/pr)^*, L2 to 1 + pr / 1 + pr^k */
2871 : void
2872 0 : sprk_get_AgL2(GEN s, GEN *A, GEN *g, GEN *L2)
2873 0 : { GEN v = gel(s,5); *A = gel(v,1); *g = gel(v,2); *L2 = gel(v,3); }
2874 : void
2875 1203090 : sprk_get_U2(GEN s, GEN *U1, GEN *U2)
2876 1203090 : { GEN v = gel(s,6); *U1 = gel(v,1); *U2 = gel(v,2); }
2877 : static int
2878 2542292 : sprk_is_prime(GEN s) { return lg(s) == 5; }
2879 :
2880 : GEN
2881 1968481 : famat_zlog_pr(GEN nf, GEN g, GEN e, GEN sprk, GEN mod)
2882 : {
2883 1968481 : GEN x, expo = sprk_get_expo(sprk);
2884 1968481 : if (mod) expo = gcdii(expo,mod);
2885 1968471 : x = famat_makecoprime(nf, g, e, sprk_get_pr(sprk), sprk_get_prk(sprk), expo);
2886 1968478 : return log_prk(nf, x, sprk, mod);
2887 : }
2888 : /* famat_zlog_pr assuming (g,sprk.pr) = 1 */
2889 : static GEN
2890 196 : famat_zlog_pr_coprime(GEN nf, GEN g, GEN e, GEN sprk, GEN MOD)
2891 : {
2892 196 : GEN x = famat_to_nf_modideal_coprime(nf, g, e, sprk_get_prk(sprk),
2893 : sprk_get_expo(sprk));
2894 196 : return log_prk(nf, x, sprk, MOD);
2895 : }
2896 :
2897 : /* o t_INT, O = [ord,fa] format for multiple of o (for Fq_log);
2898 : * return o in [ord,fa] format */
2899 : static GEN
2900 559618 : order_update(GEN o, GEN O)
2901 : {
2902 559618 : GEN p = gmael(O,2,1), z = o, P, E;
2903 559618 : long i, j, l = lg(p);
2904 559618 : P = cgetg(l, t_COL);
2905 559613 : E = cgetg(l, t_COL);
2906 616759 : for (i = j = 1; i < l; i++)
2907 : {
2908 616758 : long v = Z_pvalrem(z, gel(p,i), &z);
2909 616708 : if (v)
2910 : {
2911 603671 : gel(P,j) = gel(p,i);
2912 603671 : gel(E,j) = utoipos(v); j++;
2913 603692 : if (is_pm1(z)) break;
2914 : }
2915 : }
2916 559578 : setlg(P, j);
2917 559571 : setlg(E, j); return mkvec2(o, mkmat2(P,E));
2918 : }
2919 :
2920 : /* a in Z_K (t_COL or t_INT), pr prime ideal, sprk = sprkinit(nf,pr,k,x),
2921 : * mod positive t_INT or NULL (meaning mod=0).
2922 : * return log(a) modulo mod on SNF-generators of (Z_K/pr^k)^* */
2923 : GEN
2924 2616109 : log_prk(GEN nf, GEN a, GEN sprk, GEN mod)
2925 : {
2926 : GEN e, prk, g, U1, U2, y, ff, O, o, oN, gN, N, T, p, modpr, pr, cyc;
2927 :
2928 2616109 : if (typ(a) == t_MAT) return famat_zlog_pr(nf, gel(a,1), gel(a,2), sprk, mod);
2929 2542287 : N = NULL;
2930 2542287 : ff = sprk_get_ff(sprk);
2931 2542290 : pr = gel(ff,1); /* modpr */
2932 2542290 : g = gN = gel(ff,2);
2933 2542290 : O = gel(ff,3); /* order of g = |Fq^*|, in [ord, fa] format */
2934 2542290 : o = oN = gel(O,1); /* order as a t_INT */
2935 2542290 : prk = sprk_get_prk(sprk);
2936 2542303 : modpr = nf_to_Fq_init(nf, &pr, &T, &p);
2937 2542308 : if (mod)
2938 : {
2939 2025953 : GEN d = gcdii(o,mod);
2940 2025693 : if (!equalii(o, d))
2941 : {
2942 749835 : N = diviiexact(o,d); /* > 1, coprime to p */
2943 749756 : a = nfpowmodideal(nf, a, N, prk);
2944 749976 : oN = d; /* order of g^N mod pr */
2945 : }
2946 : }
2947 2542123 : if (equali1(oN))
2948 397550 : e = gen_0;
2949 : else
2950 : {
2951 2144657 : if (N) { O = order_update(oN, O); gN = Fq_pow(g, N, T, p); }
2952 2144651 : e = Fq_log(nf_to_Fq(nf,a,modpr), gN, O, T, p);
2953 : }
2954 : /* 0 <= e < oN is correct modulo oN */
2955 2542315 : if (sprk_is_prime(sprk)) return mkcol(e); /* k = 1 */
2956 :
2957 800310 : sprk_get_U2(sprk, &U1,&U2);
2958 800405 : cyc = sprk_get_cyc(sprk);
2959 800411 : if (mod)
2960 : {
2961 379055 : cyc = ZV_snf_gcd(cyc, mod);
2962 379039 : if (signe(remii(mod,p))) return ZV_ZV_mod(ZC_Z_mul(U1,e), cyc);
2963 : }
2964 746793 : if (signe(e))
2965 : {
2966 318199 : GEN E = N? mulii(e, N): e;
2967 318199 : a = nfmulpowmodideal(nf, a, sprk_get_gnf(sprk), Fp_neg(E, o), prk);
2968 : }
2969 : /* a = 1 mod pr */
2970 746793 : y = log_prk1(nf, a, lg(U2)-1, sprk_get_L2(sprk), prk);
2971 746830 : if (N)
2972 : { /* from DL(a^N) to DL(a) */
2973 135266 : GEN E = gel(sprk_get_cyc(sprk), 1), q = powiu(p, Z_pval(E, p));
2974 135264 : y = ZC_Z_mul(y, Fp_inv(N, q));
2975 : }
2976 746831 : y = ZC_lincomb(gen_1, e, ZM_ZC_mul(U2,y), U1);
2977 746824 : return ZV_ZV_mod(y, cyc);
2978 : }
2979 : /* true nf */
2980 : GEN
2981 90145 : log_prk_init(GEN nf, GEN pr, long k, GEN MOD)
2982 90145 : { return sprkinit(nf,pr,k,NULL,MOD);}
2983 : GEN
2984 497 : veclog_prk(GEN nf, GEN v, GEN sprk)
2985 : {
2986 497 : long l = lg(v), i;
2987 497 : GEN w = cgetg(l, t_MAT);
2988 1232 : for (i = 1; i < l; i++) gel(w,i) = log_prk(nf, gel(v,i), sprk, NULL);
2989 497 : return w;
2990 : }
2991 :
2992 : static GEN
2993 1373983 : famat_zlog(GEN nf, GEN fa, GEN sgn, zlog_S *S)
2994 : {
2995 1373983 : long i, l0, l = lg(S->U);
2996 1373983 : GEN g = gel(fa,1), e = gel(fa,2), y = cgetg(l, t_COL);
2997 1373983 : l0 = lg(S->sprk); /* = l (trivial arch. part), or l-1 */
2998 2851158 : for (i=1; i < l0; i++) gel(y,i) = famat_zlog_pr(nf, g, e, gel(S->sprk,i), S->mod);
2999 1373975 : if (l0 != l)
3000 : {
3001 190755 : if (!sgn) sgn = nfsign_arch(nf, fa, S->archp);
3002 190755 : gel(y,l0) = Flc_to_ZC(sgn);
3003 : }
3004 1373975 : return y;
3005 : }
3006 :
3007 : /* assume that cyclic factors are normalized, in particular != [1] */
3008 : static GEN
3009 257426 : split_U(GEN U, GEN Sprk)
3010 : {
3011 257426 : long t = 0, k, n, l = lg(Sprk);
3012 257426 : GEN vU = cgetg(l+1, t_VEC);
3013 592451 : for (k = 1; k < l; k++)
3014 : {
3015 335023 : n = lg(sprk_get_cyc(gel(Sprk,k))) - 1; /* > 0 */
3016 335023 : gel(vU,k) = vecslice(U, t+1, t+n);
3017 335029 : t += n;
3018 : }
3019 : /* t+1 .. lg(U)-1 */
3020 257428 : n = lg(U) - t - 1; /* can be 0 */
3021 257428 : if (!n) setlg(vU,l); else gel(vU,l) = vecslice(U, t+1, t+n);
3022 257430 : return vU;
3023 : }
3024 :
3025 : static void
3026 1990309 : init_zlog_mod(zlog_S *S, GEN bid, GEN mod)
3027 : {
3028 1990309 : GEN fa2 = bid_get_fact2(bid), MOD = bid_get_MOD(bid);
3029 1990303 : S->U = bid_get_U(bid);
3030 1990302 : S->hU = lg(bid_get_cyc(bid))-1;
3031 1990301 : S->archp = bid_get_archp(bid);
3032 1990299 : S->sprk = bid_get_sprk(bid);
3033 1990297 : S->bid = bid;
3034 1990297 : if (MOD) mod = mod? gcdii(mod, MOD): MOD;
3035 1990188 : S->mod = mod;
3036 1990188 : S->P = gel(fa2,1);
3037 1990188 : S->k = gel(fa2,2);
3038 1990188 : S->no2 = lg(S->P) == lg(gel(bid_get_fact(bid),1));
3039 1990213 : }
3040 : void
3041 380084 : init_zlog(zlog_S *S, GEN bid)
3042 : {
3043 380084 : return init_zlog_mod(S, bid, NULL);
3044 : }
3045 :
3046 : /* a a t_FRAC/t_INT, reduce mod bid */
3047 : static GEN
3048 14 : Q_mod_bid(GEN bid, GEN a)
3049 : {
3050 14 : GEN xZ = gcoeff(bid_get_ideal(bid),1,1);
3051 14 : GEN b = Rg_to_Fp(a, xZ);
3052 14 : if (gsigne(a) < 0) b = subii(b, xZ);
3053 14 : return signe(b)? b: xZ;
3054 : }
3055 : /* Return decomposition of a on the CRT generators blocks attached to the
3056 : * S->sprk and sarch; sgn = sign(a, S->arch), NULL if unknown */
3057 : static GEN
3058 381275 : zlog(GEN nf, GEN a, GEN sgn, zlog_S *S)
3059 : {
3060 : long k, l;
3061 : GEN y;
3062 381275 : a = nf_to_scalar_or_basis(nf, a);
3063 381252 : switch(typ(a))
3064 : {
3065 162479 : case t_INT: break;
3066 14 : case t_FRAC: a = Q_mod_bid(S->bid, a); break;
3067 218759 : default: /* case t_COL: */
3068 : {
3069 : GEN den;
3070 218759 : check_nfelt(a, &den);
3071 218796 : if (den)
3072 : {
3073 105 : a = Q_muli_to_int(a, den);
3074 105 : a = mkmat2(mkcol2(a, den), mkcol2(gen_1, gen_m1));
3075 105 : return famat_zlog(nf, a, sgn, S);
3076 : }
3077 : }
3078 : }
3079 381164 : if (sgn)
3080 374206 : sgn = (lg(sgn) == 1)? NULL: leafcopy(sgn);
3081 : else
3082 6958 : sgn = (lg(S->archp) == 1)? NULL: nfsign_arch(nf, a, S->archp);
3083 381173 : l = lg(S->sprk);
3084 381173 : y = cgetg(sgn? l+1: l, t_COL);
3085 922272 : for (k = 1; k < l; k++)
3086 : {
3087 541130 : GEN sprk = gel(S->sprk,k);
3088 541130 : gel(y,k) = log_prk(nf, a, sprk, S->mod);
3089 : }
3090 381142 : if (sgn) gel(y,l) = Flc_to_ZC(sgn);
3091 381153 : return y;
3092 : }
3093 :
3094 : /* true nf */
3095 : GEN
3096 43813 : pr_basis_perm(GEN nf, GEN pr)
3097 : {
3098 43813 : long f = pr_get_f(pr);
3099 : GEN perm;
3100 43813 : if (f == nf_get_degree(nf)) return identity_perm(f);
3101 38164 : perm = cgetg(f+1, t_VECSMALL);
3102 38164 : perm[1] = 1;
3103 38164 : if (f > 1)
3104 : {
3105 2912 : GEN H = pr_hnf(nf,pr);
3106 : long i, k;
3107 10808 : for (i = k = 2; k <= f; i++)
3108 7896 : if (!equali1(gcoeff(H,i,i))) perm[k++] = i;
3109 : }
3110 38164 : return perm;
3111 : }
3112 :
3113 : /* \sum U[i]*y[i], U[i] ZM, y[i] ZC. We allow lg(y) > lg(U). */
3114 : static GEN
3115 1755227 : ZMV_ZCV_mul(GEN U, GEN y)
3116 : {
3117 1755227 : long i, l = lg(U);
3118 1755227 : GEN z = NULL;
3119 1755227 : if (l == 1) return cgetg(1,t_COL);
3120 4138430 : for (i = 1; i < l; i++)
3121 : {
3122 2383310 : GEN u = ZM_ZC_mul(gel(U,i), gel(y,i));
3123 2383210 : z = z? ZC_add(z, u): u;
3124 : }
3125 1755120 : return z;
3126 : }
3127 : /* A * (U[1], ..., U[d] */
3128 : static GEN
3129 518 : ZM_ZMV_mul(GEN A, GEN U)
3130 : {
3131 : long i, l;
3132 518 : GEN V = cgetg_copy(U,&l);
3133 1057 : for (i = 1; i < l; i++) gel(V,i) = ZM_mul(A,gel(U,i));
3134 518 : return V;
3135 : }
3136 :
3137 : /* a = 1 mod pr, sprk mod pr^e, e >= 1 */
3138 : static GEN
3139 402722 : sprk_log_prk1_2(GEN nf, GEN a, GEN sprk)
3140 : {
3141 402722 : GEN U1, U2, y, L2 = sprk_get_L2(sprk);
3142 402720 : sprk_get_U2(sprk, &U1,&U2);
3143 402720 : y = ZM_ZC_mul(U2, log_prk1(nf, a, lg(U2)-1, L2, sprk_get_prk(sprk)));
3144 402712 : return ZV_ZV_mod(y, sprk_get_cyc(sprk));
3145 : }
3146 : /* true nf; assume e >= 2 */
3147 : GEN
3148 105840 : sprk_log_gen_pr2(GEN nf, GEN sprk, long e)
3149 : {
3150 105840 : GEN M, G, pr = sprk_get_pr(sprk);
3151 : long i, l;
3152 105840 : if (e == 2)
3153 : {
3154 62279 : GEN L2 = sprk_get_L2(sprk), L = gel(L2,1);
3155 62279 : G = gel(L,2); l = lg(G);
3156 : }
3157 : else
3158 : {
3159 43561 : GEN perm = pr_basis_perm(nf,pr), PI = nfpow_u(nf, pr_get_gen(pr), e-1);
3160 43561 : l = lg(perm);
3161 43561 : G = cgetg(l, t_VEC);
3162 43560 : if (typ(PI) == t_INT)
3163 : { /* zk_ei_mul doesn't allow t_INT */
3164 5642 : long N = nf_get_degree(nf);
3165 5642 : gel(G,1) = addiu(PI,1);
3166 8645 : for (i = 2; i < l; i++)
3167 : {
3168 3003 : GEN z = col_ei(N, 1);
3169 3003 : gel(G,i) = z; gel(z, perm[i]) = PI;
3170 : }
3171 : }
3172 : else
3173 : {
3174 37918 : gel(G,1) = nfadd(nf, gen_1, PI);
3175 44702 : for (i = 2; i < l; i++)
3176 6783 : gel(G,i) = nfadd(nf, gen_1, zk_ei_mul(nf, PI, perm[i]));
3177 : }
3178 : }
3179 105840 : M = cgetg(l, t_MAT);
3180 234349 : for (i = 1; i < l; i++) gel(M,i) = sprk_log_prk1_2(nf, gel(G,i), sprk);
3181 105823 : return M;
3182 : }
3183 : /* Log on bid.gen of generators of P_{1,I pr^{e-1}} / P_{1,I pr^e} (I,pr) = 1,
3184 : * defined implicitly via CRT. 'ind' is the index of pr in modulus
3185 : * factorization; true nf */
3186 : GEN
3187 413817 : log_gen_pr(zlog_S *S, long ind, GEN nf, long e)
3188 : {
3189 413817 : GEN Uind = gel(S->U, ind);
3190 413817 : if (e == 1) retmkmat( gel(Uind,1) );
3191 103143 : return ZM_mul(Uind, sprk_log_gen_pr2(nf, gel(S->sprk,ind), e));
3192 : }
3193 : /* true nf */
3194 : GEN
3195 2037 : sprk_log_gen_pr(GEN nf, GEN sprk, long e)
3196 : {
3197 2037 : if (e == 1)
3198 : {
3199 0 : long n = lg(sprk_get_cyc(sprk))-1;
3200 0 : retmkmat(col_ei(n, 1));
3201 : }
3202 2037 : return sprk_log_gen_pr2(nf, sprk, e);
3203 : }
3204 : /* a = 1 mod pr */
3205 : GEN
3206 274196 : sprk_log_prk1(GEN nf, GEN a, GEN sprk)
3207 : {
3208 274196 : if (lg(sprk) == 5) return mkcol(gen_0); /* mod pr */
3209 274196 : return sprk_log_prk1_2(nf, a, sprk);
3210 : }
3211 : /* Log on bid.gen of generator of P_{1,f} / P_{1,f v[index]}
3212 : * v = vector of r1 real places */
3213 : GEN
3214 86253 : log_gen_arch(zlog_S *S, long index) { return gel(veclast(S->U), index); }
3215 :
3216 : /* compute bid.clgp: [h,cyc] or [h,cyc,gen] */
3217 : static GEN
3218 258453 : bid_grp(GEN nf, GEN U, GEN cyc, GEN g, GEN F, GEN sarch)
3219 : {
3220 258453 : GEN G, h = ZV_prod(cyc);
3221 : long c;
3222 258466 : if (!U) return mkvec2(h,cyc);
3223 258109 : c = lg(U);
3224 258109 : G = cgetg(c,t_VEC);
3225 258119 : if (c > 1)
3226 : {
3227 228033 : GEN U0, Uoo, EX = cyc_get_expo(cyc); /* exponent of bid */
3228 228033 : long i, hU = nbrows(U), nba = lg(sarch_get_cyc(sarch))-1; /* #f_oo */
3229 228041 : if (!nba) { U0 = U; Uoo = NULL; }
3230 80421 : else if (nba == hU) { U0 = NULL; Uoo = U; }
3231 : else
3232 : {
3233 71279 : U0 = rowslice(U, 1, hU-nba);
3234 71280 : Uoo = rowslice(U, hU-nba+1, hU);
3235 : }
3236 695425 : for (i = 1; i < c; i++)
3237 : {
3238 467391 : GEN t = gen_1;
3239 467391 : if (U0) t = famat_to_nf_modideal_coprime(nf, g, gel(U0,i), F, EX);
3240 467387 : if (Uoo) t = set_sign_mod_divisor(nf, ZV_to_Flv(gel(Uoo,i),2), t, sarch);
3241 467388 : gel(G,i) = t;
3242 : }
3243 : }
3244 258120 : return mkvec3(h, cyc, G);
3245 : }
3246 :
3247 : /* remove prime ideals of norm 2 with exponent 1 from factorization */
3248 : static GEN
3249 258781 : famat_strip2(GEN fa)
3250 : {
3251 258781 : GEN P = gel(fa,1), E = gel(fa,2), Q, F;
3252 258781 : long l = lg(P), i, j;
3253 258781 : Q = cgetg(l, t_COL);
3254 258783 : F = cgetg(l, t_COL);
3255 633812 : for (i = j = 1; i < l; i++)
3256 : {
3257 375025 : GEN pr = gel(P,i), e = gel(E,i);
3258 375025 : if (!absequaliu(pr_get_p(pr), 2) || itou(e) != 1 || pr_get_f(pr) != 1)
3259 : {
3260 336390 : gel(Q,j) = pr;
3261 336390 : gel(F,j) = e; j++;
3262 : }
3263 : }
3264 258787 : setlg(Q,j);
3265 258787 : setlg(F,j); return mkmat2(Q,F);
3266 : }
3267 : static int
3268 134089 : checkarchp(GEN v, long r1)
3269 : {
3270 134089 : long i, l = lg(v);
3271 134089 : pari_sp av = avma;
3272 : GEN p;
3273 134089 : if (l == 1) return 1;
3274 47154 : if (l == 2) return v[1] > 0 && v[1] <= r1;
3275 22017 : p = zero_zv(r1);
3276 66150 : for (i = 1; i < l; i++)
3277 : {
3278 44163 : long j = v[i];
3279 44163 : if (j <= 0 || j > r1 || p[j]) return gc_long(av, 0);
3280 44128 : p[j] = 1;
3281 : }
3282 21987 : return gc_long(av, 1);
3283 : }
3284 :
3285 : /* True nf. Put ideal to form [[ideal,arch]] and set fa and fa2 to its
3286 : * factorization, archp to the indices of arch places */
3287 : GEN
3288 258805 : check_mod_factored(GEN nf, GEN ideal, GEN *fa_, GEN *fa2_, GEN *archp_, GEN MOD)
3289 : {
3290 : GEN arch, x, fa, fa2, archp;
3291 : long R1;
3292 :
3293 258805 : R1 = nf_get_r1(nf);
3294 258801 : if (typ(ideal) == t_VEC && lg(ideal) == 3)
3295 : {
3296 178686 : arch = gel(ideal,2);
3297 178686 : ideal= gel(ideal,1);
3298 178686 : switch(typ(arch))
3299 : {
3300 44597 : case t_VEC:
3301 44597 : if (lg(arch) != R1+1)
3302 7 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
3303 44590 : archp = vec01_to_indices(arch);
3304 44590 : break;
3305 134089 : case t_VECSMALL:
3306 134089 : if (!checkarchp(arch, R1))
3307 35 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
3308 134049 : archp = arch;
3309 134049 : arch = indices_to_vec01(archp, R1);
3310 134045 : break;
3311 0 : default:
3312 0 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
3313 : return NULL;/*LCOV_EXCL_LINE*/
3314 : }
3315 : }
3316 : else
3317 : {
3318 80115 : arch = zerovec(R1);
3319 80112 : archp = cgetg(1, t_VECSMALL);
3320 : }
3321 258745 : if (MOD)
3322 : {
3323 214176 : if (typ(MOD) != t_INT) pari_err_TYPE("bnrinit [incorrect cycmod]", MOD);
3324 214176 : if (mpodd(MOD) && lg(archp) != 1)
3325 231 : MOD = shifti(MOD, 1); /* ensure elements of G^MOD are >> 0 */
3326 : }
3327 258746 : if (is_nf_factor(ideal))
3328 : {
3329 40271 : fa = ideal;
3330 40271 : x = factorbackprime(nf, gel(fa,1), gel(fa,2));
3331 : }
3332 : else
3333 : {
3334 218474 : fa = idealfactor(nf, ideal);
3335 218482 : x = ideal;
3336 : }
3337 258752 : if (typ(x) != t_MAT) x = idealhnf_shallow(nf, x);
3338 258725 : if (lg(x) == 1) pari_err_DOMAIN("Idealstar", "ideal","=",gen_0,x);
3339 258725 : if (typ(gcoeff(x,1,1)) != t_INT)
3340 7 : pari_err_DOMAIN("Idealstar","denominator(ideal)", "!=",gen_1,x);
3341 :
3342 258718 : fa2 = famat_strip2(fa);
3343 258723 : if (fa_ != NULL) *fa_ = fa;
3344 258723 : if (fa2_ != NULL) *fa2_ = fa2;
3345 258723 : if (fa2_ != NULL) *archp_ = archp;
3346 258723 : return mkvec2(x, arch);
3347 : }
3348 :
3349 : /* True nf. Compute [[ideal,arch], [h,[cyc],[gen]], idealfact, [liste], U]
3350 : flag may include nf_GEN | nf_INIT */
3351 : static GEN
3352 258169 : Idealstarmod_i(GEN nf, GEN ideal, long flag, GEN MOD)
3353 : {
3354 : long i, nbp;
3355 258169 : GEN y, cyc, U, u1 = NULL, fa, fa2, sprk, x_arch, x, arch, archp, E, P, sarch, gen;
3356 :
3357 258169 : x_arch = check_mod_factored(nf, ideal, &fa, &fa2, &archp, MOD);
3358 258092 : x = gel(x_arch, 1);
3359 258092 : arch = gel(x_arch, 2);
3360 :
3361 258092 : sarch = nfarchstar(nf, x, archp);
3362 258101 : P = gel(fa2,1);
3363 258101 : E = gel(fa2,2);
3364 258101 : nbp = lg(P)-1;
3365 258101 : sprk = cgetg(nbp+1,t_VEC);
3366 258107 : if (nbp)
3367 : {
3368 218837 : GEN t = (lg(gel(fa,1))==2)? NULL: x; /* beware fa != fa2 */
3369 218837 : cyc = cgetg(nbp+2,t_VEC);
3370 218820 : gen = cgetg(nbp+1,t_VEC);
3371 554557 : for (i = 1; i <= nbp; i++)
3372 : {
3373 335709 : GEN L = sprkinit(nf, gel(P,i), itou(gel(E,i)), t, MOD);
3374 335730 : gel(sprk,i) = L;
3375 335730 : gel(cyc,i) = sprk_get_cyc(L);
3376 : /* true gens are congruent to those mod x AND positive at archp */
3377 335729 : gel(gen,i) = sprk_get_gen(L);
3378 : }
3379 218848 : gel(cyc,i) = sarch_get_cyc(sarch);
3380 218848 : cyc = shallowconcat1(cyc);
3381 218853 : gen = shallowconcat1(gen);
3382 218854 : cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
3383 : }
3384 : else
3385 : {
3386 39270 : cyc = sarch_get_cyc(sarch);
3387 39270 : gen = cgetg(1,t_VEC);
3388 39270 : U = matid(lg(cyc)-1);
3389 39270 : if (flag & nf_GEN) u1 = U;
3390 : }
3391 258109 : if (MOD) cyc = ZV_snf_gcd(cyc, MOD);
3392 258089 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
3393 258118 : if (!(flag & nf_INIT)) return y;
3394 257320 : U = split_U(U, sprk);
3395 514640 : return mkvec5(mkvec2(x, arch), y, mkvec2(fa,fa2),
3396 257321 : MOD? mkvec3(sprk, sarch, MOD): mkvec2(sprk, sarch),
3397 : U);
3398 : }
3399 :
3400 : static long
3401 63 : idealHNF_norm_pval(GEN x, GEN p)
3402 : {
3403 63 : long i, v = 0, l = lg(x);
3404 175 : for (i = 1; i < l; i++) v += Z_pval(gcoeff(x,i,i), p);
3405 63 : return v;
3406 : }
3407 : static long
3408 63 : sprk_get_k(GEN sprk)
3409 : {
3410 : GEN pr, prk;
3411 63 : if (sprk_is_prime(sprk)) return 1;
3412 63 : pr = sprk_get_pr(sprk);
3413 63 : prk = sprk_get_prk(sprk);
3414 63 : return idealHNF_norm_pval(prk, pr_get_p(pr)) / pr_get_f(pr);
3415 : }
3416 : /* true nf, L a sprk */
3417 : GEN
3418 63 : sprk_to_bid(GEN nf, GEN L, long flag)
3419 : {
3420 63 : GEN y, cyc, U, u1 = NULL, fa, fa2, arch, sarch, gen, sprk;
3421 :
3422 63 : arch = zerovec(nf_get_r1(nf));
3423 63 : fa = to_famat_shallow(sprk_get_pr(L), utoipos(sprk_get_k(L)));
3424 63 : sarch = nfarchstar(nf, NULL, cgetg(1, t_VECSMALL));
3425 63 : fa2 = famat_strip2(fa);
3426 63 : sprk = mkvec(L);
3427 63 : cyc = shallowconcat(sprk_get_cyc(L), sarch_get_cyc(sarch));
3428 63 : gen = sprk_get_gen(L);
3429 63 : cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
3430 63 : y = bid_grp(nf, u1, cyc, gen, NULL, sarch);
3431 63 : if (!(flag & nf_INIT)) return y;
3432 63 : return mkvec5(mkvec2(sprk_get_prk(L), arch), y, mkvec2(fa,fa2),
3433 : mkvec2(sprk, sarch), split_U(U, sprk));
3434 : }
3435 : GEN
3436 257897 : Idealstarmod(GEN nf, GEN ideal, long flag, GEN MOD)
3437 : {
3438 257897 : pari_sp av = avma;
3439 257897 : nf = nf? checknf(nf): nfinit(pol_x(0), DEFAULTPREC);
3440 257897 : return gerepilecopy(av, Idealstarmod_i(nf, ideal, flag, MOD));
3441 : }
3442 : GEN
3443 938 : Idealstar(GEN nf, GEN ideal, long flag) { return Idealstarmod(nf, ideal, flag, NULL); }
3444 : GEN
3445 273 : Idealstarprk(GEN nf, GEN pr, long k, long flag)
3446 : {
3447 273 : pari_sp av = avma;
3448 273 : GEN z = Idealstarmod_i(nf, mkmat2(mkcol(pr),mkcols(k)), flag, NULL);
3449 273 : return gerepilecopy(av, z);
3450 : }
3451 :
3452 : /* FIXME: obsolete */
3453 : GEN
3454 0 : zidealstarinitgen(GEN nf, GEN ideal)
3455 0 : { return Idealstar(nf,ideal, nf_INIT|nf_GEN); }
3456 : GEN
3457 0 : zidealstarinit(GEN nf, GEN ideal)
3458 0 : { return Idealstar(nf,ideal, nf_INIT); }
3459 : GEN
3460 0 : zidealstar(GEN nf, GEN ideal)
3461 0 : { return Idealstar(nf,ideal, nf_GEN); }
3462 :
3463 : GEN
3464 112 : idealstarmod(GEN nf, GEN ideal, long flag, GEN MOD)
3465 : {
3466 112 : switch(flag)
3467 : {
3468 0 : case 0: return Idealstarmod(nf,ideal, nf_GEN, MOD);
3469 98 : case 1: return Idealstarmod(nf,ideal, nf_INIT, MOD);
3470 14 : case 2: return Idealstarmod(nf,ideal, nf_INIT|nf_GEN, MOD);
3471 0 : default: pari_err_FLAG("idealstar");
3472 : }
3473 : return NULL; /* LCOV_EXCL_LINE */
3474 : }
3475 : GEN
3476 0 : idealstar0(GEN nf, GEN ideal,long flag) { return idealstarmod(nf, ideal, flag, NULL); }
3477 :
3478 : void
3479 218792 : check_nfelt(GEN x, GEN *den)
3480 : {
3481 218792 : long l = lg(x), i;
3482 218792 : GEN t, d = NULL;
3483 218792 : if (typ(x) != t_COL) pari_err_TYPE("check_nfelt", x);
3484 808362 : for (i=1; i<l; i++)
3485 : {
3486 589567 : t = gel(x,i);
3487 589567 : switch (typ(t))
3488 : {
3489 589336 : case t_INT: break;
3490 231 : case t_FRAC:
3491 231 : if (!d) d = gel(t,2); else d = lcmii(d, gel(t,2));
3492 231 : break;
3493 0 : default: pari_err_TYPE("check_nfelt", x);
3494 : }
3495 : }
3496 218795 : *den = d;
3497 218795 : }
3498 :
3499 : GEN
3500 1952427 : ZV_snf_gcd(GEN x, GEN mod)
3501 4356230 : { pari_APPLY_same(gcdii(gel(x,i), mod)); }
3502 :
3503 : /* assume a true bnf and bid */
3504 : GEN
3505 227011 : ideallog_units0(GEN bnf, GEN bid, GEN MOD)
3506 : {
3507 227011 : GEN nf = bnf_get_nf(bnf), D, y, C, cyc;
3508 227011 : long j, lU = lg(bnf_get_logfu(bnf)); /* r1+r2 */
3509 : zlog_S S;
3510 227011 : init_zlog_mod(&S, bid, MOD);
3511 226985 : if (!S.hU) return zeromat(0,lU);
3512 226985 : cyc = bid_get_cyc(bid);
3513 226983 : D = nfsign_fu(bnf, bid_get_archp(bid));
3514 227001 : y = cgetg(lU, t_MAT);
3515 227000 : if ((C = bnf_build_cheapfu(bnf)))
3516 374172 : { for (j = 1; j < lU; j++) gel(y,j) = zlog(nf, gel(C,j), gel(D,j), &S); }
3517 : else
3518 : {
3519 49 : long i, l = lg(S.U), l0 = lg(S.sprk);
3520 : GEN X, U;
3521 49 : if (!(C = bnf_compactfu_mat(bnf))) bnf_build_units(bnf); /* error */
3522 49 : X = gel(C,1); U = gel(C,2);
3523 147 : for (j = 1; j < lU; j++) gel(y,j) = cgetg(l, t_COL);
3524 126 : for (i = 1; i < l0; i++)
3525 : {
3526 77 : GEN sprk = gel(S.sprk, i);
3527 77 : GEN Xi = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
3528 231 : for (j = 1; j < lU; j++)
3529 154 : gcoeff(y,i,j) = famat_zlog_pr_coprime(nf, Xi, gel(U,j), sprk, MOD);
3530 : }
3531 49 : if (l0 != l)
3532 56 : for (j = 1; j < lU; j++) gcoeff(y,l0,j) = Flc_to_ZC(gel(D,j));
3533 : }
3534 226999 : y = vec_prepend(y, zlog(nf, bnf_get_tuU(bnf), nfsign_tu(bnf, S.archp), &S));
3535 601294 : for (j = 1; j <= lU; j++)
3536 374300 : gel(y,j) = ZV_ZV_mod(ZMV_ZCV_mul(S.U, gel(y,j)), cyc);
3537 226994 : return y;
3538 : }
3539 : GEN
3540 84 : ideallog_units(GEN bnf, GEN bid)
3541 84 : { return ideallog_units0(bnf, bid, NULL); }
3542 : GEN
3543 518 : log_prk_units(GEN nf, GEN D, GEN sprk)
3544 : {
3545 518 : GEN L, Ltu = log_prk(nf, gel(D,1), sprk, NULL);
3546 518 : D = gel(D,2);
3547 518 : if (lg(D) != 3 || typ(gel(D,2)) != t_MAT) L = veclog_prk(nf, D, sprk);
3548 : else
3549 : {
3550 21 : GEN X = gel(D,1), U = gel(D,2);
3551 21 : long j, lU = lg(U);
3552 21 : X = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
3553 21 : L = cgetg(lU, t_MAT);
3554 63 : for (j = 1; j < lU; j++)
3555 42 : gel(L,j) = famat_zlog_pr_coprime(nf, X, gel(U,j), sprk, NULL);
3556 : }
3557 518 : return vec_prepend(L, Ltu);
3558 : }
3559 :
3560 : static GEN
3561 1383236 : ideallog_i(GEN nf, GEN x, zlog_S *S)
3562 : {
3563 1383236 : pari_sp av = avma;
3564 : GEN y;
3565 1383236 : if (!S->hU) return cgetg(1, t_COL);
3566 1380940 : if (typ(x) == t_MAT)
3567 1373877 : y = famat_zlog(nf, x, NULL, S);
3568 : else
3569 7063 : y = zlog(nf, x, NULL, S);
3570 1380935 : y = ZMV_ZCV_mul(S->U, y);
3571 1380931 : return gerepileupto(av, ZV_ZV_mod(y, bid_get_cyc(S->bid)));
3572 : }
3573 : GEN
3574 1389917 : ideallogmod(GEN nf, GEN x, GEN bid, GEN mod)
3575 : {
3576 : zlog_S S;
3577 1389917 : if (!nf)
3578 : {
3579 6671 : if (mod) pari_err_IMPL("Zideallogmod");
3580 6671 : return Zideallog(bid, x);
3581 : }
3582 1383246 : checkbid(bid); init_zlog_mod(&S, bid, mod);
3583 1383236 : return ideallog_i(checknf(nf), x, &S);
3584 : }
3585 : GEN
3586 13769 : ideallog(GEN nf, GEN x, GEN bid) { return ideallogmod(nf, x, bid, NULL); }
3587 :
3588 : /*************************************************************************/
3589 : /** **/
3590 : /** JOIN BID STRUCTURES, IDEAL LISTS **/
3591 : /** **/
3592 : /*************************************************************************/
3593 : /* bid1, bid2: for coprime modules m1 and m2 (without arch. part).
3594 : * Output: bid for m1 m2 */
3595 : static GEN
3596 469 : join_bid(GEN nf, GEN bid1, GEN bid2)
3597 : {
3598 469 : pari_sp av = avma;
3599 : long nbgen, l1,l2;
3600 : GEN I1,I2, G1,G2, sprk1,sprk2, cyc1,cyc2, sarch;
3601 469 : GEN sprk, fa,fa2, U, cyc, y, u1 = NULL, x, gen;
3602 :
3603 469 : I1 = bid_get_ideal(bid1);
3604 469 : I2 = bid_get_ideal(bid2);
3605 469 : if (gequal1(gcoeff(I1,1,1))) return bid2; /* frequent trivial case */
3606 259 : G1 = bid_get_grp(bid1);
3607 259 : G2 = bid_get_grp(bid2);
3608 259 : x = idealmul(nf, I1,I2);
3609 259 : fa = famat_mul_shallow(bid_get_fact(bid1), bid_get_fact(bid2));
3610 259 : fa2= famat_mul_shallow(bid_get_fact2(bid1), bid_get_fact2(bid2));
3611 259 : sprk1 = bid_get_sprk(bid1);
3612 259 : sprk2 = bid_get_sprk(bid2);
3613 259 : sprk = shallowconcat(sprk1, sprk2);
3614 :
3615 259 : cyc1 = abgrp_get_cyc(G1); l1 = lg(cyc1);
3616 259 : cyc2 = abgrp_get_cyc(G2); l2 = lg(cyc2);
3617 259 : gen = (lg(G1)>3 && lg(G2)>3)? gen_1: NULL;
3618 259 : nbgen = l1+l2-2;
3619 259 : cyc = ZV_snf_group(shallowconcat(cyc1,cyc2), &U, gen? &u1: NULL);
3620 259 : if (nbgen)
3621 : {
3622 259 : GEN U1 = bid_get_U(bid1), U2 = bid_get_U(bid2);
3623 0 : U1 = l1==1? const_vec(lg(sprk1), cgetg(1,t_MAT))
3624 259 : : ZM_ZMV_mul(vecslice(U, 1, l1-1), U1);
3625 0 : U2 = l2==1? const_vec(lg(sprk2), cgetg(1,t_MAT))
3626 259 : : ZM_ZMV_mul(vecslice(U, l1, nbgen), U2);
3627 259 : U = shallowconcat(U1, U2);
3628 : }
3629 : else
3630 0 : U = const_vec(lg(sprk), cgetg(1,t_MAT));
3631 :
3632 259 : if (gen)
3633 : {
3634 259 : GEN uv = zkchinese1init2(nf, I2, I1, x);
3635 518 : gen = shallowconcat(zkVchinese1(gel(uv,1), abgrp_get_gen(G1)),
3636 259 : zkVchinese1(gel(uv,2), abgrp_get_gen(G2)));
3637 : }
3638 259 : sarch = bid_get_sarch(bid1); /* trivial */
3639 259 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
3640 259 : x = mkvec2(x, bid_get_arch(bid1));
3641 259 : y = mkvec5(x, y, mkvec2(fa, fa2), mkvec2(sprk, sarch), U);
3642 259 : return gerepilecopy(av,y);
3643 : }
3644 :
3645 : typedef struct _ideal_data {
3646 : GEN nf, emb, L, pr, prL, sgnU, archp;
3647 : } ideal_data;
3648 :
3649 : /* z <- ( z | f(v[i])_{i=1..#v} ) */
3650 : static void
3651 757924 : concat_join(GEN *pz, GEN v, GEN (*f)(ideal_data*,GEN), ideal_data *data)
3652 : {
3653 757924 : long i, nz, lv = lg(v);
3654 : GEN z, Z;
3655 757924 : if (lv == 1) return;
3656 222904 : z = *pz; nz = lg(z)-1;
3657 222904 : *pz = Z = cgetg(lv + nz, typ(z));
3658 371665 : for (i = 1; i <=nz; i++) gel(Z,i) = gel(z,i);
3659 223329 : Z += nz;
3660 491958 : for (i = 1; i < lv; i++) gel(Z,i) = f(data, gel(v,i));
3661 : }
3662 : static GEN
3663 469 : join_idealinit(ideal_data *D, GEN x)
3664 469 : { return join_bid(D->nf, x, D->prL); }
3665 : static GEN
3666 268461 : join_ideal(ideal_data *D, GEN x)
3667 268461 : { return idealmulpowprime(D->nf, x, D->pr, D->L); }
3668 : static GEN
3669 448 : join_unit(ideal_data *D, GEN x)
3670 : {
3671 448 : GEN bid = join_idealinit(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
3672 448 : if (lg(u) != 1) v = shallowconcat(u, v);
3673 448 : return mkvec2(bid, v);
3674 : }
3675 :
3676 : GEN
3677 49 : log_prk_units_init(GEN bnf)
3678 : {
3679 49 : GEN U = bnf_has_fu(bnf);
3680 49 : if (U) U = matalgtobasis(bnf_get_nf(bnf), U);
3681 21 : else if (!(U = bnf_compactfu_mat(bnf))) (void)bnf_build_units(bnf);
3682 49 : return mkvec2(bnf_get_tuU(bnf), U);
3683 : }
3684 : /* flag & nf_GEN : generators, otherwise no
3685 : * flag &2 : units, otherwise no
3686 : * flag &4 : ideals in HNF, otherwise bid
3687 : * flag &8 : omit ideals which cannot be conductors (pr^1 with Npr=2) */
3688 : static GEN
3689 11333 : Ideallist(GEN bnf, ulong bound, long flag)
3690 : {
3691 11333 : const long do_units = flag & 2, big_id = !(flag & 4), cond = flag & 8;
3692 11333 : const long istar_flag = (flag & nf_GEN) | nf_INIT;
3693 : pari_sp av;
3694 : long i, j;
3695 11333 : GEN nf, z, p, fa, id, BOUND, U, empty = cgetg(1,t_VEC);
3696 : forprime_t S;
3697 : ideal_data ID;
3698 : GEN (*join_z)(ideal_data*, GEN);
3699 :
3700 11333 : if (do_units)
3701 : {
3702 21 : bnf = checkbnf(bnf);
3703 21 : nf = bnf_get_nf(bnf);
3704 21 : join_z = &join_unit;
3705 : }
3706 : else
3707 : {
3708 11312 : nf = checknf(bnf);
3709 11312 : join_z = big_id? &join_idealinit: &join_ideal;
3710 : }
3711 11333 : if ((long)bound <= 0) return empty;
3712 11333 : id = matid(nf_get_degree(nf));
3713 11333 : if (big_id) id = Idealstar(nf,id, istar_flag);
3714 :
3715 : /* z[i] will contain all "objects" of norm i. Depending on flag, this means
3716 : * an ideal, a bid, or a couple [bid, log(units)]. Such objects are stored
3717 : * in vectors, computed one primary component at a time; join_z
3718 : * reconstructs the global object */
3719 11333 : BOUND = utoipos(bound);
3720 11333 : z = const_vec(bound, empty);
3721 11333 : U = do_units? log_prk_units_init(bnf): NULL;
3722 11333 : gel(z,1) = mkvec(U? mkvec2(id, empty): id);
3723 11333 : ID.nf = nf;
3724 :
3725 11333 : p = cgetipos(3);
3726 11333 : u_forprime_init(&S, 2, bound);
3727 11333 : av = avma;
3728 92834 : while ((p[2] = u_forprime_next(&S)))
3729 : {
3730 81606 : if (DEBUGLEVEL>1) err_printf("%ld ",p[2]);
3731 81606 : fa = idealprimedec_limit_norm(nf, p, BOUND);
3732 163006 : for (j = 1; j < lg(fa); j++)
3733 : {
3734 81505 : GEN pr = gel(fa,j), z2 = leafcopy(z);
3735 81514 : ulong Q, q = upr_norm(pr);
3736 : long l;
3737 81514 : ID.pr = ID.prL = pr;
3738 81514 : if (cond && q == 2) { l = 2; Q = 4; } else { l = 1; Q = q; }
3739 184436 : for (; Q <= bound; l++, Q *= q) /* add pr^l */
3740 : {
3741 : ulong iQ;
3742 103041 : ID.L = utoipos(l);
3743 103044 : if (big_id) {
3744 210 : ID.prL = Idealstarprk(nf, pr, l, istar_flag);
3745 210 : if (U)
3746 189 : ID.emb = Q == 2? empty
3747 189 : : log_prk_units(nf, U, gel(bid_get_sprk(ID.prL),1));
3748 : }
3749 860802 : for (iQ = Q,i = 1; iQ <= bound; iQ += Q,i++)
3750 757880 : concat_join(&gel(z,iQ), gel(z2,i), join_z, &ID);
3751 : }
3752 : }
3753 81501 : if (gc_needed(av,1))
3754 : {
3755 18 : if(DEBUGMEM>1) pari_warn(warnmem,"Ideallist");
3756 18 : z = gerepilecopy(av, z);
3757 : }
3758 : }
3759 11333 : return z;
3760 : }
3761 : GEN
3762 63 : gideallist(GEN bnf, GEN B, long flag)
3763 : {
3764 63 : pari_sp av = avma;
3765 63 : if (typ(B) != t_INT)
3766 : {
3767 0 : B = gfloor(B);
3768 0 : if (typ(B) != t_INT) pari_err_TYPE("ideallist", B);
3769 0 : if (signe(B) < 0) B = gen_0;
3770 : }
3771 63 : if (signe(B) < 0)
3772 : {
3773 28 : if (flag != 4) pari_err_IMPL("ideallist with bid for single norm");
3774 28 : return gerepilecopy(av, ideals_by_norm(checknf(bnf), absi(B)));
3775 : }
3776 35 : if (flag < 0 || flag > 15) pari_err_FLAG("ideallist");
3777 35 : return gerepilecopy(av, Ideallist(bnf, itou(B), flag));
3778 : }
3779 : GEN
3780 11298 : ideallist0(GEN bnf, long bound, long flag)
3781 : {
3782 11298 : pari_sp av = avma;
3783 11298 : if (flag < 0 || flag > 15) pari_err_FLAG("ideallist");
3784 11298 : return gerepilecopy(av, Ideallist(bnf, bound, flag));
3785 : }
3786 : GEN
3787 10563 : ideallist(GEN bnf,long bound) { return ideallist0(bnf,bound,4); }
3788 :
3789 : /* bid = for module m (without arch. part), arch = archimedean part.
3790 : * Output: bid for [m,arch] */
3791 : static GEN
3792 42 : join_bid_arch(GEN nf, GEN bid, GEN archp)
3793 : {
3794 42 : pari_sp av = avma;
3795 : GEN G, U;
3796 42 : GEN sprk, cyc, y, u1 = NULL, x, sarch, gen;
3797 :
3798 42 : checkbid(bid);
3799 42 : G = bid_get_grp(bid);
3800 42 : x = bid_get_ideal(bid);
3801 42 : sarch = nfarchstar(nf, bid_get_ideal(bid), archp);
3802 42 : sprk = bid_get_sprk(bid);
3803 :
3804 42 : gen = (lg(G)>3)? gel(G,3): NULL;
3805 42 : cyc = diagonal_shallow(shallowconcat(gel(G,2), sarch_get_cyc(sarch)));
3806 42 : cyc = ZM_snf_group(cyc, &U, gen? &u1: NULL);
3807 42 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
3808 42 : U = split_U(U, sprk);
3809 42 : y = mkvec5(mkvec2(x, archp), y, gel(bid,3), mkvec2(sprk, sarch), U);
3810 42 : return gerepilecopy(av,y);
3811 : }
3812 : static GEN
3813 42 : join_arch(ideal_data *D, GEN x) {
3814 42 : return join_bid_arch(D->nf, x, D->archp);
3815 : }
3816 : static GEN
3817 14 : join_archunit(ideal_data *D, GEN x) {
3818 14 : GEN bid = join_arch(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
3819 14 : if (lg(u) != 1) v = shallowconcat(u, v);
3820 14 : return mkvec2(bid, v);
3821 : }
3822 :
3823 : /* L from ideallist, add archimedean part */
3824 : GEN
3825 14 : ideallistarch(GEN bnf, GEN L, GEN arch)
3826 : {
3827 : pari_sp av;
3828 14 : long i, j, l = lg(L), lz;
3829 : GEN v, z, V, nf;
3830 : ideal_data ID;
3831 : GEN (*join_z)(ideal_data*, GEN);
3832 :
3833 14 : if (typ(L) != t_VEC) pari_err_TYPE("ideallistarch",L);
3834 14 : if (l == 1) return cgetg(1,t_VEC);
3835 14 : z = gel(L,1);
3836 14 : if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
3837 14 : z = gel(z,1); /* either a bid or [bid,U] */
3838 14 : ID.archp = vec01_to_indices(arch);
3839 14 : if (lg(z) == 3)
3840 : { /* [bid,U]: do units */
3841 7 : bnf = checkbnf(bnf); nf = bnf_get_nf(bnf);
3842 7 : if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
3843 7 : ID.emb = zm_to_ZM( rowpermute(nfsign_units(bnf,NULL,1), ID.archp) );
3844 7 : join_z = &join_archunit;
3845 : }
3846 : else
3847 : {
3848 7 : join_z = &join_arch;
3849 7 : nf = checknf(bnf);
3850 : }
3851 14 : ID.nf = nf;
3852 14 : av = avma; V = cgetg(l, t_VEC);
3853 63 : for (i = 1; i < l; i++)
3854 : {
3855 49 : z = gel(L,i); lz = lg(z);
3856 49 : gel(V,i) = v = cgetg(lz,t_VEC);
3857 91 : for (j=1; j<lz; j++) gel(v,j) = join_z(&ID, gel(z,j));
3858 : }
3859 14 : return gerepilecopy(av,V);
3860 : }
|
