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

Generated by: LCOV version 1.9