1 /* Copyright (C) 2000 The PARI group.
2
3 This file is part of the PARI/GP package.
4
5 PARI/GP is free software; you can redistribute it and/or modify it under the
6 terms of the GNU General Public License as published by the Free Software
7 Foundation; either version 2 of the License, or (at your option) any later
8 version. It is distributed in the hope that it will be useful, but WITHOUT
9 ANY WARRANTY WHATSOEVER.
10
11 Check the License for details. You should have received a copy of it, along
12 with the package; see the file 'COPYING'. If not, write to the Free Software
13 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
14
15 /*******************************************************************/
16 /* */
17 /* RNF STRUCTURE AND OPERATIONS */
18 /* */
19 /*******************************************************************/
20 #include "pari.h"
21 #include "paripriv.h"
22
23 /* eq is an rnfeq; must return a t_POL */
24 GEN
eltreltoabs(GEN eq,GEN x)25 eltreltoabs(GEN eq, GEN x)
26 {
27 GEN Pabs = gel(eq,1), a = gel(eq,2), k = gel(eq,3), T = gel(eq,4), b, s;
28 long i, v = varn(Pabs);
29 pari_sp av = avma;
30
31 if (varncmp(gvar(x), v) > 0) x = scalarpol(x,v);
32 x = RgX_nffix("eltreltoabs", T, x, 1);
33 /* Mod(X - k a, Pabs(X)) is a root of the relative polynomial */
34 if (signe(k))
35 x = RgXQX_translate(x, deg1pol_shallow(negi(k), gen_0, varn(T)), T);
36 b = pol_x(v);
37 s = gen_0;
38 for (i=lg(x)-1; i>1; i--)
39 {
40 GEN c = gel(x,i);
41 if (typ(c) == t_POL) c = RgX_RgXQ_eval(c, a, Pabs);
42 s = RgX_rem(gadd(c, gmul(b,s)), Pabs);
43 }
44 return gerepileupto(av, s);
45 }
46 GEN
rnfeltreltoabs(GEN rnf,GEN x)47 rnfeltreltoabs(GEN rnf,GEN x)
48 {
49 const char *f = "rnfeltreltoabs";
50 GEN pol;
51 checkrnf(rnf);
52 pol = rnf_get_polabs(rnf);
53 switch(typ(x))
54 {
55 case t_INT: return icopy(x);
56 case t_FRAC: return gcopy(x);
57 case t_POLMOD:
58 if (RgX_equal_var(gel(x,1), pol))
59 { /* already in 'abs' form, unless possibly if nf = Q */
60 if (rnf_get_nfdegree(rnf) == 1)
61 {
62 GEN y = gel(x,2);
63 pari_sp av = avma;
64 y = simplify_shallow(liftpol_shallow(y));
65 return gerepilecopy(av, mkpolmod(y, pol));
66 }
67 return gcopy(x);
68 }
69 x = polmod_nffix(f,rnf,x,0);
70 if (typ(x) == t_POLMOD) return rnfeltup(rnf,x);
71 retmkpolmod(eltreltoabs(rnf_get_map(rnf), x), ZX_copy(pol));
72 case t_POL:
73 if (varn(x) == rnf_get_nfvarn(rnf)) return rnfeltup(rnf,x);
74 retmkpolmod(eltreltoabs(rnf_get_map(rnf), x), ZX_copy(pol));
75 }
76 pari_err_TYPE(f,x);
77 return NULL;/*LCOV_EXCL_LINE*/
78 }
79
80 GEN
eltabstorel_lift(GEN rnfeq,GEN P)81 eltabstorel_lift(GEN rnfeq, GEN P)
82 {
83 GEN k, T = gel(rnfeq,4), R = gel(rnfeq,5);
84 if (is_scalar_t(typ(P))) return P;
85 k = gel(rnfeq,3);
86 P = lift_shallow(P);
87 if (signe(k))
88 P = RgXQX_translate(P, deg1pol_shallow(k, gen_0, varn(T)), T);
89 P = RgXQX_rem(P, R, T);
90 return QXQX_to_mod_shallow(P, T);
91 }
92 /* rnfeq = [pol,a,k,T,R], P a t_POL or scalar
93 * Return Mod(P(x + k Mod(y, T(y))), pol(x)) */
94 GEN
eltabstorel(GEN rnfeq,GEN P)95 eltabstorel(GEN rnfeq, GEN P)
96 {
97 GEN T = gel(rnfeq,4), R = gel(rnfeq,5);
98 return mkpolmod(eltabstorel_lift(rnfeq,P), QXQX_to_mod_shallow(R,T));
99 }
100 GEN
rnfeltabstorel(GEN rnf,GEN x)101 rnfeltabstorel(GEN rnf,GEN x)
102 {
103 const char *f = "rnfeltabstorel";
104 pari_sp av = avma;
105 GEN pol, T, P, NF;
106 checkrnf(rnf);
107 T = rnf_get_nfpol(rnf);
108 P = rnf_get_pol(rnf);
109 pol = rnf_get_polabs(rnf);
110 switch(typ(x))
111 {
112 case t_INT: return icopy(x);
113 case t_FRAC: return gcopy(x);
114 case t_POLMOD:
115 if (RgX_equal_var(P, gel(x,1)))
116 {
117 x = polmod_nffix(f, rnf, x, 0);
118 P = QXQX_to_mod_shallow(P,T);
119 return gerepilecopy(av, mkpolmod(x,P));
120 }
121 if (RgX_equal_var(T, gel(x,1))) { x = Rg_nffix(f, T, x, 0); goto END; }
122 if (!RgX_equal_var(pol, gel(x,1))) pari_err_MODULUS(f, gel(x,1),pol);
123 x = gel(x,2); break;
124 case t_POL: break;
125 case t_COL:
126 NF = obj_check(rnf, rnf_NFABS);
127 if (!NF) pari_err_TYPE("rnfeltabstorel, apply nfinit(rnf)",x);
128 x = nf_to_scalar_or_alg(NF,x); break;
129 default:
130 pari_err_TYPE(f,x);
131 return NULL;/*LCOV_EXCL_LINE*/
132 }
133 switch(typ(x))
134 {
135 case t_INT: return icopy(x);
136 case t_FRAC: return gcopy(x);
137 case t_POL: break;
138 default: pari_err_TYPE(f, x);
139 }
140 RgX_check_QX(x,f);
141 if (varn(x) != varn(pol))
142 {
143 if (varn(x) == varn(T)) { x = Rg_nffix(f,T,x,0); goto END; }
144 pari_err_VAR(f, x,pol);
145 }
146 switch(lg(x))
147 {
148 case 2: set_avma(av); return gen_0;
149 case 3: return gerepilecopy(av, gel(x,2));
150 }
151 END:
152 return gerepilecopy(av, eltabstorel(rnf_get_map(rnf), x));
153 }
154
155 /* x a t_VEC of rnf elements in 'alg' form (t_POL). Assume maximal rank or 0 */
156 static GEN
modulereltoabs(GEN rnf,GEN x)157 modulereltoabs(GEN rnf, GEN x)
158 {
159 GEN W=gel(x,1), I=gel(x,2), rnfeq = rnf_get_map(rnf), polabs = gel(rnfeq,1);
160 long i, j, k, m, N = lg(W)-1;
161 GEN zknf, dzknf, M;
162
163 if (!N) return cgetg(1, t_VEC);
164 zknf = rnf_get_nfzk(rnf);
165 dzknf = gel(zknf,1);
166 m = rnf_get_nfdegree(rnf);
167 M = cgetg(N*m+1, t_VEC);
168 for (k=i=1; i<=N; i++)
169 {
170 GEN c0, cid, w = gel(W,i), id = gel(I,i);
171
172 if (lg(id) == 1) continue; /* must be a t_MAT */
173 id = Q_primitive_part(id, &cid);
174 w = Q_primitive_part(eltreltoabs(rnfeq,w), &c0);
175 c0 = div_content(mul_content(c0,cid), dzknf);
176 if (typ(id) == t_INT)
177 for (j=1; j<=m; j++)
178 {
179 GEN z = RgX_rem(gmul(w, gel(zknf,j)), polabs);
180 if (c0) z = RgX_Rg_mul(z, c0);
181 gel(M,k++) = z;
182 }
183 else
184 for (j=1; j<=m; j++)
185 {
186 GEN c, z = Q_primitive_part(RgV_RgC_mul(zknf,gel(id,j)), &c);
187 z = RgX_rem(gmul(w, z), polabs);
188 c = mul_content(c, c0); if (c) z = RgX_Rg_mul(z, c);
189 gel(M,k++) = z;
190 }
191 }
192 setlg(M, k); return M;
193 }
194
195 /* Z-basis for absolute maximal order: [NF.pol, NF.zk] */
196 GEN
rnf_zkabs(GEN rnf)197 rnf_zkabs(GEN rnf)
198 {
199 GEN d, v, M = modulereltoabs(rnf, rnf_get_zk(rnf));
200 GEN T = rnf_get_polabs(rnf);
201 long n = degpol(T);
202 M = Q_remove_denom(M, &d); /* t_VEC of t_POL */
203 if (d)
204 {
205 M = RgXV_to_RgM(M,n);
206 M = ZM_hnfmodall(M, d, hnf_MODID|hnf_CENTER);
207 M = RgM_Rg_div(M, d);
208 }
209 else
210 M = matid(n);
211 v = rnf_get_ramified_primes(rnf);
212 if (lg(v) == 1)
213 {
214 GEN D = gel(rnf_get_disc(rnf),1);
215 if (!isint1(D)) pari_err_TYPE("rnf_zkabs (old style rnf)", rnf);
216 }
217 v = shallowconcat(nf_get_ramified_primes(rnf_get_nf(rnf)), v);
218 return mkvec3(T, RgM_to_RgXV(M, varn(T)), ZV_sort_uniq(v));
219 }
220
221 static GEN
mknfabs(GEN rnf,long prec)222 mknfabs(GEN rnf, long prec)
223 {
224 GEN NF;
225 if ((NF = obj_check(rnf,rnf_NFABS)))
226 { if (nf_get_prec(NF) < prec) NF = nfnewprec_shallow(NF,prec); }
227 else
228 NF = nfinit(rnf_zkabs(rnf), prec);
229 return NF;
230 }
231
232 static GEN
mkupdown(GEN rnf)233 mkupdown(GEN rnf)
234 {
235 GEN NF = obj_check(rnf, rnf_NFABS), M, zknf, dzknf;
236 long i, l;
237 zknf = rnf_get_nfzk(rnf);
238 dzknf = gel(zknf,1); if (gequal1(dzknf)) dzknf = NULL;
239 l = lg(zknf); M = cgetg(l, t_MAT);
240 gel(M,1) = vec_ei(nf_get_degree(NF), 1);
241 for (i = 2; i < l; i++)
242 {
243 GEN c = poltobasis(NF, gel(zknf,i));
244 if (dzknf) c = gdiv(c, dzknf);
245 gel(M,i) = c;
246 }
247 return Qevproj_init(M);
248 }
249 GEN
rnf_build_nfabs(GEN rnf,long prec)250 rnf_build_nfabs(GEN rnf, long prec)
251 {
252 GEN NF = obj_checkbuild_prec(rnf, rnf_NFABS, &mknfabs, &nf_get_prec, prec);
253 (void)obj_checkbuild(rnf, rnf_MAPS, &mkupdown);
254 return NF;
255 }
256
257 void
rnfcomplete(GEN rnf)258 rnfcomplete(GEN rnf)
259 { (void)rnf_build_nfabs(rnf, nf_get_prec(rnf_get_nf(rnf))); }
260
261 GEN
nf_nfzk(GEN nf,GEN rnfeq)262 nf_nfzk(GEN nf, GEN rnfeq)
263 {
264 GEN pol = gel(rnfeq,1), a = gel(rnfeq,2);
265 return Q_primpart(QXV_QXQ_eval(nf_get_zkprimpart(nf), a, pol));
266 }
267
268 static GEN
rnfdisc_get_T_i(GEN P,GEN * lim)269 rnfdisc_get_T_i(GEN P, GEN *lim)
270 {
271 *lim = NULL;
272 if (typ(P) == t_VEC && lg(P) == 3)
273 {
274 GEN L = gel(P,2);
275 long i, l;
276 *lim = L;
277 switch(typ(L))
278 {
279 case t_INT:
280 if (signe(L) <= 0) return NULL;
281 break;
282 case t_VEC: case t_COL:
283 l = lg(L);
284 for (i = 1; i < l; i++)
285 {
286 GEN p = gel(L,i);
287 if (typ(p) == t_INT)
288 { if (signe(p) <= 0) return NULL; }
289 else checkprid(p);
290 }
291 break;
292 default: return NULL;
293 }
294 P = gel(P,1);
295 }
296 return (typ(P) == t_POL)? P: NULL;
297 }
298 /* true nf */
299 GEN
rnfdisc_get_T(GEN nf,GEN P,GEN * lim)300 rnfdisc_get_T(GEN nf, GEN P, GEN *lim)
301 {
302 GEN T = rnfdisc_get_T_i(P, lim);
303 if (!T) pari_err_TYPE("rnfdisc",P);
304 return RgX_nffix("rnfdisc", nf_get_pol(nf), T, 1);
305 }
306
307 GEN
rnfpseudobasis(GEN nf,GEN pol)308 rnfpseudobasis(GEN nf, GEN pol)
309 {
310 pari_sp av = avma;
311 GEN D, z, lim;
312 nf = checknf(nf);
313 pol = rnfdisc_get_T(nf, pol, &lim);
314 z = rnfallbase(nf, pol, lim, NULL, &D, NULL, NULL);
315 return gerepilecopy(av, shallowconcat(z,D));
316 }
317
318 GEN
rnfinit0(GEN nf,GEN T,long flag)319 rnfinit0(GEN nf, GEN T, long flag)
320 {
321 pari_sp av = avma;
322 GEN lim, bas, D, f, B, DKP, rnfeq, rnf = obj_init(11, 2);
323 nf = checknf(nf);
324 T = rnfdisc_get_T(nf, T, &lim);
325 gel(rnf,11) = rnfeq = nf_rnfeq(nf,T);
326 gel(rnf,2) = nf_nfzk(nf, rnfeq);
327 bas = rnfallbase(nf, T, lim, rnf, &D, &f, &DKP);
328 B = matbasistoalg(nf,gel(bas,1));
329 gel(bas,1) = lift_if_rational( RgM_to_RgXV(B,varn(T)) );
330 gel(rnf,1) = T;
331 gel(rnf,3) = D;
332 gel(rnf,4) = f;
333 gel(rnf,5) = DKP;
334 gel(rnf,6) = cgetg(1, t_VEC); /* dummy */
335 gel(rnf,7) = bas;
336 gel(rnf,8) = lift_if_rational( RgM_inv_upper(B) );
337 gel(rnf,9) = typ(f) == t_INT? powiu(f, nf_get_degree(nf))
338 : RgM_det_triangular(f);
339 gel(rnf,10)= nf;
340 rnf = gerepilecopy(av, rnf);
341 if (flag) rnfcomplete(rnf);
342 return rnf;
343 }
344 GEN
rnfinit(GEN nf,GEN T)345 rnfinit(GEN nf, GEN T) { return rnfinit0(nf,T,0); }
346
347 GEN
rnfeltup0(GEN rnf,GEN x,long flag)348 rnfeltup0(GEN rnf, GEN x, long flag)
349 {
350 pari_sp av = avma;
351 GEN zknf, nf, NF, POL;
352 long tx = typ(x);
353 checkrnf(rnf);
354 if (flag) rnfcomplete(rnf);
355 NF = obj_check(rnf,rnf_NFABS);
356 POL = rnf_get_polabs(rnf);
357 if (tx == t_POLMOD && RgX_equal_var(gel(x,1), POL))
358 {
359 if (flag) x = nf_to_scalar_or_basis(NF,x);
360 return gerepilecopy(av, x);
361 }
362 nf = rnf_get_nf(rnf);
363 if (NF && tx == t_COL && lg(x)-1 == degpol(POL) && nf_get_degree(rnf) > 1)
364 {
365 x = flag? nf_to_scalar_or_basis(NF,x)
366 : mkpolmod(nf_to_scalar_or_alg(NF,x), POL);
367 return gerepilecopy(av, x);
368 }
369 if (NF)
370 {
371 GEN d, proj;
372 x = nf_to_scalar_or_basis(nf, x);
373 if (typ(x) != t_COL) return gerepilecopy(av, x);
374 proj = obj_check(rnf,rnf_MAPS);
375 x = Q_remove_denom(x,&d);
376 x = ZM_ZC_mul(gel(proj,1), x);
377 if (d) x = gdiv(x,d);
378 if (!flag) x = basistoalg(NF,x);
379 }
380 else
381 {
382 zknf = rnf_get_nfzk(rnf);
383 x = nfeltup(nf, x, zknf);
384 if (typ(x) == t_POL) x = mkpolmod(x, POL);
385 }
386 return gerepilecopy(av, x);
387 }
388 GEN
rnfeltup(GEN rnf,GEN x)389 rnfeltup(GEN rnf, GEN x) { return rnfeltup0(rnf,x,0); }
390
391 GEN
nfeltup(GEN nf,GEN x,GEN zknf)392 nfeltup(GEN nf, GEN x, GEN zknf)
393 {
394 GEN c, dzknf = gel(zknf,1);
395 x = nf_to_scalar_or_basis(nf, x);
396 if (typ(x) != t_COL) return x;
397 x = Q_primitive_part(x, &c);
398 if (!RgV_is_ZV(x)) pari_err_TYPE("rnfeltup", x);
399 if (gequal1(dzknf)) dzknf = NULL;
400 c = div_content(c, dzknf);
401 x = RgV_RgC_mul(zknf, x); if (c) x = RgX_Rg_mul(x, c);
402 return x;
403 }
404
405 static void
fail(const char * f,GEN x)406 fail(const char *f, GEN x)
407 { pari_err_DOMAIN(f,"element","not in", strtoGENstr("the base field"),x); }
408 /* x t_COL of length degabs */
409 static GEN
eltdown(GEN rnf,GEN x,long flag)410 eltdown(GEN rnf, GEN x, long flag)
411 {
412 GEN z,y, d, proj = obj_check(rnf,rnf_MAPS);
413 GEN M= gel(proj,1), iM=gel(proj,2), diM=gel(proj,3), perm=gel(proj,4);
414 x = Q_remove_denom(x,&d);
415 if (!RgV_is_ZV(x)) pari_err_TYPE("rnfeltdown", x);
416 y = ZM_ZC_mul(iM, vecpermute(x, perm));
417 z = ZM_ZC_mul(M,y);
418 if (!isint1(diM)) z = ZC_Z_mul(z,diM);
419 if (!ZV_equal(z,x)) fail("rnfeltdown",x);
420
421 d = mul_denom(d, diM);
422 if (d) y = gdiv(y,d);
423 if (!flag) y = basistoalg(rnf_get_nf(rnf), y);
424 return y;
425 }
426 GEN
rnfeltdown0(GEN rnf,GEN x,long flag)427 rnfeltdown0(GEN rnf, GEN x, long flag)
428 {
429 const char *f = "rnfeltdown";
430 pari_sp av = avma;
431 GEN z, T, NF, nf;
432 long v;
433
434 checkrnf(rnf);
435 NF = obj_check(rnf,rnf_NFABS);
436 nf = rnf_get_nf(rnf);
437 T = nf_get_pol(nf);
438 v = varn(T);
439 switch(typ(x))
440 { /* directly belonging to base field ? */
441 case t_INT: return icopy(x);
442 case t_FRAC:return gcopy(x);
443 case t_POLMOD:
444 if (RgX_equal_var(gel(x,1), rnf_get_polabs(rnf)))
445 {
446 if (degpol(T) == 1)
447 {
448 x = simplify_shallow(liftpol_shallow(gel(x,2)));
449 if (typ(x) != t_POL) return gerepilecopy(av,x);
450 }
451 break;
452 }
453 x = polmod_nffix(f,rnf,x,0);
454 /* x was defined mod the relative polynomial & non constant => fail */
455 if (typ(x) == t_POL) fail(f,x);
456 if (flag) x = nf_to_scalar_or_basis(nf,x);
457 return gerepilecopy(av, x);
458
459 case t_POL:
460 if (varn(x) != v) break;
461 x = Rg_nffix(f,T,x,0);
462 if (flag) x = nf_to_scalar_or_basis(nf,x);
463 return gerepilecopy(av, x);
464 case t_COL:
465 {
466 long n = lg(x)-1;
467 if (n == degpol(T) && RgV_is_QV(x))
468 {
469 if (RgV_isscalar(x)) return gcopy(gel(x,1));
470 if (!flag) return gcopy(x);
471 return basistoalg(nf,x);
472 }
473 if (NF) break;
474 }
475 default: pari_err_TYPE(f, x);
476 }
477 /* x defined mod the absolute equation */
478 if (NF)
479 {
480 x = nf_to_scalar_or_basis(NF, x);
481 if (typ(x) == t_COL) x = eltdown(rnf,x,flag);
482 return gerepilecopy(av, x);
483 }
484 z = rnfeltabstorel(rnf,x);
485 switch(typ(z))
486 {
487 case t_INT:
488 case t_FRAC: return z;
489 }
490 /* typ(z) = t_POLMOD, varn of both components is rnf_get_varn(rnf) */
491 z = gel(z,2);
492 if (typ(z) == t_POL)
493 {
494 if (lg(z) != 3) fail(f,x);
495 z = gel(z,2);
496 }
497 return gerepilecopy(av, z);
498 }
499 GEN
rnfeltdown(GEN rnf,GEN x)500 rnfeltdown(GEN rnf, GEN x) { return rnfeltdown0(rnf,x,0); }
501
502 /* vector of rnf elt -> matrix of nf elts */
503 static GEN
rnfV_to_nfM(GEN rnf,GEN x)504 rnfV_to_nfM(GEN rnf, GEN x)
505 {
506 long i, l = lg(x);
507 GEN y = cgetg(l, t_MAT);
508 for (i = 1; i < l; i++) gel(y,i) = rnfalgtobasis(rnf,gel(x,i));
509 return y;
510 }
511
512 static GEN
rnfprincipaltohnf(GEN rnf,GEN x)513 rnfprincipaltohnf(GEN rnf,GEN x)
514 {
515 pari_sp av = avma;
516 GEN bas = rnf_get_zk(rnf), nf = rnf_get_nf(rnf);
517 x = rnfbasistoalg(rnf,x);
518 x = gmul(x, gmodulo(gel(bas,1), rnf_get_pol(rnf)));
519 return gerepileupto(av, nfhnf(nf, mkvec2(rnfV_to_nfM(rnf,x), gel(bas,2))));
520 }
521
522 /* pseudo-basis for the 0 ideal */
523 static GEN
rnfideal0(void)524 rnfideal0(void) { retmkvec2(cgetg(1,t_MAT),cgetg(1,t_VEC)); }
525
526 GEN
rnfidealhnf(GEN rnf,GEN x)527 rnfidealhnf(GEN rnf, GEN x)
528 {
529 GEN z, nf, bas;
530
531 checkrnf(rnf); nf = rnf_get_nf(rnf);
532 switch(typ(x))
533 {
534 case t_INT: case t_FRAC:
535 if (isintzero(x)) return rnfideal0();
536 bas = rnf_get_zk(rnf); z = cgetg(3,t_VEC);
537 gel(z,1) = matid(rnf_get_degree(rnf));
538 gel(z,2) = gmul(x, gel(bas,2)); return z;
539
540 case t_VEC:
541 if (lg(x) == 3 && typ(gel(x,1)) == t_MAT) return nfhnf(nf, x);
542 case t_MAT:
543 return rnfidealabstorel(rnf, x);
544
545 case t_POLMOD: case t_POL: case t_COL:
546 return rnfprincipaltohnf(rnf,x);
547 }
548 pari_err_TYPE("rnfidealhnf",x);
549 return NULL; /* LCOV_EXCL_LINE */
550 }
551
552 static GEN
prodidnorm(GEN nf,GEN I)553 prodidnorm(GEN nf, GEN I)
554 {
555 long i, l = lg(I);
556 GEN z;
557 if (l == 1) return gen_1;
558 z = idealnorm(nf, gel(I,1));
559 for (i=2; i<l; i++) z = gmul(z, idealnorm(nf, gel(I,i)));
560 return z;
561 }
562
563 GEN
rnfidealnormrel(GEN rnf,GEN id)564 rnfidealnormrel(GEN rnf, GEN id)
565 {
566 pari_sp av = avma;
567 GEN nf, z = gel(rnfidealhnf(rnf,id), 2);
568 if (lg(z) == 1) return cgetg(1, t_MAT);
569 nf = rnf_get_nf(rnf); z = idealprod(nf, z);
570 return gerepileupto(av, idealmul(nf,z, rnf_get_index(rnf)));
571 }
572
573 GEN
rnfidealnormabs(GEN rnf,GEN id)574 rnfidealnormabs(GEN rnf, GEN id)
575 {
576 pari_sp av = avma;
577 GEN nf, z = gel(rnfidealhnf(rnf,id), 2);
578 if (lg(z) == 1) return gen_0;
579 nf = rnf_get_nf(rnf); z = prodidnorm(nf, z);
580 return gerepileupto(av, gmul(z, gel(rnf,9)));
581 }
582
583 static GEN
rnfidealreltoabs_i(GEN rnf,GEN x)584 rnfidealreltoabs_i(GEN rnf, GEN x)
585 {
586 long i, l;
587 GEN w;
588 x = rnfidealhnf(rnf,x);
589 w = gel(x,1); l = lg(w); settyp(w, t_VEC);
590 for (i=1; i<l; i++) gel(w,i) = lift_shallow( rnfbasistoalg(rnf, gel(w,i)) );
591 return modulereltoabs(rnf, x);
592 }
593 GEN
rnfidealreltoabs(GEN rnf,GEN x)594 rnfidealreltoabs(GEN rnf, GEN x)
595 {
596 pari_sp av = avma;
597 return gerepilecopy(av, rnfidealreltoabs_i(rnf,x));
598 }
599 GEN
rnfidealreltoabs0(GEN rnf,GEN x,long flag)600 rnfidealreltoabs0(GEN rnf, GEN x, long flag)
601 {
602 pari_sp av = avma;
603 long i, l;
604 GEN NF;
605
606 x = rnfidealreltoabs_i(rnf, x);
607 if (!flag) return gerepilecopy(av,x);
608 rnfcomplete(rnf);
609 NF = obj_check(rnf,rnf_NFABS);
610 l = lg(x); settyp(x, t_MAT);
611 for (i=1; i<l; i++) gel(x,i) = algtobasis(NF, gel(x,i));
612 return gerepileupto(av, idealhnf(NF,x));
613 }
614
615 GEN
rnfidealabstorel(GEN rnf,GEN x)616 rnfidealabstorel(GEN rnf, GEN x)
617 {
618 long n, N, j, tx = typ(x);
619 pari_sp av = avma;
620 GEN A, I, invbas;
621
622 checkrnf(rnf);
623 invbas = rnf_get_invzk(rnf);
624 if (tx != t_VEC && tx != t_MAT) pari_err_TYPE("rnfidealabstorel",x);
625 N = lg(x)-1;
626 if (N != rnf_get_absdegree(rnf))
627 {
628 if (!N) return rnfideal0();
629 pari_err_DIM("rnfidealabstorel");
630 }
631 n = rnf_get_degree(rnf);
632 A = cgetg(N+1,t_MAT);
633 I = cgetg(N+1,t_VEC);
634 for (j=1; j<=N; j++)
635 {
636 GEN t = lift_shallow( rnfeltabstorel(rnf, gel(x,j)) );
637 if (typ(t) == t_POL)
638 t = RgM_RgX_mul(invbas, t);
639 else
640 t = scalarcol_shallow(t, n);
641 gel(A,j) = t;
642 gel(I,j) = gen_1;
643 }
644 return gerepileupto(av, nfhnf(rnf_get_nf(rnf), mkvec2(A,I)));
645 }
646
647 GEN
rnfidealdown(GEN rnf,GEN x)648 rnfidealdown(GEN rnf,GEN x)
649 {
650 pari_sp av = avma;
651 GEN I;
652 if (typ(x) == t_MAT)
653 {
654 GEN d;
655 x = Q_remove_denom(x,&d);
656 if (RgM_is_ZM(x))
657 {
658 GEN NF = obj_check(rnf,rnf_NFABS);
659 if (NF)
660 {
661 GEN z, proj = obj_check(rnf,rnf_MAPS), ZK = gel(proj,1);
662 long i, lz, l;
663 x = idealhnf(NF,x);
664 if (lg(x) == 1) { set_avma(av); return cgetg(1,t_MAT); }
665 z = ZM_lll(shallowconcat(ZK,x), 0.99, LLL_KER);
666 lz = lg(z); l = lg(ZK);
667 for (i = 1; i < lz; i++) setlg(gel(z,i), l);
668 z = ZM_hnfmodid(z, gcoeff(x,1,1));
669 if (d) z = gdiv(z,d);
670 return gerepileupto(av, z);
671 }
672 }
673 }
674 x = rnfidealhnf(rnf,x); I = gel(x,2);
675 if (lg(I) == 1) { set_avma(av); return cgetg(1,t_MAT); }
676 return gerepilecopy(av, gel(I,1));
677 }
678
679 /* lift ideal x to the relative extension, returns a Z-basis */
680 GEN
rnfidealup(GEN rnf,GEN x)681 rnfidealup(GEN rnf,GEN x)
682 {
683 pari_sp av = avma;
684 long i, n;
685 GEN nf, bas, bas2, I, x2, dx;
686
687 checkrnf(rnf); nf = rnf_get_nf(rnf);
688 n = rnf_get_degree(rnf);
689 bas = rnf_get_zk(rnf); bas2 = gel(bas,2);
690
691 (void)idealtyp(&x, &I); /* I is junk */
692 x = Q_remove_denom(x, &dx);
693 x2 = idealtwoelt(nf,x);
694 I = cgetg(n+1,t_VEC);
695 for (i=1; i<=n; i++)
696 {
697 GEN c = gel(bas2,i), d;
698 if (typ(c) == t_MAT)
699 {
700 c = Q_remove_denom(c,&d);
701 d = mul_denom(d, dx);
702 c = idealHNF_mul(nf,c,x2);
703 }
704 else
705 {
706 c = idealmul(nf,c,x);
707 d = dx;
708 }
709 if (d) c = gdiv(c,d);
710 gel(I,i) = c;
711 }
712 return gerepilecopy(av, modulereltoabs(rnf, mkvec2(gel(bas,1), I)));
713 }
714 GEN
rnfidealup0(GEN rnf,GEN x,long flag)715 rnfidealup0(GEN rnf,GEN x, long flag)
716 {
717 pari_sp av = avma;
718 GEN NF, nf, proj, d, x2;
719
720 if (!flag) return rnfidealup(rnf,x);
721 checkrnf(rnf); nf = rnf_get_nf(rnf);
722 rnfcomplete(rnf);
723 proj = obj_check(rnf,rnf_MAPS);
724 NF = obj_check(rnf,rnf_NFABS);
725
726 (void)idealtyp(&x, &d); /* d is junk */
727 x2 = idealtwoelt(nf,x);
728 x2 = Q_remove_denom(x2,&d);
729 if (typ(gel(x2,2)) == t_COL) gel(x2,2) = ZM_ZC_mul(gel(proj,1),gel(x2,2));
730 x2 = idealhnf_two(NF, x2);
731 if (d) x2 = gdiv(x2,d);
732 return gerepileupto(av, x2);
733 }
734
735 /* x a relative HNF => vector of 2 generators (relative polmods) */
736 GEN
rnfidealtwoelement(GEN rnf,GEN x)737 rnfidealtwoelement(GEN rnf, GEN x)
738 {
739 pari_sp av = avma;
740 GEN y, cy, z, NF;
741
742 y = rnfidealreltoabs_i(rnf,x);
743 rnfcomplete(rnf);
744 NF = obj_check(rnf,rnf_NFABS);
745 y = matalgtobasis(NF, y); settyp(y, t_MAT);
746 y = Q_primitive_part(y, &cy);
747 y = ZM_hnf(y);
748 if (lg(y) == 1) { set_avma(av); return mkvec2(gen_0, gen_0); }
749 y = idealtwoelt(NF, y);
750 if (cy) y = RgV_Rg_mul(y, cy);
751 z = gel(y,2);
752 if (typ(z) == t_COL) z = rnfeltabstorel(rnf, nf_to_scalar_or_alg(NF, z));
753 return gerepilecopy(av, mkvec2(gel(y,1), z));
754 }
755
756 GEN
rnfidealmul(GEN rnf,GEN x,GEN y)757 rnfidealmul(GEN rnf,GEN x,GEN y)
758 {
759 pari_sp av = avma;
760 GEN nf, z, x1, x2, p1, p2, bas;
761
762 y = rnfidealtwoelement(rnf,y);
763 if (isintzero(gel(y,1))) { set_avma(av); return rnfideal0(); }
764 nf = rnf_get_nf(rnf);
765 bas = rnf_get_zk(rnf);
766 x = rnfidealhnf(rnf,x);
767 x1 = gmodulo(gmul(gel(bas,1), matbasistoalg(nf,gel(x,1))), rnf_get_pol(rnf));
768 x2 = gel(x,2);
769 p1 = gmul(gel(y,1), gel(x,1));
770 p2 = rnfV_to_nfM(rnf, gmul(gel(y,2), x1));
771 z = mkvec2(shallowconcat(p1, p2), shallowconcat(x2, x2));
772 return gerepileupto(av, nfhnf(nf,z));
773 }
774
775 /* prK wrt NF ~ Q[x]/(polabs) */
776 static GEN
rnfidealprimedec_1(GEN rnf,GEN SL,GEN prK)777 rnfidealprimedec_1(GEN rnf, GEN SL, GEN prK)
778 {
779 GEN v, piL, piK = pr_get_gen(prK);
780 long i, c, l;
781 if (pr_is_inert(prK)) return SL;
782 piL = rnfeltup0(rnf, piK, 1);
783 v = cgetg_copy(SL, &l);
784 for (i = c = 1; i < l; i++)
785 {
786 GEN P = gel(SL,i);
787 if (ZC_prdvd(piL, P)) gel(v,c++) = P;
788 }
789 setlg(v, c); return v;
790 }
791 GEN
rnfidealprimedec(GEN rnf,GEN pr)792 rnfidealprimedec(GEN rnf, GEN pr)
793 {
794 pari_sp av = avma;
795 GEN p, z, NF, nf, SL;
796 checkrnf(rnf);
797 rnfcomplete(rnf);
798 NF = obj_check(rnf,rnf_NFABS);
799 nf = rnf_get_nf(rnf);
800 if (typ(pr) == t_INT) { p = pr; pr = NULL; }
801 else { checkprid(pr); p = pr_get_p(pr); }
802 SL = idealprimedec(NF, p);
803 if (pr) z = rnfidealprimedec_1(rnf, SL, pr);
804 else
805 {
806 GEN vK = idealprimedec(nf, p), vL;
807 long l = lg(vK), i;
808 vL = cgetg(l, t_VEC);
809 for (i = 1; i < l; i++) gel(vL,i) = rnfidealprimedec_1(rnf, SL, gel(vK,i));
810 z = mkvec2(vK, vL);
811 }
812 return gerepilecopy(av, z);
813 }
814
815 GEN
rnfidealfactor(GEN rnf,GEN x)816 rnfidealfactor(GEN rnf, GEN x)
817 {
818 pari_sp av = avma;
819 GEN NF;
820 checkrnf(rnf);
821 rnfcomplete(rnf);
822 NF = obj_check(rnf,rnf_NFABS);
823 return gerepileupto(av, idealfactor(NF, rnfidealreltoabs0(rnf, x, 1)));
824 }
825
826 GEN
rnfequationall(GEN A,GEN B,long * pk,GEN * pLPRS)827 rnfequationall(GEN A, GEN B, long *pk, GEN *pLPRS)
828 {
829 long lA, lB;
830 GEN nf, C;
831
832 A = get_nfpol(A, &nf); lA = lg(A);
833 if (!nf) {
834 if (lA<=3) pari_err_CONSTPOL("rnfequation");
835 RgX_check_ZX(A,"rnfequation");
836 }
837 B = RgX_nffix("rnfequation", A,B,1); lB = lg(B);
838 if (lB<=3) pari_err_CONSTPOL("rnfequation");
839 B = Q_primpart(B);
840
841 if (!nfissquarefree(A,B))
842 pari_err_DOMAIN("rnfequation","issquarefree(B)","=",gen_0,B);
843
844 *pk = 0; C = ZX_ZXY_resultant_all(A, B, pk, pLPRS);
845 if (signe(leading_coeff(C)) < 0) C = ZX_neg(C);
846 *pk = -*pk; return Q_primpart(C);
847 }
848
849 GEN
rnfequation0(GEN A,GEN B,long flall)850 rnfequation0(GEN A, GEN B, long flall)
851 {
852 pari_sp av = avma;
853 GEN LPRS, C;
854 long k;
855
856 C = rnfequationall(A, B, &k, flall? &LPRS: NULL);
857 if (flall)
858 { /* a,b,c root of A,B,C = compositum, c = b + k a */
859 GEN a, mH0 = RgX_neg(gel(LPRS,1)), H1 = gel(LPRS,2);
860 a = QXQ_div(mH0, H1, C);
861 C = mkvec3(C, mkpolmod(a, C), stoi(k));
862 }
863 return gerepilecopy(av, C);
864 }
865 GEN
rnfequation(GEN nf,GEN pol)866 rnfequation(GEN nf, GEN pol) { return rnfequation0(nf,pol,0); }
867 GEN
rnfequation2(GEN nf,GEN pol)868 rnfequation2(GEN nf, GEN pol) { return rnfequation0(nf,pol,1); }
869 GEN
nf_rnfeq(GEN nf,GEN R)870 nf_rnfeq(GEN nf, GEN R)
871 {
872 GEN pol, a, k, junk, eq;
873 R = liftpol_shallow(R);
874 eq = rnfequation2(nf, R);
875 pol = gel(eq,1);
876 a = gel(eq,2); if (typ(a) == t_POLMOD) a = gel(a,2);
877 k = gel(eq,3);
878 return mkvec5(pol,a,k,get_nfpol(nf, &junk),R);
879 }
880 /* only allow abstorel */
881 GEN
nf_rnfeqsimple(GEN nf,GEN R)882 nf_rnfeqsimple(GEN nf, GEN R)
883 {
884 long sa;
885 GEN junk, pol;
886 R = liftpol_shallow(R);
887 pol = rnfequationall(nf, R, &sa, NULL);
888 return mkvec5(pol,gen_0/*dummy*/,stoi(sa),get_nfpol(nf, &junk),R);
889 }
890
891 /*******************************************************************/
892 /* */
893 /* RELATIVE LLL */
894 /* */
895 /*******************************************************************/
896 static GEN
nftau(long r1,GEN x)897 nftau(long r1, GEN x)
898 {
899 long i, l = lg(x);
900 GEN s = r1? gel(x,1): gmul2n(real_i(gel(x,1)),1);
901 for (i=2; i<=r1; i++) s = gadd(s, gel(x,i));
902 for ( ; i < l; i++) s = gadd(s, gmul2n(real_i(gel(x,i)),1));
903 return s;
904 }
905
906 static GEN
initmat(long l)907 initmat(long l)
908 {
909 GEN x = cgetg(l, t_MAT);
910 long i;
911 for (i = 1; i < l; i++) gel(x,i) = cgetg(l, t_COL);
912 return x;
913 }
914
915 static GEN
nftocomplex(GEN nf,GEN x)916 nftocomplex(GEN nf, GEN x)
917 {
918 GEN M = nf_get_M(nf);
919 x = nf_to_scalar_or_basis(nf,x);
920 if (typ(x) != t_COL) return const_col(nbrows(M), x);
921 return RgM_RgC_mul(M, x);
922 }
923 /* assume x a square t_MAT, return a t_VEC of embeddings of its columns */
924 static GEN
mattocomplex(GEN nf,GEN x)925 mattocomplex(GEN nf, GEN x)
926 {
927 long i,j, l = lg(x);
928 GEN v = cgetg(l, t_VEC);
929 for (j=1; j<l; j++)
930 {
931 GEN c = gel(x,j), b = cgetg(l, t_MAT);
932 for (i=1; i<l; i++) gel(b,i) = nftocomplex(nf, gel(c,i));
933 b = shallowtrans(b); settyp(b, t_COL);
934 gel(v,j) = b;
935 }
936 return v;
937 }
938
939 static GEN
nf_all_roots(GEN nf,GEN x,long prec)940 nf_all_roots(GEN nf, GEN x, long prec)
941 {
942 long i, j, l = lg(x), ru = lg(nf_get_roots(nf));
943 GEN y = cgetg(l, t_POL), v, z;
944
945 x = RgX_to_nfX(nf, x);
946 y[1] = x[1];
947 for (i=2; i<l; i++) gel(y,i) = nftocomplex(nf, gel(x,i));
948 i = gprecision(y); if (i && i <= 3) return NULL;
949
950 v = cgetg(ru, t_VEC);
951 z = cgetg(l, t_POL); z[1] = x[1];
952 for (i=1; i<ru; i++)
953 {
954 for (j = 2; j < l; j++) gel(z,j) = gmael(y,j,i);
955 gel(v,i) = cleanroots(z, prec);
956 }
957 return v;
958 }
959
960 static GEN
rnfscal(GEN m,GEN x,GEN y)961 rnfscal(GEN m, GEN x, GEN y)
962 {
963 long i, l = lg(m);
964 GEN z = cgetg(l, t_COL);
965 for (i = 1; i < l; i++)
966 gel(z,i) = gmul(conj_i(shallowtrans(gel(x,i))), gmul(gel(m,i), gel(y,i)));
967 return z;
968 }
969
970 /* x ideal in HNF */
971 static GEN
findmin(GEN nf,GEN x,GEN muf)972 findmin(GEN nf, GEN x, GEN muf)
973 {
974 pari_sp av = avma;
975 long e;
976 GEN cx, y, m, M = nf_get_M(nf);
977
978 x = Q_primitive_part(x, &cx);
979 if (gequal1(gcoeff(x,1,1))) y = M;
980 else
981 {
982 GEN G = nf_get_G(nf);
983 m = lllfp(RgM_mul(G,x), 0.75, 0);
984 if (typ(m) != t_MAT)
985 {
986 x = ZM_lll(x, 0.75, LLL_INPLACE);
987 m = lllfp(RgM_mul(G,x), 0.75, 0);
988 if (typ(m) != t_MAT) pari_err_PREC("rnflllgram");
989 }
990 x = ZM_mul(x, m);
991 y = RgM_mul(M, x);
992 }
993 m = RgM_solve_realimag(y, muf);
994 if (!m) return NULL; /* precision problem */
995 if (cx) m = RgC_Rg_div(m, cx);
996 m = grndtoi(m, &e);
997 if (e >= 0) return NULL; /* precision problem */
998 m = ZM_ZC_mul(x, m);
999 if (cx) m = ZC_Q_mul(m, cx);
1000 return gerepileupto(av, m);
1001 }
1002
1003 static int
RED(long k,long l,GEN U,GEN mu,GEN MC,GEN nf,GEN I,GEN * Ik_inv)1004 RED(long k, long l, GEN U, GEN mu, GEN MC, GEN nf, GEN I, GEN *Ik_inv)
1005 {
1006 GEN x, xc, ideal;
1007 long i;
1008
1009 if (!*Ik_inv) *Ik_inv = idealinv(nf, gel(I,k));
1010 ideal = idealmul(nf,gel(I,l), *Ik_inv);
1011 x = findmin(nf, ideal, gcoeff(mu,k,l));
1012 if (!x) return 0;
1013 if (gequal0(x)) return 1;
1014
1015 xc = nftocomplex(nf,x);
1016 gel(MC,k) = gsub(gel(MC,k), vecmul(xc,gel(MC,l)));
1017 gel(U,k) = gsub(gel(U,k), gmul(coltoalg(nf,x), gel(U,l)));
1018 gcoeff(mu,k,l) = gsub(gcoeff(mu,k,l), xc);
1019 for (i=1; i<l; i++)
1020 gcoeff(mu,k,i) = gsub(gcoeff(mu,k,i), vecmul(xc,gcoeff(mu,l,i)));
1021 return 1;
1022 }
1023
1024 static int
check_0(GEN B)1025 check_0(GEN B)
1026 {
1027 long i, l = lg(B);
1028 for (i = 1; i < l; i++)
1029 if (gsigne(gel(B,i)) <= 0) return 1;
1030 return 0;
1031 }
1032
1033 static int
do_SWAP(GEN I,GEN MC,GEN MCS,GEN h,GEN mu,GEN B,long kmax,long k,const long alpha,long r1)1034 do_SWAP(GEN I, GEN MC, GEN MCS, GEN h, GEN mu, GEN B, long kmax, long k,
1035 const long alpha, long r1)
1036 {
1037 GEN p1, p2, muf, mufc, Bf, temp;
1038 long i, j;
1039
1040 p1 = nftau(r1, gadd(gel(B,k),
1041 gmul(gnorml2(gcoeff(mu,k,k-1)), gel(B,k-1))));
1042 p2 = nftau(r1, gel(B,k-1));
1043 if (gcmp(gmulsg(alpha,p1), gmulsg(alpha-1,p2)) > 0) return 0;
1044
1045 swap(gel(MC,k-1),gel(MC,k));
1046 swap(gel(h,k-1), gel(h,k));
1047 swap(gel(I,k-1), gel(I,k));
1048 for (j=1; j<=k-2; j++) swap(gcoeff(mu,k-1,j),gcoeff(mu,k,j));
1049 muf = gcoeff(mu,k,k-1);
1050 mufc = conj_i(muf);
1051 Bf = gadd(gel(B,k), vecmul(real_i(vecmul(muf,mufc)), gel(B,k-1)));
1052 if (check_0(Bf)) return 1; /* precision problem */
1053
1054 p1 = vecdiv(gel(B,k-1),Bf);
1055 gcoeff(mu,k,k-1) = vecmul(mufc,p1);
1056 temp = gel(MCS,k-1);
1057 gel(MCS,k-1) = gadd(gel(MCS,k), vecmul(muf,gel(MCS,k-1)));
1058 gel(MCS,k) = gsub(vecmul(vecdiv(gel(B,k),Bf), temp),
1059 vecmul(gcoeff(mu,k,k-1), gel(MCS,k)));
1060 gel(B,k) = vecmul(gel(B,k),p1);
1061 gel(B,k-1) = Bf;
1062 for (i=k+1; i<=kmax; i++)
1063 {
1064 temp = gcoeff(mu,i,k);
1065 gcoeff(mu,i,k) = gsub(gcoeff(mu,i,k-1), vecmul(muf, gcoeff(mu,i,k)));
1066 gcoeff(mu,i,k-1) = gadd(temp, vecmul(gcoeff(mu,k,k-1),gcoeff(mu,i,k)));
1067 }
1068 return 1;
1069 }
1070
1071 static GEN
rel_T2(GEN nf,GEN pol,long lx,long prec)1072 rel_T2(GEN nf, GEN pol, long lx, long prec)
1073 {
1074 long ru, i, j, k, l;
1075 GEN T2, s, unro, roorder, powreorder;
1076
1077 roorder = nf_all_roots(nf, pol, prec);
1078 if (!roorder) return NULL;
1079 ru = lg(roorder);
1080 unro = cgetg(lx,t_COL); for (i=1; i<lx; i++) gel(unro,i) = gen_1;
1081 powreorder = cgetg(lx,t_MAT); gel(powreorder,1) = unro;
1082 T2 = cgetg(ru, t_VEC);
1083 for (i = 1; i < ru; i++)
1084 {
1085 GEN ro = gel(roorder,i);
1086 GEN m = initmat(lx);
1087 for (k=2; k<lx; k++)
1088 {
1089 GEN c = cgetg(lx, t_COL); gel(powreorder,k) = c;
1090 for (j=1; j < lx; j++)
1091 gel(c,j) = gmul(gel(ro,j), gmael(powreorder,k-1,j));
1092 }
1093 for (l = 1; l < lx; l++)
1094 for (k = 1; k <= l; k++)
1095 {
1096 s = gen_0;
1097 for (j = 1; j < lx; j++)
1098 s = gadd(s, gmul(conj_i(gmael(powreorder,k,j)),
1099 gmael(powreorder,l,j)));
1100 if (l == k)
1101 gcoeff(m, l, l) = real_i(s);
1102 else
1103 {
1104 gcoeff(m, k, l) = s;
1105 gcoeff(m, l, k) = conj_i(s);
1106 }
1107 }
1108 gel(T2,i) = m;
1109 }
1110 return T2;
1111 }
1112
1113 /* given a base field nf (e.g main variable y), a polynomial pol with
1114 * coefficients in nf (e.g main variable x), and an order as output
1115 * by rnfpseudobasis, outputs a reduced order. */
1116 GEN
rnflllgram(GEN nf,GEN pol,GEN order,long prec)1117 rnflllgram(GEN nf, GEN pol, GEN order,long prec)
1118 {
1119 pari_sp av = avma;
1120 long j, k, l, kmax, r1, lx, count = 0;
1121 GEN M, I, h, H, mth, MC, MPOL, MCS, B, mu;
1122 const long alpha = 10, MAX_COUNT = 4;
1123
1124 nf = checknf(nf); r1 = nf_get_r1(nf);
1125 check_ZKmodule(order, "rnflllgram");
1126 M = gel(order,1);
1127 I = gel(order,2); lx = lg(I);
1128 if (lx < 3) return gcopy(order);
1129 if (lx-1 != degpol(pol)) pari_err_DIM("rnflllgram");
1130 I = leafcopy(I);
1131 H = NULL;
1132 MPOL = matbasistoalg(nf, M);
1133 MCS = matid(lx-1); /* dummy for gerepile */
1134 PRECNF:
1135 if (count == MAX_COUNT)
1136 {
1137 prec = precdbl(prec); count = 0;
1138 if (DEBUGLEVEL) pari_warn(warnprec,"rnflllgram",prec);
1139 nf = nfnewprec_shallow(nf,prec);
1140 }
1141 mth = rel_T2(nf, pol, lx, prec);
1142 if (!mth) { count = MAX_COUNT; goto PRECNF; }
1143 h = NULL;
1144 PRECPB:
1145 if (h)
1146 { /* precision problem, recompute. If no progress, increase nf precision */
1147 if (++count == MAX_COUNT || RgM_isidentity(h)) {count = MAX_COUNT; goto PRECNF;}
1148 H = H? gmul(H, h): h;
1149 MPOL = gmul(MPOL, h);
1150 }
1151 h = matid(lx-1);
1152 MC = mattocomplex(nf, MPOL);
1153 mu = cgetg(lx,t_MAT);
1154 B = cgetg(lx,t_COL);
1155 for (j=1; j<lx; j++)
1156 {
1157 gel(mu,j) = zerocol(lx - 1);
1158 gel(B,j) = gen_0;
1159 }
1160 if (DEBUGLEVEL) err_printf("k = ");
1161 gel(B,1) = real_i(rnfscal(mth,gel(MC,1),gel(MC,1)));
1162 gel(MCS,1) = gel(MC,1);
1163 kmax = 1; k = 2;
1164 do
1165 {
1166 GEN Ik_inv = NULL;
1167 if (DEBUGLEVEL) err_printf("%ld ",k);
1168 if (k > kmax)
1169 { /* Incremental Gram-Schmidt */
1170 kmax = k; gel(MCS,k) = gel(MC,k);
1171 for (j=1; j<k; j++)
1172 {
1173 gcoeff(mu,k,j) = vecdiv(rnfscal(mth,gel(MCS,j),gel(MC,k)),
1174 gel(B,j));
1175 gel(MCS,k) = gsub(gel(MCS,k), vecmul(gcoeff(mu,k,j),gel(MCS,j)));
1176 }
1177 gel(B,k) = real_i(rnfscal(mth,gel(MCS,k),gel(MCS,k)));
1178 if (check_0(gel(B,k))) goto PRECPB;
1179 }
1180 if (!RED(k, k-1, h, mu, MC, nf, I, &Ik_inv)) goto PRECPB;
1181 if (do_SWAP(I,MC,MCS,h,mu,B,kmax,k,alpha, r1))
1182 {
1183 if (!B[k]) goto PRECPB;
1184 if (k > 2) k--;
1185 }
1186 else
1187 {
1188 for (l=k-2; l; l--)
1189 if (!RED(k, l, h, mu, MC, nf, I, &Ik_inv)) goto PRECPB;
1190 k++;
1191 }
1192 if (gc_needed(av,2))
1193 {
1194 if(DEBUGMEM>1) pari_warn(warnmem,"rnflllgram");
1195 gerepileall(av, H?10:9, &nf,&mth,&h,&MPOL,&B,&MC,&MCS,&mu,&I,&H);
1196 }
1197 }
1198 while (k < lx);
1199 MPOL = gmul(MPOL,h);
1200 if (H) h = gmul(H, h);
1201 if (DEBUGLEVEL) err_printf("\n");
1202 MPOL = RgM_to_nfM(nf,MPOL);
1203 h = RgM_to_nfM(nf,h);
1204 return gerepilecopy(av, mkvec2(mkvec2(MPOL,I), h));
1205 }
1206
1207 GEN
rnfpolred(GEN nf,GEN pol,long prec)1208 rnfpolred(GEN nf, GEN pol, long prec)
1209 {
1210 pari_sp av = avma;
1211 long i, j, n, v = varn(pol);
1212 GEN id, w, I, O, bnf, nfpol;
1213
1214 if (typ(pol)!=t_POL) pari_err_TYPE("rnfpolred",pol);
1215 bnf = nf; nf = checknf(bnf);
1216 bnf = (nf == bnf)? NULL: checkbnf(bnf);
1217 if (degpol(pol) <= 1) { w = cgetg(2, t_VEC); gel(w,1) = pol_x(v); return w; }
1218 nfpol = nf_get_pol(nf);
1219
1220 id = rnfpseudobasis(nf,pol);
1221 if (bnf && is_pm1( bnf_get_no(bnf) )) /* if bnf is principal */
1222 {
1223 GEN newI, newO;
1224 O = gel(id,1);
1225 I = gel(id,2); n = lg(I)-1;
1226 newI = cgetg(n+1,t_VEC);
1227 newO = cgetg(n+1,t_MAT);
1228 for (j=1; j<=n; j++)
1229 {
1230 GEN al = gen_if_principal(bnf,gel(I,j));
1231 gel(newI,j) = gen_1;
1232 gel(newO,j) = nfC_nf_mul(nf, gel(O,j), al);
1233 }
1234 id = mkvec2(newO, newI);
1235 }
1236
1237 id = gel(rnflllgram(nf,pol,id,prec),1);
1238 O = gel(id,1);
1239 I = gel(id,2); n = lg(I)-1;
1240 w = cgetg(n+1,t_VEC);
1241 pol = lift_shallow(pol);
1242 for (j=1; j<=n; j++)
1243 {
1244 GEN newpol, L, a, Ij = gel(I,j);
1245 a = RgC_Rg_mul(gel(O,j), (typ(Ij) == t_MAT)? gcoeff(Ij,1,1): Ij);
1246 for (i=n; i; i--) gel(a,i) = nf_to_scalar_or_alg(nf, gel(a,i));
1247 a = RgV_to_RgX(a, v);
1248 newpol = RgXQX_red(RgXQ_charpoly(a, pol, v), nfpol);
1249 newpol = Q_primpart(newpol);
1250
1251 (void)nfgcd_all(newpol, RgX_deriv(newpol), nfpol, nf_get_index(nf), &newpol);
1252 L = leading_coeff(newpol);
1253 gel(w,j) = (typ(L) == t_POL)? RgXQX_div(newpol, L, nfpol)
1254 : RgX_Rg_div(newpol, L);
1255 }
1256 return gerepilecopy(av,w);
1257 }
1258
1259 /*******************************************************************/
1260 /* */
1261 /* LINEAR ALGEBRA OVER Z_K (HNF,SNF) */
1262 /* */
1263 /*******************************************************************/
1264 /* A torsion-free module M over Z_K is given by [A,I].
1265 * I=[a_1,...,a_k] is a row vector of k fractional ideals given in HNF.
1266 * A is an n x k matrix (same k) such that if A_j is the j-th column of A then
1267 * M=a_1 A_1+...+a_k A_k. We say that [A,I] is a pseudo-basis if k=n */
1268
1269 /* Given an element x and an ideal I in HNF, gives an r such that x-r is in H
1270 * and r is small */
1271 GEN
nfreduce(GEN nf,GEN x,GEN I)1272 nfreduce(GEN nf, GEN x, GEN I)
1273 {
1274 pari_sp av = avma;
1275 GEN aI;
1276 x = nf_to_scalar_or_basis(checknf(nf), x);
1277 if (idealtyp(&I,&aI) != id_MAT || lg(I)==1) pari_err_TYPE("nfreduce",I);
1278 if (typ(x) != t_COL) x = scalarcol( gmod(x, gcoeff(I,1,1)), lg(I)-1 );
1279 else x = reducemodinvertible(x, I);
1280 return gerepileupto(av, x);
1281 }
1282 /* Given an element x and an ideal in HNF, gives an a in ideal such that
1283 * x-a is small. No checks */
1284 static GEN
element_close(GEN nf,GEN x,GEN ideal)1285 element_close(GEN nf, GEN x, GEN ideal)
1286 {
1287 pari_sp av = avma;
1288 GEN y = gcoeff(ideal,1,1);
1289 x = nf_to_scalar_or_basis(nf, x);
1290 if (typ(y) == t_INT && is_pm1(y)) return ground(x);
1291 if (typ(x) == t_COL)
1292 x = closemodinvertible(x, ideal);
1293 else
1294 x = gmul(y, gdivround(x,y));
1295 return gerepileupto(av, x);
1296 }
1297
1298 /* A + v B */
1299 static GEN
colcomb1(GEN nf,GEN v,GEN A,GEN B)1300 colcomb1(GEN nf, GEN v, GEN A, GEN B)
1301 {
1302 if (isintzero(v)) return A;
1303 return RgC_to_nfC(nf, RgC_add(A, nfC_nf_mul(nf,B,v)));
1304 }
1305 /* u A + v B */
1306 static GEN
colcomb(GEN nf,GEN u,GEN v,GEN A,GEN B)1307 colcomb(GEN nf, GEN u, GEN v, GEN A, GEN B)
1308 {
1309 if (isintzero(u)) return nfC_nf_mul(nf,B,v);
1310 if (u != gen_1) A = nfC_nf_mul(nf,A,u);
1311 return colcomb1(nf, v, A, B);
1312 }
1313
1314 /* return m[i,1..lim] * x */
1315 static GEN
element_mulvecrow(GEN nf,GEN x,GEN m,long i,long lim)1316 element_mulvecrow(GEN nf, GEN x, GEN m, long i, long lim)
1317 {
1318 long j, l = minss(lg(m), lim+1);
1319 GEN dx, y = cgetg(l, t_VEC);
1320 x = nf_to_scalar_or_basis(nf, x);
1321 if (typ(x) == t_COL)
1322 {
1323 x = zk_multable(nf, Q_remove_denom(x, &dx));
1324 for (j=1; j<l; j++)
1325 {
1326 GEN t = gcoeff(m,i,j);
1327 if (!isintzero(t))
1328 {
1329 if (typ(t) == t_COL)
1330 t = RgM_RgC_mul(x, t);
1331 else
1332 t = ZC_Q_mul(gel(x,1), t);
1333 if (dx) t = gdiv(t, dx);
1334 t = nf_to_scalar_or_basis(nf,t);
1335 }
1336 gel(y,j) = t;
1337 }
1338 }
1339 else
1340 {
1341 for (j=1; j<l; j++) gel(y,j) = gmul(x, gcoeff(m,i,j));
1342 }
1343 return y;
1344 }
1345
1346 /* u Z[s,] + v Z[t,], limitied to the first lim entries */
1347 static GEN
rowcomb(GEN nf,GEN u,GEN v,long s,long t,GEN Z,long lim)1348 rowcomb(GEN nf, GEN u, GEN v, long s, long t, GEN Z, long lim)
1349 {
1350 GEN z;
1351 if (gequal0(u))
1352 z = element_mulvecrow(nf,v,Z,t, lim);
1353 else
1354 {
1355 z = element_mulvecrow(nf,u,Z,s, lim);
1356 if (!gequal0(v)) z = gadd(z, element_mulvecrow(nf,v,Z,t, lim));
1357 }
1358 return z;
1359 }
1360
1361 /* nfbezout(0,b,A,B). Either bB = NULL or b*B */
1362 static GEN
zero_nfbezout(GEN nf,GEN bB,GEN b,GEN A,GEN B,GEN * u,GEN * v,GEN * w,GEN * di)1363 zero_nfbezout(GEN nf,GEN bB, GEN b, GEN A,GEN B,GEN *u,GEN *v,GEN *w,GEN *di)
1364 {
1365 GEN d;
1366 if (isint1(b))
1367 {
1368 *v = gen_1;
1369 *w = A;
1370 d = B;
1371 *di = idealinv(nf,d);
1372 }
1373 else
1374 {
1375 *v = nfinv(nf,b);
1376 *w = idealmul(nf,A,*v);
1377 d = bB? bB: idealmul(nf,b,B);
1378 *di = idealHNF_inv(nf,d);
1379 }
1380 *u = gen_0; return d;
1381 }
1382
1383 /* Given elements a,b and ideals A, B, outputs d = a.A+b.B and gives
1384 * di=d^-1, w=A.B.di, u, v such that au+bv=1 and u in A.di, v in B.di.
1385 * Assume A, B nonzero, but a or b can be zero (not both) */
1386 static GEN
nfbezout(GEN nf,GEN a,GEN b,GEN A,GEN B,GEN * pu,GEN * pv,GEN * pw,GEN * pdi,int red)1387 nfbezout(GEN nf,GEN a,GEN b, GEN A,GEN B, GEN *pu,GEN *pv,GEN *pw,GEN *pdi,
1388 int red)
1389 {
1390 GEN w, u, v, d, di, aA, bB;
1391
1392 if (isintzero(a)) return zero_nfbezout(nf,NULL,b,A,B,pu,pv,pw,pdi);
1393 if (isintzero(b)) return zero_nfbezout(nf,NULL,a,B,A,pv,pu,pw,pdi);
1394
1395 if (a != gen_1) /* frequently called with a = gen_1 */
1396 {
1397 a = nf_to_scalar_or_basis(nf,a);
1398 if (isint1(a)) a = gen_1;
1399 }
1400 aA = (a == gen_1)? idealhnf_shallow(nf,A): idealmul(nf,a,A);
1401 bB = idealmul(nf,b,B);
1402 d = idealadd(nf,aA,bB);
1403 if (gequal(aA, d)) return zero_nfbezout(nf,d, a,B,A,pv,pu,pw,pdi);
1404 if (gequal(bB, d)) return zero_nfbezout(nf,d, b,A,B,pu,pv,pw,pdi);
1405 /* general case is slow */
1406 di = idealHNF_inv(nf,d);
1407 aA = idealmul(nf,aA,di); /* integral */
1408 bB = idealmul(nf,bB,di); /* integral */
1409
1410 u = red? idealaddtoone_i(nf, aA, bB): idealaddtoone_raw(nf, aA, bB);
1411 w = idealmul(nf,aA,B);
1412 v = nfdiv(nf, nfsub(nf, gen_1, u), b);
1413 if (a != gen_1)
1414 {
1415 GEN inva = nfinv(nf, a);
1416 u = nfmul(nf,u,inva);
1417 w = idealmul(nf, inva, w); /* AB/d */
1418 }
1419 *pu = u; *pv = v; *pw = w; *pdi = di; return d;
1420 }
1421 /* v a vector of ideals, simplify in place the ones generated by elts of Q */
1422 static void
idV_simplify(GEN v)1423 idV_simplify(GEN v)
1424 {
1425 long i, l = lg(v);
1426 for (i = 1; i < l; i++)
1427 {
1428 GEN M = gel(v,i);
1429 if (typ(M)==t_MAT && RgM_isscalar(M,NULL))
1430 gel(v,i) = Q_abs_shallow(gcoeff(M,1,1));
1431 }
1432 }
1433 /* Given a torsion-free module x outputs a pseudo-basis for x in HNF */
1434 GEN
nfhnf0(GEN nf,GEN x,long flag)1435 nfhnf0(GEN nf, GEN x, long flag)
1436 {
1437 long i, j, def, idef, m, n;
1438 pari_sp av0 = avma, av;
1439 GEN y, A, I, J, U;
1440
1441 nf = checknf(nf);
1442 check_ZKmodule(x, "nfhnf");
1443 A = gel(x,1); RgM_dimensions(A, &m, &n);
1444 I = gel(x,2);
1445 if (!n) {
1446 if (!flag) return gcopy(x);
1447 retmkvec2(gcopy(x), cgetg(1,t_MAT));
1448 }
1449 U = flag? matid(n): NULL;
1450 idef = (n < m)? m-n : 0;
1451 av = avma;
1452 A = RgM_to_nfM(nf,A);
1453 I = leafcopy(I);
1454 J = zerovec(n); def = n;
1455 for (i=m; i>idef; i--)
1456 {
1457 GEN d, di = NULL;
1458
1459 j=def; while (j>=1 && isintzero(gcoeff(A,i,j))) j--;
1460 if (!j)
1461 { /* no pivot on line i */
1462 if (idef) idef--;
1463 continue;
1464 }
1465 if (j==def) j--;
1466 else {
1467 swap(gel(A,j), gel(A,def));
1468 swap(gel(I,j), gel(I,def));
1469 if (U) swap(gel(U,j), gel(U,def));
1470 }
1471 for ( ; j; j--)
1472 {
1473 GEN a,b, u,v,w, S, T, S0, T0 = gel(A,j);
1474 b = gel(T0,i); if (isintzero(b)) continue;
1475
1476 S0 = gel(A,def); a = gel(S0,i);
1477 d = nfbezout(nf, a,b, gel(I,def),gel(I,j), &u,&v,&w,&di,1);
1478 S = colcomb(nf, u,v, S0,T0);
1479 T = colcomb(nf, a,gneg(b), T0,S0);
1480 gel(A,def) = S; gel(A,j) = T;
1481 gel(I,def) = d; gel(I,j) = w;
1482 if (U)
1483 {
1484 S0 = gel(U,def);
1485 T0 = gel(U,j);
1486 gel(U,def) = colcomb(nf, u,v, S0,T0);
1487 gel(U,j) = colcomb(nf, a,gneg(b), T0,S0);
1488 }
1489 }
1490 y = gcoeff(A,i,def);
1491 if (!isint1(y))
1492 {
1493 GEN yi = nfinv(nf,y);
1494 gel(A,def) = nfC_nf_mul(nf, gel(A,def), yi);
1495 gel(I,def) = idealmul(nf, y, gel(I,def));
1496 if (U) gel(U,def) = nfC_nf_mul(nf, gel(U,def), yi);
1497 di = NULL;
1498 }
1499 if (!di) di = idealinv(nf,gel(I,def));
1500 d = gel(I,def);
1501 gel(J,def) = di;
1502 for (j=def+1; j<=n; j++)
1503 {
1504 GEN mc, c = gcoeff(A,i,j); if (isintzero(c)) continue;
1505 c = element_close(nf, c, idealmul(nf,d,gel(J,j)));
1506 mc = gneg(c);
1507 gel(A,j) = colcomb1(nf, mc, gel(A,j),gel(A,def));
1508 if (U) gel(U,j) = colcomb1(nf, mc, gel(U,j),gel(U,def));
1509 }
1510 def--;
1511 if (gc_needed(av,2))
1512 {
1513 if(DEBUGMEM>1) pari_warn(warnmem,"nfhnf, i = %ld", i);
1514 gerepileall(av,U?4:3, &A,&I,&J,&U);
1515 }
1516 }
1517 n -= def;
1518 A += def; A[0] = evaltyp(t_MAT)|evallg(n+1);
1519 I += def; I[0] = evaltyp(t_VEC)|evallg(n+1);
1520 idV_simplify(I);
1521 x = mkvec2(A,I);
1522 if (U) x = mkvec2(x,U);
1523 return gerepilecopy(av0, x);
1524 }
1525
1526 GEN
nfhnf(GEN nf,GEN x)1527 nfhnf(GEN nf, GEN x) { return nfhnf0(nf, x, 0); }
1528
1529 static GEN
RgV_find_denom(GEN x)1530 RgV_find_denom(GEN x)
1531 {
1532 long i, l = lg(x);
1533 for (i = 1; i < l; i++)
1534 if (Q_denom(gel(x,i)) != gen_1) return gel(x,i);
1535 return NULL;
1536 }
1537 /* A torsion module M over Z_K will be given by a row vector [A,I,J] with
1538 * three components. I=[b_1,...,b_n] is a row vector of n fractional ideals
1539 * given in HNF, J=[a_1,...,a_n] is a row vector of n fractional ideals in
1540 * HNF. A is an nxn matrix (same n) such that if A_j is the j-th column of A
1541 * and e_n is the canonical basis of K^n, then
1542 * M=(b_1e_1+...+b_ne_n)/(a_1A_1+...a_nA_n) */
1543
1544 /* x=[A,I,J] a torsion module as above. Output the
1545 * smith normal form as K=[c_1,...,c_n] such that x = Z_K/c_1+...+Z_K/c_n */
1546 GEN
nfsnf0(GEN nf,GEN x,long flag)1547 nfsnf0(GEN nf, GEN x, long flag)
1548 {
1549 long i, j, k, l, n, m;
1550 pari_sp av;
1551 GEN z,u,v,w,d,dinv,A,I,J, U,V;
1552
1553 nf = checknf(nf);
1554 if (typ(x)!=t_VEC || lg(x)!=4) pari_err_TYPE("nfsnf",x);
1555 A = gel(x,1);
1556 I = gel(x,2);
1557 J = gel(x,3);
1558 if (typ(A)!=t_MAT) pari_err_TYPE("nfsnf",A);
1559 n = lg(A)-1;
1560 if (typ(I)!=t_VEC) pari_err_TYPE("nfsnf",I);
1561 if (typ(J)!=t_VEC) pari_err_TYPE("nfsnf",J);
1562 if (lg(I)!=n+1 || lg(J)!=n+1) pari_err_DIM("nfsnf");
1563 RgM_dimensions(A, &m, &n);
1564 if (!n || n != m) pari_err_IMPL("nfsnf for empty or non square matrices");
1565
1566 av = avma;
1567 if (!flag) U = V = NULL;
1568 else
1569 {
1570 U = matid(m);
1571 V = matid(n);
1572 }
1573 A = RgM_to_nfM(nf, A);
1574 I = leafcopy(I);
1575 J = leafcopy(J);
1576 for (i = 1; i <= n; i++) gel(J,i) = idealinv(nf, gel(J,i));
1577 z = zerovec(n);
1578 for (i=n; i>=1; i--)
1579 {
1580 GEN Aii, a, b, db;
1581 long c = 0;
1582 for (j=i-1; j>=1; j--)
1583 {
1584 GEN S, T, S0, T0 = gel(A,j);
1585 b = gel(T0,i); if (gequal0(b)) continue;
1586
1587 S0 = gel(A,i); a = gel(S0,i);
1588 d = nfbezout(nf, a,b, gel(J,i),gel(J,j), &u,&v,&w,&dinv,1);
1589 S = colcomb(nf, u,v, S0,T0);
1590 T = colcomb(nf, a,gneg(b), T0,S0);
1591 gel(A,i) = S; gel(A,j) = T;
1592 gel(J,i) = d; gel(J,j) = w;
1593 if (V)
1594 {
1595 T0 = gel(V,j);
1596 S0 = gel(V,i);
1597 gel(V,i) = colcomb(nf, u,v, S0,T0);
1598 gel(V,j) = colcomb(nf, a,gneg(b), T0,S0);
1599 }
1600 }
1601 for (j=i-1; j>=1; j--)
1602 {
1603 GEN ri, rj;
1604 b = gcoeff(A,j,i); if (gequal0(b)) continue;
1605
1606 a = gcoeff(A,i,i);
1607 d = nfbezout(nf, a,b, gel(I,i),gel(I,j), &u,&v,&w,&dinv,1);
1608 ri = rowcomb(nf, u,v, i,j, A, i);
1609 rj = rowcomb(nf, a,gneg(b), j,i, A, i);
1610 for (k=1; k<=i; k++) {
1611 gcoeff(A,j,k) = gel(rj,k);
1612 gcoeff(A,i,k) = gel(ri,k);
1613 }
1614 if (U)
1615 {
1616 ri = rowcomb(nf, u,v, i,j, U, m);
1617 rj = rowcomb(nf, a,gneg(b), j,i, U, m);
1618 for (k=1; k<=m; k++) {
1619 gcoeff(U,j,k) = gel(rj,k);
1620 gcoeff(U,i,k) = gel(ri,k);
1621 }
1622 }
1623 gel(I,i) = d; gel(I,j) = w; c = 1;
1624 }
1625 if (c) { i++; continue; }
1626
1627 Aii = gcoeff(A,i,i); if (gequal0(Aii)) continue;
1628 gel(J,i) = idealmul(nf, gel(J,i), Aii);
1629 gcoeff(A,i,i) = gen_1;
1630 if (V) gel(V,i) = nfC_nf_mul(nf, gel(V,i), nfinv(nf,Aii));
1631 gel(z,i) = idealmul(nf,gel(J,i),gel(I,i));
1632 b = Q_remove_denom(gel(z,i), &db);
1633 for (k=1; k<i; k++)
1634 for (l=1; l<i; l++)
1635 {
1636 GEN d, D, p1, p2, p3, Akl = gcoeff(A,k,l);
1637 long t;
1638 if (gequal0(Akl)) continue;
1639
1640 p1 = idealmul(nf,Akl,gel(J,l));
1641 p3 = idealmul(nf, p1, gel(I,k));
1642 if (db) p3 = RgM_Rg_mul(p3, db);
1643 if (RgM_is_ZM(p3) && hnfdivide(b, p3)) continue;
1644
1645 /* find d in D = I[k]/I[i] not in J[i]/(A[k,l] J[l]) */
1646 D = idealdiv(nf,gel(I,k),gel(I,i));
1647 p2 = idealdiv(nf,gel(J,i), p1);
1648 d = RgV_find_denom(QM_gauss(p2, D));
1649 if (!d) pari_err_BUG("nfsnf");
1650 p1 = element_mulvecrow(nf,d,A,k,i);
1651 for (t=1; t<=i; t++) gcoeff(A,i,t) = gadd(gcoeff(A,i,t),gel(p1,t));
1652 if (U)
1653 {
1654 p1 = element_mulvecrow(nf,d,U,k,i);
1655 for (t=1; t<=i; t++) gcoeff(U,i,t) = gadd(gcoeff(U,i,t),gel(p1,t));
1656 }
1657
1658 k = i; c = 1; break;
1659 }
1660 if (gc_needed(av,1))
1661 {
1662 if(DEBUGMEM>1) pari_warn(warnmem,"nfsnf");
1663 gerepileall(av,U?6:4, &A,&I,&J,&z,&U,&V);
1664 }
1665 if (c) i++; /* iterate on row/column i */
1666 }
1667 if (U) z = mkvec3(z,U,V);
1668 return gerepilecopy(av, z);
1669 }
1670 GEN
nfsnf(GEN nf,GEN x)1671 nfsnf(GEN nf, GEN x) { return nfsnf0(nf,x,0); }
1672
1673 /* Given a pseudo-basis x, outputs a multiple of its ideal determinant */
1674 GEN
nfdetint(GEN nf,GEN x)1675 nfdetint(GEN nf, GEN x)
1676 {
1677 GEN pass,c,v,det1,piv,pivprec,vi,p1,A,I,id,idprod;
1678 long i, j, k, rg, n, m, m1, cm=0, N;
1679 pari_sp av = avma, av1;
1680
1681 nf = checknf(nf); N = nf_get_degree(nf);
1682 check_ZKmodule(x, "nfdetint");
1683 A = gel(x,1);
1684 I = gel(x,2);
1685 n = lg(A)-1; if (!n) return gen_1;
1686
1687 m1 = lgcols(A); m = m1-1;
1688 id = matid(N);
1689 c = new_chunk(m1); for (k=1; k<=m; k++) c[k] = 0;
1690 piv = pivprec = gen_1;
1691
1692 av1 = avma;
1693 det1 = idprod = gen_0; /* dummy for gerepileall */
1694 pass = cgetg(m1,t_MAT);
1695 v = cgetg(m1,t_COL);
1696 for (j=1; j<=m; j++)
1697 {
1698 gel(pass,j) = zerocol(m);
1699 gel(v,j) = gen_0; /* dummy */
1700 }
1701 for (rg=0,k=1; k<=n; k++)
1702 {
1703 long t = 0;
1704 for (i=1; i<=m; i++)
1705 if (!c[i])
1706 {
1707 vi=nfmul(nf,piv,gcoeff(A,i,k));
1708 for (j=1; j<=m; j++)
1709 if (c[j]) vi=gadd(vi,nfmul(nf,gcoeff(pass,i,j),gcoeff(A,j,k)));
1710 gel(v,i) = vi; if (!t && !gequal0(vi)) t=i;
1711 }
1712 if (t)
1713 {
1714 pivprec = piv;
1715 if (rg == m-1)
1716 {
1717 if (!cm)
1718 {
1719 cm=1; idprod = id;
1720 for (i=1; i<=m; i++)
1721 if (i!=t)
1722 idprod = (idprod==id)? gel(I,c[i])
1723 : idealmul(nf,idprod,gel(I,c[i]));
1724 }
1725 p1 = idealmul(nf,gel(v,t),gel(I,k)); c[t]=0;
1726 det1 = (typ(det1)==t_INT)? p1: idealadd(nf,p1,det1);
1727 }
1728 else
1729 {
1730 rg++; piv=gel(v,t); c[t]=k;
1731 for (i=1; i<=m; i++)
1732 if (!c[i])
1733 {
1734 for (j=1; j<=m; j++)
1735 if (c[j] && j!=t)
1736 {
1737 p1 = gsub(nfmul(nf,piv,gcoeff(pass,i,j)),
1738 nfmul(nf,gel(v,i),gcoeff(pass,t,j)));
1739 gcoeff(pass,i,j) = rg>1? nfdiv(nf,p1,pivprec)
1740 : p1;
1741 }
1742 gcoeff(pass,i,t) = gneg(gel(v,i));
1743 }
1744 }
1745 }
1746 if (gc_needed(av1,1))
1747 {
1748 if(DEBUGMEM>1) pari_warn(warnmem,"nfdetint");
1749 gerepileall(av1,6, &det1,&piv,&pivprec,&pass,&v,&idprod);
1750 }
1751 }
1752 if (!cm) { set_avma(av); return cgetg(1,t_MAT); }
1753 return gerepileupto(av, idealmul(nf,idprod,det1));
1754 }
1755
1756 /* reduce in place components of x[1..lim] mod D (destroy x). D in HNF */
1757 static void
nfcleanmod(GEN nf,GEN x,long lim,GEN D)1758 nfcleanmod(GEN nf, GEN x, long lim, GEN D)
1759 {
1760 GEN DZ, DZ2, dD;
1761 long i;
1762 D = Q_remove_denom(D, &dD);
1763 DZ = gcoeff(D,1,1); DZ2 = shifti(DZ, -1);
1764 for (i = 1; i <= lim; i++)
1765 {
1766 GEN c = nf_to_scalar_or_basis(nf, gel(x,i));
1767 switch(typ(c)) /* c = centermod(c, D) */
1768 {
1769 case t_INT:
1770 if (!signe(c)) break;
1771 if (dD) c = mulii(c, dD);
1772 c = centermodii(c, DZ, DZ2);
1773 if (dD) c = Qdivii(c,dD);
1774 break;
1775 case t_FRAC: {
1776 GEN dc = gel(c,2), nc = gel(c,1), N = mulii(DZ, dc);
1777 if (dD) nc = mulii(nc, dD);
1778 c = centermodii(nc, N, shifti(N,-1));
1779 c = Qdivii(c, dD ? mulii(dc,dD): dc);
1780 break;
1781 }
1782 case t_COL: {
1783 GEN dc;
1784 c = Q_remove_denom(c, &dc);
1785 if (dD) c = ZC_Z_mul(c, dD);
1786 c = ZC_hnfrem(c, dc? ZM_Z_mul(D,dc): D);
1787 dc = mul_content(dc, dD);
1788 if (ZV_isscalar(c))
1789 {
1790 c = gel(c,1);
1791 if (dc) c = Qdivii(c,dc);
1792 }
1793 else
1794 if (dc) c = RgC_Rg_div(c, dc);
1795 break;
1796 }
1797 }
1798 gel(x,i) = c;
1799 }
1800 }
1801
1802 GEN
nfhnfmod(GEN nf,GEN x,GEN D)1803 nfhnfmod(GEN nf, GEN x, GEN D)
1804 {
1805 long li, co, i, j, def, ldef;
1806 pari_sp av0=avma, av;
1807 GEN dA, dI, d0, w, p1, d, u, v, A, I, J, di;
1808
1809 nf = checknf(nf);
1810 check_ZKmodule(x, "nfhnfmod");
1811 A = gel(x,1);
1812 I = gel(x,2);
1813 co = lg(A); if (co==1) return cgetg(1,t_MAT);
1814
1815 li = lgcols(A);
1816 if (typ(D)!=t_MAT) D = idealhnf_shallow(nf, D);
1817 D = Q_remove_denom(D, NULL);
1818 RgM_check_ZM(D, "nfhnfmod");
1819
1820 av = avma;
1821 A = RgM_to_nfM(nf, A);
1822 A = Q_remove_denom(A, &dA);
1823 I = Q_remove_denom(leafcopy(I), &dI);
1824 dA = mul_denom(dA,dI);
1825 if (dA) D = ZM_Z_mul(D, powiu(dA, minss(li,co)));
1826
1827 def = co; ldef = (li>co)? li-co+1: 1;
1828 for (i=li-1; i>=ldef; i--)
1829 {
1830 def--; j=def; while (j>=1 && isintzero(gcoeff(A,i,j))) j--;
1831 if (!j) continue;
1832 if (j==def) j--;
1833 else {
1834 swap(gel(A,j), gel(A,def));
1835 swap(gel(I,j), gel(I,def));
1836 }
1837 for ( ; j; j--)
1838 {
1839 GEN a, b, S, T, S0, T0 = gel(A,j);
1840 b = gel(T0,i); if (isintzero(b)) continue;
1841
1842 S0 = gel(A,def); a = gel(S0,i);
1843 d = nfbezout(nf, a,b, gel(I,def),gel(I,j), &u,&v,&w,&di,0);
1844 S = colcomb(nf, u,v, S0,T0);
1845 T = colcomb(nf, a,gneg(b), T0,S0);
1846 if (u != gen_0 && v != gen_0) /* already reduced otherwise */
1847 nfcleanmod(nf, S, i, idealmul(nf,D,di));
1848 nfcleanmod(nf, T, i, idealdiv(nf,D,w));
1849 gel(A,def) = S; gel(A,j) = T;
1850 gel(I,def) = d; gel(I,j) = w;
1851 }
1852 if (gc_needed(av,2))
1853 {
1854 if(DEBUGMEM>1) pari_warn(warnmem,"[1]: nfhnfmod, i = %ld", i);
1855 gerepileall(av,dA? 4: 3, &A,&I,&D,&dA);
1856 }
1857 }
1858 def--; d0 = D;
1859 A += def; A[0] = evaltyp(t_MAT)|evallg(li);
1860 I += def; I[0] = evaltyp(t_VEC)|evallg(li);
1861 J = cgetg(li,t_VEC);
1862 for (i=li-1; i>=1; i--)
1863 {
1864 GEN b = gcoeff(A,i,i);
1865 d = nfbezout(nf, gen_1,b, d0,gel(I,i), &u,&v,&w,&di,0);
1866 p1 = nfC_nf_mul(nf,gel(A,i),v);
1867 if (i > 1)
1868 {
1869 d0 = idealmul(nf,d0,di);
1870 nfcleanmod(nf, p1, i, d0);
1871 }
1872 gel(A,i) = p1; gel(p1,i) = gen_1;
1873 gel(I,i) = d;
1874 gel(J,i) = di;
1875 }
1876 for (i=li-2; i>=1; i--)
1877 {
1878 d = gel(I,i);
1879 for (j=i+1; j<li; j++)
1880 {
1881 GEN c = gcoeff(A,i,j); if (isintzero(c)) continue;
1882 c = element_close(nf, c, idealmul(nf,d,gel(J,j)));
1883 gel(A,j) = colcomb1(nf, gneg(c), gel(A,j),gel(A,i));
1884 }
1885 if (gc_needed(av,2))
1886 {
1887 if(DEBUGMEM>1) pari_warn(warnmem,"[2]: nfhnfmod, i = %ld", i);
1888 gerepileall(av,dA? 4: 3, &A,&I,&J,&dA);
1889 }
1890 }
1891 idV_simplify(I);
1892 if (dA) I = gdiv(I,dA);
1893 return gerepilecopy(av0, mkvec2(A, I));
1894 }
1895