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 31565652 : get_tab(GEN nf, long *N)
33 : {
34 31565652 : GEN tab = (typ(nf) == t_MAT)? nf: gel(nf,9);
35 31565652 : *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 1060855619 : _mulii(GEN x, GEN y) {
41 1720511848 : return is_pm1(x)? (signe(x) < 0)? negi(y): y
42 1720415209 : : mulii(x, y);
43 : }
44 :
45 : GEN
46 21090 : tablemul_ei_ej(GEN M, long i, long j)
47 : {
48 : long N;
49 21090 : GEN tab = get_tab(M, &N);
50 21090 : 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 11347 : tablemul_ei(GEN M, GEN x, long i)
58 : {
59 : long j, k, N;
60 : GEN v, tab;
61 :
62 11347 : if (i==1) return gcopy(x);
63 11347 : tab = get_tab(M, &N);
64 11347 : if (typ(x) != t_COL) { v = zerocol(N); gel(v,i) = gcopy(x); return v; }
65 11347 : tab += (i-1)*N; v = cgetg(N+1,t_COL);
66 : /* wi . x = [ sum_j tab[k,j] x[j] ]_k */
67 77525 : for (k=1; k<=N; k++)
68 : {
69 66178 : pari_sp av = avma;
70 66178 : GEN s = gen_0;
71 469686 : for (j=1; j<=N; j++)
72 : {
73 403508 : GEN c = gcoeff(tab,k,j);
74 403508 : if (!gequal0(c)) s = gadd(s, gmul(c, gel(x,j)));
75 : }
76 66178 : gel(v,k) = gerepileupto(av,s);
77 : }
78 11347 : return v;
79 : }
80 : /* as tablemul_ei, assume x a ZV of correct length */
81 : GEN
82 21763368 : zk_ei_mul(GEN nf, GEN x, long i)
83 : {
84 : long j, k, N;
85 : GEN v, tab;
86 :
87 21763368 : if (i==1) return ZC_copy(x);
88 21763354 : tab = get_tab(nf, &N); tab += (i-1)*N;
89 21763234 : v = cgetg(N+1,t_COL);
90 161108965 : for (k=1; k<=N; k++)
91 : {
92 139347382 : pari_sp av = avma;
93 139347382 : GEN s = gen_0;
94 2115302742 : for (j=1; j<=N; j++)
95 : {
96 1976183085 : GEN c = gcoeff(tab,k,j);
97 1976183085 : if (signe(c)) s = addii(s, _mulii(c, gel(x,j)));
98 : }
99 139119657 : gel(v,k) = gerepileuptoint(av, s);
100 : }
101 21761583 : return v;
102 : }
103 :
104 : /* table of multiplication by wi in R[w1,..., wN] */
105 : GEN
106 39185 : ei_multable(GEN TAB, long i)
107 : {
108 : long k,N;
109 39185 : GEN m, tab = get_tab(TAB, &N);
110 39185 : tab += (i-1)*N;
111 39185 : m = cgetg(N+1,t_MAT);
112 153815 : for (k=1; k<=N; k++) gel(m,k) = gel(tab,k);
113 39185 : return m;
114 : }
115 :
116 : GEN
117 9658119 : zk_multable(GEN nf, GEN x)
118 : {
119 9658119 : long i, l = lg(x);
120 9658119 : GEN mul = cgetg(l,t_MAT);
121 9657916 : gel(mul,1) = x; /* assume w_1 = 1 */
122 31061227 : for (i=2; i<l; i++) gel(mul,i) = zk_ei_mul(nf,x,i);
123 9656024 : return mul;
124 : }
125 : GEN
126 2520 : multable(GEN M, GEN x)
127 : {
128 : long i, N;
129 : GEN mul;
130 2520 : 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 4238170 : zk_scalar_or_multable(GEN nf, GEN x)
143 : {
144 4238170 : long tx = typ(x);
145 4238170 : if (tx == t_MAT || tx == t_INT) return x;
146 4129249 : x = nf_to_scalar_or_basis(nf, x);
147 4129079 : return (typ(x) == t_COL)? zk_multable(nf, x): x;
148 : }
149 :
150 : GEN
151 21303 : nftrace(GEN nf, GEN x)
152 : {
153 21303 : pari_sp av = avma;
154 21303 : nf = checknf(nf);
155 21303 : x = nf_to_scalar_or_basis(nf, x);
156 21284 : x = (typ(x) == t_COL)? RgV_dotproduct(x, gel(nf_get_Tr(nf),1))
157 21305 : : gmulgu(x, nf_get_degree(nf));
158 21306 : 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 1008 : x = rnfeltabstorel(rnf, x);
166 665 : x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gtrace(x))
167 770 : : gmulgu(x, rnf_get_degree(rnf));
168 770 : return gerepileupto(av, x);
169 : }
170 :
171 : static GEN
172 35 : famatQ_to_famatZ(GEN fa)
173 : {
174 35 : GEN E, F, Q, P = gel(fa,1);
175 35 : long i, j, l = lg(P);
176 35 : if (l == 1 || RgV_is_ZV(P)) return fa;
177 7 : Q = cgetg(2*l, t_COL);
178 7 : F = cgetg(2*l, t_COL); E = gel(fa, 2);
179 35 : for (i = j = 1; i < l; i++)
180 : {
181 28 : GEN p = gel(P,i);
182 28 : if (typ(p) == t_INT)
183 14 : { gel(Q, j) = p; gel(F, j) = gel(E, i); j++; }
184 : else
185 : {
186 14 : gel(Q, j) = gel(p,1); gel(F, j) = gel(E, i); j++;
187 14 : gel(Q, j) = gel(p,2); gel(F, j) = negi(gel(E, i)); j++;
188 : }
189 : }
190 7 : setlg(Q, j); setlg(F, j); return mkmat2(Q, F);
191 : }
192 : static GEN
193 35 : famat_cba(GEN fa)
194 : {
195 35 : GEN Q, F, P = gel(fa, 1), E = gel(fa, 2);
196 35 : long i, j, lQ, l = lg(P);
197 35 : if (l == 1) return fa;
198 28 : Q = ZV_cba(P); lQ = lg(Q);
199 28 : F = cgetg(lQ, t_COL);
200 77 : for (j = 1; j < lQ; j++)
201 : {
202 49 : GEN v = gen_0, q = gel(Q,j);
203 49 : if (!equali1(q))
204 203 : for (i = 1; i < l; i++)
205 : {
206 161 : long e = Z_pval(gel(P,i), q);
207 161 : if (e) v = addii(v, muliu(gel(E,i), e));
208 : }
209 49 : gel(F, j) = v;
210 : }
211 28 : return mkmat2(Q, F);
212 : }
213 : static long
214 35 : famat_sign(GEN fa)
215 : {
216 35 : GEN P = gel(fa,1), E = gel(fa,2);
217 35 : long i, l = lg(P), s = 1;
218 126 : for (i = 1; i < l; i++)
219 91 : if (signe(gel(P,i)) < 0 && mpodd(gel(E,i))) s = -s;
220 35 : return s;
221 : }
222 : static GEN
223 35 : famat_abs(GEN fa)
224 : {
225 35 : GEN Q, P = gel(fa,1);
226 : long i, l;
227 35 : Q = cgetg_copy(P, &l);
228 126 : for (i = 1; i < l; i++) gel(Q,i) = absi_shallow(gel(P,i));
229 35 : return mkmat2(Q, gel(fa,2));
230 : }
231 :
232 : /* assume nf is a genuine nf, fa a famat */
233 : static GEN
234 35 : famat_norm(GEN nf, GEN fa)
235 : {
236 35 : pari_sp av = avma;
237 35 : GEN G, g = gel(fa,1);
238 : long i, l, s;
239 :
240 35 : G = cgetg_copy(g, &l);
241 112 : for (i = 1; i < l; i++) gel(G,i) = nfnorm(nf, gel(g,i));
242 35 : fa = mkmat2(G, gel(fa,2));
243 35 : fa = famatQ_to_famatZ(fa);
244 35 : s = famat_sign(fa);
245 35 : fa = famat_reduce(famat_abs(fa));
246 35 : fa = famat_cba(fa);
247 35 : g = factorback(fa);
248 35 : return gerepileupto(av, s < 0? gneg(g): g);
249 : }
250 : GEN
251 108125 : nfnorm(GEN nf, GEN x)
252 : {
253 108125 : pari_sp av = avma;
254 : GEN c, den;
255 : long n;
256 108125 : nf = checknf(nf);
257 108125 : n = nf_get_degree(nf);
258 108125 : if (typ(x) == t_MAT) return famat_norm(nf, x);
259 108090 : x = nf_to_scalar_or_basis(nf, x);
260 108091 : if (typ(x)!=t_COL)
261 26292 : return gerepileupto(av, gpowgs(x, n));
262 81799 : x = nf_to_scalar_or_alg(nf, Q_primitive_part(x, &c));
263 81799 : x = Q_remove_denom(x, &den);
264 81799 : x = ZX_resultant_all(nf_get_pol(nf), x, den, 0);
265 81798 : return gerepileupto(av, c ? gmul(x, gpowgs(c, n)): x);
266 : }
267 :
268 : static GEN
269 105 : to_RgX(GEN P, long vx)
270 : {
271 105 : return varn(P) == vx ? P: scalarpol_shallow(P, vx);
272 : }
273 :
274 : GEN
275 462 : rnfeltnorm(GEN rnf, GEN x)
276 : {
277 462 : pari_sp av = avma;
278 : GEN nf, pol;
279 : long v;
280 462 : checkrnf(rnf);
281 455 : v = rnf_get_varn(rnf);
282 455 : x = liftpol_shallow(rnfeltabstorel(rnf, x));
283 217 : nf = rnf_get_nf(rnf); pol = rnf_get_pol(rnf);
284 434 : x = (typ(x) == t_POL)
285 105 : ? rnfeltdown(rnf, nfX_resultant(nf,pol,to_RgX(x,v)))
286 217 : : gpowgs(x, rnf_get_degree(rnf));
287 217 : return gerepileupto(av, x);
288 : }
289 :
290 : /* x + y in nf */
291 : GEN
292 19021133 : nfadd(GEN nf, GEN x, GEN y)
293 : {
294 19021133 : pari_sp av = avma;
295 : GEN z;
296 :
297 19021133 : nf = checknf(nf);
298 19021133 : x = nf_to_scalar_or_basis(nf, x);
299 19021133 : y = nf_to_scalar_or_basis(nf, y);
300 19021134 : if (typ(x) != t_COL)
301 15273250 : { z = (typ(y) == t_COL)? RgC_Rg_add(y, x): gadd(x,y); }
302 : else
303 3747884 : { z = (typ(y) == t_COL)? RgC_add(x, y): RgC_Rg_add(x, y); }
304 19021131 : return gerepileupto(av, z);
305 : }
306 : /* x - y in nf */
307 : GEN
308 1324188 : nfsub(GEN nf, GEN x, GEN y)
309 : {
310 1324188 : pari_sp av = avma;
311 : GEN z;
312 :
313 1324188 : nf = checknf(nf);
314 1324188 : x = nf_to_scalar_or_basis(nf, x);
315 1324188 : y = nf_to_scalar_or_basis(nf, y);
316 1324188 : if (typ(x) != t_COL)
317 936985 : { z = (typ(y) == t_COL)? Rg_RgC_sub(x,y): gsub(x,y); }
318 : else
319 387203 : { z = (typ(y) == t_COL)? RgC_sub(x,y): RgC_Rg_sub(x,y); }
320 1324188 : return gerepileupto(av, z);
321 : }
322 :
323 : /* product of ZC x,y in (true) nf; ( sum_i x_i sum_j y_j m^{i,j}_k )_k */
324 : static GEN
325 4774211 : nfmuli_ZC(GEN nf, GEN x, GEN y)
326 : {
327 : long i, j, k, N;
328 4774211 : GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
329 :
330 27052038 : for (k = 1; k <= N; k++)
331 : {
332 22277876 : pari_sp av = avma;
333 22277876 : GEN s, TABi = TAB;
334 22277876 : if (k == 1)
335 4774208 : s = mulii(gel(x,1),gel(y,1));
336 : else
337 17503469 : s = addii(mulii(gel(x,1),gel(y,k)),
338 17503668 : mulii(gel(x,k),gel(y,1)));
339 185340229 : for (i=2; i<=N; i++)
340 : {
341 163066308 : GEN t, xi = gel(x,i);
342 163066308 : TABi += N;
343 163066308 : if (!signe(xi)) continue;
344 :
345 77456898 : t = NULL;
346 1013884525 : for (j=2; j<=N; j++)
347 : {
348 936428815 : GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
349 936428815 : if (!signe(c)) continue;
350 272095188 : p1 = _mulii(c, gel(y,j));
351 272101753 : t = t? addii(t, p1): p1;
352 : }
353 77455710 : if (t) s = addii(s, mulii(xi, t));
354 : }
355 22273921 : gel(v,k) = gerepileuptoint(av,s);
356 : }
357 4774162 : return v;
358 : }
359 : static int
360 50688355 : is_famat(GEN x) { return typ(x) == t_MAT && lg(x) == 3; }
361 : /* product of x and y in nf */
362 : GEN
363 25840951 : nfmul(GEN nf, GEN x, GEN y)
364 : {
365 : GEN z;
366 25840951 : pari_sp av = avma;
367 :
368 25840951 : if (x == y) return nfsqr(nf,x);
369 :
370 21804286 : nf = checknf(nf);
371 21804286 : if (is_famat(x) || is_famat(y)) return famat_mul(x, y);
372 21803976 : x = nf_to_scalar_or_basis(nf, x);
373 21803977 : y = nf_to_scalar_or_basis(nf, y);
374 21803976 : if (typ(x) != t_COL)
375 : {
376 15768928 : if (isintzero(x)) return gen_0;
377 10591931 : z = (typ(y) == t_COL)? RgC_Rg_mul(y, x): gmul(x,y); }
378 : else
379 : {
380 6035048 : if (typ(y) != t_COL)
381 : {
382 4178426 : if (isintzero(y)) return gen_0;
383 1334685 : z = RgC_Rg_mul(x, y);
384 : }
385 : else
386 : {
387 : GEN dx, dy;
388 1856622 : x = Q_remove_denom(x, &dx);
389 1856622 : y = Q_remove_denom(y, &dy);
390 1856622 : z = nfmuli_ZC(nf,x,y);
391 1856621 : dx = mul_denom(dx,dy);
392 1856621 : if (dx) z = ZC_Z_div(z, dx);
393 : }
394 : }
395 13783230 : return gerepileupto(av, z);
396 : }
397 : /* square of ZC x in nf */
398 : static GEN
399 4958777 : nfsqri_ZC(GEN nf, GEN x)
400 : {
401 : long i, j, k, N;
402 4958777 : GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
403 :
404 30977080 : for (k = 1; k <= N; k++)
405 : {
406 26018369 : pari_sp av = avma;
407 26018369 : GEN s, TABi = TAB;
408 26018369 : if (k == 1)
409 4958857 : s = sqri(gel(x,1));
410 : else
411 21059512 : s = shifti(mulii(gel(x,1),gel(x,k)), 1);
412 238436396 : for (i=2; i<=N; i++)
413 : {
414 212434132 : GEN p1, c, t, xi = gel(x,i);
415 212434132 : TABi += N;
416 212434132 : if (!signe(xi)) continue;
417 :
418 73571755 : c = gcoeff(TABi, k, i);
419 73571755 : t = signe(c)? _mulii(c,xi): NULL;
420 668473273 : for (j=i+1; j<=N; j++)
421 : {
422 594901548 : c = gcoeff(TABi, k, j);
423 594901548 : if (!signe(c)) continue;
424 232560595 : p1 = _mulii(c, shifti(gel(x,j),1));
425 232564531 : t = t? addii(t, p1): p1;
426 : }
427 73571725 : if (t) s = addii(s, mulii(xi, t));
428 : }
429 26002264 : gel(v,k) = gerepileuptoint(av,s);
430 : }
431 4958711 : return v;
432 : }
433 : /* square of x in nf */
434 : GEN
435 6083566 : nfsqr(GEN nf, GEN x)
436 : {
437 6083566 : pari_sp av = avma;
438 : GEN z;
439 :
440 6083566 : nf = checknf(nf);
441 6083567 : if (is_famat(x)) return famat_sqr(x);
442 6083568 : x = nf_to_scalar_or_basis(nf, x);
443 6083568 : if (typ(x) != t_COL) z = gsqr(x);
444 : else
445 : {
446 : GEN dx;
447 261374 : x = Q_remove_denom(x, &dx);
448 261373 : z = nfsqri_ZC(nf,x);
449 261372 : if (dx) z = RgC_Rg_div(z, sqri(dx));
450 : }
451 6083565 : return gerepileupto(av, z);
452 : }
453 :
454 : /* x a ZC, v a t_COL of ZC/Z */
455 : GEN
456 204316 : zkC_multable_mul(GEN v, GEN x)
457 : {
458 204316 : long i, l = lg(v);
459 204316 : GEN y = cgetg(l, t_COL);
460 793681 : for (i = 1; i < l; i++)
461 : {
462 589365 : GEN c = gel(v,i);
463 589365 : if (typ(c)!=t_COL) {
464 0 : if (!isintzero(c)) c = ZC_Z_mul(gel(x,1), c);
465 : } else {
466 589365 : c = ZM_ZC_mul(x,c);
467 589365 : if (ZV_isscalar(c)) c = gel(c,1);
468 : }
469 589365 : gel(y,i) = c;
470 : }
471 204316 : return y;
472 : }
473 :
474 : GEN
475 55071 : nfC_multable_mul(GEN v, GEN x)
476 : {
477 55071 : long i, l = lg(v);
478 55071 : GEN y = cgetg(l, t_COL);
479 375003 : for (i = 1; i < l; i++)
480 : {
481 319932 : GEN c = gel(v,i);
482 319932 : if (typ(c)!=t_COL) {
483 266687 : if (!isintzero(c)) c = RgC_Rg_mul(gel(x,1), c);
484 : } else {
485 53245 : c = RgM_RgC_mul(x,c);
486 53245 : if (QV_isscalar(c)) c = gel(c,1);
487 : }
488 319932 : gel(y,i) = c;
489 : }
490 55071 : return y;
491 : }
492 :
493 : GEN
494 192434 : nfC_nf_mul(GEN nf, GEN v, GEN x)
495 : {
496 : long tx;
497 : GEN y;
498 :
499 192434 : x = nf_to_scalar_or_basis(nf, x);
500 192434 : tx = typ(x);
501 192434 : if (tx != t_COL)
502 : {
503 : long l, i;
504 145489 : if (tx == t_INT)
505 : {
506 136501 : long s = signe(x);
507 136501 : if (!s) return zerocol(lg(v)-1);
508 129258 : if (is_pm1(x)) return s > 0? leafcopy(v): RgC_neg(v);
509 : }
510 47166 : l = lg(v); y = cgetg(l, t_COL);
511 340753 : for (i=1; i < l; i++)
512 : {
513 293587 : GEN c = gel(v,i);
514 293587 : if (typ(c) != t_COL) c = gmul(c, x); else c = RgC_Rg_mul(c, x);
515 293587 : gel(y,i) = c;
516 : }
517 47166 : return y;
518 : }
519 : else
520 : {
521 : GEN dx;
522 46945 : x = zk_multable(nf, Q_remove_denom(x,&dx));
523 46945 : y = nfC_multable_mul(v, x);
524 46945 : return dx? RgC_Rg_div(y, dx): y;
525 : }
526 : }
527 : static GEN
528 10618 : mulbytab(GEN M, GEN c)
529 10618 : { return typ(c) == t_COL? RgM_RgC_mul(M,c): RgC_Rg_mul(gel(M,1), c); }
530 : GEN
531 2520 : tablemulvec(GEN M, GEN x, GEN v)
532 : {
533 : long l, i;
534 : GEN y;
535 :
536 2520 : if (typ(x) == t_COL && RgV_isscalar(x))
537 : {
538 0 : x = gel(x,1);
539 0 : return typ(v) == t_POL? RgX_Rg_mul(v,x): RgV_Rg_mul(v,x);
540 : }
541 2520 : x = multable(M, x); /* multiplication table by x */
542 2520 : y = cgetg_copy(v, &l);
543 2520 : if (typ(v) == t_POL)
544 : {
545 2520 : y[1] = v[1];
546 13138 : for (i=2; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
547 2520 : y = normalizepol(y);
548 : }
549 : else
550 : {
551 0 : for (i=1; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
552 : }
553 2520 : return y;
554 : }
555 :
556 : GEN
557 1258075 : zkmultable_capZ(GEN mx) { return Q_denom(zkmultable_inv(mx)); }
558 : GEN
559 1419951 : zkmultable_inv(GEN mx) { return ZM_gauss(mx, col_ei(lg(mx)-1,1)); }
560 : /* nf a true nf, x a ZC */
561 : GEN
562 161876 : zk_inv(GEN nf, GEN x) { return zkmultable_inv(zk_multable(nf,x)); }
563 :
564 : /* inverse of x in nf */
565 : GEN
566 83209 : nfinv(GEN nf, GEN x)
567 : {
568 83209 : pari_sp av = avma;
569 : GEN z;
570 :
571 83209 : nf = checknf(nf);
572 83209 : if (is_famat(x)) return famat_inv(x);
573 83209 : x = nf_to_scalar_or_basis(nf, x);
574 83209 : if (typ(x) == t_COL)
575 : {
576 : GEN d;
577 35836 : x = Q_remove_denom(x, &d);
578 35836 : z = zk_inv(nf, x);
579 35836 : if (d) z = RgC_Rg_mul(z, d);
580 : }
581 : else
582 47373 : z = ginv(x);
583 83209 : return gerepileupto(av, z);
584 : }
585 :
586 : /* quotient of x and y in nf */
587 : GEN
588 36171 : nfdiv(GEN nf, GEN x, GEN y)
589 : {
590 36171 : pari_sp av = avma;
591 : GEN z;
592 :
593 36171 : nf = checknf(nf);
594 36171 : if (is_famat(x) || is_famat(y)) return famat_div(x,y);
595 36080 : y = nf_to_scalar_or_basis(nf, y);
596 36080 : if (typ(y) != t_COL)
597 : {
598 21805 : x = nf_to_scalar_or_basis(nf, x);
599 21805 : z = (typ(x) == t_COL)? RgC_Rg_div(x, y): gdiv(x,y);
600 : }
601 : else
602 : {
603 : GEN d;
604 14275 : y = Q_remove_denom(y, &d);
605 14275 : z = nfmul(nf, x, zk_inv(nf,y));
606 14275 : if (d) z = typ(z) == t_COL? RgC_Rg_mul(z, d): gmul(z, d);
607 : }
608 36080 : return gerepileupto(av, z);
609 : }
610 :
611 : /* product of INTEGERS (t_INT or ZC) x and y in (true) nf */
612 : GEN
613 4252109 : nfmuli(GEN nf, GEN x, GEN y)
614 : {
615 4252109 : if (typ(x) == t_INT) return (typ(y) == t_COL)? ZC_Z_mul(y, x): mulii(x,y);
616 3151146 : if (typ(y) == t_INT) return ZC_Z_mul(x, y);
617 2917565 : return nfmuli_ZC(nf, x, y);
618 : }
619 : GEN
620 4697349 : nfsqri(GEN nf, GEN x)
621 4697349 : { return (typ(x) == t_INT)? sqri(x): nfsqri_ZC(nf, x); }
622 :
623 : /* both x and y are RgV */
624 : GEN
625 0 : tablemul(GEN TAB, GEN x, GEN y)
626 : {
627 : long i, j, k, N;
628 : GEN s, v;
629 0 : if (typ(x) != t_COL) return gmul(x, y);
630 0 : if (typ(y) != t_COL) return gmul(y, x);
631 0 : N = lg(x)-1;
632 0 : v = cgetg(N+1,t_COL);
633 0 : for (k=1; k<=N; k++)
634 : {
635 0 : pari_sp av = avma;
636 0 : GEN TABi = TAB;
637 0 : if (k == 1)
638 0 : s = gmul(gel(x,1),gel(y,1));
639 : else
640 0 : s = gadd(gmul(gel(x,1),gel(y,k)),
641 0 : gmul(gel(x,k),gel(y,1)));
642 0 : for (i=2; i<=N; i++)
643 : {
644 0 : GEN t, xi = gel(x,i);
645 0 : TABi += N;
646 0 : if (gequal0(xi)) continue;
647 :
648 0 : t = NULL;
649 0 : for (j=2; j<=N; j++)
650 : {
651 0 : GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
652 0 : if (gequal0(c)) continue;
653 0 : p1 = gmul(c, gel(y,j));
654 0 : t = t? gadd(t, p1): p1;
655 : }
656 0 : if (t) s = gadd(s, gmul(xi, t));
657 : }
658 0 : gel(v,k) = gerepileupto(av,s);
659 : }
660 0 : return v;
661 : }
662 : GEN
663 48551 : tablesqr(GEN TAB, GEN x)
664 : {
665 : long i, j, k, N;
666 : GEN s, v;
667 :
668 48551 : if (typ(x) != t_COL) return gsqr(x);
669 48551 : N = lg(x)-1;
670 48551 : v = cgetg(N+1,t_COL);
671 :
672 347071 : for (k=1; k<=N; k++)
673 : {
674 298520 : pari_sp av = avma;
675 298520 : GEN TABi = TAB;
676 298520 : if (k == 1)
677 48551 : s = gsqr(gel(x,1));
678 : else
679 249969 : s = gmul2n(gmul(gel(x,1),gel(x,k)), 1);
680 1904206 : for (i=2; i<=N; i++)
681 : {
682 1605686 : GEN p1, c, t, xi = gel(x,i);
683 1605686 : TABi += N;
684 1605686 : if (gequal0(xi)) continue;
685 :
686 414470 : c = gcoeff(TABi, k, i);
687 414470 : t = !gequal0(c)? gmul(c,xi): NULL;
688 1662752 : for (j=i+1; j<=N; j++)
689 : {
690 1248282 : c = gcoeff(TABi, k, j);
691 1248282 : if (gequal0(c)) continue;
692 641564 : p1 = gmul(gmul2n(c,1), gel(x,j));
693 641564 : t = t? gadd(t, p1): p1;
694 : }
695 414470 : if (t) s = gadd(s, gmul(xi, t));
696 : }
697 298520 : gel(v,k) = gerepileupto(av,s);
698 : }
699 48551 : return v;
700 : }
701 :
702 : static GEN
703 204753 : _mul(void *data, GEN x, GEN y) { return nfmuli((GEN)data,x,y); }
704 : static GEN
705 734582 : _sqr(void *data, GEN x) { return nfsqri((GEN)data,x); }
706 :
707 : /* Compute z^n in nf, left-shift binary powering */
708 : GEN
709 825980 : nfpow(GEN nf, GEN z, GEN n)
710 : {
711 825980 : pari_sp av = avma;
712 : long s;
713 : GEN x, cx;
714 :
715 825980 : if (typ(n)!=t_INT) pari_err_TYPE("nfpow",n);
716 825980 : nf = checknf(nf);
717 825981 : s = signe(n); if (!s) return gen_1;
718 825981 : if (is_famat(z)) return famat_pow(z, n);
719 765339 : x = nf_to_scalar_or_basis(nf, z);
720 765341 : if (typ(x) != t_COL) return powgi(x,n);
721 647232 : if (s < 0)
722 : { /* simplified nfinv */
723 : GEN d;
724 43474 : x = Q_remove_denom(x, &d);
725 43473 : x = zk_inv(nf, x);
726 43474 : x = primitive_part(x, &cx);
727 43474 : cx = mul_content(cx, d);
728 43474 : n = negi(n);
729 : }
730 : else
731 603758 : x = primitive_part(x, &cx);
732 647223 : x = gen_pow_i(x, n, (void*)nf, _sqr, _mul);
733 647229 : if (cx)
734 53609 : x = gerepileupto(av, gmul(x, powgi(cx, n)));
735 : else
736 593620 : x = gerepilecopy(av, x);
737 647237 : return x;
738 : }
739 : /* Compute z^n in nf, left-shift binary powering */
740 : GEN
741 222072 : nfpow_u(GEN nf, GEN z, ulong n)
742 : {
743 222072 : pari_sp av = avma;
744 : GEN x, cx;
745 :
746 222072 : if (!n) return gen_1;
747 222072 : x = nf_to_scalar_or_basis(nf, z);
748 222072 : if (typ(x) != t_COL) return gpowgs(x,n);
749 185842 : x = primitive_part(x, &cx);
750 185842 : x = gen_powu_i(x, n, (void*)nf, _sqr, _mul);
751 185842 : if (cx)
752 : {
753 26332 : x = gmul(x, powgi(cx, utoipos(n)));
754 26332 : return gerepileupto(av,x);
755 : }
756 159510 : return gerepilecopy(av, x);
757 : }
758 :
759 : long
760 63 : nfissquare(GEN nf, GEN z, GEN *px)
761 : {
762 63 : pari_sp av = avma;
763 63 : long v = fetch_var_higher();
764 : GEN R;
765 63 : nf = checknf(nf);
766 63 : if (nf_get_degree(nf) == 1)
767 : {
768 21 : z = algtobasis(nf, z);
769 21 : if (!issquareall(gel(z,1), px)) return gc_long(av, 0);
770 14 : if (px) *px = gerepileupto(av, *px); else set_avma(av);
771 14 : return 1;
772 : }
773 42 : z = nf_to_scalar_or_alg(nf, z);
774 42 : R = nfroots(nf, deg2pol_shallow(gen_m1, gen_0, z, v));
775 42 : delete_var(); if (lg(R) == 1) return gc_long(av, 0);
776 35 : if (px) *px = gerepilecopy(av, nf_to_scalar_or_basis(nf, gel(R,1)));
777 14 : else set_avma(av);
778 35 : return 1;
779 : }
780 :
781 : long
782 8083 : nfispower(GEN nf, GEN z, long n, GEN *px)
783 : {
784 8083 : pari_sp av = avma;
785 8083 : long v = fetch_var_higher();
786 : GEN R;
787 8083 : nf = checknf(nf);
788 8083 : if (nf_get_degree(nf) == 1)
789 : {
790 329 : z = algtobasis(nf, z);
791 329 : if (!ispower(gel(z,1), stoi(n), px)) return gc_long(av, 0);
792 147 : if (px) *px = gerepileupto(av, *px); else set_avma(av);
793 147 : return 1;
794 : }
795 7754 : if (n <= 0)
796 0 : pari_err_DOMAIN("nfeltispower","exponent","<=",gen_0,stoi(n));
797 7754 : z = nf_to_scalar_or_alg(nf, z);
798 7754 : if (n==1)
799 : {
800 0 : if (px) *px = gerepilecopy(av, z);
801 0 : return 1;
802 : }
803 7754 : R = nfroots(nf, gsub(pol_xn(n, v), z));
804 7754 : delete_var(); if (lg(R) == 1) return gc_long(av, 0);
805 3157 : if (px) *px = gerepilecopy(av, nf_to_scalar_or_basis(nf, gel(R,1)));
806 3143 : else set_avma(av);
807 3157 : return 1;
808 : }
809 :
810 : static GEN
811 56 : idmulred(void *nf, GEN x, GEN y) { return idealmulred((GEN) nf, x, y); }
812 : static GEN
813 413 : idpowred(void *nf, GEN x, GEN n) { return idealpowred((GEN) nf, x, n); }
814 : static GEN
815 70354 : idmul(void *nf, GEN x, GEN y) { return idealmul((GEN) nf, x, y); }
816 : static GEN
817 86116 : idpow(void *nf, GEN x, GEN n) { return idealpow((GEN) nf, x, n); }
818 : GEN
819 85147 : idealfactorback(GEN nf, GEN L, GEN e, int red)
820 : {
821 85147 : nf = checknf(nf);
822 85147 : if (red) return gen_factorback(L, e, (void*)nf, &idmulred, &idpowred, NULL);
823 84790 : if (!e && typ(L) == t_MAT && lg(L) == 3) { e = gel(L,2); L = gel(L,1); }
824 84790 : if (is_vec_t(typ(L)) && RgV_is_prV(L))
825 : { /* don't use gen_factorback since *= pr^v can be done more efficiently */
826 64374 : pari_sp av = avma;
827 64374 : long i, l = lg(L);
828 : GEN a;
829 64374 : if (!e) e = const_vec(l-1, gen_1);
830 61518 : else switch(typ(e))
831 : {
832 7264 : case t_VECSMALL: e = zv_to_ZV(e); break;
833 54254 : case t_VEC: case t_COL:
834 54254 : if (!RgV_is_ZV(e))
835 0 : pari_err_TYPE("factorback [not an exponent vector]", e);
836 54254 : break;
837 0 : default: pari_err_TYPE("idealfactorback", e);
838 : }
839 64374 : if (l != lg(e))
840 0 : pari_err_TYPE("factorback [not an exponent vector]", e);
841 64374 : if (l == 1 || ZV_equal0(e)) return gc_const(av, gen_1);
842 22139 : a = idealpow(nf, gel(L,1), gel(e,1));
843 241803 : for (i = 2; i < l; i++)
844 219664 : if (signe(gel(e,i))) a = idealmulpowprime(nf, a, gel(L,i), gel(e,i));
845 22139 : return gerepileupto(av, a);
846 : }
847 20417 : return gen_factorback(L, e, (void*)nf, &idmul, &idpow, NULL);
848 : }
849 : static GEN
850 325128 : eltmul(void *nf, GEN x, GEN y) { return nfmul((GEN) nf, x, y); }
851 : static GEN
852 462144 : eltpow(void *nf, GEN x, GEN n) { return nfpow((GEN) nf, x, n); }
853 : GEN
854 265082 : nffactorback(GEN nf, GEN L, GEN e)
855 265082 : { return gen_factorback(L, e, (void*)checknf(nf), &eltmul, &eltpow, NULL); }
856 :
857 : static GEN
858 3020596 : _nf_red(void *E, GEN x) { (void)E; return gcopy(x); }
859 :
860 : static GEN
861 12502657 : _nf_add(void *E, GEN x, GEN y) { return nfadd((GEN)E,x,y); }
862 :
863 : static GEN
864 735590 : _nf_neg(void *E, GEN x) { (void)E; return gneg(x); }
865 :
866 : static GEN
867 14987696 : _nf_mul(void *E, GEN x, GEN y) { return nfmul((GEN)E,x,y); }
868 :
869 : static GEN
870 51201 : _nf_inv(void *E, GEN x) { return nfinv((GEN)E,x); }
871 :
872 : static GEN
873 10358 : _nf_s(void *E, long x) { (void)E; return stoi(x); }
874 :
875 : static const struct bb_field nf_field={_nf_red,_nf_add,_nf_mul,_nf_neg,
876 : _nf_inv,&gequal0,_nf_s };
877 :
878 219721 : const struct bb_field *get_nf_field(void **E, GEN nf)
879 219721 : { *E = (void*)nf; return &nf_field; }
880 :
881 : GEN
882 14 : nfM_det(GEN nf, GEN M)
883 : {
884 : void *E;
885 14 : const struct bb_field *S = get_nf_field(&E, nf);
886 14 : return gen_det(M, E, S);
887 : }
888 : GEN
889 10344 : nfM_inv(GEN nf, GEN M)
890 : {
891 : void *E;
892 10344 : const struct bb_field *S = get_nf_field(&E, nf);
893 10344 : return gen_Gauss(M, matid(lg(M)-1), E, S);
894 : }
895 :
896 : GEN
897 0 : nfM_ker(GEN nf, GEN M)
898 : {
899 : void *E;
900 0 : const struct bb_field *S = get_nf_field(&E, nf);
901 0 : return gen_ker(M, 0, E, S);
902 : }
903 :
904 : GEN
905 9924 : nfM_mul(GEN nf, GEN A, GEN B)
906 : {
907 : void *E;
908 9924 : const struct bb_field *S = get_nf_field(&E, nf);
909 9924 : return gen_matmul(A, B, E, S);
910 : }
911 : GEN
912 199439 : nfM_nfC_mul(GEN nf, GEN A, GEN B)
913 : {
914 : void *E;
915 199439 : const struct bb_field *S = get_nf_field(&E, nf);
916 199439 : return gen_matcolmul(A, B, E, S);
917 : }
918 :
919 : /* valuation of integral x (ZV), with resp. to prime ideal pr */
920 : long
921 25153448 : ZC_nfvalrem(GEN x, GEN pr, GEN *newx)
922 : {
923 25153448 : pari_sp av = avma;
924 : long i, v, l;
925 25153448 : GEN r, y, p = pr_get_p(pr), mul = pr_get_tau(pr);
926 :
927 : /* p inert */
928 25153759 : if (typ(mul) == t_INT) return newx? ZV_pvalrem(x, p, newx):ZV_pval(x, p);
929 24158644 : y = cgetg_copy(x, &l); /* will hold the new x */
930 24159097 : x = leafcopy(x);
931 39433952 : for(v=0;; v++)
932 : {
933 160103165 : for (i=1; i<l; i++)
934 : { /* is (x.b)[i] divisible by p ? */
935 144822544 : gel(y,i) = dvmdii(ZMrow_ZC_mul(mul,x,i),p,&r);
936 144824126 : if (r != gen_0) { if (newx) *newx = x; return v; }
937 : }
938 15280621 : swap(x, y);
939 15280621 : if (!newx && (v & 0xf) == 0xf) v += pr_get_e(pr) * ZV_pvalrem(x, p, &x);
940 15280621 : if (gc_needed(av,1))
941 : {
942 0 : if(DEBUGMEM>1) pari_warn(warnmem,"ZC_nfvalrem, v >= %ld", v);
943 0 : gerepileall(av, 2, &x, &y);
944 : }
945 : }
946 : }
947 : long
948 20892419 : ZC_nfval(GEN x, GEN P)
949 20892419 : { return ZC_nfvalrem(x, P, NULL); }
950 :
951 : /* v_P(x) != 0, x a ZV. Simpler version of ZC_nfvalrem */
952 : int
953 1243477 : ZC_prdvd(GEN x, GEN P)
954 : {
955 1243477 : pari_sp av = avma;
956 : long i, l;
957 1243477 : GEN p = pr_get_p(P), mul = pr_get_tau(P);
958 1243481 : if (typ(mul) == t_INT) return ZV_Z_dvd(x, p);
959 1243005 : l = lg(x);
960 5041664 : for (i=1; i<l; i++)
961 4527947 : if (!dvdii(ZMrow_ZC_mul(mul,x,i), p)) return gc_bool(av,0);
962 513717 : return gc_bool(av,1);
963 : }
964 :
965 : int
966 357 : pr_equal(GEN P, GEN Q)
967 : {
968 357 : GEN gQ, p = pr_get_p(P);
969 357 : long e = pr_get_e(P), f = pr_get_f(P), n;
970 357 : if (!equalii(p, pr_get_p(Q)) || e != pr_get_e(Q) || f != pr_get_f(Q))
971 336 : return 0;
972 21 : gQ = pr_get_gen(Q); n = lg(gQ)-1;
973 21 : if (2*e*f > n) return 1; /* room for only one such pr */
974 14 : return ZV_equal(pr_get_gen(P), gQ) || ZC_prdvd(gQ, P);
975 : }
976 :
977 : GEN
978 420721 : famat_nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
979 : {
980 420721 : pari_sp av = avma;
981 420721 : GEN P = gel(x,1), E = gel(x,2), V = gen_0, y = NULL;
982 420721 : long l = lg(P), simplify = 0, i;
983 420721 : if (py) { *py = gen_1; y = cgetg(l, t_COL); }
984 :
985 2270533 : for (i = 1; i < l; i++)
986 : {
987 1849812 : GEN e = gel(E,i);
988 : long v;
989 1849812 : if (!signe(e))
990 : {
991 7 : if (py) gel(y,i) = gen_1;
992 7 : simplify = 1; continue;
993 : }
994 1849805 : v = nfvalrem(nf, gel(P,i), pr, py? &gel(y,i): NULL);
995 1849805 : if (v == LONG_MAX) { set_avma(av); if (py) *py = gen_0; return mkoo(); }
996 1849805 : V = addmulii(V, stoi(v), e);
997 : }
998 420721 : if (!py) V = gerepileuptoint(av, V);
999 : else
1000 : {
1001 42 : y = mkmat2(y, gel(x,2));
1002 42 : if (simplify) y = famat_remove_trivial(y);
1003 42 : gerepileall(av, 2, &V, &y); *py = y;
1004 : }
1005 420721 : return V;
1006 : }
1007 : long
1008 3313695 : nfval(GEN nf, GEN x, GEN pr)
1009 : {
1010 3313695 : pari_sp av = avma;
1011 : long w, e;
1012 : GEN cx, p;
1013 :
1014 3313695 : if (gequal0(x)) return LONG_MAX;
1015 3311831 : nf = checknf(nf);
1016 3311830 : checkprid(pr);
1017 3311837 : p = pr_get_p(pr);
1018 3311830 : e = pr_get_e(pr);
1019 3311830 : x = nf_to_scalar_or_basis(nf, x);
1020 3311736 : if (typ(x) != t_COL) return e*Q_pval(x,p);
1021 1504922 : x = Q_primitive_part(x, &cx);
1022 1504979 : w = ZC_nfval(x,pr);
1023 1504914 : if (cx) w += e*Q_pval(cx,p);
1024 1504927 : return gc_long(av,w);
1025 : }
1026 :
1027 : /* want to write p^v = uniformizer^(e*v) * z^v, z coprime to pr */
1028 : /* z := tau^e / p^(e-1), algebraic integer coprime to pr; return z^v */
1029 : static GEN
1030 943320 : powp(GEN nf, GEN pr, long v)
1031 : {
1032 : GEN b, z;
1033 : long e;
1034 943320 : if (!v) return gen_1;
1035 416304 : b = pr_get_tau(pr);
1036 416304 : if (typ(b) == t_INT) return gen_1;
1037 115486 : e = pr_get_e(pr);
1038 115486 : z = gel(b,1);
1039 115486 : if (e != 1) z = gdiv(nfpow_u(nf, z, e), powiu(pr_get_p(pr),e-1));
1040 115486 : if (v < 0) { v = -v; z = nfinv(nf, z); }
1041 115486 : if (v != 1) z = nfpow_u(nf, z, v);
1042 115486 : return z;
1043 : }
1044 : long
1045 3642851 : nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
1046 : {
1047 3642851 : pari_sp av = avma;
1048 : long w, e;
1049 : GEN cx, p, t;
1050 :
1051 3642851 : if (!py) return nfval(nf,x,pr);
1052 1779606 : if (gequal0(x)) { *py = gen_0; return LONG_MAX; }
1053 1779548 : nf = checknf(nf);
1054 1779550 : checkprid(pr);
1055 1779549 : p = pr_get_p(pr);
1056 1779549 : e = pr_get_e(pr);
1057 1779549 : x = nf_to_scalar_or_basis(nf, x);
1058 1779551 : if (typ(x) != t_COL) {
1059 532882 : w = Q_pvalrem(x,p, py);
1060 532882 : if (!w) { *py = gerepilecopy(av, x); return 0; }
1061 324611 : *py = gerepileupto(av, gmul(powp(nf, pr, w), *py));
1062 324611 : return e*w;
1063 : }
1064 1246669 : x = Q_primitive_part(x, &cx);
1065 1246665 : w = ZC_nfvalrem(x,pr, py);
1066 1246647 : if (cx)
1067 : {
1068 618709 : long v = Q_pvalrem(cx,p, &t);
1069 618709 : *py = nfmul(nf, *py, gmul(powp(nf,pr,v), t));
1070 618709 : *py = gerepileupto(av, *py);
1071 618708 : w += e*v;
1072 : }
1073 : else
1074 627938 : *py = gerepilecopy(av, *py);
1075 1246670 : return w;
1076 : }
1077 : GEN
1078 15015 : gpnfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
1079 : {
1080 : long v;
1081 15015 : if (is_famat(x)) return famat_nfvalrem(nf, x, pr, py);
1082 15008 : v = nfvalrem(nf,x,pr,py);
1083 15008 : return v == LONG_MAX? mkoo(): stoi(v);
1084 : }
1085 :
1086 : /* true nf */
1087 : GEN
1088 108115 : coltoalg(GEN nf, GEN x)
1089 : {
1090 108115 : return mkpolmod( nf_to_scalar_or_alg(nf, x), nf_get_pol(nf) );
1091 : }
1092 :
1093 : GEN
1094 160738 : basistoalg(GEN nf, GEN x)
1095 : {
1096 : GEN T;
1097 :
1098 160738 : nf = checknf(nf);
1099 160738 : switch(typ(x))
1100 : {
1101 102004 : case t_COL: {
1102 102004 : pari_sp av = avma;
1103 102004 : return gerepilecopy(av, coltoalg(nf, x));
1104 : }
1105 33390 : case t_POLMOD:
1106 33390 : T = nf_get_pol(nf);
1107 33390 : if (!RgX_equal_var(T,gel(x,1)))
1108 0 : pari_err_MODULUS("basistoalg", T,gel(x,1));
1109 33390 : return gcopy(x);
1110 2156 : case t_POL:
1111 2156 : T = nf_get_pol(nf);
1112 2156 : if (varn(T) != varn(x)) pari_err_VAR("basistoalg",x,T);
1113 2156 : retmkpolmod(RgX_rem(x, T), ZX_copy(T));
1114 23188 : case t_INT:
1115 : case t_FRAC:
1116 23188 : T = nf_get_pol(nf);
1117 23188 : retmkpolmod(gcopy(x), ZX_copy(T));
1118 0 : default:
1119 0 : pari_err_TYPE("basistoalg",x);
1120 : return NULL; /* LCOV_EXCL_LINE */
1121 : }
1122 : }
1123 :
1124 : /* true nf, x a t_POL */
1125 : static GEN
1126 4313652 : pol_to_scalar_or_basis(GEN nf, GEN x)
1127 : {
1128 4313652 : GEN T = nf_get_pol(nf);
1129 4313654 : long l = lg(x);
1130 4313654 : if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_basis", x,T);
1131 4313550 : if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
1132 4313550 : if (l == 2) return gen_0;
1133 3300349 : if (l == 3)
1134 : {
1135 802730 : x = gel(x,2);
1136 802730 : if (!is_rational_t(typ(x))) pari_err_TYPE("nf_to_scalar_or_basis",x);
1137 802723 : return x;
1138 : }
1139 2497619 : return poltobasis(nf,x);
1140 : }
1141 : /* Assume nf is a genuine nf. */
1142 : GEN
1143 122689746 : nf_to_scalar_or_basis(GEN nf, GEN x)
1144 : {
1145 122689746 : switch(typ(x))
1146 : {
1147 80142135 : case t_INT: case t_FRAC:
1148 80142135 : return x;
1149 270651 : case t_POLMOD:
1150 270651 : x = checknfelt_mod(nf,x,"nf_to_scalar_or_basis");
1151 270527 : switch(typ(x))
1152 : {
1153 38507 : case t_INT: case t_FRAC: return x;
1154 232020 : case t_POL: return pol_to_scalar_or_basis(nf,x);
1155 : }
1156 0 : break;
1157 4081646 : case t_POL: return pol_to_scalar_or_basis(nf,x);
1158 38198241 : case t_COL:
1159 38198241 : if (lg(x)-1 != nf_get_degree(nf)) break;
1160 38198023 : return QV_isscalar(x)? gel(x,1): x;
1161 : }
1162 96 : pari_err_TYPE("nf_to_scalar_or_basis",x);
1163 : return NULL; /* LCOV_EXCL_LINE */
1164 : }
1165 : /* Let x be a polynomial with coefficients in Q or nf. Return the same
1166 : * polynomial with coefficients expressed as vectors (on the integral basis).
1167 : * No consistency checks, not memory-clean. */
1168 : GEN
1169 29034 : RgX_to_nfX(GEN nf, GEN x)
1170 : {
1171 : long i, l;
1172 29034 : GEN y = cgetg_copy(x, &l); y[1] = x[1];
1173 236353 : for (i=2; i<l; i++) gel(y,i) = nf_to_scalar_or_basis(nf, gel(x,i));
1174 29034 : return y;
1175 : }
1176 :
1177 : /* Assume nf is a genuine nf. */
1178 : GEN
1179 3465027 : nf_to_scalar_or_alg(GEN nf, GEN x)
1180 : {
1181 3465027 : switch(typ(x))
1182 : {
1183 86953 : case t_INT: case t_FRAC:
1184 86953 : return x;
1185 0 : case t_POLMOD:
1186 0 : x = checknfelt_mod(nf,x,"nf_to_scalar_or_alg");
1187 0 : if (typ(x) != t_POL) return x;
1188 : /* fall through */
1189 : case t_POL:
1190 : {
1191 4704 : GEN T = nf_get_pol(nf);
1192 4704 : long l = lg(x);
1193 4704 : if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_alg", x,T);
1194 4704 : if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
1195 4704 : if (l == 2) return gen_0;
1196 4704 : if (l == 3) return gel(x,2);
1197 3192 : return x;
1198 : }
1199 3373320 : case t_COL:
1200 : {
1201 : GEN dx;
1202 3373320 : if (lg(x)-1 != nf_get_degree(nf)) break;
1203 6689559 : if (QV_isscalar(x)) return gel(x,1);
1204 3316158 : x = Q_remove_denom(x, &dx);
1205 3316173 : x = RgV_RgC_mul(nf_get_zkprimpart(nf), x);
1206 3316311 : dx = mul_denom(dx, nf_get_zkden(nf));
1207 3316281 : return gdiv(x,dx);
1208 : }
1209 : }
1210 50 : pari_err_TYPE("nf_to_scalar_or_alg",x);
1211 : return NULL; /* LCOV_EXCL_LINE */
1212 : }
1213 :
1214 : /* gmul(A, RgX_to_RgC(x)), A t_MAT of compatible dimensions */
1215 : GEN
1216 1365 : RgM_RgX_mul(GEN A, GEN x)
1217 : {
1218 1365 : long i, l = lg(x)-1;
1219 : GEN z;
1220 1365 : if (l == 1) return zerocol(nbrows(A));
1221 1351 : z = gmul(gel(x,2), gel(A,1));
1222 2555 : for (i = 2; i < l; i++)
1223 1204 : if (!gequal0(gel(x,i+1))) z = gadd(z, gmul(gel(x,i+1), gel(A,i)));
1224 1351 : return z;
1225 : }
1226 : GEN
1227 9503802 : ZM_ZX_mul(GEN A, GEN x)
1228 : {
1229 9503802 : long i, l = lg(x)-1;
1230 : GEN z;
1231 9503802 : if (l == 1) return zerocol(nbrows(A));
1232 9502668 : z = ZC_Z_mul(gel(A,1), gel(x,2));
1233 30063895 : for (i = 2; i < l ; i++)
1234 20563999 : if (signe(gel(x,i+1))) z = ZC_add(z, ZC_Z_mul(gel(A,i), gel(x,i+1)));
1235 9499896 : return z;
1236 : }
1237 : /* x a t_POL, nf a genuine nf. No garbage collecting. No check. */
1238 : GEN
1239 8908118 : poltobasis(GEN nf, GEN x)
1240 : {
1241 8908118 : GEN d, T = nf_get_pol(nf);
1242 8908046 : if (varn(x) != varn(T)) pari_err_VAR( "poltobasis", x,T);
1243 8907913 : if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
1244 8907814 : x = Q_remove_denom(x, &d);
1245 8907924 : if (!RgX_is_ZX(x)) pari_err_TYPE("poltobasis",x);
1246 8907868 : x = ZM_ZX_mul(nf_get_invzk(nf), x);
1247 8906189 : if (d) x = RgC_Rg_div(x, d);
1248 8906295 : return x;
1249 : }
1250 :
1251 : GEN
1252 918027 : algtobasis(GEN nf, GEN x)
1253 : {
1254 : pari_sp av;
1255 :
1256 918027 : nf = checknf(nf);
1257 918024 : switch(typ(x))
1258 : {
1259 113254 : case t_POLMOD:
1260 113254 : if (!RgX_equal_var(nf_get_pol(nf),gel(x,1)))
1261 7 : pari_err_MODULUS("algtobasis", nf_get_pol(nf),gel(x,1));
1262 113247 : x = gel(x,2);
1263 113247 : switch(typ(x))
1264 : {
1265 8197 : case t_INT:
1266 8197 : case t_FRAC: return scalarcol(x, nf_get_degree(nf));
1267 105050 : case t_POL:
1268 105050 : av = avma;
1269 105050 : return gerepileupto(av,poltobasis(nf,x));
1270 : }
1271 0 : break;
1272 :
1273 249368 : case t_POL:
1274 249368 : av = avma;
1275 249368 : return gerepileupto(av,poltobasis(nf,x));
1276 :
1277 83033 : case t_COL:
1278 83033 : if (!RgV_is_QV(x)) pari_err_TYPE("nfalgtobasis",x);
1279 83027 : if (lg(x)-1 != nf_get_degree(nf)) pari_err_DIM("nfalgtobasis");
1280 83027 : return gcopy(x);
1281 :
1282 472372 : case t_INT:
1283 472372 : case t_FRAC: return scalarcol(x, nf_get_degree(nf));
1284 : }
1285 0 : pari_err_TYPE("algtobasis",x);
1286 : return NULL; /* LCOV_EXCL_LINE */
1287 : }
1288 :
1289 : GEN
1290 47488 : rnfbasistoalg(GEN rnf,GEN x)
1291 : {
1292 47488 : const char *f = "rnfbasistoalg";
1293 : long lx, i;
1294 47488 : pari_sp av = avma;
1295 : GEN z, nf, R, T;
1296 :
1297 47488 : checkrnf(rnf);
1298 47488 : nf = rnf_get_nf(rnf);
1299 47488 : T = nf_get_pol(nf);
1300 47488 : R = QXQX_to_mod_shallow(rnf_get_pol(rnf), T);
1301 47488 : switch(typ(x))
1302 : {
1303 875 : case t_COL:
1304 875 : z = cgetg_copy(x, &lx);
1305 2597 : for (i=1; i<lx; i++)
1306 : {
1307 1778 : GEN c = nf_to_scalar_or_alg(nf, gel(x,i));
1308 1722 : if (typ(c) == t_POL) c = mkpolmod(c,T);
1309 1722 : gel(z,i) = c;
1310 : }
1311 819 : z = RgV_RgC_mul(gel(rnf_get_zk(rnf),1), z);
1312 735 : return gerepileupto(av, gmodulo(z,R));
1313 :
1314 31227 : case t_POLMOD:
1315 31227 : x = polmod_nffix(f, rnf, x, 0);
1316 30954 : if (typ(x) != t_POL) break;
1317 14261 : retmkpolmod(RgX_copy(x), RgX_copy(R));
1318 1274 : case t_POL:
1319 1274 : if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = gmodulo(x,T); break; }
1320 1029 : if (varn(x) == varn(R))
1321 : {
1322 973 : x = RgX_nffix(f,nf_get_pol(nf),x,0);
1323 973 : return gmodulo(x, R);
1324 : }
1325 56 : pari_err_VAR(f, x,R);
1326 : }
1327 30994 : retmkpolmod(scalarpol(x, varn(R)), RgX_copy(R));
1328 : }
1329 :
1330 : GEN
1331 2240 : matbasistoalg(GEN nf,GEN x)
1332 : {
1333 : long i, j, li, lx;
1334 2240 : GEN z = cgetg_copy(x, &lx);
1335 :
1336 2240 : if (lx == 1) return z;
1337 2233 : switch(typ(x))
1338 : {
1339 77 : case t_VEC: case t_COL:
1340 273 : for (i=1; i<lx; i++) gel(z,i) = basistoalg(nf, gel(x,i));
1341 77 : return z;
1342 2156 : case t_MAT: break;
1343 0 : default: pari_err_TYPE("matbasistoalg",x);
1344 : }
1345 2156 : li = lgcols(x);
1346 8008 : for (j=1; j<lx; j++)
1347 : {
1348 5852 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1349 5852 : gel(z,j) = c;
1350 27377 : for (i=1; i<li; i++) gel(c,i) = basistoalg(nf, gel(xj,i));
1351 : }
1352 2156 : return z;
1353 : }
1354 :
1355 : GEN
1356 30631 : matalgtobasis(GEN nf,GEN x)
1357 : {
1358 : long i, j, li, lx;
1359 30631 : GEN z = cgetg_copy(x, &lx);
1360 :
1361 30631 : if (lx == 1) return z;
1362 30169 : switch(typ(x))
1363 : {
1364 30162 : case t_VEC: case t_COL:
1365 79024 : for (i=1; i<lx; i++) gel(z,i) = algtobasis(nf, gel(x,i));
1366 30161 : return z;
1367 7 : case t_MAT: break;
1368 0 : default: pari_err_TYPE("matalgtobasis",x);
1369 : }
1370 7 : li = lgcols(x);
1371 14 : for (j=1; j<lx; j++)
1372 : {
1373 7 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1374 7 : gel(z,j) = c;
1375 21 : for (i=1; i<li; i++) gel(c,i) = algtobasis(nf, gel(xj,i));
1376 : }
1377 7 : return z;
1378 : }
1379 : GEN
1380 10603 : RgM_to_nfM(GEN nf,GEN x)
1381 : {
1382 : long i, j, li, lx;
1383 10603 : GEN z = cgetg_copy(x, &lx);
1384 :
1385 10603 : if (lx == 1) return z;
1386 10603 : li = lgcols(x);
1387 79212 : for (j=1; j<lx; j++)
1388 : {
1389 68609 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1390 68609 : gel(z,j) = c;
1391 452164 : for (i=1; i<li; i++) gel(c,i) = nf_to_scalar_or_basis(nf, gel(xj,i));
1392 : }
1393 10603 : return z;
1394 : }
1395 : GEN
1396 88492 : RgC_to_nfC(GEN nf, GEN x)
1397 637422 : { pari_APPLY_type(t_COL, nf_to_scalar_or_basis(nf, gel(x,i))) }
1398 :
1399 : /* x a t_POLMOD, supposedly in rnf = K[z]/(T), K = Q[y]/(Tnf) */
1400 : GEN
1401 141303 : polmod_nffix(const char *f, GEN rnf, GEN x, int lift)
1402 141303 : { return polmod_nffix2(f, rnf_get_nfpol(rnf), rnf_get_pol(rnf), x,lift); }
1403 : GEN
1404 141394 : polmod_nffix2(const char *f, GEN T, GEN R, GEN x, int lift)
1405 : {
1406 141394 : if (RgX_equal_var(gel(x,1), R))
1407 : {
1408 128961 : x = gel(x,2);
1409 128961 : if (typ(x) == t_POL && varn(x) == varn(R))
1410 : {
1411 98223 : x = RgX_nffix(f, T, x, lift);
1412 98223 : switch(lg(x))
1413 : {
1414 5789 : case 2: return gen_0;
1415 12155 : case 3: return gel(x,2);
1416 : }
1417 80279 : return x;
1418 : }
1419 : }
1420 43171 : return Rg_nffix(f, T, x, lift);
1421 : }
1422 : GEN
1423 1428 : rnfalgtobasis(GEN rnf,GEN x)
1424 : {
1425 1428 : const char *f = "rnfalgtobasis";
1426 1428 : pari_sp av = avma;
1427 : GEN T, R;
1428 :
1429 1428 : checkrnf(rnf);
1430 1428 : R = rnf_get_pol(rnf);
1431 1428 : T = rnf_get_nfpol(rnf);
1432 1428 : switch(typ(x))
1433 : {
1434 98 : case t_COL:
1435 98 : if (lg(x)-1 != rnf_get_degree(rnf)) pari_err_DIM(f);
1436 49 : x = RgV_nffix(f, T, x, 0);
1437 42 : return gerepilecopy(av, x);
1438 :
1439 1162 : case t_POLMOD:
1440 1162 : x = polmod_nffix(f, rnf, x, 0);
1441 1057 : if (typ(x) != t_POL) break;
1442 714 : return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
1443 112 : case t_POL:
1444 112 : if (varn(x) == varn(T))
1445 : {
1446 42 : RgX_check_QX(x,f);
1447 28 : if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
1448 28 : x = mkpolmod(x,T); break;
1449 : }
1450 70 : x = RgX_nffix(f, T, x, 0);
1451 56 : if (degpol(x) >= degpol(R)) x = RgX_rem(x, R);
1452 56 : return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
1453 : }
1454 427 : return gerepileupto(av, scalarcol(x, rnf_get_degree(rnf)));
1455 : }
1456 :
1457 : /* Given a and b in nf, gives an algebraic integer y in nf such that a-b.y
1458 : * is "small" */
1459 : GEN
1460 259 : nfdiveuc(GEN nf, GEN a, GEN b)
1461 : {
1462 259 : pari_sp av = avma;
1463 259 : a = nfdiv(nf,a,b);
1464 259 : return gerepileupto(av, ground(a));
1465 : }
1466 :
1467 : /* Given a and b in nf, gives a "small" algebraic integer r in nf
1468 : * of the form a-b.y */
1469 : GEN
1470 259 : nfmod(GEN nf, GEN a, GEN b)
1471 : {
1472 259 : pari_sp av = avma;
1473 259 : GEN p1 = gneg_i(nfmul(nf,b,ground(nfdiv(nf,a,b))));
1474 259 : return gerepileupto(av, nfadd(nf,a,p1));
1475 : }
1476 :
1477 : /* Given a and b in nf, gives a two-component vector [y,r] in nf such
1478 : * that r=a-b.y is "small". */
1479 : GEN
1480 259 : nfdivrem(GEN nf, GEN a, GEN b)
1481 : {
1482 259 : pari_sp av = avma;
1483 259 : GEN p1,z, y = ground(nfdiv(nf,a,b));
1484 :
1485 259 : p1 = gneg_i(nfmul(nf,b,y));
1486 259 : z = cgetg(3,t_VEC);
1487 259 : gel(z,1) = gcopy(y);
1488 259 : gel(z,2) = nfadd(nf,a,p1); return gerepileupto(av, z);
1489 : }
1490 :
1491 : /*************************************************************************/
1492 : /** **/
1493 : /** LOGARITHMIC EMBEDDINGS **/
1494 : /** **/
1495 : /*************************************************************************/
1496 :
1497 : static int
1498 4444155 : low_prec(GEN x)
1499 : {
1500 4444155 : switch(typ(x))
1501 : {
1502 0 : case t_INT: return !signe(x);
1503 4444155 : case t_REAL: return !signe(x) || realprec(x) <= DEFAULTPREC;
1504 0 : default: return 0;
1505 : }
1506 : }
1507 :
1508 : static GEN
1509 23284 : cxlog_1(GEN nf) { return zerocol(lg(nf_get_roots(nf))-1); }
1510 : static GEN
1511 542 : cxlog_m1(GEN nf, long prec)
1512 : {
1513 542 : long i, l = lg(nf_get_roots(nf)), r1 = nf_get_r1(nf);
1514 542 : GEN v = cgetg(l, t_COL), p, P;
1515 542 : p = mppi(prec); P = mkcomplex(gen_0, p);
1516 1221 : for (i = 1; i <= r1; i++) gel(v,i) = P; /* IPi*/
1517 542 : if (i < l) P = gmul2n(P,1);
1518 1154 : for ( ; i < l; i++) gel(v,i) = P; /* 2IPi */
1519 542 : return v;
1520 : }
1521 : static GEN
1522 1673932 : ZC_cxlog(GEN nf, GEN x, long prec)
1523 : {
1524 : long i, l, r1;
1525 : GEN v;
1526 1673932 : x = RgM_RgC_mul(nf_get_M(nf), Q_primpart(x));
1527 1673934 : l = lg(x); r1 = nf_get_r1(nf);
1528 4292546 : for (i = 1; i <= r1; i++)
1529 2618612 : if (low_prec(gel(x,i))) return NULL;
1530 3303143 : for ( ; i < l; i++)
1531 1629211 : if (low_prec(gnorm(gel(x,i)))) return NULL;
1532 1673932 : v = cgetg(l,t_COL);
1533 4292544 : for (i = 1; i <= r1; i++) gel(v,i) = glog(gel(x,i),prec);
1534 3303141 : for ( ; i < l; i++) gel(v,i) = gmul2n(glog(gel(x,i),prec),1);
1535 1673935 : return v;
1536 : }
1537 : static GEN
1538 223321 : famat_cxlog(GEN nf, GEN fa, long prec)
1539 : {
1540 223321 : GEN G, E, y = NULL;
1541 : long i, l;
1542 :
1543 223321 : if (typ(fa) != t_MAT) pari_err_TYPE("famat_cxlog",fa);
1544 223321 : if (lg(fa) == 1) return cxlog_1(nf);
1545 223321 : G = gel(fa,1);
1546 223321 : E = gel(fa,2); l = lg(E);
1547 1084913 : for (i = 1; i < l; i++)
1548 : {
1549 861592 : GEN t, e = gel(E,i), x = nf_to_scalar_or_basis(nf, gel(G,i));
1550 : /* multiplicative arch would be better (save logs), but exponents overflow
1551 : * [ could keep track of expo separately, but not worth it ] */
1552 861592 : switch(typ(x))
1553 : { /* ignore positive rationals */
1554 16491 : case t_FRAC: x = gel(x,1); /* fall through */
1555 272231 : case t_INT: if (signe(x) > 0) continue;
1556 94 : if (!mpodd(e)) continue;
1557 38 : t = cxlog_m1(nf, prec); /* we probably should not reach this line */
1558 38 : break;
1559 589361 : default: /* t_COL */
1560 589361 : t = ZC_cxlog(nf,x,prec); if (!t) return NULL;
1561 589361 : t = RgC_Rg_mul(t, e);
1562 : }
1563 589399 : y = y? RgV_add(y,t): t;
1564 : }
1565 223321 : return y ? y: cxlog_1(nf);
1566 : }
1567 : /* Archimedean components: [e_i Log( sigma_i(X) )], where X = primpart(x),
1568 : * and e_i = 1 (resp 2.) for i <= R1 (resp. > R1) */
1569 : GEN
1570 1309041 : nf_cxlog(GEN nf, GEN x, long prec)
1571 : {
1572 1309041 : if (typ(x) == t_MAT) return famat_cxlog(nf,x,prec);
1573 1085720 : x = nf_to_scalar_or_basis(nf,x);
1574 1085719 : switch(typ(x))
1575 : {
1576 0 : case t_FRAC: x = gel(x,1); /* fall through */
1577 1148 : case t_INT:
1578 1148 : return signe(x) > 0? cxlog_1(nf): cxlog_m1(nf, prec);
1579 1084571 : default:
1580 1084571 : return ZC_cxlog(nf, x, prec);
1581 : }
1582 : }
1583 : GEN
1584 97 : nfV_cxlog(GEN nf, GEN x, long prec)
1585 : {
1586 : long i, l;
1587 97 : GEN v = cgetg_copy(x, &l);
1588 167 : for (i = 1; i < l; i++)
1589 70 : if (!(gel(v,i) = nf_cxlog(nf, gel(x,i), prec))) return NULL;
1590 97 : return v;
1591 : }
1592 :
1593 : static GEN
1594 16681 : scalar_logembed(GEN nf, GEN u, GEN *emb)
1595 : {
1596 : GEN v, logu;
1597 16681 : long i, s = signe(u), RU = lg(nf_get_roots(nf))-1, R1 = nf_get_r1(nf);
1598 :
1599 16681 : if (!s) pari_err_DOMAIN("nflogembed","argument","=",gen_0,u);
1600 16681 : v = cgetg(RU+1, t_COL); logu = logr_abs(u);
1601 18641 : for (i = 1; i <= R1; i++) gel(v,i) = logu;
1602 16681 : if (i <= RU)
1603 : {
1604 15799 : GEN logu2 = shiftr(logu,1);
1605 61635 : for ( ; i <= RU; i++) gel(v,i) = logu2;
1606 : }
1607 16681 : if (emb) *emb = const_col(RU, u);
1608 16681 : return v;
1609 : }
1610 :
1611 : static GEN
1612 1309 : famat_logembed(GEN nf,GEN x,GEN *emb,long prec)
1613 : {
1614 1309 : GEN A, M, T, a, t, g = gel(x,1), e = gel(x,2);
1615 1309 : long i, l = lg(e);
1616 :
1617 1309 : if (l == 1) return scalar_logembed(nf, real_1(prec), emb);
1618 1309 : A = NULL; T = emb? cgetg(l, t_COL): NULL;
1619 1309 : if (emb) *emb = M = mkmat2(T, e);
1620 63434 : for (i = 1; i < l; i++)
1621 : {
1622 62125 : a = nflogembed(nf, gel(g,i), &t, prec);
1623 62125 : if (!a) return NULL;
1624 62125 : a = RgC_Rg_mul(a, gel(e,i));
1625 62125 : A = A? RgC_add(A, a): a;
1626 62125 : if (emb) gel(T,i) = t;
1627 : }
1628 1309 : return A;
1629 : }
1630 :
1631 : /* Get archimedean components: [e_i log( | sigma_i(x) | )], with e_i = 1
1632 : * (resp 2.) for i <= R1 (resp. > R1) and set emb to the embeddings of x.
1633 : * Return NULL if precision problem */
1634 : GEN
1635 99952 : nflogembed(GEN nf, GEN x, GEN *emb, long prec)
1636 : {
1637 : long i, l, r1;
1638 : GEN v, t;
1639 :
1640 99952 : if (typ(x) == t_MAT) return famat_logembed(nf,x,emb,prec);
1641 98643 : x = nf_to_scalar_or_basis(nf,x);
1642 98643 : if (typ(x) != t_COL) return scalar_logembed(nf, gtofp(x,prec), emb);
1643 81962 : x = RgM_RgC_mul(nf_get_M(nf), x);
1644 81963 : l = lg(x); r1 = nf_get_r1(nf); v = cgetg(l,t_COL);
1645 108863 : for (i = 1; i <= r1; i++)
1646 : {
1647 26901 : t = gabs(gel(x,i),prec); if (low_prec(t)) return NULL;
1648 26901 : gel(v,i) = glog(t,prec);
1649 : }
1650 251397 : for ( ; i < l; i++)
1651 : {
1652 169434 : t = gnorm(gel(x,i)); if (low_prec(t)) return NULL;
1653 169435 : gel(v,i) = glog(t,prec);
1654 : }
1655 81963 : if (emb) *emb = x;
1656 81963 : return v;
1657 : }
1658 :
1659 : /*************************************************************************/
1660 : /** **/
1661 : /** REAL EMBEDDINGS **/
1662 : /** **/
1663 : /*************************************************************************/
1664 : static GEN
1665 485868 : sarch_get_cyc(GEN sarch) { return gel(sarch,1); }
1666 : static GEN
1667 664285 : sarch_get_archp(GEN sarch) { return gel(sarch,2); }
1668 : static GEN
1669 163466 : sarch_get_MI(GEN sarch) { return gel(sarch,3); }
1670 : static GEN
1671 163464 : sarch_get_lambda(GEN sarch) { return gel(sarch,4); }
1672 : static GEN
1673 163465 : sarch_get_F(GEN sarch) { return gel(sarch,5); }
1674 :
1675 : /* x not a scalar, true nf, return number of positive roots of char_x */
1676 : static long
1677 1268 : num_positive(GEN nf, GEN x)
1678 : {
1679 1268 : GEN T = nf_get_pol(nf), B, charx;
1680 : long dnf, vnf, N;
1681 1268 : x = nf_to_scalar_or_alg(nf, x); /* not a scalar */
1682 1268 : charx = ZXQ_charpoly(x, T, 0);
1683 1268 : charx = ZX_radical(charx);
1684 1268 : N = degpol(T) / degpol(charx);
1685 : /* real places are unramified ? */
1686 1268 : if (N == 1 || ZX_sturm(charx) * N == nf_get_r1(nf))
1687 1261 : return ZX_sturmpart(charx, mkvec2(gen_0,mkoo())) * N;
1688 : /* painful case, multiply by random square until primitive */
1689 7 : dnf = nf_get_degree(nf);
1690 7 : vnf = varn(T);
1691 7 : B = int2n(10);
1692 : for(;;)
1693 0 : {
1694 7 : GEN y = RgXQ_sqr(random_FpX(dnf, vnf, B), T);
1695 7 : y = RgXQ_mul(x, y, T);
1696 7 : charx = ZXQ_charpoly(y, T, 0);
1697 7 : if (ZX_is_squarefree(charx))
1698 7 : return ZX_sturmpart(charx, mkvec2(gen_0,mkoo())) * N;
1699 : }
1700 : }
1701 :
1702 : /* x a QC: return sigma_k(x) where 1 <= k <= r1+r2; correct but inefficient
1703 : * if x in Q. M = nf_get_M(nf) */
1704 : static GEN
1705 658 : nfembed_i(GEN M, GEN x, long k)
1706 : {
1707 658 : long i, l = lg(M);
1708 658 : GEN z = gel(x,1);
1709 4074 : for (i = 2; i < l; i++) z = gadd(z, gmul(gcoeff(M,k,i), gel(x,i)));
1710 658 : return z;
1711 : }
1712 : GEN
1713 0 : nfembed(GEN nf, GEN x, long k)
1714 : {
1715 0 : pari_sp av = avma;
1716 0 : nf = checknf(nf);
1717 0 : x = nf_to_scalar_or_basis(nf,x);
1718 0 : if (typ(x) != t_COL) return gerepilecopy(av, x);
1719 0 : return gerepileupto(av, nfembed_i(nf_get_M(nf),x,k));
1720 : }
1721 :
1722 : /* x a ZC */
1723 : static GEN
1724 907727 : zk_embed(GEN M, GEN x, long k)
1725 : {
1726 907727 : long i, l = lg(x);
1727 907727 : GEN z = gel(x,1); /* times M[k,1], which is 1 */
1728 2967051 : for (i = 2; i < l; i++) z = mpadd(z, mpmul(gcoeff(M,k,i), gel(x,i)));
1729 907721 : return z;
1730 : }
1731 :
1732 : /* Given floating point approximation z of sigma_k(x), decide its sign
1733 : * [0/+, 1/- and -1 for FAIL] */
1734 : static long
1735 888743 : eval_sign_embed(GEN z)
1736 : {
1737 888743 : if (typ(z) == t_REAL)
1738 : {
1739 888743 : long l = realprec(z);
1740 888743 : if (l <= LOWDEFAULTPREC
1741 888743 : || (l == LOWDEFAULTPREC + 1 && !z[l-1])) return -1; /* dubious, fail */
1742 887938 : if (expo(z) < 16 - prec2nbits(l)) return -1; /* same */
1743 : }
1744 887896 : return (signe(z) < 1)? 1: 0;
1745 : }
1746 : /* return v such that (-1)^v = sign(sigma_k(x)), x primitive ZC */
1747 : static long
1748 786030 : eval_sign(GEN M, GEN x, long k)
1749 786030 : { return eval_sign_embed( zk_embed(M, x, k) ); }
1750 :
1751 : /* check that signs[i..#signs] == s; signs = NULL encodes "totally positive" */
1752 : static int
1753 0 : oksigns(long l, GEN signs, long i, long s)
1754 : {
1755 0 : if (!signs) return s == 0;
1756 0 : for (; i < l; i++)
1757 0 : if (signs[i] != s) return 0;
1758 0 : return 1;
1759 : }
1760 : /* check that signs[i] = s and signs[i+1..#signs] = 1-s */
1761 : static int
1762 0 : oksigns2(long l, GEN signs, long i, long s)
1763 : {
1764 0 : if (!signs) return s == 0 && i == l-1;
1765 0 : return signs[i] == s && oksigns(l, signs, i+1, 1-s);
1766 : }
1767 :
1768 : /* true nf, x a ZC (primitive for efficiency) which is not a scalar; embx its
1769 : * embeddings or NULL */
1770 : static int
1771 84964 : nfchecksigns_i(GEN nf, GEN x, GEN embx, GEN signs, GEN archp)
1772 : {
1773 84964 : long l = lg(archp), i;
1774 84964 : GEN M = nf_get_M(nf), sarch = NULL;
1775 84964 : long np = -1;
1776 134855 : for (i = 1; i < l; i++)
1777 : {
1778 : long s;
1779 104083 : if (embx)
1780 102719 : s = eval_sign_embed(gel(embx,i));
1781 : else
1782 1364 : s = eval_sign(M, x, archp[i]);
1783 : /* 0 / + or 1 / -; -1 for FAIL */
1784 104083 : if (s < 0) /* failure */
1785 : {
1786 0 : long ni, r1 = nf_get_r1(nf);
1787 : GEN xi;
1788 0 : if (np < 0)
1789 : {
1790 0 : np = num_positive(nf, x);
1791 0 : if (np == 0) return oksigns(l, signs, i, 1);
1792 0 : if (np == r1) return oksigns(l, signs, i, 0);
1793 0 : sarch = nfarchstar(nf, NULL, identity_perm(r1));
1794 : }
1795 0 : xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
1796 0 : xi = Q_primpart(xi);
1797 0 : ni = num_positive(nf, nfmuli(nf,x,xi));
1798 0 : if (ni == 0) return oksigns2(l, signs, i, 0);
1799 0 : if (ni == r1) return oksigns2(l, signs, i, 1);
1800 0 : s = ni < np? 0: 1;
1801 : }
1802 104083 : if (s != (signs? signs[i]: 0)) return 0;
1803 : }
1804 30772 : return 1;
1805 : }
1806 : static void
1807 775 : pl_convert(GEN pl, GEN *psigns, GEN *parchp)
1808 : {
1809 775 : long i, j, l = lg(pl);
1810 775 : GEN signs = cgetg(l, t_VECSMALL);
1811 775 : GEN archp = cgetg(l, t_VECSMALL);
1812 2576 : for (i = j = 1; i < l; i++)
1813 : {
1814 1801 : if (!pl[i]) continue;
1815 1403 : archp[j] = i;
1816 1403 : signs[j] = (pl[i] < 0)? 1: 0;
1817 1403 : j++;
1818 : }
1819 775 : setlg(archp, j); *parchp = archp;
1820 775 : setlg(signs, j); *psigns = signs;
1821 775 : }
1822 : /* pl : requested signs for real embeddings, 0 = no sign constraint */
1823 : int
1824 14719 : nfchecksigns(GEN nf, GEN x, GEN pl)
1825 : {
1826 14719 : pari_sp av = avma;
1827 : GEN signs, archp;
1828 14719 : nf = checknf(nf);
1829 14719 : x = nf_to_scalar_or_basis(nf,x);
1830 14719 : if (typ(x) != t_COL)
1831 : {
1832 13944 : long i, l = lg(pl), s = gsigne(x);
1833 27853 : for (i = 1; i < l; i++)
1834 13909 : if (pl[i] && pl[i] != s) return gc_bool(av,0);
1835 13944 : return gc_bool(av,1);
1836 : }
1837 775 : pl_convert(pl, &signs, &archp);
1838 775 : return gc_bool(av, nfchecksigns_i(nf, x, NULL, signs, archp));
1839 : }
1840 :
1841 : /* signs = NULL: totally positive, else sign[i] = 0 (+) or 1 (-) */
1842 : static GEN
1843 163464 : get_C(GEN lambda, long l, GEN signs)
1844 : {
1845 : long i;
1846 : GEN C, mlambda;
1847 163464 : if (!signs) return const_vec(l-1, lambda);
1848 133798 : C = cgetg(l, t_COL); mlambda = gneg(lambda);
1849 343135 : for (i = 1; i < l; i++) gel(C,i) = signs[i]? mlambda: lambda;
1850 133800 : return C;
1851 : }
1852 : /* signs = NULL: totally positive at archp.
1853 : * Assume that a t_COL x is not a scalar */
1854 : static GEN
1855 274302 : nfsetsigns(GEN nf, GEN signs, GEN x, GEN sarch)
1856 : {
1857 274302 : long i, l = lg(sarch_get_archp(sarch));
1858 : GEN ex;
1859 : /* Is signature already correct ? */
1860 274297 : if (typ(x) != t_COL)
1861 : {
1862 190117 : long s = gsigne(x);
1863 190119 : if (!s) i = 1;
1864 190098 : else if (!signs)
1865 4676 : i = (s < 0)? 1: l;
1866 : else
1867 : {
1868 185422 : s = s < 0? 1: 0;
1869 313354 : for (i = 1; i < l; i++)
1870 237079 : if (signs[i] != s) break;
1871 : }
1872 190119 : ex = (i < l)? const_col(l-1, x): NULL;
1873 : }
1874 : else
1875 : { /* inefficient if x scalar, wrong if x = 0 */
1876 84180 : pari_sp av = avma;
1877 84180 : GEN cex, M = nf_get_M(nf), archp = sarch_get_archp(sarch);
1878 84188 : GEN xp = Q_primitive_part(x,&cex);
1879 84187 : ex = cgetg(l,t_COL);
1880 205888 : for (i = 1; i < l; i++) gel(ex,i) = zk_embed(M,xp,archp[i]);
1881 84189 : if (nfchecksigns_i(nf, xp, ex, signs, archp)) { ex = NULL; set_avma(av); }
1882 54155 : else if (cex) ex = RgC_Rg_mul(ex, cex); /* put back content */
1883 : }
1884 274307 : if (ex)
1885 : { /* If no, fix it */
1886 163466 : GEN MI = sarch_get_MI(sarch), F = sarch_get_F(sarch);
1887 163464 : GEN lambda = sarch_get_lambda(sarch);
1888 163464 : GEN t = RgC_sub(get_C(lambda, l, signs), ex);
1889 163467 : t = grndtoi(RgM_RgC_mul(MI,t), NULL);
1890 163465 : if (lg(F) != 1) t = ZM_ZC_mul(F, t);
1891 163457 : x = typ(x) == t_COL? RgC_add(t, x): RgC_Rg_add(t, x);
1892 : }
1893 274299 : return x;
1894 : }
1895 : /* - true nf
1896 : * - sarch = nfarchstar(nf, F);
1897 : * - x encodes a vector of signs at arch.archp: either a t_VECSMALL
1898 : * (vector of signs as {0,1}-vector), NULL (totally positive at archp),
1899 : * or a nonzero number field element (replaced by its signature at archp);
1900 : * - y is a nonzero number field element
1901 : * Return z = y (mod F) with signs(y, archp) = signs(x) (a {0,1}-vector).
1902 : * Not stack-clean */
1903 : GEN
1904 305807 : set_sign_mod_divisor(GEN nf, GEN x, GEN y, GEN sarch)
1905 : {
1906 305807 : GEN archp = sarch_get_archp(sarch);
1907 305807 : if (lg(archp) == 1) return y;
1908 272550 : if (x && typ(x) != t_VECSMALL) x = nfsign_arch(nf, x, archp);
1909 272550 : return nfsetsigns(nf, x, nf_to_scalar_or_basis(nf,y), sarch);
1910 : }
1911 :
1912 : static GEN
1913 83323 : setsigns_init(GEN nf, GEN archp, GEN F, GEN DATA)
1914 : {
1915 83323 : GEN lambda, Mr = rowpermute(nf_get_M(nf), archp), MI = F? RgM_mul(Mr,F): Mr;
1916 83325 : lambda = gmul2n(matrixnorm(MI,DEFAULTPREC), -1);
1917 83325 : if (typ(lambda) != t_REAL) lambda = gmul(lambda, uutoQ(1001,1000));
1918 83325 : if (lg(archp) < lg(MI))
1919 : {
1920 58881 : GEN perm = gel(indexrank(MI), 2);
1921 58882 : if (!F) F = matid(nf_get_degree(nf));
1922 58882 : MI = vecpermute(MI, perm);
1923 58881 : F = vecpermute(F, perm);
1924 : }
1925 83325 : if (!F) F = cgetg(1,t_MAT);
1926 83325 : MI = RgM_inv(MI);
1927 83326 : return mkvec5(DATA, archp, MI, lambda, F);
1928 : }
1929 : /* F nonzero integral ideal in HNF (or NULL: Z_K), compute elements in 1+F
1930 : * whose sign matrix at archp is identity; archp in 'indices' format */
1931 : GEN
1932 259223 : nfarchstar(GEN nf, GEN F, GEN archp)
1933 : {
1934 259223 : long nba = lg(archp) - 1;
1935 259223 : if (!nba) return mkvec2(cgetg(1,t_VEC), archp);
1936 81579 : if (F && equali1(gcoeff(F,1,1))) F = NULL;
1937 81579 : if (F) F = idealpseudored(F, nf_get_roundG(nf));
1938 81569 : return setsigns_init(nf, archp, F, const_vec(nba, gen_2));
1939 : }
1940 :
1941 : /*************************************************************************/
1942 : /** **/
1943 : /** IDEALCHINESE **/
1944 : /** **/
1945 : /*************************************************************************/
1946 : static int
1947 4192 : isprfact(GEN x)
1948 : {
1949 : long i, l;
1950 : GEN L, E;
1951 4192 : if (typ(x) != t_MAT || lg(x) != 3) return 0;
1952 4192 : L = gel(x,1); l = lg(L);
1953 4192 : E = gel(x,2);
1954 13965 : for(i=1; i<l; i++)
1955 : {
1956 9773 : checkprid(gel(L,i));
1957 9773 : if (typ(gel(E,i)) != t_INT) return 0;
1958 : }
1959 4192 : return 1;
1960 : }
1961 :
1962 : /* initialize projectors mod pr[i]^e[i] for idealchinese */
1963 : static GEN
1964 4192 : pr_init(GEN nf, GEN fa, GEN w, GEN dw)
1965 : {
1966 4192 : GEN U, E, F, FZ, L = gel(fa,1), E0 = gel(fa,2);
1967 4192 : long i, r = lg(L);
1968 :
1969 4192 : if (w && lg(w) != r) pari_err_TYPE("idealchinese", w);
1970 4192 : if (r == 1 && !dw) return cgetg(1,t_VEC);
1971 4178 : E = leafcopy(E0); /* do not destroy fa[2] */
1972 13951 : for (i = 1; i < r; i++)
1973 9773 : if (signe(gel(E,i)) < 0) gel(E,i) = gen_0;
1974 4178 : F = factorbackprime(nf, L, E);
1975 4178 : if (dw)
1976 : {
1977 679 : F = ZM_Z_mul(F, dw);
1978 1568 : for (i = 1; i < r; i++)
1979 : {
1980 889 : GEN pr = gel(L,i);
1981 889 : long e = itos(gel(E0,i)), v = idealval(nf, dw, pr);
1982 889 : if (e >= 0)
1983 882 : gel(E,i) = addiu(gel(E,i), v);
1984 7 : else if (v + e <= 0)
1985 0 : F = idealmulpowprime(nf, F, pr, stoi(-v)); /* coprime to pr */
1986 : else
1987 : {
1988 7 : F = idealmulpowprime(nf, F, pr, stoi(e));
1989 7 : gel(E,i) = stoi(v + e);
1990 : }
1991 : }
1992 : }
1993 4178 : U = cgetg(r, t_VEC);
1994 13951 : for (i = 1; i < r; i++)
1995 : {
1996 : GEN u;
1997 9773 : if (w && gequal0(gel(w,i))) u = gen_0; /* unused */
1998 : else
1999 : {
2000 9696 : GEN pr = gel(L,i), e = gel(E,i), t;
2001 9696 : t = idealdivpowprime(nf,F, pr, e);
2002 9696 : u = hnfmerge_get_1(t, idealpow(nf, pr, e));
2003 9696 : if (!u) pari_err_COPRIME("idealchinese", t,pr);
2004 : }
2005 9773 : gel(U,i) = u;
2006 : }
2007 4178 : FZ = gcoeff(F, 1, 1);
2008 4178 : F = idealpseudored(F, nf_get_roundG(nf));
2009 4178 : return mkvec2(mkvec2(F, FZ), U);
2010 : }
2011 :
2012 : static GEN
2013 2247 : pl_normalize(GEN nf, GEN pl)
2014 : {
2015 2247 : const char *fun = "idealchinese";
2016 2247 : if (lg(pl)-1 != nf_get_r1(nf)) pari_err_TYPE(fun,pl);
2017 2247 : switch(typ(pl))
2018 : {
2019 693 : case t_VEC: RgV_check_ZV(pl,fun); pl = ZV_to_zv(pl);
2020 : /* fall through */
2021 2247 : case t_VECSMALL: break;
2022 0 : default: pari_err_TYPE(fun,pl);
2023 : }
2024 2247 : return pl;
2025 : }
2026 :
2027 : static int
2028 9401 : is_chineseinit(GEN x)
2029 : {
2030 : GEN fa, pl;
2031 : long l;
2032 9401 : if (typ(x) != t_VEC || lg(x)!=3) return 0;
2033 7574 : fa = gel(x,1);
2034 7574 : pl = gel(x,2);
2035 7574 : if (typ(fa) != t_VEC || typ(pl) != t_VEC) return 0;
2036 4207 : l = lg(fa);
2037 4207 : if (l != 1)
2038 : {
2039 : GEN z;
2040 4165 : if (l != 3) return 0;
2041 4165 : z = gel(fa, 1);
2042 4165 : if (typ(z) != t_VEC || lg(z) != 3 || typ(gel(z,1)) != t_MAT
2043 4158 : || typ(gel(z,2)) != t_INT
2044 4158 : || typ(gel(fa,2)) != t_VEC)
2045 7 : return 0;
2046 : }
2047 4200 : l = lg(pl);
2048 4200 : if (l != 1)
2049 : {
2050 665 : if (l != 6 || typ(gel(pl,3)) != t_MAT || typ(gel(pl,1)) != t_VECSMALL
2051 665 : || typ(gel(pl,2)) != t_VECSMALL)
2052 0 : return 0;
2053 : }
2054 4200 : return 1;
2055 : }
2056 :
2057 : /* nf a true 'nf' */
2058 : static GEN
2059 4647 : chineseinit_i(GEN nf, GEN fa, GEN w, GEN dw)
2060 : {
2061 4647 : const char *fun = "idealchineseinit";
2062 4647 : GEN archp = NULL, pl = NULL;
2063 4647 : switch(typ(fa))
2064 : {
2065 2247 : case t_VEC:
2066 2247 : if (is_chineseinit(fa))
2067 : {
2068 0 : if (dw) pari_err_DOMAIN(fun, "denom(y)", "!=", gen_1, w);
2069 0 : return fa;
2070 : }
2071 2247 : if (lg(fa) != 3) pari_err_TYPE(fun, fa);
2072 : /* of the form [x,s] */
2073 2247 : pl = pl_normalize(nf, gel(fa,2));
2074 2247 : fa = gel(fa,1);
2075 2247 : archp = vecsmall01_to_indices(pl);
2076 : /* keep pr_init, reset pl */
2077 2247 : if (is_chineseinit(fa)) { fa = gel(fa,1); break; }
2078 : /* fall through */
2079 : case t_MAT: /* factorization? */
2080 4192 : if (isprfact(fa)) { fa = pr_init(nf, fa, w, dw); break; }
2081 0 : default: pari_err_TYPE(fun,fa);
2082 : }
2083 :
2084 4647 : if (!pl) pl = cgetg(1,t_VEC);
2085 : else
2086 : {
2087 2247 : long r = lg(archp);
2088 2247 : if (r == 1) pl = cgetg(1, t_VEC);
2089 : else
2090 : {
2091 1743 : GEN F = (lg(fa) == 1)? NULL: gmael(fa,1,1), signs = cgetg(r, t_VECSMALL);
2092 : long i;
2093 5054 : for (i = 1; i < r; i++) signs[i] = (pl[archp[i]] < 0)? 1: 0;
2094 1743 : pl = setsigns_init(nf, archp, F, signs);
2095 : }
2096 : }
2097 4647 : return mkvec2(fa, pl);
2098 : }
2099 :
2100 : /* Given a prime ideal factorization x, possibly with 0 or negative exponents,
2101 : * and a vector w of elements of nf, gives b such that
2102 : * v_p(b-w_p)>=v_p(x) for all prime ideals p in the ideal factorization
2103 : * and v_p(b)>=0 for all other p, using the standard proof given in GTM 138. */
2104 : GEN
2105 8392 : idealchinese(GEN nf, GEN x0, GEN w)
2106 : {
2107 8392 : const char *fun = "idealchinese";
2108 8392 : pari_sp av = avma;
2109 8392 : GEN x = x0, x1, x2, s, dw, F;
2110 :
2111 8392 : nf = checknf(nf);
2112 8392 : if (!w) return gerepilecopy(av, chineseinit_i(nf,x,NULL,NULL));
2113 :
2114 4907 : if (typ(w) != t_VEC) pari_err_TYPE(fun,w);
2115 4907 : w = Q_remove_denom(matalgtobasis(nf,w), &dw);
2116 4907 : if (!is_chineseinit(x)) x = chineseinit_i(nf,x,w,dw);
2117 : /* x is a 'chineseinit' */
2118 4907 : x1 = gel(x,1); s = NULL;
2119 4907 : x2 = gel(x,2);
2120 4907 : if (lg(x1) == 1) { F = NULL; dw = NULL; }
2121 : else
2122 : {
2123 4865 : GEN U = gel(x1,2), FZ;
2124 4865 : long i, r = lg(w);
2125 4865 : F = gmael(x1,1,1); FZ = gmael(x1,1,2);
2126 17596 : for (i=1; i<r; i++)
2127 12731 : if (!ZV_equal0(gel(w,i)))
2128 : {
2129 9626 : GEN t = nfmuli(nf, gel(U,i), gel(w,i));
2130 9626 : s = s? ZC_add(s,t): t;
2131 : }
2132 4865 : if (s)
2133 : {
2134 4844 : s = ZC_reducemodmatrix(s, F);
2135 4844 : if (dw && x == x0) /* input was a chineseinit */
2136 : {
2137 7 : dw = modii(dw, FZ);
2138 7 : s = FpC_Fp_mul(s, Fp_inv(dw, FZ), FZ);
2139 7 : dw = NULL;
2140 : }
2141 4844 : if (ZV_isscalar(s)) s = icopy(gel(s,1));
2142 : }
2143 : }
2144 4907 : if (lg(x2) != 1)
2145 : {
2146 1750 : s = nfsetsigns(nf, gel(x2,1), s? s: gen_0, x2);
2147 1750 : if (typ(s) == t_COL && QV_isscalar(s))
2148 : {
2149 294 : s = gel(s,1); if (!dw) s = gcopy(s);
2150 : }
2151 : }
2152 3157 : else if (!s) return gc_const(av, gen_0);
2153 4858 : return gerepileupto(av, dw? gdiv(s, dw): s);
2154 : }
2155 :
2156 : /*************************************************************************/
2157 : /** **/
2158 : /** (Z_K/I)^* **/
2159 : /** **/
2160 : /*************************************************************************/
2161 : GEN
2162 2247 : vecsmall01_to_indices(GEN v)
2163 : {
2164 2247 : long i, k, l = lg(v);
2165 2247 : GEN p = new_chunk(l) + l;
2166 6594 : for (k=1, i=l-1; i; i--)
2167 4347 : if (v[i]) { *--p = i; k++; }
2168 2247 : *--p = evallg(k) | evaltyp(t_VECSMALL);
2169 2247 : set_avma((pari_sp)p); return p;
2170 : }
2171 : GEN
2172 1040236 : vec01_to_indices(GEN v)
2173 : {
2174 : long i, k, l;
2175 : GEN p;
2176 :
2177 1040236 : switch (typ(v))
2178 : {
2179 993756 : case t_VECSMALL: return v;
2180 46480 : case t_VEC: break;
2181 0 : default: pari_err_TYPE("vec01_to_indices",v);
2182 : }
2183 46480 : l = lg(v);
2184 46480 : p = new_chunk(l) + l;
2185 140000 : for (k=1, i=l-1; i; i--)
2186 93520 : if (signe(gel(v,i))) { *--p = i; k++; }
2187 46480 : *--p = evallg(k) | evaltyp(t_VECSMALL);
2188 46480 : set_avma((pari_sp)p); return p;
2189 : }
2190 : GEN
2191 136713 : indices_to_vec01(GEN p, long r)
2192 : {
2193 136713 : long i, l = lg(p);
2194 136713 : GEN v = zerovec(r);
2195 206412 : for (i = 1; i < l; i++) gel(v, p[i]) = gen_1;
2196 136711 : return v;
2197 : }
2198 :
2199 : /* return (column) vector of R1 signatures of x (0 or 1) */
2200 : GEN
2201 993756 : nfsign_arch(GEN nf, GEN x, GEN arch)
2202 : {
2203 993756 : GEN sarch, M, V, archp = vec01_to_indices(arch);
2204 993756 : long i, s, np, n = lg(archp)-1;
2205 : pari_sp av;
2206 :
2207 993756 : if (!n) return cgetg(1,t_VECSMALL);
2208 792157 : if (typ(x) == t_MAT)
2209 : { /* factorisation */
2210 273033 : GEN g = gel(x,1), e = gel(x,2);
2211 273033 : long l = lg(g);
2212 273033 : V = zero_zv(n);
2213 720979 : for (i = 1; i < l; i++)
2214 447945 : if (mpodd(gel(e,i)))
2215 387076 : Flv_add_inplace(V, nfsign_arch(nf,gel(g,i),archp), 2);
2216 273034 : set_avma((pari_sp)V); return V;
2217 : }
2218 519124 : av = avma; V = cgetg(n+1,t_VECSMALL);
2219 519124 : x = nf_to_scalar_or_basis(nf, x);
2220 519125 : switch(typ(x))
2221 : {
2222 131500 : case t_INT:
2223 131500 : s = signe(x);
2224 131500 : if (!s) pari_err_DOMAIN("nfsign_arch","element","=",gen_0,x);
2225 131500 : set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
2226 651 : case t_FRAC:
2227 651 : s = signe(gel(x,1));
2228 651 : set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
2229 : }
2230 386974 : x = Q_primpart(x); M = nf_get_M(nf); sarch = NULL; np = -1;
2231 1170827 : for (i = 1; i <= n; i++)
2232 : {
2233 784665 : long s = eval_sign(M, x, archp[i]);
2234 784661 : if (s < 0) /* failure */
2235 : {
2236 847 : long ni, r1 = nf_get_r1(nf);
2237 : GEN xi;
2238 848 : if (np < 0)
2239 : {
2240 848 : np = num_positive(nf, x);
2241 848 : if (np == 0) { set_avma(av); return const_vecsmall(n, 1); }
2242 806 : if (np == r1){ set_avma(av); return const_vecsmall(n, 0); }
2243 420 : sarch = nfarchstar(nf, NULL, identity_perm(r1));
2244 : }
2245 420 : xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
2246 420 : xi = Q_primpart(xi);
2247 420 : ni = num_positive(nf, nfmuli(nf,x,xi));
2248 420 : if (ni == 0) { set_avma(av); V = const_vecsmall(n, 1); V[i] = 0; return V; }
2249 413 : if (ni == r1){ set_avma(av); V = const_vecsmall(n, 0); V[i] = 1; return V; }
2250 35 : s = ni < np? 0: 1;
2251 : }
2252 783849 : V[i] = s;
2253 : }
2254 386162 : set_avma((pari_sp)V); return V;
2255 : }
2256 : static void
2257 6783 : chk_ind(const char *s, long i, long r1)
2258 : {
2259 6783 : if (i <= 0) pari_err_DOMAIN(s, "index", "<=", gen_0, stoi(i));
2260 6769 : if (i > r1) pari_err_DOMAIN(s, "index", ">", utoi(r1), utoi(i));
2261 6734 : }
2262 : static GEN
2263 97244 : parse_embed(GEN ind, long r, const char *f)
2264 : {
2265 : long l, i;
2266 97244 : if (!ind) return identity_perm(r);
2267 4718 : switch(typ(ind))
2268 : {
2269 70 : case t_INT: ind = mkvecsmall(itos(ind)); break;
2270 84 : case t_VEC: case t_COL: ind = vec_to_vecsmall(ind); break;
2271 4564 : case t_VECSMALL: break;
2272 0 : default: pari_err_TYPE(f, ind);
2273 : }
2274 4718 : l = lg(ind);
2275 11452 : for (i = 1; i < l; i++) chk_ind(f, ind[i], r);
2276 4669 : return ind;
2277 : }
2278 : GEN
2279 95361 : nfeltsign(GEN nf, GEN x, GEN ind0)
2280 : {
2281 95361 : pari_sp av = avma;
2282 : long i, l;
2283 : GEN v, ind;
2284 95361 : nf = checknf(nf);
2285 95361 : ind = parse_embed(ind0, nf_get_r1(nf), "nfeltsign");
2286 95340 : l = lg(ind);
2287 95340 : if (is_rational_t(typ(x)))
2288 : { /* nfsign_arch would test this, but avoid converting t_VECSMALL -> t_VEC */
2289 : GEN s;
2290 2275 : switch(gsigne(x))
2291 : {
2292 560 : case -1:s = gen_m1; break;
2293 1708 : case 1: s = gen_1; break;
2294 7 : default: s = gen_0; break;
2295 : }
2296 2275 : set_avma(av);
2297 2275 : return (ind0 && typ(ind0) == t_INT)? s: const_vec(l-1, s);
2298 : }
2299 93065 : v = nfsign_arch(nf, x, ind);
2300 93065 : if (ind0 && typ(ind0) == t_INT) { set_avma(av); return v[1]? gen_m1: gen_1; }
2301 93051 : settyp(v, t_VEC);
2302 262311 : for (i = 1; i < l; i++) gel(v,i) = v[i]? gen_m1: gen_1;
2303 93051 : return gerepileupto(av, v);
2304 : }
2305 :
2306 : /* true nf */
2307 : GEN
2308 294 : nfeltembed_i(GEN *pnf, GEN x, GEN ind0, long prec0)
2309 : {
2310 : long i, e, l, r1, r2, prec, prec1;
2311 294 : GEN v, ind, cx, nf = *pnf;
2312 294 : nf_get_sign(nf,&r1,&r2);
2313 294 : x = nf_to_scalar_or_basis(nf, x);
2314 287 : ind = parse_embed(ind0, r1+r2, "nfeltembed");
2315 280 : l = lg(ind);
2316 280 : if (typ(x) != t_COL)
2317 : {
2318 0 : if (!(ind0 && typ(ind0) == t_INT)) x = const_vec(l-1, x);
2319 0 : return x;
2320 : }
2321 280 : x = Q_primitive_part(x, &cx);
2322 280 : prec1 = prec0; e = gexpo(x);
2323 280 : if (e > 8) prec1 += nbits2extraprec(e);
2324 280 : prec = prec1;
2325 280 : if (nf_get_prec(nf) < prec) nf = nfnewprec_shallow(nf, prec);
2326 280 : v = cgetg(l, t_VEC);
2327 : for(;;)
2328 35 : {
2329 315 : GEN M = nf_get_M(nf);
2330 938 : for (i = 1; i < l; i++)
2331 : {
2332 658 : GEN t = nfembed_i(M, x, ind[i]);
2333 658 : long e = gexpo(t);
2334 658 : if (gequal0(t) || precision(t) < prec0
2335 658 : || (e < 0 && prec < prec1 + nbits2extraprec(-e)) ) break;
2336 623 : if (cx) t = gmul(t, cx);
2337 623 : gel(v,i) = t;
2338 : }
2339 315 : if (i == l) break;
2340 35 : prec = precdbl(prec);
2341 35 : if (DEBUGLEVEL>1) pari_warn(warnprec,"eltnfembed", prec);
2342 35 : *pnf = nf = nfnewprec_shallow(nf, prec);
2343 : }
2344 280 : if (ind0 && typ(ind0) == t_INT) v = gel(v,1);
2345 280 : return v;
2346 : }
2347 : GEN
2348 294 : nfeltembed(GEN nf, GEN x, GEN ind0, long prec0)
2349 : {
2350 294 : pari_sp av = avma; nf = checknf(nf);
2351 294 : return gerepilecopy(av, nfeltembed_i(&nf, x, ind0, prec0));
2352 : }
2353 :
2354 : /* number of distinct roots of sigma(f) */
2355 : GEN
2356 1596 : nfpolsturm(GEN nf, GEN f, GEN ind0)
2357 : {
2358 1596 : pari_sp av = avma;
2359 : long d, l, r1, single;
2360 : GEN ind, u, v, vr1, T, s, t;
2361 :
2362 1596 : nf = checknf(nf); T = nf_get_pol(nf); r1 = nf_get_r1(nf);
2363 1596 : ind = parse_embed(ind0, r1, "nfpolsturm");
2364 1575 : single = ind0 && typ(ind0) == t_INT;
2365 1575 : l = lg(ind);
2366 :
2367 1575 : if (gequal0(f)) pari_err_ROOTS0("nfpolsturm");
2368 1568 : if (typ(f) == t_POL && varn(f) != varn(T))
2369 : {
2370 1547 : f = RgX_nffix("nfpolsturm", T, f,1);
2371 1547 : if (lg(f) == 3) f = NULL;
2372 : }
2373 : else
2374 : {
2375 21 : (void)Rg_nffix("nfpolsturm", T, f, 0);
2376 21 : f = NULL;
2377 : }
2378 1568 : if (!f) { set_avma(av); return single? gen_0: zerovec(l-1); }
2379 1547 : d = degpol(f);
2380 1547 : if (d == 1) { set_avma(av); return single? gen_1: const_vec(l-1,gen_1); }
2381 :
2382 1505 : vr1 = const_vecsmall(l-1, 1);
2383 1505 : u = Q_primpart(f); s = ZV_to_zv(nfeltsign(nf, gel(u,d+2), ind));
2384 1505 : v = RgX_deriv(u); t = odd(d)? leafcopy(s): zv_neg(s);
2385 : for(;;)
2386 182 : {
2387 1687 : GEN r = RgX_neg( Q_primpart(RgX_pseudorem(u, v)) ), sr;
2388 1687 : long i, dr = degpol(r);
2389 1687 : if (dr < 0) break;
2390 1687 : sr = ZV_to_zv(nfeltsign(nf, gel(r,dr+2), ind));
2391 4144 : for (i = 1; i < l; i++)
2392 2457 : if (sr[i] != s[i]) { s[i] = sr[i], vr1[i]--; }
2393 1687 : if (odd(dr)) sr = zv_neg(sr);
2394 4144 : for (i = 1; i < l; i++)
2395 2457 : if (sr[i] != t[i]) { t[i] = sr[i], vr1[i]++; }
2396 1687 : if (!dr) break;
2397 182 : u = v; v = r;
2398 : }
2399 1505 : if (single) return gc_stoi(av,vr1[1]);
2400 1498 : return gerepileupto(av, zv_to_ZV(vr1));
2401 : }
2402 :
2403 : /* True nf; return the vector of signs of x; the matrix of such if x is a vector
2404 : * of nf elements */
2405 : GEN
2406 43960 : nfsign(GEN nf, GEN x)
2407 : {
2408 : long i, l;
2409 : GEN archp, S;
2410 :
2411 43960 : archp = identity_perm( nf_get_r1(nf) );
2412 43959 : if (typ(x) != t_VEC) return nfsign_arch(nf, x, archp);
2413 35937 : l = lg(x); S = cgetg(l, t_MAT);
2414 148060 : for (i=1; i<l; i++) gel(S,i) = nfsign_arch(nf, gel(x,i), archp);
2415 35937 : return S;
2416 : }
2417 :
2418 : /* x integral elt, A integral ideal in HNF; reduce x mod A */
2419 : static GEN
2420 8490028 : zk_modHNF(GEN x, GEN A)
2421 8490028 : { return (typ(x) == t_COL)? ZC_hnfrem(x, A): modii(x, gcoeff(A,1,1)); }
2422 :
2423 : /* given an element x in Z_K and an integral ideal y in HNF, coprime with x,
2424 : outputs an element inverse of x modulo y */
2425 : GEN
2426 175 : nfinvmodideal(GEN nf, GEN x, GEN y)
2427 : {
2428 175 : pari_sp av = avma;
2429 175 : GEN a, yZ = gcoeff(y,1,1);
2430 :
2431 175 : if (equali1(yZ)) return gen_0;
2432 175 : x = nf_to_scalar_or_basis(nf, x);
2433 175 : if (typ(x) == t_INT) return gerepileupto(av, Fp_inv(x, yZ));
2434 :
2435 119 : a = hnfmerge_get_1(idealhnf_principal(nf,x), y);
2436 119 : if (!a) pari_err_INV("nfinvmodideal", x);
2437 119 : return gerepileupto(av, zk_modHNF(nfdiv(nf,a,x), y));
2438 : }
2439 :
2440 : static GEN
2441 3146385 : nfsqrmodideal(GEN nf, GEN x, GEN id)
2442 3146385 : { return zk_modHNF(nfsqri(nf,x), id); }
2443 : static GEN
2444 7744556 : nfmulmodideal(GEN nf, GEN x, GEN y, GEN id)
2445 7744556 : { return x? zk_modHNF(nfmuli(nf,x,y), id): y; }
2446 : /* assume x integral, k integer, A in HNF */
2447 : GEN
2448 5846577 : nfpowmodideal(GEN nf,GEN x,GEN k,GEN A)
2449 : {
2450 5846577 : long s = signe(k);
2451 : pari_sp av;
2452 : GEN y;
2453 :
2454 5846577 : if (!s) return gen_1;
2455 5846577 : av = avma;
2456 5846577 : x = nf_to_scalar_or_basis(nf, x);
2457 5846791 : if (typ(x) != t_COL) return Fp_pow(x, k, gcoeff(A,1,1));
2458 2656780 : if (s < 0) { k = negi(k); x = nfinvmodideal(nf, x,A); }
2459 2656780 : if (equali1(k)) return gerepileupto(av, s > 0? zk_modHNF(x, A): x);
2460 1268247 : for(y = NULL;;)
2461 : {
2462 4414733 : if (mpodd(k)) y = nfmulmodideal(nf,y,x,A);
2463 4414699 : k = shifti(k,-1); if (!signe(k)) break;
2464 3145916 : x = nfsqrmodideal(nf,x,A);
2465 : }
2466 1268201 : return gerepileupto(av, y);
2467 : }
2468 :
2469 : /* a * g^n mod id */
2470 : static GEN
2471 4699139 : nfmulpowmodideal(GEN nf, GEN a, GEN g, GEN n, GEN id)
2472 : {
2473 4699139 : return nfmulmodideal(nf, a, nfpowmodideal(nf,g,n,id), id);
2474 : }
2475 :
2476 : /* assume (num(g[i]), id) = 1 for all i. Return prod g[i]^e[i] mod id.
2477 : * EX = multiple of exponent of (O_K/id)^* */
2478 : GEN
2479 2617190 : famat_to_nf_modideal_coprime(GEN nf, GEN g, GEN e, GEN id, GEN EX)
2480 : {
2481 2617190 : GEN EXo2, plus = NULL, minus = NULL, idZ = gcoeff(id,1,1);
2482 2617190 : long i, lx = lg(g);
2483 :
2484 2617190 : if (equali1(idZ)) return gen_1; /* id = Z_K */
2485 2616904 : EXo2 = (expi(EX) > 10)? shifti(EX,-1): NULL;
2486 8198119 : for (i = 1; i < lx; i++)
2487 : {
2488 5581349 : GEN h, n = centermodii(gel(e,i), EX, EXo2);
2489 5580839 : long sn = signe(n);
2490 5580839 : if (!sn) continue;
2491 :
2492 4044278 : h = nf_to_scalar_or_basis(nf, gel(g,i));
2493 4044755 : switch(typ(h))
2494 : {
2495 2359101 : case t_INT: break;
2496 0 : case t_FRAC:
2497 0 : h = Fp_div(gel(h,1), gel(h,2), idZ); break;
2498 1685654 : default:
2499 : {
2500 : GEN dh;
2501 1685654 : h = Q_remove_denom(h, &dh);
2502 1685767 : if (dh) h = FpC_Fp_mul(h, Fp_inv(dh,idZ), idZ);
2503 : }
2504 : }
2505 4044785 : if (sn > 0)
2506 4042967 : plus = nfmulpowmodideal(nf, plus, h, n, id);
2507 : else /* sn < 0 */
2508 1818 : minus = nfmulpowmodideal(nf, minus, h, negi(n), id);
2509 : }
2510 2616770 : if (minus) plus = nfmulmodideal(nf, plus, nfinvmodideal(nf,minus,id), id);
2511 2616862 : return plus? plus: gen_1;
2512 : }
2513 :
2514 : /* given 2 integral ideals x, y in HNF s.t x | y | x^2, compute (1+x)/(1+y) in
2515 : * the form [[cyc],[gen], U], where U := ux^-1 as a pair [ZM, denom(U)] */
2516 : static GEN
2517 236778 : zidealij(GEN x, GEN y)
2518 : {
2519 236778 : GEN U, G, cyc, xp = gcoeff(x,1,1), xi = hnf_invscale(x, xp);
2520 : long j, N;
2521 :
2522 : /* x^(-1) y = relations between the 1 + x_i (HNF) */
2523 236778 : cyc = ZM_snf_group(ZM_Z_divexact(ZM_mul(xi, y), xp), &U, &G);
2524 236775 : N = lg(cyc); G = ZM_mul(x,G); settyp(G, t_VEC); /* new generators */
2525 573506 : for (j=1; j<N; j++)
2526 : {
2527 336747 : GEN c = gel(G,j);
2528 336747 : gel(c,1) = addiu(gel(c,1), 1); /* 1 + g_j */
2529 336727 : if (ZV_isscalar(c)) gel(G,j) = gel(c,1);
2530 : }
2531 236759 : return mkvec4(cyc, G, ZM_mul(U,xi), xp);
2532 : }
2533 :
2534 : /* lg(x) > 1, x + 1; shallow */
2535 : static GEN
2536 169633 : ZC_add1(GEN x)
2537 : {
2538 169633 : long i, l = lg(x);
2539 169633 : GEN y = cgetg(l, t_COL);
2540 396093 : for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
2541 169637 : gel(y,1) = addiu(gel(x,1), 1); return y;
2542 : }
2543 : /* lg(x) > 1, x - 1; shallow */
2544 : static GEN
2545 70454 : ZC_sub1(GEN x)
2546 : {
2547 70454 : long i, l = lg(x);
2548 70454 : GEN y = cgetg(l, t_COL);
2549 176853 : for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
2550 70454 : gel(y,1) = subiu(gel(x,1), 1); return y;
2551 : }
2552 :
2553 : /* x,y are t_INT or ZC */
2554 : static GEN
2555 0 : zkadd(GEN x, GEN y)
2556 : {
2557 0 : long tx = typ(x);
2558 0 : if (tx == typ(y))
2559 0 : return tx == t_INT? addii(x,y): ZC_add(x,y);
2560 : else
2561 0 : return tx == t_INT? ZC_Z_add(y,x): ZC_Z_add(x,y);
2562 : }
2563 : /* x a t_INT or ZC, x+1; shallow */
2564 : static GEN
2565 255259 : zkadd1(GEN x)
2566 : {
2567 255259 : long tx = typ(x);
2568 255259 : return tx == t_INT? addiu(x,1): ZC_add1(x);
2569 : }
2570 : /* x a t_INT or ZC, x-1; shallow */
2571 : static GEN
2572 255322 : zksub1(GEN x)
2573 : {
2574 255322 : long tx = typ(x);
2575 255322 : return tx == t_INT? subiu(x,1): ZC_sub1(x);
2576 : }
2577 : /* x,y are t_INT or ZC; x - y */
2578 : static GEN
2579 0 : zksub(GEN x, GEN y)
2580 : {
2581 0 : long tx = typ(x), ty = typ(y);
2582 0 : if (tx == ty)
2583 0 : return tx == t_INT? subii(x,y): ZC_sub(x,y);
2584 : else
2585 0 : return tx == t_INT? Z_ZC_sub(x,y): ZC_Z_sub(x,y);
2586 : }
2587 : /* x is t_INT or ZM (mult. map), y is t_INT or ZC; x * y */
2588 : static GEN
2589 255275 : zkmul(GEN x, GEN y)
2590 : {
2591 255275 : long tx = typ(x), ty = typ(y);
2592 255275 : if (ty == t_INT)
2593 184832 : return tx == t_INT? mulii(x,y): ZC_Z_mul(gel(x,1),y);
2594 : else
2595 70443 : return tx == t_INT? ZC_Z_mul(y,x): ZM_ZC_mul(x,y);
2596 : }
2597 :
2598 : /* (U,V) = 1 coprime ideals. Want z = x mod U, = y mod V; namely
2599 : * z =vx + uy = v(x-y) + y, where u + v = 1, u in U, v in V.
2600 : * zkc = [v, UV], v a t_INT or ZM (mult. by v map), UV a ZM (ideal in HNF);
2601 : * shallow */
2602 : GEN
2603 0 : zkchinese(GEN zkc, GEN x, GEN y)
2604 : {
2605 0 : GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd(zkmul(v, zksub(x,y)), y);
2606 0 : return zk_modHNF(z, UV);
2607 : }
2608 : /* special case z = x mod U, = 1 mod V; shallow */
2609 : GEN
2610 255317 : zkchinese1(GEN zkc, GEN x)
2611 : {
2612 255317 : GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd1(zkmul(v, zksub1(x)));
2613 255277 : return (typ(z) == t_INT)? z: ZC_hnfrem(z, UV);
2614 : }
2615 : static GEN
2616 237334 : zkVchinese1(GEN zkc, GEN v)
2617 : {
2618 : long i, ly;
2619 237334 : GEN y = cgetg_copy(v, &ly);
2620 492600 : for (i=1; i<ly; i++) gel(y,i) = zkchinese1(zkc, gel(v,i));
2621 237287 : return y;
2622 : }
2623 :
2624 : /* prepare to solve z = x (mod A), z = y mod (B) [zkchinese or zkchinese1] */
2625 : GEN
2626 237082 : zkchineseinit(GEN nf, GEN A, GEN B, GEN AB)
2627 : {
2628 237082 : GEN v = idealaddtoone_raw(nf, A, B);
2629 : long e;
2630 237060 : if ((e = gexpo(v)) > 5)
2631 : {
2632 83218 : GEN b = (typ(v) == t_COL)? v: scalarcol_shallow(v, nf_get_degree(nf));
2633 83218 : b= ZC_reducemodlll(b, AB);
2634 83221 : if (gexpo(b) < e) v = b;
2635 : }
2636 237073 : return mkvec2(zk_scalar_or_multable(nf,v), AB);
2637 : }
2638 : /* prepare to solve z = x (mod A), z = 1 mod (B)
2639 : * and then z = 1 (mod A), z = y mod (B) [zkchinese1 twice] */
2640 : static GEN
2641 259 : zkchinese1init2(GEN nf, GEN A, GEN B, GEN AB)
2642 : {
2643 259 : GEN zkc = zkchineseinit(nf, A, B, AB);
2644 259 : GEN mv = gel(zkc,1), mu;
2645 259 : if (typ(mv) == t_INT) return mkvec2(zkc, mkvec2(subui(1,mv),AB));
2646 35 : mu = RgM_Rg_add_shallow(ZM_neg(mv), gen_1);
2647 35 : return mkvec2(mkvec2(mv,AB), mkvec2(mu,AB));
2648 : }
2649 :
2650 : static GEN
2651 2148123 : apply_U(GEN L, GEN a)
2652 : {
2653 2148123 : GEN e, U = gel(L,3), dU = gel(L,4);
2654 2148123 : if (typ(a) == t_INT)
2655 669814 : e = ZC_Z_mul(gel(U,1), subiu(a, 1));
2656 : else
2657 : { /* t_COL */
2658 1478309 : GEN t = shallowcopy(a);
2659 1478364 : gel(t,1) = subiu(gel(t,1), 1); /* t = a - 1 */
2660 1478275 : e = ZM_ZC_mul(U, t);
2661 : }
2662 2148086 : return gdiv(e, dU);
2663 : }
2664 :
2665 : /* true nf; vectors of [[cyc],[g],U.X^-1]. Assume k > 1. */
2666 : static GEN
2667 169051 : principal_units(GEN nf, GEN pr, long k, GEN prk)
2668 : {
2669 : GEN list, prb;
2670 169051 : ulong mask = quadratic_prec_mask(k);
2671 169051 : long a = 1;
2672 :
2673 169051 : prb = pr_hnf(nf,pr);
2674 169054 : list = vectrunc_init(k);
2675 405836 : while (mask > 1)
2676 : {
2677 236779 : GEN pra = prb;
2678 236779 : long b = a << 1;
2679 :
2680 236779 : if (mask & 1) b--;
2681 236779 : mask >>= 1;
2682 : /* compute 1 + pr^a / 1 + pr^b, 2a <= b */
2683 236779 : prb = (b >= k)? prk: idealpows(nf,pr,b);
2684 236779 : vectrunc_append(list, zidealij(pra, prb));
2685 236785 : a = b;
2686 : }
2687 169057 : return list;
2688 : }
2689 : /* a = 1 mod (pr) return log(a) on local-gens of 1+pr/1+pr^k */
2690 : static GEN
2691 1327732 : log_prk1(GEN nf, GEN a, long nh, GEN L2, GEN prk)
2692 : {
2693 1327732 : GEN y = cgetg(nh+1, t_COL);
2694 1327737 : long j, iy, c = lg(L2)-1;
2695 3475821 : for (j = iy = 1; j <= c; j++)
2696 : {
2697 2147964 : GEN L = gel(L2,j), cyc = gel(L,1), gen = gel(L,2), E = apply_U(L,a);
2698 2148065 : long i, nc = lg(cyc)-1;
2699 2148065 : int last = (j == c);
2700 5806209 : for (i = 1; i <= nc; i++, iy++)
2701 : {
2702 3658125 : GEN t, e = gel(E,i);
2703 3658125 : if (typ(e) != t_INT) pari_err_COPRIME("zlog_prk1", a, prk);
2704 3658118 : t = Fp_neg(e, gel(cyc,i));
2705 3658045 : gel(y,iy) = negi(t);
2706 3658159 : if (!last && signe(t)) a = nfmulpowmodideal(nf, a, gel(gen,i), t, prk);
2707 : }
2708 : }
2709 1327857 : return y;
2710 : }
2711 : /* true nf */
2712 : static GEN
2713 56614 : principal_units_relations(GEN nf, GEN L2, GEN prk, long nh)
2714 : {
2715 56614 : GEN h = cgetg(nh+1,t_MAT);
2716 56614 : long ih, j, c = lg(L2)-1;
2717 180953 : for (j = ih = 1; j <= c; j++)
2718 : {
2719 124338 : GEN L = gel(L2,j), F = gel(L,1), G = gel(L,2);
2720 124338 : long k, lG = lg(G);
2721 303998 : for (k = 1; k < lG; k++,ih++)
2722 : { /* log(g^f) mod pr^e */
2723 179659 : GEN a = nfpowmodideal(nf,gel(G,k),gel(F,k),prk);
2724 179662 : gel(h,ih) = ZC_neg(log_prk1(nf, a, nh, L2, prk));
2725 179660 : gcoeff(h,ih,ih) = gel(F,k);
2726 : }
2727 : }
2728 56615 : return h;
2729 : }
2730 : /* true nf; k > 1; multiplicative group (1 + pr) / (1 + pr^k) */
2731 : static GEN
2732 169059 : idealprincipalunits_i(GEN nf, GEN pr, long k, GEN *pU)
2733 : {
2734 169059 : GEN cyc, gen, L2, prk = idealpows(nf, pr, k);
2735 :
2736 169051 : L2 = principal_units(nf, pr, k, prk);
2737 169058 : if (k == 2)
2738 : {
2739 112444 : GEN L = gel(L2,1);
2740 112444 : cyc = gel(L,1);
2741 112444 : gen = gel(L,2);
2742 112444 : if (pU) *pU = matid(lg(gen)-1);
2743 : }
2744 : else
2745 : {
2746 56614 : long c = lg(L2), j;
2747 56614 : GEN EX, h, Ui, vg = cgetg(c, t_VEC);
2748 180954 : for (j = 1; j < c; j++) gel(vg, j) = gmael(L2,j,2);
2749 56615 : vg = shallowconcat1(vg);
2750 56614 : h = principal_units_relations(nf, L2, prk, lg(vg)-1);
2751 56616 : h = ZM_hnfall_i(h, NULL, 0);
2752 56616 : cyc = ZM_snf_group(h, pU, &Ui);
2753 56616 : c = lg(Ui); gen = cgetg(c, t_VEC); EX = cyc_get_expo(cyc);
2754 188156 : for (j = 1; j < c; j++)
2755 131542 : gel(gen,j) = famat_to_nf_modideal_coprime(nf, vg, gel(Ui,j), prk, EX);
2756 : }
2757 169059 : return mkvec4(cyc, gen, prk, L2);
2758 : }
2759 : GEN
2760 154 : idealprincipalunits(GEN nf, GEN pr, long k)
2761 : {
2762 : pari_sp av;
2763 : GEN v;
2764 154 : nf = checknf(nf);
2765 154 : if (k == 1) { checkprid(pr); retmkvec3(gen_1,cgetg(1,t_VEC),cgetg(1,t_VEC)); }
2766 147 : av = avma; v = idealprincipalunits_i(nf, pr, k, NULL);
2767 147 : return gerepilecopy(av, mkvec3(powiu(pr_norm(pr), k-1), gel(v,1), gel(v,2)));
2768 : }
2769 :
2770 : /* true nf; given an ideal pr^k dividing an integral ideal x (in HNF form)
2771 : * compute an 'sprk', the structure of G = (Z_K/pr^k)^* [ x = NULL for x=pr^k ]
2772 : * Return a vector with at least 4 components [cyc],[gen],[HNF pr^k,pr,k],ff,
2773 : * where
2774 : * cyc : type of G as abelian group (SNF)
2775 : * gen : generators of G, coprime to x
2776 : * pr^k: in HNF
2777 : * ff : data for log_g in (Z_K/pr)^*
2778 : * Two extra components are present iff k > 1: L2, U
2779 : * L2 : list of data structures to compute local DL in (Z_K/pr)^*,
2780 : * and 1 + pr^a/ 1 + pr^b for various a < b <= min(2a, k)
2781 : * U : base change matrices to convert a vector of local DL to DL wrt gen
2782 : * If MOD is not NULL, initialize G / G^MOD instead */
2783 : static GEN
2784 425655 : sprkinit(GEN nf, GEN pr, long k, GEN x, GEN MOD)
2785 : {
2786 425655 : GEN T, p, Ld, modpr, cyc, gen, g, g0, A, prk, U, L2, ord0 = NULL;
2787 425655 : long f = pr_get_f(pr);
2788 :
2789 425654 : if(DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",k,pr);
2790 425654 : modpr = nf_to_Fq_init(nf, &pr,&T,&p);
2791 425673 : if (MOD)
2792 : {
2793 378187 : GEN A = subiu(powiu(p,f), 1), d = gcdii(A, MOD), fa = Z_factor(d);
2794 378164 : ord0 = mkvec2(A, fa); /* true order, factorization of order in G/G^MOD */
2795 378149 : Ld = gel(fa,1);
2796 378149 : if (lg(Ld) > 1 && equaliu(gel(Ld,1),2)) Ld = vecslice(Ld,2,lg(Ld)-1);
2797 : }
2798 : /* (Z_K / pr)^* */
2799 425640 : if (f == 1)
2800 : {
2801 336532 : g0 = g = MOD? pgener_Fp_local(p, Ld): pgener_Fp(p);
2802 336564 : if (!ord0) ord0 = get_arith_ZZM(subiu(p,1));
2803 : }
2804 : else
2805 : {
2806 89108 : g0 = g = MOD? gener_FpXQ_local(T, p, Ld): gener_FpXQ(T,p, &ord0);
2807 89110 : g = Fq_to_nf(g, modpr);
2808 89108 : if (typ(g) == t_POL) g = poltobasis(nf, g);
2809 : }
2810 425677 : A = gel(ord0, 1); /* Norm(pr)-1 */
2811 : /* If MOD != NULL, d = gcd(A, MOD): g^(A/d) has order d */
2812 425677 : if (k == 1)
2813 : {
2814 256766 : cyc = mkvec(A);
2815 256762 : gen = mkvec(g);
2816 256757 : prk = pr_hnf(nf,pr);
2817 256778 : L2 = U = NULL;
2818 : }
2819 : else
2820 : { /* local-gens of (1 + pr)/(1 + pr^k) = SNF-gens * U */
2821 : GEN AB, B, u, v, w;
2822 : long j, l;
2823 168911 : w = idealprincipalunits_i(nf, pr, k, &U);
2824 : /* incorporate (Z_K/pr)^*, order A coprime to B = expo(1+pr/1+pr^k)*/
2825 168912 : cyc = leafcopy(gel(w,1)); B = cyc_get_expo(cyc); AB = mulii(A,B);
2826 168897 : gen = leafcopy(gel(w,2));
2827 168904 : prk = gel(w,3);
2828 168904 : g = nfpowmodideal(nf, g, B, prk);
2829 168907 : g0 = Fq_pow(g0, modii(B,A), T, p); /* update primitive root */
2830 168910 : L2 = mkvec3(A, g, gel(w,4));
2831 168913 : gel(cyc,1) = AB;
2832 168913 : gel(gen,1) = nfmulmodideal(nf, gel(gen,1), g, prk);
2833 168901 : u = mulii(Fp_inv(A,B), A);
2834 168896 : v = subui(1, u); l = lg(U);
2835 505177 : for (j = 1; j < l; j++) gcoeff(U,1,j) = Fp_mul(u, gcoeff(U,1,j), AB);
2836 168902 : U = mkvec2(Rg_col_ei(v, lg(gen)-1, 1), U);
2837 : }
2838 : /* local-gens of (Z_K/pr^k)^* = SNF-gens * U */
2839 425690 : if (x)
2840 : {
2841 236837 : GEN uv = zkchineseinit(nf, idealmulpowprime(nf,x,pr,utoineg(k)), prk, x);
2842 236817 : gen = zkVchinese1(uv, gen);
2843 : }
2844 425628 : return mkvecn(U? 6: 4, cyc, gen, prk, mkvec3(modpr,g0,ord0), L2, U);
2845 : }
2846 : GEN
2847 3976186 : sprk_get_cyc(GEN s) { return gel(s,1); }
2848 : GEN
2849 1965809 : sprk_get_expo(GEN s) { return cyc_get_expo(sprk_get_cyc(s)); }
2850 : GEN
2851 335624 : sprk_get_gen(GEN s) { return gel(s,2); }
2852 : GEN
2853 4907326 : sprk_get_prk(GEN s) { return gel(s,3); }
2854 : GEN
2855 2539173 : sprk_get_ff(GEN s) { return gel(s,4); }
2856 : GEN
2857 2599794 : sprk_get_pr(GEN s) { GEN ff = gel(s,4); return modpr_get_pr(gel(ff,1)); }
2858 : /* L2 to 1 + pr / 1 + pr^k */
2859 : static GEN
2860 1210230 : sprk_get_L2(GEN s) { return gmael(s,5,3); }
2861 : /* lift to nf of primitive root of k(pr) */
2862 : static GEN
2863 318446 : sprk_get_gnf(GEN s) { return gmael(s,5,2); }
2864 : /* A = Npr-1, <g> = (Z_K/pr)^*, L2 to 1 + pr / 1 + pr^k */
2865 : void
2866 0 : sprk_get_AgL2(GEN s, GEN *A, GEN *g, GEN *L2)
2867 0 : { GEN v = gel(s,5); *A = gel(v,1); *g = gel(v,2); *L2 = gel(v,3); }
2868 : void
2869 1201671 : sprk_get_U2(GEN s, GEN *U1, GEN *U2)
2870 1201671 : { GEN v = gel(s,6); *U1 = gel(v,1); *U2 = gel(v,2); }
2871 : static int
2872 2539186 : sprk_is_prime(GEN s) { return lg(s) == 5; }
2873 :
2874 : GEN
2875 1965613 : famat_zlog_pr(GEN nf, GEN g, GEN e, GEN sprk, GEN mod)
2876 : {
2877 1965613 : GEN x, expo = sprk_get_expo(sprk);
2878 1965613 : if (mod) expo = gcdii(expo,mod);
2879 1965590 : x = famat_makecoprime(nf, g, e, sprk_get_pr(sprk), sprk_get_prk(sprk), expo);
2880 1965611 : return log_prk(nf, x, sprk, mod);
2881 : }
2882 : /* famat_zlog_pr assuming (g,sprk.pr) = 1 */
2883 : static GEN
2884 196 : famat_zlog_pr_coprime(GEN nf, GEN g, GEN e, GEN sprk, GEN MOD)
2885 : {
2886 196 : GEN x = famat_to_nf_modideal_coprime(nf, g, e, sprk_get_prk(sprk),
2887 : sprk_get_expo(sprk));
2888 196 : return log_prk(nf, x, sprk, MOD);
2889 : }
2890 :
2891 : /* o t_INT, O = [ord,fa] format for multiple of o (for Fq_log);
2892 : * return o in [ord,fa] format */
2893 : static GEN
2894 557736 : order_update(GEN o, GEN O)
2895 : {
2896 557736 : GEN p = gmael(O,2,1), z = o, P, E;
2897 557736 : long i, j, l = lg(p);
2898 557736 : P = cgetg(l, t_COL);
2899 557732 : E = cgetg(l, t_COL);
2900 614347 : for (i = j = 1; i < l; i++)
2901 : {
2902 614347 : long v = Z_pvalrem(z, gel(p,i), &z);
2903 614270 : if (v)
2904 : {
2905 599411 : gel(P,j) = gel(p,i);
2906 599411 : gel(E,j) = utoipos(v); j++;
2907 599443 : if (is_pm1(z)) break;
2908 : }
2909 : }
2910 557690 : setlg(P, j);
2911 557688 : setlg(E, j); return mkvec2(o, mkmat2(P,E));
2912 : }
2913 :
2914 : /* a in Z_K (t_COL or t_INT), pr prime ideal, sprk = sprkinit(nf,pr,k,x),
2915 : * mod positive t_INT or NULL (meaning mod=0).
2916 : * return log(a) modulo mod on SNF-generators of (Z_K/pr^k)^* */
2917 : GEN
2918 2612946 : log_prk(GEN nf, GEN a, GEN sprk, GEN mod)
2919 : {
2920 : GEN e, prk, g, U1, U2, y, ff, O, o, oN, gN, N, T, p, modpr, pr, cyc;
2921 :
2922 2612946 : if (typ(a) == t_MAT) return famat_zlog_pr(nf, gel(a,1), gel(a,2), sprk, mod);
2923 2539155 : N = NULL;
2924 2539155 : ff = sprk_get_ff(sprk);
2925 2539173 : pr = gel(ff,1); /* modpr */
2926 2539173 : g = gN = gel(ff,2);
2927 2539173 : O = gel(ff,3); /* order of g = |Fq^*|, in [ord, fa] format */
2928 2539173 : o = oN = gel(O,1); /* order as a t_INT */
2929 2539173 : prk = sprk_get_prk(sprk);
2930 2539184 : modpr = nf_to_Fq_init(nf, &pr, &T, &p);
2931 2539162 : if (mod)
2932 : {
2933 2022871 : GEN d = gcdii(o,mod);
2934 2022618 : if (!equalii(o, d))
2935 : {
2936 747970 : N = diviiexact(o,d); /* > 1, coprime to p */
2937 747939 : a = nfpowmodideal(nf, a, N, prk);
2938 748059 : oN = d; /* order of g^N mod pr */
2939 : }
2940 : }
2941 2539006 : if (equali1(oN))
2942 396645 : e = gen_0;
2943 : else
2944 : {
2945 2142456 : if (N) { O = order_update(oN, O); gN = Fq_pow(g, N, T, p); }
2946 2142447 : e = Fq_log(nf_to_Fq(nf,a,modpr), gN, O, T, p);
2947 : }
2948 : /* 0 <= e < oN is correct modulo oN */
2949 2539184 : if (sprk_is_prime(sprk)) return mkcol(e); /* k = 1 */
2950 :
2951 799335 : sprk_get_U2(sprk, &U1,&U2);
2952 799406 : cyc = sprk_get_cyc(sprk);
2953 799409 : if (mod)
2954 : {
2955 378052 : cyc = ZV_snf_gcd(cyc, mod);
2956 378020 : if (signe(remii(mod,p))) return ZV_ZV_mod(ZC_Z_mul(U1,e), cyc);
2957 : }
2958 745798 : if (signe(e))
2959 : {
2960 318446 : GEN E = N? mulii(e, N): e;
2961 318446 : a = nfmulpowmodideal(nf, a, sprk_get_gnf(sprk), Fp_neg(E, o), prk);
2962 : }
2963 : /* a = 1 mod pr */
2964 745798 : y = log_prk1(nf, a, lg(U2)-1, sprk_get_L2(sprk), prk);
2965 745815 : if (N)
2966 : { /* from DL(a^N) to DL(a) */
2967 135234 : GEN E = gel(sprk_get_cyc(sprk), 1), q = powiu(p, Z_pval(E, p));
2968 135234 : y = ZC_Z_mul(y, Fp_inv(N, q));
2969 : }
2970 745814 : y = ZC_lincomb(gen_1, e, ZM_ZC_mul(U2,y), U1);
2971 745812 : return ZV_ZV_mod(y, cyc);
2972 : }
2973 : /* true nf */
2974 : GEN
2975 90117 : log_prk_init(GEN nf, GEN pr, long k, GEN MOD)
2976 90117 : { return sprkinit(nf,pr,k,NULL,MOD);}
2977 : GEN
2978 497 : veclog_prk(GEN nf, GEN v, GEN sprk)
2979 : {
2980 497 : long l = lg(v), i;
2981 497 : GEN w = cgetg(l, t_MAT);
2982 1232 : for (i = 1; i < l; i++) gel(w,i) = log_prk(nf, gel(v,i), sprk, NULL);
2983 497 : return w;
2984 : }
2985 :
2986 : static GEN
2987 1370810 : famat_zlog(GEN nf, GEN fa, GEN sgn, zlog_S *S)
2988 : {
2989 1370810 : long i, l0, l = lg(S->U);
2990 1370810 : GEN g = gel(fa,1), e = gel(fa,2), y = cgetg(l, t_COL);
2991 1370809 : l0 = lg(S->sprk); /* = l (trivial arch. part), or l-1 */
2992 2845152 : for (i=1; i < l0; i++) gel(y,i) = famat_zlog_pr(nf, g, e, gel(S->sprk,i), S->mod);
2993 1370808 : if (l0 != l)
2994 : {
2995 187696 : if (!sgn) sgn = nfsign_arch(nf, fa, S->archp);
2996 187696 : gel(y,l0) = Flc_to_ZC(sgn);
2997 : }
2998 1370808 : return y;
2999 : }
3000 :
3001 : /* assume that cyclic factors are normalized, in particular != [1] */
3002 : static GEN
3003 257221 : split_U(GEN U, GEN Sprk)
3004 : {
3005 257221 : long t = 0, k, n, l = lg(Sprk);
3006 257221 : GEN vU = cgetg(l+1, t_VEC);
3007 592075 : for (k = 1; k < l; k++)
3008 : {
3009 334844 : n = lg(sprk_get_cyc(gel(Sprk,k))) - 1; /* > 0 */
3010 334847 : gel(vU,k) = vecslice(U, t+1, t+n);
3011 334856 : t += n;
3012 : }
3013 : /* t+1 .. lg(U)-1 */
3014 257231 : n = lg(U) - t - 1; /* can be 0 */
3015 257231 : if (!n) setlg(vU,l); else gel(vU,l) = vecslice(U, t+1, t+n);
3016 257239 : return vU;
3017 : }
3018 :
3019 : static void
3020 1986297 : init_zlog_mod(zlog_S *S, GEN bid, GEN mod)
3021 : {
3022 1986297 : GEN fa2 = bid_get_fact2(bid);
3023 1986287 : S->U = bid_get_U(bid);
3024 1986284 : S->hU = lg(bid_get_cyc(bid))-1;
3025 1986276 : S->archp = bid_get_archp(bid);
3026 1986272 : S->sprk = bid_get_sprk(bid);
3027 1986273 : S->bid = bid;
3028 1986273 : S->mod = mod;
3029 1986273 : S->P = gel(fa2,1);
3030 1986273 : S->k = gel(fa2,2);
3031 1986273 : S->no2 = lg(S->P) == lg(gel(bid_get_fact(bid),1));
3032 1986285 : }
3033 : void
3034 379635 : init_zlog(zlog_S *S, GEN bid)
3035 : {
3036 379635 : return init_zlog_mod(S, bid, NULL);
3037 : }
3038 :
3039 : /* a a t_FRAC/t_INT, reduce mod bid */
3040 : static GEN
3041 14 : Q_mod_bid(GEN bid, GEN a)
3042 : {
3043 14 : GEN xZ = gcoeff(bid_get_ideal(bid),1,1);
3044 14 : GEN b = Rg_to_Fp(a, xZ);
3045 14 : if (gsigne(a) < 0) b = subii(b, xZ);
3046 14 : return signe(b)? b: xZ;
3047 : }
3048 : /* Return decomposition of a on the CRT generators blocks attached to the
3049 : * S->sprk and sarch; sgn = sign(a, S->arch), NULL if unknown */
3050 : static GEN
3051 380881 : zlog(GEN nf, GEN a, GEN sgn, zlog_S *S)
3052 : {
3053 : long k, l;
3054 : GEN y;
3055 380881 : a = nf_to_scalar_or_basis(nf, a);
3056 380874 : switch(typ(a))
3057 : {
3058 162142 : case t_INT: break;
3059 14 : case t_FRAC: a = Q_mod_bid(S->bid, a); break;
3060 218718 : default: /* case t_COL: */
3061 : {
3062 : GEN den;
3063 218718 : check_nfelt(a, &den);
3064 218738 : if (den)
3065 : {
3066 105 : a = Q_muli_to_int(a, den);
3067 105 : a = mkmat2(mkcol2(a, den), mkcol2(gen_1, gen_m1));
3068 105 : return famat_zlog(nf, a, sgn, S);
3069 : }
3070 : }
3071 : }
3072 380781 : if (sgn)
3073 374033 : sgn = (lg(sgn) == 1)? NULL: leafcopy(sgn);
3074 : else
3075 6748 : sgn = (lg(S->archp) == 1)? NULL: nfsign_arch(nf, a, S->archp);
3076 380782 : l = lg(S->sprk);
3077 380782 : y = cgetg(sgn? l+1: l, t_COL);
3078 921586 : for (k = 1; k < l; k++)
3079 : {
3080 540857 : GEN sprk = gel(S->sprk,k);
3081 540857 : gel(y,k) = log_prk(nf, a, sprk, S->mod);
3082 : }
3083 380729 : if (sgn) gel(y,l) = Flc_to_ZC(sgn);
3084 380737 : return y;
3085 : }
3086 :
3087 : /* true nf */
3088 : GEN
3089 43616 : pr_basis_perm(GEN nf, GEN pr)
3090 : {
3091 43616 : long f = pr_get_f(pr);
3092 : GEN perm;
3093 43616 : if (f == nf_get_degree(nf)) return identity_perm(f);
3094 38051 : perm = cgetg(f+1, t_VECSMALL);
3095 38051 : perm[1] = 1;
3096 38051 : if (f > 1)
3097 : {
3098 2912 : GEN H = pr_hnf(nf,pr);
3099 : long i, k;
3100 10808 : for (i = k = 2; k <= f; i++)
3101 7896 : if (!equali1(gcoeff(H,i,i))) perm[k++] = i;
3102 : }
3103 38051 : return perm;
3104 : }
3105 :
3106 : /* \sum U[i]*y[i], U[i] ZM, y[i] ZC. We allow lg(y) > lg(U). */
3107 : static GEN
3108 1751663 : ZMV_ZCV_mul(GEN U, GEN y)
3109 : {
3110 1751663 : long i, l = lg(U);
3111 1751663 : GEN z = NULL;
3112 1751663 : if (l == 1) return cgetg(1,t_COL);
3113 4128460 : for (i = 1; i < l; i++)
3114 : {
3115 2376889 : GEN u = ZM_ZC_mul(gel(U,i), gel(y,i));
3116 2376842 : z = z? ZC_add(z, u): u;
3117 : }
3118 1751571 : return z;
3119 : }
3120 : /* A * (U[1], ..., U[d] */
3121 : static GEN
3122 518 : ZM_ZMV_mul(GEN A, GEN U)
3123 : {
3124 : long i, l;
3125 518 : GEN V = cgetg_copy(U,&l);
3126 1057 : for (i = 1; i < l; i++) gel(V,i) = ZM_mul(A,gel(U,i));
3127 518 : return V;
3128 : }
3129 :
3130 : /* a = 1 mod pr, sprk mod pr^e, e >= 1 */
3131 : static GEN
3132 402283 : sprk_log_prk1_2(GEN nf, GEN a, GEN sprk)
3133 : {
3134 402283 : GEN U1, U2, y, L2 = sprk_get_L2(sprk);
3135 402283 : sprk_get_U2(sprk, &U1,&U2);
3136 402284 : y = ZM_ZC_mul(U2, log_prk1(nf, a, lg(U2)-1, L2, sprk_get_prk(sprk)));
3137 402288 : return ZV_ZV_mod(y, sprk_get_cyc(sprk));
3138 : }
3139 : /* true nf; assume e >= 2 */
3140 : GEN
3141 105529 : sprk_log_gen_pr2(GEN nf, GEN sprk, long e)
3142 : {
3143 105529 : GEN M, G, pr = sprk_get_pr(sprk);
3144 : long i, l;
3145 105530 : if (e == 2)
3146 : {
3147 62166 : GEN L2 = sprk_get_L2(sprk), L = gel(L2,1);
3148 62166 : G = gel(L,2); l = lg(G);
3149 : }
3150 : else
3151 : {
3152 43364 : GEN perm = pr_basis_perm(nf,pr), PI = nfpow_u(nf, pr_get_gen(pr), e-1);
3153 43364 : l = lg(perm);
3154 43364 : G = cgetg(l, t_VEC);
3155 43363 : if (typ(PI) == t_INT)
3156 : { /* zk_ei_mul doesn't allow t_INT */
3157 5557 : long N = nf_get_degree(nf);
3158 5556 : gel(G,1) = addiu(PI,1);
3159 8531 : for (i = 2; i < l; i++)
3160 : {
3161 2975 : GEN z = col_ei(N, 1);
3162 2975 : gel(G,i) = z; gel(z, perm[i]) = PI;
3163 : }
3164 : }
3165 : else
3166 : {
3167 37806 : gel(G,1) = nfadd(nf, gen_1, PI);
3168 44588 : for (i = 2; i < l; i++)
3169 6783 : gel(G,i) = nfadd(nf, gen_1, zk_ei_mul(nf, PI, perm[i]));
3170 : }
3171 : }
3172 105527 : M = cgetg(l, t_MAT);
3173 233670 : for (i = 1; i < l; i++) gel(M,i) = sprk_log_prk1_2(nf, gel(G,i), sprk);
3174 105522 : return M;
3175 : }
3176 : /* Log on bid.gen of generators of P_{1,I pr^{e-1}} / P_{1,I pr^e} (I,pr) = 1,
3177 : * defined implicitly via CRT. 'ind' is the index of pr in modulus
3178 : * factorization; true nf */
3179 : GEN
3180 413141 : log_gen_pr(zlog_S *S, long ind, GEN nf, long e)
3181 : {
3182 413141 : GEN Uind = gel(S->U, ind);
3183 413141 : if (e == 1) retmkmat( gel(Uind,1) );
3184 102823 : return ZM_mul(Uind, sprk_log_gen_pr2(nf, gel(S->sprk,ind), e));
3185 : }
3186 : /* true nf */
3187 : GEN
3188 2037 : sprk_log_gen_pr(GEN nf, GEN sprk, long e)
3189 : {
3190 2037 : if (e == 1)
3191 : {
3192 0 : long n = lg(sprk_get_cyc(sprk))-1;
3193 0 : retmkmat(col_ei(n, 1));
3194 : }
3195 2037 : return sprk_log_gen_pr2(nf, sprk, e);
3196 : }
3197 : /* a = 1 mod pr */
3198 : GEN
3199 274136 : sprk_log_prk1(GEN nf, GEN a, GEN sprk)
3200 : {
3201 274136 : if (lg(sprk) == 5) return mkcol(gen_0); /* mod pr */
3202 274136 : return sprk_log_prk1_2(nf, a, sprk);
3203 : }
3204 : /* Log on bid.gen of generator of P_{1,f} / P_{1,f v[index]}
3205 : * v = vector of r1 real places */
3206 : GEN
3207 85915 : log_gen_arch(zlog_S *S, long index) { return gel(veclast(S->U), index); }
3208 :
3209 : /* compute bid.clgp: [h,cyc] or [h,cyc,gen] */
3210 : static GEN
3211 258251 : bid_grp(GEN nf, GEN U, GEN cyc, GEN g, GEN F, GEN sarch)
3212 : {
3213 258251 : GEN G, h = ZV_prod(cyc);
3214 : long c;
3215 258266 : if (!U) return mkvec2(h,cyc);
3216 257923 : c = lg(U);
3217 257923 : G = cgetg(c,t_VEC);
3218 257924 : if (c > 1)
3219 : {
3220 227868 : GEN U0, Uoo, EX = cyc_get_expo(cyc); /* exponent of bid */
3221 227866 : long i, hU = nbrows(U), nba = lg(sarch_get_cyc(sarch))-1; /* #f_oo */
3222 227871 : if (!nba) { U0 = U; Uoo = NULL; }
3223 80311 : else if (nba == hU) { U0 = NULL; Uoo = U; }
3224 : else
3225 : {
3226 71218 : U0 = rowslice(U, 1, hU-nba);
3227 71222 : Uoo = rowslice(U, hU-nba+1, hU);
3228 : }
3229 694982 : for (i = 1; i < c; i++)
3230 : {
3231 467109 : GEN t = gen_1;
3232 467109 : if (U0) t = famat_to_nf_modideal_coprime(nf, g, gel(U0,i), F, EX);
3233 467101 : if (Uoo) t = set_sign_mod_divisor(nf, ZV_to_Flv(gel(Uoo,i),2), t, sarch);
3234 467109 : gel(G,i) = t;
3235 : }
3236 : }
3237 257929 : return mkvec3(h, cyc, G);
3238 : }
3239 :
3240 : /* remove prime ideals of norm 2 with exponent 1 from factorization */
3241 : static GEN
3242 258578 : famat_strip2(GEN fa)
3243 : {
3244 258578 : GEN P = gel(fa,1), E = gel(fa,2), Q, F;
3245 258578 : long l = lg(P), i, j;
3246 258578 : Q = cgetg(l, t_COL);
3247 258574 : F = cgetg(l, t_COL);
3248 633456 : for (i = j = 1; i < l; i++)
3249 : {
3250 374862 : GEN pr = gel(P,i), e = gel(E,i);
3251 374862 : if (!absequaliu(pr_get_p(pr), 2) || itou(e) != 1 || pr_get_f(pr) != 1)
3252 : {
3253 336226 : gel(Q,j) = pr;
3254 336226 : gel(F,j) = e; j++;
3255 : }
3256 : }
3257 258594 : setlg(Q,j);
3258 258592 : setlg(F,j); return mkmat2(Q,F);
3259 : }
3260 : static int
3261 134050 : checkarchp(GEN v, long r1)
3262 : {
3263 134050 : long i, l = lg(v);
3264 134050 : pari_sp av = avma;
3265 : GEN p;
3266 134050 : if (l == 1) return 1;
3267 47121 : if (l == 2) return v[1] > 0 && v[1] <= r1;
3268 22014 : p = zero_zv(r1);
3269 66147 : for (i = 1; i < l; i++)
3270 : {
3271 44161 : long j = v[i];
3272 44161 : if (j <= 0 || j > r1 || p[j]) return gc_long(av, 0);
3273 44126 : p[j] = 1;
3274 : }
3275 21986 : return gc_long(av, 1);
3276 : }
3277 :
3278 : /* True nf. Put ideal to form [[ideal,arch]] and set fa and fa2 to its
3279 : * factorization, archp to the indices of arch places */
3280 : GEN
3281 258566 : check_mod_factored(GEN nf, GEN ideal, GEN *fa_, GEN *fa2_, GEN *archp_, GEN MOD)
3282 : {
3283 : GEN arch, x, fa, fa2, archp;
3284 : long R1;
3285 :
3286 258566 : R1 = nf_get_r1(nf);
3287 258563 : if (typ(ideal) == t_VEC && lg(ideal) == 3)
3288 : {
3289 178521 : arch = gel(ideal,2);
3290 178521 : ideal= gel(ideal,1);
3291 178521 : switch(typ(arch))
3292 : {
3293 44471 : case t_VEC:
3294 44471 : if (lg(arch) != R1+1)
3295 7 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
3296 44464 : archp = vec01_to_indices(arch);
3297 44464 : break;
3298 134050 : case t_VECSMALL:
3299 134050 : if (!checkarchp(arch, R1))
3300 35 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
3301 134018 : archp = arch;
3302 134018 : arch = indices_to_vec01(archp, R1);
3303 134016 : break;
3304 0 : default:
3305 0 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
3306 : return NULL;/*LCOV_EXCL_LINE*/
3307 : }
3308 : }
3309 : else
3310 : {
3311 80042 : arch = zerovec(R1);
3312 80035 : archp = cgetg(1, t_VECSMALL);
3313 : }
3314 258512 : if (MOD)
3315 : {
3316 213999 : if (typ(MOD) != t_INT) pari_err_TYPE("bnrinit [incorrect cycmod]", MOD);
3317 213999 : if (mpodd(MOD) && lg(archp) != 1)
3318 224 : MOD = shifti(MOD, 1); /* ensure elements of G^MOD are >> 0 */
3319 : }
3320 258526 : if (is_nf_factor(ideal))
3321 : {
3322 40166 : fa = ideal;
3323 40166 : x = factorbackprime(nf, gel(fa,1), gel(fa,2));
3324 : }
3325 : else
3326 : {
3327 218356 : fa = idealfactor(nf, ideal);
3328 218375 : x = ideal;
3329 : }
3330 258541 : if (typ(x) != t_MAT) x = idealhnf_shallow(nf, x);
3331 258523 : if (lg(x) == 1) pari_err_DOMAIN("Idealstar", "ideal","=",gen_0,x);
3332 258523 : if (typ(gcoeff(x,1,1)) != t_INT)
3333 7 : pari_err_DOMAIN("Idealstar","denominator(ideal)", "!=",gen_1,x);
3334 :
3335 258516 : fa2 = famat_strip2(fa);
3336 258516 : if (fa_ != NULL) *fa_ = fa;
3337 258516 : if (fa2_ != NULL) *fa2_ = fa2;
3338 258516 : if (fa2_ != NULL) *archp_ = archp;
3339 258516 : return mkvec2(x, arch);
3340 : }
3341 :
3342 : /* True nf. Compute [[ideal,arch], [h,[cyc],[gen]], idealfact, [liste], U]
3343 : flag may include nf_GEN | nf_INIT */
3344 : static GEN
3345 257930 : Idealstarmod_i(GEN nf, GEN ideal, long flag, GEN MOD)
3346 : {
3347 : long i, nbp;
3348 257930 : GEN y, cyc, U, u1 = NULL, fa, fa2, sprk, x_arch, x, arch, archp, E, P, sarch, gen;
3349 :
3350 257930 : x_arch = check_mod_factored(nf, ideal, &fa, &fa2, &archp, MOD);
3351 257884 : x = gel(x_arch, 1);
3352 257884 : arch = gel(x_arch, 2);
3353 :
3354 257884 : sarch = nfarchstar(nf, x, archp);
3355 257895 : P = gel(fa2,1);
3356 257895 : E = gel(fa2,2);
3357 257895 : nbp = lg(P)-1;
3358 257895 : sprk = cgetg(nbp+1,t_VEC);
3359 257897 : if (nbp)
3360 : {
3361 218706 : GEN t = (lg(gel(fa,1))==2)? NULL: x; /* beware fa != fa2 */
3362 218706 : cyc = cgetg(nbp+2,t_VEC);
3363 218696 : gen = cgetg(nbp+1,t_VEC);
3364 554278 : for (i = 1; i <= nbp; i++)
3365 : {
3366 335551 : GEN L = sprkinit(nf, gel(P,i), itou(gel(E,i)), t, MOD);
3367 335559 : gel(sprk,i) = L;
3368 335559 : gel(cyc,i) = sprk_get_cyc(L);
3369 : /* true gens are congruent to those mod x AND positive at archp */
3370 335561 : gel(gen,i) = sprk_get_gen(L);
3371 : }
3372 218727 : gel(cyc,i) = sarch_get_cyc(sarch);
3373 218728 : cyc = shallowconcat1(cyc);
3374 218730 : gen = shallowconcat1(gen);
3375 218735 : cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
3376 : }
3377 : else
3378 : {
3379 39191 : cyc = sarch_get_cyc(sarch);
3380 39192 : gen = cgetg(1,t_VEC);
3381 39192 : U = matid(lg(cyc)-1);
3382 39192 : if (flag & nf_GEN) u1 = U;
3383 : }
3384 257888 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
3385 257902 : if (!(flag & nf_INIT)) return y;
3386 257118 : U = split_U(U, sprk);
3387 257126 : return mkvec5(mkvec2(x, arch), y, mkvec2(fa,fa2), mkvec2(sprk, sarch), U);
3388 : }
3389 :
3390 : static long
3391 63 : idealHNF_norm_pval(GEN x, GEN p)
3392 : {
3393 63 : long i, v = 0, l = lg(x);
3394 175 : for (i = 1; i < l; i++) v += Z_pval(gcoeff(x,i,i), p);
3395 63 : return v;
3396 : }
3397 : static long
3398 63 : sprk_get_k(GEN sprk)
3399 : {
3400 : GEN pr, prk;
3401 63 : if (sprk_is_prime(sprk)) return 1;
3402 63 : pr = sprk_get_pr(sprk);
3403 63 : prk = sprk_get_prk(sprk);
3404 63 : return idealHNF_norm_pval(prk, pr_get_p(pr)) / pr_get_f(pr);
3405 : }
3406 : /* true nf, L a sprk */
3407 : GEN
3408 63 : sprk_to_bid(GEN nf, GEN L, long flag)
3409 : {
3410 63 : GEN y, cyc, U, u1 = NULL, fa, fa2, arch, sarch, gen, sprk;
3411 :
3412 63 : arch = zerovec(nf_get_r1(nf));
3413 63 : fa = to_famat_shallow(sprk_get_pr(L), utoipos(sprk_get_k(L)));
3414 63 : sarch = nfarchstar(nf, NULL, cgetg(1, t_VECSMALL));
3415 63 : fa2 = famat_strip2(fa);
3416 63 : sprk = mkvec(L);
3417 63 : cyc = shallowconcat(sprk_get_cyc(L), sarch_get_cyc(sarch));
3418 63 : gen = sprk_get_gen(L);
3419 63 : cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
3420 63 : y = bid_grp(nf, u1, cyc, gen, NULL, sarch);
3421 63 : if (!(flag & nf_INIT)) return y;
3422 63 : return mkvec5(mkvec2(sprk_get_prk(L), arch), y, mkvec2(fa,fa2),
3423 : mkvec2(sprk, sarch), split_U(U, sprk));
3424 : }
3425 : GEN
3426 257654 : Idealstarmod(GEN nf, GEN ideal, long flag, GEN MOD)
3427 : {
3428 257654 : pari_sp av = avma;
3429 257654 : nf = nf? checknf(nf): nfinit(pol_x(0), DEFAULTPREC);
3430 257656 : return gerepilecopy(av, Idealstarmod_i(nf, ideal, flag, MOD));
3431 : }
3432 : GEN
3433 924 : Idealstar(GEN nf, GEN ideal, long flag) { return Idealstarmod(nf, ideal, flag, NULL); }
3434 : GEN
3435 273 : Idealstarprk(GEN nf, GEN pr, long k, long flag)
3436 : {
3437 273 : pari_sp av = avma;
3438 273 : GEN z = Idealstarmod_i(nf, mkmat2(mkcol(pr),mkcols(k)), flag, NULL);
3439 273 : return gerepilecopy(av, z);
3440 : }
3441 :
3442 : /* FIXME: obsolete */
3443 : GEN
3444 0 : zidealstarinitgen(GEN nf, GEN ideal)
3445 0 : { return Idealstar(nf,ideal, nf_INIT|nf_GEN); }
3446 : GEN
3447 0 : zidealstarinit(GEN nf, GEN ideal)
3448 0 : { return Idealstar(nf,ideal, nf_INIT); }
3449 : GEN
3450 0 : zidealstar(GEN nf, GEN ideal)
3451 0 : { return Idealstar(nf,ideal, nf_GEN); }
3452 :
3453 : GEN
3454 98 : idealstarmod(GEN nf, GEN ideal, long flag, GEN MOD)
3455 : {
3456 98 : switch(flag)
3457 : {
3458 0 : case 0: return Idealstarmod(nf,ideal, nf_GEN, MOD);
3459 84 : case 1: return Idealstarmod(nf,ideal, nf_INIT, MOD);
3460 14 : case 2: return Idealstarmod(nf,ideal, nf_INIT|nf_GEN, MOD);
3461 0 : default: pari_err_FLAG("idealstar");
3462 : }
3463 : return NULL; /* LCOV_EXCL_LINE */
3464 : }
3465 : GEN
3466 0 : idealstar0(GEN nf, GEN ideal,long flag) { return idealstarmod(nf, ideal, flag, NULL); }
3467 :
3468 : void
3469 218733 : check_nfelt(GEN x, GEN *den)
3470 : {
3471 218733 : long l = lg(x), i;
3472 218733 : GEN t, d = NULL;
3473 218733 : if (typ(x) != t_COL) pari_err_TYPE("check_nfelt", x);
3474 808201 : for (i=1; i<l; i++)
3475 : {
3476 589461 : t = gel(x,i);
3477 589461 : switch (typ(t))
3478 : {
3479 589244 : case t_INT: break;
3480 217 : case t_FRAC:
3481 217 : if (!d) d = gel(t,2); else d = lcmii(d, gel(t,2));
3482 217 : break;
3483 0 : default: pari_err_TYPE("check_nfelt", x);
3484 : }
3485 : }
3486 218740 : *den = d;
3487 218740 : }
3488 :
3489 : GEN
3490 2141361 : ZV_snf_gcd(GEN x, GEN mod)
3491 4939102 : { pari_APPLY_same(gcdii(gel(x,i), mod)); }
3492 :
3493 : /* assume a true bnf and bid */
3494 : GEN
3495 226847 : ideallog_units0(GEN bnf, GEN bid, GEN MOD)
3496 : {
3497 226847 : GEN nf = bnf_get_nf(bnf), D, y, C, cyc;
3498 226848 : long j, lU = lg(bnf_get_logfu(bnf)); /* r1+r2 */
3499 : zlog_S S;
3500 226852 : init_zlog_mod(&S, bid, MOD);
3501 226860 : if (!S.hU) return zeromat(0,lU);
3502 226860 : cyc = bid_get_cyc(bid);
3503 226856 : if (MOD) cyc = ZV_snf_gcd(cyc, MOD);
3504 226822 : D = nfsign_fu(bnf, bid_get_archp(bid));
3505 226853 : y = cgetg(lU, t_MAT);
3506 226847 : if ((C = bnf_build_cheapfu(bnf)))
3507 373978 : { for (j = 1; j < lU; j++) gel(y,j) = zlog(nf, gel(C,j), gel(D,j), &S); }
3508 : else
3509 : {
3510 49 : long i, l = lg(S.U), l0 = lg(S.sprk);
3511 : GEN X, U;
3512 49 : if (!(C = bnf_compactfu_mat(bnf))) bnf_build_units(bnf); /* error */
3513 49 : X = gel(C,1); U = gel(C,2);
3514 147 : for (j = 1; j < lU; j++) gel(y,j) = cgetg(l, t_COL);
3515 126 : for (i = 1; i < l0; i++)
3516 : {
3517 77 : GEN sprk = gel(S.sprk, i);
3518 77 : GEN Xi = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
3519 231 : for (j = 1; j < lU; j++)
3520 154 : gcoeff(y,i,j) = famat_zlog_pr_coprime(nf, Xi, gel(U,j), sprk, MOD);
3521 : }
3522 49 : if (l0 != l)
3523 56 : for (j = 1; j < lU; j++) gcoeff(y,l0,j) = Flc_to_ZC(gel(D,j));
3524 : }
3525 226849 : y = vec_prepend(y, zlog(nf, bnf_get_tuU(bnf), nfsign_tu(bnf, S.archp), &S));
3526 600966 : for (j = 1; j <= lU; j++)
3527 374122 : gel(y,j) = ZV_ZV_mod(ZMV_ZCV_mul(S.U, gel(y,j)), cyc);
3528 226844 : return y;
3529 : }
3530 : GEN
3531 84 : ideallog_units(GEN bnf, GEN bid)
3532 84 : { return ideallog_units0(bnf, bid, NULL); }
3533 : GEN
3534 518 : log_prk_units(GEN nf, GEN D, GEN sprk)
3535 : {
3536 518 : GEN L, Ltu = log_prk(nf, gel(D,1), sprk, NULL);
3537 518 : D = gel(D,2);
3538 518 : if (lg(D) != 3 || typ(gel(D,2)) != t_MAT) L = veclog_prk(nf, D, sprk);
3539 : else
3540 : {
3541 21 : GEN X = gel(D,1), U = gel(D,2);
3542 21 : long j, lU = lg(U);
3543 21 : X = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
3544 21 : L = cgetg(lU, t_MAT);
3545 63 : for (j = 1; j < lU; j++)
3546 42 : gel(L,j) = famat_zlog_pr_coprime(nf, X, gel(U,j), sprk, NULL);
3547 : }
3548 518 : return vec_prepend(L, Ltu);
3549 : }
3550 :
3551 : static GEN
3552 1379819 : ideallog_i(GEN nf, GEN x, zlog_S *S)
3553 : {
3554 1379819 : pari_sp av = avma;
3555 : GEN y;
3556 1379819 : if (!S->hU) return cgetg(1, t_COL);
3557 1377558 : if (typ(x) == t_MAT)
3558 1370705 : y = famat_zlog(nf, x, NULL, S);
3559 : else
3560 6853 : y = zlog(nf, x, NULL, S);
3561 1377551 : y = ZMV_ZCV_mul(S->U, y);
3562 1377550 : return gerepileupto(av, ZV_ZV_mod(y, bid_get_cyc(S->bid)));
3563 : }
3564 : GEN
3565 1386498 : ideallogmod(GEN nf, GEN x, GEN bid, GEN mod)
3566 : {
3567 : zlog_S S;
3568 1386498 : if (!nf)
3569 : {
3570 6671 : if (mod) pari_err_IMPL("Zideallogmod");
3571 6671 : return Zideallog(bid, x);
3572 : }
3573 1379827 : checkbid(bid); init_zlog_mod(&S, bid, mod);
3574 1379820 : return ideallog_i(checknf(nf), x, &S);
3575 : }
3576 : GEN
3577 13755 : ideallog(GEN nf, GEN x, GEN bid) { return ideallogmod(nf, x, bid, NULL); }
3578 :
3579 : /*************************************************************************/
3580 : /** **/
3581 : /** JOIN BID STRUCTURES, IDEAL LISTS **/
3582 : /** **/
3583 : /*************************************************************************/
3584 : /* bid1, bid2: for coprime modules m1 and m2 (without arch. part).
3585 : * Output: bid for m1 m2 */
3586 : static GEN
3587 469 : join_bid(GEN nf, GEN bid1, GEN bid2)
3588 : {
3589 469 : pari_sp av = avma;
3590 : long nbgen, l1,l2;
3591 : GEN I1,I2, G1,G2, sprk1,sprk2, cyc1,cyc2, sarch;
3592 469 : GEN sprk, fa,fa2, U, cyc, y, u1 = NULL, x, gen;
3593 :
3594 469 : I1 = bid_get_ideal(bid1);
3595 469 : I2 = bid_get_ideal(bid2);
3596 469 : if (gequal1(gcoeff(I1,1,1))) return bid2; /* frequent trivial case */
3597 259 : G1 = bid_get_grp(bid1);
3598 259 : G2 = bid_get_grp(bid2);
3599 259 : x = idealmul(nf, I1,I2);
3600 259 : fa = famat_mul_shallow(bid_get_fact(bid1), bid_get_fact(bid2));
3601 259 : fa2= famat_mul_shallow(bid_get_fact2(bid1), bid_get_fact2(bid2));
3602 259 : sprk1 = bid_get_sprk(bid1);
3603 259 : sprk2 = bid_get_sprk(bid2);
3604 259 : sprk = shallowconcat(sprk1, sprk2);
3605 :
3606 259 : cyc1 = abgrp_get_cyc(G1); l1 = lg(cyc1);
3607 259 : cyc2 = abgrp_get_cyc(G2); l2 = lg(cyc2);
3608 259 : gen = (lg(G1)>3 && lg(G2)>3)? gen_1: NULL;
3609 259 : nbgen = l1+l2-2;
3610 259 : cyc = ZV_snf_group(shallowconcat(cyc1,cyc2), &U, gen? &u1: NULL);
3611 259 : if (nbgen)
3612 : {
3613 259 : GEN U1 = bid_get_U(bid1), U2 = bid_get_U(bid2);
3614 0 : U1 = l1==1? const_vec(lg(sprk1), cgetg(1,t_MAT))
3615 259 : : ZM_ZMV_mul(vecslice(U, 1, l1-1), U1);
3616 0 : U2 = l2==1? const_vec(lg(sprk2), cgetg(1,t_MAT))
3617 259 : : ZM_ZMV_mul(vecslice(U, l1, nbgen), U2);
3618 259 : U = shallowconcat(U1, U2);
3619 : }
3620 : else
3621 0 : U = const_vec(lg(sprk), cgetg(1,t_MAT));
3622 :
3623 259 : if (gen)
3624 : {
3625 259 : GEN uv = zkchinese1init2(nf, I2, I1, x);
3626 518 : gen = shallowconcat(zkVchinese1(gel(uv,1), abgrp_get_gen(G1)),
3627 259 : zkVchinese1(gel(uv,2), abgrp_get_gen(G2)));
3628 : }
3629 259 : sarch = bid_get_sarch(bid1); /* trivial */
3630 259 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
3631 259 : x = mkvec2(x, bid_get_arch(bid1));
3632 259 : y = mkvec5(x, y, mkvec2(fa, fa2), mkvec2(sprk, sarch), U);
3633 259 : return gerepilecopy(av,y);
3634 : }
3635 :
3636 : typedef struct _ideal_data {
3637 : GEN nf, emb, L, pr, prL, sgnU, archp;
3638 : } ideal_data;
3639 :
3640 : /* z <- ( z | f(v[i])_{i=1..#v} ) */
3641 : static void
3642 758808 : concat_join(GEN *pz, GEN v, GEN (*f)(ideal_data*,GEN), ideal_data *data)
3643 : {
3644 758808 : long i, nz, lv = lg(v);
3645 : GEN z, Z;
3646 758808 : if (lv == 1) return;
3647 222852 : z = *pz; nz = lg(z)-1;
3648 222852 : *pz = Z = cgetg(lv + nz, typ(z));
3649 371696 : for (i = 1; i <=nz; i++) gel(Z,i) = gel(z,i);
3650 223353 : Z += nz;
3651 492071 : for (i = 1; i < lv; i++) gel(Z,i) = f(data, gel(v,i));
3652 : }
3653 : static GEN
3654 469 : join_idealinit(ideal_data *D, GEN x)
3655 469 : { return join_bid(D->nf, x, D->prL); }
3656 : static GEN
3657 268533 : join_ideal(ideal_data *D, GEN x)
3658 268533 : { return idealmulpowprime(D->nf, x, D->pr, D->L); }
3659 : static GEN
3660 448 : join_unit(ideal_data *D, GEN x)
3661 : {
3662 448 : GEN bid = join_idealinit(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
3663 448 : if (lg(u) != 1) v = shallowconcat(u, v);
3664 448 : return mkvec2(bid, v);
3665 : }
3666 :
3667 : GEN
3668 49 : log_prk_units_init(GEN bnf)
3669 : {
3670 49 : GEN U = bnf_has_fu(bnf);
3671 49 : if (U) U = matalgtobasis(bnf_get_nf(bnf), U);
3672 21 : else if (!(U = bnf_compactfu_mat(bnf))) (void)bnf_build_units(bnf);
3673 49 : return mkvec2(bnf_get_tuU(bnf), U);
3674 : }
3675 : /* flag & nf_GEN : generators, otherwise no
3676 : * flag &2 : units, otherwise no
3677 : * flag &4 : ideals in HNF, otherwise bid
3678 : * flag &8 : omit ideals which cannot be conductors (pr^1 with Npr=2) */
3679 : static GEN
3680 11333 : Ideallist(GEN bnf, ulong bound, long flag)
3681 : {
3682 11333 : const long do_units = flag & 2, big_id = !(flag & 4), cond = flag & 8;
3683 11333 : const long istar_flag = (flag & nf_GEN) | nf_INIT;
3684 : pari_sp av;
3685 : long i, j;
3686 11333 : GEN nf, z, p, fa, id, BOUND, U, empty = cgetg(1,t_VEC);
3687 : forprime_t S;
3688 : ideal_data ID;
3689 : GEN (*join_z)(ideal_data*, GEN);
3690 :
3691 11333 : if (do_units)
3692 : {
3693 21 : bnf = checkbnf(bnf);
3694 21 : nf = bnf_get_nf(bnf);
3695 21 : join_z = &join_unit;
3696 : }
3697 : else
3698 : {
3699 11312 : nf = checknf(bnf);
3700 11312 : join_z = big_id? &join_idealinit: &join_ideal;
3701 : }
3702 11333 : if ((long)bound <= 0) return empty;
3703 11333 : id = matid(nf_get_degree(nf));
3704 11333 : if (big_id) id = Idealstar(nf,id, istar_flag);
3705 :
3706 : /* z[i] will contain all "objects" of norm i. Depending on flag, this means
3707 : * an ideal, a bid, or a couple [bid, log(units)]. Such objects are stored
3708 : * in vectors, computed one primary component at a time; join_z
3709 : * reconstructs the global object */
3710 11333 : BOUND = utoipos(bound);
3711 11333 : z = const_vec(bound, empty);
3712 11333 : U = do_units? log_prk_units_init(bnf): NULL;
3713 11333 : gel(z,1) = mkvec(U? mkvec2(id, empty): id);
3714 11333 : ID.nf = nf;
3715 :
3716 11333 : p = cgetipos(3);
3717 11333 : u_forprime_init(&S, 2, bound);
3718 11333 : av = avma;
3719 92825 : while ((p[2] = u_forprime_next(&S)))
3720 : {
3721 81613 : if (DEBUGLEVEL>1) err_printf("%ld ",p[2]);
3722 81613 : fa = idealprimedec_limit_norm(nf, p, BOUND);
3723 163002 : for (j = 1; j < lg(fa); j++)
3724 : {
3725 81510 : GEN pr = gel(fa,j), z2 = leafcopy(z);
3726 81512 : ulong Q, q = upr_norm(pr);
3727 : long l;
3728 81512 : ID.pr = ID.prL = pr;
3729 81512 : if (cond && q == 2) { l = 2; Q = 4; } else { l = 1; Q = q; }
3730 184430 : for (; Q <= bound; l++, Q *= q) /* add pr^l */
3731 : {
3732 : ulong iQ;
3733 103043 : ID.L = utoipos(l);
3734 103042 : if (big_id) {
3735 210 : ID.prL = Idealstarprk(nf, pr, l, istar_flag);
3736 210 : if (U)
3737 189 : ID.emb = Q == 2? empty
3738 189 : : log_prk_units(nf, U, gel(bid_get_sprk(ID.prL),1));
3739 : }
3740 861618 : for (iQ = Q,i = 1; iQ <= bound; iQ += Q,i++)
3741 758700 : concat_join(&gel(z,iQ), gel(z2,i), join_z, &ID);
3742 : }
3743 : }
3744 81492 : if (gc_needed(av,1))
3745 : {
3746 18 : if(DEBUGMEM>1) pari_warn(warnmem,"Ideallist");
3747 18 : z = gerepilecopy(av, z);
3748 : }
3749 : }
3750 11333 : return z;
3751 : }
3752 : GEN
3753 63 : gideallist(GEN bnf, GEN B, long flag)
3754 : {
3755 63 : pari_sp av = avma;
3756 63 : if (typ(B) != t_INT)
3757 : {
3758 0 : B = gfloor(B);
3759 0 : if (typ(B) != t_INT) pari_err_TYPE("ideallist", B);
3760 0 : if (signe(B) < 0) B = gen_0;
3761 : }
3762 63 : if (signe(B) < 0)
3763 : {
3764 28 : if (flag != 4) pari_err_IMPL("ideallist with bid for single norm");
3765 28 : return gerepilecopy(av, ideals_by_norm(checknf(bnf), absi(B)));
3766 : }
3767 35 : if (flag < 0 || flag > 15) pari_err_FLAG("ideallist");
3768 35 : return gerepilecopy(av, Ideallist(bnf, itou(B), flag));
3769 : }
3770 : GEN
3771 11298 : ideallist0(GEN bnf, long bound, long flag)
3772 : {
3773 11298 : pari_sp av = avma;
3774 11298 : if (flag < 0 || flag > 15) pari_err_FLAG("ideallist");
3775 11298 : return gerepilecopy(av, Ideallist(bnf, bound, flag));
3776 : }
3777 : GEN
3778 10563 : ideallist(GEN bnf,long bound) { return ideallist0(bnf,bound,4); }
3779 :
3780 : /* bid = for module m (without arch. part), arch = archimedean part.
3781 : * Output: bid for [m,arch] */
3782 : static GEN
3783 42 : join_bid_arch(GEN nf, GEN bid, GEN archp)
3784 : {
3785 42 : pari_sp av = avma;
3786 : GEN G, U;
3787 42 : GEN sprk, cyc, y, u1 = NULL, x, sarch, gen;
3788 :
3789 42 : checkbid(bid);
3790 42 : G = bid_get_grp(bid);
3791 42 : x = bid_get_ideal(bid);
3792 42 : sarch = nfarchstar(nf, bid_get_ideal(bid), archp);
3793 42 : sprk = bid_get_sprk(bid);
3794 :
3795 42 : gen = (lg(G)>3)? gel(G,3): NULL;
3796 42 : cyc = diagonal_shallow(shallowconcat(gel(G,2), sarch_get_cyc(sarch)));
3797 42 : cyc = ZM_snf_group(cyc, &U, gen? &u1: NULL);
3798 42 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
3799 42 : U = split_U(U, sprk);
3800 42 : y = mkvec5(mkvec2(x, archp), y, gel(bid,3), mkvec2(sprk, sarch), U);
3801 42 : return gerepilecopy(av,y);
3802 : }
3803 : static GEN
3804 42 : join_arch(ideal_data *D, GEN x) {
3805 42 : return join_bid_arch(D->nf, x, D->archp);
3806 : }
3807 : static GEN
3808 14 : join_archunit(ideal_data *D, GEN x) {
3809 14 : GEN bid = join_arch(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
3810 14 : if (lg(u) != 1) v = shallowconcat(u, v);
3811 14 : return mkvec2(bid, v);
3812 : }
3813 :
3814 : /* L from ideallist, add archimedean part */
3815 : GEN
3816 14 : ideallistarch(GEN bnf, GEN L, GEN arch)
3817 : {
3818 : pari_sp av;
3819 14 : long i, j, l = lg(L), lz;
3820 : GEN v, z, V, nf;
3821 : ideal_data ID;
3822 : GEN (*join_z)(ideal_data*, GEN);
3823 :
3824 14 : if (typ(L) != t_VEC) pari_err_TYPE("ideallistarch",L);
3825 14 : if (l == 1) return cgetg(1,t_VEC);
3826 14 : z = gel(L,1);
3827 14 : if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
3828 14 : z = gel(z,1); /* either a bid or [bid,U] */
3829 14 : ID.archp = vec01_to_indices(arch);
3830 14 : if (lg(z) == 3)
3831 : { /* [bid,U]: do units */
3832 7 : bnf = checkbnf(bnf); nf = bnf_get_nf(bnf);
3833 7 : if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
3834 7 : ID.emb = zm_to_ZM( rowpermute(nfsign_units(bnf,NULL,1), ID.archp) );
3835 7 : join_z = &join_archunit;
3836 : }
3837 : else
3838 : {
3839 7 : join_z = &join_arch;
3840 7 : nf = checknf(bnf);
3841 : }
3842 14 : ID.nf = nf;
3843 14 : av = avma; V = cgetg(l, t_VEC);
3844 63 : for (i = 1; i < l; i++)
3845 : {
3846 49 : z = gel(L,i); lz = lg(z);
3847 49 : gel(V,i) = v = cgetg(lz,t_VEC);
3848 91 : for (j=1; j<lz; j++) gel(v,j) = join_z(&ID, gel(z,j));
3849 : }
3850 14 : return gerepilecopy(av,V);
3851 : }
|
