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 - F2x.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.12.1 lcov report (development 25406-bf255ab81b) Lines: 1440 1725 83.5 %
Date: 2020-06-04 05:59:24 Functions: 175 202 86.6 %
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             : /* Not so fast arithmetic with polynomials over F_2 */
      18             : 
      19             : /***********************************************************************/
      20             : /**                                                                   **/
      21             : /**                             F2x                                   **/
      22             : /**                                                                   **/
      23             : /***********************************************************************/
      24             : /* F2x objects are defined as follows:
      25             :    An F2x is a t_VECSMALL:
      26             :    x[0] = codeword
      27             :    x[1] = evalvarn(variable number)  (signe is not stored).
      28             :    x[2] = a_0...a_31 x[3] = a_32..a_63, etc.  on 32bit
      29             :    x[2] = a_0...a_63 x[3] = a_64..a_127, etc. on 64bit
      30             : 
      31             :    where the a_i are bits.
      32             : 
      33             :    signe(x) is not valid. Use lgpol(x)!=0 instead.
      34             :    Note: pol0_F2x=pol0_Flx and pol1_F2x=pol1_Flx
      35             : */
      36             : 
      37             : INLINE long
      38         204 : F2x_degreespec(GEN x, long l)
      39             : {
      40         204 :   return (l==0)?-1:l*BITS_IN_LONG-bfffo(x[l-1])-1;
      41             : }
      42             : 
      43             : INLINE long
      44   239831110 : F2x_degree_lg(GEN x, long l)
      45             : {
      46   239831110 :   return (l==2)?-1:bit_accuracy(l)-bfffo(x[l-1])-1;
      47             : }
      48             : 
      49             : long
      50    66503870 : F2x_degree(GEN x)
      51             : {
      52    66503870 :   return F2x_degree_lg(x, lg(x));
      53             : }
      54             : 
      55             : GEN
      56          63 : monomial_F2x(long d, long vs)
      57             : {
      58          63 :   long l=nbits2lg(d+1);
      59          63 :   GEN z = zero_zv(l-1);
      60          63 :   z[1] = vs;
      61          63 :   F2x_set(z,d);
      62          63 :   return z;
      63             : }
      64             : 
      65             : GEN
      66      226148 : F2x_to_ZX(GEN x)
      67             : {
      68      226148 :   long l = 3+F2x_degree(x), lx = lg(x);
      69      226148 :   GEN z = cgetg(l,t_POL);
      70             :   long i, j ,k;
      71      453248 :   for (i=2, k=2; i<lx; i++)
      72     1007850 :     for (j=0; j<BITS_IN_LONG && k<l; j++,k++)
      73      780750 :       gel(z,k) = (x[i]&(1UL<<j))?gen_1:gen_0;
      74      226148 :   z[1]=evalsigne(l>=3)|x[1];
      75      226148 :   return z;
      76             : }
      77             : 
      78             : GEN
      79      105622 : F2x_to_Flx(GEN x)
      80             : {
      81      105622 :   long l = 3+F2x_degree(x), lx = lg(x);
      82      105621 :   GEN z = cgetg(l,t_VECSMALL);
      83             :   long i, j, k;
      84      105635 :   z[1] = x[1];
      85      234816 :   for (i=2, k=2; i<lx; i++)
      86     3264103 :     for (j=0; j<BITS_IN_LONG && k<l; j++,k++)
      87     3134922 :       z[k] = (x[i]>>j)&1UL;
      88      105635 :   return z;
      89             : }
      90             : 
      91             : GEN
      92         378 : F2x_to_F2xX(GEN z, long sv)
      93             : {
      94         378 :   long i, d = F2x_degree(z);
      95         378 :   GEN x = cgetg(d+3,t_POL);
      96        6769 :   for (i=0; i<=d; i++)
      97        6391 :     gel(x,i+2) = F2x_coeff(z,i) ? pol1_F2x(sv): pol0_F2x(sv);
      98         378 :   x[1] = evalsigne(d+1!=0)| z[1]; return x;
      99             : }
     100             : 
     101             : GEN
     102      209174 : Z_to_F2x(GEN x, long sv)
     103             : {
     104      209174 :   return mpodd(x) ? pol1_F2x(sv): pol0_F2x(sv);
     105             : }
     106             : 
     107             : GEN
     108      256612 : ZX_to_F2x(GEN x)
     109             : {
     110      256612 :   long lx = lg(x), l = nbits2lg(lx-2);
     111      256612 :   GEN z = cgetg(l,t_VECSMALL);
     112             :   long i, j, k;
     113      256612 :   z[1] = ((ulong)x[1])&VARNBITS;
     114     1186173 :   for (i=2, k=1,j=BITS_IN_LONG; i<lx; i++,j++)
     115             :   {
     116      929561 :     if (j==BITS_IN_LONG)
     117             :     {
     118      240626 :       j=0; z[++k]=0;
     119             :     }
     120      929561 :     if (mpodd(gel(x,i)))
     121      538434 :       z[k] |= 1UL<<j;
     122             :   }
     123      256612 :   return F2x_renormalize(z,l);
     124             : }
     125             : 
     126             : GEN
     127       98823 : Flx_to_F2x(GEN x)
     128             : {
     129       98823 :   long lx = lg(x), l = nbits2lg(lx-2);
     130       98820 :   GEN z = cgetg(l,t_VECSMALL);
     131             :   long i, j, k;
     132       98828 :   z[1] = x[1];
     133     3203598 :   for (i=2, k=1, j=BITS_IN_LONG; i<lx; i++,j++)
     134             :   {
     135     3104770 :     if (j==BITS_IN_LONG)
     136             :     {
     137      115639 :       j=0; z[++k] = 0;
     138             :     }
     139     3104770 :     if (x[i]&1UL)
     140     1589931 :       z[k] |= 1UL<<j;
     141             :   }
     142       98828 :   return F2x_renormalize(z,l);
     143             : }
     144             : 
     145             : GEN
     146     1507193 : F2x_to_F2v(GEN x, long N)
     147             : {
     148     1507193 :   long i, l = lg(x);
     149     1507193 :   long n = nbits2lg(N);
     150     1507189 :   GEN z = cgetg(n,t_VECSMALL);
     151     1507188 :   z[1] = N;
     152     3245256 :   for (i=2; i<l ; i++) z[i]=x[i];
     153     1597553 :   for (   ; i<n; i++) z[i]=0;
     154     1507188 :   return z;
     155             : }
     156             : 
     157             : GEN
     158        1190 : RgX_to_F2x(GEN x)
     159             : {
     160        1190 :   long l=nbits2lg(lgpol(x));
     161        1190 :   GEN z=cgetg(l,t_VECSMALL);
     162             :   long i,j,k;
     163        1190 :   z[1]=((ulong)x[1])&VARNBITS;
     164        3724 :   for(i=2, k=1,j=BITS_IN_LONG;i<lg(x);i++,j++)
     165             :   {
     166        2534 :     if (j==BITS_IN_LONG)
     167             :     {
     168        1176 :       j=0; k++; z[k]=0;
     169             :     }
     170        2534 :     if (Rg_to_F2(gel(x,i)))
     171        1253 :       z[k]|=1UL<<j;
     172             :   }
     173        1190 :   return F2x_renormalize(z,l);
     174             : }
     175             : 
     176             : /* If x is a POLMOD, assume modulus is a multiple of T. */
     177             : GEN
     178     1162387 : Rg_to_F2xq(GEN x, GEN T)
     179             : {
     180     1162387 :   long ta, tx = typ(x), v = get_F2x_var(T);
     181             :   GEN a, b;
     182     1162386 :   if (is_const_t(tx))
     183             :   {
     184     1161358 :     if (tx == t_FFELT) return FF_to_F2xq(x);
     185      244380 :     return Rg_to_F2(x)? pol1_F2x(v): pol0_F2x(v);
     186             :   }
     187        1029 :   switch(tx)
     188             :   {
     189           0 :     case t_POLMOD:
     190           0 :       b = gel(x,1);
     191           0 :       a = gel(x,2); ta = typ(a);
     192           0 :       if (is_const_t(ta)) return Rg_to_F2(a)? pol1_F2x(v): pol0_F2x(v);
     193           0 :       b = RgX_to_F2x(b); if (b[1] != v) break;
     194           0 :       a = RgX_to_F2x(a); if (F2x_equal(b,T)) return a;
     195           0 :       if (lgpol(F2x_rem(b,T))==0) return F2x_rem(a, T);
     196           0 :       break;
     197        1029 :     case t_POL:
     198        1029 :       x = RgX_to_F2x(x);
     199        1029 :       if (x[1] != v) break;
     200        1029 :       return F2x_rem(x, T);
     201           0 :     case t_RFRAC:
     202           0 :       a = Rg_to_F2xq(gel(x,1), T);
     203           0 :       b = Rg_to_F2xq(gel(x,2), T);
     204           0 :       return F2xq_div(a,b, T);
     205             :   }
     206           0 :   pari_err_TYPE("Rg_to_F2xq",x);
     207             :   return NULL; /* LCOV_EXCL_LINE */
     208             : }
     209             : 
     210             : ulong
     211         245 : F2x_eval(GEN P, ulong x)
     212             : {
     213         245 :   if (odd(x))
     214             :   {
     215         133 :     long i, lP = lg(P);
     216         133 :     ulong c = 0;
     217         266 :     for (i=2; i<lP; i++)
     218         133 :       c ^= P[i];
     219             : #ifdef LONG_IS_64BIT
     220         114 :     c ^= c >> 32;
     221             : #endif
     222         133 :     c ^= c >> 16;
     223         133 :     c ^= c >>  8;
     224         133 :     c ^= c >>  4;
     225         133 :     c ^= c >>  2;
     226         133 :     c ^= c >>  1;
     227         133 :     return c & 1;
     228             :   }
     229         112 :   else return F2x_coeff(P,0);
     230             : }
     231             : 
     232             : GEN
     233    30508868 : F2x_add(GEN x, GEN y)
     234             : {
     235             :   long i,lz;
     236             :   GEN z;
     237    30508868 :   long lx=lg(x);
     238    30508868 :   long ly=lg(y);
     239    30508868 :   if (ly>lx) swapspec(x,y, lx,ly);
     240    30508868 :   lz = lx; z = cgetg(lz, t_VECSMALL); z[1]=x[1];
     241    65324709 :   for (i=2; i<ly; i++) z[i] = x[i]^y[i];
     242    36870701 :   for (   ; i<lx; i++) z[i] = x[i];
     243    30529683 :   return F2x_renormalize(z, lz);
     244             : }
     245             : 
     246             : static GEN
     247       57750 : F2x_addspec(GEN x, GEN y, long lx, long ly)
     248             : {
     249             :   long i,lz;
     250             :   GEN z;
     251             : 
     252       57750 :   if (ly>lx) swapspec(x,y, lx,ly);
     253       57750 :   lz = lx+2; z = cgetg(lz, t_VECSMALL) + 2;
     254      167426 :   for (i=0; i<ly-3; i+=4)
     255             :   {
     256      109676 :     z[i] = x[i]^y[i];
     257      109676 :     z[i+1] = x[i+1]^y[i+1];
     258      109676 :     z[i+2] = x[i+2]^y[i+2];
     259      109676 :     z[i+3] = x[i+3]^y[i+3];
     260             :   }
     261       94682 :   for (; i<ly; i++)
     262       36932 :     z[i] = x[i]^y[i];
     263      525234 :   for (   ; i<lx; i++) z[i] = x[i];
     264       57750 :   z -= 2; z[1] = 0; return F2x_renormalize(z, lz);
     265             : }
     266             : 
     267             : GEN
     268      190258 : F2x_1_add(GEN y)
     269             : {
     270             :   GEN z;
     271             :   long lz, i;
     272      190258 :   if (!lgpol(y))
     273       90165 :     return pol1_F2x(y[1]);
     274      100093 :   lz=lg(y);
     275      100093 :   z=cgetg(lz,t_VECSMALL);
     276      100093 :   z[1] = y[1];
     277      100093 :   z[2] = y[2]^1UL;
     278      100093 :   for(i=3;i<lz;i++)
     279           0 :     z[i] = y[i];
     280      100093 :   if (lz==3) z = F2x_renormalize(z,lz);
     281      100093 :   return z;
     282             : }
     283             : 
     284             : INLINE void
     285   158083037 : F2x_addshiftipspec(GEN x, GEN y, long ny, ulong db)
     286             : {
     287             :   long i;
     288   158083037 :   if (db)
     289             :   {
     290   136045398 :     ulong dc=BITS_IN_LONG-db;
     291   136045398 :     ulong r=0, yi;
     292   216427991 :     for(i=0; i<ny-3; i+=4)
     293             :     {
     294    80382593 :       yi = uel(y,i);   x[i]   ^= (yi<<db)|r; r = yi>>dc;
     295    80382593 :       yi = uel(y,i+1); x[i+1] ^= (yi<<db)|r; r = yi>>dc;
     296    80382593 :       yi = uel(y,i+2); x[i+2] ^= (yi<<db)|r; r = yi>>dc;
     297    80382593 :       yi = uel(y,i+3); x[i+3] ^= (yi<<db)|r; r = yi>>dc;
     298             :     }
     299   261558816 :     for(  ; i<ny; i++)
     300             :     {
     301   125513418 :       ulong yi = uel(y,i); x[i] ^= (yi<<db)|r; r = yi>>dc;
     302             :     }
     303   136045398 :     if (r) x[i] ^= r;
     304             :   }
     305             :   else
     306             :   {
     307    25390746 :     for(i=0; i<ny-3; i+=4)
     308             :     {
     309     3353107 :       x[i]   ^= y[i];
     310     3353107 :       x[i+1] ^= y[i+1];
     311     3353107 :       x[i+2] ^= y[i+2];
     312     3353107 :       x[i+3] ^= y[i+3];
     313             :     }
     314    47763226 :     for(   ; i<ny; i++)
     315    25725587 :       x[i] ^= y[i];
     316             :   }
     317   158083037 : }
     318             : 
     319             : INLINE void
     320   127751832 : F2x_addshiftip(GEN x, GEN y, ulong d)
     321             : {
     322   127751832 :   ulong db, dl=dvmduBIL(d, &db);
     323   127208967 :   F2x_addshiftipspec(x+2+dl, y+2, lgpol(y), db);
     324   127439663 : }
     325             : 
     326             : static GEN
     327    18022626 : F2x_mul1(ulong x, ulong y)
     328             : {
     329    18022626 :   ulong x1=(x&HIGHMASK)>>BITS_IN_HALFULONG;
     330    18022626 :   ulong x2=x&LOWMASK;
     331    18022626 :   ulong y1=(y&HIGHMASK)>>BITS_IN_HALFULONG;
     332    18022626 :   ulong y2=y&LOWMASK;
     333             :   ulong r1,r2,rr;
     334             :   GEN z;
     335             :   ulong i;
     336    18022626 :   rr=r1=r2=0UL;
     337    18022626 :   if (x2)
     338   556090295 :     for(i=0;i<BITS_IN_HALFULONG;i++)
     339   538104149 :       if (x2&(1UL<<i))
     340             :       {
     341    48769698 :         r2^=y2<<i;
     342    48769698 :         rr^=y1<<i;
     343             :       }
     344    18022626 :   if (x1)
     345     3221972 :     for(i=0;i<BITS_IN_HALFULONG;i++)
     346     3108888 :       if (x1&(1UL<<i))
     347             :       {
     348      889587 :         r1^=y1<<i;
     349      889587 :         rr^=y2<<i;
     350             :       }
     351    18022626 :   r2^=(rr&LOWMASK)<<BITS_IN_HALFULONG;
     352    18022626 :   r1^=(rr&HIGHMASK)>>BITS_IN_HALFULONG;
     353    18022626 :   z=cgetg((r1?4:3),t_VECSMALL);
     354    18022794 :   z[2]=r2;
     355    18022794 :   if (r1) z[3]=r1;
     356    18022794 :   return z;
     357             : }
     358             : 
     359             : static GEN
     360     5074990 : F2x_mulspec_basecase(GEN x, GEN y, long nx, long ny)
     361             : {
     362             :   long l, i, j;
     363             :   GEN z;
     364     5074990 :   l = nx + ny;
     365     5074990 :   z = zero_Flv(l+1);
     366     5725117 :   for(i=0; i < ny-1; i++)
     367             :   {
     368      650487 :     GEN zi = z+2+i;
     369      650487 :     ulong yi = uel(y,i);
     370      650487 :     if (yi)
     371    35644962 :       for(j=0; j < BITS_IN_LONG; j++)
     372    34995474 :         if (yi&(1UL<<j)) F2x_addshiftipspec(zi,x,nx,j);
     373             :   }
     374             :   {
     375     5074630 :     GEN zi = z+2+i;
     376     5074630 :     ulong yi = uel(y,i);
     377     5074630 :     long c = BITS_IN_LONG-bfffo(yi);
     378    31320413 :     for(j=0; j < c; j++)
     379    26247242 :       if (yi&(1UL<<j)) F2x_addshiftipspec(zi,x,nx,j);
     380             :   }
     381     5073171 :   return F2x_renormalize(z, l+2);
     382             : }
     383             : 
     384             : static GEN
     385       38503 : F2x_addshift(GEN x, GEN y, long d)
     386             : {
     387       38503 :   GEN xd,yd,zd = (GEN)avma;
     388       38503 :   long a,lz,ny = lgpol(y), nx = lgpol(x);
     389       38503 :   long vs = x[1];
     390       38503 :   if (nx == 0) return y;
     391       38503 :   x += 2; y += 2; a = ny-d;
     392       38503 :   if (a <= 0)
     393             :   {
     394        3127 :     lz = (a>nx)? ny+2: nx+d+2;
     395        3127 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
     396       36811 :     while (xd > x) *--zd = *--xd;
     397        3127 :     x = zd + a;
     398        4332 :     while (zd > x) *--zd = 0;
     399             :   }
     400             :   else
     401             :   {
     402       35376 :     xd = new_chunk(d); yd = y+d;
     403       35376 :     x = F2x_addspec(x,yd,nx,a);
     404       35376 :     lz = (a>nx)? ny+2: lg(x)+d;
     405      742269 :     x += 2; while (xd > x) *--zd = *--xd;
     406             :   }
     407      448585 :   while (yd > y) *--zd = *--yd;
     408       38503 :   *--zd = vs;
     409       38503 :   *--zd = evaltyp(t_VECSMALL) | evallg(lz); return zd;
     410             : }
     411             : 
     412             : /* shift polynomial + gerepile */
     413             : /* Do not set evalvarn. Cf Flx_shiftip */
     414             : static GEN
     415    23124947 : F2x_shiftip(pari_sp av, GEN x, long v)
     416             : {
     417    23124947 :   long i, lx = lg(x), ly;
     418             :   GEN y;
     419    23124947 :   if (!v || lx==2) return gerepileuptoleaf(av, x);
     420       27705 :   ly = lx + v;
     421       27705 :   (void)new_chunk(ly); /* check that result fits */
     422       27704 :   x += lx; y = (GEN)av;
     423      121338 :   for (i = 2; i<lx; i++) *--y = *--x;
     424      157787 :   for (i = 0; i< v; i++) *--y = 0;
     425       27704 :   y -= 2; y[0] = evaltyp(t_VECSMALL) | evallg(ly);
     426       27705 :   set_avma((pari_sp)y); return y;
     427             : }
     428             : 
     429             : static GEN
     430         204 : F2x_to_int(GEN a, long na, long da, long bs)
     431             : {
     432         204 :   const long BIL = BITS_IN_LONG;
     433             :   GEN z, zs;
     434             :   long i,j,k,m;
     435         204 :   long lz = nbits2lg(1+da*bs);
     436         204 :   z = cgeti(lz); z[1] = evalsigne(1)|evallgefint(lz);
     437         204 :   zs = int_LSW(z); *zs = 0;
     438       13668 :   for (i=0, k=2, m=0; i<na; i++)
     439      865989 :     for (j=0; j<BIL; j++, m+=bs)
     440             :     {
     441      852729 :       if (m >= BIL)
     442             :       {
     443      173157 :         k++; if (k>=lz) break;
     444      172953 :         zs = int_nextW(zs);
     445      172953 :         *zs = 0;
     446      172953 :         m -= BIL;
     447             :       }
     448      852525 :       *zs |= ((a[i]>>j)&1UL)<<m;
     449             :     }
     450         204 :   return int_normalize(z,0);
     451             : }
     452             : 
     453             : static GEN
     454         102 : int_to_F2x(GEN x, long d, long bs)
     455             : {
     456         102 :   const long BIL = BITS_IN_LONG;
     457         102 :   long lx = lgefint(x), lz = nbits2lg(d+1);
     458             :   long i,j,k,m;
     459         102 :   GEN xs = int_LSW(x);
     460         102 :   GEN z = cgetg(lz, t_VECSMALL);
     461         102 :   z[1] = 0;
     462       13464 :   for (k=2, i=2, m=0; k < lx; i++)
     463             :   {
     464       13362 :     z[i] = 0;
     465      865647 :     for (j=0; j<BIL; j++, m+=bs)
     466             :     {
     467      852387 :       if (m >= BIL)
     468             :       {
     469      173109 :         if (++k==lx) break;
     470      173007 :         xs = int_nextW(xs);
     471      173007 :         m -= BIL;
     472             :       }
     473      852285 :       if ((*xs>>m)&1UL)
     474      385725 :         z[i]|=1UL<<j;
     475             :     }
     476             :   }
     477         102 :   return F2x_renormalize(z,lz);
     478             : }
     479             : 
     480             : static GEN
     481         102 : F2x_mulspec_mulii(GEN a, GEN b, long na, long nb)
     482             : {
     483         102 :   long da = F2x_degreespec(a,na), db = F2x_degreespec(b,nb);
     484         102 :   long bs = expu(1 + maxss(da, db))+1;
     485         102 :   GEN A = F2x_to_int(a,na,da,bs);
     486         102 :   GEN B = F2x_to_int(b,nb,db,bs);
     487         102 :   GEN z = mulii(A,B);
     488         102 :   return int_to_F2x(z,da+db,bs);
     489             : }
     490             : 
     491             : /* fast product (Karatsuba) of polynomials a,b. These are not real GENs, a+2,
     492             :  * b+2 were sent instead. na, nb = number of terms of a, b.
     493             :  * Only c, c0, c1, c2 are genuine GEN.
     494             :  */
     495             : static GEN
     496    27553249 : F2x_mulspec(GEN a, GEN b, long na, long nb)
     497             : {
     498             :   GEN a0,c,c0;
     499    27553249 :   long n0, n0a, i, v = 0;
     500             :   pari_sp av;
     501    27658590 :   while (na && !a[0]) { a++; na-=1; v++; }
     502    27578191 :   while (nb && !b[0]) { b++; nb-=1; v++; }
     503    27553249 :   if (na < nb) swapspec(a,b, na,nb);
     504    27553249 :   if (!nb) return pol0_F2x(0);
     505             : 
     506    23124551 :   av = avma;
     507    23124551 :   if (na == 1)
     508    18022650 :     return F2x_shiftip(av, F2x_mul1(*a,*b), v);
     509     5101901 :   if (na < F2x_MUL_KARATSUBA_LIMIT)
     510     5074967 :     return F2x_shiftip(av, F2x_mulspec_basecase(a, b, na, nb), v);
     511       26934 :   if (nb >= F2x_MUL_MULII_LIMIT)
     512         102 :     return F2x_shiftip(av, F2x_mulspec_mulii(a, b, na, nb), v);
     513       26832 :   i=(na>>1); n0=na-i; na=i;
     514       26832 :   a0=a+n0; n0a=n0;
     515       26908 :   while (n0a && !a[n0a-1]) n0a--;
     516             : 
     517       26832 :   if (nb > n0)
     518             :   {
     519             :     GEN b0,c1,c2;
     520             :     long n0b;
     521             : 
     522       11187 :     nb -= n0; b0 = b+n0; n0b = n0;
     523       11195 :     while (n0b && !b[n0b-1]) n0b--;
     524       11187 :     c =  F2x_mulspec(a,b,n0a,n0b);
     525       11187 :     c0 = F2x_mulspec(a0,b0,na,nb);
     526             : 
     527       11187 :     c2 = F2x_addspec(a0,a,na,n0a);
     528       11187 :     c1 = F2x_addspec(b0,b,nb,n0b);
     529             : 
     530       11187 :     c1 = F2x_mul(c1,c2);
     531       11187 :     c2 = F2x_add(c0,c);
     532             : 
     533       11187 :     c2 = F2x_add(c1,c2);
     534       11187 :     c0 = F2x_addshift(c0,c2,n0);
     535             :   }
     536             :   else
     537             :   {
     538       15645 :     c  = F2x_mulspec(a,b,n0a,nb);
     539       16129 :     c0 = F2x_mulspec(a0,b,na,nb);
     540             :   }
     541       27316 :   c0 = F2x_addshift(c0,c,n0);
     542       27316 :   return F2x_shiftip(av,c0, v);
     543             : }
     544             : 
     545             : GEN
     546    27498574 : F2x_mul(GEN x, GEN y)
     547             : {
     548    27498574 :   GEN z = F2x_mulspec(x+2,y+2, lgpol(x),lgpol(y));
     549    27499211 :   z[1] = x[1]; return z;
     550             : }
     551             : 
     552             : GEN
     553     7917273 : F2x_sqr(GEN x)
     554             : {
     555     7917273 :   const ulong sq[]={0,1,4,5,16,17,20,21,64,65,68,69,80,81,84,85};
     556             :   long i,ii,j,jj;
     557     7917273 :   long lx=lg(x), lz=2+((lx-2)<<1);
     558             :   GEN z;
     559     7917273 :   z = cgetg(lz, t_VECSMALL); z[1]=x[1];
     560    16019347 :   for (j=2,jj=2;j<lx;j++,jj++)
     561             :   {
     562     8062642 :     ulong x1=((ulong)x[j]&HIGHMASK)>>BITS_IN_HALFULONG;
     563     8062642 :     ulong x2=(ulong)x[j]&LOWMASK;
     564     8062642 :     z[jj]=0;
     565     8062642 :     if (x2)
     566    67397014 :       for(i=0,ii=0;i<BITS_IN_HALFULONG;i+=4,ii+=8)
     567    59386593 :         z[jj]|=sq[(x2>>i)&15UL]<<ii;
     568     8062642 :     z[++jj]=0;
     569     8062642 :     if (x1)
     570     5802462 :       for(i=0,ii=0;i<BITS_IN_HALFULONG;i+=4,ii+=8)
     571     5043618 :         z[jj]|=sq[(x1>>i)&15UL]<<ii;
     572             :   }
     573     7956705 :   return F2x_renormalize(z, lz);
     574             : }
     575             : 
     576             : static GEN
     577      138108 : F2x_pow2n(GEN x, long n)
     578             : {
     579             :   long i;
     580      414326 :   for(i=1;i<=n;i++)
     581      276131 :     x = F2x_sqr(x);
     582      138195 :   return x;
     583             : }
     584             : 
     585             : int
     586       37087 : F2x_issquare(GEN x)
     587             : {
     588       37087 :   const ulong mask = (ULONG_MAX/3UL)*2;
     589       37087 :   ulong i, lx = lg(x);
     590       74146 :   for (i=2; i<lx; i++)
     591       37087 :     if ((x[i]&mask)) return 0;
     592       37059 :   return 1;
     593             : }
     594             : 
     595             : /* Assume x is a perfect square */
     596             : GEN
     597      102393 : F2x_sqrt(GEN x)
     598             : {
     599      102393 :   const ulong sq[]={0,1,4,5,2,3,6,7,8,9,12,13,10,11,14,15};
     600             :   long i,ii,j,jj;
     601      102393 :   long lx=lg(x), lz=2+((lx-1)>>1);
     602             :   GEN z;
     603      102393 :   z = cgetg(lz, t_VECSMALL); z[1]=x[1];
     604      204854 :   for (j=2,jj=2;jj<lz;j++,jj++)
     605             :   {
     606      102441 :     ulong x2=x[j++];
     607      102441 :     z[jj]=0;
     608      102441 :     if (x2)
     609      864945 :       for(i=0,ii=0;ii<BITS_IN_HALFULONG;i+=8,ii+=4)
     610             :       {
     611      762509 :         ulong rl = (x2>>i)&15UL, rh = (x2>>(i+4))&15UL;
     612      762509 :         z[jj]|=sq[rl|(rh<<1)]<<ii;
     613             :       }
     614      102441 :     if (j<lx)
     615             :     {
     616         306 :       x2 = x[j];
     617         306 :       if (x2)
     618        2007 :         for(i=0,ii=0;ii<BITS_IN_HALFULONG;i+=8,ii+=4)
     619             :         {
     620        1716 :           ulong rl = (x2>>i)&15UL, rh = (x2>>(i+4))&15UL;
     621        1716 :           z[jj]|=(sq[rl|(rh<<1)]<<ii)<<BITS_IN_HALFULONG;
     622             :         }
     623             :     }
     624             :   }
     625      102413 :   return F2x_renormalize(z, lz);
     626             : }
     627             : 
     628             : static GEN
     629         673 : F2x_shiftneg(GEN y, ulong d)
     630             : {
     631         673 :   long db, dl=dvmdsBIL(d, &db);
     632         673 :   long i, ly = lg(y), lx = ly-dl;
     633         673 :   GEN x = cgetg(lx, t_VECSMALL);
     634         673 :   x[1] = y[1];
     635         673 :   if (db)
     636             :   {
     637         673 :     ulong dc=BITS_IN_LONG-db;
     638         673 :     ulong r=0;
     639        1346 :     for(i=lx-1; i>=2; i--)
     640             :     {
     641         673 :       x[i] = (((ulong)y[i+dl])>>db)|r;
     642         673 :       r = ((ulong)y[i+dl])<<dc;
     643             :     }
     644             :   }
     645             :   else
     646           0 :     for(i=2; i<lx; i++)
     647           0 :       x[i] = y[i+dl];
     648         673 :   return F2x_renormalize(x,lx);
     649             : }
     650             : 
     651             : static GEN
     652      242070 : F2x_shiftpos(GEN y, ulong d)
     653             : {
     654      242070 :   long db, dl=dvmdsBIL(d, &db);
     655      241917 :   long i, ly = lg(y), lx = ly+dl+(!!db);
     656      241917 :   GEN x = cgetg(lx, t_VECSMALL);
     657      242072 :   x[1] = y[1]; for(i=0; i<dl; i++) x[i+2] = 0;
     658      242032 :   if (db)
     659             :   {
     660      242004 :     ulong dc=BITS_IN_LONG-db;
     661      242004 :     ulong r=0;
     662      486386 :     for(i=2; i<ly; i++)
     663             :     {
     664      244382 :       x[i+dl] = (((ulong)y[i])<<db)|r;
     665      244382 :       r = ((ulong)y[i])>>dc;
     666             :     }
     667      242004 :     x[i+dl] = r;
     668             :   }
     669             :   else
     670          42 :     for(i=2; i<ly; i++)
     671          14 :       x[i+dl] = y[i];
     672      242032 :   return F2x_renormalize(x,lx);
     673             : }
     674             : 
     675             : GEN
     676      242693 : F2x_shift(GEN y, long d)
     677             : {
     678      242693 :   return d<0 ? F2x_shiftneg(y,-d): F2x_shiftpos(y,d);
     679             : }
     680             : 
     681             : #define F2x_recip2(pk,m) u = ((u&m)<<pk)|((u&(~m))>>pk);
     682             : #define F2x_recipu(pk) F2x_recip2(pk,((~0UL)/((1UL<<pk)+1)))
     683             : 
     684             : static ulong
     685           0 : F2x_recip1(ulong u)
     686             : {
     687           0 :   u = (u<<BITS_IN_HALFULONG)|(u>>BITS_IN_HALFULONG);
     688             : #ifdef LONG_IS_64BIT
     689           0 :   F2x_recipu(16);
     690             : #endif
     691           0 :   F2x_recipu(8);
     692           0 :   F2x_recipu(4);
     693           0 :   F2x_recipu(2);
     694           0 :   F2x_recipu(1);
     695           0 :   return u;
     696             : }
     697             : 
     698             : static GEN
     699           0 : F2x_recip_raw(GEN x)
     700             : {
     701             :   long i, l;
     702           0 :   GEN y = cgetg_copy(x,&l);
     703           0 :   y[1] = x[1];
     704           0 :   for (i=0; i<l-2; i++)
     705           0 :     uel(y,l-1-i) = F2x_recip1(uel(x,2+i));
     706           0 :   return y;
     707             : }
     708             : 
     709             : GEN
     710           0 : F2x_recip(GEN x)
     711             : {
     712           0 :   long lb = remsBIL(F2x_degree(x)+1);
     713           0 :   GEN y = F2x_recip_raw(x);
     714           0 :   if (lb) y = F2x_shiftneg(y,BITS_IN_LONG-lb);
     715           0 :   return y;
     716             : }
     717             : 
     718             : GEN
     719           0 : F2xn_red(GEN f, long n)
     720             : {
     721             :   GEN g;
     722             :   long i, dl, db, l;
     723           0 :   if (F2x_degree(f) < n) return zv_copy(f);
     724           0 :   dl = dvmdsBIL(n, &db); l = 2 + dl + (db>0);
     725           0 :   g = cgetg(l, t_VECSMALL);
     726           0 :   g[1] = f[1];
     727           0 :   for (i=2; i<l; i++)
     728           0 :     uel(g,i) = uel(f,i);
     729           0 :   if (db) uel(g,l-1) = uel(f,l-1)&((1UL<<db)-1);
     730           0 :   return F2x_renormalize(g, l);
     731             : }
     732             : 
     733             : static GEN
     734           0 : F2xn_mul(GEN a, GEN b, long n) { return F2xn_red(F2x_mul(a, b), n); }
     735             : 
     736             : static ulong
     737           0 : F2xn_inv_basecase1(ulong x)
     738             : {
     739             :   ulong u, v, w;
     740             :   long i;
     741           0 :   u = x>>1;
     742           0 :   v = (u&1UL)|2UL;
     743           0 :   w = u&v; w ^= w >> 1; v = (w&1UL)|(v<<1);
     744           0 :   w = u&v; w ^= w >> 2; w ^= w >> 1; v = (w&1UL)|(v<<1);
     745           0 :   w = u&v; w ^= w >> 2; w ^= w >> 1; v = (w&1UL)|(v<<1);
     746           0 :   for (i=1;i<=4;i++)
     747           0 :    { w = u&v; w ^= w >> 4; w ^= w >> 2; w ^= w >> 1; v = (w&1UL)|(v<<1); }
     748           0 :   for (i=1;i<=8;i++)
     749           0 :    { w = u&v; w ^= w >> 8; w ^= w >> 4; w ^= w >> 2; w ^= w >> 1; v = (w&1UL)|(v<<1); }
     750           0 :   for (i=1;i<=16;i++)
     751           0 :    { w = u&v; w ^= w >> 16; w ^= w >> 8; w ^= w >> 4; w ^= w >> 2; w ^= w >> 1; v = (w&1UL)|(v<<1); }
     752             : #ifdef LONG_IS_64BIT
     753           0 :   for (i=1; i<=32; i++)
     754           0 :    { w = u&v; w ^= w >> 32; w ^= w >> 16; w ^= w >> 8; w ^= w >> 4; w ^= w >> 2; w ^= w >> 1;
     755           0 :    v = (w&1UL)|(v<<1); }
     756             : #endif
     757           0 :   return (F2x_recip1(v)<<1)|1UL;
     758             : }
     759             : 
     760             : static GEN
     761           0 : F2xn_inv1(GEN v, long e)
     762             : {
     763           0 :   ulong mask = e==BITS_IN_LONG ? -1UL: ((1UL<<e)-1);
     764           0 :   return mkvecsmall2(v[1],F2xn_inv_basecase1(uel(v,2))&mask);
     765             : }
     766             : 
     767             : GEN
     768           0 : F2xn_inv(GEN f, long e)
     769             : {
     770           0 :   pari_sp av = avma, av2;
     771             :   ulong mask;
     772             :   GEN W;
     773           0 :   long n=1;
     774           0 :   if (lg(f)==2) pari_err_INV("Flxn_inv",f);
     775           0 :   if (e <= BITS_IN_LONG) return F2xn_inv1(f, e);
     776           0 :   W = F2xn_inv1(f, BITS_IN_LONG);
     777           0 :   mask = quadratic_prec_mask(divsBIL(e+BITS_IN_LONG-1));
     778           0 :   n = BITS_IN_LONG;
     779           0 :   av2 = avma;
     780           0 :   for (;mask>1;)
     781             :   {
     782             :     GEN u, fr;
     783           0 :     long n2 = n;
     784           0 :     n<<=1; if (mask & 1) n--;
     785           0 :     mask >>= 1;
     786           0 :     fr = F2xn_red(f, n);
     787           0 :     u = F2x_shift(F2xn_mul(W, fr, n), -n2);
     788           0 :     W = F2x_add(W, F2x_shift(F2xn_mul(u, W, n-n2), n2));
     789           0 :     if (gc_needed(av2,2))
     790             :     {
     791           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"F2xn_inv, e = %ld", n);
     792           0 :       W = gerepileupto(av2, W);
     793             :     }
     794             :   }
     795           0 :   return gerepileupto(av, F2xn_red(W,e));
     796             : }
     797             : 
     798             : GEN
     799      218243 : F2x_get_red(GEN T)
     800             : {
     801      218243 :   return T;
     802             : }
     803             : 
     804             : /* separate from F2x_divrem for maximal speed. */
     805             : GEN
     806    43487973 : F2x_rem(GEN x, GEN y)
     807             : {
     808             :   long dx,dy;
     809    43487973 :   long lx=lg(x);
     810    43487973 :   dy = F2x_degree(y); if (!dy) return pol0_F2x(x[1]);
     811    40929582 :   dx = F2x_degree_lg(x,lx);
     812    40951867 :   x  = F2x_copy(x);
     813   145672011 :   while (dx>=dy)
     814             :   {
     815   104593982 :     F2x_addshiftip(x,y,dx-dy);
     816   106299140 :     while (lx>2 && x[lx-1]==0) lx--;
     817   103988287 :     dx = F2x_degree_lg(x,lx);
     818             :   }
     819    41078029 :   return F2x_renormalize(x, lx);
     820             : }
     821             : 
     822             : GEN
     823    11169499 : F2x_divrem(GEN x, GEN y, GEN *pr)
     824             : {
     825    11169499 :   long dx, dy, dz, lx = lg(x), vs = x[1];
     826             :   GEN z;
     827             : 
     828    11169499 :   dy = F2x_degree(y);
     829    11169179 :   if (dy<0) pari_err_INV("F2x_divrem",y);
     830    11169162 :   if (pr == ONLY_REM) return F2x_rem(x, y);
     831    11169162 :   if (!dy)
     832             :   {
     833     1491748 :     z = F2x_copy(x);
     834     1491848 :     if (pr && pr != ONLY_DIVIDES) *pr = pol0_F2x(vs);
     835     1491848 :     return z;
     836             :   }
     837     9677414 :   dx = F2x_degree_lg(x,lx);
     838     9678041 :   dz = dx-dy;
     839     9678041 :   if (dz < 0)
     840             :   {
     841          14 :     if (pr == ONLY_DIVIDES) return dx < 0? F2x_copy(x): NULL;
     842          14 :     z = pol0_F2x(vs);
     843          14 :     if (pr) *pr = F2x_copy(x);
     844          14 :     return z;
     845             :   }
     846     9678027 :   z = zero_zv(lg(x)-lg(y)+2); z[1] = vs;
     847     9680105 :   x = F2x_copy(x);
     848    31982741 :   while (dx>=dy)
     849             :   {
     850    22305170 :     F2x_set(z,dx-dy);
     851    22273818 :     F2x_addshiftip(x,y,dx-dy);
     852    23258082 :     while (lx>2 && x[lx-1]==0) lx--;
     853    22226726 :     dx = F2x_degree_lg(x,lx);
     854             :   }
     855     9677571 :   z = F2x_renormalize(z, lg(z));
     856     9680067 :   if (!pr) { cgiv(x); return z; }
     857     9052210 :   x = F2x_renormalize(x, lx);
     858     9052113 :   if (pr == ONLY_DIVIDES) {
     859           0 :     if (lg(x) == 2) { cgiv(x); return z; }
     860           0 :     set_avma((pari_sp)(z + lg(z))); return NULL;
     861             :   }
     862     9052113 :   *pr = x; return z;
     863             : }
     864             : 
     865             : long
     866       12384 : F2x_valrem(GEN x, GEN *Z)
     867             : {
     868       12384 :   long v, v2, i, l=lg(x);
     869             :   GEN y;
     870       12384 :   if (l==2) { *Z = F2x_copy(x); return LONG_MAX; }
     871       12386 :   for (i=2; i<l && x[i]==0; i++) /*empty*/;
     872       12384 :   v = i-2;
     873       12384 :   v2 = (i < l)? vals(x[i]): 0;
     874       12384 :   if (v+v2 == 0) { *Z = x; return 0; }
     875        3546 :   l -= i-2;
     876        3546 :   y = cgetg(l, t_VECSMALL); y[1] = x[1];
     877        3546 :   if (v2 == 0)
     878           4 :     for (i=2; i<l; i++) y[i] = x[i+v];
     879        3544 :   else if (l == 3)
     880        3526 :     y[2] = ((ulong)x[2+v]) >> v2;
     881             :   else
     882             :   {
     883          18 :     const ulong sh = BITS_IN_LONG - v2;
     884          18 :     ulong r = x[2+v];
     885          36 :     for (i=3; i<l; i++)
     886             :     {
     887          18 :       y[i-1] = (x[i+v] << sh) | (r >> v2);
     888          18 :       r = x[i+v];
     889             :     }
     890          18 :     y[l-1] = r >> v2;
     891          18 :     (void)F2x_renormalize(y,l);
     892             :   }
     893        3546 :   *Z = y; return (v << TWOPOTBITS_IN_LONG) + v2;
     894             : }
     895             : 
     896             : GEN
     897           0 : F2x_deflate(GEN x, long d)
     898             : {
     899             :   GEN y;
     900           0 :   long i,id, dy, dx = F2x_degree(x);
     901           0 :   if (d <= 1) return F2x_copy(x);
     902           0 :   if (dx < 0) return F2x_copy(x);
     903           0 :   dy = dx/d; /* dy+1 coefficients + 1 extra word for variable */
     904           0 :   y = zero_zv(nbits2lg(dy+1)-1); y[1] = x[1];
     905           0 :   for (i=id=0; i<=dy; i++,id+=d)
     906           0 :     if (F2x_coeff(x,id)) F2x_set(y, i);
     907           0 :   return y;
     908             : }
     909             : 
     910             : /* write p(X) = e(X^2) + Xo(X^2), shallow function */
     911             : void
     912       49949 : F2x_even_odd(GEN p, GEN *pe, GEN *po)
     913             : {
     914       49949 :   long n = F2x_degree(p), n0, n1, i;
     915             :   GEN p0, p1;
     916             : 
     917       49951 :   if (n <= 0) { *pe = F2x_copy(p); *po = pol0_F2x(p[1]); return; }
     918             : 
     919       46304 :   n0 = (n>>1)+1; n1 = n+1 - n0; /* n1 <= n0 <= n1+1 */
     920       46304 :   p0 = zero_zv(nbits2lg(n0+1)-1); p0[1] = p[1];
     921       46305 :   p1 = zero_zv(nbits2lg(n1+1)-1); p1[1] = p[1];
     922     1466009 :   for (i=0; i<n1; i++)
     923             :   {
     924     1419690 :     if (F2x_coeff(p,i<<1)) F2x_set(p0,i);
     925     1419287 :     if (F2x_coeff(p,1+(i<<1))) F2x_set(p1,i);
     926             :   }
     927       46319 :   if (n1 != n0 && F2x_coeff(p,i<<1)) F2x_set(p0,i);
     928       46319 :   *pe = F2x_renormalize(p0,lg(p0));
     929       46318 :   *po = F2x_renormalize(p1,lg(p1));
     930             : }
     931             : 
     932             : GEN
     933      867351 : F2x_deriv(GEN z)
     934             : {
     935      867351 :   const ulong mask = ULONG_MAX/3UL;
     936      867351 :   long i,l = lg(z);
     937      867351 :   GEN x = cgetg(l, t_VECSMALL); x[1] = z[1];
     938     1740098 :   for (i=2; i<l; i++) x[i] = (((ulong)z[i])>>1)&mask;
     939      867354 :   return F2x_renormalize(x,l);
     940             : }
     941             : 
     942             : GEN
     943     2906343 : F2x_gcd(GEN a, GEN b)
     944             : {
     945     2906343 :   pari_sp av = avma;
     946     2906343 :   if (lg(b) > lg(a)) swap(a, b);
     947    20313742 :   while (lgpol(b))
     948             :   {
     949    17076026 :     GEN c = F2x_rem(a,b);
     950    17407399 :     a = b; b = c;
     951    17407399 :     if (gc_needed(av,2))
     952             :     {
     953           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"F2x_gcd (d = %ld)",F2x_degree(c));
     954           0 :       gerepileall(av,2, &a,&b);
     955             :     }
     956             :   }
     957     2897858 :   if (gc_needed(av,2)) a = gerepileuptoleaf(av, a);
     958     2897858 :   return a;
     959             : }
     960             : 
     961             : GEN
     962     1290125 : F2x_extgcd(GEN a, GEN b, GEN *ptu, GEN *ptv)
     963             : {
     964     1290125 :   pari_sp av=avma;
     965             :   GEN u,v,d,d1,v1;
     966     1290125 :   long vx = a[1];
     967     1290125 :   d = a; d1 = b;
     968     1290125 :   v = pol0_F2x(vx); v1 = pol1_F2x(vx);
     969    11623976 :   while (lgpol(d1))
     970             :   {
     971    10333852 :     GEN r, q = F2x_divrem(d,d1, &r);
     972    10333858 :     v = F2x_add(v,F2x_mul(q,v1));
     973    10333849 :     u=v; v=v1; v1=u;
     974    10333849 :     u=r; d=d1; d1=u;
     975    10333849 :     if (gc_needed(av,2))
     976             :     {
     977           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"F2x_extgcd (d = %ld)",F2x_degree(d));
     978           0 :       gerepileall(av,5, &d,&d1,&u,&v,&v1);
     979             :     }
     980             :   }
     981     1290112 :   if (ptu) *ptu = F2x_div(F2x_add(d, F2x_mul(b,v)), a);
     982     1290112 :   *ptv = v;
     983     1290112 :   if (gc_needed(av,2)) gerepileall(av,ptu?3:2,&d,ptv,ptu);
     984     1290124 :   return d;
     985             : }
     986             : 
     987             : static GEN
     988         509 : F2x_halfgcd_i(GEN a, GEN b)
     989             : {
     990         509 :   pari_sp av=avma;
     991             :   GEN u,u1,v,v1;
     992         509 :   long vx = a[1];
     993         509 :   long n = (F2x_degree(a)+1)>>1;
     994         509 :   u1 = v = pol0_F2x(vx);
     995         509 :   u = v1 = pol1_F2x(vx);
     996        8297 :   while (F2x_degree(b)>=n)
     997             :   {
     998        7788 :     GEN r, q = F2x_divrem(a,b, &r);
     999        7788 :     a = b; b = r; swap(u,u1); swap(v,v1);
    1000        7788 :     u1 = F2x_add(u1, F2x_mul(u, q));
    1001        7788 :     v1 = F2x_add(v1, F2x_mul(v, q));
    1002        7788 :     if (gc_needed(av,2))
    1003             :     {
    1004           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"F2x_halfgcd (d = %ld)",F2x_degree(b));
    1005           0 :       gerepileall(av,6, &a,&b,&u1,&v1,&u,&v);
    1006             :     }
    1007             :   }
    1008         509 :   return gerepilecopy(av, mkmat2(mkcol2(u,u1), mkcol2(v,v1)));
    1009             : }
    1010             : 
    1011             : GEN
    1012         509 : F2x_halfgcd(GEN x, GEN y)
    1013             : {
    1014             :   pari_sp av;
    1015             :   GEN M,q,r;
    1016         509 :   if (F2x_degree(y)<F2x_degree(x)) return F2x_halfgcd_i(x,y);
    1017         509 :   av = avma;
    1018         509 :   q = F2x_divrem(y,x,&r);
    1019         509 :   M = F2x_halfgcd_i(x,r);
    1020         509 :   gcoeff(M,1,1) = F2x_add(gcoeff(M,1,1), F2x_mul(q, gcoeff(M,1,2)));
    1021         509 :   gcoeff(M,2,1) = F2x_add(gcoeff(M,2,1), F2x_mul(q, gcoeff(M,2,2)));
    1022         509 :   return gerepilecopy(av, M);
    1023             : }
    1024             : 
    1025             : GEN
    1026     9173955 : F2xq_mul(GEN x,GEN y,GEN pol)
    1027             : {
    1028     9173955 :   GEN z = F2x_mul(x,y);
    1029     9174522 :   return F2x_rem(z,pol);
    1030             : }
    1031             : 
    1032             : GEN
    1033     7647151 : F2xq_sqr(GEN x,GEN pol)
    1034             : {
    1035     7647151 :   GEN z = F2x_sqr(x);
    1036     7676922 :   return F2x_rem(z,pol);
    1037             : }
    1038             : 
    1039             : GEN
    1040     1290125 : F2xq_invsafe(GEN x, GEN T)
    1041             : {
    1042     1290125 :   GEN V, z = F2x_extgcd(get_F2x_mod(T), x, NULL, &V);
    1043     1290124 :   if (F2x_degree(z)) return NULL;
    1044     1290054 :   return V;
    1045             : }
    1046             : 
    1047             : GEN
    1048     1290118 : F2xq_inv(GEN x,GEN T)
    1049             : {
    1050     1290118 :   pari_sp av=avma;
    1051     1290118 :   GEN U = F2xq_invsafe(x, T);
    1052     1290117 :   if (!U) pari_err_INV("F2xq_inv", F2x_to_ZX(x));
    1053     1290047 :   return gerepileuptoleaf(av, U);
    1054             : }
    1055             : 
    1056             : GEN
    1057      575852 : F2xq_div(GEN x,GEN y,GEN T)
    1058             : {
    1059      575852 :   pari_sp av = avma;
    1060      575852 :   return gerepileuptoleaf(av, F2xq_mul(x,F2xq_inv(y,T),T));
    1061             : }
    1062             : 
    1063             : static GEN
    1064      551672 : _F2xq_red(void *E, GEN x) { return F2x_rem(x, (GEN)E); }
    1065             : static GEN
    1066     1422854 : _F2xq_add(void *E, GEN x, GEN y) { (void)E; return F2x_add(x,y); }
    1067             : 
    1068             : static GEN
    1069     1806263 : _F2xq_cmul(void *E, GEN P, long a, GEN x)
    1070             : {
    1071     1806263 :   GEN pol = (GEN) E;
    1072     1806263 :   return F2x_coeff(P,a) ? x: pol0_F2x(pol[1]);
    1073             : }
    1074             : static GEN
    1075      656834 : _F2xq_sqr(void *E, GEN x) { return F2xq_sqr(x, (GEN) E); }
    1076             : static GEN
    1077     1318096 : _F2xq_mul(void *E, GEN x, GEN y) { return F2xq_mul(x,y, (GEN) E); }
    1078             : 
    1079             : static GEN
    1080      884381 : _F2xq_one(void *E)
    1081             : {
    1082      884381 :   GEN pol = (GEN) E;
    1083      884381 :   return pol1_F2x(pol[1]);
    1084             : }
    1085             : static GEN
    1086      259714 : _F2xq_zero(void *E)
    1087             : {
    1088      259714 :   GEN pol = (GEN) E;
    1089      259714 :   return pol0_F2x(pol[1]);
    1090             : }
    1091             : 
    1092             : GEN
    1093      129581 : F2xq_pow(GEN x, GEN n, GEN pol)
    1094             : {
    1095      129581 :   pari_sp av=avma;
    1096             :   GEN y;
    1097             : 
    1098      129581 :   if (!signe(n)) return pol1_F2x(x[1]);
    1099      121888 :   if (is_pm1(n)) /* +/- 1 */
    1100       42204 :     return (signe(n) < 0)? F2xq_inv(x,pol): F2x_copy(x);
    1101             : 
    1102       79684 :   if (signe(n) < 0) x = F2xq_inv(x,pol);
    1103       79684 :   y = gen_pow_i(x, n, (void*)pol, &_F2xq_sqr, &_F2xq_mul);
    1104       79684 :   return gerepilecopy(av, y);
    1105             : }
    1106             : 
    1107             : GEN
    1108          49 : F2xq_powu(GEN x, ulong n, GEN pol)
    1109             : {
    1110          49 :   pari_sp av=avma;
    1111             :   GEN y;
    1112          49 :   switch(n)
    1113             :   {
    1114           0 :     case 0: return pol1_F2x(x[1]);
    1115           7 :     case 1: return F2x_copy(x);
    1116          14 :     case 2: return F2xq_sqr(x,pol);
    1117             :   }
    1118          28 :   y = gen_powu_i(x, n, (void*)pol, &_F2xq_sqr, &_F2xq_mul);
    1119          28 :   return gerepilecopy(av, y);
    1120             : }
    1121             : 
    1122             : GEN
    1123          14 : F2xq_pow_init(GEN x, GEN n, long k,  GEN T)
    1124             : {
    1125          14 :   return gen_pow_init(x, n, k, (void*)T, &_F2xq_sqr, &_F2xq_mul);
    1126             : }
    1127             : 
    1128             : GEN
    1129        5246 : F2xq_pow_table(GEN R, GEN n, GEN T)
    1130             : {
    1131        5246 :   return gen_pow_table(R, n, (void*)T, &_F2xq_one, &_F2xq_mul);
    1132             : }
    1133             : 
    1134             : /* generates the list of powers of x of degree 0,1,2,...,l*/
    1135             : GEN
    1136      353265 : F2xq_powers(GEN x, long l, GEN T)
    1137             : {
    1138      353265 :   int use_sqr = 2*F2x_degree(x) >= get_F2x_degree(T);
    1139      353263 :   return gen_powers(x, l, use_sqr, (void*)T, &_F2xq_sqr, &_F2xq_mul, &_F2xq_one);
    1140             : }
    1141             : 
    1142             : GEN
    1143      223807 : F2xq_matrix_pow(GEN y, long n, long m, GEN P)
    1144             : {
    1145      223807 :   return F2xV_to_F2m(F2xq_powers(y,m-1,P),n);
    1146             : }
    1147             : 
    1148             : GEN
    1149      381979 : F2x_Frobenius(GEN T)
    1150             : {
    1151      381979 :   return F2xq_sqr(polx_F2x(get_F2x_var(T)), T);
    1152             : }
    1153             : 
    1154             : GEN
    1155      223807 : F2x_matFrobenius(GEN T)
    1156             : {
    1157      223807 :   long n = get_F2x_degree(T);
    1158      223806 :   return F2xq_matrix_pow(F2x_Frobenius(T), n, n, T);
    1159             : }
    1160             : 
    1161             : static struct bb_algebra F2xq_algebra = { _F2xq_red, _F2xq_add, _F2xq_add,
    1162             :               _F2xq_mul, _F2xq_sqr, _F2xq_one, _F2xq_zero};
    1163             : 
    1164             : GEN
    1165      609903 : F2x_F2xqV_eval(GEN Q, GEN x, GEN T)
    1166             : {
    1167      609903 :   long d = F2x_degree(Q);
    1168      609903 :   return gen_bkeval_powers(Q,d,x,(void*)T,&F2xq_algebra,_F2xq_cmul);
    1169             : }
    1170             : 
    1171             : GEN
    1172       40884 : F2x_F2xq_eval(GEN Q, GEN x, GEN T)
    1173             : {
    1174       40884 :   long d = F2x_degree(Q);
    1175       40880 :   int use_sqr = 2*F2x_degree(x) >= get_F2x_degree(T);
    1176       40882 :   return gen_bkeval(Q, d, x, use_sqr, (void*)T, &F2xq_algebra, _F2xq_cmul);
    1177             : }
    1178             : 
    1179             : static GEN
    1180       26924 : F2xq_autpow_sqr(void * T, GEN x) { return F2x_F2xq_eval(x, x, (GEN) T); }
    1181             : 
    1182             : static GEN
    1183       13232 : F2xq_autpow_mul(void * T, GEN x, GEN y) { return F2x_F2xq_eval(x, y, (GEN) T); }
    1184             : 
    1185             : GEN
    1186       18163 : F2xq_autpow(GEN x, long n, GEN T)
    1187             : {
    1188       18163 :   if (n==0) return F2x_rem(polx_F2x(x[1]), T);
    1189       18163 :   if (n==1) return F2x_rem(x, T);
    1190       18163 :   return gen_powu(x,n,(void*)T,F2xq_autpow_sqr,F2xq_autpow_mul);
    1191             : }
    1192             : 
    1193             : ulong
    1194      457826 : F2xq_trace(GEN x, GEN T)
    1195             : {
    1196      457826 :   pari_sp av = avma;
    1197      457826 :   long n = get_F2x_degree(T)-1;
    1198      457826 :   GEN z = F2xq_mul(x, F2x_deriv(get_F2x_mod(T)), T);
    1199      457826 :   return gc_ulong(av, F2x_degree(z) < n ? 0 : 1);
    1200             : }
    1201             : 
    1202             : GEN
    1203           7 : F2xq_conjvec(GEN x, GEN T)
    1204             : {
    1205           7 :   long i, l = 1+get_F2x_degree(T);
    1206           7 :   GEN z = cgetg(l,t_COL);
    1207           7 :   gel(z,1) = F2x_copy(x);
    1208         140 :   for (i=2; i<l; i++) gel(z,i) = F2xq_sqr(gel(z,i-1), T);
    1209           7 :   return z;
    1210             : }
    1211             : 
    1212             : static GEN
    1213       86343 : _F2xq_pow(void *data, GEN x, GEN n)
    1214             : {
    1215       86343 :   GEN pol = (GEN) data;
    1216       86343 :   return F2xq_pow(x,n, pol);
    1217             : }
    1218             : 
    1219             : static GEN
    1220        9037 : _F2xq_rand(void *data)
    1221             : {
    1222        9037 :   pari_sp av=avma;
    1223        9037 :   GEN pol = (GEN) data;
    1224        9037 :   long d = F2x_degree(pol);
    1225             :   GEN z;
    1226             :   do
    1227             :   {
    1228        9191 :     set_avma(av);
    1229        9191 :     z = random_F2x(d,pol[1]);
    1230        9191 :   } while (lgpol(z)==0);
    1231        9037 :   return z;
    1232             : }
    1233             : 
    1234             : static GEN F2xq_easylog(void* E, GEN a, GEN g, GEN ord);
    1235             : 
    1236             : static const struct bb_group F2xq_star={_F2xq_mul,_F2xq_pow,_F2xq_rand,hash_GEN,F2x_equal,F2x_equal1,F2xq_easylog};
    1237             : 
    1238             : GEN
    1239         378 : F2xq_order(GEN a, GEN ord, GEN T)
    1240             : {
    1241         378 :   return gen_order(a,ord,(void*)T,&F2xq_star);
    1242             : }
    1243             : 
    1244             : static long
    1245      261425 : F2x_is_smooth_squarefree(GEN f, long r)
    1246             : {
    1247      261425 :   pari_sp av = avma;
    1248      261425 :   long i, df = F2x_degree(f);
    1249             :   GEN sx, a;
    1250      261191 :   if (!df) return 1;
    1251      256128 :   a = sx = polx_F2x(f[1]);
    1252      256474 :   for(i = 1;; i++)
    1253     2110794 :   {
    1254             :     long dg;
    1255             :     GEN g;
    1256     2367268 :     a = F2xq_sqr(a, f);
    1257     2366024 :     if (F2x_equal(a, sx)) return gc_long(av,1);
    1258     2302636 :     if (i==r) return gc_long(av,0);
    1259     2137633 :     g = F2x_gcd(f, F2x_add(a,sx));
    1260     2143534 :     dg = F2x_degree(g);
    1261     2132860 :     if (dg == df) return gc_long(av,1);
    1262     2110407 :     if (dg) { f = F2x_div(f, g); df -= dg; a = F2x_rem(a, f); }
    1263             :   }
    1264             : }
    1265             : 
    1266             : static long
    1267      224660 : F2x_is_smooth(GEN g, long r)
    1268             : {
    1269      224660 :   if (lgpol(g)==0) return 0;
    1270             :   while (1)
    1271       37053 :   {
    1272      261684 :     GEN f = F2x_gcd(g, F2x_deriv(g));
    1273      261471 :     if (!F2x_is_smooth_squarefree(F2x_div(g, f), r)) return 0;
    1274       96662 :     if (F2x_degree(f)==0) return 1;
    1275       36978 :     g = F2x_issquare(f) ? F2x_sqrt(f): f;
    1276             :   }
    1277             : }
    1278             : 
    1279             : static GEN
    1280       15183 : F2x_factorel(GEN h)
    1281             : {
    1282       15183 :   GEN F = F2x_factor(h);
    1283       15184 :   GEN F1 = gel(F, 1), F2 = gel(F, 2);
    1284       15184 :   long i, l1 = lg(F1)-1;
    1285       15184 :   GEN p2 = cgetg(l1+1, t_VECSMALL);
    1286       15183 :   GEN e2 = cgetg(l1+1, t_VECSMALL);
    1287       71121 :   for (i = 1; i <= l1; ++i)
    1288             :   {
    1289       55938 :     p2[i] = mael(F1, i, 2);
    1290       55938 :     e2[i] = F2[i];
    1291             :   }
    1292       15183 :   return mkmat2(p2, e2);
    1293             : }
    1294             : 
    1295             : static GEN
    1296       25837 : mkF2(ulong x, long v) { return mkvecsmall2(v, x); }
    1297             : 
    1298             : static GEN F2xq_log_Coppersmith_d(GEN W, GEN g, long r, long n, GEN T, GEN mo);
    1299             : 
    1300             : static GEN
    1301          46 : F2xq_log_from_rel(GEN W, GEN rel, long r, long n, GEN T, GEN m)
    1302             : {
    1303          46 :   pari_sp av = avma;
    1304          46 :   long vT = get_F2x_var(T);
    1305          46 :   GEN F = gel(rel,1), E = gel(rel,2), o = gen_0;
    1306          46 :   long i, l = lg(F);
    1307         417 :   for(i=1; i<l; i++)
    1308             :   {
    1309         371 :     GEN R = gel(W, F[i]);
    1310         371 :     if (signe(R)==0) /* Already failed */
    1311           0 :       return NULL;
    1312         371 :     else if (signe(R)<0) /* Not yet tested */
    1313             :     {
    1314           1 :       setsigne(gel(W,F[i]),0);
    1315           1 :       R = F2xq_log_Coppersmith_d(W, mkF2(F[i],vT), r, n, T, m);
    1316           1 :       if (!R) return NULL;
    1317             :     }
    1318         371 :     o = Fp_add(o, mulis(R, E[i]), m);
    1319             :   }
    1320          46 :   return gerepileuptoint(av, o);
    1321             : }
    1322             : 
    1323             : static GEN
    1324          58 : F2xq_log_Coppersmith_d(GEN W, GEN g, long r, long n, GEN T, GEN mo)
    1325             : {
    1326          58 :   pari_sp av = avma, av2;
    1327          58 :   long dT = get_F2x_degree(T), vT = get_F2x_var(T);
    1328          58 :   long dg = F2x_degree(g), k = r-1, m = maxss((dg-k)/2,0);
    1329          58 :   long i, j, l = dg-m, N;
    1330          58 :   GEN v = cgetg(k+m+1,t_MAT);
    1331          58 :   long h = dT>>n, d = dT-(h<<n);
    1332          58 :   GEN p1 = pol1_F2x(vT);
    1333          58 :   GEN R = F2x_add(F2x_shift(p1, dT), T);
    1334          58 :   GEN z = F2x_rem(F2x_shift(p1, h), g);
    1335         731 :   for(i=1; i<=k+m; i++)
    1336             :   {
    1337         673 :     gel(v,i) = F2x_to_F2v(F2x_shift(z,-l),m);
    1338         673 :     z = F2x_rem(F2x_shift(z,1),g);
    1339             :   }
    1340          58 :   v = F2m_ker(v);
    1341         638 :   for(i=1; i<=k; i++)
    1342         580 :     gel(v,i) = F2v_to_F2x(gel(v,i),vT);
    1343          58 :   N = 1<<k;
    1344          58 :   av2 = avma;
    1345       17358 :   for (i=1; i<N; i++)
    1346             :   {
    1347             :     GEN p,q,qh,a,b;
    1348       17346 :     set_avma(av2);
    1349       17346 :     q = pol0_F2x(vT);
    1350      190806 :     for(j=0; j<k; j++)
    1351      173460 :       if (i&(1UL<<j))
    1352       79681 :         q = F2x_add(q, gel(v,j+1));
    1353       17346 :     qh= F2x_shift(q,h);
    1354       17346 :     p = F2x_rem(qh,g);
    1355       17346 :     b = F2x_add(F2x_mul(R, F2x_pow2n(q, n)), F2x_shift(F2x_pow2n(p, n), d));
    1356       17346 :     if (lgpol(b)==0 || !F2x_is_smooth(b, r)) continue;
    1357          70 :     a = F2x_div(F2x_add(qh,p),g);
    1358          70 :     if (F2x_degree(F2x_gcd(a,q)) &&  F2x_degree(F2x_gcd(a,p))) continue;
    1359          58 :     if (!(lgpol(a)==0 || !F2x_is_smooth(a, r)))
    1360             :     {
    1361          46 :       GEN F = F2x_factorel(b);
    1362          46 :       GEN G = F2x_factorel(a);
    1363          46 :       GEN FG = vecsmall_concat(vecsmall_append(gel(F, 1), 2), gel(G, 1));
    1364          46 :       GEN E  = vecsmall_concat(vecsmall_append(gel(F, 2), -d),
    1365          46 :                                zv_z_mul(gel(G, 2),-(1L<<n)));
    1366          46 :       GEN R  = famatsmall_reduce(mkmat2(FG, E));
    1367          46 :       GEN l  = F2xq_log_from_rel(W, R, r, n, T, mo);
    1368          46 :       if (!l) continue;
    1369          46 :       l = Fp_div(l,int2n(n),mo);
    1370          46 :       if (dg <= r)
    1371             :       {
    1372           8 :         affii(l,gel(W,g[2]));
    1373           8 :         if (DEBUGLEVEL>1) err_printf("Found %lu\n", g[2]);
    1374             :       }
    1375          46 :       return gerepileuptoint(av, l);
    1376             :     }
    1377             :   }
    1378          12 :   return gc_NULL(av);
    1379             : }
    1380             : 
    1381             : static GEN
    1382          52 : F2xq_log_find_rel(GEN b, long r, GEN T, GEN *g, ulong *e)
    1383             : {
    1384          52 :   pari_sp av = avma;
    1385             :   while (1)
    1386         457 :   {
    1387             :     GEN M;
    1388         509 :     *g = F2xq_mul(*g, b, T); (*e)++;
    1389         509 :     M = F2x_halfgcd(*g,T);
    1390         509 :     if (F2x_is_smooth(gcoeff(M,1,1), r))
    1391             :     {
    1392         266 :       GEN z = F2x_add(F2x_mul(gcoeff(M,1,1),*g), F2x_mul(gcoeff(M,1,2),T));
    1393         266 :       if (F2x_is_smooth(z, r))
    1394             :       {
    1395          52 :         GEN F = F2x_factorel(z);
    1396          52 :         GEN G = F2x_factorel(gcoeff(M,1,1));
    1397          52 :         GEN rel = mkmat2(vecsmall_concat(gel(F, 1),gel(G, 1)),
    1398          52 :                          vecsmall_concat(gel(F, 2),zv_neg(gel(G, 2))));
    1399          52 :         gerepileall(av, 2, g, &rel);
    1400          52 :         return rel;
    1401             :       }
    1402             :     }
    1403         457 :     if (gc_needed(av,2))
    1404             :     {
    1405           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"F2xq_log_find_rel");
    1406           0 :       *g = gerepileuptoleaf(av, *g);
    1407             :     }
    1408             :   }
    1409             : }
    1410             : 
    1411             : static GEN
    1412          28 : F2xq_log_Coppersmith_rec(GEN W, long r2, GEN a, long r, long n, GEN T, GEN m)
    1413             : {
    1414          28 :   long vT = get_F2x_var(T);
    1415          28 :   GEN b = polx_F2x(vT);
    1416          28 :   ulong AV = 0;
    1417          28 :   GEN g = a, bad = pol0_F2x(vT);
    1418             :   pari_timer ti;
    1419             :   while(1)
    1420          24 :   {
    1421             :     long i, l;
    1422             :     GEN V, F, E, Ao;
    1423          52 :     timer_start(&ti);
    1424          52 :     V = F2xq_log_find_rel(b, r2, T, &g, &AV);
    1425          52 :     if (DEBUGLEVEL>1) timer_printf(&ti,"%ld-smooth element",r2);
    1426          52 :     F = gel(V,1); E = gel(V,2);
    1427          52 :     l = lg(F);
    1428          52 :     Ao = gen_0;
    1429         335 :     for(i=1; i<l; i++)
    1430             :     {
    1431         307 :       GEN Fi = mkF2(F[i], vT);
    1432             :       GEN R;
    1433         307 :       if (F2x_degree(Fi) <= r)
    1434             :       {
    1435         245 :         if (signe(gel(W,F[i]))==0)
    1436           0 :           break;
    1437         245 :         else if (signe(gel(W,F[i]))<0)
    1438             :         {
    1439           7 :           setsigne(gel(W,F[i]),0);
    1440           7 :           R = F2xq_log_Coppersmith_d(W,Fi,r,n,T,m);
    1441             :         } else
    1442         238 :           R = gel(W,F[i]);
    1443             :       }
    1444             :       else
    1445             :       {
    1446          62 :         if (F2x_equal(Fi,bad)) break;
    1447          50 :         R = F2xq_log_Coppersmith_d(W,Fi,r,n,T,m);
    1448          50 :         if (!R) bad = Fi;
    1449             :       }
    1450         295 :       if (!R) break;
    1451         283 :       Ao = Fp_add(Ao, mulis(R, E[i]), m);
    1452             :     }
    1453          52 :     if (i==l) return subiu(Ao,AV);
    1454             :   }
    1455             : }
    1456             : 
    1457             : /* Coppersmith:
    1458             :  T*X^e = X^(h*2^n)-R
    1459             :  (u*x^h + v)^(2^n) = u^(2^n)*X^(h*2^n)+v^(2^n)
    1460             :  (u*x^h + v)^(2^n) = u^(2^n)*R+v^(2^n)
    1461             : */
    1462             : 
    1463             : static GEN
    1464      154816 : rel_Coppersmith(GEN u, GEN v, long h, GEN R, long r, long n, long d)
    1465             : {
    1466             :   GEN b, F, G, M;
    1467      154816 :   GEN a = F2x_add(F2x_shift(u, h), v);
    1468      154835 :   if (!F2x_is_smooth(a, r)) return NULL;
    1469       51724 :   b = F2x_add(F2x_mul(R, F2x_pow2n(u, n)), F2x_shift(F2x_pow2n(v, n),d));
    1470       51738 :   if (!F2x_is_smooth(b, r)) return NULL;
    1471        7493 :   F = F2x_factorel(a);
    1472        7494 :   G = F2x_factorel(b);
    1473       14988 :   M = mkmat2(vecsmall_concat(gel(F, 1), vecsmall_append(gel(G, 1), 2)),
    1474       14988 :              vecsmall_concat(zv_z_mul(gel(F, 2),1UL<<n), vecsmall_append(zv_neg(gel(G, 2)),d)));
    1475        7494 :   return famatsmall_reduce(M);
    1476             : }
    1477             : 
    1478             : GEN
    1479        2112 : F2xq_log_Coppersmith_worker(GEN u, long i, GEN V, GEN R)
    1480             : {
    1481        2112 :   long r = V[1], h = V[2], n = V[3], d = V[4];
    1482        2112 :   pari_sp ltop = avma;
    1483        2112 :   GEN v = mkF2(0,u[1]);
    1484        2107 :   GEN L = cgetg(2*i+1, t_VEC);
    1485        2111 :   pari_sp av = avma;
    1486             :   long j;
    1487        2111 :   long nbtest=0, rel = 1;
    1488      155451 :   for(j=1; j<=i; j++)
    1489             :   {
    1490      153330 :     v[2] = j;
    1491      153330 :     set_avma(av);
    1492      153326 :     if (F2x_degree(F2x_gcd(u,v))==0)
    1493             :     {
    1494       77427 :       GEN z = rel_Coppersmith(u, v, h, R, r, n, d);
    1495       77482 :       nbtest++;
    1496       77482 :       if (z) { gel(L, rel++) = z; av = avma; }
    1497       77482 :       if (i==j) continue;
    1498       77468 :       z = rel_Coppersmith(v, u, h, R, r, n, d);
    1499       77467 :       nbtest++;
    1500       77467 :       if (z) { gel(L, rel++) = z; av = avma; }
    1501             :     }
    1502             :   }
    1503        2121 :   setlg(L,rel);
    1504        2121 :   return gerepilecopy(ltop, mkvec2(stoi(nbtest), L));
    1505             : }
    1506             : 
    1507             : static GEN
    1508          14 : F2xq_log_Coppersmith(long nbrel, long r, long n, GEN T)
    1509             : {
    1510          14 :   long dT = get_F2x_degree(T), vT = get_F2x_var(T);
    1511          14 :   long h = dT>>n, d = dT-(h<<n);
    1512          14 :   GEN R = F2x_add(F2x_shift(pol1_F2x(vT), dT), T);
    1513          14 :   GEN u = mkF2(0,vT);
    1514          14 :   long rel = 1, nbtest = 0;
    1515          14 :   GEN M = cgetg(nbrel+1, t_VEC);
    1516          14 :   long i = 0;
    1517          14 :   GEN worker = snm_closure(is_entry("_F2xq_log_Coppersmith_worker"),
    1518             :                mkvec2(mkvecsmall4(r, h, n, d), R));
    1519             :   struct pari_mt pt;
    1520          14 :   long running, pending = 0, stop=0;
    1521          14 :   mt_queue_start(&pt, worker);
    1522          14 :   if (DEBUGLEVEL) err_printf("Coppersmith (R = %ld): ",F2x_degree(R));
    1523        2197 :   while ((running = !stop) || pending)
    1524             :   {
    1525             :     GEN L, done;
    1526             :     long j, l;
    1527        2183 :     u[2] = i;
    1528        2183 :     mt_queue_submit(&pt, 0, running ? mkvec2(u, stoi(i)): NULL);
    1529        2183 :     done = mt_queue_get(&pt, NULL, &pending);
    1530        2183 :     if (!done) continue;
    1531        2121 :     L = gel(done, 2); nbtest += itos(gel(done,1));
    1532        2121 :     l = lg(L);
    1533        2121 :     if (l > 1)
    1534             :     {
    1535        9155 :       for (j=1; j<l; j++)
    1536             :       {
    1537        7277 :         if (rel>nbrel) break;
    1538        7210 :         gel(M,rel++) = gel(L,j);
    1539        7210 :         if (DEBUGLEVEL && (rel&511UL)==0)
    1540           0 :           err_printf("%ld%%[%ld] ",rel*100/nbrel,i);
    1541             :       }
    1542             :     }
    1543        2121 :     if (rel>nbrel) stop=1;
    1544        2121 :     i++;
    1545             :   }
    1546          14 :   mt_queue_end(&pt);
    1547          14 :   if (DEBUGLEVEL) err_printf(": %ld tests\n", nbtest);
    1548          14 :   return M;
    1549             : }
    1550             : 
    1551             : static GEN
    1552          14 : smallirred_F2x(ulong n, long sv)
    1553             : {
    1554          14 :   GEN a = zero_zv(nbits2lg(n+1)-1);
    1555          14 :   a[1] = sv; F2x_set(a,n); a[2]++;
    1556         238 :   while (!F2x_is_irred(a)) a[2]+=2;
    1557          14 :   return a;
    1558             : }
    1559             : 
    1560             : static GEN
    1561          14 : check_kernel(long N, GEN M, long nbi, GEN T, GEN m)
    1562             : {
    1563          14 :   pari_sp av = avma;
    1564          14 :   long dT = get_F2x_degree(T), vT = get_F2x_var(T);
    1565          14 :   GEN K = FpMs_leftkernel_elt(M, N, m);
    1566          14 :   long i, f=0, tbs;
    1567          14 :   long l = lg(K), lm = lgefint(m);
    1568          14 :   GEN idx = diviiexact(int2um1(dT), m);
    1569          14 :   GEN g = F2xq_pow(polx_F2x(vT), idx, T);
    1570             :   GEN tab;
    1571             :   pari_timer ti;
    1572          14 :   if (DEBUGLEVEL) timer_start(&ti);
    1573          14 :   K = FpC_Fp_mul(K, Fp_inv(gel(K,2), m), m);
    1574          14 :   tbs = maxss(1, expu(nbi/expi(m)));
    1575          14 :   tab = F2xq_pow_init(g, int2n(dT), tbs, T);
    1576       57344 :   for(i=1; i<l; i++)
    1577             :   {
    1578       57330 :     GEN k = gel(K,i);
    1579       57330 :     if (signe(k)==0 || !F2x_equal(F2xq_pow_table(tab, k, T),
    1580             :                                   F2xq_pow(mkF2(i,vT), idx, T)))
    1581       52084 :       gel(K,i) = cgetineg(lm);
    1582             :     else
    1583        5246 :       f++;
    1584             :   }
    1585          14 :   if (DEBUGLEVEL) timer_printf(&ti,"found %ld/%ld logs", f, nbi);
    1586          14 :   return gerepileupto(av, K);
    1587             : }
    1588             : 
    1589             : static GEN
    1590          14 : F2xq_log_index(GEN a0, GEN b0, GEN m, GEN T0)
    1591             : {
    1592          14 :   pari_sp av = avma;
    1593          14 :   GEN  M, S, a, b, Ao=NULL, Bo=NULL, W, e;
    1594             :   pari_timer ti;
    1595          14 :   long n = F2x_degree(T0), r = (long) (sqrt((double) 2*n))-(n>100);
    1596          14 :   GEN T = smallirred_F2x(n,T0[1]);
    1597          14 :   long d = 2, r2 = 3*r/2, d2 = 2;
    1598          14 :   long N = (1UL<<(r+1))-1UL;
    1599          14 :   long nbi = itos(ffsumnbirred(gen_2, r)), nbrel=nbi*5/4;
    1600          14 :   if (DEBUGLEVEL)
    1601             :   {
    1602           0 :     err_printf("F2xq_log: Parameters r=%ld r2=%ld\n", r,r2);
    1603           0 :     err_printf("F2xq_log: Size FB=%ld rel. needed=%ld\n", nbi, nbrel);
    1604           0 :     timer_start(&ti);
    1605             :   }
    1606          14 :   S = Flx_to_F2x(Flx_ffisom(F2x_to_Flx(T0),F2x_to_Flx(T),2));
    1607          14 :   a = F2x_F2xq_eval(a0, S, T);
    1608          14 :   b = F2x_F2xq_eval(b0, S, T);
    1609          14 :   if (DEBUGLEVEL) timer_printf(&ti,"model change");
    1610          14 :   M = F2xq_log_Coppersmith(nbrel,r,d,T);
    1611          14 :   if(DEBUGLEVEL)
    1612           0 :     timer_printf(&ti,"relations");
    1613          14 :   W = check_kernel(N, M, nbi, T, m);
    1614          14 :   timer_start(&ti);
    1615          14 :   Ao = F2xq_log_Coppersmith_rec(W, r2, a, r, d2, T, m);
    1616          14 :   if (DEBUGLEVEL) timer_printf(&ti,"smooth element");
    1617          14 :   Bo = F2xq_log_Coppersmith_rec(W, r2, b, r, d2, T, m);
    1618          14 :   if (DEBUGLEVEL) timer_printf(&ti,"smooth generator");
    1619          14 :   e = Fp_div(Ao, Bo, m);
    1620          14 :   if (!F2x_equal(F2xq_pow(b0,e,T0),a0)) pari_err_BUG("F2xq_log");
    1621          14 :   return gerepileupto(av, e);
    1622             : }
    1623             : 
    1624             : static GEN
    1625       10433 : F2xq_easylog(void* E, GEN a, GEN g, GEN ord)
    1626             : {
    1627       10433 :   if (F2x_equal1(a)) return gen_0;
    1628        9831 :   if (F2x_equal(a,g)) return gen_1;
    1629        9012 :   if (typ(ord)!=t_INT) return NULL;
    1630        4856 :   if (expi(ord)<28) return NULL;
    1631          14 :   return F2xq_log_index(a,g,ord,(GEN)E);
    1632             : }
    1633             : 
    1634             : GEN
    1635        7672 : F2xq_log(GEN a, GEN g, GEN ord, GEN T)
    1636             : {
    1637        7672 :   GEN z, v = get_arith_ZZM(ord);
    1638        7672 :   ord = mkvec2(gel(v,1),ZM_famat_limit(gel(v,2),int2n(28)));
    1639        7672 :   z = gen_PH_log(a,g,ord,(void*)T,&F2xq_star);
    1640        7672 :   return z? z: cgetg(1,t_VEC);
    1641             : }
    1642             : 
    1643             : GEN
    1644      212961 : F2xq_Artin_Schreier(GEN a, GEN T)
    1645             : {
    1646      212961 :   pari_sp ltop=avma;
    1647      212961 :   long j,N = get_F2x_degree(T), vT = get_F2x_var(T);
    1648      212961 :   GEN Q = F2x_matFrobenius(T);
    1649     1478582 :   for (j=1; j<=N; j++)
    1650     1265621 :     F2m_flip(Q,j,j);
    1651      212961 :   F2v_add_inplace(gel(Q,1),a);
    1652      212961 :   Q = F2m_ker_sp(Q,0);
    1653      212961 :   if (lg(Q)!=2) return NULL;
    1654      212961 :   Q = gel(Q,1);
    1655      212961 :   Q[1] = vT;
    1656      212961 :   return gerepileuptoleaf(ltop, F2x_renormalize(Q, lg(Q)));
    1657             : }
    1658             : 
    1659             : GEN
    1660       49950 : F2xq_sqrt_fast(GEN c, GEN sqx, GEN T)
    1661             : {
    1662             :   GEN c0, c1;
    1663       49950 :   F2x_even_odd(c, &c0, &c1);
    1664       49965 :   return F2x_add(c0, F2xq_mul(c1, sqx, T));
    1665             : }
    1666             : 
    1667             : static int
    1668       18156 : F2x_is_x(GEN a) { return lg(a)==3 && a[2]==2; }
    1669             : 
    1670             : GEN
    1671       33122 : F2xq_sqrt(GEN a, GEN T)
    1672             : {
    1673       33122 :   pari_sp av = avma;
    1674       33122 :   long n = get_F2x_degree(T), vT = get_F2x_var(T);
    1675             :   GEN sqx;
    1676       33122 :   if (n==1) return F2x_copy(a);
    1677       33045 :   if (n==2) return F2xq_sqr(a,T);
    1678       18156 :   sqx = F2xq_autpow(mkF2(4, vT), n-1, T);
    1679       18156 :   return gerepileuptoleaf(av, F2x_is_x(a)? sqx: F2xq_sqrt_fast(a,sqx,T));
    1680             : }
    1681             : 
    1682             : GEN
    1683        9296 : F2xq_sqrtn(GEN a, GEN n, GEN T, GEN *zeta)
    1684             : {
    1685        9296 :   long dT = get_F2x_degree(T), vT = get_F2x_var(T);
    1686        9296 :   if (!lgpol(a))
    1687             :   {
    1688           7 :     if (signe(n) < 0) pari_err_INV("F2xq_sqrtn",a);
    1689           0 :     if (zeta)
    1690           0 :       *zeta=pol1_F2x(vT);
    1691           0 :     return pol0_F2x(vT);
    1692             :   }
    1693        9289 :   return gen_Shanks_sqrtn(a, n, int2um1(dT), zeta, (void*)T, &F2xq_star);
    1694             : }
    1695             : 
    1696             : GEN
    1697         154 : gener_F2xq(GEN T, GEN *po)
    1698             : {
    1699         154 :   long i, j, vT = get_F2x_var(T), f = get_F2x_degree(T);
    1700         154 :   pari_sp av0 = avma, av;
    1701             :   GEN g, L2, o, q;
    1702             : 
    1703         154 :   if (f == 1) {
    1704          21 :     if (po) *po = mkvec2(gen_1, trivial_fact());
    1705          21 :     return pol1_F2x(vT);
    1706             :   }
    1707         133 :   q = int2um1(f);
    1708         133 :   o = factor_pn_1(gen_2,f);
    1709         133 :   L2 = leafcopy( gel(o, 1) );
    1710         385 :   for (i = j = 1; i < lg(L2); i++)
    1711             :   {
    1712         252 :     if (absequaliu(gel(L2,i),2)) continue;
    1713         252 :     gel(L2,j++) = diviiexact(q, gel(L2,i));
    1714             :   }
    1715         133 :   setlg(L2, j);
    1716         175 :   for (av = avma;; set_avma(av))
    1717             :   {
    1718         175 :     g = random_F2x(f, vT);
    1719         175 :     if (F2x_degree(g) < 1) continue;
    1720         406 :     for (i = 1; i < j; i++)
    1721             :     {
    1722         273 :       GEN a = F2xq_pow(g, gel(L2,i), T);
    1723         273 :       if (F2x_equal1(a)) break;
    1724             :     }
    1725         154 :     if (i == j) break;
    1726             :   }
    1727         133 :   if (!po) g = gerepilecopy(av0, g);
    1728             :   else {
    1729           0 :     *po = mkvec2(int2um1(f), o);
    1730           0 :     gerepileall(av0, 2, &g, po);
    1731             :   }
    1732         133 :   return g;
    1733             : }
    1734             : 
    1735             : static GEN
    1736         798 : _F2xq_neg(void *E, GEN x) { (void) E; return F2x_copy(x); }
    1737             : static GEN
    1738       13748 : _F2xq_rmul(void *E, GEN x, GEN y) { (void) E; return F2x_mul(x,y); }
    1739             : static GEN
    1740         413 : _F2xq_inv(void *E, GEN x) { return F2xq_inv(x, (GEN) E); }
    1741             : static int
    1742         973 : _F2xq_equal0(GEN x) { return lgpol(x)==0; }
    1743             : static GEN
    1744         693 : _F2xq_s(void *E, long x)
    1745         693 : { GEN T = (GEN) E;
    1746         693 :   long v = get_F2x_var(T);
    1747         693 :   return odd(x)? pol1_F2x(v): pol0_F2x(v);
    1748             : }
    1749             : 
    1750             : static const struct bb_field F2xq_field={_F2xq_red,_F2xq_add,_F2xq_rmul,_F2xq_neg,
    1751             :                                          _F2xq_inv,_F2xq_equal0,_F2xq_s};
    1752             : 
    1753        1659 : const struct bb_field *get_F2xq_field(void **E, GEN T)
    1754             : {
    1755        1659 :   *E = (void *) T;
    1756        1659 :   return &F2xq_field;
    1757             : }
    1758             : 
    1759             : /***********************************************************************/
    1760             : /**                                                                   **/
    1761             : /**                               F2xV                                **/
    1762             : /**                                                                   **/
    1763             : /***********************************************************************/
    1764             : /* F2xV are t_VEC with F2x coefficients. */
    1765             : 
    1766             : GEN
    1767       38738 : FlxC_to_F2xC(GEN x) { pari_APPLY_type(t_COL, Flx_to_F2x(gel(x,i))) }
    1768             : GEN
    1769           0 : F2xC_to_FlxC(GEN x) { pari_APPLY_type(t_COL, F2x_to_Flx(gel(x,i))) }
    1770             : GEN
    1771       10738 : F2xC_to_ZXC(GEN x) { pari_APPLY_type(t_COL, F2x_to_ZX(gel(x,i))) }
    1772             : GEN
    1773     1730241 : F2xV_to_F2m(GEN x, long n) { pari_APPLY_type(t_MAT, F2x_to_F2v(gel(x,i), n)) }
    1774             : 
    1775             : void
    1776       20994 : F2xV_to_FlxV_inplace(GEN v)
    1777             : {
    1778       20994 :   long i, l = lg(v);
    1779       45383 :   for(i = 1; i < l;i++) gel(v,i) = F2x_to_Flx(gel(v,i));
    1780       20994 : }
    1781             : void
    1782      166839 : F2xV_to_ZXV_inplace(GEN v)
    1783             : {
    1784      166839 :   long i, l = lg(v);
    1785      366874 :   for(i = 1; i < l; i++) gel(v,i)= F2x_to_ZX(gel(v,i));
    1786      166839 : }
    1787             : 
    1788             : /***********************************************************************/
    1789             : /**                                                                   **/
    1790             : /**                             F2xX                                  **/
    1791             : /**                                                                   **/
    1792             : /***********************************************************************/
    1793             : 
    1794             : GEN
    1795     2710410 : F2xX_renormalize(GEN /*in place*/ x, long lx)
    1796     2710410 : { return FlxX_renormalize(x, lx); }
    1797             : 
    1798             : GEN
    1799      446996 : pol1_F2xX(long v, long sv) { return pol1_FlxX(v, sv); }
    1800             : 
    1801             : GEN
    1802         476 : polx_F2xX(long v, long sv) { return polx_FlxX(v, sv); }
    1803             : 
    1804             : long
    1805      206766 : F2xY_degreex(GEN b)
    1806             : {
    1807      206766 :   long deg = 0, i;
    1808      206766 :   if (!signe(b)) return -1;
    1809     1342068 :   for (i = 2; i < lg(b); ++i)
    1810     1135302 :     deg = maxss(deg, F2x_degree(gel(b, i)));
    1811      206766 :   return deg;
    1812             : }
    1813             : 
    1814             : GEN
    1815         581 : FlxX_to_F2xX(GEN B)
    1816             : {
    1817         581 :   long lb=lg(B);
    1818             :   long i;
    1819         581 :   GEN b=cgetg(lb,t_POL);
    1820         581 :   b[1]=evalsigne(1)|(((ulong)B[1])&VARNBITS);
    1821        3143 :   for (i=2; i<lb; i++)
    1822        2562 :     gel(b,i) = Flx_to_F2x(gel(B,i));
    1823         581 :   return F2xX_renormalize(b, lb);
    1824             : }
    1825             : 
    1826             : GEN
    1827        4550 : ZXX_to_F2xX(GEN B, long v)
    1828             : {
    1829        4550 :   long lb=lg(B);
    1830             :   long i;
    1831        4550 :   GEN b=cgetg(lb,t_POL);
    1832        4550 :   b[1]=evalsigne(1)|(((ulong)B[1])&VARNBITS);
    1833       21434 :   for (i=2; i<lb; i++)
    1834       16884 :     switch (typ(gel(B,i)))
    1835             :     {
    1836        6888 :     case t_INT:
    1837        6888 :       gel(b,i) = Z_to_F2x(gel(B,i), evalvarn(v));
    1838        6888 :       break;
    1839        9996 :     case t_POL:
    1840        9996 :       gel(b,i) = ZX_to_F2x(gel(B,i));
    1841        9996 :       break;
    1842             :     }
    1843       21434 :   return F2xX_renormalize(b, lb);
    1844             : }
    1845             : 
    1846             : GEN
    1847          56 : F2xX_to_FlxX(GEN B)
    1848             : {
    1849          56 :   long i, lb = lg(B);
    1850          56 :   GEN b = cgetg(lb,t_POL);
    1851         266 :   for (i=2; i<lb; i++)
    1852         210 :     gel(b,i) = F2x_to_Flx(gel(B,i));
    1853          56 :   b[1] = B[1]; return b;
    1854             : }
    1855             : 
    1856             : GEN
    1857         966 : F2xX_to_ZXX(GEN B)
    1858             : {
    1859         966 :   long i, lb = lg(B);
    1860         966 :   GEN b = cgetg(lb,t_POL);
    1861        2975 :   for (i=2; i<lb; i++)
    1862             :   {
    1863        2009 :     GEN c = gel(B,i);
    1864        2009 :     gel(b,i) = lgpol(c) ?  F2x_equal1(c) ? gen_1 : F2x_to_ZX(c) : gen_0;
    1865             :   }
    1866         966 :   b[1] = B[1]; return b;
    1867             : }
    1868             : 
    1869             : GEN
    1870       80731 : F2xX_deriv(GEN z)
    1871             : {
    1872       80731 :   long i,l = lg(z)-1;
    1873             :   GEN x;
    1874       80731 :   if (l < 2) l = 2;
    1875       80731 :   x = cgetg(l, t_POL); x[1] = z[1];
    1876      528696 :   for (i=2; i<l; i++) gel(x,i) = odd(i) ? pol0_F2x(mael(z,i+1,1)): gel(z,i+1);
    1877       80731 :   return F2xX_renormalize(x,l);
    1878             : }
    1879             : 
    1880             : static GEN
    1881        1359 : F2xX_addspec(GEN x, GEN y, long lx, long ly)
    1882             : {
    1883             :   long i,lz;
    1884             :   GEN z;
    1885        1359 :   if (ly>lx) swapspec(x,y, lx,ly);
    1886        1359 :   lz = lx+2; z = cgetg(lz, t_POL);
    1887      144699 :   for (i=0; i<ly; i++) gel(z,i+2) = F2x_add(gel(x,i), gel(y,i));
    1888        1359 :   for (   ; i<lx; i++) gel(z,i+2) = F2x_copy(gel(x,i));
    1889        1359 :   z[1] = 0; return F2xX_renormalize(z, lz);
    1890             : }
    1891             : 
    1892             : GEN
    1893      422607 : F2xX_add(GEN x, GEN y)
    1894             : {
    1895             :   long i,lz;
    1896             :   GEN z;
    1897      422607 :   long lx=lg(x);
    1898      422607 :   long ly=lg(y);
    1899      422607 :   if (ly>lx) swapspec(x,y, lx,ly);
    1900      422607 :   lz = lx; z = cgetg(lz, t_POL); z[1]=x[1];
    1901     1485462 :   for (i=2; i<ly; i++) gel(z,i) = F2x_add(gel(x,i), gel(y,i));
    1902     1677597 :   for (   ; i<lx; i++) gel(z,i) = F2x_copy(gel(x,i));
    1903      422607 :   return F2xX_renormalize(z, lz);
    1904             : }
    1905             : 
    1906             : GEN
    1907           0 : F2xX_F2x_add(GEN x, GEN y)
    1908             : {
    1909           0 :   long i, lz = lg(y);
    1910             :   GEN z;
    1911           0 :   if (signe(y) == 0) return scalarpol(x, varn(y));
    1912           0 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1913           0 :   gel(z,2) = F2x_add(gel(y,2), x);
    1914           0 :   if (lz == 3) z = F2xX_renormalize(z,lz);
    1915             :   else
    1916           0 :     for(i=3;i<lz;i++) gel(z,i) = F2x_copy(gel(y,i));
    1917           0 :   return z;
    1918             : }
    1919             : 
    1920             : GEN
    1921      480696 : F2xX_F2x_mul(GEN P, GEN U)
    1922             : {
    1923      480696 :   long i, lP = lg(P);
    1924      480696 :   GEN res = cgetg(lP,t_POL);
    1925      480696 :   res[1] = P[1];
    1926     2326459 :   for(i=2; i<lP; i++)
    1927     1845763 :     gel(res,i) = F2x_mul(U,gel(P,i));
    1928      480696 :   return F2xX_renormalize(res, lP);
    1929             : }
    1930             : 
    1931             : GEN
    1932      136374 : F2xY_F2xqV_evalx(GEN P, GEN x, GEN T)
    1933             : {
    1934      136374 :   long i, lP = lg(P);
    1935      136374 :   GEN res = cgetg(lP,t_POL);
    1936      136374 :   res[1] = P[1];
    1937      616819 :   for(i=2; i<lP; i++)
    1938      480445 :     gel(res,i) = F2x_F2xqV_eval(gel(P,i), x, T);
    1939      136374 :   return F2xX_renormalize(res, lP);
    1940             : }
    1941             : 
    1942             : GEN
    1943           0 : F2xY_F2xq_evalx(GEN P, GEN x, GEN T)
    1944             : {
    1945           0 :   pari_sp av = avma;
    1946           0 :   long n = brent_kung_optpow(get_F2x_degree(T)-1,lgpol(P),1);
    1947           0 :   GEN xp = F2xq_powers(x, n, T);
    1948           0 :   return gerepileupto(av, F2xY_F2xqV_evalx(P, xp, T));
    1949             : }
    1950             : 
    1951             : static GEN
    1952        6678 : F2xX_to_Kronecker_spec(GEN P, long n, long d)
    1953             : {
    1954        6678 :   long i, k, l, N = 2*d + 1;
    1955        6678 :   long dP = n+1;
    1956             :   GEN x;
    1957        6678 :   if (n == 0) return pol0_Flx(P[1]&VARNBITS);
    1958        6678 :   l = nbits2nlong(N*dP+d+1);
    1959        6678 :   x = zero_zv(l+1);
    1960      382369 :   for (k=i=0; i<n; i++, k+=N)
    1961             :   {
    1962      375691 :     GEN c = gel(P,i);
    1963      375691 :     F2x_addshiftip(x, c, k);
    1964             :   }
    1965        6678 :   x[1] = P[1]&VARNBITS; return F2x_renormalize(x, l+2);
    1966             : }
    1967             : 
    1968             : GEN
    1969      308188 : F2xX_to_Kronecker(GEN P, long d)
    1970             : {
    1971      308188 :   long i, k, l, N = 2*d + 1;
    1972      308188 :   long dP = degpol(P);
    1973             :   GEN x;
    1974      308188 :   if (dP < 0) return pol0_Flx(P[1]&VARNBITS);
    1975      306880 :   l = nbits2nlong(N*dP+d+1);
    1976      306880 :   x = zero_zv(l+1);
    1977     1689073 :   for (k=i=0; i<=dP; i++, k+=N)
    1978             :   {
    1979     1382193 :     GEN c = gel(P,i+2);
    1980     1382193 :     F2x_addshiftip(x, c, k);
    1981             :   }
    1982      306880 :   x[1] = P[1]&VARNBITS; return F2x_renormalize(x, l+2);
    1983             : }
    1984             : 
    1985             : static GEN
    1986     1597814 : F2x_slice(GEN x, long n, long d)
    1987             : {
    1988     1597814 :   ulong ib, il=dvmduBIL(n, &ib);
    1989     1597814 :   ulong db, dl=dvmduBIL(d, &db);
    1990     1597814 :   long lN = dl+2+(db?1:0);
    1991     1597814 :   GEN t = cgetg(lN,t_VECSMALL);
    1992     1597814 :   t[1] = x[1];
    1993     1597814 :   if (ib)
    1994             :   {
    1995     1433025 :     ulong i, ic = BITS_IN_LONG-ib;
    1996     1433025 :     ulong r = uel(x,2+il)>>ib;
    1997     1433233 :     for(i=0; i<dl; i++)
    1998             :     {
    1999         208 :       uel(t,2+i) = (uel(x,3+il+i)<<ic)|r;
    2000         208 :       r = uel(x,3+il+i)>>ib;
    2001             :     }
    2002     1433025 :     if (db)
    2003     1433025 :       uel(t,2+i) = (uel(x,3+il+i)<<ic)|r;
    2004             :   }
    2005             :   else
    2006             :   {
    2007             :     long i;
    2008      329606 :     for(i=2; i<lN; i++)
    2009      164817 :       uel(t,i) = uel(x,il+i);
    2010             :   }
    2011     1597814 :   if (db) uel(t,lN-1) &= (1UL<<db)-1;
    2012     1597814 :   return F2x_renormalize(t, lN);
    2013             : }
    2014             : 
    2015             : GEN
    2016      157433 : Kronecker_to_F2xqX(GEN z, GEN T)
    2017             : {
    2018      157433 :   long lz = F2x_degree(z)+1;
    2019      157433 :   long i, j, N = get_F2x_degree(T)*2 + 1;
    2020      157433 :   long lx = lz / (N-2);
    2021      157433 :   GEN x = cgetg(lx+3,t_POL);
    2022      157433 :   x[1] = z[1];
    2023     1755247 :   for (i=0, j=2; i<lz; i+=N, j++)
    2024             :   {
    2025     1597814 :     GEN t = F2x_slice(z,i,minss(lz-i,N));
    2026     1597814 :     t[1] = T[1];
    2027     1597814 :     gel(x,j) = F2x_rem(t, T);
    2028             :   }
    2029      157433 :   return F2xX_renormalize(x, j);
    2030             : }
    2031             : 
    2032             : GEN
    2033           0 : F2xX_to_F2xC(GEN x, long N, long sv)
    2034             : {
    2035             :   long i, l;
    2036             :   GEN z;
    2037           0 :   l = lg(x)-1; x++;
    2038           0 :   if (l > N+1) l = N+1; /* truncate higher degree terms */
    2039           0 :   z = cgetg(N+1,t_COL);
    2040           0 :   for (i=1; i<l ; i++) gel(z,i) = gel(x,i);
    2041           0 :   for (   ; i<=N; i++) gel(z,i) = pol0_F2x(sv);
    2042           0 :   return z;
    2043             : }
    2044             : 
    2045             : GEN
    2046           0 : F2xXV_to_F2xM(GEN v, long n, long sv)
    2047             : {
    2048           0 :   long j, N = lg(v);
    2049           0 :   GEN y = cgetg(N, t_MAT);
    2050           0 :   for (j=1; j<N; j++) gel(y,j) = F2xX_to_F2xC(gel(v,j), n, sv);
    2051           0 :   return y;
    2052             : }
    2053             : 
    2054             : /***********************************************************************/
    2055             : /**                                                                   **/
    2056             : /**                             F2xXV/F2xXC                           **/
    2057             : /**                                                                   **/
    2058             : /***********************************************************************/
    2059             : 
    2060             : GEN
    2061         805 : FlxXC_to_F2xXC(GEN x) { pari_APPLY_type(t_COL, FlxX_to_F2xX(gel(x,i))); }
    2062             : GEN
    2063        1799 : F2xXC_to_ZXXC(GEN x) { pari_APPLY_type(t_COL, F2xX_to_ZXX(gel(x,i))); }
    2064             : 
    2065             : /***********************************************************************/
    2066             : /**                                                                   **/
    2067             : /**                             F2xqX                                 **/
    2068             : /**                                                                   **/
    2069             : /***********************************************************************/
    2070             : 
    2071             : static GEN
    2072      784474 : get_F2xqX_red(GEN T, GEN *B)
    2073             : {
    2074      784474 :   if (typ(T)!=t_VEC) { *B=NULL; return T; }
    2075        3037 :   *B = gel(T,1); return gel(T,2);
    2076             : }
    2077             : 
    2078             : GEN
    2079        2783 : random_F2xqX(long d1, long v, GEN T)
    2080             : {
    2081        2783 :   long dT = get_F2x_degree(T), vT = get_F2x_var(T);
    2082        2783 :   long i, d = d1+2;
    2083        2783 :   GEN y = cgetg(d,t_POL); y[1] = evalsigne(1) | evalvarn(v);
    2084       18759 :   for (i=2; i<d; i++) gel(y,i) = random_F2x(dT, vT);
    2085        2783 :   return FlxX_renormalize(y,d);
    2086             : }
    2087             : 
    2088             : GEN
    2089     1006326 : F2xqX_red(GEN z, GEN T)
    2090             : {
    2091             :   GEN res;
    2092     1006326 :   long i, l = lg(z);
    2093     1006326 :   res = cgetg(l,t_POL); res[1] = z[1];
    2094     5679889 :   for(i=2;i<l;i++) gel(res,i) = F2x_rem(gel(z,i),T);
    2095     1006326 :   return F2xX_renormalize(res,l);
    2096             : }
    2097             : 
    2098             : GEN
    2099        7021 : F2xqX_F2xq_mul(GEN P, GEN U, GEN T)
    2100             : {
    2101        7021 :   long i, lP = lg(P);
    2102        7021 :   GEN res = cgetg(lP,t_POL);
    2103        7021 :   res[1] = P[1];
    2104       39249 :   for(i=2; i<lP; i++)
    2105       32228 :     gel(res,i) = F2xq_mul(U,gel(P,i), T);
    2106        7021 :   return F2xX_renormalize(res, lP);
    2107             : }
    2108             : 
    2109             : GEN
    2110      206185 : F2xqX_F2xq_mul_to_monic(GEN P, GEN U, GEN T)
    2111             : {
    2112      206185 :   long i, lP = lg(P);
    2113      206185 :   GEN res = cgetg(lP,t_POL);
    2114      206185 :   res[1] = P[1];
    2115     1130451 :   for(i=2; i<lP-1; i++) gel(res,i) = F2xq_mul(U,gel(P,i), T);
    2116      206185 :   gel(res,lP-1) = pol1_F2x(T[1]);
    2117      206185 :   return F2xX_renormalize(res, lP);
    2118             : }
    2119             : 
    2120             : GEN
    2121      206199 : F2xqX_normalize(GEN z, GEN T)
    2122             : {
    2123      206199 :   GEN p1 = leading_coeff(z);
    2124      206199 :   if (!lgpol(z) || (!degpol(p1) && p1[1] == 1)) return z;
    2125      206185 :   return F2xqX_F2xq_mul_to_monic(z, F2xq_inv(p1,T), T);
    2126             : }
    2127             : 
    2128             : GEN
    2129      154094 : F2xqX_mul(GEN x, GEN y, GEN T)
    2130             : {
    2131      154094 :   pari_sp ltop=avma;
    2132             :   GEN z, kx, ky;
    2133      154094 :   long d = get_F2x_degree(T);
    2134      154094 :   kx= F2xX_to_Kronecker(x, d);
    2135      154094 :   ky= F2xX_to_Kronecker(y, d);
    2136      154094 :   z = F2x_mul(ky, kx);
    2137      154094 :   z = Kronecker_to_F2xqX(z, T);
    2138      154094 :   return gerepileupto(ltop, z);
    2139             : }
    2140             : 
    2141             : static GEN
    2142        3339 : F2xqX_mulspec(GEN x, GEN y, GEN T, long nx, long ny)
    2143             : {
    2144        3339 :   pari_sp ltop=avma;
    2145             :   GEN z, kx, ky;
    2146        3339 :   long dT = get_F2x_degree(T);
    2147        3339 :   kx= F2xX_to_Kronecker_spec(x, nx, dT);
    2148        3339 :   ky= F2xX_to_Kronecker_spec(y, ny, dT);
    2149        3339 :   z = F2x_mul(ky, kx);
    2150        3339 :   z = Kronecker_to_F2xqX(z,T);
    2151        3339 :   return gerepileupto(ltop,z);
    2152             : }
    2153             : 
    2154             : GEN
    2155      105298 : F2xqX_sqr(GEN x, GEN T)
    2156             : {
    2157      105298 :   long d = degpol(x), l = 2*d+3;
    2158             :   GEN z;
    2159      105298 :   if (!signe(x)) return pol_0(varn(x));
    2160      105291 :   z = cgetg(l,t_POL);
    2161      105291 :   z[1] = x[1];
    2162      105291 :   if (d > 0)
    2163             :   {
    2164      105186 :     GEN u = pol0_F2x(T[1]);
    2165             :     long i,j;
    2166      269831 :     for (i=2,j=2; i<d+2; i++,j+=2)
    2167             :     {
    2168      164645 :       gel(z, j) = F2xq_sqr(gel(x, i), T);
    2169      164645 :       gel(z, j+1) = u;
    2170             :     }
    2171             :   }
    2172      105291 :   gel(z, 2+2*d) = F2xq_sqr(gel(x, 2+d), T);
    2173      105291 :   return FlxX_renormalize(z, l);
    2174             : }
    2175             : 
    2176             : static GEN
    2177      545878 : F2xqX_divrem_basecase(GEN x, GEN y, GEN T, GEN *pr)
    2178             : {
    2179             :   long vx, dx, dy, dz, i, j, sx, lr;
    2180             :   pari_sp av0, av, tetpil;
    2181             :   GEN z,p1,rem,lead;
    2182             : 
    2183      545878 :   if (!signe(y)) pari_err_INV("F2xqX_divrem",y);
    2184      545878 :   vx=varn(x); dy=degpol(y); dx=degpol(x);
    2185      545878 :   if (dx < dy)
    2186             :   {
    2187           0 :     if (pr)
    2188             :     {
    2189           0 :       av0 = avma; x = F2xqX_red(x, T);
    2190           0 :       if (pr == ONLY_DIVIDES) { set_avma(av0); return signe(x)? NULL: pol_0(vx); }
    2191           0 :       if (pr == ONLY_REM) return x;
    2192           0 :       *pr = x;
    2193             :     }
    2194           0 :     return pol_0(vx);
    2195             :   }
    2196      545878 :   lead = leading_coeff(y);
    2197      545878 :   if (!dy) /* y is constant */
    2198             :   {
    2199      115679 :     if (pr && pr != ONLY_DIVIDES)
    2200             :     {
    2201      105606 :       if (pr == ONLY_REM) return pol_0(vx);
    2202          14 :       *pr = pol_0(vx);
    2203             :     }
    2204       10087 :     if (F2x_equal1(lead)) return gcopy(x);
    2205        7014 :     av0 = avma; x = F2xqX_F2xq_mul(x,F2xq_inv(lead,T),T);
    2206        7014 :     return gerepileupto(av0,x);
    2207             :   }
    2208      430199 :   av0 = avma; dz = dx-dy;
    2209      430199 :   lead = F2x_equal1(lead)? NULL: gclone(F2xq_inv(lead,T));
    2210      430199 :   set_avma(av0);
    2211      430199 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    2212      430199 :   x += 2; y += 2; z += 2;
    2213             : 
    2214      430199 :   p1 = gel(x,dx); av = avma;
    2215      430199 :   gel(z,dz) = lead? gerepileupto(av, F2xq_mul(p1,lead, T)): leafcopy(p1);
    2216     1277002 :   for (i=dx-1; i>=dy; i--)
    2217             :   {
    2218      846803 :     av=avma; p1=gel(x,i);
    2219     2598929 :     for (j=i-dy+1; j<=i && j<=dz; j++)
    2220     1752126 :       p1 = F2x_add(p1, F2x_mul(gel(z,j),gel(y,i-j)));
    2221      846803 :     if (lead) p1 = F2x_mul(p1, lead);
    2222      846803 :     tetpil=avma; gel(z,i-dy) = gerepile(av,tetpil,F2x_rem(p1,T));
    2223             :   }
    2224      430199 :   if (!pr) { guncloneNULL(lead); return z-2; }
    2225             : 
    2226      406941 :   rem = (GEN)avma; av = (pari_sp)new_chunk(dx+3);
    2227      406941 :   for (sx=0; ; i--)
    2228             :   {
    2229      222812 :     p1 = gel(x,i);
    2230     1874262 :     for (j=0; j<=i && j<=dz; j++)
    2231     1244509 :       p1 = F2x_add(p1, F2x_mul(gel(z,j),gel(y,i-j)));
    2232      629753 :     tetpil=avma; p1 = F2x_rem(p1, T); if (lgpol(p1)) { sx = 1; break; }
    2233      285994 :     if (!i) break;
    2234      222812 :     set_avma(av);
    2235             :   }
    2236      406941 :   if (pr == ONLY_DIVIDES)
    2237             :   {
    2238           0 :     guncloneNULL(lead);
    2239           0 :     if (sx) return gc_NULL(av0);
    2240           0 :     set_avma((pari_sp)rem); return z-2;
    2241             :   }
    2242      406941 :   lr=i+3; rem -= lr;
    2243      406941 :   rem[0] = evaltyp(t_POL) | evallg(lr);
    2244      406941 :   rem[1] = z[-1];
    2245      406941 :   p1 = gerepile((pari_sp)rem,tetpil,p1);
    2246      406941 :   rem += 2; gel(rem,i) = p1;
    2247     1384921 :   for (i--; i>=0; i--)
    2248             :   {
    2249      977980 :     av=avma; p1 = gel(x,i);
    2250     3368356 :     for (j=0; j<=i && j<=dz; j++)
    2251     2390376 :       p1 = F2x_add(p1, F2x_mul(gel(z,j),gel(y,i-j)));
    2252      977980 :     tetpil=avma; gel(rem,i) = gerepile(av,tetpil, F2x_rem(p1, T));
    2253             :   }
    2254      406941 :   rem -= 2;
    2255      406941 :   guncloneNULL(lead);
    2256      406941 :   if (!sx) (void)F2xX_renormalize(rem, lr);
    2257      406941 :   if (pr == ONLY_REM) return gerepileupto(av0,rem);
    2258          35 :   *pr = rem; return z-2;
    2259             : }
    2260             : 
    2261             : static GEN
    2262        2768 : F2xX_recipspec(GEN x, long l, long n, long vs)
    2263             : {
    2264             :   long i;
    2265        2768 :   GEN z = cgetg(n+2,t_POL);
    2266        2768 :   z[1] = 0; z += 2;
    2267      146043 :   for(i=0; i<l; i++)
    2268      143275 :     gel(z,n-i-1) = F2x_copy(gel(x,i));
    2269        2768 :   for(   ; i<n; i++)
    2270           0 :     gel(z,n-i-1) = pol0_F2x(vs);
    2271        2768 :   return F2xX_renormalize(z-2,n+2);
    2272             : }
    2273             : 
    2274             : static GEN
    2275           6 : F2xqX_invBarrett_basecase(GEN T, GEN Q)
    2276             : {
    2277           6 :   long i, l=lg(T)-1, lr = l-1, k;
    2278           6 :   long sv=Q[1];
    2279           6 :   GEN r=cgetg(lr,t_POL); r[1]=T[1];
    2280           6 :   gel(r,2) = pol1_F2x(sv);
    2281         273 :   for (i=3;i<lr;i++)
    2282             :   {
    2283         267 :     pari_sp ltop=avma;
    2284         267 :     GEN u = gel(T,l-i+2);
    2285        6075 :     for (k=3;k<i;k++)
    2286        5808 :       u = F2x_add(u, F2xq_mul(gel(T,l-i+k),gel(r,k),Q));
    2287         267 :     gel(r,i) = gerepileupto(ltop, u);
    2288             :   }
    2289           6 :   r = F2xX_renormalize(r,lr);
    2290           6 :   return r;
    2291             : }
    2292             : 
    2293             : /* Return new lgpol */
    2294             : static long
    2295        3651 : F2xX_lgrenormalizespec(GEN x, long lx)
    2296             : {
    2297             :   long i;
    2298        3651 :   for (i = lx-1; i>=0; i--)
    2299        3651 :     if (lgpol(gel(x,i))) break;
    2300        3651 :   return i+1;
    2301             : }
    2302             : 
    2303             : static GEN
    2304          44 : F2xqX_invBarrett_Newton(GEN S, GEN T)
    2305             : {
    2306          44 :   pari_sp av = avma;
    2307          44 :   long nold, lx, lz, lq, l = degpol(S), i, lQ;
    2308          44 :   GEN q, y, z, x = cgetg(l+2, t_POL) + 2;
    2309          44 :   long dT = get_F2x_degree(T);
    2310          44 :   ulong mask = quadratic_prec_mask(l-2); /* assume l > 2 */
    2311        5251 :   for (i=0;i<l;i++) gel(x,i) = pol0_F2x(T[1]);
    2312          44 :   q = F2xX_recipspec(S+2,l+1,l+1,dT);
    2313          44 :   lQ = lgpol(q); q+=2;
    2314             :   /* We work on _spec_ FlxX's, all the l[xzq] below are lgpol's */
    2315             : 
    2316             :   /* initialize */
    2317          44 :   gel(x,0) = F2xq_inv(gel(q,0),T);
    2318          44 :   if (lQ>1 && F2x_degree(gel(q,1)) >= dT)
    2319           0 :     gel(q,1) = F2x_rem(gel(q,1), T);
    2320          44 :   if (lQ>1 && lgpol(gel(q,1)))
    2321          44 :   {
    2322          44 :     GEN u = gel(q, 1);
    2323          44 :     if (!F2x_equal1(gel(x,0))) u = F2xq_mul(u, F2xq_sqr(gel(x,0), T), T);
    2324          44 :     else u = F2x_copy(u);
    2325          44 :     gel(x,1) = u; lx = 2;
    2326             :   }
    2327             :   else
    2328           0 :     lx = 1;
    2329          44 :   nold = 1;
    2330         353 :   for (; mask > 1; )
    2331             :   { /* set x -= x(x*q - 1) + O(t^(nnew + 1)), knowing x*q = 1 + O(t^(nold+1)) */
    2332         309 :     long i, lnew, nnew = nold << 1;
    2333             : 
    2334         309 :     if (mask & 1) nnew--;
    2335         309 :     mask >>= 1;
    2336             : 
    2337         309 :     lnew = nnew + 1;
    2338         309 :     lq = F2xX_lgrenormalizespec(q, minss(lQ,lnew));
    2339         309 :     z = F2xqX_mulspec(x, q, T, lx, lq); /* FIXME: high product */
    2340         309 :     lz = lgpol(z); if (lz > lnew) lz = lnew;
    2341         309 :     z += 2;
    2342             :     /* subtract 1 [=>first nold words are 0]: renormalize so that z(0) != 0 */
    2343         618 :     for (i = nold; i < lz; i++) if (lgpol(gel(z,i))) break;
    2344         309 :     nold = nnew;
    2345         309 :     if (i >= lz) continue; /* z-1 = 0(t^(nnew + 1)) */
    2346             : 
    2347             :     /* z + i represents (x*q - 1) / t^i */
    2348         309 :     lz = F2xX_lgrenormalizespec (z+i, lz-i);
    2349         309 :     z = F2xqX_mulspec(x, z+i, T, lx, lz); /* FIXME: low product */
    2350         309 :     lz = lgpol(z); z += 2;
    2351         309 :     if (lz > lnew-i) lz = F2xX_lgrenormalizespec(z, lnew-i);
    2352             : 
    2353         309 :     lx = lz+ i;
    2354         309 :     y  = x + i; /* x -= z * t^i, in place */
    2355        5384 :     for (i = 0; i < lz; i++) gel(y,i) = gel(z,i);
    2356             :   }
    2357          44 :   x -= 2; setlg(x, lx + 2); x[1] = S[1];
    2358          44 :   return gerepilecopy(av, x);
    2359             : }
    2360             : 
    2361             : GEN
    2362          50 : F2xqX_invBarrett(GEN T, GEN Q)
    2363             : {
    2364          50 :   pari_sp ltop=avma;
    2365          50 :   long l=lg(T), v = varn(T);
    2366             :   GEN r;
    2367          50 :   GEN c = gel(T,l-1);
    2368          50 :   if (l<5) return pol_0(v);
    2369          50 :   if (l<=F2xqX_INVBARRETT_LIMIT)
    2370             :   {
    2371           6 :     if (!F2x_equal1(c))
    2372             :     {
    2373           0 :       GEN ci = F2xq_inv(c,Q);
    2374           0 :       T = F2xqX_F2xq_mul(T, ci, Q);
    2375           0 :       r = F2xqX_invBarrett_basecase(T,Q);
    2376           0 :       r = F2xqX_F2xq_mul(r,ci,Q);
    2377             :     } else
    2378           6 :       r = F2xqX_invBarrett_basecase(T,Q);
    2379             :   } else
    2380          44 :     r = F2xqX_invBarrett_Newton(T,Q);
    2381          50 :   return gerepileupto(ltop, r);
    2382             : }
    2383             : 
    2384             : GEN
    2385      219888 : F2xqX_get_red(GEN S, GEN T)
    2386             : {
    2387      219888 :   if (typ(S)==t_POL && lg(S)>F2xqX_BARRETT_LIMIT)
    2388          47 :     retmkvec2(F2xqX_invBarrett(S, T), S);
    2389      219841 :   return S;
    2390             : }
    2391             : 
    2392             : /* Compute x mod S where 2 <= degpol(S) <= l+1 <= 2*(degpol(S)-1)
    2393             :  *  * and mg is the Barrett inverse of S. */
    2394             : static GEN
    2395        1362 : F2xqX_divrem_Barrettspec(GEN x, long l, GEN mg, GEN S, GEN T, GEN *pr)
    2396             : {
    2397             :   GEN q, r;
    2398        1362 :   long lt = degpol(S); /*We discard the leading term*/
    2399             :   long ld, lm, lT, lmg;
    2400        1362 :   ld = l-lt;
    2401        1362 :   lm = minss(ld, lgpol(mg));
    2402        1362 :   lT  = F2xX_lgrenormalizespec(S+2,lt);
    2403        1362 :   lmg = F2xX_lgrenormalizespec(mg+2,lm);
    2404        1362 :   q = F2xX_recipspec(x+lt,ld,ld,0);               /* q = rec(x)     lq<=ld*/
    2405        1362 :   q = F2xqX_mulspec(q+2,mg+2,T,lgpol(q),lmg);   /* q = rec(x) * mg lq<=ld+lm*/
    2406        1362 :   q = F2xX_recipspec(q+2,minss(ld,lgpol(q)),ld,0);/* q = rec (rec(x) * mg) lq<=ld*/
    2407        1362 :   if (!pr) return q;
    2408        1359 :   r = F2xqX_mulspec(q+2,S+2,T,lgpol(q),lT);     /* r = q*pol        lr<=ld+lt*/
    2409        1359 :   r = F2xX_addspec(x,r+2,lt,minss(lt,lgpol(r)));/* r = x - r   lr<=lt */
    2410        1359 :   if (pr == ONLY_REM) return r;
    2411           0 :   *pr = r; return q;
    2412             : }
    2413             : 
    2414             : static GEN
    2415        1362 : F2xqX_divrem_Barrett(GEN x, GEN mg, GEN S, GEN T, GEN *pr)
    2416             : {
    2417        1362 :   GEN q = NULL, r = F2xqX_red(x, T);
    2418        1362 :   long l = lgpol(r), lt = degpol(S), lm = 2*lt-1, v = varn(S);
    2419             :   long i;
    2420        1362 :   if (l <= lt)
    2421             :   {
    2422           0 :     if (pr == ONLY_REM) return r;
    2423           0 :     if (pr == ONLY_DIVIDES) return signe(r)? NULL: pol_0(v);
    2424           0 :     if (pr) *pr = r;
    2425           0 :     return pol_0(v);
    2426             :   }
    2427        1362 :   if (lt <= 1)
    2428           0 :     return F2xqX_divrem_basecase(x,S,T,pr);
    2429        1362 :   if (pr != ONLY_REM && l>lm)
    2430             :   {
    2431           0 :     long vT = get_F2x_var(T);
    2432           0 :     q = cgetg(l-lt+2, t_POL); q[1] = S[1];
    2433           0 :     for (i=0;i<l-lt;i++) gel(q+2,i) = pol0_F2x(vT);
    2434             :   }
    2435        1362 :   while (l>lm)
    2436             :   {
    2437           0 :     GEN zr, zq = F2xqX_divrem_Barrettspec(r+2+l-lm,lm,mg,S,T,&zr);
    2438           0 :     long lz = lgpol(zr);
    2439           0 :     if (pr != ONLY_REM)
    2440             :     {
    2441           0 :       long lq = lgpol(zq);
    2442           0 :       for(i=0; i<lq; i++) gel(q+2+l-lm,i) = gel(zq,2+i);
    2443             :     }
    2444           0 :     for(i=0; i<lz; i++) gel(r+2+l-lm,i) = gel(zr,2+i);
    2445           0 :     l = l-lm+lz;
    2446             :   }
    2447        1362 :   if (pr == ONLY_REM)
    2448             :   {
    2449        1359 :     if (l > lt)
    2450        1359 :       r = F2xqX_divrem_Barrettspec(r+2,l,mg,S,T,ONLY_REM);
    2451             :     else
    2452           0 :       r = F2xX_renormalize(r, lg(r));
    2453        1359 :     setvarn(r, v); return F2xX_renormalize(r, lg(r));
    2454             :   }
    2455           3 :   if (l > lt)
    2456             :   {
    2457           3 :     GEN zq = F2xqX_divrem_Barrettspec(r+2,l,mg,S,T,pr? &r: NULL);
    2458           3 :     if (!q) q = zq;
    2459             :     else
    2460             :     {
    2461           0 :       long lq = lgpol(zq);
    2462           0 :       for(i=0; i<lq; i++) gel(q+2,i) = gel(zq,2+i);
    2463             :     }
    2464             :   }
    2465           0 :   else if (pr)
    2466           0 :     r = F2xX_renormalize(r, l+2);
    2467           3 :   setvarn(q, v); q = F2xX_renormalize(q, lg(q));
    2468           3 :   if (pr == ONLY_DIVIDES) return signe(r)? NULL: q;
    2469           3 :   if (pr) { setvarn(r, v); *pr = r; }
    2470           3 :   return q;
    2471             : }
    2472             : 
    2473             : GEN
    2474       33383 : F2xqX_divrem(GEN x, GEN S, GEN T, GEN *pr)
    2475             : {
    2476             :   GEN B, y;
    2477             :   long dy, dx, d;
    2478       33383 :   if (pr==ONLY_REM) return F2xqX_rem(x, S, T);
    2479       33383 :   y = get_F2xqX_red(S, &B);
    2480       33383 :   dy = degpol(y); dx = degpol(x); d = dx-dy;
    2481       33383 :   if (!B && d+3 < F2xqX_DIVREM_BARRETT_LIMIT)
    2482       33380 :     return F2xqX_divrem_basecase(x,y,T,pr);
    2483             :   else
    2484             :   {
    2485           3 :     pari_sp av=avma;
    2486           3 :     GEN mg = B? B: F2xqX_invBarrett(y, T);
    2487           3 :     GEN q = F2xqX_divrem_Barrett(x,mg,y,T,pr);
    2488           3 :     if (!q) return gc_NULL(av);
    2489           3 :     if (!pr || pr==ONLY_DIVIDES) return gerepilecopy(av, q);
    2490           0 :     gerepileall(av,2,&q,pr);
    2491           0 :     return q;
    2492             :   }
    2493             : }
    2494             : 
    2495             : GEN
    2496      751091 : F2xqX_rem(GEN x, GEN S, GEN T)
    2497             : {
    2498      751091 :   GEN B, y = get_F2xqX_red(S, &B);
    2499      751091 :   long dy = degpol(y), dx = degpol(x), d = dx-dy;
    2500      751091 :   if (d < 0) return F2xqX_red(x, T);
    2501      513857 :   if (!B && d+3 < F2xqX_REM_BARRETT_LIMIT)
    2502      512498 :     return F2xqX_divrem_basecase(x,y, T, ONLY_REM);
    2503             :   else
    2504             :   {
    2505        1359 :     pari_sp av=avma;
    2506        1359 :     GEN mg = B? B: F2xqX_invBarrett(y, T);
    2507        1359 :     GEN r = F2xqX_divrem_Barrett(x, mg, y, T, ONLY_REM);
    2508        1359 :     return gerepileupto(av, r);
    2509             :   }
    2510             : }
    2511             : 
    2512             : static GEN
    2513           0 : F2xqX_halfgcd_basecase(GEN a, GEN b, GEN T)
    2514             : {
    2515           0 :   pari_sp av=avma;
    2516             :   GEN u,u1,v,v1;
    2517           0 :   long vx = varn(a);
    2518           0 :   long n = lgpol(a)>>1;
    2519           0 :   u1 = v = pol_0(vx);
    2520           0 :   u = v1 = pol1_F2xX(vx, get_F2x_var(T));
    2521           0 :   while (lgpol(b)>n)
    2522             :   {
    2523           0 :     GEN r, q = F2xqX_divrem(a,b, T, &r);
    2524           0 :     a = b; b = r; swap(u,u1); swap(v,v1);
    2525           0 :     u1 = F2xX_add(u1, F2xqX_mul(u, q, T));
    2526           0 :     v1 = F2xX_add(v1, F2xqX_mul(v, q ,T));
    2527           0 :     if (gc_needed(av,2))
    2528             :     {
    2529           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"F2xqX_halfgcd (d = %ld)",degpol(b));
    2530           0 :       gerepileall(av,6, &a,&b,&u1,&v1,&u,&v);
    2531             :     }
    2532             :   }
    2533           0 :   return gerepilecopy(av, mkmat2(mkcol2(u,u1), mkcol2(v,v1)));
    2534             : }
    2535             : static GEN
    2536           0 : F2xqX_addmulmul(GEN u, GEN v, GEN x, GEN y, GEN T)
    2537             : {
    2538           0 :   return F2xX_add(F2xqX_mul(u, x, T),F2xqX_mul(v, y, T));
    2539             : }
    2540             : 
    2541             : static GEN
    2542           0 : F2xqXM_F2xqX_mul2(GEN M, GEN x, GEN y, GEN T)
    2543             : {
    2544           0 :   GEN res = cgetg(3, t_COL);
    2545           0 :   gel(res, 1) = F2xqX_addmulmul(gcoeff(M,1,1), gcoeff(M,1,2), x, y, T);
    2546           0 :   gel(res, 2) = F2xqX_addmulmul(gcoeff(M,2,1), gcoeff(M,2,2), x, y, T);
    2547           0 :   return res;
    2548             : }
    2549             : 
    2550             : static GEN
    2551           0 : F2xqXM_mul2(GEN A, GEN B, GEN T)
    2552             : {
    2553           0 :   GEN A11=gcoeff(A,1,1),A12=gcoeff(A,1,2), B11=gcoeff(B,1,1),B12=gcoeff(B,1,2);
    2554           0 :   GEN A21=gcoeff(A,2,1),A22=gcoeff(A,2,2), B21=gcoeff(B,2,1),B22=gcoeff(B,2,2);
    2555           0 :   GEN M1 = F2xqX_mul(F2xX_add(A11,A22), F2xX_add(B11,B22), T);
    2556           0 :   GEN M2 = F2xqX_mul(F2xX_add(A21,A22), B11, T);
    2557           0 :   GEN M3 = F2xqX_mul(A11, F2xX_add(B12,B22), T);
    2558           0 :   GEN M4 = F2xqX_mul(A22, F2xX_add(B21,B11), T);
    2559           0 :   GEN M5 = F2xqX_mul(F2xX_add(A11,A12), B22, T);
    2560           0 :   GEN M6 = F2xqX_mul(F2xX_add(A21,A11), F2xX_add(B11,B12), T);
    2561           0 :   GEN M7 = F2xqX_mul(F2xX_add(A12,A22), F2xX_add(B21,B22), T);
    2562           0 :   GEN T1 = F2xX_add(M1,M4), T2 = F2xX_add(M7,M5);
    2563           0 :   GEN T3 = F2xX_add(M1,M2), T4 = F2xX_add(M3,M6);
    2564           0 :   retmkmat2(mkcol2(F2xX_add(T1,T2), F2xX_add(M2,M4)),
    2565             :             mkcol2(F2xX_add(M3,M5), F2xX_add(T3,T4)));
    2566             : }
    2567             : 
    2568             : /* Return [0,1;1,-q]*M */
    2569             : static GEN
    2570           0 : F2xqX_F2xqXM_qmul(GEN q, GEN M, GEN T)
    2571             : {
    2572           0 :   GEN u, v, res = cgetg(3, t_MAT);
    2573           0 :   u = F2xX_add(gcoeff(M,1,1), F2xqX_mul(gcoeff(M,2,1), q, T));
    2574           0 :   gel(res,1) = mkcol2(gcoeff(M,2,1), u);
    2575           0 :   v = F2xX_add(gcoeff(M,1,2), F2xqX_mul(gcoeff(M,2,2), q, T));
    2576           0 :   gel(res,2) = mkcol2(gcoeff(M,2,2), v);
    2577           0 :   return res;
    2578             : }
    2579             : 
    2580             : static GEN
    2581           0 : matid2_F2xXM(long v, long sv)
    2582             : {
    2583           0 :   retmkmat2(mkcol2(pol1_F2xX(v, sv),pol_0(v)),
    2584             :             mkcol2(pol_0(v),pol1_F2xX(v, sv)));
    2585             : }
    2586             : 
    2587             : static GEN
    2588           0 : F2xqX_halfgcd_split(GEN x, GEN y, GEN T)
    2589             : {
    2590           0 :   pari_sp av=avma;
    2591             :   GEN R, S, V;
    2592             :   GEN y1, r, q;
    2593           0 :   long l = lgpol(x), n = l>>1, k;
    2594           0 :   if (lgpol(y)<=n) return matid2_F2xXM(varn(x),get_F2x_var(T));
    2595           0 :   R = F2xqX_halfgcd(RgX_shift_shallow(x,-n),RgX_shift_shallow(y,-n), T);
    2596           0 :   V = F2xqXM_F2xqX_mul2(R,x,y, T); y1 = gel(V,2);
    2597           0 :   if (lgpol(y1)<=n) return gerepilecopy(av, R);
    2598           0 :   q = F2xqX_divrem(gel(V,1), y1, T, &r);
    2599           0 :   k = 2*n-degpol(y1);
    2600           0 :   S = F2xqX_halfgcd(RgX_shift_shallow(y1,-k), RgX_shift_shallow(r,-k), T);
    2601           0 :   return gerepileupto(av, F2xqXM_mul2(S,F2xqX_F2xqXM_qmul(q,R, T), T));
    2602             : }
    2603             : 
    2604             : /* Return M in GL_2(Fp[X]) such that:
    2605             : if [a',b']~=M*[a,b]~ then degpol(a')>= (lgpol(a)>>1) >degpol(b')
    2606             : */
    2607             : 
    2608             : static GEN
    2609           0 : F2xqX_halfgcd_i(GEN x, GEN y, GEN T)
    2610             : {
    2611           0 :   if (lg(x)<=F2xqX_HALFGCD_LIMIT) return F2xqX_halfgcd_basecase(x, y, T);
    2612           0 :   return F2xqX_halfgcd_split(x, y, T);
    2613             : }
    2614             : 
    2615             : GEN
    2616           0 : F2xqX_halfgcd(GEN x, GEN y, GEN T)
    2617             : {
    2618           0 :   pari_sp av = avma;
    2619             :   GEN M,q,r;
    2620           0 :   if (!signe(x))
    2621             :   {
    2622           0 :     long v = varn(x), vT = get_F2x_var(T);
    2623           0 :     retmkmat2(mkcol2(pol_0(v),pol1_F2xX(v,vT)),
    2624             :         mkcol2(pol1_F2xX(v,vT),pol_0(v)));
    2625             :   }
    2626           0 :   if (degpol(y)<degpol(x)) return F2xqX_halfgcd_i(x, y, T);
    2627           0 :   q = F2xqX_divrem(y, x, T, &r);
    2628           0 :   M = F2xqX_halfgcd_i(x, r, T);
    2629           0 :   gcoeff(M,1,1) = F2xX_add(gcoeff(M,1,1), F2xqX_mul(q, gcoeff(M,1,2), T));
    2630           0 :   gcoeff(M,2,1) = F2xX_add(gcoeff(M,2,1), F2xqX_mul(q, gcoeff(M,2,2), T));
    2631           0 :   return gerepilecopy(av, M);
    2632             : }
    2633             : 
    2634             : static GEN
    2635      169754 : F2xqX_gcd_basecase(GEN a, GEN b, GEN T)
    2636             : {
    2637      169754 :   pari_sp av = avma, av0=avma;
    2638      698941 :   while (signe(b))
    2639             :   {
    2640             :     GEN c;
    2641      529187 :     if (gc_needed(av0,2))
    2642             :     {
    2643           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"F2xqX_gcd (d = %ld)",degpol(b));
    2644           0 :       gerepileall(av0,2, &a,&b);
    2645             :     }
    2646      529187 :     av = avma; c = F2xqX_rem(a, b, T); a=b; b=c;
    2647             :   }
    2648      169754 :   set_avma(av); return a;
    2649             : }
    2650             : 
    2651             : GEN
    2652      170125 : F2xqX_gcd(GEN x, GEN y, GEN T)
    2653             : {
    2654      170125 :   pari_sp av = avma;
    2655      170125 :   x = F2xqX_red(x, T);
    2656      170125 :   y = F2xqX_red(y, T);
    2657      170125 :   if (!signe(x)) return gerepileupto(av, y);
    2658      169754 :   while (lg(y)>F2xqX_GCD_LIMIT)
    2659             :   {
    2660             :     GEN c;
    2661           0 :     if (lgpol(y)<=(lgpol(x)>>1))
    2662             :     {
    2663           0 :       GEN r = F2xqX_rem(x, y, T);
    2664           0 :       x = y; y = r;
    2665             :     }
    2666           0 :     c = F2xqXM_F2xqX_mul2(F2xqX_halfgcd(x,y, T), x, y, T);
    2667           0 :     x = gel(c,1); y = gel(c,2);
    2668           0 :     gerepileall(av,2,&x,&y);
    2669             :   }
    2670      169754 :   return gerepileupto(av, F2xqX_gcd_basecase(x, y, T));
    2671             : }
    2672             : 
    2673             : static GEN
    2674          14 : F2xqX_extgcd_basecase(GEN a, GEN b, GEN T,  GEN *ptu, GEN *ptv)
    2675             : {
    2676          14 :   pari_sp av=avma;
    2677             :   GEN u,v,d,d1,v1;
    2678          14 :   long vx = varn(a);
    2679          14 :   d = a; d1 = b;
    2680          14 :   v = pol_0(vx); v1 = pol1_F2xX(vx, get_F2x_var(T));
    2681          63 :   while (signe(d1))
    2682             :   {
    2683          49 :     GEN r, q = F2xqX_divrem(d, d1, T, &r);
    2684          49 :     v = F2xX_add(v,F2xqX_mul(q,v1,T));
    2685          49 :     u=v; v=v1; v1=u;
    2686          49 :     u=r; d=d1; d1=u;
    2687          49 :     if (gc_needed(av,2))
    2688             :     {
    2689           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"F2xqX_extgcd (d = %ld)",degpol(d));
    2690           0 :       gerepileall(av,5, &d,&d1,&u,&v,&v1);
    2691             :     }
    2692             :   }
    2693          14 :   if (ptu) *ptu = F2xqX_div(F2xX_add(d,F2xqX_mul(b,v, T)), a, T);
    2694          14 :   *ptv = v; return d;
    2695             : }
    2696             : 
    2697             : static GEN
    2698           0 : F2xqX_extgcd_halfgcd(GEN x, GEN y, GEN T,  GEN *ptu, GEN *ptv)
    2699             : {
    2700           0 :   pari_sp av=avma;
    2701           0 :   GEN u,v,R = matid2_F2xXM(varn(x), get_F2x_var(T));
    2702           0 :   while (lg(y)>F2xqX_EXTGCD_LIMIT)
    2703             :   {
    2704             :     GEN M, c;
    2705           0 :     if (lgpol(y)<=(lgpol(x)>>1))
    2706             :     {
    2707           0 :       GEN r, q = F2xqX_divrem(x, y, T, &r);
    2708           0 :       x = y; y = r;
    2709           0 :       R = F2xqX_F2xqXM_qmul(q, R, T);
    2710             :     }
    2711           0 :     M = F2xqX_halfgcd(x,y, T);
    2712           0 :     c = F2xqXM_F2xqX_mul2(M, x,y, T);
    2713           0 :     R = F2xqXM_mul2(M, R, T);
    2714           0 :     x = gel(c,1); y = gel(c,2);
    2715           0 :     gerepileall(av,3,&x,&y,&R);
    2716             :   }
    2717           0 :   y = F2xqX_extgcd_basecase(x,y, T, &u,&v);
    2718           0 :   if (ptu) *ptu = F2xqX_addmulmul(u,v,gcoeff(R,1,1),gcoeff(R,2,1), T);
    2719           0 :   *ptv = F2xqX_addmulmul(u,v,gcoeff(R,1,2),gcoeff(R,2,2), T);
    2720           0 :   return y;
    2721             : }
    2722             : 
    2723             : /* x and y in Z[Y][X], return lift(gcd(x mod T,p, y mod T,p)). Set u and v st
    2724             :  * ux + vy = gcd (mod T,p) */
    2725             : GEN
    2726          14 : F2xqX_extgcd(GEN x, GEN y, GEN T,  GEN *ptu, GEN *ptv)
    2727             : {
    2728             :   GEN d;
    2729          14 :   pari_sp ltop=avma;
    2730          14 :   x = F2xqX_red(x, T);
    2731          14 :   y = F2xqX_red(y, T);
    2732          14 :   if (lg(y)>F2xqX_EXTGCD_LIMIT)
    2733           0 :     d = F2xqX_extgcd_halfgcd(x, y, T, ptu, ptv);
    2734             :   else
    2735          14 :     d = F2xqX_extgcd_basecase(x, y, T, ptu, ptv);
    2736          14 :   gerepileall(ltop,ptu?3:2,&d,ptv,ptu);
    2737          14 :   return d;
    2738             : }
    2739             : 
    2740             : static GEN
    2741           0 : _F2xqX_mul(void *data,GEN a,GEN b) { return F2xqX_mul(a,b, (GEN) data); }
    2742             : static GEN
    2743           0 : _F2xqX_sqr(void *data,GEN a) { return F2xqX_sqr(a, (GEN) data); }
    2744             : GEN
    2745           0 : F2xqX_powu(GEN x, ulong n, GEN T)
    2746           0 : { return gen_powu(x, n, (void*)T, &_F2xqX_sqr, &_F2xqX_mul); }
    2747             : 
    2748             : /* Res(A,B) = Res(B,R) * lc(B)^(a-r) * (-1)^(ab), with R=A%B, a=deg(A) ...*/
    2749             : GEN
    2750          21 : F2xqX_resultant(GEN a, GEN b, GEN T)
    2751             : {
    2752          21 :   long vT = get_F2x_var(T);
    2753             :   long da,db,dc;
    2754             :   pari_sp av;
    2755          21 :   GEN c,lb, res = pol1_F2x(vT);
    2756             : 
    2757          21 :   if (!signe(a) || !signe(b)) return pol0_F2x(vT);
    2758             : 
    2759          21 :   da = degpol(a);
    2760          21 :   db = degpol(b);
    2761          21 :   if (db > da)
    2762           7 :     swapspec(a,b, da,db);
    2763          21 :   if (!da) return pol1_F2x(vT); /* = res * a[2] ^ db, since 0 <= db <= da = 0 */
    2764          21 :   av = avma;
    2765          49 :   while (db)
    2766             :   {
    2767          28 :     lb = gel(b,db+2);
    2768          28 :     c = F2xqX_rem(a,b, T);
    2769          28 :     a = b; b = c; dc = degpol(c);
    2770          28 :     if (dc < 0) { set_avma(av); return pol0_F2x(vT); }
    2771             : 
    2772          28 :     if (!equali1(lb)) res = F2xq_mul(res, F2xq_powu(lb, da - dc, T), T);
    2773          28 :     if (gc_needed(av,2))
    2774             :     {
    2775           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"F2xqX_resultant (da = %ld)",da);
    2776           0 :       gerepileall(av,3, &a,&b,&res);
    2777             :     }
    2778          28 :     da = db; /* = degpol(a) */
    2779          28 :     db = dc; /* = degpol(b) */
    2780             :   }
    2781          21 :   res = F2xq_mul(res, F2xq_powu(gel(b,2), da, T), T);
    2782          21 :   return gerepileupto(av, res);
    2783             : }
    2784             : 
    2785             : /* disc P = (-1)^(n(n-1)/2) lc(P)^(n - deg P' - 2) Res(P,P'), n = deg P */
    2786             : GEN
    2787          14 : F2xqX_disc(GEN P, GEN T)
    2788             : {
    2789          14 :   pari_sp av = avma;
    2790          14 :   GEN L, dP = F2xX_deriv(P), D = F2xqX_resultant(P, dP, T);
    2791             :   long dd;
    2792          14 :   if (!lgpol(D)) return pol_0(get_F2x_var(T));
    2793          14 :   dd = degpol(P) - 2 - degpol(dP); /* >= -1; > -1 iff p | deg(P) */
    2794          14 :   L = leading_coeff(P);
    2795          14 :   if (dd && !F2x_equal1(L))
    2796           0 :     D = (dd == -1)? F2xq_div(D,L,T): F2xq_mul(D, F2xq_powu(L, dd, T), T);
    2797          14 :   return gerepileupto(av, D);
    2798             : }
    2799             : 
    2800             : /***********************************************************************/
    2801             : /**                                                                   **/
    2802             : /**                             F2xqXQ                                **/
    2803             : /**                                                                   **/
    2804             : /***********************************************************************/
    2805             : 
    2806             : GEN
    2807      114520 : F2xqXQ_mul(GEN x, GEN y, GEN S, GEN T) {
    2808      114520 :   return F2xqX_rem(F2xqX_mul(x,y,T),S,T);
    2809             : }
    2810             : 
    2811             : GEN
    2812      104332 : F2xqXQ_sqr(GEN x, GEN S, GEN T) {
    2813      104332 :   return F2xqX_rem(F2xqX_sqr(x,T),S,T);
    2814             : }
    2815             : 
    2816             : GEN
    2817           7 : F2xqXQ_invsafe(GEN x, GEN S, GEN T)
    2818             : {
    2819           7 :   GEN V, z = F2xqX_extgcd(get_F2xqX_mod(S), x, T, NULL, &V);
    2820           7 :   if (degpol(z)) return NULL;
    2821           7 :   z = F2xq_invsafe(gel(z,2),T);
    2822           7 :   if (!z) return NULL;
    2823           7 :   return F2xqX_F2xq_mul(V, z, T);
    2824             : }
    2825             : 
    2826             : GEN
    2827           7 : F2xqXQ_inv(GEN x, GEN S, GEN T)
    2828             : {
    2829           7 :   pari_sp av = avma;
    2830           7 :   GEN U = F2xqXQ_invsafe(x, S, T);
    2831           7 :   if (!U) pari_err_INV("F2xqXQ_inv",x);
    2832           7 :   return gerepileupto(av, U);
    2833             : }
    2834             : 
    2835             : struct _F2xqXQ {
    2836             :   GEN T, S;
    2837             : };
    2838             : 
    2839             : static GEN
    2840      344325 : _F2xqXQ_add(void *data, GEN x, GEN y) {
    2841             :   (void) data;
    2842      344325 :   return F2xX_add(x,y);
    2843             : }
    2844             : static GEN
    2845      480696 : _F2xqXQ_cmul(void *data, GEN P, long a, GEN x) {
    2846             :   (void) data;
    2847      480696 :   return F2xX_F2x_mul(x,gel(P,a+2));
    2848             : }
    2849             : static GEN
    2850      352258 : _F2xqXQ_red(void *data, GEN x) {
    2851      352258 :   struct _F2xqXQ *d = (struct _F2xqXQ*) data;
    2852      352258 :   return F2xqX_red(x, d->T);
    2853             : }
    2854             : static GEN
    2855      114464 : _F2xqXQ_mul(void *data, GEN x, GEN y) {
    2856      114464 :   struct _F2xqXQ *d = (struct _F2xqXQ*) data;
    2857      114464 :   return F2xqXQ_mul(x,y, d->S,d->T);
    2858             : }
    2859             : static GEN
    2860       32984 : _F2xqXQ_sqr(void *data, GEN x) {
    2861       32984 :   struct _F2xqXQ *d = (struct _F2xqXQ*) data;
    2862       32984 :   return F2xqXQ_sqr(x, d->S,d->T);
    2863             : }
    2864             : static GEN
    2865          56 : _F2xqXQ_zero(void *data) {
    2866          56 :   struct _F2xqXQ *d = (struct _F2xqXQ*) data;
    2867          56 :   return pol_0(get_F2xqX_var(d->S));
    2868             : }
    2869             : static GEN
    2870      375393 : _F2xqXQ_one(void *data) {
    2871      375393 :   struct _F2xqXQ *d = (struct _F2xqXQ*) data;
    2872      375393 :   return pol1_F2xX(get_F2xqX_var(d->S),get_F2x_var(d->T));
    2873             : }
    2874             : 
    2875             : static struct bb_algebra F2xqXQ_algebra = { _F2xqXQ_red,
    2876             :  _F2xqXQ_add, _F2xqXQ_add, _F2xqXQ_mul,_F2xqXQ_sqr,_F2xqXQ_one,_F2xqXQ_zero };
    2877             : 
    2878             : GEN
    2879           7 : F2xqXQ_pow(GEN x, GEN n, GEN S, GEN T)
    2880             : {
    2881             :   struct _F2xqXQ D;
    2882           7 :   long s = signe(n);
    2883           7 :   if (!s) return pol1_F2xX(get_F2xqX_var(S), get_F2x_var(T));
    2884           7 :   if (s < 0) x = F2xqXQ_inv(x,S,T);
    2885           7 :   if (is_pm1(n)) return s < 0 ? x : gcopy(x);
    2886           7 :   if (degpol(x) >= get_F2xqX_degree(S)) x = F2xqX_rem(x,S,T);
    2887           7 :   D.T = F2x_get_red(T);
    2888           7 :   D.S = F2xqX_get_red(S, T);
    2889           7 :   return gen_pow(x, n, (void*)&D, &_F2xqXQ_sqr, &_F2xqXQ_mul);
    2890             : }
    2891             : 
    2892             : GEN
    2893        7462 : F2xqXQ_powers(GEN x, long l, GEN S, GEN T)
    2894             : {
    2895             :   struct _F2xqXQ D;
    2896        7462 :   int use_sqr = 2*degpol(x) >= get_F2xqX_degree(S);
    2897        7462 :   D.T = F2x_get_red(T);
    2898        7462 :   D.S = F2xqX_get_red(S, T);
    2899        7462 :   return gen_powers(x, l, use_sqr, (void*)&D, &_F2xqXQ_sqr, &_F2xqXQ_mul,&_F2xqXQ_one);
    2900             : }
    2901             : 
    2902             : GEN
    2903       13885 : F2xqX_F2xqXQV_eval(GEN P, GEN V, GEN S, GEN T)
    2904             : {
    2905             :   struct _F2xqXQ D;
    2906       13885 :   D.T = F2x_get_red(T);
    2907       13885 :   D.S = F2xqX_get_red(S, T);
    2908       13885 :   return gen_bkeval_powers(P, degpol(P), V, (void*)&D, &F2xqXQ_algebra,
    2909             :                                                    _F2xqXQ_cmul);
    2910             : }
    2911             : 
    2912             : GEN
    2913      122542 : F2xqX_F2xqXQ_eval(GEN Q, GEN x, GEN S, GEN T)
    2914             : {
    2915             :   struct _F2xqXQ D;
    2916      122542 :   int use_sqr = 2*degpol(x) >= get_F2xqX_degree(S);
    2917      122542 :   D.T = F2x_get_red(T);
    2918      122542 :   D.S = F2xqX_get_red(S, T);
    2919      122542 :   return gen_bkeval(Q, degpol(Q), x, use_sqr, (void*)&D, &F2xqXQ_algebra,
    2920             :                                                     _F2xqXQ_cmul);
    2921             : }
    2922             : 
    2923             : static GEN
    2924       89061 : F2xqXQ_autpow_sqr(void * E, GEN x)
    2925             : {
    2926       89061 :   struct _F2xqXQ *D = (struct _F2xqXQ *)E;
    2927       89061 :   GEN T = D->T;
    2928       89061 :   GEN phi = gel(x,1), S = gel(x,2);
    2929       89061 :   long n = brent_kung_optpow(get_F2x_degree(T)-1,lgpol(S)+1,1);
    2930       89061 :   GEN V = F2xq_powers(phi, n, T);
    2931       89061 :   GEN phi2 = F2x_F2xqV_eval(phi, V, T);
    2932       89061 :   GEN Sphi = F2xY_F2xqV_evalx(S, V, T);
    2933       89061 :   GEN S2 = F2xqX_F2xqXQ_eval(Sphi, S, D->S, T);
    2934       89061 :   return mkvec2(phi2, S2);
    2935             : }
    2936             : 
    2937             : static GEN
    2938       33481 : F2xqXQ_autpow_mul(void * E, GEN x, GEN y)
    2939             : {
    2940       33481 :   struct _F2xqXQ *D = (struct _F2xqXQ *)E;
    2941       33481 :   GEN T = D->T;
    2942       33481 :   GEN phi1 = gel(x,1), S1 = gel(x,2);
    2943       33481 :   GEN phi2 = gel(y,1), S2 = gel(y,2);
    2944       33481 :   long n = brent_kung_optpow(get_F2x_degree(T)-1,lgpol(S1)+1,1);
    2945       33481 :   GEN V = F2xq_powers(phi2, n, T);
    2946       33481 :   GEN phi3 = F2x_F2xqV_eval(phi1,V,T);
    2947       33481 :   GEN Sphi = F2xY_F2xqV_evalx(S1,V,T);
    2948       33481 :   GEN S3 = F2xqX_F2xqXQ_eval(Sphi, S2, D->S, T);
    2949       33481 :   return mkvec2(phi3, S3);
    2950             : }
    2951             : 
    2952             : GEN
    2953       70679 : F2xqXQ_autpow(GEN aut, long n, GEN S, GEN T)
    2954             : {
    2955             :   struct _F2xqXQ D;
    2956       70679 :   D.T = F2x_get_red(T);
    2957       70679 :   D.S = F2xqX_get_red(S, T);
    2958       70679 :   return gen_powu(aut,n,&D,F2xqXQ_autpow_sqr,F2xqXQ_autpow_mul);
    2959             : }
    2960             : 
    2961             : static GEN
    2962        6916 : F2xqXQ_auttrace_mul(void *E, GEN x, GEN y)
    2963             : {
    2964        6916 :   struct _F2xqXQ *D = (struct _F2xqXQ *) E;
    2965        6916 :   GEN T = D->T;
    2966        6916 :   GEN phi1 = gel(x,1), S1 = gel(x,2), a1 = gel(x,3);
    2967        6916 :   GEN phi2 = gel(y,1), S2 = gel(y,2), a2 = gel(y,3);
    2968        6916 :   long n2 = brent_kung_optpow(get_F2x_degree(T)-1, lgpol(S1)+lgpol(a1)+1, 1);
    2969        6916 :   GEN V2 = F2xq_powers(phi2, n2, T);
    2970        6916 :   GEN phi3 = F2x_F2xqV_eval(phi1, V2, T);
    2971        6916 :   GEN Sphi = F2xY_F2xqV_evalx(S1, V2, T);
    2972        6916 :   GEN aphi = F2xY_F2xqV_evalx(a1, V2, T);
    2973        6916 :   long n = brent_kung_optpow(maxss(degpol(Sphi),degpol(aphi)),2,1);
    2974        6916 :   GEN V = F2xqXQ_powers(S2, n, D->S, T);
    2975        6916 :   GEN S3 = F2xqX_F2xqXQV_eval(Sphi, V, D->S, T);
    2976        6916 :   GEN aS = F2xqX_F2xqXQV_eval(aphi, V, D->S, T);
    2977        6916 :   GEN a3 = F2xX_add(aS, a2);
    2978        6916 :   return mkvec3(phi3, S3, a3);
    2979             : }
    2980             : 
    2981             : static GEN
    2982        4641 : F2xqXQ_auttrace_sqr(void *E, GEN x) { return F2xqXQ_auttrace_mul(E, x, x); }
    2983             : GEN
    2984        2646 : F2xqXQ_auttrace(GEN aut, long n, GEN S, GEN T)
    2985             : {
    2986             :   struct _F2xqXQ D;
    2987        2646 :   D.T = F2x_get_red(T);
    2988        2646 :   D.S = F2xqX_get_red(S, T);
    2989        2646 :   return gen_powu(aut,n,&D,F2xqXQ_auttrace_sqr,F2xqXQ_auttrace_mul);
    2990             : }
    2991             : 
    2992             : GEN
    2993          63 : F2xqXQV_red(GEN x, GEN S, GEN T)
    2994         273 : { pari_APPLY_type(t_VEC, F2xqX_rem(gel(x,i), S, T)) }

Generated by: LCOV version 1.13