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 - char.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.10.0 lcov report (development 21059-cbe0d6a) Lines: 774 797 97.1 %
Date: 2017-09-22 06:24:58 Functions: 63 63 100.0 %
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             : #include "pari.h"
      14             : #include "paripriv.h"
      15             : 
      16             : /*********************************************************************/
      17             : /**                                                                 **/
      18             : /**                  GENERIC ABELIAN CHARACTERS                     **/
      19             : /**                                                                 **/
      20             : /*********************************************************************/
      21             : /* check whether G is a znstar */
      22             : int
      23     1463959 : checkznstar_i(GEN G)
      24             : {
      25     4391681 :   return (typ(G) == t_VEC && lg(G) == 6
      26     1463763 :       && typ(znstar_get_faN(G)) == t_VEC
      27     2927722 :       && typ(gel(G,1)) == t_VEC && lg(gel(G,1)) == 3);
      28             : }
      29             : 
      30             : int
      31       40229 : char_check(GEN cyc, GEN chi)
      32       40229 : { return typ(chi) == t_VEC && lg(chi) == lg(cyc) && RgV_is_ZV(chi); }
      33             : 
      34             : /* Shallow; return [ d[1],  d[1]/d[2],...,d[1]/d[n] ] */
      35             : GEN
      36       65779 : cyc_normalize(GEN d)
      37             : {
      38       65779 :   long i, l = lg(d);
      39             :   GEN C, D;
      40       65779 :   if (l == 1) return mkvec(gen_1);
      41       65765 :   D = cgetg(l, t_VEC); gel(D,1) = C = gel(d,1);
      42       65765 :   for (i = 2; i < l; i++) gel(D,i) = diviiexact(C, gel(d,i));
      43       65765 :   return D;
      44             : }
      45             : 
      46             : /* chi character [D,C] given by chi(g_i) = \zeta_D^C[i] for all i, return
      47             :  * [d,c] such that chi(g_i) = \zeta_d^c[i] for all i and d minimal */
      48             : GEN
      49       67977 : char_simplify(GEN D, GEN C)
      50             : {
      51       67977 :   GEN d = D;
      52       67977 :   if (lg(C) == 1) d = gen_1;
      53             :   else
      54             :   {
      55       67599 :     GEN t = gcdii(d, ZV_content(C));
      56       67599 :     if (!equali1(t))
      57             :     {
      58       44163 :       C = ZC_Z_divexact(C, t);
      59       44163 :       d = diviiexact(d, t);
      60             :     }
      61             :   }
      62       67977 :   return mkvec2(d,C);
      63             : }
      64             : 
      65             : /* Shallow; ncyc from cyc_normalize(): ncyc[1] = cyc[1],
      66             :  * ncyc[i] = cyc[i]/cyc[1] for i > 1; chi character on G ~ cyc.
      67             :  * Return [d,c] such that: chi( g_i ) = e(chi[i] / cyc[i]) = e(c[i]/ d) */
      68             : GEN
      69       67256 : char_normalize(GEN chi, GEN ncyc)
      70             : {
      71       67256 :   long i, l = lg(chi);
      72       67256 :   GEN c = cgetg(l, t_VEC);
      73       67256 :   if (l > 1) {
      74       67242 :     gel(c,1) = gel(chi,1);
      75       67242 :     for (i = 2; i < l; i++) gel(c,i) = mulii(gel(chi,i), gel(ncyc,i));
      76             :   }
      77       67256 :   return char_simplify(gel(ncyc,1), c);
      78             : }
      79             : 
      80             : /* Called by function 's'. x is a group object affording ".cyc" method, and
      81             :  * chi an abelian character. Return NULL if the group is (Z/nZ)^* [special
      82             :  * case more character types allowed] and x.cyc otherwise */
      83             : static GEN
      84         483 : get_cyc(GEN x, GEN chi, const char *s)
      85             : {
      86         483 :   if (nftyp(x) == typ_BIDZ)
      87             :   {
      88         448 :     if (!zncharcheck(x, chi)) pari_err_TYPE(s, chi);
      89         448 :     return NULL;
      90             :   }
      91             :   else
      92             :   {
      93          35 :     if (typ(x) != t_VEC || !RgV_is_ZV(x)) x = member_cyc(x);
      94          35 :     if (!char_check(x, chi)) pari_err_TYPE(s, chi);
      95          35 :     return x;
      96             :   }
      97             : }
      98             : 
      99             : /* conjugate character [ZV/ZC] */
     100             : GEN
     101        3094 : charconj(GEN cyc, GEN chi)
     102             : {
     103             :   long i, l;
     104        3094 :   GEN z = cgetg_copy(chi, &l);
     105        8232 :   for (i = 1; i < l; i++)
     106             :   {
     107        5138 :     GEN c = gel(chi,i);
     108        5138 :     gel(z,i) = signe(c)? subii(gel(cyc,i), c): gen_0;
     109             :   }
     110        3094 :   return z;
     111             : }
     112             : GEN
     113          28 : charconj0(GEN x, GEN chi)
     114             : {
     115          28 :   GEN cyc = get_cyc(x, chi, "charconj");
     116          28 :   return cyc? charconj(cyc, chi): zncharconj(x, chi);
     117             : }
     118             : 
     119             : GEN
     120       22400 : charorder(GEN cyc, GEN x)
     121             : {
     122       22400 :   pari_sp av = avma;
     123       22400 :   long i, l = lg(cyc);
     124       22400 :   GEN f = gen_1;
     125       65268 :   for (i = 1; i < l; i++)
     126       42868 :     if (signe(gel(x,i)))
     127             :     {
     128       31794 :       GEN c, o = gel(cyc,i);
     129       31794 :       c = gcdii(o, gel(x,i));
     130       31794 :       if (!is_pm1(c)) o = diviiexact(o,c);
     131       31794 :       f = lcmii(f, o);
     132             :     }
     133       22400 :   return gerepileuptoint(av, f);
     134             : }
     135             : GEN
     136         168 : charorder0(GEN x, GEN chi)
     137             : {
     138         168 :   GEN cyc = get_cyc(x, chi, "charorder");
     139         168 :   return cyc? charorder(cyc, chi): zncharorder(x, chi);
     140             : }
     141             : 
     142             : /* chi character of abelian G: chi[i] = chi(z_i), where G = \oplus Z/cyc[i] z_i.
     143             :  * Return Ker chi */
     144             : GEN
     145       39886 : charker(GEN cyc, GEN chi)
     146             : {
     147       39886 :   long i, l = lg(cyc);
     148             :   GEN nchi, ncyc, m, U;
     149             : 
     150       39886 :   if (l == 1) return cgetg(1,t_MAT); /* trivial subgroup */
     151       39879 :   ncyc = cyc_normalize(cyc);
     152       39879 :   nchi = char_normalize(chi, ncyc);
     153       39879 :   m = shallowconcat(gel(nchi,2), gel(nchi,1));
     154       39879 :   U = gel(ZV_extgcd(m), 2); setlg(U,l);
     155       39879 :   for (i = 1; i < l; i++) setlg(U[i], l);
     156       39879 :   return hnfmodid(U, gel(ncyc,1));
     157             : }
     158             : GEN
     159          28 : charker0(GEN x, GEN chi)
     160             : {
     161          28 :   GEN cyc = get_cyc(x, chi, "charker");
     162          28 :   return cyc? charker(cyc, chi): zncharker(x, chi);
     163             : }
     164             : 
     165             : GEN
     166          63 : charpow(GEN cyc, GEN a, GEN N)
     167             : {
     168             :   long i, l;
     169          63 :   GEN v = cgetg_copy(a, &l);
     170          63 :   for (i = 1; i < l; i++) gel(v,i) = Fp_mul(gel(a,i), N, gel(cyc,i));
     171          63 :   return v;
     172             : }
     173             : GEN
     174        2387 : charmul(GEN cyc, GEN a, GEN b)
     175             : {
     176             :   long i, l;
     177        2387 :   GEN v = cgetg_copy(a, &l);
     178        2387 :   for (i = 1; i < l; i++) gel(v,i) = Fp_add(gel(a,i), gel(b,i), gel(cyc,i));
     179        2387 :   return v;
     180             : }
     181             : GEN
     182        8393 : chardiv(GEN cyc, GEN a, GEN b)
     183             : {
     184             :   long i, l;
     185        8393 :   GEN v = cgetg_copy(a, &l);
     186        8393 :   for (i = 1; i < l; i++) gel(v,i) = Fp_sub(gel(a,i), gel(b,i), gel(cyc,i));
     187        8393 :   return v;
     188             : }
     189             : GEN
     190          63 : charpow0(GEN x, GEN a, GEN N)
     191             : {
     192          63 :   GEN cyc = get_cyc(x, a, "charpow");
     193          63 :   return cyc? charpow(cyc, a, N): zncharpow(x, a, N);
     194             : }
     195             : GEN
     196         154 : charmul0(GEN x, GEN a, GEN b)
     197             : {
     198         154 :   const char *s = "charmul";
     199         154 :   GEN cyc = get_cyc(x, a, s);
     200         154 :   if (!cyc)
     201             :   {
     202         154 :     if (!zncharcheck(x, b)) pari_err_TYPE(s, b);
     203         154 :     return zncharmul(x, a, b);
     204             :   }
     205             :   else
     206             :   {
     207           0 :     if (!char_check(cyc, b)) pari_err_TYPE(s, b);
     208           0 :     return charmul(cyc, a, b);
     209             :   }
     210             : }
     211             : GEN
     212          42 : chardiv0(GEN x, GEN a, GEN b)
     213             : {
     214          42 :   const char *s = "chardiv";
     215          42 :   GEN cyc = get_cyc(x, a, s);
     216          42 :   if (!cyc)
     217             :   {
     218          42 :     if (!zncharcheck(x, b)) pari_err_TYPE(s, b);
     219          42 :     return znchardiv(x, a, b);
     220             :   }
     221             :   else
     222             :   {
     223           0 :     if (!char_check(cyc, b)) pari_err_TYPE(s, b);
     224           0 :     return chardiv(cyc, a, b);
     225             :   }
     226             : }
     227             : 
     228             : static GEN
     229     1143933 : chareval_i(GEN nchi, GEN dlog, GEN z)
     230             : {
     231     1143933 :   GEN o, q, r, b = gel(nchi,1);
     232     1143933 :   GEN a = FpV_dotproduct(gel(nchi,2), dlog, b);
     233             :   /* image is a/b in Q/Z */
     234     1143933 :   if (!z) return gdiv(a,b);
     235     1143555 :   if (typ(z) == t_INT)
     236             :   {
     237     1142183 :     q = dvmdii(z, b, &r);
     238     1142183 :     if (signe(r)) pari_err_TYPE("chareval", z);
     239     1142183 :     return mulii(a, q);
     240             :   }
     241             :   /* return z^(a*o/b), assuming z^o = 1 and b | o */
     242        1372 :   if (typ(z) != t_VEC || lg(z) != 3) pari_err_TYPE("chareval", z);
     243        1372 :   o = gel(z,2); if (typ(o) != t_INT) pari_err_TYPE("chareval", z);
     244        1372 :   q = dvmdii(o, b, &r); if (signe(r)) pari_err_TYPE("chareval", z);
     245        1372 :   q = mulii(a, q); /* in [0, o[ since a is reduced mod b */
     246        1372 :   z = gel(z,1);
     247        1372 :   if (typ(z) == t_VEC)
     248             :   {
     249          35 :     if (itos_or_0(o) != lg(z)-1) pari_err_TYPE("chareval", z);
     250          35 :     return gcopy(gel(z, itos(q)+1));
     251             :   }
     252             :   else
     253        1337 :     return gpow(z, q, DEFAULTPREC);
     254             : }
     255             : 
     256             : static GEN
     257        1876 : not_coprime(GEN z)
     258        1876 : { return (!z || typ(z) == t_INT)? gen_m1: gen_0; }
     259             : 
     260             : static GEN
     261          35 : get_chi(GEN cyc, GEN chi)
     262             : {
     263          35 :   if (!char_check(cyc,chi)) pari_err_TYPE("chareval", chi);
     264          35 :   return char_normalize(chi, cyc_normalize(cyc));
     265             : }
     266             : /* G a bnr.  FIXME: horribly inefficient to check that (x,N)=1, what to do ? */
     267             : static int
     268          42 : bnr_coprime(GEN G, GEN x)
     269             : {
     270          42 :   GEN t, N = gel(bnr_get_mod(G), 1);
     271          42 :   if (typ(x) == t_INT) /* shortcut */
     272             :   {
     273          14 :     t = gcdii(gcoeff(N,1,1), x);
     274          14 :     if (equali1(t)) return 1;
     275           0 :     t = idealadd(G, N, x);
     276           0 :     return equali1(gcoeff(t,1,1));
     277             :   }
     278          28 :   x = idealnumden(G, x);
     279          28 :   t = idealadd(G, N, gel(x,1));
     280          28 :   if (!equali1(gcoeff(t,1,1))) return 0;
     281          21 :   t = idealadd(G, N, gel(x,2));
     282          21 :   return equali1(gcoeff(t,1,1));
     283             : }
     284             : GEN
     285        3598 : chareval(GEN G, GEN chi, GEN x, GEN z)
     286             : {
     287        3598 :   pari_sp av = avma;
     288             :   GEN nchi, L;
     289             : 
     290        3598 :   switch(nftyp(G))
     291             :   {
     292             :     case typ_BNR:
     293          42 :       if (!bnr_coprime(G, x)) return not_coprime(z);
     294          28 :       L = isprincipalray(G, x);
     295          28 :       nchi = get_chi(bnr_get_cyc(G), chi);
     296          28 :       break;
     297             :     case typ_BNF:
     298           7 :       L = isprincipal(G, x);
     299           7 :       nchi = get_chi(bnf_get_cyc(G), chi);
     300           7 :       break;
     301             :     case typ_BIDZ:
     302        3542 :       if (checkznstar_i(G)) return gerepileupto(av, znchareval(G, chi, x, z));
     303             :       /* don't implement chars on general bid: need an nf... */
     304             :     default:
     305           7 :       pari_err_TYPE("chareval", G);
     306             :       return NULL;/* LCOV_EXCL_LINE */
     307             :   }
     308          35 :   return gerepileupto(av, chareval_i(nchi, L, z));
     309             : }
     310             : 
     311             : /* nchi = [ord,D] a quasi-normalized character (ord may be a multiple of
     312             :  * the character order); return v such that v[n] = -1 if (n,N) > 1 else
     313             :  * chi(n) = e(v[n]/ord), 1 <= n <= N */
     314             : GEN
     315         294 : ncharvecexpo(GEN G, GEN nchi)
     316             : {
     317         294 :   long N = itou(znstar_get_N(G)), ord = itou(gel(nchi,1)), i, j, l;
     318             :   GEN cyc, gen, d, t, t1, t2, t3, e, u, u1, u2, u3;
     319         294 :   GEN D = gel(nchi,2), v = const_vecsmall(N,-1);
     320         294 :   pari_sp av = avma;
     321         294 :   if (typ(D) == t_COL) {
     322         294 :     cyc = znstar_get_conreycyc(G);
     323         294 :     gen = znstar_get_conreygen(G);
     324             :   } else {
     325           0 :     cyc = znstar_get_cyc(G);
     326           0 :     gen = znstar_get_gen(G);
     327             :   }
     328         294 :   l = lg(cyc);
     329         294 :   e = u = cgetg(N+1,t_VECSMALL);
     330         294 :   d = t = cgetg(N+1,t_VECSMALL);
     331         294 :   *++d = 1;
     332         294 :   *++e = 0; v[*d] = *e;
     333         679 :   for (i = 1; i < l; i++)
     334             :   {
     335         385 :     ulong g = itou(gel(gen,i)), c = itou(gel(cyc,i)), x = itou(gel(D,i));
     336        8736 :     for (t1=t,u1=u,j=c-1; j; j--,t1=t2,u1=u2)
     337       27860 :       for (t2=d,u2=e, t3=t1,u3=u1; t3<t2; )
     338             :       {
     339       11158 :         *++d = Fl_mul(*++t3, g, N);
     340       11158 :         *++e = Fl_add(*++u3, x, ord); v[*d] = *e;
     341             :       }
     342             :   }
     343         294 :   avma = av; return v;
     344             : }
     345             : 
     346             : /*****************************************************************************/
     347             : 
     348             : static ulong
     349       24234 : lcmuu(ulong a, ulong b) { return (a/ugcd(a,b)) * b; }
     350             : static ulong
     351       11242 : zv_charorder(GEN cyc, GEN x)
     352             : {
     353       11242 :   long i, l = lg(cyc);
     354       11242 :   ulong f = 1;
     355       27482 :   for (i = 1; i < l; i++)
     356       16240 :     if (x[i])
     357             :     {
     358       12992 :       ulong o = cyc[i];
     359       12992 :       f = lcmuu(f, o / ugcd(o, x[i]));
     360             :     }
     361       11242 :   return f;
     362             : }
     363             : 
     364             : /* N > 0 */
     365             : GEN
     366       50190 : coprimes_zv(ulong N)
     367             : {
     368       50190 :   GEN v = cgetg(N+1, t_VECSMALL);
     369             :   ulong i;
     370       50190 :   v[1] = 1; for (i = 2; i <= N; i++) v[i] = (ugcd(N,i)==1);
     371       50190 :   return v;
     372             : }
     373             : /* cf zv_cyc_minimal: return k such that g*k is minimal (wrt lex) */
     374             : long
     375       39942 : zv_cyc_minimize(GEN cyc, GEN g, GEN coprime)
     376             : {
     377       39942 :   pari_sp av = avma;
     378       39942 :   long d, k, e, i, k0, bestk, l = lg(g), o = lg(coprime)-1;
     379             :   GEN best, gk, gd;
     380             :   ulong t;
     381       39942 :   if (o == 1) return 1;
     382       46284 :   for (i = 1; i < l; i++)
     383       46284 :     if (g[i]) break;
     384       39935 :   if (g[i] == 1) return 1;
     385       27650 :   k0 = Fl_invgen(g[i], cyc[i], &t);
     386       27650 :   d = cyc[i] / (long)t;
     387       27650 :   if (k0 > 1) g = vecmoduu(Flv_Fl_mul(g, k0, cyc[i]), cyc);
     388       35469 :   for (i++; i < l; i++)
     389       27986 :     if (g[i]) break;
     390       27650 :   if (i == l || d >= cyc[i]) return k0;
     391        2681 :   cyc = vecslice(cyc,i,l-1);
     392        2681 :   g   = vecslice(g,  i,l-1);
     393        2681 :   e = cyc[1];
     394        2681 :   gd = Flv_Fl_mul(g, d, e);
     395        2681 :   bestk = 1; best = g;
     396        6818 :   for (gk = g, k = d+1; k < e; k += d)
     397             :   {
     398        4137 :     long ko = k % o;
     399        4137 :     gk = Flv_add(gk, gd, e); if (!ko || !coprime[ko]) continue;
     400        3136 :     gk = vecmoduu(gk, cyc);
     401        3136 :     if (vecsmall_lexcmp(gk, best) < 0) { best = gk; bestk = k; }
     402             :   }
     403        2681 :   avma = av; return bestk == 1? k0: Fl_mul(k0, bestk, o);
     404             : }
     405             : /* g of order o in abelian group G attached to cyc. Is g a minimal generator
     406             :  * [wrt lex order] of the cyclic subgroup it generates;
     407             :  * coprime = coprimes_zv(o) */
     408             : long
     409       11214 : zv_cyc_minimal(GEN cyc, GEN g, GEN coprime)
     410             : {
     411       11214 :   pari_sp av = avma;
     412       11214 :   long d, k, e, l = lg(g), o = lg(coprime)-1; /* elt order */
     413             :   GEN gd, gk;
     414       11214 :   if (o == 1) return 1;
     415       11746 :   for (k = 1; k < l; k++)
     416       11746 :     if (g[k]) break;
     417       11214 :   if (g[k] == 1) return 1;
     418        8134 :   if (cyc[k] % g[k]) return 0;
     419        8134 :   d = cyc[k] / g[k]; /* > 1 */
     420        9149 :   for (k++; k < l; k++) /* skip following 0s */
     421        9149 :     if (g[k]) break;
     422        8134 :   if (k == l || d >= cyc[k]) return 1;
     423        1288 :   cyc = vecslice(cyc,k,l-1);
     424        1288 :   g   = vecslice(g,  k,l-1);
     425        1288 :   e = cyc[1];
     426             :   /* find k in (Z/e)^* such that g*k mod cyc is lexicographically minimal,
     427             :    * k = 1 mod d to fix the first non-zero entry */
     428        1288 :   gd = Flv_Fl_mul(g, d, e);
     429        2716 :   for (gk = g, k = d+1; k < e; k += d)
     430             :   {
     431        1743 :     long ko = k % o;
     432        1743 :     gk = Flv_add(gk, gd, e); if (!ko || !coprime[ko]) continue;
     433        1106 :     gk = vecmoduu(gk, cyc);
     434        1106 :     if (vecsmall_lexcmp(gk, g) < 0) { avma = av; return 0; }
     435             :   }
     436         973 :   avma = av; return 1;
     437             : }
     438             : 
     439             : static GEN
     440        1750 : coprime_tables(long N)
     441             : {
     442        1750 :   GEN D = divisorsu(N), v = const_vec(N, NULL);
     443        1750 :   long i, l = lg(D);
     444        1750 :   for (i = 1; i < l; i++) gel(v, D[i]) = coprimes_zv(D[i]);
     445        1750 :   return v;
     446             : }
     447             : /* enumerate all group elements, modulo (Z/cyc[1])^* */
     448             : static GEN
     449        1750 : cyc2elts_normal(GEN cyc, long maxord, GEN ORD)
     450             : {
     451        1750 :   long i, n, o, N, j = 1;
     452             :   GEN z, vcoprime;
     453             : 
     454        1750 :   if (typ(cyc) != t_VECSMALL) cyc = gtovecsmall(cyc);
     455        1750 :   n = lg(cyc)-1;
     456        1750 :   if (n == 0) return cgetg(1, t_VEC);
     457        1750 :   N = zv_prod(cyc);
     458        1750 :   z = cgetg(N+1, t_VEC);
     459        1750 :   if (1 <= maxord && (!ORD|| zv_search(ORD,1)))
     460        1393 :     gel(z,j++) = zero_zv(n);
     461        1750 :   vcoprime = coprime_tables(cyc[1]);
     462        5166 :   for (i = n; i > 0; i--)
     463             :   {
     464        3416 :     GEN cyc0 = vecslice(cyc,i+1,n), pre = zero_zv(i);
     465        3416 :     GEN D = divisorsu(cyc[i]), C = cyc2elts(cyc0);
     466        3416 :     long s, t, lD = lg(D), nC = lg(C)-1; /* remove last element */
     467       13811 :     for (s = 1; s < lD-1; s++)
     468             :     {
     469       10395 :       long o0 = D[lD-s]; /* cyc[i] / D[s] */
     470       10395 :       if (o0 > maxord) continue;
     471        8981 :       pre[i] = D[s];
     472        8981 :       if (!ORD || zv_search(ORD,o0))
     473             :       {
     474        8834 :         GEN c = vecsmall_concat(pre, zero_zv(n-i));
     475        8834 :         gel(z,j++) = c;
     476             :       }
     477       20223 :       for (t = 1; t < nC; t++)
     478             :       {
     479       11242 :         GEN chi0 = gel(C,t);
     480       11242 :         o = lcmuu(o0, zv_charorder(cyc0,chi0));
     481       11242 :         if (o <= maxord && (!ORD || zv_search(ORD,o)))
     482             :         {
     483       11214 :           GEN c = vecsmall_concat(pre, chi0);
     484       11214 :           if (zv_cyc_minimal(cyc, c, gel(vcoprime,o))) gel(z,j++) = c;
     485             :         }
     486             :       }
     487             :     }
     488             :   }
     489        1750 :   setlg(z,j); return z;
     490             : }
     491             : 
     492             : GEN
     493        1764 : chargalois(GEN G, GEN ORD)
     494             : {
     495        1764 :   pari_sp av = avma;
     496             :   long maxord, i, l;
     497        1764 :   GEN v, cyc = (typ(G) == t_VEC && RgV_is_ZVpos(G))? G: member_cyc(G);
     498        1764 :   if (lg(cyc) == 1) retmkvec(cgetg(1,t_VEC));
     499        1750 :   maxord = itou(gel(cyc,1));
     500        1750 :   if (ORD && gequal0(ORD)) ORD = NULL;
     501        1750 :   if (ORD)
     502         392 :     switch(typ(ORD))
     503             :     {
     504             :       long l;
     505             :       case t_VEC:
     506           7 :         ORD = ZV_to_zv(ORD);
     507             :       case t_VECSMALL:
     508         364 :         ORD = leafcopy(ORD);
     509         364 :         vecsmall_sort(ORD);
     510         364 :         l = lg(ORD);
     511         364 :         if (l == 1) return cgetg(1, t_VECSMALL);
     512         364 :         maxord = minss(maxord, ORD[l-1]);
     513         364 :         break;
     514             :       case t_INT:
     515          28 :         maxord = minss(maxord, itos(ORD));
     516          28 :         ORD = NULL;
     517          28 :         break;
     518           0 :       default: pari_err_TYPE("chargalois", ORD);
     519             :     }
     520        1750 :   v = cyc2elts_normal(cyc, maxord, ORD); l = lg(v);
     521        1750 :   for(i = 1; i < l; i++) gel(v,i) = zv_to_ZV(gel(v,i));
     522        1750 :   return gerepileupto(av, v);
     523             : }
     524             : 
     525             : /*********************************************************************/
     526             : /**                                                                 **/
     527             : /**                  (Z/NZ)^* AND DIRICHLET CHARACTERS              **/
     528             : /**                                                                 **/
     529             : /*********************************************************************/
     530             : 
     531             : GEN
     532       58541 : znstar0(GEN N, long flag)
     533             : {
     534       58541 :   GEN F = NULL, P, E, cyc, gen, mod, G;
     535             :   long i, i0, l, nbprimes;
     536       58541 :   pari_sp av = avma;
     537             : 
     538       58541 :   if (flag && flag != 1) pari_err_FLAG("znstar");
     539       58541 :   if ((F = check_arith_all(N,"znstar")))
     540             :   {
     541       15043 :     F = clean_Z_factor(F);
     542       15043 :     N = typ(N) == t_VEC? gel(N,1): factorback(F);
     543             :   }
     544       58541 :   if (!signe(N))
     545             :   {
     546          21 :     if (flag) pari_err_IMPL("znstar(0,1)");
     547          14 :     avma = av;
     548          14 :     retmkvec3(gen_2, mkvec(gen_2), mkvec(gen_m1));
     549             :   }
     550       58520 :   if (signe(N) < 0) N = absi(N);
     551       58520 :   if (abscmpiu(N,2) <= 0)
     552             :   {
     553       10192 :     G = mkvec3(gen_1, cgetg(1,t_VEC), cgetg(1,t_VEC));
     554       10192 :     if (flag)
     555             :     {
     556       10164 :       GEN v = const_vec(6,cgetg(1,t_VEC));
     557       10164 :       gel(v,3) = cgetg(1,t_MAT);
     558       30142 :       F = equali1(N)? mkvec2(cgetg(1,t_COL),cgetg(1,t_VECSMALL))
     559       19978 :                     : mkvec2(mkcol(gen_2), mkvecsmall(1));
     560       10164 :       G = mkvec5(mkvec2(N,mkvec(gen_0)), G, F, v, cgetg(1,t_MAT));
     561             :     }
     562       10192 :     return gerepilecopy(av,G);
     563             :   }
     564       48328 :   if (!F) F = Z_factor(N);
     565       48328 :   P = gel(F,1); nbprimes = lg(P)-1;
     566       48328 :   E = ZV_to_nv( gel(F,2) );
     567       48328 :   switch(mod8(N))
     568             :   {
     569             :     case 0:
     570       11354 :       P = shallowconcat(gen_2,P);
     571       11354 :       E = vecsmall_prepend(E, E[1]); /* add a copy of p=2 row */
     572       11354 :       i = 2; /* 2 generators at 2 */
     573       11354 :       break;
     574             :     case 4:
     575        8687 :       i = 1; /* 1 generator at 2 */
     576        8687 :       break;
     577             :     case 2: case 6:
     578        6342 :       P = vecsplice(P,1);
     579        6342 :       E = vecsplice(E,1); /* remove 2 */
     580        6342 :       i = 0; /* no generator at 2 */
     581        6342 :       break;
     582             :     default:
     583       21945 :       i = 0; /* no generator at 2 */
     584       21945 :       break;
     585             :   }
     586       48328 :   l = lg(P);
     587       48328 :   cyc = cgetg(l,t_VEC);
     588       48328 :   gen = cgetg(l,t_VEC);
     589       48328 :   mod = cgetg(l,t_VEC);
     590             :   /* treat p=2 first */
     591       48328 :   if (i == 2)
     592             :   {
     593       11354 :     long v2 = E[1];
     594       11354 :     GEN q = int2n(v2);
     595       11354 :     gel(cyc,1) = gen_2;
     596       11354 :     gel(gen,1) = subiu(q,1); /* -1 */
     597       11354 :     gel(mod,1) = q;
     598       11354 :     gel(cyc,2) = int2n(v2-2);
     599       11354 :     gel(gen,2) = utoipos(5); /* Conrey normalization */
     600       11354 :     gel(mod,2) = q;
     601       11354 :     i0 = 3;
     602             :   }
     603       36974 :   else if (i == 1)
     604             :   {
     605        8687 :     gel(cyc,1) = gen_2;
     606        8687 :     gel(gen,1) = utoipos(3);
     607        8687 :     gel(mod,1) = utoipos(4);
     608        8687 :     i0 = 2;
     609             :   }
     610             :   else
     611       28287 :     i0 = 1;
     612             :   /* odd primes, fill remaining entries */
     613      118881 :   for (i = i0; i < l; i++)
     614             :   {
     615       70553 :     long e = E[i];
     616       70553 :     GEN p = gel(P,i), q = powiu(p, e-1), Q = mulii(p, q);
     617       70553 :     gel(cyc,i) = subii(Q, q); /* phi(p^e) */
     618       70553 :     gel(gen,i) = pgener_Zp(p);/* Conrey normalization, for e = 1 also */
     619       70553 :     gel(mod,i) = Q;
     620             :   }
     621             :   /* gen[i] has order cyc[i] and generates (Z/mod[i]Z)^* */
     622       48328 :   if (nbprimes > 1) /* lift generators to (Z/NZ)^*, = 1 mod N/mod[i] */
     623      119539 :     for (i=1; i<l; i++)
     624             :     {
     625       85764 :       GEN Q = gel(mod,i), g = gel(gen,i), qinv = Fp_inv(Q, diviiexact(N,Q));
     626       85764 :       g = addii(g, mulii(mulii(subsi(1,g),qinv),Q));
     627       85764 :       gel(gen,i) = modii(g, N);
     628             :     }
     629             : 
     630             :   /* cyc[i] > 1 and remain so in the loop, gen[i] = 1 mod (N/mod[i]) */
     631       48328 :   if (!flag)
     632             :   { /* update generators in place; about twice faster */
     633        4221 :     G = gen;
     634        4487 :     for (i=l-1; i>=2; i--)
     635             :     {
     636         266 :       GEN ci = gel(cyc,i), gi = gel(G,i);
     637             :       long j;
     638         602 :       for (j=i-1; j>=1; j--) /* we want cyc[i] | cyc[j] */
     639             :       {
     640         336 :         GEN cj = gel(cyc,j), gj, qj, v, d;
     641             : 
     642         336 :         d = bezout(ci,cj,NULL,&v); /* > 1 */
     643         651 :         if (absequalii(ci, d)) continue; /* ci | cj */
     644         147 :         if (absequalii(cj, d)) { /* cj | ci */
     645         126 :           swap(gel(G,j),gel(G,i));
     646         126 :           gi = gel(G,i);
     647         126 :           swap(gel(cyc,j),gel(cyc,i));
     648         126 :           ci = gel(cyc,i); continue;
     649             :         }
     650             : 
     651          21 :         qj = diviiexact(cj,d);
     652          21 :         gel(cyc,j) = mulii(ci,qj);
     653          21 :         gel(cyc,i) = d;
     654             : 
     655             :         /* [1,v*cj/d; 0,1]*[1,0;-1,1]*diag(cj,ci)*[ci/d,-v; cj/d,u]
     656             :          * = diag(lcm,gcd), with u ci + v cj = d */
     657          21 :         gj = gel(G,j);
     658             :         /* (gj, gi) *= [1,0; -1,1]^-1 */
     659          21 :         gj = Fp_mul(gj, gi, N); /* order ci*qj = lcm(ci,cj) */
     660             :         /* (gj,gi) *= [1,v*qj; 0,1]^-1 */
     661          21 :         togglesign_safe(&v);
     662          21 :         if (signe(v) < 0) v = modii(v,ci); /* >= 0 to avoid inversions */
     663          21 :         gel(G,i) = gi = Fp_mul(gi, Fp_pow(gj, mulii(qj, v), N), N);
     664          21 :         gel(G,j) = gj;
     665          21 :         ci = d; if (absequaliu(ci, 2)) break;
     666             :       }
     667             :     }
     668        4221 :     G = mkvec3(ZV_prod(cyc), cyc, FpV_to_mod(G,N));
     669             :   }
     670             :   else
     671             :   { /* keep matrices between generators, return an 'init' structure */
     672       44107 :     GEN D, U, Ui, fao = cgetg(l, t_VEC), lo = cgetg(l, t_VEC);
     673       44107 :     F = mkvec2(P, E);
     674       44107 :     D = ZV_snf_group(cyc,&U,&Ui);
     675      141568 :     for (i = 1; i < l; i++)
     676             :     {
     677       97461 :       GEN t = gen_0, p = gel(P,i), p_1 = subiu(p,1);
     678       97461 :       long e = E[i];
     679       97461 :       gel(fao,i) = get_arith_ZZM(p_1);
     680       97461 :       if (e >= 2 && !absequaliu(p,2))
     681             :       {
     682        7595 :         GEN q = gel(mod,i), g = Fp_pow(gel(gen,i),p_1,q);
     683        7595 :         if (e == 2)
     684        6216 :           t = Fp_inv(diviiexact(subiu(g,1), p), p);
     685             :         else
     686        1379 :           t = ginv(Qp_log(cvtop(g,p,e)));
     687             :       }
     688       97461 :       gel(lo,i) = t;
     689             :     }
     690       44107 :     G = cgetg(l, t_VEC);
     691       44107 :     for (i = 1; i < l; i++) gel(G,i) = FpV_factorback(gen, gel(Ui,i), N);
     692       44107 :     G = mkvec3(ZV_prod(D), D, G);
     693       44107 :     G = mkvec5(mkvec2(N,mkvec(gen_0)), G, F,
     694             :                mkvecn(6,mod,fao,Ui,gen,cyc,lo), U);
     695             :   }
     696       48328 :   return gerepilecopy(av, G);
     697             : }
     698             : GEN
     699        4193 : znstar(GEN N) { return znstar0(N, 0); }
     700             : 
     701             : /* g has order 2^(e-2), g,h = 1 (mod 4); return x s.t. g^x = h (mod 2^e) */
     702             : static GEN
     703      289219 : Zideallog_2k(GEN h, GEN g, long e, GEN pe)
     704             : {
     705      289219 :   GEN a = Fp_log(h, g, int2n(e-2), pe);
     706      289219 :   if (typ(a) != t_INT) return NULL;
     707      289219 :   return a;
     708             : }
     709             : 
     710             : /* ord = get_arith_ZZM(p-1), simplified form of znlog_rec: g is known
     711             :  * to be a primitive root mod p^e; lo = 1/log_p(g^(p-1)) */
     712             : static GEN
     713     1104551 : Zideallog_pk(GEN h, GEN g, GEN p, long e, GEN pe, GEN ord, GEN lo)
     714             : {
     715     1104551 :   GEN gp = (e == 1)? g: modii(g, p);
     716     1104551 :   GEN hp = (e == 1)? h: modii(h, p);
     717     1104551 :   GEN a = Fp_log(hp, gp, ord, p);
     718     1104551 :   if (typ(a) != t_INT) return NULL;
     719     1104544 :   if (e > 1)
     720             :   { /* find a s.t. g^a = h (mod p^e), p odd prime, e > 0, (h,p) = 1 */
     721             :     /* use p-adic log: O(log p + e) mul*/
     722       17871 :     GEN b, p_1 = gel(ord,1);
     723       17871 :     h = Fp_mul(h, Fp_pow(g, negi(a), pe), pe);
     724             :     /* g,h = 1 mod p; compute b s.t. h = g^b */
     725       17871 :     if (e == 2) /* simpler */
     726       11081 :       b = Fp_mul(diviiexact(subiu(h,1), p), lo, p);
     727             :     else
     728        6790 :       b = padic_to_Q(gmul(Qp_log(cvtop(h, p, e)), lo));
     729       17871 :     a = addii(a, mulii(p_1, b));
     730             :   }
     731     1104544 :   return a;
     732             : }
     733             : 
     734             : int
     735      211974 : znconrey_check(GEN cyc, GEN chi)
     736      211974 : { return typ(chi) == t_COL && lg(chi) == lg(cyc) && RgV_is_ZV(chi); }
     737             : 
     738             : int
     739       94164 : zncharcheck(GEN G, GEN chi)
     740             : {
     741       94164 :   switch(typ(chi))
     742             :   {
     743         875 :     case t_INT: return 1;
     744       92995 :     case t_COL: return znconrey_check(znstar_get_conreycyc(G), chi);
     745         294 :     case t_VEC: return char_check(znstar_get_cyc(G), chi);
     746             :   }
     747           0 :   return 0;
     748             : }
     749             : 
     750             : GEN
     751       21630 : znconreyfromchar_normalized(GEN bid, GEN chi)
     752             : {
     753       21630 :   GEN nchi, U = znstar_get_U(bid);
     754       21630 :   long l = lg(chi);
     755       21630 :   if (l == 1) retmkvec2(gen_1,cgetg(1,t_VEC));
     756       21602 :   if (!RgV_is_ZV(chi) || lgcols(U) != l) pari_err_TYPE("lfunchiZ", chi);
     757       21602 :   nchi = char_normalize(chi, cyc_normalize(znstar_get_cyc(bid)));
     758       21602 :   gel(nchi,2) = ZV_ZM_mul(gel(nchi,2),U); return nchi;
     759             : }
     760             : 
     761             : GEN
     762       21539 : znconreyfromchar(GEN bid, GEN chi)
     763             : {
     764       21539 :   GEN nchi = znconreyfromchar_normalized(bid, chi);
     765       21539 :   GEN v = char_denormalize(znstar_get_conreycyc(bid), gel(nchi,1), gel(nchi,2));
     766       21539 :   settyp(v, t_COL); return v;
     767             : }
     768             : 
     769             : /* discrete log on canonical "primitive root" generators
     770             :  * Allow log(x) instead of x [usual discrete log on bid's generators */
     771             : GEN
     772     1182188 : znconreylog(GEN bid, GEN x)
     773             : {
     774     1182188 :   pari_sp av = avma;
     775             :   GEN N, L, F, P,E, y, pe, fao, gen, lo, cycg;
     776             :   long i, l;
     777     1182188 :   if (!checkznstar_i(bid)) pari_err_TYPE("znconreylog", bid);
     778     1182188 :   N = znstar_get_N(bid);
     779     1182188 :   if (typ(N) != t_INT) pari_err_TYPE("znconreylog", N);
     780     1182181 :   if (cmpiu(N, 2) <= 0) return cgetg(1, t_COL);
     781     1181985 :   cycg = znstar_get_conreycyc(bid);
     782     1181985 :   switch(typ(x))
     783             :   {
     784             :     GEN Ui;
     785             :     case t_INT:
     786     1181012 :       if (!signe(x)) pari_err_COPRIME("znconreylog", x, N);
     787     1181005 :       break;
     788             :     case t_COL: /* log_bid(x) */
     789          35 :       Ui = znstar_get_Ui(bid);
     790          35 :       if (!RgV_is_ZV(x) || lg(x) != lg(Ui)) pari_err_TYPE("znconreylog", x);
     791          35 :       return gerepileupto(av, vecmodii(ZM_ZC_mul(Ui,x), cycg));
     792             :     case t_VEC:
     793         938 :       return gerepilecopy(av, znconreyfromchar(bid, x));
     794           0 :     default: pari_err_TYPE("znconreylog", x);
     795             :   }
     796     1181005 :   F = znstar_get_faN(bid); /* factor(N) */
     797     1181005 :   P = gel(F, 1); /* prime divisors of N */
     798     1181005 :   E = gel(F, 2); /* exponents */
     799     1181005 :   L = gel(bid,4);
     800     1181005 :   pe = znstar_get_pe(bid);
     801     1181005 :   fao = gel(L,2);
     802     1181005 :   gen = znstar_get_conreygen(bid); /* local generators of (Z/p^k)^* */
     803     1181005 :   lo = gel(L,6); /* 1/log_p((g_i)^(p_i-1)) */
     804             : 
     805     1181005 :   l = lg(gen); i = 1;
     806     1181005 :   y = cgetg(l, t_COL);
     807     1181005 :   if (!mod2(N) && !mod2(x)) pari_err_COPRIME("znconreylog", x, N);
     808     1180991 :   if (absequaliu(gel(P,1), 2) && E[1] >= 2)
     809             :   {
     810      545538 :     if (E[1] == 2)
     811      256319 :       gel(y,i++) = mod4(x) == 1? gen_0: gen_1;
     812             :     else
     813             :     {
     814      289219 :       GEN a, x2, q2 = gel(pe,1);
     815      289219 :       x2 = modii(x, q2);
     816      289219 :       if (mod4(x) == 1) /* 1 or 5 mod 8*/
     817      183974 :         gel(y,i++) = gen_0;
     818             :       else /* 3 or 7 */
     819      105245 :       { gel(y,i++) = gen_1; x2 = subii(q2, x2); }
     820             :       /* x2 = 5^x mod q */
     821      289219 :       a = Zideallog_2k(x2, gel(gen,i), E[1], q2);
     822      289219 :       if (!a) pari_err_COPRIME("znconreylog", x, N);
     823      289219 :       gel(y, i++) = a;
     824             :     }
     825             :   }
     826     3466526 :   while (i < l)
     827             :   {
     828     1104551 :     GEN p = gel(P,i), q = gel(pe,i), xpe = modii(x, q);
     829     1104551 :     GEN a = Zideallog_pk(xpe, gel(gen,i), p, E[i], q, gel(fao,i), gel(lo,i));
     830     1104551 :     if (!a) pari_err_COPRIME("znconreylog", x, N);
     831     1104544 :     gel(y, i++) = a;
     832             :   }
     833     1180984 :   return gerepilecopy(av, y);
     834             : }
     835             : GEN
     836        6720 : Zideallog(GEN bid, GEN x)
     837             : {
     838        6720 :   pari_sp av = avma;
     839        6720 :   GEN y = znconreylog(bid, x), U = znstar_get_U(bid);
     840        6692 :   return gerepileupto(av, ZM_ZC_mul(U, y));
     841             : }
     842             : GEN
     843         280 : znlog0(GEN h, GEN g, GEN o)
     844             : {
     845         280 :   if (typ(g) == t_VEC)
     846             :   {
     847             :     GEN N;
     848          56 :     if (o) pari_err_TYPE("znlog [with znstar]", o);
     849          56 :     if (!checkznstar_i(g)) pari_err_TYPE("znlog", h);
     850          56 :     N = znstar_get_N(g);
     851          56 :     h = Rg_to_Fp(h,N);
     852          49 :     return Zideallog(g, h);
     853             :   }
     854         224 :   return znlog(h, g, o);
     855             : }
     856             : 
     857             : GEN
     858       42518 : znconreyexp(GEN bid, GEN x)
     859             : {
     860       42518 :   pari_sp av = avma;
     861             :   long i, l;
     862             :   GEN N, pe, gen, cycg, v, vmod;
     863             :   int e2;
     864       42518 :   if (!checkznstar_i(bid)) pari_err_TYPE("znconreyexp", bid);
     865       42518 :   cycg = znstar_get_conreycyc(bid);
     866       42518 :   switch(typ(x))
     867             :   {
     868             :     case t_VEC:
     869          21 :       x = znconreylog(bid, x);
     870          21 :       break;
     871             :     case t_COL:
     872       42497 :       if (RgV_is_ZV(x) && lg(x) == lg(cycg)) break;
     873           7 :     default: pari_err_TYPE("znconreyexp",x);
     874             :   }
     875       42511 :   pe = znstar_get_pe(bid);
     876       42511 :   gen = znstar_get_conreygen(bid); /* local generators of (Z/p^k)^* */
     877       42511 :   cycg = znstar_get_conreycyc(bid);
     878       42511 :   l = lg(x); v = cgetg(l, t_VEC);
     879       42511 :   N = znstar_get_N(bid);
     880       42511 :   e2 = !mod8(N); /* 2 generators at p = 2 */
     881      136703 :   for (i = 1; i < l; i++)
     882             :   {
     883             :     GEN q, g, m;
     884       94192 :     if (i == 1 && e2) { gel(v,1) = NULL; continue; }
     885       84161 :     q = gel(pe,i);
     886       84161 :     g = gel(gen,i);
     887       84161 :     m = modii(gel(x,i), gel(cycg,i));
     888       84161 :     m = Fp_pow(g, m, q);
     889       84161 :     if (i == 2 && e2 && signe(gel(x,1))) m = Fp_neg(m, q);
     890       84161 :     gel(v,i) = mkintmod(m, q);
     891             :   }
     892       42511 :   if (e2) v = vecsplice(v, 1);
     893       42511 :   v = chinese1_coprime_Z(v);
     894       42511 :   vmod = gel(v,1);
     895       42511 :   v = gel(v,2);
     896       42511 :   if (mpodd(v) || mpodd(N)) return gerepilecopy(av, v);
     897             :   /* handle N = 2 mod 4 */
     898        2450 :   return gerepileuptoint(av, addii(v, vmod));
     899             : }
     900             : 
     901             : /* Return Dirichlet character \chi_q(m,.), where bid = znstar(q);
     902             :  * m is either a t_INT, or a t_COL [Conrey logarithm] */
     903             : GEN
     904       40075 : znconreychar(GEN bid, GEN m)
     905             : {
     906       40075 :   pari_sp av = avma;
     907             :   GEN c, d, nchi;
     908             : 
     909       40075 :   switch(typ(m))
     910             :   {
     911             :     case t_COL:
     912             :     case t_INT:
     913       40075 :       nchi = znconrey_normalized(bid,m); /* images of primroot gens */
     914       40068 :       break;
     915             :     default:
     916           0 :       pari_err_TYPE("znconreychar",m);
     917             :       return NULL;/*LCOV_EXCL_LINE*/
     918             :   }
     919       40068 :   d = gel(nchi,1);
     920       40068 :   c = ZV_ZM_mul(gel(nchi,2), znstar_get_Ui(bid)); /* images of bid gens */
     921       40068 :   return gerepilecopy(av, char_denormalize(znstar_get_cyc(bid),d,c));
     922             : }
     923             : 
     924             : /* chi a t_INT or Conrey log describing a character. Return conductor, as an
     925             :  * integer if primitive; as a t_VEC [N,factor(N)] if not. Set *pm=m to the
     926             :  * attached primitive character: chi(g_i) = m[i]/ord(g_i)
     927             :  * Caller should use znconreylog_normalize(BID, m), once BID(conductor) is
     928             :  * computed (wasteful to do it here since BID is shared by many characters) */
     929             : GEN
     930      119672 : znconreyconductor(GEN bid, GEN chi, GEN *pm)
     931             : {
     932      119672 :   pari_sp av = avma;
     933             :   GEN q, m, F, P, E;
     934             :   long i, j, l;
     935      119672 :   int e2, primitive = 1;
     936             : 
     937      119672 :   if (!checkznstar_i(bid)) pari_err_TYPE("znconreyconductor", bid);
     938      119672 :   if (typ(chi) == t_COL)
     939             :   {
     940      118979 :     if (!znconrey_check(znstar_get_conreycyc(bid), chi))
     941           0 :       pari_err_TYPE("znconreyconductor",chi);
     942             :   }
     943             :   else
     944         693 :     chi = znconreylog(bid, chi);
     945      119672 :   l = lg(chi);
     946      119672 :   F = znstar_get_faN(bid);
     947      119672 :   P = gel(F,1);
     948      119672 :   E = gel(F,2);
     949      119672 :   if (l == 1)
     950             :   {
     951       70721 :     avma = av;
     952       70721 :     if (pm) *pm = cgetg(1,t_COL);
     953       70721 :     if (lg(P) == 1) return gen_1;
     954          14 :     retmkvec2(gen_1, trivial_fact());
     955             :   }
     956       48951 :   P = leafcopy(P);
     957       48951 :   E = leafcopy(E);
     958       48951 :   m = cgetg(l, t_COL);
     959       48951 :   e2 = (E[1] >= 3 && absequaliu(gel(P,1),2));
     960       48951 :   i = j = 1;
     961       48951 :   if (e2)
     962             :   { /* two generators at p=2 */
     963       15232 :     GEN a1 = gel(chi,1), a = gel(chi,2);
     964       15232 :     i = 3;
     965       15232 :     if (!signe(a))
     966             :     {
     967        9429 :       e2 =  primitive = 0;
     968        9429 :       if (signe(a1))
     969             :       { /* lose one generator */
     970        1939 :         E[1] = 2;
     971        1939 :         gel(m,1) = a1;
     972        1939 :         j = 2;
     973             :       }
     974             :       /* else lose both */
     975             :     }
     976             :     else
     977             :     {
     978        5803 :       long v = Z_pvalrem(a, gen_2, &a);
     979        5803 :       if (v) { E[1] -= v; E[2] = E[1]; primitive = 0; }
     980        5803 :       gel(m,1) = a1;
     981        5803 :       gel(m,2) = a;
     982        5803 :       j = 3;
     983             :     }
     984             :   }
     985       48951 :   l = lg(P);
     986      119385 :   for (; i < l; i++)
     987             :   {
     988       70434 :     GEN p = gel(P,i), a = gel(chi,i);
     989             :     /* image of g_i in Q/Z is a/cycg[i], cycg[i] = order(g_i) */
     990       70434 :     if (!signe(a)) primitive = 0;
     991             :     else
     992             :     {
     993       49133 :       long v = Z_pvalrem(a, p, &a);
     994       49133 :       E[j] = E[i]; if (v) { E[j] -= v; primitive = 0; }
     995       49133 :       gel(P,j) = gel(P,i);
     996       49133 :       gel(m,j) = a; j++;
     997             :     }
     998             :   }
     999       48951 :   setlg(m,j);
    1000       48951 :   setlg(P,j);
    1001       48951 :   setlg(E,j);
    1002       48951 :   if (pm) *pm = m; /* attached primitive  character */
    1003       48951 :   if (primitive)
    1004             :   {
    1005       26705 :     q = znstar_get_N(bid);
    1006       26705 :     if (mod4(q) == 2) primitive = 0;
    1007             :   }
    1008       48951 :   if (!primitive)
    1009             :   {
    1010       26159 :     if (e2)
    1011             :     { /* remove duplicate p=2 row from factorization */
    1012        2702 :       P = vecsplice(P,1);
    1013        2702 :       E = vecsplice(E,1);
    1014             :     }
    1015       26159 :     E = zv_to_ZV(E);
    1016       26159 :     q = mkvec2(factorback2(P,E), mkmat2(P,E));
    1017             :   }
    1018       48951 :   gerepileall(av, pm? 2: 1, &q, pm);
    1019       48951 :   return q;
    1020             : }
    1021             : 
    1022             : GEN
    1023       18676 : zncharinduce(GEN G, GEN chi, GEN N)
    1024             : {
    1025       18676 :   pari_sp av = avma;
    1026             :   GEN q, faq, P, E, Pq, Eq, CHI;
    1027             :   long i, j, l;
    1028             :   int e2;
    1029             : 
    1030       18676 :   if (!checkznstar_i(G)) pari_err_TYPE("zncharinduce", G);
    1031       18676 :   if (!zncharcheck(G, chi)) pari_err_TYPE("zncharinduce", chi);
    1032       18676 :   q = znstar_get_N(G);
    1033       18676 :   if (typ(chi) != t_COL) chi = znconreylog(G, chi);
    1034       18676 :   if (checkznstar_i(N))
    1035             :   {
    1036       18487 :     GEN faN = znstar_get_faN(N);
    1037       18487 :     P = gel(faN,1); l = lg(P);
    1038       18487 :     E = gel(faN,2);
    1039       18487 :     N = znstar_get_N(N);
    1040       18487 :     if (l > 2 && equalii(gel(P,1),gel(P,2)))
    1041             :     { /* remove duplicate 2 */
    1042        4690 :       l--;
    1043        4690 :       P = vecsplice(P,1);
    1044        4690 :       E = vecsplice(E,1);
    1045             :     }
    1046             :   }
    1047             :   else
    1048             :   {
    1049         189 :     GEN faN = check_arith_pos(N, "zncharinduce");
    1050         189 :     if (!faN) faN = Z_factor(N);
    1051             :     else
    1052           0 :       N = (typ(N) == t_VEC)? gel(N,1): factorback(faN);
    1053         189 :     P = gel(faN,1);
    1054         189 :     E = gel(faN,2);
    1055             :   }
    1056       18676 :   if (!dvdii(N,q)) pari_err_DOMAIN("zncharinduce", "N % q", "!=", gen_0, N);
    1057       18669 :   if (mod4(N) == 2)
    1058             :   { /* remove 2 */
    1059        2226 :     if (lg(P) > 1 && absequaliu(gel(P,1), 2))
    1060             :     {
    1061          28 :       P = vecsplice(P,1);
    1062          28 :       E = vecsplice(E,1);
    1063             :     }
    1064        2226 :     N = shifti(N,-1);
    1065             :   }
    1066       18669 :   l = lg(P);
    1067             :   /* q = N or q = 2N, N odd */
    1068       18669 :   if (cmpii(N,q) <= 0) return gerepilecopy(av, chi);
    1069             :   /* N > 1 => l > 1*/
    1070       16849 :   if (typ(E) != t_VECSMALL) E = ZV_to_zv(E);
    1071       16849 :   e2 = (E[1] >= 3 && absequaliu(gel(P,1),2)); /* 2 generators at 2 mod N */
    1072       16849 :   if (ZV_equal0(chi))
    1073             :   {
    1074        1589 :     avma = av;
    1075        1589 :     return equali1(N)? cgetg(1, t_COL): zerocol(l+e2 - 1);
    1076             :   }
    1077             : 
    1078       15260 :   faq = znstar_get_faN(G);
    1079       15260 :   Pq = gel(faq,1);
    1080       15260 :   Eq = gel(faq,2);
    1081       15260 :   CHI = cgetg(l+e2, t_COL);
    1082       15260 :   i = j = 1;
    1083       15260 :   if (e2)
    1084             :   {
    1085        4550 :     i = 2; j = 3;
    1086        4550 :     if (absequaliu(gel(Pq,1), 2))
    1087             :     {
    1088        2779 :       if (Eq[1] >= 3)
    1089             :       { /* 2 generators at 2 mod q */
    1090        1428 :         gel(CHI,1) = gel(chi,1);
    1091        1428 :         gel(CHI,2) = shifti(gel(chi,2), E[1]-Eq[1]);
    1092             :       }
    1093        1351 :       else if (Eq[1] == 2)
    1094             :       { /* 1 generator at 2 mod q */
    1095        1351 :         gel(CHI,1) = gel(chi,1);
    1096        1351 :         gel(CHI,2) = gen_0;
    1097             :       }
    1098             :       else
    1099           0 :         gel(CHI,1) = gel(CHI,2) = gen_0;
    1100             :     }
    1101             :     else
    1102        1771 :       gel(CHI,1) = gel(CHI,2) = gen_0;
    1103             :   }
    1104       41531 :   for (; i < l; i++,j++)
    1105             :   {
    1106       26271 :     GEN p = gel(P,i);
    1107       26271 :     long k = ZV_search(Pq, p);
    1108       26271 :     gel(CHI,j) = k? mulii(gel(chi,k), powiu(p, E[i]-Eq[k])): gen_0;
    1109             :   }
    1110       15260 :   return gerepilecopy(av, CHI);
    1111             : }
    1112             : 
    1113             : /* m a Conrey log [on the canonical primitive roots], cycg the primitive
    1114             :  * roots orders */
    1115             : GEN
    1116     1184050 : znconreylog_normalize(GEN G, GEN m)
    1117             : {
    1118     1184050 :   GEN cycg = znstar_get_conreycyc(G);
    1119             :   long i, l;
    1120     1184050 :   GEN d, M = cgetg_copy(m, &l);
    1121     1184050 :   if (typ(cycg) != t_VEC || lg(cycg) != l)
    1122           0 :     pari_err_TYPE("znconreylog_normalize",mkvec2(m,cycg));
    1123     1184050 :   for (i = 1; i < l; i++) gel(M,i) = gdiv(gel(m,i), gel(cycg,i));
    1124             :   /* m[i]: image of primroot generators g_i in Q/Z */
    1125     1184050 :   M = Q_remove_denom(M, &d);
    1126     1184050 :   return mkvec2(d? d: gen_1, M);
    1127             : }
    1128             : 
    1129             : /* return normalized character on Conrey generators attached to chi: Conrey
    1130             :  * label (t_INT), char on (SNF) G.gen* (t_VEC), or Conrey log (t_COL) */
    1131             : GEN
    1132     1183973 : znconrey_normalized(GEN G, GEN chi)
    1133             : {
    1134     1183973 :   switch(typ(chi))
    1135             :   {
    1136             :     case t_INT: /* Conrey label */
    1137         420 :       return znconreylog_normalize(G, znconreylog(G, chi));
    1138             :     case t_COL: /* Conrey log */
    1139     1183462 :       if (!RgV_is_ZV(chi)) break;
    1140     1183462 :       return znconreylog_normalize(G, chi);
    1141             :     case t_VEC: /* char on G.gen */
    1142          91 :       if (!RgV_is_ZV(chi)) break;
    1143          91 :       return znconreyfromchar_normalized(G, chi);
    1144             :   }
    1145           0 :   pari_err_TYPE("znchareval",chi);
    1146             :   return NULL;/* LCOV_EXCL_LINE */
    1147             : }
    1148             : 
    1149             : /* return 1 iff chi(-1) = -1, and 0 otherwise */
    1150             : long
    1151       74816 : zncharisodd(GEN G, GEN chi)
    1152             : {
    1153             :   long i, l, s;
    1154             :   GEN N;
    1155       74816 :   if (!checkznstar_i(G)) pari_err_TYPE("zncharisodd", G);
    1156       74816 :   if (!zncharcheck(G, chi)) pari_err_TYPE("zncharisodd", chi);
    1157       74816 :   if (typ(chi) != t_COL) chi = znconreylog(G, chi);
    1158       74816 :   N = znstar_get_N(G);
    1159       74816 :   l = lg(chi);
    1160       74816 :   s = 0;
    1161       74816 :   if (!mod8(N))
    1162             :   {
    1163       15330 :     s = mpodd(gel(chi,1));
    1164       15330 :     i = 3;
    1165             :   }
    1166             :   else
    1167       59486 :     i = 1;
    1168       74816 :   for (; i < l; i++) s += mpodd(gel(chi,i));
    1169       74816 :   return odd(s);
    1170             : }
    1171             : 
    1172             : GEN
    1173         497 : znchartokronecker(GEN G, GEN chi, long flag)
    1174             : {
    1175         497 :   pari_sp av = avma;
    1176             :   long s;
    1177             :   GEN F, o;
    1178             : 
    1179         497 :   if (flag && flag != 1) pari_err_FLAG("znchartokronecker");
    1180         497 :   s = zncharisodd(G, chi)? -1: 1;
    1181         497 :   if (typ(chi) != t_COL) chi = znconreylog(G, chi);
    1182         497 :   o = zncharorder(G, chi);
    1183         497 :   if (abscmpiu(o,2) > 0) { avma = av; return gen_0; }
    1184         231 :   F = znconreyconductor(G, chi, NULL);
    1185         231 :   if (typ(F) == t_INT)
    1186             :   {
    1187         105 :     if (s < 0) F = negi(F);
    1188         105 :     return gerepileuptoint(av, F);
    1189             :   }
    1190         126 :   F = gel(F,1);
    1191         126 :   F = (s < 0)? negi(F): icopy(F);
    1192         126 :   if (!flag)
    1193             :   {
    1194          49 :     GEN MF = znstar_get_faN(G), P = gel(MF,1);
    1195          49 :     long i, l = lg(P);
    1196         140 :     for (i = 1; i < l; i++)
    1197             :     {
    1198          91 :       GEN p = gel(P,i);
    1199          91 :       if (!dvdii(F,p)) F = mulii(F,sqri(p));
    1200             :     }
    1201             :   }
    1202         126 :   return gerepileuptoint(av, F);
    1203             : }
    1204             : 
    1205             : /* (D/.) as a character mod N; assume |D| divides N and D = 0,1 mod 4*/
    1206             : GEN
    1207        1183 : znchar_quad(GEN G, GEN D)
    1208             : {
    1209        1183 :   GEN cyc = znstar_get_conreycyc(G);
    1210        1183 :   GEN gen = znstar_get_conreygen(G);
    1211        1183 :   long i, l = lg(cyc);
    1212        1183 :   GEN chi = cgetg(l, t_COL);
    1213        1484 :   for (i = 1; i < l; i++)
    1214             :   {
    1215         301 :     long k = kronecker(D, gel(gen,i));
    1216         301 :     gel(chi,i) = (k==1)? gen_0: shifti(gel(cyc,i), -1);
    1217             :   }
    1218        1183 :   return chi;
    1219             : }
    1220             : 
    1221             : GEN
    1222        1379 : znchar(GEN D)
    1223             : {
    1224        1379 :   pari_sp av = avma;
    1225             :   GEN G, chi;
    1226        1379 :   switch(typ(D))
    1227             :   {
    1228             :     case t_INT:
    1229        1162 :       if (!signe(D) || Mod4(D) > 1) pari_err_TYPE("znchar", D);
    1230        1141 :       G = znstar0(D,1);
    1231        1141 :       chi = znchar_quad(G,D);
    1232        1141 :       break;
    1233             :     case t_INTMOD:
    1234         168 :       G = znstar0(gel(D,1), 1);
    1235         168 :       chi = znconreylog(G, gel(D,2));
    1236         168 :       break;
    1237             :     case t_VEC:
    1238          42 :       if (lg(D) != 3) pari_err_TYPE("znchar", D);
    1239          35 :       G = gel(D,1);
    1240          35 :       if (!checkznstar_i(G)) pari_err_TYPE("znchar", D);
    1241          28 :       chi = gel(D,2);
    1242          28 :       if (typ(chi) == t_VEC && lg(chi) == 3 && is_vec_t(typ(gel(chi,2))))
    1243             :       { /* normalized character */
    1244           7 :         GEN n = gel(chi,1), chic = gel(chi,2);
    1245           7 :         GEN cyc = typ(chic)==t_VEC? znstar_get_cyc(G): znstar_get_conreycyc(G);
    1246           7 :         if (!char_check(cyc, chic)) pari_err_TYPE("znchar",D);
    1247           7 :         chi = char_denormalize(cyc, n, chic);
    1248             :       }
    1249          28 :       if (!zncharcheck(G, chi)) pari_err_TYPE("znchar", D);
    1250          21 :       break;
    1251             :     default:
    1252           7 :       pari_err_TYPE("znchar", D);
    1253             :       return NULL; /*LCOV_EXCL_LINE*/
    1254             :   }
    1255        1330 :   return gerepilecopy(av, mkvec2(G, chi));
    1256             : }
    1257             : 
    1258             : /* G a znstar, not stack clean */
    1259             : GEN
    1260     1145760 : znchareval(GEN G, GEN chi, GEN n, GEN z)
    1261             : {
    1262     1145760 :   GEN nchi, N = znstar_get_N(G);
    1263             :   /* avoid division by 0 */
    1264     1145760 :   if (typ(n) == t_FRAC && !equali1(gcdii(gel(n,2), N))) return not_coprime(z);
    1265     1145753 :   n = Rg_to_Fp(n, N);
    1266     1145753 :   if (!equali1(gcdii(n, N))) return not_coprime(z);
    1267             :   /* nchi: normalized character on Conrey generators */
    1268     1143898 :   nchi = znconrey_normalized(G, chi);
    1269     1143898 :   return chareval_i(nchi, znconreylog(G,n), z);
    1270             : }
    1271             : 
    1272             : /* G is a znstar, chi a Dirichlet character */
    1273             : GEN
    1274        1036 : zncharconj(GEN G, GEN chi)
    1275             : {
    1276        1036 :   switch(typ(chi))
    1277             :   {
    1278           7 :     case t_INT: chi = znconreylog(G, chi); /* fall through */
    1279        1029 :     case t_COL: return charconj(znstar_get_conreycyc(G), chi);
    1280           7 :     case t_VEC: return charconj(znstar_get_cyc(G), chi);
    1281             :   }
    1282           0 :   pari_err_TYPE("zncharconj",chi);
    1283             :   return NULL; /*LCOV_EXCL_LINE*/
    1284             : }
    1285             : 
    1286             : /* G is a znstar, chi a Dirichlet character */
    1287             : GEN
    1288       20867 : zncharorder(GEN G,  GEN chi)
    1289             : {
    1290       20867 :   switch(typ(chi))
    1291             :   {
    1292          21 :     case t_INT: chi = znconreylog(G, chi); /*fall through*/
    1293       20790 :     case t_COL: return charorder(znstar_get_conreycyc(G), chi);
    1294          77 :     case t_VEC: return charorder(znstar_get_cyc(G), chi);
    1295           0 :     default: pari_err_TYPE("zncharorder",chi);
    1296             :              return NULL; /* LCOV_EXCL_LINE */
    1297             :   }
    1298             : }
    1299             : 
    1300             : /* G is a znstar, chi a Dirichlet character */
    1301             : GEN
    1302          21 : zncharker(GEN G, GEN chi)
    1303             : {
    1304          21 :   if (typ(chi) != t_VEC) chi = znconreychar(G, chi);
    1305          21 :   return charker(znstar_get_cyc(G), chi);
    1306             : }
    1307             : 
    1308             : /* G is a znstar, 'a' is a Dirichlet character */
    1309             : GEN
    1310          84 : zncharpow(GEN G, GEN a, GEN n)
    1311             : {
    1312          84 :   switch(typ(a))
    1313             :   {
    1314          21 :     case t_INT: return Fp_pow(a, n, znstar_get_N(G));
    1315          21 :     case t_VEC: return charpow(znstar_get_cyc(G), a, n);
    1316          42 :     case t_COL: return charpow(znstar_get_conreycyc(G), a, n);
    1317           0 :     default: pari_err_TYPE("znchapow",a);
    1318             :              return NULL; /* LCOV_EXCL_LINE */
    1319             :   }
    1320             : }
    1321             : /* G is a znstar, 'a' and 'b' are Dirichlet character */
    1322             : GEN
    1323        2394 : zncharmul(GEN G, GEN a, GEN b)
    1324             : {
    1325        2394 :   long ta = typ(a), tb = typ(b);
    1326        2394 :   if (ta == tb) switch(ta)
    1327             :   {
    1328           7 :     case t_INT: return Fp_mul(a, b, znstar_get_N(G));
    1329           7 :     case t_VEC: return charmul(znstar_get_cyc(G), a, b);
    1330        2359 :     case t_COL: return charmul(znstar_get_conreycyc(G), a, b);
    1331           0 :     default: pari_err_TYPE("zncharmul",a);
    1332             :              return NULL; /* LCOV_EXCL_LINE */
    1333             :   }
    1334          21 :   if (ta != t_COL) a = znconreylog(G, a);
    1335          21 :   if (tb != t_COL) b = znconreylog(G, b);
    1336          21 :   return charmul(znstar_get_conreycyc(G), a, b);
    1337             : }
    1338             : 
    1339             : /* G is a znstar, 'a' and 'b' are Dirichlet character */
    1340             : GEN
    1341        8400 : znchardiv(GEN G, GEN a, GEN b)
    1342             : {
    1343        8400 :   long ta = typ(a), tb = typ(b);
    1344        8400 :   if (ta == tb) switch(ta)
    1345             :   {
    1346           7 :     case t_INT: return Fp_div(a, b, znstar_get_N(G));
    1347           7 :     case t_VEC: return chardiv(znstar_get_cyc(G), a, b);
    1348        8365 :     case t_COL: return chardiv(znstar_get_conreycyc(G), a, b);
    1349           0 :     default: pari_err_TYPE("znchardiv",a);
    1350             :              return NULL; /* LCOV_EXCL_LINE */
    1351             :   }
    1352          21 :   if (ta != t_COL) a = znconreylog(G, a);
    1353          21 :   if (tb != t_COL) b = znconreylog(G, b);
    1354          21 :   return chardiv(znstar_get_conreycyc(G), a, b);
    1355             : }
    1356             : 
    1357             : /* CHI mod N = \prod_p p^e; let CHI = \prod CHI_p, CHI_p mod p^e
    1358             :  * return \prod_{p | (Q,N)} CHI_p. E.g if Q = p, return chi_p */
    1359             : GEN
    1360         224 : znchardecompose(GEN G, GEN chi, GEN Q)
    1361             : {
    1362             :   GEN c, P, E, F;
    1363             :   long l, lP, i;
    1364             : 
    1365         224 :   if (!checkznstar_i(G)) pari_err_TYPE("znchardecompose", G);
    1366         224 :   if (typ(Q) != t_INT) pari_err_TYPE("znchardecompose", Q);
    1367         224 :   if (typ(chi) == t_COL)
    1368           0 :   { if (!zncharcheck(G, chi)) pari_err_TYPE("znchardecompose", chi); }
    1369             :   else
    1370         224 :     chi = znconreylog(G, chi);
    1371         224 :   l = lg(chi);
    1372         224 :   F = znstar_get_faN(G);
    1373         224 :   c = zerocol(l-1);
    1374         224 :   P = gel(F,1); /* prime divisors of N */
    1375         224 :   lP = lg(P);
    1376         224 :   E = gel(F,2); /* exponents */
    1377         672 :   for (i = 1; i < lP; i++)
    1378             :   {
    1379         448 :     GEN p = gel(P,i);
    1380         448 :     if (i == 1 && equaliu(p,2) && E[1] >= 3)
    1381             :     {
    1382         224 :       if (!mpodd(Q))
    1383             :       {
    1384         112 :         gel(c,1) = icopy(gel(chi,1));
    1385         112 :         gel(c,2) = icopy(gel(chi,2));
    1386             :       }
    1387         224 :       i = 2; /* skip P[2] = P[1] = 2 */
    1388             :     }
    1389             :     else
    1390         224 :       if (dvdii(Q, p)) gel(c,i) = icopy(gel(chi,i));
    1391             :   }
    1392         224 :   return c;
    1393             : }
    1394             : 
    1395             : GEN
    1396          21 : zncharconductor(GEN G, GEN chi)
    1397             : {
    1398          21 :   pari_sp av = avma;
    1399          21 :   GEN F = znconreyconductor(G, chi, NULL);
    1400          21 :   if (typ(F) == t_INT) return F;
    1401          14 :   return gerepilecopy(av, gel(F,1));
    1402             : }
    1403             : GEN
    1404         119 : znchartoprimitive(GEN G, GEN chi)
    1405             : {
    1406         119 :   pari_sp av = avma;
    1407         119 :   GEN chi0, F = znconreyconductor(G, chi, &chi0);
    1408         119 :   if (typ(F) == t_INT)
    1409          56 :     chi = mkvec2(G,chi);
    1410             :   else
    1411          63 :     chi = mkvec2(znstar0(F,1), chi0);
    1412         119 :   return gerepilecopy(av, chi);
    1413             : }

Generated by: LCOV version 1.11