Code coverage tests

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

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

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

LCOV - code coverage report
Current view: top level - basemath - FF.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.12.0 lcov report (development 23332-367b47754) Lines: 1109 1170 94.8 %
Date: 2018-12-10 05:41:52 Functions: 134 140 95.7 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2000-2003  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             : /**                   Routines for handling FFELT                       **/
      20             : /**                                                                     **/
      21             : /*************************************************************************/
      22             : 
      23             : /*************************************************************************/
      24             : /**                                                                     **/
      25             : /**                   Low-level constructors                            **/
      26             : /**                                                                     **/
      27             : /*************************************************************************/
      28             : 
      29             : INLINE void
      30    28240909 : _getFF(GEN x, GEN *T, GEN *p, ulong *pp)
      31             : {
      32    28240909 :   *T=gel(x,3);
      33    28240909 :   *p=gel(x,4);
      34    28240909 :   *pp=(*p)[2];
      35    28240909 : }
      36             : 
      37             : INLINE GEN
      38    26767015 : _initFF(GEN x, GEN *T, GEN *p, ulong *pp)
      39             : {
      40    26767015 :   _getFF(x,T,p,pp);
      41    26767015 :   return cgetg(5,t_FFELT);
      42             : }
      43             : 
      44             : INLINE void
      45    19306770 : _checkFF(GEN x, GEN y, const char *s)
      46    19306770 : { if (!FF_samefield(x,y)) pari_err_OP(s,x,y); }
      47             : 
      48             : INLINE GEN
      49    26550266 : _mkFF(GEN x, GEN z, GEN r)
      50             : {
      51    26550266 :   z[1]=x[1];
      52    26550266 :   gel(z,2)=r;
      53    26550266 :   gel(z,3)=gcopy(gel(x,3));
      54    26550264 :   gel(z,4)=icopy(gel(x,4));
      55    26550271 :   return z;
      56             : }
      57             : 
      58             : INLINE GEN
      59     3052370 : _mkFF_i(GEN x, GEN z, GEN r)
      60             : {
      61     3052370 :   z[1]=x[1];
      62     3052370 :   gel(z,2)=r;
      63     3052370 :   gel(z,3)=gel(x,3);
      64     3052370 :   gel(z,4)=gel(x,4);
      65     3052370 :   return z;
      66             : }
      67             : 
      68             : INLINE GEN
      69     2835591 : mkFF_i(GEN x, GEN r)
      70             : {
      71     2835591 :   GEN z = cgetg(5,t_FFELT);
      72     2835594 :   return _mkFF_i(x,z,r);
      73             : }
      74             : 
      75             : /*************************************************************************/
      76             : /**                                                                     **/
      77             : /**                   medium-level constructors                         **/
      78             : /**                                                                     **/
      79             : /*************************************************************************/
      80             : 
      81             : static GEN
      82      431319 : Z_to_raw(GEN x, GEN ff)
      83             : {
      84             :   ulong pp;
      85             :   GEN T, p;
      86      431319 :   _getFF(ff,&T,&p,&pp);
      87      431319 :   switch(ff[1])
      88             :   {
      89             :   case t_FF_FpXQ:
      90       10914 :     return scalarpol(x, varn(T));
      91             :   case t_FF_F2xq:
      92      202174 :     return Z_to_F2x(x,varn(T));
      93             :   default:
      94      218231 :     return Z_to_Flx(x,pp,T[1]);
      95             :   }
      96             : }
      97             : 
      98             : static GEN
      99     1818172 : Rg_to_raw(GEN x, GEN ff)
     100             : {
     101     1818172 :   long tx = typ(x);
     102     1818172 :   switch(tx)
     103             :   {
     104             :   case t_INT: case t_FRAC: case t_PADIC: case t_INTMOD:
     105      431319 :     return Z_to_raw(Rg_to_Fp(x, FF_p_i(ff)), ff);
     106             :   case t_FFELT:
     107     1386853 :     if (!FF_samefield(x,ff))
     108           0 :       pari_err_MODULUS("Rg_to_raw",x,ff);
     109     1386853 :     return gel(x,2);
     110             :   }
     111           0 :   pari_err_TYPE("Rg_to_raw",x);
     112             :   return NULL; /* LCOV_EXCL_LINE */
     113             : }
     114             : 
     115             : static GEN
     116      322394 : FFX_to_raw(GEN x, GEN ff)
     117             : {
     118             :   long i, lx;
     119      322394 :   GEN y = cgetg_copy(x,&lx);
     120      322394 :   y[1] = x[1];
     121     1876190 :   for(i=2; i<lx; i++)
     122     1553796 :     gel(y, i) = Rg_to_raw(gel(x, i), ff);
     123      322394 :   switch (ff[1])
     124             :   {
     125             :   case t_FF_FpXQ:
     126        5619 :     return FpXX_renormalize(y, lx);
     127             :   case t_FF_F2xq:
     128      138956 :     return F2xX_renormalize(y, lx);
     129             :   default:
     130      177819 :     return FlxX_renormalize(y, lx);
     131             :   }
     132             : }
     133             : 
     134             : static GEN
     135       28308 : FFC_to_raw(GEN x, GEN ff)
     136       28308 : { pari_APPLY_same(Rg_to_raw(gel(x, i), ff)) }
     137             : 
     138             : static GEN
     139        8932 : FFM_to_raw(GEN x, GEN ff)
     140        8932 : { pari_APPLY_same(FFC_to_raw(gel(x, i), ff)) }
     141             : 
     142             : /* in place */
     143             : static GEN
     144      654529 : rawFq_to_FF(GEN x, GEN ff)
     145             : {
     146      654529 :   return mkFF_i(ff, typ(x)==t_INT ? scalarpol(x, varn(gel(ff,3))): x);
     147             : }
     148             : 
     149             : /* in place */
     150             : static GEN
     151      109274 : raw_to_FFX(GEN x, GEN ff)
     152             : {
     153      109274 :   long i, lx = lg(x);
     154      109274 :   for (i=2; i<lx; i++) gel(x,i) = rawFq_to_FF(gel(x,i), ff);
     155      109274 :   return x;
     156             : }
     157             : 
     158             : /* in place */
     159             : static GEN
     160      115794 : raw_to_FFC(GEN x, GEN ff)
     161             : {
     162      115794 :   long i, lx = lg(x);
     163      115794 :   for (i=1; i<lx; i++) gel(x,i) = mkFF_i(ff, gel(x,i));
     164      115794 :   return x;
     165             : }
     166             : 
     167             : /* in place */
     168             : static GEN
     169        4298 : raw_to_FFM(GEN x, GEN ff)
     170             : {
     171        4298 :   long i, lx = lg(x);
     172        4298 :   for (i=1; i<lx; i++) gel(x,i) = raw_to_FFC(gel(x,i), ff);
     173        4298 :   return x;
     174             : }
     175             : 
     176             : GEN
     177      131761 : Fq_to_FF(GEN x, GEN ff)
     178             : {
     179             :   ulong pp;
     180      131761 :   GEN r, T, p, z=_initFF(ff,&T,&p,&pp);
     181      131761 :   int is_int = typ(x)==t_INT;
     182      131761 :   switch(ff[1])
     183             :   {
     184             :   case t_FF_FpXQ:
     185        6944 :     r= is_int ? scalarpol(x, varn(T)): x;
     186        6944 :     break;
     187             :   case t_FF_F2xq:
     188       45150 :     r= is_int ? Z_to_F2x(x,T[1]): ZX_to_F2x(x);
     189       45150 :     break;
     190             :   default:
     191       79667 :     r= is_int ? Z_to_Flx(x,pp,T[1]): ZX_to_Flx(x,pp);
     192             :   }
     193      131761 :   setvarn(r, varn(T)); /* paranoia */
     194      131761 :   return _mkFF_i(ff,z,r);
     195             : }
     196             : 
     197             : /*************************************************************************/
     198             : /**                                                                     **/
     199             : /**                   Public functions                                  **/
     200             : /**                                                                     **/
     201             : /*************************************************************************/
     202             : 
     203             : /* Return true if x and y are defined in the same field */
     204             : 
     205             : static int
     206    22113564 : FF_samechar(GEN x, GEN y)
     207             : {
     208    22113564 :   return x[1] == y[1] && equalii(gel(x,4),gel(y,4));
     209             : }
     210             : 
     211             : int
     212    22113564 : FF_samefield(GEN x, GEN y)
     213             : {
     214    22113564 :   return FF_samechar(x, y) && gidentical(gel(x,3),gel(y,3));
     215             : }
     216             : 
     217             : int
     218       55412 : FF_equal(GEN x, GEN y)
     219       55412 : { return FF_samefield(x,y) && gidentical(gel(x,2),gel(y,2)); }
     220             : 
     221             : int
     222     8918915 : FF_equal0(GEN x)
     223             : {
     224     8918915 :   return lgpol(gel(x,2))==0;
     225             : }
     226             : 
     227             : int
     228       23723 : FF_equal1(GEN x)
     229             : {
     230       23723 :   GEN A = gel(x,2);
     231       23723 :   switch(x[1])
     232             :   {
     233             :   case t_FF_FpXQ:
     234        5391 :     return degpol(A)==0 && gequal1(gel(A,2));
     235             :   default:
     236       18332 :     return degpol(A)==0 && A[2]==1;
     237             :   }
     238             : }
     239             : 
     240             : static int
     241           0 : Fp_cmp_1(GEN x, GEN p)
     242           0 : { pari_sp av = avma; return gc_bool(av, equalii(x, addis(p,-1))); }
     243             : 
     244             : int
     245          42 : FF_equalm1(GEN x)
     246             : {
     247             :   ulong pp;
     248          42 :   GEN T, p, y = gel(x,2);
     249          42 :   _getFF(x,&T,&p,&pp);
     250          42 :   switch(x[1])
     251             :   {
     252             :   case t_FF_FpXQ:
     253           0 :     return (degpol(y) == 0 && Fp_cmp_1(gel(y,2), p));
     254             :   default:
     255          42 :     return (degpol(y) == 0 && uel(y,2) == pp-1);
     256             :   }
     257             : }
     258             : 
     259             : GEN
     260        5723 : FF_zero(GEN x)
     261             : {
     262             :   ulong pp;
     263        5723 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     264        5723 :   switch(x[1])
     265             :   {
     266             :   case t_FF_FpXQ:
     267         497 :     r=zeropol(varn(T));
     268         497 :     break;
     269             :   case t_FF_F2xq:
     270        1327 :     r=zero_F2x(T[1]);
     271        1327 :     break;
     272             :   default:
     273        3899 :     r=zero_Flx(T[1]);
     274             :   }
     275        5723 :   return _mkFF(x,z,r);
     276             : }
     277             : 
     278             : GEN
     279       13160 : FF_1(GEN x)
     280             : {
     281             :   ulong pp;
     282       13160 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     283       13160 :   switch(x[1])
     284             :   {
     285             :   case t_FF_FpXQ:
     286         137 :     r=pol_1(varn(T));
     287         137 :     break;
     288             :   case t_FF_F2xq:
     289        8064 :     r=pol1_F2x(T[1]);
     290        8064 :     break;
     291             :   default:
     292        4959 :     r=pol1_Flx(T[1]);
     293             :   }
     294       13160 :   return _mkFF(x,z,r);
     295             : }
     296             : 
     297             : GEN
     298         469 : FF_gen(GEN x)
     299             : {
     300             :   ulong pp;
     301         469 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     302         469 :   switch(x[1])
     303             :   {
     304             :   case t_FF_FpXQ:
     305         154 :     r = pol_x(varn(T));
     306         154 :     if (degpol(T)==1) r = FpX_rem(r, T, p);
     307         154 :     break;
     308             :   case t_FF_F2xq:
     309         154 :     r = polx_F2x(T[1]);
     310         154 :     if (F2x_degree(T)==1) r = F2x_rem(r, T);
     311         154 :     break;
     312             :   default:
     313         161 :     r = polx_Flx(T[1]);
     314         161 :     if (degpol(T)==1) r = Flx_rem(r, T, pp);
     315             :   }
     316         469 :   return _mkFF(x,z,r);
     317             : }
     318             : GEN
     319       54481 : FF_q(GEN x)
     320             : {
     321             :   ulong pp;
     322             :   GEN T, p;
     323       54481 :   _getFF(x,&T,&p,&pp);
     324       54481 :   switch(x[1])
     325             :   {
     326             :   case t_FF_FpXQ:
     327        1681 :     return powiu(p, degpol(T));
     328             :     break;
     329             :   case t_FF_F2xq:
     330        7028 :     return int2n(F2x_degree(T));
     331             :     break;
     332             :   default:
     333       45772 :     return powuu(pp,degpol(T));
     334             :   }
     335             : }
     336             : 
     337             : GEN
     338          56 : FF_p(GEN x)
     339             : {
     340          56 :   return icopy(gel(x,4));
     341             : }
     342             : 
     343             : GEN
     344     2178519 : FF_p_i(GEN x)
     345             : {
     346     2178519 :   return gel(x,4);
     347             : }
     348             : 
     349             : GEN
     350      165543 : FF_mod(GEN x)
     351             : {
     352      165543 :   switch(x[1])
     353             :   {
     354             :   case t_FF_FpXQ:
     355         365 :     return ZX_copy(gel(x,3));
     356             :   case t_FF_F2xq:
     357         371 :     return F2x_to_ZX(gel(x,3));
     358             :   default:
     359      164807 :     return Flx_to_ZX(gel(x,3));
     360             :   }
     361             : }
     362             : 
     363             : long
     364         364 : FF_f(GEN x)
     365             : {
     366         364 :   switch(x[1])
     367             :   {
     368             :   case t_FF_F2xq:
     369         147 :     return F2x_degree(gel(x,3));
     370             :   default:
     371         217 :     return degpol(gel(x,3));
     372             :   }
     373             : }
     374             : 
     375             : GEN
     376      918691 : FF_to_F2xq(GEN x)
     377             : {
     378      918691 :   switch(x[1])
     379             :   {
     380             :   case t_FF_FpXQ:
     381           0 :     return ZX_to_F2x(gel(x,2));
     382             :   case t_FF_F2xq:
     383      918691 :     return zv_copy(gel(x,2));
     384             :   default:
     385           0 :     return Flx_to_F2x(gel(x,2));
     386             :   }
     387             : }
     388             : 
     389             : GEN
     390           0 : FF_to_F2xq_i(GEN x)
     391             : {
     392           0 :   switch(x[1])
     393             :   {
     394             :   case t_FF_FpXQ:
     395           0 :     return ZX_to_F2x(gel(x,2));
     396             :   case t_FF_F2xq:
     397           0 :     return gel(x,2);
     398             :   default:
     399           0 :     return Flx_to_F2x(gel(x,2));
     400             :   }
     401             : }
     402             : 
     403             : GEN
     404      726168 : FF_to_Flxq(GEN x)
     405             : {
     406      726168 :   switch(x[1])
     407             :   {
     408             :   case t_FF_FpXQ:
     409           0 :     return ZX_to_Flx(gel(x,2),itou(gel(x,4)));
     410             :   case t_FF_F2xq:
     411           0 :     return F2x_to_Flx(gel(x,2));
     412             :   default:
     413      726168 :     return zv_copy(gel(x,2));
     414             :   }
     415             : }
     416             : 
     417             : GEN
     418           0 : FF_to_Flxq_i(GEN x)
     419             : {
     420           0 :   switch(x[1])
     421             :   {
     422             :   case t_FF_FpXQ:
     423           0 :     return ZX_to_Flx(gel(x,2),itou(gel(x,4)));
     424             :   case t_FF_F2xq:
     425           0 :     return F2x_to_Flx(gel(x,2));
     426             :   default:
     427           0 :     return gel(x,2);
     428             :   }
     429             : }
     430             : 
     431             : GEN
     432       17946 : FF_to_FpXQ(GEN x)
     433             : {
     434       17946 :   switch(x[1])
     435             :   {
     436             :   case t_FF_FpXQ:
     437       16470 :     return ZX_copy(gel(x,2));
     438             :   case t_FF_F2xq:
     439          63 :     return F2x_to_ZX(gel(x,2));
     440             :   default:
     441        1413 :     return Flx_to_ZX(gel(x,2));
     442             :   }
     443             : }
     444             : 
     445             : GEN
     446      170821 : FF_to_FpXQ_i(GEN x)
     447             : {
     448      170821 :   switch(x[1])
     449             :   {
     450             :   case t_FF_FpXQ:
     451        1214 :     return gel(x,2);
     452             :   case t_FF_F2xq:
     453        1204 :     return F2x_to_ZX(gel(x,2));
     454             :   default:
     455      168403 :     return Flx_to_ZX(gel(x,2));
     456             :   }
     457             : }
     458             : 
     459             : GEN
     460     9202144 : FF_add(GEN x, GEN y)
     461             : {
     462             :   ulong pp;
     463     9202144 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     464     9202144 :   _checkFF(x,y,"+");
     465     9202144 :   switch(x[1])
     466             :   {
     467             :   case t_FF_FpXQ:
     468      310877 :     r=FpX_add(gel(x,2),gel(y,2),p);
     469      310877 :     break;
     470             :   case t_FF_F2xq:
     471     1130457 :     r=F2x_add(gel(x,2),gel(y,2));
     472     1130456 :     break;
     473             :   default:
     474     7760810 :     r=Flx_add(gel(x,2),gel(y,2),pp);
     475             :   }
     476     9202143 :   return _mkFF(x,z,r);
     477             : }
     478             : GEN
     479      283087 : FF_sub(GEN x, GEN y)
     480             : {
     481             :   ulong pp;
     482      283087 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     483      283087 :   _checkFF(x,y,"+");
     484      283087 :   switch(x[1])
     485             :   {
     486             :   case t_FF_FpXQ:
     487        2892 :     r=FpX_sub(gel(x,2),gel(y,2),p);
     488        2892 :     break;
     489             :   case t_FF_F2xq:
     490      155101 :     r=F2x_add(gel(x,2),gel(y,2));
     491      155101 :     break;
     492             :   default:
     493      125094 :     r=Flx_sub(gel(x,2),gel(y,2),pp);
     494             :   }
     495      283087 :   return _mkFF(x,z,r);
     496             : }
     497             : 
     498             : GEN
     499     1658353 : FF_Z_add(GEN x, GEN y)
     500             : {
     501             :   ulong pp;
     502     1658353 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     503     1658353 :   switch(x[1])
     504             :   {
     505             :   case t_FF_FpXQ:
     506             :     {
     507        9234 :       pari_sp av=avma;
     508        9234 :       r=gerepileupto(av,FpX_Fp_add(gel(x,2),modii(y,p),p));
     509        9234 :       break;
     510             :     }
     511             :   case t_FF_F2xq:
     512      792285 :     r=mpodd(y)?F2x_1_add(gel(x,2)):vecsmall_copy(gel(x,2));
     513      792285 :     break;
     514             :   default:
     515      856834 :     r=Flx_Fl_add(gel(x,2),umodiu(y,pp),pp);
     516             :   }
     517     1658353 :   return _mkFF(x,z,r);
     518             : }
     519             : 
     520             : GEN
     521        1337 : FF_Q_add(GEN x, GEN y)
     522             : {
     523             :   ulong pp;
     524        1337 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     525        1337 :   switch(x[1])
     526             :   {
     527             :   case t_FF_FpXQ:
     528             :     {
     529           1 :       pari_sp av=avma;
     530           1 :       r=gerepileupto(av,FpX_Fp_add(gel(x,2),Rg_to_Fp(y,p),p));
     531           1 :       break;
     532             :     }
     533             :   case t_FF_F2xq:
     534           7 :     r=Rg_to_Fl(y,pp)?F2x_1_add(gel(x,2)):vecsmall_copy(gel(x,2));
     535           7 :     break;
     536             :   default:
     537        1329 :     r=Flx_Fl_add(gel(x,2),Rg_to_Fl(y,pp),pp);
     538             :   }
     539        1337 :   return _mkFF(x,z,r);
     540             : }
     541             : 
     542             : GEN
     543       40055 : FF_neg(GEN x)
     544             : {
     545             :   ulong pp;
     546       40055 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     547       40055 :   switch(x[1])
     548             :   {
     549             :   case t_FF_FpXQ:
     550        3979 :     r=FpX_neg(gel(x,2),p);
     551        3979 :     break;
     552             :   case t_FF_F2xq:
     553       14330 :     r=vecsmall_copy(gel(x,2));
     554       14330 :     break;
     555             :   default:
     556       21746 :     r=Flx_neg(gel(x,2),pp);
     557             :   }
     558       40055 :   return _mkFF(x,z,r);
     559             : }
     560             : 
     561             : GEN
     562       85015 : FF_neg_i(GEN x)
     563             : {
     564             :   ulong pp;
     565       85015 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     566       85015 :   switch(x[1])
     567             :   {
     568             :   case t_FF_FpXQ:
     569        1748 :     r=FpX_neg(gel(x,2),p);
     570        1748 :     break;
     571             :   case t_FF_F2xq:
     572        8477 :     r=gel(x,2);
     573        8477 :     break;
     574             :   default:
     575       74790 :     r=Flx_neg(gel(x,2),pp);
     576             :   }
     577       85015 :   return _mkFF_i(x,z,r);
     578             : }
     579             : 
     580             : GEN
     581        1687 : FF_map(GEN m, GEN x)
     582             : {
     583             :   ulong pp;
     584        1687 :   GEN r, T, p, z=_initFF(m,&T,&p,&pp);
     585        1687 :   switch(m[1])
     586             :   {
     587             :   case t_FF_FpXQ:
     588         560 :     r=FpX_FpXQ_eval(gel(x,2),gel(m,2),T,p);
     589         560 :     break;
     590             :   case t_FF_F2xq:
     591         672 :     r=F2x_F2xq_eval(gel(x,2),gel(m,2),T);
     592         672 :     break;
     593             :   default:
     594         455 :     r=Flx_Flxq_eval(gel(x,2),gel(m,2),T,pp);
     595             :   }
     596        1687 :   return _mkFF(m,z,r);
     597             : }
     598             : 
     599             : GEN
     600          42 : FF_Frobenius(GEN x, long e)
     601             : {
     602             :   ulong pp;
     603          42 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     604          42 :   ulong n = umodsu(e, FF_f(x));
     605          42 :   if (n==0) return gcopy(x);
     606          42 :   switch(x[1])
     607             :   {
     608             :   case t_FF_FpXQ:
     609          14 :     r=FpXQ_pow(gel(x,2),p,T,p);
     610          14 :     if (n>1) r=FpXQ_autpow(r,n,T,p);
     611          14 :     break;
     612             :   case t_FF_F2xq:
     613          14 :     r=F2xq_sqr(gel(x,2),T);
     614          14 :     if (n>1) r=F2xq_autpow(r,n,T);
     615          14 :     break;
     616             :   default:
     617          14 :     r=Flxq_powu(gel(x,2),pp,T,pp);
     618          14 :     if (n>1) r=Flxq_autpow(r,n,T,pp);
     619             :   }
     620          42 :   return _mkFF(x,z,r);
     621             : }
     622             : 
     623             : GEN
     624         966 : FFX_preimage(GEN x, GEN F, GEN y)
     625             : {
     626             :   GEN r, T, p, z;
     627             :   ulong pp;
     628         966 :   if (FF_equal0(x)) return FF_zero(y);
     629         875 :   z=_initFF(y,&T,&p,&pp);
     630         875 :   F = FFX_to_raw(F, y);
     631         875 :   switch(y[1])
     632             :   {
     633             :   case t_FF_FpXQ:
     634         322 :     r = FpXQX_rem(gel(x,2), F, T, p);
     635         322 :     break;
     636             :   case t_FF_F2xq:
     637         322 :     r = F2xqX_rem(F2x_to_F2xX(gel(x,2),T[1]), F, T);
     638         322 :     break;
     639             :   default:
     640         231 :     r = FlxqX_rem(Flx_to_FlxX(gel(x,2),T[1]), F, T, pp);
     641             :   }
     642         875 :   if (degpol(r) > 0) return NULL;
     643         791 :   r = (y[1] == t_FF_FpXQ)? Fq_to_FpXQ(gel(r,2),T, p): gel(r,2);
     644         791 :   return _mkFF(y,z,r);
     645             : }
     646             : 
     647             : GEN
     648     9584568 : FF_mul(GEN x, GEN y)
     649             : {
     650             :   ulong pp;
     651     9584568 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     652     9584568 :   pari_sp av=avma;
     653     9584568 :   _checkFF(x,y,"*");
     654     9584568 :   switch(x[1])
     655             :   {
     656             :   case t_FF_FpXQ:
     657      299233 :     r=FpXQ_mul(gel(x,2),gel(y,2),T,p);
     658      299233 :     break;
     659             :   case t_FF_F2xq:
     660      778625 :     r=F2xq_mul(gel(x,2),gel(y,2),T);
     661      778624 :     break;
     662             :   default:
     663     8506710 :     r=Flxq_mul(gel(x,2),gel(y,2),T,pp);
     664             :   }
     665     9584567 :   return _mkFF(x,z,gerepileupto(av, r));
     666             : }
     667             : 
     668             : GEN
     669     2242156 : FF_Z_mul(GEN x, GEN y)
     670             : {
     671             :   ulong pp;
     672     2242156 :   GEN r, T, p, A = gel(x,2), z=_initFF(x,&T,&p,&pp);
     673     2242156 :   switch(x[1])
     674             :   {
     675             :   case t_FF_FpXQ: /* modii(y,p) left on stack for efficiency */
     676       42069 :     r = FpX_Fp_mul(A, modii(y,p),p);
     677       42069 :     break;
     678             :   case t_FF_F2xq:
     679     1121155 :     r = mpodd(y)? vecsmall_copy(A): zero_Flx(A[1]);
     680     1121155 :     break;
     681             :   default:
     682     1078932 :     r = Flx_Fl_mul(A, umodiu(y,pp), pp);
     683             :   }
     684     2242156 :   return _mkFF(x,z,r);
     685             : }
     686             : 
     687             : GEN
     688        4529 : FF_Z_Z_muldiv(GEN x, GEN a, GEN b)
     689             : {
     690             :   ulong pp;
     691        4529 :   GEN r, T, p, A = gel(x,2), z=_initFF(x,&T,&p,&pp);
     692        4529 :   switch(x[1])
     693             :   {
     694             :   case t_FF_FpXQ: /* Fp_div(a,b,p) left on stack for efficiency */
     695          93 :     r = FpX_Fp_mul(A, Fp_div(a,b,p), p);
     696          93 :     break;
     697             :   case t_FF_F2xq:
     698         651 :     if (!mpodd(b)) pari_err_INV("FF_Z_Z_muldiv", b);
     699         644 :     r = mpodd(a)? vecsmall_copy(A): zero_Flx(A[1]);
     700         644 :     break;
     701             :   default:
     702        3785 :     r = Flx_Fl_mul(A, Fl_div(umodiu(a,pp),umodiu(b,pp),pp),pp);
     703             :   }
     704        4515 :   return _mkFF(x,z,r);
     705             : }
     706             : 
     707             : GEN
     708     2893703 : FF_sqr(GEN x)
     709             : {
     710             :   ulong pp;
     711     2893703 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     712     2893703 :   switch(x[1])
     713             :   {
     714             :   case t_FF_FpXQ:
     715             :     {
     716       16296 :       pari_sp av=avma;
     717       16296 :       r=gerepileupto(av,FpXQ_sqr(gel(x,2),T,p));
     718       16296 :       break;
     719             :     }
     720             :   case t_FF_F2xq:
     721      468626 :     r=F2xq_sqr(gel(x,2),T);
     722      468625 :     break;
     723             :   default:
     724     2408781 :     r=Flxq_sqr(gel(x,2),T,pp);
     725             :   }
     726     2893702 :   return _mkFF(x,z,r);
     727             : }
     728             : 
     729             : GEN
     730      216861 : FF_mul2n(GEN x, long n)
     731             : {
     732             :   ulong pp;
     733      216861 :   GEN r, T, p, A = gel(x,2), z=_initFF(x,&T,&p,&pp);
     734      216861 :   switch(x[1])
     735             :   {
     736             :   case t_FF_FpXQ:
     737             :     {
     738             :       GEN p1; /* left on stack for efficiency */
     739        2965 :       if (n>0) p1=remii(int2n(n),p);
     740          43 :       else p1=Fp_inv(remii(int2n(-n),p),p);
     741        2965 :       r = FpX_Fp_mul(A, p1, p);
     742             :     }
     743        2965 :     break;
     744             :   case t_FF_F2xq:
     745       91036 :     if (n<0) pari_err_INV("FF_mul2n", gen_2);
     746       91036 :     r = n==0? vecsmall_copy(A): zero_Flx(A[1]);
     747       91036 :     break;
     748             :   default:
     749             :     {
     750             :       ulong l1;
     751      122860 :       if (n>0) l1 = umodiu(int2n(n),pp);
     752         384 :       else l1 = Fl_inv(umodiu(int2n(-n),pp),pp);
     753      122860 :       r = Flx_Fl_mul(A,l1,pp);
     754             :     }
     755             :   }
     756      216861 :   return _mkFF(x,z,r);
     757             : }
     758             : 
     759             : GEN
     760       14078 : FF_inv(GEN x)
     761             : {
     762             :   ulong pp;
     763       14078 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     764       14078 :   pari_sp av=avma;
     765       14078 :   switch(x[1])
     766             :   {
     767             :   case t_FF_FpXQ:
     768        3298 :     r=gerepileupto(av,FpXQ_inv(gel(x,2),T,p));
     769        3298 :     break;
     770             :   case t_FF_F2xq:
     771        8471 :     r=F2xq_inv(gel(x,2),T);
     772        8471 :     break;
     773             :   default:
     774        2309 :     r=Flxq_inv(gel(x,2),T,pp);
     775             :   }
     776       14078 :   return _mkFF(x,z,r);
     777             : }
     778             : 
     779             : GEN
     780      236790 : FF_div(GEN x, GEN y)
     781             : {
     782             :   ulong pp;
     783      236790 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     784      236790 :   pari_sp av=avma;
     785      236790 :   _checkFF(x,y,"/");
     786      236790 :   switch(x[1])
     787             :   {
     788             :   case t_FF_FpXQ:
     789        3387 :     r=gerepileupto(av,FpXQ_div(gel(x,2),gel(y,2),T,p));
     790        3387 :     break;
     791             :   case t_FF_F2xq:
     792      101753 :     r=gerepileupto(av,F2xq_div(gel(x,2),gel(y,2),T));
     793      101683 :     break;
     794             :   default:
     795      131650 :     r=gerepileupto(av,Flxq_div(gel(x,2),gel(y,2),T,pp));
     796             :   }
     797      236692 :   return _mkFF(x,z,r);
     798             : }
     799             : 
     800             : GEN
     801         301 : Z_FF_div(GEN n, GEN x)
     802             : {
     803             :   ulong pp;
     804         301 :   GEN r, T, p, A = gel(x,2), z=_initFF(x,&T,&p,&pp);
     805         301 :   pari_sp av=avma;
     806         301 :   switch(x[1])
     807             :   {
     808             :   case t_FF_FpXQ:
     809           8 :     r = gerepileupto(av,FpX_Fp_mul(FpXQ_inv(A,T,p),modii(n,p),p));
     810           8 :     break;
     811             :   case t_FF_F2xq:
     812          56 :     r = F2xq_inv(A,T); /*Check for division by 0*/
     813          56 :     if(!mpodd(n)) { set_avma(av); r = zero_Flx(A[1]); }
     814          56 :     break;
     815             :   default:
     816         237 :     r = gerepileupto(av, Flx_Fl_mul(Flxq_inv(A,T,pp),umodiu(n,pp),pp));
     817             :   }
     818         301 :   return _mkFF(x,z,r);
     819             : }
     820             : 
     821             : GEN
     822         217 : FF_sqrtn(GEN x, GEN n, GEN *zetan)
     823             : {
     824             :   ulong pp;
     825         217 :   GEN r, T, p, y=_initFF(x,&T,&p,&pp);
     826         217 :   switch (x[1])
     827             :   {
     828             :   case t_FF_FpXQ:
     829          66 :     r=FpXQ_sqrtn(gel(x,2),n,T,p,zetan);
     830          59 :     break;
     831             :   case t_FF_F2xq:
     832          28 :     r=F2xq_sqrtn(gel(x,2),n,T,zetan);
     833          21 :     break;
     834             :   default:
     835         123 :     r=Flxq_sqrtn(gel(x,2),n,T,pp,zetan);
     836             :   }
     837         196 :   if (!r) pari_err_SQRTN("FF_sqrtn",x);
     838         196 :   (void)_mkFF(x, y, r);
     839         196 :   if (zetan)
     840             :   {
     841         140 :     GEN z = cgetg(lg(y),t_FFELT);
     842         140 :     *zetan=_mkFF(x, z, *zetan);
     843             :   }
     844         196 :   return y;
     845             : }
     846             : 
     847             : GEN
     848          70 : FF_sqrt(GEN x)
     849             : {
     850             :   ulong pp;
     851          70 :   GEN r, T, p, y=_initFF(x,&T,&p,&pp);
     852          70 :   switch (x[1])
     853             :   {
     854             :   case t_FF_FpXQ:
     855           1 :     r = FpXQ_sqrt(gel(x,2),T,p);
     856           1 :     break;
     857             :   case t_FF_F2xq:
     858          49 :     r = F2xq_sqrt(gel(x,2),T);
     859          49 :     break;
     860             :   default:
     861          20 :     r = Flxq_sqrt(gel(x,2),T,pp);
     862             :   }
     863          70 :   if (!r) pari_err_SQRTN("FF_sqrt",x);
     864          70 :   return _mkFF(x, y, r);
     865             : }
     866             : 
     867             : long
     868           7 : FF_issquare(GEN x)
     869             : {
     870             :   GEN T, p;
     871             :   ulong pp;
     872           7 :   _getFF(x, &T, &p, &pp);
     873           7 :   switch(x[1])
     874             :   {
     875             :   case t_FF_FpXQ:
     876           0 :     return FpXQ_issquare(gel(x,2), T, p);
     877             :   case t_FF_F2xq:
     878           7 :     return 1;
     879             :   default: /* case t_FF_Flxq: */
     880           0 :     return Flxq_issquare(gel(x,2), T, pp);
     881             :   }
     882             : }
     883             : 
     884             : long
     885         182 : FF_issquareall(GEN x, GEN *pt)
     886             : {
     887         182 :   if (!pt) return FF_issquare(x);
     888         175 :   return FF_ispower(x, gen_2, pt);
     889             : }
     890             : 
     891             : long
     892         203 : FF_ispower(GEN x, GEN K, GEN *pt)
     893             : {
     894             :   ulong pp;
     895             :   GEN r, T, p;
     896         203 :   pari_sp av = avma;
     897             : 
     898         203 :   if (FF_equal0(x)) { if (pt) *pt = gcopy(x); return 1; }
     899         203 :   _getFF(x, &T, &p, &pp);
     900         203 :   if (pt) *pt = cgetg(5,t_FFELT);
     901         203 :   switch(x[1])
     902             :   {
     903             :   case t_FF_FpXQ:
     904          71 :     r = FpXQ_sqrtn(gel(x,2),K,T,p,NULL);
     905          71 :     break;
     906             :   case t_FF_F2xq:
     907          42 :     r = F2xq_sqrtn(gel(x,2),K,T,NULL);
     908          42 :     break;
     909             :   default: /* case t_FF_Flxq: */
     910          90 :     r = Flxq_sqrtn(gel(x,2),K,T,pp,NULL);
     911          90 :     break;
     912             :   }
     913         203 :   if (!r) return gc_long(av,0);
     914         133 :   if (pt) { (void)_mkFF(x,*pt,r); }
     915         133 :   return 1;
     916             : }
     917             : 
     918             : GEN
     919          80 : FF_pow(GEN x, GEN n)
     920             : {
     921             :   ulong pp;
     922          80 :   GEN r, T, p, z=_initFF(x,&T,&p,&pp);
     923          80 :   switch(x[1])
     924             :    {
     925             :    case t_FF_FpXQ:
     926          17 :      r = FpXQ_pow(gel(x,2), n, T, p);
     927          17 :      break;
     928             :    case t_FF_F2xq:
     929           3 :      r = F2xq_pow(gel(x,2), n, T);
     930           3 :      break;
     931             :    default:
     932          60 :      r = Flxq_pow(gel(x,2), n, T, pp);
     933             :    }
     934          80 :   return _mkFF(x,z,r);
     935             : }
     936             : 
     937             : GEN
     938          28 : FF_norm(GEN x)
     939             : {
     940             :   ulong pp;
     941             :   GEN T,p;
     942          28 :   _getFF(x,&T,&p,&pp);
     943          28 :   switch (x[1])
     944             :   {
     945             :   case t_FF_FpXQ:
     946           1 :     return FpXQ_norm(gel(x,2),T,p);
     947             :   case t_FF_F2xq:
     948           7 :     return lgpol(gel(x,2))?gen_1:gen_0;
     949             :   default:
     950          20 :     return utoi(Flxq_norm(gel(x,2),T,pp));
     951             :   }
     952             : }
     953             : 
     954             : GEN
     955          28 : FF_trace(GEN x)
     956             : {
     957             :   ulong pp;
     958             :   GEN T,p;
     959          28 :   _getFF(x,&T,&p,&pp);
     960          28 :   switch(x[1])
     961             :   {
     962             :   case t_FF_FpXQ:
     963           1 :     return FpXQ_trace(gel(x,2),T,p);
     964             :   case t_FF_F2xq:
     965           7 :     return F2xq_trace(gel(x,2),T)?gen_1:gen_0;
     966             :   default:
     967          20 :     return utoi(Flxq_trace(gel(x,2),T,pp));
     968             :   }
     969             : }
     970             : 
     971             : GEN
     972          28 : FF_conjvec(GEN x)
     973             : {
     974             :   ulong pp;
     975             :   GEN r,T,p,v;
     976             :   long i,l;
     977             :   pari_sp av;
     978          28 :   _getFF(x,&T,&p,&pp);
     979          28 :   av = avma;
     980          28 :   switch(x[1])
     981             :   {
     982             :   case t_FF_FpXQ:
     983           1 :     v = FpXQ_conjvec(gel(x,2), T, p);
     984           1 :     break;
     985             :   case t_FF_F2xq:
     986           7 :     v = F2xq_conjvec(gel(x,2), T);
     987           7 :     break;
     988             :   default:
     989          20 :     v = Flxq_conjvec(gel(x,2), T, pp);
     990             :   }
     991          28 :   l = lg(v); r = cgetg(l, t_COL);
     992         252 :   for(i=1; i<l; i++)
     993         224 :     gel(r,i) = mkFF_i(x, gel(v,i));
     994          28 :   return gerepilecopy(av, r);
     995             : }
     996             : 
     997             : GEN
     998          28 : FF_charpoly(GEN x)
     999             : {
    1000             :   ulong pp;
    1001             :   GEN T,p;
    1002          28 :   pari_sp av=avma;
    1003          28 :   _getFF(x,&T,&p,&pp);
    1004          28 :   switch(x[1])
    1005             :   {
    1006             :   case t_FF_FpXQ:
    1007           1 :     return gerepileupto(av,FpXQ_charpoly(gel(x,2), T, p));
    1008             :   case t_FF_F2xq:
    1009           7 :     return gerepileupto(av,Flx_to_ZX(Flxq_charpoly(F2x_to_Flx(gel(x,2)),
    1010             :                                                    F2x_to_Flx(T), 2UL)));
    1011             :   default:
    1012          20 :     return gerepileupto(av,Flx_to_ZX(Flxq_charpoly(gel(x,2), T, pp)));
    1013             :   }
    1014             : }
    1015             : 
    1016             : GEN
    1017         154 : FF_minpoly(GEN x)
    1018             : {
    1019             :   ulong pp;
    1020             :   GEN T,p;
    1021         154 :   pari_sp av=avma;
    1022         154 :   _getFF(x,&T,&p,&pp);
    1023         154 :   switch(x[1])
    1024             :   {
    1025             :   case t_FF_FpXQ:
    1026          43 :     return gerepileupto(av,FpXQ_minpoly(gel(x,2), T, p));
    1027             :   case t_FF_F2xq:
    1028          49 :     return gerepileupto(av,Flx_to_ZX(Flxq_minpoly(F2x_to_Flx(gel(x,2)),
    1029             :                                                   F2x_to_Flx(T), 2UL)));
    1030             :   default:
    1031          62 :     return gerepileupto(av,Flx_to_ZX(Flxq_minpoly(gel(x,2), T, pp)));
    1032             :   }
    1033             : }
    1034             : 
    1035             : GEN
    1036         182 : FF_log(GEN x, GEN g, GEN ord)
    1037             : {
    1038         182 :   pari_sp av=avma;
    1039             :   ulong pp;
    1040             :   GEN r, T, p;
    1041         182 :   _getFF(x,&T,&p,&pp);
    1042         182 :   _checkFF(x,g,"log");
    1043         182 :   switch(x[1])
    1044             :   {
    1045             :   case t_FF_FpXQ:
    1046          44 :     if (!ord) ord = factor_pn_1(p,degpol(T));
    1047          44 :     r = FpXQ_log(gel(x,2), gel(g,2), ord, T, p);
    1048          16 :     break;
    1049             :   case t_FF_F2xq:
    1050          28 :     if (!ord) ord = factor_pn_1(gen_2,F2x_degree(T));
    1051          28 :     r = F2xq_log(gel(x,2), gel(g,2), ord, T);
    1052          28 :     break;
    1053             :   default:
    1054         110 :     if (!ord) ord = factor_pn_1(p,degpol(T));
    1055         110 :     r = Flxq_log(gel(x,2), gel(g,2), ord, T, pp);
    1056             :   }
    1057         154 :   return gerepileupto(av, r);
    1058             : }
    1059             : 
    1060             : GEN
    1061          28 : FF_order(GEN x, GEN o)
    1062             : {
    1063          28 :   pari_sp av=avma;
    1064             :   ulong pp;
    1065             :   GEN r, T,p;
    1066          28 :   _getFF(x,&T,&p,&pp);
    1067          28 :   switch(x[1])
    1068             :   {
    1069             :   case t_FF_FpXQ:
    1070           1 :     if (!o) o = factor_pn_1(p,degpol(T));
    1071           1 :     r = FpXQ_order(gel(x,2), o, T, p);
    1072           1 :     break;
    1073             :   case t_FF_F2xq:
    1074           7 :     if (!o) o = factor_pn_1(gen_2,F2x_degree(T));
    1075           7 :     r = F2xq_order(gel(x,2), o, T);
    1076           7 :     break;
    1077             :   default:
    1078          20 :     if (!o) o = factor_pn_1(p,degpol(T));
    1079          20 :     r = Flxq_order(gel(x,2), o, T, pp);
    1080             :   }
    1081          28 :   if (!o) r = gerepileuptoint(av,r);
    1082          28 :   return r;
    1083             : }
    1084             : 
    1085             : GEN
    1086         392 : FF_primroot(GEN x, GEN *o)
    1087             : {
    1088             :   ulong pp;
    1089         392 :   GEN r,T,p,z=_initFF(x,&T,&p,&pp);
    1090         392 :   switch(x[1])
    1091             :   {
    1092             :   case t_FF_FpXQ:
    1093          30 :     r = gener_FpXQ(T, p, o);
    1094          30 :     break;
    1095             :   case t_FF_F2xq:
    1096         154 :     r = gener_F2xq(T, o);
    1097         154 :     break;
    1098             :   default:
    1099         208 :     r = gener_Flxq(T, pp, o);
    1100             :   }
    1101         392 :   return _mkFF(x,z,r);
    1102             : }
    1103             : 
    1104             : static GEN
    1105      512575 : to_FFE(GEN P, GEN fg)
    1106             : {
    1107      512575 :   if(ell_is_inf(P))
    1108      232085 :     return ellinf();
    1109             :   else
    1110      280490 :     retmkvec2(mkFF_i(fg,gel(P,1)), mkFF_i(fg,gel(P,2)));
    1111             : }
    1112             : 
    1113             : static GEN
    1114       18109 : to_FFE_vec(GEN x, GEN ff)
    1115             : {
    1116       18109 :   long i, lx = lg(x);
    1117       18109 :   for (i=1; i<lx; i++) gel(x,i) = to_FFE(gel(x,i), ff);
    1118       18109 :   return x;
    1119             : }
    1120             : 
    1121             : static GEN
    1122        1711 : FpXQ_ell_to_a4a6(GEN E, GEN T, GEN p)
    1123             : {
    1124             :   GEN a1, a3, b2, c4, c6;
    1125        1711 :   a1 = Rg_to_FpXQ(ell_get_a1(E),T,p);
    1126        1711 :   a3 = Rg_to_FpXQ(ell_get_a3(E),T,p);
    1127        1711 :   b2 = Rg_to_FpXQ(ell_get_b2(E),T,p);
    1128        1711 :   c4 = Rg_to_FpXQ(ell_get_c4(E),T,p);
    1129        1711 :   c6 = Rg_to_FpXQ(ell_get_c6(E),T,p);
    1130        1711 :   retmkvec3(FpX_neg(FpX_mulu(c4, 27, p), p), FpX_neg(FpX_mulu(c6, 54, p), p),
    1131             :             mkvec4(Z_to_FpX(utoi(6),p,varn(T)),FpX_mulu(b2,3,p),
    1132             :                    FpX_mulu(a1,3,p),FpX_mulu(a3,108,p)));
    1133             : }
    1134             : 
    1135             : static GEN
    1136       24710 : F3xq_ell_to_a4a6(GEN E, GEN T)
    1137             : {
    1138             :   GEN a1, a3, b2, b4, b6;
    1139       24710 :   a1 = Rg_to_Flxq(ell_get_a1(E),T,3);
    1140       24710 :   a3 = Rg_to_Flxq(ell_get_a3(E),T,3);
    1141       24710 :   b2 = Rg_to_Flxq(ell_get_b2(E),T,3);
    1142       24710 :   b4 = Rg_to_Flxq(ell_get_b4(E),T,3);
    1143       24710 :   b6 = Rg_to_Flxq(ell_get_b6(E),T,3);
    1144       24710 :   if(lgpol(b2)) /* ordinary case */
    1145             :   {
    1146       17885 :     GEN b4b2 = Flxq_div(b4,b2,T,3);
    1147       17885 :     GEN a6 = Flx_sub(b6,Flxq_mul(b4b2,Flx_add(b4,Flxq_sqr(b4b2,T,3),3),T,3),3);
    1148       17885 :     retmkvec3(mkvec(b2), a6,
    1149             :        mkvec4(Fl_to_Flx(1,T[1]),b4b2,Flx_neg(a1,3),Flx_neg(a3,3)));
    1150             :   }
    1151             :   else /* super-singular case */
    1152        6825 :     retmkvec3(Flx_neg(b4, 3), b6,
    1153             :        mkvec4(Fl_to_Flx(1,T[1]),zero_Flx(T[1]), Flx_neg(a1,3), Flx_neg(a3,3)));
    1154             : }
    1155             : 
    1156             : static GEN
    1157       66812 : Flxq_ell_to_a4a6(GEN E, GEN T, ulong p)
    1158             : {
    1159             :   GEN a1, a3, b2, c4, c6;
    1160       66812 :   if(p==3) return F3xq_ell_to_a4a6(E, T);
    1161       42102 :   a1 = Rg_to_Flxq(ell_get_a1(E),T,p);
    1162       42102 :   a3 = Rg_to_Flxq(ell_get_a3(E),T,p);
    1163       42102 :   b2 = Rg_to_Flxq(ell_get_b2(E),T,p);
    1164       42102 :   c4 = Rg_to_Flxq(ell_get_c4(E),T,p);
    1165       42102 :   c6 = Rg_to_Flxq(ell_get_c6(E),T,p);
    1166       42102 :   retmkvec3(Flx_neg(Flx_mulu(c4, 27, p), p), Flx_neg(Flx_mulu(c6, 54, p), p),
    1167             :             mkvec4(Fl_to_Flx(6%p,T[1]),Flx_mulu(b2,3,p),
    1168             :                    Flx_mulu(a1,3,p),Flx_mulu(a3,108,p)));
    1169             : }
    1170             : 
    1171             : static GEN
    1172       45662 : F2xq_ell_to_a4a6(GEN E, GEN T)
    1173             : {
    1174       45662 :   long v = T[1];
    1175       45662 :   GEN a1 = Rg_to_F2xq(ell_get_a1(E),T);
    1176       45661 :   GEN a2 = Rg_to_F2xq(ell_get_a2(E),T);
    1177       45661 :   GEN a3 = Rg_to_F2xq(ell_get_a3(E),T);
    1178       45661 :   GEN a4 = Rg_to_F2xq(ell_get_a4(E),T);
    1179       45661 :   GEN a6 = Rg_to_F2xq(ell_get_a6(E),T);
    1180       45661 :   if (lgpol(a1))
    1181             :   {
    1182        7644 :     GEN a1i = F2xq_inv(a1,T);
    1183        7645 :     GEN a1i2 = F2xq_sqr(a1i,T);
    1184        7644 :     GEN a1i3 = F2xq_mul(a1i,a1i2,T);
    1185        7645 :     GEN a1i6 = F2xq_sqr(a1i3,T);
    1186        7645 :     GEN d  = F2xq_mul(a3,a1i,T);
    1187        7645 :     GEN dd = F2xq_mul(d,a1i2,T);
    1188        7645 :     GEN e  = F2xq_mul(F2x_add(a4,F2xq_sqr(d,T)),a1i,T);
    1189        7645 :     GEN ee = F2xq_mul(e,a1i3,T);
    1190        7644 :     GEN da2 = F2x_add(a2,d);
    1191        7644 :     GEN d2 = F2xq_mul(da2,a1i2,T);
    1192        7643 :     GEN d6 = F2xq_mul(F2x_add(F2x_add(F2xq_mul(F2x_add(F2xq_mul(da2,d,T),a4),d,T),a6),F2xq_sqr(e,T)),a1i6,T);
    1193        7644 :     retmkvec3(d2, d6, mkvec4(a1i,dd,pol0_F2x(v),ee));
    1194             :   }
    1195             :   else
    1196             :   { /* allow a1 = a3 = 0: singular curve */
    1197       38017 :     GEN d4 = F2x_add(F2xq_sqr(a2,T),a4);
    1198       38017 :     GEN d6 = F2x_add(F2xq_mul(a2,a4,T),a6);
    1199       38017 :     GEN a3i = lgpol(a3)? F2xq_inv(a3,T): a3;
    1200       38017 :     retmkvec3(mkvec3(a3,d4,a3i), d6, mkvec4(pol1_F2x(v),a2,pol0_F2x(T[1]),pol0_F2x(T[1])));
    1201             :   }
    1202             : }
    1203             : 
    1204             : static GEN
    1205        1744 : FqV_to_FpXQV(GEN x, GEN T)
    1206             : {
    1207        1744 :   pari_sp av = avma;
    1208        1744 :   long v = varn(T), i, s=0, l = lg(x);
    1209        1744 :   GEN y = shallowcopy(x);
    1210        8720 :   for(i=1; i<l; i++)
    1211        6976 :     if (typ(gel(x,i))==t_INT)
    1212             :     {
    1213           0 :       gel(y,i) = scalarpol(gel(x,i), v);
    1214           0 :       s = 1;
    1215             :     }
    1216        1744 :   if (!s) { set_avma(av); return x; }
    1217           0 :   return y;
    1218             : }
    1219             : 
    1220             : GEN
    1221       95670 : FF_ellcard(GEN E)
    1222             : {
    1223             :   GEN T,p;
    1224             :   ulong pp;
    1225       95670 :   GEN e = ellff_get_a4a6(E);
    1226       95670 :   GEN fg = ellff_get_field(E);
    1227       95670 :   _getFF(fg,&T,&p,&pp);
    1228       95670 :   switch(fg[1])
    1229             :   {
    1230             :   case t_FF_FpXQ:
    1231        1695 :     return FpXQ_ellcard(Fq_to_FpXQ(gel(e,1),T,p), Fq_to_FpXQ(gel(e,2),T,p),T,p);
    1232             :   case t_FF_F2xq:
    1233       45410 :     return F2xq_ellcard(gel(e,1),gel(e,2),T);
    1234             :   default:
    1235       48565 :     return Flxq_ellcard(gel(e,1),gel(e,2),T,pp);
    1236             :   }
    1237             : }
    1238             : 
    1239             : GEN
    1240          14 : FF_ellcard_SEA(GEN E, long smallfact)
    1241             : {
    1242          14 :   pari_sp av = avma;
    1243             :   GEN T,p;
    1244             :   ulong pp;
    1245          14 :   GEN e = ellff_get_a4a6(E), a4, a6, card;
    1246          14 :   GEN fg = ellff_get_field(E);
    1247          14 :   _getFF(fg,&T,&p,&pp);
    1248          14 :   switch(fg[1])
    1249             :   {
    1250             :   case t_FF_FpXQ:
    1251           0 :     a4 = Fq_to_FpXQ(gel(e,1), T, p);
    1252           0 :     a6 = Fq_to_FpXQ(gel(e,2), T, p);
    1253           0 :     card = Fq_ellcard_SEA(a4, a6, powiu(p,degpol(T)), T,p, smallfact);
    1254           0 :     break;
    1255             :   case t_FF_F2xq:
    1256           0 :     pari_err_IMPL("SEA for char 2");
    1257             :   default:
    1258          14 :     a4 = Flx_to_ZX(gel(e,1));
    1259          14 :     a6 = Flx_to_ZX(gel(e,2));
    1260          14 :     card = Fq_ellcard_SEA(a4, a6, powuu(pp,degpol(T)), Flx_to_ZX(T), p, smallfact);
    1261             :   }
    1262          14 :   return gerepileuptoint(av, card);
    1263             : }
    1264             : 
    1265             : /* return G, set m */
    1266             : GEN
    1267       19432 : FF_ellgroup(GEN E, GEN *pm)
    1268             : {
    1269             :   GEN T, p, e, fg, N;
    1270             :   ulong pp;
    1271             : 
    1272       19432 :   N = ellff_get_card(E);
    1273       19432 :   e = ellff_get_a4a6(E);
    1274       19432 :   fg = ellff_get_field(E);
    1275       19432 :   _getFF(fg,&T,&p,&pp);
    1276       19432 :   switch(fg[1])
    1277             :   {
    1278             :   case t_FF_FpXQ:
    1279          30 :     return FpXQ_ellgroup(Fq_to_FpXQ(gel(e,1),T,p),
    1280          15 :                          Fq_to_FpXQ(gel(e,2),T,p),N,T,p,pm);
    1281             :   case t_FF_F2xq:
    1282        3808 :     return F2xq_ellgroup(gel(e,1),gel(e,2),N,T,pm);
    1283             :   default:
    1284       15609 :     return Flxq_ellgroup(gel(e,1),gel(e,2),N,T,pp,pm);
    1285             :   }
    1286             : }
    1287             : 
    1288             : GEN
    1289       18109 : FF_ellgens(GEN E)
    1290             : {
    1291             :   GEN D, fg, e, m, T, p, F, e3;
    1292             :   ulong pp;
    1293             : 
    1294       18109 :   fg = ellff_get_field(E);
    1295       18109 :   e = ellff_get_a4a6(E);
    1296       18109 :   m = ellff_get_m(E);
    1297       18109 :   D = ellff_get_D(E);
    1298       18109 :   _getFF(fg,&T,&p,&pp);
    1299       18109 :   switch(fg[1])
    1300             :   {
    1301             :   case t_FF_FpXQ:
    1302           8 :     e3 = FqV_to_FpXQV(gel(e,3),T);
    1303           8 :     F = FpXQ_ellgens(Fq_to_FpXQ(gel(e,1),T,p),Fq_to_FpXQ(gel(e,2),T,p),e3,D,m,T,p);
    1304           8 :     break;
    1305             :   case t_FF_F2xq:
    1306        3738 :     F = F2xq_ellgens(gel(e,1),gel(e,2),gel(e,3),D,m,T);
    1307        3738 :     break;
    1308             :   default:
    1309       14363 :     F = Flxq_ellgens(gel(e,1),gel(e,2),gel(e,3),D,m,T, pp);
    1310       14363 :     break;
    1311             :   }
    1312       18109 :   return to_FFE_vec(F,fg);
    1313             : }
    1314             : 
    1315             : GEN
    1316         119 : FF_elldata(GEN E, GEN fg)
    1317             : {
    1318             :   GEN T,p,e;
    1319             :   ulong pp;
    1320         119 :   _getFF(fg,&T,&p,&pp);
    1321         119 :   switch(fg[1])
    1322             :   {
    1323             :   case t_FF_FpXQ:
    1324           0 :     e = FpXQ_ell_to_a4a6(E,T,p); break;
    1325             :   case t_FF_F2xq:
    1326          56 :     e = F2xq_ell_to_a4a6(E,T); break;
    1327             :   default:
    1328          63 :     e = Flxq_ell_to_a4a6(E,T,pp); break;
    1329             :   }
    1330         119 :   return mkvec2(fg,e);
    1331             : }
    1332             : 
    1333             : /* allow singular E, set E.j = 0 in this case */
    1334             : GEN
    1335      114066 : FF_ellinit(GEN E, GEN fg)
    1336             : {
    1337      114066 :   GEN T,p,e, D, c4, y = E;
    1338             :   ulong pp;
    1339             :   long i;
    1340      114066 :   _getFF(fg,&T,&p,&pp);
    1341      114066 :   switch(fg[1])
    1342             :   {
    1343             :   case t_FF_FpXQ:
    1344        1711 :     e = FpXQ_ell_to_a4a6(E,T,p);
    1345        1711 :     for(i=1;i<=12;i++) gel(y,i) = mkFF_i(fg,Rg_to_FpXQ(gel(E,i),T,p));
    1346        1711 :     break;
    1347             :   case t_FF_F2xq:
    1348       45606 :     e = F2xq_ell_to_a4a6(E,T);
    1349       45606 :     for(i=1;i<=12;i++) gel(y,i) = mkFF_i(fg,Rg_to_F2xq(gel(E,i),T));
    1350       45606 :     break;
    1351             :   default:
    1352       66749 :     e = Flxq_ell_to_a4a6(E,T,pp);
    1353       66749 :     for(i=1;i<=12;i++) gel(y,i) = mkFF_i(fg,Rg_to_Flxq(gel(E,i),T,pp));
    1354       66749 :     break;
    1355             :   }
    1356      114066 :   D = ell_get_disc(y);
    1357      114066 :   c4 = ell_get_c4(y);
    1358      114066 :   gel(y,13) = FF_equal0(D)? D: gdiv(gpowgs(c4,3), D);
    1359      114066 :   gel(y,14) = mkvecsmall(t_ELL_Fq);
    1360      114066 :   gel(y,15) = mkvec2(fg,e);
    1361      114066 :   return y;
    1362             : }
    1363             : 
    1364             : GEN
    1365       27188 : FF_elltwist(GEN E)
    1366             : {
    1367       27188 :   pari_sp av = avma;
    1368       27188 :   GEN fg = ellff_get_field(E), e = ellff_get_a4a6(E);
    1369             :   GEN T, p, a, a6, V;
    1370             :   ulong pp;
    1371       27188 :   _getFF(fg,&T,&p,&pp);
    1372       27188 :   switch (fg[1])
    1373             :   {
    1374             :   case t_FF_FpXQ:
    1375         840 :     FpXQ_elltwist(gel(e,1), gel(e,2), T, p, &a, &a6);
    1376         840 :     V = mkvec5(gen_0, gen_0, gen_0, mkFF_i(fg, a), mkFF_i(fg, a6));
    1377         840 :     break;
    1378             :   case t_FF_F2xq:
    1379        3514 :     F2xq_elltwist(gel(e,1), gel(e,2), T, &a, &a6);
    1380        7028 :     V = typ(a)==t_VECSMALL ?
    1381        3514 :           mkvec5(gen_1, mkFF_i(fg, a), gen_0, gen_0, mkFF_i(fg, a6)):
    1382           0 :           mkvec5(gen_0, gen_0, mkFF_i(fg, gel(a,1)), mkFF_i(fg, gel(a,2)), mkFF_i(fg, a6));
    1383        3514 :     break;
    1384             :   default:
    1385       22834 :     Flxq_elltwist(gel(e,1), gel(e,2), T, pp, &a, &a6);
    1386       45668 :     V = typ(a)==t_VECSMALL ?
    1387       30436 :           mkvec5(gen_0, gen_0, gen_0, mkFF_i(fg, a), mkFF_i(fg, a6)):
    1388        7602 :           mkvec5(gen_0, mkFF_i(fg, gel(a,1)), gen_0, gen_0, mkFF_i(fg, a6));
    1389             :   }
    1390       27188 :   return gerepilecopy(av, V);
    1391             : }
    1392             : 
    1393             : static long
    1394           0 : F3x_equalm1(GEN x)
    1395           0 : { return degpol(x)==0 && x[2] == 2; }
    1396             : 
    1397             : GEN
    1398      244314 : FF_ellrandom(GEN E)
    1399             : {
    1400      244314 :   pari_sp av = avma;
    1401      244314 :   GEN fg = ellff_get_field(E), e = ellff_get_a4a6(E), Q;
    1402             :   GEN T,p;
    1403             :   ulong pp;
    1404      244314 :   _getFF(fg,&T,&p,&pp);
    1405      244314 :   switch (fg[1])
    1406             :   {
    1407             :   case t_FF_FpXQ:
    1408         854 :     Q = random_FpXQE(Fq_to_FpXQ(gel(e,1),T,p), Fq_to_FpXQ(gel(e,2),T,p), T, p);
    1409         854 :     Q = FpXQE_changepoint(Q, FqV_to_FpXQV(gel(e,3), T) , T, p);
    1410         854 :     break;
    1411             :   case t_FF_F2xq:
    1412             :     {
    1413      168336 :       long d = F2x_degree(T);
    1414             :       /* if #E(Fq) = 1 return [0] */
    1415      168336 :       if (d<=2 && typ(gel(e,1)) == t_VEC)
    1416             :       { /* over F2 or F4, supersingular */
    1417        1687 :         GEN v = gel(e,1), A6 = gel(e,2), a3 = gel(v,1), A4 = gel(v,2);
    1418        1687 :         if (F2x_equal1(a3) &&
    1419         196 :              ((d==1 && F2x_equal1(A4) && F2x_equal1(A6))
    1420         595 :            || (d==2 && !lgpol(A4)     && F2x_degree(A6)==1))) return ellinf();
    1421             :       }
    1422      168252 :       Q = random_F2xqE(gel(e,1), gel(e,2), T);
    1423      168252 :       Q = F2xqE_changepoint(Q, gel(e,3), T);
    1424      168252 :       break;
    1425             :     }
    1426             :   default:
    1427             :     /* if #E(Fq) = 1 return [0] */
    1428       75124 :     if (pp==3 && degpol(T)==1 && typ(gel(e,1))==t_VECSMALL)
    1429             :     { /* over F3, supersingular */
    1430           0 :       GEN mb4 = gel(e,1), b6 = gel(e,2);
    1431           0 :       if (F3x_equalm1(mb4) && F3x_equalm1(b6)) return ellinf();
    1432             :     }
    1433       75124 :     Q = random_FlxqE(gel(e,1), gel(e,2), T, pp);
    1434       75124 :     Q = FlxqE_changepoint(Q, gel(e,3), T, pp);
    1435             :   }
    1436      244230 :   return gerepilecopy(av, to_FFE(Q, fg));
    1437             : }
    1438             : 
    1439             : GEN
    1440      247702 : FF_ellmul(GEN E, GEN P, GEN n)
    1441             : {
    1442      247702 :   pari_sp av = avma;
    1443      247702 :   GEN fg = ellff_get_field(E), e = ellff_get_a4a6(E), Q;
    1444             :   GEN T,p, Pp, Qp, e3;
    1445             :   ulong pp;
    1446      247702 :   _getFF(fg,&T,&p,&pp);
    1447      247702 :   switch (fg[1])
    1448             :   {
    1449             :   case t_FF_FpXQ:
    1450         854 :     e3 = FqV_to_FpXQV(gel(e,3),T);
    1451         854 :     Pp = FpXQE_changepointinv(RgE_to_FpXQE(P, T, p), e3, T, p);
    1452         854 :     Qp = FpXQE_mul(Pp, n, gel(e,1), T, p);
    1453         854 :     Q = FpXQE_changepoint(Qp, gel(e,3), T, p);
    1454         854 :     break;
    1455             :   case t_FF_F2xq:
    1456      184604 :     Pp = F2xqE_changepointinv(RgE_to_F2xqE(P, T), gel(e,3), T);
    1457      184604 :     Qp = F2xqE_mul(Pp, n, gel(e,1), T);
    1458      184604 :     Q = F2xqE_changepoint(Qp, gel(e,3), T);
    1459      184604 :     break;
    1460             :   default:
    1461       62244 :     Pp = FlxqE_changepointinv(RgE_to_FlxqE(P, T, pp), gel(e,3), T, pp);
    1462       62244 :     Qp = FlxqE_mul(Pp, n, gel(e,1), T, pp);
    1463       62244 :     Q = FlxqE_changepoint(Qp, gel(e,3), T, pp);
    1464             :   }
    1465      247702 :   return gerepilecopy(av, to_FFE(Q, fg));
    1466             : }
    1467             : 
    1468             : GEN
    1469        1750 : FF_ellorder(GEN E, GEN P, GEN o)
    1470             : {
    1471        1750 :   pari_sp av = avma;
    1472        1750 :   GEN fg = ellff_get_field(E), e = ellff_get_a4a6(E);
    1473             :   GEN r,T,p,Pp,e3;
    1474             :   ulong pp;
    1475        1750 :   _getFF(fg,&T,&p,&pp);
    1476        1750 :   switch (fg[1])
    1477             :   {
    1478             :   case t_FF_FpXQ:
    1479          14 :     e3 = FqV_to_FpXQV(gel(e,3),T);
    1480          14 :     Pp = FpXQE_changepointinv(RgE_to_FpXQE(P,T,p), e3, T, p);
    1481          14 :     r = FpXQE_order(Pp, o, gel(e,1), T, p);
    1482          14 :     break;
    1483             :   case t_FF_F2xq:
    1484         280 :     Pp = F2xqE_changepointinv(RgE_to_F2xqE(P,T), gel(e,3), T);
    1485         280 :     r = F2xqE_order(Pp, o, gel(e,1), T);
    1486         280 :     break;
    1487             :   default:
    1488        1456 :     Pp = FlxqE_changepointinv(RgE_to_FlxqE(P,T,pp), gel(e,3), T, pp);
    1489        1456 :     r = FlxqE_order(Pp, o, gel(e,1), T, pp);
    1490             :   }
    1491        1750 :   return gerepileupto(av, r);
    1492             : }
    1493             : 
    1494             : GEN
    1495          91 : FF_elllog(GEN E, GEN P, GEN Q, GEN o)
    1496             : {
    1497          91 :   pari_sp av = avma;
    1498          91 :   GEN fg = ellff_get_field(E), e = ellff_get_a4a6(E);
    1499             :   GEN r,T,p, Pp,Qp, e3;
    1500             :   ulong pp;
    1501          91 :   _getFF(fg,&T,&p,&pp);
    1502          91 :   switch(fg[1])
    1503             :   {
    1504             :   case t_FF_FpXQ:
    1505           0 :     e3 = FqV_to_FpXQV(gel(e,3),T);
    1506           0 :     Pp = FpXQE_changepointinv(RgE_to_FpXQE(P,T,p), e3, T, p);
    1507           0 :     Qp = FpXQE_changepointinv(RgE_to_FpXQE(Q,T,p), e3, T, p);
    1508           0 :     r = FpXQE_log(Pp, Qp, o, gel(e,1), T, p);
    1509           0 :     break;
    1510             :   case t_FF_F2xq:
    1511          42 :     Pp = F2xqE_changepointinv(RgE_to_F2xqE(P,T), gel(e,3), T);
    1512          42 :     Qp = F2xqE_changepointinv(RgE_to_F2xqE(Q,T), gel(e,3), T);
    1513          42 :     r = F2xqE_log(Pp, Qp, o, gel(e,1), T);
    1514          42 :     break;
    1515             :   default:
    1516          49 :     Pp = FlxqE_changepointinv(RgE_to_FlxqE(P,T,pp), gel(e,3), T, pp);
    1517          49 :     Qp = FlxqE_changepointinv(RgE_to_FlxqE(Q,T,pp), gel(e,3), T, pp);
    1518          49 :     r = FlxqE_log(Pp, Qp, o, gel(e,1), T, pp);
    1519             :   }
    1520          91 :   return gerepileupto(av, r);
    1521             : }
    1522             : 
    1523             : GEN
    1524        4984 : FF_ellweilpairing(GEN E, GEN P, GEN Q, GEN m)
    1525             : {
    1526        4984 :   GEN fg = ellff_get_field(E), e = ellff_get_a4a6(E);
    1527             :   GEN r,T,p, Pp,Qp, e3;
    1528             :   ulong pp;
    1529        4984 :   GEN z=_initFF(fg,&T,&p,&pp);
    1530        4984 :   pari_sp av = avma;
    1531        4984 :   switch(fg[1])
    1532             :   {
    1533             :   case t_FF_FpXQ:
    1534           7 :     e3 = FqV_to_FpXQV(gel(e,3),T);
    1535           7 :     Pp = FpXQE_changepointinv(RgE_to_FpXQE(P,T,p), e3, T, p);
    1536           7 :     Qp = FpXQE_changepointinv(RgE_to_FpXQE(Q,T,p), e3, T, p);
    1537           7 :     r = FpXQE_weilpairing(Pp, Qp, m, gel(e,1), T, p);
    1538           7 :     break;
    1539             :   case t_FF_F2xq:
    1540        4963 :     Pp = F2xqE_changepointinv(RgE_to_F2xqE(P,T), gel(e,3), T);
    1541        4963 :     Qp = F2xqE_changepointinv(RgE_to_F2xqE(Q,T), gel(e,3), T);
    1542        4963 :     r = F2xqE_weilpairing(Pp, Qp, m, gel(e,1), T);
    1543        4963 :     break;
    1544             :   default:
    1545          14 :     Pp = FlxqE_changepointinv(RgE_to_FlxqE(P,T,pp), gel(e,3), T, pp);
    1546          14 :     Qp = FlxqE_changepointinv(RgE_to_FlxqE(Q,T,pp), gel(e,3), T, pp);
    1547          14 :     r = FlxqE_weilpairing(Pp, Qp, m, gel(e,1), T, pp);
    1548             :   }
    1549        4984 :   r = gerepileupto(av, r);
    1550        4984 :   return _mkFF(fg,z,r);
    1551             : }
    1552             : 
    1553             : GEN
    1554          91 : FF_elltatepairing(GEN E, GEN P, GEN Q, GEN m)
    1555             : {
    1556          91 :   GEN fg = ellff_get_field(E), e = ellff_get_a4a6(E);
    1557             :   GEN r,T,p, Pp,Qp, e3;
    1558             :   ulong pp;
    1559          91 :   GEN z=_initFF(fg,&T,&p,&pp);
    1560          91 :   pari_sp av = avma;
    1561          91 :   switch(fg[1])
    1562             :   {
    1563             :   case t_FF_FpXQ:
    1564           7 :     e3 = FqV_to_FpXQV(gel(e,3),T);
    1565           7 :     Pp = FpXQE_changepointinv(RgE_to_FpXQE(P,T,p), e3, T, p);
    1566           7 :     Qp = FpXQE_changepointinv(RgE_to_FpXQE(Q,T,p), e3, T, p);
    1567           7 :     r = FpXQE_tatepairing(Pp, Qp, m, gel(e,1), T, p);
    1568           7 :     break;
    1569             :   case t_FF_F2xq:
    1570          28 :     Pp = F2xqE_changepointinv(RgE_to_F2xqE(P,T), gel(e,3), T);
    1571          28 :     Qp = F2xqE_changepointinv(RgE_to_F2xqE(Q,T), gel(e,3), T);
    1572          28 :     r = F2xqE_tatepairing(Pp, Qp, m, gel(e,1), T);
    1573          28 :     break;
    1574             :   default:
    1575          56 :     Pp = FlxqE_changepointinv(RgE_to_FlxqE(P,T,pp), gel(e,3), T, pp);
    1576          56 :     Qp = FlxqE_changepointinv(RgE_to_FlxqE(Q,T,pp), gel(e,3), T, pp);
    1577          56 :     r = FlxqE_tatepairing(Pp, Qp, m, gel(e,1), T, pp);
    1578             :   }
    1579          91 :   r = gerepileupto(av, r);
    1580          91 :   return _mkFF(fg,z,r);
    1581             : }
    1582             : 
    1583             : GEN
    1584      103649 : FFX_roots(GEN Pf, GEN ff)
    1585             : {
    1586      103649 :   pari_sp av = avma;
    1587             :   GEN r,T,p;
    1588             :   ulong pp;
    1589      103649 :   GEN P = FFX_to_raw(Pf, ff);
    1590      103649 :   _getFF(ff,&T,&p,&pp);
    1591      103649 :   switch(ff[1])
    1592             :   {
    1593             :   case t_FF_FpXQ:
    1594          85 :     r = FpXQX_roots(P, T, p);
    1595          85 :     break;
    1596             :   case t_FF_F2xq:
    1597       57512 :     r = F2xqX_roots(P, T);
    1598       57512 :     break;
    1599             :   default:
    1600       46052 :     r = FlxqX_roots(P, T, pp);
    1601             :   }
    1602      103649 :   return gerepilecopy(av, raw_to_FFC(r, ff));
    1603             : }
    1604             : 
    1605             : static GEN
    1606           7 : raw_to_FFXC(GEN x, GEN ff)
    1607           7 : { pari_APPLY_type(t_COL, raw_to_FFX(gel(x,i), ff)); }
    1608             : 
    1609             : static GEN
    1610         448 : raw_to_FFX_fact(GEN F, GEN ff)
    1611             : {
    1612             :   GEN y, u, v;
    1613         448 :   GEN P = gel(F,1), E = gel(F,2);
    1614         448 :   long j, l = lg(P);
    1615         448 :   y = cgetg(3,t_MAT);
    1616         448 :   u = cgetg(l,t_COL); gel(y,1) = u;
    1617         448 :   v = cgetg(l,t_COL); gel(y,2) = v;
    1618        2107 :   for (j=1; j<l; j++)
    1619             :   {
    1620        1659 :     gel(u,j) = raw_to_FFX(gel(P,j), ff);
    1621        1659 :     gel(v,j) = utoi(uel(E,j));
    1622             :   }
    1623         448 :   return y;
    1624             : }
    1625             : 
    1626             : static GEN
    1627        2244 : FFX_zero(GEN ff, long v)
    1628             : {
    1629        2244 :   GEN r = cgetg(3,t_POL);
    1630        2244 :   r[1] = evalvarn(v);
    1631        2244 :   gel(r,2) = FF_zero(ff);
    1632        2244 :   return r;
    1633             : }
    1634             : 
    1635             : static GEN
    1636      107332 : FFX_wrap2(GEN Pf, GEN Qf, GEN ff, GEN FpXQX(GEN, GEN, GEN, GEN),
    1637             :           GEN F2xqX(GEN, GEN, GEN), GEN FlxqX(GEN, GEN, GEN, ulong))
    1638             : {
    1639      107332 :   pari_sp av = avma;
    1640             :   GEN r,T,p;
    1641             :   ulong pp;
    1642      107332 :   GEN P = FFX_to_raw(Pf, ff);
    1643      107332 :   GEN Q = FFX_to_raw(Qf, ff);
    1644      107332 :   _getFF(ff,&T,&p,&pp);
    1645      107332 :   switch(ff[1])
    1646             :   {
    1647             :   case t_FF_FpXQ:
    1648        2374 :     r = FpXQX(P, Q, T, p);
    1649        2374 :     break;
    1650             :   case t_FF_F2xq:
    1651       39931 :     r = F2xqX(P, Q, T);
    1652       39931 :     break;
    1653             :   default:
    1654       65027 :     r = FlxqX(P, Q, T, pp);
    1655             :   }
    1656      107332 :   if (!lgpol(r)) { set_avma(av); return FFX_zero(ff, varn(Pf)); }
    1657      105088 :   return gerepilecopy(av, raw_to_FFX(r, ff));
    1658             : }
    1659             : 
    1660             : GEN
    1661      104161 : FFX_mul(GEN Pf, GEN Qf, GEN ff)
    1662      104161 : { return FFX_wrap2(Pf, Qf, ff, FpXQX_mul, F2xqX_mul, FlxqX_mul); }
    1663             : 
    1664             : GEN
    1665        2541 : FFX_gcd(GEN Pf, GEN Qf, GEN ff)
    1666        2541 : { return FFX_wrap2(Pf, Qf, ff, FpXQX_gcd, F2xqX_gcd, FlxqX_gcd); }
    1667             : 
    1668             : GEN
    1669        2296 : FFX_sqr(GEN Pf, GEN ff)
    1670             : {
    1671        2296 :   pari_sp av = avma;
    1672             :   GEN r,T,p;
    1673             :   ulong pp;
    1674        2296 :   GEN P = FFX_to_raw(Pf, ff);
    1675        2296 :   _getFF(ff,&T,&p,&pp);
    1676        2296 :   switch(ff[1])
    1677             :   {
    1678             :   case t_FF_FpXQ:
    1679         210 :     r = FpXQX_sqr(P, T, p);
    1680         210 :     break;
    1681             :   case t_FF_F2xq:
    1682         966 :     r = F2xqX_sqr(P, T);
    1683         966 :     break;
    1684             :   default:
    1685        1120 :     r = FlxqX_sqr(P, T, pp);
    1686             :   }
    1687        2296 :   if (!lgpol(r)) { set_avma(av); return FFX_zero(ff, varn(Pf)); }
    1688        2296 :   return gerepilecopy(av, raw_to_FFX(r, ff));
    1689             : }
    1690             : 
    1691             : GEN
    1692         616 : FFX_rem(GEN Pf, GEN Qf, GEN ff)
    1693         616 : { return FFX_wrap2(Pf, Qf, ff, FpXQX_rem, F2xqX_rem, FlxqX_rem); }
    1694             : 
    1695             : GEN
    1696          21 : FFX_resultant(GEN Pf, GEN Qf, GEN ff)
    1697             : {
    1698          21 :   pari_sp av = avma;
    1699             :   GEN r,T,p;
    1700             :   ulong pp;
    1701          21 :   GEN P = FFX_to_raw(Pf, ff);
    1702          21 :   GEN Q = FFX_to_raw(Qf, ff);
    1703          21 :   GEN z = _initFF(ff,&T,&p,&pp);
    1704          21 :   switch(ff[1])
    1705             :   {
    1706             :   case t_FF_FpXQ:
    1707           7 :     r = FpXQX_resultant(P, Q, T, p);
    1708           7 :     break;
    1709             :   case t_FF_F2xq:
    1710           7 :     r = F2xqX_resultant(P, Q, T);
    1711           7 :     break;
    1712             :   default:
    1713           7 :     r = FlxqX_resultant(P, Q, T, pp);
    1714             :   }
    1715          21 :   return gerepileupto(av, _mkFF(ff,z,r));
    1716             : }
    1717             : 
    1718             : GEN
    1719          35 : FFX_disc(GEN Pf, GEN ff)
    1720             : {
    1721          35 :   pari_sp av = avma;
    1722             :   GEN r,T,p;
    1723             :   ulong pp;
    1724          35 :   GEN P = FFX_to_raw(Pf, ff);
    1725          35 :   GEN z = _initFF(ff,&T,&p,&pp);
    1726          35 :   switch(ff[1])
    1727             :   {
    1728             :   case t_FF_FpXQ:
    1729           7 :     r = FpXQX_disc(P, T, p);
    1730           7 :     break;
    1731             :   case t_FF_F2xq:
    1732          14 :     r = F2xqX_disc(P, T);
    1733          14 :     break;
    1734             :   default:
    1735          14 :     r = FlxqX_disc(P, T, pp);
    1736             :   }
    1737          35 :   return gerepileupto(av, _mkFF(ff,z,r));
    1738             : }
    1739             : 
    1740             : static GEN
    1741          21 : gc_gcdext(pari_sp av, GEN r, GEN *u, GEN *v)
    1742             : {
    1743          21 :   if (!u && !v) return gerepilecopy(av, r);
    1744          21 :   if (u  &&  v) gerepileall(av, 3, &r, u, v);
    1745           0 :   else          gerepileall(av, 2, &r, u ? u: v);
    1746          21 :   return r;
    1747             : }
    1748             : 
    1749             : GEN
    1750          21 : FFX_extgcd(GEN Pf, GEN Qf, GEN ff, GEN *pt_Uf, GEN *pt_Vf)
    1751             : {
    1752          21 :   pari_sp av = avma;
    1753             :   GEN r,T,p;
    1754             :   ulong pp;
    1755          21 :   GEN P = FFX_to_raw(Pf, ff);
    1756          21 :   GEN Q = FFX_to_raw(Qf, ff);
    1757          21 :   _getFF(ff,&T,&p,&pp);
    1758          21 :   switch(ff[1])
    1759             :   {
    1760             :   case t_FF_FpXQ:
    1761           0 :     r = FpXQX_extgcd(P, Q, T, p, pt_Uf, pt_Vf);
    1762           0 :     break;
    1763             :   case t_FF_F2xq:
    1764           7 :     r = F2xqX_extgcd(P, Q, T, pt_Uf, pt_Vf);
    1765           7 :     break;
    1766             :   default:
    1767          14 :     r = FlxqX_extgcd(P, Q, T, pp, pt_Uf, pt_Vf);
    1768             :   }
    1769          21 :   if (pt_Uf) *pt_Uf = raw_to_FFX(*pt_Uf, ff);
    1770          21 :   if (pt_Vf) *pt_Vf = raw_to_FFX(*pt_Vf, ff);
    1771          21 :   return gc_gcdext(av, raw_to_FFX(r, ff), pt_Uf, pt_Vf);
    1772             : }
    1773             : 
    1774             : GEN
    1775           0 : FFXQ_sqr(GEN Pf, GEN Qf, GEN ff)
    1776           0 : { return FFX_wrap2(Pf, Qf, ff, FpXQXQ_sqr, F2xqXQ_sqr, FlxqXQ_sqr); }
    1777             : 
    1778             : GEN
    1779          14 : FFXQ_inv(GEN Pf, GEN Qf, GEN ff)
    1780          14 : { return FFX_wrap2(Pf, Qf, ff, FpXQXQ_inv, F2xqXQ_inv, FlxqXQ_inv); }
    1781             : 
    1782             : GEN
    1783          14 : FFXQ_mul(GEN Pf, GEN Qf, GEN Sf, GEN ff)
    1784             : {
    1785          14 :   pari_sp av = avma;
    1786             :   GEN r,T,p;
    1787             :   ulong pp;
    1788          14 :   GEN P = FFX_to_raw(Pf, ff);
    1789          14 :   GEN Q = FFX_to_raw(Qf, ff);
    1790          14 :   GEN S = FFX_to_raw(Sf, ff);
    1791          14 :   _getFF(ff,&T,&p,&pp);
    1792          14 :   switch(ff[1])
    1793             :   {
    1794             :   case t_FF_FpXQ:
    1795           0 :     r = FpXQXQ_mul(P, Q, S, T, p);
    1796           0 :     break;
    1797             :   case t_FF_F2xq:
    1798           7 :     r = F2xqXQ_mul(P, Q, S, T);
    1799           7 :     break;
    1800             :   default:
    1801           7 :     r = FlxqXQ_mul(P, Q, S, T, pp);
    1802             :   }
    1803          14 :   if (!lgpol(r)) { set_avma(av); return FFX_zero(ff, varn(Pf)); }
    1804          14 :   return gerepilecopy(av, raw_to_FFX(r, ff));
    1805             : }
    1806             : 
    1807             : long
    1808         168 : FFX_ispower(GEN Pf, long k, GEN ff, GEN *pt_r)
    1809             : {
    1810         168 :   pari_sp av = avma;
    1811             :   GEN P,T,p;
    1812             :   ulong pp;
    1813             :   long s;
    1814         168 :   if (degpol(Pf) % k) return 0;
    1815         168 :   P = FFX_to_raw(Pf, ff);
    1816         168 :   _getFF(ff,&T,&p,&pp);
    1817         168 :   switch(ff[1])
    1818             :   {
    1819             :   case t_FF_FpXQ:
    1820          56 :     s = FpXQX_ispower(P, k, T, p, pt_r);
    1821          56 :     break;
    1822             :   case t_FF_F2xq:
    1823          56 :     s = F2xqX_ispower(P, k, T, pt_r);
    1824          56 :     break;
    1825             :   default:
    1826          56 :     s = FlxqX_ispower(P, k, T, pp, pt_r);
    1827             :   }
    1828         168 :   if (s==0) return gc_long(av,0);
    1829         147 :   if (pt_r)
    1830         147 :     *pt_r = gerepilecopy(av, raw_to_FFX(*pt_r, ff));
    1831           0 :   else set_avma(av);
    1832         147 :   return 1;
    1833             : }
    1834             : 
    1835             : GEN
    1836         441 : FFX_factor(GEN Pf, GEN ff)
    1837             : {
    1838         441 :   pari_sp av = avma;
    1839             :   GEN r,T,p;
    1840             :   ulong pp;
    1841         441 :   GEN P = FFX_to_raw(Pf, ff);
    1842         441 :   _getFF(ff,&T,&p,&pp);
    1843         441 :   switch(ff[1])
    1844             :   {
    1845             :   case t_FF_FpXQ:
    1846         121 :     r = FpXQX_factor(P, T, p);
    1847         121 :     break;
    1848             :   case t_FF_F2xq:
    1849         133 :     r = F2xqX_factor(P, T);
    1850         133 :     break;
    1851             :   default:
    1852         187 :     r = FlxqX_factor(P, T, pp);
    1853             :   }
    1854         441 :   return gerepilecopy(av, raw_to_FFX_fact(r, ff));
    1855             : }
    1856             : 
    1857             : GEN
    1858           7 : FFX_factor_squarefree(GEN Pf, GEN ff)
    1859             : {
    1860           7 :   pari_sp av = avma;
    1861             :   GEN r,T,p;
    1862             :   ulong pp;
    1863           7 :   GEN P = FFX_to_raw(Pf, ff);
    1864           7 :   _getFF(ff,&T,&p,&pp);
    1865           7 :   switch(ff[1])
    1866             :   {
    1867             :   case t_FF_FpXQ:
    1868           7 :     r = FpXQX_factor_squarefree(P, T, p);
    1869           7 :     break;
    1870             :   case t_FF_F2xq:
    1871           0 :     r = F2xqX_factor_squarefree(P, T);
    1872           0 :     break;
    1873             :   default:
    1874           0 :     r = FlxqX_factor_squarefree(P, T, pp);
    1875             :   }
    1876           7 :   return gerepilecopy(av, raw_to_FFXC(r, ff));
    1877             : }
    1878             : 
    1879             : GEN
    1880           7 : FFX_ddf(GEN Pf, GEN ff)
    1881             : {
    1882           7 :   pari_sp av = avma;
    1883             :   GEN r,T,p;
    1884             :   ulong pp;
    1885           7 :   GEN P = FFX_to_raw(Pf, ff);
    1886           7 :   _getFF(ff,&T,&p,&pp);
    1887           7 :   switch(ff[1])
    1888             :   {
    1889             :   case t_FF_FpXQ:
    1890           7 :     r = FpXQX_ddf(P, T, p);
    1891           7 :     break;
    1892             :   case t_FF_F2xq:
    1893           0 :     r = F2xqX_ddf(P, T);
    1894           0 :     break;
    1895             :   default:
    1896           0 :     r = FlxqX_ddf(P, T, pp);
    1897             :   }
    1898           7 :   return gerepilecopy(av, raw_to_FFX_fact(r, ff));
    1899             : }
    1900             : 
    1901             : GEN
    1902         126 : FFX_degfact(GEN Pf, GEN ff)
    1903             : {
    1904         126 :   pari_sp av = avma;
    1905             :   GEN r,T,p;
    1906             :   ulong pp;
    1907         126 :   GEN P = FFX_to_raw(Pf, ff);
    1908         126 :   _getFF(ff, &T, &p, &pp);
    1909         126 :   switch(ff[1])
    1910             :   {
    1911             :   case t_FF_FpXQ:
    1912          42 :     r = FpXQX_degfact(P, T, p);
    1913          42 :     break;
    1914             :   case t_FF_F2xq:
    1915          42 :     r = F2xqX_degfact(P, T);
    1916          42 :     break;
    1917             :   default:
    1918          42 :     r = FlxqX_degfact(P, T, pp);
    1919             :   }
    1920         126 :   return gerepilecopy(av, r);
    1921             : }
    1922             : 
    1923             : GEN
    1924           0 : FqX_to_FFX(GEN x, GEN ff)
    1925             : {
    1926             :   long i, lx;
    1927           0 :   GEN y =  cgetg_copy(x,&lx);
    1928           0 :   y[1] = x[1];
    1929           0 :   for (i=2; i<lx; i++) gel(y,i) = Fq_to_FF(gel(x,i), ff);
    1930           0 :   return y;
    1931             : }
    1932             : 
    1933             : GEN
    1934        2412 : ffgen(GEN T, long v)
    1935             : {
    1936        2412 :   GEN A, p = NULL, ff = cgetg(5,t_FFELT);
    1937             :   long d;
    1938        2412 :   switch(typ(T))
    1939             :   {
    1940             :     case t_FFELT:
    1941          28 :       p = FF_p_i(T); T = FF_mod(T); d = degpol(T);
    1942          28 :       break;
    1943             :     case t_POL:
    1944         252 :       d = degpol(T); p = NULL;
    1945         252 :       if (d < 1 || !RgX_is_FpX(T, &p) || !p) pari_err_TYPE("ffgen",T);
    1946         252 :       T = RgX_to_FpX(T, p);
    1947             :       /* testing for irreducibility is too costly */
    1948         252 :       if (!FpX_is_squarefree(T,p)) pari_err_IRREDPOL("ffgen",T);
    1949         245 :       break;
    1950             :     case t_INT:
    1951        1720 :       d = ispseudoprimepower(T,&p);
    1952        1718 :       if (!d) pari_err_PRIME("ffgen",T);
    1953        1718 :       T = init_Fq(p, d, v);
    1954        1723 :       break;
    1955             :     case t_VEC: case t_COL:
    1956         413 :       if (lg(T) == 3) {
    1957         413 :         p = gel(T,1);
    1958         413 :         A = gel(T,2);
    1959         413 :         if (typ(p) == t_INT && typ(A) == t_INT)
    1960             :         {
    1961         413 :           d = itos(A);
    1962         413 :           T = init_Fq(p, d, v);
    1963         413 :           break;
    1964             :         }
    1965             :       }
    1966             :     default:
    1967           0 :       pari_err_TYPE("ffgen",T);
    1968           0 :       return NULL;
    1969             :   }
    1970        2409 :   if (v < 0) v = varn(T);
    1971        2409 :   if (lgefint(p)==3)
    1972             :   {
    1973        2121 :     ulong pp = p[2];
    1974        2121 :     long sv = evalvarn(v);
    1975        2121 :     if (pp==2)
    1976             :     {
    1977        1002 :       ff[1] = t_FF_F2xq;
    1978        1002 :       T = ZX_to_F2x(T); T[1] = sv;
    1979        1002 :       A = polx_F2x(sv); if (d == 1) A = F2x_rem(A, T);
    1980        1002 :       p = gen_2;
    1981             :     }
    1982             :     else
    1983             :     {
    1984        1119 :       ff[1] = t_FF_Flxq;
    1985        1119 :       T = ZX_to_Flx(T,pp); T[1] = sv;
    1986        1119 :       A = polx_Flx(sv); if (d == 1) A = Flx_rem(A, T, pp);
    1987        1119 :       p = icopy(p);
    1988             :     }
    1989             :   }
    1990             :   else
    1991             :   {
    1992         288 :     ff[1] = t_FF_FpXQ;
    1993         288 :     setvarn(T,v);
    1994         288 :     A = pol_x(v); if (d == 1) A = FpX_rem(A, T, p);
    1995         288 :     p = icopy(p);
    1996             :   }
    1997        2409 :   gel(ff,2) = A;
    1998        2409 :   gel(ff,3) = T;
    1999        2409 :   gel(ff,4) = p; return ff;
    2000             : }
    2001             : 
    2002             : GEN
    2003        2989 : p_to_FF(GEN p, long v)
    2004             : {
    2005             :   GEN A, T;
    2006        2989 :   GEN ff = cgetg(5,t_FFELT);
    2007        2989 :   if (lgefint(p)==3)
    2008             :   {
    2009        2984 :     ulong pp = p[2];
    2010        2984 :     long sv = evalvarn(v);
    2011        2984 :     if (pp==2)
    2012             :     {
    2013         329 :       ff[1] = t_FF_F2xq;
    2014         329 :       T = polx_F2x(sv);
    2015         329 :       A = pol1_F2x(sv);
    2016         329 :       p = gen_2;
    2017             :     }
    2018             :     else
    2019             :     {
    2020        2655 :       ff[1] = t_FF_Flxq;
    2021        2655 :       T = polx_Flx(sv);
    2022        2655 :       A = pol1_Flx(sv);
    2023        2655 :       p = icopy(p);
    2024             :     }
    2025             :   }
    2026             :   else
    2027             :   {
    2028           5 :     ff[1] = t_FF_FpXQ;
    2029           5 :     T = pol_x(v);
    2030           5 :     A = pol_1(v);
    2031           5 :     p = icopy(p);
    2032             :   }
    2033        2989 :   gel(ff,2) = A;
    2034        2989 :   gel(ff,3) = T;
    2035        2989 :   gel(ff,4) = p; return ff;
    2036             : }
    2037             : GEN
    2038        4179 : Tp_to_FF(GEN T, GEN p)
    2039             : {
    2040             :   GEN A, ff;
    2041             :   long v;
    2042        4179 :   if (!T) return p_to_FF(p,0);
    2043        3990 :   ff = cgetg(5,t_FFELT);
    2044        3990 :   v = varn(T);
    2045        3990 :   if (lgefint(p)==3)
    2046             :   {
    2047        3780 :     ulong pp = p[2];
    2048        3780 :     long sv = evalvarn(v);
    2049        3780 :     if (pp==2)
    2050             :     {
    2051        1526 :       ff[1] = t_FF_F2xq;
    2052        1526 :       T = ZX_to_F2x(T);
    2053        1526 :       A = pol1_F2x(sv);
    2054        1526 :       p = gen_2;
    2055             :     }
    2056             :     else
    2057             :     {
    2058        2254 :       ff[1] = t_FF_Flxq;
    2059        2254 :       T = ZX_to_Flx(T, pp);
    2060        2254 :       A = pol1_Flx(sv);
    2061        2254 :       p = icopy(p);
    2062             :     }
    2063             :   }
    2064             :   else
    2065             :   {
    2066         210 :     ff[1] = t_FF_FpXQ;
    2067         210 :     T = ZX_copy(T);
    2068         210 :     A = pol_1(v);
    2069         210 :     p = icopy(p);
    2070             :   }
    2071        3990 :   gel(ff,2) = A;
    2072        3990 :   gel(ff,3) = T;
    2073        3990 :   gel(ff,4) = p; return ff;
    2074             : }
    2075             : 
    2076             : GEN
    2077          28 : fforder(GEN x, GEN o)
    2078             : {
    2079          28 :   if (typ(x)!=t_FFELT) pari_err_TYPE("fforder",x);
    2080          28 :   return FF_order(x,o);
    2081             : }
    2082             : 
    2083             : GEN
    2084         392 : ffprimroot(GEN x, GEN *o)
    2085             : {
    2086         392 :   if (typ(x)!=t_FFELT) pari_err_TYPE("ffprimroot",x);
    2087         392 :   return FF_primroot(x,o);
    2088             : }
    2089             : 
    2090             : GEN
    2091         182 : fflog(GEN x, GEN g, GEN o)
    2092             : {
    2093         182 :   if (typ(x)!=t_FFELT) pari_err_TYPE("fflog",x);
    2094         182 :   if (typ(g)!=t_FFELT) pari_err_TYPE("fflog",g);
    2095         182 :   return FF_log(x,g,o);
    2096             : }
    2097             : 
    2098             : GEN
    2099      144438 : ffrandom(GEN ff)
    2100             : {
    2101             :   ulong pp;
    2102      144438 :   GEN r, T, p, z = _initFF(ff,&T,&p,&pp);
    2103      144438 :   switch(ff[1])
    2104             :   {
    2105             :   case t_FF_FpXQ:
    2106        6940 :     r = random_FpX(degpol(T), varn(T), p);
    2107        6940 :     break;
    2108             :   case t_FF_F2xq:
    2109      117271 :     r = random_F2x(F2x_degree(T), T[1]);
    2110      117271 :     break;
    2111             :   default:
    2112       20227 :     r = random_Flx(degpol(T), T[1], pp);
    2113             :   }
    2114      144438 :   return _mkFF(ff,z,r);
    2115             : }
    2116             : 
    2117             : int
    2118    11489033 : Rg_is_FF(GEN c, GEN *ff)
    2119             : {
    2120    11489033 :   switch(typ(c))
    2121             :   {
    2122             :   case t_FFELT:
    2123        8008 :     if (!*ff) *ff = c;
    2124        7840 :     else if (!FF_samefield(*ff, c)) return 0;
    2125        8008 :     break;
    2126             :   default:
    2127    11481025 :     return 0;
    2128             :   }
    2129        8008 :   return 1;
    2130             : }
    2131             : 
    2132             : int
    2133    11481613 : RgC_is_FFC(GEN x, GEN *ff)
    2134             : {
    2135    11481613 :   long i, lx = lg(x);
    2136    11489621 :   for (i=lx-1; i>0; i--)
    2137    11489033 :     if (!Rg_is_FF(gel(x,i), ff)) return 0;
    2138         588 :   return (*ff != NULL);
    2139             : }
    2140             : 
    2141             : int
    2142    11481102 : RgM_is_FFM(GEN x, GEN *ff)
    2143             : {
    2144    11481102 :   long j, lx = lg(x);
    2145    11481613 :   for (j=lx-1; j>0; j--)
    2146    11481536 :     if (!RgC_is_FFC(gel(x,j), ff)) return 0;
    2147          77 :   return (*ff != NULL);
    2148             : }
    2149             : 
    2150             : static GEN
    2151        2206 : FqC_to_FpXQC(GEN x, GEN T, GEN p)
    2152             : {
    2153             :   long i, lx;
    2154        2206 :   GEN y = cgetg_copy(x,&lx);
    2155       34066 :   for(i=1; i<lx; i++)
    2156       31860 :     gel(y, i) = Fq_to_FpXQ(gel(x, i), T, p);
    2157        2206 :   return y;
    2158             : }
    2159             : 
    2160             : static GEN
    2161         390 : FqM_to_FpXQM(GEN x, GEN T, GEN p)
    2162             : {
    2163             :   long i, lx;
    2164         390 :   GEN y = cgetg_copy(x,&lx);
    2165        2526 :   for(i=1; i<lx; i++)
    2166        2136 :     gel(y, i) = FqC_to_FpXQC(gel(x, i), T, p);
    2167         390 :   return y;
    2168             : }
    2169             : 
    2170             : /* for functions t_MAT -> t_MAT */
    2171             : static GEN
    2172         308 : FFM_wrap(GEN M, GEN ff, GEN (*Fq)(GEN,GEN,GEN),
    2173             :                        GEN (*Flxq)(GEN,GEN,ulong),
    2174             :                        GEN (*F2xq)(GEN,GEN))
    2175             : {
    2176         308 :   pari_sp av = avma;
    2177             :   ulong pp;
    2178             :   GEN T, p;
    2179         308 :   _getFF(ff,&T,&p,&pp); M = FFM_to_raw(M, ff);
    2180         308 :   switch(ff[1])
    2181             :   {
    2182         119 :   case t_FF_FpXQ: M = Fq(M,T,p); if (M) M = FqM_to_FpXQM(M,T,p);
    2183         119 :                   break;
    2184          70 :   case t_FF_F2xq: M = F2xq(M,T); break;
    2185         119 :   default: M = Flxq(M,T,pp); break;
    2186             :   }
    2187         308 :   if (!M) return gc_NULL(av);
    2188         273 :   return gerepilecopy(av, raw_to_FFM(M, ff));
    2189             : }
    2190             : 
    2191             : /* for functions (t_MAT, t_MAT) -> t_MAT */
    2192             : static GEN
    2193        4116 : FFM_FFM_wrap(GEN M, GEN N, GEN ff,
    2194             :              GEN (*Fq)(GEN, GEN, GEN, GEN),
    2195             :              GEN (*Flxq)(GEN, GEN, GEN, ulong),
    2196             :              GEN (*F2xq)(GEN, GEN, GEN))
    2197             : {
    2198        4116 :   pari_sp av = avma;
    2199             :   ulong pp;
    2200             :   GEN T, p;
    2201        4116 :   int is_sqr = M==N;
    2202        4116 :   _getFF(ff, &T, &p, &pp);
    2203        4116 :   M = FFM_to_raw(M, ff);
    2204        4116 :   N = is_sqr? M: FFM_to_raw(N, ff);
    2205        4116 :   switch(ff[1])
    2206             :   {
    2207         320 :   case t_FF_FpXQ: M = Fq(M, N, T, p); if (M) M = FqM_to_FpXQM(M, T, p);
    2208         320 :                   break;
    2209        1449 :   case t_FF_F2xq: M = F2xq(M, N, T); break;
    2210        2347 :   default: M = Flxq(M, N, T, pp); break;
    2211             :   }
    2212        4116 :   if (!M) return gc_NULL(av);
    2213        4025 :   return gerepilecopy(av, raw_to_FFM(M, ff));
    2214             : }
    2215             : 
    2216             : /* for functions (t_MAT, t_COL) -> t_COL */
    2217             : static GEN
    2218         196 : FFM_FFC_wrap(GEN M, GEN C, GEN ff,
    2219             :              GEN (*Fq)(GEN, GEN, GEN, GEN),
    2220             :              GEN (*Flxq)(GEN, GEN, GEN, ulong),
    2221             :              GEN (*F2xq)(GEN, GEN, GEN))
    2222             : {
    2223         196 :   pari_sp av = avma;
    2224             :   ulong pp;
    2225             :   GEN T, p;
    2226         196 :   _getFF(ff, &T, &p, &pp);
    2227         196 :   M = FFM_to_raw(M, ff);
    2228         196 :   C = FFC_to_raw(C, ff);
    2229         196 :   switch(ff[1])
    2230             :   {
    2231          63 :   case t_FF_FpXQ: C = Fq(M, C, T, p); if (C) C = FqC_to_FpXQC(C, T, p);
    2232          63 :                   break;
    2233          70 :   case t_FF_F2xq: C = F2xq(M, C, T); break;
    2234          63 :   default: C = Flxq(M, C, T, pp); break;
    2235             :   }
    2236         196 :   if (!C) return gc_NULL(av);
    2237         147 :   return gerepilecopy(av, raw_to_FFC(C, ff));
    2238             : }
    2239             : 
    2240             : GEN
    2241          77 : FFM_ker(GEN M, GEN ff)
    2242          77 : { return FFM_wrap(M,ff, &FqM_ker,&FlxqM_ker,&F2xqM_ker); }
    2243             : GEN
    2244          49 : FFM_image(GEN M, GEN ff)
    2245          49 : { return FFM_wrap(M,ff, &FqM_image,&FlxqM_image,&F2xqM_image); }
    2246             : GEN
    2247         147 : FFM_inv(GEN M, GEN ff)
    2248         147 : { return FFM_wrap(M,ff, &FqM_inv,&FlxqM_inv,&F2xqM_inv); }
    2249             : GEN
    2250          35 : FFM_suppl(GEN M, GEN ff)
    2251          35 : { return FFM_wrap(M,ff, &FqM_suppl,&FlxqM_suppl,&F2xqM_suppl); }
    2252             : 
    2253             : GEN
    2254          84 : FFM_deplin(GEN M, GEN ff)
    2255             : {
    2256          84 :   pari_sp av = avma;
    2257             :   ulong pp;
    2258             :   GEN C, T, p;
    2259          84 :   _getFF(ff, &T, &p, &pp); M = FFM_to_raw(M, ff);
    2260          84 :   switch(ff[1]) {
    2261          35 :   case t_FF_FpXQ: C = FqM_deplin(M, T, p);
    2262          35 :     if (C) C = FqC_to_FpXQC(C, T, p); break;
    2263          14 :   case t_FF_F2xq: C = F2xqM_deplin(M, T); break;
    2264          35 :   default: C = FlxqM_deplin(M, T, pp); break;
    2265             :   }
    2266          84 :   if (!C) return gc_NULL(av);
    2267          49 :   return gerepilecopy(av, raw_to_FFC(C, ff));
    2268             : }
    2269             : 
    2270             : GEN
    2271          21 : FFM_indexrank(GEN M, GEN ff)
    2272             : {
    2273          21 :   pari_sp av = avma;
    2274             :   ulong pp;
    2275             :   GEN R, T, p;
    2276          21 :   _getFF(ff,&T,&p,&pp); M = FFM_to_raw(M, ff);
    2277          21 :   switch(ff[1]) {
    2278           7 :   case t_FF_FpXQ: R = FqM_indexrank(M,T,p); break;
    2279           7 :   case t_FF_F2xq: R = F2xqM_indexrank(M,T); break;
    2280           7 :   default: R = FlxqM_indexrank(M,T,pp); break;
    2281             :   }
    2282          21 :   return gerepileupto(av, R);
    2283             : }
    2284             : 
    2285             : long
    2286          63 : FFM_rank(GEN M, GEN ff)
    2287             : {
    2288          63 :   pari_sp av = avma;
    2289             :   long r;
    2290             :   ulong pp;
    2291             :   GEN T, p;
    2292          63 :   _getFF(ff,&T,&p,&pp); M = FFM_to_raw(M, ff);
    2293          63 :   switch(ff[1])
    2294             :   {
    2295          28 :   case t_FF_FpXQ: r = FqM_rank(M,T,p); break;
    2296           7 :   case t_FF_F2xq: r = F2xqM_rank(M,T); break;
    2297          28 :   default: r = FlxqM_rank(M,T,pp); break;
    2298             :   }
    2299          63 :   return gc_long(av,r);
    2300             : }
    2301             : GEN
    2302          63 : FFM_det(GEN M, GEN ff)
    2303             : {
    2304          63 :   pari_sp av = avma;
    2305             :   ulong pp;
    2306             :   GEN d, T, p;
    2307          63 :   _getFF(ff,&T,&p,&pp); M = FFM_to_raw(M, ff);
    2308          63 :   switch(ff[1])
    2309             :   {
    2310          28 :   case t_FF_FpXQ: d = FqM_det(M,T,p); break;
    2311           7 :   case t_FF_F2xq: d = F2xqM_det(M,T); break;
    2312          28 :   default: d = FlxqM_det(M,T,pp); break;
    2313             :   }
    2314          63 :   return gerepilecopy(av, mkFF_i(ff, d));
    2315             : }
    2316             : 
    2317             : GEN
    2318          56 : FFM_FFC_gauss(GEN M, GEN C, GEN ff)
    2319             : {
    2320          56 :   return FFM_FFC_wrap(M, C, ff, FqM_FqC_gauss,
    2321             :                       FlxqM_FlxqC_gauss, F2xqM_F2xqC_gauss);
    2322             : }
    2323             : 
    2324             : GEN
    2325          63 : FFM_gauss(GEN M, GEN N, GEN ff)
    2326             : {
    2327          63 :   return FFM_FFM_wrap(M, N, ff, FqM_gauss,
    2328             :                       FlxqM_gauss, F2xqM_gauss);
    2329             : }
    2330             : 
    2331             : GEN
    2332          63 : FFM_FFC_invimage(GEN M, GEN C, GEN ff)
    2333             : {
    2334          63 :   return FFM_FFC_wrap(M, C, ff, FqM_FqC_invimage,
    2335             :                       FlxqM_FlxqC_invimage, F2xqM_F2xqC_invimage);
    2336             : }
    2337             : 
    2338             : GEN
    2339         105 : FFM_invimage(GEN M, GEN N, GEN ff)
    2340             : {
    2341         105 :   return FFM_FFM_wrap(M, N, ff, FqM_invimage,
    2342             :                       FlxqM_invimage, F2xqM_invimage);
    2343             : }
    2344             : 
    2345             : GEN
    2346          77 : FFM_FFC_mul(GEN M, GEN C, GEN ff)
    2347             : {
    2348          77 :   return FFM_FFC_wrap(M, C, ff, FqM_FqC_mul,
    2349             :                       FlxqM_FlxqC_mul, F2xqM_F2xqC_mul);
    2350             : }
    2351             : 
    2352             : GEN
    2353        3948 : FFM_mul(GEN M, GEN N, GEN ff)
    2354             : {
    2355        3948 :   return FFM_FFM_wrap(M, N, ff, FqM_mul, FlxqM_mul, F2xqM_mul);
    2356             : }

Generated by: LCOV version 1.13