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 9658193 : get_tab(GEN nf, long *N)
30 : {
31 9658193 : GEN tab = (typ(nf) == t_MAT)? nf: gel(nf,9);
32 9658193 : *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 388575046 : _mulii(GEN x, GEN y) {
38 1004748700 : return is_pm1(x)? (signe(x) < 0)? negi(y): y
39 616173654 : : mulii(x, y);
40 : }
41 :
42 : GEN
43 16919 : tablemul_ei_ej(GEN M, long i, long j)
44 : {
45 : long N;
46 16919 : GEN tab = get_tab(M, &N);
47 16919 : 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 10269 : tablemul_ei(GEN M, GEN x, long i)
55 : {
56 : long j, k, N;
57 : GEN v, tab;
58 :
59 10269 : if (i==1) return gcopy(x);
60 10269 : tab = get_tab(M, &N);
61 10269 : if (typ(x) != t_COL) { v = zerocol(N); gel(v,i) = gcopy(x); return v; }
62 10269 : tab += (i-1)*N; v = cgetg(N+1,t_COL);
63 : /* wi . x = [ sum_j tab[k,j] x[j] ]_k */
64 70007 : for (k=1; k<=N; k++)
65 : {
66 59738 : pari_sp av = avma;
67 59738 : GEN s = gen_0;
68 417312 : for (j=1; j<=N; j++)
69 : {
70 357574 : GEN c = gcoeff(tab,k,j);
71 357574 : if (!gequal0(c)) s = gadd(s, gmul(c, gel(x,j)));
72 : }
73 59738 : gel(v,k) = gerepileupto(av,s);
74 : }
75 10269 : return v;
76 : }
77 : /* as tablemul_ei, assume x a ZV of correct length */
78 : GEN
79 7370414 : zk_ei_mul(GEN nf, GEN x, long i)
80 : {
81 : long j, k, N;
82 : GEN v, tab;
83 :
84 7370414 : if (i==1) return ZC_copy(x);
85 7370400 : tab = get_tab(nf, &N); tab += (i-1)*N;
86 7370400 : v = cgetg(N+1,t_COL);
87 57265656 : for (k=1; k<=N; k++)
88 : {
89 49895256 : pari_sp av = avma;
90 49895256 : GEN s = gen_0;
91 638698584 : for (j=1; j<=N; j++)
92 : {
93 588803328 : GEN c = gcoeff(tab,k,j);
94 588803328 : if (signe(c)) s = addii(s, _mulii(c, gel(x,j)));
95 : }
96 49895256 : gel(v,k) = gerepileuptoint(av, s);
97 : }
98 7370400 : return v;
99 : }
100 :
101 : /* table of multiplication by wi in R[w1,..., wN] */
102 : GEN
103 2275 : ei_multable(GEN TAB, long i)
104 : {
105 : long k,N;
106 2275 : GEN m, tab = get_tab(TAB, &N);
107 2275 : tab += (i-1)*N;
108 2275 : m = cgetg(N+1,t_MAT);
109 2275 : for (k=1; k<=N; k++) gel(m,k) = gel(tab,k);
110 2275 : return m;
111 : }
112 :
113 : GEN
114 2989460 : zk_multable(GEN nf, GEN x)
115 : {
116 2989460 : long i, l = lg(x);
117 2989460 : GEN mul = cgetg(l,t_MAT);
118 2989460 : gel(mul,1) = x; /* assume w_1 = 1 */
119 2989460 : for (i=2; i<l; i++) gel(mul,i) = zk_ei_mul(nf,x,i);
120 2989460 : return mul;
121 : }
122 : GEN
123 1617 : multable(GEN M, GEN x)
124 : {
125 : long i, N;
126 : GEN mul;
127 1617 : 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 1593320 : zk_scalar_or_multable(GEN nf, GEN x)
140 : {
141 1593320 : long tx = typ(x);
142 1593320 : if (tx == t_MAT || tx == t_INT) return x;
143 1568894 : x = nf_to_scalar_or_basis(nf, x);
144 1568894 : return (typ(x) == t_COL)? zk_multable(nf, x): x;
145 : }
146 :
147 : GEN
148 23562 : nftrace(GEN nf, GEN x)
149 : {
150 23562 : pari_sp av = avma;
151 23562 : nf = checknf(nf);
152 23562 : x = nf_to_scalar_or_basis(nf, x);
153 70665 : x = (typ(x) == t_COL)? RgV_dotproduct(x, gel(nf_get_Tr(nf),1))
154 47103 : : gmulgs(x, nf_get_degree(nf));
155 23562 : return gerepileupto(av, x);
156 : }
157 : GEN
158 567 : rnfelttrace(GEN rnf, GEN x)
159 : {
160 567 : pari_sp av = avma;
161 567 : checkrnf(rnf);
162 567 : x = rnfeltabstorel(rnf, x);
163 1344 : x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gtrace(x))
164 868 : : gmulgs(x, rnf_get_degree(rnf));
165 476 : 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 24602 : nfnorm(GEN nf, GEN x)
181 : {
182 24602 : pari_sp av = avma;
183 24602 : nf = checknf(nf);
184 24602 : if (typ(x) == t_MAT) return famat_norm(nf, x);
185 24595 : x = nf_to_scalar_or_alg(nf, x);
186 69872 : x = (typ(x) == t_POL)? RgXQ_norm(x, nf_get_pol(nf))
187 45277 : : gpowgs(x, nf_get_degree(nf));
188 24595 : 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 15574636 : nfadd(GEN nf, GEN x, GEN y)
205 : {
206 15574636 : pari_sp av = avma;
207 : GEN z;
208 :
209 15574636 : nf = checknf(nf);
210 15574636 : x = nf_to_scalar_or_basis(nf, x);
211 15574636 : y = nf_to_scalar_or_basis(nf, y);
212 15574636 : if (typ(x) != t_COL)
213 12674572 : { z = (typ(y) == t_COL)? RgC_Rg_add(y, x): gadd(x,y); }
214 : else
215 2900064 : { z = (typ(y) == t_COL)? RgC_add(x, y): RgC_Rg_add(x, y); }
216 15574636 : return gerepileupto(av, z);
217 : }
218 : /* x - y in nf */
219 : GEN
220 1191645 : nfsub(GEN nf, GEN x, GEN y)
221 : {
222 1191645 : pari_sp av = avma;
223 : GEN z;
224 :
225 1191645 : nf = checknf(nf);
226 1191645 : x = nf_to_scalar_or_basis(nf, x);
227 1191645 : y = nf_to_scalar_or_basis(nf, y);
228 1191645 : if (typ(x) != t_COL)
229 883393 : { z = (typ(y) == t_COL)? Rg_RgC_sub(x,y): gsub(x,y); }
230 : else
231 308252 : { z = (typ(y) == t_COL)? RgC_sub(x,y): RgC_Rg_sub(x,y); }
232 1191645 : return gerepileupto(av, z);
233 : }
234 :
235 : /* product of x and y in nf */
236 : GEN
237 20574744 : nfmul(GEN nf, GEN x, GEN y)
238 : {
239 : GEN z;
240 20574744 : pari_sp av = avma;
241 :
242 20574744 : if (x == y) return nfsqr(nf,x);
243 :
244 17721166 : nf = checknf(nf);
245 17721166 : x = nf_to_scalar_or_basis(nf, x);
246 17721166 : y = nf_to_scalar_or_basis(nf, y);
247 17721166 : if (typ(x) != t_COL)
248 : {
249 13725642 : if (isintzero(x)) return gen_0;
250 9848083 : z = (typ(y) == t_COL)? RgC_Rg_mul(y, x): gmul(x,y); }
251 : else
252 : {
253 3995524 : if (typ(y) != t_COL)
254 : {
255 2839074 : if (isintzero(y)) return gen_0;
256 638106 : z = RgC_Rg_mul(x, y);
257 : }
258 : else
259 : {
260 : GEN dx, dy;
261 1156450 : x = Q_remove_denom(x, &dx);
262 1156450 : y = Q_remove_denom(y, &dy);
263 1156450 : z = nfmuli(nf,x,y);
264 1156450 : dx = mul_denom(dx,dy);
265 1156450 : if (dx) z = RgC_Rg_div(z, dx);
266 : }
267 : }
268 11642639 : return gerepileupto(av, z);
269 : }
270 : /* square of x in nf */
271 : GEN
272 4745779 : nfsqr(GEN nf, GEN x)
273 : {
274 4745779 : pari_sp av = avma;
275 : GEN z;
276 :
277 4745779 : nf = checknf(nf);
278 4745779 : x = nf_to_scalar_or_basis(nf, x);
279 4745779 : if (typ(x) != t_COL) z = gsqr(x);
280 : else
281 : {
282 : GEN dx;
283 113781 : x = Q_remove_denom(x, &dx);
284 113781 : z = nfsqri(nf,x);
285 113781 : if (dx) z = RgC_Rg_div(z, sqri(dx));
286 : }
287 4745779 : return gerepileupto(av, z);
288 : }
289 :
290 : /* x a ZC, v a t_COL of ZC/Z */
291 : GEN
292 129421 : zkC_multable_mul(GEN v, GEN x)
293 : {
294 129421 : long i, l = lg(v);
295 129421 : GEN y = cgetg(l, t_COL);
296 456688 : for (i = 1; i < l; i++)
297 : {
298 327267 : GEN c = gel(v,i);
299 327267 : if (typ(c)!=t_COL) {
300 0 : if (!isintzero(c)) c = ZC_Z_mul(gel(x,1), c);
301 : } else {
302 327267 : c = ZM_ZC_mul(x,c);
303 327267 : if (ZV_isscalar(c)) c = gel(c,1);
304 : }
305 327267 : gel(y,i) = c;
306 : }
307 129421 : return y;
308 : }
309 :
310 : GEN
311 49896 : nfC_multable_mul(GEN v, GEN x)
312 : {
313 49896 : long i, l = lg(v);
314 49896 : GEN y = cgetg(l, t_COL);
315 311752 : for (i = 1; i < l; i++)
316 : {
317 261856 : GEN c = gel(v,i);
318 261856 : if (typ(c)!=t_COL) {
319 215405 : if (!isintzero(c)) c = RgC_Rg_mul(gel(x,1), c);
320 : } else {
321 46451 : c = RgM_RgC_mul(x,c);
322 46451 : if (QV_isscalar(c)) c = gel(c,1);
323 : }
324 261856 : gel(y,i) = c;
325 : }
326 49896 : return y;
327 : }
328 :
329 : GEN
330 164507 : nfC_nf_mul(GEN nf, GEN v, GEN x)
331 : {
332 : long tx;
333 : GEN y;
334 :
335 164507 : x = nf_to_scalar_or_basis(nf, x);
336 164507 : tx = typ(x);
337 164507 : if (tx != t_COL)
338 : {
339 : long l, i;
340 121212 : if (tx == t_INT)
341 : {
342 113267 : long s = signe(x);
343 113267 : if (!s) return zerocol(lg(v)-1);
344 107097 : if (is_pm1(x)) return s > 0? leafcopy(v): RgC_neg(v);
345 : }
346 39029 : l = lg(v); y = cgetg(l, t_COL);
347 274333 : for (i=1; i < l; i++)
348 : {
349 235304 : GEN c = gel(v,i);
350 235304 : if (typ(c) != t_COL) c = gmul(c, x); else c = RgC_Rg_mul(c, x);
351 235304 : gel(y,i) = c;
352 : }
353 39029 : return y;
354 : }
355 : else
356 : {
357 : GEN dx;
358 43295 : x = zk_multable(nf, Q_remove_denom(x,&dx));
359 43295 : y = nfC_multable_mul(v, x);
360 43295 : return dx? RgC_Rg_div(y, dx): y;
361 : }
362 : }
363 : static GEN
364 7392 : mulbytab(GEN M, GEN c)
365 7392 : { return typ(c) == t_COL? RgM_RgC_mul(M,c): RgC_Rg_mul(gel(M,1), c); }
366 : GEN
367 1617 : tablemulvec(GEN M, GEN x, GEN v)
368 : {
369 : long l, i;
370 : GEN y;
371 :
372 1617 : 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 1617 : x = multable(M, x); /* multiplication table by x */
378 1617 : y = cgetg_copy(v, &l);
379 1617 : if (typ(v) == t_POL)
380 : {
381 1617 : y[1] = v[1];
382 1617 : for (i=2; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
383 1617 : 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 1617 : return y;
390 : }
391 :
392 : GEN
393 378891 : zkmultable_capZ(GEN mx) { return Q_denom(zkmultable_inv(mx)); }
394 :
395 : GEN
396 418663 : zkmultable_inv(GEN mx)
397 418663 : { return ZM_gauss(mx, col_ei(lg(mx)-1,1)); }
398 :
399 : /* nf a true nf, x a ZC */
400 : GEN
401 39772 : zk_inv(GEN nf, GEN x) { return zkmultable_inv(zk_multable(nf,x)); }
402 :
403 : /* inverse of x in nf */
404 : GEN
405 63553 : nfinv(GEN nf, GEN x)
406 : {
407 63553 : pari_sp av = avma;
408 : GEN z;
409 :
410 63553 : nf = checknf(nf);
411 63553 : x = nf_to_scalar_or_basis(nf, x);
412 63553 : if (typ(x) == t_COL)
413 : {
414 : GEN d;
415 24451 : x = Q_remove_denom(x, &d);
416 24451 : z = zk_inv(nf, x);
417 24451 : if (d) z = RgC_Rg_mul(z, d);
418 : }
419 : else
420 39102 : z = ginv(x);
421 63553 : return gerepileupto(av, z);
422 : }
423 :
424 : /* quotient of x and y in nf */
425 : GEN
426 18326 : nfdiv(GEN nf, GEN x, GEN y)
427 : {
428 18326 : pari_sp av = avma;
429 : GEN z;
430 :
431 18326 : nf = checknf(nf);
432 18326 : y = nf_to_scalar_or_basis(nf, y);
433 18326 : if (typ(y) != t_COL)
434 : {
435 10352 : x = nf_to_scalar_or_basis(nf, x);
436 10352 : z = (typ(x) == t_COL)? RgC_Rg_div(x, y): gdiv(x,y);
437 : }
438 : else
439 : {
440 : GEN d;
441 7974 : y = Q_remove_denom(y, &d);
442 7974 : z = nfmul(nf, x, zk_inv(nf,y));
443 7974 : if (d) z = RgC_Rg_mul(z, d);
444 : }
445 18326 : 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 1585310 : nfmuli(GEN nf, GEN x, GEN y)
452 : {
453 : long i, j, k, N;
454 1585310 : GEN s, v, TAB = get_tab(nf, &N);
455 :
456 1585310 : if (typ(x) == t_INT) return (typ(y) == t_COL)? ZC_Z_mul(y, x): mulii(x,y);
457 1489133 : if (typ(y) == t_INT) return ZC_Z_mul(x, y);
458 : /* both x and y are ZV */
459 1454376 : v = cgetg(N+1,t_COL);
460 5946980 : for (k=1; k<=N; k++)
461 : {
462 4492604 : pari_sp av = avma;
463 4492604 : GEN TABi = TAB;
464 4492604 : if (k == 1)
465 1454376 : s = mulii(gel(x,1),gel(y,1));
466 : else
467 6076456 : s = addii(mulii(gel(x,1),gel(y,k)),
468 6076456 : mulii(gel(x,k),gel(y,1)));
469 23144746 : for (i=2; i<=N; i++)
470 : {
471 18652142 : GEN t, xi = gel(x,i);
472 18652142 : TABi += N;
473 18652142 : if (!signe(xi)) continue;
474 :
475 12909075 : t = NULL;
476 133353841 : for (j=2; j<=N; j++)
477 : {
478 120444766 : GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
479 120444766 : if (!signe(c)) continue;
480 54236324 : p1 = _mulii(c, gel(y,j));
481 54236324 : t = t? addii(t, p1): p1;
482 : }
483 12909075 : if (t) s = addii(s, mulii(xi, t));
484 : }
485 4492604 : gel(v,k) = gerepileuptoint(av,s);
486 : }
487 1454376 : return v;
488 : }
489 : /* square of INTEGER (t_INT or ZC) x in nf */
490 : GEN
491 673020 : nfsqri(GEN nf, GEN x)
492 : {
493 : long i, j, k, N;
494 673020 : GEN s, v, TAB = get_tab(nf, &N);
495 :
496 673020 : if (typ(x) == t_INT) return sqri(x);
497 673020 : v = cgetg(N+1,t_COL);
498 5318109 : for (k=1; k<=N; k++)
499 : {
500 4645089 : pari_sp av = avma;
501 4645089 : GEN TABi = TAB;
502 4645089 : if (k == 1)
503 673020 : s = sqri(gel(x,1));
504 : else
505 3972069 : s = shifti(mulii(gel(x,1),gel(x,k)), 1);
506 56173513 : for (i=2; i<=N; i++)
507 : {
508 51528424 : GEN p1, c, t, xi = gel(x,i);
509 51528424 : TABi += N;
510 51528424 : if (!signe(xi)) continue;
511 :
512 17392361 : c = gcoeff(TABi, k, i);
513 17392361 : t = signe(c)? _mulii(c,xi): NULL;
514 250040874 : for (j=i+1; j<=N; j++)
515 : {
516 232648513 : c = gcoeff(TABi, k, j);
517 232648513 : if (!signe(c)) continue;
518 121512851 : p1 = _mulii(c, shifti(gel(x,j),1));
519 121512851 : t = t? addii(t, p1): p1;
520 : }
521 17392361 : if (t) s = addii(s, mulii(xi, t));
522 : }
523 4645089 : gel(v,k) = gerepileuptoint(av,s);
524 : }
525 673020 : 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 39900 : tablesqr(GEN TAB, GEN x)
569 : {
570 : long i, j, k, N;
571 : GEN s, v;
572 :
573 39900 : if (typ(x) != t_COL) return gsqr(x);
574 39900 : N = lg(x)-1;
575 39900 : v = cgetg(N+1,t_COL);
576 :
577 278250 : for (k=1; k<=N; k++)
578 : {
579 238350 : pari_sp av = avma;
580 238350 : GEN TABi = TAB;
581 238350 : if (k == 1)
582 39900 : s = gsqr(gel(x,1));
583 : else
584 198450 : s = gmul2n(gmul(gel(x,1),gel(x,k)), 1);
585 1453368 : for (i=2; i<=N; i++)
586 : {
587 1215018 : GEN p1, c, t, xi = gel(x,i);
588 1215018 : TABi += N;
589 1215018 : if (gequal0(xi)) continue;
590 :
591 320733 : c = gcoeff(TABi, k, i);
592 320733 : t = !gequal0(c)? gmul(c,xi): NULL;
593 1240911 : for (j=i+1; j<=N; j++)
594 : {
595 920178 : c = gcoeff(TABi, k, j);
596 920178 : if (gequal0(c)) continue;
597 482475 : p1 = gmul(gmul2n(c,1), gel(x,j));
598 482475 : t = t? gadd(t, p1): p1;
599 : }
600 320733 : if (t) s = gadd(s, gmul(xi, t));
601 : }
602 238350 : gel(v,k) = gerepileupto(av,s);
603 : }
604 39900 : return v;
605 : }
606 :
607 : static GEN
608 28429 : _mul(void *data, GEN x, GEN y) { return nfmuli((GEN)data,x,y); }
609 : static GEN
610 105405 : _sqr(void *data, GEN x) { return nfsqri((GEN)data,x); }
611 :
612 : /* Compute z^n in nf, left-shift binary powering */
613 : GEN
614 101946 : nfpow(GEN nf, GEN z, GEN n)
615 : {
616 101946 : pari_sp av = avma;
617 : long s;
618 : GEN x, cx;
619 :
620 101946 : if (typ(n)!=t_INT) pari_err_TYPE("nfpow",n);
621 101946 : nf = checknf(nf);
622 101946 : s = signe(n); if (!s) return gen_1;
623 101946 : x = nf_to_scalar_or_basis(nf, z);
624 101946 : if (typ(x) != t_COL) return powgi(x,n);
625 101470 : if (s < 0)
626 : { /* simplified nfinv */
627 : GEN d;
628 3430 : x = Q_remove_denom(x, &d);
629 3430 : x = zk_inv(nf, x);
630 3430 : x = primitive_part(x, &cx);
631 3430 : cx = mul_content(cx, d);
632 3430 : n = absi(n);
633 : }
634 : else
635 98040 : x = primitive_part(x, &cx);
636 101470 : x = gen_pow(x, n, (void*)nf, _sqr, _mul);
637 101470 : if (cx) x = gmul(x, powgi(cx, n));
638 101470 : return av==avma? gcopy(x): gerepileupto(av,x);
639 : }
640 : /* Compute z^n in nf, left-shift binary powering */
641 : GEN
642 45703 : nfpow_u(GEN nf, GEN z, ulong n)
643 : {
644 45703 : pari_sp av = avma;
645 : GEN x, cx;
646 :
647 45703 : nf = checknf(nf);
648 45703 : if (!n) return gen_1;
649 45703 : x = nf_to_scalar_or_basis(nf, z);
650 45703 : if (typ(x) != t_COL) return gpowgs(x,n);
651 17465 : x = primitive_part(x, &cx);
652 17465 : x = gen_powu(x, n, (void*)nf, _sqr, _mul);
653 17465 : if (cx) x = gmul(x, powgi(cx, utoipos(n)));
654 17465 : return av==avma? gcopy(x): gerepileupto(av,x);
655 : }
656 :
657 : static GEN
658 2393755 : _nf_red(void *E, GEN x) { (void)E; return x; }
659 :
660 : static GEN
661 9315495 : _nf_add(void *E, GEN x, GEN y) { return nfadd((GEN)E,x,y); }
662 :
663 : static GEN
664 569884 : _nf_neg(void *E, GEN x) { (void)E; return gneg(x); }
665 :
666 : static GEN
667 11276363 : _nf_mul(void *E, GEN x, GEN y) { return nfmul((GEN)E,x,y); }
668 :
669 : static GEN
670 41426 : _nf_inv(void *E, GEN x) { return nfinv((GEN)E,x); }
671 :
672 : static GEN
673 8232 : _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 177177 : const struct bb_field *get_nf_field(void **E, GEN nf)
679 177177 : { *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 8218 : nfM_inv(GEN nf, GEN M)
690 : {
691 : void *E;
692 8218 : const struct bb_field *S = get_nf_field(&E, nf);
693 8218 : return gen_Gauss(M, matid(lg(M)-1), E, S);
694 : }
695 : GEN
696 8022 : nfM_mul(GEN nf, GEN A, GEN B)
697 : {
698 : void *E;
699 8022 : const struct bb_field *S = get_nf_field(&E, nf);
700 8022 : return gen_matmul(A, B, E, S);
701 : }
702 : GEN
703 160923 : nfM_nfC_mul(GEN nf, GEN A, GEN B)
704 : {
705 : void *E;
706 160923 : const struct bb_field *S = get_nf_field(&E, nf);
707 160923 : return gen_matcolmul(A, B, E, S);
708 : }
709 :
710 : /* valuation of integral x (ZV), with resp. to prime ideal pr */
711 : long
712 5074052 : ZC_nfvalrem(GEN x, GEN pr, GEN *newx)
713 : {
714 : long i, v, l;
715 5074052 : GEN r, y, p = pr_get_p(pr), mul = pr_get_tau(pr);
716 :
717 : /* p inert */
718 5074052 : if (typ(mul) == t_INT) return newx? ZV_pvalrem(x, p, newx):ZV_pval(x, p);
719 5063069 : y = cgetg_copy(x, &l); /* will hold the new x */
720 5063069 : x = leafcopy(x);
721 7628404 : for(v=0;; v++)
722 : {
723 26845951 : for (i=1; i<l; i++)
724 : { /* is (x.b)[i] divisible by p ? */
725 24280616 : gel(y,i) = dvmdii(ZMrow_ZC_mul(mul,x,i),p,&r);
726 24280616 : if (r != gen_0) { if (newx) *newx = x; return v; }
727 : }
728 2565335 : swap(x, y);
729 2565335 : }
730 : }
731 : long
732 4836895 : ZC_nfval(GEN x, GEN P)
733 4836895 : { return ZC_nfvalrem(x, P, NULL); }
734 :
735 : /* v_P(x) != 0, x a ZV. Simpler version of ZC_nfvalrem */
736 : int
737 243740 : ZC_prdvd(GEN x, GEN P)
738 : {
739 243740 : pari_sp av = avma;
740 : long i, l;
741 243740 : GEN p = pr_get_p(P), mul = pr_get_tau(P);
742 243740 : if (typ(mul) == t_INT) return ZV_Z_dvd(x, p);
743 243523 : l = lg(x);
744 973434 : for (i=1; i<l; i++)
745 869663 : if (remii(ZMrow_ZC_mul(mul,x,i), p) != gen_0) { avma = av; return 0; }
746 103771 : 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 1301440 : nfval(GEN nf, GEN x, GEN pr)
763 : {
764 1301440 : pari_sp av = avma;
765 : long w, e;
766 : GEN cx, p;
767 :
768 1301440 : if (gequal0(x)) return LONG_MAX;
769 1299956 : nf = checknf(nf);
770 1299956 : checkprid(pr);
771 1299956 : p = pr_get_p(pr);
772 1299956 : e = pr_get_e(pr);
773 1299956 : x = nf_to_scalar_or_basis(nf, x);
774 1299956 : if (typ(x) != t_COL) return e*Q_pval(x,p);
775 106393 : x = Q_primitive_part(x, &cx);
776 106393 : w = ZC_nfval(x,pr);
777 106393 : if (cx) w += e*Q_pval(cx,p);
778 106393 : 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 19985 : b = pr_get_tau(pr);
790 19985 : if (typ(b) == t_INT) return gen_1;
791 1281 : e = pr_get_e(pr);
792 1281 : z = gel(b,1);
793 1281 : if (e != 1) z = gdiv(nfpow_u(nf, z, e), powiu(pr_get_p(pr),e-1));
794 1281 : if (v < 0) { v = -v; z = nfinv(nf, z); }
795 1281 : if (v != 1) z = nfpow_u(nf, z, v);
796 1281 : 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 52864 : w = Q_pvalrem(x,p, py);
814 52864 : if (!w) { *py = gerepilecopy(av, x); return 0; }
815 18970 : *py = gerepileupto(av, gmul(powp(nf, pr, w), *py));
816 18970 : return e*w;
817 : }
818 11886 : x = Q_primitive_part(x, &cx);
819 11886 : w = ZC_nfvalrem(x,pr, py);
820 11886 : if (cx)
821 : {
822 1134 : long v = Q_pvalrem(cx,p, &t);
823 1134 : *py = nfmul(nf, *py, gmul(powp(nf,pr,v), t));
824 1134 : *py = gerepileupto(av, *py);
825 1134 : w += e*v;
826 : }
827 : else
828 10752 : *py = gerepilecopy(av, *py);
829 11886 : 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 78407 : coltoalg(GEN nf, GEN x)
841 : {
842 78407 : return mkpolmod( nf_to_scalar_or_alg(nf, x), nf_get_pol(nf) );
843 : }
844 :
845 : GEN
846 87199 : basistoalg(GEN nf, GEN x)
847 : {
848 : GEN z, T;
849 :
850 87199 : nf = checknf(nf);
851 87199 : switch(typ(x))
852 : {
853 : case t_COL: {
854 72289 : pari_sp av = avma;
855 72289 : 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 959 : T = nf_get_pol(nf);
864 959 : if (varn(T) != varn(x)) pari_err_VAR("basistoalg",x,T);
865 959 : z = cgetg(3,t_POLMOD);
866 959 : gel(z,1) = ZX_copy(T);
867 959 : gel(z,2) = RgX_rem(x, T); return z;
868 : case t_INT:
869 : case t_FRAC:
870 13468 : T = nf_get_pol(nf);
871 13468 : z = cgetg(3,t_POLMOD);
872 13468 : gel(z,1) = ZX_copy(T);
873 13468 : 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 1440723 : pol_to_scalar_or_basis(GEN nf, GEN x)
883 : {
884 1440723 : GEN T = nf_get_pol(nf);
885 1440723 : long l = lg(x);
886 1440723 : if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_basis", x,T);
887 1440660 : if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
888 1440660 : if (l == 2) return gen_0;
889 844659 : if (l == 3)
890 : {
891 201754 : x = gel(x,2);
892 201754 : if (!is_rational_t(typ(x))) pari_err_TYPE("nf_to_scalar_or_basis",x);
893 201754 : return x;
894 : }
895 642905 : return poltobasis(nf,x);
896 : }
897 : /* Assume nf is a genuine nf. */
898 : GEN
899 81498592 : nf_to_scalar_or_basis(GEN nf, GEN x)
900 : {
901 81498592 : switch(typ(x))
902 : {
903 : case t_INT: case t_FRAC:
904 62823054 : return x;
905 : case t_POLMOD:
906 181643 : x = checknfelt_mod(nf,x,"nf_to_scalar_or_basis");
907 181573 : switch(typ(x))
908 : {
909 34454 : case t_INT: case t_FRAC: return x;
910 147119 : case t_POL: return pol_to_scalar_or_basis(nf,x);
911 : }
912 0 : break;
913 1293604 : case t_POL: return pol_to_scalar_or_basis(nf,x);
914 : case t_COL:
915 17200291 : if (lg(x)-1 != nf_get_degree(nf)) break;
916 17200228 : 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 5188 : RgX_to_nfX(GEN nf, GEN x)
926 : {
927 : long i, l;
928 5188 : GEN y = cgetg_copy(x, &l); y[1] = x[1];
929 5188 : for (i=2; i<l; i++) gel(y,i) = nf_to_scalar_or_basis(nf, gel(x,i));
930 5188 : return y;
931 : }
932 :
933 : /* Assume nf is a genuine nf. */
934 : GEN
935 182710 : nf_to_scalar_or_alg(GEN nf, GEN x)
936 : {
937 182710 : switch(typ(x))
938 : {
939 : case t_INT: case t_FRAC:
940 14881 : 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 15891 : GEN T = nf_get_pol(nf);
948 15891 : long l = lg(x);
949 15891 : if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_alg", x,T);
950 15891 : if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
951 15891 : if (l == 2) return gen_0;
952 15891 : if (l == 3) return gel(x,2);
953 15667 : return x;
954 : }
955 : case t_COL:
956 : {
957 : GEN dx;
958 151882 : if (lg(x)-1 != nf_get_degree(nf)) break;
959 303764 : if (QV_isscalar(x)) return gel(x,1);
960 123630 : x = Q_remove_denom(x, &dx);
961 123630 : x = RgV_RgC_mul(nf_get_zkprimpart(nf), x);
962 123630 : dx = mul_denom(dx, nf_get_zkden(nf));
963 123630 : 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 2318089 : ZM_ZX_mul(GEN A, GEN x)
984 : {
985 2318089 : long i, l = lg(x)-1;
986 : GEN z;
987 2318089 : if (l == 1) return zerocol(nbrows(A));
988 2317084 : z = ZC_Z_mul(gel(A,1), gel(x,2));
989 9201372 : for (i = 2; i < l ; i++)
990 6884292 : if (signe(gel(x,i+1))) z = ZC_add(z, ZC_Z_mul(gel(A,i), gel(x,i+1)));
991 2317080 : return z;
992 : }
993 : /* x a t_POL, nf a genuine nf. No garbage collecting. No check. */
994 : GEN
995 2160812 : poltobasis(GEN nf, GEN x)
996 : {
997 2160812 : GEN d, T = nf_get_pol(nf);
998 2160812 : if (varn(x) != varn(T)) pari_err_VAR( "poltobasis", x,T);
999 2160756 : if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
1000 2160756 : x = Q_remove_denom(x, &d);
1001 2160756 : if (!RgX_is_ZX(x)) pari_err_TYPE("poltobasis",x);
1002 2160735 : x = ZM_ZX_mul(nf_get_invzk(nf), x);
1003 2160735 : if (d) x = RgC_Rg_div(x, d);
1004 2160735 : return x;
1005 : }
1006 :
1007 : GEN
1008 154466 : algtobasis(GEN nf, GEN x)
1009 : {
1010 : pari_sp av;
1011 :
1012 154466 : nf = checknf(nf);
1013 154466 : switch(typ(x))
1014 : {
1015 : case t_POLMOD:
1016 38563 : 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 38556 : x = gel(x,2);
1019 38556 : switch(typ(x))
1020 : {
1021 : case t_INT:
1022 4242 : case t_FRAC: return scalarcol(x, nf_get_degree(nf));
1023 : case t_POL:
1024 34314 : av = avma;
1025 34314 : return gerepileupto(av,poltobasis(nf,x));
1026 : }
1027 0 : break;
1028 :
1029 : case t_POL:
1030 55812 : av = avma;
1031 55812 : return gerepileupto(av,poltobasis(nf,x));
1032 :
1033 : case t_COL:
1034 14913 : if (lg(x)-1 != nf_get_degree(nf)) pari_err_DIM("nfalgtobasis");
1035 14913 : return gcopy(x);
1036 :
1037 : case t_INT:
1038 45178 : 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 36001 : rnfbasistoalg(GEN rnf,GEN x)
1046 : {
1047 36001 : const char *f = "rnfbasistoalg";
1048 : long lx, i;
1049 36001 : pari_sp av = avma;
1050 : GEN z, nf, relpol, T;
1051 :
1052 36001 : checkrnf(rnf);
1053 36001 : nf = rnf_get_nf(rnf);
1054 36001 : T = nf_get_pol(nf);
1055 36001 : relpol = QXQX_to_mod_shallow(rnf_get_pol(rnf), T);
1056 36001 : 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 23884 : x = polmod_nffix(f, rnf, x, 0);
1071 23674 : if (typ(x) != t_POL) break;
1072 9639 : retmkpolmod(RgX_copy(x), RgX_copy(relpol));
1073 : case t_POL:
1074 847 : if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = gmodulo(x,T); break; }
1075 623 : if (varn(x) == varn(relpol))
1076 : {
1077 574 : x = RgX_nffix(f,nf_get_pol(nf),x,0);
1078 574 : return gmodulo(x, relpol);
1079 : }
1080 49 : pari_err_VAR(f, x,relpol);
1081 : }
1082 24682 : retmkpolmod(scalarpol(x, varn(relpol)), RgX_copy(relpol));
1083 : }
1084 :
1085 : GEN
1086 1232 : matbasistoalg(GEN nf,GEN x)
1087 : {
1088 : long i, j, li, lx;
1089 1232 : GEN z = cgetg_copy(x, &lx);
1090 :
1091 1232 : if (lx == 1) return z;
1092 1225 : switch(typ(x))
1093 : {
1094 : case t_VEC: case t_COL:
1095 28 : for (i=1; i<lx; i++) gel(z,i) = basistoalg(nf, gel(x,i));
1096 28 : return z;
1097 1197 : case t_MAT: break;
1098 0 : default: pari_err_TYPE("matbasistoalg",x);
1099 : }
1100 1197 : li = lgcols(x);
1101 4655 : for (j=1; j<lx; j++)
1102 : {
1103 3458 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1104 3458 : gel(z,j) = c;
1105 3458 : for (i=1; i<li; i++) gel(c,i) = basistoalg(nf, gel(xj,i));
1106 : }
1107 1197 : return z;
1108 : }
1109 :
1110 : GEN
1111 2604 : matalgtobasis(GEN nf,GEN x)
1112 : {
1113 : long i, j, li, lx;
1114 2604 : GEN z = cgetg_copy(x, &lx);
1115 :
1116 2604 : if (lx == 1) return z;
1117 2548 : switch(typ(x))
1118 : {
1119 : case t_VEC: case t_COL:
1120 2541 : for (i=1; i<lx; i++) gel(z,i) = algtobasis(nf, gel(x,i));
1121 2541 : 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 8337 : RgM_to_nfM(GEN nf,GEN x)
1136 : {
1137 : long i, j, li, lx;
1138 8337 : GEN z = cgetg_copy(x, &lx);
1139 :
1140 8337 : if (lx == 1) return z;
1141 8337 : li = lgcols(x);
1142 64344 : for (j=1; j<lx; j++)
1143 : {
1144 56007 : GEN c = cgetg(li,t_COL), xj = gel(x,j);
1145 56007 : gel(z,j) = c;
1146 56007 : for (i=1; i<li; i++) gel(c,i) = nf_to_scalar_or_basis(nf, gel(xj,i));
1147 : }
1148 8337 : return z;
1149 : }
1150 : GEN
1151 75873 : RgC_to_nfC(GEN nf,GEN x)
1152 : {
1153 75873 : long i, lx = lg(x);
1154 75873 : GEN z = cgetg(lx, t_COL);
1155 75873 : for (i=1; i<lx; i++) gel(z,i) = nf_to_scalar_or_basis(nf, gel(x,i));
1156 75873 : return z;
1157 : }
1158 :
1159 : /* x a t_POLMOD, supposedly in rnf = K[z]/(T), K = Q[y]/(Tnf) */
1160 : GEN
1161 61362 : polmod_nffix(const char *f, GEN rnf, GEN x, int lift)
1162 61362 : { return polmod_nffix2(f, rnf_get_nfpol(rnf), rnf_get_pol(rnf), x,lift); }
1163 : GEN
1164 61453 : polmod_nffix2(const char *f, GEN T, GEN relpol, GEN x, int lift)
1165 : {
1166 61453 : if (RgX_equal_var(gel(x,1),relpol))
1167 : {
1168 56007 : x = gel(x,2);
1169 56007 : if (typ(x) == t_POL && varn(x) == varn(relpol))
1170 : {
1171 28917 : x = RgX_nffix(f, T, x, lift);
1172 28917 : switch(lg(x))
1173 : {
1174 294 : case 2: return gen_0;
1175 3374 : case 3: return gel(x,2);
1176 : }
1177 25249 : return x;
1178 : }
1179 : }
1180 32536 : return Rg_nffix(f, T, x, lift);
1181 : }
1182 : GEN
1183 1176 : rnfalgtobasis(GEN rnf,GEN x)
1184 : {
1185 1176 : const char *f = "rnfalgtobasis";
1186 1176 : pari_sp av = avma;
1187 : GEN T, relpol;
1188 :
1189 1176 : checkrnf(rnf);
1190 1176 : relpol = rnf_get_pol(rnf);
1191 1176 : T = rnf_get_nfpol(rnf);
1192 1176 : switch(typ(x))
1193 : {
1194 : case t_COL:
1195 49 : if (lg(x)-1 != rnf_get_degree(rnf)) pari_err_DIM(f);
1196 28 : x = RgV_nffix(f, T, x, 0);
1197 21 : return gerepilecopy(av, x);
1198 :
1199 : case t_POLMOD:
1200 1043 : x = polmod_nffix(f, rnf, x, 0);
1201 1001 : if (typ(x) != t_POL) break;
1202 707 : return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
1203 : case t_POL:
1204 56 : if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = mkpolmod(x,T); break; }
1205 35 : x = RgX_nffix(f, T, x, 0);
1206 28 : if (degpol(x) >= degpol(relpol)) x = RgX_rem(x,relpol);
1207 28 : return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
1208 : }
1209 336 : return gerepileupto(av, scalarcol(x, rnf_get_degree(rnf)));
1210 : }
1211 :
1212 : /* Given a and b in nf, gives an algebraic integer y in nf such that a-b.y
1213 : * is "small" */
1214 : GEN
1215 259 : nfdiveuc(GEN nf, GEN a, GEN b)
1216 : {
1217 259 : pari_sp av = avma;
1218 259 : a = nfdiv(nf,a,b);
1219 259 : return gerepileupto(av, ground(a));
1220 : }
1221 :
1222 : /* Given a and b in nf, gives a "small" algebraic integer r in nf
1223 : * of the form a-b.y */
1224 : GEN
1225 259 : nfmod(GEN nf, GEN a, GEN b)
1226 : {
1227 259 : pari_sp av = avma;
1228 259 : GEN p1 = gneg_i(nfmul(nf,b,ground(nfdiv(nf,a,b))));
1229 259 : return gerepileupto(av, nfadd(nf,a,p1));
1230 : }
1231 :
1232 : /* Given a and b in nf, gives a two-component vector [y,r] in nf such
1233 : * that r=a-b.y is "small". */
1234 : GEN
1235 259 : nfdivrem(GEN nf, GEN a, GEN b)
1236 : {
1237 259 : pari_sp av = avma;
1238 259 : GEN p1,z, y = ground(nfdiv(nf,a,b));
1239 :
1240 259 : p1 = gneg_i(nfmul(nf,b,y));
1241 259 : z = cgetg(3,t_VEC);
1242 259 : gel(z,1) = gcopy(y);
1243 259 : gel(z,2) = nfadd(nf,a,p1); return gerepileupto(av, z);
1244 : }
1245 :
1246 : /*************************************************************************/
1247 : /** **/
1248 : /** REAL EMBEDDINGS **/
1249 : /** **/
1250 : /*************************************************************************/
1251 : static GEN
1252 49161 : sarch_get_cyc(GEN sarch) { return gel(sarch,1); }
1253 : static GEN
1254 264559 : sarch_get_archp(GEN sarch) { return gel(sarch,2); }
1255 : static GEN
1256 55248 : sarch_get_MI(GEN sarch) { return gel(sarch,3); }
1257 : static GEN
1258 55248 : sarch_get_lambda(GEN sarch) { return gel(sarch,4); }
1259 : static GEN
1260 55248 : sarch_get_F(GEN sarch) { return gel(sarch,5); }
1261 :
1262 : /* true nf, return number of positive roots of char_x */
1263 : static long
1264 1434 : num_positive(GEN nf, GEN x)
1265 : {
1266 1434 : GEN T = nf_get_pol(nf);
1267 1434 : GEN charx = ZXQ_charpoly(nf_to_scalar_or_alg(nf,x), T, 0);
1268 : long np;
1269 1434 : charx = ZX_radical(charx);
1270 1434 : np = ZX_sturmpart(charx, mkvec2(gen_0,mkoo()));
1271 1434 : return np * (degpol(T) / degpol(charx));
1272 : }
1273 :
1274 : /* x a QC: return sigma_k(x) where 1 <= k <= r1+r2; correct but inefficient
1275 : * if x in Q. M = nf_get_M(nf) */
1276 : static GEN
1277 10031 : nfembed_i(GEN M, GEN x, long k)
1278 : {
1279 10031 : long i, l = lg(M);
1280 10031 : GEN z = gel(x,1);
1281 10031 : for (i = 2; i < l; i++) z = gadd(z, gmul(gcoeff(M,k,i), gel(x,i)));
1282 10031 : return z;
1283 : }
1284 : GEN
1285 1610 : nfembed(GEN nf, GEN x, long k)
1286 : {
1287 1610 : pari_sp av = avma;
1288 1610 : nf = checknf(nf);
1289 1610 : x = nf_to_scalar_or_basis(nf,x);
1290 1610 : if (typ(x) != t_COL) return gerepilecopy(av, x);
1291 0 : return gerepileupto(av, nfembed_i(nf_get_M(nf),x,k));
1292 : }
1293 :
1294 : /* x a ZC */
1295 : static GEN
1296 393771 : zk_embed(GEN M, GEN x, long k)
1297 : {
1298 393771 : long i, l = lg(x);
1299 393771 : GEN z = gel(x,1); /* times M[k,1], which is 1 */
1300 393771 : for (i = 2; i < l; i++) z = mpadd(z, mpmul(gcoeff(M,k,i), gel(x,i)));
1301 393771 : return z;
1302 : }
1303 :
1304 : /* Given floating point approximation z of sigma_k(x), decide its sign
1305 : * [0/+, 1/- and -1 for FAIL] */
1306 : static long
1307 382480 : eval_sign_embed(GEN z)
1308 : { /* dubious, fail */
1309 382480 : if (typ(z) == t_REAL && realprec(z) <= LOWDEFAULTPREC) return -1;
1310 381540 : return (signe(z) < 1)? 1: 0;
1311 : }
1312 : /* return v such that (-1)^v = sign(sigma_k(x)), x primitive ZC */
1313 : static long
1314 310625 : eval_sign(GEN M, GEN x, long k)
1315 310625 : { return eval_sign_embed( zk_embed(M, x, k) ); }
1316 :
1317 : /* check that signs[i..#signs] == s; signs = NULL encodes "totally positive" */
1318 : static int
1319 0 : oksigns(long l, GEN signs, long i, long s)
1320 : {
1321 0 : if (!signs) return s == 0;
1322 0 : for (; i < l; i++)
1323 0 : if (signs[i] != s) return 0;
1324 0 : return 1;
1325 : }
1326 : /* check that signs[i] = s and signs[i+1..#signs] = 1-s */
1327 : static int
1328 0 : oksigns2(long l, GEN signs, long i, long s)
1329 : {
1330 0 : if (!signs) return s == 0 && i == l-1;
1331 0 : return signs[i] == s && oksigns(l, signs, i+1, 1-s);
1332 : }
1333 :
1334 : /* true nf, x a ZC (primitive for efficiency), embx its embeddings or NULL */
1335 : static int
1336 63371 : nfchecksigns_i(GEN nf, GEN x, GEN embx, GEN signs, GEN archp)
1337 : {
1338 63371 : long l = lg(archp), i;
1339 63371 : GEN M = nf_get_M(nf), sarch = NULL;
1340 63371 : long np = -1;
1341 93219 : for (i = 1; i < l; i++)
1342 : {
1343 : long s;
1344 72044 : if (embx)
1345 71855 : s = eval_sign_embed(gel(embx,i));
1346 : else
1347 189 : s = eval_sign(M, x, archp[i]);
1348 : /* 0 / + or 1 / -; -1 for FAIL */
1349 72044 : if (s < 0) /* failure */
1350 : {
1351 0 : long ni, r1 = nf_get_r1(nf);
1352 : GEN xi;
1353 0 : if (np < 0)
1354 : {
1355 0 : np = num_positive(nf, x);
1356 0 : if (np == 0) return oksigns(l, signs, i, 1);
1357 0 : if (np == r1) return oksigns(l, signs, i, 0);
1358 0 : sarch = nfarchstar(nf, NULL, identity_perm(r1));
1359 : }
1360 0 : xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
1361 0 : xi = Q_primpart(xi);
1362 0 : ni = num_positive(nf, nfmuli(nf,x,xi));
1363 0 : if (ni == 0) return oksigns2(l, signs, i, 0);
1364 0 : if (ni == r1) return oksigns2(l, signs, i, 1);
1365 0 : s = ni < np? 0: 1;
1366 : }
1367 72044 : if (s != (signs? signs[i]: 0)) return 0;
1368 : }
1369 21175 : return 1;
1370 : }
1371 : static void
1372 189 : pl_convert(GEN pl, GEN *psigns, GEN *parchp)
1373 : {
1374 189 : long i, j, l = lg(pl);
1375 189 : GEN signs = cgetg(l, t_VECSMALL);
1376 189 : GEN archp = cgetg(l, t_VECSMALL);
1377 434 : for (i = j = 1; i < l; i++)
1378 : {
1379 245 : if (!pl[i]) continue;
1380 231 : archp[j] = i;
1381 231 : signs[j] = (pl[i] < 0)? 1: 0;
1382 231 : j++;
1383 : }
1384 189 : setlg(archp, j); *parchp = archp;
1385 189 : setlg(signs, j); *psigns = signs;
1386 189 : }
1387 : /* pl : requested signs for real embeddings, 0 = no sign constraint */
1388 : int
1389 546 : nfchecksigns(GEN nf, GEN x, GEN pl)
1390 : {
1391 546 : pari_sp av = avma;
1392 : GEN signs, archp;
1393 : int res;
1394 546 : nf = checknf(nf);
1395 546 : x = nf_to_scalar_or_basis(nf,x);
1396 546 : if (typ(x) != t_COL)
1397 : {
1398 357 : long i, l = lg(pl), s = gsigne(x);
1399 791 : for (i = 1; i < l; i++)
1400 434 : if (pl[i] && pl[i] != s) { avma = av; return 0; }
1401 357 : avma = av; return 1;
1402 : }
1403 189 : pl_convert(pl, &signs, &archp);
1404 189 : res = nfchecksigns_i(nf, x, NULL, signs, archp);
1405 189 : avma = av; return res;
1406 : }
1407 :
1408 : /* signs = NULL: totally positive, else sign[i] = 0 (+) or 1 (-) */
1409 : static GEN
1410 55248 : get_C(GEN lambda, long l, GEN signs)
1411 : {
1412 : long i;
1413 : GEN C, mlambda;
1414 55248 : if (!signs) return const_vec(l-1, lambda);
1415 15754 : C = cgetg(l, t_COL); mlambda = gneg(lambda);
1416 15754 : for (i = 1; i < l; i++) gel(C,i) = signs[i]? mlambda: lambda;
1417 15754 : return C;
1418 : }
1419 : /* signs = NULL: totally positive at archp */
1420 : static GEN
1421 91676 : nfsetsigns(GEN nf, GEN signs, GEN x, GEN sarch)
1422 : {
1423 91676 : long i, l = lg(sarch_get_archp(sarch));
1424 : GEN ex;
1425 : /* Is signature already correct ? */
1426 91676 : if (typ(x) != t_COL)
1427 : {
1428 28494 : long s = gsigne(x) < 0? 1: 0;
1429 28494 : if (!signs)
1430 2674 : i = (s == 1)? 1: l;
1431 : else
1432 : {
1433 40156 : for (i = 1; i < l; i++)
1434 27402 : if (signs[i] != s) break;
1435 : }
1436 28494 : ex = (i < l)? const_col(l-1, x): NULL;
1437 : }
1438 : else
1439 : {
1440 63182 : pari_sp av = avma;
1441 63182 : GEN cex, M = nf_get_M(nf), archp = sarch_get_archp(sarch);
1442 63182 : GEN xp = Q_primitive_part(x,&cex);
1443 63182 : ex = cgetg(l,t_COL);
1444 63182 : for (i = 1; i < l; i++) gel(ex,i) = zk_embed(M,xp,archp[i]);
1445 63182 : if (nfchecksigns_i(nf, xp, ex, signs, archp)) { ex = NULL; avma = av; }
1446 42154 : else if (cex) ex = RgC_Rg_mul(ex, cex); /* put back content */
1447 : }
1448 91676 : if (ex)
1449 : { /* If no, fix it */
1450 55248 : GEN lambda = sarch_get_lambda(sarch);
1451 55248 : GEN MI = sarch_get_MI(sarch);
1452 55248 : GEN F = sarch_get_F(sarch);
1453 55248 : GEN t = RgC_sub(get_C(lambda, l, signs), ex);
1454 : long e;
1455 55248 : t = grndtoi(RgM_RgC_mul(MI,t), &e);
1456 55248 : if (lg(F) != 1) t = ZM_ZC_mul(F, t);
1457 55248 : x = typ(x) == t_COL? RgC_add(t, x): RgC_Rg_add(t, x);
1458 : }
1459 91676 : return x;
1460 : }
1461 : /* - sarch = nfarchstar(nf, F);
1462 : * - x encodes a vector of signs at arch.archp: either a t_VECSMALL
1463 : * (vector of signs as {0,1}-vector), NULL (totally positive at archp),
1464 : * or a non-zero number field element (replaced by its signature at archp);
1465 : * - y is a non-zero number field element
1466 : * Return z = y (mod F) with signs(y, archp) = signs(x) (a {0,1}-vector) */
1467 : GEN
1468 109701 : set_sign_mod_divisor(GEN nf, GEN x, GEN y, GEN sarch)
1469 : {
1470 109701 : GEN archp = sarch_get_archp(sarch);
1471 109701 : if (lg(archp) == 1) return y;
1472 90367 : nf = checknf(nf);
1473 90367 : if (x && typ(x) != t_VECSMALL) x = nfsign_arch(nf, x, archp);
1474 90367 : y = nf_to_scalar_or_basis(nf,y);
1475 90367 : return nfsetsigns(nf, x, y, sarch);
1476 : }
1477 :
1478 : static GEN
1479 14354 : setsigns_init(GEN nf, GEN archp, GEN F, GEN DATA)
1480 : {
1481 14354 : GEN lambda, Mr = rowpermute(nf_get_M(nf), archp), MI = F? RgM_mul(Mr,F): Mr;
1482 14354 : lambda = gmul2n(matrixnorm(MI,DEFAULTPREC), -1);
1483 14354 : if (lg(archp) < lg(MI))
1484 : {
1485 12719 : GEN perm = gel(indexrank(MI), 2);
1486 12719 : if (!F) F = matid(nf_get_degree(nf));
1487 12719 : MI = vecpermute(MI, perm);
1488 12719 : F = vecpermute(F, perm);
1489 : }
1490 14354 : if (!F) F = cgetg(1,t_MAT);
1491 14354 : MI = RgM_inv(MI);
1492 14354 : return mkvec5(DATA, archp, MI, lambda, F);
1493 : }
1494 : /* F non-0 integral ideal in HNF (or NULL: Z_K), compute elements in 1+F
1495 : * whose sign matrix at archp is identity; archp in 'indices' format */
1496 : GEN
1497 28074 : nfarchstar(GEN nf, GEN F, GEN archp)
1498 : {
1499 28074 : long nba = lg(archp) - 1;
1500 28074 : if (!nba) return mkvec2(cgetg(1,t_VEC), archp);
1501 13052 : if (F && equali1(gcoeff(F,1,1))) F = NULL;
1502 13052 : if (F) F = idealpseudored(F, nf_get_roundG(nf));
1503 13052 : return setsigns_init(nf, archp, F, const_vec(nba, gen_2));
1504 : }
1505 :
1506 : /*************************************************************************/
1507 : /** **/
1508 : /** IDEALCHINESE **/
1509 : /** **/
1510 : /*************************************************************************/
1511 : static int
1512 2079 : isprfact(GEN x)
1513 : {
1514 : long i, l;
1515 : GEN L, E;
1516 2079 : if (typ(x) != t_MAT || lg(x) != 3) return 0;
1517 2079 : L = gel(x,1); l = lg(L);
1518 2079 : E = gel(x,2);
1519 4942 : for(i=1; i<l; i++)
1520 : {
1521 2863 : checkprid(gel(L,i));
1522 2863 : if (typ(gel(E,i)) != t_INT) return 0;
1523 : }
1524 2079 : return 1;
1525 : }
1526 :
1527 : /* initialize projectors mod pr[i]^e[i] for idealchinese */
1528 : static GEN
1529 2079 : pr_init(GEN nf, GEN fa, GEN w, GEN dw)
1530 : {
1531 2079 : GEN U, E, F, L = gel(fa,1), E0 = gel(fa,2);
1532 2079 : long i, r = lg(L);
1533 :
1534 2079 : if (w && lg(w) != r) pari_err_TYPE("idealchinese", w);
1535 2079 : if (r == 1 && !dw) return cgetg(1,t_VEC);
1536 2072 : E = leafcopy(E0); /* do not destroy fa[2] */
1537 4935 : for (i = 1; i < r; i++)
1538 2863 : if (signe(gel(E,i)) < 0) gel(E,i) = gen_0;
1539 2072 : F = factorbackprime(nf, L, E);
1540 2072 : if (dw)
1541 : {
1542 686 : F = ZM_Z_mul(F, dw);
1543 1568 : for (i = 1; i < r; i++)
1544 : {
1545 882 : GEN pr = gel(L,i);
1546 882 : long e = itos(gel(E0,i)), v = idealval(nf, dw, pr);
1547 882 : if (e >= 0)
1548 875 : gel(E,i) = addiu(gel(E,i), v);
1549 7 : else if (v + e <= 0)
1550 0 : F = idealmulpowprime(nf, F, pr, stoi(-v)); /* coprime to pr */
1551 : else
1552 : {
1553 7 : F = idealmulpowprime(nf, F, pr, stoi(e));
1554 7 : gel(E,i) = stoi(v + e);
1555 : }
1556 : }
1557 : }
1558 2072 : U = cgetg(r, t_VEC);
1559 4935 : for (i = 1; i < r; i++)
1560 : {
1561 : GEN u;
1562 2863 : if (w && gequal0(gel(w,i))) u = gen_0; /* unused */
1563 : else
1564 : {
1565 2807 : GEN pr = gel(L,i), e = gel(E,i), t;
1566 2807 : t = idealdivpowprime(nf,F, pr, e);
1567 2807 : u = hnfmerge_get_1(t, idealpow(nf, pr, e));
1568 2807 : if (!u) pari_err_COPRIME("idealchinese", t,pr);
1569 : }
1570 2863 : gel(U,i) = u;
1571 : }
1572 2072 : F = idealpseudored(F, nf_get_roundG(nf));
1573 2072 : return mkvec2(F, U);
1574 : }
1575 :
1576 : static GEN
1577 1302 : pl_normalize(GEN nf, GEN pl)
1578 : {
1579 1302 : const char *fun = "idealchinese";
1580 1302 : if (lg(pl)-1 != nf_get_r1(nf)) pari_err_TYPE(fun,pl);
1581 1302 : switch(typ(pl))
1582 : {
1583 679 : case t_VEC: RgV_check_ZV(pl,fun); pl = ZV_to_zv(pl);
1584 : /* fall through */
1585 1302 : case t_VECSMALL: break;
1586 0 : default: pari_err_TYPE(fun,pl);
1587 : }
1588 1302 : return pl;
1589 : }
1590 :
1591 : static int
1592 4949 : is_chineseinit(GEN x)
1593 : {
1594 : GEN fa, pl;
1595 : long l;
1596 4949 : if (typ(x) != t_VEC || lg(x)!=3) return 0;
1597 3682 : fa = gel(x,1);
1598 3682 : pl = gel(x,2);
1599 3682 : if (typ(fa) != t_VEC || typ(pl) != t_VEC) return 0;
1600 1596 : l = lg(fa);
1601 1596 : if (l != 1)
1602 : {
1603 1575 : if (l != 3 || typ(gel(fa,1)) != t_MAT || typ(gel(fa,2)) != t_VEC)
1604 7 : return 0;
1605 : }
1606 1589 : l = lg(pl);
1607 1589 : if (l != 1)
1608 : {
1609 518 : if (l != 6 || typ(gel(pl,3)) != t_MAT || typ(gel(pl,1)) != t_VECSMALL
1610 518 : || typ(gel(pl,2)) != t_VECSMALL)
1611 0 : return 0;
1612 : }
1613 1589 : return 1;
1614 : }
1615 :
1616 : /* nf a true 'nf' */
1617 : static GEN
1618 2142 : chineseinit_i(GEN nf, GEN fa, GEN w, GEN dw)
1619 : {
1620 2142 : const char *fun = "idealchineseinit";
1621 2142 : GEN archp = NULL, pl = NULL;
1622 2142 : switch(typ(fa))
1623 : {
1624 : case t_VEC:
1625 1302 : if (is_chineseinit(fa))
1626 : {
1627 0 : if (dw) pari_err_DOMAIN(fun, "denom(y)", "!=", gen_1, w);
1628 0 : return fa;
1629 : }
1630 1302 : if (lg(fa) != 3) pari_err_TYPE(fun, fa);
1631 : /* of the form [x,s] */
1632 1302 : pl = pl_normalize(nf, gel(fa,2));
1633 1302 : fa = gel(fa,1);
1634 1302 : archp = vecsmall01_to_indices(pl);
1635 : /* keep pr_init, reset pl */
1636 1302 : if (is_chineseinit(fa)) { fa = gel(fa,1); break; }
1637 : /* fall through */
1638 : case t_MAT: /* factorization? */
1639 2079 : if (isprfact(fa)) { fa = pr_init(nf, fa, w, dw); break; }
1640 0 : default: pari_err_TYPE(fun,fa);
1641 : }
1642 :
1643 2142 : if (pl)
1644 : {
1645 1302 : GEN F = (lg(fa) == 1)? NULL: gel(fa,1);
1646 1302 : long i, r = lg(archp);
1647 1302 : GEN signs = cgetg(r, t_VECSMALL);
1648 1302 : for (i = 1; i < r; i++) signs[i] = (pl[archp[i]] < 0)? 1: 0;
1649 1302 : pl = setsigns_init(nf, archp, F, signs);
1650 : }
1651 : else
1652 840 : pl = cgetg(1,t_VEC);
1653 2142 : return mkvec2(fa, pl);
1654 : }
1655 :
1656 : /* Given a prime ideal factorization x, possibly with 0 or negative exponents,
1657 : * and a vector w of elements of nf, gives b such that
1658 : * v_p(b-w_p)>=v_p(x) for all prime ideals p in the ideal factorization
1659 : * and v_p(b)>=0 for all other p, using the standard proof given in GTM 138. */
1660 : GEN
1661 3668 : idealchinese(GEN nf, GEN x, GEN w)
1662 : {
1663 3668 : const char *fun = "idealchinese";
1664 3668 : pari_sp av = avma;
1665 : GEN x1, x2, s, dw, F;
1666 :
1667 3668 : nf = checknf(nf);
1668 3668 : if (!w) return gerepilecopy(av, chineseinit_i(nf,x,NULL,NULL));
1669 :
1670 2345 : if (typ(w) != t_VEC) pari_err_TYPE(fun,w);
1671 2345 : w = Q_remove_denom(matalgtobasis(nf,w), &dw);
1672 2345 : if (!is_chineseinit(x)) x = chineseinit_i(nf,x,w,dw);
1673 : /* x is a 'chineseinit' */
1674 2345 : x1 = gel(x,1); s = NULL;
1675 2345 : if (lg(x1) == 1) F = NULL;
1676 : else
1677 : {
1678 2324 : GEN U = gel(x1,2);
1679 2324 : long i, r = lg(w);
1680 2324 : F = gel(x1,1);
1681 5719 : for (i=1; i<r; i++)
1682 3395 : if (!gequal0(gel(w,i)))
1683 : {
1684 2842 : GEN t = nfmuli(nf, gel(U,i), gel(w,i));
1685 2842 : s = s? ZC_add(s,t): t;
1686 : }
1687 2324 : if (s) s = ZC_reducemodmatrix(s, F);
1688 : }
1689 2345 : if (!s) { s = zerocol(nf_get_degree(nf)); dw = NULL; }
1690 :
1691 2345 : x2 = gel(x,2);
1692 2345 : if (lg(x2) != 1) s = nfsetsigns(nf, gel(x2,1), s, x2);
1693 2345 : if (dw) s = RgC_Rg_div(s,dw);
1694 2345 : return gerepileupto(av, s);
1695 : }
1696 :
1697 : /*************************************************************************/
1698 : /** **/
1699 : /** (Z_K/I)^* **/
1700 : /** **/
1701 : /*************************************************************************/
1702 : GEN
1703 1302 : vecsmall01_to_indices(GEN v)
1704 : {
1705 1302 : long i, k, l = lg(v);
1706 1302 : GEN p = new_chunk(l) + l;
1707 3724 : for (k=1, i=l-1; i; i--)
1708 2422 : if (v[i]) { *--p = i; k++; }
1709 1302 : *--p = evallg(k) | evaltyp(t_VECSMALL);
1710 1302 : avma = (pari_sp)p; return p;
1711 : }
1712 : GEN
1713 360493 : vec01_to_indices(GEN v)
1714 : {
1715 : long i, k, l;
1716 : GEN p;
1717 :
1718 360493 : switch (typ(v))
1719 : {
1720 345793 : case t_VECSMALL: return v;
1721 14700 : case t_VEC: break;
1722 0 : default: pari_err_TYPE("vec01_to_indices",v);
1723 : }
1724 14700 : l = lg(v);
1725 14700 : p = new_chunk(l) + l;
1726 42427 : for (k=1, i=l-1; i; i--)
1727 27727 : if (signe(gel(v,i))) { *--p = i; k++; }
1728 14700 : *--p = evallg(k) | evaltyp(t_VECSMALL);
1729 14700 : avma = (pari_sp)p; return p;
1730 : }
1731 : GEN
1732 4830 : indices_to_vec01(GEN p, long r)
1733 : {
1734 4830 : long i, l = lg(p);
1735 4830 : GEN v = zerovec(r);
1736 4830 : for (i = 1; i < l; i++) gel(v, p[i]) = gen_1;
1737 4830 : return v;
1738 : }
1739 :
1740 : /* return (column) vector of R1 signatures of x (0 or 1) */
1741 : GEN
1742 345793 : nfsign_arch(GEN nf, GEN x, GEN arch)
1743 : {
1744 345793 : GEN sarch, M, V, archp = vec01_to_indices(arch);
1745 345793 : long i, s, np, n = lg(archp)-1;
1746 : pari_sp av;
1747 :
1748 345793 : if (!n) return cgetg(1,t_VECSMALL);
1749 345261 : nf = checknf(nf);
1750 345261 : if (typ(x) == t_MAT)
1751 : { /* factorisation */
1752 98207 : GEN g = gel(x,1), e = gel(x,2);
1753 98207 : V = zero_zv(n);
1754 287719 : for (i=1; i<lg(g); i++)
1755 189512 : if (mpodd(gel(e,i)))
1756 164144 : Flv_add_inplace(V, nfsign_arch(nf,gel(g,i),archp), 2);
1757 98207 : avma = (pari_sp)V; return V;
1758 : }
1759 247054 : av = avma; V = cgetg(n+1,t_VECSMALL);
1760 247054 : x = nf_to_scalar_or_basis(nf, x);
1761 247054 : switch(typ(x))
1762 : {
1763 : case t_INT:
1764 64495 : s = signe(x);
1765 64495 : if (!s) pari_err_DOMAIN("nfsign_arch","element","=",gen_0,x);
1766 64495 : avma = av; return const_vecsmall(n, (s < 0)? 1: 0);
1767 : case t_FRAC:
1768 35 : s = signe(gel(x,1));
1769 35 : avma = av; return const_vecsmall(n, (s < 0)? 1: 0);
1770 : }
1771 182524 : x = Q_primpart(x); M = nf_get_M(nf); sarch = NULL; np = -1;
1772 492020 : for (i = 1; i <= n; i++)
1773 : {
1774 310436 : long s = eval_sign(M, x, archp[i]);
1775 310436 : if (s < 0) /* failure */
1776 : {
1777 940 : long ni, r1 = nf_get_r1(nf);
1778 : GEN xi;
1779 940 : if (np < 0)
1780 : {
1781 940 : np = num_positive(nf, x);
1782 940 : if (np == 0) { avma = av; return const_vecsmall(n, 1); }
1783 802 : if (np == r1){ avma = av; return const_vecsmall(n, 0); }
1784 494 : sarch = nfarchstar(nf, NULL, identity_perm(r1));
1785 : }
1786 494 : xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
1787 494 : xi = Q_primpart(xi);
1788 494 : ni = num_positive(nf, nfmuli(nf,x,xi));
1789 494 : if (ni == 0) { avma = av; V = const_vecsmall(n, 1); V[i] = 0; return V; }
1790 356 : if (ni == r1){ avma = av; V = const_vecsmall(n, 0); V[i] = 1; return V; }
1791 0 : s = ni < np? 0: 1;
1792 : }
1793 309496 : V[i] = s;
1794 : }
1795 181584 : avma = (pari_sp)V; return V;
1796 : }
1797 : static void
1798 18557 : chk_ind(const char *s, long i, long r1)
1799 : {
1800 18557 : if (i <= 0)
1801 7 : pari_err_DOMAIN(s, "index", "<=", gen_0, stoi(i));
1802 18550 : if (i > r1)
1803 21 : pari_err_DOMAIN(s, "index", ">", stoi(r1), stoi(i));
1804 18529 : }
1805 : GEN
1806 770 : nfeltsign(GEN nf, GEN x, GEN ind0)
1807 : {
1808 770 : pari_sp av = avma;
1809 : long i, l, r1;
1810 : GEN v, ind;
1811 770 : nf = checknf(nf);
1812 770 : r1 = nf_get_r1(nf);
1813 770 : x = nf_to_scalar_or_basis(nf, x);
1814 770 : if (!ind0) ind0 = identity_perm(r1);
1815 770 : switch(typ(ind0))
1816 : {
1817 : case t_INT: case t_VEC: case t_COL:
1818 56 : ind = gtovecsmall(ind0); break;
1819 : case t_VECSMALL:
1820 714 : ind = ind0; break;
1821 : default:
1822 0 : pari_err_TYPE("nfeltsign",ind0);
1823 : return NULL; /* LCOV_EXCL_LINE */
1824 : }
1825 770 : l = lg(ind);
1826 770 : for (i = 1; i < l; i++) chk_ind("nfeltsign", ind[i], r1);
1827 749 : if (typ(x) != t_COL)
1828 : {
1829 : GEN s;
1830 21 : switch(gsigne(x))
1831 : {
1832 7 : case -1:s = gen_m1; break;
1833 7 : case 1: s = gen_1; break;
1834 7 : default: s = gen_0; break;
1835 : }
1836 21 : avma = av;
1837 21 : return typ(ind0) == t_INT? s: const_vec(l-1, s);
1838 : }
1839 728 : v = nfsign_arch(nf, x, ind);
1840 728 : if (typ(ind0) == t_INT) { avma = av; return v[1]? gen_m1: gen_1; }
1841 721 : settyp(v, t_VEC);
1842 721 : for (i = 1; i < l; i++) gel(v,i) = v[i]? gen_m1: gen_1;
1843 721 : return gerepileupto(av, v);
1844 :
1845 : }
1846 :
1847 : GEN
1848 4893 : nfeltembed(GEN nf, GEN x, GEN ind0)
1849 : {
1850 4893 : pari_sp av = avma;
1851 : long i, l, r1, r2;
1852 : GEN v, ind, cx, M;
1853 4893 : nf = checknf(nf);
1854 4893 : r1 = nf_get_r1(nf);
1855 4893 : r2 = nf_get_r2(nf);
1856 4893 : x = nf_to_scalar_or_basis(nf, x);
1857 4886 : if (!ind0) ind0 = identity_perm(r1+r2);
1858 4886 : switch(typ(ind0))
1859 : {
1860 : case t_INT: case t_VEC: case t_COL:
1861 42 : ind = gtovecsmall(ind0); break;
1862 : case t_VECSMALL:
1863 4844 : ind = ind0; break;
1864 : default:
1865 0 : pari_err_TYPE("nfeltsign",ind0);
1866 : return NULL; /* LCOV_EXCL_LINE */
1867 : }
1868 4886 : l = lg(ind);
1869 4886 : for (i = 1; i < l; i++) chk_ind("nfeltembed", ind[i], r1+r2);
1870 4879 : if (typ(x) != t_COL)
1871 : {
1872 2044 : if (typ(ind0) != t_INT) x = const_vec(l-1, x);
1873 2044 : return gerepilecopy(av, x);
1874 : }
1875 2835 : x = Q_primitive_part(x, &cx); M = nf_get_M(nf);
1876 2835 : v = cgetg(l, t_VEC);
1877 12866 : for (i = 1; i < l; i++)
1878 : {
1879 10031 : GEN t = nfembed_i(M, x, ind[i]);
1880 10031 : if (cx) t = gmul(t, cx);
1881 10031 : gel(v,i) = t;
1882 : }
1883 2835 : if (typ(ind0) == t_INT) v = gel(v,1);
1884 2835 : return gerepilecopy(av, v);
1885 : }
1886 :
1887 : /* return the vector of signs of x; the matrix of such if x is a vector
1888 : * of nf elements */
1889 : GEN
1890 1302 : nfsign(GEN nf, GEN x)
1891 : {
1892 : long i, l;
1893 : GEN archp, S;
1894 :
1895 1302 : nf = checknf(nf);
1896 1302 : archp = identity_perm( nf_get_r1(nf) );
1897 1302 : if (typ(x) != t_VEC) return nfsign_arch(nf, x, archp);
1898 182 : l = lg(x); S = cgetg(l, t_MAT);
1899 182 : for (i=1; i<l; i++) gel(S,i) = nfsign_arch(nf, gel(x,i), archp);
1900 182 : return S;
1901 : }
1902 :
1903 : /* x integral elt, A integral ideal in HNF; reduce x mod A */
1904 : static GEN
1905 602973 : zk_modHNF(GEN x, GEN A)
1906 602973 : { return (typ(x) == t_COL)? ZC_hnfrem(x, A): modii(x, gcoeff(A,1,1)); }
1907 :
1908 : /* given an element x in Z_K and an integral ideal y in HNF, coprime with x,
1909 : outputs an element inverse of x modulo y */
1910 : GEN
1911 147 : nfinvmodideal(GEN nf, GEN x, GEN y)
1912 : {
1913 147 : pari_sp av = avma;
1914 147 : GEN a, yZ = gcoeff(y,1,1);
1915 :
1916 147 : if (is_pm1(yZ)) return gen_0;
1917 147 : x = nf_to_scalar_or_basis(nf, x);
1918 147 : if (typ(x) == t_INT) return gerepileupto(av, Fp_inv(x, yZ));
1919 :
1920 77 : a = hnfmerge_get_1(idealhnf_principal(nf,x), y);
1921 77 : if (!a) pari_err_INV("nfinvmodideal", x);
1922 77 : return gerepileupto(av, zk_modHNF(nfdiv(nf,a,x), y));
1923 : }
1924 :
1925 : static GEN
1926 276029 : nfsqrmodideal(GEN nf, GEN x, GEN id)
1927 276029 : { return zk_modHNF(nfsqri(nf,x), id); }
1928 : static GEN
1929 736103 : nfmulmodideal(GEN nf, GEN x, GEN y, GEN id)
1930 736103 : { return x? zk_modHNF(nfmuli(nf,x,y), id): y; }
1931 : /* assume x integral, k integer, A in HNF */
1932 : GEN
1933 477261 : nfpowmodideal(GEN nf,GEN x,GEN k,GEN A)
1934 : {
1935 477261 : long s = signe(k);
1936 : pari_sp av;
1937 : GEN y;
1938 :
1939 477261 : if (!s) return gen_1;
1940 477261 : av = avma;
1941 477261 : x = nf_to_scalar_or_basis(nf, x);
1942 477261 : if (typ(x) != t_COL) return Fp_pow(x, k, gcoeff(A,1,1));
1943 226420 : if (s < 0) { x = nfinvmodideal(nf, x,A); k = absi(k); }
1944 226420 : for(y = NULL;;)
1945 : {
1946 502449 : if (mpodd(k)) y = nfmulmodideal(nf,y,x,A);
1947 502449 : k = shifti(k,-1); if (!signe(k)) break;
1948 276029 : x = nfsqrmodideal(nf,x,A);
1949 276029 : }
1950 226420 : return gerepileupto(av, y);
1951 : }
1952 :
1953 : /* a * g^n mod id */
1954 : static GEN
1955 424117 : nfmulpowmodideal(GEN nf, GEN a, GEN g, GEN n, GEN id)
1956 : {
1957 424117 : return nfmulmodideal(nf, a, nfpowmodideal(nf,g,n,id), id);
1958 : }
1959 :
1960 : /* assume (num(g[i]), id) = 1 for all i. Return prod g[i]^e[i] mod id.
1961 : * EX = multiple of exponent of (O_K/id)^* */
1962 : GEN
1963 192336 : famat_to_nf_modideal_coprime(GEN nf, GEN g, GEN e, GEN id, GEN EX)
1964 : {
1965 192336 : GEN EXo2, plus = NULL, minus = NULL, idZ = gcoeff(id,1,1);
1966 192336 : long i, lx = lg(g);
1967 :
1968 192336 : if (is_pm1(idZ)) return gen_1; /* id = Z_K */
1969 192336 : EXo2 = (expi(EX) > 10)? shifti(EX,-1): NULL;
1970 872282 : for (i = 1; i < lx; i++)
1971 : {
1972 679946 : GEN h, n = centermodii(gel(e,i), EX, EXo2);
1973 679946 : long sn = signe(n);
1974 679946 : if (!sn) continue;
1975 :
1976 324878 : h = nf_to_scalar_or_basis(nf, gel(g,i));
1977 324878 : switch(typ(h))
1978 : {
1979 200138 : case t_INT: break;
1980 : case t_FRAC:
1981 0 : h = Fp_div(gel(h,1), gel(h,2), idZ); break;
1982 : default:
1983 : {
1984 : GEN dh;
1985 124740 : h = Q_remove_denom(h, &dh);
1986 124740 : if (dh) h = FpC_Fp_mul(h, Fp_inv(dh,idZ), idZ);
1987 : }
1988 : }
1989 324878 : if (sn > 0)
1990 323471 : plus = nfmulpowmodideal(nf, plus, h, n, id);
1991 : else /* sn < 0 */
1992 1407 : minus = nfmulpowmodideal(nf, minus, h, absi(n), id);
1993 : }
1994 192336 : if (minus) plus = nfmulmodideal(nf, plus, nfinvmodideal(nf,minus,id), id);
1995 192336 : return plus? plus: gen_1;
1996 : }
1997 :
1998 : /* given 2 integral ideals x, y in HNF s.t x | y | x^2, compute (1+x)/(1+y) in
1999 : * the form [[cyc],[gen], U], where U := ux^-1 as a pair [ZM, denom(U)] */
2000 : static GEN
2001 20797 : zidealij(GEN x, GEN y)
2002 : {
2003 20797 : GEN U, G, cyc, xp = gcoeff(x,1,1), xi = hnf_invscale(x, xp);
2004 : long j, N;
2005 :
2006 : /* x^(-1) y = relations between the 1 + x_i (HNF) */
2007 20797 : cyc = ZM_snf_group(ZM_Z_divexact(ZM_mul(xi, y), xp), &U, &G);
2008 20797 : N = lg(cyc); G = ZM_mul(x,G); settyp(G, t_VEC); /* new generators */
2009 77623 : for (j=1; j<N; j++)
2010 : {
2011 56826 : GEN c = gel(G,j);
2012 56826 : gel(c,1) = addiu(gel(c,1), 1); /* 1 + g_j */
2013 56826 : if (ZV_isscalar(c)) gel(G,j) = gel(c,1);
2014 : }
2015 20797 : return mkvec4(cyc, G, ZM_mul(U,xi), xp);
2016 : }
2017 :
2018 : /* lg(x) > 1, x + 1; shallow */
2019 : static GEN
2020 6265 : ZC_add1(GEN x)
2021 : {
2022 6265 : long i, l = lg(x);
2023 6265 : GEN y = cgetg(l, t_COL);
2024 6265 : for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
2025 6265 : gel(y,1) = addiu(gel(x,1), 1); return y;
2026 : }
2027 : /* lg(x) > 1, x - 1; shallow */
2028 : static GEN
2029 3766 : ZC_sub1(GEN x)
2030 : {
2031 3766 : long i, l = lg(x);
2032 3766 : GEN y = cgetg(l, t_COL);
2033 3766 : for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
2034 3766 : gel(y,1) = subiu(gel(x,1), 1); return y;
2035 : }
2036 :
2037 : /* x,y are t_INT or ZC */
2038 : static GEN
2039 0 : zkadd(GEN x, GEN y)
2040 : {
2041 0 : long tx = typ(x);
2042 0 : if (tx == typ(y))
2043 0 : return tx == t_INT? addii(x,y): ZC_add(x,y);
2044 : else
2045 0 : return tx == t_INT? ZC_Z_add(y,x): ZC_Z_add(x,y);
2046 : }
2047 : /* x a t_INT or ZC, x+1; shallow */
2048 : static GEN
2049 12215 : zkadd1(GEN x)
2050 : {
2051 12215 : long tx = typ(x);
2052 12215 : return tx == t_INT? addiu(x,1): ZC_add1(x);
2053 : }
2054 : /* x a t_INT or ZC, x-1; shallow */
2055 : static GEN
2056 12215 : zksub1(GEN x)
2057 : {
2058 12215 : long tx = typ(x);
2059 12215 : return tx == t_INT? subiu(x,1): ZC_sub1(x);
2060 : }
2061 : /* x,y are t_INT or ZC; x - y */
2062 : static GEN
2063 0 : zksub(GEN x, GEN y)
2064 : {
2065 0 : long tx = typ(x), ty = typ(y);
2066 0 : if (tx == ty)
2067 0 : return tx == t_INT? subii(x,y): ZC_sub(x,y);
2068 : else
2069 0 : return tx == t_INT? Z_ZC_sub(x,y): ZC_Z_sub(x,y);
2070 : }
2071 : /* x is t_INT or ZM (mult. map), y is t_INT or ZC; x * y */
2072 : static GEN
2073 12215 : zkmul(GEN x, GEN y)
2074 : {
2075 12215 : long tx = typ(x), ty = typ(y);
2076 12215 : if (ty == t_INT)
2077 8449 : return tx == t_INT? mulii(x,y): ZC_Z_mul(gel(x,1),y);
2078 : else
2079 3766 : return tx == t_INT? ZC_Z_mul(y,x): ZM_ZC_mul(x,y);
2080 : }
2081 :
2082 : /* (U,V) = 1 coprime ideals. Want z = x mod U, = y mod V; namely
2083 : * z =vx + uy = v(x-y) + y, where u + v = 1, u in U, v in V.
2084 : * zkc = [v, UV], v a t_INT or ZM (mult. by v map), UV a ZM (ideal in HNF);
2085 : * shallow */
2086 : GEN
2087 0 : zkchinese(GEN zkc, GEN x, GEN y)
2088 : {
2089 0 : GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd(zkmul(v, zksub(x,y)), y);
2090 0 : return zk_modHNF(z, UV);
2091 : }
2092 : /* special case z = x mod U, = 1 mod V; shallow */
2093 : GEN
2094 12215 : zkchinese1(GEN zkc, GEN x)
2095 : {
2096 12215 : GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd1(zkmul(v, zksub1(x)));
2097 12215 : return (typ(z) == t_INT)? z: ZC_hnfrem(z, UV);
2098 : }
2099 : static GEN
2100 11011 : zkVchinese1(GEN zkc, GEN v)
2101 : {
2102 : long i, ly;
2103 11011 : GEN y = cgetg_copy(v, &ly);
2104 11011 : for (i=1; i<ly; i++) gel(y,i) = zkchinese1(zkc, gel(v,i));
2105 11011 : return y;
2106 : }
2107 :
2108 : /* prepare to solve z = x (mod A), z = y mod (B) [zkchinese or zkchinese1] */
2109 : GEN
2110 10752 : zkchineseinit(GEN nf, GEN A, GEN B, GEN AB)
2111 : {
2112 : GEN v;
2113 : long e;
2114 10752 : nf = checknf(nf);
2115 10752 : v = idealaddtoone_raw(nf, A, B);
2116 10752 : if ((e = gexpo(v)) > 5)
2117 : {
2118 455 : GEN b = (typ(v) == t_COL)? v: scalarcol_shallow(v, nf_get_degree(nf));
2119 455 : b= ZC_reducemodlll(b, AB);
2120 455 : if (gexpo(b) < e) v = b;
2121 : }
2122 10752 : return mkvec2(zk_scalar_or_multable(nf,v), AB);
2123 : }
2124 : /* prepare to solve z = x (mod A), z = 1 mod (B)
2125 : * and then z = 1 (mod A), z = y mod (B) [zkchinese1 twice] */
2126 : static GEN
2127 259 : zkchinese1init2(GEN nf, GEN A, GEN B, GEN AB)
2128 : {
2129 259 : GEN zkc = zkchineseinit(nf, A, B, AB);
2130 259 : GEN mv = gel(zkc,1), mu;
2131 259 : if (typ(mv) == t_INT) return mkvec2(zkc, mkvec2(subui(1,mv),AB));
2132 35 : mu = RgM_Rg_add_shallow(ZM_neg(mv), gen_1);
2133 35 : return mkvec2(mkvec2(mv,AB), mkvec2(mu,AB));
2134 : }
2135 :
2136 : static GEN
2137 356571 : apply_U(GEN L, GEN a)
2138 : {
2139 356571 : GEN e, U = gel(L,3), dU = gel(L,4);
2140 356571 : if (typ(a) == t_INT)
2141 130603 : e = ZC_Z_mul(gel(U,1), subiu(a, 1));
2142 : else
2143 : { /* t_COL */
2144 225968 : GEN t = gel(a,1);
2145 225968 : gel(a,1) = subiu(gel(a,1), 1); /* a -= 1 */
2146 225968 : e = ZM_ZC_mul(U, a);
2147 225968 : gel(a,1) = t; /* restore */
2148 : }
2149 356571 : return gdiv(e, dU);
2150 : }
2151 :
2152 : /* true nf; vectors of [[cyc],[g],U.X^-1]. Assume k > 1. */
2153 : static GEN
2154 14133 : principal_units(GEN nf, GEN pr, long k, GEN prk)
2155 : {
2156 : GEN list, prb;
2157 14133 : ulong mask = quadratic_prec_mask(k);
2158 14133 : long a = 1;
2159 :
2160 14133 : if (DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",k,pr);
2161 14133 : prb = pr_hnf(nf,pr);
2162 14133 : list = vectrunc_init(k);
2163 49063 : while (mask > 1)
2164 : {
2165 20797 : GEN pra = prb;
2166 20797 : long b = a << 1;
2167 :
2168 20797 : if (mask & 1) b--;
2169 20797 : mask >>= 1;
2170 : /* compute 1 + pr^a / 1 + pr^b, 2a <= b */
2171 20797 : if(DEBUGLEVEL>3) err_printf(" treating a = %ld, b = %ld\n",a,b);
2172 20797 : prb = (b >= k)? prk: idealpows(nf,pr,b);
2173 20797 : vectrunc_append(list, zidealij(pra, prb));
2174 20797 : a = b;
2175 : }
2176 14133 : return list;
2177 : }
2178 : /* a = 1 mod (pr) return log(a) on local-gens of 1+pr/1+pr^k */
2179 : static GEN
2180 227008 : log_prk1(GEN nf, GEN a, long nh, GEN L2, GEN prk)
2181 : {
2182 227008 : GEN y = cgetg(nh+1, t_COL);
2183 227008 : long j, iy, c = lg(L2)-1;
2184 583572 : for (j = iy = 1; j <= c; j++)
2185 : {
2186 356571 : GEN L = gel(L2,j), cyc = gel(L,1), gen = gel(L,2), E = apply_U(L,a);
2187 356571 : long i, nc = lg(cyc)-1;
2188 356571 : int last = (j == c);
2189 1245688 : for (i = 1; i <= nc; i++, iy++)
2190 : {
2191 889124 : GEN t, e = gel(E,i);
2192 889124 : if (typ(e) != t_INT) pari_err_COPRIME("zlog_prk1", a, prk);
2193 889117 : t = Fp_neg(e, gel(cyc,i));
2194 889117 : gel(y,iy) = negi(t);
2195 889117 : if (!last && signe(t)) a = nfmulpowmodideal(nf, a, gel(gen,i), t, prk);
2196 : }
2197 : }
2198 227001 : return y;
2199 : }
2200 : /* true nf */
2201 : static GEN
2202 5768 : principal_units_relations(GEN nf, GEN L2, GEN prk, long nh)
2203 : {
2204 5768 : GEN h = cgetg(nh+1,t_MAT);
2205 5768 : long ih, j, c = lg(L2)-1;
2206 18200 : for (j = ih = 1; j <= c; j++)
2207 : {
2208 12432 : GEN L = gel(L2,j), F = gel(L,1), G = gel(L,2);
2209 12432 : long k, lG = lg(G);
2210 51534 : for (k = 1; k < lG; k++,ih++)
2211 : { /* log(g^f) mod pr^e */
2212 39102 : GEN a = nfpowmodideal(nf,gel(G,k),gel(F,k),prk);
2213 39102 : gel(h,ih) = ZC_neg(log_prk1(nf, a, nh, L2, prk));
2214 39102 : gcoeff(h,ih,ih) = gel(F,k);
2215 : }
2216 : }
2217 5768 : return h;
2218 : }
2219 : /* true nf; e > 1; multiplicative group (1 + pr) / (1 + pr^k),
2220 : * prk = pr^k or NULL */
2221 : static GEN
2222 14133 : idealprincipalunits_i(GEN nf, GEN pr, long k, GEN *pU)
2223 : {
2224 14133 : GEN cyc, gen, L2, prk = idealpows(nf, pr, k);
2225 :
2226 14133 : L2 = principal_units(nf, pr, k, prk);
2227 14133 : if (k == 2)
2228 : {
2229 8365 : GEN L = gel(L2,1);
2230 8365 : cyc = gel(L,1);
2231 8365 : gen = gel(L,2);
2232 8365 : if (pU) *pU = matid(lg(gen)-1);
2233 : }
2234 : else
2235 : {
2236 5768 : long c = lg(L2), j;
2237 5768 : GEN EX, h, Ui, vg = cgetg(c, t_VEC);
2238 5768 : for (j = 1; j < c; j++) gel(vg, j) = gmael(L2,j,2);
2239 5768 : vg = shallowconcat1(vg);
2240 5768 : h = principal_units_relations(nf, L2, prk, lg(vg)-1);
2241 5768 : h = ZM_hnfall_i(h, NULL, 0);
2242 5768 : cyc = ZM_snf_group(h, pU, &Ui);
2243 5768 : c = lg(Ui); gen = cgetg(c, t_VEC); EX = gel(cyc,1);
2244 32592 : for (j = 1; j < c; j++)
2245 26824 : gel(gen,j) = famat_to_nf_modideal_coprime(nf, vg, gel(Ui,j), prk, EX);
2246 : }
2247 14133 : return mkvec4(cyc, gen, prk, L2);
2248 : }
2249 : GEN
2250 112 : idealprincipalunits(GEN nf, GEN pr, long k)
2251 : {
2252 : pari_sp av;
2253 : GEN v;
2254 112 : nf = checknf(nf);
2255 112 : if (k == 1) { checkprid(pr); retmkvec3(gen_1,cgetg(1,t_VEC),cgetg(1,t_VEC)); }
2256 105 : av = avma; v = idealprincipalunits_i(nf, pr, k, NULL);
2257 105 : return gerepilecopy(av, mkvec3(powiu(pr_norm(pr), k-1), gel(v,1), gel(v,2)));
2258 : }
2259 :
2260 : /* true nf; given an ideal pr^k dividing an integral ideal x (in HNF form)
2261 : * compute an 'sprk', the structure of G = (Z_K/pr^k)^* [ x = NULL for x=pr^k ]
2262 : * Return a vector with at least 4 components [cyc],[gen],[HNF pr^k,pr,k],ff,
2263 : * where
2264 : * cyc : type of G as abelian group (SNF)
2265 : * gen : generators of G, coprime to x
2266 : * pr^k: in HNF
2267 : * ff : data for log_g in (Z_K/pr)^*
2268 : * Two extra components are present iff k > 1: L2, U
2269 : * L2 : list of data structures to compute local DL in (Z_K/pr)^*,
2270 : * and 1 + pr^a/ 1 + pr^b for various a < b <= min(2a, k)
2271 : * U : base change matrices to convert a vector of local DL to DL wrt gen */
2272 : static GEN
2273 31654 : sprkinit(GEN nf, GEN pr, GEN gk, GEN x)
2274 : {
2275 : GEN T, p, modpr, cyc, gen, g, g0, ord0, A, prk, U, L2;
2276 31654 : long k = itos(gk), f = pr_get_f(pr);
2277 :
2278 31654 : if(DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",k,pr);
2279 31654 : modpr = nf_to_Fq_init(nf, &pr,&T,&p);
2280 : /* (Z_K / pr)^* */
2281 31654 : if (f == 1)
2282 : {
2283 22659 : g0 = g = pgener_Fp(p);
2284 22659 : ord0 = get_arith_ZZM(subiu(p,1));
2285 : }
2286 : else
2287 : {
2288 8995 : g0 = g = gener_FpXQ(T,p, &ord0);
2289 8995 : g = Fq_to_nf(g, modpr);
2290 8995 : if (typ(g) == t_POL) g = poltobasis(nf, g);
2291 : }
2292 31654 : A = gel(ord0, 1); /* Norm(pr)-1 */
2293 31654 : if (k == 1)
2294 : {
2295 17626 : cyc = mkvec(A);
2296 17626 : gen = mkvec(g);
2297 17626 : prk = pr_hnf(nf,pr);
2298 17626 : L2 = U = NULL;
2299 : }
2300 : else
2301 : { /* local-gens of (1 + pr)/(1 + pr^k) = SNF-gens * U */
2302 : GEN AB, B, u, v, w;
2303 : long j, l;
2304 14028 : w = idealprincipalunits_i(nf, pr, k, &U);
2305 : /* incorporate (Z_K/pr)^*, order A coprime to B = expo(1+pr/1+pr^k)*/
2306 14028 : cyc = leafcopy(gel(w,1)); B = gel(cyc,1); AB = mulii(A,B);
2307 14028 : gen = leafcopy(gel(w,2));
2308 14028 : prk = gel(w,3);
2309 14028 : g = nfpowmodideal(nf, g, B, prk);
2310 14028 : g0 = Fq_pow(g0, modii(B,A), T, p); /* update primitive root */
2311 14028 : L2 = mkvec3(A, g, gel(w,4));
2312 14028 : gel(cyc,1) = AB;
2313 14028 : gel(gen,1) = nfmulmodideal(nf, gel(gen,1), g, prk);
2314 14028 : u = mulii(Fp_inv(A,B), A);
2315 14028 : v = subui(1, u); l = lg(U);
2316 14028 : for (j = 1; j < l; j++) gcoeff(U,1,j) = Fp_mul(u, gcoeff(U,1,j), AB);
2317 14028 : U = mkvec2(Rg_col_ei(v, lg(gen)-1, 1), U);
2318 : }
2319 : /* local-gens of (Z_K/pr^k)^* = SNF-gens * U */
2320 31654 : if (x)
2321 : {
2322 10493 : GEN uv = zkchineseinit(nf, idealdivpowprime(nf,x,pr,gk), prk, x);
2323 10493 : gen = zkVchinese1(uv, gen);
2324 : }
2325 31654 : return mkvecn(U? 6: 4, cyc, gen, prk, mkvec3(modpr,g0,ord0), L2, U);
2326 : }
2327 : static GEN
2328 348847 : sprk_get_cyc(GEN s) { return gel(s,1); }
2329 : static GEN
2330 110534 : sprk_get_expo(GEN s)
2331 : {
2332 110534 : GEN cyc = sprk_get_cyc(s);
2333 110534 : return lg(cyc) == 1? gen_1: gel(cyc, 1);
2334 : }
2335 : static GEN
2336 25592 : sprk_get_gen(GEN s) { return gel(s,2); }
2337 : static GEN
2338 298440 : sprk_get_prk(GEN s) { return gel(s,3); }
2339 : static GEN
2340 388535 : sprk_get_ff(GEN s) { return gel(s,4); }
2341 : static GEN
2342 129532 : sprk_get_pr(GEN s) { GEN ff = gel(s,4); return modpr_get_pr(gel(ff,1)); }
2343 : /* A = Npr-1, <g> = (Z_K/pr)^*, L2 to 1 + pr / 1 + pr^k */
2344 : static void
2345 198812 : sprk_get_L2(GEN s, GEN *A, GEN *g, GEN *L2)
2346 198812 : { GEN v = gel(s,5); *A = gel(v,1); *g = gel(v,2); *L2 = gel(v,3); }
2347 : static void
2348 187906 : sprk_get_U2(GEN s, GEN *U1, GEN *U2)
2349 187906 : { GEN v = gel(s,6); *U1 = gel(v,1); *U2 = gel(v,2); }
2350 : static int
2351 388535 : sprk_is_prime(GEN s) { return lg(s) == 5; }
2352 :
2353 : static GEN
2354 110534 : famat_zlog_pr(GEN nf, GEN g, GEN e, GEN sprk)
2355 : {
2356 110534 : GEN pr = sprk_get_pr(sprk);
2357 110534 : GEN prk = sprk_get_prk(sprk);
2358 110534 : GEN x = famat_makecoprime(nf, g, e, pr, prk, sprk_get_expo(sprk));
2359 110534 : return zlog_pr(nf, x, sprk);
2360 : }
2361 : /* log_g(a) in (Z_K/pr)^* */
2362 : static GEN
2363 388535 : nf_log(GEN nf, GEN a, GEN ff)
2364 : {
2365 388535 : GEN pr = gel(ff,1), g = gel(ff,2), ord = gel(ff,3);
2366 388535 : GEN T,p, modpr = nf_to_Fq_init(nf, &pr, &T, &p);
2367 388535 : return Fq_log(nf_to_Fq(nf,a,modpr), g, ord, T, p);
2368 : }
2369 : /* a in Z_K (t_COL or t_INT), pr prime ideal, sprk = sprkinit(nf,pr,k,x).
2370 : * return log(a) on SNF-generators of (Z_K/pr^k)^**/
2371 : GEN
2372 389571 : zlog_pr(GEN nf, GEN a, GEN sprk)
2373 : {
2374 : GEN e, prk, A, g, L2, U1, U2, y;
2375 :
2376 389571 : if (typ(a) == t_MAT) return famat_zlog_pr(nf, gel(a,1), gel(a,2), sprk);
2377 :
2378 388535 : e = nf_log(nf, a, sprk_get_ff(sprk));
2379 388535 : if (sprk_is_prime(sprk)) return mkcol(e); /* k = 1 */
2380 187906 : prk = sprk_get_prk(sprk);
2381 187906 : sprk_get_L2(sprk, &A,&g,&L2);
2382 187906 : if (signe(e))
2383 : {
2384 46972 : e = Fp_neg(e, A);
2385 46972 : a = nfmulpowmodideal(nf, a, g, e, prk);
2386 : }
2387 187906 : sprk_get_U2(sprk, &U1,&U2);
2388 187906 : y = ZM_ZC_mul(U2, log_prk1(nf, a, lg(U2)-1, L2, prk));
2389 187899 : if (signe(e)) y = ZC_sub(y, ZC_Z_mul(U1,e));
2390 187899 : return vecmodii(y, sprk_get_cyc(sprk));
2391 : }
2392 : GEN
2393 6062 : zlog_pr_init(GEN nf, GEN pr, long k)
2394 6062 : { return sprkinit(checknf(nf),pr,utoipos(k),NULL);}
2395 : GEN
2396 378 : vzlog_pr(GEN nf, GEN v, GEN sprk)
2397 : {
2398 378 : long l = lg(v), i;
2399 378 : GEN w = cgetg(l, t_MAT);
2400 378 : for (i = 1; i < l; i++) gel(w,i) = zlog_pr(nf, gel(v,i), sprk);
2401 378 : return w;
2402 : }
2403 :
2404 : static GEN
2405 113474 : famat_zlog(GEN nf, GEN fa, GEN sgn, zlog_S *S)
2406 : {
2407 113474 : long i, n0, n = lg(S->U)-1;
2408 : GEN g, e, y;
2409 113474 : if (lg(fa) == 1) return zerocol(n);
2410 113474 : g = gel(fa,1);
2411 113474 : e = gel(fa,2);
2412 113474 : y = cgetg(n+1, t_COL);
2413 113474 : n0 = lg(S->sprk)-1; /* n0 = n (trivial arch. part), or n-1 */
2414 113474 : for (i=1; i <= n0; i++) gel(y,i) = famat_zlog_pr(nf, g, e, gel(S->sprk,i));
2415 113474 : if (n0 != n)
2416 : {
2417 93174 : if (!sgn) sgn = nfsign_arch(nf, fa, S->archp);
2418 93174 : gel(y,n) = Flc_to_ZC(sgn);
2419 : }
2420 113474 : return y;
2421 : }
2422 :
2423 : /* assume that cyclic factors are normalized, in particular != [1] */
2424 : static GEN
2425 25956 : split_U(GEN U, GEN Sprk)
2426 : {
2427 25956 : long t = 0, k, n, l = lg(Sprk);
2428 25956 : GEN vU = cgetg(l+1, t_VEC);
2429 50778 : for (k = 1; k < l; k++)
2430 : {
2431 24822 : n = lg(sprk_get_cyc(gel(Sprk,k))) - 1; /* > 0 */
2432 24822 : gel(vU,k) = vecslice(U, t+1, t+n);
2433 24822 : t += n;
2434 : }
2435 : /* t+1 .. lg(U)-1 */
2436 25956 : n = lg(U) - t - 1; /* can be 0 */
2437 25956 : if (!n) setlg(vU,l); else gel(vU,l) = vecslice(U, t+1, t+n);
2438 25956 : return vU;
2439 : }
2440 :
2441 : void
2442 353912 : init_zlog(zlog_S *S, GEN bid)
2443 : {
2444 353912 : GEN fa2 = bid_get_fact2(bid);
2445 353912 : S->U = bid_get_U(bid);
2446 353912 : S->hU = lg(bid_get_cyc(bid))-1;
2447 353912 : S->archp = bid_get_archp(bid);
2448 353912 : S->sprk = bid_get_sprk(bid);
2449 353912 : S->bid = bid;
2450 353912 : S->P = gel(fa2,1);
2451 353912 : S->k = gel(fa2,2);
2452 353912 : S->no2 = lg(S->P) == lg(gel(bid_get_fact(bid),1));
2453 353912 : }
2454 :
2455 : /* a a t_FRAC/t_INT, reduce mod bid */
2456 : static GEN
2457 7 : Q_mod_bid(GEN bid, GEN a)
2458 : {
2459 7 : GEN xZ = gcoeff(bid_get_ideal(bid),1,1);
2460 7 : GEN b = Rg_to_Fp(a, xZ);
2461 7 : if (gsigne(a) < 0) b = subii(b, xZ);
2462 7 : return b;
2463 : }
2464 : /* Return decomposition of a on the CRT generators blocks attached to the
2465 : * S->sprk and sarch; sgn = sign(a, S->arch), NULL if unknown */
2466 : static GEN
2467 246679 : zlog(GEN nf, GEN a, GEN sgn, zlog_S *S)
2468 : {
2469 : long k, l;
2470 : GEN y;
2471 246679 : a = nf_to_scalar_or_basis(nf, a);
2472 246679 : switch(typ(a))
2473 : {
2474 63861 : case t_INT: break;
2475 7 : case t_FRAC: a = Q_mod_bid(S->bid, a); break;
2476 : default: /* case t_COL: */
2477 : {
2478 : GEN den;
2479 182811 : check_nfelt(a, &den);
2480 182811 : if (den)
2481 : {
2482 46659 : a = Q_muli_to_int(a, den);
2483 46659 : a = mkmat2(mkcol2(a, den), mkcol2(gen_1, gen_m1));
2484 46659 : return famat_zlog(nf, a, sgn, S);
2485 : }
2486 : }
2487 : }
2488 200020 : if (sgn)
2489 34531 : sgn = (lg(sgn) == 1)? NULL: leafcopy(sgn);
2490 : else
2491 165489 : sgn = (lg(S->archp) == 1)? NULL: nfsign_arch(nf, a, S->archp);
2492 200020 : l = lg(S->sprk);
2493 200020 : y = cgetg(sgn? l+1: l, t_COL);
2494 443420 : for (k = 1; k < l; k++)
2495 : {
2496 243407 : GEN sprk = gel(S->sprk,k);
2497 243407 : gel(y,k) = zlog_pr(nf, a, sprk);
2498 : }
2499 200013 : if (sgn) gel(y,l) = Flc_to_ZC(sgn);
2500 200013 : return y;
2501 : }
2502 :
2503 : /* true nf */
2504 : GEN
2505 8344 : pr_basis_perm(GEN nf, GEN pr)
2506 : {
2507 8344 : long f = pr_get_f(pr);
2508 : GEN perm;
2509 8344 : if (f == nf_get_degree(nf)) return identity_perm(f);
2510 6944 : perm = cgetg(f+1, t_VECSMALL);
2511 6944 : perm[1] = 1;
2512 6944 : if (f > 1)
2513 : {
2514 399 : GEN H = pr_hnf(nf,pr);
2515 : long i, k;
2516 1463 : for (i = k = 2; k <= f; i++)
2517 : {
2518 1064 : if (is_pm1(gcoeff(H,i,i))) continue;
2519 840 : perm[k++] = i;
2520 : }
2521 : }
2522 6944 : return perm;
2523 : }
2524 :
2525 : /* \sum U[i]*y[i], U[i] ZM, y[i] ZC. We allow lg(y) > lg(U). */
2526 : static GEN
2527 313487 : ZMV_ZCV_mul(GEN U, GEN y)
2528 : {
2529 313487 : long i, l = lg(U);
2530 313487 : GEN z = NULL;
2531 313487 : if (l == 1) return cgetg(1,t_COL);
2532 863980 : for (i = 1; i < l; i++)
2533 : {
2534 550493 : GEN u = ZM_ZC_mul(gel(U,i), gel(y,i));
2535 550493 : z = z? ZC_add(z, u): u;
2536 : }
2537 313487 : return z;
2538 : }
2539 : /* A * (U[1], ..., U[d] */
2540 : static GEN
2541 518 : ZM_ZMV_mul(GEN A, GEN U)
2542 : {
2543 : long i, l;
2544 518 : GEN V = cgetg_copy(U,&l);
2545 518 : for (i = 1; i < l; i++) gel(V,i) = ZM_mul(A,gel(U,i));
2546 518 : return V;
2547 : }
2548 :
2549 : /* Log on bid.gen of generators of P_{1,I pr^{e-1}} / P_{1,I pr^e} (I,pr) = 1,
2550 : * defined implicitly via CRT. 'ind' is the index of pr in modulus
2551 : * factorization */
2552 : GEN
2553 50820 : log_gen_pr(zlog_S *S, long ind, GEN nf, long e)
2554 : {
2555 50820 : GEN A, sprk = gel(S->sprk,ind), Uind = gel(S->U, ind);
2556 :
2557 50820 : if (e == 1) retmkmat( gel(Uind,1) );
2558 : else
2559 : {
2560 18998 : GEN G, pr = sprk_get_pr(sprk);
2561 : long i, l;
2562 18998 : if (e == 2)
2563 : {
2564 10906 : GEN A, g, L, L2; sprk_get_L2(sprk,&A,&g,&L2); L = gel(L2,1);
2565 10906 : G = gel(L,2); l = lg(G);
2566 : }
2567 : else
2568 : {
2569 8092 : GEN perm = pr_basis_perm(nf,pr), PI = nfpow_u(nf, pr_get_gen(pr), e-1);
2570 8092 : l = lg(perm);
2571 8092 : G = cgetg(l, t_VEC);
2572 8092 : if (typ(PI) == t_INT)
2573 : { /* zk_ei_mul doesn't allow t_INT */
2574 1393 : long N = nf_get_degree(nf);
2575 1393 : gel(G,1) = addiu(PI,1);
2576 2261 : for (i = 2; i < l; i++)
2577 : {
2578 868 : GEN z = col_ei(N, 1);
2579 868 : gel(G,i) = z; gel(z, perm[i]) = PI;
2580 : }
2581 : }
2582 : else
2583 : {
2584 6699 : gel(G,1) = nfadd(nf, gen_1, PI);
2585 6909 : for (i = 2; i < l; i++)
2586 210 : gel(G,i) = nfadd(nf, gen_1, zk_ei_mul(nf, PI, perm[i]));
2587 : }
2588 : }
2589 18998 : A = cgetg(l, t_MAT);
2590 40761 : for (i = 1; i < l; i++)
2591 21763 : gel(A,i) = ZM_ZC_mul(Uind, zlog_pr(nf, gel(G,i), sprk));
2592 18998 : return A;
2593 : }
2594 : }
2595 : /* Log on bid.gen of generator of P_{1,f} / P_{1,f v[index]}
2596 : * v = vector of r1 real places */
2597 : GEN
2598 9975 : log_gen_arch(zlog_S *S, long index)
2599 : {
2600 9975 : GEN U = gel(S->U, lg(S->U)-1);
2601 9975 : return gel(U, index);
2602 : }
2603 :
2604 : /* compute bid.clgp: [h,cyc] or [h,cyc,gen] */
2605 : static GEN
2606 27013 : bid_grp(GEN nf, GEN U, GEN cyc, GEN g, GEN F, GEN sarch)
2607 : {
2608 27013 : GEN G, h = ZV_prod(cyc);
2609 : long c;
2610 27013 : if (!U) return mkvec2(h,cyc);
2611 26761 : c = lg(U);
2612 26761 : G = cgetg(c,t_VEC);
2613 26761 : if (c > 1)
2614 : {
2615 22407 : GEN U0, Uoo, EX = gel(cyc,1); /* exponent of bid */
2616 22407 : long i, hU = nbrows(U), nba = lg(sarch_get_cyc(sarch))-1; /* #f_oo */
2617 22407 : if (!nba) { U0 = U; Uoo = NULL; }
2618 11704 : else if (nba == hU) { U0 = NULL; Uoo = U; }
2619 : else
2620 : {
2621 9471 : U0 = rowslice(U, 1, hU-nba);
2622 9471 : Uoo = rowslice(U, hU-nba+1, hU);
2623 : }
2624 64288 : for (i = 1; i < c; i++)
2625 : {
2626 41881 : GEN t = gen_1;
2627 41881 : if (U0) t = famat_to_nf_modideal_coprime(nf, g, gel(U0,i), F, EX);
2628 41881 : if (Uoo) t = set_sign_mod_divisor(nf, ZV_to_Flv(gel(Uoo,i),2), t, sarch);
2629 41881 : gel(G,i) = t;
2630 : }
2631 : }
2632 26761 : return mkvec3(h, cyc, G);
2633 : }
2634 :
2635 : /* remove prime ideals of norm 2 with exponent 1 from factorization */
2636 : static GEN
2637 26698 : famat_strip2(GEN fa)
2638 : {
2639 26698 : GEN P = gel(fa,1), E = gel(fa,2), Q, F;
2640 26698 : long l = lg(P), i, j;
2641 26698 : Q = cgetg(l, t_COL);
2642 26698 : F = cgetg(l, t_COL);
2643 56315 : for (i = j = 1; i < l; i++)
2644 : {
2645 29617 : GEN pr = gel(P,i), e = gel(E,i);
2646 29617 : if (!absequaliu(pr_get_p(pr), 2) || itou(e) != 1 || pr_get_f(pr) != 1)
2647 : {
2648 25592 : gel(Q,j) = pr;
2649 25592 : gel(F,j) = e; j++;
2650 : }
2651 : }
2652 26698 : setlg(Q,j);
2653 26698 : setlg(F,j); return mkmat2(Q,F);
2654 : }
2655 :
2656 : /* Compute [[ideal,arch], [h,[cyc],[gen]], idealfact, [liste], U]
2657 : flag may include nf_GEN | nf_INIT */
2658 : static GEN
2659 26719 : Idealstar_i(GEN nf, GEN ideal, long flag)
2660 : {
2661 : long i, k, nbp, R1;
2662 26719 : GEN y, cyc, U, u1 = NULL, fa, fa2, sprk, x, arch, archp, E, P, sarch, gen;
2663 :
2664 26719 : nf = checknf(nf);
2665 26719 : R1 = nf_get_r1(nf);
2666 26719 : if (typ(ideal) == t_VEC && lg(ideal) == 3)
2667 : {
2668 12957 : arch = gel(ideal,2);
2669 12957 : ideal= gel(ideal,1);
2670 12957 : switch(typ(arch))
2671 : {
2672 : case t_VEC:
2673 12922 : if (lg(arch) != R1+1)
2674 0 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
2675 12922 : archp = vec01_to_indices(arch);
2676 12922 : break;
2677 : case t_VECSMALL:
2678 35 : archp = arch;
2679 35 : k = lg(archp)-1;
2680 35 : if (k && archp[k] > R1)
2681 7 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
2682 28 : arch = indices_to_vec01(archp, R1);
2683 28 : break;
2684 : default:
2685 0 : pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
2686 0 : return NULL;
2687 : }
2688 12950 : }
2689 : else
2690 : {
2691 13762 : arch = zerovec(R1);
2692 13762 : archp = cgetg(1, t_VECSMALL);
2693 : }
2694 26712 : if (is_nf_factor(ideal))
2695 : {
2696 350 : fa = ideal;
2697 350 : x = idealfactorback(nf, gel(fa,1), gel(fa,2), 0);
2698 : }
2699 : else
2700 : {
2701 26362 : fa = idealfactor(nf, ideal);
2702 26355 : x = ideal;
2703 : }
2704 26705 : if (typ(x) != t_MAT) x = idealhnf_shallow(nf, x);
2705 26705 : if (lg(x) == 1) pari_err_DOMAIN("Idealstar", "ideal","=",gen_0,x);
2706 26705 : if (typ(gcoeff(x,1,1)) != t_INT)
2707 7 : pari_err_DOMAIN("Idealstar","denominator(ideal)", "!=",gen_1,x);
2708 26698 : sarch = nfarchstar(nf, x, archp);
2709 26698 : fa2 = famat_strip2(fa);
2710 26698 : P = gel(fa2,1);
2711 26698 : E = gel(fa2,2);
2712 26698 : nbp = lg(P)-1;
2713 26698 : sprk = cgetg(nbp+1,t_VEC);
2714 26698 : if (nbp)
2715 : {
2716 20097 : GEN t = (nbp==1)? NULL: x;
2717 20097 : cyc = cgetg(nbp+2,t_VEC);
2718 20097 : gen = cgetg(nbp+1,t_VEC);
2719 45689 : for (i = 1; i <= nbp; i++)
2720 : {
2721 25592 : GEN L = sprkinit(nf, gel(P,i), gel(E,i), t);
2722 25592 : gel(sprk,i) = L;
2723 25592 : gel(cyc,i) = sprk_get_cyc(L);
2724 : /* true gens are congruent to those mod x AND positive at archp */
2725 25592 : gel(gen,i) = sprk_get_gen(L);
2726 : }
2727 20097 : gel(cyc,i) = sarch_get_cyc(sarch);
2728 20097 : cyc = shallowconcat1(cyc);
2729 20097 : gen = shallowconcat1(gen);
2730 20097 : cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
2731 : }
2732 : else
2733 : {
2734 6601 : cyc = sarch_get_cyc(sarch);
2735 6601 : gen = cgetg(1,t_VEC);
2736 6601 : U = matid(lg(cyc)-1);
2737 6601 : if (flag & nf_GEN) u1 = U;
2738 : }
2739 26698 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
2740 26698 : if (!(flag & nf_INIT)) return y;
2741 25900 : U = split_U(U, sprk);
2742 25900 : return mkvec5(mkvec2(x, arch), y, mkvec2(fa,fa2), mkvec2(sprk, sarch), U);
2743 : }
2744 : GEN
2745 26446 : Idealstar(GEN nf, GEN ideal, long flag)
2746 : {
2747 26446 : pari_sp av = avma;
2748 26446 : if (!nf) nf = nfinit(pol_x(0), DEFAULTPREC);
2749 26446 : return gerepilecopy(av, Idealstar_i(nf, ideal, flag));
2750 : }
2751 : GEN
2752 273 : Idealstarprk(GEN nf, GEN pr, long k, long flag)
2753 : {
2754 273 : pari_sp av = avma;
2755 273 : GEN z = Idealstar_i(nf, mkmat2(mkcol(pr),mkcols(k)), flag);
2756 273 : return gerepilecopy(av, z);
2757 : }
2758 :
2759 : /* FIXME: obsolete */
2760 : GEN
2761 0 : zidealstarinitgen(GEN nf, GEN ideal)
2762 0 : { return Idealstar(nf,ideal, nf_INIT|nf_GEN); }
2763 : GEN
2764 0 : zidealstarinit(GEN nf, GEN ideal)
2765 0 : { return Idealstar(nf,ideal, nf_INIT); }
2766 : GEN
2767 0 : zidealstar(GEN nf, GEN ideal)
2768 0 : { return Idealstar(nf,ideal, nf_GEN); }
2769 :
2770 : GEN
2771 63 : idealstar0(GEN nf, GEN ideal,long flag)
2772 : {
2773 63 : switch(flag)
2774 : {
2775 0 : case 0: return Idealstar(nf,ideal, nf_GEN);
2776 49 : case 1: return Idealstar(nf,ideal, nf_INIT);
2777 14 : case 2: return Idealstar(nf,ideal, nf_INIT|nf_GEN);
2778 0 : default: pari_err_FLAG("idealstar");
2779 : }
2780 : return NULL; /* LCOV_EXCL_LINE */
2781 : }
2782 :
2783 : void
2784 182811 : check_nfelt(GEN x, GEN *den)
2785 : {
2786 182811 : long l = lg(x), i;
2787 182811 : GEN t, d = NULL;
2788 182811 : if (typ(x) != t_COL) pari_err_TYPE("check_nfelt", x);
2789 666159 : for (i=1; i<l; i++)
2790 : {
2791 483348 : t = gel(x,i);
2792 483348 : switch (typ(t))
2793 : {
2794 387946 : case t_INT: break;
2795 : case t_FRAC:
2796 95402 : if (!d) d = gel(t,2); else d = lcmii(d, gel(t,2));
2797 95402 : break;
2798 0 : default: pari_err_TYPE("check_nfelt", x);
2799 : }
2800 : }
2801 182811 : *den = d;
2802 182811 : }
2803 :
2804 : GEN
2805 1209640 : vecmodii(GEN a, GEN b)
2806 : {
2807 : long i, l;
2808 1209640 : GEN c = cgetg_copy(a, &l);
2809 1209640 : for (i = 1; i < l; i++) gel(c,i) = modii(gel(a,i), gel(b,i));
2810 1209640 : return c;
2811 : }
2812 : GEN
2813 8253 : vecmoduu(GEN a, GEN b)
2814 : {
2815 : long i, l;
2816 8253 : GEN c = cgetg_copy(a, &l);
2817 8253 : for (i = 1; i < l; i++) c[i] = a[i] % b[i];
2818 8253 : return c;
2819 : }
2820 :
2821 : static GEN
2822 315083 : ideallog_i(GEN nf, GEN x, GEN sgn, zlog_S *S)
2823 : {
2824 315083 : pari_sp av = avma;
2825 : GEN y, cyc;
2826 315083 : if (!S->hU) return cgetg(1, t_COL);
2827 313501 : cyc = bid_get_cyc(S->bid);
2828 313501 : if (typ(x) == t_MAT)
2829 : {
2830 66822 : if (lg(x) == 1) return zerocol(lg(cyc)-1);
2831 66815 : y = famat_zlog(nf, x, sgn, S);
2832 : }
2833 : else
2834 246679 : y = zlog(nf, x, sgn, S);
2835 313487 : y = ZMV_ZCV_mul(S->U, y);
2836 313487 : return gerepileupto(av, vecmodii(y, cyc));
2837 : }
2838 :
2839 : /* Given x (not necessarily integral), and bid as output by zidealstarinit,
2840 : * compute the vector of components on the generators bid[2].
2841 : * Assume (x,bid) = 1 and sgn is either NULL or nfsign_arch(x, bid) */
2842 : GEN
2843 301909 : ideallog_sgn(GEN nf, GEN x, GEN sgn, GEN bid)
2844 : {
2845 : zlog_S S;
2846 301909 : nf = checknf(nf); checkbid(bid);
2847 301902 : init_zlog(&S, bid);
2848 301902 : if (sgn && typ(x) == t_VEC) /* vector of elements and signatures */
2849 : {
2850 21350 : long i, l = lg(x);
2851 21350 : GEN y = cgetg(l, t_MAT);
2852 21350 : for (i = 1; i < l; i++) gel(y,i) = ideallog_i(nf, gel(x,i), gel(sgn,i), &S);
2853 21350 : return y;
2854 : }
2855 280552 : return ideallog_i(nf, x, sgn, &S);
2856 : }
2857 : GEN
2858 287230 : ideallog(GEN nf, GEN x, GEN bid)
2859 : {
2860 287230 : if (!nf) return Zideallog(bid, x);
2861 280559 : return ideallog_sgn(nf, x, NULL, bid);
2862 : }
2863 :
2864 : /*************************************************************************/
2865 : /** **/
2866 : /** JOIN BID STRUCTURES, IDEAL LISTS **/
2867 : /** **/
2868 : /*************************************************************************/
2869 : /* bid1, bid2: for coprime modules m1 and m2 (without arch. part).
2870 : * Output: bid for m1 m2 */
2871 : static GEN
2872 476 : join_bid(GEN nf, GEN bid1, GEN bid2)
2873 : {
2874 476 : pari_sp av = avma;
2875 : long nbgen, l1,l2;
2876 : GEN I1,I2, G1,G2, sprk1,sprk2, cyc1,cyc2, sarch;
2877 476 : GEN sprk, fa,fa2, U, cyc, y, u1 = NULL, x, gen;
2878 :
2879 476 : I1 = bid_get_ideal(bid1);
2880 476 : I2 = bid_get_ideal(bid2);
2881 476 : if (gequal1(gcoeff(I1,1,1))) return bid2; /* frequent trivial case */
2882 259 : G1 = bid_get_grp(bid1);
2883 259 : G2 = bid_get_grp(bid2);
2884 259 : x = idealmul(nf, I1,I2);
2885 259 : fa = famat_mul_shallow(bid_get_fact(bid1), bid_get_fact(bid2));
2886 259 : fa2= famat_mul_shallow(bid_get_fact2(bid1), bid_get_fact2(bid2));
2887 259 : sprk1 = bid_get_sprk(bid1);
2888 259 : sprk2 = bid_get_sprk(bid2);
2889 259 : sprk = shallowconcat(sprk1, sprk2);
2890 :
2891 259 : cyc1 = abgrp_get_cyc(G1); l1 = lg(cyc1);
2892 259 : cyc2 = abgrp_get_cyc(G2); l2 = lg(cyc2);
2893 259 : gen = (lg(G1)>3 && lg(G2)>3)? gen_1: NULL;
2894 259 : nbgen = l1+l2-2;
2895 259 : cyc = ZV_snf_group(shallowconcat(cyc1,cyc2), &U, gen? &u1: NULL);
2896 259 : if (nbgen)
2897 : {
2898 259 : GEN U1 = bid_get_U(bid1), U2 = bid_get_U(bid2);
2899 259 : U1 = l1==1? const_vec(lg(sprk1), cgetg(1,t_MAT))
2900 259 : : ZM_ZMV_mul(vecslice(U, 1, l1-1), U1);
2901 259 : U2 = l2==1? const_vec(lg(sprk2), cgetg(1,t_MAT))
2902 259 : : ZM_ZMV_mul(vecslice(U, l1, nbgen), U2);
2903 259 : U = shallowconcat(U1, U2);
2904 : }
2905 : else
2906 0 : U = const_vec(lg(sprk), cgetg(1,t_MAT));
2907 :
2908 259 : if (gen)
2909 : {
2910 259 : GEN uv = zkchinese1init2(nf, I2, I1, x);
2911 518 : gen = shallowconcat(zkVchinese1(gel(uv,1), abgrp_get_gen(G1)),
2912 259 : zkVchinese1(gel(uv,2), abgrp_get_gen(G2)));
2913 : }
2914 259 : sarch = bid_get_sarch(bid1); /* trivial */
2915 259 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
2916 259 : x = mkvec2(x, bid_get_arch(bid1));
2917 259 : y = mkvec5(x, y, mkvec2(fa, fa2), mkvec2(sprk, sarch), U);
2918 259 : return gerepilecopy(av,y);
2919 : }
2920 :
2921 : typedef struct _ideal_data {
2922 : GEN nf, emb, L, pr, prL, sgnU, archp;
2923 : } ideal_data;
2924 :
2925 : /* z <- ( z | f(v[i])_{i=1..#v} ) */
2926 : static void
2927 86065 : concat_join(GEN *pz, GEN v, GEN (*f)(ideal_data*,GEN), ideal_data *data)
2928 : {
2929 86065 : long i, nz, lv = lg(v);
2930 : GEN z, Z;
2931 172130 : if (lv == 1) return;
2932 38143 : z = *pz; nz = lg(z)-1;
2933 38143 : *pz = Z = cgetg(lv + nz, typ(z));
2934 38143 : for (i = 1; i <=nz; i++) gel(Z,i) = gel(z,i);
2935 38143 : Z += nz;
2936 38143 : for (i = 1; i < lv; i++) gel(Z,i) = f(data, gel(v,i));
2937 : }
2938 : static GEN
2939 476 : join_idealinit(ideal_data *D, GEN x)
2940 476 : { return join_bid(D->nf, x, D->prL); }
2941 : static GEN
2942 47698 : join_ideal(ideal_data *D, GEN x)
2943 47698 : { return idealmulpowprime(D->nf, x, D->pr, D->L); }
2944 : static GEN
2945 455 : join_unit(ideal_data *D, GEN x)
2946 : {
2947 455 : GEN bid = join_idealinit(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
2948 455 : if (lg(u) != 1) v = shallowconcat(u, v);
2949 455 : return mkvec2(bid, v);
2950 : }
2951 :
2952 : /* flag & nf_GEN : generators, otherwise no
2953 : * flag &2 : units, otherwise no
2954 : * flag &4 : ideals in HNF, otherwise bid
2955 : * flag &8 : omit ideals which cannot be conductors (pr^1 with Npr=2) */
2956 : static GEN
2957 3192 : Ideallist(GEN bnf, ulong bound, long flag)
2958 : {
2959 3192 : const long cond = flag & 8;
2960 3192 : const long do_units = flag & 2, big_id = !(flag & 4);
2961 3192 : const long istar_flag = (flag & nf_GEN) | nf_INIT;
2962 3192 : pari_sp av, av0 = avma;
2963 : long i, j;
2964 3192 : GEN nf, z, p, fa, id, BOUND, U, empty = cgetg(1,t_VEC);
2965 : forprime_t S;
2966 : ideal_data ID;
2967 3192 : GEN (*join_z)(ideal_data*, GEN) =
2968 : do_units? &join_unit
2969 3192 : : (big_id? &join_idealinit: &join_ideal);
2970 :
2971 3192 : nf = checknf(bnf);
2972 3192 : if ((long)bound <= 0) return empty;
2973 3192 : id = matid(nf_get_degree(nf));
2974 3192 : if (big_id) id = Idealstar(nf,id, istar_flag);
2975 :
2976 : /* z[i] will contain all "objects" of norm i. Depending on flag, this means
2977 : * an ideal, a bid, or a couple [bid, log(units)]. Such objects are stored
2978 : * in vectors, computed one primary component at a time; join_z
2979 : * reconstructs the global object */
2980 3192 : BOUND = utoipos(bound);
2981 3192 : z = cgetg(bound+1,t_VEC);
2982 3192 : if (do_units) {
2983 14 : U = bnf_build_units(bnf);
2984 14 : gel(z,1) = mkvec( mkvec2(id, cgetg(1,t_VEC)) );
2985 : } else {
2986 3178 : U = NULL; /* -Wall */
2987 3178 : gel(z,1) = mkvec(id);
2988 : }
2989 3192 : for (i=2; i<=(long)bound; i++) gel(z,i) = empty;
2990 3192 : ID.nf = nf;
2991 :
2992 3192 : p = cgetipos(3);
2993 3192 : u_forprime_init(&S, 2, bound);
2994 3192 : av = avma;
2995 19600 : while ((p[2] = u_forprime_next(&S)))
2996 : {
2997 13216 : if (DEBUGLEVEL>1) { err_printf("%ld ",p[2]); err_flush(); }
2998 13216 : fa = idealprimedec_limit_norm(nf, p, BOUND);
2999 26859 : for (j=1; j<lg(fa); j++)
3000 : {
3001 13643 : GEN pr = gel(fa,j), z2 = leafcopy(z);
3002 13643 : ulong Q, q = upr_norm(pr);
3003 13643 : long l = (cond && q == 2)? 2: 1;
3004 :
3005 13643 : ID.pr = ID.prL = pr;
3006 33775 : for (Q = q; Q <= bound; l++, Q *= q) /* add pr^l */
3007 : {
3008 : ulong iQ;
3009 20132 : ID.L = utoipos(l);
3010 20132 : if (big_id) {
3011 217 : ID.prL = Idealstarprk(nf, pr, l, istar_flag);
3012 217 : if (do_units)
3013 : {
3014 196 : GEN sprk = bid_get_sprk(ID.prL);
3015 392 : ID.emb = lg(sprk) == 1? cgetg(1,t_VEC)
3016 196 : : vzlog_pr(nf, U, gel(sprk,1));
3017 : }
3018 : }
3019 106197 : for (iQ = Q,i = 1; iQ <= bound; iQ += Q,i++)
3020 86065 : concat_join(&gel(z,iQ), gel(z2,i), join_z, &ID);
3021 : }
3022 : }
3023 13216 : if (gc_needed(av,1))
3024 : {
3025 0 : if(DEBUGMEM>1) pari_warn(warnmem,"Ideallist");
3026 0 : z = gerepilecopy(av, z);
3027 : }
3028 : }
3029 3192 : return gerepilecopy(av0, z);
3030 : }
3031 : GEN
3032 350 : ideallist0(GEN bnf,long bound, long flag) {
3033 350 : if (flag<0 || flag>15) pari_err_FLAG("ideallist");
3034 350 : return Ideallist(bnf,bound,flag);
3035 : }
3036 : GEN
3037 2842 : ideallist(GEN bnf,long bound) { return Ideallist(bnf,bound,4); }
3038 :
3039 : /* bid = for module m (without arch. part), arch = archimedean part.
3040 : * Output: bid for [m,arch] */
3041 : static GEN
3042 56 : join_bid_arch(GEN nf, GEN bid, GEN archp)
3043 : {
3044 56 : pari_sp av = avma;
3045 : GEN G, U;
3046 56 : GEN sprk, cyc, y, u1 = NULL, x, sarch, gen;
3047 :
3048 56 : checkbid(bid);
3049 56 : G = bid_get_grp(bid);
3050 56 : x = bid_get_ideal(bid);
3051 56 : sarch = nfarchstar(nf, bid_get_ideal(bid), archp);
3052 56 : sprk = bid_get_sprk(bid);
3053 :
3054 56 : gen = (lg(G)>3)? gel(G,3): NULL;
3055 56 : cyc = diagonal_shallow(shallowconcat(gel(G,2), sarch_get_cyc(sarch)));
3056 56 : cyc = ZM_snf_group(cyc, &U, gen? &u1: NULL);
3057 56 : y = bid_grp(nf, u1, cyc, gen, x, sarch);
3058 56 : U = split_U(U, sprk);
3059 56 : y = mkvec5(mkvec2(x, archp), y, gel(bid,3), mkvec2(sprk, sarch), U);
3060 56 : return gerepilecopy(av,y);
3061 : }
3062 : static GEN
3063 56 : join_arch(ideal_data *D, GEN x) {
3064 56 : return join_bid_arch(D->nf, x, D->archp);
3065 : }
3066 : static GEN
3067 28 : join_archunit(ideal_data *D, GEN x) {
3068 28 : GEN bid = join_arch(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
3069 28 : if (lg(u) != 1) v = shallowconcat(u, v);
3070 28 : return mkvec2(bid, v);
3071 : }
3072 :
3073 : /* L from ideallist, add archimedean part */
3074 : GEN
3075 14 : ideallistarch(GEN bnf, GEN L, GEN arch)
3076 : {
3077 : pari_sp av;
3078 14 : long i, j, l = lg(L), lz;
3079 : GEN v, z, V;
3080 : ideal_data ID;
3081 : GEN (*join_z)(ideal_data*, GEN);
3082 :
3083 14 : if (typ(L) != t_VEC) pari_err_TYPE("ideallistarch",L);
3084 14 : if (l == 1) return cgetg(1,t_VEC);
3085 14 : z = gel(L,1);
3086 14 : if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
3087 14 : z = gel(z,1); /* either a bid or [bid,U] */
3088 14 : ID.nf = checknf(bnf);
3089 14 : ID.archp = vec01_to_indices(arch);
3090 14 : if (lg(z) == 3) { /* the latter: do units */
3091 7 : if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
3092 7 : ID.emb = zm_to_ZM( rowpermute(nfsign_units(bnf,NULL,1), ID.archp) );
3093 7 : join_z = &join_archunit;
3094 : } else
3095 7 : join_z = &join_arch;
3096 14 : av = avma; V = cgetg(l, t_VEC);
3097 70 : for (i = 1; i < l; i++)
3098 : {
3099 56 : z = gel(L,i); lz = lg(z);
3100 56 : gel(V,i) = v = cgetg(lz,t_VEC);
3101 56 : for (j=1; j<lz; j++) gel(v,j) = join_z(&ID, gel(z,j));
3102 : }
3103 14 : return gerepilecopy(av,V);
3104 : }
|
