Code coverage tests

This page documents the degree to which the PARI/GP source code is tested by our public test suite, distributed with the source distribution in directory src/test/. This is measured by the gcov utility; we then process gcov output using the lcov frond-end.

We test a few variants depending on Configure flags on the pari.math.u-bordeaux.fr machine (x86_64 architecture), and agregate them in the final report:

The target is to exceed 90% coverage for all mathematical modules (given that branches depending on DEBUGLEVEL or DEBUGMEM are not covered). This script is run to produce the results below.

LCOV - code coverage report
Current view: top level - basemath - ZX.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.12.0 lcov report (development 23712-7b25a218b) Lines: 559 596 93.8 %
Date: 2019-03-24 05:44:59 Functions: 77 81 95.1 %
Legend: Lines: hit not hit

          Line data    Source code
       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. 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             : #include "pari.h"
      15             : #include "paripriv.h"
      16             : 
      17             : /*******************************************************************/
      18             : /*                                                                 */
      19             : /*                                ZX                               */
      20             : /*                                                                 */
      21             : /*******************************************************************/
      22             : void
      23      245700 : RgX_check_QX(GEN x, const char *s)
      24      245700 : { if (!RgX_is_QX(x)) pari_err_TYPE(stack_strcat(s," [not in Q[X]]"), x); }
      25             : void
      26      636699 : RgX_check_ZX(GEN x, const char *s)
      27      636699 : { if (!RgX_is_ZX(x)) pari_err_TYPE(stack_strcat(s," [not in Z[X]]"), x); }
      28             : long
      29    20023724 : ZX_max_lg(GEN x)
      30             : {
      31    20023724 :   long i, prec = 0, lx = lg(x);
      32             : 
      33    20023724 :   for (i=2; i<lx; i++) { long l = lgefint(gel(x,i)); if (l > prec) prec = l; }
      34    20023724 :   return prec;
      35             : }
      36             : 
      37             : GEN
      38     8429436 : ZX_add(GEN x, GEN y)
      39             : {
      40             :   long lx,ly,i;
      41             :   GEN z;
      42     8429436 :   lx = lg(x); ly = lg(y); if (lx < ly) swapspec(x,y, lx,ly);
      43     8429436 :   z = cgetg(lx,t_POL); z[1] = x[1];
      44     8431564 :   for (i=2; i<ly; i++) gel(z,i) = addii(gel(x,i),gel(y,i));
      45     8426555 :   for (   ; i<lx; i++) gel(z,i) = icopy(gel(x,i));
      46     8429423 :   if (lx == ly) z = ZX_renormalize(z, lx);
      47     8429429 :   if (!lgpol(z)) { set_avma((pari_sp)(z + lx)); return pol_0(varn(x)); }
      48     7195483 :   return z;
      49             : }
      50             : 
      51             : GEN
      52      603317 : ZX_sub(GEN x,GEN y)
      53             : {
      54      603317 :   long i, lx = lg(x), ly = lg(y);
      55             :   GEN z;
      56      603317 :   if (lx >= ly)
      57             :   {
      58      591744 :     z = cgetg(lx,t_POL); z[1] = x[1];
      59      591862 :     for (i=2; i<ly; i++) gel(z,i) = subii(gel(x,i),gel(y,i));
      60      591742 :     if (lx == ly)
      61             :     {
      62      563204 :       z = ZX_renormalize(z, lx);
      63      563206 :       if (!lgpol(z)) { set_avma((pari_sp)(z + lx)); z = pol_0(varn(x)); }
      64             :     }
      65             :     else
      66       28538 :       for (   ; i<lx; i++) gel(z,i) = icopy(gel(x,i));
      67             :   }
      68             :   else
      69             :   {
      70       11573 :     z = cgetg(ly,t_POL); z[1] = y[1];
      71       11573 :     for (i=2; i<lx; i++) gel(z,i) = subii(gel(x,i),gel(y,i));
      72       11573 :     for (   ; i<ly; i++) gel(z,i) = negi(gel(y,i));
      73             :   }
      74      603315 :   return z;
      75             : }
      76             : 
      77             : GEN
      78      318122 : ZX_neg(GEN x)
      79             : {
      80      318122 :   long i, l = lg(x);
      81      318122 :   GEN y = cgetg(l,t_POL);
      82      318122 :   y[1] = x[1]; for(i=2; i<l; i++) gel(y,i) = negi(gel(x,i));
      83      318122 :   return y;
      84             : }
      85             : GEN
      86     4825597 : ZX_copy(GEN x)
      87             : {
      88     4825597 :   long i, l = lg(x);
      89     4825597 :   GEN y = cgetg(l, t_POL);
      90     4825663 :   y[1] = x[1];
      91    15602482 :   for (i=2; i<l; i++)
      92             :   {
      93    10776821 :     GEN c = gel(x,i);
      94    10776821 :     gel(y,i) = lgefint(c) == 2? gen_0: icopy(c);
      95             :   }
      96     4825661 :   return y;
      97             : }
      98             : 
      99             : GEN
     100       54282 : scalar_ZX(GEN x, long v)
     101             : {
     102             :   GEN z;
     103       54282 :   if (!signe(x)) return pol_0(v);
     104        4316 :   z = cgetg(3, t_POL);
     105        4316 :   z[1] = evalsigne(1) | evalvarn(v);
     106        4316 :   gel(z,2) = icopy(x); return z;
     107             : }
     108             : 
     109             : GEN
     110       19134 : scalar_ZX_shallow(GEN x, long v)
     111             : {
     112             :   GEN z;
     113       19134 :   if (!signe(x)) return pol_0(v);
     114       18854 :   z = cgetg(3, t_POL);
     115       18854 :   z[1] = evalsigne(1) | evalvarn(v);
     116       18854 :   gel(z,2) = x; return z;
     117             : }
     118             : 
     119             : GEN
     120       83509 : ZX_Z_add(GEN y, GEN x)
     121             : {
     122             :   long lz, i;
     123       83509 :   GEN z = cgetg_copy(y, &lz);
     124       83513 :   if (lz == 2) { set_avma((pari_sp)(z + 2)); return scalar_ZX(x,varn(y)); }
     125       72845 :   z[1] = y[1];
     126       72845 :   gel(z,2) = addii(gel(y,2),x);
     127       71647 :   for(i=3; i<lz; i++) gel(z,i) = icopy(gel(y,i));
     128       72846 :   if (lz==3) z = ZX_renormalize(z,lz);
     129       72846 :   return z;
     130             : }
     131             : GEN
     132       10778 : ZX_Z_add_shallow(GEN y, GEN x)
     133             : {
     134             :   long lz, i;
     135       10778 :   GEN z = cgetg_copy(y, &lz);
     136       10778 :   if (lz == 2) { set_avma((pari_sp)(z + 2)); return scalar_ZX_shallow(x,varn(y)); }
     137       10568 :   z[1] = y[1];
     138       10568 :   gel(z,2) = addii(gel(y,2),x);
     139       10568 :   for(i=3; i<lz; i++) gel(z,i) = gel(y,i);
     140       10568 :   if (lz==3) z = ZX_renormalize(z,lz);
     141       10568 :   return z;
     142             : }
     143             : 
     144             : GEN
     145       34559 : ZX_Z_sub(GEN y, GEN x)
     146             : {
     147             :   long lz, i;
     148       34559 :   GEN z = cgetg_copy(y, &lz);
     149       34559 :   if (lz == 2)
     150             :   { /* scalarpol(negi(x), v) */
     151           0 :     long v = varn(y);
     152           0 :     set_avma((pari_sp)(z + 2));
     153           0 :     if (!signe(x)) return pol_0(v);
     154           0 :     z = cgetg(3,t_POL);
     155           0 :     z[1] = evalvarn(v) | evalsigne(1);
     156           0 :     gel(z,2) = negi(x); return z;
     157             :   }
     158       34559 :   z[1] = y[1];
     159       34559 :   gel(z,2) = subii(gel(y,2),x);
     160       34559 :   for(i=3; i<lz; i++) gel(z,i) = icopy(gel(y,i));
     161       34559 :   if (lz==3) z = ZX_renormalize(z,lz);
     162       34559 :   return z;
     163             : }
     164             : 
     165             : GEN
     166      480347 : Z_ZX_sub(GEN x, GEN y)
     167             : {
     168             :   long lz, i;
     169      480347 :   GEN z = cgetg_copy(y, &lz);
     170      480347 :   if (lz == 2) { set_avma((pari_sp)(z + 2)); return scalar_ZX(x,varn(y)); }
     171      480347 :   z[1] = y[1];
     172      480347 :   gel(z,2) = subii(x, gel(y,2));
     173      480347 :   for(i=3; i<lz; i++) gel(z,i) = negi(gel(y,i));
     174      480347 :   if (lz==3) z = ZX_renormalize(z,lz);
     175      480347 :   return z;
     176             : }
     177             : 
     178             : GEN
     179     2488654 : ZX_Z_divexact(GEN y,GEN x)
     180             : {
     181     2488654 :   long i, l = lg(y);
     182     2488654 :   GEN z = cgetg(l,t_POL); z[1] = y[1];
     183     2489472 :   for(i=2; i<l; i++) gel(z,i) = diviiexact(gel(y,i),x);
     184     2488643 :   return z;
     185             : }
     186             : 
     187             : GEN
     188     8932592 : ZX_Z_mul(GEN y,GEN x)
     189             : {
     190             :   GEN z;
     191             :   long i, l;
     192     8932592 :   if (!signe(x)) return pol_0(varn(y));
     193     5878958 :   l = lg(y); z = cgetg(l,t_POL); z[1] = y[1];
     194     5880174 :   for(i=2; i<l; i++) gel(z,i) = mulii(gel(y,i),x);
     195     5878938 :   return z;
     196             : }
     197             : 
     198             : GEN
     199      120491 : ZX_mulu(GEN y, ulong x)
     200             : {
     201             :   GEN z;
     202             :   long i, l;
     203      120491 :   if (!x) return pol_0(varn(y));
     204      120491 :   l = lg(y); z = cgetg(l,t_POL); z[1] = y[1];
     205      120491 :   for(i=2; i<l; i++) gel(z,i) = mului(x,gel(y,i));
     206      120491 :   return z;
     207             : }
     208             : 
     209             : GEN
     210      181600 : ZX_shifti(GEN y, long n)
     211             : {
     212             :   GEN z;
     213             :   long i, l;
     214      181600 :   l = lg(y); z = cgetg(l,t_POL); z[1] = y[1];
     215      203120 :   for(i=2; i<l; i++) gel(z,i) = shifti(gel(y,i),n);
     216      181481 :   return ZX_renormalize(z,l);
     217             : }
     218             : 
     219             : GEN
     220       97347 : ZX_remi2n(GEN y, long n)
     221             : {
     222             :   GEN z;
     223             :   long i, l;
     224       97347 :   l = lg(y); z = cgetg(l,t_POL); z[1] = y[1];
     225      104528 :   for(i=2; i<l; i++) gel(z,i) = remi2n(gel(y,i),n);
     226       97260 :   return ZX_renormalize(z,l);
     227             : }
     228             : 
     229             : GEN
     230       39001 : ZXT_remi2n(GEN z, long n)
     231             : {
     232       39001 :   if (typ(z) == t_POL)
     233       35070 :     return ZX_remi2n(z, n);
     234             :   else
     235             :   {
     236        3931 :     long i,l = lg(z);
     237        3931 :     GEN x = cgetg(l, t_VEC);
     238        3930 :     for (i=1; i<l; i++) gel(x,i) = ZXT_remi2n(gel(z,i), n);
     239        3925 :     return x;
     240             :   }
     241             : }
     242             : 
     243             : GEN
     244         294 : zx_to_ZX(GEN z)
     245             : {
     246         294 :   long i, l = lg(z);
     247         294 :   GEN x = cgetg(l,t_POL);
     248         294 :   for (i=2; i<l; i++) gel(x,i) = stoi(z[i]);
     249         294 :   x[1] = evalsigne(l-2!=0)| z[1]; return x;
     250             : }
     251             : 
     252             : GEN
     253     2016310 : ZX_deriv(GEN x)
     254             : {
     255     2016310 :   long i,lx = lg(x)-1;
     256             :   GEN y;
     257             : 
     258     2016310 :   if (lx<3) return pol_0(varn(x));
     259     2008588 :   y = cgetg(lx,t_POL);
     260     2008628 :   for (i=2; i<lx ; i++) gel(y,i) = mului(i-1,gel(x,i+1));
     261     2008583 :   y[1] = x[1]; return y;
     262             : }
     263             : 
     264             : int
     265     1521419 : ZX_equal(GEN V, GEN W)
     266             : {
     267     1521419 :   long i, l = lg(V);
     268     1521419 :   if (lg(W) != l) return 0;
     269     4997728 :   for (i = 2; i < l; i++)
     270     3740983 :     if (!equalii(gel(V,i), gel(W,i))) return 0;
     271     1256745 :   return 1;
     272             : }
     273             : 
     274             : static long
     275   106610105 : ZX_valspec(GEN x, long nx)
     276             : {
     277             :   long vx;
     278   168248450 :   for (vx = 0; vx<nx ; vx++)
     279   168248413 :     if (signe(gel(x,vx))) break;
     280   106610105 :   return vx;
     281             : }
     282             : 
     283             : long
     284       30307 : ZX_val(GEN x)
     285             : {
     286             :   long vx;
     287       30307 :   if (!signe(x)) return LONG_MAX;
     288       31875 :   for (vx = 0;; vx++)
     289       33443 :     if (signe(gel(x,2+vx))) break;
     290       30307 :   return vx;
     291             : }
     292             : long
     293    21243463 : ZX_valrem(GEN x, GEN *Z)
     294             : {
     295             :   long vx;
     296    21243463 :   if (!signe(x)) { *Z = pol_0(varn(x)); return LONG_MAX; }
     297    41170381 :   for (vx = 0;; vx++)
     298    61097299 :     if (signe(gel(x,2+vx))) break;
     299    21243463 :   *Z = RgX_shift_shallow(x, -vx);
     300    21243463 :   return vx;
     301             : }
     302             : 
     303             : GEN
     304          21 : ZX_div_by_X_1(GEN a, GEN *r)
     305             : {
     306          21 :   long l = lg(a), i;
     307          21 :   GEN a0, z0, z = cgetg(l-1, t_POL);
     308          21 :   z[1] = a[1];
     309          21 :   a0 = a + l-1;
     310          21 :   z0 = z + l-2; *z0 = *a0--;
     311        4732 :   for (i=l-3; i>1; i--) /* z[i] = a[i+1] + z[i+1] */
     312             :   {
     313        4711 :     GEN t = addii(gel(a0--,0), gel(z0--,0));
     314        4711 :     gel(z0,0) = t;
     315             :   }
     316          21 :   if (r) *r = addii(gel(a0,0), gel(z0,0));
     317          21 :   return z;
     318             : }
     319             : 
     320             : /* Return 2^(n degpol(P))  P(x >> n) */
     321             : GEN
     322      149041 : ZX_rescale2n(GEN P, long n)
     323             : {
     324      149041 :   long i, l = lg(P), ni = n;
     325             :   GEN Q;
     326      149041 :   if (l==2) return pol_0(varn(P));
     327      149041 :   Q = cgetg(l,t_POL);
     328      149041 :   gel(Q,l-1) = icopy(gel(P,l-1));
     329     1185396 :   for (i=l-2; i>=2; i--)
     330             :   {
     331     1036355 :     gel(Q,i) = shifti(gel(P,i), ni);
     332     1036355 :     ni += n;
     333             :   }
     334      149041 :   Q[1] = P[1]; return Q;
     335             : }
     336             : 
     337             : /* Return h^deg(P) P(x / h), not memory clean. h integer, P ZX */
     338             : GEN
     339       11939 : ZX_rescale(GEN P, GEN h)
     340             : {
     341       11939 :   long l = lg(P);
     342       11939 :   GEN Q = cgetg(l,t_POL);
     343       11939 :   if (l != 2)
     344             :   {
     345       11939 :     long i = l-1;
     346       11939 :     GEN hi = h;
     347       11939 :     gel(Q,i) = gel(P,i);
     348       11939 :     if (l != 3) { i--; gel(Q,i) = mulii(gel(P,i), h); }
     349       11939 :     for (i--; i>=2; i--) { hi = mulii(hi,h); gel(Q,i) = mulii(gel(P,i), hi); }
     350             :   }
     351       11939 :   Q[1] = P[1]; return Q;
     352             : }
     353             : /* Return h^(deg(P)-1) P(x / h), P!=0, h=lt(P), memory unclean; monic result */
     354             : GEN
     355           0 : ZX_rescale_lt(GEN P)
     356             : {
     357           0 :   long l = lg(P);
     358           0 :   GEN Q = cgetg(l,t_POL);
     359           0 :   gel(Q,l-1) = gen_1;
     360           0 :   if (l != 3)
     361             :   {
     362           0 :     long i = l-1;
     363           0 :     GEN h = gel(P,i), hi = h;
     364           0 :     i--; gel(Q,i) = gel(P,i);
     365           0 :     if (l != 4) { i--; gel(Q,i) = mulii(gel(P,i), h); }
     366           0 :     for (i--; i>=2; i--) { hi = mulii(hi,h); gel(Q,i) = mulii(gel(P,i), hi); }
     367             :   }
     368           0 :   Q[1] = P[1]; return Q;
     369             : }
     370             : 
     371             : 
     372             : /*Eval x in 2^(k*BIL) in linear time*/
     373             : static GEN
     374    53434118 : ZX_eval2BILspec(GEN x, long k, long nx)
     375             : {
     376    53434118 :   pari_sp av = avma;
     377    53434118 :   long i,j, lz = k*nx, ki;
     378    53434118 :   GEN pz = cgetipos(2+lz);
     379    53432091 :   GEN nz = cgetipos(2+lz);
     380  1554243869 :   for(i=0; i < lz; i++)
     381             :   {
     382  1500808710 :     *int_W(pz,i) = 0UL;
     383  1500808710 :     *int_W(nz,i) = 0UL;
     384             :   }
     385   448710593 :   for(i=0, ki=0; i<nx; i++, ki+=k)
     386             :   {
     387   395275434 :     GEN c = gel(x,i);
     388   395275434 :     long lc = lgefint(c)-2;
     389   395275434 :     if (signe(c)==0) continue;
     390   309370286 :     if (signe(c) > 0)
     391   246713296 :       for (j=0; j<lc; j++) *int_W(pz,j+ki) = *int_W(c,j);
     392             :     else
     393    62656990 :       for (j=0; j<lc; j++) *int_W(nz,j+ki) = *int_W(c,j);
     394             :   }
     395    53435159 :   pz = int_normalize(pz,0);
     396    53435257 :   nz = int_normalize(nz,0); return gerepileuptoint(av, subii(pz,nz));
     397             : }
     398             : 
     399             : static long
     400    57585610 : ZX_expispec(GEN x, long nx)
     401             : {
     402    57585610 :   long i, m = 0;
     403   474032792 :   for(i = 0; i < nx; i++)
     404             :   {
     405   416450477 :     long e = expi(gel(x,i));
     406   416447182 :     if (e > m) m = e;
     407             :   }
     408    57582315 :   return m;
     409             : }
     410             : 
     411             : static GEN
     412    28785725 : Z_mod2BIL_ZX(GEN x, long bs, long d, long vx)
     413             : {
     414    28785725 :   long i, offset, lm = lgefint(x)-2, l = d+vx+3, sx = signe(x);
     415    28785725 :   GEN s1 = int2n(bs*BITS_IN_LONG), pol = cgetg(l, t_POL);
     416    28878219 :   int carry = 0;
     417    28878219 :   pol[1] = evalsigne(1);
     418    28878219 :   for (i=0; i<vx; i++) gel(pol,i+2) = gen_0;
     419   425193005 :   for (offset=0; i <= d+vx; i++, offset += bs)
     420             :   {
     421   396412165 :     pari_sp av = avma;
     422   396412165 :     long lz = minss(bs, lm-offset);
     423   397553395 :     GEN z = lz > 0 ?adduispec_offset(carry, x, offset, lz): utoi(carry);
     424   396441883 :     if (lgefint(z) == 3+bs) { carry = 1; z = gen_0;}
     425             :     else
     426             :     {
     427   385295336 :       carry = (lgefint(z) == 2+bs && (HIGHBIT & *int_W(z,bs-1)));
     428   385295336 :       if (carry)
     429    29557655 :         z = gerepileuptoint(av, (sx==-1)? subii(s1,z): subii(z,s1));
     430   355737681 :       else if (sx==-1) togglesign(z);
     431             :     }
     432   396314786 :     gel(pol,i+2) = z;
     433             :   }
     434    28780840 :   return ZX_renormalize(pol,l);
     435             : }
     436             : 
     437             : static GEN
     438     2491503 : ZX_sqrspec_sqri(GEN x, long nx, long ex, long v)
     439             : {
     440     2491503 :   long e = 2*ex + expu(nx) + 3;
     441     2491488 :   long N = divsBIL(e)+1;
     442     2491469 :   GEN  z = sqri(ZX_eval2BILspec(x,N,nx));
     443     2492671 :   return Z_mod2BIL_ZX(z, N, nx*2-2, v);
     444             : }
     445             : 
     446             : static GEN
     447    24628039 : ZX_mulspec_mulii(GEN x, GEN y, long nx, long ny, long ex, long ey, long v)
     448             : {
     449    24628039 :   long e = ex + ey + expu(minss(nx,ny)) + 3;
     450    24627947 :   long N = divsBIL(e)+1;
     451    24627909 :   GEN  z = mulii(ZX_eval2BILspec(x,N,nx), ZX_eval2BILspec(y,N,ny));
     452    24627607 :   return Z_mod2BIL_ZX(z, N, nx+ny-2, v);
     453             : }
     454             : 
     455             : INLINE GEN
     456    38508465 : ZX_sqrspec_basecase_limb(GEN x, long a, long i)
     457             : {
     458    38508465 :   pari_sp av = avma;
     459    38508465 :   GEN s = gen_0;
     460    38508465 :   long j, l = (i+1)>>1;
     461   146393475 :   for (j=a; j<l; j++)
     462             :   {
     463   108475974 :     GEN xj = gel(x,j), xx = gel(x,i-j);
     464   108475974 :     if (signe(xj) && signe(xx))
     465   106799527 :       s = addii(s, mulii(xj, xx));
     466             :   }
     467    37917501 :   s = shifti(s,1);
     468    37867762 :   if ((i&1) == 0)
     469             :   {
     470    21128925 :     GEN t = gel(x, i>>1);
     471    21128925 :     if (signe(t))
     472    21005292 :       s = addii(s, sqri(t));
     473             :   }
     474    37852822 :   return gerepileuptoint(av,s);
     475             : }
     476             : 
     477             : static GEN
     478     4148708 : ZX_sqrspec_basecase(GEN x, long nx, long v)
     479             : {
     480             :   long i, lz, nz;
     481             :   GEN z;
     482             : 
     483     4148708 :   lz = (nx << 1) + 1; nz = lz-2;
     484     4148708 :   lz += v;
     485     4148708 :   z = cgetg(lz,t_POL); z[1] = evalsigne(1); z += 2;
     486     4150965 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
     487    25299373 :   for (i=0; i<nx; i++)
     488    21155717 :     gel(z,i) = ZX_sqrspec_basecase_limb(x, 0, i);
     489     4143656 :   for (  ; i<nz; i++) gel(z,i) = ZX_sqrspec_basecase_limb(x, i-nx+1, i);
     490     4147887 :   z -= v+2; return z;
     491             : }
     492             : 
     493             : static GEN
     494      491086 : Z_sqrshiftspec_ZX(GEN x, long vx)
     495             : {
     496      491086 :   long i, nz = 2*vx+1;
     497      491086 :   GEN z = cgetg(2+nz, t_POL);
     498      491086 :   z[1] = evalsigne(1);
     499      491086 :   for(i=2;i<nz+1;i++) gel(z,i) = gen_0;
     500      491086 :   gel(z,nz+1) = sqri(x);
     501      491086 :   return z;
     502             : }
     503             : 
     504             : static GEN
     505    25112116 : Z_ZX_mulshiftspec(GEN x, GEN y, long ny, long vz)
     506             : {
     507    25112116 :   long i, nz = vz+ny;
     508    25112116 :   GEN z = cgetg(2+nz, t_POL);
     509    25112180 :   z[1] = evalsigne(1);
     510    25112180 :   for (i=0; i<vz; i++)   gel(z,i+2)    = gen_0;
     511    25112180 :   for (i=0; i<ny; i++) gel(z,i+vz+2) = mulii(x, gel(y,i));
     512    25112081 :   return z;
     513             : }
     514             : 
     515             : GEN
     516     7147699 : ZX_sqrspec(GEN x, long nx)
     517             : {
     518             : #ifdef PARI_KERNEL_GMP
     519     4385443 :   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};
     520     4385443 :   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};
     521             : #else
     522     2762256 :   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};
     523     2762256 :   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};
     524             : #endif
     525     7147699 :   const long nblow = numberof(low);
     526             :   pari_sp av;
     527             :   long ex, vx;
     528             :   GEN z;
     529     7147699 :   if (!nx) return pol_0(0);
     530     7132268 :   vx = ZX_valspec(x,nx); nx-=vx; x+=vx;
     531     7132341 :   if (nx==1) return Z_sqrshiftspec_ZX(gel(x, 0), vx);
     532     6641255 :   av = avma;
     533     6641255 :   ex = ZX_expispec(x,nx);
     534     6640036 :   if (nx-2 < nblow && low[nx-2]<=ex && ex<=high[nx-2])
     535     4148707 :     z = ZX_sqrspec_basecase(x, nx, 2*vx);
     536             :   else
     537     2491329 :     z = ZX_sqrspec_sqri(x, nx, ex, 2*vx);
     538     6638425 :   return gerepileupto(av, z);
     539             : }
     540             : 
     541             : GEN
     542     7145955 : ZX_sqr(GEN x)
     543             : {
     544     7145955 :   GEN z = ZX_sqrspec(x+2, lgpol(x));
     545     7144886 :   z[1] = x[1];
     546     7144886 :   return z;
     547             : }
     548             : 
     549             : GEN
     550    60124708 : ZX_mulspec(GEN x, GEN y, long nx, long ny)
     551             : {
     552             :   pari_sp av;
     553             :   long ex, ey, vx, vy;
     554             :   GEN z;
     555    60124708 :   if (!nx || !ny) return pol_0(0);
     556    49739840 :   vx = ZX_valspec(x,nx); nx-=vx; x += vx;
     557    49739791 :   vy = ZX_valspec(y,ny); ny-=vy; y += vy;
     558    49739915 :   if (nx==1) return Z_ZX_mulshiftspec(gel(x,0), y, ny, vx+vy);
     559    30291186 :   if (ny==1) return Z_ZX_mulshiftspec(gel(y,0), x, nx, vy+vx);
     560    24627801 :   av = avma;
     561    24627801 :   ex = ZX_expispec(x, nx); ey = ZX_expispec(y, ny);
     562    24627891 :   z  = ZX_mulspec_mulii(x,y,nx,ny,ex,ey,vx+vy);
     563    24626762 :   return gerepileupto(av, z);
     564             : }
     565             : GEN
     566    58905448 : ZX_mul(GEN x, GEN y)
     567             : {
     568             :   GEN z;
     569    58905448 :   if (x == y) return ZX_sqr(x);
     570    58180421 :   z = ZX_mulspec(x+2,y+2,lgpol(x),lgpol(y));
     571    58180529 :   z[1] = x[1];
     572    58180529 :   if (!signe(y)) z[1] &= VARNBITS;
     573    58180529 :   return z;
     574             : }
     575             : 
     576             : /* x,y two ZX in the same variable; assume y is monic */
     577             : GEN
     578     7429558 : ZX_rem(GEN x, GEN y)
     579             : {
     580             :   long vx, dx, dy, dz, i, j, sx, lr;
     581             :   pari_sp av0, av;
     582             :   GEN z,p1,rem;
     583             : 
     584     7429558 :   vx = varn(x);
     585     7429558 :   dy = degpol(y);
     586     7429460 :   dx = degpol(x);
     587     7429468 :   if (dx < dy) return ZX_copy(x);
     588     2794317 :   if (!dy) return pol_0(vx); /* y is constant */
     589     2794282 :   av0 = avma; dz = dx-dy;
     590     2794282 :   z=cgetg(dz+3,t_POL); z[1] = x[1];
     591     2794118 :   x += 2; y += 2; z += 2;
     592             : 
     593     2794118 :   p1 = gel(x,dx);
     594     2794118 :   gel(z,dz) = icopy(p1);
     595    13500966 :   for (i=dx-1; i>=dy; i--)
     596             :   {
     597    10707078 :     av=avma; p1=gel(x,i);
     598   106962989 :     for (j=i-dy+1; j<=i && j<=dz; j++)
     599    96267824 :       p1 = subii(p1, mulii(gel(z,j),gel(y,i-j)));
     600    10695165 :     gel(z,i-dy) = avma == av? icopy(p1): gerepileuptoint(av, p1);
     601             :   }
     602     2793888 :   rem = (GEN)avma; av = (pari_sp)new_chunk(dx+3);
     603     3887117 :   for (sx=0; ; i--)
     604             :   {
     605     4976034 :     p1 = gel(x,i);
     606    27451214 :     for (j=0; j<=i && j<=dz; j++)
     607    23568609 :       p1 = subii(p1, mulii(gel(z,j),gel(y,i-j)));
     608     3882605 :     if (signe(p1)) { sx = 1; break; }
     609     1217201 :     if (!i) break;
     610     1088919 :     set_avma(av);
     611             :   }
     612     2793686 :   lr=i+3; rem -= lr;
     613     2793686 :   rem[0] = evaltyp(t_POL) | evallg(lr);
     614     2793555 :   rem[1] = z[-1];
     615     2793555 :   p1 = gerepileuptoint((pari_sp)rem, p1);
     616     2794981 :   rem += 2; gel(rem,i) = p1;
     617    16065362 :   for (i--; i>=0; i--)
     618             :   {
     619    13271588 :     av=avma; p1 = gel(x,i);
     620   122824526 :     for (j=0; j<=i && j<=dz; j++)
     621   109580573 :       p1 = subii(p1, mulii(gel(z,j),gel(y,i-j)));
     622    13243953 :     gel(rem,i) = avma == av? icopy(p1): gerepileuptoint(av, p1);
     623             :   }
     624     2793774 :   rem -= 2;
     625     2793774 :   if (!sx) (void)ZX_renormalize(rem, lr);
     626     2793774 :   return gerepileupto(av0,rem);
     627             : }
     628             : 
     629             : /* return x(1) */
     630             : GEN
     631        4718 : ZX_eval1(GEN x)
     632             : {
     633        4718 :   pari_sp av = avma;
     634        4718 :   long i = lg(x)-1;
     635             :   GEN s;
     636        4718 :   if (i < 2) return gen_0;
     637        4193 :   s = gel(x,i); i--;
     638        4193 :   if (i == 1) return icopy(s);
     639       34608 :   for ( ; i>=2; i--)
     640             :   {
     641       31423 :     GEN c = gel(x,i);
     642       31423 :     if (signe(c)) s = addii(s, c);
     643             :   }
     644        3185 :   return gerepileuptoint(av,s);
     645             : }
     646             : 
     647             : /* reduce T mod X^n - 1. Shallow function */
     648             : GEN
     649     1479987 : ZX_mod_Xnm1(GEN T, ulong n)
     650             : {
     651     1479987 :   long i, j, L = lg(T), l = n+2;
     652             :   GEN S;
     653     1479987 :   if (L <= l) return T;
     654     1313799 :   S = cgetg(l, t_POL);
     655     1317222 :   S[1] = T[1];
     656     1317222 :   for (i = 2; i < l; i++) gel(S,i) = gel(T,i);
     657     5986114 :   for (j = 2; i < L; i++) {
     658     4674043 :     gel(S,j) = addii(gel(S,j), gel(T,i));
     659     4668892 :     if (++j == l) j = 2;
     660             :   }
     661     1312071 :   return normalizepol_lg(S, l);
     662             : }
     663             : 
     664             : /*******************************************************************/
     665             : /*                                                                 */
     666             : /*                                ZXV                              */
     667             : /*                                                                 */
     668             : /*******************************************************************/
     669             : 
     670             : int
     671          35 : ZXV_equal(GEN V, GEN W)
     672             : {
     673          35 :   long l = lg(V);
     674          35 :   if (l!=lg(W)) return 0;
     675          70 :   while (--l > 0)
     676          35 :     if (!ZX_equal(gel(V,l), gel(W,l))) return 0;
     677           0 :   return 1;
     678             : }
     679             : 
     680             : GEN
     681        7476 : ZXV_Z_mul(GEN x, GEN y)
     682        7476 : { pari_APPLY_same(ZX_Z_mul(gel(x,i), y)) }
     683             : 
     684             : GEN
     685           0 : ZXV_remi2n(GEN x, long N)
     686           0 : { pari_APPLY_same(ZX_remi2n(gel(x,i), N)) }
     687             : 
     688             : GEN
     689       82873 : ZXV_dotproduct(GEN x, GEN y)
     690             : {
     691       82873 :   pari_sp av = avma;
     692       82873 :   long i, lx = lg(x);
     693             :   GEN c;
     694       82873 :   if (lx == 1) return pol_0(varn(x));
     695       82873 :   c = ZX_mul(gel(x,1), gel(y,1));
     696      391629 :   for (i = 2; i < lx; i++)
     697             :   {
     698      308756 :     GEN t = ZX_mul(gel(x,i), gel(y,i));
     699      308756 :     if (signe(t)) c = ZX_add(c, t);
     700             :   }
     701       82873 :   return gerepileupto(av, c);
     702             : }
     703             : 
     704             : /*******************************************************************/
     705             : /*                                                                 */
     706             : /*                                ZXQM                             */
     707             : /*                                                                 */
     708             : /*******************************************************************/
     709             : 
     710             : GEN
     711     2091047 : ZXn_mul(GEN x, GEN y, long n)
     712     2091047 : { return RgXn_red_shallow(ZX_mul(x, y), n); }
     713             : 
     714             : GEN
     715        2044 : ZXn_sqr(GEN x, long n)
     716        2044 : { return RgXn_red_shallow(ZX_sqr(x), n); }
     717             : 
     718             : /*******************************************************************/
     719             : /*                                                                 */
     720             : /*                                ZXQM                             */
     721             : /*                                                                 */
     722             : /*******************************************************************/
     723             : 
     724             : static long
     725     5123564 : ZX_expi(GEN x)
     726             : {
     727     5123564 :   if (signe(x)==0) return 0;
     728     2352058 :   if (typ(x)==t_INT) return expi(x);
     729     1689739 :   return ZX_expispec(x+2, lgpol(x));
     730             : }
     731             : 
     732             : static long
     733      771701 : ZXC_expi(GEN x)
     734             : {
     735      771701 :   long i, l = lg(x), m=0;
     736     5895265 :   for(i = 1; i < l; i++)
     737             :   {
     738     5123564 :     long e = ZX_expi(gel(x,i));
     739     5123564 :     if (e > m) m = e;
     740             :   }
     741      771701 :   return m;
     742             : }
     743             : 
     744             : static long
     745      142406 : ZXM_expi(GEN x)
     746             : {
     747      142406 :   long i, l = lg(x), m=0;
     748      914107 :   for(i = 1; i < l; i++)
     749             :   {
     750      771701 :     long e = ZXC_expi(gel(x,i));
     751      771701 :     if (e > m) m = e;
     752             :   }
     753      142406 :   return m;
     754             : }
     755             : 
     756             : static GEN
     757     5123564 : ZX_eval2BIL(GEN x, long k)
     758             : {
     759     5123564 :   if (signe(x)==0) return gen_0;
     760     2352058 :   if (typ(x)==t_INT) return x;
     761     1689739 :   return ZX_eval2BILspec(x+2, k, lgpol(x));
     762             : }
     763             : 
     764             : /*Eval x in 2^(k*BIL) in linear time*/
     765             : static GEN
     766      771701 : ZXC_eval2BIL(GEN x, long k)
     767             : {
     768      771701 :   long i, lx = lg(x);
     769      771701 :   GEN A = cgetg(lx, t_COL);
     770      771701 :   for (i=1; i<lx; i++) gel(A,i) = ZX_eval2BIL(gel(x,i), k);
     771      771701 :   return A;
     772             : }
     773             : 
     774             : static GEN
     775      142406 : ZXM_eval2BIL(GEN x, long k)
     776             : {
     777      142406 :   long i, lx = lg(x);
     778      142406 :   GEN A = cgetg(lx, t_MAT);
     779      142406 :   for (i=1; i<lx; i++) gel(A,i) = ZXC_eval2BIL(gel(x,i), k);
     780      142406 :   return A;
     781             : }
     782             : 
     783             : static GEN
     784      661285 : Z_mod2BIL_ZXQ(GEN x, long bs, GEN T)
     785             : {
     786      661285 :   pari_sp av = avma;
     787      661285 :   long v = varn(T), d = 2*(degpol(T)-1);
     788      661285 :   GEN z = Z_mod2BIL_ZX(x, bs, d, 0);
     789      661285 :   setvarn(z, v);
     790      661285 :   return gerepileupto(av, ZX_rem(z, T));
     791             : }
     792             : 
     793             : static GEN
     794       26712 : ZC_mod2BIL_ZXQC(GEN x, long bs, GEN T)
     795             : {
     796       26712 :   long i, lx = lg(x);
     797       26712 :   GEN A = cgetg(lx, t_COL);
     798       26712 :   for (i=1; i<lx; i++) gel(A,i) = Z_mod2BIL_ZXQ(gel(x,i), bs, T);
     799       26712 :   return A;
     800             : }
     801             : 
     802             : static GEN
     803        8476 : ZM_mod2BIL_ZXQM(GEN x, long bs, GEN T)
     804             : {
     805        8476 :   long i, lx = lg(x);
     806        8476 :   GEN A = cgetg(lx, t_MAT);
     807        8476 :   for (i=1; i<lx; i++) gel(A,i) = ZC_mod2BIL_ZXQC(gel(x,i), bs, T);
     808        8476 :   return A;
     809             : }
     810             : 
     811             : GEN
     812        8336 : ZXQM_mul(GEN x, GEN y, GEN T)
     813             : {
     814        8336 :   long d = degpol(T);
     815             :   GEN z;
     816        8336 :   pari_sp av = avma;
     817        8336 :   if (d == 0)
     818           0 :     z = ZM_mul(simplify_shallow(x),simplify_shallow(y));
     819             :   else
     820             :   {
     821        8336 :     long ex = ZXM_expi(x), ey = ZXM_expi(y), d= degpol(T), n = lg(x)-1;
     822        8336 :     long e = ex + ey + expu(d) + expu(n) + 4;
     823        8336 :     long N = divsBIL(e)+1;
     824        8336 :     z = ZM_mul(ZXM_eval2BIL(x,N), ZXM_eval2BIL(y,N));
     825        8336 :     z = ZM_mod2BIL_ZXQM(z, N, T);
     826             :   }
     827        8336 :   return gerepileupto(av, z);
     828             : }
     829             : 
     830             : GEN
     831         140 : ZXQM_sqr(GEN x, GEN T)
     832             : {
     833         140 :   long d = degpol(T);
     834             :   GEN z;
     835         140 :   pari_sp av = avma;
     836         140 :   if (d == 0)
     837           0 :     z = ZM_sqr(simplify_shallow(x));
     838             :   else
     839             :   {
     840         140 :     long ex = ZXM_expi(x), d = degpol(T), n = lg(x)-1;
     841         140 :     long e = 2*ex + expu(d) + expu(n) + 4;
     842         140 :     long N = divsBIL(e)+1;
     843         140 :     z = ZM_sqr(ZXM_eval2BIL(x,N));
     844         140 :     z = ZM_mod2BIL_ZXQM(z, N, T);
     845             :   }
     846         140 :   return gerepileupto(av, z);
     847             : }
     848             : 
     849             : GEN
     850        4571 : QXQM_mul(GEN x, GEN y, GEN T)
     851             : {
     852        4571 :   GEN dx, nx = Q_primitive_part(x, &dx);
     853        4571 :   GEN dy, ny = Q_primitive_part(y, &dy);
     854        4571 :   GEN z = ZXQM_mul(nx, ny, T);
     855        4571 :   if (dx || dy)
     856             :   {
     857        4571 :     GEN d = dx ? dy ? gmul(dx, dy): dx : dy;
     858        4571 :     if (!gequal1(d)) z = RgM_Rg_mul(z, d);
     859             :   }
     860        4571 :   return z;
     861             : }
     862             : 
     863             : GEN
     864           7 : QXQM_sqr(GEN x, GEN T)
     865             : {
     866           7 :   GEN dx, nx = Q_primitive_part(x, &dx);
     867           7 :   GEN z = ZXQM_sqr(nx, T);
     868           7 :   if (dx) z = RgM_Rg_mul(z, gsqr(dx));
     869           7 :   return z;
     870             : }
     871             : 
     872             : static GEN
     873     1004876 : Z_mod2BIL_Fq(GEN x, long bs, GEN T, GEN p)
     874             : {
     875     1004876 :   pari_sp av = avma;
     876     1004876 :   long v = get_FpX_var(T), d = 2*(get_FpX_degree(T)-1);
     877     1004876 :   GEN z = Z_mod2BIL_ZX(x, bs, d, 0);
     878     1004876 :   setvarn(z, v);
     879     1004876 :   return gerepileupto(av, FpX_rem(z, T, p));
     880             : }
     881             : 
     882             : static GEN
     883      240729 : ZC_mod2BIL_FqC(GEN x, long bs, GEN T, GEN p)
     884             : {
     885      240729 :   long i, lx = lg(x);
     886      240729 :   GEN A = cgetg(lx, t_COL);
     887      240729 :   for (i=1; i<lx; i++) gel(A,i) = Z_mod2BIL_Fq(gel(x,i), bs, T, p);
     888      240729 :   return A;
     889             : }
     890             : 
     891             : static GEN
     892       62797 : ZM_mod2BIL_FqM(GEN x, long bs, GEN T, GEN p)
     893             : {
     894       62797 :   long i, lx = lg(x);
     895       62797 :   GEN A = cgetg(lx, t_MAT);
     896       62797 :   for (i=1; i<lx; i++) gel(A,i) = ZC_mod2BIL_FqC(gel(x,i), bs, T, p);
     897       62797 :   return A;
     898             : }
     899             : 
     900             : GEN
     901       62797 : FqM_mul_Kronecker(GEN x, GEN y, GEN T, GEN p)
     902             : {
     903       62797 :   pari_sp av = avma;
     904       62797 :   long ex = ZXM_expi(x), ey = ZXM_expi(y), d= degpol(T), n = lg(x)-1;
     905       62797 :   long e = ex + ey + expu(d) + expu(n) + 4;
     906       62797 :   long N = divsBIL(e)+1;
     907       62797 :   GEN  z = ZM_mul(ZXM_eval2BIL(x,N), ZXM_eval2BIL(y,N));
     908       62797 :   return gerepileupto(av, ZM_mod2BIL_FqM(z, N, T, p));
     909             : }
     910             : 
     911             : /*******************************************************************/
     912             : /*                                                                 */
     913             : /*                                ZXX                              */
     914             : /*                                                                 */
     915             : /*******************************************************************/
     916             : 
     917             : void
     918          42 : RgX_check_ZXX(GEN x, const char *s)
     919             : {
     920          42 :   long k = lg(x)-1;
     921        5110 :   for ( ; k>1; k--) {
     922        5068 :     GEN t = gel(x,k);
     923        5068 :     switch(typ(t)) {
     924        5054 :       case t_INT: break;
     925          14 :       case t_POL: if (RgX_is_ZX(t)) break;
     926             :       /* fall through */
     927           0 :       default: pari_err_TYPE(stack_strcat(s, " not in Z[X,Y]"),x);
     928             :     }
     929             :   }
     930          42 : }
     931             : 
     932             : /*Renormalize (in place) polynomial with t_INT or ZX coefficients.*/
     933             : GEN
     934   128563093 : ZXX_renormalize(GEN x, long lx)
     935             : {
     936             :   long i;
     937   160688054 :   for (i = lx-1; i>1; i--)
     938   151156984 :     if (signe(gel(x,i))) break;
     939   128563093 :   stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + (i+1)));
     940   128563694 :   setlg(x, i+1); setsigne(x, i!=1); return x;
     941             : }
     942             : 
     943             : long
     944        2275 : ZXX_max_lg(GEN x)
     945             : {
     946        2275 :   long i, prec = 0, lx = lg(x);
     947       13328 :   for (i=2; i<lx; i++)
     948             :   {
     949       11053 :     GEN p = gel(x,i);
     950       11053 :     long l = (typ(p) == t_INT)? lgefint(p): ZX_max_lg(p);
     951       11053 :     if (l > prec) prec = l;
     952             :   }
     953        2275 :   return prec;
     954             : }
     955             : 
     956             : GEN
     957        2163 : ZXX_Z_mul(GEN y, GEN x)
     958             : {
     959        2163 :   long i, l = lg(y);
     960        2163 :   GEN z = cgetg(l,t_POL); z[1] = y[1];
     961      151732 :   for(i=2; i<l; i++)
     962      149569 :     if(typ(gel(y,i))==t_INT)
     963           0 :       gel(z,i) = mulii(gel(y,i),x);
     964             :     else
     965      149569 :       gel(z,i) = ZX_Z_mul(gel(y,i),x);
     966        2163 :   return z;
     967             : }
     968             : 
     969             : GEN
     970          56 : ZXX_Z_add_shallow(GEN x, GEN y)
     971             : {
     972          56 :   long i, l = lg(x);
     973             :   GEN z, a;
     974          56 :   if (signe(x)==0) return scalarpol(y,varn(x));
     975          56 :   z = cgetg(l,t_POL); z[1] = x[1];
     976          56 :   a = gel(x,2);
     977          56 :   gel(z, 2) = typ(a)==t_INT? addii(a,y): ZX_Z_add(a,y);
     978         280 :   for(i=3; i<l; i++)
     979         224 :     gel(z,i) = gel(x,i);
     980          56 :   return z;
     981             : }
     982             : 
     983             : GEN
     984       53067 : ZXX_Z_divexact(GEN y, GEN x)
     985             : {
     986       53067 :   long i, l = lg(y);
     987       53067 :   GEN z = cgetg(l,t_POL); z[1] = y[1];
     988      405048 :   for(i=2; i<l; i++)
     989      351981 :     if(typ(gel(y,i))==t_INT)
     990        2856 :       gel(z,i) = diviiexact(gel(y,i),x);
     991             :     else
     992      349125 :       gel(z,i) = ZX_Z_divexact(gel(y,i),x);
     993       53067 :   return z;
     994             : }
     995             : 
     996             : /* Kronecker substitution, ZXX -> ZX:
     997             :  * P(X,Y) = sum_{0<=i<lP} P_i(X) * Y^i, where deg P_i < n.
     998             :  * Returns P(X,X^(2n-1)) */
     999             : GEN
    1000     1820514 : ZXX_to_Kronecker_spec(GEN P, long lP, long n)
    1001             : {
    1002     1820514 :   long i, k, N = (n<<1) + 1;
    1003             :   GEN y;
    1004     1820514 :   if (!lP) return pol_0(0);
    1005     1807228 :   y = cgetg((N-2)*lP + 2, t_POL) + 2;
    1006    16837480 :   for (k=i=0; i<lP; i++)
    1007             :   {
    1008             :     long j;
    1009    16837480 :     GEN c = gel(P,i);
    1010    16837480 :     if (typ(c)!=t_POL)
    1011             :     {
    1012     1152214 :       gel(y,k++) = c;
    1013     1152214 :       j = 3;
    1014             :     }
    1015             :     else
    1016             :     {
    1017    15685266 :       long l = lg(c);
    1018    15685266 :       if (l-3 >= n)
    1019           0 :         pari_err_BUG("ZXX_to_Kronecker, P is not reduced mod Q");
    1020    15685266 :       for (j=2; j < l; j++) gel(y,k++) = gel(c,j);
    1021             :     }
    1022    16837480 :     if (i == lP-1) break;
    1023    15030252 :     for (   ; j < N; j++) gel(y,k++) = gen_0;
    1024             :   }
    1025     1807228 :   y-=2; setlg(y, k+2); y[1] = evalsigne(1); return y;
    1026             : }
    1027             : 
    1028             : GEN
    1029      755614 : ZXX_to_Kronecker(GEN P, long n)
    1030             : {
    1031      755614 :   GEN z = ZXX_to_Kronecker_spec(P+2, lgpol(P), n);
    1032      755614 :   setvarn(z,varn(P)); return z;
    1033             : }
    1034             : 
    1035             : GEN
    1036           0 : ZXQX_sqr(GEN x, GEN T)
    1037             : {
    1038           0 :   pari_sp av = avma;
    1039           0 :   long n = degpol(T);
    1040           0 :   GEN z = ZXX_sqr_Kronecker(x, n);
    1041           0 :   z = Kronecker_to_ZXX(z, n, varn(T));
    1042           0 :   return gerepileupto(av, z);
    1043             : }
    1044             : 
    1045             : GEN
    1046           0 : ZXQX_mul(GEN x, GEN y, GEN T)
    1047             : {
    1048           0 :   pari_sp av = avma;
    1049           0 :   long n = degpol(T);
    1050           0 :   GEN z = ZXX_mul_Kronecker(x, y, n);
    1051           0 :   z = Kronecker_to_ZXX(z, n, varn(T));
    1052           0 :   return gerepileupto(av, z);
    1053             : }
    1054             : 
    1055             : GEN
    1056     1705358 : QX_mul(GEN x, GEN y)
    1057             : {
    1058     1705358 :   GEN dx, nx = Q_primitive_part(x, &dx);
    1059     1705358 :   GEN dy, ny = Q_primitive_part(y, &dy);
    1060     1705358 :   GEN z = ZX_mul(nx, ny);
    1061     1705358 :   if (dx || dy)
    1062             :   {
    1063     1624426 :     GEN d = dx ? dy ? gmul(dx, dy): dx : dy;
    1064     1624426 :     return ZX_Q_mul(z, d);
    1065             :   } else
    1066       80932 :     return z;
    1067             : }
    1068             : 
    1069             : GEN
    1070       45068 : QX_sqr(GEN x)
    1071             : {
    1072       45068 :   GEN dx, nx = Q_primitive_part(x, &dx);
    1073       45068 :   GEN z = ZX_sqr(nx);
    1074       45068 :   if (dx)
    1075       35022 :     return ZX_Q_mul(z, gsqr(dx));
    1076             :   else
    1077       10046 :     return z;
    1078             : }
    1079             : 
    1080             : GEN
    1081     1221884 : QX_ZX_rem(GEN x, GEN y)
    1082             : {
    1083     1221884 :   pari_sp av = avma;
    1084     1221884 :   GEN d, nx = Q_primitive_part(x, &d);
    1085     1221884 :   GEN r = ZX_rem(nx, y);
    1086     1221884 :   if (d) r = ZX_Q_mul(r, d);
    1087     1221884 :   return gerepileupto(av, r);
    1088             : }

Generated by: LCOV version 1.13