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 - buch4.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.8.0 lcov report (development 17072-bc6ca01) Lines: 380 465 81.7 %
Date: 2014-11-17 Functions: 25 31 80.6 %
Legend: Lines: hit not hit | Branches: + taken - not taken # not executed Branches: 207 335 61.8 %

           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                 :            : /*               S-CLASS GROUP AND NORM SYMBOLS                    */
      17                 :            : /*          (Denis Simon, desimon@math.u-bordeaux.fr)              */
      18                 :            : /*                                                                 */
      19                 :            : /*******************************************************************/
      20                 :            : #include "pari.h"
      21                 :            : #include "paripriv.h"
      22                 :            : 
      23                 :            : /* p > 2, T ZX, p prime, x t_INT */
      24                 :            : static long
      25                 :          0 : lemma6(GEN T, GEN p, long nu, GEN x)
      26                 :            : {
      27                 :            :   long la, mu;
      28                 :          0 :   pari_sp av = avma;
      29                 :          0 :   GEN gpx, gx = poleval(T, x);
      30                 :            : 
      31         [ #  # ]:          0 :   if (Zp_issquare(gx, p)) { avma = av; return 1; }
      32                 :            : 
      33                 :          0 :   la = Z_pval(gx, p);
      34                 :          0 :   gpx = poleval(ZX_deriv(T), x);
      35                 :          0 :   mu = signe(gpx)? Z_pval(gpx,p)
      36         [ #  # ]:          0 :                  : la+nu+1; /* mu = +oo */
      37                 :          0 :   avma = av;
      38         [ #  # ]:          0 :   if (la > mu<<1) return 1;
      39 [ #  # ][ #  # ]:          0 :   if (la >= nu<<1 && mu >= nu) return 0;
      40                 :          0 :   return -1;
      41                 :            : }
      42                 :            : /* p = 2, T ZX, x t_INT: return 1 = yes, -1 = no, 0 = inconclusive */
      43                 :            : static long
      44                 :          0 : lemma7(GEN T, long nu, GEN x)
      45                 :            : {
      46                 :            :   long odd4, la, mu;
      47                 :          0 :   pari_sp av = avma;
      48                 :          0 :   GEN gpx, oddgx, gx = poleval(T, x);
      49                 :            : 
      50         [ #  # ]:          0 :   if (Zp_issquare(gx,gen_2)) return 1;
      51                 :            : 
      52                 :          0 :   gpx = poleval(ZX_deriv(T), x);
      53                 :          0 :   la = Z_lvalrem(gx, 2, &oddgx);
      54                 :          0 :   odd4 = umodiu(oddgx,4); avma = av;
      55                 :            : 
      56                 :          0 :   mu = vali(gpx);
      57         [ #  # ]:          0 :   if (mu < 0) mu = la+nu+1; /* mu = +oo */
      58                 :            : 
      59         [ #  # ]:          0 :   if (la > mu<<1) return 1;
      60         [ #  # ]:          0 :   if (nu > mu)
      61                 :            :   {
      62                 :          0 :     long mnl = mu+nu-la;
      63         [ #  # ]:          0 :     if (odd(la)) return -1;
      64         [ #  # ]:          0 :     if (mnl==1) return 1;
      65 [ #  # ][ #  # ]:          0 :     if (mnl==2 && odd4==1) return 1;
      66                 :            :   }
      67                 :            :   else
      68                 :            :   {
      69                 :          0 :     long nu2 = nu << 1;
      70         [ #  # ]:          0 :     if (la >= nu2) return 0;
      71 [ #  # ][ #  # ]:          0 :     if (la == nu2 - 2 && odd4==1) return 0;
      72                 :            :   }
      73                 :          0 :   return -1;
      74                 :            : }
      75                 :            : 
      76                 :            : /* T a ZX, p a prime, pnu = p^nu, x0 t_INT */
      77                 :            : static long
      78                 :          0 : zpsol(GEN T, GEN p, long nu, GEN pnu, GEN x0)
      79                 :            : {
      80                 :            :   long i, res;
      81                 :          0 :   pari_sp av = avma;
      82                 :            :   GEN x, pnup;
      83                 :            : 
      84         [ #  # ]:          0 :   res = equaliu(p,2)? lemma7(T,nu,x0): lemma6(T,p,nu,x0);
      85         [ #  # ]:          0 :   if (res== 1) return 1;
      86         [ #  # ]:          0 :   if (res==-1) return 0;
      87                 :          0 :   x = x0; pnup = mulii(pnu,p);
      88         [ #  # ]:          0 :   for (i=0; i < itos(p); i++)
      89                 :            :   {
      90                 :          0 :     x = addii(x,pnu);
      91         [ #  # ]:          0 :     if (zpsol(T,p,nu+1,pnup,x)) { avma = av; return 1; }
      92                 :            :   }
      93                 :          0 :   avma = av; return 0;
      94                 :            : }
      95                 :            : 
      96                 :            : /* return 1 if equation y^2=T(x) has a rational p-adic solution (possibly
      97                 :            :  * infinite), 0 otherwise. */
      98                 :            : long
      99                 :          0 : hyperell_locally_soluble(GEN T,GEN p)
     100                 :            : {
     101                 :          0 :   pari_sp av = avma;
     102                 :            :   long res;
     103         [ #  # ]:          0 :   if (typ(T)!=t_POL) pari_err_TYPE("zpsoluble",T);
     104         [ #  # ]:          0 :   if (typ(p)!=t_INT) pari_err_TYPE("zpsoluble",p);
     105                 :          0 :   RgX_check_ZX(T, "zpsoluble");
     106 [ #  # ][ #  # ]:          0 :   res = zpsol(T,p,0,gen_1,gen_0) || zpsol(RgX_recip_shallow(T), p, 1, p, gen_0);
     107                 :          0 :   avma = av; return res;
     108                 :            : }
     109                 :            : 
     110                 :            : /* is t a square in (O_K/pr) ? Assume v_pr(t) = 0 */
     111                 :            : static long
     112                 :        140 : quad_char(GEN nf, GEN t, GEN pr)
     113                 :            : {
     114                 :        140 :   GEN ord, ordp, T, p, modpr = zk_to_Fq_init(nf, &pr,&T,&p);
     115                 :        140 :   t = nf_to_Fq(nf,t,modpr);
     116         [ +  - ]:        140 :   if (T)
     117                 :            :   {
     118                 :        140 :     ord = subis( pr_norm(pr), 1 ); /* |(O_K / pr)^*| */
     119                 :        140 :     ordp= subis( p, 1);            /* |F_p^*|        */
     120                 :        140 :     t = Fq_pow(t, diviiexact(ord, ordp), T,p); /* in F_p^* */
     121         [ +  - ]:        140 :     if (typ(t) == t_POL)
     122                 :            :     {
     123         [ -  + ]:        140 :       if (degpol(t)) pari_err_BUG("nfhilbertp");
     124                 :        140 :       t = gel(t,2);
     125                 :            :     }
     126                 :            :   }
     127                 :        140 :   return kronecker(t, p);
     128                 :            : }
     129                 :            : /* quad_char(x), x in Z, non-zero mod p */
     130                 :            : static long
     131                 :        154 : Z_quad_char(GEN x, GEN pr)
     132                 :            : {
     133                 :        154 :   long f = pr_get_f(pr);
     134         [ -  + ]:        154 :   if (!odd(f)) return 1;
     135                 :        154 :   return kronecker(x, pr_get_p(pr));
     136                 :            : }
     137                 :            : 
     138                 :            : /* (pr,2) = 1. return 1 if x in Z_K is a square in Z_{K_pr}, 0 otherwise.
     139                 :            :  * modpr = zkmodprinit(nf,pr) */
     140                 :            : static long
     141                 :          0 : psquarenf(GEN nf,GEN x,GEN pr,GEN modpr)
     142                 :            : {
     143                 :          0 :   pari_sp av = avma;
     144                 :          0 :   GEN p = pr_get_p(pr);
     145                 :            :   long v;
     146                 :            : 
     147                 :          0 :   x = nf_to_scalar_or_basis(nf, x);
     148         [ #  # ]:          0 :   if (typ(x) == t_INT) {
     149         [ #  # ]:          0 :     if (!signe(x)) return 1;
     150                 :          0 :     v = Z_pvalrem(x, p, &x) * pr_get_e(pr);
     151         [ #  # ]:          0 :     if (v&1) return 0;
     152                 :          0 :     v = (Z_quad_char(x, pr) == 1);
     153                 :            :   } else {
     154                 :          0 :     v = ZC_nfvalrem(nf, x, pr, &x);
     155         [ #  # ]:          0 :     if (v&1) return 0;
     156                 :          0 :     v = (quad_char(nf, x, modpr) == 1);
     157                 :            :   }
     158                 :          0 :   avma = av; return v;
     159                 :            : }
     160                 :            : 
     161                 :            : /* Is  x a square in (ZK / pr^(1+2e))^* ?  pr | 2 */
     162                 :            : static long
     163                 :      10983 : check2(GEN nf, GEN x, GEN zinit)
     164                 :            : {
     165                 :      10983 :   GEN zlog = ideallog(nf, x, zinit);
     166                 :      10983 :   long i, l = lg(zlog);
     167         [ +  + ]:      40726 :   for (i=1; i<l; i++) /* all elementary divisors are even (1+2e > 1) */
     168         [ +  + ]:      40278 :     if (mpodd(gel(zlog,i))) return 0;
     169                 :      10983 :   return 1;
     170                 :            : }
     171                 :            : 
     172                 :            : /* pr | 2. Return 1 if x in Z_K is square in Z_{K_pr}, 0 otherwise */
     173                 :            : static int
     174                 :       6279 : psquare2nf_i(GEN nf,GEN x,GEN pr,GEN zinit)
     175                 :            : {
     176                 :       6279 :   long v = nfvalrem(nf, x, pr, &x);
     177                 :            :   /* now (x,pr) = 1 */
     178 [ +  + ][ +  + ]:       6279 :   return v == LONG_MAX || (!odd(v) && check2(nf,x,zinit));
                 [ +  + ]
     179                 :            : }
     180                 :            : static int
     181                 :       6279 : psquare2nf(GEN nf,GEN x,GEN pr,GEN zinit)
     182                 :            : {
     183                 :       6279 :   pari_sp av = avma;
     184                 :       6279 :   long v = psquare2nf_i(nf,x,pr,zinit);
     185                 :       6279 :   avma = av; return v;
     186                 :            : }
     187                 :            : 
     188                 :            : /* pr above an odd prime */
     189                 :            : static long
     190                 :          0 : lemma6nf(GEN nf, GEN T, GEN pr, long nu, GEN x, GEN modpr)
     191                 :            : {
     192                 :          0 :   pari_sp av = avma;
     193                 :            :   long la, mu;
     194                 :          0 :   GEN gpx, gx = nfpoleval(nf, T, x);
     195                 :            : 
     196         [ #  # ]:          0 :   if (psquarenf(nf,gx,pr,modpr)) return 1;
     197                 :            : 
     198                 :          0 :   la = nfval(nf,gx,pr);
     199                 :          0 :   gpx = nfpoleval(nf, RgX_deriv(T), x);
     200         [ #  # ]:          0 :   mu = gequal0(gpx)? la+nu+1: nfval(nf,gpx,pr);
     201                 :          0 :   avma = av;
     202         [ #  # ]:          0 :   if (la > (mu<<1)) return 1;
     203 [ #  # ][ #  # ]:          0 :   if (la >= (nu<<1)  && mu >= nu) return 0;
     204                 :          0 :   return -1;
     205                 :            : }
     206                 :            : /* pr above 2 */
     207                 :            : static long
     208                 :       5775 : lemma7nf(GEN nf, GEN T, GEN pr, long nu, GEN x, GEN zinit)
     209                 :            : {
     210                 :            :   long res, la, mu, q;
     211                 :       5775 :   GEN gpx, gx = nfpoleval(nf, T, x);
     212                 :            : 
     213         [ +  + ]:       5775 :   if (psquare2nf(nf,gx,pr,zinit)) return 1;
     214                 :            : 
     215                 :       5642 :   gpx = nfpoleval(nf, RgX_deriv(T), x);
     216                 :            :   /* gx /= pi^la, pi a pr-uniformizer */
     217                 :       5642 :   la = ZC_nfvalrem(nf, gx, pr, &gx);
     218         [ +  + ]:       5642 :   mu = gequal0(gpx)? la+nu+1: nfval(nf,gpx,pr);
     219                 :            : 
     220         [ -  + ]:       5642 :   if (la > (mu<<1)) return 1;
     221         [ +  + ]:       5642 :   if (nu > mu)
     222                 :            :   {
     223         [ -  + ]:        119 :     if (la&1) return -1;
     224                 :        119 :     q = mu+nu-la; res = 1;
     225                 :            :   }
     226                 :            :   else
     227                 :            :   {
     228                 :       5523 :     long nu2 = nu<<1;
     229         [ +  + ]:       5523 :     if (la >= nu2) return 0;
     230         [ +  + ]:       5222 :     if (odd(la)) return -1;
     231                 :       5026 :     q = nu2-la; res = 0;
     232                 :            :   }
     233         [ +  + ]:       5145 :   if (q > pr_get_e(pr)<<1)  return -1;
     234         [ -  + ]:       4949 :   if (q == 1) return res;
     235                 :            : 
     236                 :            :   /* is gx a square mod pi^q ? FIXME : highly inefficient */
     237                 :       4949 :   zinit = zidealstarinit(nf, idealpows(nf,pr,q));
     238         [ +  + ]:       4949 :   if (!check2(nf, gx, zinit)) res = -1;
     239                 :       5775 :   return res;
     240                 :            : }
     241                 :            : /* zinit either a bid (pr | 2) or a modpr structure (pr | p odd).
     242                 :            :    pnu = pi^nu, pi a uniformizer */
     243                 :            : static long
     244                 :       5775 : zpsolnf(GEN nf,GEN T,GEN pr,long nu,GEN pnu,GEN x0,GEN repr,GEN zinit)
     245                 :            : {
     246                 :            :   long i, res;
     247                 :       5775 :   pari_sp av = avma;
     248                 :            :   GEN pnup;
     249                 :            : 
     250                 :      11550 :   res = typ(zinit) == t_VEC? lemma7nf(nf,T,pr,nu,x0,zinit)
     251         [ +  - ]:       5775 :                            : lemma6nf(nf,T,pr,nu,x0,zinit);
     252                 :       5775 :   avma = av;
     253         [ +  + ]:       5775 :   if (res== 1) return 1;
     254         [ +  + ]:       5642 :   if (res==-1) return 0;
     255                 :        623 :   pnup = nfmul(nf, pnu, pr_get_gen(pr));
     256                 :        623 :   nu++;
     257         [ +  + ]:       5691 :   for (i=1; i<lg(repr); i++)
     258                 :            :   {
     259                 :       5383 :     GEN x = nfadd(nf, x0, nfmul(nf,pnu,gel(repr,i)));
     260         [ +  + ]:       5383 :     if (zpsolnf(nf,T,pr,nu,pnup,x,repr,zinit)) { avma = av; return 1; }
     261                 :            :   }
     262                 :       5775 :   avma = av; return 0;
     263                 :            : }
     264                 :            : 
     265                 :            : /* Let y = copy(x); y[k] := j; return y */
     266                 :            : static GEN
     267                 :       3206 : ZC_add_coeff(GEN x, long k, long j)
     268                 :       3206 : { GEN y = shallowcopy(x); gel(y, k) = utoi(j); return y; }
     269                 :            : 
     270                 :            : /* system of representatives for Zk/pr */
     271                 :            : static GEN
     272                 :        252 : repres(GEN nf, GEN pr)
     273                 :            : {
     274                 :        252 :   long f = pr_get_f(pr), N = nf_get_degree(nf), p = itos(pr_get_p(pr));
     275                 :        252 :   long i, j, k, pi, pf = upowuu(p, f);
     276                 :        252 :   GEN rep, perm = cgetg(f+1, t_VECSMALL);
     277                 :            : 
     278                 :        252 :   perm[1] = 1;
     279         [ +  + ]:        252 :   if (f > 1) {
     280                 :        217 :     GEN H = idealhnf_two(nf,pr);
     281         [ +  + ]:        854 :     for (i = k = 2; k <= f; i++)
     282                 :            :     {
     283         [ -  + ]:        637 :       if (is_pm1(gcoeff(H,i,i))) continue;
     284                 :        637 :       perm[k++] = i;
     285                 :            :     }
     286                 :            :   }
     287                 :        252 :   rep = cgetg(pf+1,t_VEC);
     288                 :        252 :   gel(rep,1) = zerocol(N);
     289         [ +  + ]:       1141 :   for (pi=i=1; i<=f; i++,pi*=p)
     290                 :            :   {
     291                 :        889 :     long t = perm[i];
     292         [ +  + ]:       1778 :     for (j=1; j<p; j++)
     293         [ +  + ]:       4095 :       for (k=1; k<=pi; k++) gel(rep, j*pi+k) = ZC_add_coeff(gel(rep,k), t, j);
     294                 :            :   }
     295                 :        252 :   return rep;
     296                 :            : }
     297                 :            : 
     298                 :            : /* = 1 if equation y^2 = z^deg(T) * T(x/z) has a pr-adic rational solution
     299                 :            :  * (possibly (1,y,0) = oo), 0 otherwise.
     300                 :            :  * coeffs of T are algebraic integers in nf */
     301                 :            : long
     302                 :        252 : nf_hyperell_locally_soluble(GEN nf,GEN T,GEN pr)
     303                 :            : {
     304                 :            :   GEN repr, zinit, p1;
     305                 :        252 :   pari_sp av = avma;
     306                 :            : 
     307         [ -  + ]:        252 :   if (typ(T)!=t_POL) pari_err_TYPE("nf_hyperell_locally_soluble",T);
     308         [ -  + ]:        252 :   if (gequal0(T)) return 1;
     309                 :        252 :   checkprid(pr); nf = checknf(nf);
     310         [ +  - ]:        252 :   if (equaliu(pr_get_p(pr), 2))
     311                 :            :   { /* tough case */
     312                 :        252 :     zinit = Idealstar(nf, idealpows(nf,pr,1+2*pr_get_e(pr)), nf_INIT);
     313         [ -  + ]:        252 :     if (psquare2nf(nf,constant_term(T),pr,zinit)) return 1;
     314         [ -  + ]:        252 :     if (psquare2nf(nf, leading_term(T),pr,zinit)) return 1;
     315                 :            :   }
     316                 :            :   else
     317                 :            :   {
     318                 :          0 :     zinit = zkmodprinit(nf, pr);
     319         [ #  # ]:          0 :     if (psquarenf(nf,constant_term(T),pr,zinit)) return 1;
     320         [ #  # ]:          0 :     if (psquarenf(nf, leading_term(T),pr,zinit)) return 1;
     321                 :            :   }
     322                 :        252 :   repr = repres(nf,pr);
     323         [ +  + ]:        252 :   if (zpsolnf(nf,T,pr,0,gen_1,gen_0,repr,zinit)) { avma=av; return 1; }
     324                 :        140 :   p1 = pr_get_gen(pr);
     325         [ +  + ]:        140 :   if (zpsolnf(nf,RgX_recip_shallow(T),pr,1,p1,gen_0,repr,zinit)) { avma=av; return 1; }
     326                 :            : 
     327                 :        252 :   avma = av; return 0;
     328                 :            : }
     329                 :            : 
     330                 :            : /* return a * denom(a)^2, as an 'liftalg' */
     331                 :            : static GEN
     332                 :        504 : den_remove(GEN nf, GEN a)
     333                 :            : {
     334                 :            :   GEN da;
     335                 :        504 :   a = nf_to_scalar_or_basis(nf, a);
     336   [ +  -  +  - ]:        504 :   switch(typ(a))
     337                 :            :   {
     338                 :         42 :     case t_INT: return a;
     339                 :          0 :     case t_FRAC: return mulii(gel(a,1), gel(a,2));
     340                 :            :     case t_COL:
     341                 :        462 :       a = Q_remove_denom(a, &da);
     342         [ -  + ]:        462 :       if (da) a = ZC_Z_mul(a, da);
     343                 :        462 :       a = coltoliftalg(nf, a);
     344                 :        462 :       return a;
     345                 :          0 :     default: pari_err_TYPE("nfhilbert",a);
     346                 :        504 :       return NULL;/*not reached*/
     347                 :            :   }
     348                 :            : }
     349                 :            : 
     350                 :            : static long
     351                 :        252 : hilb2nf(GEN nf,GEN a,GEN b,GEN p)
     352                 :            : {
     353                 :        252 :   pari_sp av = avma;
     354                 :            :   long rep;
     355                 :            :   GEN pol;
     356                 :            : 
     357                 :        252 :   a = den_remove(nf, a);
     358                 :        252 :   b = den_remove(nf, b);
     359                 :        252 :   pol = mkpoln(3, a, gen_0, b);
     360                 :            :   /* varn(nf.pol) = 0, pol is not a valid GEN  [as in Pol([x,x], x)].
     361                 :            :    * But it is only used as a placeholder, hence it is not a problem */
     362                 :            : 
     363         [ +  + ]:        252 :   rep = nf_hyperell_locally_soluble(nf,pol,p)? 1: -1;
     364                 :        252 :   avma = av; return rep;
     365                 :            : }
     366                 :            : 
     367                 :            : /* local quadratic Hilbert symbol (a,b)_pr, for a,b (non-zero) in nf */
     368                 :            : static long
     369                 :        546 : nfhilbertp(GEN nf, GEN a, GEN b, GEN pr)
     370                 :            : {
     371                 :            :   GEN t;
     372                 :            :   long va, vb, rep;
     373                 :        546 :   pari_sp av = avma;
     374                 :            : 
     375         [ +  + ]:        546 :   if (equaliu(pr_get_p(pr), 2)) return hilb2nf(nf,a,b,pr);
     376                 :            : 
     377                 :            :   /* pr not above 2, compute t = tame symbol */
     378                 :        294 :   va = nfval(nf,a,pr);
     379                 :        294 :   vb = nfval(nf,b,pr);
     380 [ +  + ][ -  + ]:        294 :   if (!odd(va) && !odd(vb)) { avma = av; return 1; }
     381                 :            :   /* Trick: pretend the exponent is 2, result is OK up to squares ! */
     382                 :        294 :   t = famat_makecoprime(nf, mkvec2(a,b), mkvec2s(vb, -va),
     383                 :            :                         pr, idealhnf_two(nf, pr), gen_2);
     384         [ +  + ]:        294 :   if (typ(t) == t_INT) {
     385 [ +  + ][ -  + ]:         14 :     if (odd(va) && odd(vb)) t = negi(t);
     386                 :            :     /* t = (-1)^(v(a)v(b)) a^v(b) b^(-v(a)) */
     387                 :         14 :     rep = Z_quad_char(t, pr);
     388                 :            :   }
     389         [ +  + ]:        280 :   else if (ZV_isscalar(t)) {
     390                 :        140 :     t = gel(t,1);
     391 [ +  - ][ +  - ]:        140 :     if (odd(va) && odd(vb)) t = negi(t);
     392                 :            :     /* t = (-1)^(v(a)v(b)) a^v(b) b^(-v(a)) */
     393                 :        140 :     rep = Z_quad_char(t, pr);
     394                 :            :   } else {
     395 [ +  - ][ -  + ]:        140 :     if (odd(va) && odd(vb)) t = ZC_neg(t);
     396                 :            :     /* t = (-1)^(v(a)v(b)) a^v(b) b^(-v(a)) */
     397                 :        140 :     rep = quad_char(nf, t, pr);
     398                 :            :   }
     399                 :            :   /* quad. symbol is image of t by the quadratic character  */
     400                 :        546 :   avma = av; return rep;
     401                 :            : }
     402                 :            : 
     403                 :            : /* Global quadratic Hilbert symbol (a,b):
     404                 :            :  *  =  1 if X^2 - aY^2 - bZ^2 has a point in projective plane
     405                 :            :  *  = -1 otherwise
     406                 :            :  * a, b should be non-zero */
     407                 :            : long
     408                 :          7 : nfhilbert(GEN nf, GEN a, GEN b)
     409                 :            : {
     410                 :          7 :   pari_sp av = avma;
     411                 :            :   long i, l;
     412                 :            :   GEN S, S2, Sa, Sb, sa, sb;
     413                 :            : 
     414                 :          7 :   nf = checknf(nf);
     415                 :          7 :   a = nf_to_scalar_or_basis(nf, a);
     416                 :          7 :   b = nf_to_scalar_or_basis(nf, b);
     417                 :            :   /* local solutions in real completions ? [ error in nfsign if arg is 0 ]*/
     418                 :          7 :   sa = nfsign(nf, a);
     419                 :          7 :   sb = nfsign(nf, b); l = lg(sa);
     420         [ +  + ]:         14 :   for (i=1; i<l; i++)
     421 [ -  + ][ #  # ]:          7 :     if (sa[i] && sb[i])
     422                 :            :     {
     423         [ #  # ]:          0 :       if (DEBUGLEVEL>3)
     424                 :          0 :         err_printf("nfhilbert not soluble at real place %ld\n",i);
     425                 :          0 :       avma = av; return -1;
     426                 :            :     }
     427                 :            : 
     428                 :            :   /* local solutions in finite completions ? (pr | 2ab)
     429                 :            :    * primes above 2 are toughest. Try the others first */
     430                 :          7 :   Sa = idealfactor(nf, a);
     431                 :          7 :   Sb = idealfactor(nf, b);
     432                 :          7 :   S2 = idealfactor(nf, gen_2);
     433                 :          7 :   S = merge_factor(Sa, Sb, (void*)&cmp_prime_ideal, &cmp_nodata);
     434                 :          7 :   S = merge_factor(S,  S2, (void*)&cmp_prime_ideal, &cmp_nodata);
     435                 :          7 :   S = gel(S,1);
     436                 :            :   /* product of all hilbertp is 1 ==> remove one prime (above 2!) */
     437         [ +  - ]:          7 :   for (i=lg(S)-1; i>1; i--)
     438         [ +  - ]:          7 :     if (nfhilbertp(nf,a,b,gel(S,i)) < 0)
     439                 :            :     {
     440         [ -  + ]:          7 :       if (DEBUGLEVEL>3)
     441                 :          0 :         err_printf("nfhilbert not soluble at finite place %Ps\n",S[i]);
     442                 :          7 :       avma = av; return -1;
     443                 :            :     }
     444                 :          7 :   avma = av; return 1;
     445                 :            : }
     446                 :            : 
     447                 :            : long
     448                 :        560 : nfhilbert0(GEN nf,GEN a,GEN b,GEN p)
     449                 :            : {
     450                 :        560 :   nf = checknf(nf);
     451         [ +  + ]:        560 :   if (p) {
     452                 :        553 :     checkprid(p);
     453         [ +  + ]:        553 :     if (gequal0(a)) pari_err_DOMAIN("nfhilbert", "a", "=", gen_0, a);
     454         [ +  + ]:        546 :     if (gequal0(b)) pari_err_DOMAIN("nfhilbert", "b", "=", gen_0, b);
     455                 :        539 :     return nfhilbertp(nf,a,b,p);
     456                 :            :   }
     457                 :        546 :   return nfhilbert(nf,a,b);
     458                 :            : }
     459                 :            : 
     460                 :            : /* S a list of prime ideal in idealprimedec format. Return res:
     461                 :            :  * res[1] = generators of (S-units / units), as polynomials
     462                 :            :  * res[2] = [perm, HB, den], for bnfissunit
     463                 :            :  * res[3] = [] (was: log. embeddings of res[1])
     464                 :            :  * res[4] = S-regulator ( = R * det(res[2]) * \prod log(Norm(S[i])))
     465                 :            :  * res[5] = S class group
     466                 :            :  * res[6] = S */
     467                 :            : GEN
     468                 :        133 : bnfsunit(GEN bnf,GEN S,long prec)
     469                 :            : {
     470                 :        133 :   pari_sp av = avma;
     471                 :            :   long i,j,ls;
     472                 :            :   GEN p1,nf,gen,M,U,H;
     473                 :            :   GEN sunit,card,sreg,res,pow;
     474                 :            : 
     475         [ -  + ]:        133 :   if (!is_vec_t(typ(S))) pari_err_TYPE("bnfsunit",S);
     476                 :        133 :   bnf = checkbnf(bnf);
     477                 :        133 :   nf = bnf_get_nf(bnf);
     478                 :        133 :   gen = bnf_get_gen(bnf);
     479                 :            : 
     480                 :        133 :   sreg = bnf_get_reg(bnf);
     481                 :        133 :   res=cgetg(7,t_VEC);
     482                 :        133 :   gel(res,1) = gel(res,2) = gel(res,3) = cgetg(1,t_VEC);
     483                 :        133 :   gel(res,4) = sreg;
     484                 :        133 :   gel(res,5) = bnf_get_clgp(bnf);
     485                 :        133 :   gel(res,6) = S; ls=lg(S);
     486                 :            : 
     487                 :            :  /* M = relation matrix for the S class group (in terms of the class group
     488                 :            :   * generators given by gen)
     489                 :            :   * 1) ideals in S
     490                 :            :   */
     491                 :        133 :   M = cgetg(ls,t_MAT);
     492         [ +  + ]:       1869 :   for (i=1; i<ls; i++)
     493                 :            :   {
     494                 :       1736 :     p1 = gel(S,i); checkprid(p1);
     495                 :       1736 :     gel(M,i) = isprincipal(bnf,p1);
     496                 :            :   }
     497                 :            :   /* 2) relations from bnf class group */
     498                 :        133 :   M = shallowconcat(M, diagonal_shallow(bnf_get_cyc(bnf)));
     499                 :            : 
     500                 :            :   /* S class group */
     501                 :        133 :   H = ZM_hnfall(M, &U, 1);
     502                 :        133 :   card = gen_1;
     503         [ +  + ]:        133 :   if (lg(H) > 1)
     504                 :            :   { /* non trivial (rare!) */
     505                 :         56 :     GEN A, u, D = ZM_snfall_i(H, &u, NULL, 1);
     506                 :         56 :     card = detcyc(D, &i); setlg(D,i);
     507                 :         56 :     A = cgetg(i,t_VEC); pow = ZM_inv(u,gen_1);
     508         [ +  + ]:        112 :     for(i--; i; i--) gel(A,i) = idealfactorback(nf, gen, gel(pow,i), 1);
     509                 :         56 :     gel(res,5) = mkvec3(card, D, A);
     510                 :            :   }
     511                 :            : 
     512                 :            :   /* S-units */
     513         [ +  - ]:        133 :   if (ls>1)
     514                 :            :   {
     515                 :        133 :     GEN den, Sperm, perm, dep, B, A, U1 = U;
     516                 :            :     long lH, lB;
     517                 :            : 
     518                 :            :    /* U1 = upper left corner of U, invertible. S * U1 = principal ideals
     519                 :            :     * whose generators generate the S-units */
     520                 :        133 :     setlg(U1,ls); p1 = cgetg(ls, t_MAT); /* p1 is junk for mathnfspec */
     521         [ +  + ]:       1869 :     for (i=1; i<ls; i++) { setlg(U1[i],ls); gel(p1,i) = cgetg(1,t_COL); }
     522                 :        133 :     H = mathnfspec(U1,&perm,&dep,&B,&p1);
     523                 :        133 :     lH = lg(H);
     524                 :        133 :     lB = lg(B);
     525 [ +  + ][ -  + ]:        133 :     if (lg(dep) > 1 && lgcols(dep) > 1) pari_err_BUG("bnfsunit");
     526                 :            :    /*                   [ H B  ]            [ H^-1   - H^-1 B ]
     527                 :            :     * perm o HNF(U1) =  [ 0 Id ], inverse = [  0         Id   ]
     528                 :            :     * (permute the rows)
     529                 :            :     * S * HNF(U1) = _integral_ generators for S-units  = sunit */
     530                 :        133 :     Sperm = cgetg(ls, t_VEC); sunit = cgetg(ls, t_VEC);
     531         [ +  + ]:       1869 :     for (i=1; i<ls; i++) Sperm[i] = S[perm[i]]; /* S o perm */
     532                 :            : 
     533                 :        133 :     setlg(Sperm, lH);
     534         [ +  + ]:        175 :     for (i=1; i<lH; i++)
     535                 :            :     {
     536                 :         42 :       GEN v = isprincipalfact(bnf, NULL,Sperm,gel(H,i), nf_GEN|nf_FORCE);
     537                 :         42 :       gel(sunit,i) = coltoliftalg(nf, gel(v,2));
     538                 :            :     }
     539         [ +  + ]:       1827 :     for (j=1; j<lB; j++,i++)
     540                 :            :     {
     541                 :       1694 :       GEN v = isprincipalfact(bnf, gel(Sperm,i),Sperm,gel(B,j),nf_GEN|nf_FORCE);
     542                 :       1694 :       gel(sunit,i) = coltoliftalg(nf, gel(v,2));
     543                 :            :    }
     544                 :        133 :     den = ZM_det_triangular(H); H = ZM_inv(H,den);
     545                 :        133 :     A = shallowconcat(H, ZM_neg(ZM_mul(H,B))); /* top part of inverse * den */
     546                 :            :     /* HNF in split form perm + (H B) [0 Id missing] */
     547                 :        133 :     gel(res,1) = sunit;
     548                 :        133 :     gel(res,2) = mkvec3(perm,A,den);
     549                 :            :   }
     550                 :            : 
     551                 :            :   /* S-regulator */
     552                 :        133 :   sreg = mpmul(sreg,card);
     553         [ +  + ]:       1869 :   for (i=1; i<ls; i++)
     554                 :            :   {
     555                 :       1736 :     GEN p = pr_get_p( gel(S,i) );
     556                 :       1736 :     sreg = mpmul(sreg, logr_abs(itor(p,prec)));
     557                 :            :   }
     558                 :        133 :   gel(res,4) = sreg;
     559                 :        133 :   return gerepilecopy(av,res);
     560                 :            : }
     561                 :            : 
     562                 :            : static GEN
     563                 :       1365 : make_unit(GEN nf, GEN bnfS, GEN *px)
     564                 :            : {
     565                 :            :   long lB, cH, i, ls;
     566                 :            :   GEN den, gen, S, v, p1, xp, xb, N, N0, HB, perm;
     567                 :            : 
     568         [ -  + ]:       1365 :   if (gequal0(*px)) return NULL;
     569                 :       1365 :   S = gel(bnfS,6); ls = lg(S);
     570         [ -  + ]:       1365 :   if (ls==1) return cgetg(1, t_COL);
     571                 :            : 
     572                 :       1365 :   xb = nf_to_scalar_or_basis(nf,*px);
     573      [ +  -  + ]:       1365 :   switch(typ(xb))
     574                 :            :   {
     575                 :        476 :     case t_INT:  N = xb; break;
     576                 :          0 :     case t_FRAC: N = mulii(gel(xb,1),gel(xb,2)); break;
     577                 :        889 :     default: { GEN d = Q_denom(xb); N = mulii(idealnorm(nf,gmul(*px,d)), d); }
     578                 :            :   } /* relevant primes divide N */
     579         [ +  + ]:       1365 :   if (is_pm1(N)) return zerocol(ls -1);
     580                 :            : 
     581                 :       1148 :   p1 = gel(bnfS,2);
     582                 :       1148 :   perm = gel(p1,1);
     583                 :       1148 :   HB   = gel(p1,2);
     584                 :       1148 :   den  = gel(p1,3);
     585                 :       1148 :   cH = nbrows(HB);
     586                 :       1148 :   lB = lg(HB) - cH;
     587                 :       1148 :   v = zero_zv(ls-1);
     588                 :       1148 :   N0 = N;
     589         [ +  + ]:      36645 :   for (i=1; i<ls; i++)
     590                 :            :   {
     591                 :      35497 :     GEN P = gel(S,i), p = pr_get_p(P);
     592         [ +  + ]:      35497 :     if ( remii(N, p) == gen_0 )
     593                 :            :     {
     594                 :       2317 :       v[i] = nfval(nf,xb,P);
     595                 :       2317 :       (void)Z_pvalrem(N0, p, &N0);
     596                 :            :     }
     597                 :            :   }
     598         [ +  + ]:       1148 :   if (!is_pm1(N0)) return NULL;
     599                 :            :   /* here, x = S v */
     600                 :       1141 :   p1 = vecsmallpermute(v, perm);
     601                 :       1141 :   v = ZM_zc_mul(HB, p1);
     602         [ +  + ]:       1505 :   for (i=1; i<=cH; i++)
     603                 :            :   {
     604                 :        364 :     GEN r, w = dvmdii(gel(v,i), den, &r);
     605         [ -  + ]:        364 :     if (r != gen_0) return NULL;
     606                 :        364 :     gel(v,i) = w;
     607                 :            :   }
     608                 :       1141 :   p1 += cH; p1[0] = evaltyp(t_VECSMALL) | evallg(lB);
     609                 :       1141 :   v = shallowconcat(v, zc_to_ZC(p1)); /* append bottom of p1 (= [0 Id] part) */
     610                 :            : 
     611                 :       1141 :   gen = gel(bnfS,1);
     612                 :       1141 :   xp = cgetg(1, t_MAT);
     613         [ +  + ]:      36603 :   for (i=1; i<ls; i++)
     614                 :            :   {
     615                 :      35462 :     GEN e = gel(v,i);
     616         [ +  + ]:      35462 :     if (!signe(e)) continue;
     617                 :       1547 :     xp = famat_mul(xp, to_famat_shallow(gel(gen,i), negi(e)));
     618                 :            :   }
     619         [ +  - ]:       1141 :   if (lg(xp) > 1) *px = famat_mul(xp, to_famat_shallow(xb, gen_1));
     620                 :       1365 :   return v;
     621                 :            : }
     622                 :            : 
     623                 :            : /* Analog to bnfisunit, for S-units. Let v the result
     624                 :            :  * If x not an S-unit, v = []~, else
     625                 :            :  * x = \prod_{i=0}^r e_i^v[i] * prod{i=r+1}^{r+s} s_i^v[i]
     626                 :            :  * where the e_i are the field units (cf bnfisunit), and the s_i are
     627                 :            :  * the S-units computed by bnfsunit (in the same order) */
     628                 :            : GEN
     629                 :       1365 : bnfissunit(GEN bnf,GEN bnfS,GEN x)
     630                 :            : {
     631                 :       1365 :   pari_sp av = avma;
     632                 :            :   GEN v, w, nf;
     633                 :            : 
     634                 :       1365 :   bnf = checkbnf(bnf);
     635                 :       1365 :   nf = bnf_get_nf(bnf);
     636 [ +  - ][ -  + ]:       1365 :   if (typ(bnfS)!=t_VEC || lg(bnfS)!=7) pari_err_TYPE("bnfissunit",bnfS);
     637                 :       1365 :   x = nf_to_scalar_or_alg(nf,x);
     638                 :       1365 :   v = NULL;
     639         [ +  + ]:       1365 :   if ( (w = make_unit(nf, bnfS, &x)) ) v = bnfisunit(bnf, x);
     640 [ +  + ][ -  + ]:       1365 :   if (!v || lg(v) == 1) { avma = av; return cgetg(1,t_COL); }
     641                 :       1365 :   return gerepileupto(av, concat(v, w));
     642                 :            : }
     643                 :            : 
     644                 :            : static void
     645                 :        511 : pr_append(GEN nf, GEN rel, GEN p, GEN *prod, GEN *S1, GEN *S2)
     646                 :            : {
     647         [ +  + ]:        959 :   if (dvdii(*prod, p)) return;
     648                 :        448 :   *prod = mulii(*prod, p);
     649                 :        448 :   *S1 = shallowconcat(*S1, idealprimedec(nf,p));
     650                 :        448 :   *S2 = shallowconcat(*S2, idealprimedec(rel,p));
     651                 :            : }
     652                 :            : 
     653                 :            : /* N a t_INT */
     654                 :            : static void
     655                 :        119 : Zfa_pr_append(GEN nf,GEN rel,GEN N,GEN *prod,GEN *S1,GEN *S2)
     656                 :            : {
     657         [ +  + ]:        119 :   if (!is_pm1(N))
     658                 :            :   {
     659                 :        112 :     GEN v = gel(Z_factor(N),1);
     660                 :        112 :     long i, l = lg(v);
     661         [ +  + ]:        273 :     for (i=1; i<l; i++) pr_append(nf,rel,gel(v,i),prod,S1,S2);
     662                 :            :   }
     663                 :        119 : }
     664                 :            : /* N a t_INT or t_FRAC */
     665                 :            : static void
     666                 :         63 : fa_pr_append(GEN nf,GEN rel,GEN N,GEN *prod,GEN *S1,GEN *S2)
     667                 :            : {
     668         [ -  + ]:         63 :   if (typ(N) == t_FRAC)
     669                 :            :   {
     670                 :          0 :     Zfa_pr_append(nf,rel,gel(N,1),prod,S1,S2);
     671                 :          0 :     Zfa_pr_append(nf,rel,gel(N,2),prod,S1,S2);
     672                 :            :   }
     673                 :            :   else
     674                 :         63 :     Zfa_pr_append(nf,rel,N,prod,S1,S2);
     675                 :         63 : }
     676                 :            : 
     677                 :            : /* apply lift(rnfeltup) to all coeffs, without rnf structure */
     678                 :            : static GEN
     679                 :          7 : nfX_eltup(GEN nf, GEN rnfeq, GEN x)
     680                 :            : {
     681                 :            :   long i, l;
     682                 :          7 :   GEN zknf, czknf, y = cgetg_copy(x, &l);
     683                 :          7 :   y[1] = x[1]; nf_nfzk(nf, rnfeq, &zknf, &czknf);
     684         [ +  + ]:         35 :   for (i=2; i<l; i++) gel(y,i) = nfeltup(nf, gel(x,i), zknf, czknf);
     685                 :          7 :   return y;
     686                 :            : }
     687                 :            : 
     688                 :            : GEN
     689                 :         49 : rnfisnorminit(GEN T, GEN relpol, int galois)
     690                 :            : {
     691                 :         49 :   pari_sp av = avma;
     692                 :            :   long i, l, drel;
     693                 :            :   GEN prod, S1, S2, gen, cyc, bnf, nf, nfabs, rnfeq, bnfabs, k, polabs;
     694                 :         49 :   GEN y = cgetg(9, t_VEC);
     695                 :            : 
     696 [ +  - ][ -  + ]:         49 :   if (galois < 0 || galois > 2) pari_err_FLAG("rnfisnorminit");
     697                 :         49 :   T = get_bnfpol(T, &bnf, &nf);
     698 [ +  + ][ -  + ]:         49 :   if (!bnf) bnf = Buchall(nf? nf: T, nf_FORCE, DEFAULTPREC);
     699         [ +  + ]:         49 :   if (!nf) nf = bnf_get_nf(bnf);
     700                 :            : 
     701                 :         49 :   relpol = get_bnfpol(relpol, &bnfabs, &nfabs);
     702         [ -  + ]:         49 :   if (!gequal1(leading_term(relpol))) pari_err_IMPL("non monic relative equation");
     703                 :         49 :   drel = degpol(relpol);
     704         [ +  + ]:         49 :   if (drel <= 2) galois = 1;
     705                 :            : 
     706                 :         49 :   relpol = RgX_nffix("rnfisnorminit", T, relpol, 1);
     707         [ +  + ]:         49 :   if (nf_get_degree(nf) == 1) /* over Q */
     708                 :         35 :     rnfeq = mkvec5(relpol,gen_0,gen_0,T,relpol);
     709         [ +  + ]:         14 :   else if (galois == 2) /* needs eltup+abstorel */
     710                 :          7 :     rnfeq = nf_rnfeq(nf, relpol);
     711                 :            :   else /* needs abstorel */
     712                 :          7 :     rnfeq = nf_rnfeqsimple(nf, relpol);
     713                 :         49 :   polabs = gel(rnfeq,1);
     714                 :         49 :   k = gel(rnfeq,3);
     715 [ +  + ][ -  + ]:         49 :   if (!bnfabs || !gequal0(k))
     716                 :         21 :     bnfabs = Buchall(polabs, nf_FORCE, nf_get_prec(nf));
     717         [ +  + ]:         49 :   if (!nfabs) nfabs = bnf_get_nf(bnfabs);
     718                 :            : 
     719         [ +  + ]:         49 :   if (galois == 2)
     720                 :            :   {
     721         [ +  + ]:         21 :     GEN P = polabs==relpol? leafcopy(relpol): nfX_eltup(nf, rnfeq, relpol);
     722                 :         21 :     setvarn(P, fetch_var_higher());
     723                 :         21 :     galois = nfissplit(nfabs, P);
     724                 :         21 :     (void)delete_var();
     725                 :            :   }
     726                 :            : 
     727                 :         49 :   prod = gen_1; S1 = S2 = cgetg(1, t_VEC);
     728                 :         49 :   cyc = bnf_get_cyc(bnfabs);
     729                 :         49 :   gen = bnf_get_gen(bnfabs); l = lg(cyc);
     730         [ +  + ]:         77 :   for(i=1; i<l; i++)
     731                 :            :   {
     732                 :         35 :     GEN g = gel(gen,i);
     733         [ +  + ]:         35 :     if (ugcd(umodiu(gel(cyc,i), drel), drel) == 1) break;
     734                 :         28 :     Zfa_pr_append(nf,bnfabs,gcoeff(g,1,1),&prod,&S1,&S2);
     735                 :            :   }
     736         [ +  + ]:         49 :   if (!galois)
     737                 :            :   {
     738                 :         21 :     GEN Ndiscrel = diviiexact(nf_get_disc(nfabs), powiu(nf_get_disc(nf), drel));
     739                 :         21 :     Zfa_pr_append(nf,bnfabs,absi(Ndiscrel), &prod,&S1,&S2);
     740                 :            :   }
     741                 :            : 
     742                 :         49 :   gel(y,1) = bnf;
     743                 :         49 :   gel(y,2) = bnfabs;
     744                 :         49 :   gel(y,3) = relpol;
     745                 :         49 :   gel(y,4) = rnfeq;
     746                 :         49 :   gel(y,5) = prod;
     747                 :         49 :   gel(y,6) = S1;
     748                 :         49 :   gel(y,7) = S2;
     749                 :         49 :   gel(y,8) = stoi(galois); return gerepilecopy(av, y);
     750                 :            : }
     751                 :            : 
     752                 :            : /* T as output by rnfisnorminit
     753                 :            :  * if flag=0 assume extension is Galois (==> answer is unconditional)
     754                 :            :  * if flag>0 add to S all primes dividing p <= flag
     755                 :            :  * if flag<0 add to S all primes dividing abs(flag)
     756                 :            : 
     757                 :            :  * answer is a vector v = [a,b] such that
     758                 :            :  * x = N(a)*b and x is a norm iff b = 1  [assuming S large enough] */
     759                 :            : GEN
     760                 :         63 : rnfisnorm(GEN T, GEN x, long flag)
     761                 :            : {
     762                 :         63 :   pari_sp av = avma;
     763                 :            :   GEN bnf, rel, relpol, rnfeq, nfpol;
     764                 :            :   GEN nf, aux, H, U, Y, M, A, bnfS, sunitrel, futu, prod, S1, S2;
     765                 :            :   long L, i, drel, itu;
     766                 :            : 
     767 [ +  - ][ -  + ]:         63 :   if (typ(T) != t_VEC || lg(T) != 9)
     768                 :          0 :     pari_err_TYPE("rnfisnorm [please apply rnfisnorminit()]", T);
     769                 :         63 :   bnf = gel(T,1);
     770                 :         63 :   rel = gel(T,2);
     771                 :         63 :   relpol = gel(T,3);
     772                 :         63 :   rnfeq = gel(T,4);
     773                 :         63 :   drel = degpol(relpol);
     774                 :         63 :   bnf = checkbnf(bnf);
     775                 :         63 :   rel = checkbnf(rel);
     776                 :         63 :   nf = bnf_get_nf(bnf);
     777                 :         63 :   x = nf_to_scalar_or_alg(nf,x);
     778         [ -  + ]:         63 :   if (gequal0(x)) { avma = av; return mkvec2(gen_0, gen_1); }
     779         [ -  + ]:         63 :   if (gequal1(x)) { avma = av; return mkvec2(gen_1, gen_1); }
     780 [ -  + ][ #  # ]:         63 :   if (gequalm1(x) && odd(drel)) { avma = av; return mkvec2(gen_m1, gen_1); }
     781                 :            : 
     782                 :            :   /* build set T of ideals involved in the solutions */
     783                 :         63 :   nfpol = nf_get_pol(nf);
     784                 :         63 :   prod = gel(T,5);
     785                 :         63 :   S1   = gel(T,6);
     786                 :         63 :   S2   = gel(T,7);
     787 [ +  + ][ +  + ]:         63 :   if (flag && !gequal0(gel(T,8)))
     788                 :          7 :     pari_warn(warner,"useless flag in rnfisnorm: the extension is Galois");
     789         [ +  + ]:         63 :   if (flag > 0)
     790                 :            :   {
     791                 :            :     forprime_t T;
     792                 :            :     ulong p;
     793                 :         14 :     u_forprime_init(&T, 2, flag);
     794         [ +  + ]:        364 :     while ( (p = u_forprime_next(&T)) )
     795                 :        350 :       pr_append(nf,rel, utoipos(p),&prod,&S1,&S2);
     796                 :            :   }
     797         [ +  + ]:         49 :   else if (flag < 0)
     798                 :          7 :     Zfa_pr_append(nf,rel, utoipos(-flag),&prod,&S1,&S2);
     799                 :            :   /* overkill: prime ideals dividing x would be enough */
     800                 :         63 :   fa_pr_append(nf,rel,idealnorm(nf,x), &prod,&S1,&S2);
     801                 :            : 
     802                 :            :   /* computation on T-units */
     803                 :         63 :   futu = shallowconcat(bnf_get_fu(rel), bnf_get_tuU(rel));
     804                 :         63 :   bnfS = bnfsunit(bnf,S1,LOWDEFAULTPREC);
     805                 :         63 :   sunitrel = shallowconcat(futu, gel(bnfsunit(rel,S2,LOWDEFAULTPREC), 1));
     806                 :            : 
     807                 :         63 :   A = lift_intern(bnfissunit(bnf,bnfS,x));
     808                 :         63 :   L = lg(sunitrel);
     809                 :         63 :   itu = lg(nf_get_roots(nf))-1; /* index of torsion unit in bnfsunit(nf) output */
     810                 :         63 :   M = cgetg(L+1,t_MAT);
     811         [ +  + ]:       1344 :   for (i=1; i<L; i++)
     812                 :            :   {
     813                 :       1281 :     GEN u = eltabstorel(rnfeq, gel(sunitrel,i));
     814                 :       1281 :     gel(sunitrel,i) = u;
     815                 :       1281 :     u = bnfissunit(bnf,bnfS, gnorm(u));
     816         [ -  + ]:       1281 :     if (lg(u) == 1) pari_err_BUG("rnfisnorm");
     817                 :       1281 :     gel(u,itu) = lift_intern(gel(u,itu)); /* lift root of 1 part */
     818                 :       1281 :     gel(M,i) = u;
     819                 :            :   }
     820                 :         63 :   aux = zerocol(lg(A)-1); gel(aux,itu) = utoipos( bnf_get_tuN(rel) );
     821                 :         63 :   gel(M,L) = aux;
     822                 :         63 :   H = ZM_hnfall(M, &U, 2);
     823                 :         63 :   Y = RgM_RgC_mul(U, inverseimage(H,A));
     824                 :            :   /* Y: sols of MY = A over Q */
     825                 :         63 :   setlg(Y, L);
     826                 :         63 :   aux = factorback2(sunitrel, gfloor(Y));
     827                 :         63 :   x = mkpolmod(x,nfpol);
     828         [ +  + ]:         63 :   if (!gequal1(aux)) x = gdiv(x, gnorm(aux));
     829                 :         63 :   x = lift_if_rational(x);
     830 [ +  + ][ +  + ]:         63 :   if (typ(aux) == t_POLMOD && degpol(nfpol) == 1)
     831                 :         28 :     gel(aux,2) = lift_if_rational(gel(aux,2));
     832                 :         63 :   return gerepilecopy(av, mkvec2(aux, x));
     833                 :            : }
     834                 :            : 
     835                 :            : GEN
     836                 :         28 : bnfisnorm(GEN bnf, GEN x, long flag)
     837                 :            : {
     838                 :         28 :   pari_sp av = avma;
     839         [ +  + ]:         28 :   GEN T = rnfisnorminit(pol_x(fetch_var()), bnf, flag == 0? 1: 2);
     840         [ +  + ]:         28 :   GEN r = rnfisnorm(T, x, flag == 1? 0: flag);
     841                 :         28 :   (void)delete_var();
     842                 :         28 :   return gerepileupto(av,r);
     843                 :            : }

Generated by: LCOV version 1.9