1 /* Copyright (C) 2007 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 #include "pari.h"
16 #include "paripriv.h"
17
18 /*******************************************************************/
19 /* */
20 /* ZX */
21 /* */
22 /*******************************************************************/
23 void
RgX_check_QX(GEN x,const char * s)24 RgX_check_QX(GEN x, const char *s)
25 { if (!RgX_is_QX(x)) pari_err_TYPE(stack_strcat(s," [not in Q[X]]"), x); }
26 void
RgX_check_ZX(GEN x,const char * s)27 RgX_check_ZX(GEN x, const char *s)
28 { if (!RgX_is_ZX(x)) pari_err_TYPE(stack_strcat(s," [not in Z[X]]"), x); }
29 long
ZX_max_lg(GEN x)30 ZX_max_lg(GEN x)
31 {
32 long i, prec = 0, lx = lg(x);
33
34 for (i=2; i<lx; i++) { long l = lgefint(gel(x,i)); if (l > prec) prec = l; }
35 return prec;
36 }
37
38 GEN
ZX_add(GEN x,GEN y)39 ZX_add(GEN x, GEN y)
40 {
41 long lx,ly,i;
42 GEN z;
43 lx = lg(x); ly = lg(y); if (lx < ly) swapspec(x,y, lx,ly);
44 z = cgetg(lx,t_POL); z[1] = x[1];
45 for (i=2; i<ly; i++) gel(z,i) = addii(gel(x,i),gel(y,i));
46 for ( ; i<lx; i++) gel(z,i) = icopy(gel(x,i));
47 if (lx == ly) z = ZX_renormalize(z, lx);
48 if (!lgpol(z)) { set_avma((pari_sp)(z + lx)); return pol_0(varn(x)); }
49 return z;
50 }
51
52 GEN
ZX_sub(GEN x,GEN y)53 ZX_sub(GEN x,GEN y)
54 {
55 long i, lx = lg(x), ly = lg(y);
56 GEN z;
57 if (lx >= ly)
58 {
59 z = cgetg(lx,t_POL); z[1] = x[1];
60 for (i=2; i<ly; i++) gel(z,i) = subii(gel(x,i),gel(y,i));
61 if (lx == ly)
62 {
63 z = ZX_renormalize(z, lx);
64 if (!lgpol(z)) { set_avma((pari_sp)(z + lx)); z = pol_0(varn(x)); }
65 }
66 else
67 for ( ; i<lx; i++) gel(z,i) = icopy(gel(x,i));
68 }
69 else
70 {
71 z = cgetg(ly,t_POL); z[1] = y[1];
72 for (i=2; i<lx; i++) gel(z,i) = subii(gel(x,i),gel(y,i));
73 for ( ; i<ly; i++) gel(z,i) = negi(gel(y,i));
74 }
75 return z;
76 }
77
78 GEN
ZX_neg(GEN x)79 ZX_neg(GEN x)
80 {
81 long i, l = lg(x);
82 GEN y = cgetg(l,t_POL);
83 y[1] = x[1]; for(i=2; i<l; i++) gel(y,i) = negi(gel(x,i));
84 return y;
85 }
86 GEN
ZX_copy(GEN x)87 ZX_copy(GEN x)
88 {
89 long i, l = lg(x);
90 GEN y = cgetg(l, t_POL);
91 y[1] = x[1];
92 for (i=2; i<l; i++)
93 {
94 GEN c = gel(x,i);
95 gel(y,i) = lgefint(c) == 2? gen_0: icopy(c);
96 }
97 return y;
98 }
99
100 GEN
scalar_ZX(GEN x,long v)101 scalar_ZX(GEN x, long v)
102 {
103 GEN z;
104 if (!signe(x)) return pol_0(v);
105 z = cgetg(3, t_POL);
106 z[1] = evalsigne(1) | evalvarn(v);
107 gel(z,2) = icopy(x); return z;
108 }
109
110 GEN
scalar_ZX_shallow(GEN x,long v)111 scalar_ZX_shallow(GEN x, long v)
112 {
113 GEN z;
114 if (!signe(x)) return pol_0(v);
115 z = cgetg(3, t_POL);
116 z[1] = evalsigne(1) | evalvarn(v);
117 gel(z,2) = x; return z;
118 }
119
120 GEN
ZX_Z_add(GEN y,GEN x)121 ZX_Z_add(GEN y, GEN x)
122 {
123 long lz, i;
124 GEN z = cgetg_copy(y, &lz);
125 if (lz == 2) { set_avma((pari_sp)(z + 2)); return scalar_ZX(x,varn(y)); }
126 z[1] = y[1];
127 gel(z,2) = addii(gel(y,2),x);
128 for(i=3; i<lz; i++) gel(z,i) = icopy(gel(y,i));
129 if (lz==3) z = ZX_renormalize(z,lz);
130 return z;
131 }
132 GEN
ZX_Z_add_shallow(GEN y,GEN x)133 ZX_Z_add_shallow(GEN y, GEN x)
134 {
135 long lz, i;
136 GEN z = cgetg_copy(y, &lz);
137 if (lz == 2) { set_avma((pari_sp)(z + 2)); return scalar_ZX_shallow(x,varn(y)); }
138 z[1] = y[1];
139 gel(z,2) = addii(gel(y,2),x);
140 for(i=3; i<lz; i++) gel(z,i) = gel(y,i);
141 if (lz==3) z = ZX_renormalize(z,lz);
142 return z;
143 }
144
145 GEN
ZX_Z_sub(GEN y,GEN x)146 ZX_Z_sub(GEN y, GEN x)
147 {
148 long lz, i;
149 GEN z = cgetg_copy(y, &lz);
150 if (lz == 2)
151 { /* scalarpol(negi(x), v) */
152 long v = varn(y);
153 set_avma((pari_sp)(z + 2));
154 if (!signe(x)) return pol_0(v);
155 z = cgetg(3,t_POL);
156 z[1] = evalvarn(v) | evalsigne(1);
157 gel(z,2) = negi(x); return z;
158 }
159 z[1] = y[1];
160 gel(z,2) = subii(gel(y,2),x);
161 for(i=3; i<lz; i++) gel(z,i) = icopy(gel(y,i));
162 if (lz==3) z = ZX_renormalize(z,lz);
163 return z;
164 }
165
166 GEN
Z_ZX_sub(GEN x,GEN y)167 Z_ZX_sub(GEN x, GEN y)
168 {
169 long lz, i;
170 GEN z = cgetg_copy(y, &lz);
171 if (lz == 2) { set_avma((pari_sp)(z + 2)); return scalar_ZX(x,varn(y)); }
172 z[1] = y[1];
173 gel(z,2) = subii(x, gel(y,2));
174 for(i=3; i<lz; i++) gel(z,i) = negi(gel(y,i));
175 if (lz==3) z = ZX_renormalize(z,lz);
176 return z;
177 }
178
179 GEN
ZX_Z_divexact(GEN y,GEN x)180 ZX_Z_divexact(GEN y,GEN x)
181 {
182 long i, l = lg(y);
183 GEN z = cgetg(l,t_POL); z[1] = y[1];
184 for(i=2; i<l; i++) gel(z,i) = diviiexact(gel(y,i),x);
185 return z;
186 }
187
188 GEN
ZX_divuexact(GEN y,ulong x)189 ZX_divuexact(GEN y, ulong x)
190 {
191 long i, l = lg(y);
192 GEN z = cgetg(l,t_POL); z[1] = y[1];
193 for(i=2; i<l; i++) gel(z,i) = diviuexact(gel(y,i),x);
194 return z;
195 }
196
197 GEN
zx_z_divexact(GEN y,long x)198 zx_z_divexact(GEN y, long x)
199 {
200 long i, l = lg(y);
201 GEN z = cgetg(l,t_VECSMALL); z[1] = y[1];
202 for (i=2; i<l; i++) z[i] = y[i]/x;
203 return z;
204 }
205
206 GEN
ZX_Z_mul(GEN y,GEN x)207 ZX_Z_mul(GEN y,GEN x)
208 {
209 GEN z;
210 long i, l;
211 if (!signe(x)) return pol_0(varn(y));
212 l = lg(y); z = cgetg(l,t_POL); z[1] = y[1];
213 for(i=2; i<l; i++) gel(z,i) = mulii(gel(y,i),x);
214 return z;
215 }
216
217 GEN
ZX_mulu(GEN y,ulong x)218 ZX_mulu(GEN y, ulong x)
219 {
220 GEN z;
221 long i, l;
222 if (!x) return pol_0(varn(y));
223 l = lg(y); z = cgetg(l,t_POL); z[1] = y[1];
224 for(i=2; i<l; i++) gel(z,i) = mului(x,gel(y,i));
225 return z;
226 }
227
228 GEN
ZX_shifti(GEN y,long n)229 ZX_shifti(GEN y, long n)
230 {
231 GEN z;
232 long i, l;
233 l = lg(y); z = cgetg(l,t_POL); z[1] = y[1];
234 for(i=2; i<l; i++) gel(z,i) = shifti(gel(y,i),n);
235 return ZX_renormalize(z,l);
236 }
237
238 GEN
ZX_remi2n(GEN y,long n)239 ZX_remi2n(GEN y, long n)
240 {
241 GEN z;
242 long i, l;
243 l = lg(y); z = cgetg(l,t_POL); z[1] = y[1];
244 for(i=2; i<l; i++) gel(z,i) = remi2n(gel(y,i),n);
245 return ZX_renormalize(z,l);
246 }
247
248 GEN
ZXT_remi2n(GEN z,long n)249 ZXT_remi2n(GEN z, long n)
250 {
251 if (typ(z) == t_POL)
252 return ZX_remi2n(z, n);
253 else
254 {
255 long i,l = lg(z);
256 GEN x = cgetg(l, t_VEC);
257 for (i=1; i<l; i++) gel(x,i) = ZXT_remi2n(gel(z,i), n);
258 return x;
259 }
260 }
261
262 GEN
zx_to_ZX(GEN z)263 zx_to_ZX(GEN z)
264 {
265 long i, l = lg(z);
266 GEN x = cgetg(l,t_POL);
267 for (i=2; i<l; i++) gel(x,i) = stoi(z[i]);
268 x[1] = evalsigne(l-2!=0)| z[1]; return x;
269 }
270
271 GEN
ZX_deriv(GEN x)272 ZX_deriv(GEN x)
273 {
274 long i,lx = lg(x)-1;
275 GEN y;
276
277 if (lx<3) return pol_0(varn(x));
278 y = cgetg(lx,t_POL);
279 for (i=2; i<lx ; i++) gel(y,i) = mului(i-1,gel(x,i+1));
280 y[1] = x[1]; return y;
281 }
282
283 int
ZX_equal(GEN V,GEN W)284 ZX_equal(GEN V, GEN W)
285 {
286 long i, l = lg(V);
287 if (lg(W) != l) return 0;
288 for (i = 2; i < l; i++)
289 if (!equalii(gel(V,i), gel(W,i))) return 0;
290 return 1;
291 }
292
293 static long
ZX_valspec(GEN x,long nx)294 ZX_valspec(GEN x, long nx)
295 {
296 long vx;
297 for (vx = 0; vx<nx ; vx++)
298 if (signe(gel(x,vx))) break;
299 return vx;
300 }
301
302 long
ZX_val(GEN x)303 ZX_val(GEN x)
304 {
305 long vx;
306 if (!signe(x)) return LONG_MAX;
307 for (vx = 0;; vx++)
308 if (signe(gel(x,2+vx))) break;
309 return vx;
310 }
311 long
ZX_valrem(GEN x,GEN * Z)312 ZX_valrem(GEN x, GEN *Z)
313 {
314 long vx;
315 if (!signe(x)) { *Z = pol_0(varn(x)); return LONG_MAX; }
316 for (vx = 0;; vx++)
317 if (signe(gel(x,2+vx))) break;
318 *Z = RgX_shift_shallow(x, -vx);
319 return vx;
320 }
321
322 GEN
ZX_div_by_X_1(GEN a,GEN * r)323 ZX_div_by_X_1(GEN a, GEN *r)
324 {
325 long l = lg(a), i;
326 GEN a0, z0, z = cgetg(l-1, t_POL);
327 z[1] = a[1];
328 a0 = a + l-1;
329 z0 = z + l-2; *z0 = *a0--;
330 for (i=l-3; i>1; i--) /* z[i] = a[i+1] + z[i+1] */
331 {
332 GEN t = addii(gel(a0--,0), gel(z0--,0));
333 gel(z0,0) = t;
334 }
335 if (r) *r = addii(gel(a0,0), gel(z0,0));
336 return z;
337 }
338
339 /* return P(X + c) using destructive Horner, optimize for c = 1,-1 */
340 static GEN
ZX_translate_basecase(GEN P,GEN c)341 ZX_translate_basecase(GEN P, GEN c)
342 {
343 pari_sp av = avma;
344 GEN Q, R;
345 long i, k, n;
346
347 if (!signe(P) || !signe(c)) return ZX_copy(P);
348 Q = leafcopy(P);
349 R = Q+2; n = degpol(P);
350 if (equali1(c))
351 {
352 for (i=1; i<=n; i++)
353 {
354 for (k=n-i; k<n; k++) gel(R,k) = addii(gel(R,k), gel(R,k+1));
355 if (gc_needed(av,2))
356 {
357 if(DEBUGMEM>1) pari_warn(warnmem,"ZX_translate(1), i = %ld/%ld", i,n);
358 Q = gerepilecopy(av, Q); R = Q+2;
359 }
360 }
361 }
362 else if (equalim1(c))
363 {
364 for (i=1; i<=n; i++)
365 {
366 for (k=n-i; k<n; k++) gel(R,k) = subii(gel(R,k), gel(R,k+1));
367 if (gc_needed(av,2))
368 {
369 if(DEBUGMEM>1) pari_warn(warnmem,"ZX_translate(-1), i = %ld/%ld", i,n);
370 Q = gerepilecopy(av, Q); R = Q+2;
371 }
372 }
373 }
374 else
375 {
376 for (i=1; i<=n; i++)
377 {
378 for (k=n-i; k<n; k++) gel(R,k) = addmulii_inplace(gel(R,k), c, gel(R,k+1));
379 if (gc_needed(av,2))
380 {
381 if(DEBUGMEM>1) pari_warn(warnmem,"ZX_translate, i = %ld/%ld", i,n);
382 Q = gerepilecopy(av, Q); R = Q+2;
383 }
384 }
385 }
386 return gerepilecopy(av, Q);
387 }
388
389 static GEN
Z_Xpm1_powu(long n,long s,long v)390 Z_Xpm1_powu(long n, long s, long v)
391 {
392 long d, k;
393 GEN C;
394 if (!n) return pol_1(v);
395 d = (n + 1) >> 1;
396 C = cgetg(n+3, t_POL);
397 C[1] = evalsigne(1)| evalvarn(v);
398 gel(C,2) = gen_1;
399 gel(C,3) = utoipos(n);
400 for (k=2; k <= d; k++)
401 gel(C,k+2) = diviuexact(mului(n-k+1, gel(C,k+1)), k);
402 if (s < 0)
403 for (k = odd(n)? 0: 1; k <= d; k += 2)
404 togglesign_safe(&gel(C,k+2));
405 if (s > 0 || !odd(n))
406 for (k = d+1; k <= n; k++) gel(C,k+2) = gel(C,n-k+2);
407 else
408 for (k = d+1; k <= n; k++) gel(C,k+2) = negi(gel(C,n-k+2));
409 return C;
410 }
411 /* return (x+u)^n */
412 static GEN
Z_XpN_powu(GEN u,long n,long v)413 Z_XpN_powu(GEN u, long n, long v)
414 {
415 pari_sp av;
416 long k;
417 GEN B, C, V;
418 if (!n) return pol_1(v);
419 if (is_pm1(u))
420 return Z_Xpm1_powu(n, signe(u), v);
421 av = avma;
422 V = gpowers(u, n);
423 B = vecbinomial(n);
424 C = cgetg(n+3, t_POL);
425 C[1] = evalsigne(1)| evalvarn(v);
426 for (k=1; k <= n+1; k++)
427 gel(C,k+1) = mulii(gel(V,n+2-k), gel(B,k));
428 return gerepileupto(av, C);
429 }
430
431 GEN
ZX_translate(GEN P,GEN c)432 ZX_translate(GEN P, GEN c)
433 {
434 pari_sp av = avma;
435 long n = degpol(P);
436 if (n < 220)
437 return ZX_translate_basecase(P, c);
438 else
439 {
440 long d = n >> 1;
441 GEN Q = ZX_translate(RgX_shift_shallow(P, -d), c);
442 GEN R = ZX_translate(RgXn_red_shallow(P, d), c);
443 GEN S = Z_XpN_powu(c, d, varn(P));
444 return gerepileupto(av, ZX_add(ZX_mul(Q, S), R));
445 }
446 }
447
448 /* Return 2^(n degpol(P)) P(x >> n) */
449 GEN
ZX_rescale2n(GEN P,long n)450 ZX_rescale2n(GEN P, long n)
451 {
452 long i, l = lg(P), ni = n;
453 GEN Q;
454 if (l==2) return pol_0(varn(P));
455 Q = cgetg(l,t_POL);
456 gel(Q,l-1) = icopy(gel(P,l-1));
457 for (i=l-2; i>=2; i--)
458 {
459 gel(Q,i) = shifti(gel(P,i), ni);
460 ni += n;
461 }
462 Q[1] = P[1]; return Q;
463 }
464
465 /* Return h^deg(P) P(x / h), not memory clean. h integer, P ZX */
466 GEN
ZX_rescale(GEN P,GEN h)467 ZX_rescale(GEN P, GEN h)
468 {
469 long l = lg(P);
470 GEN Q = cgetg(l,t_POL);
471 if (l != 2)
472 {
473 long i = l-1;
474 GEN hi = h;
475 gel(Q,i) = gel(P,i);
476 if (l != 3) { i--; gel(Q,i) = mulii(gel(P,i), h); }
477 for (i--; i>=2; i--) { hi = mulii(hi,h); gel(Q,i) = mulii(gel(P,i), hi); }
478 }
479 Q[1] = P[1]; return Q;
480 }
481 /* Return h^(deg(P)-1) P(x / h), P!=0, h=lt(P), memory unclean; monic result */
482 GEN
ZX_rescale_lt(GEN P)483 ZX_rescale_lt(GEN P)
484 {
485 long l = lg(P);
486 GEN Q = cgetg(l,t_POL);
487 gel(Q,l-1) = gen_1;
488 if (l != 3)
489 {
490 long i = l-1;
491 GEN h = gel(P,i), hi = h;
492 i--; gel(Q,i) = gel(P,i);
493 if (l != 4) { i--; gel(Q,i) = mulii(gel(P,i), h); }
494 for (i--; i>=2; i--) { hi = mulii(hi,h); gel(Q,i) = mulii(gel(P,i), hi); }
495 }
496 Q[1] = P[1]; return Q;
497 }
498
499 /*Eval x in 2^(k*BIL) in linear time*/
500 static GEN
ZX_eval2BILspec(GEN x,long k,long nx)501 ZX_eval2BILspec(GEN x, long k, long nx)
502 {
503 pari_sp av = avma;
504 long i,j, lz = k*nx, ki;
505 GEN pz = cgetipos(2+lz);
506 GEN nz = cgetipos(2+lz);
507 for(i=0; i < lz; i++)
508 {
509 *int_W(pz,i) = 0UL;
510 *int_W(nz,i) = 0UL;
511 }
512 for(i=0, ki=0; i<nx; i++, ki+=k)
513 {
514 GEN c = gel(x,i);
515 long lc = lgefint(c)-2;
516 if (signe(c)==0) continue;
517 if (signe(c) > 0)
518 for (j=0; j<lc; j++) *int_W(pz,j+ki) = *int_W(c,j);
519 else
520 for (j=0; j<lc; j++) *int_W(nz,j+ki) = *int_W(c,j);
521 }
522 pz = int_normalize(pz,0);
523 nz = int_normalize(nz,0); return gerepileuptoint(av, subii(pz,nz));
524 }
525
526 static long
ZX_expispec(GEN x,long nx)527 ZX_expispec(GEN x, long nx)
528 {
529 long i, m = 0;
530 for(i = 0; i < nx; i++)
531 {
532 long e = expi(gel(x,i));
533 if (e > m) m = e;
534 }
535 return m;
536 }
537
538 static GEN
Z_mod2BIL_ZX(GEN x,long bs,long d,long vx)539 Z_mod2BIL_ZX(GEN x, long bs, long d, long vx)
540 {
541 long i, offset, lm = lgefint(x)-2, l = d+vx+3, sx = signe(x);
542 GEN s1 = int2n(bs*BITS_IN_LONG), pol = cgetg(l, t_POL);
543 int carry = 0;
544 pol[1] = evalsigne(1);
545 for (i=0; i<vx; i++) gel(pol,i+2) = gen_0;
546 for (offset=0; i <= d+vx; i++, offset += bs)
547 {
548 pari_sp av = avma;
549 long lz = minss(bs, lm-offset);
550 GEN z = lz > 0 ?adduispec_offset(carry, x, offset, lz): utoi(carry);
551 if (lgefint(z) == 3+bs) { carry = 1; z = gen_0;}
552 else
553 {
554 carry = (lgefint(z) == 2+bs && (HIGHBIT & *int_W(z,bs-1)));
555 if (carry)
556 z = gerepileuptoint(av, (sx==-1)? subii(s1,z): subii(z,s1));
557 else if (sx==-1) togglesign(z);
558 }
559 gel(pol,i+2) = z;
560 }
561 return ZX_renormalize(pol,l);
562 }
563
564 static GEN
ZX_sqrspec_sqri(GEN x,long nx,long ex,long v)565 ZX_sqrspec_sqri(GEN x, long nx, long ex, long v)
566 {
567 long e = 2*ex + expu(nx) + 3;
568 long N = divsBIL(e)+1;
569 GEN z = sqri(ZX_eval2BILspec(x,N,nx));
570 return Z_mod2BIL_ZX(z, N, nx*2-2, v);
571 }
572
573 static GEN
ZX_mulspec_mulii(GEN x,GEN y,long nx,long ny,long ex,long ey,long v)574 ZX_mulspec_mulii(GEN x, GEN y, long nx, long ny, long ex, long ey, long v)
575 {
576 long e = ex + ey + expu(minss(nx,ny)) + 3;
577 long N = divsBIL(e)+1;
578 GEN z = mulii(ZX_eval2BILspec(x,N,nx), ZX_eval2BILspec(y,N,ny));
579 return Z_mod2BIL_ZX(z, N, nx+ny-2, v);
580 }
581
582 INLINE GEN
ZX_sqrspec_basecase_limb(GEN x,long a,long i)583 ZX_sqrspec_basecase_limb(GEN x, long a, long i)
584 {
585 pari_sp av = avma;
586 GEN s = gen_0;
587 long j, l = (i+1)>>1;
588 for (j=a; j<l; j++)
589 {
590 GEN xj = gel(x,j), xx = gel(x,i-j);
591 if (signe(xj) && signe(xx))
592 s = addii(s, mulii(xj, xx));
593 }
594 s = shifti(s,1);
595 if ((i&1) == 0)
596 {
597 GEN t = gel(x, i>>1);
598 if (signe(t))
599 s = addii(s, sqri(t));
600 }
601 return gerepileuptoint(av,s);
602 }
603
604 static GEN
ZX_sqrspec_basecase(GEN x,long nx,long v)605 ZX_sqrspec_basecase(GEN x, long nx, long v)
606 {
607 long i, lz, nz;
608 GEN z;
609
610 lz = (nx << 1) + 1; nz = lz-2;
611 lz += v;
612 z = cgetg(lz,t_POL); z[1] = evalsigne(1); z += 2;
613 for (i=0; i<v; i++) gel(z++, 0) = gen_0;
614 for (i=0; i<nx; i++)
615 gel(z,i) = ZX_sqrspec_basecase_limb(x, 0, i);
616 for ( ; i<nz; i++) gel(z,i) = ZX_sqrspec_basecase_limb(x, i-nx+1, i);
617 z -= v+2; return z;
618 }
619
620 static GEN
Z_sqrshiftspec_ZX(GEN x,long vx)621 Z_sqrshiftspec_ZX(GEN x, long vx)
622 {
623 long i, nz = 2*vx+1;
624 GEN z = cgetg(2+nz, t_POL);
625 z[1] = evalsigne(1);
626 for(i=2;i<nz+1;i++) gel(z,i) = gen_0;
627 gel(z,nz+1) = sqri(x);
628 return z;
629 }
630
631 static GEN
Z_ZX_mulshiftspec(GEN x,GEN y,long ny,long vz)632 Z_ZX_mulshiftspec(GEN x, GEN y, long ny, long vz)
633 {
634 long i, nz = vz+ny;
635 GEN z = cgetg(2+nz, t_POL);
636 z[1] = evalsigne(1);
637 for (i=0; i<vz; i++) gel(z,i+2) = gen_0;
638 for (i=0; i<ny; i++) gel(z,i+vz+2) = mulii(x, gel(y,i));
639 return z;
640 }
641
642 GEN
ZX_sqrspec(GEN x,long nx)643 ZX_sqrspec(GEN x, long nx)
644 {
645 #ifdef PARI_KERNEL_GMP
646 const long low[]={ 17, 32, 96, 112, 160, 128, 128, 160, 160, 160, 160, 160, 176, 192, 192, 192, 192, 192, 224, 224, 224, 240, 240, 240, 272, 288, 288, 240, 288, 304, 304, 304, 304, 304, 304, 352, 352, 368, 352, 352, 352, 368, 368, 432, 432, 496, 432, 496, 496};
647 const long high[]={ 102860, 70254, 52783, 27086, 24623, 18500, 15289, 13899, 12635, 11487, 10442, 9493, 8630, 7845, 7132, 7132, 6484, 6484, 5894, 5894, 4428, 4428, 3660, 4428, 3660, 3660, 2749, 2499, 2272, 2066, 1282, 1282, 1166, 1166, 1166, 1166, 1166, 1166, 1166, 963, 963, 724, 658, 658, 658, 528, 528, 528, 528};
648 #else
649 const long low[]={ 17, 17, 32, 32, 96, 96, 96, 96, 96, 96, 96, 96, 96, 96, 112, 112, 128, 112, 112, 112, 112, 112, 128, 128, 160, 160, 112, 128, 128, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 176, 160, 160, 176, 160, 160, 176, 176, 208, 176, 176, 176, 192, 192, 176, 176, 224, 176, 224, 224, 176, 224, 224, 224, 176, 176, 176, 176, 176, 176, 176, 176, 224, 176, 176, 224, 224, 224, 224, 224, 224, 224, 240, 288, 240, 288, 288, 240, 288, 288, 240, 240, 304, 304};
650 const long high[]={ 165657, 85008, 52783, 43622, 32774, 27086, 22385, 15289, 13899, 12635, 11487, 10442, 9493, 7845, 6484, 6484, 5894, 5894, 4871, 4871, 4428, 4026, 3660, 3660, 3660, 3327, 3327, 3024, 2749, 2749, 2272, 2749, 2499, 2499, 2272, 1878, 1878, 1878, 1707, 1552, 1552, 1552, 1552, 1552, 1411, 1411, 1411, 1282, 1282, 1282, 1282, 1282, 1166, 1166, 1166, 1166, 1166, 1166, 1166, 1166, 1060, 1060, 963, 963, 963, 963, 963, 963, 963, 963, 963, 963, 963, 876, 876, 876, 876, 796, 658, 724, 658, 724, 658, 658, 658, 658, 658, 658, 658, 658, 658, 658, 658, 658, 336, 658, 658, 592, 336, 336};
651 #endif
652 const long nblow = numberof(low);
653 pari_sp av;
654 long ex, vx;
655 GEN z;
656 if (!nx) return pol_0(0);
657 vx = ZX_valspec(x,nx); nx-=vx; x+=vx;
658 if (nx==1) return Z_sqrshiftspec_ZX(gel(x, 0), vx);
659 av = avma;
660 ex = ZX_expispec(x,nx);
661 if (nx-2 < nblow && low[nx-2]<=ex && ex<=high[nx-2])
662 z = ZX_sqrspec_basecase(x, nx, 2*vx);
663 else
664 z = ZX_sqrspec_sqri(x, nx, ex, 2*vx);
665 return gerepileupto(av, z);
666 }
667
668 GEN
ZX_sqr(GEN x)669 ZX_sqr(GEN x)
670 {
671 GEN z = ZX_sqrspec(x+2, lgpol(x));
672 z[1] = x[1];
673 return z;
674 }
675
676 GEN
ZX_mulspec(GEN x,GEN y,long nx,long ny)677 ZX_mulspec(GEN x, GEN y, long nx, long ny)
678 {
679 pari_sp av;
680 long ex, ey, vx, vy;
681 GEN z;
682 if (!nx || !ny) return pol_0(0);
683 vx = ZX_valspec(x,nx); nx-=vx; x += vx;
684 vy = ZX_valspec(y,ny); ny-=vy; y += vy;
685 if (nx==1) return Z_ZX_mulshiftspec(gel(x,0), y, ny, vx+vy);
686 if (ny==1) return Z_ZX_mulshiftspec(gel(y,0), x, nx, vy+vx);
687 av = avma;
688 ex = ZX_expispec(x, nx); ey = ZX_expispec(y, ny);
689 z = ZX_mulspec_mulii(x,y,nx,ny,ex,ey,vx+vy);
690 return gerepileupto(av, z);
691 }
692 GEN
ZX_mul(GEN x,GEN y)693 ZX_mul(GEN x, GEN y)
694 {
695 GEN z;
696 if (x == y) return ZX_sqr(x);
697 z = ZX_mulspec(x+2,y+2,lgpol(x),lgpol(y));
698 z[1] = x[1];
699 if (!signe(y)) z[1] &= VARNBITS;
700 return z;
701 }
702
703 /* x,y two ZX in the same variable; assume y is monic */
704 GEN
ZX_rem(GEN x,GEN y)705 ZX_rem(GEN x, GEN y)
706 {
707 long vx, dx, dy, dz, i, j, sx, lr;
708 pari_sp av0, av;
709 GEN z,p1,rem;
710
711 vx = varn(x);
712 dy = degpol(y);
713 dx = degpol(x);
714 if (dx < dy) return ZX_copy(x);
715 if (!dy) return pol_0(vx); /* y is constant */
716 av0 = avma; dz = dx-dy;
717 z=cgetg(dz+3,t_POL); z[1] = x[1];
718 x += 2; y += 2; z += 2;
719
720 p1 = gel(x,dx);
721 gel(z,dz) = icopy(p1);
722 for (i=dx-1; i>=dy; i--)
723 {
724 av=avma; p1=gel(x,i);
725 for (j=i-dy+1; j<=i && j<=dz; j++)
726 p1 = subii(p1, mulii(gel(z,j),gel(y,i-j)));
727 gel(z,i-dy) = avma == av? icopy(p1): gerepileuptoint(av, p1);
728 }
729 rem = (GEN)avma; av = (pari_sp)new_chunk(dx+3);
730 for (sx=0; ; i--)
731 {
732 p1 = gel(x,i);
733 for (j=0; j<=i && j<=dz; j++)
734 p1 = subii(p1, mulii(gel(z,j),gel(y,i-j)));
735 if (signe(p1)) { sx = 1; break; }
736 if (!i) break;
737 set_avma(av);
738 }
739 lr=i+3; rem -= lr;
740 rem[0] = evaltyp(t_POL) | evallg(lr);
741 rem[1] = z[-1];
742 p1 = gerepileuptoint((pari_sp)rem, p1);
743 rem += 2; gel(rem,i) = p1;
744 for (i--; i>=0; i--)
745 {
746 av=avma; p1 = gel(x,i);
747 for (j=0; j<=i && j<=dz; j++)
748 p1 = subii(p1, mulii(gel(z,j),gel(y,i-j)));
749 gel(rem,i) = avma == av? icopy(p1): gerepileuptoint(av, p1);
750 }
751 rem -= 2;
752 if (!sx) (void)ZX_renormalize(rem, lr);
753 return gerepileupto(av0,rem);
754 }
755
756 /* return x(1) */
757 GEN
ZX_eval1(GEN x)758 ZX_eval1(GEN x)
759 {
760 pari_sp av = avma;
761 long i = lg(x)-1;
762 GEN s;
763 if (i < 2) return gen_0;
764 s = gel(x,i); i--;
765 if (i == 1) return icopy(s);
766 for ( ; i>=2; i--)
767 {
768 GEN c = gel(x,i);
769 if (signe(c)) s = addii(s, c);
770 }
771 return gerepileuptoint(av,s);
772 }
773
774 /* reduce T mod X^n - 1. Shallow function */
775 GEN
ZX_mod_Xnm1(GEN T,ulong n)776 ZX_mod_Xnm1(GEN T, ulong n)
777 {
778 long i, j, L = lg(T), l = n+2;
779 GEN S;
780 if (L <= l) return T;
781 S = cgetg(l, t_POL);
782 S[1] = T[1];
783 for (i = 2; i < l; i++) gel(S,i) = gel(T,i);
784 for (j = 2; i < L; i++) {
785 gel(S,j) = addii(gel(S,j), gel(T,i));
786 if (++j == l) j = 2;
787 }
788 return normalizepol_lg(S, l);
789 }
790
791 /*******************************************************************/
792 /* */
793 /* ZXV */
794 /* */
795 /*******************************************************************/
796
797 int
ZXV_equal(GEN V,GEN W)798 ZXV_equal(GEN V, GEN W)
799 {
800 long l = lg(V);
801 if (l!=lg(W)) return 0;
802 while (--l > 0)
803 if (!ZX_equal(gel(V,l), gel(W,l))) return 0;
804 return 1;
805 }
806
807 GEN
ZXV_Z_mul(GEN x,GEN y)808 ZXV_Z_mul(GEN x, GEN y)
809 { pari_APPLY_same(ZX_Z_mul(gel(x,i), y)) }
810
811 GEN
ZXV_remi2n(GEN x,long N)812 ZXV_remi2n(GEN x, long N)
813 { pari_APPLY_same(ZX_remi2n(gel(x,i), N)) }
814
815 GEN
ZXV_dotproduct(GEN x,GEN y)816 ZXV_dotproduct(GEN x, GEN y)
817 {
818 pari_sp av = avma;
819 long i, lx = lg(x);
820 GEN c;
821 if (lx == 1) return pol_0(varn(x));
822 c = ZX_mul(gel(x,1), gel(y,1));
823 for (i = 2; i < lx; i++)
824 {
825 GEN t = ZX_mul(gel(x,i), gel(y,i));
826 if (signe(t)) c = ZX_add(c, t);
827 }
828 return gerepileupto(av, c);
829 }
830
831 /*******************************************************************/
832 /* */
833 /* ZXQM */
834 /* */
835 /*******************************************************************/
836
837 GEN
ZXn_mul(GEN x,GEN y,long n)838 ZXn_mul(GEN x, GEN y, long n)
839 { return RgXn_red_shallow(ZX_mul(x, y), n); }
840
841 GEN
ZXn_sqr(GEN x,long n)842 ZXn_sqr(GEN x, long n)
843 { return RgXn_red_shallow(ZX_sqr(x), n); }
844
845 /*******************************************************************/
846 /* */
847 /* ZXQM */
848 /* */
849 /*******************************************************************/
850
851 static long
ZX_expi(GEN x)852 ZX_expi(GEN x)
853 {
854 if (signe(x)==0) return 0;
855 if (typ(x)==t_INT) return expi(x);
856 return ZX_expispec(x+2, lgpol(x));
857 }
858
859 static long
ZXC_expi(GEN x)860 ZXC_expi(GEN x)
861 {
862 long i, l = lg(x), m=0;
863 for(i = 1; i < l; i++)
864 {
865 long e = ZX_expi(gel(x,i));
866 if (e > m) m = e;
867 }
868 return m;
869 }
870
871 static long
ZXM_expi(GEN x)872 ZXM_expi(GEN x)
873 {
874 long i, l = lg(x), m=0;
875 for(i = 1; i < l; i++)
876 {
877 long e = ZXC_expi(gel(x,i));
878 if (e > m) m = e;
879 }
880 return m;
881 }
882
883 static GEN
ZX_eval2BIL(GEN x,long k)884 ZX_eval2BIL(GEN x, long k)
885 {
886 if (signe(x)==0) return gen_0;
887 if (typ(x)==t_INT) return x;
888 return ZX_eval2BILspec(x+2, k, lgpol(x));
889 }
890
891 /*Eval x in 2^(k*BIL) in linear time*/
892 static GEN
ZXC_eval2BIL(GEN x,long k)893 ZXC_eval2BIL(GEN x, long k)
894 {
895 long i, lx = lg(x);
896 GEN A = cgetg(lx, t_COL);
897 for (i=1; i<lx; i++) gel(A,i) = ZX_eval2BIL(gel(x,i), k);
898 return A;
899 }
900
901 static GEN
ZXM_eval2BIL(GEN x,long k)902 ZXM_eval2BIL(GEN x, long k)
903 {
904 long i, lx = lg(x);
905 GEN A = cgetg(lx, t_MAT);
906 for (i=1; i<lx; i++) gel(A,i) = ZXC_eval2BIL(gel(x,i), k);
907 return A;
908 }
909
910 static GEN
Z_mod2BIL_ZXQ(GEN x,long bs,GEN T)911 Z_mod2BIL_ZXQ(GEN x, long bs, GEN T)
912 {
913 pari_sp av = avma;
914 long v = varn(T), d = 2*(degpol(T)-1);
915 GEN z = Z_mod2BIL_ZX(x, bs, d, 0);
916 setvarn(z, v);
917 return gerepileupto(av, ZX_rem(z, T));
918 }
919
920 static GEN
ZC_mod2BIL_ZXQC(GEN x,long bs,GEN T)921 ZC_mod2BIL_ZXQC(GEN x, long bs, GEN T)
922 {
923 long i, lx = lg(x);
924 GEN A = cgetg(lx, t_COL);
925 for (i=1; i<lx; i++) gel(A,i) = Z_mod2BIL_ZXQ(gel(x,i), bs, T);
926 return A;
927 }
928
929 static GEN
ZM_mod2BIL_ZXQM(GEN x,long bs,GEN T)930 ZM_mod2BIL_ZXQM(GEN x, long bs, GEN T)
931 {
932 long i, lx = lg(x);
933 GEN A = cgetg(lx, t_MAT);
934 for (i=1; i<lx; i++) gel(A,i) = ZC_mod2BIL_ZXQC(gel(x,i), bs, T);
935 return A;
936 }
937
938 GEN
ZXQM_mul(GEN x,GEN y,GEN T)939 ZXQM_mul(GEN x, GEN y, GEN T)
940 {
941 long d = degpol(T);
942 GEN z;
943 pari_sp av = avma;
944 if (d == 0)
945 z = ZM_mul(simplify_shallow(x),simplify_shallow(y));
946 else
947 {
948 long ex = ZXM_expi(x), ey = ZXM_expi(y), d= degpol(T), n = lg(x)-1;
949 long e = ex + ey + expu(d) + expu(n) + 4;
950 long N = divsBIL(e)+1;
951 z = ZM_mul(ZXM_eval2BIL(x,N), ZXM_eval2BIL(y,N));
952 z = ZM_mod2BIL_ZXQM(z, N, T);
953 }
954 return gerepileupto(av, z);
955 }
956
957 GEN
ZXQM_sqr(GEN x,GEN T)958 ZXQM_sqr(GEN x, GEN T)
959 {
960 long d = degpol(T);
961 GEN z;
962 pari_sp av = avma;
963 if (d == 0)
964 z = ZM_sqr(simplify_shallow(x));
965 else
966 {
967 long ex = ZXM_expi(x), d = degpol(T), n = lg(x)-1;
968 long e = 2*ex + expu(d) + expu(n) + 4;
969 long N = divsBIL(e)+1;
970 z = ZM_sqr(ZXM_eval2BIL(x,N));
971 z = ZM_mod2BIL_ZXQM(z, N, T);
972 }
973 return gerepileupto(av, z);
974 }
975
976 GEN
QXQM_mul(GEN x,GEN y,GEN T)977 QXQM_mul(GEN x, GEN y, GEN T)
978 {
979 GEN dx, nx = Q_primitive_part(x, &dx);
980 GEN dy, ny = Q_primitive_part(y, &dy);
981 GEN z = ZXQM_mul(nx, ny, T);
982 if (dx || dy)
983 {
984 GEN d = dx ? dy ? gmul(dx, dy): dx : dy;
985 if (!gequal1(d)) z = RgM_Rg_mul(z, d);
986 }
987 return z;
988 }
989
990 GEN
QXQM_sqr(GEN x,GEN T)991 QXQM_sqr(GEN x, GEN T)
992 {
993 GEN dx, nx = Q_primitive_part(x, &dx);
994 GEN z = ZXQM_sqr(nx, T);
995 if (dx) z = RgM_Rg_mul(z, gsqr(dx));
996 return z;
997 }
998
999 static GEN
Z_mod2BIL_Fq(GEN x,long bs,GEN T,GEN p)1000 Z_mod2BIL_Fq(GEN x, long bs, GEN T, GEN p)
1001 {
1002 pari_sp av = avma;
1003 long v = get_FpX_var(T), d = 2*(get_FpX_degree(T)-1);
1004 GEN z = Z_mod2BIL_ZX(x, bs, d, 0);
1005 setvarn(z, v);
1006 return gerepileupto(av, FpX_rem(z, T, p));
1007 }
1008
1009 static GEN
ZC_mod2BIL_FqC(GEN x,long bs,GEN T,GEN p)1010 ZC_mod2BIL_FqC(GEN x, long bs, GEN T, GEN p)
1011 {
1012 long i, lx = lg(x);
1013 GEN A = cgetg(lx, t_COL);
1014 for (i=1; i<lx; i++) gel(A,i) = Z_mod2BIL_Fq(gel(x,i), bs, T, p);
1015 return A;
1016 }
1017
1018 static GEN
ZM_mod2BIL_FqM(GEN x,long bs,GEN T,GEN p)1019 ZM_mod2BIL_FqM(GEN x, long bs, GEN T, GEN p)
1020 {
1021 long i, lx = lg(x);
1022 GEN A = cgetg(lx, t_MAT);
1023 for (i=1; i<lx; i++) gel(A,i) = ZC_mod2BIL_FqC(gel(x,i), bs, T, p);
1024 return A;
1025 }
1026
1027 GEN
FqM_mul_Kronecker(GEN x,GEN y,GEN T,GEN p)1028 FqM_mul_Kronecker(GEN x, GEN y, GEN T, GEN p)
1029 {
1030 pari_sp av = avma;
1031 long ex = ZXM_expi(x), ey = ZXM_expi(y), d= degpol(T), n = lg(x)-1;
1032 long e = ex + ey + expu(d) + expu(n) + 4;
1033 long N = divsBIL(e)+1;
1034 GEN z = ZM_mul(ZXM_eval2BIL(x,N), ZXM_eval2BIL(y,N));
1035 return gerepileupto(av, ZM_mod2BIL_FqM(z, N, T, p));
1036 }
1037
1038 /*******************************************************************/
1039 /* */
1040 /* ZXX */
1041 /* */
1042 /*******************************************************************/
1043
1044 void
RgX_check_ZXX(GEN x,const char * s)1045 RgX_check_ZXX(GEN x, const char *s)
1046 {
1047 long k = lg(x)-1;
1048 for ( ; k>1; k--) {
1049 GEN t = gel(x,k);
1050 switch(typ(t)) {
1051 case t_INT: break;
1052 case t_POL: if (RgX_is_ZX(t)) break;
1053 /* fall through */
1054 default: pari_err_TYPE(stack_strcat(s, " not in Z[X,Y]"),x);
1055 }
1056 }
1057 }
1058
1059 /*Renormalize (in place) polynomial with t_INT or ZX coefficients.*/
1060 GEN
ZXX_renormalize(GEN x,long lx)1061 ZXX_renormalize(GEN x, long lx)
1062 {
1063 long i;
1064 for (i = lx-1; i>1; i--)
1065 if (signe(gel(x,i))) break;
1066 stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + (i+1)));
1067 setlg(x, i+1); setsigne(x, i!=1); return x;
1068 }
1069
1070 GEN
ZXX_evalx0(GEN y)1071 ZXX_evalx0(GEN y)
1072 {
1073 long i, l = lg(y);
1074 GEN z = cgetg(l,t_POL); z[1] = y[1];
1075 for(i=2; i<l; i++)
1076 {
1077 GEN yi = gel(y,i);
1078 gel(z,i) = typ(yi)==t_INT? yi: constant_coeff(yi);
1079 }
1080 return ZX_renormalize(z,l);
1081 }
1082
1083 long
ZXX_max_lg(GEN x)1084 ZXX_max_lg(GEN x)
1085 {
1086 long i, prec = 0, lx = lg(x);
1087 for (i=2; i<lx; i++)
1088 {
1089 GEN p = gel(x,i);
1090 long l = (typ(p) == t_INT)? lgefint(p): ZX_max_lg(p);
1091 if (l > prec) prec = l;
1092 }
1093 return prec;
1094 }
1095
1096 GEN
ZXX_Z_mul(GEN y,GEN x)1097 ZXX_Z_mul(GEN y, GEN x)
1098 {
1099 long i, l = lg(y);
1100 GEN z = cgetg(l,t_POL); z[1] = y[1];
1101 for(i=2; i<l; i++)
1102 if(typ(gel(y,i))==t_INT)
1103 gel(z,i) = mulii(gel(y,i),x);
1104 else
1105 gel(z,i) = ZX_Z_mul(gel(y,i),x);
1106 return z;
1107 }
1108
1109 GEN
ZXX_Z_add_shallow(GEN x,GEN y)1110 ZXX_Z_add_shallow(GEN x, GEN y)
1111 {
1112 long i, l = lg(x);
1113 GEN z, a;
1114 if (signe(x)==0) return scalarpol(y,varn(x));
1115 z = cgetg(l,t_POL); z[1] = x[1];
1116 a = gel(x,2);
1117 gel(z, 2) = typ(a)==t_INT? addii(a,y): ZX_Z_add(a,y);
1118 for(i=3; i<l; i++)
1119 gel(z,i) = gel(x,i);
1120 return z;
1121 }
1122
1123 GEN
ZXX_Z_divexact(GEN y,GEN x)1124 ZXX_Z_divexact(GEN y, GEN x)
1125 {
1126 long i, l = lg(y);
1127 GEN z = cgetg(l,t_POL); z[1] = y[1];
1128 for(i=2; i<l; i++)
1129 if(typ(gel(y,i))==t_INT)
1130 gel(z,i) = diviiexact(gel(y,i),x);
1131 else
1132 gel(z,i) = ZX_Z_divexact(gel(y,i),x);
1133 return z;
1134 }
1135
1136 /* Kronecker substitution, ZXX -> ZX:
1137 * P(X,Y) = sum_{0<=i<lP} P_i(X) * Y^i, where deg P_i < n.
1138 * Returns P(X,X^(2n-1)) */
1139 GEN
ZXX_to_Kronecker_spec(GEN P,long lP,long n)1140 ZXX_to_Kronecker_spec(GEN P, long lP, long n)
1141 {
1142 long i, k, N = (n<<1) + 1;
1143 GEN y;
1144 if (!lP) return pol_0(0);
1145 y = cgetg((N-2)*lP + 2, t_POL) + 2;
1146 for (k=i=0; i<lP; i++)
1147 {
1148 long j;
1149 GEN c = gel(P,i);
1150 if (typ(c)!=t_POL)
1151 {
1152 gel(y,k++) = c;
1153 j = 3;
1154 }
1155 else
1156 {
1157 long l = lg(c);
1158 if (l-3 >= n)
1159 pari_err_BUG("ZXX_to_Kronecker, P is not reduced mod Q");
1160 for (j=2; j < l; j++) gel(y,k++) = gel(c,j);
1161 }
1162 if (i == lP-1) break;
1163 for ( ; j < N; j++) gel(y,k++) = gen_0;
1164 }
1165 y-=2; setlg(y, k+2); y[1] = evalsigne(1); return y;
1166 }
1167
1168 /* shallow, n = deg(T) */
1169 GEN
Kronecker_to_ZXX(GEN z,long n,long v)1170 Kronecker_to_ZXX(GEN z, long n, long v)
1171 {
1172 long i,j,lx,l, N = (n<<1)+1;
1173 GEN x, t;
1174 l = lg(z); lx = (l-2) / (N-2);
1175 x = cgetg(lx+3,t_POL);
1176 x[1] = z[1];
1177 for (i=2; i<lx+2; i++)
1178 {
1179 t = cgetg(N,t_POL); t[1] = evalvarn(v);
1180 for (j=2; j<N; j++) gel(t,j) = gel(z,j);
1181 z += (N-2);
1182 gel(x,i) = ZX_renormalize(t,N);
1183 }
1184 N = (l-2) % (N-2) + 2;
1185 t = cgetg(N,t_POL); t[1] = evalvarn(v);
1186 for (j=2; j<N; j++) gel(t,j) = gel(z,j);
1187 gel(x,i) = ZX_renormalize(t,N);
1188 return ZXX_renormalize(x, i+1);
1189 }
1190
1191 GEN
ZXX_to_Kronecker(GEN P,long n)1192 ZXX_to_Kronecker(GEN P, long n)
1193 {
1194 GEN z = ZXX_to_Kronecker_spec(P+2, lgpol(P), n);
1195 setvarn(z,varn(P)); return z;
1196 }
1197 GEN
ZXX_mul_Kronecker(GEN x,GEN y,long n)1198 ZXX_mul_Kronecker(GEN x, GEN y, long n)
1199 { return ZX_mul(ZXX_to_Kronecker(x,n), ZXX_to_Kronecker(y,n)); }
1200
1201 GEN
ZXX_sqr_Kronecker(GEN x,long n)1202 ZXX_sqr_Kronecker(GEN x, long n)
1203 { return ZX_sqr(ZXX_to_Kronecker(x,n)); }
1204
1205 /* shallow, n = deg(T) */
1206 GEN
Kronecker_to_ZXQX(GEN z,GEN T)1207 Kronecker_to_ZXQX(GEN z, GEN T)
1208 {
1209 long i,j,lx,l, N = (degpol(T)<<1)+1;
1210 GEN x, t;
1211 l = lg(z); lx = (l-2) / (N-2);
1212 x = cgetg(lx+3,t_POL);
1213 x[1] = z[1];
1214 for (i=2; i<lx+2; i++)
1215 {
1216 t = cgetg(N,t_POL); t[1] = T[1];
1217 for (j=2; j<N; j++) gel(t,j) = gel(z,j);
1218 z += (N-2);
1219 gel(x,i) = ZX_rem(ZX_renormalize(t,N), T);
1220 }
1221 N = (l-2) % (N-2) + 2;
1222 t = cgetg(N,t_POL); t[1] = T[1];
1223 for (j=2; j<N; j++) gel(t,j) = gel(z,j);
1224 gel(x,i) = ZX_rem(ZX_renormalize(t,N), T);
1225 return ZXX_renormalize(x, i+1);
1226 }
1227
1228 GEN
ZXQX_sqr(GEN x,GEN T)1229 ZXQX_sqr(GEN x, GEN T)
1230 {
1231 pari_sp av = avma;
1232 long n = degpol(T);
1233 GEN z = ZXX_sqr_Kronecker(x, n);
1234 z = Kronecker_to_ZXQX(z, T);
1235 return gerepilecopy(av, z);
1236 }
1237
1238 GEN
ZXQX_mul(GEN x,GEN y,GEN T)1239 ZXQX_mul(GEN x, GEN y, GEN T)
1240 {
1241 pari_sp av = avma;
1242 long n = degpol(T);
1243 GEN z = ZXX_mul_Kronecker(x, y, n);
1244 z = Kronecker_to_ZXQX(z, T);
1245 return gerepilecopy(av, z);
1246 }
1247
1248 GEN
ZXQX_ZXQ_mul(GEN P,GEN U,GEN T)1249 ZXQX_ZXQ_mul(GEN P, GEN U, GEN T)
1250 {
1251 long i, lP;
1252 GEN res;
1253 res = cgetg_copy(P, &lP); res[1] = P[1];
1254 for(i=2; i<lP; i++)
1255 gel(res,i) = typ(gel(P,i))==t_POL? ZXQ_mul(U, gel(P,i), T)
1256 : gmul(U, gel(P,i));
1257 return ZXX_renormalize(res,lP);
1258 }
1259
1260 GEN
QX_mul(GEN x,GEN y)1261 QX_mul(GEN x, GEN y)
1262 {
1263 GEN dx, nx = Q_primitive_part(x, &dx);
1264 GEN dy, ny = Q_primitive_part(y, &dy);
1265 GEN z = ZX_mul(nx, ny);
1266 if (dx || dy)
1267 {
1268 GEN d = dx ? dy ? gmul(dx, dy): dx : dy;
1269 return ZX_Q_mul(z, d);
1270 } else
1271 return z;
1272 }
1273
1274 GEN
QX_sqr(GEN x)1275 QX_sqr(GEN x)
1276 {
1277 GEN dx, nx = Q_primitive_part(x, &dx);
1278 GEN z = ZX_sqr(nx);
1279 if (dx)
1280 return ZX_Q_mul(z, gsqr(dx));
1281 else
1282 return z;
1283 }
1284
1285 GEN
QX_ZX_rem(GEN x,GEN y)1286 QX_ZX_rem(GEN x, GEN y)
1287 {
1288 pari_sp av = avma;
1289 GEN d, nx = Q_primitive_part(x, &d);
1290 GEN r = ZX_rem(nx, y);
1291 if (d) r = ZX_Q_mul(r, d);
1292 return gerepileupto(av, r);
1293 }
1294
1295 GEN
QXQX_mul(GEN x,GEN y,GEN T)1296 QXQX_mul(GEN x, GEN y, GEN T)
1297 {
1298 GEN dx, nx = Q_primitive_part(x, &dx);
1299 GEN dy, ny = Q_primitive_part(y, &dy);
1300 GEN z = ZXQX_mul(nx, ny, T);
1301 if (dx || dy)
1302 {
1303 GEN d = dx ? dy ? gmul(dx, dy): dx : dy;
1304 return ZXX_Q_mul(z, d);
1305 } else
1306 return z;
1307 }
1308
1309 GEN
QXQX_sqr(GEN x,GEN T)1310 QXQX_sqr(GEN x, GEN T)
1311 {
1312 GEN dx, nx = Q_primitive_part(x, &dx);
1313 GEN z = ZXQX_sqr(nx, T);
1314 if (dx)
1315 return ZXX_Q_mul(z, gsqr(dx));
1316 else
1317 return z;
1318 }
1319
1320 GEN
QXQX_powers(GEN P,long n,GEN T)1321 QXQX_powers(GEN P, long n, GEN T)
1322 {
1323 GEN v = cgetg(n+2, t_VEC);
1324 long i;
1325 gel(v, 1) = pol_1(varn(T));
1326 if (n==0) return v;
1327 gel(v, 2) = gcopy(P);
1328 for (i = 2; i <= n; i++) gel(v,i+1) = QXQX_mul(P, gel(v,i), T);
1329 return v;
1330 }
1331
1332 GEN
QXQX_QXQ_mul(GEN P,GEN U,GEN T)1333 QXQX_QXQ_mul(GEN P, GEN U, GEN T)
1334 {
1335 long i, lP;
1336 GEN res;
1337 res = cgetg_copy(P, &lP); res[1] = P[1];
1338 for(i=2; i<lP; i++)
1339 gel(res,i) = typ(gel(P,i))==t_POL? QXQ_mul(U, gel(P,i), T)
1340 : gmul(U, gel(P,i));
1341 return ZXX_renormalize(res,lP);
1342 }
1343