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 - FpXQX_factor.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.18.1 lcov report (development 29950-285c5b69ed) Lines: 1614 1976 81.7 %
Date: 2025-02-05 09:09:51 Functions: 127 157 80.9 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2016  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             : #define DEBUGLEVEL DEBUGLEVEL_factorff
      19             : 
      20             : /*******************************************************************/
      21             : /**                                                               **/
      22             : /**           Isomorphisms between finite fields                  **/
      23             : /**                                                               **/
      24             : /*******************************************************************/
      25             : static void
      26          14 : err_Flxq(const char *s, GEN P, ulong l)
      27             : {
      28          14 :   if (!uisprime(l)) pari_err_PRIME(s, utoi(l));
      29          14 :   pari_err_IRREDPOL(s, Flx_to_ZX(get_Flx_mod(P)));
      30           0 : }
      31             : static void
      32           0 : err_FpXQ(const char *s, GEN P, GEN l)
      33             : {
      34           0 :   if (!BPSW_psp(l)) pari_err_PRIME(s, l);
      35           0 :   pari_err_IRREDPOL(s, get_FpX_mod(P));
      36           0 : }
      37             : 
      38             : /* compute the reciprocical isomorphism of S mod T,p, i.e. V such that
      39             :  * V(S)=X  mod T,p*/
      40             : static GEN
      41        2871 : Flxq_ffisom_inv_pre(GEN S, GEN T, ulong p, ulong pi)
      42             : {
      43        2871 :   pari_sp ltop = avma;
      44        2871 :   long n = get_Flx_degree(T);
      45        2871 :   GEN M = Flxq_matrix_pow_pre(S,n,n,T,p,pi);
      46        2871 :   GEN V = Flm_Flc_invimage(M, vecsmall_ei(n, 2), p);
      47        2871 :   if (!V) err_Flxq("Flxq_ffisom_inv", T, p);
      48        2871 :   return gerepileupto(ltop, Flv_to_Flx(V, get_Flx_var(T)));
      49             : }
      50             : GEN
      51           0 : Flxq_ffisom_inv(GEN S, GEN T, ulong p)
      52           0 : { return Flxq_ffisom_inv_pre(S, T, p, SMALL_ULONG(p)? 0: get_Fl_red(p)); }
      53             : 
      54             : GEN
      55         588 : FpXQ_ffisom_inv(GEN S,GEN T, GEN p)
      56             : {
      57         588 :   pari_sp ltop = avma;
      58         588 :   long n = get_FpX_degree(T);
      59         588 :   GEN M = FpXQ_matrix_pow(S,n,n,T,p);
      60         588 :   GEN V = FpM_FpC_invimage(M, col_ei(n, 2), p);
      61         588 :   if (!V) err_FpXQ("Flxq_ffisom_inv", T, p);
      62         588 :   return gerepilecopy(ltop, RgV_to_RgX(V, get_FpX_var(T)));
      63             : }
      64             : 
      65             : /* Let M the matrix of the Frobenius automorphism of Fp[X]/(T). Compute M^d
      66             :  * TODO: use left-right binary (tricky!) */
      67             : GEN
      68         399 : Flm_Frobenius_pow(GEN M, long d, GEN T, ulong p)
      69             : {
      70         399 :   pari_sp ltop=avma;
      71         399 :   long i,l = get_Flx_degree(T);
      72         399 :   GEN R, W = gel(M,2);
      73        1337 :   for (i = 2; i <= d; ++i) W = Flm_Flc_mul(M,W,p);
      74         399 :   R=Flxq_matrix_pow(Flv_to_Flx(W,get_Flx_var(T)),l,l,T,p);
      75         399 :   return gerepileupto(ltop,R);
      76             : }
      77             : 
      78             : GEN
      79          35 : FpM_Frobenius_pow(GEN M, long d, GEN T, GEN p)
      80             : {
      81          35 :   pari_sp ltop=avma;
      82          35 :   long i,l = get_FpX_degree(T);
      83          35 :   GEN R, W = gel(M,2);
      84         147 :   for (i = 2; i <= d; ++i) W = FpM_FpC_mul(M,W,p);
      85          35 :   R=FpXQ_matrix_pow(RgV_to_RgX(W, get_FpX_var(T)),l,l,T,p);
      86          35 :   return gerepilecopy(ltop,R);
      87             : }
      88             : 
      89             : /* Essentially we want to compute FqM_ker(MA-pol_x(v),U,l)
      90             :  * To avoid use of matrix in Fq we compute FpM_ker(U(MA),l) then recover the
      91             :  * eigenvalue by Galois action */
      92             : static GEN
      93       40331 : Flx_Flm_Flc_eval(GEN U, GEN MA, GEN a, ulong p)
      94             : {
      95       40331 :   long i, l = lg(U);
      96       40331 :   GEN b = Flv_Fl_mul(a, uel(U, l-1), p);
      97      163805 :   for (i=l-2; i>=2; i--)
      98      123474 :     b = Flv_add(Flm_Flc_mul(MA, b, p), Flv_Fl_mul(a, uel(U, i), p), p);
      99       40331 :   return b;
     100             : }
     101             : 
     102             : static GEN
     103       39235 : Flx_intersect_ker(GEN P, GEN MA, GEN U, ulong p)
     104             : {
     105       39235 :   pari_sp ltop = avma;
     106       39235 :   long i, vp = get_Flx_var(P), vu = get_Flx_var(U), r = get_Flx_degree(U);
     107             :   GEN V, A, R;
     108             :   ulong ib0;
     109             :   pari_timer T;
     110       39235 :   if (DEBUGLEVEL>=4) timer_start(&T);
     111       39235 :   V = Flx_div(Flx_Fl_add(monomial_Flx(1, get_Flx_degree(P), vu), p-1, p), U, p);
     112             :   do
     113             :   {
     114       40330 :     A = Flx_Flm_Flc_eval(V, MA, random_Flv(lg(MA)-1, p), p);
     115       40331 :   } while (zv_equal0(A));
     116       39236 :   if (DEBUGLEVEL>=4) timer_printf(&T,"matrix polcyclo");
     117             :   /*The formula is
     118             :    * a_{r-1} = -\phi(a_0)/b_0
     119             :    * a_{i-1} = \phi(a_i)+b_ia_{r-1}  i=r-1 to 1
     120             :    * Where a_0=A[1] and b_i=U[i+2] */
     121       39236 :   ib0 = Fl_inv(Fl_neg(U[2], p), p);
     122       39236 :   R = cgetg(r+1,t_MAT);
     123       39236 :   gel(R,1) = A;
     124       39236 :   gel(R,r) = Flm_Flc_mul(MA, Flv_Fl_mul(A,ib0, p), p);
     125       45087 :   for(i=r-1; i>1; i--)
     126             :   {
     127        5852 :     gel(R,i) = Flm_Flc_mul(MA,gel(R,i+1),p);
     128        5852 :     Flv_add_inplace(gel(R,i), Flv_Fl_mul(gel(R,r), U[i+2], p), p);
     129             :   }
     130       39235 :   return gerepileupto(ltop, Flm_to_FlxX(Flm_transpose(R),vp,vu));
     131             : }
     132             : 
     133             : static GEN
     134         182 : FpX_FpM_FpC_eval(GEN U, GEN MA, GEN a, GEN p)
     135             : {
     136         182 :   long i, l = lg(U);
     137         182 :   GEN b = FpC_Fp_mul(a, gel(U, l-1), p);
     138         966 :   for (i=l-2; i>=2; i--)
     139         784 :     b = FpC_add(FpM_FpC_mul(MA, b, p), FpC_Fp_mul(a, gel(U, i), p), p);
     140         182 :   return b;
     141             : }
     142             : 
     143             : static GEN
     144         182 : FpX_intersect_ker(GEN P, GEN MA, GEN U, GEN l)
     145             : {
     146         182 :   pari_sp ltop = avma;
     147         182 :   long i, vp = get_FpX_var(P), vu = get_FpX_var(U), r = get_FpX_degree(U);
     148             :   GEN V, A, R, ib0;
     149             :   pari_timer T;
     150         182 :   if (DEBUGLEVEL>=4) timer_start(&T);
     151         182 :   V = FpX_div(FpX_Fp_sub(pol_xn(get_FpX_degree(P), vu), gen_1, l), U, l);
     152             :   do
     153             :   {
     154         182 :     A = FpX_FpM_FpC_eval(V, MA, random_FpC(lg(MA)-1, l), l);
     155         182 :   } while (ZV_equal0(A));
     156         182 :   if (DEBUGLEVEL>=4) timer_printf(&T,"matrix polcyclo");
     157             :   /*The formula is
     158             :    * a_{r-1} = -\phi(a_0)/b_0
     159             :    * a_{i-1} = \phi(a_i)+b_ia_{r-1}  i=r-1 to 1
     160             :    * Where a_0=A[1] and b_i=U[i+2] */
     161         182 :   ib0 = Fp_inv(negi(gel(U,2)),l);
     162         182 :   R = cgetg(r+1,t_MAT);
     163         182 :   gel(R,1) = A;
     164         182 :   gel(R,r) = FpM_FpC_mul(MA, FpC_Fp_mul(A,ib0,l), l);
     165         518 :   for(i=r-1;i>1;i--)
     166         336 :     gel(R,i) = FpC_add(FpM_FpC_mul(MA,gel(R,i+1),l),
     167         336 :         FpC_Fp_mul(gel(R,r), gel(U,i+2), l),l);
     168         182 :   return gerepilecopy(ltop,RgM_to_RgXX(shallowtrans(R),vp,vu));
     169             : }
     170             : 
     171             : /* n must divide both the degree of P and Q.  Compute SP and SQ such
     172             :  * that the subfield of FF_l[X]/(P) generated by SP and the subfield of
     173             :  * FF_l[X]/(Q) generated by SQ are isomorphic of degree n.  P and Q do
     174             :  * not need to be of the same variable; if MA, resp. MB, is not NULL, must be
     175             :  * the matrix of the Frobenius map in FF_l[X]/(P), resp. FF_l[X]/(Q).
     176             :  * Implementation choice:  we assume the prime p is large so we handle
     177             :  * Frobenius as matrices. */
     178             : static void
     179       43307 : Flx_ffintersect_pre(GEN P, GEN Q, long n, ulong l, ulong li, GEN *SP, GEN *SQ, GEN MA, GEN MB)
     180             : {
     181       43307 :   pari_sp ltop = avma;
     182       43307 :   long vp = get_Flx_var(P), vq =  get_Flx_var(Q);
     183       43307 :   long np = get_Flx_degree(P), nq = get_Flx_degree(Q), e;
     184             :   ulong pg;
     185             :   GEN A, B, Ap, Bp;
     186       43308 :   if (np<=0) pari_err_IRREDPOL("FpX_ffintersect", P);
     187       43308 :   if (nq<=0) pari_err_IRREDPOL("FpX_ffintersect", Q);
     188       43308 :   if (n<=0 || np%n || nq%n)
     189           0 :     pari_err_TYPE("FpX_ffintersect [bad degrees]",stoi(n));
     190       43308 :   li = SMALL_ULONG(l)? 0: get_Fl_red(l);
     191       43308 :   e = u_lvalrem(n, l, &pg);
     192       43308 :   if(!MA) MA = Flx_matFrobenius_pre(P,l,li);
     193       43308 :   if(!MB) MB = Flx_matFrobenius_pre(Q,l,li);
     194       43308 :   A = Ap = pol0_Flx(vp);
     195       43308 :   B = Bp = pol0_Flx(vq);
     196       43308 :   if (pg > 1)
     197             :   {
     198             :     pari_timer T;
     199       40218 :     GEN ipg = utoipos(pg);
     200       40218 :     if (l%pg == 1)
     201             :     { /* more efficient special case */
     202             :       ulong L, z, An, Bn;
     203       20600 :       z = Fl_neg(rootsof1_Fl(pg, l), l);
     204       20600 :       if (DEBUGLEVEL>=4) timer_start(&T);
     205       20600 :       A = Flm_ker(Flm_Fl_add(MA, z, l),l);
     206       20600 :       if (lg(A)!=2) err_Flxq("FpX_ffintersect",P,l);
     207       20600 :       A = Flv_to_Flx(gel(A,1),vp);
     208             : 
     209       20600 :       B = Flm_ker(Flm_Fl_add(MB, z, l),l);
     210       20600 :       if (lg(B)!=2) err_Flxq("FpX_ffintersect",Q,l);
     211       20593 :       B = Flv_to_Flx(gel(B,1),vq);
     212             : 
     213       20593 :       if (DEBUGLEVEL>=4) timer_printf(&T, "FpM_ker");
     214       20593 :       An = Flxq_powu_pre(A,pg,P,l,li)[2];
     215       20593 :       Bn = Flxq_powu_pre(B,pg,Q,l,li)[2];
     216       20593 :       if (!Bn) pari_err_IRREDPOL("FpX_ffintersect", mkvec2(P,Q));
     217       20593 :       z = Fl_div(An,Bn,l);
     218       20593 :       L = Fl_sqrtn(z, pg, l, NULL);
     219       20593 :       if (L==ULONG_MAX) pari_err_IRREDPOL("FpX_ffintersect", mkvec2(P,Q));
     220       20593 :       if (DEBUGLEVEL>=4) timer_printf(&T, "Fp_sqrtn");
     221       20593 :       B = Flx_Fl_mul(B,L,l);
     222             :     }
     223             :     else
     224             :     {
     225             :       GEN L, An, Bn, z, U;
     226       19618 :       U = gmael(Flx_factor(ZX_to_Flx(polcyclo(pg, fetch_var()),l),l),1,1);
     227       19617 :       A = Flx_intersect_ker(P, MA, U, l);
     228       19618 :       B = Flx_intersect_ker(Q, MB, U, l);
     229       19617 :       if (DEBUGLEVEL>=4) timer_start(&T);
     230       19617 :       An = gel(FlxYqq_pow(A,ipg,P,U,l),2);
     231       19618 :       Bn = gel(FlxYqq_pow(B,ipg,Q,U,l),2);
     232       19618 :       if (DEBUGLEVEL>=4) timer_printf(&T,"pows [P,Q]");
     233       19618 :       z = Flxq_div_pre(An,Bn,U,l,li);
     234       19618 :       L = Flxq_sqrtn(z,ipg,U,l,NULL);
     235       19618 :       if (!L) pari_err_IRREDPOL("FpX_ffintersect", mkvec2(P,Q));
     236       19618 :       if (DEBUGLEVEL>=4) timer_printf(&T,"FpXQ_sqrtn");
     237       19618 :       B = FlxqX_Flxq_mul_pre(B,L,U,l,li);
     238       19618 :       A = FlxY_evalx_pre(A,0,l,li);
     239       19618 :       B = FlxY_evalx_pre(B,0,l,li);
     240       19618 :       (void)delete_var();
     241             :     }
     242             :   }
     243       43301 :   if (e)
     244             :   {
     245             :     GEN VP, VQ, Ay, By;
     246        3153 :     ulong lmun = l-1;
     247             :     long j;
     248        3153 :     MA = Flm_Fl_add(MA,lmun,l);
     249        3153 :     MB = Flm_Fl_add(MB,lmun,l);
     250        3153 :     Ay = pol1_Flx(vp);
     251        3153 :     By = pol1_Flx(vq);
     252        3153 :     VP = vecsmall_ei(np, 1);
     253        3153 :     VQ = np == nq? VP: vecsmall_ei(nq, 1); /* save memory */
     254        6715 :     for(j=0;j<e;j++)
     255             :     {
     256        3569 :       if (j)
     257             :       {
     258         416 :         Ay = Flxq_mul_pre(Ay,Flxq_powu_pre(Ap,lmun,P,l,li),P,l,li);
     259         416 :         VP = Flx_to_Flv(Ay,np);
     260             :       }
     261        3569 :       Ap = Flm_Flc_invimage(MA,VP,l);
     262        3569 :       if (!Ap) err_Flxq("FpX_ffintersect",P,l);
     263        3569 :       Ap = Flv_to_Flx(Ap,vp);
     264             : 
     265        3569 :       if (j)
     266             :       {
     267         416 :         By = Flxq_mul_pre(By,Flxq_powu_pre(Bp,lmun,Q,l,li),Q,l,li);
     268         416 :         VQ = Flx_to_Flv(By,nq);
     269             :       }
     270        3569 :       Bp = Flm_Flc_invimage(MB,VQ,l);
     271        3569 :       if (!Bp) err_Flxq("FpX_ffintersect",Q,l);
     272        3562 :       Bp = Flv_to_Flx(Bp,vq);
     273             :     }
     274             :   }
     275       43294 :   *SP = Flx_add(A,Ap,l);
     276       43294 :   *SQ = Flx_add(B,Bp,l);
     277       43293 :   gerepileall(ltop,2,SP,SQ);
     278       43294 : }
     279             : void
     280           0 : Flx_ffintersect(GEN P, GEN Q, long n, ulong p, GEN *SP, GEN *SQ, GEN MA, GEN MB)
     281             : {
     282           0 :   ulong pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
     283           0 :   Flx_ffintersect_pre(P, Q, n, p, pi, SP, SQ, MA, MB);
     284           0 : }
     285             : 
     286             : /* Let l be a prime number, P, Q in Z[X]; both are irreducible modulo l and
     287             :  * degree(P) divides degree(Q).  Output a monomorphism between F_l[X]/(P) and
     288             :  * F_l[X]/(Q) as a polynomial R such that Q | P(R) mod l.  If P and Q have the
     289             :  * same degree, it is of course an isomorphism.  */
     290             : GEN
     291        2871 : Flx_ffisom(GEN P,GEN Q,ulong l)
     292             : {
     293        2871 :   pari_sp av = avma;
     294             :   GEN SP, SQ, R;
     295        2871 :   ulong li = SMALL_ULONG(l)? 0: get_Fl_red(l);
     296        2871 :   Flx_ffintersect_pre(P,Q,get_Flx_degree(P),l,li,&SP,&SQ,NULL,NULL);
     297        2871 :   R = Flxq_ffisom_inv_pre(SP,P,l,li);
     298        2871 :   return gerepileupto(av, Flx_Flxq_eval_pre(R,SQ,Q,l,li));
     299             : }
     300             : 
     301             : void
     302         322 : FpX_ffintersect(GEN P, GEN Q, long n, GEN l, GEN *SP, GEN *SQ, GEN MA, GEN MB)
     303             : {
     304         322 :   pari_sp ltop = avma;
     305             :   long vp, vq, np, nq, e;
     306             :   ulong pg;
     307             :   GEN A, B, Ap, Bp;
     308         322 :   if (lgefint(l)==3)
     309             :   {
     310           0 :     ulong pp = l[2];
     311           0 :     GEN Pp = ZX_to_Flx(P,pp), Qp = ZX_to_Flx(Q,pp);
     312           0 :     GEN MAp = MA ? ZM_to_Flm(MA, pp): NULL;
     313           0 :     GEN MBp = MB ? ZM_to_Flm(MB, pp): NULL;
     314           0 :     Flx_ffintersect(Pp, Qp, n, pp, SP, SQ, MAp, MBp);
     315           0 :     *SP = Flx_to_ZX(*SP); *SQ = Flx_to_ZX(*SQ);
     316           0 :     gerepileall(ltop,2,SP,SQ);
     317           0 :     return;
     318             :   }
     319         322 :   vp = get_FpX_var(P); np = get_FpX_degree(P);
     320         322 :   vq = get_FpX_var(Q); nq = get_FpX_degree(Q);
     321         322 :   if (np<=0) pari_err_IRREDPOL("FpX_ffintersect", P);
     322         322 :   if (nq<=0) pari_err_IRREDPOL("FpX_ffintersect", Q);
     323         322 :   if (n<=0 || np%n || nq%n)
     324           0 :     pari_err_TYPE("FpX_ffintersect [bad degrees]",stoi(n));
     325         322 :   e = u_pvalrem(n, l, &pg);
     326         322 :   if(!MA) MA = FpX_matFrobenius(P, l);
     327         322 :   if(!MB) MB = FpX_matFrobenius(Q, l);
     328         322 :   A = Ap = pol_0(vp);
     329         322 :   B = Bp = pol_0(vq);
     330         322 :   if (pg > 1)
     331             :   {
     332         322 :     GEN ipg = utoipos(pg);
     333             :     pari_timer T;
     334         322 :     if (umodiu(l,pg) == 1)
     335             :     /* No need to use relative extension, so don't. (Well, now we don't
     336             :      * in the other case either, but this special case is more efficient) */
     337             :     {
     338             :       GEN L, An, Bn, z;
     339         231 :       z = negi( rootsof1u_Fp(pg, l) );
     340         231 :       if (DEBUGLEVEL>=4) timer_start(&T);
     341         231 :       A = FpM_ker(RgM_Rg_add_shallow(MA, z),l);
     342         231 :       if (lg(A)!=2) err_FpXQ("FpX_ffintersect",P,l);
     343         231 :       A = RgV_to_RgX(gel(A,1),vp);
     344             : 
     345         231 :       B = FpM_ker(RgM_Rg_add_shallow(MB, z),l);
     346         231 :       if (lg(B)!=2) err_FpXQ("FpX_ffintersect",Q,l);
     347         231 :       B = RgV_to_RgX(gel(B,1),vq);
     348             : 
     349         231 :       if (DEBUGLEVEL>=4) timer_printf(&T, "FpM_ker");
     350         231 :       An = gel(FpXQ_pow(A,ipg,P,l),2);
     351         231 :       Bn = gel(FpXQ_pow(B,ipg,Q,l),2);
     352         231 :       if (!signe(Bn)) pari_err_IRREDPOL("FpX_ffintersect", mkvec2(P,Q));
     353         231 :       z = Fp_div(An,Bn,l);
     354         231 :       L = Fp_sqrtn(z,ipg,l,NULL);
     355         231 :       if (!L) pari_err_IRREDPOL("FpX_ffintersect", mkvec2(P,Q));
     356         231 :       if (DEBUGLEVEL>=4) timer_printf(&T, "Fp_sqrtn");
     357         231 :       B = FpX_Fp_mul(B,L,l);
     358             :     }
     359             :     else
     360             :     {
     361             :       GEN L, An, Bn, z, U;
     362          91 :       U = gmael(FpX_factor(polcyclo(pg,fetch_var()),l),1,1);
     363          91 :       A = FpX_intersect_ker(P, MA, U, l);
     364          91 :       B = FpX_intersect_ker(Q, MB, U, l);
     365          91 :       if (DEBUGLEVEL>=4) timer_start(&T);
     366          91 :       An = gel(FpXYQQ_pow(A,ipg,P,U,l),2);
     367          91 :       Bn = gel(FpXYQQ_pow(B,ipg,Q,U,l),2);
     368          91 :       if (DEBUGLEVEL>=4) timer_printf(&T,"pows [P,Q]");
     369          91 :       if (!signe(Bn)) pari_err_IRREDPOL("FpX_ffintersect", mkvec2(P,Q));
     370          91 :       z = Fq_div(An,Bn,U,l);
     371          91 :       L = Fq_sqrtn(z,ipg,U,l,NULL);
     372          91 :       if (!L) pari_err_IRREDPOL("FpX_ffintersect", mkvec2(P,Q));
     373          91 :       if (DEBUGLEVEL>=4) timer_printf(&T,"FpXQ_sqrtn");
     374          91 :       B = FqX_Fq_mul(B,L,U,l);
     375          91 :       A = FpXY_evalx(A,gen_0,l);
     376          91 :       B = FpXY_evalx(B,gen_0,l);
     377          91 :       (void)delete_var();
     378             :     }
     379             :   }
     380         322 :   if (e)
     381             :   {
     382           0 :     GEN VP, VQ, Ay, By, lmun = subiu(l,1);
     383             :     long j;
     384           0 :     MA = RgM_Rg_add_shallow(MA,gen_m1);
     385           0 :     MB = RgM_Rg_add_shallow(MB,gen_m1);
     386           0 :     Ay = pol_1(vp);
     387           0 :     By = pol_1(vq);
     388           0 :     VP = col_ei(np, 1);
     389           0 :     VQ = np == nq? VP: col_ei(nq, 1); /* save memory */
     390           0 :     for(j=0;j<e;j++)
     391             :     {
     392           0 :       if (j)
     393             :       {
     394           0 :         Ay = FpXQ_mul(Ay,FpXQ_pow(Ap,lmun,P,l),P,l);
     395           0 :         VP = RgX_to_RgC(Ay,np);
     396             :       }
     397           0 :       Ap = FpM_FpC_invimage(MA,VP,l);
     398           0 :       if (!Ap) err_FpXQ("FpX_ffintersect",P,l);
     399           0 :       Ap = RgV_to_RgX(Ap,vp);
     400             : 
     401           0 :       if (j)
     402             :       {
     403           0 :         By = FpXQ_mul(By,FpXQ_pow(Bp,lmun,Q,l),Q,l);
     404           0 :         VQ = RgX_to_RgC(By,nq);
     405             :       }
     406           0 :       Bp = FpM_FpC_invimage(MB,VQ,l);
     407           0 :       if (!Bp) err_FpXQ("FpX_ffintersect",Q,l);
     408           0 :       Bp = RgV_to_RgX(Bp,vq);
     409             :     }
     410             :   }
     411         322 :   *SP = FpX_add(A,Ap,l);
     412         322 :   *SQ = FpX_add(B,Bp,l);
     413         322 :   gerepileall(ltop,2,SP,SQ);
     414             : }
     415             : /* Let l be a prime number, P, Q in Z[X]; both are irreducible modulo l and
     416             :  * degree(P) divides degree(Q).  Output a monomorphism between F_l[X]/(P) and
     417             :  * F_l[X]/(Q) as a polynomial R such that Q | P(R) mod l.  If P and Q have the
     418             :  * same degree, it is of course an isomorphism.  */
     419             : GEN
     420        2829 : FpX_ffisom(GEN P, GEN Q, GEN p)
     421             : {
     422        2829 :   pari_sp av = avma;
     423             :   GEN SP, SQ, R;
     424        2829 :   if (lgefint(p)==3)
     425             :   {
     426        2829 :     ulong pp = p[2];
     427        2829 :     GEN R = Flx_ffisom(ZX_to_Flx(P,pp), ZX_to_Flx(Q,pp), pp);
     428        2829 :     return gerepileupto(av, Flx_to_ZX(R));
     429             :   }
     430           0 :   FpX_ffintersect(P,Q,get_FpX_degree(P),p,&SP,&SQ,NULL,NULL);
     431           0 :   R = FpXQ_ffisom_inv(SP,P,p);
     432           0 :   return gerepileupto(av, FpX_FpXQ_eval(R,SQ,Q,p));
     433             : }
     434             : 
     435             : /* Let l be a prime number, P a ZX irreducible modulo l, MP the matrix of the
     436             :  * Frobenius automorphism of F_l[X]/(P).
     437             :  * Factor P over the subfield of F_l[X]/(P) of index d. */
     438             : static GEN
     439         322 : FpX_factorgalois(GEN P, GEN l, long d, long w, GEN MP)
     440             : {
     441         322 :   pari_sp ltop = avma;
     442             :   GEN R, V, Tl, z, M;
     443         322 :   long v = get_FpX_var(P), n = get_FpX_degree(P);
     444         322 :   long k, m = n/d;
     445             : 
     446             :   /* x - y */
     447         322 :   if (m == 1) return deg1pol_shallow(gen_1, deg1pol_shallow(subis(l,1), gen_0, w), v);
     448          35 :   M = FpM_Frobenius_pow(MP,d,P,l);
     449             : 
     450          35 :   Tl = leafcopy(P); setvarn(Tl,w);
     451          35 :   V = cgetg(m+1,t_VEC);
     452          35 :   gel(V,1) = pol_x(w);
     453          35 :   z = gel(M,2);
     454          35 :   gel(V,2) = RgV_to_RgX(z,w);
     455          77 :   for(k=3;k<=m;k++)
     456             :   {
     457          42 :     z = FpM_FpC_mul(M,z,l);
     458          42 :     gel(V,k) = RgV_to_RgX(z,w);
     459             :   }
     460          35 :   R = FqV_roots_to_pol(V,Tl,l,v);
     461          35 :   return gerepileupto(ltop,R);
     462             : }
     463             : /* same: P is an Flx, MP an Flm */
     464             : static GEN
     465       40423 : Flx_factorgalois(GEN P, ulong l, long d, long w, GEN MP)
     466             : {
     467       40423 :   pari_sp ltop = avma;
     468             :   GEN R, V, Tl, z, M;
     469       40423 :   long k, n = get_Flx_degree(P), m = n/d;
     470       40423 :   long v = get_Flx_var(P);
     471             : 
     472       40423 :   if (m == 1) {
     473       40024 :     R = polx_Flx(v);
     474       40024 :     gel(R,2) = z = polx_Flx(w); z[3] = l - 1; /* - y */
     475       40024 :     gel(R,3) = pol1_Flx(w);
     476       40024 :     return R; /* x - y */
     477             :   }
     478         399 :   M = Flm_Frobenius_pow(MP,d,P,l);
     479             : 
     480         399 :   Tl = leafcopy(P); Tl[1] = w;
     481         399 :   V = cgetg(m+1,t_VEC);
     482         399 :   gel(V,1) = polx_Flx(w);
     483         399 :   z = gel(M,2);
     484         399 :   gel(V,2) = Flv_to_Flx(z,w);
     485         700 :   for(k=3;k<=m;k++)
     486             :   {
     487         301 :     z = Flm_Flc_mul(M,z,l);
     488         301 :     gel(V,k) = Flv_to_Flx(z,w);
     489             :   }
     490         399 :   R = FlxqV_roots_to_pol(V,Tl,l,v);
     491         399 :   return gerepileupto(ltop,R);
     492             : }
     493             : 
     494             : GEN
     495      111517 : Flx_factorff_irred(GEN P, GEN Q, ulong p)
     496             : {
     497      111517 :   pari_sp ltop = avma, av;
     498             :   GEN SP, SQ, MP, MQ, M, FP, FQ, E, V, IR, res;
     499      111517 :   long np = get_Flx_degree(P), nq = get_Flx_degree(Q), d = ugcd(np,nq);
     500      111517 :   long i, vp = get_Flx_var(P), vq = get_Flx_var(Q);
     501             :   ulong pi, PI;
     502      111516 :   if (d==1) retmkcol(Flx_to_FlxX(P, vq));
     503       40436 :   PI = get_Fl_red(p);
     504       40436 :   pi = SMALL_ULONG(p)? 0: PI; /* PI for Fp, pi for Fp[x] */
     505       40436 :   FQ = Flx_matFrobenius_pre(Q,p,pi);
     506       40436 :   av = avma;
     507       40436 :   FP = Flx_matFrobenius_pre(P,p,pi);
     508       40436 :   Flx_ffintersect_pre(P,Q,d,p,pi,&SP,&SQ, FP, FQ);
     509       40423 :   E = Flx_factorgalois(P,p,d,vq, FP);
     510       40423 :   E = FlxX_to_Flm(E,np);
     511       40423 :   MP= Flxq_matrix_pow_pre(SP,np,d,P,p,pi);
     512       40423 :   IR= gel(Flm_indexrank(MP,p),1);
     513       40423 :   E = rowpermute(E, IR);
     514       40423 :   M = rowpermute(MP,IR);
     515       40423 :   M = Flm_inv(M,p);
     516       40423 :   MQ= Flxq_matrix_pow_pre(SQ,nq,d,Q,p,pi);
     517       40423 :   M = Flm_mul_pre(MQ,M,p,PI);
     518       40423 :   M = Flm_mul_pre(M,E,p,PI);
     519       40423 :   M = gerepileupto(av,M);
     520       40423 :   V = cgetg(d+1,t_VEC);
     521       40423 :   gel(V,1) = M;
     522      176116 :   for(i=2;i<=d;i++) gel(V,i) = Flm_mul_pre(FQ,gel(V,i-1),p,PI);
     523       40423 :   res = cgetg(d+1,t_COL);
     524      216535 :   for(i=1;i<=d;i++) gel(res,i) = Flm_to_FlxX(gel(V,i),vp,vq);
     525       40421 :   return gerepileupto(ltop,res);
     526             : }
     527             : 
     528             : /* P,Q irreducible over F_p. Factor P over FF_p[X] / Q  [factors are ordered as
     529             :  * a Frobenius cycle] */
     530             : GEN
     531       32879 : FpX_factorff_irred(GEN P, GEN Q, GEN p)
     532             : {
     533       32879 :   pari_sp ltop = avma, av;
     534             :   GEN res;
     535       32879 :   long np = get_FpX_degree(P), nq = get_FpX_degree(Q), d = ugcd(np,nq);
     536       32879 :   if (d==1) return mkcolcopy(P);
     537             : 
     538       32816 :   if (lgefint(p)==3)
     539             :   {
     540       32494 :     ulong pp = p[2];
     541       32494 :     GEN F = Flx_factorff_irred(ZX_to_Flx(P,pp), ZX_to_Flx(Q,pp), pp);
     542       32493 :     long i, lF = lg(F);
     543       32493 :     res = cgetg(lF, t_COL);
     544      182916 :     for(i=1; i<lF; i++)
     545      150424 :       gel(res,i) = FlxX_to_ZXX(gel(F,i));
     546             :   }
     547             :   else
     548             :   {
     549             :     GEN SP, SQ, MP, MQ, M, FP, FQ, E, V, IR;
     550         322 :     long i, vp = get_FpX_var(P), vq = get_FpX_var(Q);
     551         322 :     FQ = FpX_matFrobenius(Q,p);
     552         322 :     av = avma;
     553         322 :     FP = FpX_matFrobenius(P,p);
     554         322 :     FpX_ffintersect(P,Q,d,p,&SP,&SQ,FP,FQ);
     555             : 
     556         322 :     E = FpX_factorgalois(P,p,d,vq,FP);
     557         322 :     E = RgXX_to_RgM(E,np);
     558         322 :     MP= FpXQ_matrix_pow(SP,np,d,P,p);
     559         322 :     IR= gel(FpM_indexrank(MP,p),1);
     560         322 :     E = rowpermute(E, IR);
     561         322 :     M = rowpermute(MP,IR);
     562         322 :     M = FpM_inv(M,p);
     563         322 :     MQ= FpXQ_matrix_pow(SQ,nq,d,Q,p);
     564         322 :     M = FpM_mul(MQ,M,p);
     565         322 :     M = FpM_mul(M,E,p);
     566         322 :     M = gerepileupto(av,M);
     567         322 :     V = cgetg(d+1,t_VEC);
     568         322 :     gel(V,1) = M;
     569        1050 :     for(i=2;i<=d;i++)
     570         728 :       gel(V,i) = FpM_mul(FQ,gel(V,i-1),p);
     571         322 :     res = cgetg(d+1,t_COL);
     572        1372 :     for(i=1;i<=d;i++)
     573        1050 :       gel(res,i) = RgM_to_RgXX(gel(V,i),vp,vq);
     574             :   }
     575       32814 :   return gerepilecopy(ltop,res);
     576             : }
     577             : 
     578             : /* not memory-clean, as Flx_factorff_i, returning only linear factors */
     579             : static GEN
     580       29926 : Flx_rootsff_i(GEN P, GEN T, ulong p)
     581             : {
     582       29926 :   GEN V, F = gel(Flx_factor(P,p), 1);
     583       29926 :   long i, lfact = 1, nmax = lgpol(P), n = lg(F), dT = get_Flx_degree(T);
     584             : 
     585       29926 :   V = cgetg(nmax,t_COL);
     586       63701 :   for(i=1;i<n;i++)
     587             :   {
     588       33775 :     GEN R, Fi = gel(F,i);
     589       33775 :     long di = degpol(Fi), j, r;
     590       33775 :     if (dT % di) continue;
     591       32256 :     R = Flx_factorff_irred(gel(F,i),T,p);
     592       32256 :     r = lg(R);
     593       79497 :     for (j=1; j<r; j++,lfact++)
     594       47241 :       gel(V,lfact) = Flx_neg(gmael(R,j, 2), p);
     595             :   }
     596       29926 :   setlg(V,lfact);
     597       29926 :   gen_sort_inplace(V, (void*) &cmp_Flx, &cmp_nodata, NULL);
     598       29926 :   return V;
     599             : }
     600             : GEN
     601           0 : Flx_rootsff(GEN P, GEN T, ulong p)
     602             : {
     603           0 :   pari_sp av = avma;
     604           0 :   return gerepilecopy(av, Flx_rootsff_i(P, T, p));
     605             : }
     606             : 
     607             : /* dummy implementation */
     608             : static GEN
     609       16520 : F2x_rootsff_i(GEN P, GEN T)
     610             : {
     611       16520 :   return FlxC_to_F2xC(Flx_rootsff_i(F2x_to_Flx(P), F2x_to_Flx(T), 2UL));
     612             : }
     613             : 
     614             : /* not memory-clean, as FpX_factorff_i, returning only linear factors */
     615             : static GEN
     616         308 : FpX_rootsff_i(GEN P, GEN T, GEN p)
     617             : {
     618             :   GEN V, F;
     619             :   long i, lfact, nmax, n, dT;
     620         308 :   if (lgefint(p)==3)
     621             :   {
     622           0 :     ulong pp = p[2];
     623           0 :     GEN V = Flx_rootsff_i(ZX_to_Flx(P,pp), ZXT_to_FlxT(T,pp), pp);
     624           0 :     return FlxC_to_ZXC(V);
     625             :   }
     626         308 :   F = gel(FpX_factor(P,p), 1);
     627         308 :   lfact = 1; nmax = lgpol(P); n = lg(F); dT = get_FpX_degree(T);
     628             : 
     629         308 :   V = cgetg(nmax,t_COL);
     630         630 :   for(i=1;i<n;i++)
     631             :   {
     632         322 :     GEN R, Fi = gel(F,i);
     633         322 :     long di = degpol(Fi), j, r;
     634         322 :     if (dT % di) continue;
     635         322 :     R = FpX_factorff_irred(gel(F,i),T,p);
     636         322 :     r = lg(R);
     637        1260 :     for (j=1; j<r; j++,lfact++)
     638         938 :       gel(V,lfact) = Fq_to_FpXQ(Fq_neg(gmael(R,j, 2), T, p), T, p);
     639             :   }
     640         308 :   setlg(V,lfact);
     641         308 :   gen_sort_inplace(V, (void*) &cmp_RgX, &cmp_nodata, NULL);
     642         308 :   return V;
     643             : }
     644             : GEN
     645           0 : FpX_rootsff(GEN P, GEN T, GEN p)
     646             : {
     647           0 :   pari_sp av = avma;
     648           0 :   return gerepilecopy(av, FpX_rootsff_i(P, T, p));
     649             : }
     650             : 
     651             : static GEN
     652       11998 : Flx_factorff_i(GEN P, GEN T, ulong p)
     653             : {
     654       11998 :   GEN V, E, F = Flx_factor(P, p);
     655       11998 :   long i, lfact = 1, nmax = lgpol(P), n = lgcols(F);
     656             : 
     657       11998 :   V = cgetg(nmax,t_VEC);
     658       11998 :   E = cgetg(nmax,t_VECSMALL);
     659       58751 :   for(i=1;i<n;i++)
     660             :   {
     661       46767 :     GEN R = Flx_factorff_irred(gmael(F,1,i),T,p), e = gmael(F,2,i);
     662       46753 :     long j, r = lg(R);
     663       96278 :     for (j=1; j<r; j++,lfact++)
     664             :     {
     665       49525 :       gel(V,lfact) = gel(R,j);
     666       49525 :       gel(E,lfact) = e;
     667             :     }
     668             :   }
     669       11984 :   setlg(V,lfact);
     670       11984 :   setlg(E,lfact); return sort_factor_pol(mkvec2(V,E), cmp_Flx);
     671             : }
     672             : 
     673             : static long
     674        7084 : simpleff_to_nbfact(GEN F, long dT)
     675             : {
     676        7084 :   long i, l = lg(F), k = 0;
     677       90146 :   for (i = 1; i < l; i++) k += ugcd(uel(F,i), dT);
     678        7084 :   return k;
     679             : }
     680             : 
     681             : static long
     682        7084 : Flx_nbfactff(GEN P, GEN T, ulong p)
     683             : {
     684        7084 :   pari_sp av = avma;
     685        7084 :   GEN F = gel(Flx_degfact(P, p), 1);
     686        7084 :   long s = simpleff_to_nbfact(F, get_Flx_degree(T));
     687        7084 :   return gc_long(av,s);
     688             : }
     689             : 
     690             : /* dummy implementation */
     691             : static GEN
     692         420 : F2x_factorff_i(GEN P, GEN T)
     693             : {
     694         420 :   GEN M = Flx_factorff_i(F2x_to_Flx(P), F2x_to_Flx(T), 2);
     695         413 :   return mkvec2(FlxXC_to_F2xXC(gel(M,1)), gel(M,2));
     696             : }
     697             : 
     698             : /* not memory-clean */
     699             : static GEN
     700          63 : FpX_factorff_i(GEN P, GEN T, GEN p)
     701             : {
     702          63 :   GEN V, E, F = FpX_factor(P,p);
     703          63 :   long i, lfact = 1, nmax = lgpol(P), n = lgcols(F);
     704             : 
     705          63 :   V = cgetg(nmax,t_VEC);
     706          63 :   E = cgetg(nmax,t_VECSMALL);
     707         126 :   for(i=1;i<n;i++)
     708             :   {
     709          63 :     GEN R = FpX_factorff_irred(gmael(F,1,i),T,p), e = gmael(F,2,i);
     710          63 :     long j, r = lg(R);
     711         238 :     for (j=1; j<r; j++,lfact++)
     712             :     {
     713         175 :       gel(V,lfact) = gel(R,j);
     714         175 :       gel(E,lfact) = e;
     715             :     }
     716             :   }
     717          63 :   setlg(V,lfact);
     718          63 :   setlg(E,lfact); return sort_factor_pol(mkvec2(V,E), cmp_RgX);
     719             : }
     720             : 
     721             : static long
     722           0 : FpX_nbfactff(GEN P, GEN T, GEN p)
     723             : {
     724           0 :   pari_sp av = avma;
     725           0 :   GEN F = gel(FpX_degfact(P, p), 1);
     726           0 :   long s = simpleff_to_nbfact(F, get_FpX_degree(T));
     727           0 :   return gc_long(av,s);
     728             : }
     729             : 
     730             : GEN
     731           0 : FpX_factorff(GEN P, GEN T, GEN p)
     732             : {
     733           0 :   pari_sp av = avma;
     734           0 :   return gerepilecopy(av, FpX_factorff_i(P, T, p));
     735             : }
     736             : 
     737             : /***********************************************************************/
     738             : /**                                                                   **/
     739             : /**               Factorisation over finite fields                    **/
     740             : /**                                                                   **/
     741             : /***********************************************************************/
     742             : 
     743             : static GEN
     744       12259 : FlxqXQ_halfFrobenius_i(GEN a, GEN xp, GEN Xp, GEN S, GEN T, ulong p, ulong pi)
     745             : {
     746       12259 :   GEN ap2 = FlxqXQ_powu_pre(a, p>>1, S, T, p, pi);
     747       12259 :   GEN V = FlxqXQ_autsum_pre(mkvec3(xp, Xp, ap2), get_Flx_degree(T), S, T, p,pi);
     748       12259 :   return gel(V,3);
     749             : }
     750             : 
     751             : GEN
     752         470 : FlxqXQ_halfFrobenius(GEN a, GEN S, GEN T, ulong p)
     753             : {
     754         470 :   ulong pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
     755         470 :   long vT = get_Flx_var(T);
     756             :   GEN xp, Xp;
     757         470 :   T = Flx_get_red_pre(T, p, pi);
     758         470 :   S = FlxqX_get_red_pre(S, T, p, pi);
     759         470 :   xp = Flx_Frobenius_pre(T, p, pi);
     760         470 :   Xp = FlxqXQ_powu_pre(polx_FlxX(get_FlxqX_var(S), vT), p, S, T, p, pi);
     761         470 :   return FlxqXQ_halfFrobenius_i(a, xp, Xp, S, T, p, pi);
     762             : }
     763             : 
     764             : static GEN
     765        1194 : FpXQXQ_halfFrobenius_i(GEN a, GEN xp, GEN Xp, GEN S, GEN T, GEN p)
     766             : {
     767        1194 :   GEN ap2 = FpXQXQ_pow(a, shifti(p,-1), S, T, p);
     768        1194 :   GEN V = FpXQXQ_autsum(mkvec3(xp, Xp, ap2), get_FpX_degree(T), S, T, p);
     769        1194 :   return gel(V, 3);
     770             : }
     771             : 
     772             : GEN
     773         238 : FpXQXQ_halfFrobenius(GEN a, GEN S, GEN T, GEN p)
     774             : {
     775         238 :   pari_sp av = avma;
     776             :   GEN z;
     777         238 :   if (lgefint(p)==3)
     778             :   {
     779          99 :     ulong pp = p[2];
     780          99 :     long v = get_FpX_var(T);
     781          99 :     GEN Tp = ZXT_to_FlxT(T,pp), Sp = ZXXT_to_FlxXT(S, pp, v);
     782          99 :     z = FlxX_to_ZXX(FlxqXQ_halfFrobenius(ZXX_to_FlxX(a,pp,v),Sp,Tp,pp));
     783             :   }
     784             :   else
     785             :   {
     786             :     GEN xp, Xp;
     787         139 :     T = FpX_get_red(T, p);
     788         139 :     S = FpXQX_get_red(S, T, p);
     789         139 :     xp = FpX_Frobenius(T, p);
     790         139 :     Xp = FpXQXQ_pow(pol_x(get_FpXQX_var(S)), p, S, T, p);
     791         139 :     z = FpXQXQ_halfFrobenius_i(a, xp, Xp, S, T, p);
     792             :   }
     793         238 :   return gerepilecopy(av, z);
     794             : }
     795             : 
     796             : static GEN
     797       69366 : FlxqXQ_Frobenius(GEN xp, GEN Xp, GEN f, GEN T, ulong p, ulong pi)
     798             : {
     799       69366 :   ulong dT = get_Flx_degree(T), df = get_FlxqX_degree(f);
     800       69366 :   GEN q = powuu(p,dT);
     801       69366 :   if (expi(q) >= expu(dT)*(long)usqrt(df))
     802       69324 :     return gel(FlxqXQ_autpow_pre(mkvec2(xp, Xp), dT, f, T, p, pi), 2);
     803             :   else
     804          42 :     return FlxqXQ_pow_pre(pol_x(get_FlxqX_var(f)), q, f, T, p, pi);
     805             : }
     806             : 
     807             : GEN
     808        9906 : FlxqX_Frobenius_pre(GEN S, GEN T, ulong p, ulong pi)
     809             : {
     810        9906 :   pari_sp av = avma;
     811        9906 :   GEN X  = polx_FlxX(get_FlxqX_var(S), get_Flx_var(T));
     812        9906 :   GEN xp = Flx_Frobenius_pre(T, p, pi);
     813        9906 :   GEN Xp = FlxqXQ_powu_pre(X, p, S, T, p, pi);
     814        9906 :   GEN Xq = FlxqXQ_Frobenius(xp, Xp, S, T, p, pi);
     815        9906 :   return gerepilecopy(av, Xq);
     816             : }
     817             : GEN
     818         678 : FlxqX_Frobenius(GEN S, GEN T, ulong p)
     819         678 : { return FlxqX_Frobenius_pre(S, T, p, SMALL_ULONG(p)? 0: get_Fl_red(p)); }
     820             : 
     821             : static GEN
     822         493 : FpXQXQ_Frobenius(GEN xp, GEN Xp, GEN f, GEN T, GEN p)
     823             : {
     824         493 :   ulong dT = get_FpX_degree(T), df = get_FpXQX_degree(f);
     825         493 :   GEN q = powiu(p, dT);
     826         493 :   if (expi(q) >= expu(dT)*(long)usqrt(df))
     827         493 :     return gel(FpXQXQ_autpow(mkvec2(xp, Xp), dT, f, T, p), 2);
     828             :   else
     829           0 :     return FpXQXQ_pow(pol_x(get_FpXQX_var(f)), q, f, T, p);
     830             : }
     831             : 
     832             : GEN
     833         388 : FpXQX_Frobenius(GEN S, GEN T, GEN p)
     834             : {
     835         388 :   pari_sp av = avma;
     836         388 :   GEN X  = pol_x(get_FpXQX_var(S));
     837         388 :   GEN xp = FpX_Frobenius(T, p);
     838         388 :   GEN Xp = FpXQXQ_pow(X, p, S, T, p);
     839         388 :   GEN Xq = FpXQXQ_Frobenius(xp, Xp, S, T, p);
     840         388 :   return gerepilecopy(av, Xq);
     841             : }
     842             : 
     843             : static GEN
     844       70441 : F2xqXQ_Frobenius(GEN xp, GEN Xp, GEN f, GEN T)
     845             : {
     846       70441 :   ulong dT = get_F2x_degree(T), df = get_F2xqX_degree(f);
     847       70441 :   if (dT >= expu(dT)*usqrt(df))
     848       70434 :     return gel(F2xqXQ_autpow(mkvec2(xp, Xp), dT, f, T), 2);
     849             :   else
     850             :   {
     851           7 :     long v = get_F2xqX_var(f), vT = get_F2x_var(T);
     852           7 :     return F2xqXQ_pow(polx_F2xX(v,vT), int2n(dT), f, T);
     853             :   }
     854             : }
     855             : 
     856             : GEN
     857           0 : F2xqX_Frobenius(GEN S, GEN T)
     858             : {
     859           0 :   pari_sp av = avma;
     860           0 :   GEN X  = polx_F2xX(get_F2xqX_var(S), get_F2x_var(T));
     861           0 :   GEN xp = F2x_Frobenius(T);
     862           0 :   GEN Xp = F2xqXQ_sqr(X, S, T);
     863           0 :   GEN Xq = F2xqXQ_Frobenius(xp, Xp, S, T);
     864           0 :   return gerepilecopy(av, Xq);
     865             : }
     866             : 
     867             : static GEN
     868           0 : F2xqX_split_part(GEN f, GEN T)
     869             : {
     870           0 :   long n = degpol(f);
     871             :   GEN z, Xq, X;
     872           0 :   if (n <= 1) return f;
     873           0 :   X = polx_F2xX(varn(f),get_F2x_var(T));
     874           0 :   f = F2xqX_red(f, T);
     875           0 :   Xq = F2xqX_Frobenius(f, T);
     876           0 :   z = F2xX_add(Xq, X);
     877           0 :   return F2xqX_gcd(z, f, T);
     878             : }
     879             : 
     880             : static GEN
     881        9039 : FlxqX_split_part(GEN f, GEN T, ulong p)
     882             : {
     883        9039 :   long n = degpol(f);
     884             :   GEN z, Xq, X;
     885             :   ulong pi;
     886        9039 :   if (n <= 1) return f;
     887        9039 :   pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
     888        9039 :   X = polx_FlxX(varn(f),get_Flx_var(T));
     889        9039 :   f = FlxqX_red_pre(f, T, p, pi);
     890        9039 :   Xq = FlxqX_Frobenius_pre(f, T, p, pi);
     891        9039 :   z = FlxX_sub(Xq, X , p);
     892        9039 :   return FlxqX_gcd_pre(z, f, T, p, pi);
     893             : }
     894             : 
     895             : GEN
     896        1736 : FpXQX_split_part(GEN f, GEN T, GEN p)
     897             : {
     898        1736 :   if(lgefint(p)==3)
     899             :   {
     900        1720 :     ulong pp=p[2];
     901        1720 :     GEN Tp = ZXT_to_FlxT(T, pp);
     902        1720 :     GEN z = FlxqX_split_part(ZXX_to_FlxX(f, pp, get_FpX_var(T)), Tp, pp);
     903        1720 :     return FlxX_to_ZXX(z);
     904             :   } else
     905             :   {
     906          16 :     long n = degpol(f);
     907          16 :     GEN z, X = pol_x(varn(f));
     908          16 :     if (n <= 1) return f;
     909          16 :     f = FpXQX_red(f, T, p);
     910          16 :     z = FpXQX_Frobenius(f, T, p);
     911          16 :     z = FpXX_sub(z, X , p);
     912          16 :     return FpXQX_gcd(z, f, T, p);
     913             :   }
     914             : }
     915             : 
     916             : long
     917        1680 : FpXQX_nbroots(GEN f, GEN T, GEN p)
     918             : {
     919        1680 :   pari_sp av = avma;
     920        1680 :   GEN z = FpXQX_split_part(f, T, p);
     921        1680 :   return gc_long(av, degpol(z));
     922             : }
     923             : 
     924             : long
     925      123004 : FqX_nbroots(GEN f, GEN T, GEN p)
     926      123004 : { return T ? FpXQX_nbroots(f, T, p): FpX_nbroots(f, p); }
     927             : 
     928             : long
     929           0 : F2xqX_nbroots(GEN f, GEN T)
     930             : {
     931           0 :   pari_sp av = avma;
     932           0 :   GEN z = F2xqX_split_part(f, T);
     933           0 :   return gc_long(av, degpol(z));
     934             : }
     935             : 
     936             : long
     937        7319 : FlxqX_nbroots(GEN f, GEN T, ulong p)
     938             : {
     939        7319 :   pari_sp av = avma;
     940        7319 :   GEN z = FlxqX_split_part(f, T, p);
     941        7319 :   return gc_long(av, degpol(z));
     942             : }
     943             : 
     944             : static GEN
     945           0 : FlxqX_Berlekamp_ker_i(GEN Xq, GEN S, GEN T, ulong p)
     946             : {
     947           0 :   long j, N = get_FlxqX_degree(S);
     948           0 :   GEN Q  = FlxqXQ_matrix_pow(Xq,N,N,S,T,p);
     949           0 :   for (j=1; j<=N; j++)
     950           0 :     gcoeff(Q,j,j) = Flx_Fl_add(gcoeff(Q,j,j), p-1, p);
     951           0 :   return FlxqM_ker(Q,T,p);
     952             : }
     953             : 
     954             : static GEN
     955           0 : FpXQX_Berlekamp_ker_i(GEN Xq, GEN S, GEN T, GEN p)
     956             : {
     957           0 :   long j,N = get_FpXQX_degree(S);
     958           0 :   GEN Q  = FpXQXQ_matrix_pow(Xq,N,N,S,T,p);
     959           0 :   for (j=1; j<=N; j++)
     960           0 :     gcoeff(Q,j,j) = Fq_sub(gcoeff(Q,j,j), gen_1, T, p);
     961           0 :   return FqM_ker(Q,T,p);
     962             : }
     963             : 
     964             : static long
     965        2661 : isabsolutepol(GEN f)
     966             : {
     967        2661 :   long i, l = lg(f);
     968        4524 :   for(i=2; i<l; i++)
     969             :   {
     970        4153 :     GEN c = gel(f,i);
     971        4153 :     if (typ(c) == t_POL && degpol(c) > 0) return 0;
     972             :   }
     973         371 :   return 1;
     974             : }
     975             : 
     976             : #define set_irred(i) { if ((i)>ir) swap(t[i],t[ir]); ir++;}
     977             : 
     978             : static long
     979           0 : FlxqX_split_Berlekamp(GEN *t, GEN xp, GEN T, ulong p, ulong pi)
     980             : {
     981           0 :   GEN u = *t, a,b,vker,pol;
     982           0 :   long vu = varn(u), vT = get_Flx_var(T), dT = get_Flx_degree(T);
     983             :   long d, i, ir, L, la, lb;
     984             :   GEN S, X, Xp, Xq;
     985           0 :   if (degpol(u)==1) return 1;
     986           0 :   T = Flx_get_red_pre(T, p, pi);
     987           0 :   S = FlxqX_get_red_pre(u, T, p, pi);
     988           0 :   X  = polx_FlxX(get_FlxqX_var(S),get_Flx_var(T));
     989           0 :   Xp = FlxqXQ_powu_pre(X, p, S, T, p, pi);
     990           0 :   Xq = FlxqXQ_Frobenius(xp, Xp, S, T, p, pi);
     991           0 :   vker = FlxqX_Berlekamp_ker_i(Xq, S, T, p);
     992           0 :   vker = Flm_to_FlxV(vker,u[1]);
     993           0 :   d = lg(vker)-1;
     994           0 :   ir = 0;
     995             :   /* t[i] irreducible for i < ir, still to be treated for i < L */
     996           0 :   for (L=1; L<d; )
     997             :   {
     998           0 :     pol= scalarpol(random_Flx(dT,vT,p),vu);
     999           0 :     for (i=2; i<=d; i++)
    1000           0 :       pol = FlxX_add(pol, FlxqX_Flxq_mul_pre(gel(vker,i),
    1001             :                                              random_Flx(dT,vT,p), T, p, pi), p);
    1002           0 :     pol = FlxqX_red_pre(pol,T,p,pi);
    1003           0 :     for (i=ir; i<L && L<d; i++)
    1004             :     {
    1005           0 :       a = t[i]; la = degpol(a);
    1006           0 :       if (la == 1) { set_irred(i); }
    1007             :       else
    1008             :       {
    1009           0 :         pari_sp av = avma;
    1010           0 :         GEN S = FlxqX_get_red_pre(a, T, p, pi);
    1011           0 :         b = FlxqX_rem_pre(pol, S, T,p,pi);
    1012           0 :         if (degpol(b)<=0) { set_avma(av); continue; }
    1013           0 :         b = FlxqXQ_halfFrobenius_i(b, xp, FlxqX_rem_pre(Xp,S,T,p,pi), S,T,p,pi);
    1014           0 :         if (degpol(b)<=0) { set_avma(av); continue; }
    1015           0 :         gel(b,2) = Flxq_sub(gel(b,2), gen_1,T,p);
    1016           0 :         b = FlxqX_gcd_pre(a,b, T,p,pi); lb = degpol(b);
    1017           0 :         if (lb && lb < la)
    1018             :         {
    1019           0 :           b = FlxqX_normalize_pre(b, T,p,pi);
    1020           0 :           t[L] = FlxqX_div_pre(a,b,T,p,pi);
    1021           0 :           t[i]= b; L++;
    1022             :         }
    1023           0 :         else set_avma(av);
    1024             :       }
    1025             :     }
    1026             :   }
    1027           0 :   return d;
    1028             : }
    1029             : 
    1030             : static long
    1031           0 : FpXQX_split_Berlekamp(GEN *t, GEN T, GEN p)
    1032             : {
    1033           0 :   GEN u = *t, a, b, vker, pol;
    1034             :   GEN X, xp, Xp, Xq, S;
    1035           0 :   long vu = varn(u), vT = get_FpX_var(T), dT = get_FpX_degree(T);
    1036             :   long d, i, ir, L, la, lb;
    1037           0 :   if (degpol(u)==1) return 1;
    1038           0 :   T = FpX_get_red(T, p);
    1039           0 :   xp = FpX_Frobenius(T, p);
    1040           0 :   S = FpXQX_get_red(u, T, p);
    1041           0 :   X  = pol_x(get_FpXQX_var(S));
    1042           0 :   Xp = FpXQXQ_pow(X, p, S, T, p);
    1043           0 :   Xq = FpXQXQ_Frobenius(xp, Xp, S, T, p);
    1044           0 :   vker = FpXQX_Berlekamp_ker_i(Xq, S, T, p);
    1045           0 :   vker = RgM_to_RgXV(vker,vu);
    1046           0 :   d = lg(vker)-1;
    1047           0 :   ir = 0;
    1048             :   /* t[i] irreducible for i < ir, still to be treated for i < L */
    1049           0 :   for (L=1; L<d; )
    1050             :   {
    1051           0 :     pol= scalarpol(random_FpX(dT,vT,p),vu);
    1052           0 :     for (i=2; i<=d; i++)
    1053           0 :       pol = FqX_add(pol, FqX_Fq_mul(gel(vker,i),
    1054             :                                     random_FpX(dT,vT,p), T, p), T, p);
    1055           0 :     pol = FpXQX_red(pol,T,p);
    1056           0 :     for (i=ir; i<L && L<d; i++)
    1057             :     {
    1058           0 :       a = t[i]; la = degpol(a);
    1059           0 :       if (la == 1) { set_irred(i); }
    1060             :       else
    1061             :       {
    1062           0 :         pari_sp av = avma;
    1063           0 :         GEN S = FpXQX_get_red(a, T, p);
    1064           0 :         b = FqX_rem(pol, S, T,p);
    1065           0 :         if (degpol(b)<=0) { set_avma(av); continue; }
    1066           0 :         b = FpXQXQ_halfFrobenius_i(b, xp, FpXQX_rem(Xp, S, T, p), S, T, p);
    1067           0 :         if (degpol(b)<=0) { set_avma(av); continue; }
    1068           0 :         gel(b,2) = Fq_sub(gel(b,2), gen_1,T,p);
    1069           0 :         b = FqX_gcd(a,b, T,p); lb = degpol(b);
    1070           0 :         if (lb && lb < la)
    1071             :         {
    1072           0 :           b = FpXQX_normalize(b, T,p);
    1073           0 :           t[L] = FqX_div(a,b,T,p);
    1074           0 :           t[i]= b; L++;
    1075             :         }
    1076           0 :         else set_avma(av);
    1077             :       }
    1078             :     }
    1079             :   }
    1080           0 :   return d;
    1081             : }
    1082             : 
    1083             : static GEN
    1084       11501 : F2xqX_quad_roots(GEN P, GEN T)
    1085             : {
    1086       11501 :   GEN b= gel(P,3), c = gel(P,2);
    1087       11501 :   if (lgpol(b))
    1088             :   {
    1089       10227 :     GEN z, d = F2xq_div(c, F2xq_sqr(b,T),T);
    1090       10227 :     if (F2xq_trace(d,T))
    1091        1001 :       return cgetg(1, t_COL);
    1092        9226 :     z = F2xq_mul(b, F2xq_Artin_Schreier(d, T), T);
    1093        9226 :     return mkcol2(z, F2x_add(b, z));
    1094             :   }
    1095             :   else
    1096        1274 :     return mkcol(F2xq_sqrt(c, T));
    1097             : }
    1098             : 
    1099             : /* Assume p>2 and x monic */
    1100             : static GEN
    1101       14684 : FlxqX_quad_roots(GEN x, GEN T, ulong p, ulong pi)
    1102             : {
    1103       14684 :   GEN s, D, nb, b = gel(x,3), c = gel(x,2);
    1104       14684 :   D = Flx_sub(Flxq_sqr_pre(b,T,p,pi), Flx_mulu(c,4,p), p);
    1105       14684 :   nb = Flx_neg(b,p);
    1106       14684 :   if (lgpol(D)==0) return mkcol(Flx_halve(nb, p));
    1107       14649 :   s = Flxq_sqrt(D,T,p);
    1108       14649 :   if (!s) return cgetg(1, t_COL);
    1109       14187 :   s = Flx_halve(Flx_add(s,nb,p),p);
    1110       14187 :   return mkcol2(s, Flx_sub(nb,s,p));
    1111             : }
    1112             : 
    1113             : static GEN
    1114         839 : FpXQX_quad_roots(GEN x, GEN T, GEN p)
    1115             : {
    1116         839 :   GEN s, D, nb, b = gel(x,3), c = gel(x,2);
    1117         839 :   if (absequaliu(p, 2))
    1118             :   {
    1119           0 :     GEN f2 = ZXX_to_F2xX(x, get_FpX_var(T));
    1120           0 :     s = F2xqX_quad_roots(f2, ZX_to_F2x(get_FpX_mod(T)));
    1121           0 :     return F2xC_to_ZXC(s);
    1122             :   }
    1123         839 :   D = Fq_sub(Fq_sqr(b,T,p), Fq_Fp_mul(c,utoi(4),T,p), T,p);
    1124         839 :   nb = Fq_neg(b,T,p);
    1125         839 :   if (signe(D)==0)
    1126           0 :     return mkcol(Fq_to_FpXQ(Fq_halve(nb,T, p),T,p));
    1127         839 :   s = Fq_sqrt(D,T,p);
    1128         839 :   if (!s) return cgetg(1, t_COL);
    1129         825 :   s = Fq_halve(Fq_add(s,nb,T,p),T, p);
    1130         825 :   return mkcol2(Fq_to_FpXQ(s,T,p), Fq_to_FpXQ(Fq_sub(nb,s,T,p),T,p));
    1131             : }
    1132             : 
    1133             : static GEN
    1134        9247 : F2xqX_Frobenius_deflate(GEN S, GEN T)
    1135             : {
    1136        9247 :   GEN F = RgX_deflate(S, 2);
    1137        9247 :   long i, l = lg(F);
    1138       33089 :   for (i=2; i<l; i++)
    1139       23842 :     gel(F,i) = F2xq_sqrt(gel(F,i), T);
    1140        9247 :   return F;
    1141             : }
    1142             : 
    1143             : static GEN
    1144       16940 : F2xX_to_F2x(GEN x)
    1145             : {
    1146       16940 :   long l=nbits2lg(lgpol(x));
    1147       16940 :   GEN z=cgetg(l,t_VECSMALL);
    1148             :   long i,j,k;
    1149       16940 :   z[1]=x[1];
    1150       63980 :   for(i=2, k=1,j=BITS_IN_LONG;i<lg(x);i++,j++)
    1151             :   {
    1152       47040 :     if (j==BITS_IN_LONG)
    1153             :     {
    1154       16967 :       j=0; k++; z[k]=0;
    1155             :     }
    1156       47040 :     if (lgpol(gel(x,i)))
    1157       34517 :       z[k]|=1UL<<j;
    1158             :   }
    1159       16940 :   return F2x_renormalize(z,l);
    1160             : }
    1161             : 
    1162             : static GEN
    1163      205807 : F2xqX_easyroots(GEN f, GEN T)
    1164             : {
    1165      205807 :   if (F2xY_degreex(f) <= 0) return F2x_rootsff_i(F2xX_to_F2x(f), T);
    1166      189287 :   if (degpol(f)==1) return mkcol(constant_coeff(f));
    1167      154119 :   if (degpol(f)==2) return F2xqX_quad_roots(f,T);
    1168      143150 :   return NULL;
    1169             : }
    1170             : 
    1171             : /* Adapted from Shoup NTL */
    1172             : GEN
    1173       71204 : F2xqX_factor_squarefree(GEN f, GEN T)
    1174             : {
    1175       71204 :   pari_sp av = avma;
    1176             :   GEN r, t, v, tv;
    1177       71204 :   long i, q, n = degpol(f);
    1178       71204 :   GEN u = const_vec(n+1, pol1_F2xX(varn(f), get_F2x_var(T)));
    1179       71204 :   for(q = 1;;q *= 2)
    1180             :   {
    1181       80451 :     r = F2xqX_gcd(f, F2xX_deriv(f), T);
    1182       80451 :     if (degpol(r) == 0)
    1183             :     {
    1184       69531 :       gel(u, q) = F2xqX_normalize(f, T);
    1185       69531 :       break;
    1186             :     }
    1187       10920 :     t = F2xqX_div(f, r, T);
    1188       10920 :     if (degpol(t) > 0)
    1189             :     {
    1190             :       long j;
    1191       10024 :       for(j = 1;;j++)
    1192             :       {
    1193       14924 :         v = F2xqX_gcd(r, t, T);
    1194       14924 :         tv = F2xqX_div(t, v, T);
    1195       14924 :         if (degpol(tv) > 0)
    1196       11690 :           gel(u, j*q) = F2xqX_normalize(tv, T);
    1197       14924 :         if (degpol(v) <= 0) break;
    1198        4900 :         r = F2xqX_div(r, v, T);
    1199        4900 :         t = v;
    1200             :       }
    1201       10024 :       if (degpol(r) == 0) break;
    1202             :     }
    1203        9247 :     f = F2xqX_Frobenius_deflate(r, T);
    1204             :   }
    1205      417172 :   for (i = n; i; i--)
    1206      417172 :     if (degpol(gel(u,i))) break;
    1207       71204 :   setlg(u,i+1); return gerepilecopy(av, u);
    1208             : }
    1209             : 
    1210             : long
    1211          56 : F2xqX_ispower(GEN f, long k, GEN T, GEN *pt_r)
    1212             : {
    1213          56 :   pari_sp av = avma;
    1214             :   GEN lc, F;
    1215          56 :   long i, l, n = degpol(f);
    1216          56 :   if (n % k) return 0;
    1217          56 :   lc = F2xq_sqrtn(leading_coeff(f), stoi(k), T, NULL);
    1218          56 :   if (!lc) return gc_long(av,0);
    1219          56 :   F = F2xqX_factor_squarefree(f, T); l = lg(F)-1;
    1220        2030 :   for(i=1; i<=l; i++)
    1221        1981 :     if (i%k && degpol(gel(F,i))) return gc_long(av,0);
    1222          49 :   if (pt_r)
    1223             :   {
    1224          49 :     long v = varn(f);
    1225          49 :     GEN r = scalarpol(lc, v), s = pol1_F2xX(v, T[1]);
    1226        2023 :     for(i=l; i>=1; i--)
    1227             :     {
    1228        1974 :       if (i%k) continue;
    1229         406 :       s = F2xqX_mul(s, gel(F,i), T);
    1230         406 :       r = F2xqX_mul(r, s, T);
    1231             :     }
    1232          49 :     *pt_r = gerepileupto(av, r);
    1233           0 :   } else set_avma(av);
    1234          49 :   return 1;
    1235             : }
    1236             : 
    1237             : static void
    1238       49889 : F2xqX_roots_edf(GEN Sp, GEN xp, GEN Xp, GEN T, GEN V, long idx)
    1239             : {
    1240             :   pari_sp btop;
    1241       49889 :   long n = degpol(Sp);
    1242             :   GEN S, f, ff;
    1243       49889 :   long dT = get_F2x_degree(T);
    1244       49889 :   GEN R = F2xqX_easyroots(Sp, T);
    1245       49889 :   if (R)
    1246             :   {
    1247       47754 :     long i, l = lg(R)-1;
    1248      105987 :     for (i=0; i<l; i++)
    1249       58233 :       gel(V, idx+i) = gel(R,1+i);
    1250       47754 :     return;
    1251             :   }
    1252        2135 :   S = F2xqX_get_red(Sp, T);
    1253        2135 :   Xp = F2xqX_rem(Xp, S, T);
    1254        2135 :   btop = avma;
    1255             :   while (1)
    1256         588 :   {
    1257        2723 :     GEN a = random_F2xqX(degpol(Sp), varn(Sp), T);
    1258        2723 :     GEN R = gel(F2xqXQ_auttrace(mkvec3(xp, Xp, a), dT, S, T), 3);
    1259        2723 :     f = F2xqX_gcd(R, Sp, T);
    1260        2723 :     if (degpol(f) > 0 && degpol(f) < n) break;
    1261         588 :     set_avma(btop);
    1262             :   }
    1263        2135 :   f = gerepileupto(btop, F2xqX_normalize(f, T));
    1264        2135 :   ff = F2xqX_div(Sp, f, T);
    1265        2135 :   F2xqX_roots_edf(f, xp, Xp, T, V, idx);
    1266        2135 :   F2xqX_roots_edf(ff,xp, Xp, T, V, idx+degpol(f));
    1267             : }
    1268             : 
    1269             : static GEN
    1270       80787 : F2xqX_roots_ddf(GEN f, GEN xp, GEN T)
    1271             : {
    1272             :   GEN X, Xp, Xq, g, V;
    1273             :   long n;
    1274       80787 :   GEN R = F2xqX_easyroots(f, T);
    1275       80787 :   if (R) return R;
    1276       70147 :   X  = pol_x(varn(f));
    1277       70147 :   Xp = F2xqXQ_sqr(X, f, T);
    1278       70147 :   Xq = F2xqXQ_Frobenius(xp, Xp, f, T);
    1279       70147 :   g = F2xqX_gcd(F2xX_add(Xq, X), f, T);
    1280       70147 :   n = degpol(g);
    1281       70147 :   if (n==0) return cgetg(1, t_COL);
    1282       45619 :   g = F2xqX_normalize(g, T);
    1283       45619 :   V = cgetg(n+1,t_COL);
    1284       45619 :   F2xqX_roots_edf(g, xp, Xp, T, V, 1);
    1285       45619 :   return V;
    1286             : }
    1287             : static GEN
    1288       75138 : F2xqX_roots_i(GEN S, GEN T)
    1289             : {
    1290             :   GEN M;
    1291       75138 :   S = F2xqX_red(S, T);
    1292       75138 :   if (!signe(S)) pari_err_ROOTS0("F2xqX_roots");
    1293       75138 :   if (degpol(S)==0) return cgetg(1, t_COL);
    1294       75131 :   S = F2xqX_normalize(S, T);
    1295       75131 :   M = F2xqX_easyroots(S, T);
    1296       75131 :   if (!M)
    1297             :   {
    1298       70868 :     GEN xp = F2x_Frobenius(T);
    1299       70868 :     GEN F, V = F2xqX_factor_squarefree(S, T);
    1300       70868 :     long i, j, l = lg(V);
    1301       70868 :     F = cgetg(l, t_VEC);
    1302      154567 :     for (i=1, j=1; i < l; i++)
    1303       83699 :       if (degpol(gel(V,i)))
    1304       80787 :         gel(F, j++) = F2xqX_roots_ddf(gel(V,i), xp, T);
    1305       70868 :     setlg(F,j); M = shallowconcat1(F);
    1306             :   }
    1307       75131 :   gen_sort_inplace(M, (void*) &cmp_Flx, &cmp_nodata, NULL);
    1308       75131 :   return M;
    1309             : }
    1310             : 
    1311             : static GEN
    1312      183795 : FlxqX_easyroots(GEN f, GEN T, ulong p, ulong pi)
    1313             : {
    1314      183795 :   if (FlxY_degreex(f) <= 0) return Flx_rootsff_i(FlxX_to_Flx(f), T, p);
    1315      170389 :   if (degpol(f)==1) return mkcol(Flx_neg(constant_coeff(f), p));
    1316      141585 :   if (degpol(f)==2) return FlxqX_quad_roots(f,T,p,pi);
    1317      127726 :   return NULL;
    1318             : }
    1319             : 
    1320             : static GEN
    1321         574 : FlxqX_invFrobenius(GEN xp, GEN T, ulong p, ulong pi)
    1322         574 : { return Flxq_autpow_pre(xp, get_Flx_degree(T)-1, T, p, pi); }
    1323             : 
    1324             : static GEN
    1325         637 : FlxqX_Frobenius_deflate(GEN S, GEN ixp, GEN T, ulong p)
    1326             : {
    1327         637 :   GEN F = RgX_deflate(S, p);
    1328         637 :   long i, l = lg(F);
    1329         637 :   if (typ(ixp)==t_INT)
    1330           0 :     for (i=2; i<l; i++)
    1331           0 :       gel(F,i) = Flxq_pow(gel(F,i), ixp, T, p);
    1332             :   else
    1333             :   {
    1334         637 :     long n = brent_kung_optpow(get_Flx_degree(T)-1, l-2, 1);
    1335         637 :     GEN V = Flxq_powers(ixp, n, T, p);
    1336        5719 :     for (i=2; i<l; i++)
    1337        5082 :       gel(F,i) = Flx_FlxqV_eval(gel(F,i), V, T, p);
    1338             :   }
    1339         637 :   return F;
    1340             : }
    1341             : 
    1342             : /* Adapted from Shoup NTL */
    1343             : static GEN
    1344       59961 : FlxqX_factor_squarefree_i(GEN f, GEN xp, GEN T, ulong p, ulong pi)
    1345             : {
    1346       59961 :   pari_sp av = avma;
    1347       59961 :   long i, q, n = degpol(f);
    1348       59961 :   GEN u = const_vec(n+1, pol1_FlxX(varn(f),get_Flx_var(T)));
    1349       59961 :   GEN r, t, v, tv, ixp = NULL;
    1350       59961 :   for(q = 1;; q *= p)
    1351             :   {
    1352       60598 :     r = FlxqX_gcd_pre(f, FlxX_deriv(f, p), T, p, pi);
    1353       60598 :     if (degpol(r) == 0)
    1354             :     {
    1355       55743 :       gel(u, q) = FlxqX_normalize_pre(f, T, p, pi);
    1356       55743 :       break;
    1357             :     }
    1358        4855 :     t = FlxqX_div_pre(f, r, T, p, pi);
    1359        4855 :     if (degpol(t) > 0)
    1360             :     {
    1361             :       long j;
    1362        4645 :       for(j = 1;;j++)
    1363             :       {
    1364       10123 :         v = FlxqX_gcd_pre(r, t, T, p, pi);
    1365       10123 :         tv = FlxqX_div_pre(t, v, T, p, pi);
    1366       10123 :         if (degpol(tv) > 0)
    1367        8751 :           gel(u, j*q) = FlxqX_normalize_pre(tv, T, p, pi);
    1368       10123 :         if (degpol(v) <= 0) break;
    1369        5478 :         r = FlxqX_div_pre(r, v, T, p, pi);
    1370        5478 :         t = v;
    1371             :       }
    1372        4645 :       if (degpol(r) == 0) break;
    1373             :     }
    1374         637 :     if (!xp)   xp = Flx_Frobenius_pre(T, p, pi);
    1375         637 :     if (!ixp) ixp = FlxqX_invFrobenius(xp, T, p, pi);
    1376         637 :     f = FlxqX_Frobenius_deflate(r, ixp, T, p);
    1377             :   }
    1378      323906 :   for (i = n; i; i--)
    1379      323906 :     if (degpol(gel(u,i))) break;
    1380       59961 :   setlg(u,i+1); return gerepilecopy(av, u);
    1381             : }
    1382             : 
    1383             : GEN
    1384          42 : FlxqX_factor_squarefree_pre(GEN f, GEN T, ulong p, ulong pi)
    1385          42 : { return FlxqX_factor_squarefree_i(f, NULL, T, p, pi); }
    1386             : GEN
    1387           0 : FlxqX_factor_squarefree(GEN f, GEN T, ulong p)
    1388           0 : { return FlxqX_factor_squarefree_pre(f,T,p, SMALL_ULONG(p)? 0: get_Fl_red(p)); }
    1389             : 
    1390             : long
    1391          98 : FlxqX_ispower(GEN f, ulong k, GEN T, ulong p, GEN *pt_r)
    1392             : {
    1393          98 :   pari_sp av = avma;
    1394             :   GEN lc, F;
    1395          98 :   long i, l, n = degpol(f), v = varn(f);
    1396             :   ulong pi;
    1397             : 
    1398          98 :   if (n % k) return 0;
    1399          98 :   lc = Flxq_sqrtn(leading_coeff(f), stoi(k), T, p, NULL);
    1400          98 :   if (!lc) return gc_long(av,0);
    1401          98 :   pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
    1402          98 :   F = FlxqX_factor_squarefree_i(f, NULL, T, p, pi); l = lg(F)-1;
    1403        3521 :   for(i=1; i<=l; i++)
    1404        3437 :     if (i%k && degpol(gel(F,i))) return gc_long(av,0);
    1405          84 :   if (pt_r)
    1406             :   {
    1407          84 :     GEN r = scalarpol(lc, v), s = pol1_FlxX(v, T[1]);
    1408        3507 :     for(i=l; i>=1; i--)
    1409             :     {
    1410        3423 :       if (i%k) continue;
    1411         700 :       s = FlxqX_mul_pre(s, gel(F,i), T, p, pi);
    1412         700 :       r = FlxqX_mul_pre(r, s, T, p, pi);
    1413             :     }
    1414          84 :     *pt_r = gerepileupto(av, r);
    1415           0 :   } else set_avma(av);
    1416          84 :   return 1;
    1417             : }
    1418             : 
    1419             : static GEN
    1420        9442 : FlxqX_roots_split(GEN Sp, GEN xp, GEN Xp, GEN S, GEN T, ulong p, long pi)
    1421             : {
    1422        9442 :   pari_sp btop = avma;
    1423        9442 :   long n = degpol(Sp);
    1424             :   GEN f;
    1425        9442 :   long vT = get_Flx_var(T), dT = get_Flx_degree(T);
    1426             :   pari_timer ti;
    1427        9442 :   if (DEBUGLEVEL >= 7) timer_start(&ti);
    1428             :   while (1)
    1429        2102 :   {
    1430       11544 :     GEN a = deg1pol(pol1_Flx(vT), random_Flx(dT, vT, p), varn(Sp));
    1431       11544 :     GEN R = FlxqXQ_halfFrobenius_i(a, xp, Xp, S, T, p, pi);
    1432       11544 :     if (DEBUGLEVEL >= 7) timer_printf(&ti, "FlxqXQ_halfFrobenius");
    1433       11544 :     f = FlxqX_gcd_pre(FlxX_Flx_sub(R, pol1_Flx(vT), p), Sp, T, p, pi);
    1434       11544 :     if (degpol(f) > 0 && degpol(f) < n) break;
    1435        2102 :     set_avma(btop);
    1436             :   }
    1437        9442 :   return gerepileupto(btop, FlxqX_normalize_pre(f, T, p, pi));
    1438             : }
    1439             : 
    1440             : static void
    1441       53422 : FlxqX_roots_edf(GEN Sp, GEN xp, GEN Xp, GEN T, ulong p, ulong pi, GEN V, long idx)
    1442             : {
    1443             :   GEN S, f, ff;
    1444       53422 :   GEN R = FlxqX_easyroots(Sp, T, p, pi);
    1445       53422 :   if (R)
    1446             :   {
    1447       44342 :     long i, l = lg(R)-1;
    1448      102484 :     for (i=0; i<l; i++)
    1449       58142 :       gel(V, idx+i) = gel(R,1+i);
    1450       44342 :     return;
    1451             :   }
    1452        9080 :   S  = FlxqX_get_red_pre(Sp, T, p, pi);
    1453        9080 :   Xp = FlxqX_rem_pre(Xp, S, T, p, pi);
    1454        9080 :   f  = FlxqX_roots_split(Sp, xp, Xp, S, T, p, pi);
    1455        9080 :   ff = FlxqX_div_pre(Sp, f, T, p, pi);
    1456        9080 :   FlxqX_roots_edf(f, xp, Xp, T, p, pi, V, idx);
    1457        9080 :   FlxqX_roots_edf(ff,xp, Xp, T, p, pi, V, idx+degpol(f));
    1458             : }
    1459             : 
    1460             : static GEN
    1461       63886 : FlxqX_roots_ddf(GEN f, GEN xp, GEN T, ulong p, ulong pi)
    1462             : {
    1463             :   GEN X, Xp, Xq, g, V;
    1464             :   long n;
    1465       63886 :   GEN R = FlxqX_easyroots(f, T, p, pi);
    1466       63886 :   if (R) return R;
    1467       59118 :   X  = pol_x(varn(f));
    1468       59118 :   Xp = FlxqXQ_powu_pre(X, p, f, T, p, pi);
    1469       59118 :   Xq = FlxqXQ_Frobenius(xp, Xp, f, T, p, pi);
    1470       59118 :   g = FlxqX_gcd_pre(FlxX_sub(Xq, X, p), f, T, p, pi);
    1471       59118 :   n = degpol(g);
    1472       59118 :   if (n==0) return cgetg(1, t_COL);
    1473       35262 :   g = FlxqX_normalize_pre(g, T, p, pi);
    1474       35262 :   V = cgetg(n+1,t_COL);
    1475       35262 :   FlxqX_roots_edf(g, xp, Xp, T, p, pi, V, 1);
    1476       35262 :   return V;
    1477             : }
    1478             : 
    1479             : /* do not handle p==2 */
    1480             : static GEN
    1481       66494 : FlxqX_roots_i(GEN S, GEN T, ulong p)
    1482             : {
    1483       66494 :   ulong pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
    1484             :   GEN M;
    1485       66494 :   S = FlxqX_red_pre(S, T, p, pi);
    1486       66494 :   if (!signe(S)) pari_err_ROOTS0("FlxqX_roots");
    1487       66494 :   if (degpol(S)==0) return cgetg(1, t_COL);
    1488       66487 :   S = FlxqX_normalize_pre(S, T, p, pi);
    1489       66487 :   M = FlxqX_easyroots(S, T, p, pi);
    1490       66487 :   if (!M)
    1491             :   {
    1492       59528 :     GEN xp = Flx_Frobenius_pre(T, p, pi);
    1493       59528 :     GEN F, V = FlxqX_factor_squarefree_i(S, xp, T, p, pi);
    1494       59528 :     long i, j, l = lg(V);
    1495       59528 :     F = cgetg(l, t_VEC);
    1496      124051 :     for (i=1, j=1; i < l; i++)
    1497       64523 :       if (degpol(gel(V,i)))
    1498       63886 :         gel(F, j++) = FlxqX_roots_ddf(gel(V,i), xp, T, p, pi);
    1499       59528 :     setlg(F,j); M = shallowconcat1(F);
    1500             :   }
    1501       66487 :   gen_sort_inplace(M, (void*) &cmp_Flx, &cmp_nodata, NULL);
    1502       66487 :   return M;
    1503             : }
    1504             : 
    1505             : static GEN
    1506        2576 : FpXQX_easyroots(GEN f, GEN T, GEN p)
    1507             : {
    1508        2576 :   if (isabsolutepol(f)) return FpX_rootsff_i(simplify_shallow(f), T, p);
    1509        2268 :   if (degpol(f)==1) return mkcol(Fq_to_FpXQ(Fq_neg(constant_coeff(f),T,p),T,p));
    1510        1867 :   if (degpol(f)==2) return FpXQX_quad_roots(f,T,p);
    1511        1071 :   return NULL;
    1512             : }
    1513             : 
    1514             : /* Adapted from Shoup NTL */
    1515             : static GEN
    1516         210 : FpXQX_factor_Yun(GEN f, GEN T, GEN p)
    1517             : {
    1518         210 :   pari_sp av = avma;
    1519             :   GEN r, t, v, tv;
    1520         210 :   long j, n = degpol(f);
    1521         210 :   GEN u = const_vec(n+1, pol_1(varn(f)));
    1522         210 :   r = FpXQX_gcd(f, FpXX_deriv(f, p), T, p);
    1523         210 :   t = FpXQX_div(f, r, T, p);
    1524         210 :   for (j = 1;;j++)
    1525             :   {
    1526        1652 :     v = FpXQX_gcd(r, t, T, p);
    1527        1652 :     tv = FpXQX_div(t, v, T, p);
    1528        1652 :     if (degpol(tv) > 0)
    1529         252 :       gel(u, j) = FpXQX_normalize(tv, T, p);
    1530        1652 :     if (degpol(v) <= 0) break;
    1531        1442 :     r = FpXQX_div(r, v, T, p);
    1532        1442 :     t = v;
    1533             :   }
    1534         210 :   setlg(u, j+1); return gerepilecopy(av, u);
    1535             : }
    1536             : 
    1537             : GEN
    1538           7 : FpXQX_factor_squarefree(GEN f, GEN T, GEN p)
    1539             : {
    1540           7 :   if (lgefint(p)==3)
    1541             :   {
    1542           0 :     pari_sp av = avma;
    1543           0 :     ulong pp = p[2];
    1544             :     GEN M;
    1545           0 :     long vT = get_FpX_var(T);
    1546           0 :     if (pp==2)
    1547             :     {
    1548           0 :       M = F2xqX_factor_squarefree(ZXX_to_F2xX(f, vT),  ZX_to_F2x(get_FpX_mod(T)));
    1549           0 :       return gerepileupto(av, F2xXC_to_ZXXC(M));
    1550             :     }
    1551           0 :     M = FlxqX_factor_squarefree(ZXX_to_FlxX(f, pp, vT),  ZXT_to_FlxT(T, pp), pp);
    1552           0 :     return gerepileupto(av, FlxXC_to_ZXXC(M));
    1553             :   }
    1554           7 :   return FpXQX_factor_Yun(f, T, p);
    1555             : }
    1556             : 
    1557             : GEN
    1558           0 : FpXQX_roots_mult(GEN f, long n, GEN T, GEN p)
    1559             : {
    1560           0 :   pari_sp av = avma;
    1561           0 :   GEN V = FpXQX_factor_squarefree(f, T, p), W;
    1562           0 :   long l = lg(V), i;
    1563           0 :   if (l <= n) return cgetg(1,t_COL);
    1564           0 :   W = cgetg(l-n+1,t_VEC);
    1565           0 :   for (i = n; i < l; i++)
    1566           0 :     gel(W,i-n+1) = FpXQX_roots(gel(V,i), T, p);
    1567           0 :   W = shallowconcat1(W);
    1568           0 :   gen_sort_inplace(W, (void*) &cmp_RgX, &cmp_nodata, NULL);
    1569           0 :   return gerepilecopy(av, W);
    1570             : }
    1571             : 
    1572             : long
    1573          98 : FpXQX_ispower(GEN f, ulong k, GEN T, GEN p, GEN *pt_r)
    1574             : {
    1575          98 :   pari_sp av = avma;
    1576             :   GEN lc, F;
    1577          98 :   long i, l, n = degpol(f), v = varn(f);
    1578          98 :   if (n % k) return 0;
    1579          98 :   if (lgefint(p)==3)
    1580             :   {
    1581          42 :     ulong pp = p[2];
    1582          42 :     GEN fp = ZXX_to_FlxX(f, pp, varn(T));
    1583          42 :     if (!FlxqX_ispower(fp, k, ZX_to_Flx(T,pp), pp, pt_r)) return gc_long(av,0);
    1584          35 :     if (pt_r) *pt_r = gerepileupto(av, FlxX_to_ZXX(*pt_r));
    1585           0 :     else set_avma(av);
    1586          35 :     return 1;
    1587             :   }
    1588          56 :   lc = FpXQ_sqrtn(leading_coeff(f), stoi(k), T, p, NULL);
    1589          56 :   if (!lc) return gc_long(av,0);
    1590          56 :   F = FpXQX_factor_Yun(f, T, p); l = lg(F)-1;
    1591        1533 :   for(i=1; i <= l; i++)
    1592        1484 :     if (i%k && degpol(gel(F,i))) return gc_long(av,0);
    1593          49 :   if (pt_r)
    1594             :   {
    1595          49 :     GEN r = scalarpol(lc, v), s = pol_1(v);
    1596        1526 :     for(i=l; i>=1; i--)
    1597             :     {
    1598        1477 :       if (i%k) continue;
    1599         308 :       s = FpXQX_mul(s, gel(F,i), T, p);
    1600         308 :       r = FpXQX_mul(r, s, T, p);
    1601             :     }
    1602          49 :     *pt_r = gerepileupto(av, r);
    1603           0 :   } else set_avma(av);
    1604          49 :   return 1;
    1605             : }
    1606             : 
    1607             : long
    1608         210 : FqX_ispower(GEN f, ulong k, GEN T, GEN p, GEN *pt_r)
    1609         210 : { return T ? FpXQX_ispower(f, k, T, p, pt_r): FpX_ispower(f, k, p, pt_r); }
    1610             : 
    1611             : static GEN
    1612         928 : FpXQX_roots_split(GEN Sp, GEN xp, GEN Xp, GEN S, GEN T, GEN p)
    1613             : {
    1614         928 :   pari_sp btop = avma;
    1615         928 :   long n = degpol(Sp);
    1616             :   GEN f;
    1617         928 :   long vT = get_FpX_var(T), dT = get_FpX_degree(T);
    1618             :   pari_timer ti;
    1619         928 :   if (DEBUGLEVEL >= 7) timer_start(&ti);
    1620             :   while (1)
    1621         127 :   {
    1622        1055 :     GEN a = deg1pol(pol_1(vT), random_FpX(dT, vT, p), varn(Sp));
    1623        1055 :     GEN R = FpXQXQ_halfFrobenius_i(a, xp, Xp, S, T, p);
    1624        1055 :     if (DEBUGLEVEL >= 7) timer_printf(&ti, "FpXQXQ_halfFrobenius");
    1625        1055 :     f = FpXQX_gcd(FqX_Fq_sub(R, pol_1(vT), T, p), Sp, T, p);
    1626        1055 :     if (degpol(f) > 0 && degpol(f) < n) break;
    1627         127 :     set_avma(btop);
    1628             :   }
    1629         928 :   return gerepileupto(btop, FpXQX_normalize(f, T, p));
    1630             : }
    1631             : 
    1632             : static void
    1633        1893 : FpXQX_roots_edf(GEN Sp, GEN xp, GEN Xp, GEN T, GEN p, GEN V, long idx)
    1634             : {
    1635             :   GEN S, f, ff;
    1636        1893 :   GEN R = FpXQX_easyroots(Sp, T, p);
    1637        1893 :   if (R)
    1638             :   {
    1639         988 :     long i, l = lg(R)-1;
    1640        2576 :     for (i=0; i<l; i++)
    1641        1588 :       gel(V, idx+i) = gel(R,1+i);
    1642         988 :     return;
    1643             :   }
    1644         905 :   S  = FpXQX_get_red(Sp, T, p);
    1645         905 :   Xp = FpXQX_rem(Xp, S, T, p);
    1646         905 :   f  = FpXQX_roots_split(Sp, xp, Xp, S, T, p);
    1647         905 :   ff = FpXQX_div(Sp, f, T, p);
    1648         905 :   FpXQX_roots_edf(f, xp, Xp, T, p, V, idx);
    1649         905 :   FpXQX_roots_edf(ff,xp, Xp, T, p, V, idx+degpol(f));
    1650             : }
    1651             : 
    1652             : static GEN
    1653          83 : FpXQX_roots_ddf(GEN f, GEN xp, GEN T, GEN p)
    1654             : {
    1655             :   GEN X, Xp, Xq, g, V;
    1656             :   long n;
    1657          83 :   GEN R = FpXQX_easyroots(f, T, p);
    1658          83 :   if (R) return R;
    1659          83 :   X  = pol_x(varn(f));
    1660          83 :   Xp = FpXQXQ_pow(X, p, f, T, p);
    1661          83 :   Xq = FpXQXQ_Frobenius(xp, Xp, f, T, p);
    1662          83 :   g = FpXQX_gcd(FpXX_sub(Xq, X, p), f, T, p);
    1663          83 :   n = degpol(g);
    1664          83 :   if (n==0) return cgetg(1, t_COL);
    1665          83 :   g = FpXQX_normalize(g, T, p);
    1666          83 :   V = cgetg(n+1,t_COL);
    1667          83 :   FpXQX_roots_edf(g, xp, Xp, T, p, V, 1);
    1668          83 :   return V;
    1669             : }
    1670             : 
    1671             : /* do not handle small p */
    1672             : static GEN
    1673       24454 : FpXQX_roots_i(GEN S, GEN T, GEN p)
    1674             : {
    1675             :   GEN F, M;
    1676       24454 :   if (lgefint(p)==3)
    1677             :   {
    1678       23854 :     ulong pp = p[2];
    1679       23854 :     if (pp == 2)
    1680             :     {
    1681        3843 :       GEN V = F2xqX_roots_i(ZXX_to_F2xX(S,get_FpX_var(T)), ZX_to_F2x(get_FpX_mod(T)));
    1682        3843 :       return F2xC_to_ZXC(V);
    1683             :     }
    1684             :     else
    1685             :     {
    1686       20011 :       GEN V = FlxqX_roots_i(ZXX_to_FlxX(S,pp,get_FpX_var(T)),
    1687             :                             ZXT_to_FlxT(T,pp), pp);
    1688       20011 :       return FlxC_to_ZXC(V);
    1689             :     }
    1690             :   }
    1691         600 :   S = FpXQX_red(S, T, p);
    1692         600 :   if (!signe(S)) pari_err_ROOTS0("FpXQX_roots");
    1693         600 :   if (degpol(S)==0) return cgetg(1, t_COL);
    1694         600 :   S = FpXQX_normalize(S, T, p);
    1695         600 :   M = FpXQX_easyroots(S, T, p);
    1696         600 :   if (!M)
    1697             :   {
    1698          83 :     GEN xp = FpX_Frobenius(T, p);
    1699          83 :     GEN V = FpXQX_factor_Yun(S, T, p);
    1700          83 :     long i, j, l = lg(V);
    1701          83 :     F = cgetg(l, t_VEC);
    1702         166 :     for (i=1, j=1; i < l; i++)
    1703          83 :       if (degpol(gel(V,i)))
    1704          83 :         gel(F, j++) = FpXQX_roots_ddf(gel(V,i), xp, T, p);
    1705          83 :     setlg(F,j); M = shallowconcat1(F);
    1706             :   }
    1707         600 :   gen_sort_inplace(M, (void*) &cmp_RgX, &cmp_nodata, NULL);
    1708         600 :   return M;
    1709             : }
    1710             : 
    1711             : GEN
    1712       71295 : F2xqX_roots(GEN x, GEN T)
    1713             : {
    1714       71295 :   pari_sp av = avma;
    1715       71295 :   return gerepilecopy(av, F2xqX_roots_i(x, T));
    1716             : }
    1717             : 
    1718             : GEN
    1719       46483 : FlxqX_roots(GEN x, GEN T, ulong p)
    1720             : {
    1721       46483 :   pari_sp av = avma;
    1722       46483 :   if (p==2)
    1723             :   {
    1724           0 :     GEN V = F2xqX_roots_i(FlxX_to_F2xX(x), Flx_to_F2x(get_Flx_mod(T)));
    1725           0 :     return gerepileupto(av, F2xC_to_FlxC(V));
    1726             :   }
    1727       46483 :   return gerepilecopy(av, FlxqX_roots_i(x, T, p));
    1728             : }
    1729             : 
    1730             : GEN
    1731       24454 : FpXQX_roots(GEN x, GEN T, GEN p)
    1732             : {
    1733       24454 :   pari_sp av = avma;
    1734       24454 :   return gerepilecopy(av, FpXQX_roots_i(x, T, p));
    1735             : }
    1736             : 
    1737             : static GEN
    1738         553 : FE_concat(GEN F, GEN E, long l)
    1739             : {
    1740         553 :   setlg(E,l); E = shallowconcat1(E);
    1741         553 :   setlg(F,l); F = shallowconcat1(F); return mkvec2(F,E);
    1742             : }
    1743             : 
    1744             : static GEN
    1745         532 : F2xqX_factor_2(GEN f, GEN T)
    1746             : {
    1747         532 :   long v = varn(f), vT = get_F2x_var(T);
    1748         532 :   GEN r = F2xqX_quad_roots(f, T);
    1749         532 :   switch(lg(r)-1)
    1750             :   {
    1751          14 :   case 0:
    1752          14 :     return mkvec2(mkcolcopy(f), mkvecsmall(1));
    1753         511 :   case 1:
    1754         511 :     return mkvec2(mkcol(deg1pol_shallow(pol1_F2x(vT), gel(r,1), v)), mkvecsmall(2));
    1755           7 :   default: /* 2 */
    1756             :     {
    1757           7 :       GEN f1 = deg1pol_shallow(pol1_F2x(vT), gel(r,1), v);
    1758           7 :       GEN f2 = deg1pol_shallow(pol1_F2x(vT), gel(r,2), v);
    1759           7 :       GEN t = mkcol2(f1, f2), E = mkvecsmall2(1, 1);
    1760           7 :       sort_factor_pol(mkvec2(t, E), cmp_Flx);
    1761           7 :       return mkvec2(t, E);
    1762             :     }
    1763             :   }
    1764             : }
    1765             : 
    1766             : static GEN
    1767         825 : FlxqX_factor_2(GEN f, GEN T, ulong p, ulong pi)
    1768             : {
    1769         825 :   long v = varn(f), vT = get_Flx_var(T);
    1770         825 :   GEN r = FlxqX_quad_roots(f, T, p, pi);
    1771         825 :   switch(lg(r)-1)
    1772             :   {
    1773          56 :   case 0:
    1774          56 :     return mkvec2(mkcolcopy(f), mkvecsmall(1));
    1775          35 :   case 1:
    1776          70 :     return mkvec2(mkcol(deg1pol_shallow(pol1_Flx(vT),
    1777          35 :                         Flx_neg(gel(r,1), p), v)), mkvecsmall(2));
    1778         734 :   default: /* 2 */
    1779             :     {
    1780         734 :       GEN f1 = deg1pol_shallow(pol1_Flx(vT), Flx_neg(gel(r,1), p), v);
    1781         734 :       GEN f2 = deg1pol_shallow(pol1_Flx(vT), Flx_neg(gel(r,2), p), v);
    1782         734 :       GEN t = mkcol2(f1, f2), E = mkvecsmall2(1, 1);
    1783         734 :       sort_factor_pol(mkvec2(t, E), cmp_Flx);
    1784         734 :       return mkvec2(t, E);
    1785             :     }
    1786             :   }
    1787             : }
    1788             : 
    1789             : static GEN
    1790          43 : FpXQX_factor_2(GEN f, GEN T, GEN p)
    1791             : {
    1792          43 :   long v = varn(f);
    1793          43 :   GEN r = FpXQX_quad_roots(f, T, p);
    1794          43 :   switch(lg(r)-1)
    1795             :   {
    1796          14 :   case 0:
    1797          14 :     return mkvec2(mkcolcopy(f), mkvecsmall(1));
    1798           0 :   case 1:
    1799           0 :     return mkvec2(mkcol(deg1pol_shallow(gen_1, Fq_neg(gel(r,1), T, p), v)),
    1800             :         mkvecsmall(2));
    1801          29 :   default: /* 2 */
    1802             :     {
    1803          29 :       GEN f1 = deg1pol_shallow(gen_1, Fq_neg(gel(r,1), T, p), v);
    1804          29 :       GEN f2 = deg1pol_shallow(gen_1, Fq_neg(gel(r,2), T, p), v);
    1805          29 :       GEN t = mkcol2(f1, f2), E = mkvecsmall2(1, 1);
    1806          29 :       sort_factor_pol(mkvec2(t, E), cmp_RgX);
    1807          29 :       return mkvec2(t, E);
    1808             :     }
    1809             :   }
    1810             : }
    1811             : 
    1812             : static GEN
    1813         294 : F2xqX_ddf_Shoup(GEN S, GEN Xq, GEN T)
    1814             : {
    1815         294 :   pari_sp av = avma;
    1816             :   GEN b, g, h, F, f, Sr, xq;
    1817             :   long i, j, n, v, vT, dT, bo, ro;
    1818             :   long B, l, m;
    1819             :   pari_timer ti;
    1820         294 :   n = get_F2xqX_degree(S); v = get_F2xqX_var(S);
    1821         294 :   vT = get_F2x_var(T); dT = get_F2x_degree(T);
    1822         294 :   if (n == 0) return cgetg(1, t_VEC);
    1823         294 :   if (n == 1) return mkvec(get_F2xqX_mod(S));
    1824         140 :   B = n/2;
    1825         140 :   l = usqrt(B);
    1826         140 :   m = (B+l-1)/l;
    1827         140 :   S = F2xqX_get_red(S, T);
    1828         140 :   b = cgetg(l+2, t_VEC);
    1829         140 :   gel(b, 1) = polx_F2xX(v, vT);
    1830         140 :   gel(b, 2) = Xq;
    1831         140 :   bo = brent_kung_optpow(n, l-1, 1);
    1832         140 :   ro = l<=1 ? 0: (bo-1)/(l-1) + ((n-1)/bo);
    1833         140 :   if (DEBUGLEVEL>=7) timer_start(&ti);
    1834         140 :   if (dT <= ro)
    1835           0 :     for (i = 3; i <= l+1; i++)
    1836           0 :       gel(b, i) = F2xqXQ_pow(gel(b, i-1), int2n(dT), S, T);
    1837             :   else
    1838             :   {
    1839         140 :     xq = F2xqXQ_powers(gel(b, 2), bo, S, T);
    1840         140 :     if (DEBUGLEVEL>=7) timer_printf(&ti,"F2xqX_ddf_Shoup: xq baby");
    1841         140 :     for (i = 3; i <= l+1; i++)
    1842           0 :       gel(b, i) = F2xqX_F2xqXQV_eval(gel(b, i-1), xq, S, T);
    1843             :   }
    1844         140 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"F2xqX_ddf_Shoup: baby");
    1845         140 :   xq = F2xqXQ_powers(gel(b, l+1), brent_kung_optpow(n, m-1, 1), S, T);
    1846         140 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"F2xqX_ddf_Shoup: xq giant");
    1847         140 :   g = cgetg(m+1, t_VEC);
    1848         140 :   gel(g, 1) = gel(xq, 2);
    1849         168 :   for(i = 2; i <= m; i++)
    1850          28 :     gel(g, i) = F2xqX_F2xqXQV_eval(gel(g, i-1), xq, S, T);
    1851         140 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"F2xqX_ddf_Shoup: giant");
    1852         140 :   h = cgetg(m+1, t_VEC);
    1853         308 :   for (j = 1; j <= m; j++)
    1854             :   {
    1855         168 :     pari_sp av = avma;
    1856         168 :     GEN gj = gel(g, j);
    1857         168 :     GEN e = F2xX_add(gj, gel(b, 1));
    1858         168 :     for (i = 2; i <= l; i++)
    1859           0 :       e = F2xqXQ_mul(e, F2xX_add(gj, gel(b, i)), S, T);
    1860         168 :     gel(h, j) = gerepileupto(av, e);
    1861             :   }
    1862         140 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"F2xqX_ddf_Shoup: diff");
    1863         140 :   Sr = get_F2xqX_mod(S);
    1864         140 :   F = cgetg(m+1, t_VEC);
    1865         308 :   for (j = 1; j <= m; j++)
    1866             :   {
    1867         168 :     GEN u = F2xqX_gcd(Sr, gel(h,j), T);
    1868         168 :     if (degpol(u))
    1869             :     {
    1870         112 :       u = F2xqX_normalize(u, T);
    1871         112 :       Sr = F2xqX_div(Sr, u, T);
    1872             :     }
    1873         168 :     gel(F,j) = u;
    1874             :   }
    1875         140 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"F2xqX_ddf_Shoup: F");
    1876         140 :   f = const_vec(n, pol1_F2xX(v, vT));
    1877         308 :   for (j = 1; j <= m; j++)
    1878             :   {
    1879         168 :     GEN e = gel(F, j);
    1880         168 :     for (i=l-1; i >= 0; i--)
    1881             :     {
    1882         168 :       GEN u = F2xqX_gcd(e, F2xX_add(gel(g, j), gel(b, i+1)), T);
    1883         168 :       if (degpol(u))
    1884             :       {
    1885         112 :         gel(f, l*j-i) = u = F2xqX_normalize(u, T);
    1886         112 :         e = F2xqX_div(e, u, T);
    1887             :       }
    1888         168 :       if (!degpol(e)) break;
    1889             :     }
    1890             :   }
    1891         140 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"F2xqX_ddf_Shoup: f");
    1892         140 :   if (degpol(Sr)) gel(f, degpol(Sr)) = Sr;
    1893         140 :   return gerepilecopy(av, f);
    1894             : }
    1895             : 
    1896             : static GEN
    1897          91 : F2xqX_ddf_i(GEN f, GEN T, GEN X, GEN xp)
    1898             : {
    1899             :   GEN Xp, Xq;
    1900          91 :   if (!get_F2xqX_degree(f)) return cgetg(1,t_VEC);
    1901          42 :   f = F2xqX_get_red(f, T);
    1902          42 :   Xp = F2xqXQ_sqr(X, f, T);
    1903          42 :   Xq = F2xqXQ_Frobenius(xp, Xp, f, T);
    1904          42 :   return F2xqX_ddf_Shoup(f, Xq, T);
    1905             : }
    1906             : static void
    1907          42 : F2xqX_ddf_init(GEN *S, GEN *T, GEN *xp, GEN *X)
    1908             : {
    1909          42 :   *T = F2x_get_red(*T);
    1910          42 :   *S = F2xqX_normalize(get_F2xqX_mod(*S), *T);
    1911          42 :   *xp = F2x_Frobenius(*T);
    1912          42 :   *X  = polx_F2xX(get_F2xqX_var(*S), get_F2x_var(*T));
    1913          42 : }
    1914             : GEN
    1915          42 : F2xqX_degfact(GEN S, GEN T)
    1916             : {
    1917             :   GEN xp, X, V;
    1918             :   long i, l;
    1919          42 :   F2xqX_ddf_init(&S,&T,&xp,&X);
    1920          42 :   V = F2xqX_factor_squarefree(S, T); l = lg(V);
    1921         133 :   for (i=1; i < l; i++) gel(V,i) = F2xqX_ddf_i(gel(V,i), T, X, xp);
    1922          42 :   return vddf_to_simplefact(V, degpol(S));
    1923             : }
    1924             : GEN
    1925           0 : F2xqX_ddf(GEN S, GEN T)
    1926             : {
    1927             :   GEN xp, X;
    1928           0 :   F2xqX_ddf_init(&S,&T,&xp,&X);
    1929           0 :   return ddf_to_ddf2( F2xqX_ddf_i(S, T, X, xp) );
    1930             : }
    1931             : 
    1932             : static void
    1933         168 : F2xqX_edf_simple(GEN Sp, GEN xp, GEN Xp, GEN Sq, long d, GEN T, GEN V, long idx)
    1934             : {
    1935         168 :   long v = varn(Sp), n = degpol(Sp), r = n/d;
    1936             :   GEN S, f, ff;
    1937         168 :   long dT = get_F2x_degree(T);
    1938         168 :   if (r==1) { gel(V, idx) = Sp; return; }
    1939          63 :   S = F2xqX_get_red(Sp, T);
    1940          63 :   Xp = F2xqX_rem(Xp, S, T);
    1941          63 :   Sq = F2xqXQV_red(Sq, S, T);
    1942             :   while (1)
    1943          56 :   {
    1944         119 :     pari_sp btop = avma;
    1945             :     long l;
    1946         119 :     GEN w0 = random_F2xqX(n, v, T), g = w0;
    1947         140 :     for (l=1; l<d; l++) /* sum_{0<i<d} w^(q^i), result in (F_q)^r */
    1948          21 :       g = F2xX_add(w0, F2xqX_F2xqXQV_eval(g, Sq, S, T));
    1949         119 :     w0 = g;
    1950         700 :     for (l=1; l<dT; l++) /* sum_{0<i<k} w^(2^i), result in (F_2)^r */
    1951         581 :       g = F2xX_add(w0, F2xqXQ_sqr(g,S,T));
    1952         119 :     f = F2xqX_gcd(g, Sp, T);
    1953         119 :     if (degpol(f) > 0 && degpol(f) < n) break;
    1954          56 :     set_avma(btop);
    1955             :   }
    1956          63 :   f = F2xqX_normalize(f, T);
    1957          63 :   ff = F2xqX_div(Sp, f , T);
    1958          63 :   F2xqX_edf_simple(f, xp, Xp, Sq, d, T, V, idx);
    1959          63 :   F2xqX_edf_simple(ff, xp, Xp, Sq, d, T, V, idx+degpol(f)/d);
    1960             : }
    1961             : 
    1962             : static GEN
    1963         252 : F2xqX_factor_Shoup(GEN S, GEN xp, GEN T)
    1964             : {
    1965         252 :   long i, n, s = 0;
    1966             :   GEN X, Xp, Xq, Sq, D, V;
    1967         252 :   long vT = get_F2x_var(T);
    1968             :   pari_timer ti;
    1969         252 :   n = get_F2xqX_degree(S);
    1970         252 :   S = F2xqX_get_red(S, T);
    1971         252 :   if (DEBUGLEVEL>=6) timer_start(&ti);
    1972         252 :   X  = polx_F2xX(get_F2xqX_var(S), vT);
    1973         252 :   Xp = F2xqXQ_sqr(X, S, T);
    1974         252 :   Xq = F2xqXQ_Frobenius(xp, Xp, S, T);
    1975         252 :   Sq = F2xqXQ_powers(Xq, n-1, S, T);
    1976         252 :   if (DEBUGLEVEL>=6) timer_printf(&ti,"F2xqX_Frobenius");
    1977         252 :   D = F2xqX_ddf_Shoup(S, Xq, T);
    1978         252 :   if (DEBUGLEVEL>=6) timer_printf(&ti,"F2xqX_ddf_Shoup");
    1979         252 :   s = ddf_to_nbfact(D);
    1980         252 :   V = cgetg(s+1, t_COL);
    1981         763 :   for (i = 1, s = 1; i <= n; i++)
    1982             :   {
    1983         511 :     GEN Di = gel(D,i);
    1984         511 :     long ni = degpol(Di), ri = ni/i;
    1985         511 :     if (ni == 0) continue;
    1986         308 :     Di = F2xqX_normalize(Di, T);
    1987         308 :     if (ni == i) { gel(V, s++) = Di; continue; }
    1988          42 :     F2xqX_edf_simple(Di, xp, Xp, Sq, i, T, V, s);
    1989          42 :     if (DEBUGLEVEL>=6) timer_printf(&ti,"F2xqX_edf(%ld)",i);
    1990          42 :     s += ri;
    1991             :   }
    1992         252 :   return V;
    1993             : }
    1994             : 
    1995             : static GEN
    1996        1239 : F2xqX_factor_Cantor(GEN f, GEN T)
    1997             : {
    1998             :   GEN xp, E, F, V;
    1999             :   long i, j, l;
    2000        1239 :   T = F2x_get_red(T);
    2001        1239 :   f = F2xqX_normalize(f, T);
    2002        1239 :   switch(degpol(f))
    2003             :   {
    2004          14 :     case -1: retmkmat2(mkcol(f), mkvecsmall(1));
    2005          14 :     case 0: return trivial_fact();
    2006          21 :     case 1: retmkmat2(mkcol(f), mkvecsmall(1));
    2007         532 :     case 2: return F2xqX_factor_2(f, T);
    2008             :   }
    2009         658 :   if (F2xY_degreex(f) <= 0) return F2x_factorff_i(F2xX_to_F2x(f), T);
    2010         238 :   xp = F2x_Frobenius(T);
    2011         238 :   V = F2xqX_factor_squarefree(f, T);
    2012         238 :   l = lg(V);
    2013         238 :   F = cgetg(l, t_VEC);
    2014         238 :   E = cgetg(l, t_VEC);
    2015         784 :   for (i=1, j=1; i < l; i++)
    2016         546 :     if (degpol(gel(V,i)))
    2017             :     {
    2018         252 :       GEN Fj = F2xqX_factor_Shoup(gel(V,i), xp, T);
    2019         252 :       gel(F, j) = Fj;
    2020         252 :       gel(E, j) = const_vecsmall(lg(Fj)-1, i);
    2021         252 :       j++;
    2022             :     }
    2023         238 :   return sort_factor_pol(FE_concat(F,E,j), cmp_Flx);
    2024             : }
    2025             : 
    2026             : static GEN
    2027           0 : FlxqX_Berlekamp_i(GEN f, GEN T, ulong p)
    2028             : {
    2029           0 :   long lfact, d = degpol(f), j, k, lV;
    2030             :   GEN E, t, V, xp;
    2031             :   ulong pi;
    2032           0 :   switch(d)
    2033             :   {
    2034           0 :     case -1: retmkmat2(mkcolcopy(f), mkvecsmall(1));
    2035           0 :     case 0: return trivial_fact();
    2036             :   }
    2037           0 :   pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
    2038           0 :   T = Flx_get_red_pre(T, p, pi);
    2039           0 :   f = FlxqX_normalize_pre(f, T, p, pi);
    2040           0 :   if (FlxY_degreex(f) <= 0) return Flx_factorff_i(FlxX_to_Flx(f), T, p);
    2041           0 :   if (degpol(f)==2) return FlxqX_factor_2(f, T, p, pi);
    2042           0 :   xp = Flx_Frobenius_pre(T, p, pi);
    2043           0 :   V = FlxqX_factor_squarefree_i(f, xp, T, p, pi); lV = lg(V);
    2044             : 
    2045             :   /* to hold factors and exponents */
    2046           0 :   t = cgetg(d+1,t_VEC);
    2047           0 :   E = cgetg(d+1, t_VECSMALL);
    2048           0 :   lfact = 1;
    2049           0 :   for (k=1; k<lV ; k++)
    2050             :   {
    2051           0 :     if (degpol(gel(V,k))==0) continue;
    2052           0 :     gel(t,lfact) = FlxqX_normalize_pre(gel(V, k), T,p, pi);
    2053           0 :     d = FlxqX_split_Berlekamp(&gel(t,lfact), xp, T, p, pi);
    2054           0 :     for (j = 0; j < d; j++) E[lfact+j] = k;
    2055           0 :     lfact += d;
    2056             :   }
    2057           0 :   setlg(t, lfact);
    2058           0 :   setlg(E, lfact);
    2059           0 :   return sort_factor_pol(mkvec2(t, E), cmp_Flx);
    2060             : }
    2061             : 
    2062             : static GEN
    2063           0 : FpXQX_Berlekamp_i(GEN f, GEN T, GEN p)
    2064             : {
    2065           0 :   long lfact, d = degpol(f), j, k, lV;
    2066             :   GEN E, t, V;
    2067           0 :   switch(d)
    2068             :   {
    2069           0 :     case -1: retmkmat2(mkcolcopy(f), mkvecsmall(1));
    2070           0 :     case 0: return trivial_fact();
    2071             :   }
    2072           0 :   T = FpX_get_red(T, p);
    2073           0 :   f = FpXQX_normalize(f, T, p);
    2074           0 :   if (isabsolutepol(f)) return FpX_factorff_i(simplify_shallow(f), T, p);
    2075           0 :   if (degpol(f)==2) return FpXQX_factor_2(f, T, p);
    2076           0 :   V = FpXQX_factor_Yun(f, T, p); lV = lg(V);
    2077             : 
    2078             :   /* to hold factors and exponents */
    2079           0 :   t = cgetg(d+1,t_VEC);
    2080           0 :   E = cgetg(d+1, t_VECSMALL);
    2081           0 :   lfact = 1;
    2082           0 :   for (k=1; k<lV ; k++)
    2083             :   {
    2084           0 :     if (degpol(gel(V,k))==0) continue;
    2085           0 :     gel(t,lfact) = FpXQX_normalize(gel(V, k), T,p);
    2086           0 :     d = FpXQX_split_Berlekamp(&gel(t,lfact), T, p);
    2087           0 :     for (j = 0; j < d; j++) E[lfact+j] = k;
    2088           0 :     lfact += d;
    2089             :   }
    2090           0 :   setlg(t, lfact);
    2091           0 :   setlg(E, lfact);
    2092           0 :   return sort_factor_pol(mkvec2(t, E), cmp_RgX);
    2093             : }
    2094             : 
    2095             : long
    2096         268 : FlxqX_ddf_degree(GEN S, GEN XP, GEN T, ulong p)
    2097             : {
    2098         268 :   pari_sp av = avma;
    2099             :   GEN X, b, g, xq, q;
    2100             :   long i, j, n, v, B, l, m, bo, ro;
    2101             :   ulong pi;
    2102             :   pari_timer ti;
    2103             :   hashtable h;
    2104             : 
    2105         268 :   n = get_FlxqX_degree(S); v = get_FlxqX_var(S);
    2106         268 :   X = polx_FlxX(v,get_Flx_var(T));
    2107         268 :   if (gequal(X,XP)) return 1;
    2108         268 :   pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
    2109         268 :   B = n/2;
    2110         268 :   l = usqrt(B);
    2111         268 :   m = (B+l-1)/l;
    2112         268 :   T = Flx_get_red_pre(T, p, pi);
    2113         268 :   S = FlxqX_get_red_pre(S, T, p, pi);
    2114         268 :   hash_init_GEN(&h, l+2, gequal, 1);
    2115         268 :   hash_insert_long(&h, X,  0);
    2116         268 :   hash_insert_long(&h, XP, 1);
    2117         268 :   bo = brent_kung_optpow(n, l-1, 1);
    2118         268 :   ro = l<=1 ? 0: (bo-1)/(l-1) + ((n-1)/bo);
    2119         268 :   q = powuu(p, get_Flx_degree(T));
    2120         268 :   if (DEBUGLEVEL>=7) timer_start(&ti);
    2121         268 :   b = XP;
    2122         268 :   if (expi(q) > ro)
    2123             :   {
    2124         268 :     xq = FlxqXQ_powers_pre(b, brent_kung_optpow(n, l-1, 1),  S, T, p, pi);
    2125         268 :     if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_degree: xq baby");
    2126           0 :   } else xq = NULL;
    2127         635 :   for (i = 3; i <= l+1; i++)
    2128             :   {
    2129         401 :     b = xq ? FlxqX_FlxqXQV_eval_pre(b, xq, S, T, p, pi)
    2130         401 :            : FlxqXQ_pow_pre(b, q, S, T, p, pi);
    2131         401 :     if (gequal(b,X)) return gc_long(av,i-1);
    2132         367 :     hash_insert_long(&h, b, i-1);
    2133             :   }
    2134         234 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_degree: baby");
    2135         234 :   g = b;
    2136         234 :   xq = FlxqXQ_powers_pre(g, brent_kung_optpow(n, m, 1),  S, T, p, pi);
    2137         234 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_degree: xq giant");
    2138         651 :   for(i = 2; i <= m+1; i++)
    2139             :   {
    2140         562 :     g = FlxqX_FlxqXQV_eval_pre(g, xq, S, T, p, pi);
    2141         562 :     if (hash_haskey_long(&h, g, &j)) return gc_long(av, l*i-j);
    2142             :   }
    2143          89 :   return gc_long(av,n);
    2144             : }
    2145             : 
    2146             : static GEN
    2147         531 : FlxqX_ddf_Shoup(GEN S, GEN Xq, GEN T, ulong p, ulong pi)
    2148             : {
    2149         531 :   pari_sp av = avma;
    2150             :   GEN b, g, h, F, f, Sr, xq, q;
    2151             :   long i, j, n, v, vT, bo, ro;
    2152             :   long B, l, m;
    2153             :   pari_timer ti;
    2154         531 :   n = get_FlxqX_degree(S); v = get_FlxqX_var(S);
    2155         531 :   vT = get_Flx_var(T);
    2156         531 :   if (n == 0) return cgetg(1, t_VEC);
    2157         531 :   if (n == 1) return mkvec(get_FlxqX_mod(S));
    2158         377 :   B = n/2;
    2159         377 :   l = usqrt(B);
    2160         377 :   m = (B+l-1)/l;
    2161         377 :   S = FlxqX_get_red_pre(S, T, p, pi);
    2162         377 :   b = cgetg(l+2, t_VEC);
    2163         377 :   gel(b, 1) = polx_FlxX(v, vT);
    2164         377 :   gel(b, 2) = Xq;
    2165         377 :   bo = brent_kung_optpow(n, l-1, 1);
    2166         377 :   ro = l<=1 ? 0: (bo-1)/(l-1) + ((n-1)/bo);
    2167         377 :   q = powuu(p, get_Flx_degree(T));
    2168         377 :   if (DEBUGLEVEL>=7) timer_start(&ti);
    2169         377 :   if (expi(q) <= ro)
    2170          21 :     for (i = 3; i <= l+1; i++)
    2171          14 :       gel(b, i) = FlxqXQ_pow_pre(gel(b, i-1), q, S, T, p, pi);
    2172             :   else
    2173             :   {
    2174         370 :     xq = FlxqXQ_powers_pre(gel(b, 2), bo, S, T, p, pi);
    2175         370 :     if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_Shoup: xq baby");
    2176         370 :     for (i = 3; i <= l+1; i++)
    2177           0 :       gel(b, i) = FlxqX_FlxqXQV_eval_pre(gel(b, i-1), xq, S, T, p, pi);
    2178             :   }
    2179         377 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_Shoup: baby");
    2180         377 :   xq = FlxqXQ_powers_pre(gel(b, l+1), brent_kung_optpow(n, m-1, 1), S, T,p,pi);
    2181         377 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_Shoup: xq giant");
    2182         377 :   g = cgetg(m+1, t_VEC);
    2183         377 :   gel(g, 1) = gel(xq, 2);
    2184         746 :   for(i = 2; i <= m; i++)
    2185         369 :     gel(g, i) = FlxqX_FlxqXQV_eval_pre(gel(g, i-1), xq, S, T, p, pi);
    2186         377 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_Shoup: giant");
    2187         377 :   h = cgetg(m+1, t_VEC);
    2188        1123 :   for (j = 1; j <= m; j++)
    2189             :   {
    2190         746 :     pari_sp av = avma;
    2191         746 :     GEN gj = gel(g, j);
    2192         746 :     GEN e = FlxX_sub(gj, gel(b, 1), p);
    2193         816 :     for (i = 2; i <= l; i++)
    2194          70 :       e = FlxqXQ_mul_pre(e, FlxX_sub(gj, gel(b, i), p), S, T, p, pi);
    2195         746 :     gel(h, j) = gerepileupto(av, e);
    2196             :   }
    2197         377 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_Shoup: diff");
    2198         377 :   Sr = get_FlxqX_mod(S);
    2199         377 :   F = cgetg(m+1, t_VEC);
    2200        1123 :   for (j = 1; j <= m; j++)
    2201             :   {
    2202         746 :     GEN u = FlxqX_gcd_pre(Sr, gel(h, j), T, p, pi);
    2203         746 :     if (degpol(u))
    2204             :     {
    2205         300 :       u = FlxqX_normalize_pre(u, T, p, pi);
    2206         300 :       Sr = FlxqX_div_pre(Sr, u, T, p, pi);
    2207             :     }
    2208         746 :     gel(F,j) = u;
    2209             :   }
    2210         377 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_Shoup: F");
    2211         377 :   f = const_vec(n, pol1_FlxX(v, vT));
    2212        1123 :   for (j = 1; j <= m; j++)
    2213             :   {
    2214         746 :     GEN e = gel(F, j);
    2215         746 :     for (i=l-1; i >= 0; i--)
    2216             :     {
    2217         746 :       GEN u = FlxqX_gcd_pre(e, FlxX_sub(gel(g, j), gel(b, i+1), p), T, p, pi);
    2218         746 :       if (degpol(u))
    2219             :       {
    2220         300 :         gel(f, l*j-i) = u = FlxqX_normalize_pre(u, T, p, pi);
    2221         300 :         e = FlxqX_div_pre(e, u, T, p, pi);
    2222             :       }
    2223         746 :       if (!degpol(e)) break;
    2224             :     }
    2225             :   }
    2226         377 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_ddf_Shoup: f");
    2227         377 :   if (degpol(Sr)) gel(f, degpol(Sr)) = Sr;
    2228         377 :   return gerepilecopy(av, f);
    2229             : }
    2230             : 
    2231             : static GEN
    2232          42 : FlxqX_ddf_i(GEN f, GEN T, ulong p, ulong pi)
    2233             : {
    2234             :   GEN Xq;
    2235          42 :   if (!get_FlxqX_degree(f)) return cgetg(1, t_VEC);
    2236          42 :   f = FlxqX_get_red_pre(f, T, p, pi);
    2237          42 :   Xq = FlxqX_Frobenius_pre(f, T, p, pi);
    2238          42 :   return FlxqX_ddf_Shoup(f, Xq, T, p, pi);
    2239             : }
    2240             : GEN
    2241           0 : FlxqX_ddf(GEN S, GEN T, ulong p)
    2242             : {
    2243           0 :   ulong pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
    2244           0 :   T = Flx_get_red_pre(T, p, pi);
    2245           0 :   S = FlxqX_normalize_pre(get_FlxqX_mod(S), T, p, pi);
    2246           0 :   return ddf_to_ddf2( FlxqX_ddf_i(S, T, p, pi) );
    2247             : }
    2248             : GEN
    2249          42 : FlxqX_degfact(GEN S, GEN T, ulong p)
    2250             : {
    2251          42 :   ulong pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
    2252             :   GEN V;
    2253             :   long i, l;
    2254          42 :   T = Flx_get_red_pre(T, p, pi);
    2255          42 :   S = FlxqX_normalize_pre(get_FlxqX_mod(S), T, p, pi);
    2256          42 :   V = FlxqX_factor_squarefree_pre(S, T, p, pi); l = lg(V);
    2257          84 :   for (i=1; i < l; i++) gel(V,i) = FlxqX_ddf_i(gel(V,i), T, p, pi);
    2258          42 :   return vddf_to_simplefact(V, degpol(S));
    2259             : }
    2260             : 
    2261             : static void
    2262         362 : FlxqX_edf_rec(GEN S, GEN xp, GEN Xp, GEN hp, GEN t, long d, GEN T, ulong p,
    2263             :   ulong pi, GEN V, long idx)
    2264             : {
    2265         362 :   GEN Sp = get_FlxqX_mod(S);
    2266             :   GEN u1, u2, f1, f2;
    2267             :   GEN h;
    2268         362 :   h = FlxqX_get_red_pre(hp, T, p, pi);
    2269         362 :   t = FlxqX_rem_pre(t, S, T, p, pi);
    2270         362 :   Xp = FlxqX_rem_pre(Xp, h, T, p, pi);
    2271         362 :   u1 = FlxqX_roots_split(hp, xp, Xp, h, T, p, pi);
    2272         362 :   f1 = FlxqX_gcd_pre(FlxqX_FlxqXQ_eval_pre(u1, t, S, T, p, pi), Sp, T, p, pi);
    2273         362 :   f1 = FlxqX_normalize_pre(f1, T, p, pi);
    2274         362 :   u2 = FlxqX_div_pre(hp, u1, T, p, pi);
    2275         362 :   f2 = FlxqX_div_pre(Sp, f1, T, p, pi);
    2276         362 :   if (degpol(u1)==1)
    2277         262 :     gel(V, idx) = f1;
    2278             :   else
    2279         100 :     FlxqX_edf_rec(FlxqX_get_red_pre(f1, T, p, pi),
    2280             :                   xp, Xp, u1, t, d, T, p, pi, V, idx);
    2281         362 :   idx += degpol(f1)/d;
    2282         362 :   if (degpol(u2)==1)
    2283         246 :     gel(V, idx) = f2;
    2284             :   else
    2285         116 :     FlxqX_edf_rec(FlxqX_get_red_pre(f2, T, p, pi),
    2286             :                   xp, Xp, u2, t, d, T, p, pi, V, idx);
    2287         362 : }
    2288             : 
    2289             : static void
    2290         146 : FlxqX_edf(GEN Sp, GEN xp, GEN Xp, GEN Xq, long d, GEN T, ulong p, ulong pi,
    2291             :   GEN V, long idx)
    2292             : {
    2293         146 :   long n = degpol(Sp), r = n/d, vS = varn(Sp), vT = get_Flx_var(T);
    2294             :   GEN S, h, t;
    2295             :   pari_timer ti;
    2296         146 :   if (r==1) { gel(V, idx) = Sp; return; }
    2297         146 :   S = FlxqX_get_red_pre(Sp, T, p, pi);
    2298         146 :   Xp = FlxqX_rem_pre(Xp, S, T, p, pi);
    2299         146 :   Xq = FlxqX_rem_pre(Xq, S, T, p, pi);
    2300         146 :   if (DEBUGLEVEL>=7) timer_start(&ti);
    2301             :   do
    2302             :   {
    2303         162 :     GEN g = random_FlxqX(n, vS, T, p);
    2304         162 :     t = gel(FlxqXQ_auttrace_pre(mkvec2(Xq, g), d, S, T, p, pi), 2);
    2305         162 :     if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_edf: FlxqXQ_auttrace");
    2306         162 :     h = FlxqXQ_minpoly_pre(t, S, T, p, pi);
    2307         162 :     if (DEBUGLEVEL>=7) timer_printf(&ti,"FlxqX_edf: FlxqXQ_minpoly");
    2308         162 :   } while (degpol(h) != r);
    2309         146 :   Xp = FlxqXQ_powu_pre(polx_FlxX(vS, vT), p, h, T, p, pi);
    2310         146 :   FlxqX_edf_rec(S, xp, Xp, h, t, d, T, p, pi, V, idx);
    2311             : }
    2312             : 
    2313             : static void
    2314         357 : FlxqX_edf_simple(GEN Sp, GEN xp, GEN Xp, GEN Xq, long d, GEN T, ulong p,
    2315             :   ulong pi, GEN V, long idx)
    2316             : {
    2317         357 :   long v = varn(Sp), n = degpol(Sp), r = n/d;
    2318             :   GEN S, f, ff;
    2319         357 :   long vT = get_Flx_var(T), dT = get_Flx_degree(T);
    2320         357 :   if (r==1) { gel(V, idx) = Sp; return; }
    2321         175 :   S = FlxqX_get_red_pre(Sp, T, p, pi);
    2322         175 :   Xp = FlxqX_rem_pre(Xp, S, T, p, pi);
    2323         175 :   Xq = FlxqX_rem_pre(Xq, S, T, p, pi);
    2324             :   while (1)
    2325           0 :   {
    2326         175 :     pari_sp btop = avma;
    2327             :     long i;
    2328         175 :     GEN g = random_FlxqX(n, v, T, p);
    2329         175 :     GEN t = gel(FlxqXQ_auttrace_pre(mkvec2(Xq, g), d, S, T, p, pi), 2);
    2330         175 :     if (lgpol(t) == 0) continue;
    2331         245 :     for(i=1; i<=10; i++)
    2332             :     {
    2333         245 :       pari_sp btop2 = avma;
    2334         245 :       GEN r = random_Flx(dT, vT, p);
    2335         245 :       GEN R = FlxqXQ_halfFrobenius_i(FlxX_Flx_add(t, r, p), xp, Xp, S, T, p, pi);
    2336         245 :       f = FlxqX_gcd_pre(FlxX_Flx_sub(R, pol1_Flx(vT), p), Sp, T, p, pi);
    2337         245 :       if (degpol(f) > 0 && degpol(f) < n) break;
    2338          70 :       set_avma(btop2);
    2339             :     }
    2340         175 :     if (degpol(f) > 0 && degpol(f) < n) break;
    2341           0 :     set_avma(btop);
    2342             :   }
    2343         175 :   f = FlxqX_normalize_pre(f, T, p, pi);
    2344         175 :   ff = FlxqX_div_pre(Sp, f , T, p, pi);
    2345         175 :   FlxqX_edf_simple(f, xp, Xp, Xq, d, T, p, pi, V, idx);
    2346         175 :   FlxqX_edf_simple(ff, xp, Xp, Xq, d, T, p, pi, V, idx+degpol(f)/d);
    2347             : }
    2348             : 
    2349             : static GEN
    2350         342 : FlxqX_factor_Shoup(GEN S, GEN xp, GEN T, ulong p, ulong pi)
    2351             : {
    2352         342 :   long i, n, s = 0;
    2353             :   GEN X, Xp, Xq, D, V;
    2354         342 :   long dT = get_Flx_degree(T), vT = get_Flx_var(T);
    2355         342 :   long e = expi(powuu(p, dT));
    2356             :   pari_timer ti;
    2357         342 :   n = get_FlxqX_degree(S);
    2358         342 :   S = FlxqX_get_red_pre(S, T, p, pi);
    2359         342 :   if (DEBUGLEVEL>=6) timer_start(&ti);
    2360         342 :   X  = polx_FlxX(get_FlxqX_var(S), vT);
    2361         342 :   Xp = FlxqXQ_powu_pre(X, p, S, T, p, pi);
    2362         342 :   Xq = FlxqXQ_Frobenius(xp, Xp, S, T, p, pi);
    2363         342 :   if (DEBUGLEVEL>=6) timer_printf(&ti,"FlxqX_Frobenius");
    2364         342 :   D = FlxqX_ddf_Shoup(S, Xq, T, p, pi);
    2365         342 :   if (DEBUGLEVEL>=6) timer_printf(&ti,"FlxqX_ddf_Shoup");
    2366         342 :   s = ddf_to_nbfact(D);
    2367         342 :   V = cgetg(s+1, t_COL);
    2368        1386 :   for (i = 1, s = 1; i <= n; i++)
    2369             :   {
    2370        1044 :     GEN Di = gel(D,i);
    2371        1044 :     long ni = degpol(Di), ri = ni/i;
    2372        1044 :     if (ni == 0) continue;
    2373         370 :     Di = FlxqX_normalize_pre(Di, T, p, pi);
    2374         370 :     if (ni == i) { gel(V, s++) = Di; continue; }
    2375         153 :     if (ri <= e*expu(e))
    2376         146 :       FlxqX_edf(Di, xp, Xp, Xq, i, T, p, pi, V, s);
    2377             :     else
    2378           7 :       FlxqX_edf_simple(Di, xp, Xp, Xq, i, T, p, pi, V, s);
    2379         153 :     if (DEBUGLEVEL>=6) timer_printf(&ti,"FlxqX_edf(%ld)",i);
    2380         153 :     s += ri;
    2381             :   }
    2382         342 :   return V;
    2383             : }
    2384             : 
    2385             : static GEN
    2386       12759 : FlxqX_factor_Cantor(GEN f, GEN T, ulong p)
    2387             : {
    2388       12759 :   ulong pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
    2389             :   GEN xp, E, F, V;
    2390             :   long i, j, l;
    2391       12759 :   T = Flx_get_red_pre(T, p, pi);
    2392       12759 :   f = FlxqX_normalize_pre(f, T, p, pi);
    2393       12759 :   switch(degpol(f))
    2394             :   {
    2395          21 :     case -1: retmkmat2(mkcol(f), mkvecsmall(1));
    2396          21 :     case 0: return trivial_fact();
    2397          21 :     case 1: retmkmat2(mkcol(f), mkvecsmall(1));
    2398         825 :     case 2: return FlxqX_factor_2(f, T, p, pi);
    2399             :   }
    2400       11871 :   if (FlxY_degreex(f) <= 0) return Flx_factorff_i(FlxX_to_Flx(f), T, p);
    2401         293 :   xp = Flx_Frobenius_pre(T, p, pi);
    2402         293 :   V = FlxqX_factor_squarefree_i(f, xp, get_Flx_mod(T), p, pi);
    2403         293 :   l = lg(V);
    2404         293 :   F = cgetg(l, t_VEC);
    2405         293 :   E = cgetg(l, t_VEC);
    2406        1027 :   for (i=1, j=1; i < l; i++)
    2407         734 :     if (degpol(gel(V,i)))
    2408             :     {
    2409         342 :       GEN Fj = FlxqX_factor_Shoup(gel(V,i), xp, T, p, pi);
    2410         342 :       gel(F, j) = Fj;
    2411         342 :       gel(E, j) = const_vecsmall(lg(Fj)-1, i);
    2412         342 :       j++;
    2413             :     }
    2414         293 :   return sort_factor_pol(FE_concat(F,E,j), cmp_Flx);
    2415             : }
    2416             : 
    2417             : /* T must be squarefree mod p*/
    2418             : GEN
    2419         147 : FlxqX_nbfact_by_degree(GEN f, long *nb, GEN T, ulong p)
    2420             : {
    2421             :   GEN Xq, D;
    2422             :   pari_timer ti;
    2423         147 :   long i, s, n = get_FlxqX_degree(f);
    2424         147 :   ulong pi = SMALL_ULONG(p)? 0: get_Fl_red(p);
    2425         147 :   GEN V = const_vecsmall(n, 0);
    2426         147 :   pari_sp av = avma;
    2427         147 :   T = Flx_get_red_pre(T, p, pi);
    2428         147 :   f = FlxqX_get_red_pre(f, T, p, pi);
    2429         147 :   if (DEBUGLEVEL>=6) timer_start(&ti);
    2430         147 :   Xq = FlxqX_Frobenius_pre(f, T, p, pi);
    2431         147 :   if (DEBUGLEVEL>=6) timer_printf(&ti,"FlxqX_Frobenius");
    2432         147 :   D = FlxqX_ddf_Shoup(f, Xq, T, p, pi);
    2433         147 :   if (DEBUGLEVEL>=6) timer_printf(&ti,"FlxqX_ddf_Shoup");
    2434         945 :   for (i = 1, s = 0; i <= n; i++) { V[i] = degpol(gel(D,i))/i; s += V[i]; }
    2435         147 :   *nb = s; set_avma(av); return V;
    2436             : }
    2437             : 
    2438             : long
    2439           0 : FlxqX_nbfact_Frobenius(GEN S, GEN Xq, GEN T, ulong p)
    2440             : {
    2441           0 :   pari_sp av = avma;
    2442           0 :   GEN u = get_FlxqX_mod(S);
    2443             :   long s;
    2444           0 :   if (FlxY_degreex(u) <= 0)
    2445           0 :     s = Flx_nbfactff(FlxX_to_Flx(u), T, p);
    2446             :   else
    2447           0 :     s = ddf_to_nbfact(FlxqX_ddf_Shoup(S, Xq, T, p,
    2448           0 :                                       SMALL_ULONG(p)? 0: get_Fl_red(p)));
    2449           0 :   return gc_long(av,s);
    2450             : }
    2451             : 
    2452             : long
    2453        7084 : FlxqX_nbfact(GEN S, GEN T, ulong p)
    2454             : {
    2455        7084 :   pari_sp av = avma;
    2456        7084 :   GEN u = get_FlxqX_mod(S);
    2457             :   long s;
    2458        7084 :   if (FlxY_degreex(u) <= 0)
    2459        7084 :     s = Flx_nbfactff(FlxX_to_Flx(u), T, p);
    2460             :   else
    2461           0 :     s = ddf_to_nbfact(FlxqX_ddf_Shoup(S, FlxqX_Frobenius(S, T, p), T, p,
    2462           0 :                                       SMALL_ULONG(p)? 0: get_Fl_red(p)));
    2463        7084 :   return gc_long(av,s);
    2464             : }
    2465             : 
    2466             : GEN
    2467         187 : FlxqX_factor(GEN x, GEN T, ulong p)
    2468             : {
    2469         187 :   pari_sp av = avma;
    2470         187 :   return gerepilecopy(av, FlxqX_factor_Cantor(x, T, p));
    2471             : }
    2472             : 
    2473             : GEN
    2474         133 : F2xqX_factor(GEN x, GEN T)
    2475             : {
    2476         133 :   pari_sp av = avma;
    2477         133 :   return gerepilecopy(av, F2xqX_factor_Cantor(x, T));
    2478             : }
    2479             : 
    2480             : long
    2481         187 : FpXQX_ddf_degree(GEN S, GEN XP, GEN T, GEN p)
    2482             : {
    2483         187 :   pari_sp av = avma;
    2484             :   GEN X, b, g, xq, q;
    2485             :   long i, j, n, v, B, l, m, bo, ro;
    2486             :   pari_timer ti;
    2487             :   hashtable h;
    2488             : 
    2489         187 :   n = get_FpXQX_degree(S); v = get_FpXQX_var(S);
    2490         187 :   X = pol_x(v);
    2491         187 :   if (gequal(X,XP)) return 1;
    2492         187 :   B = n/2;
    2493         187 :   l = usqrt(B);
    2494         187 :   m = (B+l-1)/l;
    2495         187 :   T = FpX_get_red(T, p);
    2496         187 :   S = FpXQX_get_red(S, T, p);
    2497         187 :   hash_init_GEN(&h, l+2, gequal, 1);
    2498         187 :   hash_insert_long(&h, X,  0);
    2499         187 :   hash_insert_long(&h, simplify_shallow(XP), 1);
    2500         187 :   bo = brent_kung_optpow(n, l-1, 1);
    2501         187 :   ro = l<=1 ? 0: (bo-1)/(l-1) + ((n-1)/bo);
    2502         187 :   q = powiu(p, get_FpX_degree(T));
    2503         187 :   if (DEBUGLEVEL>=7) timer_start(&ti);
    2504         187 :   b = XP;
    2505         187 :   if (expi(q) > ro)
    2506             :   {
    2507         187 :     xq = FpXQXQ_powers(b, brent_kung_optpow(n, l-1, 1),  S, T, p);
    2508         187 :     if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_degree: xq baby");
    2509           0 :   } else xq = NULL;
    2510         695 :   for (i = 3; i <= l+1; i++)
    2511             :   {
    2512         523 :     b = xq ? FpXQX_FpXQXQV_eval(b, xq, S, T, p)
    2513         523 :            : FpXQXQ_pow(b, q, S, T, p);
    2514         523 :     if (gequal(b,X)) return gc_long(av,i-1);
    2515         508 :     hash_insert_long(&h, simplify_shallow(b), i-1);
    2516             :   }
    2517         172 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_degree: baby");
    2518         172 :   g = b;
    2519         172 :   xq = FpXQXQ_powers(g, brent_kung_optpow(n, m, 1),  S, T, p);
    2520         172 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_degree: xq giant");
    2521         756 :   for(i = 2; i <= m+1; i++)
    2522             :   {
    2523         691 :     g = FpXQX_FpXQXQV_eval(g, xq, S, T, p);
    2524         691 :     if (hash_haskey_long(&h, simplify_shallow(g), &j)) return gc_long(av,l*i-j);
    2525             :   }
    2526          65 :   return gc_long(av,n);
    2527             : }
    2528             : 
    2529             : static GEN
    2530          71 : FpXQX_ddf_Shoup(GEN S, GEN Xq, GEN T, GEN p)
    2531             : {
    2532          71 :   pari_sp av = avma;
    2533             :   GEN b, g, h, F, f, Sr, xq, q;
    2534             :   long i, j, n, v, bo, ro;
    2535             :   long B, l, m;
    2536             :   pari_timer ti;
    2537          71 :   n = get_FpXQX_degree(S); v = get_FpXQX_var(S);
    2538          71 :   if (n == 0) return cgetg(1, t_VEC);
    2539          71 :   if (n == 1) return mkvec(get_FpXQX_mod(S));
    2540          64 :   B = n/2;
    2541          64 :   l = usqrt(B);
    2542          64 :   m = (B+l-1)/l;
    2543          64 :   S = FpXQX_get_red(S, T, p);
    2544          64 :   b = cgetg(l+2, t_VEC);
    2545          64 :   gel(b, 1) = pol_x(v);
    2546          64 :   gel(b, 2) = Xq;
    2547          64 :   bo = brent_kung_optpow(n, l-1, 1);
    2548          64 :   ro = l<=1 ? 0: (bo-1)/(l-1) + ((n-1)/bo);
    2549          64 :   q = powiu(p, get_FpX_degree(T));
    2550          64 :   if (DEBUGLEVEL>=7) timer_start(&ti);
    2551          64 :   if (expi(q) <= ro)
    2552           0 :     for (i = 3; i <= l+1; i++)
    2553           0 :       gel(b, i) = FpXQXQ_pow(gel(b, i-1), q, S, T, p);
    2554             :   else
    2555             :   {
    2556          64 :     xq = FpXQXQ_powers(gel(b, 2), bo, S, T, p);
    2557          64 :     if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_Shoup: xq baby");
    2558          85 :     for (i = 3; i <= l+1; i++)
    2559          21 :       gel(b, i) = FpXQX_FpXQXQV_eval(gel(b, i-1), xq, S, T, p);
    2560             :   }
    2561          64 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_Shoup: baby");
    2562          64 :   xq = FpXQXQ_powers(gel(b, l+1), brent_kung_optpow(n, m-1, 1), S, T, p);
    2563          64 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_Shoup: xq giant");
    2564          64 :   g = cgetg(m+1, t_VEC);
    2565          64 :   gel(g, 1) = gel(xq, 2);
    2566         108 :   for(i = 2; i <= m; i++)
    2567          44 :     gel(g, i) = FpXQX_FpXQXQV_eval(gel(g, i-1), xq, S, T, p);
    2568          64 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_Shoup: giant");
    2569          64 :   h = cgetg(m+1, t_VEC);
    2570         172 :   for (j = 1; j <= m; j++)
    2571             :   {
    2572         108 :     pari_sp av = avma;
    2573         108 :     GEN gj = gel(g, j);
    2574         108 :     GEN e = FpXX_sub(gj, gel(b, 1), p);
    2575         171 :     for (i = 2; i <= l; i++)
    2576          63 :       e = FpXQXQ_mul(e, FpXX_sub(gj, gel(b, i), p), S, T, p);
    2577         108 :     gel(h, j) = gerepileupto(av, e);
    2578             :   }
    2579          64 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_Shoup: diff");
    2580          64 :   Sr = get_FpXQX_mod(S);
    2581          64 :   F = cgetg(m+1, t_VEC);
    2582         172 :   for (j = 1; j <= m; j++)
    2583             :   {
    2584         108 :     GEN u = FpXQX_gcd(Sr, gel(h,j), T, p);
    2585         108 :     if (degpol(u))
    2586             :     {
    2587          78 :       u = FpXQX_normalize(u, T, p);
    2588          78 :       Sr = FpXQX_div(Sr, u, T, p);
    2589             :     }
    2590         108 :     gel(F,j) = u;
    2591             :   }
    2592          64 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_Shoup: F");
    2593          64 :   f = const_vec(n, pol_1(v));
    2594         172 :   for (j = 1; j <= m; j++)
    2595             :   {
    2596         108 :     GEN e = gel(F, j);
    2597         150 :     for (i=l-1; i >= 0; i--)
    2598             :     {
    2599         150 :       GEN u = FpXQX_gcd(e, FpXX_sub(gel(g, j), gel(b, i+1), p), T, p);
    2600         150 :       if (degpol(u))
    2601             :       {
    2602          99 :         gel(f, l*j-i) = u = FpXQX_normalize(u, T, p);
    2603          99 :         e = FpXQX_div(e, u, T, p);
    2604             :       }
    2605         150 :       if (!degpol(e)) break;
    2606             :     }
    2607             :   }
    2608          64 :   if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_ddf_Shoup: f");
    2609          64 :   if (degpol(Sr)) gel(f, degpol(Sr)) = Sr;
    2610          64 :   return gerepilecopy(av, f);
    2611             : }
    2612             : 
    2613             : static GEN
    2614          49 : FpXQX_ddf_i(GEN f, GEN T, GEN p)
    2615             : {
    2616             :   GEN Xq;
    2617          49 :   if (!get_FpXQX_degree(f)) return cgetg(1,t_VEC);
    2618          49 :   f = FpXQX_get_red(f, T, p);
    2619          49 :   Xq = FpXQX_Frobenius(f, T, p);
    2620          49 :   return FpXQX_ddf_Shoup(f, Xq, T, p);
    2621             : }
    2622             : 
    2623             : static GEN
    2624           7 : FpXQX_ddf_raw(GEN f, GEN T, GEN p)
    2625             : {
    2626           7 :   if (lgefint(p)==3)
    2627             :   {
    2628           0 :     ulong pp = p[2];
    2629             :     GEN M;
    2630           0 :     long vT = get_FpX_var(T);
    2631           0 :     if (pp==2)
    2632             :     {
    2633           0 :       M = F2xqX_ddf(ZXX_to_F2xX(f, vT),  ZX_to_F2x(get_FpX_mod(T)));
    2634           0 :       return mkvec2(F2xXC_to_ZXXC(gel(M,1)), gel(M,2));
    2635             :     }
    2636           0 :     M = FlxqX_ddf(ZXX_to_FlxX(f, pp, vT),  ZXT_to_FlxT(T, pp), pp);
    2637           0 :     return mkvec2(FlxXC_to_ZXXC(gel(M,1)), gel(M,2));
    2638             :   }
    2639           7 :   T = FpX_get_red(T, p);
    2640           7 :   f = FpXQX_normalize(get_FpXQX_mod(f), T, p);
    2641           7 :   return ddf_to_ddf2( FpXQX_ddf_i(f, T, p) );
    2642             : }
    2643             : 
    2644             : GEN
    2645           7 : FpXQX_ddf(GEN x, GEN T, GEN p)
    2646           7 : { pari_sp av = avma; return gerepilecopy(av, FpXQX_ddf_raw(x,T,p)); }
    2647             : 
    2648             : static GEN
    2649          42 : FpXQX_degfact_raw(GEN f, GEN T, GEN p)
    2650             : {
    2651             :   GEN V;
    2652             :   long i,l;
    2653          42 :   if (lgefint(p)==3)
    2654             :   {
    2655           0 :     ulong pp = p[2];
    2656           0 :     long vT = get_FpX_var(T);
    2657           0 :     if (pp==2)
    2658           0 :       return F2xqX_degfact(ZXX_to_F2xX(f, vT),  ZX_to_F2x(get_FpX_mod(T)));
    2659             :     else
    2660           0 :       return FlxqX_degfact(ZXX_to_FlxX(f, pp, vT),  ZXT_to_FlxT(T, pp), pp);
    2661             :   }
    2662          42 :   T = FpX_get_red(T, p);
    2663          42 :   f = FpXQX_normalize(get_FpXQX_mod(f), T, p);
    2664          42 :   V = FpXQX_factor_Yun(f, T, p); l = lg(V);
    2665          84 :   for (i=1; i < l; i++) gel(V,i) = FpXQX_ddf_i(gel(V,i), T, p);
    2666          42 :   return vddf_to_simplefact(V, degpol(f));
    2667             : }
    2668             : 
    2669             : GEN
    2670          42 : FpXQX_degfact(GEN x, GEN T, GEN p)
    2671          42 : { pari_sp av = avma; return gerepilecopy(av, FpXQX_degfact_raw(x,T,p)); }
    2672             : 
    2673             : static void
    2674          23 : FpXQX_edf_rec(GEN S, GEN xp, GEN Xp, GEN hp, GEN t, long d, GEN T, GEN p, GEN V, long idx)
    2675             : {
    2676          23 :   GEN Sp = get_FpXQX_mod(S);
    2677             :   GEN u1, u2, f1, f2;
    2678             :   GEN h;
    2679          23 :   h = FpXQX_get_red(hp, T, p);
    2680          23 :   t = FpXQX_rem(t, S, T, p);
    2681          23 :   Xp = FpXQX_rem(Xp, h, T, p);
    2682          23 :   u1 = FpXQX_roots_split(hp, xp, Xp, h, T, p);
    2683          23 :   f1 = FpXQX_gcd(FpXQX_FpXQXQ_eval(u1, t, S, T, p), Sp, T, p);
    2684          23 :   f1 = FpXQX_normalize(f1, T, p);
    2685          23 :   u2 = FpXQX_div(hp, u1, T, p);
    2686          23 :   f2 = FpXQX_div(Sp, f1, T, p);
    2687          23 :   if (degpol(u1)==1)
    2688          15 :     gel(V, idx) = f1;
    2689             :   else
    2690           8 :     FpXQX_edf_rec(FpXQX_get_red(f1, T, p), xp, Xp, u1, t, d, T, p, V, idx);
    2691          23 :   idx += degpol(f1)/d;
    2692          23 :   if (degpol(u2)==1)
    2693          16 :     gel(V, idx) = f2;
    2694             :   else
    2695           7 :     FpXQX_edf_rec(FpXQX_get_red(f2, T, p), xp, Xp, u2, t, d, T, p, V, idx);
    2696          23 : }
    2697             : 
    2698             : static void
    2699           8 : FpXQX_edf(GEN Sp, GEN xp, GEN Xp, GEN Xq, long d, GEN T, GEN p, GEN V, long idx)
    2700             : {
    2701           8 :   long n = degpol(Sp), r = n/d, vS = varn(Sp);
    2702             :   GEN S, h, t;
    2703             :   pari_timer ti;
    2704           8 :   if (r==1) { gel(V, idx) = Sp; return; }
    2705           8 :   S = FpXQX_get_red(Sp, T, p);
    2706           8 :   Xp = FpXQX_rem(Xp, S, T, p);
    2707           8 :   Xq = FpXQX_rem(Xq, S, T, p);
    2708           8 :   if (DEBUGLEVEL>=7) timer_start(&ti);
    2709             :   do
    2710             :   {
    2711           8 :     GEN g = random_FpXQX(n, vS, T, p);
    2712           8 :     t = gel(FpXQXQ_auttrace(mkvec2(Xq, g), d, S, T, p), 2);
    2713           8 :     if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_edf: FpXQXQ_auttrace");
    2714           8 :     h = FpXQXQ_minpoly(t, S, T, p);
    2715           8 :     if (DEBUGLEVEL>=7) timer_printf(&ti,"FpXQX_edf: FpXQXQ_minpoly");
    2716           8 :   } while (degpol(h) != r);
    2717           8 :   Xp = FpXQXQ_pow(pol_x(vS), p, h, T, p);
    2718           8 :   FpXQX_edf_rec(S, xp, Xp, h, t, d, T, p, V, idx);
    2719             : }
    2720             : 
    2721             : static void
    2722           0 : FpXQX_edf_simple(GEN Sp, GEN xp, GEN Xp, GEN Xq, long d, GEN T, GEN p, GEN V, long idx)
    2723             : {
    2724           0 :   long v = varn(Sp), n = degpol(Sp), r = n/d;
    2725             :   GEN S, f, ff;
    2726           0 :   long vT = get_FpX_var(T), dT = get_FpX_degree(T);
    2727           0 :   if (r==1) { gel(V, idx) = Sp; return; }
    2728           0 :   S = FpXQX_get_red(Sp, T, p);
    2729           0 :   Xp = FpXQX_rem(Xp, S, T, p);
    2730           0 :   Xq = FpXQX_rem(Xq, S, T, p);
    2731             :   while (1)
    2732           0 :   {
    2733           0 :     pari_sp btop = avma;
    2734             :     long i;
    2735           0 :     GEN g = random_FpXQX(n, v, T, p);
    2736           0 :     GEN t = gel(FpXQXQ_auttrace(mkvec2(Xq, g), d, S, T, p), 2);
    2737           0 :     if (lgpol(t) == 0) continue;
    2738           0 :     for(i=1; i<=10; i++)
    2739             :     {
    2740           0 :       pari_sp btop2 = avma;
    2741           0 :       GEN r = random_FpX(dT, vT, p);
    2742           0 :       GEN R = FpXQXQ_halfFrobenius_i(FqX_Fq_add(t, r, T, p), xp, Xp, S, T, p);
    2743           0 :       f = FpXQX_gcd(FqX_Fq_add(R, gen_m1, T, p), Sp, T, p);
    2744           0 :       if (degpol(f) > 0 && degpol(f) < n) break;
    2745           0 :       set_avma(btop2);
    2746             :     }
    2747           0 :     if (degpol(f) > 0 && degpol(f) < n) break;
    2748           0 :     set_avma(btop);
    2749             :   }
    2750           0 :   f = FpXQX_normalize(f, T, p);
    2751           0 :   ff = FpXQX_div(Sp, f , T, p);
    2752           0 :   FpXQX_edf_simple(f,  xp, Xp, Xq, d, T, p, V, idx);
    2753           0 :   FpXQX_edf_simple(ff, xp, Xp, Xq, d, T, p, V, idx+degpol(f)/d);
    2754             : }
    2755             : 
    2756             : static GEN
    2757          22 : FpXQX_factor_Shoup(GEN S, GEN xp, GEN T, GEN p)
    2758             : {
    2759          22 :   long i, n, s = 0, dT = get_FpX_degree(T), e = expi(powiu(p, dT));
    2760             :   GEN X, Xp, Xq, D, V;
    2761             :   pari_timer ti;
    2762          22 :   n = get_FpXQX_degree(S);
    2763          22 :   S = FpXQX_get_red(S, T, p);
    2764          22 :   if (DEBUGLEVEL>=6) timer_start(&ti);
    2765          22 :   X  = pol_x(get_FpXQX_var(S));
    2766          22 :   Xp = FpXQXQ_pow(X, p, S, T, p);
    2767          22 :   Xq = FpXQXQ_Frobenius(xp, Xp, S, T, p);
    2768          22 :   if (DEBUGLEVEL>=6) timer_printf(&ti,"FpXQX_Frobenius");
    2769          22 :   D = FpXQX_ddf_Shoup(S, Xq, T, p);
    2770          22 :   if (DEBUGLEVEL>=6) timer_printf(&ti,"FpXQX_ddf_Shoup");
    2771          22 :   s = ddf_to_nbfact(D);
    2772          22 :   V = cgetg(s+1, t_COL);
    2773         140 :   for (i = 1, s = 1; i <= n; i++)
    2774             :   {
    2775         118 :     GEN Di = gel(D,i);
    2776         118 :     long ni = degpol(Di), ri = ni/i;
    2777         118 :     if (ni == 0) continue;
    2778          50 :     Di = FpXQX_normalize(Di, T, p);
    2779          50 :     if (ni == i) { gel(V, s++) = Di; continue; }
    2780           8 :     if (ri <= e*expu(e))
    2781           8 :       FpXQX_edf(Di, xp, Xp, Xq, i, T, p, V, s);
    2782             :     else
    2783           0 :       FpXQX_edf_simple(Di, xp, Xp, Xq, i, T, p, V, s);
    2784           8 :     if (DEBUGLEVEL>=6) timer_printf(&ti,"FpXQX_edf(%ld)",i);
    2785           8 :     s += ri;
    2786             :   }
    2787          22 :   return V;
    2788             : }
    2789             : 
    2790             : static GEN
    2791         177 : FpXQX_factor_Cantor(GEN f, GEN T, GEN p)
    2792             : {
    2793             :   GEN xp, E, F, V;
    2794             :   long i, j, l;
    2795         177 :   T = FpX_get_red(T, p);
    2796         177 :   f = FpXQX_normalize(f, T, p);
    2797         177 :   switch(degpol(f))
    2798             :   {
    2799          14 :     case -1: retmkmat2(mkcol(f), mkvecsmall(1));
    2800          14 :     case 0: return trivial_fact();
    2801          21 :     case 1: retmkmat2(mkcol(f), mkvecsmall(1));
    2802          43 :     case 2: return FpXQX_factor_2(f, T, p);
    2803             :   }
    2804          85 :   if (isabsolutepol(f)) return FpX_factorff_i(simplify_shallow(f), T, p);
    2805          22 :   xp = FpX_Frobenius(T, p);
    2806          22 :   V = FpXQX_factor_Yun(f, T, p);
    2807          22 :   l = lg(V);
    2808          22 :   F = cgetg(l, t_VEC);
    2809          22 :   E = cgetg(l, t_VEC);
    2810          44 :   for (i=1, j=1; i < l; i++)
    2811          22 :     if (degpol(gel(V,i)))
    2812             :     {
    2813          22 :       GEN Fj = FpXQX_factor_Shoup(gel(V,i), xp, T, p);
    2814          22 :       gel(E,j) = const_vecsmall(lg(Fj)-1, i);
    2815          22 :       gel(F,j) = Fj; j++;
    2816             :     }
    2817          22 :   return sort_factor_pol(FE_concat(F,E,j), cmp_RgX);
    2818             : }
    2819             : 
    2820             : long
    2821           0 : FpXQX_nbfact_Frobenius(GEN S, GEN Xq, GEN T, GEN p)
    2822             : {
    2823           0 :   pari_sp av = avma;
    2824           0 :   GEN u = get_FpXQX_mod(S);
    2825             :   long s;
    2826           0 :   if (lgefint(p)==3)
    2827             :   {
    2828           0 :     ulong pp = p[2];
    2829           0 :     long vT = get_FpX_var(T);
    2830           0 :     GEN Sp = ZXXT_to_FlxXT(S,pp,vT), Xqp = ZXX_to_FlxX(Xq,pp,vT);
    2831           0 :     s = FlxqX_nbfact_Frobenius(Sp, Xqp, ZXT_to_FlxT(T,pp), pp);
    2832             :   }
    2833           0 :   else if (isabsolutepol(u))
    2834           0 :     s = FpX_nbfactff(simplify_shallow(u), T, p);
    2835             :   else
    2836           0 :     s = ddf_to_nbfact(FpXQX_ddf_Shoup(S, Xq, T, p));
    2837           0 :   return gc_long(av,s);
    2838             : }
    2839             : 
    2840             : long
    2841           0 : FpXQX_nbfact(GEN S, GEN T, GEN p)
    2842             : {
    2843           0 :   pari_sp av = avma;
    2844           0 :   GEN u = get_FpXQX_mod(S);
    2845             :   long s;
    2846           0 :   if (lgefint(p)==3)
    2847             :   {
    2848           0 :     ulong pp = p[2];
    2849           0 :     long vT = get_FpX_var(T);
    2850           0 :     s = FlxqX_nbfact(ZXXT_to_FlxXT(S,pp,vT), ZXT_to_FlxT(T,pp), pp);
    2851             :   }
    2852           0 :   else if (isabsolutepol(u))
    2853           0 :     s = FpX_nbfactff(simplify_shallow(u), T, p);
    2854             :   else
    2855           0 :     s = ddf_to_nbfact(FpXQX_ddf_Shoup(S, FpXQX_Frobenius(S, T, p), T, p));
    2856           0 :   return gc_long(av,s);
    2857             : }
    2858             : long
    2859           0 : FqX_nbfact(GEN u, GEN T, GEN p)
    2860           0 : { return T ? FpXQX_nbfact(u, T, p): FpX_nbfact(u, p); }
    2861             : 
    2862             : static GEN
    2863           0 : FpXQX_factor_Berlekamp_i(GEN f, GEN T, GEN p)
    2864             : {
    2865           0 :   if (lgefint(p)==3)
    2866             :   {
    2867           0 :     ulong pp = p[2];
    2868             :     GEN M;
    2869           0 :     long vT = get_FpX_var(T);
    2870           0 :     if (pp==2)
    2871             :     {
    2872           0 :       M = F2xqX_factor_Cantor(ZXX_to_F2xX(f, vT),  ZX_to_F2x(get_FpX_mod(T)));
    2873           0 :       return mkvec2(F2xXC_to_ZXXC(gel(M,1)), gel(M,2));
    2874             :     }
    2875           0 :     M = FlxqX_Berlekamp_i(ZXX_to_FlxX(f, pp, vT),  ZXT_to_FlxT(T, pp), pp);
    2876           0 :     return mkvec2(FlxXC_to_ZXXC(gel(M,1)), gel(M,2));
    2877             :   }
    2878           0 :   return FpXQX_Berlekamp_i(f, T, p);
    2879             : }
    2880             : GEN
    2881           0 : FpXQX_factor_Berlekamp(GEN x, GEN T, GEN p)
    2882           0 : { pari_sp av = avma; return gerepilecopy(av, FpXQX_factor_Berlekamp_i(x,T,p)); }
    2883             : 
    2884             : static GEN
    2885       13855 : FpXQX_factor_i(GEN f, GEN T, GEN p)
    2886             : {
    2887       13855 :   if (lgefint(p)==3)
    2888             :   {
    2889       13678 :     ulong pp = p[2];
    2890             :     GEN M;
    2891       13678 :     long vT = get_FpX_var(T);
    2892       13678 :     if (pp==2)
    2893             :     {
    2894        1106 :       M = F2xqX_factor_Cantor(ZXX_to_F2xX(f, vT),  ZX_to_F2x(get_FpX_mod(T)));
    2895        1099 :       return mkvec2(F2xXC_to_ZXXC(gel(M,1)), gel(M,2));
    2896             :     }
    2897       12572 :     M = FlxqX_factor_Cantor(ZXX_to_FlxX(f, pp, vT),  ZXT_to_FlxT(T, pp), pp);
    2898       12565 :     return mkvec2(FlxXC_to_ZXXC(gel(M,1)), gel(M,2));
    2899             :   }
    2900         177 :   return FpXQX_factor_Cantor(f, T, p);
    2901             : }
    2902             : GEN
    2903       13853 : FpXQX_factor(GEN x, GEN T, GEN p)
    2904       13853 : { pari_sp av = avma; return gerepilecopy(av, FpXQX_factor_i(x,T,p)); }
    2905             : 
    2906             : long
    2907       10693 : FlxqX_is_squarefree(GEN P, GEN T, ulong p)
    2908             : {
    2909       10693 :   pari_sp av = avma;
    2910       10693 :   GEN z = FlxqX_gcd(P, FlxX_deriv(P, p), T, p);
    2911       10693 :   return gc_long(av, degpol(z)==0);
    2912             : }
    2913             : 
    2914             : long
    2915      683802 : FqX_is_squarefree(GEN P, GEN T, GEN p)
    2916             : {
    2917      683802 :   pari_sp av = avma;
    2918      683802 :   GEN z = FqX_gcd(P, FqX_deriv(P, T, p), T, p);
    2919      683802 :   return gc_long(av, degpol(z)==0);
    2920             : }
    2921             : 
    2922             : /* as RgX_to_FpXQ(FqX_to_FFX), leaving alone t_FFELT */
    2923             : static GEN
    2924         350 : RgX_to_FFX(GEN x, GEN ff)
    2925             : {
    2926             :   long i, lx;
    2927         350 :   GEN p = FF_p_i(ff), T = FF_mod(ff), y =  cgetg_copy(x,&lx);
    2928         350 :   y[1] = x[1]; if (degpol(T) == 1) T = NULL;
    2929        1120 :   for (i = 2; i < lx; i++)
    2930             :   {
    2931         770 :     GEN c = gel(x,i);
    2932         770 :     gel(y,i) = typ(c) == t_FFELT? c: Fq_to_FF(Rg_to_Fq(c,T,p), ff);
    2933             :   }
    2934         350 :   return y;
    2935             : }
    2936             : 
    2937             : #define code(t1,t2) ((t1 << 6) | t2)
    2938             : /* Check types and replace F by a monic normalized FpX having the same roots
    2939             :  * Don't bother to make constant polynomials monic */
    2940             : static GEN
    2941      152625 : factmod_init(GEN f, GEN *pD, GEN *pT, GEN *pp)
    2942             : {
    2943      152625 :   const char *s = "factormod";
    2944      152625 :   GEN T, p, D = *pD;
    2945      152625 :   if (typ(f) != t_POL) pari_err_TYPE(s,f);
    2946      152625 :   if (!D)
    2947             :   {
    2948      103439 :     long pa, t = RgX_type(f, pp, pT, &pa);
    2949      103439 :     if (t == t_FFELT) return f;
    2950           0 :     *pD = gen_0;
    2951           0 :     if (t != t_INTMOD && t != code(t_POLMOD,t_INTMOD)) pari_err_TYPE(s,f);
    2952           0 :     return RgX_to_FqX(f, *pT, *pp);
    2953             :   }
    2954       49186 :   if (typ(D) == t_FFELT) { *pD = NULL; *pT = D; return RgX_to_FFX(f,D); }
    2955       48836 :   if (!ff_parse_Tp(D, &T, &p, 1)) pari_err_TYPE(s,D);
    2956       48822 :   if (T && varncmp(varn(T), varn(f)) <= 0)
    2957           0 :     pari_err_PRIORITY(s, T, "<=", varn(f));
    2958       48822 :   *pT = T; *pp = p; return RgX_to_FqX(f, T, p);
    2959             : }
    2960             : #undef code
    2961             : 
    2962             : int
    2963       49515 : ff_parse_Tp(GEN Tp, GEN *T, GEN *p, long red)
    2964             : {
    2965       49515 :   long t = typ(Tp);
    2966       49515 :   *T = *p = NULL;
    2967       49515 :   if (t == t_INT) { *p = Tp; return cmpiu(*p, 1) > 0; }
    2968         462 :   if (t != t_VEC || lg(Tp) != 3) return 0;
    2969         455 :   *T = gel(Tp,1);
    2970         455 :   *p = gel(Tp,2);
    2971         455 :   if (typ(*p) != t_INT)
    2972             :   {
    2973         420 :     if (typ(*T) != t_INT) return 0;
    2974         420 :     swap(*T, *p); /* support both [T,p] and [p,T] */
    2975             :   }
    2976         455 :   if (red) *T = RgX_to_FpX(*T, *p);
    2977         455 :   return cmpiu(*p, 1) > 0 && typ(*T) == t_POL && RgX_is_ZX(*T);
    2978             : }
    2979             : 
    2980             : static GEN
    2981        4851 : to_Fq(GEN x, GEN T, GEN p)
    2982             : {
    2983        4851 :   long i, lx, tx = typ(x);
    2984             :   GEN y;
    2985             : 
    2986        4851 :   if (tx == t_INT)
    2987         273 :     y = mkintmod(x,p);
    2988             :   else
    2989             :   {
    2990        4578 :     if (tx != t_POL) pari_err_TYPE("to_Fq",x);
    2991        4578 :     y = cgetg_copy(x,&lx); y[1] = x[1];
    2992      134204 :     for (i=2; i<lx; i++) gel(y,i) = mkintmod(gel(x,i), p);
    2993             :   }
    2994        4851 :   return mkpolmod(y, T);
    2995             : }
    2996             : static GEN
    2997         252 : to_Fq_pol(GEN x, GEN T, GEN p)
    2998             : {
    2999         252 :   long i, lx = lg(x);
    3000         252 :   if (lx == 2)
    3001             :   {
    3002          21 :     GEN y = cgetg(3,t_POL); y[1]=x[1];
    3003          21 :     gel(y,2) = mkintmod(gen_0,p); return y;
    3004             :   }
    3005         749 :   for (i = 2; i < lx; i++) gel(x,i) = to_Fq(gel(x,i),T,p);
    3006         231 :   return x;
    3007             : }
    3008             : static GEN
    3009         189 : to_Fq_fact(GEN fa, GEN T, GEN p)
    3010             : {
    3011         189 :   GEN P = gel(fa,1);
    3012         189 :   long j, l = lg(P);
    3013         189 :   p = icopy(p); T = FpX_to_mod(T, p);
    3014         441 :   for (j=1; j<l; j++) gel(P,j) = to_Fq_pol(gel(P,j), T,p);
    3015         189 :   return fa;
    3016             : }
    3017             : static GEN
    3018         224 : to_FqC(GEN P, GEN T, GEN p)
    3019             : {
    3020         224 :   long j, l = lg(P);
    3021         224 :   p = icopy(p); T = FpX_to_mod(T, p);
    3022        4557 :   for (j=1; j<l; j++) gel(P,j) = to_Fq(gel(P,j), T,p);
    3023         224 :   return P;
    3024             : }
    3025             : 
    3026             : GEN
    3027         910 : factmod(GEN f, GEN D)
    3028             : {
    3029             :   pari_sp av;
    3030             :   GEN y, F, P, E, T, p;
    3031         910 :   f = factmod_init(f, &D, &T,&p);
    3032         903 :   if (!D) return FFX_factor(f, T);
    3033         686 :   av = avma;
    3034         686 :   F = FqX_factor(f, T, p); P = gel(F,1); E = gel(F,2);
    3035         672 :   if (!T)
    3036             :   {
    3037         483 :     y = cgetg(3, t_MAT);
    3038         483 :     gel(y,1) = FpXC_to_mod(P, p);
    3039         483 :     gel(y,2) = Flc_to_ZC(E); return gerepileupto(av, y);
    3040             :   }
    3041         189 :   F = gerepilecopy(av, mkmat2(simplify_shallow(P), Flc_to_ZC(E)));
    3042         189 :   return to_Fq_fact(F, T, p);
    3043             : }
    3044             : GEN
    3045       47975 : simplefactmod(GEN f, GEN D)
    3046             : {
    3047       47975 :   pari_sp av = avma;
    3048             :   GEN T, p;
    3049       47975 :   f = factmod_init(f, &D, &T,&p);
    3050       47975 :   if (lg(f) <= 3) { set_avma(av); return trivial_fact(); }
    3051       47912 :   f = D? FqX_degfact(f, T, p): FFX_degfact(f, T);
    3052       47912 :   return gerepileupto(av, Flm_to_ZM(f));
    3053             : }
    3054             : static GEN
    3055          14 : sqf_to_fact(GEN f)
    3056             : {
    3057          14 :   long i, j, l = lg(f);
    3058          14 :   GEN P = cgetg(l, t_COL), E = cgetg(l, t_COL);
    3059          35 :   for (i = j = 1; i < l; i++)
    3060          21 :     if (degpol(gel(f,i)))
    3061             :     {
    3062          21 :       gel(P,j) = gel(f,i);
    3063          21 :       gel(E,j) = utoi(i); j++;
    3064             :     }
    3065          14 :   setlg(P,j);
    3066          14 :   setlg(E,j); return mkvec2(P,E);
    3067             : }
    3068             : 
    3069             : GEN
    3070          35 : factormodSQF(GEN f, GEN D)
    3071             : {
    3072          35 :   pari_sp av = avma;
    3073             :   GEN F, T, p;
    3074          35 :   f = factmod_init(f, &D, &T,&p);
    3075          35 :   if (lg(f) <= 3) { set_avma(av); return trivial_fact(); }
    3076          14 :   if (!D)
    3077           7 :     F = sqf_to_fact(FFX_factor_squarefree(f, T));
    3078             :   else
    3079             :   {
    3080           7 :     F = sqf_to_fact(FqX_factor_squarefree(f,T,p));
    3081           7 :     gel(F,1) = FqXC_to_mod(gel(F,1), T,p);
    3082             :   }
    3083          14 :   settyp(F,t_MAT); return gerepilecopy(av, F);
    3084             : }
    3085             : GEN
    3086          28 : factormodDDF(GEN f, GEN D)
    3087             : {
    3088          28 :   pari_sp av = avma;
    3089             :   GEN F, T, p;
    3090          28 :   f = factmod_init(f, &D, &T,&p);
    3091          28 :   if (lg(f) <= 3) { set_avma(av); return trivial_fact(); }
    3092          14 :   if (!D) return FFX_ddf(f, T);
    3093           7 :   F = FqX_ddf(f,T,p);
    3094           7 :   gel(F,1) = FqXC_to_mod(gel(F,1), T,p);
    3095           7 :   gel(F,2) = Flc_to_ZC(gel(F,2));
    3096           7 :   settyp(F, t_MAT); return gerepilecopy(av, F);
    3097             : }
    3098             : 
    3099             : GEN
    3100       48451 : factormod0(GEN f, GEN p, long flag)
    3101             : {
    3102       48451 :   if (flag == 0) return factmod(f,p);
    3103       47975 :   if (flag != 1) pari_err_FLAG("factormod");
    3104       47975 :   return simplefactmod(f,p);
    3105             : }
    3106             : GEN
    3107      103677 : polrootsmod(GEN f, GEN D)
    3108             : {
    3109             :   pari_sp av;
    3110             :   GEN y, T, p;
    3111      103677 :   f = factmod_init(f, &D, &T,&p);
    3112      103670 :   if (!D) return FFX_roots(f, T);
    3113         308 :   av = avma; y = FqX_roots(f, T, p);
    3114         280 :   if (!T) return gerepileupto(av, FpC_to_mod(y,p));
    3115         224 :   y = gerepilecopy(av, simplify_shallow(y));
    3116         224 :   return to_FqC(y, T, p);
    3117             : }
    3118             : 
    3119             : GEN /* OBSOLETE */
    3120           0 : rootmod0(GEN f, GEN p, long flag) { (void)flag; return polrootsmod(f,p); }
    3121             : GEN /* OBSOLETE */
    3122           0 : factorff(GEN f, GEN p, GEN T)
    3123             : {
    3124           0 :   pari_sp av = avma;
    3125           0 :   GEN D = (p && T)? mkvec2(T,p): NULL;
    3126           0 :   return gerepileupto(av, factmod(f,D));
    3127             : }
    3128             : GEN /* OBSOLETE */
    3129           0 : polrootsff(GEN f, GEN p, GEN T)
    3130             : {
    3131           0 :   pari_sp av = avma;
    3132           0 :   GEN D = (p && T)? mkvec2(T,p): NULL;
    3133           0 :   return gerepileupto(av, polrootsmod(f, D));
    3134             : }

Generated by: LCOV version 1.16