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 - FpXX.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.16.2 lcov report (development 29115-f22e516b23) Lines: 951 1215 78.3 %
Date: 2024-03-19 08:06:49 Functions: 114 134 85.1 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2012  The PARI group.
       2             : 
       3             : This file is part of the PARI/GP package.
       4             : 
       5             : PARI/GP is free software; you can redistribute it and/or modify it under the
       6             : terms of the GNU General Public License as published by the Free Software
       7             : Foundation; either version 2 of the License, or (at your option) any later
       8             : version. It is distributed in the hope that it will be useful, but WITHOUT
       9             : ANY WARRANTY WHATSOEVER.
      10             : 
      11             : Check the License for details. You should have received a copy of it, along
      12             : with the package; see the file 'COPYING'. If not, write to the Free Software
      13             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      14             : 
      15             : #include "pari.h"
      16             : #include "paripriv.h"
      17             : 
      18             : /* Not so fast arithmetic with polynomials over FpX */
      19             : 
      20             : /*******************************************************************/
      21             : /*                                                                 */
      22             : /*                             FpXX                                */
      23             : /*                                                                 */
      24             : /*******************************************************************/
      25             : /*Polynomials whose coefficients are either polynomials or integers*/
      26             : 
      27             : static GEN
      28       53770 : to_ZX(GEN a, long v) { return typ(a)==t_INT? scalarpol_shallow(a,v): a; }
      29             : 
      30             : static ulong
      31     1250897 : to_FlxqX(GEN P, GEN Q, GEN T, GEN p, GEN *pt_P, GEN *pt_Q, GEN *pt_T)
      32             : {
      33     1250897 :   ulong pp = uel(p,2);
      34     1250897 :   long v = get_FpX_var(T);
      35     1250883 :   *pt_P = ZXX_to_FlxX(P, pp, v);
      36     1250789 :   if (pt_Q) *pt_Q = ZXX_to_FlxX(Q, pp, v);
      37     1250706 :   *pt_T = ZXT_to_FlxT(T, pp);
      38     1250802 :   return pp;
      39             : }
      40             : 
      41             : static GEN
      42         126 : ZXX_copy(GEN a) { return gcopy(a); }
      43             : 
      44             : GEN
      45       39934 : FpXX_red(GEN z, GEN p)
      46             : {
      47             :   GEN res;
      48       39934 :   long i, l = lg(z);
      49       39934 :   res = cgetg(l,t_POL); res[1] = z[1];
      50      281451 :   for (i=2; i<l; i++)
      51             :   {
      52      241517 :     GEN zi = gel(z,i), c;
      53      241517 :     if (typ(zi)==t_INT)
      54       14672 :       c = modii(zi,p);
      55             :     else
      56             :     {
      57      226845 :       pari_sp av = avma;
      58      226845 :       c = FpX_red(zi,p);
      59      226845 :       switch(lg(c)) {
      60           7 :         case 2: set_avma(av); c = gen_0; break;
      61       20516 :         case 3: c = gerepilecopy(av, gel(c,2)); break;
      62             :       }
      63      241517 :     }
      64      241517 :     gel(res,i) = c;
      65             :   }
      66       39934 :   return FpXX_renormalize(res,lg(res));
      67             : }
      68             : GEN
      69      446115 : FpXX_add(GEN x, GEN y, GEN p)
      70             : {
      71             :   long i,lz;
      72             :   GEN z;
      73      446115 :   long lx=lg(x);
      74      446115 :   long ly=lg(y);
      75      446115 :   if (ly>lx) swapspec(x,y, lx,ly);
      76      446115 :   lz = lx; z = cgetg(lz, t_POL); z[1]=x[1];
      77     9338123 :   for (i=2; i<ly; i++) gel(z,i) = Fq_add(gel(x,i), gel(y,i), NULL, p);
      78     1671102 :   for (   ; i<lx; i++) gel(z,i) = gcopy(gel(x,i));
      79      446115 :   return FpXX_renormalize(z, lz);
      80             : }
      81             : GEN
      82       30880 : FpXX_sub(GEN x, GEN y, GEN p)
      83             : {
      84             :   long i,lz;
      85             :   GEN z;
      86       30880 :   long lx=lg(x);
      87       30880 :   long ly=lg(y);
      88       30880 :   if (ly <= lx)
      89             :   {
      90       15521 :     lz = lx; z = cgetg(lz, t_POL); z[1]=x[1];
      91      197912 :     for (i=2; i<ly; i++) gel(z,i) = Fq_sub(gel(x,i), gel(y,i), NULL, p);
      92       33415 :     for (   ; i<lx; i++) gel(z,i) = gcopy(gel(x,i));
      93             :   }
      94             :   else
      95             :   {
      96       15359 :     lz = ly; z = cgetg(lz, t_POL); z[1]=x[1];
      97       94153 :     for (i=2; i<lx; i++) gel(z,i) = Fq_sub(gel(x,i), gel(y,i), NULL, p);
      98       56661 :     for (   ; i<ly; i++) gel(z,i) = Fq_neg(gel(y,i), NULL, p);
      99             :   }
     100       30880 :   return FpXX_renormalize(z, lz);
     101             : }
     102             : 
     103             : static GEN
     104      152221 : FpXX_subspec(GEN x, GEN y, GEN p, long nx, long ny)
     105             : {
     106             :   long i,lz;
     107             :   GEN z;
     108      152221 :   if (ny <= nx)
     109             :   {
     110      152221 :     lz = nx+2; z = cgetg(lz, t_POL);
     111     3714434 :     for (i=0; i<ny; i++) gel(z,i+2) = Fq_sub(gel(x,i), gel(y,i), NULL, p);
     112      152221 :     for (   ; i<nx; i++) gel(z,i+2) = gcopy(gel(x,i));
     113             :   }
     114             :   else
     115             :   {
     116           0 :     lz = ny+2; z = cgetg(lz, t_POL);
     117           0 :     for (i=0; i<nx; i++) gel(z,i+2) = Fq_sub(gel(x,i), gel(y,i), NULL, p);
     118           0 :     for (   ; i<ny; i++) gel(z,i+2) = Fq_neg(gel(y,i), NULL, p);
     119             :   }
     120      152221 :   z[1] = 0; return FpXX_renormalize(z, lz);
     121             : }
     122             : 
     123             : GEN
     124        2073 : FpXX_neg(GEN x, GEN p)
     125             : {
     126        2073 :   long i, lx = lg(x);
     127        2073 :   GEN y = cgetg(lx,t_POL);
     128        2073 :   y[1] = x[1];
     129       53964 :   for(i=2; i<lx; i++) gel(y,i) = Fq_neg(gel(x,i), NULL, p);
     130        2073 :   return FpXX_renormalize(y, lx);
     131             : }
     132             : 
     133             : GEN
     134       56448 : FpXX_Fp_mul(GEN P, GEN u, GEN p)
     135             : {
     136             :   long i, lP;
     137       56448 :   GEN res = cgetg_copy(P, &lP); res[1] = P[1];
     138      484271 :   for(i=2; i<lP; i++)
     139             :   {
     140      427823 :     GEN x = gel(P,i);
     141      427823 :     gel(res,i) = typ(x)==t_INT? Fp_mul(x,u,p): FpX_Fp_mul(x,u,p);
     142             :   }
     143       56448 :   return FpXX_renormalize(res,lP);
     144             : }
     145             : 
     146             : GEN
     147        7116 : FpXX_mulu(GEN P, ulong u, GEN p)
     148             : {
     149             :   long i, lP;
     150        7116 :   GEN res = cgetg_copy(P, &lP); res[1] = P[1];
     151       52335 :   for(i=2; i<lP; i++)
     152             :   {
     153       45219 :     GEN x = gel(P,i);
     154       45219 :     gel(res,i) = typ(x)==t_INT? Fp_mulu(x,u,p): FpX_mulu(x,u,p);
     155             :   }
     156        7116 :   return FpXX_renormalize(res,lP);
     157             : }
     158             : 
     159             : GEN
     160        2240 : FpXX_halve(GEN P, GEN p)
     161             : {
     162             :   long i, lP;
     163        2240 :   GEN res = cgetg_copy(P, &lP); res[1] = P[1];
     164        7770 :   for(i=2; i<lP; i++)
     165             :   {
     166        5530 :     GEN x = gel(P,i);
     167        5530 :     gel(res,i) = typ(x)==t_INT? Fp_halve(x,p): FpX_halve(x,p);
     168             :   }
     169        2240 :   return FpXX_renormalize(res,lP);
     170             : }
     171             : 
     172             : GEN
     173       13405 : FpXX_deriv(GEN P, GEN p)
     174             : {
     175       13405 :   long i, l = lg(P)-1;
     176             :   GEN res;
     177             : 
     178       13405 :   if (l < 3) return pol_0(varn(P));
     179       13055 :   res = cgetg(l, t_POL);
     180       13055 :   res[1] = P[1];
     181       80960 :   for (i=2; i<l ; i++)
     182             :   {
     183       67905 :     GEN x = gel(P,i+1);
     184       67905 :     gel(res,i) = typ(x)==t_INT? Fp_mulu(x,i-1,p): FpX_mulu(x,i-1,p);
     185             :   }
     186       13055 :   return FpXX_renormalize(res, l);
     187             : }
     188             : 
     189             : GEN
     190           0 : FpXX_integ(GEN P, GEN p)
     191             : {
     192           0 :   long i, l = lg(P);
     193             :   GEN res;
     194             : 
     195           0 :   if (l == 2) return pol_0(varn(P));
     196           0 :   res = cgetg(l+1, t_POL);
     197           0 :   res[1] = P[1];
     198           0 :   gel(res,2) = gen_0;
     199           0 :   for (i=3; i<=l ; i++)
     200             :   {
     201           0 :     GEN x = gel(P,i-1);
     202           0 :     if (signe(x))
     203             :     {
     204           0 :       GEN i1 = Fp_inv(utoi(i-2), p);
     205           0 :       gel(res,i) = typ(x)==t_INT? Fp_mul(x,i1,p): FpX_Fp_mul(x,i1,p);
     206             :     } else
     207           0 :       gel(res,i) = gen_0;
     208             :   }
     209           0 :   return FpXX_renormalize(res, l+1);
     210             : }
     211             : 
     212             : /*******************************************************************/
     213             : /*                                                                 */
     214             : /*                             (Fp[X]/(Q))[Y]                      */
     215             : /*                                                                 */
     216             : /*******************************************************************/
     217             : 
     218             : static GEN
     219     1332457 : get_FpXQX_red(GEN T, GEN *B)
     220             : {
     221     1332457 :   if (typ(T)!=t_VEC) { *B=NULL; return T; }
     222       94923 :   *B = gel(T,1); return gel(T,2);
     223             : }
     224             : 
     225             : GEN
     226          52 : random_FpXQX(long d1, long v, GEN T, GEN p)
     227             : {
     228          52 :   long dT = get_FpX_degree(T), vT = get_FpX_var(T);
     229          52 :   long i, d = d1+2;
     230          52 :   GEN y = cgetg(d,t_POL); y[1] = evalsigne(1) | evalvarn(v);
     231         284 :   for (i=2; i<d; i++) gel(y,i) = random_FpX(dT, vT, p);
     232          52 :   return FpXQX_renormalize(y,d);
     233             : }
     234             : 
     235             : /*Not stack clean*/
     236             : GEN
     237     1890802 : Kronecker_to_FpXQX(GEN Z, GEN T, GEN p)
     238             : {
     239     1890802 :   long i,j,lx,l, N = (get_FpX_degree(T)<<1) + 1;
     240     1890772 :   GEN x, t = cgetg(N,t_POL), z = FpX_red(Z, p);
     241     1889717 :   t[1] = evalvarn(get_FpX_var(T));
     242     1889781 :   l = lg(z); lx = (l-2) / (N-2);
     243     1889781 :   x = cgetg(lx+3,t_POL);
     244     1890513 :   x[1] = z[1];
     245    29935212 :   for (i=2; i<lx+2; i++)
     246             :   {
     247   238507419 :     for (j=2; j<N; j++) gel(t,j) = gel(z,j);
     248    28044691 :     z += (N-2);
     249    28044691 :     gel(x,i) = FpX_rem(FpX_renormalize(t,N), T,p);
     250             :   }
     251     1890521 :   N = (l-2) % (N-2) + 2;
     252     3381483 :   for (j=2; j<N; j++) gel(t,j) = gel(z,j);
     253     1890521 :   gel(x,i) = FpX_rem(FpX_renormalize(t,N), T,p);
     254     1890597 :   return FpXQX_renormalize(x, i+1);
     255             : }
     256             : 
     257             : GEN
     258     1936685 : FpXQX_red(GEN z, GEN T, GEN p)
     259             : {
     260     1936685 :   long i, l = lg(z);
     261     1936685 :   GEN res = cgetg(l,t_POL); res[1] = z[1];
     262    16297517 :   for(i=2;i<l;i++)
     263    14364551 :     if (typ(gel(z,i)) == t_INT)
     264      160149 :       gel(res,i) = modii(gel(z,i),p);
     265             :     else
     266    14204402 :       gel(res,i) = FpXQ_red(gel(z,i),T,p);
     267     1932966 :   return FpXQX_renormalize(res,l);
     268             : }
     269             : 
     270             : GEN
     271           0 : FpXQXV_red(GEN x, GEN T, GEN p)
     272           0 : { pari_APPLY_type(t_VEC, FpXQX_red(gel(x,i), T, p)) }
     273             : 
     274             : GEN
     275           0 : FpXQXT_red(GEN x, GEN T, GEN p)
     276             : {
     277           0 :   if (typ(x) == t_POL)
     278           0 :     return FpXQX_red(x, T, p);
     279             :   else
     280           0 :     pari_APPLY_type(t_VEC, FpXQXT_red(gel(x,i), T, p))
     281             : }
     282             : 
     283             : static GEN
     284        2191 : to_intmod(GEN x, GEN p) { retmkintmod(modii(x, p), p); }
     285             : 
     286             : GEN
     287         525 : FpXQX_to_mod(GEN z, GEN T, GEN p)
     288             : {
     289         525 :   long i,l = lg(z);
     290         525 :   GEN x = cgetg(l, t_POL);
     291         525 :   x[1] = z[1];
     292         525 :   if (l == 2) return x;
     293         525 :   p = icopy(p);
     294         525 :   T = FpX_to_mod_raw(T, p);
     295        6692 :   for (i=2; i<l; i++)
     296             :   {
     297        6167 :     GEN zi = gel(z,i);
     298        6167 :     gel(x,i) = typ(zi) == t_POL? mkpolmod(FpX_to_mod_raw(zi, p), T)
     299        6167 :                                : to_intmod(zi, p);
     300             :   }
     301         525 :   return normalizepol_lg(x,l);
     302             : }
     303             : 
     304             : static GEN
     305           0 : FpXQX_to_mod_raw(GEN z, GEN T, GEN p)
     306             : {
     307           0 :   long i,l = lg(z);
     308           0 :   GEN x = cgetg(l, t_POL);
     309           0 :   x[1] = z[1];
     310           0 :   if (l == 2) return x;
     311           0 :   for (i=2; i<l; i++)
     312             :   {
     313           0 :     GEN zi = gel(z,i);
     314           0 :     gel(x,i) = typ(zi) == t_POL? mkpolmod(FpX_to_mod_raw(zi, p), T)
     315           0 :                                : to_intmod(zi, p);
     316             :   }
     317           0 :   return normalizepol_lg(x,l);
     318             : }
     319             : 
     320             : INLINE GEN
     321           0 : FqX_to_mod_raw(GEN f, GEN T, GEN p)
     322           0 : { return T?FpXQX_to_mod_raw(f, T, p): FpX_to_mod_raw(f, p); }
     323             : 
     324             : static GEN
     325           0 : FqXC_to_mod_raw(GEN x, GEN T, GEN p)
     326           0 : { pari_APPLY_type(t_COL, FqX_to_mod_raw(gel(x,i), T, p)) }
     327             : 
     328             : GEN
     329          14 : FqXC_to_mod(GEN z, GEN T, GEN p)
     330             : {
     331             :   GEN x;
     332          14 :   long i,l = lg(z);
     333          14 :   if (!T) return FpXC_to_mod(z, p);
     334           0 :   x = cgetg(l, t_COL);
     335           0 :   if (l == 1) return x;
     336           0 :   p = icopy(p);
     337           0 :   T = FpX_to_mod_raw(T, p);
     338           0 :   for (i=1; i<l; i++)
     339           0 :     gel(x,i) = FqX_to_mod_raw(gel(z, i), T, p);
     340           0 :   return x;
     341             : }
     342             : 
     343             : GEN
     344           0 : FqXM_to_mod(GEN z, GEN T, GEN p)
     345             : {
     346             :   GEN x;
     347           0 :   long i,l = lg(z);
     348           0 :   if (!T) return FpXM_to_mod(z, p);
     349           0 :   x = cgetg(l, t_MAT);
     350           0 :   if (l == 1) return x;
     351           0 :   p = icopy(p);
     352           0 :   T = FpX_to_mod_raw(T, p);
     353           0 :   for (i=1; i<l; i++)
     354           0 :     gel(x,i) = FqXC_to_mod_raw(gel(z, i), T, p);
     355           0 :   return x;
     356             : }
     357             : 
     358             : static int
     359     3639359 : ZXX_is_ZX_spec(GEN a,long na)
     360             : {
     361             :   long i;
     362     3950837 :   for(i=0;i<na;i++)
     363     3890270 :     if(typ(gel(a,i))!=t_INT) return 0;
     364       60567 :   return 1;
     365             : }
     366             : 
     367             : static int
     368      246002 : ZXX_is_ZX(GEN a) { return ZXX_is_ZX_spec(a+2,lgpol(a)); }
     369             : 
     370             : static GEN
     371      140790 : FpXX_FpX_mulspec(GEN P, GEN U, GEN p, long v, long lU)
     372             : {
     373      140790 :   long i, lP =lg(P);
     374             :   GEN res;
     375      140790 :   res = cgetg(lP, t_POL); res[1] = P[1];
     376     7695966 :   for(i=2; i<lP; i++)
     377             :   {
     378     7555176 :     GEN Pi = gel(P,i);
     379     7555176 :     gel(res,i) = typ(Pi)==t_INT? FpX_Fp_mulspec(U, Pi, p, lU):
     380     7540499 :                                  FpX_mulspec(U, Pi+2, p, lU, lgpol(Pi));
     381     7555176 :     setvarn(gel(res,i),v);
     382             :   }
     383      140790 :   return FpXQX_renormalize(res,lP);
     384             : }
     385             : 
     386             : GEN
     387      125645 : FpXX_FpX_mul(GEN P, GEN U, GEN p)
     388      125645 : { return FpXX_FpX_mulspec(P,U+2,p,varn(U),lgpol(U)); }
     389             : 
     390             : static GEN
     391       15145 : FpXY_FpY_mulspec(GEN x, GEN y, GEN T, GEN p, long lx, long ly)
     392             : {
     393       15145 :   pari_sp av = avma;
     394       15145 :   long v = fetch_var();
     395       15145 :   GEN z = RgXY_swapspec(x,get_FpX_degree(T)-1,v,lx);
     396       15145 :   z = FpXX_FpX_mulspec(z,y,p,v,ly);
     397       15145 :   z = RgXY_swapspec(z+2,lx+ly+3,get_FpX_var(T),lgpol(z));
     398       15145 :   (void)delete_var(); return gerepilecopy(av,z);
     399             : }
     400             : 
     401             : static GEN
     402     1696519 : FpXQX_mulspec(GEN x, GEN y, GEN T, GEN p, long lx, long ly)
     403             : {
     404     1696519 :   pari_sp av = avma;
     405             :   GEN z, kx, ky;
     406             :   long n;
     407     1696519 :   if (ZXX_is_ZX_spec(y,ly))
     408             :   {
     409       15658 :     if (ZXX_is_ZX_spec(x,lx))
     410        8271 :       return FpX_mulspec(x,y,p,lx,ly);
     411             :     else
     412        7387 :       return FpXY_FpY_mulspec(x,y,T,p,lx,ly);
     413     1681112 :   } else if (ZXX_is_ZX_spec(x,lx))
     414        7758 :       return FpXY_FpY_mulspec(y,x,T,p,ly,lx);
     415     1673571 :   n = get_FpX_degree(T);
     416     1673597 :   kx = RgXX_to_Kronecker_spec(x, lx, n);
     417     1673690 :   ky = RgXX_to_Kronecker_spec(y, ly, n);
     418     1673816 :   z = Kronecker_to_FpXQX(ZX_mul(ky,kx), T, p);
     419     1673449 :   return gerepileupto(av, z);
     420             : }
     421             : 
     422             : GEN
     423     1384214 : FpXQX_mul(GEN x, GEN y, GEN T, GEN p)
     424             : {
     425     1384214 :   GEN z = FpXQX_mulspec(x+2,y+2,T,p,lgpol(x),lgpol(y));
     426     1384789 :   setvarn(z,varn(x)); return z;
     427             : }
     428             : 
     429             : GEN
     430      186666 : FpXQX_sqr(GEN x, GEN T, GEN p)
     431             : {
     432      186666 :   pari_sp av = avma;
     433             :   GEN z, kx;
     434      186666 :   if (ZXX_is_ZX(x)) return ZX_sqr(x);
     435      179492 :   kx= RgXX_to_Kronecker(x, get_FpX_degree(T));
     436      179492 :   z = Kronecker_to_FpXQX(ZX_sqr(kx), T, p);
     437      179492 :   return gerepileupto(av, z);
     438             : }
     439             : 
     440             : GEN
     441      552912 : FpXQX_FpXQ_mul(GEN P, GEN U, GEN T, GEN p)
     442             : {
     443             :   long i, lP;
     444             :   GEN res;
     445      552912 :   res = cgetg_copy(P, &lP); res[1] = P[1];
     446     2064699 :   for(i=2; i<lP; i++)
     447     2746489 :     gel(res,i) = typ(gel(P,i))==t_INT? FpX_Fp_mul(U, gel(P,i), p):
     448     1234615 :                                        FpXQ_mul(U, gel(P,i), T,p);
     449      552704 :   return FpXQX_renormalize(res,lP);
     450             : }
     451             : 
     452             : /* x and y in Z[Y][X]. Assume T irreducible mod p */
     453             : static GEN
     454      173394 : FpXQX_divrem_basecase(GEN x, GEN y, GEN T, GEN p, GEN *pr)
     455             : {
     456             :   long vx, dx, dy, dy1, dz, i, j, sx, lr;
     457             :   pari_sp av0, av, tetpil;
     458             :   GEN z,p1,rem,lead;
     459             : 
     460      173394 :   if (!signe(y)) pari_err_INV("FpX_divrem",y);
     461      173394 :   vx=varn(x); dy=degpol(y); dx=degpol(x);
     462      173394 :   if (dx < dy)
     463             :   {
     464         185 :     if (pr)
     465             :     {
     466         135 :       av0 = avma; x = FpXQX_red(x, T, p);
     467         135 :       if (pr == ONLY_DIVIDES) { set_avma(av0); return signe(x)? NULL: pol_0(vx); }
     468         135 :       if (pr == ONLY_REM) return x;
     469         135 :       *pr = x;
     470             :     }
     471         185 :     return pol_0(vx);
     472             :   }
     473      173209 :   lead = leading_coeff(y);
     474      173209 :   if (!dy) /* y is constant */
     475             :   {
     476        1408 :     if (pr && pr != ONLY_DIVIDES)
     477             :     {
     478        1044 :       if (pr == ONLY_REM) return pol_0(vx);
     479           7 :       *pr = pol_0(vx);
     480             :     }
     481         371 :     if (gequal1(lead)) return FpXQX_red(x,T,p);
     482         355 :     av0 = avma; x = FqX_Fq_mul(x, Fq_inv(lead, T,p), T,p);
     483         355 :     return gerepileupto(av0,x);
     484             :   }
     485      171801 :   av0 = avma; dz = dx-dy;
     486      171801 :   lead = gequal1(lead)? NULL: gclone(Fq_inv(lead,T,p));
     487      171801 :   set_avma(av0);
     488      171801 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
     489      171801 :   x += 2; y += 2; z += 2;
     490      178279 :   for (dy1=dy-1; dy1>=0 && !signe(gel(y, dy1)); dy1--);
     491             : 
     492      171801 :   p1 = gel(x,dx); av = avma;
     493      171801 :   gel(z,dz) = lead? gerepileupto(av, Fq_mul(p1,lead, T, p)): gcopy(p1);
     494      521194 :   for (i=dx-1; i>=dy; i--)
     495             :   {
     496      349393 :     av=avma; p1=gel(x,i);
     497     1164457 :     for (j=i-dy1; j<=i && j<=dz; j++)
     498      815064 :       p1 = Fq_sub(p1, Fq_mul(gel(z,j),gel(y,i-j),NULL,p),NULL,p);
     499      349393 :     if (lead) p1 = Fq_mul(p1, lead, NULL,p);
     500      349393 :     tetpil=avma; gel(z,i-dy) = gerepile(av,tetpil,Fq_red(p1,T,p));
     501             :   }
     502      171801 :   if (!pr) { guncloneNULL(lead); return z-2; }
     503             : 
     504      168642 :   rem = (GEN)avma; av = (pari_sp)new_chunk(dx+3);
     505      180685 :   for (sx=0; ; i--)
     506             :   {
     507      180685 :     p1 = gel(x,i);
     508      692376 :     for (j=maxss(0,i-dy1); j<=i && j<=dz; j++)
     509      511691 :       p1 = Fq_sub(p1, Fq_mul(gel(z,j),gel(y,i-j),NULL,p),NULL,p);
     510      180685 :     tetpil=avma; p1 = Fq_red(p1, T, p); if (signe(p1)) { sx = 1; break; }
     511       13968 :     if (!i) break;
     512       12043 :     set_avma(av);
     513             :   }
     514      168642 :   if (pr == ONLY_DIVIDES)
     515             :   {
     516           0 :     guncloneNULL(lead);
     517           0 :     if (sx) return gc_NULL(av0);
     518           0 :     return gc_const((pari_sp)rem, z-2);
     519             :   }
     520      168642 :   lr=i+3; rem -= lr;
     521      168642 :   rem[0] = evaltyp(t_POL) | _evallg(lr);
     522      168642 :   rem[1] = z[-1];
     523      168642 :   p1 = gerepile((pari_sp)rem,tetpil,p1);
     524      168642 :   rem += 2; gel(rem,i) = p1;
     525     1479878 :   for (i--; i>=0; i--)
     526             :   {
     527     1311236 :     av=avma; p1 = gel(x,i);
     528     4109472 :     for (j=maxss(0,i-dy1); j<=i && j<=dz; j++)
     529     2798236 :       p1 = Fq_sub(p1, Fq_mul(gel(z,j),gel(y,i-j), NULL,p), NULL,p);
     530     1311236 :     tetpil=avma; gel(rem,i) = gerepile(av,tetpil, Fq_red(p1, T, p));
     531             :   }
     532      168642 :   rem -= 2;
     533      168642 :   guncloneNULL(lead);
     534      168642 :   if (!sx) (void)FpXQX_renormalize(rem, lr);
     535      168642 :   if (pr == ONLY_REM) return gerepileupto(av0,rem);
     536       13575 :   *pr = rem; return z-2;
     537             : }
     538             : 
     539             : static GEN
     540         752 : FpXQX_addmulmul(GEN u, GEN v, GEN x, GEN y, GEN T, GEN p)
     541             : {
     542         752 :   return FpXX_add(FpXQX_mul(u, x, T, p),FpXQX_mul(v, y, T, p), p);
     543             : }
     544             : 
     545             : static GEN
     546         376 : FpXQXM_FpXQX_mul2(GEN M, GEN x, GEN y, GEN T, GEN p)
     547             : {
     548         376 :   GEN res = cgetg(3, t_COL);
     549         376 :   gel(res, 1) = FpXQX_addmulmul(gcoeff(M,1,1), gcoeff(M,1,2), x, y, T, p);
     550         376 :   gel(res, 2) = FpXQX_addmulmul(gcoeff(M,2,1), gcoeff(M,2,2), x, y, T, p);
     551         376 :   return res;
     552             : }
     553             : 
     554             : static GEN
     555         161 : FpXQXM_mul2(GEN A, GEN B, GEN T, GEN p)
     556             : {
     557         161 :   GEN A11=gcoeff(A,1,1),A12=gcoeff(A,1,2), B11=gcoeff(B,1,1),B12=gcoeff(B,1,2);
     558         161 :   GEN A21=gcoeff(A,2,1),A22=gcoeff(A,2,2), B21=gcoeff(B,2,1),B22=gcoeff(B,2,2);
     559         161 :   GEN M1 = FpXQX_mul(FpXX_add(A11,A22, p), FpXX_add(B11,B22, p), T, p);
     560         161 :   GEN M2 = FpXQX_mul(FpXX_add(A21,A22, p), B11, T, p);
     561         161 :   GEN M3 = FpXQX_mul(A11, FpXX_sub(B12,B22, p), T, p);
     562         161 :   GEN M4 = FpXQX_mul(A22, FpXX_sub(B21,B11, p), T, p);
     563         161 :   GEN M5 = FpXQX_mul(FpXX_add(A11,A12, p), B22, T, p);
     564         161 :   GEN M6 = FpXQX_mul(FpXX_sub(A21,A11, p), FpXX_add(B11,B12, p), T, p);
     565         161 :   GEN M7 = FpXQX_mul(FpXX_sub(A12,A22, p), FpXX_add(B21,B22, p), T, p);
     566         161 :   GEN T1 = FpXX_add(M1,M4, p), T2 = FpXX_sub(M7,M5, p);
     567         161 :   GEN T3 = FpXX_sub(M1,M2, p), T4 = FpXX_add(M3,M6, p);
     568         161 :   retmkmat22(FpXX_add(T1,T2, p), FpXX_add(M3,M5, p),
     569             :              FpXX_add(M2,M4, p), FpXX_add(T3,T4, p));
     570             : }
     571             : /* Return [0,1;1,-q]*M */
     572             : static GEN
     573         161 : FpXQX_FpXQXM_qmul(GEN q, GEN M, GEN T, GEN p)
     574             : {
     575         161 :   GEN u = FpXQX_mul(gcoeff(M,2,1), q, T, p);
     576         161 :   GEN v = FpXQX_mul(gcoeff(M,2,2), q, T, p);
     577         161 :   retmkmat22(gcoeff(M,2,1), gcoeff(M,2,2),
     578             :     FpXX_sub(gcoeff(M,1,1), u, p), FpXX_sub(gcoeff(M,1,2), v, p));
     579             : }
     580             : 
     581             : static GEN
     582           0 : matid2_FpXQXM(long v)
     583           0 : { retmkmat22(pol_1(v),pol_0(v),pol_0(v),pol_1(v)); }
     584             : 
     585             : static GEN
     586           0 : matJ2_FpXQXM(long v)
     587           0 : { retmkmat22(pol_0(v),pol_1(v),pol_1(v),pol_0(v)); }
     588             : 
     589             : static GEN
     590       19490 : FpXX_shift(GEN a, long n) { return RgX_shift_shallow(a, n); }
     591             : 
     592             : INLINE GEN
     593        8428 : FpXXn_red(GEN a, long n) { return RgXn_red_shallow(a, n); }
     594             : 
     595             : /* Fast resultant formula from William Hart in Flint <http://flintlib.org/> */
     596             : 
     597             : struct FpXQX_res
     598             : {
     599             :    GEN res, lc;
     600             :    long deg0, deg1, off;
     601             : };
     602             : 
     603             : INLINE void
     604           0 : FpXQX_halfres_update(long da, long db, long dr, GEN T, GEN p, struct FpXQX_res *res)
     605             : {
     606           0 :   if (dr >= 0)
     607             :   {
     608           0 :     if (!ZX_equal1(res->lc))
     609             :     {
     610           0 :       res->lc  = FpXQ_powu(res->lc, da - dr, T, p);
     611           0 :       res->res = FpXQ_mul(res->res, res->lc, T, p);
     612             :     }
     613           0 :     if (both_odd(da + res->off, db + res->off))
     614           0 :       res->res = FpX_neg(res->res, p);
     615             :   } else
     616             :   {
     617           0 :     if (db == 0)
     618             :     {
     619           0 :       if (!ZX_equal1(res->lc))
     620             :       {
     621           0 :           res->lc  = FpXQ_powu(res->lc, da, T, p);
     622           0 :           res->res = FpXQ_mul(res->res, res->lc, T, p);
     623             :       }
     624             :     } else
     625           0 :       res->res = pol_0(get_FpX_var(T));
     626             :   }
     627           0 : }
     628             : 
     629             : static GEN
     630         275 : FpXQX_halfres_basecase(GEN a, GEN b, GEN T, GEN p, GEN *pa, GEN *pb, struct FpXQX_res *res)
     631             : {
     632         275 :   pari_sp av=avma;
     633             :   GEN u,u1,v,v1, M;
     634         275 :   long vx = varn(a), vT = get_FpX_var(T), n = lgpol(a)>>1;
     635         275 :   u1 = v = pol_0(vx);
     636         275 :   u = v1 = pol_1(vx);
     637        2837 :   while (lgpol(b)>n)
     638             :   {
     639             :     GEN r, q;
     640        2562 :     q = FpXQX_divrem(a,b, T, p, &r);
     641        2562 :     if (res)
     642             :     {
     643           0 :       long da = degpol(a), db=degpol(b), dr = degpol(r);
     644           0 :       res->lc = to_ZX(gel(b,db+2),vT);
     645           0 :       if (dr >= n)
     646           0 :         FpXQX_halfres_update(da, db, dr, T, p, res);
     647             :       else
     648             :       {
     649           0 :         res->deg0 = da;
     650           0 :         res->deg1 = db;
     651             :       }
     652             :     }
     653        2562 :     a = b; b = r; swap(u,u1); swap(v,v1);
     654        2562 :     u1 = FpXX_sub(u1, FpXQX_mul(u, q, T, p), p);
     655        2562 :     v1 = FpXX_sub(v1, FpXQX_mul(v, q, T, p), p);
     656        2562 :     if (gc_needed(av,2))
     657             :     {
     658           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"FpXQX_halfgcd (d = %ld)",degpol(b));
     659           0 :       gerepileall(av,res ? 8: 6, &a,&b,&u1,&v1,&u,&v,&res->res,&res->lc);
     660             :     }
     661             :   }
     662         275 :   M = mkmat22(u,v,u1,v1); *pa = a; *pb = b;
     663         275 :   return gc_all(av, res ? 5: 3, &M, pa, pb, &res->res, &res->lc);
     664             : }
     665             : 
     666             : static GEN FpXQX_halfres_i(GEN x, GEN y, GEN T, GEN p, GEN *a, GEN *b, struct FpXQX_res *res);
     667             : 
     668             : static GEN
     669         215 : FpXQX_halfres_split(GEN x, GEN y, GEN T, GEN p, GEN *a, GEN *b, struct FpXQX_res *res)
     670             : {
     671         215 :   pari_sp av = avma;
     672             :   GEN Q, R, S, V1, V2;
     673             :   GEN x1, y1, r, q;
     674         215 :   long l = lgpol(x), n = l>>1, k, vT = get_FpX_var(T);
     675         215 :   if (lgpol(y) <= n)
     676           0 :     { *a = RgX_copy(x); *b = RgX_copy(y); return matid2_FpXQXM(varn(x)); }
     677         215 :   if (res)
     678             :   {
     679           0 :      res->lc = to_ZX(leading_coeff(y), vT);
     680           0 :      res->deg0 -= n;
     681           0 :      res->deg1 -= n;
     682           0 :      res->off += n;
     683             :   }
     684         215 :   R = FpXQX_halfres_i(FpXX_shift(x,-n),FpXX_shift(y,-n), T, p, a, b, res);
     685         215 :   if (res)
     686             :   {
     687           0 :     res->off -= n;
     688           0 :     res->deg0 += n;
     689           0 :     res->deg1 += n;
     690             :   }
     691         215 :   V1 = FpXQXM_FpXQX_mul2(R, Flxn_red(x,n), Flxn_red(y,n), T, p);
     692         215 :   x1 = FpXX_add(FpXX_shift(*a,n), gel(V1,1), p);
     693         215 :   y1 = FpXX_add(FpXX_shift(*b,n), gel(V1,2), p);
     694         215 :   if (lgpol(y1) <= n)
     695             :   {
     696          54 :     *a = x1; *b = y1;
     697          54 :     return gc_all(av, res ? 5: 3, &R, a, b, &res->res, &res->lc);
     698             :   }
     699         161 :   k = 2*n-degpol(y1);
     700         161 :   q = FpXQX_divrem(x1, y1, T, p, &r);
     701         161 :   if (res)
     702             :   {
     703           0 :     long dx1 = degpol(x1), dy1 = degpol(y1), dr = degpol(r);
     704           0 :     if (dy1 < degpol(y))
     705           0 :       FpXQX_halfres_update(res->deg0, res->deg1, dy1, T, p, res);
     706           0 :     res->lc = to_ZX(leading_coeff(y1), vT);
     707           0 :     res->deg0 = dx1;
     708           0 :     res->deg1 = dy1;
     709           0 :     if (dr >= n)
     710             :     {
     711           0 :       FpXQX_halfres_update(dx1, dy1, dr, T, p, res);
     712           0 :       res->deg0 = dy1;
     713           0 :       res->deg1 = dr;
     714             :     }
     715           0 :     res->deg0 -= k;
     716           0 :     res->deg1 -= k;
     717           0 :     res->off += k;
     718             :   }
     719         161 :   S = FpXQX_halfres_i(FpXX_shift(y1,-k), FpXX_shift(r,-k), T, p, a, b, res);
     720         161 :   if (res)
     721             :   {
     722           0 :     res->deg0 += k;
     723           0 :     res->deg1 += k;
     724           0 :     res->off -= k;
     725             :   }
     726         161 :   Q = FpXQXM_mul2(S,FpXQX_FpXQXM_qmul(q, R, T, p), T, p);
     727         161 :   V2 = FpXQXM_FpXQX_mul2(S, FpXXn_red(y1,k), FpXXn_red(r,k), T, p);
     728         161 :   *a = FpXX_add(FpXX_shift(*a,k), gel(V2,1), p);
     729         161 :   *b = FpXX_add(FpXX_shift(*b,k), gel(V2,2), p);
     730         161 :   return gc_all(av, res ? 5: 3, &Q, a, b, &res->res, &res->lc);
     731             : }
     732             : 
     733             : static GEN
     734         490 : FpXQX_halfres_i(GEN x, GEN y, GEN T, GEN p, GEN *a, GEN *b, struct FpXQX_res *res)
     735             : {
     736         490 :   if (lgpol(x) < FpXQX_HALFGCD_LIMIT)
     737         275 :     return FpXQX_halfres_basecase(x, y, T, p, a, b, res);
     738         215 :   return FpXQX_halfres_split(x, y, T, p, a, b, res);
     739             : }
     740             : 
     741             : static GEN
     742         114 : FpXQX_halfgcd_all_i(GEN x, GEN y, GEN T, GEN p, GEN *pa, GEN *pb)
     743             : {
     744             :   GEN a, b;
     745         114 :   GEN R = FpXQX_halfres_i(x, y, T, p, &a, &b, NULL);
     746         114 :   if (pa) *pa = a;
     747         114 :   if (pb) *pb = b;
     748         114 :   return R;
     749             : }
     750             : 
     751             : /* Return M in GL_2(Fp[X]/(T)[Y]) such that:
     752             : if [a',b']~=M*[a,b]~ then degpol(a')>= (lgpol(a)>>1) >degpol(b')
     753             : */
     754             : 
     755             : GEN
     756         114 : FpXQX_halfgcd_all(GEN x, GEN y, GEN T, GEN p, GEN *a, GEN *b)
     757             : {
     758         114 :   pari_sp av = avma;
     759             :   GEN R,q,r;
     760         114 :   if (lgefint(p)==3)
     761             :   {
     762           0 :     ulong pp = to_FlxqX(x, y, T, p, &x, &y, &T);
     763           0 :     R = FlxXM_to_ZXXM(FlxqX_halfgcd(x, y, T, pp));
     764           0 :     if (a) *a = Flx_to_ZX(*a);
     765           0 :     if (b) *b = Flx_to_ZX(*b);
     766           0 :     return !a && b ? gc_all(av, 2, &R, b): gc_all(av, 1+!!a+!!b, &R, a, b);
     767             :   }
     768         114 :   if (!signe(x))
     769             :   {
     770           0 :     if (a) *a = RgX_copy(y);
     771           0 :     if (b) *b = RgX_copy(x);
     772           0 :     return matJ2_FpXQXM(varn(x));
     773             :   }
     774         114 :   if (degpol(y)<degpol(x)) return FpXQX_halfgcd_all_i(x, y, T, p, a, b);
     775          26 :   q = FpXQX_divrem(y, x, T, p, &r);
     776          26 :   R = FpXQX_halfgcd_all_i(x, r, T, p, a, b);
     777          26 :   gcoeff(R,1,1) = FpXX_sub(gcoeff(R,1,1),
     778          26 :                            FpXQX_mul(q, gcoeff(R,1,2), T, p), p);
     779          26 :   gcoeff(R,2,1) = FpXX_sub(gcoeff(R,2,1),
     780          26 :                            FpXQX_mul(q, gcoeff(R,2,2), T, p), p);
     781          26 :   return !a && b ? gc_all(av, 2, &R, b): gc_all(av, 1+!!a+!!b, &R, a, b);
     782             : }
     783             : 
     784             : GEN
     785          44 : FpXQX_halfgcd(GEN x, GEN y, GEN T, GEN p)
     786          44 : { return FpXQX_halfgcd_all(x, y, T, p, NULL, NULL); }
     787             : 
     788             : static GEN
     789        3896 : FpXQX_gcd_basecase(GEN a, GEN b, GEN T, GEN p)
     790             : {
     791        3896 :   pari_sp av = avma, av0=avma;
     792       39080 :   while (signe(b))
     793             :   {
     794             :     GEN c;
     795       35184 :     if (gc_needed(av0,2))
     796             :     {
     797           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"FpXQX_gcd (d = %ld)",degpol(b));
     798           0 :       gerepileall(av0,2, &a,&b);
     799             :     }
     800       35184 :     av = avma; c = FpXQX_rem(a, b, T, p); a=b; b=c;
     801             :   }
     802        3896 :   return gc_const(av, a);
     803             : }
     804             : 
     805             : GEN
     806       14545 : FpXQX_gcd(GEN x, GEN y, GEN T, GEN p)
     807             : {
     808       14545 :   pari_sp av = avma;
     809       14545 :   if (lgefint(p) == 3)
     810             :   {
     811             :     GEN Pl, Ql, Tl, U;
     812       10540 :     ulong pp = to_FlxqX(x, y, T, p, &Pl, &Ql, &Tl);
     813       10540 :     U  = FlxqX_gcd(Pl, Ql, Tl, pp);
     814       10540 :     return gerepileupto(av, FlxX_to_ZXX(U));
     815             :   }
     816        4005 :   x = FpXQX_red(x, T, p);
     817        4005 :   y = FpXQX_red(y, T, p);
     818        4005 :   if (!signe(x)) return gerepileupto(av, y);
     819        3959 :   while (lgpol(y)>=FpXQX_GCD_LIMIT)
     820             :   {
     821          63 :     if (lgpol(y)<=(lgpol(x)>>1))
     822             :     {
     823           0 :       GEN r = FpXQX_rem(x, y, T, p);
     824           0 :       x = y; y = r;
     825             :     }
     826          63 :     (void) FpXQX_halfgcd_all(x,y, T, p, &x, &y);
     827          63 :     if (gc_needed(av,2))
     828             :     {
     829           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"FpXQX_gcd (y = %ld)",degpol(y));
     830           0 :       gerepileall(av,2,&x,&y);
     831             :     }
     832             :   }
     833        3896 :   return gerepileupto(av, FpXQX_gcd_basecase(x, y, T, p));
     834             : }
     835             : 
     836             : static GEN
     837           0 : FpXQX_extgcd_basecase(GEN a, GEN b, GEN T, GEN p, GEN *ptu, GEN *ptv)
     838             : {
     839           0 :   pari_sp av=avma;
     840             :   GEN u,v,d,d1,v1;
     841           0 :   long vx = varn(a);
     842           0 :   d = a; d1 = b;
     843           0 :   v = pol_0(vx); v1 = pol_1(vx);
     844           0 :   while (signe(d1))
     845             :   {
     846           0 :     GEN r, q = FpXQX_divrem(d, d1, T, p, &r);
     847           0 :     v = FpXX_sub(v,FpXQX_mul(q,v1,T, p),p);
     848           0 :     u=v; v=v1; v1=u;
     849           0 :     u=r; d=d1; d1=u;
     850           0 :     if (gc_needed(av,2))
     851             :     {
     852           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"FpXQX_extgcd (d = %ld)",degpol(d));
     853           0 :       gerepileall(av,5, &d,&d1,&u,&v,&v1);
     854             :     }
     855             :   }
     856           0 :   if (ptu) *ptu = FpXQX_div(FpXX_sub(d,FpXQX_mul(b,v, T, p), p), a, T, p);
     857           0 :   *ptv = v; return d;
     858             : }
     859             : 
     860             : static GEN
     861           0 : FpXQX_extgcd_halfgcd(GEN x, GEN y, GEN T, GEN p, GEN *ptu, GEN *ptv)
     862             : {
     863             :   GEN u,v;
     864           0 :   GEN V = cgetg(expu(lgpol(y))+2,t_VEC);
     865           0 :   long i, n = 0, vs = varn(x);
     866           0 :   while (lgpol(y) >= FpXQX_EXTGCD_LIMIT)
     867             :   {
     868           0 :     if (lgpol(y)<=(lgpol(x)>>1))
     869             :     {
     870           0 :       GEN r, q = FpXQX_divrem(x, y, T, p, &r);
     871           0 :       x = y; y = r;
     872           0 :       gel(V,++n) = mkmat22(pol_0(vs),pol_1(vs),pol_1(vs),FpXX_neg(q,p));
     873             :     } else
     874           0 :       gel(V,++n) = FpXQX_halfgcd_all(x, y, T, p, &x, &y);
     875             :   }
     876           0 :   y = FpXQX_extgcd_basecase(x,y, T, p, &u,&v);
     877           0 :   for (i = n; i>1; i--)
     878             :   {
     879           0 :     GEN R = gel(V,i);
     880           0 :     GEN u1 = FpXQX_addmulmul(u, v, gcoeff(R,1,1), gcoeff(R,2,1), T, p);
     881           0 :     GEN v1 = FpXQX_addmulmul(u, v, gcoeff(R,1,2), gcoeff(R,2,2), T, p);
     882           0 :     u = u1; v = v1;
     883             :   }
     884             :   {
     885           0 :     GEN R = gel(V,1);
     886           0 :     if (ptu)
     887           0 :       *ptu = FpXQX_addmulmul(u, v, gcoeff(R,1,1), gcoeff(R,2,1), T, p);
     888           0 :     *ptv   = FpXQX_addmulmul(u, v, gcoeff(R,1,2), gcoeff(R,2,2), T, p);
     889             :   }
     890           0 :   return y;
     891             : }
     892             : 
     893             : /* x and y in Z[Y][X], return lift(gcd(x mod T,p, y mod T,p)). Set u and v st
     894             :  * ux + vy = gcd (mod T,p) */
     895             : GEN
     896      133490 : FpXQX_extgcd(GEN x, GEN y, GEN T, GEN p, GEN *ptu, GEN *ptv)
     897             : {
     898      133490 :   pari_sp av = avma;
     899             :   GEN d;
     900      133490 :   if (lgefint(p) == 3)
     901             :   {
     902             :     GEN Pl, Ql, Tl, Dl;
     903      133490 :     ulong pp = to_FlxqX(x, y, T, p, &Pl, &Ql, &Tl);
     904      133489 :     Dl = FlxqX_extgcd(Pl, Ql, Tl, pp, ptu, ptv);
     905      133494 :     if (ptu) *ptu = FlxX_to_ZXX(*ptu);
     906      133490 :     *ptv = FlxX_to_ZXX(*ptv);
     907      133490 :     d = FlxX_to_ZXX(Dl);
     908             :   }
     909             :   else
     910             :   {
     911           0 :     x = FpXQX_red(x, T, p);
     912           0 :     y = FpXQX_red(y, T, p);
     913           0 :     if (lgpol(y)>=FpXQX_EXTGCD_LIMIT)
     914           0 :       d = FpXQX_extgcd_halfgcd(x, y, T, p, ptu, ptv);
     915             :     else
     916           0 :       d = FpXQX_extgcd_basecase(x, y, T, p, ptu, ptv);
     917             :   }
     918      133485 :   return gc_all(av, ptu?3:2, &d, ptv, ptu);
     919             : }
     920             : 
     921             : static GEN
     922           0 : FpXQX_halfres(GEN x, GEN y, GEN T, GEN p, GEN *a, GEN *b, GEN *r)
     923             : {
     924             :   struct FpXQX_res res;
     925             :   GEN V;
     926           0 :   long dB, vT = get_FpX_var(T);
     927             : 
     928           0 :   res.res  = *r;
     929           0 :   res.lc   = to_ZX(leading_coeff(y),vT);
     930           0 :   res.deg0 = degpol(x);
     931           0 :   res.deg1 = degpol(y);
     932           0 :   res.off = 0;
     933           0 :   V = FpXQX_halfres_i(x, y, T, p, a, b, &res);
     934           0 :   dB = degpol(*b);
     935           0 :   if (dB < degpol(y))
     936           0 :     FpXQX_halfres_update(res.deg0, res.deg1, dB, T, p, &res);
     937           0 :   *r = res.res;
     938           0 :   return V;
     939             : }
     940             : 
     941             : /* Res(A,B) = Res(B,R) * lc(B)^(a-r) * (-1)^(ab), with R=A%B, a=deg(A) ...*/
     942             : static GEN
     943          28 : FpXQX_resultant_basecase(GEN a, GEN b, GEN T, GEN p)
     944             : {
     945          28 :   pari_sp av = avma;
     946          28 :   long vT = get_FpX_var(T), da,db,dc;
     947          28 :   GEN c,lb, res = pol_1(vT);
     948             : 
     949          28 :   if (!signe(a) || !signe(b)) return pol_0(vT);
     950             : 
     951          28 :   da = degpol(a);
     952          28 :   db = degpol(b);
     953          28 :   if (db > da)
     954             :   {
     955           0 :     swapspec(a,b, da,db);
     956           0 :     if (both_odd(da,db)) res = FpX_neg(res, p);
     957             :   }
     958          28 :   if (!da) return pol_1(vT); /* = res * a[2] ^ db, since 0 <= db <= da = 0 */
     959          98 :   while (db)
     960             :   {
     961          70 :     lb = to_ZX(gel(b,db+2),vT);
     962          70 :     c = FpXQX_rem(a,b, T,p);
     963          70 :     a = b; b = c; dc = degpol(c);
     964          70 :     if (dc < 0) { set_avma(av); return pol_0(vT); }
     965             : 
     966          70 :     if (both_odd(da,db)) res = FpX_neg(res, p);
     967          70 :     if (!ZX_equal1(lb)) res = FpXQ_mul(res, FpXQ_powu(lb, da - dc, T, p), T, p);
     968          70 :     if (gc_needed(av,2))
     969             :     {
     970           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"FpXQX_resultant (da = %ld)",da);
     971           0 :       gerepileall(av,3, &a,&b,&res);
     972             :     }
     973          70 :     da = db; /* = degpol(a) */
     974          70 :     db = dc; /* = degpol(b) */
     975             :   }
     976          28 :   res = FpXQ_mul(res, FpXQ_powu(gel(b,2), da, T, p), T, p);
     977          28 :   return gerepileupto(av, res);
     978             : }
     979             : 
     980             : GEN
     981          63 : FpXQX_resultant(GEN x, GEN y, GEN T, GEN p)
     982             : {
     983          63 :   pari_sp av = avma;
     984          63 :   long dx, dy, vT = get_FpX_var(T);
     985          63 :   GEN res = pol_1(vT);
     986          63 :   if (!signe(x) || !signe(y)) return pol_0(vT);
     987          63 :   if (lgefint(p) == 3)
     988             :   {
     989          35 :     pari_sp av = avma;
     990             :     GEN Pl, Ql, Tl, R;
     991          35 :     ulong pp = to_FlxqX(x, y, T, p, &Pl, &Ql, &Tl);
     992          35 :     R = FlxqX_resultant(Pl, Ql, Tl, pp);
     993          35 :     return gerepileupto(av, Flx_to_ZX(R));
     994             :   }
     995             : 
     996          28 :   dx = degpol(x); dy = degpol(y);
     997          28 :   if (dx < dy)
     998             :   {
     999          14 :     swap(x,y);
    1000          14 :     if (both_odd(dx, dy))
    1001           0 :       res = Fp_neg(res, p);
    1002             :   }
    1003          28 :   while (lgpol(y) >= FpXQX_GCD_LIMIT)
    1004             :   {
    1005           0 :     if (lgpol(y)<=(lgpol(x)>>1))
    1006             :     {
    1007           0 :       GEN r = FpXQX_rem(x, y, T, p);
    1008           0 :       long dx = degpol(x), dy = degpol(y), dr = degpol(r);
    1009           0 :       GEN ly = FpX_red(gel(y,dy+2),p);
    1010           0 :       if (!ZX_equal1(ly)) res = FpXQ_mul(res, FpXQ_powu(ly, dx - dr, T, p), T, p);
    1011           0 :       if (both_odd(dx, dy))
    1012           0 :         res = Fp_neg(res, p);
    1013           0 :       x = y; y = r;
    1014             :     }
    1015           0 :     (void) FpXQX_halfres(x, y, T, p, &x, &y, &res);
    1016           0 :     if (gc_needed(av,2))
    1017             :     {
    1018           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"FpXQX_resultant (y = %ld)",degpol(y));
    1019           0 :       gerepileall(av,3,&x,&y,&res);
    1020             :     }
    1021             :   }
    1022          28 :   return gerepileupto(av, FpXQ_mul(res, FpXQX_resultant_basecase(x, y, T, p), T, p));
    1023             : }
    1024             : 
    1025             : /* disc P = (-1)^(n(n-1)/2) lc(P)^(n - deg P' - 2) Res(P,P'), n = deg P */
    1026             : GEN
    1027          35 : FpXQX_disc(GEN P, GEN T, GEN p)
    1028             : {
    1029          35 :   pari_sp av = avma;
    1030          35 :   GEN L, dP = FpXX_deriv(P, p), D = FpXQX_resultant(P, dP, T, p);
    1031             :   long dd;
    1032          35 :   if (!signe(D)) return pol_0(get_FpX_var(T));
    1033          35 :   dd = degpol(P) - 2 - degpol(dP); /* >= -1; > -1 iff p | deg(P) */
    1034          35 :   L = leading_coeff(P);
    1035          35 :   if (dd && !gequal1(L))
    1036           0 :     D = (dd == -1)? FpXQ_div(D,L,T,p): FpXQ_mul(D, FpXQ_powu(L, dd, T, p), T, p);
    1037          35 :   if (degpol(P) & 2) D = FpX_neg(D, p);
    1038          35 :   return gerepileupto(av, D);
    1039             : }
    1040             : 
    1041             : GEN
    1042         396 : FpXQX_dotproduct(GEN x, GEN y, GEN T, GEN p)
    1043             : {
    1044         396 :   long i, l = minss(lg(x), lg(y));
    1045             :   pari_sp av;
    1046             :   GEN c;
    1047         396 :   if (l == 2) return gen_0;
    1048         396 :   av = avma; c = gmul(gel(x,2),gel(y,2));
    1049        1642 :   for (i=3; i<l; i++) c = gadd(c, gmul(gel(x,i),gel(y,i)));
    1050         396 :   return gerepileupto(av, Fq_red(c,T,p));
    1051             : }
    1052             : 
    1053             : /***********************************************************************/
    1054             : /**                                                                   **/
    1055             : /**                       Barrett reduction                           **/
    1056             : /**                                                                   **/
    1057             : /***********************************************************************/
    1058             : 
    1059             : /* Return new lgpol */
    1060             : static long
    1061      316761 : ZXX_lgrenormalizespec(GEN x, long lx)
    1062             : {
    1063             :   long i;
    1064      317271 :   for (i = lx-1; i>=0; i--)
    1065      317271 :     if (signe(gel(x,i))) break;
    1066      316761 :   return i+1;
    1067             : }
    1068             : 
    1069             : static GEN
    1070        2918 : FpXQX_invBarrett_basecase(GEN S, GEN T, GEN p)
    1071             : {
    1072        2918 :   long i, l=lg(S)-1, lr = l-1, k;
    1073        2918 :   GEN r=cgetg(lr, t_POL); r[1]=S[1];
    1074        2918 :   gel(r,2) = gen_1;
    1075       24894 :   for (i=3; i<lr; i++)
    1076             :   {
    1077       21976 :     pari_sp av = avma;
    1078       21976 :     GEN u = gel(S,l-i+2);
    1079      221103 :     for (k=3; k<i; k++)
    1080      199127 :       u = Fq_add(u, Fq_mul(gel(S,l-i+k), gel(r,k), NULL, p), NULL, p);
    1081       21976 :     gel(r,i) = gerepileupto(av, Fq_red(Fq_neg(u, NULL, p), T, p));
    1082             :   }
    1083        2918 :   return FpXQX_renormalize(r,lr);
    1084             : }
    1085             : 
    1086             : INLINE GEN
    1087      306360 : FpXX_recipspec(GEN x, long l, long n)
    1088             : {
    1089      306360 :   return RgX_recipspec_shallow(x, l, n);
    1090             : }
    1091             : 
    1092             : static GEN
    1093         537 : FpXQX_invBarrett_Newton(GEN S, GEN T, GEN p)
    1094             : {
    1095         537 :   pari_sp av = avma;
    1096         537 :   long nold, lx, lz, lq, l = degpol(S), i, lQ;
    1097         537 :   GEN q, y, z, x = cgetg(l+2, t_POL) + 2;
    1098         537 :   ulong mask = quadratic_prec_mask(l-2); /* assume l > 2 */
    1099       38210 :   for (i=0;i<l;i++) gel(x,i) = gen_0;
    1100         537 :   q = RgX_recipspec_shallow(S+2,l+1,l+1); lQ = lgpol(q); q+=2;
    1101             :   /* We work on _spec_ FpX's, all the l[xzq] below are lgpol's */
    1102             : 
    1103             :   /* initialize */
    1104         537 :   gel(x,0) = Fq_inv(gel(q,0), T, p);
    1105         537 :   if (lQ>1) gel(q,1) = Fq_red(gel(q,1), T, p);
    1106         537 :   if (lQ>1 && signe(gel(q,1)))
    1107         537 :   {
    1108         537 :     GEN u = gel(q, 1);
    1109         537 :     if (!gequal1(gel(x,0))) u = Fq_mul(u, Fq_sqr(gel(x,0), T, p), T, p);
    1110         537 :     gel(x,1) = Fq_neg(u, T, p); lx = 2;
    1111             :   }
    1112             :   else
    1113           0 :     lx = 1;
    1114         537 :   nold = 1;
    1115        4004 :   for (; mask > 1; )
    1116             :   { /* set x -= x(x*q - 1) + O(t^(nnew + 1)), knowing x*q = 1 + O(t^(nold+1)) */
    1117        3467 :     long i, lnew, nnew = nold << 1;
    1118             : 
    1119        3467 :     if (mask & 1) nnew--;
    1120        3467 :     mask >>= 1;
    1121             : 
    1122        3467 :     lnew = nnew + 1;
    1123        3467 :     lq = ZXX_lgrenormalizespec(q, minss(lQ,lnew));
    1124        3467 :     z = FpXQX_mulspec(x, q, T, p, lx, lq); /* FIXME: high product */
    1125        3467 :     lz = lgpol(z); if (lz > lnew) lz = lnew;
    1126        3467 :     z += 2;
    1127             :     /* subtract 1 [=>first nold words are 0]: renormalize so that z(0) != 0 */
    1128        6934 :     for (i = nold; i < lz; i++) if (signe(gel(z,i))) break;
    1129        3467 :     nold = nnew;
    1130        3467 :     if (i >= lz) continue; /* z-1 = 0(t^(nnew + 1)) */
    1131             : 
    1132             :     /* z + i represents (x*q - 1) / t^i */
    1133        3467 :     lz = ZXX_lgrenormalizespec (z+i, lz-i);
    1134        3467 :     z = FpXQX_mulspec(x, z+i, T, p, lx, lz); /* FIXME: low product */
    1135        3467 :     lz = lgpol(z); z += 2;
    1136        3467 :     if (lz > lnew-i) lz = ZXX_lgrenormalizespec(z, lnew-i);
    1137             : 
    1138        3467 :     lx = lz+ i;
    1139        3467 :     y  = x + i; /* x -= z * t^i, in place */
    1140       39529 :     for (i = 0; i < lz; i++) gel(y,i) = Fq_neg(gel(z,i), T, p);
    1141             :   }
    1142         537 :   x -= 2; setlg(x, lx + 2); x[1] = S[1];
    1143         537 :   return gerepilecopy(av, x);
    1144             : }
    1145             : 
    1146             : GEN
    1147        3476 : FpXQX_invBarrett(GEN S, GEN T, GEN p)
    1148             : {
    1149        3476 :   pari_sp ltop = avma;
    1150        3476 :   long l = lg(S);
    1151             :   GEN r;
    1152        3476 :   if (l<5) return pol_0(varn(S));
    1153        3455 :   if (l<=FpXQX_INVBARRETT_LIMIT)
    1154             :   {
    1155        2918 :     GEN c = gel(S,l-1), ci=gen_1;
    1156        2918 :     if (!gequal1(c))
    1157             :     {
    1158        1870 :       ci = Fq_inv(c, T, p);
    1159        1870 :       S = FqX_Fq_mul(S, ci, T, p);
    1160        1870 :       r = FpXQX_invBarrett_basecase(S, T, p);
    1161        1870 :       r = FqX_Fq_mul(r, ci, T, p);
    1162             :     } else
    1163        1048 :       r = FpXQX_invBarrett_basecase(S, T, p);
    1164             :   }
    1165             :   else
    1166         537 :     r = FpXQX_invBarrett_Newton(S, T, p);
    1167        3455 :   return gerepileupto(ltop, r);
    1168             : }
    1169             : 
    1170             : GEN
    1171       12685 : FpXQX_get_red(GEN S, GEN T, GEN p)
    1172             : {
    1173       12685 :   if (typ(S)==t_POL && lg(S)>FpXQX_BARRETT_LIMIT)
    1174        1174 :     retmkvec2(FpXQX_invBarrett(S,T,p),S);
    1175       11511 :   return S;
    1176             : }
    1177             : 
    1178             : /* Compute x mod S where 2 <= degpol(S) <= l+1 <= 2*(degpol(S)-1)
    1179             :  * and mg is the Barrett inverse of S. */
    1180             : static GEN
    1181      153180 : FpXQX_divrem_Barrettspec(GEN x, long l, GEN mg, GEN S, GEN T, GEN p, GEN *pr)
    1182             : {
    1183             :   GEN q, r;
    1184      153180 :   long lt = degpol(S); /*We discard the leading term*/
    1185             :   long ld, lm, lT, lmg;
    1186      153180 :   ld = l-lt;
    1187      153180 :   lm = minss(ld, lgpol(mg));
    1188      153180 :   lT  = ZXX_lgrenormalizespec(S+2,lt);
    1189      153180 :   lmg = ZXX_lgrenormalizespec(mg+2,lm);
    1190      153180 :   q = FpXX_recipspec(x+lt,ld,ld);                 /* q = rec(x)     lq<=ld*/
    1191      153180 :   q = FpXQX_mulspec(q+2,mg+2,T,p,lgpol(q),lmg);    /* q = rec(x) * mg lq<=ld+lm*/
    1192      153180 :   q = FpXX_recipspec(q+2,minss(ld,lgpol(q)),ld);  /* q = rec (rec(x) * mg) lq<=ld*/
    1193      153180 :   if (!pr) return q;
    1194      152221 :   r = FpXQX_mulspec(q+2,S+2,T,p,lgpol(q),lT);      /* r = q*pol        lr<=ld+lt*/
    1195      152221 :   r = FpXX_subspec(x,r+2,p,lt,minss(lt,lgpol(r))); /* r = x - r   lr<=lt */
    1196      152221 :   if (pr == ONLY_REM) return r;
    1197       68338 :   *pr = r; return q;
    1198             : }
    1199             : 
    1200             : static GEN
    1201       85254 : FpXQX_divrem_Barrett(GEN x, GEN mg, GEN S, GEN T, GEN p, GEN *pr)
    1202             : {
    1203       85254 :   GEN q = NULL, r = FpXQX_red(x, T, p);
    1204       85254 :   long l = lgpol(r), lt = degpol(S), lm = 2*lt-1, v = varn(S);
    1205             :   long i;
    1206       85254 :   if (l <= lt)
    1207             :   {
    1208           0 :     if (pr == ONLY_REM) return r;
    1209           0 :     if (pr == ONLY_DIVIDES) return signe(r)? NULL: pol_0(v);
    1210           0 :     if (pr) *pr = r;
    1211           0 :     return pol_0(v);
    1212             :   }
    1213       85254 :   if (lt <= 1)
    1214          21 :     return FpXQX_divrem_basecase(r,S,T,p,pr);
    1215       85233 :   if (pr != ONLY_REM && l>lm)
    1216             :   {
    1217        1350 :     q = cgetg(l-lt+2, t_POL); q[1] = S[1];
    1218       87536 :     for (i=0;i<l-lt;i++) gel(q+2,i) = gen_0;
    1219             :   }
    1220      153261 :   while (l>lm)
    1221             :   {
    1222       68028 :     GEN zr, zq = FpXQX_divrem_Barrettspec(r+2+l-lm,lm,mg,S,T,p,&zr);
    1223       68028 :     long lz = lgpol(zr);
    1224       68028 :     if (pr != ONLY_REM)
    1225             :     {
    1226       34904 :       long lq = lgpol(zq);
    1227      113973 :       for(i=0; i<lq; i++) gel(q+2+l-lm,i) = gel(zq,2+i);
    1228             :     }
    1229      284012 :     for(i=0; i<lz; i++) gel(r+2+l-lm,i) = gel(zr,2+i);
    1230       68028 :     l = l-lm+lz;
    1231             :   }
    1232       85233 :   if (pr == ONLY_REM)
    1233             :   {
    1234       83883 :     if (l > lt)
    1235       83883 :       r = FpXQX_divrem_Barrettspec(r+2,l,mg,S,T,p,ONLY_REM);
    1236             :     else
    1237           0 :       r = FpXQX_renormalize(r, l+2);
    1238       83883 :     setvarn(r, v); return r;
    1239             :   }
    1240        1350 :   if (l > lt)
    1241             :   {
    1242        1269 :     GEN zq = FpXQX_divrem_Barrettspec(r+2,l,mg,S,T,p,pr ? &r: NULL);
    1243        1269 :     if (!q) q = zq;
    1244             :     else
    1245             :     {
    1246        1269 :       long lq = lgpol(zq);
    1247        8050 :       for(i=0; i<lq; i++) gel(q+2,i) = gel(zq,2+i);
    1248             :     }
    1249             :   }
    1250          81 :   else if (pr)
    1251          81 :     r = FpX_renormalize(r, l+2);
    1252        1350 :   setvarn(q, v); q = FpXQX_renormalize(q, lg(q));
    1253        1350 :   if (pr == ONLY_DIVIDES) return signe(r)? NULL: q;
    1254        1350 :   if (pr) { setvarn(r, v); *pr = r; }
    1255        1350 :   return q;
    1256             : }
    1257             : 
    1258             : GEN
    1259     1069033 : FpXQX_divrem(GEN x, GEN S, GEN T, GEN p, GEN *pr)
    1260             : {
    1261             :   GEN B, y;
    1262             :   long dy, dx, d;
    1263     1069033 :   if (pr==ONLY_REM) return FpXQX_rem(x, S, T, p);
    1264     1069033 :   y = get_FpXQX_red(S, &B);
    1265     1069025 :   dy = degpol(y); dx = degpol(x); d = dx-dy;
    1266     1069007 :   if (lgefint(p) == 3)
    1267             :   {
    1268             :     GEN a, b, t, z;
    1269     1050367 :     pari_sp tetpil, av = avma;
    1270     1050367 :     ulong pp = to_FlxqX(x, y, T, p, &a, &b, &t);
    1271     1050305 :     z = FlxqX_divrem(a, b, t, pp, pr);
    1272     1050363 :     if (pr == ONLY_DIVIDES && !z) return gc_NULL(av);
    1273     1050363 :     tetpil=avma;
    1274     1050363 :     z = FlxX_to_ZXX(z);
    1275     1050218 :     if (pr && pr != ONLY_DIVIDES && pr != ONLY_REM)
    1276     1016254 :       *pr = FlxX_to_ZXX(*pr);
    1277       33964 :     else return gerepile(av, tetpil, z);
    1278     1016240 :     gerepileallsp(av,tetpil,2, pr, &z);
    1279     1016453 :     return z;
    1280             :   }
    1281       18640 :   if (!B && d+3 < FpXQX_DIVREM_BARRETT_LIMIT)
    1282       17269 :     return FpXQX_divrem_basecase(x,y,T,p,pr);
    1283             :   else
    1284             :   {
    1285        1371 :     pari_sp av=avma;
    1286        1371 :     GEN mg = B? B: FpXQX_invBarrett(y, T, p);
    1287        1371 :     GEN q = FpXQX_divrem_Barrett(x,mg,y,T,p,pr);
    1288        1371 :     if (!q) return gc_NULL(av);
    1289        1371 :     if (!pr || pr==ONLY_DIVIDES) return gerepilecopy(av, q);
    1290         398 :     return gc_all(av, 2, &q, pr);
    1291             :   }
    1292             : }
    1293             : 
    1294             : GEN
    1295      263338 : FpXQX_rem(GEN x, GEN S, GEN T, GEN p)
    1296             : {
    1297      263338 :   GEN B, y = get_FpXQX_red(S, &B);
    1298      263338 :   long dy = degpol(y), dx = degpol(x), d = dx-dy;
    1299      263338 :   if (d < 0) return FpXQX_red(x, T, p);
    1300      242723 :   if (lgefint(p) == 3)
    1301             :   {
    1302        2736 :     pari_sp av = avma;
    1303             :     GEN a, b, t, z;
    1304        2736 :     ulong pp = to_FlxqX(x, y, T, p, &a, &b, &t);
    1305        2736 :     z = FlxqX_rem(a, b, t, pp);
    1306        2736 :     z = FlxX_to_ZXX(z);
    1307        2736 :     return gerepileupto(av, z);
    1308             :   }
    1309      239987 :   if (!B && d+3 < FpXQX_REM_BARRETT_LIMIT)
    1310      156104 :     return FpXQX_divrem_basecase(x,y, T, p, ONLY_REM);
    1311             :   else
    1312             :   {
    1313       83883 :     pari_sp av=avma;
    1314       83883 :     GEN mg = B? B: FpXQX_invBarrett(y, T, p);
    1315       83883 :     GEN r = FpXQX_divrem_Barrett(x, mg, y, T, p, ONLY_REM);
    1316       83883 :     return gerepileupto(av, r);
    1317             :   }
    1318             : }
    1319             : 
    1320             : /* x + y*z mod p */
    1321             : INLINE GEN
    1322       36478 : Fq_addmul(GEN x, GEN y, GEN z, GEN T, GEN p)
    1323             : {
    1324             :   pari_sp av;
    1325       36478 :   if (!signe(y) || !signe(z)) return Fq_red(x, T, p);
    1326       36478 :   if (!signe(x)) return Fq_mul(z,y, T, p);
    1327       36478 :   av = avma;
    1328       36478 :   return gerepileupto(av, Fq_add(x, Fq_mul(y, z, T, p), T, p));
    1329             : }
    1330             : 
    1331             : GEN
    1332       71939 : FpXQX_div_by_X_x(GEN a, GEN x, GEN T, GEN p, GEN *r)
    1333             : {
    1334       71939 :   long l = lg(a), i;
    1335             :   GEN z;
    1336       71939 :   if (lgefint(p)==3)
    1337             :   {
    1338       53700 :     pari_sp av = avma;
    1339             :     GEN ap, xp, t, z;
    1340       53700 :     ulong pp = to_FlxqX(a, NULL, T, p, &ap, NULL, &t);
    1341       53700 :     xp = ZX_to_Flx(to_ZX(x, get_FpX_var(T)), pp);
    1342       53700 :     z = FlxX_to_ZXX(FlxqX_div_by_X_x(ap, xp, t, pp, r));
    1343       53700 :     if (!r) return gerepileupto(av, z);
    1344           0 :     *r = Flx_to_ZX(*r);
    1345           0 :     return gc_all(av, 2, &z, r);
    1346             :   }
    1347       18239 :   if (l <= 3)
    1348             :   {
    1349           0 :     if (r) *r = l == 2? gen_0: gcopy(gel(a,2));
    1350           0 :     return pol_0(varn(a));
    1351             :   }
    1352       18239 :   l--; z = cgetg(l, t_POL); z[1] = a[1];
    1353       18239 :   gel(z, l-1) = gel(a,l);
    1354       54717 :   for (i=l-2; i>1; i--) /* z[i] = a[i+1] + x*z[i+1] */
    1355       36478 :     gel(z, i) = Fq_addmul(gel(a,i+1), x, gel(z,i+1), T, p);
    1356       18239 :   if (r) *r = Fq_addmul(gel(a,2), x, gel(z,2), T, p);
    1357       18239 :   return z;
    1358             : }
    1359             : 
    1360             : struct _FpXQXQ {
    1361             :   GEN T, S;
    1362             :   GEN p;
    1363             : };
    1364             : 
    1365      118528 : static GEN _FpXQX_mul(void *data, GEN a,GEN b)
    1366             : {
    1367      118528 :   struct _FpXQXQ *d=(struct _FpXQXQ*)data;
    1368      118528 :   return FpXQX_mul(a,b,d->T,d->p);
    1369             : }
    1370             : 
    1371        1729 : static GEN _FpXQX_sqr(void *data, GEN a)
    1372             : {
    1373        1729 :   struct _FpXQXQ *d=(struct _FpXQXQ*)data;
    1374        1729 :   return FpXQX_sqr(a, d->T, d->p);
    1375             : }
    1376             : 
    1377             : GEN
    1378          63 : FpXQX_powu(GEN x, ulong n, GEN T, GEN p)
    1379             : {
    1380             :   struct _FpXQXQ D;
    1381          63 :   if (n==0) return pol_1(varn(x));
    1382          63 :   D.T = T; D.p = p;
    1383          63 :   return gen_powu(x, n, (void *)&D, _FpXQX_sqr, _FpXQX_mul);
    1384             : }
    1385             : 
    1386             : GEN
    1387       16702 : FpXQXV_prod(GEN V, GEN T, GEN p)
    1388             : {
    1389       16702 :   if (lgefint(p) == 3)
    1390             :   {
    1391           0 :     pari_sp av = avma;
    1392           0 :     ulong pp = p[2];
    1393           0 :     GEN Tl = ZXT_to_FlxT(T, pp);
    1394           0 :     GEN Vl = ZXXV_to_FlxXV(V, pp, get_FpX_var(T));
    1395           0 :     Tl = FlxqXV_prod(Vl, Tl, pp);
    1396           0 :     return gerepileupto(av, FlxX_to_ZXX(Tl));
    1397             :   }
    1398             :   else
    1399             :   {
    1400             :     struct _FpXQXQ d;
    1401       16702 :     d.T=T; d.p=p;
    1402       16702 :     return gen_product(V, (void*)&d, &_FpXQX_mul);
    1403             :   }
    1404             : }
    1405             : 
    1406             : static GEN
    1407        9954 : _FpXQX_divrem(void * E, GEN x, GEN y, GEN *r)
    1408             : {
    1409        9954 :   struct _FpXQXQ *d = (struct _FpXQXQ *) E;
    1410        9954 :   return FpXQX_divrem(x, y, d->T, d->p, r);
    1411             : }
    1412             : 
    1413             : static GEN
    1414      121791 : _FpXQX_add(void * E, GEN x, GEN y)
    1415             : {
    1416      121791 :   struct _FpXQXQ *d = (struct _FpXQXQ *) E;
    1417      121791 :   return FpXX_add(x, y, d->p);
    1418             : }
    1419             : 
    1420             : static GEN
    1421        4638 : _FpXQX_sub(void * E, GEN x, GEN y) {
    1422        4638 :   struct _FpXQXQ *d = (struct _FpXQXQ*) E;
    1423        4638 :   return FpXX_sub(x,y, d->p);
    1424             : }
    1425             : 
    1426             : static struct bb_ring FpXQX_ring = { _FpXQX_add, _FpXQX_mul, _FpXQX_sqr };
    1427             : 
    1428             : GEN
    1429         623 : FpXQX_digits(GEN x, GEN B, GEN T, GEN p)
    1430             : {
    1431         623 :   long d = degpol(B), n = (lgpol(x)+d-1)/d;
    1432             :   struct _FpXQXQ D;
    1433         623 :   D.T = T; D.p = p;
    1434         623 :   return gen_digits(x, B, n, (void *)&D, &FpXQX_ring, _FpXQX_divrem);
    1435             : }
    1436             : 
    1437             : GEN
    1438         189 : FpXQXV_FpXQX_fromdigits(GEN x, GEN B, GEN T, GEN p)
    1439             : {
    1440             :   struct _FpXQXQ D;
    1441         189 :   D.T = T; D.p = p;
    1442         189 :   return gen_fromdigits(x,B,(void *)&D, &FpXQX_ring);
    1443             : }
    1444             : 
    1445             : /* Q an FpXY (t_POL with FpX coeffs), evaluate at X = x */
    1446             : GEN
    1447       51552 : FpXY_evalx(GEN Q, GEN x, GEN p)
    1448             : {
    1449       51552 :   long i, lb = lg(Q);
    1450             :   GEN z;
    1451       51552 :   z = cgetg(lb, t_POL); z[1] = Q[1];
    1452      445045 :   for (i=2; i<lb; i++)
    1453             :   {
    1454      393493 :     GEN q = gel(Q,i);
    1455      393493 :     gel(z,i) = typ(q) == t_INT? modii(q,p): FpX_eval(q, x, p);
    1456             :   }
    1457       51552 :   return FpX_renormalize(z, lb);
    1458             : }
    1459             : /* Q an FpXY, evaluate at Y = y */
    1460             : GEN
    1461       18799 : FpXY_evaly(GEN Q, GEN y, GEN p, long vx)
    1462             : {
    1463       18799 :   pari_sp av = avma;
    1464       18799 :   long i, lb = lg(Q);
    1465             :   GEN z;
    1466       18799 :   if (!signe(Q)) return pol_0(vx);
    1467       18771 :   if (lb == 3 || !signe(y)) {
    1468          84 :     z = gel(Q, 2);
    1469          84 :     return typ(z)==t_INT? scalar_ZX(z, vx): ZX_copy(z);
    1470             :   }
    1471       18687 :   z = gel(Q, lb-1);
    1472       18687 :   if (typ(z) == t_INT) z = scalar_ZX_shallow(z, vx);
    1473      242061 :   for (i=lb-2; i>=2; i--) z = Fq_add(gel(Q,i), FpX_Fp_mul(z, y, p), NULL, p);
    1474       18687 :   return gerepileupto(av, z);
    1475             : }
    1476             : /* Q an FpXY, evaluate at (X,Y) = (x,y) */
    1477             : GEN
    1478       13657 : FpXY_eval(GEN Q, GEN y, GEN x, GEN p)
    1479             : {
    1480       13657 :   pari_sp av = avma;
    1481       13657 :   return gerepileuptoint(av, FpX_eval(FpXY_evalx(Q, x, p), y, p));
    1482             : }
    1483             : 
    1484             : GEN
    1485        3852 : FpXY_FpXQV_evalx(GEN P, GEN x, GEN T, GEN p)
    1486             : {
    1487        3852 :   long i, lP = lg(P);
    1488        3852 :   GEN res = cgetg(lP,t_POL);
    1489        3852 :   res[1] = P[1];
    1490       61377 :   for(i=2; i<lP; i++)
    1491      115050 :     gel(res,i) = typ(gel(P,i))==t_INT? icopy(gel(P,i)):
    1492       57525 :                                        FpX_FpXQV_eval(gel(P,i), x, T, p);
    1493        3852 :   return FlxX_renormalize(res, lP);
    1494             : }
    1495             : 
    1496             : GEN
    1497         154 : FpXY_FpXQ_evalx(GEN P, GEN x, GEN T, GEN p)
    1498             : {
    1499         154 :   pari_sp av = avma;
    1500         154 :   long n = brent_kung_optpow(get_FpX_degree(T)-1,lgpol(P),1);
    1501         154 :   GEN xp = FpXQ_powers(x, n, T, p);
    1502         154 :   return gerepileupto(av, FpXY_FpXQV_evalx(P, xp, T, p));
    1503             : }
    1504             : 
    1505             : /*******************************************************************/
    1506             : /*                                                                 */
    1507             : /*                       (Fp[X]/T(X))[Y] / S(Y)                    */
    1508             : /*                                                                 */
    1509             : /*******************************************************************/
    1510             : 
    1511             : /*Preliminary implementation to speed up FpX_ffisom*/
    1512             : typedef struct {
    1513             :   GEN S, T, p;
    1514             : } FpXYQQ_muldata;
    1515             : 
    1516             : /* reduce x in Fp[X, Y] in the algebra Fp[X,Y]/ (S(X),T(Y)) */
    1517             : static GEN
    1518         476 : FpXYQQ_redswap(GEN x, GEN S, GEN T, GEN p)
    1519             : {
    1520         476 :   pari_sp ltop=avma;
    1521         476 :   long n = get_FpX_degree(S);
    1522         476 :   long m = get_FpX_degree(T);
    1523         476 :   long v = get_FpX_var(T);
    1524         476 :   GEN V = RgXY_swap(x,m,v);
    1525         476 :   V = FpXQX_red(V,S,p);
    1526         476 :   V = RgXY_swap(V,n,v);
    1527         476 :   return gerepilecopy(ltop,V);
    1528             : }
    1529             : static GEN
    1530         280 : FpXYQQ_sqr(void *data, GEN x)
    1531             : {
    1532         280 :   FpXYQQ_muldata *D = (FpXYQQ_muldata*)data;
    1533         280 :   return FpXYQQ_redswap(FpXQX_sqr(x, D->T, D->p),D->S,D->T,D->p);
    1534             : 
    1535             : }
    1536             : static GEN
    1537         196 : FpXYQQ_mul(void *data, GEN x, GEN y)
    1538             : {
    1539         196 :   FpXYQQ_muldata *D = (FpXYQQ_muldata*)data;
    1540         196 :   return FpXYQQ_redswap(FpXQX_mul(x,y, D->T, D->p),D->S,D->T,D->p);
    1541             : }
    1542             : 
    1543             : /* x in Z[X,Y], S in Z[X] over Fq = Z[Y]/(p,T); compute lift(x^n mod (S,T,p)) */
    1544             : GEN
    1545         182 : FpXYQQ_pow(GEN x, GEN n, GEN S, GEN T, GEN p)
    1546             : {
    1547         182 :   pari_sp av = avma;
    1548             :   FpXYQQ_muldata D;
    1549             :   GEN y;
    1550         182 :   if (lgefint(p) == 3)
    1551             :   {
    1552           0 :     ulong pp = to_FlxqX(x, NULL, T, p, &x, NULL, &T);
    1553           0 :     S = ZX_to_Flx(S, pp);
    1554           0 :     y = FlxX_to_ZXX( FlxYqq_pow(x, n, S, T, pp) );
    1555           0 :     y = gerepileupto(av, y);
    1556             :   }
    1557             :   else
    1558             :   {
    1559         182 :     D.S = S;
    1560         182 :     D.T = T;
    1561         182 :     D.p = p;
    1562         182 :     y = gen_pow(x, n, (void*)&D, &FpXYQQ_sqr, &FpXYQQ_mul);
    1563             :   }
    1564         182 :   return y;
    1565             : }
    1566             : 
    1567             : GEN
    1568       49194 : FpXQXQ_mul(GEN x, GEN y, GEN S, GEN T, GEN p) {
    1569       49194 :   return FpXQX_rem(FpXQX_mul(x, y, T, p), S, T, p);
    1570             : }
    1571             : 
    1572             : GEN
    1573      176117 : FpXQXQ_sqr(GEN x, GEN S, GEN T, GEN p) {
    1574      176117 :   return FpXQX_rem(FpXQX_sqr(x, T, p), S, T, p);
    1575             : }
    1576             : 
    1577             : /* Inverse of x in Z/pZ[X]/(pol) or NULL if inverse doesn't exist
    1578             :  * return lift(1 / (x mod (p,pol))) */
    1579             : GEN
    1580           0 : FpXQXQ_invsafe(GEN x, GEN S, GEN T, GEN p)
    1581             : {
    1582           0 :   GEN V, z = FpXQX_extgcd(get_FpXQX_mod(S), x, T, p, NULL, &V);
    1583           0 :   if (degpol(z)) return NULL;
    1584           0 :   z = gel(z,2);
    1585           0 :   z = typ(z)==t_INT ? Fp_invsafe(z,p) : FpXQ_invsafe(z,T,p);
    1586           0 :   if (!z) return NULL;
    1587           0 :   return typ(z)==t_INT ? FpXX_Fp_mul(V, z, p): FpXQX_FpXQ_mul(V, z, T, p);
    1588             : }
    1589             : 
    1590             : GEN
    1591           0 : FpXQXQ_inv(GEN x, GEN S, GEN T,GEN p)
    1592             : {
    1593           0 :   pari_sp av = avma;
    1594           0 :   GEN U = FpXQXQ_invsafe(x, S, T, p);
    1595           0 :   if (!U) pari_err_INV("FpXQXQ_inv",x);
    1596           0 :   return gerepileupto(av, U);
    1597             : }
    1598             : 
    1599             : GEN
    1600           0 : FpXQXQ_div(GEN x,GEN y,GEN S, GEN T,GEN p)
    1601             : {
    1602           0 :   pari_sp av = avma;
    1603           0 :   return gerepileupto(av, FpXQXQ_mul(x, FpXQXQ_inv(y,S,T,p),S,T,p));
    1604             : }
    1605             : 
    1606             : static GEN
    1607      125668 : _FpXQXQ_cmul(void *data, GEN P, long a, GEN x) {
    1608      125668 :   struct _FpXQXQ *d = (struct _FpXQXQ*) data;
    1609      125668 :   GEN y = gel(P,a+2);
    1610      251313 :   return typ(y)==t_INT ? FpXX_Fp_mul(x,y, d->p):
    1611      125645 :                          FpXX_FpX_mul(x,y,d->p);
    1612             : }
    1613             : static GEN
    1614       22956 : _FpXQXQ_red(void *data, GEN x) {
    1615       22956 :   struct _FpXQXQ *d = (struct _FpXQXQ*) data;
    1616       22956 :   return FpXQX_red(x, d->T, d->p);
    1617             : }
    1618             : static GEN
    1619       45926 : _FpXQXQ_mul(void *data, GEN x, GEN y) {
    1620       45926 :   struct _FpXQXQ *d = (struct _FpXQXQ*) data;
    1621       45926 :   return FpXQXQ_mul(x,y, d->S,d->T, d->p);
    1622             : }
    1623             : static GEN
    1624      176117 : _FpXQXQ_sqr(void *data, GEN x) {
    1625      176117 :   struct _FpXQXQ *d = (struct _FpXQXQ*) data;
    1626      176117 :   return FpXQXQ_sqr(x, d->S,d->T, d->p);
    1627             : }
    1628             : 
    1629             : static GEN
    1630       18963 : _FpXQXQ_one(void *data) {
    1631       18963 :   struct _FpXQXQ *d = (struct _FpXQXQ*) data;
    1632       18963 :   return pol_1(get_FpXQX_var(d->S));
    1633             : }
    1634             : 
    1635             : static GEN
    1636         132 : _FpXQXQ_zero(void *data) {
    1637         132 :   struct _FpXQXQ *d = (struct _FpXQXQ*) data;
    1638         132 :   return pol_0(get_FpXQX_var(d->S));
    1639             : }
    1640             : 
    1641             : static struct bb_algebra FpXQXQ_algebra = { _FpXQXQ_red, _FpXQX_add,
    1642             :        _FpXQX_sub, _FpXQXQ_mul, _FpXQXQ_sqr, _FpXQXQ_one, _FpXQXQ_zero };
    1643             : 
    1644             : const struct bb_algebra *
    1645         352 : get_FpXQXQ_algebra(void **E, GEN S, GEN T, GEN p)
    1646             : {
    1647         352 :   GEN z = new_chunk(sizeof(struct _FpXQXQ));
    1648         352 :   struct _FpXQXQ *e = (struct _FpXQXQ *) z;
    1649         352 :   e->T = FpX_get_red(T, p);
    1650         352 :   e->S = FpXQX_get_red(S, e->T, p);
    1651         352 :   e->p  = p; *E = (void*)e;
    1652         352 :   return &FpXQXQ_algebra;
    1653             : }
    1654             : 
    1655             : static struct bb_algebra FpXQX_algebra = { _FpXQXQ_red, _FpXQX_add,
    1656             :        _FpXQX_sub, _FpXQX_mul, _FpXQX_sqr, _FpXQXQ_one, _FpXQXQ_zero };
    1657             : 
    1658             : const struct bb_algebra *
    1659           0 : get_FpXQX_algebra(void **E, GEN T, GEN p, long v)
    1660             : {
    1661           0 :   GEN z = new_chunk(sizeof(struct _FpXQXQ));
    1662           0 :   struct _FpXQXQ *e = (struct _FpXQXQ *) z;
    1663           0 :   e->T = FpX_get_red(T, p);
    1664           0 :   e->S = pol_x(v);
    1665           0 :   e->p  = p; *E = (void*)e;
    1666           0 :   return &FpXQX_algebra;
    1667             : }
    1668             : 
    1669             : /* x over Fq, return lift(x^n) mod S */
    1670             : GEN
    1671        1908 : FpXQXQ_pow(GEN x, GEN n, GEN S, GEN T, GEN p)
    1672             : {
    1673        1908 :   pari_sp ltop = avma;
    1674             :   GEN y;
    1675             :   struct _FpXQXQ D;
    1676        1908 :   long s = signe(n);
    1677        1908 :   if (!s) return pol_1(varn(x));
    1678        1908 :   if (is_pm1(n)) /* +/- 1 */
    1679           0 :     return (s < 0)? FpXQXQ_inv(x,S,T,p): ZXX_copy(x);
    1680        1908 :   if (lgefint(p) == 3)
    1681             :   {
    1682          35 :     ulong pp = to_FlxqX(x, S, T, p, &x, &S, &T);
    1683          35 :     GEN z = FlxqXQ_pow(x, n, S, T, pp);
    1684          35 :     y = FlxX_to_ZXX(z);
    1685          35 :     return gerepileupto(ltop, y);
    1686             :   }
    1687             :   else
    1688             :   {
    1689        1873 :     T = FpX_get_red(T, p);
    1690        1873 :     S = FpXQX_get_red(S, T, p);
    1691        1873 :     D.S = S; D.T = T; D.p = p;
    1692        1873 :     if (s < 0) x = FpXQXQ_inv(x,S,T,p);
    1693        1873 :     y = gen_pow_i(x, n, (void*)&D,&_FpXQXQ_sqr,&_FpXQXQ_mul);
    1694        1873 :     return gerepilecopy(ltop, y);
    1695             :   }
    1696             : }
    1697             : 
    1698             : /* generates the list of powers of x of degree 0,1,2,...,l*/
    1699             : GEN
    1700        1965 : FpXQXQ_powers(GEN x, long l, GEN S, GEN T, GEN p)
    1701             : {
    1702             :   struct _FpXQXQ D;
    1703        1965 :   int use_sqr = 2*degpol(x) >= get_FpXQX_degree(S);
    1704        1965 :   T = FpX_get_red(T, p);
    1705        1965 :   S = FpXQX_get_red(S, T, p);
    1706        1965 :   D.S = S; D.T = T; D.p = p;
    1707        1965 :   return gen_powers(x, l, use_sqr, (void*)&D, &_FpXQXQ_sqr, &_FpXQXQ_mul,&_FpXQXQ_one);
    1708             : }
    1709             : 
    1710             : /* Let v a linear form, return the linear form z->v(tau*z)
    1711             :    that is, v*(M_tau) */
    1712             : 
    1713             : INLINE GEN
    1714         248 : FpXQX_recipspec(GEN x, long l, long n)
    1715             : {
    1716         248 :   return RgX_recipspec_shallow(x, l, n);
    1717             : }
    1718             : 
    1719             : static GEN
    1720          88 : FpXQXQ_transmul_init(GEN tau, GEN S, GEN T, GEN p)
    1721             : {
    1722             :   GEN bht;
    1723          88 :   GEN h, Sp = get_FpXQX_red(S, &h);
    1724          88 :   long n = degpol(Sp), vT = varn(Sp);
    1725          88 :   GEN ft = FpXQX_recipspec(Sp+2, n+1, n+1);
    1726          88 :   GEN bt = FpXQX_recipspec(tau+2, lgpol(tau), n);
    1727          88 :   setvarn(ft, vT); setvarn(bt, vT);
    1728          88 :   if (h)
    1729          16 :     bht = FpXQXn_mul(bt, h, n-1, T, p);
    1730             :   else
    1731             :   {
    1732          72 :     GEN bh = FpXQX_div(FpXX_shift(tau, n-1), S, T, p);
    1733          72 :     bht = FpXQX_recipspec(bh+2, lgpol(bh), n-1);
    1734          72 :     setvarn(bht, vT);
    1735             :   }
    1736          88 :   return mkvec3(bt, bht, ft);
    1737             : }
    1738             : 
    1739             : static GEN
    1740         200 : FpXQXQ_transmul(GEN tau, GEN a, long n, GEN T, GEN p)
    1741             : {
    1742         200 :   pari_sp ltop = avma;
    1743             :   GEN t1, t2, t3, vec;
    1744         200 :   GEN bt = gel(tau, 1), bht = gel(tau, 2), ft = gel(tau, 3);
    1745         200 :   if (signe(a)==0) return pol_0(varn(a));
    1746         200 :   t2 = FpXX_shift(FpXQX_mul(bt, a, T, p),1-n);
    1747         200 :   if (signe(bht)==0) return gerepilecopy(ltop, t2);
    1748         114 :   t1 = FpXX_shift(FpXQX_mul(ft, a, T, p),-n);
    1749         114 :   t3 = FpXQXn_mul(t1, bht, n-1, T, p);
    1750         114 :   vec = FpXX_sub(t2, FpXX_shift(t3, 1), p);
    1751         114 :   return gerepileupto(ltop, vec);
    1752             : }
    1753             : 
    1754             : static GEN
    1755          44 : polxn_FpXX(long n, long v, long vT)
    1756             : {
    1757          44 :   long i, a = n+2;
    1758          44 :   GEN p = cgetg(a+1, t_POL);
    1759          44 :   p[1] = evalsigne(1)|evalvarn(v);
    1760         440 :   for (i = 2; i < a; i++) gel(p,i) = pol_0(vT);
    1761          44 :   gel(p,a) = pol_1(vT); return p;
    1762             : }
    1763             : 
    1764             : GEN
    1765          44 : FpXQXQ_minpoly(GEN x, GEN S, GEN T, GEN p)
    1766             : {
    1767          44 :   pari_sp ltop = avma;
    1768             :   long vS, vT, n;
    1769             :   GEN v_x, g, tau;
    1770          44 :   vS = get_FpXQX_var(S);
    1771          44 :   vT = get_FpX_var(T);
    1772          44 :   n = get_FpXQX_degree(S);
    1773          44 :   g = pol_1(vS);
    1774          44 :   tau = pol_1(vS);
    1775          44 :   S = FpXQX_get_red(S, T, p);
    1776          44 :   v_x = FpXQXQ_powers(x, usqrt(2*n), S, T, p);
    1777          88 :   while(signe(tau) != 0)
    1778             :   {
    1779             :     long i, j, m, k1;
    1780             :     GEN M, v, tr;
    1781             :     GEN g_prime, c;
    1782          44 :     if (degpol(g) == n) { tau = pol_1(vS); g = pol_1(vS); }
    1783          44 :     v = random_FpXQX(n, vS, T, p);
    1784          44 :     tr = FpXQXQ_transmul_init(tau, S, T, p);
    1785          44 :     v = FpXQXQ_transmul(tr, v, n, T, p);
    1786          44 :     m = 2*(n-degpol(g));
    1787          44 :     k1 = usqrt(m);
    1788          44 :     tr = FpXQXQ_transmul_init(gel(v_x,k1+1), S, T, p);
    1789          44 :     c = cgetg(m+2,t_POL);
    1790          44 :     c[1] = evalsigne(1)|evalvarn(vS);
    1791         200 :     for (i=0; i<m; i+=k1)
    1792             :     {
    1793         156 :       long mj = minss(m-i, k1);
    1794         552 :       for (j=0; j<mj; j++)
    1795         396 :         gel(c,m+1-(i+j)) = FpXQX_dotproduct(v, gel(v_x,j+1), T, p);
    1796         156 :       v = FpXQXQ_transmul(tr, v, n, T, p);
    1797             :     }
    1798          44 :     c = FpXX_renormalize(c, m+2);
    1799             :     /* now c contains <v,x^i> , i = 0..m-1  */
    1800          44 :     M = FpXQX_halfgcd(polxn_FpXX(m, vS, vT), c, T, p);
    1801          44 :     g_prime = gmael(M, 2, 2);
    1802          44 :     if (degpol(g_prime) < 1) continue;
    1803          44 :     g = FpXQX_mul(g, g_prime, T, p);
    1804          44 :     tau = FpXQXQ_mul(tau, FpXQX_FpXQXQV_eval(g_prime, v_x, S, T, p), S, T, p);
    1805             :   }
    1806          44 :   g = FpXQX_normalize(g,T, p);
    1807          44 :   return gerepilecopy(ltop,g);
    1808             : }
    1809             : 
    1810             : GEN
    1811           0 : FpXQXQ_matrix_pow(GEN y, long n, long m, GEN S, GEN T, GEN p)
    1812             : {
    1813           0 :   return RgXV_to_RgM(FpXQXQ_powers(y,m-1,S,T,p),n);
    1814             : }
    1815             : 
    1816             : GEN
    1817        4191 : FpXQX_FpXQXQV_eval(GEN P, GEN V, GEN S, GEN T, GEN p)
    1818             : {
    1819             :   struct _FpXQXQ D;
    1820        4191 :   T = FpX_get_red(T, p);
    1821        4191 :   S = FpXQX_get_red(S, T, p);
    1822        4191 :   D.S=S; D.T=T; D.p=p;
    1823        4191 :   return gen_bkeval_powers(P, degpol(P), V, (void*)&D, &FpXQXQ_algebra,
    1824             :                                                    _FpXQXQ_cmul);
    1825             : }
    1826             : 
    1827             : GEN
    1828         855 : FpXQX_FpXQXQ_eval(GEN Q, GEN x, GEN S, GEN T, GEN p)
    1829             : {
    1830             :   struct _FpXQXQ D;
    1831         855 :   int use_sqr = 2*degpol(x) >= get_FpXQX_degree(S);
    1832         855 :   T = FpX_get_red(T, p);
    1833         855 :   S = FpXQX_get_red(S, T, p);
    1834         855 :   D.S=S; D.T=T; D.p=p;
    1835         855 :   return gen_bkeval(Q, degpol(Q), x, use_sqr, (void*)&D, &FpXQXQ_algebra,
    1836             :       _FpXQXQ_cmul);
    1837             : }
    1838             : 
    1839             : static GEN
    1840         507 : FpXQXQ_autpow_sqr(void * E, GEN x)
    1841             : {
    1842         507 :   struct _FpXQXQ *D = (struct _FpXQXQ *)E;
    1843         507 :   GEN S = D->S, T = D->T, p = D->p;
    1844         507 :   GEN phi = gel(x,1), S1 = gel(x,2);
    1845         507 :   long n = brent_kung_optpow(get_FpX_degree(T)-1,lgpol(S1)+1,1);
    1846         507 :   GEN V = FpXQ_powers(phi, n, T, p);
    1847         507 :   GEN phi2 = FpX_FpXQV_eval(phi, V, T, p);
    1848         507 :   GEN Sphi = FpXY_FpXQV_evalx(S1, V, T, p);
    1849         507 :   GEN S2 = FpXQX_FpXQXQ_eval(Sphi, S1, S, T, p);
    1850         507 :   return mkvec2(phi2, S2);
    1851             : }
    1852             : 
    1853             : static GEN
    1854         325 : FpXQXQ_autpow_mul(void * E, GEN x, GEN y)
    1855             : {
    1856         325 :   struct _FpXQXQ *D = (struct _FpXQXQ *)E;
    1857         325 :   GEN S = D->S, T = D->T, p = D->p;
    1858         325 :   GEN phi1 = gel(x,1), S1 = gel(x,2);
    1859         325 :   GEN phi2 = gel(y,1), S2 = gel(y,2);
    1860         325 :   long n = brent_kung_optpow(get_FpX_degree(T)-1, lgpol(S1)+1, 1);
    1861         325 :   GEN V = FpXQ_powers(phi2, n, T, p);
    1862         325 :   GEN phi3 = FpX_FpXQV_eval(phi1, V, T, p);
    1863         325 :   GEN Sphi = FpXY_FpXQV_evalx(S1, V, T, p);
    1864         325 :   GEN S3 = FpXQX_FpXQXQ_eval(Sphi, S2, S, T, p);
    1865         325 :   return mkvec2(phi3, S3);
    1866             : }
    1867             : 
    1868             : GEN
    1869         493 : FpXQXQ_autpow(GEN aut, long n, GEN S, GEN T, GEN p)
    1870             : {
    1871         493 :   pari_sp av = avma;
    1872             :   struct _FpXQXQ D;
    1873         493 :   T = FpX_get_red(T, p);
    1874         493 :   S = FpXQX_get_red(S, T, p);
    1875         493 :   D.S=S; D.T=T; D.p=p;
    1876         493 :   aut = gen_powu_i(aut,n,&D,FpXQXQ_autpow_sqr,FpXQXQ_autpow_mul);
    1877         493 :   return gerepilecopy(av, aut);
    1878             : }
    1879             : 
    1880             : static GEN
    1881           1 : FpXQXQ_auttrace_mul(void *E, GEN x, GEN y)
    1882             : {
    1883           1 :   struct _FpXQXQ *D = (struct _FpXQXQ *)E;
    1884           1 :   GEN S = D->S, T = D->T;
    1885           1 :   GEN p = D->p;
    1886           1 :   GEN S1 = gel(x,1), a1 = gel(x,2);
    1887           1 :   GEN S2 = gel(y,1), a2 = gel(y,2);
    1888           1 :   long n = brent_kung_optpow(maxss(degpol(S1),degpol(a1)),2,1);
    1889           1 :   GEN V = FpXQXQ_powers(S2, n, S, T, p);
    1890           1 :   GEN S3 = FpXQX_FpXQXQV_eval(S1, V, S, T, p);
    1891           1 :   GEN aS = FpXQX_FpXQXQV_eval(a1, V, S, T, p);
    1892           1 :   GEN a3 = FpXX_add(aS, a2, p);
    1893           1 :   return mkvec2(S3, a3);
    1894             : }
    1895             : 
    1896             : static GEN
    1897           1 : FpXQXQ_auttrace_sqr(void *E, GEN x)
    1898           1 : { return FpXQXQ_auttrace_mul(E, x, x); }
    1899             : 
    1900             : GEN
    1901           8 : FpXQXQ_auttrace(GEN aut, long n, GEN S, GEN T, GEN p)
    1902             : {
    1903           8 :   pari_sp av = avma;
    1904             :   struct _FpXQXQ D;
    1905           8 :   T = FpX_get_red(T, p);
    1906           8 :   S = FpXQX_get_red(S, T, p);
    1907           8 :   D.S=S; D.T=T; D.p=p;
    1908           8 :   aut = gen_powu_i(aut,n,&D,FpXQXQ_auttrace_sqr,FpXQXQ_auttrace_mul);
    1909           8 :   return gerepilecopy(av, aut);
    1910             : }
    1911             : 
    1912             : static GEN
    1913        1433 : FpXQXQ_autsum_mul(void *E, GEN x, GEN y)
    1914             : {
    1915        1433 :   struct _FpXQXQ *D = (struct _FpXQXQ *) E;
    1916        1433 :   GEN S = D->S, T = D->T, p = D->p;
    1917        1433 :   GEN phi1 = gel(x,1), S1 = gel(x,2), a1 = gel(x,3);
    1918        1433 :   GEN phi2 = gel(y,1), S2 = gel(y,2), a2 = gel(y,3);
    1919        1433 :   long n2 = brent_kung_optpow(get_FpX_degree(T)-1, lgpol(S1)+lgpol(a1)+1, 1);
    1920        1433 :   GEN V2 = FpXQ_powers(phi2, n2, T, p);
    1921        1433 :   GEN phi3 = FpX_FpXQV_eval(phi1, V2, T, p);
    1922        1433 :   GEN Sphi = FpXY_FpXQV_evalx(S1, V2, T, p);
    1923        1433 :   GEN aphi = FpXY_FpXQV_evalx(a1, V2, T, p);
    1924        1433 :   long n = brent_kung_optpow(maxss(degpol(Sphi),degpol(aphi)),2,1);
    1925        1433 :   GEN V = FpXQXQ_powers(S2, n, S, T, p);
    1926        1433 :   GEN S3 = FpXQX_FpXQXQV_eval(Sphi, V, S, T, p);
    1927        1433 :   GEN aS = FpXQX_FpXQXQV_eval(aphi, V, S, T, p);
    1928        1433 :   GEN a3 = FpXQXQ_mul(aS, a2, S, T, p);
    1929        1433 :   return mkvec3(phi3, S3, a3);
    1930             : }
    1931             : 
    1932             : static GEN
    1933        1275 : FpXQXQ_autsum_sqr(void * T, GEN x)
    1934        1275 : { return FpXQXQ_autsum_mul(T,x,x); }
    1935             : 
    1936             : GEN
    1937        1233 : FpXQXQ_autsum(GEN aut, long n, GEN S, GEN T, GEN p)
    1938             : {
    1939        1233 :   pari_sp av = avma;
    1940             :   struct _FpXQXQ D;
    1941        1233 :   T = FpX_get_red(T, p);
    1942        1233 :   S = FpXQX_get_red(S, T, p);
    1943        1233 :   D.S=S; D.T=T; D.p=p;
    1944        1233 :   aut = gen_powu_i(aut,n,&D,FpXQXQ_autsum_sqr,FpXQXQ_autsum_mul);
    1945        1233 :   return gerepilecopy(av, aut);
    1946             : }
    1947             : 
    1948             : GEN
    1949       44552 : FpXQXn_mul(GEN x, GEN y, long n, GEN T, GEN p)
    1950             : {
    1951       44552 :   pari_sp av = avma;
    1952             :   GEN z, kx, ky;
    1953             :   long d;
    1954       44552 :   if (ZXX_is_ZX(y) && ZXX_is_ZX(x))
    1955        6930 :     return FpXn_mul(x,y,n,p);
    1956       37622 :   d = get_FpX_degree(T);
    1957       37622 :   kx = RgXX_to_Kronecker(x, d);
    1958       37622 :   ky = RgXX_to_Kronecker(y, d);
    1959       37622 :   z = Kronecker_to_FpXQX(ZXn_mul(ky,kx,(2*d-1)*n), T, p);
    1960       37622 :   return gerepileupto(av, z);
    1961             : }
    1962             : 
    1963             : GEN
    1964           0 : FpXQXn_sqr(GEN x, long n, GEN T, GEN p)
    1965             : {
    1966           0 :   pari_sp av = avma;
    1967             :   GEN z, kx;
    1968             :   long d;
    1969           0 :   if (ZXX_is_ZX(x)) return ZXn_sqr(x, n);
    1970           0 :   d = get_FpX_degree(T);
    1971           0 :   kx= RgXX_to_Kronecker(x, d);
    1972           0 :   z = Kronecker_to_FpXQX(ZXn_sqr(kx, (2*d-1)*n), T, p);
    1973           0 :   return gerepileupto(av, z);
    1974             : }
    1975             : 
    1976             : /* (f*g) \/ x^n */
    1977             : static GEN
    1978        7392 : FpXQX_mulhigh_i(GEN f, GEN g, long n, GEN T, GEN p)
    1979             : {
    1980        7392 :   return FpXX_shift(FpXQX_mul(f,g,T, p),-n);
    1981             : }
    1982             : 
    1983             : static GEN
    1984        4690 : FpXQXn_mulhigh(GEN f, GEN g, long n2, long n, GEN T, GEN p)
    1985             : {
    1986        4690 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    1987        4690 :   return FpXX_add(FpXQX_mulhigh_i(fl, g, n2, T, p), FpXQXn_mul(fh, g, n - n2, T, p), p);
    1988             : }
    1989             : 
    1990             : /* Compute intformal(x^n*S)/x^(n+1) */
    1991             : static GEN
    1992         763 : FpXX_integXn(GEN x, long n, GEN p)
    1993             : {
    1994         763 :   long i, lx = lg(x);
    1995             :   GEN y;
    1996         763 :   if (lx == 2) return ZXX_copy(x);
    1997         763 :   y = cgetg(lx, t_POL); y[1] = x[1];
    1998        4317 :   for (i=2; i<lx; i++)
    1999             :   {
    2000        3554 :     ulong j = n+i-1;
    2001        3554 :     GEN xi = gel(x,i);
    2002        3554 :     if (!signe(xi))
    2003           0 :       gel(y,i) = gen_0;
    2004             :     else
    2005        3554 :       gel(y,i) = typ(xi)==t_INT ? Fp_divu(xi, j, p)
    2006        3554 :                                 : FpX_divu(xi, j, p);
    2007             :   }
    2008         763 :   return ZXX_renormalize(y, lx);;
    2009             : }
    2010             : 
    2011             : /* Compute intformal(x^n*S)/x^(n+1) */
    2012             : static GEN
    2013        2702 : ZlXX_integXn(GEN x, long n, GEN p, ulong pp)
    2014             : {
    2015        2702 :   long i, lx = lg(x);
    2016             :   GEN y;
    2017        2702 :   if (lx == 2) return ZXX_copy(x);
    2018        2576 :   if (!pp) return FpXX_integXn(x, n, p);
    2019        1813 :   y = cgetg(lx, t_POL); y[1] = x[1];
    2020        6883 :   for (i=2; i<lx; i++)
    2021             :   {
    2022        5070 :     GEN xi = gel(x,i);
    2023        5070 :     if (!signe(xi))
    2024          14 :       gel(y,i) = gen_0;
    2025             :     else
    2026             :     {
    2027             :       ulong j;
    2028        5056 :       long v = u_lvalrem(n+i-1, pp, &j);
    2029        5056 :       if (typ(xi)==t_INT)
    2030             :       {
    2031           0 :         if (v==0)
    2032           0 :           gel(y,i) = Fp_divu(xi, j, p);
    2033             :         else
    2034           0 :           gel(y,i) = Fp_divu(diviuexact(xi, upowuu(pp, v)), j, p);
    2035             :       } else
    2036             :       {
    2037        5056 :         if (v==0)
    2038        5056 :           gel(y,i) = FpX_divu(xi, j, p);
    2039             :         else
    2040           0 :           gel(y,i) = FpX_divu(ZX_divuexact(xi, upowuu(pp, v)), j, p);
    2041             :       }
    2042             :     }
    2043             :   }
    2044        1813 :   return ZXX_renormalize(y, lx);;
    2045             : }
    2046             : 
    2047             : GEN
    2048         714 : ZlXQXn_expint(GEN h, long e, GEN T, GEN p, ulong pp)
    2049             : {
    2050         714 :   pari_sp av = avma, av2;
    2051         714 :   long v = varn(h), n=1;
    2052         714 :   GEN f = pol_1(v), g = pol_1(v);
    2053         714 :   ulong mask = quadratic_prec_mask(e);
    2054         714 :   av2 = avma;
    2055        2702 :   for (;mask>1;)
    2056             :   {
    2057             :     GEN u, w;
    2058        2702 :     long n2 = n;
    2059        2702 :     n<<=1; if (mask & 1) n--;
    2060        2702 :     mask >>= 1;
    2061        2702 :     u = FpXQXn_mul(g, FpXQX_mulhigh_i(f, FpXXn_red(h, n2-1), n2-1, T, p), n-n2, T, p);
    2062        2702 :     u = FpXX_add(u, FpXX_shift(FpXXn_red(h, n-1), 1-n2), p);
    2063        2702 :     w = FpXQXn_mul(f, ZlXX_integXn(u, n2-1, p, pp), n-n2, T, p);
    2064        2702 :     f = FpXX_add(f, FpXX_shift(w, n2), p);
    2065        2702 :     if (mask<=1) break;
    2066        1988 :     u = FpXQXn_mul(g, FpXQXn_mulhigh(f, g, n2, n, T, p), n-n2, T, p);
    2067        1988 :     g = FpXX_sub(g, FpXX_shift(u, n2), p);
    2068        1988 :     if (gc_needed(av2,2))
    2069             :     {
    2070           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"FpXQXn_exp, e = %ld", n);
    2071           0 :       gerepileall(av2, 2, &f, &g);
    2072             :     }
    2073             :   }
    2074         714 :   return gerepileupto(av, f);
    2075             : }
    2076             : 
    2077             : GEN
    2078         178 : FpXQXn_expint(GEN h, long e, GEN T, GEN p)
    2079         178 : { return ZlXQXn_expint(h, e, T, p, 0); }
    2080             : 
    2081             : GEN
    2082           0 : FpXQXn_exp(GEN h, long e, GEN T, GEN p)
    2083             : {
    2084           0 :   if (signe(h)==0 || degpol(h)<1 || !gequal0(gel(h,2)))
    2085           0 :     pari_err_DOMAIN("FpXQXn_exp","valuation", "<", gen_1, h);
    2086           0 :   return FpXQXn_expint(FpXX_deriv(h, p), e, T, p);
    2087             : }
    2088             : 
    2089             : GEN
    2090         714 : FpXQXn_div(GEN g, GEN f, long e, GEN T, GEN p)
    2091             : {
    2092         714 :   pari_sp av = avma, av2;
    2093             :   ulong mask;
    2094             :   GEN W, a;
    2095         714 :   long v = varn(f), n = 1;
    2096             : 
    2097         714 :   if (!signe(f)) pari_err_INV("FpXXn_inv",f);
    2098         714 :   a = Fq_inv(gel(f,2), T, p);
    2099         714 :   if (e == 1 && !g) return scalarpol(a, v);
    2100         714 :   else if (e == 2 && !g)
    2101             :   {
    2102             :     GEN b;
    2103           0 :     if (degpol(f) <= 0) return scalarpol(a, v);
    2104           0 :     b = Fq_neg(gel(f,3),T,p);
    2105           0 :     if (signe(b)==0) return scalarpol(a, v);
    2106           0 :     if (!is_pm1(a)) b = Fq_mul(b, Fq_sqr(a, T, p), T, p);
    2107           0 :     W = deg1pol_shallow(b, a, v);
    2108           0 :     return gerepilecopy(av, W);
    2109             :   }
    2110         714 :   W = scalarpol_shallow(Fq_inv(gel(f,2), T, p),v);
    2111         714 :   mask = quadratic_prec_mask(e);
    2112         714 :   av2 = avma;
    2113        3416 :   for (;mask>1;)
    2114             :   {
    2115             :     GEN u, fr;
    2116        2702 :     long n2 = n;
    2117        2702 :     n<<=1; if (mask & 1) n--;
    2118        2702 :     mask >>= 1;
    2119        2702 :     fr = FpXXn_red(f, n);
    2120        2702 :     if (mask>1 || !g)
    2121             :     {
    2122        2702 :       u = FpXQXn_mul(W, FpXQXn_mulhigh(fr, W, n2, n, T, p), n-n2, T, p);
    2123        2702 :       W = FpXX_sub(W, FpXX_shift(u, n2), p);
    2124             :     }
    2125             :     else
    2126             :     {
    2127           0 :       GEN y = FpXQXn_mul(g, W, n, T, p), yt =  FpXXn_red(y, n-n2);
    2128           0 :       u = FpXQXn_mul(yt, FpXQXn_mulhigh(fr,  W, n2, n, T, p), n-n2, T, p);
    2129           0 :       W = FpXX_sub(y, FpXX_shift(u, n2), p);
    2130             :     }
    2131        2702 :     if (gc_needed(av2,2))
    2132             :     {
    2133           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"FpXQXn_inv, e = %ld", n);
    2134           0 :       W = gerepileupto(av2, W);
    2135             :     }
    2136             :   }
    2137         714 :   return gerepileupto(av, W);
    2138             : }
    2139             : 
    2140             : GEN
    2141         714 : FpXQXn_inv(GEN f, long e, GEN T, GEN p)
    2142         714 : { return FpXQXn_div(NULL, f, e, T, p); }

Generated by: LCOV version 1.14