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 - modules - krasner.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.10.0 lcov report (development 21059-cbe0d6a) Lines: 439 447 98.2 %
Date: 2017-09-22 06:24:58 Functions: 34 34 100.0 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2009  The PARI group.
       2             : 
       3             : This file is part of the PARI/GP package.
       4             : 
       5             : PARI/GP is free software; you can redistribute it and/or modify it under the
       6             : terms of the GNU General Public License as published by the Free Software
       7             : Foundation. It is distributed in the hope that it will be useful, but WITHOUT
       8             : ANY WARRANTY WHATSOEVER.
       9             : 
      10             : Check the License for details. You should have received a copy of it, along
      11             : with the package; see the file 'COPYING'. If not, write to the Free Software
      12             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      13             : 
      14             : #include "pari.h"
      15             : #include "paripriv.h"
      16             : 
      17             : /************************************************************/
      18             : /*     Computation of all the extensions of a given         */
      19             : /*               degree of a p-adic field                   */
      20             : /* Xavier Roblot                                            */
      21             : /************************************************************/
      22             : /* cf. Math. Comp, vol. 70, No. 236, pp. 1641-1659 for more details.
      23             :    Note that n is now e (since the e from the paper is always = 1) and l
      24             :    is now f */
      25             : /* Structure for a given type of extension */
      26             : typedef struct {
      27             :   GEN p;
      28             :   long e, f, j;
      29             :   long a, b, c;
      30             :   long v;     /* auxiliary variable */
      31             :   long r;     /* pr = p^r */
      32             :   GEN pr;     /* p-adic precision for poly. reduction */
      33             :   GEN q, qm1; /* p^f, q - 1 */
      34             :   GEN T;      /* polynomial defining K^ur */
      35             :   GEN frob;   /* Frobenius acting of the root of T (mod pr) */
      36             :   GEN u;      /* suitable root of unity (expressed in terms of the root of T) */
      37             :   GEN nbext;  /* number of extensions */
      38             :   GEN roottable; /* table of roots of polynomials over the residue field */
      39             :   GEN pk;     /* powers of p: [p^1, p^2, ..., p^c] */
      40             : } KRASNER_t;
      41             : 
      42             : /* Structure containing the field data (constructed with FieldData) */
      43             : typedef struct {
      44             :   GEN p;
      45             :   GEN top;  /* absolute polynomial of a = z + pi (mod pr) */
      46             :   GEN topr; /* top mod p */
      47             :   GEN z;    /* root of T in terms of a (mod pr) */
      48             :   GEN eis;  /* relative polynomial of pi in terms of z (mod pr)  */
      49             :   GEN pi;   /* prime element in terms of a */
      50             :   GEN ipi;  /* p/pi in terms of a (mod pr) (used to divide by pi) */
      51             :   GEN pik;  /* [1, pi, ..., pi^(e-1), pi^e/p] in terms of a (mod pr).
      52             :                Note the last one is different! */
      53             :   long cj;  /* number of conjugate fields */
      54             : } FAD_t;
      55             : 
      56             : #undef CHECK
      57             : 
      58             : /* Eval P(x) assuming P has positive coefficients and the result is a long */
      59             : static ulong
      60      231616 : ZX_z_eval(GEN P, ulong x)
      61             : {
      62      231616 :   long i, l = lg(P);
      63             :   ulong z;
      64             : 
      65      231616 :   if (typ(P) == t_INT) return itou(P);
      66       74361 :   if (l == 2) return 0;
      67       42770 :   z = itou(gel(P, l-1));
      68       42770 :   for (i = l-2; i > 1; i--) z = z*x + itou(gel(P,i));
      69       42770 :   return z;
      70             : }
      71             : 
      72             : /* Eval P(x, y) assuming P has positive coefficients and the result is a long */
      73             : static ulong
      74       53907 : ZXY_z_eval(GEN P, ulong x, ulong y)
      75             : {
      76       53907 :   long i, l = lg(P);
      77             :   ulong z;
      78             : 
      79       53907 :   if (l == 2) return 0;
      80       53907 :   z = ZX_z_eval(gel(P, l-1), y);
      81       53907 :   for (i = l-2; i > 1; i--) z = z*x + ZX_z_eval(gel(P, i), y);
      82       53907 :   return z;
      83             : }
      84             : 
      85             : /* P an Fq[X], where Fq = Fp[Y]/(T(Y)), a an FpX representing the automorphism
      86             :  * y -> a(y). Return a(P), applying a() coefficientwise. */
      87             : static GEN
      88       39522 : FqX_FpXQ_eval(GEN P, GEN a, GEN T, GEN p)
      89             : {
      90             :   long i, l;
      91       39522 :   GEN Q = cgetg_copy(P, &l);
      92             : 
      93       39522 :   Q[1] = P[1];
      94      179774 :   for (i = 2; i < l; i++)
      95             :   {
      96      140252 :     GEN cf = gel(P, i);
      97      140252 :     if (typ(cf) == t_POL) {
      98       94969 :       switch(degpol(cf))
      99             :       {
     100         952 :         case -1: cf = gen_0; break;
     101         511 :         case 0:  cf = gel(cf,2); break;
     102             :         default:
     103       93506 :           cf = FpX_FpXQ_eval(cf, a, T, p);
     104       93506 :           cf = simplify_shallow(cf);
     105       93506 :           break;
     106             :       }
     107             :     }
     108      140252 :     gel(Q, i) = cf;
     109             :   }
     110             : 
     111       39522 :   return Q;
     112             : }
     113             : 
     114             : /* Sieving routines */
     115             : static GEN
     116         497 : InitSieve(long l) { return zero_F2v(l); }
     117             : static long
     118         777 : NextZeroValue(GEN sv, long i)
     119             : {
     120        1890 :   for(; i <= sv[1]; i++)
     121        1393 :     if (!F2v_coeff(sv, i)) return i;
     122         497 :   return 0; /* sieve is full or out of the sieve! */
     123             : }
     124             : static void
     125        1113 : SetSieveValue(GEN sv, long i)
     126        1113 : { if (i >= 1 && i <= sv[1]) F2v_set(sv, i); }
     127             : 
     128             : /* return 1 if the data verify Ore condition and 0 otherwise */
     129             : static long
     130          21 : VerifyOre(GEN p, long e, long j)
     131             : {
     132             :   long nv, b, vb, nb;
     133             : 
     134          21 :   if (j < 0) return 0;
     135          21 :   nv = e * u_pval(e, p);
     136          21 :   b  = j%e;
     137          21 :   if (b == 0) return (j == nv);
     138          14 :   if (j > nv) return 0;
     139             :   /* j < nv */
     140          14 :   vb = u_pval(b, p);
     141          14 :   nb = vb*e;
     142          14 :   return (nb <= j);
     143             : }
     144             : 
     145             : /* Given [K:Q_p] = m and disc(K/Q_p) = p^d, return all decompositions K/K^ur/Q_p
     146             :  * as [e, f, j] with
     147             :  *   K^ur/Q_p unramified of degree f,
     148             :  *   K/K^ur totally ramified of degree e and discriminant p^(e+j-1);
     149             :  * thus d = f*(e+j-1) and j > 0 iff ramification is wild */
     150             : static GEN
     151           7 : possible_efj_by_d(GEN p, long m, long d)
     152             : {
     153             :   GEN rep, div;
     154             :   long i, ctr, l;
     155             : 
     156           7 :   if (d == 0) return mkvec(mkvecsmall3(1, m, 0)); /* unramified */
     157           7 :   div = divisorsu(ugcd(m, d));
     158           7 :   l = lg(div); ctr = 1;
     159           7 :   rep = cgetg(l, t_VEC);
     160          28 :   for (i = 1; i < l; i++)
     161             :   {
     162          21 :     long f = div[i], e = m/f, j = d/f - e + 1;
     163          21 :     if (VerifyOre(p, e, j)) gel(rep, ctr++) = mkvecsmall3(e, f, j);
     164             :   }
     165           7 :   setlg(rep, ctr); return rep;
     166             : }
     167             : 
     168             : /* return the number of extensions corresponding to (e,f,j) */
     169             : static GEN
     170        1743 : NumberExtensions(KRASNER_t *data)
     171             : {
     172             :   ulong pp, pa;
     173             :   long e, a, b;
     174             :   GEN p, q, s0, p1;
     175             : 
     176        1743 :   e = data->e;
     177        1743 :   p = data->p;
     178        1743 :   q = data->q;
     179        1743 :   a = data->a; /* floor(j/e) <= v_p(e), hence p^a | e */
     180        1743 :   b = data->b; /* j % e */
     181        1743 :   if (is_bigint(p)) /* implies a = 0 */
     182          70 :     return b == 0? utoipos(e): mulsi(e, data->qm1);
     183             : 
     184        1673 :   pp = p[2];
     185        1673 :   switch(a)
     186             :   {
     187        1617 :     case 0:  pa = 1;  s0 = utoipos(e); break;
     188          49 :     case 1:  pa = pp; s0 = mului(e, powiu(q, e / pa)); break;
     189             :     default:
     190           7 :       pa = upowuu(pp, a); /* p^a */
     191           7 :       s0 = mulsi(e, powiu(q, (e / pa) * ((pa - 1) / (pp - 1))));
     192             :   }
     193             :   /* s0 = e * q^(e*sum(p^-i, i=1...a)) */
     194        1673 :   if (b == 0) return s0;
     195             : 
     196             :   /* q^floor((b-1)/p^(a+1)) */
     197          98 :   p1 = powiu(q, sdivsi(b-1, muluu(pa, pp)));
     198          98 :   return mulii(mulii(data->qm1, s0), p1);
     199             : }
     200             : 
     201             : static GEN
     202        3213 : get_topx(KRASNER_t *data, GEN eis)
     203             : {
     204             :   GEN p1, p2, rpl, c;
     205             :   pari_sp av;
     206             :   long j;
     207             :   /* top poly. is the minimal polynomial of root(T) + root(eis) */
     208        3213 :   if (data->f == 1) return eis;
     209        1953 :   c = FpX_neg(pol_x(data->v),data->pr);
     210        1953 :   rpl = FqX_translate(eis, c, data->T, data->pr);
     211        1953 :   p1 = p2 = rpl; av = avma;
     212       21525 :   for (j = 1; j < data->f; j++)
     213             :   {
     214             :     /* compute conjugate polynomials using the Frobenius */
     215       19572 :     p1 = FqX_FpXQ_eval(p1, data->frob, data->T, data->pr);
     216       19572 :     p2 = FqX_mul(p2, p1, data->T, data->pr);
     217       19572 :     if (gc_needed(av,2)) gerepileall(av, 2, &p1,&p2);
     218             :   }
     219        1953 :   return simplify_shallow(p2); /* ZX */
     220             : }
     221             : 
     222             : /* eis (ZXY): Eisenstein polynomial over the field defined by T.
     223             :  * topx (ZX): absolute equation of root(T) + root(eis).
     224             :  * Return the struct FAD corresponding to the field it defines (GENs created
     225             :  * as clones). Assume e > 1. */
     226             : static void
     227         854 : FieldData(KRASNER_t *data, FAD_t *fdata, GEN eis, GEN topx)
     228             : {
     229             :   GEN p1, p2, p3, z, ipi, cipi, dipi, t, Q;
     230             : 
     231         854 :   fdata->p = data->p;
     232         854 :   t = leafcopy(topx); setvarn(t, data->v);
     233         854 :   fdata->top  = t;
     234         854 :   fdata->topr = FpX_red(t, data->pr);
     235             : 
     236         854 :   if (data->f == 1) z = gen_0;
     237             :   else
     238             :   { /* Compute a root of T in K(top) using Hensel's lift */
     239         238 :     z = pol_x(data->v);
     240         238 :     p1 = FpX_deriv(data->T, data->p);
     241             :     /* First lift to a root mod p */
     242             :     for (;;) {
     243         560 :       p2 = FpX_FpXQ_eval(data->T, z, fdata->top, data->p);
     244         560 :       if (gcmp0(p2)) break;
     245         322 :       p3 = FpX_FpXQ_eval(p1, z, fdata->top, data->p);
     246         322 :       z  = FpX_sub(z, FpXQ_div(p2, p3, fdata->top, data->p), data->p);
     247         322 :     }
     248             :     /* Then a root mod p^r */
     249         238 :     z = ZpX_ZpXQ_liftroot(data->T, z, fdata->top, data->p, data->r);
     250             :   }
     251             : 
     252         854 :   fdata->z  = z;
     253         854 :   fdata->eis = eis;
     254         854 :   fdata->pi  = Fq_sub(pol_x(data->v), fdata->z,
     255             :                       FpX_red(fdata->top, data->p), data->p);
     256         854 :   ipi = RgXQ_inv(fdata->pi, fdata->top);
     257         854 :   ipi = Q_remove_denom(ipi, &dipi);
     258         854 :   Q = mulii(data->pr, data->p);
     259         854 :   cipi = Fp_inv(diviiexact(dipi, data->p), Q);
     260         854 :   fdata->ipi = FpX_Fp_mul(ipi, cipi, Q); /* p/pi mod p^(pr+1) */
     261             : 
     262             :   /* Last one is set to pi^e/p (so we compute pi^e with one extra precision) */
     263         854 :   p2 = mulii(data->pr, data->p);
     264         854 :   p1 = FpXQ_powers(fdata->pi, data->e, fdata->topr, p2);
     265         854 :   gel(p1, data->e+1) = ZX_Z_divexact(gel(p1, data->e+1), data->p);
     266         854 :   fdata->pik  = p1;
     267         854 : }
     268             : 
     269             : /* return pol*ipi/p (mod top, pp) if it has integral coefficient, NULL
     270             :    otherwise. The result is reduced mod top, pp */
     271             : static GEN
     272      433335 : DivideByPi(FAD_t *fdata, GEN pp, GEN ppp, GEN pol)
     273             : {
     274             :   GEN P;
     275             :   long i, l;
     276      433335 :   pari_sp av = avma;
     277             : 
     278             :   /* "pol" is a t_INT or an FpX: signe() works for both */
     279      433335 :   if (!signe(pol)) return pol;
     280             : 
     281      426300 :   P = Fq_mul(pol, fdata->ipi, fdata->top, ppp); /* FpX */
     282      426300 :   l  = lg(P);
     283     2176174 :   for (i = 2; i < l; i++)
     284             :   {
     285             :     GEN r;
     286     1857730 :     gel(P,i) = dvmdii(gel(P,i), fdata->p, &r);
     287     1857730 :     if (r != gen_0) { avma = av; return NULL; }
     288             :   }
     289      318444 :   return FpX_red(P, pp);
     290             : }
     291             : 
     292             : /* return pol# = pol~/pi^vpi(pol~) mod pp where pol~(x) = pol(pi.x + beta)
     293             :  * has coefficients in the field defined by top and pi is a prime element */
     294             : static GEN
     295       53928 : GetSharp(FAD_t *fdata, GEN pp, GEN ppp, GEN pol, GEN beta, long *pl)
     296             : {
     297             :   GEN p1, p2;
     298       53928 :   long i, v, l, d = degpol(pol);
     299       53928 :   pari_sp av = avma;
     300             : 
     301       53928 :   if (!gequal0(beta))
     302       41636 :     p1 = FqX_translate(pol, beta, fdata->topr, pp);
     303             :   else
     304       12292 :     p1 = shallowcopy(pol);
     305       53928 :   p2 = p1;
     306      183232 :   for (v = 0;; v++)
     307             :   {
     308      475916 :     for (i = 0; i <= v; i++)
     309             :     {
     310      346612 :       GEN c = gel(p2, i+2);
     311      346612 :       c = DivideByPi(fdata, pp, ppp, c);
     312      346612 :       if (!c) break;
     313      292684 :       gel(p2, i+2) = c;
     314             :     }
     315      183232 :     if (i <= v) break;
     316      129304 :     p1 = shallowcopy(p2);
     317      129304 :   }
     318       53928 :   if (!v) pari_err_BUG("GetSharp [no division]");
     319             : 
     320       86723 :   for (l = v; l >= 0; l--)
     321             :   {
     322       86723 :     GEN c = gel(p1, l+2);
     323       86723 :     c = DivideByPi(fdata, pp, ppp, c);
     324       86723 :     if (!c) { break; }
     325             :   }
     326             : 
     327       53928 :   *pl = l; if (l < 2) return NULL;
     328             : 
     329             :   /* adjust powers */
     330       41244 :   for (i = v+1; i <= d; i++)
     331       23380 :     gel(p1, i+2) = Fq_mul(gel(p1, i+2),
     332       11690 :                           gel(fdata->pik, i-v+1), fdata->topr, pp);
     333             : 
     334       29554 :   return gerepilecopy(av, normalizepol(p1));
     335             : }
     336             : 
     337             : /* Compute roots of pol in the residue field. Use table look-up if possible */
     338             : static GEN
     339       53907 : Quick_FqX_roots(KRASNER_t *data, GEN pol)
     340             : {
     341             :   GEN rts;
     342       53907 :   ulong ind = 0;
     343             : 
     344       53907 :   if (data->f == 1)
     345       33068 :     pol = FpXY_evalx(pol, gen_0, data->p);
     346             :   else
     347       20839 :     pol = FqX_red(pol, data->T, data->p);
     348       53907 :   if (data->roottable)
     349             :   {
     350       53907 :     ind = ZXY_z_eval(pol, data->q[2], data->p[2]);
     351       53907 :     if (gel(data->roottable, ind)) return gel(data->roottable, ind);
     352             :   }
     353        2891 :   rts = FqX_roots(pol, data->T, data->p);
     354        2891 :   if (ind) gel(data->roottable, ind) = gclone(rts);
     355        2891 :   return rts;
     356             : }
     357             : 
     358             : static void
     359         133 : FreeRootTable(GEN T)
     360             : {
     361         133 :   if (T)
     362             :   {
     363         133 :     long j, l = lg(T);
     364      590373 :     for (j = 1; j < l; j++)
     365      590240 :       if (gel(T,j)) gunclone(gel(T,j));
     366             :   }
     367         133 : }
     368             : 
     369             : /* Return the number of roots of pol congruent to alpha modulo pi working
     370             :    modulo pp (all roots if alpha is NULL); if flag is non-zero, return 1
     371             :    as soon as a root is found. (Note: ppp = pp*p for DivideByPi) */
     372             : static long
     373       78281 : RootCongruents(KRASNER_t *data, FAD_t *fdata, GEN pol, GEN alpha, GEN pp, GEN ppp, long sl, long flag)
     374             : {
     375             :   GEN R;
     376             :   long s, i;
     377             : 
     378       78281 :   if (alpha)
     379             :   {
     380             :     long l;
     381       53928 :     pol = GetSharp(fdata, pp, ppp, pol, alpha, &l);
     382       53928 :     if (l <= 1) return l;
     383             :     /* decrease precision if sl gets bigger than a multiple of e */
     384       29554 :     sl += l;
     385       29554 :     if (sl >= data->e)
     386             :     {
     387       25970 :       sl -= data->e;
     388       25970 :       ppp = pp;
     389       25970 :       pp = diviiexact(pp, data->p);
     390             :     }
     391             :   }
     392             : 
     393       53907 :   R  = Quick_FqX_roots(data, pol);
     394      101108 :   for (s = 0, i = 1; i < lg(R); i++)
     395             :   {
     396       53928 :     s += RootCongruents(data, fdata, pol, gel(R, i), pp, ppp, sl, flag);
     397       53928 :     if (flag && s) return 1;
     398             :   }
     399       47180 :   return s;
     400             : }
     401             : 
     402             : /* pol is a ZXY defining a polynomial over the field defined by fdata
     403             :    If flag != 0, return 1 as soon as a root is found. Computations are done with
     404             :    a precision of pr. */
     405             : static long
     406       24353 : RootCountingAlgorithm(KRASNER_t *data, FAD_t *fdata, GEN pol, long flag)
     407             : {
     408             :   long j, l, d;
     409       24353 :   GEN P = cgetg_copy(pol, &l);
     410             : 
     411       24353 :   P[1] = pol[1];
     412       24353 :   d = l-3;
     413      109718 :   for (j = 0; j < d; j++)
     414             :   {
     415       85365 :     GEN cf = gel(pol, j+2);
     416       85365 :     if (typ(cf) == t_INT)
     417       53193 :       cf = diviiexact(cf, data->p);
     418             :     else
     419       32172 :       cf = ZX_Z_divexact(cf, data->p);
     420       85365 :     gel(P, j+2) = Fq_mul(cf, gel(fdata->pik, j+1), fdata->topr, data->pr);
     421             :   }
     422       24353 :   gel(P, d+2) = gel(fdata->pik, d+1); /* pik[d] = pi^d/p */
     423             : 
     424       24353 :   return RootCongruents(data, fdata, P, NULL, diviiexact(data->pr, data->p), data->pr, 0, flag);
     425             : }
     426             : 
     427             : /* Return non-zero if the field given by fdata defines a field isomorphic to
     428             :  * the one defined by pol */
     429             : static long
     430       15918 : IsIsomorphic(KRASNER_t *data, FAD_t *fdata, GEN pol)
     431             : {
     432             :   long j, nb;
     433       15918 :   pari_sp av = avma;
     434             : 
     435       15918 :   if (RgX_is_ZX(pol)) return RootCountingAlgorithm(data, fdata, pol, 1);
     436             : 
     437       15582 :   for (j = 1; j <= data->f; j++)
     438             :   {
     439       11725 :     GEN p1 = FqX_FpXQ_eval(pol, fdata->z, fdata->top, data->pr);
     440       11725 :     nb = RootCountingAlgorithm(data, fdata, p1, 1);
     441       11725 :     if (nb) { avma = av; return nb; }
     442       11032 :     if (j < data->f)
     443        7175 :       pol = FqX_FpXQ_eval(pol, data->frob, data->T, data->pr);
     444             :   }
     445        3857 :   avma = av; return 0;
     446             : }
     447             : 
     448             : /* Compute the number of conjugates fields of the field given by fdata */
     449             : static void
     450         854 : NbConjugateFields(KRASNER_t *data, FAD_t *fdata)
     451             : {
     452         854 :   GEN pol = fdata->eis;
     453             :   long j, nb;
     454         854 :   pari_sp av = avma;
     455             : 
     456         854 :   if (RgX_is_ZX(pol)) { /* split for efficiency; contains the case f = 1 */
     457         616 :     fdata->cj = data->e / RootCountingAlgorithm(data, fdata, pol, 0);
     458         616 :     avma = av; return;
     459             :   }
     460             : 
     461         238 :   nb = 0;
     462         882 :   for (j = 1; j <= data->f; j++)
     463             :   { /* Transform to pol. in z to pol. in a */
     464         644 :     GEN p1 = FqX_FpXQ_eval(pol, fdata->z, fdata->top, data->pr);
     465         644 :     nb += RootCountingAlgorithm(data, fdata, p1, 0);
     466             :     /* Look at the roots of conjugates polynomials */
     467         644 :     if (j < data->f)
     468         406 :       pol = FqX_FpXQ_eval(pol, data->frob, data->T, data->pr);
     469             :   }
     470         238 :   avma = av;
     471         238 :   fdata->cj = data->e * data->f / nb;
     472         238 :   return;
     473             : }
     474             : 
     475             : /* return a minimal list of polynomials generating all the totally
     476             :    ramified extensions of K^ur of degree e and discriminant
     477             :    p^{e + j - 1} in the tamely ramified case */
     478             : static GEN
     479        1582 : TamelyRamifiedCase(KRASNER_t *data)
     480             : {
     481        1582 :   long av = avma;
     482        1582 :   long g = ugcd(data->e, umodiu(data->qm1, data->e)); /* (e, q-1) */
     483        1582 :   GEN rep, eis, Xe = gpowgs(pol_x(0), data->e), m = stoi(data->e / g);
     484             : 
     485        1582 :   rep = zerovec(g);
     486        1582 :   eis = gadd(Xe, data->p);
     487        1582 :   gel(rep, 1) = mkvec2(get_topx(data,eis), m);
     488        1582 :   if (g > 1)
     489             :   {
     490         497 :     ulong pmodg = umodiu(data->p, g);
     491         497 :     long r = 1, ct = 1;
     492         497 :     GEN sv = InitSieve(g-1);
     493             :     /* let Frobenius act to get a minimal set of polynomials over Q_p */
     494        1771 :     while (r)
     495             :     {
     496             :       long gr;
     497        1554 :       GEN p1 = (typ(data->u) == t_INT)?
     498        1365 :         mulii(Fp_powu(data->u, r, data->p), data->p):
     499         588 :         ZX_Z_mul(FpXQ_powu(data->u, r, data->T, data->p), data->p);
     500         777 :       eis = gadd(Xe, p1); /* add cst coef */
     501         777 :       ct++;
     502         777 :       gel(rep, ct) = mkvec2(get_topx(data,eis), m);
     503         777 :       gr = r;
     504        1113 :       do { SetSieveValue(sv, gr); gr = Fl_mul(gr, pmodg, g); } while (gr != r);
     505         777 :       r  = NextZeroValue(sv, r);
     506             :     }
     507         497 :     setlg(rep, ct+1);
     508             :   }
     509        1582 :   return gerepilecopy(av, rep);
     510             : }
     511             : 
     512             : static long
     513         399 : function_l(GEN p, long a, long b, long i)
     514             : {
     515         399 :   long l = 1 + a - z_pval(i, p);
     516         399 :   if (i < b) l++;
     517         399 :   return (l < 1)? 1: l;
     518             : }
     519             : 
     520             : /* Structure of the coefficients set Omega. Each coefficient is [start, zr]
     521             :  * meaning all the numbers of the form:
     522             :  *   zeta_0 * p^start + ... + zeta_s * p^c (s = c - start)
     523             :  * with zeta_i roots of unity (powers of zq + zero), zeta_0 = 0 is
     524             :  * possible iff zr = 1 */
     525             : static GEN
     526         133 : StructureOmega(KRASNER_t *data, GEN *pnbpol)
     527             : {
     528         133 :   GEN nbpol, O = cgetg(data->e + 1, t_VEC);
     529             :   long i;
     530             : 
     531             :   /* i = 0 */
     532         133 :   gel(O, 1) = mkvecsmall2(1, 0);
     533         133 :   nbpol = mulii(data->qm1, powiu(data->q, data->c - 1));
     534         532 :   for (i = 1; i < data->e; i++)
     535             :   {
     536         399 :     long v_start, zero = 0;
     537             :     GEN nbcf, p1;
     538         399 :     v_start = function_l(data->p, data->a, data->b, i);
     539         399 :     p1 = powiu(data->q, data->c - v_start);
     540         399 :     if (i == data->b)
     541          98 :       nbcf = mulii(p1, data->qm1);
     542             :     else
     543             :     {
     544         301 :       nbcf = mulii(p1, data->q);
     545         301 :       zero = 1;
     546             :     }
     547         399 :     gel(O, i+1) = mkvecsmall2(v_start, zero);
     548         399 :     nbpol = mulii(nbpol, nbcf);
     549             :   }
     550         133 :   *pnbpol = nbpol; return O;
     551             : }
     552             : 
     553             : /* a random element of the finite field */
     554             : static GEN
     555       42602 : RandomFF(KRASNER_t *data)
     556             : {
     557       42602 :   long i, l = data->f + 2, p = itou(data->p);
     558       42602 :   GEN c = cgetg(l, t_POL);
     559       42602 :   c[1] = evalvarn(data->v);
     560       42602 :   for (i = 2; i < l; i++) gel(c, i) = utoi(random_Fl(p));
     561       42602 :   return ZX_renormalize(c, l);
     562             : }
     563             : static GEN
     564        2877 : RandomPol(KRASNER_t *data, GEN Omega)
     565             : {
     566        2877 :   long i, j, l = data->e + 3, end = data->c;
     567        2877 :   GEN pol = cgetg(l, t_POL);
     568        2877 :   pol[1] = evalsigne(1) | evalvarn(0);
     569       14189 :   for (i = 1; i <= data->e; i++)
     570             :   {
     571       11312 :     GEN c, cf = gel(Omega, i);
     572       11312 :     long st = cf[1], zr = cf[2];
     573             :     /* c = sum_{st <= j <= end} x_j p^j, where x_j are random Fq mod (p,upl)
     574             :      * If (!zr), insist on x_st != 0 */
     575             :     for (;;) {
     576       14154 :       c = RandomFF(data);
     577       14154 :       if (zr || signe(c)) break;
     578        2842 :     }
     579       39760 :     for (j = 1; j <= end-st; j++)
     580       28448 :       c = ZX_add(c, ZX_Z_mul(RandomFF(data), gel(data->pk, j)));
     581       11312 :     c = ZX_Z_mul(c, gel(data->pk, st));
     582       11312 :     c = FpX_red(c, data->pr);
     583       11312 :     switch(degpol(c))
     584             :     {
     585         735 :       case -1: c = gen_0; break;
     586        8022 :       case  0: c = gel(c,2); break;
     587             :     }
     588       11312 :     gel(pol, i+1) = c;
     589             :   }
     590        2877 :   gel(pol, i+1) = gen_1; /* monic */
     591        2877 :   return pol;
     592             : }
     593             : 
     594             : static void
     595         854 : CloneFieldData(FAD_t *fdata)
     596             : {
     597         854 :  fdata->top = gclone(fdata->top);
     598         854 :  fdata->topr= gclone(fdata->topr);
     599         854 :  fdata->z   = gclone(fdata->z);
     600         854 :  fdata->eis = gclone(fdata->eis);
     601         854 :  fdata->pi  = gclone(fdata->pi);
     602         854 :  fdata->ipi = gclone(fdata->ipi);
     603         854 :  fdata->pik = gclone(fdata->pik);
     604         854 : }
     605             : static void
     606         854 : FreeFieldData(FAD_t *fdata)
     607             : {
     608         854 :   gunclone(fdata->top);
     609         854 :   gunclone(fdata->topr);
     610         854 :   gunclone(fdata->z);
     611         854 :   gunclone(fdata->eis);
     612         854 :   gunclone(fdata->pi);
     613         854 :   gunclone(fdata->ipi);
     614         854 :   gunclone(fdata->pik);
     615         854 : }
     616             : 
     617             : static GEN
     618         133 : WildlyRamifiedCase(KRASNER_t *data)
     619             : {
     620         133 :   long nbext, ct, fd, nb = 0, j;
     621             :   GEN nbpol, rpl, rep, Omega;
     622             :   FAD_t **vfd;
     623             :   pari_timer T;
     624         133 :   pari_sp av = avma, av2;
     625             : 
     626             :   /* Omega = vector giving the structure of the set Omega */
     627             :   /* nbpol = number of polynomials to consider = |Omega| */
     628         133 :   Omega = StructureOmega(data, &nbpol);
     629         133 :   nbext = itos_or_0(data->nbext);
     630         133 :   if (!nbext || (nbext & ~LGBITS))
     631           0 :     pari_err_OVERFLOW("padicfields [too many extensions]");
     632             : 
     633         133 :   if (DEBUGLEVEL>0) {
     634           0 :     err_printf("There are %ld extensions to find and %Ps polynomials to consider\n", nbext, nbpol);
     635           0 :     timer_start(&T);
     636             :   }
     637             : 
     638         133 :   vfd = (FAD_t**)new_chunk(nbext);
     639         133 :   for (j = 0; j < nbext; j++) vfd[j] = (FAD_t*)new_chunk(sizeof(FAD_t));
     640             : 
     641         133 :   ct = 0;
     642         133 :   fd = 0;
     643         133 :   av2 = avma;
     644             : 
     645        3143 :   while (fd < nbext)
     646             :   { /* Jump randomly among the polynomials : seems best... */
     647        2877 :     rpl = RandomPol(data, Omega);
     648        2877 :     if (DEBUGLEVEL>3) err_printf("considering polynomial %Ps\n", rpl);
     649       16772 :     for (j = 0; j < ct; j++)
     650             :     {
     651       15918 :       nb = IsIsomorphic(data, vfd[j], rpl);
     652       15918 :       if (nb) break;
     653             :     }
     654        2877 :     if (!nb)
     655             :     {
     656         854 :       GEN topx = get_topx(data, rpl);
     657         854 :       FAD_t *f = (FAD_t*)vfd[ct];
     658         854 :       FieldData(data, f, rpl, topx);
     659         854 :       CloneFieldData(f);
     660         854 :       NbConjugateFields(data, f);
     661         854 :       nb = f->cj;
     662         854 :       fd += nb;
     663         854 :       ct++;
     664         854 :       if (DEBUGLEVEL>1)
     665           0 :         err_printf("%ld more extension%s\t(%ld/%ld, %ldms)\n",
     666             :                    nb, (nb == 1)? "": "s", fd, nbext, timer_delay(&T));
     667             :     }
     668        2877 :     avma = av2;
     669             :   }
     670             : 
     671         133 :   rep = cgetg(ct+1, t_VEC);
     672         987 :   for (j = 0; j < ct; j++)
     673             :   {
     674         854 :     FAD_t *f = (FAD_t*)vfd[j];
     675         854 :     GEN topx = ZX_copy(f->top);
     676         854 :     setvarn(topx, 0);
     677         854 :     gel(rep, j+1) = mkvec2(topx, stoi(f->cj));
     678         854 :     FreeFieldData(f);
     679             :   }
     680         133 :   FreeRootTable(data->roottable);
     681         133 :   return gerepileupto(av, rep);
     682             : }
     683             : 
     684             : /* return the minimal polynomial T of a generator of K^ur and the expression (mod pr)
     685             :  * in terms of T of a root of unity u such that u is l-maximal for all primes l
     686             :  * dividing g = (e,q-1). */
     687             : static void
     688        1715 : setUnramData(KRASNER_t *d)
     689             : {
     690        1715 :   if (d->f == 1)
     691             :   {
     692         560 :     d->T = pol_x(d->v);
     693         560 :     d->u = pgener_Fp(d->p);
     694         560 :     d->frob = pol_x(d->v);
     695             :   }
     696             :   else
     697             :   {
     698        1155 :     GEN L, z, T, p = d->p;
     699        1155 :     d->T = T = init_Fq(p, d->f, d->v);
     700        1155 :     L = gel(factoru(ugcd(d->e, umodiu(d->qm1, d->e))), 1);
     701        1155 :     z = gener_FpXQ_local(T, p, zv_to_ZV(L));
     702        1155 :     d->u = ZpXQ_sqrtnlift(pol_1(d->v), d->qm1, z, T, p, d->r);
     703        1155 :     d->frob = ZpX_ZpXQ_liftroot(T, FpXQ_pow(pol_x(d->v), p, T, p), T, p, d->r);
     704             :   }
     705        1715 : }
     706             : 
     707             : /* return [ p^1, p^2, ..., p^c ] */
     708             : static GEN
     709         133 : get_pk(KRASNER_t *data)
     710             : {
     711         133 :   long i, l = data->c + 1;
     712         133 :   GEN pk = cgetg(l, t_VEC), p = data->p;
     713         133 :   gel(pk, 1) = p;
     714         133 :   for (i = 2; i <= data->c; i++) gel(pk, i) = mulii(gel(pk, i-1), p);
     715         133 :   return pk;
     716             : }
     717             : 
     718             : /* Compute an absolute polynomial for all the totally ramified
     719             :    extensions of Q_p(z) of degree e and discriminant p^{e + j - 1}
     720             :    where z is a root of upl defining an unramified extension of Q_p */
     721             : /* See padicfields for the meaning of flag */
     722             : static GEN
     723        1743 : GetRamifiedPol(GEN p, GEN efj, long flag)
     724             : {
     725        1743 :   long e = efj[1], f = efj[2], j = efj[3], i, l;
     726        1743 :   const long v = 1;
     727             :   GEN L, pols;
     728             :   KRASNER_t data;
     729        1743 :   pari_sp av = avma;
     730             : 
     731        1743 :   if (DEBUGLEVEL>1)
     732           0 :     err_printf("  Computing extensions with decomposition [e, f, j] = [%ld, %ld, %ld]\n", e,f,j);
     733        1743 :   data.p   = p;
     734        1743 :   data.e   = e;
     735        1743 :   data.f   = f;
     736        1743 :   data.j   = j;
     737        1743 :   data.a   = j/e;
     738        1743 :   data.b   = j%e;
     739        1743 :   data.c   = (e+2*j)/e+1;
     740        1743 :   data.q   = powiu(p, f);
     741        1743 :   data.qm1 = subiu(data.q, 1);
     742        1743 :   data.v   = v;
     743        1743 :   data.r   = 1 + (long)ceildivuu(2*j + 3, e); /* enough precision */
     744        1743 :   data.pr  = powiu(p, data.r);
     745        1743 :   data.nbext = NumberExtensions(&data);
     746             : 
     747        1743 :   if (flag == 2) return data.nbext;
     748             : 
     749        1715 :   setUnramData(&data);
     750        1715 :   if (DEBUGLEVEL>1)
     751           0 :     err_printf("  Unramified part %Ps\n", data.T);
     752        1715 :   data.roottable = NULL;
     753        1715 :   if (j)
     754             :   {
     755         133 :     if (lgefint(data.q) == 3)
     756             :     {
     757         133 :       ulong npol = upowuu(data.q[2], e+1);
     758         133 :       if (npol && npol < (1<<19)) data.roottable = const_vec(npol, NULL);
     759             :     }
     760         133 :     data.pk = get_pk(&data);
     761         133 :     L = WildlyRamifiedCase(&data);
     762             :   }
     763             :   else
     764        1582 :     L = TamelyRamifiedCase(&data);
     765             : 
     766        1715 :   pols = cgetg_copy(L, &l);
     767        1715 :   if (flag == 1)
     768             :   {
     769        1631 :     GEN E = utoipos(e), F = utoipos(f), D = utoi(f*(e+j-1));
     770        4578 :     for (i = 1; i < l; i++)
     771             :     {
     772        2947 :       GEN T = gel(L,i);
     773        2947 :       gel(pols, i) = mkvec5(ZX_copy(gel(T, 1)), E,F,D, icopy(gel(T, 2)));
     774             :     }
     775             :   }
     776             :   else
     777             :   {
     778         350 :     for (i = 1; i < l; i++)
     779             :     {
     780         266 :       GEN T = gel(L,i);
     781         266 :       gel(pols, i) = ZX_copy(gel(T,1));
     782             :     }
     783             :   }
     784        1715 :   return gerepileupto(av, pols);
     785             : }
     786             : 
     787             : static GEN
     788         483 : possible_efj(GEN p, long m)
     789             : { /* maximal possible discriminant valuation d <= m * (1+v_p(m)) - 1 */
     790             :   /* 1) [j = 0, tame] d = m - f with f | m and v_p(f) = v_p(m), or
     791             :    * 2) [j > 0, wild] d >= m, j <= v_p(e)*e   */
     792         483 :   ulong m1, pve, pp = p[2]; /* pp used only if v > 0 */
     793         483 :   long ve, v = u_pvalrem(m, p, &m1);
     794         483 :   GEN L, D = divisorsu(m1);
     795         483 :   long i, taum1 = lg(D)-1, nb = 0;
     796             : 
     797         483 :   if (v) {
     798          21 :     for (pve = 1,ve = 1; ve <= v; ve++) { pve *= pp; nb += pve * ve; }
     799          21 :     nb = itos_or_0(muluu(nb, zv_sum(D)));
     800          21 :     if (!nb || is_bigint( mului(pve, sqru(v)) ) )
     801           0 :       pari_err_OVERFLOW("padicfields [too many ramification possibilities]");
     802             :   }
     803         483 :   nb += taum1; /* upper bound for the number of possible triples [e,f,j] */
     804             : 
     805         483 :   L = cgetg(nb + 1, t_VEC);
     806             :   /* 1) tame */
     807        2093 :   for (nb=1, i=1; i < lg(D); i++)
     808             :   {
     809        1610 :     long e = D[i], f = m / e;
     810        1610 :     gel(L, nb++) = mkvecsmall3(e, f, 0);
     811             :   }
     812             :   /* 2) wild */
     813             :   /* Ore's condition: either
     814             :    * 1) j = v_p(e) * e, or
     815             :    * 2) j = a e + b, with 0 < b < e and v_p(b) <= a < v_p(e) */
     816         511 :   for (pve = 1, ve = 1; ve <= v; ve++)
     817             :   {
     818          28 :     pve *= pp; /* = p^ve */
     819          56 :     for (i = 1; i < lg(D); i++)
     820             :     {
     821          28 :       long a,b, e = D[i] * pve, f = m / e;
     822          98 :       for (b = 1; b < e; b++)
     823         154 :         for (a = u_lval(b, pp); a < ve; a++)
     824          84 :           gel(L, nb++) = mkvecsmall3(e, f,  a*e + b);
     825          28 :       gel(L, nb++) = mkvecsmall3(e, f, ve*e);
     826             :     }
     827             :   }
     828         483 :   setlg(L, nb); return L;
     829             : }
     830             : 
     831             : #ifdef CHECK
     832             : static void
     833             : checkpols(GEN p, GEN EFJ, GEN pols)
     834             : {
     835             :   GEN pol, p1, p2;
     836             :   long c1, c2, e, f, j, i, l = lg(pols);
     837             : 
     838             :   if (typ(pols) == t_INT) return;
     839             : 
     840             :   e = EFJ[1];
     841             :   f = EFJ[2];
     842             :   j = EFJ[3];
     843             : 
     844             :   for (i = 1; i < l; i++)
     845             :   {
     846             :     pol = gel(pols, i);
     847             :     if (typ(pol) == t_VEC) pol = gel(pol, 1);
     848             :     if (!isirreducible(pol)) pari_err_BUG("Polynomial is reducible");
     849             :     p1 = poldisc0(pol, -1);
     850             :     if (gvaluation(p1, p) != f*(e+j-1)) pari_err_BUG("Discriminant is incorrect");
     851             :     /* only compute a p-maximal order */
     852             :     p1 = nfinitall(mkvec2(pol, mkvec(p)), 0, DEFAULTPREC);
     853             :     p2 = idealprimedec(p1, p);
     854             :     if(lg(p2) > 2) pari_err_BUG("Prime p is split");
     855             :     p2 = gel(p2, 1);
     856             :     if (cmpis(gel(p2, 3), e)) pari_err_BUG("Ramification index is incorrect");
     857             :     if (cmpis(gel(p2, 4), f)) pari_err_BUG("inertia degree is incorrect");
     858             :   }
     859             : 
     860             :   if (l == 2) return;
     861             :   if (e*f > 20) return;
     862             : 
     863             :   /* TODO: check that (random) distinct polynomials give non-isomorphic extensions */
     864             :   for (i = 1; i < 3*l; i++)
     865             :   {
     866             :     c1 = random_Fl(l-1)+1;
     867             :     c2 = random_Fl(l-1)+1;
     868             :     if (c1 == c2) continue;
     869             :     p1 = gel(pols, c1);
     870             :     if (typ(p1) == t_VEC) p1 = gel(p1, 1);
     871             :     p2 = gel(pols, c2);
     872             :     if (typ(p2) == t_VEC) p2 = gel(p2, 1);
     873             :     pol = polcompositum0(p1, p2, 0);
     874             :     pol = gel(pol, 1);
     875             :     if (poldegree(pol, -1) > 100) continue;
     876             :     p1 = factorpadic(pol, p, 2);
     877             :     p1 = gmael(p1, 1, 1);
     878             :     if (poldegree(p1, -1) == e*f) pari_err_BUG("fields are isomorphic");
     879             :     /*
     880             :       p1 = nfinitall(mkvec2(pol, mkvec(p)), 0, DEFAULTPREC);
     881             :       p2 = idealprimedec_galois(p1, p);
     882             :       if (!cmpis(mulii(gel(p2, 3), gel(p2, 4)), e*f)) pari_err_BUG("fields are isomorphic");
     883             :     */
     884             :   }
     885             : }
     886             : #endif
     887             : 
     888             : static GEN
     889         490 : pols_from_efj(pari_sp av, GEN EFJ, GEN p, long flag)
     890             : {
     891             :   long i, l;
     892         490 :   GEN L = cgetg_copy(EFJ, &l);
     893         490 :   if (l == 1) { avma = av; return flag == 2? gen_0: cgetg(1, t_VEC); }
     894        2233 :   for (i = 1; i < l; i++)
     895             :   {
     896        1743 :     gel(L,i) = GetRamifiedPol(p, gel(EFJ,i), flag);
     897             : #ifdef CHECK
     898             :     checkpols(p, gel(EFJ, i), gel(L, i));
     899             : #endif
     900             :   }
     901         490 :   if (flag == 2) return gerepileuptoint(av, ZV_sum(L));
     902         483 :   return gerepilecopy(av, shallowconcat1(L));
     903             : }
     904             : 
     905             : /* return a minimal list of polynomials generating all the extensions of
     906             :    Q_p of given degree N; if N is a t_VEC [n,d], return extensions of degree n
     907             :    and discriminant p^d. */
     908             : /* Return only the polynomials if flag = 0 (default), also the ramification
     909             :    degree, the residual degree and the discriminant if flag = 1 and only the
     910             :    number of extensions if flag = 2 */
     911             : GEN
     912         490 : padicfields0(GEN p, GEN N, long flag)
     913             : {
     914         490 :   pari_sp av = avma;
     915         490 :   long m = 0, d = -1;
     916             :   GEN L;
     917             : 
     918         490 :   if (typ(p) != t_INT) pari_err_TYPE("padicfields",p);
     919             :   /* be nice to silly users */
     920         490 :   if (!BPSW_psp(p)) pari_err_PRIME("padicfields",p);
     921         490 :   switch(typ(N))
     922             :   {
     923             :     case t_VEC:
     924           7 :       if (lg(N) != 3) pari_err_TYPE("padicfields",N);
     925           7 :       d = gtos(gel(N,2));
     926           7 :       N = gel(N,1); /* fall through */
     927             :     case t_INT:
     928         490 :       m = itos(N);
     929         490 :       if (m <= 0) pari_err_DOMAIN("padicfields", "degree", "<=", gen_0,N);
     930         490 :       break;
     931             :     default:
     932           0 :       pari_err_TYPE("padicfields",N);
     933             :   }
     934         490 :   if (d >= 0) return padicfields(p, m, d, flag);
     935         483 :   L = possible_efj(p, m);
     936         483 :   return pols_from_efj(av, L, p, flag);
     937             : }
     938             : 
     939             : GEN
     940           7 : padicfields(GEN p, long m, long d, long flag)
     941             : {
     942           7 :   long av = avma;
     943           7 :   GEN L = possible_efj_by_d(p, m, d);
     944           7 :   return pols_from_efj(av, L, p, flag);
     945             : }

Generated by: LCOV version 1.11