Code coverage tests

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

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

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

LCOV - code coverage report
Current view: top level - basemath - base1.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.8.0 lcov report (development 16962-5a32637) Lines: 1526 1625 93.9 %
Date: 2014-10-29 Functions: 110 125 88.0 %
Legend: Lines: hit not hit | Branches: + taken - not taken # not executed Branches: 964 1206 79.9 %

           Branch data     Line data    Source code
       1                 :            : /* Copyright (C) 2000  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                 :            : /**************************************************************/
      15                 :            : /*                                                            */
      16                 :            : /*                        NUMBER FIELDS                       */
      17                 :            : /*                                                            */
      18                 :            : /**************************************************************/
      19                 :            : #include "pari.h"
      20                 :            : #include "paripriv.h"
      21                 :            : 
      22                 :            : int new_galois_format = 0;
      23                 :            : 
      24                 :            : void
      25                 :       7196 : checkrnf(GEN rnf)
      26                 :            : {
      27 [ +  - ][ -  + ]:       7196 :   if (typ(rnf)!=t_VEC || lg(rnf)!=13) pari_err_TYPE("checkrnf",rnf);
      28                 :       7196 : }
      29                 :            : 
      30                 :            : GEN
      31                 :     577781 : checkbnf_i(GEN X)
      32                 :            : {
      33         [ +  + ]:     577781 :   if (typ(X) == t_VEC)
      34      [ +  +  - ]:     577669 :     switch (lg(X))
      35                 :            :     {
      36                 :     576766 :       case 11: return X;
      37                 :        903 :       case 7:  return checkbnf_i(bnr_get_bnf(X));
      38                 :            :     }
      39                 :     577781 :   return NULL;
      40                 :            : }
      41                 :            : 
      42                 :            : GEN
      43                 :    4568635 : checknf_i(GEN X)
      44                 :            : {
      45         [ +  + ]:    4568635 :   if (typ(X)==t_VEC)
      46   [ +  +  +  +  :    4568523 :     switch(lg(X))
                      + ]
      47                 :            :     {
      48                 :    4486742 :       case 10: return X;
      49                 :      76468 :       case 11: return checknf_i(bnf_get_nf(X));
      50                 :       5299 :       case 7:  return checknf_i(bnr_get_bnf(X));
      51         [ +  - ]:          7 :       case 3: if (typ(gel(X,2)) == t_POLMOD) return checknf_i(gel(X,1));
      52                 :            :     }
      53                 :    4568635 :   return NULL;
      54                 :            : }
      55                 :            : 
      56                 :            : GEN
      57                 :     576731 : checkbnf(GEN x)
      58                 :            : {
      59                 :     576731 :   GEN bnf = checkbnf_i(x);
      60         [ -  + ]:     576731 :   if (!bnf) pari_err_TYPE("checkbnf [please apply bnfinit()]",x);
      61                 :     576731 :   return bnf;
      62                 :            : }
      63                 :            : 
      64                 :            : GEN
      65                 :    4486714 : checknf(GEN x)
      66                 :            : {
      67                 :    4486714 :   GEN nf = checknf_i(x);
      68         [ +  + ]:    4486714 :   if (!nf) pari_err_TYPE("checknf [please apply nfinit()]",x);
      69                 :    4486707 :   return nf;
      70                 :            : }
      71                 :            : 
      72                 :            : void
      73                 :     168292 : checkbnr(GEN bnr)
      74                 :            : {
      75 [ +  - ][ -  + ]:     168292 :   if (typ(bnr)!=t_VEC || lg(bnr)!=7)
      76                 :          0 :     pari_err_TYPE("checkbnr [please apply bnrinit()]",bnr);
      77                 :     168292 :   (void)checkbnf(bnr_get_bnf(bnr));
      78                 :     168292 : }
      79                 :            : 
      80                 :            : void
      81                 :          0 : checkbnrgen(GEN bnr)
      82                 :            : {
      83                 :          0 :   checkbnr(bnr);
      84         [ #  # ]:          0 :   if (lg(bnr_get_clgp(bnr))<=3)
      85                 :          0 :     pari_err_TYPE("checkbnrgen [apply bnrinit(,,1), not bnrinit()]",bnr);
      86                 :          0 : }
      87                 :            : 
      88                 :            : void
      89                 :          0 : checksqmat(GEN x, long N)
      90                 :            : {
      91         [ #  # ]:          0 :   if (typ(x)!=t_MAT) pari_err_TYPE("checksqmat",x);
      92 [ #  # ][ #  # ]:          0 :   if (lg(x) == 1 || lgcols(x) != N+1) pari_err_DIM("checksqmat");
      93                 :          0 : }
      94                 :            : 
      95                 :            : GEN
      96                 :     178547 : checkbid_i(GEN bid)
      97                 :            : {
      98                 :            :   GEN f;
      99 [ +  + ][ +  + ]:     178547 :   if (typ(bid)!=t_VEC || lg(bid)!=6) return NULL;
     100                 :     173717 :   f = bid_get_mod(bid);
     101 [ +  + ][ -  + ]:     173717 :   if (typ(f)!=t_VEC || lg(f)!=3) return NULL;
     102                 :     178547 :   return bid;
     103                 :            : }
     104                 :            : void
     105                 :     173696 : checkbid(GEN bid)
     106                 :            : {
     107         [ -  + ]:     173696 :   if (!checkbid_i(bid)) pari_err_TYPE("checkbid",bid);
     108                 :     173696 : }
     109                 :            : void
     110                 :      11067 : checkabgrp(GEN v)
     111                 :            : {
     112         [ +  - ]:      11067 :   if (typ(v) == t_VEC) switch(lg(v))
              [ +  +  - ]
     113                 :            :   {
     114         [ -  + ]:      10983 :     case 4: if (typ(gel(v,3)) != t_VEC) break;
     115         [ +  + ]:      11067 :     case 3: if (typ(gel(v,2)) != t_VEC) break;
     116         [ -  + ]:      11039 :             if (typ(gel(v,1)) != t_INT) break;
     117                 :      11039 :             return;/*OK*/
     118                 :          0 :     default: break;
     119                 :            :   }
     120                 :         28 :   pari_err_TYPE("checkabgrp",v);
     121                 :            : }
     122                 :            : 
     123                 :            : GEN
     124                 :      49497 : checknfelt_mod(GEN nf, GEN x, const char *s)
     125                 :            : {
     126                 :      49497 :   GEN T = gel(x,1), a = gel(x,2), Tnf = nf_get_pol(nf);
     127         [ +  + ]:      49497 :   if (!RgX_equal_var(T, Tnf)) pari_err_MODULUS(s, T, Tnf);
     128                 :      49427 :   return a;
     129                 :            : }
     130                 :            : 
     131                 :            : void
     132                 :       2086 : check_ZKmodule(GEN x, const char *s)
     133                 :            : {
     134 [ +  - ][ -  + ]:       2086 :   if (typ(x) != t_VEC || lg(x) < 3) pari_err_TYPE(s,x);
     135         [ -  + ]:       2086 :   if (typ(gel(x,1)) != t_MAT) pari_err_TYPE(s,gel(x,1));
     136         [ -  + ]:       2086 :   if (typ(gel(x,2)) != t_VEC) pari_err_TYPE(s,gel(x,2));
     137         [ -  + ]:       2086 :   if (lg(gel(x,2)) != lgcols(x)) pari_err_DIM(s);
     138                 :       2086 : }
     139                 :            : 
     140                 :            : GEN
     141                 :      12453 : get_bnf(GEN x, long *t)
     142                 :            : {
     143   [ +  +  +  +  :      12453 :   switch(typ(x))
                      + ]
     144                 :            :   {
     145                 :         56 :     case t_POL: *t = typ_POL;  return NULL;
     146                 :         56 :     case t_QUAD: *t = typ_Q  ; return NULL;
     147                 :            :     case t_VEC:
     148   [ +  +  +  +  :      11893 :       switch(lg(x))
             +  +  +  +  
                      + ]
     149                 :            :       {
     150                 :       4375 :         case 5: *t = typ_QUA; return NULL;
     151                 :            :         case 6:
     152 [ +  + ][ +  - ]:        245 :           if (typ(gel(x,1)) == t_VEC && typ(gel(x,3)) == t_MAT)
     153                 :        112 :           { *t = typ_BID; return NULL; }
     154 [ +  + ][ +  - ]:        133 :           if (typ(gel(x,2)) == t_COL && typ(gel(x,3)) == t_INT)
     155                 :         77 :           { *t = typ_PRID; return NULL; }
     156                 :         56 :           break;
     157                 :         63 :         case 7:  *t = typ_BNR;
     158 [ +  - ][ +  - ]:         63 :           x = bnr_get_bnf(x); if (typ(x)!=t_VEC || lg(x)!=11) break;
     159                 :         63 :           return x;
     160                 :            :         case 9:
     161                 :         56 :           x = gel(x,2);
     162 [ +  - ][ +  - ]:         56 :           if (typ(x) == t_VEC && lg(x) == 4) *t = typ_GAL;
     163                 :         56 :           return NULL;
     164                 :         56 :         case 10: *t = typ_NF; return NULL;
     165                 :        140 :         case 11: *t = typ_BNF; return x;
     166                 :         56 :         case 13: *t = typ_RNF; return NULL;
     167                 :        266 :         case 17: *t = typ_ELL; return NULL;
     168                 :            :       }
     169                 :       6692 :       break;
     170                 :            :     case t_COL:
     171         [ +  + ]:        112 :       if (get_prid(x)) { *t = typ_MODPR; return NULL; }
     172                 :         56 :       break;
     173                 :            :   }
     174                 :      12453 :   *t = typ_NULL; return NULL;
     175                 :            : }
     176                 :            : 
     177                 :            : GEN
     178                 :      55489 : get_nf(GEN x, long *t)
     179                 :            : {
     180   [ +  +  +  +  :      55489 :   switch(typ(x))
                      + ]
     181                 :            :   {
     182                 :        133 :     case t_POL : *t = typ_POL; return NULL;
     183                 :        133 :     case t_QUAD: *t = typ_Q  ; return NULL;
     184                 :            :     case t_VEC:
     185   [ +  +  +  +  :      54096 :       switch(lg(x))
          +  +  +  +  +  
                      + ]
     186                 :            :       {
     187                 :            :         case 3:
     188         [ -  + ]:        133 :           if (typ(gel(x,2)) != t_POLMOD) break;
     189                 :        133 :           return get_nf(gel(x,1),t);
     190                 :        133 :         case 5: *t = typ_QUA; return NULL;
     191                 :            :         case 6:
     192 [ +  + ][ +  - ]:      49595 :           if (typ(gel(x,1)) == t_VEC && typ(gel(x,3)) == t_MAT)
     193                 :        266 :           { *t = typ_BID; return NULL; }
     194 [ +  + ][ +  - ]:      49329 :           if (typ(gel(x,2)) == t_COL && typ(gel(x,3)) == t_INT)
     195                 :      49196 :           { *t = typ_PRID; return NULL; }
     196                 :        133 :           break;
     197                 :        140 :         case 7: *t = typ_BNR;
     198 [ +  - ][ +  - ]:        140 :           x = bnr_get_bnf(x); if (typ(x)!=t_VEC || lg(x)!=11) break;
     199 [ +  - ][ +  - ]:        140 :           x = bnf_get_nf(x);  if (typ(x)!=t_VEC || lg(x)!=10) break;
     200                 :        140 :           return x;
     201                 :            :         case 9:
     202                 :        238 :           x = gel(x,2);
     203 [ +  - ][ +  - ]:        238 :           if (typ(x) == t_VEC && lg(x) == 4) *t = typ_GAL;
     204                 :        238 :           return NULL;
     205                 :        427 :         case 10: *t = typ_NF; return x;
     206                 :        154 :         case 11: *t = typ_BNF;
     207 [ +  - ][ +  - ]:        154 :           x = bnf_get_nf(x); if (typ(x)!=t_VEC || lg(x)!=10) break;
     208                 :        154 :           return x;
     209                 :        210 :         case 13: *t = typ_RNF; return NULL;
     210                 :       2933 :         case 17: *t = typ_ELL; return NULL;
     211                 :            :       }
     212                 :        266 :       break;
     213                 :            :     case t_COL:
     214         [ +  + ]:        266 :       if (get_prid(x)) { *t = typ_MODPR; return NULL; }
     215                 :        133 :       break;
     216                 :            :   }
     217                 :      55489 :   *t = typ_NULL; return NULL;
     218                 :            : }
     219                 :            : 
     220                 :            : long
     221                 :      35861 : nftyp(GEN x)
     222                 :            : {
     223   [ +  +  +  - ]:      35861 :   switch(typ(x))
     224                 :            :   {
     225                 :          7 :     case t_POL : return typ_POL;
     226                 :          7 :     case t_QUAD: return typ_Q;
     227                 :            :     case t_VEC:
     228   [ +  +  +  +  :      35847 :       switch(lg(x))
                +  +  - ]
     229                 :            :       {
     230                 :            :         case 10:
     231         [ +  + ]:      35805 :           if (typ(gel(x,1))!=t_POL) break;
     232                 :      35798 :           return typ_NF;
     233                 :            :         case 11:
     234 [ +  - ][ +  - ]:          7 :           x = bnf_get_nf(x); if (typ(x)!=t_VEC || lg(x)!=10) break;
     235                 :          7 :           return typ_BNF;
     236                 :            :         case 7:
     237 [ +  + ][ +  - ]:         14 :           x = bnr_get_bnf(x); if (typ(x)!=t_VEC || lg(x)!=11) break;
     238 [ +  - ][ +  - ]:          7 :           x = bnf_get_nf(x);  if (typ(x)!=t_VEC || lg(x)!=10) break;
     239                 :          7 :           return typ_BNR;
     240                 :            :         case 6:
     241 [ -  + ][ #  # ]:          7 :           if (typ(gel(x,1)) != t_VEC || typ(gel(x,3)) != t_MAT) break;
     242                 :          0 :           return typ_BID;
     243                 :            :         case 9:
     244                 :          7 :           x = gel(x,2);
     245 [ -  + ][ #  # ]:          7 :           if (typ(x) == t_VEC && lg(x) == 4) return typ_GAL;
     246                 :         14 :         case 17: return typ_ELL;
     247                 :            :       }
     248                 :            :   }
     249                 :      35861 :   return typ_NULL;
     250                 :            : }
     251                 :            : 
     252                 :            : /*************************************************************************/
     253                 :            : /**                                                                     **/
     254                 :            : /**                           GALOIS GROUP                              **/
     255                 :            : /**                                                                     **/
     256                 :            : /*************************************************************************/
     257                 :            : 
     258                 :            : /* exchange elements i and j in vector x */
     259                 :            : static GEN
     260                 :       5698 : transroot(GEN x, int i, int j)
     261                 :       5698 : { x = leafcopy(x); swap(gel(x,i), gel(x,j)); return x; }
     262                 :            : 
     263                 :            : GEN
     264                 :       3136 : tschirnhaus(GEN x)
     265                 :            : {
     266                 :       3136 :   pari_sp av = avma, av2;
     267                 :       3136 :   long a, v = varn(x);
     268                 :       3136 :   GEN u, y = cgetg(5,t_POL);
     269                 :            : 
     270         [ -  + ]:       3136 :   if (typ(x)!=t_POL) pari_err_TYPE("tschirnhaus",x);
     271         [ -  + ]:       3136 :   if (lg(x) < 4) pari_err_CONSTPOL("tschirnhaus");
     272         [ -  + ]:       3136 :   if (v) { u = leafcopy(x); setvarn(u,0); x=u; }
     273                 :       3136 :   y[1] = evalsigne(1)|evalvarn(0);
     274                 :            :   do
     275                 :            :   {
     276         [ +  + ]:       3248 :     a = random_bits(2); if (a==0) a  = 1; gel(y,4) = stoi(a);
     277         [ +  + ]:       3248 :     a = random_bits(3); if (a>=4) a -= 8; gel(y,3) = stoi(a);
     278         [ +  + ]:       3248 :     a = random_bits(3); if (a>=4) a -= 8; gel(y,2) = stoi(a);
     279                 :       3248 :     u = RgXQ_charpoly(y,x,v); av2 = avma;
     280                 :            :   }
     281         [ +  + ]:       3248 :   while (degpol(RgX_gcd(u,RgX_deriv(u)))); /* while u not separable */
     282         [ -  + ]:       3136 :   if (DEBUGLEVEL>1)
     283                 :          0 :     err_printf("Tschirnhaus transform. New pol: %Ps",u);
     284                 :       3136 :   avma=av2; return gerepileupto(av,u);
     285                 :            : }
     286                 :            : 
     287                 :            : /* Assume pol in Z[X], monic of degree n. Find L in Z such that
     288                 :            :  * POL = L^(-n) pol(L x) is monic in Z[X]. Return POL and set *ptk = L.
     289                 :            :  * No GC. */
     290                 :            : GEN
     291                 :       9625 : ZX_Z_normalize(GEN pol, GEN *ptk)
     292                 :            : {
     293                 :       9625 :   long i,j, sk, n = degpol(pol); /* > 0 */
     294                 :            :   GEN k, fa, P, E, a, POL;
     295                 :            : 
     296                 :       9625 :   a = pol + 2; k = gel(a,n-1); /* a[i] = coeff of degree i */
     297         [ +  + ]:      40796 :   for (i = n-2; i >= 0; i--)
     298                 :            :   {
     299                 :      37268 :     k = gcdii(k, gel(a,i));
     300 [ +  + ][ +  + ]:      37268 :     if (is_pm1(k)) { if (ptk) *ptk = gen_1; return pol; }
     301                 :            :   }
     302                 :       3528 :   sk = signe(k);
     303 [ +  + ][ +  - ]:       3528 :   if (!sk) { if (ptk) *ptk = gen_1; return pol; /* monomial! */ }
     304                 :       3388 :   fa = absi_factor_limit(k, 0); k = gen_1;
     305                 :       3388 :   P = gel(fa,1);
     306                 :       3388 :   E = gel(fa,2);
     307                 :       3388 :   POL = leafcopy(pol); a = POL+2;
     308         [ +  + ]:       7868 :   for (i = lg(P)-1; i > 0; i--)
     309                 :            :   {
     310                 :       4480 :     GEN p = gel(P,i), pv, pvj;
     311                 :       4480 :     long vmin = itos(gel(E,i));
     312                 :            :     /* find v_p(k) = min floor( v_p(a[i]) / (n-i)) */
     313         [ +  + ]:      32991 :     for (j=n-1; j>=0; j--)
     314                 :            :     {
     315                 :            :       long v;
     316         [ +  + ]:      28511 :       if (!signe(gel(a,j))) continue;
     317                 :      18347 :       v = Z_pval(gel(a,j), p) / (n - j);
     318         [ +  + ]:      18347 :       if (v < vmin) vmin = v;
     319                 :            :     }
     320         [ +  + ]:       4480 :     if (!vmin) continue;
     321                 :        952 :     pvj = pv = powiu(p,vmin); k = mulii(k, pv);
     322                 :            :     /* a[j] /= p^(v*(n-j)) */
     323         [ +  + ]:       7133 :     for (j=n-1; j>=0; j--)
     324                 :            :     {
     325         [ +  + ]:       6181 :       if (j < n-1) pvj = mulii(pvj, pv);
     326                 :       6181 :       gel(a,j) = diviiexact(gel(a,j), pvj);
     327                 :            :     }
     328                 :            :   }
     329         [ +  + ]:       9625 :   if (ptk) *ptk = k; return POL;
     330                 :            : }
     331                 :            : 
     332                 :            : /* Assume pol != 0 in Z[X]. Find C in Q, L in Z such that POL = C pol(x/L) monic
     333                 :            :  * in Z[X]. Return POL and set *pL = L. Wasteful (but correct) if pol is not
     334                 :            :  * primitive: better if caller used Q_primpart already. No GC. */
     335                 :            : GEN
     336                 :       9639 : ZX_primitive_to_monic(GEN pol, GEN *pL)
     337                 :            : {
     338                 :       9639 :   long i,j, n = degpol(pol);
     339                 :       9639 :   GEN lc = leading_term(pol), L, fa, P, E, a, POL;
     340                 :            : 
     341         [ +  + ]:       9639 :   if (is_pm1(lc))
     342                 :            :   {
     343         [ +  - ]:       9457 :     if (pL) *pL = gen_1;
     344         [ +  + ]:       9457 :     return signe(lc) < 0? ZX_neg(pol): pol;
     345                 :            :   }
     346         [ +  + ]:        182 :   if (signe(lc) < 0)
     347                 :          7 :     POL = ZX_neg(pol);
     348                 :            :   else
     349                 :        175 :     POL = leafcopy(pol);
     350                 :        182 :   a = POL+2; lc = gel(a,n);
     351                 :        182 :   fa = Z_factor_limit(lc,0); L = gen_1;
     352                 :        182 :   P = gel(fa,1);
     353                 :        182 :   E = gel(fa,2);
     354         [ +  + ]:        448 :   for (i = lg(P)-1; i > 0; i--)
     355                 :            :   {
     356                 :        266 :     GEN p = gel(P,i), pk, pku;
     357                 :        266 :     long v, j0, e = itos(gel(E,i)), k = e/n, d = k*n - e;
     358                 :            : 
     359         [ +  + ]:        266 :     if (d < 0) { k++; d += n; }
     360                 :            :     /* k = ceil(e[i] / n); find d, k such that  p^d pol(x / p^k) monic */
     361         [ +  + ]:       1106 :     for (j=n-1; j>0; j--)
     362                 :            :     {
     363         [ +  + ]:        840 :       if (!signe(gel(a,j))) continue;
     364                 :        763 :       v = Z_pval(gel(a,j), p);
     365         [ +  + ]:        833 :       while (v + d < k * j) { k++; d += n; }
     366                 :            :     }
     367                 :        266 :     pk = powiu(p,k); j0 = d/k;
     368                 :        266 :     L = mulii(L, pk);
     369                 :            : 
     370                 :        266 :     pku = powiu(p,d - k*j0);
     371                 :            :     /* a[j] *= p^(d - kj) */
     372         [ +  + ]:       1113 :     for (j=j0; j>=0; j--)
     373                 :            :     {
     374         [ +  + ]:        847 :       if (j < j0) pku = mulii(pku, pk);
     375                 :        847 :       gel(a,j) = mulii(gel(a,j), pku);
     376                 :            :     }
     377                 :        266 :     j0++;
     378                 :        266 :     pku = powiu(p,k*j0 - d);
     379                 :            :     /* a[j] /= p^(kj - d) */
     380         [ +  + ]:        791 :     for (j=j0; j<=n; j++)
     381                 :            :     {
     382         [ +  + ]:        525 :       if (j > j0) pku = mulii(pku, pk);
     383                 :        525 :       gel(a,j) = diviiexact(gel(a,j), pku);
     384                 :            :     }
     385                 :            :   }
     386         [ +  - ]:        182 :   if (pL) *pL = L;
     387                 :       9639 :   return POL;
     388                 :            : }
     389                 :            : /* Assume pol != 0 in Z[X]. Find C,L in Q such that POL = C pol(x/L)
     390                 :            :  * monic in Z[X]. Return POL and set *pL = L.
     391                 :            :  * Wasteful (but correct) if pol is not primitive: better if caller used
     392                 :            :  * Q_primpart already. No GC. */
     393                 :            : GEN
     394                 :       9429 : ZX_Q_normalize(GEN pol, GEN *pL)
     395                 :            : {
     396                 :       9429 :   GEN lc, POL = ZX_primitive_to_monic(pol, &lc);
     397                 :       9429 :   POL = ZX_Z_normalize(POL, pL);
     398         [ +  + ]:       9429 :   if (pL) *pL = gdiv(lc, *pL);
     399                 :       9429 :   return POL;
     400                 :            : }
     401                 :            : /* pol != 0 in Z[x], returns a monic polynomial POL in Z[x] generating the
     402                 :            :  * same field: there exist C in Q, L in Z such that POL(x) = C pol(x/L).
     403                 :            :  * Set *L = NULL if L = 1, and to L otherwise. No garbage collecting. */
     404                 :            : GEN
     405                 :          0 : ZX_to_monic(GEN pol, GEN *L)
     406                 :            : {
     407                 :          0 :   long n = lg(pol)-1;
     408                 :          0 :   GEN lc = gel(pol,n);
     409 [ #  # ][ #  # ]:          0 :   if (is_pm1(lc)) { *L = gen_1; return signe(lc) > 0? pol: ZX_neg(pol); }
     410                 :          0 :   return ZX_primitive_to_monic(Q_primpart(pol), L);
     411                 :            : }
     412                 :            : 
     413                 :            : /* x1*x2^2 + x2*x3^2 + x3*x4^2 + x4*x1^2 */
     414                 :            : static GEN
     415                 :       1008 : F4(GEN x)
     416                 :            : {
     417                 :       1008 :   return gadd(
     418                 :       4032 :     gmul(gel(x,1), gadd(gsqr(gel(x,2)), gmul(gel(x,4),gel(x,1)))),
     419                 :       4032 :     gmul(gel(x,3), gadd(gsqr(gel(x,4)), gmul(gel(x,2),gel(x,3)))));
     420                 :            : }
     421                 :            : 
     422                 :            : static GEN
     423                 :       1332 : roots_to_ZX(GEN z, long *e)
     424                 :            : {
     425                 :       1332 :   GEN a = roots_to_pol(z,0);
     426                 :       1332 :   GEN b = grndtoi(real_i(a),e);
     427                 :       1332 :   long e1 = gexpo(imag_i(a));
     428         [ +  + ]:       1332 :   if (e1 > *e) *e = e1;
     429                 :       1332 :   return b;
     430                 :            : }
     431                 :            : 
     432                 :            : static GEN
     433                 :       1015 : polgaloisnames(long a, long b)
     434                 :            : {
     435                 :       1015 :   const char * const t[]={"S1", "S2", "A3", "S3",
     436                 :            :        "C(4) = 4", "E(4) = 2[x]2", "D(4)", "A4", "S4",
     437                 :            :        "C(5) = 5", "D(5) = 5:2", "F(5) = 5:4", "A5", "S5",
     438                 :            :        "C(6) = 6 = 3[x]2", "D_6(6) = [3]2", "D(6) = S(3)[x]2",
     439                 :            :              "A_4(6) = [2^2]3", "F_18(6) = [3^2]2 = 3 wr 2",
     440                 :            :              "2A_4(6) = [2^3]3 = 2 wr 3", "S_4(6d) = [2^2]S(3)",
     441                 :            :              "S_4(6c) = 1/2[2^3]S(3)", "F_18(6):2 = [1/2.S(3)^2]2",
     442                 :            :              "F_36(6) = 1/2[S(3)^2]2", "2S_4(6) = [2^3]S(3) = 2 wr S(3)",
     443                 :            :              "L(6) = PSL(2,5) = A_5(6)", "F_36(6):2 = [S(3)^2]2 = S(3) wr 2",
     444                 :            :              "L(6):2 = PGL(2,5) = S_5(6)", "A6", "S6",
     445                 :            :        "C(7) = 7", "D(7) = 7:2", "F_21(7) = 7:3", "F_42(7) = 7:6",
     446                 :            :              "L(7) = L(3,2)", "A7", "S7"};
     447                 :            : 
     448                 :       1015 :    const long idx[]={0,1,2,4,9,14,30};
     449                 :       1015 :    return strtoGENstr(t[idx[a-1]+b-1]);
     450                 :            : }
     451                 :            : 
     452                 :            : static GEN
     453                 :       1015 : galois_res(long d, long n, long s, long k)
     454                 :            : {
     455                 :       1015 :   long kk = k;
     456                 :       1015 :   GEN z = cgetg(5,t_VEC);
     457         [ +  + ]:       1015 :   if (!new_galois_format)
     458                 :            :   {
     459         [ +  + ]:        511 :     switch (d) {
     460                 :            :       case 6:
     461 [ +  + ][ +  + ]:        231 :         kk = (k == 6 || k == 2)? 2: 1;
     462                 :        231 :         break;
     463                 :            :       default:
     464                 :        280 :         kk = 1;
     465                 :            :     }
     466                 :            :   }
     467                 :       1015 :   gel(z,1) = stoi(n);
     468                 :       1015 :   gel(z,2) = stoi(s);
     469                 :       1015 :   gel(z,3) = stoi(kk);
     470                 :       1015 :   gel(z,4) = polgaloisnames(d,k);
     471                 :       1015 :   return z;
     472                 :            : }
     473                 :            : 
     474                 :            : GEN
     475                 :       3003 : polgalois(GEN x, long prec)
     476                 :            : {
     477                 :       3003 :   pari_sp av = avma, av1;
     478                 :            :   long i,j,k,n,f,l,l2,e,e1,pr,ind;
     479                 :            :   GEN x1,p1,p2,p3,p4,p5,w,z,ee;
     480                 :       3003 :   const int ind5[20]={2,5,3,4, 1,3,4,5, 1,5,2,4, 1,2,3,5, 1,4,2,3};
     481                 :       3003 :   const int ind6[60]={3,5,4,6, 2,6,4,5, 2,3,5,6, 2,4,3,6, 2,5,3,4,
     482                 :            :                       1,4,5,6, 1,5,3,6, 1,6,3,4, 1,3,4,5, 1,6,2,5,
     483                 :            :                       1,2,4,6, 1,5,2,4, 1,3,2,6, 1,2,3,5, 1,4,2,3};
     484         [ -  + ]:       3003 :   if (typ(x)!=t_POL) pari_err_TYPE("galois",x);
     485                 :       3003 :   n=degpol(x);
     486         [ -  + ]:       3003 :   if (n>11) pari_err_IMPL("galois of degree higher than 11");
     487                 :       3003 :   x = Q_primpart(x);
     488                 :       3003 :   RgX_check_ZX(x, "galois");
     489         [ -  + ]:       3003 :   if (!ZX_is_irred(x)) pari_err_IRREDPOL("galois",x);
     490                 :            : 
     491         [ +  + ]:       3003 :   if (n<4)
     492                 :            :   {
     493         [ -  + ]:         84 :     if (n == 1) { avma = av; return galois_res(n,1, 1,1); }
     494         [ +  + ]:         84 :     if (n == 2) { avma = av; return galois_res(n,2,-1,1); }
     495                 :            :     /* n = 3 */
     496                 :         56 :     f = Z_issquare(ZX_disc(x));
     497                 :         56 :     avma = av;
     498         [ +  + ]:         56 :     return f? galois_res(n,3,1,1):
     499                 :            :               galois_res(n,6,-1,2);
     500                 :            :   }
     501                 :       2919 :   x1 = x = ZX_Q_normalize(x,NULL); av1=avma;
     502         [ +  + ]:       2919 :   if (n > 7) return galoisbig(x, prec);
     503                 :            :   for(;;)
     504                 :            :   {
     505                 :       1120 :     double fb = fujiwara_bound(x);
     506   [ +  +  +  +  :       1120 :     switch(n)
                      - ]
     507                 :            :     {
     508                 :        168 :       case 4: z = cgetg(7,t_VEC);
     509                 :        168 :         prec = nbits2prec((long)(fb*18.) + 64);
     510                 :            :         for(;;)
     511                 :            :         {
     512                 :        168 :           p1=QX_complex_roots(x,prec);
     513                 :        168 :           gel(z,1) = F4(p1);
     514                 :        168 :           gel(z,2) = F4(transroot(p1,1,2));
     515                 :        168 :           gel(z,3) = F4(transroot(p1,1,3));
     516                 :        168 :           gel(z,4) = F4(transroot(p1,1,4));
     517                 :        168 :           gel(z,5) = F4(transroot(p1,2,3));
     518                 :        168 :           gel(z,6) = F4(transroot(p1,3,4));
     519         [ +  - ]:        168 :           p5 = roots_to_ZX(z, &e); if (e <= -10) break;
     520                 :          0 :           prec = precdbl(prec);
     521                 :          0 :         }
     522         [ +  + ]:        168 :         if (!ZX_is_squarefree(p5)) goto tchi;
     523                 :        140 :         p2 = gel(ZX_factor(p5),1);
     524   [ +  +  +  - ]:        140 :         switch(lg(p2)-1)
     525                 :            :         {
     526                 :         56 :           case 1: f = Z_issquare(ZX_disc(x)); avma = av;
     527         [ +  + ]:         56 :             return f? galois_res(n,12,1,4): galois_res(n,24,-1,5);
     528                 :            : 
     529                 :         28 :           case 2: avma = av; return galois_res(n,8,-1,3);
     530                 :            : 
     531                 :         56 :           case 3: avma = av;
     532                 :        112 :             return (degpol(gel(p2,1))==2)? galois_res(n,4,1,2)
     533         [ +  + ]:         56 :                                          : galois_res(n,4,-1,1);
     534                 :            : 
     535                 :          0 :           default: pari_err_BUG("galois (bug1)");
     536                 :            :         }
     537                 :            : 
     538                 :        154 :       case 5: z = cgetg(7,t_VEC);
     539                 :        154 :         ee= cgetg(7,t_VECSMALL);
     540                 :        154 :         w = cgetg(7,t_VECSMALL);
     541                 :        154 :         prec = nbits2prec((long)(fb*21.) + 64);
     542                 :            :         for(;;)
     543                 :            :         {
     544                 :            :           for(;;)
     545                 :            :           {
     546                 :        154 :             p1=QX_complex_roots(x,prec);
     547         [ +  + ]:       1078 :             for (l=1; l<=6; l++)
     548                 :            :             {
     549 [ +  + ][ +  + ]:        924 :               p2=(l==1)?p1: ((l<6)?transroot(p1,1,l): transroot(p1,2,5));
     550                 :        924 :               p3=gen_0;
     551         [ +  + ]:       5544 :               for (k=0,i=1; i<=5; i++,k+=4)
     552                 :            :               {
     553                 :       4620 :                 p5 = gadd(gmul(gel(p2,ind5[k]),gel(p2,ind5[k+1])),
     554                 :       9240 :                           gmul(gel(p2,ind5[k+2]),gel(p2,ind5[k+3])));
     555                 :       4620 :                 p3 = gadd(p3, gmul(gsqr(gel(p2,i)),p5));
     556                 :            :               }
     557                 :        924 :               gel(w,l) = grndtoi(real_i(p3),&e);
     558         [ +  + ]:        924 :               e1 = gexpo(imag_i(p3)); if (e1>e) e=e1;
     559                 :        924 :               ee[l]=e; gel(z,l) = p3;
     560                 :            :             }
     561         [ +  - ]:        154 :             p5 = roots_to_ZX(z, &e); if (e <= -10) break;
     562                 :          0 :             prec = precdbl(prec);
     563                 :          0 :           }
     564         [ -  + ]:        154 :           if (!ZX_is_squarefree(p5)) goto tchi;
     565                 :        154 :           p3=gel(ZX_factor(p5),1);
     566                 :        154 :           f=Z_issquare(ZX_disc(x));
     567         [ +  + ]:        154 :           if (lg(p3)-1==1)
     568                 :            :           {
     569                 :         56 :             avma = av;
     570         [ +  + ]:         56 :             return f? galois_res(n,60,1,4): galois_res(n,120,-1,5);
     571                 :            :           }
     572         [ +  + ]:         98 :           if (!f) { avma = av; return galois_res(n,20,-1,3); }
     573                 :            : 
     574                 :         70 :           pr = - (prec2nbits(prec) >> 1);
     575         [ +  - ]:        280 :           for (l=1; l<=6; l++)
     576 [ +  + ][ +  - ]:        280 :             if (ee[l] <= pr && gequal0(poleval(p5,gel(w,l)))) break;
     577         [ -  + ]:         70 :           if (l>6) pari_err_BUG("galois (bug4)");
     578         [ +  + ]:         70 :           p2=(l==6)? transroot(p1,2,5):transroot(p1,1,l);
     579                 :         70 :           p3=gen_0;
     580         [ +  + ]:        420 :           for (i=1; i<=5; i++)
     581                 :            :           {
     582         [ +  + ]:        350 :             j = (i == 5)? 1: i+1;
     583                 :        350 :             p3 = gadd(p3,gmul(gmul(gel(p2,i),gel(p2,j)),
     584                 :        700 :                               gsub(gel(p2,j),gel(p2,i))));
     585                 :            :           }
     586                 :         70 :           p5=gsqr(p3); p4=grndtoi(real_i(p5),&e);
     587         [ +  + ]:         70 :           e1 = gexpo(imag_i(p5)); if (e1>e) e=e1;
     588         [ +  - ]:         70 :           if (e <= -10)
     589                 :            :           {
     590         [ +  + ]:         70 :             if (gequal0(p4)) goto tchi;
     591                 :         56 :             f = Z_issquare(p4); avma = av;
     592         [ +  + ]:         56 :             return f? galois_res(n,5,1,1): galois_res(n,10,1,2);
     593                 :            :           }
     594                 :          0 :           prec = precdbl(prec);
     595                 :          0 :         }
     596                 :            : 
     597                 :        602 :       case 6: z = cgetg(7, t_VEC);
     598                 :        602 :         prec = nbits2prec((long) (fb * 42) + 64);
     599                 :            :         for(;;)
     600                 :            :         {
     601                 :            :           for(;;)
     602                 :            :           {
     603                 :        602 :             p1=QX_complex_roots(x,prec);
     604         [ +  + ]:       4214 :             for (l=1; l<=6; l++)
     605                 :            :             {
     606         [ +  + ]:       3612 :               p2=(l==1)?p1:transroot(p1,1,l);
     607                 :       3612 :               p3=gen_0; k=0;
     608 [ +  + ][ +  + ]:      75852 :               for (i=1; i<=5; i++) for (j=i+1; j<=6; j++)
     609                 :            :               {
     610                 :      54180 :                 p5=gadd(gmul(gel(p2,ind6[k]),gel(p2,ind6[k+1])),
     611                 :     108360 :                         gmul(gel(p2,ind6[k+2]),gel(p2,ind6[k+3])));
     612                 :      54180 :                 p3=gadd(p3,gmul(gsqr(gmul(gel(p2,i),gel(p2,j))),p5));
     613                 :      54180 :                 k += 4;
     614                 :            :               }
     615                 :       3612 :               gel(z,l) = p3;
     616                 :            :             }
     617         [ +  - ]:        602 :             p5 = roots_to_ZX(z, &e); if (e <= -10) break;
     618                 :          0 :             prec = precdbl(prec);
     619                 :          0 :           }
     620         [ +  + ]:        602 :           if (!ZX_is_squarefree(p5)) goto tchi;
     621                 :        455 :           p2=gel(ZX_factor(p5),1);
     622   [ +  +  +  +  :        455 :           switch(lg(p2)-1)
                      - ]
     623                 :            :           {
     624                 :            :             case 1:
     625                 :        112 :               z = cgetg(11,t_VEC); ind=0;
     626                 :        112 :               p3=gadd(gmul(gmul(gel(p1,1),gel(p1,2)),gel(p1,3)),
     627                 :        336 :                       gmul(gmul(gel(p1,4),gel(p1,5)),gel(p1,6)));
     628                 :        112 :               gel(z,++ind) = p3;
     629         [ +  + ]:        448 :               for (i=1; i<=3; i++)
     630         [ +  + ]:       1344 :                 for (j=4; j<=6; j++)
     631                 :            :                 {
     632                 :       1008 :                   p2=transroot(p1,i,j);
     633                 :       1008 :                   p3=gadd(gmul(gmul(gel(p2,1),gel(p2,2)),gel(p2,3)),
     634                 :       3024 :                           gmul(gmul(gel(p2,4),gel(p2,5)),gel(p2,6)));
     635                 :       1008 :                   gel(z,++ind) = p3;
     636                 :            :                 }
     637                 :        112 :               p5 = roots_to_ZX(z, &e);
     638         [ +  - ]:        112 :               if (e <= -10)
     639                 :            :               {
     640         [ -  + ]:        112 :                 if (!ZX_is_squarefree(p5)) goto tchi;
     641                 :        112 :                 p2 = gel(ZX_factor(p5),1);
     642                 :        112 :                 f = Z_issquare(ZX_disc(x));
     643                 :        112 :                 avma = av;
     644         [ +  + ]:        112 :                 if (lg(p2)-1==1)
     645         [ +  + ]:         56 :                   return f? galois_res(n,360,1,15): galois_res(n,720,-1,16);
     646                 :            :                 else
     647         [ +  + ]:         56 :                   return f? galois_res(n,36,1,10): galois_res(n,72,-1,13);
     648                 :            :               }
     649                 :          0 :               prec = precdbl(prec); break;
     650                 :            : 
     651         [ -  + ]:        196 :             case 2: l2=degpol(gel(p2,1)); if (l2>3) l2=6-l2;
     652   [ +  +  +  - ]:        196 :               switch(l2)
     653                 :            :               {
     654                 :         56 :                 case 1: f = Z_issquare(ZX_disc(x));
     655                 :         56 :                   avma = av;
     656         [ +  + ]:         56 :                   return f? galois_res(n,60,1,12): galois_res(n,120,-1,14);
     657                 :         84 :                 case 2: f = Z_issquare(ZX_disc(x));
     658         [ +  + ]:         84 :                   if (f) { avma = av; return galois_res(n,24,1,7); }
     659         [ +  - ]:         56 :                   p3 = (degpol(gel(p2,1))==2)? gel(p2,2): gel(p2,1);
     660                 :         56 :                   f = Z_issquare(ZX_disc(p3));
     661                 :         56 :                   avma = av;
     662         [ +  + ]:         56 :                   return f? galois_res(n,24,-1,6): galois_res(n,48,-1,11);
     663         [ +  - ]:         64 :                 case 3: f = Z_issquare(ZX_disc(gel(p2,1)))
     664         [ +  - ]:         56 :                          || Z_issquare(ZX_disc(gel(p2,2)));
           [ +  +  +  + ]
     665                 :         56 :                   avma = av;
     666         [ +  + ]:         56 :                   return f? galois_res(n,18,-1,5): galois_res(n,36,-1,9);
     667                 :            :               }
     668                 :            :             case 3:
     669         [ +  + ]:        476 :               for (l2=1; l2<=3; l2++)
     670         [ +  + ]:        357 :                 if (degpol(gel(p2,l2)) >= 3) p3 = gel(p2,l2);
     671         [ +  + ]:        119 :               if (degpol(p3) == 3)
     672                 :            :               {
     673                 :         56 :                 f = Z_issquare(ZX_disc(p3)); avma = av;
     674         [ +  + ]:         56 :                 return f? galois_res(n,6,-1,1): galois_res(n,12,-1,3);
     675                 :            :               }
     676                 :            :               else
     677                 :            :               {
     678                 :         63 :                 f = Z_issquare(ZX_disc(x)); avma = av;
     679         [ +  + ]:         63 :                 return f? galois_res(n,12,1,4): galois_res(n,24,-1,8);
     680                 :            :               }
     681                 :         28 :             case 4: avma = av; return galois_res(n,6,-1,2);
     682                 :          0 :             default: pari_err_BUG("galois (bug3)");
     683                 :            :           }
     684                 :          0 :         }
     685                 :            : 
     686                 :        196 :       case 7: z = cgetg(36,t_VEC);
     687                 :        196 :         prec = nbits2prec((long)(fb*7.) + 64);
     688                 :            :         for(;;)
     689                 :            :         {
     690                 :        296 :           ind = 0; p1=QX_complex_roots(x,prec);
     691         [ +  + ]:       1776 :           for (i=1; i<=5; i++)
     692         [ +  + ]:       5920 :             for (j=i+1; j<=6; j++)
     693                 :            :             {
     694                 :       4440 :               GEN t = gadd(gel(p1,i),gel(p1,j));
     695         [ +  + ]:      14800 :               for (k=j+1; k<=7; k++) gel(z,++ind) = gadd(t, gel(p1,k));
     696                 :            :             }
     697         [ +  + ]:        296 :           p5 = roots_to_ZX(z, &e); if (e <= -10) break;
     698                 :        100 :           prec = precdbl(prec);
     699                 :        100 :         }
     700         [ -  + ]:        196 :         if (!ZX_is_squarefree(p5)) goto tchi;
     701                 :        196 :         p2=gel(ZX_factor(p5),1);
     702   [ +  +  +  +  :        196 :         switch(lg(p2)-1)
                   +  - ]
     703                 :            :         {
     704                 :         56 :           case 1: f = Z_issquare(ZX_disc(x)); avma = av;
     705         [ +  + ]:         56 :             return f? galois_res(n,2520,1,6): galois_res(n,5040,-1,7);
     706                 :         56 :           case 2: f = degpol(gel(p2,1)); avma = av;
     707 [ +  + ][ -  + ]:         56 :             return (f==7 || f==28)? galois_res(n,168,1,5): galois_res(n,42,-1,4);
     708                 :         28 :           case 3: avma = av; return galois_res(n,21,1,3);
     709                 :         28 :           case 4: avma = av; return galois_res(n,14,-1,2);
     710                 :         28 :           case 5: avma = av; return galois_res(n,7,1,1);
     711                 :          0 :           default: pari_err_BUG("galois (bug2)");
     712                 :            :         }
     713                 :            :     }
     714                 :       3192 :     tchi: avma = av1; x = tschirnhaus(x1);
     715                 :        189 :   }
     716                 :            : }
     717                 :            : 
     718                 :            : #undef _res
     719                 :            : 
     720                 :            : /* Evaluate pol in s using nfelt arithmetic and Horner rule */
     721                 :            : GEN
     722                 :      11417 : nfpoleval(GEN nf, GEN pol, GEN s)
     723                 :            : {
     724                 :      11417 :   pari_sp av=avma;
     725                 :      11417 :   long i=lg(pol)-1;
     726                 :            :   GEN res;
     727         [ -  + ]:      11417 :   if (i==1) return gen_0;
     728                 :      11417 :   res = nf_to_scalar_or_basis(nf, gel(pol,i));
     729         [ +  + ]:      28609 :   for (i-- ; i>=2; i--)
     730                 :      17192 :     res = nfadd(nf, nfmul(nf, s, res), gel(pol,i));
     731                 :      11417 :   return gerepileupto(av, res);
     732                 :            : }
     733                 :            : 
     734                 :            : static GEN
     735                 :      11914 : QX_table_nfpoleval(GEN nf, GEN pol, GEN m)
     736                 :            : {
     737                 :      11914 :   pari_sp av = avma;
     738                 :      11914 :   long i = lg(pol)-1;
     739                 :            :   GEN res, den;
     740         [ -  + ]:      11914 :   if (i==1) return gen_0;
     741                 :      11914 :   pol = Q_remove_denom(pol, &den);
     742                 :      11914 :   res = scalarcol_shallow(gel(pol,i), nf_get_degree(nf));
     743         [ +  + ]:      73731 :   for (i-- ; i>=2; i--)
     744                 :      61817 :     res = ZC_Z_add(ZM_ZC_mul(m, res), gel(pol,i));
     745         [ +  + ]:      11914 :   if (den) res = RgC_Rg_div(res, den);
     746                 :      11914 :   return gerepileupto(av, res);
     747                 :            : }
     748                 :            : 
     749                 :            : GEN
     750                 :         63 : FpX_FpC_nfpoleval(GEN nf, GEN pol, GEN a, GEN p)
     751                 :            : {
     752                 :         63 :   pari_sp av=avma;
     753                 :         63 :   long i=lg(pol)-1, n=nf_get_degree(nf);
     754                 :            :   GEN res, Ma;
     755         [ -  + ]:         63 :   if (i==1) return zerocol(n);
     756                 :         63 :   Ma = FpM_red(zk_multable(nf, a), p);
     757                 :         63 :   res = scalarcol(gel(pol,i),n);
     758         [ +  + ]:        168 :   for (i-- ; i>=2; i--)
     759                 :            :   {
     760                 :        105 :     res = FpM_FpC_mul(Ma, res, p);
     761                 :        105 :     gel(res,1) = Fp_add(gel(res,1), gel(pol,i), p);
     762                 :            :   }
     763                 :         63 :   return gerepileupto(av, res);
     764                 :            : }
     765                 :            : 
     766                 :            : /* compute s(x), not stack clean */
     767                 :            : static GEN
     768                 :       3738 : table_galoisapply(GEN nf, GEN m, GEN x)
     769                 :            : {
     770                 :       3738 :   x = nf_to_scalar_or_alg(nf, x);
     771         [ +  + ]:       3738 :   if (typ(x) != t_POL) return scalarcol(x, nf_get_degree(nf));
     772                 :       3738 :   return QX_table_nfpoleval(nf, x, m);
     773                 :            : }
     774                 :            : 
     775                 :            : /* compute s(x), not stack clean */
     776                 :            : static GEN
     777                 :       3556 : ZC_galoisapply(GEN nf, GEN s, GEN x)
     778                 :            : {
     779                 :       3556 :   x = nf_to_scalar_or_alg(nf, x);
     780         [ +  + ]:       3556 :   if (typ(x) != t_POL) return scalarcol(x, nf_get_degree(nf));
     781                 :       3556 :   return QX_table_nfpoleval(nf, x, zk_multable(nf, s));
     782                 :            : }
     783                 :            : 
     784                 :            : static GEN
     785                 :         63 : QX_galoisapplymod(GEN nf, GEN pol, GEN S, GEN p)
     786                 :            : {
     787                 :         63 :   GEN den, P = Q_remove_denom(pol,&den);
     788                 :            :   GEN pe, pe1, denpe, R;
     789         [ +  + ]:         63 :   if (den)
     790                 :            :   {
     791                 :          7 :     ulong e = Z_pval(den, p);
     792                 :          7 :     pe = powiu(p, e); pe1 = mulii(pe, p);
     793                 :          7 :     denpe = Fp_inv(diviiexact(den, pe), pe1);
     794                 :            :   } else {
     795                 :         56 :     pe = gen_1; pe1 = p; denpe = gen_1;
     796                 :            :   }
     797                 :         63 :   R = FpX_FpC_nfpoleval(nf, FpX_red(P, pe1), FpC_red(S, pe1), pe1);
     798                 :         63 :   return gdivexact(FpC_Fp_mul(R, denpe, pe1), pe);
     799                 :            : }
     800                 :            : 
     801                 :            : static GEN
     802                 :          7 : pr_galoisapply(GEN nf, GEN pr, GEN aut)
     803                 :            : {
     804                 :            :   GEN p, t, u;
     805         [ -  + ]:          7 :   if (typ(pr_get_tau(pr)) == t_INT) return pr; /* inert */
     806                 :          7 :   p = pr_get_p(pr);
     807                 :          7 :   u = QX_galoisapplymod(nf, coltoliftalg(nf, pr_get_gen(pr)), aut, p);
     808                 :          7 :   t = FpM_deplin(zk_multable(nf, u), p);
     809                 :          7 :   t = zk_scalar_or_multable(nf, t);
     810                 :          7 :   return mkvec5(p, u, gel(pr,3), gel(pr,4), t);
     811                 :            : }
     812                 :            : 
     813                 :            : static GEN
     814                 :          7 : vecgaloisapply(GEN nf, GEN aut, GEN v)
     815                 :            : {
     816                 :            :   long i, l;
     817                 :          7 :   GEN V = cgetg_copy(v, &l);
     818         [ +  + ]:         21 :   for (i = 1; i < l; i++) gel(V,i) = galoisapply(nf, aut, gel(v,i));
     819                 :          7 :   return V;
     820                 :            : }
     821                 :            : 
     822                 :            : /* x: famat or standard algebraic number, aut automorphism in ZC form
     823                 :            :  * simplified from general galoisapply */
     824                 :            : static GEN
     825                 :         49 : elt_galoisapply(GEN nf, GEN aut, GEN x)
     826                 :            : {
     827                 :         49 :   pari_sp av = avma;
     828   [ +  +  +  +  :         49 :   switch(typ(x))
                +  +  - ]
     829                 :            :   {
     830                 :          7 :     case t_INT:  return icopy(x);
     831                 :          7 :     case t_FRAC: return gcopy(x);
     832                 :          7 :     case t_POLMOD: x = gel(x,2); /* fall through */
     833                 :            :     case t_POL: {
     834                 :         14 :       GEN y = basistoalg(nf, ZC_galoisapply(nf, aut, x));
     835                 :         14 :       return gerepileupto(av,y);
     836                 :            :     }
     837                 :            :     case t_COL:
     838                 :          7 :       return gerepileupto(av, ZC_galoisapply(nf, aut, x));
     839                 :            :     case t_MAT:
     840      [ +  +  - ]:         14 :       switch(lg(x)) {
     841                 :          7 :         case 1: return cgetg(1, t_MAT);
     842                 :          7 :         case 3: retmkmat2(vecgaloisapply(nf,aut,gel(x,1)), ZC_copy(gel(x,2)));
     843                 :            :       }
     844                 :            :   }
     845                 :          0 :   pari_err_TYPE("galoisapply",x);
     846                 :         49 :   return NULL; /* not reached */
     847                 :            : }
     848                 :            : 
     849                 :            : GEN
     850                 :       4074 : galoisapply(GEN nf, GEN aut, GEN x)
     851                 :            : {
     852                 :       4074 :   pari_sp av = avma;
     853                 :            :   long lx, j;
     854                 :            :   GEN y;
     855                 :            : 
     856                 :       4074 :   nf = checknf(nf);
     857   [ +  +  +  +  :       4074 :   switch(typ(x))
             +  +  +  - ]
     858                 :            :   {
     859                 :         56 :     case t_INT:  return icopy(x);
     860                 :          7 :     case t_FRAC: return gcopy(x);
     861                 :            : 
     862                 :          7 :     case t_POLMOD: x = gel(x,2); /* fall through */
     863                 :            :     case t_POL:
     864                 :         21 :       aut = algtobasis(nf, aut);
     865                 :         21 :       y = basistoalg(nf, ZC_galoisapply(nf, aut, x));
     866                 :         21 :       return gerepileupto(av,y);
     867                 :            : 
     868                 :            :     case t_VEC:
     869                 :         56 :       aut = algtobasis(nf, aut);
     870      [ +  +  - ]:         56 :       switch(lg(x))
     871                 :            :       {
     872                 :          7 :         case 6: return gerepilecopy(av, pr_galoisapply(nf, x, aut));
     873                 :         49 :         case 3: y = cgetg(3,t_VEC);
     874                 :         49 :           gel(y,1) = galoisapply(nf, aut, gel(x,1));
     875                 :         49 :           gel(y,2) = elt_galoisapply(nf, aut, gel(x,2));
     876                 :         49 :           return gerepileupto(av, y);
     877                 :            :       }
     878                 :          0 :       break;
     879                 :            : 
     880                 :            :     case t_COL:
     881                 :       3003 :       aut = algtobasis(nf, aut);
     882                 :       3003 :       return gerepileupto(av, ZC_galoisapply(nf, aut, x));
     883                 :            : 
     884                 :            :     case t_MAT: /* ideal */
     885         [ -  + ]:        931 :       lx = lg(x); if (lx==1) return cgetg(1,t_MAT);
     886         [ -  + ]:        931 :       if (nbrows(x) != nf_get_degree(nf)) break;
     887                 :        931 :       aut = zk_multable(nf, algtobasis(nf, aut));
     888                 :        931 :       y = cgetg(lx,t_MAT);
     889         [ +  + ]:       4669 :       for (j=1; j<lx; j++) gel(y,j) = table_galoisapply(nf, aut, gel(x,j));
     890                 :        931 :       return gerepileupto(av, idealhnf_shallow(nf,y));
     891                 :            :   }
     892                 :          0 :   pari_err_TYPE("galoisapply",x);
     893                 :       4074 :   return NULL; /* not reached */
     894                 :            : }
     895                 :            : 
     896                 :            : GEN
     897                 :       1540 : nfgaloismatrix(GEN nf, GEN s)
     898                 :            : {
     899                 :            :   GEN zk, M, m;
     900                 :            :   long k, l;
     901                 :       1540 :   nf = checknf(nf);
     902                 :       1540 :   zk = nf_get_zk(nf);
     903         [ -  + ]:       1540 :   if (typ(s) != t_COL) s = algtobasis(nf, s); /* left on stack for efficiency */
     904                 :       1540 :   m = zk_multable(nf, s);
     905                 :       1540 :   l = lg(s); M = cgetg(l, t_MAT);
     906                 :       1540 :   gel(M, 1) = col_ei(l-1, 1); /* s(1) = 1 */
     907         [ +  + ]:       7098 :   for (k = 2; k < l; k++)
     908                 :       5558 :     gel(M, k) = QX_table_nfpoleval(nf, gel(zk, k), m);
     909                 :       1540 :   return M;
     910                 :            : }
     911                 :            : 
     912                 :            : static GEN
     913                 :         84 : idealquasifrob(GEN nf, GEN gal, GEN pr, GEN subg, GEN *S)
     914                 :            : {
     915                 :         84 :   pari_sp av = avma;
     916                 :         84 :   long i, n = nf_get_degree(nf), f = pr_get_f(pr);
     917                 :         84 :   GEN grp = gal_get_group(gal), pi = pr_get_gen(pr);
     918         [ +  - ]:        994 :   for (i=1; i<=n; i++)
     919                 :            :   {
     920                 :        994 :     GEN g = gel(grp,i);
     921 [ +  + ][ +  + ]:        994 :     if ((!subg && perm_order(g)==f)
     922 [ +  + ][ +  + ]:        721 :       || (subg && perm_relorder(g, subg)==f))
     923                 :            :     {
     924                 :        294 :       *S = poltobasis(nf, galoispermtopol(gal, g));
     925         [ +  + ]:        294 :       if (ZC_prdvd(nf, ZC_galoisapply(nf, *S, pi), pr)) return g;
     926                 :        210 :       avma = av;
     927                 :            :     }
     928                 :            :   }
     929                 :          0 :   pari_err_BUG("idealquasifrob [Frobenius not found]");
     930                 :         84 :   return NULL; /*NOT REACHED*/
     931                 :            : }
     932                 :            : 
     933                 :            : static void
     934                 :        112 : gal_check_pol(const char *f, GEN x, GEN y)
     935         [ -  + ]:        112 : { if (!RgX_equal_var(x,y)) pari_err_MODULUS(f,x,y); }
     936                 :            : 
     937                 :            : GEN
     938                 :         14 : idealfrobenius(GEN nf, GEN gal, GEN pr)
     939                 :            : {
     940                 :         14 :   pari_sp av = avma;
     941                 :         14 :   GEN S=NULL, g=NULL; /*-Wall*/
     942                 :            :   GEN T, p, a, b, modpr;
     943                 :            :   long f, n, s;
     944                 :         14 :   nf = checknf(nf);
     945                 :         14 :   checkgal(gal);
     946                 :         14 :   checkprid(pr);
     947                 :         14 :   gal_check_pol("idealfrobenius",nf_get_pol(nf),gal_get_pol(gal));
     948         [ +  + ]:         14 :   if (pr_get_e(pr)>1) pari_err_DOMAIN("idealfrobenius","pr.e", ">", gen_1,pr);
     949                 :          7 :   f = pr_get_f(pr); n = nf_get_degree(nf);
     950         [ -  + ]:          7 :   if (f==1) { avma = av; return identity_perm(n); }
     951                 :          7 :   modpr = zk_to_Fq_init(nf,&pr,&T,&p);
     952                 :          7 :   g = idealquasifrob(nf, gal, pr, NULL, &S);
     953                 :          7 :   a = pol_x(nf_get_varn(nf));
     954                 :          7 :   b = nf_to_Fq(nf, QX_galoisapplymod(nf, modpr_genFq(modpr), S, p), modpr);
     955         [ +  + ]:         14 :   for (s=0; !ZX_equal(a, b); s++)
     956                 :          7 :     a = Fq_pow(a, p, T, p);
     957                 :          7 :   g = perm_pow(g, Fl_inv(s, f));
     958                 :          7 :   return gerepileupto(av, g);
     959                 :            : }
     960                 :            : 
     961                 :            : static GEN
     962                 :         21 : idealinertiagroup(GEN nf, GEN gal, GEN pr)
     963                 :            : {
     964                 :         21 :   long i, n = nf_get_degree(nf);
     965                 :         21 :   GEN p, T, modpr = zk_to_Fq_init(nf,&pr,&T,&p);
     966                 :         21 :   GEN b = modpr_genFq(modpr);
     967                 :         21 :   long e = pr_get_e(pr), coprime = cgcd(e, pr_get_f(pr)) == 1;
     968                 :         21 :   GEN grp = gal_get_group(gal), pi = pr_get_gen(pr);
     969                 :         21 :   pari_sp ltop = avma;
     970         [ +  - ]:        224 :   for (i=1; i<=n; i++)
     971                 :            :   {
     972                 :        224 :     GEN iso = gel(grp,i);
     973         [ +  + ]:        224 :     if (perm_order(iso) == e)
     974                 :            :     {
     975                 :         63 :       GEN S = poltobasis(nf, galoispermtopol(gal, iso));
     976         [ +  + ]:         63 :       if (ZC_prdvd(nf, ZC_galoisapply(nf, S, pi), pr)
     977 [ -  + ][ #  # ]:         21 :           && (coprime || gequalX(nf_to_Fq(nf, galoisapply(nf,S,b), modpr))))
     978                 :         21 :           return iso;
     979                 :         42 :       avma = ltop;
     980                 :            :     }
     981                 :            :   }
     982                 :          0 :   pari_err_BUG("idealinertiagroup [no isotropic element]");
     983                 :         21 :   return NULL;
     984                 :            : }
     985                 :            : 
     986                 :            : static GEN
     987                 :         84 : idealramgroupstame(GEN nf, GEN gal, GEN pr)
     988                 :            : {
     989                 :         84 :   pari_sp av = avma;
     990                 :            :   GEN iso, frob, giso, isog, S, res;
     991                 :         84 :   long e = pr_get_e(pr), f = pr_get_f(pr);
     992         [ +  + ]:         84 :   if (e == 1)
     993                 :            :   {
     994         [ -  + ]:         63 :     if (f==1)
     995                 :          0 :       return cgetg(1,t_VEC);
     996                 :         63 :     frob = idealquasifrob(nf, gal, pr, NULL, &S);
     997                 :         63 :     avma = av;
     998                 :         63 :     res = cgetg(2, t_VEC);
     999                 :         63 :     gel(res, 1) = cyclicgroup(frob, f);
    1000                 :         63 :     return res;
    1001                 :            :   }
    1002                 :         21 :   res = cgetg(3, t_VEC);
    1003                 :         21 :   av = avma;
    1004                 :         21 :   iso = idealinertiagroup(nf, gal, pr);
    1005                 :         21 :   avma = av;
    1006                 :         21 :   giso = cyclicgroup(iso, e);
    1007                 :         21 :   gel(res, 2) = giso;
    1008         [ +  + ]:         21 :   if (f==1)
    1009                 :            :   {
    1010                 :          7 :     gel(res, 1) = giso;
    1011                 :          7 :     return res;
    1012                 :            :   }
    1013                 :         14 :   av = avma;
    1014                 :         14 :   isog = group_set(giso, nf_get_degree(nf));
    1015                 :         14 :   frob = idealquasifrob(nf, gal, pr, isog, &S);
    1016                 :         14 :   avma = av;
    1017                 :         14 :   gel(res, 1) = dicyclicgroup(iso,frob,e,f);
    1018                 :         84 :   return res;
    1019                 :            : }
    1020                 :            : 
    1021                 :            : static GEN
    1022                 :         14 : idealramgroupindex(GEN nf, GEN gal, GEN pr)
    1023                 :            : {
    1024                 :         14 :   pari_sp av = avma;
    1025                 :            :   GEN p, T, g, idx, modpr;
    1026                 :            :   long i, e, f, n;
    1027                 :            :   ulong nt,rorder;
    1028                 :         14 :   GEN grp = vecvecsmall_sort(gal_get_group(gal));
    1029                 :         14 :   e = pr_get_e(pr); f = pr_get_f(pr); n = nf_get_degree(nf);
    1030                 :         14 :   modpr = zk_to_Fq_init(nf,&pr,&T,&p);
    1031                 :         14 :   (void) u_pvalrem(n,p,&nt);
    1032                 :         14 :   rorder = e*f*(n/nt);
    1033                 :         14 :   idx = const_vecsmall(n,-1);
    1034                 :         14 :   g = modpr_genFq(modpr);
    1035         [ +  + ]:        266 :   for (i=2; i<=n; i++)
    1036                 :            :   {
    1037                 :            :     GEN iso;
    1038                 :            :     long o;
    1039         [ -  + ]:        252 :     if (idx[i]>=0) continue;
    1040                 :        252 :     iso = gel(grp,i); o = perm_order(iso);
    1041         [ +  + ]:        252 :     if (rorder%o == 0)
    1042                 :            :     {
    1043                 :        154 :       GEN piso = iso;
    1044                 :        154 :       GEN S = poltobasis(nf, galoispermtopol(gal, iso));
    1045                 :        154 :       GEN pi = pr_get_gen(pr);
    1046                 :        154 :       GEN spi = ZC_galoisapply(nf, S, pi);
    1047                 :            :       long j;
    1048                 :        154 :       idx[i] = idealval(nf, gsub(spi,pi), pr);
    1049         [ +  + ]:        154 :       if (idx[i] >=1)
    1050                 :            :       {
    1051         [ +  + ]:         56 :         if (f>1)
    1052                 :            :         {
    1053                 :         49 :           GEN b = nf_to_Fq(nf, QX_galoisapplymod(nf, g, S, p), modpr);
    1054         [ +  + ]:         49 :           if (!gequalX(b)) idx[i] = 0;
    1055                 :            :         }
    1056                 :            :       }
    1057                 :         98 :       else idx[i] = -1;
    1058         [ -  + ]:        154 :       for(j=2;j<o;j++)
    1059                 :            :       {
    1060                 :          0 :         piso = perm_mul(piso,iso);
    1061         [ #  # ]:          0 :         if(cgcd(j,o)==1) idx[piso[1]] = idx[i];
    1062                 :            :       }
    1063                 :            :     }
    1064                 :            :   }
    1065                 :         14 :   return gerepileuptoleaf(av, idx);
    1066                 :            : }
    1067                 :            : 
    1068                 :            : GEN
    1069                 :         98 : idealramgroups(GEN nf, GEN gal, GEN pr)
    1070                 :            : {
    1071                 :         98 :   pari_sp av = avma;
    1072                 :            :   GEN tbl, idx, res, set, sub;
    1073                 :            :   long i, j, e, n, maxm, p;
    1074                 :            :   ulong et;
    1075                 :         98 :   nf = checknf(nf);
    1076                 :         98 :   checkgal(gal);
    1077                 :         98 :   checkprid(pr);
    1078                 :         98 :   gal_check_pol("idealramgroups",nf_get_pol(nf),gal_get_pol(gal));
    1079                 :         98 :   e = pr_get_e(pr); n = nf_get_degree(nf);
    1080                 :         98 :   p = itos(pr_get_p(pr));
    1081         [ +  + ]:         98 :   if (e%p) return idealramgroupstame(nf, gal, pr);
    1082                 :         14 :   (void) u_lvalrem(e,p,&et);
    1083                 :         14 :   idx = idealramgroupindex(nf, gal, pr);
    1084                 :         14 :   sub = group_subgroups(galois_group(gal));
    1085                 :         14 :   tbl = subgroups_tableset(sub, n);
    1086                 :         14 :   maxm = vecsmall_max(idx)+1;
    1087                 :         14 :   res = cgetg(maxm+1,t_VEC);
    1088                 :         14 :   set = zero_F2v(n); F2v_set(set,1);
    1089         [ +  + ]:         77 :   for(i=maxm; i>0; i--)
    1090                 :            :   {
    1091         [ +  + ]:       1183 :     for(j=1;j<=n;j++)
    1092         [ +  + ]:       1120 :       if (idx[j]==i-1)
    1093                 :         56 :         F2v_set(set,j);
    1094                 :         63 :     gel(res,i) = gel(sub, tableset_find_index(tbl, set));
    1095                 :            :   }
    1096                 :         98 :   return gerepilecopy(av, res);
    1097                 :            : }
    1098                 :            : 
    1099                 :            : /* x = relative polynomial nf = absolute nf, bnf = absolute bnf */
    1100                 :            : GEN
    1101                 :        147 : get_bnfpol(GEN x, GEN *bnf, GEN *nf)
    1102                 :            : {
    1103                 :        147 :   *bnf = checkbnf_i(x);
    1104                 :        147 :   *nf  = checknf_i(x);
    1105         [ +  + ]:        147 :   if (*nf) x = nf_get_pol(*nf);
    1106         [ -  + ]:        147 :   if (typ(x) != t_POL) pari_err_TYPE("get_bnfpol",x);
    1107                 :        147 :   return x;
    1108                 :            : }
    1109                 :            : 
    1110                 :            : GEN
    1111                 :       7203 : get_nfpol(GEN x, GEN *nf)
    1112                 :            : {
    1113         [ +  + ]:       7203 :   if (typ(x) == t_POL) { *nf = NULL; return x; }
    1114                 :       7203 :   *nf = checknf(x); return nf_get_pol(*nf);
    1115                 :            : }
    1116                 :            : 
    1117                 :            : /* is isomorphism / inclusion (a \subset b) compatible with what we know about
    1118                 :            :  * basic invariants ? (degree, signature, discriminant) */
    1119                 :            : static int
    1120                 :         49 : tests_OK(GEN a, GEN nfa, GEN b, GEN nfb, long fliso)
    1121                 :            : {
    1122                 :            :   GEN da, db, fa, P, E, U;
    1123                 :         49 :   long i, nP, m = degpol(a), n = degpol(b), q = m / n; /* relative degree */
    1124                 :            : 
    1125         [ -  + ]:         49 :   if (m <= 0) pari_err_IRREDPOL("nfisincl",a);
    1126         [ -  + ]:         49 :   if (n <= 0) pari_err_IRREDPOL("nfisincl",b);
    1127 [ +  + ][ -  + ]:         49 :   if (fliso) { if (n != m) return 0; } else { if (n % m) return 0; }
                 [ -  + ]
    1128         [ +  + ]:         49 :   if (m == 1) return 1;
    1129                 :            : 
    1130 [ -  + ][ #  # ]:         42 :   if (nfa && nfb) /* both nf structures available */
    1131                 :            :   {
    1132                 :          0 :     long r1a = nf_get_r1(nfa), r1b = nf_get_r1(nfb) ;
    1133         [ #  # ]:          0 :     if (fliso)
    1134 [ #  # ][ #  # ]:          0 :       return (r1a == r1b && equalii(nf_get_disc(nfa), nf_get_disc(nfb)));
    1135                 :            :     else
    1136   [ #  #  #  # ]:          0 :       return (r1b <= r1a * q &&
    1137                 :          0 :               dvdii(nf_get_disc(nfb), powiu(nf_get_disc(nfa), q)));
    1138                 :            :   }
    1139         [ -  + ]:         42 :   da = nfa? nf_get_disc(nfa): ZX_disc(a);
    1140         [ +  + ]:         42 :   if (!signe(da)) pari_err_IRREDPOL("nfisincl",a);
    1141         [ +  + ]:         35 :   db = nfb? nf_get_disc(nfb): ZX_disc(b);
    1142         [ -  + ]:         35 :   if (!signe(db)) pari_err_IRREDPOL("nfisincl",a);
    1143         [ +  + ]:         35 :   if (fliso) return issquare(gdiv(da,db));
    1144                 :            : 
    1145 [ +  + ][ -  + ]:         21 :   if (odd(q) && signe(da) != signe(db)) return 0;
    1146                 :         21 :   fa = absi_factor_limit(da, 0);
    1147                 :         21 :   P = gel(fa,1);
    1148                 :         21 :   E = gel(fa,2); nP = lg(P) - 1;
    1149         [ +  + ]:         77 :   for (i=1; i<nP; i++)
    1150 [ +  + ][ -  + ]:         56 :     if (mod2(gel(E,i)) && !dvdii(db, powiu(gel(P,i),q))) return 0;
    1151                 :         21 :   U = gel(P,nP);
    1152 [ -  + ][ #  # ]:         21 :   if (mod2(gel(E,i)) && expi(U) < 150)
    1153                 :            :   { /* "unfactored" cofactor is small, finish */
    1154         [ #  # ]:          0 :     if (cmpiu(U, maxprime()) > 0)
    1155                 :            :     {
    1156                 :          0 :       fa = Z_factor(U);
    1157                 :          0 :       P = gel(fa,1);
    1158                 :          0 :       E = gel(fa,2);
    1159                 :            :     }
    1160                 :            :     else
    1161                 :            :     {
    1162                 :          0 :       P = mkvec(U);
    1163                 :          0 :       E = mkvec(gen_1);
    1164                 :            :     }
    1165                 :          0 :     nP = lg(P) - 1;
    1166         [ #  # ]:          0 :     for (i=1; i<=nP; i++)
    1167 [ #  # ][ #  # ]:          0 :       if (mod2(gel(E,i)) && !dvdii(db, powiu(gel(P,i),q))) return 0;
    1168                 :            :   }
    1169                 :         42 :   return 1;
    1170                 :            : }
    1171                 :            : 
    1172                 :            : /* if fliso test for isomorphism, for inclusion otherwise. */
    1173                 :            : static GEN
    1174                 :         49 : nfiso0(GEN a, GEN b, long fliso)
    1175                 :            : {
    1176                 :         49 :   pari_sp av = avma;
    1177                 :            :   long i, vb, lx;
    1178                 :            :   GEN nfa, nfb, y, la, lb;
    1179                 :            : 
    1180                 :         49 :   a = get_nfpol(a, &nfa);
    1181                 :         49 :   b = get_nfpol(b, &nfb);
    1182         [ +  + ]:         49 :   if (!nfa) { a = Q_primpart(a); RgX_check_ZX(a, "nsiso0"); }
    1183         [ +  + ]:         49 :   if (!nfb) { b = Q_primpart(b); RgX_check_ZX(b, "nsiso0"); }
    1184 [ +  + ][ -  + ]:         49 :   if (fliso && nfa && !nfb) { swap(a,b); nfb = nfa; nfa = NULL; }
                 [ #  # ]
    1185         [ -  + ]:         49 :   if (!tests_OK(a, nfa, b, nfb, fliso)) { avma = av; return gen_0; }
    1186                 :            : 
    1187         [ +  + ]:         42 :   if (nfb) lb = gen_1; else b = ZX_Q_normalize(b,&lb);
    1188         [ +  + ]:         42 :   if (nfa) la = gen_1; else a = ZX_Q_normalize(a,&la);
    1189                 :         42 :   a = leafcopy(a); setvarn(a,0);
    1190                 :         42 :   b = leafcopy(b); vb = varn(b);
    1191         [ +  + ]:         42 :   if (nfb)
    1192                 :            :   {
    1193         [ +  - ]:         14 :     if (vb == 0) nfb = gsubst(nfb, 0, pol_x(MAXVARN));
    1194                 :         14 :     y = lift_intern(nfroots(nfb,a));
    1195                 :            :   }
    1196                 :            :   else
    1197                 :            :   {
    1198         [ +  + ]:         28 :     if (vb == 0) setvarn(b, fetch_var());
    1199                 :         28 :     y = gel(polfnf(a,b),1); lx = lg(y);
    1200         [ +  + ]:        154 :     for (i=1; i<lx; i++)
    1201                 :            :     {
    1202                 :        126 :       GEN t = gel(y,i);
    1203         [ -  + ]:        126 :       if (degpol(t) != 1) { setlg(y,i); break; }
    1204                 :        126 :       gel(y,i) = gneg_i(lift_intern(gel(t,2)));
    1205                 :            :     }
    1206         [ +  + ]:         28 :     if (vb == 0) (void)delete_var();
    1207                 :         28 :     settyp(y, t_VEC);
    1208                 :         28 :     gen_sort_inplace(y, (void*)&cmp_RgX, &cmp_nodata, NULL);
    1209                 :            :   }
    1210         [ -  + ]:         42 :   lx = lg(y); if (lx==1) { avma=av; return gen_0; }
    1211         [ +  + ]:        189 :   for (i=1; i<lx; i++)
    1212                 :            :   {
    1213                 :        147 :     GEN t = gel(y,i);
    1214         [ +  + ]:        147 :     if (typ(t) == t_POL) setvarn(t, vb); else t = scalarpol(t, vb);
    1215         [ +  + ]:        147 :     if (lb != gen_1) t = RgX_unscale(t, lb);
    1216         [ +  + ]:        147 :     if (la != gen_1) t = RgX_Rg_div(t, la);
    1217                 :        147 :     gel(y,i) = t;
    1218                 :            :   }
    1219                 :         42 :   return gerepilecopy(av,y);
    1220                 :            : }
    1221                 :            : 
    1222                 :            : GEN
    1223                 :         14 : nfisisom(GEN a, GEN b) { return nfiso0(a,b,1); }
    1224                 :            : 
    1225                 :            : GEN
    1226                 :         35 : nfisincl(GEN a, GEN b) { return nfiso0(a,b,0); }
    1227                 :            : 
    1228                 :            : /*************************************************************************/
    1229                 :            : /**                                                                     **/
    1230                 :            : /**                               INITALG                               **/
    1231                 :            : /**                                                                     **/
    1232                 :            : /*************************************************************************/
    1233                 :            : 
    1234                 :            : GEN
    1235                 :       4943 : get_roots(GEN x, long r1, long prec)
    1236                 :            : {
    1237                 :            :   long i, ru;
    1238                 :            :   GEN z;
    1239         [ -  + ]:       4943 :   if (typ(x) != t_POL)
    1240                 :            :   {
    1241                 :          0 :     z = leafcopy(x);
    1242                 :          0 :     ru = (lg(z)-1 + r1) >> 1;
    1243                 :            :   }
    1244                 :            :   else
    1245                 :            :   {
    1246                 :       4943 :     long n = degpol(x);
    1247         [ +  + ]:       4943 :     z = (r1 == n)? realroots(x, NULL, prec): QX_complex_roots(x,prec);
    1248                 :       4943 :     ru = (n+r1)>>1;
    1249                 :            :   }
    1250         [ +  + ]:      13162 :   for (i=r1+1; i<=ru; i++) gel(z,i) = gel(z, (i<<1)-r1);
    1251                 :       4943 :   z[0]=evaltyp(t_VEC)|evallg(ru+1); return z;
    1252                 :            : }
    1253                 :            : 
    1254                 :            : GEN
    1255                 :          0 : nf_get_allroots(GEN nf)
    1256                 :            : {
    1257                 :          0 :   return embed_roots(nf_get_roots(nf), nf_get_r1(nf));
    1258                 :            : }
    1259                 :            : 
    1260                 :            : /* For internal use. compute trace(x mod pol), sym=polsym(pol,deg(pol)-1) */
    1261                 :            : GEN
    1262                 :      41965 : quicktrace(GEN x, GEN sym)
    1263                 :            : {
    1264                 :      41965 :   GEN p1 = gen_0;
    1265                 :            :   long i;
    1266                 :            : 
    1267         [ -  + ]:      41965 :   if (typ(x) != t_POL) return gmul(x, gel(sym,1));
    1268         [ +  - ]:      41965 :   if (signe(x))
    1269                 :            :   {
    1270                 :      41965 :     sym--;
    1271         [ +  + ]:     721798 :     for (i=lg(x)-1; i>1; i--)
    1272                 :     679833 :       p1 = gadd(p1, gmul(gel(x,i),gel(sym,i)));
    1273                 :            :   }
    1274                 :      41965 :   return p1;
    1275                 :            : }
    1276                 :            : 
    1277                 :            : static GEN
    1278                 :       2191 : get_Tr(GEN mul, GEN x, GEN basden)
    1279                 :            : {
    1280                 :       2191 :   GEN t, bas = gel(basden,1), den = gel(basden,2);
    1281                 :       2191 :   long i, j, n = lg(bas)-1;
    1282                 :       2191 :   GEN T = cgetg(n+1,t_MAT), TW = cgetg(n+1,t_COL), sym = polsym(x, n-1);
    1283                 :            : 
    1284                 :       2191 :   gel(TW,1) = utoipos(n);
    1285         [ +  + ]:      10514 :   for (i=2; i<=n; i++)
    1286                 :            :   {
    1287                 :       8323 :     t = quicktrace(gel(bas,i), sym);
    1288 [ +  + ][ +  + ]:       8323 :     if (den && den[i]) t = diviiexact(t,gel(den,i));
    1289                 :       8323 :     gel(TW,i) = t; /* tr(w[i]) */
    1290                 :            :   }
    1291                 :       2191 :   gel(T,1) = TW;
    1292         [ +  + ]:      10514 :   for (i=2; i<=n; i++)
    1293                 :            :   {
    1294                 :       8323 :     gel(T,i) = cgetg(n+1,t_COL); gcoeff(T,1,i) = gel(TW,i);
    1295         [ +  + ]:      78582 :     for (j=2; j<=i; j++) /* Tr(W[i]W[j]) */
    1296                 :      70259 :       gcoeff(T,i,j) = gcoeff(T,j,i) = ZV_dotproduct(gel(mul,j+(i-1)*n), TW);
    1297                 :            :   }
    1298                 :       2191 :   return T;
    1299                 :            : }
    1300                 :            : 
    1301                 :            : /* return [bas[i]*denom(bas[i]), denom(bas[i])], denom 1 is given as NULL */
    1302                 :            : GEN
    1303                 :       7090 : get_bas_den(GEN bas)
    1304                 :            : {
    1305                 :       7090 :   GEN b,d,den, dbas = leafcopy(bas);
    1306                 :       7090 :   long i, l = lg(bas);
    1307                 :       7090 :   int power = 1;
    1308                 :       7090 :   den = cgetg(l,t_VEC);
    1309         [ +  + ]:      42721 :   for (i=1; i<l; i++)
    1310                 :            :   {
    1311                 :      35631 :     b = Q_remove_denom(gel(bas,i), &d);
    1312                 :      35631 :     gel(dbas,i) = b;
    1313         [ +  + ]:      35631 :     gel(den,i) = d; if (d) power = 0;
    1314                 :            :   }
    1315         [ +  + ]:       7090 :   if (power) den = NULL; /* power basis */
    1316                 :       7090 :   return mkvec2(dbas, den);
    1317                 :            : }
    1318                 :            : 
    1319                 :            : /* Internal: nf partially filled. Require pol; fill zk, invzk, multable */
    1320                 :            : void
    1321                 :       2324 : nf_set_multable(GEN nf, GEN bas, GEN basden)
    1322                 :            : {
    1323                 :       2324 :   GEN T = nf_get_pol(nf), invbas, basM;
    1324                 :       2324 :   long i,j, n = degpol(T);
    1325                 :       2324 :   GEN w, den, mul = cgetg(n*n+1,t_MAT);
    1326                 :            : 
    1327         [ +  + ]:       2324 :   if (typ(bas) == t_MAT)
    1328                 :        133 :   { basM = bas; bas = RgM_to_RgXV(basM, varn(T)); }
    1329                 :            :   else
    1330                 :       2191 :     basM = RgV_to_RgM(bas, n);
    1331                 :       2324 :   gel(nf,7) = bas;
    1332                 :       2324 :   gel(nf,8) = invbas = QM_inv(basM, gen_1);
    1333                 :       2324 :   gel(nf,9) = mul;
    1334                 :            : 
    1335         [ +  + ]:       2324 :   if (!basden) basden = get_bas_den(nf_get_zk(nf)); /*integral basis*/
    1336                 :       2324 :   w   = gel(basden,1);
    1337                 :       2324 :   den = gel(basden,2);
    1338                 :            :   /* i = 1 split for efficiency, assume w[1] = 1 */
    1339         [ +  + ]:      13636 :   for (j=1; j<=n; j++)
    1340                 :      11312 :     gel(mul,j) = gel(mul,1+(j-1)*n) = col_ei(n, j);
    1341         [ +  + ]:      11312 :   for (i=2; i<=n; i++)
    1342         [ +  + ]:      81270 :     for (j=i; j<=n; j++)
    1343                 :            :     {
    1344                 :      72282 :       pari_sp av = avma;
    1345         [ +  + ]:      72282 :       GEN z = (i == j)? ZXQ_sqr(gel(w,i), T): ZXQ_mul(gel(w,i),gel(w,j), T);
    1346                 :      72282 :       z = mulmat_pol(invbas, z); /* integral column */
    1347         [ +  + ]:      72282 :       if (den)
    1348                 :            :       {
    1349                 :      64841 :         GEN d = mul_denom(gel(den,i), gel(den,j));
    1350         [ +  + ]:      64841 :         if (d) z = ZC_Z_divexact(z, d);
    1351                 :            :       }
    1352                 :      72282 :       gel(mul,j+(i-1)*n) = gel(mul,i+(j-1)*n) = gerepileupto(av,z);
    1353                 :            :     }
    1354                 :       2324 : }
    1355                 :            : 
    1356                 :            : /* as get_Tr, mul_table not precomputed */
    1357                 :            : static GEN
    1358                 :       1225 : make_Tr(GEN x, GEN basden)
    1359                 :            : {
    1360                 :       1225 :   long i,j, n = degpol(x);
    1361                 :            :   GEN c, t, d;
    1362                 :       1225 :   GEN T   = cgetg(n+1,t_MAT);
    1363                 :       1225 :   GEN sym = polsym(x, n-1);
    1364                 :       1225 :   GEN w   = gel(basden,1); /* W[i] = w[i]/den[i] */
    1365                 :       1225 :   GEN den = gel(basden,2);
    1366                 :            :   /* assume W[1] = 1, case i = 1 split for efficiency */
    1367                 :       1225 :   c = cgetg(n+1,t_COL); gel(T,1) = c;
    1368                 :       1225 :   gel(c, 1) = utoipos(n);
    1369         [ +  + ]:       5271 :   for (j=2; j<=n; j++)
    1370                 :            :   {
    1371                 :       4046 :     pari_sp av = avma;
    1372                 :       4046 :     t = quicktrace(gel(w,j), sym);
    1373         [ +  + ]:       4046 :     if (den)
    1374                 :            :     {
    1375                 :       3248 :       d = gel(den,j);
    1376         [ +  + ]:       3248 :       if (d) t = diviiexact(t, d);
    1377                 :            :     }
    1378                 :       4046 :     gel(c,j) = gerepileuptoint(av, t);
    1379                 :            :   }
    1380         [ +  + ]:       5271 :   for (i=2; i<=n; i++)
    1381                 :            :   {
    1382                 :       4046 :     c = cgetg(n+1,t_COL); gel(T,i) = c;
    1383         [ +  + ]:      31192 :     for (j=1; j<i ; j++) gel(c,j) = gcoeff(T,i,j);
    1384         [ +  + ]:      31192 :     for (   ; j<=n; j++)
    1385                 :            :     {
    1386                 :      27146 :       pari_sp av = avma;
    1387         [ +  + ]:      27146 :       t = (i == j)? ZXQ_sqr(gel(w,i), x): ZXQ_mul(gel(w,i),gel(w,j), x);
    1388                 :      27146 :       t = quicktrace(t, sym);
    1389         [ +  + ]:      27146 :       if (den)
    1390                 :            :       {
    1391                 :      25984 :         d = mul_denom(gel(den,i),gel(den,j));
    1392         [ +  + ]:      25984 :         if (d) t = diviiexact(t, d);
    1393                 :            :       }
    1394                 :      27146 :       gel(c,j) = gerepileuptoint(av, t); /* Tr (W[i]W[j]) */
    1395                 :            :     }
    1396                 :            :   }
    1397                 :       1225 :   return T;
    1398                 :            : }
    1399                 :            : 
    1400                 :            : /* compute roots so that _absolute_ accuracy of M >= prec [also holds for G] */
    1401                 :            : static void
    1402                 :       5766 : get_roots_for_M(nffp_t *F)
    1403                 :            : {
    1404                 :            :   long n, eBD, PREC;
    1405                 :            : 
    1406         [ +  - ]:       5766 :   if (F->extraprec < 0)
    1407                 :            :   { /* not initialized yet */
    1408                 :            :     double er;
    1409                 :       5766 :     n = degpol(F->x);
    1410                 :       5766 :     eBD = 1 + gexpo(gel(F->basden,1));
    1411         [ +  + ]:       5766 :     er  = F->ro? (1+gexpo(F->ro)): fujiwara_bound(F->x);
    1412         [ +  + ]:       5766 :     if (er < 0) er = 0;
    1413                 :       5766 :     F->extraprec = nbits2extraprec(n*er + eBD + log2(n));
    1414                 :            :   }
    1415                 :            : 
    1416                 :       5766 :   PREC = F->prec + F->extraprec;
    1417                 :            :   /* make sure that default accuracy is the same on 32/64bit */
    1418                 :            : #ifndef LONG_IS_64BIT
    1419         [ +  + ]:        864 :   if (odd(PREC)) PREC += EXTRAPRECWORD;
    1420                 :            : #endif
    1421 [ +  + ][ +  + ]:      10709 :   if (F->ro && gprecision(gel(F->ro,1)) >= PREC) return;
    1422                 :       4943 :   F->ro = get_roots(F->x, F->r1, PREC);
    1423                 :            : }
    1424                 :            : 
    1425                 :            : /* [bas[i]/den[i]]= integer basis. roo = real part of the roots */
    1426                 :            : static void
    1427                 :       5766 : make_M(nffp_t *F, int trunc)
    1428                 :            : {
    1429                 :       5766 :   GEN bas = gel(F->basden,1), den = gel(F->basden,2), ro = F->ro;
    1430                 :            :   GEN m, d, M;
    1431                 :       5766 :   long i, j, l = lg(ro), n = lg(bas);
    1432                 :       5766 :   M = cgetg(n,t_MAT);
    1433                 :       5766 :   gel(M,1) = const_col(l-1, gen_1); /* bas[1] = 1 */
    1434         [ +  + ]:      29658 :   for (j=2; j<n; j++)
    1435                 :            :   {
    1436                 :      23892 :     m = cgetg(l,t_COL); gel(M,j) = m;
    1437         [ +  + ]:     220494 :     for (i=1; i<l; i++) gel(m,i) = poleval(gel(bas,j), gel(ro,i));
    1438                 :            :   }
    1439         [ +  + ]:       5766 :   if (den)
    1440         [ +  + ]:      21385 :     for (j=2; j<n; j++)
    1441                 :            :     {
    1442         [ +  + ]:      18397 :       d = gel(den,j); if (!d) continue;
    1443                 :      14875 :       m = gel(M,j);
    1444         [ +  + ]:     163736 :       for (i=1; i<l; i++) gel(m,i) = gdiv(gel(m,i), d);
    1445                 :            :     }
    1446                 :            : 
    1447 [ +  + ][ +  + ]:       5766 :   if (trunc && gprecision(M) > F->prec)
    1448                 :            :   {
    1449                 :       1657 :     M     = gprec_w(M, F->prec);
    1450                 :       1657 :     F->ro = gprec_w(ro,F->prec);
    1451                 :            :   }
    1452                 :       5766 :   F->M = M;
    1453                 :       5766 : }
    1454                 :            : 
    1455                 :            : /* return G real such that G~ * G = T_2 */
    1456                 :            : static void
    1457                 :       5766 : make_G(nffp_t *F)
    1458                 :            : {
    1459                 :       5766 :   GEN G, M = F->M;
    1460                 :       5766 :   long i, j, k, r1 = F->r1, l = lg(M);
    1461                 :            : 
    1462                 :       5766 :   G = cgetg(l, t_MAT);
    1463         [ +  + ]:      35424 :   for (j=1; j<l; j++)
    1464                 :            :   {
    1465                 :      29658 :     GEN g = cgetg(l, t_COL);
    1466                 :      29658 :     GEN m = gel(M,j);
    1467                 :      29658 :     gel(G,j) = g;
    1468         [ +  + ]:     100830 :     for (k=i=1; i<=r1; i++) g[k++] = m[i];
    1469         [ +  + ]:     173870 :     for (     ; k < l; i++)
    1470                 :            :     {
    1471                 :     144212 :       GEN r = gel(m,i);
    1472         [ +  + ]:     144212 :       if (typ(r) == t_COMPLEX)
    1473                 :            :       {
    1474                 :     133336 :         gel(g,k++) = mpadd(gel(r,1), gel(r,2));
    1475                 :     133336 :         gel(g,k++) = mpsub(gel(r,1), gel(r,2));
    1476                 :            :       }
    1477                 :            :       else
    1478                 :            :       {
    1479                 :      10876 :         gel(g,k++) = r;
    1480                 :      10876 :         gel(g,k++) = r;
    1481                 :            :       }
    1482                 :            :     }
    1483                 :            :   }
    1484                 :       5766 :   F->G = G;
    1485                 :       5766 : }
    1486                 :            : 
    1487                 :            : static void
    1488                 :       5766 : make_M_G(nffp_t *F, int trunc)
    1489                 :            : {
    1490                 :       5766 :   get_roots_for_M(F);
    1491                 :       5766 :   make_M(F, trunc);
    1492                 :       5766 :   make_G(F);
    1493                 :       5766 : }
    1494                 :            : 
    1495                 :            : void
    1496                 :        671 : remake_GM(GEN nf, nffp_t *F, long prec)
    1497                 :            : {
    1498                 :        671 :   F->x  = nf_get_pol(nf);
    1499                 :        671 :   F->ro = NULL;
    1500                 :        671 :   F->r1 = nf_get_r1(nf);
    1501                 :        671 :   F->basden = get_bas_den(nf_get_zk(nf));
    1502                 :        671 :   F->extraprec = -1;
    1503                 :        671 :   F->prec = prec; make_M_G(F, 1);
    1504                 :        671 : }
    1505                 :            : 
    1506                 :            : static void
    1507                 :       5095 : nffp_init(nffp_t *F, nfbasic_t *T, GEN ro, long prec)
    1508                 :            : {
    1509                 :       5095 :   F->x  = T->x;
    1510                 :       5095 :   F->ro = ro;
    1511                 :       5095 :   F->r1 = T->r1;
    1512         [ +  + ]:       5095 :   if (!T->basden) T->basden = get_bas_den(T->bas);
    1513                 :       5095 :   F->basden = T->basden;
    1514                 :       5095 :   F->extraprec = -1;
    1515                 :       5095 :   F->prec = prec;
    1516                 :       5095 : }
    1517                 :            : 
    1518                 :            : static void
    1519                 :       3205 : get_nf_fp_compo(nfbasic_t *T, nffp_t *F, GEN ro, int trunc, long prec)
    1520                 :            : {
    1521                 :       3205 :   nffp_init(F,T,ro,prec);
    1522                 :       3205 :   make_M_G(F, trunc);
    1523                 :       3205 : }
    1524                 :            : 
    1525                 :            : static GEN
    1526                 :       2191 : get_sign(long r1, long n) { return mkvec2s(r1, (n-r1)>>1); }
    1527                 :            : 
    1528                 :            : GEN
    1529                 :       2191 : nfbasic_to_nf(nfbasic_t *T, GEN ro, long prec)
    1530                 :            : {
    1531                 :       2191 :   GEN nf = cgetg(10,t_VEC);
    1532                 :       2191 :   GEN x = T->x, absdK, Tr, D, TI, A, dA, MDI, mat = cgetg(9,t_VEC);
    1533                 :       2191 :   long n = degpol(T->x);
    1534                 :            :   nffp_t F;
    1535                 :       2191 :   get_nf_fp_compo(T, &F, ro, 0, prec);
    1536                 :            : 
    1537                 :       2191 :   gel(nf,1) = T->x;
    1538                 :       2191 :   gel(nf,2) = get_sign(T->r1, n);
    1539                 :       2191 :   gel(nf,3) = T->dK;
    1540                 :       2191 :   gel(nf,4) = T->index;
    1541                 :       2191 :   gel(nf,6) = F.ro;
    1542                 :       2191 :   gel(nf,5) = mat;
    1543                 :            : 
    1544                 :       2191 :   gel(mat,1) = F.M;
    1545                 :       2191 :   gel(mat,2) = F.G;
    1546                 :            : 
    1547                 :       2191 :   nf_set_multable(nf, T->bas, F.basden);
    1548                 :            : 
    1549                 :       2191 :   Tr = get_Tr(gel(nf,9), x, F.basden);
    1550         [ +  + ]:       2191 :   absdK = T->dK; if (signe(absdK) < 0) absdK = negi(absdK);
    1551                 :       2191 :   TI = ZM_inv(Tr, absdK); /* dK T^-1 */
    1552                 :       2191 :   A = Q_primitive_part(TI, &dA);
    1553                 :       2191 :   gel(mat,6) = A; /* primitive part of codifferent, dA its content */
    1554         [ +  + ]:       2191 :   dA = dA? diviiexact(absdK, dA): absdK;
    1555                 :       2191 :   A = ZM_hnfmodid(A, dA);
    1556                 :       2191 :   MDI = idealtwoelt(nf, A);
    1557                 :       2191 :   gel(MDI,2) = zk_scalar_or_multable(nf, gel(MDI,2));
    1558                 :       2191 :   gel(mat,7) = MDI;
    1559         [ +  + ]:       2191 :   if (is_pm1(T->index)) /* principal ideal (x'), whose norm is |dK| */
    1560                 :            :   {
    1561                 :       1365 :     D = zk_scalar_or_multable(nf, ZX_deriv(x));
    1562         [ +  + ]:       1365 :     if (typ(D) == t_MAT) D = ZM_hnfmod(D, absdK);
    1563                 :            :   }
    1564                 :            :   else
    1565                 :        826 :     D = RgM_Rg_mul(idealinv(nf, A), dA);
    1566                 :       2191 :   gel(mat,3) = RM_round_maxrank(F.G);
    1567                 :       2191 :   gel(mat,4) = Tr;
    1568                 :       2191 :   gel(mat,5) = D;
    1569         [ +  + ]:       2191 :   gel(mat,8) = T->dKP? shallowtrans(T->dKP): cgetg(1,t_VEC);
    1570                 :       2191 :   return nf;
    1571                 :            : }
    1572                 :            : 
    1573                 :            : static GEN
    1574                 :         84 : primes_certify(GEN dK, GEN dKP)
    1575                 :            : {
    1576                 :         84 :   pari_sp av = avma;
    1577                 :         84 :   long i, l = lg(dKP);
    1578                 :         84 :   GEN v, D = dK;
    1579                 :         84 :   v = vectrunc_init(l);
    1580         [ +  + ]:        413 :   for (i = 1; i < l; i++)
    1581                 :            :   {
    1582                 :        329 :     GEN p = gel(dKP,i);
    1583         [ +  + ]:        329 :     if (!isprime(p)) vectrunc_append(v, p);
    1584                 :        329 :     (void)Z_pvalrem(D, p, &D);
    1585                 :            :   }
    1586         [ -  + ]:         84 :   if (!is_pm1(D))
    1587                 :            :   {
    1588         [ #  # ]:          0 :     if (signe(D) < 0) D = negi(D);
    1589         [ #  # ]:          0 :     if (!isprime(D)) vectrunc_append(v, D);
    1590                 :            :   }
    1591                 :         84 :   fixlg(v, lg(v)); return gerepilecopy(av, v);
    1592                 :            : }
    1593                 :            : GEN
    1594                 :          7 : nfcertify(GEN nf)
    1595                 :            : {
    1596                 :          7 :   nf = checknf(nf);
    1597                 :          7 :   return primes_certify(nf_get_disc(nf),gmael(nf, 5, 8));
    1598                 :            : }
    1599                 :            : 
    1600                 :            : #if 0 /* used to check benches between HNF nf.zk and LLL-reduced nf.zk */
    1601                 :            : static GEN
    1602                 :            : hnffromLLL(GEN nf)
    1603                 :            : {
    1604                 :            :   GEN d, x;
    1605                 :            :   x = RgV_to_RgM(nf_get_zk(nf), nf_get_degree(nf));
    1606                 :            :   x = Q_remove_denom(x, &d);
    1607                 :            :   if (!d) return x; /* power basis */
    1608                 :            :   return RgM_solve(ZM_hnfmodid(x, d), x);
    1609                 :            : }
    1610                 :            : 
    1611                 :            : static GEN
    1612                 :            : nfbasechange(GEN u, GEN x)
    1613                 :            : {
    1614                 :            :   long i,lx;
    1615                 :            :   GEN y;
    1616                 :            :   switch(typ(x))
    1617                 :            :   {
    1618                 :            :     case t_COL: /* nfelt */
    1619                 :            :       return RgM_RgC_mul(u, x);
    1620                 :            : 
    1621                 :            :     case t_MAT: /* ideal */
    1622                 :            :       y = cgetg_copy(x, &lx);
    1623                 :            :       for (i=1; i<lx; i++) gel(y,i) = RgM_RgC_mul(u, gel(x,i));
    1624                 :            :       break;
    1625                 :            : 
    1626                 :            :     case t_VEC: /* pr */
    1627                 :            :       checkprid(x); y = leafcopy(x);
    1628                 :            :       gel(y,2) = RgM_RgC_mul(u, gel(y,2));
    1629                 :            :       gel(y,5) = RgM_RgC_mul(u, gel(y,5));
    1630                 :            :       break;
    1631                 :            :     default: y = x;
    1632                 :            :   }
    1633                 :            :   return y;
    1634                 :            : }
    1635                 :            : 
    1636                 :            : GEN
    1637                 :            : nffromhnfbasis(GEN nf, GEN x)
    1638                 :            : {
    1639                 :            :   long tx = typ(x);
    1640                 :            :   pari_sp av = avma;
    1641                 :            :   GEN u;
    1642                 :            :   if (!is_vec_t(tx)) return gcopy(x);
    1643                 :            :   nf = checknf(nf);
    1644                 :            :   u = hnffromLLL(nf);
    1645                 :            :   return gerepilecopy(av, nfbasechange(u, x));
    1646                 :            : }
    1647                 :            : 
    1648                 :            : GEN
    1649                 :            : nftohnfbasis(GEN nf, GEN x)
    1650                 :            : {
    1651                 :            :   long tx = typ(x);
    1652                 :            :   pari_sp av = avma;
    1653                 :            :   GEN u;
    1654                 :            :   if (!is_vec_t(tx)) return gcopy(x);
    1655                 :            :   nf = checknf(nf);
    1656                 :            :   u = ZM_inv(hnffromLLL(nf), gen_1);
    1657                 :            :   return gerepilecopy(av, nfbasechange(u, x));
    1658                 :            : }
    1659                 :            : #endif
    1660                 :            : 
    1661                 :            : /* set *pro to roots of T->x */
    1662                 :            : static GEN
    1663                 :       1890 : get_red_G(nfbasic_t *T, GEN *pro)
    1664                 :            : {
    1665                 :       1890 :   GEN G, u, u0 = NULL;
    1666                 :            :   pari_sp av;
    1667                 :       1890 :   long i, prec, n = degpol(T->x);
    1668                 :            :   nffp_t F;
    1669                 :            : 
    1670                 :       1890 :   prec = nbits2prec(n+32);
    1671                 :       1890 :   nffp_init(&F, T, NULL, prec);
    1672                 :       1890 :   av = avma;
    1673                 :       1890 :   for (i=1; ; i++)
    1674                 :            :   {
    1675                 :       1890 :     F.prec = prec; make_M_G(&F, 0); G = F.G;
    1676         [ -  + ]:       1890 :     if (u0) G = RgM_mul(G, u0);
    1677         [ -  + ]:       1890 :     if (DEBUGLEVEL)
    1678                 :          0 :       err_printf("get_red_G: starting LLL, prec = %ld (%ld + %ld)\n",
    1679                 :          0 :                   prec + F.extraprec, prec, F.extraprec);
    1680         [ +  - ]:       1890 :     if ((u = lllfp(G, 0.99, LLL_KEEP_FIRST|LLL_COMPATIBLE)))
    1681                 :            :     {
    1682         [ +  - ]:       1890 :       if (lg(u)-1 == n) break;
    1683                 :            :       /* singular ==> loss of accuracy */
    1684         [ #  # ]:          0 :       if (u0) u0 = gerepileupto(av, RgM_mul(u0,u));
    1685                 :          0 :       else    u0 = gerepilecopy(av, u);
    1686                 :            :     }
    1687                 :          0 :     prec = precdbl(prec) + nbits2extraprec(gexpo(u0));
    1688                 :          0 :     F.ro = NULL;
    1689         [ #  # ]:          0 :     if (DEBUGLEVEL) pari_warn(warnprec,"get_red_G", prec);
    1690                 :          0 :   }
    1691         [ -  + ]:       1890 :   if (u0) u = RgM_mul(u0,u);
    1692                 :       1890 :   *pro = F.ro; return u;
    1693                 :            : }
    1694                 :            : 
    1695                 :            : /* Compute an LLL-reduced basis for the integer basis of nf(T).
    1696                 :            :  * set *pro = roots of x if computed [NULL if not computed] */
    1697                 :            : static void
    1698                 :       3115 : set_LLL_basis(nfbasic_t *T, GEN *pro, double DELTA)
    1699                 :            : {
    1700                 :       3115 :   GEN B = T->bas;
    1701         [ +  + ]:       3115 :   if (T->r1 == degpol(T->x)) {
    1702                 :       1225 :     pari_sp av = avma;
    1703                 :       1225 :     GEN u, basden = T->basden;
    1704         [ +  - ]:       1225 :     if (!basden) basden = get_bas_den(B);
    1705                 :       1225 :     u = ZM_lll(make_Tr(T->x,basden), DELTA,
    1706                 :            :                LLL_GRAM|LLL_KEEP_FIRST|LLL_IM|LLL_COMPATIBLE);
    1707                 :       1225 :     B = gerepileupto(av, RgV_RgM_mul(B, u));
    1708                 :       1225 :     *pro = NULL;
    1709                 :            :   }
    1710                 :            :   else
    1711                 :       1890 :     B = RgV_RgM_mul(B, get_red_G(T, pro));
    1712                 :       3115 :   T->bas = B;
    1713                 :       3115 :   T->basden = get_bas_den(B);
    1714                 :       3115 : }
    1715                 :            : 
    1716                 :            : static int
    1717                 :       2955 : cmp_abs_ZX(GEN x, GEN y) { return gen_cmp_RgX((void*)&absi_cmp, x, y); }
    1718                 :            : /* current best: ZX x of discriminant *dx, is ZX y better than x ?
    1719                 :            :  * (if so update *dx) */
    1720                 :            : static int
    1721                 :       4362 : ZX_is_better(GEN y, GEN x, GEN *dx)
    1722                 :            : {
    1723                 :       4362 :   GEN d = ZX_disc(y);
    1724                 :            :   int cmp;
    1725         [ +  + ]:       4362 :   if (!*dx) *dx = ZX_disc(x);
    1726                 :       4362 :   cmp = absi_cmp(d, *dx);
    1727         [ +  + ]:       4362 :   if (cmp < 0) { *dx = d; return 1; }
    1728         [ +  + ]:       3564 :   if (cmp == 0) return cmp_abs_ZX(y, x) < 0;
    1729                 :       4362 :   return 0;
    1730                 :            : }
    1731                 :            : 
    1732                 :            : static void polredbest_aux(nfbasic_t *T, GEN *pro, GEN *px, GEN *pdx, GEN *pa);
    1733                 :            : /* Seek a simpler, polynomial pol defining the same number field as
    1734                 :            :  * x (assumed to be monic at this point) */
    1735                 :            : static GEN
    1736                 :         63 : nfpolred(nfbasic_t *T, GEN *pro)
    1737                 :            : {
    1738                 :         63 :   GEN x = T->x, dx, b, rev;
    1739                 :         63 :   long n = degpol(x), v = varn(x);
    1740                 :            : 
    1741         [ +  + ]:         63 :   if (n == 1) {
    1742                 :          7 :     T->x = deg1pol_shallow(gen_1, gen_m1, v);
    1743                 :          7 :     *pro = NULL; return pol_1(v);
    1744                 :            :   }
    1745                 :         56 :   polredbest_aux(T, pro, &x, &dx, &b);
    1746         [ +  + ]:         56 :   if (x == T->x) return NULL; /* no improvement */
    1747         [ -  + ]:         42 :   if (DEBUGLEVEL>1) err_printf("xbest = %Ps\n",x);
    1748                 :            : 
    1749                 :            :   /* update T */
    1750                 :         42 :   rev = QXQ_reverse(b, T->x);
    1751                 :         42 :   T->bas = QXV_QXQ_eval(T->bas, rev, x);
    1752                 :         42 :   (void)Z_issquareall(diviiexact(dx,T->dK), &(T->index));
    1753                 :         42 :   T->basden = get_bas_den(T->bas);
    1754                 :         42 :   T->dx = dx;
    1755                 :         42 :   T->x = x;
    1756                 :         42 :   *pro = NULL; /* reset */
    1757                 :         63 :   return rev;
    1758                 :            : }
    1759                 :            : 
    1760                 :            : /* let bas a t_VEC of QX giving a Z-basis of O_K. Return the index of the
    1761                 :            :  * basis. Assume bas[1] is 1 and that the leading coefficient of elements
    1762                 :            :  * of bas are of the form 1/b for a t_INT b */
    1763                 :            : GEN
    1764                 :        154 : get_nfindex(GEN bas)
    1765                 :            : {
    1766                 :        154 :   pari_sp av = avma;
    1767                 :        154 :   long n = lg(bas)-1, i;
    1768                 :            :   GEN D, d, mat;
    1769                 :            : 
    1770                 :        154 :   D = gen_1; /* assume bas[1] = 1 */
    1771         [ +  + ]:       1092 :   for (i = 2; i <= n; i++)
    1772                 :            :   { /* in most cases [e.g after nfbasis] basis is upper triangular! */
    1773                 :        945 :     GEN B = gel(bas,i), lc;
    1774         [ +  + ]:        945 :     if (degpol(B) != i-1) break;
    1775                 :        938 :     lc = gel(B, i+1);
    1776      [ +  +  - ]:        938 :     switch (typ(lc))
    1777                 :            :     {
    1778                 :        364 :       case t_INT: continue;
    1779                 :        574 :       case t_FRAC: lc = gel(lc,2); break;
    1780                 :          0 :       default: pari_err_TYPE("get_nfindex",lc);
    1781                 :            :     }
    1782                 :        574 :     D = mulii(D, lc);
    1783                 :            :   }
    1784         [ +  + ]:        154 :   if (i <= n)
    1785                 :            :   { /* not triangular after all */
    1786                 :          7 :     bas = Q_remove_denom(bas, &d);
    1787         [ -  + ]:          7 :     if (!d) { avma = av; return D; }
    1788                 :          7 :     mat = RgV_to_RgM(bas, n);
    1789                 :          7 :     d = diviiexact(powiu(d, n), ZM_det(mat));
    1790                 :          7 :     D = mulii(D,absi(d));
    1791                 :            :   }
    1792                 :        154 :   return gerepileuptoint(av, D);
    1793                 :            : }
    1794                 :            : 
    1795                 :            : /* Either nf type or ZX or [monic ZX, data], where data is either an integral
    1796                 :            :  * basis (deprecated), or listP data (nfbasis input format) to specify
    1797                 :            :  * a set of primes at with the basis order must be maximal.
    1798                 :            :  * 1) nf type (or unrecognized): return t_VEC
    1799                 :            :  * 2) ZX or [ZX, listP]: return t_POL
    1800                 :            :  * 3) [ZX, order basis]: return 0 (deprecated)
    1801                 :            :  * incorrect: return -1 */
    1802                 :            : static long
    1803                 :       3192 : nf_input_type(GEN x)
    1804                 :            : {
    1805                 :            :   GEN T, V;
    1806                 :            :   long i, d, v;
    1807      [ +  +  - ]:       3192 :   switch(typ(x))
    1808                 :            :   {
    1809                 :       3073 :     case t_POL: return t_POL;
    1810                 :            :     case t_VEC:
    1811         [ +  + ]:        119 :       if (lg(x) != 3) return t_VEC; /* nf or incorrect */
    1812                 :        112 :       T = gel(x,1); V = gel(x,2);
    1813         [ -  + ]:        112 :       if (typ(T) != t_POL) return -1;
    1814      [ +  +  - ]:        112 :       switch(typ(V))
    1815                 :            :       {
    1816                 :         35 :         case t_INT: case t_MAT: return t_POL;
    1817                 :            :         case t_VEC:
    1818         [ +  + ]:         77 :           if (RgV_is_ZV(V)) return t_POL;
    1819                 :         70 :           break;
    1820                 :          0 :         default: return -1;
    1821                 :            :       }
    1822                 :         70 :       d = degpol(T); v = varn(T);
    1823 [ +  - ][ +  - ]:         70 :       if (d<1 || !RgX_is_ZX(T) || !isint1(gel(T,d+2)) || lg(V)-1!=d) return -1;
         [ +  - ][ +  + ]
    1824         [ +  + ]:        637 :       for (i = 1; i <= d; i++)
    1825                 :            :       { /* check integer basis */
    1826                 :        588 :         GEN c = gel(V,i);
    1827      [ +  +  - ]:        588 :         switch(typ(c))
    1828                 :            :         {
    1829                 :         28 :           case t_INT: break;
    1830 [ +  + ][ +  - ]:        560 :           case t_POL: if (varn(c) == v && RgX_is_QX(c) && degpol(c) < d) break;
                 [ +  + ]
    1831                 :            :           /* fall through */
    1832                 :         14 :           default: return -1;
    1833                 :            :         }
    1834                 :            :       }
    1835                 :         49 :       return 0;
    1836                 :            :   }
    1837                 :       3192 :   return t_VEC; /* nf or incorrect */
    1838                 :            : }
    1839                 :            : 
    1840                 :            : static void
    1841                 :       2184 : nfbasic_add_disc(nfbasic_t *T)
    1842                 :            : {
    1843         [ +  + ]:       2184 :   if (!T->index) T->index = get_nfindex(T->bas);
    1844         [ +  + ]:       2184 :   if (!T->dx) T->dx = ZX_disc(T->x);
    1845         [ +  + ]:       2184 :   if (!T->dK) T->dK = diviiexact(T->dx, sqri(T->index));
    1846                 :       2184 : }
    1847                 :            : 
    1848                 :            : static void
    1849                 :       3192 : nfbasic_init(GEN x, long flag, nfbasic_t *T)
    1850                 :            : {
    1851                 :       3192 :   GEN bas, dK, dx, index, unscale = gen_1;
    1852                 :            :   long r1;
    1853                 :            : 
    1854                 :       3192 :   T->dKP = NULL;
    1855   [ +  +  +  + ]:       3192 :   switch (nf_input_type(x))
    1856                 :            :   {
    1857                 :            :     case t_POL:
    1858                 :            :     {
    1859                 :            :       nfmaxord_t S;
    1860                 :       3115 :       nfmaxord(&S, x, flag);
    1861                 :       3108 :       x = S.T;
    1862                 :       3108 :       T->x0 = S.T0;
    1863                 :       3108 :       T->dKP = S.dKP;
    1864                 :       3108 :       dK = S.dK;
    1865                 :       3108 :       index = S.index;
    1866                 :       3108 :       bas = S.basis;
    1867                 :       3108 :       dx = S.dT;
    1868                 :       3108 :       unscale = S.unscale;
    1869                 :       3108 :       r1 = ZX_sturm(x);
    1870                 :            :       break;
    1871                 :            :     }
    1872                 :            :     case t_VEC:
    1873                 :            :     { /* nf, bnf, bnr */
    1874                 :          7 :       GEN nf = checknf(x);
    1875                 :          7 :       x     = nf_get_pol(nf);
    1876                 :          7 :       dK    = nf_get_disc(nf);
    1877                 :          7 :       index = nf_get_index(nf);
    1878                 :          7 :       bas   = nf_get_zk(nf);
    1879                 :          7 :       T->x0 = x;
    1880                 :          7 :       dx = NULL;
    1881                 :          7 :       r1 = nf_get_r1(nf);
    1882                 :          7 :       break;
    1883                 :            :     }
    1884                 :            :     case 0: /* monic integral polynomial + integer basis */
    1885                 :         49 :       bas = gel(x,2); x = gel(x,1);
    1886                 :         49 :       T->x0 = x;
    1887                 :         49 :       index = NULL;
    1888                 :         49 :       dx = NULL;
    1889                 :         49 :       dK = NULL;
    1890                 :         49 :       r1 = ZX_sturm(x);
    1891                 :         49 :       break;
    1892                 :            :     default: /* -1 */
    1893                 :         21 :       pari_err_TYPE("nfbasic_init", x);
    1894                 :       3164 :       return;
    1895                 :            :   }
    1896                 :       3164 :   T->x     = x;
    1897                 :       3164 :   T->unscale = unscale;
    1898                 :       3164 :   T->r1    = r1;
    1899                 :       3164 :   T->dx    = dx;
    1900                 :       3164 :   T->dK    = dK;
    1901                 :       3164 :   T->bas   = bas;
    1902                 :       3164 :   T->basden= NULL;
    1903                 :       3164 :   T->index = index;
    1904                 :            : }
    1905                 :            : 
    1906                 :            : /* Initialize the number field defined by the polynomial x (in variable v)
    1907                 :            :  * flag & nf_RED:     try a polred first.
    1908                 :            :  * flag & nf_ORIG
    1909                 :            :  *    do a polred and return [nfinit(x), Mod(a,red)], where
    1910                 :            :  *    Mod(a,red) = Mod(v,x) (i.e return the base change). */
    1911                 :            : GEN
    1912                 :       2205 : nfinitall(GEN x, long flag, long prec)
    1913                 :            : {
    1914                 :       2205 :   const pari_sp av = avma;
    1915                 :            :   GEN nf, unscale;
    1916                 :            :   nfbasic_t T;
    1917                 :            : 
    1918                 :       2205 :   nfbasic_init(x, flag, &T);
    1919         [ -  + ]:       2184 :   if (!ZX_is_irred(T.x)) pari_err_IRREDPOL("nfinit",x);
    1920 [ +  + ][ +  + ]:       2184 :   if (!equali1(leading_term(T.x0)) && !(flag & nf_RED))
    1921                 :            :   {
    1922                 :         42 :     pari_warn(warner,"non-monic polynomial. Result of the form [nf,c]");
    1923                 :         42 :     flag |= nf_RED | nf_ORIG;
    1924                 :            :   }
    1925                 :       2184 :   unscale = T.unscale;
    1926 [ +  + ][ +  + ]:       2184 :   if (!(flag & nf_RED) && !isint1(unscale))
    1927                 :            :   { /* implies lc(x0) = 1 and L := 1/unscale is integral */
    1928                 :         63 :     long d = degpol(T.x0);
    1929                 :         63 :     GEN L = ginv(unscale); /* x = L^(-deg(x)) x0(L X) */
    1930                 :         63 :     GEN f= powiu(L, (d*(d-1)) >> 1);
    1931                 :         63 :     T.x = T.x0; /* restore original user-supplied x0, unscale data */
    1932                 :         63 :     T.unscale = gen_1;
    1933                 :         63 :     T.dx    = gmul(T.dx, sqri(f));
    1934                 :         63 :     T.bas   = RgXV_unscale(T.bas, unscale);
    1935                 :         63 :     T.index = gmul(T.index, f);
    1936                 :            :   }
    1937                 :       2184 :   nfbasic_add_disc(&T); /* more expensive after set_LLL_basis */
    1938         [ +  + ]:       2184 :   if (flag & nf_RED)
    1939                 :            :   {
    1940                 :            :     GEN ro, rev;
    1941                 :            :     /* lie to polred: more efficient to update *after* modreverse, than to
    1942                 :            :      * unscale in the polred subsystem */
    1943                 :         63 :     T.unscale = gen_1;
    1944                 :         63 :     rev = nfpolred(&T, &ro);
    1945                 :         63 :     nf = nfbasic_to_nf(&T, ro, prec);
    1946         [ +  + ]:         63 :     if (flag & nf_ORIG)
    1947                 :            :     {
    1948         [ +  + ]:         49 :       if (!rev) rev = pol_x(varn(T.x)); /* no improvement */
    1949         [ +  + ]:         49 :       if (!isint1(unscale)) rev = RgX_Rg_div(rev, unscale);
    1950                 :         49 :       nf = mkvec2(nf, mkpolmod(rev, T.x));
    1951                 :            :     }
    1952                 :         63 :     T.unscale = unscale; /* restore */
    1953                 :            :   } else {
    1954                 :       2121 :     GEN ro; set_LLL_basis(&T, &ro, 0.99);
    1955                 :       2121 :     nf = nfbasic_to_nf(&T, ro, prec);
    1956                 :            :   }
    1957                 :       2184 :   return gerepilecopy(av, nf);
    1958                 :            : }
    1959                 :            : 
    1960                 :            : GEN
    1961                 :          0 : nfinitred(GEN x, long prec)  { return nfinitall(x, nf_RED, prec); }
    1962                 :            : GEN
    1963                 :          0 : nfinitred2(GEN x, long prec) { return nfinitall(x, nf_RED|nf_ORIG, prec); }
    1964                 :            : GEN
    1965                 :       1407 : nfinit(GEN x, long prec)     { return nfinitall(x, 0, prec); }
    1966                 :            : 
    1967                 :            : GEN
    1968                 :        658 : nfinit0(GEN x, long flag,long prec)
    1969                 :            : {
    1970   [ +  +  +  - ]:        658 :   switch(flag)
    1971                 :            :   {
    1972                 :            :     case 0:
    1973                 :        637 :     case 1: return nfinitall(x,0,prec);
    1974                 :         14 :     case 2: case 4: return nfinitall(x,nf_RED,prec);
    1975                 :          7 :     case 3: case 5: return nfinitall(x,nf_RED|nf_ORIG,prec);
    1976                 :          0 :     default: pari_err_FLAG("nfinit");
    1977                 :            :   }
    1978                 :        637 :   return NULL; /* not reached */
    1979                 :            : }
    1980                 :            : 
    1981                 :            : /* assume x a bnr/bnf/nf */
    1982                 :            : long
    1983                 :      19403 : nf_get_prec(GEN x)
    1984                 :            : {
    1985                 :      19403 :   GEN nf = checknf(x), ro = nf_get_roots(nf);
    1986         [ +  - ]:      19403 :   return (typ(ro)==t_VEC)? precision(gel(ro,1)): DEFAULTPREC;
    1987                 :            : }
    1988                 :            : 
    1989                 :            : /* assume nf is an nf */
    1990                 :            : GEN
    1991                 :        671 : nfnewprec_shallow(GEN nf, long prec)
    1992                 :            : {
    1993                 :        671 :   GEN NF = leafcopy(nf);
    1994                 :            :   nffp_t F;
    1995                 :        671 :   gel(NF,5) = leafcopy(gel(NF,5));
    1996                 :        671 :   remake_GM(NF, &F, prec);
    1997                 :        671 :   gel(NF,6) = F.ro;
    1998                 :        671 :   gmael(NF,5,1) = F.M;
    1999                 :        671 :   gmael(NF,5,2) = F.G;
    2000                 :        671 :   return NF;
    2001                 :            : }
    2002                 :            : 
    2003                 :            : GEN
    2004                 :         63 : nfnewprec(GEN nf, long prec)
    2005                 :            : {
    2006                 :            :   GEN z;
    2007   [ +  +  +  - ]:         63 :   switch(nftyp(nf))
    2008                 :            :   {
    2009                 :         49 :     default: pari_err_TYPE("nfnewprec", nf);
    2010                 :          7 :     case typ_BNF: z = bnfnewprec(nf,prec); break;
    2011                 :          7 :     case typ_BNR: z = bnrnewprec(nf,prec); break;
    2012                 :            :     case typ_NF: {
    2013                 :          0 :       pari_sp av = avma;
    2014                 :          0 :       z = gerepilecopy(av, nfnewprec_shallow(checknf(nf), prec));
    2015                 :          0 :       break;
    2016                 :            :     }
    2017                 :            :   }
    2018                 :         14 :   return z;
    2019                 :            : }
    2020                 :            : 
    2021                 :            : /********************************************************************/
    2022                 :            : /**                                                                **/
    2023                 :            : /**                           POLRED                               **/
    2024                 :            : /**                                                                **/
    2025                 :            : /********************************************************************/
    2026                 :            : GEN
    2027                 :          0 : embednorm_T2(GEN x, long r1)
    2028                 :            : {
    2029                 :          0 :   pari_sp av = avma;
    2030                 :          0 :   GEN p = RgV_sumpart(x, r1);
    2031                 :          0 :   GEN q = RgV_sumpart2(x,r1+1, lg(x)-1);
    2032         [ #  # ]:          0 :   if (q != gen_0) p = gadd(p, gmul2n(q,1));
    2033         [ #  # ]:          0 :   return avma == av? gcopy(p): gerepileupto(av, p);
    2034                 :            : }
    2035                 :            : 
    2036                 :            : /* simplified version of gnorm for scalar, non-complex inputs, without GC */
    2037                 :            : static GEN
    2038                 :       6636 : real_norm(GEN x)
    2039                 :            : {
    2040   [ -  +  -  - ]:       6636 :   switch(typ(x))
    2041                 :            :   {
    2042                 :          0 :     case t_INT:  return sqri(x);
    2043                 :       6636 :     case t_REAL: return sqrr(x);
    2044                 :          0 :     case t_FRAC: return sqrfrac(x);
    2045                 :            :   }
    2046                 :          0 :   pari_err_TYPE("real_norm", x);
    2047                 :       6636 :   return NULL;
    2048                 :            : }
    2049                 :            : /* simplified version of gnorm, without GC */
    2050                 :            : static GEN
    2051                 :    1201209 : complex_norm(GEN x)
    2052                 :            : {
    2053         [ +  - ]:    1201209 :   return typ(x) == t_COMPLEX? cxnorm(x): real_norm(x);
    2054                 :            : }
    2055                 :            : /* return T2(x), argument r1 needed in case x has components whose type
    2056                 :            :  * is unexpected, e.g. all of them t_INT for embed(gen_1) */
    2057                 :            : GEN
    2058                 :       1373 : embed_T2(GEN x, long r1)
    2059                 :            : {
    2060                 :       1373 :   pari_sp av = avma;
    2061                 :       1373 :   long i, l = lg(x);
    2062                 :       1373 :   GEN c, s = NULL, t = NULL;
    2063         [ -  + ]:       1373 :   if (typ(gel(x,1)) == t_INT) return muliu(gel(x,1), 2*(l-1)-r1);
    2064         [ +  + ]:       8009 :   for (i = 1; i <= r1; i++)
    2065                 :            :   {
    2066                 :       6636 :     c = real_norm(gel(x,i));
    2067         [ +  + ]:       6636 :     s = s? gadd(s, c): c;
    2068                 :            :   }
    2069         [ +  + ]:       5603 :   for (; i < l; i++)
    2070                 :            :   {
    2071                 :       4230 :     c = complex_norm(gel(x,i));
    2072         [ +  + ]:       4230 :     t = t? gadd(t, c): c;
    2073                 :            :   }
    2074 [ +  + ][ +  + ]:       1373 :   if (t) { t = gmul2n(t,1); s = s? gadd(s,t): t; }
    2075                 :       1373 :   return gerepileupto(av, s);
    2076                 :            : }
    2077                 :            : /* return N(x) */
    2078                 :            : GEN
    2079                 :     732141 : embed_norm(GEN x, long r1)
    2080                 :            : {
    2081                 :     732141 :   pari_sp av = avma;
    2082                 :     732141 :   long i, l = lg(x);
    2083                 :     732141 :   GEN c, s = NULL, t = NULL;
    2084         [ +  + ]:     732141 :   if (typ(gel(x,1)) == t_INT) return powiu(gel(x,1), 2*(l-1)-r1);
    2085         [ +  + ]:    2083463 :   for (i = 1; i <= r1; i++)
    2086                 :            :   {
    2087                 :    1354008 :     c = gel(x,i);
    2088         [ +  + ]:    1354008 :     s = s? gmul(s, c): c;
    2089                 :            :   }
    2090         [ +  + ]:    1926434 :   for (; i < l; i++)
    2091                 :            :   {
    2092                 :    1196979 :     c = complex_norm(gel(x,i));
    2093         [ +  + ]:    1196979 :     t = t? gmul(t, c): c;
    2094                 :            :   }
    2095 [ +  + ][ +  + ]:     729455 :   if (t) s = s? gmul(s,t): t;
    2096                 :     732141 :   return gerepileupto(av, s);
    2097                 :            : }
    2098                 :            : 
    2099                 :            : typedef struct {
    2100                 :            :   long r1, v, prec;
    2101                 :            :   GEN ZKembed; /* embeddings of fincke-pohst-reduced Zk basis */
    2102                 :            :   GEN u; /* matrix giving fincke-pohst-reduced Zk basis */
    2103                 :            :   GEN M; /* embeddings of initial (LLL-reduced) Zk basis */
    2104                 :            :   GEN bound; /* T2 norm of the polynomial defining nf */
    2105                 :            :   long expo_best_disc; /* expo(disc(x)), best generator so far */
    2106                 :            : } CG_data;
    2107                 :            : 
    2108                 :            : /* characteristic pol of x (given by embeddings) */
    2109                 :            : static GEN
    2110                 :      23021 : get_pol(CG_data *d, GEN x)
    2111                 :            : {
    2112                 :            :   long e;
    2113                 :      23021 :   GEN g = grndtoi(roots_to_pol_r1(x, d->v, d->r1), &e);
    2114         [ +  + ]:      23021 :   return (e > -5)? NULL: g;
    2115                 :            : }
    2116                 :            : 
    2117                 :            : /* characteristic pol of x (given as vector on (w_i)) */
    2118                 :            : static GEN
    2119                 :      10313 : get_polchar(CG_data *d, GEN x)
    2120                 :      10313 : { return get_pol(d, RgM_RgC_mul(d->ZKembed,x)); }
    2121                 :            : 
    2122                 :            : /* Choose a canonical polynomial in the pair { z(X), (+/-)z(-X) }.
    2123                 :            :  * z a ZX with lc(z) > 0. We want to keep that property, while
    2124                 :            :  * ensuring that the leading coeff of the odd (resp. even) part of z is < 0
    2125                 :            :  * if deg z is even (resp. odd).
    2126                 :            :  * Either leave z alone (return 1) or set z <-- (-1)^deg(z) z(-X). In place. */
    2127                 :            : static int
    2128                 :      10257 : ZX_canon_neg(GEN z)
    2129                 :            : {
    2130                 :            :   long i,s;
    2131                 :            : 
    2132         [ +  + ]:      81629 :   for (i = lg(z)-2; i >= 2; i -= 2)
    2133                 :            :   { /* examine the odd (resp. even) part of z if deg(z) even (resp. odd). */
    2134                 :      76414 :     s = signe(gel(z,i));
    2135         [ +  + ]:      76414 :     if (!s) continue;
    2136                 :            :     /* non trivial */
    2137         [ +  + ]:       5042 :     if (s < 0) break; /* the condition is already satisfied */
    2138                 :            : 
    2139         [ +  + ]:       8619 :     for (; i>=2; i-=2) gel(z,i) = negi(gel(z,i));
    2140                 :       1695 :     return 1;
    2141                 :            :   }
    2142                 :      10257 :   return 0;
    2143                 :            : }
    2144                 :            : /* return a defining polynomial for Q(alpha), v = embeddings of alpha.
    2145                 :            :  * Return NULL on failure: discriminant too large or non primitive */
    2146                 :            : static GEN
    2147                 :      17773 : try_polmin(CG_data *d, nfbasic_t *T, GEN v, long flag, GEN *ai)
    2148                 :            : {
    2149                 :      17773 :   const long best = flag & nf_ABSOLUTE;
    2150                 :            :   long ed;
    2151                 :      17773 :   pari_sp av = avma;
    2152                 :            :   GEN g;
    2153         [ +  + ]:      17773 :   if (best)
    2154                 :            :   {
    2155                 :      16898 :     ed = expo(embed_disc(v, d->r1, LOWDEFAULTPREC));
    2156         [ +  + ]:      16898 :     avma = av; if (d->expo_best_disc < ed) return NULL;
    2157                 :            :   }
    2158                 :            :   else
    2159                 :        875 :     ed = 0;
    2160                 :      10472 :   g = get_pol(d, v);
    2161                 :            :   /* accuracy too low, compute algebraically */
    2162         [ +  + ]:      10472 :   if (!g) { avma = av; g = ZXQ_charpoly(*ai, T->x, varn(T->x)); }
    2163                 :      10472 :   (void)ZX_gcd_all(g, ZX_deriv(g), &g);
    2164 [ +  + ][ +  + ]:      10472 :   if (best && degpol(g) != degpol(T->x)) { avma = av; return NULL; }
    2165                 :       4459 :   g = gerepilecopy(av, g);
    2166                 :       4459 :   d->expo_best_disc = ed;
    2167         [ +  + ]:       4459 :   if (flag & nf_ORIG)
    2168                 :            :   {
    2169         [ +  + ]:        896 :     if (ZX_canon_neg(g)) *ai = RgX_neg(*ai);
    2170         [ +  + ]:        896 :     if (!isint1(T->unscale)) *ai = RgX_unscale(*ai, T->unscale);
    2171                 :            :   }
    2172                 :            :   else
    2173                 :       3563 :     (void)ZX_canon_neg(g);
    2174         [ -  + ]:       4459 :   if (DEBUGLEVEL>3) err_printf("polred: generator %Ps\n", g);
    2175                 :      17773 :   return g;
    2176                 :            : }
    2177                 :            : 
    2178                 :            : /* does x generate the correct field ? */
    2179                 :            : static GEN
    2180                 :      10313 : chk_gen(void *data, GEN x)
    2181                 :            : {
    2182                 :      10313 :   pari_sp av = avma, av1;
    2183                 :      10313 :   GEN h, g = get_polchar((CG_data*)data,x);
    2184         [ -  + ]:      10313 :   if (!g) pari_err_PREC("chk_gen");
    2185                 :      10313 :   av1 = avma;
    2186                 :      10313 :   h = ZX_gcd(g, ZX_deriv(g));
    2187         [ +  + ]:      10313 :   if (degpol(h)) { avma = av; return NULL; }
    2188         [ -  + ]:       5833 :   if (DEBUGLEVEL>3) err_printf("  generator: %Ps\n",g);
    2189                 :      10313 :   avma = av1; return gerepileupto(av, g);
    2190                 :            : }
    2191                 :            : 
    2192                 :            : static long
    2193                 :       1274 : chk_gen_prec(long N, long bit)
    2194                 :       1274 : { return nbits2prec(10 + (long)log2((double)N) + bit); }
    2195                 :            : 
    2196                 :            : /* Remove duplicate polynomials in P, updating A (same indices), in place.
    2197                 :            :  * Among elements having the same characteristic pol, choose the smallest
    2198                 :            :  * according to ZV_abscmp */
    2199                 :            : static void
    2200                 :        280 : remove_duplicates(GEN P, GEN A)
    2201                 :            : {
    2202                 :        280 :   long k, i, l = lg(P);
    2203                 :        280 :   pari_sp av = avma;
    2204                 :            :   GEN x, a;
    2205                 :            : 
    2206         [ -  + ]:        560 :   if (l < 2) return;
    2207                 :        280 :   (void)sort_factor_pol(mkmat2(P, A), cmpii);
    2208                 :        280 :   x = gel(P,1); a = gel(A,1);
    2209         [ +  + ]:       5798 :   for  (k=1,i=2; i<l; i++)
    2210         [ +  + ]:       5518 :     if (ZX_equal(gel(P,i), x))
    2211                 :            :     {
    2212         [ +  + ]:       3326 :       if (ZV_abscmp(gel(A,i), a) < 0) a = gel(A,i);
    2213                 :            :     }
    2214                 :            :     else
    2215                 :            :     {
    2216                 :       2192 :       gel(A,k) = a;
    2217                 :       2192 :       gel(P,k) = x;
    2218                 :       2192 :       k++;
    2219                 :       2192 :       x = gel(P,i); a = gel(A,i);
    2220                 :            :     }
    2221                 :        280 :   gel(A,k) = a;
    2222                 :        280 :   gel(P,k) = x;
    2223                 :        280 :   l = k+1; setlg(A,l); setlg(P,l);
    2224                 :        280 :   avma = av;
    2225                 :            : }
    2226                 :            : 
    2227                 :            : static long
    2228                 :        994 : polred_init(nfbasic_t *T, nffp_t *F, CG_data *d)
    2229                 :            : {
    2230                 :        994 :   long e, prec, n = degpol(T->x);
    2231                 :            :   double log2rho;
    2232                 :            :   GEN ro;
    2233                 :        994 :   set_LLL_basis(T, &ro, 0.9999);
    2234                 :            :   /* || polchar ||_oo < 2^e ~ 2 (n * rho)^n, rho = max modulus of root */
    2235         [ +  + ]:        994 :   log2rho = ro ? (double)gexpo(ro): fujiwara_bound(T->x);
    2236                 :        994 :   e = n * (long)(log2rho + log2((double)n)) + 1;
    2237         [ -  + ]:        994 :   if (e < 0) e = 0; /* can occur if n = 1 */
    2238                 :        994 :   prec = chk_gen_prec(n, e);
    2239                 :        994 :   get_nf_fp_compo(T, F, ro, 1, prec);
    2240                 :        994 :   d->v = varn(T->x);
    2241                 :        994 :   d->expo_best_disc = -1;
    2242                 :        994 :   d->ZKembed = NULL;
    2243                 :        994 :   d->M = NULL;
    2244                 :        994 :   d->u = NULL;
    2245                 :        994 :   d->r1= T->r1; return prec;
    2246                 :            : }
    2247                 :            : static GEN
    2248                 :        287 : findmindisc(GEN y, GEN *pa)
    2249                 :            : {
    2250                 :        287 :   GEN a = *pa, x = gel(y,1), b = gel(a,1), dx = NULL;
    2251                 :        287 :   long i, l = lg(y);
    2252         [ +  + ]:        358 :   for (i = 2; i < l; i++)
    2253                 :            :   {
    2254                 :         71 :     GEN yi = gel(y,i);
    2255         [ +  + ]:         71 :     if (ZX_is_better(yi,x,&dx)) { x = yi; b = gel(a,i); }
    2256                 :            :   }
    2257                 :        287 :   *pa = b; return x;
    2258                 :            : }
    2259                 :            : /* filter [y,b] from polred_aux: keep a single polynomial of degree n in y
    2260                 :            :  * [ the best wrt discriminant ordering ], but keep all non-primitive
    2261                 :            :  * polynomials */
    2262                 :            : static void
    2263                 :        714 : filter(GEN y, GEN b, long n)
    2264                 :            : {
    2265                 :            :   GEN x, a, dx;
    2266                 :        714 :   long i, k = 1, l = lg(y);
    2267                 :        714 :   a = x = dx = NULL;
    2268         [ +  + ]:       5236 :   for (i = 1; i < l; i++)
    2269                 :            :   {
    2270                 :       4522 :     GEN yi = gel(y,i), ai = gel(b,i);
    2271         [ +  + ]:       4522 :     if (degpol(yi) == n)
    2272                 :            :     {
    2273                 :       4354 :       pari_sp av = avma;
    2274 [ +  + ][ +  + ]:       4354 :       if (dx && !ZX_is_better(yi,x,&dx)) { avma = av; continue; }
    2275         [ +  + ]:       1015 :       if (!dx) dx = ZX_disc(yi);
    2276                 :       1015 :       x = yi; a = ai; continue;
    2277                 :            :     }
    2278                 :        168 :     gel(y,k) = yi;
    2279                 :        168 :     gel(b,k) = ai; k++;
    2280                 :            :   }
    2281         [ +  - ]:        714 :   if (dx)
    2282                 :            :   {
    2283                 :        714 :     gel(y,k) = x;
    2284                 :        714 :     gel(b,k) = a; k++;
    2285                 :            :   }
    2286                 :        714 :   setlg(y, k);
    2287                 :        714 :   setlg(b, k);
    2288                 :        714 : }
    2289                 :            : 
    2290                 :            : static GEN
    2291                 :        742 : polred_aux(nfbasic_t *T, GEN *pro, long flag)
    2292                 :            : { /* only keep polynomials of max degree and best discriminant */
    2293                 :        742 :   const long best = flag & nf_ABSOLUTE;
    2294                 :        742 :   const long orig = flag & nf_ORIG;
    2295                 :        742 :   GEN M, b, y, x = T->x;
    2296                 :        742 :   long maxi, i, j, k, v = varn(x), n = lg(T->bas)-1;
    2297                 :            :   nffp_t F;
    2298                 :            :   CG_data d;
    2299                 :            : 
    2300         [ +  + ]:        742 :   if (n == 1)
    2301                 :            :   {
    2302         [ +  + ]:         28 :     if (!best)
    2303                 :            :     {
    2304                 :         14 :       GEN ch = deg1pol_shallow(gen_1, gen_m1, v);
    2305         [ -  + ]:         14 :       return orig? mkmat2(mkcol(ch),mkcol(gen_1)): mkvec(ch);
    2306                 :            :     }
    2307                 :            :     else
    2308         [ +  + ]:         14 :       return orig? trivial_fact(): cgetg(1,t_VEC);
    2309                 :            :   }
    2310                 :            : 
    2311                 :        714 :   (void)polred_init(T, &F, &d);
    2312                 :        714 :   *pro = F.ro;
    2313                 :        714 :   M = F.M;
    2314         [ +  + ]:        714 :   if (best)
    2315                 :            :   {
    2316         [ -  + ]:        651 :     if (!T->dx) T->dx = ZX_disc(T->x);
    2317                 :        651 :     d.expo_best_disc = expi(T->dx);
    2318                 :            :   }
    2319                 :            : 
    2320                 :            :   /* n + 2 sum_{1 <= i <= n} n-i = n + n(n-1) = n*n */
    2321                 :        714 :   y = cgetg(n*n + 1, t_VEC);
    2322                 :        714 :   b = cgetg(n*n + 1, t_COL);
    2323                 :        714 :   k = 1;
    2324         [ +  + ]:        714 :   if (!best)
    2325                 :            :   {
    2326                 :         63 :     GEN ch = deg1pol_shallow(gen_1, gen_m1, v);
    2327                 :         63 :     gel(y,1) = ch; gel(b,1) = gen_1; k++;
    2328                 :            :   }
    2329         [ +  + ]:       3815 :   for (i = 2; i <= n; i++)
    2330                 :            :   {
    2331                 :            :     GEN ch, ai;
    2332                 :       3101 :     ai = gel(T->bas,i);
    2333                 :       3101 :     ch = try_polmin(&d, T, gel(M,i), flag, &ai);
    2334         [ +  + ]:       3101 :     if (ch) { gel(y,k) = ch; gel(b,k) = ai; k++; }
    2335                 :            :   }
    2336                 :        714 :   maxi = minss(n, 3);
    2337         [ +  + ]:       2681 :   for (i = 1; i <= maxi; i++)
    2338         [ +  + ]:       9303 :     for (j = i+1; j <= n; j++)
    2339                 :            :     {
    2340                 :            :       GEN ch, ai, v;
    2341                 :       7336 :       ai = gadd(gel(T->bas,i), gel(T->bas,j));
    2342                 :       7336 :       v = RgV_add(gel(M,i), gel(M,j));
    2343                 :            :       /* defining polynomial for Q(w_i+w_j) */
    2344                 :       7336 :       ch = try_polmin(&d, T, v, flag, &ai);
    2345         [ +  + ]:       7336 :       if (ch) { gel(y,k) = ch; gel(b,k) = ai; k++; }
    2346                 :            : 
    2347                 :       7336 :       ai = gsub(gel(T->bas,i), gel(T->bas,j));
    2348                 :       7336 :       v = RgV_sub(gel(M,i), gel(M,j));
    2349                 :            :       /* defining polynomial for Q(w_i-w_j) */
    2350                 :       7336 :       ch = try_polmin(&d, T, v, flag, &ai);
    2351         [ +  + ]:       7336 :       if (ch) { gel(y,k) = ch; gel(b,k) = ai; k++; }
    2352                 :            :     }
    2353                 :        714 :   setlg(y, k);
    2354                 :        714 :   setlg(b, k); filter(y, b, n);
    2355         [ +  + ]:        714 :   if (!orig) return gen_sort_uniq(y, (void*)cmpii, &gen_cmp_RgX);
    2356                 :        126 :   (void)sort_factor_pol(mkmat2(y, b), cmpii);
    2357                 :        742 :   settyp(y, t_COL); return mkmat2(b, y);
    2358                 :            : }
    2359                 :            : 
    2360                 :            : static GEN
    2361                 :         84 : Polred(GEN x, long flag, GEN fa)
    2362                 :            : {
    2363                 :         84 :   pari_sp av = avma;
    2364                 :            :   GEN ro;
    2365         [ +  + ]:         84 :   nfbasic_t T; nfbasic_init(fa? mkvec2(x,fa): x, flag & nf_PARTIALFACT, &T);
    2366                 :         77 :   return gerepilecopy(av, polred_aux(&T, &ro, flag));
    2367                 :            : }
    2368                 :            : 
    2369                 :            : /* finds "best" polynomial in polred_aux list, defaulting to T->x if none of
    2370                 :            :  * them is primitive. *px is the ZX, characteristic polynomial of Mod(*pb,T->x),
    2371                 :            :  * *pdx its discriminant. Set *pro = polroots(T->x) [ NOT *px ]. */
    2372                 :            : static void
    2373                 :        665 : polredbest_aux(nfbasic_t *T, GEN *pro, GEN *px, GEN *pdx, GEN *pb)
    2374                 :            : {
    2375                 :        665 :   GEN y, x = T->x; /* default value */
    2376                 :            :   long i, l;
    2377         [ +  + ]:        665 :   y = polred_aux(T, pro, pb? nf_ORIG|nf_ABSOLUTE: nf_ABSOLUTE);
    2378                 :        665 :   *pdx = T->dx;
    2379         [ +  + ]:        665 :   if (pb)
    2380                 :            :   {
    2381                 :        112 :     GEN a, b = deg1pol_shallow(T->unscale, gen_0, varn(x));
    2382                 :        112 :     a = gel(y,1); l = lg(a);
    2383                 :        112 :     y = gel(y,2);
    2384         [ +  + ]:        217 :     for (i=1; i<l; i++)
    2385                 :            :     {
    2386                 :        105 :       GEN yi = gel(y,i);
    2387                 :        105 :       pari_sp av = avma;
    2388         [ +  + ]:        105 :       if (ZX_is_better(yi,x,pdx)) { x = yi; b = gel(a,i); } else avma = av;
    2389                 :            :     }
    2390                 :        112 :     *pb = b;
    2391                 :            :   }
    2392                 :            :   else
    2393                 :            :   {
    2394                 :        553 :     l = lg(y);
    2395         [ +  + ]:       1099 :     for (i=1; i<l; i++)
    2396                 :            :     {
    2397                 :        546 :       GEN yi = gel(y,i);
    2398                 :        546 :       pari_sp av = avma;
    2399         [ +  + ]:        546 :       if (ZX_is_better(yi,x,pdx)) x = yi; else avma = av;
    2400                 :            :     }
    2401                 :            :   }
    2402         [ -  + ]:        665 :   if (!*pdx) *pdx = ZX_disc(x);
    2403                 :        665 :   *px = x;
    2404                 :        665 : }
    2405                 :            : GEN
    2406                 :        581 : polredbest(GEN T0, long flag)
    2407                 :            : {
    2408                 :        581 :   pari_sp av = avma;
    2409                 :            :   GEN T, dT, ro, a;
    2410                 :            :   nfbasic_t S;
    2411 [ +  - ][ -  + ]:        581 :   if (flag < 0 || flag > 1) pari_err_FLAG("polredbest");
    2412                 :        581 :   T = T0; nfbasic_init(T, nf_PARTIALFACT, &S);
    2413         [ +  + ]:        581 :   polredbest_aux(&S, &ro, &T, &dT, flag? &a: NULL);
    2414         [ +  + ]:        581 :   if (flag)
    2415                 :            :   { /* charpoly(Mod(a,T0)) = T */
    2416                 :            :     GEN b;
    2417         [ -  + ]:         28 :     if (T0 == T)
    2418                 :          0 :       b = pol_x(varn(T)); /* no improvement */
    2419                 :            :     else
    2420                 :         28 :       b = QXQ_reverse(a, T0); /* charpoly(Mod(b,T)) = S.x */
    2421         [ +  + ]:         28 :     b = (degpol(T) == 1)? gmodulo(b, T): mkpolmod(b,T);
    2422                 :         28 :     T = mkvec2(T, b);
    2423                 :            :   }
    2424                 :        581 :   return gerepilecopy(av, T);
    2425                 :            : }
    2426                 :            : /* DEPRECATED: backward compatibility */
    2427                 :            : GEN
    2428                 :         70 : polred0(GEN x, long flag, GEN fa)
    2429                 :            : {
    2430                 :         70 :   long fl = 0;
    2431         [ +  + ]:         70 :   if (flag & 1) fl |= nf_PARTIALFACT;
    2432         [ +  + ]:         70 :   if (flag & 2) fl |= nf_ORIG;
    2433                 :         70 :   return Polred(x, fl, fa);
    2434                 :            : }
    2435                 :            : 
    2436                 :            : GEN
    2437                 :         21 : polredord(GEN x)
    2438                 :            : {
    2439                 :         21 :   pari_sp av = avma;
    2440                 :            :   GEN v, lt;
    2441                 :            :   long i, n, vx;
    2442                 :            : 
    2443         [ -  + ]:         21 :   if (typ(x) != t_POL) pari_err_TYPE("polredord",x);
    2444                 :         21 :   x = Q_primpart(x); RgX_check_ZX(x,"polredord");
    2445         [ -  + ]:         21 :   n = degpol(x); if (n <= 0) pari_err_CONSTPOL("polredord");
    2446         [ +  + ]:         21 :   if (n == 1) return gerepilecopy(av, mkvec(x));
    2447                 :         14 :   lt = leading_term(x); vx = varn(x);
    2448         [ +  + ]:         14 :   if (is_pm1(lt))
    2449                 :            :   {
    2450         [ -  + ]:          7 :     if (signe(lt) < 0) x = ZX_neg(x);
    2451                 :          7 :     v = pol_x_powers(n, vx);
    2452                 :            :   }
    2453                 :            :   else
    2454                 :            :   { GEN L;
    2455                 :            :     /* basis for Dedekind order */
    2456                 :          7 :     v = cgetg(n+1, t_VEC);
    2457                 :          7 :     gel(v,1) = scalarpol_shallow(lt, vx);
    2458         [ +  + ]:         14 :     for (i = 2; i <= n; i++)
    2459                 :          7 :       gel(v,i) = RgX_Rg_add(RgX_mulXn(gel(v,i-1), 1), gel(x,n+3-i));
    2460                 :          7 :     gel(v,1) = pol_1(vx);
    2461                 :          7 :     x = ZX_Q_normalize(x, &L);
    2462                 :          7 :     v = gsubst(v, vx, monomial(ginv(L),1,vx));
    2463         [ +  + ]:         14 :     for (i=2; i <= n; i++)
    2464         [ +  - ]:          7 :       if (Q_denom(gel(v,i)) == gen_1) gel(v,i) = monomial(gen_1, i-1, vx);
    2465                 :            :   }
    2466                 :         21 :   return gerepileupto(av, polred(mkvec2(x, v)));
    2467                 :            : }
    2468                 :            : 
    2469                 :            : GEN
    2470                 :         14 : polred(GEN x) { return Polred(x, 0, NULL); }
    2471                 :            : GEN
    2472                 :          0 : smallpolred(GEN x) { return Polred(x, nf_PARTIALFACT, NULL); }
    2473                 :            : GEN
    2474                 :          0 : factoredpolred(GEN x, GEN fa) { return Polred(x, 0, fa); }
    2475                 :            : GEN
    2476                 :          0 : polred2(GEN x) { return Polred(x, nf_ORIG, NULL); }
    2477                 :            : GEN
    2478                 :          0 : smallpolred2(GEN x) { return Polred(x, nf_PARTIALFACT|nf_ORIG, NULL); }
    2479                 :            : GEN
    2480                 :          0 : factoredpolred2(GEN x, GEN fa) { return Polred(x, nf_PARTIALFACT, fa); }
    2481                 :            : 
    2482                 :            : /********************************************************************/
    2483                 :            : /**                                                                **/
    2484                 :            : /**                           POLREDABS                            **/
    2485                 :            : /**                                                                **/
    2486                 :            : /********************************************************************/
    2487                 :            : /* set V[k] := matrix of multiplication by nk.zk[k] */
    2488                 :            : static GEN
    2489                 :       1162 : set_mulid(GEN V, GEN M, GEN Mi, long r1, long r2, long N, long k)
    2490                 :            : {
    2491                 :       1162 :   GEN v, Mk = cgetg(N+1, t_MAT);
    2492                 :            :   long i, e;
    2493         [ +  + ]:       4158 :   for (i = 1; i < k; i++) gel(Mk,i) = gmael(V, i, k);
    2494         [ +  + ]:      11508 :   for (     ; i <=N; i++)
    2495                 :            :   {
    2496                 :      10346 :     v = vecmul(gel(M,k), gel(M,i));
    2497                 :      10346 :     v = RgM_RgC_mul(Mi, split_realimag(v, r1, r2));
    2498                 :      10346 :     gel(Mk,i) = grndtoi(v, &e);
    2499         [ -  + ]:      10346 :     if (e > -5) return NULL;
    2500                 :            :   }
    2501                 :       1162 :   gel(V,k) = Mk; return Mk;
    2502                 :            : }
    2503                 :            : 
    2504                 :            : static GEN
    2505                 :       1308 : ZM_image_shallow(GEN M, long *pr)
    2506                 :            : {
    2507                 :            :   long j, k, r;
    2508                 :       1308 :   GEN y, d = ZM_pivots(M, &k);
    2509                 :       1308 :   r = lg(M)-1 - k;
    2510                 :       1308 :   y = cgetg(r+1,t_MAT);
    2511         [ +  + ]:       9955 :   for (j=k=1; j<=r; k++)
    2512         [ +  + ]:       8647 :     if (d[k]) gel(y,j++) = gel(M,k);
    2513                 :       1308 :   *pr = r; return y;
    2514                 :            : }
    2515                 :            : 
    2516                 :            : /* U = base change matrix, R = Cholesky form of the quadratic form [matrix
    2517                 :            :  * Q from algo 2.7.6] */
    2518                 :            : static GEN
    2519                 :        288 : chk_gen_init(FP_chk_fun *chk, GEN R, GEN U)
    2520                 :            : {
    2521                 :        288 :   CG_data *d = (CG_data*)chk->data;
    2522                 :            :   GEN P, V, D, inv, bound, S, M;
    2523                 :        288 :   long N = lg(U)-1, r1 = d->r1, r2 = (N-r1)>>1;
    2524                 :        288 :   long i, j, prec, firstprim = 0, skipfirst = 0;
    2525                 :            :   pari_sp av;
    2526                 :            : 
    2527                 :        288 :   d->u = U;
    2528                 :        288 :   d->ZKembed = M = RgM_mul(d->M, U);
    2529                 :            : 
    2530                 :        288 :   av = avma; bound = d->bound;
    2531                 :        288 :   D = cgetg(N+1, t_VECSMALL);
    2532         [ +  + ]:       2481 :   for (i = 1; i <= N; i++)
    2533                 :            :   {
    2534                 :       2201 :     pari_sp av2 = avma;
    2535                 :       2201 :     P = get_pol(d, gel(M,i));
    2536         [ +  + ]:       2201 :     if (!P) pari_err_PREC("chk_gen_init");
    2537                 :       2193 :     (void)ZX_gcd_all(P, ZX_deriv(P), &P);
    2538                 :       2193 :     P = gerepilecopy(av2, P);
    2539                 :       2193 :     D[i] = degpol(P);
    2540         [ +  + ]:       2193 :     if (D[i] == N)
    2541                 :            :     { /* primitive element */
    2542                 :        883 :       GEN B = embed_T2(gel(M,i), r1);
    2543         [ +  + ]:        883 :       if (!firstprim) firstprim = i; /* index of first primitive element */
    2544         [ -  + ]:        883 :       if (DEBUGLEVEL>2) err_printf("chk_gen_init: generator %Ps\n",P);
    2545         [ +  + ]:        883 :       if (gcmp(B,bound) < 0) bound = gerepileuptoleaf(av2, B);
    2546                 :            :     }
    2547                 :            :     else
    2548                 :            :     {
    2549         [ -  + ]:       1310 :       if (DEBUGLEVEL>2) err_printf("chk_gen_init: subfield %Ps\n",P);
    2550         [ +  + ]:       1310 :       if (firstprim)
    2551                 :            :       { /* cycle basis vectors so that primitive elements come last */
    2552                 :        125 :         GEN u = d->u, e = M;
    2553                 :        125 :         GEN te = gel(e,i), tu = gel(u,i), tR = gel(R,i);
    2554                 :        125 :         long tS = D[i];
    2555         [ +  + ]:        349 :         for (j = i; j > firstprim; j--)
    2556                 :            :         {
    2557                 :        224 :           u[j] = u[j-1];
    2558                 :        224 :           e[j] = e[j-1];
    2559                 :        224 :           R[j] = R[j-1];
    2560                 :        224 :           D[j] = D[j-1];
    2561                 :            :         }
    2562                 :        125 :         gel(u,firstprim) = tu;
    2563                 :        125 :         gel(e,firstprim) = te;
    2564                 :        125 :         gel(R,firstprim) = tR;
    2565                 :        125 :         D[firstprim] = tS; firstprim++;
    2566                 :            :       }
    2567                 :            :     }
    2568                 :            :   }
    2569         [ +  + ]:        280 :   if (!firstprim)
    2570                 :            :   { /* try (a little) to find primitive elements to improve bound */
    2571                 :         21 :     GEN x = cgetg(N+1, t_VECSMALL);
    2572         [ -  + ]:         21 :     if (DEBUGLEVEL>1)
    2573                 :          0 :       err_printf("chk_gen_init: difficult field, trying random elements\n");
    2574         [ +  + ]:        231 :     for (i = 0; i < 10; i++)
    2575                 :            :     {
    2576                 :            :       GEN e, B;
    2577         [ +  + ]:       2310 :       for (j = 1; j <= N; j++) x[j] = (long)random_Fl(7) - 3;
    2578                 :        210 :       e = RgM_zc_mul(M, x);
    2579                 :        210 :       B = embed_T2(e, r1);
    2580         [ +  + ]:        210 :       if (gcmp(B,bound) >= 0) continue;
    2581         [ -  + ]:         35 :       P = get_pol(d, e); if (!P) pari_err_PREC( "chk_gen_init");
    2582         [ -  + ]:         35 :       if (!ZX_is_squarefree(P)) continue;
    2583         [ -  + ]:         35 :       if (DEBUGLEVEL>2) err_printf("chk_gen_init: generator %Ps\n",P);
    2584                 :         35 :       bound = B ;
    2585                 :            :     }
    2586                 :            :   }
    2587                 :            : 
    2588         [ +  - ]:        280 :   if (firstprim != 1)
    2589                 :            :   {
    2590                 :        280 :     inv = ginv( split_realimag(M, r1, r2) ); /*TODO: use QR?*/
    2591                 :        280 :     V = gel(inv,1);
    2592         [ +  + ]:       1407 :     for (i = 2; i <= r1+r2; i++) V = gadd(V, gel(inv,i));
    2593                 :            :     /* V corresponds to 1_Z */
    2594                 :        280 :     V = grndtoi(V, &j);
    2595         [ -  + ]:        280 :     if (j > -5) pari_err_BUG("precision too low in chk_gen_init");
    2596                 :        280 :     S = mkmat(V); /* 1 */
    2597                 :            : 
    2598                 :        280 :     V = cgetg(N+1, t_VEC);
    2599         [ +  - ]:       1267 :     for (i = 1; i <= N; i++,skipfirst++)
    2600                 :            :     { /* S = Q-basis of subfield generated by nf.zk[1..i-1] */
    2601                 :            :       GEN Mx, M2;
    2602                 :       1267 :       long j, k, h, rkM, dP = D[i];
    2603                 :            : 
    2604         [ +  + ]:       1267 :       if (dP == N) break; /* primitive */
    2605                 :       1162 :       Mx = set_mulid(V, M, inv, r1, r2, N, i);
    2606         [ -  + ]:       1162 :       if (!Mx) break; /* prec. problem. Stop */
    2607         [ +  + ]:       1162 :       if (dP == 1) continue;
    2608                 :        917 :       rkM = lg(S)-1;
    2609                 :        917 :       M2 = cgetg(N+1, t_MAT); /* we will add to S the elts of M2 */
    2610                 :        917 :       gel(M2,1) = col_ei(N, i); /* nf.zk[i] */
    2611                 :        917 :       k = 2;
    2612         [ +  + ]:       2225 :       for (h = 1; h < dP; h++)
    2613                 :            :       {
    2614                 :            :         long r; /* add to M2 the elts of S * nf.zk[i]  */
    2615         [ +  + ]:       7381 :         for (j = 1; j <= rkM; j++) gel(M2,k++) = ZM_ZC_mul(Mx, gel(S,j));
    2616                 :       1308 :         setlg(M2, k); k = 1;
    2617                 :       1308 :         S = ZM_image_shallow(shallowconcat(S,M2), &r);
    2618         [ +  + ]:       1308 :         if (r == rkM) break;
    2619         [ +  - ]:        818 :         if (r > rkM)
    2620                 :            :         {
    2621                 :        818 :           rkM = r;
    2622         [ +  + ]:        818 :           if (rkM == N) break;
    2623                 :            :         }
    2624                 :            :       }
    2625         [ +  + ]:        917 :       if (rkM == N) break;
    2626                 :            :       /* Q(w[1],...,w[i-1]) is a strict subfield of nf */
    2627                 :            :     }
    2628                 :            :   }
    2629                 :            :   /* x_1,...,x_skipfirst generate a strict subfield [unless N=skipfirst=1] */
    2630                 :        280 :   chk->skipfirst = skipfirst;
    2631         [ -  + ]:        280 :   if (DEBUGLEVEL>2) err_printf("chk_gen_init: skipfirst = %ld\n",skipfirst);
    2632                 :            : 
    2633                 :            :   /* should be DEF + gexpo( max_k C^n_k (bound/k)^(k/2) ) */
    2634                 :        280 :   bound = gerepileuptoleaf(av, bound);
    2635                 :        280 :   prec = chk_gen_prec(N, (gexpo(bound)*N)/2);
    2636         [ -  + ]:        280 :   if (DEBUGLEVEL)
    2637                 :          0 :     err_printf("chk_gen_init: new prec = %ld (initially %ld)\n", prec, d->prec);
    2638         [ -  + ]:        280 :   if (prec > d->prec) pari_err_BUG("polredabs (precision problem)");
    2639         [ +  + ]:        280 :   if (prec < d->prec) d->ZKembed = gprec_w(M, prec);
    2640                 :        280 :   return bound;
    2641                 :            : }
    2642                 :            : 
    2643                 :            : /* z "small" minimal polynomial of Mod(a,x), deg z = deg x */
    2644                 :            : static GEN
    2645                 :       2415 : store(GEN x, GEN z, GEN a, nfbasic_t *T, long flag, GEN u)
    2646                 :            : {
    2647                 :            :   GEN y, b;
    2648                 :            : 
    2649         [ +  + ]:       2415 :   if (u) a = RgV_RgC_mul(T->bas, ZM_ZC_mul(u, a));
    2650         [ +  + ]:       2415 :   if (flag & (nf_ORIG|nf_ADDZK))
    2651                 :            :   {
    2652                 :        189 :     b = QXQ_reverse(a, x);
    2653         [ +  + ]:        189 :     if (!isint1(T->unscale)) b = gdiv(b, T->unscale); /* not RgX_Rg_div */
    2654                 :            :   }
    2655                 :            :   else
    2656                 :       2226 :     b = NULL;
    2657                 :            : 
    2658         [ +  + ]:       2415 :   if (flag & nf_RAW)
    2659                 :         28 :     y = mkvec2(z, a);
    2660         [ +  + ]:       2387 :   else if (flag & nf_ORIG) /* store phi(b mod z). */
    2661                 :        189 :     y = mkvec2(z, mkpolmod(b,z));
    2662                 :            :   else
    2663                 :       2198 :     y = z;
    2664         [ -  + ]:       2415 :   if (flag & nf_ADDZK)
    2665                 :            :   { /* append integral basis for number field Q[X]/(z) to result */
    2666                 :          0 :     long n = degpol(x);
    2667                 :          0 :     GEN t = RgV_RgM_mul(RgXQ_powers(b, n-1, z), RgV_to_RgM(T->bas,n));
    2668                 :          0 :     y = mkvec2(y, t);
    2669                 :            :   }
    2670                 :       2415 :   return y;
    2671                 :            : }
    2672                 :            : static GEN
    2673                 :        280 : polredabs_aux(nfbasic_t *T, GEN *u)
    2674                 :            : {
    2675                 :            :   long prec;
    2676                 :            :   GEN v;
    2677                 :        280 :   FP_chk_fun chk = { &chk_gen, &chk_gen_init, NULL, NULL, 0 };
    2678                 :            :   nffp_t F;
    2679                 :        280 :   CG_data d; chk.data = (void*)&d;
    2680                 :            : 
    2681                 :        280 :   prec = polred_init(T, &F, &d);
    2682                 :        280 :   d.bound = embed_T2(F.ro, d.r1);
    2683         [ +  + ]:        280 :   if (realprec(d.bound) > prec) d.bound = rtor(d.bound, prec);
    2684                 :            :   for (;;)
    2685                 :            :   {
    2686                 :        300 :     GEN R = R_from_QR(F.G, prec);
    2687         [ +  + ]:        300 :     if (R)
    2688                 :            :     {
    2689                 :        288 :       d.prec = prec;
    2690                 :        288 :       d.M    = F.M;
    2691                 :        288 :       v = fincke_pohst(mkvec(R),NULL,-1, 0, &chk);
    2692         [ +  + ]:        288 :       if (v) break;
    2693                 :            :     }
    2694                 :         20 :     prec = precdbl(prec);
    2695                 :         20 :     get_nf_fp_compo(T, &F, NULL, 1, prec);
    2696         [ -  + ]:         20 :     if (DEBUGLEVEL) pari_warn(warnprec,"polredabs0",prec);
    2697                 :         20 :   }
    2698                 :        280 :   *u = d.u; return v;
    2699                 :            : }
    2700                 :            : 
    2701                 :            : GEN
    2702                 :        294 : polredabs0(GEN x, long flag)
    2703                 :            : {
    2704                 :        294 :   pari_sp av = avma;
    2705                 :            :   long i, l, vx;
    2706                 :            :   GEN y, a, u;
    2707                 :            :   nfbasic_t T;
    2708                 :            : 
    2709                 :        294 :   nfbasic_init(x, flag & nf_PARTIALFACT, &T);
    2710                 :        294 :   x = T.x; vx = varn(x);
    2711                 :            : 
    2712         [ +  + ]:        294 :   if (degpol(x) == 1)
    2713                 :            :   {
    2714                 :         14 :     u = NULL;
    2715                 :         14 :     y = mkvec( pol_x(vx) );
    2716                 :         14 :     a = mkvec( deg1pol_shallow(gen_1, negi(gel(x,2)), vx) );
    2717                 :         14 :     l = 2;
    2718                 :            :   }
    2719                 :            :   else
    2720                 :            :   {
    2721                 :            :     GEN v;
    2722         [ +  + ]:        280 :     if (!(flag & nf_PARTIALFACT)
    2723 [ +  + ][ -  + ]:         98 :         && T.dKP && lg(primes_certify(T.dK, T.dKP)) != 1) return gen_0;
    2724                 :        280 :     v = polredabs_aux(&T, &u);
    2725                 :        280 :     y = gel(v,1);
    2726                 :        280 :     a = gel(v,2); l = lg(a);
    2727         [ +  + ]:       6078 :     for (i=1; i<l; i++)
    2728         [ +  + ]:       5798 :       if (ZX_canon_neg(gel(y,i))) gel(a,i) = ZC_neg(gel(a,i));
    2729                 :        280 :     remove_duplicates(y,a);
    2730                 :        280 :     l = lg(a);
    2731         [ -  + ]:        280 :     if (l == 1)
    2732                 :          0 :       pari_err_BUG("polredabs (missing vector)");
    2733                 :            :   }
    2734         [ -  + ]:        294 :   if (DEBUGLEVEL) err_printf("Found %ld minimal polynomials.\n",l-1);
    2735         [ +  + ]:        294 :   if (flag & nf_ALL) {
    2736         [ +  + ]:       2135 :     for (i=1; i<l; i++) gel(y,i) = store(x, gel(y,i), gel(a,i), &T, flag, u);
    2737                 :            :   } else {
    2738                 :        287 :     GEN z = findmindisc(y, &a);
    2739                 :        287 :     y = store(x, z, a, &T, flag, u);
    2740                 :            :   }
    2741                 :        294 :   return gerepilecopy(av, y);
    2742                 :            : }
    2743                 :            : 
    2744                 :            : GEN
    2745                 :          0 : polredabsall(GEN x, long flun) { return polredabs0(x, flun | nf_ALL); }
    2746                 :            : GEN
    2747                 :          0 : polredabs(GEN x) { return polredabs0(x,0); }
    2748                 :            : GEN
    2749                 :          0 : polredabs2(GEN x) { return polredabs0(x,nf_ORIG); }
    2750                 :            : 
    2751                 :            : /* relative polredabs/best. Returns relative polynomial by default (flag = 0)
    2752                 :            :  * flag & nf_ORIG: + element (base change)
    2753                 :            :  * flag & nf_ABSOLUTE: absolute polynomial */
    2754                 :            : static GEN
    2755                 :         63 : rnfpolred_i(GEN nf, GEN relpol, long flag, long best)
    2756                 :            : {
    2757         [ +  + ]:         63 :   const char *f = best? "rnfpolredbest": "rnfpolredabs";
    2758 [ +  + ][ +  + ]:         63 :   const long abs = ((flag & nf_ORIG) && (flag & nf_ABSOLUTE));
    2759                 :            :   pari_timer ti;
    2760                 :         63 :   GEN listP = NULL, red, bas, A, P, pol, T, rnfeq;
    2761                 :         63 :   long ty = typ(relpol);
    2762                 :         63 :   pari_sp av = avma;
    2763                 :            : 
    2764         [ +  + ]:         63 :   if (ty == t_VEC) {
    2765         [ -  + ]:         14 :     if (lg(relpol) != 3) pari_err_TYPE(f,relpol);
    2766                 :         14 :     listP = gel(relpol,2);
    2767                 :         14 :     relpol = gel(relpol,1);
    2768                 :            :   }
    2769         [ -  + ]:         63 :   if (typ(relpol) != t_POL) pari_err_TYPE(f,relpol);
    2770                 :         63 :   nf = checknf(nf);
    2771         [ -  + ]:         63 :   if (DEBUGLEVEL>1) timer_start(&ti);
    2772                 :         63 :   T = nf_get_pol(nf);
    2773                 :         63 :   relpol = RgX_nffix(f, T, relpol, 0);
    2774 [ +  + ][ +  + ]:         63 :   if (best || (flag & nf_PARTIALFACT))
    2775                 :            :   {
    2776         [ +  + ]:         42 :     if (abs)
    2777                 :            :     {
    2778                 :          7 :       rnfeq = nf_rnfeq(nf, relpol);
    2779                 :          7 :       pol = gel(rnfeq,1);
    2780                 :            :     }
    2781                 :            :     else
    2782                 :            :     {
    2783                 :            :       long sa;
    2784                 :         35 :       pol = rnfequationall(nf, relpol, &sa, NULL);
    2785                 :         35 :       rnfeq = mkvec5(gen_0,gen_0,stoi(sa),T,liftpol_shallow(relpol));
    2786                 :            :     }
    2787         [ +  + ]:         42 :     bas = listP? mkvec2(pol, listP): pol;
    2788         [ +  + ]:         84 :     if (best)
    2789                 :            :     {
    2790         [ +  + ]:         35 :       if (abs) red = polredbest(bas, 1);
    2791                 :            :       else
    2792                 :            :       {
    2793                 :            :         GEN ro, x, dx, a;
    2794                 :            :         nfbasic_t T;
    2795                 :         28 :         nfbasic_init(bas, nf_PARTIALFACT, &T);
    2796                 :         28 :         polredbest_aux(&T, &ro, &x, &dx, &a);
    2797                 :         28 :         red = mkvec2(x, a);
    2798                 :            :       }
    2799                 :            :     }
    2800                 :            :     else
    2801         [ -  + ]:          7 :       red = polredabs0(bas, (abs? nf_ORIG: nf_RAW)|nf_PARTIALFACT);
    2802                 :            :   }
    2803                 :            :   else
    2804                 :            :   {
    2805                 :         21 :     GEN rnf = rnfinit(nf, relpol), M = rnf_basM(rnf);
    2806                 :         21 :     rnfeq = rnf_get_map(rnf);
    2807                 :         21 :     pol = rnf_get_polabs(rnf);
    2808                 :         21 :     bas = mkvec2(pol, RgM_to_RgXV(M, varn(pol)));
    2809         [ -  + ]:         21 :     if (DEBUGLEVEL>1) timer_printf(&ti, "absolute basis");
    2810                 :         21 :     red = polredabs0(bas, nf_RAW);
    2811                 :            :   }
    2812                 :         63 :   P = gel(red,1);
    2813                 :         63 :   A = gel(red,2);
    2814         [ -  + ]:         63 :   if (DEBUGLEVEL>1) err_printf("reduced absolute generator: %Ps\n",P);
    2815         [ +  + ]:         63 :   if (flag & nf_ABSOLUTE)
    2816                 :            :   {
    2817         [ +  + ]:         14 :     if (flag & nf_ORIG)
    2818                 :            :     {
    2819                 :          7 :       GEN a = gel(rnfeq,2); /* Mod(a,pol) root of T */
    2820                 :          7 :       GEN k = gel(rnfeq,3); /* Mod(variable(relpol),relpol) + k*a root of pol */
    2821                 :          7 :       a = RgX_RgXQ_eval(a, lift_intern(A), P); /* Mod(a, P) root of T */
    2822                 :          7 :       P = mkvec3(P, mkpolmod(a,P), gsub(A, gmul(k,a)));
    2823                 :            :     }
    2824                 :         14 :     return gerepilecopy(av, P);
    2825                 :            :   }
    2826                 :         49 :   A = eltabstorel_lift(rnfeq, A);
    2827                 :         49 :   P = RgXQ_charpoly(A, relpol, varn(relpol));
    2828                 :         49 :   P = lift_if_rational(P);
    2829         [ +  + ]:         49 :   if (flag & nf_ORIG) P = mkvec2(P, mkpolmod(RgXQ_reverse(A,relpol),P));
    2830                 :         63 :   return gerepilecopy(av, P);
    2831                 :            : }
    2832                 :            : GEN
    2833                 :         28 : rnfpolredabs(GEN nf, GEN relpol, long flag)
    2834                 :         28 : { return rnfpolred_i(nf,relpol,flag, 0); }
    2835                 :            : GEN
    2836                 :         35 : rnfpolredbest(GEN nf, GEN relpol, long flag)
    2837                 :            : {
    2838 [ +  - ][ -  + ]:         35 :   if (flag < 0 || flag > 3) pari_err_FLAG("rnfpolredbest");
    2839                 :         35 :   return rnfpolred_i(nf,relpol,flag, 1);
    2840                 :            : }

Generated by: LCOV version 1.9