Code coverage tests

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

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

The target is 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.10.0 lcov report (development 20777-d2a9243) Lines: 390 484 80.6 %
Date: 2017-06-25 05:59:24 Functions: 27 33 81.8 %
Legend: Lines: hit not hit

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

Generated by: LCOV version 1.11