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. It is distributed in the hope that it will be useful, but WITHOUT
8 : ANY WARRANTY WHATSOEVER.
9 :
10 : Check the License for details. You should have received a copy of it, along
11 : with the package; see the file 'COPYING'. If not, write to the Free Software
12 : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
13 :
14 : /*******************************************************************/
15 : /* */
16 : /* BASIC NF OPERATIONS */
17 : /* */
18 : /*******************************************************************/
19 : #include "pari.h"
20 : #include "paripriv.h"
21 :
22 : /*******************************************************************/
23 : /* */
24 : /* OPERATIONS OVER NUMBER FIELD ELEMENTS. */
25 : /* represented as column vectors over the integral basis */
26 : /* */
27 : /*******************************************************************/
28 : static GEN
29 9688329 : get_tab(GEN nf, long *N)
30 : {
31 9688329 : GEN tab = (typ(nf) == t_MAT)? nf: gel(nf,9);
32 9688329 : *N = nbrows(tab); return tab;
33 : }
34 :
35 : /* x != 0, y t_INT. Return x * y (not memory clean if x = 1) */
36 : static GEN
37 388928702 : _mulii(GEN x, GEN y) {
38 1006018049 : return is_pm1(x)? (signe(x) < 0)? negi(y): y
39 617089347 : : mulii(x, y);
40 : }
41 :
42 : GEN
43 17052 : tablemul_ei_ej(GEN M, long i, long j)
44 : {
45 : long N;
46 17052 : GEN tab = get_tab(M, &N);
47 17052 : tab += (i-1)*N; return gel(tab,j);
48 : }
49 :
50 : /* Outputs x.ei, where ei is the i-th elt of the algebra basis.
51 : * x an RgV of correct length and arbitrary content (polynomials, scalars...).
52 : * M is the multiplication table ei ej = sum_k M_k^(i,j) ek */
53 : GEN
54 10311 : tablemul_ei(GEN M, GEN x, long i)
55 : {
56 : long j, k, N;
57 : GEN v, tab;
58 :
59 10311 : if (i==1) return gcopy(x);
60 10311 : tab = get_tab(M, &N);
61 10311 : if (typ(x) != t_COL) { v = zerocol(N); gel(v,i) = gcopy(x); return v; }
62 10311 : tab += (i-1)*N; v = cgetg(N+1,t_COL);
63 : /* wi . x = [ sum_j tab[k,j] x[j] ]_k */
64 70217 : for (k=1; k<=N; k++)
65 : {
66 59906 : pari_sp av = avma;
67 59906 : GEN s = gen_0;
68 418152 : for (j=1; j<=N; j++)
69 : {
70 358246 : GEN c = gcoeff(tab,k,j);
71 358246 : if (!gequal0(c)) s = gadd(s, gmul(c, gel(x,j)));
72 : }
73 59906 : gel(v,k) = gerepileupto(av,s);
74 : }
75 10311 : return v;
76 : }
77 : /* as tablemul_ei, assume x a ZV of correct length */
78 : GEN
79 7394891 : zk_ei_mul(GEN nf, GEN x, long i)
80 : {
81 : long j, k, N;
82 : GEN v, tab;
83 :
84 7394891 : if (i==1) return ZC_copy(x);
85 7394877 : tab = get_tab(nf, &N); tab += (i-1)*N;
86 7394877 : v = cgetg(N+1,t_COL);
87 57375223 : for (k=1; k<=N; k++)
88 : {
89 49980346 : pari_sp av = avma;
90 49980346 : GEN s = gen_0;
91 636542964 : for (j=1; j<=N; j++)
92 : {
93 586562618 : GEN c = gcoeff(tab,k,j);
94 586562618 : if (signe(c)) s = addii(s, _mulii(c, gel(x,j)));
95 : }
96 49980346 : gel(v,k) = gerepileuptoint(av, s);
97 : }
98 7394877 : return v;
99 : }
100 :
101 : /* table of multiplication by wi in R[w1,..., wN] */
102 : GEN
103 2338 : ei_multable(GEN TAB, long i)
104 : {
105 : long k,N;
106 2338 : GEN m, tab = get_tab(TAB, &N);
107 2338 : tab += (i-1)*N;
108 2338 : m = cgetg(N+1,t_MAT);
109 2338 : for (k=1; k<=N; k++) gel(m,k) = gel(tab,k);
110 2338 : return m;
111 : }
112 :
113 : GEN
114 2986077 : zk_multable(GEN nf, GEN x)
115 : {
116 2986077 : long i, l = lg(x);
117 2986077 : GEN mul = cgetg(l,t_MAT);
118 2986077 : gel(mul,1) = x; /* assume w_1 = 1 */
119 2986077 : for (i=2; i<l; i++) gel(mul,i) = zk_ei_mul(nf,x,i);
120 2986077 : return mul;
121 : }
122 : GEN
123 1673 : multable(GEN M, GEN x)
124 : {
125 : long i, N;
126 : GEN mul;
127 1673 : if (typ(x) == t_MAT) return x;
128 0 : M = get_tab(M, &N);
129 0 : if (typ(x) != t_COL) return scalarmat(x, N);
130 0 : mul = cgetg(N+1,t_MAT);
131 0 : gel(mul,1) = x; /* assume w_1 = 1 */
132 0 : for (i=2; i<=N; i++) gel(mul,i) = tablemul_ei(M,x,i);
133 0 : return mul;
134 : }
135 :
136 : /* x integral in nf; table of multiplication by x in ZK = Z[w1,..., wN].
137 : * Return a t_INT if x is scalar, and a ZM otherwise */
138 : GEN
139 1593609 : zk_scalar_or_multable(GEN nf, GEN x)
140 : {
141 1593609 : long tx = typ(x);
142 1593609 : if (tx == t_MAT || tx == t_INT) return x;
143 1568793 : x = nf_to_scalar_or_basis(nf, x);
144 1568793 : return (typ(x) == t_COL)? zk_multable(nf, x): x;
145 : }
146 :
147 : GEN
148 23555 : nftrace(GEN nf, GEN x)
149 : {
150 23555 : pari_sp av = avma;
151 23555 : nf = checknf(nf);
152 23555 : x = nf_to_scalar_or_basis(nf, x);
153 70644 : x = (typ(x) == t_COL)? RgV_dotproduct(x, gel(nf_get_Tr(nf),1))
154 47089 : : gmulgs(x, nf_get_degree(nf));
155 23555 : return gerepileupto(av, x);
156 : }
157 : GEN
158 581 : rnfelttrace(GEN rnf, GEN x)
159 : {
160 581 : pari_sp av = avma;
161 581 : checkrnf(rnf);
162 581 : x = rnfeltabstorel(rnf, x);
163 1414 : x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gtrace(x))
164 924 : : gmulgs(x, rnf_get_degree(rnf));
165 490 : return gerepileupto(av, x);
166 : }
167 :
168 : /* assume nf is a genuine nf, fa a famat */
169 : static GEN
170 7 : famat_norm(GEN nf, GEN fa)
171 : {
172 7 : pari_sp av = avma;
173 7 : GEN g = gel(fa,1), e = gel(fa,2), N = gen_1;
174 7 : long i, l = lg(g);
175 21 : for (i = 1; i < l; i++)
176 14 : N = gmul(N, powgi(nfnorm(nf, gel(g,i)), gel(e,i)));
177 7 : return gerepileupto(av, N);
178 : }
179 : GEN
180 24654 : nfnorm(GEN nf, GEN x)
181 : {
182 24654 : pari_sp av = avma;
183 24654 : nf = checknf(nf);
184 24654 : if (typ(x) == t_MAT) return famat_norm(nf, x);
185 24647 : x = nf_to_scalar_or_alg(nf, x);
186 69944 : x = (typ(x) == t_POL)? RgXQ_norm(x, nf_get_pol(nf))
187 45297 : : gpowgs(x, nf_get_degree(nf));
188 24647 : return gerepileupto(av, x);
189 : }
190 :
191 : GEN
192 231 : rnfeltnorm(GEN rnf, GEN x)
193 : {
194 231 : pari_sp av = avma;
195 231 : checkrnf(rnf);
196 231 : x = rnfeltabstorel(rnf, x);
197 378 : x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gnorm(x))
198 238 : : gpowgs(x, rnf_get_degree(rnf));
199 140 : return gerepileupto(av, x);
200 : }
201 :
202 : /* x + y in nf */
203 : GEN
204 15595280 : nfadd(GEN nf, GEN x, GEN y)
205 : {
206 15595280 : pari_sp av = avma;
207 : GEN z;
208 :
209 15595280 : nf = checknf(nf);
210 15595280 : x = nf_to_scalar_or_basis(nf, x);
211 15595280 : y = nf_to_scalar_or_basis(nf, y);
212 15595280 : if (typ(x) != t_COL)
213 12689119 : { z = (typ(y) == t_COL)? RgC_Rg_add(y, x): gadd(x,y); }
214 : else
215 2906161 : { z = (typ(y) == t_COL)? RgC_add(x, y): RgC_Rg_add(x, y); }
216 15595280 : return gerepileupto(av, z);
217 : }
218 : /* x - y in nf */
219 : GEN
220 1193360 : nfsub(GEN nf, GEN x, GEN y)
221 : {
222 1193360 : pari_sp av = avma;
223 : GEN z;
224 :
225 1193360 : nf = checknf(nf);
226 1193360 : x = nf_to_scalar_or_basis(nf, x);
227 1193360 : y = nf_to_scalar_or_basis(nf, y);
228 1193360 : if (typ(x) != t_COL)
229 884072 : { z = (typ(y) == t_COL)? Rg_RgC_sub(x,y): gsub(x,y); }
230 : else
231 309288 : { z = (typ(y) == t_COL)? RgC_sub(x,y): RgC_Rg_sub(x,y); }
232 1193360 : return gerepileupto(av, z);
233 : }
234 :
235 : /* product of x and y in nf */
236 : GEN
237 20604646 : nfmul(GEN nf, GEN x, GEN y)
238 : {
239 : GEN z;
240 20604646 : pari_sp av = avma;
241 :
242 20604646 : if (x == y) return nfsqr(nf,x);
243 :
244 17745783 : nf = checknf(nf);
245 17745783 : x = nf_to_scalar_or_basis(nf, x);
246 17745783 : y = nf_to_scalar_or_basis(nf, y);
247 17745783 : if (typ(x) != t_COL)
248 : {
249 13740104 : if (isintzero(x)) return gen_0;
250 9855209 : z = (typ(y) == t_COL)? RgC_Rg_mul(y, x): gmul(x,y); }
251 : else
252 : {
253 4005679 : if (typ(y) != t_COL)
254 : {
255 2845745 : if (isintzero(y)) return gen_0;
256 640339 : z = RgC_Rg_mul(x, y);
257 : }
258 : else
259 : {
260 : GEN dx, dy;
261 1159934 : x = Q_remove_denom(x, &dx);
262 1159934 : y = Q_remove_denom(y, &dy);
263 1159934 : z = nfmuli(nf,x,y);
264 1159934 : dx = mul_denom(dx,dy);
265 1159934 : if (dx) z = ZC_Z_div(z, dx);
266 : }
267 : }
268 11655482 : return gerepileupto(av, z);
269 : }
270 : /* square of x in nf */
271 : GEN
272 4750678 : nfsqr(GEN nf, GEN x)
273 : {
274 4750678 : pari_sp av = avma;
275 : GEN z;
276 :
277 4750678 : nf = checknf(nf);
278 4750678 : x = nf_to_scalar_or_basis(nf, x);
279 4750678 : if (typ(x) != t_COL) z = gsqr(x);
280 : else
281 : {
282 : GEN dx;
283 113472 : x = Q_remove_denom(x, &dx);
284 113472 : z = nfsqri(nf,x);
285 113472 : if (dx) z = RgC_Rg_div(z, sqri(dx));
286 : }
287 4750678 : return gerepileupto(av, z);
288 : }
289 :
290 : /* x a ZC, v a t_COL of ZC/Z */
291 : GEN
292 129323 : zkC_multable_mul(GEN v, GEN x)
293 : {
294 129323 : long i, l = lg(v);
295 129323 : GEN y = cgetg(l, t_COL);
296 456842 : for (i = 1; i < l; i++)
297 : {
298 327519 : GEN c = gel(v,i);
299 327519 : if (typ(c)!=t_COL) {
300 0 : if (!isintzero(c)) c = ZC_Z_mul(gel(x,1), c);
301 : } else {
302 327519 : c = ZM_ZC_mul(x,c);
303 327519 : if (ZV_isscalar(c)) c = gel(c,1);
304 : }
305 327519 : gel(y,i) = c;
306 : }
307 129323 : return y;
308 : }
309 :
310 : GEN
311 50456 : nfC_multable_mul(GEN v, GEN x)
312 : {
313 50456 : long i, l = lg(v);
314 50456 : GEN y = cgetg(l, t_COL);
315 314622 : for (i = 1; i < l; i++)
316 : {
317 264166 : GEN c = gel(v,i);
318 264166 : if (typ(c)!=t_COL) {
319 217288 : if (!isintzero(c)) c = RgC_Rg_mul(gel(x,1), c);
320 : } else {
321 46878 : c = RgM_RgC_mul(x,c);
322 46878 : if (QV_isscalar(c)) c = gel(c,1);
323 : }
324 264166 : gel(y,i) = c;
325 : }
326 50456 : return y;
327 : }
328 :
329 : GEN
330 166663 : nfC_nf_mul(GEN nf, GEN v, GEN x)
331 : {
332 : long tx;
333 : GEN y;
334 :
335 166663 : x = nf_to_scalar_or_basis(nf, x);
336 166663 : tx = typ(x);
337 166663 : if (tx != t_COL)
338 : {
339 : long l, i;
340 122857 : if (tx == t_INT)
341 : {
342 114905 : long s = signe(x);
343 114905 : if (!s) return zerocol(lg(v)-1);
344 108644 : if (is_pm1(x)) return s > 0? leafcopy(v): RgC_neg(v);
345 : }
346 39344 : l = lg(v); y = cgetg(l, t_COL);
347 275880 : for (i=1; i < l; i++)
348 : {
349 236536 : GEN c = gel(v,i);
350 236536 : if (typ(c) != t_COL) c = gmul(c, x); else c = RgC_Rg_mul(c, x);
351 236536 : gel(y,i) = c;
352 : }
353 39344 : return y;
354 : }
355 : else
356 : {
357 : GEN dx;
358 43806 : x = zk_multable(nf, Q_remove_denom(x,&dx));
359 43806 : y = nfC_multable_mul(v, x);
360 43806 : return dx? RgC_Rg_div(y, dx): y;
361 : }
362 : }
363 : static GEN
364 7602 : mulbytab(GEN M, GEN c)
365 7602 : { return typ(c) == t_COL? RgM_RgC_mul(M,c): RgC_Rg_mul(gel(M,1), c); }
366 : GEN
367 1673 : tablemulvec(GEN M, GEN x, GEN v)
368 : {
369 : long l, i;
370 : GEN y;
371 :
372 1673 : if (typ(x) == t_COL && RgV_isscalar(x))
373 : {
374 0 : x = gel(x,1);
375 0 : return typ(v) == t_POL? RgX_Rg_mul(v,x): RgV_Rg_mul(v,x);
376 : }
377 1673 : x = multable(M, x); /* multiplication table by x */
378 1673 : y = cgetg_copy(v, &l);
379 1673 : if (typ(v) == t_POL)
380 : {
381 1673 : y[1] = v[1];
382 1673 : for (i=2; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
383 1673 : y = normalizepol(y);
384 : }
385 : else
386 : {
387 0 : for (i=1; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
388 : }
389 1673 : return y;
390 : }
391 :
392 : GEN
393 378516 : zkmultable_capZ(GEN mx) { return Q_denom(zkmultable_inv(mx)); }
394 :
395 : GEN
396 418620 : zkmultable_inv(GEN mx)
397 418620 : { return ZM_gauss(mx, col_ei(lg(mx)-1,1)); }
398 :
399 : /* nf a true nf, x a ZC */
400 : GEN
401 40104 : zk_inv(GEN nf, GEN x) { return zkmultable_inv(zk_multable(nf,x)); }
402 :
403 : /* inverse of x in nf */
404 : GEN
405 64043 : nfinv(GEN nf, GEN x)
406 : {
407 64043 : pari_sp av = avma;
408 : GEN z;
409 :
410 64043 : nf = checknf(nf);
411 64043 : x = nf_to_scalar_or_basis(nf, x);
412 64043 : if (typ(x) == t_COL)
413 : {
414 : GEN d;
415 24752 : x = Q_remove_denom(x, &d);
416 24752 : z = zk_inv(nf, x);
417 24752 : if (d) z = RgC_Rg_mul(z, d);
418 : }
419 : else
420 39291 : z = ginv(x);
421 64043 : return gerepileupto(av, z);
422 : }
423 :
424 : /* quotient of x and y in nf */
425 : GEN
426 18410 : nfdiv(GEN nf, GEN x, GEN y)
427 : {
428 18410 : pari_sp av = avma;
429 : GEN z;
430 :
431 18410 : nf = checknf(nf);
432 18410 : y = nf_to_scalar_or_basis(nf, y);
433 18410 : if (typ(y) != t_COL)
434 : {
435 10443 : x = nf_to_scalar_or_basis(nf, x);
436 10443 : z = (typ(x) == t_COL)? RgC_Rg_div(x, y): gdiv(x,y);
437 : }
438 : else
439 : {
440 : GEN d;
441 7967 : y = Q_remove_denom(y, &d);
442 7967 : z = nfmul(nf, x, zk_inv(nf,y));
443 7967 : if (d) z = RgC_Rg_mul(z, d);
444 : }
445 18410 : return gerepileupto(av, z);
446 : }
447 :
448 : /* product of INTEGERS (t_INT or ZC) x and y in nf
449 : * compute xy as ( sum_i x_i sum_j y_j m^{i,j}_k )_k */
450 : GEN
451 1589298 : nfmuli(GEN nf, GEN x, GEN y)
452 : {
453 : long i, j, k, N;
454 1589298 : GEN s, v, TAB = get_tab(nf, &N);
455 :
456 1589298 : if (typ(x) == t_INT) return (typ(y) == t_COL)? ZC_Z_mul(y, x): mulii(x,y);
457 1493513 : if (typ(y) == t_INT) return ZC_Z_mul(x, y);
458 : /* both x and y are ZV */
459 1459148 : v = cgetg(N+1,t_COL);
460 5965697 : for (k=1; k<=N; k++)
461 : {
462 4506549 : pari_sp av = avma;
463 4506549 : GEN TABi = TAB;
464 4506549 : if (k == 1)
465 1459148 : s = mulii(gel(x,1),gel(y,1));
466 : else
467 6094802 : s = addii(mulii(gel(x,1),gel(y,k)),
468 6094802 : mulii(gel(x,k),gel(y,1)));
469 23145995 : for (i=2; i<=N; i++)
470 : {
471 18639446 : GEN t, xi = gel(x,i);
472 18639446 : TABi += N;
473 18639446 : if (!signe(xi)) continue;
474 :
475 12922054 : t = NULL;
476 133095868 : for (j=2; j<=N; j++)
477 : {
478 120173814 : GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
479 120173814 : if (!signe(c)) continue;
480 54321085 : p1 = _mulii(c, gel(y,j));
481 54321085 : t = t? addii(t, p1): p1;
482 : }
483 12922054 : if (t) s = addii(s, mulii(xi, t));
484 : }
485 4506549 : gel(v,k) = gerepileuptoint(av,s);
486 : }
487 1459148 : return v;
488 : }
489 : /* square of INTEGER (t_INT or ZC) x in nf */
490 : GEN
491 674453 : nfsqri(GEN nf, GEN x)
492 : {
493 : long i, j, k, N;
494 674453 : GEN s, v, TAB = get_tab(nf, &N);
495 :
496 674453 : if (typ(x) == t_INT) return sqri(x);
497 674453 : v = cgetg(N+1,t_COL);
498 5312279 : for (k=1; k<=N; k++)
499 : {
500 4637826 : pari_sp av = avma;
501 4637826 : GEN TABi = TAB;
502 4637826 : if (k == 1)
503 674453 : s = sqri(gel(x,1));
504 : else
505 3963373 : s = shifti(mulii(gel(x,1),gel(x,k)), 1);
506 55847444 : for (i=2; i<=N; i++)
507 : {
508 51209618 : GEN p1, c, t, xi = gel(x,i);
509 51209618 : TABi += N;
510 51209618 : if (!signe(xi)) continue;
511 :
512 17360568 : c = gcoeff(TABi, k, i);
513 17360568 : t = signe(c)? _mulii(c,xi): NULL;
514 249464217 : for (j=i+1; j<=N; j++)
515 : {
516 232103649 : c = gcoeff(TABi, k, j);
517 232103649 : if (!signe(c)) continue;
518 121429408 : p1 = _mulii(c, shifti(gel(x,j),1));
519 121429408 : t = t? addii(t, p1): p1;
520 : }
521 17360568 : if (t) s = addii(s, mulii(xi, t));
522 : }
523 4637826 : gel(v,k) = gerepileuptoint(av,s);
524 : }
525 674453 : return v;
526 : }
527 :
528 : /* both x and y are RgV */
529 : GEN
530 0 : tablemul(GEN TAB, GEN x, GEN y)
531 : {
532 : long i, j, k, N;
533 : GEN s, v;
534 0 : if (typ(x) != t_COL) return gmul(x, y);
535 0 : if (typ(y) != t_COL) return gmul(y, x);
536 0 : N = lg(x)-1;
537 0 : v = cgetg(N+1,t_COL);
538 0 : for (k=1; k<=N; k++)
539 : {
540 0 : pari_sp av = avma;
541 0 : GEN TABi = TAB;
542 0 : if (k == 1)
543 0 : s = gmul(gel(x,1),gel(y,1));
544 : else
545 0 : s = gadd(gmul(gel(x,1),gel(y,k)),
546 0 : gmul(gel(x,k),gel(y,1)));
547 0 : for (i=2; i<=N; i++)
548 : {
549 0 : GEN t, xi = gel(x,i);
550 0 : TABi += N;
551 0 : if (gequal0(xi)) continue;
552 :
553 0 : t = NULL;
554 0 : for (j=2; j<=N; j++)
555 : {
556 0 : GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
557 0 : if (gequal0(c)) continue;
558 0 : p1 = gmul(c, gel(y,j));
559 0 : t = t? gadd(t, p1): p1;
560 : }
561 0 : if (t) s = gadd(s, gmul(xi, t));
562 : }
563 0 : gel(v,k) = gerepileupto(av,s);
564 : }
565 0 : return v;
566 : }
567 : GEN
568 40040 : tablesqr(GEN TAB, GEN x)
569 : {
570 : long i, j, k, N;
571 : GEN s, v;
572 :
573 40040 : if (typ(x) != t_COL) return gsqr(x);
574 40040 : N = lg(x)-1;
575 40040 : v = cgetg(N+1,t_COL);
576 :
577 278922 : for (k=1; k<=N; k++)
578 : {
579 238882 : pari_sp av = avma;
580 238882 : GEN TABi = TAB;
581 238882 : if (k == 1)
582 40040 : s = gsqr(gel(x,1));
583 : else
584 198842 : s = gmul2n(gmul(gel(x,1),gel(x,k)), 1);
585 1455440 : for (i=2; i<=N; i++)
586 : {
587 1216558 : GEN p1, c, t, xi = gel(x,i);
588 1216558 : TABi += N;
589 1216558 : if (gequal0(xi)) continue;
590 :
591 321069 : c = gcoeff(TABi, k, i);
592 321069 : t = !gequal0(c)? gmul(c,xi): NULL;
593 1241723 : for (j=i+1; j<=N; j++)
594 : {
595 920654 : c = gcoeff(TABi, k, j);
596 920654 : if (gequal0(c)) continue;
597 482706 : p1 = gmul(gmul2n(c,1), gel(x,j));
598 482706 : t = t? gadd(t, p1): p1;
599 : }
600 321069 : if (t) s = gadd(s, gmul(xi, t));
601 : }
602 238882 : gel(v,k) = gerepileupto(av,s);
603 : }
604 40040 : return v;
605 : }
606 :
607 : static GEN
608 28554 : _mul(void *data, GEN x, GEN y) { return nfmuli((GEN)data,x,y); }
609 : static GEN
610 105678 : _sqr(void *data, GEN x) { return nfsqri((GEN)data,x); }
611 :
612 : /* Compute z^n in nf, left-shift binary powering */
613 : GEN
614 102063 : nfpow(GEN nf, GEN z, GEN n)
615 : {
616 102063 : pari_sp av = avma;
617 : long s;
618 : GEN x, cx;
619 :
620 102063 : if (typ(n)!=t_INT) pari_err_TYPE("nfpow",n);
621 102063 : nf = checknf(nf);
622 102063 : s = signe(n); if (!s) return gen_1;
623 102063 : x = nf_to_scalar_or_basis(nf, z);
624 102063 : if (typ(x) != t_COL) return powgi(x,n);
625 101587 : if (s < 0)
626 : { /* simplified nfinv */
627 : GEN d;
628 3451 : x = Q_remove_denom(x, &d);
629 3451 : x = zk_inv(nf, x);
630 3451 : x = primitive_part(x, &cx);
631 3451 : cx = mul_content(cx, d);
632 3451 : n = absi(n);
633 : }
634 : else
635 98136 : x = primitive_part(x, &cx);
636 101587 : x = gen_pow(x, n, (void*)nf, _sqr, _mul);
637 101587 : if (cx) x = gmul(x, powgi(cx, n));
638 101587 : return av==avma? gcopy(x): gerepileupto(av,x);
639 : }
640 : /* Compute z^n in nf, left-shift binary powering */
641 : GEN
642 45997 : nfpow_u(GEN nf, GEN z, ulong n)
643 : {
644 45997 : pari_sp av = avma;
645 : GEN x, cx;
646 :
647 45997 : nf = checknf(nf);
648 45997 : if (!n) return gen_1;
649 45997 : x = nf_to_scalar_or_basis(nf, z);
650 45997 : if (typ(x) != t_COL) return gpowgs(x,n);
651 17759 : x = primitive_part(x, &cx);
652 17759 : x = gen_powu(x, n, (void*)nf, _sqr, _mul);
653 17759 : if (cx) x = gmul(x, powgi(cx, utoipos(n)));
654 17759 : return av==avma? gcopy(x): gerepileupto(av,x);
655 : }
656 :
657 : static GEN
658 2402169 : _nf_red(void *E, GEN x) { (void)E; return x; }
659 :
660 : static GEN
661 9332820 : _nf_add(void *E, GEN x, GEN y) { return nfadd((GEN)E,x,y); }
662 :
663 : static GEN
664 571655 : _nf_neg(void *E, GEN x) { (void)E; return gneg(x); }
665 :
666 : static GEN
667 11300121 : _nf_mul(void *E, GEN x, GEN y) { return nfmul((GEN)E,x,y); }
668 :
669 : static GEN
670 41734 : _nf_inv(void *E, GEN x) { return nfinv((GEN)E,x); }
671 :
672 : static GEN
673 8316 : _nf_s(void *E, long x) { (void)E; return stoi(x); }
674 :
675 : static const struct bb_field nf_field={_nf_red,_nf_add,_nf_mul,_nf_neg,
676 : _nf_inv,&gequal0,_nf_s };
677 :
678 178066 : const struct bb_field *get_nf_field(void **E, GEN nf)
679 178066 : { *E = (void*)nf; return &nf_field; }
680 :
681 : GEN
682 14 : nfM_det(GEN nf, GEN M)
683 : {
684 : void *E;
685 14 : const struct bb_field *S = get_nf_field(&E, nf);
686 14 : return gen_det(M, E, S);
687 : }
688 : GEN
689 8302 : nfM_inv(GEN nf, GEN M)
690 : {
691 : void *E;
692 8302 : const struct bb_field *S = get_nf_field(&E, nf);
693 8302 : return gen_Gauss(M, matid(lg(M)-1), E, S);
694 : }
695 : GEN
696 8092 : nfM_mul(GEN nf, GEN A, GEN B)
697 : {
698 : void *E;
699 8092 : const struct bb_field *S = get_nf_field(&E, nf);
700 8092 : return gen_matmul(A, B, E, S);
701 : }
702 : GEN
703 161658 : nfM_nfC_mul(GEN nf, GEN A, GEN B)
704 : {
705 : void *E;
706 161658 : const struct bb_field *S = get_nf_field(&E, nf);
707 161658 : return gen_matcolmul(A, B, E, S);
708 : }
709 :
710 : /* valuation of integral x (ZV), with resp. to prime ideal pr */
711 : long
712 5303783 : ZC_nfvalrem(GEN x, GEN pr, GEN *newx)
713 : {
714 : long i, v, l;
715 5303783 : GEN r, y, p = pr_get_p(pr), mul = pr_get_tau(pr);
716 :
717 : /* p inert */
718 5303783 : if (typ(mul) == t_INT) return newx? ZV_pvalrem(x, p, newx):ZV_pval(x, p);
719 5292793 : y = cgetg_copy(x, &l); /* will hold the new x */
720 5292793 : x = leafcopy(x);
721 7952195 : for(v=0;; v++)
722 : {
723 27929246 : for (i=1; i<l; i++)
724 : { /* is (x.b)[i] divisible by p ? */
725 25269844 : gel(y,i) = dvmdii(ZMrow_ZC_mul(mul,x,i),p,&r);
726 25269844 : if (r != gen_0) { if (newx) *newx = x; return v; }
727 : }
728 2659402 : swap(x, y);
729 2659402 : }
730 : }
731 : long
732 5067081 : ZC_nfval(GEN x, GEN P)
733 5067081 : { return ZC_nfvalrem(x, P, NULL); }
734 :
735 : /* v_P(x) != 0, x a ZV. Simpler version of ZC_nfvalrem */
736 : int
737 243866 : ZC_prdvd(GEN x, GEN P)
738 : {
739 243866 : pari_sp av = avma;
740 : long i, l;
741 243866 : GEN p = pr_get_p(P), mul = pr_get_tau(P);
742 243866 : if (typ(mul) == t_INT) return ZV_Z_dvd(x, p);
743 243621 : l = lg(x);
744 975373 : for (i=1; i<l; i++)
745 871946 : if (!dvdii(ZMrow_ZC_mul(mul,x,i), p)) { avma = av; return 0; }
746 103427 : avma = av; return 1;
747 : }
748 :
749 : int
750 28 : pr_equal(GEN P, GEN Q)
751 : {
752 28 : GEN gQ, p = pr_get_p(P);
753 28 : long e = pr_get_e(P), f = pr_get_f(P), n;
754 28 : if (!equalii(p, pr_get_p(Q)) || e != pr_get_e(Q) || f != pr_get_f(Q))
755 14 : return 0;
756 14 : gQ = pr_get_gen(Q); n = lg(gQ)-1;
757 14 : if (2*e*f > n) return 1; /* room for only one such pr */
758 7 : return ZV_equal(pr_get_gen(P), gQ) || ZC_prdvd(gQ, P);
759 : }
760 :
761 : long
762 1302189 : nfval(GEN nf, GEN x, GEN pr)
763 : {
764 1302189 : pari_sp av = avma;
765 : long w, e;
766 : GEN cx, p;
767 :
768 1302189 : if (gequal0(x)) return LONG_MAX;
769 1300698 : nf = checknf(nf);
770 1300698 : checkprid(pr);
771 1300698 : p = pr_get_p(pr);
772 1300698 : e = pr_get_e(pr);
773 1300698 : x = nf_to_scalar_or_basis(nf, x);
774 1300698 : if (typ(x) != t_COL) return e*Q_pval(x,p);
775 106813 : x = Q_primitive_part(x, &cx);
776 106813 : w = ZC_nfval(x,pr);
777 106813 : if (cx) w += e*Q_pval(cx,p);
778 106813 : avma = av; return w;
779 : }
780 :
781 : /* want to write p^v = uniformizer^(e*v) * z^v, z coprime to pr */
782 : /* z := tau^e / p^(e-1), algebraic integer coprime to pr; return z^v */
783 : static GEN
784 20104 : powp(GEN nf, GEN pr, long v)
785 : {
786 : GEN b, z;
787 : long e;
788 20104 : if (!v) return gen_1;
789 19978 : b = pr_get_tau(pr);
790 19978 : if (typ(b) == t_INT) return gen_1;
791 1274 : e = pr_get_e(pr);
792 1274 : z = gel(b,1);
793 1274 : if (e != 1) z = gdiv(nfpow_u(nf, z, e), powiu(pr_get_p(pr),e-1));
794 1274 : if (v < 0) { v = -v; z = nfinv(nf, z); }
795 1274 : if (v != 1) z = nfpow_u(nf, z, v);
796 1274 : return z;
797 : }
798 : long
799 64925 : nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
800 : {
801 64925 : pari_sp av = avma;
802 : long w, e;
803 : GEN cx, p, t;
804 :
805 64925 : if (!py) return nfval(nf,x,pr);
806 64806 : if (gequal0(x)) { *py = gcopy(x); return LONG_MAX; }
807 64750 : nf = checknf(nf);
808 64750 : checkprid(pr);
809 64750 : p = pr_get_p(pr);
810 64750 : e = pr_get_e(pr);
811 64750 : x = nf_to_scalar_or_basis(nf, x);
812 64750 : if (typ(x) != t_COL) {
813 52871 : w = Q_pvalrem(x,p, py);
814 52871 : if (!w) { *py = gerepilecopy(av, x); return 0; }
815 18977 : *py = gerepileupto(av, gmul(powp(nf, pr, w), *py));
816 18977 : return e*w;
817 : }
818 11879 : x = Q_primitive_part(x, &cx);
819 11879 : w = ZC_nfvalrem(x,pr, py);
820 11879 : if (cx)
821 : {
822 1127 : long v = Q_pvalrem(cx,p, &t);
823 1127 : *py = nfmul(nf, *py, gmul(powp(nf,pr,v), t));
824 1127 : *py = gerepileupto(av, *py);
825 1127 : w += e*v;
826 : }
827 : else
828 10752 : *py = gerepilecopy(av, *py);
829 11879 : return w;
830 : }
831 : GEN
832 147 : gpnfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
833 : {
834 147 : long v = nfvalrem(nf,x,pr,py);
835 147 : return v == LONG_MAX? mkoo(): stoi(v);
836 : }
837 :
838 : /* true nf */
839 : GEN
840 78659 : coltoalg(GEN nf, GEN x)
841 : {
842 78659 : return mkpolmod( nf_to_scalar_or_alg(nf, x), nf_get_pol(nf) );
843 : }
844 :
845 : GEN
846 87815 : basistoalg(GEN nf, GEN x)
847 : {
848 : GEN z, T;
849 :
850 87815 : nf = checknf(nf);
851 87815 : switch(typ(x))
852 : {
853 : case t_COL: {
854 72541 : pari_sp av = avma;
855 72541 : return gerepilecopy(av, coltoalg(nf, x));
856 : }
857 : case t_POLMOD:
858 483 : T = nf_get_pol(nf);
859 483 : if (!RgX_equal_var(T,gel(x,1)))
860 0 : pari_err_MODULUS("basistoalg", T,gel(x,1));
861 483 : return gcopy(x);
862 : case t_POL:
863 1015 : T = nf_get_pol(nf);
864 1015 : if (varn(T) != varn(x)) pari_err_VAR("basistoalg",x,T);
865 1015 : z = cgetg(3,t_POLMOD);
866 1015 : gel(z,1) = ZX_copy(T);
867 1015 : gel(z,2) = RgX_rem(x, T); return z;
868 : case t_INT:
869 : case t_FRAC:
870 13776 : T = nf_get_pol(nf);
871 13776 : z = cgetg(3,t_POLMOD);
872 13776 : gel(z,1) = ZX_copy(T);
873 13776 : gel(z,2) = gcopy(x); return z;
874 : default:
875 0 : pari_err_TYPE("basistoalg",x);
876 : return NULL; /* LCOV_EXCL_LINE */
877 : }
878 : }
879 :
880 : /* true nf, x a t_POL */
881 : static GEN
882 1448048 : pol_to_scalar_or_basis(GEN nf, GEN x)
883 : {
884 1448048 : GEN T = nf_get_pol(nf);
885 1448048 : long l = lg(x);
886 1448048 : if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_basis", x,T);
887 1447985 : if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
888 1447985 : if (l == 2) return gen_0;
889 852432 : if (l == 3)
890 : {
891 201026 : x = gel(x,2);
892 201026 : if (!is_rational_t(typ(x))) pari_err_TYPE("nf_to_scalar_or_basis",x);
893 201026 : return x;
894 : }
895 651406 : return poltobasis(nf,x);
896 : }
897 : /* Assume nf is a genuine nf. */
898 : GEN
899 81630391 : nf_to_scalar_or_basis(GEN nf, GEN x)
900 : {
901 81630391 : switch(typ(x))
902 : {
903 : case t_INT: case t_FRAC:
904 62906544 : return x;
905 : case t_POLMOD:
906 185297 : x = checknfelt_mod(nf,x,"nf_to_scalar_or_basis");
907 185227 : switch(typ(x))
908 : {
909 34132 : case t_INT: case t_FRAC: return x;
910 151095 : case t_POL: return pol_to_scalar_or_basis(nf,x);
911 : }
912 0 : break;
913 1296953 : case t_POL: return pol_to_scalar_or_basis(nf,x);
914 : case t_COL:
915 17241597 : if (lg(x)-1 != nf_get_degree(nf)) break;
916 17241534 : return QV_isscalar(x)? gel(x,1): x;
917 : }
918 63 : pari_err_TYPE("nf_to_scalar_or_basis",x);
919 : return NULL; /* LCOV_EXCL_LINE */
920 : }
921 : /* Let x be a polynomial with coefficients in Q or nf. Return the same
922 : * polynomial with coefficients expressed as vectors (on the integral basis).
923 : * No consistency checks, not memory-clean. */
924 : GEN
925 5853 : RgX_to_nfX(GEN nf, GEN x)
926 : {
927 : long i, l;
928 5853 : GEN y = cgetg_copy(x, &l); y[1] = x[1];
929 5853 : for (i=2; i<l; i++) gel(y,i) = nf_to_scalar_or_basis(nf, gel(x,i));
930 5853 : return y;
931 : }
932 :
933 : /* Assume nf is a genuine nf. */
934 : GEN
935 182746 : nf_to_scalar_or_alg(GEN nf, GEN x)
936 : {
937 182746 : switch(typ(x))
938 : {
939 : case t_INT: case t_FRAC:
940 14853 : return x;
941 : case t_POLMOD:
942 1365 : x = checknfelt_mod(nf,x,"nf_to_scalar_or_alg");
943 1365 : if (typ(x) != t_POL) return x;
944 : /* fall through */
945 : case t_POL:
946 : {
947 14883 : GEN T = nf_get_pol(nf);
948 14883 : long l = lg(x);
949 14883 : if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_alg", x,T);
950 14883 : if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
951 14883 : if (l == 2) return gen_0;
952 14883 : if (l == 3) return gel(x,2);
953 14659 : return x;
954 : }
955 : case t_COL:
956 : {
957 : GEN dx;
958 152954 : if (lg(x)-1 != nf_get_degree(nf)) break;
959 305908 : if (QV_isscalar(x)) return gel(x,1);
960 124548 : x = Q_remove_denom(x, &dx);
961 124548 : x = RgV_RgC_mul(nf_get_zkprimpart(nf), x);
962 124548 : dx = mul_denom(dx, nf_get_zkden(nf));
963 124548 : return gdiv(x,dx);
964 : }
965 : }
966 49 : pari_err_TYPE("nf_to_scalar_or_alg",x);
967 : return NULL; /* LCOV_EXCL_LINE */
968 : }
969 :
970 : /* gmul(A, RgX_to_RgC(x)), A t_MAT of compatible dimensions */
971 : GEN
972 1330 : RgM_RgX_mul(GEN A, GEN x)
973 : {
974 1330 : long i, l = lg(x)-1;
975 : GEN z;
976 1330 : if (l == 1) return zerocol(nbrows(A));
977 1330 : z = gmul(gel(x,2), gel(A,1));
978 2527 : for (i = 2; i < l; i++)
979 1197 : if (!gequal0(gel(x,i+1))) z = gadd(z, gmul(gel(x,i+1), gel(A,i)));
980 1330 : return z;
981 : }
982 : GEN
983 2329657 : ZM_ZX_mul(GEN A, GEN x)
984 : {
985 2329657 : long i, l = lg(x)-1;
986 : GEN z;
987 2329657 : if (l == 1) return zerocol(nbrows(A));
988 2328876 : z = ZC_Z_mul(gel(A,1), gel(x,2));
989 9234016 : for (i = 2; i < l ; i++)
990 6905139 : if (signe(gel(x,i+1))) z = ZC_add(z, ZC_Z_mul(gel(A,i), gel(x,i+1)));
991 2328877 : return z;
992 : }
993 : /* x a t_POL, nf a genuine nf. No garbage collecting. No check. */
994 : GEN
995 2166713 : poltobasis(GEN nf, GEN x)
996 : {
997 2166713 : GEN d, T = nf_get_pol(nf);
998 2166713 : if (varn(x) != varn(T)) pari_err_VAR( "poltobasis", x,T);
999 2166657 : if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
1000 2166657 : x = Q_remove_denom(x, &d);
1001 2166657 : if (!RgX_is_ZX(x)) pari_err_TYPE("poltobasis",x);
1002 2166636 : x = ZM_ZX_mul(nf_get_invzk(nf), x);
1003 2166636 : if (d) x = RgC_Rg_div(x, d);
1004 2166636 : return x;
1005 : }
1006 :
1007 : GEN
1008 154647 : algtobasis(GEN nf, GEN x)
1009 : {
1010 : pari_sp av;
1011 :
1012 154647 : nf = checknf(nf);
1013 154647 : switch(typ(x))
1014 : {
1015 : case t_POLMOD:
1016 38738 : if (!RgX_equal_var(nf_get_pol(nf),gel(x,1)))
1017 7 : pari_err_MODULUS("algtobasis", nf_get_pol(nf),gel(x,1));
1018 38731 : x = gel(x,2);
1019 38731 : switch(typ(x))
1020 : {
1021 : case t_INT:
1022 4480 : case t_FRAC: return scalarcol(x, nf_get_degree(nf));
1023 : case t_POL:
1024 34251 : av = avma;
1025 34251 : return gerepileupto(av,poltobasis(nf,x));
1026 : }
1027 0 : break;
1028 :
1029 : case t_POL:
1030 56097 : av = avma;
1031 56097 : return gerepileupto(av,poltobasis(nf,x));
1032 :
1033 : case t_COL:
1034 14753 : if (lg(x)-1 != nf_get_degree(nf)) pari_err_DIM("nfalgtobasis");
1035 14753 : return gcopy(x);
1036 :
1037 : case t_INT:
1038 45059 : case t_FRAC: return scalarcol(x, nf_get_degree(nf));
1039 : }
1040 0 : pari_err_TYPE("algtobasis",x);
1041 : return NULL; /* LCOV_EXCL_LINE */
1042 : }
1043 :
1044 : GEN
1045 36715 : rnfbasistoalg(GEN rnf,GEN x)
1046 : {
1047 36715 : const char *f = "rnfbasistoalg";
1048 : long lx, i;
1049 36715 : pari_sp av = avma;
1050 : GEN z, nf, relpol, T;
1051 :
1052 36715 : checkrnf(rnf);
1053 36715 : nf = rnf_get_nf(rnf);
1054 36715 : T = nf_get_pol(nf);
1055 36715 : relpol = QXQX_to_mod_shallow(rnf_get_pol(rnf), T);
1056 36715 : switch(typ(x))
1057 : {
1058 : case t_COL:
1059 798 : z = cgetg_copy(x, &lx);
1060 2338 : for (i=1; i<lx; i++)
1061 : {
1062 1589 : GEN c = nf_to_scalar_or_alg(nf, gel(x,i));
1063 1540 : if (typ(c) == t_POL) c = mkpolmod(c,T);
1064 1540 : gel(z,i) = c;
1065 : }
1066 749 : z = RgV_RgC_mul(gel(rnf_get_zk(rnf),1), z);
1067 686 : return gerepileupto(av, gmodulo(z,relpol));
1068 :
1069 : case t_POLMOD:
1070 24675 : x = polmod_nffix(f, rnf, x, 0);
1071 24465 : if (typ(x) != t_POL) break;
1072 9709 : retmkpolmod(RgX_copy(x), RgX_copy(relpol));
1073 : case t_POL:
1074 868 : if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = gmodulo(x,T); break; }
1075 644 : if (varn(x) == varn(relpol))
1076 : {
1077 595 : x = RgX_nffix(f,nf_get_pol(nf),x,0);
1078 595 : return gmodulo(x, relpol);
1079 : }
1080 49 : pari_err_VAR(f, x,relpol);
1081 : }
1082 25305 : retmkpolmod(scalarpol(x, varn(relpol)), RgX_copy(relpol));
1083 : }
1084 :
1085 : GEN
1086 1302 : matbasistoalg(GEN nf,GEN x)
1087 : {
1088 : long i, j, li, lx;
1089 1302 : GEN z = cgetg_copy(x, &lx);
1090 :
1091 1302 : if (lx == 1) return z;
1092 1295 : switch(typ(x))
1093 : {
1094 : case t_VEC: case t_COL:
1095 42 : for (i=1; i<lx; i++) gel(z,i) = basistoalg(nf, gel(x,i));
1096 42 : return z;
1097 1253 : case t_MAT: break;
1098 0 : default: pari_err_TYPE("matbasistoalg",x);
1099 : }
1100 1253 : li = lgcols(x);
1101 4837 : for (j=1; j<lx; j++)
1102 : {
1103 3584 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1104 3584 : gel(z,j) = c;
1105 3584 : for (i=1; i<li; i++) gel(c,i) = basistoalg(nf, gel(xj,i));
1106 : }
1107 1253 : return z;
1108 : }
1109 :
1110 : GEN
1111 2744 : matalgtobasis(GEN nf,GEN x)
1112 : {
1113 : long i, j, li, lx;
1114 2744 : GEN z = cgetg_copy(x, &lx);
1115 :
1116 2744 : if (lx == 1) return z;
1117 2688 : switch(typ(x))
1118 : {
1119 : case t_VEC: case t_COL:
1120 2681 : for (i=1; i<lx; i++) gel(z,i) = algtobasis(nf, gel(x,i));
1121 2681 : return z;
1122 7 : case t_MAT: break;
1123 0 : default: pari_err_TYPE("matalgtobasis",x);
1124 : }
1125 7 : li = lgcols(x);
1126 14 : for (j=1; j<lx; j++)
1127 : {
1128 7 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1129 7 : gel(z,j) = c;
1130 7 : for (i=1; i<li; i++) gel(c,i) = algtobasis(nf, gel(xj,i));
1131 : }
1132 7 : return z;
1133 : }
1134 : GEN
1135 8463 : RgM_to_nfM(GEN nf,GEN x)
1136 : {
1137 : long i, j, li, lx;
1138 8463 : GEN z = cgetg_copy(x, &lx);
1139 :
1140 8463 : if (lx == 1) return z;
1141 8463 : li = lgcols(x);
1142 65219 : for (j=1; j<lx; j++)
1143 : {
1144 56756 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1145 56756 : gel(z,j) = c;
1146 56756 : for (i=1; i<li; i++) gel(c,i) = nf_to_scalar_or_basis(nf, gel(xj,i));
1147 : }
1148 8463 : return z;
1149 : }
1150 : GEN
1151 77917 : RgC_to_nfC(GEN nf, GEN x)
1152 77917 : { pari_APPLY_type(t_COL, nf_to_scalar_or_basis(nf, gel(x,i))) }
1153 :
1154 : /* x a t_POLMOD, supposedly in rnf = K[z]/(T), K = Q[y]/(Tnf) */
1155 : GEN
1156 62041 : polmod_nffix(const char *f, GEN rnf, GEN x, int lift)
1157 62041 : { return polmod_nffix2(f, rnf_get_nfpol(rnf), rnf_get_pol(rnf), x,lift); }
1158 : GEN
1159 62132 : polmod_nffix2(const char *f, GEN T, GEN relpol, GEN x, int lift)
1160 : {
1161 62132 : if (RgX_equal_var(gel(x,1),relpol))
1162 : {
1163 57064 : x = gel(x,2);
1164 57064 : if (typ(x) == t_POL && varn(x) == varn(relpol))
1165 : {
1166 29568 : x = RgX_nffix(f, T, x, lift);
1167 29568 : switch(lg(x))
1168 : {
1169 287 : case 2: return gen_0;
1170 3752 : case 3: return gel(x,2);
1171 : }
1172 25529 : return x;
1173 : }
1174 : }
1175 32564 : return Rg_nffix(f, T, x, lift);
1176 : }
1177 : GEN
1178 1176 : rnfalgtobasis(GEN rnf,GEN x)
1179 : {
1180 1176 : const char *f = "rnfalgtobasis";
1181 1176 : pari_sp av = avma;
1182 : GEN T, relpol;
1183 :
1184 1176 : checkrnf(rnf);
1185 1176 : relpol = rnf_get_pol(rnf);
1186 1176 : T = rnf_get_nfpol(rnf);
1187 1176 : switch(typ(x))
1188 : {
1189 : case t_COL:
1190 49 : if (lg(x)-1 != rnf_get_degree(rnf)) pari_err_DIM(f);
1191 28 : x = RgV_nffix(f, T, x, 0);
1192 21 : return gerepilecopy(av, x);
1193 :
1194 : case t_POLMOD:
1195 1043 : x = polmod_nffix(f, rnf, x, 0);
1196 1001 : if (typ(x) != t_POL) break;
1197 707 : return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
1198 : case t_POL:
1199 56 : if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = mkpolmod(x,T); break; }
1200 35 : x = RgX_nffix(f, T, x, 0);
1201 28 : if (degpol(x) >= degpol(relpol)) x = RgX_rem(x,relpol);
1202 28 : return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
1203 : }
1204 336 : return gerepileupto(av, scalarcol(x, rnf_get_degree(rnf)));
1205 : }
1206 :
1207 : /* Given a and b in nf, gives an algebraic integer y in nf such that a-b.y
1208 : * is "small" */
1209 : GEN
1210 259 : nfdiveuc(GEN nf, GEN a, GEN b)
1211 : {
1212 259 : pari_sp av = avma;
1213 259 : a = nfdiv(nf,a,b);
1214 259 : return gerepileupto(av, ground(a));
1215 : }
1216 :
1217 : /* Given a and b in nf, gives a "small" algebraic integer r in nf
1218 : * of the form a-b.y */
1219 : GEN
1220 259 : nfmod(GEN nf, GEN a, GEN b)
1221 : {
1222 259 : pari_sp av = avma;
1223 259 : GEN p1 = gneg_i(nfmul(nf,b,ground(nfdiv(nf,a,b))));
1224 259 : return gerepileupto(av, nfadd(nf,a,p1));
1225 : }
1226 :
1227 : /* Given a and b in nf, gives a two-component vector [y,r] in nf such
1228 : * that r=a-b.y is "small". */
1229 : GEN
1230 259 : nfdivrem(GEN nf, GEN a, GEN b)
1231 : {
1232 259 : pari_sp av = avma;
1233 259 : GEN p1,z, y = ground(nfdiv(nf,a,b));
1234 :
1235 259 : p1 = gneg_i(nfmul(nf,b,y));
1236 259 : z = cgetg(3,t_VEC);
1237 259 : gel(z,1) = gcopy(y);
1238 259 : gel(z,2) = nfadd(nf,a,p1); return gerepileupto(av, z);
1239 : }
1240 :
1241 : /*************************************************************************/
1242 : /** **/
1243 : /** REAL EMBEDDINGS **/
1244 : /** **/
1245 : /*************************************************************************/
1246 : static GEN
1247 49336 : sarch_get_cyc(GEN sarch) { return gel(sarch,1); }
1248 : static GEN
1249 264685 : sarch_get_archp(GEN sarch) { return gel(sarch,2); }
1250 : static GEN
1251 54555 : sarch_get_MI(GEN sarch) { return gel(sarch,3); }
1252 : static GEN
1253 54555 : sarch_get_lambda(GEN sarch) { return gel(sarch,4); }
1254 : static GEN
1255 54555 : sarch_get_F(GEN sarch) { return gel(sarch,5); }
1256 :
1257 : /* true nf, return number of positive roots of char_x */
1258 : static long
1259 1441 : num_positive(GEN nf, GEN x)
1260 : {
1261 1441 : GEN T = nf_get_pol(nf);
1262 1441 : GEN charx = ZXQ_charpoly(nf_to_scalar_or_alg(nf,x), T, 0);
1263 : long np;
1264 1441 : charx = ZX_radical(charx);
1265 1441 : np = ZX_sturmpart(charx, mkvec2(gen_0,mkoo()));
1266 1441 : return np * (degpol(T) / degpol(charx));
1267 : }
1268 :
1269 : /* x a QC: return sigma_k(x) where 1 <= k <= r1+r2; correct but inefficient
1270 : * if x in Q. M = nf_get_M(nf) */
1271 : static GEN
1272 91 : nfembed_i(GEN M, GEN x, long k)
1273 : {
1274 91 : long i, l = lg(M);
1275 91 : GEN z = gel(x,1);
1276 91 : for (i = 2; i < l; i++) z = gadd(z, gmul(gcoeff(M,k,i), gel(x,i)));
1277 91 : return z;
1278 : }
1279 : GEN
1280 1631 : nfembed(GEN nf, GEN x, long k)
1281 : {
1282 1631 : pari_sp av = avma;
1283 1631 : nf = checknf(nf);
1284 1631 : x = nf_to_scalar_or_basis(nf,x);
1285 1631 : if (typ(x) != t_COL) return gerepilecopy(av, x);
1286 0 : return gerepileupto(av, nfembed_i(nf_get_M(nf),x,k));
1287 : }
1288 :
1289 : /* x a ZC */
1290 : static GEN
1291 394961 : zk_embed(GEN M, GEN x, long k)
1292 : {
1293 394961 : long i, l = lg(x);
1294 394961 : GEN z = gel(x,1); /* times M[k,1], which is 1 */
1295 394961 : for (i = 2; i < l; i++) z = mpadd(z, mpmul(gcoeff(M,k,i), gel(x,i)));
1296 394961 : return z;
1297 : }
1298 :
1299 : /* Given floating point approximation z of sigma_k(x), decide its sign
1300 : * [0/+, 1/- and -1 for FAIL] */
1301 : static long
1302 384384 : eval_sign_embed(GEN z)
1303 : { /* dubious, fail */
1304 384384 : if (typ(z) == t_REAL && realprec(z) <= LOWDEFAULTPREC) return -1;
1305 383437 : return (signe(z) < 1)? 1: 0;
1306 : }
1307 : /* return v such that (-1)^v = sign(sigma_k(x)), x primitive ZC */
1308 : static long
1309 311794 : eval_sign(GEN M, GEN x, long k)
1310 311794 : { return eval_sign_embed( zk_embed(M, x, k) ); }
1311 :
1312 : /* check that signs[i..#signs] == s; signs = NULL encodes "totally positive" */
1313 : static int
1314 0 : oksigns(long l, GEN signs, long i, long s)
1315 : {
1316 0 : if (!signs) return s == 0;
1317 0 : for (; i < l; i++)
1318 0 : if (signs[i] != s) return 0;
1319 0 : return 1;
1320 : }
1321 : /* check that signs[i] = s and signs[i+1..#signs] = 1-s */
1322 : static int
1323 0 : oksigns2(long l, GEN signs, long i, long s)
1324 : {
1325 0 : if (!signs) return s == 0 && i == l-1;
1326 0 : return signs[i] == s && oksigns(l, signs, i+1, 1-s);
1327 : }
1328 :
1329 : /* true nf, x a ZC (primitive for efficiency), embx its embeddings or NULL */
1330 : static int
1331 63441 : nfchecksigns_i(GEN nf, GEN x, GEN embx, GEN signs, GEN archp)
1332 : {
1333 63441 : long l = lg(archp), i;
1334 63441 : GEN M = nf_get_M(nf), sarch = NULL;
1335 63441 : long np = -1;
1336 94717 : for (i = 1; i < l; i++)
1337 : {
1338 : long s;
1339 72779 : if (embx)
1340 72590 : s = eval_sign_embed(gel(embx,i));
1341 : else
1342 189 : s = eval_sign(M, x, archp[i]);
1343 : /* 0 / + or 1 / -; -1 for FAIL */
1344 72779 : if (s < 0) /* failure */
1345 : {
1346 0 : long ni, r1 = nf_get_r1(nf);
1347 : GEN xi;
1348 0 : if (np < 0)
1349 : {
1350 0 : np = num_positive(nf, x);
1351 0 : if (np == 0) return oksigns(l, signs, i, 1);
1352 0 : if (np == r1) return oksigns(l, signs, i, 0);
1353 0 : sarch = nfarchstar(nf, NULL, identity_perm(r1));
1354 : }
1355 0 : xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
1356 0 : xi = Q_primpart(xi);
1357 0 : ni = num_positive(nf, nfmuli(nf,x,xi));
1358 0 : if (ni == 0) return oksigns2(l, signs, i, 0);
1359 0 : if (ni == r1) return oksigns2(l, signs, i, 1);
1360 0 : s = ni < np? 0: 1;
1361 : }
1362 72779 : if (s != (signs? signs[i]: 0)) return 0;
1363 : }
1364 21938 : return 1;
1365 : }
1366 : static void
1367 217 : pl_convert(GEN pl, GEN *psigns, GEN *parchp)
1368 : {
1369 217 : long i, j, l = lg(pl);
1370 217 : GEN signs = cgetg(l, t_VECSMALL);
1371 217 : GEN archp = cgetg(l, t_VECSMALL);
1372 462 : for (i = j = 1; i < l; i++)
1373 : {
1374 245 : if (!pl[i]) continue;
1375 231 : archp[j] = i;
1376 231 : signs[j] = (pl[i] < 0)? 1: 0;
1377 231 : j++;
1378 : }
1379 217 : setlg(archp, j); *parchp = archp;
1380 217 : setlg(signs, j); *psigns = signs;
1381 217 : }
1382 : /* pl : requested signs for real embeddings, 0 = no sign constraint */
1383 : int
1384 581 : nfchecksigns(GEN nf, GEN x, GEN pl)
1385 : {
1386 581 : pari_sp av = avma;
1387 : GEN signs, archp;
1388 : int res;
1389 581 : nf = checknf(nf);
1390 581 : x = nf_to_scalar_or_basis(nf,x);
1391 581 : if (typ(x) != t_COL)
1392 : {
1393 364 : long i, l = lg(pl), s = gsigne(x);
1394 805 : for (i = 1; i < l; i++)
1395 441 : if (pl[i] && pl[i] != s) { avma = av; return 0; }
1396 364 : avma = av; return 1;
1397 : }
1398 217 : pl_convert(pl, &signs, &archp);
1399 217 : res = nfchecksigns_i(nf, x, NULL, signs, archp);
1400 217 : avma = av; return res;
1401 : }
1402 :
1403 : /* signs = NULL: totally positive, else sign[i] = 0 (+) or 1 (-) */
1404 : static GEN
1405 54555 : get_C(GEN lambda, long l, GEN signs)
1406 : {
1407 : long i;
1408 : GEN C, mlambda;
1409 54555 : if (!signs) return const_vec(l-1, lambda);
1410 15754 : C = cgetg(l, t_COL); mlambda = gneg(lambda);
1411 15754 : for (i = 1; i < l; i++) gel(C,i) = signs[i]? mlambda: lambda;
1412 15754 : return C;
1413 : }
1414 : /* signs = NULL: totally positive at archp */
1415 : static GEN
1416 91711 : nfsetsigns(GEN nf, GEN signs, GEN x, GEN sarch)
1417 : {
1418 91711 : long i, l = lg(sarch_get_archp(sarch));
1419 : GEN ex;
1420 : /* Is signature already correct ? */
1421 91711 : if (typ(x) != t_COL)
1422 : {
1423 28487 : long s = gsigne(x) < 0? 1: 0;
1424 28487 : if (!signs)
1425 2667 : i = (s == 1)? 1: l;
1426 : else
1427 : {
1428 40156 : for (i = 1; i < l; i++)
1429 27402 : if (signs[i] != s) break;
1430 : }
1431 28487 : ex = (i < l)? const_col(l-1, x): NULL;
1432 : }
1433 : else
1434 : {
1435 63224 : pari_sp av = avma;
1436 63224 : GEN cex, M = nf_get_M(nf), archp = sarch_get_archp(sarch);
1437 63224 : GEN xp = Q_primitive_part(x,&cex);
1438 63224 : ex = cgetg(l,t_COL);
1439 63224 : for (i = 1; i < l; i++) gel(ex,i) = zk_embed(M,xp,archp[i]);
1440 63224 : if (nfchecksigns_i(nf, xp, ex, signs, archp)) { ex = NULL; avma = av; }
1441 41461 : else if (cex) ex = RgC_Rg_mul(ex, cex); /* put back content */
1442 : }
1443 91711 : if (ex)
1444 : { /* If no, fix it */
1445 54555 : GEN lambda = sarch_get_lambda(sarch);
1446 54555 : GEN MI = sarch_get_MI(sarch);
1447 54555 : GEN F = sarch_get_F(sarch);
1448 54555 : GEN t = RgC_sub(get_C(lambda, l, signs), ex);
1449 : long e;
1450 54555 : t = grndtoi(RgM_RgC_mul(MI,t), &e);
1451 54555 : if (lg(F) != 1) t = ZM_ZC_mul(F, t);
1452 54555 : x = typ(x) == t_COL? RgC_add(t, x): RgC_Rg_add(t, x);
1453 : }
1454 91711 : return x;
1455 : }
1456 : /* - sarch = nfarchstar(nf, F);
1457 : * - x encodes a vector of signs at arch.archp: either a t_VECSMALL
1458 : * (vector of signs as {0,1}-vector), NULL (totally positive at archp),
1459 : * or a non-zero number field element (replaced by its signature at archp);
1460 : * - y is a non-zero number field element
1461 : * Return z = y (mod F) with signs(y, archp) = signs(x) (a {0,1}-vector) */
1462 : GEN
1463 109750 : set_sign_mod_divisor(GEN nf, GEN x, GEN y, GEN sarch)
1464 : {
1465 109750 : GEN archp = sarch_get_archp(sarch);
1466 109750 : if (lg(archp) == 1) return y;
1467 90367 : nf = checknf(nf);
1468 90367 : if (x && typ(x) != t_VECSMALL) x = nfsign_arch(nf, x, archp);
1469 90367 : y = nf_to_scalar_or_basis(nf,y);
1470 90367 : return nfsetsigns(nf, x, y, sarch);
1471 : }
1472 :
1473 : static GEN
1474 14389 : setsigns_init(GEN nf, GEN archp, GEN F, GEN DATA)
1475 : {
1476 14389 : GEN lambda, Mr = rowpermute(nf_get_M(nf), archp), MI = F? RgM_mul(Mr,F): Mr;
1477 14389 : lambda = gmul2n(matrixnorm(MI,DEFAULTPREC), -1);
1478 14389 : if (lg(archp) < lg(MI))
1479 : {
1480 12747 : GEN perm = gel(indexrank(MI), 2);
1481 12747 : if (!F) F = matid(nf_get_degree(nf));
1482 12747 : MI = vecpermute(MI, perm);
1483 12747 : F = vecpermute(F, perm);
1484 : }
1485 14389 : if (!F) F = cgetg(1,t_MAT);
1486 14389 : MI = RgM_inv(MI);
1487 14389 : return mkvec5(DATA, archp, MI, lambda, F);
1488 : }
1489 : /* F non-0 integral ideal in HNF (or NULL: Z_K), compute elements in 1+F
1490 : * whose sign matrix at archp is identity; archp in 'indices' format */
1491 : GEN
1492 28200 : nfarchstar(GEN nf, GEN F, GEN archp)
1493 : {
1494 28200 : long nba = lg(archp) - 1;
1495 28200 : if (!nba) return mkvec2(cgetg(1,t_VEC), archp);
1496 13052 : if (F && equali1(gcoeff(F,1,1))) F = NULL;
1497 13052 : if (F) F = idealpseudored(F, nf_get_roundG(nf));
1498 13052 : return setsigns_init(nf, archp, F, const_vec(nba, gen_2));
1499 : }
1500 :
1501 : /*************************************************************************/
1502 : /** **/
1503 : /** IDEALCHINESE **/
1504 : /** **/
1505 : /*************************************************************************/
1506 : static int
1507 2121 : isprfact(GEN x)
1508 : {
1509 : long i, l;
1510 : GEN L, E;
1511 2121 : if (typ(x) != t_MAT || lg(x) != 3) return 0;
1512 2121 : L = gel(x,1); l = lg(L);
1513 2121 : E = gel(x,2);
1514 5054 : for(i=1; i<l; i++)
1515 : {
1516 2933 : checkprid(gel(L,i));
1517 2933 : if (typ(gel(E,i)) != t_INT) return 0;
1518 : }
1519 2121 : return 1;
1520 : }
1521 :
1522 : /* initialize projectors mod pr[i]^e[i] for idealchinese */
1523 : static GEN
1524 2121 : pr_init(GEN nf, GEN fa, GEN w, GEN dw)
1525 : {
1526 2121 : GEN U, E, F, L = gel(fa,1), E0 = gel(fa,2);
1527 2121 : long i, r = lg(L);
1528 :
1529 2121 : if (w && lg(w) != r) pari_err_TYPE("idealchinese", w);
1530 2121 : if (r == 1 && !dw) return cgetg(1,t_VEC);
1531 2114 : E = leafcopy(E0); /* do not destroy fa[2] */
1532 5047 : for (i = 1; i < r; i++)
1533 2933 : if (signe(gel(E,i)) < 0) gel(E,i) = gen_0;
1534 2114 : F = factorbackprime(nf, L, E);
1535 2114 : if (dw)
1536 : {
1537 686 : F = ZM_Z_mul(F, dw);
1538 1568 : for (i = 1; i < r; i++)
1539 : {
1540 882 : GEN pr = gel(L,i);
1541 882 : long e = itos(gel(E0,i)), v = idealval(nf, dw, pr);
1542 882 : if (e >= 0)
1543 875 : gel(E,i) = addiu(gel(E,i), v);
1544 7 : else if (v + e <= 0)
1545 0 : F = idealmulpowprime(nf, F, pr, stoi(-v)); /* coprime to pr */
1546 : else
1547 : {
1548 7 : F = idealmulpowprime(nf, F, pr, stoi(e));
1549 7 : gel(E,i) = stoi(v + e);
1550 : }
1551 : }
1552 : }
1553 2114 : U = cgetg(r, t_VEC);
1554 5047 : for (i = 1; i < r; i++)
1555 : {
1556 : GEN u;
1557 2933 : if (w && gequal0(gel(w,i))) u = gen_0; /* unused */
1558 : else
1559 : {
1560 2877 : GEN pr = gel(L,i), e = gel(E,i), t;
1561 2877 : t = idealdivpowprime(nf,F, pr, e);
1562 2877 : u = hnfmerge_get_1(t, idealpow(nf, pr, e));
1563 2877 : if (!u) pari_err_COPRIME("idealchinese", t,pr);
1564 : }
1565 2933 : gel(U,i) = u;
1566 : }
1567 2114 : F = idealpseudored(F, nf_get_roundG(nf));
1568 2114 : return mkvec2(F, U);
1569 : }
1570 :
1571 : static GEN
1572 1337 : pl_normalize(GEN nf, GEN pl)
1573 : {
1574 1337 : const char *fun = "idealchinese";
1575 1337 : if (lg(pl)-1 != nf_get_r1(nf)) pari_err_TYPE(fun,pl);
1576 1337 : switch(typ(pl))
1577 : {
1578 679 : case t_VEC: RgV_check_ZV(pl,fun); pl = ZV_to_zv(pl);
1579 : /* fall through */
1580 1337 : case t_VECSMALL: break;
1581 0 : default: pari_err_TYPE(fun,pl);
1582 : }
1583 1337 : return pl;
1584 : }
1585 :
1586 : static int
1587 5068 : is_chineseinit(GEN x)
1588 : {
1589 : GEN fa, pl;
1590 : long l;
1591 5068 : if (typ(x) != t_VEC || lg(x)!=3) return 0;
1592 3766 : fa = gel(x,1);
1593 3766 : pl = gel(x,2);
1594 3766 : if (typ(fa) != t_VEC || typ(pl) != t_VEC) return 0;
1595 1645 : l = lg(fa);
1596 1645 : if (l != 1)
1597 : {
1598 1624 : if (l != 3 || typ(gel(fa,1)) != t_MAT || typ(gel(fa,2)) != t_VEC)
1599 7 : return 0;
1600 : }
1601 1638 : l = lg(pl);
1602 1638 : if (l != 1)
1603 : {
1604 553 : if (l != 6 || typ(gel(pl,3)) != t_MAT || typ(gel(pl,1)) != t_VECSMALL
1605 553 : || typ(gel(pl,2)) != t_VECSMALL)
1606 0 : return 0;
1607 : }
1608 1638 : return 1;
1609 : }
1610 :
1611 : /* nf a true 'nf' */
1612 : static GEN
1613 2184 : chineseinit_i(GEN nf, GEN fa, GEN w, GEN dw)
1614 : {
1615 2184 : const char *fun = "idealchineseinit";
1616 2184 : GEN archp = NULL, pl = NULL;
1617 2184 : switch(typ(fa))
1618 : {
1619 : case t_VEC:
1620 1337 : if (is_chineseinit(fa))
1621 : {
1622 0 : if (dw) pari_err_DOMAIN(fun, "denom(y)", "!=", gen_1, w);
1623 0 : return fa;
1624 : }
1625 1337 : if (lg(fa) != 3) pari_err_TYPE(fun, fa);
1626 : /* of the form [x,s] */
1627 1337 : pl = pl_normalize(nf, gel(fa,2));
1628 1337 : fa = gel(fa,1);
1629 1337 : archp = vecsmall01_to_indices(pl);
1630 : /* keep pr_init, reset pl */
1631 1337 : if (is_chineseinit(fa)) { fa = gel(fa,1); break; }
1632 : /* fall through */
1633 : case t_MAT: /* factorization? */
1634 2121 : if (isprfact(fa)) { fa = pr_init(nf, fa, w, dw); break; }
1635 0 : default: pari_err_TYPE(fun,fa);
1636 : }
1637 :
1638 2184 : if (pl)
1639 : {
1640 1337 : GEN F = (lg(fa) == 1)? NULL: gel(fa,1);
1641 1337 : long i, r = lg(archp);
1642 1337 : GEN signs = cgetg(r, t_VECSMALL);
1643 1337 : for (i = 1; i < r; i++) signs[i] = (pl[archp[i]] < 0)? 1: 0;
1644 1337 : pl = setsigns_init(nf, archp, F, signs);
1645 : }
1646 : else
1647 847 : pl = cgetg(1,t_VEC);
1648 2184 : return mkvec2(fa, pl);
1649 : }
1650 :
1651 : /* Given a prime ideal factorization x, possibly with 0 or negative exponents,
1652 : * and a vector w of elements of nf, gives b such that
1653 : * v_p(b-w_p)>=v_p(x) for all prime ideals p in the ideal factorization
1654 : * and v_p(b)>=0 for all other p, using the standard proof given in GTM 138. */
1655 : GEN
1656 3759 : idealchinese(GEN nf, GEN x, GEN w)
1657 : {
1658 3759 : const char *fun = "idealchinese";
1659 3759 : pari_sp av = avma;
1660 : GEN x1, x2, s, dw, F;
1661 :
1662 3759 : nf = checknf(nf);
1663 3759 : if (!w) return gerepilecopy(av, chineseinit_i(nf,x,NULL,NULL));
1664 :
1665 2394 : if (typ(w) != t_VEC) pari_err_TYPE(fun,w);
1666 2394 : w = Q_remove_denom(matalgtobasis(nf,w), &dw);
1667 2394 : if (!is_chineseinit(x)) x = chineseinit_i(nf,x,w,dw);
1668 : /* x is a 'chineseinit' */
1669 2394 : x1 = gel(x,1); s = NULL;
1670 2394 : if (lg(x1) == 1) F = NULL;
1671 : else
1672 : {
1673 2373 : GEN U = gel(x1,2);
1674 2373 : long i, r = lg(w);
1675 2373 : F = gel(x1,1);
1676 5838 : for (i=1; i<r; i++)
1677 3465 : if (!gequal0(gel(w,i)))
1678 : {
1679 2926 : GEN t = nfmuli(nf, gel(U,i), gel(w,i));
1680 2926 : s = s? ZC_add(s,t): t;
1681 : }
1682 2373 : if (s) s = ZC_reducemodmatrix(s, F);
1683 : }
1684 2394 : if (!s) { s = zerocol(nf_get_degree(nf)); dw = NULL; }
1685 :
1686 2394 : x2 = gel(x,2);
1687 2394 : if (lg(x2) != 1) s = nfsetsigns(nf, gel(x2,1), s, x2);
1688 2394 : if (dw) s = RgC_Rg_div(s,dw);
1689 2394 : return gerepileupto(av, s);
1690 : }
1691 :
1692 : /*************************************************************************/
1693 : /** **/
1694 : /** (Z_K/I)^* **/
1695 : /** **/
1696 : /*************************************************************************/
1697 : GEN
1698 1337 : vecsmall01_to_indices(GEN v)
1699 : {
1700 1337 : long i, k, l = lg(v);
1701 1337 : GEN p = new_chunk(l) + l;
1702 3766 : for (k=1, i=l-1; i; i--)
1703 2429 : if (v[i]) { *--p = i; k++; }
1704 1337 : *--p = evallg(k) | evaltyp(t_VECSMALL);
1705 1337 : avma = (pari_sp)p; return p;
1706 : }
1707 : GEN
1708 361389 : vec01_to_indices(GEN v)
1709 : {
1710 : long i, k, l;
1711 : GEN p;
1712 :
1713 361389 : switch (typ(v))
1714 : {
1715 346689 : case t_VECSMALL: return v;
1716 14700 : case t_VEC: break;
1717 0 : default: pari_err_TYPE("vec01_to_indices",v);
1718 : }
1719 14700 : l = lg(v);
1720 14700 : p = new_chunk(l) + l;
1721 42427 : for (k=1, i=l-1; i; i--)
1722 27727 : if (signe(gel(v,i))) { *--p = i; k++; }
1723 14700 : *--p = evallg(k) | evaltyp(t_VECSMALL);
1724 14700 : avma = (pari_sp)p; return p;
1725 : }
1726 : GEN
1727 4830 : indices_to_vec01(GEN p, long r)
1728 : {
1729 4830 : long i, l = lg(p);
1730 4830 : GEN v = zerovec(r);
1731 4830 : for (i = 1; i < l; i++) gel(v, p[i]) = gen_1;
1732 4830 : return v;
1733 : }
1734 :
1735 : /* return (column) vector of R1 signatures of x (0 or 1) */
1736 : GEN
1737 346689 : nfsign_arch(GEN nf, GEN x, GEN arch)
1738 : {
1739 346689 : GEN sarch, M, V, archp = vec01_to_indices(arch);
1740 346689 : long i, s, np, n = lg(archp)-1;
1741 : pari_sp av;
1742 :
1743 346689 : if (!n) return cgetg(1,t_VECSMALL);
1744 345660 : nf = checknf(nf);
1745 345660 : if (typ(x) == t_MAT)
1746 : { /* factorisation */
1747 98207 : GEN g = gel(x,1), e = gel(x,2);
1748 98207 : V = zero_zv(n);
1749 287719 : for (i=1; i<lg(g); i++)
1750 189512 : if (mpodd(gel(e,i)))
1751 164144 : Flv_add_inplace(V, nfsign_arch(nf,gel(g,i),archp), 2);
1752 98207 : avma = (pari_sp)V; return V;
1753 : }
1754 247453 : av = avma; V = cgetg(n+1,t_VECSMALL);
1755 247453 : x = nf_to_scalar_or_basis(nf, x);
1756 247453 : switch(typ(x))
1757 : {
1758 : case t_INT:
1759 64502 : s = signe(x);
1760 64502 : if (!s) pari_err_DOMAIN("nfsign_arch","element","=",gen_0,x);
1761 64502 : avma = av; return const_vecsmall(n, (s < 0)? 1: 0);
1762 : case t_FRAC:
1763 35 : s = signe(gel(x,1));
1764 35 : avma = av; return const_vecsmall(n, (s < 0)? 1: 0);
1765 : }
1766 182916 : x = Q_primpart(x); M = nf_get_M(nf); sarch = NULL; np = -1;
1767 493574 : for (i = 1; i <= n; i++)
1768 : {
1769 311605 : long s = eval_sign(M, x, archp[i]);
1770 311605 : if (s < 0) /* failure */
1771 : {
1772 947 : long ni, r1 = nf_get_r1(nf);
1773 : GEN xi;
1774 947 : if (np < 0)
1775 : {
1776 947 : np = num_positive(nf, x);
1777 947 : if (np == 0) { avma = av; return const_vecsmall(n, 1); }
1778 809 : if (np == r1){ avma = av; return const_vecsmall(n, 0); }
1779 494 : sarch = nfarchstar(nf, NULL, identity_perm(r1));
1780 : }
1781 494 : xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
1782 494 : xi = Q_primpart(xi);
1783 494 : ni = num_positive(nf, nfmuli(nf,x,xi));
1784 494 : if (ni == 0) { avma = av; V = const_vecsmall(n, 1); V[i] = 0; return V; }
1785 356 : if (ni == r1){ avma = av; V = const_vecsmall(n, 0); V[i] = 1; return V; }
1786 0 : s = ni < np? 0: 1;
1787 : }
1788 310658 : V[i] = s;
1789 : }
1790 181969 : avma = (pari_sp)V; return V;
1791 : }
1792 : static void
1793 2926 : chk_ind(const char *s, long i, long r1)
1794 : {
1795 2926 : if (i <= 0) pari_err_DOMAIN(s, "index", "<=", gen_0, stoi(i));
1796 2919 : if (i > r1) pari_err_DOMAIN(s, "index", ">", utoi(r1), utoi(i));
1797 2898 : }
1798 : GEN
1799 1715 : nfeltsign(GEN nf, GEN x, GEN ind0)
1800 : {
1801 1715 : pari_sp av = avma;
1802 : long i, l, r1;
1803 : GEN v, ind;
1804 1715 : nf = checknf(nf);
1805 1715 : r1 = nf_get_r1(nf);
1806 1715 : x = nf_to_scalar_or_basis(nf, x);
1807 1715 : if (!ind0) ind0 = identity_perm(r1);
1808 1715 : switch(typ(ind0))
1809 : {
1810 : case t_INT: case t_VEC: case t_COL:
1811 56 : ind = gtovecsmall(ind0); break;
1812 : case t_VECSMALL:
1813 1659 : ind = ind0; break;
1814 : default:
1815 0 : pari_err_TYPE("nfeltsign",ind0);
1816 : return NULL; /* LCOV_EXCL_LINE */
1817 : }
1818 1715 : l = lg(ind);
1819 1715 : for (i = 1; i < l; i++) chk_ind("nfeltsign", ind[i], r1);
1820 1694 : if (typ(x) != t_COL)
1821 : {
1822 : GEN s;
1823 35 : switch(gsigne(x))
1824 : {
1825 14 : case -1:s = gen_m1; break;
1826 14 : case 1: s = gen_1; break;
1827 7 : default: s = gen_0; break;
1828 : }
1829 35 : avma = av;
1830 35 : return typ(ind0) == t_INT? s: const_vec(l-1, s);
1831 : }
1832 1659 : v = nfsign_arch(nf, x, ind);
1833 1659 : if (typ(ind0) == t_INT) { avma = av; return v[1]? gen_m1: gen_1; }
1834 1652 : settyp(v, t_VEC);
1835 1652 : for (i = 1; i < l; i++) gel(v,i) = v[i]? gen_m1: gen_1;
1836 1652 : return gerepileupto(av, v);
1837 : }
1838 :
1839 : GEN
1840 63 : nfeltembed(GEN nf, GEN x, GEN ind0, long prec0)
1841 : {
1842 63 : pari_sp av = avma;
1843 : long i, e, l, r1, r2, prec, prec1;
1844 : GEN v, ind, cx;
1845 63 : nf = checknf(nf); nf_get_sign(nf,&r1,&r2);
1846 63 : x = nf_to_scalar_or_basis(nf, x);
1847 56 : if (!ind0) ind0 = identity_perm(r1+r2);
1848 56 : switch(typ(ind0))
1849 : {
1850 : case t_INT: case t_VEC: case t_COL:
1851 42 : ind = gtovecsmall(ind0); break;
1852 : case t_VECSMALL:
1853 14 : ind = ind0; break;
1854 : default:
1855 0 : pari_err_TYPE("nfeltsign",ind0);
1856 : return NULL; /* LCOV_EXCL_LINE */
1857 : }
1858 56 : l = lg(ind);
1859 56 : for (i = 1; i < l; i++) chk_ind("nfeltembed", ind[i], r1+r2);
1860 49 : if (typ(x) != t_COL)
1861 : {
1862 0 : if (typ(ind0) != t_INT) x = const_vec(l-1, x);
1863 0 : return gerepilecopy(av, x);
1864 : }
1865 49 : x = Q_primitive_part(x, &cx);
1866 49 : prec1 = prec0; e = gexpo(x);
1867 49 : if (e > 8) prec1 += nbits2extraprec(e);
1868 49 : prec = prec1;
1869 49 : if (nf_get_prec(nf) < prec) nf = nfnewprec_shallow(nf, prec);
1870 49 : v = cgetg(l, t_VEC);
1871 : for(;;)
1872 : {
1873 56 : GEN M = nf_get_M(nf);
1874 140 : for (i = 1; i < l; i++)
1875 : {
1876 91 : GEN t = nfembed_i(M, x, ind[i]);
1877 91 : long e = gexpo(t);
1878 91 : if (gequal0(t) || precision(t) < prec0
1879 91 : || (e < 0 && prec < prec1 + nbits2extraprec(-e)) ) break;
1880 84 : if (cx) t = gmul(t, cx);
1881 84 : gel(v,i) = t;
1882 : }
1883 56 : if (i == l) break;
1884 7 : prec = precdbl(prec);
1885 7 : if (DEBUGLEVEL>1) pari_warn(warnprec,"eltnfembed", prec);
1886 7 : nf = nfnewprec_shallow(nf, prec);
1887 7 : }
1888 49 : if (typ(ind0) == t_INT) v = gel(v,1);
1889 49 : return gerepilecopy(av, v);
1890 : }
1891 :
1892 : /* return the vector of signs of x; the matrix of such if x is a vector
1893 : * of nf elements */
1894 : GEN
1895 1323 : nfsign(GEN nf, GEN x)
1896 : {
1897 : long i, l;
1898 : GEN archp, S;
1899 :
1900 1323 : nf = checknf(nf);
1901 1323 : archp = identity_perm( nf_get_r1(nf) );
1902 1323 : if (typ(x) != t_VEC) return nfsign_arch(nf, x, archp);
1903 182 : l = lg(x); S = cgetg(l, t_MAT);
1904 182 : for (i=1; i<l; i++) gel(S,i) = nfsign_arch(nf, gel(x,i), archp);
1905 182 : return S;
1906 : }
1907 :
1908 : /* x integral elt, A integral ideal in HNF; reduce x mod A */
1909 : static GEN
1910 604086 : zk_modHNF(GEN x, GEN A)
1911 604086 : { return (typ(x) == t_COL)? ZC_hnfrem(x, A): modii(x, gcoeff(A,1,1)); }
1912 :
1913 : /* given an element x in Z_K and an integral ideal y in HNF, coprime with x,
1914 : outputs an element inverse of x modulo y */
1915 : GEN
1916 140 : nfinvmodideal(GEN nf, GEN x, GEN y)
1917 : {
1918 140 : pari_sp av = avma;
1919 140 : GEN a, yZ = gcoeff(y,1,1);
1920 :
1921 140 : if (is_pm1(yZ)) return gen_0;
1922 140 : x = nf_to_scalar_or_basis(nf, x);
1923 140 : if (typ(x) == t_INT) return gerepileupto(av, Fp_inv(x, yZ));
1924 :
1925 70 : a = hnfmerge_get_1(idealhnf_principal(nf,x), y);
1926 70 : if (!a) pari_err_INV("nfinvmodideal", x);
1927 70 : return gerepileupto(av, zk_modHNF(nfdiv(nf,a,x), y));
1928 : }
1929 :
1930 : static GEN
1931 277124 : nfsqrmodideal(GEN nf, GEN x, GEN id)
1932 277124 : { return zk_modHNF(nfsqri(nf,x), id); }
1933 : static GEN
1934 735860 : nfmulmodideal(GEN nf, GEN x, GEN y, GEN id)
1935 735860 : { return x? zk_modHNF(nfmuli(nf,x,y), id): y; }
1936 : /* assume x integral, k integer, A in HNF */
1937 : GEN
1938 476237 : nfpowmodideal(GEN nf,GEN x,GEN k,GEN A)
1939 : {
1940 476237 : long s = signe(k);
1941 : pari_sp av;
1942 : GEN y;
1943 :
1944 476237 : if (!s) return gen_1;
1945 476237 : av = avma;
1946 476237 : x = nf_to_scalar_or_basis(nf, x);
1947 476237 : if (typ(x) != t_COL) return Fp_pow(x, k, gcoeff(A,1,1));
1948 226327 : if (s < 0) { x = nfinvmodideal(nf, x,A); k = absi(k); }
1949 226327 : for(y = NULL;;)
1950 : {
1951 503451 : if (mpodd(k)) y = nfmulmodideal(nf,y,x,A);
1952 503451 : k = shifti(k,-1); if (!signe(k)) break;
1953 277124 : x = nfsqrmodideal(nf,x,A);
1954 277124 : }
1955 226327 : return gerepileupto(av, y);
1956 : }
1957 :
1958 : /* a * g^n mod id */
1959 : static GEN
1960 423100 : nfmulpowmodideal(GEN nf, GEN a, GEN g, GEN n, GEN id)
1961 : {
1962 423100 : return nfmulmodideal(nf, a, nfpowmodideal(nf,g,n,id), id);
1963 : }
1964 :
1965 : /* assume (num(g[i]), id) = 1 for all i. Return prod g[i]^e[i] mod id.
1966 : * EX = multiple of exponent of (O_K/id)^* */
1967 : GEN
1968 192168 : famat_to_nf_modideal_coprime(GEN nf, GEN g, GEN e, GEN id, GEN EX)
1969 : {
1970 192168 : GEN EXo2, plus = NULL, minus = NULL, idZ = gcoeff(id,1,1);
1971 192168 : long i, lx = lg(g);
1972 :
1973 192168 : if (is_pm1(idZ)) return gen_1; /* id = Z_K */
1974 192168 : EXo2 = (expi(EX) > 10)? shifti(EX,-1): NULL;
1975 871484 : for (i = 1; i < lx; i++)
1976 : {
1977 679316 : GEN h, n = centermodii(gel(e,i), EX, EXo2);
1978 679316 : long sn = signe(n);
1979 679316 : if (!sn) continue;
1980 :
1981 324248 : h = nf_to_scalar_or_basis(nf, gel(g,i));
1982 324248 : switch(typ(h))
1983 : {
1984 199620 : case t_INT: break;
1985 : case t_FRAC:
1986 0 : h = Fp_div(gel(h,1), gel(h,2), idZ); break;
1987 : default:
1988 : {
1989 : GEN dh;
1990 124628 : h = Q_remove_denom(h, &dh);
1991 124628 : if (dh) h = FpC_Fp_mul(h, Fp_inv(dh,idZ), idZ);
1992 : }
1993 : }
1994 324248 : if (sn > 0)
1995 322848 : plus = nfmulpowmodideal(nf, plus, h, n, id);
1996 : else /* sn < 0 */
1997 1400 : minus = nfmulpowmodideal(nf, minus, h, absi(n), id);
1998 : }
1999 192168 : if (minus) plus = nfmulmodideal(nf, plus, nfinvmodideal(nf,minus,id), id);
2000 192168 : return plus? plus: gen_1;
2001 : }
2002 :
2003 : /* given 2 integral ideals x, y in HNF s.t x | y | x^2, compute (1+x)/(1+y) in
2004 : * the form [[cyc],[gen], U], where U := ux^-1 as a pair [ZM, denom(U)] */
2005 : static GEN
2006 20790 : zidealij(GEN x, GEN y)
2007 : {
2008 20790 : GEN U, G, cyc, xp = gcoeff(x,1,1), xi = hnf_invscale(x, xp);
2009 : long j, N;
2010 :
2011 : /* x^(-1) y = relations between the 1 + x_i (HNF) */
2012 20790 : cyc = ZM_snf_group(ZM_Z_divexact(ZM_mul(xi, y), xp), &U, &G);
2013 20790 : N = lg(cyc); G = ZM_mul(x,G); settyp(G, t_VEC); /* new generators */
2014 77609 : for (j=1; j<N; j++)
2015 : {
2016 56819 : GEN c = gel(G,j);
2017 56819 : gel(c,1) = addiu(gel(c,1), 1); /* 1 + g_j */
2018 56819 : if (ZV_isscalar(c)) gel(G,j) = gel(c,1);
2019 : }
2020 20790 : return mkvec4(cyc, G, ZM_mul(U,xi), xp);
2021 : }
2022 :
2023 : /* lg(x) > 1, x + 1; shallow */
2024 : static GEN
2025 6265 : ZC_add1(GEN x)
2026 : {
2027 6265 : long i, l = lg(x);
2028 6265 : GEN y = cgetg(l, t_COL);
2029 6265 : for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
2030 6265 : gel(y,1) = addiu(gel(x,1), 1); return y;
2031 : }
2032 : /* lg(x) > 1, x - 1; shallow */
2033 : static GEN
2034 3766 : ZC_sub1(GEN x)
2035 : {
2036 3766 : long i, l = lg(x);
2037 3766 : GEN y = cgetg(l, t_COL);
2038 3766 : for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
2039 3766 : gel(y,1) = subiu(gel(x,1), 1); return y;
2040 : }
2041 :
2042 : /* x,y are t_INT or ZC */
2043 : static GEN
2044 0 : zkadd(GEN x, GEN y)
2045 : {
2046 0 : long tx = typ(x);
2047 0 : if (tx == typ(y))
2048 0 : return tx == t_INT? addii(x,y): ZC_add(x,y);
2049 : else
2050 0 : return tx == t_INT? ZC_Z_add(y,x): ZC_Z_add(x,y);
2051 : }
2052 : /* x a t_INT or ZC, x+1; shallow */
2053 : static GEN
2054 12215 : zkadd1(GEN x)
2055 : {
2056 12215 : long tx = typ(x);
2057 12215 : return tx == t_INT? addiu(x,1): ZC_add1(x);
2058 : }
2059 : /* x a t_INT or ZC, x-1; shallow */
2060 : static GEN
2061 12215 : zksub1(GEN x)
2062 : {
2063 12215 : long tx = typ(x);
2064 12215 : return tx == t_INT? subiu(x,1): ZC_sub1(x);
2065 : }
2066 : /* x,y are t_INT or ZC; x - y */
2067 : static GEN
2068 0 : zksub(GEN x, GEN y)
2069 : {
2070 0 : long tx = typ(x), ty = typ(y);
2071 0 : if (tx == ty)
2072 0 : return tx == t_INT? subii(x,y): ZC_sub(x,y);
2073 : else
2074 0 : return tx == t_INT? Z_ZC_sub(x,y): ZC_Z_sub(x,y);
2075 : }
2076 : /* x is t_INT or ZM (mult. map), y is t_INT or ZC; x * y */
2077 : static GEN
2078 12215 : zkmul(GEN x, GEN y)
2079 : {
2080 12215 : long tx = typ(x), ty = typ(y);
2081 12215 : if (ty == t_INT)
2082 8449 : return tx == t_INT? mulii(x,y): ZC_Z_mul(gel(x,1),y);
2083 : else
2084 3766 : return tx == t_INT? ZC_Z_mul(y,x): ZM_ZC_mul(x,y);
2085 : }
2086 :
2087 : /* (U,V) = 1 coprime ideals. Want z = x mod U, = y mod V; namely
2088 : * z =vx + uy = v(x-y) + y, where u + v = 1, u in U, v in V.
2089 : * zkc = [v, UV], v a t_INT or ZM (mult. by v map), UV a ZM (ideal in HNF);
2090 : * shallow */
2091 : GEN
2092 0 : zkchinese(GEN zkc, GEN x, GEN y)
2093 : {
2094 0 : GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd(zkmul(v, zksub(x,y)), y);
2095 0 : return zk_modHNF(z, UV);
2096 : }
2097 : /* special case z = x mod U, = 1 mod V; shallow */
2098 : GEN
2099 12215 : zkchinese1(GEN zkc, GEN x)
2100 : {
2101 12215 : GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd1(zkmul(v, zksub1(x)));
2102 12215 : return (typ(z) == t_INT)? z: ZC_hnfrem(z, UV);
2103 : }
2104 : static GEN
2105 11011 : zkVchinese1(GEN zkc, GEN v)
2106 : {
2107 : long i, ly;
2108 11011 : GEN y = cgetg_copy(v, &ly);
2109 11011 : for (i=1; i<ly; i++) gel(y,i) = zkchinese1(zkc, gel(v,i));
2110 11011 : return y;
2111 : }
2112 :
2113 : /* prepare to solve z = x (mod A), z = y mod (B) [zkchinese or zkchinese1] */
2114 : GEN
2115 10752 : zkchineseinit(GEN nf, GEN A, GEN B, GEN AB)
2116 : {
2117 : GEN v;
2118 : long e;
2119 10752 : nf = checknf(nf);
2120 10752 : v = idealaddtoone_raw(nf, A, B);
2121 10752 : if ((e = gexpo(v)) > 5)
2122 : {
2123 455 : GEN b = (typ(v) == t_COL)? v: scalarcol_shallow(v, nf_get_degree(nf));
2124 455 : b= ZC_reducemodlll(b, AB);
2125 455 : if (gexpo(b) < e) v = b;
2126 : }
2127 10752 : return mkvec2(zk_scalar_or_multable(nf,v), AB);
2128 : }
2129 : /* prepare to solve z = x (mod A), z = 1 mod (B)
2130 : * and then z = 1 (mod A), z = y mod (B) [zkchinese1 twice] */
2131 : static GEN
2132 259 : zkchinese1init2(GEN nf, GEN A, GEN B, GEN AB)
2133 : {
2134 259 : GEN zkc = zkchineseinit(nf, A, B, AB);
2135 259 : GEN mv = gel(zkc,1), mu;
2136 259 : if (typ(mv) == t_INT) return mkvec2(zkc, mkvec2(subui(1,mv),AB));
2137 35 : mu = RgM_Rg_add_shallow(ZM_neg(mv), gen_1);
2138 35 : return mkvec2(mkvec2(mv,AB), mkvec2(mu,AB));
2139 : }
2140 :
2141 : static GEN
2142 355605 : apply_U(GEN L, GEN a)
2143 : {
2144 355605 : GEN e, U = gel(L,3), dU = gel(L,4);
2145 355605 : if (typ(a) == t_INT)
2146 130421 : e = ZC_Z_mul(gel(U,1), subiu(a, 1));
2147 : else
2148 : { /* t_COL */
2149 225184 : GEN t = gel(a,1);
2150 225184 : gel(a,1) = subiu(gel(a,1), 1); /* a -= 1 */
2151 225184 : e = ZM_ZC_mul(U, a);
2152 225184 : gel(a,1) = t; /* restore */
2153 : }
2154 355605 : return gdiv(e, dU);
2155 : }
2156 :
2157 : /* true nf; vectors of [[cyc],[g],U.X^-1]. Assume k > 1. */
2158 : static GEN
2159 14126 : principal_units(GEN nf, GEN pr, long k, GEN prk)
2160 : {
2161 : GEN list, prb;
2162 14126 : ulong mask = quadratic_prec_mask(k);
2163 14126 : long a = 1;
2164 :
2165 14126 : if (DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",k,pr);
2166 14126 : prb = pr_hnf(nf,pr);
2167 14126 : list = vectrunc_init(k);
2168 49042 : while (mask > 1)
2169 : {
2170 20790 : GEN pra = prb;
2171 20790 : long b = a << 1;
2172 :
2173 20790 : if (mask & 1) b--;
2174 20790 : mask >>= 1;
2175 : /* compute 1 + pr^a / 1 + pr^b, 2a <= b */
2176 20790 : if(DEBUGLEVEL>3) err_printf(" treating a = %ld, b = %ld\n",a,b);
2177 20790 : prb = (b >= k)? prk: idealpows(nf,pr,b);
2178 20790 : vectrunc_append(list, zidealij(pra, prb));
2179 20790 : a = b;
2180 : }
2181 14126 : return list;
2182 : }
2183 : /* a = 1 mod (pr) return log(a) on local-gens of 1+pr/1+pr^k */
2184 : static GEN
2185 226112 : log_prk1(GEN nf, GEN a, long nh, GEN L2, GEN prk)
2186 : {
2187 226112 : GEN y = cgetg(nh+1, t_COL);
2188 226112 : long j, iy, c = lg(L2)-1;
2189 581710 : for (j = iy = 1; j <= c; j++)
2190 : {
2191 355605 : GEN L = gel(L2,j), cyc = gel(L,1), gen = gel(L,2), E = apply_U(L,a);
2192 355605 : long i, nc = lg(cyc)-1;
2193 355605 : int last = (j == c);
2194 1243679 : for (i = 1; i <= nc; i++, iy++)
2195 : {
2196 888081 : GEN t, e = gel(E,i);
2197 888081 : if (typ(e) != t_INT) pari_err_COPRIME("zlog_prk1", a, prk);
2198 888074 : t = Fp_neg(e, gel(cyc,i));
2199 888074 : gel(y,iy) = negi(t);
2200 888074 : if (!last && signe(t)) a = nfmulpowmodideal(nf, a, gel(gen,i), t, prk);
2201 : }
2202 : }
2203 226105 : return y;
2204 : }
2205 : /* true nf */
2206 : static GEN
2207 5768 : principal_units_relations(GEN nf, GEN L2, GEN prk, long nh)
2208 : {
2209 5768 : GEN h = cgetg(nh+1,t_MAT);
2210 5768 : long ih, j, c = lg(L2)-1;
2211 18200 : for (j = ih = 1; j <= c; j++)
2212 : {
2213 12432 : GEN L = gel(L2,j), F = gel(L,1), G = gel(L,2);
2214 12432 : long k, lG = lg(G);
2215 51534 : for (k = 1; k < lG; k++,ih++)
2216 : { /* log(g^f) mod pr^e */
2217 39102 : GEN a = nfpowmodideal(nf,gel(G,k),gel(F,k),prk);
2218 39102 : gel(h,ih) = ZC_neg(log_prk1(nf, a, nh, L2, prk));
2219 39102 : gcoeff(h,ih,ih) = gel(F,k);
2220 : }
2221 : }
2222 5768 : return h;
2223 : }
2224 : /* true nf; e > 1; multiplicative group (1 + pr) / (1 + pr^k),
2225 : * prk = pr^k or NULL */
2226 : static GEN
2227 14126 : idealprincipalunits_i(GEN nf, GEN pr, long k, GEN *pU)
2228 : {
2229 14126 : GEN cyc, gen, L2, prk = idealpows(nf, pr, k);
2230 :
2231 14126 : L2 = principal_units(nf, pr, k, prk);
2232 14126 : if (k == 2)
2233 : {
2234 8358 : GEN L = gel(L2,1);
2235 8358 : cyc = gel(L,1);
2236 8358 : gen = gel(L,2);
2237 8358 : if (pU) *pU = matid(lg(gen)-1);
2238 : }
2239 : else
2240 : {
2241 5768 : long c = lg(L2), j;
2242 5768 : GEN EX, h, Ui, vg = cgetg(c, t_VEC);
2243 5768 : for (j = 1; j < c; j++) gel(vg, j) = gmael(L2,j,2);
2244 5768 : vg = shallowconcat1(vg);
2245 5768 : h = principal_units_relations(nf, L2, prk, lg(vg)-1);
2246 5768 : h = ZM_hnfall_i(h, NULL, 0);
2247 5768 : cyc = ZM_snf_group(h, pU, &Ui);
2248 5768 : c = lg(Ui); gen = cgetg(c, t_VEC); EX = gel(cyc,1);
2249 32592 : for (j = 1; j < c; j++)
2250 26824 : gel(gen,j) = famat_to_nf_modideal_coprime(nf, vg, gel(Ui,j), prk, EX);
2251 : }
2252 14126 : return mkvec4(cyc, gen, prk, L2);
2253 : }
2254 : GEN
2255 112 : idealprincipalunits(GEN nf, GEN pr, long k)
2256 : {
2257 : pari_sp av;
2258 : GEN v;
2259 112 : nf = checknf(nf);
2260 112 : if (k == 1) { checkprid(pr); retmkvec3(gen_1,cgetg(1,t_VEC),cgetg(1,t_VEC)); }
2261 105 : av = avma; v = idealprincipalunits_i(nf, pr, k, NULL);
2262 105 : return gerepilecopy(av, mkvec3(powiu(pr_norm(pr), k-1), gel(v,1), gel(v,2)));
2263 : }
2264 :
2265 : /* true nf; given an ideal pr^k dividing an integral ideal x (in HNF form)
2266 : * compute an 'sprk', the structure of G = (Z_K/pr^k)^* [ x = NULL for x=pr^k ]
2267 : * Return a vector with at least 4 components [cyc],[gen],[HNF pr^k,pr,k],ff,
2268 : * where
2269 : * cyc : type of G as abelian group (SNF)
2270 : * gen : generators of G, coprime to x
2271 : * pr^k: in HNF
2272 : * ff : data for log_g in (Z_K/pr)^*
2273 : * Two extra components are present iff k > 1: L2, U
2274 : * L2 : list of data structures to compute local DL in (Z_K/pr)^*,
2275 : * and 1 + pr^a/ 1 + pr^b for various a < b <= min(2a, k)
2276 : * U : base change matrices to convert a vector of local DL to DL wrt gen */
2277 : static GEN
2278 31696 : sprkinit(GEN nf, GEN pr, GEN gk, GEN x)
2279 : {
2280 : GEN T, p, modpr, cyc, gen, g, g0, ord0, A, prk, U, L2;
2281 31696 : long k = itos(gk), f = pr_get_f(pr);
2282 :
2283 31696 : if(DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",k,pr);
2284 31696 : modpr = nf_to_Fq_init(nf, &pr,&T,&p);
2285 : /* (Z_K / pr)^* */
2286 31696 : if (f == 1)
2287 : {
2288 22673 : g0 = g = pgener_Fp(p);
2289 22673 : ord0 = get_arith_ZZM(subiu(p,1));
2290 : }
2291 : else
2292 : {
2293 9023 : g0 = g = gener_FpXQ(T,p, &ord0);
2294 9023 : g = Fq_to_nf(g, modpr);
2295 9023 : if (typ(g) == t_POL) g = poltobasis(nf, g);
2296 : }
2297 31696 : A = gel(ord0, 1); /* Norm(pr)-1 */
2298 31696 : if (k == 1)
2299 : {
2300 17675 : cyc = mkvec(A);
2301 17675 : gen = mkvec(g);
2302 17675 : prk = pr_hnf(nf,pr);
2303 17675 : L2 = U = NULL;
2304 : }
2305 : else
2306 : { /* local-gens of (1 + pr)/(1 + pr^k) = SNF-gens * U */
2307 : GEN AB, B, u, v, w;
2308 : long j, l;
2309 14021 : w = idealprincipalunits_i(nf, pr, k, &U);
2310 : /* incorporate (Z_K/pr)^*, order A coprime to B = expo(1+pr/1+pr^k)*/
2311 14021 : cyc = leafcopy(gel(w,1)); B = gel(cyc,1); AB = mulii(A,B);
2312 14021 : gen = leafcopy(gel(w,2));
2313 14021 : prk = gel(w,3);
2314 14021 : g = nfpowmodideal(nf, g, B, prk);
2315 14021 : g0 = Fq_pow(g0, modii(B,A), T, p); /* update primitive root */
2316 14021 : L2 = mkvec3(A, g, gel(w,4));
2317 14021 : gel(cyc,1) = AB;
2318 14021 : gel(gen,1) = nfmulmodideal(nf, gel(gen,1), g, prk);
2319 14021 : u = mulii(Fp_inv(A,B), A);
2320 14021 : v = subui(1, u); l = lg(U);
2321 14021 : for (j = 1; j < l; j++) gcoeff(U,1,j) = Fp_mul(u, gcoeff(U,1,j), AB);
2322 14021 : U = mkvec2(Rg_col_ei(v, lg(gen)-1, 1), U);
2323 : }
2324 : /* local-gens of (Z_K/pr^k)^* = SNF-gens * U */
2325 31696 : if (x)
2326 : {
2327 10493 : GEN uv = zkchineseinit(nf, idealdivpowprime(nf,x,pr,gk), prk, x);
2328 10493 : gen = zkVchinese1(uv, gen);
2329 : }
2330 31696 : return mkvecn(U? 6: 4, cyc, gen, prk, mkvec3(modpr,g0,ord0), L2, U);
2331 : }
2332 : static GEN
2333 347769 : sprk_get_cyc(GEN s) { return gel(s,1); }
2334 : static GEN
2335 110268 : sprk_get_expo(GEN s)
2336 : {
2337 110268 : GEN cyc = sprk_get_cyc(s);
2338 110268 : return lg(cyc) == 1? gen_1: gel(cyc, 1);
2339 : }
2340 : static GEN
2341 25634 : sprk_get_gen(GEN s) { return gel(s,2); }
2342 : static GEN
2343 297278 : sprk_get_prk(GEN s) { return gel(s,3); }
2344 : static GEN
2345 387737 : sprk_get_ff(GEN s) { return gel(s,4); }
2346 : static GEN
2347 129266 : sprk_get_pr(GEN s) { GEN ff = gel(s,4); return modpr_get_pr(gel(ff,1)); }
2348 : /* A = Npr-1, <g> = (Z_K/pr)^*, L2 to 1 + pr / 1 + pr^k */
2349 : static void
2350 197916 : sprk_get_L2(GEN s, GEN *A, GEN *g, GEN *L2)
2351 197916 : { GEN v = gel(s,5); *A = gel(v,1); *g = gel(v,2); *L2 = gel(v,3); }
2352 : static void
2353 187010 : sprk_get_U2(GEN s, GEN *U1, GEN *U2)
2354 187010 : { GEN v = gel(s,6); *U1 = gel(v,1); *U2 = gel(v,2); }
2355 : static int
2356 387737 : sprk_is_prime(GEN s) { return lg(s) == 5; }
2357 :
2358 : static GEN
2359 110268 : famat_zlog_pr(GEN nf, GEN g, GEN e, GEN sprk)
2360 : {
2361 110268 : GEN pr = sprk_get_pr(sprk);
2362 110268 : GEN prk = sprk_get_prk(sprk);
2363 110268 : GEN x = famat_makecoprime(nf, g, e, pr, prk, sprk_get_expo(sprk));
2364 110268 : return zlog_pr(nf, x, sprk);
2365 : }
2366 : /* log_g(a) in (Z_K/pr)^* */
2367 : static GEN
2368 387737 : nf_log(GEN nf, GEN a, GEN ff)
2369 : {
2370 387737 : GEN pr = gel(ff,1), g = gel(ff,2), ord = gel(ff,3);
2371 387737 : GEN T,p, modpr = nf_to_Fq_init(nf, &pr, &T, &p);
2372 387737 : return Fq_log(nf_to_Fq(nf,a,modpr), g, ord, T, p);
2373 : }
2374 : /* a in Z_K (t_COL or t_INT), pr prime ideal, sprk = sprkinit(nf,pr,k,x).
2375 : * return log(a) on SNF-generators of (Z_K/pr^k)^**/
2376 : GEN
2377 388773 : zlog_pr(GEN nf, GEN a, GEN sprk)
2378 : {
2379 : GEN e, prk, A, g, L2, U1, U2, y;
2380 :
2381 388773 : if (typ(a) == t_MAT) return famat_zlog_pr(nf, gel(a,1), gel(a,2), sprk);
2382 :
2383 387737 : e = nf_log(nf, a, sprk_get_ff(sprk));
2384 387737 : if (sprk_is_prime(sprk)) return mkcol(e); /* k = 1 */
2385 187010 : prk = sprk_get_prk(sprk);
2386 187010 : sprk_get_L2(sprk, &A,&g,&L2);
2387 187010 : if (signe(e))
2388 : {
2389 46587 : e = Fp_neg(e, A);
2390 46587 : a = nfmulpowmodideal(nf, a, g, e, prk);
2391 : }
2392 187010 : sprk_get_U2(sprk, &U1,&U2);
2393 187010 : y = ZM_ZC_mul(U2, log_prk1(nf, a, lg(U2)-1, L2, prk));
2394 187003 : if (signe(e)) y = ZC_sub(y, ZC_Z_mul(U1,e));
2395 187003 : return vecmodii(y, sprk_get_cyc(sprk));
2396 : }
2397 : GEN
2398 6062 : zlog_pr_init(GEN nf, GEN pr, long k)
2399 6062 : { return sprkinit(checknf(nf),pr,utoipos(k),NULL);}
2400 : GEN
2401 378 : vzlog_pr(GEN nf, GEN v, GEN sprk)
2402 : {
2403 378 : long l = lg(v), i;
2404 378 : GEN w = cgetg(l, t_MAT);
2405 378 : for (i = 1; i < l; i++) gel(w,i) = zlog_pr(nf, gel(v,i), sprk);
2406 378 : return w;
2407 : }
2408 :
2409 : static GEN
2410 113208 : famat_zlog(GEN nf, GEN fa, GEN sgn, zlog_S *S)
2411 : {
2412 113208 : long i, n0, n = lg(S->U)-1;
2413 : GEN g, e, y;
2414 113208 : if (lg(fa) == 1) return zerocol(n);
2415 113208 : g = gel(fa,1);
2416 113208 : e = gel(fa,2);
2417 113208 : y = cgetg(n+1, t_COL);
2418 113208 : n0 = lg(S->sprk)-1; /* n0 = n (trivial arch. part), or n-1 */
2419 113208 : for (i=1; i <= n0; i++) gel(y,i) = famat_zlog_pr(nf, g, e, gel(S->sprk,i));
2420 113208 : if (n0 != n)
2421 : {
2422 93174 : if (!sgn) sgn = nfsign_arch(nf, fa, S->archp);
2423 93174 : gel(y,n) = Flc_to_ZC(sgn);
2424 : }
2425 113208 : return y;
2426 : }
2427 :
2428 : /* assume that cyclic factors are normalized, in particular != [1] */
2429 : static GEN
2430 26082 : split_U(GEN U, GEN Sprk)
2431 : {
2432 26082 : long t = 0, k, n, l = lg(Sprk);
2433 26082 : GEN vU = cgetg(l+1, t_VEC);
2434 50946 : for (k = 1; k < l; k++)
2435 : {
2436 24864 : n = lg(sprk_get_cyc(gel(Sprk,k))) - 1; /* > 0 */
2437 24864 : gel(vU,k) = vecslice(U, t+1, t+n);
2438 24864 : t += n;
2439 : }
2440 : /* t+1 .. lg(U)-1 */
2441 26082 : n = lg(U) - t - 1; /* can be 0 */
2442 26082 : if (!n) setlg(vU,l); else gel(vU,l) = vecslice(U, t+1, t+n);
2443 26082 : return vU;
2444 : }
2445 :
2446 : void
2447 353324 : init_zlog(zlog_S *S, GEN bid)
2448 : {
2449 353324 : GEN fa2 = bid_get_fact2(bid);
2450 353324 : S->U = bid_get_U(bid);
2451 353324 : S->hU = lg(bid_get_cyc(bid))-1;
2452 353324 : S->archp = bid_get_archp(bid);
2453 353324 : S->sprk = bid_get_sprk(bid);
2454 353324 : S->bid = bid;
2455 353324 : S->P = gel(fa2,1);
2456 353324 : S->k = gel(fa2,2);
2457 353324 : S->no2 = lg(S->P) == lg(gel(bid_get_fact(bid),1));
2458 353324 : }
2459 :
2460 : /* a a t_FRAC/t_INT, reduce mod bid */
2461 : static GEN
2462 7 : Q_mod_bid(GEN bid, GEN a)
2463 : {
2464 7 : GEN xZ = gcoeff(bid_get_ideal(bid),1,1);
2465 7 : GEN b = Rg_to_Fp(a, xZ);
2466 7 : if (gsigne(a) < 0) b = subii(b, xZ);
2467 7 : return b;
2468 : }
2469 : /* Return decomposition of a on the CRT generators blocks attached to the
2470 : * S->sprk and sarch; sgn = sign(a, S->arch), NULL if unknown */
2471 : static GEN
2472 246091 : zlog(GEN nf, GEN a, GEN sgn, zlog_S *S)
2473 : {
2474 : long k, l;
2475 : GEN y;
2476 246091 : a = nf_to_scalar_or_basis(nf, a);
2477 246091 : switch(typ(a))
2478 : {
2479 63749 : case t_INT: break;
2480 7 : case t_FRAC: a = Q_mod_bid(S->bid, a); break;
2481 : default: /* case t_COL: */
2482 : {
2483 : GEN den;
2484 182335 : check_nfelt(a, &den);
2485 182335 : if (den)
2486 : {
2487 46407 : a = Q_muli_to_int(a, den);
2488 46407 : a = mkmat2(mkcol2(a, den), mkcol2(gen_1, gen_m1));
2489 46407 : return famat_zlog(nf, a, sgn, S);
2490 : }
2491 : }
2492 : }
2493 199684 : if (sgn)
2494 34559 : sgn = (lg(sgn) == 1)? NULL: leafcopy(sgn);
2495 : else
2496 165125 : sgn = (lg(S->archp) == 1)? NULL: nfsign_arch(nf, a, S->archp);
2497 199684 : l = lg(S->sprk);
2498 199684 : y = cgetg(sgn? l+1: l, t_COL);
2499 442552 : for (k = 1; k < l; k++)
2500 : {
2501 242875 : GEN sprk = gel(S->sprk,k);
2502 242875 : gel(y,k) = zlog_pr(nf, a, sprk);
2503 : }
2504 199677 : if (sgn) gel(y,l) = Flc_to_ZC(sgn);
2505 199677 : return y;
2506 : }
2507 :
2508 : /* true nf */
2509 : GEN
2510 8344 : pr_basis_perm(GEN nf, GEN pr)
2511 : {
2512 8344 : long f = pr_get_f(pr);
2513 : GEN perm;
2514 8344 : if (f == nf_get_degree(nf)) return identity_perm(f);
2515 6944 : perm = cgetg(f+1, t_VECSMALL);
2516 6944 : perm[1] = 1;
2517 6944 : if (f > 1)
2518 : {
2519 399 : GEN H = pr_hnf(nf,pr);
2520 : long i, k;
2521 1463 : for (i = k = 2; k <= f; i++)
2522 : {
2523 1064 : if (is_pm1(gcoeff(H,i,i))) continue;
2524 840 : perm[k++] = i;
2525 : }
2526 : }
2527 6944 : return perm;
2528 : }
2529 :
2530 : /* \sum U[i]*y[i], U[i] ZM, y[i] ZC. We allow lg(y) > lg(U). */
2531 : static GEN
2532 312885 : ZMV_ZCV_mul(GEN U, GEN y)
2533 : {
2534 312885 : long i, l = lg(U);
2535 312885 : GEN z = NULL;
2536 312885 : if (l == 1) return cgetg(1,t_COL);
2537 862524 : for (i = 1; i < l; i++)
2538 : {
2539 549639 : GEN u = ZM_ZC_mul(gel(U,i), gel(y,i));
2540 549639 : z = z? ZC_add(z, u): u;
2541 : }
2542 312885 : return z;
2543 : }
2544 : /* A * (U[1], ..., U[d] */
2545 : static GEN
2546 518 : ZM_ZMV_mul(GEN A, GEN U)
2547 : {
2548 : long i, l;
2549 518 : GEN V = cgetg_copy(U,&l);
2550 518 : for (i = 1; i < l; i++) gel(V,i) = ZM_mul(A,gel(U,i));
2551 518 : return V;
2552 : }
2553 :
2554 : /* Log on bid.gen of generators of P_{1,I pr^{e-1}} / P_{1,I pr^e} (I,pr) = 1,
2555 : * defined implicitly via CRT. 'ind' is the index of pr in modulus
2556 : * factorization */
2557 : GEN
2558 50820 : log_gen_pr(zlog_S *S, long ind, GEN nf, long e)
2559 : {
2560 50820 : GEN A, sprk = gel(S->sprk,ind), Uind = gel(S->U, ind);
2561 :
2562 50820 : if (e == 1) retmkmat( gel(Uind,1) );
2563 : else
2564 : {
2565 18998 : GEN G, pr = sprk_get_pr(sprk);
2566 : long i, l;
2567 18998 : if (e == 2)
2568 : {
2569 10906 : GEN A, g, L, L2; sprk_get_L2(sprk,&A,&g,&L2); L = gel(L2,1);
2570 10906 : G = gel(L,2); l = lg(G);
2571 : }
2572 : else
2573 : {
2574 8092 : GEN perm = pr_basis_perm(nf,pr), PI = nfpow_u(nf, pr_get_gen(pr), e-1);
2575 8092 : l = lg(perm);
2576 8092 : G = cgetg(l, t_VEC);
2577 8092 : if (typ(PI) == t_INT)
2578 : { /* zk_ei_mul doesn't allow t_INT */
2579 1393 : long N = nf_get_degree(nf);
2580 1393 : gel(G,1) = addiu(PI,1);
2581 2261 : for (i = 2; i < l; i++)
2582 : {
2583 868 : GEN z = col_ei(N, 1);
2584 868 : gel(G,i) = z; gel(z, perm[i]) = PI;
2585 : }
2586 : }
2587 : else
2588 : {
2589 6699 : gel(G,1) = nfadd(nf, gen_1, PI);
2590 6909 : for (i = 2; i < l; i++)
2591 210 : gel(G,i) = nfadd(nf, gen_1, zk_ei_mul(nf, PI, perm[i]));
2592 : }
2593 : }
2594 18998 : A = cgetg(l, t_MAT);
2595 40761 : for (i = 1; i < l; i++)
2596 21763 : gel(A,i) = ZM_ZC_mul(Uind, zlog_pr(nf, gel(G,i), sprk));
2597 18998 : return A;
2598 : }
2599 : }
2600 : /* Log on bid.gen of generator of P_{1,f} / P_{1,f v[index]}
2601 : * v = vector of r1 real places */
2602 : GEN
2603 9975 : log_gen_arch(zlog_S *S, long index)
2604 : {
2605 9975 : GEN U = gel(S->U, lg(S->U)-1);
2606 9975 : return gel(U, index);
2607 : }
2608 :
2609 : /* compute bid.clgp: [h,cyc] or [h,cyc,gen] */
2610 : static GEN
2611 27139 : bid_grp(GEN nf, GEN U, GEN cyc, GEN g, GEN F, GEN sarch)
2612 : {
2613 27139 : GEN G, h = ZV_prod(cyc);
2614 : long c;
2615 27139 : if (!U) return mkvec2(h,cyc);
2616 26887 : c = lg(U);
2617 26887 : G = cgetg(c,t_VEC);
2618 26887 : if (c > 1)
2619 : {
2620 22456 : GEN U0, Uoo, EX = gel(cyc,1); /* exponent of bid */
2621 22456 : long i, hU = nbrows(U), nba = lg(sarch_get_cyc(sarch))-1; /* #f_oo */
2622 22456 : if (!nba) { U0 = U; Uoo = NULL; }
2623 11704 : else if (nba == hU) { U0 = NULL; Uoo = U; }
2624 : else
2625 : {
2626 9471 : U0 = rowslice(U, 1, hU-nba);
2627 9471 : Uoo = rowslice(U, hU-nba+1, hU);
2628 : }
2629 64386 : for (i = 1; i < c; i++)
2630 : {
2631 41930 : GEN t = gen_1;
2632 41930 : if (U0) t = famat_to_nf_modideal_coprime(nf, g, gel(U0,i), F, EX);
2633 41930 : if (Uoo) t = set_sign_mod_divisor(nf, ZV_to_Flv(gel(Uoo,i),2), t, sarch);
2634 41930 : gel(G,i) = t;
2635 : }
2636 : }
2637 26887 : return mkvec3(h, cyc, G);
2638 : }
2639 :
2640 : /* remove prime ideals of norm 2 with exponent 1 from factorization */
2641 : static GEN
2642 26824 : famat_strip2(GEN fa)
2643 : {
2644 26824 : GEN P = gel(fa,1), E = gel(fa,2), Q, F;
2645 26824 : long l = lg(P), i, j;
2646 26824 : Q = cgetg(l, t_COL);
2647 26824 : F = cgetg(l, t_COL);
2648 56490 : for (i = j = 1; i < l; i++)
2649 : {
2650 29666 : GEN pr = gel(P,i), e = gel(E,i);
2651 29666 : if (!absequaliu(pr_get_p(pr), 2) || itou(e) != 1 || pr_get_f(pr) != 1)
2652 : {
2653 25634 : gel(Q,j) = pr;
2654 25634 : gel(F,j) = e; j++;
2655 : }
2656 : }
2657 26824 : setlg(Q,j);
2658 26824 : setlg(F,j); return mkmat2(Q,F);
2659 : }
2660 :
2661 : /* Compute [[ideal,arch], [h,[cyc],[gen]], idealfact, [liste], U]
2662 : flag may include nf_GEN | nf_INIT */
2663 : static GEN
2664 26845 : Idealstar_i(GEN nf, GEN ideal, long flag)
2665 : {
2666 : long i, k, nbp, R1;
2667 26845 : GEN y, cyc, U, u1 = NULL, fa, fa2, sprk, x, arch, archp, E, P, sarch, gen;
2668 :
2669 26845 : nf = checknf(nf);
2670 26845 : R1 = nf_get_r1(nf);
2671 26845 : if (typ(ideal) == t_VEC && lg(ideal) == 3)
2672 : {
2673 12957 : arch = gel(ideal,2);
2674 12957 : ideal= gel(ideal,1);
2675 12957 : switch(typ(arch))
2676 : {
2677 : case t_VEC:
2678 12922 : if (lg(arch) != R1+1)
2679 0 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
2680 12922 : archp = vec01_to_indices(arch);
2681 12922 : break;
2682 : case t_VECSMALL:
2683 35 : archp = arch;
2684 35 : k = lg(archp)-1;
2685 35 : if (k && archp[k] > R1)
2686 7 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
2687 28 : arch = indices_to_vec01(archp, R1);
2688 28 : break;
2689 : default:
2690 0 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
2691 0 : return NULL;
2692 : }
2693 12950 : }
2694 : else
2695 : {
2696 13888 : arch = zerovec(R1);
2697 13888 : archp = cgetg(1, t_VECSMALL);
2698 : }
2699 26838 : if (is_nf_factor(ideal))
2700 : {
2701 350 : fa = ideal;
2702 350 : x = idealfactorback(nf, gel(fa,1), gel(fa,2), 0);
2703 : }
2704 : else
2705 : {
2706 26488 : fa = idealfactor(nf, ideal);
2707 26481 : x = ideal;
2708 : }
2709 26831 : if (typ(x) != t_MAT) x = idealhnf_shallow(nf, x);
2710 26831 : if (lg(x) == 1) pari_err_DOMAIN("Idealstar", "ideal","=",gen_0,x);
2711 26831 : if (typ(gcoeff(x,1,1)) != t_INT)
2712 7 : pari_err_DOMAIN("Idealstar","denominator(ideal)", "!=",gen_1,x);
2713 26824 : sarch = nfarchstar(nf, x, archp);
2714 26824 : fa2 = famat_strip2(fa);
2715 26824 : P = gel(fa2,1);
2716 26824 : E = gel(fa2,2);
2717 26824 : nbp = lg(P)-1;
2718 26824 : sprk = cgetg(nbp+1,t_VEC);
2719 26824 : if (nbp)
2720 : {
2721 20139 : GEN t = (nbp==1)? NULL: x;
2722 20139 : cyc = cgetg(nbp+2,t_VEC);
2723 20139 : gen = cgetg(nbp+1,t_VEC);
2724 45773 : for (i = 1; i <= nbp; i++)
2725 : {
2726 25634 : GEN L = sprkinit(nf, gel(P,i), gel(E,i), t);
2727 25634 : gel(sprk,i) = L;
2728 25634 : gel(cyc,i) = sprk_get_cyc(L);
2729 : /* true gens are congruent to those mod x AND positive at archp */
2730 25634 : gel(gen,i) = sprk_get_gen(L);
2731 : }
2732 20139 : gel(cyc,i) = sarch_get_cyc(sarch);
2733 20139 : cyc = shallowconcat1(cyc);
2734 20139 : gen = shallowconcat1(gen);
2735 20139 : cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
2736 : }
2737 : else
2738 : {
2739 6685 : cyc = sarch_get_cyc(sarch);
2740 6685 : gen = cgetg(1,t_VEC);
2741 6685 : U = matid(lg(cyc)-1);
2742 6685 : if (flag & nf_GEN) u1 = U;
2743 : }
2744 26824 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
2745 26824 : if (!(flag & nf_INIT)) return y;
2746 26026 : U = split_U(U, sprk);
2747 26026 : return mkvec5(mkvec2(x, arch), y, mkvec2(fa,fa2), mkvec2(sprk, sarch), U);
2748 : }
2749 : GEN
2750 26572 : Idealstar(GEN nf, GEN ideal, long flag)
2751 : {
2752 26572 : pari_sp av = avma;
2753 26572 : if (!nf) nf = nfinit(pol_x(0), DEFAULTPREC);
2754 26572 : return gerepilecopy(av, Idealstar_i(nf, ideal, flag));
2755 : }
2756 : GEN
2757 273 : Idealstarprk(GEN nf, GEN pr, long k, long flag)
2758 : {
2759 273 : pari_sp av = avma;
2760 273 : GEN z = Idealstar_i(nf, mkmat2(mkcol(pr),mkcols(k)), flag);
2761 273 : return gerepilecopy(av, z);
2762 : }
2763 :
2764 : /* FIXME: obsolete */
2765 : GEN
2766 0 : zidealstarinitgen(GEN nf, GEN ideal)
2767 0 : { return Idealstar(nf,ideal, nf_INIT|nf_GEN); }
2768 : GEN
2769 0 : zidealstarinit(GEN nf, GEN ideal)
2770 0 : { return Idealstar(nf,ideal, nf_INIT); }
2771 : GEN
2772 0 : zidealstar(GEN nf, GEN ideal)
2773 0 : { return Idealstar(nf,ideal, nf_GEN); }
2774 :
2775 : GEN
2776 63 : idealstar0(GEN nf, GEN ideal,long flag)
2777 : {
2778 63 : switch(flag)
2779 : {
2780 0 : case 0: return Idealstar(nf,ideal, nf_GEN);
2781 49 : case 1: return Idealstar(nf,ideal, nf_INIT);
2782 14 : case 2: return Idealstar(nf,ideal, nf_INIT|nf_GEN);
2783 0 : default: pari_err_FLAG("idealstar");
2784 : }
2785 : return NULL; /* LCOV_EXCL_LINE */
2786 : }
2787 :
2788 : void
2789 182335 : check_nfelt(GEN x, GEN *den)
2790 : {
2791 182335 : long l = lg(x), i;
2792 182335 : GEN t, d = NULL;
2793 182335 : if (typ(x) != t_COL) pari_err_TYPE("check_nfelt", x);
2794 664668 : for (i=1; i<l; i++)
2795 : {
2796 482333 : t = gel(x,i);
2797 482333 : switch (typ(t))
2798 : {
2799 387239 : case t_INT: break;
2800 : case t_FRAC:
2801 95094 : if (!d) d = gel(t,2); else d = lcmii(d, gel(t,2));
2802 95094 : break;
2803 0 : default: pari_err_TYPE("check_nfelt", x);
2804 : }
2805 : }
2806 182335 : *den = d;
2807 182335 : }
2808 :
2809 : GEN
2810 1205829 : vecmodii(GEN x, GEN y)
2811 1205829 : { pari_APPLY_same(modii(gel(x,i), gel(y,i))) }
2812 :
2813 : GEN
2814 95004 : vecmoduu(GEN x, GEN y)
2815 95004 : { pari_APPLY_ulong(uel(x,i) % uel(y,i)) }
2816 :
2817 : static GEN
2818 314481 : ideallog_i(GEN nf, GEN x, GEN sgn, zlog_S *S)
2819 : {
2820 314481 : pari_sp av = avma;
2821 : GEN y, cyc;
2822 314481 : if (!S->hU) return cgetg(1, t_COL);
2823 312899 : cyc = bid_get_cyc(S->bid);
2824 312899 : if (typ(x) == t_MAT)
2825 : {
2826 66808 : if (lg(x) == 1) return zerocol(lg(cyc)-1);
2827 66801 : y = famat_zlog(nf, x, sgn, S);
2828 : }
2829 : else
2830 246091 : y = zlog(nf, x, sgn, S);
2831 312885 : y = ZMV_ZCV_mul(S->U, y);
2832 312885 : return gerepileupto(av, vecmodii(y, cyc));
2833 : }
2834 :
2835 : /* Given x (not necessarily integral), and bid as output by zidealstarinit,
2836 : * compute the vector of components on the generators bid[2].
2837 : * Assume (x,bid) = 1 and sgn is either NULL or nfsign_arch(x, bid) */
2838 : GEN
2839 301321 : ideallog_sgn(GEN nf, GEN x, GEN sgn, GEN bid)
2840 : {
2841 : zlog_S S;
2842 301321 : nf = checknf(nf); checkbid(bid);
2843 301314 : init_zlog(&S, bid);
2844 301314 : if (sgn && typ(x) == t_VEC) /* vector of elements and signatures */
2845 : {
2846 21392 : long i, l = lg(x);
2847 21392 : GEN y = cgetg(l, t_MAT);
2848 21392 : for (i = 1; i < l; i++) gel(y,i) = ideallog_i(nf, gel(x,i), gel(sgn,i), &S);
2849 21392 : return y;
2850 : }
2851 279922 : return ideallog_i(nf, x, sgn, &S);
2852 : }
2853 : GEN
2854 286600 : ideallog(GEN nf, GEN x, GEN bid)
2855 : {
2856 286600 : if (!nf) return Zideallog(bid, x);
2857 279929 : return ideallog_sgn(nf, x, NULL, bid);
2858 : }
2859 :
2860 : /*************************************************************************/
2861 : /** **/
2862 : /** JOIN BID STRUCTURES, IDEAL LISTS **/
2863 : /** **/
2864 : /*************************************************************************/
2865 : /* bid1, bid2: for coprime modules m1 and m2 (without arch. part).
2866 : * Output: bid for m1 m2 */
2867 : static GEN
2868 476 : join_bid(GEN nf, GEN bid1, GEN bid2)
2869 : {
2870 476 : pari_sp av = avma;
2871 : long nbgen, l1,l2;
2872 : GEN I1,I2, G1,G2, sprk1,sprk2, cyc1,cyc2, sarch;
2873 476 : GEN sprk, fa,fa2, U, cyc, y, u1 = NULL, x, gen;
2874 :
2875 476 : I1 = bid_get_ideal(bid1);
2876 476 : I2 = bid_get_ideal(bid2);
2877 476 : if (gequal1(gcoeff(I1,1,1))) return bid2; /* frequent trivial case */
2878 259 : G1 = bid_get_grp(bid1);
2879 259 : G2 = bid_get_grp(bid2);
2880 259 : x = idealmul(nf, I1,I2);
2881 259 : fa = famat_mul_shallow(bid_get_fact(bid1), bid_get_fact(bid2));
2882 259 : fa2= famat_mul_shallow(bid_get_fact2(bid1), bid_get_fact2(bid2));
2883 259 : sprk1 = bid_get_sprk(bid1);
2884 259 : sprk2 = bid_get_sprk(bid2);
2885 259 : sprk = shallowconcat(sprk1, sprk2);
2886 :
2887 259 : cyc1 = abgrp_get_cyc(G1); l1 = lg(cyc1);
2888 259 : cyc2 = abgrp_get_cyc(G2); l2 = lg(cyc2);
2889 259 : gen = (lg(G1)>3 && lg(G2)>3)? gen_1: NULL;
2890 259 : nbgen = l1+l2-2;
2891 259 : cyc = ZV_snf_group(shallowconcat(cyc1,cyc2), &U, gen? &u1: NULL);
2892 259 : if (nbgen)
2893 : {
2894 259 : GEN U1 = bid_get_U(bid1), U2 = bid_get_U(bid2);
2895 259 : U1 = l1==1? const_vec(lg(sprk1), cgetg(1,t_MAT))
2896 259 : : ZM_ZMV_mul(vecslice(U, 1, l1-1), U1);
2897 259 : U2 = l2==1? const_vec(lg(sprk2), cgetg(1,t_MAT))
2898 259 : : ZM_ZMV_mul(vecslice(U, l1, nbgen), U2);
2899 259 : U = shallowconcat(U1, U2);
2900 : }
2901 : else
2902 0 : U = const_vec(lg(sprk), cgetg(1,t_MAT));
2903 :
2904 259 : if (gen)
2905 : {
2906 259 : GEN uv = zkchinese1init2(nf, I2, I1, x);
2907 518 : gen = shallowconcat(zkVchinese1(gel(uv,1), abgrp_get_gen(G1)),
2908 259 : zkVchinese1(gel(uv,2), abgrp_get_gen(G2)));
2909 : }
2910 259 : sarch = bid_get_sarch(bid1); /* trivial */
2911 259 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
2912 259 : x = mkvec2(x, bid_get_arch(bid1));
2913 259 : y = mkvec5(x, y, mkvec2(fa, fa2), mkvec2(sprk, sarch), U);
2914 259 : return gerepilecopy(av,y);
2915 : }
2916 :
2917 : typedef struct _ideal_data {
2918 : GEN nf, emb, L, pr, prL, sgnU, archp;
2919 : } ideal_data;
2920 :
2921 : /* z <- ( z | f(v[i])_{i=1..#v} ) */
2922 : static void
2923 86065 : concat_join(GEN *pz, GEN v, GEN (*f)(ideal_data*,GEN), ideal_data *data)
2924 : {
2925 86065 : long i, nz, lv = lg(v);
2926 : GEN z, Z;
2927 172130 : if (lv == 1) return;
2928 38143 : z = *pz; nz = lg(z)-1;
2929 38143 : *pz = Z = cgetg(lv + nz, typ(z));
2930 38143 : for (i = 1; i <=nz; i++) gel(Z,i) = gel(z,i);
2931 38143 : Z += nz;
2932 38143 : for (i = 1; i < lv; i++) gel(Z,i) = f(data, gel(v,i));
2933 : }
2934 : static GEN
2935 476 : join_idealinit(ideal_data *D, GEN x)
2936 476 : { return join_bid(D->nf, x, D->prL); }
2937 : static GEN
2938 47698 : join_ideal(ideal_data *D, GEN x)
2939 47698 : { return idealmulpowprime(D->nf, x, D->pr, D->L); }
2940 : static GEN
2941 455 : join_unit(ideal_data *D, GEN x)
2942 : {
2943 455 : GEN bid = join_idealinit(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
2944 455 : if (lg(u) != 1) v = shallowconcat(u, v);
2945 455 : return mkvec2(bid, v);
2946 : }
2947 :
2948 : /* flag & nf_GEN : generators, otherwise no
2949 : * flag &2 : units, otherwise no
2950 : * flag &4 : ideals in HNF, otherwise bid
2951 : * flag &8 : omit ideals which cannot be conductors (pr^1 with Npr=2) */
2952 : static GEN
2953 3192 : Ideallist(GEN bnf, ulong bound, long flag)
2954 : {
2955 3192 : const long cond = flag & 8;
2956 3192 : const long do_units = flag & 2, big_id = !(flag & 4);
2957 3192 : const long istar_flag = (flag & nf_GEN) | nf_INIT;
2958 3192 : pari_sp av, av0 = avma;
2959 : long i, j;
2960 3192 : GEN nf, z, p, fa, id, BOUND, U, empty = cgetg(1,t_VEC);
2961 : forprime_t S;
2962 : ideal_data ID;
2963 3192 : GEN (*join_z)(ideal_data*, GEN) =
2964 : do_units? &join_unit
2965 3192 : : (big_id? &join_idealinit: &join_ideal);
2966 :
2967 3192 : nf = checknf(bnf);
2968 3192 : if ((long)bound <= 0) return empty;
2969 3192 : id = matid(nf_get_degree(nf));
2970 3192 : if (big_id) id = Idealstar(nf,id, istar_flag);
2971 :
2972 : /* z[i] will contain all "objects" of norm i. Depending on flag, this means
2973 : * an ideal, a bid, or a couple [bid, log(units)]. Such objects are stored
2974 : * in vectors, computed one primary component at a time; join_z
2975 : * reconstructs the global object */
2976 3192 : BOUND = utoipos(bound);
2977 3192 : z = cgetg(bound+1,t_VEC);
2978 3192 : if (do_units) {
2979 14 : U = bnf_build_units(bnf);
2980 14 : gel(z,1) = mkvec( mkvec2(id, cgetg(1,t_VEC)) );
2981 : } else {
2982 3178 : U = NULL; /* -Wall */
2983 3178 : gel(z,1) = mkvec(id);
2984 : }
2985 3192 : for (i=2; i<=(long)bound; i++) gel(z,i) = empty;
2986 3192 : ID.nf = nf;
2987 :
2988 3192 : p = cgetipos(3);
2989 3192 : u_forprime_init(&S, 2, bound);
2990 3192 : av = avma;
2991 19600 : while ((p[2] = u_forprime_next(&S)))
2992 : {
2993 13216 : if (DEBUGLEVEL>1) { err_printf("%ld ",p[2]); err_flush(); }
2994 13216 : fa = idealprimedec_limit_norm(nf, p, BOUND);
2995 26859 : for (j=1; j<lg(fa); j++)
2996 : {
2997 13643 : GEN pr = gel(fa,j), z2 = leafcopy(z);
2998 13643 : ulong Q, q = upr_norm(pr);
2999 13643 : long l = (cond && q == 2)? 2: 1;
3000 :
3001 13643 : ID.pr = ID.prL = pr;
3002 33775 : for (Q = q; Q <= bound; l++, Q *= q) /* add pr^l */
3003 : {
3004 : ulong iQ;
3005 20132 : ID.L = utoipos(l);
3006 20132 : if (big_id) {
3007 217 : ID.prL = Idealstarprk(nf, pr, l, istar_flag);
3008 217 : if (do_units)
3009 : {
3010 196 : GEN sprk = bid_get_sprk(ID.prL);
3011 392 : ID.emb = lg(sprk) == 1? cgetg(1,t_VEC)
3012 196 : : vzlog_pr(nf, U, gel(sprk,1));
3013 : }
3014 : }
3015 106197 : for (iQ = Q,i = 1; iQ <= bound; iQ += Q,i++)
3016 86065 : concat_join(&gel(z,iQ), gel(z2,i), join_z, &ID);
3017 : }
3018 : }
3019 13216 : if (gc_needed(av,1))
3020 : {
3021 0 : if(DEBUGMEM>1) pari_warn(warnmem,"Ideallist");
3022 0 : z = gerepilecopy(av, z);
3023 : }
3024 : }
3025 3192 : return gerepilecopy(av0, z);
3026 : }
3027 : GEN
3028 350 : ideallist0(GEN bnf,long bound, long flag) {
3029 350 : if (flag<0 || flag>15) pari_err_FLAG("ideallist");
3030 350 : return Ideallist(bnf,bound,flag);
3031 : }
3032 : GEN
3033 2842 : ideallist(GEN bnf,long bound) { return Ideallist(bnf,bound,4); }
3034 :
3035 : /* bid = for module m (without arch. part), arch = archimedean part.
3036 : * Output: bid for [m,arch] */
3037 : static GEN
3038 56 : join_bid_arch(GEN nf, GEN bid, GEN archp)
3039 : {
3040 56 : pari_sp av = avma;
3041 : GEN G, U;
3042 56 : GEN sprk, cyc, y, u1 = NULL, x, sarch, gen;
3043 :
3044 56 : checkbid(bid);
3045 56 : G = bid_get_grp(bid);
3046 56 : x = bid_get_ideal(bid);
3047 56 : sarch = nfarchstar(nf, bid_get_ideal(bid), archp);
3048 56 : sprk = bid_get_sprk(bid);
3049 :
3050 56 : gen = (lg(G)>3)? gel(G,3): NULL;
3051 56 : cyc = diagonal_shallow(shallowconcat(gel(G,2), sarch_get_cyc(sarch)));
3052 56 : cyc = ZM_snf_group(cyc, &U, gen? &u1: NULL);
3053 56 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
3054 56 : U = split_U(U, sprk);
3055 56 : y = mkvec5(mkvec2(x, archp), y, gel(bid,3), mkvec2(sprk, sarch), U);
3056 56 : return gerepilecopy(av,y);
3057 : }
3058 : static GEN
3059 56 : join_arch(ideal_data *D, GEN x) {
3060 56 : return join_bid_arch(D->nf, x, D->archp);
3061 : }
3062 : static GEN
3063 28 : join_archunit(ideal_data *D, GEN x) {
3064 28 : GEN bid = join_arch(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
3065 28 : if (lg(u) != 1) v = shallowconcat(u, v);
3066 28 : return mkvec2(bid, v);
3067 : }
3068 :
3069 : /* L from ideallist, add archimedean part */
3070 : GEN
3071 14 : ideallistarch(GEN bnf, GEN L, GEN arch)
3072 : {
3073 : pari_sp av;
3074 14 : long i, j, l = lg(L), lz;
3075 : GEN v, z, V;
3076 : ideal_data ID;
3077 : GEN (*join_z)(ideal_data*, GEN);
3078 :
3079 14 : if (typ(L) != t_VEC) pari_err_TYPE("ideallistarch",L);
3080 14 : if (l == 1) return cgetg(1,t_VEC);
3081 14 : z = gel(L,1);
3082 14 : if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
3083 14 : z = gel(z,1); /* either a bid or [bid,U] */
3084 14 : ID.nf = checknf(bnf);
3085 14 : ID.archp = vec01_to_indices(arch);
3086 14 : if (lg(z) == 3) { /* the latter: do units */
3087 7 : if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
3088 7 : ID.emb = zm_to_ZM( rowpermute(nfsign_units(bnf,NULL,1), ID.archp) );
3089 7 : join_z = &join_archunit;
3090 : } else
3091 7 : join_z = &join_arch;
3092 14 : av = avma; V = cgetg(l, t_VEC);
3093 70 : for (i = 1; i < l; i++)
3094 : {
3095 56 : z = gel(L,i); lz = lg(z);
3096 56 : gel(V,i) = v = cgetg(lz,t_VEC);
3097 56 : for (j=1; j<lz; j++) gel(v,j) = join_z(&ID, gel(z,j));
3098 : }
3099 14 : return gerepilecopy(av,V);
3100 : }
|
