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 to exceed 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 - base3.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.16.2 lcov report (development 29115-f22e516b23) Lines: 2040 2174 93.8 %
Date: 2024-03-28 08:06:56 Functions: 228 241 94.6 %
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; either version 2 of the License, or (at your option) any later
       8             : version. It is distributed in the hope that it will be useful, but WITHOUT
       9             : ANY WARRANTY WHATSOEVER.
      10             : 
      11             : Check the License for details. You should have received a copy of it, along
      12             : with the package; see the file 'COPYING'. If not, write to the Free Software
      13             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      14             : 
      15             : /*******************************************************************/
      16             : /*                                                                 */
      17             : /*                       BASIC NF OPERATIONS                       */
      18             : /*                                                                 */
      19             : /*******************************************************************/
      20             : #include "pari.h"
      21             : #include "paripriv.h"
      22             : 
      23             : #define DEBUGLEVEL DEBUGLEVEL_nf
      24             : 
      25             : /*******************************************************************/
      26             : /*                                                                 */
      27             : /*                OPERATIONS OVER NUMBER FIELD ELEMENTS.           */
      28             : /*     represented as column vectors over the integral basis       */
      29             : /*                                                                 */
      30             : /*******************************************************************/
      31             : static GEN
      32    40035882 : get_tab(GEN nf, long *N)
      33             : {
      34    40035882 :   GEN tab = (typ(nf) == t_MAT)? nf: gel(nf,9);
      35    40035882 :   *N = nbrows(tab); return tab;
      36             : }
      37             : 
      38             : /* x != 0, y t_INT. Return x * y (not memory clean if x = 1) */
      39             : static GEN
      40  1082884877 : _mulii(GEN x, GEN y) {
      41  1749750823 :   return is_pm1(x)? (signe(x) < 0)? negi(y): y
      42  1749606354 :                   : mulii(x, y);
      43             : }
      44             : 
      45             : GEN
      46       21090 : tablemul_ei_ej(GEN M, long i, long j)
      47             : {
      48             :   long N;
      49       21090 :   GEN tab = get_tab(M, &N);
      50       21090 :   tab += (i-1)*N; return gel(tab,j);
      51             : }
      52             : 
      53             : /* Outputs x.ei, where ei is the i-th elt of the algebra basis.
      54             :  * x an RgV of correct length and arbitrary content (polynomials, scalars...).
      55             :  * M is the multiplication table ei ej = sum_k M_k^(i,j) ek */
      56             : GEN
      57       11347 : tablemul_ei(GEN M, GEN x, long i)
      58             : {
      59             :   long j, k, N;
      60             :   GEN v, tab;
      61             : 
      62       11347 :   if (i==1) return gcopy(x);
      63       11347 :   tab = get_tab(M, &N);
      64       11347 :   if (typ(x) != t_COL) { v = zerocol(N); gel(v,i) = gcopy(x); return v; }
      65       11347 :   tab += (i-1)*N; v = cgetg(N+1,t_COL);
      66             :   /* wi . x = [ sum_j tab[k,j] x[j] ]_k */
      67       77525 :   for (k=1; k<=N; k++)
      68             :   {
      69       66178 :     pari_sp av = avma;
      70       66178 :     GEN s = gen_0;
      71      469686 :     for (j=1; j<=N; j++)
      72             :     {
      73      403508 :       GEN c = gcoeff(tab,k,j);
      74      403508 :       if (!gequal0(c)) s = gadd(s, gmul(c, gel(x,j)));
      75             :     }
      76       66178 :     gel(v,k) = gerepileupto(av,s);
      77             :   }
      78       11347 :   return v;
      79             : }
      80             : /* as tablemul_ei, assume x a ZV of correct length */
      81             : GEN
      82    23830437 : zk_ei_mul(GEN nf, GEN x, long i)
      83             : {
      84             :   long j, k, N;
      85             :   GEN v, tab;
      86             : 
      87    23830437 :   if (i==1) return ZC_copy(x);
      88    23830423 :   tab = get_tab(nf, &N); tab += (i-1)*N;
      89    23830285 :   v = cgetg(N+1,t_COL);
      90   169367839 :   for (k=1; k<=N; k++)
      91             :   {
      92   145541467 :     pari_sp av = avma;
      93   145541467 :     GEN s = gen_0;
      94  2139905985 :     for (j=1; j<=N; j++)
      95             :     {
      96  1994508417 :       GEN c = gcoeff(tab,k,j);
      97  1994508417 :       if (signe(c)) s = addii(s, _mulii(c, gel(x,j)));
      98             :     }
      99   145397568 :     gel(v,k) = gerepileuptoint(av, s);
     100             :   }
     101    23826372 :   return v;
     102             : }
     103             : 
     104             : /* table of multiplication by wi in R[w1,..., wN] */
     105             : GEN
     106       39186 : ei_multable(GEN TAB, long i)
     107             : {
     108             :   long k,N;
     109       39186 :   GEN m, tab = get_tab(TAB, &N);
     110       39186 :   tab += (i-1)*N;
     111       39186 :   m = cgetg(N+1,t_MAT);
     112      153818 :   for (k=1; k<=N; k++) gel(m,k) = gel(tab,k);
     113       39186 :   return m;
     114             : }
     115             : 
     116             : GEN
     117    10689733 : zk_multable(GEN nf, GEN x)
     118             : {
     119    10689733 :   long i, l = lg(x);
     120    10689733 :   GEN mul = cgetg(l,t_MAT);
     121    10689567 :   gel(mul,1) = x; /* assume w_1 = 1 */
     122    34157901 :   for (i=2; i<l; i++) gel(mul,i) = zk_ei_mul(nf,x,i);
     123    10685860 :   return mul;
     124             : }
     125             : GEN
     126        2520 : multable(GEN M, GEN x)
     127             : {
     128             :   long i, N;
     129             :   GEN mul;
     130        2520 :   if (typ(x) == t_MAT) return x;
     131           0 :   M = get_tab(M, &N);
     132           0 :   if (typ(x) != t_COL) return scalarmat(x, N);
     133           0 :   mul = cgetg(N+1,t_MAT);
     134           0 :   gel(mul,1) = x; /* assume w_1 = 1 */
     135           0 :   for (i=2; i<=N; i++) gel(mul,i) = tablemul_ei(M,x,i);
     136           0 :   return mul;
     137             : }
     138             : 
     139             : /* x integral in nf; table of multiplication by x in ZK = Z[w1,..., wN].
     140             :  * Return a t_INT if x is scalar, and a ZM otherwise */
     141             : GEN
     142     4886345 : zk_scalar_or_multable(GEN nf, GEN x)
     143             : {
     144     4886345 :   long tx = typ(x);
     145     4886345 :   if (tx == t_MAT || tx == t_INT) return x;
     146     4720075 :   x = nf_to_scalar_or_basis(nf, x);
     147     4719983 :   return (typ(x) == t_COL)? zk_multable(nf, x): x;
     148             : }
     149             : 
     150             : GEN
     151       21305 : nftrace(GEN nf, GEN x)
     152             : {
     153       21305 :   pari_sp av = avma;
     154       21305 :   nf = checknf(nf);
     155       21306 :   x = nf_to_scalar_or_basis(nf, x);
     156       21284 :   x = (typ(x) == t_COL)? RgV_dotproduct(x, gel(nf_get_Tr(nf),1))
     157       21305 :                        : gmulgu(x, nf_get_degree(nf));
     158       21307 :   return gerepileupto(av, x);
     159             : }
     160             : GEN
     161        1015 : rnfelttrace(GEN rnf, GEN x)
     162             : {
     163        1015 :   pari_sp av = avma;
     164        1015 :   checkrnf(rnf);
     165        1008 :   x = rnfeltabstorel(rnf, x);
     166         665 :   x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gtrace(x))
     167         770 :                           : gmulgu(x, rnf_get_degree(rnf));
     168         770 :   return gerepileupto(av, x);
     169             : }
     170             : 
     171             : static GEN
     172          35 : famatQ_to_famatZ(GEN fa)
     173             : {
     174          35 :   GEN E, F, Q, P = gel(fa,1);
     175          35 :   long i, j, l = lg(P);
     176          35 :   if (l == 1 || RgV_is_ZV(P)) return fa;
     177           7 :   Q = cgetg(2*l, t_COL);
     178           7 :   F = cgetg(2*l, t_COL); E = gel(fa, 2);
     179          35 :   for (i = j = 1; i < l; i++)
     180             :   {
     181          28 :     GEN p = gel(P,i);
     182          28 :     if (typ(p) == t_INT)
     183          14 :     { gel(Q, j) = p; gel(F, j) = gel(E, i); j++; }
     184             :     else
     185             :     {
     186          14 :       gel(Q, j) = gel(p,1); gel(F, j) = gel(E, i); j++;
     187          14 :       gel(Q, j) = gel(p,2); gel(F, j) = negi(gel(E, i)); j++;
     188             :     }
     189             :   }
     190           7 :   setlg(Q, j); setlg(F, j); return mkmat2(Q, F);
     191             : }
     192             : static GEN
     193          35 : famat_cba(GEN fa)
     194             : {
     195          35 :   GEN Q, F, P = gel(fa, 1), E = gel(fa, 2);
     196          35 :   long i, j, lQ, l = lg(P);
     197          35 :   if (l == 1) return fa;
     198          28 :   Q = ZV_cba(P); lQ = lg(Q);
     199          28 :   F = cgetg(lQ, t_COL);
     200          77 :   for (j = 1; j < lQ; j++)
     201             :   {
     202          49 :     GEN v = gen_0, q = gel(Q,j);
     203          49 :     if (!equali1(q))
     204         203 :       for (i = 1; i < l; i++)
     205             :       {
     206         161 :         long e = Z_pval(gel(P,i), q);
     207         161 :         if (e) v = addii(v, muliu(gel(E,i), e));
     208             :       }
     209          49 :     gel(F, j) = v;
     210             :   }
     211          28 :   return mkmat2(Q, F);
     212             : }
     213             : static long
     214          35 : famat_sign(GEN fa)
     215             : {
     216          35 :   GEN P = gel(fa,1), E = gel(fa,2);
     217          35 :   long i, l = lg(P), s = 1;
     218         126 :   for (i = 1; i < l; i++)
     219          91 :     if (signe(gel(P,i)) < 0 && mpodd(gel(E,i))) s = -s;
     220          35 :   return s;
     221             : }
     222             : static GEN
     223          35 : famat_abs(GEN fa)
     224             : {
     225          35 :   GEN Q, P = gel(fa,1);
     226             :   long i, l;
     227          35 :   Q = cgetg_copy(P, &l);
     228         126 :   for (i = 1; i < l; i++) gel(Q,i) = absi_shallow(gel(P,i));
     229          35 :   return mkmat2(Q, gel(fa,2));
     230             : }
     231             : 
     232             : /* assume nf is a genuine nf, fa a famat */
     233             : static GEN
     234          35 : famat_norm(GEN nf, GEN fa)
     235             : {
     236          35 :   pari_sp av = avma;
     237          35 :   GEN G, g = gel(fa,1);
     238             :   long i, l, s;
     239             : 
     240          35 :   G = cgetg_copy(g, &l);
     241         112 :   for (i = 1; i < l; i++) gel(G,i) = nfnorm(nf, gel(g,i));
     242          35 :   fa = mkmat2(G, gel(fa,2));
     243          35 :   fa = famatQ_to_famatZ(fa);
     244          35 :   s = famat_sign(fa);
     245          35 :   fa = famat_reduce(famat_abs(fa));
     246          35 :   fa = famat_cba(fa);
     247          35 :   g = factorback(fa);
     248          35 :   return gerepileupto(av, s < 0? gneg(g): g);
     249             : }
     250             : GEN
     251      195226 : nfnorm(GEN nf, GEN x)
     252             : {
     253      195226 :   pari_sp av = avma;
     254             :   GEN c, den;
     255             :   long n;
     256      195226 :   nf = checknf(nf);
     257      195226 :   n = nf_get_degree(nf);
     258      195226 :   if (typ(x) == t_MAT) return famat_norm(nf, x);
     259      195191 :   x = nf_to_scalar_or_basis(nf, x);
     260      195191 :   if (typ(x)!=t_COL)
     261      112406 :     return gerepileupto(av, gpowgs(x, n));
     262       82785 :   x = nf_to_scalar_or_alg(nf, Q_primitive_part(x, &c));
     263       82785 :   x = Q_remove_denom(x, &den);
     264       82785 :   x = ZX_resultant_all(nf_get_pol(nf), x, den, 0);
     265       82785 :   return gerepileupto(av, c ? gmul(x, gpowgs(c, n)): x);
     266             : }
     267             : 
     268             : static GEN
     269         105 : to_RgX(GEN P, long vx)
     270             : {
     271         105 :   return varn(P) == vx ? P: scalarpol_shallow(P, vx);
     272             : }
     273             : 
     274             : GEN
     275         462 : rnfeltnorm(GEN rnf, GEN x)
     276             : {
     277         462 :   pari_sp av = avma;
     278             :   GEN nf, pol;
     279             :   long v;
     280         462 :   checkrnf(rnf);
     281         455 :   v = rnf_get_varn(rnf);
     282         455 :   x = liftpol_shallow(rnfeltabstorel(rnf, x));
     283         217 :   nf = rnf_get_nf(rnf); pol = rnf_get_pol(rnf);
     284         434 :   x = (typ(x) == t_POL)
     285         105 :     ? rnfeltdown(rnf, nfX_resultant(nf,pol,to_RgX(x,v)))
     286         217 :     : gpowgs(x, rnf_get_degree(rnf));
     287         217 :   return gerepileupto(av, x);
     288             : }
     289             : 
     290             : /* x + y in nf */
     291             : GEN
     292    23295159 : nfadd(GEN nf, GEN x, GEN y)
     293             : {
     294    23295159 :   pari_sp av = avma;
     295             :   GEN z;
     296             : 
     297    23295159 :   nf = checknf(nf);
     298    23295159 :   x = nf_to_scalar_or_basis(nf, x);
     299    23295159 :   y = nf_to_scalar_or_basis(nf, y);
     300    23295159 :   if (typ(x) != t_COL)
     301    17567689 :   { z = (typ(y) == t_COL)? RgC_Rg_add(y, x): gadd(x,y); }
     302             :   else
     303     5727470 :   { z = (typ(y) == t_COL)? RgC_add(x, y): RgC_Rg_add(x, y); }
     304    23295159 :   return gerepileupto(av, z);
     305             : }
     306             : /* x - y in nf */
     307             : GEN
     308     1803268 : nfsub(GEN nf, GEN x, GEN y)
     309             : {
     310     1803268 :   pari_sp av = avma;
     311             :   GEN z;
     312             : 
     313     1803268 :   nf = checknf(nf);
     314     1803268 :   x = nf_to_scalar_or_basis(nf, x);
     315     1803268 :   y = nf_to_scalar_or_basis(nf, y);
     316     1803268 :   if (typ(x) != t_COL)
     317     1276527 :   { z = (typ(y) == t_COL)? Rg_RgC_sub(x,y): gsub(x,y); }
     318             :   else
     319      526741 :   { z = (typ(y) == t_COL)? RgC_sub(x,y): RgC_Rg_sub(x,y); }
     320     1803268 :   return gerepileupto(av, z);
     321             : }
     322             : 
     323             : /* product of ZC x,y in (true) nf; ( sum_i x_i sum_j y_j m^{i,j}_k )_k */
     324             : static GEN
     325     8804142 : nfmuli_ZC(GEN nf, GEN x, GEN y)
     326             : {
     327             :   long i, j, k, N;
     328     8804142 :   GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
     329             : 
     330    43181499 :   for (k = 1; k <= N; k++)
     331             :   {
     332    34377401 :     pari_sp av = avma;
     333    34377401 :     GEN s, TABi = TAB;
     334    34377401 :     if (k == 1)
     335     8804143 :       s = mulii(gel(x,1),gel(y,1));
     336             :     else
     337    25573119 :       s = addii(mulii(gel(x,1),gel(y,k)),
     338    25573258 :                 mulii(gel(x,k),gel(y,1)));
     339   221631978 :     for (i=2; i<=N; i++)
     340             :     {
     341   187258081 :       GEN t, xi = gel(x,i);
     342   187258081 :       TABi += N;
     343   187258081 :       if (!signe(xi)) continue;
     344             : 
     345    94780662 :       t = NULL;
     346  1065788210 :       for (j=2; j<=N; j++)
     347             :       {
     348   971009305 :         GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
     349   971009305 :         if (!signe(c)) continue;
     350   283542504 :         p1 = _mulii(c, gel(y,j));
     351   283546899 :         t = t? addii(t, p1): p1;
     352             :       }
     353    94778905 :       if (t) s = addii(s, mulii(xi, t));
     354             :     }
     355    34373897 :     gel(v,k) = gerepileuptoint(av,s);
     356             :   }
     357     8804098 :   return v;
     358             : }
     359             : static int
     360    74184785 : is_famat(GEN x) { return typ(x) == t_MAT && lg(x) == 3; }
     361             : /* product of x and y in nf */
     362             : GEN
     363    36100201 : nfmul(GEN nf, GEN x, GEN y)
     364             : {
     365             :   GEN z;
     366    36100201 :   pari_sp av = avma;
     367             : 
     368    36100201 :   if (x == y) return nfsqr(nf,x);
     369             : 
     370    32059161 :   nf = checknf(nf);
     371    32059161 :   if (is_famat(x) || is_famat(y)) return famat_mul(x, y);
     372    32058852 :   x = nf_to_scalar_or_basis(nf, x);
     373    32058853 :   y = nf_to_scalar_or_basis(nf, y);
     374    32058849 :   if (typ(x) != t_COL)
     375             :   {
     376    21710117 :     if (isintzero(x)) return gen_0;
     377    15700925 :     z = (typ(y) == t_COL)? RgC_Rg_mul(y, x): gmul(x,y); }
     378             :   else
     379             :   {
     380    10348732 :     if (typ(y) != t_COL)
     381             :     {
     382     4493349 :       if (isintzero(y)) return gen_0;
     383     1598270 :       z = RgC_Rg_mul(x, y);
     384             :     }
     385             :     else
     386             :     {
     387             :       GEN dx, dy;
     388     5855383 :       x = Q_remove_denom(x, &dx);
     389     5855384 :       y = Q_remove_denom(y, &dy);
     390     5855384 :       z = nfmuli_ZC(nf,x,y);
     391     5855384 :       dx = mul_denom(dx,dy);
     392     5855384 :       if (dx) z = ZC_Z_div(z, dx);
     393             :     }
     394             :   }
     395    23154569 :   return gerepileupto(av, z);
     396             : }
     397             : /* square of ZC x in nf */
     398             : static GEN
     399     7332327 : nfsqri_ZC(GEN nf, GEN x)
     400             : {
     401             :   long i, j, k, N;
     402     7332327 :   GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
     403             : 
     404    40498709 :   for (k = 1; k <= N; k++)
     405             :   {
     406    33166458 :     pari_sp av = avma;
     407    33166458 :     GEN s, TABi = TAB;
     408    33166458 :     if (k == 1)
     409     7332453 :       s = sqri(gel(x,1));
     410             :     else
     411    25834005 :       s = shifti(mulii(gel(x,1),gel(x,k)), 1);
     412   259896950 :     for (i=2; i<=N; i++)
     413             :     {
     414   226747514 :       GEN p1, c, t, xi = gel(x,i);
     415   226747514 :       TABi += N;
     416   226747514 :       if (!signe(xi)) continue;
     417             : 
     418    82463459 :       c = gcoeff(TABi, k, i);
     419    82463459 :       t = signe(c)? _mulii(c,xi): NULL;
     420   682111075 :       for (j=i+1; j<=N; j++)
     421             :       {
     422   599647479 :         c = gcoeff(TABi, k, j);
     423   599647479 :         if (!signe(c)) continue;
     424   234108462 :         p1 = _mulii(c, shifti(gel(x,j),1));
     425   234112990 :         t = t? addii(t, p1): p1;
     426             :       }
     427    82463596 :       if (t) s = addii(s, mulii(xi, t));
     428             :     }
     429    33149436 :     gel(v,k) = gerepileuptoint(av,s);
     430             :   }
     431     7332251 :   return v;
     432             : }
     433             : /* square of x in nf */
     434             : GEN
     435     8858015 : nfsqr(GEN nf, GEN x)
     436             : {
     437     8858015 :   pari_sp av = avma;
     438             :   GEN z;
     439             : 
     440     8858015 :   nf = checknf(nf);
     441     8858014 :   if (is_famat(x)) return famat_sqr(x);
     442     8858016 :   x = nf_to_scalar_or_basis(nf, x);
     443     8858018 :   if (typ(x) != t_COL) z = gsqr(x);
     444             :   else
     445             :   {
     446             :     GEN dx;
     447     2629159 :     x = Q_remove_denom(x, &dx);
     448     2629162 :     z = nfsqri_ZC(nf,x);
     449     2629166 :     if (dx) z = RgC_Rg_div(z, sqri(dx));
     450             :   }
     451     8858026 :   return gerepileupto(av, z);
     452             : }
     453             : 
     454             : /* x a ZC, v a t_COL of ZC/Z */
     455             : GEN
     456      204316 : zkC_multable_mul(GEN v, GEN x)
     457             : {
     458      204316 :   long i, l = lg(v);
     459      204316 :   GEN y = cgetg(l, t_COL);
     460      793681 :   for (i = 1; i < l; i++)
     461             :   {
     462      589365 :     GEN c = gel(v,i);
     463      589365 :     if (typ(c)!=t_COL) {
     464           0 :       if (!isintzero(c)) c = ZC_Z_mul(gel(x,1), c);
     465             :     } else {
     466      589365 :       c = ZM_ZC_mul(x,c);
     467      589365 :       if (ZV_isscalar(c)) c = gel(c,1);
     468             :     }
     469      589365 :     gel(y,i) = c;
     470             :   }
     471      204316 :   return y;
     472             : }
     473             : 
     474             : GEN
     475       55071 : nfC_multable_mul(GEN v, GEN x)
     476             : {
     477       55071 :   long i, l = lg(v);
     478       55071 :   GEN y = cgetg(l, t_COL);
     479      375003 :   for (i = 1; i < l; i++)
     480             :   {
     481      319932 :     GEN c = gel(v,i);
     482      319932 :     if (typ(c)!=t_COL) {
     483      266687 :       if (!isintzero(c)) c = RgC_Rg_mul(gel(x,1), c);
     484             :     } else {
     485       53245 :       c = RgM_RgC_mul(x,c);
     486       53245 :       if (QV_isscalar(c)) c = gel(c,1);
     487             :     }
     488      319932 :     gel(y,i) = c;
     489             :   }
     490       55071 :   return y;
     491             : }
     492             : 
     493             : GEN
     494      192434 : nfC_nf_mul(GEN nf, GEN v, GEN x)
     495             : {
     496             :   long tx;
     497             :   GEN y;
     498             : 
     499      192434 :   x = nf_to_scalar_or_basis(nf, x);
     500      192434 :   tx = typ(x);
     501      192434 :   if (tx != t_COL)
     502             :   {
     503             :     long l, i;
     504      145489 :     if (tx == t_INT)
     505             :     {
     506      136501 :       long s = signe(x);
     507      136501 :       if (!s) return zerocol(lg(v)-1);
     508      129258 :       if (is_pm1(x)) return s > 0? leafcopy(v): RgC_neg(v);
     509             :     }
     510       47166 :     l = lg(v); y = cgetg(l, t_COL);
     511      340753 :     for (i=1; i < l; i++)
     512             :     {
     513      293587 :       GEN c = gel(v,i);
     514      293587 :       if (typ(c) != t_COL) c = gmul(c, x); else c = RgC_Rg_mul(c, x);
     515      293587 :       gel(y,i) = c;
     516             :     }
     517       47166 :     return y;
     518             :   }
     519             :   else
     520             :   {
     521             :     GEN dx;
     522       46945 :     x = zk_multable(nf, Q_remove_denom(x,&dx));
     523       46945 :     y = nfC_multable_mul(v, x);
     524       46945 :     return dx? RgC_Rg_div(y, dx): y;
     525             :   }
     526             : }
     527             : static GEN
     528       10618 : mulbytab(GEN M, GEN c)
     529       10618 : { return typ(c) == t_COL? RgM_RgC_mul(M,c): RgC_Rg_mul(gel(M,1), c); }
     530             : GEN
     531        2520 : tablemulvec(GEN M, GEN x, GEN v)
     532             : {
     533             :   long l, i;
     534             :   GEN y;
     535             : 
     536        2520 :   if (typ(x) == t_COL && RgV_isscalar(x))
     537             :   {
     538           0 :     x = gel(x,1);
     539           0 :     return typ(v) == t_POL? RgX_Rg_mul(v,x): RgV_Rg_mul(v,x);
     540             :   }
     541        2520 :   x = multable(M, x); /* multiplication table by x */
     542        2520 :   y = cgetg_copy(v, &l);
     543        2520 :   if (typ(v) == t_POL)
     544             :   {
     545        2520 :     y[1] = v[1];
     546       13138 :     for (i=2; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
     547        2520 :     y = normalizepol(y);
     548             :   }
     549             :   else
     550             :   {
     551           0 :     for (i=1; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
     552             :   }
     553        2520 :   return y;
     554             : }
     555             : 
     556             : GEN
     557     1258314 : zkmultable_capZ(GEN mx) { return Q_denom(zkmultable_inv(mx)); }
     558             : GEN
     559     1573528 : zkmultable_inv(GEN mx) { return ZM_gauss(mx, col_ei(lg(mx)-1,1)); }
     560             : /* nf a true nf, x a ZC */
     561             : GEN
     562      315217 : zk_inv(GEN nf, GEN x) { return zkmultable_inv(zk_multable(nf,x)); }
     563             : 
     564             : /* inverse of x in nf */
     565             : GEN
     566      236551 : nfinv(GEN nf, GEN x)
     567             : {
     568      236551 :   pari_sp av = avma;
     569             :   GEN z;
     570             : 
     571      236551 :   nf = checknf(nf);
     572      236551 :   if (is_famat(x)) return famat_inv(x);
     573      236551 :   x = nf_to_scalar_or_basis(nf, x);
     574      236551 :   if (typ(x) == t_COL)
     575             :   {
     576             :     GEN d;
     577      189178 :     x = Q_remove_denom(x, &d);
     578      189178 :     z = zk_inv(nf, x);
     579      189178 :     if (d) z = RgC_Rg_mul(z, d);
     580             :   }
     581             :   else
     582       47373 :     z = ginv(x);
     583      236551 :   return gerepileupto(av, z);
     584             : }
     585             : 
     586             : /* quotient of x and y in nf */
     587             : GEN
     588       36169 : nfdiv(GEN nf, GEN x, GEN y)
     589             : {
     590       36169 :   pari_sp av = avma;
     591             :   GEN z;
     592             : 
     593       36169 :   nf = checknf(nf);
     594       36169 :   if (is_famat(x) || is_famat(y)) return famat_div(x,y);
     595       36078 :   y = nf_to_scalar_or_basis(nf, y);
     596       36078 :   if (typ(y) != t_COL)
     597             :   {
     598       21805 :     x = nf_to_scalar_or_basis(nf, x);
     599       21805 :     z = (typ(x) == t_COL)? RgC_Rg_div(x, y): gdiv(x,y);
     600             :   }
     601             :   else
     602             :   {
     603             :     GEN d;
     604       14273 :     y = Q_remove_denom(y, &d);
     605       14273 :     z = nfmul(nf, x, zk_inv(nf,y));
     606       14273 :     if (d) z = typ(z) == t_COL? RgC_Rg_mul(z, d): gmul(z, d);
     607             :   }
     608       36078 :   return gerepileupto(av, z);
     609             : }
     610             : 
     611             : /* product of INTEGERS (t_INT or ZC) x and y in (true) nf */
     612             : GEN
     613     4248683 : nfmuli(GEN nf, GEN x, GEN y)
     614             : {
     615     4248683 :   if (typ(x) == t_INT) return (typ(y) == t_COL)? ZC_Z_mul(y, x): mulii(x,y);
     616     3167148 :   if (typ(y) == t_INT) return ZC_Z_mul(x, y);
     617     2948719 :   return nfmuli_ZC(nf, x, y);
     618             : }
     619             : GEN
     620     4703215 : nfsqri(GEN nf, GEN x)
     621     4703215 : { return (typ(x) == t_INT)? sqri(x): nfsqri_ZC(nf, x); }
     622             : 
     623             : /* both x and y are RgV */
     624             : GEN
     625           0 : tablemul(GEN TAB, GEN x, GEN y)
     626             : {
     627             :   long i, j, k, N;
     628             :   GEN s, v;
     629           0 :   if (typ(x) != t_COL) return gmul(x, y);
     630           0 :   if (typ(y) != t_COL) return gmul(y, x);
     631           0 :   N = lg(x)-1;
     632           0 :   v = cgetg(N+1,t_COL);
     633           0 :   for (k=1; k<=N; k++)
     634             :   {
     635           0 :     pari_sp av = avma;
     636           0 :     GEN TABi = TAB;
     637           0 :     if (k == 1)
     638           0 :       s = gmul(gel(x,1),gel(y,1));
     639             :     else
     640           0 :       s = gadd(gmul(gel(x,1),gel(y,k)),
     641           0 :                gmul(gel(x,k),gel(y,1)));
     642           0 :     for (i=2; i<=N; i++)
     643             :     {
     644           0 :       GEN t, xi = gel(x,i);
     645           0 :       TABi += N;
     646           0 :       if (gequal0(xi)) continue;
     647             : 
     648           0 :       t = NULL;
     649           0 :       for (j=2; j<=N; j++)
     650             :       {
     651           0 :         GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
     652           0 :         if (gequal0(c)) continue;
     653           0 :         p1 = gmul(c, gel(y,j));
     654           0 :         t = t? gadd(t, p1): p1;
     655             :       }
     656           0 :       if (t) s = gadd(s, gmul(xi, t));
     657             :     }
     658           0 :     gel(v,k) = gerepileupto(av,s);
     659             :   }
     660           0 :   return v;
     661             : }
     662             : GEN
     663       48551 : tablesqr(GEN TAB, GEN x)
     664             : {
     665             :   long i, j, k, N;
     666             :   GEN s, v;
     667             : 
     668       48551 :   if (typ(x) != t_COL) return gsqr(x);
     669       48551 :   N = lg(x)-1;
     670       48551 :   v = cgetg(N+1,t_COL);
     671             : 
     672      347071 :   for (k=1; k<=N; k++)
     673             :   {
     674      298520 :     pari_sp av = avma;
     675      298520 :     GEN TABi = TAB;
     676      298520 :     if (k == 1)
     677       48551 :       s = gsqr(gel(x,1));
     678             :     else
     679      249969 :       s = gmul2n(gmul(gel(x,1),gel(x,k)), 1);
     680     1904206 :     for (i=2; i<=N; i++)
     681             :     {
     682     1605686 :       GEN p1, c, t, xi = gel(x,i);
     683     1605686 :       TABi += N;
     684     1605686 :       if (gequal0(xi)) continue;
     685             : 
     686      414470 :       c = gcoeff(TABi, k, i);
     687      414470 :       t = !gequal0(c)? gmul(c,xi): NULL;
     688     1662752 :       for (j=i+1; j<=N; j++)
     689             :       {
     690     1248282 :         c = gcoeff(TABi, k, j);
     691     1248282 :         if (gequal0(c)) continue;
     692      641564 :         p1 = gmul(gmul2n(c,1), gel(x,j));
     693      641564 :         t = t? gadd(t, p1): p1;
     694             :       }
     695      414470 :       if (t) s = gadd(s, gmul(xi, t));
     696             :     }
     697      298520 :     gel(v,k) = gerepileupto(av,s);
     698             :   }
     699       48551 :   return v;
     700             : }
     701             : 
     702             : static GEN
     703      271127 : _mul(void *data, GEN x, GEN y) { return nfmuli((GEN)data,x,y); }
     704             : static GEN
     705      809108 : _sqr(void *data, GEN x) { return nfsqri((GEN)data,x); }
     706             : 
     707             : /* Compute z^n in nf, left-shift binary powering */
     708             : GEN
     709      884876 : nfpow(GEN nf, GEN z, GEN n)
     710             : {
     711      884876 :   pari_sp av = avma;
     712             :   long s;
     713             :   GEN x, cx;
     714             : 
     715      884876 :   if (typ(n)!=t_INT) pari_err_TYPE("nfpow",n);
     716      884876 :   nf = checknf(nf);
     717      884878 :   s = signe(n); if (!s) return gen_1;
     718      884878 :   if (is_famat(z)) return famat_pow(z, n);
     719      824236 :   x = nf_to_scalar_or_basis(nf, z);
     720      824240 :   if (typ(x) != t_COL) return powgi(x,n);
     721      706130 :   if (s < 0)
     722             :   { /* simplified nfinv */
     723             :     GEN d;
     724       43473 :     x = Q_remove_denom(x, &d);
     725       43473 :     x = zk_inv(nf, x);
     726       43473 :     x = primitive_part(x, &cx);
     727       43473 :     cx = mul_content(cx, d);
     728       43473 :     n = negi(n);
     729             :   }
     730             :   else
     731      662657 :     x = primitive_part(x, &cx);
     732      706122 :   x = gen_pow_i(x, n, (void*)nf, _sqr, _mul);
     733      706122 :   if (cx)
     734       54006 :     x = gerepileupto(av, gmul(x, powgi(cx, n)));
     735             :   else
     736      652116 :     x = gerepilecopy(av, x);
     737      706129 :   return x;
     738             : }
     739             : /* Compute z^n in nf, left-shift binary powering */
     740             : GEN
     741      343733 : nfpow_u(GEN nf, GEN z, ulong n)
     742             : {
     743      343733 :   pari_sp av = avma;
     744             :   GEN x, cx;
     745             : 
     746      343733 :   if (!n) return gen_1;
     747      343733 :   x = nf_to_scalar_or_basis(nf, z);
     748      343733 :   if (typ(x) != t_COL) return gpowgs(x,n);
     749      307293 :   x = primitive_part(x, &cx);
     750      307293 :   x = gen_powu_i(x, n, (void*)nf, _sqr, _mul);
     751      307293 :   if (cx)
     752             :   {
     753      112754 :     x = gmul(x, powgi(cx, utoipos(n)));
     754      112754 :     return gerepileupto(av,x);
     755             :   }
     756      194539 :   return gerepilecopy(av, x);
     757             : }
     758             : 
     759             : long
     760          63 : nfissquare(GEN nf, GEN z, GEN *px)
     761             : {
     762          63 :   pari_sp av = avma;
     763          63 :   long v = fetch_var_higher();
     764             :   GEN R;
     765          63 :   nf = checknf(nf);
     766          63 :   if (nf_get_degree(nf) == 1)
     767             :   {
     768          21 :     z = algtobasis(nf, z);
     769          21 :     if (!issquareall(gel(z,1), px)) return gc_long(av, 0);
     770          14 :     if (px) *px = gerepileupto(av, *px); else set_avma(av);
     771          14 :     return 1;
     772             :   }
     773          42 :   z = nf_to_scalar_or_alg(nf, z);
     774          42 :   R = nfroots(nf, deg2pol_shallow(gen_m1, gen_0, z, v));
     775          42 :   delete_var(); if (lg(R) == 1) return gc_long(av, 0);
     776          35 :   if (px) *px = gerepilecopy(av, nf_to_scalar_or_basis(nf, gel(R,1)));
     777          14 :   else set_avma(av);
     778          35 :   return 1;
     779             : }
     780             : 
     781             : long
     782        8083 : nfispower(GEN nf, GEN z, long n, GEN *px)
     783             : {
     784        8083 :   pari_sp av = avma;
     785        8083 :   long v = fetch_var_higher();
     786             :   GEN R;
     787        8083 :   nf = checknf(nf);
     788        8083 :   if (nf_get_degree(nf) == 1)
     789             :   {
     790         329 :     z = algtobasis(nf, z);
     791         329 :     if (!ispower(gel(z,1), stoi(n), px)) return gc_long(av, 0);
     792         147 :     if (px) *px = gerepileupto(av, *px); else set_avma(av);
     793         147 :     return 1;
     794             :   }
     795        7754 :   if (n <= 0)
     796           0 :     pari_err_DOMAIN("nfeltispower","exponent","<=",gen_0,stoi(n));
     797        7754 :   z = nf_to_scalar_or_alg(nf, z);
     798        7754 :   if (n==1)
     799             :   {
     800           0 :     if (px) *px = gerepilecopy(av, z);
     801           0 :     return 1;
     802             :   }
     803        7754 :   R = nfroots(nf, gsub(pol_xn(n, v), z));
     804        7754 :   delete_var(); if (lg(R) == 1) return gc_long(av, 0);
     805        3157 :   if (px) *px = gerepilecopy(av, nf_to_scalar_or_basis(nf, gel(R,1)));
     806        3143 :   else set_avma(av);
     807        3157 :   return 1;
     808             : }
     809             : 
     810             : static GEN
     811          56 : idmulred(void *nf, GEN x, GEN y) { return idealmulred((GEN) nf, x, y); }
     812             : static GEN
     813         413 : idpowred(void *nf, GEN x, GEN n) { return idealpowred((GEN) nf, x, n); }
     814             : static GEN
     815       70354 : idmul(void *nf, GEN x, GEN y) { return idealmul((GEN) nf, x, y); }
     816             : static GEN
     817       86116 : idpow(void *nf, GEN x, GEN n) { return idealpow((GEN) nf, x, n); }
     818             : GEN
     819       85149 : idealfactorback(GEN nf, GEN L, GEN e, int red)
     820             : {
     821       85149 :   nf = checknf(nf);
     822       85149 :   if (red) return gen_factorback(L, e, (void*)nf, &idmulred, &idpowred, NULL);
     823       84792 :   if (!e && typ(L) == t_MAT && lg(L) == 3) { e = gel(L,2); L = gel(L,1); }
     824       84792 :   if (is_vec_t(typ(L)) && RgV_is_prV(L))
     825             :   { /* don't use gen_factorback since *= pr^v can be done more efficiently */
     826       64375 :     pari_sp av = avma;
     827       64375 :     long i, l = lg(L);
     828             :     GEN a;
     829       64375 :     if (!e) e = const_vec(l-1, gen_1);
     830       61519 :     else switch(typ(e))
     831             :     {
     832        7264 :       case t_VECSMALL: e = zv_to_ZV(e); break;
     833       54255 :       case t_VEC: case t_COL:
     834       54255 :         if (!RgV_is_ZV(e))
     835           0 :           pari_err_TYPE("factorback [not an exponent vector]", e);
     836       54255 :         break;
     837           0 :       default: pari_err_TYPE("idealfactorback", e);
     838             :     }
     839       64375 :     if (l != lg(e))
     840           0 :       pari_err_TYPE("factorback [not an exponent vector]", e);
     841       64375 :     if (l == 1 || ZV_equal0(e)) return gc_const(av, gen_1);
     842       22139 :     a = idealpow(nf, gel(L,1), gel(e,1));
     843      241817 :     for (i = 2; i < l; i++)
     844      219678 :       if (signe(gel(e,i))) a = idealmulpowprime(nf, a, gel(L,i), gel(e,i));
     845       22139 :     return gerepileupto(av, a);
     846             :   }
     847       20417 :   return gen_factorback(L, e, (void*)nf, &idmul, &idpow, NULL);
     848             : }
     849             : static GEN
     850      325136 : eltmul(void *nf, GEN x, GEN y) { return nfmul((GEN) nf, x, y); }
     851             : static GEN
     852      462155 : eltpow(void *nf, GEN x, GEN n) { return nfpow((GEN) nf, x, n); }
     853             : GEN
     854      265083 : nffactorback(GEN nf, GEN L, GEN e)
     855      265083 : { return gen_factorback(L, e, (void*)checknf(nf), &eltmul, &eltpow, NULL); }
     856             : 
     857             : static GEN
     858     3020596 : _nf_red(void *E, GEN x) { (void)E; return gcopy(x); }
     859             : 
     860             : static GEN
     861    12502657 : _nf_add(void *E, GEN x, GEN y) { return nfadd((GEN)E,x,y); }
     862             : 
     863             : static GEN
     864      735590 : _nf_neg(void *E, GEN x) { (void)E; return gneg(x); }
     865             : 
     866             : static GEN
     867    14987696 : _nf_mul(void *E, GEN x, GEN y) { return nfmul((GEN)E,x,y); }
     868             : 
     869             : static GEN
     870       51201 : _nf_inv(void *E, GEN x) { return nfinv((GEN)E,x); }
     871             : 
     872             : static GEN
     873       10358 : _nf_s(void *E, long x) { (void)E; return stoi(x); }
     874             : 
     875             : static const struct bb_field nf_field={_nf_red,_nf_add,_nf_mul,_nf_neg,
     876             :                                         _nf_inv,&gequal0,_nf_s };
     877             : 
     878      219721 : const struct bb_field *get_nf_field(void **E, GEN nf)
     879      219721 : { *E = (void*)nf; return &nf_field; }
     880             : 
     881             : GEN
     882          14 : nfM_det(GEN nf, GEN M)
     883             : {
     884             :   void *E;
     885          14 :   const struct bb_field *S = get_nf_field(&E, nf);
     886          14 :   return gen_det(M, E, S);
     887             : }
     888             : GEN
     889       10344 : nfM_inv(GEN nf, GEN M)
     890             : {
     891             :   void *E;
     892       10344 :   const struct bb_field *S = get_nf_field(&E, nf);
     893       10344 :   return gen_Gauss(M, matid(lg(M)-1), E, S);
     894             : }
     895             : 
     896             : GEN
     897           0 : nfM_ker(GEN nf, GEN M)
     898             : {
     899             :    void *E;
     900           0 :    const struct bb_field *S = get_nf_field(&E, nf);
     901           0 :    return gen_ker(M, 0, E, S);
     902             : }
     903             : 
     904             : GEN
     905        9924 : nfM_mul(GEN nf, GEN A, GEN B)
     906             : {
     907             :   void *E;
     908        9924 :   const struct bb_field *S = get_nf_field(&E, nf);
     909        9924 :   return gen_matmul(A, B, E, S);
     910             : }
     911             : GEN
     912      199439 : nfM_nfC_mul(GEN nf, GEN A, GEN B)
     913             : {
     914             :   void *E;
     915      199439 :   const struct bb_field *S = get_nf_field(&E, nf);
     916      199439 :   return gen_matcolmul(A, B, E, S);
     917             : }
     918             : 
     919             : /* valuation of integral x (ZV), with resp. to prime ideal pr */
     920             : long
     921    26030034 : ZC_nfvalrem(GEN x, GEN pr, GEN *newx)
     922             : {
     923    26030034 :   pari_sp av = avma;
     924             :   long i, v, l;
     925    26030034 :   GEN r, y, p = pr_get_p(pr), mul = pr_get_tau(pr);
     926             : 
     927             :   /* p inert */
     928    26030174 :   if (typ(mul) == t_INT) return newx? ZV_pvalrem(x, p, newx):ZV_pval(x, p);
     929    25035288 :   y = cgetg_copy(x, &l); /* will hold the new x */
     930    25036246 :   x = leafcopy(x);
     931    41041691 :   for(v=0;; v++)
     932             :   {
     933   164527125 :     for (i=1; i<l; i++)
     934             :     { /* is (x.b)[i] divisible by p ? */
     935   148511624 :       gel(y,i) = dvmdii(ZMrow_ZC_mul(mul,x,i),p,&r);
     936   148518259 :       if (r != gen_0) { if (newx) *newx = x; return v; }
     937             :     }
     938    16015501 :     swap(x, y);
     939    16015501 :     if (!newx && (v & 0xf) == 0xf) v += pr_get_e(pr) * ZV_pvalrem(x, p, &x);
     940    16015501 :     if (gc_needed(av,1))
     941             :     {
     942           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"ZC_nfvalrem, v >= %ld", v);
     943           0 :       gerepileall(av, 2, &x, &y);
     944             :     }
     945             :   }
     946             : }
     947             : long
     948    21766683 : ZC_nfval(GEN x, GEN P)
     949    21766683 : { return ZC_nfvalrem(x, P, NULL); }
     950             : 
     951             : /* v_P(x) != 0, x a ZV. Simpler version of ZC_nfvalrem */
     952             : int
     953     1243679 : ZC_prdvd(GEN x, GEN P)
     954             : {
     955     1243679 :   pari_sp av = avma;
     956             :   long i, l;
     957     1243679 :   GEN p = pr_get_p(P), mul = pr_get_tau(P);
     958     1243704 :   if (typ(mul) == t_INT) return ZV_Z_dvd(x, p);
     959     1243165 :   l = lg(x);
     960     5041477 :   for (i=1; i<l; i++)
     961     4527443 :     if (!dvdii(ZMrow_ZC_mul(mul,x,i), p)) return gc_bool(av,0);
     962      514034 :   return gc_bool(av,1);
     963             : }
     964             : 
     965             : int
     966         357 : pr_equal(GEN P, GEN Q)
     967             : {
     968         357 :   GEN gQ, p = pr_get_p(P);
     969         357 :   long e = pr_get_e(P), f = pr_get_f(P), n;
     970         357 :   if (!equalii(p, pr_get_p(Q)) || e != pr_get_e(Q) || f != pr_get_f(Q))
     971         336 :     return 0;
     972          21 :   gQ = pr_get_gen(Q); n = lg(gQ)-1;
     973          21 :   if (2*e*f > n) return 1; /* room for only one such pr */
     974          14 :   return ZV_equal(pr_get_gen(P), gQ) || ZC_prdvd(gQ, P);
     975             : }
     976             : 
     977             : GEN
     978      420721 : famat_nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
     979             : {
     980      420721 :   pari_sp av = avma;
     981      420721 :   GEN P = gel(x,1), E = gel(x,2), V = gen_0, y = NULL;
     982      420721 :   long l = lg(P), simplify = 0, i;
     983      420721 :   if (py) { *py = gen_1; y = cgetg(l, t_COL); }
     984             : 
     985     2270533 :   for (i = 1; i < l; i++)
     986             :   {
     987     1849812 :     GEN e = gel(E,i);
     988             :     long v;
     989     1849812 :     if (!signe(e))
     990             :     {
     991           7 :       if (py) gel(y,i) = gen_1;
     992           7 :       simplify = 1; continue;
     993             :     }
     994     1849805 :     v = nfvalrem(nf, gel(P,i), pr, py? &gel(y,i): NULL);
     995     1849805 :     if (v == LONG_MAX) { set_avma(av); if (py) *py = gen_0; return mkoo(); }
     996     1849805 :     V = addmulii(V, stoi(v), e);
     997             :   }
     998      420721 :   if (!py) V = gerepileuptoint(av, V);
     999             :   else
    1000             :   {
    1001          42 :     y = mkmat2(y, gel(x,2));
    1002          42 :     if (simplify) y = famat_remove_trivial(y);
    1003          42 :     gerepileall(av, 2, &V, &y); *py = y;
    1004             :   }
    1005      420721 :   return V;
    1006             : }
    1007             : long
    1008     5632785 : nfval(GEN nf, GEN x, GEN pr)
    1009             : {
    1010     5632785 :   pari_sp av = avma;
    1011             :   long w, e;
    1012             :   GEN cx, p;
    1013             : 
    1014     5632785 :   if (gequal0(x)) return LONG_MAX;
    1015     5619566 :   nf = checknf(nf);
    1016     5619563 :   checkprid(pr);
    1017     5619563 :   p = pr_get_p(pr);
    1018     5619562 :   e = pr_get_e(pr);
    1019     5619561 :   x = nf_to_scalar_or_basis(nf, x);
    1020     5619457 :   if (typ(x) != t_COL) return e*Q_pval(x,p);
    1021     2375736 :   x = Q_primitive_part(x, &cx);
    1022     2375784 :   w = ZC_nfval(x,pr);
    1023     2375733 :   if (cx) w += e*Q_pval(cx,p);
    1024     2375732 :   return gc_long(av,w);
    1025             : }
    1026             : 
    1027             : /* want to write p^v = uniformizer^(e*v) * z^v, z coprime to pr */
    1028             : /* z := tau^e / p^(e-1), algebraic integer coprime to pr; return z^v */
    1029             : static GEN
    1030      951496 : powp(GEN nf, GEN pr, long v)
    1031             : {
    1032             :   GEN b, z;
    1033             :   long e;
    1034      951496 :   if (!v) return gen_1;
    1035      424480 :   b = pr_get_tau(pr);
    1036      424480 :   if (typ(b) == t_INT) return gen_1;
    1037      121940 :   e = pr_get_e(pr);
    1038      121940 :   z = gel(b,1);
    1039      121940 :   if (e != 1) z = gdiv(nfpow_u(nf, z, e), powiu(pr_get_p(pr),e-1));
    1040      121940 :   if (v < 0) { v = -v; z = nfinv(nf, z); }
    1041      121940 :   if (v != 1) z = nfpow_u(nf, z, v);
    1042      121940 :   return z;
    1043             : }
    1044             : long
    1045     3651028 : nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
    1046             : {
    1047     3651028 :   pari_sp av = avma;
    1048             :   long w, e;
    1049             :   GEN cx, p, t;
    1050             : 
    1051     3651028 :   if (!py) return nfval(nf,x,pr);
    1052     1787783 :   if (gequal0(x)) { *py = gen_0; return LONG_MAX; }
    1053     1787726 :   nf = checknf(nf);
    1054     1787726 :   checkprid(pr);
    1055     1787726 :   p = pr_get_p(pr);
    1056     1787726 :   e = pr_get_e(pr);
    1057     1787726 :   x = nf_to_scalar_or_basis(nf, x);
    1058     1787727 :   if (typ(x) != t_COL) {
    1059      538510 :     w = Q_pvalrem(x,p, py);
    1060      538510 :     if (!w) { *py = gerepilecopy(av, x); return 0; }
    1061      330239 :     *py = gerepileupto(av, gmul(powp(nf, pr, w), *py));
    1062      330239 :     return e*w;
    1063             :   }
    1064     1249217 :   x = Q_primitive_part(x, &cx);
    1065     1249212 :   w = ZC_nfvalrem(x,pr, py);
    1066     1249203 :   if (cx)
    1067             :   {
    1068      621257 :     long v = Q_pvalrem(cx,p, &t);
    1069      621257 :     *py = nfmul(nf, *py, gmul(powp(nf,pr,v), t));
    1070      621257 :     *py = gerepileupto(av, *py);
    1071      621257 :     w += e*v;
    1072             :   }
    1073             :   else
    1074      627946 :     *py = gerepilecopy(av, *py);
    1075     1249218 :   return w;
    1076             : }
    1077             : GEN
    1078       15015 : gpnfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
    1079             : {
    1080             :   long v;
    1081       15015 :   if (is_famat(x)) return famat_nfvalrem(nf, x, pr, py);
    1082       15008 :   v = nfvalrem(nf,x,pr,py);
    1083       15008 :   return v == LONG_MAX? mkoo(): stoi(v);
    1084             : }
    1085             : 
    1086             : /* true nf */
    1087             : GEN
    1088      304906 : coltoalg(GEN nf, GEN x)
    1089             : {
    1090      304906 :   return mkpolmod( nf_to_scalar_or_alg(nf, x), nf_get_pol(nf) );
    1091             : }
    1092             : 
    1093             : GEN
    1094      357732 : basistoalg(GEN nf, GEN x)
    1095             : {
    1096             :   GEN T;
    1097             : 
    1098      357732 :   nf = checknf(nf);
    1099      357732 :   switch(typ(x))
    1100             :   {
    1101      298795 :     case t_COL: {
    1102      298795 :       pari_sp av = avma;
    1103      298795 :       return gerepilecopy(av, coltoalg(nf, x));
    1104             :     }
    1105       33390 :     case t_POLMOD:
    1106       33390 :       T = nf_get_pol(nf);
    1107       33390 :       if (!RgX_equal_var(T,gel(x,1)))
    1108           0 :         pari_err_MODULUS("basistoalg", T,gel(x,1));
    1109       33390 :       return gcopy(x);
    1110        2359 :     case t_POL:
    1111        2359 :       T = nf_get_pol(nf);
    1112        2359 :       if (varn(T) != varn(x)) pari_err_VAR("basistoalg",x,T);
    1113        2359 :       retmkpolmod(RgX_rem(x, T), ZX_copy(T));
    1114       23188 :     case t_INT:
    1115             :     case t_FRAC:
    1116       23188 :       T = nf_get_pol(nf);
    1117       23188 :       retmkpolmod(gcopy(x), ZX_copy(T));
    1118           0 :     default:
    1119           0 :       pari_err_TYPE("basistoalg",x);
    1120             :       return NULL; /* LCOV_EXCL_LINE */
    1121             :   }
    1122             : }
    1123             : 
    1124             : /* true nf, x a t_POL */
    1125             : static GEN
    1126     4551657 : pol_to_scalar_or_basis(GEN nf, GEN x)
    1127             : {
    1128     4551657 :   GEN T = nf_get_pol(nf);
    1129     4551658 :   long l = lg(x);
    1130     4551658 :   if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_basis", x,T);
    1131     4551555 :   if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
    1132     4551555 :   if (l == 2) return gen_0;
    1133     3538354 :   if (l == 3)
    1134             :   {
    1135      818413 :     x = gel(x,2);
    1136      818413 :     if (!is_rational_t(typ(x))) pari_err_TYPE("nf_to_scalar_or_basis",x);
    1137      818406 :     return x;
    1138             :   }
    1139     2719941 :   return poltobasis(nf,x);
    1140             : }
    1141             : /* Assume nf is a genuine nf. */
    1142             : GEN
    1143   160259462 : nf_to_scalar_or_basis(GEN nf, GEN x)
    1144             : {
    1145   160259462 :   switch(typ(x))
    1146             :   {
    1147    96238641 :     case t_INT: case t_FRAC:
    1148    96238641 :       return x;
    1149      555841 :     case t_POLMOD:
    1150      555841 :       x = checknfelt_mod(nf,x,"nf_to_scalar_or_basis");
    1151      555703 :       switch(typ(x))
    1152             :       {
    1153       85428 :         case t_INT: case t_FRAC: return x;
    1154      470275 :         case t_POL: return pol_to_scalar_or_basis(nf,x);
    1155             :       }
    1156           0 :       break;
    1157     4081383 :     case t_POL: return pol_to_scalar_or_basis(nf,x);
    1158    59387440 :     case t_COL:
    1159    59387440 :       if (lg(x)-1 != nf_get_degree(nf)) break;
    1160    59387135 :       return QV_isscalar(x)? gel(x,1): x;
    1161             :   }
    1162          96 :   pari_err_TYPE("nf_to_scalar_or_basis",x);
    1163             :   return NULL; /* LCOV_EXCL_LINE */
    1164             : }
    1165             : /* Let x be a polynomial with coefficients in Q or nf. Return the same
    1166             :  * polynomial with coefficients expressed as vectors (on the integral basis).
    1167             :  * No consistency checks, not memory-clean. */
    1168             : GEN
    1169       29090 : RgX_to_nfX(GEN nf, GEN x)
    1170             : {
    1171             :   long i, l;
    1172       29090 :   GEN y = cgetg_copy(x, &l); y[1] = x[1];
    1173      236661 :   for (i=2; i<l; i++) gel(y,i) = nf_to_scalar_or_basis(nf, gel(x,i));
    1174       29090 :   return y;
    1175             : }
    1176             : 
    1177             : /* Assume nf is a genuine nf. */
    1178             : GEN
    1179     3868475 : nf_to_scalar_or_alg(GEN nf, GEN x)
    1180             : {
    1181     3868475 :   switch(typ(x))
    1182             :   {
    1183       87206 :     case t_INT: case t_FRAC:
    1184       87206 :       return x;
    1185           0 :     case t_POLMOD:
    1186           0 :       x = checknfelt_mod(nf,x,"nf_to_scalar_or_alg");
    1187           0 :       if (typ(x) != t_POL) return x;
    1188             :       /* fall through */
    1189             :     case t_POL:
    1190             :     {
    1191        4704 :       GEN T = nf_get_pol(nf);
    1192        4704 :       long l = lg(x);
    1193        4704 :       if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_alg", x,T);
    1194        4704 :       if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
    1195        4704 :       if (l == 2) return gen_0;
    1196        4704 :       if (l == 3) return gel(x,2);
    1197        3192 :       return x;
    1198             :     }
    1199     3776515 :     case t_COL:
    1200             :     {
    1201             :       GEN dx;
    1202     3776515 :       if (lg(x)-1 != nf_get_degree(nf)) break;
    1203     7475452 :       if (QV_isscalar(x)) return gel(x,1);
    1204     3698769 :       x = Q_remove_denom(x, &dx);
    1205     3698814 :       x = RgV_RgC_mul(nf_get_zkprimpart(nf), x);
    1206     3698918 :       dx = mul_denom(dx, nf_get_zkden(nf));
    1207     3698910 :       return gdiv(x,dx);
    1208             :     }
    1209             :   }
    1210          50 :   pari_err_TYPE("nf_to_scalar_or_alg",x);
    1211             :   return NULL; /* LCOV_EXCL_LINE */
    1212             : }
    1213             : 
    1214             : /* gmul(A, RgX_to_RgC(x)), A t_MAT of compatible dimensions */
    1215             : GEN
    1216        1365 : RgM_RgX_mul(GEN A, GEN x)
    1217             : {
    1218        1365 :   long i, l = lg(x)-1;
    1219             :   GEN z;
    1220        1365 :   if (l == 1) return zerocol(nbrows(A));
    1221        1351 :   z = gmul(gel(x,2), gel(A,1));
    1222        2555 :   for (i = 2; i < l; i++)
    1223        1204 :     if (!gequal0(gel(x,i+1))) z = gadd(z, gmul(gel(x,i+1), gel(A,i)));
    1224        1351 :   return z;
    1225             : }
    1226             : GEN
    1227    10301133 : ZM_ZX_mul(GEN A, GEN x)
    1228             : {
    1229    10301133 :   long i, l = lg(x)-1;
    1230             :   GEN z;
    1231    10301133 :   if (l == 1) return zerocol(nbrows(A));
    1232    10299999 :   z = ZC_Z_mul(gel(A,1), gel(x,2));
    1233    32148173 :   for (i = 2; i < l ; i++)
    1234    21850900 :     if (signe(gel(x,i+1))) z = ZC_add(z, ZC_Z_mul(gel(A,i), gel(x,i+1)));
    1235    10297273 :   return z;
    1236             : }
    1237             : /* x a t_POL, nf a genuine nf. No garbage collecting. No check.  */
    1238             : GEN
    1239     9704634 : poltobasis(GEN nf, GEN x)
    1240             : {
    1241     9704634 :   GEN d, T = nf_get_pol(nf);
    1242     9704612 :   if (varn(x) != varn(T)) pari_err_VAR( "poltobasis", x,T);
    1243     9704479 :   if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
    1244     9704452 :   x = Q_remove_denom(x, &d);
    1245     9704746 :   if (!RgX_is_ZX(x)) pari_err_TYPE("poltobasis",x);
    1246     9704673 :   x = ZM_ZX_mul(nf_get_invzk(nf), x);
    1247     9702730 :   if (d) x = RgC_Rg_div(x, d);
    1248     9702827 :   return x;
    1249             : }
    1250             : 
    1251             : GEN
    1252      919011 : algtobasis(GEN nf, GEN x)
    1253             : {
    1254             :   pari_sp av;
    1255             : 
    1256      919011 :   nf = checknf(nf);
    1257      919011 :   switch(typ(x))
    1258             :   {
    1259      113254 :     case t_POLMOD:
    1260      113254 :       if (!RgX_equal_var(nf_get_pol(nf),gel(x,1)))
    1261           7 :         pari_err_MODULUS("algtobasis", nf_get_pol(nf),gel(x,1));
    1262      113247 :       x = gel(x,2);
    1263      113247 :       switch(typ(x))
    1264             :       {
    1265        8197 :         case t_INT:
    1266        8197 :         case t_FRAC: return scalarcol(x, nf_get_degree(nf));
    1267      105050 :         case t_POL:
    1268      105050 :           av = avma;
    1269      105050 :           return gerepileupto(av,poltobasis(nf,x));
    1270             :       }
    1271           0 :       break;
    1272             : 
    1273      249368 :     case t_POL:
    1274      249368 :       av = avma;
    1275      249368 :       return gerepileupto(av,poltobasis(nf,x));
    1276             : 
    1277       83036 :     case t_COL:
    1278       83036 :       if (!RgV_is_QV(x)) pari_err_TYPE("nfalgtobasis",x);
    1279       83029 :       if (lg(x)-1 != nf_get_degree(nf)) pari_err_DIM("nfalgtobasis");
    1280       83029 :       return gcopy(x);
    1281             : 
    1282      473355 :     case t_INT:
    1283      473355 :     case t_FRAC: return scalarcol(x, nf_get_degree(nf));
    1284             :   }
    1285           0 :   pari_err_TYPE("algtobasis",x);
    1286             :   return NULL; /* LCOV_EXCL_LINE */
    1287             : }
    1288             : 
    1289             : GEN
    1290       47488 : rnfbasistoalg(GEN rnf,GEN x)
    1291             : {
    1292       47488 :   const char *f = "rnfbasistoalg";
    1293             :   long lx, i;
    1294       47488 :   pari_sp av = avma;
    1295             :   GEN z, nf, R, T;
    1296             : 
    1297       47488 :   checkrnf(rnf);
    1298       47488 :   nf = rnf_get_nf(rnf);
    1299       47488 :   T = nf_get_pol(nf);
    1300       47488 :   R = QXQX_to_mod_shallow(rnf_get_pol(rnf), T);
    1301       47488 :   switch(typ(x))
    1302             :   {
    1303         875 :     case t_COL:
    1304         875 :       z = cgetg_copy(x, &lx);
    1305        2597 :       for (i=1; i<lx; i++)
    1306             :       {
    1307        1778 :         GEN c = nf_to_scalar_or_alg(nf, gel(x,i));
    1308        1722 :         if (typ(c) == t_POL) c = mkpolmod(c,T);
    1309        1722 :         gel(z,i) = c;
    1310             :       }
    1311         819 :       z = RgV_RgC_mul(gel(rnf_get_zk(rnf),1), z);
    1312         735 :       return gerepileupto(av, gmodulo(z,R));
    1313             : 
    1314       31227 :     case t_POLMOD:
    1315       31227 :       x = polmod_nffix(f, rnf, x, 0);
    1316       30954 :       if (typ(x) != t_POL) break;
    1317       14261 :       retmkpolmod(RgX_copy(x), RgX_copy(R));
    1318        1274 :     case t_POL:
    1319        1274 :       if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = gmodulo(x,T); break; }
    1320        1029 :       if (varn(x) == varn(R))
    1321             :       {
    1322         973 :         x = RgX_nffix(f,nf_get_pol(nf),x,0);
    1323         973 :         return gmodulo(x, R);
    1324             :       }
    1325          56 :       pari_err_VAR(f, x,R);
    1326             :   }
    1327       30994 :   retmkpolmod(scalarpol(x, varn(R)), RgX_copy(R));
    1328             : }
    1329             : 
    1330             : GEN
    1331        2240 : matbasistoalg(GEN nf,GEN x)
    1332             : {
    1333             :   long i, j, li, lx;
    1334        2240 :   GEN z = cgetg_copy(x, &lx);
    1335             : 
    1336        2240 :   if (lx == 1) return z;
    1337        2233 :   switch(typ(x))
    1338             :   {
    1339          77 :     case t_VEC: case t_COL:
    1340         273 :       for (i=1; i<lx; i++) gel(z,i) = basistoalg(nf, gel(x,i));
    1341          77 :       return z;
    1342        2156 :     case t_MAT: break;
    1343           0 :     default: pari_err_TYPE("matbasistoalg",x);
    1344             :   }
    1345        2156 :   li = lgcols(x);
    1346        8008 :   for (j=1; j<lx; j++)
    1347             :   {
    1348        5852 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1349        5852 :     gel(z,j) = c;
    1350       27377 :     for (i=1; i<li; i++) gel(c,i) = basistoalg(nf, gel(xj,i));
    1351             :   }
    1352        2156 :   return z;
    1353             : }
    1354             : 
    1355             : GEN
    1356       30629 : matalgtobasis(GEN nf,GEN x)
    1357             : {
    1358             :   long i, j, li, lx;
    1359       30629 :   GEN z = cgetg_copy(x, &lx);
    1360             : 
    1361       30631 :   if (lx == 1) return z;
    1362       30169 :   switch(typ(x))
    1363             :   {
    1364       30162 :     case t_VEC: case t_COL:
    1365       79025 :       for (i=1; i<lx; i++) gel(z,i) = algtobasis(nf, gel(x,i));
    1366       30162 :       return z;
    1367           7 :     case t_MAT: break;
    1368           0 :     default: pari_err_TYPE("matalgtobasis",x);
    1369             :   }
    1370           7 :   li = lgcols(x);
    1371          14 :   for (j=1; j<lx; j++)
    1372             :   {
    1373           7 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1374           7 :     gel(z,j) = c;
    1375          21 :     for (i=1; i<li; i++) gel(c,i) = algtobasis(nf, gel(xj,i));
    1376             :   }
    1377           7 :   return z;
    1378             : }
    1379             : GEN
    1380       10603 : RgM_to_nfM(GEN nf,GEN x)
    1381             : {
    1382             :   long i, j, li, lx;
    1383       10603 :   GEN z = cgetg_copy(x, &lx);
    1384             : 
    1385       10603 :   if (lx == 1) return z;
    1386       10603 :   li = lgcols(x);
    1387       79212 :   for (j=1; j<lx; j++)
    1388             :   {
    1389       68609 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1390       68609 :     gel(z,j) = c;
    1391      452164 :     for (i=1; i<li; i++) gel(c,i) = nf_to_scalar_or_basis(nf, gel(xj,i));
    1392             :   }
    1393       10603 :   return z;
    1394             : }
    1395             : GEN
    1396      145892 : RgC_to_nfC(GEN nf, GEN x)
    1397      895722 : { pari_APPLY_type(t_COL, nf_to_scalar_or_basis(nf, gel(x,i))) }
    1398             : 
    1399             : /* x a t_POLMOD, supposedly in rnf = K[z]/(T), K = Q[y]/(Tnf) */
    1400             : GEN
    1401      141303 : polmod_nffix(const char *f, GEN rnf, GEN x, int lift)
    1402      141303 : { return polmod_nffix2(f, rnf_get_nfpol(rnf), rnf_get_pol(rnf), x,lift); }
    1403             : GEN
    1404      141394 : polmod_nffix2(const char *f, GEN T, GEN R, GEN x, int lift)
    1405             : {
    1406      141394 :   if (RgX_equal_var(gel(x,1), R))
    1407             :   {
    1408      128961 :     x = gel(x,2);
    1409      128961 :     if (typ(x) == t_POL && varn(x) == varn(R))
    1410             :     {
    1411       98223 :       x = RgX_nffix(f, T, x, lift);
    1412       98223 :       switch(lg(x))
    1413             :       {
    1414        5789 :         case 2: return gen_0;
    1415       12155 :         case 3: return gel(x,2);
    1416             :       }
    1417       80279 :       return x;
    1418             :     }
    1419             :   }
    1420       43171 :   return Rg_nffix(f, T, x, lift);
    1421             : }
    1422             : GEN
    1423        1428 : rnfalgtobasis(GEN rnf,GEN x)
    1424             : {
    1425        1428 :   const char *f = "rnfalgtobasis";
    1426        1428 :   pari_sp av = avma;
    1427             :   GEN T, R;
    1428             : 
    1429        1428 :   checkrnf(rnf);
    1430        1428 :   R = rnf_get_pol(rnf);
    1431        1428 :   T = rnf_get_nfpol(rnf);
    1432        1428 :   switch(typ(x))
    1433             :   {
    1434          98 :     case t_COL:
    1435          98 :       if (lg(x)-1 != rnf_get_degree(rnf)) pari_err_DIM(f);
    1436          49 :       x = RgV_nffix(f, T, x, 0);
    1437          42 :       return gerepilecopy(av, x);
    1438             : 
    1439        1162 :     case t_POLMOD:
    1440        1162 :       x = polmod_nffix(f, rnf, x, 0);
    1441        1057 :       if (typ(x) != t_POL) break;
    1442         714 :       return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
    1443         112 :     case t_POL:
    1444         112 :       if (varn(x) == varn(T))
    1445             :       {
    1446          42 :         RgX_check_QX(x,f);
    1447          28 :         if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
    1448          28 :         x = mkpolmod(x,T); break;
    1449             :       }
    1450          70 :       x = RgX_nffix(f, T, x, 0);
    1451          56 :       if (degpol(x) >= degpol(R)) x = RgX_rem(x, R);
    1452          56 :       return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
    1453             :   }
    1454         427 :   return gerepileupto(av, scalarcol(x, rnf_get_degree(rnf)));
    1455             : }
    1456             : 
    1457             : /* Given a and b in nf, gives an algebraic integer y in nf such that a-b.y
    1458             :  * is "small" */
    1459             : GEN
    1460         259 : nfdiveuc(GEN nf, GEN a, GEN b)
    1461             : {
    1462         259 :   pari_sp av = avma;
    1463         259 :   a = nfdiv(nf,a,b);
    1464         259 :   return gerepileupto(av, ground(a));
    1465             : }
    1466             : 
    1467             : /* Given a and b in nf, gives a "small" algebraic integer r in nf
    1468             :  * of the form a-b.y */
    1469             : GEN
    1470         259 : nfmod(GEN nf, GEN a, GEN b)
    1471             : {
    1472         259 :   pari_sp av = avma;
    1473         259 :   GEN p1 = gneg_i(nfmul(nf,b,ground(nfdiv(nf,a,b))));
    1474         259 :   return gerepileupto(av, nfadd(nf,a,p1));
    1475             : }
    1476             : 
    1477             : /* Given a and b in nf, gives a two-component vector [y,r] in nf such
    1478             :  * that r=a-b.y is "small". */
    1479             : GEN
    1480         259 : nfdivrem(GEN nf, GEN a, GEN b)
    1481             : {
    1482         259 :   pari_sp av = avma;
    1483         259 :   GEN p1,z, y = ground(nfdiv(nf,a,b));
    1484             : 
    1485         259 :   p1 = gneg_i(nfmul(nf,b,y));
    1486         259 :   z = cgetg(3,t_VEC);
    1487         259 :   gel(z,1) = gcopy(y);
    1488         259 :   gel(z,2) = nfadd(nf,a,p1); return gerepileupto(av, z);
    1489             : }
    1490             : 
    1491             : /*************************************************************************/
    1492             : /**                                                                     **/
    1493             : /**                   LOGARITHMIC EMBEDDINGS                            **/
    1494             : /**                                                                     **/
    1495             : /*************************************************************************/
    1496             : 
    1497             : static int
    1498     4444157 : low_prec(GEN x)
    1499             : {
    1500     4444157 :   switch(typ(x))
    1501             :   {
    1502           0 :     case t_INT: return !signe(x);
    1503     4444157 :     case t_REAL: return !signe(x) || realprec(x) <= DEFAULTPREC;
    1504           0 :     default: return 0;
    1505             :   }
    1506             : }
    1507             : 
    1508             : static GEN
    1509       23284 : cxlog_1(GEN nf) { return zerocol(lg(nf_get_roots(nf))-1); }
    1510             : static GEN
    1511         542 : cxlog_m1(GEN nf, long prec)
    1512             : {
    1513         542 :   long i, l = lg(nf_get_roots(nf)), r1 = nf_get_r1(nf);
    1514         542 :   GEN v = cgetg(l, t_COL), p,  P;
    1515         542 :   p = mppi(prec); P = mkcomplex(gen_0, p);
    1516        1221 :   for (i = 1; i <= r1; i++) gel(v,i) = P; /* IPi*/
    1517         542 :   if (i < l) P = gmul2n(P,1);
    1518        1154 :   for (     ; i < l; i++) gel(v,i) = P; /* 2IPi */
    1519         542 :   return v;
    1520             : }
    1521             : static GEN
    1522     1673935 : ZC_cxlog(GEN nf, GEN x, long prec)
    1523             : {
    1524             :   long i, l, r1;
    1525             :   GEN v;
    1526     1673935 :   x = RgM_RgC_mul(nf_get_M(nf), Q_primpart(x));
    1527     1673934 :   l = lg(x); r1 = nf_get_r1(nf);
    1528     4292546 :   for (i = 1; i <= r1; i++)
    1529     2618612 :     if (low_prec(gel(x,i))) return NULL;
    1530     3303144 :   for (     ; i <  l;  i++)
    1531     1629212 :     if (low_prec(gnorm(gel(x,i)))) return NULL;
    1532     1673932 :   v = cgetg(l,t_COL);
    1533     4292545 :   for (i = 1; i <= r1; i++) gel(v,i) = glog(gel(x,i),prec);
    1534     3303143 :   for (     ; i <  l;  i++) gel(v,i) = gmul2n(glog(gel(x,i),prec),1);
    1535     1673934 :   return v;
    1536             : }
    1537             : static GEN
    1538      223321 : famat_cxlog(GEN nf, GEN fa, long prec)
    1539             : {
    1540      223321 :   GEN G, E, y = NULL;
    1541             :   long i, l;
    1542             : 
    1543      223321 :   if (typ(fa) != t_MAT) pari_err_TYPE("famat_cxlog",fa);
    1544      223321 :   if (lg(fa) == 1) return cxlog_1(nf);
    1545      223321 :   G = gel(fa,1);
    1546      223321 :   E = gel(fa,2); l = lg(E);
    1547     1084913 :   for (i = 1; i < l; i++)
    1548             :   {
    1549      861592 :     GEN t, e = gel(E,i), x = nf_to_scalar_or_basis(nf, gel(G,i));
    1550             :     /* multiplicative arch would be better (save logs), but exponents overflow
    1551             :      * [ could keep track of expo separately, but not worth it ] */
    1552      861592 :     switch(typ(x))
    1553             :     { /* ignore positive rationals */
    1554       16491 :       case t_FRAC: x = gel(x,1); /* fall through */
    1555      272231 :       case t_INT: if (signe(x) > 0) continue;
    1556          94 :         if (!mpodd(e)) continue;
    1557          38 :         t = cxlog_m1(nf, prec); /* we probably should not reach this line */
    1558          38 :         break;
    1559      589361 :       default: /* t_COL */
    1560      589361 :         t = ZC_cxlog(nf,x,prec); if (!t) return NULL;
    1561      589361 :         t = RgC_Rg_mul(t, e);
    1562             :     }
    1563      589399 :     y = y? RgV_add(y,t): t;
    1564             :   }
    1565      223321 :   return y ? y: cxlog_1(nf);
    1566             : }
    1567             : /* Archimedean components: [e_i Log( sigma_i(X) )], where X = primpart(x),
    1568             :  * and e_i = 1 (resp 2.) for i <= R1 (resp. > R1) */
    1569             : GEN
    1570     1309043 : nf_cxlog(GEN nf, GEN x, long prec)
    1571             : {
    1572     1309043 :   if (typ(x) == t_MAT) return famat_cxlog(nf,x,prec);
    1573     1085722 :   x = nf_to_scalar_or_basis(nf,x);
    1574     1085722 :   switch(typ(x))
    1575             :   {
    1576           0 :     case t_FRAC: x = gel(x,1); /* fall through */
    1577        1148 :     case t_INT:
    1578        1148 :       return signe(x) > 0? cxlog_1(nf): cxlog_m1(nf, prec);
    1579     1084574 :     default:
    1580     1084574 :       return ZC_cxlog(nf, x, prec);
    1581             :   }
    1582             : }
    1583             : GEN
    1584          97 : nfV_cxlog(GEN nf, GEN x, long prec)
    1585             : {
    1586             :   long i, l;
    1587          97 :   GEN v = cgetg_copy(x, &l);
    1588         167 :   for (i = 1; i < l; i++)
    1589          70 :     if (!(gel(v,i) = nf_cxlog(nf, gel(x,i), prec))) return NULL;
    1590          97 :   return v;
    1591             : }
    1592             : 
    1593             : static GEN
    1594       16681 : scalar_logembed(GEN nf, GEN u, GEN *emb)
    1595             : {
    1596             :   GEN v, logu;
    1597       16681 :   long i, s = signe(u), RU = lg(nf_get_roots(nf))-1, R1 = nf_get_r1(nf);
    1598             : 
    1599       16681 :   if (!s) pari_err_DOMAIN("nflogembed","argument","=",gen_0,u);
    1600       16681 :   v = cgetg(RU+1, t_COL); logu = logr_abs(u);
    1601       18641 :   for (i = 1; i <= R1; i++) gel(v,i) = logu;
    1602       16681 :   if (i <= RU)
    1603             :   {
    1604       15799 :     GEN logu2 = shiftr(logu,1);
    1605       61635 :     for (   ; i <= RU; i++) gel(v,i) = logu2;
    1606             :   }
    1607       16681 :   if (emb) *emb = const_col(RU, u);
    1608       16681 :   return v;
    1609             : }
    1610             : 
    1611             : static GEN
    1612        1309 : famat_logembed(GEN nf,GEN x,GEN *emb,long prec)
    1613             : {
    1614        1309 :   GEN A, M, T, a, t, g = gel(x,1), e = gel(x,2);
    1615        1309 :   long i, l = lg(e);
    1616             : 
    1617        1309 :   if (l == 1) return scalar_logembed(nf, real_1(prec), emb);
    1618        1309 :   A = NULL; T = emb? cgetg(l, t_COL): NULL;
    1619        1309 :   if (emb) *emb = M = mkmat2(T, e);
    1620       63434 :   for (i = 1; i < l; i++)
    1621             :   {
    1622       62125 :     a = nflogembed(nf, gel(g,i), &t, prec);
    1623       62125 :     if (!a) return NULL;
    1624       62125 :     a = RgC_Rg_mul(a, gel(e,i));
    1625       62125 :     A = A? RgC_add(A, a): a;
    1626       62125 :     if (emb) gel(T,i) = t;
    1627             :   }
    1628        1309 :   return A;
    1629             : }
    1630             : 
    1631             : /* Get archimedean components: [e_i log( | sigma_i(x) | )], with e_i = 1
    1632             :  * (resp 2.) for i <= R1 (resp. > R1) and set emb to the embeddings of x.
    1633             :  * Return NULL if precision problem */
    1634             : GEN
    1635       99953 : nflogembed(GEN nf, GEN x, GEN *emb, long prec)
    1636             : {
    1637             :   long i, l, r1;
    1638             :   GEN v, t;
    1639             : 
    1640       99953 :   if (typ(x) == t_MAT) return famat_logembed(nf,x,emb,prec);
    1641       98644 :   x = nf_to_scalar_or_basis(nf,x);
    1642       98644 :   if (typ(x) != t_COL) return scalar_logembed(nf, gtofp(x,prec), emb);
    1643       81963 :   x = RgM_RgC_mul(nf_get_M(nf), x);
    1644       81963 :   l = lg(x); r1 = nf_get_r1(nf); v = cgetg(l,t_COL);
    1645      108863 :   for (i = 1; i <= r1; i++)
    1646             :   {
    1647       26901 :     t = gabs(gel(x,i),prec); if (low_prec(t)) return NULL;
    1648       26901 :     gel(v,i) = glog(t,prec);
    1649             :   }
    1650      251397 :   for (   ; i < l; i++)
    1651             :   {
    1652      169435 :     t = gnorm(gel(x,i)); if (low_prec(t)) return NULL;
    1653      169435 :     gel(v,i) = glog(t,prec);
    1654             :   }
    1655       81962 :   if (emb) *emb = x;
    1656       81962 :   return v;
    1657             : }
    1658             : 
    1659             : /*************************************************************************/
    1660             : /**                                                                     **/
    1661             : /**                        REAL EMBEDDINGS                              **/
    1662             : /**                                                                     **/
    1663             : /*************************************************************************/
    1664             : static GEN
    1665      486122 : sarch_get_cyc(GEN sarch) { return gel(sarch,1); }
    1666             : static GEN
    1667      659836 : sarch_get_archp(GEN sarch) { return gel(sarch,2); }
    1668             : static GEN
    1669      163591 : sarch_get_MI(GEN sarch) { return gel(sarch,3); }
    1670             : static GEN
    1671      163592 : sarch_get_lambda(GEN sarch) { return gel(sarch,4); }
    1672             : static GEN
    1673      163591 : sarch_get_F(GEN sarch) { return gel(sarch,5); }
    1674             : 
    1675             : /* x not a scalar, true nf, return number of positive roots of char_x */
    1676             : static long
    1677        1268 : num_positive(GEN nf, GEN x)
    1678             : {
    1679        1268 :   GEN T = nf_get_pol(nf), B, charx;
    1680             :   long dnf, vnf, N;
    1681        1268 :   x = nf_to_scalar_or_alg(nf, x); /* not a scalar */
    1682        1268 :   charx = ZXQ_charpoly(x, T, 0);
    1683        1268 :   charx = ZX_radical(charx);
    1684        1268 :   N = degpol(T) / degpol(charx);
    1685             :   /* real places are unramified ? */
    1686        1268 :   if (N == 1 || ZX_sturm(charx) * N == nf_get_r1(nf))
    1687        1261 :     return ZX_sturmpart(charx, mkvec2(gen_0,mkoo())) * N;
    1688             :   /* painful case, multiply by random square until primitive */
    1689           7 :   dnf = nf_get_degree(nf);
    1690           7 :   vnf = varn(T);
    1691           7 :   B = int2n(10);
    1692             :   for(;;)
    1693           0 :   {
    1694           7 :     GEN y = RgXQ_sqr(random_FpX(dnf, vnf, B), T);
    1695           7 :     y = RgXQ_mul(x, y, T);
    1696           7 :     charx = ZXQ_charpoly(y, T, 0);
    1697           7 :     if (ZX_is_squarefree(charx))
    1698           7 :       return ZX_sturmpart(charx, mkvec2(gen_0,mkoo())) * N;
    1699             :   }
    1700             : }
    1701             : 
    1702             : /* x a QC: return sigma_k(x) where 1 <= k <= r1+r2; correct but inefficient
    1703             :  * if x in Q. M = nf_get_M(nf) */
    1704             : static GEN
    1705         658 : nfembed_i(GEN M, GEN x, long k)
    1706             : {
    1707         658 :   long i, l = lg(M);
    1708         658 :   GEN z = gel(x,1);
    1709        4074 :   for (i = 2; i < l; i++) z = gadd(z, gmul(gcoeff(M,k,i), gel(x,i)));
    1710         658 :   return z;
    1711             : }
    1712             : GEN
    1713           0 : nfembed(GEN nf, GEN x, long k)
    1714             : {
    1715           0 :   pari_sp av = avma;
    1716           0 :   nf = checknf(nf);
    1717           0 :   x = nf_to_scalar_or_basis(nf,x);
    1718           0 :   if (typ(x) != t_COL) return gerepilecopy(av, x);
    1719           0 :   return gerepileupto(av, nfembed_i(nf_get_M(nf),x,k));
    1720             : }
    1721             : 
    1722             : /* x a ZC */
    1723             : static GEN
    1724      900816 : zk_embed(GEN M, GEN x, long k)
    1725             : {
    1726      900816 :   long i, l = lg(x);
    1727      900816 :   GEN z = gel(x,1); /* times M[k,1], which is 1 */
    1728     2951605 :   for (i = 2; i < l; i++) z = mpadd(z, mpmul(gcoeff(M,k,i), gel(x,i)));
    1729      900805 :   return z;
    1730             : }
    1731             : 
    1732             : /* Given floating point approximation z of sigma_k(x), decide its sign
    1733             :  * [0/+, 1/- and -1 for FAIL] */
    1734             : static long
    1735      882541 : eval_sign_embed(GEN z)
    1736             : {
    1737      882541 :   if (typ(z) == t_REAL)
    1738             :   {
    1739      882541 :     long l = realprec(z);
    1740      882541 :     if (l <= LOWDEFAULTPREC
    1741      882541 :       || (l == LOWDEFAULTPREC + 1 && !z[l-1])) return -1; /* dubious, fail */
    1742      881736 :     if (expo(z) < 16 - l) return -1; /* same */
    1743             :   }
    1744      881693 :   return (signe(z) < 1)? 1: 0;
    1745             : }
    1746             : /* return v such that (-1)^v = sign(sigma_k(x)), x primitive ZC */
    1747             : static long
    1748      786019 : eval_sign(GEN M, GEN x, long k)
    1749      786019 : { return eval_sign_embed( zk_embed(M, x, k) ); }
    1750             : 
    1751             : /* check that signs[i..#signs] == s; signs = NULL encodes "totally positive" */
    1752             : static int
    1753           0 : oksigns(long l, GEN signs, long i, long s)
    1754             : {
    1755           0 :   if (!signs) return s == 0;
    1756           0 :   for (; i < l; i++)
    1757           0 :     if (signs[i] != s) return 0;
    1758           0 :   return 1;
    1759             : }
    1760             : /* check that signs[i] = s and signs[i+1..#signs] = 1-s */
    1761             : static int
    1762           0 : oksigns2(long l, GEN signs, long i, long s)
    1763             : {
    1764           0 :   if (!signs) return s == 0 && i == l-1;
    1765           0 :   return signs[i] == s && oksigns(l, signs, i+1, 1-s);
    1766             : }
    1767             : 
    1768             : /* true nf, x a ZC (primitive for efficiency) which is not a scalar; embx its
    1769             :  * embeddings or NULL */
    1770             : static int
    1771       80204 : nfchecksigns_i(GEN nf, GEN x, GEN embx, GEN signs, GEN archp)
    1772             : {
    1773       80204 :   long l = lg(archp), i;
    1774       80204 :   GEN M = nf_get_M(nf), sarch = NULL;
    1775       80204 :   long np = -1;
    1776      126340 :   for (i = 1; i < l; i++)
    1777             :   {
    1778             :     long s;
    1779       97898 :     if (embx)
    1780       96534 :       s = eval_sign_embed(gel(embx,i));
    1781             :     else
    1782        1364 :       s = eval_sign(M, x, archp[i]);
    1783             :     /* 0 / + or 1 / -; -1 for FAIL */
    1784       97898 :     if (s < 0) /* failure */
    1785             :     {
    1786           0 :       long ni, r1 = nf_get_r1(nf);
    1787             :       GEN xi;
    1788           0 :       if (np < 0)
    1789             :       {
    1790           0 :         np = num_positive(nf, x);
    1791           0 :         if (np == 0)  return oksigns(l, signs, i, 1);
    1792           0 :         if (np == r1) return oksigns(l, signs, i, 0);
    1793           0 :         sarch = nfarchstar(nf, NULL, identity_perm(r1));
    1794             :       }
    1795           0 :       xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
    1796           0 :       xi = Q_primpart(xi);
    1797           0 :       ni = num_positive(nf, nfmuli(nf,x,xi));
    1798           0 :       if (ni == 0)  return oksigns2(l, signs, i, 0);
    1799           0 :       if (ni == r1) return oksigns2(l, signs, i, 1);
    1800           0 :       s = ni < np? 0: 1;
    1801             :     }
    1802       97898 :     if (s != (signs? signs[i]: 0)) return 0;
    1803             :   }
    1804       28442 :   return 1;
    1805             : }
    1806             : static void
    1807         775 : pl_convert(GEN pl, GEN *psigns, GEN *parchp)
    1808             : {
    1809         775 :   long i, j, l = lg(pl);
    1810         775 :   GEN signs = cgetg(l, t_VECSMALL);
    1811         775 :   GEN archp = cgetg(l, t_VECSMALL);
    1812        2576 :   for (i = j = 1; i < l; i++)
    1813             :   {
    1814        1801 :     if (!pl[i]) continue;
    1815        1403 :     archp[j] = i;
    1816        1403 :     signs[j] = (pl[i] < 0)? 1: 0;
    1817        1403 :     j++;
    1818             :   }
    1819         775 :   setlg(archp, j); *parchp = archp;
    1820         775 :   setlg(signs, j); *psigns = signs;
    1821         775 : }
    1822             : /* pl : requested signs for real embeddings, 0 = no sign constraint */
    1823             : int
    1824       14719 : nfchecksigns(GEN nf, GEN x, GEN pl)
    1825             : {
    1826       14719 :   pari_sp av = avma;
    1827             :   GEN signs, archp;
    1828       14719 :   nf = checknf(nf);
    1829       14719 :   x = nf_to_scalar_or_basis(nf,x);
    1830       14719 :   if (typ(x) != t_COL)
    1831             :   {
    1832       13944 :     long i, l = lg(pl), s = gsigne(x);
    1833       27853 :     for (i = 1; i < l; i++)
    1834       13909 :       if (pl[i] && pl[i] != s) return gc_bool(av,0);
    1835       13944 :     return gc_bool(av,1);
    1836             :   }
    1837         775 :   pl_convert(pl, &signs, &archp);
    1838         775 :   return gc_bool(av, nfchecksigns_i(nf, x, NULL, signs, archp));
    1839             : }
    1840             : 
    1841             : /* signs = NULL: totally positive, else sign[i] = 0 (+) or 1 (-) */
    1842             : static GEN
    1843      163592 : get_C(GEN lambda, long l, GEN signs)
    1844             : {
    1845             :   long i;
    1846             :   GEN C, mlambda;
    1847      163592 :   if (!signs) return const_vec(l-1, lambda);
    1848      133926 :   C = cgetg(l, t_COL); mlambda = gneg(lambda);
    1849      342929 :   for (i = 1; i < l; i++) gel(C,i) = signs[i]? mlambda: lambda;
    1850      133926 :   return C;
    1851             : }
    1852             : /* signs = NULL: totally positive at archp.
    1853             :  * Assume that a t_COL x is not a scalar */
    1854             : static GEN
    1855      274442 : nfsetsigns(GEN nf, GEN signs, GEN x, GEN sarch)
    1856             : {
    1857      274442 :   long i, l = lg(sarch_get_archp(sarch));
    1858             :   GEN ex;
    1859             :   /* Is signature already correct ? */
    1860      274443 :   if (typ(x) != t_COL)
    1861             :   {
    1862      195018 :     long s = gsigne(x);
    1863      195018 :     if (!s) i = 1;
    1864      194997 :     else if (!signs)
    1865        4676 :       i = (s < 0)? 1: l;
    1866             :     else
    1867             :     {
    1868      190321 :       s = s < 0? 1: 0;
    1869      323928 :       for (i = 1; i < l; i++)
    1870      245306 :         if (signs[i] != s) break;
    1871             :     }
    1872      195018 :     ex = (i < l)? const_col(l-1, x): NULL;
    1873             :   }
    1874             :   else
    1875             :   { /* inefficient if x scalar, wrong if x = 0 */
    1876       79425 :     pari_sp av = avma;
    1877       79425 :     GEN cex, M = nf_get_M(nf), archp = sarch_get_archp(sarch);
    1878       79428 :     GEN xp = Q_primitive_part(x,&cex);
    1879       79429 :     ex = cgetg(l,t_COL);
    1880      194229 :     for (i = 1; i < l; i++) gel(ex,i) = zk_embed(M,xp,archp[i]);
    1881       79429 :     if (nfchecksigns_i(nf, xp, ex, signs, archp)) { ex = NULL; set_avma(av); }
    1882       51725 :     else if (cex) ex = RgC_Rg_mul(ex, cex); /* put back content */
    1883             :   }
    1884      274451 :   if (ex)
    1885             :   { /* If no, fix it */
    1886      163591 :     GEN MI = sarch_get_MI(sarch), F = sarch_get_F(sarch);
    1887      163592 :     GEN lambda = sarch_get_lambda(sarch);
    1888      163592 :     GEN t = RgC_sub(get_C(lambda, l, signs), ex);
    1889      163589 :     t = grndtoi(RgM_RgC_mul(MI,t), NULL);
    1890      163586 :     if (lg(F) != 1) t = ZM_ZC_mul(F, t);
    1891      163583 :     x = typ(x) == t_COL? RgC_add(t, x): RgC_Rg_add(t, x);
    1892             :   }
    1893      274445 :   return x;
    1894             : }
    1895             : /* - true nf
    1896             :  * - sarch = nfarchstar(nf, F);
    1897             :  * - x encodes a vector of signs at arch.archp: either a t_VECSMALL
    1898             :  *   (vector of signs as {0,1}-vector), NULL (totally positive at archp),
    1899             :  *   or a nonzero number field element (replaced by its signature at archp);
    1900             :  * - y is a nonzero number field element
    1901             :  * Return z = y (mod F) with signs(y, archp) = signs(x) (a {0,1}-vector).
    1902             :  * Not stack-clean */
    1903             : GEN
    1904      305978 : set_sign_mod_divisor(GEN nf, GEN x, GEN y, GEN sarch)
    1905             : {
    1906      305978 :   GEN archp = sarch_get_archp(sarch);
    1907      305978 :   if (lg(archp) == 1) return y;
    1908      272693 :   if (x && typ(x) != t_VECSMALL) x = nfsign_arch(nf, x, archp);
    1909      272693 :   return nfsetsigns(nf, x, nf_to_scalar_or_basis(nf,y), sarch);
    1910             : }
    1911             : 
    1912             : static GEN
    1913       83405 : setsigns_init(GEN nf, GEN archp, GEN F, GEN DATA)
    1914             : {
    1915       83405 :   GEN lambda, Mr = rowpermute(nf_get_M(nf), archp), MI = F? RgM_mul(Mr,F): Mr;
    1916       83409 :   lambda = gmul2n(matrixnorm(MI,DEFAULTPREC), -1);
    1917       83411 :   if (typ(lambda) != t_REAL) lambda = gmul(lambda, uutoQ(1001,1000));
    1918       83411 :   if (lg(archp) < lg(MI))
    1919             :   {
    1920       58891 :     GEN perm = gel(indexrank(MI), 2);
    1921       58891 :     if (!F) F = matid(nf_get_degree(nf));
    1922       58891 :     MI = vecpermute(MI, perm);
    1923       58891 :     F = vecpermute(F, perm);
    1924             :   }
    1925       83412 :   if (!F) F = cgetg(1,t_MAT);
    1926       83412 :   MI = RgM_inv(MI);
    1927       83412 :   return mkvec5(DATA, archp, MI, lambda, F);
    1928             : }
    1929             : /* F nonzero integral ideal in HNF (or NULL: Z_K), compute elements in 1+F
    1930             :  * whose sign matrix at archp is identity; archp in 'indices' format */
    1931             : GEN
    1932      259373 : nfarchstar(GEN nf, GEN F, GEN archp)
    1933             : {
    1934      259373 :   long nba = lg(archp) - 1;
    1935      259373 :   if (!nba) return mkvec2(cgetg(1,t_VEC), archp);
    1936       81653 :   if (F && equali1(gcoeff(F,1,1))) F = NULL;
    1937       81653 :   if (F) F = idealpseudored(F, nf_get_roundG(nf));
    1938       81637 :   return setsigns_init(nf, archp, F, const_vec(nba, gen_2));
    1939             : }
    1940             : 
    1941             : /*************************************************************************/
    1942             : /**                                                                     **/
    1943             : /**                         IDEALCHINESE                                **/
    1944             : /**                                                                     **/
    1945             : /*************************************************************************/
    1946             : static int
    1947        4192 : isprfact(GEN x)
    1948             : {
    1949             :   long i, l;
    1950             :   GEN L, E;
    1951        4192 :   if (typ(x) != t_MAT || lg(x) != 3) return 0;
    1952        4192 :   L = gel(x,1); l = lg(L);
    1953        4192 :   E = gel(x,2);
    1954       13965 :   for(i=1; i<l; i++)
    1955             :   {
    1956        9773 :     checkprid(gel(L,i));
    1957        9773 :     if (typ(gel(E,i)) != t_INT) return 0;
    1958             :   }
    1959        4192 :   return 1;
    1960             : }
    1961             : 
    1962             : /* initialize projectors mod pr[i]^e[i] for idealchinese */
    1963             : static GEN
    1964        4192 : pr_init(GEN nf, GEN fa, GEN w, GEN dw)
    1965             : {
    1966        4192 :   GEN U, E, F, FZ, L = gel(fa,1), E0 = gel(fa,2);
    1967        4192 :   long i, r = lg(L);
    1968             : 
    1969        4192 :   if (w && lg(w) != r) pari_err_TYPE("idealchinese", w);
    1970        4192 :   if (r == 1 && !dw) return cgetg(1,t_VEC);
    1971        4178 :   E = leafcopy(E0); /* do not destroy fa[2] */
    1972       13951 :   for (i = 1; i < r; i++)
    1973        9773 :     if (signe(gel(E,i)) < 0) gel(E,i) = gen_0;
    1974        4178 :   F = factorbackprime(nf, L, E);
    1975        4178 :   if (dw)
    1976             :   {
    1977         679 :     F = ZM_Z_mul(F, dw);
    1978        1568 :     for (i = 1; i < r; i++)
    1979             :     {
    1980         889 :       GEN pr = gel(L,i);
    1981         889 :       long e = itos(gel(E0,i)), v = idealval(nf, dw, pr);
    1982         889 :       if (e >= 0)
    1983         882 :         gel(E,i) = addiu(gel(E,i), v);
    1984           7 :       else if (v + e <= 0)
    1985           0 :         F = idealmulpowprime(nf, F, pr, stoi(-v)); /* coprime to pr */
    1986             :       else
    1987             :       {
    1988           7 :         F = idealmulpowprime(nf, F, pr, stoi(e));
    1989           7 :         gel(E,i) = stoi(v + e);
    1990             :       }
    1991             :     }
    1992             :   }
    1993        4178 :   U = cgetg(r, t_VEC);
    1994       13951 :   for (i = 1; i < r; i++)
    1995             :   {
    1996             :     GEN u;
    1997        9773 :     if (w && gequal0(gel(w,i))) u = gen_0; /* unused */
    1998             :     else
    1999             :     {
    2000        9696 :       GEN pr = gel(L,i), e = gel(E,i), t;
    2001        9696 :       t = idealdivpowprime(nf,F, pr, e);
    2002        9696 :       u = hnfmerge_get_1(t, idealpow(nf, pr, e));
    2003        9696 :       if (!u) pari_err_COPRIME("idealchinese", t,pr);
    2004             :     }
    2005        9773 :     gel(U,i) = u;
    2006             :   }
    2007        4178 :   FZ = gcoeff(F, 1, 1);
    2008        4178 :   F = idealpseudored(F, nf_get_roundG(nf));
    2009        4178 :   return mkvec2(mkvec2(F, FZ), U);
    2010             : }
    2011             : 
    2012             : static GEN
    2013        2247 : pl_normalize(GEN nf, GEN pl)
    2014             : {
    2015        2247 :   const char *fun = "idealchinese";
    2016        2247 :   if (lg(pl)-1 != nf_get_r1(nf)) pari_err_TYPE(fun,pl);
    2017        2247 :   switch(typ(pl))
    2018             :   {
    2019         693 :     case t_VEC: RgV_check_ZV(pl,fun); pl = ZV_to_zv(pl);
    2020             :       /* fall through */
    2021        2247 :     case t_VECSMALL: break;
    2022           0 :     default: pari_err_TYPE(fun,pl);
    2023             :   }
    2024        2247 :   return pl;
    2025             : }
    2026             : 
    2027             : static int
    2028        9401 : is_chineseinit(GEN x)
    2029             : {
    2030             :   GEN fa, pl;
    2031             :   long l;
    2032        9401 :   if (typ(x) != t_VEC || lg(x)!=3) return 0;
    2033        7574 :   fa = gel(x,1);
    2034        7574 :   pl = gel(x,2);
    2035        7574 :   if (typ(fa) != t_VEC || typ(pl) != t_VEC) return 0;
    2036        4207 :   l = lg(fa);
    2037        4207 :   if (l != 1)
    2038             :   {
    2039             :     GEN z;
    2040        4165 :     if (l != 3) return 0;
    2041        4165 :     z = gel(fa, 1);
    2042        4165 :     if (typ(z) != t_VEC || lg(z) != 3 || typ(gel(z,1)) != t_MAT
    2043        4158 :                         || typ(gel(z,2)) != t_INT
    2044        4158 :                         || typ(gel(fa,2)) != t_VEC)
    2045           7 :       return 0;
    2046             :   }
    2047        4200 :   l = lg(pl);
    2048        4200 :   if (l != 1)
    2049             :   {
    2050         665 :     if (l != 6 || typ(gel(pl,3)) != t_MAT || typ(gel(pl,1)) != t_VECSMALL
    2051         665 :                                           || typ(gel(pl,2)) != t_VECSMALL)
    2052           0 :       return 0;
    2053             :   }
    2054        4200 :   return 1;
    2055             : }
    2056             : 
    2057             : /* nf a true 'nf' */
    2058             : static GEN
    2059        4647 : chineseinit_i(GEN nf, GEN fa, GEN w, GEN dw)
    2060             : {
    2061        4647 :   const char *fun = "idealchineseinit";
    2062        4647 :   GEN archp = NULL, pl = NULL;
    2063        4647 :   switch(typ(fa))
    2064             :   {
    2065        2247 :     case t_VEC:
    2066        2247 :       if (is_chineseinit(fa))
    2067             :       {
    2068           0 :         if (dw) pari_err_DOMAIN(fun, "denom(y)", "!=", gen_1, w);
    2069           0 :         return fa;
    2070             :       }
    2071        2247 :       if (lg(fa) != 3) pari_err_TYPE(fun, fa);
    2072             :       /* of the form [x,s] */
    2073        2247 :       pl = pl_normalize(nf, gel(fa,2));
    2074        2247 :       fa = gel(fa,1);
    2075        2247 :       archp = vecsmall01_to_indices(pl);
    2076             :       /* keep pr_init, reset pl */
    2077        2247 :       if (is_chineseinit(fa)) { fa = gel(fa,1); break; }
    2078             :       /* fall through */
    2079             :     case t_MAT: /* factorization? */
    2080        4192 :       if (isprfact(fa)) { fa = pr_init(nf, fa, w, dw); break; }
    2081           0 :     default: pari_err_TYPE(fun,fa);
    2082             :   }
    2083             : 
    2084        4647 :   if (!pl) pl = cgetg(1,t_VEC);
    2085             :   else
    2086             :   {
    2087        2247 :     long r = lg(archp);
    2088        2247 :     if (r == 1) pl = cgetg(1, t_VEC);
    2089             :     else
    2090             :     {
    2091        1743 :       GEN F = (lg(fa) == 1)? NULL: gmael(fa,1,1), signs = cgetg(r, t_VECSMALL);
    2092             :       long i;
    2093        5054 :       for (i = 1; i < r; i++) signs[i] = (pl[archp[i]] < 0)? 1: 0;
    2094        1743 :       pl = setsigns_init(nf, archp, F, signs);
    2095             :     }
    2096             :   }
    2097        4647 :   return mkvec2(fa, pl);
    2098             : }
    2099             : 
    2100             : /* Given a prime ideal factorization x, possibly with 0 or negative exponents,
    2101             :  * and a vector w of elements of nf, gives b such that
    2102             :  * v_p(b-w_p)>=v_p(x) for all prime ideals p in the ideal factorization
    2103             :  * and v_p(b)>=0 for all other p, using the standard proof given in GTM 138. */
    2104             : GEN
    2105        8392 : idealchinese(GEN nf, GEN x0, GEN w)
    2106             : {
    2107        8392 :   const char *fun = "idealchinese";
    2108        8392 :   pari_sp av = avma;
    2109        8392 :   GEN x = x0, x1, x2, s, dw, F;
    2110             : 
    2111        8392 :   nf = checknf(nf);
    2112        8392 :   if (!w) return gerepilecopy(av, chineseinit_i(nf,x,NULL,NULL));
    2113             : 
    2114        4907 :   if (typ(w) != t_VEC) pari_err_TYPE(fun,w);
    2115        4907 :   w = Q_remove_denom(matalgtobasis(nf,w), &dw);
    2116        4907 :   if (!is_chineseinit(x)) x = chineseinit_i(nf,x,w,dw);
    2117             :   /* x is a 'chineseinit' */
    2118        4907 :   x1 = gel(x,1); s = NULL;
    2119        4907 :   x2 = gel(x,2);
    2120        4907 :   if (lg(x1) == 1) { F = NULL; dw = NULL; }
    2121             :   else
    2122             :   {
    2123        4865 :     GEN  U = gel(x1,2), FZ;
    2124        4865 :     long i, r = lg(w);
    2125        4865 :     F = gmael(x1,1,1); FZ = gmael(x1,1,2);
    2126       17596 :     for (i=1; i<r; i++)
    2127       12731 :       if (!ZV_equal0(gel(w,i)))
    2128             :       {
    2129        9626 :         GEN t = nfmuli(nf, gel(U,i), gel(w,i));
    2130        9626 :         s = s? ZC_add(s,t): t;
    2131             :       }
    2132        4865 :     if (s)
    2133             :     {
    2134        4844 :       s = ZC_reducemodmatrix(s, F);
    2135        4844 :       if (dw && x == x0) /* input was a chineseinit */
    2136             :       {
    2137           7 :         dw = modii(dw, FZ);
    2138           7 :         s = FpC_Fp_mul(s, Fp_inv(dw, FZ), FZ);
    2139           7 :         dw = NULL;
    2140             :       }
    2141        4844 :       if (ZV_isscalar(s)) s = icopy(gel(s,1));
    2142             :     }
    2143             :   }
    2144        4907 :   if (lg(x2) != 1)
    2145             :   {
    2146        1750 :     s = nfsetsigns(nf, gel(x2,1), s? s: gen_0, x2);
    2147        1750 :     if (typ(s) == t_COL && QV_isscalar(s))
    2148             :     {
    2149         294 :       s = gel(s,1); if (!dw) s = gcopy(s);
    2150             :     }
    2151             :   }
    2152        3157 :   else if (!s) return gc_const(av, gen_0);
    2153        4858 :   return gerepileupto(av, dw? gdiv(s, dw): s);
    2154             : }
    2155             : 
    2156             : /*************************************************************************/
    2157             : /**                                                                     **/
    2158             : /**                           (Z_K/I)^*                                 **/
    2159             : /**                                                                     **/
    2160             : /*************************************************************************/
    2161             : GEN
    2162        2247 : vecsmall01_to_indices(GEN v)
    2163             : {
    2164        2247 :   long i, k, l = lg(v);
    2165        2247 :   GEN p = new_chunk(l) + l;
    2166        6594 :   for (k=1, i=l-1; i; i--)
    2167        4347 :     if (v[i]) { *--p = i; k++; }
    2168        2247 :   *--p = _evallg(k) | evaltyp(t_VECSMALL);
    2169        2247 :   set_avma((pari_sp)p); return p;
    2170             : }
    2171             : GEN
    2172     1042241 : vec01_to_indices(GEN v)
    2173             : {
    2174             :   long i, k, l;
    2175             :   GEN p;
    2176             : 
    2177     1042241 :   switch (typ(v))
    2178             :   {
    2179      995544 :    case t_VECSMALL: return v;
    2180       46697 :    case t_VEC: break;
    2181           0 :    default: pari_err_TYPE("vec01_to_indices",v);
    2182             :   }
    2183       46697 :   l = lg(v);
    2184       46697 :   p = new_chunk(l) + l;
    2185      140448 :   for (k=1, i=l-1; i; i--)
    2186       93751 :     if (signe(gel(v,i))) { *--p = i; k++; }
    2187       46697 :   *--p = _evallg(k) | evaltyp(t_VECSMALL);
    2188       46697 :   set_avma((pari_sp)p); return p;
    2189             : }
    2190             : GEN
    2191      136810 : indices_to_vec01(GEN p, long r)
    2192             : {
    2193      136810 :   long i, l = lg(p);
    2194      136810 :   GEN v = zerovec(r);
    2195      206543 :   for (i = 1; i < l; i++) gel(v, p[i]) = gen_1;
    2196      136809 :   return v;
    2197             : }
    2198             : 
    2199             : /* return (column) vector of R1 signatures of x (0 or 1) */
    2200             : GEN
    2201      995544 : nfsign_arch(GEN nf, GEN x, GEN arch)
    2202             : {
    2203      995544 :   GEN sarch, M, V, archp = vec01_to_indices(arch);
    2204      995544 :   long i, s, np, n = lg(archp)-1;
    2205             :   pari_sp av;
    2206             : 
    2207      995544 :   if (!n) return cgetg(1,t_VECSMALL);
    2208      793939 :   if (typ(x) == t_MAT)
    2209             :   { /* factorisation */
    2210      273861 :     GEN g = gel(x,1), e = gel(x,2);
    2211      273861 :     long l = lg(g);
    2212      273861 :     V = zero_zv(n);
    2213      722635 :     for (i = 1; i < l; i++)
    2214      448775 :       if (mpodd(gel(e,i)))
    2215      387906 :         Flv_add_inplace(V, nfsign_arch(nf,gel(g,i),archp), 2);
    2216      273860 :     set_avma((pari_sp)V); return V;
    2217             :   }
    2218      520078 :   av = avma; V = cgetg(n+1,t_VECSMALL);
    2219      520079 :   x = nf_to_scalar_or_basis(nf, x);
    2220      520079 :   switch(typ(x))
    2221             :   {
    2222      132454 :     case t_INT:
    2223      132454 :       s = signe(x);
    2224      132454 :       if (!s) pari_err_DOMAIN("nfsign_arch","element","=",gen_0,x);
    2225      132454 :       set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
    2226         651 :     case t_FRAC:
    2227         651 :       s = signe(gel(x,1));
    2228         651 :       set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
    2229             :   }
    2230      386974 :   x = Q_primpart(x); M = nf_get_M(nf); sarch = NULL; np = -1;
    2231     1170813 :   for (i = 1; i <= n; i++)
    2232             :   {
    2233      784653 :     long s = eval_sign(M, x, archp[i]);
    2234      784644 :     if (s < 0) /* failure */
    2235             :     {
    2236         849 :       long ni, r1 = nf_get_r1(nf);
    2237             :       GEN xi;
    2238         848 :       if (np < 0)
    2239             :       {
    2240         848 :         np = num_positive(nf, x);
    2241         848 :         if (np == 0) { set_avma(av); return const_vecsmall(n, 1); }
    2242         806 :         if (np == r1){ set_avma(av); return const_vecsmall(n, 0); }
    2243         420 :         sarch = nfarchstar(nf, NULL, identity_perm(r1));
    2244             :       }
    2245         420 :       xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
    2246         420 :       xi = Q_primpart(xi);
    2247         420 :       ni = num_positive(nf, nfmuli(nf,x,xi));
    2248         420 :       if (ni == 0) { set_avma(av); V = const_vecsmall(n, 1); V[i] = 0; return V; }
    2249         413 :       if (ni == r1){ set_avma(av); V = const_vecsmall(n, 0); V[i] = 1; return V; }
    2250          35 :       s = ni < np? 0: 1;
    2251             :     }
    2252      783830 :     V[i] = s;
    2253             :   }
    2254      386160 :   set_avma((pari_sp)V); return V;
    2255             : }
    2256             : static void
    2257       35483 : chk_ind(const char *s, long i, long r1)
    2258             : {
    2259       35483 :   if (i <= 0) pari_err_DOMAIN(s, "index", "<=", gen_0, stoi(i));
    2260       35469 :   if (i > r1) pari_err_DOMAIN(s, "index", ">", utoi(r1), utoi(i));
    2261       35434 : }
    2262             : static GEN
    2263      125944 : parse_embed(GEN ind, long r, const char *f)
    2264             : {
    2265             :   long l, i;
    2266      125944 :   if (!ind) return identity_perm(r);
    2267       33418 :   switch(typ(ind))
    2268             :   {
    2269          70 :     case t_INT: ind = mkvecsmall(itos(ind)); break;
    2270          84 :     case t_VEC: case t_COL: ind = vec_to_vecsmall(ind); break;
    2271       33264 :     case t_VECSMALL: break;
    2272           0 :     default: pari_err_TYPE(f, ind);
    2273             :   }
    2274       33418 :   l = lg(ind);
    2275       68852 :   for (i = 1; i < l; i++) chk_ind(f, ind[i], r);
    2276       33369 :   return ind;
    2277             : }
    2278             : GEN
    2279      124061 : nfeltsign(GEN nf, GEN x, GEN ind0)
    2280             : {
    2281      124061 :   pari_sp av = avma;
    2282             :   long i, l;
    2283             :   GEN v, ind;
    2284      124061 :   nf = checknf(nf);
    2285      124061 :   ind = parse_embed(ind0, nf_get_r1(nf), "nfeltsign");
    2286      124040 :   l = lg(ind);
    2287      124040 :   if (is_rational_t(typ(x)))
    2288             :   { /* nfsign_arch would test this, but avoid converting t_VECSMALL -> t_VEC */
    2289             :     GEN s;
    2290       30975 :     switch(gsigne(x))
    2291             :     {
    2292       16366 :       case -1:s = gen_m1; break;
    2293       14602 :       case 1: s = gen_1; break;
    2294           7 :       default: s = gen_0; break;
    2295             :     }
    2296       30975 :     set_avma(av);
    2297       30975 :     return (ind0 && typ(ind0) == t_INT)? s: const_vec(l-1, s);
    2298             :   }
    2299       93065 :   v = nfsign_arch(nf, x, ind);
    2300       93065 :   if (ind0 && typ(ind0) == t_INT) { set_avma(av); return v[1]? gen_m1: gen_1; }
    2301       93051 :   settyp(v, t_VEC);
    2302      262311 :   for (i = 1; i < l; i++) gel(v,i) = v[i]? gen_m1: gen_1;
    2303       93051 :   return gerepileupto(av, v);
    2304             : }
    2305             : 
    2306             : /* true nf */
    2307             : GEN
    2308         294 : nfeltembed_i(GEN *pnf, GEN x, GEN ind0, long prec0)
    2309             : {
    2310             :   long i, e, l, r1, r2, prec, prec1;
    2311         294 :   GEN v, ind, cx, nf = *pnf;
    2312         294 :   nf_get_sign(nf,&r1,&r2);
    2313         294 :   x = nf_to_scalar_or_basis(nf, x);
    2314         287 :   ind = parse_embed(ind0, r1+r2, "nfeltembed");
    2315         280 :   l = lg(ind);
    2316         280 :   if (typ(x) != t_COL)
    2317             :   {
    2318           0 :     if (!(ind0 && typ(ind0) == t_INT)) x = const_vec(l-1, x);
    2319           0 :     return x;
    2320             :   }
    2321         280 :   x = Q_primitive_part(x, &cx);
    2322         280 :   prec1 = prec0; e = gexpo(x);
    2323         280 :   if (e > 8) prec1 += nbits2extraprec(e);
    2324         280 :   prec = prec1;
    2325         280 :   if (nf_get_prec(nf) < prec) nf = nfnewprec_shallow(nf, prec);
    2326         280 :   v = cgetg(l, t_VEC);
    2327             :   for(;;)
    2328          35 :   {
    2329         315 :     GEN M = nf_get_M(nf);
    2330         938 :     for (i = 1; i < l; i++)
    2331             :     {
    2332         658 :       GEN t = nfembed_i(M, x, ind[i]);
    2333         658 :       long e = gexpo(t);
    2334         658 :       if (gequal0(t) || precision(t) < prec0
    2335         658 :                      || (e < 0 && prec < prec1 + nbits2extraprec(-e)) ) break;
    2336         623 :       if (cx) t = gmul(t, cx);
    2337         623 :       gel(v,i) = t;
    2338             :     }
    2339         315 :     if (i == l) break;
    2340          35 :     prec = precdbl(prec);
    2341          35 :     if (DEBUGLEVEL>1) pari_warn(warnprec,"eltnfembed", prec);
    2342          35 :     *pnf = nf = nfnewprec_shallow(nf, prec);
    2343             :   }
    2344         280 :   if (ind0 && typ(ind0) == t_INT) v = gel(v,1);
    2345         280 :   return v;
    2346             : }
    2347             : GEN
    2348         294 : nfeltembed(GEN nf, GEN x, GEN ind0, long prec0)
    2349             : {
    2350         294 :   pari_sp av = avma; nf = checknf(nf);
    2351         294 :   return gerepilecopy(av, nfeltembed_i(&nf, x, ind0, prec0));
    2352             : }
    2353             : 
    2354             : /* number of distinct roots of sigma(f) */
    2355             : GEN
    2356        1596 : nfpolsturm(GEN nf, GEN f, GEN ind0)
    2357             : {
    2358        1596 :   pari_sp av = avma;
    2359             :   long d, l, r1, single;
    2360             :   GEN ind, u, v, vr1, T, s, t;
    2361             : 
    2362        1596 :   nf = checknf(nf); T = nf_get_pol(nf); r1 = nf_get_r1(nf);
    2363        1596 :   ind = parse_embed(ind0, r1, "nfpolsturm");
    2364        1575 :   single = ind0 && typ(ind0) == t_INT;
    2365        1575 :   l = lg(ind);
    2366             : 
    2367        1575 :   if (gequal0(f)) pari_err_ROOTS0("nfpolsturm");
    2368        1568 :   if (typ(f) == t_POL && varn(f) != varn(T))
    2369             :   {
    2370        1547 :     f = RgX_nffix("nfpolsturm", T, f,1);
    2371        1547 :     if (lg(f) == 3) f = NULL;
    2372             :   }
    2373             :   else
    2374             :   {
    2375          21 :     (void)Rg_nffix("nfpolsturm", T, f, 0);
    2376          21 :     f = NULL;
    2377             :   }
    2378        1568 :   if (!f) { set_avma(av); return single? gen_0: zerovec(l-1); }
    2379        1547 :   d = degpol(f);
    2380        1547 :   if (d == 1) { set_avma(av); return single? gen_1: const_vec(l-1,gen_1); }
    2381             : 
    2382        1505 :   vr1 = const_vecsmall(l-1, 1);
    2383        1505 :   u = Q_primpart(f); s = ZV_to_zv(nfeltsign(nf, gel(u,d+2), ind));
    2384        1505 :   v = RgX_deriv(u); t = odd(d)? leafcopy(s): zv_neg(s);
    2385             :   for(;;)
    2386         182 :   {
    2387        1687 :     GEN r = RgX_neg( Q_primpart(RgX_pseudorem(u, v)) ), sr;
    2388        1687 :     long i, dr = degpol(r);
    2389        1687 :     if (dr < 0) break;
    2390        1687 :     sr = ZV_to_zv(nfeltsign(nf, gel(r,dr+2), ind));
    2391        4144 :     for (i = 1; i < l; i++)
    2392        2457 :       if (sr[i] != s[i]) { s[i] = sr[i], vr1[i]--; }
    2393        1687 :     if (odd(dr)) sr = zv_neg(sr);
    2394        4144 :     for (i = 1; i < l; i++)
    2395        2457 :       if (sr[i] != t[i]) { t[i] = sr[i], vr1[i]++; }
    2396        1687 :     if (!dr) break;
    2397         182 :     u = v; v = r;
    2398             :   }
    2399        1505 :   if (single) return gc_stoi(av,vr1[1]);
    2400        1498 :   return gerepileupto(av, zv_to_ZV(vr1));
    2401             : }
    2402             : 
    2403             : /* True nf; return the vector of signs of x; the matrix of such if x is a vector
    2404             :  * of nf elements */
    2405             : GEN
    2406       43960 : nfsign(GEN nf, GEN x)
    2407             : {
    2408             :   long i, l;
    2409             :   GEN archp, S;
    2410             : 
    2411       43960 :   archp = identity_perm( nf_get_r1(nf) );
    2412       43960 :   if (typ(x) != t_VEC) return nfsign_arch(nf, x, archp);
    2413       35938 :   l = lg(x); S = cgetg(l, t_MAT);
    2414      148063 :   for (i=1; i<l; i++) gel(S,i) = nfsign_arch(nf, gel(x,i), archp);
    2415       35938 :   return S;
    2416             : }
    2417             : 
    2418             : /* x integral elt, A integral ideal in HNF; reduce x mod A */
    2419             : static GEN
    2420     8353558 : zk_modHNF(GEN x, GEN A)
    2421     8353558 : { return (typ(x) == t_COL)?  ZC_hnfrem(x, A): modii(x, gcoeff(A,1,1)); }
    2422             : 
    2423             : /* given an element x in Z_K and an integral ideal y in HNF, coprime with x,
    2424             :    outputs an element inverse of x modulo y */
    2425             : GEN
    2426         175 : nfinvmodideal(GEN nf, GEN x, GEN y)
    2427             : {
    2428         175 :   pari_sp av = avma;
    2429         175 :   GEN a, yZ = gcoeff(y,1,1);
    2430             : 
    2431         175 :   if (equali1(yZ)) return gen_0;
    2432         175 :   x = nf_to_scalar_or_basis(nf, x);
    2433         175 :   if (typ(x) == t_INT) return gerepileupto(av, Fp_inv(x, yZ));
    2434             : 
    2435         119 :   a = hnfmerge_get_1(idealhnf_principal(nf,x), y);
    2436         119 :   if (!a) pari_err_INV("nfinvmodideal", x);
    2437         119 :   return gerepileupto(av, zk_modHNF(nfdiv(nf,a,x), y));
    2438             : }
    2439             : 
    2440             : static GEN
    2441     3077656 : nfsqrmodideal(GEN nf, GEN x, GEN id)
    2442     3077656 : { return zk_modHNF(nfsqri(nf,x), id); }
    2443             : static GEN
    2444     7612845 : nfmulmodideal(GEN nf, GEN x, GEN y, GEN id)
    2445     7612845 : { return x? zk_modHNF(nfmuli(nf,x,y), id): y; }
    2446             : /* assume x integral, k integer, A in HNF */
    2447             : GEN
    2448     5775401 : nfpowmodideal(GEN nf,GEN x,GEN k,GEN A)
    2449             : {
    2450     5775401 :   long s = signe(k);
    2451             :   pari_sp av;
    2452             :   GEN y;
    2453             : 
    2454     5775401 :   if (!s) return gen_1;
    2455     5775401 :   av = avma;
    2456     5775401 :   x = nf_to_scalar_or_basis(nf, x);
    2457     5775633 :   if (typ(x) != t_COL) return Fp_pow(x, k, gcoeff(A,1,1));
    2458     2629529 :   if (s < 0) { k = negi(k); x = nfinvmodideal(nf, x,A); }
    2459     2629529 :   if (equali1(k)) return gerepileupto(av, s > 0? zk_modHNF(x, A): x);
    2460     1239066 :   for(y = NULL;;)
    2461             :   {
    2462     4316750 :     if (mpodd(k)) y = nfmulmodideal(nf,y,x,A);
    2463     4316755 :     k = shifti(k,-1); if (!signe(k)) break;
    2464     3077299 :     x = nfsqrmodideal(nf,x,A);
    2465             :   }
    2466     1239044 :   return gerepileupto(av, y);
    2467             : }
    2468             : 
    2469             : /* a * g^n mod id */
    2470             : static GEN
    2471     4627902 : nfmulpowmodideal(GEN nf, GEN a, GEN g, GEN n, GEN id)
    2472             : {
    2473     4627902 :   return nfmulmodideal(nf, a, nfpowmodideal(nf,g,n,id), id);
    2474             : }
    2475             : 
    2476             : /* assume (num(g[i]), id) = 1 for all i. Return prod g[i]^e[i] mod id.
    2477             :  * EX = multiple of exponent of (O_K/id)^* */
    2478             : GEN
    2479     2618519 : famat_to_nf_modideal_coprime(GEN nf, GEN g, GEN e, GEN id, GEN EX)
    2480             : {
    2481     2618519 :   GEN EXo2, plus = NULL, minus = NULL, idZ = gcoeff(id,1,1);
    2482     2618519 :   long i, lx = lg(g);
    2483             : 
    2484     2618519 :   if (equali1(idZ)) return gen_1; /* id = Z_K */
    2485     2618031 :   EXo2 = (expi(EX) > 10)? shifti(EX,-1): NULL;
    2486     8200458 :   for (i = 1; i < lx; i++)
    2487             :   {
    2488     5582530 :     GEN h, n = centermodii(gel(e,i), EX, EXo2);
    2489     5582063 :     long sn = signe(n);
    2490     5582063 :     if (!sn) continue;
    2491             : 
    2492     3973161 :     h = nf_to_scalar_or_basis(nf, gel(g,i));
    2493     3973533 :     switch(typ(h))
    2494             :     {
    2495     2315087 :       case t_INT: break;
    2496           0 :       case t_FRAC:
    2497           0 :         h = Fp_div(gel(h,1), gel(h,2), idZ); break;
    2498     1658446 :       default:
    2499             :       {
    2500             :         GEN dh;
    2501     1658446 :         h = Q_remove_denom(h, &dh);
    2502     1658598 :         if (dh) h = FpC_Fp_mul(h, Fp_inv(dh,idZ), idZ);
    2503             :       }
    2504             :     }
    2505     3973575 :     if (sn > 0)
    2506     3971757 :       plus = nfmulpowmodideal(nf, plus, h, n, id);
    2507             :     else /* sn < 0 */
    2508        1818 :       minus = nfmulpowmodideal(nf, minus, h, negi(n), id);
    2509             :   }
    2510     2617928 :   if (minus) plus = nfmulmodideal(nf, plus, nfinvmodideal(nf,minus,id), id);
    2511     2618008 :   return plus? plus: gen_1;
    2512             : }
    2513             : 
    2514             : /* given 2 integral ideals x, y in HNF s.t x | y | x^2, compute (1+x)/(1+y) in
    2515             :  * the form [[cyc],[gen], U], where U := ux^-1 as a pair [ZM, denom(U)] */
    2516             : static GEN
    2517      236789 : zidealij(GEN x, GEN y)
    2518             : {
    2519      236789 :   GEN U, G, cyc, xp = gcoeff(x,1,1), xi = hnf_invscale(x, xp);
    2520             :   long j, N;
    2521             : 
    2522             :   /* x^(-1) y = relations between the 1 + x_i (HNF) */
    2523      236779 :   cyc = ZM_snf_group(ZM_Z_divexact(ZM_mul(xi, y), xp), &U, &G);
    2524      236771 :   N = lg(cyc); G = ZM_mul(x,G); settyp(G, t_VEC); /* new generators */
    2525      573482 :   for (j=1; j<N; j++)
    2526             :   {
    2527      336733 :     GEN c = gel(G,j);
    2528      336733 :     gel(c,1) = addiu(gel(c,1), 1); /* 1 + g_j */
    2529      336714 :     if (ZV_isscalar(c)) gel(G,j) = gel(c,1);
    2530             :   }
    2531      236749 :   return mkvec4(cyc, G, ZM_mul(U,xi), xp);
    2532             : }
    2533             : 
    2534             : /* lg(x) > 1, x + 1; shallow */
    2535             : static GEN
    2536      169659 : ZC_add1(GEN x)
    2537             : {
    2538      169659 :   long i, l = lg(x);
    2539      169659 :   GEN y = cgetg(l, t_COL);
    2540      396152 :   for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
    2541      169666 :   gel(y,1) = addiu(gel(x,1), 1); return y;
    2542             : }
    2543             : /* lg(x) > 1, x - 1; shallow */
    2544             : static GEN
    2545       70452 : ZC_sub1(GEN x)
    2546             : {
    2547       70452 :   long i, l = lg(x);
    2548       70452 :   GEN y = cgetg(l, t_COL);
    2549      176849 :   for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
    2550       70452 :   gel(y,1) = subiu(gel(x,1), 1); return y;
    2551             : }
    2552             : 
    2553             : /* x,y are t_INT or ZC */
    2554             : static GEN
    2555           0 : zkadd(GEN x, GEN y)
    2556             : {
    2557           0 :   long tx = typ(x);
    2558           0 :   if (tx == typ(y))
    2559           0 :     return tx == t_INT? addii(x,y): ZC_add(x,y);
    2560             :   else
    2561           0 :     return tx == t_INT? ZC_Z_add(y,x): ZC_Z_add(x,y);
    2562             : }
    2563             : /* x a t_INT or ZC, x+1; shallow */
    2564             : static GEN
    2565      255298 : zkadd1(GEN x)
    2566             : {
    2567      255298 :   long tx = typ(x);
    2568      255298 :   return tx == t_INT? addiu(x,1): ZC_add1(x);
    2569             : }
    2570             : /* x a t_INT or ZC, x-1; shallow */
    2571             : static GEN
    2572      255345 : zksub1(GEN x)
    2573             : {
    2574      255345 :   long tx = typ(x);
    2575      255345 :   return tx == t_INT? subiu(x,1): ZC_sub1(x);
    2576             : }
    2577             : /* x,y are t_INT or ZC; x - y */
    2578             : static GEN
    2579           0 : zksub(GEN x, GEN y)
    2580             : {
    2581           0 :   long tx = typ(x), ty = typ(y);
    2582           0 :   if (tx == ty)
    2583           0 :     return tx == t_INT? subii(x,y): ZC_sub(x,y);
    2584             :   else
    2585           0 :     return tx == t_INT? Z_ZC_sub(x,y): ZC_Z_sub(x,y);
    2586             : }
    2587             : /* x is t_INT or ZM (mult. map), y is t_INT or ZC; x * y */
    2588             : static GEN
    2589      255307 : zkmul(GEN x, GEN y)
    2590             : {
    2591      255307 :   long tx = typ(x), ty = typ(y);
    2592      255307 :   if (ty == t_INT)
    2593      184865 :     return tx == t_INT? mulii(x,y): ZC_Z_mul(gel(x,1),y);
    2594             :   else
    2595       70442 :     return tx == t_INT? ZC_Z_mul(y,x): ZM_ZC_mul(x,y);
    2596             : }
    2597             : 
    2598             : /* (U,V) = 1 coprime ideals. Want z = x mod U, = y mod V; namely
    2599             :  * z =vx + uy = v(x-y) + y, where u + v = 1, u in U, v in V.
    2600             :  * zkc = [v, UV], v a t_INT or ZM (mult. by v map), UV a ZM (ideal in HNF);
    2601             :  * shallow */
    2602             : GEN
    2603           0 : zkchinese(GEN zkc, GEN x, GEN y)
    2604             : {
    2605           0 :   GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd(zkmul(v, zksub(x,y)), y);
    2606           0 :   return zk_modHNF(z, UV);
    2607             : }
    2608             : /* special case z = x mod U, = 1 mod V; shallow */
    2609             : GEN
    2610      255342 : zkchinese1(GEN zkc, GEN x)
    2611             : {
    2612      255342 :   GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd1(zkmul(v, zksub1(x)));
    2613      255302 :   return (typ(z) == t_INT)? z: ZC_hnfrem(z, UV);
    2614             : }
    2615             : static GEN
    2616      237371 : zkVchinese1(GEN zkc, GEN v)
    2617             : {
    2618             :   long i, ly;
    2619      237371 :   GEN y = cgetg_copy(v, &ly);
    2620      492656 :   for (i=1; i<ly; i++) gel(y,i) = zkchinese1(zkc, gel(v,i));
    2621      237314 :   return y;
    2622             : }
    2623             : 
    2624             : /* prepare to solve z = x (mod A), z = y mod (B) [zkchinese or zkchinese1] */
    2625             : GEN
    2626      237131 : zkchineseinit(GEN nf, GEN A, GEN B, GEN AB)
    2627             : {
    2628      237131 :   GEN v = idealaddtoone_raw(nf, A, B);
    2629             :   long e;
    2630      237120 :   if ((e = gexpo(v)) > 5)
    2631             :   {
    2632       83246 :     GEN b = (typ(v) == t_COL)? v: scalarcol_shallow(v, nf_get_degree(nf));
    2633       83246 :     b= ZC_reducemodlll(b, AB);
    2634       83251 :     if (gexpo(b) < e) v = b;
    2635             :   }
    2636      237119 :   return mkvec2(zk_scalar_or_multable(nf,v), AB);
    2637             : }
    2638             : /* prepare to solve z = x (mod A), z = 1 mod (B)
    2639             :  * and then         z = 1 (mod A), z = y mod (B) [zkchinese1 twice] */
    2640             : static GEN
    2641         259 : zkchinese1init2(GEN nf, GEN A, GEN B, GEN AB)
    2642             : {
    2643         259 :   GEN zkc = zkchineseinit(nf, A, B, AB);
    2644         259 :   GEN mv = gel(zkc,1), mu;
    2645         259 :   if (typ(mv) == t_INT) return mkvec2(zkc, mkvec2(subui(1,mv),AB));
    2646          35 :   mu = RgM_Rg_add_shallow(ZM_neg(mv), gen_1);
    2647          35 :   return mkvec2(mkvec2(mv,AB), mkvec2(mu,AB));
    2648             : }
    2649             : 
    2650             : static GEN
    2651     2148966 : apply_U(GEN L, GEN a)
    2652             : {
    2653     2148966 :   GEN e, U = gel(L,3), dU = gel(L,4);
    2654     2148966 :   if (typ(a) == t_INT)
    2655      670645 :     e = ZC_Z_mul(gel(U,1), subiu(a, 1));
    2656             :   else
    2657             :   { /* t_COL */
    2658     1478321 :     GEN t = shallowcopy(a);
    2659     1478378 :     gel(t,1) = subiu(gel(t,1), 1); /* t = a - 1 */
    2660     1478296 :     e = ZM_ZC_mul(U, t);
    2661             :   }
    2662     2148919 :   return gdiv(e, dU);
    2663             : }
    2664             : 
    2665             : /* true nf; vectors of [[cyc],[g],U.X^-1]. Assume k > 1. */
    2666             : static GEN
    2667      169057 : principal_units(GEN nf, GEN pr, long k, GEN prk)
    2668             : {
    2669             :   GEN list, prb;
    2670      169057 :   ulong mask = quadratic_prec_mask(k);
    2671      169057 :   long a = 1;
    2672             : 
    2673      169057 :   prb = pr_hnf(nf,pr);
    2674      169065 :   list = vectrunc_init(k);
    2675      405843 :   while (mask > 1)
    2676             :   {
    2677      236789 :     GEN pra = prb;
    2678      236789 :     long b = a << 1;
    2679             : 
    2680      236789 :     if (mask & 1) b--;
    2681      236789 :     mask >>= 1;
    2682             :     /* compute 1 + pr^a / 1 + pr^b, 2a <= b */
    2683      236789 :     prb = (b >= k)? prk: idealpows(nf,pr,b);
    2684      236788 :     vectrunc_append(list, zidealij(pra, prb));
    2685      236780 :     a = b;
    2686             :   }
    2687      169054 :   return list;
    2688             : }
    2689             : /* a = 1 mod (pr) return log(a) on local-gens of 1+pr/1+pr^k */
    2690             : static GEN
    2691     1328580 : log_prk1(GEN nf, GEN a, long nh, GEN L2, GEN prk)
    2692             : {
    2693     1328580 :   GEN y = cgetg(nh+1, t_COL);
    2694     1328589 :   long j, iy, c = lg(L2)-1;
    2695     3477521 :   for (j = iy = 1; j <= c; j++)
    2696             :   {
    2697     2148957 :     GEN L = gel(L2,j), cyc = gel(L,1), gen = gel(L,2), E = apply_U(L,a);
    2698     2148800 :     long i, nc = lg(cyc)-1;
    2699     2148800 :     int last = (j == c);
    2700     5807839 :     for (i = 1; i <= nc; i++, iy++)
    2701             :     {
    2702     3658907 :       GEN t, e = gel(E,i);
    2703     3658907 :       if (typ(e) != t_INT) pari_err_COPRIME("zlog_prk1", a, prk);
    2704     3658900 :       t = Fp_neg(e, gel(cyc,i));
    2705     3658881 :       gel(y,iy) = negi(t);
    2706     3659047 :       if (!last && signe(t)) a = nfmulpowmodideal(nf, a, gel(gen,i), t, prk);
    2707             :     }
    2708             :   }
    2709     1328564 :   return y;
    2710             : }
    2711             : /* true nf */
    2712             : static GEN
    2713       56615 : principal_units_relations(GEN nf, GEN L2, GEN prk, long nh)
    2714             : {
    2715       56615 :   GEN h = cgetg(nh+1,t_MAT);
    2716       56615 :   long ih, j, c = lg(L2)-1;
    2717      180955 :   for (j = ih = 1; j <= c; j++)
    2718             :   {
    2719      124340 :     GEN L = gel(L2,j), F = gel(L,1), G = gel(L,2);
    2720      124340 :     long k, lG = lg(G);
    2721      304001 :     for (k = 1; k < lG; k++,ih++)
    2722             :     { /* log(g^f) mod pr^e */
    2723      179661 :       GEN a = nfpowmodideal(nf,gel(G,k),gel(F,k),prk);
    2724      179660 :       gel(h,ih) = ZC_neg(log_prk1(nf, a, nh, L2, prk));
    2725      179661 :       gcoeff(h,ih,ih) = gel(F,k);
    2726             :     }
    2727             :   }
    2728       56615 :   return h;
    2729             : }
    2730             : /* true nf; k > 1; multiplicative group (1 + pr) / (1 + pr^k) */
    2731             : static GEN
    2732      169060 : idealprincipalunits_i(GEN nf, GEN pr, long k, GEN *pU)
    2733             : {
    2734      169060 :   GEN cyc, gen, L2, prk = idealpows(nf, pr, k);
    2735             : 
    2736      169057 :   L2 = principal_units(nf, pr, k, prk);
    2737      169061 :   if (k == 2)
    2738             :   {
    2739      112448 :     GEN L = gel(L2,1);
    2740      112448 :     cyc = gel(L,1);
    2741      112448 :     gen = gel(L,2);
    2742      112448 :     if (pU) *pU = matid(lg(gen)-1);
    2743             :   }
    2744             :   else
    2745             :   {
    2746       56613 :     long c = lg(L2), j;
    2747       56613 :     GEN EX, h, Ui, vg = cgetg(c, t_VEC);
    2748      180948 :     for (j = 1; j < c; j++) gel(vg, j) = gmael(L2,j,2);
    2749       56613 :     vg = shallowconcat1(vg);
    2750       56615 :     h = principal_units_relations(nf, L2, prk, lg(vg)-1);
    2751       56616 :     h = ZM_hnfall_i(h, NULL, 0);
    2752       56616 :     cyc = ZM_snf_group(h, pU, &Ui);
    2753       56616 :     c = lg(Ui); gen = cgetg(c, t_VEC); EX = cyc_get_expo(cyc);
    2754      188160 :     for (j = 1; j < c; j++)
    2755      131544 :       gel(gen,j) = famat_to_nf_modideal_coprime(nf, vg, gel(Ui,j), prk, EX);
    2756             :   }
    2757      169063 :   return mkvec4(cyc, gen, prk, L2);
    2758             : }
    2759             : GEN
    2760         154 : idealprincipalunits(GEN nf, GEN pr, long k)
    2761             : {
    2762             :   pari_sp av;
    2763             :   GEN v;
    2764         154 :   nf = checknf(nf);
    2765         154 :   if (k == 1) { checkprid(pr); retmkvec3(gen_1,cgetg(1,t_VEC),cgetg(1,t_VEC)); }
    2766         147 :   av = avma; v = idealprincipalunits_i(nf, pr, k, NULL);
    2767         147 :   return gerepilecopy(av, mkvec3(powiu(pr_norm(pr), k-1), gel(v,1), gel(v,2)));
    2768             : }
    2769             : 
    2770             : /* true nf; given an ideal pr^k dividing an integral ideal x (in HNF form)
    2771             :  * compute an 'sprk', the structure of G = (Z_K/pr^k)^* [ x = NULL for x=pr^k ]
    2772             :  * Return a vector with at least 4 components [cyc],[gen],[HNF pr^k,pr,k],ff,
    2773             :  * where
    2774             :  * cyc : type of G as abelian group (SNF)
    2775             :  * gen : generators of G, coprime to x
    2776             :  * pr^k: in HNF
    2777             :  * ff  : data for log_g in (Z_K/pr)^*
    2778             :  * Two extra components are present iff k > 1: L2, U
    2779             :  * L2  : list of data structures to compute local DL in (Z_K/pr)^*,
    2780             :  *       and 1 + pr^a/ 1 + pr^b for various a < b <= min(2a, k)
    2781             :  * U   : base change matrices to convert a vector of local DL to DL wrt gen
    2782             :  * If MOD is not NULL, initialize G / G^MOD instead */
    2783             : static GEN
    2784      425762 : sprkinit(GEN nf, GEN pr, long k, GEN x, GEN MOD)
    2785             : {
    2786      425762 :   GEN T, p, Ld, modpr, cyc, gen, g, g0, A, prk, U, L2, ord0 = NULL;
    2787      425762 :   long f = pr_get_f(pr);
    2788             : 
    2789      425761 :   if(DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",k,pr);
    2790      425761 :   modpr = nf_to_Fq_init(nf, &pr,&T,&p);
    2791      425782 :   if (MOD)
    2792             :   {
    2793      378273 :     GEN o = subiu(powiu(p,f), 1), d = gcdii(o, MOD), fa = Z_factor(d);
    2794      378241 :     ord0 = mkvec2(o, fa); /* true order, factorization of order in G/G^MOD */
    2795      378235 :     Ld = gel(fa,1);
    2796      378235 :     if (lg(Ld) > 1 && equaliu(gel(Ld,1),2)) Ld = vecslice(Ld,2,lg(Ld)-1);
    2797             :   }
    2798             :   /* (Z_K / pr)^* */
    2799      425750 :   if (f == 1)
    2800             :   {
    2801      336641 :     g0 = g = MOD? pgener_Fp_local(p, Ld): pgener_Fp(p);
    2802      336663 :     if (!ord0) ord0 = get_arith_ZZM(subiu(p,1));
    2803             :   }
    2804             :   else
    2805             :   {
    2806       89109 :     g0 = g = MOD? gener_FpXQ_local(T, p, Ld): gener_FpXQ(T,p, &ord0);
    2807       89109 :     g = Fq_to_nf(g, modpr);
    2808       89109 :     if (typ(g) == t_POL) g = poltobasis(nf, g);
    2809             :   }
    2810      425776 :   A = gel(ord0, 1); /* Norm(pr)-1 */
    2811             :   /* If MOD != NULL, d = gcd(A, MOD): g^(A/d) has order d */
    2812      425776 :   if (k == 1)
    2813             :   {
    2814      256863 :     cyc = mkvec(A);
    2815      256860 :     gen = mkvec(g);
    2816      256856 :     prk = pr_hnf(nf,pr);
    2817      256859 :     L2 = U = NULL;
    2818             :   }
    2819             :   else
    2820             :   { /* local-gens of (1 + pr)/(1 + pr^k) = SNF-gens * U */
    2821             :     GEN AB, B, u, v, w;
    2822             :     long j, l;
    2823      168913 :     w = idealprincipalunits_i(nf, pr, k, &U);
    2824             :     /* incorporate (Z_K/pr)^*, order A coprime to B = expo(1+pr/1+pr^k)*/
    2825      168913 :     cyc = leafcopy(gel(w,1)); B = cyc_get_expo(cyc); AB = mulii(A,B);
    2826      168896 :     gen = leafcopy(gel(w,2));
    2827      168895 :     prk = gel(w,3);
    2828      168895 :     g = nfpowmodideal(nf, g, B, prk);
    2829      168915 :     g0 = Fq_pow(g0, modii(B,A), T, p); /* update primitive root */
    2830      168913 :     L2 = mkvec3(A, g, gel(w,4));
    2831      168915 :     gel(cyc,1) = AB;
    2832      168915 :     gel(gen,1) = nfmulmodideal(nf, gel(gen,1), g, prk);
    2833      168898 :     u = mulii(Fp_inv(A,B), A);
    2834      168902 :     v = subui(1, u); l = lg(U);
    2835      505154 :     for (j = 1; j < l; j++) gcoeff(U,1,j) = Fp_mul(u, gcoeff(U,1,j), AB);
    2836      168903 :     U = mkvec2(Rg_col_ei(v, lg(gen)-1, 1), U);
    2837             :   }
    2838             :   /* local-gens of (Z_K/pr^k)^* = SNF-gens * U */
    2839      425774 :   if (x)
    2840             :   {
    2841      236863 :     GEN uv = zkchineseinit(nf, idealmulpowprime(nf,x,pr,utoineg(k)), prk, x);
    2842      236854 :     gen = zkVchinese1(uv, gen);
    2843             :   }
    2844      425709 :   return mkvecn(U? 6: 4, cyc, gen, prk, mkvec3(modpr,g0,ord0), L2, U);
    2845             : }
    2846             : GEN
    2847     3978086 : sprk_get_cyc(GEN s) { return gel(s,1); }
    2848             : GEN
    2849     1966788 : sprk_get_expo(GEN s) { return cyc_get_expo(sprk_get_cyc(s)); }
    2850             : GEN
    2851      335715 : sprk_get_gen(GEN s) { return gel(s,2); }
    2852             : GEN
    2853     4909514 : sprk_get_prk(GEN s) { return gel(s,3); }
    2854             : GEN
    2855     2540295 : sprk_get_ff(GEN s) { return gel(s,4); }
    2856             : GEN
    2857     2600854 : sprk_get_pr(GEN s) { GEN ff = gel(s,4); return modpr_get_pr(gel(ff,1)); }
    2858             : /* L2 to 1 + pr / 1 + pr^k */
    2859             : static GEN
    2860     1211169 : sprk_get_L2(GEN s) { return gmael(s,5,3); }
    2861             : /* lift to nf of primitive root of k(pr) */
    2862             : static GEN
    2863      318446 : sprk_get_gnf(GEN s) { return gmael(s,5,2); }
    2864             : /* A = Npr-1, <g> = (Z_K/pr)^*, L2 to 1 + pr / 1 + pr^k */
    2865             : void
    2866           0 : sprk_get_AgL2(GEN s, GEN *A, GEN *g, GEN *L2)
    2867           0 : { GEN v = gel(s,5); *A = gel(v,1); *g = gel(v,2); *L2 = gel(v,3); }
    2868             : void
    2869     1202506 : sprk_get_U2(GEN s, GEN *U1, GEN *U2)
    2870     1202506 : { GEN v = gel(s,6); *U1 = gel(v,1); *U2 = gel(v,2); }
    2871             : static int
    2872     2540291 : sprk_is_prime(GEN s) { return lg(s) == 5; }
    2873             : 
    2874             : GEN
    2875     1966591 : famat_zlog_pr(GEN nf, GEN g, GEN e, GEN sprk, GEN mod)
    2876             : {
    2877     1966591 :   GEN x, expo = sprk_get_expo(sprk);
    2878     1966592 :   if (mod) expo = gcdii(expo,mod);
    2879     1966567 :   x = famat_makecoprime(nf, g, e, sprk_get_pr(sprk), sprk_get_prk(sprk), expo);
    2880     1966589 :   return log_prk(nf, x, sprk, mod);
    2881             : }
    2882             : /* famat_zlog_pr assuming (g,sprk.pr) = 1 */
    2883             : static GEN
    2884         196 : famat_zlog_pr_coprime(GEN nf, GEN g, GEN e, GEN sprk, GEN MOD)
    2885             : {
    2886         196 :   GEN x = famat_to_nf_modideal_coprime(nf, g, e, sprk_get_prk(sprk),
    2887             :                                        sprk_get_expo(sprk));
    2888         196 :   return log_prk(nf, x, sprk, MOD);
    2889             : }
    2890             : 
    2891             : /* o t_INT, O = [ord,fa] format for multiple of o (for Fq_log);
    2892             :  * return o in [ord,fa] format */
    2893             : static GEN
    2894      558031 : order_update(GEN o, GEN O)
    2895             : {
    2896      558031 :   GEN p = gmael(O,2,1), z = o, P, E;
    2897      558031 :   long i, j, l = lg(p);
    2898      558031 :   P = cgetg(l, t_COL);
    2899      558024 :   E = cgetg(l, t_COL);
    2900      614775 :   for (i = j = 1; i < l; i++)
    2901             :   {
    2902      614775 :     long v = Z_pvalrem(z, gel(p,i), &z);
    2903      614719 :     if (v)
    2904             :     {
    2905      599863 :       gel(P,j) = gel(p,i);
    2906      599863 :       gel(E,j) = utoipos(v); j++;
    2907      599883 :       if (is_pm1(z)) break;
    2908             :     }
    2909             :   }
    2910      557990 :   setlg(P, j);
    2911      557987 :   setlg(E, j); return mkvec2(o, mkmat2(P,E));
    2912             : }
    2913             : 
    2914             : /* a in Z_K (t_COL or t_INT), pr prime ideal, sprk = sprkinit(nf,pr,k,x),
    2915             :  * mod positive t_INT or NULL (meaning mod=0).
    2916             :  * return log(a) modulo mod on SNF-generators of (Z_K/pr^k)^* */
    2917             : GEN
    2918     2614074 : log_prk(GEN nf, GEN a, GEN sprk, GEN mod)
    2919             : {
    2920             :   GEN e, prk, g, U1, U2, y, ff, O, o, oN, gN,  N, T, p, modpr, pr, cyc;
    2921             : 
    2922     2614074 :   if (typ(a) == t_MAT) return famat_zlog_pr(nf, gel(a,1), gel(a,2), sprk, mod);
    2923     2540280 :   N = NULL;
    2924     2540280 :   ff = sprk_get_ff(sprk);
    2925     2540291 :   pr = gel(ff,1); /* modpr */
    2926     2540291 :   g = gN = gel(ff,2);
    2927     2540291 :   O = gel(ff,3); /* order of g = |Fq^*|, in [ord, fa] format */
    2928     2540291 :   o = oN = gel(O,1); /* order as a t_INT */
    2929     2540291 :   prk = sprk_get_prk(sprk);
    2930     2540308 :   modpr = nf_to_Fq_init(nf, &pr, &T, &p);
    2931     2540325 :   if (mod)
    2932             :   {
    2933     2023984 :     GEN d = gcdii(o,mod);
    2934     2023702 :     if (!equalii(o, d))
    2935             :     {
    2936      748221 :       N = diviiexact(o,d); /* > 1, coprime to p */
    2937      748172 :       a = nfpowmodideal(nf, a, N, prk);
    2938      748367 :       oN = d; /* order of g^N mod pr */
    2939             :     }
    2940             :   }
    2941     2540121 :   if (equali1(oN))
    2942      397374 :     e = gen_0;
    2943             :   else
    2944             :   {
    2945     2142849 :     if (N) { O = order_update(oN, O); gN = Fq_pow(g, N, T, p); }
    2946     2142850 :     e = Fq_log(nf_to_Fq(nf,a,modpr), gN, O, T, p);
    2947             :   }
    2948             :   /* 0 <= e < oN is correct modulo oN */
    2949     2540325 :   if (sprk_is_prime(sprk)) return mkcol(e); /* k = 1 */
    2950             : 
    2951      800013 :   sprk_get_U2(sprk, &U1,&U2);
    2952      800119 :   cyc = sprk_get_cyc(sprk);
    2953      800123 :   if (mod)
    2954             :   {
    2955      378769 :     cyc = ZV_snf_gcd(cyc, mod);
    2956      378733 :     if (signe(remii(mod,p))) return ZV_ZV_mod(ZC_Z_mul(U1,e), cyc);
    2957             :   }
    2958      746498 :   if (signe(e))
    2959             :   {
    2960      318446 :     GEN E = N? mulii(e, N): e;
    2961      318446 :     a = nfmulpowmodideal(nf, a, sprk_get_gnf(sprk), Fp_neg(E, o), prk);
    2962             :   }
    2963             :   /* a = 1 mod pr */
    2964      746498 :   y = log_prk1(nf, a, lg(U2)-1, sprk_get_L2(sprk), prk);
    2965      746552 :   if (N)
    2966             :   { /* from DL(a^N) to DL(a) */
    2967      135235 :     GEN E = gel(sprk_get_cyc(sprk), 1), q = powiu(p, Z_pval(E, p));
    2968      135232 :     y = ZC_Z_mul(y, Fp_inv(N, q));
    2969             :   }
    2970      746552 :   y = ZC_lincomb(gen_1, e, ZM_ZC_mul(U2,y), U1);
    2971      746550 :   return ZV_ZV_mod(y, cyc);
    2972             : }
    2973             : /* true nf */
    2974             : GEN
    2975       90118 : log_prk_init(GEN nf, GEN pr, long k, GEN MOD)
    2976       90118 : { return sprkinit(nf,pr,k,NULL,MOD);}
    2977             : GEN
    2978         497 : veclog_prk(GEN nf, GEN v, GEN sprk)
    2979             : {
    2980         497 :   long l = lg(v), i;
    2981         497 :   GEN w = cgetg(l, t_MAT);
    2982        1232 :   for (i = 1; i < l; i++) gel(w,i) = log_prk(nf, gel(v,i), sprk, NULL);
    2983         497 :   return w;
    2984             : }
    2985             : 
    2986             : static GEN
    2987     1371792 : famat_zlog(GEN nf, GEN fa, GEN sgn, zlog_S *S)
    2988             : {
    2989     1371792 :   long i, l0, l = lg(S->U);
    2990     1371792 :   GEN g = gel(fa,1), e = gel(fa,2), y = cgetg(l, t_COL);
    2991     1371791 :   l0 = lg(S->sprk); /* = l (trivial arch. part), or l-1 */
    2992     2847111 :   for (i=1; i < l0; i++) gel(y,i) = famat_zlog_pr(nf, g, e, gel(S->sprk,i), S->mod);
    2993     1371790 :   if (l0 != l)
    2994             :   {
    2995      188496 :     if (!sgn) sgn = nfsign_arch(nf, fa, S->archp);
    2996      188496 :     gel(y,l0) = Flc_to_ZC(sgn);
    2997             :   }
    2998     1371790 :   return y;
    2999             : }
    3000             : 
    3001             : /* assume that cyclic factors are normalized, in particular != [1] */
    3002             : static GEN
    3003      257375 : split_U(GEN U, GEN Sprk)
    3004             : {
    3005      257375 :   long t = 0, k, n, l = lg(Sprk);
    3006      257375 :   GEN vU = cgetg(l+1, t_VEC);
    3007      592333 :   for (k = 1; k < l; k++)
    3008             :   {
    3009      334956 :     n = lg(sprk_get_cyc(gel(Sprk,k))) - 1; /* > 0 */
    3010      334961 :     gel(vU,k) = vecslice(U, t+1, t+n);
    3011      334961 :     t += n;
    3012             :   }
    3013             :   /* t+1 .. lg(U)-1 */
    3014      257377 :   n = lg(U) - t - 1; /* can be 0 */
    3015      257377 :   if (!n) setlg(vU,l); else gel(vU,l) = vecslice(U, t+1, t+n);
    3016      257380 :   return vU;
    3017             : }
    3018             : 
    3019             : static void
    3020     1987890 : init_zlog_mod(zlog_S *S, GEN bid, GEN mod)
    3021             : {
    3022     1987890 :   GEN fa2 = bid_get_fact2(bid);
    3023     1987884 :   S->U = bid_get_U(bid);
    3024     1987885 :   S->hU = lg(bid_get_cyc(bid))-1;
    3025     1987884 :   S->archp = bid_get_archp(bid);
    3026     1987885 :   S->sprk = bid_get_sprk(bid);
    3027     1987885 :   S->bid = bid;
    3028     1987885 :   S->mod = mod;
    3029     1987885 :   S->P = gel(fa2,1);
    3030     1987885 :   S->k = gel(fa2,2);
    3031     1987885 :   S->no2 = lg(S->P) == lg(gel(bid_get_fact(bid),1));
    3032     1987884 : }
    3033             : void
    3034      379924 : init_zlog(zlog_S *S, GEN bid)
    3035             : {
    3036      379924 :   return init_zlog_mod(S, bid, NULL);
    3037             : }
    3038             : 
    3039             : /* a a t_FRAC/t_INT, reduce mod bid */
    3040             : static GEN
    3041          14 : Q_mod_bid(GEN bid, GEN a)
    3042             : {
    3043          14 :   GEN xZ = gcoeff(bid_get_ideal(bid),1,1);
    3044          14 :   GEN b = Rg_to_Fp(a, xZ);
    3045          14 :   if (gsigne(a) < 0) b = subii(b, xZ);
    3046          14 :   return signe(b)? b: xZ;
    3047             : }
    3048             : /* Return decomposition of a on the CRT generators blocks attached to the
    3049             :  * S->sprk and sarch; sgn = sign(a, S->arch), NULL if unknown */
    3050             : static GEN
    3051      381180 : zlog(GEN nf, GEN a, GEN sgn, zlog_S *S)
    3052             : {
    3053             :   long k, l;
    3054             :   GEN y;
    3055      381180 :   a = nf_to_scalar_or_basis(nf, a);
    3056      381181 :   switch(typ(a))
    3057             :   {
    3058      162434 :     case t_INT: break;
    3059          14 :     case t_FRAC: a = Q_mod_bid(S->bid, a); break;
    3060      218733 :     default: /* case t_COL: */
    3061             :     {
    3062             :       GEN den;
    3063      218733 :       check_nfelt(a, &den);
    3064      218746 :       if (den)
    3065             :       {
    3066         104 :         a = Q_muli_to_int(a, den);
    3067         105 :         a = mkmat2(mkcol2(a, den), mkcol2(gen_1, gen_m1));
    3068         105 :         return famat_zlog(nf, a, sgn, S);
    3069             :       }
    3070             :     }
    3071             :   }
    3072      381083 :   if (sgn)
    3073      374160 :     sgn = (lg(sgn) == 1)? NULL: leafcopy(sgn);
    3074             :   else
    3075        6923 :     sgn = (lg(S->archp) == 1)? NULL: nfsign_arch(nf, a, S->archp);
    3076      381087 :   l = lg(S->sprk);
    3077      381087 :   y = cgetg(sgn? l+1: l, t_COL);
    3078      922067 :   for (k = 1; k < l; k++)
    3079             :   {
    3080      541031 :     GEN sprk = gel(S->sprk,k);
    3081      541031 :     gel(y,k) = log_prk(nf, a, sprk, S->mod);
    3082             :   }
    3083      381036 :   if (sgn) gel(y,l) = Flc_to_ZC(sgn);
    3084      381048 :   return y;
    3085             : }
    3086             : 
    3087             : /* true nf */
    3088             : GEN
    3089       43617 : pr_basis_perm(GEN nf, GEN pr)
    3090             : {
    3091       43617 :   long f = pr_get_f(pr);
    3092             :   GEN perm;
    3093       43617 :   if (f == nf_get_degree(nf)) return identity_perm(f);
    3094       38052 :   perm = cgetg(f+1, t_VECSMALL);
    3095       38052 :   perm[1] = 1;
    3096       38052 :   if (f > 1)
    3097             :   {
    3098        2912 :     GEN H = pr_hnf(nf,pr);
    3099             :     long i, k;
    3100       10808 :     for (i = k = 2; k <= f; i++)
    3101        7896 :       if (!equali1(gcoeff(H,i,i))) perm[k++] = i;
    3102             :   }
    3103       38052 :   return perm;
    3104             : }
    3105             : 
    3106             : /* \sum U[i]*y[i], U[i] ZM, y[i] ZC. We allow lg(y) > lg(U). */
    3107             : static GEN
    3108     1752933 : ZMV_ZCV_mul(GEN U, GEN y)
    3109             : {
    3110     1752933 :   long i, l = lg(U);
    3111     1752933 :   GEN z = NULL;
    3112     1752933 :   if (l == 1) return cgetg(1,t_COL);
    3113     4131843 :   for (i = 1; i < l; i++)
    3114             :   {
    3115     2379014 :     GEN u = ZM_ZC_mul(gel(U,i), gel(y,i));
    3116     2378950 :     z = z? ZC_add(z, u): u;
    3117             :   }
    3118     1752829 :   return z;
    3119             : }
    3120             : /* A * (U[1], ..., U[d] */
    3121             : static GEN
    3122         518 : ZM_ZMV_mul(GEN A, GEN U)
    3123             : {
    3124             :   long i, l;
    3125         518 :   GEN V = cgetg_copy(U,&l);
    3126        1057 :   for (i = 1; i < l; i++) gel(V,i) = ZM_mul(A,gel(U,i));
    3127         518 :   return V;
    3128             : }
    3129             : 
    3130             : /* a = 1 mod pr, sprk mod pr^e, e >= 1 */
    3131             : static GEN
    3132      402408 : sprk_log_prk1_2(GEN nf, GEN a, GEN sprk)
    3133             : {
    3134      402408 :   GEN U1, U2, y, L2 = sprk_get_L2(sprk);
    3135      402407 :   sprk_get_U2(sprk, &U1,&U2);
    3136      402407 :   y = ZM_ZC_mul(U2, log_prk1(nf, a, lg(U2)-1, L2, sprk_get_prk(sprk)));
    3137      402398 :   return ZV_ZV_mod(y, sprk_get_cyc(sprk));
    3138             : }
    3139             : /* true nf; assume e >= 2 */
    3140             : GEN
    3141      105620 : sprk_log_gen_pr2(GEN nf, GEN sprk, long e)
    3142             : {
    3143      105620 :   GEN M, G, pr = sprk_get_pr(sprk);
    3144             :   long i, l;
    3145      105620 :   if (e == 2)
    3146             :   {
    3147       62258 :     GEN L2 = sprk_get_L2(sprk), L = gel(L2,1);
    3148       62258 :     G = gel(L,2); l = lg(G);
    3149             :   }
    3150             :   else
    3151             :   {
    3152       43362 :     GEN perm = pr_basis_perm(nf,pr), PI = nfpow_u(nf, pr_get_gen(pr), e-1);
    3153       43365 :     l = lg(perm);
    3154       43365 :     G = cgetg(l, t_VEC);
    3155       43365 :     if (typ(PI) == t_INT)
    3156             :     { /* zk_ei_mul doesn't allow t_INT */
    3157        5558 :       long N = nf_get_degree(nf);
    3158        5558 :       gel(G,1) = addiu(PI,1);
    3159        8533 :       for (i = 2; i < l; i++)
    3160             :       {
    3161        2975 :         GEN z = col_ei(N, 1);
    3162        2975 :         gel(G,i) = z; gel(z, perm[i]) = PI;
    3163             :       }
    3164             :     }
    3165             :     else
    3166             :     {
    3167       37807 :       gel(G,1) = nfadd(nf, gen_1, PI);
    3168       44590 :       for (i = 2; i < l; i++)
    3169        6783 :         gel(G,i) = nfadd(nf, gen_1, zk_ei_mul(nf, PI, perm[i]));
    3170             :     }
    3171             :   }
    3172      105623 :   M = cgetg(l, t_MAT);
    3173      233884 :   for (i = 1; i < l; i++) gel(M,i) = sprk_log_prk1_2(nf, gel(G,i), sprk);
    3174      105611 :   return M;
    3175             : }
    3176             : /* Log on bid.gen of generators of P_{1,I pr^{e-1}} / P_{1,I pr^e} (I,pr) = 1,
    3177             :  * defined implicitly via CRT. 'ind' is the index of pr in modulus
    3178             :  * factorization; true nf */
    3179             : GEN
    3180      413439 : log_gen_pr(zlog_S *S, long ind, GEN nf, long e)
    3181             : {
    3182      413439 :   GEN Uind = gel(S->U, ind);
    3183      413439 :   if (e == 1) retmkmat( gel(Uind,1) );
    3184      102925 :   return ZM_mul(Uind, sprk_log_gen_pr2(nf, gel(S->sprk,ind), e));
    3185             : }
    3186             : /* true nf */
    3187             : GEN
    3188        2037 : sprk_log_gen_pr(GEN nf, GEN sprk, long e)
    3189             : {
    3190        2037 :   if (e == 1)
    3191             :   {
    3192           0 :     long n = lg(sprk_get_cyc(sprk))-1;
    3193           0 :     retmkmat(col_ei(n, 1));
    3194             :   }
    3195        2037 :   return sprk_log_gen_pr2(nf, sprk, e);
    3196             : }
    3197             : /* a = 1 mod pr */
    3198             : GEN
    3199      274135 : sprk_log_prk1(GEN nf, GEN a, GEN sprk)
    3200             : {
    3201      274135 :   if (lg(sprk) == 5) return mkcol(gen_0); /* mod pr */
    3202      274135 :   return sprk_log_prk1_2(nf, a, sprk);
    3203             : }
    3204             : /* Log on bid.gen of generator of P_{1,f} / P_{1,f v[index]}
    3205             :  * v = vector of r1 real places */
    3206             : GEN
    3207       86104 : log_gen_arch(zlog_S *S, long index) { return gel(veclast(S->U), index); }
    3208             : 
    3209             : /* compute bid.clgp: [h,cyc] or [h,cyc,gen] */
    3210             : static GEN
    3211      258372 : bid_grp(GEN nf, GEN U, GEN cyc, GEN g, GEN F, GEN sarch)
    3212             : {
    3213      258372 :   GEN G, h = ZV_prod(cyc);
    3214             :   long c;
    3215      258389 :   if (!U) return mkvec2(h,cyc);
    3216      258039 :   c = lg(U);
    3217      258039 :   G = cgetg(c,t_VEC);
    3218      258051 :   if (c > 1)
    3219             :   {
    3220      227974 :     GEN U0, Uoo, EX = cyc_get_expo(cyc); /* exponent of bid */
    3221      227974 :     long i, hU = nbrows(U), nba = lg(sarch_get_cyc(sarch))-1; /* #f_oo */
    3222      227989 :     if (!nba) { U0 = U; Uoo = NULL; }
    3223       80397 :     else if (nba == hU) { U0 = NULL; Uoo = U; }
    3224             :     else
    3225             :     {
    3226       71255 :       U0 = rowslice(U, 1, hU-nba);
    3227       71260 :       Uoo = rowslice(U, hU-nba+1, hU);
    3228             :     }
    3229      695244 :     for (i = 1; i < c; i++)
    3230             :     {
    3231      467270 :       GEN t = gen_1;
    3232      467270 :       if (U0) t = famat_to_nf_modideal_coprime(nf, g, gel(U0,i), F, EX);
    3233      467257 :       if (Uoo) t = set_sign_mod_divisor(nf, ZV_to_Flv(gel(Uoo,i),2), t, sarch);
    3234      467252 :       gel(G,i) = t;
    3235             :     }
    3236             :   }
    3237      258051 :   return mkvec3(h, cyc, G);
    3238             : }
    3239             : 
    3240             : /* remove prime ideals of norm 2 with exponent 1 from factorization */
    3241             : static GEN
    3242      258732 : famat_strip2(GEN fa)
    3243             : {
    3244      258732 :   GEN P = gel(fa,1), E = gel(fa,2), Q, F;
    3245      258732 :   long l = lg(P), i, j;
    3246      258732 :   Q = cgetg(l, t_COL);
    3247      258729 :   F = cgetg(l, t_COL);
    3248      633697 :   for (i = j = 1; i < l; i++)
    3249             :   {
    3250      374969 :     GEN pr = gel(P,i), e = gel(E,i);
    3251      374969 :     if (!absequaliu(pr_get_p(pr), 2) || itou(e) != 1 || pr_get_f(pr) != 1)
    3252             :     {
    3253      336327 :       gel(Q,j) = pr;
    3254      336327 :       gel(F,j) = e; j++;
    3255             :     }
    3256             :   }
    3257      258728 :   setlg(Q,j);
    3258      258727 :   setlg(F,j); return mkmat2(Q,F);
    3259             : }
    3260             : static int
    3261      134086 : checkarchp(GEN v, long r1)
    3262             : {
    3263      134086 :   long i, l = lg(v);
    3264      134086 :   pari_sp av = avma;
    3265             :   GEN p;
    3266      134086 :   if (l == 1) return 1;
    3267       47157 :   if (l == 2) return v[1] > 0 && v[1] <= r1;
    3268       22020 :   p = zero_zv(r1);
    3269       66150 :   for (i = 1; i < l; i++)
    3270             :   {
    3271       44163 :     long j = v[i];
    3272       44163 :     if (j <= 0 || j > r1 || p[j]) return gc_long(av, 0);
    3273       44128 :     p[j] = 1;
    3274             :   }
    3275       21987 :   return gc_long(av, 1);
    3276             : }
    3277             : 
    3278             : /* True nf. Put ideal to form [[ideal,arch]] and set fa and fa2 to its
    3279             :  * factorization, archp to the indices of arch places */
    3280             : GEN
    3281      258712 : check_mod_factored(GEN nf, GEN ideal, GEN *fa_, GEN *fa2_, GEN *archp_, GEN MOD)
    3282             : {
    3283             :   GEN arch, x, fa, fa2, archp;
    3284             :   long R1;
    3285             : 
    3286      258712 :   R1 = nf_get_r1(nf);
    3287      258711 :   if (typ(ideal) == t_VEC && lg(ideal) == 3)
    3288             :   {
    3289      178662 :     arch = gel(ideal,2);
    3290      178662 :     ideal= gel(ideal,1);
    3291      178662 :     switch(typ(arch))
    3292             :     {
    3293       44576 :       case t_VEC:
    3294       44576 :         if (lg(arch) != R1+1)
    3295           7 :           pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    3296       44569 :         archp = vec01_to_indices(arch);
    3297       44569 :         break;
    3298      134086 :       case t_VECSMALL:
    3299      134086 :         if (!checkarchp(arch, R1))
    3300          35 :           pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    3301      134052 :         archp = arch;
    3302      134052 :         arch = indices_to_vec01(archp, R1);
    3303      134051 :         break;
    3304           0 :       default:
    3305           0 :         pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    3306             :         return NULL;/*LCOV_EXCL_LINE*/
    3307             :     }
    3308             :   }
    3309             :   else
    3310             :   {
    3311       80049 :     arch = zerovec(R1);
    3312       80041 :     archp = cgetg(1, t_VECSMALL);
    3313             :   }
    3314      258659 :   if (MOD)
    3315             :   {
    3316      214125 :     if (typ(MOD) != t_INT) pari_err_TYPE("bnrinit [incorrect cycmod]", MOD);
    3317      214125 :     if (mpodd(MOD) && lg(archp) != 1)
    3318         231 :       MOD = shifti(MOD, 1); /* ensure elements of G^MOD are >> 0 */
    3319             :   }
    3320      258659 :   if (is_nf_factor(ideal))
    3321             :   {
    3322       40250 :     fa = ideal;
    3323       40250 :     x = factorbackprime(nf, gel(fa,1), gel(fa,2));
    3324             :   }
    3325             :   else
    3326             :   {
    3327      218410 :     fa = idealfactor(nf, ideal);
    3328      218436 :     x = ideal;
    3329             :   }
    3330      258686 :   if (typ(x) != t_MAT) x = idealhnf_shallow(nf, x);
    3331      258674 :   if (lg(x) == 1) pari_err_DOMAIN("Idealstar", "ideal","=",gen_0,x);
    3332      258674 :   if (typ(gcoeff(x,1,1)) != t_INT)
    3333           7 :     pari_err_DOMAIN("Idealstar","denominator(ideal)", "!=",gen_1,x);
    3334             : 
    3335      258667 :   fa2 = famat_strip2(fa);
    3336      258664 :   if (fa_ != NULL) *fa_ = fa;
    3337      258664 :   if (fa2_ != NULL) *fa2_ = fa2;
    3338      258664 :   if (fa2_ != NULL) *archp_ = archp;
    3339      258664 :   return mkvec2(x, arch);
    3340             : }
    3341             : 
    3342             : /* True nf. Compute [[ideal,arch], [h,[cyc],[gen]], idealfact, [liste], U]
    3343             :    flag may include nf_GEN | nf_INIT */
    3344             : static GEN
    3345      258073 : Idealstarmod_i(GEN nf, GEN ideal, long flag, GEN MOD)
    3346             : {
    3347             :   long i, nbp;
    3348      258073 :   GEN y, cyc, U, u1 = NULL, fa, fa2, sprk, x_arch, x, arch, archp, E, P, sarch, gen;
    3349             : 
    3350      258073 :   x_arch = check_mod_factored(nf, ideal, &fa, &fa2, &archp, MOD);
    3351      258034 :   x = gel(x_arch, 1);
    3352      258034 :   arch = gel(x_arch, 2);
    3353             : 
    3354      258034 :   sarch = nfarchstar(nf, x, archp);
    3355      258050 :   P = gel(fa2,1);
    3356      258050 :   E = gel(fa2,2);
    3357      258050 :   nbp = lg(P)-1;
    3358      258050 :   sprk = cgetg(nbp+1,t_VEC);
    3359      258058 :   if (nbp)
    3360             :   {
    3361      218795 :     GEN t = (lg(gel(fa,1))==2)? NULL: x; /* beware fa != fa2 */
    3362      218795 :     cyc = cgetg(nbp+2,t_VEC);
    3363      218784 :     gen = cgetg(nbp+1,t_VEC);
    3364      554445 :     for (i = 1; i <= nbp; i++)
    3365             :     {
    3366      335646 :       GEN L = sprkinit(nf, gel(P,i), itou(gel(E,i)), t, MOD);
    3367      335654 :       gel(sprk,i) = L;
    3368      335654 :       gel(cyc,i) = sprk_get_cyc(L);
    3369             :       /* true gens are congruent to those mod x AND positive at archp */
    3370      335652 :       gel(gen,i) = sprk_get_gen(L);
    3371             :     }
    3372      218799 :     gel(cyc,i) = sarch_get_cyc(sarch);
    3373      218799 :     cyc = shallowconcat1(cyc);
    3374      218805 :     gen = shallowconcat1(gen);
    3375      218804 :     cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
    3376             :   }
    3377             :   else
    3378             :   {
    3379       39263 :     cyc = sarch_get_cyc(sarch);
    3380       39263 :     gen = cgetg(1,t_VEC);
    3381       39262 :     U = matid(lg(cyc)-1);
    3382       39262 :     if (flag & nf_GEN) u1 = U;
    3383             :   }
    3384      258030 :   if (MOD) cyc = ZV_snf_gcd(cyc, MOD);
    3385      258009 :   y = bid_grp(nf, u1, cyc, gen, x, sarch);
    3386      258052 :   if (!(flag & nf_INIT)) return y;
    3387      257268 :   U = split_U(U, sprk);
    3388      257269 :   return mkvec5(mkvec2(x, arch), y, mkvec2(fa,fa2), mkvec2(sprk, sarch), U);
    3389             : }
    3390             : 
    3391             : static long
    3392          63 : idealHNF_norm_pval(GEN x, GEN p)
    3393             : {
    3394          63 :   long i, v = 0, l = lg(x);
    3395         175 :   for (i = 1; i < l; i++) v += Z_pval(gcoeff(x,i,i), p);
    3396          63 :   return v;
    3397             : }
    3398             : static long
    3399          63 : sprk_get_k(GEN sprk)
    3400             : {
    3401             :   GEN pr, prk;
    3402          63 :   if (sprk_is_prime(sprk)) return 1;
    3403          63 :   pr = sprk_get_pr(sprk);
    3404          63 :   prk = sprk_get_prk(sprk);
    3405          63 :   return idealHNF_norm_pval(prk, pr_get_p(pr)) / pr_get_f(pr);
    3406             : }
    3407             : /* true nf, L a sprk */
    3408             : GEN
    3409          63 : sprk_to_bid(GEN nf, GEN L, long flag)
    3410             : {
    3411          63 :   GEN y, cyc, U, u1 = NULL, fa, fa2, arch, sarch, gen, sprk;
    3412             : 
    3413          63 :   arch = zerovec(nf_get_r1(nf));
    3414          63 :   fa = to_famat_shallow(sprk_get_pr(L), utoipos(sprk_get_k(L)));
    3415          63 :   sarch = nfarchstar(nf, NULL, cgetg(1, t_VECSMALL));
    3416          63 :   fa2 = famat_strip2(fa);
    3417          63 :   sprk = mkvec(L);
    3418          63 :   cyc = shallowconcat(sprk_get_cyc(L), sarch_get_cyc(sarch));
    3419          63 :   gen = sprk_get_gen(L);
    3420          63 :   cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
    3421          63 :   y = bid_grp(nf, u1, cyc, gen, NULL, sarch);
    3422          63 :   if (!(flag & nf_INIT)) return y;
    3423          63 :   return mkvec5(mkvec2(sprk_get_prk(L), arch), y, mkvec2(fa,fa2),
    3424             :                 mkvec2(sprk, sarch), split_U(U, sprk));
    3425             : }
    3426             : GEN
    3427      257801 : Idealstarmod(GEN nf, GEN ideal, long flag, GEN MOD)
    3428             : {
    3429      257801 :   pari_sp av = avma;
    3430      257801 :   nf = nf? checknf(nf): nfinit(pol_x(0), DEFAULTPREC);
    3431      257801 :   return gerepilecopy(av, Idealstarmod_i(nf, ideal, flag, MOD));
    3432             : }
    3433             : GEN
    3434         924 : Idealstar(GEN nf, GEN ideal, long flag) { return Idealstarmod(nf, ideal, flag, NULL); }
    3435             : GEN
    3436         273 : Idealstarprk(GEN nf, GEN pr, long k, long flag)
    3437             : {
    3438         273 :   pari_sp av = avma;
    3439         273 :   GEN z = Idealstarmod_i(nf, mkmat2(mkcol(pr),mkcols(k)), flag, NULL);
    3440         273 :   return gerepilecopy(av, z);
    3441             : }
    3442             : 
    3443             : /* FIXME: obsolete */
    3444             : GEN
    3445           0 : zidealstarinitgen(GEN nf, GEN ideal)
    3446           0 : { return Idealstar(nf,ideal, nf_INIT|nf_GEN); }
    3447             : GEN
    3448           0 : zidealstarinit(GEN nf, GEN ideal)
    3449           0 : { return Idealstar(nf,ideal, nf_INIT); }
    3450             : GEN
    3451           0 : zidealstar(GEN nf, GEN ideal)
    3452           0 : { return Idealstar(nf,ideal, nf_GEN); }
    3453             : 
    3454             : GEN
    3455         105 : idealstarmod(GEN nf, GEN ideal, long flag, GEN MOD)
    3456             : {
    3457         105 :   switch(flag)
    3458             :   {
    3459           0 :     case 0: return Idealstarmod(nf,ideal, nf_GEN, MOD);
    3460          91 :     case 1: return Idealstarmod(nf,ideal, nf_INIT, MOD);
    3461          14 :     case 2: return Idealstarmod(nf,ideal, nf_INIT|nf_GEN, MOD);
    3462           0 :     default: pari_err_FLAG("idealstar");
    3463             :   }
    3464             :   return NULL; /* LCOV_EXCL_LINE */
    3465             : }
    3466             : GEN
    3467           0 : idealstar0(GEN nf, GEN ideal,long flag) { return idealstarmod(nf, ideal, flag, NULL); }
    3468             : 
    3469             : void
    3470      218749 : check_nfelt(GEN x, GEN *den)
    3471             : {
    3472      218749 :   long l = lg(x), i;
    3473      218749 :   GEN t, d = NULL;
    3474      218749 :   if (typ(x) != t_COL) pari_err_TYPE("check_nfelt", x);
    3475      808220 :   for (i=1; i<l; i++)
    3476             :   {
    3477      589472 :     t = gel(x,i);
    3478      589472 :     switch (typ(t))
    3479             :     {
    3480      589255 :       case t_INT: break;
    3481         217 :       case t_FRAC:
    3482         217 :         if (!d) d = gel(t,2); else d = lcmii(d, gel(t,2));
    3483         217 :         break;
    3484           0 :       default: pari_err_TYPE("check_nfelt", x);
    3485             :     }
    3486             :   }
    3487      218748 :   *den = d;
    3488      218748 : }
    3489             : 
    3490             : GEN
    3491     1952135 : ZV_snf_gcd(GEN x, GEN mod)
    3492     4355198 : { pari_APPLY_same(gcdii(gel(x,i), mod)); }
    3493             : 
    3494             : /* assume a true bnf and bid */
    3495             : GEN
    3496      226975 : ideallog_units0(GEN bnf, GEN bid, GEN MOD)
    3497             : {
    3498      226975 :   GEN nf = bnf_get_nf(bnf), D, y, C, cyc;
    3499      226974 :   long j, lU = lg(bnf_get_logfu(bnf)); /* r1+r2 */
    3500             :   zlog_S S;
    3501      226973 :   init_zlog_mod(&S, bid, MOD);
    3502      226972 :   if (!S.hU) return zeromat(0,lU);
    3503      226972 :   cyc = bid_get_cyc(bid);
    3504      226971 :   D = nfsign_fu(bnf, bid_get_archp(bid));
    3505      226963 :   y = cgetg(lU, t_MAT);
    3506      226958 :   if ((C = bnf_build_cheapfu(bnf)))
    3507      374100 :   { for (j = 1; j < lU; j++) gel(y,j) = zlog(nf, gel(C,j), gel(D,j), &S); }
    3508             :   else
    3509             :   {
    3510          49 :     long i, l = lg(S.U), l0 = lg(S.sprk);
    3511             :     GEN X, U;
    3512          49 :     if (!(C = bnf_compactfu_mat(bnf))) bnf_build_units(bnf); /* error */
    3513          49 :     X = gel(C,1); U = gel(C,2);
    3514         147 :     for (j = 1; j < lU; j++) gel(y,j) = cgetg(l, t_COL);
    3515         126 :     for (i = 1; i < l0; i++)
    3516             :     {
    3517          77 :       GEN sprk = gel(S.sprk, i);
    3518          77 :       GEN Xi = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
    3519         231 :       for (j = 1; j < lU; j++)
    3520         154 :         gcoeff(y,i,j) = famat_zlog_pr_coprime(nf, Xi, gel(U,j), sprk, MOD);
    3521             :     }
    3522          49 :     if (l0 != l)
    3523          56 :       for (j = 1; j < lU; j++) gcoeff(y,l0,j) = Flc_to_ZC(gel(D,j));
    3524             :   }
    3525      226958 :   y = vec_prepend(y, zlog(nf, bnf_get_tuU(bnf), nfsign_tu(bnf, S.archp), &S));
    3526      601203 :   for (j = 1; j <= lU; j++)
    3527      374244 :     gel(y,j) = ZV_ZV_mod(ZMV_ZCV_mul(S.U, gel(y,j)), cyc);
    3528      226959 :   return y;
    3529             : }
    3530             : GEN
    3531          84 : ideallog_units(GEN bnf, GEN bid)
    3532          84 : { return ideallog_units0(bnf, bid, NULL); }
    3533             : GEN
    3534         518 : log_prk_units(GEN nf, GEN D, GEN sprk)
    3535             : {
    3536         518 :   GEN L, Ltu = log_prk(nf, gel(D,1), sprk, NULL);
    3537         518 :   D = gel(D,2);
    3538         518 :   if (lg(D) != 3 || typ(gel(D,2)) != t_MAT) L = veclog_prk(nf, D, sprk);
    3539             :   else
    3540             :   {
    3541          21 :     GEN X = gel(D,1), U = gel(D,2);
    3542          21 :     long j, lU = lg(U);
    3543          21 :     X = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
    3544          21 :     L = cgetg(lU, t_MAT);
    3545          63 :     for (j = 1; j < lU; j++)
    3546          42 :       gel(L,j) = famat_zlog_pr_coprime(nf, X, gel(U,j), sprk, NULL);
    3547             :   }
    3548         518 :   return vec_prepend(L, Ltu);
    3549             : }
    3550             : 
    3551             : static GEN
    3552     1381011 : ideallog_i(GEN nf, GEN x, zlog_S *S)
    3553             : {
    3554     1381011 :   pari_sp av = avma;
    3555             :   GEN y;
    3556     1381011 :   if (!S->hU) return cgetg(1, t_COL);
    3557     1378715 :   if (typ(x) == t_MAT)
    3558     1371687 :     y = famat_zlog(nf, x, NULL, S);
    3559             :   else
    3560        7028 :     y = zlog(nf, x, NULL, S);
    3561     1378707 :   y = ZMV_ZCV_mul(S->U, y);
    3562     1378706 :   return gerepileupto(av, ZV_ZV_mod(y, bid_get_cyc(S->bid)));
    3563             : }
    3564             : GEN
    3565     1387691 : ideallogmod(GEN nf, GEN x, GEN bid, GEN mod)
    3566             : {
    3567             :   zlog_S S;
    3568     1387691 :   if (!nf)
    3569             :   {
    3570        6671 :     if (mod) pari_err_IMPL("Zideallogmod");
    3571        6671 :     return Zideallog(bid, x);
    3572             :   }
    3573     1381020 :   checkbid(bid); init_zlog_mod(&S, bid, mod);
    3574     1381011 :   return ideallog_i(checknf(nf), x, &S);
    3575             : }
    3576             : GEN
    3577       13762 : ideallog(GEN nf, GEN x, GEN bid) { return ideallogmod(nf, x, bid, NULL); }
    3578             : 
    3579             : /*************************************************************************/
    3580             : /**                                                                     **/
    3581             : /**               JOIN BID STRUCTURES, IDEAL LISTS                      **/
    3582             : /**                                                                     **/
    3583             : /*************************************************************************/
    3584             : /* bid1, bid2: for coprime modules m1 and m2 (without arch. part).
    3585             :  * Output: bid for m1 m2 */
    3586             : static GEN
    3587         469 : join_bid(GEN nf, GEN bid1, GEN bid2)
    3588             : {
    3589         469 :   pari_sp av = avma;
    3590             :   long nbgen, l1,l2;
    3591             :   GEN I1,I2, G1,G2, sprk1,sprk2, cyc1,cyc2, sarch;
    3592         469 :   GEN sprk, fa,fa2, U, cyc, y, u1 = NULL, x, gen;
    3593             : 
    3594         469 :   I1 = bid_get_ideal(bid1);
    3595         469 :   I2 = bid_get_ideal(bid2);
    3596         469 :   if (gequal1(gcoeff(I1,1,1))) return bid2; /* frequent trivial case */
    3597         259 :   G1 = bid_get_grp(bid1);
    3598         259 :   G2 = bid_get_grp(bid2);
    3599         259 :   x = idealmul(nf, I1,I2);
    3600         259 :   fa = famat_mul_shallow(bid_get_fact(bid1), bid_get_fact(bid2));
    3601         259 :   fa2= famat_mul_shallow(bid_get_fact2(bid1), bid_get_fact2(bid2));
    3602         259 :   sprk1 = bid_get_sprk(bid1);
    3603         259 :   sprk2 = bid_get_sprk(bid2);
    3604         259 :   sprk = shallowconcat(sprk1, sprk2);
    3605             : 
    3606         259 :   cyc1 = abgrp_get_cyc(G1); l1 = lg(cyc1);
    3607         259 :   cyc2 = abgrp_get_cyc(G2); l2 = lg(cyc2);
    3608         259 :   gen = (lg(G1)>3 && lg(G2)>3)? gen_1: NULL;
    3609         259 :   nbgen = l1+l2-2;
    3610         259 :   cyc = ZV_snf_group(shallowconcat(cyc1,cyc2), &U, gen? &u1: NULL);
    3611         259 :   if (nbgen)
    3612             :   {
    3613         259 :     GEN U1 = bid_get_U(bid1), U2 = bid_get_U(bid2);
    3614           0 :     U1 = l1==1? const_vec(lg(sprk1), cgetg(1,t_MAT))
    3615         259 :               : ZM_ZMV_mul(vecslice(U, 1, l1-1),   U1);
    3616           0 :     U2 = l2==1? const_vec(lg(sprk2), cgetg(1,t_MAT))
    3617         259 :               : ZM_ZMV_mul(vecslice(U, l1, nbgen), U2);
    3618         259 :     U = shallowconcat(U1, U2);
    3619             :   }
    3620             :   else
    3621           0 :     U = const_vec(lg(sprk), cgetg(1,t_MAT));
    3622             : 
    3623         259 :   if (gen)
    3624             :   {
    3625         259 :     GEN uv = zkchinese1init2(nf, I2, I1, x);
    3626         518 :     gen = shallowconcat(zkVchinese1(gel(uv,1), abgrp_get_gen(G1)),
    3627         259 :                         zkVchinese1(gel(uv,2), abgrp_get_gen(G2)));
    3628             :   }
    3629         259 :   sarch = bid_get_sarch(bid1); /* trivial */
    3630         259 :   y = bid_grp(nf, u1, cyc, gen, x, sarch);
    3631         259 :   x = mkvec2(x, bid_get_arch(bid1));
    3632         259 :   y = mkvec5(x, y, mkvec2(fa, fa2), mkvec2(sprk, sarch), U);
    3633         259 :   return gerepilecopy(av,y);
    3634             : }
    3635             : 
    3636             : typedef struct _ideal_data {
    3637             :   GEN nf, emb, L, pr, prL, sgnU, archp;
    3638             : } ideal_data;
    3639             : 
    3640             : /* z <- ( z | f(v[i])_{i=1..#v} ) */
    3641             : static void
    3642      758017 : concat_join(GEN *pz, GEN v, GEN (*f)(ideal_data*,GEN), ideal_data *data)
    3643             : {
    3644      758017 :   long i, nz, lv = lg(v);
    3645             :   GEN z, Z;
    3646      758017 :   if (lv == 1) return;
    3647      223147 :   z = *pz; nz = lg(z)-1;
    3648      223147 :   *pz = Z = cgetg(lv + nz, typ(z));
    3649      371665 :   for (i = 1; i <=nz; i++) gel(Z,i) = gel(z,i);
    3650      223319 :   Z += nz;
    3651      491937 :   for (i = 1; i < lv; i++) gel(Z,i) = f(data, gel(v,i));
    3652             : }
    3653             : static GEN
    3654         469 : join_idealinit(ideal_data *D, GEN x)
    3655         469 : { return join_bid(D->nf, x, D->prL); }
    3656             : static GEN
    3657      268445 : join_ideal(ideal_data *D, GEN x)
    3658      268445 : { return idealmulpowprime(D->nf, x, D->pr, D->L); }
    3659             : static GEN
    3660         448 : join_unit(ideal_data *D, GEN x)
    3661             : {
    3662         448 :   GEN bid = join_idealinit(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
    3663         448 :   if (lg(u) != 1) v = shallowconcat(u, v);
    3664         448 :   return mkvec2(bid, v);
    3665             : }
    3666             : 
    3667             : GEN
    3668          49 : log_prk_units_init(GEN bnf)
    3669             : {
    3670          49 :   GEN U = bnf_has_fu(bnf);
    3671          49 :   if (U) U = matalgtobasis(bnf_get_nf(bnf), U);
    3672          21 :   else if (!(U = bnf_compactfu_mat(bnf))) (void)bnf_build_units(bnf);
    3673          49 :   return mkvec2(bnf_get_tuU(bnf), U);
    3674             : }
    3675             : /*  flag & nf_GEN : generators, otherwise no
    3676             :  *  flag &2 : units, otherwise no
    3677             :  *  flag &4 : ideals in HNF, otherwise bid
    3678             :  *  flag &8 : omit ideals which cannot be conductors (pr^1 with Npr=2) */
    3679             : static GEN
    3680       11333 : Ideallist(GEN bnf, ulong bound, long flag)
    3681             : {
    3682       11333 :   const long do_units = flag & 2, big_id = !(flag & 4), cond = flag & 8;
    3683       11333 :   const long istar_flag = (flag & nf_GEN) | nf_INIT;
    3684             :   pari_sp av;
    3685             :   long i, j;
    3686       11333 :   GEN nf, z, p, fa, id, BOUND, U, empty = cgetg(1,t_VEC);
    3687             :   forprime_t S;
    3688             :   ideal_data ID;
    3689             :   GEN (*join_z)(ideal_data*, GEN);
    3690             : 
    3691       11333 :   if (do_units)
    3692             :   {
    3693          21 :     bnf = checkbnf(bnf);
    3694          21 :     nf = bnf_get_nf(bnf);
    3695          21 :     join_z = &join_unit;
    3696             :   }
    3697             :   else
    3698             :   {
    3699       11312 :     nf = checknf(bnf);
    3700       11312 :     join_z = big_id? &join_idealinit: &join_ideal;
    3701             :   }
    3702       11333 :   if ((long)bound <= 0) return empty;
    3703       11333 :   id = matid(nf_get_degree(nf));
    3704       11333 :   if (big_id) id = Idealstar(nf,id, istar_flag);
    3705             : 
    3706             :   /* z[i] will contain all "objects" of norm i. Depending on flag, this means
    3707             :    * an ideal, a bid, or a couple [bid, log(units)]. Such objects are stored
    3708             :    * in vectors, computed one primary component at a time; join_z
    3709             :    * reconstructs the global object */
    3710       11333 :   BOUND = utoipos(bound);
    3711       11333 :   z = const_vec(bound, empty);
    3712       11333 :   U = do_units? log_prk_units_init(bnf): NULL;
    3713       11333 :   gel(z,1) = mkvec(U? mkvec2(id, empty): id);
    3714       11333 :   ID.nf = nf;
    3715             : 
    3716       11333 :   p = cgetipos(3);
    3717       11333 :   u_forprime_init(&S, 2, bound);
    3718       11333 :   av = avma;
    3719       92469 :   while ((p[2] = u_forprime_next(&S)))
    3720             :   {
    3721       81603 :     if (DEBUGLEVEL>1) err_printf("%ld ",p[2]);
    3722       81603 :     fa = idealprimedec_limit_norm(nf, p, BOUND);
    3723      162650 :     for (j = 1; j < lg(fa); j++)
    3724             :     {
    3725       81514 :       GEN pr = gel(fa,j), z2 = leafcopy(z);
    3726       81513 :       ulong Q, q = upr_norm(pr);
    3727             :       long l;
    3728       81513 :       ID.pr = ID.prL = pr;
    3729       81513 :       if (cond && q == 2) { l = 2; Q = 4; } else { l = 1; Q = q; }
    3730      184073 :       for (; Q <= bound; l++, Q *= q) /* add pr^l */
    3731             :       {
    3732             :         ulong iQ;
    3733      103038 :         ID.L = utoipos(l);
    3734      103039 :         if (big_id) {
    3735         210 :           ID.prL = Idealstarprk(nf, pr, l, istar_flag);
    3736         210 :           if (U)
    3737         189 :             ID.emb = Q == 2? empty
    3738         189 :                            : log_prk_units(nf, U, gel(bid_get_sprk(ID.prL),1));
    3739             :         }
    3740      860565 :         for (iQ = Q,i = 1; iQ <= bound; iQ += Q,i++)
    3741      758005 :           concat_join(&gel(z,iQ), gel(z2,i), join_z, &ID);
    3742             :       }
    3743             :     }
    3744       81136 :     if (gc_needed(av,1))
    3745             :     {
    3746          18 :       if(DEBUGMEM>1) pari_warn(warnmem,"Ideallist");
    3747          18 :       z = gerepilecopy(av, z);
    3748             :     }
    3749             :   }
    3750       11333 :   return z;
    3751             : }
    3752             : GEN
    3753          63 : gideallist(GEN bnf, GEN B, long flag)
    3754             : {
    3755          63 :   pari_sp av = avma;
    3756          63 :   if (typ(B) != t_INT)
    3757             :   {
    3758           0 :     B = gfloor(B);
    3759           0 :     if (typ(B) != t_INT) pari_err_TYPE("ideallist", B);
    3760           0 :     if (signe(B) < 0) B = gen_0;
    3761             :   }
    3762          63 :   if (signe(B) < 0)
    3763             :   {
    3764          28 :     if (flag != 4) pari_err_IMPL("ideallist with bid for single norm");
    3765          28 :     return gerepilecopy(av, ideals_by_norm(checknf(bnf), absi(B)));
    3766             :   }
    3767          35 :   if (flag < 0 || flag > 15) pari_err_FLAG("ideallist");
    3768          35 :   return gerepilecopy(av, Ideallist(bnf, itou(B), flag));
    3769             : }
    3770             : GEN
    3771       11298 : ideallist0(GEN bnf, long bound, long flag)
    3772             : {
    3773       11298 :   pari_sp av = avma;
    3774       11298 :   if (flag < 0 || flag > 15) pari_err_FLAG("ideallist");
    3775       11298 :   return gerepilecopy(av, Ideallist(bnf, bound, flag));
    3776             : }
    3777             : GEN
    3778       10563 : ideallist(GEN bnf,long bound) { return ideallist0(bnf,bound,4); }
    3779             : 
    3780             : /* bid = for module m (without arch. part), arch = archimedean part.
    3781             :  * Output: bid for [m,arch] */
    3782             : static GEN
    3783          42 : join_bid_arch(GEN nf, GEN bid, GEN archp)
    3784             : {
    3785          42 :   pari_sp av = avma;
    3786             :   GEN G, U;
    3787          42 :   GEN sprk, cyc, y, u1 = NULL, x, sarch, gen;
    3788             : 
    3789          42 :   checkbid(bid);
    3790          42 :   G = bid_get_grp(bid);
    3791          42 :   x = bid_get_ideal(bid);
    3792          42 :   sarch = nfarchstar(nf, bid_get_ideal(bid), archp);
    3793          42 :   sprk = bid_get_sprk(bid);
    3794             : 
    3795          42 :   gen = (lg(G)>3)? gel(G,3): NULL;
    3796          42 :   cyc = diagonal_shallow(shallowconcat(gel(G,2), sarch_get_cyc(sarch)));
    3797          42 :   cyc = ZM_snf_group(cyc, &U, gen? &u1: NULL);
    3798          42 :   y = bid_grp(nf, u1, cyc, gen, x, sarch);
    3799          42 :   U = split_U(U, sprk);
    3800          42 :   y = mkvec5(mkvec2(x, archp), y, gel(bid,3), mkvec2(sprk, sarch), U);
    3801          42 :   return gerepilecopy(av,y);
    3802             : }
    3803             : static GEN
    3804          42 : join_arch(ideal_data *D, GEN x) {
    3805          42 :   return join_bid_arch(D->nf, x, D->archp);
    3806             : }
    3807             : static GEN
    3808          14 : join_archunit(ideal_data *D, GEN x) {
    3809          14 :   GEN bid = join_arch(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
    3810          14 :   if (lg(u) != 1) v = shallowconcat(u, v);
    3811          14 :   return mkvec2(bid, v);
    3812             : }
    3813             : 
    3814             : /* L from ideallist, add archimedean part */
    3815             : GEN
    3816          14 : ideallistarch(GEN bnf, GEN L, GEN arch)
    3817             : {
    3818             :   pari_sp av;
    3819          14 :   long i, j, l = lg(L), lz;
    3820             :   GEN v, z, V, nf;
    3821             :   ideal_data ID;
    3822             :   GEN (*join_z)(ideal_data*, GEN);
    3823             : 
    3824          14 :   if (typ(L) != t_VEC) pari_err_TYPE("ideallistarch",L);
    3825          14 :   if (l == 1) return cgetg(1,t_VEC);
    3826          14 :   z = gel(L,1);
    3827          14 :   if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
    3828          14 :   z = gel(z,1); /* either a bid or [bid,U] */
    3829          14 :   ID.archp = vec01_to_indices(arch);
    3830          14 :   if (lg(z) == 3)
    3831             :   { /* [bid,U]: do units */
    3832           7 :     bnf = checkbnf(bnf); nf = bnf_get_nf(bnf);
    3833           7 :     if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
    3834           7 :     ID.emb = zm_to_ZM( rowpermute(nfsign_units(bnf,NULL,1), ID.archp) );
    3835           7 :     join_z = &join_archunit;
    3836             :   }
    3837             :   else
    3838             :   {
    3839           7 :     join_z = &join_arch;
    3840           7 :     nf = checknf(bnf);
    3841             :   }
    3842          14 :   ID.nf = nf;
    3843          14 :   av = avma; V = cgetg(l, t_VEC);
    3844          63 :   for (i = 1; i < l; i++)
    3845             :   {
    3846          49 :     z = gel(L,i); lz = lg(z);
    3847          49 :     gel(V,i) = v = cgetg(lz,t_VEC);
    3848          91 :     for (j=1; j<lz; j++) gel(v,j) = join_z(&ID, gel(z,j));
    3849             :   }
    3850          14 :   return gerepilecopy(av,V);
    3851             : }

Generated by: LCOV version 1.14