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-bordeaux1.fr machine (x86_64 architecture), and agregate them in the final report:

The target is 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 - FpX_factor.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.8.0 lcov report (development 16624-25b9976) Lines: 1591 1712 92.9 %
Date: 2014-06-24 Functions: 142 152 93.4 %
Legend: Lines: hit not hit | Branches: + taken - not taken # not executed Branches: 861 1132 76.1 %

           Branch data     Line data    Source code
       1                 :            : /* Copyright (C) 2012  The PARI group.
       2                 :            : 
       3                 :            : This file is part of the PARI/GP package.
       4                 :            : 
       5                 :            : PARI/GP is free software; you can redistribute it and/or modify it under the
       6                 :            : terms of the GNU General Public License as published by the Free Software
       7                 :            : Foundation. It is distributed in the hope that it will be useful, but WITHOUT
       8                 :            : ANY WARRANTY WHATSOEVER.
       9                 :            : 
      10                 :            : Check the License for details. You should have received a copy of it, along
      11                 :            : with the package; see the file 'COPYING'. If not, write to the Free Software
      12                 :            : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      13                 :            : 
      14                 :            : #include "pari.h"
      15                 :            : #include "paripriv.h"
      16                 :            : 
      17                 :            : /***********************************************************************/
      18                 :            : /**                                                                   **/
      19                 :            : /**               Factorisation over finite field                     **/
      20                 :            : /**                                                                   **/
      21                 :            : /***********************************************************************/
      22                 :            : 
      23                 :            : /*******************************************************************/
      24                 :            : /*                                                                 */
      25                 :            : /*           ROOTS MODULO a prime p (no multiplicities)            */
      26                 :            : /*                                                                 */
      27                 :            : /*******************************************************************/
      28                 :            : /* Check types and replace F by a monic normalized FpX having the same roots
      29                 :            :  * Don't bother to make constant polynomials monic */
      30                 :            : static void
      31                 :        175 : factmod_init(GEN *F, GEN p)
      32                 :            : {
      33         [ -  + ]:        175 :   if (typ(p)!=t_INT) pari_err_TYPE("factmod",p);
      34         [ -  + ]:        175 :   if (signe(p) < 0) pari_err_PRIME("factmod",p);
      35         [ -  + ]:        175 :   if (typ(*F)!=t_POL) pari_err_TYPE("factmod",*F);
      36         [ +  + ]:        175 :   if (lgefint(p) == 3)
      37                 :            :   {
      38                 :        110 :     ulong pp = p[2];
      39         [ -  + ]:        110 :     if (pp < 2) pari_err_PRIME("factmod", p);
      40                 :        110 :     *F = RgX_to_Flx(*F, pp);
      41         [ +  + ]:        110 :     if (lg(*F) > 3) *F = Flx_normalize(*F, pp);
      42                 :            :   }
      43                 :            :   else
      44                 :            :   {
      45                 :         65 :     *F = RgX_to_FpX(*F, p);
      46         [ +  + ]:         65 :     if (lg(*F) > 3) *F = FpX_normalize(*F, p);
      47                 :            :   }
      48                 :        175 : }
      49                 :            : /* as above, assume p prime and *F a ZX */
      50                 :            : static void
      51                 :     395779 : ZX_factmod_init(GEN *F, GEN p)
      52                 :            : {
      53         [ +  + ]:     395779 :   if (lgefint(p) == 3)
      54                 :            :   {
      55                 :     392901 :     ulong pp = p[2];
      56                 :     392901 :     *F = ZX_to_Flx(*F, pp);
      57         [ +  + ]:     392901 :     if (lg(*F) > 3) *F = Flx_normalize(*F, pp);
      58                 :            :   }
      59                 :            :   else
      60                 :            :   {
      61                 :       2878 :     *F = FpX_red(*F, p);
      62         [ +  + ]:       2878 :     if (lg(*F) > 3) *F = FpX_normalize(*F, p);
      63                 :            :   }
      64                 :     395779 : }
      65                 :            : 
      66                 :            : /* return 1,...,p-1 [not_0 = 1] or 0,...,p [not_0 = 0] */
      67                 :            : static GEN
      68                 :         25 : all_roots_mod_p(ulong p, int not_0)
      69                 :            : {
      70                 :            :   GEN r;
      71                 :            :   ulong i;
      72         [ +  + ]:         25 :   if (not_0) {
      73                 :         15 :     r = cgetg(p, t_VECSMALL);
      74         [ +  + ]:         45 :     for (i = 1; i < p; i++) r[i] = i;
      75                 :            :   } else {
      76                 :         10 :     r = cgetg(p+1, t_VECSMALL);
      77         [ +  + ]:         60 :     for (i = 0; i < p; i++) r[i+1] = i;
      78                 :            :   }
      79                 :         25 :   return r;
      80                 :            : }
      81                 :            : 
      82                 :            : /* X^n - 1 */
      83                 :            : static GEN
      84                 :         14 : Flx_Xnm1(long sv, long n, ulong p)
      85                 :            : {
      86                 :         14 :   GEN t = cgetg(n+3, t_VECSMALL);
      87                 :            :   long i;
      88                 :         14 :   t[1] = sv;
      89                 :         14 :   t[2] = p - 1;
      90         [ +  + ]:         57 :   for (i = 3; i <= n+1; i++) t[i] = 0;
      91                 :         14 :   t[i] = 1; return t;
      92                 :            : }
      93                 :            : /* X^n + 1 */
      94                 :            : static GEN
      95                 :          9 : Flx_Xn1(long sv, long n, ulong p)
      96                 :            : {
      97                 :          9 :   GEN t = cgetg(n+3, t_VECSMALL);
      98                 :            :   long i;
      99                 :            :   (void) p;
     100                 :          9 :   t[1] = sv;
     101                 :          9 :   t[2] = 1;
     102         [ +  + ]:         37 :   for (i = 3; i <= n+1; i++) t[i] = 0;
     103                 :          9 :   t[i] = 1; return t;
     104                 :            : }
     105                 :            : 
     106                 :            : static ulong
     107                 :          4 : Fl_nonsquare(ulong p)
     108                 :            : {
     109                 :          4 :   long k = 2;
     110                 :          4 :   for (;; k++)
     111                 :            :   {
     112                 :          8 :     long i = krouu(k, p);
     113         [ -  + ]:          8 :     if (!i) pari_err_PRIME("Fl_nonsquare",utoipos(p));
     114         [ +  + ]:          8 :     if (i < 0) return k;
     115                 :          4 :   }
     116                 :            : }
     117                 :            : 
     118                 :            : /* f monic Flx, f(0) != 0. Return a monic squarefree g with the same
     119                 :            :  * roots as f */
     120                 :            : static GEN
     121                 :         10 : Flx_cut_out_roots(GEN f, ulong p)
     122                 :            : {
     123                 :         10 :   GEN g = Flx_mod_Xnm1(f, p-1, p); /* f mod x^(p-1) - 1 */
     124 [ +  + ][ +  - ]:         10 :   if (g != f && degpol(g) >= 0) {
     125                 :          5 :     (void)Flx_valrem(g, &g); /* reduction may introduce 0 root */
     126                 :          5 :     g = Flx_gcd(g, Flx_Xnm1(g[1], p-1, p), p);
     127                 :          5 :     g = Flx_normalize(g, p);
     128                 :            :   }
     129                 :         10 :   return g;
     130                 :            : }
     131                 :            : 
     132                 :            : /* by checking f(0..p-1) */
     133                 :            : GEN
     134                 :         10 : Flx_roots_naive(GEN f, ulong p)
     135                 :            : {
     136                 :         10 :   long d, n = 0;
     137                 :         10 :   ulong s = 1UL, r;
     138                 :         10 :   GEN q, y = cgetg(degpol(f) + 1, t_VECSMALL);
     139                 :         10 :   pari_sp av2, av = avma;
     140                 :            : 
     141         [ -  + ]:         10 :   if (Flx_valrem(f, &f)) y[++n] = 0;
     142                 :         10 :   f = Flx_cut_out_roots(f, p);
     143                 :         10 :   d = degpol(f);
     144         [ -  + ]:         10 :   if (d < 0) return all_roots_mod_p(p, n == 0);
     145                 :         10 :   av2 = avma;
     146         [ +  + ]:        210 :   while (d > 1) /* d = current degree of f */
     147                 :            :   {
     148                 :        205 :     q = Flx_div_by_X_x(f, s, p, &r); /* TODO: FFT-type multi-evaluation */
     149         [ +  + ]:        205 :     if (r) avma = av2; else { y[++n] = s; d--; f = q; av2 = avma; }
     150         [ +  + ]:        205 :     if (++s == p) break;
     151                 :            :   }
     152         [ +  + ]:         10 :   if (d == 1)
     153                 :            :   { /* -f[2]/f[3], root of deg 1 polynomial */
     154                 :          5 :     r = Fl_mul(p - Fl_inv(f[3], p), f[2], p);
     155         [ +  - ]:          5 :     if (r >= s) y[++n] = r; /* otherwise double root */
     156                 :            :   }
     157                 :         10 :   avma = av; fixlg(y, n+1); return y;
     158                 :            : }
     159                 :            : static GEN
     160                 :       4092 : Flx_root_mod_2(GEN f)
     161                 :            : {
     162                 :       4092 :   int z1, z0 = !(f[2] & 1);
     163                 :            :   long i,n;
     164                 :            :   GEN y;
     165                 :            : 
     166         [ +  + ]:      15908 :   for (i=2, n=1; i < lg(f); i++) n += f[i];
     167                 :       4092 :   z1 = n & 1;
     168                 :       4092 :   y = cgetg(z0+z1+1, t_VECSMALL); i = 1;
     169         [ +  + ]:       4092 :   if (z0) y[i++] = 0;
     170         [ +  + ]:       4092 :   if (z1) y[i  ] = 1;
     171                 :       4092 :   return y;
     172                 :            : }
     173                 :            : static ulong
     174                 :         10 : Flx_oneroot_mod_2(GEN f)
     175                 :            : {
     176                 :            :   long i,n;
     177         [ -  + ]:         10 :   if (!(f[2] & 1)) return 0;
     178         [ +  + ]:         40 :   for (i=2, n=1; i < lg(f); i++) n += f[i];
     179         [ +  + ]:         10 :   if (n & 1) return 1;
     180                 :         10 :   return 2;
     181                 :            : }
     182                 :            : 
     183                 :            : static GEN FpX_roots_i(GEN f, GEN p);
     184                 :            : static GEN Flx_roots_i(GEN f, ulong p);
     185                 :            : static GEN FpX_Berlekamp_i(GEN f, GEN pp, long flag);
     186                 :            : 
     187                 :            : /* Generic driver to computes the roots of f modulo pp, using 'Roots' when
     188                 :            :  * pp is a small prime.
     189                 :            :  * if (gpwrap), check types thoroughly and return t_INTMODs, otherwise
     190                 :            :  * assume that f is an FpX, pp a prime and return t_INTs */
     191                 :            : static GEN
     192                 :      40575 : rootmod_aux(GEN f, GEN pp, GEN (*Roots)(GEN,ulong), int gpwrap)
     193                 :            : {
     194                 :      40575 :   pari_sp av = avma;
     195                 :            :   GEN y;
     196         [ +  + ]:      40575 :   if (gpwrap)
     197                 :         55 :     factmod_init(&f, pp);
     198                 :            :   else
     199                 :      40520 :     ZX_factmod_init(&f, pp);
     200      [ +  +  + ]:      40575 :   switch(lg(f))
     201                 :            :   {
     202                 :         10 :     case 2: pari_err_ROOTS0("rootmod");
     203                 :         25 :     case 3: avma = av; return cgetg(1,t_COL);
     204                 :            :   }
     205         [ +  + ]:      40540 :   if (typ(f) == t_VECSMALL)
     206                 :            :   {
     207                 :      39241 :     ulong p = pp[2];
     208         [ +  + ]:      39241 :     if (p == 2)
     209                 :       4092 :       y = Flx_root_mod_2(f);
     210                 :            :     else
     211                 :            :     {
     212         [ -  + ]:      35149 :       if (!odd(p)) pari_err_PRIME("rootmod",utoi(p));
     213                 :      35149 :       y = Roots(f, p);
     214                 :            :     }
     215                 :      39241 :     y = Flc_to_ZC(y);
     216                 :            :   }
     217                 :            :   else
     218                 :       1299 :     y = FpX_roots_i(f, pp);
     219         [ +  + ]:      40540 :   if (gpwrap) y = FpC_to_mod(y, pp);
     220                 :      40565 :   return gerepileupto(av, y);
     221                 :            : }
     222                 :            : /* assume that f is a ZX an pp a prime */
     223                 :            : GEN
     224                 :      40520 : FpX_roots(GEN f, GEN pp)
     225                 :      40520 : { return rootmod_aux(f, pp, Flx_roots_i, 0); }
     226                 :            : /* no assumptions on f and pp */
     227                 :            : GEN
     228                 :         10 : rootmod2(GEN f, GEN pp) { return rootmod_aux(f, pp, &Flx_roots_naive, 1); }
     229                 :            : GEN
     230                 :         45 : rootmod(GEN f, GEN pp) { return rootmod_aux(f, pp, &Flx_roots_i, 1); }
     231                 :            : GEN
     232                 :         50 : rootmod0(GEN f, GEN p, long flag)
     233                 :            : {
     234      [ +  +  - ]:         50 :   switch(flag)
     235                 :            :   {
     236                 :         40 :     case 0: return rootmod(f,p);
     237                 :         10 :     case 1: return rootmod2(f,p);
     238                 :          0 :     default: pari_err_FLAG("polrootsmod");
     239                 :            :   }
     240                 :         40 :   return NULL; /* not reached */
     241                 :            : }
     242                 :            : 
     243                 :            : /* assume x reduced mod p > 2, monic. */
     244                 :            : static int
     245                 :         15 : FpX_quad_factortype(GEN x, GEN p)
     246                 :            : {
     247                 :         15 :   GEN b = gel(x,3), c = gel(x,2);
     248                 :         15 :   GEN D = subii(sqri(b), shifti(c,2));
     249                 :         15 :   return kronecker(D,p);
     250                 :            : }
     251                 :            : /* assume x reduced mod p, monic. Return one root, or NULL if irreducible */
     252                 :            : GEN
     253                 :       2885 : FpX_quad_root(GEN x, GEN p, int unknown)
     254                 :            : {
     255                 :       2885 :   GEN s, D, b = gel(x,3), c = gel(x,2);
     256                 :            : 
     257         [ -  + ]:       2885 :   if (equaliu(p, 2)) {
     258         [ #  # ]:          0 :     if (!signe(b)) return c;
     259         [ #  # ]:          0 :     return signe(c)? NULL: gen_1;
     260                 :            :   }
     261                 :       2885 :   D = subii(sqri(b), shifti(c,2));
     262                 :       2885 :   D = remii(D,p);
     263 [ +  + ][ +  + ]:       2885 :   if (unknown && kronecker(D,p) == -1) return NULL;
     264                 :            : 
     265                 :       2879 :   s = Fp_sqrt(D,p);
     266                 :            :   /* p is not prime, go on and give e.g. maxord a chance to recover */
     267         [ +  + ]:       2879 :   if (!s) return NULL;
     268                 :       2885 :   return Fp_halve(subii(s,b), p);
     269                 :            : }
     270                 :            : static GEN
     271                 :       1356 : FpX_otherroot(GEN x, GEN r, GEN p)
     272                 :            : {
     273                 :       1356 :   GEN s = addii(gel(x,3), r);
     274         [ -  + ]:       1356 :   if (!signe(s)) return s;
     275         [ +  + ]:       1356 :   s = subii(p, s); if (signe(s) < 0) s = addii(s,p);
     276                 :       1356 :   return s;
     277                 :            : }
     278                 :            : 
     279                 :            : /* disc(x^2+bx+c) = b^2 - 4c */
     280                 :            : static ulong
     281                 :     164218 : Fl_disc_bc(ulong b, ulong c, ulong p)
     282                 :     164218 : { return Fl_sub(Fl_sqr(b,p), Fl_double(Fl_double(c,p),p), p); }
     283                 :            : /* p > 2 */
     284                 :            : static ulong
     285                 :     143848 : Flx_quad_root(GEN x, ulong p, int unknown)
     286                 :            : {
     287                 :     143848 :   ulong s, b = x[3], c = x[2];
     288                 :     143848 :   ulong D = Fl_disc_bc(b, c, p);
     289 [ +  + ][ +  + ]:     143848 :   if (unknown && krouu(D,p) == -1) return p;
     290                 :     108380 :   s = Fl_sqrt(D,p);
     291         [ +  + ]:     108380 :   if (s==~0UL) return p;
     292                 :     143848 :   return Fl_halve(Fl_sub(s,b, p), p);
     293                 :            : }
     294                 :            : static ulong
     295                 :     107983 : Flx_otherroot(GEN x, ulong r, ulong p)
     296                 :     107983 : { return Fl_neg(Fl_add(x[3], r, p), p); }
     297                 :            : 
     298                 :            : 
     299                 :            : /* 'todo' contains the list of factors to be split.
     300                 :            :  * 'done' the list of finished factors, no longer touched */
     301                 :            : struct split_t { GEN todo, done; };
     302                 :            : static void
     303                 :      13534 : split_init(struct split_t *S, long max)
     304                 :            : {
     305                 :      13534 :   S->done = vectrunc_init(max);
     306                 :      13534 :   S->todo = vectrunc_init(max);
     307                 :      13534 : }
     308                 :            : #if 0
     309                 :            : /* move todo[i] to done */
     310                 :            : static void
     311                 :            : split_convert(struct split_t *S, long i)
     312                 :            : {
     313                 :            :   long n = lg(S->todo)-1;
     314                 :            :   vectrunc_append(S->done, gel(S->todo,i));
     315                 :            :   if (n) gel(S->todo,i) = gel(S->todo, n);
     316                 :            :   setlg(S->todo, n);
     317                 :            : }
     318                 :            : #endif
     319                 :            : /* append t to todo */
     320                 :            : static void
     321                 :       7929 : split_add(struct split_t *S, GEN t) { vectrunc_append(S->todo, t); }
     322                 :            : /* delete todo[i], add t to done */
     323                 :            : static void
     324                 :       7929 : split_moveto_done(struct split_t *S, long i, GEN t)
     325                 :            : {
     326                 :       7929 :   long n = lg(S->todo)-1;
     327                 :       7929 :   vectrunc_append(S->done, t);
     328         [ +  - ]:       7929 :   if (n) gel(S->todo,i) = gel(S->todo, n);
     329                 :       7929 :   setlg(S->todo, n);
     330                 :            : 
     331                 :       7929 : }
     332                 :            : /* append t to done */
     333                 :            : static void
     334                 :      17629 : split_add_done(struct split_t *S, GEN t)
     335                 :      17629 : { vectrunc_append(S->done, t); }
     336                 :            : /* split todo[i] into a and b */
     337                 :            : static void
     338                 :       6164 : split_todo(struct split_t *S, long i, GEN a, GEN b)
     339                 :            : {
     340                 :       6164 :   gel(S->todo, i) = a;
     341                 :       6164 :   split_add(S, b);
     342                 :       6164 : }
     343                 :            : /* split todo[i] into a and b, moved to done */
     344                 :            : static void
     345                 :       5127 : split_done(struct split_t *S, long i, GEN a, GEN b)
     346                 :            : {
     347                 :       5127 :   split_moveto_done(S, i, a);
     348                 :       5127 :   split_add_done(S, b);
     349                 :       5127 : }
     350                 :            : 
     351                 :            : /* by splitting, assume p > 2 prime, deg(f) > 0, and f monic */
     352                 :            : static GEN
     353                 :       1299 : FpX_roots_i(GEN f, GEN p)
     354                 :            : {
     355                 :            :   GEN pol, pol0, a, q;
     356                 :            :   struct split_t S;
     357                 :            : 
     358                 :       1299 :   split_init(&S, lg(f)-1);
     359                 :       1299 :   settyp(S.done, t_COL);
     360         [ +  + ]:       1299 :   if (ZX_valrem(f, &f)) split_add_done(&S, gen_0);
     361   [ +  +  +  + ]:       1299 :   switch(degpol(f))
     362                 :            :   {
     363                 :          5 :     case 0: return S.done;
     364                 :          5 :     case 1: split_add_done(&S, subii(p, gel(f,2))); return S.done;
     365                 :            :     case 2: {
     366                 :        684 :       GEN s, r = FpX_quad_root(f, p, 1);
     367         [ +  - ]:        684 :       if (r) {
     368                 :        684 :         split_add_done(&S, r);
     369                 :        684 :         s = FpX_otherroot(f,r, p);
     370                 :            :         /* f not known to be square free yet */
     371         [ +  - ]:        684 :         if (!equalii(r, s)) split_add_done(&S, s);
     372                 :            :       }
     373                 :        684 :       return sort(S.done);
     374                 :            :     }
     375                 :            :   }
     376                 :            : 
     377                 :        605 :   a = FpXQ_pow(pol_x(varn(f)), subiu(p,1), f,p);
     378         [ -  + ]:        605 :   if (lg(a) < 3) pari_err_PRIME("rootmod",p);
     379                 :        605 :   a = FpX_Fp_sub_shallow(a, gen_1, p); /* a = x^(p-1) - 1 mod f */
     380                 :        605 :   a = FpX_gcd(f,a, p);
     381         [ -  + ]:        605 :   if (!degpol(a)) return S.done;
     382                 :        605 :   split_add(&S, FpX_normalize(a,p));
     383                 :            : 
     384                 :        605 :   q = shifti(p,-1);
     385                 :        605 :   pol0 = icopy(gen_1); /* constant term, will vary in place */
     386                 :        605 :   pol = deg1pol_shallow(gen_1, pol0, varn(f));
     387                 :        605 :   for (pol0[2] = 1;; pol0[2]++)
     388                 :            :   {
     389                 :       1280 :     long j, l = lg(S.todo);
     390         [ +  + ]:       1280 :     if (l == 1) return sort(S.done);
     391 [ -  + ][ #  # ]:        675 :     if (pol0[2] == 100 && !BPSW_psp(p)) pari_err_PRIME("polrootsmod",p);
     392         [ +  + ]:       1415 :     for (j = 1; j < l; j++)
     393                 :            :     {
     394                 :        740 :       GEN c = gel(S.todo,j);
     395      [ +  +  + ]:        740 :       switch(degpol(c))
     396                 :            :       { /* convert linear and quadratics to roots, try to split the rest */
     397                 :            :         case 1:
     398                 :         40 :           split_moveto_done(&S, j, subii(p, gel(c,2)));
     399                 :         40 :           j--; l--; break;
     400                 :            :         case 2: {
     401                 :        630 :           GEN r = FpX_quad_root(c, p, 0), s = FpX_otherroot(c,r, p);
     402                 :        630 :           split_done(&S, j, r, s);
     403                 :        630 :           j--; l--; break;
     404                 :            :         }
     405                 :            :         default: {
     406                 :            :           /* b = pol^(p-1)/2 - 1 */
     407                 :         70 :           GEN b = FpX_Fp_sub_shallow(FpXQ_pow(pol,q, c,p), gen_1, p);
     408                 :            :           long db;
     409                 :         70 :           b = FpX_gcd(c,b, p); db = degpol(b);
     410 [ +  - ][ +  + ]:         70 :           if (db && db < degpol(c))
     411                 :            :           {
     412                 :         65 :             b = FpX_normalize(b, p);
     413                 :         65 :             c = FpX_div(c,b, p);
     414                 :         65 :             split_todo(&S, j, b, c);
     415                 :            :           }
     416                 :            :         }
     417                 :            :       }
     418                 :            :     }
     419                 :       1974 :   }
     420                 :            : }
     421                 :            : 
     422                 :            : /* Assume f is normalized */
     423                 :            : static ulong
     424                 :         10 : Flx_cubic_root(GEN ff, ulong p)
     425                 :            : {
     426                 :         10 :   GEN f = Flx_normalize(ff,p);
     427                 :         10 :   ulong pi = get_Fl_red(p);
     428         [ -  + ]:         10 :   ulong a = f[4], b=f[3], c=f[2], p3 = p%3==1 ? (2*p+1)/3 :(p+1)/3;
     429                 :         10 :   ulong t = Fl_mul_pre(a, p3, p, pi), t2 = Fl_sqr_pre(t, p, pi);
     430                 :         10 :   ulong A = Fl_sub(b, Fl_triple(t2, p), p);
     431                 :         10 :   ulong B = Fl_add(Fl_mul_pre(t, Fl_sub(Fl_double(t2, p), b, p), p ,pi), c, p);
     432                 :         10 :   ulong A3 =  Fl_mul_pre(A, p3, p, pi);
     433                 :         10 :   ulong A32 = Fl_sqr_pre(A3, p, pi), A33 = Fl_mul_pre(A3, A32, p, pi);
     434                 :         10 :   ulong S = Fl_neg(B,p), P = Fl_neg(A3,p);
     435                 :         10 :   ulong D = Fl_add(Fl_sqr_pre(S, p, pi), Fl_double(Fl_double(A33, p), p), p);
     436                 :         10 :   ulong s = Fl_sqrt_pre(D, p, pi), vS1, vS2;
     437         [ -  + ]:         10 :   if (s!=~0UL)
     438                 :            :   {
     439         [ #  # ]:          0 :     ulong S1 = S==s ? S: Fl_halve(Fl_sub(S, s, p), p);
     440         [ #  # ]:          0 :     if (p%3==2) /* 1 solutions */
     441                 :          0 :       vS1 = Fl_powu_pre(S1, (2*p-1)/3, p, pi);
     442                 :            :     else
     443                 :            :     {
     444                 :          0 :       vS1 = Fl_sqrtl_pre(S1, 3, p, pi);
     445         [ #  # ]:          0 :       if (vS1==~0UL) return p; /*0 solutions*/
     446                 :            :       /*3 solutions*/
     447                 :            :     }
     448         [ #  # ]:          0 :     vS2 = P? Fl_mul_pre(P, Fl_inv(vS1, p), p, pi): 0;
     449                 :          0 :     return Fl_sub(Fl_add(vS1,vS2, p), t, p);
     450                 :            :   }
     451                 :            :   else
     452                 :            :   {
     453                 :         10 :     pari_sp av = avma;
     454                 :         10 :     GEN S1 = mkvecsmall2(Fl_halve(S, p), Fl_halve(1UL, p));
     455                 :         10 :     GEN vS1 = Fl2_sqrtn_pre(S1, utoi(3), D, p, pi, NULL);
     456                 :            :     ulong Sa;
     457         [ -  + ]:         10 :     if (!vS1) return p; /*0 solutions, p%3==2*/
     458                 :         10 :     Sa = vS1[1];
     459         [ -  + ]:         10 :     if (p%3==1) /*1 solutions*/
     460                 :            :     {
     461                 :          0 :       ulong Fa = Fl2_norm_pre(vS1, D, p, pi);
     462         [ #  # ]:          0 :       if (Fa!=P)
     463                 :          0 :         Sa = Fl_mul(Sa, Fl_div(Fa, P, p),p);
     464                 :            :     }
     465                 :         10 :     avma = av;
     466                 :         10 :     return Fl_sub(Fl_double(Sa,p),t,p);
     467                 :            :   }
     468                 :            : }
     469                 :            : 
     470                 :            : /* assume p > 2 prime */
     471                 :            : static ulong
     472                 :        988 : Flx_oneroot_i(GEN f, ulong p)
     473                 :            : {
     474                 :            :   GEN pol, a;
     475                 :            :   ulong q;
     476                 :            :   long da;
     477                 :            : 
     478         [ +  + ]:        988 :   if (Flx_val(f)) return 0;
     479   [ +  +  -  + ]:        983 :   switch(degpol(f))
     480                 :            :   {
     481                 :        108 :     case 1: return Fl_neg(f[2], p);
     482                 :        425 :     case 2: return Flx_quad_root(f, p, 1);
     483         [ #  # ]:          0 :     case 3: if (p>3) return Flx_cubic_root(f, p); /*FALL THROUGH*/
     484                 :            :   }
     485                 :            : 
     486                 :        450 :   a = Flxq_powu(polx_Flx(f[1]), p - 1, f,p);
     487         [ -  + ]:        450 :   if (lg(a) < 3) pari_err_PRIME("rootmod",utoipos(p));
     488                 :        450 :   a = Flx_Fl_add(a, p-1, p); /* a = x^(p-1) - 1 mod f */
     489                 :        450 :   a = Flx_gcd(f,a, p);
     490                 :        450 :   da = degpol(a);
     491         [ +  + ]:        450 :   if (!da) return p;
     492                 :         10 :   a = Flx_normalize(a,p);
     493                 :            : 
     494                 :         10 :   q = p >> 1;
     495                 :         10 :   pol = polx_Flx(f[1]);
     496                 :         10 :   for(pol[2] = 1;; pol[2]++)
     497                 :            :   {
     498 [ -  + ][ #  # ]:         20 :     if (pol[2] == 1000 && !uisprime(p)) pari_err_PRIME("Flx_oneroot",utoipos(p));
     499   [ -  -  +  + ]:         20 :     switch(da)
     500                 :            :     {
     501                 :          0 :       case 1: return Fl_neg(a[2], p);
     502                 :          0 :       case 2: return Flx_quad_root(a, p, 0);
     503         [ +  - ]:         10 :       case 3: if (p>3) return Flx_cubic_root(a, p); /*FALL THROUGH*/
     504                 :            :       default: {
     505                 :         10 :         GEN b = Flx_Fl_add(Flxq_powu(pol,q, a,p), p-1, p);
     506                 :            :         long db;
     507                 :         10 :         b = Flx_gcd(a,b, p); db = degpol(b);
     508 [ +  - ][ +  - ]:         10 :         if (db && db < da)
     509                 :            :         {
     510                 :         10 :           b = Flx_normalize(b, p);
     511         [ +  - ]:         10 :           if (db <= (da >> 1)) {
     512                 :         10 :             a = b;
     513                 :         10 :             da = db;
     514                 :            :           } else {
     515                 :          0 :             a = Flx_div(a,b, p);
     516                 :          0 :             da -= db;
     517                 :            :           }
     518                 :            :         }
     519                 :            :       }
     520                 :            :     }
     521                 :        998 :   }
     522                 :            : }
     523                 :            : /* assume p > 2 prime */
     524                 :            : static GEN
     525                 :       1522 : FpX_oneroot_i(GEN f, GEN p)
     526                 :            : {
     527                 :            :   GEN pol, pol0, a, q;
     528                 :            :   long da;
     529                 :            : 
     530         [ -  + ]:       1522 :   if (ZX_val(f)) return gen_0;
     531      [ +  +  + ]:       1522 :   switch(degpol(f))
     532                 :            :   {
     533                 :         82 :     case 1: return subii(p, gel(f,2));
     534                 :       1425 :     case 2: return FpX_quad_root(f, p, 1);
     535                 :            :   }
     536                 :            : 
     537                 :         15 :   a = FpXQ_pow(pol_x(varn(f)), subiu(p,1), f,p);
     538         [ -  + ]:         15 :   if (lg(a) < 3) pari_err_PRIME("rootmod",p);
     539                 :         15 :   a = FpX_Fp_sub_shallow(a, gen_1, p); /* a = x^(p-1) - 1 mod f */
     540                 :         15 :   a = FpX_gcd(f,a, p);
     541                 :         15 :   da = degpol(a);
     542         [ -  + ]:         15 :   if (!da) return NULL;
     543                 :         15 :   a = FpX_normalize(a,p);
     544                 :            : 
     545                 :         15 :   q = shifti(p,-1);
     546                 :         15 :   pol0 = icopy(gen_1); /* constant term, will vary in place */
     547                 :         15 :   pol = deg1pol_shallow(gen_1, pol0, varn(f));
     548                 :         15 :   for (pol0[2]=1; ; pol0[2]++)
     549                 :            :   {
     550 [ -  + ][ #  # ]:         50 :     if (pol0[2] == 1000 && !BPSW_psp(p)) pari_err_PRIME("FpX_oneroot",p);
     551      [ +  +  + ]:         50 :     switch(da)
     552                 :            :     {
     553                 :         10 :       case 1: return subii(p, gel(a,2));
     554                 :          5 :       case 2: return FpX_quad_root(a, p, 0);
     555                 :            :       default: {
     556                 :         35 :         GEN b = FpX_Fp_sub_shallow(FpXQ_pow(pol,q, a,p), gen_1, p);
     557                 :            :         long db;
     558                 :         35 :         b = FpX_gcd(a,b, p); db = degpol(b);
     559 [ +  - ][ +  - ]:         35 :         if (db && db < da)
     560                 :            :         {
     561                 :         35 :           b = FpX_normalize(b, p);
     562         [ +  + ]:         35 :           if (db <= (da >> 1)) {
     563                 :         25 :             a = b;
     564                 :         25 :             da = db;
     565                 :            :           } else {
     566                 :         10 :             a = FpX_div(a,b, p);
     567                 :         10 :             da -= db;
     568                 :            :           }
     569                 :            :         }
     570                 :            :       }
     571                 :            :     }
     572                 :       1557 :   }
     573                 :            : }
     574                 :            : 
     575                 :            : ulong
     576                 :        450 : Flx_oneroot(GEN f, ulong p)
     577                 :            : {
     578                 :        450 :   pari_sp av = avma;
     579                 :            :   ulong r;
     580      [ -  -  + ]:        450 :   switch(lg(f))
     581                 :            :   {
     582                 :          0 :     case 2: return 0;
     583                 :          0 :     case 3: avma = av; return p;
     584                 :            :   }
     585         [ -  + ]:        450 :   if (p == 2) return Flx_oneroot_mod_2(f);
     586                 :        450 :   r = Flx_oneroot_i(Flx_normalize(f, p), p);
     587                 :        450 :   avma = av; return r;
     588                 :            : }
     589                 :            : 
     590                 :            : /* assume that p is prime */
     591                 :            : GEN
     592                 :       2070 : FpX_oneroot(GEN f, GEN pp) {
     593                 :       2070 :   pari_sp av = avma;
     594                 :       2070 :   ZX_factmod_init(&f, pp);
     595      [ -  -  + ]:       2070 :   switch(lg(f))
     596                 :            :   {
     597                 :          0 :     case 2: avma = av; return gen_0;
     598                 :          0 :     case 3: avma = av; return NULL;
     599                 :            :   }
     600         [ +  + ]:       2070 :   if (typ(f) == t_VECSMALL)
     601                 :            :   {
     602                 :        548 :     ulong r, p = pp[2];
     603         [ +  + ]:        548 :     if (p == 2)
     604                 :         10 :       r = Flx_oneroot_mod_2(f);
     605                 :            :     else
     606                 :        538 :       r = Flx_oneroot_i(f, p);
     607                 :        548 :     avma = av;
     608         [ +  + ]:        548 :     return (r == p)? NULL: utoi(r);
     609                 :            :   }
     610                 :       1522 :   f = FpX_oneroot_i(f, pp);
     611         [ -  + ]:       1522 :   if (!f) { avma = av; return NULL; }
     612                 :       2070 :   return gerepileuptoint(av, f);
     613                 :            : }
     614                 :            : 
     615                 :            : /*******************************************************************/
     616                 :            : /*                                                                 */
     617                 :            : /*                     FACTORISATION MODULO p                      */
     618                 :            : /*                                                                 */
     619                 :            : /*******************************************************************/
     620                 :            : 
     621                 :            : /* Functions giving information on the factorisation. */
     622                 :            : 
     623                 :            : /* u in Z[X], return kernel of (Frob - Id) over Fp[X] / u */
     624                 :            : GEN
     625                 :         58 : FpX_Berlekamp_ker(GEN u, GEN p)
     626                 :            : {
     627                 :         58 :   pari_sp ltop=avma;
     628                 :         58 :   long j,N = degpol(u);
     629                 :         58 :   GEN XP = FpXQ_pow(pol_x(varn(u)),p,u,p);
     630                 :         58 :   GEN Q  = FpXQ_matrix_pow(XP,N,N,u,p);
     631         [ +  + ]:       1600 :   for (j=1; j<=N; j++)
     632                 :       1542 :     gcoeff(Q,j,j) = Fp_sub(gcoeff(Q,j,j), gen_1, p);
     633                 :         58 :   return gerepileupto(ltop, FpM_ker(Q,p));
     634                 :            : }
     635                 :            : 
     636                 :            : GEN
     637                 :      14710 : F2x_Berlekamp_ker(GEN u)
     638                 :            : {
     639                 :      14710 :   pari_sp ltop=avma;
     640                 :      14710 :   long j,N = F2x_degree(u);
     641                 :            :   GEN Q, XP;
     642                 :            :   pari_timer T;
     643                 :      14710 :   timer_start(&T);
     644                 :      14710 :   XP = F2xq_sqr(polx_F2x(u[1]),u);
     645                 :      14710 :   Q  = F2xq_matrix_pow(XP,N,N,u);
     646         [ +  + ]:      76965 :   for (j=1; j<=N; j++)
     647                 :      62255 :     F2m_flip(Q,j,j);
     648         [ -  + ]:      14710 :   if(DEBUGLEVEL>=9) timer_printf(&T,"Berlekamp matrix");
     649                 :      14710 :   Q = F2m_ker_sp(Q,0);
     650         [ -  + ]:      14710 :   if(DEBUGLEVEL>=9) timer_printf(&T,"kernel");
     651                 :      14710 :   return gerepileupto(ltop,Q);
     652                 :            : }
     653                 :            : 
     654                 :            : GEN
     655                 :     147398 : Flx_Berlekamp_ker(GEN u, ulong l)
     656                 :            : {
     657                 :     147398 :   pari_sp ltop=avma;
     658                 :     147398 :   long j,N = degpol(u);
     659                 :            :   GEN Q, XP;
     660                 :            :   pari_timer T;
     661                 :     147398 :   timer_start(&T);
     662                 :     147398 :   XP = Flxq_powu(polx_Flx(u[1]),l,u,l);
     663                 :     147398 :   Q  = Flxq_matrix_pow(XP,N,N,u,l);
     664         [ +  + ]:     661328 :   for (j=1; j<=N; j++)
     665                 :     513930 :     coeff(Q,j,j) = Fl_sub(coeff(Q,j,j),1,l);
     666         [ -  + ]:     147398 :   if(DEBUGLEVEL>=9) timer_printf(&T,"Berlekamp matrix");
     667                 :     147398 :   Q = Flm_ker_sp(Q,l,0);
     668         [ -  + ]:     147398 :   if(DEBUGLEVEL>=9) timer_printf(&T,"kernel");
     669                 :     147398 :   return gerepileupto(ltop,Q);
     670                 :            : }
     671                 :            : 
     672                 :            : /* product of terms of degree 1 in factorization of f */
     673                 :            : GEN
     674                 :      41257 : FpX_split_part(GEN f, GEN p)
     675                 :            : {
     676                 :      41257 :   long n = degpol(f);
     677                 :      41257 :   GEN z, X = pol_x(varn(f));
     678         [ +  + ]:      41257 :   if (n <= 1) return f;
     679                 :      41247 :   f = FpX_red(f, p);
     680                 :      41247 :   z = FpXQ_pow(X, p, f, p);
     681                 :      41247 :   z = FpX_sub(z, X, p);
     682                 :      41257 :   return FpX_gcd(z,f,p);
     683                 :            : }
     684                 :            : 
     685                 :            : /* Compute the number of roots in Fp without counting multiplicity
     686                 :            :  * return -1 for 0 polynomial. lc(f) must be prime to p. */
     687                 :            : long
     688                 :      12647 : FpX_nbroots(GEN f, GEN p)
     689                 :            : {
     690                 :      12647 :   pari_sp av = avma;
     691                 :      12647 :   GEN z = FpX_split_part(f, p);
     692                 :      12647 :   avma = av; return degpol(z);
     693                 :            : }
     694                 :            : 
     695                 :            : int
     696                 :          0 : FpX_is_totally_split(GEN f, GEN p)
     697                 :            : {
     698                 :          0 :   long n=degpol(f);
     699                 :          0 :   pari_sp av = avma;
     700                 :            :   GEN z;
     701         [ #  # ]:          0 :   if (n <= 1) return 1;
     702         [ #  # ]:          0 :   if (cmpui(n, p) > 0) return 0;
     703                 :          0 :   f = FpX_red(f, p);
     704                 :          0 :   z = FpXQ_pow(pol_x(varn(f)), p, f, p);
     705                 :          0 :   avma = av; return gequalX(z);
     706                 :            : }
     707                 :            : 
     708                 :            : /* Flv_Flx( Flm_Flc_mul(x, Flx_Flv(y), p) ) */
     709                 :            : static GEN
     710                 :     153637 : Flm_Flx_mul(GEN x, GEN y, ulong p)
     711                 :            : {
     712                 :     153637 :   long i,k,l, ly = lg(y)-1;
     713                 :            :   GEN z;
     714                 :     153637 :   long vs=y[1];
     715         [ -  + ]:     153637 :   if (ly==1) return zero_Flx(vs);
     716                 :     153637 :   l = lgcols(x);
     717                 :     153637 :   y++;
     718                 :     153637 :   z = zero_zv(l) + 1;
     719         [ +  - ]:     153637 :   if (SMALL_ULONG(p))
     720                 :            :   {
     721         [ +  + ]:    3103385 :     for (k=1; k<ly; k++)
     722                 :            :     {
     723                 :            :       GEN c;
     724         [ +  + ]:    2949748 :       if (!y[k]) continue;
     725                 :    1251988 :       c = gel(x,k);
     726         [ +  + ]:    1251988 :       if (y[k] == 1)
     727         [ +  + ]:    3164204 :         for (i=1; i<l; i++)
     728                 :            :         {
     729                 :    3047424 :           z[i] += c[i];
     730         [ -  + ]:    3047424 :           if (z[i] & HIGHBIT) z[i] %= p;
     731                 :            :         }
     732                 :            :       else
     733         [ +  + ]:   58131682 :         for (i=1; i<l; i++)
     734                 :            :         {
     735                 :   56996474 :           z[i] += c[i] * y[k];
     736         [ -  + ]:   56996474 :           if (z[i] & HIGHBIT) z[i] %= p;
     737                 :            :         }
     738                 :            :     }
     739         [ +  + ]:    3668923 :     for (i=1; i<l; i++) z[i] %= p;
     740                 :            :   }
     741                 :            :   else
     742                 :            :   {
     743         [ #  # ]:          0 :     for (k=1; k<ly; k++)
     744                 :            :     {
     745                 :            :       GEN c;
     746         [ #  # ]:          0 :       if (!y[k]) continue;
     747                 :          0 :       c = gel(x,k);
     748         [ #  # ]:          0 :       if (y[k] == 1)
     749         [ #  # ]:          0 :         for (i=1; i<l; i++)
     750                 :          0 :           z[i] = Fl_add(z[i], c[i], p);
     751                 :            :       else
     752         [ #  # ]:          0 :         for (i=1; i<l; i++)
     753                 :          0 :           z[i] = Fl_add(z[i], Fl_mul(c[i],y[k],p), p);
     754                 :            :     }
     755                 :            :   }
     756 [ +  - ][ +  + ]:     973705 :   while (--l && !z[l]);
     757         [ -  + ]:     153637 :   if (!l) return zero_Flx(vs);
     758                 :     153637 :   *z-- = vs; return z;
     759                 :            : }
     760                 :            : 
     761                 :            : /* z must be squarefree mod p*/
     762                 :            : GEN
     763                 :      54932 : Flx_nbfact_by_degree(GEN z, long *nb, ulong p)
     764                 :            : {
     765                 :      54932 :   long lgg, d = 0, e = degpol(z);
     766                 :      54932 :   GEN D = zero_zv(e);
     767                 :      54932 :   pari_sp av = avma;
     768                 :      54932 :   GEN g, w,  PolX = polx_Flx(z[1]);
     769                 :      54932 :   GEN XP = Flxq_powu(PolX,p,z,p);
     770                 :      54932 :   GEN MP = Flxq_matrix_pow(XP,e,e,z,p);
     771                 :            : 
     772                 :      54932 :   w = PolX; *nb = 0;
     773         [ +  + ]:     172581 :   while (d < (e>>1))
     774                 :            :   { /* here e = degpol(z) */
     775                 :     153637 :     d++;
     776                 :     153637 :     w = Flm_Flx_mul(MP, w, p); /* w^p mod (z,p) */
     777                 :     153637 :     g = Flx_gcd(z, Flx_sub(w, PolX, p), p);
     778         [ +  + ]:     153637 :     lgg = degpol(g); if (!lgg) continue;
     779                 :            : 
     780                 :      49518 :     e -= lgg;
     781                 :      49518 :     D[d] = lgg/d; *nb += D[d];
     782         [ -  + ]:      49518 :     if (DEBUGLEVEL>5) err_printf("   %3ld fact. of degree %3ld\n", D[d], d);
     783         [ +  + ]:      49518 :     if (!e) break;
     784                 :      13530 :     z = Flx_div(z, g, p);
     785                 :      13530 :     w = Flx_rem(w, z, p);
     786                 :            :   }
     787         [ +  + ]:      54932 :   if (e)
     788                 :            :   {
     789         [ -  + ]:      18944 :     if (DEBUGLEVEL>5) err_printf("   %3ld fact. of degree %3ld\n",1,e);
     790                 :      18944 :     D[e] = 1; (*nb)++;
     791                 :            :   }
     792                 :      54932 :   avma = av; return D;
     793                 :            : }
     794                 :            : 
     795                 :            : /* z must be squarefree mod p*/
     796                 :            : long
     797                 :      25092 : Flx_nbfact(GEN z, ulong p)
     798                 :            : {
     799                 :      25092 :   pari_sp av = avma;
     800                 :      25092 :   long nb; (void)Flx_nbfact_by_degree(z, &nb, p);
     801                 :      25092 :   avma = av; return nb;
     802                 :            : }
     803                 :            : 
     804                 :            : long
     805                 :    1570688 : Flx_nbroots(GEN f, ulong p)
     806                 :            : {
     807                 :    1570688 :   long n = degpol(f);
     808                 :    1570688 :   pari_sp av = avma;
     809                 :            :   GEN z, X;
     810         [ +  + ]:    1570688 :   if (n <= 1) return n;
     811                 :    1570633 :   X = polx_Flx(f[1]);
     812                 :    1570633 :   z = Flxq_powu(X, p, f, p);
     813                 :    1570633 :   z = Flx_sub(z, X, p);
     814                 :    1570633 :   z = Flx_gcd(z, f, p);
     815                 :    1570688 :   avma = av; return degpol(z);
     816                 :            : }
     817                 :            : 
     818                 :            : long
     819                 :         20 : FpX_nbfact(GEN u, GEN p)
     820                 :            : {
     821                 :         20 :   pari_sp av = avma;
     822                 :         20 :   GEN vker = FpX_Berlekamp_ker(u, p);
     823                 :         20 :   avma = av; return lg(vker)-1;
     824                 :            : }
     825                 :            : 
     826                 :            : static GEN
     827                 :        202 : try_pow(GEN w0, GEN pol, GEN p, GEN q, long r)
     828                 :            : {
     829                 :        202 :   GEN w2, w = FpXQ_pow(w0,q, pol,p);
     830                 :            :   long s;
     831         [ +  + ]:        202 :   if (gequal1(w)) return w0;
     832         [ +  + ]:        438 :   for (s=1; s<r; s++,w=w2)
     833                 :            :   {
     834                 :        253 :     w2 = gsqr(w);
     835                 :        253 :     w2 = FpX_rem(w2, pol, p);
     836         [ +  + ]:        253 :     if (gequal1(w2)) break;
     837                 :            :   }
     838         [ +  - ]:        202 :   return gequalm1(w)? NULL: w;
     839                 :            : }
     840                 :            : 
     841                 :            : static GEN
     842                 :      39402 : Flx_try_pow(GEN w0, GEN pol, ulong p, GEN q, long r)
     843                 :            : {
     844                 :      39402 :   GEN w2, w = Flxq_pow(w0,q, pol,p);
     845                 :            :   long s;
     846         [ +  + ]:      39402 :   if (Flx_equal1(w)) return w0;
     847         [ +  + ]:     165151 :   for (s=1; s<r; s++,w=w2)
     848                 :            :   {
     849                 :     133045 :     w2 = Flxq_sqr(w,pol,p);
     850         [ +  + ]:     133045 :     if (Flx_equal1(w2)) break;
     851                 :            :   }
     852 [ +  + ][ -  + ]:      39402 :   return degpol(w)==0 && uel(w,2) == p-1 ? NULL: w;
     853                 :            : }
     854                 :            : 
     855                 :            : /* INPUT:
     856                 :            :  *  m integer (converted to polynomial w in Z[X] by stoFpX)
     857                 :            :  *  p prime; q = (p^d-1) / 2^r
     858                 :            :  *  t[0] polynomial of degree k*d product of k irreducible factors of degree d
     859                 :            :  *  t[0] is expected to be normalized (leading coeff = 1)
     860                 :            :  * OUTPUT:
     861                 :            :  *  t[0],t[1]...t[k-1] the k factors, normalized */
     862                 :            : 
     863                 :            : static void
     864                 :      14073 : F2x_split(ulong m, GEN *t, long d)
     865                 :            : {
     866                 :            :   long l, v, dv;
     867                 :            :   pari_sp av0, av;
     868                 :            :   GEN w;
     869                 :            : 
     870         [ +  + ]:      15583 :   dv = F2x_degree(*t); if (dv==d) return;
     871                 :       1510 :   v=(*t)[1]; av0=avma;
     872                 :       1510 :   for(av=avma;;avma=av)
     873                 :            :   {
     874                 :       2612 :     GEN w0 = w = F2xq_powu(polx_F2x(v), m-1, *t); m += 2;
     875         [ +  + ]:      17935 :     for (l=1; l<d; l++) w = F2x_add(w0, F2xq_sqr(w, *t));
     876                 :       2612 :     w = F2x_gcd(*t,w);
     877 [ +  + ][ +  + ]:       2612 :     l = F2x_degree(w); if (l && l!=dv) break;
     878                 :       1102 :   }
     879                 :       1510 :   w = gerepileupto(av0, w);
     880                 :       1510 :   l /= d; t[l]=F2x_div(*t,w); *t=w;
     881                 :       1510 :   F2x_split(m,t+l,d);
     882                 :       1510 :   F2x_split(m,t,  d);
     883                 :            : }
     884                 :            : 
     885                 :            : static void
     886                 :      74434 : Flx_split(GEN *t, long d, ulong p, GEN q, long r)
     887                 :            : {
     888                 :            :   long l, v, dv;
     889                 :            :   pari_sp av0, av;
     890                 :            :   GEN w;
     891                 :            : 
     892         [ +  + ]:     102571 :   dv=degpol(*t); if (dv==d) return;
     893                 :      28137 :   v=(*t)[1]; av0=avma;
     894                 :      28137 :   for(av=avma;;avma=av)
     895                 :            :   {
     896                 :      39402 :     w = random_Flx(dv,v,p);
     897                 :      39402 :     w = Flx_try_pow(w,*t,p,q,r);
     898         [ +  + ]:      39402 :     if (!w) continue;
     899                 :      30568 :     w = Flx_Fl_add(w, p-1, p);
     900                 :      30568 :     w = Flx_gcd(*t,w, p);
     901 [ +  + ][ +  + ]:      30568 :     l = degpol(w); if (l && l!=dv) break;
     902                 :      11265 :   }
     903                 :      28137 :   w = Flx_normalize(w, p);
     904                 :      28137 :   w = gerepileupto(av0, w);
     905                 :      28137 :   l /= d; t[l]=Flx_div(*t,w,p); *t=w;
     906                 :      28137 :   Flx_split(t+l,d,p,q,r);
     907                 :      28137 :   Flx_split(t,  d,p,q,r);
     908                 :            : }
     909                 :            : 
     910                 :            : 
     911                 :            : static void
     912                 :        310 : FpX_split(GEN *t, long d, GEN  p, GEN q, long r)
     913                 :            : {
     914                 :        310 :   long l, v, dv = degpol(*t);
     915                 :            :   pari_sp av;
     916                 :            :   GEN w;
     917                 :            : 
     918         [ +  + ]:        460 :   if (dv==d) return;
     919                 :        150 :   v = varn(*t);
     920                 :        150 :   av = avma;
     921                 :         52 :   for(;; avma = av)
     922                 :            :   {
     923                 :        202 :     w = random_FpX(dv, v, p);
     924                 :        202 :     w = try_pow(w,*t,p,q,r);
     925         [ -  + ]:        202 :     if (!w) continue;
     926                 :        202 :     w = FpX_Fp_sub_shallow(w, gen_1, p);
     927                 :        202 :     w = FpX_gcd(*t,w, p); l=degpol(w);
     928 [ +  + ][ +  - ]:        202 :     if (l && l!=dv) break;
     929                 :         52 :   }
     930                 :        150 :   w = FpX_normalize(w, p);
     931                 :        150 :   w = gerepileupto(av, w);
     932                 :        150 :   l /= d; t[l]=FpX_div(*t,w,p); *t=w;
     933                 :        150 :   FpX_split(t+l,d,p,q,r);
     934                 :        150 :   FpX_split(t,  d,p,q,r);
     935                 :            : }
     936                 :            : 
     937                 :            : static int
     938         [ +  + ]:     999425 : cmpGuGu(GEN a, GEN b) { return (ulong)a < (ulong)b? -1: (a == b? 0: 1); }
     939                 :            : 
     940                 :            : /* p > 2 */
     941                 :            : static GEN
     942                 :          5 : FpX_is_irred_2(GEN f, GEN p, long d)
     943                 :            : {
     944      [ -  -  + ]:          5 :   switch(d)
     945                 :            :   {
     946                 :            :     case -1:
     947                 :          0 :     case 0: return NULL;
     948                 :          0 :     case 1: return gen_1;
     949                 :            :   }
     950         [ -  + ]:          5 :   return FpX_quad_factortype(f, p) == -1? gen_1: NULL;
     951                 :            : }
     952                 :            : /* p > 2 */
     953                 :            : static GEN
     954                 :         10 : FpX_degfact_2(GEN f, GEN p, long d)
     955                 :            : {
     956   [ -  -  -  + ]:         10 :   switch(d)
     957                 :            :   {
     958                 :          0 :     case -1:retmkvec2(mkvecsmall(-1),mkvecsmall(1));
     959                 :          0 :     case 0: return trivial_fact();
     960                 :          0 :     case 1: retmkvec2(mkvecsmall(1), mkvecsmall(1));
     961                 :            :   }
     962      [ +  +  - ]:         10 :   switch(FpX_quad_factortype(f, p)) {
     963                 :          5 :     case  1: retmkvec2(mkvecsmall2(1,1), mkvecsmall2(1,1));
     964                 :          5 :     case -1: retmkvec2(mkvecsmall(2), mkvecsmall(1));
     965                 :         10 :     default: retmkvec2(mkvecsmall(1), mkvecsmall(2));
     966                 :            :   }
     967                 :            : }
     968                 :            : 
     969                 :            : GEN
     970                 :         25 : prime_fact(GEN x) { retmkmat2(mkcolcopy(x), mkcol(gen_1)); }
     971                 :            : GEN
     972                 :     154900 : trivial_fact(void) { retmkmat2(cgetg(1,t_COL), cgetg(1,t_COL)); }
     973                 :            : 
     974                 :            : /* Mod(0,p) * x, where x is f's main variable */
     975                 :            : static GEN
     976                 :         10 : Mod0pX(GEN f, GEN p)
     977                 :         10 : { return scalarpol(mkintmod(gen_0, p), varn(f)); }
     978                 :            : static GEN
     979                 :         10 : zero_fact_intmod(GEN f, GEN p) { return prime_fact(Mod0pX(f,p)); }
     980                 :            : 
     981                 :            : /* not gerepile safe */
     982                 :            : static GEN
     983                 :         42 : FpX_factor_2(GEN f, GEN p, long d)
     984                 :            : {
     985                 :            :   GEN r, s, R, S;
     986                 :            :   long v;
     987                 :            :   int sgn;
     988   [ -  +  -  + ]:         42 :   switch(d)
     989                 :            :   {
     990                 :          0 :     case -1: retmkvec2(mkcol(pol_0(varn(f))), mkvecsmall(1));
     991                 :          2 :     case  0: retmkvec2(cgetg(1,t_COL), cgetg(1,t_VECSMALL));
     992                 :          0 :     case  1: retmkvec2(mkcol(f), mkvecsmall(1));
     993                 :            :   }
     994                 :         40 :   r = FpX_quad_root(f, p, 1);
     995         [ +  + ]:         40 :   if (!r) return mkvec2(mkcol(f), mkvecsmall(1));
     996                 :         33 :   v = varn(f);
     997                 :         33 :   s = FpX_otherroot(f, r, p);
     998         [ +  - ]:         33 :   if (signe(r)) r = subii(p, r);
     999         [ +  - ]:         33 :   if (signe(s)) s = subii(p, s);
    1000         [ +  - ]:         33 :   sgn = cmpii(s, r); if (sgn < 0) swap(s,r);
    1001                 :         33 :   R = deg1pol_shallow(gen_1, r, v);
    1002         [ -  + ]:         33 :   if (!sgn) return mkvec2(mkcol(R), mkvecsmall(2));
    1003                 :         33 :   S = deg1pol_shallow(gen_1, s, v);
    1004                 :         42 :   return mkvec2(mkcol2(R,S), mkvecsmall2(1,1));
    1005                 :            : }
    1006                 :            : static GEN
    1007                 :         57 : FpX_factor_deg2(GEN f, GEN p, long d, long flag)
    1008                 :            : {
    1009      [ +  +  + ]:         57 :   switch(flag) {
    1010                 :          5 :     case 2: return FpX_is_irred_2(f, p, d);
    1011                 :         10 :     case 1: return FpX_degfact_2(f, p, d);
    1012                 :         57 :     default: return FpX_factor_2(f, p, d);
    1013                 :            :   }
    1014                 :            : }
    1015                 :            : 
    1016                 :            : static int
    1017                 :      17251 : F2x_quad_factortype(GEN x)
    1018         [ +  + ]:      17251 : { return x[2] == 7 ? -1: x[2] == 6 ? 1 :0; }
    1019                 :            : 
    1020                 :            : static GEN
    1021                 :          5 : F2x_is_irred_2(GEN f, long d)
    1022 [ +  - ][ +  - ]:          5 : { return d == 1 || (d==2 && F2x_quad_factortype(f) == -1)? gen_1: NULL; }
                 [ +  - ]
    1023                 :            : 
    1024                 :            : static GEN
    1025                 :        420 : F2x_degfact_2(GEN f, long d)
    1026                 :            : {
    1027         [ -  + ]:        420 :   if (!d) return trivial_fact();
    1028         [ +  + ]:        420 :   if (d == 1) return mkvec2(mkvecsmall(1), mkvecsmall(1));
    1029      [ +  +  + ]:        415 :   switch(F2x_quad_factortype(f)) {
    1030                 :         95 :     case 1: return mkvec2(mkvecsmall2(1,1), mkvecsmall2(1,1));
    1031                 :         90 :     case -1:return mkvec2(mkvecsmall(2), mkvecsmall(1));
    1032                 :        420 :     default: return mkvec2(mkvecsmall(1), mkvecsmall(2));
    1033                 :            :   }
    1034                 :            : }
    1035                 :            : 
    1036                 :            : static GEN
    1037                 :      21515 : F2x_factor_2(GEN f, long d)
    1038                 :            : {
    1039                 :      21515 :   long v = f[1];
    1040         [ -  + ]:      21515 :   if (d < 0) pari_err(e_ROOTS0,"Flx_factor_2");
    1041         [ +  + ]:      21515 :   if (!d) return mkvec2(cgetg(1,t_COL), cgetg(1,t_VECSMALL));
    1042         [ +  + ]:      16290 :   if (d == 1) return mkvec2(mkcol(f), mkvecsmall(1));
    1043      [ +  +  + ]:      11750 :   switch(F2x_quad_factortype(f))
    1044                 :            :   {
    1045                 :       2070 :   case -1: return mkvec2(mkcol(f), mkvecsmall(1));
    1046                 :       8614 :   case 0:  return mkvec2(mkcol(mkvecsmall2(v,2+F2x_coeff(f,0))), mkvecsmall(2));
    1047                 :      21515 :   default: return mkvec2(mkcol2(mkvecsmall2(v,2),mkvecsmall2(v,3)), mkvecsmall2(1,1));
    1048                 :            :   }
    1049                 :            : }
    1050                 :            : static GEN
    1051                 :      21940 : F2x_factor_deg2(GEN f, long d, long flag)
    1052                 :            : {
    1053      [ +  +  + ]:      21940 :   switch(flag) {
    1054                 :          5 :     case 2: return F2x_is_irred_2(f, d);
    1055                 :        420 :     case 1: return F2x_degfact_2(f, d);
    1056                 :      21940 :     default: return F2x_factor_2(f, d);
    1057                 :            :   }
    1058                 :            : }
    1059                 :            : 
    1060                 :            : static void
    1061                 :         13 : split_squares(struct split_t *S, GEN g, ulong p)
    1062                 :            : {
    1063                 :         13 :   ulong q = p >> 1;
    1064                 :         13 :   GEN a = Flx_mod_Xnm1(g, q, p); /* mod x^(p-1)/2 - 1 */
    1065         [ +  + ]:         13 :   if (degpol(a) < 0)
    1066                 :            :   {
    1067                 :            :     ulong i;
    1068                 :          4 :     split_add_done(S, (GEN)1);
    1069         [ +  + ]:         12 :     for (i = 2; i <= q; i++) split_add_done(S, (GEN)Fl_sqr(i,p));
    1070                 :            :   } else {
    1071                 :          9 :     (void)Flx_valrem(a, &a);
    1072         [ +  - ]:          9 :     if (degpol(a))
    1073                 :            :     {
    1074                 :          9 :       a = Flx_gcd(a, Flx_Xnm1(g[1], q, p), p);
    1075         [ +  - ]:          9 :       if (degpol(a)) split_add(S, Flx_normalize(a, p));
    1076                 :            :     }
    1077                 :            :   }
    1078                 :         13 : }
    1079                 :            : static void
    1080                 :         13 : split_nonsquares(struct split_t *S, GEN g, ulong p)
    1081                 :            : {
    1082                 :         13 :   ulong q = p >> 1;
    1083                 :         13 :   GEN a = Flx_mod_Xn1(g, q, p); /* mod x^(p-1)/2 + 1 */
    1084         [ +  + ]:         13 :   if (degpol(a) < 0)
    1085                 :            :   {
    1086                 :          4 :     ulong i, z = Fl_nonsquare(p);
    1087                 :          4 :     split_add_done(S, (GEN)z);
    1088         [ +  + ]:         12 :     for (i = 2; i <= q; i++) split_add_done(S, (GEN)Fl_mul(z, Fl_sqr(i,p), p));
    1089                 :            :   } else {
    1090                 :          9 :     (void)Flx_valrem(a, &a);
    1091         [ +  - ]:          9 :     if (degpol(a))
    1092                 :            :     {
    1093                 :          9 :       a = Flx_gcd(a, Flx_Xn1(g[1], q, p), p);
    1094         [ +  - ]:          9 :       if (degpol(a)) split_add(S, Flx_normalize(a, p));
    1095                 :            :     }
    1096                 :            :   }
    1097                 :         13 : }
    1098                 :            : /* p > 2. f monic Flx, f(0) != 0. Add to split_t structs coprime factors
    1099                 :            :  * of g = \prod_{f(a) = 0} (X - a). Return 0 when f(x) = 0 for all x in Fp* */
    1100                 :            : static int
    1101                 :      12235 : split_Flx_cut_out_roots(struct split_t *S, GEN f, ulong p)
    1102                 :            : {
    1103                 :      12235 :   GEN a, g = Flx_mod_Xnm1(f, p-1, p); /* f mod x^(p-1) - 1 */
    1104                 :      12235 :   long d = degpol(g);
    1105         [ +  + ]:      12235 :   if (d < 0) return 0;
    1106         [ -  + ]:      12210 :   if (g != f) { (void)Flx_valrem(g, &g); d = degpol(g); } /*kill powers of x*/
    1107         [ +  + ]:      12210 :   if (!d) return 1;
    1108         [ +  + ]:       1155 :   if (p <= 1.4 * (ulong)d) {
    1109                 :            :     /* small p; split further using x^((p-1)/2) +/- 1.
    1110                 :            :      * 30% degree drop makes the extra gcd worth it. */
    1111                 :         13 :     split_squares(S, g, p);
    1112                 :         13 :     split_nonsquares(S, g, p);
    1113                 :            :   } else { /* large p; use x^(p-1) - 1 directly */
    1114                 :       1142 :     a = Flxq_powu(polx_Flx(f[1]), p-1, g,p);
    1115         [ -  + ]:       1142 :     if (lg(a) < 3) pari_err_PRIME("rootmod",utoipos(p));
    1116                 :       1142 :     a = Flx_Fl_add(a, p-1, p); /* a = x^(p-1) - 1 mod g */
    1117                 :       1142 :     g = Flx_gcd(g,a, p);
    1118         [ +  - ]:       1142 :     if (degpol(g)) split_add(S, Flx_normalize(g,p));
    1119                 :            :   }
    1120                 :      12235 :   return 1;
    1121                 :            : }
    1122                 :            : 
    1123                 :            : /* by splitting, assume p > 2 prime, deg(f) > 0, and f monic */
    1124                 :            : static GEN
    1125                 :      35139 : Flx_roots_i(GEN f, ulong p)
    1126                 :            : {
    1127                 :            :   GEN pol, g;
    1128                 :      35139 :   long v = Flx_valrem(f, &g);
    1129                 :            :   ulong q;
    1130                 :            :   struct split_t S;
    1131                 :            : 
    1132                 :            :   /* optimization: test for small degree first */
    1133      [ +  +  + ]:      35139 :   switch(degpol(g))
    1134                 :            :   {
    1135                 :            :     case 1: {
    1136                 :      16289 :       ulong r = p - g[2];
    1137         [ +  + ]:      16289 :       return v? mkvecsmall2(0, r): mkvecsmall(r);
    1138                 :            :     }
    1139                 :            :     case 2: {
    1140                 :       6615 :       ulong r = Flx_quad_root(g, p, 1), s;
    1141 [ +  + ][ -  + ]:       6615 :       if (r == p) return v? mkvecsmall(0): cgetg(1,t_VECSMALL);
    1142                 :       6600 :       s = Flx_otherroot(g,r, p);
    1143         [ +  + ]:       6600 :       if (r < s)
    1144         [ +  + ]:       3956 :         return v? mkvecsmall3(0, r, s): mkvecsmall2(r, s);
    1145         [ +  + ]:       2644 :       else if (r > s)
    1146         [ +  + ]:       2599 :         return v? mkvecsmall3(0, s, r): mkvecsmall2(s, r);
    1147                 :            :       else
    1148         [ -  + ]:         45 :         return v? mkvecsmall2(0, s): mkvecsmall(s);
    1149                 :            :     }
    1150                 :            :   }
    1151                 :      12235 :   q = p >> 1;
    1152                 :      12235 :   split_init(&S, lg(f)-1);
    1153                 :      12235 :   settyp(S.done, t_VECSMALL);
    1154         [ +  + ]:      12235 :   if (v) split_add_done(&S, (GEN)0);
    1155         [ +  + ]:      12235 :   if (! split_Flx_cut_out_roots(&S, g, p))
    1156                 :         25 :     return all_roots_mod_p(p, lg(S.done) == 1);
    1157                 :      12210 :   pol = polx_Flx(f[1]);
    1158                 :      12210 :   for (pol[2]=1; ; pol[2]++)
    1159                 :            :   {
    1160                 :      16409 :     long j, l = lg(S.todo);
    1161         [ +  + ]:      16409 :     if (l == 1) { vecsmall_sort(S.done); return S.done; }
    1162 [ -  + ][ #  # ]:       4199 :     if (pol[2] == 100 && !uisprime(p)) pari_err_PRIME("polrootsmod",utoipos(p));
    1163         [ +  + ]:      18762 :     for (j = 1; j < l; j++)
    1164                 :            :     {
    1165                 :      14563 :       GEN c = gel(S.todo,j);
    1166      [ +  +  + ]:      14563 :       switch(degpol(c))
    1167                 :            :       {
    1168                 :            :         case 1:
    1169                 :       2762 :           split_moveto_done(&S, j, (GEN)(p - c[2]));
    1170                 :       2762 :           j--; l--; break;
    1171                 :            :         case 2: {
    1172                 :       4497 :           ulong r = Flx_quad_root(c, p, 0), s = Flx_otherroot(c,r, p);
    1173                 :       4497 :           split_done(&S, j, (GEN)r, (GEN)s);
    1174                 :       4497 :           j--; l--; break;
    1175                 :            :         }
    1176                 :            :         default: {
    1177                 :       7304 :           GEN b = Flx_Fl_add(Flxq_powu(pol,q, c,p), p-1, p); /* pol^(p-1)/2 */
    1178                 :            :           long db;
    1179                 :       7304 :           b = Flx_gcd(c,b, p); db = degpol(b);
    1180 [ +  + ][ +  + ]:       7304 :           if (db && db < degpol(c))
    1181                 :            :           {
    1182                 :       6099 :             b = Flx_normalize(b, p);
    1183                 :       6099 :             c = Flx_div(c,b, p);
    1184                 :       6099 :             split_todo(&S, j, b, c);
    1185                 :            :           }
    1186                 :            :         }
    1187                 :            :       }
    1188                 :            :     }
    1189                 :      39338 :   }
    1190                 :            : }
    1191                 :            : 
    1192                 :            : GEN
    1193                 :          0 : Flx_roots(GEN f, ulong p)
    1194                 :            : {
    1195                 :          0 :   pari_sp av = avma;
    1196      [ #  #  # ]:          0 :   switch(lg(f))
    1197                 :            :   {
    1198                 :          0 :     case 2: pari_err_ROOTS0("Flx_roots");
    1199                 :          0 :     case 3: avma = av; return cgetg(1, t_VECSMALL);
    1200                 :            :   }
    1201         [ #  # ]:          0 :   if (p == 2) return Flx_root_mod_2(f);
    1202                 :          0 :   return gerepileuptoleaf(av, Flx_roots_i(Flx_normalize(f, p), p));
    1203                 :            : }
    1204                 :            : 
    1205                 :            : /* assume x reduced mod p, monic. */
    1206                 :            : static int
    1207                 :      20370 : Flx_quad_factortype(GEN x, ulong p)
    1208                 :            : {
    1209                 :      20370 :   ulong b = x[3], c = x[2];
    1210                 :      20370 :   return krouu(Fl_disc_bc(b, c, p), p);
    1211                 :            : }
    1212                 :            : static GEN
    1213                 :          0 : Flx_is_irred_2(GEN f, ulong p, long d)
    1214                 :            : {
    1215         [ #  # ]:          0 :   if (!d) return NULL;
    1216         [ #  # ]:          0 :   if (d == 1) return gen_1;
    1217         [ #  # ]:          0 :   return Flx_quad_factortype(f, p) == -1? gen_1: NULL;
    1218                 :            : }
    1219                 :            : static GEN
    1220                 :      20380 : Flx_degfact_2(GEN f, ulong p, long d)
    1221                 :            : {
    1222         [ -  + ]:      20380 :   if (!d) return trivial_fact();
    1223         [ +  + ]:      20380 :   if (d == 1) return mkvec2(mkvecsmall(1), mkvecsmall(1));
    1224      [ +  +  + ]:      20370 :   switch(Flx_quad_factortype(f, p)) {
    1225                 :       9530 :     case 1: return mkvec2(mkvecsmall2(1,1), mkvecsmall2(1,1));
    1226                 :      10405 :     case -1:return mkvec2(mkvecsmall(2), mkvecsmall(1));
    1227                 :      20380 :     default: return mkvec2(mkvecsmall(1), mkvecsmall(2));
    1228                 :            :   }
    1229                 :            : }
    1230                 :            : /* p > 2 */
    1231                 :            : static GEN
    1232                 :     141043 : Flx_factor_2(GEN f, ulong p, long d)
    1233                 :            : {
    1234                 :            :   ulong r, s;
    1235                 :            :   GEN R,S;
    1236                 :     141043 :   long v = f[1];
    1237         [ -  + ]:     141043 :   if (d < 0) pari_err(e_ROOTS0,"Flx_factor_2");
    1238         [ +  + ]:     141043 :   if (!d) return mkvec2(cgetg(1,t_COL), cgetg(1,t_VECSMALL));
    1239         [ +  + ]:     131360 :   if (d == 1) return mkvec2(mkcol(f), mkvecsmall(1));
    1240                 :      82935 :   r = Flx_quad_root(f, p, 1);
    1241         [ +  + ]:      82935 :   if (r==p) return mkvec2(mkcol(f), mkvecsmall(1));
    1242                 :      56697 :   s = Flx_otherroot(f, r, p);
    1243                 :      56697 :   r = Fl_neg(r, p);
    1244                 :      56697 :   s = Fl_neg(s, p);
    1245         [ +  + ]:      56697 :   if (s < r) lswap(s,r);
    1246                 :      56697 :   R = mkvecsmall3(v,r,1);
    1247         [ +  + ]:      56697 :   if (s == r) return mkvec2(mkcol(R), mkvecsmall(2));
    1248                 :      47242 :   S = mkvecsmall3(v,s,1);
    1249                 :     141043 :   return mkvec2(mkcol2(R,S), mkvecsmall2(1,1));
    1250                 :            : }
    1251                 :            : static GEN
    1252                 :     161423 : Flx_factor_deg2(GEN f, ulong p, long d, long flag)
    1253                 :            : {
    1254      [ -  +  + ]:     161423 :   switch(flag) {
    1255                 :          0 :     case 2: return Flx_is_irred_2(f, p, d);
    1256                 :      20380 :     case 1: return Flx_degfact_2(f, p, d);
    1257                 :     161423 :     default: return Flx_factor_2(f, p, d);
    1258                 :            :   }
    1259                 :            : }
    1260                 :            : 
    1261                 :            : static GEN
    1262                 :       6886 : F2x_factcantor_i(GEN f, long flag)
    1263                 :            : {
    1264                 :       6886 :   long j, e, nbfact, d = F2x_degree(f);
    1265                 :            :   ulong k;
    1266                 :            :   GEN X, E, f2, g, g1, u, v, y, t;
    1267                 :            : 
    1268         [ +  + ]:       6886 :   if (d <= 2) return F2x_factor_deg2(f, d, flag);
    1269                 :            :   /* to hold factors and exponents */
    1270         [ +  + ]:       6431 :   t = flag ? cgetg(d+1,t_VECSMALL): cgetg(d+1,t_VEC);
    1271                 :       6431 :   E = cgetg(d+1, t_VECSMALL);
    1272                 :       6431 :   X = polx_F2x(f[1]);
    1273                 :       6431 :   e = nbfact = 1;
    1274                 :            :   for(;;)
    1275                 :            :   {
    1276                 :       8714 :     f2 = F2x_gcd(f,F2x_deriv(f));
    1277 [ +  + ][ +  + ]:       8714 :     if (flag == 2 && F2x_degree(f2) > 0) return NULL;
    1278                 :       8684 :     g1 = F2x_div(f,f2);
    1279                 :       8684 :     k = 0;
    1280         [ +  + ]:      17462 :     while (F2x_degree(g1)>0)
    1281                 :            :     {
    1282                 :            :       pari_sp av;
    1283                 :            :       long du, dg, dg1;
    1284         [ +  + ]:       8843 :       k++; if (k%2==0) { k++; f2 = F2x_div(f2,g1); }
    1285                 :       8843 :       u = g1; g1 = F2x_gcd(f2,g1);
    1286                 :       8843 :       du = F2x_degree(u);
    1287                 :       8843 :       dg1 = F2x_degree(g1);
    1288         [ +  + ]:       8843 :       if (dg1>0)
    1289                 :            :       {
    1290                 :        538 :         f2= F2x_div(f2,g1);
    1291         [ +  + ]:        538 :         if (du == dg1) continue;
    1292                 :        327 :         u = F2x_div( u,g1);
    1293                 :        327 :         du -= dg1;
    1294                 :            :       }
    1295                 :            :       /* here u is square-free (product of irred. of multiplicity e * k) */
    1296                 :       8632 :       v = X;
    1297                 :       8632 :       av = avma;
    1298         [ +  + ]:      46983 :       for (d=1; d <= du>>1; d++)
    1299                 :            :       {
    1300                 :      38416 :         v = F2xq_sqr(v, u);
    1301                 :      38416 :         g = F2x_gcd(F2x_add(v, X), u);
    1302                 :      38416 :         dg = F2x_degree(g);
    1303         [ +  + ]:      38416 :         if (dg <= 0) {avma = (pari_sp)v; v = gerepileuptoleaf(av,v); continue;}
    1304                 :            :         /* g is a product of irred. pols, all of which have degree d */
    1305                 :      11455 :         j = nbfact+dg/d;
    1306         [ +  + ]:      11455 :         if (flag)
    1307                 :            :         {
    1308         [ +  + ]:        402 :           if (flag == 2) return NULL;
    1309         [ +  + ]:        762 :           for ( ; nbfact<j; nbfact++) { t[nbfact]=d; E[nbfact]=e*k; }
    1310                 :            :         }
    1311                 :            :         else
    1312                 :            :         {
    1313                 :      11053 :           gel(t,nbfact) = g;
    1314                 :      11053 :           F2x_split(2,&gel(t,nbfact),d);
    1315         [ +  + ]:      23616 :           for (; nbfact<j; nbfact++) E[nbfact]=e*k;
    1316                 :            :         }
    1317                 :      11390 :         du -= dg;
    1318                 :      11390 :         u = F2x_div(u,g);
    1319                 :      11390 :         v = F2x_rem(v,u);
    1320                 :      11390 :         av = avma;
    1321                 :            :       }
    1322         [ +  + ]:       8567 :       if (du)
    1323                 :            :       {
    1324         [ +  + ]:       7680 :         if (flag) t[nbfact] = du;
    1325                 :       7021 :         else  gel(t,nbfact) = u;
    1326                 :       7680 :         E[nbfact++]=e*k;
    1327                 :            :       }
    1328                 :            :     }
    1329         [ +  + ]:       8619 :     if (F2x_degree(f2)==0) break;
    1330                 :       2283 :     e <<=1; f = F2x_sqrt(f2);
    1331                 :       2283 :   }
    1332         [ +  + ]:       6336 :   if (flag == 2) return gen_1; /* irreducible */
    1333                 :       6331 :   setlg(t, nbfact);
    1334                 :       6331 :   setlg(E, nbfact); y = mkvec2(t, E);
    1335                 :       6886 :   return flag ? sort_factor(y, (void*)cmpGuGu, cmp_nodata)
    1336         [ +  + ]:       6331 :               : sort_factor_pol(y, cmpGuGu);
    1337                 :            : }
    1338                 :            : GEN
    1339                 :       5248 : F2x_factcantor(GEN f, long flag)
    1340                 :            : {
    1341                 :       5248 :   pari_sp av = avma;
    1342                 :       5248 :   GEN z = F2x_factcantor_i(f, flag);
    1343         [ -  + ]:       5248 :   if (flag == 2) { avma = av; return z; }
    1344                 :       5248 :   return gerepilecopy(av, z);
    1345                 :            : }
    1346                 :            : 
    1347                 :            : int
    1348                 :        100 : F2x_is_irred(GEN f) { return !!F2x_factcantor_i(f,2); }
    1349                 :            : 
    1350                 :            : void
    1351                 :        396 : F2xV_to_FlxV_inplace(GEN v)
    1352                 :            : {
    1353                 :            :   long i;
    1354         [ +  + ]:       1273 :   for(i=1;i<lg(v);i++) gel(v,i)= F2x_to_Flx(gel(v,i));
    1355                 :        396 : }
    1356                 :            : void
    1357                 :     301134 : FlxV_to_ZXV_inplace(GEN v)
    1358                 :            : {
    1359                 :            :   long i;
    1360         [ +  + ]:     799633 :   for(i=1;i<lg(v);i++) gel(v,i)= Flx_to_ZX(gel(v,i));
    1361                 :     301134 : }
    1362                 :            : void
    1363                 :      53268 : F2xV_to_ZXV_inplace(GEN v)
    1364                 :            : {
    1365                 :            :   long i;
    1366         [ +  + ]:     111871 :   for(i=1;i<lg(v);i++) gel(v,i)= F2x_to_ZX(gel(v,i));
    1367                 :      53268 : }
    1368                 :            : 
    1369                 :            : /* factor f mod pp.
    1370                 :            :  * flag = 1: return the degrees, not the factors
    1371                 :            :  * flag = 2: return NULL if f is not irreducible */
    1372                 :            : static GEN
    1373                 :     123900 : Flx_factcantor_i(GEN f, ulong p, long flag)
    1374                 :            : {
    1375                 :     123900 :   long j, e, nbfact, d = degpol(f);
    1376                 :            :   ulong k;
    1377                 :            :   GEN X, E, f2, g, g1, u, v, q, y, t;
    1378         [ +  + ]:     123900 :   if (p==2) { /*We need to handle 2 specially */
    1379                 :       1523 :     GEN F = F2x_factcantor_i(Flx_to_F2x(f),flag);
    1380         [ +  + ]:       1523 :     if (flag==0) F2xV_to_FlxV_inplace(gel(F,1));
    1381                 :       1523 :     return F;
    1382                 :            :   }
    1383                 :            :   /* Now we assume p odd */
    1384         [ +  + ]:     122377 :   if (d <= 2) return Flx_factor_deg2(f, p, d, flag);
    1385                 :            : 
    1386                 :            :   /* to hold factors and exponents */
    1387         [ +  + ]:     101951 :   t = cgetg(d+1, flag? t_VECSMALL: t_VEC);
    1388                 :     101951 :   E = cgetg(d+1, t_VECSMALL);
    1389                 :     101951 :   X = polx_Flx(f[1]);
    1390                 :     101951 :   e = nbfact = 1;
    1391                 :            :   for(;;)
    1392                 :            :   {
    1393                 :     102186 :     f2 = Flx_gcd(f,Flx_deriv(f,p), p);
    1394 [ +  + ][ +  + ]:     102186 :     if (flag == 2 && degpol(f2) > 0) return NULL;
    1395                 :     102086 :     g1 = Flx_div(f,f2,p);
    1396                 :     102086 :     k = 0;
    1397         [ +  + ]:     207177 :     while (degpol(g1)>0)
    1398                 :            :     {
    1399                 :            :       pari_sp av;
    1400                 :            :       long du,dg;
    1401         [ +  + ]:     105391 :       k++; if (k%p==0) { k++; f2 = Flx_div(f2,g1,p); }
    1402                 :     105391 :       u = g1; g1 = Flx_gcd(f2,g1, p);
    1403         [ +  + ]:     105391 :       if (degpol(g1)>0)
    1404                 :            :       {
    1405                 :       3440 :         u = Flx_div( u,g1,p);
    1406                 :       3440 :         f2= Flx_div(f2,g1,p);
    1407                 :            :       }
    1408                 :     105391 :       du = degpol(u);
    1409         [ +  + ]:     105391 :       if (du <= 0) continue;
    1410                 :            : 
    1411                 :            :       /* here u is square-free (product of irred. of multiplicity e * k) */
    1412                 :     104190 :       v=X;
    1413                 :     104190 :       av = avma;
    1414         [ +  + ]:     348520 :       for (d=1; d <= du>>1; d++)
    1415                 :            :       {
    1416                 :     244630 :         v = Flxq_powu(v, p, u, p);
    1417                 :     244630 :         g = Flx_gcd(Flx_sub(v, X, p), u, p);
    1418                 :     244630 :         dg = degpol(g);
    1419         [ +  + ]:     244630 :         if (dg <= 0) {avma = (pari_sp)v; v = gerepileuptoleaf(av,v); continue;}
    1420                 :            :         /* g is a product of irred. pols, all of which have degree d */
    1421                 :     101493 :         j = nbfact+dg/d;
    1422         [ +  + ]:     101493 :         if (flag)
    1423                 :            :         {
    1424         [ +  + ]:      83333 :           if (flag == 2) return NULL;
    1425         [ +  + ]:     277187 :           for ( ; nbfact<j; nbfact++) { t[nbfact]=d; E[nbfact]=e*k; }
    1426                 :            :         }
    1427                 :            :         else
    1428                 :            :         {
    1429                 :      18160 :           GEN pd = powuu(p, d);
    1430                 :            :           long r;
    1431                 :      18160 :           g = Flx_normalize(g, p);
    1432                 :      18160 :           gel(t,nbfact) = g; q = subis(pd,1);
    1433                 :      18160 :           r = vali(q); q = shifti(pd,-r);
    1434                 :      18160 :           Flx_split(&gel(t,nbfact),d,p,q,r);
    1435         [ +  + ]:      64457 :           for (; nbfact<j; nbfact++) E[nbfact]=e*k;
    1436                 :            :         }
    1437                 :     101193 :         du -= dg;
    1438                 :     101193 :         u = Flx_div(u,g,p);
    1439                 :     101193 :         v = Flx_rem(v,u,p);
    1440                 :     101193 :         av = avma;
    1441                 :            :       }
    1442         [ +  + ]:     103890 :       if (du)
    1443                 :            :       {
    1444         [ +  + ]:      46354 :         if (flag) t[nbfact] = du;
    1445                 :       7379 :         else  gel(t,nbfact) = Flx_normalize(u, p);
    1446                 :      46354 :         E[nbfact++]=e*k;
    1447                 :            :       }
    1448                 :            :     }
    1449         [ +  + ]:     101786 :     if (degpol(f2)==0) break;
    1450                 :        235 :     e *= p; f = Flx_deflate(f2, p);
    1451                 :        235 :   }
    1452         [ +  + ]:     101551 :   if (flag == 2) return gen_1; /* irreducible */
    1453                 :     101536 :   setlg(t, nbfact);
    1454                 :     101536 :   setlg(E, nbfact); y = mkvec2(t, E);
    1455                 :     123900 :   return flag ? sort_factor(y, (void*)cmpGuGu, cmp_nodata)
    1456         [ +  + ]:     101536 :               : sort_factor_pol(y, cmpGuGu);
    1457                 :            : }
    1458                 :            : GEN
    1459                 :       5156 : Flx_factcantor(GEN f, ulong p, long flag)
    1460                 :            : {
    1461                 :       5156 :   pari_sp av = avma;
    1462                 :       5156 :   GEN z = Flx_factcantor_i(Flx_normalize(f,p),p,flag);
    1463         [ -  + ]:       5156 :   if (flag == 2) { avma = av; return z; }
    1464                 :       5156 :   return gerepilecopy(av, z);
    1465                 :            : }
    1466                 :            : GEN
    1467                 :     109477 : Flx_degfact(GEN f, ulong p)
    1468                 :            : {
    1469                 :     109477 :   pari_sp av = avma;
    1470                 :     109477 :   GEN z = Flx_factcantor_i(Flx_normalize(f,p),p,1);
    1471                 :     109477 :   return gerepilecopy(av, z);
    1472                 :            : }
    1473                 :            : int
    1474                 :        415 : Flx_is_irred(GEN f, ulong p) { return !!Flx_factcantor_i(f,p,2); }
    1475                 :            : 
    1476                 :            : /* factor f (FpX or Flx) mod pp.
    1477                 :            :  * flag = 1: return the degrees, not the factors
    1478                 :            :  * flag = 2: return NULL if f is not irreducible.
    1479                 :            :  * Not gerepile-safe */
    1480                 :            : static GEN
    1481                 :         65 : FpX_factcantor_i(GEN f, GEN pp, long flag)
    1482                 :            : {
    1483                 :         65 :   long j, nbfact, d = degpol(f);
    1484                 :            :   ulong k;
    1485                 :            :   GEN X, E,y,f2,g,g1,v,q;
    1486                 :            :   GEN t;
    1487         [ +  + ]:         65 :   if (typ(f) == t_VECSMALL)
    1488                 :            :   { /* lgefint(pp) = 3 */
    1489                 :            :     GEN F;
    1490                 :         40 :     ulong p = pp[2];
    1491         [ +  + ]:         40 :     if (p==2) {
    1492                 :         15 :       F = F2x_factcantor_i(Flx_to_F2x(f),flag);
    1493         [ +  + ]:         15 :       if (flag==0) F2xV_to_ZXV_inplace(gel(F,1));
    1494                 :            :     } else {
    1495                 :         25 :       F = Flx_factcantor_i(f,p,flag);
    1496         [ +  + ]:         25 :       if (flag==0) FlxV_to_ZXV_inplace(gel(F,1));
    1497                 :            :     }
    1498                 :         40 :     return F;
    1499                 :            :   }
    1500                 :            :   /*Now, we can assume that p is large and odd*/
    1501         [ +  + ]:         25 :   if (d <= 2) return FpX_factor_deg2(f,pp,d,flag);
    1502                 :            : 
    1503                 :            :   /* to hold factors and exponents */
    1504         [ -  + ]:          5 :   t = cgetg(d+1, flag? t_VECSMALL: t_VEC);
    1505                 :          5 :   E = cgetg(d+1, t_VECSMALL);
    1506                 :          5 :   X = pol_x(varn(f)); nbfact = 1;
    1507                 :          5 :   f2 = FpX_gcd(f,ZX_deriv(f), pp);
    1508 [ -  + ][ #  # ]:          5 :   if (flag == 2 && degpol(f2)) return NULL;
    1509         [ +  - ]:          5 :   g1 = degpol(f2)? FpX_div(f,f2,pp): f; /* squarefree */
    1510                 :          5 :   k = 0;
    1511         [ +  + ]:         15 :   while (degpol(g1)>0)
    1512                 :            :   {
    1513                 :            :     long du,dg;
    1514                 :            :     GEN u, S;
    1515                 :         10 :     k++;
    1516                 :         10 :     u = g1; g1 = FpX_gcd(f2,g1, pp);
    1517         [ +  + ]:         10 :     if (degpol(g1)>0)
    1518                 :            :     {
    1519                 :          5 :       u = FpX_div( u,g1,pp);
    1520                 :          5 :       f2= FpX_div(f2,g1,pp);
    1521                 :            :     }
    1522                 :         10 :     du = degpol(u);
    1523         [ -  + ]:         10 :     if (du <= 0) continue;
    1524                 :            : 
    1525                 :            :     /* here u is square-free (product of irred. of multiplicity e * k) */
    1526                 :         10 :     v=X;
    1527         [ +  + ]:         10 :     S = du==1 ?  cgetg(1, t_VEC): FpXQ_powers(FpXQ_pow(v, pp, u, pp), du-1, u, pp);
    1528         [ +  + ]:         20 :     for (d=1; d <= du>>1; d++)
    1529                 :            :     {
    1530                 :         10 :       v = FpX_FpXQV_eval(v, S, u, pp);
    1531                 :         10 :       g = FpX_gcd(ZX_sub(v, X), u, pp);
    1532                 :         10 :       dg = degpol(g);
    1533         [ -  + ]:         10 :       if (dg <= 0) continue;
    1534                 :            : 
    1535                 :            :       /* g is a product of irred. pols, all of which have degree d */
    1536                 :         10 :       j = nbfact+dg/d;
    1537         [ -  + ]:         10 :       if (flag)
    1538                 :            :       {
    1539         [ #  # ]:          0 :         if (flag == 2) return NULL;
    1540         [ #  # ]:          0 :         for ( ; nbfact<j; nbfact++) { t[nbfact]=d; E[nbfact]=k; }
    1541                 :            :       }
    1542                 :            :       else
    1543                 :            :       {
    1544                 :         10 :         GEN pd = powiu(pp,d);
    1545                 :            :         long r;
    1546                 :         10 :         g = FpX_normalize(g, pp);
    1547                 :         10 :         gel(t,nbfact) = g; q = subis(pd,1);
    1548                 :         10 :         r = vali(q); q = shifti(q,-r);
    1549                 :         10 :         FpX_split(&gel(t,nbfact),d,pp,q,r);
    1550         [ +  + ]:        170 :         for (; nbfact<j; nbfact++) E[nbfact]=k;
    1551                 :            :       }
    1552                 :         10 :       du -= dg;
    1553                 :         10 :       u = FpX_div(u,g,pp);
    1554                 :         10 :       v = FpX_rem(v,u,pp);
    1555                 :            :     }
    1556         [ +  + ]:         10 :     if (du)
    1557                 :            :     {
    1558         [ -  + ]:          5 :       if (flag) t[nbfact] = du;
    1559                 :          5 :       else gel(t,nbfact) = FpX_normalize(u, pp);
    1560                 :          5 :       E[nbfact++]=k;
    1561                 :            :     }
    1562                 :            :   }
    1563         [ -  + ]:          5 :   if (flag == 2) return gen_1; /* irreducible */
    1564                 :          5 :   setlg(t, nbfact);
    1565                 :          5 :   setlg(E, nbfact); y = mkvec2(t, E);
    1566                 :         65 :   return flag ? sort_factor(y, (void*)&cmpGuGu, cmp_nodata)
    1567         [ -  + ]:          5 :               : sort_factor_pol(y, cmpii);
    1568                 :            : }
    1569                 :            : GEN
    1570                 :          0 : FpX_factcantor(GEN f, GEN pp, long flag)
    1571                 :            : {
    1572                 :          0 :   pari_sp av = avma;
    1573                 :            :   GEN z;
    1574                 :          0 :   ZX_factmod_init(&f,pp);
    1575                 :          0 :   z = FpX_factcantor_i(f,pp,flag);
    1576         [ #  # ]:          0 :   if (flag == 2) { avma = av; return z; }
    1577                 :          0 :   return gerepilecopy(av, z);
    1578                 :            : }
    1579                 :            : 
    1580                 :            : static GEN
    1581                 :        120 : factmod_aux(GEN f, GEN p, GEN (*Factor)(GEN,GEN,long), long flag)
    1582                 :            : {
    1583                 :        120 :   pari_sp av = avma;
    1584                 :            :   long j, nbfact;
    1585                 :            :   GEN z, y, t, E, u, v;
    1586                 :            : 
    1587                 :        120 :   factmod_init(&f, p);
    1588      [ +  +  + ]:        120 :   switch(lg(f))
    1589                 :            :   {
    1590                 :         10 :     case 3: avma = av; return trivial_fact();
    1591                 :         10 :     case 2: return gerepileupto(av, zero_fact_intmod(f, p));
    1592                 :            :   }
    1593                 :        100 :   z = Factor(f,p,flag); t = gel(z,1); E = gel(z,2);
    1594                 :        100 :   y = cgetg(3, t_MAT); nbfact = lg(t);
    1595                 :        100 :   u = cgetg(nbfact,t_COL); gel(y,1) = u;
    1596                 :        100 :   v = cgetg(nbfact,t_COL); gel(y,2) = v;
    1597         [ +  + ]:        100 :   if (flag)
    1598         [ +  + ]:         75 :     for (j=1; j<nbfact; j++)
    1599                 :            :     {
    1600                 :         50 :       gel(u,j) = utoi(uel(t,j));
    1601                 :         50 :       gel(v,j) = utoi(uel(E,j));
    1602                 :            :     }
    1603                 :            :   else
    1604         [ +  + ]:       5510 :     for (j=1; j<nbfact; j++)
    1605                 :            :     {
    1606                 :       5435 :       gel(u,j) = FpX_to_mod(gel(t,j), p);
    1607                 :       5435 :       gel(v,j) = utoi(uel(E,j));
    1608                 :            :     }
    1609                 :        120 :   return gerepileupto(av, y);
    1610                 :            : }
    1611                 :            : GEN
    1612                 :         65 : factcantor0(GEN f, GEN p, long flag)
    1613                 :         65 : { return factmod_aux(f, p, &FpX_factcantor_i, flag); }
    1614                 :            : GEN
    1615                 :         55 : factmod(GEN f, GEN p)
    1616                 :         55 : { return factmod_aux(f, p, &FpX_Berlekamp_i, 0); }
    1617                 :            : 
    1618                 :            : /* Use this function when you think f is reducible, and that there are lots of
    1619                 :            :  * factors. If you believe f has few factors, use FpX_nbfact(f,p)==1 instead */
    1620                 :            : int
    1621                 :         10 : FpX_is_irred(GEN f, GEN p) {
    1622                 :         10 :   ZX_factmod_init(&f,p);
    1623                 :         10 :   return !!FpX_factcantor_i(f,p,2);
    1624                 :            : }
    1625                 :            : GEN
    1626                 :          0 : FpX_degfact(GEN f, GEN p) {
    1627                 :          0 :   pari_sp av = avma;
    1628                 :            :   GEN z;
    1629                 :          0 :   ZX_factmod_init(&f,p);
    1630                 :          0 :   z = FpX_factcantor_i(f,p,1);
    1631                 :          0 :   return gerepilecopy(av, z);
    1632                 :            : }
    1633                 :            : GEN
    1634                 :         40 : factcantor(GEN f, GEN p) { return factcantor0(f,p,0); }
    1635                 :            : GEN
    1636                 :         25 : simplefactmod(GEN f, GEN p) { return factcantor0(f,p,1); }
    1637                 :            : 
    1638                 :            : /* set x <-- x + c*y mod p */
    1639                 :            : /* x is not required to be normalized.*/
    1640                 :            : static void
    1641                 :     371609 : Flx_addmul_inplace(GEN gx, GEN gy, ulong c, ulong p)
    1642                 :            : {
    1643                 :            :   long i, lx, ly;
    1644                 :     371609 :   ulong *x=(ulong *)gx;
    1645                 :     371609 :   ulong *y=(ulong *)gy;
    1646         [ +  + ]:     731983 :   if (!c) return;
    1647                 :     360374 :   lx = lg(gx);
    1648                 :     360374 :   ly = lg(gy);
    1649         [ -  + ]:     360374 :   if (lx<ly) pari_err_BUG("lx<ly in Flx_addmul_inplace");
    1650         [ +  + ]:     360374 :   if (SMALL_ULONG(p))
    1651         [ +  + ]:    2212718 :     for (i=2; i<ly;  i++) x[i] = (x[i] + c*y[i]) % p;
    1652                 :            :   else
    1653         [ +  + ]:      10805 :     for (i=2; i<ly;  i++) x[i] = Fl_add(x[i], Fl_mul(c,y[i],p),p);
    1654                 :            : }
    1655                 :            : 
    1656                 :            : /* return a random polynomial in F_q[v], degree < d1 */
    1657                 :            : GEN
    1658                 :       1865 : FqX_rand(long d1, long v, GEN T, GEN p)
    1659                 :            : {
    1660                 :       1865 :   long i, d = d1+2, k = get_FpX_degree(T), w = get_FpX_var(T);
    1661                 :       1865 :   GEN y = cgetg(d,t_POL); y[1] = evalsigne(1) | evalvarn(v);
    1662         [ +  + ]:      12473 :   for (i=2; i<d; i++) gel(y,i) = random_FpX(k, w, p);
    1663                 :       1865 :   (void)normalizepol_lg(y,d); return y;
    1664                 :            : }
    1665                 :            : 
    1666                 :            : #define set_irred(i) { if ((i)>ir) swap(t[i],t[ir]); ir++;}
    1667                 :            : /* assume x1 != 0 */
    1668                 :            : static GEN
    1669                 :      80378 : deg1_Flx(ulong x1, ulong x0, ulong sv)
    1670                 :            : {
    1671                 :      80378 :   return mkvecsmall3(sv, x0, x1);
    1672                 :            : }
    1673                 :            : 
    1674                 :            : long
    1675                 :      30766 : F2x_split_Berlekamp(GEN *t)
    1676                 :            : {
    1677                 :      30766 :   GEN u = *t, a, b, vker;
    1678                 :      30766 :   long lb, d, i, ir, L, la, sv = u[1], du = F2x_degree(u);
    1679                 :            : 
    1680         [ +  + ]:      30766 :   if (du == 1) return 1;
    1681         [ +  + ]:      19791 :   if (du == 2)
    1682                 :            :   {
    1683         [ -  + ]:       5081 :     if (F2x_quad_factortype(u) == 1) /* 0 is a root: shouldn't occur */
    1684                 :            :     {
    1685                 :          0 :       t[0] = mkvecsmall2(sv, 2);
    1686                 :          0 :       t[1] = mkvecsmall2(sv, 3);
    1687                 :          0 :       return 2;
    1688                 :            :     }
    1689                 :       5081 :     return 1;
    1690                 :            :   }
    1691                 :            : 
    1692                 :      14710 :   vker = F2x_Berlekamp_ker(u);
    1693                 :      14710 :   lb = lgcols(vker);
    1694                 :      14710 :   d = lg(vker)-1;
    1695                 :      14710 :   ir = 0;
    1696                 :            :   /* t[i] irreducible for i < ir, still to be treated for i < L */
    1697         [ +  + ]:      18704 :   for (L=1; L<d; )
    1698                 :            :   {
    1699                 :            :     GEN pol;
    1700         [ +  + ]:       3994 :     if (d == 2)
    1701                 :       3772 :       pol = F2v_to_F2x(gel(vker,2), sv);
    1702                 :            :     else
    1703                 :            :     {
    1704                 :        222 :       GEN v = zero_zv(lb);
    1705                 :        222 :       v[1] = du;
    1706                 :        222 :       v[2] = random_Fl(2); /*Assume vker[1]=1*/
    1707         [ +  + ]:        692 :       for (i=2; i<=d; i++)
    1708         [ +  + ]:        470 :         if (random_Fl(2)) F2v_add_inplace(v, gel(vker,i));
    1709                 :        222 :       pol = F2v_to_F2x(v, sv);
    1710                 :            :     }
    1711 [ +  + ][ +  + ]:       8185 :     for (i=ir; i<L && L<d; i++)
    1712                 :            :     {
    1713                 :       4191 :       a = t[i]; la = F2x_degree(a);
    1714 [ +  + ][ -  + ]:       4191 :       if (la == 1) { set_irred(i); }
    1715         [ -  + ]:       4184 :       else if (la == 2)
    1716                 :            :       {
    1717         [ #  # ]:          0 :         if (F2x_quad_factortype(a) == 1) /* 0 is a root: shouldn't occur */
    1718                 :            :         {
    1719                 :          0 :           t[i] = mkvecsmall2(sv, 2);
    1720                 :          0 :           t[L] = mkvecsmall2(sv, 3); L++;
    1721                 :            :         }
    1722         [ #  # ]:          0 :         set_irred(i);
    1723                 :            :       }
    1724                 :            :       else
    1725                 :            :       {
    1726                 :       4184 :         pari_sp av = avma;
    1727                 :            :         long lb;
    1728                 :       4184 :         b = F2x_rem(pol, a);
    1729         [ +  + ]:       4184 :         if (F2x_degree(b) <= 0) { avma=av; continue; }
    1730                 :       3887 :         b = F2x_gcd(a,b); lb = F2x_degree(b);
    1731 [ +  - ][ +  - ]:       3887 :         if (lb && lb < la)
    1732                 :            :         {
    1733                 :       3887 :           t[L] = F2x_div(a,b);
    1734                 :       3887 :           t[i]= b; L++;
    1735                 :            :         }
    1736                 :          0 :         else avma = av;
    1737                 :            :       }
    1738                 :            :     }
    1739                 :            :   }
    1740                 :      30766 :   return d;
    1741                 :            : }
    1742                 :            : 
    1743                 :            : /* p != 2 */
    1744                 :            : long
    1745                 :     164495 : Flx_split_Berlekamp(GEN *t, ulong p)
    1746                 :            : {
    1747                 :     164495 :   GEN u = *t, a,b,vker;
    1748                 :     164495 :   long d, i, ir, L, la, lb, sv = u[1];
    1749                 :     164495 :   long l = lg(u);
    1750                 :            :   ulong po2;
    1751                 :            : 
    1752         [ -  + ]:     164495 :   if (p == 2)
    1753                 :            :   {
    1754                 :          0 :     *t = Flx_to_F2x(*t);
    1755                 :          0 :     d = F2x_split_Berlekamp(t);
    1756         [ #  # ]:          0 :     for (i = 1; i <= d; i++) t[i] = F2x_to_Flx(t[i]);
    1757                 :          0 :     return d;
    1758                 :            :   }
    1759                 :     164495 :   la = degpol(u);
    1760         [ +  + ]:     164495 :   if (la == 1) return 1;
    1761         [ +  + ]:     154709 :   if (la == 2)
    1762                 :            :   {
    1763                 :       7311 :     ulong r = Flx_quad_root(u,p,1);
    1764         [ +  + ]:       7311 :     if (r != p)
    1765                 :            :     {
    1766                 :       3005 :       t[0] = deg1_Flx(1, p - r, sv); r = Flx_otherroot(u,r,p);
    1767                 :       3005 :       t[1] = deg1_Flx(1, p - r, sv);
    1768                 :       3005 :       return 2;
    1769                 :            :     }
    1770                 :       4306 :     return 1;
    1771                 :            :   }
    1772                 :            : 
    1773                 :     147398 :   vker = Flx_Berlekamp_ker(u,p);
    1774                 :     147398 :   vker = Flm_to_FlxV(vker, sv);
    1775                 :     147398 :   d = lg(vker)-1;
    1776                 :     147398 :   po2 = p >> 1; /* (p-1) / 2 */
    1777                 :     147398 :   ir = 0;
    1778                 :            :   /* t[i] irreducible for i < ir, still to be treated for i < L */
    1779         [ +  + ]:     341411 :   for (L=1; L<d; )
    1780                 :            :   {
    1781                 :     194013 :     GEN pol = zero_zv(l-2);
    1782                 :     194013 :     pol[1] = sv;
    1783                 :     194013 :     pol[2] = random_Fl(p); /*Assume vker[1]=1*/
    1784         [ +  + ]:     565622 :     for (i=2; i<=d; i++)
    1785                 :     371609 :       Flx_addmul_inplace(pol, gel(vker,i), random_Fl(p), p);
    1786                 :     194013 :     (void)Flx_renormalize(pol,l-1);
    1787                 :            : 
    1788 [ +  + ][ +  + ]:     461400 :     for (i=ir; i<L && L<d; i++)
    1789                 :            :     {
    1790                 :     267387 :       a = t[i]; la = degpol(a);
    1791 [ +  + ][ +  + ]:     267387 :       if (la == 1) { set_irred(i); }
    1792         [ +  + ]:     240260 :       else if (la == 2)
    1793                 :            :       {
    1794                 :      42065 :         ulong r = Flx_quad_root(a,p,1);
    1795         [ +  + ]:      42065 :         if (r != p)
    1796                 :            :         {
    1797                 :      37184 :           t[i] = deg1_Flx(1, p - r, sv); r = Flx_otherroot(a,r,p);
    1798                 :      37184 :           t[L] = deg1_Flx(1, p - r, sv); L++;
    1799                 :            :         }
    1800         [ +  + ]:      42065 :         set_irred(i);
    1801                 :            :       }
    1802                 :            :       else
    1803                 :            :       {
    1804                 :     198195 :         pari_sp av = avma;
    1805                 :     198195 :         b = Flx_rem(pol, a, p);
    1806         [ +  + ]:     198195 :         if (degpol(b) <= 0) { avma=av; continue; }
    1807                 :     188965 :         b = Flx_Fl_add(Flxq_powu(b,po2, a,p), p-1, p);
    1808                 :     188965 :         b = Flx_gcd(a,b, p); lb = degpol(b);
    1809 [ +  + ][ +  + ]:     188965 :         if (lb && lb < la)
    1810                 :            :         {
    1811                 :     109319 :           b = Flx_normalize(b, p);
    1812                 :     109319 :           t[L] = Flx_div(a,b,p);
    1813                 :     109319 :           t[i]= b; L++;
    1814                 :            :         }
    1815                 :      79646 :         else avma = av;
    1816                 :            :       }
    1817                 :            :     }
    1818                 :            :   }
    1819                 :     164495 :   return d;
    1820                 :            : }
    1821                 :            : 
    1822                 :            : long
    1823                 :        138 : FpX_split_Berlekamp(GEN *t, GEN p)
    1824                 :            : {
    1825                 :        138 :   GEN u = *t, a,b,po2,vker;
    1826                 :        138 :   long d, i, ir, L, la, lb, vu = varn(u);
    1827         [ +  + ]:        138 :   if (lgefint(p) == 3)
    1828                 :            :   {
    1829                 :         95 :     ulong up = p[2];
    1830         [ -  + ]:         95 :     if (up == 2)
    1831                 :            :     {
    1832                 :          0 :       *t = ZX_to_F2x(*t);
    1833                 :          0 :       d = F2x_split_Berlekamp(t);
    1834         [ #  # ]:          0 :       for (i = 0; i < d; i++) t[i] = F2x_to_ZX(t[i]);
    1835                 :            :     }
    1836                 :            :     else
    1837                 :            :     {
    1838                 :         95 :       *t = ZX_to_Flx(*t, up);
    1839                 :         95 :       d = Flx_split_Berlekamp(t, up);
    1840         [ +  + ]:        765 :       for (i = 0; i < d; i++) t[i] = Flx_to_ZX(t[i]);
    1841                 :            :     }
    1842                 :         95 :     return d;
    1843                 :            :   }
    1844                 :         43 :   la = degpol(u);
    1845         [ +  + ]:         43 :   if (la == 1) return 1;
    1846         [ +  + ]:         40 :   if (la == 2)
    1847                 :            :   {
    1848                 :          2 :     GEN r = FpX_quad_root(u,p,1);
    1849         [ -  + ]:          2 :     if (r)
    1850                 :            :     {
    1851                 :          0 :       t[0] = deg1pol_shallow(gen_1, subii(p,r), vu); r = FpX_otherroot(u,r,p);
    1852                 :          0 :       t[1] = deg1pol_shallow(gen_1, subii(p,r), vu);
    1853                 :          0 :       return 2;
    1854                 :            :     }
    1855                 :          2 :     return 1;
    1856                 :            :   }
    1857                 :         38 :   vker = FpX_Berlekamp_ker(u,p);
    1858                 :         38 :   vker = RgM_to_RgXV(vker,vu);
    1859                 :         38 :   d = lg(vker)-1;
    1860                 :         38 :   po2 = shifti(p, -1); /* (p-1) / 2 */
    1861                 :         38 :   ir = 0;
    1862                 :            :   /* t[i] irreducible for i < ir, still to be treated for i < L */
    1863         [ +  + ]:         88 :   for (L=1; L<d; )
    1864                 :            :   {
    1865                 :         50 :     GEN pol = scalar_ZX_shallow(randomi(p), vu);
    1866         [ +  + ]:        216 :     for (i=2; i<=d; i++)
    1867                 :        166 :       pol = ZX_add(pol, ZX_Z_mul(gel(vker,i), randomi(p)));
    1868                 :         50 :     pol = FpX_red(pol,p);
    1869 [ +  + ][ +  + ]:        140 :     for (i=ir; i<L && L<d; i++)
    1870                 :            :     {
    1871                 :         90 :       a = t[i]; la = degpol(a);
    1872 [ +  + ][ +  + ]:         90 :       if (la == 1) { set_irred(i); }
    1873         [ +  + ]:         78 :       else if (la == 2)
    1874                 :            :       {
    1875                 :          9 :         GEN r = FpX_quad_root(a,p,1);
    1876         [ +  - ]:          9 :         if (r)
    1877                 :            :         {
    1878                 :          9 :           t[i] = deg1pol_shallow(gen_1, subii(p,r), vu); r = FpX_otherroot(a,r,p);
    1879                 :          9 :           t[L] = deg1pol_shallow(gen_1, subii(p,r), vu); L++;
    1880                 :            :         }
    1881         [ +  + ]:          9 :         set_irred(i);
    1882                 :            :       }
    1883                 :            :       else
    1884                 :            :       {
    1885                 :         69 :         pari_sp av = avma;
    1886                 :         69 :         b = FpX_rem(pol, a, p);
    1887         [ +  + ]:         69 :         if (degpol(b) <= 0) { avma=av; continue; }
    1888                 :         68 :         b = FpX_Fp_sub_shallow(FpXQ_pow(b,po2, a,p), gen_1, p);
    1889                 :         68 :         b = FpX_gcd(a,b, p); lb = degpol(b);
    1890 [ +  + ][ +  + ]:         68 :         if (lb && lb < la)
    1891                 :            :         {
    1892                 :         33 :           b = FpX_normalize(b, p);
    1893                 :         33 :           t[L] = FpX_div(a,b,p);
    1894                 :         33 :           t[i]= b; L++;
    1895                 :            :         }
    1896                 :         35 :         else avma = av;
    1897                 :            :       }
    1898                 :            :     }
    1899                 :            :   }
    1900                 :        138 :   return d;
    1901                 :            : }
    1902                 :            : 
    1903                 :            : GEN
    1904                 :      10715 : FqX_deriv(GEN f, /*unused*/GEN T, GEN p) {
    1905                 :      10715 :   (void)T; return FpXX_red(RgX_deriv(f), p);
    1906                 :            : }
    1907                 :            : 
    1908                 :            : long
    1909                 :        330 : FqX_is_squarefree(GEN P, GEN T, GEN p)
    1910                 :            : {
    1911                 :        330 :   pari_sp av = avma;
    1912                 :        330 :   GEN z = FqX_gcd(P, FqX_deriv(P, T, p), T, p);
    1913                 :        330 :   avma = av;
    1914                 :        330 :   return degpol(z)==0;
    1915                 :            : }
    1916                 :            : 
    1917                 :            : static GEN
    1918                 :      53263 : F2x_Berlekamp_i(GEN f, long flag)
    1919                 :            : {
    1920                 :      53263 :   long e, nbfact, val, d = F2x_degree(f);
    1921                 :            :   ulong k, j;
    1922                 :            :   GEN y, E, f2, g1, t;
    1923                 :            : 
    1924         [ +  + ]:      53263 :   if (d <= 2) return F2x_factor_deg2(f, d, flag);
    1925                 :            : 
    1926                 :            :   /* to hold factors and exponents */
    1927         [ -  + ]:      31778 :   t = cgetg(d+1, flag? t_VECSMALL: t_VEC);
    1928                 :      31778 :   E = cgetg(d+1,t_VECSMALL);
    1929                 :      31778 :   val = F2x_valrem(f, &f);
    1930                 :      31778 :   e = nbfact = 1;
    1931         [ +  + ]:      31778 :   if (val) {
    1932 [ -  + ][ #  # ]:       6619 :     if (flag == 2 && val > 1) return NULL;
    1933         [ -  + ]:       6619 :     if (flag == 1)
    1934                 :          0 :       t[1] = 1;
    1935                 :            :     else
    1936                 :       6619 :       gel(t,1) = polx_F2x(f[1]);
    1937                 :       6619 :     E[1] = val; nbfact++;
    1938                 :            :   }
    1939                 :            : 
    1940                 :            :   for(;;)
    1941                 :            :   {
    1942                 :      57649 :     f2 = F2x_gcd(f,F2x_deriv(f));
    1943 [ -  + ][ #  # ]:      57649 :     if (flag == 2 && F2x_degree(f2)) return NULL;
    1944         [ +  + ]:      57649 :     g1 = F2x_degree(f2)? F2x_div(f,f2): f; /* squarefree */
    1945                 :      57649 :     k = 0;
    1946         [ +  + ]:      96433 :     while (F2x_degree(g1)>0)
    1947                 :            :     {
    1948                 :            :       GEN u;
    1949         [ +  + ]:      38784 :       k++; if (k%2 == 0) { k++; f2 = F2x_div(f2,g1); }
    1950                 :      38784 :       u = g1; /* deg(u) > 0 */
    1951         [ +  + ]:      38784 :       if (!F2x_degree(f2)) g1 = pol1_F2x(0); /* only its degree (= 0) matters */
    1952                 :            :       else
    1953                 :            :       {
    1954                 :            :         long dg1;
    1955                 :      10044 :         g1 = F2x_gcd(f2,g1);
    1956                 :      10044 :         dg1 = F2x_degree(g1);
    1957         [ +  + ]:      10044 :         if (dg1)
    1958                 :            :         {
    1959                 :       8788 :           f2= F2x_div(f2,g1);
    1960         [ +  + ]:       8788 :           if (F2x_degree(u) == dg1) continue;
    1961                 :        770 :           u = F2x_div( u,g1);
    1962                 :            :         }
    1963                 :            :       }
    1964                 :            :       /* u is square-free (product of irred. of multiplicity e * k) */
    1965                 :      30766 :       gel(t,nbfact) = u;
    1966                 :      30766 :       d = F2x_split_Berlekamp(&gel(t,nbfact));
    1967 [ -  + ][ #  # ]:      30766 :       if (flag == 2 && d != 1) return NULL;
    1968         [ -  + ]:      30766 :       if (flag == 1)
    1969         [ #  # ]:          0 :         for (j=0; j<(ulong)d; j++) t[nbfact+j] = F2x_degree(gel(t,nbfact+j));
    1970         [ +  + ]:      65419 :       for (j=0; j<(ulong)d; j++) E[nbfact+j] = e*k;
    1971                 :      30766 :       nbfact += d;
    1972                 :            :     }
    1973         [ +  + ]:      57649 :     j = F2x_degree(f2); if (!j) break;
    1974                 :      25871 :     e *= 2; f = F2x_deflate(f2, 2);
    1975                 :      25871 :   }
    1976         [ -  + ]:      31778 :   if (flag == 2) return gen_1; /* irreducible */
    1977                 :      31778 :   setlg(t, nbfact);
    1978                 :      31778 :   setlg(E, nbfact); y = mkvec2(t,E);
    1979                 :      53263 :   return flag ? sort_factor(y, (void*)&cmpGuGu, cmp_nodata)
    1980         [ -  + ]:      31778 :               : sort_factor_pol(y, cmpGuGu);
    1981                 :            : }
    1982                 :            : 
    1983                 :            : static GEN
    1984                 :     303346 : Flx_Berlekamp_i(GEN f, ulong p, long flag)
    1985                 :            : {
    1986                 :     303346 :   long e, nbfact, val, d = degpol(f);
    1987                 :            :   ulong k, j;
    1988                 :            :   GEN y, E, f2, g1, t;
    1989                 :            : 
    1990         [ -  + ]:     303346 :   if (p == 2)
    1991                 :            :   {
    1992                 :          0 :     GEN F = F2x_Berlekamp_i(Flx_to_F2x(f),flag);
    1993         [ #  # ]:          0 :     if (flag==0) F2xV_to_FlxV_inplace(gel(F,1));
    1994                 :          0 :     return F;
    1995                 :            :   }
    1996         [ +  + ]:     303346 :   if (d <= 2) return Flx_factor_deg2(f,p,d,flag);
    1997                 :            : 
    1998                 :            :   /* to hold factors and exponents */
    1999         [ -  + ]:     162349 :   t = cgetg(d+1, flag? t_VECSMALL: t_VEC);
    2000                 :     162349 :   E = cgetg(d+1,t_VECSMALL);
    2001                 :     162349 :   val = Flx_valrem(f, &f);
    2002                 :     162349 :   e = nbfact = 1;
    2003         [ +  + ]:     162349 :   if (val) {
    2004 [ -  + ][ #  # ]:       4459 :     if (flag == 2 && val > 1) return NULL;
    2005         [ -  + ]:       4459 :     if (flag == 1)
    2006                 :          0 :       t[1] = 1;
    2007                 :            :     else
    2008                 :       4459 :       gel(t,1) = polx_Flx(f[1]);
    2009                 :       4459 :     E[1] = val; nbfact++;
    2010                 :            :   }
    2011                 :            : 
    2012                 :            :   for(;;)
    2013                 :            :   {
    2014                 :     166207 :     f2 = Flx_gcd(f,Flx_deriv(f,p), p);
    2015 [ -  + ][ #  # ]:     166207 :     if (flag == 2 && degpol(f2)) return NULL;
    2016         [ +  + ]:     166207 :     g1 = degpol(f2)? Flx_div(f,f2,p): f; /* squarefree */
    2017                 :     166207 :     k = 0;
    2018         [ +  + ]:     356741 :     while (degpol(g1)>0)
    2019                 :            :     {
    2020                 :            :       GEN u;
    2021         [ +  + ]:     190534 :       k++; if (k%p == 0) { k++; f2 = Flx_div(f2,g1,p); }
    2022                 :     190534 :       u = g1; /* deg(u) > 0 */
    2023         [ +  + ]:     190534 :       if (!degpol(f2)) g1 = pol1_Flx(0); /* only its degree (= 0) matters */
    2024                 :            :       else
    2025                 :            :       {
    2026                 :      29641 :         g1 = Flx_gcd(f2,g1, p);
    2027         [ +  + ]:      29641 :         if (degpol(g1))
    2028                 :            :         {
    2029                 :      29540 :           f2= Flx_div(f2,g1,p);
    2030         [ +  + ]:      29540 :           if (lg(u) == lg(g1)) continue;
    2031                 :       3406 :           u = Flx_div( u,g1,p);
    2032                 :            :         }
    2033                 :            :       }
    2034                 :            :       /* u is square-free (product of irred. of multiplicity e * k) */
    2035                 :     164400 :       u = Flx_normalize(u,p);
    2036                 :     164400 :       gel(t,nbfact) = u;
    2037                 :     164400 :       d = Flx_split_Berlekamp(&gel(t,nbfact), p);
    2038 [ -  + ][ #  # ]:     164400 :       if (flag == 2 && d != 1) return NULL;
    2039         [ -  + ]:     164400 :       if (flag == 1)
    2040         [ #  # ]:          0 :         for (j=0; j<(ulong)d; j++) t[nbfact+j] = degpol(gel(t,nbfact+j));
    2041         [ +  + ]:     477733 :       for (j=0; j<(ulong)d; j++) E[nbfact+j] = e*k;
    2042                 :     164400 :       nbfact += d;
    2043                 :            :     }
    2044         [ -  + ]:     166207 :     if (!p) break;
    2045         [ +  + ]:     166207 :     j = degpol(f2); if (!j) break;
    2046         [ -  + ]:       3858 :     if (j % p) pari_err_PRIME("factmod",utoi(p));
    2047                 :       3858 :     e *= p; f = Flx_deflate(f2, p);
    2048                 :       3858 :   }
    2049         [ -  + ]:     162349 :   if (flag == 2) return gen_1; /* irreducible */
    2050                 :     162349 :   setlg(t, nbfact);
    2051                 :     162349 :   setlg(E, nbfact); y = mkvec2(t,E);
    2052                 :     303346 :   return flag ? sort_factor(y, (void*)&cmpGuGu, cmp_nodata)
    2053         [ -  + ]:     162349 :               : sort_factor_pol(y, cmpGuGu);
    2054                 :            : }
    2055                 :            : 
    2056                 :            : /* f an FpX or an Flx */
    2057                 :            : static GEN
    2058                 :     353224 : FpX_Berlekamp_i(GEN f, GEN pp, long flag)
    2059                 :            : {
    2060                 :     353224 :   long e, nbfact, val, d = degpol(f);
    2061                 :            :   ulong k, j;
    2062                 :            :   GEN y, E, f2, g1, t;
    2063                 :            : 
    2064         [ +  + ]:     353224 :   if (typ(f) == t_VECSMALL)
    2065                 :            :   {/* lgefint(pp) == 3 */
    2066                 :     353147 :     ulong p = pp[2];
    2067                 :            :     GEN F;
    2068         [ +  + ]:     353147 :     if (p == 2) {
    2069                 :      53263 :       F = F2x_Berlekamp_i(Flx_to_F2x(f), flag);
    2070         [ +  - ]:      53263 :       if (flag==0) F2xV_to_ZXV_inplace(gel(F,1));
    2071                 :            :     } else {
    2072                 :     299884 :       F = Flx_Berlekamp_i(f, p, flag);
    2073         [ +  - ]:     299884 :       if (flag==0) FlxV_to_ZXV_inplace(gel(F,1));
    2074                 :            :     }
    2075                 :     353147 :     return F;
    2076                 :            :   }
    2077                 :            :   /* p is large (and odd) */
    2078         [ +  + ]:         77 :   if (d <= 2) return FpX_factor_deg2(f,pp,d,flag);
    2079                 :            : 
    2080                 :            :   /* to hold factors and exponents */
    2081         [ -  + ]:         40 :   t = cgetg(d+1, flag? t_VECSMALL: t_VEC);
    2082                 :         40 :   E = cgetg(d+1,t_VECSMALL);
    2083                 :         40 :   val = ZX_valrem(f, &f);
    2084                 :         40 :   e = nbfact = 1;
    2085         [ +  + ]:         40 :   if (val) {
    2086 [ -  + ][ #  # ]:          2 :     if (flag == 2 && val > 1) return NULL;
    2087         [ -  + ]:          2 :     if (flag == 1)
    2088                 :          0 :       t[1] = 1;
    2089                 :            :     else
    2090                 :          2 :       gel(t,1) = pol_x(varn(f));
    2091                 :          2 :     E[1] = val; nbfact++;
    2092                 :            :   }
    2093                 :            : 
    2094                 :         40 :   f2 = FpX_gcd(f,ZX_deriv(f), pp);
    2095 [ -  + ][ #  # ]:         40 :   if (flag == 2 && degpol(f2)) return NULL;
    2096         [ +  + ]:         40 :   g1 = degpol(f2)? FpX_div(f,f2,pp): f; /* squarefree */
    2097                 :         40 :   k = 0;
    2098         [ +  + ]:         84 :   while (degpol(g1)>0)
    2099                 :            :   {
    2100                 :            :     GEN u;
    2101                 :         44 :     k++;
    2102                 :         44 :     u = g1; /* deg(u) > 0 */
    2103         [ +  + ]:         44 :     if (!degpol(f2)) g1 = pol_1(0); /* only its degree (= 0) matters */
    2104                 :            :     else
    2105                 :            :     {
    2106                 :          4 :       g1 = FpX_gcd(f2,g1, pp);
    2107         [ +  - ]:          4 :       if (degpol(g1))
    2108                 :            :       {
    2109                 :          4 :         f2= FpX_div(f2,g1,pp);
    2110         [ +  + ]:          4 :         if (lg(u) == lg(g1)) continue;
    2111                 :          3 :         u = FpX_div( u,g1,pp);
    2112                 :            :       }
    2113                 :            :     }
    2114                 :            :     /* u is square-free (product of irred. of multiplicity e * k) */
    2115                 :         43 :     u = FpX_normalize(u,pp);
    2116                 :         43 :     gel(t,nbfact) = u;
    2117                 :         43 :     d = FpX_split_Berlekamp(&gel(t,nbfact), pp);
    2118 [ -  + ][ #  # ]:         43 :     if (flag == 2 && d != 1) return NULL;
    2119         [ -  + ]:         43 :     if (flag == 1)
    2120         [ #  # ]:          0 :       for (j=0; j<(ulong)d; j++) t[nbfact+j] = degpol(gel(t,nbfact+j));
    2121         [ +  + ]:        128 :     for (j=0; j<(ulong)d; j++) E[nbfact+j] = e*k;
    2122                 :         43 :     nbfact += d;
    2123                 :            :   }
    2124         [ -  + ]:         40 :   if (flag == 2) return gen_1; /* irreducible */
    2125                 :         40 :   setlg(t, nbfact);
    2126                 :         40 :   setlg(E, nbfact); y = mkvec2(t,E);
    2127                 :     353224 :   return flag ? sort_factor(y, (void*)&cmpGuGu, cmp_nodata)
    2128         [ -  + ]:         40 :               : sort_factor_pol(y, cmpii);
    2129                 :            : }
    2130                 :            : GEN
    2131                 :     353179 : FpX_factor(GEN f, GEN p)
    2132                 :            : {
    2133                 :     353179 :   pari_sp av = avma;
    2134                 :     353179 :   ZX_factmod_init(&f, p);
    2135                 :     353179 :   return gerepilecopy(av, FpX_Berlekamp_i(f, p, 0));
    2136                 :            : }
    2137                 :            : GEN
    2138                 :      12289 : Flx_factor(GEN f, ulong p)
    2139                 :            : {
    2140                 :      12289 :   pari_sp av = avma;
    2141         [ +  + ]:      12289 :   GEN F = (degpol(f)>log2(p))? Flx_factcantor_i(f,p,0): Flx_Berlekamp_i(f,p,0);
    2142                 :      12289 :   return gerepilecopy(av, F);
    2143                 :            : }
    2144                 :            : GEN
    2145                 :          0 : F2x_factor(GEN f)
    2146                 :            : {
    2147                 :          0 :   pari_sp av = avma;
    2148                 :          0 :   return gerepilecopy(av, F2x_Berlekamp_i(f, 0));
    2149                 :            : }
    2150                 :            : 
    2151                 :            : GEN
    2152                 :         65 : factormod0(GEN f, GEN p, long flag)
    2153                 :            : {
    2154      [ +  +  - ]:         65 :   switch(flag)
    2155                 :            :   {
    2156                 :         45 :     case 0: return factmod(f,p);
    2157                 :         20 :     case 1: return simplefactmod(f,p);
    2158                 :          0 :     default: pari_err_FLAG("factormod");
    2159                 :            :   }
    2160                 :         65 :   return NULL; /* not reached */
    2161                 :            : }
    2162                 :            : 
    2163                 :            : /*******************************************************************/
    2164                 :            : /*                                                                 */
    2165                 :            : /*                     FACTORIZATION IN F_q                        */
    2166                 :            : /*                                                                 */
    2167                 :            : /*******************************************************************/
    2168                 :            : static GEN
    2169                 :       6830 : to_Fq(GEN x, GEN T, GEN p)
    2170                 :            : {
    2171                 :       6830 :   long i, lx, tx = typ(x);
    2172                 :            :   GEN y;
    2173                 :            : 
    2174         [ +  + ]:       6830 :   if (tx == t_INT)
    2175                 :       3145 :     y = mkintmod(x,p);
    2176                 :            :   else
    2177                 :            :   {
    2178         [ -  + ]:       3685 :     if (tx != t_POL) pari_err_TYPE("to_Fq",x);
    2179                 :       3685 :     lx = lg(x);
    2180                 :       3685 :     y = cgetg(lx,t_POL); y[1] = x[1];
    2181         [ +  + ]:     152900 :     for (i=2; i<lx; i++) gel(y,i) = mkintmod(gel(x,i), p);
    2182                 :            :   }
    2183                 :       6830 :   return mkpolmod(y, T);
    2184                 :            : }
    2185                 :            : 
    2186                 :            : static GEN
    2187                 :       3125 : to_Fq_pol(GEN x, GEN T, GEN p)
    2188                 :            : {
    2189                 :       3125 :   long i, lx = lg(x);
    2190         [ +  + ]:       9385 :   for (i=2; i<lx; i++) gel(x,i) = to_Fq(gel(x,i),T,p);
    2191                 :       3125 :   return x;
    2192                 :            : }
    2193                 :            : 
    2194                 :            : static GEN
    2195                 :        180 : to_Fq_fact(GEN P, GEN E, GEN T, GEN p, pari_sp av)
    2196                 :            : {
    2197                 :            :   GEN y, u, v;
    2198                 :        180 :   long j, l = lg(P), nbf = lg(P);
    2199                 :            : 
    2200                 :        180 :   u = cgetg(nbf,t_COL);
    2201                 :        180 :   v = cgetg(nbf,t_COL);
    2202         [ +  + ]:       3305 :   for (j=1; j<l; j++)
    2203                 :            :   {
    2204                 :       3125 :     gel(u,j) = simplify_shallow(gel(P,j)); /* may contain pols of degree 0 */
    2205                 :       3125 :     gel(v,j) = utoi(uel(E,j));
    2206                 :            :   }
    2207                 :        180 :   y = gerepilecopy(av, mkmat2(u, v));
    2208                 :        180 :   u = gel(y,1);
    2209                 :        180 :   p = icopy(p);
    2210                 :        180 :   T = FpX_to_mod(T, p);
    2211         [ +  + ]:       3305 :   for (j=1; j<nbf; j++) gel(u,j) = to_Fq_pol(gel(u,j), T,p);
    2212                 :        180 :   return y;
    2213                 :            : }
    2214                 :            : static GEN
    2215                 :         25 : to_FqC(GEN P, GEN T, GEN p, pari_sp av)
    2216                 :            : {
    2217                 :            :   GEN u;
    2218                 :         25 :   long j, l = lg(P), nbf = lg(P);
    2219                 :            : 
    2220                 :         25 :   u = cgetg(nbf,t_COL);
    2221         [ +  + ]:        595 :   for (j=1; j<l; j++)
    2222                 :        570 :     gel(u,j) = simplify_shallow(gel(P,j)); /* may contain pols of degree 0 */
    2223                 :         25 :   u = gerepilecopy(av, u);
    2224                 :         25 :   p = icopy(p);
    2225                 :         25 :   T = FpX_to_mod(T, p);
    2226         [ +  + ]:        595 :   for (j=1; j<nbf; j++) gel(u,j) = to_Fq(gel(u,j), T,p);
    2227                 :         25 :   return u;
    2228                 :            : }
    2229                 :            : 
    2230                 :            : GEN
    2231                 :       1881 : FlxqXQ_halfFrobenius(GEN a, GEN S, GEN T, ulong p)
    2232                 :            : {
    2233                 :       1881 :   long vT = get_Flx_var(T);
    2234                 :       1881 :   GEN xp = Flxq_powu(polx_Flx(vT), p, T, p);
    2235                 :       1881 :   GEN Xp = FlxqXQ_pow(polx_FlxX(varn(S), vT), utoi(p), S, T, p);
    2236                 :       1881 :   GEN ap2 = FlxqXQ_pow(a,utoi(p>>1), S, T, p);
    2237                 :       1881 :   GEN V = FlxqXQV_autsum(mkvec3(xp, Xp, ap2), get_Flx_degree(T), S, T, p);
    2238                 :       1881 :   return gel(V,3);
    2239                 :            : }
    2240                 :            : 
    2241                 :            : GEN
    2242                 :       1840 : FpXQXQ_halfFrobenius(GEN a, GEN S, GEN T, GEN p)
    2243                 :            : {
    2244         [ +  + ]:       1840 :   if (lgefint(p)==3)
    2245                 :            :   {
    2246                 :       1791 :     ulong pp = p[2];
    2247                 :       1791 :     long v = get_FpX_var(T);
    2248                 :       1791 :     GEN Tp = ZXT_to_FlxT(T,pp), Sp = ZXX_to_FlxX(S, pp, v);
    2249                 :       1791 :     return FlxX_to_ZXX(FlxqXQ_halfFrobenius(ZXX_to_FlxX(a,pp,v),Sp,Tp,pp));
    2250                 :            :   }
    2251                 :            :   else
    2252                 :            :   {
    2253                 :         49 :     GEN xp = FpXQ_pow(pol_x(get_FpX_var(T)), p, T, p);
    2254                 :         49 :     GEN Xp = FpXQXQ_pow(pol_x(varn(S)), p, S, T, p);
    2255                 :         49 :     GEN ap2 = FpXQXQ_pow(a,shifti(p,-1), S, T, p);
    2256                 :         49 :     GEN V = FpXQXQV_autsum(mkvec3(xp,Xp,ap2), get_FpX_degree(T), S, T, p);
    2257                 :       1840 :     return gel(V,3);
    2258                 :            :   }
    2259                 :            : }
    2260                 :            : 
    2261                 :            : GEN
    2262                 :        755 : FlxqX_Frobenius(GEN S, GEN T, ulong p)
    2263                 :            : {
    2264                 :        755 :   pari_sp av = avma;
    2265                 :        755 :   long n = get_Flx_degree(T), vT = get_Flx_var(T);
    2266                 :        755 :   GEN X  = polx_FlxX(varn(S), vT);
    2267                 :        755 :   GEN xp = Flxq_powu(polx_Flx(vT), p, T, p);
    2268                 :        755 :   GEN Xp = FlxqXQ_pow(X, utoi(p), S, T, p);
    2269                 :        755 :   GEN Xq = gel(FlxqXQV_autpow(mkvec2(xp,Xp), n, S, T, p), 2);
    2270                 :        755 :   return gerepilecopy(av, Xq);
    2271                 :            : }
    2272                 :            : 
    2273                 :            : GEN
    2274                 :         87 : FpXQX_Frobenius(GEN S, GEN T, GEN p)
    2275                 :            : {
    2276                 :         87 :   pari_sp av = avma;
    2277                 :         87 :   long n = get_FpX_degree(T);
    2278                 :         87 :   GEN X  = pol_x(varn(S));
    2279                 :         87 :   GEN xp = FpXQ_pow(pol_x(get_FpX_var(T)), p, T, p);
    2280                 :         87 :   GEN Xp = FpXQXQ_pow(X, p, S, T, p);
    2281                 :         87 :   GEN Xq = gel(FpXQXQV_autpow(mkvec2(xp,Xp), n, S, T, p), 2);
    2282                 :         87 :   return gerepilecopy(av, Xq);
    2283                 :            : }
    2284                 :            : 
    2285                 :            : static GEN
    2286                 :        562 : FqX_Frobenius_powers(GEN S, GEN T, GEN p)
    2287                 :            : {
    2288                 :        562 :   long N = degpol(S);
    2289         [ +  - ]:        562 :   if (lgefint(p)==3)
    2290                 :            :   {
    2291                 :        562 :     ulong pp = p[2];
    2292                 :        562 :     GEN Tp = ZXT_to_FlxT(T, pp), Sp = ZXX_to_FlxX(S, pp, get_FpX_var(T));
    2293                 :        562 :     GEN Xq = FlxqX_Frobenius(Sp, Tp, pp);
    2294                 :        562 :     return FlxqXQ_powers(Xq, N-1, Sp, Tp, pp);
    2295                 :            :   } else
    2296                 :            :   {
    2297                 :          0 :     GEN Xq = FpXQX_Frobenius(S, T, p);
    2298                 :        562 :     return FpXQXQ_powers(Xq, N-1, S, T, p);
    2299                 :            :   }
    2300                 :            : }
    2301                 :            : 
    2302                 :            : static GEN
    2303                 :       1463 : FqX_Frobenius_eval(GEN x, GEN V, GEN S, GEN T, GEN p)
    2304                 :            : {
    2305         [ +  - ]:       1463 :   if (lgefint(p)==3)
    2306                 :            :   {
    2307                 :       1463 :     ulong pp = p[2];
    2308                 :       1463 :     long v = get_FpX_var(T);
    2309                 :       1463 :     GEN Tp = ZXT_to_FlxT(T, pp), Sp = ZXX_to_FlxX(S, pp, v);
    2310                 :       1463 :     GEN xp = ZXX_to_FlxX(x, pp, v);
    2311                 :       1463 :     return FlxX_to_ZXX(FlxqX_FlxqXQV_eval(xp, V, Sp, Tp, pp));
    2312                 :            :   }
    2313                 :            :   else
    2314                 :       1463 :     return FpXQX_FpXQXQV_eval(x, V, S, T, p);
    2315                 :            : }
    2316                 :            : 
    2317                 :            : static GEN
    2318                 :          5 : FpXQX_split_part(GEN f, GEN T, GEN p)
    2319                 :            : {
    2320                 :          5 :   long n = degpol(f);
    2321                 :          5 :   GEN z, X = pol_x(varn(f));
    2322         [ -  + ]:          5 :   if (n <= 1) return f;
    2323                 :          5 :   f = FpXQX_red(f, T, p);
    2324                 :          5 :   z = FpXQX_Frobenius(f, T, p);
    2325                 :          5 :   z = FpXX_sub(z, X , p);
    2326                 :          5 :   return FpXQX_gcd(z, f, T, p);
    2327                 :            : }
    2328                 :            : 
    2329                 :            : static GEN
    2330                 :         25 : FlxqX_split_part(GEN f, GEN T, ulong p)
    2331                 :            : {
    2332                 :         25 :   long n = degpol(f);
    2333                 :         25 :   GEN z, Xq, X = polx_FlxX(varn(f),get_Flx_var(T));
    2334         [ -  + ]:         25 :   if (n <= 1) return f;
    2335                 :         25 :   f = FlxqX_red(f, T, p);
    2336                 :         25 :   Xq = FlxqX_Frobenius(f, T, p);
    2337                 :         25 :   z = FlxX_sub(Xq, X , p);
    2338                 :         25 :   return FlxqX_gcd(z, f, T, p);
    2339                 :            : }
    2340                 :            : 
    2341                 :            : long
    2342                 :         15 : FpXQX_nbroots(GEN f, GEN T, GEN p)
    2343                 :            : {
    2344                 :         15 :   pari_sp av = avma;
    2345                 :            :   GEN z;
    2346         [ +  + ]:         15 :   if(lgefint(p)==3)
    2347                 :            :   {
    2348                 :         10 :     ulong pp=p[2];
    2349                 :         10 :     z = FlxqX_split_part(ZXX_to_FlxX(f,pp,varn(T)),ZXT_to_FlxT(T,pp),pp);
    2350                 :            :   }
    2351                 :            :   else
    2352                 :          5 :     z = FpXQX_split_part(f, T, p);
    2353                 :         15 :   avma = av; return degpol(z);
    2354                 :            : }
    2355                 :            : 
    2356                 :            : long
    2357                 :        830 : FqX_nbroots(GEN f, GEN T, GEN p)
    2358         [ +  + ]:        830 : { return T ? FpXQX_nbroots(f, T, p): FpX_nbroots(f, p); }
    2359                 :            : 
    2360                 :            : long
    2361                 :         15 : FlxqX_nbroots(GEN f, GEN T, ulong p)
    2362                 :            : {
    2363                 :         15 :   pari_sp av = avma;
    2364                 :         15 :   GEN z = FlxqX_split_part(f, T, p);
    2365                 :         15 :   avma = av; return degpol(z);
    2366                 :            : }
    2367                 :            : 
    2368                 :            : GEN
    2369                 :         48 : FlxqX_Berlekamp_ker(GEN u, GEN T, ulong p)
    2370                 :            : {
    2371                 :         48 :   pari_sp ltop=avma;
    2372                 :         48 :   long j,N = degpol(u);
    2373                 :         48 :   GEN Xq = FlxqX_Frobenius(u,T,p);
    2374                 :         48 :   GEN Q  = FlxqXQ_matrix_pow(Xq,N,N,u,T,p);
    2375         [ +  + ]:        210 :   for (j=1; j<=N; j++)
    2376                 :        162 :     gcoeff(Q,j,j) = Flx_Fl_add(gcoeff(Q,j,j), p-1, p);
    2377                 :         48 :   return gerepileupto(ltop, FlxqM_ker(Q,T,p));
    2378                 :            : }
    2379                 :            : 
    2380                 :            : GEN
    2381                 :         50 : FpXQX_Berlekamp_ker(GEN u, GEN T, GEN p)
    2382                 :            : {
    2383         [ +  + ]:         50 :   if (lgefint(p)==3)
    2384                 :            :   {
    2385                 :         48 :     ulong pp=p[2];
    2386                 :         48 :     long v = get_FpX_var(T);
    2387                 :         48 :     GEN Tp = ZXT_to_FlxT(T,pp), Sp = ZXX_to_FlxX(u,pp,v);
    2388                 :         48 :     return FlxM_to_ZXM(FlxqX_Berlekamp_ker(Sp, Tp, pp));
    2389                 :            :   } else
    2390                 :            :   {
    2391                 :          2 :     pari_sp ltop=avma;
    2392                 :          2 :     long j,N = degpol(u);
    2393                 :          2 :     GEN Xq = FpXQX_Frobenius(u,T,p);
    2394                 :          2 :     GEN Q  = FpXQXQ_matrix_pow(Xq,N,N,u,T,p);
    2395         [ +  + ]:         10 :     for (j=1; j<=N; j++)
    2396                 :          8 :       gcoeff(Q,j,j) = Fq_sub(gcoeff(Q,j,j), gen_1, T, p);
    2397                 :         50 :     return gerepileupto(ltop, FqM_ker(Q,T,p));
    2398                 :            :   }
    2399                 :            : }
    2400                 :            : 
    2401                 :            : long
    2402                 :          0 : FpXQX_nbfact(GEN u, GEN T, GEN p)
    2403                 :            : {
    2404                 :          0 :   pari_sp av = avma;
    2405                 :          0 :   GEN vker = FpXQX_Berlekamp_ker(u, T, p);
    2406                 :          0 :   avma = av; return lg(vker)-1;
    2407                 :            : }
    2408                 :            : 
    2409                 :            : long
    2410                 :          0 : FqX_nbfact(GEN u, GEN T, GEN p)
    2411                 :            : {
    2412         [ #  # ]:          0 :   return T ? FpX_nbfact(u, p): FpXQX_nbfact(u, T, p);
    2413                 :            : }
    2414                 :            : 
    2415                 :            : long
    2416                 :         50 : FqX_split_Berlekamp(GEN *t, GEN T, GEN p)
    2417                 :            : {
    2418                 :         50 :   GEN u = *t, a,b,vker,pol;
    2419                 :         50 :   long vu = varn(u), vT = varn(T), dT = degpol(T);
    2420                 :            :   long d, i, ir, L, la, lb;
    2421                 :         50 :   T = FpX_get_red(T, p);
    2422                 :         50 :   vker = FpXQX_Berlekamp_ker(u,T,p);
    2423                 :         50 :   vker = RgM_to_RgXV(vker,vu);
    2424                 :         50 :   d = lg(vker)-1;
    2425                 :         50 :   ir = 0;
    2426                 :            :   /* t[i] irreducible for i < ir, still to be treated for i < L */
    2427         [ +  + ]:        131 :   for (L=1; L<d; )
    2428                 :            :   {
    2429                 :         81 :     pol= scalarpol(random_FpX(dT,vT,p),vu);
    2430         [ +  + ]:        225 :     for (i=2; i<=d; i++)
    2431                 :        144 :       pol = FqX_add(pol, FqX_Fq_mul(gel(vker,i),
    2432                 :            :                                     random_FpX(dT,vT,p), T, p), T, p);
    2433                 :         81 :     pol = FpXQX_red(pol,T,p);
    2434 [ +  + ][ +  + ]:        201 :     for (i=ir; i<L && L<d; i++)
    2435                 :            :     {
    2436                 :        120 :       a = t[i]; la = degpol(a);
    2437 [ +  + ][ +  + ]:        120 :       if (la == 1) { set_irred(i); }
    2438                 :            :       else
    2439                 :            :       {
    2440                 :        115 :         pari_sp av = avma;
    2441                 :        115 :         b = FqX_rem(pol, a, T,p);
    2442         [ +  + ]:        115 :         if (degpol(b)<=0) { avma=av; continue; }
    2443                 :         97 :         b = FpXQXQ_halfFrobenius(b, a,T,p);
    2444         [ +  + ]:         97 :         if (degpol(b)<=0) { avma=av; continue; }
    2445                 :         60 :         gel(b,2) = Fq_sub(gel(b,2), gen_1,T,p);
    2446                 :         60 :         b = FqX_gcd(a,b, T,p); lb = degpol(b);
    2447 [ +  - ][ +  - ]:         60 :         if (lb && lb < la)
    2448                 :            :         {
    2449                 :         60 :           b = FqX_normalize(b, T,p);
    2450                 :         60 :           t[L] = FqX_div(a,b,T,p);
    2451                 :         60 :           t[i]= b; L++;
    2452                 :            :         }
    2453                 :          0 :         else avma = av;
    2454                 :            :       }
    2455                 :            :     }
    2456                 :            :   }
    2457                 :         50 :   return d;
    2458                 :            : }
    2459                 :            : 
    2460                 :            : /* split into r factors of degree d */
    2461                 :            : static void
    2462                 :       2820 : FqX_split(GEN *t, long d, GEN q, GEN S, GEN T, GEN p)
    2463                 :            : {
    2464                 :       2820 :   GEN u = *t;
    2465                 :       2820 :   long l, v, is2, cnt, dt = degpol(u), dT = degpol(T);
    2466                 :            :   pari_sp av;
    2467                 :            :   pari_timer ti;
    2468                 :            :   GEN w,w0;
    2469                 :            : 
    2470         [ +  + ]:       5640 :   if (dt == d) return;
    2471                 :       1263 :   v = varn(*t);
    2472         [ -  + ]:       1263 :   if (DEBUGLEVEL > 6) timer_start(&ti);
    2473                 :       1263 :   av = avma; is2 = equaliu(p, 2);
    2474                 :       1263 :   for(cnt = 1;;cnt++, avma = av)
    2475                 :            :   { /* splits *t with probability ~ 1 - 2^(1-r) */
    2476                 :       1865 :     w = w0 = FqX_rand(dt,v, T,p);
    2477         [ +  + ]:       1865 :     if (degpol(w) <= 0) continue;
    2478         [ +  + ]:       2634 :     for (l=1; l<d; l++) /* sum_{0<i<d} w^(q^i), result in (F_q)^r */
    2479                 :        821 :       w = RgX_add(w0, FqX_Frobenius_eval(w, S, u, T, p));
    2480                 :       1813 :     w = FpXQX_red(w, T,p);
    2481         [ +  + ]:       1813 :     if (is2)
    2482                 :            :     {
    2483                 :        115 :       w0 = w;
    2484         [ +  + ]:        230 :       for (l=1; l<dT; l++) /* sum_{0<i<k} w^(2^i), result in (F_2)^r */
    2485                 :            :       {
    2486                 :        115 :         w = FqX_rem(FqX_sqr(w,T,p), *t, T,p);
    2487                 :        115 :         w = FpXX_red(RgX_add(w0,w), p);
    2488                 :            :       }
    2489                 :            :     }
    2490                 :            :     else
    2491                 :            :     {
    2492                 :       1698 :       w = FpXQXQ_halfFrobenius(w, *t, T, p);
    2493                 :            :       /* w in {-1,0,1}^r */
    2494         [ +  + ]:       1698 :       if (degpol(w) <= 0) continue;
    2495                 :       1181 :       gel(w,2) = gadd(gel(w,2), gen_1);
    2496                 :            :     }
    2497                 :       1296 :     w = FqX_gcd(*t,w, T,p); l = degpol(w);
    2498 [ +  + ][ +  + ]:       1296 :     if (l && l != dt) break;
    2499                 :        602 :   }
    2500                 :       1263 :   w = gerepileupto(av,FqX_normalize(w,T,p));
    2501         [ -  + ]:       1263 :   if (DEBUGLEVEL > 6)
    2502                 :          0 :     err_printf("[FqX_split] splitting time: %ld (%ld trials)\n",
    2503                 :            :                timer_delay(&ti),cnt);
    2504                 :       1263 :   l /= d; t[l] = FqX_div(*t,w, T,p); *t = w;
    2505                 :       1263 :   FqX_split(t+l,d,q,S,T,p);
    2506                 :       2820 :   FqX_split(t  ,d,q,S,T,p);
    2507                 :            : }
    2508                 :            : 
    2509                 :            : /*******************************************************************/
    2510                 :            : /*                                                                 */
    2511                 :            : /*                  FACTOR USING TRAGER'S TRICK                    */
    2512                 :            : /*                                                                 */
    2513                 :            : /*******************************************************************/
    2514                 :            : static GEN
    2515                 :         45 : FqX_frob_deflate(GEN f, GEN T, GEN p)
    2516                 :            : {
    2517                 :         45 :   GEN F = RgX_deflate(f, itos(p)), frobinv = powiu(p, degpol(T)-1);
    2518                 :         45 :   long i, l = lg(F);
    2519         [ +  + ]:        170 :   for (i=2; i<l; i++) gel(F,i) = Fq_pow(gel(F,i), frobinv, T,p);
    2520                 :         45 :   return F;
    2521                 :            : }
    2522                 :            : /* Factor _sqfree_ polynomial a on the finite field Fp[X]/(T). Assumes
    2523                 :            :  * varncmp (varn(T), varn(A)) > 0 */
    2524                 :            : static GEN
    2525                 :         25 : FqX_split_Trager(GEN A, GEN T, GEN p)
    2526                 :            : {
    2527                 :         25 :   GEN c = NULL, P, u, fa, n;
    2528                 :            :   long lx, i, k;
    2529                 :            : 
    2530                 :         25 :   u = A;
    2531                 :         25 :   n = NULL;
    2532         [ +  + ]:         50 :   for (k = 0; cmpui(k, p) < 0; k++)
    2533                 :            :   {
    2534                 :            :     GEN U;
    2535                 :         40 :     c = deg1pol_shallow(stoi(k) , gen_0, varn(T));
    2536                 :         40 :     U = FqX_translate(u, c, T, p);
    2537                 :         40 :     n = FpX_FpXY_resultant(T, U, p);
    2538         [ +  + ]:         40 :     if (FpX_is_squarefree(n, p)) break;
    2539                 :         25 :     n = NULL;
    2540                 :            :   }
    2541         [ +  + ]:         25 :   if (!n) return NULL;
    2542         [ -  + ]:         15 :   if (DEBUGLEVEL>4) err_printf("FqX_split_Trager: choosing k = %ld\n",k);
    2543                 :            :   /* n guaranteed to be squarefree */
    2544                 :         15 :   fa = FpX_factor(n, p); fa = gel(fa,1); lx = lg(fa);
    2545         [ -  + ]:         15 :   if (lx == 2) return mkcol(A); /* P^k, P irreducible */
    2546                 :            : 
    2547                 :         15 :   P = cgetg(lx,t_COL);
    2548                 :         15 :   c = FpX_neg(c,p);
    2549         [ +  + ]:         30 :   for (i=lx-1; i>1; i--)
    2550                 :            :   {
    2551                 :         15 :     GEN F = FqX_translate(gel(fa,i), c, T, p);
    2552                 :         15 :     F = FqX_normalize(FqX_gcd(u, F, T, p), T, p);
    2553 [ +  - ][ -  + ]:         15 :     if (typ(F) != t_POL || degpol(F) == 0)
    2554                 :          0 :       pari_err_IRREDPOL("FqX_split_Trager [modulus]",T);
    2555                 :         15 :     u = FqX_div(u, F, T, p);
    2556                 :         15 :     gel(P,i) = F;
    2557                 :            :   }
    2558                 :         25 :   gel(P,1) = u; return P;
    2559                 :            : }
    2560                 :            : 
    2561                 :            : static long
    2562                 :       2425 : isabsolutepol(GEN f)
    2563                 :            : {
    2564                 :       2425 :   long i, l = lg(f);
    2565         [ +  + ]:       8469 :   for(i=2; i<l; i++)
    2566                 :            :   {
    2567                 :       7865 :     GEN c = gel(f,i);
    2568 [ +  + ][ +  + ]:       7865 :     if (typ(c) == t_POL && degpol(c) > 0) return 0;
    2569                 :            :   }
    2570                 :       2425 :   return 1;
    2571                 :            : }
    2572                 :            : 
    2573                 :            : static void
    2574                 :        524 : add(GEN z, GEN g, long d) { vectrunc_append(z, mkvec2(utoipos(d), g)); }
    2575                 :            : /* return number of roots of u; assume deg u >= 0 */
    2576                 :            : long
    2577                 :        377 : FqX_split_deg1(GEN *pz, GEN u, GEN T, GEN p)
    2578                 :            : {
    2579                 :        377 :   long dg, N = degpol(u);
    2580                 :        377 :   GEN v, S, g, X, z = vectrunc_init(N+1);
    2581                 :            : 
    2582                 :        377 :   *pz = z;
    2583         [ -  + ]:        377 :   if (N == 0) return 0;
    2584         [ +  + ]:        377 :   if (N == 1) return 1;
    2585                 :        352 :   v = X = pol_x(varn(u));
    2586                 :        352 :   S = FqX_Frobenius_powers(u, T, p);
    2587                 :        352 :   vectrunc_append(z, S);
    2588                 :        352 :   v = FqX_Frobenius_eval(v, S, u, T, p);
    2589                 :        352 :   g = FqX_gcd(FpXX_sub(v,X,p),u, T,p);
    2590                 :        352 :   dg = degpol(g);
    2591         [ +  + ]:        352 :   if (dg > 0) add(z, FqX_normalize(g,T,p), dg);
    2592                 :        377 :   return dg;
    2593                 :            : }
    2594                 :            : 
    2595                 :            : /* return number of factors; z not properly initialized if deg(u) <= 1 */
    2596                 :            : long
    2597                 :        200 : FqX_split_by_degree(GEN *pz, GEN u, GEN T, GEN p)
    2598                 :            : {
    2599                 :        200 :   long nb = 0, d, dg, N = degpol(u);
    2600                 :        200 :   GEN v, S, g, X, z = vectrunc_init(N+1);
    2601                 :            : 
    2602                 :        200 :   *pz = z;
    2603         [ -  + ]:        200 :   if (N <= 1) return 1;
    2604                 :        200 :   v = X = pol_x(varn(u));
    2605                 :        200 :   S = FqX_Frobenius_powers(u, T, p);
    2606                 :        200 :   vectrunc_append(z, S);
    2607         [ +  + ]:        480 :   for (d=1; d <= N>>1; d++)
    2608                 :            :   {
    2609                 :        280 :     v = FqX_Frobenius_eval(v, S, u, T, p);
    2610                 :        280 :     g = FqX_gcd(FpXX_sub(v,X,p),u, T,p);
    2611         [ +  + ]:        280 :     dg = degpol(g); if (dg <= 0) continue;
    2612                 :            :     /* all factors of g have degree d */
    2613                 :        225 :     add(z, FqX_normalize(g, T,p), dg / d); nb += dg / d;
    2614                 :        225 :     N -= dg;
    2615         [ +  + ]:        225 :     if (N)
    2616                 :            :     {
    2617                 :         25 :       u = FqX_div(u,g, T,p);
    2618                 :         25 :       v = FqX_rem(v,u, T,p);
    2619                 :            :     }
    2620                 :            :   }
    2621         [ -  + ]:        200 :   if (N) { add(z, FqX_normalize(u, T,p), 1); nb++; }
    2622                 :        200 :   return nb;
    2623                 :            : }
    2624                 :            : 
    2625                 :            : /* see roots_from_deg1() */
    2626                 :            : static GEN
    2627                 :        219 : FqXC_roots_from_deg1(GEN x, GEN T, GEN p)
    2628                 :            : {
    2629                 :        219 :   long i,l = lg(x);
    2630                 :        219 :   GEN r = cgetg(l,t_COL);
    2631         [ +  + ]:        536 :   for (i=1; i<l; i++) { GEN P = gel(x,i); gel(r,i) = Fq_neg(gel(P,2), T, p); }
    2632                 :        219 :   return r;
    2633                 :            : }
    2634                 :            : 
    2635                 :            : static GEN
    2636                 :        284 : FqX_split_equal(GEN L, GEN S, GEN T, GEN p)
    2637                 :            : {
    2638                 :        284 :   long n = itos(gel(L,1));
    2639                 :        284 :   GEN u = gel(L,2), z = cgetg(n + 1, t_COL);
    2640                 :        284 :   gel(z,1) = u;
    2641                 :        284 :   FqX_split((GEN*)(z+1), degpol(u) / n, powiu(p, degpol(T)), S, T, p);
    2642                 :        284 :   return z;
    2643                 :            : }
    2644                 :            : GEN
    2645                 :        239 : FqX_split_roots(GEN z, GEN T, GEN p, GEN pol)
    2646                 :            : {
    2647                 :        239 :   GEN S = gel(z,1), L = gel(z,2), rep = FqX_split_equal(L, S, T, p);
    2648         [ +  + ]:        239 :   if (pol) rep = shallowconcat(rep, FqX_div(pol, gel(L,2), T,p));
    2649                 :        239 :   return rep;
    2650                 :            : }
    2651                 :            : GEN
    2652                 :         40 : FqX_split_all(GEN z, GEN T, GEN p)
    2653                 :            : {
    2654                 :         40 :   GEN S = gel(z,1), rep = cgetg(1, t_VEC);
    2655                 :         40 :   long i, l = lg(z);
    2656         [ +  + ]:         85 :   for (i = 2; i < l; i++)
    2657                 :         45 :     rep = shallowconcat(rep, FqX_split_equal(gel(z,i), S, T, p));
    2658                 :         40 :   return rep;
    2659                 :            : }
    2660                 :            : 
    2661                 :            : /* not memory-clean, as FpX_factorff_i, returning only linear factors */
    2662                 :            : static GEN
    2663                 :        419 : FpX_rootsff_i(GEN P, GEN T, GEN p)
    2664                 :            : {
    2665                 :        419 :   GEN V, F = gel(FpX_factor(P,p), 1);
    2666                 :        419 :   long i, lfact = 1, nmax = lgpol(P), n = lg(F), dT = degpol(T);
    2667                 :            : 
    2668                 :        419 :   V = cgetg(nmax,t_COL);
    2669         [ +  + ]:        969 :   for(i=1;i<n;i++)
    2670                 :            :   {
    2671                 :        550 :     GEN R, Fi = gel(F,i);
    2672                 :        550 :     long di = degpol(Fi), j, r;
    2673         [ +  + ]:        550 :     if (dT % di) continue;
    2674                 :        464 :     R = FpX_factorff_irred(gel(F,i),T,p);
    2675                 :        464 :     r = lg(R);
    2676         [ +  + ]:       1511 :     for (j=1; j<r; j++,lfact++)
    2677                 :       1047 :       gel(V,lfact) = Fq_neg(gmael(R,j, 2), T, p);
    2678                 :            :   }
    2679                 :        419 :   setlg(V,lfact);
    2680                 :        419 :   gen_sort_inplace(V, (void*) &cmp_RgX, &cmp_nodata, NULL);
    2681                 :        419 :   return V;
    2682                 :            : }
    2683                 :            : GEN
    2684                 :          0 : FpX_rootsff(GEN P, GEN T, GEN p)
    2685                 :            : {
    2686                 :          0 :   pari_sp av = avma;
    2687                 :          0 :   return gerepilecopy(av, FpX_rootsff_i(P, T, p));
    2688                 :            : }
    2689                 :            : 
    2690                 :            : static GEN
    2691                 :       1234 : F2xqX_quad_roots(GEN P, GEN T)
    2692                 :            : {
    2693                 :       1234 :   GEN b= gel(P,3), c = gel(P,2);
    2694         [ +  + ]:       1234 :   if (degpol(b)>=0)
    2695                 :            :   {
    2696                 :        703 :     GEN z, d = F2xq_div(c, F2xq_sqr(b,T),T);
    2697         [ +  + ]:        703 :     if (F2xq_trace(d,T))
    2698                 :        278 :       return cgetg(1, t_COL);
    2699                 :        425 :     z = F2xq_mul(b, F2xq_Artin_Schreier(d, T), T);
    2700                 :        425 :     return mkcol2(z, F2x_add(b, z));
    2701                 :            :   }
    2702                 :            :   else
    2703                 :       1234 :     return mkcol(F2xq_sqrt(c, T));
    2704                 :            : }
    2705                 :            : 
    2706                 :            : static GEN
    2707                 :       1449 : FqX_quad_roots(GEN x, GEN T, GEN p)
    2708                 :            : {
    2709                 :       1449 :   GEN s, D, nb, b = gel(x,3), c = gel(x,2);
    2710         [ +  + ]:       1449 :   if (equaliu(p, 2))
    2711                 :            :   {
    2712                 :       1234 :     GEN f2 = ZXX_to_F2xX(x, get_FpX_var(T));
    2713                 :       1234 :     s = F2xqX_quad_roots(f2, ZX_to_F2x(get_FpX_mod(T)));
    2714                 :       1234 :     return F2xC_to_ZXC(s);
    2715                 :            :   }
    2716                 :        215 :   D = Fq_sub(Fq_sqr(b,T,p), Fq_Fp_mul(c,utoi(4),T,p), T,p);
    2717                 :        215 :   nb = Fq_neg(b,T,p);
    2718         [ -  + ]:        215 :   if (signe(D)==0)
    2719                 :          0 :     return mkcol(Fq_halve(nb,T, p));
    2720                 :        215 :   s = Fq_sqrt(D,T,p);
    2721         [ +  + ]:        215 :   if (!s) return cgetg(1, t_COL);
    2722                 :        210 :   s = Fq_halve(Fq_add(s,nb,T,p),T, p);
    2723                 :       1449 :   return mkcol2(s,Fq_sub(nb,s,T,p));
    2724                 :            : }
    2725                 :            : 
    2726                 :            : static GEN
    2727                 :       2145 : FqX_roots_i(GEN f, GEN T, GEN p)
    2728                 :            : {
    2729                 :            :   GEN R;
    2730                 :       2145 :   f = FqX_normalize(f, T, p);
    2731         [ -  + ]:       2145 :   if (!signe(f)) pari_err_ROOTS0("FqX_roots");
    2732         [ +  + ]:       2145 :   if (isabsolutepol(f))
    2733                 :            :   {
    2734                 :        424 :     f = simplify_shallow(f);
    2735         [ +  + ]:        424 :     if (typ(f) == t_INT) return cgetg(1, t_COL);
    2736                 :        419 :     return FpX_rootsff_i(f, T, p);
    2737                 :            :   }
    2738         [ +  + ]:       1721 :   if (degpol(f)==2)
    2739                 :       1449 :     return gen_sort(FqX_quad_roots(f,T,p), (void*) &cmp_RgX, &cmp_nodata);
    2740      [ +  +  + ]:        272 :   switch( FqX_split_deg1(&R, f, T, p) )
    2741                 :            :   {
    2742                 :         53 :   case 0: return cgetg(1, t_COL);
    2743         [ +  + ]:        131 :   case 1: if (lg(R) == 1) { R = mkvec(f); break; }
    2744                 :            :             /* fall through */
    2745                 :        194 :   default: R = FqX_split_roots(R, T, p, NULL);
    2746                 :            :   }
    2747                 :        219 :   R = FqXC_roots_from_deg1(R, T, p);
    2748                 :        219 :   gen_sort_inplace(R, (void*) &cmp_RgX, &cmp_nodata, NULL);
    2749                 :       2145 :   return R;
    2750                 :            : }
    2751                 :            : 
    2752                 :            : GEN
    2753                 :       3745 : FqX_roots(GEN x, GEN T, GEN p)
    2754                 :            : {
    2755                 :       3745 :   pari_sp av = avma;
    2756         [ +  + ]:       3745 :   if (!T) return FpX_roots(x, p);
    2757                 :       3745 :   return gerepileupto(av, FqX_roots_i(x, T, p));
    2758                 :            : }
    2759                 :            : 
    2760                 :            : static long
    2761                 :         10 : FqX_sqf_split(GEN *t0, GEN q, GEN T, GEN p)
    2762                 :            : {
    2763                 :         10 :   GEN *t = t0, u = *t, v, S, g, X;
    2764                 :         10 :   long d, dg, N = degpol(u);
    2765                 :            : 
    2766         [ -  + ]:         10 :   if (N == 1) return 1;
    2767                 :         10 :   v = X = pol_x(varn(u));
    2768                 :         10 :   S = FqX_Frobenius_powers(u, T, p);
    2769         [ +  + ]:         20 :   for (d=1; d <= N>>1; d++)
    2770                 :            :   {
    2771                 :         10 :     v = FqX_Frobenius_eval(v, S, u, T, p);
    2772                 :         10 :     g = FqX_normalize(FqX_gcd(FpXX_sub(v,X,p),u, T,p),T,p);
    2773         [ -  + ]:         10 :     dg = degpol(g); if (dg <= 0) continue;
    2774                 :            : 
    2775                 :            :     /* all factors of g have degree d */
    2776                 :         10 :     *t = g;
    2777                 :         10 :     FqX_split(t, d, q, S, T, p);
    2778                 :         10 :     t += dg / d;
    2779                 :         10 :     N -= dg;
    2780         [ -  + ]:         10 :     if (N)
    2781                 :            :     {
    2782                 :          0 :       u = FqX_div(u,g, T,p);
    2783                 :          0 :       v = FqX_rem(v,u, T,p);
    2784                 :            :     }
    2785                 :            :   }
    2786         [ -  + ]:         10 :   if (N) *t++ = u;
    2787                 :         10 :   return t - t0;
    2788                 :            : }
    2789                 :            : 
    2790                 :            : /* not memory-clean */
    2791                 :            : static GEN
    2792                 :        180 : FpX_factorff_i(GEN P, GEN T, GEN p)
    2793                 :            : {
    2794                 :        180 :   GEN V, E, F = FpX_factor(P,p);
    2795                 :        180 :   long i, lfact = 1, nmax = lgpol(P), n = lgcols(F);
    2796                 :            : 
    2797                 :        180 :   V = cgetg(nmax,t_VEC);
    2798                 :        180 :   E = cgetg(nmax,t_VECSMALL);
    2799         [ +  + ]:        495 :   for(i=1;i<n;i++)
    2800                 :            :   {
    2801                 :        315 :     GEN R = FpX_factorff_irred(gmael(F,1,i),T,p), e = gmael(F,2,i);
    2802                 :        315 :     long j, r = lg(R);
    2803         [ +  + ]:       3670 :     for (j=1; j<r; j++,lfact++)
    2804                 :            :     {
    2805                 :       3355 :       gel(V,lfact) = gel(R,j);
    2806                 :       3355 :       gel(E,lfact) = e;
    2807                 :            :     }
    2808                 :            :   }
    2809                 :        180 :   setlg(V,lfact);
    2810                 :        180 :   setlg(E,lfact); return sort_factor_pol(mkvec2(V,E), cmp_RgX);
    2811                 :            : }
    2812                 :            : GEN
    2813                 :          0 : FpX_factorff(GEN P, GEN T, GEN p)
    2814                 :            : {
    2815                 :          0 :   pari_sp av = avma;
    2816                 :          0 :   return gerepilecopy(av, FpX_factorff_i(P, T, p));
    2817                 :            : }
    2818                 :            : 
    2819                 :            : /* assumes varncmp (varn(T), varn(f)) > 0 */
    2820                 :            : static GEN
    2821                 :        290 : FqX_factor_i(GEN f, GEN T, GEN p)
    2822                 :            : {
    2823                 :        290 :   long pg, j, k, e, N, lfact, pk, d = degpol(f);
    2824                 :            :   GEN E, f2, f3, df1, df2, g1, u, q, t;
    2825                 :            : 
    2826      [ +  +  + ]:        290 :   switch(d)
    2827                 :            :   {
    2828                 :          5 :     case -1: retmkmat2(mkcolcopy(f), mkvecsmall(1));
    2829                 :          5 :     case 0: return trivial_fact();
    2830                 :            :   }
    2831                 :        280 :   T = FpX_normalize(T, p);
    2832                 :        280 :   f = FqX_normalize(f, T, p);
    2833         [ +  + ]:        280 :   if (isabsolutepol(f)) return FpX_factorff_i(simplify_shallow(f), T, p);
    2834                 :            : 
    2835                 :        100 :   pg = itos_or_0(p);
    2836                 :        100 :   df2  = NULL; /* gcc -Wall */
    2837                 :        100 :   t = cgetg(d+1,t_VEC);
    2838                 :        100 :   E = cgetg(d+1, t_VECSMALL);
    2839                 :            : 
    2840                 :        100 :   q = powiu(p, degpol(T));
    2841                 :        100 :   e = lfact = 1;
    2842                 :        100 :   pk = 1;
    2843                 :        100 :   f3 = NULL;
    2844                 :        100 :   df1 = FqX_deriv(f, T, p);
    2845                 :            :   for(;;)
    2846                 :            :   {
    2847                 :            :     long nb0;
    2848         [ +  + ]:        190 :     while (!signe(df1))
    2849                 :            :     { /* needs d >= p: pg = 0 can't happen  */
    2850                 :         45 :       pk *= pg; e = pk;
    2851                 :         45 :       f = FqX_frob_deflate(f, T, p);
    2852                 :         45 :       df1 = FqX_deriv(f, T, p); f3 = NULL;
    2853                 :            :     }
    2854         [ +  + ]:        145 :     f2 = f3? f3: FqX_gcd(f,df1, T,p);
    2855         [ +  + ]:        145 :     if (!degpol(f2)) u = f;
    2856                 :            :     else
    2857                 :            :     {
    2858                 :         45 :       g1 = FqX_div(f,f2, T,p);
    2859                 :         45 :       df2 = FqX_deriv(f2, T,p);
    2860         [ +  + ]:         45 :       if (gequal0(df2)) { u = g1; f3 = f2; }
    2861                 :            :       else
    2862                 :            :       {
    2863                 :         35 :         f3 = FqX_gcd(f2,df2, T,p);
    2864         [ +  + ]:         35 :         u = degpol(f3)? FqX_div(f2, f3, T,p): f2;
    2865                 :         35 :         u = FqX_div(g1, u, T,p);
    2866                 :            :       }
    2867                 :            :     }
    2868                 :            :     /* u is square-free (product of irreducibles of multiplicity e) */
    2869                 :        145 :     N = degpol(u);
    2870         [ +  + ]:        145 :     if (N) {
    2871                 :        110 :       nb0 = lfact;
    2872                 :        110 :       gel(t,lfact) = FqX_normalize(u, T,p);
    2873         [ +  + ]:        110 :       if (N == 1) lfact++;
    2874                 :            :       else
    2875                 :            :       {
    2876         [ +  + ]:         75 :         if (!equaliu(p,2))
    2877                 :         50 :           lfact += FqX_split_Berlekamp(&gel(t,lfact), T, p);
    2878                 :            :         else
    2879                 :            :         {
    2880                 :         25 :           GEN P = FqX_split_Trager(gel(t,lfact), T, p);
    2881         [ +  + ]:         25 :           if (P) {
    2882         [ +  + ]:         45 :             for (j = 1; j < lg(P); j++) gel(t,lfact++) = gel(P,j);
    2883                 :            :           } else {
    2884         [ -  + ]:         10 :             if (DEBUGLEVEL) pari_warn(warner, "FqX_split_Trager failed!");
    2885                 :         10 :             lfact += FqX_sqf_split(&gel(t,lfact), q, T, p);
    2886                 :            :           }
    2887                 :            :         }
    2888                 :            :       }
    2889         [ +  + ]:        305 :       for (j = nb0; j < lfact; j++) E[j] = e;
    2890                 :            :     }
    2891                 :            : 
    2892         [ +  + ]:        145 :     if (!degpol(f2)) break;
    2893                 :         45 :     f = f2; df1 = df2; e += pk;
    2894                 :         45 :   }
    2895                 :        100 :   setlg(t, lfact);
    2896                 :        100 :   setlg(E, lfact);
    2897         [ +  + ]:        295 :   for (j=1; j<lfact; j++) gel(t,j) = FqX_normalize(gel(t,j), T,p);
    2898                 :        100 :   (void)sort_factor_pol(mkvec2(t, E), cmp_RgX);
    2899                 :        100 :   k = 1;
    2900         [ +  + ]:        195 :   for (j = 2; j < lfact; j++)
    2901                 :            :   {
    2902         [ +  + ]:         95 :     if (RgX_equal(gel(t,j), gel(t,k)))
    2903                 :         10 :       E[k] += E[j];
    2904                 :            :     else
    2905                 :            :     { /* new factor */
    2906                 :         85 :       k++;
    2907                 :         85 :       E[k] = E[j];
    2908                 :         85 :       gel(t,k) = gel(t,j);
    2909                 :            :     }
    2910                 :            :   }
    2911                 :        100 :   setlg(t, k+1);
    2912                 :        290 :   setlg(E, k+1); return mkvec2(t, E);
    2913                 :            : }
    2914                 :            : 
    2915                 :            : static void
    2916                 :        205 : ffcheck(pari_sp *av, GEN *f, GEN *T, GEN p)
    2917                 :            : {
    2918                 :            :   long v;
    2919         [ -  + ]:        205 :   if (typ(*T)!=t_POL) pari_err_TYPE("factorff",*T);
    2920         [ -  + ]:        205 :   if (typ(*f)!=t_POL) pari_err_TYPE("factorff",*f);
    2921         [ -  + ]:        205 :   if (typ(p)!=t_INT) pari_err_TYPE("factorff",p);
    2922                 :        205 :   v = varn(*T);
    2923         [ -  + ]:        205 :   if (varncmp(v, varn(*f)) <= 0)
    2924                 :          0 :     pari_err_PRIORITY("factorff", *T, "<=", varn(*f));
    2925                 :        205 :   *T = RgX_to_FpX(*T, p); *av = avma;
    2926                 :        205 :   *f = RgX_to_FqX(*f, *T,p);
    2927                 :        205 : }
    2928                 :            : GEN
    2929                 :        205 : factorff(GEN f, GEN p, GEN T)
    2930                 :            : {
    2931                 :            :   pari_sp av;
    2932                 :            :   GEN z;
    2933 [ +  + ][ -  + ]:        205 :   if (!p || !T)
    2934                 :            :   {
    2935                 :            :     long pa, t;
    2936         [ -  + ]:         25 :     if (typ(f) != t_POL) pari_err_TYPE("factorff",f);
    2937                 :         25 :     T = p = NULL;
    2938                 :         25 :     t = RgX_type(f, &p, &T, &pa);
    2939         [ -  + ]:         25 :     if (t != t_FFELT) pari_err_TYPE("factorff",f);
    2940                 :         25 :     return FFX_factor(f,T);
    2941                 :            :   }
    2942                 :        180 :   ffcheck(&av, &f, &T, p); z = FqX_factor_i(f, T, p);
    2943                 :        205 :   return to_Fq_fact(gel(z,1),gel(z,2), T,p, av);
    2944                 :            : }
    2945                 :            : GEN
    2946                 :         50 : polrootsff(GEN f, GEN p, GEN T)
    2947                 :            : {
    2948                 :            :   pari_sp av;
    2949                 :            :   GEN z;
    2950 [ +  + ][ -  + ]:         50 :   if (!p || !T)
    2951                 :            :   {
    2952                 :            :     long pa, t;
    2953         [ -  + ]:         25 :     if (typ(f) != t_POL) pari_err_TYPE("polrootsff",f);
    2954                 :         25 :     T = p = NULL;
    2955                 :         25 :     t = RgX_type(f, &p, &T, &pa);
    2956         [ -  + ]:         25 :     if (t != t_FFELT) pari_err_TYPE("polrootsff",f);
    2957                 :         25 :     return FFX_roots(f,T);
    2958                 :            :   }
    2959                 :         25 :   ffcheck(&av, &f, &T, p); z = FqX_roots_i(f, T, p);
    2960                 :         50 :   return to_FqC(z, T,p, av);
    2961                 :            : }
    2962                 :            : 
    2963                 :            : /* factorization of x modulo (T,p). Assume x already reduced */
    2964                 :            : GEN
    2965                 :        455 : FqX_factor(GEN x, GEN T, GEN p)
    2966                 :            : {
    2967                 :        455 :   pari_sp av = avma;
    2968         [ +  + ]:        455 :   if (!T) return FpX_factor(x, p);
    2969                 :        455 :   return gerepilecopy(av, FqX_factor_i(x, T, p));
    2970                 :            : }
    2971                 :            : /* See also: Isomorphisms between finite field and relative
    2972                 :            :  * factorization in polarit3.c */

Generated by: LCOV version 1.9